Qual é a maneira mais rápida de calcular o maior autovalor de uma matriz geral?

27

EDIT: Estou testando se algum autovalor tem uma magnitude de um ou mais.

Preciso encontrar o maior autovalor absoluto de uma grande matriz esparsa e não simétrica.

Eu tenho usado a eigen()função de R , que usa o QR algo do EISPACK ou LAPACK para encontrar todos os autovalores e, em seguida, utilizo abs()para obter os valores absolutos. No entanto, eu preciso fazer isso mais rápido.

Eu também tentei usar a interface ARPACK no igraphpacote R. No entanto, deu um erro para uma das minhas matrizes.

A implementação final deve ser acessível a partir de R.

Provavelmente haverá vários autovalores da mesma magnitude.

Você tem alguma sugestão?

EDIT: Precisão só precisa ser 1e-11. Até agora, uma matriz "típica" foi . Consegui fazer uma fatoração de QR nele. No entanto, também é possível ter muito maiores. Atualmente, estou começando a ler sobre o algoritmo de Arnoldi. Eu entendo que isso está relacionado ao Lanczsos.386×386

EDIT2: Se eu tiver várias matrizes que estou "testando" e sei que existe uma submatriz grande que não varia. É possível ignorá-lo / descartá-lo?

poder
fonte
Veja minha resposta aqui: scicomp.stackexchange.com/a/1679/979 . Este é um tópico de pesquisa atual e os métodos atuais podem fazer melhor do que Lanczos. O problema de calcular valores singulares é equivalente ao problema de calcular valores próprios.
dranxo
2
Matriz de 400 x 400! = Grande. Além disso, o que significa maior se "Provavelmente haverá vários valores próprios da mesma magnitude". Em terras numpy: linalg.eig (random.normal (size = (400.400))) leva cerca de meio segundo. Isso é muito lento?
precisa saber é o seguinte
@meawoppl sim meio segundo é muito lento. Isso ocorre porque faz parte de outro algo que executa esse cálculo várias vezes.
poder
1
@power gotcah. Você tem uma aproximação ao vetor próprio. ou seja, é provavelmente semelhante à última solução, ou você pode adivinhar qual é a sua estrutura?
27513

Respostas:

14

Depende muito do tamanho da sua matriz, no caso em larga escala, também se é escassa e na precisão que você deseja obter.

Se sua matriz é muito grande para permitir uma única fatoração e você precisa de alta precisão, o algoritmo de Lanczsos é provavelmente o caminho mais rápido. No caso não simétrico, é necessário o algoritmo de Arnoldi, que é numericamente instável; portanto, uma implementação precisa abordar isso (é um pouco difícil de curar).

Se esse não for o seu problema, forneça informações mais específicas em sua pergunta. Em seguida, adicione um comentário a esta resposta e eu a atualizarei.

Edit: [Isso era para a versão antiga da pergunta, buscando o maior valor próprio.] Como sua matriz é pequena e aparentemente densa, eu faria a iteração de Arnoldi em B = (IA) ^ {- 1}, usando uma inicial permitiu que a fatoração triangular de IA tivesse multiplicação barata por B. (Ou calcule uma inversa explícita, mas isso custa três vezes mais que a fatoração.) Você deseja testar se B tem um valor próprio negativo. Trabalhando com B no lugar de A, os autovalores negativos são muito melhores separados; portanto, se houver um, você deverá convergir rapidamente.

Mas estou curioso para saber de onde vem o seu problema. Matrizes não simétricas geralmente têm autovalores complexos, portanto '' maior '' nem é bem definido. Portanto, você deve saber mais sobre o seu problema, o que pode ajudar a sugerir como resolvê-lo ainda mais rápido e / ou com mais confiabilidade.

Edit2: É difícil conseguir com Arnoldi um subconjunto particular de interesse. Para obter os valores autônomos absolutamente maiores, faça a iteração do subespaço usando a matriz original, com um tamanho de subespaço correspondendo ou excedendo o número de valores próprios que se espera que sejam próximos a 1 ou maior em magnitude. Em matrizes pequenas, isso será mais lento que o algoritmo QR, mas em matrizes grandes será muito mais rápido.

Arnold Neumaier
fonte
Preciso testar se o maior valor próprio é maior que 1. A precisão precisa ser apenas 1e-11. Até agora, uma matriz "típica" tem 386 x 386. Consegui fazer uma fatoração QR nela. No entanto, também é possível ter muito maiores. Atualmente, estou começando a ler sobre o algoritmo de Arnoldi. Eu entendo que isso está relacionado ao Lanczsos.
poder
Esta informação pertence à sua pergunta - por isso edite-a e adicione mais informações (por que os autovalores são reais? Ou o que significa maior?) - veja a edição da minha resposta.
Arnold Neumaier 21/03
desculpe por não me explicar claramente. Também não expliquei claramente que os autovalores são complexos. Estou testando se algum autovalor tem magnitude igual ou superior a um.
poder
1
Isso faz mais sentido, mas agora minha receita com funciona bem apenas se o autovalor ruim for realmente real> 1. Por outro lado, as novas informações provavelmente implicam que você tem pouca escolha além de calcular todos os autovalores. - Atualize sua pergunta para transmitir as informações extras! (IA)1
Arnold Neumaier 23/03
1
veja editar 2 na minha resposta
Arnold Neumaier 30/03
7

A Iteração de Potência (ou Método de Potência), por exemplo, o que Dan está descrevendo, deve sempre convergir, embora à taxa .|λn1/λn|

Se estiver próximo de λ n , será lento, mas você pode usar a extrapolação para contornar isso. Pode parecer complicado, mas uma implementação em pseudo-código é apresentada no artigo.λn1λn

Pedro
fonte
1
e se | λ (n − 1) | = | λ (n) | ?
poder
@power, a Iteração de energia normal não convergirá. Não sei até que ponto os métodos de extrapolação farão a distinção entre os diferentes autovalores. Você precisará ler o artigo para isso.
Pedro Pedro
2
|λn1|=|λn|λnλn1
você tem uma referência a um artigo ou livro acadêmico que apóie isso? Além disso, e se \ lambda_ {n} for complexo?
potência
5
Se houver vários autovalores diferentes do módulo máximo, a iteração de potência converge apenas em circunstâncias excepcionais. Geralmente oscila de maneira um tanto imprevisível.
Arnold Neumaier
5

Houve algumas boas pesquisas sobre isso recentemente. As novas abordagens usam "algoritmos aleatórios", que exigem apenas algumas leituras de sua matriz para obter boa precisão nos maiores valores próprios. Isso contrasta com as iterações de potência que requerem várias multiplicações de vetores matriciais para alcançar alta precisão.

Você pode ler mais sobre a nova pesquisa aqui:

http://math.berkeley.edu/~strain/273.F10/martinsson.tygert.rokhlin.randomized.decomposition.pdf

http://arxiv.org/abs/0909.4061

Este código fará isso por você:

http://cims.nyu.edu/~tygert/software.html

https://bitbucket.org/rcompton/pca_hgdp/raw/be45a1d9a7077b60219f7017af0130c7f43d7b52/pca.m

http://code.google.com/p/redsvd/

https://cwiki.apache.org/MAHOUT/stochastic-singular-value-decomposition.html

Se o seu idioma de escolha não estiver lá, você poderá rolar seu próprio SVD aleatório com bastante facilidade; requer apenas uma multiplicação de vetores de matriz seguida de uma chamada para um SVD pronto para uso.

dranxo
fonte
4

Aqui você encontrará uma introdução algorítmica ao algoritmo Jacobi-Davidson, que calcula o valor próprio máximo.

Neste artigo, os aspectos matemáticos são explorados. O JD permite matrizes gerais (reais ou complexas) e pode ser usado para calcular intervalos de valores próprios.

Aqui você pode encontrar várias implementações de bibliotecas JDQR e JDQZ (incluindo uma interface C, à qual você deve poder vincular a partir do R).

Último folego
fonte
Não pude encontrar nenhuma literatura que declare explicitamente que o método de Jacobi-Davidson funciona para uma matriz geral real.
poder
A menos que todo artigo indique explicitamente uma restrição e o argumento de convergência se baseie na restrição que não importa.
Death Breath
Aqui está outra explicação do JD. As matrizes consideradas são completamente gerais. Nenhuma estrutura especial é explorada e os resultados específicos das matrizes hermitianas são comparados e contrastados; por exemplo, a convergência para matrizes gerais é quadrática, mas cúbica para matrizes hermitianas.
Death Breath
obrigado por isso. Não encontrei nenhum código C para uma matriz geral, portanto terei que escrever o meu. Os links para os algoritmos parecem ser apenas para matrizes hermetianas.
poder
1
@power, você também não encontrará na literatura um resultado que declare que as implementações padrão do QR convergem para uma matriz geral real - que é um problema em aberto e, de fato, há pouco tempo foi encontrado um contra-exemplo para o código QR no LAPACK.
Federico Poloni 22/10
2

Na sua postagem original, você diz:

"Eu também tentei usar a interface ARPACK no pacote igraph R. No entanto, ocorreu um erro para uma das minhas matrizes."

Eu estaria interessado em saber mais sobre o erro. Se você puder disponibilizar essa matriz publicamente em algum lugar, eu estaria interessado em experimentar o ARPACK nela.

Com base no que li acima, eu esperaria que o ARPACK fizesse um trabalho muito bom ao extrair os maiores (ou alguns dos maiores) autovalores de uma matriz esparsa. Para ser mais específico, eu esperaria que os métodos da Arnoldi funcionassem bem nesse caso e, é claro, é nisso que o ARPACK se baseia.

A lenta convergência do método de potência quando existem autovalores próximos na região de interesse foi mencionada acima. Arnoldi aprimora isso iterando com vários vetores em vez do método power.

Bill Greene
fonte
Vou ver se consigo encontrar meu trabalho naquela época. Eu trabalhei nisso há um ano.
potência
0

Não é a maneira mais rápida, mas uma maneira razoavelmente rápida é simplesmente atingir um vetor (inicialmente aleatório) com a matriz repetidamente e normalizar a cada poucos passos. Eventualmente, convergirá para o maior vetor próprio, e o ganho na norma para uma única etapa é o valor próprio associado.

Isso funciona melhor quando o maior autovalor é substancialmente maior que qualquer outro autovalor. Se outro autovalor está perto em magnitude ao maior, isso vai levar um tempo para convergir, e isso pode ser difícil determinar se ele tem convergido.

Dan
fonte
1
Obrigado Dan, no entanto: Nas minhas matrizes, alguns dos outros autovalores terão uma magnitude semelhante (se não a mesma) que a maior. O seu método é semelhante à Iteração de energia e Iteração de quociente de Rayleigh? Batterson e Smillie (1990) escrevem que, para algumas matrizes não simétricas, a iteração de quociente de Rayleigh não convergirá. Batterson, S., Smillie, J (1990) "Quociente de Rayleigh Iteração para não simétrica Matrizes", Mathematics of Computation, vol 55, núm 191, P 169-178
poder
Se outros valores próprios tiverem a mesma magnitude que o maior ... então esses valores também não são "o maior" também?
Ely
@EMS: Eles ainda seriam "maiores valores próprios", mas a presença de mais de um maior ainda mataria a convergência.
Dan
Só estou me perguntando para qual valor próprio você deseja que ele converja. Coisas como o quociente de Rayleigh / método de energia são entendidas quando existe um maior valor próprio distinto. Sua pergunta pede para encontrar o maior valor próprio, mas parece que isso não está realmente definido para o seu problema. Estou apenas enganado pelo título do post.
Ely
-1

O pacote RARPACK funciona para mim. E parece ser muito rápido, pois é apenas uma interface para o ARPACK, o pacote padrão para álgebra linear esparsa (o que significa calcular alguns autovalores e autovetores).

HoangDT
fonte
Bem-vindo ao SciComp! Como a pergunta indica, o ARPACK não funciona para o OP, portanto, essa resposta não é realmente útil.
Christian Clason
@HoangDT Esta questão é anterior ao rARPACK
power