I. Introduction
En 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.
638302
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 :
638459
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 :
638460
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 :
638461
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 :
638371
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].
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 :
638372
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).
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).
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 :
Valeurs retournées (T, pvalue) :
# 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. :
638350
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 :
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 :
# 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 :
def fonction_proba(self, x, cumulative = False):
...
return ..
Arguments des fonctions de probabilité :
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 :
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)/(len(self.valeurs)-1)
# intervalle total des valeurs [a, b] : [x0-u/2, xn + u/2]
a = self.valeurs - 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 :
# 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 :
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 :
# 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 :
# 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 :
Elle renvoie un objet Test défini à partir d'une 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
Exemple de résultat renvoyé :
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 :
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,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 = + effectifs_est +
# 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
self.test_actif.pvalue = 1 - stats.chi2.cdf(self.test_actif.T, 1) # ou self.test_actif.pvalue = result_test
# 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 :
# 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 :
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*valeurs # 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*(valeurs - 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)/(len(self.valeurs)-1)
# intervalle total des valeurs [a, b] : [x0-u/2, xn + u/2]
a = self.valeurs - 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)/(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 = + self.effectifs[indice_debut+1:indice_fin] +
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,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 = + effectifs_est +
# 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
self.test_actif.pvalue = 1 - stats.chi2.cdf(self.test_actif.T, 1) # ou self.test_actif.pvalue = result_test
# 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 :
638440
Puis, génère le graphique à partir de la série statistique et de la fonction de probabilité :
638398
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%C3%A9_de_libert%C3%A9_(statistiques)
https://fr.wikipedia.org/wiki/Loi_uniforme_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
Soutenez le club developpez.com en souscrivant un abonnement pour que nous puissions continuer à vous proposer des publications.