Next Article in Journal
Pharmacognostic Study, Diuretic Activity and Acute Oral Toxicity of the Leaves of Xiphidium caeruleum Aubl. Collected in Two Different Phenological Stages
Next Article in Special Issue
Changes in Pomological and Physical Parameters in Rosehips during Ripening
Previous Article in Journal
Detergent-Resistant Membranes in Chloroplasts and Mitochondria of the Halophyte Salicornia perennans under Salt Stress
Previous Article in Special Issue
Effects of High Voltage Electrical Discharge (HVED) on Endogenous Hormone and Polyphenol Profile in Wheat
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Plants
Volume 12
Issue 6
10.3390/plants12061266
plants-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Population Structure, Genetic Diversity and Candidate Genes for the Adaptation to Environmental Stress in Picea koraiensis

by
Ya Wang
1,
Zeping Jiang
2,
Aili Qin
2,
Fude Wang
3,
Ermei Chang
1,
Yifu Liu
2,
Wen Nie
2,
Cancan Tan
2,
Yanchao Yuan
2,
Yao Dong
2,
Ruizhi Huang
2,
Zirui Jia
1,* and
Junhui Wang
1,*
1
State Key Laboratory of Tree Genetics and Breeding, Research Institute of Forestry, Chinese Academy of Forestry, Beijing 100091, China
2
Key Laboratory of Forest Ecology and Environment of National Forestry and Grassland Administration, Ecology and Nature Conservation Institute, Chinese Academy of Forestry, Beijing 100091, China
3
Forestry Research Institute in Heilongjiang Province, Harbin 150081, China
*
Authors to whom correspondence should be addressed.
Plants 2023, 12(6), 1266; https://doi.org/10.3390/plants12061266
Submission received: 20 January 2023 / Revised: 7 March 2023 / Accepted: 8 March 2023 / Published: 10 March 2023
(This article belongs to the Special Issue Plant Science Research in Nanjing Agricultural University: Commemorating the University’s 120th Anniversary)
Browse Figures
Review Reports Versions Notes

Abstract

:
Picea koraiensis is major silvicultural and timber species in northeast China, and its distribution area is an important transition zone for genus spruce migration. The degree of intraspecific differentiation of P. koraiensis is high, but population structure and differentiation mechanisms are not clear. In this study, 523,761 single nucleotide polymorphisms (SNPs) were identified in 113 individuals from 9 populations of P. koraiensis by genotyping-by-sequencing (GBS). Population genomic analysis showed that P. koraiensis was divided into three geoclimatic regions: Great Khingan Mountains climatic region, Lesser Khingan Mountains climatic region, and Changbai Mountain climatic region. Mengkeshan (MKS) population on the northern edge of the distribution area and Wuyiling (WYL) population located in the mining area are two highly differentiated groups. Selective sweep analysis showed that MKS and WYL populations had 645 and 1126 selected genes, respectively. Genes selected in the MKS population were associated with flowering and photomorphogenesis, cellular response to water deficit, and glycerophospholipid metabolism; genes selected in the WYL population were associated with metal ion transport, biosynthesis of macromolecules, and DNA repair. Climatic factors and heavy metal stress drives divergence in MKS and WYL populations, respectively. Our findings provide insights into adaptive divergence mechanisms in Picea and will contribute to molecular breeding studies.

1. Introduction

Picea is the third-largest genus in Pinaceae and a major component of coniferous forests in the northern hemisphere. Most conifers have the characteristics of long generations, large population size [1,2,3], and a low level of interspecific reproductive isolation [4]. Spruce underwent rapid differentiation during the Pliocene, Qinghai–Tibet Plateau uplift, and Quaternary glacial oscillation [5]. The continuous gene flow among sympatric and parapatric species makes the evolutionary relationship more complex [6,7]. Geographic isolation and gene introgression are two possible driving factors of species differentiation [8,9]. In conifers with long generations, a long period of time can elapse from the onset of differentiation to the cessation of gene flow. Studies have suggested that gene flow between different species continues to decline over time when the species is under either isolation or selective pressure [10]. Understanding gene introgression and selection pressure is of great significance for further understanding spruce phylogeny.
Picea koraiensis, belonging to the genus Picea, is distributed in the northeast of the Great Khingan Mountains, Lesser Khingan Mountains, Zhang Guangcai Mountains, Wanda Mountains, and Changbai Mountains in China, with elevations ranging from 300 m to 1600 m. As an independent species, P. koraiensis has generally been accepted by the academic community, but its population genetic structure and geographical migration and spread history are not clear. There are different hypotheses about the origin of the genus spruce, the “North American Origin Theory” [11,12,13] and the “Asian Origin Theory” [14,15,16], which are supported by different DNA evidence. Whether spruce originated in North America or southeastern Asia, its migration route is through the northeast plain of China, from North America to the Qinghai–Tibet Plateau, or from east Asia or the Qinghai–Tibet Plateau to North America. Therefore, the distribution area of P. koraiensis, the Chinese northeast plain, is an important transition zone for spruce migration. Many fossils of the genus Picea from the Jurassic to Early Cretaceous, which is the earliest fossil of the genus Picea found thus far, have been found in the distribution area of P. koraiensis [17,18,19]. Some studies have been performed on the natural population division and genetic diversity of P. koraiensis. For example, Zhang et al. [20] divided the natural population of P. koraiensis into five climatic zones. The genetic diversity of the P. koraiensis population was analyzed by phenotypic biosystematics, seedling morphology, growth, biomass, allozyme variation [20], mineral element content [21], and random amplified polymorphic DNA (RAPD) molecular markers [22]. Variations within and among populations of P. koraiensis have been found. The variation within populations was greater than the variation among populations, and the genetic distance of P. koraiensis was positively correlated with geographical location to a certain extent. The level of genetic differentiation among different geographical populations of P. koraiensis is very high, but the population genetic structure, driving factors, and molecular mechanism of differentiation are not clear, and the phylogeography of P. koraiensis at the whole genome level has not been reported.
In this study, P. koraiensis was used as the research object to study the genetic structure of the natural population and the adaptive differentiation mechanism of different geographical populations to environmental stress by genotyping-by-sequencing. Based on the above research, we put forward the following questions: (i) What is the genetic relationship, genetic diversity, and differentiation of P. koraiensis in the natural distribution area? (ii) Is there any obvious differentiation among populations in different geographical distribution areas of P. koraiensis? (iii) If there is obvious differentiation, what is the molecular mechanism of differentiation motivation and adaptability?

2. Results

2.1. Characterization of GBS-Seq Data and SNPs

The genotyping-by-sequencing (GBS) of 127 samples produced the original reads of 563.64 Gb. After quality filtering, the clean reads of 532.685 Gb were retained and the average value of each individual was 4.19 Gb. The average efficiency was 94.5%. The sequencing results showed that the average Q20 value, Q30 value, and GC content were 97.38%, 91.63%, and 40.21%, respectively. After comparing the filtered data with P. abies, the average comparison rate was 95.33%. The number of SNP (single nucleotide polymorphism) bits detected was 10,375,590. Phylogenetic analysis was conducted with 127 samples, and a total of 523,761 SNPs with missing ratio ≤20% minor allele frequency (maf) ≥0.01 were identified.

2.2. Phylogenetic and Population Genetic Structure Analysis

According to the phylogenetic tree constructed by the neighbor-joining (NJ) method (Figure 1), with P. jezoensis, P. likiangensis, and P. pungens as outgroups, 112 individuals of P. koraiensis were divided into 3 groups: Great Khingan Mountains group (Group I), Changbai Mountain group (Group II), and Lesser Khingan Mountains group (Group III). Among these groups, Group I is composed of the Mengkeshan (MKS) population in the northern Great Khingan Mountains, which is located at the base. Group II is composed of four populations: Hailin (HL), Muling (ML), Linjiang (LJ), and Tianqiaoling (TQL) populations, and Group III is composed of the Wuyiling (WYL), Aihui (AH), Gaofeng (GF) and Zhanhe (ZH) populations. Group II and Group III are sister clades. The results of the phylogenetic tree showed that the genetic relationship of different populations of P. koraiensis was related to geographical distribution. The population genetic structure was further analyzed by ADMIXTURE v1.23 software [23]. The number of clusters (K) was set from 2 to 10, and the best number of clusters was selected according to the minimum cross-validation rate (Figure 2). The optimal number of clusters was three. When K = 2 to 4, the error rate of cross-validation was small, indicating that 2 to 4 clusters is more reasonable. When K = 2, 127 samples were divided into P. koraiensis and the outgroup, and the MKS population of P. koraiensis had a small amount of outgroup gene infiltration. When K = 3, P. koraiensis was divided into 2 groups, the WYL population and other populations. When K = 4, other populations had extensive gene flow except for the WYL and GF populations. The results showed that gene introgression occurred in the MKS population and outgroups, and the WYL population had greater differentiation from other populations. To further analyze the population genetic structure, 112 individuals of P. koraiensis were analyzed by principal component analysis (PCA) (Figure 3). The results showed that P. koraiensis could be divided into three main clusters: MKS, WYL, and other populations. The results of PCA and admixture analyses indicated that there is an obvious differentiation between the two populations (MKS and WYL) and others.

2.3. Population Diversity and Selective Sweep Analysis

The data are shown in Table 1; the π(nucleotide diversity) of the P. koraiensis population was 0.00194 to 0.00231. Among the populations, the π value of the WYL population was the lowest (0.001943), while the π value of LJ was the highest (0.00235). The Fst (fixation index) among P. koraiensis was 0.0254 to 0.1569, among which the mean Fst between the WYL population and other populations was the highest (0.1211), followed by the MKS population (0.0983), and the mean Fst between the TQL population and other populations was the lowest (0.0538). This result is consistent with the phylogenetic analysis, indicating that the WYL and MKS populations have a high degree of differentiation from other populations, and these two populations may be isolated from other populations and have less gene flow, resulting in a significant decrease in nucleotide diversity. That is, these two populations experienced a high degree of selection in the process of evolution. The degree of differentiation of the LJ, TQL, and ZH populations was low, and the three populations had frequent gene flow with other populations and high nucleotide diversity. The Tajima’s D values of the nine populations of P. koraiensis were all positive (mean 0.1554), indicating that these nine populations may have experienced balanced selection.
The main purpose of selective sweep analysis was to find the selected genes of the population, which further reveals the adaptive mechanism of population evolution. In our study, four populations of Changbai Mountain group were used as background group. The MKS and WYL populations were used as test groups, respectively. Fixation index (Fst) and nucleotide diversity ratios (π ratio) were used in combination to identify regions where selective sweeps occurred during species differentiation. A total of 645 and 1126 selected genes were identified in the MKS and WYL populations, respectively (Tables S1 and S2).
The results showed that the selected regions of the MKS population were Fst (MKS) ≥ 0.19 and π ratio ≥ 1.57 (Figure 4a). Based on GO (gene ontology) enrichment analysis, the molecular functions of the genes were concentrated mainly in the response to radiation, response to light stimulus, macromolecule methylation, response to reactive oxygen species, nitrogen compound metabolic process, ubiquitin-dependent protein catabolic process, and nucleic acid metabolic process (Figure 5a; Table S3). As shown in Figure 6a, the KEGG (Kyoto Encyclopedia of Genes and Genomes) pathways were concentrated mainly in ABC transporters, RNA transport, inositol phosphate metabolism, glycerophospholipid metabolism, RNA polymerase, nucleotide excision repair, mismatch repair, flavonoid biosynthesis, cysteine and methionine metabolism, peroxisome, endocytosis, and biosynthesis of amino acids. Among these genes, we found genes related to photoperiod, flowering, and photomorphogenesis (e.g., CUL4, GATA9, SBH2, and MBD9), and genes involved in the cell response to water scarcity (such as RPK2, CYT1, and MGL), which may be related to the climate adaptation of the MKS population to the northern Great Khingan Mountains. In the selected region, we also found genes related to metal ion transport (e.g., ABCB family genes), flavonoid biosynthesis-related genes (e.g., ANS, DFR), cysteine metabolism-related genes (e.g., DNMT1, DNMT3A), amino acid biosynthesis genes (e.g., GlyA, thrB, cysK, RPE, and TYRAAT), peroxidase metabolism-related genes (e.g., AGXT, SOD1, and PEX1), genes related to response to reactive oxygen species (e.g., HSP, SODCC, CPN60-2, and CYT1), signal transduction-related genes (e.g., CTR1), and oxidative stress defense-related genes (e.g., MYB family genes), showing that the MKS population may have been subjected to stress. In addition, genes related to DNA repair (e.g., CUL4, MNAT, RFA1, MLH3, and POLD2) were also found, suggesting that the MKS population may have DNA damage under adverse conditions.
In the WYL population, the selected regions were Fst (WYL) ≥ 0.32 and π ratio ≥ 2.5 (Figure 4b). Based on the GO enrichment analysis, the molecular functions of the selected genes were concentrated mainly in chromosome condensation, regulation of macromolecule metabolic process, gene expression, cation transport, RNA metabolic process, metal ion transport, ATP hydrolysis coupled proton transport, etc. (Figure 5b, Table S4). As shown in Figure 6b, the KEGG pathways were concentrated mainly in phosphonate and phosphinate metabolism, phagosome, glycerophospholipid metabolism, oxidative phosphorylation, ribosome, RNA polymerase, peroxisome, glutathione metabolism, proteasome, etc. Among these genes, we found genes related to plant–pathogen interactions (e.g., CNGC, CPK, FLS2, NOA1, CERK1, CML, RPS2, MPK6, and KCS), response to salt-stress-related genes (e.g., CYT1, FAR4, MYB78, and PPR40), bivalent metal ion transport-related genes (e.g., TL17, VHA_D, SAUR32, AGT1, TIP2-1, and GLR 3.3), univalent metal ion transport-related genes (e.g., HAK family genes, YMF19, VHA, AVP1, and OEP64), plant hormone signal-transduction-related genes (e.g., IAA, SAUR, ARR-B, GID1, MPK6, and EBF1_2), flavonoid-biosynthesis-related genes (e.g., ANS, FLS, ANR, and LAR), and DNA-repair-related genes (e.g., TOP3, RFA1, BRCA2, MLH1, PMS2, and RBX1), showing that the WYL population may differentiate under the stress of environmental factors, and the metal ions in the soil of its distribution area may be higher than the metal ions in the soil of other population areas.

3. Discussion

3.1. Geographical Isolation, Climate Heterogeneity, and Gene Introgression

In our study, whole-genome SNPs were obtained by the GBS-seq technique, and the phylogeny, population structure, and genetic differentiation of P. koraiensis distributed in northeast China were analyzed. The results showed that the genetic relationship among different populations of P. koraiensis was consistent with the division of geographical and climatic regions: there was gene introgression between the MKS population and outgroup, and there was obvious genetic differentiation among the MKS and WYL populations and other populations. Climate heterogeneity, gene introgression, and geographical isolation are the main reasons for the differentiation of P. koraiensis.
Climatic factors such as temperature and precipitation are likely to be among the important factors driving the differentiation of P. koraiensis due to different habitats. The phylogenetic results divided P. koraiensis into three groups corresponding to three geoclimatic regions: the Great Khingan Mountains climatic region, the Lesser Khingan Mountains climatic region, and the Changbai Mountain climatic region. The climate of the northern Great Khingan Mountains area is characterized by long sunshine, low precipitation, and low temperature; the climate of the Changbai Mountain region is characterized by medium sunshine, medium precipitation, and medium temperature; and the climate of the Lesser Khingan Mountains area is characterized by short sunshine, low precipitation, and low temperature [20]. Zhang et al. [20] divided the distribution of P. koraiensis in China into five sections: the Changbai Mountains area, the Laoye Mountain–Wanda Mountains area, the Lesser Khingan Mountains area, the north part of Great Khingan Mountains, and the southwest part of Great Khingan Mountains. Among these sections, the spruce species distributed in the southwestern part of the Great Khingan Mountains is P. meyeri var. mongolica. The Laoye Mountain–Wanda Mountains area belongs to the branches of Changbai Mountain, with little difference in climate factors. The results of the phylogenetic and genetic diversity analyses showed that the populations of the Laoye Mountain–Wanda Mountains and Changbai Mountain areas were clustered into monophyletic groups. Therefore, we classified the two populations in the Laoye Mountain–Wanda Mountains (HL and ML populations) in the Changbai Mountain climatic zone. Among other species, there are also examples of temperature and moisture influencing population differentiation, such as Liaodong oak [24], eucalyptus in Australia [25], and alpine conifers in the western Alps [26].
Gene flow and introgression are also important causes of intraspecific differentiation in P. koraiensis. On the one hand, the size of gene flow between different groups within the species and the degree of differentiation were inversely proportional. These two groups, WYL and MKS, have little gene exchange with other groups, low nucleotide diversity, and greater differentiation; the ZH, TQL, and LJ populations have extensive gene flow with other groups, high nucleotide diversity, and low differentiation. On the other hand, gene introgression between species promotes the differentiation of different populations within P. koraiensis. The results of phylogenetic and population differentiation analysis showed that the MKS population was more differentiated from other populations and had a small genetic component of outgroups. Because of the large overlap in distribution between the MKS population of P. koraiensis and its sympatric species P. jezoensis, we speculate that there is a small amount of gene introgression from P. jezoensis to the MKS population. Many closely related species of P. koraiensis have been successfully crossed with closely related species of P. jezoensis in hybridization experiments [15]. Fossil evidence suggests that during the Pleistocene, P. jezoensis and its extinct relatives were widespread in northeastern Russia [27], that the ancestors of P. koraiensis and its relatives may have appeared in east Asia during the same period [28], and that the ancestors of P. jezoensis and P. koraiensis underwent infiltration before P. koraiensis and its relatives’ migration to Eurasia [29]. Gene introgression between other species of the genus spruce has been reported repeatedly: P. glauca, P. engelmannii, and P. sitchensis in western North America [30,31,32,33,34,35] and P. likiangensis, P. purpurea, and P. wilsonii in the Qinghai–Tibet Plateau [6,7,36,37,38]
Geographical isolation is an important reason for the intraspecific differentiation of P. koraiensis. At the northern edge of the P. koraiensis distribution area, the MKS population has the characteristics of geographical and genetic isolation because of spatial isolation and limited gene communication; that is, the differentiation characteristics of the MKS population accord with the central–marginal hypothesis, in which marginal populations show lower genetic variation and higher differentiation than central populations. There are similar examples in other species, such as Abies sachalinensis in Hokkaido and Pinus cembra L. in the Alps [39,40]. The WYL population is distributed in the Lesser Khingan Mountains gold mining area. Due to the influence of mining and different geological conditions, the natural environment shows fragmentation and islanding, forming an ecologically isolated island, which is different from other distribution areas of P. koraiensis. Furthermore, there is geographical isolation. There are similar examples of other species distributed in mining areas, such as Caryophyllaceae [41,42], Scopelophila cataractae [43], and Colocasia esculenta [44]. In summary, geographical isolation leads to the differentiation of species and the formation of obvious geographical differences to a great extent; the same is true of Lamiaceae and Circaeasteraceae [45,46].

3.2. Molecular Mechanisms of Population Differentiation

To further detect natural selection signals in different populations of P. koraiensis and explore the differentiation mechanism, we performed selection sweep analysis to identify candidate genes based on KEGG and GO enrichment analyses involving various biological functions, using populations of the Changbai Mountain group as background populations and the more differentiated MKS and WYL populations as control populations. The key functional genes for which differentiation between different populations of P. koraiensis was subject to selection were to improve the resistance of the population, to promote flowering and fruiting ability for adaptation to different habitats, and to continuously reproduce offspring.
The MKS population is distributed in the northern part of the Great Khingan Mountains and is geographically distant from other populations, with large climatic differences. Climatic heterogeneity drives significant differentiation of the MKS Mountain population, and the function of genes shows the molecular adaptation mechanism of the MKS population to its habitat. The northern Great Khingan Mountains have long sunshine hours (approximately 2600 h per year and 2300 h per year in the Changbai Mountain region), and light affects the photoperiod, circadian rhythm, and photomorphic building of plants. The results of our selection sweep analysis showed that 15 genes in the MKS population were enriched in response to radiation, and 14 genes were enriched in response to light stimulus. Among these genes, CUL4 in Arabidopsis can regulate the flowering time in a short period of time during the day and is involved in photomorphogenesis [47,48]; GATA9 is involved in the regulation of circadian rhythms in Arabidopsis [49]; SBH2 is involved in photomorphogenesis [50]; and MBD9 regulates photoperiod and flowering in Arabidopsis through DNA methylation [51,52]. The northern Great Khingan Mountains are characterized by low annual precipitation (463 mm average annual rainfall and 680 mm average annual rainfall in Changbai Mountain), and our results show that 5 of the genes are enriched in response to water deprivation, with RPK2 being involved in pollen germination, response to cold, and response to water deprivation [53]; CYT1 being involved in defense responses to bacteria, hormone homeostasis, and salt stress [54,55]; and MGL being involved in cellular response to sulfate starvation and cellular response to water deprivation [56]. The northern Great Khingan Mountains are characterized by low temperatures (−2.4 °C annual mean temperature and 3.2 °C annual mean temperature in Changbai Mountain), and organisms living in this environment face various growth-related challenges from low temperatures, such as reduced lipid membrane fluidity. The lipid bilayer transmits external signals to its interior and protects the integrity of various organelle membranes in unfavorable environments. The results showed that 86 genes were associated with membrane-bound organelles and that 6 genes were involved in glycerophospholipid metabolism. Among these genes, PLC is involved in leaf membrane lipid remodeling [57] and regulation of stomatal movement and root development under abiotic stress [58,59]. Remodeling of membrane lipids may protect membranes from low temperatures. We hypothesize that these genes may play a major role in maintaining membrane integrity under low-temperature conditions. The functions of these genes suggest that MKS populations have evolved a range of genes adapted to the long light hours, low annual rainfall, and low temperatures of the habitat.
Among the Lesser Khingan Mountains group, the WYL population is distributed in the gold mining area, and there are also large amounts of iron ore and agate in its distribution area. The WYL population is distributed at a low altitude (average elevation approximately 350 m), and due to tailings pollution caused by mining, monsoon effects, and rainfall leaching, many heavy metals, such as Cd, Cr, Cu, Ni, Pb, and Zn, in the soil exceed background values, especially Pb, Cd, Zn, and Cu [60,61]. Heavy metals cause oxidative damage, disruption of cellular homeostasis, DNA strand breaks, fragmentation of proteins or cell membranes, and damage to light and pigments [62]. Plants undergo a series of physiological changes in response to heavy metal stress to adapt to their environment. The first line of defense of plants against heavy metal stress is to block the transport of metal ions in the roots through the plasma membrane [63]. The enrichment analysis showed that 19 genes are involved in the transport of metal ions, of which 5 genes were enriched in ABC transporter, where the ABC transporter family (e.g., ABCB1, ABCA15, ABC1, ABCG28, and ABC1K8) is involved in the efflux of metal ions from the plasma membrane and in the resistance of plants to heavy metals and aluminum toxicity [64,65]. For example, the ABCB family genes enhance resistance to cadmium and lead in Arabidopsis [65], and ABCB1 is upregulated in mango when it is subjected to lead stress [66]. After heavy metals penetrate the first line of defense and enter the cell, plants tolerate or neutralize the toxicity produced by heavy metals through the synthesis of biomolecules such as metal chelators and specific amino acids [63]. Forty-eight genes were identified as regulators of cellular macromolecular biosynthetic processes, among which PCYT2, PCYT1, EPT1, and PHO1 are involved in the synthesis of phosphate and hypophosphate [67,68,69,70,71,72,73]; ANS, FLS, and ANR are involved in flavonoid biosynthesis [74,75,76]; CysCK is involved in the synthesis of cysteine [77]; and G6PD, SPE3, GST, and RRM1 are involved in glutathione synthesis [78,79,80,81]. Phosphates, flavonoids, cysteines, glutathione, and other phytochelators are markers of heavy metal stress in plants [82,83]. GST genes are upregulated in horse rushes, tomatoes, and radishes when they are subjected to heavy metal stresses such as cadmium [84,85,86]. When both of the above mechanisms against heavy metal stress are exhausted, plants activate oxidative stress defense mechanisms and the synthesis of stress-related proteins and signaling molecules [63]. Seven genes were associated with oxidative stress defense in plants, among which MPK6 and FLS2 were enriched in the MAPK pathway and IAA, SAUR, and GID1 were enriched in the hormone signaling pathway. In Arabidopsis, MPK6 can be induced by Cd ions and Cu ions [87,88,89,90], and the IAA (AUX) concentration increases to promote lateral root growth and protect plants from Cd [91]. In addition, the regulation of transcription factor family expression by heavy metal stress has been reported [92,93,94]. The transcription factors identified in this study are members of the MYB family (MYB4, MYB52, MYB58, MYB78), WRKY family (WRKY7, WRKY35, and WRKY42), NAC, and BZIP68. In Arabidopsis, MYB family genes such as MYB4 are highly induced by Cd and Zn metal stress [95]. The expression of WRKY family members was significantly elevated under Cu and Cd metal exposure [96]. BZIP transcription factor expression was also induced by Cd stress [97]. For DNA damage caused by heavy metals, plants initiate DNA damage repair mechanisms with six genes involved in DNA repair, of which RBX1 and RFA1 are involved in nucleotide excision repair [98,99], MLH1 and PMS2 are involved in mismatch repair [69,100,101], and TOP3 and BRCA2 are involved in homologous recombination [102,103], suggesting that the WYL population evolved multiple DNA repair pathways. In addition, compared with Changbai Mountain, the sunshine in the Lesser Khingan Mountains area is short and the temperature is low. The WYL population also identified genes related to light morphogenesis, circadian rhythm, and cold tolerance. For example, in Arabidopsis thaliana, the CSN4 and CSN8 participate in the development of seedlings in the dark and promote the development of light morphology in the dark [104,105]; ABC1K8 and BLH1 are involved in regulating circadian rhythms [106,107]; and LSM1A, CIPK9, and PSRP2 are involved in plant adaptation to cold environment [108,109,110]. Our results show that the genome of the WYL population not only responds to the climatic environment of the Lesser Khingan Mountains area, but also may be stressed by heavy metals in mining areas, and the genome has undergone adaptive evolution.

4. Materials and Methods

4.1. Sample Collection and GBS Analysis

We collected samples according to the 4 climatic zones of P. koraiensis and identified a total of 9 populations of 10 to 15 individuals each, for a total of 112 individuals (Figure 7, Table 2). A total of 2 populations of the Changbai Mountain climate zone (Region I): Tianqiaoling (TQL) and Linjiang (LJ); 2 populations of the Laoye Mountain, Zhangguangcai Mountain, and Wanda Mountains climate zone (Region II): Muling (ML) and Hailin (HL); 4 populations of the Lesser Khingan Mountains climate zone (Region III): Wuyiling (WYL), Zhenhe (ZH), Aihui (AH), and Gaofeng (GF); and 1 population of the Great Khingan Mountains climate zone (Region IV): Mengkeshan (MKS). Three species, P. jezoensis (sympatric species of P. koraiensis), P. likiangensis, and P. pungens, were used as outgroups, with five individuals per species. Plants from each population were separated by at least 100 m at the time of sampling, and needles were dried directly using silica gel and kept at low temperature. DNA extraction from plant needles was performed using a modified bromocetyltrimethylamine (CTAB) method, followed by a Nanodrop-1000 spectrophotometer (Nanodrop, MA, USA) to detect nucleic acid concentration and DNA purity on 1% agarose gel electrophoresis.
DNA samples were qualified, sequenced, and constructed in GBS libraries according to the methods reported by Poland et al. [111]. Each sample extracted from 1.5 μg of DNA was digested in 96-well plates with EcoRI and NiaIII restriction enzymes. The digestion reaction was performed in 96-well plates using 100 ng of DNA for each individual reaction digested with EcoRI and NIaIII (New England Biolabs, Ipswich, MA, USA). The digested product was mixed with A1 and A2 linkers at 25 pmoL per well, adding linkers at both ends of the DNA. Library pooling, selection of size (400–600 bp) on a 1% agarose gel, column purification using PCR purification kit (NEB), and amplification of 12 cycles using Phusion DNA polymerase (NEB) were carried out. After the second column wash, the average fragment size was estimated on a BioAnalyzer 2100 (Agilent, Santa Clara, CA, USA) using a DNA1000 chip, and library quantification was performed using PicoGreen (Invitrogen, Carlsbad, CA, USA). Sequencing with PE125 on HiSeq4000 (Illumina, San Diego, CA, USA) adjusted the common library to 10 nmol. Numerous raw sequencing sequences (raw reads) were obtained through high-throughput sequencing. To ensure data quality, the quality control of the original data should be carried out before information analysis, and data filtering should be used to reduce data interference information. We filtered the dismounted raw reads more strictly to obtain high-quality clean reads for subsequent information analysis. The conditions for filtering were as follows: (1) reads containing adapter sequences were removed; (2) reads containing unknown base N ratios greater than 10 were removed; and (3) low-quality reads (the number of bases of mass value Q ≤ 10 accounted for more than 50% of the entire read) were removed for genomic electronic digestion assessment. To check the efficiency of digestion and analyze all markers, we digested the reference genome according to the site of restriction enzymes and counted the corresponding conditions for later evaluation. Sequencing evaluation: statistics are made of output data, including sequencing data yield, sequencing error rate, Q20 content, Q30 content, N content, GC content, etc., to evaluate whether the library construction and sequencing are successful.
Using the comparison software Burrows-Wheeler Aligner (BWA) ver.0.7.12 (Li, Cambridge, UK) [112], the filtered reads were compared with the reference genome (European spruce) by the mem algorithm, and the alignment parameter was-k 32-M. After the comparison, the results were marked by Picard (1.129) software. The Unified Genotyper module of the software GATK (3.446) is used to detect the variant of multiple samples of the processed comparison files, and the detected variations are filtered by VariantFiltration. The filtering parameters are -Window 4, -filter “QD < 4.0 || FS > 60.0 || MQ < 40.0”, -G_filter “GQ < 20”. ANNOVAR was used to make functional comments on the detected variant. Using VCFTOOL v.0.1.11 (Danecek, Cambridge, UK) [113], the SNP loci (127 individuals) obtained by GBS sequencing were filtered if the missing ratio ≤20% dint maf ≥0.01 and all indels were removed.

4.2. Phylogenetic and Population Structure Analyses

The phylogenetic tree of 127 individuals was constructed by the Neighbor-Joining method and the P-diatance model [114] based on MEGA version X software [115] and corrected by Jukes-Cantor, Tajima-Nei, and the Maximum Composite Likelihood model. A bootstrap consensus tree was formed from a sample of 1000 replicates. iTOL version 6.4 was used to edit and beautify the phylogenetic tree [116].
To carry out the cluster analysis of the population structure, the parameters of the ancestry coefficient matrix Q and the population allele frequency matrix F and K were calculated by admixture software [23] and finally plotted by the R software package (http://www.r-project.org/, accessed on 12 January 2023). The K value range was set from 1 to 100, and the default settings were used for the other parameters. The optimal clustering number was determined according to the minimum value of the cross-validation (CV) error, and the R language package was used to draw the results.
The genetic structure of the spruce population was evaluated by principal component analysis (PCA). We used GCTA software (http://cnsgenomics.com/software/gcta/pca.html, accessed on 12 January 2023) to reduce the dimension of the selected SNP loci and calculate the eigenvectors and eigenvalues. Three principal components were obtained by using GCTA, and the PCA distribution map was drawn by using the R language package. The values in the axis label parentheses represent the percentage of the principal component that explains the overall variance.

4.3. Genetic Diversity and Selective Sweep Analysis

Nucleotide diversity (π) (also known as θ π, Pi), genetic diversity analysis (Tajima’s D), and genetic differentiation among populations (Fst) were analyzed by the PopGenome module in the R software package. According to the phylogenetic results, P. koraiensis was divided into three groups: the Changbai Mountain group, the Great Khingan Mountains group, and the Lesser Khingan Mountains group. Because of the moderate climate and frequent gene communication with other populations in the Changbai Mountain group, we took four populations (HL, ML, LJ, and TQL) in this group as background groups. Combined with the results of Fst, the MKS population of the Great Khingan Mountains group and WYL populations of the Lesser Khingan Mountains group were used as control groups to carry out selective sweep analysis. We scanned the genome in 100 kb sliding windows with a step size of 10 kb to identify significant regions based on a fixed index (Fst) and nucleotide diversity ratio (π ratio). We respectively screened the significant regions according to the top 5%, and the intersection region is the selected region. Genes in the selected region were defined as the selected gene, namely, the candidate gene. After obtaining the candidate genes in the selected region, the genes were analyzed by GO (gene ontology) (www.geneontology.org, accessed on 13 January 2023) and KEGG enrichment analysis. Then, the variations in SNPs in the coding region of the candidate genes were annotated, and the core functional mutations were speculated to provide clues for further verification.

5. Conclusions

In our study, phylogenetic, population structure and selective sweep analyses of different populations of P. koraiensis in northeastern China were performed by GBS-seq sequencing technology. The results of the phylogenetic analysis showed that P. koraiensis was divided into three geoclimatic zones: the Great Khingan Mountains climatic region, the Lesser Khingan Mountains climatic region, and the Changbai Mountain climatic region, with introgression in the MKS population and outgroups. Geographic isolation, climatic factors, and introgression are the driving factors for the differentiation of P. koraiensis. The MKS population at the northern edge of the distribution area and the WYL population in the mining area are the two populations with greater differentiation. The differentiation of the MKS population is consistent with the central–marginal hypothesis, and the results of selective sweep analysis show that the genes subject to selection removal in the MKS population are related mainly to flowering and photomorphogenesis, cellular response to water deficit, response to reactive oxygen species, flavonoid synthesis, and DNA repair. Climatic heterogeneity factors (e.g., sunlight, precipitation, and temperature) drive the differentiation of MKS; genes selected in the WYL population are related mainly to metal ion transport and macromolecule biosynthesis, such as flavonoids, phytohormone signaling, and DNA repair, and the WYL population may be subject to heavy metal stress to produce differentiation. Our results provide insights into adaptive evolution molecular mechanisms of genus spruce and contribute to spruce resistance improvement via genomic-assisted breeding.

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/plants12061266/s1, Table S1: List of selected genes in P. koraiensis MKS (Mengkeshan) population; Table S2: List of selected genes in P. koraiensis WYL (Wuyiling) population; Table S3: GO categories of selected genes in P. koraiensis MKS (Mengkeshan) population; Table S4: GO categories of selected genes in P. koraiensis WYL (Wuyiling) population.

Author Contributions

Conceptualization, Y.W., Z.J. (Zirui Jia) and J.W.; data curation, Z.J. (Zirui Jia) and J.W.; formal analysis, Y.W. and Y.L.; funding acquisition, Z.J. (Zirui Jia); investigation, F.W., Y.W. and J.W.; project administration, Z.J. (Zirui Jia); resources, F.W., Y.W. and J.W.; software, Y.W., Y.L., W.N., C.T., Y.D., R.H. and Y.Y.; supervision, Z.J. (Zirui Jia), J.W. and Z.J. (Zeping Jiang); writing—original draft, Y.W.; writing—review and editing, Y.W., A.Q., Z.J. (Zirui Jia), Z.J. (Zeping Jiang) and E.C. All authors have read and agreed to the published version of the manuscript.

Funding

This study was supported by the “Study of Phylogeography and Quaternary glacial refuge in Picea koraiensis (31500540)” and the “study on Phylogeny and Phylogeography in Keteleeria (32271694)” funding by the National Natural Science Foundation of China.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The sequence reads associated with the paper are stored at the NCBI Sequence Read Archive in the BioProject PRJNA664562.

Acknowledgments

We thank Genedenovo Bioinformatics Technology Co., Ltd., for technical support on sequencing.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Schuster, W.S.F.; Mitton, J.B. Paternity and gene dispersal in limber pine (Pinus flexilis James). Heredity 2000, 84, 348–361. [Google Scholar] [CrossRef] [PubMed]
  2. Brown, G.R.; Gill, G.P.; Kuntz, R.J.; Langley, C.H.; Neale, D.B. Nucleotide diversity and linkage disequilibrium in loblolly pine. Proc. Natl. Acad. Sci. USA 2004, 101, 15255–15260. [Google Scholar] [CrossRef] [Green Version]
  3. Syring, J.; Farrell, K.; Businský, R.; Cronn, R.; Liston, A. Widespread genealogical nonmonophyly in species of pinus subgenus strobus. Syst. Biol. 2007, 56, 163–181. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  4. Ledig, F.T.; Hodgskiss, P.D.; Krutovskii, K.V.; Neale, D.B.; Eguiluz-Piedra, T. Relationships among the spruces (Picea, Pinaceae) of southwestern North America. Syst. Bot. 2004, 29, 275–295. [Google Scholar] [CrossRef] [Green Version]
  5. Ran, J.-H.; Shen, T.-T.; Wang, M.-M.; Wang, X.-Q. Phylogenomics resolves the deep phylogeny of seed plants and indicates partial convergent or homoplastic evolution between Gnetales and angiosperms. Proc. R. Soc. B Biol. Sci. 2018, 285, 20181012. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  6. Sun, Y.; Abbott, R.J.; Lu, Z.; Mao, K.; Zhang, L.; Wang, X.; Ru, D.; Liu, J. Reticulate evolution within a spruce (Picea) species complex revealed by population genomic analysis. Evolution 2018, 72, 2669–2681. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  7. Du, F.K.; Petit, R.J.; Liu, J.Q. More introgression with less gene flow: Chloroplast vs. mitochondrial DNA in the Picea asperata complex in China, and comparison with other Conifers. Mol. Ecol. 2009, 18, 1396–1407. [Google Scholar] [CrossRef]
  8. Wu, S.; Wang, Y.; Wang, Z.; Shrestha, N.; Liu, J. Species divergence with gene flow and hybrid speciation on the Qinghai–Tibet Plateau. New Phytol. 2022, 234, 392–404. [Google Scholar] [CrossRef]
  9. Ding, W.-N.; Ree, R.H.; Spicer, R.A.; Xing, Y.-W. Ancient orogenic and monsoon-driven assembly of the world’s richest temperate alpine flora. Science 2020, 369, 578–581. [Google Scholar] [CrossRef]
  10. Owens, G.L.; Rieseberg, L.H. Hybrid incompatibility is acquired faster in annual than in perennial species of sunflower and tarweed. Evolution 2014, 68, 893–900. [Google Scholar] [CrossRef]
  11. Feng, S.; Ru, D.; Sun, Y.; Mao, K.; Milne, R.; Liu, J. Trans-lineage polymorphism and nonbifurcating diversification of the genus Picea. New Phytol. 2019, 222, 576–587. [Google Scholar] [CrossRef] [Green Version]
  12. Ran, J.-H.; Wei, X.-X.; Wang, X.-Q. Molecular phylogeny and biogeography of Picea (Pinaceae): Implications for phylogeographical studies using cytoplasmic haplotypes. Mol. Phylogenetics Evol. 2006, 41, 405–419. [Google Scholar] [CrossRef] [PubMed]
  13. Sigurgeirsson, A.; Szmidt, A.E. Phylogenetic and biogeographic implications of chloroplast DNA variation in Picea. Nord. J. Bot. 1993, 13, 233–246. [Google Scholar] [CrossRef]
  14. Lockwood, J.D.; Aleksić, J.M.; Zou, J.; Wang, J.; Liu, J.; Renner, S.S. A new phylogeny for the genus Picea from plastid, mitochondrial, and nuclear sequences. Mol. Phylogenet. Evol. 2013, 69, 717–727. [Google Scholar] [CrossRef] [PubMed]
  15. Wright, J.W. Species crossability in spruce in relation to distribution and taxonomy. For. Sci. 1955, 1, 319–349. [Google Scholar]
  16. Nienstaedt, H. Genetics of White Spruce; US Department of Agriculture: Washington, DC, USA, 1972.
  17. Ren, Y.G.; Wan, C.B.; Qiao, X.Y.; Jiao, Y.G.; Chi, H.Y.; Zhao, C.B. Palynological assemblages of the Shahezi formation in Changwu area, north of Songliao Basin. J. Jilin Univ. 2014, 33, 406–412. [Google Scholar]
  18. Wang, R.; Wang, Y.; Chen, Y. Fossil woods from late cretaceous of Heilong-Jiang Province, Northeast China, and their palaeoenviromental implications. J. Integr. Plant Biol. 1997, 39, 10. [Google Scholar]
  19. Ding, Q.H.; Zhang, W.; Zheng, S.L. Research on fossil woods from the Fuxin formation in West Liaoning. Liaoning Geol. 2000, 17, 284–291. [Google Scholar] [CrossRef]
  20. Zhang, H.; Sun, L.; Han, J.; Feng, B. Study on the genetic diversity of Picea koraiensis. Zhiwu Yanjiu 2003, 23, 224–229. [Google Scholar]
  21. Jiang, X.; Gao, X.; Li, Y. Character diversity of mineral nutrition of koyama spruce populations. J. Northeast For. Univ. 2002, 30, 24–29. [Google Scholar]
  22. Wang, Q.; Ren, X.; Jiang, J. Genetic diversity for the provenance of Picea koraiensis by RAPD markers. J. Northeast For. Univ. 2004, 32, 1–3. [Google Scholar]
  23. Alexander, D.H.; Novembre, J.; Lange, K. Fast model-based estimation of ancestry in unrelated individuals. Genome Res. 2009, 19, 1655–1664. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  24. Yang, J.; Vázquez, L.; Feng, L.; Liu, Z.; Zhao, G. Climatic and soil factors shape the demographical history and genetic diversity of a deciduous oak (Quercus liaotungensis) in Northern China. Front. Plant Sci. 2018, 9, 1534. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  25. Costa e Silva, J.; Potts, B.; Harrison, P.A.; Bailey, T. Temperature and rainfall are separate agents of selection shaping population differentiation in a forest tree. Forests 2019, 10, 1145. [Google Scholar] [CrossRef] [Green Version]
  26. Đurović, S.Z.; Temunović, M.; Niketić, M.; Tomović, G.; Schönswetter, P.; Frajman, B.; Lavergne, S. Impact of Quaternary climatic oscillations on phylogeographic patterns of three habitat-segregated Cerastium taxa endemic to the Dinaric Alps. J. Biogeogr. 2021, 48, 2022–2036. [Google Scholar] [CrossRef]
  27. Aizawa, M.; Yoshimaru, H.; Saito, H.; Katsuki, T.; Kawahara, T.; Kitamura, K.; Shi, F.; Kaji, M. Phylogeography of a northeast Asian spruce, Picea jezoensis, inferred from genetic variation observed in organelle DNA markers. Mol. Ecol. 2007, 16, 3393–3405. [Google Scholar] [CrossRef]
  28. Ran, J.-H.; Shen, T.-T.; Liu, W.-J.; Wang, P.-P.; Wang, X.-Q. Mitochondrial introgression and complex biogeographic history of the genus Picea. Mol. Phylogenet. Evol. 2015, 93, 63–76. [Google Scholar] [CrossRef]
  29. Sullivan, A.R.; Eldfjell, Y.; Schiffthaler, B.; Delhomme, N.; Asp, T.; Hebelstrup, K.H.; Keech, O.; Öberg, L.; Møller, I.M.; Arvestad, L.; et al. The mitogenome of norway spruce and a reappraisal of mitochondrial recombination in plants. Genome Biol. Evol. 2019, 12, 3586–3598. [Google Scholar] [CrossRef]
  30. Sutton, B.C.S.; Flanagan, D.J.; Gawley, J.R.; Newton, C.H.; Lester, D.T.; El-Kassaby, Y.A. Inheritance of chloroplast and mitochondrial DNA in Picea and composition of hybrids from introgression zones. Theor. Appl. Genet. 1991, 82, 242–248. [Google Scholar] [CrossRef] [PubMed]
  31. Rajora, O.P.; Dancik, B.P. Population genetic variation, structure, and evolution in Engelmann spruce, white spruce, and their natural hybrid complex in Alberta. Can. J. Bot. 2000, 78, 768–780. [Google Scholar] [CrossRef]
  32. Hamilton, J.A.; Lexer, C.; Aitken, S.N. Genomic and phenotypic architecture of a spruce hybrid zone (Picea sitchensis × P. glauca). Mol. Ecol. 2013, 22, 827–841. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  33. De La Torre, A.R.; Lin, Y.-C.; Van de Peer, Y.; Ingvarsson, P.K. Genome-wide analysis reveals diverged patterns of codon bias, gene expression, and rates of sequence evolution in Picea gene families. Genome Biol. Evol. 2015, 7, 1002–1015. [Google Scholar] [CrossRef] [Green Version]
  34. de Lafontaine, G.; Prunier, J.; Gérardi, S.; Bousquet, J. Tracking the progression of speciation: Variable patterns of introgression across the genome provide insights on the species delimitation between progenitor-derivative spruces (Picea mariana × P. rubens). Mol. Ecol. 2015, 24, 5229–5247. [Google Scholar] [CrossRef] [PubMed]
  35. Tsuda, Y.; Chen, J.; Stocks, M.; Källman, T.; Sønstebø, J.H.; Parducci, L.; Semerikov, V.; Sperisen, C.; Politov, D.; Ronkainen, T.; et al. The extent and meaning of hybridization and introgression between Siberian spruce (Picea obovata) and Norway spruce (Picea abies): Cryptic refugia as stepping stones to the west? Mol. Ecol. 2016, 25, 2773–2789. [Google Scholar] [CrossRef]
  36. Du, F.K.; Peng, X.L.; Liu, J.Q.; Lascoux, M.; Hu, F.S.; Petit, R.J. Direction and extent of organelle DNA introgression between two spruce species in the Qinghai-Tibetan Plateau. New Phytol. 2011, 192, 1024–1033. [Google Scholar] [CrossRef]
  37. Li, Y.; Stocks, M.; Hemmilä, S.; Källman, T.; Zhu, H.; Zhou, Y.; Chen, J.; Liu, J.; Lascoux, M. Demographic histories of four spruce (Picea) species of the Qinghai-Tibetan Plateau and neighboring areas inferred from multiple nuclear loci. Mol. Biol. Evol. 2010, 27, 1001–1014. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  38. Zou, J.; Sun, Y.; Li, L.; Wang, G.; Yue, W.; Lu, Z.; Wang, Q.; Liu, J. Population genetic evidence for speciation pattern and gene flow between Picea wilsonii, P. morrisonicola and P. neoveitchii. Ann. Bot. 2013, 112, 1829–1844. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  39. Kitamura, K.; Uchiyama, K.; Ueno, S.; Ishizuka, W.; Tsuyama, I.; Goto, S. Geographical gradients of genetic diversity and differentiation among the southernmost marginal populations of Abies sachalinensis revealed by EST-SSR polymorphism. Forests 2020, 11, 233. [Google Scholar] [CrossRef] [Green Version]
  40. Tóth, E.G.; Tremblay, F.; Housset, J.M.; Bergeron, Y.; Carcaillet, C. Geographic isolation and climatic variability contribute to genetic differentiation in fragmented populations of the long-lived subalpine conifer Pinus cembra L. in the western Alps. BMC Evol. Biol. 2019, 19, 190. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  41. Mengoni, A.; Barabesi, C.; Gonnelli, C.; Galardi, F.; Gabbrielli, R.; Bazzicalupo, M. Genetic diversity of heavy metal-tolerant populations in Silene paradoxa L. (Caryophyllaceae): A chloroplast microsatellite analysis. Mol. Ecol. 2001, 10, 1909–1916. [Google Scholar] [CrossRef]
  42. Mengoni, A.; Gonnelli, C.; Galardi, F.; Gabbrielli, R.; Bazzicalupo, M. Genetic diversity and heavy metal tolerance in populations of Silene paradoxa L. (Caryophyllaceae): A random amplified polymorphic DNA analysis. Mol. Ecol. 2000, 9, 1319–1324. [Google Scholar] [CrossRef]
  43. Boquete, M.T.; Schmid, M.W.; Wagemaker, N.C.A.M.; Carey, S.B.; McDaniel, S.F.; Richards, C.; Alonso, C. Molecular basis of intraspecific differentiation for heavy metal tolerance in the copper moss Scopelophila cataractae. Environ. Exp. Bot. 2022, 201, 104970. [Google Scholar] [CrossRef]
  44. Boonmee, S.; Neeratanaphan, L.; Tanee, T.; Khamon, P. The genetic differentiation of Colocasia esculenta growing in gold mining areas with arsenic contamination. Environ. Monit. Assess. 2015, 187, 227–235. [Google Scholar] [CrossRef] [PubMed]
  45. Luo, D.; Xu, B.; Li, Z.-M.; Sun, H. The ‘Ward Line–Mekong–Salween Divide’ is an important floristic boundary between the eastern Himalaya and Hengduan mountains: Evidence from the phylogeographical structure of subnival herbs Marmoritis complanatum (Lamiaceae). Bot. J. Linn. Soc. 2017, 185, 482–496. [Google Scholar] [CrossRef]
  46. Zhang, X.; Sun, Y.; Landis, J.B.; Zhang, J.; Yang, L.; Lin, N.; Zhang, H.; Guo, R.; Li, L.; Zhang, Y.; et al. Genomic insights into adaptation to heterogeneous environments for the ancient relictual Circaeaster agrestis (Circaeasteraceae, Ranunculales). New Phytol. 2020, 228, 285–301. [Google Scholar] [CrossRef]
  47. Zhu, L.; Bu, Q.; Xu, X.; Paik, I.; Huang, X.; Hoecker, U.; Deng, X.W.; Huq, E. CUL4 forms an E3 ligase with COP1 and SPA to promote light-induced degradation of PIF1. Nat. Commun. 2015, 6, 7245. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  48. Chen, Y.; Hu, X.; Liu, S.; Su, T.; Huang, H.; Ren, H.; Gao, Z.; Wang, X.; Lin, D.; Wang, Q.; et al. Regulation of Arabidopsis photoreceptor CRY2 by two distinct E3 ubiquitin ligases. Nat. Commun. 2021, 12, 2155. [Google Scholar] [CrossRef] [PubMed]
  49. Manfield, I.W.; Devlin, P.F.; Jen, C.-H.; Westhead, D.R.; Gilmartin, P.M. Conservation, convergence, and divergence of light-responsive, circadian-regulated, and tissue-specific expression patterns during evolution of the Arabidopsis GATA gene family. Plant Physiol. 2007, 143, 941–958. [Google Scholar] [CrossRef] [Green Version]
  50. Jones, A.M.; Xuan, Y.; Xu, M.; Wang, R.-S.; Ho, C.-H.; Lalonde, S.; You, C.H.; Sardi, M.I.; Parsa, S.A.; Smith-Valle, E.; et al. Border control—a membrane-linked interactome of Arabidopsis. Science 2014, 344, 711–716. [Google Scholar] [CrossRef]
  51. Yaish, M.W.F.; Peng, M.; Rothstein, S.J. AtMBD9 modulates Arabidopsis development through the dual epigenetic pathways of DNA methylation and histone acetylation. Plant J. 2009, 59, 123–135. [Google Scholar] [CrossRef]
  52. Peng, M.; Cui, Y.; Bi, Y.-M.; Rothstein, S.J. AtMBD9: A protein with a methyl-CpG-binding domain regulates flowering time and shoot branching in Arabidopsis. Plant J. 2006, 46, 282–296. [Google Scholar] [CrossRef]
  53. Shimizu, N.; Ishida, T.; Yamada, M.; Shigenobu, S.; Tabata, R.; Kinoshita, A.; Yamaguchi, K.; Hasebe, M.; Mitsumasu, K.; Sawa, S. BAM 1 and RECEPTOR-LIKE PROTEIN KINASE 2 constitute a signaling pathway and modulate CLE peptide-triggered growth inhibition in Arabidopsis root. New Phytol. 2015, 208, 1104–1113. [Google Scholar] [CrossRef]
  54. Barth, C.; Gouzd, Z.A.; Steele, H.P.; Imperio, R.M. A mutation in GDP-mannose pyrophosphorylase causes conditional hypersensitivity to ammonium, resulting in Arabidopsis root growth inhibition, altered ammonium metabolism, and hormone homeostasis. J. Exp. Bot. 2010, 61, 379–394. [Google Scholar] [CrossRef] [Green Version]
  55. Huang, C.; He, W.; Guo, J.; Chang, X.; Su, P.; Zhang, L. Increased sensitivity to salt stress in an ascorbate-deficient Arabidopsis mutant. J. Exp. Bot. 2005, 56, 3041–3049. [Google Scholar] [CrossRef] [Green Version]
  56. Joshi, V.; Jander, G. Arabidopsis methionine gamma-lyase is regulated according to isoleucine biosynthesis needs but plays a subordinate role to threonine deaminase. Plant Physiol. 2009, 151, 367–378. [Google Scholar] [CrossRef] [Green Version]
  57. Su, Y.; Li, M.; Guo, L.; Wang, X. Different effects of phospholipase Dζ2 and non-specific phospholipase C4 on lipid remodeling and root hair growth in Arabidopsis response to phosphate deficiency. Plant J. 2018, 94, 315–326. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  58. Sagar, S.; Singh, A. Emerging role of phospholipase C mediated lipid signaling in abiotic stress tolerance and development in plants. Plant Cell Rep. 2021, 40, 2123–2133. [Google Scholar] [CrossRef] [PubMed]
  59. Zhao, A.; Wang, X.; Dong, J.; Hou, Z.; Zhao, W.; Gao, X.; Yang, M. Advances in the Study of Phospholipase C Response to Stress in Plants. Biotechnol. Bull. 2021, 37, 154. [Google Scholar] [CrossRef]
  60. Xia, X.; Yang, Z.; Cui, Y.; Li, Y.; Hou, Q.; Yu, T. Soil heavy metal concentrations and their typical input and output fluxes on the southern Song-nen Plain, Heilongjiang Province, China. J. Geochem. Explor. 2014, 139, 85–96. [Google Scholar] [CrossRef]
  61. Zhao, W.; Shan, Y.; Zhao, L.; Nie, L. Spatial distribution characteristics of five heavy metals elements in rice: A case study for different basins in Heilongjiang Province. In Advances in Applied Chemistry and Industrial Catalysis; CRC Press: Boca Raton, FL, USA, 2022; pp. 390–397. [Google Scholar]
  62. Krämer, U.; Clemens, S. Functions and homeostasis of zinc, copper, and nickel in plants. In Topics in Current Genetics; Springer Verlag: New York, NY, USA, 2005; pp. 216–271. [Google Scholar]
  63. Manara, A. Plant responses to heavy metal toxicity. In Plants and Heavy Metals; Furini, A., Ed.; Springer: Dordrecht, The Netherlands, 2012; pp. 27–53. [Google Scholar]
  64. Tommasini, R.; Vogt, E.; Fromenteau, M.; Hörtensteiner, S.; Matile, P.; Amrhein, N.; Martinoia, E. An ABC-transporter of Arabidopsis thaliana has both glutathione-conjugate and chlorophyll catabolite transport activity. Plant J. 1998, 13, 773–780. [Google Scholar] [CrossRef] [PubMed]
  65. Kim, D.Y.; Bovet, L.; Maeshima, M.; Martinoia, E.; Lee, Y. The ABC transporter AtPDR8 is a cadmium extrusion pump conferring heavy metal resistance. Plant J. 2007, 50, 207–218. [Google Scholar] [CrossRef] [PubMed]
  66. He, G.; Tian, W.; Qin, L.; Meng, L.; Wu, D.; Huang, Y.; Li, D.; Zhao, D.; He, T. Identification of novel heavy metal detoxification proteins in Solanum tuberosum: Insights to improve food security protection from metal ion stress. Sci. Total Environ. 2021, 779, 146197. [Google Scholar] [CrossRef]
  67. Tilley, D.M.; Evans, C.R.; Larson, T.M.; Edwards, K.A.; Friesen, J.A. Identification and characterization of the nuclear isoform of Drosophila melanogaster CTP: Phosphocholine cytidylyltransferase. Biochemistry 2008, 47, 11838–11846. [Google Scholar] [CrossRef]
  68. Fullerton, M.D.; Hakimuddin, F.; Bakovic, M. Developmental and metabolic effects of disruption of the mouse CTP: Phosphoethanolamine cytidylyltransferase gene (Pcyt2). Mol. Cell. Biol. 2007, 27, 3327–3336. [Google Scholar] [CrossRef] [Green Version]
  69. Gaudet, P.; Livstone, M.S.; Lewis, S.E.; Thomas, P.D. Phylogenetic-based propagation of functional annotations within the Gene Ontology consortium. Brief. Bioinform. 2011, 12, 449–462. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  70. Morash, S.C.; McMaster, C.R.; Hjelmstad, R.H.; Bell, R.M. Studies employing Saccharomyces cerevisiae cpt1 and ept1 null mutants implicate the CPT1 gene in coordinate regulation of phospholipid biosynthesis. J. Biol. Chem. 1994, 269, 28769–28776. [Google Scholar] [CrossRef]
  71. Yap, W.S.; Shyu, P.; Gaspar, M.L.; Jesch, S.A.; Marvalim, C.; Prinz, W.A.; Henry, S.A.; Thibault, G. Yeast FIT2 homolog is necessary to maintain cellular proteostasis by regulating lipid homeostasis. J. Cell Sci. 2020, 133, jcs248526. [Google Scholar] [CrossRef] [PubMed]
  72. Wild, R.; Gerasimaite, R.; Jung, J.-Y.; Truffault, V.; Pavlovic, I.; Schmidt, A.; Saiardi, A.; Jessen, H.J.; Poirier, Y.; Hothorn, M.; et al. Control of eukaryotic phosphate homeostasis by inositol polyphosphate sensor domains. Science 2016, 352, 986–990. [Google Scholar] [CrossRef] [Green Version]
  73. Wege, S.; Khan, G.A.; Jung, J.-Y.; Vogiatzaki, E.; Pradervand, S.; Aller, I.; Meyer, A.J.; Poirier, Y. The EXS domain of PHO1 participates in the response of shoots to phosphate deficiency via a root-to-shoot signal. Plant Physiol. 2016, 170, 385–400. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  74. Stracke, R.; De Vos, R.C.H.; Bartelniewoehner, L.; Ishihara, H.; Sagasser, M.; Martens, S.; Weisshaar, B. Metabolomic and genetic analyses of flavonol synthesis in Arabidopsis thaliana support the in vivo involvement of leucoanthocyanidin dioxygenase. Planta 2009, 229, 427–445. [Google Scholar] [CrossRef]
  75. Prescott, A.G.; Stamford, N.P.J.; Wheeler, G.; Firmin, J.L. In vitro properties of a recombinant flavonol synthase from Arabidopsis thaliana. Phytochemistry 2002, 60, 589–593. [Google Scholar] [CrossRef]
  76. Gutha, L.R.; Casassa, L.F.; Harbertson, J.F.; Naidu, R.A. Modulation of flavonoid biosynthetic pathway genes and anthocyanins due to virus infection in grapevine (Vitis vinifera L.) leaves. BMC Plant Biol. 2010, 10, 1–18. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  77. Satishchandran, C.; Markham, G.D. Mechanistic studies of Escherichia coli adenosine-5′-phosphosulfate kinase. Arch. Biochem. Biophys. 2000, 378, 210–215. [Google Scholar] [CrossRef]
  78. Ganea, E.; Harding, J.J. Glutathione-related enzymes and the eye. Curr. Eye Res. 2006, 31, 1–11. [Google Scholar] [CrossRef]
  79. Peng, C.; Xia, Y.; Zhang, W.; Wu, Y.; Shi, J. Proteomic analysis unravels response and antioxidation defense mechanism of rice plants to copper oxide nanoparticles: Comparison with bulk particles and dissolved Cu ions. ACS Agric. Sci. Technol. 2022, 2, 671–683. [Google Scholar] [CrossRef]
  80. Nebert, D.W.; Vasiliou, V. Analysis of the glutathione S-transferase (GST) gene family. Hum. Genom. 2004, 1, 1–5. [Google Scholar] [CrossRef]
  81. Teng, Z.; Zheng, W.; Jiang, S.; Hong, S.-B.; Zhu, Z.; Zang, Y. Role of melatonin in promoting plant growth by regulating carbon assimilation and ATP accumulation. Plant Sci. 2022, 319, 111276. [Google Scholar] [CrossRef] [PubMed]
  82. Emamverdian, A.; Ding, Y.; Mokhberdoran, F.; Xie, Y. Heavy metal stress and some mechanisms of plant defense response. Sci. World J. 2015, 2015, 1–18. [Google Scholar] [CrossRef] [Green Version]
  83. Saba, H.; Jyoti, P.; Neha, S. Mycorrhizae and phytochelators as remedy in heavy metal contaminated land remediation. Int. Res. J. Environ. Sci. 2013, 2, 74–78. [Google Scholar]
  84. Liu, Q.; Zhang, Y.; Wang, Y.; Wang, W.; Gu, C.; Huang, S.; Yuan, H.; Dhankher, O.P. Quantitative proteomic analysis reveals complex regulatory and metabolic response of Iris lactea Pall. var. chinensis to cadmium toxicity. J. Hazard. Mater. 2020, 400, 123165. [Google Scholar] [CrossRef]
  85. Borges, K.L.R.; Salvato, F.; Alcântara, B.K.; Nalin, R.S.; Piotto, F.Â.; Azevedo, R.A. Temporal dynamic responses of roots in contrasting tomato genotypes to cadmium tolerance. Ecotoxicology 2018, 27, 245–258. [Google Scholar] [CrossRef]
  86. Xu, L.; Wang, Y.; Zhang, F.; Tang, M.; Chen, Y.; Wang, J.; Karanja, B.K.; Luo, X.; Zhang, W.; Liu, L. Dissecting root proteome changes reveals new insight into cadmium stress response in radish (Raphanus sativus L.). Plant Cell Physiol. 2017, 58, 1901–1913. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  87. Asai, T.; Tena, G.; Plotnikova, J.; Willmann, M.R.; Chiu, W.-L.; Gomez-Gomez, L.; Boller, T.; Ausubel, F.M.; Sheen, J. MAP kinase signalling cascade in Arabidopsis innate immunity. Nature 2002, 415, 977–983. [Google Scholar] [CrossRef] [PubMed]
  88. Pitzschke, A.; Schikora, A.; Hirt, H. MAPK cascade signalling networks in plant defence. Curr. Opin. Plant Biol. 2009, 12, 421–426. [Google Scholar] [CrossRef]
  89. Takahashi, F.; Mizoguchi, T.; Yoshida, R.; Ichimura, K.; Shinozaki, K. Calmodulin-dependent activation of MAP kinase for ROS homeostasis in Arabidopsis. Mol. Cell 2011, 41, 649–660. [Google Scholar] [CrossRef]
  90. Sethi, V.; Raghuram, B.; Sinha, A.K.; Chattopadhyay, S. A mitogen-activated protein kinase cascade module, MKK3-MPK6 and MYC2, is involved in blue light-mediated seedling development in Arabidopsis. Plant Cell 2014, 26, 3343–3357. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  91. Vitti, A.; Nuzzaci, M.; Scopa, A.; Tataranni, G.; Remans, T.; Vangronsveld, J.; Sofo, A. Auxin and cytokinin metabolism and root morphological modifications in Arabidopsis thaliana seedlings infected with Cucumber mosaic virus (CMV) or exposed to cadmium. Int. J. Mol. Sci. 2013, 14, 6889–6902. [Google Scholar] [CrossRef] [Green Version]
  92. Farinati, S.; DalCorso, G.; Varotto, S.; Furini, A. The Brassica juncea BjCdR15, an ortholog of Arabidopsis TGA3, is a regulator of cadmium uptake, transport and accumulation in shoots and confers cadmium tolerance in transgenic plants. New Phytol. 2010, 185, 964–978. [Google Scholar] [CrossRef]
  93. Wang, Y.; Gao, C.; Liang, Y.; Wang, C.; Yang, C.; Liu, G. A novel bZIP gene from Tamarix hispida mediates physiological responses to salt stress in tobacco plants. J. Plant Physiol. 2010, 167, 222–230. [Google Scholar] [CrossRef] [PubMed]
  94. Smeets, K.; Opdenakker, K.; Remans, T.; Forzani, C.; Hirt, H.; Vangronsveld, J.; Cuypers, A.N.N. The role of the kinase OXI1 in cadmium- and copper-induced molecular responses in Arabidopsis thaliana. Plant Cell Environ. 2013, 36, 1228–1238. [Google Scholar] [CrossRef]
  95. van De Mortel, J.E.; Schat, H.; Moerland, P.D.; Van Themaat, E.V.L.; Van Der Ent, S.; Blankestijn, H.; Ghandilyan, A.; Tsiatsiani, S.; Aarts, M.G.M. Expression differences for genes involved in lignin, glutathione and sulphate metabolism in response to cadmium in Arabidopsis thaliana and the related Zn/Cd-hyperaccumulator Thlaspi caerulescens. Plant Cell Environ. 2008, 31, 301–324. [Google Scholar] [CrossRef]
  96. Opdenakker, K.; Remans, T.; Vangronsveld, J.; Cuypers, A. Mitogen-Activated Protein (MAP) kinases in plant metal stress: Regulation and responses in comparison to other biotic and abiotic stresses. Int. J. Mol. Sci. 2012, 13, 7828–7853. [Google Scholar] [CrossRef]
  97. Ramos, J.; Clemente, M.R.; Naya, L.; Loscos, J.; Pérez-Rontomé, C.; Sato, S.; Tabata, S.; Becana, M. Phytochelatin synthases of the model legume Lotus japonicus. A small multigene family with differential response to cadmium and alternatively spliced variants. Plant Physiol. 2007, 143, 1110–1118. [Google Scholar] [CrossRef] [Green Version]
  98. Gracheva, E.; Chitale, S.; Wilhelm, T.; Rapp, A.; Byrne, J.; Stadler, J.; Medina, R.; Cardoso, M.C.; Richly, H. ZRF1 mediates remodeling of E3 ligases at DNA lesion sites during nucleotide excision repair. J. Cell Biol. 2016, 213, 185–200. [Google Scholar] [CrossRef] [Green Version]
  99. Bélanger, F.; Angers, J.-P.; Fortier, É.; Hammond-Martel, I.; Costantino, S.; Drobetsky, E.; Wurtele, H. Mutations in replicative stress response pathways are associated with s phase-specific defects in nucleotide excision repair. J. Biol. Chem. 2016, 291, 522–537. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  100. Lipkin, S.M.; Wang, V.; Jacoby, R.; Banerjee-Basu, S.; Baxevanis, A.D.; Lynch, H.T.; Elliott, R.M.; Collins, F.S. MLH3: A DNA mismatch repair gene associated with mammalian microsatellite instability. Nat. Genet. 2000, 24, 27–35. [Google Scholar] [CrossRef] [PubMed]
  101. de Barros, A.C.; Takeda, A.A.S.; Dreyer, T.R.; Velazquez-Campoy, A.; Kobe, B.; Fontes, M.R.M. DNA mismatch repair proteins MLH1 and PMS2 can be imported to the nucleus by a classical nuclear import pathway. Biochimie 2018, 146, 87–96. [Google Scholar] [CrossRef] [PubMed]
  102. Hartung, F.; Suer, S.; Puchta, H. Two closely related RecQ helicases have antagonistic roles in homologous recombination and DNA repair in Arabidopsis thaliana. Proc. Natl. Acad. Sci. USA 2007, 104, 18836–18841. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  103. Holloman, W.K. Unraveling the mechanism of BRCA2 in homologous recombination. Nat. Struct. Mol. Biol. 2011, 18, 748–754. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  104. Wei, N.; Kwok, S.F.; von Arnim, A.G.; Lee, A.; McNellis, T.W.; Piekos, B.; Deng, X.W. Arabidopsis COP8, COP10, and COP11 genes are involved in repression of photomorphogenic development in darkness. Plant Cell 1994, 6, 629–643. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  105. Franciosini, A.; Moubayidin, L.; Du, K.; Matari, N.H.; Boccaccini, A.; Butera, S.; Vittorioso, P.; Sabatini, S.; Jenik, P.D.; Costantino, P.; et al. The COP9 SIGNALOSOME is required for postembryonic meristem maintenance in Arabidopsis thaliana. Mol. Plant. 2015, 8, 1623–1634. [Google Scholar] [CrossRef] [Green Version]
  106. Manara, A.; DalCorso, G.; Guzzo, F.; Furini, A. Loss of the Atypical Kinases ABC1K7 and ABC1K8 Changes the Lipid Composition of the Chloroplast Membrane. Plant. Cell. Physiol. 2015, 56, 1193–1204. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  107. Staneloni, R.J.; Rodriguez-Batiller, M.J.; Legisa, D.; Scarpin, M.R.; Agalou, A.; Cerdán, P.D.; Meijer, A.H.; Ouwerkerk, P.B.; Casal, J.J. Bell-like homeodomain selectively regulates the high-irradiance response of phytochrome A. Proc. Natl. Acad. Sci. USA 2009, 106, 13624–13629. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  108. Perea-Resa, C.; Carrasco-López, C.; Catalá, R.; Turečková, V.; Novak, O.; Zhang, W.; Sieburth, L.; Jiménez-Gómez, J.M.; Salinas, J. The LSM1-7 Complex Differentially Regulates Arabidopsis Tolerance to Abiotic Stress Conditions by Promoting Selective mRNA Decapping. Plant Cell 2016, 28, 505–520. [Google Scholar] [CrossRef] [Green Version]
  109. Kanwar, P.; Sanyal, S.K.; Mahiwal, S.; Ravi, B.; Kaur, K.; Fernandes, J.L.; Yadav, A.K.; Tokas, I.; Srivastava, A.K.; Suprasanna, P.; et al. CIPK9 targets VDAC3 and modulates oxidative stress responses in Arabidopsis. Plant J. 2022, 109, 241–260. [Google Scholar] [CrossRef]
  110. Xu, T.; Lee, K.; Gu, L.; Kim, J.I.; Kang, H. Functional characterization of a plastid-specific ribosomal protein PSRP2 in Arabidopsis thaliana under abiotic stress conditions. Plant Physiol. Biochem. 2013, 73, 405–411. [Google Scholar] [CrossRef]
  111. Poland, J.A.; Brown, P.J.; Sorrells, M.E.; Jannink, J.-L. Development of high-density genetic maps for barley and wheat using a novel two-enzyme genotyping-by-sequencing approach. PLoS ONE 2012, 7, e32253. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  112. Li, H.; Durbin, R. Fast and accurate short read alignment with Burrows–Wheeler transform. Bioinformatics 2009, 25, 1754–1760. [Google Scholar] [CrossRef] [Green Version]
  113. Danecek, P.; Auton, A.; Abecasis, G.; Albers, C.A.; Banks, E.; DePristo, M.A.; Handsaker, R.E.; Lunter, G.; Marth, G.T.; Sherry, S.T.; et al. The variant call format and VCFtools. Bioinformatics 2011, 27, 2156–2158. [Google Scholar] [CrossRef] [PubMed]
  114. Jin, L.N.; Nei, M. Limitations of the evolutionary parsimony method of phylogenetic analysis. Mol. Biol. Evol. 1990, 7, 82–102. [Google Scholar] [CrossRef]
  115. Kumar, S.; Stecher, G.; Li, M.; Knyaz, C.; Tamura, K. MEGA X: Molecular evolutionary genetics analysis across computing platforms. Mol. Biol. Evol. 2018, 35, 1547–1549. [Google Scholar] [CrossRef] [PubMed]
  116. Letunic, I.; Bork, P. Interactive tree of life (iTOL) v5: An online tool for phylogenetic tree display and annotation. Nucleic Acids Res. 2021, 49, W293–W296. [Google Scholar] [CrossRef] [PubMed]
Plants 12 01266 g001 550
Figure 1. Phylogenetic tree and structure of 127 samples of P. koraiensis and outgroups with bootstrap value. PUN: Picea pungens; LIK: Picea likiangensis; JEZ: Picea jezoensis. I: Great Khingan Mountains group; II: Changbai Mountain group; III: Lesser Khingan Mountains group.
Figure 1. Phylogenetic tree and structure of 127 samples of P. koraiensis and outgroups with bootstrap value. PUN: Picea pungens; LIK: Picea likiangensis; JEZ: Picea jezoensis. I: Great Khingan Mountains group; II: Changbai Mountain group; III: Lesser Khingan Mountains group.
Plants 12 01266 g001
Plants 12 01266 g002 550
Figure 2. Cross-validation error value based on ADMIXTURE v1.23 software. The best number of clusters was 3 according to the minimum cross-validation rate.
Figure 2. Cross-validation error value based on ADMIXTURE v1.23 software. The best number of clusters was 3 according to the minimum cross-validation rate.
Plants 12 01266 g002
Plants 12 01266 g003 550
Figure 3. Principal component analysis of 9 populations of P. koraiensis. (a) PC1 and pPC2; (b) PC2 and PC3. PC: principal component.
Figure 3. Principal component analysis of 9 populations of P. koraiensis. (a) PC1 and pPC2; (b) PC2 and PC3. PC: principal component.
Plants 12 01266 g003
Plants 12 01266 g004 550
Figure 4. Genomic regions with selective sweep in MKS and WYL population: (a) MKS population; (b) WYL population.
Figure 4. Genomic regions with selective sweep in MKS and WYL population: (a) MKS population; (b) WYL population.
Plants 12 01266 g004
Plants 12 01266 g005 550
Figure 5. Annotation based on the GO database. The Y-axis indicates the top 30 GO terms, and the X-axis indicates the numbers of genes annotated to the terms and their P values: (a) MKS population; (b) WYL population.
Figure 5. Annotation based on the GO database. The Y-axis indicates the top 30 GO terms, and the X-axis indicates the numbers of genes annotated to the terms and their P values: (a) MKS population; (b) WYL population.
Plants 12 01266 g005
Plants 12 01266 g006 550
Figure 6. KEGG enrichment analysis of selected genes. The color scale indicates the p value: (a) MKS population; (b) WYL population.
Figure 6. KEGG enrichment analysis of selected genes. The color scale indicates the p value: (a) MKS population; (b) WYL population.
Plants 12 01266 g006
Plants 12 01266 g007 550
Figure 7. Geographical distributions of 9 populations of P. koraiensis with altitude as the base map.
Figure 7. Geographical distributions of 9 populations of P. koraiensis with altitude as the base map.
Plants 12 01266 g007
Table
Table 1. Fst between 9 populations, π and Tajima. D neutral test in P. koraiensis.
Table 1. Fst between 9 populations, π and Tajima. D neutral test in P. koraiensis.
PopulationGFMKSMLHLTQLAHWYLLJZHπ × 10−2Tajima’s D
GF0 0.2181 0.3892
MKS0.1018 0 0.2285 0.1205
ML0.0643 0.0878 0 0.2289 0.0526
HL0.0628 0.0909 0.0399 0 0.2272 0.1307
TQL0.0599 0.0838 0.0294 0.0371 0 0.2319 0.0297
AH0.0529 0.0938 0.0500 0.0544 0.0462 0 0.2212 0.1657
WYL0.1282 0.1569 0.1166 0.1182 0.1115 0.1191 0 0.19430.5110
LJ0.0623 0.0873 0.0322 0.0399 0.0254 0.0488 0.1137 0 0.2351 0.0507
ZH0.0491 0.0842 0.0408 0.0416 0.0371 0.0450 0.1042 0.0404 0 0.2314 0.1875
Table
Table 2. Locations of 9 P. koraiensis populations.
Table 2. Locations of 9 P. koraiensis populations.
PopulationNO.Collection SiteLatitudeLongitude
GF10Mengke Mountain, Heilongjiang52.63124.31
MKS13Hailin, Heilongjiang44.42128.54
ML15Muling, Heilongjiang44.21130.21
HL15Linjiang, Jilin 41.99127.24
TQL13Tianqiaoling, Jilin 43.36129.37
AH14Wuyinling, Heilongjiang 48.67129.42
WYL14Zhanhe, Heilongjiang 48.33128.07
LJ11Aihui, Heilongjiang 49.62126.81
ZH15Gaofeng, Heilongjiang 49.43125.23
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Wang, Y.; Jiang, Z.; Qin, A.; Wang, F.; Chang, E.; Liu, Y.; Nie, W.; Tan, C.; Yuan, Y.; Dong, Y.; et al. Population Structure, Genetic Diversity and Candidate Genes for the Adaptation to Environmental Stress in Picea koraiensis. Plants 2023, 12, 1266. https://doi.org/10.3390/plants12061266

AMA Style

Wang Y, Jiang Z, Qin A, Wang F, Chang E, Liu Y, Nie W, Tan C, Yuan Y, Dong Y, et al. Population Structure, Genetic Diversity and Candidate Genes for the Adaptation to Environmental Stress in Picea koraiensis. Plants. 2023; 12(6):1266. https://doi.org/10.3390/plants12061266

Chicago/Turabian Style

Wang, Ya, Zeping Jiang, Aili Qin, Fude Wang, Ermei Chang, Yifu Liu, Wen Nie, Cancan Tan, Yanchao Yuan, Yao Dong, and et al. 2023. "Population Structure, Genetic Diversity and Candidate Genes for the Adaptation to Environmental Stress in Picea koraiensis" Plants 12, no. 6: 1266. https://doi.org/10.3390/plants12061266

APA Style

Wang, Y., Jiang, Z., Qin, A., Wang, F., Chang, E., Liu, Y., Nie, W., Tan, C., Yuan, Y., Dong, Y., Huang, R., Jia, Z., & Wang, J. (2023). Population Structure, Genetic Diversity and Candidate Genes for the Adaptation to Environmental Stress in Picea koraiensis. Plants, 12(6), 1266. https://doi.org/10.3390/plants12061266

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop