SVD pour trouver la plus grande valeur propre d'une matrice 50x50 - est-ce que je perds beaucoup de temps?

J'ai un programme qui calcule la plus grande valeur propre de nombreuses matrices 50x50 symétriques réelles en effectuant des décompositions de valeurs singulières sur chacune d'entre elles. Le SVD est un goulot d'étranglement dans le programme.

Existe-t-il des algorithmes beaucoup plus rapides pour trouver la plus grande valeur propre, ou l'optimisation de cette partie ne donnerait-elle pas beaucoup de retour sur investissement?

linear-algebra eigensystem Anna
la source

Pourriez-vous donner plus d'informations sur vos matrices, par exemple si vous savez quoi que ce soit sur leur structure, la plage de leurs valeurs propres ou leur similitude?

Pedro

C'est une matrice de covariance (

). Les tests montrent que toutes les valeurs propres sauf les 5 ou plus grandes sont proches de zéro et que la plus grande valeur propre est au moins ~ 20% plus grande que la deuxième plus grande. Puisqu'il y a beaucoup de valeurs propres proches de zéro, je suppose que la plage n'est pas importante? Il peut être redimensionné dans n'importe quelle plage. L'échelle que j'utilise actuellement me donne une plage de 150 à 200.

X X^{T}

$XX^T$

Anna

De plus, la matrice n'est pas très singulière, donc le problème SVD est bien conditionné.

Anna

Puisque

est symétrique et positif (semi) défini, vous pouvez utiliser la factorisation de Cholesky au lieu du SVD. La factorisation de Cholesky prend beaucoup moins de flops à calculer que la SVD, mais étant une méthode exacte, il faut toujours

flops.

X X^{T}

$XX^T$

O (n^{3})

$O(n^3)$

Ken

@Anna: Avez-vous essayé l'une des nombreuses approches proposées ici? Je serais assez curieux de savoir ce qui a le mieux fonctionné pour vous ...

Pedro

Réponses:

Selon la précision dont vous avez besoin pour la plus grande valeur propre, vous pouvez essayer d'utiliser l' itération de puissance .

Pour votre exemple spécifique, j'irais jusqu'à ne pas former explicitement, mais calculer à chaque itération. Le calcul de nécessiterait des opérations alors que le produit matrice-vecteur ne nécessite que . $A=XX^\mathsf{T}$ $x \leftarrow X(X^\mathsf{T}x)$ $A$ $\mathcal O(n^3)$ $\mathcal O(n^2)$

Le taux de convergence dépend de la séparation entre les deux plus grandes valeurs propres, donc cela peut ne pas être une bonne solution dans tous les cas,

Pedro
la source

Si la plus grande valeur propre est 20% plus grande que la suivante, l'itération de puissance devrait converger assez rapidement (toutes les autres valeurs propres sont amorties par un facteur 5/6 à chaque itération, vous obtenez donc un chiffre pour 13 itérations).

Wolfgang Bangerth

Les méthodes du sous-espace de Krylov sont strictement meilleures que les méthodes de puissance, car elles contiennent le vecteur de l'itération de puissance avec le même nombre d'itérations.

Jack Poulson

@JackPoulson: Oui, mais chaque itération est plus coûteuse à calculer ... Est-ce que cela en vaudrait vraiment la peine pour un si petit problème?

Pedro

@Pedro: bien sûr, les matvecs nécessitent un travail quadratique et le quotient de Rayleigh eigensolve et l'expansion ultérieure sont triviaux en comparaison.

Jack Poulson

Frais de code? Depuis que @JackPoulson a abordé la question, B. Parlett et al (1982) («On Estimating the Largest Eigenvalue with the Lanczos Algorithm») comparent la méthode de puissance, la méthode de puissance + l'accélération Aitken, et une application de Lanczos ciblant la plus grande valeur propre d'un réel symétrique (ou hermitien) pos. def. matrice. Ils concluent que la méthode de Lanczos est plus efficace si une précision même modeste (de la première valeur propre par rapport à la seconde) est nécessaire, et mieux pour éviter les erreurs de convergence.

hardmath

Si seulement 5 valeurs propres sont très significatives, l'algorithme de Lanczsos avec comme multiplication matrice-vecteur devrait donner une convergence linéaire rapide après 5 étapes initiales, d'où une valeur propre plus grande et assez précise avec peu d'itérations. $X(X^Tx)$

Arnold Neumaier
la source

Envisagez-vous (@ArnoldNeumaier) quelque chose comme ça , convenablement simplifié (

)? Il est intéressant de noter qu'il donne une approximation différente de Lanczos si un troisième vecteur est conservé, sur le même sous-espace de Krylov.

B = T = I

$B = T = I$

hardmath

Non; Je voulais dire l'algorithme Lanczsos standard, mais j'avais pressé écrit CG. Maintenant corrigé.

Arnold Neumaier

Pour une matrice semi-définie positive telle que cela peut valoir la peine d'accélérer la convergence avec un décalage de spectre . Autrement dit, un scalaire adapté est choisi et la méthode de la puissance est appliquée à au lieu d' . $A = XX^T$ $\mu$ $A - \mu I$ $A$

Quelques itérations de la méthode de puissance de base devraient vous donner une estimation approximative de la plus grande valeur propre . En supposant que la valeur propre dominante a la multiplicité 1 et que toutes les autres sont dans $||Ax||/||x||$ $\lambda_1$ $[0,\frac{5}{6} \lambda_1]$ $A - \frac{5}{12} \lambda_1 I$ $\frac{7}{12} \lambda_1$ $[\frac{-5}{12} \lambda_1, \frac{5}{12} \lambda_1]$

En d'autres termes, vous augmenteriez la dominance de la plus grande valeur propre de 20% sur la prochaine plus grande à 40% sur la prochaine plus grande (valeur absolue d'une) valeur propre. La convergence géométrique de la méthode de puissance s'accélérerait en conséquence. Une fois que la plus grande valeur propre de $A - \mu I$ est trouvé avec une précision suffisante, $\lambda_1$ est estimé en ajoutant le décalage $\mu$ qui avait été enlevé.

Notez que vous n'avez pas besoin de former explicitement $A - \mu I$ car $(A - \mu I)x = X(X^Tx) - \mu x$ peut toujours être calculé avec $O(n^2)$ effort.

hardmath
la source

This would seem to require having a good idea of what the magnitude of the second-largest eigenvalue is. How would you approximate it in such a case?

Pedro

@Pedro: Applying the shifted power iteration only requires an estimate of the largest eigenvalue

λ_{1}

$\lambda_1$ , but as with the unshifted power method, analysis shows the rate of convergence depends on

| λ_{2} | / | λ_{1} |

$|\lambda_2|/|\lambda_1|$ (throwing in abs. values that are unnecessary for the pos. semi-definite case at hand). In turn observed rates of convergence can be used to estimate

| λ_{2} | / | λ_{1} |

$|\lambda_2|/|\lambda_1|$ , and hence the size of

λ_{2}

$\lambda_2$ relative to

λ_{1}

$\lambda_1$ if desired. I was suggesting what benefit you'd see in a case such as Anna describes in her comments below the Question.

hardmath