Table des matières

1 Position du problème 6
1.1 Rappels sur les options européennes . . . . . . . . . . . . . . . . . . . . . . . . . . . 6

1.2 Objectifs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8

2 Pricing des options européennes 9

2.1 Partie théorique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9

2.1.1 Formule de Black & Scholes . . . . . . . . . . . . . . . . . . . . . . . . . . . 9

2.1.2 Méthode de Monte Carlo sans réduction . . . . . . . . . . . . . . . . . . . . 11

2.1.3 Échantillonnage préférentiel . . . . . . . . . . . . . . . . . . . . . . . . . . . 13

2.1.4 Variable de contrôle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14

2.1.5 Variables antithétiques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14

2.2 Partie numérique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15

2.2.1 Listing des programmes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16

2.2.2 Résultats numériques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17

2
Table des gures
1 Méthode de Monte Carlo sans réduction - Call. Valeur exacte en pointillés noirs . . . . . . . 18

2 Méthode de Monte Carlo sans réduction - Put. Valeur exacte en pointillés noirs . . . . . . . 20

3 Echantillonnage préférentiel - Call. Valeur exacte en pointillés noirs . . . . . . . . . . . . . . 22

4 Variable de contrôle - Call. Valeur exacte en pointillés noirs . . . . . . . . . . . . . . . . . . 24

5 Variables antithétiques - Call. Valeur exacte en pointillés noirs . . . . . . . . . . . . . . . . . 26

3
Liste des tableaux
1 Caractéristiques d'une option européenne. On notera que C0 = C et P0 = P . . . . . . . . . 6

2 Méthode sans réduction. Prélèvement de quelques valeurs approchées du Call en

fonction du nombre de simulations . . . . . . . . . . . . . . . . . . . . . . . . . . . 17

3 Méthode Sans réduction. Prélèvement de quelques valeurs approchées du Put en

fonction du nombre de simulations . . . . . . . . . . . . . . . . . . . . . . . . . . . 19

4 Echantillonnage préférentiel. Prélèvement de quelques valeurs approchées du Call

en fonction du nombre de simulations . . . . . . . . . . . . . . . . . . . . . . . . . . 21

5 Variable de contrôle. Prélèvement de quelques valeurs approchées du Call en fonction

du nombre de simulations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23

6 Variables antithétiques. Prélèvement de quelques valeurs approchées du Call en

fonction du nombre de simulations . . . . . . . . . . . . . . . . . . . . . . . . . . . 25

7 Temps de calcul du Call pour 106 simulations . . . . . . . . . . . . . . . . . . . . . 27

8 Précisions pour 106 simulations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27

4
Introduction
Le but de ce travail pratique est d'étudier les options européennes d'un point de vue numérique

en calculant leurs primes par méthode de Monte Carlo. Nous implémenterons plusieurs variantes

de la méthode an de les comparer en terme de temps de calcul aussi bien qu'en terme de précision.

Des rappels théoriques seront d'abord exposés avant d'étayer les programmes écrits sous Scilab

et de comparer les résultats obtenus.

5
1 Position du problème

1.1 Rappels sur les options européennes

Une option européenne est un produit dérivé où l'une des parties acquiert le droit, mais non

l'obligation, d'acheter (Call) ou de vendre (Put) à l'autre partie un actif , dit sous-jacent, à un

prix et une date xés à l'avance.

En général, ce droit est acquis contre une prime payée immédiatement au vendeur de l'option

et qui peut être comparée, abusivement, à des frais d'assurance.

Le prix d'exercice (Strike) est celui auquel sera échangé l'actif sous-jacent si l'option est exercée.

De même, la maturité, ou l'échéance, est la date à laquelle le sous-jacent sera échangé en cas

d'exercice de l'option.

Dans toute la suite, on adoptera les notations suivantes pour désigner les caractéristiques d'une

option européenne :

Call Put

Prime C P
Strike K K
Maturité T T
Valeur à la date t ∈ [0, T ] Ct Pt
Prix du sous-jacent à la date t ∈ [0, T ] St St

Table 1  Caractéristiques d'une option européenne. On notera que C 0 = C et P0 = P

La rationnalité de l'investisseur lui impose de n'exercer l'option que s'il réalise un bénéce à

la date de l'exercice. Un raisonnement simple basé sur la comparaison du Strike et du prix du

sous-jacent conduit aux valeurs à échéances suivantes :

CT = (ST − K)+

6
 PT = (K − ST )+

On remarque en particulier que les valeurs à échéance d'un Call et d'un Put de mêmes carac-

téristiques sont liées par la relation de parité :

CT − PT = ST − K (1)

1
Si l'on note ρ0,T le facteur d'actualisation sans risque, c'est-à-dire la somme à placer à t=0

pour obtenir, à coup sûr, une unité monétaire en t=T alors le gain net réalisé par l'acheteur du

Call européen est :

C0 C
CT − = (ST − K)+ −
ρ0,T ρ0,T

De même, l'acheteur du Put européen réalise un gain net de :

P0 P
PT − = (K − ST )+ −
ρ0,T ρ0,T

La parité de la relation 1 peut être étendue à toute date t ∈ [0, T ]. En eet, en l'absence

d'opportunité d'arbitrage (i.e : il est impossible de gagner de l'argent à coup sûr avec une mise

nette initiale nulle), on montre que :

ρ0,T
Ct − Pt = St − K ·
ρ0,t

On se place maintenant dans le cas d'un taux exponentiel r constant sur [0, T ].

Le modèle introduit par Black & Scholes valorise le Call européen, à la date de signature du

contrat, comme suit :

1. Dans le cas particulier d'un taux exponentiel r constant sur [0, T ] : ρ0,T = e−rT

7
 h
i
C=E eβ·G − K +
(2)

Où G ∼ N (0, 1) et β un coecient de l'ordre de l'unité.

De même, pour le Put européen :

h
i
P =E K − eβ·G +
(3)

1.2 Objectifs

2
Les formules 2 et 3 ramènent le pricing des options européennes à un calcul d'espérance où la

loi de la variable aléatoire est non seulement connue mais aussi usuelle.

L'objet de ce travail pratique réside dans ce calcul qui sera conduit de diérentes manières :

 Par intégration directe, en introduisant la densité de la loi normale centrée réduite. Ceci

donnera les valeurs exactes (d'après ce modèle) du Call et du Put à la date initiale.

 Par la méthode de Monte Carlo, pour obtenir des valeurs approchées des options européennes

à comparer aux valeurs exactes. On implémentera plusieurs versions de la la méthode de

Monte Carlo :

1. D'abord sans réduction de la variance

2. Ensuite avec un échantillonnage préférentiel

3. Puis par recours à une variable de contrôle

4. Enn, en introduisant des variables antithétiques

On veillera à comparer les diérentes méthodes d'approximation, notamment en terme de précision

et de temps de calcul.

2. Bien sûr, lorsque les hypothèses du modèle qui y a conduit sont vériées...

8
2 Pricing des options européennes

2.1 Partie théorique

2.1.1 Formule de Black & Scholes

On souhaite calculer sous forme explicite, en fonction de K et de β, les valeurs à t=0 des

options européennes. Pour cela, il sut, par exemple, de valoriser le Call et d'en déduire ensuite

le prix du Put.

La formule 2 de Black & Scholes pour le Call permet d'écrire :

ˆ g2
βg

e− 2
C= e − K + √ dg
R 2π
log (K)
Pour alléger l'écriture, introduisons le réel α tel que K = eβα , autrement dit α= .
β



On a de cette façon : eβg − K +
= eβg − eβα 1[α,∞[ (g) , d'où :

ˆ g2
 1 2

2(
∞ β2 − g −2βg+β 2 ) −
e e 2
C= e 2 · √ − eβα · √  dg
α 2π 2π

Or :

ˆ 1 2 ˆ ∞ − 1 (g−β)2
e− 2 (g −2βg+β )
2
∞ β2 β2 e 2
·
e2 √ dg = e 2 √ dg
α 2π α 2π
β2
ˆ ∞ − u2
e 2
(Changement de variable ane u = g − β ) = e 2 √ du
α−β 2π
β2
= e 2 φ (α − β)

1
´ ∞ e− 2 u
2

Avec φ la fonction dénie sur R par φ : x 7−→ x √ du .

2π

9
 Il vient que :
β2
C=e 2 φ (α − β) − eβα φ (α)

Soit encore :
β2
   
log (K) log (K)
C= e2φ −β −K ·φ
β β

Ce qu'on peut encore exprimer en faisant apparaître la fonction de répartition N de la loi

normale centrée réduite. En eet, cette loi étant symétrique, on a φ (x) = N (−x) pour tout x ∈ R̄.

Par conséquent :
β2
   
log (K) log (K)
C= e2N β− −K ·N − (4)
β β
βG



Pour le Put, il sut de remarquer que C − P = E e − K = E eβG − K et d'utiliser la

formule donnant la transformée de Laplace d'une loi normale centrée réduite :

β2
βG


E e = e2

Ensuite, compte-tenu de la symétrie de la loi de G et tout calcul fait :

β2
   
log (K) log (K)
P =K ·φ − −e 2 φ β−
β β

Soit encore :

β2
   
log (K) log (K)
P =K ·N − e2N −β
β β
 
log (K) 1
Application numérique : En prenant β=1 et K = 1, N = N (0) = et :
β 2

1 1
C = e 2 N (1) −
2
1 1
P = − e 2 (1 − N (1))
2

10
 En exécutant les commandes suivantes dans la console de Scilab :

[ N1 , N2 ] = c d f n o r ( "PQ" , 1 , 0 , 1 ) ;

C = e x p ( 0 . 5 ) ∗ N1 − 0.5;

P = 0.5 − e x p ( 0 . 5 ) ∗ N2 ;

On obtient :

C = 0.8871430

P = 0.2384217

En particulier, le Call coûte plus cher que le Put. On peut le concevoir par le fait que l'acheteur

du Call anticipe une hausse du prix du sous-jacent, ce qui ramené à la valeur actuelle signie que

S0 ≥ K · ρ0,T . Associée à la parité Call-Put, cette inégalité induit C ≥ P. On peut éventuelle-

ment le justier autrement en disant que le payo du Put est toujours borné par K alors que,

théoriquement, celui du Call est illimité, il semble normal qu'il soit plus cher.

Comme annoncé plus haut, on va maintenant s'intéresser à l'approximation du Call et du Put

par plusieurs variantes de la méthode de Monte Carlo.

On adoptera dans toute la suite des valeurs unitaires pour K et β.

2.1.2 Méthode de Monte Carlo sans réduction

Le but de cette section est d'écrire une méthode de Monte Carlo pour le calcul approché de C

et de P qui sera alors comparé au calcul exact. On étudiera la précision de la méthode en fonction

du nombre de simulations.

On a d'après la relation 2 :

C = E (g(G))

11
 Où G ∼ N (0, 1) et g est la fonction dénie sur R par g(x) = (ex − 1)+ .


Puisque : g(G) = |g(G)| < eG et E eG = e0.5 < ∞ alors g(G) ∈ L1 .

Soient (Gi )i≥1 des variables aléatoires i.i.d de même loi que G. D'après ce qui précède, les g(Gi )

sont toutes intégrables ; la loi forte des grands nombres assure que :

n
1X p.s
Iˆn = g(Gi ) ∼ C
n i=1 n→∞

On peut donc approcher C par la moyenne empirique des g(Gi ).

E g (G)2 < E e2G = e2 < ∞,

De plus, le T.C.L est applicable et fournit un intervalle de

conance de niveau 95% sous la forme :

σ
IC95% = C ± 1.96 √
n

Où σ 2 = V ar (E [g (G)]). Bien sûr, dans la pratique, ce paramètre est souvent inaccessible, on

remplacera :

 C par son estimateur Iˆn

n
1 X 2
 σ2 par son estimateur non biaisé : Sn2 = g (Gi ) − Iˆn
n − 1 i=1
Ceci est justié par le résultat suivant :

√ 
n ˆ 
loi
In − C −→ G ∼ N (0, 1)
Sn n→∞

Ainsi, tout compte fait, l'intervalle de conance à 95% sera de la forme :

Sn
IC95% = Iˆn ± 1.96 √
n

On procède de la même façon pour approcher le Put en se basant sur la formule :

12
 h
i
P =E 1 − eG +

Remarquons toutefois que la variable aléatoire 1 − eG +
présente la particularité d'être bornée

par 1 et que donc sa variance l'est aussi.

2.1.3 Échantillonnage préférentiel

La variable aléatoire G est gaussienne centrée, cela implique que les valeurs prises par G sont,

dans un sens à préciser, majoritairement situées dans un intervalle pas trop grand centré en 0. Plus

précisément :

P (|G| ≤ 1.96) = 0.95

Cette remarque suggère l'utilisation d'un échantillonnage préférentiel basé sur l'approximation

eg − 1 ' g valable seulement au proche voisinage de 0.

En eet :

ˆ ∞ 2
eg − 1 ge−0.5g
C = · √ dg
0 g 2π
ˆ ∞ √2g
1 e − 1 −g
(Changement de variable 0.5g 2 = y ) = √ √ · e dg
2π 0 2g
√ !
1 e 2Y − 1
= √ E √
2π 2Y
= E (f (Y ))

Où Y est une variable aléatoire qui suit la loi exponentielle de paramètre 1, i.e : Y ∼ E (1), et
√
2y
e −1
f : y 7−→ √ (dénie sur ]0, +∞[).
2 πy
Cette nouvelle formule va nous permettre d'approcher C par la méthode de Monte Carlo.

13
2.1.4 Variable de contrôle

La relation de parité Call-Put s'écrit à la date initiale :

12
C − P = E eG − 1 = e 2 − 1

√
Par conséquent : C = P + e−1 et l'accès à une valeur approchée du Put fournit immédiatement

une valeur approchée de C. De même, un intervalle de conance sur P conduit par une simple

translation à un intervalle de conance sur C.

2.1.5 Variables antithétiques

Parce que la variable gaussienne G est centrée, −G a la même loi que G. On a alors :

n
1X p.s
Iˆn = g(Gi ) ∼ C
n i=1 n→∞

n
1X p.s
Jˆn = g(−Gi ) ∼ C
n i=1 n→∞

Or, on sait que si deux suites de variables aléatoires réelles convergent presque sûrement alors

toute combinaison linéaire de ces suites converge également presque sûrement (ce qui est généra-

lement faux pour la convergence en loi, sauf si elle se fait, comme ici, vers des constantes). Par

conséquent :

p.s
Iˆn + Jˆn ∼ 2C
n→∞

Il est donc possible d'approcher la valeur du Call par :

14
 n
1 X
C ' [g (Gi ) + g(−Gi )]
2n i=1

1
La variance de cet estimateur est V2 = [Var (g(G) + Cov (g(G), g(−G))] alors que celle de
2n
1
Iˆn est V1 = Var (g (Gi )).
n
Or :

Cov (g (G) , g (−G)) = E [g (G) g (−G)] − E [g (G)]2

= E g (G) g(G) − eG − 1 − E [g (G)]2





= −E g (G) eG − 1 + Var (g(G))




De sorte que :

1 1 h G
G
i
V2 = Var (g (G)) − E e − 1 e − 1 + < V1
n 2n

On est donc théoriquement certain de réduire la variance de l'estimation, il ne reste plus qu'à

savoir, numériquement, jusqu'à quel point.

2.2 Partie numérique

Les diérentes méthodes ci-dessus ont été implémentées sous Scilab. Dans un premier temps, on

expose les programmes, fonctions et scripts, avant de discuter et comparer les résultats numériques

obtenus.

15
2.2.1 Listing des programmes

Fonction principale : Il est apparu opportun de créer une fonction principale nommée monte-

Carlo() qui puisse être utilisée pour toutes les méthodes. Les diérences entre ces dernières sont

prises en compte dans un argument de la fonction principale.

Plus précisément, il y a un noyau commun à toutes les variantes de Monte Carlo puisque

l'espérance à estimer s'écrit toujours sous la forme I = E (f (X)), seule la fonction f varie d'une

méthode à l'autre. Il convient donc de dénir f comme une entrée de la fonction monteCarlo().

On trouvera la fonction principale en Annexe 1.

Fonctions génératrices d'échantillons :

MC sans réduction, Call : La fonction est nommée sansReducCall().Voir Annexe 2.

MC sans réduction, Put : La fonction est nommée sansReducPut(). Voir Annexe 3

Echantillonnage préférentiel : La fonction est nommée preferCall() .Voir Annexe 4

Variable de contrôle : Il s'agit plutôt d'un script présenté en annexe 5.

Variables antithétiques : La fonction est nommée varAnti(). Voir Annexe 6.

Autres programmes :
 titre() : c'est une fonction appelée dans la fonction principale et qui se charge de déterminer

les titres des graphiques en fonction de la méthode de Monte Carlo choisie. Voir Annexe 7.

 Un script qui appelle la fonction monteCarlo() pour répondre aux questions du travail pra-

tique. Voir Annexe 8.

16
2.2.2 Résultats numériques

Dans tous les cas, le nombre de simulations est xé à 106 . Les diérentes valeurs visualisées sur

les graphiques sont alors prélevées à des indices multiples de 1000.

Méthode de Monte Carlo Sans réduction, approximation du Call : On constate que dès

la 10000ème simulation, au millième près, la méthode donne une valeur approchée qui correspond

à la valeur exacte. Porter le nombre de simulations jusqu'à 106 n'apporte rien de très remarquable

pour l'estimation, il est même possible qu'elle en soit détériorée comme le montre le tableau 2. En

revanche, l'augmentation du nombre de simulations peut avoir un eetpsychologique en rassurant

sur l'intervalle de conance de l'estimation. En eet, comme le montre la gure 1, la largeur de ce

1
dernier est d'autant plus petite que la taille N de l'échantillon est élevée (elle décroît en √ ).
N

Enn, soulignons que le temps de calcul de la 106 ème estimation est de : 0.113.

Nombre de simulations Valeur approchée du call

3
10 0.9147517
5 · 103 0.9001346
105 0.8872050
5 · 105 0.8859757
106 − 1 0.8842927

Table 2  Méthode sans réduction. Prélèvement de quelques valeurs approchées du Call en fonction
du nombre de simulations

17
 Figure 1  Méthode de Monte Carlo sans réduction - Call. Valeur exacte en pointillés noirs

Méthode de Monte Carlo Sans réduction, approximation du Put : On obtient cette fois-

ci une bonne approximation au centième près dès la 1000ème simulation. L'estimation s'améliore

nettement quand on augmente le nombre de simulations (jusqu'à être satisfaisante à 10−4 près).

Ceci est dû au fait que l'intervalle de conance sur le Put se resserre beaucoup plutôt que pour le

Call et d'une façon plus importante. Puisque le nombre de simulations est identique dans les deux

cas, ceci ne doit être attribué qu'à la variance du Put qui est plus petite que celle du Call.

Enn, soulignons que le temps de calcul de la 106 ème estimation est de : 0.114.

18
 Nombre de simulations Valeur approchée du call

103 0.2304315
5 · 103 0.2455560
105 0.2379595
5 · 105 0.2385309
106 − 1 0.2384556

Table 3  Méthode Sans réduction. Prélèvement de quelques valeurs approchées du Put en fonction
du nombre de simulations

19
 Figure 2  Méthode de Monte Carlo sans réduction - Put. Valeur exacte en pointillés noirs

Méthode de Monte Carlo avec échantillonnage préférentiel, approximation du Call :

Contrairement à la méthode de l'exercice 2, on gagne ici en précision dès la 1000ème simulation.

Celle-ci donne une estimation d'environ 0.88 quand l'estimation correspondante de l'exercice 2 est

d'environ 0.91. L'intervalle de conance se resserre plus vite que pour la méthode sans réduction,

20
ce qui conrme ce qui a été prévu dans la partie théorique (variance inférieure).

Le temps de calcul de la 106 ème estimation est de : 0.134, ce qui est dans le même ordre de

grandeur que la méthode de l'exercice 2.

Nombre de simulations Valeur approchée du call

3
10 0.8833481
5 · 103 0.8921928
105 0.8889665
5 · 105 0.8879513
106 − 1 0.8881118

Table 4  Echantillonnage préférentiel. Prélèvement de quelques valeurs approchées du Call en

fonction du nombre de simulations

21
 Figure 3  Echantillonnage préférentiel - Call. Valeur exacte en pointillés noirs

Variable de contrôle, approximation du Call : Cette méthode fait gagner en précision, ce

gain est en réalité hérité de l'estimation du Put qui a une variance inférieure à celle du Call.

Le temps de calcul de la 106 ème estimation est de : 0.113, qui est le même que dans l'exercice

2 puisque le Put est estimé de la même manière.

22
 Nombre de simulations Valeur approchée du call

103 0.8832064
5 · 103 0.8771537
105 0.8872630
5 · 105 0.8873483
106 − 1 0.8872287

Table 5  Variable de contrôle. Prélèvement de quelques valeurs approchées du Call en fonction

du nombre de simulations

23
 Figure 4  Variable de contrôle - Call. Valeur exacte en pointillés noirs

Variable de contrôle, approximation du Call : Aux plus bas nombres de simulations, la pré-

cision de cette méthode est plus grande que celle de l'exercice 2. Elles restent toutefois comparables

lorsque la taille d'échantillon dépasse 105 .

Le temps de calcul de la 106 ème estimation est de : 0.313, qui est plus que le double du temps

24
de calcul de l'exercice 2. Ce résultat était attendu car le nombre d'opérations arithmétiques est

plus grand pour cette méthode.

Nombre de simulations Valeur approchée du call

3
10 0.8699393
5 · 103 0.8555398
105 0.8866082
5 · 105 0.8892706
106 − 1 0.8884238

Table 6  Variables antithétiques. Prélèvement de quelques valeurs approchées du Call en fonction

du nombre de simulations

25
 Figure 5  Variables antithétiques - Call. Valeur exacte en pointillés noirs

Comparaison des diérentes méthodes Récapitulons les temps de calcul et les précisions

des diérentes méthodes pour 106 simulations.

26
 Méthode Temps de calcul

Sans réduction 0.114

Echantillonnage préférentiel 0.134
Variable de contrôle 0.113
Variables antithétiques 0.313

Table 7  Temps de calcul du Call pour 106 simulations

σ
Méthode Précision sous forme √
N
Sans réduction 0.0039859
Echantillonnage préférentiel 0.0008883
Variable de contrôle 0.0005822
Variables antithétiques 0.0025510

Table 8  Précisions pour 106 simulations

Il apparaît alors évident que la méthode la plus précise est celle de la variable de contrôle tout

en ayant le temps de calcul minimal. Ceci s'hérite de la variance du Put qui est inférieure à celle

du Call.

L'échantillonnage préférentiel fournit une précision comparable à celle de la variable de contrôle

mais avec un temps de calcul supérieur.

La méthode des variables antithétiques augmente considérablement le temps de calcul sans

apporter une précision notablement inférieure à la méthode sans réduction.

27
Conclusion
Ce travail pratique a été l'occasion d'implémenter et de comparer diérentes méthodes de Monte

Carlo appliquées au calcul des primes des options européennes. Ce qui complète les connaissances

théoriques acquises en cours de mathématiques nancières et ne pourra qu'être très utile dans le

monde professionnel où les méthodes de Monte Carlo sont très utilisées.

28