Documentation
TENAX class
- class pyTENAX.tenax.TENAX(return_period: list[int | float], durations: list[int], time_resolution: int, beta: float | int = 4, temp_time_hour: int = 24, alpha=0.05, n_monte_carlo=20000, tolerance=0.1, min_event_duration=30, storm_separation_time=24, left_censoring: list = [0, 1], niter_smev=100, niter_tenax=100, temp_res_monte_carlo=0.001, temp_delta=10, init_param_guess=[0.7, 0, 2, 0], min_rain: float | int = 0)[source]
Bases:
object- TNX_tenax_bootstrap_uncertainty(P, T, blocks_id, Ts, temp_method='norm', method_root_scalar='brentq', minimize_method='Nelder-Mead', parallel=False, n_jobs=-1)[source]
Bootstrap uncertainty estimation for the TENAX model.
- Parameters:
P (np.ndarray) – Precipitation ordinary events data.
T (np.ndarray) – Temperature ordinary events data.
blocks_id (np.ndarray) – Block identifiers (e.g., years) for each event.
Ts (np.ndarray) – Array of temperature values for the Monte Carlo integration.
temp_method (str, optional) – Distribution used for the temperature model. Defaults to
"norm".method_root_scalar (str, optional) – Root-finding method for model inversion. Defaults to
"brentq".minimize_method (str, optional) – Optimisation method passed to
magnitude_model. Defaults to"Nelder-Mead".parallel (bool, optional) – Run bootstrap iterations in parallel using
ProcessPoolExecutor. Warm start is disabled whenparallel=Truebecause iterations run in separate processes and cannot share state. Defaults toFalse.n_jobs (int, optional) – Number of worker processes.
-1uses all available cores. Ignored whenparallel=False. Defaults to-1.
Notes
Warm start (sequential mode only): each iteration reuses the fitted parameters from the previous iteration as the initial guess for
magnitude_model. Since consecutive bootstrap samples are drawn from the same dataset, their optimal parameters are typically close to each other, so the optimiser converges in fewer steps.self.init_param_guessis restored to its original value after the loop so that subsequent calls are not affected.- Returns:
F_phat_unc (np.ndarray) – Magnitude model parameters from each bootstrap sample, shape
(niter, 4).g_phat_unc (np.ndarray) – Temperature model parameters from each bootstrap sample, shape
(niter, 2)for"norm"or(niter, 3)for"skewnorm".RL_unc (np.ndarray) – Return levels from each bootstrap sample, shape
(niter, len(return_period)).n_unc (np.ndarray) – Mean number of events per block from each bootstrap sample, shape
(niter,).n_err (int) – Number of iterations where model fitting failed.
- associate_vars(dict_ordinary, data_temperature, dates_temperature, method='vectorized')[source]
Associate temperature with each ordinary event.
The associated temperature is the mean over the past
temp_time_hourhours before the event. Events for which no temperature can be found are dropped.- Parameters:
dict_ordinary (dict) – Dictionary of ordinary events as returned by get_ordinary_events_values.
data_temperature (np.ndarray) – Full temperature time series.
dates_temperature (np.ndarray) – Timestamps of the full temperature dataset.
method (str, optional) –
Backend used for the temperature association. Defaults to
"iterrows". One of:"iterrows"— pandas iterrows loop (original behaviour)."vectorized"— fully vectorised numpy using cumulative-sum trick; O(1) per event after one O(n) precomputation pass.
- Returns:
dict_ordinary (dict) – Ordinary events with an added
Tcolumn, keyed by duration. Example: ``{“10”: pd.DataFrame(columns=[‘year’, ‘oe_time’, ‘ordinary’, ‘T’])}``.
dict_dropped_oe (dict) – Ordinary events dropped because no temperature was found, keyed by duration.
n_ordinary_per_year_new (pd.DataFrame) – DataFrame with the count of ordinary events per year after dropping events without temperature.
- get_ordinary_events(data: DataFrame | ndarray, dates: ndarray, name_col: str = 'value', check_gaps=True) list[source]
Extract ordinary precipitation events from a time series.
Groups timesteps at or above
self.min_raininto independent storm events separated by at leastself.storm_separation_timehours. Optionally removes events too close to dataset boundaries or data gaps.- Parameters:
data (Union[pd.DataFrame, np.ndarray]) – Precipitation values.
dates (np.ndarray) – Timestamps of the precipitation data.
name_col (str, optional) – Column name to use when
datais a DataFrame. Defaults to “value”.check_gaps (bool, optional) – Remove events that fall within
storm_separation_timeof the dataset boundaries or internal data gaps. Defaults to True.
- Returns:
List of np.ndarray, each containing the timestamps of one ordinary event (values >=
self.min_rainseparated by more thanself.storm_separation_timehours).- Return type:
list
- get_ordinary_events_values(data: ndarray, dates: ndarray, arr_dates_oe: ndarray, method: str = 'vectorized') Tuple[Dict[str, DataFrame], Dict[str, DataFrame]][source]
Extract ordinary events and annual maxima from precipitation data.
- Parameters:
data (np.ndarray) – Full precipitation time series.
dates (np.ndarray) – Timestamps of the full precipitation dataset.
arr_dates_oe (np.ndarray) – End and start times of ordinary events as returned by remove_short.
method (str, optional) –
Backend used for the sliding-window maximum search. Defaults to
"vectorized". One of:"vectorized"— pure numpy,np.convolveper event."njit"— numba JIT-compiled loop, single-threaded."njit_parallel"— numba JIT-compiled loop, parallelised over events. Requiresnumbato be installed.
- Returns:
dict_ordinary (dict) – Key is duration (str), value is a
pd.DataFramewith columnsyear,oe_time,ordinary(event depth/intensity).dict_AMS (dict) – Key is duration (str), value is a
pd.DataFramewith columnsyearandAMS(annual maximum value).
- magnitude_model(data_oe_prec, data_oe_temp, thr, b_set=None, b_exp=False, minimize_method='Nelder-Mead')[source]
Fits the data to the magnitude model of TENAX.
- Parameters:
data_oe_prec (numpy.ndarray) – Array of precipitation ordinary events data.
data_oe_temp (numpy.ndarray) – Array of temperature ordinary events data.
thr (numpy.float64) – Magnitude of precipitation threshold.
b_set (NoneType or float) – Set value of b. fits magnitude model with a specified value for b.
b_exp (bool) – If True, uses the exponential rather than linear fit for b.
minimize_method (str, optional) – Optimisation method passed to
scipy.optimize.minimize. Defaults to'Nelder-Mead'. Tested also'L-BFGS-B', though it can lead to some instability
- Returns:
phat (np.ndarray) – Fitted parameters
[kappa_0, b, lambda_0, a].loglik (float) – Log-likelihood of the selected model.
loglik_H1 (float or None) – Log-likelihood of the alternative hypothesis (H1).
Nonewhenb_setis provided.loglik_H0shape (float or None) – Log-likelihood of the null hypothesis (constant shape).
Nonewhenb_setis provided.
- model_inversion(F_phat, g_phat, n, Ts, gen_P_mc=False, gen_RL=True, temp_method='norm', method_root_scalar='brentq', b_exp=False)[source]
Inversion of the TENAX model to predict return levels or plot model.
- Parameters:
F_phat (numpy.ndarray) – distribution values. F_phat = [kappa_0,b,lambda_0,a].
g_phat (numpy.ndarray) – [mu, sigma] of temperature distribution.
n (float) – Mean number of ordinary events per year.
Ts (numpy.ndarray) – Array of T values to use in the Monte Carlo.
gen_P_mc (bool, optional) – Specify whether to generate Monte Carlo values for precipitation. The default is False.
gen_RL (bool, optional) – Specify whether to generate return levels. The default is True.
temp_method (str, optional) – Type of fit used for the temperature model. The default is “norm”.
method_root_scalar (str, optional) – method used for inversion. The default is “brentq”.
b_exp (bool) – If True, uses the exponential rather than linear fit for b.
- Returns:
ret_lev (np.ndarray or list) – Return levels at periods specified in
self.return_period. Empty list[]ifgen_RL=False.T_mc (np.ndarray) – Monte Carlo generated temperature values, shape
(n_monte_carlo, 1).P_mc (np.ndarray or list) – Monte Carlo generated precipitation values. Empty list
[]ifgen_P_mc=False.
- remove_incomplete_years(data_pr: DataFrame, name_col='value', nan_to_zero=True) DataFrame[source]
Delete incomplete years in precipitation data.
An incomplete year is defined as a year where observations are missing above a given threshold.
- Parameters:
data_pr (pd.DataFrame) – Dataframe containing (hourly) precipitation values.
name_col (str, optional) – Column name in data_pr with precipitation values. Defaults to “value”.
nan_to_zero (bool, optional) – Set nan to zero. Defaults to True.
- Returns:
Dataframe containing (hourly) precipitation values with incomplete years removed.
- Return type:
pd.DataFrame
- remove_short(list_ordinary: list) Tuple[ndarray, ndarray, DataFrame][source]
Remove ordinary events that are too short.
- Parameters:
list_ordinary (list) – List of ordinary events as returned by get_ordinary_events(). Each event may contain pd.Timestamp or np.datetime64 values.
- Returns:
arr_vals (np.ndarray) – Boolean array (all True) of length equal to the number of kept events, one entry per event that passed the duration filter.
arr_dates (np.ndarray) – Array of (end, start) date tuples for each kept event.
n_ordinary_per_year (pd.DataFrame) – DataFrame with the count of ordinary events per year.
- temperature_model(data_oe_temp, beta=0, method='norm')[source]
Fit temperature data to the TENAX temperature model.
- Parameters:
data_oe_temp (np.ndarray) – Temperature data.
beta (float, optional) – Shape parameter of the generalised normal distribution. If 0, uses
self.beta. Defaults to 0.method (str, optional) – Distribution to fit.
"norm"uses the generalised normal;"skewnorm"uses a skewed normal. Defaults to"norm".
- Returns:
g_phat – Fitted parameters.
[mu, sigma]for"norm";[alpha, loc, scale]for"skewnorm".- Return type:
np.ndarray
Supporting utils
- pyTENAX.tenax.MC_tSMEV_cdf(y, wbl_phat, n)[source]
Evaluate the Monte Carlo SMEV CDF at given values.
- Parameters:
y (float or array-like) – Value(s) at which to evaluate the CDF.
wbl_phat (np.ndarray) – Weibull parameters, shape
(N, 2): columns are[scale, shape].n (float) – Power applied to the average probability.
- Returns:
CDF value(s) for input
y.- Return type:
np.ndarray
- pyTENAX.tenax.SMEV_Mc_inversion(wbl_phat, n, target_return_periods, vguess, method_root_scalar)[source]
Invert the MC-SMEV CDF to find quantiles for target return periods.
- Parameters:
wbl_phat (np.ndarray) – Weibull parameters, shape
(N, 2): columns are[scale, shape].n (float) – Power applied to the average probability.
target_return_periods (list or array-like) – Desired return periods.
vguess (np.ndarray) – Initial value grid for root-finding.
method_root_scalar (str) – Root-finding method passed to
scipy.optimize.root_scalar.
- Returns:
Quantiles corresponding to each target return period.
- Return type:
np.ndarray
- pyTENAX.tenax.gen_norm_loglik(x: ndarray, par: list, beta: float) float[source]
Compute the log-likelihood for the Generalized normal distribution.
- Parameters:
x (np.ndarray) – Data points.
par (list) – Parameters
[mu, sigma].beta (float) – Shape parameter.
- Returns:
Log-likelihood value.
- Return type:
float
- pyTENAX.tenax.gen_norm_pdf(x: ndarray, mu: float, sigma: float, beta: float) ndarray[source]
Compute the Generalized normal distribution PDF.
- Parameters:
x (np.ndarray) – Data points.
mu (float) – Location parameter.
sigma (float) – Scale parameter.
beta (float) – Shape parameter.
- Returns:
Generalized normal distribution PDF values.
- Return type:
np.ndarray
- pyTENAX.tenax.inverse_magnitude_model(F_phat, eT, qs, b_exp=False)[source]
Calculate percentiles from the Weibell magnitude model
- Parameters:
F_phat (numpy.ndarray) – distribution values. F_phat = [kappa_0,b,lambda_0,a].
eT (numpy.ndarray) – Temperature values from which to produce distribution.
qs (list) – list of percentiles to calculate (between 0 and 1). e.g. [0.85,0.95,0.99].
b_exp (bool) – If True, uses the exponential rather than linear fit for b.
- Returns:
percentile_lines – array with shape length(qs) by length(eT) giving the magnitudes for each eT. percentile_lines[0] are the values for qs[0].
- Return type:
numpy.ndarray
- pyTENAX.tenax.randdf(size, df, flag)[source]
Generate random numbers from a user-defined PDF or CDF.
Pythonised version of MATLAB’s randdf coded by halleyhit on Aug. 15th, 2018. Email: halleyhit@sjtu.edu.cn or halleyhit@163.com
- Parameters:
size (int or tuple) – Output size.
10→ 1-D array of 10;(10, 2)→ 10×2 matrix.df (np.ndarray) – 2-row matrix: first row is function values (PDF or CDF), second row is the corresponding sampling points.
flag (str) –
"pdf"or"cdf".
- Returns:
Random samples drawn according to the defined distribution.
- Return type:
np.ndarray
- pyTENAX.tenax.wbl_leftcensor_loglik(theta, t0, x1, t1, thr)[source]
Compute log-likelihood for a left-censored Weibull distribution.
Shape and scale parameters depend linearly on temperature. Observations below the threshold are left-censored.
- Parameters:
theta (array-like) – Parameter vector
[kappa_0, b, lambda_0, a].t0 (np.ndarray) – Temperature values for censored events (precipitation below
thr).x1 (np.ndarray) – Precipitation values at or above
thr.t1 (np.ndarray) – Temperature values corresponding to
x1.thr (float) – Left-censoring threshold.
- Returns:
Log-likelihood value.
- Return type:
float
- pyTENAX.tenax.wbl_leftcensor_loglik_H0shape(theta, t0, x1, t1, thr)[source]
Compute log-likelihood for a left-censored Weibull with constant shape (H0).
Same as wbl_leftcensor_loglik but with
b=0, i.e. shape does not depend on temperature. Used as the null hypothesis in the likelihood-ratio test.- Parameters:
theta (array-like) – Parameter vector
[kappa_0, b, lambda_0, a](bis ignored).t0 (np.ndarray) – Temperature values for censored events (precipitation below
thr).x1 (np.ndarray) – Precipitation values at or above
thr.t1 (np.ndarray) – Temperature values corresponding to
x1.thr (float) – Left-censoring threshold.
- Returns:
Log-likelihood value.
- Return type:
float
- pyTENAX.tenax.wbl_leftcensor_loglik_bset(theta, t0, x1, t1, thr, b_set)[source]
Compute log-likelihood for a left-censored Weibull with fixed b.
Same as wbl_leftcensor_loglik but
bis fixed tob_setand not optimised.- Parameters:
theta (array-like) – Parameter vector
[kappa_0, b, lambda_0, a](bis overridden byb_set).t0 (np.ndarray) – Temperature values for censored events (precipitation below
thr).x1 (np.ndarray) – Precipitation values at or above
thr.t1 (np.ndarray) – Temperature values corresponding to
x1.thr (float) – Left-censoring threshold.
b_set (float) – Fixed value of
b(shape temperature-dependence parameter).
- Returns:
Log-likelihood value.
- Return type:
float
Plotting utils
- pyTENAX.plotting.TNX_FIG_magn_model(P: ndarray, T: ndarray, F_phat: ndarray, thr: float, eT: ndarray, qs: list, obscol='r', valcol='b', xlimits: list = [-12, 30], ylimits: list = [0.1, 1000], b_exp=False) None[source]
Plots Figure 2a of Marra et al. (2024), the observed T-P pairs and the W model percentiles.
- Parameters:
P (np.ndarray) – Precipitation values.
T (np.ndarray) – Temperature values.
F_phat (np.ndarray) – Distribution values. F_phat = [kappa_0,b,lambda_0,a]..
thr (float) – Precipitation threshold for left-censorig.
eT (np.ndarray) – x values to plot W model.
qs (list) – Percentiles to calculate W.
obscol (str, optional) – Color code to plot observations. Defaults to “r”.
valcol (str, optional) – Color code to plot model. Defaults to “b”.
xlimits (list, optional) – x limits of plot. Defaults to [-12, 30].
ylimits (list, optional) – y limits of plot. Defaults to [0.1, 1000].
- pyTENAX.plotting.TNX_FIG_scaling(P, T, P_mc, T_mc, F_phat, niter_smev, eT, iTs, qs=[0.99], obscol='r', valcol='b', xlimits=[-15, 30], ylimits=[0.4, 1000])[source]
Plots figure 5.
- Parameters:
P (numpy.ndarray) – precipitation values
T (numpy.ndarray) – temperature values
P_mc (numpy.ndarray) – Monte Carlo generated precipitation values.
T_mc (numpy.ndarray) – Monte Carlo generated temperature values.
F_phat (numpy.ndarray) – distribution values. F_phat = [kappa_0,b,lambda_0,a].
niter_smev (int) – Number of iterations for uncertainty for the SMEV model .
eT (numpy.ndarray) – x (temperature) values to produce distribution for magnitude model.
iTs (numpy.ndarray) – x (temperature) values to produce distribution for quantile regression, binning, and TENAX.
qs (list) – percentiles to calculate W.
obscol (string, optional) – color code to plot observations. The default is ‘r’.
valcol (string, optional) – color code to plot model. The default is ‘b’.
xlimits (list, optional) – [min_x,max_x]. x limits to plot. The default is [-15,30].
ylimits (list, optional) – [min_y,max_y]. y limits to plot. The default is [0.4,1000].
- Returns:
scaling_rate – scaling rate of
- Return type:
float
- pyTENAX.plotting.TNX_FIG_temp_model(T, g_phat, beta, eT, obscol='r', valcol='b', obslabel='Observations', vallabel='Temperature model g(T)', xlimits=[-15, 30], ylimits=[0, 0.06], method='norm')[source]
Plots the observational and model temperature pdf
- Parameters:
T (numpy.ndarray) – Array of observed temperatures.
g_phat (numpy.ndarray) – [mu, sigma] of temperature distribution.
beta (float) – value of beta in generalised normal distribution.
eT (numpy.ndarray) – x (temperature) values to produce distribution.
obscol (string, optional) – color to plot observations. The default is ‘r’.
valcol (string, optional) – color to plot magnitude model. The default is ‘b’.
obslabel (string, optional) – Label for observations. The default is ‘observations’.
vallabel (string, optional) – Label for model plot. The default is ‘temperature model g(T)’.
xlimits (list, optional) – limits for the x axis [lower_x_limit, upper_x_limit]. The default is [-15,30].
ylimits (list, optional) – limits for the y axis [lower_y_limit, upper_y_limit]. The default is [0,0.06].
- Returns:
hist (numpy.ndarray) – pdf values of observed distribution.
pdf_values (numpy.ndarray) – pdf values of fitted model.
- pyTENAX.plotting.TNX_FIG_valid(AMS: DataFrame, RP: list, RL: ndarray, smev_RL: ndarray | list = [], RL_unc: ndarray | list = [], smev_RL_unc=0, TENAXcol='b', obscol_shape='g+', smev_colshape='--r', TENAXlabel='The TENAX model', obslabel='Observed annual maxima', smevlabel='The SMEV model', alpha=0.2, xlimits: list = [1, 200], ylimits: list = [0, 50]) None[source]
Plots figure 4 of Marra et al. (2024).
- Parameters:
AMS (pd.DataFrame) – Dataframe containing annual maxima.
RP (list) – Return periods to plot.
RL (np.ndarray) – Return levels calculated by TENAX.
smev_RL (Union[np.ndarray, list], optional) – Return levels calculated by SMEV. Defaults to [].
RL_unc (int, optional) – Uncertainty of return levels calculated by TENAX. Only relevant if smev_RL is provided. Defaults to 0.
smev_RL_unc (int, optional) – Uncertainty of return levels calculated by SMEV. Only relevant if smev_RL is provided. Defaults to 0.
TENAXcol (str, optional) – Linestyle for TENAX data to use in plot. Defaults to “b”.
obscol_shape (str, optional) – Linestyle for annual maxima data to use in plot. Defaults to “g+”.
smev_colshape (str, optional) – Linestyle for SMEV data to use in plot. Defaults to “–r”.
TENAXlabel (str, optional) – Label for TENAX data to use in plot. Defaults to “The TENAX model”.
obslabel (str, optional) – Label for annual maxima observation data to use in plot. Defaults to “Observed annual maxima”.
smevlabel (str, optional) – Label for SMEV data to use in plot. Defaults to “The SMEV model”.
alpha (float, optional) – Transparency to use in plot. Defaults to 0.2.
xlimits (list, optional) – x limits of plot. Defaults to [1, 200].
ylimits (list, optional) – y limits of plot. Defaults to [0, 50].
Supporting SMEV class
The SMEV implementation in pyTENAX builds upon the concepts from the pysmev library.
- class pyTENAX.smev.SMEV(return_period: list[int | float], durations: list[int], time_resolution: int, tolerance: float = 0.1, min_event_duration: int = 30, storm_separation_time: int = 24, left_censoring: list = [0, 1], min_rain: float | int = 0)[source]
Bases:
object- do_smev_all(dict_ordinary: Dict[str, DataFrame], n: float) Dict[str, dict][source]
Run SMEV parameter estimation and return level computation for all durations.
- Parameters:
dict_ordinary (Dict[str, pd.DataFrame]) – Dictionary of ordinary events per duration, as returned by get_ordinary_events_values.
n (float) – Mean number of ordinary events per year.
- Returns:
Keys are duration strings (e.g.
"10"). Each value is a dict with keys'SMEV_phat'(list[float]of length 2:[shape, scale]) and'RLs'(floatornp.ndarrayof return levels, one per return period).- Return type:
Dict[str, dict]
- estimate_smev_parameters(ordinary_events: ndarray | Series | list, data_portion: list[Tuple[int, float]]) list[float][source]
Estimate shape and scale parameters of the Weibull distribution.
- Parameters:
ordinary_events (np.ndarray or pd.Series or list) – Values of ordinary events.
data_portion (list) – Lower and upper limits of the probabilities of data to be used for the parameters estimation.
- Returns:
Shape and scale parameters of the Weibull distribution.
- Return type:
list[float]
- get_ordinary_events(data: DataFrame | ndarray, dates: ndarray, name_col: str = 'value', check_gaps=True) list[source]
Extract ordinary precipitation events from a time series.
Groups timesteps at or above
self.min_raininto independent storm events separated by at leastself.storm_separation_timehours. Optionally removes events too close to dataset boundaries or data gaps.- Parameters:
data (Union[pd.DataFrame, np.ndarray]) – Precipitation values.
dates (np.ndarray) – Timestamps of the precipitation data.
name_col (str, optional) – Column name to use when
datais a DataFrame. Defaults to “value”.check_gaps (bool, optional) – Remove events that fall within
storm_separation_timeof the dataset boundaries or internal data gaps. Defaults to True.
- Returns:
List of np.ndarray, each containing the timestamps of one ordinary event (values >=
self.min_rainseparated by more thanself.storm_separation_timehours).- Return type:
list
- get_ordinary_events_values(data: ndarray, dates: ndarray, arr_dates_oe: ndarray, method: str = 'vectorized') Tuple[Dict[str, DataFrame], Dict[str, DataFrame]][source]
Extract ordinary events and annual maxima from precipitation data.
- Parameters:
data (np.ndarray) – Full precipitation time series.
dates (np.ndarray) – Timestamps of the full precipitation dataset.
arr_dates_oe (np.ndarray) – End and start times of ordinary events as returned by remove_short.
method (str, optional) –
Backend used for the sliding-window maximum search. Defaults to
"vectorized". One of:"vectorized"— pure numpy,np.convolveper event."njit"— numba JIT-compiled loop, single-threaded."njit_parallel"— numba JIT-compiled loop, parallelised over events. Requiresnumbato be installed.
Notes
When using
method="njit"ormethod="njit_parallel", the first call in a Python session triggers JIT compilation and can take several seconds. Run a warmup call before any timed or production code:# Warmup — compile both kernels once at session start S_SMEV.get_ordinary_events_values( data=df_arr, dates=df_dates, arr_dates_oe=arr_dates, method="njit" ) S_SMEV.get_ordinary_events_values( data=df_arr, dates=df_dates, arr_dates_oe=arr_dates, method="njit_parallel" ) # Subsequent calls are fast dict_ordinary, dict_AMS = S_SMEV.get_ordinary_events_values( data=df_arr, dates=df_dates, arr_dates_oe=arr_dates, method="njit_parallel" )
- Returns:
dict_ordinary (dict) – Key is duration (str), value is a
pd.DataFramewith columnsyear,oe_time,ordinary(event depth/intensity). Example: ``{“10”: pd.DataFrame(columns=[‘year’, ‘oe_time’, ‘ordinary’])}``.
dict_AMS (dict) – Key is duration (str), value is a
pd.DataFramewith columnsyearandAMS(annual maximum value).
- static get_stats(df: DataFrame) Tuple[Series, Series, Series, Series][source]
Compute statistics of precipitation values.
Statistics are total precipitation per year, mean precipitation per year, standard deviation of precipitation per year, and count of precipitation events per year.
- Parameters:
df (pd.DataFrame) – Dataframe with precipitation values.
- Returns:
total_prec (pd.Series) – Total precipitation per year.
mean_prec (pd.Series) – Mean precipitation per year.
sd_prec (pd.Series) – Standard deviation of precipitation per year.
count_prec (pd.Series) – Count of precipitation events per year.
- remove_incomplete_years(data_pr: DataFrame, name_col='value', nan_to_zero=True) DataFrame[source]
Delete incomplete years in precipitation data.
An incomplete year is defined as a year where observations are missing above a given threshold.
- Parameters:
data_pr (pd.DataFrame) – Dataframe containing (hourly) precipitation values.
name_col (str, optional) – Column name in data_pr with precipitation values. Defaults to “value”.
nan_to_zero (bool, optional) – Set nan to zero. Defaults to True.
- Returns:
Dataframe containing (hourly) precipitation values with incomplete years removed.
- Return type:
pd.DataFrame
- remove_short(list_ordinary: list) Tuple[ndarray, ndarray, DataFrame][source]
Remove ordinary events that are too short.
- Parameters:
list_ordinary (list) – List of ordinary events as returned by get_ordinary_events(). Each event may contain pd.Timestamp or np.datetime64 values.
- Returns:
arr_vals (np.ndarray) – Boolean array (all True) of length equal to the number of kept events, one entry per event that passed the duration filter.
arr_dates (np.ndarray) – Array of (end, start) date tuples for each kept event.
n_ordinary_per_year (pd.DataFrame) – DataFrame with the count of ordinary events per year.
- smev_bootstrap_uncertainty(P: ndarray, blocks_id: ndarray, niter: int, n: float)[source]
Bootstrap uncertainty of SMEV return values.
- Parameters:
P (np.ndarray) – Array of precipitation data.
blocks_id (np.ndarray) – Array of block identifiers (e.g., years).
niter (int) – Number of bootstrap iterations.
n (float) – SMEV parameter n.
- Returns:
Array with bootstrapped return value uncertainty.
- Return type:
np.ndarray
- smev_return_values(return_period: int | float | list | ndarray, shape: float, scale: float, n: float) float | ndarray[source]
Calculate return values (rainfall intensity) from Weibull parameters.
- Parameters:
return_period (Union[int, float, list, np.ndarray]) – Return period(s) of interest. Scalar returns a float, array-like returns np.ndarray.
shape (float) – Shape parameter value.
scale (float) – Scale parameter value.
n (float) – SMEV parameter n.
- Returns:
Rainfall intensity value(s).
- Return type:
Union[float, np.ndarray]
Weibull tail test
Original code written in MATLAB is available at: https://zenodo.org/records/7234708 For now, users can refer to the file tests/test_weibull_test.py to see how to use this code.
Created on Wed Sep 10 15:53:57 2025
@author1: Yaniv yaniv.goldschmidt@unipd.it @author2: PetrVey
The test is described in: - Marra F, W Amponsah, SM Papalexiou, 2023. Non-asymptotic Weibull tails explain the statistics of extreme daily precipitation. Adv. Water Resour., 173, 104388, https://doi.org/10.1016/j.advwatres.2023.104388
Original code written in MATLAB is available at: https://zenodo.org/records/7234708
- pyTENAX.wbl_tail_test.check_confidence_interval(annual_max_indexes, records_df, p_confidence, annual_max, censor_value, p_out_dicts_lst)[source]
Arguments: - annual_max_indexes (list): List of indexes in the record of the annual/block maxima - records_df (dataframe): df with all the synthetic records - p_confidence (float): Probability to be used for the test. confidence interval = 1-p_confidence - annual_max (list): List of values over which the hypothesis is tested, i.e. block maxima - censor_value (float): The threshold for left censoring the record - p_out_dicts_lst (list): List of dicts - Each censor value tested gets a dict as follow: {censor_value:p_out}
Returns: - p_out_dicts_lst (list): Same list as in the arguments, after appending dict for the tested censor_value —————————————————————————–
- pyTENAX.wbl_tail_test.create_synthetic_records(seed_random: int, synthetic_records_amount: int, record_size: int, shape: float, scale: float) DataFrame[source]
The synthetic records contain random ordinary events sampled uniformly from the Weibull distribution. These synthetic records use as the basis for extracting the confidence interval.
- Parameters:
seed_random (-) – Value that determines the starting point for the pseudorandom number generator => due to reproducibility
synthetic_records_amount (-) – Value that determines how many synthetic records to generate => number of stochastic realizations
record_size (-) – The number of ordinary events in the record
shape (-) – Weibull distribution parameter
scale (-) – Weibull distribution parameter
- Returns:
pd.DataFrame – DataFrame with all the synthetic records. Each row represents separate synthetic record.
—————————————————————————–
- pyTENAX.wbl_tail_test.estimate_smev_param_without_AM(ordinary_events: ndarray | Series | list, censor_value, annual_max_indexes)[source]
Arguments: - ordinary_events ([np.ndarray, pd.Series, list): values of ordinary events - without zeros!!! - censor_value (float): The threshold for left censoring the record - annual_max_indexes (list): List of indexes in the record of the annual/block maxima, COMING ALREADY FROM SORTED ARRAY!
Returns: - shape, scale (floats): Weibull distribution parameters —————————————————————————–
- pyTENAX.wbl_tail_test.find_optimal_threshold(p_out_dicts_lst, p_confidence)[source]
If all threshold rejected - it will return 1
Arguments: - p_out_dicts_lst (list): List of dicts for each of the censor values tested, as follow: {censor_value:p_out} - p_confidence (float): Probability to be used for the test. confidence interval = 1-p_confidence
Returns: - optimal_threshold (float): The minimal threshold from which p_out <= p_confidence for all bigger thresholds.
If all threshold rejected - it will return 1. If not all thresholds rejected, (1-optimal_threshold) is the portion of the record that can be assumed to be distributed Weibull.
- range_of_optimal: list
List of thresholds where p_out < p_confidence.
- pyTENAX.wbl_tail_test.plot_curve(p_out_dicts_lst, p_confidence, optimal_threshold)[source]
Arguments: - p_out_dicts_lst (list): List of dicts for each of the censor values tested, as follow: {censor_value:p_out} - p_confidence (float): Probability to be used for the test. confidence interval = 1-p_confidence - optimal_threshold (float): The optimal left censoring threshold - csv_filename (str): Name of the input CSV file (without extension) to use for the output plot name
Returns: - Saves figure to monte_carlo/monte_carlo_output directory —————————————————————————–
- pyTENAX.wbl_tail_test.weibul_test_MC(ordinary_events_df: DataFrame, pr_field: str, hydro_year_field: str, seed_random: int = 42, synthetic_records_amount: int = 500, p_confidence: float = 0.1, make_plot: bool = True, censor_AM: bool = True, censor_values: ndarray | list = array([0., 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95]))[source]
This function will return the optimal left censoring threshold. If all threshold rejected - it will return 1. If not all thresholds rejected, (1-optimal_threshold) is the portion of the record that can be assumed to be distributed Weibull.
Warning: If the returned optimal_threshold == 1.0, this does not necessarily mean that all thresholds failed. It may indicate that your tested threshold range extended too high (e.g., up to 0.95 or higher), where results are dominated by stochastic sampling noise of extremes. In such cases, you should check whether some lower thresholds worked. As a rule of thumb, thresholds above ~0.95 are usually unreliable, though the exact cutoff is case dependent.
- Parameters:
ordinary_events_df (-) – One column pandas dataframe of the ordinary events - without zeros!!!
pr_field (-) – The name of the column with the precipitation values
hydro_year_field (-) – The name of the column with the hydrological-years / blocks values
seed_random (-) – Value that determines the starting point for the pseudorandom number generator => due to reproducibility
synthetic_records_amount (-) – Value that determines how many synthetic records to generate. => number of stochastic realizations
p_confidence (-) – Probability to be used for the test. confidence interval = 1-p_confidence
make_plot (-) – Choose whether or not to include the plot
censor_AM (-) – Choose whether or not the annual maximas should be included in ordinary events and test
censor_values_range (-) – The censoring thresholds which should be tested, nativally range from 0 to 1 in 0.05 step
- Returns:
- optimal_threshold (Union[float, int]) – The optimal left censoring threshold, or 1 if all rejected, or 1111 if there is a problem with Weibull parameters fit
- estimated_params (list) – Estimated weibull parameters of the optimal threshold (None if optimal==1)
- range_of_optimal (list) – List of thresholds where p_out < p_confidence.
- p_out_dicts_lst (list) – Fraction of block maxima outside of the Y = 1-p_out confidence interval
—————————————————————————–