Comment créer des graphiques de variables ajoutées dans R

L’analyse de régression linéaire multiple constitue l’un des piliers méthodologiques les plus omniprésents au sein des sciences empiriques contemporaines, qu’il s’agisse de l’économétrie, de l’épidémiologie ou de la psychologie quantitative. Toutefois, l’interprétation des modèles statistiques complexes se heurte fréquemment à un écueil épistémologique fondamental : la confusion récurrente entre les corrélations bivariées brutes et les contributions marginales réelles de chaque variable prédictive. Dans une matrice de données multidimensionnelle où les covariables partagent inévitablement des portions substantielles de variance, l’examen isolé des nuages de points conventionnels s’avère au mieux insuffisant, et au pire trompeur, masquant des phénomènes d’inversion de signe, des effets d’entraînement ou des dynamiques de suppression statistique.

Pour dépasser les limites inhérentes aux tableaux de coefficients purement numériques, les statisticiens ont formalisé des outils de diagnostic graphique sophistiqués capables d’isoler visuellement la trajectoire conditionnelle d’un régresseur donné. Parmi ces innovations, le graphique de variable ajoutée, désigné dans la littérature anglo-saxonne sous le terme d’Added Variable Plot (AVP) ou de graphique de régression partielle, s’impose comme un instrument d’une puissance analytique singulière. En projetant les résidus d’une double orthogonalisation, cet outil graphique permet au praticien d’observer l’association nette entre un prédicteur et la variable réponse, comme si l’ensemble des autres covariables du modèle avait été préalablement neutralisé par contrôle statistique direct.

Le présent article propose une exploration encyclopédique, théorique et opérationnelle de la création et de l’interprétation des graphiques de variables ajoutées dans l’environnement de programmation statistique R. À travers l’étude approfondie des fondements algébriques issus du théorème de Frisch-Waugh-Lovell, l’analyse d’applications concrètes en psychométrie et le déploiement de solutions logicielles automatisées ou manuelles, ce guide offre aux chercheurs, analystes de données et étudiants avancés l’ensemble des clés méthodologiques nécessaires pour diagnostiquer avec rigueur la linéarité, la présence d’observations influentes et la validité structurelle de leurs modèles de régression multivariée.

1. Introduction conceptuelle aux graphiques de variables ajoutées

1.1 Définition et terminologie statistique

Dans l’arsenal du diagnostic de modélisation statistique, le concept de graphique de variable ajoutée renvoie à une représentation bidimensionnelle spécifique conçue pour illustrer l’impact marginal d’un prédicteur continu ou catégoriel sur une variable dépendante, une fois prise en compte l’influence de l’ensemble des autres variables indépendantes incluses dans le modèle. Communément désigné sous les appellations synonymes de graphique de régression partielle, de graphique de résidu partiel ajusté ou d’Added Variable Plot, ce dispositif graphique ne doit pas être confondu avec un simple tracé de résidus ordinaires ni avec un nuage de points standardisé confrontant directement la réponse observée au régresseur d’intérêt.

La distinction formelle entre une corrélation bivariée brute et une relation conditionnelle partielle repose sur la notion d’espace vectoriel et de décomposition de la variance. Alors que le tracé bivarié traditionnel ignore superbement l’architecture de covariance sous-jacente unissant le prédicteur examiné aux autres éléments du système d’équations, l’AVP opère une projection orthogonale rigoureuse. Cette opération isole la variance résiduelle propre au régresseur d’intérêt d’une part, et la variance résiduelle propre à la réponse d’autre part. Dans la modélisation statistique linéaire multivariée contemporaine, l’AVP sert de pont indispensable entre l’abstraction algébrique des matrices de projection et la perception géométrique de l’ajustement des moindres carrés ordinaires.

Il convient de souligner que cette terminologie s’ancre dans la théorie de l’identification causale et de l’ajustement statistique développée dans des manuels de référence tels que ceux documentés par la comunauté R officielle sur le CRAN. Le terme de variable ajoutée illustre précisément la question contrefactuelle posée par l’analyste : quelle est la trajectoire d’information additionnelle apportée par cette covariable spécifique lorsqu’elle vient s’ajouter à un modèle qui comprendrait déjà tous les autres facteurs explicatifs ?

1.2 Objectifs analytiques fondamentaux

L’objectif analytique premier des graphiques de variables ajoutées réside dans la visualisation explicite de l’effet net d’un prédicteur sur la variable réponse. Dans la majorité des jeux de données observationnelles, les prédicteurs présentent un degré non négligeable d’intercorrélation, phénomène souvent qualifié de multicolinéarité quasi-stochastique. Dès lors, le coefficient d’une variable au sein d’une régression multiple ne quantifie nullement son lien direct et isolé avec la réponse, mais bien sa contribution résiduelle conditionnelle. L’AVP matérialise graphiquement cette dynamique en neutralisant de manière bidirectionnelle la variance partagée avec les autres covariables, offrant ainsi une lecture directe de la relation épurée de toute contamination parasite.

Au-delà de la simple visualisation de la pente, ces graphiques jouent un rôle prépondérant dans l’évaluation intuitive de la linéarité résiduelle. Un tableau d’estimation paramétrique classique se contente d’imposer une contrainte de linéarité globale et de fournir un estimateur ponctuel assorti d’une erreur-type. Or, cet estimateur peut masquer une forme fonctionnelle curviligne, exponentielle ou discontinue qui invaliderait formellement l’hypothèse de spécification du modèle linéaire. L’examen attentif d’un AVP permet d’inspecter visuellement si la relation entre les résidus demeure strictement rectiligne ou si elle requiert des transformations mathématiques complémentaires, telles que l’incorporation de polynômes ou de fonctions splines.

Enfin, l’apport qualitatif des AVP surpasse largement les indices numériques d’un tableau d’estimation récapitulatif. Tandis qu’une statistique de Student ou une valeur de p synthétise l’effet sous une forme univariée condensée sensible aux artéfacts d’échantillonnage, le tracé de variable ajoutée contextualise chaque observation individuelle, rendant immédiatement perceptibles l’hétérogénéité locale des données et les anomalies de dispersion que les métriques agrégées ont tendance à occulter.

1.3 Historique et évolution dans l’analyse de régression

L’ancrage historique des représentations de régression partielle remonte aux fondations mêmes de l’économétrie moderne et de l’analyse macro-structurelle au cours du premier tiers du vingtième siècle. Les premières formulations conceptuelles de l’orthogonalisation séquentielle trouvent leur origine dans les travaux séminaux de Ragnar Frisch et Frederick V. Waugh publiés en 1933, ultérieurement généralisés et complétés par Michael C. Lovell en 1963. Ce corpus théorique, immortalisé sous le nom de théorème de Frisch-Waugh-Lovell, démontra mathématiquement qu’une régression multivariée complexe pouvait être scindée sans perte d’information en une succession de régressions univariées appliquées à des variables résidualisées.

Néanmoins, la formalisation moderne et la dénomination explicite des graphiques de variables ajoutées sous une perspective de diagnostic statistique opérationnel sont principalement attribuables aux contributions majeures de David A. Belsley, Edwin Kuh et Roy E. Welsch dans leur ouvrage fondateur de 1980, intitulé Regression Diagnostics: Identifying Influential Data and Sources of Collinearity, ainsi qu’aux travaux parallèles de R. Dennis Cook et Sanford Weisberg. Ces auteurs ont déplacé la focale depuis la simple équivalence calculatoire des estimateurs vers l’évaluation de la stabilité du modèle et l’identification des points atypiques exerçant une influence démesurée sur les surfaces de réponse.

L’avènement de la micro-informatique et l’évolution subséquente des logiciels d’analyse ont profondément bouleversé l’accessibilité de ces outils. Historiquement confinés à des calculs matriciels manuels fastidieux ou à des scripts personnalisés réservés aux spécialistes d’algorithmes numériques, les AVP ont été intégrés nativement au sein des bibliothèques de calcul moderne. Dans l’écosystème R, cette transition s’est matérialisée par l’élaboration de modules spécialisés permettant de générer ces graphiques en une seule ligne de commande, démocratisant ainsi une démarche méthodologique autrefois prohibitive.

2. Fondements mathématiques et dérivation statistique

2.1 Le théorème de Frisch-Waugh-Lovell (FWL)

Pour appréhender rigoureusement les propriétés des graphiques de variables ajoutées, il est impératif d’examiner la formulation algébrique du théorème de Frisch-Waugh-Lovell au sein du cadre canonique des moindres carrés ordinaires. Considérons un modèle de régression linéaire formulé sous notation matricielle où le vecteur de la variable dépendante Y, de dimension n par 1, est régressé sur une matrice de prédicteurs partitionnée en deux sous-ensembles distincts : X1 de dimension n par k1, et X2 de dimension n par k2. Le modèle complet s’énonce sous la forme matricielle Y = X1 * Beta1 + X2 * Beta2 + Epsilon, où Epsilon désigne le vecteur des perturbations stochastiques non observées.

Le théorème de Frisch-Waugh-Lovell stipule que l’estimateur des moindres carrés ordinaires associé au bloc de coefficients Beta2 peut être obtenu de manière rigoureusement identique via une procédure de projection orthogonale séquentielle. Définissons la matrice de projection orthogonale, ou matrice annihilatrice associée à l’espace vectoriel engendré par les colonnes de X1, par l’expression M1 = I – X1 * (X1′ * X1)^(-1) * X1′, où I représente la matrice identité de dimension n par n, et X1′ désigne la transposée de la matrice X1. Cette matrice de projection M1 est symétrique et idempotente, ce qui implique les égalités fondamentales M1′ = M1 et M1 * M1 = M1.

En pré-multipliant le modèle global par M1, nous obtenons l’équation transformée M1 * Y = M1 * X1 * Beta1 + M1 * X2 * Beta2 + M1 * Epsilon. Or, en vertu de la définition de la matrice résiduelle, le produit M1 * X1 s’annule identiquement, car toute colonne de X1 est strictement orthogonale à l’espace des résidus généré par M1. Dès lors, le système se simplifie en M1 * Y = M1 * X2 * Beta2 + M1 * Epsilon. L’estimateur ponctuel de Beta2 issu de cette régression projetée s’écrit algébriquement ((M1 * X2)’ * (M1 * X2))^(-1) * (M1 * X2)’ * (M1 * Y). Grâce à l’idempotence de M1, cette expression se réduit immédiatement à (X2′ * M1 * X2)^(-1) * X2′ * M1 * Y, ce qui correspond exactement à l’estimateur de Beta2 obtenu lors de l’inversion globale de la matrice partitionnée du modèle initial.

2.2 Mécanisme de double résidualisation

L’implémentation opérationnelle du théorème FWL dans le cadre d’un graphique de variable ajoutée se traduit par un mécanisme séquentiel qualifié de double résidualisation. Supposons que l’analyste souhaite isoler la contribution spécifique d’un unique prédicteur scalaire Xk au sein d’un modèle comprenant un ensemble d’autres covariables regroupées sous la matrice notée X_(-k), incluant le vecteur colonne unitaire correspondant à la constante d’ajustement. La procédure se décompose en deux régressions auxiliaires parfaitement distinctes.

La première régression auxiliaire consiste à régresser la variable dépendante d’origine Y sur l’ensemble des prédicteurs résiduels X_(-k). Cette estimation par les moindres carrés génère un vecteur de valeurs ajustées, à partir duquel on extrait le premier vecteur de résidus, noté e(Y|X_(-k)) = Y – Y_ajusté_sans_Xk. Ce résidu capture la totalité des variations de la variable réponse qui ne peuvent pas être expliquées linéairement par les covariables d’ajustement. Il représente ce qu’il reste à expliquer dans la réponse après avoir épuisé le pouvoir prédictif de toutes les autres variables du système.

La seconde régression auxiliaire prend le prédicteur cible Xk comme variable dépendante temporaire et le régresse sur la même matrice de covariables d’ajustement X_(-k). On en extrait le second vecteur de résidus, noté e(Xk|X_(-k)) = Xk – Xk_ajusté_sans_Xk. Ce second vecteur résiduel matérialise la variation endogène et idiosyncratique de la variable Xk qui est totalement indépendante et non corrélée aux autres prédicteurs. Le graphique de variable ajoutée est alors rigoureusement défini comme le tracé en coordonnées cartésiennes de e(Y|X_(-k)) en ordonnée contre e(Xk|X_(-k)) en abscisse. La pente de la droite des moindres carrés ajustée sur ce nuage résiduel bivarié correspond précisément au coefficient partiel de Xk dans la régression multivariée initiale.

2.3 Propriétés statistiques de la pente partielle

L’intérêt méthodologique fondamental de l’AVP repose sur un ensemble de propriétés de conservation mathématique exceptionnelles. La première propriété, issue directement des développements précédents, est l’identité numérique absolue entre la pente de la droite de régression simple ajustée sur les résidus et le coefficient partiel Beta_k issu du modèle de régression multiple global. Si l’on estime la régression univariée e(Y|X_(-k)) = Beta_partiel * e(Xk|X_(-k)) + u, l’estimateur OLS de Beta_partiel vérifie sans la moindre approximation l’égalité Beta_partiel = Beta_k.

La deuxième propriété concerne la préservation intégrale des résidus finaux du modèle. Les résidus issus de la régression simple des résidus sur les résidus sont strictement équivalents aux résidus Epsilon du modèle multiple initial complet. En effet, u = e(Y|X_(-k)) – Beta_k * e(Xk|X_(-k)) = M_(-k) * Y – Beta_k * M_(-k) * Xk = M_(-k) * (Y – Beta_k * Xk), ce qui correspond point par point au vecteur d’erreurs du modèle multivarié originel.

Cependant, une nuance critique d’ordre inférentiel doit être impérativement soulignée quant aux degrés de liberté. Bien que la somme des carrés des résidus de la régression univariée auxiliaire soit identique à celle du modèle complet, la régression univariée simple ne prend en compte par défaut que n – 2 degrés de liberté pour l’estimation de la variance de l’erreur. Or, le modèle multivarié d’origine consomme n – p – 1 degrés de liberté, où p représente le nombre total de paramètres estimés. Par conséquent, pour que l’erreur-type du coefficient, la statistique t de Student et la valeur de p associées à la pente de l’AVP coïncident rigoureusement avec les résultats du modèle multivarié, le calcul de la variance de l’erreur doit impérativement intégrer la correction relative aux degrés de liberté globaux du système.

3. Pertinence et applications en psychologie quantitative

3.1 Résolution des problèmes de confusion en sciences du comportement

En psychologie quantitative et dans l’ensemble des sciences comportementales, les devis de recherche sont très fréquemment quasi-expérimentaux ou observationnels. Les chercheurs sont confrontés à des matrices de corrélations denses au sein desquelles des variables parasites telles que l’âge chronologique, le statut socio-économique, le niveau de scolarité ou l’état de fatigue cognitive viennent systématiquement polluer les relations fonctionnelles théoriques. L’omission ou la mauvaise prise en compte de ces facteurs confondants expose l’analyste à des inférences totalement erronées, pouvant aller jusqu’à l’observation empirique d’associations fallacieuses.

Les graphiques de variables ajoutées constituent un outil visuel de premier ordre pour démasquer des mécanismes complexes de suppression statistique ou d’inversion d’effet. Un exemple canonique réside dans le cas où une corrélation brute bivariée entre un marqueur psychologique et une performance cognitive apparaît initialement positive, alors que l’effet conditionnel ajusté sur l’âge se révèle fortement négatif. L’AVP illustre avec une clarté remarquable comment l’élimination de la covariance due à l’âge redresse la trajectoire de la relation et fait émerger la réalité de l’effet clinique sous-jacent. L’inspection graphique offre au chercheur l’assurance visuelle que l’effet rapporté ne constitue pas un artefact provoqué par une distribution déséquilibrée des variables sociodémographiques.

De surcroît, la visualisation des résidus partiels permet d’identifier si l’effet d’une intervention psychologique ou d’une variable différentielle persiste de manière homogène à travers les différents strates des covariables de contrôle, fournissant ainsi une validation empirique robuste des hypothèses causales formulées dans le cadre de protocoles non randomisés.

3.2 Diagnostic des modèles de personnalité et d’attitudes

L’évaluation des structures psychométriques, notamment dans le cadre des modèles de personnalité tels que le modèle en cinq facteurs (Big Five) ou l’évaluation des attitudes par échelles de Likert, se heurte inévitablement au problème de la corrélation modérée à élevée entre sous-dimensions théoriques. Par exemple, l’anxiété et la dépression partagent une part substantielle de variance commune attribuable à un facteur de second ordre d’affectivité négative générale ou de névrosisme. Lorsqu’un psychologue régresse un comportement d’évitement sur ces deux construits simultanément, l’interprétation des coefficients peut devenir extrêmement délicate en raison de cette colinéarité intrinsèque.

Dans ce contexte, l’AVP permet d’évaluer directement l’incrément de validité prédictive d’une sous-échelle spécifique au-delà des autres facettes psychométriques. En observant le graphique de variable ajoutée dédié à l’anxiété après contrôle de la dépression, le praticien peut immédiatement apprécier la part de variance unique qui contribue à la prédiction du critère. Si le nuage de points résiduels s’affaisse le long d’une bande horizontale amorphe, cela indique visuellement que la variable examinée n’apporte aucune information discriminante supplémentaire, quand bien même sa corrélation simple avec la variable dépendante serait hautement significative sur le plan statistique.

Cette approche graphique protège les psychologues contre l’écueil de la redondance conceptuelle des instruments de mesure. Elle facilite l’épuration des batteries de tests en mettant en lumière les composantes qui n’expliquent aucune variation spécifique du comportement une fois les dimensions centrales contrôlées.

3.3 Éthique méthodologique et reproductibilité

À l’ère des crises de réplicabilité qui ont traversé la psychologie et les sciences biomédicales, la transparence analytique s’est muée en un impératif méthodologique incontournable. Les publications scientifiques ont trop longtemps reposé de manière exclusive sur la présentation de matrices de corrélation et de tableaux de régression formatés selon les normes de l’American Psychological Association (APA), sans fournir d’éléments tangibles quant à la validité structurelle sous-jacente des modèles ajustés. Les tableaux numériques présentent en effet l’inconvénient majeur de masquer les distorsions locales causées par des sous-échantillons particuliers ou des valeurs aberrantes de levier.

L’intégration systématique des graphiques de variables ajoutées dans les protocoles de diagnostic pré-publication répond directement aux exigences accrues de rigueur méthodologique. Selon les préconisations du groupe de travail sur l’inférence statistique de l’APA, les chercheurs sont vivement encouragés à démontrer que leurs conclusions ne dépendent pas d’une poignée d’observations atypiques logées dans des régions marginales de l’espace multidimensionnel. L’AVP offre une attestation graphique transparente montrant que l’effet linéaire proclamé traverse de part en part la distribution échantillonnale de manière robuste.

Adopter une telle démarche de diagnostic visuel participe activement à la lutte contre les pratiques de recherche douteuses, telles que le dragage de données (p-hacking), en forçant l’expérimentateur à confronter ses hypothèses théoriques à la géométrie réelle de ses résidus conditionnels avant de revendiquer la découverte d’un effet psychologique statistiquement robuste.

4. Préparation de l’environnement de calcul dans R

4.1 Installation et chargement de l’écosystème car

L’implémentation la plus emblématique, la plus robuste et la plus largement documentée des graphiques de variables ajoutées au sein de l’écosystème open source R est hébergée au sein du package spécialisé nommé car, acronyme de Companion to Applied Regression, développé et maintenu sous la direction de John Fox et Sanford Weisberg. Cette bibliothèque logicielle constitue le prolongement computationnel direct des manuels méthodologiques de référence en sociologie quantitative et en économétrie appliquée, fournissant des suites de fonctions diagnostiques introuvables dans la distribution de base de R.

Pour exploiter les potentialités offertes par ce package, l’utilisateur doit procéder à son installation depuis les serveurs officiels du Comprehensive R Archive Network en exécutant l’instruction conventionnelle install.packages(« car »). Il est pertinent de noter que le package car s’appuie sur une hiérarchie de dépendances logicielles rigoureuses, incluant notamment carData, qui héberge les jeux de données classiques d’illustration, ainsi que des modules matriciels de calcul d’erreurs-types robustes tels que sandwich ou lme4.

Une fois l’installation menée à son terme, le chargement du module s’effectue via l’appel library(car). Lors de cette phase, l’analyste veillera à surveiller les messages d’avertissement relatifs aux conflits d’espaces de noms. En effet, certaines fonctions génériques de car peuvent entrer en collision lexicale avec des modules de manipulation de données appartenant au tidyverse, à l’image des fonctions select ou recode. Une saine pratique de programmation consistera, en cas d’ambiguïté au sein de scripts industriels, à spécifier explicitement le préfixe du package par l’opérateur de double deux-points car::avPlots.

4.2 Structure et exploration des jeux de données d’exemple

Afin d’illustrer la mécanique concrète des AVP, le recours à des jeux de données historiques intégrés permet de garantir la reproductibilité immédiate des résultats. Le jeu de données Duncan, inclus dans le package carData, représente un cas d’école universellement plébiscité. Ce recueil de données documente le prestige professionnel de quarante-cinq professions aux États-Unis au cours des années 1950, en fonction du niveau d’éducation moyen des professionnels et du pourcentage de ces derniers percevant un revenu supérieur à un seuil monétaire déterminé.

Un autre cadre expérimental d’un intérêt remarquable réside dans le jeu de données Prestige, documentant le statut socio-économique de 102 professions canadiennes en relation avec le revenu moyen, le niveau d’instruction et la proportion de femmes exerçant dans ces domaines. Dans une optique psychométrique alternative, l’analyste peut générer des structures synthétiques simulant des passations de tests psychologiques standardisés comprenant des scores de mémoire de travail, de vitesse de traitement, de régulation émotionnelle et d’âge chronologique.

Avant d’ajuster le moindre modèle de régression, une phase exploratoire préliminaire s’impose. Il convient de vérifier scrupuleusement la typologie des variables via les fonctions structurelles de R (telles que str() ou summary()). L’analyste doit s’assurer que les variables catégorielles sont formellement encodées sous forme de facteurs explicites et que les covariables continues ne présentent pas de singularités d’encodage, comme des valeurs aberrantes impossibles résultant d’erreurs de saisie manuelle.

4.3 Vérification des prérequis computationnels

L’exactitude mathématique des résidus projetés par les fonctions de variables ajoutées dépend en amont de la conformité de la matrice de conception. Le premier prérequis concerne la gestion des observations incomplètes ou valeurs manquantes. Par défaut, les fonctions d’estimation de modèles linéaires sous R appliquent l’option na.action = na.omit, qui exclut de manière listwise toute observation comportant au moins une donnée manquante sur l’une des variables mobilisées dans la formule.

Ce mécanisme d’exclusion automatique peut altérer considérablement la composition de l’échantillon analysé si les données manquantes ne sont pas réparties de façon totalement aléatoire. De surcroît, si l’on tente d’extraire manuellement des résidus à partir de sous-modèles qui n’intègrent pas exactement les mêmes observations en raison d’exclusions différentielles, les vecteurs résiduels ne posséderont plus la même longueur dimensionnelle n, provoquant l’interruption immédiate des calculs matriciels. Il est donc recommandé d’assainir le tableau de données en amont en figeant un sous-échantillon complet dédié à la régression.

Le second prérequis computationnel touche à la vérification du rang de la matrice de covariance. En présence d’une colinéarité parfaite entre deux ou plusieurs prédicteurs, la matrice X’X devient singulière et non inversible, forçant R à éliminer unilatéralement des colonnes d’estimation. Enfin, sur le plan ergonomique, l’utilisateur configurera sa fenêtre de visualisation graphique au sein de son environnement de développement (tel que RStudio) pour accueillir des affichages multipanneaux matriciels sans distorsion de ratio hauteur-largeur.

5. Ajustement du modèle de régression linéaire multiple de base

5.1 Spécification formelle de la syntaxe lm()

L’ancrage opératoire de l’analyse débute par l’ajustement formel du modèle de régression linéaire multiple à l’aide de la fonction native lm() de R. Sur le plan mathématique, nous formulons l’hypothèse que la variable dépendante continue Y se structure comme une combinaison linéaire de k prédicteurs continus ou catégoriels, augmentée d’un terme résiduel d’erreur non corrélé de moyenne nulle et de variance constante.

Au niveau de la syntaxe R, la relation s’énonce selon le formalisme symbolique de Wilkinson et Rogers au moyen de l’opérateur tilde. Par exemple, si l’on modélise le prestige professionnel à partir du revenu (income) et de l’éducation (education) au sein du tableau Duncan, l’instruction s’écrit de la manière suivante : modele_duncan <- lm(prestige ~ income + education, data = Duncan). Cette écriture instancie un objet statistique de classe formelle « lm », qui encapsule au sein d’une liste structurée les coefficients de régression estimés par l’algorithme de décomposition QR, le vecteur des valeurs ajustées, les résidus ordinaires, ainsi que la structure de la matrice de conception modèle.

Le stockage de cet objet ajusté dans une variable dédiée constitue une étape impérative. En effet, toutes les fonctions de diagnostic ultérieures, qu’il s’agisse des procédures natives d’extraction de résidus ou des algorithmes spécialisés du package car, prendront cet objet modélisé comme argument d’entrée fondamental pour en extraire la géométrie interne sans qu’il soit nécessaire de respécifier les formules initiales.

5.2 Analyse paramétrique des sorties de summary()

L’exécution de la commande standard summary(modele_duncan) produit un condensé statistique textuel divisé en plusieurs compartiments analytiques hautement normalisés. La première section fournit la distribution empirique des résidus ordinaires (minimum, quartiles, médiane, maximum), permettant de juger de prime abord de la symétrie de l’erreur. Vient ensuite le tableau fondamental des coefficients paramétriques, structuré en quatre colonnes : l’estimation ponctuelle du paramètre (Estimate), l’erreur-type conventionnelle (Std. Error), la statistique de test t de Student (t value) et la probabilité critique bilatérale associée (Pr(>|t|)).

Bien que ce tableau indique la significativité marginale apparente de chaque prédicteur, sa lecture isolée comporte des risques d’interprétation majeurs. Par exemple, une valeur de p extrêmement basse associée au prédicteur du revenu atteste que l’on peut rejeter l’hypothèse nulle d’un coefficient partiel nul dans la population, sous l’hypothèse stricte que le modèle est correctement spécifié. Toutefois, cette métrique demeure totalement aveugle quant à l’uniformité de l’effet à travers l’échantillon. Le coefficient de détermination global R-carré et sa version corrigée (Adjusted R-squared), situés en bas de sortie textuelle, quantifient la proportion globale de variance de la réponse expliquée conjointement par les prédicteurs, mais n’apportent aucun éclairage sur la validité structurelle locale des relations individuelles.

C’est précisément face à ces lacunes inhérentes aux métriques scalaires agrégées que la transition vers une analyse graphique diagnostique s’avère méthodologiquement indispensable pour tout analyste soucieux de robustesse.

5.3 Validation des hypothèses classiques de Gauss-Markov

Pour que les estimateurs des moindres carrés ordinaires soient déclarés optimaux au sens du théorème de Gauss-Markov (garantissant qu’ils sont les meilleurs estimateurs linéaires sans biais, ou BLUE), un ensemble d’hypothèses probabilistes fondamentales doit être vérifié : l’espérance conditionnelle des erreurs doit être nulle, la variance des perturbations doit être strictement constante (homoscédasticité), les erreurs doivent être mutuellement non corrélées, et dans une perspective inférentielle exacte, les résidus doivent suivre une loi normale multivariée.

L’approche diagnostique traditionnelle consiste à inspecter les graphiques par défaut générés par l’instruction plot(modele_duncan). Ces tracés comprennent généralement le graphique des résidus ordinaires en fonction des valeurs prédites, le tracé quantile-quantile normal (Q-Q plot), le graphique de l’échelle de dispersion standardisée et le tracé des résidus contre les leviers statistiques avec les frontières de distance de Cook. Bien que cruciaux, ces diagnostics globaux présentent une faiblesse intrinsèque : ils collapsed l’ensemble de la structure multivariée sur une dimension unique, celle de la prédiction synthétique Y_chapeau.

Par conséquent, si une violation d’hypothèse (telle qu’une courbure non linéaire flagrante ou une hétéroscédasticité aiguë) n’affecte qu’un seul des prédicteurs inclus dans le modèle, elle sera souvent diluée ou partiellement masquée lors de l’examen de ces résidus globaux. L’analyste se trouve alors dans l’incapacité d’imputer l’anomalie à un paramètre déterminé. Il devient donc indispensable de décomposer la structure d’ajustement prédicteur par prédicteur, ce qui constitue la raison d’être exacte des graphiques de variables ajoutées.

6. Génération des graphiques avec la fonction avPlots()

6.1 Syntaxe canonique et exécution basique

L’extraction et l’affichage direct de l’ensemble des graphiques de variables ajoutées pour un modèle de régression linéaire ajusté s’exécutent avec une simplicité déconcertante grâce à la fonction avPlots() fournie par la bibliothèque car. Dès lors que l’objet résultant de l’ajustement du modèle linéaire a été stocké sous l’identifiant modele_duncan, l’analyste n’a qu’à soumettre au terminal l’appel canonique suivant : avPlots(modele_duncan). Cette seule instruction déclenche une série d’opérations matricielles sous-jacentes totalement automatisées.

La fonction procède séquentiellement à la double résidualisation pour chaque variable indépendante incluse dans le modèle, extrait les vecteurs résiduels correspondants et divise automatiquement la fenêtre graphique active de R en une grille multi-panneaux symétrique. Dans le cas d’un modèle comportant deux covariables explicatives, avPlots() génère une matrice graphique de deux panneaux distincts disposés côte à côte, chacun dédié à un prédicteur unique.

Le décodage sémiotique des axes est fondamental pour éviter toute confusion lors de la lecture immédiate. L’axe des ordonnées affiche systématiquement les résidus de la variable réponse conditionnés par les autres prédicteurs, généralement étiquetés sous la nomenclature prestige | others (ce qui signifie littéralement : la part de prestige résiduelle non expliquée par les autres variables du modèle). L’axe des abscisses présente les résidus du prédicteur concerné après neutralisation de ses partenaires de covariance, par exemple income | others. Chaque point tracé dans l’espace graphique correspond strictement à une unité statistique de l’échantillon d’origine, projetée selon ses deux coordonnées résiduelles.

6.2 Sélection ciblée de prédicteurs spécifiques

Lorsque le modèle statistique investigué intègre un volume important de régresseurs (par exemple dans le cadre de modèles de contrôle épidémiologiques ou psychométriques comprenant dix, quinze ou vingt covariables), l’appel indistinct de la commande avPlots() aboutit à la génération d’une grille graphique surchargée. Dans un tel agencement, la taille individuelle de chaque sous-panneau se trouve drastiquement réduite, rendant l’inspection visuelle des nuages de points fastidieuse et imprécise.

Pour parer à cette contrainte d’espace, le package met à disposition deux stratégies opérationnelles d’une grande souplesse. La première consiste à employer l’argument dédié terms au sein de la fonction collective, sous la forme avPlots(modele_duncan, terms = ~ income). L’utilisateur peut ainsi transmettre une formule partielle listant le ou les seuls prédicteurs qu’il souhaite soumettre à l’examen visuel, préservant ainsi la lisibilité de la fenêtre de rendu.

La seconde alternative repose sur l’usage de la fonction au singulier avPlot(), conçue spécifiquement pour cibler un régresseur isolé à la fois. Sa syntaxe requiert la transmission de l’objet du modèle suivie de la désignation sous forme de chaîne de caractères du prédicteur d’intérêt : avPlot(modele_duncan, variable = « income »). Cette commande unitaire permet de générer un graphique d’une haute résolution occupant l’intégralité du périphérique graphique, idéal pour l’exportation vers des supports de diffusion scientifique ou des présentations académiques ciblées.

6.3 Gestion des paramètres graphiques natifs de avPlots

La fonction avPlots() intègre un ensemble d’arguments formels permettant d’adapter le rendu esthétique et informatif des tracés sans nécessiter de post-traitement complexe. Par défaut, la fonction superpose sur le nuage de points résiduels une droite de régression linéaire ajustée par les moindres carrés ordinaires, représentant visuellement la pente partielle du coefficient estimé. Les attributs graphiques de cette droite peuvent être ajustés par le truchement des paramètres col.lines (pour modifier la couleur de la droite d’ajustement) et lwd (définissant la largeur du trait vectoriel).

Les marqueurs ponctuels représentant les observations individuelles peuvent également être modulés en couleur (col), en géométrie de symbole (pch) et en échelle d’agrandissement (cex). De surcroît, la fonction insère par défaut un quadrillage en arrière-plan afin de matérialiser visuellement les axes de référence d’ordonnée et d’abscisse nulles (passant par le point d’intersection central c(0,0)). Ce comportement peut être contrôlé par l’argument booléen grid = TRUE ou FALSE.

Enfin, les titres principaux et les intitulés des axes peuvent être entièrement personnalisés en surchargeant les valeurs assignées par défaut au moyen des arguments textuels usuels main, xlab et ylab. Cette flexibilité permet au chercheur francophone de substituer aux libellés automatiques en langue anglaise des dénominations épistémologiquement adaptées à la diffusion de ses travaux au sein de revues scientifiques francophones.

7. Construction manuelle d’un graphique de variable ajoutée

7.1 Algorithme de calcul étape par étape sans package tiers

Bien que l’automatisation proposée par les bibliothèques spécialisées soit d’une redoutable efficacité, la programmation manuelle étape par étape d’un graphique de variable ajoutée à l’aide des seules fonctions fondamentales de l’interpréteur R constitue un exercice d’une valeur pédagogique inestimable. Cet exercice permet de lever toute opacité quant à la nature des calculs matriciels effectués sous le capot. Considérons à titre d’exemple le modèle reliant le statut social (Y) aux variables d’éducation (X1) et de revenu (X2).

L’algorithme de décomposition s’amorce par la mise en œuvre de la première régression auxiliaire visant à résidualiser la variable réponse Y vis-à-vis du revenu X2. L’analyste programme l’instruction : reg_y_sans_x1 <- lm(prestige ~ income, data = Duncan). À l’issue de cet ajustement partiel, nous extrayons le vecteur résiduel de la réponse, qui renferme la variance du prestige non imputable au revenu, via l’assignation : res_y <- residuals(reg_y_sans_x1).

La seconde phase consiste à résidualiser le prédicteur d’intérêt éducation (X1) vis-à-vis de l’autre covariable revenu (X2). Nous formulons l’ajustement : reg_x1_sans_x2 <- lm(education ~ income, data = Duncan), puis nous procédons à l’extraction de son vecteur d’erreurs : res_x1 <- residuals(reg_x1_sans_x2). Ayant isolé les deux composantes orthogonales, nous assemblons ces deux vecteurs au sein d’une structure de tableau de données élémentaire data.frame(res_x1, res_y). Enfin, la construction visuelle s’exécute par le recours aux primitives graphiques de base : plot(res_x1, res_y), complétée immédiatement par l’adjonction de la droite de tendance linéaire via l’instruction abline(lm(res_y ~ res_x1), col = « red », lwd = 2).

7.2 Vérification de la conformité empirique

Une fois le tracé manuel déployé, il est indispensable de procéder au contrôle formel de son identité statistique avec les estimateurs du modèle multivarié global. Ajustons la régression simple des résidus manuels : modele_verif <- lm(res_y ~ res_x1). L’interrogation de son coefficient par l’expression coef(modele_verif)[« res_x1 »] révèle une équivalence numérique rigoureusement parfaite, jusqu’à la dernière décimale de précision flottante de la machine, avec le coefficient obtenu dans le modèle global modele_complet <- lm(prestige ~ education + income, data = Duncan) interrogé par coef(modele_complet)[« education »].

De plus, l’examen de l’ordonnée à l’origine (l’intercept) de ce modèle de vérification démontre une valeur strictement égale à zéro (aux erreurs d’arrondi numérique près, de l’ordre de 10 à la puissance -16). Cela découle directement de la propriété mathématique des moindres carrés voulant que la moyenne de tout vecteur de résidus issu d’un modèle incluant une constante soit identiquement nulle.

Enfin, si l’on superpose graphiquement le tracé manuel ainsi généré avec la sortie produite par car::avPlot(modele_complet, « education »), on constate une coïncidence géométrique absolue de l’ensemble des points et de la pente de la droite. Cette démonstration empirique confirme sans équivoque que l’AVP n’est rien d’autre que la projection directe de la relation épurée par le théorème de Frisch-Waugh-Lovell.

7.3 Avantages pédagogiques de la programmation manuelle

Développer soi-même l’algorithme d’orthogonalisation présente des vertus heuristiques indéniables pour les chercheurs et les étudiants en méthodologie quantitative. Cette démarche démystifie complètement la notion parfois opaque de contrôle statistique. Alors que l’expression verbale contrôler pour la variable Z peut sembler abstraite dans un manuel de méthodologie, le fait de constater concrètement que ce contrôle consiste simplement à retirer des variables Y et X leurs prédictions linéaires respectives fondées sur Z confère une clarté conceptuelle définitive à l’opération.

Par ailleurs, la maîtrise de cette mécanique manuelle libère le chercheur des contraintes d’implémentation imposées par les packages tiers. Par exemple, si l’analyste souhaite appliquer des estimateurs de régression non linéaires ou hautement robustes (tels que des estimateurs de type M, MM ou S implémentés via des algorithmes de régression quantile ou des moindres carrés pondérés) pour lesquels aucune fonction avPlots automatisée n’a été implémentée par les développeurs, la logique de double résidualisation manuelle demeure intégralement transposable et opérationnelle.

L’expérimentateur peut ainsi concevoir des représentations de variables ajoutées semi-paramétriques ou robustes totalement sur-mesure, en appliquant des algorithmes de filtrage adaptatif en amont de la mise en relation graphique des deux résidus finaux.

8. Interprétation diagnostique avancée des tracés

8.1 Lecture de la direction et de la pente partielle

L’interprétation analytique d’un graphique de variable ajoutée requiert une attention aiguë portée à plusieurs caractéristiques géométriques fondamentales du nuage résiduel. Le premier élément diagnostique réside dans l’orientation et l’inclinaison de la droite d’ajustement partielle. Une pente ascendante robuste indique que, à niveau d’éducation constant, une augmentation du revenu est associée de manière conditionnelle à une progression du prestige perçu.

À l’inverse, une pente descendante matérialiserait une association négative nette, révélant parfois des phénomènes insoupçonnés d’inversion d’effet par rapport à la corrélation brute bivariée initiale. L’absence notable d’inclinaison (une droite parfaitement horizontale confondue avec l’axe des abscisses) traduit visuellement le constat que la variable d’intérêt n’exerce aucun effet explicatif incrémental sur la réponse une fois les autres prédicteurs mobilisés dans l’architecture du modèle.

Il importe également d’évaluer visuellement la dispersion des points autour de cette droite. Plus le nuage de points résiduels est étroitement resserré autour de la trajectoire linéaire, plus l’erreur-type du coefficient partiel est restreinte, traduisant une grande précision de l’estimateur. Si au contraire le nuage adopte une conformation diffuse et amorphe sans structure directionnelle dominante, la statistique t sera mécaniquement faible et la valeur de p élevée, matérialisant l’incapacité de la variable à réduire substantiellement la variance inexpliquée du système.

8.2 Détection des relations non linéaires sous-jacentes

L’un des apports diagnostiques les plus spectaculaires des AVP réside dans leur propension à révéler des distorsions de spécification fonctionnelle que les indicateurs numériques conventionnels dissimulent totalement. L’estimateur OLS contraint par construction la relation à s’exprimer sous une forme rectiligne stricte, matérialisée par la droite superposée. Cependant, l’inspection visuelle du nuage résiduel peut mettre en évidence des configurations géométriques curvilignes manifestes.

Si l’on observe que les résidus forment une parabole en U (ou en U inversé), avec des points situés systématiquement au-dessus de la droite aux deux extrémités de l’axe horizontal et sous la droite dans la région centrale, cette structure incurvée signale sans ambiguïté une relation quadratique sous-jacente entre le prédicteur et la variable réponse. Le modèle linéaire additif initial souffre alors d’un biais de sous-spécification. L’analyste se voit contraint de corriger sa formule en incorporant un terme d’ordre supérieur (tel qu’un polynôme de degré 2) ou en recourant à des fonctions non linéaires ou splines de lissage.

D’autres morphologies fonctionnelles peuvent également émerger, telles que des dynamiques de saturation asymptotique ou des effets de seuil brutal, invitant le modélisateur à appliquer des transformations logarithmiques ou logistiques sur la variable originale avant toute réestimation.

8.3 Identification de l’hétéroscédasticité conditionnelle

La variance des erreurs conditionnelles constitue l’un des piliers de validité des inférences paramétriques basées sur la loi de Student. L’examen des graphiques de variables ajoutées offre une modalité diagnostique particulièrement fine pour traquer l’hétéroscédasticité conditionnelle spécifiquement associée à un prédicteur déterminé. Dans une situation conforme aux hypothèses théoriques d’homoscédasticité, la dispersion verticale des points autour de la droite de régression doit demeurer rigoureusement uniforme et stationnaire d’un bout à l’autre de l’axe des abscisses.

Si l’analyste identifie visuellement une morphologie en éventail ou en entonnoir (par exemple, des résidus très resserrés pour les faibles valeurs de e(Xk|others) mais qui s’évasent considérablement à mesure que l’on progresse vers les valeurs positives élevées), cela démontre que la variance de la perturbation dépend intimement de l’amplitude de ce régresseur particulier. Une telle configuration invalide formellement le calcul classique des erreurs-types par la méthode de Gauss-Markov.

Face à cette découverte graphique, le chercheur devra obligatoirement substituer aux erreurs-types conventionnelles des estimateurs d’erreurs-types robustes à l’hétéroscédasticité (comme les estimateurs sandwich de type HC3 ou HC4 théorisés par Halbert White et Long & Ervin) ou envisager des transformations stabilisatrices de variance sur la variable dépendante d’origine par le truchement de la famille de transformations de Box-Cox.

9. Détection et traitement des observations influentes

9.1 Identification automatique des points atypiques avec id=list()

L’un des atouts opérationnels majeurs de l’implémentation avPlots() au sein de la bibliothèque car réside dans son protocole automatisé de détection et d’étiquetage textuel des observations potentiellement problématiques. Dans un espace multidimensionnel, une observation peut présenter un profil d’atypicité tout à fait invisible lors de l’examen univarié de ses composantes brutes, mais devenir extraordinairement saillante une fois projetée sur le sous-espace des résidus partiels.

Par défaut, la fonction avPlots() applique un algorithme de marquage statistique qui affiche le numéro d’identifiant (ou l’étiquette de ligne, rownames) des points les plus extrêmes directement à proximité immédiate de leurs marqueurs graphiques respectifs. Ce comportement d’identification est paramétrable via l’argument formel id. Cet argument accepte une liste de sous-paramètres d’une grande précision technique, tels que id = list(method = « mahal », n = 3, cex = 0.8, col = « blue »).

Le sous-paramètre n régit le nombre quantitatif d’observations devant être obligatoirement étiquetées sur chaque panneau graphique. L’argument method définit le critère statistique présidant à la désignation des cas suspects : la méthode « mahal » s’appuie sur la distance multidimensionnelle de Mahalanobis calculée dans l’espace résiduel bidimensionnel du graphique, tandis que la méthode « identify » permet une sélection interactive à la souris directement sur l’écran d’affichage, offrant à l’utilisateur une immersion exploratoire immédiate dans ses données.

9.2 Distinction cruciale : Levantage vs Influence vs Outlier

Pour exploiter pleinement la portée diagnostique d’un graphique de variable ajoutée, le modélisateur doit impérativement maîtriser la distinction conceptuelle et géométrique unissant trois notions statistiques souvent amalgamées : le point aberrant (outlier), le levier statistique (leverage) et l’observation influente (influential point).

Un point aberrant résiduel se manifeste dans l’AVP par une divergence verticale majeure le long de l’axe des ordonnées. Cette observation possède une valeur de réponse résiduelle e(Y|others) anormalement éloignée de la tendance générale, sans pour autant présenter nécessairement une valeur extrême sur le prédicteur examiné. Si ce point se situe au centre de l’axe des abscisses (près de zéro), il augmentera l’erreur quadratique moyenne et gonflera artificiellement l’erreur-type du coefficient partiel, diminuant la puissance du test statistique sans pour autant modifier substantiellement l’angle d’inclinaison de la droite.

À l’inverse, un point à fort levier partiel se positionne de façon isolée aux confins extrêmes de l’axe horizontal, c’est-à-dire loin vers la droite ou vers la gauche sur l’axe e(Xk|others). Il correspond à un individu dont la combinaison de covariables rend la valeur de Xk extraordinairement inhabituelle au regard du reste de l’échantillon. Si ce point de fort levier s’aligne fidèlement sur la trajectoire rectiligne dessinée par les autres observations, il stabilise le modèle sans introduire de distorsion.

En revanche, le péril méthodologique suprême survient lorsqu’un point conjugue simultanément un fort levier partiel et un résidu partiel disjoint. Situé dans un coin marginal du graphique hors de l’alignement général, ce point agit comme un point d’appui mécanique qui exerce un véritable pouvoir de bascule (tilt) sur la droite des moindres carrés. Une seule observation de cette nature possède le pouvoir délétère de renverser le signe d’un coefficient ou de fabriquer de toutes pièces une significativité statistique totalement factice.

9.3 Mesures correctives et analyses de sensibilité

La découverte visuelle d’un cas suspect étiqueté sur un graphique de variable ajoutée ne doit en aucun cas conduire à sa suppression aveugle et irréfléchie du fichier de données. L’analyste doit en premier lieu procéder à une vérification minutieuse de l’intégrité de la saisie numérique : s’agit-il d’une erreur d’encodage typographique ou d’une mauvaise affectation d’unité de mesure ? Si la valeur s’avère parfaitement valide sur le plan empirique, elle incarne une réalité observationnelle légitime qu’il convient de comprendre et non de censurer arbitrairement.

La stratégie analytique recommandée réside dans la conduite systématique d’une analyse de sensibilité documentée. Le chercheur ajuste le modèle linéaire sur l’échantillon complet, puis ré-estime ce même modèle en excluant formellement la ou les observations étiquetées comme hautement influentes dans les AVP. Par exemple, au sein du jeu de données Duncan, les professions de ministre du culte (minister) et de chef de train (conductor) sont notoirement connues pour exercer une influence démesurée sur les coefficients partiels de l’éducation et du revenu.

En comparant rigoureusement les sorties paramétriques avec et sans ces cas critiques, l’analyste évalue la stabilité structurelle de ses inférences. Si les conclusions scientifiques se trouvent radicalement modifiées par le retrait d’une unique unité statistique, la plus grande réserve s’impose dans l’interprétation théorique. Dans de telles circonstances, l’analyste délaissera les moindres carrés ordinaires au profit d’approches de régression robuste, telles que la fonction rlm() de la bibliothèque MASS, qui procède à une pondération adaptative décroissante des observations en fonction de l’amplitude de leurs résidus.

10. Personnalisation esthétique et publication scientifique avec ggplot2

10.1 Extraction des données résiduelles pour ggplot2

Bien que les graphiques générés par les fonctions du moteur graphique de base de R via car::avPlots() soient techniquement irréprochables sur le plan du diagnostic exploratoire, leurs qualités esthétiques demeurent limitées et s’accordent difficilement avec les exigences typographiques et stylistiques des publications scientifiques internationales contemporaines. L’écosystème graphique unifié centré sur la grammaire des graphiques et implémenté par la prestigieuse bibliothèque ggplot2 développée par Hadley Wickham s’impose aujourd’hui comme la norme graphique de référence.

Pour transposer l’architecture d’un AVP au sein de cet univers visuel, il est nécessaire d’extraire préalablement les résidus partiels sous la forme d’une structure rectangulaire ordonnée (tidy data frame ou tibble). En nous appuyant sur la syntaxe matricielle ou sur des boucles automatisées, nous pouvons itérer à travers l’ensemble des régresseurs du modèle pour calculer les double résidualisations respectives.

Une démarche élégante consiste à structurer ces résultats au format long. Le tableau résultant comportera trois colonnes primordiales : une variable catégorielle désignant le prédicteur considéré (Predictor), une colonne numérique regroupant les résidus horizontaux correspondants (Partial_Predictor_Residuals) et une colonne numérique consignant les résidus de la variable réponse (Partial_Response_Residuals). L’adoption de ce format long conditionne la capacité du moteur graphique à générer des visualisations facettées d’une rigueur formelle absolue.

10.2 Construction d’un graphique moderne multi-facettes

Disposant d’un jeu de données résiduelles correctement formaté, nous pouvons initialiser la construction d’une visualisation multi-facettes combinant l’ensemble des prédicteurs au sein d’une composition graphique harmonisée. L’instruction fondamentale s’articule autour de l’appel ggplot(df_residus, aes(x = Partial_Predictor_Residuals, y = Partial_Response_Residuals)), enrichi immédiatement de la couche géométrique des marqueurs ponctuels via geom_point().

Afin de valoriser la lisibilité des distributions à forte densité, il est vivement judicieux de moduler l’opacité des points par le truchement du paramètre alpha (par exemple alpha = 0.6) tout en calibrant délicatement leur taille par size = 2. L’étape cruciale consiste ensuite à superposer simultanément deux trajectoires modélisées distinctes : d’une part, la droite de régression linéaire stricte via geom_smooth(method = « lm », color = « darkblue », se = TRUE), affichant son intervalle de confiance à quatre-vingt-quinze pour cent ; et d’autre part, une courbe de lissage local non paramétrique de type Loess via geom_smooth(method = « loess », color = « firebrick », linetype = « dashed », se = FALSE).

L’adjonction de l’opérateur facet_wrap(~ Predictor, scales = « free_x ») constitue l’aboutissement de cette architecture. L’argument scales = « free_x » est ici capital, dans la mesure où les échelles d’abscisses reflètent les résidus de prédicteurs distincts s’exprimant dans des unités métriques hétérogènes (par exemple des pourcentages de revenus confrontés à des années d’études), tandis que l’axe des ordonnées demeure strictement invariant et synchronisé.

10.3 Formatage typographique selon les normes de publication

L’optimisation finale de la planche visuelle destinée à l’insertion dans un manuscrit scientifique requiert l’épuration totale des éléments d’arrière-plan superflus. L’application d’un thème minimaliste tel que theme_bw() ou theme_classic() permet de supprimer le fond grisâtre standard et d’instaurer des lignes d’axes noires d’une netteté parfaite. La typographie textuelle des étiquettes d’axes, des titres de facettes et des légendes sera minutieusement calibrée à l’aide des fonctions de personnalisation élémentaire element_text(), en optant pour des polices vectorielles neutres à l’image d’Arial ou d’Helvetica.

Pour assurer un étiquetage textuel irréprochable des points atypiques sans risquer d’enchevêtrement ou de superposition illisible des étiquettes sur les marqueurs, l’intégration de l’extension spécialisée ggrepel constitue une solution de choix. Grâce à la géométrie geom_text_repel(), les identifiants des observations suspectes se déploient harmonieusement autour des résidus extrêmes en s’écartant dynamiquement des points denses par un algorithme de répulsion physique.

Enfin, l’exportation de l’objet graphique final s’exécute par l’entremise de la fonction vectorielle ggsave(). L’analyste prendra soin de consigner la planche au format PDF pour une conservation vectorielle intégrale sans tramage matriciel, ou au format TIFF haute résolution à une densité minimale de 300 à 600 points par pouce (DPI), en stipulant formellement les dimensions physiques en centimètres requises par les directives éditoriales de la revue académique ciblée.

11. Diagnostics différentiels et comparaison graphique

11.1 AVP versus Graphiques de résidus partiels (Component-plus-Residual)

L’une des confusions conceptuelles les plus persistantes au sein de la littérature statistique appliquée concerne la distinction primordiale entre les graphiques de variables ajoutées (AVP) et les graphiques de résidus partiels conventionnels, désignés couramment sous le vocable de graphiques composante-plus-résidu (Component-plus-Residual plots ou CR plots), accessibles dans R via la fonction crPlots() du package car.

Bien que ces deux outils partagent la vocation générale de sonder les relations conditionnelles, leurs formules de calcul diffèrent fondamentalement. Dans un graphique CR, l’axe vertical ne trace pas les résidus de la réponse, mais la quantité construite Beta_k * X_ik + e_i, c’est-à-dire la composante linéaire estimée du prédicteur k augmentée du vecteur des résidus complets du modèle, tandis que l’axe horizontal conserve la valeur brute du prédicteur X_ik sans la moindre résidualisation. En conséquence, les graphiques CR s’avèrent infiniment supérieurs et spécifiquement optimisés pour révéler la forme fonctionnelle précise (courbure, non-linéarité) que prend la relation le long de l’échelle d’origine de la variable indépendante.

En contrepartie absolue, les graphiques CR sont structurellement inaptes à diagnostiquer l’effet de levier et l’influence des points. En effet, l’axe horizontal d’un CR plot ne tenant aucun compte des intercorrélations unissant Xk aux autres prédicteurs, un point peut paraître tout à fait banal sur l’échelle brute de Xk tout en constituant une anomalie multidimensionnelle colossale dans l’espace des covariables. Les graphiques de variables ajoutées (AVP) demeurent ainsi les seuls et uniques instruments graphiques capables d’isoler simultanément le levier partiel réel, l’influence par bascule et l’impact marginal net sur le coefficient de régression.

11.2 AVP versus Graphiques d’effets marginaux (marginal effects)

Une autre distinction cruciale doit être opérée vis-à-vis des graphiques d’effets marginaux ou de prédictions marginales ajustées, devenus extraordinairement populaires dans les sciences sociales grâce à des modules R modernes tels que ggeffects ou marginaleffects. Les graphiques d’effets marginaux ont pour mission exclusive de traduire les prédictions d’un modèle ajusté sous une forme substantielle intuitive. Ils projettent la valeur attendue de la variable réponse (par exemple une probabilité prédite ou un score théorique) en faisant varier un prédicteur d’intérêt tout en figeant les autres covariables à des valeurs constantes de référence, telles que leur moyenne ou leur médiane.

De fait, les tracés d’effets marginaux illustrent des prédictions idéalisées issues du modèle, et non des résidus observés. Ils appartiennent au registre de la communication des résultats finaux et de la vulgarisation théorique des effets statistiques. En revanche, les graphiques de variables ajoutées appartiennent rigoureusement au registre du diagnostic d’invalidation du modèle. L’AVP met en scène les données brutes résidualisées contre la modélisation, exposant impitoyablement les failles de l’ajustement, la dispersion non expliquée et les observations qui déstabilisent l’estimation.

En résumé méthodologique, le chercheur prudent mobilisera les AVP en amont lors de sa phase de laboratoire afin d’éprouver la robustesse de sa spécification et de traquer les anomalies distributionnelles, puis emploiera les graphiques d’effets marginaux en aval au sein du corps de son article pour communiquer de manière pédagogique l’ampleur sociologique ou clinique des associations mises en évidence.

11.3 AVP versus Nuages de points bivariés bruts

La confrontation directe entre un nuage de points bivarié standard (Y en fonction de Xk) et un graphique de variable ajoutée constitue sans doute la démonstration la plus spectaculaire du paradoxe de Simpson et des illusions d’optique générées par la confusion statistique en milieu observationnel. Dans de multiples architectures de données, la trajectoire bivariée brute observée entre deux grandeurs s’avère diamétralement opposée à la réalité de leur lien conditionnel.

Imaginons une étude où le temps passé sur une tâche cognitive (X) est positivement corrélé au nombre d’erreurs commises (Y) à travers un nuage bivarié brut ascendant bien visible. Un observateur superficiel conclurait que la délibération prolongée nuit à la performance cognitive. Cependant, l’introduction au sein du modèle d’une covariable d’ajustement fondamentale, telle que la difficulté intrinsèque de l’exercice proposé, renverse totalement la perspective.

En projetant le graphique de variable ajoutée reliant les résidus du temps aux résidus des erreurs après extraction de l’effet de la difficulté, la droite des moindres carrés affiche désormais une pente fortement négative. À complexité de problème rigoureusement égale, allouer davantage de temps de réflexion réduit significativement le taux d’erreur. L’AVP a littéralement opéré une purge géométrique de la variable de confusion masquée, substituant à une corrélation bivariée trompeuse la vérité de l’association conditionnelle partielle.

12. Bonnes pratiques, pièges courants et recommandations de reporting

12.1 Erreurs fréquentes d’interprétation à éviter

L’utilisation des graphiques de variables ajoutées est fréquemment entachée d’erreurs d’interprétation récurrentes chez les praticiens non statisticiens. Le premier écueil réside dans la tentative erronée d’interpréter les valeurs numériques lues sur l’axe des abscisses comme s’il s’agissait des mesures réelles de la variable indépendante. Une valeur de zéro sur l’abscisse d’un AVP ne correspond aucunement à une absence de la caractéristique physique chez le sujet, mais signifie simplement que l’observation se situe exactement sur la valeur moyenne prédite pour cette variable par les autres prédicteurs.

Une seconde méprise classique consiste à tenter d’assigner une signification sociologique ou théorique à l’ordonnée à l’origine du tracé partiel. Comme démontré dans les fondements algébriques du théorème de Frisch-Waugh-Lovell, la droite d’un AVP passe nécessairement et immuablement par le barycentre théorique c(0,0). Vouloir interpréter l’intercept d’une telle régression résiduelle n’a absolument aucun sens méthodologique.

Enfin, une vigilance extrême s’impose quant à l’impact de la multicolinéarité sévère sur la géométrie même de l’AVP. Lorsqu’une variable prédictive se trouve presque entièrement expliquée par une combinaison linéaire des autres régresseurs (associée à un facteur d’inflation de la variance, ou VIF, supérieur à 10), la variance résiduelle e(Xk|others) s’effondre dramatiquement. Sur le plan visuel, le nuage de points se trouve alors violemment compressé le long d’une bande verticale étroite centrée sur zéro. Cette compression artificielle réduit le bras de levier des observations, instabilise de façon spectaculaire la pente de la droite et expose l’estimation paramétrique à une volatilité d’échantillonnage extrême.

12.2 Gestion des variables catégorielles et termes d’interaction

Lorsque le modèle statistique intègre des variables qualitatives nominales encodées sous la forme de variables indicatrices binaires (variables muettes ou dummies), l’interprétation des sorties de avPlots() requiert une grille de lecture adaptée. Pour un facteur binaire découpé en une indicatrice 0/1, l’axe horizontal de l’AVP présentera deux sous-grappes distinctes de résidus, reflétant les écarts résiduels des deux groupes d’individus par rapport à leur propension conditionnelle moyenne.

Dans ce contexte, la pente de la droite continue de matérialiser fidèlement l’écart moyen conditionnel ajusté séparant la catégorie cible de la catégorie de référence. Néanmoins, pour des facteurs polytomiques comprenant un grand nombre de modalités (par exemple une variable de région géographique découpée en cinq indicatrices distinctes), la fragmentation du facteur en cinq graphiques partiels disjoints complique l’appréciation globale de l’effet du construit. L’analyste se tournera alors vers des tests d’hypothèses linéaires collectives de type Wald ou analyses de variance multivariées d’Anova de type II ou III.

Concernant les termes d’interaction multiplicative (modération statistique), la prudence s’accentue. Un terme d’interaction noté X1:X2 générera son propre panneau avPlot dédié. Les résidus de ce terme représentent la fraction de l’effet modérateur qui ne peut être prédite ni par les effets principaux individuels, ni par les autres régresseurs du système. La lecture de ce graphique permet de vérifier visuellement si l’effet de modération proclamé n’est pas le fruit exclusif d’une poignée de données atypiques logées dans les zones marginales du produit croisé des deux distributions.

12.3 Normes de rédaction et communication des résultats académiques

La reddition de comptes relative au diagnostic par variables ajoutées au sein d’un article scientifique évalué par les pairs doit obéir à des standards de rédaction formels d’une grande rigueur. Il est proscrit d’affirmer sommairement dans le texte que les hypothèses de linéarité et d’influence ont été vérifiées sans expliciter le protocole diagnostique concret mis en œuvre par l’équipe de recherche.

Dans la section méthodologique du manuscrit, les auteurs spécifieront clairement le recours aux graphiques de régression partielle : « L’absence de distorsion non linéaire et la stabilité des coefficients partiels face aux observations influentes ont été formellement diagnostiquées par l’inspection systématique des graphiques de variables ajoutées (Added Variable Plots) générés via le package car sous R. » En cas de détection d’observations à fort levier ayant requis une analyse de sensibilité, les conclusions comparatives doivent être obligatoirement consignées sous la forme d’un tableau comparatif des modèles inséré dans le corps du texte ou versé au sein des annexes électroniques (Supplementary Materials).

Enfin, dans l’esprit du mouvement pour la science ouverte et la reproductibilité computationnelle promu par les directives du Center for Open Science, les chercheurs fourniront en libre accès sur des plateformes pérennes (telles que l’Open Science Framework ou GitHub) l’intégralité du script R ayant présidé à la confection des graphiques, assorti des instructions de compilation permettant à n’importe quel réviseur ou lecteur de régénérer à l’identique l’ensemble des tracés résiduels présentés.

Références

• Belsley, D. A., Kuh, E., & Welsch, R. E. (1980). Regression diagnostics: Identifying influential data and sources of collinearity. John Wiley & Sons. https://doi.org/10.1002/0471725153
• Cook, R. D., & Weisberg, S. (1982). Residuals and influence in regression. Chapman & Hall.
• Fox, J., & Weisberg, S. (2019). An R companion to applied regression (3rd ed.). SAGE Publications. https://socialsciences.mcmaster.ca/jfox/Books/Companion/
• Frisch, R., & Waugh, F. V. (1933). Partial time regressions as compared with individual trends. Econometrica, 1(4), 387–401. https://doi.org/10.2307/1907330
• Long, J. S., & Ervin, L. H. (2000). Using heteroscedasticity consistent standard errors in the linear regression model. The American Statistician, 54(3), 217–224. https://doi.org/10.1080/00031305.2000.10474549
• Lovell, M. C. (1963). Seasonal adjustment and other data transformations: An analysis of alternative methods. Econometrica, 31(4), 993–1010.
• R Core Team. (2023). R: A language and environment for statistical computing. R Foundation for Statistical Computing. https://www.R-project.org/
• Wickham, H. (2016). ggplot2: Elegant graphics for data analysis. Springer-Verlag. https://ggplot2.tidyverse.org/