# How do I propagate error values through a matrix diagonalization?

In the general case, if I have a symmetric, invertible (real-valued) matrix where each entry in the matrix has a separate error associated with it (or, in the more general case, where each value has a separate positive and negative error associated with it), how can I rigorously propagate these error values through the eigendecomposition to get error values for each of the eigenvalues and each entry in the eigenvector matrix? For my current application, it is as critical that I avoid overestimating the error as it is that I avoid underestimating it.

-

The propagation will depend on the diagonalization algorithm--which might be a black box--as well as the multivariate distribution of the errors. Pursuing an analytical solution therefore looks unpromising. Why not just compute an empirical distribution? That is, draw a large number of variants of the original matrix from the hypothesized error distribution and diagonalize them. Study the output distribution of the eigenvectors and eigenvalues.

There are some subtleties, because there will not be a definite matching among the lists of eigenvalues. For instance, in one iteration the sorted eigenvalues might be $(1.0, 0.99, 0.17)$ and in the next they might be $(1.01, 0.98, 0.17)$. Is the $1.01$ in the latter a slight variation of the $1.0$ in the former, or perhaps has the $0.99$ been perturbed into $1.01$ and the $1.0$ into $0.98$? It is impossible to know. Thus, you need to characterize the multivariate distribution of multisets of eigenvalues rather than $n$-tuples of eigenvalues.

The same problem attaches to the eigenvectors, but it gets worse, because there is no unique normalization of the eigenvectors. (They are determined only up to sign.) However, these problems are no different in nature than the ambiguities present in other geometric problems such as characterizing the directions of linear features in a plane (which can be given only up to a multiple of 180 degrees) and so should not present any additional conceptual challenge; they are just going to be a nuisance.

Here is an example of the empirical distributions of the sorted eigenvalues of a 4 by 4 matrix, using 2500 draws from the error distribution.

This scatterplot matrix also shows the lines y=x on each plot to emphasize the constraints imposed by sorting the four eigenvalues.

-
Thanks for the great answer. I thought of this already, and I was afraid it might be the only way. As an outsider to the field of mathematics, it just seems to me that someone would have figured out by now how to do this analytically and independent of the eigendecomposition algorithm, but perhaps I'm just nieve as to the intricacies involved in the problem. –  David Hollman Oct 21 '11 at 19:35
@David Consider even the simplest non-trivial setting with 2 by 2 matrices. The eigenvalues are roots of a quadratic equation whose coefficients are the trace and determinant of the matrix. You would first have to propagate a 3-variate error distribution through these operations (which is easy for the trace and not easy at all for the determinant) but then you would have to propagate the resulting bivariate distribution through the root-finding step, which can be done only approximately, not analytically. Even in this case, simulation is the way to go. –  whuber Oct 21 '11 at 19:41