The annualised volatility¶
A variance per period as a volatility a year, sigma = sqrt(v x n): the standard deviation of the log return over a year, from the variance over one sampling period and the number of periods a year, the returns of consecutive periods independent, so their variances add.
Both dimensionless: a variance per period (an EWMA's, a realised one's) and a count of periods a year (Hull's 252 trading days; a game's samples a day times its days a year). The volatility is the spread lusine.finance.d1 and the option kernels take over a year. Source: Hull, Options, Futures, and Other Derivatives, ch. 15 (the Black-Scholes-Merton model), Estimating volatility from historical data: sigma_annual = sigma_period sqrt(periods a year). Checked against a reference in numpy, which has its own tests.
- Source
- Hull, Options, Futures, and Other Derivatives, ch. 15 (the Black-Scholes-Merton model), Estimating volatility from historical data.
- Proof
- reference, invariant, limit.
annualised_volatility.trait
namespace lusine.finance;
float annualised_volatility(float variance, float periods)
{
return sqrt(variance * periods);
}
A period's variance as a volatility a year, sqrt(v n), the periods' returns independent so their variances add (Hull, ch. 15, Estimating volatility from historical data).
import os, sys
root = os.getcwd()
while not os.path.isfile(os.path.join(root, 'docsite.json')):
root = os.path.dirname(root)
sys.path.insert(0, os.path.join(root, 'tools'))
import numpy as np
import matplotlib.pyplot as plt
from lusine import trait_run
F32 = np.float32
U = float(np.finfo(F32).eps) / 2 # binary32's unit roundoff, 2^-24
The kernel¶
The law as the library states it, and as this page runs it.
print(trait_run.source('finance/annualised_volatility.trait'))
// lusine.finance.annualised_volatility: a variance per period as a volatility a year, sigma = sqrt(v x n): the standard
// deviation of the log return over a year, from the variance over one sampling period and the number of periods a
// year, the returns of consecutive periods independent, so their variances add. Both dimensionless: a variance per
// period (an EWMA's, a realised one's) and a count of periods a year (Hull's 252 trading days; a game's samples a day
// times its days a year). The volatility is the spread lusine.finance.d1 and the option kernels take over a year.
// Source: Hull, Options, Futures, and Other Derivatives, ch. 15 (the Black-Scholes-Merton model), Estimating
// volatility from historical data: sigma_annual = sigma_period sqrt(periods a year).
// Checked against a reference in numpy, which has its own tests.
namespace lusine.finance;
float annualised_volatility(float variance, float periods)
{
return sqrt(variance * periods);
}
The program¶
A lattice of estimates, each a variance per period and the periods a year, both the host's; main() binds the kernel over them.
VOL = """namespace vol;
lusine.version version = { { number = 1 } };
struct estimate { host float variance; host float periods; float sigma; };
lattice<estimate, 1> estimates(spacing = 1 m, boundary = wall);
void main() { estimates.sigma = lusine.finance.annualised_volatility(estimates.variance, estimates.periods); }
"""
print(VOL)
def annualised(variance, periods):
"""Each estimate's volatility a year."""
n = max(np.size(variance), np.size(periods))
run = trait_run.Program(['finance/annualised_volatility.trait'], VOL, {'vol.estimates': [n]})
run.store('vol.estimates', 'variance', np.broadcast_to(np.asarray(variance, dtype=float), (n,)))
run.store('vol.estimates', 'periods', np.broadcast_to(np.asarray(periods, dtype=float), (n,)))
return run.steps('vol.main', 1, [('vol.estimates', 'sigma')])[0][0]
namespace vol;
lusine.version version = { { number = 1 } };
struct estimate { host float variance; host float periods; float sigma; };
lattice<estimate, 1> estimates(spacing = 1 m, boundary = wall);
void main() { estimates.sigma = lusine.finance.annualised_volatility(estimates.variance, estimates.periods); }
Reference: Hull's trading year¶
A daily standard deviation of 1% over 252 trading days is an annual volatility of 15.87% (Hull, Options, Futures, and Other Derivatives, ch. 15): the kernel's, from the daily variance, to the book's four digits.
sigma = float(annualised(0.01 ** 2, 252.0)[0])
assert abs(sigma - 0.1587) < 0.0001, sigma
print(f'{sigma * 100:.4f}% a year')
15.8745% a year
Reference: the float64 reference, within the square root's rounding¶
Three hundred variances from 1e-6 to 1e-2 over 1 to 100,000 periods a year: each within 2 ulps of the float64 reference's, since the product rounds once and the square root is correctly rounded.
from lusine.kernels.finance.annualised_volatility import annualised_volatility as reference
rng = np.random.default_rng(417)
v = rng.uniform(1e-6, 1e-2, 300)
n = rng.integers(1, 100000, 300).astype(float)
got = annualised(v, n)
want = reference(F32(v).astype(float), F32(n).astype(float))
assert np.all(np.abs(got - want) <= 4 * U * want), np.max(np.abs(got / want - 1))
print('worst relative difference', float(np.max(np.abs(got / want - 1))))
worst relative difference 7.340632024455118e-08
Invariant: variances add over the year¶
Squared, the volatility is the variance over the year, v n, within the root's rounding; and four times the periods doubles it.
sq = got.astype(float) ** 2
assert np.all(np.abs(sq / (F32(v).astype(float) * n) - 1) <= 8 * U)
four = annualised(v, 4.0 * n)
assert np.all(np.abs(four / (2.0 * got) - 1) <= 8 * U)
print('sigma^2 = v n and sigma(4 n) = 2 sigma(n), for all 300')
fig, ax = plt.subplots(figsize=(7, 3))
periods = np.array([1, 12, 52, 252, 365, 252 * 390], dtype=float)
ax.plot(periods, annualised(0.01 ** 2, periods) * 100, 'o-')
ax.set_xscale('log')
ax.set_xlabel('periods a year')
ax.set_ylabel('volatility a year (%) at 1% a period')
plt.show()
sigma^2 = v n and sigma(4 n) = 2 sigma(n), for all 300
Limit: nothing, and one period¶
No variance is no volatility, and one period a year is the period's own standard deviation, bit for bit.
out = annualised([0.0, 0.01 ** 2, 2.5e-5], [252.0, 1.0, 1.0])
assert out[0] == 0.0 and out[1] == F32(np.sqrt(F32(0.01 ** 2))) and out[2] == F32(np.sqrt(F32(2.5e-5)))
print(out)
[0. 0.01 0.005]