Собственные векторы и разложения

Направления, которые не поворачиваются

Посмотрите на сдвиг $\begin{pmatrix} 1 & 1 \\ 0 & 1 \end{pmatrix}$ или растяжение в виджете из прошлых уроков. Почти все векторы после преобразования поворачиваются. Но некоторые остаются на своей прямой — только растягиваются, сжимаются или разворачиваются. Это и есть собственные векторы.

Собственный вектор $v \ne 0$ матрицы $A$ и его собственное значение $\lambda$ :

Av = \lambda v

Собственные векторы — «оси» преобразования, а собственные значения — коэффициенты растяжения вдоль этих осей.

import numpy as np
import vis

A = np.array([[3.0, 1], [0, 2]])
vals, vecs = np.linalg.eig(A)
vals, vecs = vals.real, vecs.real      # eig может вернуть комплексный тип; здесь мнимые части нулевые
print("собственные значения:", vals)
print("собственные векторы (в столбцах):")
print(vecs.round(3))

for i in range(2):
    v = vecs[:, i]
    print(f"A·v{i + 1} = {A @ v} = {vals[i]:.0f}·v{i + 1}")

vis.transform(A, title="Векторы на собственных прямых только растягиваются", vecs=[vecs[:, 0], vecs[:, 1]])

На графике розовые стрелки — образы собственных векторов: они лежат на тех же прямых, что и исходные.

Что говорят собственные значения

$\lambda > 1$ — растяжение вдоль этого направления, $0 < \lambda < 1$ — сжатие.
$\lambda < 0$ — разворот в противоположную сторону.
$\lambda = 0$ — направление схлопывается: матрица вырождена.
Произведение всех собственных значений равно определителю, а сумма — следу (сумме диагонали).

У поворота на 90° вещественных собственных векторов нет: каждое направление поворачивается.

Многократное применение

Собственные значения объясняют, что происходит при многократном умножении на матрицу: компонента вдоль направления с $|\lambda| > 1$ растёт как $\lambda^n$ , с $|\lambda| < 1$ — затухает. Это прямо связано с проблемой взрывающихся и затухающих градиентов в глубоких сетях: градиент на каждом слое умножается на матрицу, и если масштабы не контролировать, через 50 слоёв он становится огромным или исчезает.

import numpy as np
import vis

A = np.array([[1.1, 0], [0, 0.9]])
v = np.array([1.0, 1.0])
xs, ys = [], []
for _ in range(40):
    v = A @ v
    xs.append(v[0]); ys.append(v[1])
vis.line({"компонента с λ = 1.1": xs, "компонента с λ = 0.9": ys}, title="40 умножений на матрицу", xlabel="шаг")
print(f"1.1⁴⁰ = {1.1 ** 40:.1f},  0.9⁴⁰ = {0.9 ** 40:.4f}")

Степенной метод

Из этого же наблюдения получается простой алгоритм: если много раз умножать любой вектор на $A$ и нормировать, он повернётся к собственному вектору с наибольшим $|\lambda|$ — остальные компоненты затухнут относительно него. Так (в усложнённом виде) Google считал PageRank. Вы реализуете степенной метод в задании.

Симметричные матрицы

Матрица симметрична, если $A = A^T$ . Такие матрицы в ML повсюду: ковариационная матрица, $X^TX$ , матрица расстояний, гессиан. И у них отличные свойства:

все собственные значения вещественные;
собственные векторы можно выбрать взаимно перпендикулярными (ортонормированными).

Значит, симметричная матрица — это просто растяжение вдоль перпендикулярных осей: $A = Q \Lambda Q^T$ , где столбцы $Q$ — собственные векторы, $\Lambda$ — диагональная матрица собственных значений. Для симметричных матриц используйте np.linalg.eigh — она быстрее и возвращает значения по возрастанию.

import numpy as np

A = np.array([[2.0, 1], [1, 3]])
vals, Q = np.linalg.eigh(A)
print("λ =", vals.round(4))
print("Qᵀ·Q =\n", (Q.T @ Q).round(10))                 # единичная: векторы ортонормированы
print("Q·Λ·Qᵀ = A:", np.allclose(Q @ np.diag(vals) @ Q.T, A))

PCA: метод главных компонент

Данные часто имеют много признаков, но «живут» вблизи пространства меньшей размерности. PCA находит направления наибольшего разброса данных и проецирует данные на них.

Центрируем данные: вычитаем среднее каждого признака.
Считаем ковариационную матрицу $C = \frac{1}{N-1} X_c^T X_c$ (d×d). Её элемент $c_{ij}$ показывает, как признаки $i$ и $j$ меняются вместе.
Собственные векторы $C$ — главные компоненты, собственные значения — дисперсия данных вдоль них.
Проецируем данные на $k$ компонент с наибольшими значениями.

import numpy as np
import vis

rng = np.random.default_rng(1)
X = rng.normal(size=(300, 2)) @ np.array([[3.0, 1.2], [0, 0.6]])   # вытянутое облако

Xc = X - X.mean(axis=0)
C = Xc.T @ Xc / (len(X) - 1)
vals, vecs = np.linalg.eigh(C)
print("ковариационная матрица:\n", C.round(2))
print("дисперсия вдоль компонент:", vals.round(2))

vis.scatter(Xc, title="Облако точек", equal=True)
vis.vectors([vecs[:, 1] * np.sqrt(vals[1]) * 2, vecs[:, 0] * np.sqrt(vals[0]) * 2],
            labels=["1-я компонента", "2-я"], title="Главные компоненты (длина ~ разброс)", lim=7)

Сжимаем цифры из 64 измерений в 2

Набор рукописных цифр: каждая картинка 8×8 — вектор из 64 чисел. Спроецируем все 1797 картинок на две главные компоненты и раскрасим по настоящей цифре. Метки PCA не видит — и всё равно похожие цифры собираются вместе.

import numpy as np
import ml_data
import vis

X, y = ml_data.load_digits()
vis.images(X[:10].reshape(-1, 8, 8), titles=[str(int(t)) for t in y[:10]])

Xc = X - X.mean(axis=0)
vals, vecs = np.linalg.eigh(Xc.T @ Xc / (len(X) - 1))
order = np.argsort(vals)[::-1]
vals, vecs = vals[order], vecs[:, order]

explained = vals / vals.sum()
print(f"первые 2 компоненты объясняют {explained[:2].sum():.0%} разброса, первые 10 — {explained[:10].sum():.0%}")
vis.line(np.cumsum(explained)[:30], title="Доля объяснённого разброса", xlabel="число компонент")

Z = Xc @ vecs[:, :2]
keep = np.isin(y, [0, 1, 4, 6])        # четыре цифры, чтобы было нагляднее
colors = np.searchsorted([0, 1, 4, 6], y[keep])   # цвет: 0 — синий, 1 — оранжевый, 4 — зелёный, 6 — розовый
vis.scatter(Z[keep], labels=colors, title="Цифры 0, 1, 4, 6 в двух главных компонентах")
vis.images(vecs[:, :6].T.reshape(-1, 8, 8), titles=[f"PC{i + 1}" for i in range(6)], diverging=True)

Последняя строка картинок — сами главные компоненты, «базовые узоры», из которых складываются цифры.

SVD — разложение любой матрицы

Собственное разложение есть не у всех матриц, а сингулярное разложение (SVD) — у любой, даже прямоугольной:

A = U \Sigma V^T

Любое линейное преобразование — это поворот ( $V^T$ ), растяжение вдоль осей ( $\Sigma$ ) и ещё один поворот ( $U$ ). Числа на диагонали $\Sigma$ — сингулярные значения. Отбросив маленькие, получаем лучшее приближение матрицы меньшего ранга — так сжимают данные и ускоряют большие модели (метод LoRA при дообучении языковых моделей основан на той же идее малого ранга).

import numpy as np
import ml_data
import vis

X, y = ml_data.load_digits()
img = X[0].reshape(8, 8)
U, S, Vt = np.linalg.svd(img)
print("сингулярные значения:", S.round(1))
approx = [U[:, :r] @ np.diag(S[:r]) @ Vt[:r] for r in (1, 2, 3, 8)]
vis.images([img] + approx, titles=["оригинал", "ранг 1", "ранг 2", "ранг 3", "ранг 8"])

Итог

$Av = \lambda v$ : собственный вектор лишь растягивается в $\lambda$ раз, не поворачиваясь.
Многократное умножение усиливает направления с $|\lambda| > 1$ и гасит с $|\lambda| < 1$ — отсюда взрыв и затухание градиентов и степенной метод.
Симметричные матрицы: вещественные $\lambda$ , ортогональные собственные векторы, $A = Q\Lambda Q^T$ .
PCA: собственные векторы ковариационной матрицы — направления наибольшего разброса; позволяет сжать и визуализировать данные.
SVD $A = U\Sigma V^T$ есть у любой матрицы и даёт лучшие приближения малого ранга.