Méthode des moindres carrés généralisés
Jean Debord
30 avril 2003
1 Théorie
Soit à ajuster une fonction à une série de m points expérimentaux. On
suppose que, pour chaque observation k, la valeur observée de la variable
dépendante (yk ) est distribuée selon une loi normale de moyenne ŷk (ŷk étant
la valeur calculée de la variable dépendante) et de variance vk .
La valeur calculée ŷk est une fonction de la variable indépendante xk et
du vecteur des paramètres de la régression, B = [b0 , b1 , ...] :
ŷk = f (xk , B) (1)
vk est une fonction de ŷk et du vecteur des paramètres de la variance,
Θ = [θ0 , θ1 , ...], telle que :
vk = θ0 · g(ŷk , θ1 , ...) (2)
où g est une fonction spécifiée par l’utilisateur, par exemple une fonction
puissance :
g(ŷk , θ1 , ...) = (ŷk )θ1 (3)
Si l’on dispose de N courbes, correspondant chacune à un vecteur de
paramètres Bi (i = 1...N ), et si l’on considère que la fonction de variance
est la même pour toutes les courbes, le modèle général est alors le modèle de
régression non linéaire hétéroscédastique :
yik = f (xik , Bi ) + θ0 · g[f (xik , Bi ), θ1 , ...] · ik (4)
Les ik sont des variables aléatoires supposées indépendantes et identi-
quement distribuées selon une loi normale de moyenne nulle et de variance
unité. Dans ces conditions, la méthode qui s’impose de façon théorique est
1
la méthode du maximum de vraisemblance, qui consiste à minimiser simul-
tanément en tous les Bi et Θ la pseudo-vraisemblance (pseudo-likelihood,
P L) : " #
N N X mi
X X (yik − ŷik )2
PL = P Li = + ln vik (5)
i=1 i=1 k=1 vik
soit :
mi
" #
X [yik − f (xik , Bi )]2
P Li = + ln θ0 · g[f (xik , Bi ), θ1 , ...] (6)
k=1 θ0 · g[f (xik , Bi ), θ1 , ...]
et :
mi
N X
" #
X [yik − f (xik , Bi )]2
PL = + ln θ0 · g[f (xik , Bi ), θ1 , ...] (7)
i=1 k=1 θ0 · g[f (xik , Bi ), θ1 , ...]
où mi désigne le nombre d’observations pour la courbe i
Cependant, la mise en oeuvre pratique de ce calcul pour une pseudo-
vraisemblance un peu compliquée peut présenter des difficultés. Il vient alors
naturellement à l’esprit de choisir la méthode des moindres carrés pondérés.
Malheureusement il est possible parfois que l’estimateur ainsi construit ne soit
pas consistant, c’est-à-dire n’approche jamais la vraie valeur du paramètre.
Pour nous tirer de ce mauvais pas, nous allons faire une adaptation simple
de ce procédé en utilisant une méthode à deux étapes, c’est-à-dire que les
paramètres de la régression sont déterminés pour chaque courbe et les pa-
ramètres de la variance pour toutes les courbes en même temps [réf. 1,2].
Plus précisément, le processus itératif est tel qu’à l’itération courante :
1. Les paramètres de la variance sont estimés en minimisant la pseudo-
vraisemblance globale P L par rapport à Θ, avec les valeurs courantes
des paramètres de régression Bi
2. Les paramètres de la régression sont déterminés, pour chaque courbe,
en minimisant la pseudo-vraisemblance individuelle P Li par rapport à
Bi , en utilisant les paramètres de variance obtenus à l’étape précédente.
Ces deux étapes sont répétées jusqu’à convergence des paramètres de la
variance.
Les estimations initiales des paramètres de régression Bi sont obtenues
par régression non linéaire pondérée, les variances étant approchées par :
vik = g(yik , θ1 , ...) (8)
2
où θ1 , .. ont des valeurs fixées par l’utilisateur (p. ex. 1). Ceci revient
à calculer, pour chaque courbe i, le vecteur Bi qui minimise la somme des
carrés des écarts pondérée :
mi
X [yik − f (xik , Bi )]2
Si = (9)
k=1 g[f (xik , Bi ), θ1 , ...]
Le paramètre linéaire θ0 peut alors être estimé à partir des résultats de
cette régression pondérée, en minimisant la pseudo-vraisemblance totale P L
par rapport à θ0 seul, ce qui donne :
N
1 X
θ0 = Si (10)
Mt i=1
PN
où Mt désigne le nombre total d’observations, Mt = i=1 mi
2 Programmation
Les calculs précédents ont fait l’objet d’une programmation en Pascal [réf.
3] selon l’algorithme suivant :
1. Initialiser θj (j > 0)
2. Pour chaque courbe i : Trouver Bi qui minimise (9)
3. Estimer θ0 par (10)
4. Iter ← 0
5. Répéter :
(a) Iter ← Iter + 1
(b) Θ0 ← Θ
(c) Trouver Θ qui minimise (7)
(d) Pour chaque courbe i : Trouver Bi qui minimise (6)
jusqu’à (|θj − θj0 | < T OL · |θj0 | ∀j) ou (Iter > M AXIT ER)
où Iter désigne le numéro de l’itération, M AXIT ER le nombre maximal
d’itérations et T OL la précision requise sur les paramètres de la variance.
Le coeur du programme est constitué par la procédure de minimisation.
En principe, n’importe lequel des nombreux algorithmes existants peut conve-
nir. Nous avons pour notre part utilisé l’algorithme BFGS avec détermination
numérique du gradient.
3
3 Références
1. M. Davidian, D. M. Giltinan, Biometrics (1993) 49, 59-73
2. J. Debord, T. Dantoine, K. Suchaud, M. Harel, B. Verneuil, L. Merle.
Analusis (1997) 25, 293-297
3. J. Debord, K. Suchaud. Programme [Link] ([Link]
fr/pages_perso/[Link]/tpmath/[Link])