Précision de l'inversion de matrice
Théorie
L’inversion de matrices de grande taille est difficile. Newman(1974) propose une méthode pour tester la précision de l’inversion.
Soit \(A\) la matrice (carrée, non singulière) à inverser. \(X\) est l’approximation de \(A^{-1}\) obtenue par calcul. \(N(Y)=\sqrt{\sum_{i,j}{y_{i,j}^2}}\) est la norme euclidienne d’une matrice \(Y=(y_{i,j})\).
\(R = I - AX\) est l’écart (idéalement nul) entre la matrice identité et le produit de \(A\) et son inverse calculé.
La norme de \(R\) n’est pas un bon indicateur de la précision de l’inversion (Newman, 1974, section 3). En revanche, la norme de \(A^{-1} - X\) en est un, mais \(A^{-1}\) est inconnu. Newman montre que \[\frac{N(XR)}{1 + N(XR)} \leq N(A^{-1} - X) \leq \frac{N(XR)}{1 - N(XR)},\] ce qui permet d’encadrer la valeur de \(N(A^{-1} - X)\) avec les valeurs calculées de \(X\) et \(R\).
L’objectif est que cette valeur soit petite, idéalement nulle si \(A^{-1} = X\), par rapport à la norme de \(X\).
Exemple
Le paramètre size
fixe la taille (nombre de ligne et de colonnes) des matrices utilisées, ici une matrice de 1000 x 1000.
size <- 1000
Une matrice est crée:
A <- matrix(rnorm(size^2), nrow = size)
Elle est inversée numériquement (calcul de \(X\)) puis \(R\) et \(XR\) sont calculées:
X <- solve(A)
AX <- A %*% X
R <- diag(size) - AX
XR <- X %*% R
La fonction norme euclidienne est définie:
norm.matrix <- function(X) {
sqrt(sum(X^2))
}
Les normes des matrices sont calculées:
# Matrice à inverser
norm.matrix(A)
## [1] 1000.815
# Inversion numérique
norm.matrix(X)
## [1] 27.63621
# Ecart d'inversion
norm_R <- norm.matrix(R)
norm_XR <- norm.matrix(XR)
Finalement, l’encadrement de la norme de \(A^{-1} - X\) est fourni:
norm_XR / (1 + norm_R)
## [1] 3.514881e-12
norm_XR / (1 - norm_R)
## [1] 3.514881e-12
La norme de \(A^{-1} - X\) est ici très proche de 0 : l’inversion de matrice dans R fonctionne bien.
Référence
Newman, M. (1974). How to Determine the Accuracy of the Output of a Matrix Inversion Program. JOURNAL OF RESEARCH of the National Bureau of Standards - B. Mathematical Sciences, 788(2), 65–68.