Mouvement brownien géométrique

Les prix varient en pourcentage

Le mouvement brownien arithmétique ajoute le même bruit en dollars qu'une action coûte 5 ou 500, et il peut passer sous zéro. Une action se comporte autrement : un mouvement de 1 % reste 1 % à n'importe quel prix. Le mouvement brownien géométrique place l'aléa dans la variation en pourcentage : dS = mu * S * dt + sigma * S * dW. C'est le modèle de Black–Scholes. mu est le rendement espéré annuel et sigma la volatilité, la grandeur cotée sur le marché.

import math

def gbm_step(s, mu, sigma, dt, z):
    return s * math.exp((mu - 0.5 * sigma ** 2) * dt + sigma * math.sqrt(dt) * z)

print(round(gbm_step(100.0, 0.08, 0.2, 1 / 252, 0.0), 4))     # a calm day: only the (mu - sigma^2/2) drift acts
print(round(gbm_step(100.0, 0.08, 0.2, 1 / 252, 1.0), 4))     # a one-sigma up day
print(round(gbm_step(100.0, 0.08, 0.2, 1 / 252, -50.0), 4))   # an absurd fall: still positive

La solution exacte et le moins un demi sigma carré

La résolution de l'équation donne S(t) = S(0) * exp((mu - sigma**2 / 2) * t + sigma * W(t)). Pourquoi ce curieux - sigma**2 / 2 ? Parce que l'exponentielle est convexe : d'après l'inégalité de Jensen, le bruit augmente la moyenne de exp(bruit), et soustraire sigma**2 / 2 l'annule exactement, de sorte que le prix *espéré* croît au taux mu. Les *trajectoires* typiques croissent plus lentement, à mu - sigma**2 / 2, d'où une médiane inférieure à la moyenne. Le sujet suivant, le lemme d'Itô, établit ce terme.

import numpy as np

def gbm_path(s0, mu, sigma, dt, z):
    log_moves = (mu - 0.5 * sigma ** 2) * dt + sigma * np.sqrt(dt) * z
    return s0 * np.exp(np.concatenate([[0.0], np.cumsum(log_moves)]))

rng = np.random.default_rng(1)
z = rng.standard_normal(252)
path = gbm_path(100.0, 0.08, 0.2, 1 / 252, z)
print(path[:4].round(3), round(float(path[-1]), 3))

Beaucoup de trajectoires : la moyenne et la médiane

Simulez des milliers de trajectoires d'un an et regardez où elles finissent. Le prix final *moyen* est proche de S0 * exp(mu * T), la formule du prix espéré. La *médiane* est plus basse, S0 * exp((mu - sigma**2/2) * T), et avec une forte volatilité l'écart est grand : la plupart des trajectoires finissent sous la moyenne, tirée vers le haut par quelques-unes spectaculaires. Cette asymétrie à droite, une loi log-normale, est la forme des prix d'actions sur longs horizons.

import numpy as np

rng = np.random.default_rng(5)
s0, mu, sigma, T, steps, paths = 100.0, 0.08, 0.4, 1.0, 252, 20000
dt = T / steps
z = rng.standard_normal((paths, steps))
final = s0 * np.exp(np.sum((mu - 0.5 * sigma ** 2) * dt + sigma * np.sqrt(dt) * z, axis=1))

print(round(float(final.mean()), 1), round(s0 * np.exp(mu * T), 1))                      # close to each other: the expected price
print(round(float(np.median(final)), 1), round(s0 * np.exp((mu - 0.5 * sigma ** 2) * T), 1))   # the median is lower
print(round(float(np.mean(final < s0)), 3))                                              # the share of paths that lose money

Estimer mu et sigma

À partir de prix échantillonnés tous les dt ans, les rendements logarithmiques r = diff(log(prices)) sont des lois normales indépendantes de moyenne (mu - sigma**2/2) * dt et d'écart-type sigma * sqrt(dt). Donc sigma = std(r) / sqrt(dt) et mu = mean(r) / dt + sigma**2 / 2. L'estimation de la volatilité est bonne avec quelques centaines de points ; la dérive est notoirement bruitée, car son erreur standard ne diminue qu'avec le *temps* total, pas avec le nombre d'observations.

import numpy as np

def estimate_gbm(prices, dt):
    r = np.diff(np.log(prices))
    sigma = r.std(ddof=1) / np.sqrt(dt)
    return r.mean() / dt + 0.5 * sigma ** 2, sigma

rng = np.random.default_rng(8)
dt = 1 / 252
for years in (1, 10, 100):
    n = int(years * 252)
    log_moves = (0.10 - 0.5 * 0.25 ** 2) * dt + 0.25 * np.sqrt(dt) * rng.standard_normal(n)
    prices = 100 * np.exp(np.concatenate([[0.0], np.cumsum(log_moves)]))
    mu, sigma = estimate_gbm(prices, dt)
    print(years, "years:", round(float(mu), 3), round(float(sigma), 3))           # true mu 0.10 and sigma 0.25

• MBG : S(t) = S(0) * exp((mu - sigma**2 / 2) * t + sigma * W(t)) ; les variations en pourcentage sont aléatoires et les prix restent positifs.
• Le - sigma**2 / 2 annule la convexité de exp : le prix espéré est S0 * exp(mu * t), mais la trajectoire médiane croît à mu - sigma**2 / 2.
• Estimez à partir des rendements log : sigma = std(r) / sqrt(dt), mu = mean(r) / dt + sigma**2 / 2 ; la volatilité est bien estimée, la dérive mal.

À suivre : le lemme d'Itô, la règle de dérivation composée qui explique d'où vient le - sigma**2 / 2.

Lab : ce qui est noté

• Faire un pas exact avec gbm_step25 pts
• Simuler une trajectoire de prix avec gbm_path25 pts
• Moyenner tous les futurs avec gbm_expected25 pts
• Relire les paramètres avec estimate_gbm25 pts