## 2 Materials and Methods
### Methodology  
#### 2.3.1 Data Preprocessing  
To obtain accurate albedo, emissivity, land surface temperature, and latent heat flux data, we used the quality control data of MODIS to extract the data first. In other words, we extracted the pixel data of the value that represented good quality as a quality control mask and filtered the origin data to obtain surface parameter data with good quality. The wide-band emissivity was needed for the calculation of the surface longwave radiation. However, the MODIS products only provided some narrow-band emissivity \begin{table}
\begin{tabular}{c c c c c} \hline \hline
**Data** & **Data Source** & **MAE 1** & \begin{tabular}{c} **Spatial** \\ **Resolution** \\ \end{tabular} & \
\begin{tabular}{c} **Temporal** \\ **Resolution** \\ \end{tabular} \\ \hline Land use &
\begin{tabular}{c} Land-use database \\ from CAS \\ \end{tabular} & \(\geq\)90\% 2 & 1 km & - \\ Albedo & MCD43B3 & \(-\)0.008 & 1 km & 8 d \\ Emissivity & MOD11A2 & 0.001–0.005 & 1 km & 8 d \\ LST & MOD11A2 & \(\leq\)1 \({}^{\circ}\)C & 1 km & 8 d \\ Latent heat flux & MOD16A2 & 0.31 to 0.33 mm.d\({}^{-1}\) & 1 km & 8 d \\ Snow & MOD10A2 & \(\geq\)93\% 3 & 500 m & 8 d \\ Radiation data & ERA-Interim & - & 0.75\({}^{\circ}\) & Monthly \\ \hline \hline \end{tabular}  
* MAE here means ‘mean absolute error’. 2 Data here means the accuracy of land use data. 3 Data here means the accuracy of snow cover data.  
\end{table} 
Table 2: Accuracy range and resolution of the main dataset. information. Hence, we used a nonlinear fitting algorithm of high precision to fit the wide-band emissivity from the 2 narrow-band datasets offered by the MODIS products [36]. The fitting algorithm was as follows:  
\[\varepsilon_{b}=0.273+1.778\varepsilon_{31}-1.807\varepsilon_{31}\varepsilon_{32 }-1.037\varepsilon_{32}+1.774\varepsilon_{32}^{2} \tag{1}\]  
where \(\varepsilon_{b}\) is the wide-band emissivity and \(\varepsilon_{31}\) and \(\varepsilon_{31}\) are the narrow-band emissivity of band 31 and band 32 of MOD11A2, respectively. The MODIS products only provide \"white sky\" and \"black sky\" albedo data, but we needed to use the actual surface albedo in this study. Therefore, according to sky scattering factors, we calculated the actual surface albedo based on the \"white sky\" and \"black sky\" albedo data [37]. The actual albedo can be calculated as follows:  
\[\alpha_{blue}=S\alpha_{white-sky}+(1-S)\alpha_{black-sky}. \tag{2}\]  
where \(\alpha_{blue}\) is the actual albedo (i.e., the blue sky albedo), \(\alpha_{white-sky}\) is the white sky albedo, \(\alpha_{black-sky}\) is the black sky albedo, and \(S\) is the sky scattering factor. \(S\) can be calculated as follows:  
\[S=a\cdot(cos(SZA))^{b} \tag{3}\]  
where \(a\) and \(b\) are the regression coefficients, which are usually equal to 0.1 and \(-0.8\), respectively, in the mid-latitude region; \(SZA\) is the solar zenith angle, which can be extracted with the MODIS albedo product. Finally, because of the significant impact of snow on the surface parameters, we improved the surface parameter data with the snow data. To ensure consistent spatial resolution, we used the nearest neighbor sampling method to resample the maximum extent of snow cover data to a 1 km spatial resolution. Then, we extracted the pixel data of value 25, which means the land without snow cover as a mask, and filtered the surface parameter data to obtain snow-free surface parameter data. 2.3.2 Calculation of Radiative Forcing and the Heat Flux to the Atmosphere Due to Cropland Expansion  
Surface radiative forcing is commonly used to describe the net radiation energy at the Earth's surface after climate factor changes over a certain time scale (e.g., inter-annual scale); it is expressed in units of W/m\({}^{2}\). Because the technique avoids the complex feedback processes of the Earth's gas system, surface radiative forcing can be used to make relatively simple comparisons of the contributions of several climatic factors to climate change.