Matrix Query + Block Krylov SVD

Идеи

• Быстрые алгоритмы для матриц расстояний в различных метриках
• Randomized Block Krylov Method для частичной задачи на сингулярные числа матрицы расстояний

Запросы к матрице расстояний

• Пусть имеется датасет $X\in {\mathbb{R}}^{n\times d}$.
• Матрица расстояний $A\in {\mathbb{R}}^{n\times n}$ содержит расстояния $A_{ij}=\rho (x_i,x_j)$, где $\rho (\cdot ,\cdot )$ - некоторая метрика.
• Запросом к матрице $A$ будем называть матрично-векторное умножение $Ay,y\in {\mathbb{R}}^n$.
• Запросы позволяют решать численные задачи линейной алгебры без хранения матрицы $A$

Быстрые запросы к матрице расстояний в ${\ell }_2^2$ метрике

k-ая координата для query $Ay$, где A - матрица расстояний, может быть найдена следующим образом:

$$(Ay{)}_k=\sum _{j=1}^n{y}_j\parallel x_k-x_j{\parallel }_2^2=||x_k|{|}_2^2\sum _{j=1}^n{y}_j+\sum _{j=1}^n{y}_j||x_j|{|}_2^2-2⟨x_k,\sum _{j=1}^n{y}_j{x}_j⟩$$

Составим простой query-алгоритм

1. $v=\sum _{i=1}^n{y}_i{x}_i$

2. $S_1=\sum _{i=1}^n{y}_i$

3. $S_2=\sum _{i=1}^n{y}_i||x_i|{|}_2^2$

4. $\text{ans}(k)=S_1||x_k|{|}_2^2+S_2-2⟨x_k,v⟩$

Сложность $O(nd)$

Сравнение наивного алгоритма запроса с быстрым

Основная особенность — отсутствие preprocessing

Быстрые запросы к матрице расстояний в ${\ell }_1$ метрике

$k$-ая координата для query $Ay$, где A - матрица расстояний, может быть найдена следующим образом:

$$(Ay{)}_k=\sum _{j=1}^n{y}_j\parallel x_k-x_j\parallel =\sum _{j=1}^n\sum _{i=1}^d{y}_j|x_{k,i}-x_{j,i}|=\sum _{i=1}^d\sum _{j=1}^n{y}_j|x_{k,i}-x_{j,i}|$$

Для каждого признака $i=\overline{1,d}$ введем перестановку ${\pi }^i$, которая соответствует отсортированному в порядке возрастания массиву значений $x_{j,i}$ (по столбцам). Тогда:

$$(Ay{)}_k=\sum _{i=1}^d\sum _{j=1}^n{y}_j|x_{k,i}-x_{j,i}|=\sum _{i=1}^d(\sum _{j:{\pi }^i(k)⩽{\pi }^i(j)}{y}_j(x_{j,i}-x_{k,i})+\sum _{j:{\pi }^i(k)>{\pi }^i(j)}{y}_j(x_{k,i}-x_{j,i}))$$

Перегруппируем значения в скобках:

$$(Ay{)}_k=\sum _{i=1}^d(x_{k,i}(\sum _{j:{\pi }^i(k)>{\pi }^i(j)}{y}_j-\sum _{j:{\pi }^i(k)⩽{\pi }^i(j)}{y}_j)+\sum _{j:{\pi }^i(k)⩽{\pi }^i(j)}{y}_j{x}_{j,i}-\sum _{j:{\pi }^i(k)>{\pi }^i(j)}{y}_j{x}_{j,i})$$ $$(Ay{)}_k=\sum _{i=1}^d(x_{k,i}(S_3-S_4)+S_2-S_1)$$

Сложность алгоритма

• Препроцессинг - поиск перестановок ${\pi }^i$

• Cложность препроцессинга $O(nd\log n)$

• Cложность query = $O(nd)$

Сравнение наивного алгоритма запроса с быстрым

Частичная задача на сингулярные числа

Randomized Block Krylov Method [MM15]:

Вход: $A\in {\mathbb{R}}^{n\times d}$, ошибка $ϵ\in (0,1)$, ранг $k\le n,d$

Выход: ${\hat{U}}_k,{\hat{\mathrm{\Sigma }}}_k,{\hat{V}}_k$

1. $q=\mathrm{\Theta }(\log (d)/\sqrt{ϵ})$
2. $\mathrm{\Pi }\sim \mathcal{N}(0,1{)}^{d\times k}$
3. $K:=[A^{T}A\mathrm{\Pi },(A^{T}A{)}^2\mathrm{\Pi },\dots ,(A^{T}A{)}^q\mathrm{\Pi }]\in {\mathbb{R}}^{n\times qk}$
4. $K=QR,Q\in {\mathbb{R}}^{n\times qk}$
5. $\hat{U},\hat{\mathrm{\Sigma }},\hat{V}\leftarrow \text{SVD}(AQ)$
6. Вернуть ${\hat{U}}_k,{\hat{\mathrm{\Sigma }}}_k,Q{\hat{V}}_k$

Подробности в bksvd.py.

Сходимость

Теорема (10, 11, 12 из [MM15]). С вероятностью $0.99$ алгоритм BKSVD возвращает $k$-ранговую аппроксимацию $\widehat{A_k}={\hat{U}}_k{\hat{\mathrm{\Sigma }}}_k{\hat{V}}_k$ такую, что

$$\begin{array}{rl}\parallel A-{\hat{A}}_k{\parallel }_2& \le (1+ϵ)\parallel A-A_k{\parallel }_2\\ \parallel A-{\hat{A}}_k{\parallel }_F& \le (1+ϵ)\parallel A-A_k{\parallel }_F\\ |{\sigma }_i-{\hat{\sigma }}_i|& \le ϵ\cdot {\sigma }_{k+1}^2\end{array}$$

Для этого требуется $q=\mathrm{\Theta }\left(k/\sqrt{ϵ}\right)$ запросов к матрице $A$.

Теорема (1.3 из [BCW22]). Для $ϵ>0$ и константы $p\ge 1$ найдётся распределение $\mathcal{D}$ матриц $n\times n$ таких, что для $A\sim \mathcal{D}$ любой алгоритм, который хотя с константной вероятностью возвращает вектор $v:$

$$\parallel A-Av{v}^T{\parallel }_{S_p}^p\le (1+ϵ)\underset{\parallel u{\parallel }_1=1}{min}\parallel A-Au{u}^T{\parallel }_{S_p}^p,$$

требует $\mathrm{\Omega }(1/{ϵ}^{1/3})$ запросов к матрице $A$.

Gaussian Mixture (24K объектов)

Время svds	120.7 сек	77.3 сек
Время постр. dist matrix	10.55 сек	17.46 сек
Время preproc	0.227 сек	-

MNIST (15K объектов)

Время svds	13.95 сек	11.68 сек
Время постр. dist matrix	525.2 сек	819.41 сек
Время preproc	3.1 сек	-

Препятствия для ускорения других методов

Просто заменить matvec на query не получится:

• Предобуславливатели

• Уравнение коррекции Якоби

Возможное решение: использовать и запросы, и матрицу расстояний.

Литература

• [IS22] Indyk, P., Silwal, S.. (2022). Faster Linear Algebra for Distance Matrices.
• [MM15] Musco, C., & Musco, C.. (2015). Randomized Block Krylov Methods for Stronger and Faster Approximate Singular Value Decomposition.
• [BCW22] Bakshi, A., Clarkson, K., & Woodruff, D.. (2022). Low-Rank Approximation with $1/{\epsilon }^{1/3}$ Matrix-Vector Products.