Moindres carrés ordinaires et régularisation

Pour une matrice de conception $X$, des cibles $y$ et un vecteur de coefficients $w$, les moindres carrés ordinaires résolvent

$$\hat{w}=\arg\min_w\|Xw-y\|_2^2.$$

La pseudo-inverse fournit une solution aux moindres carrés, $\hat{w}=X^{+}y$. Lorsque les colonnes de $X$ sont presque linéairement dépendantes, de faibles variations des données peuvent entraîner de fortes variations des coefficients. Les prédictions peuvent rester acceptables alors que les coefficients individuels deviennent instables.

Du résidu à la solution

$X$ comporte $n$ lignes (observations) et $p$ colonnes (coefficients), avec $y\in\mathbb{R}^n$ et $w\in\mathbb{R}^p$. Ajoutez une colonne de 1 pour ajuster une constante. Pour une droite $\hat y_i=b+ax_i$, le résidu $r_i=y_i-\hat y_i$ est un écart vertical dans l’unité de la cible, pas la distance minimale à la droite. Dans des coordonnées cartésiennes non redimensionnées, la distance perpendiculaire vaut $|r_i|/\sqrt{1+a^2}$ ; la minimiser résoudrait un autre problème, autorisant aussi un déplacement en $x$.

L’objectif $J(w)=\|Xw-y\|_2^2/(2n)$ conserve les mêmes minimiseurs OLS et donne

$$\nabla J(w)=\frac{1}{n}X^\top(Xw-y),\qquad \nabla^2J(w)=\frac{1}{n}X^\top X.$$

La Hessienne est semi-définie positive : un gradient nul correspond donc à un minimum global. L’annuler donne les équations normales :

$$X^\top X\hat w=X^\top y.$$

Si $\operatorname{rank}(X)=p$ (donc $n\geq p$), la solution est unique et vaut $(X^\top X)^{-1}X^\top y$. Sinon, ajouter un vecteur du noyau de $X$ aux coefficients ne change pas les valeurs ajustées : les coefficients ne sont pas uniques. La pseudo-inverse de Moore–Penrose choisit le vecteur de norme euclidienne minimale ; les valeurs ajustées sur l’entraînement restent uniques. Par exemple, des colonnes identiques n’identifient que la somme de leurs coefficients. Les prédictions peuvent diverger si de nouvelles entrées ne respectent plus cette relation entre colonnes.

Utilisez un solveur de moindres carrés fondé sur QR ou SVD plutôt qu’une inversion explicite ; former $X^\top X$ élève au carré le conditionnement spectral lorsque $X$ est de rang colonne plein. Le guide des modèles linéaires de scikit-learn présente OLS et sa résolution par SVD.

Trois observations à la main

Prenons $(x,y)=(0,1),(1,2),(2,2)$. La colonne constante et la colonne $x$ sont indépendantes. Avec $\bar x=1$ et $\bar y=5/3$,

$$a=\frac{\sum_i(x_i-\bar x)(y_i-\bar y)}{\sum_i(x_i-\bar x)^2} =\frac{1}{2},\qquad b=\bar y-a\bar x=\frac{7}{6}.$$

Les prédictions sont $(7/6,10/6,13/6)$, les résidus $(-1/6,1/3,-1/6)$ et leur somme de carrés vaut $1/6$. Les égalités $\sum_i r_i=0$ et $\sum_i x_i r_i=0$ vérifient les équations normales. La référence constante, moyenne d’entraînement $5/3$, obtient une somme de carrés de $2/3$. Cette comparaison porte sur l’entraînement, pas sur la généralisation.

from math import isclose
from statistics import mean

x, y = [0, 1, 2], [1, 2, 2]
xbar, ybar = mean(x), mean(y)
sxx = sum((v - xbar) ** 2 for v in x)
if sxx == 0:
    raise ValueError("A constant x cannot identify slope and intercept separately")
a = sum((u - xbar) * (v - ybar) for u, v in zip(x, y)) / sxx
b = ybar - a * xbar
r = [v - (b + a * u) for u, v in zip(x, y)]
assert isclose(a, 0.5) and isclose(b, 7 / 6)
assert isclose(sum(v * v for v in r), 1 / 6)
assert isclose(sum(r), 0, abs_tol=1e-12)
assert isclose(sum(u * v for u, v in zip(x, r)), 0, abs_tol=1e-12)

Ajustement et inférence statistique

Calculer le minimiseur ne nécessite aucune hypothèse gaussienne. Pour l’usage statistique, la documentation de régression statsmodels explicite le modèle de covariance des erreurs ; choisir une formule d’incertitude exige ce modèle supplémentaire. Sous un modèle de moyenne conditionnelle correctement spécifié $y=X\beta+\varepsilon$, avec $\mathbb{E}[\varepsilon\mid X]=0$ et un rang colonne plein, OLS est conditionnellement sans biais. Si, de plus, $\operatorname{Var}(\varepsilon\mid X)=\sigma^2 I$, alors

$$\operatorname{Var}(\hat w\mid X)=\sigma^2(X^\top X)^{-1}.$$

Pour $n>p$, la somme des carrés résiduels divisée par $n-p$ estime $\sigma^2$ sous ces hypothèses. Des erreurs conditionnellement gaussiennes justifient en outre l’inférence t/F exacte habituelle en échantillon fini. Des erreurs hétéroscédastiques ou corrélées nécessitent des estimateurs d’incertitude adaptés ; les erreurs-types robustes ne corrigent ni une mauvaise moyenne conditionnelle ni la confusion causale. Un intervalle pour la réponse moyenne diffère d’un intervalle de prédiction pour une nouvelle observation bruitée.

Les grands résidus peuvent dominer la perte quadratique, et les observations à fort levier (valeurs de variables inhabituelles) déplacer fortement la droite. Examinez les résidus selon les valeurs ajustées et le temps, puis comparez les erreurs sur données réservées à une référence avant d’interpréter les coefficients.

Objectifs pénalisés

Ci-dessous, $\lambda\geq0$ et $0\leq\rho\leq1$. Pour ces formules pénalisées, centrez $X,y$ avec les moyennes d’entraînement et ne mettez que les pentes dans $w$ ; retrouvez la constante non pénalisée par $\bar y-\bar x^\top w$. Mettez les variables à l’échelle dans les plis d’entraînement, car les pénalités dépendent des unités. La normalisation fait partie de la définition : Ridge dans scikit-learn utilise $\|Xw-y\|^2+\alpha\|w\|^2$, donc la formule Ridge ici correspond à $\alpha=2n\lambda$. Les conventions Lasso et Elastic Net correspondent à $\alpha=\lambda$, avec l1_ratio égal à $\rho$. À $\lambda$ identique, le terme Ridge d’Elastic Net avec $\rho=0$ a la moitié de l’intensité de celui de la formule Ridge ci-dessous.

Ridge

$$\hat{w}_{\text{ridge}} = \arg\min_w \frac{1}{2n}\|Xw-y\|_2^2 + \lambda\|w\|_2^2.$$

Ridge réduit les coefficients corrélés ou mal identifiés et les conserve généralement tous non nuls.

Lasso

$$\hat{w}_{\text{lasso}} = \arg\min_w \frac{1}{2n}\|Xw-y\|_2^2 + \lambda\|w\|_1.$$

La pénalité $\ell_1$ peut amener des coefficients exactement à zéro. Cette parcimonie est une propriété de l’objectif ajusté ; elle ne prouve pas que les variables retenues sont les seules importantes ni qu’elles sont causales. Avec des variables fortement corrélées, le membre sélectionné peut être instable.

Elastic Net

$$\hat{w}_{\text{EN}} = \arg\min_w \frac{1}{2n}\|Xw-y\|_2^2 + \lambda\left[ \rho\|w\|_1 + \frac{1-\rho}{2}\|w\|_2^2 \right].$$

Elastic Net associe la parcimonie à la stabilisation de Ridge et se révèle souvent utile lorsque les prédicteurs forment des groupes corrélés.

Règles pratiques

• Ajuster la mise à l’échelle et les transformations de base uniquement sur les données d’entraînement ; appliquer ensuite la transformation enregistrée aux données de validation et de test.
• Ne pas pénaliser l’ordonnée à l’origine, sauf si la formulation le prévoit explicitement.
• Choisir $\lambda$ et $\rho$ au moyen de la validation ou de la validation croisée au sein du processus d’entraînement.
• Comparer l’erreur de prédiction, la stabilité des coefficients et la simplicité opérationnelle, pas seulement la perte d’entraînement.
• Employer des objectifs robustes ou quantiles lorsque l’erreur quadratique ne représente pas la cible ou le coût des erreurs souhaité.
• Séparer la modélisation prédictive des affirmations inférentielles ; l’estimation de l’incertitude exige des hypothèses explicites sur le processus de génération des données.

Les variables polynomiales et les interactions peuvent rendre la relation avec les entrées non linéaire tout en conservant des coefficients ajustés linéaires. Le modèle obtenu hérite néanmoins des risques d’extrapolation et de surajustement liés à la base choisie.

Consultez le guide des modèles linéaires de scikit-learn pour les solveurs actuels et les API des estimateurs.

Mes lectures · Les notes ouvertes sont mémorisées dans ce navigateur.