On Friday, August 26, 1983, Bilbao was celebrating its Aste Nagusia or Great Week, the main annual festivity in the city, when it and other municipalities of the Basque Country, Burgos, and Cantabria suffered devastating flooding due to heavy rains. In 24 hours, the volume of water registered 600 liters per square meter. Across all the affected areas, the weather service recorded 1.5 billion tons of water. In areas of Bilbao, the water reached a height of 5 meters (15 feet). Transportation, electricity and gas services, drinking water, food, telephone, and many other basic services were severely affected. 32 people died in Biscay, 4 people died in Cantabria, 2 people died in Alava, and 2 people died Burgos. 5 more people went missing.
7.4 How often will such rainfall happen?
How often does it rain 50 mm in 1 day? What about 100 mm in 1 day? How big is a “once-in-a-century event”?
Let’s examine Bilbao’s daily rainfall (mm), between 1947 to 2021
import stuff
import matplotlib.pyplot as pltimport numpy as npimport pandas as pdimport seaborn as snssns.set_theme(style="ticks", font_scale=1.5)import urllib.requestfrom matplotlib.dates import DateFormatterimport matplotlib.dates as mdatesimport altair as altalt.data_transformers.disable_max_rows()from scipy.stats import genextreme
define function and download data
def download_data(station_name, station_code): url_daily ='https://www.ncei.noaa.gov/data/global-historical-climatology-network-daily/access/' url_monthly ='https://www.ncei.noaa.gov/data/gsom/access/'# download daily data - uncomment to make this work urllib.request.urlretrieve(url_daily + station_code +'.csv', station_name +'_daily.csv')# download monthly data urllib.request.urlretrieve(url_monthly + station_code +'.csv', station_name +'_monthly.csv')download_data('BILBAO', 'SPE00120611')
load data and plot precipitation
df = pd.read_csv('BILBAO_daily.csv', sep=",", parse_dates=['DATE'], index_col='DATE')# IMPORTANT!! daily precipitation data is in tenths of mm, divide by 10 to get it in mm.df['PRCP'] = df['PRCP'] /10fig, ax = plt.subplots(figsize=(10,7))ax.plot(df['PRCP'])ax.set(xlabel="date", ylabel="daily rainfall (mm)", title="Bilbao, Spain, 1947--2021" )ax.annotate("26 August 1983", xy=(pd.to_datetime('1983-08-26'), 2500), xycoords='data', xytext=(0.7, 0.95), textcoords='axes fraction', fontsize=16, va="center", arrowprops=dict(facecolor='black', shrink=0.05));
On the week of 22-28 August 1983, Bilbao’s weather station measured 45 cm of rainfall!
Show the code
fig, ax = plt.subplots(figsize=(10,7))one_week = df.loc['1983-08-22':'1983-08-28', 'PRCP']bars = ax.bar(one_week.index, one_week)ax.set_xlabel("date")ax.set_ylabel("daily rainfall (mm)")ax.set_title("Bilbao, Spain, August 1983")# write daily rainfallfor i inrange(len(one_week)): ax.text(one_week.index[i], one_week.iloc[i], f"{one_week.iloc[i]:.0f}", ha="center", fontsize=16)# ax.text(one_week.index[i], one_week[i], f"{one_week[i]:.0f}", ha="center", fontsize=16);ax.text(0.1, 0.8, f"Total rainfall during this week:\n{one_week.sum():.0f} mm", transform=ax.transAxes, fontsize=16)# Define the date format# https://docs.python.org/2/library/datetime.html#strftime-and-strptime-behaviordate_form = DateFormatter("%b-%d")ax.xaxis.set_major_formatter(date_form)# Ensure a major tick for each day using (interval=1)# https://matplotlib.org/stable/api/dates_api.html#date-tickersax.xaxis.set_major_locator(mdates.DayLocator(interval=1))
Let’s analyze this data and find out how rare such events are. First we need to find the annual maximum for each hydrological year. Do we see any seasonal patterns with our eyes? Play with the widget below.
widget
# Altair only recognizes column data; it ignores index values.# You can plot the index data by first resetting the index# I know that I've just made 'DATE' the index, but I want to have this here nonetheless so I can refer to this in the futuredf_new = df.reset_index()source = df_new[['DATE', 'PRCP']]brush = alt.selection_interval(encodings=['x'])base = alt.Chart(source).mark_line().encode( x ='DATE:T', y ='PRCP:Q').properties( width=600, height=200)upper = base.encode( alt.X('DATE:T', scale=alt.Scale(domain=brush)), alt.Y('PRCP:Q', scale=alt.Scale(domain=(0,100))))lower = base.properties( height=60).add_params(brush)alt.vconcat(upper, lower)