Description Usage Arguments Details Value Acknowledgments Note Author(s) References See Also Examples
View source: R/flognormgpdcon.r
Maximum likelihood estimation for fitting the extreme value mixture model with log-normal for bulk distribution upto the threshold and conditional GPD above threshold with continuity at threshold. With options for profile likelihood estimation for threshold and fixed threshold approach.
1 2 3 4 5 6 7 8 9 10 11 12 13 | flognormgpdcon(x, phiu = TRUE, useq = NULL, fixedu = FALSE,
pvector = NULL, std.err = TRUE, method = "BFGS",
control = list(maxit = 10000), finitelik = TRUE, ...)
llognormgpdcon(x, lnmean = 0, lnsd = 1, u = qlnorm(0.9, lnmean,
lnsd), xi = 0, phiu = TRUE, log = TRUE)
nllognormgpdcon(pvector, x, phiu = TRUE, finitelik = FALSE)
proflulognormgpdcon(u, pvector, x, phiu = TRUE, method = "BFGS",
control = list(maxit = 10000), finitelik = TRUE, ...)
nlulognormgpdcon(pvector, u, x, phiu = TRUE, finitelik = FALSE)
|
x |
vector of sample data |
phiu |
probability of being above threshold (0, 1) or logical, see Details in
help for |
useq |
vector of thresholds (or scalar) to be considered in profile likelihood or
|
fixedu |
logical, should threshold be fixed (at either scalar value in |
pvector |
vector of initial values of parameters or |
std.err |
logical, should standard errors be calculated |
method |
optimisation method (see |
control |
optimisation control list (see |
finitelik |
logical, should log-likelihood return finite value for invalid parameters |
... |
optional inputs passed to |
lnmean |
scalar mean on log scale |
lnsd |
scalar standard deviation on log scale (positive) |
u |
scalar threshold value |
xi |
scalar shape parameter |
log |
logical, if |
The extreme value mixture model with log-normal bulk and GPD tail with continuity at threshold is fitted to the entire dataset using maximum likelihood estimation. The estimated parameters, variance-covariance matrix and their standard errors are automatically output.
See help for fnormgpd
for details, type help fnormgpd
.
Only the different features are outlined below for brevity.
The GPD sigmau
parameter is now specified as function of other parameters, see
help for dlognormgpdcon
for details, type help lognormgpdcon
.
Therefore, sigmau
should not be included in the parameter vector if initial values
are provided, making the full parameter vector
(lnmean
, lnsd
, u
, xi
) if threshold is also estimated and
(lnmean
, lnsd
, xi
) for profile likelihood or fixed threshold approach.
Non-positive data are ignored.
Log-likelihood is given by llognormgpdcon
and it's
wrappers for negative log-likelihood from nllognormgpdcon
and nlulognormgpdcon
. Profile likelihood for single
threshold given by proflulognormgpdcon
. Fitting function
flognormgpdcon
returns a simple list with the
following elements
call : | optim call |
x : | data vector x |
init : | pvector |
fixedu : | fixed threshold, logical |
useq : | threshold vector for profile likelihood or scalar for fixed threshold |
nllhuseq : | profile negative log-likelihood at each threshold in useq |
optim : | complete optim output |
mle : | vector of MLE of parameters |
cov : | variance-covariance matrix of MLE of parameters |
se : | vector of standard errors of MLE of parameters |
rate : | phiu to be consistent with evd |
nllh : | minimum negative log-likelihood |
n : | total sample size |
lnmean : | MLE of log-normal mean |
lnsd : | MLE of log-normal standard deviation |
u : | threshold (fixed or MLE) |
sigmau : | MLE of GPD scale (estimated from other parameters) |
xi : | MLE of GPD shape |
phiu : | MLE of tail fraction (bulk model or parameterised approach) |
se.phiu : | standard error of MLE of tail fraction |
See Acknowledgments in
fnormgpd
, type help fnormgpd
.
When pvector=NULL
then the initial values are:
MLE of log-normal parameters assuming entire population is log-normal; and
threshold 90% quantile (not relevant for profile likelihood for threshold or fixed threshold approaches);
MLE of GPD shape parameter above threshold.
Yang Hu and Carl Scarrott carl.scarrott@canterbury.ac.nz
http://www.math.canterbury.ac.nz/~c.scarrott/evmix
http://en.wikipedia.org/wiki/Lognormal_distribution
http://en.wikipedia.org/wiki/Generalized_Pareto_distribution
Scarrott, C.J. and MacDonald, A. (2012). A review of extreme value threshold estimation and uncertainty quantification. REVSTAT - Statistical Journal 10(1), 33-59. Available from http://www.ine.pt/revstat/pdf/rs120102.pdf
Hu, Y. (2013). Extreme value mixture modelling: An R package and simulation study. MSc (Hons) thesis, University of Canterbury, New Zealand. http://ir.canterbury.ac.nz/simple-search?query=extreme&submit=Go
Solari, S. and Losada, M.A. (2004). A unified statistical model for hydrological variables including the selection of threshold for the peak over threshold method. Water Resources Research. 48, W10541.
Other lognormgpd: flognormgpd
,
lognormgpdcon
, lognormgpd
Other lognormgpdcon: flognormgpd
,
lognormgpdcon
, lognormgpd
Other normgpdcon: fgngcon
,
fhpdcon
, fnormgpdcon
,
fnormgpd
, gngcon
,
gng
, hpdcon
,
hpd
, normgpdcon
,
normgpd
Other flognormgpdcon: lognormgpdcon
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 | ## Not run:
set.seed(1)
par(mfrow = c(2, 1))
x = rlnorm(1000)
xx = seq(-0.1, 10, 0.01)
y = dlnorm(xx)
# Continuity constraint
fit = flognormgpdcon(x)
hist(x, breaks = 100, freq = FALSE, xlim = c(-0.1, 10), ylim = c(0, 0.8))
lines(xx, y)
with(fit, lines(xx, dlognormgpdcon(xx, lnmean, lnsd, u, xi), col="red"))
abline(v = fit$u, col = "red")
# No continuity constraint
fit2 = flognormgpd(x, phiu = FALSE)
with(fit2, lines(xx, dlognormgpd(xx, lnmean, lnsd, u, sigmau, xi, phiu), col="blue"))
abline(v = fit2$u, col = "blue")
legend("topright", c("True Density","No continuity constraint","With continuty constraint"),
col=c("black", "blue", "red"), lty = 1)
# Profile likelihood for initial value of threshold and fixed threshold approach
fitu = flognormgpdcon(x, useq = seq(1, 5, length = 20))
fitfix = flognormgpdcon(x, useq = seq(1, 5, length = 20), fixedu = TRUE)
hist(x, breaks = 100, freq = FALSE, xlim = c(-0.1, 10), ylim = c(0, 0.8))
lines(xx, y)
with(fit, lines(xx, dlognormgpdcon(xx, lnmean, lnsd, u, xi), col="red"))
abline(v = fit$u, col = "red")
with(fitu, lines(xx, dlognormgpdcon(xx, lnmean, lnsd, u, xi), col="purple"))
abline(v = fitu$u, col = "purple")
with(fitfix, lines(xx, dlognormgpdcon(xx, lnmean, lnsd, u, xi), col="darkgreen"))
abline(v = fitfix$u, col = "darkgreen")
legend("topright", c("True Density","Default initial value (90% quantile)",
"Prof. lik. for initial value", "Prof. lik. for fixed threshold"),
col=c("black", "red", "purple", "darkgreen"), lty = 1)
## End(Not run)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.