Régression linéaire multiple

Lorsqu'on dispose de n observations (y_i, x_i1, ..., x_ip), i ∈ {1, n} , qui sont des réalisations des variables aléatoires (Y_i, X_i1, ..., X_ip), l'équation de régression s'écrit

${\displaystyle y_{i}=a_{0}+a_{1}x_{i1}+\ldots +a_{p}x_{ip}+\varepsilon _{i}\qquad i=1,\ldots ,n\,}$

La problématique reste la même que pour la régression simple :

• estimer les paramètres a₀, a₁, ..., a_p en exploitant les observations ;
• évaluer la précision de ces estimateurs ;
• mesurer le pouvoir explicatif du modèle ;
• évaluer l'influence des variables dans le modèle :
  • globalement (les p variables en bloc) et,
  • individuellement (chaque variable) ;
• évaluer la qualité du modèle lors de la prédiction (intervalle de prédiction) ;
• détecter les observations qui peuvent influencer exagérément les résultats (points atypiques).

On peut adopter une écriture condensée qui rend la lecture et la manipulation de l'ensemble plus facile. Les équations suivantes

${\displaystyle {\begin{cases}y_{1}=a_{0}+a_{1}x_{1,1}+\ldots +a_{p}x_{1,p}+\varepsilon _{1}\\y_{2}=a_{0}+a_{1}x_{2,1}+\ldots +a_{p}x_{2,p}+\varepsilon _{2}\\\cdots \\y_{n}=a_{0}+a_{1}x_{n,1}+\ldots +a_{p}x_{n,p}+\varepsilon _{n}\end{cases}}}$

peuvent être résumées avec la notation matricielle

${\displaystyle {\begin{pmatrix}y_{1}\\\vdots \\y_{n}\end{pmatrix}}={\begin{pmatrix}1&x_{1,1}&\cdots &x_{1,p}\\\vdots &\vdots &\ddots &\vdots \\1&x_{n,1}&\cdots &x_{n,p}\end{pmatrix}}{\begin{pmatrix}a_{0}\\a_{1}\\\vdots \\a_{p}\\\end{pmatrix}}+{\begin{pmatrix}\varepsilon _{1}\\\vdots \\\varepsilon _{n}\\\end{pmatrix}}}$

Soit de manière compacte :

${\displaystyle y=Xa+\varepsilon \,}$

avec

• y est de dimension (n, 1)
• X est de dimension (n, p+1)
• a est de dimension (p+1, 1)
• ε est de dimension (n, 1)

La première colonne de la matrice X sert à indiquer que la régression est effectuée avec constante (ici ${\displaystyle a_{0}}$).

Comme en régression simple, les hypothèses permettent de déterminer : les propriétés des estimateurs (biais, convergence) ; et leurs lois de distributions (pour les estimations par intervalle et les tests d'hypothèses).

Il existe principalement deux catégories d'hypothèses :

Hypothèses stochastiques

• H₁ : Les X_j sont déterminées sans erreurs, j = 1, …, p ;
• H₂ : ${\displaystyle \mathbb {E} (\varepsilon _{i})=0\,}$ Le modèle est bien spécifié en moyenne ;
• H₃ : ${\displaystyle {\text{Var}}(\varepsilon _{i})=\sigma ^{2}\ \forall {i}\,}$ Homoscédasticité des erreurs (variance constante)
• H₄ : ${\displaystyle \mathrm {cov} (\varepsilon _{i},\varepsilon _{j})=0\ \forall {i\neq j}\,}$ Pas d'autocorrélation des erreurs.
• H₅ : ${\displaystyle \mathrm {cov} (X_{i},\varepsilon _{j})=0\ \forall {i\neq j}\,}$ Les erreurs sont linéairement indépendantes des variables exogènes.
• H₆ : ${\displaystyle \varepsilon \sim {\mathcal {N}}_{n}(0,\sigma ^{2}I_{n})\,}$ Les erreurs suivent une loi normale multidimensionnelle (H₆ implique les hypothèses H₂, H₃ et H₄, la réciproque étant fausse car les trois hypothèses réunies n'impliquent pas que ε soit un vecteur gaussien).

Hypothèses structurelles

• H₇ : absence de colinéarité entre les variables explicatives, i.e. X ^TX est régulière, det(X ^TX) ≠ 0 et (X ^TX)⁻¹ existe (remarque : c'est équivalent à rang(X) = rang(X ^TX) = p + 1) ;
• H₈ : ${\displaystyle {\frac {1}{n}}X^{T}X}$ tend vers une matrice finie non singulière Q lorsque n → +∞ ;
• H₉ : ${\displaystyle n>p+1\,}$ Le nombre d'observations est supérieur au nombre de variables + 1 (la constante). S'il y avait égalité, le nombre d'équations serait égal au nombre d'inconnues a_j, la droite de régression passerait par tous les points, nous serions face à un problème d'interpolation linéaire (voir Interpolation numérique).

Écriture matricielle de l'hypothèse H₆

${\displaystyle \mathrm {H_{2}:} \mathbb {E} (\varepsilon )=\mathbb {E} {\begin{pmatrix}\varepsilon _{1}\\\vdots \\\varepsilon _{n}\end{pmatrix}}={\begin{pmatrix}0\\\vdots \\0\end{pmatrix}}}$

Sous l'hypothèse d'homoscedasticité et d'absence d'auto-corrélation, la matrice de variance-covariance du vecteur des erreurs peut s'écrire :

${\displaystyle \mathrm {H_{3}\ {\mbox{et}}\ H_{4}:} \ \mathrm {cov} (\varepsilon )=\sigma ^{2}I_{n}=\sigma ^{2}{\begin{pmatrix}1&0&\cdots &0\\0&1&\cdots &0\\\vdots &&\ddots &\vdots \\0&\cdots &\cdots &1\end{pmatrix}}={\begin{pmatrix}\sigma ^{2}&0&\cdots &0\\0&\sigma ^{2}&\cdots &0\\\vdots &&\ddots &\vdots \\0&\cdots &\cdots &\sigma ^{2}\end{pmatrix}}}$

Dans certains cas, l'hypothèse (H₁) est intenable : les régresseurs X sont supposés aléatoires. Mais dans ce cas, on suppose que X est aléatoire mais est indépendant de l'aléa ε. On remplace alors l'hypothèse (H₂) par une hypothèse sur l'espérance conditionnelle :

${\displaystyle \mathrm {H_{2}:} \mathbb {E} (\varepsilon _{i}\mid X)=0\,}$

De même, il faudrait changer en conséquence les hypothèses (H₃), (H₄) et aussi (H₅).

Du modèle complet :

${\displaystyle y_{i}=a_{0}+a_{1}x_{i,1}+\cdots +a_{p}x_{i,p}+\epsilon _{i}\,}$

On va estimer les paramètres et obtenir:

${\displaystyle {\hat {y_{i}}}={\hat {a}}_{0}+{\hat {a}}_{1}x_{i,1}+\cdots +{\hat {a}}_{p}{x}_{i,p}\,}$

Les résidus estimés sont la différence entre la valeur de y observée et estimée. Soit :

Le principe des moindres carrés consiste à rechercher les valeurs des paramètres qui minimisent la somme des carrés des résidus.

${\displaystyle \min \sum _{i=1}^{n}{\hat {\epsilon }}_{i}^{2}=\min _{{\hat {a}}_{0},.,{\hat {a}}_{p}}\sum _{i=1}^{n}(y_{i}-{\hat {a}}_{0}-{\hat {a}}_{1}x_{i,1}-\cdots -{\hat {a}}_{p}x_{i,p})^{2}}$.

Ce qui revient à rechercher les solutions de ${\displaystyle {\frac {\partial (\sum {\hat {\epsilon }}_{i}^{2})}{\partial {\hat {a}}_{j}}}=0\,}$. Nous avons j =p + 1 équations, dites équations normales, à résoudre.

La solution obtenue est l'estimateur des moindres carrés ordinaires, il s'écrit :

avec X^T la transposée de X

Remarques:

• Pourquoi minimiser la somme des carrés plutôt que la simple somme? Cela tient, en partie, au fait que la moyenne de ces résidus sera 0, et donc que nous disposerons de résidus positifs et négatifs. Une simple somme les annulerait, ce qui n'est pas le cas avec les carrés.
• si les x_j sont centrés, 1/n (X ^TX) correspond à la matrice de variance-covariance des variables exogènes ; s'ils sont centrés et réduits, 1/n (X ^TX) correspond à la matrice de corrélation.

Interprétation géométrique, algébrique et statistique de l'estimateur MCO (Moindres Carrés Ordinaires)

• L'estimateur MCO correspond à une projection orthogonale du vecteur Y sur l'espace formé par les vecteurs X.
• L'estimateur MCO correspond à une matrice inverse généralisée du système Y = Xa pour mettre a en évidence. En effet, si on multiplie à gauche par l'inverse généralisée ${\displaystyle (X^{T}X)^{-1}X^{T}}$ on a :

${\displaystyle (X^{T}X)^{-1}X^{T}Y=(X^{T}X)^{-1}X^{T}Xa=a}$

• L'estimateur MCO est identique à l'estimateur obtenu par le principe du maximum de vraisemblance.

Si les hypothèses initiales sont respectées, l'estimateur des MCO possède d'excellentes propriétés.

Cette propriété se base seulement sur les hypothèses d'espérance nulle des résidus. La présence d'autocorrélation ou d'hétéroscédasticité n'affecte pas ce résultat.

Ceci signifie qu'il n'existe pas d'estimateur linéaire sans biais de a qui ait une variance plus petite. Cette propriété en anglais est désignée par BLUE, pour best linear unbiased estimator. La preuve est donnée par le théorème de Gauss-Markov.

Ce résultat est obtenu sans l'hypothèse de normalité des résidus (H₆).

La matrice de variance-covariance des coefficients est importante car elle renseigne sur la variance de chaque coefficient estimé, et permet de faire des tests d'hypothèse, notamment de voir si chaque coefficient est significativement différent de zéro. Elle est définie par :

${\displaystyle \operatorname {Var} ({\hat {a}})\equiv \Sigma =\operatorname {E} [({\hat {a}}-a)({\hat {a}}-a)^{T}]}$

Sous les hypothèses d'espérance nulle, d'absence d'autocorrélation et d'homoscédasticité des résidus (H₁ à H₅), on a :

${\displaystyle \operatorname {Var} ({\hat {a}})=\sigma ^{2}(X'X)^{-1}}$

Cette formule ne s'applique cependant que dans le cas où les résidus sont homoscédastiques et sans auto-corrélation, ce qui permet d'écrire la matrice des erreurs comme :

${\displaystyle {\textrm {Cov}}[\varepsilon ]=\sigma ^{2}I_{n}\,}$

S'il y a de l'hétéroscédasticité ou de l'auto-corrélation, et donc ${\displaystyle {\textrm {Cov}}[\varepsilon ]\neq \sigma ^{2}I_{n}}$, il est possible de rectifier la matrice de variance-covariance estimée par :

• la matrice de variance-covariance de White (ou Eicker-White (1967, 1980)), consistante en cas d'hétéroscédasticité (en anglais HC pour Heteroskedasticity Consistent).
• la matrice de variance-covariance de Newey-West (1987), consistante en cas d'hétéroscédasticité et d'auto-corrélation (en anglais HAC pour Heteroskedasticity and Autocorrelation Consistent).

Pour la variance du résidu ${\displaystyle \sigma ^{2}\equiv \operatorname {Var} [\varepsilon ]}$, on peut utiliser l'estimateur sans biais construit à partir de la variance des résidus observés :

${\displaystyle s^{2}\equiv {\hat {\sigma }}^{2}={\frac {1}{n-p-1}}\sum _{i=1}^{N}{\hat {\varepsilon }}_{i}^{2}}$

Les ${\displaystyle {\hat {\varepsilon }}}$ correspondent aux résidus observés : ${\displaystyle {\hat {\varepsilon }}=Y-{\hat {Y}}}$.

On remarque deux choses par rapport à l'estimateur classique de la variance :

${\displaystyle s_{n-1}^{2}\equiv {\hat {\sigma }}^{2}={\frac {1}{n-1}}\sum _{i=1}^{n}\left(y_{i}-{\overline {y}}\right)^{2}}$,

• on n'inclut pas l'espérance des résidus, car celle-ci est supposée être de zéro (selon ${\displaystyle H_{2}}$). Surtout, les résidus du modèle ont exactement une moyenne de zéro lorsqu'une constante est introduite dans le modèle.
• La somme des carrés est divisée par n - p - 1 = n - (p + 1) et non par n-1. En fait, n-p-1 correspond aux degrés de liberté du modèle (le nombre d'observations moins le nombre de coefficients à estimer). On remarque effectivement que ${\displaystyle \operatorname {E} ({\hat {\varepsilon }}'{\hat {\varepsilon }})=\sigma ^{2}(n-p-1)}$.

Il existe également un autre estimateur, obtenu par la méthode du maximum de vraisemblance, qui est cependant biaisé :

${\displaystyle s^{2}\equiv {\hat {\sigma }}^{2}={\frac {1}{n}}\sum _{i=1}^{N}{\hat {\varepsilon }}_{i}^{2}}$

Il suffit de remplacer la variance théorique des résidus, σ², par son estimateur sans biais des moindres carrés : ${\displaystyle s^{2}\equiv {\hat {\sigma }}^{2}={\frac {1}{n-p-1}}\sum _{i=1}^{N}{\hat {\varepsilon }}_{i}^{2}}$

L'estimateur de la matrice de variance-covariance des résidus devient :

${\displaystyle {\widehat {\operatorname {Var} }}[{\hat {a}}]\equiv {\hat {\Sigma }}_{\hat {a}}={\hat {\sigma }}^{2}(X^{T}X)^{-1}}$

La variance estimée ${\displaystyle {\hat {\sigma }}_{{\hat {a}}_{j}}^{2}}$ de l'estimation du paramètre â_j est lue sur la diagonale principale de cette matrice.

Après avoir obtenu l'estimateur, son espérance et une estimation de sa variance, il ne reste plus qu'à calculer sa loi de distribution pour produire une estimation par intervalle et réaliser des tests d'hypothèses.

En partant de l'hypothèse

${\displaystyle \epsilon _{i}\sim {\mathcal {N}}(0,\sigma ^{2})\,}$,

on peut montrer

1. ${\displaystyle {\frac {{\hat {a}}_{j}-a_{j}}{\sigma _{{\hat {a}}_{j}}}}\sim {\mathcal {N}}(0,1)}$
2. ${\displaystyle (n-p-1){\frac {{\hat {\sigma }}_{{\hat {a}}_{j}}^{2}}{\sigma _{{\hat {a}}_{j}}^{2}}}\sim \chi _{n-p-1}^{2}}$

Le rapport d'une loi normale et de la racine carrée d'une loi du χ² normalisée par ses degrés de liberté aboutit à une loi de Student. On en déduit donc la statistique :

${\displaystyle t={\frac {{\hat {a}}_{j}-a_{j}}{{\hat {\sigma }}_{{\hat {a}}_{j}}}}\sim \mathrm {T} (n-p-1)}$

elle suit une loi de Student à (n - p - 1) degrés de liberté.

À partir de ces informations, il est possible de calculer les intervalles de confiance des estimations des coefficients.

Il est également possible de procéder à des tests d'hypothèses, notamment les tests d'hypothèses de conformité à un standard. Parmi les différents tests possibles, le test de nullité du coefficient (H₀ : a_j = 0, contre H₁ : a_j ≠ 0) tient un rôle particulier : il permet de déterminer si la variable x_j joue un rôle significatif dans le modèle. Il faut néanmoins être prudent quant à ce test. L'acceptation de l'hypothèse nulle peut effectivement indiquer une absence de corrélation entre la variable incriminée et la variable endogène ; mais il peut également résulter de la forte corrélation de x_j avec une autre variable exogène, son rôle est masqué dans ce cas, laissant à croire une absence d'explication de la part de la variable.

L'évaluation globale de la pertinence du modèle de prédiction s'appuie sur l'équation d'analyse de variance SCT = SCE + SCR, où

• SCT, somme des carrés totaux, traduit la variabilité totale de l'endogène ;
• SCE, somme des carrés expliqués, traduit la variabilité expliquée par le modèle ;
• SCR, somme des carrés résiduels correspond à la variabilité non-expliquée par le modèle.

Toutes ces informations sont résumées dans un tableau, le tableau d'analyse de variance.

Source de variation	Somme des carrés	Degrés de liberté	Carrés moyens
Expliquée	${\displaystyle SCE=\sum _{i}({\hat {y}}_{i}-{\bar {y}})^{2}}$	p	${\displaystyle CME={\frac {SCE}{p}}}$
Résiduelle	${\displaystyle SCR=\sum _{i}(y_{i}-{\hat {y}}_{i})^{2}}$	n – p – 1	${\displaystyle CMR={\frac {SCR}{n-p-1}}}$
Totale	${\displaystyle SCT=\sum _{i}(y_{i}-{\bar {y}})^{2}}$	n – 1

Dans le meilleur des cas, SCR = 0, le modèle arrive à prédire exactement toutes les valeurs de y à partir des valeurs des x_j. Dans le pire des cas, SCE = 0, le meilleur prédicteur de y est sa moyenne y.

Un indicateur spécifique permet de traduire la variance expliquée par le modèle, il s'agit du coefficient de détermination. Sa formule est la suivante :

${\displaystyle R^{2}={\frac {SCR}{SCT}}=1-{\frac {SCE}{SCT}}\,}$

${\displaystyle R={\sqrt {R^{2}}}\,}$ est le coefficient de corrélation multiple.

Dans une régression avec constante, nous avons forcément

${\displaystyle 0\leqslant R^{2}\leqslant 1}$.

Enfin, si le R ² est certes un indicateur pertinent, il présente un défaut parfois ennuyeux, il a tendance à mécaniquement augmenter à mesure que l'on ajoute des variables dans le modèle. De ce fait, il est inopérant si l'on veut comparer des modèles comportant un nombre différent de variables. Il est conseillé dans ce cas d'utiliser le coefficient de détermination ajusté qui est corrigé des degrés de libertés. Le R ² ajusté est toujours inférieur au R ².

Le R ² est un indicateur simple, on comprend aisément que plus il s'approche de la valeur 1, plus le modèle est intéressant. En revanche, il ne permet pas de savoir si le modèle est statistiquement pertinent pour expliquer les valeurs de y.

Nous devons nous tourner vers les tests d'hypothèses pour vérifier si la liaison mise en évidence avec la régression n'est pas un simple artefact.

La formulation du test d'hypothèse qui permet d'évaluer globalement le modèle est la suivante :

• H₀ : a₁ = a₂ = … = a_p = 0 ;
• H₁ : un des coefficients au moins est non nul.

La statistique dédiée à ce test s'appuie (parmi les différentes formulations possibles) sur le R ², il s'écrit :

${\displaystyle F_{calc}={\frac {\frac {R^{2}}{p}}{\frac {1-R^{2}}{n-p-1}}}}$,

et suit une loi de Fisher à (p, n - p - 1) degrés de liberté.

La région critique du test est donc : rejet de H₀ si et seulement si F_calc > F_{1 - α}(p, n - p - 1), où α est le risque de première espèce.

Une autre manière de lire le test est de comparer la p-value (probabilité critique du test) avec α : si elle est inférieure, l'hypothèse nulle est rejetée.

• Régis Bourbonnais, Économétrie, Dunod, 1998 (ISBN 2100038605)
• Yadolah Dodge et Valentin Rousson, Analyse de régression appliquée, Dunod, 2004 (ISBN 2100486594)
• R. Giraud, N. Chaix, Économétrie, Puf, 1994
• C. Labrousse, Introduction à l'économétrie -- Maîtrise d'économétrie, Dunod, 1983
• J. Confais, M. Le Guen, Premiers pas en régression linéaire, La Revue Modulad, N°35, 2006, pp220–363[1],

• Ajustement de profil
• Modèles de régression multiple postulés et non postulés
• Régression linéaire
• Régression polynomiale
• Régression elliptique