key: cord-329618-kywhulpc
authors: Xu, Cheng; Evensen, Øystein; Munang’andu, Hetron Mweemba
title: A de novo transcriptome analysis shows that modulation of the JAK-STAT signaling pathway by salmonid alphavirus subtype 3 favors virus replication in macrophage/dendritic-like TO-cells
date: 2016-05-23
journal: BMC Genomics
DOI: 10.1186/s12864-016-2739-6
sha: 
doc_id: 329618
cord_uid: kywhulpc

BACKGROUND: The Janus kinase (Jak) and signaling transducer activator of transcription (Stat) pathway mediates the signaling of genes required for cellular development and homeostasis. To elucidate the effect of type I IFN on the Jak/stat pathway in salmonid alphavirus subtype 3 (SAV3) infected macrophage/dendritic like TO-cells derived from Atlantic salmon (Salmo salar L) headkidney leukocytes, we used a differential transcriptome analysis by RNA-seq and the Kyoto encyclopedia of genes and genomes (KEGGs) pathway analysis to generate a repertoire of de novo assembled genes from type I IFN treated and non-treated TO-cells infected with SAV3. RESULTS: Concurrent SAV3 infection with type I IFN treatment of TO-cells suppressed SAV3 structural protein (SP) expression by 2log(10) at 2 days post infection compared to SAV3 infection without IFN treatment which paved way to evaluating the impact of type I IFN on expression of Jak/stat pathway genes in SAV3 infected TO-cells. In the absence of type I IFN treatment, SAV3 downregulated several Jak/stat pathway genes that included type I and II receptor genes, Jak2, tyrosine kinase 2 (Tyk2), Stat3 and Stat5 pointing to possible failure to activate the Jak/stat signaling pathway and inhibition of signal transducers caused by SAV3 infection. Although the suppressor of cytokine signaling (SOCS) genes 1 and 3 were upregulated in the IFN treated cells, only SOCS3 was downregulated in the SAV3 infected cells which points to inhibition of SOCS3 by SAV3 infection in TO-cells. CONCLUSION: Data presented in this study shows that SAV3 infection downregulates several genes of the Jak/stat pathway, which could be an immune evasion strategy, used to block the transcription of antiviral genes that would interfere with SAV3 replication in TO-cells. Overall, we have shown that combining de novo assembly with pathway based transcriptome analyses provides a contextual approach to understanding the molecular networks of genes that form the Jak/stat pathway in TO-cells infected by SAV3.

The Jak/stat pathway has been shown to contribute positively to innate immune responses against viral infections from drosophila to vertebrates [1] [2] [3] . Viruses have evolved counter measures that block different components of the Jak/stat pathway to prevent production of antiviral compounds [3, 4] . While very little information is available for salmon pancreas disease virus (SPDV also referred to as salmon alphavirus -SAV) causing pancreas disease of Atlantic salmon (Salmo salar L) and trout (Oncorhynchus mykiss), it has been shown that other alphaviruses like chikungunya virus (ChikV) inhibit the stimulation type I and II IFN by blocking Jak/stat signaling, likely through NSP2 [1, 3] . Lin et al. [5] showed that Japanese encephalitis virus (JEV) block IFN induced Jak/stat signaling through a protein tyrosine phosphatase mediated mechanism. As pointed out by Fujii et al. [6] that strategies of counteracting Jak/stat signaling vary among viruses which include suppressing the phosphorylation of the Jaks, inhibition of nuclear translocation of activated Stats, degradation of different Jak/stat components as well as inducing the expression of negative regulators such as the SOCS proteins by host cells. This shows that viruses alter the Jak/stat signaling at different stages of the pathway.

Janus kinase (Jak) and signaling transducer and activator of transcription (Stat) is an intracellular signaling pathway that mediates the effects of a large number of cytokines and growth factors [7] . It is a pleiotropic cascade used to transduce a network of signals for the development and homeostasis of cells from vertebrates to insects. Jak activation stimulates cell proliferation, differentiation, migration and apoptosis. In general, the Jak/stat pathway is based on a few principle components involving a variety of ligands and their receptors [8] . Intracellular activation of the Jaks takes places when ligand binding induces dimerization of receptor subunits on cell surfaces [8] . For signal transduction through dimerization, the cytoplasmic domain of the receptor subunits must be associated with Jak tyrosine kinases, which juxtaposes the Jaks to allow for their transphosphorylation. The transphophorylated Jaks phosphorylate additional targets that include the Stats. Stats are latent transcriptional factors found in the cytoplasm in an inactive state until they are activated. Phosphorylated Stats migrate to the nucleus where they modulate the transcription of different target genes. In general, the Jak/stat pathway provides a direct mechanism in which extracellular signals are translated into transcriptional responses that modulate the expression of different target genes.

Apart from Stats, other effector molecules that contribute to Jak/stat signaling include the signal-transducing adaptor molecules (STAMs) that facilitate the transcriptional activity of specific target genes such as transcriptional regulator Myc (Myc) [9] and the stat-interacting protein (StlP) which serves as a scaffold to phosphorylate the Stats by the Jaks [8] . On the other hand, three major classes of negative regulators have been identified which include the suppressor of cytokine signaling (SOCS), protein inhibitors of activated stats (PIAS) and protein tyrosine phosphates (PTPs) [10] . The SOCS family serves as a negative feedback loop on the Jak/stat circuit in which activated Stats stimulate the transcription of SOCS genes that bind to phosphorylated Jaks and their receptors to turn off the signaling pathway [11, 12] . The PTPs includes the SHP-2 phosphatases, which has two SH2 (src-homology-2) domains that bind to the Jaks to initiate the dephosphorylation of the Jak and Stat proteins [13] . The PIAS proteins bind to activated Stat dimers to prevent their signal transducing activities although the exact mechanism in which this is carried out has not been established [13, 14] . Put together, the Jak/stat cascade is a tightly controlled pathway regulated by different effector proteins and negative regulators that prevent overstimulation of the paracrine pathway.

Therefore, to identify genes that alter Jak/stat signaling induced by viral infections calls for a comprehensive profiling of genes induced by virus replication during infection progression in host cells. In the present study we used RNA-seq and put together a de novo assembly of genes belonging to the Jak/stat pathway induced by salmonid alphavirus subtype-3 (SAV-3) infection in TOcells. TO-cells are a continuous cell-line derived from Atlantic salmon headkidney leukocytes characterized to possess dendritic cell-like properties [15] . SAV3 is a potent inducer of type I IFN responses in TO-cells while at the same time, the virus replicates vividly under high IFN I levels [15] [16] [17] . Pretreatment of permissive cells with recombinant salmon IFN I will protect the cells against CPE but the effect of this treatment will depend on time of treatment relative to infection [17] . The mechanisms employed by SAV3 to circumvent or block the antiviral responses of the host cells are not known. Hence, in this study we wanted to identify the profile of genes linked to the Jak/stat pathway induced by SAV3 infection in TO-cells and to compare their expression patterns in cells simultaneously infected with SAV3 and treated with type I IFN to SAV3-only infected cells. The aim was to understand any modulation of the antiviral responses elicited by SAV3 infection compared to partly protected (IFN I/SAV3) cells. To do this, we used the Kyoto encyclopedia of genes and genomes (KEGG) database analysis to generate pathways of differentially expressed genes generated by de novo assembly using RNA-seq from SAV3 infected TO-cells. In general, we anticipate that data presented here shall serve as a roadmap to elucidating the cellular networks linked to the Jak/stat pathway used by SAV3 to evade the host defensins in infected fish.

TO-cells derived from Atlantic salmon (Salmo salar L) head kidney leukocytes characterized to possess macrophage/dendritic-like properties were use in this study [15, 16] . The TO-cells, were propagated in HMEM (Eagle's minimal essential medium [MEM] with Hanks' balanced salt solution [BSS] ) supplemented with L-glutamine, MEM nonessential amino acids, gentamicin sulfate, 10 % FBS and were incubated at 20°C. The Salmon alphavirus subtype 3 (Genebank accession JQ799139) used to inoculate the TO-cells has previously been described [18] . One batch of TO-cells was treated with 500 ng/ml of Atlantic salmon Type I IFN in triplicates and was simultaneously infected with SAV3. Another batch of TO-cells was only infected with SAV-3, without type I IFN treatment. The virus was inoculated at an MOI 1 when the cells were 80 % confluent in both groups. Thereafter, both the SAV-3 infected cells with and without IFN treatment were incubated at 15°C in maintenance media using HMEM growth media supplemented with 2 % FBS. The mock group was only treated with maintenance media without SAV3 infection and no type I IFN treatment. Cells were harvested after 48 h and were used for RNA extraction used for RNA-seq analysis. Another set of three independent replicates subjected to the same treatment was used for qRT-PCR analysis. All studies in TOcells were carried out in triplicates.

Extraction of total RNA was carried out using the RNeasy mini Kit with on-column DNase treatment according to the manufacturer's instructions (Qiagen, Hilden, Germany). The concentration and quality of RNA was analyzed using a Nanodrop ND1000 (Nanodrop Technologies, Wilmington, USA) and Agilent 2100 Bioanalyzer (Agilent Technologies, USA), respectively.

Triplicates of equal quantities of total RNA from IFNtreated and non-treated-TO-cells infected with SAV3 together with Mock cells were mixed to prepare the pooled RNA sample for each group. The pooled samples were treated with DNase I for the degradation of any possible DNA contamination followed by enrichment using oligo(dT) magnetic beads. Thereafter, the mRNA was fragmented into short fragments of approximately 200 bp by mixing with the fragmentation buffer. This was followed by first strand cDNA synthesis using random hexamer-primer, which was followed by adding a buffer containing dNTPs, RNase H and DNA polymerase I to synthesize the second strand. The generated double strand cDNA was then purified with magnetic beads and the end reparation and 3'-end single nucleotide A (adenine) addition was then performed. Thereafter, sequence adaptors were ligated to the fragments, which were later enriched by PCR amplification. The Agilent 2100 Bioanaylzer and ABI StepOnePlus Real-Time PCR System (Biorad.com) were used to qualify and quantify the sample library during quality check (QC step). After the QC step, library products were ready for RNA-sequencing using Illumina HiSeqTM 2000, BGI-Hong Kong. For RNA-seq analysis, clean reads were obtained after removal of adaptor sequences, reads having greater than 10 % unknown bases as well as removal of low quality bases (base with quality value ≤ 5) greater than 50 % in a read.

Raw reads generated by illumina sequencing were used to generate a library of clean reads after quality check, which was used for transcriptome de novo assembly using the Trinity software (http://trinityrnaseq.github.io [19] . Thereafter, assembled unigenes were used for annotation for the classification of gene functions by searching different protein databasesusing BlastX version 2.2.23 as previously described [18] . Data from BlastX was used to extract the coding regions (CDS) from unigene sequences and translate them into peptide sequences. As for unigenes without hits in BlastX, the ESTScan was used to predict their CDS and to decide their sequence direction. Unigenes having NR annotation were further analyzed with Blast2go (https://www.blast2go.com/) to obtain their gene ontology (GO) annotations for classification according to GO functions using the Web Gene Ontology (WEGO) annotation software.

Identification of differentially expressed genes (DEGs) was carried out as previously described [18] . In this study, DEGs were generated based on comparing the RPKM mapped reads from TO cells treated with type I IFN simultaneously infected with SAV3 versus mock TO-cells that were designated as TO-VS-IFN/SAV3, while the second comparison was based on TO-cell infected with SAV-3 only, without IFN treatment, versus mock TO-cells that were designated as TO-VS-SAV3. Only genes with a threshold of false discovery rate (FDR) <0.001 and an absolute value log 2 ratio > 1 were considered as differentially expressed. Thereafter, all DEGs were assigned KEGG orthologue (KO) identifiers which were later used for pathway analysis using the KEGG database (http://www.genome.jp/kegg/). In addition, we used the TIGR gene indices clustering tools (TGICL) (http://www.tigr.org/tdb/tgi) [19, 20] software to remove redundancy.

To confirm the validity of the DEGs generated by RNA-Seq, 13 DEGs obtained from the TO-VS-IFN/SAV3 and TO-VS-SAV3 were randomly selected for the qRT-PCR validation test. Gene expression of three independent replicates of TO-cells simultaneously treated with type I IFN and infected with SAV-3 (TO-VS-IFN/SAV3) was compared with none-type I IFN treated cells infected with SAV3 (TO-VS-SAV3) relative to gene expression of untreated TO-cells before normalized to β-actin. Cells used for independent replicate biological testing were independent of the cells used for RNA-seq data analysis. The qRT-PCR was carried out as previously described [18] using 100 ng total RNA as a template for each gene in a master-Mix final volume of 25 μl according to manufacturer's recommendations (Roche Diagnostics, Mannheim, Germany). Primers used in the study are shown in Table 1 and the primer-master-mix solutions were first incubated for reverse transcription at 50°C for 10 min followed by PCR initial activation at 95°C for 5 min and 40 amplification cycles (10 s at 95°C and 30 s at 60°C). The specificity of each PCR products was confirmed by melting-curve analysis and agarose gel electrophoresis. The delta C t method was used to calculate the fold-increase in gene expression relative to the β-actin control group [21] . All quantifications were normalized to β-actin control group. The qRT-PCR was done thrice, with each repeat study having three independent replicates of IFN treated (TO-VS-IFN/SAV3) and non-treated cells (TO-VS-SAV3). The mean of each gene from the three independent replicates after three repeats was used to compute the correlation coefficient test between RNA-seq and qRT-PCR data.

In order to determine the effect of type I IFN treatment on SAV3 replication in TO-cells, virus quantification was determined after 48 hrs replication at 15°C for the IFNα treated and non-treated groups using a guided sequence format for the structural proteins (SPs) and non-structural proteins (NSPs) of SAV3 (GeneBank accession nos. AFJ15605.1 and AF15604.1) based on the RNA-seq data. The Fragment per kilobase per million (FPRM) method was used to read the number of fragments of NSPs and SPs expressed by SAV3 after 48 hrs replication at 15°C in the IFNα treated and non-treated cells. In order to verify the validity of the virus quantification data generated by RNA-Seq, quantitative RT-PCR (qRT-PCR) was used to determine the SAV3 replication level in the IFNα treated and non-treated cells using the E2 SP expressed during SAV3 replication. The primer sequences used for E2 quantification are shown in Table 1 and the qRT-PCR method used has previously been described by Xu et al. [17] .

First we documented the effect of IFN I treatment on the SAV3 SP (structural protein) expression. We found that the non-IFNα treated cells had 10 2 higher expression levels of SAV3 SPs than the IFN I treated cells (Table 2 ). This shows that simultaneous treatment with IFN I reduced the level of synthesis of mRNA SPs from subgenomic viral RNA. Similarly, the NSP mRNAs (non-structural proteins mRNAs) expressed in the non-IFNα treated cells were 10 0.8 higher than in the IFNα treated group (IFN/SAV3). The log2 ratio of the IFN-plusSAV3/SAV3 for the SPs and NSPs showed a downregulation of -7.4092 and -2.3219, respectively, showing that type I IFN added concurrently with SAV-3 infection significantly reduced virus replication in TO cells ( Table 2) . Consistent with RNA-seq data shown in Table 2 , Fig. 1 shows a 20 fold-reduction in SAV3 replication in the IFN I treated cells (IFN/SAV3) relative to the untreated cells based on the E2 SPs detected by qPCR using the Cp-value test. Thus, qPCR confirms the validity of our RNA-seq data analysis by showing that a downregulation of SAV3 SPs in IFN I treated cells (IFN/ SAV3) versus SAV3 infected cells only.

Based on data generated by illumina sequencing, we constructed three libraries that yielded a total of 23 331 

Bcl-X-R AGTCCTAGATACCTCCCTCC 788, 24 608 338 and 23 421 631 raw reads for the IFNtreated and non-IFN treated TO-cells infected with SAV3 together with mock TO-cells, respectively. The total unique match in total mapped reads was >90 % after removal of redundancy using the TGICL software (http://www.tigr.org/tdb/tgi) suggesting low redundancy. After quality check for clean reads and filtration [18] , 20,115 unigenes from the TO-VS-IFN/SAV3 and TO-SV-SAV3 groups assigned to the KEGG orthologues (KO) were identified using the KEGG database (http:// www.genome.ad.jp/kegg/), which has been widely used to generate pathway networks for model organisms such as the Zebra fish (Danio rerio) [22, 23] and mouse (Mus musculus) [24] as well as several non-model organisms [25, 26] . The profile of genes identified using KO-based on the RPKM method in which only genes with an FDR < 0.001 and an absolute value log 2 ratio > 1 were considered differentially expressed showed that the TO-VS-IFN/SAV3 had 1629 DEGs whose KOs were mapped to 237 pathways of which 49 DEGs belonged the Jak/stat pathway (Table 3) . On the other hand, the TO-VS-SAV3 had a total of 9135 DEGs assigned to KO-identifiers, which was 7506 DEGs more than the TO-VS-IFN/SAV3. As a result, the TO-VS-SAV3 had 252 pathways identified by the KEGG database being 15 pathways more than those identified in the TO-VS-IFN/SAV3 (Table 3 ). In terms of Jak/stat genes, the TO-VS-SAV3 had a total of 135 DEGs expressed being 86 DEGs more than the TO-VS-IFN/SAV3. In general, the TO-VS-SAV3 had more DEGs and pathways identified using the KEGG database analysis than the TO-VS-IFN/SAV3 although the TO-VS-IFN/SAV3 had a higher enrichment Qvalue (1.136E-06, p value = 8.07E-08) than the TO-VS-SAV3 (5.899E-01, p value = 8.07E-01) ( Table 3) .

To gain more insights on the impact of type I IFN treatment on the Jak/stat pathway in TO-cells infected with SAV-3, we compared the DEGs expressed from the TO-VS-IFN/SAV3 with DEGs from the TO-VS-SAV3. Figure 2 shows that of the 135 DEGs expressed in the TO-VS-SAV3, 119 DEGs were downregulated while 16 were upregulated. In the TO-VS-IFN/SAV3 group 40 DEGs were upregulated while 9 DEGs were downregulated. Put together, these findings show that the proportion of upregulated DEGs in the TO-VS-IFN/SAV3 (81.6 %, n = 49) was higher than in the TO-SV-SAV3 group (11.9 %, n = 135) while the proportion of downregulated DEGs for the TO-SV-SAV3 (88.1 %, n = 135) was higher than the TO-VS-IFN/SAV3 group (18.4 %, n = 49). This is strongly suggesting that type I IFNs upregulate several genes of the Jak/stat pathway that were downregulated by SAV3 infection alone in TO-cells. Hence, to further elucidate the type of genes suppressed by SAV3 infection in the absence of type I IFN treatment, we broadly classified the DEGs in the Jak/stat pathways for both the TO-VS-IFN/SAV3 and TO-VS-SAV3 into type I and II receptor genes, Jaks, Stats and negative regulator genes as shown below.

Differentially expressed genes belonging to receptor ligands that activate the Jak/stat pathway were classified into type I and II receptor-genes. As shown in Fig. 2 and Table 3 , the type I receptor family was subdivided into the gamma chain (γC), gp130 and single chain (homodimer) receptor genes while the type II receptor family comprised of type I and II IFN receptor genes. The major difference Table 4 ). The type II receptor genes were subdivided into two groups namely the interleukin (IL) and IFN receptor genes. As shown in Table 3 , the IL-receptor genes comprised of IL-10RB and IL-20RA of which both were upregulated in the TO-VS-IFN/SAV3 and TO-VS-SAV3. The IFN-receptor genes comprised of IFN-α/β and IFNγ receptor genes ( Table 4 ). The IFN-α/β receptor genes comprised of IFNR1 and IFNR2 which were upregulated in both the TO-VS-SAV3 and TO-VS-IFN/SAV3 while the IFNγ receptor genes comprised of IFNGR1 and IFNGR2 that were only expressed in the TO-VS-IFN/ SAV3 in which both were upregulated (Table 4 ). Fig. 4 shows a schematic layout of the common IFN-receptor signaling pathways based on mammalian model species [8, 27] of which applying data presented in Table 3 shows that the IFN-receptor-genes expressed by the TO-VS-IFN/SAV3 were linked to both pathways A and B while the IFN-receptor-genes expressed by the TO-VS-SAV3 was only linked to pathway A.

Three genes belonging to the Jaks were differentially expressed namely Jak1, Jak2 and tyrosine kinase (Tyk2) (Fig. 3 , Table 5 ). All three Jaks were differentially expressed in the TO-VS-SAV3 while only Jak1 was differentially expressed in the TO-VS-IFN/SAV3. Jak1 was upregulated in both the TO-VS-SAV3 and TO-VS-IFN/ SAV3 while Jak2 and Tyk2 were downregulated in the TO-VS-SAV3 (Table 5 ). Network pathways linking the Jaks to other genes in the Jak/stat pathway are shown in Fig. 5A for the TO-VS-IFN/SAV3 and in Fig. 4B for TO-VS-SAV3. Note that Fig. 5A shows upregulation of all Jaks being in conformity with observations made in Table 5 in which only Jak1 was differentially expressed in the TO-VS-IFN/SAV3 and it was upregulated. On the other hand, Fig. 5B shows a mixed representation of upregulated and downregulated Jak-genes, which is in agreement with observations shown in Table 4 in which Jak1 was upregulated whiles Jak2 and Tyk2 were downregulated.

Five Stat genes were differentially expressed in this study namely Stats 1, 2, 3, 5B and 6 ( Table 5 ). As shown in Fig. 2A , Stats 1, 2, 3 and 5B were differentially expressed in the TO-VS-SAV3 of which Stats 1 and 2 were upregulated while Stats 3 and 5B were downregulated (Table 5 ). In the TO-VS-IFN/SAV3, Stats 1, 2 and 6 were upregulated (Table 5 , Fig. 3B ). The network links of the Stats to other genes in the Jak/stat pathway for the TO-VS-SAV3 and TO-VS-IFN/SAV3 are shown in Fig. 5A and B, respectively. Note that all Stats in the TO-VS-IFN/SAV3 were upregulated in Fig. 5A conforming to observations shown in Table 5 in which all Stats were upregulated while Fig. 5B shows a mixed representation of upregulated and downregulated Stats for the TO-VS-SAV3 being in line with observation shown in Table 4 in which Stats 1 and 2 were upregulated while Stats 3 and 5B were downregulated.

The negative regulatory genes for the Jak/stat pathway differentially expressed in this study were classified into Table 5 ). Table 5 shows that SOCS1 and SOCS3 were both upregulated in the TO-VS-IFN/SAV3 while SOCS3 was downregulated and SOCS1 was upregulated in the TO-VS-SAV3. Figure 5A and B shows network pathways linking the SOCS to other genes in the Jak/stat pathway in which all SOCSs were upregulated in the TO-VS-IFN/SAV3 (Fig. 5A) , which conforms to observations shown in Table 4 . The network pathway in the TO-VS-SAV3 shows a mixed representation of upregulated and downregulated SOCSs (Fig. 5B) , which conforms to observations shown in Table 5 where SOCS1 was upregulated while SOCS3 was downregulated in the TO-VS-SAV3. As for PIAS, it was downregulated in the TO-VS-SAV3 and upregulated in the TO-VS-IFN/SAV3 ( Fig. 5A and B ).

Genes that regulate apoptosis and cell-cycle activities were broadly classified into the Stat, MAPK and PI3K-AKT pathways as shown in Table 6 based on network pathways shown in Fig. 5A and B. In the Stat-pathway, the IFN regulatory factor 9 (IRF), creb binding protein (CBP) and threonine-protein kinase pim-1 (Pim-1) were upregulated in the TO-VS-IFN/SAV3 while in the TO-VS-SAV3, the signal transducer and adaptor molecule 2 (STAM2), CBP, Pim-1 and Pim-3 were downregulated (Table 6 ). In the MAPK pathway shown in Fig. 5A and B, the transcriptional regulator Myc-2 (Myc-2) and cyclin D (CycD1) were upregulated in the TO-VS-IFN/ SAV3 while the Src homology region 2 domaincontaining phosphatase-2 (SHP2) which is represented by the protein tyrosine phosphatase non-receptor 11 (PTPN11) in Table 5 , son of sevenless (SOS), Myc-2, CycD, B-cell lymphoma extra-large (BclXL) and Sprouty-related EVH1 domain-containing protein 1 (Spred) were all downregulated in the TO-VS-SAV3 (Table 6 and Fig. 5B ). The phosphatidylinositol-4,5bisphosphate 3-kinase (PI3K) and V-akt viral oncogene 2-like protein (AKT2) that regulate the PI3K-AKTpathway for apoptosis were downregulated in TO-VS-SAV3 (Table 6 ) while only PI3K was downregulated in the TO-VS-IFN/SAV3. Network pathways linking the apoptosis and cell cycle genes for the TO-VS-SAV3 and TO-VS-IFN/SAV3 are shown in Fig. 5A and B. In general, all apoptosis and cell proliferation genes differentially expressed in the TO-VS-SAV3 were downregulated while the majority differentially expressed in the TO-VS-IFN/SAV3 was upregulated indicating that in the absence of type I IFN treatment SAV3 suppressed the expression of most of anti-apoptosis and cell proliferation genes.

Validation of RNA-Seq data by qRT-PCR Figure 6 shows validation of RNA-seq data based correlation coefficient test between RNA-seq and qRT-PCR data. Our findings show a significantly high linear correlation (r 2 = 0.8832, P < 0.0001) between RNA-seq and qRT-PCR for the TO-VS-IFN/SAV3 and TO-VS-SAV3. Put together, qRT-PCR confirms the validity of our RNA-seq data analysis by showing a high correlation between RNA-seq and qRT-PCR.

De novo and reference gene guided transcriptome assembly Current advances in next generation sequencing show that RNA-seq is poised to enhance our understanding of molecular networks of genes that form the innate and adaptive immune system in teleosts fish. In general, there are two approaches used to assemble a transcriptome from RNA-seq data which is by use of a de novo assembly or a reference gene guided approach. In situations where no reference gene exists as is the case for non-model organisms like Atlantic salmon, a de novo assembly is the only option for sequence assembly by which raw reads from RNA-seq are assembled into contigs without the guidance of a reference gene. Hence, in the present study we used a de novo approach to assemble a repertoire of genes that form the Jak/stat pathway from Atlantic salmon headkidney/macrophage like TOcells with an aim to understand SAV3 modulation of Fig. 4 Shows the schematic presentation of the IFN-α/β and IFNγ receptors linked to the Jaks for the TO-VS-IFN/SAV3 and TO-VS-SAV3. In pathway A, the IFN-α/β binds Jak1 and Tyk2 receptor subunits linked to Stats 1 and 2, followed by a pathway that translocates the phophorylated Stats into nucleus. In pathway B, the IFNγ receptor is shown to bind to Jaks 1 and 2, which are linked to Stat1 and the phosphorylated Stats are shown to translocate to the nucleus. Based on data shown in Table 2 , the TO-VS-IFN/SAV3 uses both pathways A and B using the IFN-α/β and IFNγ ligands while the TO-VS-SAV3 only uses pathway A using the IFN-α/β ligand. This model is adapted from IFN-receptor signaling pathways for the Jak/stat pathway in higher vertebrates [41, 42, 70] responses post infection. Unlike microarray that uses predefined probes, de novo assembly provides a unique opportunity to identify novel genes which could be missed by none-coded micro-RNA probes. After de novo assembly, we used the TGICL software, which is a highly refined protocol used for the analysis of genes and EST sequences needed for transcriptome analysis [20] . Using this software, sequences are first cleaned to remove adaptors, low quality sequences without quality score, vector and adaptor contaminants and redundancies. Thereafter, sequences are searched pairwise and grouped into clusters that are assembled at high stringency to produce tentative concesus (TC) sequences. And as pointed out by Pertea et al. [19] , individual assembly of each cluster to produce TCs has the advantage of producing large, more concensus sequences while eliminating potentially misclustered sequences. The TCs produced by the TGICL software have been used to construct species specific databases and to identify novel genes [20] . For example, several genes such as STAM2, SOS, and PIAS not previously reported from Atlantic salmon were matched to their orthologues in model species such as zebrafish and mammalia in this study suggesting that the Jak/stat pathway induced in macrophage/dendritic cell like TO-cells derived from Atlantic salmon headkidney cells expresses a repertoire of genes comparable to the Jak/stat signaling pathway observed in higher vertebrate species. On the other hand, we used a reference gene guided approach to detect and quantify the E2 SP of SAV3 based on the fact that the genomic sequence of SAV3 has been fully sequenced and characterized in salmonids [17, 28, 29] . Hence, both the de novo and reference gene guided assemblies were used in this study.

Virus replication in TO-cells of simultaneous IFN I/ SAV3 treatment/infection (IFN/SAV3) was two logs (10 2 ) lower than levels detected in SAV3-only infected cells. For the latter, it shows that SAV3 is a potent inducer of IFN signaling pathways in TO-cells post infection and also that the virus had the capacity to replicate in the presence of high levels of type I IFNs and timing of IFN I treatment is crucial [18] . These findings are consistent with observations made for other alphavirus studies in which it has been shown that timing in type I IFN treatment has a significant influence on replication of alphaviruses in infected cells [1, 30, 31] .

RNA-seq is perfect for gene network construction and pathway analysis [32] . And as pointed out by Khatri et al. [33] , pathway analysis has become the first choice to understand the underlying biological mechanisms of DEGs because it reduces the complexicity of RNA-seq data analysis whose vast array of sequence data is often difficult to interpret but instead combining RNA-seq with pathway analysis increases the explanatory power of transcriptome data analysis. In pathway analysis the first step is to generate a transcriptome of DEGs as shown in this study followed by searching for functional networks of genes expressed in the transcriptome based on the fact that there is a tendency for co-regulated genes to have functional correlations in the pathway. Khatri et al. [33] shown that pathway software development has progress through three different generations in the last decade. The first generation of pathway analytical softwares comprise of the over-representation analysis (ORA) softwares that include GenMAPP (http:// www.genmapp.org), GOstat (http://gostat.wehi.edu.au) and GOToolBox (http://genome.crg.es/GOToolBox/). These software establish pathway analysis by putting together over expressed genes that show significant changes in expression based on gene level statistic [34, 35] . The second generation of softwares consist of functional class scoring (FCS) approaches such as the Gene Set Enrichment Analysis (GSEA) software (http://www.broadinstitute.org/ gsea/) and GeneTrail (http://genetrail.bioinf.uni-sb.de) that compute both the statistical significance of individual genes that form a pathway (gene level statistic) and the statistical significance of the entire pathway (pathway level statistic) [36] [37] [38] . The third generation of softwares comprises of pathway topology (PT) based approaches such as the KEGG [39] and Panther DB [40] that do not only encorporate individual gene and pathway level statistics computed by the the ORA and FCS softwares, but they also provide information about gene products that interact with each other in a given pathway, how they interact (e.g. activation, inhibition, etc), and where thay interact (e.g. cytoplasm, nucleus, endoplasmic reticulum, etc). Hence, in this study we used the KEGG pathway analysis software which is a third generation PT based software to identify networks of genes that form the Jak/stat signaling pathway from a de novo assembled transcriptome with the aim to better understand the modulation of genes expressed by TO cells infected with SAV3. Consistent with all TP-based softwares, the KEGG pathways used in this study shows gene level statistical analyses for all DEGs and in terms of pathway level statistic analysis, our findings show that concurrent treatment of TO-cells with type I IFN together with SAV3 infection significantly enriched the Jak/stat pathway from a Qvalue of 5.889E-01 (P-value = 8.07E-01) to 1.136E-06 (Pvalue = 8.07E-08). This was supported by upregulation of several activator and signal transducer genes in the TO-VS-IFN/SAV3 group compared to the TO-VS-SAV3 group. As for genes that interact with each other and their interaction sites, the modulation is complex and below we pinpoint some of the key areas in the pathway where SAV3 infection interferes.

Receptor genes that activate Jak/stat signaling Type I receptor genes generated in response to SAV3 infection and type I IFN treatment were classified into the γC, gp130 and homodimer receptor families based on receptor classification used for higher vertebrates [41, 42] . Wang et al. [43] reviewed the γC receptor system for teleost fish in which they reported of IL-2RB, IL-4RA, IL-13RA and IL-15RA in salmonids which were also differentially expressed in this study. SAV3 infected cells (TO-VS-SAV3) had more type I receptor genes downregulated than type I IFN treated cells (TO-SV-IFN/SAV3) suggesting that SAV3 like other alphaviruses in higher vertebrates, inhibits the activation of several type I receptor genes post infection. The type II receptor family generated in this study was classfied into IL and IFN receptor genes.

As for genes of the IFN receptor genes, IFNR1 and IFNR2 were expressed in both the TO-VS-IFN/SAV3 and TO-VS-SAV3 suggesting that activation of IFN-α/β receptors occurs as consequence of infection as well as IFN I treatment (as would be expected). However, it is important to point out that although our findings show the presence of both IFNR1 and IFNR2, which is in agreement with observations made by Ali et al. [44] who reported the presence of both IFNR1 and IFNR2 based on a transcriptome analysis of rainbow trout spleen. Studies carried out by Sun et al. [45] suggest that Atlantic salmon possess two clusters type I IFN receptor genes located on different chromosomes which allow for a large repertoire of IFN receptors than mammals and zebrafish. Hence, there is need for further studies to compare the IFNR1 and IFNR2 genes detected in this study with the IFN receptor genes repertoire reported by Sun et al. [45] and to determine their functional roles in the Jakstat pathway signaling. The expression of IFNGR1 and IFNGR2 was only found in the TO-VS-IFN/SAV3 cells. This suggests that SAV3 infection does not activate the IFNγ-receptors. The importance of IFN-γ in eliciting an antiviral effect in vitro remains non-conclusive. We have shown that IFN-γ exhibits marginal antiviral effect against SAV3 in vitro [16] . In contrast Sun et al. [46] suggested that recombinant trout IFN-γ had antiviral effect, although they also concluded that the observed effect of IFN-γ is partly dependent on IFN-α induction [47] . Additional studies are needed to understand the importance of IFN-γ in relation to SAV3 infection.

Four Jak-family members have been identified namely Jak1-Jak3 and Tyk2 in mammalia [48] and fish [49] [50] [51] [52] of which Jak1, Jak2 and Tyk2 transcripts generated in this study were matched to their mammalian and fish orthologues. The mechanisms underlying activation, regulation and pleiotropic signaling in fish have not been studied in detail. Here we find that Jak1 was upregulated in SAV3 infected TO cells, with or without IFN-I treatment, while Jak2 and Tyk2 were only differentially expressed in the TO-VS-SAV3 groups, and importantly, they were both downregulated. Thus, these findings suggest that SAV3 infection activates Jak1 signaling while it suppresses activation of Tyk2 and Jak2. These observations are in line with what is seen for alphavirus infections in higher vertebrates in which blockage of Tyk2 was shown for Japanese encephalitis virus (JEV) [5] , while Dengue virus type 2 was found to antagonize the Jak/stat pathway by downregulating Tyk2-stat signaling in human dendritic cells [53] . Lin et al. [5] showed that JEV infection blocked tyrosine phosphorylation of IFN associated Tyk2 without affecting the expression of IFNα/β on the cell surface. We found that IFNα was upregulated when Tyk2 was downregulated, and while it is possible that SAV3 uses similar mechanisms to what is seen for other alphaviruses to block the phosphorylation of Tyk2 thereby inhibiting Jak/stat signaling in TO-cells, this remains to be shown.

Type I IFN mediates tyrosine phosphorylation by binding to receptors associated with Jak1 and Tyk2 which activates Stats 1 and 2, while IFNγ mediates its signaling through Jaks 1 and 2 to activate Stat1 in mammalian cells [27, 42] . The TO-SV-IFN/SAV3 showed upregulation of IFNγ receptor genes, corresponding to the observed upregulation of Jak1 and Stat1, while downregulation of Tyk2 in SAV3 infected cells showed a corresponding downregulation of Stat3 and Stat5. Hence, it is possible that mechanisms in salmon and human cells are comparable, however, to identify the different Jaks that activate specific Stats in fish cells, there is need for detailed investigations. Studies in higher vertebrates show that alphaviruses block the Jak/stat signaling pathway by inhibiting the phosphorylation of different Stats to prevent the transcription of antiviral genes in infected cells. For example, JEV infections block tyrosine phosporylation of Stats 1, 2 and 3 [46] while CHIKV nsPs were potent inhibitors of IFN-induced Jak/stat signaling by blocking stat1 phosphorylation of the Jak-receptors either by type I or II IFNs [1] . The observed downregulation of several type I and II receptor genes in SAV3 infected cells could resulted in inhibition of Jak2 and Tyk2 phosphorylation, which would align with the downregulation of most Stat genes. Although these data suggest that SAV3 could be using mechanisms similar to those used by mammalian alphaviruses to block the Jak/stat pathway, there is need for detailed investigations to consolidate these observations in TO-cells.

The negative regulatory genes linked to inhibition of the Jak/STAT pathway include genes that form a negative feedback loop comprising of members of the SOCS family and genes that directly inhibit the Jak/stat pathway such as PTP1b, SHP and PIAS. In salmonids, several SOCS genes have been cloned and characterized [53, 54] and studies carried out by Skjesol et al. [55] show that SOCS1 supresseses type I and II signaling of the Jak/stat pathway in Atlantic salmon. In mammalia, SOCS1 and SOCS3 but not SOCS2 inhibit IFN-mediated antiviral activities [52] . In this study, SOCS1 was upregulated in both TO-VS-IFN/SAV3 and TO-VS-SAV3 while SOCS3 was only upregulated in the TO-VS-IFN/SAV3 and downregulated in the TO-VS-SAV3 suggesting that SAV3 replication could have suppressed the expression of SOCS3.

Among the genes that directly inhibit Jak signaling, only PIAS was differentially expressed and it was upregulated in the TO-VS-IFN/SAV3 and downregulated in the TO-VS-SAV3, suggesting that SAV3 suppressed the expression of PIAS in TO-cells. Although PIAS has not been characterized in salmonids, its expression at transcript level in this study suggests that it exists in Atlantic salmon and that it likely to play an important role in regulating the Jak/stat pathway in TO-cells. Overall, these negative regulators were all downregulated in the TO-VS-SAV3 cells, except for SOCS1, suggesting that SAV-3 upregulates SOCS1 to prolong its survival in infected cells.

Here we find that apoptosis related genes Pim-1, c-Myc and CycD are upregulated in the TO-VS-IFN/SAV3 and downregulated in the TO-VS-SAV3 cells. Apoptosis following virus infection is an important defense mechanisms given that cell autodestruction is one of the safest ways of limiting virus spread to neighboring cells in an infected host [56] . However, viruses have evolved defense mechanisms that evade or delay apoptosis [56] [57] [58] by altering the expression of anti-apoptosis target genes that regulate cell growth and proliferation [59] [60] [61] [62] . Our finding can be interpreted as a survival strategy of SAV3 to prolong its replication cycle in TO-cells. Further, it is also possible that the expression of PI3K and AKT could be regulating SAV3 induced apoptosis in TO-cells as shown for other alphavirus where PI3K-AKT pathway regulates cell survival by preventing caspase-8 activation to prevent apoptosis in virus infected cells [63] [64] [65] . Mastrangelo et al. [66, 67] have shown that Bcl-2 and BcI-xL limit apoptosis upon infection with alphaviruses while CycD has been shown to exert the same effects in neuronal cells infected by Sindbis virus [68] . Similarly, Pim-1 acts synergestically with c-Myc [61] , to inhibit apoptosis related mitochondrial dysfunction through the bcl-2 dependent pathway [62, 69] . It remains to be shown whether downregulation of apoptosis target genes in the TO-VS-SAV3 cells could be an evasion strategy to shutdown cell-death or to slowdown apoptosis to prolong replication time in infected TO-cells.

This study provides a global overview of genes that form the Jak/stat pathway induced by SAV3 infection and the combined IFN I/SAV3 infection allow us to pinpoint modulating effects of SAV3 infection on the Jak/stat signaling pathway. It is important to note that the data presented here only provides a repertoire of downregulated and upregulated genes interlinked in a de novo assembled network pathway analysis, and does not show the exact mechanism by which the modulatory processes take place. Thus a lot of work remains to be done to elucidate the exact interaction points and modulation of the Jak/stat encoded proteins. genes generated by de novo assembly in this study. Further, the study shows that combining de novo assembly with pathway based analysis increases the explanatory power of transcriptome data analysis by putting together a network of genes that regulate hostpathogen interaction processes at cellular level. We anticipate that the Jak/stat pathway genes put forth in this study shall serve as a roadmap for elucidating molecular mechanisms used by SAV3 infection to evade host defensins at cellular level.

The RNA-sequencing data generated in this study has been deposited in the National Center for Biotechnology Information (NCBI) Gene Expression Omnibus (GEO) database accession number GSE64095 (www.ncbi.nih.gov/ geo Accession number GSE64095).

Chikungunya Virus Nonstructural Protein 2 Inhibits Type I/II Interferon-Stimulated JAK-STAT Signaling

JAK/STAT Pathway in Drosophila Immunity

Interferons and viruses: an interplay between induction, signalling, antiviral responses and virus countermeasures

STAT-5, Activates the ATM DNA Damage Pathway to Induce HPV 31 Genome Amplification upon Epithelial Differentiation

Blocking of interferon-induced Jak-Stat signaling by Japanese encephalitis virus NS5 through a protein tyrosine phosphatase-mediated mechanism

Molecular mechanisms for suppression of interferon signal transduction pathways caused by viral infections

The Jak/STAT pathway

The JAK/STAT signaling pathway

STAM/EAST/Hbp adapter proteins -integrators of signalling pathways

Negative regulation of cytokine signaling

Suppressors of cytokine signaling (SOCS): negative regulators of signal transduction

Negative regulation of cytokine signaling by the SOCS proteins

Negative regulation of the JAK/STAT pathway

PIASx is a transcriptional co-repressor of signal transducer and activator of transcription 4

A highly phagocytic cell line TO from Atlantic salmon is CD83 positive and M-CSFR negative, indicating a dendritic-like cell type

A salmonid cell line (TO) for production of infectious salmon anaemia virus (ISAV)

Alpha interferon and not gamma interferon inhibits salmonid alphavirus subtype 3 replication in vitro

De Novo assembly and transcriptome analysis of Atlantic salmon macrophage/dendritic-like TO cells following type I IFN treatment and Salmonid alphavirus subtype-3 infection

TIGR Gene Indices clustering tools (TGICL): a software system for fast clustering of large EST datasets

The TIGR Gene Indices: clustering and assembling EST and known genes and integration with eukaryotic genomes

Real-time, fluorescence-based quantitative PCR: a snapshot of current procedures and preferences

Transcriptome profiling of zebrafish infected with Streptococcus suis

Behavioral and neurogenomic transcriptome changes in wild-derived zebrafish with fluoxetine treatment

Differential Hippocampal Gene Expression and Pathway Analysis in an Etiology-Based Mouse Model of Major Depressive Disorder

Transcriptome analysis of chicken kidney tissues following coronavirus avian infectious bronchitis virus infection

De Novo Assembly and Characterization of Two Transcriptomes Reveal Multiple Light-Mediated Functions in the Scallop Eye (Bivalvia: Pectinidae)

Detection and antigenic characterization of salmonid alphavirus isolates from sera obtained from farmed Atlantic salmon, Salmo salar L., and farmed rainbow trout, Oncorhynchus mykiss (Walbaum)

Characterization of untranslated regions of the salmonid alphavirus 3 (SAV3) genome and construction of a SAV3 based replicon

Antiviral Activity of Alpha-Interferon in Sindbis Virus-Infected Cells Is Restored by Anti-E2 Monoclonal-Antibody Treatment

Early Events in Alphavirus Replication Determine the Outcome of Infection

behavior and next-generation RNA sequencing

Ten Years of Pathway Analysis: Current Approaches and Outstanding Challenges

Bioinformatics enrichment tools: paths toward the comprehensive functional analysis of large gene lists

Ontological analysis of gene expression data: current tools, limitations, and open problems

GlobalANCOVA: exploration and assessment of gene group effects

A multivariate extension of the gene set enrichment analysis

Non-linear tests for identifying differentially expressed genes or genetic networks

KEGG: kyoto encyclopedia of genes and genomes

PANTHER: a library of protein families and subfamilies indexed by function

The Jak-STAT pathway

Signaling through the JAK/ STAT pathway, recent advances and future challenges

The gamma-chain cytokine/ receptor system in fish: More ligands and receptors

Characterization of the rainbow trout spleen transcriptome and identification of immune-related genes

Atlantic salmon possesses two clusters of type I interferon receptor genes on different chromosomes, which allows for a larger repertoire of interferon receptors than in zebrafish and mammals

Antiviral activity of salmonid gamma interferon against infectious pancreatic necrosis virus and salmonid alphavirus and its dependency on type I interferon

A direct cross-talk between interferon-gamma and sonic hedgehog signaling that leads to the proliferation of neuronal precursor cells

The Janus kinases (Jaks)

The JAK and STAT family members of the mandarin fish Siniperca chuatsi: Molecular cloning, tissues distribution and immunobiological activity

Genomic organization and characterization of the promoter region of the roundspotted pufferfish (Tetraodon fluviatilis) JAK1 kinase gene

Complete genomic organization and promoter analysis of the round-spotted pufferfish JAK1, JAK2, JAK3, and TYK2 genes

The Atlantic salmon protein tyrosine kinase Tyk2: Molecular cloning, modulation of expression and function

Rainbow trout suppressor of cytokine signalling (SOCS)-1, 2 and 3: Molecular identification, expression and modulation

Identification of suppressor of cytokine signalling (SOCS) 6, 7, 9 and CISH in rainbow trout Oncorhynchus mykiss and analysis of their expression in relation to other known trout SOCS

Functional conservation of suppressors of cytokine signaling proteins between teleosts and mammals: Atlantic salmon SOCS1 binds to JAK/STAT family members and suppresses type I and II IFN signaling

Regulators of apoptosis on the road to persistent alphavirus infection

Viral control of mitochondrial apoptosis

Regulation of apoptosis by viral gene products

c-myc target genes involved in cell growth, apoptosis, and metabolism

Hepatitis B virus X protein induces apoptosis in hepatoma cells through inhibiting Bcl-xL expression

Synergistic roles for Pim-1 and c-Myc in STAT3-mediated cell cycle progression and antiapoptosis

Pim-1: a serine/threonine kinase with a role in cell survival, proliferation, differentiation and tumorigenesis

Inhibition of the PI3K-Akt signaling pathway enhances the sensitivity of Fas-mediated apoptosis in human gastric carcinoma cell line, MKN-45

PI3-kinase regulates survival of chronic lymphocytic leukemia B-cells by preventing caspase-8 activation

Inhibition of phosphatidylinositol 3 '-kinase/AKT signaling promotes apoptosis of primary effusion lymphoma cells

Bcl-2 inhibits apoptosis and extends recombinant protein production in cells infected with Sindbis viral vectors

Part I. Bcl-2 and Bcl-x(L) limit apoptosis upon infection with alphavirus vectors

Cyclin-dependent kinases participate in death of neurons evoked by DNA-damaging agents

The PIM-1 serine kinase prolongs survival and inhibits apoptosis-related mitochondrial dysfunction in part through a bcl-2-dependent pathway

GAS elements: A few nucleotides with a major impact on cytokine-induced gene expression

This study was funded by the Research Council of Norway, project "The Atlantic salmon genome sequence as a tool for precision breeding" project number 226275, and the Targeted disease prophylaxis in European Fish farming project, grant agreement 311993. We thank Jiang Li, at the Beijing Genomics Institute (BGI) for assistance with some bioinformatics analysis.

Authors' contribution CX = Carried out laboratory experiments and data analysis; HMM = carried out data analysis and preparation of the manuscript; ØE = mobilizing of resources, data analysis, preparation of the manuscript and overall supervision of the project. All authors read and approved submission of the manuscript for publication.

The authors declare that they have no competing interests. Submit your next manuscript to BioMed Central and we will help you at every step: