Impact Factor 6.429

The 5th most cited journal in Immunology

Original Research ARTICLE

Front. Immunol., 22 November 2017 | https://doi.org/10.3389/fimmu.2017.01634

Dynamic Changes in Host Gene Expression following In Vitro Viral Mimic Stimulation in Crocodile Cells

imageSubir Sarker1, imageYinan Wang2, imageBrenden Warren-Smith1 and imageKarla J. Helbig1*
  • 1Department of Physiology, Anatomy and Microbiology, School of Life Sciences, La Trobe University, Melbourne, VIC, Australia
  • 2Genomics Research Platform, La Trobe University, Melbourne, VIC, Australia

The initial control of viral infection in a host is dominated by a very well orchestrated early innate immune system; however, very little is known about the ability of a host to control viral infection outside of mammals. The reptiles offer an evolutionary bridge between the fish and mammals, with the crocodile having evolved from the archosauria clade that included the dinosaurs, and being the largest living reptile species. Using an RNA-seq approach, we have defined the dynamic changes of a passaged primary crocodile cell line to stimulation with both RNA and DNA viral mimics. Cells displayed a marked upregulation of many genes known to be involved in the mammalian response to viral infection, including viperin, Mx1, IRF7, IRF1, and RIG-I with approximately 10% of the genes being uncharacterized transcripts. Both pathway and genome analysis suggested that the crocodile may utilize the main known mammalian TLR and cytosolic antiviral RNA signaling pathways, with the pathways being responsible for sensing DNA viruses less clear. Viral mimic stimulation upregulated the type I interferon, IFN-Omega, with many known antiviral interferon-stimulated genes also being upregulated. This work demonstrates for the first time that reptiles show functional regulation of many known and unknown antiviral pathways and effector genes. An enhanced knowledge of these ancient antiviral pathways will not only add to our understanding of the host antiviral innate response in non-mammalian species, but is critical to fully comprehend the complexity of the mammalian innate immune response to viral infection.

Introduction

The innate immune system carries a substantial burden of defense against viral pathogens. The study of this response across animal species in recent years, as well as the examination of the phylogenetic conservation of these responses has changed our concept of innate immunity; however, much of this work has been performed in mammals. Major explorations of antiviral innate immune responses in non-mammalian species remains very limited, but may be critical to fully comprehend the complexity of the mammalian innate immune response to viral infection, and the discovery of novel antiviral therapeutics (1).

Mammalian cells have been shown to orchestrate elaborate defense mechanisms to detect and inhibit viral replication. Immediately after viral sensing by the host cells, the innate immune response is initiated by germline-encoded molecules termed pattern recognition receptors (PRRs). PRPs recognize conserved features of viruses and other microorganisms, known as pathogen-associated molecular patterns (PAMPs), which are small molecular motifs recognized as non-self, such as microbial nucleic acids, proteins, and carbohydrates (2). Following the recognitions of PAMPs, PRRs initiate a set of signaling cascades, which ultimately result in the production of interferon (IFN) and the upregulation of hundreds of interferon-stimulated genes (ISGs) (3). The expression of these ISGs is known to limit pathogens, in particularly viral pathogens, although the exact role of the majority of these ISGs remains unknown (4), more specifically we have very little understanding of how this response is orchestrated in non-mammalian vertebrate species.

There is an enormous lack of information surrounding antiviral innate immunity in the Reptilia class, which represents a bridge between fish and mammals. Both type I and type III IFNs are known to be central cytokines in the antiviral response in mammals, inducing the upregulation of hundreds of antiviral effector genes (5). To date, both functional type I and III IFNs have been found in the genomes of amphibians (6, 7), with type I IFNs also being found in fish (8). In reptiles, type III IFN has recently been found in the genome of lizards (9); however, the pathways that upregulate their production have not been described to date, and there is very little information on the downstream antiviral effector genes that may be responsible for viral control in reptiles. Recent transcriptomics work performed on non-infected tissue in the lizard has been able to identify the presence of multiple known PRRs in the reptile (10) with recent studies by our group also showing a number of known ISGs to be upregulated in the presence of viral infection in a reptile in vitro (11). The ISGs, viperin, 2′-5′-oligoadenylate synthetase (OASL) and IFN-induced GTP-binding protein Mx1 were demonstrated to be upregulated in saltwater crocodile, C. porosus LV-1 cells in the presence of the viral mimics, dsRNA, and dsDNA, as well as in the presence of replicating dengue virus (11). This study also demonstrated that crocodile viperin retained its antiviral activity and was able to inhibit dengue viral replication in vitro. Given the large number of viruses described to infect reptiles (12), these studies only give a very small insight into the induction of antiviral pathways in the reptile.

The world’s largest living reptile species, the saltwater crocodile, is a member of the prehistoric order Crocodylia, evolved from the archosauria clade that includes the dinosaurs, pterosaurs, crocodilians, and birds, the latter two being the only extant members of the clade (13). In recent years, two novel herpesviruses, crocodyline herpesvirus 1 and crocodyline herpesvirus 2 (CrHV-1 and -2, respectively) (14, 15), adenovirus, and poxvirus have emerged as significant viral infections of the crocodile, with multiple bacterial and fungal pathogens also being detected in this saltwater crocodile (16). However, very little is known about the ability to control these pathogens by the host immune system. The recent advancement of the crocodilian genome sequence (17) has given the opportunity to unveil the innate immune pathways in this ancient species, in particular its response to viral pathogens.

In recent years, transcriptomes profiling using high-throughput RNA sequencing (RNA-seq) technology has provided unprecedented opportunities to study the host response to infection against a wide range of viral and bacterial infections (1822). RNA-seq technology is both an efficient and accurate tool to reveal the systemic changes in host gene expression in response to infectious pathogens, which could help to unearth a better understanding of innate immune pathways in Reptilia. In the present study, we have used RNA-seq technology to comprehensively study the host transcriptomic profile following viral mimic stimulation of crocodile cell lines using both dsRNA and dsDNA. This study provides a global view for the first time, of nucleic acid-specific and post-stimulation time-specific mRNA profiles in the saltwater crocodile (C. porosus), adding significantly to the body of knowledge surrounding the early innate host response of reptiles to a virus.

Materials and Methods

Cells Stock and Stimulation Using dsRNA and dsDNA

The Crocodylus porosus liver cell line, LV-1, used in this experiment was generated by the Berrimah Veterinary Laboratory (Department of Primary Industry and Resources, Government of Northern Territory), and was maintained at 28°C without CO2, in M199 media containing 10% FCS and antibiotics as outlined previously (11). The materials and methods used to stimulate LV-1 with poly dA:dT and low molecular weight (LMW) poly I:C (Invivogen, CA, USA) have been described previously (11). Three replicates per group were stimulated with poly dA:dT (dsDNA) and LMW poly I:C (dsRNA) at a concentration of 2 µg/mL, and cells were harvested at either 8, 24, 48, 72 or 96 h to examine the dynamic changes in host gene expression.

RNA Extraction and RT-PCR

Total cellular RNA was isolated from each individual replicate and purified using an Isolate II RNA Mini Kit incorporating an on-column DNase treatment step (Bioline). The quality and quantity of the RNA was assessed using a NanoDrop spectrophotometer, and an Agilent 2100 Bioanalyzer (Agilent Technologies, USA). A260/280 ratios >1.8 and RNA integrity numbers >9.0 were standard for all total RNA samples purified across the stimulation time course. Initial synthesis of cDNA and subsequent real-time PCR (qPCR) were performed on an ABI 7000 as previously described using the primers and qPCR conditions, including the control genes GAPDH for the LV-1 cells, and primer sets for viperin, OASL, MX1, and ISG20L (11). Primers for CXCL10 were 5′-tgtgagcgccttgagatcat and 5′-gctgccacgtttagacttgtt; RTP2 5′-gtgacttcagcgagccagta and 5′-tccacggactctccatagca; ACOD1 5′-agtgggactactgggtagca and 5′-agaccatgcctagctgcatt.

RNA-seq Library Construction and Sequencing

The protocol for RNA-seq library preparation was adapted from the Illumina TruSeq® RNA sample preparation v2 Kit. Twelve strand-specific Illumina® RNA-seq libraries were generated (three libraries, for each group of control; three replicates for dsRNA at the 24 h stimulation; and six replicates for each group of dsRNA and dsDNA at the 48 h stimulation) using 0.5 µg of total RNA. Total RNA was heated at 65°C for 5 min to denature any secondary structure and facilitate binding of the poly(A) RNA to the oligo-dT beads. Purification of poly(A) RNA was performed using Illumina TruSeq® RNA sample preparation v2 Kit according to the manufacturer’s instructions (Illumina® Inc., San Diego, CA, USA). Purified poly(A) RNA was then fragmented and primed using 19.5 µL of Elute, Prime and Fragment Mix containing random hexamers (Illumina® Inc., San Diego, CA, USA) for 8 min at 94°C.

Synthesis of first-strand cDNA was performed immediately by incubating fragmented and primed mRNA with first strand master mix (Illumina, USA) and SuperScript® II reverse transcriptase (Invitrogen™) at a ratio of 9:1. The reaction mixtures were run in a thermal cycler at 25°C for 10 min, at 42°C for 50 min, and 70°C for 15 min. Second-strand cDNA synthesis was initiated immediately, by adding 25 µL of second-strand master mix (Illumina® Inc., San Diego, CA, USA) to each well, when the thermal cycle reached 4°C. After gentle and thorough mixing, the plate was incubated in a pre-heated thermal cycle at 16°C for 1 h. The double-stranded cDNA (ds cDNA) was subsequently purified by using AMPure XP beads (Invitrogen™, USA) according to the manufacturer’s instructions and eluted in 50 µL of the resuspension buffer (Illumina® Inc., San Diego, CA, USA).

Blunt-end repair of ds cDNA was performed in a 100 µL reaction containing 10 µL diluted End Repair Control to 1/100 in resuspension buffer and 40 µL of End Repair Mix to each well (Illumina, USA). Reactions were incubated on a pre-heated thermal cycler at 30°C for 30 min and the ds cDNA was cleaned up using AMPure XP beads (Invitrogen™, USA) according to the protocol described in Illumina TruSeq® RNA sample preparation v2 Kit, and eluted in 15 µL of the resuspension buffer (Illumina® Inc., San Diego, CA, USA).

To facilitate Illumina® adaptor ligation, a single “A” nucleotide was added to the 3′ ends of the blunt-end-repaired cDNA samples. Fifteen microliters of purified phosphorylated blunt-end-repaired cDNA was included in a final 30 µL reaction mixture containing 2.5 µL diluted A-tailing control to 1/100 in resuspension buffer and 12.5 µL of A-tailing mix to each well (Illumina® Inc., San Diego, CA, USA). The reaction mixtures were run in a thermal cycler at 37°C for 30 min followed by 70°C for 5 min.

Illumina® RNA-seq adaptor ligations were performed in a reaction volume of 37.5 µL containing 30 µL of phosphorylated blunt-ended cDNA in addition to 2.5 µL diluted ligation control, 2.5 µL ligation mix, and 2.5 of custom indexed adaptors to each well (see Table S1 in Supplementary Material for barcode index sequences). Reaction mixtures were incubated for 30°C for 10 min followed by inactivation of the ligation with 5 µL of stop ligation buffer into each well (Illumina® Inc., San Diego, CA, USA). Adaptor-ligated cDNA was purified using AMPure XP beads (Invitrogen™, USA) according to the protocol described in Illumina TruSeq® RNA sample preparation v2 Kit and eluted in resuspension buffer in a final volume of 20 µL.

PCR amplification of selectively enriched DNA fragments (50 µL) was performed using 20 µL of adaptor-ligated cDNA, 5 µL of Illumina® PCR primer cocktail and 25 µL of Illumina® PCR master mix to each well (Illumina® Inc., San Diego, CA, USA). PCR amplification reactions were performed with the following temperature cycling profile: 98°C initial denaturation for 30 s; 15 cycles of 98°C for 10 s, 60°C for 30 s, and 72°C for 30 s; and 72°C final extension step for 5 min. PCR products were purified to remove PCR-generated adaptor-dimers using AMPure XP beads (Invitrogen™, USA) according to the protocol described in Illumina TruSeq® RNA sample preparation v2 Kit with final elution in 30 µL of resuspension buffer.

All RNA-seq libraries were quantified and assessed using an Agilent Tape Station (Agilent Technologies) by the Australian Genomic Research Facility (AGRF, Melbourne) and confirmed insert sizes of 100–125 bp for all individual libraries. Individual RNA-seq libraries were standardized and pooled in equimolar quantities. The quantity and quality of the final pooled library was assessed as described above prior to sequencing by the facility. Cluster generation and sequencing of the pooled RNA-seq libraries were sequenced as paired-end using Illumina® HiSeq HT chemistry according to the manufacturer’s instructions. RNA-seq data from this study have been deposited in the NCBI Sequence Read Achieve (SRA) under the accession number of PRJNA399550 (SUB2572918, SUB2982437) (http://www.ncbi.nlm.nih.gov/sra/).

Transcriptome Data Analysis

An initial quality check was performed on each of the raw read files using FastQC (version 0.11.5) (23) to determine the best sequence read quality for filtering strategy. Low-quality reads were filtered considering the following criteria: (i) reads containing more than 25% bases with a phred score <20; (ii) reads with an average quality score <20; (iii) reads containing more than 10% of skipped bases (marked as “N”). The trimmomatic (version 0.32) (24) was used to remove the adapter sequences. The HiSat2 RNA-seq strand-specific aligner software package (version 2.0.5) with default parameters (25) was used to align filtered sequence reads to the most recent version of the Crocodylus porosus reference genome (assembly GCA_001723895.1 CroPor_comp1). Aligned sequence reads in individual BAM files were then used for a final quality check such as quality control of gene body coverage, duplication level, splice junction assessment, deletion, and mismatch profiles using RSeQC (version 2.6.4) (26), and all samples successfully passed.

The featureCounts tool, which is part of Subread software package (version 1.4.6p5) (27, 28), was used to quantify gene expression as reads count against Crocodylus porosus reference genome annotation (GCA_001723895.1 CroPor_comp1). Raw counts were loaded into the Bioconductor edgeR package (v 3.18.1) (29) and was used for the differentially expressed genes (DEGs) analysis. Expressed genes with a set threshold were defined as log2-count per million (logCPM) >0.5 in at least three different replicates. For each library, a trimmed mean of M values based normalization factor was calculated to eliminate composition biases between libraries (30). Furthermore, we also performed Voom transforms of RNA-seq data (31) for linear modeling, and it transformed count data to logCPM, which was used to estimate the mean variance relationships to compute the observation-level weights. The Bioconductor limma package was used to test for differential expression. Empirical Bayes method (ebayes) shrinkage was performed on the variances, and estimated moderated t-statistics and the associated p values. DEGs were defined by setting a false discovery rate (FDR) ≤0.01. Top 200 and 500 most variable genes were used to calculate a matrix of Euclidean distance for examining hierarchical clustering of samples in heat maps, using the heatmap.2 function of gplots (v 3.0.1).

Functional Analysis of DEGs

Considering the lack of gene ontology (GO) for this non-model saltwater crocodile, we annotate GO term associated with crocodile protein coding genes by first comparing ortholog genes of Gallus gallus (chicken) as a reference organism and identified the molecular functions, biological processes, and cellular components for 8,399 genes, and then using human ortholog genes annotated another 2,928. TopGO (v 2.28.0) was used to perform the enrichment analysis on the DEGs with the total 11,327 genes as background. To characterize the identified genes from DEG analysis, a GO-based trend test was performed using Fisher’s exact test; p values <0.001 as a cut-off for statistically significant enrichment. To retrieve biological pathways related to the crocodile, we obtained KEGG orthologous gene information using the KEGG Automatic Annotation Server (32) with a eukaryote representative set as reference and default parameters, and DEG genes were linked to KEGG pathways using an in-house Perl script and edited based on KEGG Mapper results.

Results

Preliminary Analysis of Dynamic Changes in Selected Host Genes

To better understand the ability of the crocodilian family members to elicit an early innate immune response to the recognition of dsRNA and dsDNA, we conducted a time-course stimulation of the C. porosus liver cell line, LV-1. LV-1 cells were stimulated for 8, 24, and 48 h to assess the ability of the C. porosus cells to activate PRRs, and initiate early innate signaling pathways to induce ISGs. A thorough examination of the C. porosus genome database for the presence of ISGs revealed assigned gene products with close homologies to other species of ISG20L, Mx1, and viperin. The selected genes were screened using primers and already established qPCR protocols on LV-1 cells (11), following 8, 24, and 48 h of stimulation with dsDNA and dsRNA. As can be seen in Figure 1, the ISG viperin responded strongest to both stimuli, with maximal expression following dsRNA plateauing at 24–48 h, and maximal expression following dsDNA stimulation seen at 48 h, albeit considerably lower than was observed for dsRNA stimulation (6-fold upregulation versus 387-fold upregulation, respectively). No observable response in ISG20L was seen for either stimuli, and a maximal increase in Mx1 mRNA to dsRNA was observed at 48 h (13-fold change) (Figure 1C). Stimulation of LV-1 cells with dsRNA induced the strongest response in ISG mRNA expression, which appeared to peak around 24–48 h, with a much weaker ISG response seen against dsDNA, which was maximal at 48 h post stimulation (hps).

FIGURE 1
www.frontiersin.org

Figure 1. Time-course stimulation of LV-1 cells with dsRNA and dsDNA. Quantitative real-time PCR was performed to determine the changes in gene regulation of the ISGs. Relative expression levels of viperin, Mx1, and ISG20L were determined in comparison to GAPDH at 8 hrs (A), 24 hrs (B) and 48 hrs (C). The graph shows mean ± SD of values from three replicates per group. **Significant at the p < 0.001 level, *significant at the p < 0.05 level.

Preliminary Analysis of the High-Throughput RNA-seq Data

To further understand the dynamic gene expression profiles and their role in the immune pathway, RNA-seq was performed to explore the transcriptomes from crocodile LV-1 cells stimulated with the viral mimics, dsRNA and dsDNA at 24 and 48 hps, the times of maximal ISG expression in Figure 1. More than 63.21 million reads per library were generated, of which 58.87 million reads per library (93.21%) remained after adapter sequence and poor quality reads trimming (Table S1 in Supplementary Material). Alignment of the clean RNA-seq reads to the C. porosus reference genome yielded an average of 46.45 million reads (73.57%) per library uniquely mapped to the crocodile genome (Figure S1 in Supplementary Material). The RNA-seq read depth was distributed evenly along the whole body of the genes (Figures S2A,B in Supplementary Material), reflecting no obvious bias being introduced during randomly primed reverse transcription and subsequent RNA-seq.

Analysis of Differential Gene Expression from RNA-seq Transcriptomes

Following the preliminary RNA-seq analysis, the sequence reads that mapped to unique locations in the C. porosus reference genome were used to quantify gene expression, compared between the control and the stimulated groups; and to list the DEGs with a FDR threshold of ≤0.05. Using a fold change threshold (±1.5), the total numbers of both up- and downregulated DEGs for dsRNA and dsDNA were calculated. A significantly greater number of upregulated genes than the number of downregulated genes was observed in both stimulated groups at all time points (Figure 2A). Comparison between the DEG profiles of dsRNA stimulated LV-1 cells revealed that 768 genes were co-upregulated in both the 24 and 48 h time points (Table S5 in Supplementary Material), indicating that the majority of genes induced by dsRNA viral mimics continue to increase in expression 24 hps (Figure 2B; Tables S5–S7 in Supplementary Material for gene lists corresponding to individual subgroupings within Figure 2B), with a total of 895 more DEGs being upregulated at 48 h compared to 24 hps (Table S11 in Supplementary Material for gene lists corresponding to individual subgroupings within Figure 2B). This pattern was also observed for downregulated genes following dsRNA viral mimic stimulation, with the highest number of genes being downregulated at the 48 h time point only (Figures 2C, 572 DEGs). Considerably lower numbers of genes (634 DEGs) were upregulated at 48 hps in dsDNA-stimulated LV-1 cells compared to stimulation with dsRNA at any time point; however, comparable levels of genes were downregulated at least 1.5-fold following stimulation with either dsDNA or dsRNA for 48 h (406 DEGs versus 572 DEGs, respectively, Figure 2C). Interestingly, 218 DEGs were uniquely upregulated and 192 DEGs uniquely downregulated in response to dsDNA viral mimics in the LV-1 cells, with 201 and 15 commonly DEGs up- or downregulated, respectively, between dsDNA and dsRNA viral mimics (Figures 2B,C, respectively) (all gene tables can be seen in Tables S2–S4 in Supplementary Material and Tables S12–S17 in Supplementary Material for gene lists corresponding to individual subgroupings within Figures 2B,C).

FIGURE 2
www.frontiersin.org

Figure 2. Overview of the RNA sequencing data under a time-course stimulation of LV-1 cells with dsRNA and dsDNA. LV-1 cells were harvested at 24 and 48 h post stimulation (hps). (A) Using a fold change threshold of ≥1.5 up or downregulated differentially expressed genes (DEGs) were identified from the comparison among control with dsRNA- and dsDNA-stimulated cells (DEGs were identified based on a false discovery rate q-value threshold of less than 0.05). (B,C) Venn diagrams showing of overlapping comparison of DEG profiles for dsRNA and dsDNA. The mRNA differential expressions in dsRNA- and dsDNA-stimulated LV-1 cells were depicted in three overlapping circles for 1.5-fold up (B) and down (C) regulation at 24 and 48 hps. (D) MA plots showing DEGs for dsRNA- and dsDNA-stimulated LV-1 cells. The y-axis represents the log fold change observed for each mRNA transcript and the y-axis represents the average-log expression values for each transcript. The upregulated genes are in red color, whereas downregulated genes are in green. Data for the genes that were not classified as differentially expressed are plotted as black. (E) Volcano plots showing DEGs for dsRNA- and dsDNA-stimulated LV-1 cells. The x-axis represents the log2 fold change observed for each mRNA transcript and the y-axis represents the log10 value of p values of the significant test between replicates for each transcript. Top 100 DEGs are highlighted in volcano plots.

To illustrate the gene expression pattern across all individual biological samples at various time points within the stimulated groups, volcano plots (Figure 2D), and MA (mean average) plots (Figure 2E) were generated. There was a very good correlation between the fold change differences and p-value and/or average-log expression of the DEGs. In addition, DEGs were further visualized by a hierarchical clustering and an MDS (multi-dimensional scaling) plot was generated (Figures S3A,B in Supplementary Material). These data indicate a good clustering of samples according to the levels of similarities in the gene expression patterns, with a clear distinction between the dsRNA- and dsDNA-stimulated LV-1 cells, and importantly the gene expression profiles were able to differentiate between the two main groups of stimulated C. porosus LV-1 cells (Figure 3A in Supplementary Material). Furthermore, dsRNA- and dsDNA-stimulated samples at 24 and 48 hps were distinctly distributed and closely clustered according to the stimuli and time point (Figure 3B in Supplementary Material).

Distinct and Dynamic Changes in Host DEGs in Response to dsRNA and dsDNA Stimulation

To further investigate the role of DEGs in host cells, we focused on the top-ranked 15 up- and downregulated genes from dsRNA- and dsDNA-stimulated cells (Tables 13). There was a significant increase in upregulation of immune-related genes in dsRNA-stimulated LV-1 cells at 48 h in comparison to 24 h (Table 2). Chemokine pathway-associated genes, such as C-X-C motif chemokine 10 (CXCL-10) and CXCL-11 were among the 15 most upregulated genes for dsRNA-stimulated LV-1 cells and were also top-ranked upregulated genes in dsDNA-stimulated cells at 48 hps (Table 3). Importantly, there was remarkable upregulation of innate immune regulatory and antimicrobial genes in both dsRNA- and dsDNA-stimulated LV-1 cells at 48 h, including the Radical S-adenosyl methionine domain-containing protein 2 (RSAD2), a well-documented broad antiviral gene, with immunomodulatory properties [reviewed in Helbig and Beard (33)], IFN-induced protein with tetratricopeptide repeats 5 (IFIT5) and probable ATP-dependent RNA helicase DDX60 (DDX60), both involved in the detection or enhancement of cytosolic RNA (34, 35) and aconitate decarboxylase 1 [ACOD1 or immune responsive gene 1], a relatively unexplored molecule that has been shown to have antimicrobial effects (36). Furthermore, programmed cell death 1 ligand 2 (PDCD1LG2) was found to be one of the most upregulated genes in the dsRNA- and dsDNA-stimulated LV-1 cells, although its expression was higher in dsRNA-stimulated cells than in dsDNA-stimulated cells (11.6 versus 6.58 log2 fold change); PDCD1LG2 is thought to limit viral-induced cellular damage due to overactive T-cells (37). Outside of the top 15 upregulated genes across data sets, many other molecules known to be pivotal to the early innate immune response to viral infection in mammals were commonly upregulated, as can be seen in the TLR and RIG-I pathways shown in Figures S4–S7 in Supplementary Material. Other important immune regulatory genes commonly upregulated were the interferon regulatory factors (IRF)1, 7, and 8; the type I interferon receptors IFNAR1/2; and the type I IFN, IFN-omega. The top 15 downregulated genes differed exclusively between dsRNA and dsDNA viral mimic simulation, with only one common DEG seen between both the stimulated LV-1 cells at the 24 and 48 h time points; a Rho-mediated GTP-binding protein Rho, a protein known to regulate intracellular actin dynamics. As might be expected, the up- and downregulation of gene transcripts stimulated at two different time points by dsRNA fell into three main categories, those that were differentially regulated initially presenting at 48 h only, and those that initially presented at 24 hps and then either continued to increase/decrease in expression, or remained unchanged in their expression between 24 and 48 h (genes demonstrating a significantly increased or decreased expression level between 24 and 48 h of dsRNA stimulation can be found in Tables S18 and S19 in Supplementary Material).

TABLE 1
www.frontiersin.org

Table 1. The top 15 upregulated and downregulated genes for dsRNA stimulated versus control LV-1 cells at 24 h post stimulation as ranked by fold change and false discovery rate (FDR)-adjusted P values of ≤0.05.

TABLE 2
www.frontiersin.org

Table 2. The top 15 upregulated and downregulated genes for dsRNA stimulated versus control LV-1 cells at 48 h post stimulation as ranked by fold change and false discovery rate (FDR)-adjusted P values of ≤0.05.

TABLE 3
www.frontiersin.org

Table 3. The top 15 upregulated and downregulated genes for dsDNA-stimulated versus control LV-1 cells at 48 h post stimulation as ranked by fold change and false discovery rate (FDR)-adjusted P values of ≤0.05.

In order to validate the RNAseq data and to determine the longevity of the host innate response to viral mimic stimulation, we performed a longer time-course experiment, examining genes known to be involved in effective immune responses, and shown to be highly upregulated within the RNAseq data set for both dsRNA and dsDNA stimulation, including ACOD1, CXCL10, RTP2, and OASL. As can be seen in Figures 3A,B, both dsRNA and dsDNA viral mimic stimulation was able to significantly upregulate all selected candidate genes. Interestingly, dsRNA stimulation of LV-1 cells induced maximal expression of candidate genes between 48 and 72 h, with high expression still observed at 96 hps. In comparison, dsDNA stimulation induced a much slower response, with candidate genes showing maximal regulation at 96 h, excepting for CXCL10 (Figure 3B).

FIGURE 3
www.frontiersin.org

Figure 3. Post RNAseq validation with longer timecourse. Quantitative real-time PCR was performed to determine the changes in gene regulation of selected candidate genes. (A,B) Relative expression levels of aconitate decarboxylase 1 (ACOD1), CXCL10, RTP2, and OASL were determined following 48, 72, and 96 h of viral mimic stimulation. The graph shows mean ± SD of values from three replicates per group.

Functional Categorization and Canonical Pathways of DEGs Detected with RNA-seq

To categorize functional networks and identify enriched biological process of GO terms, the DEGs at 24 and 48 hps were imported into the Bioconductor topGO package (38). There were 49 and 50 Biological Processes at 24 and 48 h activated significantly by dsRNA (Classic fisher’s test P-value ≤0.001) (Tables S20 and S21 in Supplementary Material). Figures 4A,B represents the list of activated Biological Processes involved in immunological function after LV-1 cells were stimulated with dsRNA at 24 and 48 hps, respectively. Among the top-ranked (Classic fisher’s test P-value ≤0.001) Biological Processes, cytokine-mediated signaling pathway, CD4/CD8-positive alpha-beta lineage commitment, positive regulation of cytokine production, positive regulation of lymphocyte activation, immune response, and defense response to other organism were activated both at 24 and 48 hps with dsRNA (Figures 4A,B). In addition, at 24 hps with dsRNA, there were other Biological Processes that were significantly over-represented including I-kappaB kinase/NF-kappaB (IKK/NF-kB) signaling, inflammatory response, response to interferon-beta, negative regulation of viral life cycle, and positive regulation of defense response to virus (Figure 4A). The significantly over-represented Biological Processes involved in immunological functions were relatively similar at both time points following dsRNA stimulation, excepting the presence of gene sets within the negative regulation of I-kappaB kinase/NF-kappaB (IKK/NF-kB) signaling at 48 hps with dsRNA (Figure 4B).

FIGURE 4
www.frontiersin.org

Figure 4. The top-ranked biological processes gene ontology functions identified by topGO package. Pie charts based on the DE genes involved in immunological function of the enriched biological processes generated using differentially expressed genes (DEGs) at 24 hps (A) and 48 hps (B). The values for each function represent the ratio of DEGs versus the total annotated genes for each functional category.

Analysis of the dsDNA viral mimic stimulation of LV-1 cells showed only 28 Biological Processes activated significantly (Classic fisher’s test P-value ≤0.001) at 48 hps with dsDNA (Table S22 in Supplementary Material). Among these, only three Biological Processes were related to immunological function, including, cellular response to virus, defense response to virus, and negative regulation of viral process, all known to have direct and/or indirect influence on the immune system.

The enriched DEGs at both post-stimulated time points were also analyzed using KEGG Mapper (version 2.8) using the American Alligator (Alligator mississippiensis; amj) as a search model to identify canonical pathways. DEGs due to dsRNA and dsDNA stimulation were indicated as either red (upregulated) or blue (downregulated), with at least 30 hits required to consider the canonical pathways enriched. In this study, we identified 54 and 72 canonical pathways that were significantly enriched at 24 and 48 hps with dsRNA, respectively (Tables S23 and S24 in Supplementary Material), whereas 66 canonical pathways were identified at 48 hps with dsDNA (Table S25 in Supplementary Material). It is notable that a large number of DEGs were not found using the amj search model, although there were six separate canonical pathways identified at both time-course stimulations with dsRNA or dsDNA, which have immunological functions; five of these were common to all stimulation time points (highlighted as red in Tables S23–S25 in Supplementary Material). These canonical pathways included the Toll-like receptor signaling pathway which is shown to overlay with the DE gene expression results in Figure 5 (Figures S4 and S5 in Supplementary Material) and the RIG-I-like receptor signaling pathway which is presented in Figure 6 (Figures S6 and S7 in Supplementary Material). The notable differences included higher upregulation of key downstream signaling molecules such as activators of receptors (JAK-STAT signaling pathway) and activator of the TLR3 viral sensing (MyD88-independent pathway) in dsRNA- than dsDNA-stimulated LV-1 cells at 48 h; as well as the downstream signaling for the indirect P13K–Akt pathway being upregulated at 24 h following dsDNA stimulation (Figures S5 in Supplementary Material), whereas it was downregulated in the presence of dsRNA at 48 h (Figure 5).

FIGURE 5
www.frontiersin.org

Figure 5. The enriched canonical pathway for toll-like receptor signaling at 48 hps with dsRNA. Pathway analysis using KEGG Mapper allowed us to identify the pathways that were differentially expressed between dsRNA-stimulated and non-stimulated LV-1 cells. Red and blue shading indicates increased and decreased expression, respectively, in dsRNA-stimulated LV-1 cells relative to the non-stimulated control cells. White and green shading indicates non-expression and non-differential expression, respectively. Solid and dashed lines represent direct and indirect interactions, respectively.

FIGURE 6
www.frontiersin.org

Figure 6. The enriched canonical pathway for RIG-I like receptor signaling at 48 hps with dsRNA. Pathway analysis using KEGG Mapper allowed us to identify the pathways that were differentially expressed between dsRNA-stimulated and non-stimulated LV-1 cells. Red and blue shading indicates increased and decreased expression, respectively, in dsRNA-stimulated LV-1 cells relative to the non-stimulated control cells. White and green shading indicates non-expression and non-differential expression, respectively. Solid and dashed lines represent direct and indirect interactions, respectively.

We also further investigated the RIG-I-like receptor signaling pathway for dsRNA- and dsDNA-stimulated LV-1 cells at various time points. The RIG-I- and Mda5-mediated downstream signaling pathways, leading to the production of inflammatory cytokines and Type I IFN were activated to a greater extent in dsRNA stimulated LV-1 cells at 48 h (Figure 6) than at 24 h (Figures S6 in Supplementary Material), although the adaptor protein, IPS-1, which is known to be involved in RIG-I- and Mda5-mediated antiviral immune responses was downregulated at 48 hps.

Discussion

In recent years, transcriptome profiling has enabled us to generate an unprecedented global view of the extent and complexity of gene transcription for a number of eukaryotic species, and has revealed a better understanding of the host genes and cellular pathways that are activated and perturbed in response to viral pathogens (3943). Recent viral infections reported in the crocodile are novel, and do not have tissue culture models as yet. This study was aimed at describing the underlying host immune response against viral pathogens in this ancient animal using viral mimics, and the results have been demonstrated for the first time that stimulation of a primary passaged cell line from a crocodile (LV-1) is able to functionally upregulate multiple innate immune pathways, using RNA-seq analysis.

Innate immunity is a hosts first line of defense against viral infection, and is induced rapidly following detection of the pathogen. This response can include the production of antimicrobial peptides (44), complement, and lectins (45), however, the recognition of foreign viral products by PRRs, and the subsequent production of IFN remains the most potent host response to viral infection (46). Our knowledge of the signaling pathways following activation of the PRRs, as well as the downstream antiviral effector molecules, ISGs, induced by these pathways is based mostly on mammals, mainly humans and mice. Major explorations of antiviral innate immune responses in non-mammalian species remain very limited. There have been a number of studies in non-mammalian vertebrate species, such as the amphibians and the birds, describing the presence of various PRRs and IFNs in the genome, as well the potential role of varying IFN types; however, there is very little work analyzing whether these gene products are expressed in a functional context, or describing their antiviral effector genes for most species [reviewed in Chen et al. (47)]; with the exception of the chicken, which is an important agriculturally species (47, 48). The reptiles represent a bridge between the fish and mammals, and are the only ectothermic amniotes; this positions them as a pivotal species to not only enhance our evolutionary knowledge of the early innate systems controlling viral infection but also to provide insight into the complicated and intricate mammalian system that controls early viral infections.

The saltwater crocodile is the world’s largest living reptile species, and is evolved from the archosauria clade, which also includes the dinosaurs, pterosaurs, crocodilians, and birds; the latter two being the only living members (13). There a number of reported DNA and RNA viral infections of the crocodile, and this work describes the use of both DNA and RNA viral mimics at optimal time points, to stimulate the early innate responses in a recently derived crocodile cell line LV-1 (49). Many of the viruses infecting crocodiles are novel, and do not have tissue culture models as yet, likewise, the use of viral mimics ensured that we did not have to account for the ability of viruses to evade the innate immune system and disable varying arms of host protection (50). Using an RNA-seq approach and the recently annotated C. porosus genome as a reference, we were able to determine that the C. porosus host response to viral DNA and RNA mimics induced the expression of 634 and 1,663 genes respectively, with 406 and 572 genes being downregulated at 48 hps. There was considerable overlap between the two gene sets at 48 h, with 25 and 37% of genes up- or downregulated following dsRNA stimulation, being regulated similarly with dsDNA stimulation. Of these genes, approximately 10% belonged to uncharacterized open reading frames in the crocodile genome, currently unannotated, demonstrating the potential novelty of host proteins that may be involved in viral control in these animals.

RNA sequencing analysis was performed at two time points following RNA viral mimic stimulation, an early time point of 24 h, as well as the 48 h time point described above. These time points were chosen due to the early upregulation of initial ISGs tested prior to the RNA-seq experiment, in comparison to a distinct lack of ISG expression in the DNA viral mimic-stimulated cells at this point (Figure 1); in addition, ISGs have been demonstrated in mammalian cells previously to have distinct expression kinetics, with many genes being up- and downregulated within a 24 h period (51). In the context of the highest upregulated DEG’s, there was little difference between the 24 and 48 h stimulation in the LV-1 cells (Tables 1 and 2); however, there was only one common gene found to be downregulated in both time points following RNA viral mimic stimulation, Rho-related GTP-binding protein RhoB, a gene involved in intracellular protein trafficking, that has also previously been shown to be involved in enhanced entry of a select group of viruses (52).

The C. porosus genome annotation is an on-going project, as such much of the transcriptome remains unannotated, and currently there is no information regarding GO. Consequently, to understand the critical role of the DEGs in the biological processes, we annotated GO terms using ortholog genes of Gallus gallus (chicken) for 8,399 genes and human ortholog genes for another 2,928. As expected, GO-term analysis revealed a marked increase in biological processes associated with immunological cell functions in all groups stimulated with nucleic acid viral mimics, including “cellular response to virus,” “defence response to virus,” and “immune response,” with general immunological processes being over-represented in all sample sets (Tables S20–22 in Supplementary Material). Interestingly, the biological processes involved in the set of DEGs from the DNA viral mimic stimulated cells demonstrated significantly less pathways common to immune processes in comparison to the RNA viral mimic-stimulated cells. However, in terms of pathway analysis, utilizing KEGG mapper, the canonical pathways found to be significantly represented in both DNA and RNA viral mimic-stimulated cells were remarkably similar (Tables S23–25 in Supplementary Material).

Many common innate immune pathways known to be involved in the recognition of both viral dsDNA and dsRNA in mammals in the cytosol were all shown to be significantly represented in the crocodile cell lines stimulated with the viral replication mimics. Canonical pathway analysis demonstrated the MAPK signaling pathway to be significantly represented in all three datasets, being in the top five canonical pathways (Tables S23–25 in Supplementary Material). This signaling pathway has previously been shown to be pivotal in the upregulation of inflammatory genes following activation of both dsRNA and dsDNA signaling pathways after host viral recognition (53, 54). RNA viral mimic stimulation was shown to activate gene sets present in both the RIG-I and TLR3 signaling pathways, the predominant pathways known to sense viral RNA (54). All pathway members required to induce IFN production were positively regulated, with the exception of IRF3 and IPS-1, which was downregulated at the 48 h RNA time point. Both of these proteins are known to be present at high basal levels, and are activated via phosphorylation alone, or phosphorylation and ubiquitination, respectively, following viral sensing by RIG and TLR3 (5557). Stimulation of LV-1 cells with the DNA viral mimics also induced gene sets common to the RIG-I and TLR signaling pathways, although as many of these genes are inducible via type I IFN production (discussed below), a likely result of DNA stimulation, it is not surprising that they are differentially regulated. Although essentially an RNA sensor of 5′ tri-phosphate ssRNA, RIG-I has previously been reported to detect viral DNA following modifications via RNA polymerase III (58) in human and murine cell lines, and it is possible that non-mammalian host cells may also perform this function, but further work will need to be done to confirm this.

No prior work has been performed to identify functional dsDNA receptors in non-mammalian vertebrates, however the C. porosus genome contains both cGAS and DDX41, two known mammalian receptors for dsDNA sensing in the cytosol [reviewed in Paludan and Bowie (59)], as well as TLR21, which is known to act as a functional homolog of TLR9 in the chicken (60). Interestingly, the crocodile genome appears to lack an open reading frame for STING (or MITA), the signaling adaptor that has been shown to be central to the upregulation of the antiviral cytokines, type I and III IFN, following viral DNA sensing by all dsDNA receptors, even though it is present in the very closely related genomes of the American alligator and the gharial [reviewed in Burdette and Vance (61)]. However, presumably TLR21 signals through a similar pathway to the other endosomal TLRs, which appear to be represented in the significant canonical pathways following both RNA and DNA viral mimic stimulation (Tables S23–25 in Supplementary Material). Further detailed studies will be required to elucidate the functional source of ISG induction following DNA stimulation in the crocodile, and the specific pathways involved.

The activation of PRR pathways via the sensing of viral nucleic acid, results in the production of the antiviral cytokines, termed IFNs. It is these IFNs that are known to induce the upregulation of hundreds of ISGs via the JAK–STAT pathway, in response to viral infection, many of which are known to be potently antiviral in the context of mammals (3, 4). The C. porosus genome contains one documented open reading frame with a similarity to the type I IFN-omega. This gene product was significantly upregulated during stimulation with RNA viral mimics at both the early and late time points (Tables S2 and S3 in Supplementary Material), and is a probable candidate for activation of the JAK/STAT pathway in the crocodile, via the type I INFA receptors. There were also multiple highly upregulated known antiviral ISGs present in both RNA and DNA viral mimic stimulated LV-1 cells, including viperin (RSAD-2), MX1, MOV10, CH25H, ADAR, DDX58/60, IRF1/7, ISG15 and TRIM5/25 [reviewed in Schoggins and Rice (4)].

This work is demonstrated for the first time that stimulation of a primary passaged cell line from a crocodile, or perhaps any reptile, is able to functionally upregulate multiple PPR pathway members, and, subsequently, induce a large sub-set of genes with potential antiviral function in response to viral mimics. We have also shown that multiple conserved canonical pathways are likely in play in the host response to viral infection in this reptile, with many of the viral-responsive genes being novel (approximately 10%). Further work needs to be performed to determine the functionality of potential important players in the antiviral innate host response of the reptile to a viral infection. This work will underpin these studies, and offers a better understanding of these pathways and the effector genes responsible for control of viral infection in non-mammalian species. An enhanced knowledge of these ancient antiviral pathways will not only add to our understanding of the host antiviral innate response in non-mammalian species, but is critical to fully comprehend the complexity of the mammalian innate immune response to viral infection.

Author Contributions

SS and KH conceived and designed the experiments. SS performed the experiments. BW-S performed the post validation of RNAseq data using some selected candidate genes. SS and YW analyzed the data. SS and KH contributed reagents and materials. SS and KH wrote the paper. All authors read and approved the final manuscript.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

The authors would like to thank the Berrimah Veterinary Laboratory (Department of Primary Industry and Resources, Government of Northern Territory, Australia) for providing the C. porosus liver cell line (LV-1) used in this experiment. The authors also like to thank La Trobe University for the financial support under Start Up Grant to SS and KH.

Supplementary Material

The Supplementary Material for this article can be found online at http://www.frontiersin.org/article/10.3389/fimmu.2017.01634/full#supplementary-material.

References

1. Riera Romo M, Perez-Martinez D, Castillo Ferrer C. Innate immunity in vertebrates: an overview. Immunology (2016) 148(2):125–39. doi:10.1111/imm.12597

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Chan YK, Gack MU. Viral evasion of intracellular DNA and RNA sensing. Nat Rev Microbiol (2016) 14(6):360–73. doi:10.1038/nrmicro.2016.45

CrossRef Full Text | Google Scholar

3. Sen GC, Sarkar SN. The interferon-stimulated genes: targets of direct signaling by interferons, double-stranded RNA, and viruses. In: Pitha PM, editor. Interferon: The 50th Anniversary. Berlin, Heidelberg: Springer (2007). p. 233–50.

Google Scholar

4. Schoggins JW, Rice CM. Interferon-stimulated genes and their antiviral effector functions. Curr Opin Virol (2011) 1(6):519–25. doi:10.1016/j.coviro.2011.10.008

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Hoffmann H-H, Schneider WM, Rice CM. Interferons and viruses: an evolutionary arms race of molecular interactions. Trends Immunol (2015) 36(3):124–38. doi:10.1016/j.it.2015.01.004

CrossRef Full Text | Google Scholar

6. Qi Z, Nie P, Secombes CJ, Zou J. Intron-containing type I and type III IFN coexist in amphibians: refuting the concept that a retroposition event gave rise to type I IFNs. J Immunol (2010) 184(9):5038–46. doi:10.4049/jimmunol.0903374

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Grayfer L, De Jesús Andino F, Robert J. The amphibian (Xenopus laevis) type I interferon response to frog virus 3: new insight into ranavirus pathogenicity. J Virol (2014) 88(10):5766–77. doi:10.1128/JVI.00223-14

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Ding Y, Ao J, Huang X, Chen X. Identification of two subgroups of type I IFNs in perciforme fish large yellow croaker Larimichthys crocea provides novel insights into function and regulation of fish type I IFNs. Front Immunol (2016) 7:343. doi:10.3389/fimmu.2016.00343

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Chen SN, Zhang XW, Li L, Ruan BY, Huang B, Huang WS, et al. Evolution of IFN-λ in tetrapod vertebrates and its functional characterization in green anole lizard (Anolis carolinensis). Dev Comp Immunol (2016) 61:208–24. doi:10.1016/j.dci.2016.04.004

CrossRef Full Text | Google Scholar

10. Priyam M, Tripathy M, Rai U, Ghorai SM. Tracing the evolutionary lineage of pattern recognition receptor homologues in vertebrates: an insight into reptilian immunity via de novo sequencing of the wall lizard splenic transcriptome. Vet Immunol Immunopathol (2016) 172:26–37. doi:10.1016/j.vetimm.2016.03.002

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Milic NL, Davis S, Carr JM, Isberg S, Beard MR, Helbig KJ. Sequence analysis and characterisation of virally induced viperin in the saltwater crocodile (Crocodylus porosus). Dev Comp Immunol (2015) 51(1):108–15. doi:10.1016/j.dci.2015.03.001

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Marschang RE. Viruses infecting reptiles. Viruses (2011) 3(11):2087. doi:10.3390/v3112087

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Christopher AB. Phylogenetic approaches toward Crocodylian history. Annu Rev Earth Planet Sci (2003) 31:357–97. doi:10.1146/annurev.earth.31.100901.141308

CrossRef Full Text | Google Scholar

14. McCowan C, Shepherdley C, Slocombe RF. Herpesvirus-like particles in the skin of a saltwater crocodile (Crocodylus porosus). Aust Vet J (2004) 82(6):375–7. doi:10.1111/j.1751-0813.2004.tb11109.x

CrossRef Full Text | Google Scholar

15. Shilton CM, Jerrett IV, Davis S, Walsh S, Benedict S, Isberg SR, et al. Diagnostic investigation of new disease syndromes in farmed Australian saltwater crocodiles (Crocodylus porosus) reveals associations with herpesviral infection. J Vet Diagn Invest (2016) 28(3):279–90. doi:10.1177/1040638716642268

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Buenviaje GN, Ladds PW, Martin Y. Pathology of skin diseases in crocodiles. Aust Vet J (1998) 76(5):357–63. doi:10.1111/j.1751-0813.1998.tb12368.x

PubMed Abstract | CrossRef Full Text | Google Scholar

17. St John JA, Braun EL, Isberg SR, Miles LG, Chong AY, Gongora J, et al. Sequencing three crocodilian genomes to illuminate the evolution of archosaurs and amniotes. Genome Biol (2012) 13(1):415. doi:10.1186/gb-2012-13-1-415

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Juranic Lisnic V, Babic Cac M, Lisnic B, Trsan T, Mefferd A, Das Mukhopadhyay C, et al. Dual analysis of the murine cytomegalovirus and host cell transcriptomes reveal new aspects of the virus-host cell interface. PLoS Pathog (2013) 9(9):e1003611. doi:10.1371/journal.ppat.1003611

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Ertl R, Klein D. Transcriptional profiling of the host cell response to feline immunodeficiency virus infection. Virol J (2014) 11:52. doi:10.1186/1743-422x-11-52

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Jones M, Dry IR, Frampton D, Singh M, Kanda RK, Yee MB, et al. RNA-seq analysis of host and viral gene expression highlights interaction between varicella zoster virus and keratinocyte differentiation. PLoS Pathog (2014) 10(1):e1003896. doi:10.1371/journal.ppat.1003896

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Mei B, Ding X, Xu HZ, Wang MT. Global gene expression changes in human peripheral blood after H7N9 infection. Gene (2014) 551(2):255–60. doi:10.1016/j.gene.2014.08.062

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Veneman WJ, de Sonneville J, van der Kolk K-J, Ordas A, Al-Ars Z, Meijer AH, et al. Analysis of RNAseq datasets from a comparative infectious disease zebrafish model using GeneTiles bioinformatics. Immunogenetics (2015) 67(3):135–47. doi:10.1007/s00251-014-0820-3

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Andrews S. FastQC: A Quality Control Tool for High Throughput Sequence. (2016). Available: http://www.bioinformatics.babraham.ac.uk/projects/fastqc

Google Scholar

24. Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics (2014) 30(15):2114–20. doi:10.1093/bioinformatics/btu170

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Kim D, Langmead B, Salzberg SL. HISAT: a fast spliced aligner with low memory requirements. Nat Methods (2015) 12(4):357–60. doi:10.1038/nmeth.3317

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Wang L, Wang S, Li W. RSeQC: quality control of RNA-seq experiments. Bioinformatics (2012) 28(16):2184–5. doi:10.1093/bioinformatics/bts356

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Liao Y, Smyth GK, Shi W. The Subread aligner: fast, accurate and scalable read mapping by seed-and-vote. Nucleic Acids Res (2013) 41(10):e108. doi:10.1093/nar/gkt214

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Liao Y, Smyth GK, Shi W. featureCounts: an efficient general purpose program for assigning sequence reads to genomic features. Bioinformatics (2014) 30(7):923–30. doi:10.1093/bioinformatics/btt656

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Robinson MD, McCarthy DJ, Smyth GK. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics (2010) 26(1):139–40. doi:10.1093/bioinformatics/btp616

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Robinson MD, Oshlack A. A scaling normalization method for differential expression analysis of RNA-seq data. Genome Biol (2010) 11(3):R25. doi:10.1186/gb-2010-11-3-r25

CrossRef Full Text | Google Scholar

31. Law CW, Chen Y, Shi W, Smyth GK. voom: precision weights unlock linear model analysis tools for RNA-seq read counts. Genome Biol (2014) 15(2):R29. doi:10.1186/gb-2014-15-2-r29

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Moriya Y, Itoh M, Okuda S, Yoshizawa AC, Kanehisa M. KAAS: an automatic genome annotation and pathway reconstruction server. Nucleic Acids Res (2007) 35(Web Server issue):W182–5. doi:10.1093/nar/gkm321

CrossRef Full Text | Google Scholar

33. Helbig KJ, Beard MR. The role of viperin in the innate antiviral response. J Mol Biol (2014) 426(6):1210–9. doi:10.1016/j.jmb.2013.10.019

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Miyashita M, Oshiumi H, Matsumoto M, Seya T. DDX60, a DEXD/H box helicase, is a novel antiviral factor promoting RIG-I-like receptor-mediated signaling. Mol Cell Biol (2011) 31(18):3802–19. doi:10.1128/mcb.01368-10

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Katibah GE, Qin Y, Sidote DJ, Yao J, Lambowitz AM, Collins K. Broad and adaptable RNA structure recognition by the human interferon-induced tetratricopeptide repeat protein IFIT5. Proc Natl Acad Sci U S A (2014) 111(33):12025–30. doi:10.1073/pnas.1412842111

CrossRef Full Text | Google Scholar

36. El-Zaatari M, Kao JY. Role of dietary metabolites in regulating the host immune response in gastrointestinal disease. Front Immunol (2017) 8(51). doi:10.3389/fimmu.2017.00051

CrossRef Full Text | Google Scholar

37. Sharpe AH, Wherry EJ, Ahmed R, Freeman GJ. The function of programmed cell death 1 and its ligands in regulating autoimmunity and infection. Nat Immunol (2007) 8(3):239–45. doi:10.1038/ni1443

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Alexa A, Rahnenfuhrer J. topGO: Enrichment Analysis for Gene Ontology. (2016).

Google Scholar

39. Wang Z, Gerstein M, Snyder M. RNA-Seq: a revolutionary tool for transcriptomics. Nat Rev Genet (2009) 10:57–63. doi:10.1038/nrg2484

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Glennon NB, Jabado O, Lo MK, Shaw ML. Transcriptome profiling of the virus-induced innate immune response in pteropus vampyrus and its attenuation by Nipah virus interferon antagonist functions. J Virol (2015) 89(15):7550–66. doi:10.1128/jvi.00302-15

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Mehrbod P, Harun MS, Shuid AN, Omar AR. Transcriptome analysis of feline infectious peritonitis virus infection. Methods Mol Biol (2015) 1282:241–50. doi:10.1007/978-1-4939-2438-7_20

CrossRef Full Text | Google Scholar

42. Park SJ, Kumar M, Kwon HI, Seong RK, Han K, Song JM, et al. Dynamic changes in host gene expression associated with H5N8 avian influenza virus infection in mice. Sci Rep (2015) 5:16512. doi:10.1038/srep16512

CrossRef Full Text | Google Scholar

43. Lamontagne J, Mell JC, Bouchard MJ. Transcriptome-wide analysis of hepatitis B virus-mediated changes to normal hepatocyte gene expression. PLoS Pathog (2016) 12(2):e1005438. doi:10.1371/journal.ppat.1005438

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Klotman ME, Chang TL. Defensins in innate antiviral immunity. Nat Rev Immunol (2006) 6(6):447–56. doi:10.1038/nri1860

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Mason CP, Tarr AW. Human lectins and their roles in viral infections. Molecules (2015) 20(2):2229–71. doi:10.3390/molecules20022229

CrossRef Full Text | Google Scholar

46. Sen GC. Viruses and interferons. Annu Rev Microbiol (2001) 55(1):255–81. doi:10.1146/annurev.micro.55.1.255

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Chen S, Cheng A, Wang M. Innate sensing of viruses by pattern recognition receptors in birds. Vet Res (2013) 44:82. doi:10.1186/1297-9716-44-82

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Santhakumar D, Rubbenstroth D, Martinez-Sobrido L, Munir M. Avian interferons and their antiviral effectors. Front Immunol (2017) 8:49. doi:10.3389/fimmu.2017.00049

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Lorna M, Steven D, Cathy S, Sally I, Amanda C, Jaime G. Hunting Viruses in Crocodiles. Barton, ACT: RIRDC Publication (2012).

Google Scholar

50. Bowie AG, Unterholzner L. Viral evasion and subversion of pattern-recognition receptor signalling. Nat Rev Immunol (2008) 8(12):911–22. doi:10.1038/nri2436

PubMed Abstract | CrossRef Full Text | Google Scholar

51. Bolen CR, Ding S, Robek MD, Kleinstein SH. Dynamic expression profiling of type I and type III interferon-stimulated hepatocytes reveals a stable hierarchy of gene expression. Hepatology (2014) 59(4):1262–72. doi:10.1002/hep.26657

CrossRef Full Text | Google Scholar

52. Quinn K, Brindley MA, Weller ML, Kaludov N, Kondratowicz A, Hunt CL, et al. Rho GTPases modulate entry of Ebola virus and vesicular stomatitis virus pseudotyped vectors. J Virol (2009) 83(19):10176–86. doi:10.1128/jvi.00422-09

PubMed Abstract | CrossRef Full Text | Google Scholar

53. Barber GN. Innate immune DNA sensing pathways: STING, AIMII and the regulation of interferon production and inflammatory responses. Curr Opin Immunol (2011) 23(1):10–20. doi:10.1016/j.coi.2010.12.015

PubMed Abstract | CrossRef Full Text | Google Scholar

54. Jensen S, Thomsen AR. Sensing of RNA viruses: a review of innate immune receptors involved in recognizing RNA virus invasion. J Virol (2012) 86(6):2900–10. doi:10.1128/jvi.05738-11

PubMed Abstract | CrossRef Full Text | Google Scholar

55. Panne D, McWhirter SM, Maniatis T, Harrison SC. Interferon regulatory factor 3 is regulated by a dual phosphorylation-dependent switch. J Biol Chem (2007) 282(31):22816–22. doi:10.1074/jbc.M703019200

PubMed Abstract | CrossRef Full Text | Google Scholar

56. Castanier C, Zemirli N, Portier A, Garcin D, Bidere N, Vazquez A, et al. MAVS ubiquitination by the E3 ligase TRIM25 and degradation by the proteasome is involved in type I interferon production after activation of the antiviral RIG-I-like receptors. BMC Biol (2012) 10:44. doi:10.1186/1741-7007-10-44

PubMed Abstract | CrossRef Full Text | Google Scholar

57. Liu S, Cai X, Wu J, Cong Q, Chen X, Li T, et al. Phosphorylation of innate immune adaptor proteins MAVS, STING, and TRIF induces IRF3 activation. Science (2015) 347(6227):aaa2630. doi:10.1126/science.aaa2630

PubMed Abstract | CrossRef Full Text | Google Scholar

58. Ablasser A, Bauernfeind F, Hartmann G, Latz E, Fitzgerald KA, Hornung V. RIG-I dependent sensing of poly(dA-dT) via the induction of an RNA polymerase III transcribed RNA intermediate. Nat Immunol (2009) 10:1065–72. doi:10.1038/ni.1779

CrossRef Full Text | Google Scholar

59. Paludan SR, Bowie AG. Immune sensing of DNA. Immunity (2013) 38(5):870–80. doi:10.1016/j.immuni.2013.05.004

CrossRef Full Text | Google Scholar

60. Brownlie R, Zhu J, Allan B, Mutwiri GK, Babiuk LA, Potter A, et al. Chicken TLR21 acts as a functional homologue to mammalian TLR9 in the recognition of CpG oligodeoxynucleotides. Mol Immunol (2009) 46(15):3163–70. doi:10.1016/j.molimm.2009.06.002

PubMed Abstract | CrossRef Full Text | Google Scholar

61. Burdette DL, Vance RE. STING and the innate immune response to nucleic acids in the cytosol. Nat Immunol (2013) 14(1):19–26. doi:10.1038/ni.2491

PubMed Abstract | CrossRef Full Text | Google Scholar

Keywords: crocodile, RNA-sequencing, transcriptome analysis, virus, immune response, innate immunity, interferon, reptile

Citation: Sarker S, Wang Y, Warren-Smith B and Helbig KJ (2017) Dynamic Changes in Host Gene Expression following In Vitro Viral Mimic Stimulation in Crocodile Cells. Front. Immunol. 8:1634. doi: 10.3389/fimmu.2017.01634

Received: 29 August 2017; Accepted: 09 November 2017;
Published: 22 November 2017

Edited by:

Thomas A. Kufer, University of Hohenheim, Germany

Reviewed by:

Baubak Bajoghli, Universitätsklinikum Tübingen, Germany
Junji Xing, Houston Methodist Research Institute, United States

Copyright: © 2017 Sarker, Wang, Warren-Smith and Helbig. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) or licensor are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Karla J. Helbig, k.helbig@latrobe.edu.au