  We consider sample covariance matrices $S_N=\frac{1}{p}\Sigma_N^{1/2}X_NX_N^*
\Sigma_N^{1/2}$ where $X_N$ is a $N \times p$ real or complex matrix with
i.i.d. entries with finite $12^{\rm th}$ moment and $\Sigma_N$ is a $N \times
N$ positive definite matrix. In addition we assume that the spectral measure of
$\Sigma_N$ almost surely converges to some limiting probability distribution as
$N \to \infty$ and $p/N \to \gamma >0.$ We quantify the relationship between
sample and population eigenvectors by studying the asymptotics of functionals
of the type $\frac{1}{N} \text{Tr} (g(\Sigma_N) (S_N-zI)^{-1})),$ where $I$ is
the identity matrix, $g$ is a bounded function and $z$ is a complex number.
This is then used to compute the asymptotically optimal bias correction for
sample eigenvalues, paving the way for a new generation of improved estimators
of the covariance matrix and its inverse.
