The exponentially decaying behavior of the momentum-momentum and momentumforce time correlation functions of Brownian motion at large times has been extensively used for the numerical evaluation of the friction coefficient γ from molecular dynamics simulations. We perform numerical analysis on these methods and address issues related to the appropriate choice of large time and the rate of convergence of these methods. To this end, we obtain asymptotic expansions of the time correlation functions with respect to the reduced mass μ of the Brownian particle. For two important limit procedures of achieving the Brownian limit, certain forms of the asymptotic expansion of the Mori memory function K(t) are introduced by physical arguments, and then the asymptotic expansions of the time correlation functions are expressed in terms of the limit of K(t), i.e., K 0 (t) = limμ→∞ K(t), and the next-order correction term K 1 (t). For the infinite mass limit case, we show that the numerical methods of estimating γ from the exponential decay of the time correlation functions produce γ + O(μ −1), where the first-order correction depends on the microscopic nature of K 0 (t) (i.e., deviation of K 0 (t) from the Dirac delta function γδ(t)) as well as the contribution of K 1 (t) (i.e., deviation of K(t) from K 0 (t)). We also analyze the ratio of the momentum-force correlation function to the momentum-momentum correlation function, which gives instantaneous exponential decay rate. For the thermodynamic limit case, we consider the Rayleigh gas system to investigate the finite-volume effect due to the boundary conditions and to demonstrate that the lowest-order terms of the asymptotic expansions may fail to describe some characteristic behavior of the time correlation functions. We perform systematic molecular dynamics simulations of the Rayleigh gas system to confirm the theoretical predictions.