Introduction 

Tropical montane forests are among the world’s most biodiverse ecosystems and constitute one of the largest pools of terrestrial carbon on Earth, thereby contributing to the modulation of global biogeochemical cycles and climate ([4]). These tropical montane forests are currently at risk as a result of intense anthropogenic resource exploitation. Threats include forest clearing for pasture and crops, illegal road constructions as well as settlement and habitat fragmentation ([18]). This situation is expected to be aggravated by ongoing climate change and human pressure, but data confirming this assumption remain scarce. In addition, the potential influence of climate (i.e., temperature and precipitation) on forest productivity is still poorly known in these latitudes, even though climate-growth linkages are seen as crucial to implement adequate adaptation targets to mitigate the negative effects of climate change. The impacts of the global climate system on terrestrial forest ecosystems are manifold, in part due to the dependence of both temperature and precipitation on complex atmospheric and oceanic circulation patterns ([33]). Therefore, a proper attribution of the respective role of different sources of climate variability on forest productivity is critically needed so as to better understand and anticipate the future evolution of the forests in a context of climate change ([4], [27]).

The ways by which climate change can impact productivity in forest ecosystems are complex. For instance, while evidence from tropical montane forests suggests a net negative impact of warming on tree photosynthetic capacity ([11]), studies realized in temperate and boreal latitudes support seemingly opposite results. In addition, scaling responses from leaf photosynthesis to stem growth are not necessarily direct; information regarding climatic drivers of growth in tropical trees thus remains a big challenge. As a result, studies based on empirical growth evidence in Neotropical forests remain scarce, mainly because of (i) the widespread lack of long-term forest data and (ii) the absence of a clear seasonality, which hampers the formation of well-defined growth rings ([48], [51]). Yet, at higher elevations in South America, tree growth is limited by temperature and so tree species can form annual rings, even at low latitudes, and thereby provide valuable insights into growth-climate relations, natural forests dynamics, and responses of trees to climate change ([51], [5]).

Natural climatic variability controls temperature and precipitation fluctuations in the Neotropics at inter-annual and decadal scales ([13]), and this variability therefore indirectly controls plant growth. Regulation of natural climate variability in South America is primarily driven by the superposition of several large-scale phenomena including the El Niño Southern Oscillation (ENSO), Pacific Decadal Oscillation (PDO) and the Antartic Oscillation (AAO - [13]). The two phases of ENSO, El Niño and La Niña, have been demonstrated to exert major control on interannual climate variability over South America, and thereby also on rainfall and temperature anomalies ([1], [33]). The PDO ([30]), on the other hand, is responsible for much of the multi-decadal precipitation variability over South America ([13]). Several studies in South America have found clear links between large-scale climatic modes and tree growth ([5]). Villalba et al. ([48]) for instance demonstrated that tree-rings from montane forests in north-west Argentina are strongly connected to large-scale circulation patterns. Mundo et al. ([31]) showed a connection between Araucaria araucana K. Koch growth and the Antarctic Oscillation (AAO) in Argentina.

Besides climate, tree growth also responds to factors at the level of individual trees and forest stands, such as tree size, tree-to-tree competition, or altitude ([43]). Among these factors, tree size is decisive when it comes to improve our understanding of ontogenetic growth patterns associated to resource use ([52], [28]). For instance, it has been demonstrated that competition for light has a strong influence on the growth of small trees, whereas competition for nutrients affects trees of all sizes ([9]).

A correct identification of the relative influence of each source of variability (i.e., climate, endogenous factors, site factors) on tree growth in Neotropical montane forests as well as a correct identification of the complex network of relationships between variables are thus crucially needed as they will allow implementation of adequate climate change adaptation and preservation strategies for these forest ecosystems that are currently at risk ([27], [28]). Structural equation modeling (SEM) has proven to be a valuable tool to explore complex data in which variables are hierarchically structured ([37], [14]). SEM allows for testing theoretical hypotheses on ecological issues and concepts in the presence of measurement errors and uncontrolled variation, where both direct and indirect effects are considered and where relationships can be simultaneously evaluated ([29]). This method was defined to analyze multiple regressions simultaneously, following a sequence specified by a causal diagram ([42], [15]). This diagram is designed by the user and based on theoretical concepts that represent the variables and their connections that are hypothesized to exist for a situation being modeled ([34]).

Here we explore the relationship between tree radial growth and climatic variability at annual timescale in a mountain forest of the Peruvian Andes. The focus of the study is on a natural and well-preserved Podocarpus glomeratus Don. forest located in the National Sanctuary of Ampay, Cordillera Vilcabamba (Peru). We hypothesize that P. glomeratus growth is responsive to local climate conditions (i.e., precipitation and temperature), which in turn are driven by modes of natural climatic variability (ENSO and PDO) and global climate change trends. To this end, we analyze tree-ring widths using a hierarchical framework in the form of Structural Equation Models (SEM) to unravel the relative contributions of each of the climate components as well as the role of other growth determinants such as tree size.

  Materials and methods 

Study site

This study has been conducted in a protected P. glomeratus forest in the province of Abancay (Apurimac region, Peru). This forest is the main protected feature of the Ampay National Sanctuary (central coordinates: 13° 33′ S, 72° 51′ W; extension: 38.5 km²; altitudinal range: 2900 to 5235 m a.s.l.) located in the Cordillera Vilcabamba, Peruvian central Andes (Fig. 1). The protected part of the forest covers approximately 6 km², representing 41% of the woodlands and 17% of the sanctuary’s territory, and extends from 3252 to 3539 m a.s.l.

Fig. 1 - The study site is located in the Central Peruvian Andes, in the Apurimac region. The P. glomeratus forest range from 3252 to 3539 m a.s.l.

  Enlarge/Shrink   Download   Full Width  Open in Viewer

The region is dominated by an Oceanic climate typically found in tropical mountains and defined as “subtropical highland variety” (Cwb) according to the Köppen-Geiger classification ([35]). The study area is characterized by a cold and dry season from April to September with some frosts between June and September; and a mild and rainy season extending from October to March. Mean temperature ranges from below 0 °C above 5000 m a.s.l. to 10-18 °C in the lower parts of the protected area. Average annual precipitation in Abancay is around 618 mm with maximum rainfall from December to March ([19]).

P. glomeratus represents one of the least studied podocarp species and its presence is restricted to small forested areas in Bolivia, Ecuador, and Peru ([54]). It appears mostly as single trees or in small patches within cloud forests and montane rainforests and can reach altitudes of around 4000 m a.s.l. ([10]). The species has been classified as “near threatened” by the IUCN Red List of Threatened Species ([12]) because the remaining P. glomeratus forests are highly fragmented and exploited for timber, and the patches are eventually converted into agricultural land.

Individuals of P. glomeratus exhibit different stem forms that can range from simple in young individuals (Fig. 2a) to very complicated, lobate growth patterns in older specimens (Fig. 2b, Fig. 2c, Fig. 2d). In other species of the Podocarpus genus, the lobate stem growth has been shown to cause wedging and missing rings ([21]). Furthermore, intra-annual density fluctuations within the earlywood have been reported to form “false rings”, which are sometimes difficult to differentiate from regular growth ring boundaries ([21]).

Fig. 2 - Tree details showing the lobate form of the P. glomeratus stems in several localities within the protected forest. Canopy cover depends on land-use types that prevailed before the strict protection of the Ampay Sanctuary. Red points indicate the localization of sampled trees.

  Enlarge/Shrink   Download   Full Width  Open in Viewer

Wood sample collection and measurements

Fieldwork was realized between September and October 2015 within the Ampay Sanctuary. A total of 80 Podocarpus glomeratus Don. trees were sampled, including specimens of different size classes. Trees were selected randomly to include as much variability as possible and to guarantee that the sample was representative of the local tree population. Wood samples were collected with Pressler increment borers (inner diameter 5.5 mm) at breast height. We extracted several increment cores per tree to take account of the multi-stem nature of Podocarpus spp. Samples were mounted on woody supports with glue and air-dried for at least one week. Sample surface preparation included polishing with sandpaper of increasing gradation. The age of trees was first estimated by counting the number of visible tree rings under a stereomicroscope. Ring widths were then measured using a digital LINTAB positioning table connected to a Leica stereomicroscope and TSAPWin scientific software at a resolution of 0.01 mm. From the sampled trees, we measured tree-ring widths from 40 suitable samples from 33 different trees. By contrast, many other samples had to be excluded from analyses as they did not exhibit continuous growth-ring records but instead contained discontinuities or were of poor wood quality.

To facilitate cross-dating of series from individual trees we used pointer years, i.e., characteristic growth signals that were common to all tree individuals ([53]). Finally, annual trunk diameter increment (d) was reconstructed by discounting tree-ring widths backwards. We then computed the corresponding Basal Area (BA) for every year by applying the formula (eqn. 1):

We express annual radial growth as a Basal Area Increment (BAI) in cm² and define it as the difference between basal areas of consecutive years ([3]). BAI values were standardized for the period 1967-2014 following procedures commonly used to obtain z-scores (i.e., subtracting the average of the data and dividing the result by the standard deviation; [8], [17]).

Climate information

In this paper, we analyzed relationships between BAI and the Southern Oscillation Index (SOI - [38]), a measure indicative of ENSO activity. Prolonged periods of negative (positive) SOI values coincide with abnormally warm (cold) ocean waters across the eastern tropical Pacific typical of El Niño (La Niña) episodes. SOI information and values were retrieved from ⇒ https:/­/­www.­ncdc.­noaa.­gov/­teleconnections/­enso/­indicators/­soi/­. In this study, we also considered the Pacific Decadal Oscillation (PDO - [30]) as it is associated with low-frequency climate variability over South America ([13]). PDO information and values were retrieved from ⇒ https:/­/­www.­ncdc.­noaa.­gov/­teleconnections/­pdo/­.

Monthly temperature and precipitation data were retrieved from the CRU TS 4.01 web site (⇒ https:/­/­crudata.­uea.­ac.­uk/­cru/­data/­hrg/­cru_ts_4.­01/­). We used yearly climatic data as; (i) differences between seasons are generally weak in terms of temperature and precipitation in tropical mountain forests; (ii) conifers have been shown to sometimes develop several growth activations within the same year ([7]). We thus calculated average values of the climatic variables from September to August (i.e., hydrological year) to adjust climatic data to the annual resolution of tree growth values. Drought years were defined by SENAMHI ([40]) using the Palmer Drought Index defined in terms of deficit of precipitation.

Statistical analysis

We applied a Structural Equation Model (SEM - [42]) to explain the hierarchical effects of (i) the main climatic modes affecting the Peruvian Andes (ENSO and PDO), (ii) climatic variables limiting growth of trees in mountain areas (precipitation and minimum temperature - [46]), (iii) forest factors (tree size), and (iv) the calendar year (as a proxy of global change trends of atmospheric CO₂ and soil deposits of N - [36]) on P. glomeratus growth (i.e., BAI) over the past 50 years. SEM analyses were realized for the period 1967-2014 for which we had the largest possible number of individuals (27 trees) covering a representative time span of 48 years.

SEM analyses were conducted in R using the “piecewiseSEM” package which allows implementation of mixed linear regressions in a SEM framework ([23]). First, we constructed a hypothetical model - represented in a causal graph - in which variables (i.e., predictors and response variables) were related by direct causal links or paths (Fig. 3). The relations between variables used to construct the causal model were based on findings from previous work (summarized in Tab. 1) and are also illustrated in Fig. S1 (Supplementary material). In this model, we assumed that climatic variability (i.e., T_min and precipitation) can be modeled as a function of climatic modes (ENSO and PDO, paths 7, 8, 9, and 10 in Fig. 3) and the calendar year (paths 1 and 2 in Fig. 3). Specifically, climate change trends were represented by linear trends in both minimum temperatures (T_min) and precipitation, which are denoted by the relationship between them and the calendar years ([36], [27], [28]). Likewise, BAI appears here as a function of climatic variability (paths 5 and 6 in Fig. 3), tree ontogenetic development (path 4 in Fig. 3), as well as the calendar year in a direct relationship (path 3 in Fig. 3). Tree-ring widths define sequential growth measurements within individual trees. Thus, we considered the tree ID as the random term, thereby assuming that measurements taken in the same individual are correlated. Moreover, tree rings represent a temporal series and we thus included autoregressive parameters in the regression model to consider potential temporal autocorrelation. We tested models differing only in the degree of autocorrelation (0 to 4 degrees of temporal autocorrelation) to account for the influence of previous years on current tree growth ([32]). We selected the best model by using the Akaike Information Criterion corrected for small sample size (AICc). In order to identify missing paths (i.e., potential direct causal relations without theoretical support in the initial model) in the hypothetical model we run a χ²-test on the Fisher’s C statistic ([23]).

Fig. 3 - Causal model used to illustrate the structural equation modelling (SEM) process. This model evaluates the hypothesis that P. glomeratus growth is correlated with local climatic conditions (i.e., precipitation and temperature), which in turn are driven by modes of natural climatic variability (ENSO and PDO) and tree endogenous factors (i.e., tree size).

  Enlarge/Shrink   Download   Full Width  Open in Viewer

  Results 

Wood samples and tree description

We measured tree-ring widths from 40 tree-ring series from 33 trees exhibiting visibly continuous growth-ring records and lacking obvious discontinuities or poor wood quality (Fig. 4). The core samples used exhibited between 48 and 127 annual growth rings, after the exclusion of false rings and the insertion of missing rings.

Fig. 4 - P. glomeratus chronology from the Ampay Sanctuary. Graph showing (i) the number of samples available to calculate the Basal Area Increment (BAI) (in green), (ii) reconstructed BAI average for the 30 samples (from 27 trees used in the SEM (in black); (iii) standardized values of BAI for the period 1967-2014 are represented with a blue line; (iv) PDO Index variability. Asterisks indicate droughts (red: severe to extreme drought; Palmer Drought Index < -3; green: moderate drought, Palmer Drought Index -2 to -3; [40]). Red arrows indicate the strong and very strong El Niño events in 1972-73, 1982-83 and 1997-98. The blue arrows indicate the strong La Niña events in 1973-74, 1975-76, 1988-89, 1998-99, 1999-00, 2007-08, and 2010-11. El Niño and La Niña events are based on the Oceanic Niño Index (ONI - ⇒ http:/­/­origin.­cpc.­ncep.­noaa.­gov/­products/­analysis_monitoring/­ensostuff/­ONI_v5.­php ).

  Enlarge/Shrink   Download   Full Width  Open in Viewer

The P. glomeratus trees sampled in the Ampay sanctuary showed clear tree-ring boundaries with rows of earlywood cells separated from one another by latewood tracheid with slightly thickened cell walls (Fig. 5a, Fig. 5b). Tree rings occasionally disappeared on samples or merged with other growth rings in other parts of the stem (Fig. 5c), a phenomenon known as wedging or incomplete rings. Furthermore, we also identified intra-annual wood density variations within the earlywood zone of some rings, often close to the previous growth ring boundary, a growth pattern known as false rings (Fig. 5b).

Fig. 5 - Micro-section of Podocarpus glomeratus stem wood. a) Clear ring boundaries forming annual bands. b) Tree ring (white arrows) with a density variation band inside (i.e. false ring, red arrows). c) Detail of a wedging tree ring.

  Enlarge/Shrink   Download   Full Width  Open in Viewer

SEM analyses were limited to the period 1967-2014 for which the largest possible number of individuals (27 trees) was available with 30 increment core samples. The average BAI of these 30 samples exhibits strong negative anomalies in 1976, 1994, and 2000, with this negative becoming a consistently negative trend since 2003 (Fig. 4).

Influence of climate, global change, and forest factors on P. glomeratus growth

Fig. 6 displays the results of the SEM. The hypothesized paths between the variables that were validated by the model are paths 2, 4, 5, 7, and 9. The climate modes ENSO and PDO are significantly correlated with T_min (Pearson’s correlation coefficients r = -0.54 and r = -0.31, respectively; see p-values in Tab. 2), which in turn had a direct effect on BAI (Pearson’s r = -0.12). Our analysis also shows that the calendar year, as a proxy of the temporal lineal trends of temperature related with the increase of atmospheric greenhouse gases, exhibits a positive significant correlation with T_min (Pearson’s r = 0.59), indicating a rise of T_min in the area which in turn has a negative indirect influence on BAI.

Fig. 6 - Structural Equation Models (SEMs). Arrows represent relationships between variables. Only significant paths (P ≥ 0.05) among all possible paths are shown. Standardized regression coefficients are shown in red besides the arrows.

  Enlarge/Shrink   Download   Full Width  Open in Viewer

On the other hand, at the level of individual trees, mean tree size shows a positive direct influence (Pearson’s r = 0.7) on growth. Neither ENSO, PDO, nor calendar years had any significant effect on precipitation. In turn, the parameter estimate associated to precipitation in the basal area increment model was not significant (path 6 in Tab. 2).

  Discussion and conclusions 

In this paper, we provide evidence indicating that the ENSO and PDO control BAI of the investigated P. glomeratus forest indirectly through their influence on minimum temperatures, and that ongoing climate warming, through its direct impact on T_min increase¸ also has an indirect negative effect on BAI. This study thereby shows how global warming is affecting this sensitive high-mountain area and ecosystem. We also provide evidence on how tree growth can be severely affected in tropical montane forest in the future.

Our results suggest that the calendar year, as a proxy of global change trends ([36]), has an indirect negative effect on BAI through the recent and ongoing increase of minimum air temperatures (T_min). A likely explanation of this outcome is the change in growth phenology in the sense that early xylem reactivation can render tissues vulnerable to freeze injury by late frost ([25]) due to problems in hydraulic transport. This idea is supported by (i) a clear and positive trend of T_min in the Peruvian Central Andes over the last decades ([22]), and (ii) an increase in frost frequency ([44]), thus suggesting a new scenario of late frost vulnerability for P. glomeratus populations in this Peruvian mountain environment.

Results also suggest a lack of influence of precipitation on BAI. From a statistical perspective, an interpretation of these outcomes can be related with the short temporary extension of the available precipitation data, therefore hampering a correct identification of the rainfall signal. At the same time, however, our tree-ring chronology is indeed sensitive to extreme meteorological droughts in some years. For instance, the strong decrease in average BAI values recorded in 1976, 1994, and 2000 are likely related to the severe droughts observed in the region in 1974-75, 1991-93 and during the drought associated with the strong El Niño in 1997-98 ([40]). In a similar way, the La Niña event that occurred between 1999 and 2001 could have influenced tree growth as well with more frequent frost events resulting from a significant temperature decrease of around 2.5 °C ([41]). Our findings indicate that the ENSO and PDO climatic modes indirectly control BAI through their influence on T_min. Indeed, in the Peruvian Andes, temperature is negatively correlated with the South Oscillation Index (SOI - [22]) at inter-annual scales, therefore establishing a link with the main ENSO phases: as such, temperature is higher during El Niño and lower during La Niña ([1], [49]). This correlation is statistically significant in our datasets with a standardized value of -0.54 (p < 0.0001), a finding that agrees with previous studies showing that the increase of observed T_mean over the region may be associated with an intensification of El Niño events since the late 1970s ([49]), which in turn would lead to tree-growth anomalies in South America ([47]). Our results also exhibit a significant correlation between PDO and T_min (standardized value of -0.32, p = 0.02). According to Garreaud et al. ([13]), PDO-related anomalies of temperature over South America are spatially similar to ENSO, but their amplitude is lower. During the period covered by this study, a climatic shift in the Pacific can be observed in 1976-77, coincident with a phase change in the PDO from negative to positive (Fig. 4). This fluctuation was associated with changes in ENSO ([45]) and can be discerned over South America in precipitation and surface temperature changes ([30]). Since 1999, the PDO has remained mostly in its negative phase, a phenomenon that was accompanied by an increased frequency of La Niña occurrences in the Pacific basin, and by only few, weak to moderate El Niño events until 2015.

Previous research has demonstrated a direct relation between rising air temperatures and the upward shift of treelines in other parts of the world ([20]). This elevational change can expose trees and plants in general to physiological drought conditions as a result of seasonal water stress and more frequent frost events ([25]). In the case of P. glomeratus, previous research in the Bolivian Andes ([2]) indicates that the species is indeed sensitive to altitudinal shifts. This sensitivity is related with the exposure of trees to warmer (colder) conditions as well as to hydric stress (through differences in solar radiation intensity) at lower (higher) altitudes. Thus, we hypothesize that an eventual shift of the P. glomeratus treeline due to rising air temperatures in the Peruvian Central Andes could have direct implications on P. glomeratus population health. Furthermore, and because of differences in the magnitude of optimum elevational shifts among coexisting species, this upward move could also trigger a disruption of biotic interactions and the ecological networks in the Ampay Sanctuary ([24]).

We conclude that over the last two decades, tree growth of P. glomeratus of the Ampay Sanctuary in the Central Andes of Peru has declined progressively, and that this decline can be related statistically to the rise of T_min. This study therefore adds new information on tree growth in Neotropical mountain forests and thereby highlights that these environments are indeed particularly vulnerable to climate change and changes in modes of natural climatic variability - as both processes can alter growth of plants by controlling local climatic conditions. Should the warming trend persist in the future, mountain environments such as Ampay Sanctuary with its severely threatened P. glomeratus forest could be seriously and irreversibly damaged, which would also have consequences for other coexisting animals and plants as well as for the human communities that are dependent upon their functions and services. Besides, a reduction of tree productivity as a consequence of climate warming would also slow down carbon mitigation, thus affecting the carbon biogeochemical cycle and ensuing influences on climate and ocean processes. In view of the endangered status of P. glomeratus and the effects of ongoing climate warming on its vitality and growth, we call for more research on its regeneration, phenology, ecology, and possible adaptation to climate change, so as to improve the protection and management of forests like the one in the Ampay Sanctuary and in a context of future climate change.

  Acknowledegments 

Authors thank the PACC project funded by the Swiss Agency for Development and Cooperation (SDC) for partial funding of this study. We thank the Peruvian National Service of Protected Areas (SERNANP), and especially the director of the Ampay National Sanctuary, Amilcar Osorio and the biologist Jaime Valenzuela Trujillo for their help during the field work. We also acknowledge the logistic assistance provided by Ljubika Indira Ruiz and her team belonging to the Peruvian Institute of Civil Defence (INDECI), as well as the important support of Sandra Paula Villacorta Chambi from the Peruvian Geological Institute (INGEMMET).

CR, JB and JM conceived and designed the research. CR and JB carried out surveys. CR measured and analysed the samples. JM helped with methodology and data analysis. CR developed and wrote the manuscript. MS performed the language revision and helped to draft the manuscript.

  References

  Supplementary Material

Fig. S1 - Correlation coefficients plots between (a-e) the Basal Area Increment (BAI) and the different variables used to construct the causal model; (f-l) the climatic modes (i.e., Pacific Decadal Oscillation and El Niño-South Oscillation) represented by their respective indices, PDO and SOI, and the other predictors; and (o-q) the calendar year and precipitation, temperature and the tree size to illustrate temporal evolution of these variables.