## 3 Results
### Random Forests  
In total, four random forests were tuned: for step 1, the RF regression with MODIS FSC (used in A1 and A2) and the regression with MODIS NDSI (used in A3), and for step 2, the RF classification with SCE S2 (in A1) and the RF regression with NDSI S2 (in A2 and A3). The tuning results of the other RFs are given in the Supplementary Materials. The two random forests in the workflow were tuned to extensively explore the predictive behavior of RF. The results reveal that this step is superfluous, and no significant improvement in performance is obtained when varying the hyperparameters considered, as highlighted by other studies [92]. The RF of the first step, which was used to complete the time series of the MODIS snow product, was fitted with a training sample of 25% of each MODIS scene available for the study area in the time frame under consideration, equal to 1.9 \(\times\) 10\({}^{4}\) pixels. No other training sizes were tested. First, the optimal number of trees to grow, between 100 and 1000, was investigated as a function of the explained variance (R\({}^{2}\) score) and the computational time. Based on the RF tuning results, we used an RF with 500 trees, _mtry_ equal to 7, and a minimum _node size_ of 3. In the second phase, tuning tests of the _ntree_ hyperparameter were carried out with 100, 200, and 300 trees, the results of which showed that, with larger forests, the calculation time considerably increases with almost no improvement. An RF of 100 trees was then grown, with _mtry_ equal to 3 and a minimum _node size_ of 2. The RF tuning results are showed in the Supplementary Materials. In terms of computational resources (Table 7), A1 is the most efficient, taking an average of 1.5 min per scene, while A2 and A3 reconstruct each scene in an average of 2 up to 3 min. Furthermore, the random forest training in the second step is faster in A1 compared to the other two approaches, considering a training dataset of 2.3 \(\times\) 10\({}^{7}\) observations and eight variables; above all, when reconstructing daily images, A1 requires 20 GB RAM to process, while A2 and A3 occupy more than 110 GB, which is not available on most machines. \begin{table}
\begin{tabular}{c c c c} \hline \hline  & **A1 vs. A2** & **A1 vs. A3** & **A2 vs. A3** \\ \hline \(\chi^{2}\) & 604.4 & 464.5 & 7.6 \\ \(p\)-value & \textless{}2.2 \(\times\) 10\({}^{-16}\) & \textless{}2.2 \(\times\) 10\({}^{-16}\) & 0.0058 \\ \hline \hline \end{tabular}
\end{table} 
Table 6: McNemar test between the three approaches’ performances against weather stations. \begin{table}
\begin{tabular}{c c c c} \hline \hline  & **A1** & **A2** & **A3** \\ \hline RF training (step 2) time (min) & 43 & 56 & 70 \\ Prediction time per image (min) & 1.5 & 2 & 3 \\ RF prediction (step 2) resources (GB of RAM) & 20 & \textgreater{}110 & \textgreater{}110 \\ \hline \hline \end{tabular}
\end{table} 
Table 7: Computational and time resources required by the processing. In the _ranger_ package, is it possible to compute the 'impurity_corrected' measure Gini index for classification and the variance reduction of the responses for regression [93]. The Gini index is a splitting measure that determines the information gained at each split node. Specifically, at each split, the Gini, or impurity index, measures the variable ability to separate observations into two distinct classes, lowering entropy. The average of all Gini scores obtained in the forest for variable yields the Gini variable importance measure [94]. The variance reduction is a method of calculating the importance of features on out-of-bag (OOB) samples based on the average decrease in accuracy [95]. In both cases, the higher the score, the better the predictive ability of the variable. Figure 10 shows the variable importance of the RFs used in the first (a) and second steps (b). The most important features in the RF of step 1 are the day of the year, elevation, and year, followed by latitude (Y) and longitude (X), and, lastly, aspect and slope.