The combination of genomic offset and niche modelling provides insights into climate change-driven vulnerability

Anthropogenic climate change is a primary driver of environmental change and biodiversity loss, leading to range shifts, population declines, and extinctions. Traditional modelling of species’ responses to climate change has focused on abiotic and biotic environmental changes using ecological niche models, with limited integration of ecological genomics that captures genotype–climate associations and potential genetic changes required under future climates. Recent ecological genomic studies reveal population-specific genotype–climate associations, implying that intraspecific variation may shape differential responses to climate change. However, ecological niche models often assume uniform responses across populations. There is a pressing need to incorporate intraspecific genomic variation into niche modelling for finer-scale assessments of climate vulnerability, especially in mountainous systems characterized by strong environmental heterogeneity and local adaptation. The Sino-Himalayan Mountains (including the Southwest Mountains and East Himalayas) are biodiversity hotspots that face rapid climatic change. Two montane bird species, the Green-backed tit (Parus monticolus) and Elliot's laughingthrush (Trochalopteron elliotii), occupy different elevational and ecological zones across this region and show evidence of local adaptation. This study integrates ecological genomics and niche modelling to evaluate population-specific responses to future climate, jointly considering genomic offset (mismatch between current genotype–climate association and future climates) and niche suitability change. A genome–niche index identifies populations least disrupted by climate change as potential donors for evolutionary rescue, and landscape genetic analyses evaluate rescue routes and connectivity. The framework underscores the value of combining genomic offset, niche suitability modelling, and landscape connectivity to guide conservation under climate change.

Previous work has shown that climate change drives widespread ecological impacts, including local extinctions and range shifts, and that ecological niche models (ENMs) can project future distributions and climatic suitability. However, ENM projections can face challenges under novel climates and often ignore intraspecific adaptive genetic variation. Ecological genomics has begun to bridge this gap by mapping genotype–environment associations and predicting genomic vulnerability (genomic offset) to future climates. Studies have demonstrated that incorporating adaptive genetic variation can reduce projected range losses and predict climate-driven population declines. In mountainous regions, strong topographic and climatic gradients drive local adaptation and diversification in birds, suggesting potential for population-specific responses to climate change. Prior phylogeographic studies of the focal species indicate population isolation in heterogeneous landscapes and long-term refugia in parts of the Sino-Himalayan Mountains, providing context for assessing population-level genomic adaptation and climate vulnerability.

Study area and species: The ranges of Trochalopteron elliotii and Parus monticolus span the Sino-Himalayan Mountains (Southwest Mountains and East Himalayas) and mid-elevation mountains in Central China. Ecological zones include the southern Tibetan zone (cold-dry alpine meadow/shrubland), the southwest mountainous zone (vertical zonation from subtropical broadleaf to alpine meadows), the eastern Himalayan zone (humid evergreen forests to meadows), and the western mountainous plateau zone. Sampling and genome assembly: A de novo 10X Genomics linked-read assembly was generated for T. elliotii (1.12 Gb; contig N50 95 Kb; scaffold N50 2.702 Mb; BUSCO completeness: 90% complete single-copy, 2% duplicated, 4.5% fragmented). Protein-coding genes (17,585) were annotated using homology-based methods. Whole-genome resequencing was conducted for 55 T. elliotii and 58 P. monticolus individuals (average coverage 19.14× and 15.48×, respectively). Reads were aligned to the de novo T. elliotii genome and the great tit (Parus major) genome. Variant calling used GATK HaplotypeCaller and Samtools with stringent filters (biallelic SNPs; minQ > 30; depth 7–1000; max-missing 5; SNPs ≥5 bp from indels). SNP discovery yielded ~10.3 million SNPs in T. elliotii and ~3.9 million in P. monticolus. Genomic offset with gradientForest: From WorldClim, 19 bioclimatic variables (2.5 arc-min) were used. A random subset of 50,000 SNPs with MAF > 0.05 was analysed, building 500 regression trees per SNP. SNPs with R² > 0 were considered predictive and used to estimate variable importance and transform climatic space into genomic space (allele frequency turnover). Top, uncorrelated variables (Pearson r < 0.7) were selected. Genomic offset was computed as the Euclidean distance between current and projected genomic compositions for future climates. Future climate scenarios used four GCMs (CCSM4, MIROC5, MPI-ESM-LR, CNRM-CM5-2) under RCP4.5/8.5 for 2050 and 2070. Offsets were averaged across GCMs. Differences among climate-tolerant groups were tested via Wilcoxon rank-sum with FDR correction. Genomic offset with GDM: To address small sample sizes at some localities, generalized dissimilarity modelling (GDM) related pairwise FST (from 100,000 randomly sampled SNPs) to selected, uncorrelated climatic variables (e.g., BIO1, BIO2, BIO4, BIO7, BIO14 for T. elliotii; BIO1, BIO2, BIO4, BIO5, BIO12, BIO14 for P. monticolus). Variable importance used scaled I-splines; projections to future climates followed gradientForest procedures. Genotype–climate association and outlier SNPs: Three complementary approaches (LFMM with K=3; RDA; dbRDA) identified outlier SNPs associated with the top climate variables. P values were FDR-adjusted for LFMM; RDA/dbRDA outliers were SNPs with loadings > 3 SDs. SNPs identified by all three methods were retained. SNP annotation used SnpEff and pathway enrichment with KOBAS. Genetic clustering and climate-tolerant groups: Admixture (v1.3) identified K=3 as optimal (smallest CV error). Individuals with >60% ancestry were assigned to cold-dry, warm-dry, or warm-humid tolerant groups. Ecological niche modelling (ENM): Ensemble modelling in Biomod2 combined MaxEnt, GBM, GAM, and MARS. Occurrence records came from museums and GBIF; spatial thinning at 10 km reduced sampling bias (Moran’s I). After filtering and removing admixed localities for T. elliotii, 292 T. elliotii and 619 P. monticolus records remained, partitioned by tolerant groups (T. elliotii: 175 cold-dry, 49 warm-humid, 68 warm-dry; P. monticolus: 284 cold-dry, 155 warm-humid, 180 warm-dry). Autocorrelated variables were removed (r > 0.7). Model tuning used ENMeval (feature classes L/Q/H; RM 0.5–5). For each group: 10,000 background points; prevalence 0.5; 5-fold CV with 70/30 split; models with AUC < 0.8 or TSS < 0.6 were excluded; remaining models were ensemble-weighted by TSS. Projections were restricted to analogous climates. Future projections used the same four GCMs and RCPs (2050, 2070). Niche suitability change (NSC) = future suitability − current suitability; Wilcoxon tests with FDR compared groups. Genome–niche index (GNI): To integrate genomic offset (go) and NSC (nsc), only areas with increasing future suitability (nsc > 0) were considered. A weighted index gni = nsc^α × go^(1−α) was derived by minimizing deviation via an artificial bee colony algorithm (ABCoptim), normalizing nsc and go to 0.1–0.9. Optimal α under RCP8.5 2050 was 0.52725 for T. elliotii and 0.617253 for P. monticolus. Landscape genetic analysis: Resistance surfaces for elevation, slope, elevation SD, habitat suitability, and land cover were generated in Circuitscape. Resistance for continuous variables was based on absolute deviation from the group’s average; land cover assigned low resistance to forest/shrubland, higher to savanna, cropland, barren, and water (sensitivity tests validated assignments). Pairwise FST (100,000 SNPs) among localities was modeled with MLPE linear mixed-effects models (lme4), using z-transformed predictors, AICc-based model selection, Akaike weights, and marginal/conditional R² (MuMIn). Demographic inference: Fastsimcoal2 compared three models (no gene flow; secondary contact; continuous gene flow) using 2D folded SFS per group pair, with 100 replicates of 100,000 simulations each, mutation rate 3.3e−9 per site per generation. Best models were selected by AIC, with bootstrap to estimate confidence intervals.

• De novo genome assembly for T. elliotii: 1.12 Gb, contig N50 95 Kb, scaffold N50 2.702 Mb; 17,585 protein-coding genes; BUSCO completeness 90% (single-copy), 2% duplicated, 4.5% fragmented. • Whole-genome resequencing: 55 T. elliotii (19.14×) and 58 P. monticolus (15.48×), identifying ~10.3M and ~3.9M SNPs, respectively. • Climate variables most associated with genomic variation: P. monticolus top five were BIO3, BIO18, BIO9, BIO19, BIO5; T. elliotii top five were BIO2, BIO10, BIO7, BIO19, BIO4. • Outlier SNPs associated with climate (identified by LFMM, RDA, dbRDA): 72 in T. elliotii and 798 in P. monticolus; 25 and 204 of these located in coding/promoter regions across 23 and 147 genes, respectively. Functional enrichment highlighted catalytic and metabolic processes; shared climate-associated gene CRB1 showed strong genotype–climate correlations (T. elliotii: r = 0.795 with BIO2; P. monticolus: r = 0.835 with BIO3), suggesting roles in temperature adaptation. • Genetic clustering (Admixture, K = 3) defined climate-tolerant groups: T. elliotii — 33% cold-dry, 15% warm-humid, 16% warm-dry; P. monticolus — 10% cold-dry, 40% warm-humid, 48% warm-dry. • Genomic offset (gradientForest and GDM) under multiple GCMs and RCPs consistently showed higher offsets in western range populations (southern Tibetan zone for T. elliotii; eastern Himalayan zone for P. monticolus) than in eastern/southern regions. Cold-dry tolerant groups had significantly greater genomic offsets than warm-humid and warm-dry groups (Wilcoxon tests, FDR-adjusted P < 0.001 across datasets and methods). • Ecological niche modelling (ensemble of MaxEnt, GBM, GAM, MARS) performed well (TSS 0.71–0.88; AUC 0.90–0.96). Niche suitability change (NSC) under future climates decreased markedly more for warm-tolerant (warm-dry, warm-humid) groups than for cold-dry groups (Wilcoxon, FDR < 0.001), indicating greater future habitat decline for warm-tolerant populations. Cold-tolerant populations showed potential niche expansion. • Genome–niche index (integrating NSC and genomic offset) identified central Southwest Mountains populations as having minimal genome–niche interruption and as potential donor populations for evolutionary rescue (α ≈ 0.53 for T. elliotii; α ≈ 0.62 for P. monticolus under RCP8.5 2050). • Landscape connectivity analyses indicated viable dispersal routes through central Southwest Mountains enabling southward, westward, eastward, and northward movement. Evolutionary rescue is plausible via westward migration toward expanding niches in the southern Tibetan zone (T. elliotii) and eastern Himalayan zone (P. monticolus).

Integrating genomic offset with ecological niche modelling reveals that population-level adaptation mediates heterogeneous climate vulnerability across complex mountain landscapes. Although cold-tolerant populations exhibit higher genomic offset (indicating more genetic change needed to track future climates), their niche suitability is projected to be relatively stable or even expand, likely due to steep elevational gradients facilitating upslope/downslope movements. In contrast, warm-tolerant populations show smaller genomic offsets but face substantial declines in niche suitability, especially in patchy mid-elevational mountains with limited room for uphill shifts, heightening risk of habitat loss. The combined genome–niche index identifies populations least disrupted by climate change—primarily in the central Southwest Mountains—as potential donors for evolutionary rescue. Landscape genetic modelling supports connectivity corridors allowing multi-directional dispersal through these central regions, thereby enabling rescue to areas where suitable niches are projected to expand. These findings underscore that vulnerability assessments should consider intraspecific adaptive variation, ecological suitability dynamics, and landscape permeability simultaneously to prioritize conservation interventions and potential assisted gene flow or habitat connectivity enhancements.

This study demonstrates a framework that integrates genomic offset, population-specific ecological niche modelling, and landscape connectivity to assess climate change-driven vulnerability at intraspecific resolution. In two Sino-Himalayan montane birds, cold-tolerant populations have higher genomic offsets but are buffered by stable or improving niche suitability and connectivity, whereas warm-tolerant populations face pronounced suitability declines. A genome–niche index identifies central Southwest Mountains populations as likely donors for evolutionary rescue, and landscape analyses indicate feasible dispersal routes. The approach provides actionable insights for conservation planning under climate change. Future research should incorporate additional predictors of adaptive capacity and demographic resilience—such as phenotypic plasticity, genetic diversity, effective population size, and biotic interactions—and employ experimental validations (e.g., common garden, transplant experiments) and functional assays of climate-sensitive genes to refine forecasts and management strategies.

• Sample size constraints at certain localities (some with only 2–3 individuals) necessitated reliance on large SNP datasets and complementary modelling (GDM) to mitigate biases; gradientForest ideally benefits from ≥4 individuals per site. • ENM projections were limited to analogous climates to reduce extrapolation uncertainty; nevertheless, future novel climates can challenge transferability. • Genomic offset is a proxy for adaptive challenge and does not directly measure phenotypic plasticity, standing adaptive variation, or evolutionary rates. • Actual species responses depend on additional factors not explicitly modeled here, including phenotypic plasticity, evolutionary potential, effective population size, dispersal capacity, and interspecific interactions. • Arbitrary components in resistance surface parameterization for land cover were addressed via sensitivity tests, but residual uncertainty remains. • The equal contribution and potential unlisted affiliations do not affect analyses but may limit attribution clarity; nonetheless, analytical reproducibility is supported by publicly available data and code.