Abstract: Elastic reverse time migration (RTM) uses the elastic wave equation to extrapolate multicomponent seismic data to the subsurface and separate the elastic wavefield into P- and S-waves. P- and S-wave separation is a necessary step in elastic RTM to avoid crosstalk between coupled wavefields. However, the current curl-divergence operator-based separation method has a polarity reversal problem in PS imaging, and vector separation methods often have separation artifacts at the interface, which affects the quality of the imaging stack. We propose a non-artifact P- and S-wave separation method based on the first-order velocity-strain equation. This equation is used for wavefield extrapolation and separation in the first-order staggered-grid finite-difference scheme, and the storage and calculation amounts are consistent with the classical first-order velocity-stress equation. The separation equation does not calculate the partial derivatives of the elastic parameters, and thus, there is no artifact in the separated P- and S-waves. During wavefield extrapolation, the dynamic characteristics of the reflected wave undergo some changes, but the transmitted wavefield is accurate; therefore, it does not affect the dynamic characteristics of the final migration imaging. Through numerical examples of 2D simple models, part SEAM model, BP model, and 3D 4-layer model, different wavefield separation methods and corresponding elastic RTM imaging results are analyzed. We found that the velocity-strain based elastic RTM can image subsurface structures well, without spike artifacts caused by separation artifacts, and without polarity reversal phenomenon of the PS imaging.