Light-absorbing Particles in Snow and Ice: Measurement and Modeling of Climatic and Hydrological impact

Light absorbing particles (LAPs, e.g., black carbon, brown carbon, and dust) influence water and energy budgets of the atmosphere and snowpack in multiple ways (e.g., Twomey et al., 1984; Albrecht, 1989; Hansen et al., 1997; Ramanathan et al., 2001; Lau and Kim, 2006; Qian et al., 2009, 2011; Bond et al., 2013). In addition to their effects associated with atmospheric heating by absorption of solar radiation and interactions with clouds, LAP in snow reduces snow albedo and absorbed more solar radiation (a.k.a., snow darkening), further accelerating the snow aging process and the speed of snowpack melt (Warren and Wiscombe, 1980; Hansen and Nazarenko, 2004; Ming et al., 2009; Xu et al., 2009a). LAP in snow was identified as one of major forcing agents affecting climate change in the fourth assessment report (AR4) of Intergovernmental Panel on Climate Change (IPCC) (IPCC, 2007). However, the uncertainty level in quantifying this effect remains very high even in the most recent fifth assessment report (AR5) by IPCC (IPCC, 2013).

(1) Light absorbing particles

Black carbon (BC, also referred to as elemental carbon, EC) and dust are two primary types of LAP that have been focused on. In this review article, we mainly focus on the BC, with the other LAPs briefly mentioned. BC is produced by incomplete combustion of carbonaceous material, mainly fossil fuels and biomass. For example, vehicular emissions contribute 25% of the global average BC production (Cooke et al., 1999). Carbon combustion products are usually classified as BC (soot) and organic carbon (OC). The BC component of the combustion particles have a typical size of 0.1μm, but very quickly after emission combine with other components to form mixed particles of larger size. Particles with BC and other light-absorbing components heat the air by absorbing solar radiation, converting the solar radiation into internal energy (raising the temperature of the particle), and emitting, at the higher temperature, thermal-infrared radiation, which is absorbed selectively by air molecules (Jacobson, 2004). BC particles are hydrophobic upon emission but can mix internally or externally with other more hygroscopic aerosol species, such as sulfate and OC (e.g., Moteki et al., 2007; Adachi et al., 2010). They are removed from the atmosphere within days to weeks by dry deposition and/or precipitation processes (e.g., Textor et al., 2006).

(2) Snow albedo

Snow albedo is important in determining the surface energy transfer within the cryosphere. For example, 100 ng g-1 of BC in snow of grain radius 1000 microns will reduce the visible-wavelength albedo by 10% (Fig. 1b of Warren, 2013). Major perturbations to the albedo of land surface can be caused by snow over very short time frames. Snow cover can increase the albedo of grassland by a factor of 3-4 and forested regions by a factor of 2-3 (Betts and Ball, 1997; Thomas and Rowntree, 1992). This has been emphasized by the results of many albedo sensitivity inter-comparison studies (Cess et al., 1991; Randall et al., 1994). Other studies (e.g. Barnett et al., 1988; Walland and Simmonds, 1996) have shown that the long-term average of snow accumulation or melt patterns may significantly alter regional climate and have a strong impact on the atmospheric general circulation.

(3) LAP in snow and ice (LAPSI)

Perceptions persist about the purity of fresh snow, but measurements tell more insightful story. Optical and electron microscopes showed that a typical snow crystal contains thousands of particles, including BC. It is hypothesized that wet deposition (via snow and rain) is the primary removal mechanism for airborne LAP. Dry deposition can also be significant, accounting for several tens of percent of deposition (Davidson et al., 1985). Surface darkening caused by dirty snow is likely to accelerate the snow aging process, which further lowers snow albedo and increases the speed of snowpack melt (Flanner et al., 2007). Once snow starts to melt, LAPs may accumulate at the snow surface, which increases LAP concentrations through the post-depositional enrichment, further darkening snow, warming the snowpack and accelerating the melting (Flanner et al., 2007; Doherty et al., 2013; Qian et al., 2014).

(4) Multiple effects of LAPs

In contrast to the straightforward radiative effect of greenhouse gases (GHGs) in the atmosphere, the mechanisms by which atmospheric LAPs influence the lower troposphere and snowpack energy budget are complex (Flanner et al., 2009). First, airborne LAPs can warm the troposphere via solar heating but cool the surface, which could potentially change the atmospheric stability and general circulation. Second, LAPs may serve as a source of cloud condensation nuclei and/or ice nuclei and affect microphysical properties of clouds, thus impacting precipitation. Third, deposited LAPs reduce snow reflectance (snow darkening) and upwelling shortwave radiative flux, and thus warm the surface (Painter et al., 2007).

(5) Snow albedo feedbacks

Modeling studies have suggested that this snow darkening mechanism has greater warming and snow-melting efficacy than any other anthropogenic agent (Hansen et al., 2005, Flanner et al., 2007; Qian et al., 2011, Skiles et al., 2012). This large impact results from a series of positive feedback mechanisms: (i) as melt commences some of the LAPs are not washed away by melt-water but accumulate at the surface or inside of snowpack (e.g., Conway et al., 1996; Flanner et al., 2007; Doherty et al., 2010); (ii) warming of the snow via reduced albedo increases the snow grain sizes, which further lowers snow albedo (e.g., Warren and Wiscombe, 1980; Hadley and Kirchstetter, 2012); (iii) With sufficient snow melt more of the darker underlying surface is exposed, leading to the well-known "snow albedo feedback" (e.g., Warren and Wiscombe, 1980; Hansen and Nazarenko, 2004; Flanner et al., 2007; Brandt et al., 2011; Hadley and Kirchstetter, 2012). The magnitude of the positive feedback through the LAP enrichment depends on the scavenging efficiency (SE) of LAP by snowmelt water (Qian et al., 2014). (Doherty et al., 2013) analyzed field measurements of the vertical distribution of LAP in the Arctic snow during the melt season and found a significant melt-induced amplification of BC concentration in surface snow, by up to a factor of five. (Xu et al., 2012) also revealed a post-depositional enrichment of BC in surface snow by measuring the redistribution of BC in snow-pits sampled monthly on a Tianshan glacier in Northern China.

(6) Importance of LAP effect

Disentangling the influences of aerosols from natural and anthropogenic sources is important for understanding how changes in biomass and fuel-burning might help mitigate glacial melt. More generally, understanding controls on snow cover timing and extent is important because most of the inter-annual variability in mid- and high-latitude planetary albedo is caused by changes in snow and sea-ice cover (Qu and Hall, 2006). Shifts in the timing and amount of runoff from snowmelt and glacier melt due to the combined effects of climate change and LAPs drive uncertainty in predicting runoff and stresses local water supplies. The relative contributions of absorption of sunlight by LAP and climate change to the loss of glacier mass and modification of the glacier runoff are unknown, due to the poor understanding and model uncertainty associated with LAPSI process and to the lack of sufficient in situ and distributed observations of spectral albedo, concentrations of LAP, and coincident changes in snow depth, ice height, snow water equivalent, and ice mass.

In this paper, we review the progress made in measuring the LAPSI in Arctic, Tibetan Plateau and other mid-latitude regions (Section 2) and in modeling the mass concentrations, radiative effects, and climatic and hydrological impact of LAPSI (Section 3). A few recent studies focus on BC deposition in Antarctica (e.g. Bisiaux et al., 2012a, 2012b), but we didn't include them in this review because BC concentrations in Antarctica are too low to be climatically important. We close by identifying research directions that would help reduce uncertainties in quantifying the impact of LAPSI on global and regional climate and hydrological cycles (Section 4).

Four primary measurement techniques have been used to measure LAPs in snow and ice. None of these are direct chemical measurements; instead each leverages specific properties of the particles (e.g. volatility, light-absorbing properties, solubility) to constrain mass mixing ratios of relevant components (BC/EC, OC, mineral dust). Measured mass mixing ratios of BC in particular may differ using different techniques, depending on what properties are used to define its presence (Bond et al., 2013). These measurements can be accompanied by chemical analyses, such as of particulate organic components (e.g. Hagler et al., 2007a, 2007b; Legrand et al., 2007; McConnell et al., 2007a, 2007b; Ming et al., 2008; Dang and Hegg, 2014), dust composition (e.g. Painter et al., 2007; Hadley et al., 2010), or snow water composition, or by source-receptor analysis, by back-trajectory analyses and factor analysis (e.g. Hegg et al., 2009, 2010) to elucidate absorbing particulate sources.

(1) OC and EC mass mixing ratios via thermo-optical (TO) analysis

Particles are collected by filtering snow water through quartz-fiber filters. OC and BC masses are determined by heating the sample and detecting CO₂ produced as a function of temperature, with BC/OC separation determined by optically monitoring changes in transmittance or reflectance as the particles are volatilized then combusted (e.g., Chow et al., 1993; Birch and Cary, 1996; Watson et al., 2005).

(2) Refractory BC mass size distributions and total mass mixing ratio via laser-induced incandescence

A nebulization system is used to inject snow water particles into a sample air stream. In the Single Particle Soot Photometer (SP2), individual particles are heated with lasers to incandescence temperatures. The incandescence emission is measured and used to quantify individual BC particle mass (Schwarz et al., 2006, 2008; Moteki and Kondo, 2010). In addition to BC mass size distributions, the SP2 also detects the coating thickness on BC.

(3) Spectrally-resolved total particulate absorption and estimate of BC mass mixing ratio

Particles are collected by filtering snow water through nuclepore filters. Filters are optically analyzed for spectrally-resolved absorption. Spectral absorption properties are used to apportion absorption to BC and non-BC components, and combined with volume of snow water filtered to estimate BC mass mixing ratios (e.g. Clarke and Noone, 1985; Grenfell et al., 2011; Doherty et al., 2010, 2014b). Below we will refer to this as the nuclepore-light-absorption (NLA) method.

(4) Gravimetrically-determined particulate mass mixing ratio

Particles are collected onto filters, followed by gravimetric analysis of total particulate mass. Combined with a measure of total snow mass, this yields snow mass mixing ratios of all particles. In regions where dust strongly dominates snow particulate mass, this is considered an estimate of snow dust mass (e.g. Aoki et al., 2006; Painter et al., 2012b; Painter et al., 2013a, 2013b).

Just as each of these techniques measures a different set of properties each has limitations and potential biases. Ideally, all four would be used in conjunction to provide complementary information but this is not generally practical. The method(s) used for a given study should be selected based on the properties and processes of interest and with an awareness of each techniques strengths and weaknesses.

The thermo-optical technique has been widely used for atmospheric measurements, nominally making snow and atmospheric measurements comparable, and it quantifies two components of interest (BC and OC). However, tests of atmospheric aerosol samples show that different operational protocols can lead to up to an order of magnitude range in estimated BC concentrations (Schmid et al., 2001). Optimized protocols have been developed to help minimize artifacts and improve OC/BC separation (Fung, 1990; Conny et al., 2003; Subramanian et al., 2006; Cavalli et al., 2010) but are not universally used. An additional issue for snow samples is that the capture efficiency of quartz-fiber filters for water samples is not well quantified and can be quite low. Capture efficiency of individual filters have been reported at 30% and lower (Hadley et al., 2010; Aamaas et al., 2011; Torres et al., 2013), with high variability from sample to sample (Aamaas et al., 2011; Torres et al., 2013) and limited understanding of controlling mechanisms.

The SP2 provides the most direct measurement of BC, the component primarily focused on in studies of climate forcing by snow albedo reduction. The SP2 is increasingly being used for atmospheric measurements, again allowing for comparability of atmospheric and snow data sets. However, the SP2 has challenges specific to analyzing snow water samples. First, different nebulization systems have widely different efficiencies in getting particles from snow water into a sample air stream (Schwarz et al., 2012; Ohata et al., 2013; Wendl et al., 2014). In addition, if BC in a sample is attached to larger particles (e.g. dust, sand), it may not be efficiently aerosolized and/or it may stick to sample tubing, leading to low biases. Second, the SP2 efficiently measures BC of ∼80-700 nm mass-equivalent diameter, with some variation depending on operational configuration. This size range encompasses the vast majority of atmospheric BC, but it appears that BC particles shift to larger sizes when incorporated in snow, especially with melting and refreezing either in the ambient or after snow collection (Schwarz et al., 2012, 2013), again possibly leading to low biases in measured BC.

Analysis of spectral absorption by snow particles allows direct measurement of the quantity of interest for snow albedo reduction, and the simplicity of the measurement allows for analysis of a large number of samples, but this technique does not measure what components lead to that absorption. Estimates of BC mass mixing ratios and apportionment of absorption to BC and non-BC components require making assumptions about the mass absorption cross-sections by these components (Clarke and Noone, 1985; Grenfell et al., 2011). While this property is relatively well constrained for bare BC, it is highly uncertain for the non-BC components or for coated BC. Thus, estimated BC mixing ratios can be highly uncertain (e.g. Fig. 16 of Doherty et al., 2010). Further, it appears that this technique over-estimates BC mass mixing ratios by up to a factor of 1.5-2.0 or more when a large (>50%) fraction of particulate absorption is due to non-BC components (Schwarz et al., 2012), possibly due to optical effects of particle/filter interactions.

Gravimetric measurements of snow particulates are technologically straightforward, nominally making it amenable to processing larger numbers of samples, and less prone to poorly understood biases. However, it requires collecting, melting and filtering sufficient snow that the mass of particles on the filter is greater than the uncertainty in the mass of the filter itself, which in most locations is not practical. In addition, the mass of a specific component——e.g. dust——can only be estimated where it is clear that component strongly dominates the total mass or where independent measurements such as those described above can be used to estimate, e.g., mass contributions by combustion aerosol versus mineral dust.

The first measurements of LAPs in Arctic snow were by (Clarke and Noone, 1985), who collected the particles on nuclepore filters and measured their spectral absorption. They also attempted to separate the contributions of BC and soil dust to the measured absorption, by use of the absorption Ångstr\"om exponent (AAE). That study was updated by (Doherty et al., 2010), who conducted large-area surveys of Arctic snow using essentially the same method for samples collected in 2006-09. Estimated BC mixing ratios (ng g^-1) were as follows: Greenland 3, Arctic Ocean 7, Arctic Canada 8, subarctic Canada 14, Svalbard 13, northern Norway 21, western Arctic Russia 26, and northeastern Siberia 17. In regions which (Clarke and Noone, 1985) had also sampled, the new measurements were somewhat lower than those of 25 years prior, consistent with the observed decline of BC in Arctic near-surface air (Sharma et al., 2013).

The AAE indicated that typically 40% of the light absorption by LAP is due to non-BC constituents. Chemical analyses of filters and meltwater, input to a receptor model, indicate that the major source of BC in most parts of the Arctic is biomass burning, but industrial sources dominate in Svalbard and the central Arctic Ocean (Hegg et al., 2010). When the snow surface layer melts, much of the BC is left at the top of the snowpack rather than carried away in meltwater, thus causing a positive feedback on snowmelt. This process was studied in Greenland, Alaska, and Norway (Doherty et al., 2013). In the percolation zone of South Greenland at the end of July, the subsurface snow had 2 ng g^-1 but the top 5 cm had ∼ 20 ng g^-1.

TO analyses for BC mixing ratio often obtain smaller values than the NLA method, but not always. In Svalbard, (Forsstr\"om et al., 2009) obtained a median BC mixing ratio of 4 ng g^-1, by comparison with 13 ng g^-1 for the NLA method. At the Summit camp on Greenland, TO measurements gave 0.2-1.2 ng g^-1 (Hagler et al., 2007b), compared to 1-4 for the NLA method. Values similar to those of Hagler et al. were obtained in northwest Greenland by (Aoki et al., 2014) using the TO method.

(Forsstr\"om et al., 2013) surveyed snow in the European Arctic using the TO method, obtaining BC mixing ratios of 40-80 ng g^-1 in Scandinavia, 11-14 in Svalbard, and 7-42 on the sea ice of Fram Strait. (Aamaas et al., 2011) measured BC in snow in the vicinity of coal-mining villages in Svalbard, showing the decline of coal-dust pollution with distance, with an e-folding distance of a few km.

Trends of Arctic BC for the past 200 years have been obtained from ice cores. To obtain monthly time-resolution in an ice core the SP2 method is appropriate, because it requires only a small sample size (a few ml of water). (McConnell et al., 2007b) analyzed two ice cores in southwest Greenland, where the air trajectories come mainly from North America. BC from forest fires gave a pre-industrial mixing ratio of 2 ng g^-1; then industrialization using coal resulted in a peak of ∼10 ng g^-1 about 1910. Improvements of combustion efficiency and replacement of coal by cleaner fuels caused the mixing ratio to decline to ∼ 2 ng g^-1 by 1950, similar to the preindustrial value.

Tibetan Plateau (TP) is the home of major components of the glaciers in High Asia, which are very sensitive to climate change and fluctuate with climatic cooling and warming. Concentrations of atmospheric LAPs over the TP highly depend on elevations and distance to emission sources (e.g. Cao et al., 2010), which also influence the deposition of LAPs onto snow/glacier and the induced snow darkening effect. Modeling studies have suggested that BC is the second most important factor (after GHGs) causing the fast melting of Himalayan and Tibet glaciers (Ramanathan and Carmichael, 2008; Menon et al., 2010), although uncertainties in the model results can be over one order of magnitude (Kopacz et al., 2011; Qian et al., 2011).

Investigators making field measurements first noticed the deposition of light-absorbing dust on snow in the Tianshan Mountains and their effects on reducing snow albedo (Zeng et al., 1984). (Fujita, 2007) found that the effect of dust deposition to snow on the runoff of a TP glacier depends on the surface albedo. When the surface albedo is above 0.7, the effect of dust can be neglected, but it becomes important when dust lowers snow albedo below this value.

BC came into the sights of researchers investigating the background concentrations of BC in TP glaciers many years later. (Xu et al., 2006) reported on BC concentrations in the snow samples from seven TP glaciers, and showed the spatial distribution of BC (in the range of 4 to 80 ng g^-1) in the highly elevated glaciers. The investigation of BC deposition in glaciers was later expanded to the broad area of High Asia including the Tianshan Mountains in the Xinjiang autonomous region (Ming et al., 2009; Ming et al., 2013b). BC concentrations in High-Asia glaciers are similar to those in the Arctic and western American mountains (order of 10 ng g^-1; see Section 2.4), but are much less than in heavily industrialized areas such as northern China (Huang et al., 2011; Ming et al., 2013b; Wang et al., 2013b). Up to 2012, 18 High-Asia glaciers had been sampled for BC in snow and ice cores by Chinese researchers. These studies found that BC concentrations in High-Asia glaciers primarily depend on elevation (i.e. higher sites have lower BC concentrations) and secondarily on regional emission intensities, precipitation and snow melting conditions. Figure 1 shows the relationship between BC concentrations in snow/glaciers and the elevation of the sampling site for samples gathered on the TP, based on observational data published in the literature. There is a decreasing trend in BC with increasing elevation. Within the same glacier area, BC concentrations tend to increase with decreasing elevation (red and blue dots in the figure), which is likely due to the important role of post-depositional enhancement of BC in surface snow at lower elevations, where snow melt is more frequent and/or intense. Variations in BC for samples at the same elevation likely, to first order, reflect the impact of differences in emissions and precipitation on the atmospheric deposition of BC at different geographical locations.

Figure 1.  Scatter plots for BC mixing ratios (ng g^-1) in snow/ glacier and the elevation (m) of sampling sites over TP based on observations from Ming et al. (2008, 2009, 2013a,b), Wang et al. (2012) and Xu et al. (2006, 2009a, 2009b, 2012). The dashed line is a linear fit to all points. Dots in the same color (red or blue) indicate samples from the same glaciers.

The strong and fast melting of a TP glacier can cause BC particles to become concentrated at the snow surface. This process will result in a dark layer in the snow-pit profile of a glacier surface at the level where melting occurred, and then was buried by subsequent snowfall. The surface enrichment of BC has been reported in previous studies focusing on the Tianshan mountains and central TP (Ming et al., 2009; Xu et al., 2012). This process probably explains why high BC concentrations were observed in some glaciers far away from major emission sources.

Ice cores have been used to study the deposition history of BC in the TP region. The first BC-deposition record retrieved from Tibetan glaciers was reported by (Ming et al., 2008), who found a mean BC mixing ratio in snow of 20 ng g^-1 throughout the second half of 20 ^th century, and a fast increasing trend to 50 ng g^-1 after 1990. In that study, variations in the mixing ratio of BC in snow was primarily attributed to variations in atmospheric transport by summer monsoons. Further, they estimated the historical mean radiative forcing due to absorption by BC in snow at about 1 W m^-2. A more comprehensive record showed that mixing ratios of BC in snow were generally increasing in Pamir, a hinterland of TP, and in the Himalaya during 1990-2000 (Xu et al., 2009a).

In some sampled glaciers, mineral dust particles are the dominant insoluble impurities other than BC in terms of particle mass (Chen et al., 2013; Qu et al., 2014). The measurements of snow albedo and other physical parameters associated with LAPs in snow were first conducted on a winter snowpack on Mountain Nyainqentanglha, in the central TP, in 2011 (Ming et al., 2013a). This study suggested that the impact of dust on snow albedo can exceed that of BC in the TP, mainly due to much larger concentrations of dust. (Wang et al., 2014b) reported an increasing trend in OC concentration and the ratio of OC to BC in a southeastern TP glacier since early 1990s based on ice-core measurements (Xu et al., 2009a), and they suggested that more attention to snow darkening effect of OC is merited because of its non-negligible light absorption and the recent rapid increases in its emissions.

LAPs in seasonal snow can significantly reduce its albedo and increase snowmelt rates in mid-latitudes. Here, the snowpack is often exposed to more intense insolation and generally has higher mixing ratios of BC and dust in snow than at high latitudes; thus, the surface albedo at mid-latitudes can be more strongly affected by darkened snow, despite a shorter snow-cover season (Bond et al., 2013). While much of the focus has been on BC in Arctic snow, there may be seasonally and regionally significant radiative forcing due to BC in snow at mid-latitudes, with impacts on snowmelt rates and timing. Several field campaigns have been dedicated to measuring BC concentrations in seasonal snow and glacier ice at mid-latitudes.

2.4.1. North America

The earliest observations of BC in snow in North America were carried out in spring of 1980 in the Cascade Mountains of Washington (∼ 47^°N, 121^°W, at ∼ 1000\;m MSL) by (Grenfell et al., 1981). Samples were taken from the top 10 cm of the snow. For these samples, BC mixing ratios were estimated using an earlier version of the NLA method, where quartz-fiber filters were used instead of nuclepore filters. (Lin et al., 1973). Estimated BC mixing ratios were in the range of 22-59 ng g^-1.

(Clarke and Noone, 1985) collected fresh and 2-week old snow samples at Hurricane Hill on the Olympic Peninsula, Washington (48^°N, 123.5^°W, at about 1500m AMSL) in March 1984. These samples were analyzed using the NLA method (Lin et al., 1973). Here, BC mixing ratios were found to be larger in fresh snow samples (18.5 ng g^-1) than in the aged samples (10.1-15.4 ng g^-1).

During the winters of 1982-85, (Ch\'ylek et al., 1987) collected 29 snow samples at collection sites in southern New Mexico and west Texas (32^°N, 106^°W), and analyzed them using the TO method. They found the average mixing ratio of BC in snow to vary between 4.9 and 15.9 ng g^-1, which is of the same order of magnitude as previous estimates of BC in snow in Washington State's Cascade Mountains (Grenfell et al., 1981) and Olympic Peninsula (Clarke and Noone, 1985).

Fresh snow samples were also collected during several snow storms in the winter of 1995-96 in Halifax, Canada (45^°N, 64^°W) by (Ch\'ylek et al., 1999) and analyzed using the thermo-optical method. The estimated BC mixing ratios in the Halifax area varied between 4.3 and 32 ng g^-1, with an average value of 11 ng g^-1.

The melting of the snow packs in the Sierra Nevada and Southern Cascade mountain ranges is one of the crucial contributions to the source of fresh water for California's agriculture and population. (Hadley et al., 2010) measured BC mixing ratios in snow at the three remote sites in the Sierra Nevada snowpack between late February and mid-April of 2006 using a modified version of thermo-optical analysis, as described in detail in (Hadley et al., 2008). The estimated BC mixing ratios in snow ranged from 1.7 to 12.9 ng g^-1; (Hadley et al., 2008) argue these mixing ratios are sufficient to perturb both snowmelt and surface temperature. The mixing ratios of BC in snow in the Sierra Nevada were similar that in falling snow near Halifax, Nova Scotia as measured by (Ch\'ylek et al., 1999), but a factor of 2-4 lower than predicted by (Flanner et al., 2007) and (Qian et al., 2009).

While (Hadley et al., 2008) focused on BC in snow in the Sierra Nevada, (Painter et al., 2012b) show that dust plays a very significant role in snow particulate light absorption in the San Juan range of the Sierra Nevada mountains. Here, the loading of dust in snow is so heavy in late Spring that it is visually obvious from photographs or satellite (e.g. see Fig. 2 of Painter et al., 2012b). Spring coincides with both large dust deposition events (see their Fig. 7) and the onset of melt, which removes snow water from the surface of the snowpack but appears to leave most dust particles in the surface snow. Assuming that the snow particulate mass is dominated by dust, Painter et al. gravimetrically determined the mass mixing ratio of dust in snow. At the end of the melt season, mixing ratios are on the order of <1 μg g^-1 to almost 5 μg g^-1 depending on the year (2005-10), compared to ∼ 1-10 ng g^-1 BC mixing ratios of (Hadley et al., 2008), for Sierra Nevada snow samples). While the mass absorption efficiency of BC is much greater than that of dust, this does not compensate for the three-orders of magnitude difference in mass mixing ratios.

More recently, (Doherty et al., 2014b) reported on a survey of seasonal snow conducted in January, February and March 2013 across 67 sites in the northwest U.S., the North American Great Plains and the Canadian boreal forest. They used NLA method to measure snow particulate light absorption and estimate mass mixing ratios of BC. Average mixing ratios of BC in snow across sites in the Pacific Northwest, Intra-mountain Northwest and Canadian sub-regions of their study were 20-37 ng g^-1; in the Great Plains region, regional-average mixing ratios were higher: 44 ng g^-1 for surface snow, and 78 ng g^-1 for sub-surface snow. In the Northern U.S. Plains, on average about 50% of the visible-wavelength light absorption was due to non-BC constituents, which (Doherty et al., 2014b) conclude is mostly local soil. At some sites in this region, they estimated nearly all (>85%) of absorption was due to soil in the snow, likely from local sources. A source apportionment analysis attributes BC in the snow to a mix of pollution (fossil fuel burning) and biomass burning (including anthropogenic biomass buring, e.g. wood stove emissions), with the relative rolls of each varying from site to site. As noted above, the data from all of these field measurements are useful for evaluating aerosol-climate models in simulating BC in snow and further investigating the climate impacts of BC and other light-absorbing impurities in snow in mid-latitude.

2.4.2. North & East Asia

Developing countries in South and East Asia have become major BC emission source regions, starting in the 1980s (Novakov et al., 2003; Bond et al., 2007). Dust emissions are also significant in parts of East Asia, specifically northern China. Being downwind of China, Japan is a receptor region for dust storms coming off of the deserts of North China. Accordingly, (Aoki et al., 2006) studied dust in snow at a site in Sapporo, Japan, December, 2003-March, 2004. As with (Painter et al., 2012b), (Aoki et al., 2006) assumed that snow particulate mass was dominated by dust and used gravimetric analysis to estimate dust mass mixing ratios in snow. For most of the study period, LAP mass mixing ratios in the snow were ∼ 1-10 μg g^-1, but in March they rapidly rose to ∼ 100 μg g^-1, and remained elevated through the end of the month. This high-dust period was also the period when snow was melting, leading to amplified surface snow mixing ratios of dust. Coincident with the snow sampling, (Aoki et al., 2006) measured broadband albedo, and show that the visible-wavelength broadband snow albedo varies from ∼ 0.8-0.9 preceding the heavy dust event. Following the event, albedo declined to below 0.8, reaching a low of <0.4 at the end of the study period.

The first large-area survey of LAPs in seasonal snow at 46 sites spanning six provinces of northern China was conducted in January and February of 2010 (Huang et al., 2011; Wang et al., 2013b). About 400 snow samples were collected and analyzed using the NLA techniques pioneered by (Clarke and Noone, 1985), and later modified and used for the Arctic snow survey by (Doherty et al., 2010) and North American survey of (Doherty et al., 2014b). A portion of snow melt-water was used for a source-attribution study, reported by (Zhang et al., 2013), analogous to that done for the (Doherty et al., 2010) Arctic samples (Hegg et al., 2009, 2010) and North American samples (included in Doherty et al., 2014b).

(Wang et al., 2013b) found that BC in snow over northern China had large spatial variability. The lowest concentrations were in the remote northeast China, along the south border of Siberia, where mixing ratios of BC in surface snow were in the range of 50-150 ng g^-1, with a median of 117 ng g^-1. In contrast, in industrial northeast China (south of 48^°N), the BC concentrations were typically in the range of 1000-2000 ng g^-1. Across the grassland of Inner Mongolia most sites had BC concentrations in the range of 100-600 ng g^-1. In the northeast China, LAPs in snow were dominated by BC. However, in the arid Qilian Mountains, at the northern border of TP, local soil and desert dust were the dominant sources of LAPs in snow, as indicated by the yellow-brown color of the sampling filters, so that Wang and colleagues did not report the values of BC concentrations in this region. Dust also significantly contributed to snow particulate absorption in the Inner Mongolia study sites. (Zhang et al., 2013) did a source-attribution study to explore the sources of LAPs in the snow. They found that three major sources contributed to the measured light absorption in snow over the entire sampling region: soil dust (53% of absorption), industrial pollution (27%), and biomass and/or biofuel burning (20%). In the arid Qilian Mountains, the soil dust accounted for ∼ 90% of absorption by LAPs on average.

A second large-area field campaign was conducted in January and February of 2012, to measure the LAPs in seasonal snow in western China. Snow was sampled at 36 sites covering the northern part of Xinjiang and the southeast of Qinghai provinces. Using visual estimates of absorption by particles on nuclepore filters (sampled/filtered using the same techniques as Wang et al., 2013b), (Ye et al., 2012) estimated that LAPs in snow in Xinjiang were dominated by BC particles, but with much lower BC concentrations than found in northeast China by (Wang et al., 2013b).

In global models, LAPs in snow are typically represented in a snow model within the land surface component of the model (e.g., Oleson et al., 2010; Watanabe et al., 2010; Yasunari et al., 2014). Models can use either prescribed or prognostically-determined (calculate by the model interactively) atmospheric aerosol concentrations and aerosol deposition fluxes, and they can use either prescribed meteorology (e.g. reanalysis data sets) or meteorology as simulated by the model to determine snowfall rates. Snow BC mixing ratios, snow albedo, radiative forcing, and climate impacts can be calculated for a single location or over regional or global domains. Depending on the model, feedbacks that affect snow BC mixing ratios and snow properties that determine albedo (e.g. snow grain size) and snow cover may or may not be included (e.g., Fig. 29 and Table 20 of Bond et al., 2013). These and other model configuration choices and——for prognostic runs——the emissions data set used as input to the model will affect calculated snow BC mixing ratios (e.g., Bond et al., 2013; M\'en\'egoz et al., 2014; Yasunari et al., 2014).

(Bond et al., 2013) recently reviewed BC in the climate system and they summarized model studies of BC in snow and ice. Here we review modeling of BC and the other LAPs such as dust and OC in snow and ice with our perspective that includes more recent model studies of BC and other LAPs in snow. Therefore, this section's role should be the combination of the complementation of (Bond et al., 2013), recent updates, our own perspectives, and possible uncertainties on the modeling. Here we also summarize the global models that currently treat LAPs in snow (see Table 1).

3.1.1. Comparison of model simulated BC in snow to observations

Estimates of how LAPs in snow affects snow albedo reach back to as early at 1960-80s (e.g., Onuma et al., 1967; Higuchi and Nagoshi, 1977; Warren and Wiscombe, 1980; Grenfell et al., 1981; Clarke and Noone, 1985). However, a few modeling studies in 2004-2005 (e.g., Hansen and Nazarenko, 2004; Jacobson, 2004; Hansen et al., 2005) drew renewed attention to the BC snow albedo effect (SAE) and its influence on climate. Following this, forcing by this effect was included in the AR4 of the Intergovernmental Panel on Climate Change (IPCC, 2007; Forster et al., 2007). In these and later global model studies, the BC SAE has been the primary focus. However, as noted in Section 2 (see references therein), field measurement have shown that dust is a significant, and sometimes dominant, source of snow particulate light absorption in regions such as the U.S. Great Plains, the western Great Plains of China and the southwest U.S., as well as likely some glacial regions.

Since publication of IPCC AR4 (IPCC, 2007), the number of global modeling studies of the BC SAE has increased rapidly (e.g., Flanner et al., 2007, 2009; Koch et al., 2009; Menon et al., 2010; Skeie et al., 2011). The early study of (Hansen and Nazarenko, 2004) used a simplified representation of the effects of BC on snow albedo, with fixed albedo reductions across broad regions, essentially as a sensitivity study. Following (Hansen and Nazarenko, 2004), (Hansen et al., 2005) also used a simple estimate on snow albedo reduction caused by BC, which was proportional to local BC depositions calculated by the method of (Koch, 2001).

(Jacobson, 2004) prognostically calculated snow BC mixing ratios and showed that the globally modeled BC in snow and sea ice in his model were in general agreement with the observations available at the time (i.e. Clarke and Noone, 1985; Ch\'ylek et al., 1987; Warren and Clarke, 1990; Ch\'ylek et al., 1999; Grenfell et al., 2002). However, as (Jacobson, 2004) pointed out, a problem with this comparison was that the BC emission inventory used in the model was for the year 1996 but some BC measurements were carried out in the 1980s.

(Flanner et al., 2007) simulated BC in snow and calculated the radiative transfer and heating by the Snow, Ice, and Aerosol Radiative model (SNICAR) model (Flanner and Zender, 2005), and more comprehensively compared the simulated BC with the available observations (see their Table 2). Again, however, as for the comparison of (Jacobson, 2004), the observed and modeled years were different. (Flanner et al., 2007) reported a wider range of differences in the comparison of modeled snow BC mixing ratios, in this case as modeled and measured in top 2-cm of the snowpack (see their Fig. 4). With the difficulty of the comparisons caused by the time mismatch, the discussion beyond the climatologically possible variations (CPV) at each location was impossible. (Qian et al., 2009) used a regional Weather Research and Forecasting (WRF-Chem) model to track BC in snow with the algorithm of Jacobson (2004) and calculated the snow albedo reductions caused by the BC in snow over Western United States. Their simulations showed the range of BC in snow of 10-120 ng g^-1, which was close to the reported values from the observations (Grenfell et al., 1981; Clarke and Noone, 1985; Ch\'ylek et al., 1987; Hansen and Nazarenko, 2004, and references therein) and global model simulations (Jacobson, 2004; Flanner et al., 2007). However, the comparison by (Qian et al., 2009) again suffered the problem of mismatch in the time periods of observations and simulations. Using the GEOS-Chem CTM with the assimilated meteorological fields by the Goddard Earth Observing System, version 5 (GEOS-5; e.g., Rienecker et al., 2008, 2011), (Wang et al., 2011a) applied a very simplified method to calculate BC mass concentrations in snow during 2007-09, which was just calculated as the ratio of deposited BC to water fluxes. They showed a similar range of the modeled BC concentrations in comparison to the observations by (Doherty et al., 2010).

(Kopacz et al., 2011) also applied the similar method with the GEOS-Chem model to estimate BC content in snow over the TP and Himalayan (HTP) region. However, because these studies didn't consider the many processes that can affect BC in snow after deposition (so-called post-depositional process) (Conway et al., 1996; Doherty et al., 2010, 2013; Aamaas et al., 2011; Sterle et al., 2013), these comparisons can only serve to check the possible range of simulated BC concentrations in snow if all the deposited BC were well mixed in the snow, with no post-depositional effect. (Qian et al., 2011) used the same model by the National Center for Atmospheric Research (NCAR) Community Atmosphere Model (NCAR-CAM 3.1), as (Flanner et al., 2009) did, to examine the impact of carbonaceous aerosols in the atmosphere and BC and dust in snow on the Asian monsoon climate and hydrological cycle over the Tibetan Plateau and Himalaya (TPH) region. Their simulated BC concentrations in snow over the southern slope of the TP were greater than 100 ng g^-1 (exceeding 200 ng g^-1 in some areas, see their Fig. 4). This was larger than the observed concentrations over wide areas of the TP and Himalayas (Xu et al., 2006; Ming et al., 2009). (Qian et al., 2011) point out that comparisons between their simulations and observations were only possible at the order-of-magnitude level due to the following issues: (1) the representativeness of one grid point in the model and the local observations was different; (2) the time period discrepancies between the multi-year mean from the simulations and observations made in specific years; (3) discrepancies in the depth of the snow layer where BC mixing ratios are estimated in the model versus the depth where it is measured in the observations. The simulated snow cover in CAM3 was much larger than that shown in Moderate-resolution Imaging Spectroradiometer (MODIS) satellite data. It's worth noting that the bias in snow cover fractions has been largely corrected in the later version of NCAR-CAM5.1 model, as shown in (Qian et al., 2014).

(Dou et al., 2012) carried out simulations using the NASA GISS composition-climate model (GISS-E2-PUCCHINI), driven by a linear relaxation of winds from two reanalysis data (NCEP: Kalnay et al., 1996; MERRA: Rienecker et al., 2011). They analyzed model results for 2006-09, using annually-repeating year-2000 monthly anthropogenic emissions (Lamarque et al., 2010) and monthly biomass burning emissions within each modeled year from the Global Fire Emission Database (GFED) version 3 (van der Werf et al., 2010). They found that the ratio of observed-to-modeled snow BC mixing ratios from two model runs were in the range of 0.73-1.81 for the Arctic Ocean, Canada and Alaska, Russia, Svalbard, and Greenland (see their Table 4). The models underestimated BC in snow over Russian Arctic in spring (see their Fig. 3), and they attributed this mainly to bias in the biomass burning emissions in the simulations.

(Bond et al., 2013) compared global modeling results from (Flanner et al., 2009) and (Koch et al., 2009) to observations made in a number of regions (Clarke and Noone, 1985; Ch\'ylek et al., 1987, 1995; Cachier and Pertuisot, 1994; Hagler et al., 2007a, 2007b; McConnell et al., 2007b; Forsstr\"om et al., 2009; Doherty et al., 2010). In Table 21 of (Bond et al., 2013), the range of ratios in modeled-to-observed snow BC mixing ratios varied from 0.6 to14 in the comparison to (Flanner et al., 2009) and 0.3 to 2.3 in the comparison to (Koch et al., 2009). In both cases, the largest biases were for Greenland in summer. For the study by (Flanner et al., 2009), they attributed the very large (factor of 14) summertime high bias in Greenland to the difficulty in representing snow melt processes on BC in snow.

More recently, (Zhao et al., 2014) coupled the SNICAR model (Flanner et al., 2007) with the WRF-Chem model and performed a case study of BC and dust in seasonal snow over Northern China. They showed that the model-simulated spatial variability of LAPs in snow is generally consistent with observations (Wang et al., 2013b), resolving local maximum mixing ratios of BC in snow of up to 5000 ng g^-1 near sources, although some quantitative uncertainties remain. Their study represents a significant effort in using a regional modeling framework at high resolution and better-resolved physics to study the effect of LAPs in snow.

Two recent studies (Lee et al., 2013; Jiao et al., 2014) carried out multi-model comparisons of globally modeled BC in snow, using results from global models that participated in the Atmospheric Chemistry and Climate Model Intercomparison Project (ACCMIP; Lamarque et al., 2013; also see at: http://www.giss.nasa.gov/projects/accmip/) and the Aerosol Comparisons between Observations and Models (AeroCom) project (e.g., Kinne et al., 2006; Schulz et al., 2006, 2009; also see at: http://aerocom.met.no/aerocomhome.html), respectively. The two studies applied similar methodology, in which monthly aerosol deposition fluxes from the global models and meteorology (i.e. snowfall rates) were prescribed model inputs, and these were used to calculate off-line mixing ratios of BC in snow. Both studies used the Community Land Model version 4 (CLM4; Oleson et al., 2010; Lawrence et al., 2011) and Community Ice CodE version 4 (CICE4) (Holland et al., 2012). CLM4 includes in-snow processes that affect snow BC mixing ratios. (Lee et al., 2013) used the year-2000 aerosol deposition fluxes given by the 8 ACCMIP models and prescribed meteorology (e.g. snowfall) from the years 1994-2000 as inputs to the model, then reported averaged quantities 1996-2000. (Jiao et al., 2014) used dust and BC deposition fluxes given by the set of AeroCom models (see their Table 1 and references therein), and coupled this with six-hour resolution reanalysis data for years 2004-09 for meteorology. Here the reanalysis data used was a blended data set from the NCEP/NCAR reanalysis (Kistler et al., 2001) and from the Climatic Research Unit (CRU; see at: ftp://nacp.ornl.gov/synthesis/2009/frescati/model_driver/cru_ ncep/analysis/ readme.htm). The two studies both compared their simulated snow BC mixing ratios to the observations from the field measurements of (Doherty et al., 2010), covering a wide area in the Arctic. (Lee et al., 2013) also discussed the historical comparisons of BC deposition and concentrations in snow derived from the available ice-core measurements (McConnell et al., 2007b; McConnell and Edwards, 2008; McConnell, 2010; Ming et al., 2008; Thevenon et al., 2009; Xu et al., 2009a; Bisiaux et al., 2012a, 2012b). An issue common to the two studies is that they used monthly-mean aerosol deposition in models to compare with the sporadic point observations. In (Lee et al., 2013), the differences in the mixing ratio of BC in snow between the off-line simulations and observations were mostly within a factor of 2-3 on average, except for the Arctic Ocean (model underestimate by a factor of 2-5 in NCAR-CAM5.1) and Greenland (model overestimate by a factor of 4-8). Similar comparisons between the off-line simulations and the observations by (Jiao et al., 2014) showed that the differences were mostly within one order of magnitude, with some exceptions for both the phase I and II AeroCom cases under their settings of meltwater SE (see their Figs. 2 and 3).

According to what we have learned from these two studies on the multi-model comparisons with prescribed deposition fluxes and meteorology (Lee et al., 2013; Jiao et al., 2014), we still have much difficulty in validating model simulated BC in snow beyond understanding CPV at the grid points where the observations were carried out. The temporal mismatch also causes difficulty in the validation because BC mixing ratios can have large monthly fluctuations even within a single winter season at a given location, caused by (for example) the migration of BC in snow with melt (e.g., Aoki et al., 2011; Skeie et al., 2011; Forsstr\"om et al., 2013). In addition, (Doherty et al., 2014a) point out that the method of calculating BC mass mixing ratios in snow when CLM4.0 is run using prescribed aerosol deposition fluxes, such as done in the (Lee et al., 2013) and (Jiao et al., 2014) studies, introduces a high bias, due to the decoupling of snowfall and BC wet deposition rates.

3.1.2. Comparison of model simulated non-BC LAP in snow to observations

As discussed in previous studies (e.g., Higuchi and Nagoshi, 1977; Warren and Wiscombe, 1980; Painter et al., 2007, 2010, 2012a, 2012b; Gautam et al., 2013; Wang et al., 2013b; Doherty et al., 2014b), dust and soil in snow also have an important role in snow albedo reduction. Most of the global modeling studies as reviewed in the previous section have mainly focused on BC SAE. A limited number of global models have considered the effect of both the BC and dust (or even OC or tar balls) in snow and/or ice (Table 1).

(Aoki and Tanaka, 2008) started considering both the dust and BC SAE in a land surface model (Hosaka et al., 2005) coupled with a chemical transport model (Model of Aerosol Species IN the Global AtmospheRe, MASINGAR: Tanaka et al., 2003, 2007) at Meteorological Research Institute (MRI) in Japan. They considered the Snow Impurity Factor (SIF), which is the product of mass mixing ratio in snow and mass absorption coefficients. The SIF was calculated for both dust and BC, and was used to calculate snow albedo reductions (Aoki et al., 1999, 2000), using the Physically Based Snow Albedo Model (PBSAM) (Aoki et al., 2011). Recently MRI has developed a global climate model, MRI-CGCM3 (Yukimoto et al., 2012), which is a part of their earth system model (MRI-ESM1; Yukimoto et al., 2011) and considers the aerosol mixing in snow, based on (Aoki et al., 2003). The PBSAM above is coupled to the Snow Metamorphism and Albedo Process (SMAP) model (Niwano et al., 2012) in the Hydrology, Atmosphere, and Land (HAL) model. Using these models, Aoki and Tanaka (2008, 2011) calculate the radiative forcing due to BC and/or dust in snow (see Section 3.2.1) but they did not report any model-observation comparisons of the mixing ratios of BC and dust in snow.

(Flanner et al., 2009) carried out global simulations of that accounted for darkening of snow by both dust and BC using the NCAR-CAM3.1 model (Collins et al., 2006). They found that observed declines in springtime snow cover over Eurasia were better captured in the simulation when forcing by BC in snow was included. The treatment of both BC and dust SAE has been included in more recent versions of the NCAR Community Earth System Model CESM/CAM5 model (version 4 or later), but an earlier version of the model used by Flanner et al. (2007, 2009) and (Qian et al., 2011) that included forcing by dust and BC in snow was not publicly available at that time (Table 1).

More recently, several global models have been developed to account for the SAE due to both BC and dust (Watanabe et al., 2010; Watanabe et al., 2011; Bentsen et al., 2013; M\'en\'egoz et al., 2013, 2014; Yasunari et al., 2014). (Jacobson, 2012) updated an earlier version of his model that included forcing by BC in sea ice and snow to also include dust and tar balls (mentioned later) depositions to snow. The Japanese global models, the Model for Interdisciplinary Research on Climate (MIROC) for climate (MIROC5; Watanabe et al., 2010) and the earth system modeling (MIROC-ESM; Watanabe et al., 2011), now include the depositions of BC and dust to snow in their land surface model, Minimal Advanced Treatments of Surface Interaction and Runoff (MATSIRO; Takata et al., 2003). As mentioned by (Watanabe et al., 2010), the MATSIRO treats BC and dust in snow as a combined mass (i.e., no separation on these constituents, calculated with total deposition-adjusted mass to the snow so that it reflects the relative difference of dust and BC absorptions in snow). The snow albedo scheme, considering snow aging, was based on (Yang et al., 1997) and their dirt component was modified to incorporate the snow impurity concentrations (Watanabe et al., 2010; Watanabe et al., 2011; Nitta et al., 2014). However, MATSIRO also calculates the combined mass concentration without the deposition adjustment, which was, for example, provided as a variable to the Coupled Model Intercomparison Project Phase 5 (CMIP5) project (e.g., Taylor et al., 2012; also see at: http://cmip-pcmdi.\break llnl.gov/index.html) (Kumiko Takata, Kengo Sudo, Shingo Watanabe, and Masahiro Watanabe, 2014, personal communication). Hence, the simulated LAP mass concentrations from this model can only be compared to any observations on the combined mass of BC and dust.

The Norwegian Earth System Model (NorESM) (Bentsen et al., 2013) uses the main components of the Community Climate System Model version 4 (CCSM4; Vertenstein et al., 2010; Gent et al., 2011) including CLM4 (Oleson et al., 2010; Lawrence et al., 2011) with SNICAR (e.g., Flanner and Zender, 2005, 2006; Flanner et al., 2007, 2009) and CICE4 (Holland et al., 2012), but also including their own components (Bentsen et al., 2013; also see Table 1). The Laboratoire de M\'et\'eorologie Dynamique (LMDZ) model with Interaction of chemistry and aerosol (LMDZ-INCA) (Hauglustaine et al., 2004; Krinner et al., 2005, 2006; Hourdin et al., 2006; Balkanski et al., 2010; M\'en\'egoz et al., 2013; Szopa et al., 2013), uses a parameterization to account for particle flushing rates through the snowpack during melt, and this has been used to investigate how dust affected the Asian glaciers during the last glacial maximum (Krinner et al., 2006). Using this model, (M\'en\'egoz et al., 2014) simulated BC and dust mass concentrations in snow over the Himalayas at two different model horizontal resolutions (a coarse resolution of ∼ 350 km and a stretched fine one of ∼ 50 km), and compared this to ice core data obtained from the Mera Glacier (Ginot et al., 2014). Due to a lack of reliable observations, M\'en\'egoz et al. (2013, 2014) chose to neglect the flushing effect both for dust and BC, assuming that the most less-hydrophilic particles stay at the surface of the snow cover. They showed a good model-observation comparison for dust in snow, but large differences for BC in snow (up to a factor of 60 in annual mean and a factor of 30 in simulated inter-monsoon and monsoon values).

(M\'en\'egoz et al., 2014) further discussed several points as follows to explain the large discrepancies in BC: (1) the model resolution induces discrepancies in the altitude in the model gridbox versus at the sampling site; (2) differences in the vertical resolution of the modeled snow layer versus that in the ice core sample, and the lack of consideration of the flushing effect of BC by liquid water in the model; (3) uncertainty in the winter snow layer in the ice core data due to the effects of strong winds (Wagnon et al., 2013; Ginot et al., 2014) were not reflected in their model; (4) larger uncertainties in the actual BC concentrations in the atmosphere and snow result from the methodologies and sample-treatment differences in these measurements (Kaspari et al., 2014; Petzold et al., 2013; Lim et al., 2014), for example, SP2 (Kaspari et al., 2011; Ginot et al., 2014) vs. thermal optical method (Ming et al., 2008); (5) an additional point after (M\'en\'egoz et al., 2014) to explain the discrepancies is the horizontal heterogeneity of BC in snow, which could not be reproduced in the coarse gridded models (Martin M\'en\'egoz and Gerhard Krinner, 2014, personal communication).

In addition, the ice-core dating itself by (Ginot et al., 2014) is uncertain, and the lengths of monsoon and non-monsoon seasons vary year by year. That is why the quantitative comparisons at the seasonal scale (inter-monsoon and monsoon) are not reliable. To reduce this sort of uncertainty, annual-scale comparison, like the ACCMIP comparison with the annual ice core data (Lee et al., 2013), is more appropriate if the ice core data do not have very accurate time separations in dating (i.e., monthly or shorter separation). As also learned from (M\'en\'egoz et al., 2014), over the high-elevation regions like the Himalayas, there are lots of difficulties to carry out reliable comparisons on LAP in snow, even just for the discussions on CPV (Lee et al., 2013; Jiao et al., 2014), because of all the possible uncertainties.

(Flanner et al., 2009) and (Aoki et al., 2011) pointed out the importance of including SAE of light-absorbing OC in future studies of snow-albedo effect. Indeed, in the Arctic, (Doherty et al., 2010) found that on average about 40% of light absorption is due to non-BC components. Since almost all the particulate light absorptions in snow in this pan-Arctic study were attributed to combustion particles (Hegg et al., 2009, 2010), the non-BC absorption is likely OC (or "brown carbon"; Moosm\"uller et al., 2009, and references therein). Several global models currently have the capability of considering OC in snow. This includes CESM, which use the CLM4 land model (Oleson et al., 2010; Lawrence et al., 2011). However, the effect of OC is not considered in snow over sea ice in CICE4 (Holland et al., 2012), and in NorESM (Bentsen et al., 2013). The model of (Jacobson, 2012) accounts for tar balls as well as soil dust in snow and sea ice. Tar balls, which can be clearly distinguished from OC and BC, are generated in smoke from biomass and biofuel burning, and have a typical particle size range of 30-500 nm and low hygroscopicity (P\'osfai et al., 2004; Adachi and Buseck, 2011). The NASA GEOS-5 model recently incorporated a snow darkening module, called GOddard SnoW Impurity Module (GOSWIM), in the land surface model (Yasunari et al., 2014), which considers OC in snow in addition to dust and BC treatments in land snow. To the best of our knowledge, (Yasunari et al., 2014) is the first study to validate simulated mass concentrations of dust, BC, and OC from off-line and on-line experiments against snow impurity observations, as given by (Aoki et al., 2011) for bi-weekly observations covering one winter at Sapporo in Japan. (Skeie et al., 2011) carried out a similar validation, but only for BC in snow. The limitation of (Yasunari et al., 2014) is that the validation was only done for the LAPs in the surface snow at a single location, and the effects on the observations of changes between the snow sampling timings, such as snow impurity re-distributions and flushing amount, were still uncertain. (Yasunari et al., 2014) concluded that, in their model, the deposition rates of dust and BC needed to be increased by factors of 4.3 and 3.06, respectively, to explain the mass mixing ratios of dust and BC in snow at Sapporo in Japan. A recent global modeling study by Lin et al. (2014) carried out global off-line SDE simulations for organic aerosol radiative effect over the snow and ice (see Section 3.2.1) with CLM4 (Oleson et al., 2010; Lawrence et al., 2011) and CICE4 (e.g., Holland et al., 2012), using the aerosol depositions from the Integrated Massively Parallel Atmospheric Chemical Transport (IMPACT) model simulations, though they did not report the simulated mass concentrations of the organic aerosols in snow and ice.

3.1.3. Possible uncertainties in globally modeled snow impurities

Uncertainties in model representation of each of the following processes may propagate to the simulated LAPs in snow: (1) Emissions; (2) Transport; (3) Deposition; and (4) Post-depositional processes. Here we briefly explain these with some examples.

(1) Emissions: Some global modeling studies (e.g., Flanner et al., 2009; Skeie et al., 2011; Forsstr\"om et al., 2013) used the GFED biomass burning emissions (van der Werf et al., 2006, 2010), while NASA GEOS-5 used the Quick Fire Emissions Dataset (QFED) [e.g., Petrenko et al., 2012; Darmenov and da Silva, (2014)\footnote[1]Darmenov, A. S., and A. da Silva, 2014: The Quick Fire Emissions Dataset (QFED)——Documentation of versions 2.1, 2.2 and 2.4. Technical Report Series on Global Modeling and Data Assimilation, NASA/TM-2014-104606. [Will be available at: http://gmao.gsfc.nasa.gov/pubs/tm/]]. The total emissions in QFED are, in general, higher than in GFED [Darmenov and da Silva, (2014) ^a; Fig. 12 in Ichoku and Ellison, 2014]. (Dou et al., 2012) argued that the underestimates of BC in snow over the Russian Arctic in their simulations were due to biases in biomass burning emissions. Therefore, comparisons of model runs with different up-to-date biomass burning inventories are needed to test how this influences modeled mixing ratios of LAPSI. Emission uncertainties may also come from other source sectors, both in terms of overall magnitude of emissions and their temporal and spatial distributions, which will impact on the seasonal variations of airborne LAP deposition to snow/ice. Recently, Wang et al. (2014a) implemented an explicit source tagging technique in the CAM5 model to provide a detailed characterization of the destiny of BC particles emitted from different geographical regions and source sectors (e.g., biomass burning vs. fossil fuels), and quantified the sensitivity of BC deposition in the Arctic to regional BC emission uncertainties. This modeling framework has also been applied to quantify contributions from regional emission sources to the observed BC in a TP glacier and assist in explaining the observed long-term trend (Wang et al., 2014b). For dust emissions, (Textor et al., 2006) reported an annual mean dust emission rate and its defined diversity within the available AeroCom global models of 1840 Tg yr^-1 and 49%, respectively. As summarized in Table 1, none of the global models used in (Textor et al., 2006) considered SAE of dust at the time of their publication. However the large uncertainty in dust emissions among the models shown in the Textor et al. study likely reflects similarly large uncertainties in model-simulated dust in snow. Inter-model comparisons of dust in snow and sea ice should also be a focus of future studies.

(2) and (3), Transport and Deposition: (Qian et al., 2014) recently carried out sensitivity studies using global simulations where they changed flushing (or scavenging) efficiencies by liquid water in the snowpack and snow-aging factors from the default settings with the two versions of the CAM5 (Neale et al., 2010) with the improvements in aerosol transport and atmospheric wet removal processes made by Wang et al. (2013a; called IMPRV simulation). In (Qian et al., 2014), the improvements of the atmospheric transport and wet deposition processes in the model by (Wang et al., 2013a) reduced BC deposition in mid-latitudes and thus allowed the transport of more BC to the Arctic, followed by more deposition into the Arctic snow, which largely removed both the high bias in BC in snow over Northern China and the low bias in the Arctic, in comparison to the observations (Doherty et al., 2010; Wang et al., 2013b) (see Fig. 1 of Qian et al., 2014). This suggests that model performance on deposition processes has an important role in determining the global distribution of snow LAPs in different regions. Over the Himalayan region, (Yasunari et al., 2013) compared several estimates of BC dry deposition amounts from the off-line simulations and the output results from two global models, and reported the impacts on snow albedo reductions, in which the lowest bound of deposition amount was estimated with a fixed slow dry deposition velocity as given by (Yasunari et al., 2010). (Yasunari et al., 2013) suggested that snow surface roughness and winds are the keys to reduce the uncertainties of in dry deposition of BC over snow surface, which is likely applicable to the other non-BC LAPs as well. In addition, as discussed in (M\'en\'egoz et al., 2014), large spatial variability of both dry and wet depositions may also affect the model-estimated LAPs in snow.

(4) Post-depositional Processes: An important in-snow process that can affect surface snow LAP mixing ratios——and therefore snow albedo——is the efficient or inefficient flushing of particles through the snowpack with snow melt water. This appears to lead to either significant reduction or enhancement of LAP at the snow surface during melt depending on the LAP properties such as particle size and hydrophilicity, but the magnitude of this effect is highly uncertain because of the very limited number of observations studying this processes (Higuchi and Nagoshi, 1977; Conway et al., 1996; Doherty et al., 2010, 2013; Aamaas et al., 2011; Sterle et al., 2013). As summarized in Table 1, some global model accounts for the flushing effect of the snow impurity (e.g., Jacobson, 2004; Flanner et al., 2007; Koch et al., 2009; Yasunari et al., 2014) but some models do not (e.g., Watanabe et al., 2010; Watanabe et al., 2011; Skeie et al., 2011; Forsstr\"om et al., 2013; M\'en\'egoz et al., 2013, 2014). Some, such as LMDZ-ORCHIDEE-INCA, can account for this effect depending on its settings (Krinner et al., 2006; M\'en\'egoz et al., 2013, 2014; Table 1).

(Qian et al., 2014) conducted a sensitivity study on how BC in snow responds to perturbations to two important parameters controlling post-depositional processes in snow, meltwater SE and snow-aging factor, using the version of CAM5 that includes the improvements in wet deposition discussed by (Wang et al., 2013a). Although the magnitude of perturbations was somewhat arbitrary (factors of 0.1 and 10 were applied to the default values), their results showed that the influence of changing model parameters on BC in snow can be very significant and has regional/seasonal dependence.

Figure 2 is an example of comparisons on BC mass concentration in the top of snow layer simulated by three global models for March-May mean in 2008: Oslo-CTM (Forsstr\"om et al., 2013; Ragnhild Skeie, 2014, personal communication), GEOS-5 and CAM5 (Qian et al., 2014). The results clearly showed large differences in the simulated BC mass mixing ratios in snow among the different global simulations with and without the BC flushing effect in the. The sensitivity simulation with a low scavenging efficiency of BC by melt water (SEl from Qian et al., 2014) showed BC mass mixing ratios more similar to those simulated by OsloCTM2 (Forsstr\"om et al., 2013) than did the simulation with a high scavenging efficiency (SEh from Qian et al., 2014). The GEOS-5 simulated snow BC mixing ratios were roughly within the range between the SEl and SEh CAM5 simulations. As also discussed by (Qian et al., 2014), the comparisons here strongly support the idea that uncertainties in the flushing effect have a large impact on the simulated LAP mass concentration in snow. These large uncertainties on the flushing effect propagate to uncertainties on snow albedo, radiative forcing, and related feedbacks between the land and atmosphere.

Figure 2.  Comparisons of simulated mean BC mass concentrations in the top snow layers for March-May (MAM) in 2008 by three global models: (a) A replay run GEOS-5 simulation, for which the MERRA re-analysis atmospheric fields by Rienecker et al. (2011) were used. The mass concentration was calculated from the outputs of the total mass of hydrophobic and hydrophilic BC and snow water equivalent (SWE) in the top snow layer. The main settings were mostly similar to the Exp GEOS-5 or TE1 experiment in Yasunari et al. (2014), but started from a different year (June 29, 2001) and the saved monthly mean data for MAM in 2008 were used. The snow darkening module in GEOS-5, GOSWIM, only calculated snow albedo effect (SAE) over the land surface components (i.e., excluding the land ice and snow over sea ice) when SWE was of equal to or more than 0.01 kg m^-2 at each time step (Yasunari et al., 2014), and the shaded grids in gray probably may indicate complete loss of the BC mass from the top snow layers or non SAE calculations because of considered smaller amount of SWE; (b) An OsloCTM2 simulation, for which the data were used in Forsstr\"om et al. (2013) and provided by Ragnhild Skeie (2014, personal communication). Then, this panel was adapted from Fig. 6 of Forsstr\"om et al. (2013) with the provided simulated data; (c) and (d) NCAR CAM5 simulations for the SEh and SEl (high and low SE by melt-water) experiments in Qian et al. (2014), for which the outputs were obtained from their study and the mass concentrations were calculated with the mass of BC and SWE in the top snow layer. See the differences on SAE treatment among the global models in Table 1. The selected map orientation, time period, and scale of the color bar for this figure was based on Fig. 6 of Forsstr\"om et al. (2013).

Overall, we can conclude that the flushing (or scavenging) effect has an important role in determining the concentration of LAP in snow and its uncertainty during the snow-melting season, as also discussed by (Doherty et al., 2013). However, reliable parameterizations of this flushing effect and other in-snow process in models, within the limited number of global models that consider SAE (Table 1), are still very difficult because of the poor understanding of the post-depositional processes due to the limited number of observations on these processes (Higuchi and Nagoshi, 1977; Conway et al., 1996; Doherty et al., 2010, 2013; Aamaas et al., 2011; Sterle et al., 2013).

3.2.1. Radiative forcing

A small initial snow albedo reduction may have a large net forcing because the induced warming affects the snow grain size, sublimation rates, and melt rates of snow, all of which enhance LAP-induced snowpack albedo reduction, resulting in an amplification of the radiative forcing (e.g. Flanner et al., 2007; Qian et al., 2009; Bond et al., 2013). (Hansen and Nazarenko, 2004) published the first global model study estimating the radiative forcing of BC in snow/ice with a simplified approach that did not prognostically calculate BC concentrations in snow/ice. (Hansen and Nazarenko, 2004) used some BC-in-snow measurements available at the time to estimate the adjusted radiative forcing induced by the assumed BC due to visible snow albedo reductions of 2.5%, 1%, 5%, and 1%, respectively (see their Table 2), in the Arctic, Greenland, other snow-covered areas in the northern hemisphere, and the snow-covered areas in the southern hemisphere except for Antarctica. (Hansen and Nazarenko, 2004) reported the adjusted an radiative forcing value of 0.16 W m^-2 from the one case above (Case 1) but this was revised to 0.08 W m^-2 by (Hansen et al., 2005) because they discussed that BC contamination in snow at the period of (Clarke and Noone, 1985) should be larger than that at the present, which was similarly reviewed in (Aoki and Tanaka, 2011). Based on these two studies, IPCC AR4 (IPCC, 2007; Forster et al., 2007) adopted a forcing range of 0.10 0.10 W m^-2. Although (Jacobson, 2004) did not report forcing values for BC in snow (only climate impacts), (Bond et al., 2013) estimated the present-day adjusted forcing as 0.08 W m^-2 based on his modeled data (see their Table 19), in good agreement with the estimate of (Hansen et al., 2005). Several studies either followed the approach of (Hansen and Nazarenko, 2004) by applying regionally uniform snow albedo changes (Wang et al., 2011b), scaled albedo changes based on model-predicted BC deposition rates (Hansen et al., 2005; Shindell and Faluvegi, 2009), or used on-line predicted BC content in snow to calculate the snow albedo change and its relevant radiative forcing (Flanner et al., 2007, 2009; Aoki and Tanaka, 2008, 2011; Rypdal et al., 2009; Koch et al., 2009; Skeie et al., 2011). (Aoki and Tanaka, 2011) separately calculated the radiative forcing values of BC, dust, and BC+dust in snow with a 15-yr model integration in their Table 7.2 of +0.09 0.27, +0.04 0.30, and +0.42 0.35 W m^-2, respectively (+0.22, +0.20, and +0.7 W m^-2, respectively from the 3-yr integration simulations in Aoki and Tanaka, 2008; see their Table 1). These forcing values were calculated as the annual mean difference of net shortwave radiation at the top of atmosphere (TOA) under whole sky conditions, with and without snow impurities, including feedback processes between the land and atmosphere (Teruo Aoki and Taichu Y. Tanaka, 2013, personal communication). Lin et al. (2014) recently reported adjusted radiative forcing for organic aerosols in snow and sea ice with off-line global simulations (see Section 3.1.2), in which cases with assumed high (low) absorption by organic aerosols produced +0.0025 (+0.0009) and +0.00055 (+0.00016) W m^-2 over the land snow and sea ice, respectively (see their Table 10). These estimates included forcing by both primary and secondary sources of organic aerosol, and high forcing estimates in this model were as large as 24% of the forcing by BC.

(Bond et al., 2013) provided a comprehensive overview on the radiative forcing of BC, both in the atmosphere and in snow. Here we do not intend to summarize the main results again from the literature regarding the forcing of LAPSI, instead, we present a few highlights from (Bond et al., 2013). The best estimates of forcing by them were +0.04 [+0.01 to +0.09] W m^-2 and +0.035 [+0.008 to +0.078] W m^-2, respectively, for all-source and industrial-era adjusted forcing by BC in snow (see their Section 8.1). They used scaling methods, which were applied to the modelled values, and obtained the adjusted forcing estimates consistently from all sources (fossil fuel, biofuel, and biomass burning) at the present day. For their forcing in the industrial era, they used a scaling method, based on the ratio of industrial-era to all-source BC emissions, so that the fractional differences between pre-industrial (1750) and present day on the BC contributions on snow from fossil fuel-plus-biofuel and open biomass burning emissions were accounted.

Recent high-profile summaries of BC forcing on snow or ice include those by (Bond et al., 2013) and the IPCC AR5 (IPCC, 2013; Boucher et al., 2013). The IPCC AR5 adopted a radiative forcing of BC on snow and sea-ice of +0.04 W m-2, based on (Bond et al., 2013). (Bond et al., 2013) reported the adjusted industrial-era forcings (1750-2010) from BC in snow (+0.035 W m^-2) and sea-ice (+0.011 W m^-2). The sum of these forcings (+0.046 W m^-2) was subjectively rounded down to +0.04 W m^-2 in AR5 because many models, such as the versions of CAM applied in previous derivations of the forcing (Flanner et al., 2007; 2009), tend to overestimate snow cover and, therefore, BC-in-snow forcing on the TPH (Qian et al., 2011). Radiative forcing has several definitions (e.g., Bond et al., 2013, Table 2) and here we distill some of the essential differences, for many readers who are not in the field of radiative forcing, in commonly-used metrics: (1) "instantaneous", (2) "adjusted", and (3) "effective" radiative forcings. (1) The instantaneous forcing includes no feedback, and simply represents the difference (usually calculated every time-step) in net TOA or tropopause radiative flux caused by the immediate presence of LAPs in snow and sea-ice. (2) The adjusted forcing includes both the instantaneous direct influence of LAP in snow and also changes in net flux due to "fast feedbacks" such as cloud response, snow grain size, and the surface enrichment of LAPs in snow and sea-ice, which happen on timescales on the order of days. (3) The effective forcing includes, in addition to fast feedbacks encompassed by the adjusted forcing, the radiative effects of longer-term feedbacks such as changes in snow and sea-ice cover, which may generally happen on monthly or seasonal time-scales. Such changes in snow cover, caused by earlier melt induced by LAPs, likely represent the largest amplifier of forcing from LAPs in snow. The effective forcing can be often derived by multiplying the instantaneous or adjusted forcing by an "efficacy" factor, such as those derived by (Hansen et al., 2005), (Flanner et al., 2007), or (Hansen et al., 2007) from equilibrium or transient climate change simulations. For the case of BC-in-snow, this efficacy is roughly 3 (Flanner et al., 2007; Hansen et al., 2007; Bond et al., 2013), implying that a unit of instantaneous or adjusted forcing from BC in snow will trigger three-fold greater equilibrium temperature response than the same adjusted forcing by CO₂. However, one can also obtain the effective forcing from annual mean or climatological mean differences in the net TOA radiation budget if the long-term feedbacks are reflected in the difference (Myhre et al., 2013, AR5 Chapter 8). If the computational method applies prescribed sea surface temperature and neglects sea-ice feedbacks, this technique would likely provide a smaller forcing estimate than the effective forcing reported in (Bond et al., 2013), which was derived from the product of the efficacy and adjusted forcing.

As pointed out in previous sections, TPH are one of hottest spots in investigating the climate impacts of LAPSI. Indeed, model simulations of (Qian et al., 2011) using CAM3 estimate global maximum LAPSI mixing ratios over the Tibetan Plateau. Because of the high LAPSI content and large incident solar radiation at the low latitudes and high elevations of the Tibetan Plateau, it exhibits the largest surface radiative flux changes induced by LAPSI of any snow-covered region in the world. In their study, LAPSI-induced snow albedo perturbations generate surface radiative flux changes of 5-25 W m^-2 during spring, with a maximum in April or May. (Flanner et al., 2007) also estimated that the largest annual average surface forcing induced by BC-in-snow is over the TP, namely, 1.5 W m^-2 averaged over the whole region. Their forcing averaged over the only snow-covered areas reached a maximum of about 10 to 20 W m^-2 in spring. These are adjusted forcings, because they include some of the short-term feedback process in the snow. The modeling study by (Flanner et al., 2009) also included snow albedo reduction due to dust. However, Flanner et al. (2007, 2009) and (Qian et al., 2011) studies used a model that significantly over-estimates snow cover across the Tibetan Plateau, and therefore over-estimates the forcing due to reduced snow albedo in this region. (M\'en\'egoz et al., 2014) estimated a surface net adjusted forcing of 1-3 W m^-2 over the Himalayas, where the adjusted forcing again includes the effects of snow aging processes (see their Fig. 3c).

(Kopacz et al., 2011) calculated instantaneous radiative forcing (without and with the presence of BC particles in snow) of +3.78 to +15.6 W m^-2 at the five grid points with snow covers over TPH where are the locations of the five glaciers, with minima in the winter (approximately in the range of 3 to 11 W m^-2 across the different sites) and maxima in summer (approximately 7 to 16 W m^-2) (see their Fig. 4). They used the GEOS-Chem simulated BC deposition fluxes to calculate the BC mass concentrations in snow and the relationship between BC content in snow and snow albedo given by (Ming et al., 2009). The reduction in albedo due to dust or soil in the snow is not accounted for in their study. The recent offline simulations of LAPSI forcing and in situ radiation flux measurements in the Himalayan and TP glaciers have presented more reliable results (Ming et al., 2013b). They found the impacts of LAPSI explained less than 8% and 2% of the net radiation flux on the TP and Himalayan glaciers (at one site from each), respectively, and they concluded that the influence of BC on snow was not significant.

While limited measurements and modeling studies have been conducted with a focus on TPH a more robust estimate of forcing by BC-in-snow in this region is still hindered by (1) too sparse and short-term measurements of LAPSI, which have large spatial and seasonal variations, (2) limited data on the mixing ratios of dust in snow together with the mixing ratios of BC in snow, in a region frequently subject to dust storms, and (3) the inherent difficulty of using a global model with coarse spatial resolution to accurately represent snow cover, atmospheric transport, and the deposition of BC and dust in a mountainous region with complex terrain (also see section 3.1.1 and in section 3.1.2).

3.2.2. Climatic and hydrological impact

Modeling studies investigating the climate impacts of BC in snow at the global scale indicate that BC in snow produces surface warming both in the Arctic and across the northern hemisphere (Hansen and Nazarenko, 2004; Jacobson, 2004; Hansen et al., 2005; Flanner et al., 2007, 2009; Koch et al., 2009; Shindell and Faluvegi, 2009; Rypdal et al., 2009; Goldenson et al., 2012). Those studies suggested that the climate efficacy of BC-in-snow is about 2-4 times that of CO₂, implying that LAP has a larger climate impact than its relatively smaller averaged direct radiative forcing would suggest.

Most of modeling studies on climatic and hydrological impact of LAPSI focus on the Tibetan Plateau and a few other high-mountain regions. Indeed, the Tibetan Plateau has long been identified as critical in regulating the Asian hydrological cycle and monsoon climate (e.g. Yanai et al., 1992; Wu and Zhang, 1998). On the one hand, the glaciers in the Himalaya and on the Tibetan Plateau act as a water storage tower for many Asian countries. Long-term trends and/or seasonal shifts in water supply provided by TPH may significantly affect agriculture, hydropower and even national security for the developing countries in the region. On the another hand, the Tibetan Plateau also exerts significant mechanical and thermal forcing that influences the South Asian and East Asian monsoon systems (Manabe and Terpstra, 1974; Yeh et al., 1979). Anomalous snow cover can influence the energy and water exchange between the land surfaces and the lower troposphere by modulating radiation and water and heat flux (Cohen and Rind, 1991), which in turn could affect rainfall in China and India in the subsequent summer (e.g. Wu and Qian, 2003).

Himalayan glaciers have been in general retreat since the mid-1800s, in particular the central to eastern Himalayan glaciers have been experiencing the most rapid retreat (Qin et al., 2006; NRC report, 2012). Generally the attribution of this glacier retreat has been to climate change, at times with erroneous and detrimental claims (e.g. the IPCC claim of year 2035 loss of all Himalayan glaciers) (Cogley et al., 2010). Observational evidence indicates that the surface temperatures on the Tibetan Plateau have increased by about 1.8^°C over the past 50 years (Wang et al., 2008). The IPCC AR4 (Meehl et al., 2007) reported that under the A1B emissions scenario, a 4^°C warming will likely occur over the Tibetan Plateau during the next 100 years. Both observed and projected warming over the Tibetan Plateau are much larger than for the global average.

Accelerated melting of the global snowpack and glaciers is generally driven by warming due to increasing greenhouse gas concentrations (Barnett et al., 2005), but the larger rate of warming and speed of glacier retreat in this region suggests that additional mechanisms may be involved (Xu et al., 2009a). For example, heating of the atmosphere by light absorbing aerosols and melting induced by LAP in glacial snow and ice could be substantial (e.g., Ramanathan and Carmichael, 2008, Lau et al., 2010, Qian et al., 2011). The TPH are located in close proximity to some of the most intense sources of strongly light-absorbing BC in the world (Bond and Bergstrom, 2006; Bond et al., 2013). For example, South Asia, especially the Indo-Gangetic Plain (IGP), is one of the largest BC emission sources in the world (Ramanathan et al., 2007). The southern side of the Himalaya is directly exposed to Indian emissions and more likely to be impacted by BC than the northern side (see Fig. 3). Himalayan ice core records indicate a significant increase in the deposition of both BC and particulate organic carbon to snow on TP, especially since 1990 (Ming et al., 2008; Xu et al., 2009a; Kaspari et al., 2011). The prevailing westerly air flow carries to the TPH a large flux of mineral dust aerosol from the arid regions of Southwest Asia, the Arabian Peninsula and the Thar Desert.

Figure 3.  Global distribution of annual mean surface forcing induced by BC + Dust in snow, averaged only when snow is present (Unit: W m^-2). [Reprinted from Qian et al. (2011)]

(Qian et al., 2011) conducted a series of numerical experiments with a global model (CAM3) to assess the relative impacts of anthropogenic CO₂ and carbonaceous particles in the atmosphere and snow on the snowpack over the Tibetan Plateau and subsequent impacts on the Asian monsoon climate and hydrological cycle. They found that BC-in-snow increases the surface air temperature by around 1.0^°C and reduces spring snowpack over the Tibetan Plateau more efficiently than the increase of CO₂ and carbonaceous particles in the atmosphere. As a result, runoff shows an earlier melt trend, i.e. increasing during late winter and early spring but decreasing during late spring and early summer. They also defined a so-called snowmelt efficacy, i.e. the snowpack reduction per unit degree of warming induced by the forcing agent. They found the snowmelt efficiency is 1-4 times larger for LAPSI than CO₂ increase during April-July. LAPSI directly increases the net solar radiation by reducing surface albedo, while CO₂ first warms the air by absorbing more longwave radiation then warms the surface and melts the snow by air-surface heat transfer. This may be a reason why LAP can accelerate snowmelt more efficiently than CO₂.

The simulations of (Qian et al., 2011) also show that during boreal spring aerosols are transported by southwesterly flow, causing some LAPs to reach higher altitudes and deposit to the snowpack/glacier over the Tibetan Plateau. While LAP in the atmosphere directly absorbs sunlight and warms the air, the darkened snow surface polluted by LAP absorbs more sunlight and warms the surface directly. Both effects enhance the upward motion of the atmosphere and spur deep convection along the Tibetan Plateau during the pre-monsoon season, resulting in a possible earlier onset of the South Asian Monsoon and increase of moisture, cloudiness and convective precipitation over northern India. In East Asia, LAPSI have a more significant impact on monsoon circulation in July than the CO₂ and LAP in the atmosphere. The role of the Tibetan Plateau as a heat pump is elevated from spring through summer as the land-sea thermal contrast increases to strengthen the East Asian Monsoon, probably due to the increase of both sensible heat flux associated with the warm skin temperature and latent heat flux associated with increased soil moisture. As a result, both southern China and northern China become wetter while central China (i.e. Yangtze River Basin) becomes drier, which is a near-zonal anomaly pattern consistent with the dominant mode of precipitation variability in East Asia. However, the snow impurity effects reported in (Qian et al., 2011) likely represent some upper limits as snowpack is remarkably overestimated over the Tibetan Plateau due to excessive precipitation. It will be critical to improve the precipitation and snowpack simulation for improving the estimates of the climatic and hydrologic effects of snowpack LAPs.

(Lau and Kim, 2006) proposed the so-called Elevated Heat Pump (EHP) effect, whereby heating of the atmosphere by elevated absorbing aerosols over northern India and Tibetan Plateau strengthens local atmospheric circulation, leading to a northward shift of the South Asian monsoon rain belt in the late boreal spring and early summer season. (Lau et al., 2010) found from climate model experiments that atmospheric heating and feedback induced by LAP (mainly dust and BC) could result in a surface air warming exceeding 1- 1.5^°C over the Tibetan Plateau, and the western Himalayas. Recent studies have also attributed the change of the South Asian monsoon system and accelerated warming of the troposphere over the Tibetan Plateau to atmospheric heating and surface cooling by aerosols (Gautam et al., 2009; Prasad et al., 2009; Ganguly et al., 2012). From GCM experiment, (Lau et al., 2010) found that as part of the EHP effect due to heating by LAP in the atmosphere, the tropospheric and surface heating over the Tibetan Plateau during the pre-monsoon months (April-May) are amplified by reduced cloudiness and adiabatic warming of subsiding air over the Tibetan Plateau, stemming from circulation changes induced by high LAP concentration over the Indo-Gangetic Plain and the Himalaya foothills. The authors also found that feedback processes associated with changes in surface heat and latent heat fluxes over the Tibetan Plateau can play an important role in accelerating the upper tropospheric warming over the Tibetan Plateau during the pre-monsoon season, leading to an early monsoon onset over northern India. The roles of LAPSI in affecting the atmospheric and terrestrial water cycles over the Tibetan Plateau, and their possible interaction with LAP heating and induced water cycle feedback processes in affecting the Asian monsoon are open questions, and subjects for future investigations.

Outside HTP works on LAPSI on hydrology are just emerging. (Painter et al., 2010) used the Variable Infiltration Capacity model to study the impacts of dust on snow albedo and runoff from the Upper Colorado River Basin, and found that the peak runoff at key hydrologic sites occurs on average three weeks earlier in 2003 versus in 1916, due to increased snow dust loading resulting from anthropogenic activities. They also found that increases in evapotranspiration from earlier exposure of vegetation and soils decreases annual runoff by more than 1.0 billion cubic meters or 5% of the annual average. (Qian et al., 2009) simulated the deposition of BC on snow and the resulting impact on snowpack and the hydrological cycle in the western United States (see Fig. 4), based on the chemistry version of the Weather Research and Forecasting model (WRF-Chem). With a higher spatial resolution of a regional model, WRF-Chem can better simulate the snowpack over mountainous region and capture the large spatial variability in BC deposition that reflects the localized emissions and the influence of the complex terrain than can a global model. The BC-induced snow albedo perturbations increase the surface air temperature, and reduce the snow net accumulation and spring snowmelt. These effects are strongest over the central Rockies and southern Alberta, where BC deposition and snowpack overlap the most. The change to surface radiation flux and temperature is around 50%-80% under a doubled snow albedo scenario against control simulation, but snowpack reduction is nonlinearly accelerated.

Figure 4.  Accumulated BC deposition (dry + wet) for particles less than 10 mm in diameter for the periods (a) DJF and (b) MAM. Depositions units are μg m^-2. Areas in WRF-Chem with seasonal mean snow cover greater than 1 cm in depth are overlaid with white hatching. [Reprinted from Qian et al. (2009)]

Future studies quantifying the climatic effects of LAPSI are needed to address some challenges. First, Accurate estimates of this climatic effect requires accurate representation of LAP deposition fluxes, snow and ice accumulation rates, snow aging processes, and the melt-water scavenging efficiency of LAP in the snowpack. All of these are difficult to constrain given the scarce observations and large variability in representing these processes from model to model. Second, climate forcing calculations must account for the feedback processes, including those associated with LAPSI effect, as discussed in Introduction Section. Not all of these feedback processes are included in current climate models. Third, interpretation across different model studies is complicated by the fact that the climate response to the radiative forcing and rapid adjustments differ substantially depending on the types of LAPSI, and model representation of snow and aerosol physics. Fourth, studies need to include the impact of all types of LAPSI. To date, the focus has largely been on the snow albedo impacts of BC, but mineral dust and organic carbon have both been shown to be important——and sometimes the dominant——particulate absorbers in snowpack in many regions. The role of dust/soil on snow albedo in glaciers is particularly under-studies, aside from the work in the U.S. southwest (i.e. Painter et al., 2010and related work).

As reviewed in previous sections,\;light-absorbing\;aerosols deposited on snow and glaciers can alter their surface albedo and accelerate melting. However, there are a number of critical details that hinder our understanding of the role of aerosols in the climate system. For example, the scientific community is not yet able to quantitatively characterize how changes in aerosols change particular components of the climate system (e.g. particular cloud system types, particular circulation features such as monsoons), and in particular the cryosphere such as glaciers on the Himalayas and Greenland. These challenges have broad ramifications for our understanding of the climate system. We do not know whether changes in aerosols have influenced important feedback mechanisms in the climate system, e.g. would adding aerosols or removing them significantly change the snow albedo response to CO₂ forcing?

What we suggest for the future global or regional modeling studies of the effects of LAP in snow are as follows:

(1) Currently, most of the previous comparisons, which were reviewed in section 3.1, were carried out with discrepancies between modeled time periods and in situ observation time periods (e.g., Lee et al., 2013; Jiao et al., 2014; Qian et al., 2014). For example, (Qian et al., 2014) compared their modeled climatological values with in situ observations (Doherty et al., 2010; Wang et al., 2013b). However, these kinds of comparisons cannot fully assess the bias between models and observations because of the mismatch in time. Therefore, we encourage more comparisons to be made during "the same time periods" from diurnal to seasonal time scales.

In addition, we should also pay attention to the calculation method of LAP mass concentrations in snow in off-line simulations with prescribed LAP deposition fluxes, since the relative rates of snow and LAP deposition fluxes are important as suggested by (Doherty et al., 2014a). They pointed out that high biases are possible if the prescribed inputs of LAP deposition fluxes (i.e., BC in their study) and snow deposition flux from independent data and/or from different times and sources are used, because of the inconsistency in the calculation of spatial/temporal averages between the inputs.

Furthermore, for the purpose of better comparisons between off-line models, the availability of observed data on vertical snow information such as snow depth for sampling (i.e., sampling intervals) together with LAP mass concentration (e.g., Doherty et al., 2010; Ye et al., 2012; Wang et al., 2013b) is essential to carry out proper comparisons of models to observations. Without the observed snow information, the uncertainty such as point (2) discussed in (M\'en\'egoz et al., 2014) (see Section 3.1.2) would increase.

(2) Comparisons between observations and the multi-model off-line simulations (Lee et al., 2013; Jiao et al., 2014), and single global model simulations (e.g., Bond et al., 2013; Forsstr\"om et al., 2013; Qian et al., 2014) showed that most of the differences fell in the range of one order of magnitude. However, currently we do not know how large a difference in LAPSI mixing ratios between simulations and observations is sufficient for present and future climate simulations and for meaningful discussion of the many feedbacks in global models. i.e.: What is the actual impact of an error in the mass mixing ratio of LAPSI of a factor of 2, 5, or 10 on local, regional, and global climate? In future studies, this point should be addressed in multi-model calculations with different treatments of LAPSI effects.

(3) Currently only a limited number of global models can simulate both BC and the other LAPs in snow (Table 1 and references therein). Progress is being made on better understanding of snow impurity processes, such as the mixture conditions of LAP and snow particles as discussed by (Flanner et al., 2012), observed flushing/scavenging efficiency of LAPs (e.g., Conway et al., 1996; Doherty et al., 2013), and the enrichment in snow-melting process (e.g., Aamaas et al., 2011; Xu et al., 2012; Painter et al., 2012b). Global model that currently don't represent these processes are strongly encouraged to have LAP-in-snow processes incorporated. Further, the global models listed in Table 1 should also be improved in the LAP treatments with better approximations based on future process-oriented studies. Then, we can better understand the possible uncertainties among many global models not only from off-line comparisons like (Lee et al., 2013) and (Jiao et al., 2014), but also from comparison of on-line global simulations.

(4) As suggested by (Forsstr\"om et al., 2013) and (Yasunari et al., 2014), more frequent snow impurity measurements at, e.g., daily resolution are necessary for model validation in terms of different years and locations. (Skeie et al., 2011) compared modeled BC in snow to the available observations at several locations including the Arctic region (Forsstr\"om et al., 2009). (Yasunari et al., 2014) compared modeled dust, BC, and OC concentrations in snow to available relevant observations (Aoki et al., 2011) at a single location in a time window from pre-melting to melting. In addition, weekly BC measurements from the snow surface were carried out in (Forsstr\"om et al., 2013). However, currently, vertically well-resolved daily or shorter-time observations of LAP in snow covering at least one winter are not available at any location. That is why we still do not know in detail the temporally frequent fluctuations of the impurities in snow. To monitor the behavior of the LAP movements in snow in detail, we need at least day-by-day sampling to capture the transition from pre-melting to melting time. Then, we can discuss more on the effect of flushing efficiency and movement of LAP in snow. In addition, this kind of day-by-day sampling should also be extended to several sites from mid-latitudes to polar regions for validating not only individual models but also multiple model inter-comparisons.

(5) In addition, we should consider the role of snow algae (i.e. in cryoconite) on glaciers and the ice sheets in the parameterizations of SAE in future studies, as also suggested by (Aoki et al., 2011). The snow algae affect snow albedo significantly in some locations/seasons (e.g., Takeuchi et al., 2001; Takeuchi, 2009; Stibal et al., 2012, and references therein). For example, (Takeuchi et al., 2001) showed the cryoconite obtained from Yala Glacier in the Himalayas that had granular like shape in darker color, including significant amount of humic acid and reducing spectral albedo (a flat shape in the visible range) on the glacier to the range of 10%-20% (see their Fig. 9b). They blamed the albedo reduction to the inclusion of humic material in the cryoconite, probably induced from the decomposed snow algae or/and other organics via biological activities. A lack of the consideration of the biological activities on snow may generate further biases in estimating the radiative forcing and climatic impact of LAPSI, in addition to the current uncertainties of BC and other LAPs deposited to snow on glaciers in particular.

Increasing measurements (in-situ or remote sensing) is the most important and urgent task for the near future. For example in the Himalayas and TP, although the uniqueness of the region and the importance of comprehensive observations in the atmosphere and at the surface have been recognized, to date only limited small-scale field campaigns have been conducted in this high-altitude region. A lack of observations in the Himalayas region has led to a wide range of explanations for observed glacier retreat, as well as to widely-criticized, erroneous projections of glacier retreat (IPCC, 2007; Cogley et al., 2010). Because of the critical dependence of the region's population on water from these glaciers and climate modulation by the snow and ice cover, and due to our poor understanding of important physical processes related to the glaciers, the U.S. National Research Council appointed a committee to analyze the scientific understanding of Himalayan glaciers, their impact on regional water cycle, and the impact of glacier change on the population of South Asia (NRC, 2012). They concluded that the lacking of observational data in the TPH region has led to misunderstandings about the effects of climate change on glacier retreat rates. To fill those established data gaps, a suite of in-situ and remote sensing measurements and a hierarchy of numerical models are critical.

The important science questions that remain to be answered are as follows:

(1) What are the optical/physical properties, spatial gradients, radiative effect, source and transport pathways, and deposition rates of LAP in atmosphere? What are the optical and physical properties and concentrations of LAP in the snow/ice on the ground, and their linkages with their counterparts in the atmosphere?

(2) How do temperature, precipitation amount and frequency, and the deposition of LAP influence surface albedo (and its potential feedbacks between land and atmosphere), surface melt, and runoff?

(3) What are the relative roles of atmospheric warming induced by CO₂, LAP or other warming agents, and snow surface darkening induced by LAPSI in accelerating the snowpack and glacier melting, changing water resources, and modifying global and regional climate in the next 20-30 years?

Successfully addressing these questions will lead to improved understanding of key processes and our ability to better model atmosphere-cryosphere interactions and aerosol-cloud-precipitation interactions over the Arctic and the mid-latitudes, leading to improved predictions of the impact of aerosols on climate, the cryosphere and the hydrologic cycle.