Transcriptome Analysis Reveals Comprehensive Insights into the Early Immune Response of Large Yellow Croaker (Larimichthys crocea) Induced by Trivalent Bacterial Vaccine

Vaccination is an effective and safe strategy for combating bacterial diseases in fish, but the mechanisms underlying the early immune response after vaccination remain to be elucidated. In the present study, we used RNA-seq technology to perform transcriptome analysis of spleens from large yellow croaker (Larimichthys crocea) induced by inactivated trivalent bacterial vaccine (Vibrio parahaemolyticus, Vibrio alginolyticus and Aeromonas hydrophila). A total of 2,789 or 1,511 differentially expressed genes (DEGs) were obtained at 24 or 72 h after vaccination, including 1,132 or 842 remarkably up-regulated genes and 1,657 or 669 remarkably down-regulated genes, respectively. Gene ontology and Kyoto Encyclopedia of Genes and Genomes enrichments revealed that numerous DEGs belong to immune-relevant genes, involved in many immune-relevant pathways. Most of the strongly up-regulated DEGs are innate defense molecules, such as antimicrobial peptides, complement components, lectins, and transferrins. Trivalent bacterial vaccine affected the expressions of many components associated with bacterial ligand–depending Toll-like receptor signaling pathways and inflammasome formation, indicating that multiple innate immune processes were activated at the early period of vaccination in large yellow croaker. Moreover, the expression levels of genes involved in antigen processing were also up-regulated by bacterial vaccine. However, the expression levels of several T cell receptors and related CD molecules and signal transducers were down-regulated, suggesting that the T cell receptor signaling pathway was rapidly suppressed after vaccination. These results provide the comprehensive insights into the early immune response of large yellow croaker to vaccination and valuable information for developing a highly immunogenic vaccine against bacterial infection in teleosts.


Introduction
Large yellow croaker (Larimichthys crocea) is one of the most important aquaculture species in China, with the annual yield exceeding any other net-cage-farmed marine fish species [1].However, the diseases caused by bacteria, such as Vibrio parahaemolyticus [2], Vibrio alginolyticus [3], and Aeromonas hydrophila [4], outbreak frequently and result in tremendous economic losses in large yellow croaker aquaculture industry.Chemotherapeutants and antibiotics are useful for preventing bacterial infection, but long term treatment often leads to resistance and environmental pollution [5].Vaccination is an effective and safe strategy for combating bacterial diseases in fish [6].A number of vaccines are commercially available for use in the aquaculture industry, such as formalin-inactivated Yersinia ruckeri [7], Aeromonas salmonicida [8], and Vibrio anguillarum [9].The effects of vaccines were mainly evaluated by serum or mucus antibody titers and survival rate [10].Recently, several studies on the gene expression changes have been performed to understand the immune response to bacterial vaccines in fish.In Atlantic salmon (Salmo salar), the increased expressions of antibacterial proteins and proteases were observed after vaccination with live Aeromonas salmonicida [11].In zebrafish (Danio rerio), multiple pathways including acute phase response, complement activation, and antigen processing and presentation were remarkably affected during the early stage of vaccination [12].The genes involved in inflammation and antioxidant defense were also upregulated in Atlantic cod (Gadus morhua) following vaccination with heat-killed Vibrio anguillarum [13].
An inactivated trivalent bacterial vaccine, consisting of V. alginolyticus, V. parahaemolyticus, and A. hydrophila, was developed and provided an effective protection against bacterial infection in large yellow croaker, with a relative percent survival (RPS) of 88.9% [14].This trivalent bacterial vaccine has been used to study the expression patterns of immune-relevant genes in large yellow croaker, such as goose-type lysozyme [15], caspase 9 [16], C-type lectinlike receptor [17], CXCL8 [18], CXCL13 [19], and stefin [20], indicating that this bacterial vaccine could induce the expression of innate immunity molecules in large yellow croaker.Currently, high-throughput RNA sequencing (RNA-seq) technology provides a powerful and cost-effective approach to understand the immune response evoked by different stimulus.For instance, RNA-seq was applied to identify immune-candidate genes, infection markers, and putative signaling pathways in rainbow trout (Oncorhynchus mykiss) [21], tilapia (Oreochromis niloticus) [22], and grass carp (Ctenopharyngodon idellus) [23] after bacterial infection.
Here, we used RNA-seq to examine the transcriptional profiles of the large yellow croaker spleens at early time points (0, 24, and 72 h) following immunization with inactivated trivalent vaccine (V.alginolyticus, V. parahaemolyticus and A. hydrophila).A series of differentially expressed genes (DEGs) were identified and found to be involved in the innate immunity and acquired immunity processes.Moreover, quantitative real-time PCR was used to verify expression changes of partial immune-relevant genes.These data provide valuable information to better understand the early immune response after vaccination in large yellow croaker.

Ethics statement
All handling of fish was approved by the Institutional Animal Care and Use Committee of Fujian Province, and performed in strict accordance with the Regulations for the Administration of Affairs Concerning Experimental Animals, under protocol license number: SYXK (MIN) 2007-0004.All surgeries were carried out under Tricaine-S anesthesia, and all efforts were made to minimize suffering.

Inactivated trivalent bacterial vaccine and fish samples
The V. parahaemolyticus, V. alginolyticus and A. hydrophila were isolated from the diseased large yellow croaker previously [14].For vaccine preparation, the three bacteria were cultured in Luria-Bertain (LB) for 16 h in a shaker incubator (180 rpm) at 30˚C, respectively.After collected and washed with the sterilized phosphate-buffered saline (PBS, pH 7.4), V. parahaemolyticus, V. alginolyticus and A. hydrophila were mixed at same proportion (1.0×10 9 CFU/mL each) and inactivated with 0.5% formalin (v/v) at 4˚C for 4 h.The mixed inactivated trivalent vaccine was washed with the sterilized PBS.To evaluate the effectiveness of the formalin-inactivation, 50 μL of formalin-inactivated suspensions were spread on LB agar (1.5%, w/v) and incubated at 37˚C for 24 h.No bacterial growth was observed, indicating that the bacterial cells were completely inactivated.
Large yellow croakers (length: 20 ± 1.76 cm; weight: 100 ± 21.5 g) were obtained from the Mari-culture farm in Lianjiang, Fuzhou, China.The fish were maintained at 20˚C in aerated water tanks with a flow-through seawater supply.After one week of acclimatization, a group of 30 fish were injected intraperitoneally with inactivated trivalent bacterial vaccine (1.0×10 9 CFU/mL each) at a dose of 0.2 mL/100 g fish.The spleen tissues were harvested from 10 large yellow croakers at each time point (0, 24, and 72 h) after vaccination.

RNA extraction, library preparation, and sequencing
Total RNA was isolated from pooled spleen tissues of 6 fish sampled at each time point for RNA-seq analysis using TRIZOL1 Reagent (Invitrogen, USA), respectively, and treated with RNase-free DNase I (Takara, China) at 37˚C for 1 h to remove residual genomic DNA.The RNA samples (5 μg each) were used to construct RNA-seq libraries by Illumina mRNA-Seq Prep Kit, and the libraries were sequenced by paired-end sequencing on the Illumina HiSeq 2000 sequencing platform (Illumina, USA).

Identification of the differentially expressed genes
Clean reads were picked out from the raw reads of each library following removal of the reads containing adaptor sequences, reads with an N (unknown bases in read) percentage higher than 10%, and low quality reads (> 50% of bases with a quantity score Q-value 10) using an in-house C ++ program.The clean reads were aligned to the large yellow croaker genome (JRPU00000000) [24] using the software Burrows-Wheeler Aligner (BWA, parameter: -o 1 -e 63 -i 90 -L -k 2 -l 31 -t 4 -q 10) [25], and mapped to the coding sequences (CDS) of large yellow croaker genome using Bowtie2 (Parameter: -q -phred64 -sensitive -dpad 0 -gbar 99999999 -mp 1,1 -np 1 -score-min L,0,-0.1 -I 1 -X 1000 -no-mixed -no-discordant -p 1 -k 200) [26].The gene expression levels were calculated based on FPKM (fragments per kilobase of transcripts per million fragments mapped) values by using RSEM (RNA Seq by Expectation Maximization) with default setting [27].Comparisons of gene expression difference between control (0 h) and other time points after vaccination (24 and 72 h) were performed by strict Poisson distribution algorithm.Genes with the absolute value of the log 2 Ratio ! 1 and FDR (false discovery rate) 0.001 were defined as DEGs [28].

Gene ontology and kyoto encyclopedia of genes and genomes (KEGG) analysis
All the DEGs were mapped to the gene ontology (GO) database (http://www.geneontology.org/) using Blast2GO software [29].KEGG Orthology (KO) was carried out by BLASTx based on the KEGG database (http://www.genome.jp/kegg/).Finally, GO and KO enrichment analysis were performed according to the hypergeometric distribution test.For GO enrichment analysis, all of the P-values were performed with Bonferroni correction.We have selected a corrected P-value of 0.05 as a threshold to determine significant enrichment of the gene sets.In contrast, for KO enrichment analysis, we have used a Q-value of 0.05 as the threshold to determine significant enrichment of the gene sets.

Quantitative real-time PCR
Quantitative real-time PCR was performed using Mastercycler epgradient realplex4 (Eppendorf, Germany) with SYBR Green (Takara, China).Primer sets were designed based on the coding sequences of identified genes in large yellow croaker genome and specific PCR product of each gene was observed on the agarose gels (S1 Table and S1 Fig) .Total RNA was extracted from pooled spleen tissues of three fish sampled at each time point (0, 24, and 72 h) after vaccination for Real-time PCR analysis.Real-time PCR was performed in a total volume of 20 μL, and cycling conditions were 95˚C for 1 min, followed by 40 cycles of 95˚C for 15 s, 58˚C for 15 s, 72˚C for 20 s.Each real-time PCR assay was repeated three times with different batches of fish.The expression level of target gene was normalized by that of β-actin using the 2 -ΔΔCT method [30].The fold change was calculated as the average expression level of target gene in trivalent vaccine-challenged samples at each time point (24 and 72 h) divided by that in the control sample (0 h).The data of real-time PCR was analyzed using GraphPad Prism 5 software and expressed as the standard error of the mean (SEM).Two-tailed Student's t test was used for the significance test between the experimental group and the control group.A P-value <0.05 was considered to be statistically significant.

Sequence analysis of the transcriptomes
Three transcriptomes were sequenced from pooled spleens of large yellow croaker after vaccination.More than 5.4 gigabase (Gb) of data were generated in each library (Table 1).After eliminating adaptor sequences and low-quality sequences, a total of 63,568,408, 60,346,074, and 60,445,278 clean reads were obtained from the three libraries (0, 24, and 72 h), respectively.Average 70% of the clean reads could be mapped to the genome of large yellow croaker with 19,648, 19,106, and 19,145 genes, respectively.Raw sequencing reads data has been deposited in SRA (NCBI) under the accession number of SRP092778.

Analysis of differential expression genes
The gene expression levels were calculated based on FPKM values, and gene expression difference was analyzed between each time point (24 and 72 h) after vaccination and control (0 h).At 24 h, 2,789 genes were differentially expressed (log 2 Ratio ! 1, FDR 0.001), including 1,132 up-regulated genes and 1,657 down-regulated genes (Fig 1).Meanwhile, a total of 1,511 differentially expressed genes were found at 72 h after vaccination, including 842 up-regulated genes and 669 down-regulated genes.All the differentially expressed genes were provided in S2 Table.

Annotation of immune-relevant genes
To better understand the early immune response of large yellow croaker immunized with bacterial vaccine, the immune-relevant DEGs were collected and clustered (Table 3).Lots of  3), and some immunoglobulins (Ig heavy chain V region, Ig heavy chain V-III region CAM, Ig heavy chain V-III region HIL, Ig kappa chain V-IV region JI) were down-regulated (S2 Table ).

Verification of the expression changes of partial immune-relevant genes
Partial DEGs were validated by quantitative real-time PCR.The mRNA transcripts of C3, CLEC4E, F-lectin 1, LysG, Hep-1, TLR5M, AP-1, IL-1β, IL-12, TNF-α, TAP1, PSMA1, CTSL, and CTSO were observably increased at 24 or 72 h after vaccination, while the expression levels of TCRβ were down-regulated (Fig 5).Fold changes of those genes from quantitative real-time PCR were well corresponding to the results of RNA-seq expression analysis, supporting the reliability of the transcriptome data.

Discussion
To date RNA-seq based transcriptome profiling has been widely used in fish for identifying host determinants of response to bacterial infection [22,23,[31][32][33].In previous study, we have demonstrated that A. hydrophila infection could activate Toll-like receptor, JAK-STAT, and MAPK pathways, while inhibit TCR signaling pathway at the early period in large yellow croaker [34].In this study, we analyzed the spleen transcriptome profiles of large yellow croaker after vaccination and found abundant DEGs were involved in innate defense, ligand recognition, antigen processing, and TCR signaling processes.Innate immune response is considered as the first line of host defense in opposing pathogenic infection in fish.The components of the innate defense are evolutionarily conserved in organisms, including antimicrobial peptides, lysozymes, complement, lectins, and transferrins  [35].In this study, the expression of several innate defense molecules was increased obviously after vaccination (Table 3).Their up-regulations were further confirmed by real-time PCR, which were well corresponding to the results from transcriptome analysis (Fig 5A).The similar findings were also found in the studies on Atlantic salmon [11]and Atlantic cod [13], indicating that bacterial vaccine immunization could significantly increase the expression of innate immunity molecules in fish.The activation of complement system usually initiates a cascade of biochemical reactions, resulting in antigen elimination via cell membrane lysis [36].C3 and C7 of large yellow croaker have been proved to participate in the response to bacterial challenge [37,38].Lectins possess binding activity towards carbohydrate residues located on bacterial surfaces and aid neutralization of pathogens [39].The transferrins, lysozymes, and antibacterial peptides function as the growth inhibitors of bacteria against the infection of bacterial pathogens [15,40,41].These data suggest that immunization with trivalent vaccine intensely enhances the innate defense of large yellow croaker.
The activation of innate immune system against invading microbes mainly relies on the recognition of microbial components by pattern-recognition receptors (PRRs) [42].Several PRRs, including TLR5M, NLRP1, NLRP3, MMR1, PGRP-SC2, and CTLR were up-regulated significantly after vaccination (Table 3 and Fig 5B).TLR5M is known to recognize bacterial flagellin, signal through MyD88 and TRAF6, trigger the activation of MAPK and NF-κB, and induce inflammatory factors production [43].Due to V. alginolyticus, V. parahaemolyticus and A. hydrophila are all flagellated bacteria, TLR5M was up-regulated by 7.8-fold at 24 h after vaccination with trivalent bacterial vaccine.Eventually, its downstream signal transducers (TRAF6, IκBα, ERK), transcriptional factors (RelB, AP-1), and effectors (TNF-α, IL-6, IL-12, CCL4) were also up-regulated, indicating TLR5M mediated signaling pathway was activated by trivalent vaccine (Fig 3).Other signal transducers, such as IRAK1/4, TAK1, and IKKα/β/γ were changed slightly after vaccination, perhaps due to these genes play their roles through the changes in phosphorylation levels [44].NLRP1 or NLRP3 can bind with their adaptor proteins ASC and capase-1 to form the inflammasome complex after recognizing corresponding ligands (LPS, Muramyl Dipeptide, or bacterial RNA), and then promote cleavage and secretion of biologically active pro-inflammatory cytokine IL-1β [45][46][47].In this study, NLRP1, NLRP3, ASC, and IL-1β were significantly induced by vaccination (Fig 3), which may result in increase of the biologically active IL-1β to recruit inflammatory cells for eliminating pathogens.These results demonstrate that multiple innate immune processes were activated at the early stage of vaccination in large yellow croaker.
Spleen is a major secondary lymphoid organ in the fish, containing a large number of B cells, T cells, neutrophils, and macrophages.The Neutrophils and macrophages as antigen-presenting cells, possess heterogeneous intracellular pathways responsible for generating complexes of MHC class I and class II molecules with peptide antigens, for presentation to T cells [48][49][50][51].Bacterial antigens are internalized into host cells via phagosomes or endosomes.Specifically, the antigenic peptides are degraded in the cytoplasm by proteasome, and then transported to the endoplasmic reticulum by TAP and loaded onto MHC-I molecules with the help of ERAP, TAP, TAPBP, and CLAR [52].Antigenic peptides also can be captured and degraded in endosome and lysosome with the help of CTS, HSP70, and PPT2, and then loaded onto MHC-II molecules [53].In our study, CTS, PPT2, CLAR, TAP, HSP70, and TAPBP were significantly up-regulated after vaccination (Fig 4), while slightly changes were found in the MHC-I and and TCR signaling pathway (C).The expression levels of target gene were expressed relative to those of βactin in each sample using the 2 -ΔΔCT method.Two-tailed Student's t test was used for the significance test analysis, * P < 0.05, ** P < 0.01.doi:10.1371/journal.pone.0170958.g005MHC-II expression levels (Table 3).These results suggest that antigen processing may be activated at the early stage of vaccination, but the antigen presentation was in the process of being started.
Moreover, the expression levels of several genes involved in TCR signaling pathway, including TCRα, TCRβ, related CD molecules (CD3γ, CD4, and CD8β), and downstream signaling molecules (FYN, ZAP-70, PLC-γ1, PDK, and PKCθ), were down-regulated by trivalent bacterial vaccine (Table 3), suggesting that the TCR signaling pathway may be suppressed at the early stage of vaccination.The similar results were also found in other teleost fish [22,34,54,55].In zebrafish, expression of TCR pathway-related genes was found to decrease at 24 h following immunization with attenuated Edwardsiella tarda vaccine [54].In half-smooth tongue sole (Cynoglossus semilaevis), TCR pathway was suppressed at 20 h after Vibrio anguillarum infection [55].Expression levels of TCRα, β, and δ in tilapia (Oreochromis niloticus) were down-regulated by Streptococcus iniae at 24 h post-infection [22].In our previous study, TCR pathway was also suppressed at 24 h after Aeromonas hydrophila infection in large yellow croaker [34].However, the mechanism by which TCR pathway was suppressed at the early period of immunization with bacterial vaccine or bacterial infection remains unknown.
In conclusion, transcriptome analysis was utilized to investigate the early immune response of large yellow croaker to bacterial vaccine immunization.The most of strongly up-regulated genes are innate defense molecules, bacterial ligand-depending PRRs, inflammatory factors, and antigen processing related molecules, indicating that multiple innate immune processes and acquired immune initiation were activated at the early stage of vaccination in large yellow croaker.However, TCR signaling pathway was suppressed after vaccination.The findings observed in this study are similar to immune responses of large yellow croaker infected by live A. hydrophila, implying that inactivated vaccine has the similar effects on induction of host defense by bacterial infection.These results comprehensively reveal the early immune response in the large yellow croaker induced by trivalent vaccine and provide valuable information for developing a highly immunogenic vaccine in teleosts.

Fig 1 .
Fig 1. Statistic of differentially expressed genes (DEGs) at different time points after vaccination.The comparisons of gene expression difference between control (0 h) and each time point after vaccination (24 and 72 h) were performed based on FPKM values.The significant DEGs are defined as log 2 Ratio ! 1 and FDR 0.001.doi:10.1371/journal.pone.0170958.g001

Fig 2 .
Fig 2. Enrichment of Gene Ontology classifications for the DEGs.The DEGs were annotated into three major functional categories: biological process, cellular component, and molecular function.The X-axis indicates the names of GO terms (2nd level) and the Y-axis indicates the number of DEGs.A, Gene Ontology classifications of DEGs at 24 h after vaccination; B, Gene Ontology classifications of DEGs at 72 h after vaccination.doi:10.1371/journal.pone.0170958.g002

Fig 3 .
Fig 3. Pattern-recognition receptors and putative signaling pathways.The predicted pathway was constructed based on the references [43, 56]and the knowledge of mammalian species.Red background indicates significantly up-regulated expression, and yellow background indicates no significant differential expression.The arrows represent promotion, the solid lines indicate direct relationships between genes, and the dashed lines indicate that more than one step is involved in the process.After binding to flagellin, TLR5M triggers sequential recruitment of MyD88, members of IRAK (interleukin-1 receptor-associated kinase) family and TRAF6 that resulted in activation of nuclear factor NF-κB and AP-1, thus inducing production of inflammatory cytokines such as IL-1β, IL-6, IL-12, and TNF-α.After recognizing corresponding ligands (LPS, Muramyl Dipeptide, or bacterial RNA), NLRP1 or NLRP3 can bind with their adaptor proteins ASC and capase-1 to form the inflammasome complex, and the latter promotes cleavage and secretion of biologically active pro-inflammatory cytokine IL-1β.doi:10.1371/journal.pone.0170958.g003

Fig 4 .Fig 5 .
Fig 4. Putative antigen processing and presentation pathway.The predicted antigen processing and presentation pathway of large yellow croaker was constructed based on the knowledge of mammalian species.Red background indicates significantly up-regulated expression, and yellow background indicates no significant differential expression.The arrows represent promotion, the solid lines indicate direct relationships between genes, and the dashed lines indicate that more than one step is involved in the process.Bacterial antigens are internalized into host cells via phagosomes or endosomes.Specifically, the antigenic peptides are degraded in the cytoplasm by proteasome, and then transported to the endoplasmic reticulum by TAP and loaded onto MHC-I molecules with the help of ERAP, TAP, TAPBP, and CLAR.Antigenic peptides can also be captured and degraded in endosome and lysosome with the help of CTS, HSP70, and PPT2, and then loaded onto MHC-II molecules.doi:10.1371/journal.pone.0170958.g004