Non-Stationary GEV.
NSGEV(
formulas,
data,
response = NULL,
psi = NULL,
est = c("none", "optim", "nloptr"),
trace = 0
)
A named list with three formulas.
The data frame to be used.
Name of the column in data
that gives the
response. If NULL
, the model will not be estimated but still
can have a limited use.
Named vector of coefficients. If the model is (re)-estimated, this vector is used a initial value for the coefficients.
Character or logical. The default is to use the parameter
values given in psi
.
Integer level of verbosity.
An object with class "NSGEV"
The model involves a vector \(\boldsymbol{\psi}\) of \(p\) parameters which are called the coefficients of the model. The GEV parameters \(\mu\) (loc), \(\sigma\) (scale) and \(\xi\) (shape) are expressed as formulas using the parameters in \(\boldsymbol{\psi}\) and the model variables. The three GEV parameters form a vector \(\boldsymbol{\theta}(\mathbf{x}) = [\mu(\mathbf{x}),\,\sigma(\mathbf{x}),\,xi(\mathbf{x})]^\top\) which depends on the covariates hence varies across the observations. The observations are assumed to be independent conditional on the covariates.
Note that the GEV parameters are always assumed to be given in the
order loc
, scale
, shape
. The data frame data
can not for now use these names for its columns.
df <- data.frame(t = 1:10)
## built a model with given coefficients
psi <- c("alpha" = 1, "beta" = 0.01, "delta" = 0.6, "xi" = 0.06)
ns0 <- NSGEV(formulas = list("loc" = ~ alpha + beta * t, "scale" = ~ delta, "shape" = ~ xi),
data = df, psi = psi)
## simulate a path
set.seed(1234)
ysim <- simulate(ns0, nsim = 1, psi = psi)
df2 <- cbind(df, y = ysim[ , 1L])
ns1 <- NSGEV(formulas = list("loc" = ~ alpha + beta * t, "scale" = ~ delta, "shape" = ~ xi),
data = df2, response = "y", psi = psi, est = "optim")
## try an exponential link
ns2 <- NSGEV(formulas = list("loc" = ~ exp(alpha + beta * t), "scale" = ~ delta, "shape" = ~ xi),
data = df2, response = "y", psi = psi, est = "optim")
## compare the estimation with that of ismev::gev.fit
require(ismev)
ns1.ismev <- gev.fit(xdat = df2$y, ydat = as.matrix(df), mul = 1, show = FALSE)
rbind("NSGEV" = c(ns1$estimate, "negLogL" = ns1$negLogL),
"ismev" = c(ns1.ismev$mle, "negLogL" = ns1.ismev$nllh))
#> alpha beta delta xi negLogL
#> NSGEV 1.143777 -0.002918431 0.6033124 -0.4359059 8.22994
#> ismev 1.143704 -0.002906355 0.6032869 -0.4358915 8.22994
## Try an expoential link
ns2.ismev <- gev.fit(xdat = df2$y, ydat = as.matrix(df), mul = 1, mulink = exp, show = FALSE)
rbind("NSGEV" = c(ns2$estimate, "negLogL" = ns2$negLogL),
"ismev" = c(ns2.ismev$mle, "negLogL" = ns2.ismev$nllh))
#> alpha beta delta xi negLogL
#> NSGEV 0.3976449 -0.11871098 1.5582157 -1.0861830 8.106294
#> ismev 0.1284920 -0.00155176 0.6035595 -0.4351343 8.230150