So I have read that using pairwise deletion can lead to a non-positive definite matrix. However, I tried to find a source which showed this mathematically. Can someone point me to that?
My intuition thus far is that the factor $ N - 1 $ by which I divide the covariance between two variables differs. Also that using different sample sizes leads to different levels of accuracy in $ \hat{\sigma}_{ij} $. But I would appreciate a more detailed explanation.
The simple explanation is that the mathematical proof that covariance matrices are positive semi-definite breaks down when using pairwise deletion. But to see with more details what is happening, let us look at a simple example. The simplest example possible needs 3 variables, since we need different pairs of observations being left out for different pairs of variables. Let us make a simple data example which has (approximately) a singular (that is semi-definite) covariance matrix, so that a small perturbation could make it non-definite, for instance, with a small, negative, determinant. Let the covariance matrix be $$ \Sigma = \begin{pmatrix} 1 & -\frac12 & -\frac12 \\ -\frac12 & 1 & -\frac12 \\ -\frac12 & -\frac12 & 1\end{pmatrix} $$ which you can check has determinant 0.
An example data set with approximately this covariance matrix is
XX
[,1] [,2] [,3]
1 9.29 10.86 9.85
2 9.10 11.03 9.87
3 12.23 9.17 8.60
4 10.91 9.63 9.46
5 9.66 9.90 10.43
6 10.03 9.73 10.23
7 9.61 10.62 9.77
8 9.54 7.97 12.48
9 10.63 9.76 9.61
10 9.00 11.32 9.68
with the following pairs plot:
Let us decide for a few variables we can declare as missing, NA. We need the covariances to get a little bit more negative, but without changing the variances too much. So look in the plot for some concordant pairs, that is, pairs contributing positively to the empirical covariance. That way the covariance will be more negative when this pairs is left out. But at the same time, choose pairs with values close to their variables means, so that the variance will not change too much. Looking at the plot shows a possible solution:
XXX <- XX
XXX[9, 1] <- NA
XXX[4, 3] <- NA
XXX[7, 2] <- NA
and then check that this works:
cov(XXX, use="pairwise")
[,1] [,2] [,3]
[1,] 1.0689000 -0.5853143 -0.5133893
[2,] -0.5853143 1.0732500 -0.6492607
[3,] -0.5133893 -0.6492607 1.0831194
> det(cov(XXX, use="pairwise"))
[1] -0.2521741