key: cord-0788978-ewbtb3yc
authors: Davy, Christina M.; Donaldson, Michael E.; Subudhi, Sonu; Rapin, Noreen; Warnecke, Lisa; Turner, James M.; Bollinger, Trent K.; Kyle, Christopher J.; Dorville, Nicole A. S.-Y.; Kunkel, Emma L.; Norquay, Kaleigh J. O.; Dzal, Yvonne A.; Willis, Craig K. R.; Misra, Vikram
title: White-nose syndrome is associated with increased replication of a naturally persisting coronaviruses in bats
date: 2018-10-19
journal: Sci Rep
DOI: 10.1038/s41598-018-33975-x
sha: d495854289e2ffd80ba9ac86688abe6c7577b545
doc_id: 788978
cord_uid: ewbtb3yc

Spillover of viruses from bats to other animals may be associated with increased contact between them, as well as increased shedding of viruses by bats. Here, we tested the prediction that little brown bats (Myotis lucifugus) co-infected with the M. lucifugus coronavirus (Myl-CoV) and with Pseudogymnoascus destructans (Pd), the fungus that causes bat white-nose syndrome (WNS), exhibit different disease severity, viral shedding and molecular responses than bats infected with only Myl-CoV or only P. destructans. We took advantage of the natural persistence of Myl-CoV in bats that were experimentally inoculated with P. destructans in a previous study. Here, we show that the intestines of virus-infected bats that were also infected with fungus contained on average 60-fold more viral RNA than bats with virus alone. Increased viral RNA in the intestines correlated with the severity of fungus-related pathology. Additionally, the intestines of bats infected with fungus exhibited different expression of mitogen-activated protein kinase pathway and cytokine related transcripts, irrespective of viral presence. Levels of coronavirus antibodies were also higher in fungal-infected bats. Our results suggest that the systemic effects of WNS may down-regulate anti-viral responses in bats persistently infected with M. lucifugus coronavirus and increase the potential of virus shedding.

Quantitation of Myl-CoV and M. lucifugus RNA through reverse transcription quantitative PCR (RT-qPCR) and dual-RNA-sequencing indicated that co-infected bats had significantly higher levels of Myl-CoV RNA than bats infected with virus alone. The amount of Myl-CoV RNA correlated with the severity of WNS pathology in co-infected bats. This phenomenon was associated with specific molecular responses to co-infection, even in the intestines of bats where only one of the two pathogens was directly interacting with the host tissue. The levels of antibodies against Myl-CoV nucleocapsid (N) protein were also higher in co-infected bats. Each key result is discussed in detail below.

Bats co-infected with the fungus P. destructans and the virus Myl-CoV contained higher levels of Myl-CoV RNA. Myl-CoV genomic RNA was detected in bats infected only with Myl-CoV (virus-infected; 7/18), co-infected bats (European P. destructans (3/13) , or with North American P. destructans (7/16) 37 ). There was no difference in the frequency of Myl-CoV detected among these treatments (p-value = 0.801). We pooled bats infected with the two P. destructans isolates for all further analyses and tested whether co-infection with P. destructans and Myl-CoV correlated with an increase in viral replication. Our RT-qPCR data showed that the co-infected bats contained 60-fold more Myl-CoV RNA on average than the virus-infected bats (Mann Whitney test; p-value = 0.014; Fig. 1 ).

Relative quantities of Myl-CoV RNA detected in the ileum of the virus-infected bats were low and showed low variation (Standard Deviation of 1/ΔCT = 0.005), compared to the relative quantities of Myl-CoV RNA in the co-infected bats (Standard Deviation of 1/ΔCT = 0.108; Fig. 1 ). The severity of WNS fungal pathology varied in co-infected bats, and we therefore tested whether relative quantities of viral RNA in the ileum correlated with the severity of WNS symptoms. Levels of WNS severity were scored based on fungal hyphae on the wings, secondary bacteria in wing lesions, oedema, necrosis and inflammation in wing lesions, and levels of neutrophils in lung, spleen and liver interstitium. Severity scores for wing tissue, secondary bacteria in lesions, and neutrophils in the lung interstitium positively correlated with relative amounts of coronavirus RNA in hibernating bats (Table 1) .

To determine the extent to which Myl-CoV and P. destructans infection interact to influence gene expression in bat intestines, we performed a transcriptomic analysis on bat intestines comparing gene expression among the SCIENTIFIC RepoRTs | (2018) 8:15508 | DOI: 10 .1038/s41598-018-33975-x uninfected, virus-infected, fungus-infected, and co-infected treatments ( Fig. 2(A) ). RNA sequencing resulted in ~700 million paired-end reads passing filters, 65% of which aligned to the M. lucifugus genome (Table S1 ). Pairwise differential gene expression varied widely among the four treatments with relatively low overlap in differentially expressed transcripts (Fig. 2B ,C), Supplementary Fig. 1 ). Similar transcript expression occurred between the uninfected and virus-infected bats, and between the fungus-infected and co-infected bats.

The fungus-infected bats exhibited a much stronger response, differentially expressing 324 transcripts compared to the uninfected bats (Table S3 ). These transcripts were enriched for only two gene ontology (GO) terms (cell-cell junction and plasma membrane part; Table S8 ). The co-infected bats differentially expressed 634 transcripts relative to the uninfected bats (Table S4) . These transcripts showed significant enrichment for 16 GO terms (Table S8 ). The co-infected and fungus-infected bats shared 108 similar differentially expressed transcripts and overlapped in one enriched GO term relative to the uninfected bats (plasma membrane part; Table S8 ).

When we directly compared responses of bats among the four treatments, response of the virus-infected bats differed strongly from the responses of fungus-infected or co-infected bats (virus-infected vs. fungus-infected: 461 differentially expressed transcripts and 9 significantly enriched GO terms; virus-infected vs. co-infected: 473 transcripts and 43 enriched GO terms; Tables S5, S6, S7; Supplementary Fig. 1 ). These differences in gene expression patterns included genes that clustered in two processes relevant to host-pathogen interactions -the mitogen-activated protein kinase (MAPK) pathways and cytokine and innate immune responses. Table 2 lists genes from the two processes that were significantly either up or down-regulated when virus-infected bats were compared to co-infected bats. For the MAPK pathway-related transcripts, genes such as RSU1 and RERG were up-regulated while those, such as STYK1, RRAD, MAP3K and SRC were down-regulated. For cytokine-related genes several transcripts were suppressed. When we compared the expression of the same genes for bats with WNS (combining fungus-infected and co-infected bats) and all bats without WNS (combining uninfected and virus-infected bats), we found similar differences (last two columns of Table 2 ). This suggested that superficial infection with fungus, P. destructans, was the driving factor for altered gene expression in the bat intestines. To confirm the results of the RNA-seq analysis, we selected 4 genes from Table 2 , namely IRF1, RERG, SRC and IL22RA1, to be verified by RT-qPCR. We also included interleukin 10 (IL10) due to its biological relevance to immune regulation and because we had previously observed an increase in its expression related to fungal infection 43 . As we wanted to confirm whether WNS was driving gene expression in the intestines of bats, we performed a two-group analysis for the RT-qPCR data. We combined all the bats without WNS into a single group (uninfected + virus-infected) and all the bats with WNS into the other group (fungus-infected + co-infected; Fig. 3 (A)). Expression of Ras-like estrogen regulated growth inhibitor (RERG) increased while expression of Interleukin 22 receptor subunit alpha 1 (IL22 RA1) genes decreased in bats with WNS, irrespective of viral infection (Fig. 3C,E) . Expression of the immune modulatory cytokine IL10 tended to be higher in bats with WNS than in bats without WNS, but the difference was not statistically significant (p-value = 0.07) (Fig. 3(F) ).

White-nose syndrome is associated with increased coronavirus antibody levels in the co-infected bats. In 2017, we performed a similar study, experimentally exposing 63 M. lucifugus to P. destructans as described in Warnecke et al. 37 . We performed IgG ELISA on blood plasma to detect Myl-CoV (coronavirus) N protein antibodies and found that 21/63 were positive for antibodies against the coronavirus. Of those 21 bats, 7 had detectable coronavirus RNA in their intestines suggesting an active infection, and 3 out of the 7 had been experimentally infected with P. destructans during the course of the study. We compared the ELISA optical density (O.D.) values of these virus-infected bats to co-infected bats ( Fig. 4(A) ) and found that the presence of P. destructans was associated with increased levels of coronavirus antibodies (Mann Whitney test, p value = 0.03; Fig. 4(B) ).

Our findings suggest that systemic responses of bats to WNS results in increased coronavirus replication and consequently, increased viral shedding, which may lead to subsequent infection of susceptible animals. Coronavirus infection may in turn increase the severity of WNS pathology. This is the first study to examine the systemic effects of co-infection on either bat coronavirus or WNS, and our results raise important questions in regard to zoonotic spillover events. Although events of successful viral spillover to distantly related species are thought to be extremely rare, in recent years several coronaviruses have spilled over, including SARS-CoV 1 , MERS-CoV 2-5 , PEDV-CoV 6 and SADS-CoV 7 . These viruses are thought to have originated in bats. In addition, circumstantial evidence suggests that most alpha and beta coronaviruses that parasitize other mammals may have originated in bats as well 46 . If so, then understanding host-pathogen interactions between bats and coronaviruses could inform our ability to predict or manage the risk of spillover. In this study, we showed that a coronavirus exhibits low activity in its natural host, M. lucifugus, but that co-infection with a fungus increases the quantity of viral RNA

Virus-infected vs. Coinfected Table 2 . RNA-sequencing identified differential expression of transcripts related to the MAPK pathway and to cytokine-related processes, comparing gene expression in the ileum of little brown bats (Myotis lucifugus) infected only with the M. lucifugus coronavirus (Myl-CoV; virus-infected) or co-infected with Myl-CoV and Pseudogymnoascus destructans. The last 2 columns show the same comparisons made after grouping bats that were not exposed to the fungus, and bats that were exposed to the fungus and exhibiting symptoms of WNS (irrespective of their viral infection status). a (Uninfected + virus-infected) vs. (fungus-infected + co-infected) b Positive log2 fold-change values indicate higher expression in the second listed treatments relative to the first. with a red arrow indicating the two groups ("with fungus" and "without fungus") that were compared. (B-F) The relative transcript levels of each gene for bats with and without WNS, depicted as reciprocal of Cycle threshold (Ct) normalized separately (ΔCt) for levels of transcripts for GAPDH in each sample. Statistical significance was calculated based on the independent Mann Whitney test. The difference in the two groups was significant for RERG and IL22RA1 genes.

in the intestines. We have no reason to expect zoonotic transmission of the coronavirus i.e. Myl-CoV, but similar co-infection mechanisms may operate in tropical bat species harbouring potentially zoonotic viruses.

Our results suggest that secondary skin infection with the fungus, P. destructans, substantially increases the level of viral RNA in the intestine of hibernating bats. We showed that infection of the skin with P. destructans can cause profound changes in gene expression in the intestines, despite a lack of direct contact between intestinal tissue and the fungus. Infection with P. destructans causes modulation of a number of immune responses, including down-regulation of interleukin and cell proliferation genes which may compromise bats' ability to suppress viral activity (Fig. 5) . Taken together, our results have implications for epidemiological studies of P. destructans, the WNS fungus and for research into viral spillovers, which should consider the potential implications of co-infections that increase viral shedding.

Complex strategies allow viruses to remain endemic in populations. These include a continuously replenished source of susceptible hosts for viruses that cause short-term acute infections with long-lasting immunity (e.g. measles virus), antigenic drift of virus (e.g. influenza virus) or waning immunity (e.g. respiratory syncytial virus) that allows reinfection, and long-lasting latent (e.g. herpesviruses) or persistent infections (e.g. pestiviruses) with sustained or periodic shedding. It is not yet clear how bat viruses are maintained in their natural host populations, or how they avoid extinction as host populations become immune and less susceptible. Persistent infections can be established in cultured cells with viruses that may have originated in bats, including Ebola virus 47 and SARS-CoV 48-50 , but whether these viruses persist in their primary hosts is not known. Studies of persistence of bat viruses in infected bats have produced equivocal results. The lack of direct evidence supporting specific models of persistence or transmission dynamics represents a major knowledge gap in bat-virus ecology 21 .

We maintained M. lucifugus in controlled laboratory hibernation chambers for four months during these experiments, and we detected the coronavirus i.e. Myl-CoV, at the end of hibernation. These data imply that the coronavirus can persist in its host for at least the duration of hibernation, particularly as nucleotide variability among the detected coronavirus isolates showed that spread of coronavirus among bats within a chamber was unlikely 32 . In an extensive study of New World alphacoronaviruses, no target viruses were detected in the rectal swabs of individual bats sampled over time 46 , suggesting that persistence and intensity of shedding varies among species or viruses. The authors concluded that the targeted coronaviruses do not persist in their hosts but are maintained in populations by the introduction of new susceptible individuals. However, their results could also reflect viral persistence in individual animals, with low baseline levels of virus replication and undetectable shedding interspersed with periods of increased replication and shedding that did not occur during the sampling period.

Periodic or seasonal increases in virus shedding associated with parturition, lactation, nutritional deprivation or environmental stress 21, 29 suggest persistent or latent viruses may be activated by hormonal or other systemic cues. Direct evidence linking a specific trigger to increased shedding has not yet been found. However, viral replication in rodent and bat cells persistently infected with Ebola virus increased greatly following modulation of the Ras/MAPK pathway with lipopolysaccharides or phorbol esters, and with the resulting suppression of the cells' interferon response [47] [48] [49] . In experimental systems, the immune modulatory cytokine IL10 also influences viral persistence and replication [51] [52] [53] , although more study is required to clarify the effects of circulating cytokines on the replication of persistently infecting viruses. Nevertheless, these results suggest that circumstances which induce anti-inflammatory cytokines or suppress anti-viral innate responses, may provide a trigger for increased shedding of persistently infecting virus.

We discovered that bats with WNS (fungus-infected and co-infected) had significantly lower intestinal levels of transcripts for IL22RA1 and other interferon-related genes as compared to uninfected bats, and we observed the same trend in IL10 (although it was not significant; p value = 0.07). IL22RA1 is the receptor present on host cells, including intestinal cells, which help in initiating cellular signalling in response to IL22 produced by T-cells 54 . IL-22 leads to an increase in anti-microbial peptide production, cellular protection against damage and increases cellular proliferation 55 . Therefore, reduced IL-22 signalling in the intestines of bats with WNS, might suppress the bat defences that control the coronavirus infection. Additionally, previous studies have shown that the anti-inflammatory gene, IL10, is expressed more in the lungs of bats with WNS than in bats without it 43 . We saw a similar trend with the levels of IL10 in the intestines which might play a role in suppressing the immune response against the coronavirus. Another altered cytokine gene which was of interest was the suppressor of cytokine signalling-6 (SOCS6) gene. Fungal-infected bats showed lower levels of SOCS6 transcripts, lack of which has been implicated in mild growth retardation in mice 56 . Overall, our results suggest that WNS triggers changes in gene expression in the ileum (Fig. 5) . These may influence expression of interferon-stimulated-genes (ISGs), thereby leading to increased viral replication at the site of viral persistence. Interferon-related transcripts were more highly expressed in the ileum of virus-infected bats that did not have WNS, suggesting that the bat's response to WNS causes down-regulation of interferon activity. Interferons may control coronavirus replication, as seen in cases of SARS-CoV 57 and MERS-CoV 58 . Therefore, a decrease in interferon activity might cause an increase in coronavirus (Myl-CoV) replication. In addition to interferon-related genes, we also found that RERG, which is related to growth inhibition, was upregulated in the fungus-infected bats when compared to virus-infected bats. Upregulation of RERG could affect the rate of cell proliferation in the intestines 59 . Finally, this cascade of responses is associated with increased severity of WNS symptoms. Bats with WNS experience a range of systemic disturbances including dehydration, hypovolemia, metabolic acidosis and fat depletion 40, 41 , neutrophil infiltration of the lung interstitium, and increased expression of transcripts related to anti-microbial and pro-and anti-inflammatory cytokines 43 . Taken together, this evidence suggests that hibernating bats respond systemically to superficial fungal infection, and this hypothesis is further supported by our observations of altered gene expression in the ileum of fungus-infected bats.

Based on our results, we propose a model for how secondary infections may increase the replication and subsequent shedding of persistently infecting virus (Fig. 5) . The establishment of WNS (or other secondary infection) impacts the tissue with which that pathogen interacts (in the case of P. destructans, the skin). Direct interactions between the host and the secondary pathogen are limited to the affected tissue, but the systemic response to the disease triggers a cascade of immune responses, including increased release of cytokines or neutrophils. Affected cells such as intra-alveolar macrophages in the lungs or cells lining the intestine, may produce pro-or anti-inflammatory molecules and influence cells that harbour viral genomes. This cascade of host responses disrupts the equilibrium between the persistently infecting virus and the cell's innate immune response, leading to a dramatic increase in the expression of coronavirus (Myl-CoV) replication. Our assays were unfortunately limited to analysing viral and cytokine transcripts rather than protein, because reagents for detecting bat viral and host proteins are not yet available. We were not able to perform serial dilutions of the plasma to precisely quantify anti-viral titre due to the limitation in the amount of plasma obtained from each bat. The sample size was also small for this assay, because only 7 of the sampled bats had detectable levels of coronavirus in their intestines and were positive for viral antibodies. Despite these limitations, we demonstrated higher antibodies against the coronavirus in the plasma of co-infected bats when compared with virus-infected bats. This increased antibody level in co-infected bats might reflect an adaptive immune response to increased coronavirus replication in the intestines.

Our proposed hypothesis for the mechanism driving increased viral replication following pathogenic co-infection was worth testing, but our results are also consistent with an alternative hypothesis. Increased viral replication or viral load may affect the severity and population-level impacts of WNS. Bat mortality following the arrival of WNS varies widely from site to site, with populations decreasing from 30% to 99% 35 . Variation in the microclimate, and other ecological factors may drive some of this variation 60 , but our data suggest that cryptic viral infections may also play a role in determining survival rates for bats hibernating in sites colonized by P. destructans. We recommend that future studies on population-wide impact of WNS incorporate viral sampling to help better understand the role of co-infections on bat populations in the wild.

Sample acquisition. Fifty-four male M. lucifugus were collected from a WNS-free cave in Manitoba, Canada in November 2010. Details of the experimental design as well as protocols for collecting and transporting bats, infection with P. destructans, maintenance of bats in hibernation and sample collection have been described previously 37, 43 . Briefly, bats in groups of 18 were either sham-inoculated or inoculated with North American or European isolates of P. destructans. Bats were housed at 7 °C and >97% relative humidity with ad libitum water. All bats were equipped with data loggers to monitor skin temperatures. Bats were euthanized during the experiment when humanely required or at the termination of the experiment 120 days after inoculation. Immediately following euthanasia samples from segments of wing as well as various tissues were preserved in RNAlater (Qiagen, 76016) or in formalin. Samples in RNAlater were kept at −20 °C until they were processed. North American and European isolates of P. destructans caused similar disease outcomes 37 , so we did not differentiate between the strains in subsequent analysis. The procedures for care, handling and euthanasia of bats were approved by the University Committee on Animal Care and Supply of the University of Saskatchewan (Protocol #20100120). Bats were collected under the province of Manitoba Wildlife Scientific Permit WB11145.

In 2017, a further 129 M. lucifugus were collected from a WNS-free cave in Manitoba, Canada in January under the Manitoba Sustainable Development Wildlife Scientific Permit No. SAR16009. Bats were euthanized during the experiment when humanely required or at the termination of the experiment 70 days after infection and a similar experiment was performed at the University of Winnipeg as described above (Protocol #AE08399).

Histological classification. During necropsy, we collected representative samples for histopathology from all major organ systems. In addition, representative samples were taken from all areas of the wing and rolled on dental wax before placing in 10% neutral buffered formalin. Tissues were processed routinely for histology. Five µm sections were cut and stained with periodic acid-Schiff stain to highlight fungal hyphae. Liver and other tissues were processed routinely and stained with hematoxylin and eosin. Wings were scored on a scale of 0 to 5 with 5 being very severe with >50% of wing covered in fungal hyphae. We used a bacterial score from 0 to 5, with 5 indicating wide-spread and abundant bacteria being present in many areas within the dermis and underlying connective tissues. Average scores from 5 sections of wing were used for analysis. Interstitial lung neutrophil assessment was similarly evaluated on a scale of 0 to 5, with 5 being very severe. Average scores from the 5 sections were used for analysis.

Tissues preserved in RNAlater were homogenized in 2 ml sealed vials with a 5 mm stainless steel bead, 0.1 g of 0.1 mm zirconia/silica beads and 350 μl Buffer RLT Plus (with β-mercaptoethanol, RNeasy Plus Mini Kit, Qiagen, 74136) using a Retsch MM400 Oscillating Mill at 30 Hz for 4 min. Total RNA was extracted following the manufacturers protocol. RNA integrity was assessed using RNA 6000 Nano Kit (Agilent, 5067-1511) with the Agilent 2100 Bioanalyzer. cDNA Synthesis. cDNA was synthesized from 1 μg of RNA (or less if concentrations were too low) per reaction using QuantiTect Reverse Transcription Kit (Qiagen 205313). cDNA samples were stored at −80 °C until they were used for PCR.

Polymerase Chain Reaction (PCR). Tissue samples were identified by their submission numbers with no reference to treatment class prior to analysis with PCR, so that evaluation of the results could not be inadvertently biased by knowledge of the treatment. We used semi-nested PCR to detect Myl-CoV. Primers were designed from the partial sequence of Rocky Mountain bat coronavirus replicase (accession number EF544563) ( Table S9) purified using MinElute PCR Purification Kit (Qiagen, 28006) and verified by sequencing (Macrogen, Korea). If more than one DNA band was present, the appropriate size band was cut out and purified using QIAquick Gel Extraction Kit (Qiagen, 28706) before sequencing.

in conjunction with QuantiFast SYBR Green PCR Kit (Qiagen 204056). We quantified coronavirus with RNA primers MyCVF2 and MyCV R1 (Table S9) . For initial experiments data were normalized to two transcriptsglyceraldehyde-3-phosphate dehydrogenase (GAPDH) and beta-actin 43 . As there were no differences in results, all subsequent experiments used only GAPDH as a normalizer using primers GAPDH US and GAPDH DS (designed for use in humans but also amplify M. lucifugus transcripts - Table S9 ). As well, a no-template (negative) control was included with every set of primers. Each 25 μl reaction contained: 1 μM of each primer set, 12.5 μl SYBR Green Master Mix and 8.5 μl of diluted cDNA.

To verify the RNAseq data, cDNA from ileum samples in which coronavirus RNA had been detected via RT-qPCR were analysed using the following primers, IL22RA1, IRF1, RERG and SRC (for sequence of primers see Table S9 ). Primers were designed by aligning primers described for quantitating human cytokines (PrimerBank) with annotated transcripts of M. lucifugus genes: c-jun (Accession number: XM_006096110.1), cyclin D1 (XM_006098046.1), IL10 (XM_006094865.1) and TNF alpha (XM_006104644.1). The interferon beta primers were designed using the annotated transcript for the E. fuscus gene (XM_008145044.1), which also amplify transcripts from M. lucifugus. Primer efficiencies were determined from cycle threshold (Ct) values of purified PCR products serially diluted and re-amplified. Primers amplified targets with an efficiency of about 100% and in all cases the identities of the PCR products were confirmed by their specific dissociation temperature, specific sizes on agarose gels and by sequencing.

We observed primer-dimers in some reactions in addition to the PCR product. The dimers dissociated at 77 °C, while the specific coronavirus polymerase product dissociated at 83 °C. To avoid false positives due to primer-dimers, the thermocycler was programmed to read at 80 °C (in the cycle after the primer-dimer had dissociated, and before dissociation of the target product). The thermal profile used was 95 °C for 5 min followed by 40 cycles of 95 °C for 10 sec, 60 °C for 30 sec (readings taken at 80 °C), and a final cycle of dissociation of product 95 °C for 1 min, 55 °C for 30 sec and 95 °C for 30 sec (readings taken at every degree between 55 °C and 95 °C).

Only results from reactions that yielded unambiguous results were used for analysis.

To explore the mechanisms driving high virus load in bats with WNS, we performed RNA-seq analysis which could potentially screen all targets in the bat intestinal cells. We targeted the ileum transcriptome because this is the tissue in which Myl-CoV is present in detectable concentrations 32 . Extraction of RNA from ileum tissue, which includes the ileum and potential gut contents have been described in previous sections. Bats were screened for Myl-CoV using RT-qPCR, and bats were assigned post hoc to treatment groups representing four infection histories ( Fig. 1(A) ): 1) Uninfected (bats were not infected with virus or fungus; n = 5), 2) Virus-infected (bats were naïve to the fungus but had a persistent Myl-CoV infection; n = 4), 3) Fungus-infected (bats were experimentally infected with P. destructans and no virus was detected; n = 3), or 4) Co-infected (bats with persistent Myl-CoV infections that were also experimentally infected with P. destructans; n = 4). All samples had adequate RNA quality for sequencing (i.e. RIN value >7).

Tissues were homogenized in 2 ml sealed vials with a 5 mm steel bead, 0.1 g of 0.1 mm zirconium silica beads, 350 µL of RLT buffer (with β-mercaptoethanol) (RNeasy Plus Kit, Qiagen) using a Retsch MM400 tissue homogenizer at 30 Hz twice for 2 minutes each. Total RNA from tissues was extracted using the procedure provided with the RNeasy Plus Kit. cDNA library preparation and RNA-sequencing. Total RNA was sent to The Centre for Applied Genomics at The Hospital for Sick Children (Toronto, Canada). RNA quality was assessed using a Bioanalyzer (Agilent Technologies). We retained all samples with a DV200 (percentage of RNA fragments greater than 200 nt) greater than 85% (Table S1 ), discarding one Co-infected sample with a DV200 = 42%. Poly(A) mRNA was enriched using oligo dT-beads, and cDNA libraries were prepared using the NEBNext Ultra Directional RNA Library Prep Kit for Illumina (New England BioLabs). Barcoded libraries were pooled in equimolar quantities, and the sixteen libraries were sequenced on three lanes of a HiSeq. 2500 System (Illumina Inc.), which generated 126 bp paired-end reads.

RNA-sequencing read alignment and analysis. We used FastQC v0.11.5 61 to assess sequence quality and Trimmomatic v0.36 62 to remove the adapter sequences and low-quality bases from reads with the following settings: Illumina clop:2:30:10, leading:3, tailing:3, slidingwindow:4:15, minlength:36. We used TopHat v2.1.1 63 to align the trimmed paired-end reads from each library, separately, to the Ensembl M. lucifugus genome sequence (Myoluc2.0 64 ) in strand-specific mode (fr-firststrand) with mate-inner-dist values specific for the insert size of each library. We used featureCounts 65 to count reads mapped to the Myoluc2.0 genome annotation in strand-specific mode (reversely stranded), counting paired-end reads as fragments, counting only those fragments where both reads aligned successfully, counting multi-mapping fragments, and excluding chimeric fragments. We assessed the variability within and between the treatments using the R package SARTools v.1.3.0 66 . The featureCount-estimated gene counts were transformed by a variance stabilizing method (VST) using SARTools.

We identified differentially expressed genes between each of the treatments using DESeq2 v.1.12.3, run in SARTools. Custom SARTools-based DESeq settings included: cooksCutoff = TRUE (perform outliers detection), independentFiltering = TRUE, alpha = 0.05 (threshold of statistical significance), pAdjustMethod = BH (benjamini hochberg p-value adjustment method), and locfunc = median (estimate size factors). Differentially expressed genes were identified as having a fold-change >2 and false discovery rate (FDR)-corrected p-values < 0.05 67 . We SCIENTIFIC RepoRTs | (2018) 8:15508 | DOI:10.1038/s41598-018-33975-x produced volcano plots representing the differential expression comparisons by plotting the log of the adjusted p value as a function of the log ratio of differential expression. We used the Ensembl gene IDs identified by DESeq2 as input for the web-based g:Profiler 68 to test for gene ontology (GO) term enrichment among the differentially expressed genes, using a FDR significance threshold <0.05. These GO-terms and their corresponding p-values were used in REViGO 69 to visualize significant enrichment of biological processes.

Myl-CoV N protein expressed in infected E. coli BL21 cells was used as positive antigen, and GST-tagged protein expressed in uninfected BL21 cells was used as negative antigen. 96-well Costar high-binding round-bottom assay plates were coated with 0.05 µg/well of either antigen diluted in 0.1 M phosphate buffered saline (pH 7.4) in a total volume of 100 µl. Plates were covered and incubated overnight at 4 °C and washed three times with 300 µl of PBS-Tween 20 (0.1%) immediately prior to use. Serum samples were diluted to 1:100 in PBS-Tween 20 (0.2%) supplemented with 5% fetal bovine serum (Gibco, Thermofisher). 100 µl of each sample was added in parallel to a positive and negative antigen plate and incubated at 37 °C for one hour and washed as above. A peroxidase-labelled goat anti-bat IgG secondary antibody (0.05 µg in 100 µl per well, Bethyl labs) was added, incubated for one hour at 37 °C and washed as above. Peroxidase substrate (2,2'-azino-bis (3-ethylbenzthiazoline-6sulfonic acid)) was added to each well and colour development was quantified 30 minutes later by measuring the optical density at 405 nm using an ELISA microplate reader. The ELISA cut-off value (0.39) was calculated as the Statistical analysis. Data from RT-qPCR and histopathological scores were analysed with SPSS Statistics version 23. The relative levels of a transcript for each bat were calculated as RT-qPCR Cycle threshold (Ct) normalized separately (ΔCt) to the "house-keeping" gene GAPDH. A ΔCt reduction of one (1) indicates an approximately two-fold higher concentration of RNA. The significance of differences of mean values of ΔCt between co-infected bats and virus-infected bats were determined using an independent-samples Mann-Whitney U test. We calculated Pearson's coefficients to test the correlation between ΔCt levels for coronavirus polymerase cDNA for bats in each treatment class, and average scores for fungal hyphae, secondary bacteria, oedema, necrosis and inflammation in wing lesions, as well as bacteremia and levels of neutrophils in lung, spleen and liver interstitium.

Bat studies were carried out in strict compliance with Canadian Council on Animal Care guidelines and the procedure for care, handling, and euthanasia of bats were approved by the University Committee on Animal Care and Supply of the University of Saskatchewan (protocol s#20100120).

All RNA-seq fastq files have been submitted to the NCBI Sequence Read Archive database (accession number SRX3752319-SRX3752333). These data files will be released to the public upon acceptance of this manuscript for publication.

Isolation and characterization of a bat SARS-like coronavirus that uses the ACE2 receptor

Rooting the phylogenetic tree of middle East respiratory syndrome coronavirus by characterization of a conspecific virus from an African bat

Close relative of human Middle East respiratory syndrome coronavirus in bat

MERS-related betacoronavirus in Vespertilio superans bats

Further Evidence for Bats as the Evolutionary Source of Middle East Respiratory Syndrome Coronavirus

Origin, evolution, and genotyping of emergent porcine epidemic diarrhea virus strains in the United States

Fatal swine acute diarrhoea syndrome caused by an HKU2-related coronavirus of bat origin

Isolation of Hendra virus from pteropid bats: a natural reservoir of Hendra virus

Nipah virus infection in bats (order Chiroptera) in peninsular Malaysia

Isolation of genetically diverse Marburg viruses from Egyptian fruit bats

Fruit bats as reservoirs of Ebola virus

Evolutionary insights into the ecology of coronaviruses

Bats host major mammalian paramyxoviruses

Bat rabies-a Gordian knot

Ebolavirus and other filoviruses

Coronaviruses in bent-winged bats (Miniopterus spp

Pteropid bats are confirmed as the reservoir hosts of henipaviruses: a comprehensive experimental study of virus transmission

Experimental Nipah virus infection in pteropid bats (Pteropus poliocephalus)

Antiviral immune responses of bats: a review

Bat flight and zoonotic viruses

Ecological dynamics of emerging bat virus spillover

Urban habituation, ecological connectivity and epidemic dampening: the emergence of Hendra virus from flying foxes (Pteropus spp.)

Social network characteristics and predicted pathogen transmission in summer colonies of female big brown bats (Eptesicus fuscus)

Amplification of emerging viruses in a bat colony

Hendra virus infection dynamics in Australian fruit bats

Reproduction and nutritional stress are risk factors for Hendra virus infection in little red flying foxes (Pteropus scapulatus)

Spatial and temporal patterns of Zaire ebolavirus antibody prevalence in the possible reservoir bat species

Characterization of Nipah virus from naturally infected Pteropus vampyrus bats

Physiological stress and Hendra virus in flying-foxes (Pteropus spp

Seasonal pulses of Marburg virus circulation in juvenile Rousettus aegyptiacus bats coincide with periods of increased risk of human infection

A persistently infecting coronavirus in hibernating Myotis lucifugus, the North American little brown bat

1 coronaviruses in bats in North America

Bat white-nose syndrome: an emerging fungal pathogen?

An emerging disease causes regional population collapse of a common North American bat species

Experimental infection of bats with Geomyces destructans causes white-nose syndrome

Inoculation of bats with European Geomyces destructans supports the novel pathogen hypothesis for the origin of white-nose syndrome

The White-Nose Syndrome Transcriptome: Activation of Anti-fungal Host Responses in Wing Tissue of Hibernating Little Brown Myotis

The other white-nose syndrome transcriptome: Tolerant and susceptible hosts respond differently to the pathogen

White-nose syndrome initiates a cascade of physiologic disturbances in the hibernating bat host

Pathophysiology of white-nose syndrome in bats: a mechanistic model linking wing damage to mortality

White-nose syndrome increases torpid metabolic rate and evaporative water loss in hibernating bats

Activation of innate immune-response genes in little brown bats (Myotis lucifugus) infected with the fungus Pseudogymnoascus destructans

Within-Host Competition between Two Entomopathogenic Fungi and a Granulovirus In

Infection and co-infection by the amphibian chytrid fungus and ranavirus in wild Costa Rican frogs

Alphacoronaviruses in New World bats: prevalence, persistence, phylogeny, and potential for interaction with humans

Stimulation of Ebola virus production from persistent infection through activation of the Ras/MAPK pathway

Mechanisms of establishment of persistent SARS-CoV-infected cells

Characterization of persistent SARS-CoV infection in Vero E6 cells

Severe acute respiratory syndrome coronavirus persistence in Vero cells

IL-10: achieving balance during persistent viral infection

Shifting hierarchies of interleukin-10-producing T cell populations in the central nervous system during acute and persistent viral encephalomyelitis

The role of IL-10 in regulating immunity to persistent viral infections

Interferon-λ and interleukin 22 act synergistically for the induction of interferon-stimulated genes and control of rotavirus infection

Interleukin-22: a novel T-and NK-cell derived cytokine that regulates the biology of tissue cells

SOCS-6 binds to insulin receptor substrate 4, and mice lacking the SOCS-6 gene exhibit mild growth retardation

Interferon-beta 1a and SARS coronavirus replication

Inhibition of novel beta coronavirus replication by a combination of interferon-alpha2b and ribavirin

RERG is a novel ras-related, estrogen-regulated and growth-inhibitory gene in breast cancer

Sociality, density-dependence and microclimates determine the persistence of populations suffering from a novel fungal disease, white-nose syndrome

Quality Control tool for High Throughput Sequence Data

Trimmomatic: a flexible trimmer for Illumina sequence data

TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions

featureCounts: an efficient general purpose program for assigning sequence reads to genomic features

SARTools: A DESeq2-and EdgeR-Based R Pipeline for Comprehensive Differential Analysis of RNA-Seq Data

Controlling the False Discovery Rate: A Practical and Powerful Approach to Multiple Testing

Profiler-a web server for functional interpretation of gene lists (2016 update)

REVIGO summarizes and visualizes long lists of gene ontology terms

Funding for this study was provided by a United States Fish and Wildlife Service grant to Craig K.R. Willis, Trent K. Bollinger and Vikram Misra, Natural Sciences and Engineering Research Council (Discovery) grants to Christopher J. Kyle, Vikram Misra and Craig K.R. Willis, a University of Saskatchewan One Health Initiative grant to Vikram Misra and a fellowship within the Postdoc Programme of the Deutscher Akademischer Austausch Dienst (DAAD), German Academic Exchange Service to Lisa Warnecke. The funders had no role in study design, data collection and analysis, decision to publish or preparation of the manuscript.

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.