Precipitation variability, vegetation turnover, and anthropogenic disturbance over the last millennium in the Atacama highlands of northern Chile (19°S)

The Late-Holocene history of hydroclimatic variability in the Atacama Desert offers insights into the effects of precipitation and humans on ecosystems in one of the most extremely arid regions of the world. However, understanding the effects of regional precipitation variability in relation to local ecological stressors remains to be fully resolved. Here, we present a pollen-based qualitative precipitation reconstruction derived from fossil rodent middens recovered from two sites near Laguna Roja (LRO; n = 25) and Isluga (ISL; n = 15) in the Atacama highlands (19°S) of northern Chile. At LRO, the fossil pollen record shows multi-centennial hydroclimatic anomalies during the last millennium, with wetter than present phases at 1155–1130, 865–670, and 215–80 cal yrs BP, and similar to present conditions between 1005 and 880 cal yrs BP. In contrast, the ISL record shows a wet phase during 1115–840 cal yrs BP, suggesting that meso-ecological processes were as important in vegetation turnover as regional hydroclimate anomalies. Wetter conditions derived from LRO partially overlap with the Medieval Climate Anomaly (865–670 cal yrs BP) and with the latest part of the Little Ice Age (215–80 cal yrs BP). Furthermore, no strong anthropogenic signal was identified possibly related to the remote location of the records. Palynological diversity analyses evidence increasing diversification of plant communities during wet events at both sites. In correlation to existing regional hydroclimatic records from the Western Andes, our precipitation reconstruction verifies that centennial-scale changes in the strength of the South American Summer Monsoon (SASM) and partial influence of El Niño-like (ENSO) conditions drove vegetation turnover in the Atacama Desert during the last millennium.


Introduction
The drylands of the South America, expanding from southeastern Patagonia to the Peruvian northern coast, are among the most vulnerable regions to the impacts of aridification and land-use intensification worldwide (Schickhoff et al., 2022;Vuille et al., 2018). Meteorological data shows significant changes in precipitation patterns toward increasing aridity over the past several decades (Satgé et al., 2019;Torres-Batlló and Martí-Cardona, 2020;Zubieta et al., 2021). Furthermore, high-resolution tree-ring precipitation reconstructions provide evidence of widespread and unprecedented regional decline in rainfall since the 1930s and continuing until present times (Morales et al., 2012(Morales et al., , 2020. The incidence of extreme drought events and seasonal precipitation anomalies have direct impacts on the biota and socioeconomic sectors across the region, such as increased degradation of native ecosystems and decreased agricultural productivity (Lima et al., 2016;Satgé et al., 2019). Given that future climate projections suggest increasing drought frequency and intensity in the region including the Andean and Atacama highlands (Minvielle and Garreaud, 2011;Thibeault et al., 2012), paleoclimate reconstructions could help to assess the past impacts of hydroclimate variability on vegetation and human societies as well as provide relevant information for the design of effective strategies for climate change adaptation. Understanding long-term precipitation variability depends on the availability of extended records. However, in the Central Andes, instrumental data are only available since the mid-20th century.
Over the last decade, paleoclimatic reconstructions have elucidated the climatic variability in the Atacama Desert and the western Andean highlands (de Porras et al., 2021;Gayo et al., 2012a; Precipitation variability, vegetation turnover, and anthropogenic disturbance over the last millennium in the Atacama highlands of northern Chile (19°S) González-Pinilla et al., 2021;Jara et al., 2019Jara et al., , 2020Mujica et al., 2015). The Medieval Climate Anomaly (MCA) and the Little Ice Age (LIA) were two climatic events that produced important landscape transformations in the region. During the MCA (1200-650 cal yrs BP), wetter than modern conditions favored the expansion of riparian wetlands and agrarian societies flourished in the core of the desert, today an otherwise hyper-arid and barren landscape (Gayo et al., 2012a;Nester et al., 2007). In the Atacama Desert, the first part of the LIA (450-230 yrs BP) was characterized by drier than modern conditions, whereas in the Andean highlands, colder conditions favored the advance of glaciers between 400 and 100 cal yrs BP (Rabatel et al., 2008;Thompson et al., 1995, but see Gayo et al., 2012b;Nester et al., 2007). Glacier expansion started at 300 cal yrs BP and continued until 215 cal yrs BP (Rabatel et al., 2008). Following the LIA (100 cal yrs BP), dry conditions prevailed, and the region experienced an increase in temperatures and a persistent decrease in precipitation associated with the rise of greenhouse gas emissions (Morales et al., 2012(Morales et al., , 2020. In addition to climate variability, human activities have progressively affected portions of the Atacama highlands, valleys, and coast for millennia. The region was initially occupied by hunter-gatherers during the Late Pleistocene (~13,000 cal yrs BP) (Gayo et al., 2019;Santoro et al., 2017). Mobile groups likely hunted large artiodactyls such as guanacos (Lama guanicoe), vicuñas (Vicugna vicugna), and taruka deer (Hippocamelus antisensis) and set fires on natural vegetation, including grass steppes and woodlands (Santoro and Núñez, 1987). With the domestication of camelids (Lama glama and Vicugna pacos) and crop plants (Solanum tuberosum and Chenopodium quinoa), camelid herding in the highlands and irrigation farming in the lowland valleys eventually became the major economic strategies in the region starting ca. 3200 cal yrs BP (Agüero, 2005;Agüero and Uribe, 2011;Núñez and Santoro, 2011). With the European conquest of the Andes during the 16th century, landscape transformations increased through a surge in mining activities, herding of introduced animals (sheep, cattle, donkeys) and farming of introduced cultigens (onions, olives, and grapes) (Gavira Márquez, 2005;Marquet, 1998;Munizaga et al., 1975). At the onset of the Anthropocene, land-use has intensified leading to overgrazing, soil erosion, and ecosystem degradation (Pizarro et al., 2011). While the human footprint on the Atacama highlands is discernable, its extent and intensity have yet remained to be fully explored (Gayo et al., 2019).
In this study, we present a new reconstruction of precipitation variability based on pollen recovered from 40 rodent middens for the northern Atacama highlands during the last ~1100 years. We employed rodent middens, because they constitute a key environmental archive as they are abundant and rich repositories of plant remains (e.g. seeds, flowers, and pollen) (Betancourt and Saavedra, 2002;Holmgren et al., 2008;Latorre et al., 2002Latorre et al., , 2003. Previous research has demonstrated that pollen from rodent middens preserves excellent records of vegetation changes through time (Late Pleistocene to present) and space (local and regional scales) in the Atacama Desert (e.g. de Porras et al., 2017de Porras et al., , 2021Maldonado et al., 2005;Mujica et al., 2015). Furthermore, modern vegetation-pollen calibrations in rodent middens have demonstrated that these pollen assemblages provide balanced and robust representations of vegetation composition and distribution along environmental gradients that may also be used to understand past climate change (de Porras et al., 2015;Zamora-Allendes, 2015). Building on this baseline work, our pollen-based reconstruction was assisted by a modern pollen study from unconsolidated rodent middens recovered along an elevation transect at 19°S to determine the association between the current pollen spectra and the modern vegetation composition (Figure 1). Precipitation changes in time were conducted on the pollen assemblages and then validated against other regional records to assess the extent of past precipitation anomalies and teleconnections between the Atacama highlands, Altiplano, and the Pacific Ocean.
Note. Please refer to the online version of the article to view this figure in color. transport, and rainfall over the central Andes and is the product of the southward shift of the upper-level anticyclone which is directly related to enhanced summer easterly flow (Garreaud et al., 2003). Interannual summer precipitation variability in the Atacama is controlled by two modes of tropical precipitation: (a) the northerly mode, which transports moisture from the Amazon basin, located north and northeast and (b) the southerly mode, which brings moisture mainly associated with monsoonal precipitation from the Gran Chaco and the South Atlantic Ocean, located toward the southeast (Chaves and Nobre, 2004;Vuille and Keimig, 2004). The northerly mode is influenced by the El Niño-Southern Oscillation (ENSO) (Arnaud et al., 2001). ENSO is a quasiperiodic variation in sea surface temperature over the tropical Pacific Ocean. During El Niño years, warmer and drier conditions are experienced in the northern Atacama highlands, while during La Niña years, cooler and wetter conditions prevail (Garreaud, 2009;Garreaud and Aceituno, 2001;Vuille et al., 2000).
In northern Chile's Atacama Desert, three major west-east ecological zones are broadly distinguished: the Pacific coast and Coastal Cordillera; the Central Depression located between the coast and the Andes; and the Andean highlands and Altiplano. Within the Andean highlands and Altiplano, several distinct vegetation belts occur along an elevational gradient, which provide distinctive habitats for a diversity of plants and animals, including endemic species adapted to the extreme environmental conditions. Regional patterns of vegetation distribution are mostly controlled by rainfall, excluding elevations above 4800 m asl, where low temperatures and the snowline restrict plant growth (Betancourt and Saavedra, 2002;Díaz et al., 2019;Latorre et al., 2002). Four distinctive vegetation belts are clearly distinguishable on the western slope of the Central Andes: Prepuna, Puna, High Andean Steppe, and Subnival (Arroyo et al., 1988;Villagrán et al., 1981Villagrán et al., , 1983. Plant richness increases with elevation, reaching its maximum at the Puna, followed by a decrease in the High Andean Steppe and Subnival (Villagrán et al., 1983).
At 19°S, the Prepuna occurs above the absolute desert, between 2600 and 3300 m asl, and is inhabited by sporadic xerophytic shrubs of Ambrosia artemisioides, intermixed with Atriplex microphylla, Browningia candelaris and Corryocactus brevistylus. Annuals herbs include Descurainia stricta, Cistanthe amarantoides, C. celosioides, Nolana tarapacana, Reyesia juniperoides, Tetragonia macrocarpa, and Tagetes multiflora (García et al., 2018;Trivelli-Jolly and Valdivia-Ríos, 2009;Villagrán et al., 1999). The Puna, traditionally described as a single unit (3300-4000 m asl), is composed mostly by shrubs, and two vegetation sub-belts levels can be distinguished (Luebert and Pliscoff, 2006;Villagrán et al., 1999). The Low Puna (3300-3600 m asl; commonly known as traditional Puna) is dominated by shrubs of Baccharis boliviensis, B. tola, and Adesmia spinossisima, and accompanied by Fabiana ramulosa, Balbisia microphylla, and Diplostephium meyennii (García et al., 2018;Trivelli-Jolly and Valdivia-Ríos, 2009 Today, Aymara and Atacameño/Kunza indigenous communities practice a mixed economic strategy that complements animal husbandry with migration to Chilean cities to sell their labor. Native herds of llamas and alpacas share pastures with introduced sheep, all of which rely on natural fodder from steppe grasslands and wetlands (Dransart and Ortega Perrier, 2021). Herders periodically burn high-altitude vegetation, particularly steppes, to promote production of more palatable, young grass shoots for livestock and to clear land for agriculture. Small-scale agriculture is restricted to localized microenvironments in the Puna and characterized by extensive rotational cultivation of frost-resistant varieties of potatoes, quinoa, and barley as well as intensive agriculture of tomatoes and other leafy vegetables facilitated by artificial irrigation in the alluvial plains of the western lowland valleys draining into the Pacific (García et al., 2018).
Laguna Roja (LRO; 3700 m asl; 19.05°S, 69.25°W) is a small shallow lake, located in the High Puna vegetation belt and characterized by an intense red color and high concentrations of arsenic and boron (Cornejo et al., 2019;Tapia et al., 2021). Isluga (ISL; 4100 m asl; 19.19°S, 68.90° W) is a stratovolcano located near the border with Bolivia, situated within the High Andean Steppe. Both Laguna Roja and Volcan Isluga are considered religious sites by indigenous communities. Volcan Isluga, in particular, is one of the most sacred mountains in the region, and rain-calling ceremonies are celebrated at the mountain summit (Dransart and Ortega Perrier, 2021;Reinhard, 2002).

Materials and methods
We retrieved 25 fossil rodent middens near LRO and within the uppermost part of the High Puna (3500-3700 m asl) and 15 fossil rodent middens near ISL inside a P. tarapacana woodland in the High Andean Steppe (~3900 m asl). Additionally, in 2012, we collected 16 modern unconsolidated rodent middens along a westeast altitudinal transect from Tana to Isluga (TA-IS) between 2900 to 4600 m asl in the western Andes of northern Chile (19°S). Middens were recovered every 100 m in elevation to establish the representation of vegetation belts in modern rodent midden pollen assemblages.
Fossil middens were recovered using a hammer and chisel and cleaned in the field; whereas, unconsolidated rodent middens were collected using a spatula according to established procedures (Spaulding et al., 1990). A sample of each fossil rodent midden was soaked in distilled water for 24-48 h to dissolve the enclosing urine and subsequently screened through a 120 μm sieve to separate macro-(seeds, leaves, stems, etc.) from microremains (midden matrix). Modern rodent middens were processed following the protocol outlined in de Porras et al. (2015). A subsample of 1-cm 3 of the micro-remains was processed following standard methods for pollen extraction, including acetolysis, hydrochloric and hydrofluoric acid to remove cellulose and to dissolve carbonates and silicates (Faegri and Iversen, 1989). A tablet of exotic Lycopodium spores was added to each sample to calculate pollen concentration. In each sample, we identified a minimum of 300 pollen grains whenever it was possible, excluding Cyperaceae and other aquatic taxa. Identification of pollen grains and palynomorphs relied on the pollen reference collection at the Laboratorio de Paleoecología y Paleoclima (Centro de Estudios Avanzados en Zonas Áridas) and palynological atlases and keys (Heusser, 1971;Markgraf and D'Antoni, 1978;Miesen et al., 2015;Sandoval et al., 2011).
To determine the age of each midden and generate a chronological sequence of the midden series, a sample of three to four fecal pellets from each fossil midden was dated by accelerator mass spectrometry (AMS) radiocarbon dating. We controlled for contamination by sampling from the interior of the middens and avoiding loose pellets from the outer area of the midden. One set of samples were AMS 14 C dated at Direct AMS Lab (D-AMS) and another at the Pennsylvania State University AMS Radiocarbon Lab (PSUAMS). Rodent midden 14 C ages were calibrated to cal yrs BP with the SHCal20 Southern Hemisphere curve (Hogg et al., 2020) using Oxcal 4.4 (Bronk Ramsey, 2021) (Table 1).

Data analysis
Fossil middens were grouped according to their provenance and chronology as determined by their calibrated radiocarbon ages. A Constrained Incremental Sum of Squares (CONISS) cluster analysis was further performed to divide the fossil midden sequences into zones of similar pollen composition, considering all local pollen taxa contributing at least 2% (Grimm, 1987). To calculate the relative abundance of each taxon, we used the total terrestrial pollen sum per midden sample and plotted as percentages using the software Tilia (Grimm, 1991).
The reconstruction of past precipitation changes was based on the documented relationships between vegetation and summer precipitation along the west-east altitudinal gradient in the western Andean Andes of Chile ( middens are interpreted as altitudinal shifts in vegetation and used to estimate past relative changes in precipitation as well as the result of anthropogenic effects as suggested by the relative frequency of introduced taxa and indicators of disturbance (e.g. Amaranthaceae, Trifolium). This is supported by the modern pollen-vegetation relationship presented here and other studies (Collao-Alvarado et al., 2015;de Porras et al., 2015de Porras et al., , 2017Jara et al., 2019;Ortuño et al., 2011). We compared each fossil pollen zone to its modern pollen analog to reconstruct past hydroclimatic conditions in comparison to present climate within that modern vegetation zone. For instance, at LRO, we interpreted phases of relative high abundance of High Andean Steppe (Puna) pollen types as a downward (upward) shift of this vegetation belt in altitude regarding the present thus indicating wetter (drier) than modern conditions. We additionally reconstructed past precipitation variability by qualitatively comparing the fossil pollen assemblages to the modern pollen rain from each rodent midden and assigning a code to each fossil pollen zone as M (similar than modern); D-(slightly drier than present); W+ (slightly wetter than present) and W++ (moderately wetter than present) (de Porras et al., 2017, 2021; Zamora-Allendes, 2015). The reconstruction of past precipitation was compared to published precipitation reconstructions from a diversity of paleoclimate studies, including rodent middens, ice cores, tree-ring chronologies, and wetland deposits. We estimated palynological variability using several diversity indices. Due to differences in taxonomic precision, differences in pollen production by each taxa and pollen dispersal strategies, it is important to note that there is a nonlinear relationship between palynological diversity and biological diversity (Meltsov et al., 2011;Matthias et al., 2015;Sugita et al., 1999). Therefore, to compare vegetation diversity among samples with varying number of pollen grains, we computed Fisher's alpha diversity (α) as S = a × ln(1+n/a), where S is number of pollen taxa, n is number of pollen grains, and a is Fisher's alpha (Fisher et al., 1943). Pielou's evenness index (J) was used to estimate variability in taxa abundances within a pollen assemblage, calculated as J = H′/ log (S), where H' is the Shannon index and S the number of taxa for each sample (Pielou, 1969). The Shannon index was calculated as H′ = -∑[(p i ) × log(p i )], where p i is the percentage of pollen grains for each taxon in a sample. Given that estimates of palynological richness and evenness can be potentially biased due to different count size between samples, a rarefaction method (J(S n )) was employed to estimate the expected number of taxa in a sample (Birks and Line, 1992;Simberloff, 1972). In the case of LRO middens, a pollen sum of 310 pollen grains was used and in ISL middens a pollen sum of 111 pollen grains. This estimate allows standardization of count-size and thus comparisons of richness between samples. All analysis were conducted with the vegan package (Oksanen et al., 2007) in R environment version 4.0.3 (R Core Team, 2021).

Modern (unconsolidated) rodent middens
The pollen assemblages from modern rodent middens recovered from the TA-IS transect reflect modern vegetation belts along an elevational gradient. The pollen assemblages cluster into five groups, which were collected from the same modern vegetation belts (Figure 2). A detailed description of the pollen composition recovered from modern rodent middens associated with each vegetation belt is presented as Supplemental Information S1, available online. From 2900 to 3500 m asl, the pollen assemblage resembles the Prepuna vegetation, characterized by high abundances of Ambrosia type (5-30%) and Brassicaceae pollen (25-45%), while Amaranthaceae appears sporadically (<5%), except in one sample (15%). Between 3600 and 3900 m asl, the pollen assemblage is characteristic of the Low Puna vegetation belt. This midden series show high abundances of Baccharis type pollen (~50%) and Senecio/Parastrephia type pollen (~20%), accompanied by Brassicaceae pollen (~15%) and Fabaceae (Adesmia type) pollen, whereas Amaranthaceae pollen is found occasionally (<1%). From 4000 to 4200 m asl, the pollen assemblage is distinctive of the High Puna showing a sharp increase in Poaceae pollen (~40%) and a decrease of Asteraceae pollen (Senecio/ Parastrephia type, 40-20% and Baccharis type, 20-10%), abundance of Amaranthaceae pollen remains relatively low (<3%). In the High Andean Steppe (4300-4500 m asl), the pollen assemblage shows the dominance of Poaceae (<40%), followed by Senecio/Parastrephia type (25-30%), whilst Amaranthaceae pollen is also barely found (~1%). Above 4500 m asl, the pollen assemblage reflects the Subnival vegetation belt and is characterized by the high abundances of Caryophyllaceae pollen (40-10%) and Apiaceae pollen (20-5%) and a decrease of Poaceae pollen (15-25%), Amaranthaceae pollen is found occasionally (1%). Average annual precipitation increased from <30 mm yr -1 at 2900 m asl to >90 mm yr −1 above 4600 m asl. In contrast, mean annual temperature, associated with the altitudinal transect, exhibited a constant decrease from 9°C at 2900 m asl to less than 3°C at elevations above 4600 m asl.

Laguna Roja rodent midden record
Pollen concentrations were relatively high (average = 51,952 grains/ gram midden). Moreover, pollen preservation was good as the abundance of broken grains was characteristically low (<5%). We identified 44 pollen taxa in the rodent midden series. Most of the LRO rodent middens dated between 85 and 990 cal yrs BP (n = 21), while a few dated between 1005 and 1155 cal yrs BP (n = 3), and one dated 1821 cal yrs BP (Table 1; Supplemental Information S2a, available online). The CONISS cluster analysis identified three major zones and five subzones that we describe below (Figure 3a). Asteraceae  (Baccharis and Senecio/Parastrephia) pollen dominated the pollen assemblages, followed by Poaceae pollen, suggesting that shrubs, intermixed with grasses dominated the landscape for at the last ~1155 years.

Isluga rodent midden record
Pollen preservation of the ISL rodent middens was good as the percentage of broken grains was relatively low (<5%), and in most samples the pollen concentration was relatively high (average = 14,600 grains/gram midden). A total of 29 taxa were identified, from which Asteraceae and Poaceae pollen were the most abundant, followed by Polylepis pollen. Other persistent elements include Azorella, Amaranthaceae, Gentianaceae, and Mutisia. The ISL rodent middens showed sparse temporal coverage concentrated mostly during the last 160 years (n = 10), with four of the middens dating between 1115 and 840 cal yrs BP, and a single sample dating to 4417 cal yrs BP (Table 1; Supplemental Information S2b, available online). The CONISS cluster analysis identified three major pollen zones that we characterize below (Figure 3b).
ISL 3 (160 cal yrs BP-modern) was characterized by an increase in the abundance of Poaceae pollen (30-57%). Polylepis (11-38%) and Baccharis type (9-24%) pollen were relatively common but fluctuated considerably between samples. Azorella, Adesmia, and Senecio/Parastrephia type also varied substantially reaching up to 15%, while exotic taxa, such as Myrtaceae and Pinus, were occasionally present. The pollen assemblage suggests modern conditions (M). Diversity indexes show similar values to

Modern pollen-vegetation relationships
The pollen spectra recovered from modern rodent middens spanning a 1703 m elevation gradient (equivalent to 91 mm in precipitation gradient and 7°C in mean annual temperature) in the Atacama highlands of northern Chile shows clear differences in Pollen diversity. This study also elucidates changes in the diversity of plant communities during the last millennia. Diversity indices showed that vegetation turnover took place during dry and wet events. At LRO, the highest palynological diversity values were recorded during three wet events (1155-1130, 865-670, and215-180 cal yr BP) and during the transition from wet to dry conditions between 600 and 425 cal yr BP. In contrast, periods of lower palynological diversity were observed during a wet event (660-605 cal yr BP) and a period with similar to modern conditions (1005-880 cal yr BP). General trends show that during periods of increasing precipitation, plant communities diversified due to an increase in abundance of Puna taxa and the retraction of grasses and other plant elements from the Prepuna. The retraction of grasses into LRO could have favored the establishment of local species (e.g. Deyeuxia, Festuca, and Anatherostipa), whichAt ISL, increasing palynological diversity coincided with one wet event (1115-840 cal yr BP) and during the last 160 years. In this site, relatively diverse plant communities are associated with the expansion of Polylepis woodlands, which form stable microclimates and constitute the habitat for several plant and animal species (Fjeldså and Kessler, 1996;Tapella et al., 2021). The results highlight the role of precipitation as a strong driver of vegetation composition in the region, but also show complex ecological dynamics that drive community assembly across sites. The pollen of Polylepis recovered from the rodent midden series likely indicates the presence P. tarapacana. This is the only species of the genus that currently inhabits the semiarid highlands (above 4000 m asl) in the Central Andes (Boza Espinoza and Kessler, 2022). The low abundance of Polylepis pollen at 4417 cal yrs BP suggests either the presence of scattered individuals due to dry and cold conditions or the early formation of a Polylepis woodland. Between 1115 and 840 cal yrs BP, increased abundance of Polylepis pollen indicates that these woodlands were already stablished in Isluga and constituted a conspicuous element of the landscape. Today, Polylepis woodlands persist as small fragments across the region, partially due to wood collection that reduced forest cover and increased habitat degradation. Our study also highlights the importance of using rodent middens to conduct local and regional reconstructions of Polylepis over time, as they contained abundant and well-preserved Polylepis pollen.

Late-Holocene vegetation change in the Atacama highlands (19°S)
Laguna Roja. The pollen record from fossil rodent middens at Laguna Roja provides a reconstruction of the vegetation changes during the last ~1820 cal yrs BP within the High Puna. Around 1820 years ago (LRO 1), shrublands of Baccharis and Senecio/Parastrephia dominated the region, while grasses were relatively scarce. Palynological diversity and palynological evenness were comparatively high suggesting a diverse but homogeneous vegetation community.
Between 1155 and 1130 cal yrs BP (LRO 2), the pollen assemblage suggests the dominance of High Puna shrubs (Baccharis and Senecio/Parastrephia) and the intrusion of High Andean grasses into the vicinity of Laguna Roja, which is ~200 m lower than their modern distribution. Palynological and palynological diversity were relatively high during this interval, suggesting a diverse but rather homogenous plant community. The midden series also shows a slight increase in abundances of Amaranthaceae during this period, which could indicate quinoa cultivation and/or soil disturbance by domestic camelids (García and Ajata, 2016;García et al., 2018;Ledru et al., 2013;Medina et al., 2017). It is also possible that an increase in abundance of Amaranthaceae is associated microenvironmental conditions, such as salt accumulation in soils that could have favored their local establishment.
Between 1005 and 880 cal yrs BP (LRO 3a), High Puna shrubs (Baccharis and Senecio/Parastrephia) increased in abundance, while Prepuna taxa (e.g. Amaranthaceae, Ambrosia, Ephedra) and High Andean Steppe grasses withdrew, resulting in a vegetation quite similar to present at the site. During this period, high values of Fisher's alpha and Pielou's evenness point to a relatively diverse vegetation.
Following this trend, by 865-670 cal yrs BP (LRO 3b), Baccharis and Senecio/Parastrephia shrubs remained relatively dominant. Adesmia and Poaceae increased in abundance, but Ambrosia declined. Moreover, this period exhibited low values of palynological diversity (Fisher's alpha), but high values of Pielou's evenness. The pollen assemblages indicate a turnover of the vegetation, associated with the incursion of grasses and the expansion of High Puna vegetation.
The interval between 665 and 605 cal yrs BP (LRO 3c) includes a pollen assemblage relatively similar to the one before, with the exception that grasses expanded. In tandem, palynological evenness suggests that the landscape was more homogenous with a slight increase in palynological diversity. This could be associated with a temporary retraction of Puna vegetation adapted to more xeric conditions (e.g. Chuquiraga, Ophryosporus) as well as the persistence of grasses that could form dense tussocks (e.g. Deyeuxia, Festuca, and Anatherostipa).
The period between 600 and 425 cal yrs BP (LRO-4a) was characterized by the contraction of grasses and the expansion of Adesmia shrubs. During this time, palynological diversity values slightly increase and evenness values remain low, suggesting a more diverse plant community.
The pollen assemblages between 215 and 180 cal yrs BP (LRO-4b) show the relative dominance of Puna shrubs (Baccharis, Senecio/Parastrephia, and Adesmia) and an expansion of grasses. Palynological diversity also increased, but evenness values remain low, associated with an increase of Puna elements (e.g.

Chuquiraga, Mutisia, Fabiana).
Isluga. The rodent midden record from Isluga shows that ~4415 cal yrs BP (ISL 1), grass steppe with Baccharis shrubs dominated the landscape, while Polylepis woodlands was scarce. The composition of the pollen assemblage also shows that Subnival vegetation, particularly Azorella, was common. This pollen assemblage points to the upslope migration of Puna vegetation and downslope migration of Subnival vegetation. Palynological diversity and evenness show low values, suggesting a relatively species-poor but heterogenous plant community, where grasses were the dominant element.
Between 1115 and 840 cal yrs BP (ISL 2), grasses, Polylepis woodlands, and Azorella cushions became more abundant, while Puna shrubs (Bacharis and Senecio/Parastrephia) declined, implying a relative expansion of the High Andean Steppe and the retraction of Puna vegetation. This period also shows higher values of palynological diversity and evenness, pointing to a more diverse vegetation.
The last 160 years (ISL 3) were marked by a relative dominance of grasses, followed by Polylepis woodlands. Baccharis shrubs and Azorella showed a minor increase at 110 cal yrs BP. This interval also featured the presence of Pinus and Myrtaceae (likely Eucalyptus), which reflects the introduction of exotic tree species during the 19th century (Luzar, 2007), and is evidence of increasing landscape transformation by humans (Domic et al., 2018). Palynological and evenness also increase indicated a more plant community richer in taxa but relatively homogenous.
Modern samples (last ~70 years) show high variation of Asteraceae shrubs (Baccharis and Senecio/Parastrephia), Polylepis, and Azorella between samples. This could be associated with higher sample resolution of this interval, spatial heterogeneity, and internal restructuring of plant communities due to the interacting effects of land-use intensification and ongoing climate change. Differences in the abundance of Polylepis and Azorella between middens could be associated with extraction of individuals for firewood and the presence of relict plants in areas of difficult access. Variability of Puna shrubs could be connected with their gradual upslope migration due to increasing aridification experienced in the region during the last century (Morales et al., 2012). However, given that dated samples show large standard deviation values, results should be taken with consideration.

Late-Holocene vegetation and hydroclimate
Our findings complement previously published climate and vegetation reconstructions from fossil rodent middens along the northern Atacama Desert and provide a picture of vegetation dynamics during the last 12 centuries (de Porras et al., 2017;Gayo et al., 2012a;González-Pinilla et al., 2021;Jara et al., 2022;Olson et al., 2020). This new information in combination with other paleoclimatological reconstructions allows us to assess vegetation dynamics and the role of climatic variability and human disturbance during the last millennia. Specifically, the LRO rodent midden series shows evidence of periods of multi-centennial changes in summer rainfall with increased precipitation between 1155-1130, 865-605, and 215-180 cal yrs BP, and a period of similar to present precipitation conditions between 1005 and 880 cal yrs BP (Figure 4a). Additionally, the ISL rodent midden series shows a wet event between 1115 and 840 cal yrs BP (Figure 4b). The discrepancies in precipitation between LRO and (e) pollen percentage of aff. Parastrephia from La Higuera midden series (Mujica et al., 2015). (f ) Precipitation reconstruction based on tree-ring chronologies from the Altiplano region (Morales et al., 2020). (g) Southern Oscillation Index reconstruction (Yan et al., 2011). M (similar than modern); W+ (slightly wetter than present) and W++ (moderately wetter than present).
ISL between ~1005 and 840 cal yrs BP could be due to the relative position of both sites with respect to the vegetation gradient and the sensitivity of vegetation in those sites to a given precipitation change.
The centennial wet anomalies observed at LRO during the last ~1155 years correspond with wet events observed in several paleohydrological records from the western Andean slope and the Atacama Desert (Figure 1). Humid conditions at LRO between 1155 and 1130 cal yrs BP agree with a wet interval between 2000-1000 yrs reconstructed from Prosopis tree-rings at Pampa del Tamarugal (21°S) (Olson et al., 2020). The similar to present conditions inferred from the LRO rodent midden series between 1005 and 880 cal yrs BP, while a reduction in precipitation at Laguna Ceusis (21°S) took place between 1400 and 900 cal yrs BP (Figure 4c; Jara et al., 2020).
In the northern Atacama Desert, wet conditions are associated with an increase in precipitation during the Medieval Climate Anomaly (MCA; 1100-700 cal yrs BP) due to the strengthening of the SASM in association with the southward displacement of the ITCZ (Gayo et al., 2012a;González-Pinilla et al., 2021;Nester et al., 2007). The wet event recorded at the LRO midden series between ~1000 and 500 years ago, coincided with elevated water tables within canyons and springs across the northern and central Atacama Desert (18-23°S;Rech, 2001;Rech et al., 2002Rech et al., , 2003Tully et al., 2019). Although several paleoclimate reconstructions have demonstrated that the period from 1100 and 600 yr BP was relatively humid, our study also shows greater precipitation variability, which could be due to the greater sensitivity of rodent middens in comparison with other types of paleoclimate archives (e.g. wetland sediment records) that possess long response times.
The LRO record covers the last portion of the Little Ice Age (215-180 cal yrs BP), during which the high abundance of Asteraceae shrubs suggests wet conditions. Given that glacier studies in the tropical Andes (16°S) suggest cooler conditions in comparison with modern conditions (2.2-1.1°C cooler) between 600 and 100 cal yr BP, we consider that a decrease in temperatures could have also driven the expansion of shrubs and the retraction of Prepuna vegetation during this period (Jomelli et al., 2011;Rabatel et al., 2008). This is supported by the multivariate analysis of pollen spectra, which shows that abundance of pollen of taxa from the Prepuna are positively associated with temperature (Supplemental Information S3, available online).
The LRO pollen record shows some similarities with the precipitation reconstruction derived from fossil rodent middens at Toconce (~22°S), that signals dry conditions between 930-550 and wet conditions between 200 and 90 cal yr BP (Jara et al., 2022). Wet conditions were also documented between 855-840 and 520-450 cal yrs BP in fossil rodent middens from Cuesta Chita (~22°S; Figure 4d) (de Porras et al., 2021). Positive precipitation anomalies were also recorded at La Higuera (18°S), where an increase in Parastrephia pollen and bigger fecal pellets (a proxy for rodent body size and precipitation) suggest increased rainfall during a pluvial event between 190 and 100 cal yrs BP (Figure 4e) (Mujica et al., 2015). This later event has also been detected in a composite tree-ring chronology from the western Altiplano (Figure 4f; Morales et al., 2020) and Lago Chungará (18°S) water level reconstructions at ~150 cal yrs BP (Valero-Garcés et al., 2003). These relative differences suggest spatial variation in drought intensity from north to south in the region.
Hydroclimate anomalies recorded at LRO mostly differ with precipitation reconstructions from the Altiplano. A reconstruction based on δD wax and δ 18 O calcite suggest dry conditions at Lake Orurillo (14°S, 100 km north of Lake Titicaca) between 1035 and 925 cal yrs BP (Arnold et al., 2021). Lake Titicaca (16°S) experienced a moderate lowering in lake levels between 830 and 680 cal yrs BP (Weide et al., 2017). Lake levels at Lago Chungará (18°S) decreased between 700 and 45 cal yrs BP (Valero-Garcés et al., 2003). A reduction in precipitation was also recorded at Marcacocha (13°S) and at Cerro Llamoca (14°S) in Peru between 1070 and 550 cal yrs BP (Chepstow-Lusty et al., 2009;Schittek et al., 2015). The calcite and δ 18 O records from Laguna Pumacocha (10°S) show a dry event between 1050 and 850 cal yrs BP, followed by a wet event between 500 and 300 cal yrs BP (Bird et al., 2011).
The overall pollen assemblages from the ISL rodent midden series agrees with positive precipitation anomalies observed at Pampa del Tamarugal (21°S) based Radiocarbon dating of macrobotanical remains between 1100 and 680 cal yrs BP (Gayo et al., 2012a(Gayo et al., , 2012bNester et al., 2007). In contrast, the ISL record disagrees with precipitation reconstructions from Andean tropical ice-cores that suggest decreased monsoonal precipitation between 1000 and 680 cal yrs BP. Specifically, the ice-core oxygen isotopic record from Nevado Sajama (18°S; located ~100 km northeast from Isluga) indicates a dry period of between 1050 and 450 cal yrs BP, suggesting decreased SASM activity (Thompson et al., 1998). The Nevado Huascarán (9°S) ice-core oxygen isotopic record also shows a trend toward enriched values of δ 18 O between 1050 and 650 cal yrs BP (Thompson et al., 1995).
Over the last 2000 years, centennial precipitation changes have been attributed to ENSO-like activity in the tropical Pacific SST-gradient (Yan et al., 2011;Figure 4g). In the northern Atacama Desert, ENSO activity constitutes as a primary factor driving modern inter-annual rainfall variations (Houston, 2006) by influencing the northerly mode (Vuille and Keimig, 2004). During La Niña years, the Atacama highlands receive increased summer rainfall, promoting recharge of ground-water reservoirs and increasing water runoff of perennial streams. Paleoclimate records indicate that during most of the MCA (~1050-600 cal yrs BP), the tropical Pacific experienced century scale changes ENSO variation with an intensification toward the end (Gayo et al., 2012b;Jiang et al., 2021;Lawman et al., 2020;Yan et al., 2011). This evidence suggests that increased precipitation at Laguna Roja might have been related to La Niña-like conditions toward the end of the LIA, between 215 and 180 cal yrs BP.
The hydrological reconstructions based on fossil rodent middens from the Atacama highlands indicate a close correlation with precipitation reconstructions from the Atacama Desert and western Andean slope. The synchronicities of wet events suggest that the reconstructed centennial precipitation anomalies were regional climatic events driven by the strengthening of the SASM during the last ~1000 years. Furthermore, discrepancies in the paleoclimate reconstructions from the Atacama highlands and Altiplano suggest that the precipitation anomalies recorded at LRO and ISL were largely decoupled from north-south shifts of the ITCZ.
Pollen assemblages from both rodent midden series do not show significant changes in vegetation composition associated to human activities (e.g. an increase of pollen from disturbance-tolerant taxa). The ISL record only shows the presence of pollen of exotic trees (Myrtaceae and Pinus). The scarcity of archeological studies in the region further limits our capacity to address how humans managed and transform the region over the last millennium. Further paleoclimatic and paleoecological reconstructions, particularly from the western highlands of Bolivia where midden records are presently absent, can contribute additional information toward building a more comprehensive understanding of past hydroclimate dynamics and anthropogenic impact on vegetation turnover in the northern Atacama highlands.

Conclusions
The pollen analysis from fossil rodent middens shows that significant changes in the altitudinal distribution of vegetation communities and plant diversity occurred in the northern Atacama highlands over the last ~1160 years. At Laguna Roja, positive multi-centennial hydroclimatic anomalies were identified between 1155-1130, 865-670, and 215-180 cal yrs BP and similar to present conditions between 1005 and 880 cal yrs BP. At Isluga, a centennial wet event was detected between 1115 and 840 cal yrs BP. The coeval timing of these events with other paleoclimate records from the Atacama Desert and western Andean slope support the existence of important pluvial events on a regional scale during the last millennium. Furthermore, the synchronous response of the vegetation from the Atacama lowlands and highlands implies that rainfall anomalies impacted the hydrology, vegetation composition, and diversity at a regional scale. The positive hydroclimatic anomalies observed at LRO reflect partially a strengthening of the SASM. In particular, the positive precipitation anomaly recorded at LRO during part of the MCA (860-605 cal yrs BP) was likely caused by a prevalence of La Niña-like conditions at the centennial scale. Both rodent midden series suggest that changes in the strength of the South American Summer Monsoon drove precipitation variability in the Atacama highlands during the last millennium.

Author contributions
Individual contributions of each co-author to the article are: MEdP and AM designed the research. MEdP, AZA, AM collected the samples. AID, MEdP, AZA collected and analyzed the data and wrote the draft of the manuscript. JMC, AM and SJI assisted with the data analysis. AID, MEdP, SIJ, JMC and AM obtained funding. All authors commented, provided input to the manuscript, and approved the submitted version.

Data availability
The pollen data presented in this article will be freely available in the Neotoma Paleoecological Database (https://www.neotomadb. org/) once it is accepted for publication.

Funding
The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: We thank the indigenous communities of Laguna Roja and Isluga for permission to work in their territories. This study was funded by FONDECYT postdoctoral grant #3160443, IFS I-1-D-4839-2, FONDECYT #3130511, FONDECYT #1181829, REDES-180099-CONICYT, Millennium Nucleus UPWELL (ANID NCN19_153), and NSF DEB #2208411.

Supplemental material
Supplemental material for this article is available online.  (3)