# NOT RUN {
### Example for estimating the parameters of the spliced Weibull-GPD
# }
# NOT RUN {
data("lossdat")
data=lossdat[[1]]$Loss
## starting values - method of moments
# Weibull distribution (shape, scale)
model <- function(x) c(F1 = x[2]*gamma(1+1/x[1])-mean(data),
F2 = x[2]^2*(gamma(1+2/x[1])-(gamma(1+1/x[1]))^2)-sd(data)^2)
weibullpar <- rootSolve::multiroot(f = model,start = c((sd(data)/mean(data))^(-1.086),
mean(data)/(gamma(1+1/((sd(data)/mean(data))^(-1.086))))))$root
# generalized Pareto distribution (tresh, beta, xi)
thresh=quantile(data,.9) # threshold
data2=data[which(data>thresh)]
model=function(x) c(F1=thresh+x[1]/(1-x[2])-mean(data2),
F2=x[1]^2/(1-x[2])^2/(1-2*x[2])-sd(data2)^2)
gpdpar=c(thresh,rootSolve::multiroot(f=model,start=c(mean(data2),1/mean(data2)))$root)
## parameters of the prior distribution
prior <- list(xi = c(gpdpar[3],sd(data)),tau = c(quantile(data,.9),sd(data)/10),
beta = c(gpdpar[2],sd(data)),wscale = c(weibullpar[2],sd(data)),
wshape = c(weibullpar[1],sd(data)))
## estimation of (shape, scale, thresh, beta, xi)
fitSplicedBayesWeibullGPD(data, prior = prior, proposal_scale = evmix::fweibullgpd(data,
method = "Nelder-Mead", pvector = c(weibullpar[1], weibullpar[2], gpdpar[1],
gpdpar[2], gpdpar[3]))$se, start = evmix::fweibullgpd(data,
method = "Nelder-Mead", pvector = c(weibullpar[1], weibullpar[2], gpdpar[1],
gpdpar[2], gpdpar[3]))$optim$par)
# }
Run the code above in your browser using DataLab