Travail 2 : Partie simulations.
A d´eposer sur le portail de cours le 30 octobre avant 19h00 `
Probl`eme 1 : M´ethode de rejet pour la loi normal
Nous avons ´etabli dans le cours que la m´ethode de rejet permet de simuler des donn´ees
de la loi normale N (0, 1) selon l’algorithme suivant :
Algorithme
a) G´en´erer U1, U2 et U3 uniform´ement sur [0, 1]
b) Y = − ln(U1)
c) Si U2 ≤ exp [−(Y − 1)2/2] aller `a d) et e), sinon aller `a a).
d) |Z| = Y
e) Z = Y si U3 ≤ 1/2 sinon X = −Y .
Remarquons que nous pouvons simuler des donn´ees de X ∼ N (µ, σ2
) `a partir de Z ∼ N (0, 1)
en utilisant la relation :
X = σZ + µ.
Questions :
1. Simulez `a partir de cet algorithme 10000 observations. Repr´esentez ces donn´ees par
un histogramme. Comparer ce dernier avec la densit´e th´eorique.
2. On veut utiliser cette algorithme pour estimer les moments d’une variable al´eatoire
de loi normale X ∼ N (µ = 24, σ2 = 4). On sait que les deux premiers moments sont :
E(X) = µ = 24 et E(X
2
) = µ
2 + σ
2 = 580.
Utiliser le fait que
(a) Utilisez l’algorithme pr´ec´edent pour approcher E(X) et E(X2
) par simulation
(avec k = 100 r´ep´etitions). Comparez avec les valeurs exactes. Indication : rappelons que
E(X
k
) ≈
1
N
X
N
i=1
X
k
i k = 1, 2, . . .
o`u X1, . . . , XN des valeurs simul´ees `a partir de la normale N N (µ = 24, σ2 = 4).
(b) Utilisez l’algorithme pr´ec´edent pour approcher E(X5
) par simulation (avec k =
100 r´ep´etitions).
Probl`eme 2 : M´ethode de rejet pour la loi Beta
Rappelons que la densit´e de la loi uniforme sur [0, 1] est :
g(x) =
1 si 0 ≤ x ≤ 1,
0 sinon.
1
On veut utiliser simuler des donn´ees d’une variable al´eatoire de loi beta dont la densit´e est :
f(x) =
1260x
4
(1 − x)
5
si 0 ≤ x ≤ 1,
0 sinon.
Ceci peut faire en utilisant la m´ethode de rejet bas´ee sur la densit´e g pour ´etablir un algorithme permettant de g´en´erer des donn´ees issues de la densit´e f.
1. Ecrire en d´etail un algorithme permettant de simuler des donn´ees `a partir de la va- ´
riable al´eatoire X de densit´e f(x).
2. Simulez `a partir de cet algorithme 10000 observations. Repr´esentez ces donn´ees par
un histogramme. Comparer ce dernier avec la densit´e th´eorique (Repr´esentez ces deux
quantit´es dans le mˆeme graphe).
3. En utilisant la simulation pr´ec´edente, estimez P(0 ≤ X ≤ 0.70) avec une taille
d’´echantillon N = 1000 et k = 100 r´ep´etitions. Indication utiliser le fait que :
P(0 ≤ X ≤ 0.70) = card
i ∈ 1, . . . , N} : 0 ≤ Xi ≤ 0.70
N
o`u X1, . . . , XN un ´echantillon simul´e `a partir de la loi de X.
4. En ´ecrivant P(0 ≤ X ≤ 0.45) sous forme d’int´egrale `a savoir :
P(0 ≤ X ≤ 0.70) = Z 0.70
0
f(x)dx.
Estimez P(0 ≤ X ≤ 0.70) en utilisant la m´ethode de Monte-Carlo vue dans le cours
N = 1000 et k = 100 r´ep´etitions.
5. Estimez cette probabilit´e en utilisant la m´ethode du Trap`eze rappel´ee dans le cours.
6. Finalement, calculez exactement cette probabilit´e. Comparez les r´esultats des m´ethodes pr´ec´edentes `a la valeur exacte de cette probabilit´e. Quelle est votre conclusion ?
Probl`eme 3 : Calcul d’int´egrale impropre
La m´ethode de Monte Carlo peut ˆetre utilis´ee pour approximer une int´egrale impropre
(c-`a-d, une int´egrale dont au moins une de ces bornes est infinie). Pour ce faire, consid´erons
l’int´egrale impropre suivant
I =
Z ∞
0
f(x)dx <>
L’id´ee principale consiste `a ´ecrire l’int´egrale I sous forme d’une esp´erance math´ematique
d’une variable al´eatoire X `a valeur sur l’intervalle [0,∞[. Ainsi, si g d´esigne la densit´e de X,
alors on peut d´eduire que
I =
Z ∞
0
f(x)
g(x)
g(x)dx = E
f(X)
g(X)
.
2
Par cons´equent :
I ≈
1
N
X
N
i=1
f(Xi)
g(Xi)
o`u X1, . . . , XN un ´echantillon al´eatoire g´en´er´e `a partir de la densit´e g. Pour illustrer cette
m´ethode, consid´erons la fameuse int´egrale de Gauss :
I =
Z ∞
0
e
−x
2
dx.
Gauss a montr´e que cette l’int´egrale I =
√
π/2 = 0.8862269. Notre objectif est d’´evaluer
cette int´egrale par la m´ethode de Monte Carlo. Ceci peut ˆetre r´ealis´e comme suit :
1. Montrer
I = E
e
−X(X−1)
o`u X une varaible al´eatoire de loi exponentielle de param`etre λ = 1. Rappelons que
la densit´e d’une loi exponentielle de param`etre λ = 1 est :
g(x) =
e
−x
si x ≥ 0,
0 sinon.
2. En d´eduire que
I ≈
1
N
X
N
i=1
e
− ln(Ui)(ln(Ui)+1)
o`u U1, . . . , UN un ´echantillon g´en´er´e de la loi uniforme [0, 1].
3. En prenant N = 1000 avec un nombre de r´epitition k = 100, calculer `a l’aide de la
formule pr´ec´edente une approximation de I. Comparez avec la valeur exacte de I.
4. Donner une approximation de I en utilisant la formule de trap`eze vu au cours en
prenant le nombre de subdivision n = 1000 et a = 0 et b = 2000. Comparer avec la
valeur exacte de I.
5. Comparez les approximations obtenues dans les questions 3) et 4).
6. En applicant la mˆeme d´emarche utilis´ee aux points 1), 2) et 3), calculer par la m´ethode
de Monte Carlo les int´egrales suivantes :
6.1) : J =
R ∞
0
e
−x
3
dx.
6.2) : K =
R π
0
sin(x
2
)dx.
Probl`eme 4 : Calcul d’int´egrale double
La m´ethode de Monte Carlo peut ˆetre utilis´ee pour approximer une int´egrale double :
D =
Z b
a
Z d
c
f(x, y)dxdy.
3
On peut montrer que (vous pouvez ajouter la preuve comme question facultative)
D = (b − a)(d − c)E
f