Effect of aerosol-induced snow darkening on the direct radiative effect of aerosols over the Himalayan region

Regional heterogeneity in direct and snow albedo forcing of aerosols over the Himalayan cryosphere was investigated using a regional climate model coupled with the community land model having snow, ice and aerosol radiation module. Deposition of absorbing aerosols like dust (natural) and black carbon (BC) (anthropogenic) decreases the snow albedo (snow darkening) over the Himalayas. Western Himalayas experiences a large reduction in the snow albedo (0.037) despite having lower BC mass concentration compared to central (0.014) and eastern (0.005) Himalayas. The contribution of BC and dust to the snow albedo reduction is comparable over the western and eastern Himalayas. The inclusion of aerosol-induced snow darkening in to the model reduces its bias with respect to the satellite derived surface albedo by 59%, 53% and 35% over western, central and eastern Himalayas respectively during the spring season. Since surface albedo decides the sign and magnitude of aerosol direct radiative forcing, aerosol induced snow darkening significantly affects the direct radiative effects of aerosols. Hence, the aerosol-induced decrease in snow albedo causes an early reversal in the sign of aerosol direct radiative forcing at the top of the atmosphere from warming to cooling over the western and central Himalayas, which can have implications in the radiation balance and water security over the region.


Introduction
Himalayan-Hindu-Kush mountains and Tibetan plateau encompasses the largest area of frozen water outside the north and south polar regions and thus called the 'Third Pole' . Despite its geographical position over the subtropics, the region maintains a significant amount of snow and ice due to its elevation (Bolch et al 2012). The seasonal snow and glacier melt support major rivers of Asia, which are the source of fresh water for the diverse ecosystem and billions of people inhabited in this region (Bookhagen and Burbank 2010). Himalayan mountains also play a significant role in regional hydroclimate by the dynamic and thermodynamic forcing (Sato andKimura 2007, Boos andKuang 2010).
The trends in the glacier mass budget are highly variable across the Himalayan region, with a decreasing trend over the central and eastern Himalayas and an increasing or no trend over the Karakoram region (Bolch et al 2012). Several studies attributed the accelerated melting of snow and glaciers to the warming due to greenhouse gases (Kulkarni et al 2002). Recent studies have shown that absorbing aerosols like black carbon (BC) and dust warms the upper troposphere by directly absorbing the solar radiation (Lau et al 2006) and reduces the snow albedo when deposited over the snowcovered region, thereby traps more solar energy within snowpack (Flanner et al 2007, Nair et al 2013, Yasunari et al 2015, Usha et al 2020, Li et al 2021. The snow darkening due to the deposition of absorbing aerosols (hereafter snow albedo effect or SAE) plays a major role in the radiative forcing, regional temperature, snow runoff and hydrological cycle over the Himalayan-Tibetan region (Nair et al 2013, Ma et al 2019, Sarangi et al 2019, Usha et al 2020. On the other hand, the direct interaction of atmospheric aerosols with radiation through scattering and absorption (aerosol direct effect or DRE) leads to surface cooling and atmospheric warming, which affect the snow over the region , Nair et al 2017. Since the magnitude and sign of DRE depend on the surface albedo (Haywood andShine 1995, Satheesh 2002), the effect of SAE on DRE needs to be investigated to reduce the uncertainty in aerosol radiative forcing.
Using in situ measurements and climate models, several studies were carried out over the Himalayan-Hindu-Kush-Tibetan region to understand the sources and concentration of absorbing aerosols in the atmosphere/snow and its implications on snow albedo and regional hydroclimate (Yasunari et al 2010, Ming et al 2012, Nair et al 2013, Sarangi et al 2019, Das et al 2020, Gul et al 2021. Observations by several investigators showed a high concentration of absorbing aerosols in the atmosphere and on snow over the entire Himalayan region (Ming et al 2008, Yasunari et al 2010, Nair et al 2013, Kaspari et al 2014. Modelling studies have shown a significant influence of aerosol deposition on snow albedo over the Himalayas. However, the importance of SAE in simulating the seasonal snow cover has not been investigated over the Himalayan region. Climate models have intrinsic difficulties in simulating aerosol fields over South Asia (Nair et al 2012) and seasonal snow cover over the Himalayan-Tibetan region (Qian et al 2011, Sarangi et al 2019. So regional climate model capable of simulating the aerosol and snow properties over the Himalayas has an added-value in climate impact assessments. By considering the large variability in the meteorology, snow properties and aerosol distribution along the Himalayan mountain stretch, the SAE and DRE over western, central and eastern Himalayas during the melting period (March-September) are examined in the present study. The high aerosol loading and solar insolation make the spring and summer monsoon period more conducive for strong aerosol-radiation and aerosol-snow interactions over this region. The SAE and DRE over the Himalayas during the melting period and the implications of SAE on the estimation of DRE are investigated for the first time over the region.

Model, experiment design and validation
Regional climate model (RegCM4.6) coupled with the online aerosol scheme and community land model having snow, ice and aerosol radiative (SNI-CAR) scheme is used for the simulation of aerosol DRE and SAE over the Himalayan region. The domain and model configurations are analogous to Usha et al (2020). Briefly, the model configuration has a horizontal resolution of 50 km with 18 vertical sigma pressure levels keeping model top pressure at 50 hPa. The parameterization schemes used are: the radiative transfer module from the NCAR model CCM3, the parametrization of convection is Grell scheme over land and Tiedtke scheme over the ocean, the sub-grid explicit moisture scheme for microphysics, the University of Washington turbulence closure model for boundary layer physics and the community land model version 4.5 (CLM4.5) for land surface processes. The radiative transfer across the snowpack is calculated using the SNICAR module of CLM4.5 (Flanner et al 2007(Flanner et al , 2009, which estimates the snow albedo based on the snow microphysics, ageing and the concentration of impurities. The SNICAR module assumes spherical snow grains and external mixing of aerosols and snow. The aerosol module comprises BC, organic carbon, sulfate, sea salt and dust (Solmon et al 2006). The dust scheme has four size bins (0.1-1 µm, 1-2.5 µm, 2.5-5 µm, and 5-10 µm) and the emission flux is estimated using CLM4.5 dust mobilization scheme. The emission flux of BC and organic carbon are taken from MaCCity inventory (Solmon et al 2006). The indirect effect of aerosols is not included in this study since RegCM4.6 does not have an explicit representation of the aerosol-cloud interaction. More details of the model dynamics and the parameterization schemes are already discussed in Giorgi et al (2012). Lateral and initial meteorological conditions are provided by ERA-Interim reanalysis dataset and sea surface temperature is taken from National Oceanic and Atmospheric Administration Optimum Interpolation dataset.
Three sets of experiments were carried out for the period 2010-2015 with the first year (2010) as spinup. These experiments include (a) control experiment (CNTRL) without the effect of aerosols, (b) simulation with the effects of aerosols in the atmosphere and within the snowpack (AERO) and (c) the third experiment (AeroATM), with the same configuration of AERO run but without impurities in the SNICAR module to exclude SAE. Thus, the AERO simulation includes both DRE and SAE, whereas AeroATM has DRE only. Clear-sky DRE values are reported in this study to avoid the influence of change in cloud cover on DRE while estimating the effect of SAE on DRE using AERO (with SAE) and AeroATM (without SAE) simulations.
The meteorological and aerosol fields simulated using RegCM4.6 have been validated with groundbased as well as satellite observations. Usha et al (2020) have carried out a detailed evaluation of the model performance with the same configuration. A brief description of the evaluation experiments is given in the supplementary file (available online at stacks.iop.org/ERL/16/064004/mmedia). It is found that the model is able to capture the spring and monsoon precipitation reasonably well (23% and −11% bias over the central and eastern Himalayas respectively). Most importantly, the model could simulate the spatial distribution and seasonal evolution of the snow cover fraction from spring to summer (June to September) monsoon season. The model simulation of the west to east gradient of the snow cover which persisted throughout the year (although there exists an overestimation over the western part of the Himalayas during spring) is shown in figure S1. Further, the model capabilities in simulating the aerosol parameters (BC mass concentration and aerosol optical depth (AOD)) have been evaluated using the measurements from the ground-based network stations (figure S2). AOD values are within 25% bias with respect to the measurements over different parts of the Indian region (Usha et al 2020) and the ratio of simulated to measured values of BC ranges from 1.38 to 0.39 during spring. The BC and dust deposited on the top layer of snow are also evaluated with respect to the in-situ data available in the literature, where the bias ranges from −4.5 to 50 µg kg −1 (RMSE of 26.8 µg kg −1 ) for BC and 15 to 29 mg kg −1 (RMSE of 12.6 mg kg −1 ) for dust. On the backdrop of the confidence gained from the model evaluation, the radiative effect of the aerosol-induced snow albedo feedback across the Himalayan region is examined in the following section.

Snow darkening due to aerosols
The topography of the Himalayan region and its subdivisions (western, central and eastern Himalayas, following Bolch et al (2012)) are shown in figure 1.
The mean values for all the variables discussed in this manuscript are selected above the 2500 m within each subdivision since it is observed that the snow cover extends up to a minimum altitude of about 2500 m over the Himalayas (Kulkarni et al 2010). Even though the snow cover and snowline elevation across the Himalayas is highly heterogeneous in space and time and it is challenging to fix a particular altitude as the cut off for the entire Himalayas, we have considered a cut-off altitude of 2500 m to minimize the direct influence of local sources on DRE and SAE values. The central and eastern Himalayas experience higher aerosol loading compared to the western Himalayas during spring as shown in figure 1(b). The seasonal variation of AOD depicts spring and summer high over central Himalayas (figure S3). The concentration of BC aerosols deposited on the top layer of snow shows that central Himalayas has the higher BC content up to 116 µg kg −1 than the other parts of Himalayas. This could be due to its proximity to the Indo-Gangetic Plains, which is one of the major aerosol hotspots over Asia. Several studies have reported the transport of BC during spring from the Indo-Gangetic Plains and eastern low-lying regions to the Himalayas through synoptic as well as local circulations (Cong et al 2015. The prevailing large-scale westerlies also transport dust from the major desert regions over central and west Asia to the Himalayas and Indian region. The concentration of dust aerosols deposited on the top layer on snow decreases from west (43 mg kg −1 ) to east (28 mg kg −1 ).
The aerosol-induced changes in the snow albedo (figure 1(c)) show a drastic decrease from west to east. The fraction of the area covered with snow is larger over the western Himalayas and is persistent over a longer time than the eastern Himalayas (figures S1 and 1(b)). During summer monsoon season, AOD increases by about 0.06-0.1 with respect to spring over central and eastern Himalayas respectively and the fractional snow-covered area reduced drastically over the Himalayas. The impact of BC and dust on the snow albedo are comparable over western and eastern Himalayas, whereas the impact of BC significantly dominates over the central Himalayas. The albedo reduction due to BC over the western Himalayas is higher than the central and eastern Himalayas, even though the BC on snow is lower over the western Himalayas. The albedo change due to dust over the western Himalayas is 4.7 times higher than the eastern Himalayas even though the dust loading on snow over the western Himalayas is only 1.5 times higher than the eastern Himalayas. Hence, this contrasting pattern of change in snow albedo and the amount of BC and dust deposited on the snow surface clearly shows that snow albedo reduction not only depends on the amount of aerosols deposited but also on the snow microphysics and snow thickness. For example, 100 µg kg −1 of BC on snow reduces the broadband albedo of the fresh snow by 1.6% and aged snow by 3.4% for a constant snow depth of 1 m. Since the western Himalayas have high snow thickness (1.39 m) compared to the eastern Himalayas (depth 0.08 m), 81.2 µg kg −1 of BC on western Himalaya results in higher albedo reduction than that of 100.6 µg kg −1 of BC on eastern Himalaya. Even though the impurity in the shallow snowpack reduces the albedo, this reduction is compensated by masking the influence of the underlying ground. Thus, albedo change is reduced over the shallower snowpack. It is interesting to note that, though the dust (mg kg −1 ) loading on snow is significantly higher than the BC (µg kg −1 ), the change in albedo due to BC and dust are comparable, which is due to the high absorption crosssection of BC compared to dust. Usha et al (2020) and Sarangi et al (2020) have also reported comparable contribution of BC and dust on snow albedo change over the Himalayas. The basic assumption of spherical snow grain and external mixing of aerosols and snow could lead to the overestimation of albedo change. Dang et al (2016) estimated that for a spherical snow grain of effective radius of 100 µm, the BC of 100 µg kg −1 reduced the albedo by 0.019 whereas it was 0.012 for non-spherical equidimensional grain. The mean albedo effect decreases by 31% for nonspherical snow grains relative to spherical grain while the internal mixing of BC enhances the mean albedo changes by 30%-60% with respect to external mixing (He et al 2018).

Radiative effects
Monthly evolution of the regional mean values of DRE (at the top of the atmosphere (TOA) and surface) and SAE for western, central and eastern Himalayas are shown in figure 2. The SAE values decreased from spring to summer monsoon whereas the mag- nitude of DRE at the surface increased (more cooling) from March to September. The high values of SAE over the western Himalayas gradually decreased towards the east, while DRE values are high over the eastern and central Himalayas compared to that of the western Himalaya. The reduction in SAE from spring to summer is associated with the rapid decrease in snow cover during summer due to the high solar insolation and aerosol-induced snow darkening. The high solar insolation at the surface (which is modulated by the cloud cover and aerosol loading in the atmosphere) directly melts the snow and enhances the SAE. For a constant snow albedo change of 0.05, the regional variation of solar flux at the surface would results in SAEs of 11.4, 13 and 11 Wm −2 over western, central and eastern Himalayas as respectively during spring. It should be noted that clear-sky DRE values are reported in this study. In the case of cloudysky conditions, DRE values decrease depending on the cloud fraction and hence result in lesser implications on the surface temperature and snow cover compared to SAE. Though the clear-sky DRE at surface is comparable to SAE over the Himalayas, the allsky DRE is smaller than the SAE due to the masking effects of clouds. Usha et al (2020) have reported that warming due to SAE is higher than cooling due to DRE over the western Himalayas during spring.
To investigate the effect of SAE on the simulation of snow cover over the Himalayas, the monthly evolution of surface albedo at the visible band (0.3-0.7 µm) simulated using RegCM4.6 is compared with MODIS observations during the melting period (March-September) (figure 3). The surface albedo decreases from winter to summer and from west to east. The control simulation overestimates the snow albedo during spring over central and eastern Himalayas and throughout the melting season over the western Himalayas. The snow cover reduced significantly due to melting during summer and simulated snow albedo matches with observation, especially over central and eastern Himalayas. The overestimation of the snow albedo by the control simulation is significantly reduced by the SAE over the entire Himalayan region during the spring season. It is interesting to note that, the model simulations are closer to the observations only when the SAE is included, which confirms that the inclusion of the aerosolinduced snow darkening in the climate models is very crucial for the proper representation of the radiation balance over the Himalayas. Hence, the SAE-induced snow melting is very important for the model to capture the reduction in the albedo during spring, especially over the central and eastern Himalayas. Though the positive bias in albedo still persists in both simulations over the western and eastern Himalayas during spring, the magnitude has decreased significantly with SAE compared to the control run. The inclusion of SAE reduced the model bias in albedo by 59% over the western Himalayas, 53% over the central Himalayas and 35% over the eastern Himalayas respectively.
Several global modelling studies on the aerosol induced snow albedo feedback especially due to BC deposition have suggested that the Himalayan-Tibetan region experiences the largest warming effects (Flanner et al 2007, Qian et al 2011, Yasunari et al 2015, Xu et al 2016, which is much higher than the cooling due to atmospheric aerosols (DRE) (Flanner et al 2009, Yasunari et al 2010, Ming et al 2012, Nair et al 2013, Usha et al 2020. Sarangi et al (2019) using WRF coupled with SNICAR estimated a mean surface radiative forcing of ∼40 Wm −2 due to SAE in spring over the Himalayas and suggested that the western Himalayas are more vulnerable to this phenomenon. As the intensity of solar insolation increases from winter to spring and summer, there exists an increase in the dust activity over the subtropics and its transport to the Himalayas and also the uplifting of BC to Himalayan foothills from the south Asian region (Cong et al 2015. It also coincides with the commencement of the seasonal melting of snow over this region which is very crucial for the hydrology of the surrounding regions. Overestimation of snow cover over the Himalayan-Tibetan region, especially over the western Himalayas, has been present in many climate models like community atmosphere model (CAM  (Das et al 2020). These authors attributed this bias to the (a) coarse resolution of the climate models which smoothen the complex terrain of Himalayan mountains, (b) excessive wintertime snowfall and (c) uncertainties in the snow cover retrieval algorithm of MODIS. Flanner et al (2009) showed that the inclusion of SAE within the climate models improved the warming trend over the springtime Eurasian region. Xu et al (2016) suggested that the proper incorporation of the aerosolinduced warming of the atmosphere, snow darkening as well as feedbacks are very important for the accurate estimation of the temperature trend over the cryosphere at high altitudes. This suggests a complex interplay between several feedback processes involved in the aerosol-induced SAE. The multiple feedback processes result in accelerated melting of snow and eventually the exposure of the underlying surface with much lower albedo, which causes further surface warming (Flanner et al 2007).

Effects of snow darkening on direct effects of aerosols
The DRE of aerosols has a strong dependence on surface albedo. So, to understand the effect of the decrease in snow albedo due to aerosol deposition on DRE, we have estimated the critical single scattering albedo. An aerosol system would have a warming (cooling) effect if the aerosol albedo is less (greater) than the critical single scattering albedo (Haywood and Shine 1995). For a reflecting surface like snow, the critical single scattering albedo is high and a small contribution of aerosol absorption can result in the warming of the earth-atmosphere system (Nair et al 2013). For the present study, critical single scattering albedo is estimated using surface albedo of AERO (DRE & SAE) and AeroATM (DRE) experiments (figure S3) following Haywood and Shine (1995). The effect of change in critical single scattering albedo due to snow darkening is very clear in the aerosol DRE at TOA as shown in figure 4. Over the western Himalayas, the DRE at TOA estimated without snow darkening effect (DRE no SAE ) is positive (warming) from March to September, whereas the magnitude of DRE decreases with the inclusion of SAE (DRE with SAE ) and the sign of DRE reverse from positive (warming) to negative (cooling) during July to September period. A similar phenomenon occurs over the central Himalayas as well, where the transition of radiative forcing at TOA from warming to cooling happen much early in spring ( figure 4(b)). Because of the low snow cover fraction and surface albedo over the central and eastern Himalayas, DRE of aerosols results in cooling at TOA (critical albedo < aerosol albedo) during most of the melting period, which is further strengthened (more cooling) by SAE. Over the eastern Himalayas, the DRE in presence of SAE (DRE with SAE ) is more negative from March to May compared to DRE no SAE simulation. Hence, the forcing due to aerosols at TOA is overestimated by 3.8, 2.4 and 1 Wm −2 over western, central and eastern Himalayan regions respectively when the aerosol-induced snow darkening is not considered in the simulation. The DRE no SAE and DRE with SAE are comparable (the difference between them is less than 1 Wm −2 ) over central and eastern Himalayas during summer monsoon due to the low SAE (figure 2), which is associated with lesser snow cover over these regions compared to the western Himalayas. The small changes in DRE due to SAE over the central and eastern Himalayas during the summer monsoon may not have much physical significance, since SAE itself is close to zero during this period. In addition, all the inherent uncertainties associated with climate modelling make the assessment of change in DRE due to SAE less reliable during the summer monsoon period over the eastern and central Himalayas.
The change in surface albedo due to aerosol deposition significantly affects the DRE at TOA, which was not considered in the earlier studies. The spatial variation of change in DRE at TOA due to SAE is shown in figure S5, where the DRE reduces by 6 Wm −2 over the western Himalayas, with larger reductions mostly over the low altitude regions during the spring. The atmospheric warming by the aerosols has not changed significantly (<1 Wm −2 ) due to snow darkening, showing more surface cooling because of DRE over a relatively darker surface (when the snow darkening is present). This SAE induced increase in surface cooling due to DRE partly compensate for the warming of the SAE. The quantification of the feedback processes associated with the SAE and its implications on DRE is challenging. To address this aspect, we have estimated the forcing efficiency of aerosols (DRE/AOD) for AERO and AeroATM simulations. The difference in forcing efficiency for DRE simulations with and without SAE are −27.5, −11 and −3.3 Wm −2 τ −1 over western, central and eastern Himalayas respectively during spring. This is mostly attributed to the snow darkening effects since the changes in the AOD due to meteorological feedbacks is less than 9% with respect to AeroATM experiment. The aerosol DRE of about 2 Wm −2 over the snow-covered surface for Hanle (western Himalayas) was reported by Nair et al (2013), which decreased to −1.5 Wm −2 over the non-snow-covered surface during spring. There have been several studies emphasizing the influence of the underlying surface on the DRE of aerosols at TOA (Haywood andShine 1995, Satheesh 2002). These authors reported that as surface albedo reduces, the DRE at TOA decreases and the surface experience more cooling.
The DRE and SAE of aerosols estimated using various datasets (in situ, satellite, and modelling) showed large variabilities due to the heterogeneities in aerosol sources and surface properties (albedo, snow cover, snow depth and snow ageing) (Kang et al 2020). Relative dominance of DRE or SAE on total aerosol forcing not only depends on the aerosol concentration (in the atmosphere or snow) but also on the surface and snow properties. This study demonstrates that the sign and magnitude of the DRE are closely associated with SAE. The western and central Himalayas are more vulnerable to SAE, whereas the eastern Himalayas is mostly dominated by the DRE of aerosols. The aerosol induced snow darkening over the snow-covered Himalayan region during the melting season (March-September) plays a crucial role in deciding the radiative perturbation caused by the aerosols. This study highlights the need for more field experiments and modelling efforts to understand the complex interaction and feedback of aerosols on regional climate over the high-mountain Himalayas.

Conclusions
The regional variability in the direct and SAE of aerosols within the high-altitude Himalayan mountains was studied using RegCM-4.6 coupled with aerosol and snow module. The important findings are (a) Western Himalayas experiences a large reduction in the snow albedo despite having lower AOD and BC mass concentration compared to central and eastern Himalayas. (b) Even though the absorption efficiency of airborne dust is low (compared to BC), its impact on snow darkening is comparable to BC due to the large mass loading of dust. (c) The positive bias between the simulated and satellite measured surface albedo reduces by 47%, 55% and 33% over western, central, and eastern Himalayas respectively when the aerosol induced snow albedo feedback is included in the model simulation. (d) The decline in surface albedo due to aerosol deposition causes a reversal (from warming to cooling) of the direct effect of aerosols at TOA over the western and central Himalayas. The forcing due to aerosols at TOA is overestimated by 3.8, 2.4 and 1 Wm −2 over western, central and eastern Himalayas respectively when SAE is not considered.
The results presented in this study need to be interpreted by keeping the following uncertainties and limitations into consideration, (a) The RegCM4.6 overestimate the snow cover fraction during the winter season, especially over the western Himalayas, which is a common bias seen in most of the climate models. Therefore, there exists an immediate need to improve the parameterization of the snow cover fraction over high altitude Himalayas for the accurate estimation of SAE. (b) The bias in the simulated precipitation alters the scavenging of the aerosols especially during the monsoon season, which leads to the bias in the AOD and the radiative forcing estimations. (c) The SNICAR module used in the present study assumes spherical snow grains and external mixing of aerosols and snow. Recent studies have shown that the inclusion of non-spherical grains decreases the SAE whereas the internal mixing of the impurities enhances SAE (Dang et al 2016, He et al 2017He et al , 2018. This could lead to uncertainties in the present SAE estimation.
The results of this study imply that the aerosol induced snow albedo feedback process plays an important role in radiation balance over the Himalayan region. It is crucial to incorporate snow darkening effects of aerosols in the climate models for the accurate estimation of the surface energy balance as well as aerosol-radiation interaction over the high altitude Himalayas. The implications of SAE on radiation balance and hydrological cycles need to be investigated further for the water security of the south Asian region.

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.