Линейная алгебра
Матричные функции: eᴬ и её родственники
torch.matrix_exp - реальная функция в PyTorch, используемая в Neural ODE. Матричная экспонента eᴬ - не экзотика из учебника, а инструмент непрерывных RNN и нормализующих потоков. В scipy она используется в медицинской визуализации МРТ и системах управления ежедневно.
- Neural ODE (Chen et al. 2018): решение dh/dt = f(h,t) через torch.matrix_exp - параметров в 10x меньше ResNet
- Управление: e^(At) - матрица перехода системы dx/dt = Ax, аналитическое решение
- Квантовые вычисления: эволюция состояния кубита - e^(-iHt/ℏ), H - гамильтониан
- Диффузионная МРТ: матричный логарифм тензора диффузии для интерполяции между вокселями
- Нормализующие потоки: параметризация ортогональных матриц через e^(A-Aᵀ)
**PyTorch содержит функцию `torch.matrix_exp`.** Она вычисляет матричную экспоненту eᴬ - и именно она лежит в основе Neural ODE, непрерывных RNN и нескольких семейств нормализующих потоков. Матричная экспонента - это не экзотика из учебника, это инструмент из production-кода. `scipy.linalg.expm` вызывается в коде медицинских систем анализа МРТ, рекуррентных сетей и генеративных моделей ежедневно.
**О чём этот урок на самом деле**: как функция числа f(x) превращается в функцию матрицы f(A). Это не метафора - это строгая конструкция через степенные ряды и спектральное разложение. И у неё есть три прямых применения в ML: матричная экспонента в Neural ODE, softmax как нормировка exp, логарифм матрицы в геометрии ковариационных матриц.
Как функция числа становится функцией матрицы
Для числа x функция eˣ определяется рядом Тейлора: eˣ = 1 + x + x²/2! + x³/3! + ... Ключевое наблюдение: в этом ряду используются только **сложение** и **умножение** - операции, которые определены и для матриц. Подставляем матрицу A вместо числа x:
**Контр-интуиция**: для чисел eᵃ · eᵇ = eᵃ⁺ᵇ всегда. Для матриц eᴬ · eᴮ = eᴬ⁺ᴮ **только если AB = BA**. Если A и B не коммутируют - формула неверна. Это источник половины ошибок при работе с матричными экспонентами.
Диагонализация: быстрый путь к f(A)
Если матрица A диагонализируема, то A = PΛP⁻¹, где Λ - диагональная матрица собственных значений. Тогда A² = PΛ²P⁻¹, A³ = PΛ³P⁻¹, и ряд сворачивается красиво:
Neural ODE: непрерывная глубина сети
Обычная нейросеть - это дискретная последовательность слоёв: h₁ = f₁(h₀), h₂ = f₂(h₁), ... Neural ODE (Chen et al., 2018, NeurIPS) делает глубину **непрерывной**: вместо дискретных слоёв задаётся дифференциальное уравнение dh/dt = f(h, t, θ). Тогда выход сети - это h(T), которое вычисляется ODE-солвером.
**Производные в Neural ODE** считаются adjoint-методом - без хранения промежуточных активаций. Это даёт O(1) памяти по глубине - главное преимущество перед классическими ResNet. Статья «Neural Ordinary Differential Equations» набрала 8000+ цитирований за 5 лет.
Softmax: нормировка экспоненты в действии
Softmax - это **нормировка** функции exp применённой к вектору. Связь с матричной экспонентой прямая: softmax(x)ᵢ = eˣⁱ / Σⱼ eˣʲ. Это не матричная экспонента - но та же операция exp, только над вектором, нормированная до суммы 1. В трансформерах attention score = softmax(QᵀK / √d) - именно так, через softmax каждого строчного вектора произведения матриц.
Логарифм матрицы: геодезика в пространстве SPD
~Матричный логарифм~{log(A) - матрица X такая, что eˣ = A} - обратная операция к матричной экспоненте. Это не просто математическая абстракция: **в медицинской визуализации (DTI МРТ)** диффузионные тензоры - это симметричные положительно определённые матрицы (SPD). Обычная средняя арифметическая двух SPD матриц не является SPD. Геодезическое среднее - является.
**Зачем это в ML**: в метрическом обучении (metric learning) с ковариационными дескрипторами - например, для анализа текстур и лиц (Riemannian SPD networks) - именно логарифм матрицы используется для «выпрямления» многообразия SPD в евклидово пространство. Пакет `geomstats` в Python реализует это для нейросетей.
Матричный квадратный корень и whitening
~Матричный квадратный корень~{A^(1/2): матрица X такая, что X·X = A} применяется в **whitening** - процедуре, которая делает признаки некоррелированными и равномерно масштабированными. Это стандартный шаг перед PCA и в некоторых архитектурах нормализации.
Где встречаются матричные функции
Практика: матричная экспонента
Как функция числа становится функцией матрицы
Для диагонализируемой матрицы $A = V\Lambda V^{-1}$ функция определяется покоординатно: $f(A) = V \cdot \text{diag}(f(\lambda_1), \ldots, f(\lambda_n)) \cdot V^{-1}$. Для неквадратных или несимметричных матриц используют разложение Жордана или определение через степенной ряд: $f(A) = \sum_{k=0}^{\infty} \frac{f^{(k)}(0)}{k!} A^k$.
**Матричная экспонента** $e^A = I + A + A^2/2! + A^3/3! + \ldots$ - наиболее важная матричная функция. Для диагонализируемой $A$: $e^A = Ve^\Lambda V^{-1}$, где $e^\Lambda = \text{diag}(e^{\lambda_1}, \ldots, e^{\lambda_n})$. Ключевое свойство: $e^{tA}$ - решение системы $\dot{x} = Ax$: $x(t) = e^{tA} x(0)$.
Neural ODE: слой нейросети как ODE
torch.matrix_exp в продакшне
Neural ODE (Chen et al., 2018): вместо дискретных слоёв $h_{t+1} = f(h_t, W)$ используют непрерывную динамику $dh/dt = f(h(t), W)$. Решение: $h(T) = e^{WT} h(0)$ для линейного случая. В PyTorch: `torch.matrix_exp(W * T) @ h_0`. Параметров в 10x меньше чем в ResNet при сопоставимом качестве на MNIST.
Как вычисляется $e^A$ для диагонализируемой матрицы $A = V\Lambda V^{-1}$?
Для $A = V\Lambda V^{-1}$: $e^A = Ve^\Lambda V^{-1}$, где $e^\Lambda$ - диагональная матрица с $e^{\lambda_i}$ на диагонали. Это следует из степенного ряда: $A^k = V\Lambda^k V^{-1}$, поэтому $\sum A^k/k! = V(\sum \Lambda^k/k!)V^{-1} = Ve^\Lambda V^{-1}$.
Матричный логарифм и геометрия ковариационных матриц
Логарифм матрицы $\log A$ - обратная к экспоненте: $e^{\log A} = A$. Для $A \succ 0$: $\log A = V \cdot \text{diag}(\log \lambda_1, \ldots, \log \lambda_n) \cdot V^\top$. Существует только при положительных собственных значениях.
**Геодезика между ковариационными матрицами**: пространство SPD-матриц ($A \succ 0$) - риманово многообразие. Среднее двух матриц $\Sigma_1, \Sigma_2$ по кратчайшему пути (геодезике) - не $\frac{1}{2}(\Sigma_1 + \Sigma_2)$, а $\Sigma_1^{1/2}(\Sigma_1^{-1/2}\Sigma_2\Sigma_1^{-1/2})^{1/2}\Sigma_1^{1/2}$. Это используется в анализе МРТ (диффузионная тензорная томография) и ковариационных нейросетях.
Softmax как нормировка матричной экспоненты
Связь с attention
Softmax: $\text{softmax}(z)_i = e^{z_i} / \sum_j e^{z_j}$. Для вектора $z$ это нормировка покоординатной экспоненты. В attention: $\text{softmax}(QK^\top/\sqrt{d})$ применяет это к каждой строке матрицы. Матричная экспонента $e^A$ применяет полноценную матричную версию - учитывает взаимодействие между компонентами через $A^2, A^3, \ldots$.
При каком условии матричный логарифм $\log A$ существует и единственен?
$\log A = V \cdot \text{diag}(\log \lambda_i) \cdot V^{-1}$. Логарифм вещественного числа определён только для положительных чисел, поэтому требуем $\lambda_i > 0$. Для матриц с отрицательными или комплексными собственными значениями логарифм либо комплексен, либо неединственен.
Для скаляра: eˣ = I·1 + x + x²/2! + x³/3! + ... Для матрицы: eᴬ = I + A + A²/2! + A³/3! + ... Где I - единичная матрица (аналог числа 1). Ряд сходится для ЛЮБОЙ квадратной матрицы A. Пример: A = [[0, -1], [1, 0]] (поворот на t=π/4) A¹ = [[0, -1], [1, 0]] A² = [[-1, 0], [0, -1]] = -I A³ = [[0, 1], [-1, 0]] A⁴ = [[1, 0], [0, 1]] = +I (цикл длины 4) eᴬᵗ = [[cos t, -sin t], - матрица поворота на угол t! [sin t, cos t]]
Если A = P·diag(λ₁, ..., λₙ)·P⁻¹, то: f(A) = P · diag(f(λ₁), ..., f(λₙ)) · P⁻¹ Пример: eᴬ для A с собственными значениями -1, -2 eᴬ = P · [[e⁻¹, 0 ], · P⁻¹ [ 0, e⁻²]] Интуиция: функция от матрицы = применить функцию к каждому собственному значению, оставив собственные векторы нетронутыми. Для недиагонализируемых: разложение Шура A = QTQ*, где T - верхнетреугольная. scipy.linalg.expm использует именно это + аппроксимацию Паде.
Для линейного случая f(h, t) = A·h: dh/dt = A·h, h(0) = h₀ Решение: h(t) = eᴬᵗ · h₀ Один «слой» Neural ODE с глубиной T: output = eᴬᵀ · input Параметры: матрица A (n² чисел вместо n² на каждый слой). Глубина T - гиперпараметр, не влияющий на число параметров. → Непрерывная глубина сети без взрыва числа параметров.
Пространство SPD матриц - риманово многообразие. Геодезическая кривая между A и B: γ(t) = A^(1/2) · (A^(-1/2) · B · A^(-1/2))^t · A^(1/2) Среднее (t = 0.5): M = A^(1/2) · exp(0.5 · log(A^(-1/2) · B · A^(-1/2))) · A^(1/2) Пример: A и B - ковариационные матрицы двух групп пациентов. Междиагностическое среднее M - более корректная оценка, чем (A + B)/2, которое может не быть SPD. Применяется: FSL, FreeSurfer, DiPy - стандартные пакеты нейровизуализации.
- **Neural ODE / непрерывные RNN**: torch.matrix_exp(A * t) как один «слой» любой глубины
- **Softmax / attention в трансформерах**: exp(x) + нормировка = softmax; одна операция во всех LLM
- **Riemannian SPD networks**: logm(A) для «выпрямления» SPD многообразия
- **Whitening / ZCA**: C^(-1/2) через sqrtm для декорреляции признаков
- **Normalizing flows (FFJORD)**: Непрерывный поток через Neural ODE: log det вычисляется через trace
Упражнения
- Чем вычисление eᴬ через ряд Тейлора отличается от scipy.linalg.expm? Почему прямое суммирование ряда численно неустойчиво? — Прямое суммирование: при больших ||A|| члены ряда сначала растут, затем падают - catastrofic cancellation; Expm использует Pade approximation + scaling and squaring: A → A/2^s, вычислить eᴬ/²ˢ, затем возвести в квадрат s раз; scipy.linalg.expm реализует алгоритм Al-Mohy-Higham 2009 - оптимальный по числу умножений матриц
- Почему eᴬ · eᴮ ≠ eᴬ⁺ᴮ в общем случае? Когда равенство выполняется? — eᴬ⁺ᴮ = eᴬ · eᴮ только если AB = BA (матрицы коммутируют); Общий случай: eᴬ⁺ᴮ = eᴬ · eᴮ · e^(-[A,B]/2) · ... (формула BCH); Именно из-за некоммутативности квантовая механика так сложна - операторы H и p не коммутируют
- Softmax - это частный случай матричной экспоненты или отдельная операция? — Softmax - это exp поэлементно к вектору + нормировка, не матричная экспонента; Матричная экспонента - нелинейная операция над матрицей как целым объектом; Формальная связь: softmax(x) = диагональ нормированного exp от диагональной матрицы с x на диагонали
Что унести из урока
- **eᴬ = I + A + A²/2! + ...** - ряд Тейлора, сходится для любой матрицы
- **Для диагонализируемой A = PΛP⁻¹**: f(A) = P·diag(f(λᵢ))·P⁻¹ - функция от матрицы = функция от собственных значений
- **Neural ODE**: dh/dt = A·h → h(t) = eᴬᵗ·h₀; непрерывная глубина за счёт одной матрицы A
- **Softmax** = exp поэлементно + нормировка; одна операция в каждом attention head
- **log(A)** - геодезика в SPD пространстве; стандарт в нейровизуализации (DTI) и Riemannian DL
- **eᴬ · eᴮ = eᴬ⁺ᴮ только если AB = BA** - ловушка, которая ломает интуицию из скалярного случая
Куда дальше
Матричные функции - мост между линейной алгеброй и квантовой механикой, где каждая эволюция системы - это eᴬ.
- Линейная алгебра в квантовой механике — Квантовая эволюция U(t) = e^(-iHt/ℏ) - матричная экспонента антиэрмитова оператора
- Собственные векторы и диагонализация — f(A) через спектральное разложение требует диагонализации или разложения Шура
- SVD — Полярное разложение A = UP через sqrtm(AᵀA) - ещё одно применение матричного корня