Large-Scale Screening of Growth-Related Variants in Chinese Tongue Sole (Cynoglossus semilaevis)

SONG Weihao, ZHU He, WANG Yujue, ZHANG Kai, ZHANG Quanqi, and HE Yan

Large-Scale Screening of Growth-Related Variants in Chinese Tongue Sole ()

SONG Weihao, ZHU He, WANG Yujue, ZHANG Kai, ZHANG Quanqi, and HE Yan*

,,,266003,

Growth is the most valuable economic trait for improving aquaculture fish species, since fast growth can reduce labor cost and make more economic benefits. However, the knowledge about how many and which genes are related to the growth of Chi- nese tongue sole is limited. High-throughput sequencing screening of variants is a fast, economical and accurate assay to identify genes related to growth in crops, livestock and other aquaculture products. In the present study, genome-wide resequencing of 30 female Chinese tongue sole individuals from large and small groups to screen growth trait-related variations (SNPs, InDels) was carried out. In total, 6545735 SNPs and 1016745 InDels were detected, while 31 genes related to growth traits were identified. Their functions were mainly involved in muscle architecture, post-embryonic development, neurosensory development and hormone regulation. Fur- thermore, 18 of the 31 genes related to growth trait are located on W chromosome, indicating that W chromosome plays an important role in regulating the body size of female tongue soles. The markers and genes identified in our study can be applied to selective breeding of tongue sole and thus will promote the aquaculture industry and increase economic efficiency.

; growth; resequencing; SNP; indel

1 Introduction

Growth-related traits are of great importance for aqua- tic species as they influence the production directly. As a breeding target, growth-related traits are controlled by many genes and multiple pathways. Some previous studies have identified some candidate genes related to growth, reveal- ing the genetic basis of growth in aquatic species.() gene was one of the earliest genes in growth, development, reproduction, immune function, food conver- sion and appetite of fish (Quik, 2010).gene is an ideal candidate for fish growth trait selection, because many growth-relatedgenes’ single nucleotide polymor- phisms (SNPs) have been reported in fish species, such as Atlantic salmon (Gross and Nilsson, 1999), Chinook sal- mon and yellow croaker (Docker and Heath, 2002; Ni, 2012). Additionally, GnRH peptides stimulate GH release in goldfish bothand, and in common carp(Lin, 1993; Lin, 1995). In the process of muscle production and regulation, myogenic regulatory factors activate the transcription of genes involved in ske- letal muscle production through multiple pathways, there- by promoting the transformation of quiescent muscle satellite cells and other types of cells into myoblasts and pro- moting the differentiation of myoblasts into mature cells (Sabourin and Rudnicki, 2000). The paired box protein 3 and 7 () and severalgenes also contribute to thecontrol of muscle differentiation (Rescan and Ralliere, 2010).

Variations in DNA sequences are known to be the basis for trait variation. For marker-assisted selection (MAS), both genetic markers and linked genes can be used for se- lecting fish with specific characteristics. In the initial study, growth traits were mainly localized to the somatotropic axis genes, the myogenic regulatory factor genes, and the transforming growth factor genes (Valente, 2013). Some studies suggest that,(insulin-like growth factor), and() are candidate genes for mar- ker-assisted selection (De-Santis and Jerry, 2007). Single- nucleotide polymorphisms (SNPs) (Sánchez-Ramos, 2012), simple sequence repeats (SSRs) (Almuly, 2005), variable number of tandem repeats (VNTRs) and restriction fragment length polymorphisms (RFLPs) have been identified in large numbers in fish growth-related traits (Al- muly, 2005; Almuly, 2008; Sánchez-Ramos, 2012). However, identifying growth-related genes one by one is inefficient. With the advent of high-throughput se- quencing technology, it becomes accessible for exploring and genotyping a large number of variants in a genome- scale to locate quantitative trait locus (QTLs) and trait-re- lated SNPs. A QTL for growth in linkage group 1 was de- tected in European sea bass, after which two different QTLs were identified for body weight in the same species (Chatziplis, 2007; Massault, 2010). In Asian seabass, 43 QTLs were detected for growth traits during the development and two QTLs intervals at two stages were overlapped, while the others were mapped onto different positions (Xia, 2013). Twenty-two SNP markers and a mitochondrial haplotype were found to be significantly associated with growth traits in rainbow trout with the full transcriptome analysis of RNA sequences (Salem,2012). Moreover, 435 growth traits-related SNPs were iden- tified through transcriptome SNP analysis between fast- and slow-growing turbot (Robledo, 2017). Several miRNAs have been known to regulate muscle growth and development.For example, 168 small RNAs are express- ed in the muscles of fast myoblasts, and the expression re- gulation of several small RNAs is related to the transition from hyperplasia to hypertrophy during development (Johns- ton., 2009).

Tongue sole has become a very important maricultural fish in China. At present, although families with fast growth and high female ratio have been screened out, a breeding generation of Chinese tongue sole needs 2–3 years, and se- lection breeding faces huge time cost. Using 1007 SSR mar- kers, four weight-and body-width-related quantitative trait loci (QTLs) were detected in Chinese tongue sole (Song, 2012).andgenes of Chines tongue sole were cloned (Ma, 2011, 2012). Furthermore, the ex- pression levels of these two genes in females were signi- ficantly higher than those in males. There was no signifi- cant difference in GHRH expression between males and females before the age of 8 months, but GHRH expres- sion in females was significantly higher than that in males during the age of 9–12 months (Ji, 2011). These stu- dies are all at the one-gene level to study the characteris- tics of tongue sole, which is of great significance to the study of the growth of tongue sole and other fish species. However, functional studies related to growth traits at the genome-wide level has not been carried out.

To gain knowledge on the key genes related to growth in Chinese tongue sole, genome-wide resequencing studies were conducted in large and small groups in this study. Two methods, Fisher’s exact test combined with FSTstatistics were used to identify genetic markers associated with growth. GO and KEGG enrichment analysis of genes were further carried out for the candidate genes. The markers and genes identified in our study can be applied to selective breed- ing of tongue sole and thus will promote the aquaculture industry and increase economic efficiency.

2 Materials and Methods

2.1 DNA Extraction

Extreme-sized individuals of Chinese tongue soles were selectedfrom mixed cultures. Muscle samples were taken from female individuals (both physiological and genetic sex). Body length, width, and weight were measured. To- tal DNA was extracted using the Phenol-chloroform me- thod following the protocol of our laboratory and dissolved with TE to preserve genomic DNA. The quality and the concentration of extracted DNA were detected using agarose gel electrophoresis and Nanodrop 2000 spec- trophotometers.Then DNA sampleswere frozen at −20℃.

2.2 WGS and Quality Control (QC)

The original image data generated by the sequencing ma- chine were converted into sequence database calling (Illumina pipeline CASAVA v1.8.2) and then subjected to quality control (QC) procedure to remove unusable reads:

1) The reads contain the Illumina library construction adapters;

2) The reads contain more than 10% unknown bases (N bases);

3) One end of the read contains more than 50% of low- quality bases (sequencing quality value≤5).

2.3 Variant Detection and Annotation

Sequencing reads were aligned to the reference genome using BWA-0.7.17 with default parameters (Li and Dur- bin, 2009). Duplicate removal was performed using samtools and PICARD (http://picard.sourceforge.net).

The raw SNP/InDel sets are called by samtools-1.9 with the parameters as ‘-q 1 -C 50 -m 2 -F 0.002 -d 1000’ (Li., 2009). Then we filtered these sets using the follow- ing criteria:

1) The mapping quality>20;

2) The depth of the variate position>4.

SnpEff v.4.266 was used for the functional annotation of variants with thereference geno- me (Cingolani, 2012; Chen, 2014).

The TBtools software was used to plot the statistical re- sults of the variants (Chen, 2018).

2.4 Statistical Difference Analysis

Differential allele frequency and genotype frequency be-tween large group (LG) and small group (SG) samples were obtained by script following the formula ‘allele/genotype frequency=(the count of a type of allele/genotype in one locus)/(total count of allele/genotype in one locus)’. Thesig- nificance of the allele frequency and genotype frequency be- tween the two groups was tested by Fisher’s exact test. ThoseSNPs showing-values<0.001 and further MAF>0.10 were considered significant. GO term and KEGG pathway enrich- ment analyses of the annotated genes were performed for theselected genes using the OmicShare tools (www.omicshare. com/tools), a free online platform for data analysis.

2.5 Selective Scanning by F Statistics

An FSTvalue was calculated for each window to compare the difference between any two pools. Before the cal- culation, the selected windows and defining sweep windows were identified. VCFtools v0.1.13 was applied to cal- culated Weir & Cockerham’s FSTvalues on a 10kb window size and with 5kb step size between windows (Danecek, 2011).

2.6 Kompetitive Allele Specific PCR(KASP) Genotyping

KASP high-throughput genotyping is a SNPline genotyping assay based on KASP technology. KASP is a com- petitive allele-specific PCR that performs genotyping of SNPs by performing precise biallelic discrimination in a wide range of genomic DNA samples based on the specific matching of primer-terminal bases.

1) Design 2 SNP PCR primers, each corresponding to an SNP allele. The 3’ ends are respectively two alleles, and the 5’ ends are respectively FAM/HEX tag sequences.

2) Fluorescent probes were designed to be identical to the tag sequence, with a FAM/HEX fluorophore at each of the 5’ end. Besides, a quenching probe with a quenching group at the 3’ end was designed.

3) Perform specific PCR amplification and perform SNP genotyping by detecting fluorescence.

2.7 Ethics Statement

The study was conducted following the guidelines and regulations established by the Chinese Government Principles for the Utilization and Care of Animals Used in Tes- ting, Research, and Training. All the experimental proto-cols were permitted and approved by College of Marine Life Science, Ocean University of China (Qingdao, Chi- na).

3 Results

3.1 DNA Resequencing Output

We conducted a resequencing analysis of two groups of Chinese tongue sole among 15 large individuals and 15 small individuals (Fig.1). The exclusion of pseudo males was carried out following a previously described method (Liu, 2014). A total of 30 libraries were constructed using total DNAs from muscle tissue. Each library was sequenced with an Illumina Hiseq™ 4000. As a result, a total of 229.86G raw data were generated, for an average of 7.66G per sample, ranging from 5.43G to 12.70G. After filtering, a total of 228.78G clean data was retained, for an average of 7.63G per sample, ranging from 5.39G to 12.67G. The average effective rate is 99.50% (ranging from 98.64% to 99.77%). The high quality of resequenc- ing (Q20≥94.13%, Q30≥87.86%) and the normal distribution of GC content (40.47%–43.69%) means the next step mapping is effective (Table 1). On average, 96.44% of the clean reads align to the Chinese tongue sole geno- me, ranging from 94.12% to 97.06%. The average depth is 14.07X, ranging from 9.93X to 22.57X (Table 2). The mapping results are qualified, which can be used in subsequent mutation detection and related analysis.

Fig.1 Statistics of two groups of Chinese tongue sole samples for resequencing. A, Body weight of the LG and SG groups; B, Body length and width of the LG and SG groups.

Table 1 The statistics of re-sequencing data

Table 2 The mapping rate of resequencing data to the genome

3.2 Variants Calling in Two Groups of Chinese Tongue Sole

A total of 6545735 SNPs and 1016745 InDels were call- ed based on resequencing data. The sequences containing these SNPs were annotated according to the Chinese ton- gue sole genome information. Because several effects can be annotated to one locus, a total of 34132325 effected SNPs and 5506066 effected InDels were detected at the level of the entire genome (including repeated annotation). A total of 11846123 effected SNPs (34.71%) and 1796516effected InDels (32.63%) were found in exonic regions.Among these SNPs, 244456 correspond to non-synony- mous variants and 551545 correspond to synonymous va- riants. There are 34158680 transitions and 21419741 trans- versions among all SNPs, and the Ts/Tv ratio is 1.59.

3.3 Association Analysis

Allele frequencies and genotype frequencies were estimated for LG (large group) and SG (small group) of Chi- nese tongue sole, and the significance of each frequency was tested. Significant difference in allele frequencies be- tween LG and SG were detected in 1206 SNPs and 428 Indels,while significant difference in genotype frequencies were found in 830 SNPs and 149 InDels. These 2151 significant variants were located in 1324 genes.Among these variants, 414 significant variants were located in up- stream of genes; 170 significant variants were located in exonic regions, which include 58 synonymous variants and 23 non-synonymous variants; 1077 significant variants were located in intron regions; 1027 significant variants were located in non-genetic regions. (A particular type of variant of one gene may be another type of variant of another gene).

Fig.2 SNP/Indel calling and annotating results. a, Chromosome scale; b, SNP density, the light to dark indicates low to high density; c, Indel density, the light to dark indicates low to high density; d, FST between LG and SG groups; e, Fisher P-value of SNPs; f, Fisher P-value of Indels; g, Chromosome numbers. In e and f, red color stands for allele frequencies; green and blue represent frequencies of different genotypes.

3.4 The GO and KEGG Enrichment of Growth-Related Markers

To better understand the function of these significantly divergent genes between LG and SG groups of Chinese tongue sole, GO term enrichment was conducted. As the results, 669 GO terms with<0.05 were significantly en- riched in whole divergent genes, and 28 GO terms with FDR<0.05 were significantly enriched (Fig.3). Among the 28 most significantly enriched GO terms, there were 17 GO categories within the ‘Biological process’ and the three most abundant terms were ‘cell differentiation’, ‘single-or- ganism cellular process’ and ‘cellular developmental pro- cess’. There were also a large number of genes being involved in ‘single-organism process’, ‘response to stimu- lus’, ‘developmental process’, ‘single-organism develop- mental process’ and ‘single-multicellular organism process’. In the ‘Cellular component’ division, there were 11 categories and the top three categories were ‘cell junction’, ‘synapse’, and ‘membrane’. Though there were no most significantly enriched terms in the ‘Molecular function’ division, the top three categories were ‘protein tyrosine kinase activity’, ‘metalloaminopeptidase activity’ and ‘transmembrane receptor protein tyrosine kinase activity’.

Subsequently, KEGG pathway enrichment revealed that 16 pathways were significantly enriched with<0.05 (Fig.4). There were 7 pathways (N-Glycan biosynthesis, other types of O-glycan biosynthesis, glycosphingolipid biosynthesis-lacto and neolacto series, fatty acid biosynthesis, citrate cycle, biotin metabolism, and glycosphingo-lipid biosynthesis-ganglio series) significantly enriched in metabolism division. In cellular community division, focal adhesion was significantly enriched. In signal transduction division, the Rap1 signaling pathway and Ras sig- naling pathway were significantly enriched. Among the pathways at the top of the non-significant ranking, salivary secretion, pancreatic secretion, steroid hormone bio- synthesis, and GnRH signaling pathway were enriched. There were also a large number of genes being enriched in the PI3K-Akt signaling pathway and the actin cytoske- letonregulation.

Fig.3 GO enrichment of candidate genes. Enriched GO terms of selected marker genes (P<0.05).

Fig.4 KEGG pathway enrichment of candidate genes. A, Count of enriched KEGG terms of selected marker genes of KEGG A and B class; B, Enriched KEGG pathways of selected marker genes (P<0.05). The size bar indicates enriched gene numbers of each KEGG term. The color bar indicates P value from low (green) to high (red).

3.5 Selective Scanning by F Statistics

Genome scanning by comparing the single locus estimates of FSTbetween different populations might identify regions of the genome that have been subjected to diversifying selection (Holsinger and Weir, 2009). Among all 91551 windows (10kb in length, sliding in 5kb steps across the Chinese tongue sole genome), 87939 windows contained >10 variants while covering 96.05% of the genome and 99.77% of variants (Fig.5), which were used to detect signatures of selective sweeps. 4397 windows were selected, and all of them have significantly higher Weir & Cockerham’s FSTvalues (5% right tail). Consequently, we identified a total of 42.93Mb genomic regions (4.80% of the genome) with strong selective sweep signals between LG and SG groups (Fig.2).

Fig.5 SNP distribution in windows.

Combining the significant different frequencies, we de- tected 571 variants. Among these variants, 130 significant variants were located in upstream of genes, 46 significant variants were located in exonic regions, which include 10 synonymous variants and 7 non-synonymous variants, 275 significant variants were located in intron regions, 347 sig- nificant variants were located in non-genetic regions. We detected that non-synonymous variants were located at the coding region of Pro-interleukin-16 (), Ephrin type- A receptor6 (), Interleukin-8(), Alpha-1A ad- renergic receptor(), tRNA (cytosine(34)-C(5))-me thyltransferase (, mitochondrial) and Uncharacte- rized protein C11orf65 homolog. And we also detected sy- nonymous variants in the coding region of Acetoacetyl- CoA synthetase (), Chromodomain-helicase-DNA- binding protein 9 (), Ellis-van Creveld syndrome pro- tein homolog (), F-box/WD repeat-containing protein 8 (), Ankyrin repeat and KH domain-containing protein 1 (), E3 ubiquitin-protein ligase MYCBP2(), tRNA (cytosine(34)-C(5))-methyltransferase (, mitochondrial), PR domain zinc finger protein 8 () and Ral guanine nucleotide dissociation stimulator ().

In addition to those variants located in coding regions, we also scan the untranslated regions (3’, 5’ UTR and up- stream), and growth-related genes such as Bone morphogenetic protein 7 (), Vascular endothelial growth fac- tor receptor 1 (), Desmin (), Muscleblind-like pro- tein 1 (), Solute carrier family 2 and facilitated glu- cose transporter member 8 () were detected.

3.6 Resequencing Results Verified by KASP Genotyping

Among the candidate growth-related genes, we random-ly selected 9 extremely significant loci and performed KASP genotyping on the sequenced samples to verify the reliability of the resequencing results. The genotyping re- sults of 8 loci were consistent with the resequencing re- sults, and there was a significant difference between large and small groups (Table 3), which proved that the rese- quencing results were reliable.

4 Discussion

In present study, we re-sequenced two groups of female Chinese tongue sole, a large group and a small group,from the same population with different growth rates, and more than two thousand variants were detected. Most of these variants were related to different aspects of growth. A large number of densely distributed divergent variants were clustered at autosome 11 and sex chromosome W. In- terestingly, several selected genes are evenly distributed on the sex chromosome W. Quantitative traits located by genome-wide resequencing has been widely used in crop and animal husbandry breeding. And it is very effective and greatly shortens the breeding cycle (Yano, 2016).

Growth is coordinated and regulated by different aspect of body systems,including muscle architecture, muscle regeneration and muscle fusion, embryo development, neu- rological regulation and endocrine regulation, energy metabolism, protein synthesis and degradation. In our work, we detected genes from several aspects related to Chinese tongue sole growth. For instance, there are lots of genes involving in different kinds of metabolism pathways, digestive, excretory, endocrine and immune systems, neuro- genesis and sensory formation, regulation of actin cyto- skeleton and steroid hormone biosynthesis.

4.1 Muscle Architecture

Unlike other vertebrates, fish’s muscles are separated by discrete layers with different fiber types. Vertebrates have two main types of striated muscle fibers: red and white, which specialize in low-speed swimming and burstsof maximum in fish, respectively (Sanger and Stoiber, 2001). The fat content of fish is generally low, though the speed of muscle growth and the size of muscle fibers are the main indicators to measure the growth of fish. We de- tected six myosin genes and one intermediate filament genewas detected in the highly differentiated region by F-statistics. Desmin() is a muscle-specific type III intermediate filament, which is essential to the structureand function of muscle. It plays an important role in main- taining sarcomere structure, connecting Z disc and form- ing myofibrils. These genes are not only connected with the cytoskeleton of sarcomere cells, but also with the nu- cleus and mitochondria, thus providing strength for mu- scle fibers during activity (Hnia, 2015).

Table 3 The comparison of P-values between resequencing and KASP genotyping

Notes: Fre., frequency; *<0.05, **<0.01.

4.2 Post-Embryonic Growth

Teleosts often exhibit an indeterminate growth pattern, with body size and muscle mass increasing until mortality or senescence occurs (Johnston, 2011). In addition, the body size increases dramatically between embryo and adult in most teleosts. Therefore, there must be certain me- chanisms that control the rapid development of the body and the sustained growth of muscles. Since muscle cells are multinucleated, terminal differentiated tissue, post- embryonic growth requires a source of proliferative MPC. The adaxial cells of zebrafish differentiate into slow mus- cle fibers of adult fish, while the lateral presomitic cells remain deep in the muscle and differentiate into fast muscle fibers (Devoto, 1996).

Interestingly, we found that a series of selected genes were associated with cellular development, including Chro- modomain-helicase-DNA-binding protein 9 (), Ellis- van Creveld syndrome protein homolog () and PR do- main zinc finger protein 8 (). Each of them had a synonymous SNP located in FSTselective regions.participates in the differentiation of progenitor cells du- ring osteogenesis and participates in the relaxed chroma- tin structure of mouse oocyte, which is a prerequisite for post-fertilization versatility (Shur, 2006; Ooga, 2018).is a component of EvC complex that regulates ciliary Hedgehog (Hh) signaling, and participates in endo- chondral growth and skeletal development (Ruiz-Perez, 2007).is a possible histone methyltransfe- rase that preferentially acts on the ninth Lys of Histone H3 (Eom, 2009). Steroid-producing marker genes such asandare involved in the control of steroid production through transcriptional inhibition (Ross, 2012). It forms a transcriptional inhibitory complex with bhlhe22, which controls genes involved in nerve deve- lopment and differentiation. In the retina, rod bipolar and type 2 non-conical bipolar cells are required to survive (Jung, 2015). Beyond these,bone morphogenetic pro- tein 7() can induce cartilage and bone formation, which is a bone-inducing factor leading to epithelial os- teogenesis and plays a role in calcium regulation and bonebalance. In addition,plays an important role in the de- velopment of mammalian skeleton, kidney, eyes and brown adipose tissue (Dudley, 1995; Katagiri, 1998; Boon, 2013). During the development of zebrafish,is strongly expressed in the eyes, ears, pronephros and gastrointestinal system (Shawi and Serluca, 2008). Prostaglandin E synthase 2 () catalyzes the conver- sion of PGH2 into prostaglandin E2 () which is more stable. On one hand, PGE2 production and transport are essential for the formation of cilia and photoreceptor de- velopment in zebrafish (Li, 2019). On the other hand, PGE2 signal participates in medaka ovulation follicle rupture and results in the cytoskeleton rearrangement of GC actin during ovulation (Takahashi, 2018). There is a basement membrane on the outer surface of Schwann cells, which plays an important role in peripheral nerve re- generation. In zebrafish embryos, the lack of, a mes- senger of gene output, can damage the development of Schwann cells (Seytanoglu, 2016). Nucleoporinis necessary for the nucleus to export poly (a) tail-con- taining RNA into the cytoplasm. It may be involved in the terminal step of the nuclear pore complex (NPC) in theprocess of RNA transfer (Bolger, 2008), and it is also one of the genes we screened.

In addition to a large number of genes selected for cel- lular development, we have also screened a series of genes for cell proliferation and differentiation, which may be re- lated to post-embryonic growth. In interleukin-8 (), two non-synonymous SNPs were detected and were consistent with the FSTselection regions.is an autocrine factor that promotes the growth of fibroblasts. In humans,plays a role in growth and maintenance of ectopic endo- metrial tissues by chemically attracting and stimulating leucocyte to secrete growth factors and cytokines, and by directly affecting the proliferation of endometrial cells (Ari- ci, 2002). In alpha-1a adrenergic receptor (), one non-synonymous SNP in FSTselection regions was detect- ed.is a member of the G protein-coupled receptor superfamily. They activate mitotic responses and regulate the growth and proliferation of many cells (Luttrell, 2008). In ral guanine nucleotide dissociation stimulator (), three upstream SNPs, two 3’UTR variants and one syn- onymous SNP were detected.stimulates the sepa- ration of ras-related ral-a and Ral-b GTPases from GDP, thus allowing GTP to bind and activate GTPases. ral-a is a multifunctional GTPase involved in many cell processes, including gene expression, cell migration, cell prolifera- tion, carcinogenic transformation and membrane transport. Its multiple functions are realized by interacting with dif- ferent downstream effectors. At the late stage of cell divi- sion, after the formation of bridges between dividing cells is completed, ral-b mediates the recruitment of external cysts to the midbody to drive separation (Cascone, 2008). Vascular endothelial growth factor receptor 1 () is a tyrosine protein kinase.It is the cell surface receptor of vascular endothelial growth factor A, vascular endothe- lial growth factor B and prostaglandin F, and plays an im- portant role in the development of embryonic vascular system, angiogenesis, cell survival, cell migration, macro- phage function, chemotaxis and regulation of cancer cell invasion. It may also play an important role in negative regulation of embryonic angiogenesis by inhibiting endo- thelial cell proliferation (Kendall and Thomas, 1993). Fur- thermore,can promote the expression of GFP in zebra- fish arteries and veins as a special enhancer (Bussmann, 2010).

4.3 Neurological and Sensory Development

The development of nerves and senses affects the individual’s speed of movement, the ability to find food, the ability to avoid harm, and the regulation of digestion and absorption, so it is inseparable from the individual’s growth. A transcription factor for brain development, a vesicle transporter and a secreted glycoprotein involved in regulating chemosensory cilia in olfactory neurons and optic nerve extension were screened by differential expression analysis in turbot (Robledo, 2017). A non-synony- mous SNP of ephrin type-a receptor 6 () was de- tected. It is a tyrosine kinase of EPH receptor. Its mem- brane-binding ligand ephrins plays a key role in pattern formation and morphogenesis. An important role of EPH receptor and renin is to mediate contact-dependent rejection, which is related to the pathological findings of axons and crest cells and the restriction of cell mixing between the posterior brain segments. Biochemical studies have shown that the degree of EPH receptor polymerization can regulate cell response, and actin cytoskeleton is the main target of EPH receptor-activated intracellular path- way (Wilkinson, 2000). In mice and zebrafish,was found to be associated with retinal axon guidance, which is one of the eye development receptors, and inhibited by dimer HMX1 (Marcelli, 2014). Lipoxygenase ho- mology domain-containing protein 1 () is required for normal function of inner ear hair cells (Grillet, 2009). Acetoacetyl-CoA synthetase () was detected as a synonymous SNP in selected regions.It is a ketone body- utilizing enzyme for the synthesis of cholesterol and fatty acids and is highly expressed in the brain.

4.4 Hormonal Regulation

Pituitary growth hormone () is generally considered to be essential for postnatal growth of mammals and is one of the most widely studied genes in the field of fish growth. In fishes, growth hormone induces muscle growth by regulating the expression of several genes of(),,and insulin-like growthfactor system, and myogenic regulatory factors (MRFs). IGFs stimulates myogenic cell proliferation, differentia- tion and protein synthesis through MAPK/ERK and PI3K/ Akt/Tor signaling pathways and eliminates protein degra- dation and atrophy through PI3K/Akt/FoxO signaling pa- thway (Fuentes, 2013). On one hand, molecular va- riation markers related to growth hormone receptor gene were detected in several animals between fast-growing spe- cies and slow-growing species (Zhang, 2017). On the other hand, growth hormone transgenic strains may have higher growth rate (Chen, 2015). In our research, we detected a few genes associated with some other growth- related hormones and their receptors,such as solute car- rier family 2,(),which is an insulin-regulated facilitative glucose transporter. Bindingin a glucose-inhibi- table manner (Ibberson, 2000), protein FAM102A () plays a role in estrogen action (Wang, 2004).(),(), parathyroid hormone-related protein (),() and acidic fibroblast growth factor intracellular-binding protein () were all detected in di-vergent variants analysis, but not in selective scanning ana- lysis.

4.5 Protein Degradation

Proteins in the body are constantly decomposed and re- placed. In many teleosts, skeletal muscle also undergoes a process of accelerating protein decomposition during the seasonal stages of fasting and gonadal maturation.When protein degradation far exceeds protein synthesis, atrophy occurs. Protein degradation is complex, including ubiqui- tin-proteasome system, calpain, NF-KB pathway and ly- sosome. Damaged proteins and proteins with short half- lives are modified by ubiquitination, resulting in their de- gradation by a multi-catalytic protease complex called pro- teasome (Johnston, 2011). We detected two genes related to protein degradation. One is tumor necrosis factor receptor type 1-associated DEATH domain protein(), which contains binding molecules interacting with TNFRSF1a/TNFR1 and mediating programmed cell death signals and activation of NF-kappab. The other is the F- BOX/WD repeat containing protein 8 () that mediates the degradation of MAP4K1/hpk1, through which affecting cell proliferation and differentiation (Wang, 2014).

4.6 Sexual Dimorphism

The genetic sex of Chinese tongue sole is determined by W and Z chromosomes and exhibits female-biased se- xual size dimorphism. 64.35% of the module genes were found on the W chromosome in the WGCNA analysis of the body size sexual dimorphism (Chen, 2014). In the present study, 58.06% (18/31) genes associated with growth were detected on the W chromosome. Present and previous studies have shown that the W chromosome of Chinese tongue sole not only has genes affecting the male and female body size sexual dimorphism, but also has genes controling the size of female tongue sole.The results were consistent with the results of the body size se- xual dimorphism research that a large number of genes that led to male and female body size dimorphism bias were located on the W chromosome (Weir and Cockerham, 1984).

Although we have screened a large number of potential growth-related loci, further validation work needs to be done. This is probably the result of polygenic interaction, and its mechanism needs further research. This study pro- vides a large number of potential research options for the study of genes related to growth traits in Chinese tongue sole, whereas further validation work needs to be done to investigate the gene functions.

5 Conclusions

Instead of random sampling, we selected individuals with extreme weight in population for sampling and used sta- tistical difference analysis to screen out a large number of growth-related mutation sites. Assuming that extreme groups come from different populations, FSTmethod was used to further screen regions with significant differences between two groups in the genome.The results of statistical difference analysis were intersected to further locate genes and loci significantly related to growth traits. It was found that there were a large number of sex-independent but growth-related genes on the W chromosome of the Chi- nese tongue sole. It was the presence of these genes that led to the larger body size of female (ZW) than male (ZZ) at the same age.

Acknowledgement

This work was supported by the National Natural Sci- ence Foundation of China (No. 31402292).

Almuly, R., Poleg-Danin, Y., Gorshkov, S., Gorshkova, G., Rapo- port, B., Soller, M., Kashi, Y., and Funkenstein, B., 2005. Cha- racterization of the 5’ flanking region of the growth hormone gene of the marine teleost, gilthead sea bream: analysis of a polymorphic microsatellite in the proximal pro- moter., 71 (3): 479-490, DOI: 10.1111/j.1444-2906.2005.00991.x.

Almuly, R., Skopal, T., and Funkenstein, B., 2008. Regulatory re- gions in the promoter and first intron ofgrowth hormone gene: Repression of gene activity by a polymorphic minisatellite., 3 (1): 43-50, DOI: 10.1016/j.cbd.2006.12.006.

Arici, A., 2002. Local cytokines in endometrial tissue: The role of interleukin-8 in the pathogenesis of endometriosis.,955 (1): 101-109, DOI: 10.1111/j.1749-6632.2002.tb02770.x.

Bolger, T.A., Folkmann, A.W., Tran, E.J., and Wente, S.R., 2008. The mRNA export factor Gle1 and inositol hexakisphosphate regulate distinct stages of translation.,134 (4): 624-633, DOI: 10.1016/j.cell.2008.06.027.

Boon, M.R., van den Berg, S.A., Wang, Y., van den Bossche, J., Karkampouna, S., Bauwens, M., De Saint-Hubert, M., van der Horst, G., Vukicevic, S., and de Winther, M.P., 2013. BMP7 activates brown adipose tissue and reduces diet-induced obe- sity only at subthermoneutrality.,8 (9): e74083, DOI: 10.1371/journal.pone.0074083.

Bussmann, J., Bos, F.L., Urasaki, A., Kawakami, K., Duckers, H.J., and Schulte-Merker, S., 2010. Arteries provide essential gui-dance cues for lymphatic endothelial cells in the zebrafish trunk.,137 (16): 2653-2657, DOI: 10.1242/dev.048207.

Cascone, I., Selimoglu, R., Ozdemir, C., Del Nery, E., Yeaman, C., White, M., and Camonis, J., 2008. Distinct roles of RalA and RalB in the progression of cytokinesis are supported by dis- tinct RalGEFs.,27 (18): 2375-2387, DOI: 10.1038/emboj.2008.166.

Chatziplis, D., Batargias, C., Tsigenopoulos, C.S., Magoulas, A., Kollias, S., Kotoulas, G., Volckaert, F.A., and Haley, C.S., 2007. Mapping quantitative trait loci in European sea bass (): The BASSMAP pilot study.,272: S172-S182, DOI: 10.1016/j.aquaculture.2007.08.022.

Chen, C., Xia, R., Chen, H., and He, Y., 2018. TBtools, a Toolkit for Biologists integrating various biological data handling toolswith a user-friendly interface.: 289660, DOI: 10.1101/289660.

Chen, J., Luo, Q., Bao, H., Liao, L., Li, Y., Zhu, Z., Wang, Y., and Hu, W., 2015. The integration characteristics of the exo- genous growth hormone gene in a transgenic common carp(L.) with fast-growth performance.,60 (19): 1654-1660, DOI: 10.1007/s11434-015-0893-x.

Chen, S., Zhang, G., Shao, C., Huang, Q., Liu, G., Zhang, P., Song, W., An, N., Chalopin, D., and Volff, J.N., 2014. Whole-genome sequence of a flatfish provides insights into ZW sex chromosome evolution and adaptation to a benthic lifestyle.,46 (3): 253, DOI: https://doi.org/10.1038/ng.2890.

Cingolani, P., Platts, A., Wang, L.L., Coon, M., Nguyen, T., Wang, L., Land, S.J., Lu, X., and Ruden, D.M., 2012. A pro- gram for annotating and predicting the effects of single nu- cleotide polymorphisms, SnpEff: SNPs in the genome ofstrain w1118; iso-2; iso-3.,6 (2): 80-92, DOI: https://doi.org/10.4161/fly.19695.

Danecek, P., Auton, A., Abecasis, G., Albers, C.A., Banks, E., De- Pristo, M.A., Handsaker, R.E., Lunter, G., Marth, G.T., and Sherry, S.T., 2011. The variant call format and VCFtools., 27 (15): 2156-2158, DOI: 10.1093/bioinformatics/btr330.

De-Santis, C., and Jerry, D.R., 2007. Candidate growth genes in finfish–Where should we be looking?,272 (1-4): 22-38, DOI: 10.1016/j.aquaculture.2007.08.036.

Devoto, S.H., Melançon, E., Eisen, J.S., and Westerfield, M., 1996. Identification of separate slow and fast muscle precursor cells, prior to somite formation.,122 (11): 3371-3380, DOI: 10.1101/gad.10.22.2935.

Docker, M.F., and Heath, D.D., 2002. PCR-based markers de- tect genetic variation at growth and immune function-related loci in chinook salmon ()., 2 (4): 606-609, DOI: 10.1046/j.1471-8286.2002.00315.x.

Dudley, A.T., Lyons, K.M., and Robertson, E.J., 1995. A re- quirement for bone morphogenetic protein-7 during develop- ment of the mammalian kidney and eye.,9 (22): 2795-2807, DOI: 10.1101/gad.9.22.2795.

Eom, G.H., Kim, K., Kim, S.M., Kee, H.J., Kim, J.Y., Jin, H.M., Kim, J.R., Kim, J.H., Choe, N., and Kim, K.B., 2009. Histone methyltransferase PRDM8 regulates mouse testis ste- roidogenesis.,388 (1): 131-136, DOI: 10.1016/j.bbrc.2009.07.134.

Fuentes, E.N., Valdés, J.A., Molina, A., and Björnsson, B.T.,2013. Regulation of skeletal muscle growth in fish by the growth hormone-insulin-like growth factor system.,192: 136-148, DOI: 10.1016/j.ygcen.2013.06.009.

Grillet, N., Schwander, M., Hildebrand, M.S., Sczaniecka, A., Kolatkar, A., Velasco, J., Webster, J.A., Kahrizi, K., Najma- badi, H., and Kimberling, W.J., 2009. Mutations in LOXHD1, an evolutionarily conserved stereociliary protein, disrupt hair cell function in mice and cause progressive hearing loss in humans.,85 (3): 328-337, DOI: 10.1016/j.ajhg.2009.07.017.

Gross, R., and Nilsson, J., 1999. Restriction fragment length poly- morphism at the growth hormone 1 gene in Atlantic salmon (L.) and its association with weight among the offspring of a hatchery stock.,173 (1-4): 73-80, DOI: 10.1016/S0044-8486(98)00470-0.

Hnia, K., Ramspacher, C., Vermot, J., and Laporte, J., 2015. Des- tumin in muscle and associated diseases: Beyond the struc- ral function.,360 (3): 591-608, DOI: 10.1007/s00441-014-2016-4.

Holsinger, K.E., and Weir, B.S., 2009. Genetics in geographi- cally structured populations: Defining, estimating and interpre- ting F(ST).,10 (9): 639, DOI: 10.1038/nrg2611.

Ibberson, M., Uldry, M., and Thorens, B., 2000. GLUTX1, a no- vel mammalian glucose transporter expressed in the central nervous system and insulin-sensitive tissues.,275 (7): 4607-4612, DOI: 10.1074/jbc.275.7.4607.

Ji, X.S., Liu, H.W., Chen, S.L., Jiang, Y.L., and Tian, Y.S., 2011. Growth differences and dimorphic expression of growthhormone (GH) in female and maleafter male sexual maturation.,4 (1): 9-16, DOI: 10.1016/j.margen.2010.11.002.

Johnston, I.A., Bower, N.I., and Macqueen, D.J., 2011. Growth and the regulation of myotomal muscle mass in teleost fish.,214 (10): 1617-1628, DOI: 10.1242/jeb.038620.

Johnston, I.A., Lee, H.T., Macqueen, D.J., Paranthaman, K., Ka- washima, C., Anwar, A., Kinghorn, J.R., and Dalmay, T., 2009.Embryonic temperature affects muscle fibre recruitment in adult zebrafish: Genome-wide changes in gene and microRNA expression associated with the transition from hyperplastic to hypertrophic growth phenotypes.,212 (12): 1781-1793, DOI: 10.1242/jeb.029918.

Jung, C.C., Atan, D., Ng, D., Ploder, L., Ross, S.E., Klein, M., Birch, D.G., Diez, E., and McInnes, R.R., 2015. Transcrip- tion factor PRDM8 is required for rod bipolar and type 2 OFF-cone bipolar cell survival and amacrine subtype identity.,112 (23): E3010-E3019, DOI: 10.1073/pnas.1505870112.

Katagiri, T., Boorla, S., Frendo, J.L., Hogan, B.L., and Karsenty, G., 1998. Skeletal abnormalities in doubly heterozygous Bmp4 and Bmp7 mice.,22 (4): 340-348, DOI: 10.1002/(SICI)1520-6408(1998)22:4<340::AID-DVG4>3.0.CO;2-6.

Kendall, R.L., and Thomas, K.A., 1993. Inhibition of vascular endothelial cell growth factor activity by an endogenously en- coded soluble receptor.,90 (22): 10705-10709, DOI: 10.1073/pnas.90.22.10705.

Li, H., and Durbin, R., 2009. Fast and accurate short read align- ment with Burrows-Wheeler transform.,25 (14): 1754-1760, DOI: 10.1093/bioinformatics/btp324.

Li, H., Handsaker, B., Wysoker, A., Fennell, T., Ruan, J., Homer, N., Marth, G., Abecasis, G., and Durbin, R., 2009. The sequence alignment/map format and SAMtools.,25 (16): 2078-2079, DOI: 10.1093/bioinformatics/btp352.

Li, W., Jin, D., and Zhong, T.P., 2019. Photoreceptor cell de- velopment requires prostaglandin signaling in the zebrafish re-tina.,510 (2): 230-235, DOI: 10.1016/j.bbrc.2019.01.073.

Lin, H., Lu, M., Lin, X., Zhang, W., Sun, Y., and Chen, L., 1995. Effects of gonadotropin-releasing hormone (GnRH) analogsand sex steroids on growth hormone (GH) secretion and growth in common carp () and grass carp ().,135 (1-3): 173-184, DOI: 10.1016/0044-8486(95)01012-2.

Lin, X.W., Lin, H.R., and Peter, R.E., 1993. Growth hormone and gonadotropin secretion in the common carp (L.):interactions of gonadotropin-releasing hormone, somatostatin, and the dopamine agonist apomorphine.,89 (1): 62-71, DOI: 10.1006/gcen.1993.1009.

Liu, Y., Chen, S., Gao, F., Meng, L., Hu, Q., Song, W., Shao, C., and Lv, W., 2014. SCAR-transformation of sex-specific SSR marker and its application in half-smooth tongue sole ().,22 (6): 787-792.

Luttrell, L.M., 2008. Reviews in molecular biology and bio-technology: Transmembrane signaling by G protein-coupled receptors.,39 (3): 239-264, DOI: 10.1007/s12033-008-9031-1.

Ma, Q., Liu, S.F., Zhuang, Z.M., Lin, L., Sun, Z. Z., Liu, C. L., Ma, H., Su, Y. Q., and Tang, Q. S., 2012. Genomic structure, polymorphism and expression analysis of the growth hormone (GH) gene in female and male half-smooth tongue sole ().,493 (1): 92-104, DOI: 10.1016/j.gene.2011.11.015.

Ma, Q., Liu, S.F., Zhuang, Z.M., Sun, Z.Z., Liu, C.L., Su, Y.Q., and Tang, Q.S., 2011. Molecular cloning, expression ana- lysis of insulin-like growth factor I (IGF-I) gene and IGF-I se- rum concentration in female and male Tongue sole ().,160 (4): 208-214, DOI: 10.1016/j.cbpb.2011.08.008.

Marcelli, F., Boisset, G., and Schorderet, D.F., 2014. A dime- rized HMX1 inhibits EPHA6/epha4b in mouse and zebrafish retinas.,9 (6): e100096, DOI: 10.1371/journal.pone.0100096.

Massault, C., Hellemans, B., Louro, B., Batargias, C., Van Houdt, J., Canario, A., Volckaert, F., Bovenhuis, H., Haley, C., and De Koning, D., 2010. QTL for body weight, morphometric traits and stress response in European sea bass., 41 (4): 337-345, DOI: 10.1111/j.1365-2052.2009.02010.x.

Ni, J., You, F., Xu, J., Xu, D., Wen, A., Wu, Z., Xu, Y., and Zhang, P., 2012. Single nucleotide polymorphisms in intron 1 and intron 2 ofgrowth hormone gene are correlated with growth traits.,30 (2): 279-285, DOI: 10.1007/s00343-012-1078-y.

Ooga, M., Funaya, S., Hashioka, Y., Fujii, W., Naito, K., Suzuki, M.G., and Aoki, F., 2018. Chd9 mediates highly loosened chro- matin structure in growing mouse oocytes.,500 (3): 583-588, DOI: 10.1016/j.bbrc.2018.04.105.

Quik, E.H., van Dam, P.S., and Kenemans, J.L., 2010. Growth hormone and selective attention: A review.,34 (8): 1137-1143, DOI: 10.1016/j.neubiorev.2010.01.001.

Rescan, P.Y., and Ralliere, C., 2010. Agene is expressed in the myogenic lineage during trout embryonic development.,54 (5): 913-918, DOI: 10.1387/ijdb.092893pr.

Robledo, D., Rubiolo, J.A., Cabaleiro, S., Martínez, P., and Bou- za, C., 2017. Differential gene expression and SNP associa- tion between fast-and slow-growing turbot ().,7 (1): 12105, DOI: 10.1038/s41598-017-12459-4.

Ross, S.E., McCord, A.E., Jung, C., Atan, D., Mok, S.I., Hem- berg, M., Kim, T.K., Salogiannis, J., Hu, L., and Cohen, S., 2012. Bhlhb5 and Prdm8 form a repressor complex involved in neuronal circuit assembly.,73 (2): 292-303, DOI: 10.1016/j.neuron.2011.09.035.

Ruiz-Perez, V.L., Blair, H.J., Rodriguez-Andres, M.E., Blanco, M.J., Wilson, A., Liu, Y.N., Miles, C., Peters, H., and Good- ship, J.A., 2007. Evc is a positive mediator of Ihh-regulated bone growth that localises at the base of chondrocyte cilia.,134 (16): 2903-2912, DOI: 10.1242/dev.007542.

Sabourin, L.A., and Rudnicki, M.A., 2000. The molecular re- gulation of myogenesis.,57 (1): 16-25, DOI: 10.1034/j.1399-0004.2000.570103.x.

Salem, M., Vallejo, R.L., Leeds, T.D., Palti, Y., Liu, S., Sab- bagh, A., Rexroad III, C.E., and Yao, J., 2012. RNA-Seq iden-tifies SNP markers for growth traits in rainbow trout., 7 (5): e36264, DOI: 10.1371/journal.pone.0036264.

Sánchez-Ramos, I., Cross, I., Mácha, J., Martínez-Rodríguez, G., Krylov, V., and Rebordinos, L., 2012. Assessment of tools for marker-assisted selection in a marine commercial species: Sig- nificant association between MSTN-1 gene polymorphism and growth traits.,2012: 369802, DOI: 10.1100/2012/369802.

Sanger, A., and Stoiber, W., 2001. Muscle fiber diversity and plas- ticity.,18: 187-250, DOI: 10.1016/S1546-5098(01)18008-8.

Seytanoglu, A., Alsomali, N.I., Valori, C.F., McGown, A., Kim, H.R., Ning, K., Ramesh, T., Sharrack, B., Wood, J.D., and Azzouz, M., 2016. Deficiency in the mRNA export mediator Gle1 impairs Schwann cell development in the zebrafish em- bryo.,322: 287-297, DOI: 10.1016/j.neurosci- ence.2016.02.039.

Shawi, M., and Serluca, F.C., 2008. Identification of a BMP7 homolog in zebrafish expressed in developing organ systems.,8 (6): 369-375, DOI: 10.1016/j.gep.2008.05.004.

Shur, I., Socher, R., and Benayahu, D., 2006.association of CReMM/CHD9 with promoters in osteogenic cells.,207 (2): 374-378, DOI: 10.1002/jcp.20586.

Song, W., Li, Y., Zhao, Y., Liu, Y., Niu, Y., Pang, R., Miao, G., Liao, X., Shao, C., and Gao, F., 2012. Construction of a high-density microsatellite genetic linkage map and mapping of sexual and growth-related traits in half-smooth tongue sole ()., 7 (12): e52097, DOI: 10.1371/journal.pone.0052097.

Takahashi, T., Hagiwara, A., and Ogiwara, K., 2018. Prostaglan- dins in teleost ovulation: A review of the roles with a view to comparison with prostaglandins in mammalian ovulation., 461: 236-247, DOI: 10.1016/j.mce.2017.09.019.

Valente, L.M., Moutou, K.A., Conceição, L.E., Engrola, S., Fernandes, J.M., and Johnston, I.A., 2013. What determines growth potential and juvenile quality of farmed fish species?,5: S168-S193, DOI: 10.1111/raq.12020.

Wang, D.Y., Fulthorpe, R., Liss, S.N., and Edwards, E.A., 2004. Identification of estrogen-responsive genes by complemen- tary deoxyribonucleic acid microarray and characterization of a novel early estrogen-induced gene:., 18 (2): 402-411, DOI: 10.1210/me.2003-0202.

Wang, H., Chen, Y., Lin, P., Li, L., Zhou, G., Liu, G., Logsdon, C., Jin, J., Abbruzzese, J.L., and Tan, T.H., 2014. The CUL7/F-box and WD repeat domain containing 8 (CUL7/Fbxw8) ubiquitin ligase promotes degradation of hematopoietic pro- genitor kinase 1.,289 (7): 4009-4017, DOI: 10.1074/jbc.M113.520106.

Weir, B.S., and Cockerham, C.C., 1984. Estimating F-statistics for the analysis of population structure.,38 (6): 1358-1370, DOI: 10.1111/j.1558-5646.1984.tb05657.x.

Wilkinson, D.G., 2000. Eph receptors and ephrins: Regulators of guidance and assembly., 196: 177-244, DOI: 10.1016/S0074-7696(00)96005-4.

Xia, J.H., Lin, G., He, X., Liu, P., Liu, F., Sun, F., Tu, R., and Yue, G.H., 2013. Whole genome scanning and association mapping identified a significant association between growth and a SNP in the IFABP-a gene of the Asian seabass.,14 (1): 295, DOI: 10.1186/1471-2164-14-295.

Yano, K., Yamamoto, E., Aya, K., Takeuchi, H., Lo, P. C., Hu, L., Yamasaki, M., Yoshida, S., Kitano, H., and Hirano, K., 2016. Genome-wide association study using whole-genome sequen- cing rapidly identifies new genes influencing agronomic traits in rice.,48 (8): 927, DOI: 10.1038/ng.3596.

Zhang, B., Shang, P., Tao, Z., Qiangba, Y., Wang, Z., and Zhang, H., 2017. Effect of a single nucleotide polymorphism in the growth hormone secretagogue receptor () gene on growth rate in pigs.,634: 68-73, DOI: 10.1016/j.gene.2017.09.007.

June 3, 2020;

July 22, 2020;

December 1, 2020

© Ocean University of China, Science Press and Springer-Verlag GmbH Germany 2021

. E-mail: yanhe@ouc.edu.cn

(Edited by Qiu Yantao)