I. IntroductionEn statistique, le test du khi carré, aussi dit du khi-deux, d’après sa désignation symbolique χ2, est un test statistique où la statistique de test suit une loi du χ2 sous l'hypothèse nulle.
Par exemple, il permet de tester l'adéquation d'une série de données à une famille de lois de probabilité ou de tester l'indépendance entre deux variables aléatoires.
On souhaite dans notre cas implémenter le test du khi carré à l'aide de la librairie Python scipy.stats, pour nous permettre ensuite de vérifier qu'une série statistique suit une certaine loi de probabilité (Loi uniforme, loi de poisson, etc.).
Note : les exemples proposés dans ce billet ainsi que les définitions sont tous issus de la page wikipedia Test du χ2
II. Principe du test d'adéquation
Le test du χ2 d'adéquation ou de conformité permet de vérifier si un échantillon d'une variable aléatoire Y donne des observations comparables à celles d'une loi de probabilité P définie a priori dont on pense, pour des raisons théoriques ou pratiques, qu'elle devrait être la loi de Y.
L’hypothèse nulle (notée H0) est donc la suivante : « la variable aléatoire Y suit la loi de probabilité P ».
En termes de p-valeur, l'hypothèse nulle (l'observation est suffisamment proche de la théorie) est généralement rejetée lorsque p ≤ 0.05.
Il s'agit donc de vérifier qu'une série de données statistiques suit bien une loi de probabilité définie a priori, comme une loi de Poisson pour une variable discrète, ou une loi normale pour une variable continue.
II-A. Adéquation à une loi uniforme
On va maintenant tester l'hypothèse selon laquelle un dé à six faces n'est pas truqué, avec un risque α = 0.05.
L'hypothèse nulle H0 que l'on souhaite rejeter est donc ici : « Le dé est équilibré ».
Si le dé est lancé 600 fois de suite et s'il est équilibré, on s'attend donc que sur ces 600 jets, chaque chiffre tombe près de 100 fois.
Autrement dit, si le dé est équilibré, alors la variable k, représentant le numéro de la face supérieure, suit une loi uniforme. Elle peut ainsi prendre les 6 valeurs (1, 2, 3,.., 6) avec la même probabilité 1/6 pour chaque valeur.
Vérifier l'hypothèse « Le dé est équilibré », revient donc à vérifier que la variable k suit bien une loi uniforme.
Supposons que notre expérience donne les résultats suivants :
En considérant l'hypothèse nulle vraie, la valeur de la variable T du khi carré est donnée par la formule :
T = (88 − 100)2/100 + (109 − 100)2/100 + (107 − 1002/100 + (94 − 100)2/100 + (105 − 100)2/100 + (97 − 100)2/100 = 3.44
Le nombre de degrés de liberté est ici de 6-1 = 5. En effet, 88 + 109 + 107 + 94 + 105 + 97 = 600 et si l'on connaît par exemple les nombres de fois où l'on obtient les chiffres 1 à 5, on connaît le nombre de fois où l'on obtient le chiffre 6 : 600 - (88 + 109 + 107 + 94 + 105) = 97.
Ainsi, la statistique T suit la loi du χ2 à cinq degrés de liberté.
Cette loi du χ2 donne la valeur en deçà de laquelle on considère le tirage comme conforme avec un risque α = 0.05 :
P(T < 11.07) = 0.95
Puisque 3.50 < 11.07, on ne peut pas rejeter l'hypothèse nulle : ces données statistiques ne permettent pas de considérer que le dé est truqué, et on ne peut pas non plus rejeter l'hypothèse que la variable k suive une loi uniforme.
II-B. Conformité à une loi de Poisson
On considère une variable aléatoire Y prenant des valeurs entières positives ou nulles. Un échantillonnage de 100 valeurs de cette variable se répartit comme suit :
On souhaite tester l'hypothèse selon laquelle Y suit une loi de Poisson, avec un risque α = 0.05.
La valeur du paramètre de cette loi de Poisson est obtenue en calculant la moyenne pondérée de la série statistique, ce qui donne ici λ = 1,02.
Comme il s'agit d'une estimation on diminuera le nombre de degrés de liberté d'une unité.
Les effectifs attendus pour une loi de Poisson de paramètre λ sont :
On regroupe les effectifs supérieurs ou égaux à 3 dans une même classe, ceux supérieurs à 4 étant trop petits. La variable T prend alors la valeur 2.97. Or, la loi du χ2 à deux degrés de liberté donne :
P(T < 5.99) = 0.95
Donc, on ne rejette pas l'hypothèse que la variable aléatoire Y suive une loi de Poisson, au risque d'erreur de 5 %.
Si vous souhaitez avoir plus de précisions sur ce test vous pouvez consulter cette page wikipedia.
III. Module scipy.stats
Cette librairie Python contient un grand nombre de distributions de probabilités (distribution uniforme, de Poisson, etc..) ainsi que des fonctions permettant d'effectuer des tests statistiques comme le test du khi-carré.
III-A. Loi uniforme
Un objet uniform est utilisé pour une variable aléatoire continue suivant une loi uniforme :
Dans sa forme standard, la distribution est uniforme sur [0, 1]. En utilisant les paramètres loc=a et scale=b-a, on obtient une distribution uniforme sur [loc, loc + scale].
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 | from scipy import stats # fonction de répartition de la loi uniforme continue sur l'intervalle [0, 1] p1 = stats.uniform.cdf(0.5) # Fonction de répartition cdf - calcul de P(X ≤ 0.5) = 0.5 # fonction de répartition de la loi uniforme sur l'intervalle [1, 10] p2 = stats.uniform.cdf(x=5, loc=1,scale=10-1) # Fonction de répartition cdf - calcul de P(X ≤ 5) |
III-B. Loi de Poisson
Un objet poisson est utilisé pour une variable aléatoire discrète suivant une loi de Poisson :
La fonction de masse en k de cette loi est appelée ainsi :
p1 = stats.poisson.pmf(k, mu).
où k désigne la variable aléatoire discrète et le paramètre mu la moyenne λ de la série statistique.
La fonction de répartition en k :
p2 = stats.poisson.cdf(k, mu).
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 | from scipy import stats # fonction de masse en k=2 de la loi de Poisson de paramètre 1.02 p1 = stats.poisson.pmf(2, 1.02) # fonction de répartition en k=2 de la loi de Poisson de paramètre 1.02 p2 = stats.poisson.cdf(2, 1.02) |
III-C. Loi normale
Un objet norm est utilisé pour une variable aléatoire continue suivant une loi normale.
Le paramètre loc désigne la moyenne et l'argument scale l'écart type.
La fonction de répartition de la loi normale centrée réduite en x est appelée ainsi :
p1 = stats.norm.cdf(x).
La fonction de répartition de la loi normale pour une variable aléatoire de moyenne m et d'écart type s :
p2 = stats.norm.cdf(x, loc=m, scale=s).
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 | from scipy import stats # fonction de répartition en x=0.2 de la loi normale centrée et réduite N(0,1) p1 = stats.norm.cdf(0.2) # calcul de la probabilité P(X ≤ 0.2) # fonction de répartition en x=3.5 de la loi normale N(2.5, 4) p2 = stats.norm.cdf(3.5, loc=2.5, scale=4) # calcul de la probabilité P(X ≤ 3.5) |
III-D. Test du khi carré
La fonction scipy.stats.chisquare(f_obs, f_esp, ddof=degres_liberte) est utilisée pour vérifier la conformité d'une série statistique à une loi de Probabilité.
Elle a comme arguments :
- f_obs : tableau des effectifs observés
- f_exp : tableau des effectifs attendus (optionnel). Ils sont évalués à l'aide d'une fonction de densité ou de répartition suivant une probabilité.
- ddof : valeur entière représentant le nombre de degrés de liberté.
Valeurs retournées (T, pvalue) :
- T : valeur de type float résultat du test du khi carré.
- pvalue : valeur numérique comprise entre 0 et 1 indiquant si on peut ou pas accepter l'hypothèse nulle (elle est généralement rejetée lorsque p ≤ 0,05).
| Code Python : | Sélectionner tout |
1 2 | # test si les effectifs observées sont conformes aux effectifs estimées à l'aide d'une fonction de probabilité resultat_test = stats.chisquare(effectifs_obs, effectifs_est, ddof=degres_liberte) |
IV. Implémentation du test du khi carré
Comme on l'a vu précédemment, le test du khi carré s'effectue sur une série de valeurs.
Nous avons donc besoin de mémoriser dans des listes les valeurs et les effectifs de la série statistique, mais aussi de définir les fonctions Python pour chacune des lois de probabilités ainsi que pour le test du chi carré.
On va donc naturellement créer une classe Python intégrant ces attributs et ces méthodes. :
Notre classe comportera un constructeur, c'est-à-dire une méthode particulière __init__() dont le code est exécuté quand la classe est instanciée.
Elle va nous permettre d'initialiser les listes des valeurs et des effectifs avec les arguments de cette méthode au moment de la création de l'objet :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 | class SerieStat(): def __init__(self, valeurs, effectifs): # méthode constructeur de la classe self.valeurs = valeurs # valeurs de la série statistique self.effectifs = effectifs # effectifs des valeurs self.test_actif = None # pas de test actif au moment de la création de l'objet #... |
Pour tester cette méthode, nous ajoutons ces lignes au module :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 | # série de valeurs avec leurs effectifs valeurs = [1, 2, 3, 4, 5, 6] effectifs = [88, 109, 107, 94, 105, 97] # création d'un objet SerieStat avec la série de valeurs et les effectifs serie_stat = SerieStat(valeurs, effectifs) print("valeurs = " + str(serie_stat.valeurs)) print("effectifs = " + str(serie_stat.effectifs)) |
Le code affiche :
valeurs = [1, 2, 3, 4, 5, 6]
effectifs = [88, 109, 107, 94, 105, 97]
IV-A. Fonctions de probabilité
Ensuite, il nous faut ajouter à notre classe les différentes fonctions de probabilité qui doivent toutes avoir la même structure :
| Code Python : | Sélectionner tout |
1 2 3 4 |
def fonction_proba(self, x, cumulative = False):
...
return .. |
Arguments des fonctions de probabilité :
- x ou k: valeur pour laquelle on souhaite calculer la probabilité ;
- cumulative: indique si la fonction est cumulative ou pas, valeur False par défaut.
Les paramètres comme la moyenne, l'écart type, la valeur mini ou la valeur maxi, seront tous obtenus à partir de la série de valeurs.
IV-A-1. Loi uniforme
Nous ajoutons maintenant des méthodes à notre classe afin de déterminer les effectifs des valeurs selon une loi uniforme discrète et selon une loi uniforme continue.
Ces fonctions vont donc nous permettre de tester l'adéquation de la série de valeurs à ces lois :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 | class SerieStat(): ... def uniforme_discrete(self, k, cumulative = False): # loi uniforme discrète pour k ϵ {1, 2, ..., n} : fonction de masse ou de répartition en k # nombre n de valeurs différentes n = len(self.valeurs) if not cumulative: # si c'est une fonction de masse. # retourne la probabilité P(X = k) return 1/n else: # sinon c'est une fonction de répartition # retourne la probabilité P(X ≤ k) return k/n def uniforme_continue(self, x, cumulative = False): # loi uniforme sur [a, b] : fonction cumulative ou non cumulative # intervalle entre 2 valeurs consécutives : u = (self.valeurs[-1]-self.valeurs[0])/(len(self.valeurs)-1) # intervalle total des valeurs [a, b] : [x0-u/2, xn + u/2] a = self.valeurs[0] - u/2; b = self.valeurs[-1] + u/2 if not cumulative: # si la fonction n'est pas cumulative. # valeurs de x1 et x2, avec [x1, x2] centré sur x x1 = x - u/2; x2 = x + u/2 # renvoie la probabilité P(x1 ≤ X ≤ x2) donnée par la valeur de la fonction de répartition entre x1 et x2 return stats.uniform.cdf(x2,loc=a,scale=b-a) - stats.uniform.cdf(x1,loc=a,scale=b-a) else: # sinon # renvoie la probabilité P(a ≤ X ≤ x+u/2) donnée par la valeur de la fonction de répartition entre a et x+u/2 return (stats.uniform.cdf(x+u/2,loc=a,scale=b-a) - stats.uniform.cdf(a,loc=a,scale=b-a)) |
La méthode uniforme_continue utilise la fonction de répartition stats.uniform.cdf(x,loc=a,scale=b-a) du module scipy.stats.
Pour tester une des méthodes, nous ajoutons simplement ces lignes de code :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 10 11 12 | # série de valeurs avec les effectifs valeurs = [1, 2, 3, 4, 5, 6] effectifs = [88, 109, 107, 94, 105, 97] # création d'un objet SerieStat avec la série de valeurs et les effectifs serie_stat = SerieStat(valeurs, effectifs) # probabilité d'avoir un 6 : P(k=3) = 1/6 = 0.16666666666666.. p = serie_stat.uniforme_discrete(3) # affiche le résultat print("P(k=3) = " + str(p)) |
Le code affiche :
P(k=3) = 0.16666666666666666
IV-A-2. Loi de Poisson
Nous ajoutons maintenant une autre méthode pour obtenir l'effectif de chacune des valeurs selon une loi de Poisson et ainsi pouvoir tester la conformité de la série de valeurs à cette loi :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 | class SerieStat(): ... def poisson(self, k, cumulative=False): # loi de Poisson de paramètre m : fonction de masse ou de répartition en k # moyenne pondérée de la série statistique m = self.moyenne_ponderee() if not cumulative: # si pas cumulative # retourne la valeur en k de la fonction de masse de paramètre m return stats.poisson.pmf(k,m) else: # sinon # retourne la valeur en k de la fonction de répartition de paramètre m return stats.poisson.cdf(k,m) |
Cette méthode utilise la fonction de masse stats.poisson.pmf(x,m) et de répartition stats.poisson.cdf(x,m) du module scipy.stats. La moyenne pondérée de la série statistique est évaluée à l'aide de la méthode moyenne_ponderee() présente dans cette même classe.
Pour effectuer les tests nous ajoutons simplement ces lignes de code :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 10 11 | # série de valeurs avec leurs effectifs valeurs = [0, 1, 2, 3, 4] effectifs = [31, 45, 16, 7, 1] # création d'un objet SerieStat avec la série de valeurs et les effectifs serie_stat = SerieStat(valeurs, effectifs) # calcul de la probabilité d'avoir le 2 : P(k=2) = 0.18758148787803533 p = serie_stat.poisson(2) # affiche le résultat print("P(k=2) = " + str(p)) |
Le code affiche :
P(k=2) = 0.18758148787803533
La classe contient également une méthode donnant la valeur de la fonction de répartition de la loi normale.
IV-B. Test du khi carré
Nous définissons enfin une méthode test_khi2 permettant de vérifier qu'une série de valeurs suit bien une loi de probabilité.
Cette fonction a en particulier un paramètre fonction_proba qui est un objet de type fonction correspondant à la fonction de probabilité à utiliser pour le test :
| Code Python : | Sélectionner tout |
1 2 | # test si la série de valeurs suit une loi uniforme discrète resultat_test = serie_stat.test_khi2(serie_stat.uniforme_discrete, "Loi uniforme discrète") |
Arguments de la fonction :
- fonction_proba : objet de type fonction représentant la fonction de probabilité à utiliser pour le test : serie_stat.uniforme_discrete ;
- nom_loi : nom de la loi de probabilité prise en compte pour le test : "Loi uniforme discrète", "Loi de Poisson", etc. ;
- indice_debut et indice_fin : indices dans la liste de valeurs de la 1re et de la dernière classe prise en compte pour le test (arguments optionnels) ;
- degres_liberte: nombre de degrés de liberté (argument optionnel).
Elle renvoie un objet Test défini à partir d'une dataclass :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 10 | # dataclass utilisée pour définir un objet Test @dataclass class Test: type_test: str # type de test : "khi2", etc. loi_probabilite: str # loi de probabilité : "Loi de Poisson" moyenne: float # moyenne de la série statistique ecart_type: float # écart type de la série statistique degres_liberte: int # nombre de degrés de liberté pour le test T: float # valeur de T résultat du test du khi2 pvalue: float # p-value |
Exemple de résultat renvoyé :
| Code Python : | Sélectionner tout |
1 2 | Test(type_test='khi2', loi_probabilite='Loi de Poisson', moyenne=1.02, ecart_type=0.9162968951164245, degres_liberte=2, T=2.9714065378260637, pvalue=0.0847481400891209) |
On donne maintenant le code complet de la méthode test_khi2 à ajouter à la classe SerieStat :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 | class SerieStat(): ... def test_khi2(self, fonction_proba, nom_loi, indice_debut = None, indice_fin = None, degres_liberte=None): # si le nombre de degrés de liberté est omis alors on prend le nombre de valeurs différentes - 1 if degres_liberte is None: degres_liberte = len(valeurs)-1 # la série ne change pas : pas de réduction de la série de valeurs classes, effectifs_obs = (self.valeurs, self.effectifs) else: # sinon # regroupement des effectifs des valeurs mini et maxi de la série : classes = valeurs[indice_debut:indice_fin+1] classes, effectifs_obs = self.regrouper(indice_debut, indice_fin) # initialisation de la liste des effectifs attendus ou estimés effectifs_est=[] # moyenne pondérée de la série statistique m = self.moyenne_ponderee() # écart type de la série statistique s = self.ecart_type() # copie des paramètres du test actif dans l'attribut test_actif de la classe self.test_actif = Test('khi2', nom_loi, m, s, degres_liberte, None, None) # parcours des classes du milieu : classe d'indice 1 -> classe d'indice (n-1) for ci in classes[1:-1]: # calcul de l'effectif estimé pour xi : P(X = xi)*effectif_total ei = fonction_proba(ci)*self.effectif_total() effectifs_est.append(ei) # ajout de l'effetif à la liste # calcul de l'effectif estimé pour la 1re classe : P(X ≤ xi)*effectif_total e0 = fonction_proba(classes[0],cumulative=True)*self.effectif_total() # calcul de l'effectif estimé pour la dernière classe : (1 - P(X < xi))*effectif_total en = (1-fonction_proba(classes[-2],cumulative=True))*self.effectif_total() effectifs_est = [e0] + effectifs_est + [en] # résultat du test du khi carré : (T, pvalue) result_test = stats.chisquare(effectifs_obs, effectifs_est, ddof=degres_liberte) # copie de la valeur de T et de p-value dans self.test_actif self.test_actif.T = result_test[0] self.test_actif.pvalue = 1 - stats.chi2.cdf(self.test_actif.T, 1) # ou self.test_actif.pvalue = result_test[1] # copie du contenu de self.test_actif dans une variable resultat_test = self.test_actif self.test_actif = None # désactive le test # retourne le résultat du test sous la forme : Test(type_test='khi2', loi_probabilite='..', moyenne=.., ecart_type=.., degres_liberte=.., T=.., pvalue=..) return resultat_test |
La méthode test_khi2 utilise la fonction stats.chisquare(f_obs, f_esp, ddof) du module scipy.stats.
Le choix des classes et le regroupement des effectifs des valeurs mini et maxi de la série est réalisée à l'aide de la méthode regrouper() présente dans cette même classe.
Pour effectuer les tests nous ajoutons enfin ces lignes de code :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 10 11 12 | # série de valeurs avec les effectifs valeurs = [0, 1, 2, 3, 4] effectifs = [31, 45, 16, 7, 1] # création d'un objet SerieStat avec la série de valeurs et les effectifs serie_stat = SerieStat(valeurs, effectifs) # test si la série de valeurs suit une loi de Poisson avec regroupement des effectifs supérieurs ou égaux à 3 dans une même classe resultat_test = serie_stat.test_khi2(serie_stat.poisson, "Loi de Poisson", indice_debut=0, indice_fin=3, degres_liberte = degres_liberte) # affiche le résultat du test du khi carré print(resultat_test) |
Le code affiche :
Test(type_test='khi2', loi_probabilite='Loi de Poisson', moyenne=1.02, ecart_type=0.9162968951164245, degres_liberte=2, T=2.9714065378260637, pvalue=0.0847481400891209)
V. Module complet
Le module contenant la classe SerieStat avec en plus une méthode pour tracer sur un même graphique les points de la série statistique et ceux de la fonction de masse choisie :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 | from scipy import stats # librairie contenant les distributions de probabilité et la fonction de test du khi carré import matplotlib.pyplot as plt from dataclasses import dataclass # dataclass utilisée pour définir un objet Test @dataclass class Test: type_test: str # type de test : "khi2", etc. loi_probabilite: str # loi de probabilité : "Loi de Poisson" moyenne: float # moyenne de la série statistique ecart_type: float # écart type de la série statistique degres_liberte: int # nombre de degrés de liberté pour le test T: float # valeur de T résultat du test du khi2 pvalue: float # p-value class SerieStat(): def __init__(self, valeurs, effectifs): # méthode constructeur de la classe self.valeurs = valeurs # valeurs de la série statistique self.effectifs = effectifs # effectifs des valeurs self.test_actif = None # pas de test actif au moment de la création de l'objet def effectif_total(self): # renvoie l'effectif total des valeurs return sum(self.effectifs) def moyenne_ponderee(self): # calcul de la moyenne pondérée de la série de valeurs if self.test_actif: # si 1 test est actif : moyenne déjà calculée lors du test # renvoie la moyenne des valeurs du test actif return self.test_actif.moyenne else: # sinon # initialisation de la somme pondérée somme_ponderee = 0 # parcours des indices des valeurs for i in range(len(valeurs)): somme_ponderee += effectifs[i]*valeurs[i] # somme_ponderee = somme_ponderee + ni*xi # renvoie la moyenne pondérée : somme_ponderee / effectif_total return somme_ponderee/self.effectif_total() def ecart_type(self): # calcul de l'écart type de la série de valeurs if self.test_actif: # si 1 test est actif : écart type déjà calculé lors du test # renvoie l'écart type des valeurs du test actif return self.test_actif.ecart_type else: # moyenne pondérée m = self.moyenne_ponderee() # initialisation de la somme somme_ecarts = 0 # parcours des indices des valeurs for i in range(len(valeurs)): somme_ecarts += effectifs[i]*(valeurs[i] - m)**2 # somme_ecarts = somme_ecarts + ni*(xi-m)^2 # renvoie l'écart type : racine carrée de la somme obtenue divisée par l'effectif total return pow(somme_ecarts/self.effectif_total(),0.5) def uniforme_discrete(self, k, cumulative = False): # loi uniforme discrète pour k ϵ {1, 2, ..., n} : fonction de masse ou de répartition en k # nombre n de valeurs différentes n = len(self.valeurs) if not cumulative: # si c'est une fonction de masse. # retourne la probabilité P(X = k) return 1/n else: # sinon c'est une fonction de répartition # retourne la probabilité P(X ≤ k) return k/n def uniforme_continue(self, x, cumulative = False): # loi uniforme sur [a, b] : fonction cumulative ou non cumulative # intervalle entre 2 valeurs consécutives : u = (self.valeurs[-1]-self.valeurs[0])/(len(self.valeurs)-1) # intervalle total des valeurs [a, b] : [x0-u/2, xn + u/2] a = self.valeurs[0] - u/2; b = self.valeurs[-1] + u/2 if not cumulative: # si la fonction n'est pas cumulative. # valeurs de x1 et x2, avec [x1, x2] centré sur x x1 = x - u/2; x2 = x + u/2 # renvoie la probabilité P(x1 ≤ X ≤ x2) donnée par la valeur de la fonction de répartition entre x1 et x2 return stats.uniform.cdf(x2,loc=a,scale=b-a) - stats.uniform.cdf(x1,loc=a,scale=b-a) else: # sinon # renvoie la probabilité P(a ≤ X ≤ x+u/2) donnée par la valeur de la fonction de répartition entre a et x+u/2 return (stats.uniform.cdf(x+u/2,loc=a,scale=b-a) - stats.uniform.cdf(a,loc=a,scale=b-a)) def poisson(self, k, cumulative=False): # loi de Poisson de paramètre m : fonction de masse ou de répartition en k # moyenne pondérée de la série statistique m = self.moyenne_ponderee() if not cumulative: # si pas cumulative # retourne la valeur en k de la fonction de masse de paramètre m return stats.poisson.pmf(k,m) else: # sinon # retourne la valeur en k de la fonction de répartition de paramètre m return stats.poisson.cdf(k,m) def normale(self, x, cumulative=False): # loi Normale de moyenne m et d'écart type s : fonction de répartition en x # intervalle entre 2 valeurs consécutives : u = (self.valeurs[-1]-self.valeurs[0])/(len(self.valeurs)-1) # moyenne pondérée de la série statistique m = self.moyenne_ponderee() # écart type de la série statistique s = self.ecart_type() if not cumulative: # si la fonction n'est pas cumulative. # valeurs de x1 et x2, avec [x1, x2] centré sur x x1 = x - u/2; x2 = x + u/2 # retourne la valeur de la fonction de répartition de la loi normale N(m, e) entre x1 et x2 return (stats.norm.cdf(x2, loc=m, scale=s) - stats.norm.cdf(x1, loc=m, scale=s)) else: # sinon # retourne la valeur de la fonction de répartition de la loi normale N(m, e) en x+u/2 return stats.norm.cdf(x+u/2, loc=m, scale=s) def regrouper(self, indice_debut, indice_fin): # regroupement des valeurs mini et maxi de la série dans 2 classes (première et dernière classe de la série) # regroupement des effectifs des valeurs mini et maxi de la série # sélection des classes prise en compte classes = self.valeurs[indice_debut:indice_fin+1] # calculs des effectifs pour la première et la dernière classe e0 = sum(self.effectifs[0:indice_debut+1]) en = sum(self.effectifs[indice_fin:]) # ajout des effectifs en tête et en queue de liste effectifs = [e0] + self.effectifs[indice_debut+1:indice_fin] + [en] return (classes, effectifs) def test_khi2(self, fonction_proba, nom_loi, indice_debut = None, indice_fin = None, degres_liberte=None): # si le nombre de degrés de liberté est omis alors on prend le nombre de valeurs différentes - 1 if degres_liberte is None: degres_liberte = len(valeurs)-1 # la série ne change pas : pas de réduction de la série de valeurs classes, effectifs_obs = (self.valeurs, self.effectifs) else: # sinon # regroupement des effectifs des valeurs mini et maxi de la série : classes = valeurs[indice_debut:indice_fin+1] classes, effectifs_obs = self.regrouper(indice_debut, indice_fin) # initialisation de la liste des effectifs attendus ou estimés effectifs_est=[] # moyenne pondérée de la série statistique m = self.moyenne_ponderee() # écart type de la série statistique s = self.ecart_type() # copie des paramètres du test actif dans l'attribut test_actif de la classe self.test_actif = Test('khi2', nom_loi, m, s, degres_liberte, None, None) # parcours des classes du milieu : classe d'indice 1 -> classe d'indice (n-1) for ci in classes[1:-1]: # calcul de l'effectif estimé pour xi : P(X = xi)*effectif_total ei = fonction_proba(ci)*self.effectif_total() effectifs_est.append(ei) # ajout de l'effetif à la liste # calcul de l'effectif estimé pour la 1re classe : P(X ≤ xi)*effectif_total e0 = fonction_proba(classes[0],cumulative=True)*self.effectif_total() # calcul de l'effectif estimé pour la dernière classe : (1 - P(X < xi))*effectif_total en = (1-fonction_proba(classes[-2],cumulative=True))*self.effectif_total() effectifs_est = [e0] + effectifs_est + [en] # résultat du test du khi carré : (T, pvalue) result_test = stats.chisquare(effectifs_obs, effectifs_est, ddof=degres_liberte) # copie de la valeur de T et de p-value dans self.test_actif self.test_actif.T = result_test[0] self.test_actif.pvalue = 1 - stats.chi2.cdf(self.test_actif.T, 1) # ou self.test_actif.pvalue = result_test[1] # copie du contenu de self.test_actif dans une variable resultat_test = self.test_actif self.test_actif = None # désactive le test # retourne le résultat du test sous la forme : Test(type_test='khi2', loi_probabilite='..', moyenne=.., ecart_type=.., degres_liberte=.., T=.., pvalue=..) return resultat_test def khi2(self, p, degres_liberte): # nombre de degrés de liberté du test df = degres_liberte # valeur du khi-carré pour p T = stats.chi2.ppf(p, df) return T def tracer(self, fonction_proba, nom_loi): # moyenne pondérée de la série statistique m = self.moyenne_ponderee() # écart type de la série statistique s = self.ecart_type() # copie des paramètres du test actif dans l'attribut test_actif de la classe self.test_actif = Test('Comparaison graphique', nom_loi, m, s, None, None, None) # liste des valeurs x = self.valeurs # liste des effectifs y1 = self.effectifs # initialisation de la liste des effectifs attendus y2=[] for xi in x: # parcours des valeurs # calcul de l'effectif de xi : P(X = xi)*effectif_total yi = fonction_proba(xi)*self.effectif_total() y2.append(yi) # ajout de l'effectif à la liste # trace les points de la série statistique (en bleu) et de la fonction de densité (en rouge) plt.scatter(x, y1, c = 'blue', marker = 'o' , s = 50, label = 'Série') plt.scatter(x, y2, c = 'red', marker = 'o' , s = 50, label = nom_loi) # relie les points de la série statistique (en bleu) et de la fonction de densité (en rouge) plt.plot(x, y1, c = 'grey') plt.plot(x, y2, c = 'grey') # ajout d'étiquettes sur l'axe des x et des y plt.xlabel('Valeurs') plt.ylabel('Effectifs') # ajoute une légende en haut à droite du graphique plt.legend(loc="upper right") # ajoute un titre plt.title("Ajustement d'une série de points à une " + nom_loi) self.test_actif = None # désactive le test # affiche le graphique obtenu plt.show() # Test d'adéquation d'une série statistique à une loi uniforme print("I. Test d'adéquation d'une série statistique à une loi uniforme : \n") # série de valeurs avec les effectifs valeurs = [1, 2, 3, 4, 5, 6] effectifs = [88, 109, 107, 94, 105, 97] # création d'un objet SerieStat avec la série de valeurs et les effectifs serie_stat = SerieStat(valeurs, effectifs) # nombre de degrés de liberté degres_liberte = 5 # on teste si la série de valeurs suit une loi uniforme discrète resultat_test = serie_stat.test_khi2(serie_stat.uniforme_discrete, "Loi uniforme discrète") # affiche le résultat du test du khi carré print(resultat_test) # prise en compte du risque α α = 0.05 # valeur complémentaire p = 1 - α # calcul du Khi2, valeur en deçà de laquelle on considère la série statistique comme conforme avec un risque α=(1-p). t = serie_stat.khi2(p, degres_liberte) print() # affiche la valeur de t telle que : P(T < t) = 0.95 print("P(T < " + str(t) + ") = " + str(p)) # Test d'adéquation d'une série statistique à une loi de Poisson print();print() print("II. Test d'adéquation d'une série statistique à une loi de Poisson : \n") # série de valeurs avec les effectifs valeurs = [0, 1, 2, 3, 4] effectifs = [31, 45, 16, 7, 1] # création d'un objet SerieStat avec la série de valeurs et les effectifs serie_stat = SerieStat(valeurs, effectifs) # nombre de degrés de liberté degres_liberte = 2 # on teste si la série de valeurs suit une loi de Poisson avec regroupement des effectifs des valeurs supérieurs ou égales à 3 dans une même classe resultat_test = serie_stat.test_khi2(serie_stat.poisson, "Loi de Poisson", indice_debut=0, indice_fin=3, degres_liberte = degres_liberte) # affiche le résultat du test du khi carré print(resultat_test) # prise en compte du risque α = 0.05 # valeur complémentaire p = 1 - α # calcul du Khi2, valeur en deçà de laquelle on considère la série statistique comme conforme avec un risque α. t = serie_stat.khi2(p, degres_liberte) print() # affiche la valeur de t telle que : P(T < t) = 0.95 print("P(T < " + str(t) + ") = " + str(p)) # Trace les points de la série de valeurs comparés à ceux donnés par la fonction de masse de la loi de Poisson serie_stat.tracer(serie_stat.poisson,'Loi de Poisson') |
Le code affiche :
Puis, génère le graphique à partir de la série statistique et de la fonction de probabilité :
VI. Conclusion
Après avoir montré comment utiliser le test du khi carré pour vérifier qu'une série de données statistiques suit bien une loi de probabilité, nous avons pu implémenter ce test à l'aide de la classe SerieStat et de ses différentes méthodes.
Chacun pourra ensuite ajouter de nouvelles méthodes à cette classe pour par exemple effectuer des tests en utilisant d'autres lois de probabilité.
Sources :
https://fr.wikipedia.org/wiki/Test_du_%CF%87%C2%B2
https://fr.wikipedia.org/wiki/Valeur_p
https://fr.wikipedia.org/wiki/Degr%C...(statistiques)
https://fr.wikipedia.org/wiki/Loi_un..._discr%C3%A8te
https://fr.wikipedia.org/wiki/Loi_uniforme_continue
https://fr.wikipedia.org/wiki/Loi_de_Poisson
https://docs.scipy.org/doc/scipy/reference/stats.html
https://docs.python.org/fr/3/tutorial/classes.html