MR 2 – Non Parametric Approaches

A non-parametric risk estimate is read straight off a set of observed (or occasionally simulated) profit and loss outcomes, without first committing to a normal, Student t or any other named distribution. The data are allowed to describe the shape of the distribution by themselves, so features such as a heavy left tail or a lopsided body show up in the estimate automatically rather than being assumed away.

The price of that freedom is one large assumption that never goes away: the recent past must be a reasonable guide to the near future. If the coming weeks look nothing like the sample window, a non-parametric estimate will describe a world that no longer exists. Before relying on one of these methods, a risk manager therefore has to judge whether the historical sample is representative of the horizon that matters.

Most of the discussion focuses on value at risk (VaR) and expected shortfall (ES), but nothing is lost by that focus. A more general coherent risk measure is a weighted average of quantiles, so once a method can produce a quantile at any confidence level, the same machinery can be reused: choose a weighting function, divide the tail into a number of probability slices, estimate the quantile for each slice and average them with the chosen weights. One practical wrinkle is that the number of slices is not always a free choice here, because a fixed historical sample only offers so many distinct quantiles; bootstrapping, covered later, is a convenient way to enlarge the effective sample when that limit binds.

Where the methods sit on the parametric spectrum

Plain historical simulation is the purest member of the family: it simply re-prices today’s positions on past returns and ranks the outcomes. Bootstrapping and kernel smoothing keep that character but squeeze more out of the same sample, either by resampling it many times or by fitting a smooth density through it. The weighted versions go a step further and bolt parametric pieces onto the HS engine, such as a decay rate for old data or a GARCH volatility forecast. Once a formula of that kind is involved, the method is better described as semi-parametric, because it still lets the data decide the shape of the distribution while a model adjusts the data to current conditions. The appeal of these hybrids is that they keep the broad HS framework while recognising that the coming period may differ from the sample in predictable ways, for instance because volatility today is higher than it was a year ago.

Building the historical simulation P/L series

Historical simulation (HS) starts from the positions held today and asks what they would have earned on each day of a past sample. Suppose the portfolio contains n assets, the current amount invested in asset i is w_i, and R_i,t is the return (possibly after mapping) on asset i in subperiod t. The simulated P/L for subperiod t is simply the position-weighted sum of returns.

This series is a hypothetical one. It is generally different from the P/L the firm actually booked, because the real portfolio changed composition over the sample and mapping introduces approximations. Strictly, it is the P/L the current portfolio would have produced if it had been rebalanced at every close so that the amount held in each asset stayed at w_i, with gains taken out and losses topped up each day.

A useful by-product of this construction is that any number of positions collapses into a single P/L number per day. There is no covariance matrix to estimate and no curse of dimensionality to fight, which is why non-parametric methods are often the natural choice when a portfolio contains hundreds or thousands of risk factors.

Once the HS P/L series exists, basic HS amounts to sorting it. Plot the losses as a histogram, find the point beyond which the worst (1 − α) share of outcomes sits, and that point is the VaR at confidence level α. ES is then the average of the losses that lie beyond the VaR.

Take a sample of 1000 daily observations, which is about four years of data at 250 trading days a year. At 95% confidence the tail holds 5% of 1000, or 50 observations. The convention adopted here treats the VaR as the largest loss that is not in the tail, so the VaR is the 51st highest loss. On the illustrative data set this is 1.704, and the ES, the mean of the 50 biggest losses, is 2.196. Any point between the 50th and 51st largest losses (their mid-point, for instance) would also be defensible, but with a sample this large the choice hardly matters.

The data behind these figures were in fact drawn from a standard normal distribution, for which the true 95% VaR and ES are 1.645 and 2.063. The HS estimates miss by roughly 0.06 and 0.13 despite a sample of 1000, and the histogram is visibly ragged rather than bell-shaped. Sampling error of this size is the normal state of affairs for HS, which is why large samples are worth having wherever they are available.

In practice nobody needs the histogram itself. Sorting the losses from largest to smallest and picking the 51st entry, or using a spreadsheet function that returns the k-th largest value, gives the same answer.

General rule for the VaR observation

With n observations and confidence level α, the VaR is the negative of the r-th lowest P/L, where r = n(1 − α) + 1. For n = 1000 this gives the 51st lowest P/L at 95% and the 11th lowest at 99%. For n = 100 the 95% VaR is the sixth largest loss and the 96% VaR is the fifth largest.

Why HS estimates are imprecise

The VaR from basic HS is a single order statistic, so it depends entirely on where one particular observation happened to land. Shift a handful of tail draws and the estimate moves with them. ES is a little more stable because it averages several tail values, but with only 50 tail points in 1000 observations it still carries visible noise. Both estimates are also discrete: they can only take values that actually appear in the sample. These features motivate the refinements in the next sections, which either average across many resamples, smooth the empirical distribution, or attach an explicit measure of precision to the estimate.

A bootstrap builds many new data sets out of the one we already have. Each new sample has the same size as the original and is filled by drawing observations at random from the original, putting each one back after it is drawn. Because of the replacement step, a typical resample contains some observations twice or more and omits others altogether, so every resample tells a slightly different story about the tail.

Applied to HS, the recipe is straightforward. Compute the VaR from each resample, then take the mean of those resample VaRs as the bootstrapped point estimate. For ES, each resample gives its own ES (the average of the resample losses beyond the resample VaR), and the best ES estimate is the mean of these.

On the standard normal test data with 1000 resamples, the bootstrapped estimates come out at about 1.669 for VaR and 2.114 for ES. Both sit noticeably closer to the true values of 1.645 and 2.063 than the raw estimates of 1.704 and 2.196, which hints that averaging across resamples smooths away some of the noise in a single sample. The exact figures move slightly from run to run because the resampling itself is random.

The bootstrap has two further uses that recur through this lesson. First, the spread of the resample estimates measures how precise the estimate is, which leads directly to confidence intervals. Second, it lets the analyst generate as many quantile slices as needed when estimating a coherent risk measure, which removes the limit that a fixed sample places on the number of slices.

A small resampling illustration

Suppose an original sample has just five P/L values, labelled a to e. One resample, drawn with replacement, might come out as b, d, d, a, e, so c is left out and d appears twice. Another might be e, e, c, b, b. Each of these five-value sets is a legitimate alternative history consistent with the original data, and a VaR computed from each will differ. A real application uses the full HS sample (1000 values in the test data) and a large number of resamples, and the computer does the drawing with a uniform random number generator that picks an index between 1 and the sample size for every slot.

Basic HS wastes information and offers only a coarse grid of confidence levels. With 100 observations the 95% VaR is the sixth largest loss and the 96% VaR the fifth, but a 95.1% VaR has no observation of its own, even though it obviously lies between those two values. In general, n observations support VaRs at no more than n distinct confidence levels.

Non-parametric density estimation addresses both problems by treating the sample as draws from some unknown distribution and building a smooth stand-in for its density. Doing so forces some decisions that basic HS hides: how wide the bins should be, where they should be centred, and whether to summarise the data with a histogram, a naive estimator or a kernel. Theory favours kernels, and good choices on these questions can improve the VaR and ES estimates.

Turning a histogram into a surrogate density

The simplest version takes the histogram (with bar heights measuring relative frequency) and joins the mid-points at the top of neighbouring bars with straight lines. Using the bars themselves as a density would give an odd, step-shaped curve with sharp corners; the line through the mid-points is a more sensible approximation. Once it is drawn, the bars can be forgotten and the area under the line treated as a probability density.

Every confidence level now corresponds to a tail area under the curve, and the VaR is the point that splits the area into α on one side and 1 − α on the other. When that point falls between two loss observations, the computation amounts to a weighted average of the two neighbouring losses, with weights chosen so that the vertical line lands on the correct area split. A kernel goes one step further, replacing the straight segments with smooth curves.

Does the kernel matter in practice?

Empirical work by Butler and Schachter (1998) on real trading portfolios found that kernel-based VaRs differ only a little from basic HS figures, that different kernels give similar answers, and that adaptive Epanechnikov and adaptive Gaussian kernels did best where differences existed. Applying normal, box, triangular and Epanechnikov kernels to the standard normal test data gave the same VaR, 1.735, for all four. That figure is a little higher than the basic HS estimate of 1.704 and, since the true value is 1.645, slightly less accurate in this instance.

Two practical points follow. Kernels are better in theory, but crude straight-line interpolation is more transparent and easier to check. And because the P/L data already carry errors from mapping and data quality, there is a natural limit to how much precision any refinement can buy.

Plotting HS VaR and ES against the confidence level, say from 0.90 to 0.99, is easy and instructive. The VaR curve rises in an uneven, jumpy way because each point on it is a single loss observation and inherits that observation’s randomness. The ES curve rises more smoothly, because each point is an average of many tail losses and the averaging damps the noise. On the standard normal test data, VaR runs from a little above 1 at 90% to roughly 1.9 at 99%, while ES climbs from about 1.55 to about 2.4.

Changing the holding period is much harder. The methods so far give risk estimates for a holding period equal to the data frequency, so daily P/L produces a daily VaR. Nothing in theory stops us building a weekly or monthly HS P/L series and reading a VaR off that, but the data run out quickly.

The last column is added here to show how few tail observations remain.

With weekly P/L formed as the sum of five daily P/Ls, 1000 days yield only 200 weekly observations, and a monthly series has just 50. A 95% monthly ES would then rest on two or three numbers. Even if a much longer history existed, its older observations might say little about present conditions, so the data set itself puts a practical ceiling on how long a holding period HS can handle.

Putting the two dimensions together, a risk manager can in principle draw a surface showing VaR or ES against both confidence level and holding period. The confidence-level direction of that surface is cheap to fill in from one data set, whereas the holding-period direction is quickly starved of observations. That asymmetry is one reason why longer-horizon risk is often handled with scaling assumptions or with simulation methods such as filtered historical simulation, which builds multi-day paths from daily data rather than needing a long history of multi-day returns.

A point estimate of 1.704 says nothing about how far the true VaR might be from it. Two methods give a full distribution for the estimate and hence a confidence interval, and both are generally preferable to intervals built from quantile standard errors.

Order statistics approach

The theory of order statistics (OS), set out in detail in a later section, gives the whole distribution function of the VaR estimate rather than a single number. Percentiles of that distribution give the interval bounds, and its median, mean or mode offer alternative point estimates if wanted. Applied to the test data, the 5th and 95th percentiles of the distribution of the 95% VaR are 1.552 and 1.797, so there is 90% confidence that the true VaR lies in [1.552, 1.797].

For ES, the method maps from VaR to ES: each percentile point on the VaR distribution is matched with the corresponding percentile of the ES distribution. This gives a 90% interval for the 95% ES of [2.021, 2.224]. The OS approach can also be paired with a kernel rather than a histogram to squeeze more information from the data, as Butler and Schachter (1998) develop.

Bootstrap approach

The bootstrap route is more direct. Draw 1000 resamples, compute the 95% VaR of each and read the interval off the percentiles of those 1000 values. The resulting 90% interval is [1.554, 1.797]. The histogram of resample VaRs looks surprisingly broken up, since VaR can only take values equal to sample losses, yet repeating the exercise gives fairly stable bounds. Doing the same for ES, where each resample ES is the mean of its losses beyond the resample VaR, gives a histogram that is far better behaved (each ES is an average of tail values) and a 90% interval of [1.986, 2.271].

Bootstrap intervals use 1000 resamples. The width column is derived here as upper minus lower bound.

The two methods agree closely for VaR and reasonably well for ES, so either is a sensible choice in practice. Note that the true VaR of 1.645 and true ES of 2.063 both fall inside every interval in the table.

Basic HS uses a rolling window of the last n observations, and inside that window every observation counts equally. A return observed j days ago enters today’s P/L series and the next n series after it with exactly the same influence, then drops out completely once it is more than n days old. The weight is 1/n for any observation younger than n and zero for anything older, with nothing in between.

That rule is hard to defend on several grounds.

Equal weights ignore changing market conditions

Treating a calm day and a turbulent day as equally informative makes little sense when conditions differ. Shimko et al. (1998) point to natural gas, whose prices swing more in winter than in summer. An equal-weighted sample mixes the two seasons, so winter risk is understated and summer risk overstated; one fix is to give more weight to same-season observations when estimating a winter VaR.

Equal weights also make VaR slow to react. A sudden market crash may change a 95% VaR not at all, because a single extreme day merely joins the tail without moving the cut-off point. Only if losses continue over the following days does the VaR start to climb, which is rather like locking the stable door after the horse has gone. ES does respond to the first shock, because the new loss enters the tail average immediately, which illustrates why ES can be the more informative measure.

Pritsker (2001) notes a related blind spot. HS looks only at the loss tail, so a short equity position that enjoys large profits during a sharp sell-off registers no rise in risk, even though the higher volatility clearly raises the chance of a damaging rally. Pritsker also shows that, if VaR is set at confidence level α, an HS VaR rises on any given day only with probability 1 − α (the chance of a new tail loss), whereas under a GARCH(1,1) process at its long-run mean volatility the true VaR should rise about 32% of the time. Most of the time risk goes up, HS fails to register it.

The iid assumption and the ghost effect

Equal weighting implicitly treats each observation as equally likely and independent through time. Volatility clustering, the tendency of calm and stormy periods to bunch together, contradicts that assumption directly.

The abrupt drop from 1/n to zero creates a further problem. An observation that is n − 1 days old is treated as just as relevant as yesterday, yet one day later it is worth nothing, which throws information away for no good reason. The visible symptom is the ghost effect: a single large loss, or a small cluster of them, keeps the VaR elevated for n days and then, when the old loss leaves the window, the VaR drops sharply even though nothing in today’s market has changed. The fall is an artefact of window length and weighting, not a change in risk.

The weighted approaches that follow all keep the HS framework but bring in some parametric structure, which makes them semi-parametric. The first, proposed by Boudoukh, Richardson and Whitelaw (BRW, 1998), replaces equal weights with weights that shrink as observations age.

The weight on an observation one day old is w(1). The weight on a two-day-old observation is λw(1), on a three-day-old observation λ²w(1), and so on, where the decay parameter λ lies between 0 and 1. A value near 1 means slow decay; a value well below 1 means old data lose relevance quickly. Setting w(1) = (1 − λ)/(1 − λⁿ) makes the n weights sum to one, which gives the general weight.

To use it, sort the P/L observations from worst to best, place each weight beside its observation, and accumulate weights from the worst loss downwards until the running total passes 1 − α. The VaR is the loss at which that happens, and when the target percentile falls between two observations the loss can be interpolated.

Why age weighting helps

The approach has four main attractions. It generalises HS, since equal weighting is the limiting case λ → 1; if HS is like driving while watching only the rear-view mirror, equal weights are safe on a straight road and age weights cope with gentle bends. It makes VaR and ES respond faster to large recent losses, and handles clusters of them better, as Example 6 shows. It softens ghost effects, because an observation leaving the sample drops from λⁿw(1), which is far below 1/n for sensible values (here about 0.0000661 against 0.002), rather than from 1/n. Finally, since old observations fade anyway, the sample can be allowed to grow without limit, so no data are ever discarded and the jumps that cause ghost effects disappear entirely.

Where it falls short

Down-weighting reduces the effective sample size, and a run of large profits or losses can distort the implied risk profile. Pritsker shows that age-weighted VaR can still under-react to rising risk. Most troubling, the method is ad hoc: apart from the special case λ = 1, no asset return process is known for which BRW weighting is theoretically correct.

Hull and White (HW, 1998b) rescale past returns so that they reflect today’s volatility instead of the volatility that prevailed when they occurred. If a market moves 1.5% a day now but moved only 1% a day a month ago, raw month-old returns understate the moves to expect tomorrow, and HS based on them understates risk. If last month’s volatility was 2% a day, the same raw data overstate tomorrow’s moves and HS overstates risk.

To forecast VaR for day T, let r_t,i be the historical return on asset i on day t, σ_t,i the GARCH or EWMA volatility forecast for that day made at the close of day t − 1, and σ_T,i the latest forecast of volatility for asset i. Each historical return is replaced by its volatility-adjusted counterpart.

Returns from calmer periods are scaled up when current volatility is higher, and returns from stormier periods are scaled down when it is lower. The HS P/L series is then rebuilt from the adjusted returns and VaR or ES is estimated in the usual equal-weighted way. Volatility weighting obviously requires estimates of both current and past volatility to be available.

The adjustment is applied to every day in the sample, so a single historical series can be stretched in some places and compressed in others. The table shows four past days for one asset, rescaled to a current volatility forecast of 1.2% a day.

Original illustrative figures. Scaling factor = 1.2 divided by the day’s volatility forecast.

The worst raw return, −2.6% on day 3, shrinks to −1.56% because it came from a stormy period, while the milder −1.8% on day 1 grows into the worst adjusted loss, −2.7%, because it came from a calm one. The ranking of losses, and therefore the VaR, can change once volatility is taken into account.

Advantages over equal and age weighting

The HW method handles volatility changes directly and naturally, whereas equal-weighted HS ignores them and age weighting treats them in a crude, restrictive way. Its risk estimates track current volatility and so bring GARCH forecast information into HS. Because adjusted returns can be larger than any return actually observed, VaR and ES can exceed the largest historical loss, something plain HS can never do. HW also report empirical evidence that their VaRs beat the BRW approach.

The method extends easily. It can be combined with age weighting to sharpen responses to large losses and reduce ghost effects, and confidence intervals can be produced by applying order statistics or the bootstrap to the HW-adjusted P/L rather than the raw HS P/L.

Volatility weighting corrects each asset separately, but co-movements can shift as well. Correlation weighting, based on a suggestion by Duffie and Pan (1997), adjusts historical returns so that they carry today’s correlation structure. It is a little more involved than volatility weighting and is easiest to follow assuming any Hull-White volatility adjustment has already been made.

Suppose there are m positions and R is the (possibly volatility-adjusted) 1 × m vector of historical returns for some period. Those returns reflect an m × m covariance matrix Σ, which factorises as σCσ^T, where σ is the diagonal matrix of volatilities and C is the historical correlation matrix. The goal is to turn R into R̄, which carries the current correlation matrix C̄.

If both correlation matrices are positive definite, each has a matrix square root obtainable from a Choleski (Cholesky) decomposition: A for the historical matrix and Ā for the current one. Both return vectors can then be written as the relevant square root multiplying the same uncorrelated noise vector ε.

Inverting the first line gives ε = A⁻¹R, and substituting into the second yields the adjusted returns.

The adjusted returns carry the current correlation matrix and, combined with volatility weighting, the full current covariance matrix. This makes correlation weighting a broad generalisation of the HW approach.

Filtered historical simulation, developed by Barone-Adesi and co-authors in a series of papers from 1998 onwards, is a semi-parametric bootstrap. It bootstraps returns inside a conditional volatility model such as GARCH, so the bootstrap keeps HS free of distributional assumptions while the volatility model supplies a careful treatment of changing volatility.

The four steps for a one-day VaR on a single asset

Step 1, fit a volatility model. The model must be rich enough to capture the main features of the data. Barone-Adesi and colleagues favour an asymmetric GARCH (AGARCH), which captures time-varying volatility and volatility clustering and also lets negative and positive surprises affect volatility differently, the so-called leverage effect.

The mean daily return μ can often be ignored for volatility purposes. The term γ is what produces asymmetry, since it makes the effect of last period’s surprise depend on whether that surprise was positive or negative.

Step 2, standardise. Use the model to forecast volatility for each day in the sample, then divide each realised return by its volatility forecast. The resulting standardised returns should be close to iid and so suitable for HS.

Step 3, bootstrap and rescale. Draw a large number M of standardised returns at random, with replacement, and multiply each by the AGARCH forecast of tomorrow’s volatility. This gives M simulated returns, each scaled to current market conditions.

Step 4, read off the VaR. Each simulated return implies a possible end-of-day portfolio value and loss, and the VaR is the loss at the chosen confidence level. ES can be estimated from the same simulated losses.

Multiple assets and longer horizons

For a multi-asset portfolio, fit a multivariate GARCH or AGARCH, standardise the whole vector of asset returns, and let each bootstrap draw pick an entire day’s vector rather than a single portfolio return. Because assets are drawn together from the same day, existing correlations are preserved without specifying any multivariate density.

For a longer horizon, draw a standardised return to simulate tomorrow, feed it back into the model to update the volatility forecast for the following day, draw again, and continue until the horizon is reached. That path produces one simulated P/L; repeating the whole path many times builds the P/L histogram from which VaR is read.

Attractions and caveats

FHS joins non-parametric HS with a sophisticated volatility model, so it reflects changing conditions. It runs quickly even for large portfolios. Like HW, it can produce VaR and ES larger than the biggest historical loss. It keeps the correlation structure without needing a covariance matrix or a conditional return distribution, can be adapted to autocorrelation or cross-correlation, and can give confidence intervals when combined with OS or the bootstrap. Evidence suggests it performs well.

The evidence is not one-sided, though. Pritsker’s (2001) simulation study tentatively finds that FHS may under-weight extreme observations and time-varying correlations, and that FHS VaR might understate true 10-day VaR by about 10%. Results from Barone-Adesi et al. on real data conflict with the last finding, so the overall verdict on FHS is mixed. Switching from VaR to ES partly addresses the first concern; a richer cross-correlation structure addresses the second.

Confidence intervals from FHS

Because FHS already contains a bootstrap, extending it to produce a confidence interval is easy. If M estimates of VaR are wanted, simulate M batches of filtered P/L values, for example 100M or 1000M values split into M sets of 100 or 1000. Each set yields one HS VaR, and the histogram of the M VaRs gives the interval bounds in the usual percentile way. The OS route also works, since it only needs a set of P/L values paired with probabilities, which a P/L series passed through the FHS filter provides.

Advantages

• They are intuitive and conceptually simple, and most can be built in a spreadsheet with modest effort.
• Making no parametric assumption about P/L, they accommodate fat tails, skewness and other non-normal features that trouble parametric models.
• In principle they handle any kind of position, derivatives included.
• Practitioners widely regard HS as working fairly well, even though formal empirical evidence is mixed.
• They avoid the burdens parametric methods face in high dimensions: no covariance matrices and no curse of dimensionality.
• The data are often readily available, from public vendors such as Bloomberg or as a by-product of marking positions to market.
• Results are easy to explain to senior management, supervisors and rating agencies.
• Confidence intervals for VaR and ES are easy to obtain.
• They can be refined with parametric add-ons (age, volatility and correlation weighting, and FHS) that turn them into semi-parametric methods.

Disadvantages

The central weakness is near-total dependence on the historical data set, along with the practical challenge of obtaining reliable, timely series for every current position (often hard for emerging-market holdings). Several problems follow:

• A quiet sample period yields VaR and ES that are too low for current risks; a volatile one yields estimates that are too high.
• Structural shifts, such as a permanent change in exchange rate risk, take time to show up, and sudden turbulence registers slowly.
• Extreme losses that are unlikely to recur can dominate the estimates.
• Most methods suffer, to some extent, from ghost or shadow effects.
• Plausible events that simply did not happen in the sample receive no weight.
• Estimates are constrained by the largest historical loss. Simple HS cannot extrapolate beyond it at all; more refined versions loosen the constraint but do not remove it, so non-parametric methods are weak at handling extremes, especially with small or medium samples.

Many of these can be eased: semi-parametric adjustments deal with volatility, turbulence and correlation shifts, and age weighting with a growing sample deals with ghost effects.

Summary table compiled for this lesson from the discussion above.

How long should the sample window be?

Precision needs a reasonably large sample, and a common rule of thumb is at least two years of daily data, roughly 500 observations at 250 trading days a year, and often more. A very long window brings its own costs, because the longer it is:

• the more stale data it contains;
• the longer one-off past events distort results and the longer ghost effects take to clear;
• the more current market news is diluted by old observations, making estimates less responsive; and
• the harder data collection becomes, particularly for new or emerging-market instruments with short histories that are hard to proxy.

Overall, non-parametric methods have a reasonable track record, often beat parametric methods built on simple assumptions such as normality, and work fairly well when markets are reasonably stable. They should still be supplemented with other approaches and, wherever possible, with stress testing of what-if scenarios; they should never be relied on alone.

Each observation in a sample of n P/L values can be viewed as a VaR estimate at some implied confidence level. With n = 1000, the negative of the 51st smallest P/L is a 95% VaR and the negative of the 11th smallest is a 99% VaR. The general rule takes the α VaR as the negative of the r-th lowest observation with r = n(1 − α) + 1. The r-th lowest value in a sample is the r-th order statistic, and the distribution of order statistics is well understood.

Suppose the observations x₁, …, x_n come from a known distribution function F(x), with order statistics x₍₁₎ ≤ x₍₂₎ ≤ … ≤ x_(n). The number of observations that do not exceed a fixed value x follows a binomial distribution.

The r-th order statistic is at or below x exactly when at least r observations are, so its distribution function is a binomial tail sum.

G_r(x) delivers both a point estimate and an interval. Its median (the 50th percentile) is a natural best estimate of VaR, and its lower and upper percentiles bound the confidence interval. The calculation is accurate, easy in a spreadsheet, and works for any F(x), whether normal, t or an empirical distribution. To apply it, specify F(x), set the parameters and evaluate G_r.

Results for a standard normal

Setting F(x) to the standard normal gives OS estimates of normal VaR. The table shows the 5th, 50th and 95th percentiles of the VaR distribution alongside the conventional formula VaR. The final row divides the interval width by the median; widths are reported to the precision of the published bounds.

Bounds are the 5th and 95th percentiles of the VaR distribution from the order statistic formula.

Three patterns stand out. The interval is wide for small samples (42.2% of the median at n = 100) and narrows steadily as n grows (4.2% at n = 10 000). The median converges to the conventional 1.645 as n increases. And the interval is a little wider, relative to the median, at the more extreme 99% level than at 90% or 95%. Recomputing the bounds from the binomial formula reproduces every entry to within 0.002.

The same method works with a t, Gumbel or any other parametric distribution in place of the normal, and with an empirical distribution, which gives intervals for HS VaR. It is fully general across parametric and non-parametric VaR and ES, reasonable even in small samples because it does not lean on asymptotic theory (although small samples are naturally less precise), and easy to implement. It is also better than intervals built from quantile standard errors, which rely on asymptotics and force symmetric intervals, a real problem for extreme VaRs and ESs.

The bootstrap replaces mathematical derivation with simulated resampling from the data at hand, which makes it possible to judge the accuracy of an estimator without strong parametric assumptions or closed-form interval formulas. The idea has old roots but was developed and popularised by Bradley Efron, who also named it after pulling oneself up by one’s own bootstraps. It is generally easier, more powerful and more reliable than traditional interval methods, and it can also supply alternative point estimates. It is preferred to the older jackknife, which re-estimates the parameter leaving out one observation at a time and can approximate poorly when the estimator is a non-smooth function of the data.

Why closed-form intervals fall short

Suppose we estimate a parameter θ (a mean, variance, quantile and so on) with a sample estimator θ̂. Getting θ̂ is easy; getting an interval for it usually requires theory that exists only in special cases. For a variance, if the data are normal then (n − 1)s²/σ² follows a chi-squared distribution with n − 1 degrees of freedom. Writing χ²_α,n−1 for its α point, the 90% interval for that ratio and the implied interval for σ² are as follows.

For instance, with n = 30 and s² = 4, the chi-squared points with 29 degrees of freedom are 17.708 and 42.557, so the interval for σ² is [116 / 42.557, 116 / 17.708] = [2.726, 6.551]. Drop the normality assumption and the distribution of s² is usually unknown, so no closed-form interval is available. Medians, correlations and tail probabilities pose similar difficulties (quantiles are an exception, since order statistics give their distribution).

The basic procedure and the percentile interval

Start with an original sample of size n. Draw a resample of the same size using a uniform random number between 1 and n to pick each observation, returning it to the pool after each draw. Estimate θ from the resample, giving θ̂^B(i). Repeat B times. The B resample estimates form a bootstrapped sample of the estimator, and if the interval should cover the central 1 − 2α of probability, its ends are the α and 1 − α quantiles of that sample.

This is the method behind the [1.554, 1.797] VaR interval earlier. It needs no parametric theory, asymptotic or otherwise. Bootstrapping straight from the data in this way is the non-parametric version; a parametric bootstrap resamples from an assumed distribution, and FHS is a semi-parametric example.

The bias-corrected and accelerated (BCa) interval

The plain percentile interval struggles when the estimator is biased. The bias-corrected and accelerated approach improves on it in both theory and practice by swapping the α and 1 − α quantiles for adjusted levels α₁ and α₂.

When â and ẑ₀ are both zero, BCa collapses to the percentile interval. Otherwise they shift the end-points. The bias correction ẑ₀ is the standard normal inverse of the proportion of bootstrap estimates below the original θ̂. The acceleration â measures how fast the standard error of θ̂ changes with the true θ, acting as a skewness correction, and can be estimated from an initial bootstrap or jackknife run.

BCa is more accurate than the percentile interval and is transformation-respecting: if the parameter is transformed, say from a variance to a standard deviation, the interval end-points transform correctly as well.

Bootstrap point estimates and bias

The mean of the B resample estimates is an alternative point estimator that is often better than the raw θ̂, and it gives a direct estimate of bias.

A bias estimate can be used to correct the estimator or prediction errors, but it can have a large standard error, in which case the corrected estimate may be noisier than the biased one. Variance-reduction methods such as importance sampling can improve accuracy for a given computational effort, and with modern computers the extra computation is no longer a serious drawback.

Bootstrap estimates usually have little bias but can have substantial variance, which comes from two places: sampling variability (we hold a sample of size n, not the population) and resampling variability (we take B resamples, not infinitely many). The estimated standard error of θ̂ is the standard deviation of the resample estimates.

The standard error estimate is itself random. Its variance splits into the variance of its expected value plus the expected value of its variance, and Efron and Tibshirani (1993) rewrite this in terms of the moments of the bootstrap distribution. Closed-form results exist only when θ is the mean, and they simplify further only under normality, so they are of limited use. Two brute-force alternatives work generally: a jackknife-after-bootstrap, which bootstraps first and then jackknifes the bootstrapped data, and a double bootstrap, which bootstraps a sample of standard error estimates and measures their spread. Trying different values of B then shows how many resamples bring the variance down to an acceptable level.

The Andrews and Buchinsky rule for choosing B

If the main question is simply how large B must be, Andrews and Buchinsky (1997) offer a neater route. Treat the standard error from infinitely many resamples, s_∞, as the ideal. Choose a target probability τ close to 1 and a tolerance, bound, on the percentage deviation of s_B from s_∞. B is chosen so that s_B lies within that tolerance with probability τ.

For large B the required number of resamples is approximately proportional to 2500(κ − 1) divided by bound squared, scaled by a critical value c_τ that depends on the target probability, where κ is the kurtosis of the bootstrapped estimates.

Because κ is unknown, Andrews and Buchinsky propose a three-step procedure:

• Assume κ = 3, which turns 2500(κ − 1) into 5000, and compute a preliminary B₀ = int(5000c_τ/bound²), where int(a) is the smallest integer at least as large as a.
• Run B₀ resamples and estimate the kurtosis κ̂ of the bootstrapped estimates.
• Compute B₁ = 2500(κ̂ − 1)c_τ/bound² and use max(B₀, B₁) resamples.

The procedure does not estimate the variance of s_B directly, but it does something more practical: it sets B to achieve a target precision, and it works for any parameter and any distribution of the bootstrapped estimates, unlike the mean-only and normal-only formulas. It can be refined further, for example by correcting the kurtosis estimate for bias, and an analogous method targets the accuracy of confidence intervals.

A fat-tailed bootstrap distribution needs more resamples. With κ̂ = 6, the final count is B₁ = 2500 × 5 × c_τ/bound², which is 2.5 times the preliminary B₀, because 12 500 / 5000 = 2.5.

Bootstrapping when observations are not independent

The main weakness of the standard bootstrap is that it assumes observations are independent through time, and it can mislead when they are not. There are three common fixes:

• Parametric filtering. Model the dependence, for example with GARCH, and bootstrap the residuals, which should be independent. This requires identifying and estimating the right model, which brings model and parameter risk.
• Block bootstrap. Split the data into non-overlapping blocks of equal length and resample whole blocks. This can whiten the data, since observations either side of a block boundary are treated as independent; refinements such as random block lengths help but add complexity.
• Dependent selection probabilities. Let the chance of picking an observation depend on what was just picked, so that if observation i has just been drawn from chronologically ordered data, observation i + 1 is more likely than most others to come next.