You can run this notebook in a live session or view it on Github.

Diagnosing Potential Predictability¶

This demo demonstrates climpred’s capabilities to diagnose areas containing potentially predictable variations from a control/verification alone without requiring multi-member, multi-initialization simulations. This notebook identifies the slow components of internal variability that indicate potential predictability. Here, we showcase a set of methods to show regions indicating probabilities for decadal predictability.

# linting
%matplotlib inline
import climpred
# Sea surface temperature
varname = "tos"

Diagnostic Potential Predictability (DPP)¶

We can first use the Resplandy et al.  and Séférian et al.  method for computing the unbiased climpred.stats.dpp() by not chunking the time dimension.

# calculate DPP with m=10
DPP10 = climpred.stats.dpp(control3d, m=10, chunk=False)
# calculate a threshold by random shuffling (based on bootstrapping with replacement at 95% significance level)
threshold = climpred.bootstrap.dpp_threshold(
control3d, m=10, chunk=False, iterations=10, sig=95
)
# plot grid cells where DPP above threshold
DPP10.where(DPP10 > threshold).plot(
yincrease=False, vmin=-0.1, vmax=0.6, cmap="viridis"
) Now, we can turn on chunking (default) to use the Boer  method.

# chunk = True signals the Boer 2004 method
DPP10 = climpred.stats.dpp(control3d, m=10, chunk=True)
threshold = climpred.bootstrap.dpp_threshold(
control3d, m=10, chunk=True, iterations=10, sig=95
)
DPP10.where(DPP10 > 0).plot(yincrease=False, vmin=-0.1, vmax=0.6, cmap="viridis")
/Users/aaron.spring/anaconda3/envs/climpred-dev/lib/python3.9/site-packages/xarray/core/nputils.py:155: RuntimeWarning: Degrees of freedom <= 0 for slice.
result = getattr(npmodule, name)(values, axis=axis, **kwargs)
/Users/aaron.spring/anaconda3/envs/climpred-dev/lib/python3.9/site-packages/xarray/core/nputils.py:155: RuntimeWarning: Degrees of freedom <= 0 for slice.
result = getattr(npmodule, name)(values, axis=axis, **kwargs) Variance-Weighted Mean Period¶

climpred.stats.varweighted_mean_period() uses a periodogram based on a control simulation to extract the mean period of variations, which are weighted by the respective variance. Regions with a high mean period value indicate low-frequency variations with are potentially predictable .

vwmp = climpred.stats.varweighted_mean_period(control3d, dim="time")
threshold = climpred.bootstrap.varweighted_mean_period_threshold(
control3d, iterations=10
)
vwmp.where(vwmp > threshold).plot(yincrease=False, robust=True) Lag-1 Autocorrelation¶

The lag-1 autocorrelation also indicates where slower modes of variability occur by identifying regions with high temporal correlation .

from esmtools.stats import autocorr

# use climpred.bootstrap._bootstrap_func to wrap any stats function. esmtools.stats.autocorr computes the autocorrelation
# coefficient out to N lags. The first lag is at lag 0, so we select lead=1).
threshold = climpred.bootstrap._bootstrap_func(
autocorr, control3d, nlags=2, resample_dim="time", iterations=10
corr_ef.where(corr_ef > threshold).plot(yincrease=False, robust=False) Decorrelation time¶

Taking the lagged correlation further over all lags, climpred.stats.decorrelation_time() shows the time after which the autocorrelation fell beyond its e-folding .

threshold = climpred.bootstrap._bootstrap_func(
climpred.stats.decorrelation_time, control3d, "time", iterations=10
)
decorr_time = climpred.stats.decorrelation_time(control3d)
decorr_time.where(decorr_time > threshold).plot(yincrease=False, robust=False) Verify diagnostic potential predictability in initialized simulations¶

Do we find predictability in the areas highlighted above also in perfect-model experiments?

/Users/aaron.spring/Coding/climpred/climpred/utils.py:191: UserWarning: Assuming annual resolution starting Jan 1st due to numeric inits. Please change init to a datetime if it is another resolution. We recommend using xr.CFTimeIndex as init, see https://climpred.readthedocs.io/en/stable/setting-up-data.html.
warnings.warn(
/Users/aaron.spring/Coding/climpred/climpred/utils.py:191: UserWarning: Assuming annual resolution starting Jan 1st due to numeric inits. Please change init to a datetime if it is another resolution. We recommend using xr.CFTimeIndex as init, see https://climpred.readthedocs.io/en/stable/setting-up-data.html.
warnings.warn(
bootstrap_skill = pm.bootstrap(
metric="rmse",
comparison="m2e",
dim=["init", "member"],
reference=["uninitialized"],
iterations=10,
)
/Users/aaron.spring/Coding/climpred/climpred/checks.py:229: UserWarning: Consider chunking input ds along other dimensions than needed by algorithm, e.g. spatial dimensions, for parallelized performance increase.
warnings.warn(
init_skill = bootstrap_skill.sel(results="verify skill", skill="initialized")
# p value: probability that random uninitialized forecasts perform better than initialized
p = bootstrap_skill.sel(results="p", skill="uninitialized")
init_skill.where(p <= 0.05)[varname].plot( 