This study investigated the effects of spatial heterogeneity in atmospheric forcing, land use and land cover, soil properties, and topography on the modeling of evapotranspiration, runoff, and their components (i.e., canopy evaporation, ground evaporation, transpiration, surface and subsurface runoff).  The design of the numerical experiments is reasonable, the methods are innovative, and the results are insightful. With the above considerations, I suggests a publication with several clarifications.

Detailed comments

- L121--L124: Better put them in the Methodology section.

- L135, 2013a Li et al., 2013: typos?

- L136--L140: What is the purpose of these statements? They duplicate the discussion. On the other hand, more descriptions of ELM on the used parameterizations are desirable.

- L166--L167: How did you resample soil properties? How soil texture and organic matter are upscaled, and how soil color is downscaled?

- L199--L209: The calculation of Sobol's sensitivity index is still a bit confused to me. Assuming that X is the 30-year monthly value, SENSITIVITY(X) is the Sobol's sensitivity index of X, and 30-YEAR-SEASONAL-AVERAGE(X) is the seasonal average of X, did you calculate the index as SENSITIVITY(30-year-seasonal-average(X)) or 30-year-seasonal-average(SENSITIVITY(X))?

- Table 4: Did you first propose this approach?

- L260: It would be to compare the spatial variability with the temporal variability or the mean value of ET/P.

- L305--L328 and Figure 5: I did not get the objective of these contents. Since ATM and LULC are dominant, it is quite natural that their contributions to the spatial variability are complementary. Deletion seems fine and would make the paper concise.

- L392-L395: Do the interplay between the spatial variabilities in ATM, LULC, SOIL, and TOPO increase or decrease the overall variability?

- Section 3.5 and Figure 8: It is interesting to see the comparison between ERA5-Land and the 0.125-degree ELM simulations. Since the ERA5-Land atmospheric forcing is interpolated from the ~31-km ERA5 data with the consideration of elevation dependency (doi:10.5194/essd-13-4349-2021) rather than upscaled from the finer-scale observations, the spatial variability from ERA5-Land would be smaller than that from ELM. In Figure 8, it would be better to show the difference rather than absolute difference to check this.

This manuscript examines the relative importance of four heterogeneity sources on the spatially aggregated variability of some modeled components of the water and energy partitioning over CONUS. The study is appropriately motivated and, overall, well written. The presented methods are adequate and innovative, and the results reasonable. I suggest publication after some minor revisions and clarifications.

L136 – L140: Since this study did not use the described model developments, I recommend removing this paragraph. Instead, details on the methods used to characterize the land surface heterogeneity within ELM for this study would be appreciated.

L166 – L167: Please clarify the methods used to resample the original datasets to 0.125°.

L251 – L252: Please clarify the methods used to resample.

Section 2.4: Have you also considered analyzing the first-order Sobol sensitivity index? The first-order index would measure just the direct effect of each heterogeneity source on the variance of the model. Compared to the total sensitivity, the first order index would explain the importance of the interactions between heterogeneity sources (e.g., topography and atmospheric forcing; soils and topography).

L361 – L370: A more comprehensive analysis of the factors explaining the seasonal variation of the importance of heterogeneity sources would be appreciated here. For instance, what drives the increased relevance of topography in zones of the East in Summer and Fall.

Table 5: I recommend moving this table to Supplementary material and including Figure S5 in the main text (primarily, since L361 to L370 mainly focus on analyzing maps).