fgev(x, start, ..., nsloc = NULL, prob = NULL, std.err = TRUE,
corr = FALSE, method = "BFGS", warn.inf = TRUE)start is omitted the routine attempts to find good
starting values using moment estimators.optim. If parameters
of the model are included they will be held fixed at the
values given (see Examples).x, for linear modelling of the location
parameter.
The data frame is treated as a covariate matrix (excluding the
intercept).
A numeric vector can be given as an aNULL (the default),
or a probability in the closed interval [0,1].TRUE (the default), the standard
errors are returned.TRUE, the correlation matrix is
returned.optim for
details).TRUE (the default), a warning is
given if the negative log-likelihood is infinite when evaluated at
the starting values.c("gev","evd"). The generic accessor functions fitted (or
fitted.values), std.errors,
deviance, logLik and
AIC extract various features of the
returned object.
The functions profile and profile2d are
used to obtain deviance profiles for the model parameters.
In particular, profiles of the quantile $z_p$ can be
calculated and plotted when $\code{prob} = p$.
The function anova compares nested models.
The function plot produces diagnostic plots.
An object of class c("gev","evd") is a list containing
at most the following components
optim.x.nsloc.x.prob.prob is NULL
(the default), this will also be an element of param.prob is NULL (the default):
For stationary models the parameter names are loc, scale
and shape, for the location, scale and shape parameters
respectively.
For non-stationary models, the parameter names are loc,
locx1, ..., locxn, scale and
shape, where x1, ..., xn are the column names
of nsloc, so that loc is the intercept of the
linear model, and locx1, ..., locxn
are the ncol(nsloc) coefficients.
If nsloc is a vector it is converted into a single column
data frame with column name trend, and hence the associated
trend parameter is named loctrend. If $\code{prob} = p$ is a probability:
The fit is performed using a different parameterization.
Let $a$, $b$ and $s$ denote the location, scale
and shape parameters of the GEV distribution.
For stationary models, the distribution is parameterized
using $(z_p,b,s)$, where
$$z_p = a - b/s (1 - (-\log(1 - p))^s)$$
is such that $G(z_p) = 1 - p$, where $G$ is the
GEV distribution function.
$\code{prob} = p$ is therefore the probability in the upper
tail corresponding to the quantile $z_p$.
If prob is zero, then $z_p$ is the upper end point
$a - b/s$, and $s$ is restricted to the negative
(Weibull) axis.
If prob is one, then $z_p$ is the lower end point
$a - b/s$, and $s$ is restricted to the positive
(Frechet) axis.
The parameter names are quantile, scale
and shape, for $z_p$, $b$ and $s$
respectively.
For non-stationary models the parameter $z_p$ is again given by
the equation above, but $a$ becomes the intercept of the linear
model for the location parameter, so that quantile replaces
(the intercept) loc, and hence the parameter names are
quantile, locx1, ..., locxn,
scale and shape, where x1, ..., xn are
the column names of nsloc.
In either case:
For non-stationary fitting it is recommended that the covariates
within the linear model for the location parameter are (at least
approximately) centered and scaled (i.e. that the columns of
nsloc are centered and scaled), particularly if automatic
starting values are used, since the starting values for the
associated parameters are then zero.
anova.evd, optim,
plot.gev, profile.evd,
profile2d.evduvdata <- rgev(100, loc = 0.13, scale = 1.1, shape = 0.2)
trend <- (-49:50)/100
M1 <- fgev(uvdata, nsloc = trend, control = list(trace = 1))
M2 <- fgev(uvdata)
M3 <- fgev(uvdata, shape = 0)
M4 <- fgev(uvdata, scale = 1, shape = 0)
anova(M1, M2, M3, M4)
plot(M2)
M2P <- profile(M2)
plot(M2P)
rnd <- runif(100, min = -.5, max = .5)
fgev(uvdata, nsloc = data.frame(trend = trend, random = rnd))
fgev(uvdata, nsloc = data.frame(trend = trend, random = rnd), locrandom = 0)
uvdata <- rgev(100, loc = 0.13, scale = 1.1, shape = 0.2)
M1 <- fgev(uvdata, prob = 0.1)
M2 <- fgev(uvdata, prob = 0.01)
M1P <- profile(M1, which = "quantile")
M2P <- profile(M2, which = "quantile")
plot(M1P)
plot(M2P)Run the code above in your browser using DataLab