Virus genomes interact with hundreds of host RNAs

We have previously performed proximity ligation using biotinylated psoralen (SPLASH) in human neuronal progenitor cells cell (hNPCs) and Huh7 cells to study intramolecular RNA-RNA interactions within DENV and ZIKV genomes [12]. Besides capturing interactions within an RNA molecule, SPLASH also crosslinks intermolecular virus-host RNA interactions to provide information on the identity of host RNAs and the locations on which they bind to DENV and ZIKV RNA genomes. We performed analysis on four different strains of ZIKV in hNPCs and DENV-1 in Huh7 cells (Fig. 1a). The virus genomes show highly reproducible patterns of binding to host RNAs across biological replicates (Additional file 1: Fig. S1, R > 0.9). After filtering for highly consistent interactions across biological replicates, we observed that hundreds of human RNAs, including coding RNAs and long and short non-coding RNAs, interact with DENV and ZIKV genomes (Fig. 1b, Additional file 2: Table S1). These RNAs include Malat1, SND1, SSR2, SSR3, and SEC61A1, several of which were previously found in other siRNA and CRISPR host factor screens, confirming their importance in virus pathogenicity. Most of the viral interacting RNAs are localized in the cytosol (Additional file 1: Fig. S2a,b, Additional file 3: Table S2). GO term analysis of these factors showed that DENV and ZIKV binding RNAs are enriched for virus transcription and SRP signaling, agreeing with the importance of translation at the ER membrane for DENV and ZIKV replication [6] (Fig. 1c).

Fig. 1

DENV and ZIKV genomes interact with hundreds of host RNAs in cells. a Experimental workflow to identify virus-host RNA interactions globally using biotinylated psoralen and proximity ligation sequencing. b Pie chart showing the number of RNAs from different RNA classes that interact with DENV and ZIKV viruses. c GO term enrichments of the host RNAs that interact with different strains of DENV and ZIKV viruses. The Y-axis indicates log (p-value) and the X-axis indicates the number of folds that the genes are enriched in binding to the virus genome as compared to the background

To better understand the nature of these virus-host interactions, we calculated the locations and frequency of virus-host interactions along with host and virus genomes. We observed that the entire genome of ZIKV interacts with host RNAs inside cells, although different classes of RNAs such as snoRNAs or mRNAs can bind to different regions along the virus genome (Fig. 2a). Across the different virus domains, we observed that the NS2B-coding region tends to be particularly enriched for binding to host RNAs (Fig. 2b). We observed that some of the virus-host RNA interactions are extremely specific, whereby only one region along the entire virus genome binds to the host RNA (Fig. 2c), or only one region along the host RNA binds to the virus genome (Fig. 2d). We also observed other RNAs that bind to the virus genome promiscuously by interacting with multiple regions along the viral genome (Fig. 2c), or by using multiple regions along the host RNA to interact with the genome (Fig. 2d). Along the host RNAs, we did not observe a preference for using either the 5′UTR, CDS, or 3′UTR to interact with the virus (Additional file 1: Fig. S2c), although bases that interact with ZIKV genomes are more evolutionarily conserved than non-interacting bases (Additional file 1: Fig. S2d). Additionally, abundant RNAs are also not significantly associated with high total interaction counts along the virus or high interaction counts at a specific site (Additional file 1: Fig. S2e), indicating that our observed virus-host interactions are not just a result of random interactions of the virus genome with abundant RNAs.

Fig. 2

General properties of Zika-host RNA-RNA interactions. a Locations of interacting host RNAs on the viral genomes. Aggregate data and data classified by type of host RNA are shown. We observed differences in binding patterns by different classes of RNA along the genome. b The number of interactors in each coding region normalized by the length of the respective region. We observed the highest abundance of interactions in the NS2B-coding region. c Number of Zika genome binding sites for any specific host RNA (“Zika to host”). Most interactions bind at a specific, unique site. Promiscuous interactions of host RNA are present but considerably rarer. Inset, top: An example of a highly specific interaction between host RNA and virus genome, whereby SNORD27 only binds to a single region in ZIKV. Bottom: An example of promiscuous interaction between host RNA and virus genome whereby 7SK binds to many locations along ZIKV. d Number of interaction sites on host RNAs that bind to ZIKV genome (“host to Zika”). Data shows mostly unique interactions, but a higher propensity for 2 or more interaction sites. Inset, top: An example showing only 1 host region in ZNF485 interacting with the virus genome. Bottom: An example showing that many regions along SLC25A6 bind to the virus genome. e Sites along the Zika genome that show high numbers of interactions (99th percentile of interactions) show a significant preference for single-stranded, high SHAPE reactivity segments inside virion particles. f Host RNAs interacting at ZIKV position 1217 with low predicted interaction energy corresponding to low observed read counts and high interaction energies coinciding with high observed read counts, indicating interaction energetics is a significant factor driving specific host-virus interactions. g Aggregate interaction energy statistics for percentiles of observed interaction counts, again with high interaction counts corresponding to high predicted interaction energy and significant differences between all classes. h Volcano plot showing the fold change in gene expression of host RNAs after 24 h of ZIKV infection. ZIKV interactors are colored as red dots, while the non-interactors are in blue. The dotted lines indicate a 1.5-fold change in gene expression

We observed that host RNAs frequently interact with initially single-stranded regions along the virus (Figs. 3c and 5b). As these regions could serve as intermolecular interaction hubs due to their relative structural accessibility, we tested whether this observation holds globally by determining the single-strand propensity of high versus low intermolecular interaction sites along ZIKV. We used SHAPE-MaP reactivities of ZIKV genome inside virion particles as it represents the structural state of the genome before being exposed to host RNAs inside cells [12, 13]. SHAPE-MaP reactivities show that regions that participate in high intermolecular interactions inside cells tend to be more single-stranded in virion particles (Fig. 2e), indicating that highly accessible regions along the genome could act as nucleating sites for host RNA interactions. We also observed that highly interacting ZIKV regions could bind to multiple host RNAs (Fig. 2f), suggesting that host RNAs could compete to interact with ZIKV and DENV genomes to result in different gene regulation outcomes. Interestingly, we observed that the strength of RNA binding to the virus genome (indicated by the number of SPLASH interactions) is correlated with predicted energetics of binding (Fig. 2g), and length of the interactions (Additional file 1: Fig. S2f), but not the GC content (Additional file 1: Fig. S2g). These analyses indicate that binding energetics is key to a host RNA’s ability to outcompete other RNAs when multiple RNAs can bind to the same site on the virus genome.

Fig. 3

DENV and ZIKV interact with host miRNAs. a, b Line plot showing the locations along ZIKV (a) and DENV (b) that interact with host miRNAs. Y-axis indicates the number of interaction counts between the virus genome and a specific miRNA. The X-axis indicates the position along the virus genome. c Top, secondary structure modeling of ZIKV genome before and after miR-19a-3p binding. Bottom: Predicted ZIKV: miR-19a-3p interactions using the program RNAcofold. The ZIKV interacting sequence is in red and the miRNA sequence is in blue. d, e Bar charts showing the amount of Zika virus detected inside Huh7 cells using qPCR, 2, 16, 24, and 48 h post-infection in cells that are transfected with miR-19a-3p or control (d), or in cells that are transfected with miR-19a-3p inhibitor or control (e). The data in cells with over-expression of miR-19a-3p or inhibitor is normalized to its control

Previously, in our study of identifying SARS-CoV-2 interactome, we had observed that SARS-CoV-2-interacting RNAs are preferentially stabilized during infection even though many transcripts are downregulated [14]. To examine how the host RNA changes in gene expression upon ZIKV infection, we performed RNA sequencing experiments at 0 and 24 h post-ZIKV or upon mock infection in Huh7 cells. We observed that a similar number of host RNAs are up- (908) or down (911)-regulated by at least 1.5-fold in gene expression upon ZIKV infection (Fig. 2h). In contrast to SARS-CoV-2 interacting RNAs, ZIKV-interacting RNAs show a slight downregulation in gene expression upon ZIKV infection (Additional file 1: Fig. S2h). Out of 15 ZIKV-binding transcripts that show gene expression changes, 14 of them are downregulated, while only 1 (Y-RNA) RNA is upregulated (Additional file 4: Table S3, Fig. 2h). However, the functional mechanisms of how host RNAs are downregulated upon virus infection and how their downregulation impacts virus fitness remain to be studied.

Virus genomes can interact with small non-coding RNAs to promote viral replication

In addition to mRNAs, we also observed that many short non-coding RNAs, including snoRNAs and miRNAs, interact with the viral genomes (Fig. 1b). We detected the binding of 5 miRNAs each to DENV and ZIKV genomes respectively (Fig. 3a, b). To determine whether any of our detected miRNAs could regulate virus pathogenicity, we focused on miR-19a-3p which was previously identified in an independent virus-host interaction study. Still, its impact on virus replication was unknown [15]. miR19a-3p is an important regulator in cancer, whereby it is known to target the key tumor suppressor PTEN [16]. Using the program RNA22 [17], we confirmed that miR19a-3p is predicted to hybridize strongly with the ZIKV genome at the base 9691, the same location where we detected miR-19a-3p:ZIKV interaction experimentally using SPLASH (Fig. 3c). To determine the impact of miR-19a-3p:ZIKV interaction on viral replication, we performed miRNA over-expression experiments by introducing a mimic of the miR19a-3p into Huh7 cells prior to ZIKV infection. We showed that the over-expression of miR-19a-3p mimic increased the relative amount of Zika RNA inside cells (Fig. 3d), suggesting that the miRNA acts as a virus-permissive factor. Additionally, the downregulation of miR-19a-3p in Huh7 cells by introducing an inhibitor before ZIKV infection resulted in a decrease in ZIKV genome abundance by about 50% over 48 h (Fig. 3e), confirming that miR-19a-3p promotes ZIKV replication.

To further investigate how miR-19a-3p can promote ZIKV growth, we studied the underlying RNA structures and alternative structures along the ZIKV genome in the absence of miR19a-3p binding. We observed two regions along the Zika genome can form alternative structures with upstream and downstream elements (Additional file 1: Fig. S3a), including a region (9860–9890 bases) that was previously identified to be highly structured and structurally conserved across different ZIKV (Additional file 1: Fig. S3a) [12]. As such, we hypothesize that miR-19a-3p binding to the ZIKV genome could potentially release this region to allow it to fold properly. Additionally, we also observed that miR-19a-3p levels do not significantly change upon ZIKV infection (Additional file 1: Fig. S3b). We hypothesized that the ZIKV genome could also serve as a sponge, whereby the binding of miR-19a-3p to ZIKV genome could result in having less miRNA being available to interact with its endogenous targets. Indeed, we observed a modest increase in PTEN mRNA levels, which is a target of miR-19a-3p, by an average of 20% in the presence of ZIKV (Additional file 1: Fig. S3c). Additionally, we also observed that the gene expression of other miR-19a-3p predicted targets is increased upon ZIKV infection (Additional file 1: Fig. S3d). These experiments indicate that miR-19a-3p could act as a pro-viral factor partially by rearranging ZIKV genome structure and by modulating the regulation of its downstream targets.

Host RNAs that bind consistently across viruses are enriched for non-coding RNAs

Highly abundant RNAs could interact with virus RNAs in a non-specific manner. As such, we identified RNAs that interact with viral RNAs with an interaction count that was higher than would be expected from their abundance (Additional file 1: Fig. S2e, Additional file 5: Table S4). Additionally, while specific host RNAs can bind to the different DENV and ZIKV viruses, we were curious to know whether there are host RNAs that can bind to at least 4 out of 5 viruses, and hence serve as general DENV/ZIKV interactors upon virus infection. We overlapped highly interacting RNAs that bind to the different viruses and identified 32 host RNAs that bind to all viruses, as well as 66 RNAs that bind to at least 4 out of 5 viruses (Fig. 4a, Additional file 1: Fig. S4a, Additional file 5: Table S4). Eleven out of 32 of the common interactors are noncoding RNAs, including snoRNAs. Interestingly, the location of the top binding peaks of these RNAs is highly conserved across the viruses (Fig. 4b), suggesting that their binding on the viral genomes is specific. To confirm that the virus-host interactions that we captured inside cells are indeed present, we performed pulldown experiments using biotinylated oligos against one of the non-coding RNAs that bind to multiple viruses (7SK). Pulldown and qPCR experiments showed that these host RNAs indeed interact with ZIKV as expected (Fig. 4c).

Fig. 4

Noncoding RNAs interact with DENV and ZIKV genomes and impact virus fitness. a Heatmap showing the host RNAs that interact with DENV-1 and 4 strains of Zika. Host RNAs that interact with more than 3 viruses are listed on the right of the heatmap. b Locations of top 50 virus-host interaction sites along DENV-1 and the 4 ZIKV strains. Many of the top virus-host interaction sites are conserved across the viruses (shown in red). c qPCR analysis of ZIKV, U2, U31 and 28S rRNA pulled down by 7SK and tetrahymena RNA. Tetrahymena RNA is used as a negative control. **, *** indicate p-values ≤ 0.01 and 0.001 respectively, using Student’s T-test. d Locations along ZIKV that bind to 7SK RNA. The Y-axis indicates the number of SPLASH interactions between ZIKV and 7SK at that position. e Top, schematic of the strongest ZIKV-7SK interactions along the ZIKV genome. Bottom, predicted pairing interactions between ZIKV and 7SK using the program RNAcofold. 7SK sequences are shown in blue and the two ZIKV interaction sequences are in red and green respectively. f, g Bar charts showing the amount of ZIKV inside Huh7 cells at 2,16, 24, and 48 h post-infection, using qPCR analysis. ZIKV amount inside cells is decreased upon knockdown of 7SK using ASO99 (f) and after its interaction with 7SK is blocked using a 2’O-methylated anti-sense oligo to 7SK (ASO88) (g)

7SK is a member of the 7SK ribonucleoprotein (snRNP) complex that regulates polymerase II transcription and is also known to control the transcription of HIV [18, 19]. We observed extensive binding of the 5′end of 7SK along DENV and ZIKV genomes with the strongest interactions occurring within the NS3 and NS5 regions (bases 4760 and 9250, Fig. 4d, e). To determine whether 7SK RNA binding affects ZIKV replication, we designed two independent antisense oligos (ASOs) to knock down 7SK RNAs inside cells [20]. Transfection of ASOs against 7SK did not cause cellular toxicity (Additional file 1: Fig. S4b) and resulted in more than 90% knockdown of 7SK RNA for each ASO (Additional file 1: Fig. S4c,d). We observed that the downregulation of 7SK RNA resulted in a decrease in ZIKV levels inside cells (Fig. 4f, Additional file 1: Fig. S4e), suggesting that 7SK RNA promotes ZIKV replication. To confirm that it is indeed the interaction between 7SK-ZIKV that impacts ZIKV replication, we blocked 7SK and ZIKV interaction by using a 2’O-methylated antisense oligo that binds to 7SK. 2’O-methylated oligos have been used to block RNA-RNA interactions without triggering the RNase H cleavage mechanism [21]. We designed the 2’O-methylated oligo to target the strongest 7SK-ZIKV interactions in the cell and observed that blocking 7SK-ZIKV interactions resulted in a decrease in ZIKV expression without significantly impacting endogenous 7SK level (Fig. 4g, Additional file 1: Fig. S4f), consistent with the 7SK knockdown experiments. These experiments support our hypothesis that 7SK is a virus-permissive factor that facilitates ZIKV replication.

The binding of DYNLT1 to DENV and ZIKV viruses inhibits viral replication

In addition to non-coding RNAs, we also performed knockdown experiments of coding RNAs that are found to bind to DENV and ZIKV viruses, including RPL37A and EEF1A1 (Additional file 1: Fig. S5a). Knockdown of these RNAs resulted in a decrease in the amount of Zika being produced in the cells, indicating that they may promote viral replication (Additional file 1: Fig. S5b). Additionally, we observed that an mRNA coding for Dynein Light Chain Tctex-Type 1 (DYNLT1), which was previously associated with DENV infection, binds to ZAFR, ZBRA, and DENV-1 strains. SPLASH interaction reads indicated that the 3′UTR of the DYNLT1 binds to the dumbbell region in the viral 3′UTR of the viral genome (Fig. 5a, b). To validate this interaction, we performed an electrophoretic mobility shift assay (EMSA [22]) and showed that ZIKV was directly associated with DYNLT1 3′UTR (Fig. 5d). Furthermore, introducing mutations in DYNTL1 3′UTR reduced this physical interaction (Fig. 5c, d). To determine the biological significance of this interaction, we performed a knockdown of DYNLT1 using specific siRNAs and monitored virus replication upon lower levels of DYNLT1. RT-qPCR confirmed the DYNLT1 mRNA was reduced by > 60% from 24 to 72 h after siRNA transfection (Additional file 1: Fig. S5c) and that no cytotoxicity was detected upon DYNLT1 siRNA transfection (Additional file 1: Fig. S5d). Importantly, the replication of different ZIKV strains was significantly increased on Huh7 cells upon DYNLT1 knockdown (Fig. 5e, f). These data suggest that in contrast to miR19a-3p and 7SK, DYNLT1 may play a role in repressing viral replication.

Fig. 5

DYNLT1 is an anti-ZIKV factor that binds to ZIKV 3′UTRs. a Top, Location along DYNLT1 that binds to ZIKV. Bottom, Location along ZIKV that binds to DYNLT1. b Schematic and predicted RNA structures on ZIKV before and after DYNLT1 binding to the 3′UTR. The DYNLT1 sequence is in blue and the ZIKV sequence is in black and red. c Structure model of DYNLT1-ZIKV interaction using RNAcofold program. DYNLT1 RNA is in blue while the ZIKV genome is in red. The mutations on DYNLT1 3′UTR RNA are indicated in purple above the original bases. d EMSA data showing ZIKV and DYNLT1 RNA-RNA interactions using different amounts of WT and MT DYNLT1 3′UTR RNA. e Barplots showing the effect of DYNLT1 knockdown on ZIKV replication. Huh7 cells were transfected with 20 nM of DYNLT1 siRNA or a control siRNA. At 24 h post-transfection, cells were infected with recombinant ZIKV strain PRVABC59 at a multiplicity of infection (MOI) 0.5. The infectivities of the supernatants 1 day after infection were determined by plaque assay. Means and standard deviations (SDs) from six independent experiments are presented. f DYNLT1 knockdown reduces the replication of the ZIKV African lineage. Huh7 cells were transfected with 20 nM of DYNLT1 siRNA or a control siRNA. At 24 h post-transfection, cells were infected with a ZIKV strain Dakar containing a nanoluciferase gene. At 24 h post-infection, luciferase signals were measured. The relative luciferase signals were obtained by normalizing the luciferase readouts of the DYNLT1 siRNA-treated groups to those of the siRNA control. Means and SDs from five independent experiments are presented. g Overexpression of DYNLT1 3′UTR suppresses ZIKV replication. Huh7 cells were transfected with DYNLT 3′UTR, DYNLT1 3′UTR with mutation, or pXJ vector (100, 200, or 400 ng per well in a 24-well plate). At 24 h post-transfection, cells were infected with recombinant ZIKV strain PRVABC59 containing a nanoluciferase gene at an MOI of 0.5. At 24 h post-infection, luciferase signals were measured. The relative luciferase signals were obtained by normalizing the luciferase readouts of each group to those of the pXJ vector control. Means and SDs from four independent experiments are presented. h DYNLT1 3′UTR inhibits ZIKV replication. Huh7 cells were transfected with 200 ng of DYNLT 3′UTR, DYNLT1 3′UTR with mutation, or pXJ vector. At 24 h post-transfection, cells were infected with recombinant ZIKV strain PRVABC59 containing a nanoluciferase gene at an MOI of 0.5. At given time points, cells were harvested. The relative luciferase signals were obtained by normalizing the luciferase readouts of each group to those of control at 2 h post-infection. Means and SDs from four independent experiments are presented. i Venn diagram showing the amount of overlap between RNAs that interact with all 5 DENV and ZIKV viruses and RNAs that interact with SARS-CoV-2

As the SPLASH data showed that it is the 3′UTR of DYNLT1 that is primarily responsible for binding to ZIKV, we cloned the 3′UTR of DYNLT1 into a mammalian expression vector and monitored ZIKV replication upon the overexpression of the 3′UTR of DYNLT1 using the ZIKV strain PRVABC-59 containing a nanoluciferase gene (ZIKV PRV-Nluc). We observed that luciferase signals were inhibited by the DYNLT1 3′UTR in a dose-dependent manner 1 day post-infection (Fig. 5g). Additionally, the over-expression of DYNLT1 3′UTR resulted in a decrease in ZIKV amounts at early time points post-infection (9 h post-transfection, Fig. 5h), suggesting that DYNLT1 3′UTR may act on early stages of the virus life cycle. Interestingly, mutations that reduce ZIKV-DYNLT1 interactions inhibited ZIKV replication to a lesser extent (Fig. 5g), suggesting that the pairing between DYNLT1 and ZIKV is important. We also confirmed using RT-qPCR and Western blot analysis that overexpression of DYNLT1 3′UTR or 3′UTR with mutations did not alter the endogenous DYNLT1 mRNA and protein levels (Additional file 1: Fig. S5e,f). These results supported that DYNLT1 3′UTR suppresses ZIKV replication through the pairing between DYNLT1 3′UTR and viral 3′UTR. In addition to ZIKV 3′UTR, SPLASH results showed that DYNLT1 3′UTR can bind to DENV-1 3′UTR (Additional file 1: Fig. S5g). Additionally, overexpression of DYNLT1 3′UTR also inhibits the replication of both DENV-1 and DENV-2, while overexpression of mutant DYNLT1 3′UTR failed to suppress viral replication (Additional file 1: Fig. S5h,i), confirming our observation that DYNLT1 acts generally on both ZIKV and DENV viruses. However, how DYNLT1 binding suppresses ZIKV/DENV replication remains to be explored. By performing motif searches and filtering for RBPs with experimental binding evidence, we observed that 24 RBPs could bind near the DYNLT1 binding site (Additional file 6: Table S5). As such, more experiments are needed to confirm whether DYNLT1 3′UTR restricts virus growth by directly disrupting the dumbbell structure of its 3′UTR or by blocking one of the RBPs from binding to the virus.

Lastly, to determine whether the host RNAs that bind to DENV/ZIKV are unique to these flaviviruses or whether they can also interact with other RNA viruses such as SARS-CoV-2, we overlapped the conserved DENV/ZIKV interactors with SARS-CoV-2 interactors. Out of 32 RNAs that bind to DENV/ZIKV, we observed that 22% (7) of them, including Y-RNA, EIF1, RPS3, RPLP1, SNORA62, SNORA80E, SNORD100, bind to all three—DENV, ZIKV, and SARS-CoV-2—genomes (Fig. 5i). Some of these host factors, such as Y-RNA and RPLP1 are already known to be key pro-viral host factors involved in RIG-I activation and translation respectively, and our study suggests that there could be additional regulation at the level of RNA between the host and viral genomes.