Warning: Undefined array key "DOKU_PREFS" in /home/clients/a4e6fc1ce1761b72982b805de0f418c4/web/python/mesrecettespython/inc/common.php on line 2082
loi_normale [Les recettes Python de Tyrtamos]

Les recettes Python de Tyrtamos

Outils pour utilisateurs

Outils du site

Piste:
loi_normale

Différences

Ci-dessous, les différences entre deux révisions de la page.

Lien vers cette vue comparative

Les deux révisions précédentes Révision précédente
Prochaine révision 		Révision précédente
Dernière révision Les deux révisions suivantes
loi_normale [2008/04/16 09:22]
tyrtamos 		loi_normale [2014/03/17 09:21]
tyrtamos
Ligne 1: Ligne 1:
 ====== Loi normale (ou loi de Gauss) ====== ====== Loi normale (ou loi de Gauss) ======
 +
 +[modification le 17/3/2014: fonction de répartition de la loi normale centrée réduite]
  
 ===== Rappel de quelques définitions concernant la loi normale ou loi de Gauss ===== ===== Rappel de quelques définitions concernant la loi normale ou loi de Gauss =====
Ligne 37: Ligne 39:
  
 <m>Prob(<u)~=~int{-infty}{u}{1/{sqrt{2pi}}e^{-1/2{u}^2}}du</m> <m>Prob(<u)~=~int{-infty}{u}{1/{sqrt{2pi}}e^{-1/2{u}^2}}du</m>
 +
 +=== Calcul par intégration numérique ===
  
 Pour calculer cette surface, on va faire quelque chose qui ressemble à l'intégrale: on va décomposer la surface en "tranches" verticales toute petites (2000 tranches par écart-type!). Pour chaque tranche, on peut en calculer la hauteur puisqu'on a l'équation de la courbe, et on peut donc en calculer la surface puisqu'on a la largeur. On fera la somme de tous ces petits rectangles et on aura la surface recherchée, donc la probabilité recherchée.  Pour calculer cette surface, on va faire quelque chose qui ressemble à l'intégrale: on va décomposer la surface en "tranches" verticales toute petites (2000 tranches par écart-type!). Pour chaque tranche, on peut en calculer la hauteur puisqu'on a l'équation de la courbe, et on peut donc en calculer la surface puisqu'on a la largeur. On fera la somme de tous ces petits rectangles et on aura la surface recherchée, donc la probabilité recherchée. 
Ligne 44: Ligne 48:
 Ce qui fait que nous calculerons toujours la surface pour abs(u), et nous corrigerons si le u initial était négatif. Ce qui fait que nous calculerons toujours la surface pour abs(u), et nous corrigerons si le u initial était négatif.
  
-On va considérer aussi que pour %%u>=7.56%% (7.56 x écarts-types), la probabilité est égale à 1. On peut en effet vérifier que les limites du calcul en double précision sont atteintes au delà. De même pour %%u<=-7.56%%, valeur pour laquelle la fonction renverra une probabilité nulle.+On va arrondir le résultat à 7 chiffres après la virgule: il est rare que les tables donnent plus de 5 chiffres après la virgule, et, compte tenu de la grande quantité de calculs, l'accumulation des petites erreurs de calcul liées à la précision des flottants, nous obligent à limiter la précision obtenue. Pour une plus grande précision, il faudrait faire les calculs avec le module decimal de Python.
  
 \\ \\
-<code python> +<code python>from math import *
-#!/usr/bin/python +
-# -*- coding: utf-8 -*- +
-from __future__ import division +
-  +
-from math import *+
  
 def pgaussred(x): def pgaussred(x):
-    """pgaussred(x): probabilité qu'une variable aléatoire distribuée selon une loi normale réduite soit inférieure à x"""+    """fonction de répartition de la loi normale centrée réduite 
 +       (= probabilité qu'une variable aléatoire distribuée selon 
 +       cette loi soit inférieure à x) 
 +       calcul par intégration numérique: 2000 tranches par écart-type 
 +       résultat arrondi à 7 chiffres après la virgule 
 +    """
     if x==0:     if x==0:
         return 0.5         return 0.5
-    if x>=7.56: 
-        return 1.0 
-    if x<=-7.56: 
-        return 0.0 
     u=abs(x)     u=abs(x)
     n=int(u*2000)     n=int(u*2000)
Ligne 69: Ligne 69:
     f1=k     f1=k
     p=0.5     p=0.5
-    for i in xrange(0,n):+    for i in range(0,n):
         u2=u1+du         u2=u1+du
         f2=k*exp(-0.5*u2*u2)         f2=k*exp(-0.5*u2*u2)
Ligne 75: Ligne 75:
         u1=u2         u1=u2
         f1=f2         f1=f2
-    if x>0: +    if x<0: 
-        return p +        p = 1.0-p 
-    else: +    return round(p, 7) 
-        return 1.0-p+</code>
  
-# exemples d'utilisation +Exemples d'utilisation 
-print pgaussred(-8) # affiche: 0.0  + 
-print pgaussred(-7) # affiche: 1.26010313295e-012 +<code python> 
-print pgaussred(-6) # affiche: 9.86570047878e-010 +for x in range(-6, 7): 
-print pgaussred(-5) # affiche: 2.86651707593e-007 +    P = pgaussred(x) 
-print pgaussred(-4) # affiche: 3.1671252966e-005 +    print(x, P)
-print pgaussred(-3) # affiche: 0.0013498983086 +
-print pgaussred(-2) # affiche: 0.0227501341978 +
-print pgaussred(-1) # affiche: 0.158655258973 +
-print +
-print pgaussred(0) # affiche: 0.5 +
-print +
-print pgaussred(1) # affiche: 0.841344741027 +
-print pgaussred(2) # affiche: 0.977249865802 +
-print pgaussred(3) # affiche: 0.998650101691 +
-print pgaussred(4) # affiche: 0.999968328747 +
-print pgaussred(5) # affiche: 0.999999713348 +
-print pgaussred(6) # affiche: 0.999999999013 +
-print pgaussred(7) # affiche: 0.999999999999 +
-print pgaussred(8) # affiche: 1.0+
 </code> </code>
 +
 +Ce qui donne:
 +
 +<code python>
 +-6 0.0
 +-5 3e-07
 +-4 3.17e-05
 +-3 0.0013499
 +-2 0.0227501
 +-1 0.1586553
 +0 0.5
 +1 0.8413447
 +2 0.9772499
 +3 0.9986501
 +4 0.9999683
 +5 0.9999997
 +6 1.0
 +</code>
 +
 +
 +=== Calcul par une formule simplifiée ===
 +
 +Certains livres de mathématiques propose des formules approchées qui donnent de bons résultats plus rapidement. J'ai retenu une des formules proposées par l'excellent livre "Handbook of Mathematical Functions" de Abramovitz & Stegun. Voilà une proposition de codage Python:
 +
 +<code python>from math import *
 +
 +def pgaussred(x):
 +    """fonction de répartition de la loi normale centrée réduite
 +       (= probabilité qu'une variable aléatoire distribuée selon
 +       cette loi soit inférieure à x)
 +       formule simplifiée proposée par Abramovitz & Stegun dans le livre
 +       "Handbook of Mathematical Functions" (erreur < 7.5e-8)
 +    """
 +    u = abs(x) # car la formule n'est valable que pour x>=0
 +
 +    Z = 1/(sqrt(2*pi))*exp(-u*u/2) # ordonnée de la LNCR pour l'absisse u
 +
 +    b1 = 0.319381530
 +    b2 = -0.356563782
 +    b3 = 1.781477937
 +    b4 = -1.821255978
 +    b5 = 1.330274429
 +
 +    t = 1/(1+0.2316419*u)
 +    t2 = t*t
 +    t4 = t2*t2
 +
 +    P = 1-Z*(b1*t + b2*t2 + b3*t2*t + b4*t4 + b5*t4*t)
 +
 +    if x<0:
 +        P = 1.0-P # traitement des valeurs x<0
 +
 +    return round(P, 7) # retourne une valeur arrondie à 7 chiffres
 +</code>
 +
 +Cette formule donne de bons résultats avec une erreur < 7.5e-8. Il est d'ailleurs très rare qu'une table donne plus que 5 chiffres après la virgule.
 +
 +Exemple d'utilisation:
 +
 +<code python>
 +for x in range(-6, 7):
 +    P = pgaussred(x)
 +    print(x, P)
 +</code>
 +
 +Ce qui donne:
 +
 +<code python>
 +-6 0.0
 +-5 3e-07
 +-4 3.17e-05
 +-3 0.0013499
 +-2 0.0227501
 +-1 0.1586553
 +0 0.5
 +1 0.8413447
 +2 0.9772499
 +3 0.99865
 +4 0.9999683
 +5 0.9999997
 +6 1.0
 +</code>
 +
 +Résultats quasi identiques au code précédent
  
 ==== Cas de la loi normale ==== ==== Cas de la loi normale ====
loi_normale.txt · Dernière modification: 2014/03/17 09:22 de tyrtamos

Outils de la page

Sauf mention contraire, le contenu de ce wiki est placé sous les termes de la licence suivante : CC Attribution-Noncommercial-Share Alike 4.0 International
CC Attribution-Noncommercial-Share Alike 4.0 International Donate Powered by PHP Valid HTML5 Valid CSS Driven by DokuWiki