Haplotype-resolved genome of a papeda provides insights into the geographical origin and evolution of Citrus

Abstract

The publication of several high-quality genomes has contributed greatly to clarifying the evolution of citrus. However, due to their complex genetic backgrounds, the origins and evolution of many citrus species remain unclear. We assembled de novo the 294-Mbp chromosome-level genome of a more than 200-year-old primitive papeda (DYC002). Comparison between the two sets of homologous chromosomes of the haplotype-resolved genome revealed 1.2% intragenomic variations, including 1.75 million SNPs, 149,471 insertions and 154,215 deletions. Using this genome as a reference, we resequenced and performed population and phylogenetic analyses of 378 representative citrus accessions. Our study confirmed that the primary origin center of core Citrus species is in South China, particularly in the Himalaya-Hengduan Mountains. Papeda species are an ancient Citrus type compared with C. ichangensis. We found that the evolution of the Citrus genus followed two radiations through two routes (to East China and Southeast Asia) along river systems. Evidence for the origin and evolution of some individual citrus species was provided. Papeda probably played an important role in the origins of Australian finger lime, citrons, Honghe papeda and pummelos; Ichang papeda originated from Yuanjiang city of Yunnan Province, China, and C. mangshanensis has a close relationship with kumquat and Ichang papeda. Moreover, the Hunan and Guangdong Provinces of China are predicted to be the origin center of mandarin, sweet orange and sour orange. Additionally, our study revealed that fruit bitterness was significantly selected against during citrus domestication. Taken together, this study provides new insight into the origin and evolution of citrus species and may serve as a valuable genomic resource for citrus breeding and improvement.

Keywords: citrus; domestication; evolution; haplotype‐resolved genome; papeda.

Figure 1

Genome analysis of a wild papeda DYC002 and representative citrus species (A) Phylogenetic analysis of DYC002 and nine other representative citrus species using genome‐wide SNPs on the degenerated sites after the resequencing data were mapped to the assembled DYC002 monoploid genome. Pie charts associated with each species show the proportions of gene families that underwent expansion (red) or contraction (green). (B) GO enrichment analysis of the genes in rapidly expanded families in DYC002. (C) K _s value distribution of the paralogous gene pairs in the DYC002, P. trifoliata, C. sinensis, and C. grandis genomes. The black vertical line indicates the whole‐genome duplication event that occurred in the Citrus–Poncirus lineage (K _s = 1.25).

Figure 2

Comparative analysis of the two haplotype genomes of DYC002 (A) Large structural variations between the two haplotype genomes. (B) Venn diagram displaying the numbers of allele‐specifically expressed genes (ASEGs) identified in leaves, flower, small fruit, root, and stem, including 448 genes that consistently showed allele‐specific expression in all five tissues. (C) An apparent translocation event occurred between haplotype1‐chr3 and haplotype2‐chr6.

Figure 3

Evolutionary analyses of citrus species using genome‐wide single base variations (A) A phylogenetic tree of 378 representative citrus species using genome‐wide single nucleotide polymorphisms (SNPs) on the nuclear genome. Chinese box orange (Atalantia buxifolia (Poir.) Oliv.) was used as the outgroup. Colored names indicate the type of citrus to which each species belongs. The barplots on the outside circles represent the mean value of bitterness in the juice sacs (JS; blue) and segment membrane (SM; red) of each group of species. The color of the node showed the values of bootstrap percentage (red: 70 < bp ≤ 100, blue: BP ≤ 70). Scale bar represents the bitterness value. (B) Chronogram of citrus speciation. The reconstructed divergence time tree matches the phylogenetic tree in (A) and includes two distinct and well separated clades. The divergence times of P. trifoliata (9.5 Mya) and DYC002 (7.5 Mya) are used as constraints for the Bayesian estimation of species divergence times using the approximate likelihood method. (C) A working model for the origin and dispersal of citrus. All Citrus species originated in the Himalayas–Hengduan Mountains and were dispersed through two radiations. One radiation was directed toward Southeast Asia via human migration and rivers. Another radiation was toward the direction of Southeast China via rivers and human migration. The red star shows the location of the citrus fossil specimen Citrus linczangensis. The colored lines with an arrow showed the migration and evolution routes of different Citrus species. Papeda species originated in Northern Myanmar and Western Yunnan and was dispersed through the river basin of Dulong River, Nu River, Longchuan River, Lancang River and Red River. Pummelos could also be dispersed via those rivers. Citrons originated from the Himalayas–Hengduan Mountains and migrated through human activity. C. ichangensis originated in Yuanjiang in Yunnan Province and dispersed via the Zhujiang River, the Yangtze River and the Xijiang River. Mandarins originated in Nanling Mountains in Central‐south China and dispersed through Yangtze River, Zhujiang River and human migrations.

Figure 4

Demographic and migration analyses of citrus species (A) A TreeMix phylogeny with six migration edges that are indicated by colored lines. Each population contains between one and five individuals. We used Atalantia buxifolia (Poir.) Oliv. (AB001) as an outgroup. (B) Demographic history of six citrus lineages including the populations of Hybrid Kumquat, Primitive Mandarin, Rong'an Kumquat, Wild Kumquat, Wild Mandarin, and Zhejiang Kumquat inferred by SMC++.

Figure 5

Selective sweep analysis using the genome of representative accessions XP‐CLR scores and Fst are plotted across the nine citrus chromosomes. The Fst and XP‐CLR scores was calculated for each 10‐kb window. Fst values were shown on the top as gray bars. XP‐CLR scores were shown in gray points. The candidate regions for selective sweep are indicated by red points and bars, respectively. Gene IDs are shown for bitterness‐related genes near the selective sweep regions. (A) Comparison between the high bitterness samples of Ichang papeda (YCC1) and Yuanjiang Ichang papeda (YCC_YJ) versus low bitterness samples of Ichang papeda (YCC) and papeda (DYC). These two groups displayed different levels of bitterness in juice sacs. (B) Comparison between group one that includes high bitterness Poncirus (Z1) versus low bitterness Poncirus (Z) and kumquat (JG) that display different levels of bitterness in juice sacs. CCR2: Cinnamoyl‐CoA reductase; 4CL: 4‐coumarate‐CoA ligase‐like; ANR: Anthocyanidin reductase; UGTs: UDP‐glycosyltransferase family; CHI: Chalcone isomerase family.