Coping with alpine habitats: genomic insights into the adaptation strategies of Triplostegia glandulifera (Caprifoliaceae)

Abstract How plants find a way to thrive in alpine habitats remains largely unknown. Here we present a chromosome-level genome assembly for an alpine medicinal herb, Triplostegia glandulifera (Caprifoliaceae), and 13 transcriptomes from other species of Dipsacales. We detected a whole-genome duplication event in T. glandulifera that occurred prior to the diversification of Dipsacales. Preferential gene retention after whole-genome duplication was found to contribute to increasing cold-related genes in T. glandulifera. A series of genes putatively associated with alpine adaptation (e.g. CBFs, ERF-VIIs, and RAD51C) exhibited higher expression levels in T. glandulifera than in its low-elevation relative, Lonicera japonica. Comparative genomic analysis among five pairs of high- vs low-elevation species, including a comparison of T. glandulifera and L. japonica, indicated that the gene families related to disease resistance experienced a significantly convergent contraction in alpine plants compared with their lowland relatives. The reduction in gene repertory size was largely concentrated in clades of genes for pathogen recognition (e.g. CNLs, prRLPs, and XII RLKs), while the clades for signal transduction and development remained nearly unchanged. This finding reflects an energy-saving strategy for survival in hostile alpine areas, where there is a tradeoff with less challenge from pathogens and limited resources for growth. We also identified candidate genes for alpine adaptation (e.g. RAD1, DMC1, and MSH3) that were under convergent positive selection or that exhibited a convergent acceleration in evolutionary rate in the investigated alpine plants. Overall, our study provides novel insights into the high-elevation adaptation strategies of this and other alpine plants.


Introduction
Plants have adapted to a wide array of habitats, including many extreme environments, including high elevations, deserts, and polar regions.Among these, alpine ecosystems are a matter of great concern because they are usually highly biodiverse and highly sensitive to climate change [1,2].High-elevation habitats are characterized by lower temperature, lower oxygen, and higher ultraviolet (UV) radiation compared with most low-elevation habitats.In response to these environmental challenges, alpine plants generally develop distinctive morphological traits, such as a dwarf stature, a cushion form, and a reduction in leaf size, which could result from convergent adaptive evolution [3,4].Recently, various omics technologies, such as genome sequencing, transcriptomics, and proteomics, have advanced research on the adaptation of alpine plants [5][6][7][8].The molecular mechanisms of alpine adaptation have been gradually revealed via the detection of positively selected genes (PSGs), fast-evolving genes (FEGs), expanded gene families, and differentially expressed genes (DEGs) [5][6][7][8].
Candidate genes discovered to date for high-elevation adaptation are involved in different functional pathways, especially DNA repair, UV-B tolerance, and cold tolerance [5,9].Approximately 100 f lowering plant families have alpine representatives [3], and recent research into their adaptations has only revealed the tip of the iceberg.Comparative genomics studies for highvs low-elevation species that explore the underlying genomic convergence of alpine adaptation among different families are very rare [5].
Due to global warming, the geographic range of some alpine plants is decreasing [10,11].Triplostegia glandulifera Wall.ex DC. is a typical alpine plant, of which the distribution is continuously shrinking, because it is highly susceptible to changes in habitat and climate [12].Triplostegia glandulifera is mainly distributed in the Himalaya-Hengduan Mountain region and Taiwan region, Central Sulawesi, and New Guinea, with the highest elevation recorded at 4000 m [13].Its root is a key component of traditional Chinese medicine used for liver protection and the treatment of kidney diseases.Triplostegia glandulifera belongs to Caprifoliaceae, the honeysuckle family (Dipsacales), which includes many coldadapted species from high elevations and high latitudes [14,15].The genome of the low-elevation species Lonicera japonica (Caprifoliaceae) has recently been published, making an omics comparison between high-and low-elevation Caprifoliaceae species possible [16].Therefore, T. glandulifera presents a promising candidate for investigating the adaptive evolution and environmental sensitivity of high-elevation plants.Understanding its adaptation strategies will be helpful for its cultivation and conservation.However, as of now, no relevant research has been published on T. glandulifera or any other Caprifoliaceae.
In this study we report a chromosome-level genome assembly of T. glandulifera and present 13 transcriptomes from other species of Dipsacales.Through comparative genomic analyses involving Dipsacales species or pairs of high-vs low-elevation species from different angiosperm families, we characterize genomic features, whole-genome duplications (WGDs), differential gene expression, gene family evolution, and adaptive signals in gene sequences.We find that preferential gene retention after WGD increased the copy number of cold-related genes in T. glandulifera.We also detect many genes with higher expression levels in T. glandulifera compared with its low-elevation relative, L. japonica, a series of which are putatively associated with alpine adaptation.We document the convergent contraction of gene families related to disease resistance in T. glandulifera and other investigated alpine plants when compared with their lowland relatives.We identify a set of FEGs and PSGs in the investigated alpine plants, many of which might be involved in alpine adaptation, such as DNA repair, cold response, and hypoxia response.Our findings provide novel insights into the adaptation strategies of alpine plants and also have implications for the conservation of these environmentally sensitive plants.

Genome sequencing, assembly, and annotation
The genome size of T. glandulifera (Fig. 1A) was estimated to be 719.94 Mb in the k-mer analysis (Supplementary Data Fig.S1) and 668.05 Mb in the f low cytometry test (Supplementary Data Fig.S2).We used a total of 124.17 Gb Illumina, 98.30 Gb PacBio, 95.89 Gb 10x Genomics, and 149.61Gb Hi-C sequencing reads for the de novo assembly of the T. glandulifera genome (Supplementary Data Table S1).After the primary assembly, correction, polishing, and scaffolding (see Materials and methods section), we obtained a final assembly of 680.38 Mb (contig N50 = 1.81 Mb, scaffold N50 = 66.68 Mb) within 1233 scaffolds (Supplementary Data Tables S2 and S3).The majority (∼91.25%) of the genome was anchored to nine pseudomolecules corresponding to the chromosomes (n = 9, Fig. 1B), with the sizes ranging from 36.35 to 110.58 Mb (Supplementary Data Fig.S3, Supplementary Data Table S4).The mapping rate of the assembled genome reached 98.14% based on the Illumina short reads (Supplementary Data Table S5).The genome and protein BUSCO [17] scores were 97.4 and 95.4%, respectively (Supplementary Data Table S6).Approximately 94.47% of the T. glandulifera transcriptome data can be mapped back to the assembled genome (Supplementary Data Table S7).
Repeat sequences comprised 62.88% of the T. glandulifera genome (Supplementary Data Table S8), with transposable elements (TEs) being the major component (Supplementary Data Table S9).Long terminal repeats (LTRs) were the most abundant type of TEs, covering 52.09% of the genome (Supplementary Data Table S9).Using a combination of de novo, homology-based, and RNA sequence-based predictions, we identified 32 123 proteincoding genes, with an average length of 3351 bp (Supplementary Data Table S10).Notably, 95.53% of the protein-coding genes were annotated in at least one of the six functional databases (Supplementary Data Table S11).Additionally, we identified 1427 miRNA, 755 tRNA, 6905 rRNA, and 552 snRNA genes (see Materials and methods section and Supplementary Data Table S12).

Phylogenomics and whole-genome duplication
We conducted the phylogenomic and molecular dating analyses based on 106 single-copy genes derived from the genome sequences of T. glandulifera, L. japonica, Daucus carota (Apiaceae), and the transcriptome sequences from 15 other species of Dipsacales covering Viburnaceae and all subfamilies in Caprifoliaceae, including 13 newly sequenced species and two species for which transcriptomes were previously published (see Materials and methods section, Supplementary Data Fig.S4 and Supplementary Data Tables S13 and S14).Triplostegia glandulifera was included in subfamily Dipsacoideae of Caprifoliaceae with strong support for the topology [T.glandulifera (Scabiosa tschiliensis, Dipsacus fullonum)] (Supplementary Data Fig.S5).
Our various analyses including synonymous substitutions per synonymous site (K s ) age distribution, phylogenetic approaches, and synteny assessments (see Materials and methods section) suggested at least five WGDs across the Dipsacales phylogeny (Fig. 2, Supplementary Data Figs S6-S8).One of these WGDs is shared by all investigated species of Dipsacales, another is shared at least by D. fullonum and S. tschiliensis, and the remaining three are specific to Valeriana officinalis, Zabelia bif lora, and Morina kokonorica, respectively (Fig. 2, Supplementary Data Figs S6-S8).The K s distribution for the T. glandulifera paralogs showed a signature peak for the WGD event (K s 0.964) following the γ-WGD (K s 2.11) shared by all core eudicots [20] (Fig. 2C).The K s value of paralogs in T. glandulifera was higher than that of orthologs shared by T. glandulifera and other species of Dipsacales but lower than that shared by T. glandulifera and D. carota (Fig. 2C).Combined with K s distributions for 16 other species of Dipsacales, this duplication in T. glandulifera seems to be shared by all Dipsacales species, herein named D-WGD (Fig. 2C, Supplementary Data Fig.S6).Additional signature peaks for recent WGDs were also observed in D. fullonum, S. tschiliensis, V. officinalis, Z. bif lora, and M. kokonorica (Supplementary Data Fig.S6).Through the Multi-tAxon Paleopolyploidy Search (MAPS) analysis, the Dipsacalesspecific D-WGD was further confirmed, while the recent WGDs in D. fullonum and S. tschiliensis were found to be shared (Fig. 2B, Supplementary Data Table S15).Intra-and intergenomic syntenic analyses revealed a 2:1 syntenic depth ratio between T. glandulifera and Vitis vinifera [21], as well as between L. japonica and V. vinifera, and a 1:1 syntenic relationship between T. glandulifera and L. japonica (Fig. 2D, Supplementary Data Figs S7-S8).These results provide the syntenic evidence for the D-WGD event shared by T. glandulifera and L. japonica.The age of the D-WGD was estimated to be ∼91.97 Ma, just before the diversification of Dipsacales during the Late Cretaceous (Fig. 2A, Supplementary Data Table S16).
In the T. glandulifera genome we identified 24 840 (77.33%) duplicated genes derived from five different modes of duplication (WGD; DSD, dispersed duplication; TRD, transposed duplication; TD, tandem duplication; PD, proximal duplication) using Dup-Gen_finder [22].WGD was the most prevalent type of gene duplication in T. glandulifera, as well as in L. japonica (Supplementary Data Fig.S9, Supplementary Data Table S22).The T. glandulifera WGD gene sets showed a significant enrichment in biological processes linked to signal response (such as response to cold, response to auxin, and response to gibberellin) and plant development (such as f loral organ development, leaf development, and shoot system development) in the functional enrichment analyses (Supplementary Data Figs S10 and S11, Supplementary Data Tables S23-S32).

Gene evolution and differential expression
Cold stress is a major abiotic stress that significantly impacts plant growth, survival, and geographical distribution, especially in alpine plants [3,23].Among various cold-tolerance mechanisms, the CBF (C-repeat binding transcription factor)-dependent cold signaling pathway is considered to be the most well-known; it is highly conserved among various plant species [23,24].In this pathway, CBFs are the key transcription factors [24].We found that WGD-derived genes in T. glandulifera enriched the number of genes that function in response to cold (Supplementary Data Fig.S10); we further surveyed the CBFs and 10 well-known positive regulators for cold tolerance in the CBF-dependent cold signaling pathway in Dipsacales and investigated the contribution of WGD in the evolution of these genes [23,24] (Fig. 3A, Supplementary Data Figs S12-S22).
Phylogenetic analysis of CBFs suggested that the D-WGD event resulted in two clades (DipsCBFa and DipsCBFb) in Dipsacales (Supplementary Data Fig.S12).In the DipsCBFa and DipsCBFb clades, CBFs further expanded by TD events, leading to five homologs in T. glandulifera and nine in L. japonica (Fig. 3A, Supplementary Data Fig.S12).Gene expression boxplots revealed that the average expression levels of CBFs in T. glandulifera were significantly higher than those in L. japonica (P = 0.0051) (Fig. 3C and Supplemental Data Fig.S12).Phylogenetic analyses also suggest that genes encoding RGAs (rice G-protein α subunits), SUMO E3 ligase SIZ1 (SAP and Miz1), ICEs (inducers of CBF expression), and CAMTAs (calmodulin-binding transcription activators) were duplicated into multiple copies through the D-WGD event in Dipsacales (Supplementary Data Figs S15 and S17-S19).Cold-responsive genes (CORs) are regulated under CBFs and directly cope with cold stress, with most of these genes encoding late embryogenesis abundant (LEA) proteins (Fig. 3A) [25].Triplostegia glandulifera had 77 LEA proteins, while L. japonica had 64 (Fig. 3A and B, Supplementary Data Table S33), and the average expression levels of LEAs in T. glandulifera were significantly higher than those in L. japonica (P = 0.0058) (Fig. 3D).
Furthermore, we conducted pairwise comparisons between T. glandulifera and L. japonica leaf transcriptomes (Supplementary Data Table S17) to identify DEGs (Fig. 3E and F, Supplementary   S18).F The significantly enriched GO terms for upregulated genes in T. glandulifera compared with L. japonica.The size of the circles shows the number of genes in one GO term.The color of the circles displays the statistical significance of enriched GO terms.Data Table S18).We identified a total of 3758 DEGs (2143 upregulated and 1615 downregulated genes) in T. glandulifera compared with L. japonica (Fig. 3E).Among the upregulated genes, the functional enrichment terms were significantly overrepresented in categories such as 'response to nitrogen compound', 'response to wounding', and 'response to water deprivation' (Fig. 3F, Supplementary Data Tables S34 and S35).More importantly, we discovered a series of upregulated genes in T. glandulifera with potential functions in alpine adaptation, including genes related to cold tolerance (e.g.CBFs, LEAs, CZF1) [23][24][25], hypoxia response (ERF-VIIs and ADH1s) [26], and DNA repair (e.g.RAD1, RAD51C, UNG) [27,28] (Supplementary Data Table S18).

Convergent contraction of gene families related to disease resistance
Gene family expansion and contraction analysis using CAFÉ [29] (see Materials and methods section and Supplementary Data Tables S19 and S20) indicated that the expanded gene families in T. glandulifera were significantly enriched in metabolic functions, while contracted ones were significantly enriched in plantpathogen interaction pathways (Supplementary Data Fig.S23, Supplementary Data Tables S36-S39).To investigate whether convergent changes in gene copy number occurred in distantly related alpine plants, we detected the convergent expansion or contraction of gene families in high-elevation plants compared with their low-elevation relatives from five different families: Caprifoliaceae (T.glandulifera vs L. japonica), Asteraceae (Erigeron breviscapus vs E. canadensis), Ericaceae (Rhododendron williamsianum vs R. ovatum), Brassicaceae (Crucihimalaya himalaica vs Capsella rubella), and Salicaceae (Salix brachista vs S. viminalis) (see Materials and methods section and Supplementary Data Table S21).All five pairs of plants underwent independent lineage-specific WGDs, without the effect of species-specific WGD (Fig. 2C, Supplementary Data Fig.S24).Only one gene family, belonging to the histidine phosphatase superfamily, displayed a striking expansion in all of these alpine plant genomes (Supplementary Data Table S40).In contrast, 28 gene families showed a significantly convergent contraction, including several genes associated with plant immunity (Fig. 4A, Supplementary Data Table S40).Examples include intracellular nucleotide-binding leucine-rich repeat receptors (NBS-LRRs, NLRs for short), the largest gene family for plant disease resistance [30], the L-type lectin receptorlike kinases (LecRKs) for plant innate immunity [31], and the cysteine-rich receptor-like kinases (CRKs) for disease resistance and programmed cell death [32].
Plants have evolved a two-tier immune system to recognize and activate defense against pathogen infections, including pattern-triggered immunity (PTI) and effector-triggered immunity (ETI) [33].The initiation of PTI and ETI needs cell-surface and intracellular immune receptors, respectively (Supplementary Data Fig.S25).Leucine-rich repeat-domain-containing receptorlike kinases (LRR-RLKs, RLKs for short) and receptor-like proteins (LRR-RLPs, RLPs for short) without the protein kinase domain are two major components of cell-surface immune receptors that recognize pathogen-associated molecular patterns (PAMPs) and activate PTI [34].NLRs serve as intracellular immune receptors that recognize pathogen-secreted effectors, activate ETI, and induce a hypersensitive cell death response [34].

Convergent evolution of fast-evolving or positively selected genes in alpine plants
Convergent adaptation is often associated with genes under positive selection (PSGs) [40].Similarly, genes with accelerated evolutionary rates (FEGs) are often associated with positive selection [41].To identify genes that might have evolved due to specific adaptations of high-elevation plants, we detected FEGs with a significantly higher K a /K s ratio in the high-elevation plants compared with the low-elevation ones and PSGs with a few sites under positive selection in high-elevation species (see Materials and methods section, Fig. 5, and Supplementary Data Table S55).
We identified 290 FEGs and six PSGs in the investigated alpine plants, and there was no overlap of FEGs and PSGs (Supplementary Data Table S55).Using the functional annotation information from the best BLAST hits of A. thaliana through UniPort and AceView, we found that 12 FEGs and one PSG were involved in DNA repair (e.g.DMC1, RAD1, ARP); four FEGs were involved in the response to cold (e.g.FIE, FBP7, ANNAT1); one FEG had a role in the response to light (e.g.HPR, MYB27, COP9); two FEGs were involved in the response to hypoxia (PCO3); and seven FEGs had a role in the response to water deprivation (e.g.ABI5, SDIR1, CDSP32) (Fig. 5, Supplementary Data Table S55).A series of genes associated with other stress responses as well as plant development were also found here, such as CFL1 for cuticle development [42], NAP1 for trichome morphogenesis [43], and VCL1 for vacuole formation [44] (Fig. 5, Supplementary Data Table S55).In addition, eight genes were found to be upregulated in the alpine plant T. glandulifera compared with its lowland relative L. japonica, such as RAD1 (TgChr02G08093.1),SDIR1 (TgChr05G20036.1),COP9 (TgChr01G01573.1),and CFL1 (TgChr08G26115.1)(Supplementary Data Table S18).

Discussion
In response to selective pressures posed by the alpine environment (e.g.low temperature, low oxygen, high UV-B radiation), high-elevation plants have evolved a series of adaptations [3,4].Here, we generated a chromosome-level genome assembly of the high-elevation species T. glandulifera (Caprifoliaceae) (Fig. 1B), and performed comparative genomic/transcriptomic analyses with its low-elevation relative L. japonica; we also included high-vs lowelevation comparisons of other phylogenetically distant alpine plants, to explore the potential genetic basis for the alpine adaptations.
Our genomic analysis revealed one relatively recent WGD event in the T. glandulifera genome, which was estimated to occur just before the diversification of Dipsacales during the Late Cretaceous (Fig. 2A).The preferential gene retention after WGD increased the copy number of cold-related genes in T. glandulifera, especially the positive regulators in the CBF-dependent cold signaling pathway (Supplementary Data Figs S10, S12, S15, S17, and S19).This finding highlights the potentially important role of WGD events for plants to colonize cold environments, such as high elevations or high latitudes [45][46][47][48].
By comparing differences in gene expression between T. glandulifera and L. japonica (Fig. 3E), we found that several highly Figure 5. Fast-evolving or positively selected genes with potential functions for alpine adaptation in T. glandulifera and other four alpine plants.The gene labeled with an asterisk is an PSG and others are FEGs.The genes in bold black exhibited higher expression levels in T. glandulifera than its lowland relative L. japonica.
expressed genes in T. glandulifera have potential functions for highelevation adaptation, such as CBFs and LEAs for cold tolerance, ERF-VIIs and ADH1s for hypoxia response, and RAD1 and RAD51C for DNA repair [23][24][25][26][27][28] (Supplementary Data Table S18).This implies that increasing the expression levels of some important genes is also an important adaptation strategy in alpine plants.
Through comparative genomic analysis among five pairs of high-vs low-elevation species, including T. glandulifera and L. japonica, we revealed that the gene families related to disease resistance have been contracted in a convergent fashion across five alpine plants compared with their lowland relatives (Fig. 4).Here, not only NLRs, as found in the previous study [5], but also two other key components of the plant immune receptors (RLPs and RLKs families) were observed to be reduced in the gene repertory size in alpine plants (Fig. 4).The observed streamlining of the alpine plant immune system is likely due to the relaxed selection pressure from pathogens, as microbial richness and diversity are lower in alpine regions than in the corresponding lowlands [49,50].Similarly, contraction of plant immune receptors has also been found in aquatic, parasitic, and carnivorous plants, where reduced pressure from pathogens is evident [51].The loss of genes related to immune response was also reported in some animals, and has similarly been suggested to be another adaptation to alpine living [52,53].
The loss of immunity genes makes sense energetically.For example, the immunity conferred by NLRs requires extensive energy, given that the activation of NLRs requires ATP binding [54].
In addition, disease resistance responses (PTI and/or ETI) could negatively regulate the expression of development-related genes [55].Therefore, unnecessary immune induction or a high copy number of genes used for disease resistance when there is limited challenge from pathogens would be selected against [56].In alpine environments, low oxygen dramatically reduces the efficiency of cellular ATP production [57], and low temperature also inhibits the rate of ATP synthesis [58].So, the co-contraction of the gene repertoire size for plant immune receptors in alpine plants ref lects an energy-saving strategy for survival in hostile alpine environments [59].In each family of plant immune receptors, the gene repertory size of the clades for pathogen recognition were found to exhibit a convergent reduction in alpine plants, such as sensor NLRs (CNLs) [34], prRLPs [38], and XII RLKs [39] (Fig. 4B and C), whereas helper NLRs (RNLs) for the transduction of immune signals [34], bRLPs and non-XII RLKs related to growth and development [38,39] showed nearly no change in gene repertory size in high-vs low-elevation species (Fig. 4B and C).As a whole, the genes for pathogen recognition showed relatively lower expression levels than the ones for survival and growth in the same gene family (Supplementary Data Fig.S31).This fits the idea that immune receptors would be perfectly 'off' in the absence of a trigger pathogen [60].
In addition to the molecular convergence in gene copy number, we also found a group of genes with a convergent change of gene evolutionary rate in the investigated high-elevation plants (Supplementary Data Table S55).The putative functions of these genes included DNA repair, response to cold, and response to hypoxia, etc. (Fig. 5, Supplementary Data Table S55), which are all related to the major environmental stresses inhibiting the survival and growth of alpine plants.Three genes, HPR, MYB27, and ARP, were also found to be PSGs in other studies of alpine adaptation [5,61].Strong UV radiation in alpine environments causes indirect damage to DNA.Here, 12 FEGs and 1 PSG for DNA repair were detected in the alpine plants, of which potential functions cover various modes of DNA repair (Fig. 5): DNA damage response (RAD1, PDCD5, REX1, and BRU1), base excision repair (ARP, FEN1, and LIG6), mismatch repair (MSH3 and MSH7), direct reversal of damage (AlkB), homologous recombination repair (DMC1 and RAD54), and nucleotide excision repair (RFC2) [27,28].These will be good candidates for future studies of gene function, to better understand the adaptation to DNA damage in alpine plants.
All of these findings suggest that alpine plants, including T. glandulifera, have developed diverse strategies to adapt to alpine habitats.Our study provides valuable genomic resources and important candidate genes for further in-depth exploration of high-elevation adaptation combining multi-omics data with an evolutionary developmental approach.Importantly, cold-adapted alpine plants with their streamlined immune systems are facing dual pressures from elevated temperatures as well as increased pathogen activity as a result of global climate change.We propose that conservation efforts for T. glandulifera and other alpine plants must prioritize the development of in situ and ex situ protection strategies that consider both biotic and abiotic factors.

Identification of fast-evolving genes and positively selected genes
The one-to-one orthologs among the above-mentioned five pairs of high-vs low-elevation plants were extracted using OrthoFinder (v.2.3.3)[86].Then each orthologous gene set was aligned using MAFFT (v.7.310) [87] and trimmed using GBlocks [116].To detect the common FEGs of high-elevation species, we used the branch model in CODEML from PAML [91].High-elevation species were set as the foreground branches and low-elevation species as background branches.The null model assumed that all branches evolved at the same rate, and the alternative model allowed the foreground branch to evolve under a different rate.To detect positive selection on a few codons along high-elevation species, we used the optimized branch-site model in CODEML from PAML [91], with one model allowing sites to be under positive selection on the foreground branch and another model assuming sites to evolve neutrally and under purifying selection.The likelihood ratio test was employed to determine significant differences between different models for each orthologs.The false discovery rate (FDR) in multiple testing was used to correct significance levels of P-values [117].The genes with an FDR-adjusted P-value <0.05 and a higher dN/dS in the foreground branch than in the background branch were considered as FEGs.The genes with an FDR-adjusted P-value <0.05 and an ω >1 of the foreground branch were considered as PSGs.

Figure 1 .
Figure 1.Triplostegia glandulifera plant and its genome.A Chromosomes (2n = 18) (right panel) and morphology of T. glandulifera with the f lower (left panel), leaf, stem and root displayed.B Landscape of the T. glandulifera genome.Tracks i-v represent nine pseudochromosomes (Mb), GC content, repeat element density, LTR density, and gene density, respectively.The inner curve lines indicate collinear gene blocks.

Figure 2 .
Figure 2. Time tree of Dipsacales and genome evolution of T. glandulifera.A Chronogram showing divergence times and WGD events in Dipsacales (including Viburnaceae and Caprifoliaceae).The asterisks show predicated WGDs in Dipsacales, and different colors correspond to different K s peaks in Supplementary Data Fig.S6.The blue dots indicate the placement of the calibrations for molecular dating.The blue curve on the top left of the tree shows the f luctuations of global average temperature over the past 100 Ma [19].The species photos were provided by Liu Bing, Kai-Lin Dong, Jian-Fei Ye, and Li-Xin Zhou.B Multi-tAxon Paleopolyploidy Search (MAPS) on the phylogeny.The red line indicates the percentage of subtrees that contained a gene duplication shared by descendant species at each node.The null simulations and positive simulations are displayed by black lines and gray lines, respectively.The red and purple asterisks represent the D-WGD in Dipsacales and another recent WGD shared by D. fullonum and S. tschiliensis, respectively.C K s distribution from paralogs of T. glandulifera (gray line) and orthologs of T. glandulifera vs each of the four species (colorful lines) (P.scabiosifolia, W. coraeensis, S. williamsii, and D. carota).D-WGD specific to Dipsacales and γ-WGD shared by all core eudicots are highlighted.D Synteny blocks among chromosomes of T. glandulifera, L. japonica, and V. vinifera.Synteny analysis suggested the 2:1 syntenic relationship between T. glandulifera and V. vinifera, and 1:1 between T. glandulifera and L. japonica (blue highlighted).

Figure 3 .
Figure 3. Gene evolution and differential expression.A Overview of the regulation of the CBF-dependent cold signaling pathway.The numbers of corresponding genes in T. glandulifera and L. japonica are shown in parentheses in the format of [T.glandulifera, L. japonica].The red asterisks indicate the retention of duplicated copies after the D-WGD event in T. glandulifera.B Phylogenetic tree of LEA genes from T. glandulifera, L. japonica, D. carota, and A. thaliana.The red asterisks indicate the well-known Arabidopsis CORs.C Gene expression boxplot for CBF genes in T. glandulifera and L. japonica ( * * P < 0.01).D Gene expression boxplot for LEA genes in T. glandulifera and L. japonica ( * * P < 0.01).E DEGs in T. glandulifera compared with L. japonica.Red dots represent significantly upregulated genes, blue dots represent significantly downregulated genes, and gray dots represent non-significant ones.Some upregulated genes with potential functions in alpine adaptation are shown (Supplementary Data TableS18).F The significantly enriched GO terms for upregulated genes in T. glandulifera compared with L. japonica.The size of the circles shows the number of genes in one GO term.The color of the circles displays the statistical significance of enriched GO terms.

Figure 4 .
Figure 4. Convergent contraction of gene families related to disease resistance in T. glandulifera and other alpine plants.A Scatter plot showing significantly convergently expanded (yellow-brown)/contracted (green) gene families in high-elevation plants compared with their low-elevation relatives.The circle sizes stand for log 2 (P-adjust), where P-adjust is the P-value of the binomial test adjusted for multiple testing.Boxplots for gene percentages of plant immune receptors between high-elevation plants (yellow-brown) and low-elevation ones (green): NBS-LRRs (NLRs), LRR-RLPs (RLPs), LRR-RLKs (RLKs) in the upper panel, and CNL, TNL, and RNL clades of NLRs in the lower panel (B); the pathogen-responsive RLPs (prRLPs) and basal RLPs (bRLPs) of RLPs in the upper panel, and the XII-and non-XII clades of RLKs in the lower panel (C).The percentage is the number of identified genes/number of searched genes.A paired sample t-test was used to analyze significant differences between the groups ( * P < 0.05).D Phylogenetic tree of NLR genes from T. glandulifera, L. japonica, four other pairs of investigated plants, A. thaliana, and species with known function NLRs.E Phylogenetic tree of RLP genes from T. glandulifera, L. japonica, four other pairs of investigated plants, A. thaliana, and species with known function RLPs.F Phylogenetic tree of RLK genes from T. glandulifera, L. japonica, four other pairs of investigated plants, and A. thaliana.E-F The yellow-brown lines label genes from high-elevation plants, the green lines label genes from low-elevation plants, and the black lines label known function genes and genes from A. thaliana.The well-known NLRs, RLPs, and RLKs are shown with names.