Using the replica method, we develop an analytical approach to compute the characteristic function for the probability $mathcal{P}_N(K,lambda)$ that a large $N times N$ adjacency matrix of sparse random graphs has $K$ eigenvalues below a threshold $lambda$. The method allows to determine, in principle, all moments of $mathcal{P}_N(K,lambda)$, from which the typical sample to sample fluctuations can be fully characterized. For random graph models with localized eigenvectors, we show that the index variance scales linearly with $N gg 1$ for $|lambda| > 0$, with a model-dependent prefactor that can be exactly calculated. Explicit results are discussed for Erdos-Renyi and regular random graphs, both exhibiting a prefactor with a non-monotonic behavior as a function of $lambda$. These results contrast with rotationally invariant random matrices, where the index variance scales only as $ln N$, with an universal prefactor that is independent of $lambda$. Numerical diagonalization results confirm the exactness of our approach and, in addition, strongly support the Gaussian nature of the index fluctuations.