Moindres carrés : approximation avec un polynôme du second degré

Supposons que l’on veuille approximer un nuage de points avec un polynôme du second degré : y(x)=ax^2+bx+c.

SecondDegreePolynomial

Le nuage de points est donné par un ensemble de n points ayant pour coordonnées {x_i,y_i}. L’objectif est d’estimer \hat{a} et \hat{b} de façon à ce que le polynome d’équation y(x)=\hat{a}x^2+\hat{b}x + \hat{c} approche au mieux le nuage de points. Pour chaque point x_i nous allons minimiser la différence entre y_i et y(x_i), c-a-d. nous voulons minimiser \sum\limits_{i=1}^n{(y_i-y(x_i))^2}. La formulation matricielle du système est donnée par:

 \left[ \begin{matrix} {x_1}^2 & x_1 & 1 \\{x_2}^2 & x_2 & 1 \\... & ... & ... \\{x_n}^2 & x_n & 1 \\\end{matrix} \right].\left[ \begin{matrix} \hat{a} \\\hat{b} \\\hat{c} \end{matrix} \right]=\left[ \begin{matrix} y_1 \\y_2 \\... \\y_n \\\end{matrix} \right]

Définissons A, B et \hat{x}:

\begin{matrix}  A=\left[ \begin{matrix} {x_1}^2 & x_1 & 1 \\{x_2}^2 & x_2 & 1 \\... & ... & ... \\{x_n}^2 & x_n & 1 \\\end{matrix} \right] & B=\left[ \begin{matrix} y_1 \\y_2 \\...  \\y_n \\\end{matrix} \right] &\hat{x}=\left[ \begin{matrix} \hat{a} \\ \hat{b} \\ \hat{c} \end{matrix} \right] \end{matrix}

Le système est maintenant donné par:

 A.\hat{x}=B

La solution optimale est :

 \hat{x}=A^{+}.B = A^{T}(A.A^{T})^{-1}.B

A^{+} est la pseudo-inverse de A donnée par A^{+}=A^{T}(A.A^{T})^{-1}.

Téléchargement

6 réflexions au sujet de « Moindres carrés : approximation avec un polynôme du second degré »

  1. Bonjour,
    Je connais bien cette problématique.
    Je ne contredis pas la conclusion, mais je n’ai pas vu l’argumentation et surtout je doute qu’un utilisateur lambda puisse en tirer profit.
    La mise en forme matricielle est, certes, très attrayante, mais pourquoi utiliser le calcul matriciel, alors qu’il ne s’agit que de résoudre un système linéaire à trois inconnues ?
    Enfin pour ce type de problème, il me semble préférable de chercher un polynôme de degré quatre.
    Cordialement.

    1. toi ne pas connaitre cette problematique :
      * “il ne s’agit que de résoudre un système linéaire à trois inconnues ”
      -> comme tu peux le voir il y a un nuage de points donc NON

      * “il me semble préférable de chercher un polynôme de degré quatre”
      -> …… désespérant….

      conclusion : l’utilisateur lambda que tu es n’en a effectivement pas tirer parti

  2. Bonjour,

    J’ai essayer de mon côté avec l’exemple suivant :
    soit les points x y suivants : p1 = -3, 2 ; p2 = -2, 0.5 ; p3 = 1,-1 ; p4 = 2,-1 ; p5 = 3,0 ; p6 = 4,2 ; p7 = 5,4

    A x A_transposé serait du coup une matrice 7 x 7 :
    91 43 7 31 73 133 211
    43 21 3 13 31 57 91
    7 3 3 7 13 21 31
    31 13 7 21 43 73 111
    73 31 13 43 91 157 241
    133 57 21 73 157 273 421
    211 91 31 111 241 421 651

    si l’on veut avoir l’inverse de A x A_transposé, le déterminant est nul.
    Du coup je suis bloqué, il y a t’il une solution?
    Peut être un paramètre que j’ai oublié?

    1. Ça fonctionne pour moi. Le code suivant fait l’approximation :

      %% Matlab script for least square linear approximation
      %% This script shows the approximation with a second degree polynom
      %% Written by Philippe Lucidarme
      %% https://www.lucidarme.me
      close all;
      clear all;
      clc;

      % Set the coefficients a, b and c (y=ax²+bx+c)
      a=0.8;
      b=-5;
      c=-20;
      Noise=10;

      % Create a set of points with normal distribution
      %X=[-10:0.2:10];
      %Y=a*X.*X + b*X + c + Noise*randn(1,size(X,2));

      %p1 = -3, 2 ; p2 = -2, 0.5 ; p3 = 1,-1 ; p4 = 2,-1 ; p5 = 3,0 ; p6 = 4,2 ; p7 = 5,4

      X = [-3, -2, 1, 2, 3, 4, 5];
      Y = [ 2, 0.5, -1, -1, 0, 2, 4];

      % Prepare matrices
      A=[ (X.*X)’ , X’ , ones(size(X,2),1) ];
      B=Y’;

      % Least square approximation
      x=pinv(A)*B;

      % Display result
      plot (X,Y,’.’);
      hold on;
      plot (X,A*x,’g’);
      grid on;
      xlabel (‘x’);
      ylabel (‘y’);
      legend (‘Input data’,’Approximation’);
      title (‘Least square approximation’);

  3. Bonjour,

    merci de votre réponse,

    Je n’ai pas essayer avec matlab, mais de faire le calcul “à la main” dans un premier temps.
    J’imagine que matlab utilise quelques astuces pour résoudre certains blocages de calculs courants.

    J’ai dans premier temps essayer de mon côté de réaliser mon programme en php avec plusieurs librairie dédiées aux calculs de matrices. Ce fut un échec …

    J’ai donc opter pour une autre méthode utilisant notamment la résolution par méthode du pivot de Gauss … qui fonctionne!

    Ceci dit pour vérifier ma méthodologie avec les librairies php, si vous en avez le temps, pourriez vous décomposer la ligne : “x=pinv(A)*B;”. Je pense que l’ordre dans lequel je fait mes calculs de matrice n’est pas bon.
    J’aimerai tester plusieurs méthodes pour voir laquelle est la plus rapide (plusieurs milliers de calculs par jour) et de ne pas rester sur un échec!

  4. Matlab utilise la pseudoinverse de Moore-Penrose qui s’appuie la décomposition en valeur singulières.
    https://www.mathworks.com/help/matlab/ref/pinv.html

    En revanche, je pense que pour la méthode des moindres carrés, ça devrait fonctionner avec cette formule :
    https://ans.wiki/q/quelle-est-la-formule-de-la-matrice-pseudo-inverse/

    Elle est aussi utilisée sur cette page pour la régression linéaire :
    https://www.lucidarme.me/approximation-lineaire-par-les-moindres-carres/

Laisser un commentaire

Votre adresse de messagerie ne sera pas publiée. Les champs obligatoires sont indiqués avec *