Numerically determining rank of a matrix? I'm looking for a way to numerically determine rank of a covariance matrix. For instance, computing eigenvalues a gram matrix with rank 100 using scipy implementation, I see plot like below. Red dots represent negative eigenvalues. Basically there's a big drop in log of eigenvalue after rank 100 and half of the eigenvalues are negative. Is there a more robust way to turn this intuition into an algorithm?

BTW, this plot is for a covariance matrix of data file readable as np.genfromtxt(fname, delimiter= ",")
 A: Given $A\in\mathbb{R}^{m\times n}$ and its computed singular values $\sigma_1\geq\sigma_2\geq\cdots\geq\sigma_k\geq 0$, $k=\min\{m,n\}$, the rank is determined through the notion of a numerical rank. 
For an $\epsilon>0$ we define the $\epsilon$-numerical rank $r_\epsilon$ of $A$ as the largest integer $r\in[0,k]$ for which $\sigma_r\geq\epsilon\sigma_1$. Usual choice of $\epsilon$, supported by sensitivity and rounding error analyses, is $\epsilon=\max\{m,n\}\varepsilon$, where $\varepsilon$ is the machine precision of the given floating point arithmetic implementation (e.g., for double precision arithmetic, we have $\varepsilon=2^{-52}\approx 2.22\cdot 10^{-16}$.
In your case, since the given matrix is symmetric, you can use the absolute values of the computed eigenvalues in place of the singular values. 
If you won't need to keep the computed spectrum of your matrix for further calculations, you can also use the function numpy.linalg.matrix_rank with the same default behavior as described above (which is also used by the Matlab's implementation of the rank function).
