This function displays traceplots of the scaling parameter from the proposal distribution of the adaptive MCMC scheme and the associated acceptance probability.
diagnostics(mcmc)
a graph of traceplots of the scaling parameter from the proposal distribution of the adaptive MCMC scheme and the associated acceptance probability.
An output of the fGEV
or fExtDep.np
function with method="Bayesian"
.
Simone Padoan, simone.padoan@unibocconi.it, https://faculty.unibocconi.it/simonepadoan/; Boris Beranger, borisberanger@gmail.com https://www.borisberanger.com; Giulia Marcon, giuliamarcongm@gmail.com
When mcmc
is the output of fGEV
then this corresponds to a marginal estimation and therefore diagnostics
will display in a first plot the value of
When mcmc
is the output of fExtDep.np
, then this corresponds to an estimation of the dependence structure following the procedure given in Algorithm 1 of Beranger et al. (2021). If the margins are jointly estimated with the dependence (step 1 and 2 of the algorithm) then diagnostics
provides trace plots of the corresponding scaling parameters (
Beranger, B., Padoan, S. A. and Sisson, S. A. (2021). Estimation and uncertainty quantification for extreme quantile regions. Extremes, 24, 349-375.
fExtDep.np
.
##################################################
### Example - Pollution levels in Milan, Italy ###
##################################################
if (FALSE) {
### Here we will only model the dependence structure
data(MilanPollution)
data <- Milan.winter[,c("NO2","SO2")]
data <- as.matrix(data[complete.cases(data),])
# Thereshold
u <- apply(data, 2, function(x) quantile(x, prob=0.9, type=3))
# Hyperparameters
hyperparam <- list(mu.nbinom = 6, var.nbinom = 8, a.unif=0, b.unif=0.2)
### Standardise data to univariate Frechet margins
f1 <- fGEV(data=data[,1], method="Bayesian", sig0 = 0.1, nsim = 5e+4)
diagnostics(f1)
burn1 <- 1:30000
gev.pars1 <- apply(f1$param_post[-burn1,],2,mean)
sdata1 <- trans2UFrechet(data=data[,1], pars=gev.pars1, type="GEV")
f2 <- fGEV(data=data[,2], method="Bayesian", sig0 = 0.1, nsim = 5e+4)
diagnostics(f2)
burn2 <- 1:30000
gev.pars2 <- apply(f2$param_post[-burn2,],2,mean)
sdata2 <- trans2UFrechet(data=data[,2], pars=gev.pars2, type="GEV")
sdata <- cbind(sdata1,sdata2)
### Bayesian estimation using Bernstein polynomials
pollut1 <- fExtDep.np(method="Bayesian", data=sdata, u=TRUE,
mar.fit=FALSE, k0=5, hyperparam = hyperparam, nsim=5e+4)
diagnostics(pollut1)
}
Run the code above in your browser using DataLab