MIT 18.06SC, часть 4: SVD, смена базиса и псевдообратная матрица

Это перевод официальных конспектов (lecture summaries) к трём лекциям Гилберта Стрэнга из курса MIT 18.06SC. Карта всего курса — в путеводителе. Перевод распространяется под лицензией CC BY-NC-SA 4.0, © MIT OpenCourseWare.

Как читать вместе с курсом ChinaAI

Лекция 29 — к разделу «SVD: рентген матрицы» урока 1.1.3; перед ней прочитайте часть 3 о симметричных матрицах — SVD строится из собственных векторов $A^T A$. Лекция 31 — к разделу «Низкоранговое приближение»: там та же идея сжатия, только с SVD. Лекция 33 о псевдообратной матрице пригодится в уроке о линейной регрессии: решение методом наименьших квадратов — это $A^{+} b$. Блоки «Проверьте в NumPy» и примечания переводчика — наше дополнение.

Лекция 29. Сингулярное разложение

Сингулярное разложение матрицы (singular value decomposition) обычно называют SVD. Это последнее и лучшее разложение матрицы:

$$A = U \Sigma V^T,$$

где $U$ — ортогональная матрица, $\Sigma$ — диагональная, а $V$ — ортогональная.

В разложении $A = U \Sigma V^T$ матрица $A$ может быть любой. Мы знаем, что если $A$ симметрична и положительно определена, то её собственные векторы ортогональны и можно записать $A = Q \Lambda Q^T$. Это частный случай SVD, где $U = V = Q$. Для более общей $A$ в SVD нужны две разные матрицы $U$ и $V$.

Мы также научились записывать $A = S \Lambda S^{-1}$, где $S$ — матрица из $n$ различных собственных векторов $A$. Однако $S$ может не быть ортогональной; матрицы $U$ и $V$ в SVD ортогональны всегда.

Как это работает

Можно думать об $A$ как о линейном преобразовании, которое переводит вектор $v_1$ из её пространства строк в вектор $u_1 = A v_1$ в её пространстве столбцов. SVD возникает, когда мы ищем ортогональный базис пространства строк, который переходит в ортогональный базис пространства столбцов: $A v_i = \sigma_i u_i$.

Найти ортогональный базис пространства строк нетрудно — процесс Грама — Шмидта сразу его даёт. Но в общем случае нет причин ожидать, что $A$ переведёт этот базис в другой ортогональный базис.

Возможно, вы спросите о векторах в нуль-пространствах $A$ и $A^T$. С ними проблем нет — о них позаботятся нули на диагонали $\Sigma$.

На языке матриц

Суть задачи — найти ортонормированный базис $v_1, v_2, \dots, v_r$ пространства строк $A$, для которого

$$A \begin{bmatrix} v_1 & v_2 & \cdots & v_r \end{bmatrix} = \begin{bmatrix} \sigma_1 u_1 & \sigma_2 u_2 & \cdots & \sigma_r u_r \end{bmatrix} = \begin{bmatrix} u_1 & u_2 & \cdots & u_r \end{bmatrix} \begin{bmatrix} \sigma_1 & & & \\ & \sigma_2 & & \\ & & \ddots & \\ & & & \sigma_r \end{bmatrix},$$

где $u_1, u_2, \dots, u_r$ — ортонормированный базис пространства столбцов $A$. Когда мы добавим нуль-пространства, это равенство превратится в $AV = U\Sigma$. (Ортонормированные базисы $v_1, \dots, v_r$ и $u_1, \dots, u_r$ можно дополнить до ортонормированных базисов всего пространства как угодно. Поскольку $v_{r+1}, \dots, v_n$ лежат в нуль-пространстве $A$, диагональные элементы $\sigma_{r+1}, \dots, \sigma_n$ будут равны 0.)

Столбцы $V$ и $U$ — базисы пространства строк и пространства столбцов соответственно. Обычно $U \ne V$, но если $A$ положительно определена, для её пространства строк и пространства столбцов можно взять один и тот же базис!

Примечание переводчика: в оригинале порядок перепутан («столбцы $U$ и $V$ — базисы пространства строк и столбцов соответственно»). Верно так, как в переводе: $v_i$ живут в пространстве строк, $u_i$ — в пространстве столбцов.

Вычисление

Пусть $A$ — обратимая матрица $\begin{bmatrix} 4 & 4 \\ -3 & 3 \end{bmatrix}$. Мы хотим найти векторы $v_1$ и $v_2$ в пространстве строк $\mathbb{R}^2$, $u_1$ и $u_2$ в пространстве столбцов $\mathbb{R}^2$ и положительные числа $\sigma_1$ и $\sigma_2$ так, чтобы $v_i$ были ортонормированы, $u_i$ были ортонормированы, а $\sigma_i$ были коэффициентами растяжения, для которых $A v_i = \sigma_i u_i$.

Это большой шаг к тому, чтобы найти ортонормированные матрицы $V$ и $U$ и диагональную матрицу $\Sigma$, для которых

$$AV = U\Sigma.$$

Поскольку $V$ ортогональна, обе части можно умножить на $V^{-1} = V^T$ и получить

$$A = U \Sigma V^T.$$

Вместо того чтобы искать $U$, $V$ и $\Sigma$ одновременно, умножим обе части на $A^T = V \Sigma^T U^T$:

$$A^T A = V \Sigma^T U^T U \Sigma V^T = V \Sigma^2 V^T = V \begin{bmatrix} \sigma_1^2 & & & \\ & \sigma_2^2 & & \\ & & \ddots & \\ & & & \sigma_n^2 \end{bmatrix} V^T.$$

Это форма $Q \Lambda Q^T$; теперь $V$ можно найти, диагонализовав симметричную положительно определённую (или полуопределённую) матрицу $A^T A$. Столбцы $V$ — собственные векторы $A^T A$, а собственные числа $A^T A$ — это $\sigma_i^2$. (Берём $\sigma_i$ равным положительному квадратному корню из $\lambda_i$.)

Чтобы найти $U$, делаем то же самое с $A A^T$.

Примечание переводчика: в оригинале в средней части формулы опечатка — «$V \Sigma U^{-1} U \Sigma V^T$». Сокращение работает потому, что $U^T U = I$ для ортогональной $U$.

Пример SVD

Вернёмся к нашей матрице $A = \begin{bmatrix} 4 & 4 \\ -3 & 3 \end{bmatrix}$. Начнём с вычисления

$$A^T A = \begin{bmatrix} 4 & -3 \\ 4 & 3 \end{bmatrix} \begin{bmatrix} 4 & 4 \\ -3 & 3 \end{bmatrix} = \begin{bmatrix} 25 & 7 \\ 7 & 25 \end{bmatrix}.$$

Собственные векторы этой матрицы дадут нам векторы $v_i$, а собственные числа — числа $\sigma_i$.

Два ортогональных собственных вектора $A^T A$ — это $\begin{bmatrix} 1 \\ 1 \end{bmatrix}$ и $\begin{bmatrix} 1 \\ -1 \end{bmatrix}$. Чтобы получить ортонормированный базис, положим $v_1 = \begin{bmatrix} 1/\sqrt{2} \\ 1/\sqrt{2} \end{bmatrix}$ и $v_2 = \begin{bmatrix} 1/\sqrt{2} \\ -1/\sqrt{2} \end{bmatrix}$. Их собственные числа $\sigma_1^2 = 32$ и $\sigma_2^2 = 18$. Теперь у нас есть:

$$\underbrace{\begin{bmatrix} 4 & 4 \\ -3 & 3 \end{bmatrix}}_{A} = \underbrace{\begin{bmatrix} \ & \ \\ \ & \ \end{bmatrix}}_{U} \underbrace{\begin{bmatrix} 4\sqrt{2} & 0 \\ 0 & 3\sqrt{2} \end{bmatrix}}_{\Sigma} \underbrace{\begin{bmatrix} 1/\sqrt{2} & 1/\sqrt{2} \\ 1/\sqrt{2} & -1/\sqrt{2} \end{bmatrix}}_{V^T}.$$

Можно было бы найти $U$ отсюда, но для практики найдём $U$ через ортонормированные собственные векторы $u_1$ и $u_2$ матрицы $A A^T = U \Sigma^2 U^T$:

$$A A^T = \begin{bmatrix} 4 & 4 \\ -3 & 3 \end{bmatrix} \begin{bmatrix} 4 & -3 \\ 4 & 3 \end{bmatrix} = \begin{bmatrix} 32 & 0 \\ 0 & 18 \end{bmatrix}.$$

К счастью, $A A^T$ оказалась диагональной. Хочется взять $u_1 = \begin{bmatrix} 1 \\ 0 \end{bmatrix}$ и $u_2 = \begin{bmatrix} 0 \\ 1 \end{bmatrix}$, как профессор Стрэнг сделал на лекции, но поскольку $A v_2 = \begin{bmatrix} 0 \\ -3\sqrt{2} \end{bmatrix}$, на самом деле $u_2 = \begin{bmatrix} 0 \\ -1 \end{bmatrix}$ и $U = \begin{bmatrix} 1 & 0 \\ 0 & -1 \end{bmatrix}$. Заодно это даёт возможность перепроверить вычисление $\sigma_1$ и $\sigma_2$.

Итак, SVD матрицы $A$:

$$\underbrace{\begin{bmatrix} 4 & 4 \\ -3 & 3 \end{bmatrix}}_{A} = \underbrace{\begin{bmatrix} 1 & 0 \\ 0 & -1 \end{bmatrix}}_{U} \underbrace{\begin{bmatrix} 4\sqrt{2} & 0 \\ 0 & 3\sqrt{2} \end{bmatrix}}_{\Sigma} \underbrace{\begin{bmatrix} 1/\sqrt{2} & 1/\sqrt{2} \\ 1/\sqrt{2} & -1/\sqrt{2} \end{bmatrix}}_{V^T}.$$

Почему знак важен: собственный вектор определён с точностью до знака, поэтому $u_2$ «по $A A^T$» мог выйти и $(0, 1)$, и $(0, -1)$. Но $u_i$ и $v_i$ связаны условием $A v_i = \sigma_i u_i$ с положительным $\sigma_i$, и знак $u_i$ надо согласовать с выбранным $v_i$. Надёжнее всего считать $u_i = A v_i / \sigma_i$.

Пример с нуль-пространством

Теперь пусть $A = \begin{bmatrix} 4 & 3 \\ 8 & 6 \end{bmatrix}$. У неё одномерное нуль-пространство и одномерные пространства строк и столбцов.

Пространство строк $A$ состоит из кратных $\begin{bmatrix} 4 \\ 3 \end{bmatrix}$. Пространство столбцов $A$ состоит из кратных $\begin{bmatrix} 4 \\ 8 \end{bmatrix}$. Нуль-пространство и левое нуль-пространство перпендикулярны пространству строк и пространству столбцов соответственно.

Единичные базисные векторы пространства строк и пространства столбцов: $v_1 = \begin{bmatrix} 0{,}8 \\ 0{,}6 \end{bmatrix}$ и $u_1 = \begin{bmatrix} 1/\sqrt{5} \\ 2/\sqrt{5} \end{bmatrix}$. Чтобы вычислить $\sigma_1$, найдём ненулевое собственное число $A^T A$:

$$A^T A = \begin{bmatrix} 4 & 8 \\ 3 & 6 \end{bmatrix} \begin{bmatrix} 4 & 3 \\ 8 & 6 \end{bmatrix} = \begin{bmatrix} 80 & 60 \\ 60 & 45 \end{bmatrix}.$$

Это матрица ранга 1, так что одно её собственное число должно быть 0. Другое тогда равно следу, так что $\sigma_1^2 = 125$. Найдя единичные векторы, перпендикулярные $u_1$ и $v_1$ (базисные векторы левого нуль-пространства и нуль-пространства соответственно), видим, что SVD матрицы $A$:

$$\underbrace{\begin{bmatrix} 4 & 3 \\ 8 & 6 \end{bmatrix}}_{A} = \underbrace{\frac{1}{\sqrt{5}} \begin{bmatrix} 1 & 2 \\ 2 & -1 \end{bmatrix}}_{U} \underbrace{\begin{bmatrix} \sqrt{125} & 0 \\ 0 & 0 \end{bmatrix}}_{\Sigma} \underbrace{\begin{bmatrix} 0{,}8 & 0{,}6 \\ 0{,}6 & -0{,}8 \end{bmatrix}}_{V^T}.$$

Проверка числами: $\sigma_1 u_1 v_1^T = \sqrt{125} \cdot \frac{1}{\sqrt{5}} \begin{bmatrix} 1 \\ 2 \end{bmatrix} \begin{bmatrix} 0{,}8 & 0{,}6 \end{bmatrix} = 5 \begin{bmatrix} 0{,}8 & 0{,}6 \\ 1{,}6 & 1{,}2 \end{bmatrix} = \begin{bmatrix} 4 & 3 \\ 8 & 6 \end{bmatrix}$ — вся матрица собрана из одного «слоя» ранга 1.

Сингулярное разложение объединяет темы линейной алгебры от положительно определённых матриц до четырёх фундаментальных подпространств:

• $v_1, v_2, \dots, v_r$ — ортонормированный базис пространства строк;
• $u_1, u_2, \dots, u_r$ — ортонормированный базис пространства столбцов;
• $v_{r+1}, \dots, v_n$ — ортонормированный базис нуль-пространства;
• $u_{r+1}, \dots, u_m$ — ортонормированный базис левого нуль-пространства.

Это «правильные» базисы, потому что $A v_i = \sigma_i u_i$.

Проверьте в NumPy (дополнение переводчика)

import numpy as np

A = np.array([[ 4., 4.],
              [-3., 3.]])
U, s, Vt = np.linalg.svd(A)     # NumPy возвращает Vᵀ, а σ — одномерным массивом
print(s)                        # [5.6569 4.2426] = 4√2 и 3√2
print(U)                        # [[-1.  0.]
                                #  [ 0.  1.]]
print(Vt)                       # [[-0.7071 -0.7071]
                                #  [-0.7071  0.7071]]
# Знаки у NumPy другие, чем у Стрэнга: пары (u_i, v_i) можно одновременно умножить на -1.
print(np.allclose(U @ np.diag(s) @ Vt, A))   # True

B = np.array([[4., 3.],
              [8., 6.]])
U, s, Vt = np.linalg.svd(B)
print(s)                        # [11.1803  0.    ] — σ₁ = √125, σ₂ = 0 (ранг 1)
print(s[0] * np.outer(U[:, 0], Vt[0]))       # [[4. 3.] [8. 6.]] — один слой ранга 1

Лекция 31. Смена базиса; сжатие изображений

Мы уже знаем, что удачный выбор базиса упрощает вычисления. Одно из применений этого принципа — сжатие изображений. Видео лекций, музыка и другие источники данных несут много информации; хранить и передавать её эффективно можно, только сменив базис, в котором она записана.

Сжатие изображений

Пусть один кадр видео лекции — это 512 на 512 пикселей и видео чёрно-белое. Камера записывает уровень яркости каждого из $512^2$ пикселей; в этом смысле каждый кадр — вектор в векторном пространстве размерности $512^2$.

В стандартном базисе этого пространства на каждый пиксель приходится один вектор. Передавать все $512^2$ компонент каждого кадра в стандартном базисе — слишком широкий канал, но если сменить базис по стандарту сжатия изображений JPEG, вполне приличную копию видео можно передать очень экономно.

Например, если сообщать уровень освещённости пиксель за пикселем, нет экономного способа передать информацию «весь кадр чёрный». Но если один из базисных векторов соответствует тому, что у всех пикселей одинаковая яркость (скажем, 1), то запись пустой доски передаётся очень экономно.

Вместе с вектором из одних единиц можно взять базисный вектор, в котором чередуются 1 и $-1$, или вектор, у которого половина единиц и половина $-1$, — он соответствует картинке, светлой слева и тёмной справа. От выбора базиса напрямую зависит, сколько данных придётся скачать, чтобы посмотреть видео, и лучший базис для лекций по алгебре может отличаться от лучшего базиса для боевиков!

Базисные векторы Фурье

$$\begin{bmatrix} 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \end{bmatrix},\ \begin{bmatrix} 1 \\ \omega \\ \omega^2 \\ \omega^3 \\ \omega^4 \\ \omega^5 \\ \omega^6 \\ \omega^7 \end{bmatrix},\ \begin{bmatrix} 1 \\ \omega^2 \\ \omega^4 \\ \omega^6 \\ \omega^8 \\ \omega^{10} \\ \omega^{12} \\ \omega^{14} \end{bmatrix},\ \cdots,\ \begin{bmatrix} 1 \\ \omega^7 \\ \omega^{14} \\ \omega^{21} \\ \omega^{28} \\ \omega^{35} \\ \omega^{42} \\ \omega^{49} \end{bmatrix}$$

Самый известный базис — базис Фурье, тесно связанный с матрицами Фурье, которые мы изучали раньше. Базис JPEG состоит из косинусов — вещественных частей $\omega^{jk}$.

Примечание переводчика: здесь $\omega = e^{2\pi i/8}$ — комплексный корень восьмой степени из единицы (матрицы Фурье разбирались на занятии 3.2, лекция 26).

Метод разбивает квадрат $512 \times 512$ пикселей на блоки со стороной 8 пикселей, по 64 пикселя в каждом. Затем информация о яркости этих пикселей сжимается — например, отбрасываются все коэффициенты ниже некоторого порога, выбранного так, чтобы после их удаления разница была почти незаметна.

$$\text{сигнал } x \ \xrightarrow{\text{без потерь}}\ 64 \text{ коэффициента } c \ \xrightarrow{\text{сжатие с потерями}}\ \hat{c} \text{ (много нулей)} \ \longrightarrow\ \hat{x} = \sum \hat{c}_i v_i$$

В видео можно сжимать не только каждый кадр, но и последовательности кадров. Соседние кадры почти не отличаются. Если сделать всё правильно, кодировать и сжимать нужно лишь разницу между кадрами, а не каждый кадр целиком.

Вейвлет-базис Хаара

$$\begin{bmatrix} 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \\ 1 \end{bmatrix},\ \begin{bmatrix} 1 \\ 1 \\ 1 \\ 1 \\ -1 \\ -1 \\ -1 \\ -1 \end{bmatrix},\ \begin{bmatrix} 1 \\ 1 \\ -1 \\ -1 \\ 0 \\ 0 \\ 0 \\ 0 \end{bmatrix},\ \begin{bmatrix} 0 \\ 0 \\ 0 \\ 0 \\ 1 \\ 1 \\ -1 \\ -1 \end{bmatrix},\ \begin{bmatrix} 1 \\ -1 \\ 0 \\ 0 \\ 0 \\ 0 \\ 0 \\ 0 \end{bmatrix},\ \cdots,\ \begin{bmatrix} 0 \\ 0 \\ 0 \\ 0 \\ 0 \\ 0 \\ 1 \\ -1 \end{bmatrix}$$

Ближайший конкурент метода JPEG использует вейвлет-базис. (JPEG2000 улучшает вейвлеты Хаара, показанные выше.) В вейвлет-базисе Хаара для $\mathbb{R}^8$ ненулевые элементы — наполовину 1 и наполовину $-1$ (кроме вектора из одних единиц). При этом половина или даже три четверти элементов базисного вектора могут быть нулями. Эти векторы выбраны ортогональными, и их можно нормировать до ортонормированных.

Сжатие и матрицы

Линейная алгебра нужна, чтобы найти коэффициенты $c_i$ при переходе от стандартного базиса (яркость каждого пикселя) к базису Фурье или вейвлет-базису. Например, мы хотим записать:

$$x = c_1 w_1 + \dots + c_8 w_8.$$

Но это просто линейная комбинация векторов вейвлет-базиса. Если $W$ — матрица, столбцы которой — векторы вейвлет-базиса, то наша задача — всего лишь найти $c$:

$$x = W \begin{bmatrix} c_1 \\ c_2 \\ \vdots \\ c_8 \end{bmatrix}.$$

Итак, $c = W^{-1} x$.

Вычисления будут быстрее и проще, если не тратить много времени на обращение матрицы (например, если $W^{-1} = W^T$) или на умножение на обратную. Поэтому в сжатии изображений критерии хорошего базиса такие:

• умножение на матрицу базиса и на обратную к ней выполняется быстро (как в БПФ или в вейвлет-базисе);
• хорошее сжатие — изображение можно приблизить всего несколькими базисными векторами. Большинство компонент $c_i$ малы, и их можно смело обнулить.

Смена базиса

Векторы. Пусть столбцы матрицы $W$ — векторы нового базиса. Тогда, если $x$ — вектор в старом базисе, его можно перевести в вектор $c$ в новом базисе по соотношению

$$x = Wc.$$

Матрицы преобразований. Пусть у нас есть линейное преобразование $T$. Если в базисе $v_1, v_2, \dots, v_8$ у $T$ матрица $A$, а в базисе $w_1, w_2, \dots, w_8$ — матрица $B$, то оказывается, что $A$ и $B$ обязаны быть подобными матрицами. Иначе говоря, $B = M^{-1} A M$ для некоторой матрицы смены базиса $M$.

Напоминание. Если у нас есть базис $v_1, v_2, \dots, v_8$ и мы знаем $T(v_i)$ для каждого $i$, то, пользуясь линейностью $T$, можно найти $T(v) = c_1 T(v_1) + c_2 T(v_2) + \dots + c_8 T(v_8)$ для любого вектора $v = c_1 v_1 + c_2 v_2 + \dots + c_8 v_8$ пространства. Элементы столбца $i$ матрицы $A$ — это коэффициенты выходного вектора $T(v_i)$.

Если наш базис состоит из собственных векторов преобразования, то есть $T(v_i) = \lambda_i v_i$, то $A = \Lambda$ — (диагональная) матрица собственных чисел. Было бы замечательно использовать для обработки изображений базис из собственных векторов, но найти такой базис намного дороже, чем просто взять базис Фурье или вейвлет-базис.

Итог

При смене базиса координаты векторов меняются по правилу $x = Wc$. Элементы матриц меняются по правилу $B = M^{-1} A M$.

Проверьте в NumPy (дополнение переводчика)

import numpy as np

# Столбцы W — восемь векторов Хаара (первые пять показаны в формуле выше)
W = np.array([[1,  1,  1,  0,  1,  0,  0,  0],
              [1,  1,  1,  0, -1,  0,  0,  0],
              [1,  1, -1,  0,  0,  1,  0,  0],
              [1,  1, -1,  0,  0, -1,  0,  0],
              [1, -1,  0,  1,  0,  0,  1,  0],
              [1, -1,  0,  1,  0,  0, -1,  0],
              [1, -1,  0, -1,  0,  0,  0,  1],
              [1, -1,  0, -1,  0,  0,  0, -1]], dtype=float)
print(np.diag(W.T @ W))   # [8. 8. 4. 4. 2. 2. 2. 2.] — вне диагонали нули: столбцы ортогональны

x = np.array([5, 5, 5, 5, 1, 1, 1, 1.])   # «слева светло, справа темно»
c = np.linalg.solve(W, x)                  # c = W⁻¹ x
print(np.round(c, 6))     # [3. 2. 0. 0. 0. 0. 0. 0.] — 8 пикселей описаны двумя числами

Лекция 33. Левая и правая обратные; псевдообратная матрица

Псевдообратных матриц не будет на экзамене, но эта лекция поможет к нему подготовиться.

Двусторонняя обратная

Двусторонняя обратная матрицы $A$ — это матрица $A^{-1}$, для которой $A A^{-1} = I = A^{-1} A$. Это то, что мы называли обратной к $A$. Здесь $r = n = m$; у матрицы $A$ полный ранг.

Левая обратная

Напомним: у $A$ полный ранг по столбцам, если её столбцы независимы, то есть $r = n$. В этом случае нуль-пространство $A$ состоит только из нулевого вектора. Уравнение $Ax = b$ либо имеет ровно одно решение $x$, либо не имеет решений.

Матрица $A^T A$ — обратимая симметричная матрица $n \times n$, так что $(A^T A)^{-1} A^T A = I$. Мы говорим, что $A^{-1}_{\text{left}} = (A^T A)^{-1} A^T$ — левая обратная к $A$. (Могут быть и другие левые обратные, но эта — наша любимая.) Тот факт, что $A^T A$ обратима, когда у $A$ полный ранг по столбцам, был центральным в нашем разговоре о методе наименьших квадратов.

Заметьте, что $A A^{-1}_{\text{left}}$ — матрица $m \times m$, которая равна единичной, только если $m = n$. У прямоугольной матрицы не может быть двусторонней обратной, потому что у неё самой или у её транспонированной есть ненулевое нуль-пространство.

Правая обратная

Если у $A$ полный ранг по строкам, то $r = m$. Нуль-пространство $A^T$ содержит только нулевой вектор; строки $A$ независимы. Уравнение $Ax = b$ всегда имеет хотя бы одно решение; размерность нуль-пространства $A$ равна $n - m$, так что будет $n - m$ свободных переменных и (если $n > m$) бесконечно много решений!

У матриц полного ранга по строкам есть правые обратные $A^{-1}_{\text{right}}$, для которых $A A^{-1}_{\text{right}} = I$. Самая удобная из них — $A^T (A A^T)^{-1}$. Проверка: $A$, умноженная на $A^T (A A^T)^{-1}$, — это $I$.

Псевдообратная

У обратимой матрицы ($r = m = n$) нуль-пространство и левое нуль-пространство содержат только нулевой вектор. У матрицы полного ранга по столбцам $r = n$ только нулевой вектор в нуль-пространстве. У матрицы полного ранга по строкам $r = m$ только нулевой вектор в левом нуль-пространстве. Остаётся рассмотреть матрицу $A$, для которой $r < n$ и $r < m$.

Если у $A$ полный ранг по столбцам и $A^{-1}_{\text{left}} = (A^T A)^{-1} A^T$, то

$$A A^{-1}_{\text{left}} = A (A^T A)^{-1} A^T = P$$

— матрица проекции $\mathbb{R}^m$ на пространство столбцов $A$. Ближе к равенству $AM = I$ подойти нельзя.

Аналогично, если у $A$ полный ранг по строкам, то $A^{-1}_{\text{right}} A = A^T (A A^T)^{-1} A$ — матрица проекции $\mathbb{R}^n$ на пространство строк $A$.

Обращать матрицы мешают нетривиальные нуль-пространства. Если $Ax = 0$ для какого-то ненулевого $x$, то нет никакой надежды найти матрицу $A^{-1}$, которая обратит этот процесс и даст $A^{-1} 0 = x$.

Вектор $Ax$ всегда лежит в пространстве столбцов $A$. Более того, соответствие между векторами $x$ в ($r$-мерном) пространстве строк и векторами $Ax$ в ($r$-мерном) пространстве столбцов взаимно однозначно. Иначе говоря, если $x \ne y$ — векторы пространства строк $A$, то $Ax \ne Ay$ в пространстве столбцов $A$. (Доказательство этого стало бы хорошей экзаменационной задачей.)

Доказательство: если x ≠ y, то Ax ≠ Ay

Предположим, что утверждение неверно. Тогда найдутся $x \ne y$ в пространстве строк $A$, для которых $Ax = Ay$. Но тогда $A(x - y) = 0$, так что $x - y$ лежит в нуль-пространстве $A$. Но пространство строк $A$ замкнуто относительно линейных комбинаций (например, вычитания), так что $x - y$ лежит и в пространстве строк. Единственный вектор, который лежит и в нуль-пространстве, и в пространстве строк, — нулевой, так что $x - y = 0$. Это противоречит предположению, что $x$ и $y$ не равны.

Значит, отображение $x \mapsto Ax$ из пространства строк в пространство столбцов обратимо. Обратное к нему отображение называется псевдообратным (pseudoinverse), и оно очень полезно статистикам в работе с линейной регрессией — они не всегда могут гарантировать, что у их матриц полный ранг по столбцам $r = n$.

Как найти псевдообратную A⁺

Псевдообратная $A^{+}$ к $A$ — это матрица, для которой $x = A^{+} A x$ для всех $x$ в пространстве строк $A$. Нуль-пространство $A^{+}$ — это нуль-пространство $A^T$.

Начнём с сингулярного разложения $A = U \Sigma V^T$. Напомним, что $\Sigma$ — матрица $m \times n$, элементы которой равны нулю, кроме сингулярных чисел $\sigma_1, \sigma_2, \dots, \sigma_r$, стоящих на диагонали её первых $r$ строк. Матрицы $U$ и $V$ ортогональны, поэтому обратить их легко. Остаётся найти псевдообратную только для $\Sigma$.

Ближе всего к обратной для $\Sigma$ матрица $\Sigma^{+}$ размера $n \times m$, у которой в первых $r$ строках на диагонали стоят $1/\sigma_1, 1/\sigma_2, \dots, 1/\sigma_r$. Если $r = n = m$, то $\Sigma^{+} = \Sigma^{-1}$. Произведение $\Sigma$ и $\Sigma^{+}$ — всегда квадратная матрица, у которой первые $r$ диагональных элементов равны 1, а все остальные — 0.

Если $A = U \Sigma V^T$, то её псевдообратная

$$A^{+} = V \Sigma^{+} U^T.$$

(Вспомните, что $Q^T = Q^{-1}$ для ортогональных матриц $U$, $V$ или $Q$.)

Похожий результат получился бы, если бы мы поставили ненулевые элементы в правый нижний угол $\Sigma^{+}$, но лишние ненулевые элементы нам ни к чему.

Заключение

Хотя псевдообратных не будет на экзамене, многие темы, которые мы затронули, обсуждая их (четыре подпространства, SVD, ортогональные матрицы), там, скорее всего, появятся.

Проверьте в NumPy (дополнение переводчика)

import numpy as np

# Левая обратная для «высокой» матрицы полного ранга по столбцам (3 точки, 2 параметра прямой)
A = np.array([[1., 0.],
              [1., 1.],
              [1., 2.]])
A_left = np.linalg.inv(A.T @ A) @ A.T
print(np.round(A_left @ A, 6))                 # [[1. 0.] [0. 1.]] — слева получаем I
print(np.allclose(A_left, np.linalg.pinv(A)))  # True: здесь левая обратная = псевдообратная
b = np.array([1., 2., 2.])
print(A_left @ b)       # [1.1667 0.5] — прямая y = 1.1667 + 0.5·t по методу наименьших квадратов

# Псевдообратная для вырожденной матрицы ранга 1 из лекции 29
B = np.array([[4., 3.],
              [8., 6.]])
U, s, Vt = np.linalg.svd(B)
S_plus = np.diag([1 / s[0], 0.0])              # обращаем только ненулевое σ
B_plus = Vt.T @ S_plus @ U.T                   # A⁺ = V Σ⁺ Uᵀ
print(np.allclose(B_plus, np.linalg.pinv(B)))  # True
print(B_plus @ B)       # [[0.64 0.48] [0.48 0.36]] — проекция на пространство строк, а не I

Где это встретится в курсе (комментарий переводчика)

• Низкоранговое приближение и LoRA. Оставив в SVD только $k$ самых больших $\sigma_i$, получаем лучшее приближение матрицы ранга $k$. На этой идее построены сжатие весов и LoRA: обновление весов ищут сразу в виде произведения двух узких матриц.
• PCA (урок 2.3.2): главные компоненты — это правые сингулярные векторы $v_i$ центрированной матрицы данных.
• Линейная регрессия (урок 2.2.1): np.linalg.lstsq и np.linalg.pinv внутри считают именно SVD — так решение не ломается, даже если признаки линейно зависимы.

Словарик

EN	RU
singular value decomposition (SVD)	сингулярное разложение
singular value, singular vector	сингулярное число, сингулярный вектор
orthogonal matrix, orthonormal basis	ортогональная матрица, ортонормированный базис
row space / column space	пространство строк / пространство столбцов
nullspace / left nullspace	нуль-пространство / левое нуль-пространство
positive semidefinite	положительно полуопределённая
change of basis	смена (замена) базиса
standard basis	стандартный базис
Fourier basis, wavelet basis	базис Фурье, вейвлет-базис
lossless / lossy compression	сжатие без потерь / с потерями
similar matrices	подобные матрицы
two-sided / left / right inverse	двусторонняя / левая / правая обратная
full column rank / full row rank	полный ранг по столбцам / по строкам
pseudoinverse	псевдообратная матрица (Мура — Пенроуза)
least squares	метод наименьших квадратов

Дальше: производные по матрицам — путеводитель по MIT 18.S096 «Matrix Calculus»; назад к карте курса 18.06SC.