Как найти общую обратную матрицу в Python с помощью SVD?

Итак, я пытаюсь найти обобщенную обратную A+ матрицы A, используя svd

это мой код для SVD

A4 = np.array([[3, 1,2,3,4],
               [5,2,3,4,5],
               [7,4,5,6,7]])
print(A4)
U4, s4, V4 = np.linalg.svd(A4)
D4 = sp.linalg.diagsvd(s4, 3, 5)

print('UDVt =\n', U4@D4@V4 )

Который отлично работает и дает мне правильную матрицу

И это код, который я пытаюсь использовать, чтобы получить обобщенное обратное

D4p = np.array([[D4[0,0]**-1,0,0],
                [0,D4[1,1]**-1,0],
                [0,0,D4[2,2]**-1],
                [0,0,0],
                [0,0,0]])
A4plus = np.transpose(V4) @ D4p @ np.transpose(U4) 
print(np.round(A4plus @ A4,2))
print(np.round(A4 @ A4plus,2))

Я попытался изменить транспонированную матрицу.

Так? Это странный способ сделать это (со всем этим ручным расширением). Но, кажется, это работает, не так ли? Вы получаете именно то, что дает np.linalg.pinv(A4) (я не знаю, делает ли он то же самое или нет. SVD — это один из способов вычисления обратного Мура-Пенроуза, но не единственный). Или именно то, что A4.T @ np.linalg.inv([email protected]) дает. Это другой путь. А A4@A4plus — это более или менее идентичность.

chrslg 26.04.2024 23:21

Разве не этого вы ожидали?

chrslg 26.04.2024 23:22

@chrslg Да, это то, чего я ожидал от A4@A4plus, но A4plus@A4 неправильный. Но я согласен с вами, что он дает то же самое, что и np.linalg.pinv(A4). SVD для меня новичок, поэтому мне было интересно, почему результаты такие разные, и я подумал, что допустил ошибку в коде. Но, кажется, с кодом все в порядке.

Ginto 27.04.2024 00:02

Не ждите, что A4plus @ A4 будет идентификационной матрицей. Обычно этого достичь невозможно, так как у вас недостаточно ранга строки. А A4 @ A4plus — это идентичность, как и должно быть.

Dr. V 27.04.2024 00:36
Почему в Python есть оператор "pass"?
Почему в Python есть оператор "pass"?
Оператор pass в Python - это простая концепция, которую могут быстро освоить даже новички без опыта программирования.
Некоторые методы, о которых вы не знали, что они существуют в Python
Некоторые методы, о которых вы не знали, что они существуют в Python
Python - самый известный и самый простой в изучении язык в наши дни. Имея широкий спектр применения в области машинного обучения, Data Science,...
Основы Python Часть I
Основы Python Часть I
Вы когда-нибудь задумывались, почему в программах на Python вы видите приведенный ниже код?
LeetCode - 1579. Удаление максимального числа ребер для сохранения полной проходимости графа
LeetCode - 1579. Удаление максимального числа ребер для сохранения полной проходимости графа
Алиса и Боб имеют неориентированный граф из n узлов и трех типов ребер:
Оптимизация кода с помощью тернарного оператора Python
Оптимизация кода с помощью тернарного оператора Python
И последнее, что мы хотели бы показать вам, прежде чем двигаться дальше, это
Советы по эффективной веб-разработке с помощью Python
Советы по эффективной веб-разработке с помощью Python
Как веб-разработчик, Python может стать мощным инструментом для создания эффективных и масштабируемых веб-приложений.
0
4
59
2
Перейти к ответу Данный вопрос помечен как решенный

Ответы 2

Ответ принят как подходящий

Предполагая, что строки не являются линейно зависимыми (иначе S-матрица имеет элементы нулевого значения и я не смогу их инвертировать), можно сказать

   U * S * Vt = A
=> A * (Vt.T * 1/S * U.T) = I

В коде:

from numpy import array, round
from numpy.linalg import svd
from scipy.linalg import diagsvd

A = array([[3, 1, 2, 3, 4],
           [5, 2, 3, 4, 5],
           [7, 4, 5, 6, 7]])

U, S, Vt = svd(A)

S_inv = diagsvd(1/S, A.shape[1], A.shape[0])
A_inv = Vt.T @ S_inv @ U.T

print(round(A @ A_inv, 5))

В результате это дает псевдообратное A

[[ 1.  0. -0.]  [ 0.  1.  0.]  [-0. -0.  1.]]

Вот соответствующий абзац по этому поводу: https://inst.eecs.berkeley.edu/~ee127/sp21/livebook/def_pseudo_inv.html

Хорошо, последний комментарий проясняет, в чем заключается настоящий вопрос: не о SVD как таковом, а о псевдообратном. Почему, если A⁺ является псевдоинверсией, A⁺A не является идентифицирующей (хотя AA⁺ является)?

Не имеет значения, как вычисляется A⁺.

Ap = np.linalg.pinv(A)
Ap = [email protected]([email protected])
U,s,V = np.linalg.svd(A); Ap = ([email protected](1/s, 3, 5)@V).T

так как результат тот же.

Итак, ответ на реальный вопрос прост: псевдоинверсия не является инверсией. Мы используем его, потому что обратного не существует.

Если бы A⁺A также было тождественным, то A⁺ было бы инверсией A. ∀u,v, Au=v ⇔ A⁺v=u, будет означать, что A обратимо, а его обратным является A⁺.

Чего не может быть. Так как это матрица 3×5.

Думайте об этом как о теории информации. Если бы A⁺A было тождеством, это означало бы, что если у вас есть 5 чисел для хранения в файле, вы могли бы поместить их в вектор U=[u₁,u₂,u₃,u₄,u₅]. Вычислите V=AV. Это дает вам только 3 числа V=[v₁,v₂,v₃]. Сохраните эти три числа в файле. Затем, когда вы захотите получить файл, просто вычислите U'=A⁺V, что в точности равно U.

Идеальный алгоритм сжатия. Это будет работать для всех данных U. ∀U, V=AU, тогда U'=A⁺V ⇒ U'=A⁺AU = IU = U. Не может быть правдой. Вы не можете представить весь возможный набор из 5 чисел, из которых 3 числа. Вы не можете вместить все 5D-пространство в 3D-пространство.

Алгоритм сжатия, конечно, существует. Но они работают только для подмножества возможных данных. Они работают, потому что мы знаем, что фактические данные, которые люди на самом деле используют, представляют собой лишь очень-очень небольшое подмножество всех возможных данных. Сжатие изображений работает только потому, что 99,99999999% изображений выглядят для нас как белый шум, и только оставшиеся 0,00000001% выглядят как изображение чего-либо. В противном случае это было бы невозможно.

Возвращаясь в линейный мир, вот что у вас здесь есть. «Сжатие» 5 чисел до 3, которое работает только для 3D-подмножества вашего 5D-пространства. То есть, если ваши данные представляют собой линейную комбинацию строк A.

Например, для U=[1,2,3,4,5] AU = [51,69,103]. И A⁺AU = [1,2,3,4,5]. Но это потому, что [1,2,3,4,5] — это 4*A₁-5A₂+2A₃ (3 строки A).

Итак, если вы случайно знаете, что U, с которым вы сталкиваетесь, лежат не во всем пятимерном пространстве, а только в небольшом его трехмерном подпространстве, то вы действительно можете ожидать, что эта псевдообратная работа будет работать (то есть, кстати, СВД - это нахождение действительно актуального подпространства).

Но если U может быть чем угодно, а не просто линейной комбинацией строк A, например [1,0,4,8,2]. И AU — это [43,59,89], и тогда A⁺AU — это [1,2,3,4,5], а не [1,0,4,8,2] (конечно, я выбрал свой Например, A@[1,2,3,4,5] = A@[1,0,4,8,2]. Но дело в том, что оно не может быть уникальным. Вы не можете иметь биекцию между 3D-пространствами. и 5D-пространство).

Кстати, [1,2,3,4,5] = A⁺@A@[1,0,4,8,2] не является каким-либо другим u, например Au=A[1,0,4,8,2 ]; а не просто любой другой вектор, который разделяет результат Au с [1,0,4,8,2]. Это тот, который представляет собой линейную комбинацию строк А. А это означает другое. Это означает, что A⁺A[1,0,4,8,2] — это линейная комбинация строк A, ближайшая (с точки зрения квадратической разности) к [1,0,4,8,2].

Тот факт, что A⁺A@[1,0,4,8,2] = [1,2,3,4,5] означает, что [1,2,3,4,5] является линейной комбинацией строк A, которые минимизируют ||AU-[1,0,4,8,2]||

Это еще один способ увидеть псевдоинверсию.

Если бы A была матрицей ранга 5. Тогда линейные комбинации строк A будут составлять все 5D-пространство. Другими словами, любой вектор X будет линейной комбинацией строк A. Итак, вы ожидаете, что A⁺AX будет X. Потому что X является линейной комбинацией строк A, и среди тех, к которым можно быть ближе X, чем сам X?
Итак, A⁺AX=X ∀X. Итак, A⁺A=Id. Таким образом, A⁺ — это не просто псевдоинверсия. Это обратное.

Но поскольку A не имеет ранга 5, то не все векторы X представляют собой линейную комбинацию строк A. И затем, если вы хотите найти линейную комбинацию строк A, напоминающих X, вам нужно решить задачу наименьших квадратов и найти Y, например, Y представляет собой линейную комбинацию строк A и ||X-Y|| является минимальным. Оно не может быть 0 каждый раз. Но вы можете подойти как можно ближе. Сказав Y=A⁺AX

Или, используя комбинации строк: A⁺ᵀ X — это коэффициенты строк A, наиболее близких к X. Например, для X=[1,0,4,8,2] Ap.T@X = [4 ,-5,2], что означает, что вы не можете получить αA₁+βA₂+γA₃ ближе к X, чем при α=4, β=-5 γ=2. И тогда αA₁+βA₂+γA₃ = A⁺AX =[1,2,3,4,5].

Но, конечно, если бы это было [1,0,4,8,2] вместо [1,2,3,4,5], это означало бы, что любой вектор X представляет собой линейную комбинацию строк A. И это не может быть, да? Никакой 5D-вектор не может быть линейной комбинацией трех векторов.

вр; доктор

A⁺A не является тождеством, потому что в противном случае оно было бы обратным, а не псевдоинверсным.

AA⁺ является тождественным, потому что AA⁺Y (Y представляет собой 3D-вектор) является ближайшей к Y линейной комбинацией столбцов A. Если A имеет ранг 3, существует одна линейная комбинация столбцов A, которая в точности соответствует Y для всех Y. Следовательно, AA⁺=Id

A⁺A не является тождественным, потому что A⁺AX (X представляет собой 5D-вектор) является ближайшей к X линейной комбинацией строк A. И вы не можете иметь точную комбинацию трех строк A, которая является X, для всех X.

Другие вопросы по теме