Td simulation

En général, ce générateur n’est que « pseudo-aléatoire » au sens où il est construit sur un algorithme déterministe. Mais on suppose qu’il est parfait, et que différents appels à cet algorithme fournissent des variables aléatoires de loi uniforme sur [0, 1], Indépendantes. A partir de cette brique de base, on peut obtenir toutes les lois, par des considérations assez simples. 1. Utilisation de propriétés de certaines lois (1) La loi de Bernoulli de paramètre p. C’est aussi une brique de base.

L’algorithme le plus simple consiste à couper l’intervalle [O, 1] n deux sous-intervalles, l’un de taille p, l’autre de taille 1 – p : U = rand(l); if (U < p) thenX= 1; else X = O end Ou encore sous scilab En effet on vérifie bien que X est à valeurs dans {0; 1} et que (2) Loi discrète (Par exemple) Loi sur {1, 34, 48} de loi P(l) = 32 , 16 , P(48) 16 On partitionne l'intervalle [O, 1] en trois sous-intervalles de taille respective 23 , 16 et 16 . Une manière de faire cela consiste prendre [0, 321. [23 , 23 + 61 1, [23 + 16, 1].

L’algorithme est alors U rand(l); if (U 2/3) then X = 1; elseif (U < 5/6) then X = 34; else X = 48 Ou encore sous scilab : X = (U < 2/3) + 34 ((1J > 2/3) & (IJ < 5/6)) +48 * (IJ > 5/6) (3) Loi discrète uniforme sur {1 2 On va maintenant partitionner l’intervalle [0, 1] dans lequel LI prend ses valeurs en k sousintervalles de taille kl . Il faut ensuite déterminer dans lequel de ces sous-intervalles se trouve IJ Une manière de faire cela est d’utiliser la partie entière supérieure : X = ceil(k 2 Chapitre Ill.

Simulation de variables aléatoires (4) Loi binomiale de param 2 Bernoulli jusqu’à obtenir un succès : while(lJ > 1; U = rand(l); end x (6) Loi uniforme sur un intervalle On utilise le fait qu’une transformation affine d’une loi uniforme sur un intervalle est une loi uniforme sur l’image de l’intervalle. Ainsi il suffit de trouver a et b tel que l’image par x ax+b de [O, 1] soit puis : 2. Inversion de la fonction de répartition Rappel du théorème fondamental • Proposition 1 . Soit F une fonction de répartition.

On appelle fonction quantile la fonction F – (u) inf {x u}, Alors si U est une variable de loi uniforme sur [O, 1], la variable X = F – (U ) a pour fonction de répartition Ainsi, à partir d’un générateur de loi uniforme sur [O, 1], on peut n principe simuler n’importe quelle loi, à condition d’avoir un algorithme pour calculer F – . Rappelons que lorsque X est à densité, F- est la fonction inverse (au sens fonction réciproque) de la fonction de répartition. Exemple : simulation d’une loi exponentielle de paramètre X.

La fonction de répartition est x 1 – exp(-Xx), sa réciproque est x —M ln(l —x). D’où l’algorithme (car U et 1 —U ont la même loi 3 simple consiste à donner une méthode pour simuler des points suivant une distribution uniforme sur un borélien de R2 de mesure finie. 3. Méthode du rejet 3 Théorème 2. Méthode du rejet 1 Soient B c A deux boréliens de R2 de mesure de Lebesgue firmes O s MA). soit (xn ) une suite de vecteurs aléatoires indépendants de loi uniforme sur A.

On définit pour tout entier n Yn IB (Xn Alors – La suite (Yn ) est une suite de variables aléatoires indépendantes suivant une loi de Bernoulli de paramètre p = À(B) – T = min{n 0, Xn E B} suit une loi géométrique de paramètre p et est en particulier P. S. finie – les variables T et XT sont indépendantes – XT suit une loi uniforme sur B. Ainsi à partir d’une loi uniforme sur un borélien A il est assez facile e simuler une loi uniforme sur un borélien B c A. Cela prendra par cette méthode d’autant plus de temps que A est plus gros que B.

Pour A on prend souvent un pavé [a, b] x [c, d] : en effet simuler une loi uniforme sur ce pavé se fait en simulant les deux coordonnées de manière indépendante, l’abscisse uniforme sur [a, b] et l’ordonnée uniforme sur Exemple : l’algorithme scilab X = rand(l) Y = rand(l) 4 unité, au bout d’un nombre géométrique d’itérations, de paramètre Tt/4. On peut aussi sous scilab obtenir un échantillon d’un coup : l’algorithme U = rand(1000, 1) V = rand(1000, 1) Jef ind(U 2+V2<1) simule un échantillon de points uniformes indépendants sur un quart de disque unité.

Le nombre de points de l’échantillon est aléatoire : c’est ici une binomiale de paramètres (1000, rV4). Démonstration. — Les Yn sont indépendantes et de Bernoulli puisqu’à valeurs dans {0, 1} ; leur paramètre est P (xn E B) = MB) Se déduit du point précédent. – Soit cc B un borélien de R2 et n à 1. n P(XI E /B,X2E B)n-l P(Xn e C) (1 – p)n-l MC) On obtient bien un produit d’une fonction de n par une fonction de C : XT et T sont donc bien indépendantes Comme on connait la 101 de T on en déduit immédiatement P(XT : XT suit donc bien une loi uniforme sur B.

Cette méthode s’utilise facilement comme on l’a vu pour simuler par exemple une loi uniforme sur un disque ; mais on peut a our simuler une variable S 8 densité de probabilité sur R. Si le vecteur (X Y) E R2 suit une loi uniforme sur {Cx, y) e R2 , O < y < f (x)}, alors X a pour densité Ainsi pour simuler une v. a. ayant une densité donnée, il peut suffire en vertu des propositions précédentes de trouver un pavé qui contienne {Cx, y) e R2 , O < y < f et d'appliquer la méthode du rejet (ce qui écessite d'avoir un algorithme pour calculer f On peut noter que le premier point est impossible si f a un support non borné...

Démonstration. On commence par calculer y) E R2 , O < y f (application directe de Fubini). Ensuite si x R, on obtient encore par Fubini f (x)dx 1 Comment faire lorsque la loi qu'on souhaite simuler n'a pas un support borné ? On va généraliser les principes ci-dessus. Le lemme ci-dessous est la réciproque du théorème précédent : il dit que si on sait simuler une variable aléatoire de densité g (par exemple par inversion de la fonction de épartition), alors on peut pour tout c > O facilement simuler un vecteur aléatoire de loi uniforme sur {(x,y) ER2, O < y cg(x)}.

Lemme 4. Soit g une densité de probabilité sur R. Soient X et IJ deux variables aléatoires indé endantes, telles que Xa pour densité e loi uniforme sur [O, 1]. S 1]. Alors pour tout c > O, le vecteur (X, ) suit une loi uniforme sur y) e R2 , O < y < cg(x)}. Démonstration. Soit c 0. Il est clair que (X, cg(X)U ) est à valeurs dans A = y) R2 , O < y < cg(x)}, et on vérifie facilement que MA) = c. Soit tp une fonction borélienne bornée de R2 . cg(x) cp(x, cg(x)u)dxdu - R ( 0 dxdy = A NA) Enfin on obtient une méthode assez générale de simulation par rejet .

Théorème 5. Méthode du rejet 3 Soit f une densité de probabilités sur R. On suppose qu’il existe c > 0 et une densité de probabilités g telle que pour tout réel x f (x) cg(x). Soient (Un ) une suite de variables aléatoires indépendantes de loi uniforme sur [O, 1] et (Xn ) une suite de variables aléatoires indépendantes de densité g, indépendantes de (Un ). Alors T – min{n > O, cg(Xn )lJn < f (Xn )} suit une loi géométrique de paramètre Ic et cg(XT a pour densité g. Ainsi, si on veut simuler une v. a. e densité f , il suffit de trouver un réel c > O et une densité g que l’on sait simuler telles que f (x) cg(x). Démonstration. D’après le lemme précédent le vecteur (Xn , cg(Xn )lJn ) suit une loi uniforme sur A = y) R2 , O < y < cg(x)}. Il suffit donc d'appliquer le premier théorème concernant la méthode du reiet avec B = V) E rejet avec B = Y) E , o q < f MB) = = c. 4. Box Muller C'est une méthode pour simuler un couple de variables aléatoires indépendantes de loi normale centrée réduite. Elle est basée sur le théorème suivant : 5 Théorème 6.