The uncertainty of subsurface resistivity can be determined using Bayesian inversion during electromagnetic exploration. However, the extensive sampling necessary for Bayesian inversion may result in high computational costs and prevent its application in complex, large-scale 2D and 3D electromagnetic inversions. For the first attempt, we implemented a 2.5D (i.e., 3D source with 2D geological model) Bayesian inversion to analyze ground transient electromagnetics method (TEM) and made some adaptations and improvements. To be specific, we employed geostatistical-based random modeling as a prior constraint strategy and incorporated additional prior information through model mapping to improve its performance in adapting to complex terrains. To perform pre-inversion of prior parameters, we utilized a spatially-constrained trans-dimensional Bayesian inversion approach, extracting priors for the formal inversion from the posterior of the pre-inversion. This could not only reduce computation time but also lessen the formal Bayesian inversion's reliance on prior selection. More importantly, to obtain more realistic uncertainty assessments, we designed a likelihood suited to the noise characteristics of TEM. The extensive statistical analysis of posterior results provides richer insights than previous methods, showing Bayesian inversion's superiority in both uncertainty reliability and resistivity recovery effectiveness over deterministic inversion and continuous 1D Bayesian inversion. Comparative studies with multiple simulated models and applications to complex terrain goaf areas validate the effectiveness and superiority of the proposed algorithm.