Bonnes pratiques et reproductibilité

Bonnes pratiques et reproductibilité#

Ce dernier chapitre prend du recul sur l’ensemble du livre pour aborder les questions qui conditionnent la valeur réelle d’une analyse statistique : est-ce que je teste ce que je crois tester ? Mon résultat sera-t-il reproductible ? Mes conclusions sont-elles solides ou fragiles ? La statistique moderne traverse une crise profonde — la crise de la réplication — dont les causes sont largement liées à de mauvaises pratiques que ce chapitre cherche à identifier et corriger.

La crise de la réplication#

Des exemples marquants#

Psychologie sociale (2011-2015). L’étude de « priming » de Bargh et al. (1996), qui prétendait que lire des mots liés à la vieillesse ralentissait la marche, a échoué à se reproduire à grande échelle. Le projet Reproducibility Project (OSC, 2015) a tenté de reproduire 100 expériences de psychologie publiées : seulement 39% des effets ont été répliqués avec un résultat significatif.

Médecine. Ioannidis (2005) dans le célèbre article « Why Most Published Research Findings Are False » démontre que dans de nombreux domaines biomédicaux, avec des puissances typiques de 20-40%, des biais de publication et des comparaisons multiples, la majorité des p-valeurs significatives publiées sont des faux positifs.

Économie. Camerer et al. (2016) ont tenté de reproduire 18 expériences en économie publiées dans Science et Nature : 61% se sont répliquées.

Les causes structurelles#

_images/7fdad297e43e54c6948cd0d5420643ad3701af1adf80e37aaf9aeab897d68d89.png

La p-valeur : ce qu’elle est, ce qu’elle n’est pas#

Définition correcte#

La p-valeur est la probabilité d’observer une statistique de test aussi ou plus extrême que celle obtenue, si H₀ est vraie :

\[p = P(|T| \geq |t_{\text{obs}}| \mid H_0)\]

Ce que la p-valeur N’EST PAS#

Mythes courants sur la p-valeur

FAUX : « p < 0.05 signifie que mon résultat est vrai avec 95% de probabilité »
FAUX : « p = 0.04 est plus important que p = 0.06 »
FAUX : « p > 0.05 prouve que l’effet est nul »
FAUX : « La p-valeur mesure la taille de l’effet »
FAUX : « p < 0.05 signifie que le résultat est important pratiquement »

VRAI : La p-valeur mesure la compatibilité des données avec H₀. Une p-valeur faible dit « ces données seraient rares sous H₀ », rien de plus.

L’ASA Statement de 2016#

L’American Statistical Association a publié en 2016 un communiqué officiel listant six principes sur les p-valeurs :

Les p-valeurs peuvent indiquer l’incompatibilité des données avec un modèle statistique spécifié
Les p-valeurs ne mesurent pas la probabilité que H₀ soit vraie, ni la probabilité que les résultats soient dus au hasard
Les conclusions scientifiques et décisions ne devraient pas être basées uniquement sur le franchissement d’un seuil de p-valeur
L’inférence correcte requiert un reporting complet et transparent
Une p-valeur n’est pas une mesure de la taille de l’effet ou de l’importance d’un résultat
En elle-même, une p-valeur ne mesure pas l’évidence pour un modèle ou une hypothèse

Show code cell source

Hide code cell source

# Simulation : distribution des p-valeurs sous H0 (doit être uniforme)
# et sous H1 (doit être concentrée près de 0)
n_sim = 5000
n_obs_sim = 30

p_sous_H0 = []
p_sous_H1_petit = []   # Cohen's d = 0.2 (petit effet)
p_sous_H1_moyen = []   # Cohen's d = 0.5 (effet moyen)
p_sous_H1_grand = []   # Cohen's d = 0.8 (grand effet)

for _ in range(n_sim):
    a = rng.normal(0, 1, n_obs_sim)
    b0 = rng.normal(0, 1, n_obs_sim)
    b1 = rng.normal(0.2, 1, n_obs_sim)
    b2 = rng.normal(0.5, 1, n_obs_sim)
    b3 = rng.normal(0.8, 1, n_obs_sim)
    p_sous_H0.append(stats.ttest_ind(a, b0).pvalue)
    p_sous_H1_petit.append(stats.ttest_ind(a, b1).pvalue)
    p_sous_H1_moyen.append(stats.ttest_ind(a, b2).pvalue)
    p_sous_H1_grand.append(stats.ttest_ind(a, b3).pvalue)

fig, axes = plt.subplots(1, 4, figsize=(16, 4))
configs = [
    (p_sous_H0, 'Sous H₀\n(d = 0)', 'lightgray', '#555555'),
    (p_sous_H1_petit, 'Sous H₁\n(d = 0.2, petit)', '#aed6f1', 'steelblue'),
    (p_sous_H1_moyen, 'Sous H₁\n(d = 0.5, moyen)', '#f9e79f', '#d4a017'),
    (p_sous_H1_grand, 'Sous H₁\n(d = 0.8, grand)', '#fadbd8', '#e74c3c'),
]

for ax, (pvals, titre, couleur_fond, couleur_bord) in zip(axes, configs):
    ax.hist(pvals, bins=20, range=(0, 1), color=couleur_fond,
            edgecolor=couleur_bord, linewidth=0.5, density=True)
    ax.axhline(1.0, color='gray', linestyle='--', linewidth=1, label='Uniforme U[0,1]')
    ax.axvline(0.05, color='orangered', linewidth=1.5, linestyle='--', label='α = 0.05')
    puissance = np.mean(np.array(pvals) < 0.05)
    ax.set_title(f'{titre}\nPuissance = {puissance:.1%}')
    ax.set_xlabel('p-valeur')
    ax.set_ylabel('Densité')
    if 'd = 0' in titre:
        ax.text(0.5, 0.8, 'Uniforme\n(attendu)', ha='center', fontsize=9,
                transform=ax.transAxes, color='gray')

plt.suptitle(f'Distribution des p-valeurs — {n_sim} simulations, n={n_obs_sim} par groupe',
             fontsize=13, y=1.01)
plt.tight_layout()
plt.show()

_images/ba918586a3aef250fff0570799dc10ef9c22dbe4b52c02ff0ea707441b9bb400.png

Taille d’effet : pourquoi c’est plus informatif#

Les mesures de taille d’effet#

Cohen’s d (deux groupes, variable continue) : $$d = \frac{\mu_1 - \mu_2}{\sigma_{\text{poolé}}} \quad \text{avec } \sigma_{\text{poolé}} = \sqrt{\frac{(n_1-1)s_1^2 + (n_2-1)s_2^2}{n_1+n_2-2}}$$

Convention de Cohen : $d = 0.2$ (petit), $d = 0.5$ (moyen), $d = 0.8$ (grand).

r de Pearson (mesure d’association ou d’effet pour le test t) : $$r = \frac{t}{\sqrt{t^2 + df}}$$

η² (eta carré) et ω² (omega carré) pour l’ANOVA : $$\eta^2 = \frac{SS_{\text{effet}}}{SS_{\text{total}}}, \quad \omega^2 = \frac{SS_{\text{effet}} - df_{\text{effet}} \cdot MS_{\text{résidu}}}{SS_{\text{total}} + MS_{\text{résidu}}}$$

η² surestime la taille d’effet (biais positif) ; ω² est un estimateur sans biais, recommandé.

Show code cell source

Hide code cell source

def cohen_d(groupe1, groupe2):
    n1, n2 = len(groupe1), len(groupe2)
    s_pool = np.sqrt(((n1-1)*groupe1.std(ddof=1)**2 + (n2-1)*groupe2.std(ddof=1)**2) / (n1+n2-2))
    return (groupe1.mean() - groupe2.mean()) / s_pool

def r_depuis_t(t, df):
    return t / np.sqrt(t**2 + df)

def omega_carre(ss_effet, df_effet, ms_residuel, ss_total):
    return (ss_effet - df_effet * ms_residuel) / (ss_total + ms_residuel)

# Illustration : même effet (d=0.5), mais n très différents → p-valeurs très différentes
np.random.seed(2024)
effets_demo = [0.5, 0.5, 0.5, 0.5]
ns_demo = [10, 30, 100, 500]

print("Impact de n sur la p-valeur — même taille d'effet Cohen's d = 0.5\n")
print(f"{'n/groupe':>10} | {'Cohen\'s d':>10} | {'p-valeur':>10} | {'Signif.':>8}")
print("-" * 46)
for n_d in ns_demo:
    g1 = rng.normal(0, 1, n_d)
    g2 = rng.normal(0.5, 1, n_d)
    d  = cohen_d(g1, g2)
    t_s, p = stats.ttest_ind(g1, g2)
    print(f"{n_d:>10} | {d:>10.3f} | {p:>10.4f} | {'Oui ✓' if p < 0.05 else 'Non':>8}")

# Landscape des tailles d'effet
fig, axes = plt.subplots(1, 2, figsize=(13, 5))

# Visualisation de Cohen's d
ax = axes[0]
d_vals = [0.2, 0.5, 0.8, 1.2]
x_range = np.linspace(-4, 5, 500)
couleurs_d = ['steelblue', 'seagreen', 'orangered', 'purple']
labels_d  = ['Petit (d=0.2)', 'Moyen (d=0.5)', 'Grand (d=0.8)', 'Très grand (d=1.2)']

# Distribution contrôle
ax.fill_between(x_range, norm.pdf(x_range, 0, 1), alpha=0.3, color='gray', label='Contrôle')
ax.plot(x_range, norm.pdf(x_range, 0, 1), color='gray', linewidth=1.5)

for d, col, lab in zip(d_vals, couleurs_d, labels_d):
    ax.plot(x_range, norm.pdf(x_range, d, 1), color=col, linewidth=2, label=lab)
    # Overlap coefficient (OVL)
    ovl = 2 * norm.cdf(-abs(d)/2)
    ax.annotate(f'OVL={ovl:.0%}', xy=(d, norm.pdf(d, d, 1)), color=col,
                fontsize=8, ha='center', va='bottom')

ax.set_xlabel('Valeur')
ax.set_ylabel('Densité')
ax.set_title('Distributions avec différentes tailles d\'effet (Cohen\'s d)\npar rapport à un groupe contrôle σ=1')
ax.legend(fontsize=8)

# p-valeur vs taille d'effet pour différents n
ax = axes[1]
d_range = np.linspace(0.05, 1.5, 100)
for n_i, col_i in [(15, 'lightblue'), (30, 'steelblue'), (100, 'orangered'), (300, 'seagreen')]:
    # Puissance du test t bilatéral
    nc = d_range * np.sqrt(n_i / 2)  # paramètre de non-centralité
    t_crit = t_dist.ppf(0.975, df=2*n_i-2)
    puissance_vals = 1 - nct.cdf(t_crit, df=2*n_i-2, nc=nc) + \
                     nct.cdf(-t_crit, df=2*n_i-2, nc=nc)
    ax.plot(d_range, puissance_vals, color=col_i, linewidth=2, label=f'n={n_i}')

ax.axhline(0.80, color='gray', linestyle='--', linewidth=1, label='80% puissance')
ax.axhline(0.05, color='lightgray', linestyle=':', linewidth=1)
ax.axvline(0.2, color='gray', alpha=0.4, linewidth=1)
ax.axvline(0.5, color='gray', alpha=0.4, linewidth=1)
ax.axvline(0.8, color='gray', alpha=0.4, linewidth=1)
ax.text(0.2, 0.02, 'petit', ha='center', fontsize=8, color='gray')
ax.text(0.5, 0.02, 'moyen', ha='center', fontsize=8, color='gray')
ax.text(0.8, 0.02, 'grand', ha='center', fontsize=8, color='gray')
ax.set_xlabel('Cohen\'s d')
ax.set_ylabel('Puissance')
ax.set_title('Puissance du test t selon n et taille d\'effet (α=5%)')
ax.legend(fontsize=9)

plt.suptitle('Taille d\'effet : une information que la p-valeur seule ne donne pas', fontsize=13, y=1.01)
plt.tight_layout()
plt.show()

Impact de n sur la p-valeur — même taille d'effet Cohen's d = 0.5

  n/groupe |  Cohen's d |   p-valeur |  Signif.
----------------------------------------------
        10 |     -1.690 |     0.0014 |    Oui ✓
        30 |     -0.797 |     0.0031 |    Oui ✓
       100 |     -0.600 |     0.0000 |    Oui ✓
       500 |     -0.430 |     0.0000 |    Oui ✓

_images/0f7dad7e187370b529e500f690eb1c5138f0bc53b4e4606d47630b2f0866509e.png

Puissance a priori vs a posteriori#

Puissance a priori#

Calculée avant la collecte des données, elle détermine la taille d’échantillon nécessaire pour détecter un effet minimal d’intérêt (SESOI, Smallest Effect Size of Interest) avec une probabilité suffisante (généralement 80% ou 90%).

Puissance a posteriori : inutile et trompeuse#

La puissance a posteriori (ou observed power) consiste à calculer la puissance en utilisant l’effet observé dans les données. C’est une erreur conceptuelle grave :

Show code cell source

Hide code cell source

# Démonstration : la puissance a posteriori est une transformation monotone de la p-valeur
# Pour un test t à deux groupes, si p = 0.05, la puissance obs. = 50% toujours !

n_demo = 40
n_sim = 3000
p_vals_demo = []
pwr_obs_demo = []

for _ in range(n_sim):
    g1 = rng.normal(0, 1, n_demo)
    g2 = rng.normal(0.3, 1, n_demo)  # H1 vraie, d=0.3
    t_s, p = stats.ttest_ind(g1, g2)
    p_vals_demo.append(p)
    # Puissance observée = puissance à l'effet observé
    d_obs = cohen_d(g1, g2)
    nc_obs = abs(d_obs) * np.sqrt(n_demo / 2)
    df_t = 2*n_demo - 2
    t_crit_obs = t_dist.ppf(0.975, df=df_t)
    pwr_a_post = 1 - nct.cdf(t_crit_obs, df=df_t, nc=nc_obs) + \
                 nct.cdf(-t_crit_obs, df=df_t, nc=nc_obs)
    pwr_obs_demo.append(pwr_a_post)

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

ax = axes[0]
ax.scatter(p_vals_demo, pwr_obs_demo, s=4, alpha=0.2, color='steelblue')
ax.axhline(0.5, color='orangered', linestyle='--', linewidth=1.5,
           label='Puissance obs. = 50% ↔ p = α')
ax.axvline(0.05, color='orangered', linestyle='--', linewidth=1.5)
ax.set_xlabel('p-valeur')
ax.set_ylabel('Puissance a posteriori')
ax.set_title('Puissance a posteriori vs p-valeur\n(relation biunivoque : l\'une détermine l\'autre)')
ax.legend(fontsize=9)

ax = axes[1]
# Bonne pratique : IC plutôt que test binaire
effets_obs = [cohen_d(
    rng.normal(0, 1, n_demo),
    rng.normal(0.3, 1, n_demo)
) for _ in range(50)]

se_d_approx = np.sqrt(2/n_demo + np.array(effets_obs)**2 / (2*n_demo))
ic_low_d  = np.array(effets_obs) - 1.96 * se_d_approx
ic_high_d = np.array(effets_obs) + 1.96 * se_d_approx
couleurs_ic = ['steelblue' if ic_l > 0 else 'orangered' if ic_h < 0 else 'gray'
               for ic_l, ic_h in zip(ic_low_d, ic_high_d)]

for i, (d, il, ih, col) in enumerate(sorted(zip(effets_obs, ic_low_d, ic_high_d, couleurs_ic),
                                             key=lambda x: x[0])):
    ax.errorbar(d, i, xerr=[[d-il], [ih-d]], fmt='o', color=col, capsize=3, elinewidth=1, markersize=3)

ax.axvline(0, color='gray', linewidth=1, linestyle='--', label='H₀ : d = 0')
ax.axvline(0.3, color='seagreen', linewidth=1.5, linestyle='--', label='Vrai d = 0.3')
ax.set_xlabel('Cohen\'s d estimé (+ IC 95%)')
ax.set_ylabel('Simulation (triées par d)')
ax.set_title('50 réplications : IC sur d vs décision binaire\n(bleu = sig. positif, gris = non sig.)')
ax.legend(fontsize=8)

plt.suptitle('Puissance a posteriori : inutile. Intervalles de confiance : informatifs.', fontsize=13, y=1.01)
plt.tight_layout()
plt.show()

_images/10ff36b2e2a22ea9cf21a673bf953fa5ad44209a78f5be8370f06217f6301aff.png

P-hacking et ses conséquences#

Simulation du p-hacking#

Le p-hacking désigne l’ensemble des pratiques qui augmentent artificiellement la probabilité d’obtenir $p < 0.05$ : chercher plusieurs analyses, supprimer des données gênantes, ajouter des covariables jusqu’à la significativité, etc.

Show code cell source

Hide code cell source

def analyse_mauvaise_pratique(graine):
    """
    Mauvaise pratique : tenter 10 analyses différentes sur les mêmes données
    et ne rapporter que la plus significative.
    """
    rng_p = np.random.default_rng(graine)
    n = 50
    # Données sous H0 : aucun effet réel
    x = rng_p.normal(0, 1, n)
    y = rng_p.normal(0, 1, n)
    w = rng_p.normal(0, 1, n)  # covariable

    analyses = []
    # 1. Test t simple
    _, p = stats.ttest_ind(x[:25], x[25:])
    analyses.append(p)
    # 2. Test sur sous-échantillon
    _, p = stats.ttest_ind(x[:20], x[20:40])
    analyses.append(p)
    # 3. Corrélation x-y
    _, p = stats.pearsonr(x, y)
    analyses.append(p)
    # 4. Test sur valeurs extrêmes exclues
    mask = np.abs(x) < 1.5
    if mask.sum() > 10:
        _, p = stats.ttest_1samp(x[mask], 0)
        analyses.append(p)
    # 5-7. Corrélations partielles simulées
    for shift in [0.1, 0.2, -0.1]:
        _, p = stats.pearsonr(x + shift*w, y)
        analyses.append(p)
    # 8-10. Sous-groupes
    for q1, q2 in [(25, 75), (33, 67), (10, 90)]:
        mask_q = (x > np.percentile(x, q1)) & (x < np.percentile(x, q2))
        if mask_q.sum() > 5:
            _, p = stats.ttest_1samp(x[mask_q], 0)
            analyses.append(p)

    return min(analyses)  # on ne rapporte que la plus petite p-valeur

# Comparaison : analyse honnête vs p-hacking
N_SIMS = 3000
p_honnetes = [stats.ttest_ind(rng.normal(0,1,25), rng.normal(0,1,25)).pvalue
              for _ in range(N_SIMS)]
p_hackees = [analyse_mauvaise_pratique(i) for i in range(N_SIMS)]

faux_pos_honnete = np.mean(np.array(p_honnetes) < 0.05)
faux_pos_hack    = np.mean(np.array(p_hackees) < 0.05)

print(f"Taux de faux positifs (analyse honnête)  : {faux_pos_honnete:.1%}")
print(f"Taux de faux positifs (p-hacking × 10)  : {faux_pos_hack:.1%}")
print(f"Inflation du taux de faux positifs : ×{faux_pos_hack/faux_pos_honnete:.1f}")

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

bins_p = np.linspace(0, 1, 21)
ax = axes[0]
ax.hist(p_honnetes, bins=bins_p, color='steelblue', alpha=0.7, density=True,
        label=f'Honnête (FP = {faux_pos_honnete:.1%})')
ax.hist(p_hackees, bins=bins_p, color='orangered', alpha=0.5, density=True,
        label=f'P-hacking (FP = {faux_pos_hack:.1%})')
ax.axhline(1.0, color='gray', linestyle='--', linewidth=1, label='Uniforme U[0,1]')
ax.axvline(0.05, color='black', linewidth=1.5)
ax.set_xlabel('p-valeur'); ax.set_ylabel('Densité')
ax.set_title('Distribution des p-valeurs sous H₀\n(aucun effet réel)')
ax.legend(fontsize=9)

# Biais de publication : funnel plot asymétrique
ax = axes[1]
# Simulation de résultats publiés (biais vers les positifs)
n_etudes = 200
effets_vrais_nuls = rng.normal(0, 0.1, n_etudes)  # tous nuls
ses_vars = rng.uniform(0.1, 0.5, n_etudes)
t_stats  = effets_vrais_nuls / ses_vars
p_pub    = 2 * norm.cdf(-np.abs(t_stats))
# Publication biaisée : publier si p < 0.05 avec prob 90%, sinon 10%
publie = (p_pub < 0.05) | (rng.random(n_etudes) < 0.10)

ax.scatter(effets_vrais_nuls[publie], 1/ses_vars[publie],
           c=['orangered' if p < 0.05 else 'steelblue' for p in p_pub[publie]],
           alpha=0.6, s=30)
ax.set_xlabel('Taille d\'effet estimée')
ax.set_ylabel('Précision (1/SE)')
ax.set_title('Funnel plot asymétrique (biais de publication)\n(rouge = p<0.05, bleu = p≥0.05)')
ax.axvline(0, color='gray', linewidth=1, linestyle='--')

legend_elems = [
    mpatches.Patch(color='orangered', alpha=0.6, label='p < 0.05 (sur-représentés)'),
    mpatches.Patch(color='steelblue', alpha=0.6, label='p ≥ 0.05 (sous-représentés)'),
]
ax.legend(handles=legend_elems, fontsize=9)

plt.suptitle('P-hacking et biais de publication', fontsize=13, y=1.01)
plt.tight_layout()
plt.show()

Taux de faux positifs (analyse honnête)  : 4.7%
Taux de faux positifs (p-hacking × 10)  : 51.9%
Inflation du taux de faux positifs : ×11.1

_images/c2641b5082663556e44f42b061d0e2296460ada48da44e76ed1445dc9ba5fdbf.png

Rapport statistique complet#

Ce qu’un rapport doit contenir#

Les recommandations APA (7e édition), CONSORT (essais cliniques) et STROBE (études observationnelles) convergent vers les éléments suivants :

Statistiques descriptives : moyenne ± écart-type (ou médiane [IQR]) selon la distribution
Statistique de test : valeur de la statistique et degrés de liberté (ex : t(58) = 2.34)
P-valeur exacte : ne pas écrire « p < 0.05 » mais « p = 0.023 »
Taille d’effet avec IC 95% : Cohen’s d, r, η²…
Puissance a priori ou taille d’échantillon justifiée
Hypothèse exacte testée (bilatérale ? unilatérale ?)
Violations des hypothèses vérifiées et signalées

# Exemple de rapport complet avec pingouin (ou manual si non disponible)
rng_r = np.random.default_rng(2025)
groupe_trt = rng_r.normal(8.2, 1.5, 35)  # score qualité de vie après traitement
groupe_ctl = rng_r.normal(7.1, 1.6, 38)  # score qualité de vie groupe contrôle

# Test t
t_s, p_t = stats.ttest_ind(groupe_trt, groupe_ctl)
df_t = len(groupe_trt) + len(groupe_ctl) - 2
d_obs = cohen_d(groupe_trt, groupe_ctl)
r_obs = r_depuis_t(t_s, df_t)

# IC sur Cohen's d (approximation)
se_d = np.sqrt((len(groupe_trt)+len(groupe_ctl))/(len(groupe_trt)*len(groupe_ctl))
               + d_obs**2/(2*(len(groupe_trt)+len(groupe_ctl)-2)))
ic_d_low  = d_obs - 1.96 * se_d
ic_d_high = d_obs + 1.96 * se_d

print("═" * 60)
print("  RAPPORT STATISTIQUE COMPLET")
print("═" * 60)
print(f"\n  STATISTIQUES DESCRIPTIVES")
print(f"  Traitement : M = {groupe_trt.mean():.2f}, SD = {groupe_trt.std(ddof=1):.2f}, n = {len(groupe_trt)}")
print(f"  Contrôle   : M = {groupe_ctl.mean():.2f}, SD = {groupe_ctl.std(ddof=1):.2f}, n = {len(groupe_ctl)}")
print()
print(f"  HYPOTHÈSE TESTÉE")
print(f"  H₀ : μ_traitement = μ_contrôle (test bilatéral, α = 0.05)")
print()
print(f"  RÉSULTAT DU TEST")
print(f"  Test t de Student pour groupes indépendants")
print(f"  t({df_t}) = {t_s:.3f}, p = {p_t:.4f}")
print()
print(f"  TAILLE D'EFFET")
print(f"  Cohen's d = {d_obs:.3f} (IC 95% : [{ic_d_low:.3f}, {ic_d_high:.3f}])")
print(f"  r = {r_obs:.3f} ({'petit' if abs(r_obs) < 0.3 else 'moyen' if abs(r_obs) < 0.5 else 'grand'})")
print()
print(f"  VÉRIFICATION DES HYPOTHÈSES")
levene_s, levene_p = stats.levene(groupe_trt, groupe_ctl)
shapiro_t, shapiro_pt = stats.shapiro(groupe_trt)
shapiro_c, shapiro_pc = stats.shapiro(groupe_ctl)
print(f"  Égalité des variances (Levene) : F = {levene_s:.3f}, p = {levene_p:.4f}")
print(f"  Normalité traitement (Shapiro) : W = {shapiro_t:.3f}, p = {shapiro_pt:.4f}")
print(f"  Normalité contrôle   (Shapiro) : W = {shapiro_c:.3f}, p = {shapiro_pc:.4f}")
print()
print(f"  INTERPRÉTATION")
verdict = 'Rejet de H₀' if p_t < 0.05 else 'Non rejet de H₀'
print(f"  {verdict} : {'un' if p_t < 0.05 else 'aucun'} effet significatif du traitement sur la qualité de vie")
print(f"  (différence de {groupe_trt.mean()-groupe_ctl.mean():.2f} points, {d_obs:.2f} écart-types)")
print("═" * 60)

if PINGOUIN_DISPO:
    print("\n  Résultat pingouin :")
    df_pg = pd.DataFrame({
        'score': np.concatenate([groupe_trt, groupe_ctl]),
        'groupe': ['traitement']*len(groupe_trt) + ['contrôle']*len(groupe_ctl)
    })
    result_pg = pg.ttest(groupe_trt, groupe_ctl)
    print(result_pg.to_string())

════════════════════════════════════════════════════════════
  RAPPORT STATISTIQUE COMPLET
════════════════════════════════════════════════════════════

  STATISTIQUES DESCRIPTIVES
  Traitement : M = 8.12, SD = 1.64, n = 35
  Contrôle   : M = 6.96, SD = 1.79, n = 38

  HYPOTHÈSE TESTÉE
  H₀ : μ_traitement = μ_contrôle (test bilatéral, α = 0.05)

  RÉSULTAT DU TEST
  Test t de Student pour groupes indépendants
  t(71) = 2.879, p = 0.0053

  TAILLE D'EFFET
  Cohen's d = 0.674 (IC 95% : [0.202, 1.147])
  r = 0.323 (moyen)

  VÉRIFICATION DES HYPOTHÈSES
  Égalité des variances (Levene) : F = 0.407, p = 0.5254
  Normalité traitement (Shapiro) : W = 0.970, p = 0.4484
  Normalité contrôle   (Shapiro) : W = 0.980, p = 0.7302

  INTERPRÉTATION
  Rejet de H₀ : un effet significatif du traitement sur la qualité de vie
  (différence de 1.16 points, 0.67 écart-types)
════════════════════════════════════════════════════════════

  Résultat pingouin :
               T        dof alternative     p_val          CI95   cohen_d     power   BF10
T_test  2.888947  70.999602   two-sided  0.005122  [0.36, 1.96]  0.674414  0.810452  7.835

Reproductibilité computationnelle#

Bonnes pratiques Python#

# Rapport d'environnement reproductible
import sys
import importlib

def rapport_environnement():
    """Génère un rapport des versions des packages utilisés."""
    packages = ['numpy', 'scipy', 'pandas', 'matplotlib', 'seaborn',
                'sklearn', 'statsmodels']
    print("ENVIRONNEMENT D'EXÉCUTION")
    print(f"  Python : {sys.version.split()[0]}")
    for pkg in packages:
        try:
            mod = importlib.import_module(pkg)
            version = getattr(mod, '__version__', '?')
            print(f"  {pkg:<15} : {version}")
        except ImportError:
            print(f"  {pkg:<15} : non installé")

rapport_environnement()

ENVIRONNEMENT D'EXÉCUTION
  Python : 3.13.5
  numpy           : 2.3.5
  scipy           : 1.17.0
  pandas          : 3.0.0
  matplotlib      : 3.10.8
  seaborn         : 0.13.2
  sklearn         : 1.8.0
  statsmodels     : 0.14.6

Pré-enregistrement#

Le pré-enregistrement consiste à déposer publiquement, avant la collecte des données, le protocole complet de l’étude : hypothèses, plan d’analyse, taille d’échantillon, méthodes statistiques. Cela sépare :

Phase confirmatoire : tester les hypothèses pré-enregistrées → les p-valeurs sont interprétables
Phase exploratoire : générer de nouvelles hypothèses → les p-valeurs sont exploratoires

Plateformes : Open Science Framework (osf.io), AsPredicted (aspredicted.org).

Principes FAIR#

Les données et le code doivent être :

Findable (trouvables) : identifiants permanents (DOI), métadonnées
Accessible : accessible en ligne, format ouvert
Interoperable (interopérable) : formats standards (CSV, JSON)
Reusable : licence claire, documentation

Checklist finale#

Show code cell source

Hide code cell source

checklist_data = {
    "Design": [
        ("Hypothèse formulée a priori (pré-enregistrée si possible)", True),
        ("Taille d'échantillon justifiée par une analyse de puissance", True),
        ("Niveau α fixé avant l'analyse (généralement 0.05)", True),
        ("Métrique principale unique définie a priori", True),
    ],
    "Collecte": [
        ("Randomisation vérifiée (si applicable)", True),
        ("Protocole suivi sans déviation non documentée", True),
        ("Données manquantes documentées et mécanisme discuté", True),
    ],
    "Analyse": [
        ("Hypothèses du test vérifiées (normalité, homoscédasticité...)", True),
        ("Méthode d'analyse conforme au pré-enregistrement", True),
        ("Comparaisons multiples corrigées (Bonferroni, FDR...)", True),
        ("P-valeur exacte rapportée (pas seulement 'p < 0.05')", True),
        ("Taille d'effet + IC 95% calculés et rapportés", True),
        ("Analyses exploratoires clairement identifiées comme telles", True),
    ],
    "Rapport": [
        ("Statistiques descriptives complètes par groupe", True),
        ("Statistique de test + degrés de liberté rapportés", True),
        ("Code et données archivés (OSF, GitHub, Zenodo...)", True),
        ("Version des logiciels et packages documentée", True),
        ("Résultats négatifs rapportés également", True),
    ],
}

fig = plt.figure(figsize=(14, 9))
ax = fig.add_subplot(111)
ax.axis('off')

y_pos = 0.97
couleurs_cat = ['#2980b9', '#27ae60', '#8e44ad', '#e67e22']
x_col2 = 0.52

categs = list(checklist_data.items())
mid = len(categs) // 2
colonne1 = categs[:mid]
colonne2 = categs[mid:]

for col_idx, categorie_items in enumerate([colonne1, colonne2]):
    x0 = 0.01 if col_idx == 0 else x_col2
    y_c = 0.97
    for cat_i, (categorie, items) in enumerate(categorie_items):
        col = couleurs_cat[cat_i + col_idx * mid] if (cat_i + col_idx * mid) < len(couleurs_cat) else '#555555'
        ax.text(x0, y_c, f'■ {categorie}', fontsize=11, fontweight='bold', color=col,
                transform=ax.transAxes, va='top')
        y_c -= 0.055
        for (item, ok) in items:
            symbole = '☐'
            ax.text(x0 + 0.02, y_c, f'{symbole}  {item}', fontsize=9,
                    color='#333333', transform=ax.transAxes, va='top')
            y_c -= 0.048
        y_c -= 0.02

ax.set_title('Checklist : avant de conclure "effet statistiquement significatif"',
             fontsize=13, pad=15)
plt.tight_layout()
plt.show()

_images/d0768a9e0e427bd42fb6f7978c2a6b40a4f9068cabbde07be408cf60dfec7df9.png

Récapitulatif visuel du livre#

fig, ax = plt.subplots(figsize=(14, 7))
ax.axis('off')

etapes_livre = [
    ("Partie I\nExplorer", ["Statistiques\ndescriptives", "Distributions\nen pratique",
                             "Corrélation", "Visualisation"], '#2980b9'),
    ("Partie II\nInférence", ["Estimation\n& IC", "Tests\nparamétriques",
                              "Tests non-\nparamétriques", "Comparaisons\nmultiples"], '#27ae60'),
    ("Partie III\nModélisation", ["Régression\nlinéaire", "Régression\nmultiple",
                                   "Régression\nlogistique", "GLM"], '#8e44ad'),
    ("Partie IV\nAvancé", ["Bayésien", "MCMC", "Séries\ntemp.", "Multivarié"], '#e67e22'),
    ("Partie V\nPratique", ["A/B testing", "Données\nmanquantes",
                             "Monte Carlo", "Bonnes\npratiques"], '#e74c3c'),
]

x_start = 0.01
y_main = 0.75
y_sub  = 0.35

for i, (titre, sous_titres, col) in enumerate(etapes_livre):
    x = x_start + i * 0.20
    # Boîte principale
    fancy = mpatches.FancyBboxPatch((x, y_main - 0.08), 0.18, 0.22,
                                     boxstyle='round,pad=0.02',
                                     facecolor=col, alpha=0.2, edgecolor=col, linewidth=2,
                                     transform=ax.transAxes)
    ax.add_patch(fancy)
    ax.text(x + 0.09, y_main + 0.08, titre, ha='center', va='center',
            fontsize=10, fontweight='bold', color=col, transform=ax.transAxes)

    # Flèche vers sous-éléments
    ax.annotate('', xy=(x + 0.09, y_sub + 0.15), xytext=(x + 0.09, y_main - 0.08),
                xycoords='axes fraction', textcoords='axes fraction',
                arrowprops=dict(arrowstyle='->', color=col, lw=1.5))

    # Sous-éléments
    for j, st in enumerate(sous_titres):
        xj = x + (j % 2) * 0.09
        yj = y_sub + 0.10 - (j // 2) * 0.12
        ax.text(xj + 0.04, yj, st, ha='center', va='center',
                fontsize=8, color=col, transform=ax.transAxes,
                bbox=dict(boxstyle='round,pad=0.2', facecolor='white',
                          edgecolor=col, alpha=0.8))

# Flèche de progression
ax.annotate('', xy=(0.99, y_main + 0.03), xytext=(0.01, y_main + 0.03),
            xycoords='axes fraction', textcoords='axes fraction',
            arrowprops=dict(arrowstyle='->', color='gray', lw=2))

ax.text(0.5, 0.05, 'De l\'exploration à la reproductibilité : un parcours complet en statistiques appliquées',
        ha='center', va='center', fontsize=11, style='italic', color='#555555',
        transform=ax.transAxes)

ax.set_title('Statistiques appliquées — Lôc Cosnier\nRécapitulatif du livre',
             fontsize=14, pad=10)
plt.tight_layout()
plt.show()

_images/a1f81c8b75a2040cee242f1cf0954e8ede58794fc71824602f402cd3900da366.png

Conclusion#

La statistique n’est pas une boîte noire qui transforme des données en vérités. C’est un langage de quantification de l’incertitude qui, utilisé correctement, permet de tirer des conclusions robustes et reproductibles depuis des observations imparfaites.

Les meilleures pratiques ne sont pas un fardeau bureaucratique — elles sont la condition nécessaire pour que les résultats soient utiles au-delà du bureau de celui qui les a produits. Pré-enregistrer, rapporter les tailles d’effet, archiver le code et les données, distinguer exploration et confirmation : ces habitudes transforment une analyse anecdotique en contribution scientifique durable.

Les cinq leçons de ce livre

Décrire avant de tester : les statistiques descriptives et les visualisations révèlent ce que les tests ne peuvent pas.
La p-valeur seule ne suffit pas : toujours accompagner d’un intervalle de confiance et d’une taille d’effet.
La puissance est une décision de design : se poser la question avant, pas après.
Exploration ≠ Confirmation : les hypothèses nées des données ne peuvent pas être testées sur les mêmes données.
La reproductibilité est une responsabilité collective : code, données, protocoles — tout documenter, tout partager.

Bonnes pratiques et reproductibilité

Contenu

Bonnes pratiques et reproductibilité#

La crise de la réplication#

Des exemples marquants#

Les causes structurelles#

La p-valeur : ce qu’elle est, ce qu’elle n’est pas#

Définition correcte#

Ce que la p-valeur N’EST PAS#

L’ASA Statement de 2016#

Taille d’effet : pourquoi c’est plus informatif#

Les mesures de taille d’effet#

Puissance a priori vs a posteriori#

Puissance a priori#

Puissance a posteriori : inutile et trompeuse#

P-hacking et ses conséquences#

Simulation du p-hacking#

Rapport statistique complet#

Ce qu’un rapport doit contenir#

Reproductibilité computationnelle#

Bonnes pratiques Python#

Pré-enregistrement#

Principes FAIR#

Checklist finale#

Récapitulatif visuel du livre#

Conclusion#