## 1 Introduction
Hence, suitable environmental variables should be selected to determine the spatial variation characteristics of SOC in flat areas and utilize them for SOC mapping and precision agriculture. Several scholars have identified many other alternative factors to respond to the spatial variation of soil properties and solve the difficulty in selecting environmental factors in plain areas or flat terrain regions. [PERSON] et al. [29] and [PERSON] et al. [30] presented a new land surface dynamic feedback (LSDF)model by comparing the temporal responses to a rainfall event to map soil texture, SOC, and other properties. The LSDF model combined with land surface spectral or temperature variations uses short-time remote sensing images to predict soil properties [31]. The contradiction between return time and spatial resolution limits the development of high-precision soil maps. Hyperspectral images have been used to quantitatively predict soil properties through the spectral reflectance of the surface soil [32; 33]. However, surface vegetation and scant hyperspectral images hinder their use in large areas. Agricultural land occupies approximately 38.18% of the world area based on the data of World Bank in 2016 and the main land use type among all land use types. Agricultural production and activity constantly influence the change of SOC storage [34; 35]. Thus, the spatial and temporal distribution rules of SOC in agricultural lands should be investigated. SOC as the main soil fertilizer influences the soil structure and crop growth [36]. With the improvement of remote sensing technology, increasing studies have focused on recording the growth status of crops in different phenological periods [37; 38; 39]. VI is designed to enhance the contribution of vegetation properties and allow the reliable spatial and temporal intercomparisons of terrestrial photosynthetic activity and canopy structural variations [40]. [PERSON] et al. [41] concluded that an obvious exponential relationship exists between broadband NDVI and gross primary productivity (GPP), while a linear relationship occurs between broadband NDVI and the fraction of absorbed photosynthetically active radiation (fAPAR). This condition indicated that the VIs of remote sensing images can record crop variation. Thus, from this perspective, remote sensing VIs may be used to reflect the spatial variation of soil properties (e.g., SOC). However, at present, many studies have only considered environmental variables (e.g., NDVI time series characteristics) at a certain time when modeling with environmental elements and ignored the variability with time. [PERSON] et al. [42] used NDVI data from April 1987 as the parameter in modeling. [PERSON] and [PERSON] [43] adopted an NDVI map from a July 26, 2004 image when mapping the spatial distribution of SOC. [PERSON] et al. [44] mapped SOC using data mining techniques with some ancillary data that included ratio VI (RVI), soil-adjusted VI (SAVI), and NDVI on March 28, 2013. However, the temporal characteristics of environmental variables should be considered because the variation of soil physical and chemical properties requires time. In previous studies, NDVI was used as an important index to monitor the growth status and cover the vegetation degree. NDVI played an important role in remote sensing applications. However, considering only NDVI at a certain point in time may cause unrealistic results because the sequential feature of NDVI data can reveal additional information about the study object for a time and lead a highly comprehensive approach. Numerous studies based on VI time series have been applied in various fields. [PERSON] et al. [45] extracted winter wheat information on the basis of time series NDVI in the Guanzhong area. [PERSON] et al. [46] analyzed the land damage and recovery process in a rare earth mining area using multisource sequential NDVI.