Non-Gaussian diffusion of nanoparticles in complex media disrupts Einstein's picture of Brownian motion, and non-Gaussianity is thought to be closely related to diffusing diffusivity generated by spatiotemporal heterogeneities. However, the correlation between non-Gaussianity and the dynamics of heterogeneous environments in anomalous diffusion remains uncertain. Inspired by a recent study by Alexandre et al. [Phys. Rev. Lett. 130, 077101 (2023)10.1103/PhysRevLett.130.077101], we demonstrate that non-Gaussianity can be deciphered through the spatiotemporal evolution of heterogeneity-dependent diffusivity distribution. Using diffusion experiments in a linear temperature field and Brownian dynamics simulations, we found that short-time non-Gaussianity can be predicted based on the boundary ratio of the diffusivity distribution; the long-time non-Gaussianity either approaches an asymptotic value of -2 or scales with 1/t, depending on the dominance of particle migration. The temporal variation of non-Gaussianity is determined by an effective Péclet number, which represents a competition between the varying rate of diffusivity and the diffusivity of diffusivity and reveals whether the tail distribution expands or contracts. The tail is more Gaussian than exponential over long times, with exceptions significantly dependent on the diffusivity distribution. Our findings provide a versatile framework for understanding non-Gaussian diffusion in probability space, and shed light on establishing a diffusion spectrum in cells and characterizing nanomedicine transport in biological microenvironment using non-Gaussian statistics.