Reverse time migration (RTM) has become an extremely important imaging method due to its ability of handling complicated velocity models. In recent years, there has been an increasing need to perform RTM with anisotropic velocity models such as transversely isotropic media with a vertical symmetry axis (VTI). Finite difference (FD) is one of the most widely used numerical methods in RTM. However, it often suffers from serious numerical dispersion, resulting in significant loss of imaging resolution. Over the years, a nearly analytic discrete operator has been developed to approximate the high-order spatial partial differential operators, which is effective in suppressing numerical dispersion. Using this operator, we developed a nearly analytic symplectic (NAS) method to numerically solve acoustic VTI wave equations, which preserves the symplectic structure. We have investigated its properties including the stability condition, numerical dispersion relation, and computational efficiency. Also, we applied the NAS method to solve wave equations in RTM in acoustic VTI media. This new method not only reduced numerical dispersion and improved computational efficiency, compared to the conventional numerical methods such as the Lax-Wendroff correction (LWC) method, but it also was symplectic. Our forward-modeling results showed that the computational efficiency of NAS is about four times and its overall memory requirement is about 75% of those of LWC. The RTM results of the 2D prestack Hess acoustic VTI model showed that RTM by using NAS obtained better image quality than that by using LWC when coarse grids were used, which demonstrated that the NAS method can be promising in large-scale anisotropic RTM.