Upload
phamkhanh
View
213
Download
0
Embed Size (px)
Citation preview
Metodo do gradiente projetado espectral para minimizacao comrestricoes de ortogonalidade
Tiara MartiniInstituto de Matematica, Estatıstica e Computacao Cientıfica (IMECC), Unicamp
13083-859, Campinas, SP
E-mail: [email protected]
Juliano de Bem FranciscoDepartamento de Matematica, UFSC
88040-900, Florianopolis, SC
E-mail: [email protected]
Resumo: Apresentamos um metodo globalmente convergente e nao monotono para minimizacaocom restricoes de ortogonalidade. Desenvolvido recentemente por Francisco e Viloche Bazan [12]esse metodo esta baseado nas ideias dos metodos de regiao de confianca e Levenberg-Marquardt.Dessa maneira, os subproblemas consistem em minimizar um modelo quadratico da funcao obje-tivo sujeito ao conjunto de restricoes. Outros metodos viaveis sao citados, entre eles um metodode busca curvilinear e um de minimizacao ao longo de geodesicas. O desempenho dos metodosquando aplicados ao Problema de Procrustes Ortogonal e ilustrado numericamente.
Palavras-chave: Restricoes de ortogonalidade, algoritmo nao monotono, Metodo do GradienteProjetado Espectral, Problema de Procrustes Ortogonal, minimizacao em variedades de Stiefel
1 Introducao
O problema de minimizacao restrito ao conjunto das matrizes com colunas ortonormais aparecefrequentemente em importantes classes de problemas de otimizacao e, assim como as restricoeslinear e quadratica, possui uma estrutura especial a ser explorada. Como exemplo, citamosaplicacoes em otimizacao polinomial [8], otimizacao combinatorial [7], problemas de autovalorese analise de componentes principais (PCA) esparsos [2]. De acordo com Wen e Yin [15] e, emgeral, um problema do tipo NP-hard. Alem disso as restricoes de ortogonalidade podem geraralgumas dificuldades: conduzir a muitos mınimos locais; nao garantir a convergencia para omınimo global, exceto em casos mais simples; e manter a viabilidade da sequencia de pontos,visto que preservar as restricoes de ortogonalidade numericamente e muito caro, em geral, envolveuma decomposicao SVD.
Metodos usais para a resolucao dessa classe de problemas constituem em regiao de confianca[11], busca curvilinear [15], relaxacao [5], metodos de Newton e gradiente conjugado ao longode geodesicas [9]. O metodo estudado utiliza um parametro de regularizacao que controla otamanho do passo, a exemplo dos metodos de regiao de confianca.
2 Minimizacao com restricoes de ortogonalidade
Seja h : Rm×n → Rn×n, m > n, definida por h(X) = XTX − I, e considere Ω = X ∈ Rm×n :h(X) = 0.
Aplicaremos o metodo para minimizacao em conjuntos fechados desenvolvido em [12], pararesolver o seguinte problema de minimizacao com restricoes de ortogonalidade,
353
min f(X)s.a. X ∈ Ω,
(1)
em que f : Rm×n → R e uma funcao de classe C2.Note que,
g(X) =(
∂f∂X1
(X) · · · ∂f∂Xn
(X))∈ Rm×n,
em que Xi e a i-esima coluna de X. Alem disso, como a norma de Frobenius ||X||F =√n, para
todo X ∈ Ω, segue que Ω e compacto e g e Lipschitz contınua em Ω.O proximo teorema estabelece uma equivalencia entre as Condicoes KKT do problema (1)
e a projecao sobre o espaco tangente de Ω, o que simplifica consideravelmente a verificacao deponto estacionario.
Teorema 1 Seja X∗ ∈ Ω. Sao equivalentes:
(i) X∗(XT∗ g(X∗) + g(X∗)
TX∗)− 2g(X∗) = 0;
(ii) X∗ satisfaz as Condicoes KKT de (1).
2.1 Algoritmo
Considere as constantes 0 < σmin < σmax < ∞, Lu > Lf1, para cada k = 0, 1, 2..., escolha
σkspc ∈ [σmin, σmax], e defina
σkρ =
σkspc
2 , se 0 < ρ < Lu
Lu, caso contrario
, (2)
em que, ρ > 0 e um parametro de regularizacao.Assim, o modelo quadratico Qk
ρ : Rm×n → R utilizado sera dado por:
Qkρ(X) = traco(g(Xk)
T (X −Xk)) +σkρ + ρ
2||X −Xk||2F .
Para o parametro σkspc vamos escolher, sempre que possıvel, o parametro espectral de Barzilai-
Borwein [3]. Dessa forma, para dois iterados consecutivos, Xk−1 e Xk, temos
σkbb =
traco((g(Xk)T − g(Xk−1)
T )(Xk −Xk−1))
||Xk −Xk−1||2F,
e definimos
σkspc =
1, para k = 0minmaxσk
bb, σmin, σmax, para k > 1. (3)
Com estas escolhas, o algoritmo estudado e uma variacao do Metodo Gradiente ProjetadoEspectral nao monotono de [4], para minimizacao com restricoes de ortogonalidade, e a solucaodo subproblema
min Qkρ(X)
s.a. XTX = IX ∈ Rm×n.
fica dada por X∗ = UV T , em que UΣV T e a decomposicao SVD reduzida da matriz
Wk = Xk −1
ρ+ σkρ
g(Xk).
A nao monotonia do metodo esta baseada na tecnica de busca linear nao monotona propostaem [13] e apresentada na definicao a seguir.
1Lf e a constante de Lipschitz associada ao gradiente de f
354
Definicao 1 Sejam xkk∈N uma sequencia definida por xk+1 = xk + tkdk, dk = 0, tk ∈ (0, 1),γ ∈ (0, 1), M ∈ Z+ e m(0) = 0. Para cada k > 1, escolha m(k) recursivamente por
0 6 m(k) 6 minm(k − 1) + 1,M,
e definaf(xl(k)) = maxf(xk−j) : 0 6 j 6 m(k).
Dizemos que xk+1 satisfaz a Condicao de Armijo Nao Monotona se,
f(xk+1) 6 f(xl(k)) + γtk∇f(xk)Tdk.
Assim, podemos agora estabelecer o algoritmo para problemas com restricoes de ortogo-nalidade, que denominamos Gradiente Projetado para minimizacao em variedades de Stiefel(GPST).
Algoritmo 1: Algoritmo para problemas com restricoes de ortogonalidade (GPST)
Tome X0 ∈ Ω, M ∈ Z+, η ∈ (0, 1], β1 ∈ (0, 12 ], 0 < σmin < σmax < ∞ e 1 < ξ1 6 ξ2 < ∞.
Faca k = 0 e m(0) = 0.
Passo 1. Calcule g(Xk), f(Xl(k)) e faca ρ =σkspc
2 como em (3);
Passo 2. Tome σkρ definido em (2) e calcule a SVD reduzida de Wk: Wk = UΣV T ;
Passo 3. Defina Xkρ = UV T e encontre Xk
ρ ∈ Ω tal que Qkρ(X
kρ ) 6 ηQk
ρ(Xkρ). Se
Qkρ(X
kρ ) = 0, declare Xk ponto estacionario;
Passo 4. Defina
Ψk(X) = traco(g(Xk)T (X −Xk)) +
σkρ
2||X −Xk||2F .
Sef(Xk
ρ ) 6 f(Xl(k)) + β1Ψk(Xkρ ),
defina ρk = ρ, Xk+1 = Xkρ , faca k = k + 1 e volte ao Passo 1.
Senao, escolha ρnovo ∈ [ξ1ρ, ξ2ρ], faca ρ = ρnovo e volte ao Passo 2.
O teorema a seguir garante a convergencia global da sequencia de iterados Xk.
Teorema 2 Seja X um ponto de acumulacao da sequencia gerada pelo Algoritmo 1. Entao, Xe um ponto estacionario do problema (1) que nao e ponto de maximo local. Alem disso, se oconjunto dos pontos estacionarios do problema for finito, entao a sequencia converge.
3 Resultados
Implementamos o algoritmo apresentado e comparamos com metodos tradicionais para resol-ver problemas de minimizacao com restricoes de ortogonalidade, a saber um metodo de buscacurvilinear (OptStiefel) [15] e um metodo de minimizacao em geodesicas (sg min) [1].
Os testes foram executados utilizando Matlab 7.10 em um processador intel CORE i3 comHD de 500 Gb e 4 Gb de memoria Ram, sendo utilizadas os seguintes valores para as constantesdo Algoritmo 1: β1 = 10−4, ξ1 = ξ2 = 5, σmin = 10−10, Lu = Lf , σmax = Lu e M = 7.Assumiremos convergencia quando
||Xk(XTk g(Xk) + g(Xk)
TXk)− 2g(Xk)||F 6 10−6,
limitando o tempo de execucao (TMAX) em 600 segundos.
355
3.1 Testes numericos e problemas
Um problema do tipo Procrustes consiste em analisar o processo de preservacao de uma trans-formacao a um conjunto de formas. Tal processo envolve rotacao e mudanca de escala em umcerto conjunto de dados, a fim de move-los para um quadro comum de referencia onde possamser comparados. Os dados conhecidos sao dispostos em matrizes e a transformacao de rotacaoa ser determinada e uma matriz com colunas ortonormais. Em um problema de otimizacaoconsideramos a definicao a seguir.
Definicao 2 Seja X ∈ Rm×n, > n. Dadas as matrizes A ∈ Rp×m, C ∈ Rn×q e B ∈ Rp×q, oproblema denominado Problema de Procrustes Ortogonal Ponderado (WOPP) consisteem resolver
min ||AXC −B||2Fs.a. XTX = I.
(4)
Quando C = I diremos que (4) e o Problema de Procrustes Ortogonal (OPP) .
De acordo com Elden e Park [10], o problema OPP com m = n e dito balanceado e o caso emque n < m e chamado de nao balanceado. Para o caso balanceado, temos uma formula fechadapara a solucao, a saber X∗ = UV T , em que ATB = UΣV T e a decomposicao SVD reduzida deATB. Para o caso nao balanceado e para WOPP, a solucao e encontrada atraves de metodositerativos e, portanto, nao e necessariamente o mınimo global.
Para os testes numericos vamos considerar n = q, p = m, A = PSRT e C = QΛQT , com Pe R matrizes ortogonais com entradas randomicas, Q ∈ Rn×n matriz de Householder, Λ ∈ Rn×n
diagonal com entradas uniformemente distribuıdas no intervalo [12 , 2] e S diagonal definida demodo diferente para cada um dos problemas. Como ponto inicial X0 vamos gerar uma matrizrandomica satisfazendo X0 ∈ Ω. Todos os valores randomicos foram gerados com a funcao randndo Matlab.
Diante da solucao exata do problema, vamos criar uma ja conhecida solucao Q∗ tomandoB = AQ∗C, em que Q∗ ∈ Rm×n e matriz de Householder randomica, para monitorar o compor-tamento dos iterados Xk.
Os problemas testados foram retirados de [16] e sao descritos abaixo.
Problema 1. Os elementos da diagonal de S sao gerados atraves de uma distribuicao normalno intervalo [10, 12].Problema 2. A diagonal de S e dada por: Sii = i + 2ri, em que ri sao numeros aleatoriosdistribuıdos uniformemente no intervalo [0, 1].
Problema 3. Cada entrada da diagonal de S e gerada da seguinte maneira: Sii = 1+ 99(i−1)m+1 +
2ri, com ri distribuıdos uniformemente no intervalo [0, 1].A Tabela 1 expoe os resultados alcancados para os Problemas 1 e 3. Denotamos o numero
de iteracoes por Ni e o numero de avaliacoees da funcao por Na. O resıduo ||AXkC − B||2F , anorma do erro na solucao, ||Xk − Q∗||2F , e o tempo (em segundos) tambem sao apresentados.Vale ressaltar que o uso da nao monotonia foi eficaz, conforme observado na Figura 1.
Buscando reduzir o esforco computacional exigido no Algoritmo 1 durante o calculo da SVD,propomos encontrar essa decomposicao por meio de um algorimo alternativo [6]. Como resultadoconseguimos uma reducao media de 21% no tempo computacional.
4 Consideracoes finais
Apresentamos um metodo para resolver problemas com restricoes de ortogonalidade com con-vergencia global para pontos estacionarios. Com intuito de verificar os resultados teoricos eanalisar o desempenho do algoritmo apresentado, aplicamos o mesmo no problema WOPP.Apesar da programacao ter sido realizada sem rigor numerico, os resultados alcancados foram
356
Problema 1 com m = 500 e q = 70
Metodo Ni Na Tempo Resıduo Erro
GPST 20 25 2.230 3.52e-08 2.59e-09OptStiefel 31 31 1.716 3.67e-08 2.96e-09sg min 33 134 3.631 1.88e-08 1.41e-09
Problema 3 com m = 500 e q = 20
GPST 1080 1505 10.85 1.42e-07 5.19e-08OptStiefel 662 715 4.882 4.23e-07 1.53e-07sg min - - TMAX - -
Tabela 1: Desempenho dos metodos
0 100 200 300 400 50010
−8
10−6
10−4
10−2
100
102
Res
íduo
Número de iterações0 100 200 300 400 500
10−10
10−8
10−6
10−4
10−2
100
Nor
ma
do e
rro
na s
oluç
ão
Número de iterações
Figura 1: Comportamento do Algoritmo 1 para o Problema 2 com m = 50 e q = 10
promissores quando comparados com os obtidos via metodos usuais. Logo, o algoritmo estudadopode ser considerado competitivo diante dos metodos ja consolidados.
Ao algoritmo, incorporamos a ideia do calculo da SVD de uma maneira alternativa. O mesmoapresentou uma melhora significativa no desempenho quando utilizado para problemas em queo numero de linhas e muito maior do que o numero de colunas.
Adicionalmente, ressaltamos que a abordagem apresentada tambem pode ser utilizada pararesolver o Problema de Procrustes Simetrico [14].
Referencias
[1] T. A. Arias, A. Edelman, S. T. Smith. “The geometry of algorithms with ortogonalityconstraints”. SIAM Journal on Matrix Analysis and Applications, v. 20, n. 2, pp. 303-353,1998.
[2] A. Aspremont, L. E. Ghaoui, M. Jordan and G. R. Lanckriet. “A direct formulation forsparse pca using semidefinite programming”. SIAM Review, n. 49, pp. 434-448, 2007.
[3] J. Barzilai, J. M. Borwein. “Two point step size gradient methods”. IMA Journal of Nu-merical Analysis, v. 8, pp. 141-148, 1988.
[4] E. G. Birgin, J. M. Martınez, M. Raydan. “Nonmonotone spectral projected gradientmethods on convex sets”. SIAM Journal on Optimization, v. 10, pp. 1196-1211, 2000.
[5] A. W. Bojanczyk and A. Lutoborski. “The Procrustes Problem for Orthogonal Stiefel Ma-trices”. SIAM Journal on Scientific Computing, v.21, pp. 1291-1304, 1998.
357
[6] Y. Chahlaoui, K. Gallivan, P. Van Dooren. “Recursive calculation of dominant singularsubspaces”. SIAM Journal on Matrix Analysis and Applications, v. 25, n. 2, pp. 445-463,1999.
[7] W. J. Cook, W. H. Cunningham, W. R. Pulleyblank, A. Schrijver. “Combinatorial optimi-zation”. John Wiley & Sons, Canada, 1997.
[8] E. de Klerk, M. Laurent, P. A. Parrilo. “A PTAS for the minimization of polynomials offixed degree over the simplex”. Theoretical Computer Science, pp. 210-225, 2006.
[9] A. Edelman, R. Lippert. “Nonlinear Eigenvalue Problems with Orthogonality Constraints”.In Templates for the Solution of Algebraic Eigenvalue Problems, A Practical Guide, pp. 290-314, Philadelphia, 2000.
[10] L. Elden, H. Park. “A procrustes problem on Stiefel manifolds”. Numerishe Mathematik,v. 82, pp. 599-619, 1999.
[11] J. B. Francisco and J. M. Martınez and L. Martınez. “Globally convergent trust-regionmethods for self-consistent fiel electronic structure calculations”. Journal of Chemical Phy-sics, v.121, pp. 10863-10878, 2004.
[12] J. B. Francisco and F. S. Viloche Bazan. “Nonmonotone algorithm for minimization onclosed sets with application to minimization on Stiefel manifolds”. Journal of Computationaland Applied Mathematics, v. 236, n. 10, pp. 2717-2727, 2012.
[13] L. Grippo, F. Lampariello, S. Lucidi. “A nonmonotone line search technique for Newton’smethod”. SIAM Journal on Numerical Analysis, v. 23, n. 4, pp. 707-716, 1986.
[14] Nicholas J. Higham. “The symmetric procrustes problem”. BIT Numerical Mathematics,v. 28, n. 1, pp. 133-143, 1988.
[15] Z. Wen, W. Yin. “A feasible method for optimization with orthogonality constraints”.Optimization Online, pp. 1-30, 2010.
[16] Z. Zhang, K. Du. “Successive projection method for solving the unbalanced procrustesproblem”. Science in China: Series A Mathematics, v. 49, n. 7, pp. 971-986, 2006.
358