ExtremalDep (version 0.0.4-2)

summary_ExtDep: Summary of MCMC algorithm.

Description

This function computes summaries on the posterior sample obtained from the adaptive MCMC scheme for the non-parametric estimation of a bivariate dependence structure.

Usage

summary_ExtDep(object, mcmc, burn, cred=0.95, plot=FALSE, ...)

Value

The function returns a list with the following objects:

k.median, k.up, k.low:

Posterior median, upper and lower bounds of the CI for the estimated Bernstein polynomial degree \(\kappa\);

h.mean, h.up, h.low:

Posterior mean, upper and lower bounds of the CI for the estimated angular density \(h\);

A.mean, A.up, A.low:

Posterior mean, upper and lower bounds of the CI for the estimated Pickands dependence function \(A\);

p0.mean, p0.up, p0.low:

Posterior mean, upper and lower bounds of the CI for the estimated point mass \(p_0\);

p1.mean, p1.up, p1.low:

Posterior mean, upper and lower bounds of the CI for the estimated point mass \(p_1\);

A_post:

Posterior sample for Pickands dependence function;

h_post:

Posterior sample for angular density;

eta.diff_post:

Posterior sample for the Bernstein polynomial coefficients (\(\eta\) parametrisation);

beta_post:

Posterior sample for the Bernstein polynomial coefficients (\(\beta\) parametrisation);

p0_post, p1_post:

Posterior sample for point masses \(p_0\) and \(p_1\);

w:

A vector of values on the bivariate simplex where the angular density and Pickands dependence function were evaluated;

burn:

The argument provided;

If the margins were also fitted, the list given as object would contain mar1 and mar2 and the function would also output:

mar1.mean, mar1.up, mar1.low:

Posterior mean, upper and lower bounds of the CI for the estimated marginal parameter on the first component;

mar2.mean, mar2.up, mar2.low:

Posterior mean, upper and lower bounds of the CI for the estimated marginal parameter on the second component;

mar1_post:

Posterior sample for the estimated marginal parameter on the first component;

mar2_post:

Posterior sample for the estimated marginal parameter on the second component;

Arguments

object

A vector of values on \([0,1]\). If missing, a regular grid of length \(100\) is considered.

mcmc

An output of the fExtDep.np function with method="Bayesian".

burn

A positive integer indicating the burn-in period.

cred

A value in \([0,1]\) indicating the level of the credibility intervals to be computed.

plot

A logical value; if TRUE a summary of the estimated dependence is displayed by calling the plot_ExtDep.np function with type="summary".

...

Additional graphical parameters for plot_ExtDep.np when plot=TRUE.

Author

Simone Padoan, , https://faculty.unibocconi.it/simonepadoan/; Boris Beranger, https://www.borisberanger.com

Details

For each value say \(\omega \in [0,1]\) given, the complement \(1-\omega \) is automatically computed to define the observation \((\omega,1-\omega)\) on the bivariate unit simplex.

It is obvious that the value of burn must be greater than the number of iterations in the mcmc algorithm. This can be found in mcmc.

See Also

fExtDep.np.

Examples

####################################################
### 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.0001, 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.0001, 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)
pollut1_sum <- summary_ExtDep(mcmc=pollut1, burn=3e+4, plot=TRUE)
}