De l'uniforme à n'importe quoi
Un générateur fournit des nombres répartis uniformément entre 0 et 1. Un modèle a besoin d'autres formes : une pièce qui gagne 30 % du temps, un temps d'attente entre deux transactions, un rendement quotidien normal. L'astuce est toujours la même : appliquer une transformation à des nombres uniformes. La plus générale est la transformation inverse : si F est la fonction de répartition de ce que vous voulez, alors F⁻¹(u) pour u uniforme a exactement cette loi. Des tableaux entiers traversent la transformation d'un coup.
import math
import numpy as np
# Exponential waiting times with rate 2 (a mean of 0.5): F(x) = 1 - exp(-rate * x), so x = -ln(1 - u) / rate
u = np.array([0.1, 0.5, 0.9])
print((-np.log(1 - u) / 2).round(4).tolist()) # [0.0527, 0.3466, 1.1513] small u, short wait
rng = np.random.default_rng(42)
waits = -np.log(1 - rng.random(100000)) / 2
print(round(waits.mean(), 2)) # 0.5 the mean is 1 / rateRésultats discrets et loi normale
Pour une table discrète de résultats, rng.choice(values, size=n, p=probs) fait la transformation inverse pour vous : il tire chaque valeur avec sa propre probabilité, comme une journée à -2 %, 0 % ou +3 % avec des probabilités de 20 %, 50 % et 30 %. Pour la loi normale, l'inverse de sa fonction de répartition n'a pas de formule simple ; une astuce classique est la transformation de Box-Muller : à partir de deux nombres uniformes indépendants u1 et u2, z = sqrt(-2 ln u1) * cos(2 pi u2) est un tirage normal centré réduit. En pratique, rng.standard_normal(n) le fait plus vite, mais l'écrire une fois montre qu'il n'y a pas de magie.
rng = np.random.default_rng(3)
days = rng.choice(np.array([-0.02, 0.0, 0.03]), size=100000, p=[0.2, 0.5, 0.3])
print(round(days.mean(), 3), round(0.2 * -0.02 + 0.3 * 0.03, 3)) # 0.005 0.005 the sample mean matches the table
rng = np.random.default_rng(42)
u1, u2 = 1 - rng.random(100000), rng.random(100000) # 1 - random() is in (0, 1]: never log(0)
z = np.sqrt(-2 * np.log(u1)) * np.cos(2 * math.pi * u2)
print(abs(z.mean()) < 0.02, round(z.std(), 2)) # True 1.0 mean 0 and standard deviation 1Des tirages normaux aux prix
Un tirage normal centré réduit z est un choc ; c'est le modèle qui décide comment il fait bouger un prix. Le modèle usuel d'une action, le mouvement brownien géométrique, dit que le prix après t années vaut S = S0 * exp((mu - 0,5 sigma²) t + sigma sqrt(t) z). Il est log-normal : toujours positif, avec la longue queue à droite des vrais prix. Le terme -0,5 sigma² corrige la moyenne pour que le prix espéré croisse exactement au taux mu : E[S] = S0 * exp(mu t). La médiane est inférieure à la moyenne, c'est pourquoi un résultat typique est moins bon que le résultat moyen.
s0, mu, sigma, t = 100, 0.08, 0.2, 1
z = np.array([-1.0, 0.0, 1.0]) # a one-sigma loss, no shock, a one-sigma gain
prices = s0 * np.exp((mu - 0.5 * sigma ** 2) * t + sigma * math.sqrt(t) * z)
print(prices.round(2).tolist()) # [86.94, 106.18, 129.69] all positive, the upside is larger
draws = np.random.default_rng(0).standard_normal(100000)
simulated = s0 * np.exp((mu - 0.5 * sigma ** 2) * t + sigma * math.sqrt(t) * draws)
print(round(simulated.mean(), 1), round(s0 * math.exp(mu * t), 1)) # 108.3 108.3 the mean grows at mu
print(round(float(np.median(simulated)), 1)) # 106.1 the median is lower# NumPy ships ready-made samplers too; the transformations above show what they do inside
z = np.random.default_rng(0).standard_normal(100000)
print(abs(z.mean()) < 0.02, round(z.std(), 2)) # True 1.0 standard_normal: mean 0, standard deviation 1S0 * exp(mu t) pour le prix. Une erreur de signe ou une racine carrée oubliée donne des nombres d'apparence plausible mais faux, et seul un tel contrôle le révèle.- Transformation inverse : appliquer la fonction de répartition inverse à des nombres uniformes ; exponentielle :
-ln(1 - u) / rate. rng.choice(values, size=n, p=probs)échantillonne une table discrète ; Box-Mullersqrt(-2 ln u1) * cos(2 pi u2)fabrique des tirages normaux.- Prix après
tannées :S0 * exp((mu - 0,5 sigma²) t + sigma sqrt(t) z), log-normal ; vérifiez toujours un grand échantillon contre la moyenne exacte.
À suivre : l'estimation de Monte Carlo, qui utilise des tirages aléatoires pour calculer des moyennes et des probabilités sans formule.
Lab : ce qui est noté
- Tirer dans une table de résultats avec
sample_discrete25 pts - Inverser la fonction de répartition avec
exponential_from_uniform25 pts - Fabriquer des nombres en cloche avec
box_muller25 pts - Transformer des chocs en prix futurs avec
terminal_prices25 pts