An assessment of human impacts on endangered red pandas (Ailurus fulgens) living in the Himalaya

Abstract Anthropogenic factors play an important role in shaping the distribution of wildlife species and their habitats, and understanding the influence of human activities on endangered species can be key to improving conservation efforts as well as the implementation of national strategies for sustainable development. Here, we used species distribution modeling to assess human impacts on the endangered red panda (Ailurus fulgens) in high‐altitude regions of Nepal. We found that the distance to paths (tracks used by people and animals), livestock density, human population density, and annual mean temperature were the most important factors determining the habitat suitability for red pandas in Nepal. This is the first study that attempts to use comprehensive environmental and anthropogenic variables to predict habitat suitability for the red pandas at a national level. The suitable habitat identified by this study is important and could serve as a baseline for the development of conservation strategies for the red panda in Nepal.

To protect red panda habitat, managers need broad-scale geographic information. While numerous studies have been conducted to assess habitats, conservation threats, and diets of red pandas at local scales in Nepal (Bista et al., 2017;Bista, Panthi, & Weiskopf, 2018;Panthi et al., 2012Panthi et al., , 2015Panthi et al., , 2017Thapa & Basnet, 2015), few studies have investigated the species distribution and threats to their habitat at national and regional scales (Acharya et al., 2018;Kandel et al., 2015;Thapa et al., 2018). Anthropogenic factors play an important role in shaping the distribution of wildlife species and their habitats (Lewis et al., 2017), and understanding the influence of human activities on endangered species can be key to improving conservation efforts as well as the implementation of national strategies for sustainable development. Although the red panda is facing serious anthropogenic pressure (Acharya et al., 2018;Glatston et al., 2015;Panthi et al., 2017;, previous studies did not thoroughly consider anthropogenic factors when modeling the habitat of this species (Kandel et al., 2015;Thapa et al., 2018). Consequently, anthropogenic impacts on the red panda and its habitat remain unclear, and a comprehensive assessment of the suitable habitat for red pandas in Nepal is not available. Due to insufficient information on red panda habitat at large spatial scales, conservation partners such as the government of Nepal, World Wildlife Fund, National Trust for Nature Conservation, and Red Panda Network have been unable to prepare effective policies, plans, and strategies for red panda conservation in Nepal. In this study, we aim to assess human impact on endangered species living in high-altitude regions in Nepal by using the red panda as an example. Our specific objectives are to (a) quantify suitable habitat for red pandas across Nepal; (b) determine the role of anthropogenic factors to predict suitable habitat for red pandas. The information from this study will be useful for the government of Nepal and conservation partners to prepare and implement policies, plans, and strategies for immediate and long-term conservation of red panda in Nepal.

| Study area
Nepal is situated in the central part of the Himalaya and covers an area of 147,181 km 2 . Nepal has diverse climates due to the large variation in elevation, varying from tropical lowlands in the south to alpine cold semi-desert in the trans-Himalayan zone (Ohsawa, Shakya, & Numata, 1986). The average annual rainfall is around 1,000-2,000 mm, but sometimes it exceeds 3,000 mm in some lower parts of the country (Ichiyanagi, Yamanaka, Murajic, & Vaidyad, 2007). Nepal has diverse geography ranging from very rugged and permanently snow and ice-covered Himalayan Mountains in the north to tropical alluvial plains in the south. Due to variation in climate and topography, Nepal is classified into five physiographic zones (i.e., Terai, Siwalik, middle Mountain, high Mountain, and Himalaya; Barnekow Lillesø, Shrestha, Dhakal, Nayaju, & Shrestha, 2005;Shrestha, Shrestha, Chaudhary, & Chaudhary, 2010). In spite of economic obstacles, the government of Nepal has established 20 protected areas that cover more than 23% of the total land area of the country: 12 national parks, six conservation areas, one wildlife reserve, and one hunting reserve ( Figure 1) (DNPWC, 2017). These protected areas provide natural habitat for elephant, musk deer, red panda, rhino, snow leopard, tiger, wild buffalo, wild dog, and other threatened wildlife (DNPWC, 2017).

| Red panda occurrence data
We compiled two datasets including 30 first-hand and 295 secondhand red panda occurrence records ( Figure 1). The second-hand occurrence records were obtained from published research articles and unpublished government reports of Nepal. All second-hand data were collected between 2009 and 2016 using a Global Positioning System (GPS). The sources of these second-hand data are listed in Appendix 1. Based on the spatial distribution of the second-hand data, we interviewed a number of red panda experts and local park rangers to identify other potential red panda habitats for primary data collection. We carried out fieldwork in September and October 2017 in Langtang National Park, Ilam, Panchthar, and Dhading districts of Nepal. In the field, the direct and indirect signs of red pandas (i.e., droppings) were recorded using a GPS by adopting the purposive sampling.

| Topographical variables
A digital elevation model (DEM) with a spatial resolution of 1 km was downloaded from the USGS website (https ://earth explo rer.usgs. gov/; USGS/EarthExplorer, 2017), and the slope and aspect were derived from the DEM using ArcGIS software (ESRI, 2017).

| Vegetation-related variables
Satellite-derived normalized difference vegetation index (NDVI) is a commonly used vegetation index for ecological research. In this study, we used the NDVI time series to model red panda habitat.
Since most of the secondary red panda occurrence data were collected between 2009 and 2013, we downloaded atmospherically corrected 10-day composite NDVI images with a spatial resolution of 1 km over the same period (180 images, three images per month) acquired by SPOT4 and SPOT5 Vegetation (VGT) sensor from the  (Bjørn, 2009) European Space Agency product distribution portal (http://www. vito-eodata.be; Vito, 2017). We smoothed these NDVI images using an adaptive Savitzky-Golay filter in TIMESAT (Jönsson & Eklundh, 2004). The seasonal characteristics of five full phonological cycles were constructed based on the five years' time series NDVI data and statistical products (i.e., maximum, mean, minimum, standard deviation, and amplitude). The resulting smoothed data were used as environmental variables in our model. The forest cover data for the region were obtained from Advance Land Observing Satellite (http://www.eorc.jaxa.jp/ALOS/en; JAXA EORC, 2017). In addition, forest canopy height data with a 1-km spatial resolution was obtained from the Spatial Data Access Tool (see https ://webmap.

| Human population density
Human population density with a spatial resolution of 1 km was downloaded from the socio-economic data and application center (http://sedac.ciesin.colum bia.edu; CIESIN, 2000).

| Livestock density
Livestock (cattle, goat, and sheep) density with a spatial resolution of 1 km was obtained from the Center for Earth Observation and Citizen Science (see https ://www.geo-wiki.org)" (Robinson et al., 2014).

| Distance to roads
Road networks were downloaded from the Geofabrik website (http://downl oad.geofa brik.de/asia/nepal.html; OpenStreetMap Contributors, 2017). We then generated a raster file of the distance to roads with a spatial resolution of 1 km using ArcGIS (ESRI, 2017).

| Distance to paths
Path (tracks used by people and animals) networks were downloaded from the Geofabrik website (http://downl oad.geofa brik.de/ asia/nepal.html; OpenStreetMap Contributors, 2017). We then generated a raster file of the distance to paths with a spatial resolution of 1 km using ArcGIS (ESRI, 2017).

| Distance to human settlements
Settlement points throughout Nepal were obtained from the Department of Survey, Nepal. A raster layer of distance to human settlements with a spatial resolution of 1 km was created using ArcGIS (ESRI, 2017).

| Land cover and land use
Land use and land cover with a 1-km spatial resolution were obtained from the Fine Resolution Observation and Monitoring Global Land Cover website (FROM-GLC) (http://data.ess.tsing hua.edu.cn; Li et al., 2016).
F I G U R E 2 Correlation matrix of environmental and anthropogenic variables. Cool colored (blue) squares indicate a positive correlation and warm colored (red) squares indicate a negative correlation; darker colored squares indicate stronger correlation and paler colored squares indicate a weaker correlation

| Multicollinearity analysis
Removing the highly correlated (|r| > .70) variables for species distribution models is recommended for reliable and unbiased output (Braunisch et al., 2013;Dormann et al., 2013). We used ArcGIS to extract the values of these variables at species presence points (ESRI, 2017) and conducted a multicollinearity analysis between these variables using the 'mctest' package in R (R Core Team, 2018) ( Figure 2). Finally, 18 highly correlated variables were removed from the dataset, and the remaining 17 variables were used for habitat modeling (Table 1).

| Ecological niche model
The maximum entropy (MaxEnt) model is one of the most reliable and robust model for species distribution and habitat suitability modeling (Phillips, Anderson, & Schapire, 2006). In addition, built-in jackknife tests in the program allow users to estimate the significance of individual variables in computing the habitat suitability (Elith et al., 2006). We used the MaxEnt program version 3.4.0 (https ://github.com/mrmax ent/Maxent) to develop environmental niche models. In this study, no primary and secondary data of red panda occurrence points were reported from two physiographical regions of Nepal: Terai and Siwalik. Therefore, these two regions were excluded from the current study to reduce modeling bias. The

| Model scenarios, evaluation, and statistical analysis
We ran the model with two different scenarios to assess the impact of anthropogenic variables on red panda habitat prediction. First, we ran the model using only environmental variables.
Next, we ran the model using both environmental and anthropogenic variables. Assessment of prediction accuracy is essential to validate the models and to understand model performance. We randomly selected fifty percent of the species occurrence points for training and used the other fifty percent to test both models.
To evaluate the accuracy of the model predictions, we used both threshold-independent and threshold-dependent methods. For the threshold-independent method, the area under the receiveroperator curve (AUC) of models was reported (Phillips et al., 2006;Wiley, McNyset, Peterson, Robins, & Stewart, 2003). The higher the AUC, the higher the model performance was. An AUC < 0.7 indicates poor model performance, 0.7-0.9 indicates moderate performance, and >0.9 indicates excellent performance (Pearce & Ferrier, 2000). Although AUC is a commonly used model evaluation parameter, it is influenced by the geographic extent of the models (Lobo, Jiménez-valverde, & Real, 2008). Therefore, we also used the threshold-dependent method, that is, true skill statistic (TSS) to evaluate the accuracy of the model predictions (Allouche, Tsoar, & Kadmon, 2006;Merow, Smith, & Silander, 2013). True skill statistic was calculated for all model outputs (0-9 replications), and the final TSS was averaged from all 10 replicates. We tested the accuracy of the 10 replicates and found that they were normally distributed for all models (Shapiro-Wilk test, p = .05).
TA B L E 1 Environmental and anthropogenic variables used for modeling the red panda habitat suitability

| Predicted suitable habitat with and without the use of anthropogenic variables
The model based on the environmental variables identified a total of 18,193 km 2 of suitable habitat for red pandas in Nepal (Figure 3).
The model based on both environmental and anthropogenic variables identified a total of 13,781 km 2 of suitable habitat for red pandas throughout Nepal (Figure 4)  Note: For each model scenario, the AUC and TSS were given as the average values of ten replicates. Superscript letters indicate significant differences among the means of AUC and TSS. Different superscript letters indicate significant differences at p < .05 (T test).

| Variables affecting red panda habitat suitability at a national level
Analysis of the contribution of environmental and anthropogenic variables to the predictive model indicated that distance to paths, annual mean temperature (Bio1), livestock density, and human population density were the most important variables contributing to the prediction of suitable red panda habitat in Nepal ( Figure 5). It is notable that among these top four variables, three of them are anthropogenic variables. We also found that variables such as the canopy height, land use and land cover, standard deviation of NDVI, distance to roads, slope, aspect, temperature seasonality (Bio4), and the precipitation of driest month (Bio14) barely contributed to the prediction of suitable habitat for red pandas at a large spatial scale in Nepal. The remaining five variables, including forest cover, distance to settlements, NDVI minimum, mean diurnal range (Bio2), and annual precipitation, had a moderate contribution to the model prediction.
The response curves of the top four variables contributing to the prediction of red panda habitat ( Figure 6) indicate that the optimal habitat for red pandas occurred in areas where the mean annual temperature (Bio1) was between 5°C and 10°C (Figure 6a). The probability of suitable habitat for red pandas increased with increasing distance to the nearest paths, but decreased dramatically after approximately 2 km from the paths (Figure 6b). The relationships between red panda habitat suitability and livestock density and human population density were negative (Figure 6c,d). An increase in livestock density, as well as human population density, significantly reduced habitat suitability for red pandas.

| D ISCUSS I ON
We successfully predicted suitable habitat for red pandas in Nepal using both environmental and anthropogenic variables. Our results show that three out of the four top predictor variables are anthropogenic factors, that is, the distance to paths, livestock density, and human population density, which all have a negative impact on red panda habitat suitability. Nepal is famous for tourism and several tourist routes, and paths for recreational trekking have been constructed in the high-altitude regions near red panda habitat. In Nepal, many local people also live in very high mountains (Chidi, 2009). These people manage the facilities for tourists and use these paths for their daily livelihood such as fuelwood and forest products collection. If the flow of local people and tourists increase, the negative impact of human paths may increase significantly in the near future. In addition to tourism, livestock is an important source of cash income for farm households in the high mountains of Nepal.
However, a number of local-level studies have reported that livestock grazing has a negative impact on red pandas (Acharya et al., 2018;Yonzon & Hunter, 1991b). This is part of a larger trend of livestock grazing contributing to biodiversity loss around the world (Alkemade et al., 2013). For example, in China, free-ranging livestock consumes considerable amounts of bamboo,  (Hull et al., 2014;Li, Pimm, Li, Zhao, & Luo, 2017). Livestock grazing also has a negative impact on grouse populations worldwide (Dettenmaier, Messmer, Hovick, & Dahlgren, 2017). Similarly, our study revealed that the high livestock density has a significant negative impact on red panda habitat at a large spatial scale in Nepal.
Biodiversity is facing serious anthropogenic impacts and is declining rapidly throughout the world (Maxwell et al., 2016;Tittensor et al., 2014). There is growing evidence that human population growth is a major cause of wildlife loss (WWF, 2018 We also recommend promulgating legislation to allow livestock in meadows but not the forest with understory bamboo, and to prohibit the collection of fodder and fuelwood from core habitat of red panda to manage the local people and wildlife in a win-win situation.
In our study, we used both environmental and anthropogenic variables to achieve a more accurate and reliable prediction of suitable habitat for red panda. We estimated that approximately 13,800 km 2 of suitable habitats are available for red pandas in Nepal, which is significantly lower than the previous studies conducted by Kandel et al. (2015) and Thapa et al. (2018), who Reserve (Panthi et al., 2012(Panthi et al., , 2015(Panthi et al., , 2017, and Rara National Park, Nepal .
In this study, we predicted that there is suitable red panda habitat has inside 13 protected areas of Nepal, but we found that of suitable red panda habitat in Shivapuri Nagarjun National Park.
As this park is the closest protected area to Kathmandu, the capital city of Nepal, this could also help attract wildlife tourists.
This study identified suitable habitat for red panda in patches of varying size. In addition to conserving large habitat patches, restoring the unsuitable area around small habitat patches and improving habitat quality is recommend for long-term conservation of the red panda. Similar to the recommendation of Bista et al. (2019), we recommend preparing and implementing site-specific conservation plans to conserve this species and its habitat. Although this study only considered a single species, we showed that wildlife of the Himalayan region faces anthropogenic pressure. Conservationists should pay more attention to this region for the conservation of specific species and overall biodiversity. In the future, researchers should also identify the impacts of other factors like climate and land use change on red pandas.
The modeling was done with presence only data, so this study couldnot account the imperfect detection of the species. We are not modeling the probability of occurrence of red pandas but rather an index of their habitat suitability, due to the lack of absence data. We used only one sample (presence point of red panda) from one grid having one-km resolution to lessen spatial autocorrelation.

ACK N OWLED G M ENTS
We thank the Netherlands Fellowship Programme (NFP) for provid-

CO N FLI C T O F I NTE R E S T
None declared.

AUTH O R CO NTR I B UTI O N S
S.P. and T.W. conceived the project and designed the study. S.P., T.W., and A.T. collected the occurrence points. S.P., T.W., and Y.S.
analyzed data interpreted the results. S.P. and T.W. wrote the manuscript. All authors critically reviewed the manuscript.