Combien de régularisation ajouter pour rendre SVD stable?

10

J'ai utilisé le SVD d'Intel MKL ( dgesvdvia SciPy) et j'ai remarqué que les résultats sont considérablement différents lorsque je change de précision entre float32et float64lorsque ma matrice est mal conditionnée / n'est pas pleine. Existe-t-il un guide sur le montant minimum de régularisation à ajouter pour rendre les résultats insensibles aux float32-> float64changements?

En particulier, en faisant , je vois que la norme L de V T X se déplace d'environ 1 lorsque je change de précision entre et . La norme L 2 de A est 10 5 et elle a environ 200 valeurs propres nulles sur 784 au total.A=UDVTLVTXfloat32float64L2A105

Faire SVD sur avec λ = 10 - 3 a fait disparaître la différence.λI+Aλ=103

Yaroslav Bulatov
la source
Quelle est la taille d'une matrice N × N A pour cet exemple (est-ce même une matrice carrée)? 200 valeurs propres nulles ou valeurs singulières? Une norme Frobenius | | A | | F pour un exemple représentatif serait également utile. NN×NA||A||F
Anton Menshov
Dans ce cas, la matrice 784 x 784, mais je suis plus intéressé par la technique générale pour trouver une bonne valeur de lambda
Yaroslav Bulatov
Alors, la différence de uniquement dans les dernières colonnes correspond-elle aux valeurs singulières nulles? V
Nick Alger
2
S'il y a plusieurs valeurs singulières égales, le svd n'est pas unique. Dans votre exemple, je suppose que le problème vient des multiples valeurs singulières nulles et qu'une précision différente conduit à un choix différent de la base de l'espace singulier respectif. Je ne sais pas pourquoi ça change quand on régularise ...
Dirk
1
... qu'est-ce que ? X
Federico Poloni

Réponses:

1

Bien que la question ait une excellente réponse, voici une règle de base pour les petites valeurs singulières, avec un tracé.

Si une valeur singulière est non nulle mais très petite, vous devez définir sa réciproque à zéro, car sa valeur apparente est probablement un artefact d'erreur d'arrondi, pas un nombre significatif. Une réponse plausible à la question "comment petit est petit?" consiste à éditer de cette façon toutes les valeurs singulières dont le rapport à la plus grande est inférieur à fois la précision de la machine ϵ .Nϵ

- Recettes numériques p. 795

Ajouté: les deux lignes suivantes calculent cette règle empirique.

#!/usr/bin/env python2

from __future__ import division
import numpy as np
from scipy.sparse.linalg import svds  # sparse, dense or LinOp

#...............................................................................
def howsmall( A, singmax=None ):
    """ singular values < N float_eps sing_max  may be iffy, questionable
        "How small is small ?"
        [Numerical Recipes p. 795](http://apps.nrbook.com/empanel/index.html?pg=795)
    """
        # print "%d singular values are small, iffy" % (sing < howsmall(A)).sum()
        # small |eigenvalues| too ?
    if singmax is None:
        singmax = svds( A, 1, return_singular_vectors=False )[0]  # v0=random

    return max( A.shape ) * np.finfo( A.dtype ).eps * singmax


La matrice de Hilbert semble être largement utilisée comme cas de test pour l'erreur d'arrondi:

entrez la description de l'image ici

Ici, les bits de poids faible dans les mantisses de la matrice de Hilbert sont mis à zéro A.astype(np.float__).astype(np.float64), puis np.linalg.svdexécutés float64. (Les résultats avec svdtous float32sont à peu près les mêmes.)

Une simple troncature à float32pourrait même être utile pour débruiter des données de grande dimension, par exemple pour la classification des trains / tests.

De vrais cas de test seraient les bienvenus.

denis
la source
btw, scipy semble ajouter un facteur de 1e3 pour float32 et 1e6 pour float64, curieux d'où ils viennent
Yaroslav Bulatov
@Yaroslav Bulatov, numpyet scipy.linalg.svdappelez LAPACK gesdd , voir le paramètre JOBRdans dgejsv: "Spécifie la plage pour les valeurs singulières. Délivre la licence pour définir zéro petites valeurs singulières positives si elles sont en dehors ..." ( scipy.sparse.linalg.svdsencapsule ARPACK et a un paramètre tol, Tolérance pour les valeurs singulières.)
denis
13

A=ATM=UΣVT

H=[0MMT0]=[U00V][0ΣΣ0][U00V]T

ϵ>0

Aϵ=[1ϵϵ1]=VΛϵVT,Bϵ=[1+ϵ001ϵ]=UΛϵUT
Λϵ=diag(1+ϵ,1ϵ)
V=12[1111],U=[1001].
AϵBϵVUϵ>0U,VUV

M0=U0Σ0V0Tfloat64Mϵ=UϵΣϵVϵTfloat32Σ0,Σϵϵ107U0,UϵV0,Vϵ

Richard Zhang
la source
Cet exemple provient-il de: users.math.msu.edu/users/markiwen/Teaching/MTH995/Papers/… ?
Memming
1
Voilà une excellente référence. Je ne sais pas, j'ai appris cet exemple particulier il y a de nombreuses années en cours de mathématiques :-)
Richard Zhang