Moindres carrés sur une équation matricielle

[sylvain231]

je pense que j'arriverais à le programmer mais il me faudra une bibliothèque d'algèbre linéaire comme Numpy pour Python ou Lapack pour C++ car programmer à la main la résolution d'un système à neuf équations et neuf inconnues est complexe sinon (j'utilise les déterminants pour la résolution des systèmes)

[GaBuZoMeu]

j'utilise les déterminants pour la résolution des systèmes)

Pas top !

[tournesol]

Si tu as les coef , tu peux utiliser wims .

[sylvain231]

j'ai écrit le code suivant en Python, ça se lit bien tu peux me dire si c'est bon ?

Code: Tout sélectionner

import numpy as np
import random as r
N=100
A=np.zeros((3,N));#tableau de connues 3xN
B=np.zeros((3,N));#tableau de connues 3xN

#on remplit aléatoirement A et B
for i in range(3):
   for j in range(N):
      A[i,j]=r.random()
      B[i,j]=r.random()
      
M=np.zeros((3,3));#tableau d'inconnues 3x3
S=np.zeros((9,9));#coefficients du système 9x9
R=np.zeros(9);#coefficients résultats du système 9

for j in range(3):
   for l in range(3):
      for k in range(3):
         for i in range(N):
            S[l+3*j,k+3*j]=S[l+3*j,k+3*j]+B[k,i]*B[l,i];

for j in range(3):
   for l in range(3):
      for i in range(N):
         R[l+3*j]=R[l+3*j]+A[j,i]*B[l,i];
         
X = np.linalg.solve(S, R)

for j in range(3):
   for k in range(3):
      M[j,k]=X[k+3*j];
      
print(M)

[GaBuZoMeu]

J'aurais utilisé la fonction sum, mais ça me semble bien correspondre au système que je t'ai donné.

[sylvain231]

OK merci ma seule hésitation était pour les +3*j

[sylvain231]

Comment vérifier que c’est bon ? En fixant N à 1 et en faisant la multiplication matricielle ?

[sylvain231]

J’ai vérifié et obtient une erreur de 0 pour N=1 c’est bon signe !

[GaBuZoMeu]

Si tu prends tes et au hasard, tu ne verras pas grand chose. Mais si tu fixes , choisis des au hasard et prends des qui sont des perturbations pas trop grandes, au hasard, des , alors tu devrais retrouver .

[sylvain231]

Ok merci

[tournesol]

Pourquoi ne recherches tu pas une approximation affine en minimisant , avec M et P à déterminer ?

[Ben314]

Salut,
Sinon, concernant la partie programmation, tu peut pas mal raccourcir ton code vu que le système linéaire à résoudre correspond à chercher l'inverse d'une matrice 3x3 et que numpy sait le faire directement.
Si tu note la matrice dont les colonnes sont les et celle dont les colonnes sont les alors ta fonction à minimiser est où la norme sur les matrices est la "norme 2".
Et comme cette norme est issue du produit scalaire la différentielle de en est donc les points critique de sont les tels que c'est à dire et c'est ça ton système linéaire de 9 équations à 9 inconnues.
Ce système admet évidement une unique solution (qui est évidement ton minimum) lorsque la matrice est inversible et dans ce cas, la solution est .

Remarque : si et inversible, on retrouve bien évidement que la solution est .

[sylvain231]

quand je fait je trouve une matrice 1x1 !

[sylvain231]

voici mon code :

Code: Tout sélectionner

print("M2=",np.matmul(A.T,np.matmul(B,np.linalg.inv(np.matmul(B.T,B)))))

[GaBuZoMeu]

Ben314, j'ai l'impression que tu t'es fichu dedans en transposant du mauvais côté. Me trompé-je ?

[sylvain231]

la solution est donc ?

[GaBuZoMeu]

Ben314, j'ai l'impression que tu t'es fichu dedans en transposant du mauvais côté. Me trompé-je ?

Non, OK, ça revient bien exactement au système que j'avais écrit :

[tournesol]

Je trouve pour différentielle -2trace (t(MB-A)HB) , ce qui interdit d'isoler H .
J'y ai cru , ça aurait été trop simple. En general XY=0 ssi Im(Y) inclus dans ker(X) .
Quid de l'approximation affine ?

[sylvain231]

alors pourquoi je trouve une matrice 1x1 ?

[GaBuZoMeu]

Non non, tournesol, il faut y croire !

Soit une matrice de colonnes . Alors



Avec on a et . L'application linéaire tangente de en est donc