Multi-Feature Driver Variable Fusion Downscaling TROPOMI Solar-Induced Chlorophyll Fluorescence Approach

Abstract

Solar-induced chlorophyll fluorescence (SIF), as a direct indicator of vegetation photosynthesis, offers a more accurate measure of plant photosynthetic dynamics than traditional vegetation indices. However, the current SIF satellite products have low spatial resolution, limiting their application in fine-scale agricultural research. To address this, we leveraged MODIS data at a 1 km resolution, including bands b1, b2, b3, and b4, alongside indices such as the NDVI, EVI, NIRv, OSAVI, SAVI, LAI, FPAR, and LST, covering October 2018 to May 2020 for Shandong Province, China. Using the Random Forest (RF) model, we downscaled SIF data from 0.05° to 1 km based on invariant spatial scaling theory, focusing on the winter wheat growth cycle. Various machine learning models, including CNN, Stacking, Extreme Random Trees, AdaBoost, and GBDT, were compared, with Random Forest yielding the best performance, achieving R² = 0.931, RMSE = 0.052 mW/m²/nm/sr, and MAE = 0.031 mW/m²/nm/sr for 2018–2019 and R² = 0.926, RMSE = 0.058 mW/m²/nm/sr, and MAE = 0.034 mW/m²/nm/sr for 2019–2020. The downscaled SIF products showed a strong correlation with TanSIF and GOSIF products (R² > 0.8), and consistent trends with GPP further confirmed the reliability of the 1 km SIF product. Additionally, a time series analysis of Shandong Province’s wheat-growing areas revealed a strong correlation (R² > 0.8) between SIF and multiple vegetation indices, underscoring its utility for regional crop monitoring.

1. Introduction

Solar-induced chlorophyll fluorescence (SIF) offers a more direct method for assessing vegetation photosynthesis compared to traditional vegetation indices [1,2,3,4], particularly in estimating gross primary productivity (GPP) and monitoring vegetation health status [5,6,7]. However, satellite-derived SIF data often suffer from low low spatial resolution and sparse sampling, which hinders the detection of subtle variations in small-scale regions. This issue is especially critical in agricultural and ecological research [8,9]. Therefore, conducting downscaling studies on SIF data is of great importance to enhance the applicability in these fields [10,11,12].

Satellite-based SIF products are widely utilized for large-scale vegetation monitoring, productivity estimation, and spatial–temporal pattern analysis [13,14,15,16]. Nevertheless, these products are limited by their coarse resolution and fragmented data coverage. Key satellites involved in SIF observation include the Scanning Imaging Absorption Spectrometer for the Atmosphere (SCIAMACHY), the Greenhouse Gas Observing Satellite (GOSAT), the Global Ozone Monitoring-2 (GOME-2), the Orbiting Carbon Observatory-2 (OCO-2), the China Carbon Dioxide Observing Satellite (TanSat), and the Tropospheric Monitoring Instrument (TROPOMI). GOSAT was the first satellite to provide global SIF data products [17], with a spatial resolution of 10.5 km and monthly scale products at a resolution of 2° × 2°. SCIAMACHY continuously samples data with an average resolution of 1.5° × 1.5° for daily products, while GOME-2 offers SIF products at 0.5° × 0.5° [18]. On the other hand, the SIF products from OCO-2 and TanSat have smaller sampling footprints (1.3 × 2.25 km for OCO-2 and 2 × 2 km for TanSat), but their spatial discretization makes it challenging to achieve global coverage. Although TROPOMI provides continuous SIF data, which is advantageous compared to the other satellites, its spatial resolution remains relatively coarse at 3.5 km × 7 km. This still poses challenges for applications that require high spatial resolution, such as monitoring crop growth, estimating crop yields, and tracking pests and diseases at small regional scales [19].

Reconstructing SIF can address the challenges of low spatial resolution and data discretization. For example, Ma Yan used MODIS reflectance, NDVI, COS (SZA) and reanalyzed meteorological temperature data to develop a downscaling model for the 0.5° × 0.5° GOME-2 SIF product through the Random Forest model. This model enabled the generation of a global DSIF product with a spatial resolution of 0.05° × 0.05° [20]. Similarly, the discrete strip TanSat SIF (2 km × 2 km) was also generated into the global coverage spatial scale extension product TanSIF with a spatial resolution of 0.05° × 0.05° based on the same driving variables and algorithms [21]. Li and Xiao applied the Cubist model, using EVI, PAR, VPD, and T to reconstruct the discontinuous OCO-2 SIF product into the long-term global GOSIF product at a spatial resolution of 0.05° × 0.05° [22]. Yu et al. combined the OCO-2 and MODIS MCD43A4 products, incorporating one to seven reflectance bands, to reclassify land cover types through Frankenberg’s ideas [23], producing 16-day SIF products with a resolution of 0.05° [24], and it was found that they have significant advantages in plant growth analysis. Jeong et al. utilized OCO-3 based on the continuous GK-2A near-infrared, blue, rate, and red bands, as well as shortwave radiation and VPD as the driving variables to study the reconstruction of the SIF for the continents of East Asia and Oceania using the XGBoost model and to produce the dataset GEOSIF [25]. Despite these advancements, most SIF downscaling studies have been applied at large spatial scales (e.g., global), making it challenging to apply them to smaller regional areas.

As one of the world’s major food crops, the widespread cultivation of wheat has a significant impact on the global agricultural economy [26]. The growth cycle of winter wheat is closely related to climate change, and accurate monitoring of its growth status is crucial for improving yield and quality [27]. In recent years, SIF has demonstrated a unique potential for application in crop yield estimation [28,29]. Therefore, real-time monitoring and prediction of winter wheat growth can be effectively improved by downscaled SIF data, which can provide a scientific basis for farmers to make more accurate decisions in the planting process, thus optimizing agricultural management and improving yields [30]. To solve the above problems, this study uses MODIS data with 1 km spatial resolution as driving variables, combines these variables with TROPOMI enhanced SIF (eSIF) with a resolution of 0.05° through resampling, utilizes a variety of models to train the training sample dataset, compares the accuracy metrics of the model calculations, and selects the optimal algorithmic model to construct the downscale model to improve the accuracy of downscaled SIF and invert the 1 km resolution SIF products in Shandong Province. In addition, we verified the consistency of 1 km resolution SIF with other SIF satellite products from the perspective of a time series, verified the reliability of 1 km SIF from the standpoint of the SIF-GPP vegetation mechanism relationship, and finally, applied 1 km SIF to analyze the growth of winter wheat in Shandong Province. In this study, we aimed to improve the spatial resolution of SIF data to make it more applicable to small-scale agricultural areas, thereby improving the ability to monitor and manage winter wheat cultivation. This improvement is particularly important in areas such as pest and disease control, yield estimation, and crop growth monitoring. The downscaled SIF data provide more accurate and regionalized SIF information, which can help to better manage winter wheat growth in Shandong Province [31].

2. Materials and Methods

2.1. Study Area

The research area for this study is Shandong Province, located in Eastern China (Figure 1). It spans approximately 158,100 km², situated between latitudes 34°23′-38°17′ N and longitudes 114°48′–122°42′ E [32]. With a coastal strip in the northeastern portion of the province and a high elevation in the center, Shandong Province boasts a diverse and complex topography that includes plains, hills, basins, mountains, lakes, and more. One of the most notable features of Shandong Province’s terrain is the high level of naturalization within it. Winter wheat and summer corn are the primary crops farmed in China’s Shandong Province, a significant agricultural region [33].

2.2. Data and Preprocessing

2.2.1. SIF Data

The Tropospheric Observatory Programmed Observatory (TROPOMI) is carried on the Sentinel-5P (Sentinel-5P) satellite, which has a spatial resolution of 3.5 km × 7 km and a revisit period of 17 days, enabling it to acquire global daily images [34]. Liu Liangyun et al. reanalyzed the TROPOMI SIF products in the 743–758 nm band range [35], synthesizing 0.05° resolution enhanced SIF products (eSIF) with improved signal-to-noise ratios compared to the original TROPOMI SIF data. This was achieved by performing an 8-day back-performing without altering the underlying SIF mechanism, thus preserving the physiological characteristics and data structure of original SIF product. The main genesis of SIF is shown in Equation (1) [4,36,37]. To improve the reliability of the downscaled SIF product, the experiments selected the eSIF dataset of 2018–2020 and chose the vegetation index data related to the chlorophyll content of plants, along with meteorological data representing vegetation growth environments, as the driving variables to carry out the downscaling study [38,39,40].

$$SIF=PAR\times fpar\times {\Phi }_{SIF}$$

(1)

In Equation (1), PAR stands for light and effective radiation, fpar is the proportion of photosynthetically active radiation absorbed by the vegetation, and Φ_SIF is the fluorescence quantum yield at the photosystem level.

2.2.2. MODIS Data

The bidirectional reflectance distribution function adjustment (BRDF) product of the Moderate Resolution Imaging Spectroradiometer (MODIS) MCD43A4 version 6 for 2018–2020 was selected as the main driving variable for constructing the downscaling model. The BRDF includes reflectance data in four bands: red (b1), near-infrared (b2), blue (b3), and green (b4) at a resolution of 500 m. Additionally, the Normalized Difference Vegetation Index (NDVI), Enhanced Vegetation Index (EVI), Soil-Adjusted Vegetation Index (SAVI), Near-Infrared Vegetation Index (NIRv), Optimized Soil-Adjusted Vegetation Index (OSAVI), Leaf Area Index (LAI), and Vegetation Percentage of Cover (FPAR) data of MOD15A2 with 500 m spatial resolution were selected as the important driving variables, which can represent the characteristics of the data structure of the SIF both empirically and mechanically and have been widely used in the inverse of the ecological remote sensing [41]. To ensure the consistency of spatial resolution of the above feature variables, they were resampled to 1 km, and the surface temperature (LST) data of the MOD11A2 product were introduced in this study with a resolution of 1 km to capture the physiological features of SIF more accurately [42].

2.2.3. Validation Data

To verify the accuracy of the downscaled SIF products, we chose two types of SIF data for comparison. The first dataset is the SIF data from the TanSat satellite, generated by Shanshan Du et al. [43], and the TanSIF product produced by Yan Ma et al. [44] based on the satellite’s data in the 758 nm band, with a temporal resolution of 4 days, a spatial resolution of 0.05°, and a coverage period of 2017–2019, which is the dataset chosen for this experiment for the 2018–2019 range. The second dataset is the GOSIF product for the 2018–2020 period, which was derived by Li and Xiao based on the discrete OCO-2 SIF reconstruction with a resolution of 0.05° × 0.05°. The accuracy of the downscaled SIF product was validated by comparing it with both datasets. Furthermore, to assess the reliability of downscaled SIF based on the GPP-SIF mechanism relationship, we used the GOSIF GPP dataset estimated by Li and Xiao et al. from 2018 to 2020 as an additional validation [45].

2.3. Downscaling Method and Model Construction

Currently, empirical statistical modeling is widely used in SIF downscaling studies. The basic principle of this approach is to first select driving variables that are highly correlated with the downscaled variables and then aggregate these variables to match the spatial resolution of the target variables to establish spatial correspondence. Based on the theory of invariant spatial-scale relationships, high-resolution driving variables are used to predict new high-resolution data [40,46,47,48,49,50,51]. To improve the accuracy of downscaling products, this study selected six models for downscaling, including:

(1)

Convolutional Neural Network (CNN): CNN is a fundamental deep learning algorithm, composed of convolutional, pooling, and fully connected layers. The convolutional layer, as the core component, extracts features from the input data. The pooling layer reduces the spatial dimensions of the feature map, which not only decreases the model parameters but also helps mitigate overfitting. Through the weight connections, CNNs can effectively learn and represent complex data features [52].

(2)

Random Forest Regression (RF): Proposed by Leo Breiman, RF is a regression model that integrates multiple decision trees. Each tree is trained on a different subset of samples and features, with predictions made based on the majority vote of the trees. After constructing the model, the out-of-bag (OOB) error is computed. which provides a reliable estimation of the model’s future performance [53,54].

(3)

Extreme Random Tree Regression (ET): ET is similar to the RF model in that it is also a regression model composed of multiple decision trees. However, ET differs from RF in that it uses the entire training set without random sampling and trains all the trees on this full dataset [55]. After dividing the features, a feature value is randomly selected to divide the decision tree [56].

(4)

Adaptive Augmented Regression (ADA): The AdaBoost model iteratively trains a series of weak regressors and combines them into a strong regressor, each weak regressor is trained on the weighted training set, with the weights adjusted based on the errors from the previous iteration, allowing the model to focus on the samples that were previously misclassified [57,58].

(5)

Gradient boosted tree regression (GBDT): GBDT builds a stronger model by combining multiple decision trees, with each decision tree learning the residuals of the previous one. This iterative process helps to improve the accuracy of the overall model [59,60] by minimizing the difference between the original data and the predicted values.

(6)

Stacking integrated learning: Stacking is a powerful machine learning technique that combines multiple base learning models to improve the overall model prediction performance. In this study, we use RF, ET, ADA, and GBDT to construct the Stacking integrated learning model. The accuracy evaluation indices of each model are compared, and the model with better accuracy is selected to construct the SIF downscaling model [61,62].

2.4. Accuracy Assessment Indicators

In this study, the coefficient of determination (R² ), mean absolute error (MAE), and root mean square error (RMSE) were chosen as the accuracy indices to assess the downscaling model [16,63]; in Equations (2)–(4), SIF_t denotes the true value of SIF, SIF_p denotes the predicted value of SIF, SIF denotes the mean value of SIF, and n denotes the total number of samples.

$$MAE={\displaystyle \frac{{\sum }_{i=1}^n\left|{SIF}_i-{SIF}_p\right|}{n}}$$

(2)

$$RMSE=\sqrt{{\displaystyle \frac{{\sum }_{i=1}^n{({SIF}_i-{SIF}_p)}^2}{n}}}$$

(3)

$$R^2=1-{\displaystyle \frac{{\sum }_{i=1}^n{\left({SIF}_i-{SIF}_p\right)}^2}{{\sum }_{i=1}^n{\left({SIF}_i-\overline{SIF}\right)}^2}}$$

(4)

2.5. Downscaling Process

Experimentally, b1, b2, b3, b4, NDVI, EVI, NIRv, OSAVI, SAVI, LAI, FPAR, and LST were selected as the driving variables for constructing the downscaling model. During the model construction process, the resolution of the driving variables was first adjusted to match that of the eSIF data. Then, the corresponding image values were extracted according to the image points of the eSIF, to maintain data accuracy and consistency. The experiments were selected from October 2018 to May 2019 and from October 2019 to May 2020 for modeling, respectively. The downscaling model construction process is illustrated in Figure 2.

(1)

Data pre-preprocessing: The downscaled sample dataset was constructed by analyzing the causes of SIF resampling. The selected MODIS data were resampled to the same spatial resolution as the eSIF data, and the driving variables were extracted based on eSIF.

(2)

Data downscaling: The SIF downscaling model is constructed based on the CNN model, Stacking integrated learning model, Random Forest model, Extreme Random Tree model, Adaptive Enhanced Regression model, and Gradient Boosting Tree model, respectively, and R², RMSE, and MAE are selected as the accuracy evaluation indices, and the model with higher accuracy is selected as the SIF downscaling model. Then, each driving variable is resampled to 1 km, input into downscaling, and based on the assumption that the spatial relationship before and after downscaling is unchanged, the 1 km resolution eSIF inversion is realized.

(3)

Accuracy Verification: In order to determine the accuracy of 1 km SIF, a month-by-month verification was conducted based on a time series. Additionally, the downscaled 1 km SIF products were compared with other SIF satellite products, including TanSIF and GOSIF, to verify their consistency. The reliability of the 1 km SIF was further tested by performing correlation analysis with GPP data, based on the GPP-SIF mechanism. Finally, 1 km SIF was utilized to monitor the growth of wheat in Shandong Province, and its applicability was tested by comparing it with the vegetation index and its trend.

3. Results

3.1. Characterization Variable Correlation

In this paper, reflectance data (b1, b2, b3, and b4); their calculated vegetation indices (e.g., NDVI, EVI, NIRv, OSAVI, and SAVI); LAI; FPAR; and LST were selected as the driving variables for the downscaling model. A correlation analysis was conducted between the driving variables and SIF, with the Pearson correlation coefficient value (R) shown in Figure 3. Among the variables, b2, NDVI, EVI, NIRv, OSAVI, SAVI, and LST exhibited strong correlations with SIF. Based on this, dataset 1 was constructed using these seven variables (b2, NDVI, EVI, NIRv, OSAVI, SAVI, and LST), along with SIF data from October 2018 to May 2019. Additionally, dataset 2 was created using 12 driving variables (b1, b2, b3, b4, NDVI, EVI, NIRv, OSAVI, SAVI, LAI, FPAR, and LST), along with the SIF data.

3.2. Comparison of Downscaling Model Accuracy

The models selected for this study include CNN, Stacking integrated learning model, RF, ET, ADA, and GBDT. These models were optimized, and the parameters were adjusted to the appropriate values for training. The dataset 1 was input into each of the six models for training, and the experiment extracted a total of 138,845 pieces of sample data, of which 70% was used as the training set and 30% as the test set. Figure 4 displays the density plots of the training true values and the predicted eSIF for each model. The plots show that, despite using the same downscaled sample dataset, the models exhibit good precision, indicating high accuracy. Although the accuracy metrics are satisfactory, the fitted lines in the plots are somewhat dispersed from the ideal 1:1 line.

Table 1. Results of the accuracy evaluation metrics computed in each model for dataset 1.

Models	R2	RMSE	MAE
CNN	0.814	0.091	0.057
Stacking	0.886	0.061	0.0372
RF	0.896	0.060	0.0352
GBDT	0.883	0.065	0.042
ET	0.841	0.075	0.0472
ADA	0.76	0.083	0.059

presents the accuracy evaluation metrics for the training results from each model using dataset 1.

The experiment inputs dataset 2 (from October 2018 to May 2019) into six models for training, respectively, and it is found that, after adding the driving variables, each model not only improves the prediction accuracy of the model at a higher level, but according to Figure 5, it can also be seen that, in the density plot of each model, the 1:1 line is more closely fit to the fitted line, and the density of the dots is more centralized. Comparing Figure 4 and Figure 5, the scatter density of RF, Stacking, and GBDT is more centralized, while CNN, ET, and ADA exhibit more dispersed scatter points. Based on these observations, the RF model, which demonstrated the highest accuracy and the best prediction performance, was chosen to construct the downscaling model and invert the 1 km SIF.

Table 2. Results of the accuracy evaluation metrics computed in each model for dataset 2.

Models	R2	RMSE	MAE
CNN	0.882	0.083	0.055
Stacking	0.924	0.055	0.034
RF	0.931	0.052	0.031
GBDT	0.928	0.053	0.033
ET	0.895	0.065	0.041
ADA	0.832	0.088	0.066

shows the accuracy evaluation of the training results of dataset 2 in each model indicator.

For the dataset from October 2019 to May 2020, a downscaled model based on the RF algorithm was trained using the same driving variables. Figure 6 demonstrates the fitted density plot of the real value and the predicted eSIF based on the RF algorithm, and the results show that the model has a great fitting effect. The specific accuracy metrics are R² = 0.926, RMSE = 0.058, and MAE = 0.034. The fitted line is very close to the 1:1 line, which shows a high degree of consistency between the model predictions and the actual observations.

The driving variables with a 1 km resolution were input into the downscaling model to generate the 1 km resolution SIF product (denoted as SIF 1 km) for Shandong Province from October 2018 to May 2019. Figure 7 demonstrates the eSIF maps synthesized for each month from October 2018 to May 2019, while Figure 8 demonstrates the downscaled SIF for the same period, with a resolution of 1 km. The comparison reveals that the downscaled SIF has a finer spatial resolution. From October 2018 to November 2018, wheat was mainly planted in the northwest and southwest of Shandong Province. During this period, as the temperatures gradually declined and the vegetation cover decreased, the SIF values showed a decreasing trend. As Shandong Province is located in Northern China, the temperature is low in winter, and the winter wheat after sowing in the fall is affected by the temperature and grows slowly. With temperatures increased, the winter wheat entered the greening period, and the SIF values rose significantly. By May, when the average temperatures rose to the point where vegetation was in the optimal growth conditions and overall vegetation growth was vigorous. Geographically, the central part of Shandong Province is mountainous, and the northeastern coastal region has complex topography. Higher altitudes have lower SIF values, because long-wave radiation becomes a heat source instead of direct sunlight and vegetation is less stimulated by light. As the terrain rises, the air is thinner, the temperature rises slowly, and the vegetation has a slower growth period relative to the plains. In the northeastern coastal region, the cooler temperatures and slower vegetation recovering are influenced by the nearby ocean. The experimental results show that the downscaled 1 km resolution SIF data more accurately reflects local fine-scale changes compared to the eSIF.

From October 2019 to May 2020, the same methodology was applied to generate a 1 km resolution SIF product through monthly-scale SIF inversion for Shandong Province. Figure 9 and Figure 10 present the inversion results, showing both the synthesized eSIF and the downscaled 1 km SIF. Compared to the previous year, the SIF trends are the same. By comparing the SIF data products before and after downscaling, the distribution of data is unchanged, and the downscaling model constructed based on the Random Forest (RF) algorithm can accurately invert the 1 km SIF, which can more accurately observe the changes in vegetation cover and provide strong support for studying agriculture-related information through SIF changes.

3.3. Timing Verification 1 km SIF

In order to verify the applicability of the downscaled 1 km SIF, a correlation analysis with eSIF was performed to verify the stability of the downscaled model. The downscaled 1 km SIF products were resampled to match the 0.05° spatial resolution of the original eSIF, ensuring consistency in the spatial and temporal data. Figure 11 and Figure 12 show the scatter density plots of SIF correlations before and after downscaling, for the periods from October 2018 to May 2019 and October 2019 to May 2020. These plots demonstrate a strong correlation (R² > 0.8) between the downscaled SIF products and the original eSIF data, indicating that the downscaling model is stable and its predictions are highly accurate. At the same time, this also verifies the validity of the theoretical assumption based on the invariance of the scaling relationship. Overall, the downscaling method effectively preserves the characteristics of the original data, providing detailed SIF information at a higher spatial resolution, which provides a reliable database for further vegetation ecology research and agricultural production monitoring.

3.4. Verification of the Accuracy and Reliability of the 1 km SIF

The 2018–2019 TanSIF and 2018–2020 GOSIF time scales were aggregated, and to maintain the spatial and temporal consistency of the SIF data, the downscaled 1 km SIF products were resampled to 0.05° for correlation validation analysis, to ensure the accuracy of the 1 km SIF. Figure 13 presents the density plots of the validated correlation, where (a) and (b) show the fit between the downscaled SIF with TanSIF/GOSIF for the first year, respectively, while (c) shows the fit between the downscaled SIF and GOSIF for the second year. The upper part and the right sections of the plot represent the data trends of the downscaled 1 km SIF and other SIF satellite products, respectively. The results show that the R² values exceed 0.8, with RMSE and MAE both below 0.1. Figure 13 shows that the downscaled 1 km SIF products are consistent with TanSIF and GOSIF products in data trends, indicating that the downscaled SIF products have good consistency with TanSIF and GOSIF products.

GPP (ecosystem primary productivity) is one of the key parameters of terrestrial ecosystems, which directly reflects the efficiency of vegetation in utilizing light energy. Fu et al. highlighted that total primary productivity is a major component of carbon flux in the Earth’s ecosystems [64]. Koehler et al. were the first to explore the coupling between satellite-observed SIF data and the GPP, and subsequent studies have shown a good linear relationship with GPP from SIF data obtained from ground-based observations, tower-based, and satellite platforms [65]. GPP can be expressed as the product of photosynthetically active radiation absorbed by the canopy (APAR) and the true light energy use efficiency (LUEp) by the light energy utilization rate model, namely

$$GPP=APAR\times LU{E}_P=PAR\times fPAR\times LU{E}_{\mathrm{pmax}}\times f(stress)$$

(5)

In addition, LUE_pmax is the maximum light energy utilization, and f(stress) is the environmental limiting factor. The linear relationship between SIF and GPP can be obtained by means of the linear relationship between SIF and GPP in Equation (6) by association with Equation (5). The equation can be expressed as

$$GPP=SIF\times {\displaystyle \frac{LU{E}_p}{{\Phi }_{SIF}}}\times {\displaystyle \frac{1}{\epsilon }}=SIF\times {\displaystyle \frac{LU{E}_{\mathrm{p}}}{SI{F}_{yield}}}$$

(6)

A linear relationship between SIF and GPP can be seen from Equation (6). To verify the reliability of the downscaled SIF, the coupling relationship of SIF-GPP is examined. For this purpose, the GOSIF GPP dataset, estimated based on GOSIF products, was selected for this experiment. This dataset utilizes eight different forms of SIF-GPP relationships to estimate the GPP at 91 eddy covariance flux sites worldwide, and the ensemble-averaged eight-day GPP shows a strong correlation with the flux tower GPP. The dataset has a temporal resolution of 8 days and a spatial resolution of 0.05°, providing a continuous dataset. In this experiment, two years of downscaled 1 km SIF products were synthesized and resampled to match the spatial resolution of the GOSIF GPP. The GPP data were also synthesized into one year for consistency. The results obtained by correlation analysis are shown in Figure 14. Figure 14a,b show the fitted plots of SIF versus GPP for the first and second years, respectively. It can be observed that the two-year SIF and GPP show a strong correlation and the trend distribution above and to the right of the scatter plot is the same, which verifies the reliability of the downscaled 1 km SIF product.

3.5. Growth Trend Analysis of Wheat by 1 km SIF

Previous applications of satellite SIF products in crop growth monitoring and yield estimation have been somewhat more limited. To further assess the applicability of downscaled SIF products in crop growth analysis, spatial SIF values were extracted from wheat-growing areas in Shandong Province. Using the GLASS LAI remote sensing products and long-term observation of climatic records from agro-meteorological stations, we defined the inflection points and peaks of the key climatic periods of each crop. These key nodes of the LAI characteristic curves were identified by the inflection point and threshold methods, thus accurately capturing the growth dynamics of wheat in a year. In addition, the SIF, LAI, NDVI, and EVI data at the same spatial and temporal scales were selected to analyze the response of their changes to the growth status of wheat over eight days. The results showed that the trends of the SIF values closely followed those of the LAI, NDVI, and EVI, with SIF changes exhibiting greater stability than the other vegetation indices (Figure 15). Specifically, October marked the sowing period of wheat. From November to January, during winter, wheat’s vegetation index fluctuated significantly due to temperature effects, and the SIF values were also influenced by temperature and snow cover, which reduced sunlight absorption, reaching their lowest point in January. As temperatures increased from February to April, the vegetation index gradually rose, and in April, wheat growth reached its peak, coinciding with the highest SIF values. By May, as the wheat was harvested, both the SIF and vegetation index values began to decrease gradually. The wheat growth cycle showed similar trends in both years.

Furthermore, we conducted a correlation analysis between downscaled 1 km SIF and three vegetation indices. Figure 16 shows that the correlations in both years were high (R² > 0.8), with the strongest correlation observed with the LAI (R² > 0.95). This indicates that downscaled SIF not only possesses similar structural characteristics to traditional vegetation indices but also can outperform them to a certain extent due to its mechanism of spontaneous luminescence based on sunlight stimulation. Therefore, the accuracy of downscaled SIF and its potential and effectiveness in this application area was verified through crop growth.

4. Discussion

In previous studies of SIF-scale reconstruction, most have focused on using single models, with limited comparisons among multiple models and a lack of analysis on the driving variables. Furthermore, many studies did not verify the accuracy of reconstructed SIF in terms of the time series and the GPP-SIF mechanism, making it difficult to apply the results in smaller provincial areas [41,66]. To address these gaps, we chose various driving variables, such as the 1 km scale MODIS reflectance bands (b1, b2, b3, and b4); NDVI; EVI; NIRv; OSAVI; SAVI; LAI; FPAR; and LST for downscaling modeling. Vegetation indices calculated based on reflectance bands reflect the structural properties of SIF, which are associated with the chlorophyll content of vegetation, explaining their relationship with SIF mechanistically, and different vegetation indices provide different vegetation information, reflecting the various data structures of vegetation; among which, the NDVI and EVI can accurately monitor the changes in vegetation [67]. NIRv, as a new development of vegetation index, can reduce the noise effect caused by surface vegetation pollution and retrieve more details of vegetation dynamics compared to the NDVI [68]. The SAVI and OSAVI are particularly useful for monitoring sparsely vegetated areas, as they mitigate the influence of the soil background [69,70,71]. Additionally, LST is effective in capturing the growing environment information of SIF, better reflecting the physiological characteristics of vegetation. While the LAI and FPAR show a weaker correlation with SIF, they remain ecologically significant, particularly in agroforestry systems, as they reflect the plant canopy structure and vitality. These indices can improve the characterization of vegetation structure in the model training process [37,72].

By incorporating additional variables such as b1, b3, b4, LAI, and FPAR, the inclusion of multiple feature variables significantly improved the model’s prediction accuracy. Meanwhile, the accuracy of each downscaled SIF model was significantly improved, to more accurately invert the 1 km SIF, and a month-by-month comparison between the inverted 1 km SIF and the original eSIF image after inputting the high-resolution dataset into the model revealed that the SIF data structure did not change and there was a good correlation between the two. These findings align with similar studies that have demonstrated improvements in downscaling accuracy by adding more comprehensive vegetation data [73,74,75]. However, as shown in previous research, the inversion of 1 km SIF is still incomplete due to the lack of some driving variables, which can be overcome by adding more data or interpolating missing values [76,77]. The comparison of the monthly 1 km SIF with the original eSIF dataset indicated that, while the data structures remained consistent, a complete 1 km SIF reconstruction can be hindered by insufficient driving variables. This observation is consistent with the findings of Jacobson et al. (2023), who noted that missing variables during reconstruction can lead to inaccuracies in SIF products at a finer scale [78].

To verify the accuracy of the 1 km SIF, traditional validation methods using other SIF products were employed, as well as a novel approach based on the GOSIF-GPP mechanism. Although temporal and spatial resolution discrepancies between products caused some errors in data synthesis and aggregation, the overall correlation between the two was strong [67,79]. However, due to the lower resolution of the GPP data, it failed to capture finer details of the vegetation dynamics, as shown in the study by Wang et al. (2021) [80]. Higher-resolution GPP data would provide more accurate validation, as confirmed by recent work in the field [81].

5. Conclusions

Due to the low resolution of SIF data, its application in small-scale agricultural monitoring is challenging. This study applied a downscaling approach, utilizing the Random Forest (RF) model with MODIS products (b1, b2, b3, b4, NDVI, EVI, NIRv, OSAVI, SAVI, LST, LAI, and FPAR) as the driving variables. The model achieved high accuracy for the period from October 2018 to May 2019, with R² = 0.931, RMSE = 0.052 mW/m²/nm/sr, and MAE = 0.031 mW/m²/nm/sr. From October 2019 to May 2020, the model’s accuracy remained high, with R² = 0.926, RMSE = 0.058 mW/m²/nm/sr, and MAE = 0.034 mW/m²/nm/sr. The 0.05° resolution eSIF products were downscaled to 1 km resolution, and validation using time series data confirmed that the downscaled SIF maintained a high accuracy (R² > 0.8). Comparison with TanSIF and GOSIF products showed strong correlations (R² > 0.8), and the data distribution trend was consistent with GOSIF-GPP, demonstrating a strong SIF-GPP linear relationship.

This study demonstrates the successful downscaling of SIF products for improved crop growth monitoring in precision agriculture. It provides valuable data and methodologies for future research in remote sensing applications. Future studies may explore further enhancements in model accuracy and the integration of higher-resolution datasets to extend the applicability of downscaled SIF in ecological and agricultural monitoring.