Рандомизированная линейная алгебра

Матрица данных 1 миллион строк на 10 000 столбцов не помещается в оперативную память для полного SVD. Рандомизированный SVD (Halko et al., 2011) аппроксимирует первые k сингулярных векторов за O(mnk) вместо O(mn min(m,n)). Именно его использует sklearn.decomposition.TruncatedSVD.

• sklearn TruncatedSVD: внутри randomized SVD Halko 2011 - стандарт для разреженных матриц
• Обработка геномных данных: матрица клетка-ген 100K x 30K, rSVD для PCA за минуты
• Рекомендательные системы: аппроксимация матрицы рейтингов через random projection + SVD
• Johnson-Lindenstrauss: случайная проекция сохраняет расстояния - используется в LSH
• Нейросети: sketching матрицы Гессе для анализа кривизны без полного вычисления

Лемма Джонсона-Линденштраусса

**sklearn.utils.extmath.randomized_svd** - эта функция с 2011 года работает внутри TruncatedSVD, PCA, NMF в scikit-learn. На матрице 50 000 × 10 000 она даёт результат за 3 секунды; точное SVD - за 5 минут. Погрешность - меньше машинной точности для первых k компонент. Как случайность стала точным инструментом - суть этого урока.

**О чём этот урок**: точное SVD матрицы 100K×50K занимает часы. Но часто нужны только первые k сингулярных векторов - для PCA, рекомендаций, NLP. Рандомизированная линейная алгебра даёт 99.9% точности за O(mn·k) вместо O(mn²), превращая случайность из проблемы в вычислительный ресурс.

Лемма Джонсона-Линденштраусса: размерность не важна

Интуиция говорит: чтобы сохранить расстояния между точками при проекции, нужно пространство почти такой же размерности. **Лемма Джонсона-Линденштраусса (JL) говорит обратное**: n точек из R^d можно вложить в R^k при k = O(log n / ε²) с сохранением всех попарных расстояний с точностью ε - и k не зависит от d.

ДЛЯ n точек в R^d и ε > 0 существует линейное отображение f: R^d → R^k, где k = O(log n / ε²), такое что для всех пар (u, v): (1 - ε)·||u - v||² ≤ ||f(u) - f(v)||² ≤ (1 + ε)·||u - v||² ПОРАЗИТЕЛЬНЫЙ ФАКТ: k не зависит от d! ЧИСЛА: n = 10 000 документов, ε = 0.1 (10% погрешность) k ≈ 2·ln(10 000) / 0.01 ≈ 1 840 Исходная размерность d = 100 000 (TF-IDF вектор) Целевая размерность k = 1 840 (в 54 раза меньше) Все попарные расстояния сохранены с точностью 10%

**Где JL применяется прямо сейчас**: Google использует random projections для Locality Sensitive Hashing (LSH) - мгновенный поиск похожих документов в миллиардах. Random Fourier Features (RFF) аппроксимируют ядро RBF через проекцию - ускоряют kernel SVM от O(n³) до O(n·k). Compressed sensing в МРТ восстанавливает полный снимок из в 10 раз меньшего числа измерений.

Рандомизированный SVD

Randomized SVD: 100x быстрее за те же результаты

Точное SVD матрицы m×n стоит O(m·n·min(m,n)) операций. Если нужны только k первых сингулярных векторов (k << min(m,n)), алгоритм Halko-Martinsson-Tropp (HMT, 2011) вычисляет их за O(m·n·k) - в min(m,n)/k раз быстрее. Это стандартный алгоритм внутри sklearn и SciPy.

ВХОД: матрица A (m×n), целевой ранг k, oversampling p (обычно 10) ШАГ 1: Случайная проекция Ω ~ N(0,1) размера n×(k+p) Y = A·Ω [m×(k+p)] ШАГ 2: QR-разложение Y Y = Q·R Q - ортонормированный базис ШАГ 3: Проекция A на малое пространство B = Qᵀ·A [(k+p)×n] ШАГ 4: Точное SVD малой матрицы B = Û·Σ·Vᵀ дёшево: O((k+p)²·n) ШАГ 5: Восстановление левых сингулярных векторов U = Q·Û СЛОЖНОСТЬ: O(m·n·(k+p)) vs O(m·n²) для полного SVD OVERSAMPLING p=10: уменьшает вероятность промаха без oversampling - ошибка может быть большой с p=10 - ошибка экспоненциально убывает

**Почему oversampling работает**: случайная матрица Ω с вероятностью почти 1 проецирует k главных направлений A в Y. Но маленький промах в направлении может дать плохую аппроксимацию. Дополнительные p=10 столбцов резко снижают эту вероятность - ошибка убывает как exp(-p). Итерационные степенные итерации (параметр n_iter) дополнительно улучшают качество для матриц с медленным убыванием сингулярных значений.

Sketching и Random Fourier Features

Sketching: решение задач на части данных

Линейная регрессия на 10 миллионах точек - решать её точно стоит O(10⁷ · d²). **Sketch-and-solve**: умножить матрицу A и вектор b на случайную матрицу S размера s×m (s << m), решить меньшую систему за O(s·d²). С правильным выбором S получается (1+ε)-приближение к оптимуму.

**Leverage scores - не случайная субвыборка**: равномерная субвыборка строк матрицы даёт плохие результаты когда строки неоднородны по 'важности'. Leverage score строки i = ||U[i,:]||², где U - левые сингулярные векторы. Строки с высоким leverage score 'держат' задачу - их нужно брать с большей вероятностью. Оптимальная субвыборка по leverage scores даёт (1+ε)-приближение с O(d log d / ε²) строками.

Random Fourier Features: ядра без матриц Грама

Kernel SVM с RBF-ядром строит матрицу Грама n×n - при n=100 000 это 10 ГБ памяти и O(n²) вычислений. **Random Fourier Features** (Rahimi & Recht, 2007): аппроксимировать ядро k(x,y) = exp(-||x-y||²/2σ²) через явное низкоразмерное вложение φ(x) ∈ R^D такое что k(x,y) ≈ φ(x)ᵀ·φ(y). Это переводит задачу из O(n²) к O(n·D) - и SVM превращается в обычную линейную регрессию.

**Рандомизированные методы в ML-инфраструктуре** Где случайность ускоряет production-системы - **sklearn TruncatedSVD / PCA**: randomized_svd (HMT) по умолчанию. Для матриц > 500x500 sklearn автоматически переключается на рандомизированный алгоритм. 10-100x ускорение при k << min(m,n). - **LSH поиск (Google, Pinecone)**: JL проекции для locality-sensitive hashing. Миллиарды документов - мгновенный поиск похожих через случайные гиперплоскости. Основа SimHash для веб-дедупликации. - **Kernel SVM → RFF → линейная модель**: Random Fourier Features для RBF ядра. sklearn.kernel_approximation.RBFSampler. O(n*D) вместо O(n²). Стандарт для больших датасетов без GPU. - **Compressed Sensing (МРТ, CS camera)**: Восстановление разреженного сигнала из случайных измерений. МРТ-снимок из в 10 раз меньшего числа измерений через JL + L1-минимизацию. Сокращает время сканирования в 10x.

Практика: случайная проекция

**Почему randomized SVD быстрее точного при k << min(m,n), но не при k = min(m,n)?** Hints: Какова сложность каждого шага HMT?; При k = min(m,n) что происходит с матрицей B? - HMT стоит O(m·n·k) - при k = min(m,n) это то же самое что O(m·n²) для точного - Ускорение пропорционально min(m,n)/k - только при k << min(m,n) выгодно - Для k = 50 из матрицы 10000×5000: HMT в 100x быстрее - Для полного SVD рандомизация не помогает - нет сокращения размерности --- **Лемма JL сохраняет расстояния. В каких ML-задачах это достаточно, а в каких - нет?** Hints: Какие алгоритмы используют только расстояния?; Что теряется при случайной проекции? - Достаточно: k-NN, k-means, SVM (если ядро через расстояния), LSH поиск - все используют только попарные расстояния - Недостаточно: PCA (нужны направления, не расстояния), интерпретация признаков, регрессия с интерпретируемыми коэффициентами - Случайная проекция уничтожает семантику признаков - после JL нельзя сказать 'признак 3 важнее признака 7' - Для кластеризации в высоких размерностях JL часто применяют как препроцессинг --- **Почему leverage scores лучше равномерной субвыборки для sketch-and-solve?** Hints: Что происходит если в данных есть выбросы?; Что leverage score измеряет? - Leverage score строки i = ||U[i,:]||², где U - левые сингулярные векторы A - Строки с высоким leverage score 'определяют' пространство строк матрицы - без них задача меняется кардинально - Равномерная субвыборка с вероятностью 1-k/n пропустит такую строку - катастрофа для точности - Оптимальная субвыборка по leverage scores даёт (1+ε)-приближение с O(d log d / ε²) строками вместо m

Что унести из урока

• **JL лемма**: n точек из R^d → R^k, k = O(log n / ε²), все расстояния с точностью ε; k не зависит от d
• **Randomized SVD (HMT)**: O(m·n·k) вместо O(m·n²); 10-100x быстрее; погрешность ~ машинная точность для первых k векторов
• **sklearn.utils.extmath.randomized_svd** - стандартный инструмент для PCA, NMF, TruncatedSVD
• **Sketching**: умножение на S s×m → задача в s раз меньше; CountSketch/OSNAP работают за O(nnz · log n)
• **Random Fourier Features**: k(x,y) ≈ φ(x)ᵀφ(y); kernel SVM → линейная модель за O(n·D)
• **Leverage scores**: не равномерная выборка, а пропорциональная важности строк - оптимальная гарантия точности
• Случайность здесь - вычислительный ресурс: среднее по случайным проекциям сходится к точному ответу

Связи

Рандомизированные методы строятся поверх классической линейной алгебры и применяются везде где размерность велика

• SVD и сингулярные числа — Рандомизированный HMT аппроксимирует те же U, Σ, V - нужно понимать что именно ускоряется
• Тензорные разложения — Рандомизированное SVD - строительный блок быстрых Tucker и TT разложений
• Разреженная линейная алгебра — Sketching и итерационные методы решают похожие задачи для разреженных систем

Связанные уроки

• alg-33-randomized