Why are eigen and svd decompositions of a covariance matrix based on sparse data yielding different results?

Question

I am trying to decompose a covariance matrix based on a sparse / gappy data set. I'm noticing that the sum of lambda (explained variance), as calculated with svd, is being amplified with increasingly gappy data. Without gaps, svd and eigen yeild the same results.

This does not seem to happen with an eigen decomposition. I had been leaning towards using svd because the lambda values are always positive, but this tendency is worrying. Is there some sort of correction that needs to be applied, or should I avoid svd altogether for such a problem.

###Make complete and gappy data set
set.seed(1)
x <- 1:100
y <- 1:100
grd <- expand.grid(x=x, y=y)

#complete data
z <- matrix(runif(dim(grd)[1]), length(x), length(y))
image(x,y,z, col=rainbow(100))

#gappy data
zg <- replace(z, sample(seq(z), length(z)*0.5), NaN)
image(x,y,zg, col=rainbow(100))


###Covariance matrix decomposition
#complete data
C <- cov(z, use="pair")
E <- eigen(C)
S <- svd(C)

sum(E$values)
sum(S$d)
sum(diag(C))


#gappy data (50%)
Cg <- cov(zg, use="pair")
Eg <- eigen(Cg)
Sg <- svd(Cg)

sum(Eg$values)
sum(Sg$d)
sum(diag(Cg))



###Illustration of amplification of Lambda
set.seed(1)
frac <- seq(0,0.5,0.1)
E.lambda <- list()
S.lambda <- list()
for(i in seq(frac)){
    zi <- z
    NA.pos <- sample(seq(z), length(z)*frac[i])
    if(length(NA.pos) > 0){
        zi <- replace(z, NA.pos, NaN)
    }
    Ci <- cov(zi, use="pair")
    E.lambda[[i]] <- eigen(Ci)$values
	S.lambda[[i]] <- svd(Ci)$d
}


x11(width=10, height=5)
par(mfcol=c(1,2))
YLIM <- range(c(sapply(E.lambda, range), sapply(S.lambda, range)))

#eigen
for(i in seq(E.lambda)){
    if(i == 1) plot(E.lambda[[i]], t="n", ylim=YLIM, ylab="lambda", xlab="", main="Eigen Decomposition")
    lines(E.lambda[[i]], col=i, lty=1)
}
abline(h=0, col=8, lty=2)
legend("topright", legend=frac, lty=1, col=1:length(frac), title="fraction gaps")

    #svd
for(i in seq(S.lambda)){
    if(i == 1) plot(S.lambda[[i]], t="n", ylim=YLIM, ylab="lambda", xlab="", main="Singular Value Decomposition")
    lines(S.lambda[[i]], col=i, lty=1)
}
abline(h=0, col=8, lty=2)
legend("topright", legend=frac, lty=1, col=1:length(frac), title="fraction gaps")

Answer 1

You need to do the sum of the absolute value of eigen values i.e, sum(abs(Eg$values)) and compare it with the sum of the singular values. They would be equal.

The reason is that if you multiply the rows or columns that correspond to the negative eigenvalues with $-1$, then the eigen-value of the new matrix becomes positive and the orthogonality of the eigen-vectors is not disturbed.

The proof of the converse of this beautiful theorem appeared in The algebra of hyperboloids of revolution, Javier F. Cabrera, Linear Algebra and its Applications, Princeton University (now at Rutgers).

Another way to reason this is by the fact that sqrt(eigen(t(Cg) %*% Cg)) are equal to the singular values of Cg. But when the eigenvalues are negative, the data must be represented in a hermitian form with the complex plane taken into account, which is what was missed in the original formulation, i.e the data formed by the symmetric square-root of the matrix with negative eigen values would have complex entries.