Introduction

The ranges of mountain plant communities are sensitive to climate change (Telwala et al 2013). As treeline ecotones invade high elevations, alpine tundra biomes shrink or shift upslope. Previous work has shown treelines are distributionally linked to temperature (Gehrig-Fasel et al 2007; Körner 2012), with new trees prevented from establishment above the upper treeline limit by unfavorable conditions. As climate has warmed, many treelines around the planet have shifted upward (Harsch et al 2009; Gaire et al 2014; Greenwood and Jump 2014). Nevertheless, not all sites have responded in the same way, with some species and regions remaining stable or actually retreating (Harsch et al 2009). Within regions already threatened by habitat change, estimating the degree to which future climatic shifts may influence habitat suitability for treeline-forming species is essential for both scientific understanding and conservation planning. Loss of alpine tundra to treeline advance (Moen et al 2004) reduces the habitat for endangered tundra species. For instance, the endangered snow leopard in the Himalayas could lose as much as 30% of its habitat under future warming scenarios and experience an increase in competition from other cat species (Forrest et al 2012; Lovari et al 2013). This combination of factors (loss of habitat and increased competition) has the potential to severely reduce Himalayan biodiversity.

The Himalayas contain an expansive area of understudied treeline habitat that is likely to respond to ongoing climate change. They also contain a steep elevational gradient, complex topography, high warming rate, biodiversity hotspots, and proximity to indigenous communities (Xu and Grumbine 2014). These features make the Himalayas a complex environment where treeline-related studies carried out elsewhere cannot be generalized. Therefore, studies focusing on the distribution of treeline species in the Himalayas are needed to understand and conserve Himalayan ecosystems (Singh et al 2013). Previous studies of treelines in this area have indicated mixed treeline response to climate change over a limited area (Gaire et al 2014; Shrestha et al 2014; Chhetri and Cairns 2015; Suwal et al 2016). The degree to which these trends scale up to broader patterns remains unresolved. Few modeling efforts have been conducted to identify contemporary climatic limits of Himalayan treeline-species ranges and predict distributions under potential future climate regimes (Schickhoff et al 2015). Recently, modeling studies of the environmental niches of Betula utilis in Uttarakhand, India (Singh et al 2013) and the Himalayan range (Schickhoff et al 2015; Bobrowski et al 2017) have been carried out, both indicating potential habitat shift. We wish to add to this knowledge by examining 3 representative Himalayan treeline-forming species from Nepal to gauge the future of this region.

In this study, we used species distribution modeling to determine if the distributions of these 3 species are defined by climate and predict whether their ranges are likely to expand or contract under projected climate conditions through an increase or decrease in suitable habitat. The study area is the country of Nepal, which covers the majority of the Himalayan range and has undergone significant shifts in climate within the last 50 years (Shrestha and Devakota 2010). The 3 study species, Abies spectabilis (D. Don), Betula utilis (D. Don), and Pinus wallichiana (A. B. Jacks), are important components of the Nepalese subalpine ecosystem and are dominant within the Nepalese Himalayan treeline (Chhetri et al 2017). We posed 2 questions:

Material and methods

Study area

Nepal is a mountainous country and occupies the central part of the Himalayas (Figure 1). Analysis of climate records has indicated that temperature has increased at a higher rate in the Nepalese Himalayas than in other mountain areas, especially since the 1950s (Shrestha and Devkota 2010). This warming has progressed at a steady rate since the mid-1970s and is more pronounced at higher elevations. In addition, precipitation (rain, snow, sleet, or hail) has increased by 13 mm per year on average, while the number of rainy days has decreased by 0.8 days per year, suggesting more intense rainfall (Shrestha et al 2000). Several climate models have predicted a steady increase in temperature throughout Nepal accompanied by a decrease in monsoon rainfall in the north and an increase in the south (NCVST 2009; Shrestha and Aryal 2011).

FIGURE 1 

Distribution of the 3 study species.

Results

We developed SDMs for 3 treeline species on the basis of 240 presence records. Omission and predicted area curves indicated moderate to good model fit for the 3 study species (Supplemental material, Figure S1:  http://dx.doi.org/10.1659/MRD-JOURNAL-D-17-00071.S1). All of the observed omission rates on training samples and test samples were close to the predicted omission rates. The models for all 3 species performed better than random, as indicated by the average test AUC (area under the curve) values for the 15 replicate runs of A. spectabilis (0.89 ± 0.05), B. utilis (0.87 ± 0.05), and P. wallichiana (0.87 ± 0.06) (Supplemental material, Figure S2:  http://dx.doi.org/10.1659/MRD-JOURNAL-D-17-00071.S1). The results indicated that current suitable habitat for A. spectabilis, B. utilis, and P. wallichiana cover 5130 km², 9822 km², and 9764 km², respectively (Figure 2; Table 1). The model also predicted that the most suitable habitat for A. spectabilis and B. utilis was in the eastern and central part of Nepal, and the most suitable habitat for P. wallichiana was in the central part of Nepal.

FIGURE 2 

Suitable habitat for the 3 study species under current and projected climate scenarios.

TABLE 1 

Area of suitable habitat for the 3 treeline species under the current climate and 3 future climate scenarios. Numbers in parentheses indicate percent of Nepal's total land area; numbers outside parentheses represent blocks of about 100 km².

A Jackknife test (Supplemental material, Figure S3:  http://dx.doi.org/10.1659/MRD-JOURNAL-D-17-00071.S1) indicated that for all 3 species elevation is the most important predictor of habitat suitability. The climatic variables isothermality (Bio3) and mean temperature of the coldest quarter (Bio11) were consistently important predictors for all 3 species. The percentage of each variable's contribution to species distribution in the model is presented in Table 2. Elevation, isothermality, mean temperature of the coldest quarter, precipitation of the coldest quarter (Bio19), and topographic position index (TPI) had the most influence on the distribution of A. spectabilis. Elevation, mean temperature of the coldest quarter, isothermality, eastness, and precipitation of the warmest quarter (Bio18) had the most influence on the distribution of B. utilis. The distribution of P. wallichiana depended primarily on elevation, isothermality, topographic position index, mean temperature of the coldest quarter, and temperature annual range (Bio7). Response curves showing relationships between the top 3 topographic-climatic variables and species occurrence probability are presented in Figures S4–S6 (Supplemental material,  http://dx.doi.org/10.1659/MRD-JOURNAL-D-17-00071.S1).

TABLE 2 

Contribution of different topographic-climatic variables to the distribution of the 3 study species.^a)

The proportion of suitable habitat for the 3 study species varied under the projected climate scenarios (Figure 2). The model predicted that suitable habitat area for A. spectabilis will increase under all 3 scenarios in 2050 and 2070. The increase in area is projected to be greater in 2070 than in 2050, except under the RCP 6.0 scenario, where a dip in suitable habitat is predicted in 2050 (Table 1). The average elevation of A. spectabilis was higher under all warming scenarios (Table 3; Figure 3). The model predicted that the suitable habitat for B. utilis will increase in both area (Table 1) and elevation (Table 3; Figure 4) under all scenarios in 2050 and 2070. The model predicted that the suitable habitat for P. wallichiana will not significantly decrease or increase except under the RCP 6.0 scenario in 2070 (Table 1), but the average elevation of P. wallichiana was higher in all warming scenarios (Table 3; Figure 5).

TABLE 3 

Elevation (masl, average ± standard deviation) of suitable habitat for the 3 study species under the current climate and 3 future climate scenarios. We compared all future elevation estimates to current values using an independent sample t-test; all comparisons were significant (P < 0.001).

FIGURE 3 

Change in average elevation of suitable habitat for A. spectabilis under projected climate scenarios.

FIGURE 4 

Change in average elevation of suitable habitat for B. utilis under projected climate scenarios.

FIGURE 5 

Change in average elevation of suitable habitat for P. wallichiana under projected climate scenarios.

Discussion

The distribution of suitable habitat for treeline species in the Nepalese Himalaya is likely to shift in response to climate change. The distribution of the 3 species examined here was mostly explained by climatic variables, indicating that climatic shifts will likely influence future habitat suitability. Explicit analysis confirmed this conclusion. Our simulations indicated that the suitable habitat for the 3 treeline species will shift toward higher elevations under predicted future climate conditions and may undergo an overall range expansion. These predicted shifts have serious consequences for both forest and tundra species in this area.

The contemporary SDMs performed well in describing current patterns when examined through observational and statistical analysis. All model AUC values were close to 0.90, which is considered a moderately good fit (Kramer-Schadt et al 2013) and comparable to values obtained in similar studies in this region (Rhododendron spp. AUC = 0.78, Kumar 2012; B. utilis AUC = 0.92, Schickhoff et al 2015). Our model results matched the actual existing range of these species-based fieldwork data and existing land cover maps and ecological maps.

The primary topographic-climatic variables describing the contemporary models for our 3 study species were meaningful given existing knowledge of the ecology of treeline species. Elevation was the most important predictive variable in the model for all 3 of our species. Elevation may not be directly associated with plant physiology, but it plays an important role in controlling atmospheric pressure, solar radiation, precipitation, and cloud cover (Oke and Thompson 2015), and thus it is a strong indicator of climatic variables that influence physiology. Mean temperature of the coldest quarter was one of the most important climatic variables explaining the contemporary distribution models for all 3 study species. This finding is consistent with ecological studies that have found that the upper elevational distribution limit of many treeline species is often determined by low temperatures (Körner 2012). Additionally, mean temperature of the coldest quarter and isothermality are related to growing season length, which has been identified as a secondary factor limiting treeline distribution (Körner 2012).

We should remain conservative in interpreting these findings, as we removed variables that were highly correlated with those included in our model, and simple spatial autocorrelation is an inherent effect in any climate study (Kumar 2012). Nevertheless, the strong climatic associations we observed do suggest, in the absence of other limitations, that ongoing climate change will influence the distribution of these treeline species.

Although the shift in suitable habitat distribution varied under the projected climate scenarios, future treeline advance is likely for the 3 study species. In terms of change in the area of suitable habitat, our models predicted minimal levels of expansion for A. spectabilis and P. wallichiana relative to B. utilis, which will double or triple under RCP 6.0 and RCP 8.5 scenarios. This is likely due to differences in species-specific requirements that define their distributions. P. wallichiana prefers dry valleys and foothills and is an early successional species (Ghimire et al 2011). In contrast, both A. spectabilis and B. utilis grow on moist north-facing slopes (Ghimire and Lekhak 2007). A. spectabilis has a broader distribution than B. utilis, often also occurring on south-facing slopes. B. utilis is often found in areas that receive snowmelt from mountain peaks above the treeline (Shrestha et al 2007).

Future distribution range will depend on how these species-specific microsites will change. General circulation models project warmer days and nights in the future, which will result in more snowmelt, which benefits B. utilis more than A. spectabilis and P. wallichiana. Western Nepal is drier than eastern Nepal, and monsoon precipitation is projected to increase in western Nepal, which will make western parts more suitable for tree establishment. Our model projected that B. utilis will occupy most of western Nepal under both RCP 6.0 and RCP 8.5 scenarios. Nevertheless, we saw average treeline extent advancing under all future climatic scenarios. Many climate models—including the Long Ashton Research Station Weather Generator (Agarwal et al 2014), Hadley Centre Coupled Model Version 3, and Regional Climate Model Version 3 (NCVST 2009)—project an increase in Himalayan winter temperatures, which will prolong the growing season, making the area above the treeline more favorable for treeline species.

These results are similar to those of previous studies examining species in this area. Schickhoff et al (2015) predicted that B. utilis habitat will shift northward throughout all of the Himalayas under future climate scenarios. Zomer et al (2014) predicted that subalpine conifer forest zones will shift upward by over 400 m between 2000 and 2050. On the other hand, Kumar (2012) predicted that the habitat of Rhododendron spp. will shrink under future climate scenarios in the Sikkim region of the Indian Himalayas, so not all species will necessarily follow this same path.

These results should be considered within the limitations of using a topographic-climatic approach to predict future species distributions (Macias-Fauria and Johnson 2013; Zong et al 2014). Limitations not addressed in this analysis may prevent these species from occupying the entirety of their fundamental niche in the future. These could include dispersal limitations of each species (Boisvert-Marsh et al 2014), which may limit their ability to establish in newly available suitable habitat. In addition, our projections do not factor in human disturbance (timber and fuelwood harvesting), herbivory (cf. Cairns and Moen 2004), plant physiology, soil type, snow cover, or land cover. Human disturbance and topographical factors can limit upslope treeline advance (Leonelli et al 2009).

Shrestha and Bawa (2014) emphasized the need for high-resolution environmental data that capture microclimates, edaphic conditions, vegetation dynamics, and landscape heterogeneity in SDMs. Kollas et al (2014) also emphasized the use of high-spatial-resolution temperature data for predictive modeling of temperature-based niche envelopes, recommending the use of topographic-climatic variables with a resolution of less than 100 m. They also suggested the use in modeling of absolute minimum temperatures instead of long-term means, because minimum temperatures determine the phenology of tree species at the cold limit. Building on higher-resolution datasets in future studies will improve the accuracy of these results. High-resolution climate data are not available for Nepal. Nevertheless, this study provides results based on the best available data, isolating climatic influences from other ecological limitations, and providing essential information for future environmental management.

Conclusions

Our model produced a good fit of contemporary species distributions, identifying suitable habitat for 3 dominant treeline-forming species of the Nepalese Himalayas under present and potential future climates. This approach could be usefully applied to other treeline species in the region. Our results indicate that the treeline ecotone is likely to transition throughout this region, which will likely have significant impacts on the associated plant and animal species. Our models predicted that an area above the existing treeline will become suitable for tree establishment. Nevertheless, this establishment will be controlled by factors like natural or human disturbances and ecological interactions with the surrounding shrub communities. Future work examining disturbance factors and species interaction, and incorporating high-resolution satellite imagery and a DEM, will improve the accuracy of this research.