2011-08-23 22 views
11

Próbuję konwertować kod, który zawiera \ operator z Matlab (Octave) na język Python. Przykładowy kodLeft Matrix Division i Numpy Solve

B = [2;4] 
b = [4;4] 
B \ b 

To działa i daje 1,2 jako odpowiedź. Korzystanie z tej strony internetowej

http://mathesaurus.sourceforge.net/matlab-numpy.html

Tłumaczyłem, że jako:

import numpy as np 
import numpy.linalg as lin 
B = np.array([[2],[4]]) 
b = np.array([[4],[4]]) 
print lin.solve(B,b) 

To dało mi błąd:

numpy.linalg.linalg.LinAlgError: Array must be square 

Dlaczego Matlab \ pracuje non macierzy kwadratowej do B?

Jakieś rozwiązanie tego?

Odpowiedz

14

Z MathWorks documentation do podziału matrycy lewej:

If A is an m-by-n matrix with m ~= n and B is a column vector with m components, or a matrix with several such columns, then X = A\B is the solution in the least squares sense to the under- or overdetermined system of equations AX = B. In other words, X minimizes norm(A*X - B), the length of the vector AX - B.

Równowartość w numpy jest np.linalg.lstsq:

In [15]: B = np.array([[2],[4]]) 

In [16]: b = np.array([[4],[4]]) 

In [18]: x,resid,rank,s = np.linalg.lstsq(B,b) 

In [19]: x 
Out[19]: array([[ 1.2]]) 
8

Matlab będzie faktycznie wiele różnych czynności w czasie, stosuje się operator \, zależnie kształt użytych matryc (więcej szczegółów w artykule here). Na przykład, Matlab zwraca rozwiązanie najmniejszych kwadratów, zamiast bezpośrednio rozwiązywać równanie liniowe, jak w przypadku macierzy kwadratowej. Aby uzyskać takie samo zachowanie w trybie numpy, wykonaj następujące czynności:

import numpy as np 
import numpy.linalg as lin 
B = np.array([[2],[4]]) 
b = np.array([[4],[4]]) 
print np.linalg.lstsq(B,b)[0] 

które powinno dać takie samo rozwiązanie jak Matlab.

1

Można tworzyć lewy odwrotność:

import numpy as np 
import numpy.linalg as lin 
B = np.array([[2],[4]]) 
b = np.array([[4],[4]]) 

B_linv = lin.solve(B.T.dot(B), B.T) 
c = B_linv.dot(b) 
print('c\n', c) 

Wynik:

c 
[[ 1.2]] 

Faktycznie, możemy po prostu uruchomić solver raz, bez tworząc odwrócony, tak:

c = lin.solve(B.T.dot(B), B.T.dot(b)) 
print('c\n', c) 

Wynik:

c 
[[ 1.2]] 

.... tak jak poprzednio

Dlaczego? Ponieważ:

Mamy:

enter image description here

pomnożyć przez przez B.T, daje nam:

enter image description here

Teraz B.T.dot(B) jest kwadratowy, pełen ranking, ma odwrotną. I dlatego możemy pomnożyć przez odwrotność B.T.dot(B) lub użyć rozwiązania, jak wyżej, aby uzyskać c.