Giant dust particles at Nevado Illimani: a proxy of summertime deep convection over the Bolivian Altiplano

A deeper understanding of past atmospheric circulation variability in the Central Andes is a high-priority topic in paleoclimatology, mainly because of the necessity to validate climate models used to predict future precipitation trends and to develop mitigation and/or adaptation strategies for future climate change scenarios in this region. Within this context, we here investigate an 18-years firn core drilled at the Nevado Illimani in order to interpret its mineral dust record in relation to seasonal processes, in particular atmospheric circulation and deep convection. The core was dated by annual layer counting based on seasonal oscillations of dust, calcium and stable isotopes. Geochemical and mineralogical data show that dust is regionally-sourced in winter and summer. During austral summer (wet season) an increase in the relative proportion of giant dust particles (ø > 20 μm) is observed, in association to oscillations of stable isotope records (δD, δD, δD, δ O). It seems that at Nevado Illimani both the deposition of dust and the isotopic signature of precipitation are influenced by atmospheric deep convection, which is also related to the total amount of precipitation in the area. This hypothesis is corroborated by local meteorological data. The interpretation of giant particle and stable isotope records suggests that during La Niña years, summer convection activity is enhanced, in agreement with atmospheric circulation studies. Giant particles and stable isotopes, when considered together, can be therefore used as a new proxy for obtaining information about deep convective activity in the past, which is ultimately related to paleo-ENSO variability.

proxies. Dust particles entrapped in firn samples seem to originate from regional sources during both winter and summer, despite minor mineralogical differences between the two seasons are observed. An interesting result concerns the presence of giant dust particles (presenting a diameter larger than 20 µm), whose variability is correlated to the stable isotope record. Our data and local meteorological observations suggest that this relationship could be explained as a consequence of atmospheric deep convection over the Bolivian Altiplano during summer. This study shows for the first time that the presence of giant dust particles in Andean firn and ice is strictly controlled by climatic processes. We found clear evidence that the convective activity over the Altiplano, reconstructed through the analysis of giant particles, is enhanced during summer periods and in particular during La Niña years, in agreement with observations concerning atmospheric circulation anomalies in the area (Vuille, 1999). From this perspective, this study demonstrates the great potential of giant particles records, which is strongly influenced by climatic and meteorological processes at high altitude continental glaciers. This a first exploratory work, analysis of a longer ice core would be desirable in the future to investigate the relationships between giant dust particles deposition, atmospheric deep convection and periodic climatic phenomena (ENSO).
Nevado Illimani consists of a granodiorite pluton of Late Oligocene age, with a short belt composed by a coeval dacitic flow located near the south-western border of the pluton (McBride et al., 1983;Jiménez and López-Velásquez, 2008). In June 2017, a 23.8 m firn core (corresponding to 13.75 m water equivalent) was drilled at an altitude of 6350 m a.s.l. on the saddle between the two Nevado Illimani summits, approximately where two deep ice cores were recovered in June 1999 (Knüsel et al., 2003). The expedition was coordinated by a French, Russian, Bolivian and Brazilian team and integrated the Ice Memory project (Université Grenoble Alpes Foundation). An EM-100-1000 electromechanical ice core drill (Cryosphere Research Solutions, Columbus, Ohio, USA) was used for the drilling and three cores were extracted, two down the bedrock (136 and 134 m) and the core for this study (23.8 m).
The cores (diameter of 10 cm), consisting of 24 sections of approximately 1 m length, were transported by mountain porters from the drilling site to the base camp during the night in order to prevent melting. Once at the base camp, they were immediately transported to a refrigerated container located in La Paz, where temperature was set at -20 °C. After the drilling campaign, the container was shipped to the Institut des Géosciences de l'Environnement (IGE, Université Grenoble Alpes, France) where the core sections were weighted and cut longitudinally using a vertical band saw in a cold room (at -15 °C).
Ice layers (less than 5 cm thick) occurred frequently along the firn core (Fig. S1), generally in closely spaced groups of 2 or 3 individual layers, indicating occasional melting events at the Nevado Illimani site. In addition, less than 1 cm thick ice layers (or possibly wind crusts, hardly distinguishable from the former at greater depths) commonly occurred along the core. One quarter of the original core was dedicated to dust analyses and transported for this purpose to the EUROCOLD facility of the University of Milano-Bicocca (Italy). There, firn sections were transversely cut at 5 cm using a horizontal band saw with a cobalt steel blade and 464 samples were obtained. These were manually decontaminated by mechanical scraping with a clean ceramic knife inside a laminar flow high-efficiency particle air (HEPA) ISO 5 Class bench located in an ISO 6 Class cold room. Once decontaminated, the samples were put into clean Corning centrifuge tubes and kept frozen until the measurements.

Coulter Counter Analysis
Samples were melted at room temperature, and a ~10 mL aliquot from each was transferred to an Accuvette Beckman Coulter vial, previously washed with Millipore Q-POD ® Element ultra-pure water (in an ISO 5 Class laminar flow bench located inside an ISO 6 Class clean laboratory). Each sample was treated following standardized protocols (Delmonte et al., 2002). A Beckman Coulter Multisizer 4 equipped with a 100 µm orifice was used to measure dust concentration and grain size (400 size channels within the 2-60 µm interval of spherical equivalent diameter). Each sample was measured twice (consuming 0.5 mL per measurement), the mean standard deviation between the two measurements being lower than 3%.
Systematic analysis of ultra-pure water blanks allows estimating a mean signal to noise ratio around 97. Particle mass and number as well as the proportion (%) of giant particles having diameter >20 µm (coarse silt), were calculated from the particle volume and number size distributions (GPPms and GPPnb, respectively) assuming particle sphericity and a mean density of 2.5 g cm -3 .

Instrumental Neutron Activation Analysis
A set of 10 samples was dedicated to Instrumental Neutron Activation Analysis (INAA). Samples were selected from different depth intervals along the core (see Table 1, "N" series, and Table S1 for precise depths), in order to be representative of both the dry and the wet seasons. The samples were filtered using PTFE Millipore © membranes (0.45 µm pore size, 11.3 mm diameter) previously rinsed in an ultra-pure 5% solution of bi-distilled HNO3-according to the procedures adopted by Baccolo et al. (2015). The filtration took place in an ISO 5 Class laminar flow bench. Two blank membranes, that underwent the same cleaning procedures, were prepared by filtering 300 mL of MilliQ (Millipore © ) water.
For calibration and quality control, we used certified solid standards: USGS AGV2 (ground andesite), USGS BCR2 (ground basalt), NIST SRM 2709a (San Joaquin soil) and NIST SRM 2710a (Montana soil). In addition, standard acid solutions for each analyzed element were prepared with concentrations in the order of µg g -1 . Blanks for the empty flask and for the ultrapure acid solution used to prepare the liquid standards were also measured. Samples, standards and blanks were irradiated at the Applied Nuclear Energy Laboratory (LENA, University of Pavia, Italy) by a Triga Mark II reactor of 250 kW. The "Lazy Susan" channel, neutron flux equal to 2.40 ± 0.24 x 10 12 s -1 cm -2 , was used to identify Ce, Cs, Eu, Hf, La, Sc, Sm, Th, and Yb. Samples were successively transferred at the Radioactivity Laboratory of the University of Milano-Bicocca, in order to 4 100 105 110 acquire gamma spectra by means of a high-purity Germanium detector HpGe (ORTEC, GWL series), following the standardized procedure developed for low-background INAA (Baccolo et al., 2016).
The masses of the elements in each sample were determined by comparing spectra related to standards and samples (Baccolo et al., 2016). In order to compare different spectra, the time of acquisition, the radionuclide decay constant, the cooling time and a factor considering radioactive decay during the acquisition were kept into account. The detection limits were calculated considering three times the standard deviation of the blank signal. The uncertainties for each element were calculated based on the mass measurements, the adjustment for the spectrum, the subtraction of the blanks and the standard concentration uncertainties. Errors for the elemental concentrations in our samples ranged from 3% for La to 17% for Cs, and the detection limits ranged from 0.1 µg g -1 for Sm to 7 µg g -1 for Ce (Table S2). Full analytical details can be found in Baccolo et al. (2016).
The Enrichment Factor (EF) normalization was calculated for each element considering as a reference the mean composition of the upper continental crust (UCC) (Rudnick and Gao, 2003). Scandium (Sc) was chosen as the crustal reference element following Eq. (1): Scandium was chosen as reference element because it is poorly affected by processes altering its mobility in hosting minerals, and its biogeochemical cycle is almost unaffected by anthropogenic activities (Sen and Peucker-Ehrenbrink, 2012).
In addition, Sc is highly correlated with other lithogenic elements, such as Ce (r = 0.997), used by Eichler et al. (2015) as a crustal reference for Nevado Illimani samples, and La (r = 0.989). The choice of Sc has been also determined by its easy and precise determination through INAA.

Micro-Raman Spectroscopy
We used single-grain Raman spectroscopy to identify mineralogy of dust particles having a diameter smaller than 5 µm. This because this kind of analysis was carried out for provenance purposes, thus excluding the large and the giant particles of obvious local origin. A set of four samples (see Table 1, "R" series, and Table S1) was prepared, following the procedure described in previous studies (Delmonte et al., 2017;Paleari et al., 2019) specifically developed for small dust grains. Two samples are representative of mineral dust deposited in the dry season (high dust concentration) whereas two represent dust from the wet season (low dust concentration, or "background"). Measurements were performed by using an InVia Renishaw micro-Raman spectrometer (Nd YAG laser source, λ = 532 nm) available at the Laboratory for Provenance Studies (UNIMIB). We identified the mineralogy of more than 630 grains, excluding organic particles possibly related to contamination and particles with an undetermined spectrum or with no signal.

Stable Isotope and Ion Chromatography Analysis
The dust analyses described above used one quarter of the longitudinally-cut firn core. A second quarter was shipped in a frozen state to the Climate Change Institute (CCI, University of Maine, USA) for ion chromatography (IC) and water stable isotope analysis.
At the CCI, in a cold room set at -20 °C, we cut longitudinal sections of the core with a vertical band saw to separate an inner and an outer part. The latter was sampled by transverse cuts approximately every 12 cm using a stainless-steel hand saw (resulting in 190 samples) and stored in plastic bottles for stable isotope ratio determination. Decontamination of the inner part was performed by scraping with a clean ceramic knife under a laminar flow HEPA bench inside the cold room. Then, the decontaminated inner part was sampled for IC analysis by a continuous melter system (Osterberg et al., 2006) also in an ISO 6 Class clean room. The mean sample resolution was 3 cm, resulting in 767 samples. We measured Ca² + concentration using a ThermoScientific TM Dionex TM Ion Chromatograph ICS-6000 analytical system fitted with suppressed conductivity detectors with a Dionex AS-HV autosampler. The method detection limit (MDL) was defined as 3 times the standard deviation of the blank samples (MilliQ water, 10 blank samples). The detection limit for Ca 2+ was 21.05 µg L -1 .
The D and the δD, δ 18 O were determined using a Picarro L2130-i wavelength-scanned cavity ring-down spectroscopy instrument (Picarro Inc., USA) with a precision of 0.1 ‰.

Seasonal variability of proxies and firn core chronology
We established a chronology for the Nevado Illimani firn core based on annual layer counting (ALC), considering the pronounced seasonal oscillation of dust concentration, calcium and water stable isotopes (Fig. 2). Dust concentration variations, which are recognized for being useful for ALC in tropical and continental ice cores Kutuzov et al., 2019), span about two-orders of magnitude between the summer and the winter. Dust concentration varies from ~850 ppb during the wetter season, to ~4,800 ppb during the dryer season (median values). When considering extreme values, the variation range exceeds three orders of magnitude, being the lowest concentration during the wet season 70 ppb and the highest one during the dry season 20-30 ppm. The well-defined oscillatory pattern of dust concentration variability reflects the extreme seasonality of precipitation over the Altiplano and the succession of dry and wet conditions. Sublimation has a limited influence to this seasonality (Ginot et al., 2002). Dust concentration is in accordance with the Ca 2+ record and in agreement with literature studies (Knüsel et al., 2003). Ionic calcium can be primarily associated to calcium sulphate (CaSO4) and/or calcium carbonate (CaCO3) (Kutuzov et al., 2019). Because scarcity of calcium carbonates was revealed by mineralogical analyses (Fig. 3, see below), we argue that most of the ionic calcium observed in firn samples is present as a soluble species, probably CaSO4, and not detectable through Raman spectroscopy on single insoluble particles. However, we 6 160 165 170 175 180 cannot exclude that Ca-bearing aerosols might have been initially a mixture of pure gypsum and calcium carbonates that successively reacted with atmospheric H2SO4 (Röthlisberger et al., 2000).
The regular succession of dry dusty periods and wet periods can be associated with the seasonal onset and decay of the Bolivian High, a high pressure system which is well developed and centered over Bolivia (Lenters and Cook, 1997). When the Bolivian High is particularly strong and displaced southward of its climatological position, easterly flow in the high troposphere is enhanced as well as moisture advection from the interior of the continent to the Altiplano. This moisture transport from the Amazon Basin toward the Altiplano induces a notable amount of precipitation over the Altiplano (wet season), associated with strong summer convection. The relatively low dust concentration found in the Illimani snow during the summer period is therefore related to particle dilution in the snowpack because of increased precipitation and reduced regional dust mobilization deriving from wetter soil conditions. Conversely, during winter (JJA) months, conditions over the Altiplano are typically dry, leading to higher dust availability. At that time of the year, the winter westerly flow over the entire region promotes eastward dust transport towards the Nevado Illimani, leading to massive dust deposition in the firn layers representing -on average -about 85% of the total annual dust mass there deposited.
Seasonal variations of the water stable isotopes in snow precipitated over the Andes are also useful for dating. However, the Andean isotopic signal led to divergent interpretations (Vimeux et al., 2009). Whereas in polar ice cores the water isotopic signature is chiefly related to temperature (Uemura et al., 2012), the isotopic composition of tropical precipitation can be affected by a larger number of factors . It is well known that the so called "amount effect" leads to an anti-correlation between the amount of precipitation and the proportion of heavier isotopes in the precipitation. This effect is in turn related to an ensemble of physical and microphysical processes producing a robust signal on the isotopic composition of precipitation (Dansgaard, 1964;Risi et al., 2008;Vuille et al., 2003). In this context, deep atmospheric convection also plays a role on stable isotope composition (Vimeux et al., 2005). Modelling studies (e.g. Bony et al., 2008;Risi et al., 2008) reveal that the stronger the convective activity during a particular event, the higher the total amount of precipitation and thus the more depleted the isotopic composition of precipitation. Satellite data (Samuels-Crow et al., 2014) reveal that during the summer season the isotopic composition of water vapor strongly depends on convective activity. These observations lead the authors to conclude that the isotopic composition of snow from the tropical Andes mainly reflects tropical convection.
Convective precipitation over the Bolivian Altiplano is enhanced during the wet summer season, leading to the emergence of clear seasonal oscillations in the stable isotope records of Nevado Illimani firn core which can be used for ALC and to develop a chronology.
Interestingly, we note a close correspondence between the variability of stable isotopes and the proportion of giant particles in firn (Fig. 2a): oscillations of the stable isotope record (δD, δ D in Fig. 2a)  Considering the pronounced seasonal changes in dust concentration, Ca 2+ and stable isotopes, it is possible to assign to the base of the Illimani firn core an age corresponding to the beginning of 1999 AD. The firn record thus covers the 18-years period from early 2017 to early 1999 and the average accumulation rate can be estimated in the order of approximately 750 mm of water equivalent per year, slightly higher than the one inferred by Knüsel et al. (2003).

Dust provenance: mineralogy and geochemistry
The mineralogical composition of fine dust (<5 µm) deposited onto the Illimani firn layers reveals that the most abundant mineral phases are quartz, feldspars (alkali feldspars and plagioclase) and phyllosilicates (Fig. 3, Table S3). Phyllosilicates are mainly represented by muscovite-illite (and/or smectite, hardly distinguishable from illite by their Raman spectra) and secondarily by kaolinite (representing 7.5% in the dry season and 1.5% in the wet season). Altogether, quartz, feldspars and phyllosilicates account for 75-78% of mineral particles during both the wet and the dry season, but phyllosilicates are particularly abundant during the wet season, when they represent approximately 44% of minerals (Fig. 3a). We believe that the increased abundance of muscovite-illite during the wet season is related to different depositional regimes. During the dry winter, aerodynamic platy-like phyllosilicates can remain in the atmosphere for longer periods and are only partially deposited. On the contrary, during the wet summer strong scavenging is associated to heavy precipitations at Nevado Illimani (Bonnaveira, 2004), enhancing the removal of mineral particles from the atmosphere, including phyllosilicates.
Titanium oxides and iron oxides/hydroxides are present in all samples. Hematite is twice as abundance as goethite. This is typical in regions dominated by arid conditions or where a prolonged warm and dry season is followed by a shorter and wetter period (Journet et al., 2014). Accessory minerals include carbonate and tourmaline, and very rare pyroxenes (Table   S3). Such a mineralogical composition is coherent with the felsic to intermediate plutonic volcanic source rocks, suggesting that most of the dust deposited at Nevado Illimani has a local/regional provenance, both in the wet and in the dry seasons.
Low-background INAA analyses allowed determining the EFs for different rare earth elements (REEs), which are nonmobile and therefore widely used as provenance tracers (McLennan, 1989;Moreno et al., 2006;Gabrielli et al., 2010). In Fig. 4a the Yb/La and the Eu/Sm elemental ratios are used to compare dust samples retrieved from Nevado Illimani firn samples and literature data concerning geological samples from the Altiplano-Puna Volcanic Complex (Lindsay et al., 2001;Ort et al., 1996) and potential source areas (PSAs) in South America (Gaiero et al., 2004(Gaiero et al., , 2013. The Yb/La ratio can be used to appreciate wether heavy and light REEs are enriched or depleted with respect to each other, whereas the Eu/Sm ratio is a proxy for the europium anomaly, usually calculated considering Gd, not detected by our analytical method. The comparison reveals that Nevado Illimani dust has a composition similar to APVC crystal-rich ignimbrites, pointing to a correspondence with samples from the Northern Puna region, and not with samples from the salt lakes present in the Altiplano (Uyuni and Coipasa salars). These pieces of evidence agree with previous analyses of strontium and neodymium isotopes on Nevado Illimani ice core dust  and with the geochemical signature of local sources in the Altiplano (Gili et al., 2017), supporting the hypothesis that dust deposited at Nevado Illimani is sourced by sediments present in Southern Altiplano-Northen Puna areas.
The EFs of Ce, La, Sm, Eu and Yb is similar between samples with higher percentage of giant particles (GPPms > 35%, characteristic from the wet season) and samples with background GPPms (Fig. 4b), supporting what has been observed in relation to the mineralogical composition of samples from the wet and dry seasons. However, two important exceptions to this pattern occur in relation to Hf and Cs. High GPPms samples show an anomalous Hf enrichment when compared to background samples, a feature that was observed also in Saharan dust samples (Castillo et al., 2008) and attributed to the presence of detrital zircon ((Zr, Hf)SiO4) in samples showing the coarsest grain size. As discussed by Vlastelic et al. (2015) the Hf enrichment observed in samples where the percentage of giant particles is high might be related to the presence of a few silt-sized zircon grains, which would in turn require high energy (turbulence) to lift them and keep them in suspension into the atmosphere. In this study, however, zircon grains were not detected by Raman Spectroscopy, which indicated greater abundance of phyllosilicates in dust deposited during the wet summer season. The Hf enrichment in these samples may thus be related instead to tiny zircon inclusions within phyllosilicate particles.
Another deviation from background is related to cesium. Because of the very common Cs-K exchange into interlayer sites in illite-muscovite, the Cs enrichment in samples from the wet season may also be related to the greater abundance of phyllosilicates (Cremers et al., 1988;Rosso et al., 2001).

Relationship between the giant particles and deep convection
The absolute concentration of dust in firn and ice cores depends on many factors including the snow accumulation rate, the dust source strength (which includes soil aridity/wetness, vegetation cover and any other factor influencing the quantity of particles available for deflation), as well as transport processes, which also affect the residence time of particles in the atmosphere (mainly in the case of long-range transport). Conversely, the dust size distribution and the relative proportion of particles within a given grain size merely depends on transport conditions (Delmonte et al., 2004(Delmonte et al., , 2017. In the case of Nevado Illimani, where dust is mostly regionally-sourced, we believe that the dust concentration is primarily modulated by the seasonally-varying source strength, mainly depending on source aridity and humidity, and by accumulation rate. Conversely, the dust size distribution, in particular the relative mass (or number) of giant particles is related to the relative strength of local turbulence and deep convection. Indeed, we observe an interesting correlation (Fig. 2a) between the relative mass of giant particle percentage (GPPms) and δD, both related to convective activity.
Seasonal mean values (Fig. 5) show that the correlation between these two proxies is significant only in the summer wet isotopically-depleted precipitation (e.g. Risi et al., 2008). Therefore, giant particles on the Nevado Illimani glacier can be reasonably used as proxy for deep summer convective precipitation. Given their size and geochemical/mineralogical fingerprint, we confidently associate them to local/regional convective activity near the Bolivian Altiplano.
In order to test the hypothesis of a relationship between giant particles and convective precipitation, we analyzed monthly precipitation records from five meteorological stations located in the central Andes (Fig. 1). Data was provided by SENAMHI, Bolivia (www.senamhi.gob.bo/sismet), whereas monthly outgoing longwave radiation (OLR) data on a 2.5° x 2.5° grid box (Liebmann and Smith, 1996) was obtained from NOAA/OAR/ESRL PSD, Boulder, Colorado, USA (https://www.esrl.noaa.gov/psd/). OLR data centered at 17.5 °S, 70 °W was used as an index of the convective precipitation over the Altiplano, as it presents strong negative correlations with regional rainfall observations (Garreaud and Aceituno, 2001). These datasets were re-sampled into DJF (December to February) and JJA (June to August) time series and compared with our GPPms series. Results reported in  (Aceituno, 1996) it is expected that only precipitation data from stations in the closest vicinity of the Nevado Illimani show good agreement with glaciological data.
Also, GPPms is negatively correlated with the DJF OLR centered over the Altiplano (Table 2), indicating that deep convection increases giant particle entrainment and suspension, humidity and precipitation over the region.
During dry seasons, conversely, GPPms is not significantly correlated with δD, as shown by seasonal mean values reported in Fig. 5. We conclude that the more intense is summer convection, the higher is the relative mass of giant dust particles suspended in the atmosphere and the more depleted is the δD.
Interestingly, Fig. 5 shows that over the 18-years period spanned by the Nevado Illimani firn core analyzed in this work, the summer seasons showing the most intense levels of convection were 1999-2000, 2000-2001, 2005-2006, 2007-2008, 2010-2011 and 2011-2012. These all correspond to La Niña years, as it is also clear from Fig. S2, where seasonal PC1 (proxy for convection) is compared to the Oceanic Niño Index. Some months after cooling of surface waters in the Pacific Ocean along the western coast of South America, we observe intensification of convective activity over the Altiplano. It is well known that the El Niño-Southern Oscillation phenomenon has a significant impact on climate over the Altiplano, especially during the summer season. In particular, meteorological data show that La Niña conditions intensify the meridional pressure gradient on the northern side of the Bolivian High, leading to stronger high troposphere easterly winds, increased eastward upslope flow and enhanced moisture transport (Garreaud, 1999;Vuille, 1999). Our data show that during summer seasons characterized by La Niña events, the atmospheric patterns typical of wet summer seasons are emphasized. Thus, we propose a new approach for future studies aimed at the reconstruction of paleo-ENSO at high elevation glaciological sites in the Andes based on giant particles and stable isotopes of snow, which can be used as a complement to a number of other climate proxies and modeling experiments providing insights into past ENSO variability.