Investigation of Aerosol Climatology and Long-Range Transport of Aerosols over Pokhara, Nepal

This study presents the spectral monthly and seasonal variation of aerosol optical depth (τAOD), single scattering albedo (SSA), and aerosol absorption optical depth (AAOD) between 2010 and 2018 obtained from the Aerosol Robotic Network (AERONET) over Pokhara, Nepal. The analysis of these column-integrated aerosol optical data suggests significant monthly and seasonal variability of aerosol physical and optical properties. The pre-monsoon season (March to May) has the highest observed τAOD(0.75 ± 0.15), followed by winter (December to February, 0.47 ± 0.12), post-monsoon (October and November, 0.39 ± 0.08), and monsoon seasons (June to September, 0.27 ± 0.13), indicating seasonal aerosol loading over Pokhara. The variability of Ångström parameters, α, and β, were computed from the linear fit line in the logarithmic scale of spectral τAOD, and used to analyze the aerosol physical characteristics such as particle size and aerosol loading. The curvature of spectral τAOD, α’, computed from the second-order polynomial fit, reveals the domination by fine mode aerosol particles in the post-monsoon and winter seasons, with coarse mode dominating in monsoon, and both modes contributing in the pre-monsoon. Analysis of air mass back trajectories and observation of fire spots along with aerosol optical data and aerosol size spectra suggest the presence of mixed types of transboundary aerosols, such as biomass, urban-industrial, and dust aerosols in the atmospheric column over Pokhara.


Introduction
Atmospheric aerosols have a significant impact on the Earth's atmospheric radiation budget, due to their direct scattering and absorption characteristics, as well as an indirect impact on microphysics and clouds' formation [1][2][3]. In recent decades, there has been increasing concern about aerosols' impact on melting snow and ice in the high Himalaya and over the Tibetan Plateau [4][5][6], with indications that a significant portion of the aerosols arrived from the Indo-Gangetic Plains [7,8]. To date, though, there

Cimel Sun Photometer
An automatic sun-and-sky scanning Cimel Sun Photometer is located on the roof of the Shangrila Village Resort in Pokhara's south-western suburbs, as a part of AERONET. Descriptions of this network and methods for retrieving aerosol optical data have been published before [12,14,20].

Cimel Sun Photometer
An automatic sun-and-sky scanning Cimel Sun Photometer is located on the roof of the Shangrila Village Resort in Pokhara's south-western suburbs, as a part of AERONET. Descriptions of this network and methods for retrieving aerosol optical data have been published before [12,14,20].
The sun photometer measures direct sun radiances at 0.34 µm, 0.38 µm, 0.44 µm, 0.50 µm, 0.675 µm, 0.87 µm, 1.02 µm and 1.64 µm wavelengths during the daytime. The aerosol optical depth (τ AOD ) is derived by correcting attenuation due to Rayleigh scattering, absorption by ozone, and gaseous Atmosphere 2020, 11, 874 4 of 16 components in direct spectral measurements [12,14,15]. In addition, the sun photometer also measures diffuse sky radiance at four wavelengths 0.44 µm, 0.675 µm, 0.87 µm, and 1.02 µm [14,16,20]. These solar extinction measurements are used to retrieve aerosol columnar inversion products, such as volume size distribution and single scattering albedo. Three data quality levels, Level 1.0 (unscreened), Level 1.5 (cloud screened), and Level 2.0 (cloud-screened and quality-assured), are provided for analysis of aerosol data [16,20], with Level 2.0 data made available after the instrument is returned to NASA during annual swap outs. All products of AERONET are automatically computed and made available within the AERONET website [15]. The estimated uncertainty in computed AOD is reported to range from ±0.01 to ±0.02, which is spectrally dependent, and is found higher on the UV region [12,16,23].

Basic Equations and Definitions
The spectral τ AOD is guided by Ångström exponent (α), as given by power-law equation [34].
In logarithmic format can be written as: This gives the values of Ångström parameters α and β The Ångström exponent α can be further defined from the spectral AOD(λ) as The deviation of observed data from the linear fit line can be tested by using a second-order polynomial fit, with coefficients of the polynomial fit α 2 ,α 1 , and α 0 , as defined in previous studies [17][18][19].
The aerosol absorption characteristics can be exhibited by using SSA, aerosol absorption optical depth (AAOD) and AAE, and are related as below [2,3,35], and,

Tools Used for Backward Trajectories and Fire Spots
Four clusters were generated for each season by calculating five days backward air mass trajectories, starting at 500 m over the receptor site Pokhara (28.19 • N, 83.19 • E), for every day at 0:00, 6:00, 12:00, and 18:00 UTC, based on the HYSPLIT model [13]. A free software plugin called TrajStat was used from MeteoInfo for the calculations [36]. Fire spots were obtained from NASA's Fire Information for Resource Management System (FIRMS),which allows visualizing large scale bio-mass burning activities around the region. These active fire data were observed from NASA's Moderate Resolution Imaging Spectroradiometer (MODIS) and NASA's Visible Infrared Imaging Radiometer Suite (VIIRS), however, these satellite images will not provide the household level of biomass burning. A combination of satellite images of fire spots with air masses back trajectory will be presented to demonstrate a significant picture of aerosol climatology over the region [37][38][39].
In this study, we have used the data of the year 2017 to present a cluster analysis for indicating the seasonal aerosol sources' characterization over the observation site Pokhara.

Variability of Spectral Columnar AOD and Precipitable Water
In this section, we present the monthly and seasonal variability of spectral AOD and the temporal variability of AOD, by comparing with columnar precipitable water vapor (PW). Figure 2 shows the monthly mean spectral AOD values, τ AOD, over the 9 years from 2010 to 2018, at seven different wavelengths-0.34 µm, 0.38 µm, 0.44 µm, 0.50 µm, 0.675 µm, 0.87 µm, and 1.02 µm. Our statistical analysis excluded monthly averaged aerosol data from any months that had τ AOD data for less than ten days.
Atmosphere 2020, 11, x FOR PEER REVIEW 5 of 15 In this study, we have used the data of the year 2017 to present a cluster analysis for indicating the seasonal aerosol sources' characterization over the observation site Pokhara.

Variability of Spectral Columnar AOD and Precipitable Water
In this section, we present the monthly and seasonal variability of spectral AOD and the temporal variability of AOD, by comparing with columnar precipitable water vapor (PW). Figure 2 shows the monthly mean spectral AOD values, τAOD, over the 9 years from 2010 to 2018, at seven different wavelengths-0.34 μm, 0.38 μm, 0.44 μm, 0.50 μm, 0.675 μm, 0.87 μm, and 1.02 μm. Our statistical analysis excluded monthly averaged aerosol data from any months that had τAOD data for less than ten days. Due to the above constraint on data availability, the mean values for each calendar month across the nine years were computed from a varying number of months, as given in the parenthesis: January (6), February (8), March (8), April (7), May (8), June (5), July (4), August (3), September (4), October (6), November (6), and December (7). Level 2 aerosol data were very scarce on the AERONET website for the rainy season months of June, July, and August, due to cloud screening of the data and the requirement for the sun photometer to remain parked in a protected position whenever the rain sensor was wet.
The general trend of spectral variations shows that the τAOD is higher at a shorter wavelength and decreases at longer wavelengths. We found the spectral τAOD highest in April, followed by March, May, February, June (January, November), (December, October) September and (August and July). Months placed in parentheses have almost identical spectral τAOD at all wavelengths. In the monsoon season, the gradient of spectral τAOD in the longer wavelength decreases. Distinct features of seasonal τAOD can be observed in Figure 2, with the highest aerosol loading in the months of pre-monsoon season followed by winter, post-monsoon, and monsoon season. In the months of two seasons, winter and post-monsoon, the aerosol loadings are not significantly deviated from each other.
The temporal variation of light attenuation is also observed using a monthly mean AOD at 0.50 μm(τ0.50), along with a column-averaged precipitable water level (PW) in centimeters for different months and seasons ( Figure 3). The τ0.50, ranges from 0.31 to 0.56 in winter season with average value Due to the above constraint on data availability, the mean values for each calendar month across the nine years were computed from a varying number of months, as given in the parenthesis: January (6), February (8), March (8), April (7), May (8), June (5), July (4), August (3), September (4), October (6), November (6), and December (7). Level 2 aerosol data were very scarce on the AERONET website for the rainy season months of June, July, and August, due to cloud screening of the data and the requirement for the sun photometer to remain parked in a protected position whenever the rain sensor was wet.
The general trend of spectral variations shows that the τ AOD is higher at a shorter wavelength and decreases at longer wavelengths. We found the spectral τ AOD highest in April, followed by March, May, February, June (January, November), (December, October) September and (August and July). Months placed in parentheses have almost identical spectral τ AOD at all wavelengths. In the monsoon season, the gradient of spectral τ AOD in the longer wavelength decreases. Distinct features of seasonal τ AOD can be observed in Figure 2, with the highest aerosol loading in the months of pre-monsoon Atmosphere 2020, 11, 874 6 of 16 season followed by winter, post-monsoon, and monsoon season. In the months of two seasons, winter and post-monsoon, the aerosol loadings are not significantly deviated from each other.
The temporal variation of light attenuation is also observed using a monthly mean AOD at 0.50 µm(τ 0.50 ), along with a column-averaged precipitable water level (PW) in centimeters for different months and seasons ( Figure 3). The τ 0.50 , ranges from 0.31 to 0.56 in winter season with average value of 0.43 ± 0.12, 0.63 to 0.92 in pre-monsoon season with average value of 0.75 ± 0.15, 0.18 to 0.47 in monsoon season with average value of 0.27 ± 0.13, and 0.33 to 0.45, in post-monsoon season with average value of 0.39 ± 0.08. Variation of PW shows a similar trend of rainfall in Pokhara with an increase from months of winter to monsoon seasons, and then decreasing in post-monsoon season [40]. In monsoon season, the variation on τ 0.50 is significantly affected by the rainfall. 08. Variation of PW shows a similar trend of rainfall in Pokhara with an increase from months of winter to monsoon seasons, and then decreasing in post-monsoon season [40]. In monsoon season, the variation on τ0.50 is significantly affected by the rainfall. The variations of τ0.50 and PW show two different characters. Between December and April, these two parameters correlate significantly, with high R 2 (0.91) and low p-values (0.01), and while including data of May R 2 drops to 0.35 and p-value increases to 0.21. However, in monsoon and the post-monsoon months, they vary inversely with R 2 equal to 0.43 and p-value 0.15. Previous studies have shown that atmospheric water vapor can serve as a medium for igniting multiphase reaction to the formation of gas to particle transformation, and play a key role for hygroscopic growth of aerosols, which ultimately affects the aerosol optical properties [41,42]. In this study, we argue that τ0.50, which has a significant correlation with PW, from December to April, is mainly associated with the actual aerosol loading, and might not be linked with aerosol particles' hygroscopic effect growth. However, a chemical analysis of aerosols over the observation site, such as the observation of black carbon, dust, sulfate, and organic carbon concentrations due to their hygroscopic characters, can give a picture of the aerosol hygroscopic impact on τ0.50, which we have not done in this study. A study of aerosol loading based on the Ångström turbidity coefficient (β) is also presented in the next section, that will support the effect of aerosol loading for higher AOD. A previous study on the latitudinal variation of aerosol optical properties over the IGP region to the Central Himalayas during the premonsoon season has also explored the AOD due to aerosol loading over the observation site, rather than depicting hygroscopic growth of aerosol particles [43].

Angstrom Exponents, Curvature of AOD Spectra, SSA, AAOD and AAE
This section uses linear and second-order polynomial fit to study monthly and seasonal variation of spectral AOD in the wavelength ranges of 0.34 to 1.02 μm. Ångström parameters have been investigated in the past, by using spectral τAOD by distinguishing several AERONET sites with a variety of individual aerosol types such as biomass burning, urban and industrial, and desert dust [17]. Figure 4 shows monthly Ångström parameters (α and β) and curvature of the spectral AOD curve (α ). We observed that the average seasonal α as 1.14 ± 0.01 in winter, 1.17 ± 0.07 in the pre-monsoon season, 1.10 ± 0.08 in monsoon, and 1.22 ± 0.08 in the post-monsoon season, based on the monthly The variations of τ 0.50 and PW show two different characters. Between December and April, these two parameters correlate significantly, with high R 2 (0.91) and low p-values (0.01), and while including data of May R 2 drops to 0.35 and p-value increases to 0.21. However, in monsoon and the post-monsoon months, they vary inversely with R 2 equal to 0.43 and p-value 0.15. Previous studies have shown that atmospheric water vapor can serve as a medium for igniting multiphase reaction to the formation of gas to particle transformation, and play a key role for hygroscopic growth of aerosols, which ultimately affects the aerosol optical properties [41,42]. In this study, we argue that τ 0.50 , which has a significant correlation with PW, from December to April, is mainly associated with the actual aerosol loading, and might not be linked with aerosol particles' hygroscopic effect growth. However, a chemical analysis of aerosols over the observation site, such as the observation of black carbon, dust, sulfate, and organic carbon concentrations due to their hygroscopic characters, can give a picture of the aerosol hygroscopic impact on τ 0.50 , which we have not done in this study. A study of aerosol loading based on the Ångström turbidity coefficient (β) is also presented in the next section, that will support the effect of aerosol loading for higher AOD. A previous study on the latitudinal variation of aerosol optical properties over the IGP region to the Central Himalayas during the pre-monsoon season has also explored the AOD due to aerosol loading over the observation site, rather than depicting hygroscopic growth of aerosol particles [43].

Angstrom Exponents, Curvature of AOD Spectra, SSA, AAOD and AAE
This section uses linear and second-order polynomial fit to study monthly and seasonal variation of spectral AOD in the wavelength ranges of 0.34 to 1.02 µm. Ångström parameters have been investigated in the past, by using spectral τ AOD by distinguishing several AERONET sites with a variety of individual aerosol types such as biomass burning, urban and industrial, and desert dust [17]. Figure 4 shows monthly Ångström parameters (α and β) and curvature of the spectral AOD curve (α 2 ). We observed that the average seasonal α as 1.14 ± 0.01 in winter, 1.17 ± 0.07 in the pre-monsoon season, 1.10 ± 0.08 in monsoon, and 1.22 ± 0.08 in the post-monsoon season, based on the monthly mean data. Similarly, α (= −2α 2 ) were obtained 1.03 ± 0.23 for winter, 0.55 ± 0.14 for pre-monsoon, 0.23 ± 0.38 for monsoon, and 1.22 ± 0.08 for post-monsoon seasons. The seasonal turbidity parameters, β, were also found 0.18 ± 0.05 in winter, 0.31 ± 0.04 in pre-monsoon, 0.11 ± 0.06 in monsoon, and 0.16 ± 0.04 in post-monsoon seasons. Moreover, β gives a picture of aerosol loading with the higher the value, the higher the aerosol loading, and similarly, smaller β shows lower aerosol loading. We observe from Figure 4 that α shows more or less the same from one season to another, even though with a significant variation of τ 0.50 , by indicating the seasonal differences of the aerosol size spectrum.
Atmosphere 2020, 11, x FOR PEER REVIEW 7 of 15 were also found 0.18 ± 0.05 in winter, 0.31 ± 0.04 in pre-monsoon, 0.11 ± 0.06 in monsoon, and 0.16 ± 0.04 in post-monsoon seasons. Moreover, β gives a picture of aerosol loading with the higher the value, the higher the aerosol loading, and similarly, smaller β shows lower aerosol loading. We observe from Figure 4 that α shows more or less the same from one season to another, even though with a significant variation of τ0.50, by indicating the seasonal differences of the aerosol size spectrum. However, α can be seen to be the highest in post-monsoon, followed by pre-monsoon, winter, and monsoon. We found different values of α , even for a similar value of α, which indicates that the variation of aerosol microphysical properties can be presented by a climatological pattern of α rather than α. Figure 5a,b show that the significant variation on α are observed for the months even with the similar α. It was observed mainly in the months of a shift from one season to another season, in which prevailing air masses will also be in transition. Higher values of α (greater or close to 1) in the post-monsoon season followed by winter season provide a picture of columnar aerosol size distribution showing strong contribution by fine mode particles, which are mainly originated from anthropogenic, biomass burning, urban and industrial sources. The medium values (close to 0.5) in the pre-monsoon indicate a bimodal distribution of particles, and while the lowest values of α (close to 0.2) in monsoon season indicate a dominance by coarse mode particles.
The monthly mean variation of τ0.50 and β is shown in Figure 5a. We found that τ0.50 and β show a significant correlation (R 2 = 0.98 and p-value = 4.6 10 , scatter plot is not shown). This result shows that higher values of β correlating with τ0.50 was associated with aerosol loading, mainly dominated by fine mode aerosols. In July and August (peak monsoon) of monsoon, we observed However, α can be seen to be the highest in post-monsoon, followed by pre-monsoon, winter, and monsoon. We found different values of α , even for a similar value of α, which indicates that the variation of aerosol microphysical properties can be presented by a climatological pattern of α rather than α. Figure 5a,b show that the significant variation on α are observed for the months even with the similar α. It was observed mainly in the months of a shift from one season to another season, in which prevailing air masses will also be in transition. Higher values of α (greater or close to 1) in the post-monsoon season followed by winter season provide a picture of columnar aerosol size distribution showing strong contribution by fine mode particles, which are mainly originated from anthropogenic, biomass burning, urban and industrial sources. The medium values (close to 0.5) in the pre-monsoon indicate a bimodal distribution of particles, and while the lowest values of α (close to 0.2) in monsoon season indicate a dominance by coarse mode particles.
The monthly mean variation of τ 0.50 and β is shown in Figure 5a. We found that τ 0.50 and β show a significant correlation (R 2 = 0.98 and p-value = 4.6 ×10 −11 , scatter plot is not shown). This result shows that higher values of β correlating with τ 0.50 was associated with aerosol loading, mainly dominated by fine mode aerosols. In July and August (peak monsoon) of monsoon, we observed lowest β along with low τ 0.50 , and smallest values of α(∼ 1), and α (smallest, closest to 0.2) supports the contribution of coarse mode particles compared to fine mode particles on overall AOD. The significantly low τ 0.50 in monsoon compared to other seasons is associated with significant rainfall in this season [25,44].
We also analyzed α computed at different spectral bands, as shown in Figure 5b to identify the aerosol types. The difference in α at spectral bands (α (0.675−0.87) − α (0.34−0.38) ) was significantly higher in post-monsoon and winter seasons, indicating a size distribution dominated by fine modes (Figure 5b). The difference lowers in months of pre-monsoon seasons, supporting the company of coarse mode as well. This difference is negative in the monsoon season, which indicated domination by coarse mode aerosols [17,19].
Atmosphere 2020, 11, x FOR PEER REVIEW 8 of 15 We also analyzed α computed at different spectral bands, as shown in Figure 5b to identify the aerosol types. The difference in α at spectral bands (α ( . . ) − α ( . . ) ) was significantly higher in post-monsoon and winter seasons, indicating a size distribution dominated by fine modes (Figure 5b). The difference lowers in months of pre-monsoon seasons, supporting the company of coarse mode as well. This difference is negative in the monsoon season, which indicated domination by coarse mode aerosols [17,19]. The aspect of change in aerosol size spectrum, based on Angstrom parameters, is also examined by using aerosol volume size distribution data obtained from AERONET over Pokhara site. July to September data are not analyzed due to the limitation of availability of data in AERONET. Figure 6 shows the monthly mean of the aerosol volume size distribution for different seasons, along with α and α, and reveals a bimodal structure of aerosol sizes for each month. It is observed from Figure 6 that the volume size distribution for the months are significantly dominated by accumulation mode of aerosols (with higher α ), both modes of aerosol sizes (with the medium α ), and coarse modes (with smaller α ). For this study, July to September volume size distribution data are not available for comparison. The aspect of change in aerosol size spectrum, based on Angstrom parameters, is also examined by using aerosol volume size distribution data obtained from AERONET over Pokhara site. July to September data are not analyzed due to the limitation of availability of data in AERONET. Figure 6 shows the monthly mean of the aerosol volume size distribution for different seasons, along with α and α, and reveals a bimodal structure of aerosol sizes for each month. It is observed from Figure 6 that the volume size distribution for the months are significantly dominated by accumulation mode of aerosols (with higher α ), both modes of aerosol sizes (with the medium α ), and coarse modes (with smaller α ). For this study, July to September volume size distribution data are not available for comparison.   Figure 7, bottom two figures, shows AAOD spectra. Previous studies reported that SSA spectra of different AERONET locations with dust containing aerosols have increased SSA with increasing wavelength, while locations dominated by urban industrial or biomass burning decrease with increasing wavelength [2,4,45,46]. This trend was also used in the past studies to differentiate between carbonaceous aerosols and dust in different locations [2,10,47,48].
The behaviors of increasing of SSA spectra up to 0.675 μm resemble aerosol components from dust aerosols, and decreasing from 0.675 μm to 1.02 μm shows similar to that of aerosol products of biomass burning, urban and industrial activities. Similar results can also be found in previous studies to define a mixed type of aerosols [45,46].   Figure 7, bottom two figures, shows AAOD spectra. Previous studies reported that SSA spectra of different AERONET locations with dust containing aerosols have increased SSA with increasing wavelength, while locations dominated by urban industrial or biomass burning decrease with increasing wavelength [2,4,45,46]. This trend was also used in the past studies to differentiate between carbonaceous aerosols and dust in different locations [2,10,47,48].
The behaviors of increasing of SSA spectra up to 0.675 µm resemble aerosol components from dust aerosols, and decreasing from 0.675 µm to 1.02 µm shows similar to that of aerosol products of biomass burning, urban and industrial activities. Similar results can also be found in previous studies to define a mixed type of aerosols [45,46]. The spectral dependence of AAOD was also used to compute the absorption Ångström exponent (AAE), using a linear regression fit on the logarithmic scale plot [2,45,46]. AAE data were found for January (1.  [2,[49][50][51]. These previous investigations have shown that AAE values vary from 1.2 to 3 for dust, 0.75 to 1.3 for urban and industrial aerosols, and 1.2 to 2 for biomass burning. Similarly, in different studies on the AERONET sites dominated by an optical mixture of smoke, dust, and industrial and urban pollution, have reported AAE in the ranges of 1.2 to 1.8 [52]. AAE observed in this study also lies in the range, indicating the absorbing behavior of aerosol components obtained in the mixed type of aerosols. Figure 8 shows the seasonal clusters of five days air mass back trajectories arriving at the observation site Pokhara at an altitude of 500 m and active fire spots (red dots) for 2017. The percentage contribution of each cluster is also shown in the figure for each season. The air masses reaching Pokhara valley follow two distinctive pathways from the Indo Gangetic Plain (IGP) region, during winter and pre-monsoon (Western Nepal, West India, and Pakistan), when the influence of strong western disturbances occur. During monsoon, majority of air masses arrived from the eastern IGP region and Bay of Bengal. Dense fire spots are observed during winter and pre-monsoon period over the region, which can enhance the emission of aerosols from biomass burning, that could be transported to Pokhara and influence the air quality and enhance AOD during that period [53]. In addition, during the post-monsoon, widespread crop residue burning occurs in the north-west part of India [54]. Comparatively few active fire spots are detected during the monsoon season, which might be due to cloud cover and heavy rainfall in South Asia, caused by moist air from the Arabian Sea and Bay of Bengal. The spectral dependence of AAOD was also used to compute the absorption Ångström exponent (AAE), using a linear regression fit on the logarithmic scale plot [2,45,46]. AAE data were found for January (1.71 ± 0.04), February (1.50 ± 0.03), March (1.42 ± 0.03), April (1.38 ± 0.03) and May (1.47 ± 0.03), with ranges 1.38 to 1.71. Various studies have been done in the IGP region and different locations of AERONET to investigate the AAE based on types of aerosols [2,[49][50][51]. These previous investigations have shown that AAE values vary from 1.2 to 3 for dust, 0.75 to 1.3 for urban and industrial aerosols, and 1.2 to 2 for biomass burning. Similarly, in different studies on the AERONET sites dominated by an optical mixture of smoke, dust, and industrial and urban pollution, have reported AAE in the ranges of 1.2 to 1.8 [52]. AAE observed in this study also lies in the range, indicating the absorbing behavior of aerosol components obtained in the mixed type of aerosols. Figure 8 shows the seasonal clusters of five days air mass back trajectories arriving at the observation site Pokhara at an altitude of 500 m and active fire spots (red dots) for 2017. The percentage contribution of each cluster is also shown in the figure for each season. The air masses reaching Pokhara valley follow two distinctive pathways from the Indo Gangetic Plain (IGP) region, during winter and pre-monsoon (Western Nepal, West India, and Pakistan), when the influence of strong western disturbances occur. During monsoon, majority of air masses arrived from the eastern IGP region and Bay of Bengal. Dense fire spots are observed during winter and pre-monsoon period over the region, which can enhance the emission of aerosols from biomass burning, that could be transported to Pokhara and influence the air quality and enhance AOD during that period [53]. In addition, during the post-monsoon, widespread crop residue burning occurs in the north-west part of India [54]. Comparatively few active fire spots are detected during the monsoon season, which might be due to cloud cover and heavy rainfall in South Asia, caused by moist air from the Arabian Sea and Bay of Bengal. It is evident from previous analyses of aerosol chemical composition over the Himalayan foothill in Nepal that the emission from crop-residue burning over the IGP has a significant impact on the air quality over the regions [24,28,46,[55][56][57][58][59]. The cluster analysis of trajectories in this study also indicates the high likelihood that the aerosol population's physical and optical properties over the Pokhara valley could be influenced significantly by the regional transport of air masses from polluted regions of South Asia. This could be further confirmed by looking at MODIS visible imagery, showing continuous aerosol haze layers extending up from the IGP into the Himalayan valleys (Figure 9a,b). The MODIS satellite image taken on 27 October 2017 reasonably shows emissions from biomass burning. The biomass burning smoke funnels through the densely populated and industrialized areas on IGP, and after mixing with anthropogenic pollution background, it comfortably transports to Nepal's Pokhara valley, and is also suggested by the HYSPLIT back trajectory analysis (Figure 9b).  It is evident from previous analyses of aerosol chemical composition over the Himalayan foothill in Nepal that the emission from crop-residue burning over the IGP has a significant impact on the air quality over the regions [24,28,46,[55][56][57][58][59]. The cluster analysis of trajectories in this study also indicates the high likelihood that the aerosol population's physical and optical properties over the Pokhara valley could be influenced significantly by the regional transport of air masses from polluted regions of South Asia. This could be further confirmed by looking at MODIS visible imagery, showing continuous aerosol haze layers extending up from the IGP into the Himalayan valleys (Figure 9a,b). The MODIS satellite image taken on 27 October 2017 reasonably shows emissions from biomass burning. The biomass burning smoke funnels through the densely populated and industrialized areas on IGP, and after mixing with anthropogenic pollution background, it comfortably transports to Nepal's Pokhara valley, and is also suggested by the HYSPLIT back trajectory analysis (Figure 9b).

Investigation of Aerosol Sources and Types
of South Asia. This could be further confirmed by looking at MODIS visible imagery, showing continuous aerosol haze layers extending up from the IGP into the Himalayan valleys (Figure 9a,b). The MODIS satellite image taken on 27 October 2017 reasonably shows emissions from biomass burning. The biomass burning smoke funnels through the densely populated and industrialized areas on IGP, and after mixing with anthropogenic pollution background, it comfortably transports to Nepal's Pokhara valley, and is also suggested by the HYSPLIT back trajectory analysis (Figure 9b).

Conclusions
We have been able to study a cluster analysis of the aerosol climatology over Pokhara, an AERONET site located on a Himalayan foothill, based on long-term aerosol optical properties, size spectra and regional aerosol sources, and the long-range transport of aerosols over the observation site. The variation of AOD spectra and the magnitude of AOD are strongly associated with the change in seasons that brings different air masses over the observation site.
The maximum τ 0.50 in the pre-monsoon season (0.75 ± 0.15), followed by winter (0.43 ± 0.12), post-monsoon (0.39 ± 0.08), and monsoon (0.27 ± 0.13) seasons, show different aerosol loading, by a varying amount of wind-driven long-range transport of aerosols (mainly dust and aerosols from biomass burning, urban-industrial activities). The strong correlation between τ 0.50 and β also shows the association of seasonal aerosol loadings to influence the τ 0.50 . The air masses back trajectory and local and regional fire spots were supported, to investigate the effect of transboundary air pollution, from different sectors over Pokhara. We find that Pokhara receives dominating westerly air masses during the post-monsoon season, indicating a significant effect of biomass burning, such as crop residue burning over the northeast part of Pakistan. Both westerly and southwesterly air masses crossing over the Thar desert of India affect the area during the pre-monsoon and winter seasons. The pre-monsoon season (April and May) is mainly crop harvesting time, such as wheat, in the northwest area of Pakistan, and the smoke of burning crop residue easily diffuses to the atmosphere to be transported to Pokhara. The monsoon season is significantly associated with the transport of aerosols from the side of the Bay of Bengal. The trajectory analysis provides a unique background knowledge of aerosol components, which arrives over Pokhara by mixing with anthropogenic pollution background and natural aerosols, while passing over the heavily polluted IGP or coming from the IGP.
The study of spectral variation of AOD by using α ,the first derivative of the angstrom exponent (α) with the wavelength in log scale distinguishes the aerosol size distribution, even with a similar value of α. Furthermore, α is found to be more sensitive for months when there is a transition from one season to another. This investigation suggests that the analysis of α gives a more sensible complement of α to characterize more fully wavelength dependence of AOD and the comparative influence of aerosol size spectra of two modes in aerosol loading. Intermediate AAE values range from 1.38 to 1.78, and the increasing and decreasing nature of SSA with wavelengths (analyzed only for winter and pre-monsoon seasons due to limitations of availability of data) proved beneficial for identifying the aerosols over Pokhara as a mixed type of aerosols, such as dust, aerosols from biomass burning and urban-industrial activities.