Левое Матричное Деление и Numpy Solve - PullRequest
12 голосов
/ 23 августа 2011

Я пытаюсь преобразовать код, содержащий оператор \, из Matlab (Octave) в Python. Пример кода

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

Это работает и выдает 1.2 в качестве ответа. Использование этой веб-страницы

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

Я перевел это как:

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

Это дало мне ошибку:

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

Почему Matlab \ работает с неквадратной матрицей для B?

Какие-нибудь решения для этого?

Ответы [ 3 ]

15 голосов
/ 23 августа 2011

С Документация MathWorks для деления левой матрицы:

Если A - матрица размером m на n с m ~ = n, а B - вектор столбцов с m компоненты или матрица с несколькими такими столбцами, то X = A \ B является решение в смысле наименьших квадратов для недооцененных или переопределенных система уравнений AX = B. Другими словами, X минимизирует норму (A * X - B), длина вектора AX - B.

Эквивалент numpy равен 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 голосов
/ 23 августа 2011

Matlab фактически выполняет ряд различных операций, когда используется оператор \, в зависимости от формы используемых матриц (см. здесь для получения более подробной информации). В вашем примере, Matlab возвращает решение наименьших квадратов, а не решение линейного уравнения напрямую, как это происходит с квадратной матрицей. Чтобы получить такое же поведение в numpy, сделайте следующее:

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]

, который должен дать вам то же решение, что и Matlab.

2 голосов
/ 06 февраля 2017

Вы можете сформировать левое обратное:

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)

Результат:

c
 [[ 1.2]]

На самом деле, мы можем просто запустить решатель один раз, не формируя обратное, например:

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

Результат:

c
 [[ 1.2]]

.... как и раньше

Почему? Потому что:

Имеем:

enter image description here

Умножить на B.T, дает нам:

enter image description here

Теперь B.T.dot(B) - квадрат, полный ранг, имеет обратное. И поэтому мы можем умножить на обратное значение B.T.dot(B) или использовать решатель, как указано выше, чтобы получить c.

...