Approximation linéaire

Supposons que l'on dispose de \(n\) données expérimentales \((x_i, y_i)\) et que l'on cherche à déduire une loi \(y=f(x)\) dépendant linéairement de \(m\) coefficients :

\[y = a_1 f_1(x) + a_2 f_2(x) + \dots + a_m f_m(x)\]

Idéalement, les \((x_i, y_i)\) vérifieraient cette relation, auquel cas on aurait le système de \(n\) équations linéaires à \(m\) inconnues \((a_i, ~\dots ,~a_m)\) :

\[ \left[ \begin{array}{c} y_1 \\ y_2 \\ \vdots \\ y_n \end{array} \right] = \left[ \begin{array}{cccc} f_1(x_1) & f_2(x_1) & \cdots & f_m(x_1) \\ f_1(x_2) & f_2(x_2) & \cdots & f_m(x_2) \\ \vdots & \vdots & & \vdots \\ f_1(x_n) & f_2(x_n) & \cdots & f_m(x_n) \end{array} \right] \left[ \begin{array}{c} a_1 \\ a_2 \\ \vdots \\ a_m \end{array} \right]\]

soit sous forme matricielle \(Y=F.~A\) auquel cas on trouverait simplement les coefficients \(a_i\) par inversion :

\[A = F^{-1} . Y\]

En MATLAB, ce genre d'instruction peut se programmer simplement à l'aide de l'opérateur \, et on écrirait alors

1
A = F\Y;

Cela dit le nombre d'équations n'est pas égal au nombre d'inconnues dans le système linéaire ci-dessus, il est (normalement) supérieur, et le système est surcontraint. Cela ne fait rien, l'opérateur \ le résout alors «au mieux» c'est-à-dire qu'il trouve les \(a_i\) qui minimise la somme des carrés des résidus :

\[ \min_{a_1\dots a_m} \sum_{i=1}^n \left[y_i - \sum_{j=1}^m a_j f_j \left( x_i \right)\right]^2\]

Exemple

On cherche à fitter (Pardons pour l'anglicisme, mais tout le monde l'utilise, et il n'y a pas d'équivalent en français ...) des points expérimentaux \((x_i, y_i)\) par une fonction de la forme :

\[f(x) = a_1 + \frac{a_2}{t} + \frac{a_3}{t^2}\]

La séquence d'instructions est la suivante :

1
2
3
4
5
6
7
8
9
10
11
12
xi = [1/4 1/2 1 2 3]; xi = xi';     % on définit les abscisses et ordonnees
yi = [18 7 4 2 1.1]; yi = yi';      %                       en colonnes
 
F = [ones(size(xi)) 1./xi 1./(xi.^2)];  % on peut écrire aussi
                                                % xi./xi a la place de
                                                % ones(size(xi))
 
 
A = F \ yi
 
x = 1/4:1/16:3;                  % On trace la courbe obtenue
plot(xi,yi,'o', x, A(1) + A(2)./x + A(3)./(x.^2))

Ce programme trace la courbe fittée ainsi que les points expérimentaux :

Approximation linéaire par la commande \
Approximation linéaire par la commande \

Cette méthode est utilisable dès que la fonction cherchée dépend linéairement des coefficients inconnus (le cas opposé est présenté dans la section suivante). En particulier, cela donne de bon résultat pour faire de l'approximation polynomiale (les \(f_j\) sont des monômes \(x^j\) ), mais dans ce cas on peut utiliser directement la fonction polyfit, qui en plus des coefficients du polynôme peut renvoyer des informations sur la qualité de l'approximation réalisée (voir la doc). polyfit attend en entrée simplement les \(x_i\) , les \(y_i\) et l'ordre du polynôme recherché. La fonction polyval peut être utilisée pour recalculer le polynôme d'approximation en tout point.

Exemple

L'exemple suivant calcule deux polynômes d'approximation, un d'ordre 3, l'autre d'ordre 4 et trace le résultat.

1
2
3
4
5
6
7
8
xi = [-2 -1.5 -1 0 1 1.5 2];
yi = [3.5  2 1.5 0 1.3 2.5 3.9];
 
A4 = polyfit (xi, yi, 4)
A3 = polyfit (xi, yi, 3)
 
x = -2:0.1:2;
plot(xi, yi, 'o', x, polyval(A4,x), x, polyval(A3,x), '--' )
Approximation linéaire par la commande polyfit
Approximation linéaire par la commande polyfit