SEDIMENT CORING AND PREVIOUS PROXY ANALYSES

Parallel sediment cores (core 1 and core 2) were collected from LM1 in 1989 from an anchored platform using a square-rod piston corer for the deeper sediments (> 1 m) and a PVC pipe fitted with a rubber piston for the near-surface sediments. The sediments of core 1 are the focus of this study and of the bulk and compoundspecific carbon isotope analyses carried out by Lane et al. (2011). All other proxy analyses, including pollen, microscopic charcoal, macroscopic charcoal, diatom, and basic geochemical analyses, were conducted using sediments from core 2 (Horn, 1993; Haberyan and Horn, 1999; League and Horn, 2000). The chronology for LM1 core 1 is based on six AMS radiocarbon dates on bulk sediment (Lane et al., 2011). Radiocarbon dates were converted to calibrated years B.P. (where P. = A.D. 1950) using the CALIB 6.0 computer program (Stuiver and Reimer, 1993) and the data set of Reimer et al. (2004). Sedimentation rates were calculated using the weighted means of the probability distributions of the calibrated ages (Telford et al., 2004), and used to linearly interpolate ages for lake sediment horizons located between the positions of radiocarbon-dated materials.

n-ALKANE EXTRACTION AND ISOTOPIC ANALYSIS

Alkanes were extracted from LMl core 1 sediment samples (1–16 g dry weight) using an ASE 300 accelerated solvent extraction system (ASE 300, Dionex, California, U.S.A.). Samples were flushed with a mixture of 50% dichloromethane and 50% methanol at 125 °C at a pressure of 1500 psi for 10 min. This procedure was repeated three times, yielding a total organic extract that was dried under high purity nitrogen. The aliphatic fraction was then isolated from the total organic extract by solid phase extraction through a silica-gel column with hexane. Straight chain monomers were further isolated from branched and cyclic compounds using urea adduction.

Compound-specific hydrogen isotope analyses were performed on a continuous flow Finnigan MAT Delta-plus XL mass spectrometer interfaced with an Agilent 6890 GC fitted with a DB5-MS silica column (60 m, 0.25 mm i.d., 0.25 µm film thickness) via a Finnigan GC/C III interface equipped with a pyrolysis reactor at 1430 °C. Injection temperatures were 300 °C and the oven temperature program was 60 °C isothermal for 2 min, 15 °C/min to 170 °C, 4 °C/min to 320 °C, and 320 °C isothermal for 30 min. Hydrogen gas and alkane mixtures (provided by A. Schimmelmann, Indiana University) with known hydrogen isotope compositions were used as standards. The Indiana University n-alkane mixture was injected following every fourth sample to monitor precision of the instrument. The standard deviation of these standard analyses over the course of this study was less than 5‰. ³H factors were calculated daily using the Isodat software with pulses of increasing reference gas amount. All samples were run in duplicate. Hydrogen isotopic compositions are reported in standard δ-per mil notation relative to Vienna Standard Mean Ocean Water (VSMOW) where:

n-ALKANE SOURCEAND δD INTERPRETATION

While the C₂₉, C₃₁, and C₃₃ n-alkanes are typically interpreted as the most reliable biomarkers for terrestrial vegetation, analyses of the n-alkane distribution of leaf waxes in modern specimens of Chusquea bamboos in the Colombian Andes indicated a disproportionate predominance of C₂₅ and C₂₇ n-alkanes (Boom et al., 2002). This finding, together with the close correspondence of C₂₅ and C₂₇ δD values through the LM1 sediment profile, lead us to consider that these two n-alkanes are originating from the same source and should be representative of the páramo environment surrounding LM1, where Chusquea subtessellata dominates in both cover and biomass. The comparatively different δD pattern observed for the longer-chain C₂₉, C₃₁, and C₃₃ n-alkanes in the LM1 sediments, and their relatively high abundance despite the dominance of C. subtessellata in the LM1 watershed (Lane et al., 2011) suggests that these alkanes reflect a different vegetation source. We conjecture that they may derive from montane forest ecosystems at lower elevations, where the waxes are abrasively removed by winds and transported upslope (Schefuβ et al., 2003).

Our original interpretation, based on carbon isotopes, was that upslope transport was likely an insignificant source of n-alkanes to the LM1 sediments, and probably did not strongly affect stable isotope data (Lane et al., 2011). However, the opposite trends in δD values of the shorter-chain (C₂₅ and C₂₇) and longer-chain (C₂₉, C₃₁, and C₃₃) n-alkanes in our new data set indicate that the plants producing these biomarkers are not accessing the same water sources or originating from the same locale. Garcin et al. (2012) found a similar discrepancy in δD values of sedimentary C₂₉ and C₃₁ alkanes collected along an elevational and latitudinal transect in Cameroon. Garcin et al. were unable to definitively explain this discrepancy, but proposed the possibility that plants producing an abundance of C₃₁ alkanes were accessing evaporated lake waters while plants producing an abundance of C₂₉ alkanes were more directly accessing meteoric water prior to significant evaporation. We similarly consider the most parsimonious explanation for the diverging temporal trends in LM1 δD data to be that plants are accessing water of differing δD composition, especially during the early to middle Holocene. However, it is unlikely that such divergent trends in δD values could be produced wholly within the páramo ecosystem surrounding LM1 because of the predominance of shallow, poorly developed soils and high daily solar radiation under cloud-free conditions, which combine to make the highaltitude páramo ecosystem particularly prone to soil water evaporation prior to incorporation by plants (Polissar and Freeman, 2010). Thus, we propose that the differing trends in δD values of shorterchain and longer-chain n-alkanes are the result of significant elevational differences in precipitation and subsequent evaporative enrichment of these meteoric waters prior to uptake by plants. We propose that it is these elevational differences that ultimately produced differing δD temporal trends of n-alkanes sourced from lower-elevation forest communities versus those produced in the high-elevation páramo ecosystem. With this in mind, we consider temporal variations in the δD of C₂₅ and C₂₇ n-alkanes isolated from the sediments of LM1 to primarily reflect either changes in the δD composition of meteoric waters entering the Chirripó páramo ecosystem or variations in soil water evaporation in the páramo. Conversely, the δD values of C₂₉, C₃₁, and C₃₃ n-alkanes more likely represent the δD composition of meteoric waters entering lower-elevation montane forests and any subsequent evaporative enrichment of deuterium prior to uptake by vegetation.

Prior studies have established an important life form effect on ecosystem n-alkane δD signatures (Liu et al., 2006; Smith and Freeman, 2006; Polissar and Freeman, 2010), but available paleoecological data indicate that the páramo ecosystem surrounding LM1 was established soon after deglaciation of the valley ∼10,000 cal yr B.P. (Horn, 1993), making life form an unlikely driver of the millennial-scale δD dynamics. Compound-specific δ¹³C data from LM1 indicate a potentially higher predominance of C₄ photosynthetic plants in the late Pleistocene and early Holocene, but the Chirripó páramo ecosystem likely included C. subtessellata as a significant component of the landscape at that time and throughout the postglacial period (Horn, 1993; Lane et al., 2011).

The most likely causes of millennial-scale variations in the δD composition of meteoric waters entering the Chirripó páramo ecosystem prior to soil water evaporation are amount or altitude effects. The amount effect has been documented as the primary driver of meteoric water isotopic composition along the Central American isthmus, with higher rainfall totals yielding more negative δ¹⁸O and δD values (Dansgaard, 1964; Lachniet and Patterson, 2002). δD values of meteoric waters also show an inverse relationship with elevation (Lachniet and Patterson, 2002), a phenomenon primarily driven by changes in atmospheric temperature. While we can assume no significant change in the elevation of LM1 over the time span of this study, we must consider the possibility that isotherms may have shifted up- or downslope in the Cordillera de Talamanca during the Holocene. Available paleoclimate data indicate temperature depressions of 7–8 °C during the Last Glacial Maximum in the highlands of Costa Rica (Hooghiemstra et al., 1992; Orvis and Horn, 2000), but similar reconstructions have not been possible during the Holocene due to a lack of suitable sites and proxies. δD values display a positive relationship with temperature change, with decreasing temperatures causing a decrease in δD values of meteoric waters and vice versa (Dansgaard, 1964). Fossil pollen data from LM1 show no clear evidence of shifting treelines during the Holocene (Horn, 1993), leading us to assume that any temperature effects on the δD composition of precipitation are relatively minor through the Holocene.

Evaporative enrichment of deuterium in soils is considered to be the primary driver of variations in the apparent fractionation (ε_a) between meteoric water δD and terrestrial n-alkane δD in páramo ecosystems (Polissar and Freeman, 2010). The term ε_a accounts for the fractionation effects on meteoric waters during biosynthesis and soil water évapotranspiration. Thus, ε_a increases significantly with increased aridity as soil water is evaporatively enriched in deuterium pushing it closer to the δD composition of meteoric waters (Polissar and Freeman, 2010). Polissar and Freeman (2010) reported ε_a values from approximately -90 to -110‰ in modern páramo ecosystems in Venezuela. Discounting the anomalous δD value of the uppermost sample in the LM1 record, and assuming an average C₂₇ δD value of -175‰ for C₂₇ n-alkanes in the LM1 profile (Fig. 2; Appendix) and a meteoric water δD composition of about -90‰ (Lachniet and Patterson, 2002), our estimated ε_a (-85‰) falls close to the high end of this range. With these variables in mind, we consider the most likely drivers of Holocene variation in C₂₅ and C₂₇ δD in the LM1 record to be amount effects and/or evaporative enrichment of soil water that increase ε_a values. In either case, increased aridity should increase C₂₅ and C₂₇ δD values, while an increase in precipitation should decrease C₂₅ and C₂₇ δD values.

THE LM1 TERRESTRIAL n-ALKANE δD RECORD

The anomalously high δD values in all of the measured n-alkanes in the uppermost sample are difficult to explain (Fig. 2). We find it unlikely that such a high δD value could have resulted from hydrological changes alone, as no other proxies indicate such severe shifts. Lichtfouse et al. (1994) documented the potential algal origin of some long-chain alkanes following low-temperature diagenesis in sedimentary systems that could explain this anomaly if algal organisms were imparting a very different fractionation on meteoric waters. A wide range of fractionation factors has been reported for algae under similar growth conditions (Zhang and Sachs, 2007). Alternatively, low-temperature diagenetic processes could explain the decreased relative abundance of C₃₁ and C₃₃ n-alkanes in the uppermost sample, although we find it puzzling that this effect would only be obvious in our uppermost sample, when the older sediments have had more time to undergo diagenesis.

Discounting the anomalous upper sample, the pattern of C₂₅ and C₂₇ δD variability through the Holocene corresponds well with the terrestrial n-alkane δ¹³C record presented by Lane et al. (2011; Fig. 3). This result supports the hypothesis that the observed variations in terrestrial n-alkane δ¹³C are a signal of shifting ecosystem drought stress during the Holocene. The overall pattern of relatively high δD values during the early Holocene, trending toward more negative values during the middle Holocene, and back to more positive values during the late Holocene, indicates millennial-scale shifts in Holocene aridity around LM1—specifically, relatively arid conditions during the early and late Holocene, interrupted by a relatively moist middle Holocene. Our reconstruction from isotopic analyses of the LM1 sediments shows similarities to the Holocene titanium record from Cariaco Basin sediments, which reflects precipitation in adjacent northern South America (Fig. 3) and has been interpreted as a record of millennial-scale ITCZ dynamics in response to the orbital precession cycle (Haug et al., 2001). The Cariaco Basin titanium record indicates a more rapid transition to mesic conditions during the early Holocene than the LM1 δD record, but the general trend toward more mesic conditions during the middle Holocene is apparent in both records. The slower response of the LM1 record to changes in ITCZ dynamics may be related to more rapid early Holocene warming of northern South America (Cariaco Basin sediment source area) compared to the Central American isthmus, which is surrounded by large water bodies. The high thermal inertia of the Pacific Ocean and the Caribbean Sea relative to the South American landmass may have restricted migration of the ITCZ over the isthmus until sea surface temperatures warmed sufficiently. Similar restrictions in northward ITCZ mobility occur today during El Niño events due to relatively cool sea surface temperatures near the Central American isthmus (Giannini et al., 2001). The tight linkage between modern ITCZ migration and precipitation around LM1 supports the hypothesis that shifting ITCZ migrational dynamics are the primary mechanism responsible for shifting paleohydrology in the circum-Caribbean region during the Holocene.

On Cerro Chirripó, precipitation is tightly linked to ITCZ dynamics via the trade wind inversion (TWI; Lane et al., 2011). When the ITCZ is displaced southwards, convective activity in the region decreases, and the TWI intensifies and lowers, often to an elevation below the highest peaks of the Cordillera de Talamanca (Dohrenwend, 1972; Stadtmüller, 1987; Hastenrath, 1991). This inhibits cloud formation and precipitation on the high mountain tops, thereby increasing E/P ratios in the lakes and surrounding watersheds. Conversely, northward migration of the ITCZ promotes convective activity in the region, and weakens and raises the elevation of the TWI, allowing for cloud formation and precipitation on Cerro Chirripó and the other high peaks of the Cordillera de Talamanca.

FIGURE 3.

Comparison of Lago de las Morrenas 1 core 1 stable isotope data (bulk sedimentary δ¹³C, C₂₇ n-alkane δ¹³C, and C₂₇ n-alkane δD) to Lago de las Morrenas 1 core 2 macroscopic charcoal influx data (League and Horn, 2000), and the Cariaco Basin titanium record (Ocean Drilling Program Site 1002; Haug et al., 2001). The age model for the Cariaco Basin data has been adjusted slightly such that B.P. refers to years before A.D. 1950. Dates extrapolated from the lowest date (11,230 cal yr B.P.) are italicized and marked with asterisks; these should be regarded as estimates for graphing purposes only (see Lane et al., 2011). V-PDB = Vienna-Peedee Belemnite, VSMOW = Vienna Standard Mean Ocean Water.

If ITCZ and TWI dynamics are the primary controls on millennial-scale precipitation dynamics on Cerro Chirripó, this could explain the predominantly inverse relationship between the LMI C₂₅ and C₂₇ δD record and the C₂₉, C₃₁, and C₃₃ δD record (Figs. 2 and 4). A more southerly mean annual position of the ITCZ and a related lowering of the mean annual height of the TWI during the early and late Holocene could lower mean annual cloud heights and lifting condensation level (LCL). This would increase aridity in pâramo ecosystems, but could leave subalpine and montane forests at lower elevations relatively unaffected, as these forests would still lie within, or below, the cloud belt. Limited available climate data from the Monteverde Cloud Forest of Costa Rica indicate relatively stable monthly average relative humidity (always >77%; Pounds et al., 1999; Johnson et al., 2005), even during periods of intensified TWIs that can cause relative humidity on the high peaks to drop well below 50% (Dohrenwend, 1972).

This altitudinal humidity and precipitation pattern would explain the relative stability of the C₂₉ and C₃₁ biomarker δD records that we interpret as having a primarily montane or subalpine source when compared to the C₂₅ and C₂₇ biomarker δD records that we interpret as primarily originating from the high-elevation páramo ecosystem (Fig. 2). This could also explain why C₂₉, C₃₁, and C₃₃ δD values increase between ∼7700 and 4200 cal yr B.P. while C₂₅ and C₂₇ δD values are becoming more negative. If the ITCZ were located in a more northerly mean annual position during this time, the TWI would have been positioned, on average, at a higher elevation, promoting cloud formation and precipitation in the páramos and perhaps increasing water stress at lower elevations via an upward vertical displacement of the LCL and associated clouds that typically bathe montane and subalpine cloud forests year-round. Pollen evidence from the La Chonta bog at 2600 m elevation in the Cordillera de Talamanca (Hooghiemstra et al., 1992) suggests shifts in forest composition that may be associated with early to middle Holocene aridity, perhaps resulting from an increase in mean annual TWI altitude. Existing climate models support this hypothesis as they indicate higher LCLs and reduced cloud contact with tropical cloud forests with increased temperatures (Still et al., 1999). This was likely the situation in southern Central America during the early to middle Holocene, when northern hemisphere solar insolation was significantly higher than southern hemisphere solar insolation (Voorhies and Metcalfe, 2007) and the ITCZ was located in a more northerly mean annual position (Haug et al., 2001). Simultaneous evidence of increased aridity in pollen records from relatively low- (∼1360 m) and high-elevation (∼3280 m) sites located near the lower and upper boundary of modern cloud forest distributions in the Peruvian Andes (Bush et al., 2004, 2005) may indicate a narrowed altitudinal range of these specialized ecosystems during the early to middle Holocene in South America. This pattern would be consistent with decreased ITCZ-sourced rainfall and related evapotranspiration rates at low elevations that can significantly impact cloud formation (Lawton et al., 2001) and lowered LCLs as a result of a lower mean annual TWI altitude.

FIGURE 4.

Comparison of δ¹³C and C₂₇ n-alkane δD isotope data with Lago de las Morrenas 1 core 2 macroscopic charcoal data (League and Horn, 2000), C₂₇ n-alkane δD values, and C₃₁ n-alkane δD values. Inferred climate conditions on Cerro Chirripó are based on relative changes in C₂₇ n-alkane δD values. Inferred climate conditions at lower elevations (montane and subalpine forests below the páramo) are based on relative changes in C₃₁ n-alkane δD values. The terms “mesic” and “arid” are relative and refer to temporal changes in nalkane δD values compared to the preceding time period. VSMOW = Vienna Standard Mean Ocean Water.

MILLENNIAL-SCALE PRECIPITATION VARIABILITY AND FIRE DYNAMICS IN THE LM1 WATERSHED

The C₂₇ δD record supports the interpretations of League and Horn (2000), who speculated that millennial-scale patterns of macroscopic charcoal deposition in LM1 core 2 were related to Holocene precipitation variability, with wetter conditions leading to fewer fires and lower charcoal influx. The lowest C₂₅ and C₂₇ δD values in LM1, between ∼7700 and 4200 cal yr B.P., coincide with very low macroscopic charcoal influx, supporting the hypothesis that relatively moist conditions during the middle Holocene suppressed fires on the Chirripó massif (Fig. 4). Mayle and Power (2012) documented a coincident mid-Holocene rise in fire activity in the Amazon basin that they associate with increased aridity. This pattern of a wetter Central American isthmus and a drier Amazon basin would be expected if the ITCZ were located in a more northerly mean annual position during this time. Conversely, the sudden rise in C₂₇ δD values after ∼4200 cal yr B.P. coincides with a steep increase in macroscopic charcoal influx, suggesting that increased regional aridity led to increased fire activity around LM1. This rise in charcoal ∼4200 cal yr B.P. also coincides with a sudden decrease in mean Ti concentrations in the Cariaco Basin, interpreted as decreased ITCZ-related precipitation in northern South America, further indicating a tight linkage between ITCZ dynamics and fire dynamics in the Chirripó páramo.

The documentation of the relationship between millennialscale precipitation and fire dynamics in the Chirripó páramo helps fill a gap in the understanding of neotropical páramo environments. While modern fire regimes are strongly shaped by humans, charcoal in ancient sediments shows that fire was present in neotropical páramos long before human settlement (Horn and Rappelle, 2009). Fires require dry fuels to spread, and, whether set by humans or nature, depend on short-term climate and weather. The relationships between human activity, climate, and fire are complex and are only beginning to be documented in neotropical highlands (Asbjornsen and Wickel, 2009; Horn and Rappelle, 2009). Resolving links between climate and fire over millennial time scales, as we have done for the Chirripó páramo with hydrogen isotope and macroscopic charcoal analysis, may provide key information for modeling and managing potential changes in fire susceptibility or frequency in páramo habitats under changing climate.