La Méthode de Monte Carlo : MMC

La méthode de Monte-Carlo est une méthode probabiliste permettant le calcul approche d'intégrales (simples ou multiples) de fonctions quelle que soit leur régularité. C'est cette propriété qui explique son intérêt par rapport aux méthodes déterministes classiques.

Pour simplifier cette présentation, on suppose que l'on cherche a calculer : \(p=\displaystyle\int_0^1 f(t) \, dt\) pour une fonction continue dur [0 ;1] à valeurs dans [0 ;1].

Il existe deux méthodes

Méthode dite du 'rejet'

Rappelons que p est l'aire du domaine : \(D = {(x; y) \in [0; 1]^2 \quad y \leq f(x) }\)

Une première méthode possible est de tirer aléatoirement un grand nombre de points du carre \([0; 1]^2\) et de faire le quotient entre le nombre de points situés dans le domaine \(D\) et le nombre total de points.

Calcul de \(\displaystyle\int_0^1 e^{-x^2} \, dx\) par la méthode du rejet

1
nbsimul<-1000; compteur<-0
2
for (i in 1:nbsimul)
3
{
4
x<-runif(1); y<-runif(1)
5
if (y<exp(-x^2))(compteur<-sum(compteur,1))
6
}
7
aire<-compteur/nbsimul
8
aire

Méthode d'espérance

On admet le résultat suivant :

Si \(X\) est une variable aléatoire suivant une loi uniforme sur [0,1] et si \(f\) est une fonction continue sur [0,1] alors la variable aléatoire \(Y = f(X)\) possède une espérance égale à \(p=\displaystyle\int_0^1 f(t) \, dt\)

Calcul de \(\displaystyle\int_0^1 e^{-x^2} \, dx\) par la méthode d'espérance

1
nbsimul<-1000; y<-NULL
2
for (i in 1:nbsimul){x<-runif(1);y[i]<-exp(-x^2)}
3
aire<-sum(y)/nbsimul
4
aire

Exercice

Déterminer, par la méthode de Monte Carlo, l'aire de la figure ci-dessous sous forme d'une feuille, c'est à dire, l'aire du domaine \(D = {(x; y) \in [0; 1]^2 \quad x^3 \leq y \leq \sqrt{x} }\).

On utilisera les deux fonctions : \(x \rightarrow x^3\) et \(x \rightarrow \sqrt{x}\)

Rejet

1
nbsimul<-1000; compteur<-0
2
for (i in 1:nbsimul)
3
{
4
x<-runif(1); y<-runif(1)
5
if (y<sqrt(x)&&y>x^3)(compteur<-sum(compteur,1))
6
}
7
aire<-compteur/nbsimul
8
aire
1
b_inf<-aire-1.96*sqrt(aire*(1-aire))/sqrt(nbsimul)
2
b_sup<-aire+1.96*sqrt(aire*(1-aire))/sqrt(nbsimul)
3
print ('intervalle de confiance  95%')
1
cat(paste('[',round(b_inf,3),',',round(b_sup,3),']'))

Espérance

1
nbsimul<-1000; y<-NULL
2
for (i in 1:nbsimul){x<-runif(1);y[i]<-sqrt(x)-x^3}
3
aire<-sum(y)/nbsimul
4
aire
1
b_inf<-aire-1.96*sqrt(aire*(1-aire))/sqrt(nbsimul)
2
b_sup<-aire+1.96*sqrt(aire*(1-aire))/sqrt(nbsimul)
3
print ('intervalle de confiance  95%')
1
at(paste('[',round(b_inf,3),',',round(b_sup,3),']'))

Autre version

1
par(mfrow=c(2,2))# quatre figure dans une fen^etre
2
#-------------------------------- Fonction Mont Carlo-------------------------------#
3
MonteCarlo<-function(nbsimul=100){
4
curve(x^3, 0,1,n=1000,col='blue')# bordure basse de la feuille
5
curve(sqrt(x),0,1,1000,col='blue',add=T)# bordure haute
6
abline(h=0,v=0,lty=4)# les axes du repee
7
compteur<-NULL
8
for (i in 1:nbsimul){
9
x<-runif(1)# tirage d'un nombre aleaoire entre 0 et 1
10
y<-runif(1)
11
if (y<=sqrt(x) && y>=x^3){
12
compteur<-sum(compteur,1)
13
points(x,y,col='green',pch='•')# place dans le reperele point (x,y)
14
}
15
else (points(x,y,col='red',pch='-'))
16
}
17
aire<-compteur/nbsimul# estimation frequetiste de l'aire
18
text(x=0.5,y=0.5,paste(" estimation" ,round(aire,4),"\n taille",nbsimul),col=rgb(1,0,1),font=15,cex=0.8)
19
}
20
#----------------------------------- fin de la fonction----------------------------------------#
21
#---------programme principal -------------#
22
for( n in c(100,250,500,1000))(MonteCarlo(n))