VaR_ES_bounds {qrmtools} | R Documentation |
Compute the worst and best Value-at-Risk (VaR) and the best expected shortfall (ES) for given marginal distributions.
## Homogeneous case crude_VaR_bounds(level, qF, d = NULL, ...) VaR_bounds_hom(level, d, method = c("Wang", "Wang.Par", "dual"), interval = NULL, tol = NULL, ...) dual_bound(s, d, pF, tol = .Machine$double.eps^0.25, ...) ## Inomogeneous case ## Workhorses rearrange(X, tol = 0, tol.type = c("relative", "absolute"), n.lookback = ncol(X), max.ra = Inf, method = c("worst.VaR", "best.VaR", "best.ES"), sample = TRUE, is.sorted = FALSE, trace = FALSE, ...) block_rearrange(X, tol = 0, tol.type = c("absolute", "relative"), n.lookback = ncol(X), max.ra = Inf, method = c("worst.VaR", "best.VaR", "best.ES"), sample = TRUE, trace = FALSE, ...) ## User interfaces ## Rearrangement Algorithm RA(level, qF, N, abstol = 0, n.lookback = length(qF), max.ra = Inf, method = c("worst.VaR", "best.VaR", "best.ES"), sample = TRUE) ## Adaptive Rearrangement Algorithm ARA(level, qF, N.exp = seq(8, 19, by = 1), reltol = c(0, 0.01), n.lookback = length(qF), max.ra = 10*length(qF), method = c("worst.VaR", "best.VaR", "best.ES"), sample = TRUE) ## Adaptive Block Rearrangement Algorithm ABRA(level, qF, N.exp = seq(8, 19, by = 1), absreltol = c(0, 0.01), n.lookback = NULL, max.ra = Inf, method = c("worst.VaR", "best.VaR", "best.ES"), sample = TRUE)
level |
confidence level alpha for VaR and ES (e.g., 0.99). |
d |
dimension (number of risk factors; >= 2). For
|
qF |
|
method |
|
interval |
initial interval (a
|
tol |
|
tol.type |
|
n.lookback |
number of rearrangements to look back for deciding about numerical convergence. Use this option with care. |
s |
dual bound evaluation point. |
pF |
marginal loss distribution function. |
X |
( |
max.ra |
maximal number of (considered) column rearrangements
of the underlying matrix of quantiles (can be set to |
N |
number of discretization points. |
N.exp |
exponents of the number of discretization points
(a |
abstol |
absolute convergence tolerance epsilon
to determine the individual convergence, i.e., the change in the computed
minimal row sums (for |
reltol |
|
absreltol |
|
sample |
|
is.sorted |
|
trace |
|
... |
|
For d = 2
, VaR_bounds_hom()
uses the method of
Embrechts et al. (2013,
Proposition 2). For method = "Wang"
and method = "Wang.Par"
the method presented in McNeil et al. (2015, Prop. 8.32) is
implemented; this goes back to Embrechts et al. (2014, Prop. 3.1; note that
the published version of this paper contains typos for both bounds).
This requires
one uniroot()
and, for the generic method = "Wang"
,
one integrate()
. The critical part for the
generic method = "Wang"
is the lower endpoint of the initial
interval for uniroot()
. If the (marginal)
distribution function has finite first moment, this can be taken as
0. However, if it has infinite first moment, the lower endpoint has to
be positive (but must lie below the unknown root). Note that the upper
endpoint (1-alpha)/d also happens to be a
root and thus one needs a proper initional interval containing the
root and being stricticly contained in
(1-alpha)/d).
In the case of Pareto margins, Hofert et al. (2015) have
derived such an initial (which is used by
method = "Wang.Par"
).
Also note that the chosen smaller default tolerances for
uniroot()
in case of method = "Wang"
and
method = "Wang.Par"
are crucial for obtaining reliable
VaR values; see Hofert et al. (2015).
For method = "dual"
for computing worst VaR, the method
presented of Embrechts et al. (2013, Proposition 4) is implemented.
This requires two (nested) uniroot()
, and an
integrate()
. For the inner root-finding procedure to
find a root, the lower endpoint of the provided initial
interval
has to be “sufficiently large”.
Note that these approaches for computing the
VaR bounds in the homogeneous case are numerically non-trivial;
see the source code and vignette("VaR_bounds",
package = "qrmtools")
for more details. As a
rule of thumb, use method = "Wang"
if you have to (i.e., if the
margins are not Pareto) and method = "Wang.Par"
if you can (i.e.,
if the margins are Pareto). It is not recommended to use
(the numerically even more challenging) method = "dual"
.
Concerning the inhomogeneous case, rearrange()
is an auxiliary
function (workhorse). It is called by RA()
and ARA()
.
After a column rearrangement of X
, the tolerance between the
minimal row sum (for the worst VaR) or maximal row sum (for the best
VaR) or expected shortfall (obtained from the row sums; for the best
ES) after this rearrangement and the one of n.lookback
rearrangement steps before is computed and convergence determined.
For performance reasons, no input checking is done for
rearrange()
and it can change in future versions to (futher)
improve run time. Overall it should only be used by experts.
block_rearrange()
, the workhorse underlying ABRA()
,
is similar to rearrange()
in that it
checks whether convergence has occurred after every rearrangement by
comparing the change to the row sum variance from n.lookback
rearrangement steps back. block_rearrange()
differs from
rearrange
in the following ways. First, instead of single columns,
whole (randomly chosen) blocks (two at a time) are chosen and
oppositely ordered. Since some of the ideas for improving the speed of
rearrange()
do not carry over to block_rearrange()
, the
latter should in general not be as fast as the former.
Second, instead of using minimal or maximal row
sums or expected shortfall to determine numerical convergence,
block_rearrange()
uses the variance of the vector of row sums
to determine numerical convergence. By default, it targets a variance
of 0 (which is also why the default tol.type
is "absolute"
).
For the Rearrangement Algorithm RA()
, convergence of
\underline{s}_N and \overline{s}_N is determined if the
minimal row sum (for the worst VaR) or maximal row sum (for the best
VaR) or expected shortfall (obtained from the row sums; for the best ES)
satisfies the specified abstol
(so <= eps)
after at most max.ra
-many column rearrangements. This is different
from Embrechts et al. (2013) who use < eps and
only check for convergence after an iteration through all
columns of the underlying matrix of quantiles has been completed.
For the Adaptive Rearrangement Algorithm ARA()
and the Adaptive Block Rearrangement Algorithm ABRA()
,
convergence of \underline{s}_N and \overline{s}_N
is determined if, after at most max.ra
-many column
rearrangements, the (the individual relative tolerance)
reltol[1]
is satisfied and the
relative (joint) tolerance between both bounds is at most reltol[2]
.
Note that RA()
, ARA()
and ABRA()
need to evalute the
0-quantile (for the lower bound for the best VaR) and
the 1-quantile (for the upper bound for the
worst VaR). As the algorithms, due to performance reasons, can only
handle finite values, the 0-quantile and the 1-quantile need to be
adjusted if infinite. Instead of the 0-quantile,
the alpha/(2N)-quantile is
computed and instead of the 1-quantile the
alpha+(1-alpha)(1-1/(2N))-quantile
is computed for such margins (if the 0-quantile or the 1-quantile is
finite, no adjustment is made).
rearrange()
, block_rearrange()
, RA()
, ARA()
and ABRA()
compute \underline{s}_N and
\overline{s}_N which are, from a practical
point of view, treated as bounds for the worst (i.e., largest) or the
best (i.e., smallest) VaR or the best (i.e., smallest ES), but which are
not known to be such bounds from a theoretical point of view; see also above.
Calling them “bounds” for worst/best VaR or best ES is thus
theoretically not correct (unless proven) but “practical”.
The literature thus speaks of (\underline{s}_N, \overline{s}_N) as
the rearrangement gap.
More details not provided here can be found in the references listed below.
crude_VaR_bounds()
returns crude lower and upper bounds for
VaR at confidence level alpha for any
d-dimensional model with marginal quantile functions
specified by qF
.
VaR_bounds_hom()
returns the best and worst VaR at
confidence level alpha for d risks with equal
distribution function specified by the ellipsis ...
.
dual_bound()
returns the value of the dual bound D(s) as
given in Embrechts, Puccetti, Rüschendorf
(2013, Eq. (12)).
rearrange()
and block_rearrange()
return a
list
containing
bound
:computed \underline{s}_N or \overline{s}_N.
tol
:reached tolerance (i.e., the (absolute or
relative) change of the minimal row sum (for
method = "worst.VaR"
) or maximal row sum
(for method = "best.VaR"
) or expected shortfall (for
method = "best.ES"
) after the last rearrangement).
converged
:logical
indicating whether
the desired (absolute or relative) tolerance tol
has been
reached.
opt.row.sums
:vector
containing the
computed optima (minima for method = "worst.VaR"
; maxima
for method = "best.VaR"
; expected shortfalls for
method = "best.ES"
) for the row sums after each (considered)
rearrangement.
X.rearranged
:(N
, d
)-matrix
containing the rearranged X
.
RA()
returns a list
containing
bounds
:bivariate vector containing the computed \underline{s}_N and \overline{s}_N (the so-called rearrangement range) which are typically treated as bounds for worst/best VaR or best ES; see also above.
rel.ra.gap
:reached relative tolerance (also known as relative rearrangement gap) between \underline{s}_N and \overline{s}_N computed with respect to \overline{s}_N.
ind.abs.tol
:bivariate vector
containing
the reached individual absolute tolerances (i.e., the absolute change
of the minimal row sums
(for method = "worst.VaR"
) or maximal row sums
(for method = "best.VaR"
) or expected shortfalls
(for mehtod = "best.ES"
) for computing \underline{s}_N
and \overline{s}_N;
see also tol
returned by rearrange()
above).
converged
:bivariate logical
vector
indicating convergence of the computed \underline{s}_N and
\overline{s}_N (i.e., whether the desired tolerances were
reached).
num.ra
:bivariate vector containing the number of column rearrangments of the underlying matrices of quantiles for \underline{s}_N and \overline{s}_N.
opt.row.sums
:list
of length two containing
the computed optima (minima for method = "worst.VaR"
; maxima
for method = "best.VaR"
; expected shortfalls for
method = "best.ES"
) for the row sums after each
(considered) column rearrangement for the computed
\underline{s}_N and \overline{s}_N; see also
rearrange()
.
X
:initially constructed (N
, d
)-matrices
of quantiles for computing
\underline{s}_N and \overline{s}_N.
X.rearranged
:rearranged matrices X
(for
\underline{s}_N and \overline{s}_N).
ARA()
and ABRA()
return a list
containing
bounds
:see RA()
.
rel.ra.gap
:see RA()
.
tol
:trivariate vector
containing
the reached individual (relative for ARA()
; absolute for
ABRA()
) tolerances and the reached joint
relative tolerance (computed with respect to \overline{s}_N).
converged
:trivariate logical
vector
indicating individual convergence of the computed
\underline{s}_N (first entry) and \overline{s}_N
(second entry) and indicating joint convergence of the two bounds
according to the attained joint relative tolerance (third entry).
N.used
:actual N
used for computing
the (final) \underline{s}_N and \overline{s}_N.
num.ra
:see RA()
; computed for N.used
.
opt.row.sums
:see RA()
; computed for N.used
.
X
:see RA()
; computed for N.used
.
X.rearranged
:see RA()
; computed for N.used
.
Marius Hofert
Embrechts, P., Puccetti, G., Rüschendorf, L., Wang, R. and Beleraj, A. (2014). An Academic Response to Basel 3.5. Risks 2(1), 25–48.
Embrechts, P., Puccetti, G. and Rüschendorf, L. (2013). Model uncertainty and VaR aggregation. Journal of Banking \& Finance 37, 2750–2764.
McNeil, A. J., Frey, R. and Embrechts, P. (2015). Quantitative Risk Management: Concepts, Techniques, Tools. Princeton University Press.
Hofert, M., Memartoluie, A., Saunders, D. and Wirjanto, T. (2017). Improved Algorithms for Computing Worst Value-at-Risk. Statistics \& Risk Modeling or, for an earlier version, http://arxiv.org/abs/1505.02281.
Bernard, C., Rüschendorf, L. and Vanduffel, S. (2013). Value-at-Risk bounds with variance constraints. See http://papers.ssrn.com/sol3/papers.cfm?abstract_id=2342068.
Bernard, C. and McLeish, D. (2014). Algorithms for Finding Copulas Minimizing Convex Functions of Sums. See http://arxiv.org/abs/1502.02130v3.
vignette("VaR_bounds", package = "qrmtools")
for more example calls, numerical challenges
encoutered and a comparison of the different methods for computing
the worst (i.e., largest) Value-at-Risk.
### 1 Reproducing selected examples of McNeil et al. (2015; Table 8.1) ######### ## Setup alpha <- 0.95 d <- 8 theta <- 3 qF <- rep(list(function(p) qPar(p, shape = theta)), d) ## Worst VaR N <- 5e4 set.seed(271) system.time(RA.worst.VaR <- RA(alpha, qF = qF, N = N, method = "worst.VaR")) RA.worst.VaR$bounds stopifnot(RA.worst.VaR$converged, all.equal(RA.worst.VaR$bounds[["low"]], RA.worst.VaR$bounds[["up"]], tol = 1e-4)) ## Best VaR N <- 5e4 set.seed(271) system.time(RA.best.VaR <- RA(alpha, qF = qF, N = N, method = "best.VaR")) RA.best.VaR$bounds stopifnot(RA.best.VaR$converged, all.equal(RA.best.VaR$bounds[["low"]], RA.best.VaR$bounds[["up"]], tol = 1e-4)) ## Best ES N <- 5e4 # actually, we need a (much larger) N here (but that's time consuming) set.seed(271) system.time(RA.best.ES <- RA(alpha, qF = qF, N = N, method = "best.ES")) RA.best.ES$bounds stopifnot(RA.best.ES$converged, all.equal(RA.best.ES$bounds[["low"]], RA.best.ES$bounds[["up"]], tol = 5e-1)) ### 2 More Pareto examples (d = 2, d = 8; hom./inhom. case; explicit/RA/ARA) ### alpha <- 0.99 # VaR confidence level th <- 2 # Pareto parameter theta qF <- function(p, theta = th) qPar(p, shape = theta) # Pareto quantile function pF <- function(q, theta = th) pPar(q, shape = theta) # Pareto distribution function ### 2.1 The case d = 2 ######################################################### d <- 2 # dimension ## Explicit VaRbounds <- VaR_bounds_hom(alpha, d = d, qF = qF) # (best VaR, worst VaR) ## Adaptive Rearrangement Algorithm (ARA) set.seed(271) # set seed (for reproducibility) ARAbest <- ARA(alpha, qF = rep(list(qF), d), method = "best.VaR") ARAworst <- ARA(alpha, qF = rep(list(qF), d)) ## Rearrangement Algorithm (RA) with N as in ARA() RAbest <- RA(alpha, qF = rep(list(qF), d), N = ARAbest$N.used, method = "best.VaR") RAworst <- RA(alpha, qF = rep(list(qF), d), N = ARAworst$N.used) ## Compare stopifnot(all.equal(c(ARAbest$bounds[1], ARAbest$bounds[2], RAbest$bounds[1], RAbest$bounds[2]), rep(VaRbounds[1], 4), tolerance = 0.004, check.names = FALSE)) stopifnot(all.equal(c(ARAworst$bounds[1], ARAworst$bounds[2], RAworst$bounds[1], RAworst$bounds[2]), rep(VaRbounds[2], 4), tolerance = 0.003, check.names = FALSE)) ### 2.2 The case d = 8 ######################################################### d <- 8 # dimension ## Compute VaR bounds with various methods I <- crude_VaR_bounds(alpha, qF = qF, d = d) # crude bound VaR.W <- VaR_bounds_hom(alpha, d = d, method = "Wang", qF = qF) VaR.W.Par <- VaR_bounds_hom(alpha, d = d, method = "Wang.Par", shape = th) VaR.dual <- VaR_bounds_hom(alpha, d = d, method = "dual", interval = I, pF = pF) ## Adaptive Rearrangement Algorithm (ARA) (with different relative tolerances) set.seed(271) # set seed (for reproducibility) ARAbest <- ARA(alpha, qF = rep(list(qF), d), reltol = c(0.001, 0.01), method = "best.VaR") ARAworst <- ARA(alpha, qF = rep(list(qF), d), reltol = c(0.001, 0.01)) ## Rearrangement Algorithm (RA) with N as in ARA and abstol (roughly) chosen as in ARA RAbest <- RA(alpha, qF = rep(list(qF), d), N = ARAbest$N.used, abstol = mean(tail(abs(diff(ARAbest$opt.row.sums$low)), n = 1), tail(abs(diff(ARAbest$opt.row.sums$up)), n = 1)), method = "best.VaR") RAworst <- RA(alpha, qF = rep(list(qF), d), N = ARAworst$N.used, abstol = mean(tail(abs(diff(ARAworst$opt.row.sums$low)), n = 1), tail(abs(diff(ARAworst$opt.row.sums$up)), n = 1))) ## Compare stopifnot(all.equal(c(VaR.W[1], ARAbest$bounds, RAbest$bounds), rep(VaR.W.Par[1],5), tolerance = 0.004, check.names = FALSE)) stopifnot(all.equal(c(VaR.W[2], VaR.dual[2], ARAworst$bounds, RAworst$bounds), rep(VaR.W.Par[2],6), tolerance = 0.003, check.names = FALSE)) ## Using (some of) the additional results computed by (A)RA() xlim <- c(1, max(sapply(RAworst$opt.row.sums, length))) ylim <- range(RAworst$opt.row.sums) plot(RAworst$opt.row.sums[[2]], type = "l", xlim = xlim, ylim = ylim, xlab = "Number or rearranged columns", ylab = paste0("Minimal row sum per rearranged column"), main = substitute("Worst VaR minimal row sums ("*alpha==a.*","~d==d.*" and Par("* th.*"))", list(a. = alpha, d. = d, th. = th))) lines(1:length(RAworst$opt.row.sums[[1]]), RAworst$opt.row.sums[[1]], col = "royalblue3") legend("bottomright", bty = "n", lty = rep(1,2), col = c("black", "royalblue3"), legend = c("upper bound", "lower bound")) ## => One should use ARA() instead of RA() ### 3 "Reproducing" examples from Embrechts et al. (2013) ###################### ### 3.1 "Reproducing" Table 1 (but seed and eps are unknown) ################### ## Left-hand side of Table 1 N <- 50 d <- 3 qPar <- rep(list(qF), d) p <- alpha + (1-alpha)*(0:(N-1))/N # for 'worst' (= largest) VaR X <- sapply(qPar, function(qF) qF(p)) cbind(X, rowSums(X)) ## Right-hand side of Table 1 set.seed(271) res <- RA(alpha, qF = qPar, N = N) row.sum <- rowSums(res$X.rearranged$low) cbind(res$X.rearranged$low, row.sum)[order(row.sum),] ### 3.2 "Reproducing" Table 3 for alpha = 0.99 ################################# ## Note: The seed for obtaining the exact results as in Table 3 is unknown N <- 2e4 # we use a smaller N here to save run time eps <- 0.1 # absolute tolerance xi <- c(1.19, 1.17, 1.01, 1.39, 1.23, 1.22, 0.85, 0.98) beta <- c(774, 254, 233, 412, 107, 243, 314, 124) qF.lst <- lapply(1:8, function(j){ function(p) qGPD(p, shape = xi[j], scale = beta[j])}) set.seed(271) res.best <- RA(0.99, qF = qF.lst, N = N, abstol = eps, method = "best.VaR") print(format(res.best$bounds, scientific = TRUE), quote = FALSE) # close to first value of 1st row res.worst <- RA(0.99, qF = qF.lst, N = N, abstol = eps) print(format(res.worst$bounds, scientific = TRUE), quote = FALSE) # close to last value of 1st row ### 4 Further checks ########################################################### ## Calling the workhorses directly set.seed(271) ra <- rearrange(X) bra <- block_rearrange(X) stopifnot(ra$converged, bra$converged, all.equal(ra$bound, bra$bound, tolerance = 6e-3)) ## Checking ABRA against ARA set.seed(271) ara <- ARA (alpha, qF = qPar) abra <- ABRA(alpha, qF = qPar) stopifnot(ara$converged, abra$converged, all.equal(ara$bound[["low"]], abra$bound[["low"]], tolerance = 2e-3), all.equal(ara$bound[["up"]], abra$bound[["up"]], tolerance = 6e-3))