required_pkgs <- c("tidyverse", "glmmTMB", "marginaleffects", "patchwork", "GGally", "distributional", "tinytable", "performance", "TMB", "Matrix")
# Reinstalling Matrix/TMB/glmmTMB together (rather than one at a time) avoids the
# most common cause of glmmTMB load failures: an ABI mismatch between glmmTMB's
# compiled code and the Matrix/TMB versions it was linked against.
missing_pkgs <- required_pkgs[!vapply(required_pkgs, requireNamespace, logical(1), quietly = TRUE)]
if (length(missing_pkgs) > 0) {
install.packages(missing_pkgs)
}
# Verify glmmTMB actually loads and can fit a trivial model - this is the step
# that catches ABI mismatches that "successful" installation can still leave behind.
glmmtmb_ok <- tryCatch({
library(glmmTMB)
test_fit <- glmmTMB(mpg ~ wt, data = mtcars)
TRUE
}, error = function(e) {
message("glmmTMB failed to load/fit: ", conditionMessage(e))
FALSE
})
if (!glmmtmb_ok) {
message("Reinstalling TMB and glmmTMB from source to resolve a likely ABI mismatch...")
install.packages(c("TMB", "glmmTMB"), type = "source")
# Re-check after the source reinstall
glmmtmb_ok <- tryCatch({
library(glmmTMB)
glmmTMB(mpg ~ wt, data = mtcars)
TRUE
}, error = function(e) {
message("Still failing: ", conditionMessage(e))
message("If this persists, check that you have a working C++ compiler installed ",
"(Rtools on Windows, Xcode command line tools on macOS), then restart R and re-run this chunk.")
FALSE
})
}
stopifnot("glmmTMB did not install/load successfully - see messages above" = glmmtmb_ok)À quoi sert ce document
Ceci est un court document préparatoire — veuillez le parcourir avant la séance. Il sert à deux choses seulement :
- Préparer votre installation de
Ravec toutes les librairies dont nous aurons besoin, afin que le temps de l’atelier soit consacré aux interactions, et non au dépannage d’installations. - Présenter les conventions de code et les fonctions que vous verrez revenir tout au long de l’atelier, pour qu’elles vous soient familières plutôt que nouvelles.
Vous n’avez pas besoin de tout maîtriser ici. L’objectif de l’atelier est de vous apprendre à penser, analyser et rapporter les interactions entre prédicteurs dans un modèle de régression — c’est la compétence que vous venez développer. Tout ce qui suit n’est qu’un échafaudage pour vous y amener plus rapidement.
Une remarque sur les données simulées
Chaque exemple de cet atelier utilise des données simulées : nous inventons une « vraie » relation (une ordonnée à l’origine, des pentes, parfois une interaction), nous générons des données fictives à partir de celle-ci avec des fonctions comme rnorm() ou rbinom(), puis nous ajustons un modèle pour voir s’il retrouve ce que nous y avons mis.
Nous faisons cela délibérément — non pas parce que de vrais jeux de données écologiques ne sont pas disponibles, mais parce que la simulation nous permet de vérifier notre travail : nous connaissons la vérité terrain, donc nous pouvons confirmer que le modèle et les résultats de marginaleffects se comportent comme prévu, avant de leur faire confiance sur de vraies données où la vérité est inconnue.
Vous n’êtes pas tenus de comprendre chaque ligne du code de simulation. Survolez-le, saisissez l’idée générale (« d’accord, on invente des poids de mouflons fictifs qui dépendent du sexe et de la profondeur de neige »), et passez à la suite. Ce qui compte pour les objectifs d’apprentissage de l’atelier, c’est ce qui vient après la simulation : comment ajuster le modèle, comment en extraire la bonne quantité, et comment rapporter cette quantité correctement. Si vous arrivez à suivre la logique d’une simulation sans retracer chaque appel de fonction, vous êtes en bonne posture.
Avant de commencer : installer les librairies
Exécutez le bloc ci-dessous une seule fois pour installer tout ce dont vous avez besoin. Il est écrit pour éviter les échecs d’installation les plus fréquents de glmmTMB :
- Incompatibilité binaire vs source :
glmmTMBest lié àTMBetMatrix, qui doivent être compilés avec des versions ABI compatibles. Installer les trois ensemble à partir du même instantané de dépôt (plutôt que de mettre à jour un seul à la fois) évite les erreurs de chargement du type « package was installed with a different version ». - Objets compilés périmés après une mise à jour de R ou de
Matrix: si vous avez récemment mis à jour R et que vous voyez une erreur du typeError: package or namespace load failed for 'glmmTMB' ... unable to load shared object, réinstallerTMBetglmmTMBà partir des sources avec--no-multiarch(viainstall.packages(type = "source")) règle généralement le problème. - Chaîne de compilation manquante : les installations à partir des sources de
TMB/glmmTMBnécessitent un compilateur C++ fonctionnel. Sous Windows, cela signifie Rtools; sous macOS, les outils en ligne de commande de Xcode. Le test ci-dessous vous indique si la compilation est indisponible avant que vous ne rencontriez une erreur cryptique.
Si le bloc signale un échec et que les solutions suggérées ne le résolvent pas, écrivez-moi avant l’atelier (pas le matin même) afin que nous ayons le temps de régler la situation.
Préparation de la session
Exécutez ceci au début de chaque session (une fois que le bloc d’installation ci-dessus a réussi une première fois). Nous réutiliserons exactement ce bloc en début d’atelier.
library(tidyverse) # Data manipulation and plotting
library(glmmTMB) # Generalized linear (mixed) models
library(marginaleffects) # Compute marginal/conditional effects
library(patchwork) # Combine multiple plots
library(GGally) # Pairwise plots
library(distributional) # Simulate mixture/gamma/normal distributions (Part 3)
library(tinytable) # Format results tables (Part 3)
library(performance) # Compute R2/pseudo-R2 for model fit (Part 3)
theme_set(theme_minimal(base_size = 14))
set.seed(4127)Vérification rapide : est-ce que glmmTMB fonctionne bien ?
Le bloc d’installation ci-dessus ajuste déjà un modèle trivial en coulisses pour vérifier que glmmTMB se charge correctement, mais il le fait silencieusement. Exécutez le bloc ci-dessous pour le voir fonctionner de bout en bout sur un petit jeu de données intégré (mtcars) — un résumé de modèle avec coefficients, erreurs-types et valeurs-p devrait s’afficher sous le code.
test_mod <- glmmTMB(mpg ~ wt, data = mtcars, family = gaussian())
summary(test_mod) Family: gaussian ( identity )
Formula: mpg ~ wt
Data: mtcars
AIC BIC logLik -2*log(L) df.resid
166.0 170.4 -80.0 160.0 29
Dispersion estimate for gaussian family (sigma^2): 8.7
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 37.2851 1.8180 20.509 <2e-16 ***
wt -5.3445 0.5413 -9.873 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Si cela produit une erreur, cela signifie qu’un aspect de l’installation n’a pas fonctionné — réexécutez le bloc d’installation ci-dessus, redémarrez R, et réessayez avant l’atelier. Si l’erreur persiste, écrivez-moi avec le message exact.
Conventions de code et fonctions utilisées tout au long
Voici un aide-mémoire rapide des conventions et des fonctions incontournables utilisées à répétition pendant l’atelier. Considérez ceci comme un glossaire à survoler maintenant et à consulter plus tard — pas quelque chose à mémoriser.
Le pipe natif, |>. Partout où vous voyez data |> mutate(...) |> ggplot(...), lisez de gauche à droite comme « prendre data, et ensuite appliquer mutate(), et ensuite transmettre le résultat à ggplot() ». C’est équivalent à imbriquer les appels (ggplot(mutate(data, ...))), mais plus facile à lire et à modifier séquentiellement. Nous utilisons le pipe natif de R (|>), pas l’ancien pipe de magrittr (%>%) — pour nos besoins ici, ils se comportent de la même façon.
Simuler des données avec tibble() et des fonctions de tirage aléatoire. Chaque étude de cas commence par simuler des données plutôt que de charger un jeu de données réel (voir « Une remarque sur les données simulées » ci-dessus). Nous choisissons d’abord des valeurs de paramètres « vraies » (une ordonnée à l’origine, une ou plusieurs pentes), puis nous générons des observations à partir d’un processus aléatoire connu :
rnorm(n, mean, sd)— tire d’une distribution Normale (variables réponses continues, p. ex. poids, longueur de pétale)rbinom(n, size = 1, prob)— tire d’une distribution de Bernoulli/Binomiale (variables réponses de présence/absence)rpois(n, lambda)— tire d’une distribution de Poisson (variables réponses de comptage, p. ex. abondance)runif(n, min, max)— tire d’une distribution Uniforme (souvent utilisée pour simuler un prédicteur, pas une variable réponse)
Simuler de cette façon signifie que nous connaissons la vérité terrain (l’ordonnée à l’origine et les pentes exactes utilisées pour générer les données), ce qui nous permet de vérifier directement si un modèle ajusté la retrouve — une habitude qui vaut la peine d’être conservée dans votre propre travail chaque fois que vous voulez valider une approche de modélisation avant de l’appliquer à de vraies données.
Ajuster des modèles avec glmmTMB(). Notre fonction de référence pour l’ajustement de modèles tout au long de l’atelier est glmmTMB(reponse ~ predicteurs, data = ..., family = ...). Points clés à retenir :
y ~ x1 * x2se développe eny ~ x1 + x2 + x1:x2— les effets principaux des deux prédicteurs et leur interaction, le tout en un raccourcifamily = gaussian()pour les variables réponses continues,binomial(link = "logit")pour la présence/absence,poisson(link = "log")pour les comptagesdispformula = ~ xpermet à des prédicteurs d’affecter la variance de la variable réponse, pas seulement sa moyenne
La librairie marginaleffects. C’est l’outil central pour l’objectif principal de l’atelier — extraire et rapporter la bonne quantité à partir d’un modèle comportant une interaction. Les coefficients bruts d’un modèle (issus de summary()) ne sont souvent pas la quantité que vous voulez réellement rapporter — particulièrement en présence d’interactions ou de fonctions de lien non linéaires. marginaleffects convertit les coefficients en réponses directes et interprétables :
predictions()— la valeur prédite par le modèle pour la variable réponse à des valeurs précises des covariables (une prédiction conditionnelle)avg_predictions()— la même chose, mais moyennée sur une distribution de valeurs des covariables (une prédiction marginale)slopes()— l’effet prédit par le modèle (la dérivée/pente) d’un prédicteur sur la variable réponse, à des valeurs précises des covariablesavg_slopes()— la même chose, moyennée sur une distribution de valeurs des covariablesdatagrid(...)— construit la grille de covariables à laquellepredictions()/slopes()sont évaluées; les arguments peuvent être des valeurs fixes, des fonctions commemeanoufivenum(le résumé en cinq nombres), ou des vecteurs de valeursby = "group"— plutôt que de regrouper (moyenner) sur une variable, garder les résultats répartis selon ses niveauxhypotheses(hypothesis = "b1 - b2 = 0")— teste des combinaisons linéaires d’estimations, p. ex. « la pente du groupe A est-elle différente de celle du groupe B ? »
set.seed(). Appelé une fois près du début de chaque script (ou avant un exercice pour lequel vous voulez pouvoir comparer vos réponses) afin que la simulation « aléatoire » soit réellement reproductible — exécuter le même code deux fois donne des nombres identiques. Si vos résultats diffèrent de ceux d’un collègue malgré un code identique, vérifiez si vous avez tous les deux utilisé la même graine (seed) au même endroit dans le script.
La librairie distributional. Utilisée dans la Partie 3 (fonctions de lien) pour simuler des prédicteurs, pas seulement des variables réponses. Plutôt que de tirer directement avec des fonctions de base de R comme rnorm()/rgamma(), nous décrivons d’abord une distribution comme un objet, puis nous en générons des valeurs :
dist_normal(mu, sd),dist_gamma(shape, rate),dist_beta(shape1, shape2),dist_poisson(lambda)— chacune construit un objet de distribution avec ces paramètres, sans tirer de nombres pour l’instantdist_mixture(dist_1, dist_2, weights = c(w1, w2))— combine deux ou plusieurs objets de distribution en un mélange (p. ex. une variable « couvert forestier » bimodale construite à partir de deux distributions bêta, l’une asymétrique vers le bas et l’autre vers le haut)generate(dist_object, n)— tirenvaleurs aléatoires d’un objet de distribution (l’équivalent d’appeler directementrnorm()/rgamma()/etc., mais qui fonctionne de façon uniforme pour n’importe quelle distribution construite avec les fonctions ci-dessus)
Pourquoi s’embêter avec ceci plutôt que d’appeler directement rbeta()/rgamma() ? Principalement par commodité lorsque la distribution d’un prédicteur est plus complexe qu’une seule famille nommée — dist_mixture() en particulier facilite la construction d’un prédicteur bimodal ou multimodal (comme des parcelles regroupées à faible ou à fort couvert forestier) en combinant deux distributions standards, ce qui est malaisé à faire avec les seules fonctions r*() de base de R. Comme pour le reste du code de simulation, vous n’avez pas besoin de pouvoir en construire une à partir de zéro — reconnaissez simplement dist_*() comme « définir une distribution » et generate() comme « en tirer des nombres ».
Quelques autres fonctions que vous verrez utilisées pour le rapport de résultats, et non pour la simulation :
plogis(x)/qlogis(p)— les fonctions logit inverse et logit, c’est-à-dire la conversion manuelle entre l’échelle log-cote (log-odds) et l’échelle de probabilité (utilisées en complément, ou à la place, de l’argumenttype = "response"vstype = "link"demarginaleffects)fixef(model)— extrait directement les coefficients à effets fixes d’un objetglmmTMBajusté, utile lorsque vous voulez faire un calcul manuel (p. ex. vérifier un résultat deslopes()) plutôt que de vous fier uniquement àmarginaleffectspredict(model, newdata, type = "link")— prédictions de modèle en R de base sur l’échelle du lien (p. ex. log-cote);marginaleffects::predictions()est généralement préférable dès qu’il s’agit d’incertitude, maispredict()apparaît dans les dérivations manuellesperformance::r2()(et les fonctions apparentéesr2_tjur(),r2_nagelkerke()) — calcule le R² d’un modèle (ou un pseudo-R² approprié pour les modèles linéaires généralisés où le R² ordinaire ne s’applique pas), utilisé dans les tableaux de résultats aux côtés des estimations d’effetstinytable::tt()— transforme un tableau de données en un tableau bien mis en forme pour un rapport ou un manuscrit; vous le verrez clore la plupart des blocs de code produisant des tableaux de résultatsglue::glue("[{low}, {high}]")— une commodité pour construire des chaînes de caractères formatées (p. ex. assembler une étiquette d’intervalle de confiance « [borne inférieure, borne supérieure] » à partir de deux colonnes numériques), plus facile à lire quepaste0()
Où concentrer votre attention avant l’atelier
Compte tenu de tout ce qui précède, voici où votre temps de lecture préparatoire est le mieux investi :
- Faire fonctionner le bloc d’installation avec succès. C’est le seul élément non négociable — faites-le maintenant, pas le matin de l’atelier.
- Relire deux fois les puces sur
marginaleffects. Ces cinq fonctions (predictions,avg_predictions,slopes,avg_slopes,hypotheses) constituent la véritable boîte à outils analytique que nous utiliserons pour interpréter et rapporter les interactions. Tout le reste de ce document existe pour appuyer leur utilisation correcte. - Ne vous souciez pas de la mécanique de la simulation. Vous verrez
rnorm(),tibble(),ifelse(), etc. faire le travail d’inventer de faux mouflons, fleurs et parcelles fictives — le but de ces lignes est de produire des données avec une interaction connue, pas de vous enseigner la simulation elle-même.
À bientôt à l’atelier — venez avec des questions sur l’interprétation et le rapport des effets d’interaction, puisque c’est là-dessus que nous passerons la majeure partie de notre temps.