In this paper, the transition joint probability density function of the solution of geometric Brownian motion (GBM) equation is obtained via Lie group theory of differential equations (DEs). Lie symmetry analysis is applied to find new solutions for time-fractional Fokker–Planck–Kolmogorov equation of GBM. This analysis classifies the forms of the solutions for the equation by the similarity variables arising from the symmetry operators. Finally, an analytic method called invariant subspace method is applied in order to find another exact solution.