Mitochondrial Genome Architecture and Evolutionary Origin of the Yao Silkworm, a Living Fossil of the Domestic Silkworm Bombyx mori (Lepidoptera: Bombycidae)

Abstract The Yao silkworm is a unique silkworm resource producing yellow flat plate silk that has only been reared by the Baiku Yao ethnic group in Nandan County, Guangxi Province, China for a thousand years. Here, we report the mitochondrial genomes (mitogenomes) of five Yao silkworm strains and 10 local Guangxi strains of the domestic silkworm (Bombyx mori) L. (Lepidoptera: Bombycidae), and use the resulting mitogenomes and the available Bombyx mitogenomes to characterize their genome architecture and trace the evolutionary origin of the Yao silkworm. The five Yao silkworm mitogenomes exhibited genome architectures identical to typical set of 37 mitochondrial genes (13 protein-coding genes, 22 transfer RNAs, and two ribosomal RNAs) and a high level of genome sequence similarity with the domestic silkworm. Mitogenome-based phylogenetic reconstruction provided solid evidence that the Yao silkworm shares a common ancestor with the domestic silkworm. Sliding window analysis uncovered a distinct variation pattern in the mitogenome between the Yao silkworm and the other domestic silkworm strains. The phylogenetic analyses revealed a basal placement of the Yao silkworm among all available domestic silkworm strains, indicating that the Yao silkworm is an ancient population of the domestic silkworm. Our data indicated that the Yao silkworm (B. mori) is a lineage of the domestic silkworm, which for the first time provides insights into the origin of the Yao silkworm.

The Baiku Yao ethnic group in Nandan County, Guangxi Province, China has been regarded as a 'living fossil of human civilization' by the United Nations Educational, Scientific and Cultural Organization (Tan andJiang 2012, Guan et al. 2014). Baiku Yao, a branch of the Yao ethnic group with a population size of approximately 30,000, is so named because the men wear white knee-length trousers and scarves, while the women wear pleated skirts and black scarves (Wen 2006). Since directly entering modern society from a primitive society, the Baiku Yao ethnic group still carries much cultural information on the transition between a matrilineal society and a patrilineal society. Because the Baiku Yao ethnic group seldom interacts with other nationalities, their unique customs and cultures such as clothing, intraethnic marriages, funerals, bronze drums, and dietary habits are still completely preserved (Yan et al. 2012).
The present study focused on the Yao silkworm, which has just been known by silkworm researchers since 1986. The Yao silkworm is reared only in 29 villages in Nandan County (N 24°25′ ~ 25°37′; E 107°01′ ~ 107°53′) (Fig. 1A) and has been reared to produce flat plate silk by the Baiku Yao ethnic group (Fig. 1B) for a thousand years; thus, it has been regarded as a 'living fossil of Chinese silk culture' (Lu et al. 2013). Similar to its ethnic culture, the Yao silkworm is a highly pure local population (Guan et al. 2014). The legend has it that the Yao silkworm was granted by the mountain deities Yaojia and Lanla from its wild counterpart, which suggests that the Yao silkworm might have its own origin history (Lu et al. 2013) that is different from the domestic silkworm (Bombyx mori) L. (Lepidoptera: Bombycidae), whose birthplace might be in the Yellow River basin (Jiang 1982, Luo 1993. Generally, the Yao silkworm eggs are incubated by women who carry them close to their skin to hasten their hatching. The larvae are reared at home for 40-50 d, and the mature larvae are placed on a board to produce flat plate silk that can then be used for clothing, weddings, and funerals (Tan and Jiang 2012). The Yao silkworm is a unique silkworm resource characterized by the production of univoltine, flat plate silk and natural golden yellow silk, and is a natural trimolter (Guan et al. 2014, Bai et al. 2015; however, genetic information about the Yao silkworm remains severely limited. Chinese has utilized domestic silkworms for over 5000 years (Goldsmith et al. 2005, Xiang et al. 2005. Historical and archaeological evidence has suggested that the domestication of domestic silkworms might have occurred independently at multiple centers, including the Yellow River basin, the Yangtse River basin, and ethnic minority areas (Jiang 1982). In the past decade, enormous efforts on understanding the origin of the domestic silkworm have been made at the molecular level and have confirmed that the wild silkworm Chinese Bombyx mandarina Moore (Lepidoptera: Bombycidae) is the wild ancestor of the domestic silkworm, instead of the Japanese B. mandarina (Lu et al. 2002, Li et al. 2005, Xiang et al. 2005, Arunkumar et al. 2006, Li et al. 2010. Genome resequencing of a large number of the domestic silkworm strains has also suggested that the domestic silkworm might be involved in a short domestication event (Xia et al. 2009) and that silkworms might have been initially domesticated as trimolting lines (Xiang et al. 2018  distributed throughout Asia and exhibit great diversity (Astaurov and Rovinskaya 1959) makes this issue more complex. A recent study of silkworms revealed that Chinese B. mandarina populations can be divided into two genetically distinctive subtypes, northern Chinese B. mandarina and southern Chinese B. mandarina, and that the former is the wild ancestor of the domestic silkworm (Chen et al. 2019). The Yao silkworm is a case in point (Luo et al. 1993, Lu et al. 2013); however, it was not included in these previous studies.
Does the Yao silkworm have its own origin history? The fact that no wild counterparts have been recorded in Guangxi opposes the theory of the independent domestication event of the Yao silkworm. In this context, we undertook the present work on the origin of the Yao silkworm. Our primary objective was to explore the evolutionary origin of the Yao silkworm through mitochondrial genome (mitogenome) comparison. The animal mitogenome encodes a typical set of 37 mitochondrial genes, including 13 protein-coding genes (PCGs), 22 transfer RNA (tRNA) genes, two ribosomal RNA (rRNA) genes, and a noncoding A+T-rich region. Due to maternal inheritance, small genome size, and relatively high evolutionary rate, the animal mitogenome has played an increasingly important role in phylogenetic studies (Cameron 2014). A recent study has suggested that the whole mitogenome will provide good access to the maternal origin of the domestic silkworm (Chen et al. 2019).
In the present study, complete mitogenomes from five strains of the Yao silkworm and 10 local Guangxi strains of the domestic silkworm were determined. By combining our new data with previously reported sequences, we obtained a dataset of 71 mitogenome sequences for all the available Bombyx data. Then, these sequences were compared to characterize their genome architecture, and mitogenome-based phylogenetic analyses were performed to trace the origin of the Yao silkworm.

Silkworm Collection
In the present study, five inbred strains of the Yao silkworm and 10 local Guangxi strains of the domestic silkworm were used (Table 1). All these strains were preserved normally at the Guangxi Institute of Sericulture Science, China or the Silkworm Gene Bank of Southwest University, China. The five Yao silkworm strains exhibiting different agronomic characteristics were initially collected from Nandan County, Guangxi, China.

Mitogenome Sequencing, Assembly and Annotation
A single pupa of each strain was used to extract the total genomic DNA using a Wizard Genomic DNA Purification Kit (Promega, Madison, WI). Standard Illumina sequencing was performed in the paired-end mode by Novogene Biotech Co., Ltd. (Beijing, China). We generated ~15× coverage for each sequenced sample. The raw short read data of 150 bp were then used to assemble the mitogenome, with strain C108 of the domestic silkworm and its annotation (AB070264; Yukuhiro et al. 2002) as the reference. Raw short reads were mapped onto the reference mitogenome by program SOAP v1.09 (Li et al. 2008). SOAPsnp ) was used to calculate the posterior probability of each possible genotype at every genome position from the alignment results. The consensus sequence was then structured based on the highest probability. New mitogenome sequences were annotated in MITOS (Bernt et al. 2013) and then refined manually by aligning against the reference sequence.

Data Collection and Sequence Comparison
We downloaded 38 other mitogenome sequences of the domestic silkworm available to date from the National Centre for Biotechnology Information. By combining our new data with previously reported sequences, we generated a dataset of 71 whole mitogenome sequences of the domestic silkworm and its wild relatives. Nucleotide sequences were aligned by using Clustal Omega with the default settings (https://www.ebi.ac.uk/Tools/msa/clustalo/) (Sievers et al. 2011) and revised by eye. DnaSP v6 (Rozas et al. 2017) was used for sliding window analysis to determine the nucleotide diversity of the mitogenomes of silkworms, with a step size of 200 bp and a window length of 800 bp.

Phylogenetic Reconstruction
The entire mitogenomes of the Yao silkworm and the domestic silkworm available were used to infer their phylogenetic relationships with maximum likelihood analysis and Bayesian inference algorithm. Two non-Bombyx species, Antheraea pernyi Guérin-Méneville (Lepidoptera: Saturniidae) (GenBank Accession no. AY242996), belonging to Saturniidae and Rondotia menciana Moore (Lepidoptera: Bombycidae) (Kong and Yang 2015), belonging to Bombycidae served as the outgroups. The GTR+G+I model was chosen according to the lowest Bayesian Information Criterion scores with MEGA X (Kumar et al. 2018). Maximum likelihood analysis was conducted with IQ-TREE 1.6.12 (Nguyen et al. 2015), and the bootstrap values were evaluated with 1,000 iterations. Bayesian inference analysis was performed using MrBayes v3.1.2 (Huelsenbeck and Ronquist 2001). The posterior distributions were estimated with Markov Chain Monte Carlo (MCMC) sampling. In this study, the MCMC search was conducted for 11,500,000 generations, with sampling every 100 generations, and the average standard deviation of split frequencies reached 0.0076. The first 25% of trees were then discarded as 'burn-in'. Two independent runs were performed for each method. The Interactive Tree Of Life (https://itol.embl.de) was used to visualize the obtained phylogenetic trees (Letunic and Bork 2021).

Characterization of the Mitogenome of the Yao Silkworm
In this study, the mitogenomes of five strains of the Yao silkworm and 10 local Guangxi strains of the domestic silkworm were first determined to characterize their genome architecture. These five Yao silkworm strains exhibited a low cocoon shell ratio of approximately 7.8%, which was comparable to the wild silkworm B. mandarina (7-8%), but smaller than that of the domestic silkworm B. mori (15%; Li et al. 2000). All 15 mitogenomes presented here were double-stranded circular molecules and contained 37 typical mitochondrial genes (13 PCGs, 22 tRNAs, and 2 rRNAs) and an A+Trich region, with identical genome architectures, gene orientation, and gene order. The mitogenome of strain Yao_2B2 was 15,657 bp in size, while that of the other four strains (2B1, 17W, 17B, and 2W) was 15,656 bp. The genome size of the Yao silkworm was well within the range detected in the sequenced B. mori strains but shorter than those found in B. mandarina. The composition of the five Yao silkworm mitogenomes was highly biased toward the A+T nucleotide (81.4%), which was in line with that of the sequenced B. mori strains. In the mitogenomes of the Yao silkworm, the 12 PCGs started with ATN (ATG/ATA/ATT), while COI started with CGA. Ten PCGs ended with TAA, while COI, COII, and NAD4L ended with the incomplete stop codons T or TA. Twenty-two typical tRNA genes were identified in each of the mitogenomes of the Yao silkworm, with sizes ranging between 65 bp and 73 bp. All tRNAs in each of the mitogenomes of the Yao silkworm, except for tRNA Ser (AGN), which lacked a stable dihydrouridine arm, showed the typical cloverleaf secondary structure observed in the mitochondrial tRNA genes of previously reported mitogenomes of silkworms. The A+T-rich region in strain Yao_2B2 was 495 bp, while those of the other four strains were 494 bp. The size of the A+T-rich region of the Yao silkworm was highly similar to those observed in the sequenced B. mori strains. A 126 bp repeat element has been well characterized in Bombyx species (Yukuhiro et al. 2002). Only one copy of this repeat element was found in each of the mitogenomes of the five Yao silkworm strains, as has been observed in 48 sequenced B. mori strains and 15 Chinese B. mandarina collections (Pan et al. 2008, Chen et al. 2019, whereas three copies of the tandem replication of this repeat element have been observed in the Japanese B. mandarina (Yukuhiro et al. 2002). The length differences in silkworm mitogenomes are primarily the result of the size variation in the A+T-rich region (Pan et al. 2008).
Among these five mitogenome sequences, we identified 13 single nucleotide variants (SNVs) across the 15,657 bp consensus aligned sequence (Fig. 2), thus indicating that they were very closely related, consistent with their recent selection history (Zhang et al. 2017). Twelve of 13 SNVs occurred in the coding region, and one occurred in the A+T-rich region. Only one variable site was found between Yao_2W and Yao_17B, while 11 variant sites occurred between Yao_2W and Yao_2B1. When compared to the first reference mitogenome sequence of B. mori strain C108 (AB070264), the number of variant sites between each pair ranged from 23 (Yao_17B and Yao_17W) to 25 (Yao_2B1 and Yao_2B2).

The Yao Silkworm Has Highly Sequence Similarity with the Domestic Silkworm, but Exhibits a Different Variation Pattern
Within the Bombyx silkworms, four distinct haplotype clades have been identified based on the mitogenome sequence, including the B. mori clade, Japanese B. mandarina clade, northern Chinese B. mandarina clade, and southern Chinese B. mandarina clade (Chen et al. 2019). The five mitogenome sequences of the Yao silkworm, together with 10 mitogenome sequences of local Guangxi strains of the domestic silkworm determined in this study, were aligned with 38 known mitogenome sequences of the domestic silkworm, and the result indicated a high similarity of the nucleotide sequences between the Yao silkworm and B. mori (data not shown). Analysis of the genetic distance also indicated that the Yao silkworm had the smallest genetic distance from strains of the B. mori clade, a value of 0.002 comparable to that within strains of the domestic silkworm (Table 2). Both sequence alignment and genetic distance showed that the Yao silkworm is a lineage of the B. mori clade rather than the B. mandarina clade.
Sliding window analysis was conducted to investigate the difference in variation patterns between the Yao silkworm and the domestic silkworm strains. The average value of nucleotide diversity (Pi) in the five Yao silkworm strains was 0.00042, while that of the 48 domestic silkworm strains was 0.00158. Through the sliding analysis, we found two regions with distinct variation patterns between the two groups (Fig. 3). One was the region between COIII and tRNA Phe , a fragment of approximately 2.3 kb. Within this region, the Pi values of 3 mutation hotspots for the 48 domestic silkworm strains were significantly higher (>0.003). The other region, a fragment of approximately 1.8 kb, was positioned between tRNA His and tRNA Thr : one mutation hotspot for the 48 domestic silkworm strains was found with a Pi value of 0.0028 (tRNA His -ND4).

The Yao Silkworm Shares a Common Ancestor with the Domestic Silkworm, but Possesses a Basal Position in the Phylogenetic Tree
To evaluate the phylogenetic relationship of the Yao silkworm and the other domestic silkworm strains, we performed a breakthrough phylogenetic analysis based on 71 whole mitogenome sequences of Bombyx silkworms, including five Yao silkworm sequences, 48 domestic silkworm sequences, three northern Chinese wild silkworm Fig. 2. Variable sites in the mitogenomes of five Yao silkworm strains (2B1, 2B2, 2W, 17B, and 17W) and the domestic silkworm strain (C108-AB070264). In total, 35 variable sites were identified among these six strains.
sequences, 13 southern Chinese wild silkworm sequences, and two Japanese wild silkworm sequences. The mitogenome-based phylogenetic tree provided us with a well-resolved phylogeny of Bombyx silkworms (Fig. 4). In the phylogenetic tree, the five Yao silkworm strains and 48 domestic silkworm strains clustered together with a great value support of 100%, suggesting that they share a common ancestor; that is, the Yao silkworm is a lineage of the domestic silkworm. The fact that the Yao silkworm and 48 domestic silkworm strains collected from all over the world formed a monophyletic group further provided molecular evidence for a single-domestication origin of the domestic silkworm. Accordingly, the phylogenetic tree implied that the legendary early independent domestication event might not have occurred in the Yao silkworm, which is in good agreement with the fact that no wild silkworms have been recorded in Guangxi, China. Moreover, phylogenetic analyses revealed the basal placement of the Yao silkworm among all available domestic silkworm strains.

Discussion
In this study, the mitogenomes of five Yao silkworm strains were characterized and then used to trace their evolutionary relationship with the domestic silkworm. The Yao silkworm is a unique silkworm resource only found in 29 villages in Nandan County located in South China that has been utilized to produce silk by the Baiku Yao ethnic group for a thousand years. Due to the high morphological similarity, the Yao silkworm has been suggested as a relatively pure local population of the domestic silkworm (Guan et al. 2014). On the origin of the Yao silkworm, no historical records are available, although the legend has it that the Yao silkworm was developed by the mountain deities (Tan andJiang 2012, Lu et al. 2013). Some researchers hold that the Yao silkworm might have its own origin history, which is different from the domestic silkworm, in which domestication might have occurred within the area adjacent to the Yellow River (Jiang 1982, Luo 1993. The fact that no wild counterparts have been obtained or recorded in Guangxi Province makes the independent domestication event of the Yao silkworm impossible. Our analysis of the mitogenomes of five Yao silkworm strains revealed a great similarity with the domestic silkworm strains in genome size, gene composition, and genomic organization and a large difference with the wild silkworm B. mandarina. Our phylogenetic analysis based on the mitogenome further clustered the Yao silkworm together with the domestic silkworm strains. This is the first study to show molecular evidence that the Yao silkworm shares a common ancestor with the domestic silkworm strains. Sliding window analysis uncovered a different variation pattern in the mitogenome between the Yao silkworm and the other domestic silkworm strains. Further analysis identified two regions with distinct variation patterns between the two groups: one was located between COIII and tRNA Phe , and the other was located between tRNA His and tRNA Thr . The distinct variation pattern between the two groups may be derived from their different selection pressures or genetic drift. The domestic silkworm was successfully domesticated from the wild silkworm through artificial selection for thousands of years, with improvements in the cocoon shell ratio from 7-8% in the wild silkworm to 15% in the domestic silkworm (Li et al. 2000). Our previous investigation indicated that the Yao silkworm exhibits a low cocoon shell ratio (7.8%), which is in line with the fact that the Yao silkworm is reared in a primitive state (Zhang et al. 2017). The findings from molecular and economic trait analyses suggest that the Yao silkworm is a unique silkworm resource that should be given more attention for its perseveration. In the near future, comparative analysis of the nucleic genomes between the Yao silkworm and the domestic silkworm will provide novel insights into the genetic improvement of the domestic silkworm.
Our mitogenome-based phylogenetic analysis revealed the basal placement of the Yao silkworm within all available domestic  Fig. 3. Sliding window analysis of the mitogenomes of the Yao silkworm (n = 5) and the domestic silkworm (n = 48). Window length: 800 bp; step size: 200 bp. X-axis: position of the middle point of the window. Y-axis: nucleotide diversity per window. Region A is located between COIII and tRNA Phe , and Region B is positioned between tRNA His and tRNA Thr . silkworm strains. The basal placement of Yao silkworm that exhibits natural trimolter further supports the view that the domestic silkworm might have been initially domesticated as trimoulting lines according to genome resequencing (Li et al. 2010). The findings presented here suggest that the Yao silkworm is a relatively primitive lineage of the domestic silkworm that can be regarded as the living fossil of the domestic silkworm.
In summary, our study presented 15 complete mitogenomes from five Yao silkworm strains and 10 local Guangxi strains of the domestic silkworm. By comparative analysis of 71 complete mitogenomes including our new data and previously reported sequences of Bombyx, we found that the Yao silkworm is a lineage of the domestic silkworm despite reproductive isolation for a long period of time between them. We also uncovered a distinct variation pattern in the mitogenome between the Yao silkworm and the other domestic silkworm strains. The basal placement of the Yao silkworm in the phylogenetic tree among all available domestic silkworm strains indicated that the Yao silkworm is an ancient population of the domestic silkworm (B. mori) that can be regarded as the living fossil of the domestic silkworm.