We consider sample covariance matrices constructed from real or complex i.i.d. variates with finite 12th moment. We assume that the population covariance matrix is positive definite and its spectral measure almost surely converges to some limiting probability distribution as the number of variables and the number of observations go to infinity together, with their ratio converging to a finite positive limit. We quantify the relationship between sample and population eigenvectors, by studying the asymptotics of a broad family of functionals that generalizes the Stieltjes transform of the spectral measure. 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.