With the rapid proliferation of big data and high-dimensional datasets, distributed estimation methods have become indispensable for large-scale statistical inference. However, traditional centralized distributed computing frameworks often suffer from communication bottlenecks and privacy risks. To address these challenges, we propose a novel distributed stochastic optimization algorithm, termed gradient-based Markov subsampling Gradient Tracking with Variance Reduction (GMS-GT-VR). This algorithm seamlessly integrates gradient tracking and variance reduction techniques with a data-driven gradient-based Markov subsampling (GMS) strategy to solve large-scale distributed optimization problems over multi-agent networks efficiently. By leveraging a lightweight coordination mechanism for adaptive sampling, GMS-GT-VR effectively reduces both communication overhead and computational complexity, while achieving accelerated convergence. We rigorously establish its theoretical convergence rate of $\mathcal{O}(1/k)$ and conduct a detailed complexity analysis. To validate the theoretical findings, we perform extensive experiments on large-scale datasets. The results consistently demonstrate that GMS-GT-VR outperforms existing variance-reduction methods in terms of communication efficiency and convergence speed, achieving a significant reduction in the number of iterations required.