Case study: return moments and value at risk

The data

The series is 1,000 daily log returns $r_t = \ln(P_t / P_{t-1})$, as in case study 1.4.7. They are invented for illustration, generated with a fixed seed so that every cell gets the same numbers. Case studies 5.4.6 and 5.9.4 return to the same series.

The recipe gives them a mean of 0.04% a day, the mean daily log return that case study 1.4.7 annualised. It draws from a Student $t$ distribution with 4 degrees of freedom (lesson 5.4.4), which produces occasional very large days. Scaling by $\frac{0.01}{\sqrt{2}}$ sets the volatility near 1%, because that $t$ distribution has variance 2.

A position worth $V$ at the start of a day is worth $Ve^{r}$ at the end, so it loses $V(1 - e^{r})$. That is $V$ times the simple return with its sign flipped. Risk totals are built from such amounts, as case study 1.4.7 advised, not from log returns.

Four moments of a sample

Treat each of the $n = 1{,}000$ days as one equally likely outcome, with probability $\frac{1}{n}$. Then the sample is a distribution, and its moments are the expectations of lessons 5.3.1 to 5.3.3.

The volatility is close to the 1% of the recipe. The skewness is mildly negative, and the excess kurtosis is far from a normal's 0: the tails are fat.

The mean is the least certain number

The recipe's mean is $+0.04\%$ a day, yet the sample's is $-0.0061\%$. That is not an error: 1,000 days cannot pin a mean down closely.

Lesson 5.3.2's rules say how closely. The sample mean is $\frac{1}{n}(r_1 + \cdots + r_n)$. If the days are independent, the variances of the $n$ terms add, and the factor $\frac{1}{n}$ is squared. So the variance of the mean is $\frac{\sigma^2}{n}$.

With the recipe's volatility $\sigma = 1\%$ and $n = 1{,}000$, the mean's standard deviation is $\frac{1\%}{\sqrt{1000}} \approx 0.0316\%$. The sample mean sits 1.46 of these below 0.04%, an ordinary miss. Lessons 5.7.2 and 5.8.2 build on this.

What the kurtosis says

Excess kurtosis of 5.29 is easier to read as counts. Standardise each day to $z_t$ and count the days far from the mean. A normal distribution has $P(|Z| > 3) = 2N(-3)$, with $N$ the area of lesson 2.6.3, so it expects $1000 \times 2N(-3) \approx 2.7$ such days.

The middle is crowded too: 76.2% of days lie within one volatility of the mean, against 68.3% for a normal. Fat tails and a crowded middle come together, with fewer days in the shoulders between.

That fragility is lesson 5.3.3's warning made concrete. Worse, the recipe's Student $t$ with 4 degrees of freedom has an infinite kurtosis. So 5.29 does not estimate any true value, and a longer sample would give a different number.

Treat a sample kurtosis as a flag for fat tails, not as a precise input.

Value at risk as a quantile

Value at risk (VaR) at the 99% level over one day is the loss that is exceeded on only 1% of days. In return terms it starts from the 1% quantile of the daily returns, $q_{0.01}$: the return with 1% of days below it.

From the data

Sort the 1,000 returns. The worst ten are $-5.47\%$, $-5.37\%$, $-5.06\%$, $-4.65\%$, $-3.82\%$, $-3.70\%$, $-3.66\%$, $-3.57\%$, $-3.08\%$ and $-2.99\%$; the eleventh is $-2.65\%$. The 1% quantile lies between the tenth and the eleventh.

NumPy's np.quantile(r, 0.01) gives $-2.6562\%$. The historical VaR is then $1{,}000{,}000 \times (1 - e^{-0.026562})$, which is $26,212.

From a normal model

A normal model keeps only the mean $\bar{r}$ and the volatility $s$. If returns were normal, the chance of a day below $\bar{r} + zs$ would be $N(z)$ (lesson 2.6.3; lesson 5.4.3 develops this). So the 1% quantile is $\bar{r} + zs$ with $N(z) = 0.01$.

That $z$ is $-2.3263$. SciPy finds it as norm.ppf(0.01), the inverse of norm.cdf. The quantile is $\bar{r} - 2.3263s = -2.3495\%$, and the normal VaR is $23,221.

The normal figure is 11.4% below the historical one. It also fails its own test: it should be exceeded on about 10 of 1,000 days, and it was exceeded on 18.

Why kurtosis makes the normal VaR too small

The normal model matches the mean and the volatility exactly. What it cannot match is the shape. The data puts more weight far out in the tails and in the middle, and less in the shoulders.

VaR at 99% sits out in the tail, where the normal has too little weight.

So the error depends on the level. The table compares the two VaRs per $1,000,000 at six confidence levels:

At 90% and 95% the normal VaR is the larger, because the shoulders it imagines are emptier in the data. The two cross near 96%. Beyond that the normal is too small, and the gap widens: at 99.9% the data says 70% more.

A model checked only at 95% can look cautious there and still understate the 99% loss. That is why fat tails matter most for exactly the numbers risk managers report.

Expected shortfall

VaR says where the worst 1% of days begins, not how bad those days are. The expected shortfall (ES) at 99% is the average loss on the days beyond the VaR: here the 10 worst days.

For this series the 99% ES is $40,490, against a VaR of $26,212. A normal model with the same mean and volatility puts it near $26,573 (from a seeded simulation), far below the data. The Basel rules for banks' market risk use a 97.5% expected shortfall in place of 99% VaR, partly because it looks into the tail.

import numpy as np

rng = np.random.default_rng(50)
R = 0.0004 + 0.01 * rng.standard_t(4, size=1000) / np.sqrt(2)
SMALL = np.array([0.02, -0.1, 0.04, -0.05, 0.0, 0.1, -0.08, 0.01, -0.02, 0.06])


def _raises(f, *args):
    try:
        f(*args)
    except ValueError:
        return True
    return False


def test_historical_var():
    """The shared returns give a 99% historical VaR of about $26,212 per $1,000,000"""
    got = historical_var(R)
    if got is not None and abs(got - 26561.89) < 1:
        raise AssertionError("got 26,562: that treats the log return as a simple one. The loss is V(1 - e^q), np.expm1 helps")
    expect_close(got, 26212.2246829505, what="historical_var(r)", rel=1e-9)


def test_normal_var():
    """A normal with the same mean and volatility gives about $23,221"""
    expect_close(normal_var(R), 23221.176628094614, what="normal_var(r)", rel=1e-9)


def test_expected_shortfall():
    """The 99% expected shortfall, the average loss on the 10 worst days, is about $40,490"""
    expect_close(expected_shortfall(R), 40489.60362640658, what="expected_shortfall(r)", rel=1e-9)


def test_lower_level():
    """At 95% the normal VaR, about $16,493, is above the historical one, about $15,503"""
    expect_close(historical_var(R, 0.95), 15503.09856173464, what="historical_var(r, 0.95)", rel=1e-9)
    expect_close(normal_var(R, 0.95), 16492.55296317398, what="normal_var(r, 0.95)", rel=1e-9)


def test_small_sample():
    """Ten returns at the 80% level: VaR about $54,461, shortfall about $86,023 from the two worst days"""
    expect_close(historical_var(SMALL, 0.8), 54460.864109603746, what="historical_var(small, 0.8)", rel=1e-9)
    expect_close(expected_shortfall(SMALL, 0.8), 86023.11778870232, what="expected_shortfall(small, 0.8)", rel=1e-9)
    expect_close(normal_var(SMALL, 0.8), 50445.66311356407, what="normal_var(small, 0.8)", rel=1e-9)


def test_position_size():
    """A $2,500,000 position has 2.5 times the VaR: about $65,531"""
    expect_close(historical_var(R, 0.99, 2_500_000), 65530.56170737624, what="historical_var(r, 0.99, 2_500_000)", rel=1e-9)


def test_bad_level():
    """A level of 0, 1 or 1.5 raises ValueError in all three functions"""
    for f in (historical_var, normal_var, expected_shortfall):
        for level in (0, 1, 1.5):
            assert _raises(f, R, level), f"{f.__name__}(r, {level}) should raise ValueError"


def test_too_few_days():
    """50 days are too few for a 99% level: ValueError"""
    for f in (historical_var, normal_var, expected_shortfall):
        assert _raises(f, R[:50], 0.99), f"{f.__name__} with 50 days at 99% should raise ValueError"