15 minutes de Python

Méthode de Monte-Carlo : calculer des valeurs approchées à l’aide de simulation

La méthode de Monte-Carlo repose sur une idée centrale des probabilités : une moyenne empirique obtenue à partir de nombreuses simulations se rapproche d’une espérance. Elle permet ainsi d’approximer numériquement des probabilités, des intégrales ou des sommes de séries lorsqu’un calcul exact est difficile ou impossible.

Fondement mathématique : loi faible des grands nombres et théorème de transfert

• Soit \((X_n)_{n \in \mathbb{N}^*}\) une suite de variables aléatoires indépendantes, de même loi qu’une variable aléatoire \(X\), admettant une espérance et une variance. On pose, pour tout \(n \in \mathbb{N}^*\) :

  \(\displaystyle \overline{X}_n = \frac{1}{n}\sum_{k=1}^n X_k\)

  La loi faible des grands nombres permet d’affirmer que la suite \((\overline{X}_n)\) converge en probabilité vers \(\mathbb{E}(X)\), c’est-à-dire :

  \(\displaystyle \forall \varepsilon > 0,\ \lim_{n \to +\infty}\mathbb{P}\!\left(|\overline{X}_n - \mathbb{E}(X)| \ge \varepsilon\right)=0\)

  Autrement dit, lorsque \(n\) est grand, il est fortement probable que la moyenne empirique \(\overline{X}_n\) soit proche de l’espérance \(\mathbb{E}(X)\). Cela fournit un moyen de proposer une valeur approchée de \(\mathbb{E}(X)\) à l’aide de simulations.

• Le théorème de transfert permet d’utiliser ce résultat pour calculer une valeur approchée d’une somme ou d’une intégrale, en l’écrivant comme l’espérance d’une variable aléatoire de la forme \( \varphi(X) \) où \( X\) est de loi connue.
  • Cas discret. Si \(X\) est une variable aléatoire discrète et si \(\varphi\) est une fonction définie sur \( X(\Omega)\), alors

    \(\displaystyle \mathbb{E}(\varphi(X)) = \sum_{x \in X(\Omega)} \varphi(x)\,\mathbb{P}(X=x)\)

    l’existence de l’espérance étant équivalente à l’absolue convergence de la série.

    Ainsi, une somme peut être interprétée comme une espérance et approchée par une moyenne empirique.

  • Cas à densité. Si \(X\) admet une densité \(f\) nulle en dehors d’un intervalle \( \left]a,b\right[ \) et si \(\varphi\) est une fonction continue sur \( ]a,b [\), éventuellement privé d’un nombre fini de points, alors :

    \(\displaystyle \mathbb{E}(\varphi(X)) = \int_{a}^{b} \varphi(x)\,f(x)\,\mathrm{d}x,\)

    l’existence de l’espérance étant équivalente à l’absolue convergence de l’intégrale.

    Ainsi, une intégrale peut être interprétée comme une espérance et approchée par simulation.

Principe de la méthode de Monte-Carlo

Quand on cherche à calculer une valeur approchée d’une somme ou d’une intégrale :

• on commence par l’écrire comme l’espérance d’une variable aléatoire de la forme \( \varphi(X) \) où \( X\) est une variable aléatoire de loi usuelle (uniforme, exponentielle, normale par exemple).,
• on simule alors un grand nombre de variables aléatoires \( X_1,\dots,X_n \) indépendantes et toutes de même loi que \( X\),
• on calcule la moyenne des valeurs prises par \( \varphi(X_1),\dots,\varphi(X_n) \) : \[ \frac{1}{n} \sum_{k=1}^n \varphi(X_k)\]

Lorsque le nombre de simulations est grand, la valeur obtenue fournit une bonne approximation de la quantité cherchée.

À retenir

• le résultat varie d’une exécution à l’autre, mais se stabilise lorsque le nombre de simulations augmente ;
• la méthode repose toujours sur une interprétation en termes d’espérance;
• la première étape est toujours : identifier clairement la variable aléatoire simulée.

import numpy as np
import numpy.random as rd

N = 10**5

U = rd.random(N)
Y = np.exp(-U**2)

I = np.mean(Y) 
print(I)