## 4 Conclusions, Limitations, and Future Research  
Satellite observations of SIF provide a viable proxy for GPP at the global scale [3; 6]. However, these observational data are noisy and spatially incomplete so scientists need denoised and gap-filled Level 3 data products to better characterize the relationship between SIF and GPP. Here, we use a geostatistical methodology called cokriging to produce coSIF, a Level 3 SIF data product at a 0.05-degree, monthly resolution with uncertainty quantification. Our highly parallelizable approach leverages OCO-2 observations of both SIF and XCO\({}_{2}\), which we show helps to reduce the uncertainty of the Level 3 data product in months of high GPP. It achieves this by modeling the cross-dependence between SIF and XCO\({}_{2}\) one month later and, moreover, it accounts for non-stationary measurement errors in the observations. Due to its fine resolution, coSIF can be aggregated to coarser resolutions, along with coherent uncertainty quantification, which enables statistical comparison with other SIF products at varying resolutions. \begin{table}
\begin{tabular}{l l c c c c} \hline \hline
**Block** & **Method** & **BIAS** & **RASPE** & **INT** & **DSS** & **MDSS** \\ \hline Corn Belt & Cokriging & \(-0.06\) & \(\mathbf{0.58}\) & \(\mathbf{2.85}\) & \(-\mathbf{0.13}\) & \(-\mathbf{311.94}\) \\  & Kriging & \(-\mathbf{0.01}\) & \(0.59\) & \(2.93\) & \(-0.12\) & \(-308.63\) \\  & Trend surface & \(0.03\) & \(0.62\) & \(3.03\) & \(0.03\) & \(21.32\) \\ \hline Cropland & Cokriging & \(\mathbf{0.01}\) & \(\mathbf{0.56}\) & \(2.95\) & \(-\mathbf{0.22}\) & \(-\mathbf{356.18}\) \\  & Kriging & \(0.03\) & \(\mathbf{0.56}\) & \(\mathbf{2.93}\) & \(-0.20\) & \(-355.44\) \\  & Trend surface & \(0.07\) & \(0.60\) & \(3.05\) & \(-0.08\) & \(-57.77\) \\ \hline \hline \end{tabular}
\end{table} 
Table 1: Validation metrics for the Corn Belt and Cropland \(5^{\circ}\times 5^{\circ}\) validation blocks in July 2021. Bold typeface indicates the best performing method according to each metric. Locations of the validation blocks are shown in Figure 3.  
We limit our analysis of coSIF to a large region over North America, but by splitting global land into contiguous regions where geostatistical stationarity assumptions and the use of chordal distances are appropriate, our approach can yield a global data product, albeit with different estimates of spatial dependence parameters in different regions. Though we do not formally incorporate temporal dependence in the bivariate covariance model, this is the subject of future research. Incorporating dependence with additional environmental processes in a higher-dimensional multivariate covariance model is another promising extension. The cokriging method could also be used to produce improved predictions of \(\text{XCO}_{2}\); one would simply swap the roles of \(Y_{1}(\cdot)\) and \(Y_{2}(\cdot)\) to produce a de-noised, gap-filled data product for \(\text{XCO}_{2}\) with SIF as the secondary variable, and with uncertainty quantification of \(\text{XCO}_{2}\). The following supporting information can be downloaded at: [[https://www.mdpi.com/article/10.3390/rs15164038/s1](https://www.mdpi.com/article/10.3390/rs15164038/s1)]([https://www.mdpi.com/article/10.3390/rs15164038/s1](https://www.mdpi.com/article/10.3390/rs15164038/s1)), Dataset S1: Fitted parameters of the bivariate spatial models for February, April, July, and October 2021; Table S1: Average ratio of the coSIF RMSPE to the measurement error standard deviation for February, April, July, and October 2021; Figures S1-S3: The coSIF data product (predictions and corresponding RMSPEs) for February, April, and October 2021; Tables S2-S4: Validation metrics for cokriging, kriging, and trend-surface-only prediction in February, April, and October 2021.