Низкоранговое приближение

В математике приближение низкого ранга относится к процессу приближения заданной матрицы матрицей более низкого ранга. Точнее, это задача минимизации , в которой функция стоимости измеряет соответствие между заданной матрицей (данными) и приближающей матрицей (переменной оптимизации) при условии ограничения, что приближающая матрица имеет пониженный ранг . Задача используется для математического моделирования и сжатия данных . Ограничение ранга связано с ограничением на сложность модели, которая подходит к данным. В приложениях часто существуют другие ограничения на приближающую матрицу помимо ограничения ранга, например, неотрицательность и структура Ганкеля .

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

Данный

• спецификация структуры ,${\displaystyle {\mathcal {S}}:\mathbb {R} ^{n_{p}}\to \mathbb {R} ^{m\times n}}$
• вектор параметров структуры ,${\displaystyle p\in \mathbb {R} ^{n_{p}}}$
• норма , и${\displaystyle \|\cdot \|}$
• желаемый ранг ,${\displaystyle r}$

${\displaystyle {\text{minimize}}\quad {\text{over}}{\widehat {p}}\quad \|p-{\widehat {p}}\|\quad {\text{subject to}}\quad \operatorname {rank} {\big (}{\mathcal {S}}({\widehat {p}}){\big )}\leq r.}$

• Идентификация линейной системы , в этом случае аппроксимирующая матрица имеет структуру Ганкеля .
• Машинное обучение , в этом случае аппроксимирующая матрица имеет нелинейную структуру.
• Рекомендательные системы , в которых матрица данных имеет пропущенные значения , а аппроксимация является категориальной .
• Завершение матрицы расстояний , в этом случае существует ограничение положительной определенности.
• Обработка естественного языка , в этом случае приближение неотрицательно .
• Компьютерная алгебра , в этом случае приближение структурировано по Сильвестру .

Неструктурированная задача с соответствием, измеряемым нормой Фробениуса , т.е.

${\displaystyle {\text{minimize}}\quad {\text{over}}{\widehat {D}}\quad \|D-{\widehat {D}}\|_{\text{F}}\quad {\text{subject to}}\quad \operatorname {rank} {\big (}{\widehat {D}}{\big )}\leq r}$

имеет аналитическое решение в терминах сингулярного разложения матрицы данных. Результат называется леммой аппроксимации матриц или теоремой Эккарта–Янга–Мирского . Эта проблема была первоначально решена Эрхардом Шмидтом ^[1] в бесконечномерном контексте интегральных операторов (хотя его методы легко обобщаются на произвольные компактные операторы в гильбертовых пространствах) и позднее переоткрыта К. Эккартом и Г. Янгом . ^[2] Л. Мирский обобщил результат на произвольные унитарно инвариантные нормы. ^[3] Пусть

${\displaystyle D=U\Sigma V^{\top }\in \mathbb {R} ^{m\times n},\quad m\geq n}$

будет разложением сингулярных значений , где — прямоугольная диагональная матрица с сингулярными значениями . Для заданного , разбиения , , и следующим образом:${\displaystyle D}$${\displaystyle \Sigma =:\operatorname {diag} (\sigma _{1},\ldots ,\sigma _{m})}$${\displaystyle m\times m}$${\displaystyle \sigma _{1}\geq \ldots \geq \sigma _{m}}$${\displaystyle r\in \{1,\dots ,m-1\}}$${\displaystyle U}$${\displaystyle \Сигма}$${\displaystyle V}$

${\displaystyle U=:{\begin{bmatrix}U_{1}&U_{2}\end{bmatrix}},\quad \Sigma =:{\begin{bmatrix}\Sigma _{1}&0\\0&\Sigma _{2}\end{bmatrix}},\quad {\text{и}}\quad V=:{\begin{bmatrix}V_{1}&V_{2}\end{bmatrix}},}$

где есть , есть и есть . Тогда ранговая матрица, полученная из усеченного сингулярного разложения${\displaystyle U_{1}}$${\displaystyle m\times r}$${\displaystyle \Sigma _{1}}$${\displaystyle r\times r}$${\displaystyle V_{1}}$${\displaystyle r\times n}$${\displaystyle r}$

${\displaystyle {\widehat {D}}^{*}=U_{1}\Sigma _{1}V_{1}^{\top },}$

таково, что

${\displaystyle \|D-{\widehat {D}}^{*}\|_{\text{F}}=\min _{\operatorname {rank} ({\widehat {D}})\leq r}\|D-{\widehat {D}}\|_{\text{F}}={\sqrt {\sigma _{r+1}^{2}+\cdots +\sigma _{m}^{2}}}.}$

Минимизатор уникален тогда и только тогда, когда .${\displaystyle {\widehat {D}}^{*}}$${\displaystyle \sigma _{r+1}\neq \sigma _{r}}$

Пусть будет вещественной (возможно, прямоугольной) матрицей с . Предположим, что ${\displaystyle A\in \mathbb {R} ^{m\times n}}$${\displaystyle m\leq n}$

${\displaystyle A=U\Sigma V^{\top }}$

— разложение по сингулярным значениям . Напомним, что и — ортогональные матрицы, а — диагональная матрица с элементами, такими что .${\displaystyle A}$${\displaystyle U}$${\displaystyle V}$${\displaystyle \Sigma }$${\displaystyle m\times n}$ ${\displaystyle (\sigma _{1},\sigma _{2},\cdots ,\sigma _{m})}$${\displaystyle \sigma _{1}\geq \sigma _{2}\geq \cdots \geq \sigma _{m}\geq 0}$

Мы утверждаем, что наилучшее ранговое приближение к в спектральной норме, обозначенной как , дается выражением${\displaystyle k}$${\displaystyle A}$${\displaystyle \|\cdot \|_{2}}$

${\displaystyle A_{k}:=\sum _{i=1}^{k}\sigma _{i}u_{i}v_{i}^{\top }}$

где и обозначают -й столбец и , соответственно.${\displaystyle u_{i}}$${\displaystyle v_{i}}$${\displaystyle i}$${\displaystyle U}$${\displaystyle V}$

Во-первых, обратите внимание, что у нас есть

${\displaystyle \|A-A_{k}\|_{2}=\left\|\sum _{i=1}^{\color {red}{n}}\sigma _{i}u_{i}v_{i}^{\top }-\sum _{i=1}^{\color {red}{k}}\sigma _{i}u_{i}v_{i}^{\top }\right\|_{2}=\left\|\sum _{i=\color {red}{k+1}}^{n}\sigma _{i}u_{i}v_{i}^{\top }\right\|_{2}=\sigma _{k+1}}$

Поэтому нам нужно показать, что если , где и имеют столбцы, то .${\displaystyle B_{k}=XY^{\top }}$${\displaystyle X}$${\displaystyle Y}$${\displaystyle k}$${\displaystyle \|A-A_{k}\|_{2}=\sigma _{k+1}\leq \|A-B_{k}\|_{2}}$

Так как имеет столбцы, то должна быть нетривиальная линейная комбинация первых столбцов , т.е.${\displaystyle Y}$${\displaystyle k}$${\displaystyle k+1}$${\displaystyle V}$

${\displaystyle w=\gamma _{1}v_{1}+\cdots +\gamma _{k+1}v_{k+1},}$

такой, что . Без потери общности, мы можем масштабировать так, что или (эквивалентно) . Следовательно,${\displaystyle Y^{\top }w=0}$${\displaystyle w}$${\displaystyle \|w\|_{2}=1}$${\displaystyle \gamma _{1}^{2}+\cdots +\gamma _{k+1}^{2}=1}$

${\displaystyle \|A-B_{k}\|_{2}^{2}\geq \|(A-B_{k})w\|_{2}^{2}=\|Aw\|_{2}^{2}=\gamma _{1}^{2}\sigma _{1}^{2}+\cdots +\gamma _{k+1}^{2}\sigma _{k+1}^{2}\geq \sigma _{k+1}^{2}.}$

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

Пусть будет вещественной (возможно, прямоугольной) матрицей с . Предположим, что${\displaystyle A\in \mathbb {R} ^{m\times n}}$${\displaystyle m\leq n}$

${\displaystyle A=U\Sigma V^{\top }}$

является разложением по сингулярным значениям .${\displaystyle A}$

Мы утверждаем, что наилучшее ранговое приближение к в норме Фробениуса, обозначаемое как , определяется выражением${\displaystyle k}$${\displaystyle A}$${\displaystyle \|\cdot \|_{F}}$

${\displaystyle A_{k}=\sum _{i=1}^{k}\sigma _{i}u_{i}v_{i}^{\top }}$

где и обозначают -й столбец и , соответственно.${\displaystyle u_{i}}$${\displaystyle v_{i}}$${\displaystyle i}$${\displaystyle U}$${\displaystyle V}$

Во-первых, обратите внимание, что у нас есть

${\displaystyle \|A-A_{k}\|_{F}^{2}=\left\|\sum _{i=k+1}^{n}\sigma _{i}u_{i}v_{i}^{\top }\right\|_{F}^{2}=\sum _{i=k+1}^{n}\sigma _{i}^{2}}$

Поэтому нам нужно показать, что если , где и имеют столбцы, то${\displaystyle B_{k}=XY^{\top }}$${\displaystyle X}$${\displaystyle Y}$${\displaystyle k}$

${\displaystyle \|A-A_{k}\|_{F}^{2}=\sum _{i=k+1}^{n}\sigma _{i}^{2}\leq \|A-B_{k}\|_{F}^{2}.}$

По неравенству треугольника со спектральной нормой, если то . Предположим и соответственно обозначают ранговое приближение к и методом SVD, описанным выше. Тогда для любого${\displaystyle A=A'+A''}$${\displaystyle \sigma _{1}(A)\leq \sigma _{1}(A')+\sigma _{1}(A'')}$${\displaystyle A'_{k}}$${\displaystyle A''_{k}}$${\displaystyle k}$${\displaystyle A'}$${\displaystyle A''}$${\displaystyle i,j\geq 1}$

${\displaystyle {\begin{aligned}\sigma _{i}(A')+\sigma _{j}(A'')&=\sigma _{1}(A'-A'_{i-1})+\sigma _{1}(A''-A''_{j-1})\\&\geq \sigma _{1}(A-A'_{i-1}-A''_{j-1})\\&\geq \sigma _{1}(A-A_{i+j-2})\qquad ({\text{since }}{\rm {rank}}(A'_{i-1}+A''_{j-1})\leq i+j-2))\\&=\sigma _{i+j-1}(A).\end{aligned}}}$

Так как , когда и заключаем, что для${\displaystyle \sigma _{k+1}(B_{k})=0}$${\displaystyle A'=A-B_{k}}$${\displaystyle A''=B_{k}}$${\displaystyle i\geq 1,j=k+1}$

${\displaystyle \sigma _{i}(A-B_{k})\geq \sigma _{k+i}(A).}$

Поэтому,

${\displaystyle \|A-B_{k}\|_{F}^{2}=\sum _{i=1}^{n}\sigma _{i}(A-B_{k})^{2}\geq \sum _{i=k+1}^{n}\sigma _{i}(A)^{2}=\|A-A_{k}\|_{F}^{2},}$

по мере необходимости.

Норма Фробениуса равномерно взвешивает все элементы ошибки аппроксимации . Априорные знания о распределении ошибок могут быть учтены путем рассмотрения задачи взвешенной аппроксимации низкого ранга${\displaystyle D-{\widehat {D}}}$

${\displaystyle {\text{minimize}}\quad {\text{over }}{\widehat {D}}\quad \operatorname {vec} (D-{\widehat {D}})^{\top }W\operatorname {vec} (D-{\widehat {D}})\quad {\text{subject to}}\quad \operatorname {rank} ({\widehat {D}})\leq r,}$

где векторизует матрицу по столбцам, а — заданная положительно (полу)определенная весовая матрица.${\displaystyle {\text{vec}}(A)}$ ${\displaystyle A}$${\displaystyle W}$

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

В случае некоррелированных весов задача взвешенной аппроксимации низкого ранга также может быть сформулирована следующим образом: ^{[4] [5]} для неотрицательной матрицы и матрицы мы хотим минимизировать по матрицам, , ранга которых не превышает .${\displaystyle W}$${\displaystyle A}$${\displaystyle \sum _{i,j}(W_{i,j}(A_{i,j}-B_{i,j}))^{2}}$${\displaystyle B}$${\displaystyle r}$

Пусть . Для , самый быстрый алгоритм работает за время. ^{[6] [7]} Одна из важных использованных идей называется Oblivious Subspace Embedding (OSE), она впервые предложена Сарлосом. ^[8]${\displaystyle \|A\|_{p}=\left(\sum _{i,j}|A_{i,j}^{p}|\right)^{1/p}}$${\displaystyle p=2}$${\displaystyle nnz(A)+n\cdot poly(k/\epsilon )}$

Для известно, что эта норма L1 по входу более надежна, чем норма Фробениуса при наличии выбросов и указывается в моделях, где гауссовские предположения о шуме могут не применяться. Естественно стремиться минимизировать . ^[9] Для и существуют некоторые алгоритмы с доказуемыми гарантиями. ^{[10] [11]}${\displaystyle p=1}$${\displaystyle \|B-A\|_{1}}$${\displaystyle p=0}$${\displaystyle p\geq 1}$

Пусть и будут двумя точечными множествами в произвольном метрическом пространстве. Пусть представляют матрицу, где . Такие матрицы расстояний обычно вычисляются в программных пакетах и ​​имеют приложения к изучению многообразий изображений, распознаванию рукописного ввода и многомерной развертке. В попытке уменьшить размер их описания ^{[12] [13]} можно изучить низкоранговую аппроксимацию таких матриц.${\displaystyle P=\{p_{1},\ldots ,p_{m}\}}$${\displaystyle Q=\{q_{1},\ldots ,q_{n}\}}$${\displaystyle A}$${\displaystyle m\times n}$${\displaystyle A_{i,j}=dist(p_{i},q_{i})}$

Задачи аппроксимации низкого ранга в распределенной и потоковой среде рассматривались в ^{[14] .}

Используя эквивалентности

${\displaystyle \operatorname {rank} ({\widehat {D}})\leq r\quad \iff \quad {\text{there are }}P\in \mathbb {R} ^{m\times r}{\text{ and }}L\in \mathbb {R} ^{r\times n}{\text{ such that }}{\widehat {D}}=PL}$

и

${\displaystyle \operatorname {rank} ({\widehat {D}})\leq r\quad \iff \quad {\text{there is full row rank }}R\in \mathbb {R} ^{m-r\times m}{\text{ such that }}R{\widehat {D}}=0}$

задача взвешенной аппроксимации низкого ранга становится эквивалентной задачам оптимизации параметров

${\displaystyle {\text{minimize}}\quad {\text{over }}{\widehat {D}},P{\text{ and }}L\quad \operatorname {vec} ^{\top }(D-{\widehat {D}})W\operatorname {vec} (D-{\widehat {D}})\quad {\text{subject to}}\quad {\widehat {D}}=PL}$

и

${\displaystyle {\text{minimize}}\quad {\text{over }}{\widehat {D}}{\text{ and }}R\quad \operatorname {vec} ^{\top }(D-{\widehat {D}})W\operatorname {vec} (D-{\widehat {D}})\quad {\text{subject to}}\quad R{\widehat {D}}=0\quad {\text{and}}\quad RR^{\top }=I_{r},}$

где — единичная матрица размера .${\displaystyle I_{r}}$${\displaystyle r}$

Представление изображения ограничения ранга предполагает метод оптимизации параметров, в котором функция стоимости минимизируется поочередно по одной из переменных ( или ) при фиксированной другой. Хотя одновременная минимизация по обеим и является сложной двояковыпуклой задачей оптимизации, минимизация по одной из переменных является линейной задачей наименьших квадратов и может быть решена глобально и эффективно.${\displaystyle P}$${\displaystyle L}$${\displaystyle P}$${\displaystyle L}$

Полученный алгоритм оптимизации (называемый чередующимися проекциями) глобально сходится с линейной скоростью сходимости к локально оптимальному решению задачи взвешенной аппроксимации низкого ранга. Начальное значение для параметра (или ) должно быть указано. Итерация останавливается, когда выполняется определенное пользователем условие сходимости.${\displaystyle P}$${\displaystyle L}$

Реализация алгоритма чередующихся проекций для взвешенной аппроксимации низкого ранга в Matlab :

function [dh, f] = wlra_ap ( d, w, p, tol, maxiter ) [ m , n ] = size ( d ); r = size ( p , 2 ); f = inf ; for i = 2 : maxiter % минимизация по L bp = kron ( eye ( n ), p ); vl = ( bp ' * w * bp ) \ bp ' * w * d (:); l = reshape ( vl , r , n ); % минимизация по P bl = kron ( l ' , eye ( m )); vp = ( bl ' * w * bl ) \ bl ' * w * d (:); p = reshape ( vp , m , r ); % проверка условия выхода dh = p * l ; dd = d - dh ; f ( i ) = dd (:) ' * w * dd (:); если abs ( f ( i - 1 ) - f ( i )) < tol , break , end endfor                                                                                        

Алгоритм чередующихся проекций использует тот факт, что задача аппроксимации низкого ранга, параметризованная в форме изображения, является билинейной по переменным или . Билинейная природа задачи эффективно используется в альтернативном подходе, называемом переменными проекциями. ^[15]${\displaystyle P}$${\displaystyle L}$

Рассмотрим снова задачу взвешенной аппроксимации низкого ранга, параметризованную в форме изображения. Минимизация по переменной (линейная задача наименьших квадратов) приводит к замкнутому выражению ошибки аппроксимации как функции${\displaystyle L}$${\displaystyle P}$

${\displaystyle f(P)={\sqrt {\operatorname {vec} ^{\top }(D){\Big (}W-W(I_{n}\otimes P){\big (}(I_{n}\otimes P)^{\top }W(I_{n}\otimes P){\big )}^{-1}(I_{n}\otimes P)^{\top }W{\Big )}\operatorname {vec} (D)}}.}$

Исходная задача, таким образом, эквивалентна нелинейной задаче наименьших квадратов минимизации относительно . Для этой цели могут быть использованы стандартные методы оптимизации, например, алгоритм Левенберга-Марквардта .${\displaystyle f(P)}$${\displaystyle P}$

Реализация алгоритма переменных проекций для взвешенной аппроксимации низкого ранга в Matlab :

function [dh, f] = wlra_varpro ( d, w, p, tol, maxiter ) prob = optimset(); prob.solver = ' lsqnonlin ' ; prob.options = optimset ( ' MaxIter ' , maxiter , ' TolFun ' , tol ) ; prob.x0 = p ; prob.objection = @ ( p ) cost_fun ( p , d , w ) ; [ p , f ] = lsqnonlin ( prob ) ; [ f , vl ] = cost_fun ( p , d , w ); dh = p * reshape ( vl , size ( p , 2 ) , size ( d , 2 ) ) ;                                       функция [f, vl] = cost_fun ( p, d, w ) bp = kron ( eye ( size ( d , 2 )), p ); vl = ( bp ' * w * bp ) \ bp ' * w * d (:); f = d (:) ' * w * ( d (:) - bp * vl );                           

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

Обычно мы хотим, чтобы наше новое решение не только имело низкий ранг, но и удовлетворяло другим выпуклым ограничениям из-за требований приложения. Наша интересующая нас проблема будет выглядеть следующим образом:

${\displaystyle {\text{minimize}}\quad {\text{over }}{\widehat {p}}\quad \|p-{\widehat {p}}\|\quad {\text{subject to}}\quad \operatorname {rank} {\big (}{\mathcal {S}}({\widehat {p}}){\big )}\leq r{\text{ and }}g({\widehat {p}})\leq 0}$

Эта проблема имеет множество приложений в реальном мире, включая восстановление хорошего решения из неточной (полуопределенного программирования) релаксации. Если дополнительное ограничение линейно, например, мы требуем, чтобы все элементы были неотрицательными, то проблема называется структурированной аппроксимацией низкого ранга. ^[16] Более общая форма называется выпукло-ограниченной аппроксимацией низкого ранга.${\displaystyle g({\widehat {p}})\leq 0}$

Эта задача полезна при решении многих задач. Однако она сложна из-за сочетания выпуклых и невыпуклых (низкого ранга) ограничений. Были разработаны различные методы на основе различных реализаций . Однако метод переменного направления множителей (ADMM) может быть применен для решения невыпуклой задачи с выпуклой целевой функцией, ограничениями ранга и другими выпуклыми ограничениями ^[17] и, таким образом, подходит для решения нашей вышеуказанной задачи. Более того, в отличие от общих невыпуклых задач, ADMM гарантирует сходимость допустимого решения, пока его двойственная переменная сходится в итерациях.${\displaystyle g({\widehat {p}})\leq 0}$

• Аппроксимация матрицы CUR производится из строк и столбцов исходной матрицы.

• MT Chu, RE Funderlic, RJ Plemmons, Структурированная аппроксимация низкого ранга, Линейная алгебра и ее приложения, том 366, 1 июня 2003 г., страницы 157–172 doi :10.1016/S0024-3795(02)00505-0

• Пакет C++ для структурированной аппроксимации низкого ранга