Calculer un chi 2 sur un tableau de contingence#

$\chi_2$ et tableau de contingence, avec numpy, avec scipy ou sans.

formule#

Le test du $\chi_2$ (wikipedia) sert à comparer deux distributions. Il peut être appliqué sur un tableau de contingence pour comparer la distributions observée avec la distribution qu’on observerait si les deux facteurs du tableau étaient indépendants. On note $M=(m_{ij})$ une matrice de dimension $I \times J$ . Le test du $\chi_2$ se calcule comme suit :

$M = \sum_{ij} m_{ij}$
$\forall i, \; m_{i \bullet} = \sum_j m_{ij}$
$\forall j, \; m_{\bullet j} = \sum_i m_{ij}$
$\forall i,j \; n_{ij} = \frac{m_{i \bullet} m_{\bullet j}}{N}$

Avec ces notations :

$T = \sum_{ij} \frac{ (m_{ij} - n_{ij})^2}{n_{ij}}$

La variable aléatoire $T$ suit asymptotiquement une loi du $\chi_2$ à $(I-1)(J-1)$ degrés de liberté (table). Comment le calculer avec numpy ?

tableau au hasard#

On prend un petit tableau qu’on choisit au hasard, de préférence non carré pour détecter des erreurs de calculs.

[1]:

import numpy

M = numpy.array([[4, 5, 2, 1], [6, 3, 1, 7], [10, 14, 6, 9]])
M

[1]:

array([[ 4,  5,  2,  1],
       [ 6,  3,  1,  7],
       [10, 14,  6,  9]])

calcul avec scipy#

Evidemment, il existe une fonction en python qui permet de calculer la statistique $T$ : chi2_contingency.

[2]:

from scipy.stats import chi2_contingency

chi2, pvalue, degrees, expected = chi2_contingency(M)
chi2, degrees, pvalue

[2]:

(6.168598503892621, 6, 0.4045712090580829)

calcul avec numpy#

[3]:

N = M.sum()
ni = numpy.array([M[i, :].sum() for i in range(M.shape[0])])
nj = numpy.array([M[:, j].sum() for j in range(M.shape[1])])
ni, nj, N

[3]:

(array([12, 17, 39]), array([20, 22,  9, 17]), 68)

Et comme c’est un usage courant, numpy propose une façon de faire sans écrire une boucle avec la fonction sum :

[4]:

ni = M.sum(axis=1)
nj = M.sum(axis=0)
ni, nj, N

[4]:

(array([12, 17, 39]), array([20, 22,  9, 17]), 68)

[5]:

nij = ni.reshape(M.shape[0], 1) * nj / N
nij

[5]:

array([[ 3.52941176,  3.88235294,  1.58823529,  3.        ],
       [ 5.        ,  5.5       ,  2.25      ,  4.25      ],
       [11.47058824, 12.61764706,  5.16176471,  9.75      ]])

[6]:

d = (M - nij) ** 2 / nij
d.sum()

[6]:

6.168598503892621

[ ]:

Notebook on github