Learn R Programming

ExtremalDep (version 1.0.0)

returns: Compute return values

Description

Predicts the probability of future simultaneous exceedances.

Usage

returns(x, summary.mcmc, y, plot = FALSE, data = NULL, ...)

Value

A numeric vector whose length equals the number of rows of the input y.

Arguments

x

An object of class ExtDep_npBayes, the output of the fExtDep.np function with method = "Bayesian".

summary.mcmc

The output of the summary.ExtDep_npBayes function.

y

A 2-column matrix of unobserved thresholds.

plot

Logical. If TRUE, then the returns are plotted using plot.ExtDep_npBayes.

data

As in plot.ExtDep_npBayes. Required if plot = TRUE.

...

Additional graphical parameters when plot = TRUE. See plot.ExtDep_npBayes.

Details

Computes for a range of unobserved extremes (larger than those observed in a sample) the pointwise mean from the posterior predictive distribution of such predictive values. The probabilities are calculated through

$$ P(Y_1 > y_1, Y_2 > y_2) = \frac{2}{k} \sum_{j=0}^{k-2} (\eta_{j+1} - \eta_j) \times \left( \frac{(j+1) B\!\left(\frac{y_1}{y_1+y_2} \mid j+2, k-j-1\right)}{y_1} - \frac{(k-j-1) B\!\left(\frac{y_2}{y_1+y_2} \mid k-j, j+1\right)}{y_2} \right) $$

where \(B(x \mid a,b)\) denotes the cumulative distribution function of a Beta random variable with shape parameters \(a, b > 0\). See Marcon et al. (2016, p. 3323) for details.

References

Marcon, G., Padoan, S. A. and Antoniano-Villalobos, I. (2016). Bayesian inference for the extremal dependence. Electronic Journal of Statistics, 10, 3310--3337.

Examples

Run this code
#########################################################
### Example 1 - daily log-returns between the GBP/USD ###
###             and GBP/JPY exchange rates            ###
#########################################################

if(interactive()){

  data(logReturns)

  mm_gbp_usd <- ts(logReturns$USD, start = c(1991, 3), end = c(2014, 12), frequency = 12)
  mm_gbp_jpy <- ts(logReturns$JPY, start = c(1991, 3), end = c(2014, 12), frequency = 12)

  # Detect seasonality and trend in the time series of maxima:
  seas_usd <- stl(mm_gbp_usd, s.window = "period")
  seas_jpy <- stl(mm_gbp_jpy, s.window = "period")

  # Remove the seasonality and trend:
  mm_gbp_usd_filt <- mm_gbp_usd - rowSums(seas_usd$time.series[,-3])
  mm_gbp_jpy_filt <- mm_gbp_jpy - rowSums(seas_jpy$time.series[,-3])

  # Estimation of margins and dependence
  mm_gbp <- cbind(as.vector(mm_gbp_usd_filt), as.vector(mm_gbp_jpy_filt))

  hyperparam <- list(mu.nbinom = 3.2, var.nbinom = 4.48)
  gbp_mar <- fExtDep.np(x = mm_gbp, method = "Bayesian", par10 = rep(0.1, 3),
                        par20 = rep(0.1, 3), sig10 = 0.0001, sig20 = 0.0001, k0 = 5,
                        hyperparam = hyperparam, nsim = 5e4)

  gbp_mar_sum <- summary(object = gbp_mar, burn = 3e4, plot = TRUE)

  mm_gbp_range <- apply(mm_gbp, 2, quantile, c(0.9, 0.995))

  y_gbp_usd <- seq(from = mm_gbp_range[1, 1], to = mm_gbp_range[2, 1], length = 20)
  y_gbp_jpy <- seq(from = mm_gbp_range[1, 2], to = mm_gbp_range[2, 2], length = 20)
  y <- as.matrix(expand.grid(y_gbp_usd, y_gbp_jpy, KEEP.OUT.ATTRS = FALSE))

  ret_marg <- returns(x = gbp_mar, summary.mcmc = gbp_mar_sum, y = y, plot = TRUE,
                      data = mm_gbp, xlab = "GBP/USD exchange rate", ylab = "GBP/JPY exchange rate")
}

#########################################################
### Example 2 - reproducing results from Marcon et al. ###
#########################################################

if (FALSE) {
  set.seed(1890)
  data <- evd::rbvevd(n = 100, dep = 0.6, asy = c(0.8, 0.3),
                      model = "alog", mar1 = c(1, 1, 1))

  hyperparam <- list(a.unif = 0, b.unif = .5, mu.nbinom = 3.2, var.nbinom = 4.48)
  pm0 <- list(p0 = 0.06573614, p1 = 0.3752118)

  mcmc <- fExtDep.np(method = "Bayesian", data = data, mar.fit = FALSE, k0 = 5,
                     pm0 = pm0, prior.k = "nbinom", prior.pm = "unif",
                     hyperparam = hyperparam, nsim = 5e5)

  w <- seq(0.001, 0.999, length = 100)
  summary.mcmc <- summary(object = mcmc, w = w, burn = 4e5, plot = TRUE)

  plot(x = mcmc, type = "A", summary.mcmc = summary.mcmc)
  plot(x = mcmc, type = "h", summary.mcmc = summary.mcmc)
  plot(x = mcmc, type = "pm", summary.mcmc = summary.mcmc)
  plot(x = mcmc, type = "k", summary.mcmc = summary.mcmc)

  y <- seq(10, 100, 2)
  y <- as.matrix(expand.grid(y, y))
  probs <- returns(out = mcmc, summary.mcmc = summary.mcmc, y = y, plot = TRUE)
}

Run the code above in your browser using DataLab