Generation of autotriploid Carassius auratus
3nRR were generated by crossing female RCC and male 4nRR during the breeding season (Fig. 1a, b, c, Fig. 2). Chromosome number analysis showed that 92.5% of metaphases had 100 chromosomes in RCC, implying that they were diploids with 100 chromosomes (2n = 100) (Fig. 3a, Table 1). Analysis of 4nRR showed that 87.0% of metaphases had 200 chromosomes, implying that they were tetraploids with 200 chromosomes (4n = 200) (Fig. 3b, Table 1). Analysis of 3nRR showed that 91.5% of metaphases had 150 chromosomes, suggesting that they were triploids with 150 chromosomes (3n = 150) (Fig. 3c, Table 1).
Fertility of autotriploid Carassius auratus
Testes of RCC and 4nRR (Fig. 4a, b) contained spermatogonia (SG), spermatocytes (SC) and a large number of mature spermatid (ST), whereas the mature sperm was not observed in 3nRR (Fig. 4c). Ovaries of RCC, 4nRR and 3nRR contained second, third and fourth phase oocytes (Fig. 4d, e, f). These results indicated that all ovaries, and RCC and 4nRR testes were fertile whereas 3nRR testes were sterile.
Different sizes of eggs and water-like semen were collected during the reproductive season from two years old males and females of 3nRR, respectively (Fig. 5). Ploidy levels of the offspring resulting from a cross of female 3nRR and male RCC (Figure 1c, d, e, f, g) were determined by measuring the chromosome number (Figure 3d, e, f, Table 1). Analysis showed that female 3nRR produced haploid, diploid and triploid eggs.
Transcriptome sequencing and sequence alignment
Optical density (OD) ratio A260/A280 and RNA integrity numbers (RINs) of the RNA in 12 samples (Table 2) were 2.1 and 8.0-8.8, respectively (Additional file 1). These results indicate that all samples were free from contamination and their quality met the requirements for transcriptome sequencing.
RNA-seq from gonadal tissue samples of autotriploid fish and their parents was performed by Illumina. RNA-Seq results are presented in Table 3 and Table 4. Number of clean reads from the 12 RNA-seq libraries ranged from 39,624,312 to 50,588,484. All clean reads were then aligned to the RCC genome sequences using HISAT2 software. Mapped genome reads ranged from 24,237,536 to 42,474,296, genome map rates ranged from 59.79% to 91.55%, and unique match rates ranged from 57.84% to 85.93%.
Identification of differentially expressed genes (DEGs)
Analysis of F3nRR and FRCC showed that a total of 13,467 DEGs were downregulated whereas 6,688 DEGs were up-regulated (Fig. 6a). DEGs between F3nRR and FRCC included forkhead box L2 (FOXL2), LIM homeobox 8 (LHX8), lysine acetyltransferase 8 (KAT8), BCL2 apoptosis regulator (BCL2), doublesex and mab-3 related transcription factor 1 (DMRT1), ovarian serine protease (OSP) and CCM2 scaffold protein (CCM2). Analysis of M3nRR and F4nRR showed that a total of 1,886 DEGs were downregulated and 2,221 DEGs were up-regulated (Fig. 6b). DEGs between M3nRR and F4nRR included septin 12 (SEPT12), ATPase copper transporting beta (ATP7B), CF transmembrane conductance regulator (CFTR), cAMP responsive element modulator (CREM), cytochrome P450 family 26 subfamily B member 1 (CYP26B1), EF-hand calcium binding domain 2 (EFCAB2) and inhibitor of kappa light polypeptide gene enhancer in B-cells and kinase complex-associated protein (IKBKAP).
GO and KEGG enrichment analysis of DEGs
GO enrichment analysis of the biological process, cellular component and molecular function categories yielded 242, 38 and 51 terms, respectively, for F3nRR vs FRCC, and 223, 28 and 29 for M3nRR vs M4nRR group. (Additional file 2, 3). The most-enriched GO-terms for F3nRR vs FRCC group were “induction of programmed cell death” in the biological process category, “neuron projection” in the cellular component category, and “channel activity” and “passive transmembrane transporter activity” in the molecular function category. The most-enriched GO-terms for M3nRR vs F4nRR group were “extracellular region part” in the cellular component category; “kinase activity” and “transferase activity, transferring phosphorus-containing groups” in the molecular function category; and “response to osmotic stress” in the biological process category (Fig. 7).
KEGG analysis of all DEGs showed that 44 and 62 signaling pathways were enriched in the F3nRR vs FRCC group and M3nRR vs M4nRR group, respectively (Additional file 4, 5). The top 20 most enriched KEGG pathways are shown in Figure 8. The five most-enriched pathways in the F3nRR vs FRCC group were “ion channels” (ko04040), “cAMP signaling pathway” (ko04024), “focal adhesion” (ko04510), “glycosaminoglycan binding proteins” (ko00536) and “glycosyltransferases” (ko01003). Moreover, several pathways implicated in female fertility of F3nRR were identified, including “MAPK signaling pathway - plant” (ko04016), and “p53 signaling pathway” (ko04115). The five most enriched pathways for the M3nRR vs M4nRR group were “Ion channels” (ko04040), “rap1 signaling pathway” (ko04015), “ras signaling pathway” (ko04014), “alcoholism” (ko05034) and “axon guidance” (ko04360). Notably, four of the top 20 most-enriched pathways, “regulation of actin cytoskeleton” (ko04810), “calcium signaling pathway” (ko04020), “tight junction” (ko04530) and “cytokines and growth factors” (ko04052), play important roles in cellular processes such as differentiation, proliferation, migration and apoptosis, implying that they are potentially involved in male sterility of M3nRR.
Hub genes related to the fertility in 3nRR were identified
GO and KEGG enrichment analyses and published literature showed that a total of 80 genes out of the DEGs identified in the F3nRR vs FRCC group were female fertility-related genes. On the other hand, 25 genes out of the DEGs in the M3nRR vs M4nRR group are implicated in male sterility. To further identify hub genes associated with 3nRR fertility, PPI of the fertility-related genes was constructed using STRING tool and analysis was carried out using Cytoscape software. Analysis of PPI network of female fertility-related genes showed that 6 hub genes with degrees more than 15 showed strong interaction with other node proteins (Fig. 9a, Additional file 6). Furthermore, PPI of male sterility-related genes showed that 2 hub genes, with degrees more than 5 showed strong interaction with other node proteins (Fig. 9b, Additional file 7).
RT-qPCR verification of the DEGs
To verify RNA-Seq results, twenty DEGs were chosen for validation by RT-qPCR. Among the 20 DEGs, 5 DEGs were up-regulated in the F3nRR vs FRCC and M3nRR vs M4nRR groups, respectively; whereas 5 DEGs were down-regulated in the F3nRR vs FRCC and M3nRR vs M4nRR groups, respectively. Expression profiles of the twenty DEGs obtained by RT-qPCR and RNA-Seq were similar, implying that RNA-Seq results were reliable (Fig. 10).