Sockeye salmon demonstrate robust yet distinct transcriptomic kidney responses to rhabdovirus (IHNV) exposure and infection

doi:10.21203/rs.2.10333/v1

Background

Aquatic Rhabdoviruses are globally significant pathogens associated with disease in both wild and cultured fish. Infectious hematopoietic necrosis virus (IHNV) is a Rhabdovirus that causes the internationally regulated disease infectious hematopoietic necrosis (IHN) in most species of salmon. Yet not all naïve salmon exposed to IHNV become diseased, and the mechanisms by which some individuals evade or rapidly clear infection following exposure are poorly understood. Here we used RNA-sequencing to evaluate transcriptomic changes in Sockeye salmon, a keystone species in the North Pacific and natural host for IHNV, to evaluate the consequences of IHNV exposure and/or infection on host cell transcriptional pathways.

Results

Immersion challenge of Sockeye salmon smolts with IHNV resulted in approximately 33% infection prevalence, where both prevalence and viral kidney load peaked at 7 days post challenge (dpc). De novo assembly of kidney transcriptomes at 7 dpc revealed that both infected and exposed but noninfected individuals experienced substantial transcriptomic modification; however, stark variation in gene expression patterns were observed between exposed but noninfected, infected, and unexposed populations. GO and KEGG pathway enrichment in concert with differential expression analysis identified that kidney responses in exposed but noninfected fish emphasized a global pattern of transcriptional down-regulation, particularly for pathways involved in DNA transcription, protein biosynthesis and macromolecule metabolism. In contrast, transcriptomes of infected fish demonstrated a global emphasis of transcriptional up-regulation highlighting pathways involved in antiviral response, inflammation, apoptosis, and RNA processing. Quantitative PCR was subsequently used to highlight differential and time-specific regulation of acute phase, antiviral, inflammatory, cell boundary, and metabolic responsive transcripts in both infected and exposed but noninfected groups.

Conclusions

Waterborne exposure with IHNV has a dramatic effect on the Sockeye salmon kidney transcriptome that is discrete between resistant and acutely susceptible individuals. Here we identify that metabolic, acute phase and cell boundary pathways are transcriptionally affected by IHNV and tissue responses to local infection are highly divergent from those generated as part of a disseminated response. These data suggest that primary resistance of naïve fish to IHNV may involve global responses that encourage reduced cellular signaling rather than promoting classical innate antiviral responses.

Epigenetics & Genomics

Sockeye salmon

Infectious hematopoietic necrosis virus

RNA-seq transcriptome

Differential gene expression

metabolism

acute phase response

extracellular matrix

Infectious hematopoietic necrosis virus (IHNV) is an ecologically and economically significant pathogen of salmon and of specific concern to the World Organisation for Animal Health (Office International des Epizooties; OIE) due to its high virulence and potential for transmission within global salmon fisheries [ 1]. Large losses of revenue and fish have occurred in salmon aquaculture as a result of this virus [ 2, 3], and it has been a management concern to both wild and cultured Pacific salmon in North America for almost 60 years [ 4]. IHNV also represents one of at least 11 highly virulent aquatic Rhabdoviruses which together make up one of the most destructive and costly group of infectious disease-causing agents in freshwater. Like all Rhabdoviruses, the IHNV genome (11 KB) consists of a negative-sense single-stranded RNA transcript which encodes a nucleoprotein (N), phosphoprotein (P), matrix protein (M) glycoprotein (G) and a large RNA-dependent RNA polymerase protein (L) [ 5]. Specific to the aquatic Novirhabdoviruses, the IHNV genome additionally encodes a sixth non-virion protein (Nv) believed to function as a host immune suppressor [ 6, 7]. IHNV enters fish through the fin bases, gill or oral/gastrointestinal tract [ 8, 9], disseminates through the circulatory system to virtually all bodily organs with highest loads developing in the anterior kidney (also called the head kidney or pronephros; hereafter referred to simply as kidney) and spleen [ 8, 10]. It is shed into the environment in high quantities during primary amplification, most likely via external mucus [ 11, 12].

Sockeye salmon (Oncorhynchus nerka) is a keystone species in the North Pacific with important ecological, social and economic significance [ 13, 14]; they are also a natural host of IHNV [ 15]. Early developmental stages (alevins, fry and parr) of Sockeye salmon are highly susceptible to the co-evolved U genogroup of IHNV and epizootics are often associated with high mortality [ 15-18]. Older fish appear to be naturally less susceptible to infection; however, IHNV can be transmitted to and cause disease in Sockeye salmon within seawater environs [ 19-21]. Nevertheless, waterborne exposure of Sockeye salmon to IHNV in seawater does not typically result in 100% infection prevalence in heterogeneous populations even at a relatively high exposure dose [ 19, 22, 23]. The mechanisms by which some fish are able to evade or rapidly clear IHNV infection following exposure are of great interest in mitigating disease associated with this and other fish rhabdoviruses.

A commonality for rhabdoviral resistance in naïve fish is that it appears to be decided very early during the infection process and is at least partially determined by heritable genetic factors [ 16, 24-28]. Multiple studies have indicated that this resistance (i.e., the ability of an individual to avoid disease and clear or avoid systemic infection), likely involves elements associated with epidermal tissues [ 25-28]. However, at least one study has also noted that resistance can be maintained even when epidermal tissues are mostly subverted such as following intraperitoneal injection [ 24]. This suggests that similar or alternative resistance factors may also be present in internal tissues. Currently, the best known innate cellular mechanism to inhibit aquatic rhabdoviral replication and promote host resistance in salmon is the activation of the classic interferon response pathway [ 29-32], and non-specific activation of interferon prior to IHNV exposure has been demonstrated to provide significant protection against IHNV infection [ 32-34]. Nevertheless, interferon signaling does not appear to correlate with early resistance in naturally resistant fish [ 24, 35], and the mechanism(s) by which some individuals are able to evade or rapidly clear infection following rhabdoviral exposure remains unclear.

We recently explored the transcriptomic responses of the Sockeye salmon head kidney to Piscine orthoreovirus (PRV) infection and the consequences associated with an IHNV superinfection of PRV-infected individuals [ 23]. Although that study’s main focus was on elucidating responses specific to PRV, it also identified substantial transcriptome modification at 7 dpc with IHNV. Interestingly, global transcriptomic changes within the IHNV exposed group were not only observed in fish which had become infected, but also in individuals which had evaded or rapidly cleared the infection. Further, it appeared that many (>50%) of the transcriptionally altered genes were either unique or had reversed expression in fish which had become IHNV infected versus those which were refractory to the detectable virus replication. We hypothesised that further exploration into these transcriptomic differences may provide insight into the innate cellular responses associated with the ability to evade or rapidly clear IHNV infection following exposure and offer a more complete picture of the pathways involved in host antiviral responsiveness to rhabdoviral infection.

Here, we aimed to revisit and expand on our previous work exploring the transcriptomic responses in the Sockeye salmon kidney following IHNV exposure with a specific focus on the transcriptional differences between resistant (IHNV-) and susceptible (IHNV+) individuals. We also aimed to better capture and assess a global snapshot of the molecular pathways that may be involved in innate antiviral responsiveness to IHNV.

IHNV challenge identified both resistant and susceptible individuals to systemic infection with highest detection in kidney relative to blood or brain

Not all Sockeye salmon exposed to IHNV in this challenge experiment developed kidney infections [ 23]. Nevertheless, kidney samples proved the most likely relative to blood or brain to identify IHNV during any phase of infection and typically had the highest relative loads; particularly early during the infection process (Fig. 1). IHNV reached comparable prevalence within blood and kidney samples (33%; n=8/24 in each tissue) and was highest at 7 dpc in both sample types. In brain, highest prevalence was not observed until 20 dpc and was lower (17%; n=4/24) than observed in either blood or kidney (Fig 1). IHNV loads were highest in kidney relative to blood or brain in most instances; however, there were a few cases late during infection (at 20 and 48 dpc) where blood and/or brain loads were higher than those observed in the kidney (Fig 1a, b & c; Additional file 1). However, no IHNV could be recovered from blood or brain samples in fish which did not also have kidney infections during the first 20 dpc – a period which in our experience encompasses the acute infection phase of IHNV under these experimental conditions. This indicated that fish which did not develop kidney infections during the first 20 dpc were most likely refractory to internalization and systemic dissemination of IHNV.

Fig. 1 Detection of IHNV in experimentally challenged Sockeye salmon. The quantity of IHNV nucleoprotein (N) transcripts for each fish (n=24) as well as relative prevalence (insert) at five discrete time points is presented for kidney (a), brain (b) and blood (c) as detected by qPCR. Kidney samples used for RNA-seq are highlighted in yellow. Agreement for the positive/negative (+/-) detection of IHNV transcripts between kidney and brain (d) as well as kidney and blood (e) is also presented. Hk: head kidney; Br: brain; Bl: blood.

RNA-seq identified significant quantitative differences between IHNV mRNA transcripts with strong correlative agreement to qPCR

Four kidney samples identified to have IHNV as well as four kidney samples identified to have both IHNV and PRV were incorporated into our original RNA-seq transcriptome assembly [ 23]. For the current study, we retained the PRV/IHNV kidney libraries solely in order to correlate RNA-seq and real time quantitative PCR (qPCR) quantitative estimates of IHNV (Fig. 2b); these libraries were not used in transcriptome analysis. In the four non-PRV infected libraries, IHNV was observed to span four base-ten orders of magnitude as estimated by qPCR: one with low (~3,000), two with mid (50,000-100,000), and one with high (~10,000,000 copies IHNV N transcripts per µg total RNA) viral loads.

All six IHNV mRNA protein-coding transcripts (N, P, M, G, Nv, L) were identified during RNA-seq de novo assembly, although G and Nv coding sequences appeared to be combined into a single putative transcript (Fig. 2a). It is unknown if this combined G/Nv transcript was a possible error in Trinity assembly due to read trimming or if the mRNA termination signal between these genomic sections was indeed completely (or partially) ignored during transcription. Nevertheless, the quantitative pattern of expression for each of the five IHNV transcripts identified by RNA-seq was strikingly similar in all four IHNV infected kidney libraries (Fig. 2d) where the relative load of N and P transcripts was significantly higher than M and L at 7 dpc (Fig. 2c). Accuracy for RNA-seq by expectation maximization (RSEM) quantitative estimates was highly supported by correlative agreement to standard-curve-derived quantities identified by qPCR (r=0.974; p < 0.001) (Fig. 2b). Although qPCR was more sensitive than RNA-seq for detecting low levels of IHNV RNA, this was somewhat expected due to the higher dilution of template used in RNA-seq relative to IHNV-targeted qPCR assessments, as well as due to the use of polyadenylation-targeted oligo(dT) priming (specific to mRNAs only) for RNA-seq compared to random hexamer priming used for qPCR cDNA synthesis which detects both genomic and mRNA of IHNV [ 23].

Fig. 2 IHNV detection in the kidney of waterborne challenged Sockeye salmon using RNA-seq. (a) RNA-seq assembled isoform transcripts with UniProt Blastx (E < 10^-5) annotation aligned against the published IHNV genome (yellow). Dark areas represent sequence variation within an overall nucleotide agreement of 97.4%. (b) Correlation plot of IHNV N transcript detection by qPCR relative to RNA-seq RSEM quantities. Pearson correlation coefficient (r) and associated p-value are presented for log₁₀-transformed data. (c) Mean (± SE; n=4) proportion of total IHNV transcripts for each of the five structural viral mRNAs identified by RNA-seq. Letters indicate significant groupings (p < 0.05). (d) IHNV mRNA transcript abundance (RSEM normalized count) identified by RNA-seq in each of the four IHNV infected kidney libraries. N: nucleocapsid; P: polymerase-associated phosphoprotein; M: matrix protein; G: surface glycoprotein; and L: polymerase coding transcripts.

Global transcriptomic expression in Sockeye salmon kidneys revealed stark variation between infected, exposed but noninfected, and mock control fish

Using multidimensional scaling (MDS) and principal component analysis (PCA), we previously reported that seven of the eight IHNV-infected fish oriented separately from all other groups (both PRV and controls), and that differential expression analysis identified 32-42% of all component unigenes to have significant differences in expression following IHNV exposure [ 23]. However, both MDS and PCA assessments in that study used only the top 500 unigenes contributing variance to calculate pairwise distances between libraries and included samples collected at two time points relative to PRV exposure. Although this provided strong representation within the context of this broader PRV/IHNV study, it was noted that temporal changes contributed significant variance in global expression independent of viral exposure. Therefore, we removed non-IHNV relevant libraries from our current analysis thereby removing temporal variation (only one time point for RNA-seq sampling was pertinent to IHNV), and also expanded the analysis to include all unigenes which contributed variance (not just the top 500). This refined our focus for observing IHNV-associated global expression patterns. Further, we expanded the depth of our PCA global expression by assessing isoform transcripts in addition to component unigenes. This identified that two linear modes of variation contributed to the majority of variance in both component unigene and isoform transcript libraries (70% and 59%, respectively; Fig 3a & b). Library ordination within these two linear modes of variation somewhat surprisingly revealed that exposure to IHNV, independent of whether fish had detectable infections, appeared to be the main driving factor for the highest source of variation between libraries (44% in unigenes, 41% in isoform transcripts; Fig. 3c & d). The second largest variation (26% in unigenes, 18% in isoform transcripts) appeared mainly compelled by the internal prevalence (detectable infection) of IHNV in the kidney. Viewed together, these variations exposed starkly separate groupings in global expression patterns between mock non-exposed controls, exposed but noninfected, and infected individuals at 7 dpc (Fig. 3c & d). These analyses also identified that the relative ordination for component unigene expression closely mirrored that of the larger isoform transcript libraries they represent.

Fig. 3 Global variance in gene expression between RNA-seq libraries in relation to IHNV exposure. The 10 most significant linear modes of variation (dimensions) identified during principal component analysis of (a) component unigene or (b) isoform transcript normalized RSEM counts are presented following regularized-logarithmic (rlog) transformation. Orientation of each RNA-seq library within the two highest linear modes of variance is also provided based on (c) component unigene and (d) isoform transcript RSEM counts.

Transcriptome responses favored down-regulation in exposed but noninfected fish and up-regulation in fish infected with IHNV

Given the clear visual separation in transcriptomic expression between noninfected, infected and mock control libraries afforded by PCA, we sought to identify if any global trends in transcriptional regulation might be contributing to the overall variance in expression observed between these treatment groups. We therefore considered the total number of significantly regulated component unigenes and isoform transcripts from noninfected and infected libraries relative to mock non-exposed controls based on three criteria: (i) two probability stringency cut-offs (p < 0.05 and p < 0.001) based on a single false discovery rate (α = 0.05), (ii) the relative proportion for up- and down-regulation in each instance, as well as (iii) how the composition of regulation in these treatment groups compared based on fold change. This revealed that at a stringency of p < 0.05, noninfected fish had approximately an equal number of significantly up and downregulated component unigenes (11,883 up; 11,613 down) or isoform transcripts (17,511 up; 17,117 down) (Table 1), although down-regulation in both instances possessed more genes/isoforms associated with higher fold changes (Fig. 4a, b, e, f). In contrast, IHNV infected fish had similar quantities of upregulated unigenes and isoform transcripts to noninfected fish (11,442 unigenes; 16,966, isoforms); however, the number of downregulated unigenes and isoforms was substantially reduced: 7,960 downregulated unigenes and 9,589 downregulated isoforms representing 70% and 56% of upregulated quantities, respectively (Fig. 4c-f). Furthermore, a larger number of unigenes/isoforms with higher fold changes were associated with upregulation in infected fish, and these general patterns were accentuated when regulation was considered at an even higher statistical stringency. At p < 0.001, a clear emphasis in down-regulation was observed in noninfected fish whereas libraries of infected fish were highlighted with a strong emphasis for upregulation (Fig 4). This was further made evident by directly comparing the regulation of infected libraries to those of noninfected fish rather than mock controls (Table 4).

Similar to PCA analysis, both cumulative count data and relative fold induction patterns were highly mirrored between component unigene and isoform transcript libraries. This indicated that the drivers for mRNA transcriptional responses to IHNV in the Sockeye salmon kidney were mainly regulated at the gene rather than isoform level. As such, our subsequent investigations focused mainly on targeting component unigene rather than isoform transcript libraries.

Fig. 4 Differential gene expression for IHNV infected and noninfected Sockeye salmon. The numbers of differentially regulated (a) component unigenes and (b) isoform transcripts of exposed but noninfected fish (n=4) at FDR adjusted p < 0.05 and p < 0.001 are presented in relation to fold change of time-matched mock-challenged control RSEM quantities 7 dpc. Differentially regulated (c) component unigenes and (d) isoform transcripts of IHNV infected fish (n=4) are also provided, as well as comparisons for the total number of up- and downregulated (e) component unigenes and (f) isoform transcripts for both noninfected and infected groups. The number of differentially regulated unigenes and transcripts common to both infected and noninfected groups is indicated (shaded).

Insight into biological pathways regulated by IHNV exposure and/or infection provided by contig functional annotation

A major challenge for understanding the consequences associated with differential regulation in a transcriptome dataset is in appropriating correct interpretation for contig function, particularly for non-model species such as Sockeye salmon. Here, we first annotated our de novo assembled contigs based on homology to gene sequences within the nucleotide (nt) database hosted by the National Center for Biotechnology Information (NCBI). This resulted in 89% of component unigenes considered for differential expression in this study to have high similarity (E < 1e^-5) with at least one previously identified genomic/transcriptomic sequence from another eukaryotic organism; the vast majority of which (85 of the 89%) having closest similarity to transcripts identified or predicted in Atlantic salmon (Salmo salar) (Fig 5a).

The high similarity of nucleotide sequences assembled here for Sockeye salmon to those identified/predicted in Atlantic salmon provided further evidence as to the overall high quality of our de novo assembly [ 23]; however, the relationship between nucleotide homology and functional orthology can be hard to predict and current strategies for annotating contig function typically focus on comparing putative amino acid rather than nucleotide sequences [ 36]. We therefore used the NCBI blastx tool to predict the amino acid sequence from putative protein-coding unigenes/isoforms and searched for putative orthologs (E < 1e^-5) within the UniProt 2017_1 protein database. In contrast to the relatively large number of annotated salmon nucleotide sequences now available in NCBI, comparatively few proteins have been functionally verified from salmon, other fish species, or non-model organisms in general within protein databases (Swiss-Prot, TrEMBL or UniProt). As a result, these databases are mostly populated with sequences from model mammalian species, namely human and mouse. This was clearly evident following blastx analysis, where less than 33% of Sockeye component unigenes in this study had high similarity to protein sequences available in the UniProt 2017_1 database, and the highest proportion (approximately 10 of the 33%) was associated with human protein sequences (Fig. 5b). Nevertheless, the relative proportion for protein annotation was highly conserved across all differentially expressed gene (DEG) derived datasets (data not shown) and thus allowed for comparison across treatment groups without eliciting bias in comparing protein-associated functional enrichment.

To identify cellular molecular pathways differentially regulated by IHNV exposure and/or infection, we performed an enrichment analysis for KEGG terms associated with significantly differentially regulated (p < 0.05) protein-coding unigenes using KEGG Mapper software [ 37]. As KEGG interactions and pathways are verified at a single species level and are not meant to accommodate multi-species comparison, we were forced to ignore terms associated with proteins of non-model species and only consider those mapped in the model species of either human or mouse which together represented 56% of all KEGG term annotations. This revealed that the majority of terms encompassing both up and downregulation in fish exposed to IHNV were likely associated with pathways relating to host-virus interactions (e.g., pathways regulated during infection with Human T-lymphotropic virus- type 1 and Epstein-Barr virus as well as those involved in viral carcinogenesis) and cell proliferation (e.g., pathways of cancer, PI3K-Akt signaling, Ras/MAPK signaling, and others). Further, these viral-associated responses appeared to strongly influence pathways relating to metabolism (Table 2). Refinement of the enrichment analysis to focus only on KEGG terms associated with downregulated unigenes in exposed but noninfected fish and upregulated unigenes in IHNV infected fish (the two conditions emphasised during global transcript expression profiling) identified that a major portion of enriched terms downregulated in resistant fish appeared to correspond with cell proliferation and macromolecule internalization (e.g., pathways of cancer, PI3K-Akt signaling, MAPK signaling and endocytosis) whereas terms associated with up-regulation in infected fish highlighted pathways relating to protein production (e.g., pathways in cancer, RNA transport, ribosome and spliceosome functioning, as well as protein processing; Table 2). This was further supported by GO term enrichment analysis performed using GOrilla [ 38], which identified that unigenes uniquely downregulated in noninfected fish generally corresponded to metabolic, DNA transcription, and protein biosynthesis cellular processes, whereas unigenes uniquely upregulated in IHNV infected fish appeared to be involved with virus recognition, inflammation, apoptosis, and RNA processing (Fig. 5c & d). Indication for the dichotomous transcriptional responses relating to GO and KEGG pathways was yet further supported by comparing the most significantly differentially regulated transcripts with either the highest fold changes or relative quantitative abundance between noninfected and infected libraries as presented in Additional file 2.

Fig. 5. Transcriptome annotation and Select GO enrichment. The quantity of component unigenes of Sockeye salmon identified with high similarity (Expect-value < 10e-5) to previously published nucleotide sequences of eukaryotic organisms (a) or the amino acid sequence for which the component unigenes putatively translate (b). The percentage of component unigenes annotated within the filtered (with inconsistent and lowly-expressed unigenes removed) transcriptome is provided. Significantly enriched (p <0.05) Gene Ontology (GO) terms related to biological processes, their interactions, and their specificity of function in association with component unigenes uniquely downregulated in exposed but non-IHNV infected fish (c) or uniquely upregulated in IHNV infected fish (d) are also provided. General groupings of terms relating to an over-arching biological process are circled.

Gene expression by qPCR verified differential enlistment of antiviral, acute phase, cell boundary, and metabolic transcripts in IHNV infected versus exposed but noninfected individuals.

RNA-seq transcriptomic analyses were generated for only one time point relative to IHNV (7 dpc) and included challenge groups of no more than four individuals (four exposed but noninfected, four infected and four mock controls). However, at least four additional samples not submitted for RNA-seq were also available in each IHNV-associated challenge group at 7 dpc, as were sets of samples collected at four other time points following IHNV exposure. We therefore utilized these samples for qPCR investigations to both validate RNA-seq DEG analyses within the larger dataset available at 7 dpc, as well as to expand our temporal understanding of the transcriptional expression patterns of select unigenes (Table 3). Gene targets were selected for qPCR following examination of the most significant and highly differentially regulated transcripts (Additional file 2) and/or by identifying unigenes which were suspected to contribute to cellular processes and pathways highlighted as being IHNV responsive during functional enrichment analysis. The selected unigenes were also chosen to represent both classic (e.g., interferon) and potentially novel (e.g., acute phase response) contributions to host defense against IHNV. This selection resulted in three genes associated with the acute phase response (SAA, FTH1 and HAMP), two with extracellular matrix turnover (MMP9 and TIMP2), two relating to inflammation (TNFA and IL1B), four relating to antiviral defense (IFNA, IFNG, MX1 and RSAD2), and four involved in the metabolic pathways specific to either mitochondrial generation of ATP (MT-CO3), glycolysis (PKM), pentose phosphate (TKTL2), or arachidonic acid breakdown (PTGS1) (Table 4).

Of the 15 gene targets selected for qPCR analysis, three showed different patterns of expression by qPCR compared to what had been observed by RNA-seq: MMP9, TKTL2, and MT-C03 (Additional file 3). For each of these three targets, significant differences identified by RNA-seq were not confirmed by qPCR; however, although MMP9 was not verified as significantly downregulated using qPCR for the four samples suggested as significantly downregulated by RNA-seq, analysis using an expanded (n=8 rather than n=4) sample set at this 7 dpc time point did prove MMP9 to be significantly downregulated by qPCR. Further, of the 12 genes identified as differentially regulated by RNA-seq which were attempted to be verified by qPCR, only one (PKM) was not confirmed to be significantly differentially regulated within the expanded 7 dpc dataset (Additional file 3). All three acute phase response-associated genes (SAA, FTH1 and HAMP) were determined to be significantly differentially regulated in response to IHNV at 7 dpc. Both SAA and FTH1 were upregulated only in IHNV infected individuals, whereas HAMP was upregulated in both noninfected and infected individuals at 7dpc (Additional file 3). Inflammatory (IL1B and TNFA) and antiviral defense-associated genes (IFNG, IFNA, MX1 and RSAD2) were also only significantly differentially regulated in IHNV infected individuals at 7dpc, whereas genes representing metabolic pathways associated with energy production (PKM, MT-C03 and TKTL2) were not significantly impacted. Expression of PTGS1 and TIMP2 (associated with arachidonic acid and extracellular matrix turnover, respectively) was significantly reduced in infected fish relative to mock controls (Additional file 3).

In considering gene expression over a period of 1 to 48 dpc, we observed that time, independent of treatment, significantly affected the transcription of all 15 genes targeted by qPCR (Table 4). This, in concert with our previous observations for high temporal changes in transcriptome quantitative estimates by RNA-seq independent of treatment, further highlights the absolute requirement for time-matched controls in live fish experiments such as this if treatment-specific effects are to be determined. Treatment, either dependently or independently considered in relation to time, was also significant by two-way ANOVA for all genes except IL1B (Table 4) which had high variability in expression and only appeared differentially regulated in individuals with correspondingly high blood/kidney loads of IHNV (data not shown).

Analysis for transcriptional regulation at each time point during this study identified that up-regulation of acute phase, inflammatory and antiviral associated target genes peaked at 7dpc whereas down-regulation of metabolic and extracellular matrix associated target genes peaked later during the recovery phase of infection (Fig. 6). None of the qPCR targeted genes associated with inflammation or viral recognition and antiviral defense were significantly elevated in exposed but noninfected fish at any of the time points assessed during this study and were only responsive in fish with detectable IHNV kidney infections. Although acute phase genes appeared mainly responsive to the presence of IHNV, HAMP was also noted to be significantly upregulated in exposed but noninfected fish at 7 dpc; although to a lesser extent than seen in infected fish (Fig. 6; Additional file 3). MMP9 was significantly downregulated in infected fish at 20 dpc as was TIMP2; however TIMP2 was also significantly downregulated at 7 and 48 dpc in infected fish as well as in noninfected fish at 20 dpc. A similar pattern was observed in metabolic associated genes where PKM, TKTL2 and PTGS1 were all significantly downregulated for at least one later time point (20 or 48 dpc) in IHNV infected fish with TKTL2 and PTGS1 also being downregulated in noninfected fish at 20 dpc (Fig.6).

Fig. 6. Temporal expression patterns of selected gene targets in the kidney of Sockeye salmon as assessed by qPCR. Significant (*: p < 0.05; **: p < 0.01; ***: p < 0.001) differential regulation relative to time-matched mock controls is presented. For sample sizes specific to each treatment and time point see Table 3.

Successful internalization and dissemination of IHNV from epithelial surfaces is needed to establish kidney infections and IHN disease [ 39]. In this study, Sockeye salmon which did not develop kidney infections during the first 20 dpc also appeared refractory to systemic IHNV dissemination (Fig. 1). These refractive fish were not only resistant to IHNV during an initial 1 h waterborne exposure but also to putative IHNV shed by infected cohabitants – a presumably substantial quantity during the first week post challenge assuming comparable shedding occurs from Sockeye salmon as has been documented for juvenile Rainbow trout [ 12]. Thus, the kidney transcriptome profiles of IHNV-infected and IHNV-exposed but noninfected Sockeye salmon at 7 dpc in this study likely represent individuals which were either generally susceptible to or resistant to IHNV internalization and dissemination, respectively.

Mucosal epithelia of fish have gained considerable attention in recent years with regards to their immunological functions [ 40] and it is tempting to speculate that recognition of IHNV at these sites contributed to or possibly dictated the disseminated signaling in IHNV resistant fish. However, given that naïve rhabdoviral resistance can be maintained even when epidermal tissues are subverted [ 24] suggests that additional factors may be involved. In this study we did not sample mucosal surfaces and therefore do not know what type of immunological role mucosal tissues may or may not have contributed to IHNV resistance in this study. What was strikingly apparent from this work, however, was that some Sockeye salmon have the capacity to strongly respond to an environmental exposure with IHNV without developing detectable systemic infections. Moreover, host responsiveness was not restricted to epithelial sites but rather disseminated into uninfected tissues such as the kidney to promote global transcriptional changes. What signaling molecules or pathways are involved in this recognition were not determined in this study, but given the global transcriptomic modifications observed in the kidney, it stands to reason that the signaling pathways involved would also affect other internal organs such as the liver and spleen.

Interferons are potent pleiotropic cytokines that can broadly alter cellular functions [ 41], and we initially hypothesized that dissemination of interferons or their effectors from epithelial surfaces following IHNV exposure may be at least partially responsible for the altered kidney transcriptomes observed in noninfected fish in this study. However, classical innate antiviral response pathways including effectors regulated by type I or type II interferon were not transcriptionally induced as part of the kidney response in exposed but noninfected fish. Indeed, exposed but noninfected fish experienced a general down-regulation of many interferon responsive genes involved in antiviral defense [e.g., PI3K-Akt signaling, MAPK signaling, human papillomavirus infection and Human T-Lymphotropic Virus Type 1 (HTLV-1) infection KEGG pathways; see Table 2] at 7 dpc and IFNA, IFNG, MX1 and RSAD2 were not significantly induced in exposed but noninfected fish relative to mock controls at any time point in this study as analysed by qPCR (Fig. 6). This supports previous indications that classic innate antiviral signaling via interferons are not a necessary part of early rhabdoviral resistance in naturally resistant naïve fish [ 24, 25] and were also probably not responsible for the global transcriptional changes observed in kidneys of exposed but noninfected fish in this study.

Conversely, both type I and type II interferons were transcriptionally activated in Sockeye salmon kidneys which became acutely infected with IHNV as were other transcripts associated with innate antiviral defense (Table 2; Fig 5c). Specifically, the robust transcriptional activation of type II interferon (IFNG) (Figure 6) in accompaniment with increased transcriptional processes associated with apoptosis (Fig 5c) are suggestive of robust T-cell responses that are a likely important component of the host defensive strategy to acute Rhabdoviral kidney infections in salmon [ 42, 43].

Acute phase response proteins (APP) are also central constituents in the vertebrate innate immune response [ 44, 45], and we identified differential transcription of multiple acute phase proteins in the kidneys of both exposed but non-infect and infected salmon via RNA-seq that were at least partially confirmed (2 of 3 selected genes) by qPCR. Specifically, SAA, a gene that codes a type I APP, was one of the most highly upregulated transcripts responsive to acute IHNV kidney infection (200-600 fold between 7 and 20 dpc). In mammals, SAA protein is produced mainly by hepatocytes in the liver [ 44]; but in teleost fish, transcription of the SAA gene has been document in response to pathogenic stimulation in multiple organs including the skin, intestine, gills, kidney and spleen [ 46, 47]. It is questionable as to whether transcription of SAA in fish outside the liver leads to increases in serum SAA protein [ 47]; however, our findings confirm that strong transcriptional regulation of SAA is an integrated part of the localized innate pathogen recognition in more than just the liver in teleost fish [ 46, 47]. Transcription of SAA was not, however, significantly induced in the kidneys of exposed but noninfected fish at any time point analyzed during this study. Given that type I APP including SAA are thought to be induced primarily through the inflammatory cytokines IL-1 and TNF [ 45], our results suggest that neither IL-1 or TNF proteins are likely the primary mediators for the disseminated defensive strategy of fish refractive to IHNV in this study. In contrast, HAMP, a gene that codes a type II APP (Hepcidin) involved in the regulation of iron absorption (enterocytes), recycling (macrophages) and storage (hepatocytes) [ 48, 49], was upregulated in kidneys of both infected and noninfected fish during peak infection pressure (Fig. 6). HAMP transcription is considered to be mainly regulated by IL-6 like cytokines [ 45, 48], and it is therefore possible that IL-6 like cytokines could be involved in both the local responses to acute rhabdoviral infection as well as part of a disseminated systemic response to Rhabdovirus exposure. However, SAA upregulation is also tied with responses to stress via the hypothalamic–pituitary–interrenal axis in Rainbow trout (Oncorhynchus mykiss) [ 50] that would likely involve IL-6 like cytokine signaling; and the lack of upregulation of SAA suggests that general stress-associated signaling was not the main driver of the response in exposed but noninfected fish in this study. It is possible that atypical pathways initiating HAMP and SAA transcription are encountered in salmonid fish relative to mammals, or that there are uncharacterized independent transcriptional manipulators of these genes in addition to IL-6 like cytokines. In either case, acute phase response elements appear differentially involved with the innate immune response to Rhabdovirus exposure and infection in salmonid fish which warrants further investigation.

Although we were unable to definitively implicate a specific signaling molecule (or set of molecules) for initiating kidney transcriptomic changes in IHNV exposed but noninfected fish in this study, it is evident that one of the main consequences of this response involves the down-regulation of protein production and metabolism pathways in the kidney [e.g., metabolic, cancer, biosynthesis of secondary metabolites, proteoglycans in cancer, viral carcinogenesis KEGG pathways (Table 2); as well as similar GO processes (Fig. 5c)]. Given that viral replication is dependent on active cellular metabolism and protein production, we hypothesize that the transcriptional regulation of these pathways within the kidney is, at least in part, a defense employed by the host to limit resources required for viral replication thereby reducing the potential for the virus to establish a foothold for robust replication in target tissues. These findings align with integrated human-virus metabolic stoichiometric modelling predictions for a host rather than viral-optimal state [ 51] and point to the possibility that rapid global down-regulation of pathways involving protein production and metabolism may be key components of natural Rhabdoviral resistance in fish.

In contrast, increased energetic consequences appear to be involved in host kidney infections actively responding to local IHNV infections. This transcriptional alteration in IHNV infected fish involving a number of cellular metabolic pathways that involved in protein synthesis (Table 2; Fig 5c) suggests that this antiviral defensiveness and/or cellular damage control comes at a cost of increased protein production. This is further supported by the qPCR monitoring data of three important metabolism-associated genes (PKM, TKTL2 and PTGS1) beyond the peak period of acute kidney infection that showed that all three became significantly downregulated in IHNV infected fish relative to earlier time points or time-matched controls. This indicates that cellular defenses against acute IHNV kidney infections and tissue repair not only alters cellular energetic pathways, but may ultimately manifest in a delayed energetic cost to the organ or organism.

Cell boundary processes, particularly those involved in endocytosis, are key components involved in cellular homeostasis, metabolism, and both intra and intercellular signaling [ 52]. It is therefore not surprising that we see transcriptional changes relating to these processes significantly altered in the Sockeye salmon kidney in response to IHNV; both as a result of acute infection or during a disseminated response following recognition in an alternative organ or tissue. For exposed but noninfected fish, a general downregulation of genes associated with endocytosis was contrasted by a general upregulation in fish that developed kidney infections. The role of endocytosis beyond a vehicle for Rhabdovirus entry [ 53] is not well characterized in fish, but as a defensive strategy, it can be imagined that reducing cellular endocytosis (i.e., cellular ‘eating and drinking’) during viral attack may be beneficial. Similarly, virus infection and tissue damage could be envisioned to create a need for elevated endocytosis and the recruitment of professional phagocytic cells to clean up and rebuild damaged tissues. In either case, the differential expression of these cell boundary process are closely linked to cellular energetics and further point to potential energetic manipulations associated with viral resistance and/or clearance.

Similar differences in immunological and metabolic gene expression have been reported between susceptible and resistant families of Atlantic salmon challenged with infectious pancreatic necrosis virus (IPNV) [ 54, 55]. Gene expression at 1 and 5 dpc of resistant and susceptible families of Atlantic salmon fry challenged with IPNV has shown resistant families to have upregulation of genes related to endocrine function and down regulation of genes involved in tissue differentiation, protein degradation, metabolism and immunological functions when compared to individuals from susceptible families [ 54]. Similarly, individuals in susceptible families (15 to 25% cumulative mortality) have had greater immune responses (higher number of genes and higher expression levels) when compared to individuals from resistant families (<4% cumulative mortality) at 1, 7 and 20 dpc [ 55]. This trend in expression becomes generally reversed during extreme challenge with IPNV (i.e., resulting in cumulative mortality of 76%; [ 56]), further supporting that innate antiviral and inflammatory responses are important in eliminating acute viral infections once established but may not be as important for protecting against initial viral dissemination and infection.

For Sockeye salmon that became systemically infected with IHNV, some (8/24) developed brain infections during the later (resolving) phase of kidney infections which presumably occurred as a consequence of secondary neurological infection following primary infection of the hematopoietic tissues [ 57]. The declining infection prevalence of IHNV in brain of Sockeye at the end of this study suggests IHNV could be cleared by the brain in some instances; however, at least one fish sampled at the end of this study (48 dpc) had a substantial IHNV load in the brain indicating a potential for IHNV to establish persistence in the brain of saltwater-adapted Sockeye salmon juveniles similar to what has been reported in freshwater fry and parr [ 22]. Nonetheless, additional studies are required to understand whether there is link and/or transition between acute vs persistent IHNV infections and the interplay of host responses associated with each.

Through transcriptome profiling of kidney tissues, we confirmed the transcriptional induction of classic antiviral response pathways in Sockeye salmon acutely infected with IHNV. This response was closely tied with the initiation of inflammatory processes, heightened apoptotic signaling and enhancing RNA degradation. Moreover, we identified these responses included robust transcriptional activation of the acute phase response associated gene SAA that as yet has an unknown role for against rhabdovirus infection. In fish that were exposed but did not become infected with IHNV, we also saw global transcriptional changes in the kidney but in a strikingly different manner than that of infected fish; centering on reducing metabolism and protein synthesis processes that linked with reduced endocytosis and extracellular matrix manipulations. The divergence and breadth of these differential transcriptional responses to Rhabdovirus exposure and infection suggest complex and far reaching implications for initiating immunological pathways in fish. Although these energy manipulations are almost certainly part of compensatory allostasis, it is as yet unclear if or how much an overall energetic cost is incurred to the organism for mounting these responses, particularly in fish which are able to avoid or rapidly clear infection. Therefore, future exploration into consequences associated with the transcriptionally affected pathways identified in this study is warranted to establish the energetic costs associated with IHNV exposure and infection. Linking the transcriptional pathways identified herein into physiological consequence will provide greater insight into how innate resistance to rhabdoviral infection is established in fish.

IHNV challenge, sampling, viral screening, and de novo transcriptome assembly

This study incorporated and expanded upon data and samples collected as part of a larger study which challenged Sockeye salmon juveniles with PRV, IHNV, and PRV followed by IHNV [ 23]. As such, the methods associated with generating the experimental treatments, subsequent sampling procedures, viral screening techniques, and the de novo RNA-seq transcriptome assemblies of selected kidney samples have already been described in detail within our previously published report [ 23]. Nevertheless, a brief summary of these procedures as well as additional considerations pertinent to our current investigations are presented herein.

Sockeye salmon juveniles (approximately 40 g each) obtained from Inch Creek Hatchery, British Columbia, Canada, were exposed to 10³ plaque forming units (PFU) of IHNV per mL for 1 h in a static 100 L seawater bath and subsequently cultured in a 360 L system with 5-6 L per min flow-through 11°C sand-filtered UV-treated seawater. No fish exposed to PRV were considered for gene expression explorations within our current analyses; however, as all the fish challenged with IHNV came from within a broader study design incorporating PRV, it should be noted that all fish considered here were administered an intra-peritoneal (i.p.) injection of either Leibovitz-15 (L15) cell culture media or clarified blood lysate of PRV-free Atlantic salmon (diluted 1:10 in L15) exactly 14 days prior to IHNV exposure so as to provide time-matched vehicular control for PRV investigations (see Fig 1a and Table 1 of Polinski et al. [ 23]).

Tissue samples of anterior (head) kidney, blood and brain were collected from a subset of fish (n=4) in each replicate treatment tank (two replicates per treatment) at 1, 3, 20, and 48 dpc and immediately frozen in liquid nitrogen. RNA was extracted and a subset reverse transcribed for the screening of IHNV as previously described [ 23]. Following viral screening, a separate subsample of RNA from anterior kidneys of 8 IHNV challenged fish (4 infected and 4 noninfected) as well as 4 mock control fish at 7 dpc were independently prepared for high-throughput RNA sequencing (RNA-seq) along with 24 PRV-related samples [ 23]. RNA-seq was performed at the Genome Quebec Innovation Centre on an Illumina HiSeq 2000 platform and resulting sequences (45-78 million reads per library) were assembled de novo using Trinity [ 23, 58, 59]. This yielded two quantitative datasets relating to each sample library: a collection of “isoforms” that putatively represent individually unique mRNA transcripts, and a collection of “components” that consist of all isoforms assigned to a single de Bruijn graph during assembly which putatively represent all variant transcripts coded from a single gene sequence [ 58]. In this manuscript, we refer to these as “isoform transcripts” and “component unigenes”, respectively. The relative abundance for each isoform transcript and component unigene was then estimated using the RNA-Seq by Expectation-Maximization (RSEM) method [ 60] within the Trinity software framework and used for differential expression analysis.

Transcriptome annotation and differential expression analysis

We previously reported Trinotate and blastx functional annotation for all component unigenes represented in the Sockeye kidney transcriptome based on data available in the UniProt 2013_11 database [ 23]. However, as the quantity of new entries within the UniProt database has more than doubled in the last four years and many previous entries have been updated, we performed these blastx explorations again using the UniProt 2017_1 database to improve the quality and yield of our transcriptome annotations. Further, we expanded these assessments to include all isoform transcripts and not just those representatives of component unigenes (i.e. the longest isoform transcript variant per unigene). Analyses were conducted at the Canadian Centre for Computational Genomics at McGill University. As in Similar to our previous study, functional annotation was only reported for putative protein-coding isoform transcripts identified by Trinotate which had blastx alignments in the UniProt 2017_1 database with an Expect value (E) < 1e^-5. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) associations were also assigned where available.

Differential expression of both component unigenes and isoform transcripts in libraries of unexposed (n=4), noninfected (n=4), and IHNV infected (n=4) at 7 dpc were compared using DESeq2 [ 61] based on RSEM quantities calculated during de novo assembly. To reduce noise in global transcriptome comparisons and enhance focus in downstream annotation analyses, non-informative transcripts and unigenes without an RSEM quantity of at least 0.5 per million mapped reads in at least four libraries were removed prior to analysis using edgeR software [ 62]. Differential expression was conducted based on previously described protocols [ 23, 63] and a complete set of R-language commands along with session information applied for these analyses are provided in Additional file 4. Significance for differential expression was assessed for IHNV infected and noninfected libraries relative to unexposed mock controls using two cut-off values (p < 0.05 and p < 0.001) following a Benjamini and Hochberg false discovery rate (FDR) adjustment of 5%. Differentially expressed component unigenes were then evaluated for GO and KEGG pathway enrichment using GOrilla [ 38] and KEGG Mapper [ 37] software, respectively, based on pathways and relationships verified in mammalian (human) systems. Cytoscape 3.4.0 [ 64] was then used to visualize relationships between significantly enriched GO terms.

Gene expression by real-time qPCR

Sample and target gene selection

All fish were screened for the presence of IHNV prior to gene expression analyses. A temporal summary for the quantity of samples used for qPCR gene expression analyses is provided in Table 3. Selected samples in the context of a complete sample inventory can also be seen in Additional file 1. As the number of IHNV infected fish was never higher than eight for any sampling time point, all IHNV infected fish were considered for qPCR gene expression assessment. Eight noninfected fish and eight control fish were also selected at each time point, as well as an additional eight noninfected fish at 1 dpc when none of the fish had yet developed detectable internal IHNV infections.

Fifteen unigene targets were selected for qPCR explorations and were chosen to represent cell molecular pathways implicated by GO/KEGG enrichment of RNA-seq DEG at 7 dpc and in many instances represented highly dysregulated transcripts with either an extremely high fold change or overall relative quantity as noted in Additional file 3. A summary of selected gene targets in association with the specific qPCR primer sequences used in this study is provided in Additional file 5. Target-specific primers were designed using Primer3 [ 65] within the Geneious 9.1.7 software platform [ 66] to amplify all highly expressed protein-coding isoform transcripts within the targeted component unigene. All targeted transcripts were verified to code for proteins which had appropriate conserved functional domains [ 67] identified during blastx sequencing alignment using the National Center for Biotechnology Information (NCBI) non-redundant (nr) protein sequence database.

cDNA synthesis and qPCR

A portion of RNA from each sample was treated with DNase I (Life Technologies), cleaned and concentrated using RNeasy MinElute Cleanup (Qiagen), and 1 µg of cleaned RNA was reverse transcribed using a High Capacity cDNA Reverse Transcription Kit (Life Technologies) as previously described [ 23]. Real-time qPCR analyses were conducted on StepOne-Plus real-time detection system (Applied Biosystems) using SYBR green chemistry. Each PCR reaction consisted of 2X SYBR mastermix (Life Technologies), forward and reverse primers (500 nM each; Table 3), and 1 μL cDNA template (diluted 1:5 in nuclease-free water) to a final volume of 15 μL. Samples were assayed in duplicate with a five-step, fourfold dilution series of pooled cDNA included in each run to calculate amplification efficiency, linearity, and provide inter-run calibration. No-template controls as well as no-reverse transcriptase controls were also included on each run. Cycling conditions consisted of an initial activation of DNA polymerase at 95 °C for 10 min, followed by 40 cycles of 5 s at 95 °C, 20 s at 60 °C, and 10 s at 72 °C with florescence measured at the end of each 72°C step. Melt curve analyses were subsequently performed to ensure amplification specificity.

Quantitative estimation and normalization

Relative quantities were calculated from the qPCR raw fluorescence data using the global fitting mechanistic model of Carr and Moore [ 68] within the qpcR R-statistical software package [ 69]. Relative quantities of target genes were then normalized using geNorm [ 70] to the three most stable (ACTB, XRCC6 and CAPN1) of eight putative reference genes as determined by relative stability within RNA-seq libraries at 7dpc (Additional file 6) which was confirmed within the qPCR dataset using geNorm (M-value 0.654).

Statistical comparisons

Corrected normalized relative quantities were compared at 7 dpc by one-way ANOVA and Tukey post hoc test to provide comparison to results obtained by DESeq2 from RNA-seq libraries. A two-way ANOVA followed by Bonferroni post-tests was used to compare temporal expression of each target gene relative to time matched controls across the temporal dataset. As both L15 and L15 + Atlantic salmon blood lysate controls were used for this evaluation (Table 1), one-way ANOVAs were used to ensure no significant difference (p < 0.05) was identified between these control types at any sampling time point prior to pooling. Analyses were performed using Graphpad Prism 6.0. Sum of squares calculations were adjusted for unequal sample sizes where appropriate and a Log₁₀ transformation was applied to all data prior to analysis based on the standard assumption for a log-normal distribution of expression data [ 71].

IHNV: Infectious hematopoietic necrosis virus IHN: Infectious hematopoietic necrosis (disease) VHSV: Viral haemorrhagic septicaemia virus PRV: Piscine orthoreovirus HBSS: Hank’s balanced salt solution PFU: plaque forming units L15: Leibovitz-15 media dpc: days post challenge RSEM: RNA-seq by expectation-maximization FDR: False discovery rate SAA: Serum amyloid A gene FTH1: Ferritin 1 gene HAMP: Hepcidin gene MMP9: Matrix metalloprotease 9 gene TIMP2: Tissue inhibitor of metalloproteinase 2 gene MX1: Myxovirus resistance protein gene RSAD2: Radical S-adenosyl methionine domain-containing protein 2 (a.k.a., Viperin) gene IFNA: Interferon alpha gene IFNG: Interferon gamma-2 gene IL1B: Interleukin 1 beta gene TNFA: Tumor necrosis factor alpha gene PKM: Pyruvate kinase PKM gene MT-C03: Mitochondrial cytochrome oxidase subunit 3 gene TKTL2: Transketolase-like 2 gene PTGS1: Cyclooxygenase 1 gene ACTB: Beta actin gene XRCC6: X-ray repair cross-complementing protein 6 gene CAPN1: Calpain 1 gene DEG: differentially expressed unigene

Ethics and consent to participate

All work with animals was performed in strict accordance with the recommendations set out in the Canadian Council on Animal Care (CCAC) Guide to the Care and Use of Experimental Animals. The protocols were approved by the Pacific Region Animal Care Committee (Animal Use Protocol Number: 13-013 and 15-003). All fish handling was performed under Aquacalm® (Syndel Laboratories Ltd) or MS222 anaesthesia. Euthanasia of all subjects was conducted via a lethal dose of MS222.

Consent for publication

Not Applicable

Availability of data and material

The RNA-seq datasets supporting the conclusions of this article are available in the NCBI Sequence Read Archive repository, study SRP0740078. All other datasets are included within the article, its additional files, or as part of a previous study published in BMC Genomics with open access [ 23].

Competing interests

The authors declare they have no competing interests.

Funding

Financial support was provided through the Program for Aquaculture Regulatory Research project PARR-2014-P-13 awarded to KAG and SCJ. Additional funding was provided by a Canada Research Chair (Tier 2) in marine biotechnology and a Natural Sciences and Engineering Research Council (NSERC) discovery grant to MLR.

Author’s contributions

MPP participated in study design, carried out RNA-seq data interpretation, performed differential expression, qPCR data analyses, and drafted the manuscript. JCB performed RNA extractions, reverse transcription, viral qPCR diagnostic screening, and all gene-targeted qPCR. MLR, KAG and SCJ conceived the study, participated in its design, and aided in the collection of samples. All authors read, contributed to, and approved the final manuscript.

Acknowledgements

The authors would like to extend their sincere gratitude to the staff at Génome Québec Innovation Centre and the Canadian Centre for Computational Genomics who conducted the RNA-seq sequencing and initial bioinformatic analyses, respectively, in generating the reference transcriptome and its annotations used in this study. The authors would also like to thank Holly Hicklin and Elizabeth Shemming for fish care, Jon Richard for technical assistance, and Marlies Rise for sampling assistance.

1. World Organization for Animal Health (OIE): Aquatic Animal Health Code. 2017. www.oie.int

2. Congleton JL: Getting in control of deadly IHN. Focus on Renewable Resources 1988, 14:18-19.

3. Garver KA, Mahony AA, Stucchi D, Richard J, Van Woensel C, Foreman M: Estimation of parameters influencing waterborne transmission of infectious Hematopoietic Necrosis Virus (IHNV) in Atlantic Salmon (Salmo salar). 2013.

4. Amend DF, Yasutake WT, Mead RW: A hematopoietic virus disease of rainbow trout and sockeye salmon. Transactions of the American Fisheries Society 1969, 98(4):796-804.

5. Schu H, Mundt E, Mettenleiter TC: Complete genomic sequence of viral hemorrhagic septicemia virus, a fish rhabdovirus. Virus genes 1999, 19(1):59-65.

6. Thoulouze M-I, Bouguyon E, Carpentier C, Brémont M: Essential role of the NV protein of Novirhabdovirus for pathogenicity in rainbow trout. Journal of virology 2004, 78(8):4098-4107.

7. Wu Y, Guo M, Hua X, Duan K, Lian G, Sun L, Tang L, Xu Y, Liu M, Li Y: The role of infectious hematopoietic necrosis virus (IHNV) proteins in the modulation of NF-κB pathway during IHNV infection. Fish & shellfish immunology 2017, 63:500-506.

8. Drolet B, Rohovec J, Leong J: The route of entry and progression of infectious haematopoietic necrosis virus in Oncorhynchus mykiss (Walbaum): a sequential immunohistochemical study. Journal of Fish Diseases 1994, 17(4):337-344.

9. Harmache A, LeBerre M, Droineau S, Giovannini M, Brémont M: Bioluminescence imaging of live infected salmonids reveals that the fin bases are the major portal of entry for Novirhabdovirus. Journal of virology 2006, 80(7):3655-3659.

10. Chilmonczyk S, Winton J: Involvement of rainbow trout leucocytes in the pathogenesis of infectious hematopoietic necrosis. Diseases of Aquatic Organisms 1994:89-94.

11. LaPatra S, Rohovec J, Fryer J: Detection of infectious hematopoietic necrosis virus in fish mucus. Fish Pathology 1989, 24(4):197-202.

12. Wargo AR, Scott RJ, Kerr B, Kurath G: Replication and shedding kinetics of infectious hematopoietic necrosis virus in juvenile rainbow trout. Virus research 2017, 227:200-211.

13. Gende SM, Edwards RT, Willson MF, Wipfli MS: Pacific Salmon in Aquatic and Terrestrial Ecosystems: Pacific salmon subsidize freshwater and terrestrial ecosystems through several pathways, which generates unique management and conservation issues but also provides valuable research opportunities. AIBS Bulletin 2002, 52(10):917-928.

14. Knapp G, Guettabi M, Goldsmith OS: The economic importance of the Bristol Bay salmon industry. In. Anchorage, AK: Institute of Social and Economic Research; 2013.

15. Williams IV, Amend DF: A natural epizootic of infectious hematopoietic necrosis in fry of sockeye salmon (Oncorhynchus nerka) at Chilko Lake, British Columbia. Journal of the Fisheries Board of Canada 1976, 33(7):1564-1567.

16. Amend DF, Nelson JR: Variation in the susceptibility of sockeye salmon Oncorhynchus nerka to infectious haemopoietic necrosis virus. Journal of Fish Biology 1977, 11(6):567-573.

17. Garver KA, LaPatra SE, Kurath G: Efficacy of an infectious hematopoietic necrosis (IHN) virus DNA vaccine in Chinook Oncorhynchus tshawytscha and sockeye O. nerka salmon. Diseases of Aquatic Organisms 2005, 64(1):13-22.

18. Purcell M, Garver K, Conway C, Elliott D, Kurath G: Infectious haematopoietic necrosis virus genogroup‐specific virulence mechanisms in sockeye salmon, Oncorhynchus nerka (Walbaum), from Redfish Lake, Idaho. Journal of Fish Diseases 2009, 32(7):619-631.

19. Traxler G, Roome J, Kent M: Transmission of infectious hematopoietic necrosis virus in seawater. Diseases of Aquatic Organisms 1993, 16:111-111.

20. Long A, Richard J, Hawley L, LaPatra SE, Garver KA: Transmission potential of infectious hematopoietic necrosis virus in APEX-IHN®-vaccinated Atlantic salmon. Diseases of aquatic organisms 2017, 122(3):213-221.

21. Traxler G, Roome J, Lauda K, LaPatra S: Appearance of infectious hematopoietic necrosis virus (IHNV) and neutralizing antibodies in sockeye salmon Onchorynchus nerka during their migration and maturation period. Diseases of Aquatic Organisms 1997, 28(1):31-38.

22. Müller A, Sutherland BJ, Koop BF, Johnson SC, Garver KA: Infectious hematopoietic necrosis virus (IHNV) persistence in Sockeye Salmon: influence on brain transcriptome and subsequent response to the viral mimic poly (I: C). BMC genomics 2015, 16(1):634.

23. Polinski MP, Bradshaw JC, Inkpen SM, Richard J, Fritsvold C, Poppe TT, Rise ML, Garver KA, Johnson SC: De novo assembly of Sockeye salmon kidney transcriptomes reveal a limited early response to piscine reovirus with or without infectious hematopoietic necrosis virus superinfection. BMC genomics 2016, 17(1):848.

24. Purcell MK, LaPatra SE, Woodson JC, Kurath G, Winton JR: Early viral replication and induced or constitutive immunity in rainbow trout families with differential resistance to Infectious hematopoietic necrosis virus (IHNV). Fish & shellfish immunology 2010, 28(1):98-105.

25. Quillet E, Dorson M, Aubard G, Torhy C: In vitro assay to select rainbow trout with variable resistance/susceptibility to viral haemorrhagic septicaemia virus. Diseases of Aquatic Organisms 2007, 76(1):7-16.

26. Dorson M, Torhy C: Viral haemorrhagic septicaemia virus replication in external tissue excised from rainbow trout, Oncorhynchus mykiss (Walbaum), and hybrids of different susceptibilities. Journal of Fish Diseases 1993, 16(4):403-408.

27. Quillet E, Dorson M, Aubard G, Torhy C: In vitro viral haemorrhagic septicaemia virus replication in excised fins of rainbow trout: correlation with resistance to waterborne challenge and genetic variation. Diseases of Aquatic Organisms 2001, 45(3):171-182.

28. Quillet E, Dorson M, Le Guillou S, Benmansour A, Boudinot P: Wide range of susceptibility to rhabdoviruses in homozygous clones of rainbow trout. Fish & shellfish immunology 2007, 22(5):510-519.

29. Eaton W: Anti-viral activity in four species of salmonids following exposure to poly inosinic: cytidylic acid. Dis Aquat Organ 1990, 9(3):193-198.

30. Saint-Jean SR, Pérez-Prieto SI: Interferon mediated antiviral activity against salmonid fish viruses in BF-2 and other cell lines. Veterinary immunology and immunopathology 2006, 110(1):1-10.

31. Langevin C, van der Aa L, Houel A, Torhy C, Briolat V, Lunazzi A, Harmache A, Bremont M, Levraud J-P, Boudinot P: Zebrafish ISG15 exerts a strong antiviral activity against RNA and DNA viruses and regulates the interferon response. Journal of virology 2013, 87(18):10025-10036.

32. Purcell MK, Kurath G, Garver KA, Herwig RP, Winton JR: Quantitative expression profiling of immune response genes in rainbow trout following infectious haematopoietic necrosis virus (IHNV) infection or DNA vaccination. Fish & shellfish immunology 2004, 17(5):447-462.

33. Kim HJ, Oseko N, Nishizawa T, Yoshimizu M: Protection of rainbow trout from infectious hematopoietic necrosis (IHN) by injection of infectious pancreatic necrosis virus (IPNV) or Poly (I: C). Diseases of Aquatic Organisms 2009, 83(2):105-113.

34. Ooi EL, Verjan N, Haraguchi I, Oshima T, Kondo H, Hirono I, Aoki T, Kiyono H, Yuki Y: Innate immunomodulation with recombinant interferon-α enhances resistance of rainbow trout (Oncorhynchus mykiss) to infectious hematopoietic necrosis virus. Developmental & Comparative Immunology 2008, 32(10):1211-1220.

35. Trobridge G, LaPatra S, Kim C, Leong J: Mx mRNA expression and RFLP analysis of rainbow trout Oncorhynchus mykiss genetic crosses selected for susceptibility or resistance to IHNV. Diseases of Aquatic Organisms 2000, 40(1):1-7.

36. Ng PC, Henikoff S: SIFT: Predicting amino acid changes that affect protein function. Nucleic acids research 2003, 31(13):3812-3814.

37. Kanehisa M, Goto S, Sato Y, Furumichi M, Tanabe M: KEGG for integration and interpretation of large-scale molecular data sets. Nucleic Acids Research 2011, 40(D1):D109-D114.

38. Eden E, Navon R, Steinfeld I, Lipson D, Yakhini Z: GOrilla: a tool for discovery and visualization of enriched GO terms in ranked gene lists. BMC Bioinformatics 2009, 10(1):48.

39. Cain KD, Polinski MP: Infectious diseases of coldwater fish in fresh water. Diseases and disorders of finfish in cage culture 2014:60.

40. Gomez D, Sunyer JO, Salinas I: The mucosal immune system of fish: the evolution of tolerating commensals while fighting pathogens. Fish & shellfish immunology 2013, 35(6):1729-1739.

41. Fritsch SD, Weichhart T: Effects of interferons and viruses on metabolism. Frontiers in immunology 2016, 7:630.

42. Barber GN: Host defense, viruses and apoptosis. Cell death and differentiation 2001, 8(2):113.

43. Dalton DK, Pitts-Meek S, Keshav S, Figari IS, Bradley A, Stewart TA: Multiple defects of immune cell function in mice with disrupted interferon-gamma genes. Science 1993, 259(5102):1739-1742.

44. Cray C, Zaias J, Altman NH: Acute phase response in animals: a review. Comparative medicine 2009, 59(6):517-526.

45. Ebersole JL, Cappelli D: Acute‐phase reactants in infections and inflammatory diseases. Periodontology 2000 2000, 23(1):19-49.

46. Revathy KS, Umasuthan N, Whang I, Lee Y, Lee S, Oh M-J, Jung S-J, Choi CY, Park C-J, Park H-C: A novel acute phase reactant, serum amyloid A-like 1, from Oplegnathus fasciatus: genomic and molecular characterization and transcriptional expression analysis. Developmental & Comparative Immunology 2012, 37(2):294-305.

47. Villarroel F, Casado A, Vásquez J, Matamala E, Araneda B, Amthauer R, Enriquez R, Concha MI: Serum amyloid A: a typical acute-phase reactant in rainbow trout? Developmental & Comparative Immunology 2008, 32(10):1160-1169.

48. Nemeth E, Valore EV, Territo M, Schiller G, Lichtenstein A, Ganz T: Hepcidin, a putative mediator of anemia of inflammation, is a type II acute-phase protein. Blood 2003, 101(7):2461-2463.

49. Rossi E: Hepcidin-the iron regulatory hormone. Clinical Biochemist Reviews 2005, 26(3):47.

50. Talbot AT, Pottinger TG, Smith TJ, Cairns MT: Acute phase gene expression in rainbow trout (Oncorhynchus mykiss) after exposure to a confinement stressor: a comparison of pooled and individual data. Fish & shellfish immunology 2009, 27(2):309-317.

51. Aller S, Scott A, Sarkar-Tyson M, Soyer OS: Integrated human-virus metabolic stoichiometric modelling predicts host-based antiviral targets against Chikungunya, Dengue and Zika viruses. Journal of The Royal Society Interface 2018, 15(146):20180125.

52. Sigismund S, Confalonieri S, Ciliberto A, Polo S, Scita G, Di Fiore PP: Endocytosis and signaling: cell logistics shape the eukaryotic cell plan. Physiological reviews 2012, 92(1):273-366.

53. Albertini AA, Baquero E, Ferlin A, Gaudin Y: Molecular and cellular aspects of rhabdovirus entry. Viruses 2012, 4(1):117-139.

54. Reyes-López FE, Romeo JS, Vallejos-Vidal E, Reyes-Cerpa S, Sandino AM, Tort L, Mackenzie S, Imarai M: Differential immune gene expression profiles in susceptible and resistant full-sibling families of Atlantic salmon (Salmo salar) challenged with infectious pancreatic necrosis virus (IPNV). Developmental & Comparative Immunology 2015, 53(1):210-221.

55. Robledo D, Taggart JB, Ireland JH, McAndrew BJ, Starkey WG, Haley CS, Hamilton A, Guy DR, Mota-Velasco JC, Gheyas AA: Gene expression comparison of resistant and susceptible Atlantic salmon fry challenged with Infectious Pancreatic Necrosis virus reveals a marked contrast in immune response. BMC genomics 2016, 17(1):279.

56. Cofre C, Gonzalez R, Moya J, Vidal R: Phenotype gene expression differences between resistant and susceptible salmon families to IPNV. Fish physiology and biochemistry 2014, 40(3):887-896.

57. LaPatra S, Lauda K, Jones G, Walker S, Shewmaker B, Morton A: Characterization of IHNV isolates associated with neurotropism. Veterinary research 1995, 26(5-6):433-437.

58. Grabherr MG, Haas BJ, Yassour M, Levin JZ, Thompson DA, Amit I, Adiconis X, Fan L, Raychowdhury R, Zeng Q: Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nature biotechnology 2011, 29(7):644-652.

59. Haas BJ, Papanicolaou A, Yassour M, Grabherr M, Blood PD, Bowden J, Couger MB, Eccles D, Li B, Lieber M: De novo transcript sequence reconstruction from RNA-seq using the Trinity platform for reference generation and analysis. Nature protocols 2013, 8(8):1494-1512.

60. Li B, Dewey CN: RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome. BMC Bioinformatics 2011, 12(1):323.

61. Love MI, Huber W, Anders S: Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol 2014, 15(12):550.

62. Robinson MD, McCarthy DJ, Smyth GK: edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics 2010, 26.

63. Love MI, Anders S, Kim V, Huber W: RNA-Seq workflow: gene-level exploratory analysis and differential expression. F1000Research 2015, 4.

64. Shannon P, Markiel A, Ozier O, Baliga NS, Wang JT, Ramage D, Amin N, Schwikowski B, Ideker T: Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Research 2003, 13(11):2498-2504.

65. Rozen S, Skaletsky H: Primer3 on the WWW for general users and for biologist programmers. In: Bioinformatics methods and protocols. Springer; 2000: 365-386.

66. Kearse M, Moir R, Wilson A, Stones-Havas S, Cheung M, Sturrock S, Buxton S, Cooper A, Markowitz S, Duran C: Geneious Basic: an integrated and extendable desktop software platform for the organization and analysis of sequence data. Bioinformatics 2012, 28(12):1647-1649.

67. Marchler-Bauer A, Lu S, Anderson JB, Chitsaz F, Derbyshire MK, DeWeese-Scott C, Fong JH, Geer LY, Geer RC, Gonzales NR: CDD: a Conserved Domain Database for the functional annotation of proteins. Nucleic Acids Research 2010, 39(suppl_1):D225-D229.

68. Carr AC, Moore SD: Robust quantification of polymerase chain reactions using global fitting. 2012.

69. Ritz C, Spiess A-N: qpcR: an R package for sigmoidal model selection in quantitative real-time polymerase chain reaction analysis. Bioinformatics 2008, 24(13):1549-1551.

70. Vandesompele J, De Preter K, Pattyn F, Poppe B, Van Roy N, De Paepe A, Speleman F: Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. Genome Biol 2002, 3(7):18.

71. Beal J: Biochemical complexity drives log-normal variation in genetic expression. Engineering Biology 2017, 1(1):55-60.

Table 1. Kidney transcriptome differential expression summary at 7 dpc IHNV

	Component unigenes		Isoform transcripts
Comparison	Down	Up	Down	Up
p < 0.05
Mock vs. all challenged	10,491 (20%)	12,526 (24%)	15,144 (15%)	19,859 (20%)
Mock vs. noninfected	11,613 (22%)	11,883 (23%)	17,117 (17%)	17,511 (18%)
Mock vs. IHNV infected	7,960 (15%)	11,442 (22%)	9,589 (10%)	16,966 (17%)
Noninfected vs. IHNV infected	2,982 (6%)	5,239 (10%)	2659 (3%)	6,505 (7%)
p < 0.001
Mock vs. all challenged	3656 (7%)	3,853 (7%)	4,719 (5%)	6,064 (6%)
Mock vs. noninfected	3,383 (7%)	2,032 (4%)	3,706 (4%)	2,957 (3%)
Mock vs. IHNV infected	2,215 (4%)	4,820 (9%)	2,392 (2%)	6,779 (7%)
Noninfected vs. IHNV infected	273 (0.5%)	1,859 (4%)	296 (0.3%)	2,492 (3%)

Table 2. KEGG pathway enrichment in differentially expressed unigenes in the kidney of Sockeye salmon following IHNV exposure.

Pathway	Name	Assigned unigenes
Regulated by IHNV exposure (exposed but noninfected & infected)
ko01100	Metabolic pathways	600
ko05200	Pathways in cancer	197
ko01110	Biosynthesis of secondary metabolites	183
ko04151	PI3K-Akt signaling pathway	129
ko04144	Endocytosis	125
ko01130	Biosynthesis of antibiotics	123
ko05166	HTLV-I infection	121
ko04010	MAPK signaling pathway	112
ko04014	Ras signaling pathway	102
ko05205	Proteoglycans in cancer	99
ko04015	Rap1 signaling pathway	99
ko05169	Epstein-Barr virus infection	98
ko04810	Regulation of actin cytoskeleton	97
ko04141	Protein processing in endoplasmic reticulum	96
ko04510	Focal adhesion	95
ko03013	RNA transport	94
ko05016	Huntington's disease	93
ko01120	Microbial metabolism in diverse environments	88
ko05203	Viral carcinogenesis	87
ko00230	Purine metabolism	87
Upregulated in infected fish
Ko01100	Metabolic pathways	156
Ko05200	Pathways in cancer	58
Ko03010	Ribosome	53
ko05016	Huntington's disease	52
Ko03040	Spliceosome	52
Ko04714	Thermogenesis	51
Ko00190	Oxidative phosphorylation	45
Ko05010	Alzheimer’s disease	45
Ko04141	Protein processing in endoplasmic reticulum	44
Ko03013	RNA transport	43
Downregulated in exposed but noninfected fish
Ko01100	Metabolic pathways	180
Ko05200	Pathways in cancer	110
Ko05165	Human papillomavirus infection	74
ko04151	PI3K-Akt signaling pathway	69
ko04144	Endocytosis	68
ko04010	MAPK signaling pathway	66
Ko01110	Biosynthesis of secondary metabolites	64
ko05166	HTLV-I infection	60
Ko05205	Proteoglycans in cancer	59
Ko05203	Viral carcinogenesis	57

KEGG pathways with greatest number of assigned unigenes significantly (p < 0.05) regulated in response to exposure and/or infection with IHNV at 7 dpc as mapped in mammalian (human) systems.

Table 3. Temporal sampling used for differential expression by qPCR/RNA-seq

Days post challenge	Number of fish screened
	Mock exposure	Exposed, noninfected	Exposed, infected
1	8	16	0
4	8	8	5
7	8 (4)	8 (4)	8 (4)
20	8	8	6
48	8	8	2

Number of Sockeye salmon kidney samples used in qPCR (or RNA-seq where applicable) gene expression analysis at each of five time points following IHNV challenge.

Table 4. Differential expression associated with qPCR gene target genes in kidney of Sockeye salmon following IHNV Exposure.

Gene	Protein Encoded	Two-way ANOVA result
		Treatment	Time	Interaction
SAA	Serum amyloid A	***	***	***
FTH1	Ferritin	ns	**	*
HAMP	Hepcidin	***	***	***
MMP9	Matrix metalloprotease 9	***	***	*
TIMP2	Tissue inhibitor of metalloproteinase 2	***	***	*
IFNA	Interferon alpha	***	***	***
IFNG	Interferon gamma 2	***	***	***
MX1	Myxovirus resistance	***	***	***
RSAD2	Viperin	***	***	***
IL1B	Interleukin 1 beta	ns	***	ns
TNFA	Tumor necrosis factor alpha	**	***	**
PKM	Pyruvate kinase PKM	**	***	*
MT-CO3	Mitochondrial cytochrome oxidase subunit 3	ns	***	*
TKTL2	Transketolase-like 2	***	***	***
PTGS1	Cyclooxygenase 1	***	***	**

Significance (*: p < 0.05; **: p < 0.01; ***: p < 0.001; ns: not significant) for differential expression is presented in association with treatment (IHNV infected, exposed but noninfected, or unexposed mock control), time (1, 4, 7, 20, 48 dpc) or the interaction between treatment and time.

Sockeye salmon demonstrate robust yet distinct transcriptomic kidney responses to rhabdovirus (IHNV) exposure and infection

Status:

Journal Publication

Version 1

Abstract

Figures