| Title: | Generalized L-Moments Estimation for Extreme Value Distributions |
|---|---|
| Description: | Provides generalized L-moments estimation methods for the generalized extreme value ('GEV') distribution. Implements both stationary 'GEV' and non-stationary 'GEV11' models where location and scale parameters vary with time. Includes various penalty functions (Martins-Stedinger, Park, Cannon, Coles-Dixon) for shape parameter regularization. Also provides model averaging estimation ('ma.gev') that combines MLE and L-moment methods with multiple weighting schemes for robust high quantile estimation. The 'GLME' methodology is described in Shin et al. (2025a) <doi:10.48550/arXiv.2512.20385>. The non-stationary L-moment method is based on Shin et al. (2025b) <doi:10.1007/s42952-025-00325-3>. The model averaging method is described in Shin et al. (2026) <doi:10.1007/s00477-025-03167-x>. See also Hosking (1990) <doi:10.1111/j.2517-6161.1990.tb01775.x> for L-moments theory and Martins and Stedinger (2000) <doi:10.1029/1999WR900330> for penalized likelihood methods. |
| Authors: | Yonggwan Shin [aut, cre] (ORCID: <https://orcid.org/0000-0001-6966-6511>), Seokkap Ko [aut, ctb] (ORCID: <https://orcid.org/0000-0001-9218-7236>), Jihong Park [ctb] (ORCID: <https://orcid.org/0009-0003-3191-9968>), Yire Shin [aut, dtc] (ORCID: <https://orcid.org/0000-0003-1297-5430>), Jeong-Soo Park [aut, ths] (ORCID: <https://orcid.org/0000-0002-8460-4869>) |
| Maintainer: | Yonggwan Shin <[email protected]> |
| License: | GPL (>= 3) |
| Version: | 2.0.1 |
| Built: | 2026-10-26 07:05:29 UTC |
| Source: | https://github.com/sygstat/glmom |
Annual maximum daily rainfall data from Bangkok, Thailand. This dataset is used for demonstrating model averaging methods for high quantile estimation in extreme value analysis.
bangkokbangkok
A data frame with 58 rows and 5 columns:
Annual maximum daily rainfall in mm (numeric)
2nd largest annual daily rainfall in mm (numeric)
3rd largest annual daily rainfall in mm (numeric)
4th largest annual daily rainfall in mm (numeric)
5th largest annual daily rainfall in mm (numeric)
Thai Meteorological Department (TMD; https://www.tmd.go.th)
Shin, Y., Shin, Y., & Park, J. S. (2026). Model averaging with mixed criteria for estimating high quantiles of extreme values: Application to heavy rainfall. Stochastic Environmental Research and Risk Assessment, 40(2), 47. doi:10.1007/s00477-025-03167-x
data(bangkok) head(bangkok) # Estimate high quantiles using model averaging result <- ma.gev(bangkok$X1, quant = c(0.99, 0.995)) print(result$qua.ma)data(bangkok) head(bangkok) # Estimate high quantiles using model averaging result <- ma.gev(bangkok$X1, quant = c(0.99, 0.995)) print(result$qua.ma)
Estimates parameters of a non-stationary GEV distribution using multiple methods: Weighted Least Squares (WLS), GN16 method, and the proposed L-moment method from Shin et al. (2025, J. Korean Stat. Soc.).
As of GLmom v2.0.0 this function is deprecated: it was renamed to
lme.gev11 (for the proposed method), and the auxiliary
methods are available as strup.gev11 and
GN16.gev11. This wrapper reassembles the v1.x output
format from those functions and will be removed in a future release.
gado.prop_11(xdat, ntry = 20, ftol = 1e-06)gado.prop_11(xdat, ntry = 20, ftol = 1e-06)
xdat |
A numeric vector of data to be fitted. |
ntry |
Number of attempts for optimization (default 20). |
ftol |
Function tolerance for optimization (default 1e-6). |
A list containing:
para.prop - L-moment based estimates (proposed method).
para.gado - GN16 method estimates.
para.wls - Weighted least squares estimates.
strup.org - Original non-stationary WLSE by Strup method.
lme.sta - Stationary L-moment estimates.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y. & Park, J.-S. (2025). Building nonstationary extreme value model using L-moments. Journal of the Korean Statistical Society, 54, 947-970. doi:10.1007/s42952-025-00325-3
lme.gev11 for the recommended interface,
strup.gev11, GN16.gev11,
glme.gev11 for the full GLME method,
nsgev for the simple interface.
data(PhliuAgromet) # deprecated wrapper; suppressWarnings() silences the deprecation notice result <- suppressWarnings(gado.prop_11(PhliuAgromet$prec, ntry = 5)) print(result$para.prop) print(result$lme.sta)data(PhliuAgromet) # deprecated wrapper; suppressWarnings() silences the deprecation notice result <- suppressWarnings(gado.prop_11(PhliuAgromet$prec, ntry = 5)) print(result$para.prop) print(result$lme.sta)
This function estimates the parameters of the Generalized Extreme Value distribution by the generalized method of L-moment estimation (GLME), which minimizes the generalized L-moment distance plus a penalty (preference) function on the shape parameter.
glme.gev( xdat, ntry = 5, pen = "beta", pen.choice = 1, mu = -0.5, std = 0.2, p = 6, c1 = 3, c2 = 1, c0 = 0.35, q = NULL, show = FALSE, method = "BFGS", maxit = 70, abstol = 1e-05 )glme.gev( xdat, ntry = 5, pen = "beta", pen.choice = 1, mu = -0.5, std = 0.2, p = 6, c1 = 3, c2 = 1, c0 = 0.35, q = NULL, show = FALSE, method = "BFGS", maxit = 70, abstol = 1e-05 )
xdat |
A numeric vector of data to be fitted. |
ntry |
Number of attempts for parameter estimation (default 5). Higher values increase the chance of finding the global optimum by trying different initial conditions. |
pen |
Type of penalty function: Choose among "norm", "beta" (default), "ms", "park", "cannon", "cd", and "no" (without penalty function). |
pen.choice |
Choice number of penalty function specifying hyperparameters
(default 1). For "beta": 1-6 correspond to (p, c1, c2) =
(6,3,1), (6,5,2), (6,7,3), (2,3,0.5), (2,5,1), (2,7,1.5).
For "norm": 1-4 correspond to (mu, std) =
(-0.5,0.25), (-0.5,0.15), (-0.6,0.25), (-0.6,0.15).
Set |
mu |
Mean hyperparameter for "norm" penalty function (default -0.5). |
std |
Standard deviation hyperparameter for "norm" penalty function (default 0.2). |
p |
Shape hyperparameter for "beta" penalty function (default 6). |
c1 |
Scaling hyperparameter for "beta" penalty function (default 3). |
c2 |
Upper limit hyperparameter for "beta" penalty function (default 1). |
c0 |
Half-width of the adaptive beta penalty support (default 0.35). |
q |
Optional fixed second shape parameter for the beta penalty. If given, the beta penalty support is fixed instead of data-adaptive. |
show |
If TRUE, prints the objective value and parameters of each try. |
method |
Optimization method passed to |
maxit |
Maximum number of iterations for |
abstol |
Absolute and relative convergence tolerance for |
The equations for the L-moments for LME of the GEVD are
where and .
Next, we define the generalized L-moments distance (GLD) as;
where is the variance-covariance matrix of the sample L-moments up to the third order.
The glme.gev function returns a list containing the following elements:
para.glme - The estimated parameters of the Generalized Extreme Value distribution.
para.lme - The L-moment estimates of the parameters.
nllh.glme - The penalized negative log-likelihood of the GLME solution.
convergence - 0 if the optimization converged, 5 otherwise.
pen, pen_pen.choice - The penalization method (and choice) used.
p_q - (for beta penalty) The p and q values used.
c1_c2, c0_c1_c2 - (for beta penalty) The hyperparameters used.
mu_std - (for norm penalty) The mu and std values used.
covinv.lmom - The inverse of the covariance matrix of the L-moments.
lcovdet - The log determinant of the covariance matrix.
data - The input data (used by the plot method).
The object has class "glme"; see GLmom-methods.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y., Park, J. & Park, J.-S. (2025). Generalized method of L-moment estimation for stationary and nonstationary extreme value models. arXiv preprint arXiv:2512.20385. doi:10.48550/arXiv.2512.20385
glme.gev11 for non-stationary GEV estimation,
ma.gev for model averaging estimation,
glme.like for the objective function,
quagev.NS for quantile computation.
# Load example heavy-rainfall data data(haenam) x <- haenam$X1 # Estimate GEV parameters using beta penalty (default) result <- glme.gev(x, ntry = 5) print(result$para.glme) # Using Martins-Stedinger penalty result_ms <- glme.gev(x, ntry = 5, pen = "ms") print(result_ms$para.glme)# Load example heavy-rainfall data data(haenam) x <- haenam$X1 # Estimate GEV parameters using beta penalty (default) result <- glme.gev(x, ntry = 5) print(result$para.glme) # Using Martins-Stedinger penalty result_ms <- glme.gev(x, ntry = 5, pen = "ms") print(result_ms$para.glme)
This function estimates parameters of the non-stationary GEV11 model where mu(t) = mu0 + mu1*t and sigma(t) = exp(sigma0 + sigma1*t), by the generalized method of L-moment estimation (GLME).
Starting from the WLS pre-estimate of strup.gev11,
the penalized objective (generalized L-moment distance based on the
asymptotic covariance of the sample L-moments, plus penalty) is
minimized over (mu0, sigma0, xi) with mu1 and sigma1 held at their
pre-estimated values. With pen="no" the result equals the pure
L-moment method lme.gev11 of Shin et al. (2025b).
glme.gev11( xdat, ntry = 5, ftol = 1e-06, opt.choose = "nllh", pen = "beta", pen.choice = 1, mu = -0.55, std = 0.2, p = 6, c1 = 5, c2 = 2, c0 = 0.35, q = NULL, show = FALSE, init.rob = TRUE, glme.pre = NULL )glme.gev11( xdat, ntry = 5, ftol = 1e-06, opt.choose = "nllh", pen = "beta", pen.choice = 1, mu = -0.55, std = 0.2, p = 6, c1 = 5, c2 = 2, c0 = 0.35, q = NULL, show = FALSE, init.rob = TRUE, glme.pre = NULL )
xdat |
A numeric vector of data to be fitted. |
ntry |
Number of attempts for parameter estimation (default 5). |
ftol |
Tolerance for convergence (default 1e-6). |
opt.choose |
Selection criterion among the multi-start solutions: "nllh" (default, penalized negative log-likelihood) or "gof" (goodness-of-fit by the number-of-exceedances measure). |
pen |
Type of penalty function: "norm", "beta" (default), "ms", "park", "cannon", "cd", or "no". |
pen.choice |
Choice number for penalty hyperparameters (default 1).
See |
mu |
Mean for normal penalty (default -0.55). |
std |
Std for normal penalty (default 0.2). |
p |
Shape for beta penalty (default 6). |
c1 |
Scaling for beta penalty (default 5). |
c2 |
Limit for beta penalty (default 2). |
c0 |
Half-width of the adaptive beta penalty support (default 0.35). |
q |
Optional fixed second shape parameter for the beta penalty. |
show |
If TRUE, prints the objective value of each try. |
init.rob |
Use robust regression for the WLS initialization (default TRUE). |
glme.pre |
Deprecated (ignored). The WLS pre-estimate is always used. |
A list containing:
para.glme - Proposed GLME estimates (mu0, mu1, sigma0, sigma1, xi).
nllh.glme - Penalized negative log-likelihood of the solution.
convergence - 0 if converged, 5 otherwise.
pen, pen_pen.choice - Penalty method (and choice) used.
p_q, c1_c2, c0_c1_c2 - (for beta) hyperparameters used.
mu_std - (for norm) hyperparameters used.
strup.org - WLSE by the Strupczewski method (strup.para).
para.wls - Modified WLS estimates (strup.mdfy).
para.gado - GN16 method estimates.
lme.sta - Stationary L-moment estimates.
para.lme - Pure L-moment estimates; returned only when
pen="no" (then identical to para.glme). For other penalties
use lme.gev11 directly.
precis - (only when pen="no") precision of the L-moment
equation solution.
data - The input data (used by the plot method).
The object has class "glme11"; see GLmom-methods.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y., Park, J. & Park, J.-S. (2025). Generalized method of L-moment estimation for stationary and nonstationary extreme value models. arXiv preprint arXiv:2512.20385. doi:10.48550/arXiv.2512.20385
Shin, Y., Shin, Y. & Park, J.-S. (2025). Building nonstationary extreme value model using L-moments. Journal of the Korean Statistical Society, 54, 947-970. doi:10.1007/s42952-025-00325-3
glme.gev for stationary GEV estimation,
lme.gev11 for the pure L-moment method (no penalty),
strup.gev11 for the WLS pre-estimation,
GN16.gev11 for the GN16 method,
quagev.NS for non-stationary quantile computation.
data(PhliuAgromet) x <- PhliuAgromet$prec # Estimate non-stationary GEV11 parameters result <- glme.gev11(x, ntry = 5) print(result$para.glme) # Proposed GLME estimatesdata(PhliuAgromet) x <- PhliuAgromet$prec # Estimate non-stationary GEV11 parameters result <- glme.gev11(x, ntry = 5) print(result$para.glme) # Proposed GLME estimates
This function calculates the likelihood (or more precisely, a penalized negative log-likelihood) for the Generalized L-moments estimation of the Generalized Extreme Value (GEV) distribution.
glme.like( par, xdat = xdat, slmgev = slmgev, covinv = covinv, lcovdet = lcovdet, mu = mu, std = std, lme = lme, pen = pen, p = p, c1 = c1, c2 = c2, c0 = 0.35, q = NULL )glme.like( par, xdat = xdat, slmgev = slmgev, covinv = covinv, lcovdet = lcovdet, mu = mu, std = std, lme = lme, pen = pen, p = p, c1 = c1, c2 = c2, c0 = 0.35, q = NULL )
par |
A vector of GEV parameters (location, scale, shape). |
xdat |
A numeric vector of data. |
slmgev |
Sample L-moments of the data. |
covinv |
Inverse of the covariance matrix of the sample L-moments. |
lcovdet |
Log determinant of the covariance matrix. |
mu |
Mean for the normal penalization (used when pen='norm'). |
std |
Standard deviation for the normal penalization (used when pen='norm'). |
lme |
L-moment estimates of the parameters. |
pen |
Penalization method: 'norm', 'beta', 'ms', 'park', 'cannon', 'cd', or 'no'. |
p |
Shape parameter for beta penalty. |
c1 |
Scaling parameter for beta penalty. |
c2 |
Upper limit parameter for beta penalty. |
c0 |
Half-width of the adaptive beta penalty support (default 0.35). |
q |
Optional fixed second shape parameter for the beta penalty. |
The function performs the following steps: 1. Checks if the parameters are within valid ranges. 2. Calculates the expected L-moments based on the current parameters. 3. Computes the difference between expected and sample L-moments. 4. Calculates the generalized L-moments distance. 5. Applies a penalization term based on the specified method. 6. Returns the sum of the L-moments distance and the penalization term.
A numeric value representing the penalized negative log-likelihood. A lower value indicates a better fit.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y., Park, J. & Park, J.-S. (2025). Generalized method of L-moment estimation for stationary and nonstationary extreme value models. arXiv preprint arXiv:2512.20385. doi:10.48550/arXiv.2512.20385
glme.gev which calls this function for optimization.
data(haenam) x <- haenam$X1 slm <- lmomco::lmoms(x, nmom = 3) cov_mat <- lmomco::lmoms.cov(x, nmom = 3) lme_par <- lmomco::pargev(slm)$para glme.like(par = lme_par, xdat = x, slmgev = slm, covinv = solve(cov_mat), lcovdet = log(det(cov_mat)), mu = -0.5, std = 0.2, lme = lme_par, pen = "beta", p = 6, c1 = 3, c2 = 1)data(haenam) x <- haenam$X1 slm <- lmomco::lmoms(x, nmom = 3) cov_mat <- lmomco::lmoms.cov(x, nmom = 3) lme_par <- lmomco::pargev(slm)$para glme.like(par = lme_par, xdat = x, slmgev = slm, covinv = solve(cov_mat), lcovdet = log(det(cov_mat)), mu = -0.5, std = 0.2, lme = lme_par, pen = "beta", p = 6, c1 = 3, c2 = 1)
The estimation functions return classed objects: glme.gev
returns an object of class "glme", glme.gev11 of
class "glme11", lme.gev11 of class "lme11",
and ma.gev of class "magev". The objects remain
ordinary lists, so all documented fields (e.g., $para.glme) stay
accessible as before.
print shows the parameter estimates and the penalty settings;
summary additionally reports the auxiliary estimates stored in
the object; plot draws a quantile-quantile diagnostic: against
the fitted GEV distribution for "glme", on the standard Gumbel
scale (after the transformation of Eq. 8 in the manuscript) for
"glme11" and "lme11", and via magev.qqplot
for "magev".
## S3 method for class 'glme' print(x, ...) ## S3 method for class 'glme' summary(object, ...) ## S3 method for class 'glme' plot(x, ...) ## S3 method for class 'glme11' print(x, ...) ## S3 method for class 'glme11' summary(object, ...) ## S3 method for class 'glme11' plot(x, ...) ## S3 method for class 'lme11' print(x, ...) ## S3 method for class 'lme11' summary(object, ...) ## S3 method for class 'lme11' plot(x, ...) ## S3 method for class 'magev' print(x, ...) ## S3 method for class 'magev' summary(object, ...) ## S3 method for class 'magev' plot(x, ...)## S3 method for class 'glme' print(x, ...) ## S3 method for class 'glme' summary(object, ...) ## S3 method for class 'glme' plot(x, ...) ## S3 method for class 'glme11' print(x, ...) ## S3 method for class 'glme11' summary(object, ...) ## S3 method for class 'glme11' plot(x, ...) ## S3 method for class 'lme11' print(x, ...) ## S3 method for class 'lme11' summary(object, ...) ## S3 method for class 'lme11' plot(x, ...) ## S3 method for class 'magev' print(x, ...) ## S3 method for class 'magev' summary(object, ...) ## S3 method for class 'magev' plot(x, ...)
x, object
|
A fitted object. |
... |
Further arguments passed to plotting functions; ignored by
|
x (or object), invisibly.
data(haenam) fit <- glme.gev(haenam$X1) fit summary(fit) plot(fit)data(haenam) fit <- glme.gev(haenam$X1) fit summary(fit) plot(fit)
Estimates the parameters of the non-stationary GEV11 model by the quantile-based GN16 method (Gilleland and Katz, 2016; as adapted in Shin et al., 2025b), with both the original and the modified specification returned.
GN16.gev11(xdat, rob = FALSE)GN16.gev11(xdat, rob = FALSE)
xdat |
A numeric vector of data to be fitted. |
rob |
If TRUE, uses robust regression in the internal steps (default FALSE). |
A list containing:
para.gado.org - Original GN16 estimates (mu0, mu1, sigma0, sigma1, xi).
para.gado.mdfy - Modified GN16 estimates.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Gilleland, E. & Katz, R. W. (2016). extRemes 2.0: An extreme value analysis package in R. Journal of Statistical Software, 72(8), 1-39.
Shin, Y., Shin, Y. & Park, J.-S. (2025). Building nonstationary extreme value model using L-moments. Journal of the Korean Statistical Society, 54, 947-970. doi:10.1007/s42952-025-00325-3
glme.gev11, lme.gev11,
strup.gev11.
data(PhliuAgromet) result <- GN16.gev11(PhliuAgromet$prec) print(result$para.gado.org)data(PhliuAgromet) result <- GN16.gev11(PhliuAgromet$prec) print(result$para.gado.org)
Annual maximum daily rainfall data from Haenam, South Korea. This dataset is used for demonstrating model averaging methods for high quantile estimation in extreme value analysis.
haenamhaenam
A data frame with 52 rows and 2 columns:
Year of observation (integer, 1971-2022)
Annual maximum daily rainfall in mm (numeric)
Korea Meteorological Administration (KMA; https://www.weather.go.kr)
Shin, Y., Shin, Y., & Park, J. S. (2026). Model averaging with mixed criteria for estimating high quantiles of extreme values: Application to heavy rainfall. Stochastic Environmental Research and Risk Assessment, 40(2), 47. doi:10.1007/s00477-025-03167-x
data(haenam) head(haenam) # Estimate high quantiles using model averaging result <- ma.gev(haenam$X1, quant = c(0.98, 0.99, 0.995)) print(result$qua.ma)data(haenam) head(haenam) # Estimate high quantiles using model averaging result <- ma.gev(haenam$X1, quant = c(0.98, 0.99, 0.995)) print(result$qua.ma)
Generates multiple starting parameter sets for multi-start optimization
in GLME estimation, for both the stationary GEV (model="gev00")
and the non-stationary GEV11 (model="gev11"). Uses L-moment
estimates as a base and adds random perturbations.
For model="gev11" the second column is returned on the log scale
(matching sigma0 of the GEV11 parameterization), and when pretheta
is supplied the first candidate is a perturbation of the pre-estimate.
init.glme( data = NULL, ntry = NULL, model = "gev00", pretheta = NULL, xdat = NULL )init.glme( data = NULL, ntry = NULL, model = "gev00", pretheta = NULL, xdat = NULL )
data |
A numeric vector of data to be fitted. |
ntry |
Number of initial parameter sets to generate (minimum 2). |
model |
Either "gev00" (stationary, default) or "gev11". |
pretheta |
Optional pre-estimate (mu0, mu1, sigma0, sigma1, xi)
used for the first candidate when |
xdat |
Deprecated alias of |
A matrix with ntry rows and 3 columns,
where each row is a candidate starting point for optimization.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
glme.gev, glme.gev11 which use this
function internally.
data(haenam) inits <- init.glme(haenam$X1, ntry = 5) print(inits)data(haenam) inits <- init.glme(haenam$X1, ntry = 5) print(inits)
Estimates the parameters of the non-stationary GEV11 model
(mu(t) = mu0 + mu1*t, sigma(t) = exp(sigma0 + sigma1*t), constant xi)
by the pure L-moment method of Shin et al. (2025b, JKSS).
Starting from the WLS pre-estimate of strup.gev11,
the three L-moment equations (matching the sample L-moments of the
Gumbel-transformed data with the population L-moments of the standard
Gumbel distribution) are solved by Broyden's method for
(mu0, sigma0, xi), and the best solution is chosen by a goodness-of-fit
measure.
This function replaces gado.prop_11 of GLmom v1.x
(renamed in v2.0.0, with slightly changed input/output).
It is also equivalent to glme.gev11(xdat, pen="no").
lme.gev11(xdat, ntry = 5, ftol = 1e-05, show = FALSE)lme.gev11(xdat, ntry = 5, ftol = 1e-05, show = FALSE)
xdat |
A numeric vector of data to be fitted. |
ntry |
Number of attempts (initial points) for the nonlinear solver (default 5). |
ftol |
Function tolerance for the solver (default 1e-5). |
show |
If TRUE, prints the precision and parameters of each try. |
A list containing:
lme.gev11 - The L-moment estimates (mu0, mu1, sigma0, sigma1, xi). If no solution is found, the WLS pre-estimate is returned.
precis - Precision (mean absolute residual of the L-moment equations) of the selected solution.
data - The input data (used by the plot method).
The object has class "lme11"; see GLmom-methods.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y. & Park, J.-S. (2025). Building nonstationary extreme value model using L-moments. Journal of the Korean Statistical Society, 54, 947-970. doi:10.1007/s42952-025-00325-3
glme.gev11 for the penalized (GLME) version,
strup.gev11 for the WLS pre-estimation,
GN16.gev11 for the GN16 method,
nsgev and gado.prop_11 for the v1.x
compatibility wrappers.
data(PhliuAgromet) result <- lme.gev11(PhliuAgromet$prec, ntry = 5) print(result$lme.gev11)data(PhliuAgromet) result <- lme.gev11(PhliuAgromet$prec, ntry = 5) print(result$lme.gev11)
This function estimates high quantiles of the Generalized Extreme Value (GEV) distribution using model averaging with mixed criteria. It combines Maximum Likelihood Estimation (MLE) and L-moment Estimation (LME) to construct candidate submodels and assign weights effectively.
ma.gev( data = NULL, quant = c(0.98, 0.99, 0.995), weight = "like1", numk = 12, B = 200, varcom = TRUE, trim = 0, fig = FALSE, bma = FALSE, pen = "norm", CD = FALSE, remle = FALSE )ma.gev( data = NULL, quant = c(0.98, 0.99, 0.995), weight = "like1", numk = 12, B = 200, varcom = TRUE, trim = 0, fig = FALSE, bma = FALSE, pen = "norm", CD = FALSE, remle = FALSE )
data |
A numeric vector of data to be fitted (e.g., annual maxima). |
quant |
The probabilities corresponding to high quantiles to be estimated. Default is c(0.98, 0.99, 0.995). |
weight |
The weighting method name. Options are:
Variants with numbers indicate left trimming level (0, 1, or 2). |
numk |
The number of candidate submodels K. Default is 12. |
B |
The number of bootstrap samples. Default is 200. |
varcom |
Logical. Whether to compute variance of quantile estimates. Default is TRUE. |
trim |
The number of left trimming for L-moments. Usually 0 (default), 1, or 2. |
fig |
Logical. Whether to produce diagnostic plots. Default is FALSE. |
bma |
Logical. Whether to use Bayesian Model Averaging. Default is FALSE. |
pen |
Penalty type for BMA prior: 'norm' (normal, default) or 'beta'. |
CD |
Logical. Whether to compute Coles-Dixon penalized MLE. Default is FALSE. |
remle |
Logical. Whether to compute restricted MLE. Default is FALSE. |
The model averaging approach works as follows:
MLE and LME of GEV parameters are computed
K candidate shape parameters (xi) are selected from profile likelihood CI
For each candidate xi, MLE with fixed xi is computed
Weights are assigned based on the selected method
Final quantile estimates are weighted averages across submodels
The weighting schemes include:
'like': AIC-based weights using likelihood with fixed xi
'gLd': Weights based on generalized L-moment distance
'med': Weights based on median and L-moment distance
'cvt': Conventional AIC weights
When bma=TRUE, Bayesian model averaging is applied with prior specified by pen.
A list containing:
mle.hosking - MLE estimates in Hosking style (mu, sigma, xi)
qua.mle - Quantile estimates from MLE
mle.cov3 - Covariance matrix of MLE (3x3)
qua.se.mle.delta - Standard errors of MLE quantiles (delta method)
lme - L-moment estimates (mu, sigma, xi)
lme.cov3 - Covariance matrix of LME (bootstrap)
qua.lme - Quantile estimates from LME
qua.se.lme.boots - Standard errors of LME quantiles (bootstrap)
qua.ma - Model-averaged quantile estimates
w.ma - Weights used for model averaging
fixw.se.ma - Asymptotic SE under fixed weights
ranw.se.ma - Asymptotic SE under random weights
surr - Surrogate model parameters (mu, sigma, xi)
pick_xi - Selected xi values for K submodels
qua.bma - (if bma=TRUE) BMA quantile estimates
w.bma - (if bma=TRUE) BMA weights
mle.CD - (if CD=TRUE) Coles-Dixon penalized MLE
qua.CD - (if CD=TRUE) Quantile estimates from CD-MLE
remle1 - (if remle=TRUE) Restricted MLE (first constraint)
qua.remle1 - (if remle=TRUE) Quantile estimates from remle1
remle2 - (if remle=TRUE) Restricted MLE (second constraint)
qua.remle2 - (if remle=TRUE) Quantile estimates from remle2
quant - The quantile probabilities used
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y., & Park, J. S. (2026). Model averaging with mixed criteria for estimating high quantiles of extreme values: Application to heavy rainfall. Stochastic Environmental Research and Risk Assessment, 40(2), 47. doi:10.1007/s00477-025-03167-x
glme.gev for stationary GLME estimation,
magev.ksensplot for K sensitivity analysis,
magev.qqplot for Q-Q diagnostic plots,
magev.rlplot for return level plots.
# Load example data data(haenam) x <- haenam$X1 # Basic usage with likelihood weights result <- ma.gev(x, quant = c(0.95, 0.99), weight = 'like1', B = 100) print(result$qua.ma) # Model-averaged quantiles print(result$qua.mle) # MLE quantiles for comparison print(result$qua.lme) # LME quantiles for comparison # Using generalized L-moment distance weights result2 <- ma.gev(x, quant = c(0.95, 0.99), weight = 'gLd', B = 100) print(result2$w.ma) # Model weights# Load example data data(haenam) x <- haenam$X1 # Basic usage with likelihood weights result <- ma.gev(x, quant = c(0.95, 0.99), weight = 'like1', B = 100) print(result$qua.ma) # Model-averaged quantiles print(result$qua.mle) # MLE quantiles for comparison print(result$qua.lme) # LME quantiles for comparison # Using generalized L-moment distance weights result2 <- ma.gev(x, quant = c(0.95, 0.99), weight = 'gLd', B = 100) print(result2$w.ma) # Model weights
Plots return level estimates, standard errors, and first-order differences across different numbers of candidate submodels K. This helps identify stable regions where estimates converge and select an optimal K value.
magev.ksensplot( data = NULL, q.cut = 0.6, mink = 4, maxk = 20, quant = c(0.99, 0.995) )magev.ksensplot( data = NULL, q.cut = 0.6, mink = 4, maxk = 20, quant = c(0.99, 0.995) )
data |
A numeric vector of data to be fitted (e.g., annual maxima). |
q.cut |
Quantile cutoff for determining stability (default 0.6). |
mink |
Minimum number of candidate submodels to test (default 4). |
maxk |
Maximum number of candidate submodels to test (default 20). |
quant |
The probabilities for high quantile estimation. Default is c(0.99, 0.995). |
The function computes MAGEV estimates for K ranging from mink to
maxk. For each K, it calculates:
Return level estimates (black points)
Normalized standard errors (blue line)
First-order differences (red line with triangles)
The optimal K is selected as the smallest value where both the normalized
standard error and first-order difference are below their respective
q.cut quantile cutoffs. The selected K is indicated by a purple
vertical line.
The optimal K value (integer) selected by the algorithm.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y., & Park, J. S. (2026). Model averaging with mixed criteria for estimating high quantiles of extreme values: Application to heavy rainfall. Stochastic Environmental Research and Risk Assessment, 40(2), 47. doi:10.1007/s00477-025-03167-x
ma.gev for the main model averaging function.
data(haenam) optimal_k <- magev.ksensplot(haenam$X1) print(optimal_k)data(haenam) optimal_k <- magev.ksensplot(haenam$X1) print(optimal_k)
Creates a 2x2 panel of Q-Q plots comparing observed vs. fitted quantiles for different estimation methods: MLE, LME, surrogate (MA), and REMLE.
magev.qqplot(data = NULL, zx = NULL)magev.qqplot(data = NULL, zx = NULL)
data |
A numeric vector of observed data. |
zx |
A list object returned by |
The function creates four Q-Q plots:
Upper left: MLE (Maximum Likelihood Estimation)
Upper right: LME (L-moment Estimation)
Lower left: Surrogate MA (Model Averaging surrogate)
Lower right: REMLE (Restricted MLE, if available)
Points close to the 45-degree diagonal line indicate good model fit.
NULL. The function produces a plot as a side effect.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y., & Park, J. S. (2026). Model averaging with mixed criteria for estimating high quantiles of extreme values: Application to heavy rainfall. Stochastic Environmental Research and Risk Assessment, 40(2), 47. doi:10.1007/s00477-025-03167-x
ma.gev for the main model averaging function,
magev.rlplot for return level plots.
data(haenam) qq <- c(seq(0.01, 0.99, by = 0.01), 0.995, 0.999) zx <- ma.gev(haenam$X1, quant = qq, weight = 'like1', numk = 9, varcom = FALSE, remle = TRUE) magev.qqplot(data = haenam$X1, zx = zx)data(haenam) qq <- c(seq(0.01, 0.99, by = 0.01), 0.995, 0.999) zx <- ma.gev(haenam$X1, quant = qq, weight = 'like1', numk = 9, varcom = FALSE, remle = TRUE) magev.qqplot(data = haenam$X1, zx = zx)
Displays fitted return levels with 95 return period on a log scale.
magev.rlplot(par = NULL, se.vec = NULL, data = NULL)magev.rlplot(par = NULL, se.vec = NULL, data = NULL)
par |
A numeric vector of GEV parameters (mu, sigma, xi) in Hosking style.
Typically from |
se.vec |
A numeric vector of standard errors for the quantile estimates
corresponding to the plotting positions. Typically from |
data |
A numeric vector of observed data (annual maxima). |
The plot shows:
Black line: Fitted return level curve
Blue lines: 95
Black points: Observed data at empirical return periods
The x-axis (return period) is on a log scale, ranging from 0.1 to 900 years.
NULL. The function produces a plot as a side effect.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y., & Park, J. S. (2026). Model averaging with mixed criteria for estimating high quantiles of extreme values: Application to heavy rainfall. Stochastic Environmental Research and Risk Assessment, 40(2), 47. doi:10.1007/s00477-025-03167-x
ma.gev for the main model averaging function,
magev.qqplot for Q-Q diagnostic plots.
data(haenam) ff <- c(seq(0.01, 0.09, by = 0.01), 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 0.93, 0.95, 0.98, 0.99, 0.993, 0.995, 0.998, 0.999) zx <- ma.gev(haenam$X1, quant = ff, weight = 'like1', numk = 9, varcom = TRUE) magev.rlplot(par = zx$surr$par, se.vec = zx$ranw.se.ma, data = haenam$X1)data(haenam) ff <- c(seq(0.01, 0.09, by = 0.01), 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 0.93, 0.95, 0.98, 0.99, 0.993, 0.995, 0.998, 0.999) zx <- ma.gev(haenam$X1, quant = ff, weight = 'like1', numk = 9, varcom = TRUE) magev.rlplot(par = zx$surr$par, se.vec = zx$ranw.se.ma, data = haenam$X1)
Computes the Martins-Stedinger Beta(6,9) prior probability for the GEV
shape parameter on the interval .
MS_pk(para = para, p = 6, q = 9)MS_pk(para = para, p = 6, q = 9)
para |
A vector of GEV parameters (location, scale, shape). |
p |
Shape parameter for beta distribution (default 6). |
q |
Shape parameter for beta distribution (default 9). |
Prior probability value (scalar).
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Martins, E. S. & Stedinger, J. R. (2000). Generalized maximum-likelihood generalized extreme-value quantile estimators for hydrologic data. Water Resources Research, 36(3), 737-744. doi:10.1029/1999WR900330
pk.beta for the adaptive beta penalty,
glme.gev which uses these penalty functions.
# Evaluate MS prior at xi = -0.2 MS_pk(para = c(100, 20, -0.2))# Evaluate MS prior at xi = -0.2 MS_pk(para = c(100, 20, -0.2))
Estimates parameters of a non-stationary Generalized Extreme Value (GEV) distribution using the L-moment-based algorithm from Shin et al. (2025, J. Korean Stat. Soc.).
This is a convenience wrapper around lme.gev11, providing
compatibility with the original nsgev package interface (and with
GLmom v1.x).
nsgev(xdat, ntry = 20, ftol = 1e-06)nsgev(xdat, ntry = 20, ftol = 1e-06)
xdat |
A numeric vector of data to be fitted (e.g., annual maximum values). |
ntry |
Number of attempts for optimization (default 20). |
ftol |
Function tolerance for optimization (default 1e-6). |
A list containing:
para.prop - The proposed L-moment based estimates
(mu0, mu1, sigma0, sigma1, xi).
precis - Precision of the optimization.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y. & Park, J.-S. (2025). Building nonstationary extreme value model using L-moments. Journal of the Korean Statistical Society, 54, 947-970. doi:10.1007/s42952-025-00325-3
lme.gev11 for the recommended interface,
glme.gev11 for the full GLME method with penalty functions,
gado.prop_11 for detailed estimation results.
data(PhliuAgromet) result <- nsgev(PhliuAgromet$prec, ntry = 5) print(result$para.prop)data(PhliuAgromet) result <- nsgev(PhliuAgromet$prec, ntry = 5) print(result$para.prop)
Estimates GEV location and scale parameters from L-moments while keeping
the shape parameter fixed at a user-specified value. Modified from
lmomco::pargev().
pargev.kfix(lmom, kfix = 0.1, checklmom = TRUE, ...)pargev.kfix(lmom, kfix = 0.1, checklmom = TRUE, ...)
lmom |
L-moments object. |
kfix |
Fixed shape parameter value. |
checklmom |
Whether to check L-moment validity. |
... |
Additional arguments. |
A list with components:
Character "gev"
Numeric vector of GEV parameters (xi=location, alpha=scale, kappa=shape)
Character "pargev"
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Hosking, J. R. M. (1990). L-moments: Analysis and estimation of distributions using linear combinations of order statistics. Journal of the Royal Statistical Society, Series B, 52(1), 105-124. doi:10.1111/j.2517-6161.1990.tb01775.x
glme.gev for GLME estimation,
pargev for the original L-moment GEV fitting.
data(haenam) lmom <- lmomco::lmoms(haenam$X1, nmom = 3) pargev.kfix(lmom, kfix = -0.1)data(haenam) lmom <- lmomco::lmoms(haenam$X1, nmom = 3) pargev.kfix(lmom, kfix = -0.1)
Annual maximum daily precipitation and station metadata from the Phliu Agrometeorological Station, Chanthaburi province, Thailand (1984-2023). A Mann-Kendall test indicates a significant increasing trend in the annual maxima (tau = 0.235, p = 0.033), so this dataset is used for demonstrating non-stationary (GEV11) estimation; it is the heavy rainfall application of Shin et al. (2025).
PhliuAgrometPhliuAgromet
A data frame with 40 rows and 9 columns:
Station identifier (character)
Year of observation (numeric, 1984-2023)
Annual maximum daily precipitation in mm (numeric)
Station name (character)
Climate zone code (character)
Station latitude in degrees (numeric)
Station longitude in degrees (numeric)
Record start year (integer)
Record end year (numeric)
Phliu Agrometeorological Station, Thailand.
Shin, Y., Shin, Y., Park, J., & Park, J. S. (2025). Generalized method of L-moment estimation for stationary and nonstationary extreme value models. arXiv preprint arXiv:2512.20385. doi:10.48550/arXiv.2512.20385
data(PhliuAgromet) head(PhliuAgromet) # Fit non-stationary GEV11 model by GLME result <- glme.gev11(PhliuAgromet$prec, ntry = 5) print(result$para.glme)data(PhliuAgromet) head(PhliuAgromet) # Fit non-stationary GEV11 model by GLME result <- glme.gev11(PhliuAgromet$prec, ntry = 5) print(result$para.glme)
Computes the data-adaptive Beta preference (penalty) function for the GEV
shape parameter, used by both the stationary (glme.gev) and
the non-stationary (glme.gev11) GLME methods (v2.0.0 unified
implementation).
When q is NULL (adaptive mode), the support is
where
is lme.center[3], and
for
(otherwise ). When q is given (fixed mode), the support
is (min.xi, max.xi); this mode is used for the fixed literature
penalties ("ms", "park", "cannon").
pk.beta.stnary is kept as a backward-compatible alias. Note that
since v2.0.0 both names use this unified implementation (the v1.x default
hyperparameters c0=0.3, c2=5 changed to c0=0.35, c2=2, and
out-of-support values now return 1e-100 instead of 1).
pk.beta( para = NULL, lme.center = NULL, p = NULL, q = NULL, min.xi = -1, max.xi = 0.5, c0 = 0.35, c1 = 10, c2 = 2 ) pk.beta.stnary( para = NULL, lme.center = NULL, p = NULL, q = NULL, min.xi = -1, max.xi = 0.5, c0 = 0.35, c1 = 10, c2 = 2 )pk.beta( para = NULL, lme.center = NULL, p = NULL, q = NULL, min.xi = -1, max.xi = 0.5, c0 = 0.35, c1 = 10, c2 = 2 ) pk.beta.stnary( para = NULL, lme.center = NULL, p = NULL, q = NULL, min.xi = -1, max.xi = 0.5, c0 = 0.35, c1 = 10, c2 = 2 )
para |
A vector of GEV parameters; only |
lme.center |
L-moment estimates used as the center (adaptive mode). |
p |
First shape parameter of the beta distribution. |
q |
Second shape parameter (optional; if provided, fixed mode is used). |
min.xi |
Lower support limit in fixed mode (default -1). |
max.xi |
Upper support limit in fixed mode (default 0.5). |
c0 |
Half-width of the adaptive support (default 0.35). |
c1 |
Scaling hyperparameter for q (default 10). |
c2 |
Upper limit hyperparameter for q (default 2). |
A list containing:
Preference function value (scalar)
Shape parameter p used
Shape parameter q used
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y., Park, J. & Park, J.-S. (2025). Generalized method of L-moment estimation for stationary and nonstationary extreme value models. arXiv preprint arXiv:2512.20385. doi:10.48550/arXiv.2512.20385
pk.norm.stnary for the normal penalty,
MS_pk for the Martins-Stedinger penalty,
glme.gev, glme.gev11.
# Adaptive beta preference for xi = -0.2 centered at LME xi = -0.15 pk.beta(para = c(100, 20, -0.2), lme.center = c(100, 20, -0.15), p = 6, c1 = 3, c2 = 1)# Adaptive beta preference for xi = -0.2 centered at LME xi = -0.15 pk.beta(para = c(100, 20, -0.2), lme.center = c(100, 20, -0.15), p = 6, c1 = 3, c2 = 1)
Computes a normal distribution-based preference (penalty) function value
for the GEV shape parameter. The alias new_pf_norm is provided
for backward compatibility.
pk.norm.stnary(para = NULL, mu = NULL, std = NULL) new_pf_norm(para = NULL, mu = NULL, std = NULL)pk.norm.stnary(para = NULL, mu = NULL, std = NULL) new_pf_norm(para = NULL, mu = NULL, std = NULL)
para |
A vector of GEV parameters. |
mu |
Mean for normal distribution. |
std |
Standard deviation for normal distribution. |
Preference function value (scalar).
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y., Park, J. & Park, J.-S. (2025). Generalized method of L-moment estimation for stationary and nonstationary extreme value models. arXiv preprint arXiv:2512.20385. doi:10.48550/arXiv.2512.20385
pk.beta for the beta penalty,
MS_pk for the Martins-Stedinger penalty,
glme.gev which uses these penalty functions.
# Normal preference with mean=-0.5, sd=0.2 at xi=-0.2 pk.norm.stnary(para = c(100, 20, -0.2), mu = -0.5, std = 0.2)# Normal preference with mean=-0.5, sd=0.2 at xi=-0.2 pk.norm.stnary(para = c(100, 20, -0.2), mu = -0.5, std = 0.2)
Calculates quantiles for non-stationary GEV models including GEV11, GEV10, GEV20, and stationary GEV00.
quagev.NS(f = NULL, para = NULL, nsample = NULL, model = NULL)quagev.NS(f = NULL, para = NULL, nsample = NULL, model = NULL)
f |
Probability (or vector of probabilities) for quantile calculation. |
para |
Parameter vector. For GEV11: (mu0, mu1, sigma0, sigma1, xi). |
nsample |
Number of time points (sample size). |
model |
Model type: "gev11", "gev10", "gev20", "gev01", or "gev00"/"gev". |
A matrix of quantiles (nsample x length(f)) or a vector if f is scalar.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y., Park, J. & Park, J.-S. (2025). Generalized method of L-moment estimation for stationary and nonstationary extreme value models. arXiv preprint arXiv:2512.20385. doi:10.48550/arXiv.2512.20385
glme.gev11 for non-stationary GEV estimation,
glme.gev for stationary GEV estimation,
ran.gev_all for random number generation.
# GEV11 model: time-varying quantiles para <- c(84.55, 1.03, 2.91, 0.009, -0.08) # mu0, mu1, sigma0, sigma1, xi q99 <- quagev.NS(f = 0.99, para = para, nsample = 53, model = "gev11") print(q99)# GEV11 model: time-varying quantiles para <- c(84.55, 1.03, 2.91, 0.009, -0.08) # mu0, mu1, sigma0, sigma1, xi q99 <- quagev.NS(f = 0.99, para = para, nsample = 53, model = "gev11") print(q99)
Generates a random sample from a stationary or non-stationary GEV model. For non-stationary models the location and/or scale vary with time t = 1, ..., nsample: mu(t) = mu0 + mu1*t + mu2*t^2 and sigma(t) = exp(sigma0 + sigma1*t).
ran.gev_all(nsample, para, model = NULL)ran.gev_all(nsample, para, model = NULL)
nsample |
Sample size (number of time points). |
para |
Parameter vector; its length and meaning depend on
|
model |
Model type: "gev", "gev00", "gev01", "gev10", "gev11", or "gev20". |
A numeric vector of length nsample.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Shin, Y., Shin, Y. & Park, J.-S. (2025). Building nonstationary extreme value model using L-moments. Journal of the Korean Statistical Society, 54, 947-970. doi:10.1007/s42952-025-00325-3
quagev.NS for the corresponding quantile function,
glme.gev11, lme.gev11.
set.seed(1) # Non-stationary GEV11 sample of size 50 x <- ran.gev_all(50, para = c(100, 0.5, 3, 0.005, -0.2), model = "gev11") plot(x, type = "l")set.seed(1) # Non-stationary GEV11 sample of size 50 x <- ran.gev_all(50, para = c(100, 0.5, 3, 0.005, -0.2), model = "gev11") plot(x, type = "l")
Estimates the parameters of the non-stationary GEV11 model by the
weighted least squares (WLS) method of Strupczewski and Kaczmarek (2001),
with the final-specification modification described in Shin et al.
(2025b, JKSS). The result is used as the pre-estimate (starting value)
by glme.gev11 and lme.gev11.
strup.gev11(xdat, init.rob = TRUE, wls.rob = FALSE, strup.fin.rob = TRUE)strup.gev11(xdat, init.rob = TRUE, wls.rob = FALSE, strup.fin.rob = TRUE)
xdat |
A numeric vector of data to be fitted. |
init.rob |
If TRUE (default), the initial location trend is fitted
by robust regression ( |
wls.rob |
If TRUE, the log-scale regression uses robust regression (default FALSE). |
strup.fin.rob |
If TRUE (default), the final location trend re-specification uses robust regression. |
A list containing:
strup.para - Original WLS estimates (mu0, mu1, sigma0, sigma1, xi).
strup.mdfy - Modified (re-specified) WLS estimates.
Yonggwan Shin, Seokkap Ko, Jihong Park, Yire Shin, Jeong-Soo Park
Strupczewski, W. G. & Kaczmarek, Z. (2001). Non-stationary approach to at-site flood frequency modelling II. Weighted least squares estimation. Journal of Hydrology, 248(1-4), 143-151.
Shin, Y., Shin, Y. & Park, J.-S. (2025). Building nonstationary extreme value model using L-moments. Journal of the Korean Statistical Society, 54, 947-970. doi:10.1007/s42952-025-00325-3
glme.gev11, lme.gev11,
GN16.gev11.
data(PhliuAgromet) result <- strup.gev11(PhliuAgromet$prec) print(result$strup.mdfy)data(PhliuAgromet) result <- strup.gev11(PhliuAgromet$prec) print(result$strup.mdfy)