Activité

Comparer deux proportions par permutation

La proportion d’accidents avec victimes diffère-t-elle entre semaine et fin de semaine ?

60-90 minutes, indicatives · Intermédiaire

À produire : Tableau des deux proportions, distribution nulle simulée, p-valeur Monte-Carlo et référence exacte, puis six réponses commentées.

Prérequis : Comprendre proportion et tableau croisé. Lire un histogramme et distinguer différence relative et point de pourcentage. Exécuter un script R ou utiliser l’outil interactif.

Préparation enseignante

Prévoir 60 à 90 minutes, durée indicative. Extraire la trousse entière, ouvrir Donnees-bleues.Rproj, puis exécuter installer-packages.R avant la séance avec Internet. Les rapports 2022 sont inclus; le script fonctionne ensuite hors ligne. L’utilisation en classe n’est pas documentée.

Sans programmation, suivre l’outil interactif. La lecture courte prépare l’interprétation d’une p-valeur.

Télécharger la trousse avec les données

121 Ko. Préparation : 03/10/2026 à 09:51 (heure du Québec).

Les fichiers de données sont inclus dans la trousse.

Conditions de la source · Provenance et empreinte du ZIP

Fichiers inclus

rapports_accident_2022_prepares.csv : 108 186 lignes, 6 colonnes.

Contexte : une question et une population

La question porte sur la présence de victimes parmi les accidents rapportés en semaine et la fin de semaine. Le CSV annuel publié par la SAAQ contient 108 186 rapports pour 2022. Aucun tirage ni retrait de lignes n’est effectué. « Avec victimes » regroupe les gravités « Léger » et « Mortel ou grave ». Les codes SEM et FDS désignent lundi au vendredi et samedi ou dimanche.

Le dénominateur de chaque proportion est le nombre de rapports du groupe. Sans nombre de trajets ou distance parcourue, ces données ne permettent pas d’estimer le risque d’avoir un accident. Les variables de jour et de gravité sont complètes dans cette version. Les variables et la méthode précisent les contrôles.

L’écart du fichier est descriptif. Pour apprendre un test, on ajoute un modèle nul : à marges fixes, les statuts de victime seraient échangeables entre tous les rapports. Sous ce modèle, il n’y a pas d’association entre groupe de jours et présence de victimes. L’alternative permet une différence dans les deux directions. Cette échangeabilité globale n’est pas établie pour les rapports administratifs; régions, météo, heures et regroupements peuvent intervenir.

Consignes

Exécuter datasets/rapports-accident/activite-permutation.R. Il utilise 5 000 permutations, graine 20261003.
Identifier les numérateurs et dénominateurs, puis calculer l’écart fin de semaine moins semaine en points de pourcentage.
Examiner le premier mélange. Vérifier les effectifs de groupes et le total des accidents avec victimes, avant et après.
Lire l’histogramme des écarts simulés sous le modèle nul. Compter les valeurs dont l’écart absolu atteint ou dépasse celui observé.
Calculer (k + 1)/(B + 1). Comparer les 500, 2 000 et 5 000 premiers mélanges et la référence exacte.
Répondre aux six exercices, puis consulter le corrigé. Les tableaux et figures sont enregistrés dans outputs/.

Résultats calculés

Groupe	Rapports	Avec victimes	Proportion (%)
Lundi au vendredi	82 168	16 740	20,37
Samedi ou dimanche	26 018	5 758	22,13

L’écart est de 1,76 points de pourcentage. Il est calculé avec les proportions non arrondies, puis arrondi pour l’affichage.

Deux proportions observées parmi les rapports de 2022, calculées séparément pour la semaine et la fin de semaine.

Histogramme des écarts simulés sous le modèle nul, avec les deux seuils du test bilatéral en orange.

Le trait orange plein donne l’écart observé. Le trait orange pointillé donne le seuil de même amplitude dans l’autre direction. Les barres comptent des permutations. Ici, 0 mélanges sur 5000 atteignent l’un des seuils. La p-valeur Monte-Carlo vaut 0,0002.

B	p Monte-Carlo	Résolution minimale	p de référence exacte
500	0,0020	0,0020	1,20e-09
2000	0,0005	0,0005	1,20e-09
5000	0,0002	0,0002	1,20e-09

La référence exacte vaut 1,20e-09 pour le même critère bilatéral. Elle est calculée par la loi hypergéométrique, sans simulation; « exacte » qualifie le calcul sous le modèle, pas sa pertinence pour ces données. Elle conserve les mêmes hypothèses. Le détail avancé figure dans la méthode.

Six exercices

Quels sont les numérateurs et dénominateurs des deux proportions ? Que mesure leur différence ?
Dans le premier mélange, quelles quantités restent identiques et quelles quantités peuvent changer ?
Pourquoi compter les valeurs positives et négatives dont l’écart absolu atteint l’écart observé ?
Aucun mélange aussi extrême a-t-il été obtenu ? Peut-on écrire p = 0 ? Que vaut la p-valeur avec B = 5 000 ?
Si B augmente, change-t-on le nombre de rapports ou l’écart observé ? Pourquoi la référence exacte est-elle beaucoup plus petite que la p-valeur simulée ?
Évaluer la phrase : « La fin de semaine cause plus d’accidents et la probabilité que cette conclusion soit fausse est la p-valeur. » Formuler une conclusion acceptable et nommer une hypothèse à discuter.

Corrigé commenté

Semaine : 16740 accidents avec victimes parmi 82168 rapports, soit 20,37 %. Fin de semaine : 5758 parmi 26018, soit 22,13 %. Leur différence de 1,76 points décrit ces rapports, sans mesurer un risque par trajet.
Les groupes restent de tailles 82 168 et 26 018, pour 108 186 rapports. Le total reste 22 498 accidents avec victimes. La répartition des statuts entre groupes peut changer, donc leurs proportions et leur écart aussi.
L’alternative est bilatérale, définie avant le calcul. On utilise la valeur absolue et on inclut les égalités. Le signe indique le sens de la différence observée; il ne transforme pas le test en test unilatéral.
Avec cette graine, 0 mélanges sont aussi extrêmes. Le calcul est (0 + 1)/(5 000 + 1), soit 0,0002. La formule empêche une p-valeur nulle. Zéro mélange extrême renseigne sur la résolution du calcul, sans donner une probabilité exactement nulle.
Les rapports, les proportions et leur différence restent fixes. B contrôle les simulations. Quand les événements extrêmes sont très rares sous le modèle, quelques milliers de mélanges ne les résolvent pas : la p-valeur simulée peut atteindre son plancher. La référence exacte somme directement leurs probabilités et n’a pas ce plancher.
Dans le fichier, la proportion avec victimes est supérieure de 1,76 points la fin de semaine. Cet écart est peu compatible avec le modèle nul de mélange global. Une généralisation suppose de justifier l’échangeabilité; les effets de contexte et les regroupements ne sont pas contrôlés ici. Les jours ne sont pas affectés aléatoirement et les données d’exposition manquent. La p-valeur n’est ni la probabilité que l’hypothèse nulle soit vraie ni celle que la conclusion soit fausse.

Limites à faire nommer

Pour décrire uniquement le CSV publié, le tableau suffit; on n’a pas besoin d’un test pour montrer que ses deux proportions diffèrent. Le test apporte un modèle de référence hypothétique. Une faible p-valeur ne valide pas ce modèle, ne mesure pas l’importance de l’écart et ne supprime pas les différences de contexte entre groupes.

Le rapport d’accident est l’unité. Un accident peut avoir plusieurs victimes; aucune somme du code de dénombrement des victimes n’est utilisée. Pour reproduire les résultats ci-dessus dans l’application, choisir 5 000 permutations et la graine 20261003, puis calculer.

Code complet et sources

Afficher le script R

# Comparer deux proportions par permutation. Aurélien Nicosia (2026).
# Code MIT; données SAAQ et contenus originaux CC BY 4.0.
# Rapports publiés pour 2022. Le modèle ne mesure pas un risque par trajet.
suppressPackageStartupMessages({
  library(readr)
  library(dplyr)
  library(ggplot2)
  library(scales)
})

accident_groups <- c(SEM = "Lundi au vendredi", FDS = "Samedi ou dimanche")
accident_gravity <- c("Dommages matériels seulement",
  "Dommages matériels inférieurs au seuil de rapportage", "Léger", "Mortel ou grave")

read_permutation_accidents <- function(path) {
  data <- read_csv(path, col_types = cols(.default = col_character()),
                   show_col_types = FALSE)
  stopifnot(all(c("annee", "jour_semaine_code", "gravite") %in% names(data)),
    nrow(data) == 108186L, !anyNA(data[c("annee", "jour_semaine_code", "gravite")]),
    all(data$annee == "2022"), all(data$jour_semaine_code %in% names(accident_groups)),
    all(data$gravite %in% accident_gravity))
  data |> transmute(group = jour_semaine_code,
    victim = gravite %in% c("Léger", "Mortel ou grave"))
}

permutation_settings <- function(B = 2000L, seed = 20261003L) {
  stopifnot(length(B) == 1L, is.finite(B), B == as.integer(B), B >= 100L, B <= 20000L,
    length(seed) == 1L, is.finite(seed), seed >= 0, seed <= 2147483646,
    seed == as.integer(seed))
  list(B = as.integer(B), seed = as.integer(seed))
}

permutation_counts <- function(data) {
  stopifnot(is.logical(data$victim), !anyNA(data),
    all(data$group %in% names(accident_groups)), all(names(accident_groups) %in% data$group))
  tibble(group = names(accident_groups), label = unname(accident_groups),
    n = vapply(names(accident_groups), function(g) sum(data$group == g), integer(1)),
    victims = vapply(names(accident_groups), function(g)
      sum(data$victim[data$group == g]), integer(1))) |>
    mutate(proportion = victims / n)
}

# Exact signifie ici : somme de la loi conditionnelle, sans erreur Monte-Carlo.
# Le critère bilatéral est |p_FDS - p_SEM|, avec les égalités incluses.
# Ce n’est pas la convention bilatérale par probabilités de fisher.test().
permutation_exact <- function(n_weekend, n_weekday, victims, observed_weekend) {
  total <- as.double(n_weekend) + n_weekday
  support <- seq.int(max(0, victims - n_weekday), min(victims, n_weekend))
  numerator <- support * total - victims * n_weekend
  observed <- observed_weekend * total - victims * n_weekend
  extreme <- abs(numerator) >= abs(observed)
  min(1, sum(dhyper(support[extreme], m = victims, n = total - victims, k = n_weekend)))
}

permutation_accidents <- function(data, settings) {
  counts <- permutation_counts(data)
  n <- as.double(nrow(data))
  nw <- counts$n[counts$group == "FDS"]
  ns <- n - nw
  total_victims <- sum(data$victim)
  observed_w <- counts$victims[counts$group == "FDS"]
  # Le numérateur entier évite une ambiguïté d’arrondi pour les ex aequo.
  observed_numerator <- observed_w * n - total_victims * nw
  difference <- function(x) 100 * (x * n - total_victims * nw) / (nw * ns)

  old_kind <- RNGkind()
  had_seed <- exists(".Random.seed", envir = .GlobalEnv, inherits = FALSE)
  if (had_seed) old_seed <- get(".Random.seed", envir = .GlobalEnv)
  on.exit({
    do.call(RNGkind, as.list(old_kind))
    if (had_seed) assign(".Random.seed", old_seed, envir = .GlobalEnv)
    else if (exists(".Random.seed", envir = .GlobalEnv, inherits = FALSE))
      rm(".Random.seed", envir = .GlobalEnv)
  })
  RNGkind("Mersenne-Twister", "Inversion", "Rejection")
  set.seed(settings$seed)
  # Premier mélange explicite : chaque valeur binaire est réaffectée sans remise.
  first <- sample(data$victim, size = n, replace = FALSE)
  first_w <- sum(first[data$group == "FDS"])
  # Les autres mélanges sont simulés par leur décompte hypergéométrique.
  # Même loi que la permutation de tous les labels, sans B tableaux de 108186 lignes.
  simulated_w <- c(first_w, rhyper(settings$B - 1L,
    m = total_victims, n = n - total_victims, k = nw))
  numerators <- simulated_w * n - total_victims * nw
  replicates <- tibble(repetition = seq_len(settings$B),
    victims_weekend = simulated_w, difference_pp = difference(simulated_w),
    extreme = abs(numerators) >= abs(observed_numerator))
  list(settings = settings, counts = counts, n = n, victims = total_victims,
    observed_pp = difference(observed_w), replicates = replicates,
    exact_p = permutation_exact(nw, ns, total_victims, observed_w),
    first = tibble(position = seq_len(n), group = data$group,
      observed_victim = data$victim, permuted_victim = first))
}

permutation_summary <- function(result) {
  counts <- result$counts
  extreme <- sum(result$replicates$extreme)
  tibble(n_weekday = counts$n[1L], victims_weekday = counts$victims[1L],
    proportion_weekday = counts$proportion[1L], n_weekend = counts$n[2L],
    victims_weekend = counts$victims[2L], proportion_weekend = counts$proportion[2L],
    difference_pp = result$observed_pp, extreme = extreme,
    B = nrow(result$replicates), p_mc = (extreme + 1) / (nrow(result$replicates) + 1),
    resolution = 1 / (nrow(result$replicates) + 1), p_exact = result$exact_p,
    seed = result$settings$seed)
}

permutation_format <- function(x, digits = 2L) {
  formatC(x, format = "f", digits = digits, big.mark = " ", decimal.mark = ",")
}
permutation_p_format <- function(x) {
  if (x < 0.0001) formatC(x, format = "e", digits = 2, decimal.mark = ",")
  else permutation_format(x, 4L)
}

permutation_plot <- function(result, observed = FALSE) {
  if (observed) {
    return(ggplot(result$counts, aes(factor(label, levels = unname(accident_groups)), proportion)) +
      geom_col(fill = "#185b83", width = 0.55) +
      geom_text(aes(label = paste0(permutation_format(100 * proportion), " %")), vjust = -0.5) +
      scale_x_discrete(labels = c("Lundi au\nvendredi", "Samedi ou\ndimanche")) +
      scale_y_continuous(labels = label_percent(decimal.mark = ","), limits = c(0, 0.3)) +
      labs(x = NULL, y = "Proportion d’accidents avec victimes",
        title = "Les proportions observées dans les rapports 2022",
        caption = "Dénominateur : les accidents rapportés de chaque groupe, sans nombre de trajets.") +
      theme_minimal(base_size = 13))
  }
  summary <- permutation_summary(result)
  ggplot(result$replicates, aes(difference_pp)) +
    geom_histogram(bins = 35, fill = "#185b83", colour = "white", linewidth = 0.3) +
    geom_vline(xintercept = 0, colour = "#536779", linetype = "dotted") +
    geom_vline(xintercept = result$observed_pp, colour = "#c36b25", linewidth = 1) +
    geom_vline(xintercept = -result$observed_pp, colour = "#c36b25", linetype = "dashed") +
    scale_x_continuous(labels = label_number(decimal.mark = ",")) +
    labs(x = "Écart FDS - semaine\n(points de pourcentage)", y = "Nombre de permutations",
      title = "Que produirait le mélange des statuts de victime ?",
      subtitle = "Trait orange plein : écart observé; pointillé orange : seuil opposé du test bilatéral.",
      caption = paste0(summary$extreme, " mélanges au moins aussi extrêmes sur ", summary$B,
        "; p Monte-Carlo = ", permutation_p_format(summary$p_mc), ". Modèle conditionnel d’échangeabilité.")) +
    theme_minimal(base_size = 13) + theme(panel.grid.minor = element_blank())
}

# ACTIVITÉ : exécution autonome depuis le projet de la trousse.
fichier <- "data/processed/rapports-accident/rapports_accident_2022_prepares.csv"
if (!file.exists(fichier)) stop("Extraire la trousse entière et ouvrir Donnees-bleues.Rproj.")
accidents <- read_permutation_accidents(fichier)
reglages <- permutation_settings(B = 5000L)
resultats <- permutation_accidents(accidents, reglages)
resume <- permutation_summary(resultats)
print(resultats$counts)
print(resume)

# Une permutation conserve les groupes et le nombre total de victimes.
premier_melange <- resultats$first
print(head(premier_melange, 12))
comptes_melanges <- permutation_counts(tibble(group = premier_melange$group,
  victim = premier_melange$permuted_victim))
print(comptes_melanges)

# Préfixes des mêmes simulations : comparer uniquement B.
stabilite <- bind_rows(lapply(c(500L, 2000L, 5000L), function(B) {
  prefixe <- resultats
  prefixe$settings$B <- B
  prefixe$replicates <- head(resultats$replicates, B)
  permutation_summary(prefixe)
}))
print(stabilite)
graphique_observe <- permutation_plot(resultats, observed = TRUE)
graphique_permutations <- permutation_plot(resultats)
print(graphique_observe)
print(graphique_permutations)

dir.create("outputs", showWarnings = FALSE)
write_csv(resume, "outputs/permutations-resume.csv")
write_csv(resultats$counts, "outputs/permutations-proportions.csv")
write_csv(stabilite, "outputs/permutations-stabilite.csv")
write_csv(premier_melange, "outputs/permutations-premier-melange.csv")
write_csv(resultats$replicates, "outputs/permutations-repetitions.csv")
ggsave("outputs/permutations-observe.png", graphique_observe,
  width = 9, height = 5, dpi = 160, bg = "white")
ggsave("outputs/permutations-distribution.png", graphique_permutations,
  width = 9, height = 5, dpi = 160, bg = "white")

SAAQ, Rapports d’accident, CSV 2022, acquis le 3 octobre 2026; attribution et licence. Phipson et Smyth (2010), Permutation P-values Should Never Be Zero, pour le calcul Monte-Carlo. R Core Team, loi hypergéométrique, consultée le 3 octobre 2026, pour la référence à marges fixes. Contenus originaux CC BY 4.0; code MIT.

Vérifier le travail

Les dénominateurs sont les 82168 rapports en semaine et les 26018 rapports de fin de semaine.
Le mélange conserve le nombre total de rapports et d’accidents avec victimes.
La p-valeur Monte-Carlo utilise (k + 1) / (B + 1), avec les égalités incluses.
La conclusion nomme l’échangeabilité et l’absence de données d’exposition.

Objectifs et adaptations

Objectifs

Définir le dénominateur d’une proportion et calculer un écart en points de pourcentage.
Expliquer les quantités conservées dans une permutation.
Calculer une p-valeur bilatérale et distinguer B du nombre de rapports.
Formuler une conclusion conditionnelle aux hypothèses sans causalité ni risque par trajet.

Adaptations

Suivre le parcours sans programmation dans l’application.
Comparer 500, 2000 et 5000 mélanges d’une même suite.
Vérifier la loi hypergéométrique sur un petit tableau et distinguer les conventions bilatérales.

Statut documenté : Calculs et script exécutés; utilisation en classe non documentée.

Contribution pédagogique : Aurélien Nicosia. Usage en cours non documenté.

Ajout: 3 octobre 2026
Mise à jour: 3 octobre 2026