Sequencing assembly and annotation
A total of 263.21M raw reads for reference database was generated based on BGISEQ-500 sequencer and 262.05M clean reads were remained after filtering out the low quality reads. For the clean reads, more than 97.89% of the bases had a phred value > 20, and more than 89.35% of them had a phred value > 30. An overview of the assembly results was provided in Additional file 2 Table S2. The raw data were uploaded to the NCBI sequence read archive database with the accession number PRJNA632458.
In total, 48,339 unigenes were generated and annotated by searching the sequences against the NR, KEGG, GO and TFs databases, which produced 17,737 (36.69%), 13,336 (27.59%), 6,128 (12.68%) and 1,756 (3.63%) hits respectively (Table 1). There were only 1,752 (3.62%) were annotated in all of the databases, and 17,847 (36.92%) were annotated in at least one database. Among these, GO classification analysis showed that 48,339 unigenes were categorized into 21 biological processes, 16 cellular components and 9 molecular functions. The KEGG classification analysis showed that total unigenes were assigned to 6 KEGG categories in the top level: 12 metabolisms, 11 human diseases, 10 organismal systems, 4 cellular processes, 4 genetic information processing and 3 environmental information processing.
DEGs analysis and immune related pathway identification
DEGs analysis among the four groups (Pre_F, Post_F, Pre_M and Post_M)was applied to identify gene changes pre/post-spawning. All DEGs with the absolute value of fold change > 2, adjust P-value < 0.01 are showed in Fig.1. A total of 1,805 DEGs were detected between the Pre_F vs. Post_F, of which 722 and 715 unigenes were up-regulated and down-regulated respectively. For the Pre_M vs. Post_M comparison, 2,488 DEGs were detected, of which 772 and 1716 unigenes were up-regulated and down-regulated respectively. As shown in Fig. 1, the number of up-regulated genes were basically equal between male groups (Pre_M vs. Post_M) and female groups (Pre_F vs. Post_F), but the number of down-regulated genes in male was more than double in female after spawning.
To further explore the function of the identified DEGs, GO classification analyses were conducted and 7 categories attracted our attention, which were cellular process, metabolic process, biological regulation, cell, organelle, binding and catalytic activity. In these categories, while the number of DEGs in male was at least 40% more than that of DEGs in female, the proportion of down-regulated genes in male was at least 50% more than that of down-regulated genes in female. These categories were marked in Figure 2, and a list of DEGs numbers was provided in Table 2. Additionally, in the majority of KEGG classification categories, the number of DEGs in both genders were approximately in the same, except “human diseases” that the number of DEGs in male was 80%(256/141)more than that of DEGs in female, and the proportion of down-regulated genes in male was 171% (192/39) more than that of down-regulated genes in female. The function of the more down-regulated DEGs in male is to participate in metabolism, catalytic activity, cell function and disease, which implied male reduced more these competences after spawning.
GO term and KEGG pathway enrichment analyses were performed to explore the DEGs related to immune regulation. We screened differential genes on GO term by keywords “immune” in Pre_F vs. Post_F and Pre_M vs. Post_M comparisons, and enriched 11 and 4 differential genes, respectively. Here, both quantity of DEGs and the type of GO enrichment were more in females than in males (Additional file 3, Fig.S1). To better understand the immune-related pathways changes pre/post-spawning, we evaluated significant enriched pathways about “Immune system” of secondary level of KEGG pathway enrichment, and enriched 19 and 20 differential genes in female and male respectively (Fig. 3). Interestingly, common pathways with high value enriched in different genders included IL-17 signaling pathway, NOD-like receptor signaling pathway, TLR signal pathway, RIG-I-like receptor signaling pathway, which played an important role in innate immunity. The results here implied that immune-related cytokines and inflammation were indeed affected by spawning.
Differences expression of immune-related genes between reproductive states
Based on the DEGs enriched by KEGG pathway in the Section 3.2 (19 and 20 genes in female and male respectively), 22 genes were further screened by Nr annotation. Based on their function in pathway, we classified them into seven categories: complement and coagulation cascades, apoptosis, TLR signal pathway, NFκB signal pathway, leukocyte transendothelial migration, MAPK signal pathway and others (Fig. 4). The expressions of genes in complement and coagulation cascades are almost 2/3 significantly down-regulated after spawning in males, whereas the analysis showed the expressions of genes in apoptosis, TLR signal pathway, NFκB signal pathway and HSP family are significantly up-regulated after spawning in females. It means that male decreases immune-competent of complement, while female improves immune function in apoptosis, antimicrobial peptides and heat shock.
To validate the reliability of DEGs identified by RNA-Seq, 10 out of these 22 genes were selected for further relative expression quantified by qRT-PCR. These genes belong to NFκB signal pathway (IKK, DTHD1), TLR signal pathway (TLR4), apoptosis (BIRC7, IAP and CREM), MAPK signal pathway (ARAF-like, MAPKK3), complement and coagulation cascades (CTRP2) and others immune-related genes (PSME3). The expressions of these genes pre/post-spawning showed in Fig.5. Compared their expression in RNA-seq (Fig. 4), 8/10 genes had a common expression pattern, only BIRC7 and IAP are not consistent in male.
Gene expression changes pre/post-spawning under Vibrio challenge
The transcriptome analysis revealed that immune-related genes changes pre/post-spawning in clam M. petechialis in natural condition. Here, we want to know how these genes expression pre/post-pawning under the pathogen challenge. Four immune-related DEGs (ARAF-like, CREM, TLR4 and TBK1) were selected for further examination. As shown in Fig. 6, without Vibrio challenge (0h), the expression patterns of these DEGs were consistent with the changes revealed by RNA-seq under natural condition (Fig. 4), that was, after spawning their expressions were significantly increased in females but there was no significant difference in males. Under Vibrio challenge, the four genes also showed up-regulated expression in females at 12h or 24h; in male, the non-differentially expressed genes showed a significant up-regulated expression after spawning at 24h. These results indicated that the response of these immune-related genes became more sensitive and their expression changes more drastically after spawning under the pathogen challenge. It means pathogen reinforced the expression differences of immune-related genes between pre-spawning groups and post-spawning groups.
Assessment of immune status of clam under Vibrio challenge
To further verify the effects of spawning, we chose two previously reported genes, DOUX and MITF both involving the immune defense of M. petechialis [24, 37], to detect the immune changes of the clams pre/post- spawning. As shown in Fig. 7, under bacterial immersion challenge, the mRNA expression of MpDOUX and MpMITF were significantly up-regulated after spawning in both genders. We also tested the protein changes of MpDOUX and MpMITF in female clams by western blot because we perceived that many immune genes had more significant differences after spawning in female based on transcriptome data (Fig. 4). As shown on Fig. 8 and Additional file 4, the protein expressions increased at 24 h for MpDOUX and at 12 h for MpMITF in the post-spawning individuals under Vibrio challenge. The up-regulation at the transcription and protein levels of these maker genes further confirmed our conclusion that pathogen reinforced the expression differences of immune-related genes between pre-spawning and post-spawning groups.