Методы ортогонализации матриц

Материал из MachineLearning.

Статья написана с использованием LLM GPT-5.3-mini и проверена участником ~~Dovlat Demin~~

Ортогонализация матрицы — это процесс построения по заданной матрице $A \in \mathbb{R}^{m \times n}$ (обычно $m \ge n$ ) такой матрицы $Q$ с ортонормированными столбцами, что линейная оболочка столбцов $Q$ совпадает с линейной оболочкой столбцов $A$ . Иными словами, столбцы $Q$ образуют ортонормированный базис подпространства, натянутого на столбцы исходной матрицы. Если матрица $A$ квадратная и невырожденная, задачу часто понимают как поиск ближайшей (в смысле нормы Фробениуса) ортогональной матрицы $Q$ , такой что $A = Q H$ с симметричной положительно полуопределённой матрицей $H$ — полярное разложение.

Практически важным инструментом является QR-разложение $A = QR$ , где $Q$ — матрица с ортонормированными столбцами, а $R$ — верхняя треугольная. Оно лежит в основе решения задач наименьших квадратов, вычисления собственных значений и многих алгоритмов машинного обучения.

Статья систематизирует основные методы ортогонализации, их численные свойства и применения в вычислительной математике и анализе данных.

Содержание

1 Ортогональные матрицы и их свойства
2 Постановка задачи ортогонализации
3 Основные методы ортогонализации
4 Сравнение методов
5 QR-разложение как приложение ортогонализации
6 Ортогонализация в машинном обучении
7 Современные направления и GPU-реализации
8 Заключение
9 Ссылки

Ортогональные матрицы и их свойства

Матрица $Q \in \mathbb{R}^{m \times n}$ с $m \ge n$ называется ортогональной (точнее, имеющей ортонормированные столбцы), если $Q^T Q = I_n.$ В квадратном случае ( $m=n$ ) дополнительно выполняется $Q Q^T = I$ , то есть $Q^{-1} = Q^T$ .

Ключевые свойства:

Сохранение евклидовой нормы: $\|Qx\|_2 = \|x\|_2$ для любого вектора $x$ .
Сохранение углов и скалярных произведений: $(Qx)^T (Qy) = x^T y$ .
Сохранение расстояний: $\|Qx - Qy\|_2 = \|x - y\|_2$ .
Число обусловленности: $\kappa_2(Q) = 1$ , что означает, что умножение на ортогональную матрицу не усиливает относительные ошибки.
Все сингулярные числа ортогональной матрицы равны 1.

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

Постановка задачи ортогонализации

Дана матрица $A \in \mathbb{R}^{m \times n}$ с линейно независимыми столбцами. Требуется найти:

матрицу $Q \in \mathbb{R}^{m \times n}$ с ортонормированными столбцами, такую что $\operatorname{span}(Q) = \operatorname{span}(A)$ ;
или квадратную ортогональную матрицу $Q \in \mathbb{R}^{n \times n}$ , минимизирующую $\|A - Q\|_F$ при дополнительном условии, что $A$ близка к ортогональной.

Эти задачи эквивалентны построению QR-разложения, полярного разложения или сингулярного разложения. Выбор метода диктуется требованиями к точности, вычислительной сложности и архитектуре вычислителя.

Основные методы ортогонализации

Процесс Грама–Шмидта

Классический процесс Грама–Шмидта (CGS)

Для столбцов $a_1,\ldots,a_n$ матрицы $A$ ортонормированные векторы $q_1,\ldots,q_n$ вычисляются последовательно.

Первый вектор:

$v_1=a_1,$

$q_1=\frac{v_1}{\|v_1\|_2}.$

Для $k=2,\ldots,n$ :

$v_k=a_k-\sum_{i=1}^{k-1}(q_i^Ta_k)\,q_i,$

$q_k=\frac{v_k}{\|v_k\|_2}.$

Здесь проекция вектора $a$ на направление $q_i$ определяется как

$\operatorname{proj}_{q_i}(a)=(q_i^Ta)\,q_i.$

CGS крайне чувствителен к ошибкам округления: при почти линейно зависимых столбцах вычисленные векторы $q_k$ быстро теряют ортогональность. На практике классический алгоритм используется редко.

Модифицированный процесс Грама–Шмидта (MGS)

В MGS на шаге $k$ из всех ещё не обработанных векторов $a_j\ (j>k)$ немедленно вычитается проекция на только что полученный $q_k$ : $\begin{aligned} v_k^{(0)} &= a_k,\\ v_k^{(i)} &= v_k^{(i-1)} - (q_i^T v_k^{(i-1)})\, q_i,\quad i=1,\dots,k-1,\\ q_k &= v_k^{(k-1)} / \|v_k^{(k-1)}\|_2. \end{aligned}$

MGS численно устойчивее CGS и гарантирует малое отклонение от ортогональности, сравнимое с машинным эпсилон, если матрица хорошо обусловлена. Для плохо обусловленных матриц потеря ортогональности всё же происходит, но значительно медленнее. MGS лежит в основе некоторых реализаций QR-разложения, например, в алгоритме Арнольди.

Отражения Хаусхолдера

Отражение Хаусхолдера задаётся матрицей $H = I - 2 \frac{v v^T}{v^T v},$ где $v$ — некоторый ненулевой вектор. Матрица $H$ ортогональна и симметрична. С её помощью можно обнулить все компоненты вектора, кроме первой: для заданного вектора $x$ подбирают $v = x + \operatorname{sign}(x_1) \|x\|_2 \, e_1$ , тогда $Hx = \mp \|x\|_2 \, e_1$ .

Для построения QR-разложения матрицы $A$ последовательно применяют отражения слева: на $k$ -м шаге строят матрицу $H_k$ , обнуляющую поддиагональные элементы $k$ -го столбца. После $n$ шагов $H_n \cdots H_1 A = \begin{pmatrix} R \\ 0 \end{pmatrix}, \quad Q = H_1 \cdots H_n.$

Метод Хаусхолдера обладает превосходной численной устойчивостью: вычисленный $Q$ ортогонален с точностью до машинного эпсилон, а $R$ является точным для слегка возмущённой матрицы. Это стандартный выбор для плотных матриц в библиотеках LAPACK.

Вращения Гивенса

Вращение Гивенса $G(i,j,\theta)$ действует в плоскости $(i,j)$ и имеет вид единичной матрицы с четырьмя изменёнными элементами: $G_{ii}=G_{jj}=\cos\theta$ , $G_{ij}=-G_{ji}=\sin\theta$ . Угол $\theta$ выбирается так, чтобы обнулить конкретный элемент вектора или матрицы.

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

матрица разрежена и нужно занулить лишь немногие элементы;
требуется высокая степень параллелизма (блочные алгоритмы);
алгоритм реализуется на архитектурах с ограниченной памятью.

Их численная устойчивость также высока. Вычислительная сложность примерно вдвое выше, чем у отражений Хаусхолдера, для плотных матриц, но для специальных структур она может быть значительно снижена.

Полярное разложение

Полярное разложение матрицы $A \in \mathbb{R}^{m \times n}$ ( $m \ge n$ ) — представление $A = Q H,$ где $Q \in \mathbb{R}^{m \times n}$ имеет ортонормированные столбцы, а $H \in \mathbb{R}^{n \times n}$ симметрична и положительно полуопределена. Если $A$ квадратная невырожденная, то $H$ положительно определена, а $Q$ — ортогональная матрица, являющаяся ближайшей к $A$ в норме Фробениуса (задача ортогонального Прокруста).

Полярное разложение можно вычислить:

через SVD: если $A = U \Sigma V^T$ , то $Q = U V^T$ , $H = V \Sigma V^T$ ;
итерационно, с помощью метода Ньютона–Шульца (см. ниже).

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

Ортогонализация на основе сингулярного разложения

Любую матрицу $A$ можно разложить как $A = U \Sigma V^T$ . Столбцы $U$ образуют ортонормированный базис столбцового пространства $A$ . Поэтому матрица левых сингулярных векторов $U$ непосредственно даёт искомую ортогонализацию. Для квадратной матрицы $A$ ближайшая ортогональная в норме Фробениуса равна $Q = U V^T$ .

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

Итерации Ньютона–Шульца

Для квадратной матрицы $A$ полярный фактор $Q$ можно получить итерационным методом Ньютона–Шульца: $X_{k+1} = \frac{1}{2} X_k \bigl(3I - X_k^T X_k\bigr).$ Если начальное приближение $X_0 = A$ достаточно близко к ортогональной матрице, последовательность квадратично сходится к $Q$ . На практике для ускорения сходимости и гарантии устойчивости применяют масштабирование: $X_{k+1} = \frac{1}{2} \bigl( \mu_k X_k + \mu_k^{-1} (X_k^\dagger)^T \bigr),$ или адаптивно выбирают параметр релаксации.

Главное преимущество метода — операции сводятся к матричным умножениям, которые превосходно параллелятся на GPU. Поэтому итерации Ньютона–Шульца нашли применение в современных оптимизаторах глубоких нейросетей (см. раздел про Muon). Недостаток — необходимость в хорошем начальном приближении и возможная расходимость при сильном отклонении от ортогональности.

Сравнение методов

Метод	Численная устойчивость	Вычислительная сложность	Параллельная реализация	Пригодность для GPU	Используется в QR	Используется в Глубоком обучении	Достоинства	Недостатки
Классический Грам–Шмидт (CGS)	Низкая	$2mn^2$	Плохая (сильная последовательность)	Низкая	Да (редко)	Нет	Простота реализации	Катастрофическая потеря ортогональности
Модифицированный Грам–Шмидт (MGS)	Средняя	$2mn^2$	Средняя (можно частично векторизовать)	Удовлетворительная	Да (алгоритм Арнольди)	Редко	Хорошая устойчивость для хорошо обусловленных матриц	Потеря ортогональности при плохой обусловленности
Отражения Хаусхолдера	Высокая	$2mn^2 - \frac{2}{3}n^3$	Хорошая (блочные версии)	Хорошая (LAPACK)	Стандартный алгоритм QR	Редко (инициализация весов)	Обратная устойчивость, оптимален для плотных матриц	Избыточен для разреженных структур
Вращения Гивенса	Высокая	$\approx 3mn^2 - n^3$	Хорошая (можно параллелить)	Хорошая (для разреженных/ленточных)	Да (для специальных структур)	Нет	Гибкость, избирательное обнуление	Дороже Хаусхолдера для плотных матриц
Полярное разложение (через SVD)	Высокая	$O(mn^2)$ (SVD)	Сложная	Ограниченная	Нет	Да (тонкая настройка)	Наилучшая ортогональная аппроксимация	Высокая стоимость SVD
SVD-ортогонализация	Очень высокая	$O(mn^2)$ (практически больше)	Сложная	Умеренная (существуют GPU-реализации)	Нет (но даёт QR)	Да (PCA, инициализация)	Максимальная точность, выявление ранга	Высокая стоимость, не всегда дифференцируема
Ньютон–Шульц	Средняя (локальная сходимость)	$O(n^3)$ за итерацию, общее $\sim 10 n^3$	Отличная (матричные умножения)	Превосходная	Нет	Да (Muon, ортогонализация градиентов)	Безусловная параллельность, эффективен на GPU	Требует близкого начального приближения, возможна расходимость

QR-разложение как приложение ортогонализации

QR-разложение представляет матрицу $A \in \mathbb{R}^{m \times n}$ в виде $A = QR,$ где $Q \in \mathbb{R}^{m \times m}$ ортогональная, а $R \in \mathbb{R}^{m \times n}$ верхняя треугольная (в экономичной версии $Q \in \mathbb{R}^{m \times n}$ , $R \in \mathbb{R}^{n \times n}$ ). Оно напрямую строится методами ортогонализации столбцов: Хаусхолдера, Гивенса или (модифицированным) Грамом–Шмидтом.

Основные применения:

решение линейных систем и задач наименьших квадратов $\min_x \|Ax - b\|_2$ ;
вычисление собственных значений (QR-алгоритм);
построение ортогональных базисов в подпространствах Крылова.

Численные свойства QR-разложения полностью определяются выбранным методом ортогонализации. Отражения Хаусхолдера гарантируют обратную устойчивость, поэтому именно они реализованы в большинстве стандартных библиотек.

Ортогонализация в машинном обучении

Инициализация весов и ортогональные ограничения

Ортогональная инициализация весовых матриц (Saxe et al., 2014) глубоких линейных сетей позволяет избежать проблемы затухающих/взрывающихся градиентов. Для нелинейных сетей строгая ортогональность поддерживается с помощью регуляризации $\|W^T W - I\|$ или применения специальных параметризаций (например, через экспоненту кососимметрической матрицы).

Стабилизация рекуррентных сетей

В рекуррентных нейронных сетях (RNN) ортогональные или унитарные матрицы скрытого состояния предотвращают экспоненциальный рост или затухание градиентов. Модели uRNN (Arjovsky et al., 2016), expRNN, ортогональные LSTM используют параметризацию ортогональной матрицы через произведение отражений Хаусхолдера или матричную экспоненту, а обновление весов может включать шаг ортогонализации (например, полярное разложение).

Спектральная нормализация

Спектральная нормализация (Miyato et al., 2018) ограничивает спектральную норму весовой матрицы единицей, что стабилизирует обучение генеративных состязательных сетей (GAN). Хотя сама по себе она не делает матрицу ортогональной, она тесно связана с оценкой максимального сингулярного числа, и в комбинации с другими методами может способствовать близости весов к ортогональным.

Ортогонализация в оптимизаторах: пример Muon

Современный оптимизатор Muon (Bernstein et al., 2024) применяет ортогонализацию матриц обновления весов с помощью итераций Ньютона–Шульца. Для параметра-матрицы $W \in \mathbb{R}^{m \times n}$ градиентное обновление сначала выравнивается по норме, затем к нему применяется несколько итераций $X_{k+1} = \frac{1}{2} X_k (3I - X_k^T X_k)$ до достижения почти ортогональной матрицы, которая и прибавляется к весам. Это позволяет использовать полную матричную структуру градиента и значительно ускоряет обучение больших моделей.

PCA и сингулярное разложение

Метод главных компонент (PCA) основан на SVD матрицы данных и даёт ортонормированный базис (главные компоненты), в котором дисперсия данных максимальна. Ортогонализация здесь ключевой этап, выполняемый обычно через SVD, который на больших данных аппроксимируется рандомизированными алгоритмами.

Современные направления и GPU-реализации

С ростом размеров данных и моделей акцент смещается в сторону методов, максимально использующих параллелизм GPU:

Блочные алгоритмы Хаусхолдера и Гивенса, реализованные в MAGMA, cuSOLVER, минимизируют коммуникации и эффективны для задач средней размерности.
Рандомизированная ортогонализация: сначала строится случайная проекция, а затем применяется QR/SVD малого размера; используется для приближённого PCA и ускорения обучения.
Итерации Ньютона–Шульца с адаптивным шагом активно развиваются для обучения нейросетей благодаря исключительной производительности на тензорных ядрах.
Дифференцируемая ортогонализация через полярное разложение или матричную экспоненту позволяет встраивать ортогональные ограничения непосредственно в граф вычислений, делая обучение end-to-end.

Ведутся исследования быстрой ортогонализации на сверхбольших разреженных матрицах, гибридных методов (стохастический Грам–Шмидт) и квантово-инспирированных тензорных разложений.

Заключение

Выбор метода ортогонализации определяется компромиссом между точностью, вычислительной сложностью и аппаратными возможностями. Отражения Хаусхолдера остаются золотым стандартом для плотных матриц в научных вычислениях, тогда как итерации Ньютона–Шульца находят новую жизнь в глубоком обучении. Понимание численных свойств каждого алгоритма позволяет строить надёжные и быстрые программные системы.

Ссылки

[1] Golub G. H., Van Loan C. F. *Matrix Computations*. 4th ed. Johns Hopkins University Press, 2013. [2] Trefethen L. N., Bau D. *Numerical Linear Algebra*. SIAM, 1997. [3] Higham N. J. *Functions of Matrices: Theory and Computation*. SIAM, 2008. [4] Demmel J. W. *Applied Numerical Linear Algebra*. SIAM, 1997. [5] Higham N. J. *Accuracy and Stability of Numerical Algorithms*. 2nd ed. SIAM, 2002. [6] Strang G. *Linear Algebra and Learning from Data*. Wellesley-Cambridge Press, 2019. [7] Bernstein J., Vahdat A., Yue Y., Liu M.-Y. *Muon: An optimizer for matrix parameters based on the matrix sign function*. arXiv:2406.19169, 2024. [8] Saxe A. M., McClelland J. L., Ganguli S. *Exact solutions to the nonlinear dynamics of learning in deep linear neural networks*. ICLR 2014. [9] Arjovsky M., Shah A., Bengio Y. *Unitary Evolution Recurrent Neural Networks*. ICML 2016. [10] Miyato T., Kataoka T., Koyama M., Yoshida Y. *Spectral Normalization for Generative Adversarial Networks*. ICLR 2018. [11] Higham N. J. *Computing the polar decomposition—with applications*. SIAM J. Sci. Stat. Comput., 7(4), 1160–1174, 1986. [12] Strang G. *Introduction to Linear Algebra*. 5th ed. Wellesley-Cambridge Press, 2016. [13] Horn R. A., Johnson C. R. *Matrix Analysis*. 2nd ed. Cambridge University Press, 2012.

Источник — «http://www.machinelearning.ru/wiki/index.php?title=%D0%9C%D0%B5%D1%82%D0%BE%D0%B4%D1%8B_%D0%BE%D1%80%D1%82%D0%BE%D0%B3%D0%BE%D0%BD%D0%B0%D0%BB%D0%B8%D0%B7%D0%B0%D1%86%D0%B8%D0%B8_%D0%BC%D0%B0%D1%82%D1%80%D0%B8%D1%86»

Категории: Линейная алгебра | Численные методы | Машинное обучение