The electromagnetic scattering of buried gyrotropic anisotropic media is crucial for resource exploration and environmental monitoring. However, the existing analytical solutions for electromagnetic scattering of a gyrotropic anisotropic sphere are primarily limited to free-space cases due to computational complexity. To address this limitation, an analytical solution that combines spherical vector wave functions (SVWFs), the T-matrix method, the image method, and the addition theorem of SVWFs is proposed in this work. The proposed method is detailed as follows. The transmitted field of a vertically incident plane wave transmitting through the ground serves as the first incident field on the gyrotropic anisotropic sphere, which can be expanded in terms of SVWFs. Using the analytical solution for a gyrotropic anisotropic sphere in free space, expressions for the internal electromagnetic field are derived. Based on the orthogonality of the SVWFs in the surface of the buried gyrotropic anisotropic sphere, the first scattered field is obtained. This scattered field then acts as the incident field on the ground, and its reflection is calculated using the image method. The reflected field can then serve as the secondary incident field for the dielectric sphere, and this process is repeated iteratively until the field components on the ground converge.Unlike the existing methods of computing the field at a fixed point for buried homogeneous cylinder or isotropic sphere, the proposed method computes the electric field distribution along a line L on the ground, which is parallel to both the Y-axis and the sphere’s central axis. The comparison of the results from the proposed method with FEKO simulation results shows their excellent agreement with each other, with an average relative error below 0.1%, thereby validating the correctness of the proposed analytical solution. Moreover, compared with FEKO simulation method, the proposed analytical method indicates a significant advantage in computational efficiency. Using the analytical model established in this work, the influence of incident wave frequency, buried depth and other parameters on the distribution of electric field along the Y-axis is also analyzed in detail. These findings provide practical value for enhancing the accuracy of geological exploration and the reliability of environmental monitoring.