Tests paramétriques

Les tests paramétriques constituent l’arsenal classique de l’inférence statistique : ils permettent de répondre à des questions concrètes — « ce traitement est-il efficace ? », « ces deux machines produisent-elles des pièces de même calibre ? », « les trois méthodes d’enseignement donnent-elles des résultats équivalents ? » — avec une rigueur formelle. Ce chapitre présente les tests les plus utilisés, leurs hypothèses, leur implémentation en Python et leur interprétation correcte.

Show code cell source

Hide code cell source

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec
import seaborn as sns
from scipy import stats
import statsmodels.api as sm
import statsmodels.formula.api as smf
from statsmodels.stats.anova import AnovaRM
from statsmodels.stats import power as smp
import pingouin as pg
import warnings
warnings.filterwarnings('ignore')

rng = np.random.default_rng(42)

plt.rcParams.update({
    'figure.dpi': 110,
    'axes.spines.top': False,
    'axes.spines.right': False,
    'font.size': 11,
})
sns.set_palette('husl')

Rappels : cadre des tests d’hypothèses

Un test statistique repose sur une logique de réfutation : on formule une hypothèse nulle \(H_0\) (statu quo, absence d’effet), on calcule la probabilité d’observer des données au moins aussi extrêmes si \(H_0\) était vraie, et on décide de rejeter ou non \(H_0\).

Éléments fondamentaux :

• \(H_0\) (hypothèse nulle) : hypothèse par défaut, que l’on cherche à réfuter. Exemple : \(\mu = \mu_0\), ou \(\mu_1 = \mu_2\).

• \(H_1\) (hypothèse alternative) : ce que l’on cherche à établir. Peut être bilatérale (\(\mu \neq \mu_0\)) ou unilatérale (\(\mu > \mu_0\)).

• p-valeur : \(P(\text{statistique} \geq t_{\text{obs}} \mid H_0)\). Ce n’est pas la probabilité que \(H_0\) soit vraie — c’est la probabilité d’obtenir ce résultat ou plus extrême si \(H_0\) était vraie.

• Niveau \(\alpha\) : seuil de rejet fixé avant le test (souvent 0,05). Si \(p < \alpha\), on rejette \(H_0\).

Les deux types d’erreurs :

Décision \ Réalité	\(H_0\) vraie	\(H_0\) fausse
Ne pas rejeter \(H_0\)	Décision correcte	Erreur de type II (\(\beta\))
Rejeter \(H_0\)	Erreur de type I (\(\alpha\))	Décision correcte

La puissance du test est \(1 - \beta\) : la probabilité de détecter un effet réel.

p-valeur : ce qu’elle n’est pas

La p-valeur n’est pas :

• La probabilité que \(H_0\) soit vraie

• La probabilité que les résultats soient dus au hasard

• Une mesure de l’importance pratique de l’effet

Une p-valeur < 0,05 signifie seulement que les données sont incompatibles avec \(H_0\) au seuil \(\alpha = 5\%\). La taille d’effet mesure l’importance pratique.

Test t à un échantillon

Contexte et hypothèses

On dispose d’un échantillon \(x_1, \ldots, x_n\) et on veut tester si la moyenne de la population dont il est issu est égale à une valeur de référence \(\mu_0\).

\[H_0 : \mu = \mu_0 \qquad H_1 : \mu \neq \mu_0\]

La statistique de test est : $\(t = \frac{\bar{x} - \mu_0}{s / \sqrt{n}}\)$

Sous \(H_0\), \(t \sim \mathcal{T}(n-1)\) (loi de Student à \(n-1\) degrés de liberté), à condition que les observations soient indépendantes et que la population soit normale (ou \(n\) suffisamment grand pour le TCL).

Exemple : calibrage d’une machine

Une machine est censée produire des pièces de 50 mm. On mesure 30 pièces. La machine est-elle bien calibrée ?

# Simulation : machine légèrement déréglée
mu_reel = 50.8  # vrai écart
sigma = 2.0
n = 30
mesures = rng.normal(mu_reel, sigma, n)

# Test t à un échantillon
t_stat, p_val = stats.ttest_1samp(mesures, popmean=50.0)

print(f"Moyenne observée : {mesures.mean():.3f} mm")
print(f"Statistique t    : {t_stat:.3f}")
print(f"p-valeur         : {p_val:.4f}")
print(f"Degrés de liberté: {n - 1}")
print()
if p_val < 0.05:
    print("→ Rejet de H0 : la machine est probablement déréglée.")
else:
    print("→ Non-rejet de H0 : pas de preuve de dérèglage.")

Moyenne observée : 50.834 mm
Statistique t    : 2.939
p-valeur         : 0.0064
Degrés de liberté: 29

→ Rejet de H0 : la machine est probablement déréglée.

Taille d’effet : d de Cohen

La p-valeur seule est insuffisante. Le d de Cohen mesure l’ampleur de l’écart en unités d’écarts-types :

\[d = \frac{\bar{x} - \mu_0}{s}\]

Valeur absolue de \(d\)	Interprétation
0,2	Petit effet
0,5	Effet modéré
0,8	Grand effet

# Taille d'effet
cohens_d = (mesures.mean() - 50.0) / mesures.std(ddof=1)
print(f"d de Cohen : {cohens_d:.3f}")

# Avec pingouin (résultat complet)
result = pg.ttest(mesures, 50.0)
print("\nRésultat complet avec pingouin :")
print(result[['T', 'dof', 'alternative', 'p_val', 'CI95', 'cohen_d', 'power']].to_string())

d de Cohen : 0.537

Résultat complet avec pingouin :
               T  dof alternative     p_val            CI95   cohen_d     power
T_test  2.939154   29   two-sided  0.006397  [50.25, 51.41]  0.536614  0.810742

Visualisation : distribution sous H₀ et statistique observée

Show code cell source

Hide code cell source

fig, axes = plt.subplots(1, 2, figsize=(13, 4))

# Axe gauche : données brutes
ax = axes[0]
ax.hist(mesures, bins=10, color='steelblue', alpha=0.7, edgecolor='white')
ax.axvline(50.0, color='crimson', lw=2, linestyle='--', label='Référence μ₀ = 50')
ax.axvline(mesures.mean(), color='navy', lw=2, label=f'Moyenne obs. = {mesures.mean():.2f}')
ax.set_xlabel('Longueur (mm)')
ax.set_ylabel('Effectif')
ax.set_title('Distribution des mesures')
ax.legend()

# Axe droit : distribution de Student sous H0 avec région de rejet
ax = axes[1]
df_val = n - 1
x = np.linspace(-4, 4, 400)
y = stats.t.pdf(x, df=df_val)
ax.plot(x, y, 'k-', lw=2, label=f'Student({df_val})')

# Régions de rejet bilatérales (alpha = 0.05)
t_crit = stats.t.ppf(0.975, df=df_val)
mask_left = x < -t_crit
mask_right = x > t_crit
ax.fill_between(x[mask_left], y[mask_left], alpha=0.4, color='crimson', label='Région de rejet (5%)')
ax.fill_between(x[mask_right], y[mask_right], alpha=0.4, color='crimson')

# Statistique observée
ax.axvline(t_stat, color='navy', lw=2.5, linestyle='--', label=f't obs. = {t_stat:.2f}')
ax.set_xlabel('t')
ax.set_ylabel('Densité')
ax.set_title('Distribution sous H₀ et statistique de test')
ax.legend(fontsize=9)

plt.tight_layout()
plt.savefig('_static/06_ttest1.png', bbox_inches='tight')
plt.show()

Test t à deux échantillons indépendants

Student vs Welch

Deux groupes indépendants \(x_1, \ldots, x_{n_1}\) et \(y_1, \ldots, y_{n_2}\). On teste \(H_0 : \mu_X = \mu_Y\).

Test de Student (variantes égales) : \(t = \frac{\bar{x} - \bar{y}}{s_p \sqrt{1/n_1 + 1/n_2}}\) avec \(s_p^2\) variance poolée.

Test de Welch (variances inégales, recommandé par défaut) : ajuste les degrés de liberté selon l’approximation de Welch-Satterthwaite.

Quelle variante utiliser ?

En pratique, utilisez toujours le test de Welch (equal_var=False dans scipy). Il est robuste même lorsque les variances sont égales, alors que le test de Student peut être anti-conservateur lorsque les variances diffèrent. Vérifier l’égalité des variances avec un test de Levene et adapter en conséquence est une bonne pratique.

Exemple : deux méthodes pédagogiques

# Deux groupes d'élèves (scores sur 100)
n_A, n_B = 35, 40
groupe_A = rng.normal(68, 12, n_A)  # méthode traditionnelle
groupe_B = rng.normal(74, 15, n_B)  # méthode active

# Test de Levene (égalité des variances)
lev_stat, lev_p = stats.levene(groupe_A, groupe_B)
print(f"Test de Levene : F = {lev_stat:.3f}, p = {lev_p:.4f}")
print(f"→ Variances {'probablement égales' if lev_p > 0.05 else 'probablement inégales'}")
print()

# Test de Welch (recommandé)
t_welch, p_welch = stats.ttest_ind(groupe_A, groupe_B, equal_var=False)
# Test de Student (pour comparaison)
t_student, p_student = stats.ttest_ind(groupe_A, groupe_B, equal_var=True)

print(f"Test de Welch   : t = {t_welch:.3f}, p = {p_welch:.4f}")
print(f"Test de Student : t = {t_student:.3f}, p = {p_student:.4f}")
print()

# Avec pingouin
res = pg.ttest(groupe_A, groupe_B, correction=True)
print("Résultat pingouin (Welch) :")
print(res[['T', 'dof', 'p_val', 'CI95', 'cohen_d', 'power']].to_string())

Test de Levene : F = 0.760, p = 0.3860
→ Variances probablement égales

Test de Welch   : t = -0.835, p = 0.4064
Test de Student : t = -0.828, p = 0.4101

Résultat pingouin (Welch) :
               T        dof     p_val           CI95   cohen_d     power
T_test -0.835156  72.983305  0.406356  [-6.82, 2.79]  0.191745  0.129384

Show code cell source

Hide code cell source

fig, axes = plt.subplots(1, 2, figsize=(13, 4))

# Boxplot + swarmplot
ax = axes[0]
df_scores = pd.DataFrame({
    'Score': np.concatenate([groupe_A, groupe_B]),
    'Méthode': ['Traditionnelle'] * n_A + ['Active'] * n_B
})
sns.boxplot(data=df_scores, x='Méthode', y='Score', ax=ax,
            palette=['#4C72B0', '#DD8452'], width=0.5)
sns.stripplot(data=df_scores, x='Méthode', y='Score', ax=ax,
              color='black', alpha=0.4, size=4, jitter=True)
ax.set_title('Distribution des scores par méthode')
ax.set_ylabel('Score (/100)')

# Intervalles de confiance sur les moyennes
ax = axes[1]
for i, (groupe, label, color) in enumerate(zip(
    [groupe_A, groupe_B], ['Traditionnelle', 'Active'],
    ['#4C72B0', '#DD8452']
)):
    mean = groupe.mean()
    sem = stats.sem(groupe)
    ci = stats.t.interval(0.95, len(groupe)-1, loc=mean, scale=sem)
    ax.errorbar(i, mean, yerr=[[mean - ci[0]], [ci[1] - mean]],
                fmt='o', color=color, capsize=8, capthick=2,
                markersize=10, linewidth=2, label=label)

ax.set_xticks([0, 1])
ax.set_xticklabels(['Traditionnelle', 'Active'])
ax.set_ylabel('Score moyen ± IC 95%')
ax.set_title('Comparaison des moyennes avec IC 95%')
ax.set_xlim(-0.5, 1.5)
ax.legend()

plt.tight_layout()
plt.savefig('_static/06_ttest2.png', bbox_inches='tight')
plt.show()

Test t apparié

Données appariées : avant/après

Lorsque les deux mesures sont liées (même individu mesuré deux fois, jumeaux, paires matchées), les observations ne sont pas indépendantes. Le test t apparié exploite cette structure en travaillant sur les différences \(d_i = x_{i,\text{après}} - x_{i,\text{avant}}\).

\[H_0 : \mu_d = 0 \qquad t = \frac{\bar{d}}{s_d / \sqrt{n}}\]

Le test t apparié est beaucoup plus puissant que le test à deux échantillons indépendants quand la corrélation intra-paire est élevée, car il élimine la variabilité inter-individus.

# Exemple médical : pression artérielle avant/après traitement
n_patients = 25
avant = rng.normal(145, 15, n_patients)  # mmHg, avant traitement
# Effet du traitement : réduction de ~10 mmHg + variabilité individuelle
effet_traitement = rng.normal(10, 5, n_patients)
apres = avant - effet_traitement

differences = apres - avant

# Test t apparié
t_app, p_app = stats.ttest_rel(avant, apres)

# Test t indépendant (incorrect ici, pour illustration)
t_ind, p_ind = stats.ttest_ind(avant, apres)

print("Test t apparié (correct) :")
print(f"  t = {t_app:.3f}, p = {p_app:.6f}")
print()
print("Test t indépendant (incorrect pour ces données) :")
print(f"  t = {t_ind:.3f}, p = {p_ind:.4f}")
print()
print(f"Différence moyenne : {differences.mean():.2f} mmHg")
print(f"IC 95% sur la différence : {np.percentile(differences, [2.5, 97.5])}")
print()

res_app = pg.ttest(avant, apres, paired=True)
print("Résultat pingouin (apparié) :")
print(res_app[['T', 'dof', 'p_val', 'CI95', 'cohen_d', 'power']].to_string())

Test t apparié (correct) :
  t = 8.092, p = 0.000000

Test t indépendant (incorrect pour ces données) :
  t = 2.370, p = 0.0219

Différence moyenne : -9.12 mmHg
IC 95% sur la différence : [-21.81791761  -0.7118586 ]

Résultat pingouin (apparié) :
               T  dof         p_val          CI95   cohen_d     power
T_test  8.092432   24  2.573756e-08  [6.8, 11.45]  0.670261  0.895176

Show code cell source

Hide code cell source

fig, axes = plt.subplots(1, 3, figsize=(15, 4))

# Paires avant/après
ax = axes[0]
for i in range(n_patients):
    ax.plot([0, 1], [avant[i], apres[i]], 'o-', color='steelblue', alpha=0.4, linewidth=1)
ax.plot([0, 1], [avant.mean(), apres.mean()], 'o-', color='crimson',
        linewidth=3, markersize=10, label='Moyenne')
ax.set_xticks([0, 1])
ax.set_xticklabels(['Avant', 'Après'])
ax.set_ylabel('Pression artérielle (mmHg)')
ax.set_title('Pression artérielle par patient')
ax.legend()

# Distribution des différences
ax = axes[1]
ax.hist(differences, bins=10, color='steelblue', alpha=0.7, edgecolor='white')
ax.axvline(0, color='crimson', lw=2, linestyle='--', label='H₀ : diff = 0')
ax.axvline(differences.mean(), color='navy', lw=2, label=f'Moyenne = {differences.mean():.1f}')
ax.set_xlabel('Différence (Après − Avant) en mmHg')
ax.set_ylabel('Effectif')
ax.set_title('Distribution des différences')
ax.legend()

# Comparaison puissance : apparié vs indépendant
ax = axes[2]
n_range = np.arange(5, 50)
# Approximation : effet moyen / std différences vs effet moyen / std poolée
d_apparie = abs(differences.mean()) / differences.std(ddof=1)
sigma_pool = np.sqrt((avant.var(ddof=1) + apres.var(ddof=1)) / 2)
d_indep = abs(differences.mean()) / sigma_pool

pow_app = [smp.TTestPower().power(effect_size=d_apparie, nobs=n, alpha=0.05) for n in n_range]
pow_ind = [smp.TTestIndPower().power(effect_size=d_indep, nobs1=n, alpha=0.05, ratio=1) for n in n_range]

ax.plot(n_range, pow_app, 'steelblue', lw=2, label='T apparié')
ax.plot(n_range, pow_ind, 'crimson', lw=2, linestyle='--', label='T indépendant')
ax.axhline(0.8, color='gray', linestyle=':', lw=1.5, label='Puissance 80%')
ax.set_xlabel('Taille d\'échantillon n')
ax.set_ylabel('Puissance')
ax.set_title('Puissance : apparié vs indépendant')
ax.legend(fontsize=9)
ax.set_ylim(0, 1.05)

plt.tight_layout()
plt.savefig('_static/06_apparie.png', bbox_inches='tight')
plt.show()

ANOVA à un facteur

Principe

L’ANOVA (ANalysis Of VAriance) étend le test t à plus de deux groupes. On décompose la variance totale en variance inter-groupes (due aux traitements) et variance intra-groupes (résiduelle).

\[F = \frac{SS_{\text{inter}} / (k-1)}{SS_{\text{intra}} / (N-k)} = \frac{MS_{\text{entre}}}{MS_{\text{intra}}}\]

où \(k\) est le nombre de groupes et \(N\) le nombre total d’observations.

Hypothèses de l’ANOVA :

1. Les observations sont indépendantes

2. Les résidus sont normalement distribués

3. Les variances sont homogènes (homoscédasticité) — vérifiée par Levene ou Bartlett

Pourquoi ne pas faire plusieurs tests t ?

Avec \(k\) groupes, faire tous les tests t par paires donne \(\binom{k}{2}\) tests. Au seuil \(\alpha = 0,05\), la probabilité d’au moins une fausse alarme est \(1 - (1-\alpha)^m\) où \(m\) est le nombre de tests. Pour \(k = 5\) groupes (10 tests), ce taux atteint ~40%. L’ANOVA contrôle globalement \(\alpha\).

Exemple : trois régimes alimentaires

# Trois régimes, effet sur la glycémie à jeun (mg/dL)
regime_A = rng.normal(95, 10, 30)   # régime standard
regime_B = rng.normal(88, 10, 30)   # régime faible glucides
regime_C = rng.normal(92, 10, 30)   # régime méditerranéen

# ANOVA avec scipy
f_stat, p_anova = stats.f_oneway(regime_A, regime_B, regime_C)
print(f"ANOVA à un facteur :")
print(f"  F = {f_stat:.3f}")
print(f"  p = {p_anova:.4f}")
print()

# Table ANOVA complète avec pingouin
df_glyc = pd.DataFrame({
    'glycemie': np.concatenate([regime_A, regime_B, regime_C]),
    'regime': ['Standard'] * 30 + ['Faible glucides'] * 30 + ['Méditerranéen'] * 30
})

aov = pg.anova(data=df_glyc, dv='glycemie', between='regime', detailed=True)
print("Table ANOVA (pingouin) :")
print(aov.to_string(index=False))
print()
print(f"η² (eta-carré, taille d'effet) : {aov.loc[0, 'np2']:.3f}")

ANOVA à un facteur :
  F = 5.911
  p = 0.0039

Table ANOVA (pingouin) :
Source          SS  DF         MS        F    p_unc      np2
regime 1227.941587   2 613.970794 5.911488 0.003915 0.119638
Within 9035.874239  87 103.860623      NaN      NaN      NaN

η² (eta-carré, taille d'effet) : 0.120

Tests post-hoc : Tukey HSD et Bonferroni

L’ANOVA indique si au moins une paire de moyennes diffère — les tests post-hoc précisent lesquelles.

# Tests post-hoc avec pingouin
print("Post-hoc Tukey HSD :")
tukey = pg.pairwise_tukey(data=df_glyc, dv='glycemie', between='regime')
print(tukey[['A', 'B', 'mean_A', 'mean_B', 'diff', 'se', 'T', 'p_tukey']].to_string(index=False))
print()

print("Post-hoc Bonferroni :")
bonf = pg.pairwise_tests(data=df_glyc, dv='glycemie', between='regime', padjust='bonf')
print(bonf[['A', 'B', 'T', 'dof', 'p_unc', 'p_corr']].to_string(index=False))

Post-hoc Tukey HSD :
              A             B    mean_A    mean_B      diff       se         T  p_tukey
Faible glucides Méditerranéen 87.130916 90.942226 -3.811310 2.631357 -1.448420 0.320880
Faible glucides      Standard 87.130916 96.143083 -9.012167 2.631357 -3.424912 0.002681
  Méditerranéen      Standard 90.942226 96.143083 -5.200857 2.631357 -1.976492 0.124127

Post-hoc Bonferroni :
              A             B         T  dof    p_unc   p_corr
Faible glucides Méditerranéen -1.385644 58.0 0.171161 0.513483
Faible glucides      Standard -3.635859 58.0 0.000590 0.001769
  Méditerranéen      Standard -1.957010 58.0 0.055165 0.165496

Show code cell source

Hide code cell source

fig, axes = plt.subplots(1, 3, figsize=(16, 4.5))

# Boxplot comparatif
ax = axes[0]
palette = {'Standard': '#4C72B0', 'Faible glucides': '#55A868', 'Méditerranéen': '#C44E52'}
sns.boxplot(data=df_glyc, x='regime', y='glycemie', ax=ax, palette=palette)
sns.stripplot(data=df_glyc, x='regime', y='glycemie', ax=ax,
              color='black', alpha=0.3, size=3)
ax.set_title('Glycémie par régime')
ax.set_xlabel('Régime')
ax.set_ylabel('Glycémie à jeun (mg/dL)')
ax.tick_params(axis='x', rotation=10)

# Décomposition de la variance (barres)
ax = axes[1]
ss_inter = aov.loc[0, 'SS']
ss_intra = aov.loc[1, 'SS']
ss_total = ss_inter + ss_intra
bars = ax.bar(['SS inter\n(régimes)', 'SS intra\n(résidus)'],
              [ss_inter, ss_intra],
              color=['#C44E52', '#4C72B0'], alpha=0.8)
ax.bar_label(bars, fmt='%.1f', padding=3)
ax.set_ylabel('Somme des carrés')
ax.set_title(f'Décomposition de la variance\n(η² = {ss_inter/ss_total:.2f})')

# Distribution F sous H0 avec statistique observée
ax = axes[2]
df1, df2 = 2, len(df_glyc) - 3
x_f = np.linspace(0, 8, 400)
y_f = stats.f.pdf(x_f, df1, df2)
ax.plot(x_f, y_f, 'k-', lw=2, label=f'F({df1}, {df2})')
f_crit = stats.f.ppf(0.95, df1, df2)
mask_rejet = x_f > f_crit
ax.fill_between(x_f[mask_rejet], y_f[mask_rejet], alpha=0.4, color='crimson', label='Région de rejet (5%)')
ax.axvline(f_stat, color='navy', lw=2.5, linestyle='--', label=f'F obs. = {f_stat:.2f}')
ax.set_xlabel('F')
ax.set_ylabel('Densité')
ax.set_title('Distribution de Fisher sous H₀')
ax.legend(fontsize=9)

plt.tight_layout()
plt.savefig('_static/06_anova1.png', bbox_inches='tight')
plt.show()

ANOVA à deux facteurs

L’ANOVA à deux facteurs analyse simultanément deux variables catégorielles et leur éventuelle interaction.

Exemple : médicament × dose

# Plan factoriel 2×3 : 2 médicaments × 3 doses
rng2 = np.random.default_rng(123)
n_par_cellule = 15

medicaments = ['Méd. A', 'Méd. B']
doses = ['Faible', 'Moyenne', 'Forte']

data_rows = []
for med in medicaments:
    for dose in doses:
        # Effets principaux + interaction
        base = 80
        effet_med = 5 if med == 'Méd. B' else 0
        effet_dose = {'Faible': 0, 'Moyenne': 8, 'Forte': 14}[dose]
        # Interaction : Méd B est plus efficace à forte dose
        interaction = 6 if (med == 'Méd. B' and dose == 'Forte') else 0
        mu = base - effet_med - effet_dose - interaction
        vals = rng2.normal(mu, 8, n_par_cellule)
        for v in vals:
            data_rows.append({'medicament': med, 'dose': dose, 'reponse': v})

df_fact = pd.DataFrame(data_rows)

# ANOVA à deux facteurs avec pingouin
aov2 = pg.anova(data=df_fact, dv='reponse', between=['medicament', 'dose'], detailed=True)
print("ANOVA à deux facteurs :")
print(aov2[['Source', 'SS', 'DF', 'MS', 'F', 'p_unc', 'np2']].to_string(index=False))

ANOVA à deux facteurs :
           Source          SS  DF          MS         F        p_unc      np2
       medicament 1626.149100   1 1626.149100 30.469331 3.696024e-07 0.266179
             dose 4518.597961   2 2259.298981 42.332728 1.925892e-13 0.501973
medicament * dose  150.941926   2   75.470963  1.414108 2.488704e-01 0.032573
         Residual 4483.082558  84   53.370030       NaN          NaN      NaN

Show code cell source

Hide code cell source

fig, axes = plt.subplots(1, 2, figsize=(13, 4.5))

# Interaction plot
ax = axes[0]
means = df_fact.groupby(['medicament', 'dose'])['reponse'].mean().reset_index()
dose_order = ['Faible', 'Moyenne', 'Forte']
for med, color in zip(medicaments, ['#4C72B0', '#DD8452']):
    sub = means[means['medicament'] == med].copy()
    sub['dose'] = pd.Categorical(sub['dose'], categories=dose_order, ordered=True)
    sub = sub.sort_values('dose')
    ax.plot(sub['dose'], sub['reponse'], 'o-', color=color, lw=2, markersize=8, label=med)
ax.set_xlabel('Dose')
ax.set_ylabel('Réponse moyenne')
ax.set_title('Graphique d\'interaction\n(croisement = interaction significative)')
ax.legend()

# Boxplot groupé
ax = axes[1]
sns.boxplot(data=df_fact, x='dose', y='reponse', hue='medicament',
            order=dose_order, ax=ax, palette=['#4C72B0', '#DD8452'])
ax.set_xlabel('Dose')
ax.set_ylabel('Réponse')
ax.set_title('Distribution par dose et médicament')
ax.legend(title='Médicament')

plt.tight_layout()
plt.savefig('_static/06_anova2.png', bbox_inches='tight')
plt.show()

Test de Fisher : comparer deux variances

Le test F de Fisher teste \(H_0 : \sigma_1^2 = \sigma_2^2\) via la statistique \(F = s_1^2 / s_2^2\).

# Exemple : variabilité de deux procédés industriels
proc_A = rng.normal(100, 5, 30)
proc_B = rng.normal(100, 8, 30)

var_A = proc_A.var(ddof=1)
var_B = proc_B.var(ddof=1)
F_obs = var_A / var_B

# p-valeur bilatérale
df1, df2 = len(proc_A) - 1, len(proc_B) - 1
p_droite = stats.f.sf(F_obs, df1, df2)
p_gauche = stats.f.cdf(F_obs, df1, df2)
p_fisher = 2 * min(p_droite, p_gauche)

print(f"Variance procédé A : {var_A:.2f}")
print(f"Variance procédé B : {var_B:.2f}")
print(f"F = s²A / s²B = {F_obs:.3f}")
print(f"p-valeur (bilatéral) = {p_fisher:.4f}")
print()
print("Note : Levene est plus robuste que Fisher pour tester l'égalité des variances.")
lev_stat, lev_p = stats.levene(proc_A, proc_B)
print(f"Test de Levene : F = {lev_stat:.3f}, p = {lev_p:.4f}")

Variance procédé A : 27.22
Variance procédé B : 60.72
F = s²A / s²B = 0.448
p-valeur (bilatéral) = 0.0345

Note : Levene est plus robuste que Fisher pour tester l'égalité des variances.
Test de Levene : F = 2.848, p = 0.0968

Vérification des hypothèses

Normalité : test de Shapiro-Wilk

Show code cell source

Hide code cell source

fig, axes = plt.subplots(2, 4, figsize=(16, 7))

# Générer différentes distributions pour illustration
distributions = {
    'Normal': rng.normal(0, 1, 60),
    'Légèrement asymétrique': rng.exponential(1, 60) - 1,
    'Fortement asymétrique': rng.exponential(0.3, 60),
    'Bimodale': np.concatenate([rng.normal(-2, 0.8, 30), rng.normal(2, 0.8, 30)]),
}

for i, (name, data) in enumerate(distributions.items()):
    # Histogramme + courbe normale ajustée
    ax = axes[0, i]
    ax.hist(data, bins=15, density=True, alpha=0.6, color='steelblue', edgecolor='white')
    x_range = np.linspace(data.min(), data.max(), 200)
    ax.plot(x_range, stats.norm.pdf(x_range, data.mean(), data.std()), 'crimson', lw=2)
    stat, p = stats.shapiro(data)
    ax.set_title(f'{name}\nShapiro p = {p:.3f}', fontsize=9)
    ax.set_xlabel('Valeur')

    # QQ-plot
    ax = axes[1, i]
    (osm, osr), (slope, intercept, r) = stats.probplot(data, dist='norm')
    ax.plot(osm, osr, 'o', color='steelblue', alpha=0.6, markersize=4)
    ax.plot(osm, slope * np.array(osm) + intercept, 'crimson', lw=2)
    ax.set_title(f'QQ-plot', fontsize=9)
    ax.set_xlabel('Quantiles théoriques')
    ax.set_ylabel('Quantiles observés')

axes[0, 0].set_ylabel('Densité')
axes[1, 0].set_ylabel('Quantiles observés')
plt.suptitle('Vérification de normalité : histogrammes et QQ-plots', y=1.01, fontsize=12)
plt.tight_layout()
plt.savefig('_static/06_normalite.png', bbox_inches='tight')
plt.show()

Homoscédasticité : Levene et Bartlett

# Trois groupes, dont un avec variance différente
g1 = rng.normal(50, 5, 30)
g2 = rng.normal(52, 5, 30)
g3 = rng.normal(54, 15, 30)  # variance beaucoup plus grande

lev_stat, lev_p = stats.levene(g1, g2, g3)
bart_stat, bart_p = stats.bartlett(g1, g2, g3)

print("Test de Levene (robuste aux non-normalités) :")
print(f"  F = {lev_stat:.3f}, p = {lev_p:.4f}")
print()
print("Test de Bartlett (sensible aux non-normalités) :")
print(f"  χ² = {bart_stat:.3f}, p = {bart_p:.4f}")
print()
print(f"Variance G1: {g1.var(ddof=1):.1f}")
print(f"Variance G2: {g2.var(ddof=1):.1f}")
print(f"Variance G3: {g3.var(ddof=1):.1f}")

Test de Levene (robuste aux non-normalités) :
  F = 26.480, p = 0.0000

Test de Bartlett (sensible aux non-normalités) :
  χ² = 64.849, p = 0.0000

Variance G1: 14.3
Variance G2: 33.6
Variance G3: 279.1

Levene vs Bartlett

Le test de Levene est préféré en pratique : il est robuste aux écarts à la normalité, contrairement au test de Bartlett qui suppose que les groupes sont normalement distribués. Si Levene est significatif, envisagez le test de Welch (deux groupes) ou le test de Brown-Forsythe (plusieurs groupes).

Puissance statistique

La puissance \(1 - \beta\) est la probabilité de détecter un effet réel. Elle dépend de :

• La taille de l’effet (\(d\) de Cohen, \(f\) de Cohen, \(\eta^2\)…)

• La taille de l’échantillon \(n\)

• Le seuil \(\alpha\)

Show code cell source

Hide code cell source

# Courbes de puissance pour le test t (deux échantillons indépendants)
fig, axes = plt.subplots(1, 2, figsize=(13, 4.5))

# Puissance vs n pour différentes tailles d'effet
ax = axes[0]
n_range = np.arange(5, 150)
for d, label, color in zip([0.2, 0.5, 0.8, 1.2],
                            ['Petit (d=0.2)', 'Modéré (d=0.5)', 'Grand (d=0.8)', 'Très grand (d=1.2)'],
                            ['#e41a1c', '#ff7f00', '#4daf4a', '#377eb8']):
    power_vals = [smp.TTestIndPower().power(effect_size=d, nobs1=n, alpha=0.05, ratio=1)
                  for n in n_range]
    ax.plot(n_range, power_vals, color=color, lw=2, label=label)

ax.axhline(0.80, color='black', linestyle=':', lw=1.5, label='Puissance 80%')
ax.axhline(0.90, color='gray', linestyle=':', lw=1.5, label='Puissance 90%')
ax.set_xlabel('Taille d\'échantillon (par groupe)')
ax.set_ylabel('Puissance (1 - β)')
ax.set_title('Courbes de puissance — Test t indépendant')
ax.legend(fontsize=8)
ax.set_ylim(0, 1.05)

# Calcul du n nécessaire
ax = axes[1]
alphas = [0.01, 0.05, 0.10]
d_range = np.linspace(0.2, 1.5, 100)
for alpha, color in zip(alphas, ['#e41a1c', '#4daf4a', '#377eb8']):
    n_needed = [smp.TTestIndPower().solve_power(effect_size=d, alpha=alpha, power=0.80, ratio=1)
                for d in d_range]
    ax.plot(d_range, n_needed, color=color, lw=2, label=f'α = {alpha}')

ax.set_xlabel('Taille d\'effet (d de Cohen)')
ax.set_ylabel('N nécessaire par groupe')
ax.set_title('Taille d\'échantillon nécessaire (puissance = 80%)')
ax.legend()
ax.set_ylim(0, 400)
ax.axvline(0.5, color='gray', linestyle=':', lw=1)
ax.text(0.5, 380, 'Effet modéré', ha='center', fontsize=9, color='gray')

plt.tight_layout()
plt.savefig('_static/06_puissance.png', bbox_inches='tight')
plt.show()

# Calcul de puissance et de n nécessaire
print("=== Calcul de puissance ===")
analysis = smp.TTestIndPower()

# Puissance pour n donné
n_obs = 30
d_cible = 0.5
pow_calc = analysis.power(effect_size=d_cible, nobs1=n_obs, alpha=0.05, ratio=1)
print(f"Puissance avec n={n_obs} par groupe, d={d_cible} : {pow_calc:.1%}")

# N nécessaire pour puissance cible
pow_cible = 0.80
n_nec = analysis.solve_power(effect_size=d_cible, alpha=0.05, power=pow_cible, ratio=1)
print(f"\nN nécessaire pour 80% de puissance, d={d_cible} : {n_nec:.0f} par groupe")

# ANOVA : puissance
print("\n=== ANOVA — puissance ===")
# f de Cohen pour ANOVA (f = 0.10 petit, 0.25 modéré, 0.40 grand)
f_cohen = 0.25
k_groupes = 3
anova_power = smp.FTestAnovaPower()
n_anova = anova_power.solve_power(effect_size=f_cohen, alpha=0.05, power=0.80, k_groups=k_groupes)
print(f"N par groupe pour ANOVA ({k_groupes} groupes), f={f_cohen}, puissance=80% : {n_anova:.0f}")

=== Calcul de puissance ===
Puissance avec n=30 par groupe, d=0.5 : 47.8%

N nécessaire pour 80% de puissance, d=0.5 : 64 par groupe

=== ANOVA — puissance ===
N par groupe pour ANOVA (3 groupes), f=0.25, puissance=80% : 157

Résumé : choisir le bon test

Show code cell source

Hide code cell source

# Tableau récapitulatif
recap = pd.DataFrame({
    'Situation': [
        '1 groupe, μ connu',
        '2 groupes indépendants',
        '2 groupes indépendants (var ≠)',
        '2 groupes appariés',
        'k groupes indépendants',
        'k groupes, 2 facteurs',
        '2 variances',
    ],
    'Test': [
        't à 1 échantillon',
        't de Student',
        't de Welch (recommandé)',
        't apparié',
        'ANOVA à 1 facteur',
        'ANOVA à 2 facteurs',
        'Test F de Fisher',
    ],
    'Fonction scipy': [
        'ttest_1samp()',
        'ttest_ind(equal_var=True)',
        'ttest_ind(equal_var=False)',
        'ttest_rel()',
        'f_oneway()',
        'smf.ols() + anova_lm()',
        'f.ppf() / levene()',
    ],
    'Taille d\'effet': [
        'd de Cohen',
        'd de Cohen',
        'd de Cohen',
        'd de Cohen',
        'η² (eta-carré)',
        'η² partiel',
        'R² de variance',
    ]
})

print(recap.to_string(index=False))

                     Situation                    Test             Fonction scipy Taille d'effet
             1 groupe, μ connu       t à 1 échantillon              ttest_1samp()     d de Cohen
        2 groupes indépendants            t de Student  ttest_ind(equal_var=True)     d de Cohen
2 groupes indépendants (var ≠) t de Welch (recommandé) ttest_ind(equal_var=False)     d de Cohen
            2 groupes appariés               t apparié                ttest_rel()     d de Cohen
        k groupes indépendants       ANOVA à 1 facteur                 f_oneway() η² (eta-carré)
         k groupes, 2 facteurs      ANOVA à 2 facteurs     smf.ols() + anova_lm()     η² partiel
                   2 variances        Test F de Fisher         f.ppf() / levene() R² de variance

Bonnes pratiques des tests paramétriques

1. Vérifier les hypothèses avant le test : normalité (Shapiro-Wilk, QQ-plot), homoscédasticité (Levene).

2. Rapporter la taille d’effet (d de Cohen, η²) en plus de la p-valeur.

3. Calculer la puissance a priori pour déterminer la taille d’échantillon, et a posteriori pour interpréter les résultats non significatifs.

4. Utiliser Welch par défaut pour les comparaisons à deux groupes — il est robuste.

5. Ne pas confondre significativité statistique et importance pratique : un effet très faible peut être significatif avec un grand \(n\).

6. Fixer \(\alpha\) avant de voir les données.