В прошлом уроке про матрицы мы видели: у симметричной матрицы есть ортогональные собственные направления. Сейчас окажется, что именно они отвечают на важнейший вопрос анализа данных — вдоль каких осей облако точек разбросано сильнее всего. Так рождается метод главных компонент: не магия, а прямое следствие ортогональности и множителей Лагранжа.

Проекция на единичный вектор

Пусть у нас nn объектов, у каждого — вектор признаков xix_i. Спроецируем их на единичное направление uu: проекция ii-го объекта равна числу zi=xiuz_i=x_i^\top u. Насколько «широко» точки разъезжаются вдоль этой оси, измеряет выборочная дисперсия проекций. Если данные центрированы (среднее вычтено), она равна

1nizi2=u ⁣(1nXcXc) ⁣u=uSu,S=1nXcXc.\frac1n\sum_i z_i^2=u^\top\!\left(\frac1n X_c^\top X_c\right)\!u=u^\top S u, \qquad S=\frac1n X_c^\top X_c.

Здесь SS — ковариационная матрица. Первая главная компонента — это направление, вдоль которого разброс максимален:

u1=argmaxu=1uSu.u_1=\arg\max_{\lVert u\rVert=1} u^\top S u.

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

Рисунок шире экрана — проведите по немуОткрыть целиком ↗
Слева облако центрированных ирисов по длине и ширине лепестка, красная ось первой компоненты вдоль вытянутого облака. Справа кривая дисперсии проекции от угла оси: она поднимается к максимуму 3,66 около 23 градусов и падает почти до нуля перпендикулярно
Рис. 38.1. Дисперсия проекции зависит от направления

Каждому углу оси проекции отвечает своя дисперсия uSuu^\top S u. Максимум 3,663{,}66 достигается вдоль вытянутой оси облака — это и есть первая главная компонента.

Центрировать обязательно

Вычитание среднего — не формальность. Без него первая ось может просто указывать из начала координат к далёкому облаку, описывая его положение, а не внутренний разброс.

Вращайте ось над реальным облаком ирисов и следите за столбиком дисперсии справа: он максимален, когда ось совпадает с длинной осью облака. Переключатель центрирования по-настоящему пересчитывает данные — выключите его и увидите, как ось теряет смысл, уводя к началу координат.

От Лагранжа к собственному вектору

Откуда берётся оптимальное направление? Это задача на условный экстремум: максимизировать uSuu^\top S u при ограничении uu=1u^\top u=1. Составим лагранжиан

L(u,λ)=uSuλ(uu1).\mathcal L(u,\lambda)=u^\top S u-\lambda\,(u^\top u-1).

Приравняв производную к нулю, получаем 2Su2λu=02Su-2\lambda u=0, то есть

Su=λu.S u=\lambda u.

Оптимальное направление — собственный вектор ковариационной матрицы, а само λ\lambda равно дисперсии проекции вдоль него. Максимум достигается на собственном векторе с наибольшим собственным значением. Никакого отдельного алгоритма: PCA — это задача на собственные значения.

Собственные векторы ковариации ортогональны

Ковариационная матрица симметрична, а у симметричной матрицы собственные векторы, отвечающие разным собственным значениям, взаимно перпендикулярны. Поэтому главные компоненты образуют ортонормированный базис: вторая ось перпендикулярна первой, третья — первым двум, и так далее.

Рисунок шире экрана — проведите по немуОткрыть целиком ↗
Облако ирисов трёх видов с эллипсом рассеяния. Красная длинная стрелка — первая компонента с собственным значением 3,66, фиолетовая короткая перпендикулярная стрелка — вторая компонента с собственным значением 0,04; оси строго перпендикулярны
Рис. 38.2. Главные оси — ортогональные собственные векторы ковариации

Эллипс рассеяния и две его оси. Длинная — первая компонента (λ1=3,66\lambda_1=3{,}66), короткая перпендикулярная — вторая (λ2=0,04\lambda_2=0{,}04). Длины полуосей равны λ\sqrt{\lambda}.

Русская линия здесь ведёт к самому понятию ортогональности в статистике. П. Л. Чебышёв и его школа развили теорию ортогональных многочленов и метод наименьших квадратов, где взаимно перпендикулярные направления позволяют раскладывать величину на независимые вклады. Главные компоненты — это ровно такой ортогональный базис, только построенный по разбросу самих данных.

Следующие компоненты и SVD

Вторая компонента максимизирует дисперсию среди направлений, перпендикулярных первой, третья — перпендикулярных первым двум. На практике компоненты почти никогда не ищут через явную ковариацию: вместо этого раскладывают саму центрированную матрицу данных сингулярным разложением

Xc=UΣV.X_c=U\Sigma V^\top.

Столбцы VV — это и есть главные компоненты, а собственные значения ковариации связаны с сингулярными числами: λj=σj2/n\lambda_j=\sigma_j^2/n. Такой путь численно устойчивее, потому что не требует возводить данные в квадрат при построении SS.

Спектр и объяснённая дисперсия

Как понять, сколько компонент оставить? Собственные значения, выстроенные по убыванию, образуют спектр (scree plot): первые несколько велики, дальше — быстрый спад. Накопленная доля объяснённой дисперсии

Rm=λ1++λmλ1++λdR_m=\frac{\lambda_1+\dots+\lambda_m}{\lambda_1+\dots+\lambda_d}

говорит, какую долю разброса удерживают первые mm осей. На реальных рукописных цифрах 8×88\times8 уже 1313 компонент из 6464 удерживают 80%80\% дисперсии.

Рисунок шире экрана — проведите по немуОткрыть целиком ↗
Слева спектр собственных значений цифр в логарифмической шкале круто падает, синяя накопленная дисперсия растёт до 100 процентов и достигает 80 процентов при 13 компонентах. Справа реконструкция цифры четыре по 2, 8, 16 и 40 компонентам, от размытой к чёткой
Рис. 38.3. Спектр компонент и реконструкция цифры

Слева: спектр (красный) круто падает, накопленная дисперсия (синяя) достигает 80%80\% уже при 1313 компонентах. Справа: цифра, восстановленная по mm компонентам, от размытой к чёткой.