I am calculating the linearly independent modes (and their uncertainties) of a large set of correlated variables by diagonalizing their covariance matrix, $C$, using GSL's eigen_symmv, which finds the eigenvalues, which I use as the diagonal of $D$.
https://www.gnu.org/software/gsl/manual/html_node/Real-Symmetric-Matrices.html
If I calculate the covariance elements over the entire data set, the algorithm is stable, and I always get positive eigenvalues in $D$.
However, I need to look at smaller sub-sets of the data set. I can apply the same method to smaller samples from the whole distribution, but as I approach the size of the slice I need to look at, I start to get more negative eigenvalues in $D$. (I also get more anomalously small eigenvalues).
I have checked that my input matrix, $C$, is always real and symmetric. And the elements of $C$ over a small subset, qualitatively, look similar to the $C$ I get when using the whole data set. I am not exceeding the range of the double data type, and there are no zero or undefined elements anywhere in the matrix.
Shouldn't I always be getting non-negative eigenvalues if $C$ is real and symmetric? Is this a math problem or a numerical problem?