Ajustement de courbe

Un article de Wikipédia, l'encyclopédie libre.
Aller à : navigation, rechercher
Ajustement par itérations d'une courbe bruitée par un modèle de pic asymétrique (méthode de Gauss-Newton avec facteur d'amortissement variable).

L'ajustement de courbe est une technique d'analyse d'une courbe expérimentale, consistant à construire une courbe à partir de fonctions mathématiques et d'ajuster les paramètres de ces fonctions pour se rapprocher de la courbe mesurée — on parle donc aussi d'ajustement de paramètres. On utilise souvent le terme anglais curve fitting, profile fitting ou simplement fitting, pour désigner cette méthode ; on utilise souvent le franglais « fitter une courbe » pour dire « ajuster une courbe ».

On utilise des méthodes de régression. Dans les cas simples, il s'agit de régression multilinéaire si la loi est linéaire pour tous les paramètres, ou de régression polynomiale lorsque l'on utilise un polynôme pour simuler le phénomène (les paramètres physiques pouvant être déduits des coefficients du polynôme).

Les méthodes de régression classiques permettent de déterminer les paramètres à partir de calculs sur les données, mais sont inapplicables si la fonction est trop complexe. Il faut alors travailler par essai-erreur pour se rapprocher d'une solution, au sens de la méthode des moindres carrés. La solution n'est pas nécessairement unique.

Ajustement algébrique ou géométrique[modifier | modifier le code]

On distingue deux types d'ajustement :

  • les ajustement algébriques consistent à minimiser l'écart vertical (en y) entre la courbe modèle et les points expérimentaux ;
  • les ajustements géométriques consistent à minimiser la distance perpendiculairement à la courbe modèle ; c'est le cas par exemple de la régression circulaire ou elliptique, qui consiste à trouver le cercle, resp. l'ellipse, la plus proche des points

Dans le cas d'un ajustement géométrique, on parle de la méthode des moindres carrés totaux (total least square, TLS) : en effet, on prend en compte les deux coordonnées x et y pour la détermination de l'écart quadratique.

Fonction modèle utilisée[modifier | modifier le code]

Régression sur un nuage de points par un polynôme de degré croissant.

Dans certains cas, on a un modèle théorique permettant de prévoir la forme de la courbe ; la méthode d'ajustement permet de déterminer les paramètres de l'échantillon. Dans d'autres cas, on utilise une fonction empirique ; on s'intéresse alors en général à la surface, la largeur ou à la position du maximum de la fonction.

On utilise souvent des polynômes : ce sont des modèles qui se calculent facilement. Pour les formes de type « pic » (distributions unimodales), on utilise fréquemment des fonctions gaussiennes, lorentziennes ou bien des combinaisons (fonctions ou pseudo-fonctions de Voigt), ou encore des fonctions de Pearson. les courbes présentant un amortissement comportent fréquemment une composante exponentielle négative (fonction en e-x) ; les courbes en S peuvent être modélisées par une fonction sigmoïde..

De manière générale, on a une fonction ƒ modèle ayant n paramètres p1, p2, …, pn qui relie l'abscisse x à l'ordonnée y :

y = ƒ(p1, p2,…, pn, x)

et l'on compare cette fonction avec les m points expérimentaux

[(x1, y1exp), (x2, y2exp),…, (xm, ymexp)].

La fonction ƒ peut être parfois décomposée en plusieurs fonctions ƒ1, ƒ2… soit qu'elle en est la somme, le produit, le produit de convolution

ƒ = ƒ1 + ƒ2 +…
ou
ƒ = ƒ1׃2×…
ou
ƒ = ƒ1* ƒ2* …

Démarche de régression[modifier | modifier le code]

Considérations générales[modifier | modifier le code]

On calcule les points yical de la courbe simulée :

yical = ƒ(xi ).

Habituellement, on utilise un algorithme de régression visant à minimiser l'écart quadratique entre la courbe simulée et la courbe expérimentale ; on parle d'algorithme de minimisation de l'erreur.

On définit donc un facteur de fiabilité (reliability factor) R :

\mathrm{R} = \sqrt{\frac{\sum_i (y_i^{\mathrm{exp}}-y_i^{\mathrm{cal}})^2}{\sum_i {y_i^{\mathrm{exp}}}^2}}

yiexp est le i-ème point mesuré et yical est le i-ème point calculé. Le facteur R est similaire dans son expression au coefficient de corrélation multiple (R 2 étant alors le coefficient de détermination). On utilise plus couramment le facteur de fiabilité pondérée (weighted reliability factor) Rwp :

\mathrm{R_{wp}} = \sqrt{
\frac{\sum_i w_i \cdot (y_i^{\mathrm{exp}}-y_i^{\mathrm{cal}})^2}
{\sum_i w_i \cdot {y_i^{\mathrm{exp}}}^2}
}

wi est le poids attribué au point i ; ce poids représente l'incertitude associée au point i.

La régression consiste à trouver le jeu de paramètres (p1, p2, …, pn ) tel que Rwp est minimal. On a donc notamment

\forall i, \frac{\partial \mathrm{R_{wp}}}{\partial p^i} = 0.

Dans les cas simples, le système d'équations que l'on obtient peut se résoudre « simplement » ; en particulier pour les régressions multilinéaire, on peut le résoudre par inversion de matrice.

Régression non linéaire[modifier | modifier le code]

Article détaillé : Régression non linéaire.

Lorsqu'il n'est pas possible de résoudre simplement le système d'équations, on peut utiliser la méthode de Newton-Raphson.

C'est un algorithme itératif : à chaque étape, on regarde de quelle manière une petite variation d'un paramètre fait varier la fonction de fiabilité R. Le jeu de paramètre retenu pour l'étape suivante est obtenue par extrapolation de ces petites variations, par linéarisation : on cherche les valeurs des paramètres qui rendraient R nul si les variations de R étaient proportionnelles aux variations des paramètres. Le calcul s'arrête lorsque l'on n'arrive plus à diminuer le facteur de fiabilité, on parle de convergence.

Si l'on note (pk) l'ensemble des paramètres, l'intensité calculée en chaque point i à l'étape j s'exprime par

y_i^{\mathrm{cal}} = f(p_j^1, p_j^2, \ldots,p_j^n, x_i).

Pour simplifier les calculs, on peut faire un développement limité du premier ordre de cette fonction ƒ, alors en appliquant des incréments Δpk aux paramètres, on peut écrire

f(p_j^1 + \Delta p_j^1, \ldots,p_j^n + \Delta p^n, x_i) \simeq f(p_j^1, \ldots,p_j^n, x_i) + \sum_k \frac{\partial f}{\partial p^k}(x_i) \cdot \Delta p^k

en imposant yical = yiexp, on obtient ainsi un système d'équations linéaires (nous omettons l'indice j pour plus de clarté) :

\begin{pmatrix}
\frac{\partial f}{\partial p^1}(x_1) & \frac{\partial f}{\partial p^2}(x_1) & \cdots
   & \frac{\partial f}{\partial p^n}(x_1) \\
\frac{\partial f}{\partial p^1}(x_2) & \frac{\partial f}{\partial p^2}(x_2) & \cdots
   & \frac{\partial f}{\partial p^n}(x_2) \\
\vdots & \vdots & & \vdots \\
\frac{\partial f}{\partial p^1}(x_m) & \frac{\partial f}{\partial p^2}(x_m) & \cdots
   & \frac{\partial f}{\partial p^n}(x_m) \\
\end{pmatrix}
\begin{pmatrix}
\Delta p^1 \\ \Delta p^2 \\ \vdots \\ \Delta p^n
\end{pmatrix}
+
\begin{pmatrix}
f(p^1, \ldots, p^n, x_1) \\
 f(p^1, \ldots, p^n, x_2) \\
 \vdots \\
 f(p^1,\ldots, p^n, x_m) \\
\end{pmatrix}
=
\begin{pmatrix}y^{\mathrm{exp}}_1 \\
y^{\mathrm{exp}}_2 \\
\vdots \\
y^{\mathrm{exp}}_m \\
\end{pmatrix}

qui peut se résoudre par une pseudo-inversion de matrice, ce qui permet de calculer les Δpk et donc les valeurs des paramètres à l'étape suivante j + 1 :

si l'on appelle A la matrice \left ( \frac{\partial f}{\partial p^k}(x_i) \right )_{i, k}, le système d'équations s'écrit
A×(Δpkj)k + (yical)i = (yiexp)i
La matrice pseudo-inverse B a les mêmes dimensions que tA (la matrice transposée de A) et vérifie
  • A×B×A = A ;
  • B×A×B = B ;
  • A×B et B×A sont hermitiennes.

La matrice pseudo-inverse peut se calculer par décomposition en valeurs singulières.

On prend alors
pkj)k = B×[(yiexp)i - (yical)i]
et
(pkj + 1)k = (pkj)k + (Δpkj)k.

La pseudo-inversion est justifiée par le fait que le premier membre du système d'équations s'écrit alors

A×(Δpkj)k = A×B×[(yiexp)i - (yical)i] ;

si on le multiplie par B, il devient :

B×A×B×[(yiexp)i + (yical)i] = B×(yiexp)i + B×(yical)i

et donc on a bien

B×[A×(Δpkj)k - (yical)i] = B×(yiexp)i.

Notons qu'il s'agit d'une condition suffisante, mais pas nécessaire.

Il s'agit bien d'une méthode de Newton-Raphson, puisque l'on linéarise la fonction et que l'on cherche le point qui annule R (R = 0 si yical = yiexp pour tout i).

En remarquant qu'A est en fait la matrice jacobienne de ƒ (A = J), on voit quer le système d'équations est constitué des équations normales. Les incréments des paramètres peuvent alors se calculer par

pk)k = -(tJ×J)−1×tJ×[(yiexp)i - (yical)i].

On retrouve alors l'algorithme de Gauss-Newton.

Parfois, la linéarisation est « trop violente » : si la fonction présente une « courbure » importante lorsque les paramètres varient, l'application de la « loi tangente » peut éloigner de la solution plutôt que rapprocher. On peut donc utiliser un coefficient d'atténuation d ≤ 1, on a alors :

(pkj + 1)k = (pkj)k + d×(Δpkj)k.

On peut aussi utiliser un algorithme plus complexe utilisant les dérivées secondes.

Il peut arriver que le programme converge vers un minimum local de R qui n'est pas le minimum absolu. Pour éviter cette situation, on utilise l'algorithme de Metropolis-Hastings : lorsque l'on a convergé, on fait varier les paramètres d'une certaine valeur et on recommence le processus pour voir si l'on arrive à converger vers un autre jeu de paramètres ayant un facteur de fiabilité plus faible.

Si la simulation était parfaite, le facteur de fiabilité aurait une valeur dépendant du rapport signal sur bruit. Si l'on sait calculer ce rapport signal sur bruit, c'est-à-dire si l'on connaît la loi de probabilité régissant les fluctuations du signal, on peut alors déterminer un facteur de fiabilité incompressible R0. Si la simulation est parfaite, on a alors

R = R0

de fait, la qualité de la simulation est souvent exprimée par le rapport R/R0, qui doit tendre vers 1 au fur et à mesure des étapes itératives.

Régression circulaire[modifier | modifier le code]

La régression circulaire consiste à trouver le « meilleur cercle », au sens des moindres carrés, décrivant un ensemble de points. C'est un problème de régression géométrique non-linéaire. Cependant, un choix astucieux de la fonction d'erreur permet de se ramener à un problème linéaire.

Article détaillé : Régression circulaire.

Régression elliptique[modifier | modifier le code]

Régression elliptique

La régression elliptique consiste à trouver la « meilleure ellipse », au sens des moindres carrés, décrivant un ensemble de points.

Applications[modifier | modifier le code]

Voir aussi[modifier | modifier le code]

Articles connexes[modifier | modifier le code]

Liens externes[modifier | modifier le code]

Notes et références[modifier | modifier le code]