After optimization using all possible combinations of resistance layers, we generated 21 models for each genetic distance. Model comparisons using AICc did not converge on the same results for the DNei and DPCA distances. Elevation+roads was identified as the best multi-surface model using DPCA distance, whereas elevation+roughness was the best-supported model using DNei (Table S8). Considering the small delta AICc (< 2) and comparable marginal R2 of the second-best model, we believe the performance of the elevation+roughness+roads model is as good as the top model using DNei. The bootstrap analysis using DNei supports the elevation+roughness model as being the best, followed by the elevation+roads model. In contrast, bootstrap analysis using DPCA distance indicated the elevation model as the best, followed by the elevation+roads model (Table S9). Finally, model comparison using RCM also generated inconsistent results between genetic distances. In this case, the elevation+roughness+roads model was superior to all alternative models using DNei, whereas the Landuse+human density model was selected as the best using DPCA (Fig. S10). Overall, elevation was the most significant single variable, regardless of the genetic distance used. Models including elevation and roads exhibited the highest support according to AICc and the second highest support based on bootstrap analysis, using both types of genetic distance.   
 
Prediction of suitable habitats and diversity loss under climate change
After eliminating highly correlated factors, we implemented six factors to construct species distribution models: three precipitation-related (BIO6, BIO8, BIO9) and three temperature-related (BIO13, BIO14, BIO19) variables. Sixty models were generated with high predictive accuracies with an AUC of 0.97 (SD = 0.005). The response curves and relative importance (Fig. 5a-b) showed that the occurrence of P. bengalensis was fundamentally limited by the precipitation of the driest month (BIO14), the precipitation of the coldest quarter (BIO19), and the mean temperature of the wettest quarter (BIO8). Using the consensus predictions (TSS = 0.91, SD = 0.032), suitable habitat for P. bengalensis will have shrunk by 2070 under both emission models, with the predicted contraction being more pronounced according to the model of increased global warming (RCP 8.5) (Fig. 5c). In 2070, habitats at lower altitudes and southern mountainous regions are less suitable. A reduction in highly suitable habitats further to the north and in rugged mountainous areas is also predicted. There were no significant differences in any of the diversity indices among groups (all with P > 0.05) according to the two values of habitat suitability we used (Fig. 5d).