tp3.utf8.md

Introduction

Pour commencer, lisez la partie du cours relative à ces deux théorèmes fondamentaux.

Quelques distributions usuelles:

Loi Notation Esperance Variance Support Densité/masse sur \(\texttt{R}\)
Binomiale \(\mathcal B(n, p)\) \(np\) \(np(1-p)\) \(k\in\{0, ..., n\}\) \(C_n^k p^k(1-p)^{n-k}\) \(\texttt{rbinom, dbinom, pbinom, qbinom}\)
Poisson \(\mathcal P(\lambda)\) \(\lambda\) \(\lambda\) \(k\in\mathbb N\) \(\frac{e^{-\lambda}\lambda^k}{k!}\) \(\texttt{rpois, dpois, ppois, qpois}\)
Normale \(\mathcal N(\mu, \sigma^2)\) \(\mu\) \(\sigma^2\) \(x\in\mathbb R\) \(\frac{1}{\sigma\sqrt{2\pi}}\exp\left(-\frac 1 2\left(\frac{x- \mu}{\sigma}\right)^2\right)\) \(\texttt{rnorm, dnorm, pnorm, qnorm}\)
Exponentielle \(\mathcal E(\lambda)\) \(\frac 1\lambda\) \(\frac{1}{\lambda^2}\) \(x\in[0, +\infty[\) \(\lambda e^{-\lambda x}, x\geq 0\) \(\texttt{rexp, dexp, pexp, qexp}\)

Loi forte des grands nombres et Théorème Centrale Limite

Pour chacune des lois suivantes :

  • \(\mathcal N(0,1)\): Loi normale centrée réduite.

  • \(\mathcal N(4,9)\): Loi normale de moyenne 4 et de variance 9.

  • \(\mathcal B(5,1/4)\) : Loi binomiale de taille 5 et de paramètre 1/4 (\(\color{blue}{\texttt{?rbinom}}\))..

  • \(\mathcal B(1/3)\): Loi de Bernoulli de paramètre 1/3 (\(\color{blue}{\texttt{?rbinom, size=1}}\)).

  • \(\mathcal E(2)\): Loi exponentielle de paramètre 2 (\(\color{blue}{\texttt{?rexp}}\)).

  • \(\mathcal P(2)\): Loi de Poisson de paramètre 2 (\(\color{blue}{\texttt{?rpois}}\)).

  1. Ecrire une fonction en R permettant de générer ces lois de taille \(N\) (échantillon de taille \(N\)).

  2. Vérifiez pour chaque distribution la loi forte des grands nombres. \[\bar X_n = \dfrac{1}{n} \sum_{i=1}^n X_i \xrightarrow[n \to + \infty]{} \mathbb E[X_1] \]

  3. On veut maintenant simuler \(M\) échantillons de taille \(N\). Ecrire une fonction permettant de le faire.

  4. Vérifiez pour chaque distribution le théorème centrale limite. \[\bar X_n=\frac 1 n\sum_{i=1}^nX_i \xrightarrow[n\to +\infty]{} \mathcal N(\mu,\sigma^2/n)\]

Méthodes de Monte Carlo

Ces méthodes visent à calculer une intégrale par des simulations de variables aléatoires : \[ \rm I = \int g(x)f(x)\,dx \]\(f\) est la fonction densité d’une variable aléaoire \(X\). Ainsi, on a \[ \rm I=\mathbb E[g(X)] \] D’après la loi forte des grands nombres, il suffit de simuler un échantillon i.i.d \(X_1,\ldots,X_n\), de loi à densité \(f\) et de calculer \[ \frac 1n \sum_{i=1}^n g(X_i), \] cette quantité tendant vers \(I=\mathbb E[g(X)]\) quand \(n\to +\infty\).

exemple :

On cherche à calculer \[ \rm I=\int_{0}^1 \frac{\ln(1+x)}{\sqrt{1-x^2}}\,dx \]

En R, cela donne

n<-1000000
ech<-runif(n)
mean(log(1+ech)/sqrt(1-ech^2))
## [1] 0.7424894
integrate(function(x) log(1+x)/sqrt(1-x^2),0,1)
## 0.7431381 with absolute error < 2.2e-06

On peut également utiliser une autre loi que la loi uniforme : par exemple la loi beta, notée \(\mathcal B(a,b)\). La densité de cette loi est \[ \forall x\in [0,1],\, f(x) = \frac{1}{B(a,b)} x^{a-1}(1-x)^{b-1} \]\(B(a,b)=\dfrac{\Gamma(a)\Gamma(b)}{\Gamma(a+b)}\). Pour rappel, la fonction gamma, notée \(\Gamma(x)\) vaut \[ \Gamma(x) = \int_0^{+\infty} t^{x-1}\exp(-t)\,dt \]

On peut donc écrire \(I\) sous la forme \[ \rm I = \int_0^1 \frac{g(x)}{f(x)}f(x)\,dx \]\(f\) est la fonction densité de la loi Beta.

En R, cela donne

n<-1000000
a<-5 # 1er paramètre de la loi Beta
b<-1 # 2ème paramètre de la loi Beta
ech<-rbeta(n,a,b)
g<-function(x) log(1+x)/sqrt(1-x^2)
mean(g(ech)/dbeta(ech,a,b))
## [1] 0.7403863
integrate(g,0,1)
## 0.7431381 with absolute error < 2.2e-06

TP

Approximation de \(\pi\)

Soit \(f\) la fonction définie sur \([0,1]\) par \[ \forall x\in [0,1],\, f(x) = \sqrt{1-x^2} \] Faire une procédure permettant de calculer \[ I=\int_0^1 f(x)\,dx \]

En déduire une approximation de \(\pi\).

N<- 10000  # nombre de simulations
x<-runif(N) # échantillon iid de loi Uniforme sur [0,1]
res<-mean(sqrt(1-x^2))
4*res # Approximation de Pi
## [1] 3.130869

Une autre façon d’approcher \(\pi\) est de simuler des points dans le carré unité et de compter le nombre de points qui sont à une distance de l’origine plus petite que \(1\).

N<- 1000000  # nombre de simulations
x<-runif(N) # échantillon iid de loi Uniforme sur [0,1]
y<-runif(N)
res<-mean(sqrt(x^2+y^2) <= 1)
4 * res # Approximation de Pi
## [1] 3.144816

Calculer également par Monte Carlo l’intégrale suivante : \[ \int_0^1 \dfrac{1}{1+x^2}\,dx \]

Analytiquement avec \(f(x) =\arctan x, f'(x) = \frac{1}{1+x^2}\) \[ \int_0^1 \dfrac{1}{1+x^2}\,dx = \left[\arctan(x)\right]_0^1 = \pi/4 \]

Calcul d’intégrales

Exercice 1

Nous considérons une variable aléatoire \(X\) de loi \(\Gamma(a,b)\) de densité de probabilité : \[ f(x) = \frac{b^a}{\Gamma(a)}x^{a-1}\exp(-bx) \] Dans la suite, on prend \(a= 4\) et \(b= 1\).

  1. Simuler un échantillon de 1000 réalisations de \(X\) (\(\texttt{?rgamma}\)).

  2. Donner, par des simulations, des approximations de l’espérance et de la variance de \(X\).

  3. Calculer l’intégrale suivante par Monte Carlo \[ \int_0^{+\infty} \cos(2x)x^{7/3}\exp(-3.2x)\,dx \]

Exercice 2

On cherche à calculer \[ \rm I = \int_0^{+\infty} \exp(-x^2/2)\,dx \]

  1. Montrer grâce à la loi normale et par le calcul que \(\rm I=\sqrt{\dfrac{\pi}{2}}\).

  2. En effectuant un changement de variables \(y=x^2/2\) dans l’intégrale, prouvez que \[ \rm I= \int_0^{+\infty} \frac{1}{\sqrt{2y}}\exp(-y)\,dy. \]

  1. Montrer que cette intégrale peut être vue comme l’espérance d’une fonction d’une variable aléatoire \(X\) qui suit une loi exponentielle. Donner le paramètre de cette loi exponentielle.
  1. Mettre en oeuvre un programme permettant de calculer cette intégrale par la Monte Carlo. Comparer avec la fonction integrate.

Exercice 3 : En dimension 2

On s’intéresse à l’intégrale \[ \rm I = \iint_{\mathcal D} \sqrt{r^2-x^2-y^2}\,dxdy \]\(r>0\) et \(\mathcal D=\{(x,y)\in \mathbb R^2,\, x^2+y^2\leq r^2\}\).

Mettre en oeuvre un programme permettant de calculer cette intégrale double.