Aproximação de posto baixo

Em matemática, a aproximação de posto baixo é um problema de minimização, no qual a função de custo mede o ajuste entre uma dada matriz (os dados) e uma matriz de aproximação (a variável de otimização), sujeita a uma restrição de que a matriz de aproximação tenha posto reduzido. O problema é usado para modelagem matemática e compressão de dados. A restrição de classificação está relacionada a uma restrição na complexidade de um modelo que se ajusta aos dados. Em aplicações, muitas vezes há outras restrições na matriz de aproximação além da restrição de posto, por exemplo, não negatividade e estrutura de Hankel.

A aproximação de posto baixo está intimamente relacionada com:

• análise de componentes principais,
• análise fatorial,
• mínimos quadrados totais,
• análise semântica latente,
• regressão ortogonal, e
• decomposição em modos dinâmicos

Dados

• uma especificação de estrutura ${\displaystyle {\mathcal {S}}:\mathbb {R} ^{n_{p}}\to \mathbb {R} ^{m\times n}}$ ,
• um vetor de parâmetros de estrutura ${\displaystyle p\in \mathbb {R} ^{n_{p}}}$ ,
• uma norma ${\displaystyle \|\cdot \|}$, e
• o posto desejado ${\displaystyle r}$ ,

${\displaystyle {\text{minimizar}}\quad {\text{sobre }}{\widehat {p}}\quad \|p-{\widehat {p}}\|\quad {\text{sujeito a}}\quad \operatorname {rank} {\big (}{\mathcal {S}}({\widehat {p}}){\big )}\leq r.}$

• Identificação de sistema linear, caso em que a matriz de aproximação é estruturada por Hankel.^[1][2]
• Aprendizado de máquina, caso em que a matriz de aproximação é estruturada de forma não linear.^[3]
• Sistemas de recomendação, em que a matriz de dados tem valores faltantes e a aproximação é categórica.
• Completamento de matriz de distâncias, caso em que há uma restrição de que as matrizes sejam definidas positivas.
• Processamento de linguagem natural, caso em que a aproximação é não negativa.
• Álgebra computacional, caso em que a aproximação é estruturada por Sylvester.

O problema não estruturado com ajuste medido pela norma de Frobenius, ou seja,

${\displaystyle {\text{minimizar}}\quad {\text{sobre }}{\widehat {D}}\quad \|D-{\widehat {D}}\|_{\text{F}}\quad {\text{sujeito a }}\quad \operatorname {rank} {\big (}{\widehat {D}}{\big )}\leq r,}$

tem solução analítica em termos da decomposição em valores singulares da matriz de dados. O resultado é referido como o lema de aproximação de matrizes ou teorema de Eckart–Young–Mirsky. Este problema foi originalmente resolvido por Erhard Schmidt^[4] no contexto de dimensão infinita de operadores integrais (embora seus métodos facilmente se generalizem para operadores compactos arbitrários em espaços de Hilbert) e posteriormente redescoberto por C. Eckart e G. Young.^[5] L. Mirsky generalizou o resultado para normas arbitrárias unitariamente invariantes.^[6] Sejam

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

a decomposição em valores singulares de ${\displaystyle D}$ e partição ${\displaystyle U}$, ${\displaystyle \Sigma =:\operatorname {diag} (\sigma _{1},\ldots ,\sigma _{m})}$, e ${\displaystyle V}$ como segue:

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

em que ${\displaystyle U_{1}}$ é ${\displaystyle m\times r}$, ${\displaystyle \Sigma _{1}}$ é ${\displaystyle r\times r}$, e ${\displaystyle V_{1}}$ é ${\displaystyle n\times r}$ . Então a matriz de posto-${\displaystyle r}$, obtida a partir da decomposição em valores singulares truncada

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

é tal que

${\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}}}.}$

O minimizador ${\displaystyle {\widehat {D}}^{*}}$ é único se, e somente se, ${\displaystyle \sigma _{r+1}\neq \sigma _{r}}$.

Seja ${\displaystyle A\in \mathbb {R} ^{m\times n}}$ uma matriz real (possivelmente retangular) com ${\displaystyle m\geq n}$. Suponha que

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

é a decomposição em valores singulares de ${\displaystyle A}$. Lembre-se que ${\displaystyle U}$ e ${\displaystyle V}$ são matrizes ortogonais, e ${\displaystyle \Sigma }$ é uma matriz diagonal ${\displaystyle m\times n}$ com entradas ${\displaystyle (\sigma _{1},\sigma _{2},\cdots ,\sigma _{n})}$ tais que ${\displaystyle \sigma _{1}\geq \sigma _{2}\geq \cdots \geq \sigma _{n}\geq 0}$.

Afirmamos que a melhor aproximação de posto ${\displaystyle k}$ de ${\displaystyle A}$ na norma espectral, denotada por ${\displaystyle \|\cdot \|_{2}}$, é dada por

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

em que ${\displaystyle u_{i}}$ e ${\displaystyle v_{i}}$ denotam as ${\displaystyle i}$-ésimas colunas de ${\displaystyle U}$ e ${\displaystyle V}$, respectivamente.

Primeiro, note que temos

${\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}}$

Portanto, precisamos mostrar que se ${\displaystyle B_{k}=XY^{\top }}$, em que ${\displaystyle X}$ e ${\displaystyle Y}$ têm ${\displaystyle k}$ colunas, então ${\displaystyle \|A-A_{k}\|_{2}=\sigma _{k+1}\leq \|A-B_{k}\|_{2}}$.

Como ${\displaystyle Y}$ tem ${\displaystyle k}$ colunas, então deve haver uma combinação linear não trivial das primeiras ${\displaystyle k+1}$ colunas de ${\displaystyle V}$, ou seja,

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

tal que ${\displaystyle Y^{\top }w=0}$. Sem perda de generalidade, podemos escalar ${\displaystyle w}$ de modo que ${\displaystyle \|w\|_{2}=1}$ ou (equivalentemente) ${\displaystyle \gamma _{1}^{2}+\cdots +\gamma _{k+1}^{2}=1}$. Portanto,

${\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}.}$

O resultado segue tomando a raiz quadrada de ambos os lados da desigualdade acima.

Seja ${\displaystyle A\in \mathbb {R} ^{m\times n}}$ uma matriz real (possivelmente retangular) com ${\displaystyle m\geq n}$. Suponha que

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

é a decomposição em valores singulares de ${\displaystyle A}$.

Afirmamos que a melhor aproximação de posto ${\displaystyle k}$ de ${\displaystyle A}$ na norma de Frobenius, denotada por ${\displaystyle \|\cdot \|_{F}}$, é dada por

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

em que ${\displaystyle u_{i}}$ e ${\displaystyle v_{i}}$ denotam as ${\displaystyle i}$-ésimas colunas de ${\displaystyle U}$ e ${\displaystyle V}$, respectivamente.

Primeiro, note que temos

${\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}}$

Portanto, precisamos mostrar que se ${\displaystyle B_{k}=XY^{\top }}$, com ${\displaystyle X}$ e ${\displaystyle Y}$ tendo ${\displaystyle k}$ colunas, então

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

Pela desigualdade triangular com a norma espectral, se ${\displaystyle A=A'+A''}$ então ${\displaystyle \sigma _{1}(A)\leq \sigma _{1}(A')+\sigma _{1}(A'')}$. Suponha que ${\displaystyle A'_{k}}$ e ${\displaystyle A''_{k}}$ denotam respectivamente as aproximações de posto ${\displaystyle k}$ de ${\displaystyle A'}$ e ${\displaystyle A''}$ pelo método SVD descrito acima. Então, para qualquer ${\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{como }}{\rm {rank}}(A'_{i-1}+A''_{j-1})\leq {\rm {rank\,}}(A_{i+j-2}))\\&=\sigma _{i+j-1}(A).\end{aligned}}}$

Como ${\displaystyle \sigma _{k+1}(B_{k})=0}$, quando ${\displaystyle A'=A-B_{k}}$ e ${\displaystyle A''=B_{k}}$ concluímos que para ${\displaystyle i\geq 1,j=k+1}$

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

Portanto,

${\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},}$

como desejado.

A norma de Frobenius pondera uniformemente todos os elementos do erro de aproximação ${\displaystyle D-{\widehat {D}}}$. O conhecimento prévio sobre a distribuição dos erros pode ser levado em consideração considerando o problema de aproximação de posto baixo ponderada

${\displaystyle {\text{minimizar}}\quad {\text{sobre }}{\widehat {D}}\quad \operatorname {vec} (D-{\widehat {D}})^{\top }W\operatorname {vec} (D-{\widehat {D}})\quad {\text{sujeito a}}\quad \operatorname {rank} ({\widehat {D}})\leq r,}$

em que ${\displaystyle {\text{vec}}(A)}$ vetoriza a matriz ${\displaystyle A}$ por colunas e ${\displaystyle W}$ é uma matriz de peso positiva (semi-)definida dada.

O problema geral de aproximação de posto baixo ponderada não admite uma solução analítica em termos de decomposição de valores singulares e é resolvido por métodos de otimização local, que não garantem que uma solução global ótima seja encontrada.

No caso de pesos não correlacionados, o problema de aproximação de posto baixo ponderada também pode ser formulado desta forma:^[7][8] para uma matriz não negativa ${\displaystyle W}$ e uma matriz ${\displaystyle A}$ queremos minimizar ${\displaystyle \sum _{i,j}(W_{i,j}(A_{i,j}-B_{i,j}))^{2}}$ sobre matrizes ${\displaystyle B}$, de posto no máximo ${\displaystyle r}$.

Seja ${\displaystyle \|A\|_{p}=\left(\sum _{i,j}|A_{i,j}^{p}|\right)^{1/p}}$ . Para ${\displaystyle p=2}$, o algoritmo mais rápido é executado em tempo ${\displaystyle nnz(A)+n\cdot poly(k/\epsilon )}$.^[9][10] Uma das ideias importantes usadas é chamada de Oblivious Subspace Embedding (OSE), proposta pela primeira vez por Sarlos.^[11]

Para ${\displaystyle p=1}$, sabe-se que esta norma L1 por entradas é mais robusta do que a norma de Frobenius na presença de outliers e é indicada em modelos para os quais as suposições gaussianas sobre o ruído podem não se aplicar. É natural procurar minimizar ${\displaystyle \|B-A\|_{1}}$.^[12] Para ${\displaystyle p=0}$ e ${\displaystyle p\geq 1}$, existem alguns algoritmos com garantias prováveis.^[13][14]

Sejam ${\displaystyle P=\{p_{1},\ldots ,p_{m}\}}$ e ${\displaystyle Q=\{q_{1},\ldots ,q_{n}\}}$ dois conjuntos de pontos em um espaço métrico arbitrário. Seja ${\displaystyle A}$ uma matriz ${\displaystyle m\times n}$ em que ${\displaystyle A_{i,j}=dist(p_{i},q_{i})}$. Tais matrizes de distâncias são comumente calculadas em pacotes de software e têm aplicações para aprendizado de variedades de imagens, reconhecimento de escrita manual e desdobramento multidimensional. Na tentativa de reduzir seu tamanho de descrição,^[15][16] pode-se estudar uma aproximação de posto baixo de tais matrizes.

Os problemas de aproximação de posto baixo nos modelos distribuídos e de streaming foram considerados por Boutsidis et al.^[17]

Usando as equivalências

${\displaystyle \operatorname {rank} ({\widehat {D}})\leq r\quad \iff \quad {\text{existem }}P\in \mathbb {R} ^{m\times r}{\text{ e }}L\in \mathbb {R} ^{r\times n}{\text{ tais que }}{\widehat {D}}=PL}$

e

${\displaystyle \operatorname {rank} ({\widehat {D}})\leq r\quad \iff \quad {\text{existe }}R\in \mathbb {R} ^{m-r\times m}{\text{ de posto completo tal que }}R{\widehat {D}}=0}$

o problema de aproximação de posto baixo ponderada torna-se equivalente aos problemas de otimização de parâmetros

${\displaystyle {\text{minimizar}}\quad {\text{sobre }}{\widehat {D}},P{\text{ e }}L\quad \operatorname {vec} ^{\top }(D-{\widehat {D}})W\operatorname {vec} (D-{\widehat {D}})\quad {\text{sujeito a}}\quad {\widehat {D}}=PL}$

e

${\displaystyle {\text{minimizar}}\quad {\text{sobre }}{\widehat {D}}{\text{ e }}R\quad \operatorname {vec} ^{\top }(D-{\widehat {D}})W\operatorname {vec} (D-{\widehat {D}})\quad {\text{sujeito a }}\quad R{\widehat {D}}=0\quad {\text{e}}\quad RR^{\top }=I_{r},}$

em que ${\displaystyle I_{r}}$ é a matriz identidade de tamanho ${\displaystyle r}$.

A representação por imagem da restrição de posto sugere um método de otimização de parâmetros no qual a função de custo é minimizada alternativamente sobre uma das variáveis (${\displaystyle P}$ ou ${\displaystyle L}$) com a outra fixa. Embora a minimização simultânea de ${\displaystyle P}$ e ${\displaystyle L}$ seja um problema de otimização biconvexo difícil, a minimização sobre uma das variáveis sozinha é um problema linear de mínimos quadrados e pode ser resolvida globalmente e eficientemente.

O algoritmo de otimização resultante (chamado de projeções alternadas) é globalmente convergente com uma taxa de convergência linear para uma solução localmente ótima do problema de aproximação de posto baixo ponderada. O valor inicial para ${\displaystyle P}$ (ou ${\displaystyle L}$) deve ser fornecido. A iteração é interrompida quando uma condição de convergência definida pelo usuário é satisfeita.

O algoritmo de projeções alternadas para aproximação de posto baixo ponderada pode ser implementado em Matlab da seguinte forma:

function [dh, f] = wlra_ap(d, w, p, tol, maxiter)
[m, n] = size(d); r = size(p, 2); f = inf;
for i = 2:maxiter
    % minimização sobre L
    bp = kron(eye(n), p);
    vl = (bp' * w * bp) \ bp' * w * d(:);
    l  = reshape(vl, r, n);
    % minimização sobre P
    bl = kron(l', eye(m));
    vp = (bl' * w * bl) \ bl' * w * d(:);
    p  = reshape(vp, m, r);
    % verificação da condição de parada
    dh = p * l; dd = d - dh;
    f(i) = dd(:)' * w * dd(:);
    if abs(f(i - 1) - f(i)) < tol, break, end
endfor

O algoritmo de projeções alternadas explora o fato de que o problema de aproximação de posto baixo, parametrizado na forma da imagem, é bilinear nas variáveis ${\displaystyle P}$ ou ${\displaystyle L}$. A natureza bilinear do problema é efetivamente utilizada em uma abordagem alternativa, chamada de projeções variáveis.^[18]

Considere novamente o problema de aproximação de posto baixo ponderada, parametrizado na forma da imagem. A minimização em relação à variável ${\displaystyle L}$ (um problema linear de mínimos quadrados) leva à expressão de forma fechada do erro de aproximação em função de ${\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)}}.}$

O problema original é, portanto, equivalente ao problema não linear de mínimos quadrados de minimizar ${\displaystyle f(P)}$ em relação a ${\displaystyle P}$. Para este propósito, métodos de otimização padrão, por exemplo, o algoritmo de Levenberg-Marquardt podem ser usados.

Implementação Matlab do algoritmo de projeções variáveis para aproximação ponderada de baixa classificação:

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.objective = @(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));

function [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);

A abordagem de projeções variáveis também pode ser aplicada a problemas de aproximação de posto baixo parametrizados na forma de kernel. O método é eficaz quando o número de variáveis eliminadas é muito maior do que o número de variáveis de otimização deixadas no estágio de minimização não linear por mínimos quadrados. Tais problemas ocorrem na identificação do sistema, parametrizado na forma de kernel, onde as variáveis eliminadas são a trajetória de aproximação e as variáveis restantes são os parâmetros do modelo. No contexto de sistemas lineares invariantes no tempo, a etapa de eliminação é equivalente à suavização de Kalman.

Normalmente, queremos que nossa nova solução não seja apenas de posto baixo, mas também satisfaça outras restrições convexas devido aos requisitos da aplicação. Nosso problema de interesse seria o seguinte,

${\displaystyle {\text{minimizar}}\quad {\text{sobre }}{\widehat {p}}\quad \|p-{\widehat {p}}\|\quad {\text{sujeito a}}\quad \operatorname {rank} {\big (}{\mathcal {S}}({\widehat {p}}){\big )}\leq r{\text{ e }}g({\widehat {p}})\leq 0}$

Este problema tem muitas aplicações no mundo real, inclusive para recuperar uma boa solução de um relaxamento inexato (programação semidefinida). Se restrição adicional ${\displaystyle g({\widehat {p}})\leq 0}$ é linear, tal como quando se exige que todos os elementos sejam não-negativos, o problema é chamado de aproximação estruturada de posto baixo.^[19] A forma mais geral é chamada de aproximação de posto baixo restrita convexa.

Este problema é útil para resolver muitos problemas. No entanto, é um desafio devido à combinação das restrições convexas e não convexas (posto baixo). Diferentes técnicas foram desenvolvidas com base em diferentes realizações de ${\displaystyle g({\widehat {p}})\leq 0}$. No entanto, o Método de Multiplicadores de Direção Alternada (ADMM) pode ser aplicado para resolver o problema não convexo com função objetivo convexa, restrições de posto e outras restrições convexas,^[20] e, portanto, é adequado para resolver nosso problema acima. Além disso, diferentemente dos problemas gerais não convexos, o ADMM garantirá a convergência de uma solução viável desde que sua variável dual convirja nas iterações.

• A aproximação matricial CUR é feita a partir das linhas e colunas da matriz original

• MT Chu, RE Funderlic, RJ Plemmons, Aproximação estruturada de baixo escalão, Álgebra Linear e suas Aplicações, Volume 366, 1º de junho de 2003, Páginas 157–172doi:10.1016/S0024-3795(02)00505-0

• Pacote C++ para aproximação de posto baixo estruturada