Non-Stationary UGEV Modeling

• Non-stationary UGEV is a modeling framework where GEV parameters vary with covariates, capturing time- or space-varying extremes.
• It employs regression techniques like GAMs and linear models to link covariates to location and scale, with shape often held constant for stability.
• The method facilitates dynamic return level estimation and risk analysis in fields such as climate science and traffic safety.

Non-stationary univariate generalized extreme value (UGEV) modelling studies block maxima with a generalized extreme value law whose parameters are not fixed, but vary with time, space, or other covariates. In applied form, annual or seasonal maxima are represented as $Y(x)\sim \mathrm{GEV}(\mu(x), \psi(x), \xi(x))$, with cumulative distribution function

$$F_\text{GEV}(y; \mu, \psi, \xi) = \exp\left\{-\left[1 + \xi\left(\frac{y - \mu}{\psi}\right)\right]^{-1/\xi}\right\}$$

where $\mu$ is the location parameter, $\psi>0$ the scale parameter, and $\xi$ the shape parameter. In the broader extreme-value literature, the same topic also includes limit laws for maxima of non-stationary stochastic processes with non-identical marginals, non-uniform dependence, and possible clustering, so that the non-stationary setting is both a regression problem and a limit-theoretic one (Youngman, 2020, Freitas et al., 2015).

1. Core formulation and meaning of non-stationarity

The classical stationary setup assumes a single GEV distribution for all block maxima. The non-stationary formulation replaces constant parameters by functions of observed covariates. In the simplest notation, $\mu=\mu(x)$, $\psi=\psi(x)$, and $\xi=\xi(x)$, so that the distribution of extremes can change across season, elevation, longitude and latitude, traffic state, roadway class, or time. This allows the analysis to capture changes in the distribution of extremes across time, space, or other covariates.

In many environmental applications, the block-maxima interpretation is explicit. One study used maximum 24-hour precipitation during each period per year, analyzed separately for early summer and later summer, and treated location and scale as functions of time while keeping shape constant. The candidate trend structures for location and scale were “Linear (L),” “Parabolic (P),” and “No trend (N),” with model selection based on the minimum AIC. That formulation expresses a central feature of non-stationary UGEV: return levels are no longer constant but time-varying because $\mu(t)$ and $\sigma(t)$ evolve (Meng et al., 2023).

A recurring practical convention is to let the shape parameter remain constant while allowing location and scale to vary. The supplied studies give two reasons for this convention. First, the scale parameter is naturally handled on a log scale to preserve positivity. Second, covariate-varying shape is described as inherently hard to estimate, and fixing $$F_\text{GEV}(y; \mu, \psi, \xi) = \exp\left\{-\left[1 + \xi\left(\frac{y - \mu}{\psi}\right)\right]^{-1/\xi}\right\}$$0 aids stability in practice (Meng et al., 2023, André et al., 2023).

2. Non-stationary extreme value laws for stochastic processes

A second meaning of non-stationary UGEV concerns extreme value laws for processes $$F_\text{GEV}(y; \mu, \psi, \xi) = \exp\left\{-\left[1 + \xi\left(\frac{y - \mu}{\psi}\right)\right]^{-1/\xi}\right\}$$1 whose joint distributions are not invariant under time shift and whose marginals $$F_\text{GEV}(y; \mu, \psi, \xi) = \exp\left\{-\left[1 + \xi\left(\frac{y - \mu}{\psi}\right)\right]^{-1/\xi}\right\}$$2 may all differ. In that setting, the central object is the probability of no exceedance over a triangular array of thresholds $$F_\text{GEV}(y; \mu, \psi, \xi) = \exp\left\{-\left[1 + \xi\left(\frac{y - \mu}{\psi}\right)\right]^{-1/\xi}\right\}$$3. The basic exceedance-balance assumption is

$$F_\text{GEV}(y; \mu, \psi, \xi) = \exp\left\{-\left[1 + \xi\left(\frac{y - \mu}{\psi}\right)\right]^{-1/\xi}\right\}$$4

for some $$F_\text{GEV}(y; \mu, \psi, \xi) = \exp\left\{-\left[1 + \xi\left(\frac{y - \mu}{\psi}\right)\right]^{-1/\xi}\right\}$$5.

The main theorem in this framework states that if the sum of exceedance probabilities converges, the dependence condition $$F_\text{GEV}(y; \mu, \psi, \xi) = \exp\left\{-\left[1 + \xi\left(\frac{y - \mu}{\psi}\right)\right]^{-1/\xi}\right\}$$6 and its primed version $$F_\text{GEV}(y; \mu, \psi, \xi) = \exp\left\{-\left[1 + \xi\left(\frac{y - \mu}{\psi}\right)\right]^{-1/\xi}\right\}$$7 hold, and an extremal index $$F_\text{GEV}(y; \mu, \psi, \xi) = \exp\left\{-\left[1 + \xi\left(\frac{y - \mu}{\psi}\right)\right]^{-1/\xi}\right\}$$8 exists, then

$$F_\text{GEV}(y; \mu, \psi, \xi) = \exp\left\{-\left[1 + \xi\left(\frac{y - \mu}{\psi}\right)\right]^{-1/\xi}\right\}$$9

Here $\mu$0 measures clustering: $\mu$1 corresponds to effectively independent extremes, whereas $\mu$2 indicates clustering. The theory replaces stronger stationary-style mixing assumptions by weaker conditions tailored to non-autonomous dynamical systems, sequential dynamical systems, and random dynamical systems (Freitas et al., 2015).

This framework yields phenomena absent in the stationary case. For sequential $\mu$3-transformations, if $\mu$4 is periodic then $\mu$5, where $\mu$6 is the prime period, whereas if $\mu$7 is not periodic then $\mu$8. The same paper also reports that if the sequence $\mu$9 approaches $\psi>0$0 slowly, the extremal index at periodic points may anomalously be $\psi>0$1, highlighting a fundamental difference from stationary theory (Freitas et al., 2015).

3. Regression structures for non-stationary UGEV

A dominant modern formulation uses generalized additive models (GAMs). In this approach, each parameter is linked to a linear predictor built from basis expansions: $\psi>0$2 with

$\psi>0$3

This produces additive smooth, possibly non-linear dependence on covariates and supports spatial smooths, cyclic splines for seasonality, and tensor-product interactions. The R package evgam was designed precisely to let extreme value distribution parameters have generalized additive model forms, estimated objectively using Laplace’s method, and the paper illustrates both spatial models for daily precipitation accumulations in Colorado and temporal models for daily maximum temperatures in Fort Collins (Youngman, 2020).

A simpler but widely used specification treats parameters as linear or log-linear in covariates. The Python package nsEVDx supports

$\psi>0$4

$\psi>0$5

$\psi>0$6

The configuration vector controls which parameters are stationary and which are covariate-dependent, making it possible to fit, for example, a model with non-stationary location and stationary scale and shape (Kafle et al., 8 Sep 2025).

Other formulations encode non-stationarity through change-points or hierarchical structure. In non-stationary wireless channel modelling, parameters are piecewise constant over groups determined by measurable external factors, written as

$\psi>0$7

In corridor-level crash occurrence risk, a Hierarchical Bayesian spatiotemporal grouped random parameter framework embeds a non-stationary UGEV in which fixed effects depend on group-level and sample-level covariates and selected coefficients vary by group, with partial pooling across locations to address unobserved heterogeneity (Mehrnia et al., 2021, Anis et al., 2 Sep 2025).

A further Bayesian development reparametrizes the GEV in terms of an $\psi>0$8-quantile $\psi>0$9, a quantile range $\xi$0, and shape $\xi$1, so that regression targets remain interpretable even when mean or variance may not exist. This reparametrization was introduced together with a blended GEV distribution that combines the left tail of a Gumbel distribution with the right tail of a Fréchet distribution to obtain unbounded support (Castro-Camilo et al., 2021).

4. Estimation, smoothing, uncertainty, and model selection

For GAM-based non-stationary UGEV, estimation is based on a penalized log-likelihood,

$\xi$2

where $\xi$3 is a penalty matrix indexed by smoothing parameters. A central issue is smoothing-parameter selection. In evgam, this is handled by REML, integrating out coefficients using Laplace’s method; model fitting proceeds through inner iterations for $\xi$4 and outer iterations for $\xi$5. The approach is described as closely linked with the methodology of Wood (2011, 2016) for smoothing parameter estimation in GAMs beyond exponential families (Youngman, 2020).

Frequentist and Bayesian estimation coexist in current software and methodological work. nsEVDx uses maximum likelihood via scipy.optimize and also supports Metropolis-Hastings, Metropolis-adjusted Langevin Algorithm (MALA), and Hamiltonian Monte Carlo (HMC). In the grouped random-parameter crash-risk model, the full hierarchical model is estimated via MCMC in MultiBUGS, with model comparison by DIC, WAIC, and LOOIC, and predictive validation by ROC-AUC. In environmental challenge settings, penalized likelihood was combined with cross-validation, CRPS, AIC, and BIC, and a stationary, semiparametric bootstrap based on the stationary bootstrap was used to account for temporal dependence in covariates and their relationship with the response (Kafle et al., 8 Sep 2025, Anis et al., 2 Sep 2025, André et al., 2023).

Robust alternatives to likelihood-based estimation have been developed for non-stationary GEV models using L-moments. One line of work combines L-moments and robust regression, transforms non-stationary data to standardized residuals that should follow a standard Gumbel distribution, and solves for parameters by matching the first three L-moments of the transformed data to the standard Gumbel values $\xi$6, $\xi$7, and $\xi$8. Another development, the generalized method of L-moment estimation (GLME), defines a generalized L-moment distance

$\xi$9

interprets it through a multivariate normal likelihood approximation, and augments it with penalty functions or prior information, especially on the shape parameter. The same literature proposes cross-validated generalized L-moment distance as an alternative to likelihood-based criteria for model selection (Shin et al., 1 Jun 2025, Shin et al., 23 Dec 2025).

Uncertainty quantification is similarly plural. evgam provides delta-method-based standard errors for predicted parameters or return levels and simulation from the approximate posterior distribution of the parameters. Challenge-oriented POT/GAM work used bootstrap-based confidence intervals, while Bayesian implementations rely on posterior predictive exceedance probabilities and posterior visualization (Youngman, 2020, André et al., 2023, Kafle et al., 8 Sep 2025).

5. Return levels, diagnostics, and empirical uses

In non-stationary UGEV, the $\mu=\mu(x)$0-year return level is itself time-varying. A standard formula is

$\mu=\mu(x)$1

for $\mu=\mu(x)$2. Because $\mu=\mu(x)$3 and $\mu=\mu(x)$4 change with time, return levels also change with time, reflecting adaptation to climate change or other driving factors. A related proposal introduces a redefined return level under nonstationarity, obtained by solving

$\mu=\mu(x)$5

which emphasizes that stationary return-period interpretations are altered when the marginal law evolves (Meng et al., 2023, Shin et al., 1 Jun 2025).

Diagnostics are correspondingly varied. GAM-based implementations examine summaries and fitted smooths; challenge settings relied on QQ-plots for transformed exceedances and coverage measurements; non-stationary threshold models employed mean residual life plots, parameter stability plots, PP-plots, and QQ-plots. The shared objective is to determine whether the fitted model adequately captures the tail while controlling model complexity and preserving extrapolative credibility (Youngman, 2020, André et al., 2023, Mehrnia et al., 2021).

Empirical studies show why non-stationary UGEV is used. In Northeast China, analysis of 107 stations during 1959–2017 found that negative trends dominate extreme precipitation in early summer but positive trends prevail in later summer, and that early and later summer should be discussed separately. The same study reported that in early summer all return levels (2-, 20-, 50-, 100-year) show significant increasing trends, while in later summer the 2-year return level declines slightly and the 20-, 50-, and 100-year return levels show only slight increases. It also reported that normal EP events more frequent in Liaoning and extreme events more likely Jilin/Heilongjiang (Meng et al., 2023).

In traffic safety, the HBSGRP-UGEV framework was applied to geometry-aware two-dimensional time-to-collision indicators extracted from the Argoverse-2 dataset at 28 locations along Miami’s Biscayne Boulevard. Results show that the HBSGRP-UGEV framework outperforms fixed-parameter alternatives, reducing DIC by up to 7.5% for V-V and 3.1% for V-I near-misses, with ROC-AUC values of 0.89 for V-V segments, 0.82 for V-V intersections, 0.79 for V-I segments, and 0.75 for V-I intersections (Anis et al., 2 Sep 2025).

The Lancopula Utopiversity team used a flexible modelling technique, based on generalized additive models, for non-stationary univariate time series in the EVA (2023) Conference Data Challenge. Their diagnostics indicated generally good performance for the observed data, while also documenting that out-of-sample coverage of true quantiles was fair but imperfect and that confidence intervals may be somewhat narrow (André et al., 2023).

6. Limitations, misconceptions, and related directions

A persistent misconception is that non-stationarity is fully handled by adding a linear time trend to the location parameter. The supplied literature suggests a broader picture. Non-stationarity may reside in location, scale, threshold, or dependence; may be smooth, linear, parabolic, cyclic, piecewise constant, or hierarchical; and may arise from time, space, environmental factors, or jointly from roadway geometry and dynamic variables. A plausible implication is that model form is part of the scientific question rather than a purely technical choice (Youngman, 2020, Meng et al., 2023, Anis et al., 2 Sep 2025).

Another misconception is that allowing every GEV parameter to vary is always preferable. Several papers state that the shape parameter is difficult to estimate reliably when covariate-dependent, so it is often held constant for tractability and robustness. Flexible semi-parametric models can overfit when too much flexibility is allowed, high $\mu=\mu(x)$6 values can be problematic, and predictions outside the range of the data can be unreliable unless the underlying trend is well supported. In addition, coverage of confidence intervals can be lower than ideal because of over-fitting, underestimation of uncertainty by the bootstrap, or misspecification (Youngman, 2020, André et al., 2023, Shin et al., 1 Jun 2025).

A separate issue concerns the support of the classical GEV in regression settings. When $\mu=\mu(x)$7, the standard GEV has a finite lower endpoint, and with covariates this endpoint becomes covariate-dependent. The blended GEV was introduced specifically to address this usually overlooked issue by smoothly combining the left tail of a Gumbel distribution with the right tail of a Fréchet distribution, thereby producing unbounded support. The same work introduced property-preserving penalized complexity priors for the shape parameter, restricting $\mu=\mu(x)$8 to preserve first and second moments (Castro-Camilo et al., 2021).

Related developments extend beyond univariate marginal regression while remaining closely connected to the non-stationary UGEV problem. One approach models non-stationarity in the tail through a scedasis function $\mu=\mu(x)$9, defined by

$\psi=\psi(x)$0

and develops tests for spatial homogeneity and temporal homogeneity together with pseudo-maximum likelihood estimation of a common extreme value index under pooling of non-stationary and dependent observations. Another line of work treats cases where both response and covariates are non-stationary block maxima, using a transformed-stationary approach and a regression manifold derived from a bivariate extreme value copula rather than a univariate GEV margin (Einmahl et al., 2020, Bernoussi et al., 25 Jun 2025).

Taken together, these developments show that non-stationary UGEV is not a single model but a family of formulations for evolving extremal behaviour. Its modern practice spans limit theory for non-stationary processes, covariate-dependent GEV regression, smooth and hierarchical parameterizations, penalized and Bayesian inference, robust L-moment alternatives, and application-specific diagnostics directed at time-varying risk.