Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Move functions from general_utils to eva #59

Merged
merged 2 commits into from
Nov 17, 2023
Merged

Move functions from general_utils to eva #59

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
0709341
Create eva.py
stellema Nov 17, 2023
b29a647
Update eva.py
stellema Nov 17, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions docs/reference/index.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,7 @@ by running ``bias_correction`` at the command line
unseen.bias_correction
unseen.bootstrap
unseen.dask_setup
unseen.eva
unseen.fileio
unseen.general_utils
unseen.independence
Expand Down
10 changes: 5 additions & 5 deletions docs/user_guide/worked_example-HadGEM3-GC31-MM.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,7 @@
"from unseen import independence\n",
"from unseen import bias_correction\n",
"from unseen import similarity\n",
"from unseen import general_utils\n",
"from unseen import eva\n",
"from unseen import stability\n",
"from unseen import moments"
]
Expand Down Expand Up @@ -1026,7 +1026,7 @@
"metadata": {},
"outputs": [],
"source": [
"agcd_shape, agcd_loc, agcd_scale = general_utils.fit_gev(agcd_ds['pr'].values)"
"agcd_shape, agcd_loc, agcd_scale = eva.fit_gev(agcd_ds['pr'].values)"
]
},
{
Expand Down Expand Up @@ -3101,12 +3101,12 @@
],
"source": [
"model_da_indep.plot.hist(bins=50, density=True, alpha=0.7, facecolor='tab:blue')\n",
"model_raw_shape, model_raw_loc, model_raw_scale = general_utils.fit_gev(model_da_indep.values, generate_estimates=True)\n",
"model_raw_shape, model_raw_loc, model_raw_scale = eva.fit_gev(model_da_indep.values, generate_estimates=True)\n",
"model_raw_pdf = gev.pdf(xvals, model_raw_shape, model_raw_loc, model_raw_scale)\n",
"plt.plot(xvals, model_raw_pdf, color='tab:blue', linewidth=4.0, label='model')\n",
"\n",
"model_da_bc.plot.hist(bins=50, density=True, alpha=0.7, facecolor='tab:orange')\n",
"model_bc_shape, model_bc_loc, model_bc_scale = general_utils.fit_gev(model_da_bc.values, generate_estimates=True)\n",
"model_bc_shape, model_bc_loc, model_bc_scale = eva.fit_gev(model_da_bc.values, generate_estimates=True)\n",
"model_bc_pdf = gev.pdf(xvals, model_bc_shape, model_bc_loc, model_bc_scale)\n",
"plt.plot(xvals, model_bc_pdf, color='tab:orange', linewidth=4.0, label='model (corrected)')\n",
"\n",
Expand Down Expand Up @@ -3335,7 +3335,7 @@
"source": [
"fig = plt.figure(figsize=[6, 4])\n",
"ax = fig.add_subplot()\n",
"general_utils.plot_gev_return_curve(\n",
"eva.plot_gev_return_curve(\n",
" ax,\n",
" model_da_bc_stacked,\n",
" rx5day_max,\n",
Expand Down
33 changes: 16 additions & 17 deletions docs/user_guide/worked_example.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,7 +7,7 @@ Hobart extreme rainfall
During April of 1960,
a multi-day rainfall event caused severe flooding in Hobart
(see the `Tasmanian flood history <http://www.bom.gov.au/tas/flood/flood_history/flood_history.shtml#yr1960_1969>`__)
page for details).
page for details).

In this worked example,
we'll estimate the likelihood of this record rainfall event by applying the UNSEEN approach to
Expand Down Expand Up @@ -41,7 +41,7 @@ and apply simple temporal and spatial aggregation.

We can use it extract the Hobart grid point from the AGCD data files
and apply the necessary temporal aggregation to calculate the commonly used Rx5day metric
(the highest annual 5-day rainfall total).
(the highest annual 5-day rainfall total).

.. code-block:: python

Expand Down Expand Up @@ -101,7 +101,7 @@ selection and aggregation.
date_issued: 2023-05-19 06:19:17
attribution: Data should be cited as : Australian Bureau of...
copyright: (C) Copyright Commonwealth of Australia 2023, ...
history:
history:


It can be a good idea to compute the Dataset before going too much further with the analysis,
Expand Down Expand Up @@ -167,10 +167,10 @@ to get an estimate of the likelihood of the 1960 event.
.. code-block:: python

from scipy.stats import genextreme as gev
from unseen import general_utils
from unseen import eva

agcd_shape, agcd_loc, agcd_scale = eva.fit_gev(agcd_ds['pr'].values)

agcd_shape, agcd_loc, agcd_scale = general_utils.fit_gev(agcd_ds['pr'].values)

event_probability = gev.sf(rx5day_max, agcd_shape, loc=agcd_loc, scale=agcd_scale)
event_return_period = 1. / event_probability
event_percentile = (1 - event_probability) * 100
Expand Down Expand Up @@ -199,7 +199,7 @@ For example:
cat HadGEM3-GC31-MM_dcppA-hindcast_pr_files.txt


.. code-block:: none
.. code-block:: none

/g/data/oi10/replicas/CMIP6/DCPP/MOHC/HadGEM3-GC31-MM/dcppA-hindcast/s1960-r1i1p1f2/day/pr/gn/v20200417/pr_day_HadGEM3-GC31-MM_dcppA-hindcast_s1960-r1i1p1f2_gn_19601101-19601230.nc
/g/data/oi10/replicas/CMIP6/DCPP/MOHC/HadGEM3-GC31-MM/dcppA-hindcast/s1960-r1i1p1f2/day/pr/gn/v20200417/pr_day_HadGEM3-GC31-MM_dcppA-hindcast_s1960-r1i1p1f2_gn_19610101-19611230.nc
Expand Down Expand Up @@ -231,14 +231,14 @@ We've used similar keyword arguments as for the AGCD data
(``open_mfforecast`` uses ``open_dataset`` to open each individual file)
with a couple of additions:

- The ``n_ensemble_members`` and ``n_time_files`` arguments help the function sort the contents of the input file list
- The ``reset_times`` option ensures that after resampling (e.g. here we calculate the annual mean from daily data) the month assigned to each time axis value matches the initialisation month
- The ``complete_time_agg_periods`` argument makes sure that incomplete calendar years (e.g. the first year for a forecast that starts in November) aren't included
- The ``n_ensemble_members`` and ``n_time_files`` arguments help the function sort the contents of the input file list
- The ``reset_times`` option ensures that after resampling (e.g. here we calculate the annual mean from daily data) the month assigned to each time axis value matches the initialisation month
- The ``complete_time_agg_periods`` argument makes sure that incomplete calendar years (e.g. the first year for a forecast that starts in November) aren't included

.. code-block:: python

print(model_ds)


.. code-block:: none

Expand Down Expand Up @@ -330,7 +330,7 @@ To perform this test, we can use the ``independence`` module:

from unseen import independence

mean_correlations, null_correlation_bounds = independence.run_tests(model_ds['pr'])
mean_correlations, null_correlation_bounds = independence.run_tests(model_ds['pr'])
independence.create_plot(
mean_correlations,
null_correlation_bounds,
Expand Down Expand Up @@ -392,12 +392,12 @@ to see the effect of the bias correction.
import matplotlib.pyplot as plt

model_da_indep.plot.hist(bins=50, density=True, alpha=0.7, facecolor='tab:blue')
model_raw_shape, model_raw_loc, model_raw_scale = general_utils.fit_gev(model_da_indep.values, generate_estimates=True)
model_raw_shape, model_raw_loc, model_raw_scale = eva.fit_gev(model_da_indep.values, generate_estimates=True)
model_raw_pdf = gev.pdf(xvals, model_raw_shape, model_raw_loc, model_raw_scale)
plt.plot(xvals, model_raw_pdf, color='tab:blue', linewidth=4.0, label='model')

model_da_bc.plot.hist(bins=50, density=True, alpha=0.7, facecolor='tab:orange')
model_bc_shape, model_bc_loc, model_bc_scale = general_utils.fit_gev(model_da_bc.values, generate_estimates=True)
model_bc_shape, model_bc_loc, model_bc_scale = eva.fit_gev(model_da_bc.values, generate_estimates=True)
model_bc_pdf = gev.pdf(xvals, model_bc_shape, model_bc_loc, model_bc_scale)
plt.plot(xvals, model_bc_pdf, color='tab:orange', linewidth=4.0, label='model (corrected)')

Expand Down Expand Up @@ -529,13 +529,13 @@ Once we've stacked our model data so it's one dimensional,
units: mm d-1
bias_correction_method: multiplicative
bias_correction_period: 1970-01-01-2018-12-30


.. code-block:: python

fig = plt.figure(figsize=[6, 4])
ax = fig.add_subplot()
general_utils.plot_gev_return_curve(
eva.plot_gev_return_curve(
ax,
model_da_bc_stacked,
rx5day_max,
Expand All @@ -555,4 +555,3 @@ Once we've stacked our model data so it's one dimensional,

.. image:: return_curve.png
:width: 700

Loading
Loading