TP 5

Simulation de variables aléatoires à densité

Exercice 1
1. Écrire une fonction norcentreereduite() simulant la loi N (0, 1) à l’aide de la fonction rd.random.

2. Écrire une fonction normale(m,s) simulant la loi N (m, s2 ) à partir de la fonction

norcentreereduite().
3. Écrire une fonction Normale(m,s,N) donnant un vecteur contenant N réalisations de la loi
N (m, s2 ).

1. En utilisant le premier point de la propriété précédente :

1 def norcentreereduite():
2 u = rd.random(12)
3 x = np.sum(u)-6
4 return x

2. En utilisant le deuxième point de la propriété précédente :

1 def normale(m,s):
2 x = norcentreereduite()
3 y = s*x+m
4 return y

3. On peut procéder ainsi :

1 def Normale (m,s,N):

2 X = np.zeros(N)
3 for k in range(N):
4 X[k] = normale(m,s)
5 return X

Exercice 2
1. Rappeler l’expression de la fonction de répartition de U ,→ U (]0, 1[).
2. Démontrer le théorème précédent.


0 si x ≤ 0


1. On a FU : x ∈ R 7→ x si 0 < x < 1 .

1

si x ≥ 1

2. Posons Y = F −1 (U ). On a Y (Ω) = F −1 (]0, 1[) =]a, b[= X(Ω). Et pour tout x ∈]a, b[, on
a:

FY (x) = P (Y ≤ x) = P (F −1 (U ) ≤ x) = P (U ≤ F (x)) = F (x)

F (x)∈]0,1[
ECG2 - Maths approfondies Lycée Louis Pergaud

Les fonctions de répartitions de X et Y sont égales : elles suivent donc la même loi.

Exercice 3
1. Rappeler l’expression de la fonction de répartition d’une loi exponentielle, et montrer qu’elle
réalise une bijection de R∗+ sur ]0, 1[.
Déterminer sa bijection réciproque.
2. (a) Écrire une fonction exponentielle(lambda) simulant une loi E (λ) à partir de la fonction
rd.random().
(b) Écrire une fonction Exponentielle(lambda,N) donnant un échantillon de taille N de la
loi E (λ).
3. (a) Créer un vecteur de taille 10000 contenant 10000 simulations d’une variable aléatoire suivant
la loi E (1/2).
(b) En utilisant les commandes np.mean et np.std de la librairie numpy, vérifier que la moyenne
et l’écart-type empiriques (c’est-à-dire de ce vecteur) sont bien conformes à ce qu’on attend.
4. (a) Écrire une fonction gamma(n) simulant la loi γ(n) pour tout n ∈ N∗ .
(b) Écrire une fonction Gamma(n,N) renvoyant un vecteur contenant N réalisations de la loi
γ(n).

1. Notons F la fonction de répartition de la loi E (λ). On a :

0 si x ≤ 0
(
F : x ∈ R 7→ .
1 − e−λx si x > 0

F est continue sur R∗+ , strictement croissante sur cet intervalle (car F est dérivable et
F 0 (x) = λe−λx > 0) et on a lim F = 0 et lim F = 1. Par le théorème de la bijection, F
0+ +∞
réalise une bijection de R∗+ dans ]0, 1[.

Déterminons F −1 :]0, 1[→ R∗+ . Soit pour cela x ∈ R∗+ et y ∈]0, 1[. On résout :

y = F (x) ⇔ y = 1 − e−λx ⇔ e−λx = 1 − y

1
⇔ x = − ln(1 − y)
λ
1
Ainsi on a F −1 : y ∈]0, 1[→ − ln(1 − y).
λ
2. (a) En utilisant le principe d’inversion :

1 def exponentielle(lbd):
2 u = rd.random()
3 x = -(1/lbd)*np.log(1-u)
4 return x

(b) En adaptant la fonction Normale ici :

1 def Exponentielle(lbd,N):
2 X = np.zeros(N)
3 for k in range(N):
4 X[k] = exponentielle(lbd)
5 return X

2
ECG2 - Maths approfondies Lycée Louis Pergaud

3. (a) On utilise les lignes de codes suivantes :

1 E = Exponentielle(0.5,10000)
2 print(np.mean(E),np.std(E))

(b) Rappelons que les commandes np.mean et np.std calculent respectivement la moyenne
et l’écart-type d’une série statistique. On obtient ici 2.0344615 et 2.0258905. C’est
bien conforme à ce qu’on attend, puisque l’espérance et la variance d’une loi E (λ) sont
1
égales à , c’est-à-dire 2 ici.
λ
4. (a) Rappelons pour commencer que si X1 , . . . , Xn sont indépendantes et de même loi
E (1) = γ(1), alors X1 + · · · + Xn suit une loi γ(n) par stabilité de la loi gamma. En
utilisant ce résultat, on en déduit la fonction suivante pour simuler une loi γ(n) :

1 def gamma(n):
2 X = Exponentielle(1,n)
3 y = np.sum(X)
4 return y

(b) En adaptant toujours le même programme que dans Normale ou Exponentielle :

1 def Gamma(n,N):
2 X = np.zeros(N)
3 for k in range(N):
4 X[k] = gamma(n)
5 return X

Exercice 4
On rappelle qu’une variable X suit une loi de Cauchy si elle admet pour fonction de répartition la
fonction :
1 π
 
F : x ∈ R 7→ arctan(x) + .
π 2
1 1
Une densité de X est alors la fonction f : t 7→ .
π 1 + t2
1. Montrer que F réalise une bijection de R sur ]0, 1[, et déterminer F −1 .

2. (a) Écrire une fonction cauchy() simulant une loi de Cauchy à partir de la fonction rd.random.
(b) Écrire une fonction Cauchy(N) donnant un échantillon de taille N de la loi de Cauchy.

3. (a) Créer un vecteur de taille 10000 contenant 10000 simulations d’une variable aléatoire suivant
la loi de Cauchy.
(b) Utiliser la commande np.mean avec ce vecteur pour évaluer l’espérance d’une loi de Cauchy.
Recommencer avec plusieurs échantillons. Que constatez-vous ? Une variable suivant la loi
de Cauchy admet-elle une espérance ?

1. F est continue et strictement croissante car la fonction arctangente l’est, et on a lim F = 0

−∞
et lim F = 1. Par le théorème de la bijection, F réalise une bijection de R dans ]0, 1[.
+∞

3
ECG2 - Maths approfondies Lycée Louis Pergaud

Déterminons F −1 : ]0, 1[→ R. Soit pour cela x ∈ R et y ∈]0, 1[. On résout :

1 π π
 
y = F (x) ⇔ y = arctan(x) + ⇔ πy − = arctan(x)
π 2 2
π
 
⇔ x = tan πy −
2
π
 
Ainsi on a F −1 : y ∈]0, 1[→ tan πy − .
2
2. (a) En utilisant la méthode d’inversion, on obtient :

1 def cauchy():
2 u = rd.random()
3 c = np.tan(np.pi*u-np.pi/4)
4 return c

(b) Toujours sur le même principe :

1 def Cauchy(N):
2 C = np.zeros(N)
3 for k in range(N):
4 C[k] = cauchy()
5 return C

3. (a) On utilise la commande :

1 C = Cauchy(10000)

(b) On utilise :

1 print(np.mean(C))

On obtient 0.3452672. En exécutant plusieurs fois ce code, on obtient - 1.7578416,

13.079537, . . . Les valeurs obtenues sont très dispersées, alors qu’on pourrait s’attendre
à ce qu’elles soient voisine d’une valeur qui serait l’espérance d’une loi de Cauchy. Cela
suggère qu’une variable suivant une loi de Cauchy n’admet pas d’espérance.
Vérifions que c’est effectivement le cas. Une variable
Z X suivant une loi de Cauchy
+∞
admet une espérance si et seulement si l’intégrale |t|f (t) dt converge, soit par
Z +∞ −∞
t
parité si et seulement si I = dt converge. La fonction t ∈ R+ 7→
0 π(1 + t2 )
t
étant continue sur [0, +∞[, l’intégrale I est généralisée en +∞. Et on a :
π(1 + t2 )
t 1
• ∼ .
π(1 + t ) t→+∞ πt
2

• 1t ≥ 0 pour tout t > 0.

dt
Z +∞
• diverge (intégrale de Riemann en +∞ d’exposant 1 ≤ 1).
1 t
Par théorème de comparaison, l’intégrale I diverge, et E(X) n’existe pas.

Exercice 5
1. Simuler avec la fonction Exponentielle N = 10000 valeurs de la loi E (0.5).

4
ECG2 - Maths approfondies Lycée Louis Pergaud

2. Tracer la courbe représentative de la densité f de la loi E (0.5).

3. Tracer l’histogramme des fréquences de l’échantillon obtenu (on prendra pour cela une subdivi-
sion c de l’intervalle [0, 10] en p = 100 intervalles de même longueur).
Comparer l’histogramme des fréquences de l’échantillon à la courbe représentative de f . Qu’en
pensez vous ?

4. Procéder de même pour la loi γ(3).

1. On procède ainsi :

1 # Echantillon
2 N = 10000
3 lbd = 1/2
4 x = Exponentielle(lbd,N)

2. Traçons la courbe représentative de la fonction f : t ∈ R+ 7→ (1/2) exp(−t/2).

5 # densit é
6 def f(t):
7 y = lbd*np.exp(-lbd*t)
8 return y
9 t = np.linspace(0,10,100)
10 y = f(t)
11 plt.plot(t,y)

3. On utilise la commande plt.hist pour tracer l’histogramme des fréquences de l’échantillon

x.

12 # histogramme des fré quences

13 c = np.linspace(0,10,100)
14 plt.hist(x,c,density='True',edgecolor = 'k')
15 plt.show()

On obtient :

On observe graphiquement que les aires des rectangles formant l’histogramme des fréquences
sont proches des aires délimitées par la courbe représentative de la densité. On peut donc
conclure que Exponentielle(0.5) renvoie bien des simulations suivant la loi E (0.5).

5
ECG2 - Maths approfondies Lycée Louis Pergaud

4. Procédons de même avec la loi γ(3).

1 # Echantillon
2 N = 10000
3 n = 3
4 x = Gamma(n,N)
5

6 # densit é
7 def f(t):
8 y = (t**2)*np.exp(-t)/2
9 return y
10 t = np.linspace(0,10,100)
11 y = f(t)
12 plt.plot(t,y)
13

14 # histogramme des fré quences

15 c = np.linspace(0,10,100)
16 plt.hist(x,c,density='True',edgecolor = 'k')
17

18 plt.show()

On obtient :

On observe là aussi que les aires des rectangles formant l’histogramme des fréquences sont
proches des aires délimitées par la courbe représentative de la densité : la simulation de la
loi γ(3) par la commande gamma(3) est donc pertinente.

Exercice 6
1. Expliquer les lignes de commandes proposées pour tracer la fonction de répartition empirique.

2. Simuler avec la fonction Exponentielle N = 10000 réalisations de la loi E (1).

3. Tracer la fonction de répartition empirique de l’échantillon obtenu et la fonction de répartition

théorique de cette loi. Comparer.

1. Expliquons les lignes de commandes suivantes permettant de tracer la fonction de réparti-

tion empirique associée à l’échantillon x de taille N.

1 plt.step(np.sort(x),np.arange(0,1,1/N))
2 plt.show()

6
ECG2 - Maths approfondies Lycée Louis Pergaud

Le vecteur np.sort(x) est de même longueur N que x et contient les valeurs (ou modalités)
de x triées dans l’ordre croissant. Notons que l’échantillon x provenant d’une simula-
tion d’une loi continue, aucune valeur n’apparait plusieurs fois (toutes les modalités sont
d’effectif 1).
Le vecteur np.arange(0,1,1/N) est lui aussi de longueur N, et correspond à une progression
arithmétique de 0 à 1 (exclu) de raison 1/N.
La commande plt.step(np.sort(x),np.arange(0,1,1/N)) trace alors la fonction G en
escalier d’abscisse np.sort(x) et d’ordonnée np.arange(0,1,1/N). Pour tout x ∈ R, il
existe k ∈ J0, N − 2K tel que x soit compris entre la k-ème modalité yk et la (k + 1)-ème
modalité yk+1 de np.sort(x). Et G(x) est alors égale à k/N , qui est la fréquence des
modalités inférieures ou égales à x. G est donc bien la fonction de répartition empirique
associée à l’échantillon x.

2. On le fait avec la commande suivante :

1 # Echantillon
2 N = 10000
3 lbd = 1
4 x = Exponentielle(lbd,N)

3. On utilise le code suivant :

5 # Fct de rép thé orique

6 def F(x):
7 y = 1-np.exp(-lbd*x)
8 return y
9 u = np.linspace(0,6,100)
10 v = F(u)
11 plt.plot(u,v)
12

13 # Fct de rép empirique

14 plt.step(np.sort(x),np.arange(0,1,1/N))
15

16 plt.show()

On obtient :

Ces fonctions sont très proches l’une de l’autre, notre simulation semble donc pertinente.

7
ECG2 - Maths approfondies Lycée Louis Pergaud

Exercice 7
1. Calculer Φ(0), Φ(1), Φ(1.96).

2. Soit X une variable aléatoire suivant la loi N (3, 22 ). Calculer P (X > 10), P (0 ≤ X < 3).

3. Tester les simulations obtenues des lois N (0, 1) et N (2, 32 ) en traçant les fonctions de répartition
théoriques et empiriques.

4. Comparer ces résultats avec ceux obtenus par la fonction rd.normal.

1. On utilise le code suivant :

1 import scipy.special as sp
2 p1 = sp.ndtr(0)
3 p2 = sp.ndtr(1)
4 p1 = sp.ndtr(1.96)
5 print(p1,p2,p3)

On obtient respectivement les valeurs 0.5, 0.8413447 et 0.9750021.

2. On a :
X −3 7
 
P (X > 10) = 1 − P (X ≤ 10) = 1 − P (X − 3 ≤ 7) = 1 − P ≤ = 1 − Φ(3.5)
2 2

À l’aide de la commande :

1 p = sp.ndtr(3.5)
2 disp(1-p)

on obtient que cette probabilité vaut 0.0002326290790355401 . On procède de même

pour l’autre probabilité.

3. Pour la loi N (0, 1), avec un échantillon de taille N = 10000 :

1 # Echantillon
2 N = 10000
3 m = 0
4 s = 1
5 x = Normale(m,s,N)
6

7 # Fct de rép thé orique

8 import scipy.special as sp
9 u = np.linspace(-4,4,100)
10 v = sp.ndtr(u)
11 plt.plot(u,v)
12

13 # Fct de rép empirique

14 plt.step(np.sort(x),np.arange(0,1,1/N))
15

16 plt.show()

On obtient :

8
ECG2 - Maths approfondies Lycée Louis Pergaud

Ainsi la simulation de la loi N (0, 1) par la fonction Normale semble pertinente.

De même pour la loi N (2, 9) :

1 # Echantillon
2 N = 10000
3 m = 2
4 s = 3
5 x = Normale(m,s,N)
6

7 # Fct de rép thé orique

8 import scipy.special as sp
9 u = np.linspace(-6,10,100)
10 v = sp.ndtr((u-m)/s)
11 plt.plot(u,v)
12

13 # Fct de rép empirique

14 plt.step(np.sort(x),np.arange(0,1,1/N))
15

16 plt.show()

On obtient :

Là aussi, notre simulation semble pertinente.

4. Avec la fonction rd.normal :

1 # Echantillon

9
ECG2 - Maths approfondies Lycée Louis Pergaud

2 N = 10000
3 m = 0
4 s = 1
5 x = rd.normal(m,s,N)
6

7 # Fct de rép thé orique

8 import scipy.special as sp
9 u = np.linspace(-4,4,100)
10 v = sp.ndtr(u)
11 plt.plot(u,v)
12

13 # Fct de rép empirique

14 plt.step(np.sort(x),np.arange(0,1,1/N))
15

16 plt.show()

On obtient :

Les résultats obtenus avec rd.normal sont similaires de ceux obtenus avec la fonction
Normale.

Exercice 8 (FF)
Soient a ∈ R∗+ et (X1 , . . . , Xn ) une famille de variables aléatoires indépendantes, identiquement dis-
tribuées suivant une loi uniforme sur [0, a]. On pose :

U = min(X1 , . . . , Xn ) et V = max(X1 , . . . , Xn ).

1. Déterminer les lois de U et V .

2. On considère la fonction suivante :

1 def simulation(a,n):
2 nb_sim = 10000
3 R = a*rd.random([nb_sim,n])
4 couple = np.array([np.min(R,1), np.max(R,1)])
5 return couple

Expliquer la fonction simulation.

3. Prenons n = 2 et a = 1. Comparer graphiquement la qualité de cette simulation (fonction de

répartition empirique/théorique).

10
ECG2 - Maths approfondies Lycée Louis Pergaud

1. Nous l’avons traité en TD. On obtient les fonctions de répartitions suivantes pour U et V
respectivement :
 
0 si x < 0 0 si x < 0

 

FU : x ∈ R 7→ x n et FV : x ∈ R 7→ x n
1− 1− si x ∈ [0, a] si x ∈ [0, a] .



 a  a
1

si x > a 1

si x > a

2. La commande a*rd.random([nb_sim,n]) crée une matrice à np_sim = 10000 lignes et

n colonnes dont chaque coefficient a été pris aléatoirement comme réalisation d’une loi
uniforme continue sur [0, a].
La commande np.min(R,1) renvoie le minimum de chaque ligne de la matrice R, c’est-à-
dire le minimum de n réalisations d’une loi uniforme sur [0, a]. Ainsi, np.min(R,1) est une
matrice colonne contenant 10000 réalisations de la variable U . De même, np.max(R,1)
renvoie une matrice colonne contenant 10000 réalisations de la variable V .
La fonction simulation renvoie donc un échantillon de taille 10000 pour chacune des lois
des variables U et V .

3. On peut procéder comme suit pour U (on procèderait de même pour V ).

1 # é chantillon
2 n = 2
3 a = 1
4 U = simulation(a,n)[0]
5

6 # Fonction de ré partion thé orique

7 def F(x):
8 return 1-(1-x/a)**n
9 x = np.linspace(0,1,100)
10 y = F(x)
11 plt.plot(x,y)
12

13 # Fonction de ré partition empirique

14 plt.step(np.sort(U), np.arange(0,1,1/10000))
15

16 plt.show()

On obtient la représentation graphique suivante :

11
ECG2 - Maths approfondies Lycée Louis Pergaud

Les courbes étant pratiquement confondues, cela valide cette simulation de la variable U .

Exercice 9 (FFF - Différentes simulations de la loi de Pareto)

Soient k ∈ N∗ et λ > 0. 
0 si x ≤ λ
1. Déterminer la valeur de r pour laquelle fλ : x → r est une densité de probabilité.

k+1
sinon
x
Si X admet pour densité fλ , on dit que X suit la loi de Pareto de paramètre λ et k.
2. Déterminer la fonction de répartition d’une variable suivant la loi de Pareto de paramètres λ et
k.
3. En utilisant la méthode d’inversion, simuler une variable aléatoire suivant une loi de Pareto de
paramètres λ et k.
4. Soient X1 , · · · , Xk k variables aléatoires indépendantes suivant toutes la loi uniforme sur ]0, 1].
λ
On pose alors Y = .
max(X1 , · · · , Xk )
(a) Montrer que Y suit une loi de Pareto de paramètres λ et k.
(b) En déduire une autre méthode pour simuler la loi de Pareto.
5. Comparons à présent les deux méthodes obtenues de simulation d’une loi de Pareto.
(a) En simulant N = 100000 réalisations d’une loi de Pareto à l’aide de chacune de ces méth-
odes, déterminer laquelle des deux est la plus rapide.
On pourra à cet effet utiliser la commande time.clock() de la librairie time. Pour cela,
on exécute la commande tps1 = time.clock() juste avant le début de la simulation, puis
la commande tps2 = time.clock() juste après. La différence entre ces différents temps
donnera le temps d’exécution de la portion de code encadrée (en secondes).
(b) Prenons λ = 1. Comparer graphiquement la qualité de chacune des deux simulations d’une
loi de Pareto (histogramme des fréquences/densité théorique).

1. La fonction fλ est bien continue sur R sauf éventuellement en λ, et positive. De plus :

Z +∞ Z +∞
r r
fλ (x) dx = dx = .
−∞ λ xk+1 kλk

Et donc cette intégrale vaut 1 si, et seulement si, r = kλk .

2. Pour x ≤ λ, on a P (X ≤ x) = 0, et si x > λ :

1
Z x Z x  k
λ
P (X ≤ x) = fλ (t) dt = kλ dt = 1 − .
−∞ λ tk+1 x

3. FX réalise une bijection croissante de ]λ, +∞[ sur ]0, 1[. Pour trouver sa bijection, prenons
x ∈]λ, +∞[ et y ∈]0, 1[ et résolvons :
 k
λ λ
FX (x) = y ⇔ 1 − =y⇔x= 1 .
x (1 − y) k

λ
Ainsi, FX−1 (y) = 1 .
(1 − y) k
En appliquant la méthode d’inversion, on peut simuler une loi de Pareto de paramètres λ
et k comme suit.

12
ECG2 - Maths approfondies Lycée Louis Pergaud

1 def pareto1(lbd,k):
2 y = rd.random()
3 return lbd/((1-y)**(1/k))

4. (a) On vérifie que Z = max(X1 , . . . , Xk ) est une variable à densité, avec pour fonction de
répartition : 
0 si x ≤ 0


FZ : x ∈ R 7→ x n si x ∈]0, 1[ .

1

si x ≥ 0

Cherchons la fonction de répartition de Y = Zλ . Puisque Z(Ω) =]0, 1], alors Y (Ω) =

]λ, +∞[ et FY (x) = 0 pour tout x ≤ λ. Et pour x > λ :

λ λ λ
   
P (Y ≤ x) = P ≤x =P ≤ Z = 1 − P (Z < )
Z x x
λ
= 1 − P (Z ≤ ) car Z est à densité
x
 k
λ
=1− .
x
Donc Y suit une loi de Pareto de paramètres λ et k.
(b) On peut procéder comme suit, en suivant la définition de la variable Y .

1 def pareto2(lbd,k):
2 Z = np.max(rd.random(k))
3 return lbd/Z

Exercice 10 (FFF - Méthode de simulation par rejet)

Soient X une variable à densité et f une densité de X. On suppose que f est bornée et à support
borné, c’est-à-dire qu’il existe K ∈ R+ , a, b ∈ R avec a < b tels que :

∀x ∈ R, 0 ≤ f (x) ≤ K et f (x) = 0 si x ∈
/ [a, b].

Rappelons que pour tout α < β, l’aire Aα,β comprise entre la courbe, l’axe des abscisses et les droites
d’équations x = α et x = β correspond exactement à la probabilité que la variable aléatoire X prenne
les valeurs comprises entre α et β :
Z β
Aα,β = f (t) dt = P (α ≤ X ≤ β).
α

Pour simuler la variable X, on peut tirer un

point au hasard sous la courbe représentative
de f de manière uniforme. On prend alors pour
réalisation l’abscisse de ce point.
Précisons que pour tirer un point sous la courbe,
on tire un point au hasard dans le rectangle
[a, b] × [0, K] à l’aide de deux lois uniformes. Si
ce dernier est sous la courbe, on le garde. Sinon,
on recommence.

1. Quelle est la loi du nombre d’essais avant de tirer un point qui est bien situé sous la courbe de
f ?

13
ECG2 - Maths approfondies Lycée Louis Pergaud

2. Soit X une variable aléatoire dont une densité f est définie sur R par :

6x(1 − x) si x ∈ [0, 1],

(
f (x) =
0 sinon.

(a) Tracer le graphe de f à l’aide de Python. Quelles valeurs choisir pour a, b et K ?

(b) Écrire une fonction Python d’en-tête def simulation() qui renvoie une réalisation de la
variable X.
(c) Créer un vecteur contenant N = 10000 réalisations de la loi de X. Comparer l’histogramme
des fréquences de cet échantillon à la courbe représentative de f . Conclure.

1. On répète une même épreuve à deux issues (succès ou échec) de manière identique et
indépendante, et on compte le nombre d’essais pour obtenir le premier succès. Il s’agit
donc d’une loi géométrique de paramètre p la probabilité de succès.
La probabilité de succès est égale au rapport de l’aire sous la courbe représentative de f ,
qui est égale à 1 puisque f est une densité, avec l’aire du rectangle [a, b] × [0, K], égale à
(b − a) × K. Ainsi, p = (b−a)×K
1
.

2. (a) On peut procéder comme suit.

1 def f(x):
2 return 6*x*(1-x)
3

4 x = np.linspace(0,1,100)
5 y = f(x)
6 plt.plot(x,y)
7 plt.show()

On obtient la représentation graphique suivante.

On peut prendre a = 0, b = 1, K = 3/2.

(b) On prend un point (x, y) au hasard dans le rectangle [0, 1] × [0, 3/2] à l’aide de la
commande rd.random. Puis on teste si y > f (x), c’est-à-dire si le point (x, y) est sous
la courbe de f ou non. Tant que ce point est au dessus, on recommence en prenant de
nouveau un point au hasard dans le rectangle. Et lorsque la boucle s’arrête, le point
(x, y) est en dessous de la courbe, et on renvoie x.
On peut procéder comme suit :

1 def simulation():
2 x = rd.random()

14
ECG2 - Maths approfondies Lycée Louis Pergaud

3 y = (3/2)*rd.random()
4 while y>f(x) :
5 x = rd.random()
6 y = (3/2)*rd.random()
7 return x

(c) On procède comme d’habitude.

1 N = 10000
2 x = np.zeros(N)
3 for k in range(N):
4 x[k] = simulation()
5 plt.hist(x,100,density='True',edgecolor='k')
6 plt.show()

Voici la représentation graphique obtenue :

En prenant N = 1000000, on obtient :

Cette simulation de la variable X semble donc tout à fait pertinente.

15