User

I. Introduction

En statistique, le test du khi carré, aussi dit du khi-deux, d’après sa désignation symbolique χ², est un test statistique où la statistique de test suit une loi du χ² 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 χ²



II. Principe du test d'adéquation

Le test du χ² 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 H₀) 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 H₀ 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)²/100 + (109 − 100)²/100 + (107 − 100²/100 + (94 − 100)²/100 + (105 − 100)²/100 + (97 − 100)²/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 χ² à cinq degrés de liberté.

Cette loi du χ² 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 χ² à 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 - Visualiser dans une fenêtre à part

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 - Visualiser dans une fenêtre à part

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 - Visualiser dans une fenêtre à part

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 :

1. f_obs : tableau des effectifs observés
2. 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é.
3. ddof : valeur entière représentant le nombre de degrés de liberté.

Valeurs retournées (T, pvalue) :

1. T : valeur de type float résultat du test du khi carré.
2. 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 - Visualiser dans une fenêtre à part

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 - Visualiser dans une fenêtre à part

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 - Visualiser dans une fenêtre à part

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 - Visualiser dans une fenêtre à part

1
2
3
4
 
    def fonction_proba(self, x, cumulative = False):
	...
        return ..

Arguments des fonctions de probabilité :

1. x ou k: valeur pour laquelle on souhaite calculer la probabilité ;
2. 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 - Visualiser dans une fenêtre à part

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 - Visualiser dans une fenêtre à part

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 - Visualiser dans une fenêtre à part

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 - Visualiser dans une fenêtre à part

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 - Visualiser dans une fenêtre à part

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 :

1. fonction_proba : objet de type fonction représentant la fonction de probabilité à utiliser pour le test : serie_stat.uniforme_discrete ;
2. nom_loi : nom de la loi de probabilité prise en compte pour le test : "Loi uniforme discrète", "Loi de Poisson", etc. ;
3. 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) ;
4. 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 - Visualiser dans une fenêtre à part

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 - Visualiser dans une fenêtre à part

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 - Visualiser dans une fenêtre à part

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 - Visualiser dans une fenêtre à part

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 - Visualiser dans une fenêtre à part

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