Snow Ensemble Uncertainty Project (SEUP): quantiﬁcation of snow water equivalent uncertainty across North America via ensemble land surface modeling

obtained from the operational, global analysis products. The model assimilates early in-situ surface observations as well as modern high-resolution satellite datasets. In this study, the operational real-time data from the ECMWF-Integrated Forecast System (IFS) are used; the meteorological fields are provided on a 0.25 ° grid (roughly 25 km) at a 3-hour interval, and generated by 95 assimilating available atmospheric observations every 12 hours into a forecast model with surface meteorological fields (e.g., precipitation and radiation), which are diagnosed from the model output (Dee Flemming


Introduction
The following supplemental information provides additional detail on the SEUP ensemble, which is comprised of 12 ensemble members, created by the combination of four different land surface models (LSMs, Text S1) and three different forcing datasets 5 (Text S2).
Text S1 Land Surface Models (LSMs) S1.1 Noah version 2.7.1 (Noah2.7.1) Noah2.7.1 (Chen et al., 1996;Ek et al., 2003;Koren et al., 1999) is a LSM that has evolved through community efforts to simulate land surface temperature, snow depth, SWE, canopy water content, surface energy fluxes, water balance, soil 10 temperature, and soil moisture based on the Oregon State University LSM (Mahrt and Pan, 1984). The Noah model has been widely used in operational and research applications. For example, the Noah model has been adopted as the regional and global operational numerical weather prediction model at the NCEP and the 557 th Weather Wing of the U.S. Air Force (formerly the Air Force Weather Agency). It also has been extensively evaluated in both the off-line mode (e.g., Chen et al., 1996;Mitchell et al., 2004) and the coupled mode (e.g., Chen et al., 1997;Ek et al., 2003). More detailed descriptions of the Noah physics 15 and development are presented by Ek et al., (2003) and Koren et al. (1999). This study employs the Noah version 2.7.1, which has a single canopy, and a single snow layer. The one-layer snow physical model simulates SWE as the residual of snowfall minus the sum of snowmelt and sublimation. Snowfall occurs whenever there is nonzero precipitation and the near surface air temperature is less than 0℃. Snow melting occurs when the temperature for canopy, ground, and snow is greater than the freezing point of water when SWE is greater than zero. Refreezing occurs when the layer temperature is less than the freezing 20 point, and the liquid water content is greater than zero. Snow albedo representation uses a parameter --maximum snow albedo --which is the upper bound of the fraction of reflected incoming shortwave radiation for deep snowpacks. Within Noah2.7.1, the maximum snow albedo can vary spatially but is temporally constant at given locations; these values were derived from Defense Meteorological Satellite Program (DMSP) imagery (Robinson and Kukla, 1985).

S1.2 Noah-MP version 3.6 (Noah-MP3.6) 25
The Noah-MP3.6 LSM (Niu et al., 2011) is built upon the legacy of the Noah model, but with new and multiple options for selected processes. Most notably, Noah-MP3.6 performs separate energy budget calculations for the bare and vegetated portions of the grid cell; shortwave radiation is computed over the entire grid considering canopy gap probabilities, but longwave radiation, latent heat, sensible heat, and ground heat fluxes are computed separately over the vegetated area and bare ground. Noah-MP3.6 includes a four-layer soil model and up to three layers for snow, depending on snow depth. The revised snow scheme allows for snow compaction from the weight of overlying layers and melt metamorphism. Snow cover fraction is computed as a function of snow depth, ground roughness length, and snow density (Niu & Yang, 2007). Noah-MP3.6 has two options for calculating snow surface albedo: the Biosphere-Atmosphere Transfer Scheme (BATS; (Yang et al., 1997) and a Canadian Land Surface Scheme (CLASS; (Verseghy, 1991). CLASS is used here, which computes snow albedo based on the albedo of fresh snow and the age of the snow. For partitioning between rain and snow, Noah-MP3.6 has three options. Two 35 options are simple comparisons of air temperature versus the freezing temperature (Tair < Tfrz + 2.2 K and Tair < Tfrz). The third option uses a function based on air temperature where the fraction of snow decreases between 1° and 2℃ (Jordan, 1991), which is used in this study. Niu et al. (2011) performed offline comparisons for Noah and Noah-MP3.6 for two separate field sites and showed that Noah-MP3.6 improves the simulation of SWE, snow depth, density, and snow skin temperature. Compared with Noah, Noah-MP3.6, allows for the retention of meltwater within the snowpack and the refreezing of melt water, which 40 improves the simulation of snow evolution during the melting season.

S1.3 Catchment LSM version Fortuna 2.5 (CLSM-F2.5)
The Catchment version used in this study is the Fortuna 2.5 version used in the Modern-Era Retrospective Analysis for Research and Applications (MERRA) with improved land surface variables (MERRA-Land; Reichle et al., 2011). The Catchment LSM (CLSM) shares the model development legacy with the Mosaic LSM (Koster and Suarez, 1996), which has 45 been enhanced by adding a shallow groundwater module and by using a catchment-based approach (Koster et al., 2000). The basic computational unit of CLSM-F2.5 is the hydrological catchment (a.k.a., watershed) rather than a uniform latitudelongitude grid cell. It uses hydrological catchments as tiles with boundaries defined by topography. Each catchment is separated into three distinct and dynamically varying subareas: a saturated region, an unsaturated region, and a wilting region, which allows the sub-grid scale simulation of differing energy and water budget partitioning based on the moisture state. For ease of 50 model intercomparison, however, the CLSM-F2.5 is used on a regular latitude-longitude grid in this study. CLSM-F2.5 includes a three-layer snowpack model that incorporates snow physics including densification, snowmelt, refreeze, and snow insulating properties (Lynch-Stieglitz, 1994). The prognostic state variables for each of the three layers include SWE, snow cold content, and snow depth. The model uses an important parameter, SWE_min, which describes the minimum SWE that must be present per unit surface area before the model considers the surface to be snow covered. In the 55 SEUP runs, SWE_min is set to 0.013 / 2 . When there is surface melt or when rain falls on existing snow, water can percolate into the lower snow layers, where it may refreeze. Each layer of newly-fallen snow is initially given a density of 150 kg-m -3 , which increases as the snow ages. Snow surface albedo is treated with different reflectance values for the visible (VIS) and near-infrared (NIR) radiation bands, with reductions in albedo implemented by both vegetation masking and fractional snow cover (Stieglitz et al., 2001).

S1.4 Joint UK Land Environment Simulator (JULES)
The Joint UK Land Environment Simulator (JULES; Best et al., 2011;Clark et al., 2011) is based on the Met Office Surface Exchange Scheme (MOSES; Cox et al., 1999;Essery et al., 2003). This land surface model has been coupled to an atmospheric global circulation model, as done at the Hadley Center for Climate Prediction and Research (Blyth et al., 2006), or run in standalone mode driven by forcing data. In this study, the multi-layer snow model is set at a maximum of 3-layers. The actual 65 number of layers depends on the snow depth, and not all the layers necessarily exist at any one time. Snow first accumulates until it reaches twice its prescribed thickness in the lowest layer, then splits into 2 layers with the upper part staying fixed in thickness and continuing to accumulate in the lower layer, then splits again. For each layer in the snowpack, snow thickness, density, temperature, grain size, and ice/liquid content are dynamically calculated and these parameters control the heat capacity and conductance. The density of each layer changes with time, liquid water forms, percolates down, and freezes again, 70 sublimation happens, and the albedo is changed via snow grain size. JULES has two types of snow albedo schemes: the diagnostic albedo scheme and the prognostic albedo scheme. In this study, the diagnostic snow albedo parameterization option is selected (consistent with the operational configuration at the UK Met Office, Parallel Suite 41 (PS41)) which only uses snow depth and cannot represent the impacts of snow aging on the surface albedo. A detailed explanation can be found in Best et al. (2011). 75

S2.1 Modern-Era Retrospective Analysis for Research and Applications, version 2 (MERRA2)
The Modern-Era Retrospective Analysis for Research and Applications, version 2 (MERRA2; (Gelaro et al., 2017;Molod et al., 2015) is the long-term global reanalysis produced by NASA's Global Modeling and Assimilation Office (GMAO) that assimilates in situ and satellite observations into a global atmospheric model. It spans the satellite observing era (from 1980 to 80 the present) with a regularly-gridded, homogeneous record of the global atmosphere, and improved land surface representation (i.e., the use of observation-based precipitation to force the land surface). MERRA2 has a native spatial resolution of 0.5° latitude by 0.625° longitude (roughly 50 km) and hourly temporal resolution.

S2.3 European Centre for Medium-Range Weather Forecasts (ECMWF)
The European Centre for Medium-Range Weather Forecasts (ECMWF; Molteni et al., 1996) data is obtained from the operational, global analysis products. The model assimilates early in-situ surface observations as well as modern highresolution satellite datasets. In this study, the operational real-time data from the ECMWF-Integrated Forecast System (IFS) are used; the meteorological fields are provided on a 0.25° grid (roughly 25 km) at a 3-hour interval, and generated by 95 assimilating available atmospheric observations every 12 hours into a forecast model with surface meteorological fields (e.g., precipitation and radiation), which are diagnosed from the model output (Dee et al., 2011;Flemming et al., 2015).

Rockies, USA
Cascades