Least-squares reverse time migration (LSRTM) can eliminate imaging artifacts in an iterative way based on the concept of inversion, and it can restore imaging amplitude step by step. LSRTM can provide a high-resolution migration section and can be applied to irregular and poor-quality seismic data and achieve good results. Steeply dipping reflectors and complex faults are imaged by using wavefield extrapolation based on a two-way wave equation. However, the high computational cost limits the method\'s application in practice. A fast approach to realize LSRTM in the imaging domain is provided in this paper to reduce the computational cost significantly and enhance its computational e? ciency. The method uses the Kronecker decomposition algorithm to estimate the Hessian matrix. A low-rank matrix can be used to calculate the Kronecker factor, which involves the calculation of Green\'s function at the source and receiver point. The approach also avoids the direct construction of the whole Hessian matrix. Factorization-based LSRTM calculates the production of low-rank matrices instead of repeatedly calculating migration and demigration. Unlike traditional LSRTM, factorization-based LSRTM can reduce calculation costs considerably while maintaining comparable imaging quality. While having the same imaging effect, factorization-based LSRTM consumes half the running time of conventional LSRTM. In this regard, the application of factorization-based LSRTM has a promising advantage in reducing the computational cost. Ambient noise caused by this method can be removed by applying a commonly used filtering method without significantly degrading the imaging quality.