Cours 6: tests

import numpy as np
import matplotlib.pyplot as plt
from scipy.special import binom, erf, erfinv
import scipy.stats

exercice

Test d'ajustement du chi2 pour le dé

comptes = np.array([7, 8, 12, 15,10,8])
n = sum(comptes)
n

60

# statistiques de test donnee par la propriete du chi2
z_obs = (((comptes - 10)**2)/10).sum()
# p valeur
deg_liberte = 6-1
p_valeur = 1 - scipy.stats.chi2.cdf(z_obs, deg_liberte)
print(f'p_valeur = {100*p_valeur:.2f}%')

p_valeur = 46.66%

# densite et fonction de repartition du chi2
x = np.linspace(0, 10, 1000)
d = 3
f_chi2 = scipy.stats.chi2.pdf(x, d)
F_chi2 = scipy.stats.chi2.cdf(x, d)
plt.figure()
plt.plot(x, f_chi2 , label=f'densite')
plt.plot(x, F_chi2 , label=f'fonc. de rep.')
plt.title(f'densite de fonction de rep de loi chi2 a {d} deg. de lib.')
plt.legend()
plt.show()

# echantilloner une multinomiale
n = 60
pi = [1/6, 1/6, 1/6, 1/6, 1/6, 1/6]
m = np.random.multinomial(n, pi)

Correspondance entre le niveau alpha d'un test et le taux de faux negatifs en pratique

pi = np.ones(6) / 6 # probabilites de la loi
n = 100
N = 1000000
for alpha in [0.001, 0.01, 0.05]:
    false_negative_rate = (scipy.stats.chisquare(np.random.multinomial(n, pi, N), axis=1)[1] < alpha).mean()
    print(f'alpha={100*alpha:.3f} %, false_negative_rate={100*false_negative_rate:.3f} %')

alpha=0.100 %, false_negative_rate=0.103 %
alpha=1.000 %, false_negative_rate=0.967 %
alpha=5.000 %, false_negative_rate=4.905 %

exercice 2

Test d'ajustement du chi2 pour loi de poisson

# les donnees
values = np.arange(6)
effectifs = np.array([16, 20, 8 , 5, 2, 1])

n = effectifs.sum()
# on estime lambda
lambd = (values* effectifs).sum() / n
# on calcule les effectifs moyens
effectifs_moyens_h_0 = np.zeros(6,dtype='float')

# n * P(X=i) pour i = 0 ... 4
for i in range(5):
    effectifs_moyens_h_0[i] = 52 * np.exp(-lambd) * lambd**i / np.math.factorial(i) 
# n * P(X>=5) 
for i in range(5, 100):
    effectifs_moyens_h_0[5] += 52 * np.exp(-lambd) * lambd**i / np.math.factorial(i)

print(f'n = {n}, lambda = {lambd}, effectifs theoriques {effectifs_moyens_h_0}')

n = 52, lambda = 1.2307692307692308, effectifs theoriques [15.18752683 18.69234072 11.5029789   4.71917083  1.45205256  0.44593015]

# on regroupe les trois derniers groupes pour avoir un effectif theorique d'au moins 5
observed = np.concatenate([effectifs[:-3], [effectifs[-3:].sum()]])
expected =  np.concatenate([effectifs_moyens_h_0[:-3], [effectifs_moyens_h_0[-3:].sum()]])
with np.printoptions(precision=1):
    print(observed, expected)

[16 20  8  8] [15.2 18.7 11.5  6.6]

# on fait le test de chi2 en precisant qu'il y a un deg. de librete en moins
scipy.stats.chisquare(observed, expected, ddof=1)

Power_divergenceResult(statistic=1.4906850113518044, pvalue=0.47457172859003194)