On rappelle le théorème fondamental de la simulation.
Pour toute variable aléatoire d-dimensionnelle X, il existe une fonction borélienne L de ℝn dans
ℝd telle que X ait la même loi que L(U1,...,Un), où U1,...,Un sont n variables aléatoires
indépendantes et de loi uniforme sur [0, 1].
Ainsi on a par exemple :
pour la loi binômiale B(p,n),
pour la loi géométrique G(p) sur 1, 2, 3...,
pour la loi normale N(0, 1),
2 La méthode de la fonction de répartition inverse ou fonction quantile
On peut trouver L à l’aide de la fonction de répartition de la variable aléatoire réelle à simuler
lorsque celle-ci possède une inverse simple à calculer.
On note FX(x) = P(X ≤ x) la fonction de répartion de la variable aléatoire X.
La méthode des quantiles nous dit alors que si on note GX(y) = inf (x,FX(x) ≥ y),
alors GX(U) a la même loi que X, si U suit une loi uniforme sur [0, 1]. On prend donc
L(U) = GX(U).
Ainsi pour une loi exponentielle E(λ), on :
Question 1Simulation de variables aléatoires de loi binômiale et géométrique, de loiexponentielle et normale.
En utilisant les rappels ci-dessus, réaliser pour chacune des ces quatre lois :
l’écriture d’une fonction Scilab de simulation d’un tirage aléatoire ,
un algorithme de simulation de s tirages (ici encore une fonction Scilab dépendant desparamètres de la loi et du nombre s de tirages),
un histogramme pour vérifier l’adéquation à la loi originale.
Question 2Simulation de loi normale.
On donne une procédure Scilab de simulation de variables aléatoires de loi normaleN(0, 1) qu’on écrit dans un fichier s_normal.sci:
function [x]=s_normal(s,a,b)
x=[];
for i=1:s,
z=rand(1,2);
y=a+b*sqrt(-2*log(y(1)))*cos(2*%pi*y(2))
x=[x y],
Écrire une autre procédure s_normal.scine faisant appel qu’une seule fois à la fonctionrand.
3 La méthode du rejet
Si X est une variable aléatoire réelle bornée admettant la densité p sur R, si (Uk(1),Uk(2))k>0 est
une suite de variables aléatoires indépendantes de loi uniforme sur un rectangle de R2 contenant le
graphe de p, alors
définit une variable aléatoire de même loi que X.
Question 3Simulation de la loi beta.
En utilisant le rappel ci-dessus, écrire un algorithme de simulation d’une loi beta deparamètres r,s > 0, par exemple r = 2,s = 2
4 Simulation de loi discrètes
La fonction grand permet la simulation des lois discrètes usuelles. Voici pour le compléter deux
fonctions. La première engendre un échantillon d’une loi discrète, par la méthode d’inversion, la
seconde calcule les fréquences empiriques d’un tel échantillon. Ces deux fonctions peuvent être
placées à la suite l’une de l’autre dans un fichier frequences.sci.
La fonction ech_dist(x,d,m,n) retourne une matrice de taille m*n dont les coefficients sont
des réalisations indépendantes de la loi sur x spécifiée par le vecteur d, normalise a
1.
function e = ech_dist(x,d,m,n)
taille=min(max(size(x)),length(d));// ajuster les longueurs
x=x(1:taille);
d=d(1:taille);
loi = d/sum(d); // normaliser par la somme
loi=cumsum(loi); // calculer la fonction de repartition
for i=1:m,
for j=1:n,
k=1;
r=rand(1,1); // appel de random
while r>loi(k), // simulation par inversion
k=k+1,
end;
e(i,j)=x(k);
end;
end
La fonction freq_emp(ech) calcule les fréquences empiriques des valeurs différentes de ech et
retourne un vecteur des valeurs différentes de ech (v) et un vecteur des fréquences correspondantes
(f)
function [v,f] = freq_emp(ech)
taille=max(size(ech)); // taille de l'echantillon
v=[ech(1)]; // valeurs differentes
for k = 2:taille, // parcourir l'echantillon
if ech(k)~=v then, // la valeur v(k) est nouvelle
v = [v,ech(k)]; // la rajouter
end;
end;
v = mtlb_fliplr(sort(v)); // trier les valeurs trouvees
nbval=max(size(v)); // nombre de valeurs differentes
effectifs = []; // effectifs des valeurs
for k = 1:nbval, // parcourir les valeurs
e = length(find(ech==v(k)));// calculer l'effectif de la valeur k
effectifs = [effectifs,e]; // le rajouter
end;
f=effectifs/sum(effectifs); // calculer les frequences
Question 4Charger ces deux fonctions par getf et réaliser les deux exemples suivants.