With increasing globalization, trade, and human movement, the rate of alien species introduction has increased all around the globe. In addition, climate change is thought to exacerbate the situation by allowing range expansion of invasive species into new areas. Predicting the distribution of invasive species under conditions of climate change is important for identifying susceptible areas of invasion and developing strategies for limiting their expansion. We used Maxent modeling to predict the distribution of one of the world's most aggressive invasive weeds, Ageratina adenophora (Sprengel) R. King and H. Robinson, in the Chitwan–Annapurna Landscape (CHAL) of Nepal under current conditions and 3 future climate change trajectories based on 3 representative concentration pathways (RCPs 2.6, 4.5, and 8.5) in 2 different time periods (2050 and 2070) using species occurrence data, and bioclimatic and topographic variables. Minimum temperature in the coldest month was the most important variable affecting the distribution of A. adenophora. About 38% (12,215 km2) of the CHAL area is climatically suitable for A. adenophora, with the Middle Mountain physiographic region being the most suitable one. A predicted increase in current suitable areas ranges from 1 to 2% under future climate scenarios (RCP 2.6 and RCP 8.5). All protected areas and 3 physiographic regions (Siwaliks, High Mountain, High Himalaya) are likely to gain climatically suitable areas in future climate scenarios. The upper elevational distribution limit of the weed is expected to expand by 31–48 m in future climate scenarios, suggesting that the weed will colonize additional areas at higher elevations in the future. In conclusion, our results showed that a vast area of CHAL is climatically suitable for A. adenophora. Expected further range expansion and upslope migration in the future make it essential to initiate effective management measures to prevent further negative impacts of this invasive plant.
Biological invasions, a major driver of global environmental changes, are posing serious threats to global biodiversity and ecosystem functioning (IPBES 2019). Climate change is likely to further amplify the risks of biological invasions (Walther et al 2009; Bradley, Wilcove, et al 2010). Biological invasions and climate change act synergistically, and this synergistic relation between the 2 parameters of global change has been identified as a major threat to biodiversity (Dukes and Mooney 1999; Walther et al 2009; Mainka and Howard 2010). The profound negative impacts of invasive plant species on the diversity of native species, soil dynamics, and ecosystem processes, which cause ecological and economic losses, are well known (Marbuah et al 2014; Villa and Hulme 2017; Castro-Diez et al 2019). These impacts are likely to be further exacerbated by climate change, enhancing traits that promote invasiveness and creating a more hospitable climate for invasive species to cross geographic barriers, thereby facilitating range expansion in new areas (Dukes and Mooney 1999; Stacowicz et al 2002; Walther et al 2009; Bradley, Blumenthal, et al 2010; Bellard et al 2013). In comparison to native plant species, parameters of global change, such as increased temperature and CO2 enrichment, enhance the performance of invasive species, imposing a threat of further spread (Liu et al 2017). Therefore, predicting the distribution of invasive weeds under climate change scenarios and identifying the areas potentially at risk are urgent needs for effective management planning to minimize ecological and economic impacts.
A first step to identify the risk of invasions is to use ecological niche models (ENMs) to predict suitable ecological niches for a species across a landscape. These relate documented presence records of the focal species with the environmental or spatial characteristics of the potential sites (Elith and Leathwick 2009; Franklin 2009). The niche concept is central to ENMs and is based on Hutchinson's (1957) concept of fundamental and realized niches (Araujo and Guisan 2006). It is highly likely that an invasive species at an early stage of invasion occupies only a small fraction (ie realized niche) of the fundamental niche in the introduced range, and there is always a risk of invasion in the unoccupied part of the fundamental niche (Soberon and Nakamura 2009). However, controversies persist over which facets of the niche are projected by ENMs (Araujo and Guisan 2006; McInerny and Etienne 2012). ENMs have been gaining popularity and are widely used by ecologists in invasive species risk assessments (Qin et al 2016; Suarez-Mota et al 2016; Wan et al 2017; Shrestha et al 2018; Thapa et al 2018). Among different ENMs, Maxent is one of the most popular species distribution modeling tools. This model uses presence-only records and has been commonly used in building habitat suitability maps for invasive species (Phillips et al 2006; Merow et al 2013; West et al 2016; Lamsal et al 2018).
Among 124 countries, Nepal has the third highest threat to agriculture sectors from invasive species (Paini et al 2016). To date, 179 species of flowering plants are naturalized, and, among them, 26 species are reported to be invasive in Nepal (Shrestha 2019; Shrestha, Budha, et al 2019). The crofton weed, Ageratina adenophora (Sprengel) R. King and H. Robinson (Asteraceae), is one of the most noxious invasive weeds in many parts of Asia, Oceania, and Africa. It has had serious ecological impacts on native biodiversity and caused enormous economic losses (Poudel et al 2019). It is ranked as the most problematic invasive weed in Nepal (Tiwari et al 2005).
Though many studies have been carried out to investigate the potential distribution of A. adenophora on broader spatial scales in Nepal (Shrestha and Shrestha 2019), China (Wang and Wang 2006; Wang et al 2017), the Himalayas (Lamsal et al 2018; Thapa et al 2018), and South Africa (Tererai and Wood 2014), there is a lack of such studies on smaller scales where management strategies are implemented. One of the most important landscapes in Nepal is the Chitwan–Annapurna Landscape (CHAL), located in central Nepal. CHAL harbors rich biodiversity due to its wide elevation gradient (200–8091 m above sea level [masl]), diverse topography, and climatic variations (subtropical to alpine) (WWF 2013; MFSC 2016). Invasion by alien species has already been recognized as a major threat to biodiversity in CHAL (WWF 2013), and A. adenophora is the invasive weed most prioritized by the local communities for management in natural ecosystems due to its negative impacts on biodiversity and livelihoods (Shrestha, Shrestha, et al 2019). Therefore, there is an urgent need to recognize potential areas of distribution of A. adenophora in CHAL under the current climate and identify areas at risk of being invaded by this weed under future climate scenarios.
In this study, we used the Maxent modeling tool to predict the current and future potential distribution of A. adenophora in CHAL using occurrence records from different sources. The objective was to prepare habitat suitability maps for the weed under current climatic conditions and future climate scenarios (RCP 2.6, RCP 4.5, and RCP 8.5 in the years 2050 and 2070) to identify the key environmental factors influencing its distribution and areas at risk of invasion. Information on its potential distribution will be very useful for the scientific community and managers in developing future monitoring and management strategies to prevent further expansion of the weed in this landscape.
CHAL is located in central Nepal and covers 19 districts. This landscape has a wide elevation gradient, ranging from 200 to 8091 masl, and covers an area of 32,057 km2 (WWF 2013). It spans 4 physiographic regions, namely, Siwalik, Middle Mountain, High Mountain, and High Himalaya. They have diverse climatic conditions, from subtropical in Siwalik to alpine in the High Himalaya, and a cold and dry climate in Trans-Himalayan regions. This geographic and climatic diversity shapes the habitat and environmental conditions for CHAL's rich biodiversity, which includes more than 104 species of mammals (Bhuju et al 2007), 500 species of birds (Baral and Inskipp 2005; Bhuju et al 2007), and 3430 species of plants, with high levels of endemism and genetic diversity (BPP 1995). Forests and grasslands are the main natural ecosystems, occupying 35.5% and 8.6% of the landscape, respectively, whereas 21.1% of the area is under agriculture (WWF 2013). The region has a population of 4.5 million people (CBS 2013). The average minimum and maximum temperatures are 5°C and 40°C, and the average annual rainfall ranges from 165 to 5244 mm (MFSC 2016). The landscape includes portions of 4 globally recognized ecoregions and comprises 3 national parks (Chitwan, Parsa, and Langtang) and 2 conservation areas (Annapurna and Manaslu) (MFSC 2016). Annapurna Conservation Area and Chitwan National Park are among the sites with a high number of visiting tourists (DNPWC 2018). Most of the lower and mid-hill forests of this region are at risk of fragmentation and conversion to other vegetation types due to climate change (Thapa et al 2015). In comparison to the eastern and western regions of Nepal, central Nepal, where CHAL is located, hosts higher numbers of naturalized plant species (Bhattarai et al 2014). The combination of diverse natural environment along with anthropogenic disturbances has made the region vulnerable to biological invasions (WWF 2013).
A. adenophora is native to Mexico but is established in 40 countries outside its native range (Poudel et al 2019). In Nepal, 26 invasive plant species are reported, 20 of which have been documented in the CHAL region, where A. adenophora is the most problematic weed in natural ecosystems, and its management is highly prioritized by local people (Tiwari et al 2005; Shrestha 2019; Shrestha, Shrestha, et al 2019). It was first reported in 1958, having been accidentally introduced from the eastern border of India to Nepal (Tiwari et al 2005). It covers a wide elevation gradient of 400–3280 masl in Nepal (Siwakoti et al 2016). The ability of A. adenophora to occupy a wide range of climatic habitats and spread rapidly can be attributed to its phenotypic plasticity, allelopathy, and ability to alter the soil microbial community to favor its further invasion (Poudel et al 2019). It has reduced forage supply, displaced native plant species, causing loss of biodiversity, and prevented forest regeneration in CHAL (WWF 2013; Shrestha, Shrestha, et al 2019).
Species occurrence data
In total, 686 occurrence points were collected from different sources. We noted 245 occurrence points from secondary sources (Siwakoti et al 2016; Shrestha and Shrestha 2019). The remaining 441 occurrence points were collected by the first author during field visits in 2016–2019. Road networks are the major conduit for dispersal of A. adenophora (Dong et al 2008); therefore, occurrence data were mainly collected through field surveys along roadsides and trekking routes (Figure 1). This sampling bias was addressed by spatial filtering of the data. Duplicate records of occurrence points were deleted and spatially thinned using the spThin package (Aiello-Lammens et al 2015) in the R software (version 3.4.4) (R Core Team 2017), so that only 1 single location occurred in each 1 km2 grid cell. Spatial filtering makes it possible to reduce overfitting to sampling bias in ENMs (Boria et al 2014). In total, 403 occurrence records obtained after filtering were used to build the models (Figure 1).
Environmental and bioclimatic variables
In December 2018, we downloaded 19 grid-based bioclimatic variables that represent annual trends, seasonality, and extreme climatic conditions from the WorldClim database (version 1.4) ( www.worldclim.org; Hijmans et al 2005) at a spatial resolution of 30 arc-seconds (∼1 km2) ( Appendix S1 (mred-40-02-15_s01.pdf), Supplemental material, https://doi.org/10.1659/MRD-JOURNAL-D-19-00069.1.S1). Elevation was obtained from the Shuttle Radar Topographic Mission (SRTM) at 90 m spatial resolution. This was then resampled into 30 arc-second spatial resolution by using the nearest neighbor resampling technique in ArcGIS (version 10.3). Slope and aspect rasters of the study area were derived from the elevation data.
To predict climatically suitable areas in future climate scenarios, we chose projections from the Community Climate System Model (CCSM4) under b1 emission scenarios, which are based on the fifth phase of the Coupled Model Intercomparison Project5 (CMIP5) (Gent et al 2011). We selected 3 greenhouse gas (GHG) emission scenarios, also known as Representative Carbon Pathways (RCP 2.6, RCP 4.5, and RCP 8.5), for 2 different time periods (2050 and 2070) as adopted by the Intergovernmental Panel on Climate Change in its Fifth Assessment Report (AR5) (IPCC 2013). RCP 2.6, RCP 4.5, and RCP 8.5 represent the lowest, medium, and highest emission scenarios, corresponding to a 1.0°C, 1.4–1.8°C, and 2.0–3.7°C projected increase in global mean surface temperature, respectively (van Vuuren et al 2011; IPCC 2013).
The datasets were extracted for the study area (CHAL) using the Spatial Analyst Tool and the Extraction Tool in ArcGIS (version 10.3). These datasets were converted from raster format to ASCII files in ArcGIS. We repeated the procedure to prepare the predicted distribution maps for the 2 future climate scenarios (2050 and 2070).
Maxent (version 3.3.3) software was used to build the ENMs (Phillips and Dudik 2008). Maxent is one of the most commonly used habitat suitability modeling techniques. It uses presence-only data and is widely used for invasive species (Phillips et al 2006; Shrestha et al 2018; Maharjan et al 2019). Maxent works well with incomplete or limited data, so it can provide robust estimates of potentially suitable habitats for invasive species at small spatial scales (Jarnevich et al 2006; Jarnevich and Reynolds 2011; West et al 2016).
To reduce multicollinearity among predictor variables (19 bioclimatic and 3 topographic variables) and overfitting of the model, pairwise correlation analyses were performed in R, and highly correlated variables with a Pearson's correlation coefficient ≥0.8 were removed ( Appendix S2 (mred-40-02-15_s01.pdf), Supplemental material, https://doi.org/10.1659/MRD-JOURNAL-D-19-00069.1.S1) (Merow et al 2013). While selecting one variable from each pair of highly correlated variables (r2 > 0.8), special attention was paid to existing biological and ecological insights into the species. Ultimately, 7 bioclimatic variables—isothermality (Bio 3), minimum temperature of the coldest month (Bio 6), temperature annual range (Bio 7), precipitation of the driest month (Bio 14), precipitation seasonality (Bio 15), precipitation of the warmest quarter (Bio 18), and precipitation of the coldest quarter (Bio 19)—and 2 topographic variables—aspect and slope—were used as predictors to build the habitat suitability model.
The Maxent model used 75% of the data for training and the remaining 25% for testing. We used a logistic format because it improves model calibration by estimating the probability of a species being present depending on environmental variables (Phillips and Dudik 2008). A convergence threshold of 10–5, a maximum iteration value of 5000, 15 replications with a replicated run type subsample, and 10,000 random background points were used to build the model. As a threshold rule, we chose tenth percentile training presence on the basis of the area under the curve (AUC) and true skill statistic (TSS). Tenth percentile training presence omits the 10% of localities or training presence records with the lowest predicted values (Radosavljevic and Anderson 2014) and is highly conservative in estimating species tolerance with respect to each climatic variable (Svenning et al 2008). The remaining parameters were kept at their default values.
We imported the Maxent output, which is continuous data with values ranging from 0 to 1, into ArcGIS (version 10.3) and classified the map using the Reclassify Tool into 2 classes, suitable habitat and unsuitable habitat, on the basis of a tenth percentile training presence logistic threshold. In this way, a binary habitat suitability map was created for the current and all future climate scenarios. We also calculated the climatically suitable area for present and future climate scenarios, as well as changes in suitable areas in terms of gain, loss, and stable areas in the future under all scenarios in ArcGIS. Data for physiographic regions and protected areas were then clipped to projected maps with suitable areas for current and future climate scenarios to calculate the suitable areas in these physiographic regions and protected areas. Changes in the upper and lower elevational distribution range under future climate scenarios in comparison to current were quantified using the Extraction Tool and digital elevation model (DEM) raster for all maps, current and future, in ArcGIS.
Threshold-independent (area under the receiver operating characteristic [ROC] curve [AUC]) and threshold-dependent (TSS) measures of model accuracy were used to evaluate model performance (Fielding and Bell 1997; Allouche et al 2006; Franklin 2009). AUC values range from 0–1.0, with 0.5–0.7 considered low, 0.7–0.9 moderate, and >0.9 high (Swets 1988; Manel et al 2001). The TSS value ranges from –1 to +1, where +1 indicates a perfect agreement, and 0 or less indicates a performance no better than random (Allouche et al 2006). Marginal response curves were used to visually investigate the relationship between environmental variables (predictors) and the predicted index of habitat suitability of A. adenophora. The relative contribution of different predictor variables to the Maxent model was assessed by the variable percentage contribution and jackknife procedures (Elith et al 2011). The jackknife test of variable importance helps to identify those variables with important individual effects (Elith et al 2011). Two jackknife tests were taken into account. The jackknife test of regularized training gain shows the training gain of a variable when used in isolation and the training gain of a variable when omitted, and it compares these values to the training gain of all variables. Similarly, the jackknife test of AUC based on the AUC of test data shows the predictive performance of the variable when used in isolation and the predictive performance of the variable when omitted, and it compares these values with the AUC value when all variables are used (Phillips 2017).
Model performance and variable contribution
The current model for A. adenophora performed better than random, with a mean training AUC value of 0.85, mean test AUC value of 0.80, and a TSS value of 0.52 (Table 1). Responses of each predictor variable are shown in the response curve ( Appendix S3 (mred-40-02-15_s01.pdf), Supplemental material, https://doi.org/10.1659/MRD-JOURNAL-D-19-00069.1.S1). Out of 9 predictor variables used for model building, minimum temperature of the coldest month (Bio 6) contributed the most (48.7%), followed by precipitation of the warmest quarter (Bio 18). Aspect had the lowest contribution (2.9%) (Table 2). Maxent's jackknife test of variable importance also showed that minimum temperature of the coldest month (Bio 6) had the highest training gain and AUC, followed by precipitation of the warmest quarter (Bio 18) when used in isolation (Figure 2). The response curve of the variable minimum temperature of the coldest month (Bio 6) showed that the probability that the weed will occur below 1°C was the lowest, increasing with increasing minimum temperature ( Appendix S3 (mred-40-02-15_s01.pdf), Supplemental material, https://doi.org/10.1659/MRD-JOURNAL-D-19-00069.1.S1). In addition, the minimum temperature of the coldest month of about 96% of the occurrence points used in model building was above 0°C ( Appendix S4 (mred-40-02-15_s01.pdf), Supplemental material, https://doi.org/10.1659/MRD-JOURNAL-D-19-00069.1.S1). Similarly, the response curve of precipitation of the warmest quarter (Bio 18) indicated that the probability the weed will occur increased with increasing precipitation in the warmest quarter above 500 mm ( Appendix S3 (mred-40-02-15_s01.pdf), Supplemental material, https://doi.org/10.1659/MRD-JOURNAL-D-19-00069.1.S1).
Model evaluation matrices.
Relative contribution of the environmental variables to the Maxent model built for current climatic conditions.
Current potential distribution
Currently, 38% (12,215 km2) of the total area of CHAL is suitable for A. adenophora (Figure 3; Table 3). All districts of CHAL within the elevational range of 119–2824 masl had climatically suitable areas for this weed. The most suitable areas for A. adenophora were found to be in the Middle Mountain physiographic region (75%), followed by Siwalik (37%), High Mountain (29%), and High Himalaya (0.1%) (Table 4). High Himalaya is hardly suitable for the weed, with the fewest suitable areas found in the districts of Manang and Mustang. Nonetheless, all districts of the CHAL region were found to have climatically suitable areas. Because the Middle Mountain and High Mountain regions had the most suitable areas for the weed, we tried to observe the change in climatically suitable areas in 3 protected areas situated in these physiographic regions. Among the 3 protected areas—Annapurna Conservation Area, Langtang National Park, and Manaslu Conservation Area—Langtang National Park was predicted to have the highest percentage of area suitable for this weed (Figure 4).
Predicted climatically suitable area for A. adenophora under current and future climate scenarios in the Chitwan–Annapurna Landscape, Nepal.
Predicted suitable area for A. adenophora in different physiographic regions of the Chitwan–Annapurna Landscape, Nepal.
Future invasion risk and change in habitat suitability
The predicted climatically suitable areas for A. adenophora would increase under RCP 2.6 for the year 2070 and RCP 4.5 for both the years 2050 and 2070 (Figure 3; Table 3). The highest increase (2%) in area of suitable habitat was predicted for the year 2070 under RCP 2.6 and 4.5. This gain in suitable areas was more prominent in districts like Lamjung, Gorkha, Dhading, Makwanpur, Chitwan, and Tanahun. However, an increase in radiative force (from RCP 4.5 to 8.5) would decrease climatically suitable areas for the weed in both the years 2050 and 2070 (Figure 3; Table 3).
Though the suitable area was predicted to decrease in extreme climate scenarios (RCP 8.5), the upper elevational distribution limit would expand by 31 m and 42 m for the years 2050 and 2070, respectively (Figure 5). Though a small decrease (24 m) in upper elevation limit was predicted under RCP 2.6 in 2050, the highest increase of 48 m was expected under RCP 2.6 for the year 2070 in comparison to current climatic conditions. However, under the medium emission scenario of RCP 4.5, the model predicted a contraction in the upper elevational limit for both years 2050 and 2070. In contrast, the lower elevation limit of A. adenophora would either remain stable or contract in future climate scenarios (Figure 5).
With climate change, all physiographic regions except the Middle Mountain region were expected to gain climatically suitable areas. Though a minimal loss in suitable areas was predicted for Middle Mountain, this region will still contain the most suitable areas for the weed, followed by High Mountain, Siwalik, and High Himalaya, in all future climate scenarios (Table 4). In 4 of the future climate scenarios—RCP 2.6 in 2050 and 2070, and RCP 8.5 in 2050 and 2070—Siwalik will gain climatically suitable areas. Except under RCP 8.5 for the year 2050, High Mountain will also gain suitable areas in all future climate scenarios. Furthermore, under RCP 4.5 and 8.5 for the year 2070, the climatically suitable area is also predicted to increase in High Himalaya. Among all physiographic regions, the percentage gain in suitable area was highest (2.4%) in Siwalik under RCP 8.5 for the year 2050 (Table 4).
Like under current climatic conditions, in future climate scenarios, Langtang National Park will have more suitable areas than the other 2 protected areas (Figure 4). Under RCP 2.6 and 4.5 for both 2050 and 2070, climatically suitable areas are predicted to increase in Annapurna Conservation Area, whereas under extreme climate scenarios, it will lose some suitable areas. For Langtang National Park, climatically suitable areas will increase under RCP 2.6 and 8.5 for both years, but not under RCP 4.5 (likewise for both years). For Manaslu Conservation Area, a remarkable gain in suitable areas is predicted only under RCP 8.5 for the year 2070.
This study is the first to predict current and future suitable habitat for A. adenophora in CHAL, Nepal. The model evaluation parameters (AUC and TSS) obtained for our models both lie within a range that confirms the robustness of the models (Table 1): AUC values above 0.8 and TSS values closer to 1 are considered to be acceptable (Swets 1988; Manel 2001; Allouche et al 2006). Though our study did not include a model transferability assessment, this would provide valuable information for model validation as well as model selection (Wenger and Olden 2012).
Climatic factors, such as temperature and precipitation, play a pivotal role in determining the pathways and success of plant invasions (Kathiresan and Gualbert 2016; Wang et al 2017). Our predictions suggest that the minimum temperature in the coldest month is the most influential factor for the distribution of A. adenophora. This is in line with the findings of other researchers (Wang and Wang 2006; Zhu et al 2007; Wang et al 2017; Lamsal et al 2018; Thapa et al 2018; Datta et al 2019). According to our model, the minimum temperature in winter (Bio 6) alone contributes about 49% to the habitat suitability model (Table 2), and the probability of occurrence of the weed decreased to almost 0 below the minimum winter temperature of 1°C ( Appendix S3 (mred-40-02-15_s01.pdf), Supplemental material, https://doi.org/10.1659/MRD-JOURNAL-D-19-00069.1.S1). This indicates that the occurrence of this weed in high mountains is constrained by the lower minimum winter temperature, which is supported by observations in India, where the low temperature in winter limits the uppermost distribution range of A. adenophora (Datta et al 2017).
Our ENMs were built based on the occurrence data of invaded regions only, which indicates that the potential niche that we estimated is part of the fundamental niche of the species (Elith and Leathwick 2009; Soberon and Nakamura 2009). Moreover, factors such as biotic interactions and dispersal limitations also restrict the species from occupying its full potential niche (Soberon and Nakamura 2009). Thus, use of occurrence data from both native and invaded ranges would provide a more accurate estimate of the potential niche of the species (Jimenez-Valverde et al 2011).
In current climatic conditions, as well as in future climate scenarios, Middle Mountain is found to have more climatically suitable areas for A. adenophora compared to other physiographic regions. The elevation of the Middle Mountain region ranges from 1000 to 2500 masl (DHM 2017), which lies within the suitable range for the distribution of this weed (Wang and Wang 2006; Zhu et al 2007). A recent study that modeled the distribution of 24 invasive alien plants in Nepal also identified the Middle Mountain region as having a particularly large area of invasion hotspots, with suitable areas for the greatest number of species studied (Shrestha and Shrestha 2019). Displacement of native species, such as Artemisia indica and Urtica dioca, and reduction of the ground vegetation layer of Digitaria sp., Eragostris sp., and Imperata cylindrica by A. adenophora have already been observed in the Middle Mountain region (Tiwari et al 2005; Baral et al 2017). Therefore, this region requires the urgent attention of policymakers and land resource managers to implement effective management plans to prevent further spread of this weed. In contrast, the High Himalaya region is unsuitable for the weed because this region has a low minimum winter temperature (–10°C to 5°C) and low annual precipitation (400–1000 mm) (DHM 2017), which limit the distribution of the weed (Datta et al 2017). Most of the current areas predicted to be suitable for A. adenophora were found along road verges and river networks, which is in accordance with a study conducted in China (Wang and Wang 2006). These river and road networks enhance the rapid spread of A. adenophora (Dong et al 2008; Sang et al 2010). Therefore, control and monitoring efforts for the management of this weed should be focused on areas near rivers and roads.
There is growing evidence that climate change is likely to increase the risk of plant invasions, creating more suitable areas in the future (Bradley, Wilcove, et al 2010). Our projections also showed an increase in climatically suitable areas for A. adenophora in the future. The weed will gain suitable areas in 3 future climate scenarios: RCP 2.6 (2050) and RCP 4.5 (2050 and 2070). A recent study also reported that climatically suitable niches for A. adenophora would expand by 5.3% under RCP 6.0 in Nepal (Shrestha and Shrestha 2019). Studies in other parts of the world have also predicted the expansion of climatically suitable areas for A. adenophora in future climate scenarios (Wang and Wang 2006; Zhu et al 2007; Wang et al 2017; Lamsal et al 2018; Thapa et al 2018). However, a study conducted in Bhutan predicted a contraction of suitable areas for A. adenophora by 0.22% in 2050 (Thiney et al 2019). An experimental warming study conducted in China revealed that warming (2°C rise) increased biomass allocation and canopy cover of the weed, making it more stress tolerant (He et al 2012). Similarly, CO2 enrichment was also reported to increase the relative growth rate and biomass allocation of the weed (Lei et al 2012). Thus, an increase in suitable areas in future climate scenarios with increased temperature and CO2 concentration can be attributed to the innate and evolutionarily increased ecophysiological tolerances of A. adenophora favoring its growth (Blossey and Notzold 1995; He et al 2012; Lei et al 2012). Although, under extreme climate scenarios (RCP 8.5), a loss in total climatically suitable areas is predicted, the weed will still thrive and gain suitable areas in protected habitat and physiographic regions. Districts like Lamjung, Gorkha, Dhading, Makwanpur, Chitwan, and Tanahun are predicted to be vulnerable to further invasion by the weed due to climate change, so it is crucial to implement scientifically informed site-specific management policies, with the participation of local communities.
Despite its apomictic nature, with the associated evolutionary constraints, niche expansion has been observed in A. adenophora (Datta et al 2019). In addition, the weed is found to exhibit phenotypic plasticity that helps it to occupy a broader climatic niche (Zhao et al 2012). Our future climate models also demonstrate an expansion of the upper elevational distribution limit of the weed. Our results confirm findings in the Western Himalaya, where the weed was predicted to expand its upper elevational limit by 981 m compared to current climatic conditions (Thapa et al 2018). This indicates that the weed will spread toward cooler and drier regions in future. A similar trend was observed in a study of spatiotemporal patterns in China (Zhu et al 2007). Increased cold tolerance due to epigenetic modifications might help the weed to gain more suitable habitats in cooler and drier places at high elevations (Xie et al 2015). With climate change, all 3 protected areas will gain areas of potential suitable habitat. Though the weed has already been identified as the most problematic weed in Annapurna Conservation Area, impacting native diversity and livelihoods (Thapa and Maharjan 2014), no such studies have been carried out in Langtang National Park and Manaslu Conservation Area. Upward movement and colonization of A. adenophora due to recent climate changes have already been observed in Langtang National Park (Lamsal et al 2017). Our model also predicted that Langtang National Park has the highest proportion of climatically suitable areas in current climatic conditions and will continue to do so in future climate scenarios. A range shift of A. adenophora might threaten the habitat of 2 endangered animals, Ailurus fulgens (red panda) and Moschus chrysogaster (Himalayan musk deer), in Langtang National Park (Lamsal et al 2017). Thus, this information should act as a prompt for land managers, the scientific community, conservationists, and policymakers to adopt precautionary measures and formulate effective policies to prevent the further spread of this weed into new regions.
Taking climatic and topographic factors into consideration, our model predicted elevational range expansion, as well as an increase in suitable areas in future climate scenarios. However, other factors, such as biotic interactions, dispersal ability, demography, evolution, adaptation, and land-use change, also play key roles in determining the species range shift with climate change (Sinclair et al 2010; Urban et al 2016). Furthermore, long-term temporal predictions in climate change scenarios are associated with 2 other main errors. The first is the extrapolation of data beyond the training range to nonanalogue environmental conditions (climate change scenarios), which might make predictions unreliable (Fitzpatrick and Hargrove 2009). Another risk is that, as the relationship between climatic variables could change with time, the correlation structure of current and future climate variables might also change, thus reducing the certainty of models (Dormann et al 2013). To overcome these uncertainties and achieve more realistic predictions of species distribution, all factors, abiotic and biotic, that shape the distribution of invasive species should be incorporated in the species distribution models (Gonzalez-Salazar et al 2013; Leach et al 2016). Nevertheless, ENMs provide predictive information on species distribution required by vegetation managers and conservation practitioners for developing effective strategies to prevent further invasion (Peterson 2003). Thus, our study provides useful information about the current distribution of A. adenophora and identifies areas that may be at risk in the future on a local scale, demonstrating the urgent need for formulating effective management strategies to mitigate the impact of the weed.
Conclusions and management recommendations
Our study suggests that A. adenophora could spread further under future climate scenarios while retaining most of the currently suitable areas. Among the 4 physiographic regions in CHAL, the Middle Mountain region currently has the highest proportion and the High Himalaya region currently has the lowest proportion of climatically suitable areas for A. adenophora, and this is projected to continue under future projected climate scenarios. Similarly, Langtang National Park contains a higher percentage of areas suitable for the weed compared to the other 2 protected areas under current and future climate scenarios. All protected areas will gain additional suitable areas in future climate scenarios. Furthermore, it is predicted that the weed will expand its distribution range to higher elevations in future climate scenarios; this will amplify the consequences of climate change, which is already impacting these areas. For physiographic regions (Siwalik and High Himalaya) and protected areas (Manaslu Conservation Area) that have few suitable areas for invasion, regular inspection of habitats is needed to allow effective action to be taken in time to prevent further expansion of the weed. Management strategies for smaller and accessible invaded areas could either be mechanical control, for example, by hand pulling, or chemical control, by using herbicides (Parsons and Cuthbertson 2001; Di Tomaso et al 2013). Local communities in Nepal use A. adenophora for cattle bedding, composting, and also for making bio-briquettes (Tiwari et al 2005; Baral et al 2017; Shrestha, Shrestha, et al 2019); these cultural control methods can also be employed to manage the weeds, as long as precautions are taken to prevent seed dispersal. However, for regions like Middle Mountain and protected areas like Langtang National Park that have large areas at risk of invasion, effective management options, such as an integrated weed management approach, must be adopted and implemented along with regular monitoring of suitable habitats. By identifying areas that are potentially at risk in the future, our study constitutes a helpful resource for managers and policymakers to take appropriate and timely action to minimize the risk of invasion by A. adenophora associated with climate change.
This work was funded by the U.S. Agency for International Development (USAID) Bureau of Food Security under the Cooperative Agreement No. AID-OAA-L-15-00001 as part of the Feed the Future Innovation Lab for Integrated Pest Management. Any opinions, findings, conclusions, or recommendations expressed here are those of the authors alone. Some of the species location data used in this study were gathered by BBS during fieldwork supported by the International Foundation for Science (IFS, Sweden) and National Trust for Nature Conservation (Nepal). We express our gratitude to Dr DN Shah for his support during climatic data preparation. Thanks go to Ms Sara Hendery for English language editing in the manuscript. Comments and suggestions by anonymous reviewers greatly helped to improve the manuscript.
APPENDIX S1 (mred-40-02-15_s01.pdf) Bioclimatic variables used for modeling suitable habitat of Ageratina adenophora.
APPENDIX S2 (mred-40-02-15_s01.pdf) Correlation matrix of 19 bioclimatic and 3 topographic variables.
APPENDIX S3 (mred-40-02-15_s01.pdf) Marginal response curves of 2 predictor variables.
APPENDIX S4 (mred-40-02-15_s01.pdf) Scatter plot of occurrence points and minimum temperature of the coldest month for current climate.