Transcriptome Analysis Reveals the Immunosuppression in Tiger Pufferfish (Takifugu rubripes) under Cryptocaryon irritans Infection

Abstract

The tiger pufferfish (Takifugu rubripes), also known as fugu, has recently suffered from severe C. irritans infections under aquaculture environment, yet the underlying immune mechanisms against the parasite remain poorly understood. In this study, we conducted a comprehensive transcriptome analysis of the gill tissue from infected and uninfected fish using PacBio long-read (one pooled sample each for seriously infected and healthy individuals, respectively) and Illumina short-read (three pools for mildly infected, seriously infected, and healthy individuals, respectively) RNA sequencing technologies. After aligning sequence data to fugu's reference genome, 47,307 and 34,413 known full-length transcripts were identified and profiled in healthy and infected fish, respectively. Similarly, we identified and profiled 1126 and 803 novel genes that were obtained from healthy and infected fish, respectively. Interestingly, we found a decrease in the number of alternative splicing (AS) events and long non-coding RNAs (lncRNAs) after infection with C. irritans, suggesting that they may be involved in the regulation of the immune response in fugu. There were 687 and 1535 differentially expressed genes (DEGs) in moderately and heavily infected fish, respectively, compared to uninfected fish. Kyoto Encyclopedia of Genes and Genomes (KEGG) analyses showed that immune-related DEGs in the two comparison groups were mainly enriched in cytokine-cytokine receptor interactions, ECM-receptor interactions, T-cell receptor signaling pathways, Th1 and Th2 cell differentiation, and Th17 cell differentiation pathways. Further analysis revealed that a large number of immune-related genes were downregulated in infected fish relative to uninfected ones, such as CCR7, IL7R, TNFRSF21, CD4, COL2A1, FOXP3B, and ITGA8. Our study suggests that C. irritans is potentially a highly efficient parasite that may disrupt the defense mechanisms of fugu against it. In addition, in combination of short-read RNA sequencing and previous genome-wide association analyses, we identified five key genes (NDUFB6, PRELID1, SMOX, SLC25A4, and DENND1B) that might be closely associated with C. irritans resistance. This study not only provides valuable resources of novel genic transcripts for further research, but also provides new insights into the immune mechanisms underlying C. irritans infection response in farmed fugu.

Keywords: RNA-seq; aquaculture; fugu; immune response; parasite.

Figure 1

Characteristics of transcripts and novel gene annotation. (A) Classification statistics of full-length transcripts in the FL-HG. (B) Classification statistics of full-length transcripts in the FL-DG. (C) Statistic results of the NR database, Swiss–Prot protein database, KEGG database, and GO database. The abscissa represents different databases, and the ordinate represents the number of annotations.

Figure 2

Alternative splicing (AS) analysis based on the full-length transcripts. (A) Seven types of AS forms. Blue blocks represent constitutive exons, orange blocks represent alternatively spliced exons. (B) Statistic results of AS events: SE (skipping exon), MX (mutually exclusive exon), A5 (alternative 5′ splice site), A3 (alternative 3′ splice site), RI (retained intron), AF (alternative first exon), and AL (alternative last exon).

Figure 3

Long non-coding RNA (lncRNA) analysis based on the full-length transcripts. (A) Venn diagram of lncRNA prediction results in the FL-HG. (B) Venn diagram of lncRNA prediction results in the FL-DG. (C) Statistical graph of lncRNA classification results in the FL-HG. (D) Statistical graph of lncRNA classification results in the FL-DG.

Figure 4

Sample relationship analysis. (A) Principal component analysis (PCA) of the genes in terms of variance across samples. (B) Sample correlation heat map.

Figure 5

DEGs expression analysis. (A) Volcano plot of DEGs in the HG vs. MG group. (B) Volcano plot of DEGs in the HG vs. SG group. (C) Hierarchical clustering analysis DEGs in the HG vs. MG group. (D) Hierarchical clustering analysis of DEGs in the HG vs. SG group.

Figure 6

GO and KEGG function enrichment analysis of DEGs. (A) GO function enrichment analysis of the HG vs. MG group (top 20 enriched terms). (B) KEGG function enrichment analysis of the HG vs. MG group (20 enriched terms). (C) GO function enrichment analysis of the HG vs. SG group (top 20 enriched terms). (D) KEGG function enrichment analysis of the HG vs. SG group (top 20 enriched terms).

Figure 7

Pathway based GSEA results of the HG vs. MG (A) and HG vs. SG (B) groups.

Figure 8

PPI networks of selected key DEGs. The red color indicates the hub genes.