Simulation d'une intégrale de variable aléatoire

Darksasukee

Bonjour,
dans un projet que je dois faire, je bloque sur quelque chose. J'aurais donc besoin d'aide s'il vous plaît.
Soit $S_{t}$ le prix d'un actif financier tel que :
$S_{t}=S_{T_{0}} \exp \left(\sigma W_{t}-\frac{\sigma^{2}}{2} t+r t\right)$

On sait qu'on peut calculer cette intégrale :
$A_{S}\left(0, T\right)=\frac{1}{T} \int_{0}^{T} S_{u} d u$ en utilisant les sommes de Rienman. Ainsi on approche $A_{S}\left(0, T\right)$ à l'aide de 

$\mathbf{Y}_{\mathbf{T}}^{\mathbf{r}, \mathbf{n}}=\mathbf{h} \sum_{\mathbf{k}=\mathbf{0}}^{\mathbf{n}-1} \mathbf{S}_{\mathbf{t}_{\mathrm{k}}}$ qu'on appellera (2).

À l'aide d'un dvlpmt de Taylor, on peut approcher $A_{S}\left(0, T\right)$ plus précisément. On l'approche ainsi grâce à
 $\mathbf{Y}_{\mathbf{T}}^{\mathbf{p}, \mathbf{n}}=\frac{\mathbf{1}}{\mathbf{T}} \sum_{\mathbf{k}=\mathbf{0}}^{\mathbf{n}-\mathbf{1}} \mathbf{S}_{\mathbf{t}_{\mathbf{k}}}\left(\mathbf{h}+\frac{\mathbf{r h}^{2}}{\mathbf{2}}+\sigma \int_{\mathbf{t}_{\mathbf{k}}}^{\mathbf{t}_{\mathbf{k}+\mathbf{1}}}\left(\mathbf{W}_{\mathbf{u}}-\mathbf{W}_{\mathbf{t}_{\mathbf{k}}}\right) \mathbf{d u}\right)$ avec $W_{t}$ un mouvement brownien donc il suit une $N(0,t)$.

Je me suis dit que pour approcher l'intégrale $\int_{\mathbf{t}_{\mathbf{k}}}^{\mathbf{t}_{\mathbf{k}+\mathbf{1}}}\left(\mathbf{W}_{\mathbf{u}}-\mathbf{W}_{\mathbf{t}_{\mathbf{k}}}\right) \mathbf{d u}$ on dit déjà qu'elle est égale à $\int_{\mathbf{t}_{\mathbf{k}}}^{\mathbf{t}_{\mathbf{k}+\mathbf{1}}}\mathbf{W}_{\mathbf{u}} \mathbf{d u}-({\mathbf{t}_{\mathbf{k+1}}}-{\mathbf{t}_{\mathbf{k}}})\mathbf{W}_{\mathbf{t}_{\mathbf{k}}}$ et qu'on peut utiliser la méthode de (2), c'est à dire les sommes de Rienman, pour approcher l'intégrale $\int_{\mathbf{t}_{\mathbf{k}}}^{\mathbf{t}_{\mathbf{k}+\mathbf{1}}}\mathbf{W}_{\mathbf{u}} \mathbf{d u}$.

Cependant le document qui me présente la méthode ne me dit pas de faire comme cela, voici sa méthode: "Remarque 2.2. en pratique, pour simuler ce schéma, nous devons, à chaque temps, simuler $W_{t_{k+1}}$ sachant $W_{t_{k}}$ et $\left(\int_{t_{k}}^{t_{k+1}} W_{u} d u \mid W_{t_{k}}, W_{t_{k+1}}\right)$. Pour la deuxième v.a, on utilise la loi (2) et pour la première le fait que $\left(W_{t_{k+1}}-W_{t_{k}}, k=\right.$ $0, \ldots, N-1)$ sont des variables Gaussiennes i.i.d."

Je ne comprends pas du tout cette remarque 2.2. Pour parle t-elle de $W_{t_{k+1}}$ sachant $W_{t_{k}}$ et de la loi du couple $\left(\int_{t_{k}}^{t_{k+1}} W_{u} d u \mid W_{t_{k}}, W_{t_{k+1}}\right)$ ??

Voici le lien de l'article si ça peut aider, il faut télécharger le document, c'est la page 5 remarque 2.2 : https://www.researchgate.net/publication/2317162_Competitive_Monte_Carlo_methods_for_the_Pricing_of_Asian_Options

Merci d'avance pour votre d'aide.

Barjovrille

Bonjour, je ne suis pas un expert mais voici ce que j'ai compris : 
Dans la pratique il faut se poser la question qu'est ce que ton ordi sait simuler, je pars de l'hypothèse qu'il sait simuler des gaussiennes.
A partir de ça tu veux simuler $\mathbf{Y}_{\mathbf{T}}^{\mathbf{p}, \mathbf{n}}$. Il faut donc simuler les $S_{t_k}$, mais pour simuler les $S_{t_k}$ il faut simuler les $W_{t_k}$, mais ton ordi ne sait pas simuler ça directement, on y va étape par étape, on suppose que tu connais $W_{t_0}$ (une simulation) (sinon je ne sais pas comment faire )  

Tu veux simuler $W_{t_1}$ comme ça tu as une simulation de $S_{t_1}$, mais ton ordi ne sait pas simuler $W_{t_1}$ directement. 
Si tu écris $W_{t_1}=W_{t_1} - W_{t_0} + W_{t_0}$ la différence tu sais la simuler c'est une gaussienne et le $W_{t_0}$ tu le connais, donc tu simules ta gaussienne et tu ajoute le $W_{t_0}$, tu obtiens une simulation de $W_{t_1}$ ( cette façon de simuler c'est bien une simulation de $W_{t_1}$ sachant $W_{t_0}$), maintenant tu connais $W_{t_1}$, tu veux une simulation de $W_{t_2}$ afin d'avoir une simulation de $S_{t_2}$ tu refais la même chose en écrivant $W_{t_2}=W_{t_2} - W_{t_1} + W_{t_1}$, et tu remarques que la différence c'est une gaussienne ... et tu continues.

Pour l'intégrale (ce n'est pas un couple je crois que c'est une notation pour  la loi de la v.a intégrale sachant...) tu peux simuler avec Riemann mais, pour simuler avec Riemann il te faut $W_{t_k}$ et $W_{t_{k+1}}$ (toujours à cause du problème qu'est ce que ton ordi peut simuler directement), et c'est bien une approximation de l'intégrale sachant $W_{t_k}$ et $W_{t_{k+1}}$

Bon je ne suis pas sûr mais ça me parait être la voie à emprunter.

Darksasukee

Punaise mais oui c'est pas la loi couple, merci ! @Barjovrille

Je travaille sur python ainsi pour simuler $W_{t_{0}}$, j'ai juste à simuler une $N(0,{t_{0}})$, la fonction numpy.random.normal le fait très bien.

Le but de mon projet est de calculer l'espérance de $\mathbf{Y}_{\mathbf{T}}^{\mathbf{p}, \mathbf{n}}=\frac{\mathbf{1}}{\mathbf{T}} \sum_{\mathbf{k}=\mathbf{0}}^{\mathbf{n}-\mathbf{1}} \mathbf{S}_{\mathbf{t}_{\mathbf{k}}}\left(\mathbf{h}+\frac{\mathbf{r h}^{2}}{\mathbf{2}}+\sigma \int_{\mathbf{t}_{\mathbf{k}}}^{\mathbf{t}_{\mathbf{k}+\mathbf{1}}}\left(\mathbf{W}_{\mathbf{u}}-\mathbf{W}_{\mathbf{t}_{\mathbf{k}}}\right) \mathbf{d u}\right)$. Je vais donc utiliser la méthode de Monte Carlo. Peux-tu me dire si mon raisonnement te semble correcte s'il te plaît?

Supposons que $t_{0}=0$, $t_{n}=T$ et $t_{k}=\frac{kT}{n}$ ainsi on a $t_{k}-t_{k-1}=\frac{T}{n}$
Ainsi pour le code, je vais simuler par exemple 50 000 vecteurs qui seront chacun une trajectoire différente d'un mouvement brownien.  La première coordonnée sera 0 (car par convention $W_{0}=N(0,0)=0$) et la k ième coordonnée est la somme de la k-1 ième coordonnée + une $N(0,\frac{T}{n})$.
 Évidemment chaque vecteur a des coordonnées différentes, car il représente une trajectoire différente de mon mouvement brownien.

Si je prends un vecteur parmi les 50 000 vecteurs que j'ai créés, je peux simuler une réalisation de $\mathbf{Y}_{\mathbf{T}}^{\mathbf{p}, \mathbf{n}}=\frac{\mathbf{1}}{\mathbf{T}} \sum_{\mathbf{k}=\mathbf{0}}^{\mathbf{n}-\mathbf{1}} \mathbf{S}_{\mathbf{t}_{\mathbf{k}}}\left(\mathbf{h}+\frac{\mathbf{r h}^{2}}{\mathbf{2}}+\sigma \int_{\mathbf{t}_{\mathbf{k}}}^{\mathbf{t}_{\mathbf{k}+\mathbf{1}}}\left(\mathbf{W}_{\mathbf{u}}-\mathbf{W}_{\mathbf{t}_{\mathbf{k}}}\right) \mathbf{d u}\right)$ en implémentant pour chaque terme de la somme de $\mathbf{Y}_{\mathbf{T}}^{\mathbf{p}, \mathbf{n}}$ l'intégrale $\int_{\mathbf{t}_{\mathbf{k}}}^{\mathbf{t}_{\mathbf{k}+\mathbf{1}}}\mathbf{W}_{\mathbf{u}}\mathbf{d u}$.

Pour calculer cette intégrale pour tous les ${k}$, je discrétise en $M$ intervalles mon intervalle $\left[\mathbf{t}_{\mathbf{k}};\mathbf{t}_{\mathbf{k}+\mathbf{1}}\right]$. Ainsi on pose $v_{a}=\mathbf{t}_{\mathbf{k}}+(\mathbf{t}_{\mathbf{k}+\mathbf{1}}-\mathbf{t}_{\mathbf{k}})\frac{a}{M}$ pour $a$ allant de 0 à $M$

Je dois donc créer un nouveau vecteur de taille $M-1$ tel que la première coordonnée sera $W_{t_{k}}$ et la p ième coordonnée est la somme de la p-1 ième coordonnée + une $N(0,\frac{T}{n})$. ainsi je pourrais approcher pour chaque terme de la somme, l'intégrale $\int_{\mathbf{t}_{\mathbf{k}}}^{\mathbf{t}_{\mathbf{k}+\mathbf{1}}}\mathbf{W}_{\mathbf{t}_{\mathbf{k}}}\mathbf{d u}$ à l'aide de la série de Riemann 
$\frac{\mathbf{t}_{\mathbf{k+1}}-\mathbf{t}_{\mathbf{k}}}{M} \sum_{\mathbf{u}=\mathbf{0}}^{\mathbf{M}-1} \mathbf{W}_{\mathbf{v}_{\mathrm{u}}}$

Je fais tout ce raisonnement pour les 50 000 vecteurs que j'ai créés au début, je fais la moyenne des 50 000 valeurs de $\mathbf{Y}_{\mathbf{T}}^{\mathbf{p}, \mathbf{n}}$ que j'ai obtenues et ceci me donnera une approximation de l'espérance de $\mathbf{Y}_{\mathbf{T}}^{\mathbf{p}, \mathbf{n}}$ d'après la méthode de Monte Carlo.

Barjovrille

Ah tu es avec $W_0=0$ c'est un peu plus simple alors, ta pas besoin de faire de sommer avec la $k-1$ième coordonnée (pour ton vecteur de normale) car pour tout $k \in \N$  $W_{t_k}$ suis la loi $N_(0,t_k)$ que tu peux simuler directement le truc dont je te parlais c'est justement le cas où $W_0 \neq 0$.

Sinon je pense que c'est bon. On dirait qu'il y a beaucoup de calculs, (mais comme je ne pratique pas beaucoup je ne sais pas si c'est normal ou il y a un moyen de faire plus optimisé, mais ce n'est pas le but de ton projet de savoir ça je pense).

Darksasukee

Désolé pour la réponse tardive j'étais au ski.

Par contre, pour simuler mon mouvement brownien, je suis obligé de sommer les coordonnées non ? Ce que je veux dire je suis obligé de faire cela : "La première coordonnée sera 0 (car par convention $W_{0}=N(0,0)=0$) et la k ième coordonnée est la somme de la k-1 ième coordonnée + une $N(0,\frac{T}{n})$." Car $W_{t_{k+1}}$ dépend de $W_{t_{k}}$, en conséquence, je ne peux pas juste simuler une $N(0,{t_{k+1}})$ et une $N(0,{t_{k}})$.

Barjovrille

Oui, pardon je suis aller trop vite, tu dois simuler $(W_{t_1},W_{t_2}-W_{t_1},...,W_{t_n}-W_{t_{n-1}})$ (ça ce sont des normales indépendantes) puis tu sommes la première coordonné avec la deuxième coordonnée (on va appeler le résultat la nouvelle deuxième coordonnée), tu prends la  "nouvelle deuxième coordonnée" et tu la somme avec la troisième coordonnée, tu prends "la nouvelle troisième coordonnée" tu sommes avec la coordonnée suivante et .... Sur python je crois que tu peux faire ça rapidement avec la fonction "cumsum" de numpy (somme cumulée).

Et à la fin tu as bien le vecteur $(W_{t_1},W_{t_2},...,W_{t_n})$.

Darksasukee

Oui voila exactement !