Effect of the SARS-CoV-2 Delta-associated G15U mutation on the s2m element dimerization and its interactions with miR-1307-3p

The SARS-CoV-2 Delta variant s2m acquired the G > U mutation at position 15

We performed a bioinformatics analysis to monitor sequence mutants of the SARS-CoV-2 s2m in cases worldwide, specifically targeting India, the United Kingdom, and the United States, three countries which reported a significant number of SARS-CoV-2 sequences to the Global Initiative on Sharing All Influenza Data (GISAID) database (Elbe and Buckland-Merrett 2017). For this mutational analysis, we used a custom R script that extracted and aligned the s2m sequences in all collected SARS-CoV-2 cases (Frye et al. 2023). Our bioinformatics analysis, which sampled SARS-CoV-2 sequences from April 2020 to October 2021, aligned and analyzed a total of 2,034,376 reported sequences to the GISAID database. A total 1980 sequences (0.097%) were removed due to insertions which disrupt the alignment. From this data set, we identified a G to U mutation at position 15 of the s2m, denoted as s2m G15U (Fig. 1, top). The earliest observation of this s2m mutation was in April 2020, within the United Kingdom and the United States, at 0.23% prevalence. While the s2m G15U prevalence remained low for the majority of 2020, we observed an increase in s2m G15U prevalence in April 2021 and May 2021 in the United Kingdom and the United States, respectively (Fig. 3A, orange triangle and blue circle, dashed lines). Interestingly, this increase was concurrent with the rise of the SARS-CoV-2 Delta variant (PANGO lineages B.1.617.2 + AY.*) (Fig. 3A, orange and blue, solid lines). The presence of the s2m G15U rose from 0.23% in March 2021 to a 95.97% plateau in June 2021 in the United Kingdom, and from 0.82% in April 2021 to a 92.40% plateau in July 2021 in the United States. This was simultaneous with the increase in the Delta variant, which rose from 2.03% in March 2021 to 95.68% in June 2021 in the United Kingdom, and from 0.62% in April 2021 to 92.86% in July 2021 in the United States. We determined that the s2m G15U and Delta prevalence maintained an average Pearson's correlation coefficient (R) of 0.9999 for both the United Kingdom and the United States from April 2020 through October 2021, spanning the range of our entire bioinformatics analysis (Supplemental Fig. 1A,B). Additionally, as a control, we studied the Alpha variant (PANGO lineages B.1.1.7 + Q.*) of SARS-CoV-2 and found that only 0.30% of Alpha sequences from November 2020 to August 2021, when this variant was dominant, contained the s2m G15U (Supplemental Fig. 2). This shows that the s2m G15U mutant, while occurring in Alpha at low prevalence, did not rise in prevalence as it had in the Delta variant.

View larger version:

• In this window
• In a new window

• Download as PowerPoint Slide

FIGURE 3.

(A) Prevalence of the s2m G15U mutation (dashed) and Delta variant (solid) within total reported SARS-CoV-2 cases in the United States (blue), United Kingdom (orange), and India (purple). (B) The s2m G15U mutant (purple dashed line, square markers) rose in prevalence before the rise of the Delta variant (purple solid line). Kappa (pink solid line) was also found to maintain the s2m G15U mutation and comprised of most cases within India that were not accounted for with Delta. The combination of the Kappa and Delta variants’ prevalence (green solid line) is tightly correlated with the s2m G15U mutant, mirroring trends later observed as Delta spread worldwide.

While the s2m G15U mutation and Delta variant show a strong positive correlation in prevalence throughout the spread of Delta worldwide, we noted a chronological discrepancy in the s2m G15U and Delta prevalence in India from November 2020 to March 2021, in which the s2m G15U had a noticeable difference in correlation from Delta (Fig. 3A, purple dashed square and solid lines). This difference is of particular interest given that the Delta variant originated in India. Here, the s2m G15U mutant increased significantly from 1.74% in December 2020 to 97.46% by May 2021, with the largest deviation between s2m G15U and Delta in March 2021, as the Delta variant comprised of only 12.06% of submitted cases compared to the 49.87% of s2m G15U. However, through May 2021, both data sets converged as they reached 94.86% and 97.46% for Delta and s2m G15U, respectively. Analysis of these two data sets revealed R = 0.96619 (Supplemental Fig. 1C, green), indicating a positive correlation with greater deviation compared to prior data sets. In efforts to determine the cause of this discrepancy between s2m G15U and Delta, we screened all Delta cases for the percentage of s2m G15U containing sequences and found that the Delta variant was comprised of 99.32% s2m G15U (Supplemental Fig. 3, blue dashed line). We also analyzed India for other dominant VOCs which may contain the s2m G15U and noted that the Delta variant arose concurrently with a sister variant, Kappa (PANGO lineage B.1.617.1), and found that the s2m G15U is present in 99.05% of all submitted Kappa cases (Supplemental Fig. 3). This finding indicated that Kappa had accounted for a majority of the s2m G15U cases that are not observed within Delta sequences. Combination of the Kappa and Delta prevalence, when plotted against s2m G15U, revealed R = 0.9997 (Supplemental Fig. 1C, purple), indicating a strong positive correlation between the rise in s2m G15U and overall combined incidence of the Kappa and Delta variants. Thus, we attributed the rise of the s2m G15U to a concurrent increase in Kappa and Delta cases. It is to be noted, however, that the Kappa variant was subsequently outcompeted by Delta in March 2021 (Fig. 3B, purple and pink lines), attributed to Delta's increased ACE2 binding affinity and antibody escape (Kannan et al. 2021; Tian et al. 2021; Liu et al. 2022). This trend continues following the drop in Kappa prevalence, in which the combined Kappa and Delta variant prevalence maintains its high correlation with the s2m G15U mutant (Fig. 3B, purple dashed line and green solid line).

The strong positive correlation between the s2m G15U and the Delta variant through its rise in India, and eventual spread worldwide, indicates that the s2m G15U rose in tandem with Delta and was found to be in an average of 99.32% of all Delta sequences. Thus, the presence of s2m G15U in such vast number of Delta SARS-CoV-2 cases and the high conservancy and predicted roles of the s2m in the viral life cycle, warranted its further characterization.

S2m G15U secondary structure characterization

To evaluate the effect of the G15U mutation on the secondary structure of the s2m, we utilized ¹H NMR spectroscopy, comparing the imino proton resonance region for the s2m and s2m G15U at 19°C. The SARS-CoV-2 s2m imino proton resonances were previously assigned and thus served as a baseline for our work (Wacker et al. 2020; Imperatore et al. 2022). The previously assigned imino protons corresponding to base pairs in the s2m lower stem (13.19 ppm—G37, 13.42 ppm—U38, and 11.62 ppm—G39) are identical between the two spectra (Fig. 4A). This indicates that the s2m G15U lower stem is not perturbed, consistent with the predicted secondary structures of s2m reference and s2m G15U (Fig. 1). Similarly, the resonance previously assigned to the U31 imino proton has the same chemical shift (13.93 ppm) in the spectrum for the s2m and s2m G15U. We used ¹H–¹H NOESY NMR spectroscopy on s2m G15U at 19°C to confirm these assignments and assign the remaining imino proton resonances which differ in chemical shift between the s2m reference and s2m G15U. As expected, the resonances assigned to U31 and U38 imino protons (13.93 ppm and 13.42 ppm) have a strong NOE with the H2 of the corresponding adenines they are base paired with (A13 H2 at 7.23 ppm for U31 NH and A4 H2 at 7.18 ppm for U38 NH) (Fig. 4B). In addition, the resonance at 13.66 ppm gives a single strong NOE at 7.39 ppm, a signature of an A–U base pair, allowing us to assign it to the U15 imino proton (Fig. 4B). This result confirms that the G15U mutation closes the base pair U15–A29, as predicted for the s2m G15U secondary structure (Fig. 1, top). Based on its neighboring base pairs, the U31 imino should give rise to NOEs with the G30 and G32 imino protons, and indeed, two cross peaks are observable, at 13.5 ppm and at 11.76 ppm, for the U31 imino proton resonance at 13.93 ppm. Based on its neighboring base pairs, the U15 imino proton should give rise to two NOEs: with the imino protons of G30 and G28. We only observe one clear NOE at 13.5 ppm from the U15 imino proton at 13.66 ppm. Since the imino at 13.5 ppm gives NOEs to both the U31 imino and the U15 imino, we unambiguously assign it to the G30 imino, as this is the only imino proton in the proximity of both U31 and U15 imino protons. Then, by exclusion, the second NOE observed for the U31 imino (13.93 ppm; 11.76 ppm) leads to the assignment of the resonance at 11.76 ppm to the G32 imino proton. This assignment is also supported by the observed NOE between the G32 imino at 11.76 ppm and the A13 H2 at 7.23 ppm. We also observed NOEs between the imino protons of U38 and G37 (13.42 ppm; 13.19 ppm), as well as between the imino of G37 and the A4 H2 at 7.18 ppm, confirming their previous assignment (Fig. 4C; Wacker et al. 2020; Imperatore et al. 2022). A weak NOE is also observed between the imino of G39 at 11.62 ppm and the A4 H2 at 7.18 ppm, confirming the G39 imino assignment (Fig. 4B). Finally, there are two remaining imino proton resonances at 12.03 ppm and 12.7 ppm, which do not give rise to NOEs that would allow us to unambiguously assign them to G28 and G17 iminos. The imino proton resonances of U1 and U2 are not visible in either spectrum due to the dynamic nature of these end base pairs (Imperatore et al. 2022). Based on these assignments we conclude that aside from the obvious mutation which closes the U15–A29 base pair, the s2m G15U secondary structure does not deviate far from that of the s2m.

View larger version:

• In this window
• In a new window

• Download as PowerPoint Slide

FIGURE 4.

(A) One-dimensional ¹H NMR spectroscopy of the SARS-CoV-2 s2m G15U (top) and SARS-CoV-2 s2m reference (bottom). The labeled resonances were assigned previously for the s2m reference and in the current work for s2m G15U (Wacker et al. 2020; Imperatore et al. 2022). (B) ¹H–¹H NOESY spectra of the imino proton resonances of s2m G15U with labeled NOEs between base paired U imino and A H2 and between G imino protons and neighboring A H2. (C) ¹H–¹H NOESY spectra of the imino proton resonances for the s2m G15U with labeled NOE between neighboring intrastrand and interstrand Gs and Us.

To further characterize the effect this newly formed base pair has on the stability of the s2m, we performed UV thermal denaturation experiments on s2m G15U and s2m, monitoring the absorbance changes at 260 nm (Fig. 5). The melting temperature, T_m, of each structure was determined from the first derivative plots (Supplemental Fig. 4) to be 44°C and 47°C for the s2m and s2m G15U, respectively. This indicated that, as expected, the closing of the U15–A29 base pair in the s2m G15U upper stem increases the overall stability of the s2m G15U (Figs. 1, 5), which may have implications for the virus’ ability to dimerize and form ED structures. In our previous work, we established that the s2m elements of SARS-CoV and SARS-CoV-2 have different dimerization abilities and determined that this is due to the presence of a single nucleotide mutation in s2m G31U located in the s2m upper stem (Imperatore et al. 2022). Such dimerization events in viruses have been shown as critical for genome packaging, inducing conformational changes for translational regulation, reverse transcription processes, as well as being important for genetic recombination (Berkhout and van Wamel 1996; Mikkelsen et al. 2000; Aagaard et al. 2004; Moore et al. 2009; Tengs et al. 2013; Terada et al. 2014). Therefore, next we evaluated the effect of the single nucleotide mutation, the G15U, on the Delta SARS-CoV-2 s2m dimerization.

View larger version:

• In this window
• In a new window

• Download as PowerPoint Slide

FIGURE 5.

UV thermal denaturation of the s2m and s2m G15U at 1 µM in 10 mM cacodylic acid pH 6.5 and absence of MgCl₂. The melting temperatures, T_m, were 44°C and 47°C for the s2m reference and s2m G15U, respectively.

Characterization of the SARS-CoV-2 s2m G15U dimerization

Our previous work established the presence of a 4 nt palindromic “GUAC” sequence within the terminal loop of the s2m element that participates in loop–loop kissing dimerization (Imperatore et al. 2022). These KD structures are formed through base-pairing between two RNA hairpin loops that contain complementary unpaired nucleotides. Kissing dimers forming with palindromic sequences have been previously reported in other viral systems to be essential for template switching and recombination, a process important for repairing genetic mistakes and increasing diversity within the virus, which is critical for the evolution of new viral strains (Mikkelsen et al. 2000; Andersen et al. 2003). These KDs can then be converted to a thermodynamically stable ED structure through the chaperone activity of the capsid protein when in the context of the entire viral genome (Fig. 2). However, in vitro, short isolated KDs have been reported to form ED conformations spontaneously, such as in hepatitis C virus (HCV) and human immunodeficiency virus (HIV-1) (Muriaux et al. 1996; Mihailescu and Marino 2004; Shetty et al. 2010). In the case of HIV-1, recombination events have led to high rates of immune escape and drug resistance. Thus, given the potential of s2m to be involved in recombination events it is important to determine if the G15U mutation affects its dimerization, as SARS-CoV-2 is still developing new variants (Morris et al. 1999; Rambaut et al. 2004; Dubois et al. 2018).

Therefore, to analyze the G15U mutation's effect upon dimerization of the SARS-CoV-2 s2m, we used TBE/TBM nondenaturing gel electrophoresis, a method used previously in various viral systems to distinguish between KD versus duplex conformations of dimer initiation sites (Berkhout et al. 2002; Li et al. 2006; Sun et al. 2007; Kenyon et al. 2013). Conditions which maintain MgCl₂ concentrations (TBM) stabilize the formation of the KD; however, upon chelation of the Mg²⁺ ions by EDTA in TBE, the KD dissociates, whereas the duplex conformation, which does not require Mg²⁺ ions for stability, is retained. When MgCl₂ is retained in the TBM gel, a prominent monomer band (indicated by M) is present for both the s2m and s2m G15U (Fig. 6, right, arrows 1 and 1′). There are two dimer bands present for the s2m reference (Fig. 6, right, arrows 2 and 3), and a single faint dimer band is present for s2m G15U (Fig. 6, right, arrow 3). Our previous report has established that these two upper bands in the s2m reference (Fig. 6, right, arrows 2 and 3) correspond to a KD and ED, respectively (Imperatore et al. 2022). Although both the KD and ED conformations contain the same number of nucleotides, their migration patterns can differ in native gels as these separate RNAs based on both size and shape. Duplex structures are linear, whereas KDs can be present in either a bent or linear state. This would give rise to a different shape that may affect the migration pattern, explaining the difference in the two bands (Ennifar et al. 1999; Rist and Marino 2002). Although unlikely in the context of the entire genome, when the s2m is isolated in its 41 nt state, spontaneous conversion to an ED structure is possible and has been reported in studies involving isolated dimerization initiation sites in other viruses (Muriaux et al. 1996; Mihailescu and Marino 2004; Shetty et al. 2010). Similarly, it is interesting to note that, despite having the same size, the monomer bands for the s2m reference and s2m G15U migrate differently in the TBM gel (Fig. 6, right, arrows 1 and 1′), suggesting that they have different monomer conformations in the presence of MgCl₂. Our MD simulations for the s2m revealed that it has a kinked 3D structure; however, the additional U15–A29 base pair in the upper stem of s2m G15U results in a linear structure that lacks the L-shaped kink, rationalizing the increased migration of s2m G15U through the TBM gel as compared to the original s2m form (Kensinger et al. 2022; Makowski et al. 2023). Thus, we cannot determine if the single dimer band observed for the s2m G15U in the TBM gel (Fig. 6, right, arrow 3) originates from a mixture of linear KD and duplex conformations which will migrate at the same position, or from the duplex conformation formed spontaneously from an unstable KD during incubation and migration through the gel (Shetty et al. 2010).

View larger version:

• In this window
• In a new window

• Download as PowerPoint Slide

FIGURE 6.

TBE (left) and TBM (right) gel nondenaturing electrophoresis analysis of the SARS-CoV-2 s2m reference and the s2m G15U. The two s2m elements were compared using Mg²⁺ dependence where 1 µM RNA was incubated in the presence of 1, 5, or 10 mM MgCl₂ for 1 h, then split between two gels—TBE which chelated the Mg²⁺ ions, and TBM that maintained the Mg content. In the TBE, both s2m elements prefer their monomeric (M) states (left, arrow 1). In TBM, the monomeric states are still preferred (right, arrows 1 and 1′). Different migration patterns may be observed due to different shapes of the s2ms, as described in the text. Arrows 2 and 3 indicate the KD and ED, respectively. Degradation products are noted with (*).

In the TBE gel, after an hour incubation with MgCl₂, a prominent monomer s2m band exists for both the s2m reference and the s2m G15U mutant (Fig. 6, left, arrow 1). We also observed faint monomer bands marked by (*) originating from degradation products. There are also faint dimer bands present in the TBE gel, which likely originate from the duplex conformations of the intact s2m elements and of its minor degradation products (see below). It is important to note that the TBM gel was run for 4 h in the presence of MgCl₂, as opposed to the TBE gel that was run for 2 h, allowing more time for spontaneous conversion of the s2m KD to the duplex conformation, and hence, there is a more prominent duplex band in the TBM gel (Fig. 6, right, arrow 3), as compared to the TBE gel.

However, in the context of the full virus, this conversion from the KD to the ED conformation would not occur spontaneously and would instead require the assistance of the N protein's chaperone activity to facilitate the transition. This process is dependent on the overall stability of the s2m element as base pairs in each monomer will have to be broken, disrupting the monomeric hairpin structure. In coronaviruses, the N protein is responsible for packaging the RNA into ribonucleoprotein complexes to form the mature virion prior to being exported from the host cell (Zúñiga et al. 2007; Aduri et al. 2013). Our previous work has shown that the SARS-CoV-2 N protein selectively acts as a molecular chaperone to convert the s2m KD structure to the stable duplex by a similar mechanism observed for the dimer initiation sites of HIV-1 and HCV (Rist and Marino 2002; Shetty et al. 2010; Shimakami et al. 2012; Imperatore et al. 2022). In HIV-1, a retrovirus, dimerization of the RNA is important for the packaging of two copies of the genome in the mature virion (Dubois et al. 2018). In HCV, a single-stranded RNA virus, this dimerization is suggested to be a regulator of RNA replication and viral genome packaging (Shetty et al. 2010). Additionally, the dimer initiation sites in HIV-1, HCV, and in the murine leukemia virus have been proposed to be implicated in recombination events (Torrent et al. 1994; Greatorex 2004; Masante et al. 2015). We proposed previously that the s2m dimerization mediated by the viral N protein may facilitate SARS-CoV-2 recombination; thus, here we investigated if the G15U mutation affects the s2m dimerization properties (Imperatore et al. 2022).

To test the chaperone ability of the N protein to convert the s2m G15U to the duplex structure, as compared to the s2m reference, a 1 µM s2m RNA sample was incubated with 1 mM MgCl₂ for 30 min, followed by the addition of 2 µM N protein for an additional 30 min. After the incubation with the N or control proteins, proteinase K was added to degrade the proteins prior to loading the samples on the TBE and TBM gels. As demonstrated previously, in TBM both the s2m reference and s2m G15U preferred their monomeric states (Fig. 7, right, arrows 1 and 1′) in the absence of the N protein; however, in the presence of the N protein the two upper bands in the s2m reference became more prominent (Fig. 7, right, arrows 2 and 3; Imperatore et al. 2022). Similarly, the dimer band of s2m G15U, which was barely visible on its own, becomes clear in the presence of the N protein (Fig. 7, right, arrow 3). As discussed above, we cannot distinguish if the s2m G15U KD is unstable (hence the dimer band indicated by arrow 2 is absent in s2m G15U) or if it is a linear structure that migrates at the same position with the duplex conformation. To try to distinguish between these two conformations, we have performed similar experiments incubating the samples with the N protein for times ranging between 1 and 30 min but did not observe any differences in the formation of the dimer band in s2m G15U (data not shown). The TBE/TBM native PAGE method is limited in that it does not allow us to access shorter incubation times with the N protein where we could potentially observe differences in the dimerization of s2m G15U by being able to distinguish between the KD and the ED conformations. When the Mg²⁺ ions are chelated in the TBE gel, a prominent dimer band remains for both s2m and s2m G15U that were incubated with the N protein (Fig. 7, left, arrow 3), indicating that the N protein facilitates the conversion to the duplex conformation in both samples. Note that the faint dimer bands marked by (*) do not line up with the main duplex bands (Fig. 7, left, arrow 3), indicating that they originate from the dimerization of the degradation products. The addition of the P protein or bovine serum albumin (BSA) to the sample had no effect, validating the specificity of the N protein as a molecular chaperone.

View larger version:

• In this window
• In a new window

• Download as PowerPoint Slide

FIGURE 7.

Nondenaturing gel electrophoresis of the s2m reference and s2m G15U converted to an ED structure by the viral N protein. Both s2m elements were incubated with 1 µM MgCl₂ for 30 min, an additional 30 min with the N protein, and 5 min with proteinase K to digest the proteins prior to the samples being split in TBE (left) and TBM (right) conditions for electrophoresis. In both cases, the s2m remains preferentially monomeric (arrows 1 and 1′). The N protein converted both elements to an ED structure (arrow 3), and in the TBM gel, the N protein also stabilized the KD of the s2m reference (right, arrow 2). However, the band corresponding to the KD of the s2m G15U sample is still not observed. The P protein and BSA were used as controls and had no effect on the conversion of either s2m KD to the duplex conformation.

Based on these results, we conclude that while the s2m and s2m G15U are converted by the N protein to the duplex conformation, the s2m G15U KD conformation is either not as stable or forms a different conformation than the s2m KD. Specifically, the absence of an s2m G15U KD band at the same position observed for the s2m KD (Fig. 7, right, arrow 2) indicates that the s2m G15U KD could be linear, migrating at the same position with the linear duplex conformation, different than the KD structure of the s2m. Alternately, this could indicate that the s2m G15U KD is able to convert spontaneously to the duplex conformation in the absence of the N protein more efficiently than the s2m; however, considering that the s2m G15U upper stem is more stable and must be opened during this conversion, this scenario is less likely. A difference in KD conformation is more probable as the monomers of these two s2m elements have a different shape based on their different migration in the TBM gel (Fig. 6 right, arrows 1 and 1′). Moreover, our MD simulation findings show that the monomer s2m G15U is linear, while the original s2m monomer is kinked; this is consistent with this hypothesis (Makowski et al. 2023). Despite these observed differences in the s2m G15U KD versus the s2m reference KD, we further demonstrate that this mutation does not affect the ability of the viral N protein to convert the s2m G15U KD to a stable duplex structure, which is the in vivo process of duplex formation. It has been suggested that the N protein preferentially binds to double-stranded RNA segments in the central region of the s2m, specifically nucleotides G28, G30, and U31, making a stable protein–RNA complex for virion packaging (Padroni et al. 2022; Morse et al. 2023). Therefore, the closing of the U15–A29 base pair, despite making the s2m G15U upper stem more stable, which would be unfavorable for the conversion of the s2m G15U KD to the duplex conformation, may compensate for this unfavorable effect by further assisting in N protein binding.

The findings that the G15U mutation does not affect the ability of the N protein to convert the s2m G15U KD to the duplex conformation are contrasting the effect another single nucleotide mutation in the s2m upper stem (G31U) has on this conversion. In our previous work, we showed that the G31U mutation that exists between the SARS-CoV and SARS-CoV-2 s2m was responsible for drastic differences in the s2m dimerization, even in the presence of the N protein (Imperatore et al. 2022). The SARS-CoV s2m readily forms KDs that could not be converted to the ED by the N protein. Conversely, the SARS-CoV-2 s2m which differs by only 2 nt (only one of which, G31U, we showed is responsible for the dimerization differences) from SARS-CoV s2m, exists primarily as a monomer that forms KDs which are subsequently converted to the duplex conformation by the N protein (Imperatore et al. 2022). While the single nucleotide G31U mutation was shown to be responsible for this difference in s2m dimerization properties between SARS-CoV s2m and the SARS-CoV-2 s2m, we show here that the s2m G15U single nucleotide mutation present in the Delta variant has limited effects on the element's dimerization (Imperatore et al. 2022).

SARS-CoV-2 s2m G15U interactions with miR-1307-3p

We and others previously identified two miR-1307-3p binding sites on the s2m element in the SARS-CoV-2 3′-UTR (Fig. 8A), and we revealed its specificity for this miR (Arisan et al. 2020; Balmeh et al. 2020; Chan et al. 2020; Alam and Lipovich 2021; Imperatore et al. 2022). The exact role of miR-1307-3p in the SARS-CoV-2 infection is not known; however, it has been shown that miR-1307-3p is the highest expressed microRNA in Vero cells infected by SARS-CoV-2, and the levels of miR-1307-3p are elevated in severe cases of COVID-19 as compared with mild cases (Arisan et al. 2022; Zarei Ghobadi et al. 2022). Thus, we investigated if the G15U mutation within the s2m has any effect upon its interactions with miR-1307-3p. The binding interactions between miR-1307-3p and s2m G15U as compared to the s2m reference were assessed by TBE/TBM native PAGE. The s2m and s2m G15U were boiled and snap-cooled, then incubated with 1 mM MgCl₂ for 30 min, after which miR-1307-3p was added in 1:1 and 1:2, and 1:4 s2m:miR ratios and incubated for an additional 30 min at 25°C. The samples were then split and electrophoresed in TBE and TBM conditions, respectively. In the TBE gel, two bands are present for miR-1307-3p, a monomer and dimer (Fig. 8B, left, arrows 1 and 3), while both the s2m and s2m G15U are present as monomers (Fig. 8B, left, arrow 2). In the presence of miR-1307-3p, two upper bands become apparent for the s2m, which were assigned to the complex of s2m with one miR-1307-3p molecule (Fig. 8B, left, arrow 4) and two miR-1307-3p molecules, respectively (Fig. 8B, left, arrow 5), given the relative size of the bands compared to the molecular weight marker (Imperatore et al. 2022). However, the intensities of these upper bands are significantly decreased for s2m G15U, suggesting that miR-1307-3p binds less to this mutant. Note that there is a faint band at almost the same position, indicated by arrow 6, in the free s2m and s2m G15U lanes that corresponds to the duplex conformation (82 nt). This is very close in size to the complex formed by s2m with two miR-1307-3p molecules (85 nt), so it can no longer be distinguished as miR is titrated into the sample at higher concentration ratios. When the Mg²⁺ ions are retained in the TBM gel, the complex of the s2m with two miR-1307-3p molecules is stabilized (Fig. 8B, right, arrow 5) and a new higher-order complex appears, which we assigned to an s2m dimer interacting with two miR-1307-3p molecules, as described in our previous work (Fig. 8B, right, arrow 6; Imperatore et al. 2022). Further, arrows 7 and 8, which appear faintly, may be the result of end-to-end RNA LEGO complexes, tertiary RNA assemblies that are Mg²⁺ dependent, mediated by the palindromic sequence CGUG at the 3′ end of miR-1307-3p, (Fig. 8A; Horiya et al. 2003). In contrast, while Mg²⁺ ions stabilize the complex formed by s2m G15U with two miR-1307-3p molecules (Fig. 8B, right, arrow 5), the band corresponding to the complex formed with s2m G15U is less prominent, and the higher-order complexes seen with the s2m are not apparent. This result suggests that the closing of the U15–A29 base pair creates a more stable s2m upper stem, making it harder for the miR-1307-3p to invade it.

View larger version:

• In this window
• In a new window

• Download as PowerPoint Slide

FIGURE 8.

(A) Predicted miR-1307-3p binding sites to the s2m. The two miR-1307-3p molecules are shown in purple and green, respectively, whereas s2m is shown in gray. (B) Nondenaturing TBE (left) and TBM (right) gels of s2m and s2m G15U binding to miR-1307-3p in 1:1 and 1:2 and 1:4 s2m:miR ratios in the presence of 1 mM MgCl₂. In TBE, both s2m elements prefer their monomeric states (arrow 2), whereas the miR-1307-3p exists as a monomer and dimer (arrows 1 and 3). As miR is added, upper complexes are formed between the two (arrows 4 and 5). Arrow 6 reveals the presence of the ED of the s2ms. Overall, the band intensity of the complexes formed by s2m G15U with miR-1307-3p is reduced. In TBM conditions, the s2m reference can form a stable complex with two miR-1307-3p (arrow 5) and a complex formed by two s2m each bound to one miR-1307-3p is also apparent (arrow 6). Arrows 7 and 8 indicate higher molecular weight complex that could be potentially formed by RNA LEGO type interactions mediated by the palindromic sequence GUGC at the 3′ end of miR-1307-3p interacting with another miR, as it is bound to a s2m. Arrows 4 and 4′ show the KD and duplex formed by the s2m reference, respectively. The monomers of s2m and s2m G15U migrate differently (arrows 2 and 2′), as described in the text.

To determine the relative binding affinity difference between the s2m and s2m G15U, a fluorescently tagged miR-1307-3p-DY547 was used in native TBM gels monitoring the fluorescence of DY547 (excitation 558 and emission at 574 nm) (Fig. 9A,B, left), followed by staining in SYBR Gold to observe all RNA bands and overlay with the fluorescent image (Fig. 9A,B, right). The free miR-1307-3p-DY547, exists as a mixture of monomer and dimer, as monitored by the DY547 fluorescence (Fig. 9A, arrows 1 and 3). However, since the DY547 tag adds bulkiness, miR-1307-3p-DY547 has an altered migration pattern as compared to the unlabeled miR-1307-3p. Three upper bands are observable (Fig. 9, arrows 6, 7, and 8) when increasing ratios of miR-1307-3p-DY547 were incubated with s2m. Similar to the TBM binding gel with the unlabeled miR, we assign the bound complex of one s2m to two miR (arrow 6), a s2m dimer with each s2m bound by one miR (arrow 7), and finally a magnesium-stabilized RNA LEGO interaction complex formed through interactions of the palindromic sequence found at the end of the miR-1307-3p as it is bound to a s2m (arrow 8) (Horiya et al. 2003). To quantify these results, the intensities of the bands corresponding to the miR-1307-3p-DY547 complexes with s2m (arrows 6, 7, and 8) were measured using the ImageJ gel analysis software (Gallagher 2010) and normalized to the relative intensity of the free miR-1307-3p-DY547 lane (arrows 1 and 3). These experiments were performed in triplicate, and the averaged, normalized relative band intensities were in the range of 0.19–1.32 for the miR:s2m ratios of 0.25–4 miR (Fig. 9A, left, bottom of the gel and Fig. 9C). In the s2m G15U gel, like in the TBM binding gel with the unlabeled miR, only the complex formed by one s2m G15U with two copies of miR-1307-3p-DY547 is observed (Fig. 9B, arrow 5). After normalizing the intensity of the band corresponding to the miR-1307-3p-DY547 complex with s2m G15U (arrow 5) to that of the free miR-1307-3p-DY547 (arrows 1 and 3), the average, normalized relative band intensity across the triplicate experiments was in the range of 0.07–0.48 for the miR:s2m G15U ratios of 0.25–4 miR (Fig. 9B, left, bottom of the gel and Fig. 9C). Taken together these results indicate that miR-1307-3p-DY547 has reduced binding to s2m G15U, as compared to the s2m reference (Fig. 9C).

View larger version:

• In this window
• In a new window

• Download as PowerPoint Slide

FIGURE 9.

Nondenaturing TBM gel of s2m and s2m G15U binding to a fluorescently tagged DY547 miR-1307-3p in ratios ranging from 1:0.25 to 1:4 s2m:miR in the presence of 1 mM MgCl₂. The gels were visualized with a fluorescent imager at 600 nm wavelength (left), stained with SYBR Gold, then overlaid (right). (A) In the s2m binding gels, the blue arrows 1 and 3 correspond to the free DY547 miR-1307-3p in its monomeric and dimeric state, respectively. Arrow 2 indicates the monomeric s2m, whereas arrows 4 and 5 represent the KD and ED, respectively. Bound complex arrows that were used for quantification are indicated in red and show the s2m bound to two miR (arrow 6), the s2m dimer bound to two miRs (one miR to one s2m) (arrow 7), and complexes formed by end-to-end palindromic sequence interactions of the miR bound to the s2m (arrow 8). (B) When the miR is bound to the s2m G15U in TBM conditions, the miR exists as both a monomer and dimer, arrows 1 and 4, respectively. Arrow 2 indicates the monomeric s2m G15U, whereas arrow 3 represents a dimeric form, and these are only present in the right image of the gel stained in SYBR Gold due to it lacking a fluorescent tag. A bound complex between the s2m and two DY547 miR-1307-3p is indicated by arrow 5. Relative band intensities of the bound complexes (red arrows) were calculated using ImageJ as described in the text, and the intensities of the shown gels are listed. (C) The relative band intensities were averaged across the triplicate experiments (representative images shown in A and B). Across all miR ratios, the s2m consistently had a higher relative band intensity of its complexes than the s2m G15U. The standard deviation is included.

Disease severity in patients with SARS-CoV-2 has fluctuated throughout the pandemic, ranging from asymptomatic to mild to severe. In severe cases of COVID-19, hyperactive immune responses lead to the overproduction and release of cytokines, termed the “cytokine storm.” These cytokines are proteins critical for the regulation of the host's defense mechanisms (Vial and Descotes 1995). An example of cytokines includes interleukins (ILs), such as IL18, and their receptors (e.g., IL6R, IL10RA, IL10RB, IL12RB2, IL17RA), and it has been shown that their overexpression is linked to the onset of cytokine storms in COVID-19 patients (Agarwal et al. 2015). A suggested regulator of these ILs and their receptors is host miR-1307-3p (Mehta et al. 2020; Yang et al. 2020; Imperatore et al. 2022). MiR-1307-3p has also been proposed to reduce glucose regulating protein 78 (GRP78), a protein thought to aid in binding of the SARS-CoV-2 S protein to the host cell's receptors (Elfiky 2020; Ibrahim et al. 2020; Carlos et al. 2021). By allowing the virus to hijack this miR, GRP78 abundance in the cells would be increased, aiding viral entry into the cell. Viruses have been shown to hijack host miRNAs for their own benefit (Guo and Steitz 2014; Trobaugh and Klimstra 2017). One example of this is HCV, which uses miR-122 to enhance its replication (Jopling et al. 2005; Shimakami et al. 2012). Previously, upon finding that the s2m binds two copies of miR-1307-3p, we proposed that this miR is hijacked by SARS-CoV-2 for the benefit of the virus (Imperatore et al. 2022). However, considering these new results, our previous hypothesis must be amended as it seems that the virus has evolved to reduce its interactions with this miR in the Delta variant. We show that the s2m G15U binds to less miR-1307-3p as compared to the s2m reference and suggest the newly formed U15–A29 base pair formation as causative for these diminished interactions, as the s2m upper stem must be opened for new base pairs to be formed with miR-1307-3p. This implies that miR-1307-3p may instead have a detrimental effect on the virus, as our data supports that the virus preferentially evolved to reduce this interaction. Recent work has identified that host miRs may target the viral genome as part of the immune response, in which case, the virus would preferentially evolve to evade this attack (Lecellier et al. 2005; Huang et al. 2007; Otsuka et al. 2007; Skalsky and Cullen 2010). MiRs normally bind to the target mRNA 3′-UTR to signal its degradation and/or to regulate its translation (O'Brien et al. 2018). Often, the miR binding sites on mammalian mRNAs are highly conserved, preserving the regulatory role of these miRs (Friedman et al. 2009). Given the high conservation of the s2m in the 3′-UTR of several viral families, it is possible that the host immune response could use miRs to target this element, and several miRs are dysregulated upon SARS-CoV-2 infection (Farr et al. 2021). Thus, by targeting the s2m of SARS-CoV-2, the host may be utilizing miR-1307-3p to alter the translation of viral proteins (Panda et al. 2022). It is interesting to note that since the SARS-CoV-2 subgenomic RNAs share the same 3′-UTR with the genomic RNA, miR-1307-3p could target the translation of both the viral structural and nonstructural proteins. The host miR-1307-3p has been shown to be up-regulated in SARS-CoV-2 infected cells, up to a 38-fold increase, and is implicated in the silencing of SARS-CoV-2 post-transcriptionally (Arisan et al. 2022; Panda et al. 2022; Zarei Ghobadi et al. 2022). Additionally, increased expression of miR-1307-3p has been shown to result in decreased replication of the SARS-CoV-2 genome (Soremekun et al. 2020; Zarei Ghobadi et al. 2022). Therefore, mutations which decrease the binding interactions of miR-1307-3p with the s2m, such as those we observe in the s2m G15U mutant, would be optimal for overall viral fitness of Delta compared to the previous variants that contain the s2m.

Interestingly, our recent work has shown that the Omicron variant, the most prominent viral strain to emerge after Delta, contains a 26 nt deletion within the s2m element that eliminates the miR-1307-3p binding sites (Colson et al. 2022; Frye et al. 2023). Omicron has also been shown to have further reduced disease severity and increased infectivity when compared to Delta and other prior variants (Chaguza et al. 2022). We emphasize that our results could be evidence of a stepwise reduction of the miR-1307-3p binding sites as the s2m has evolved to reduce miR-1307-3p binding in Delta and subsequently removed the binding sites with the 26 nt deletion in Omicron (Frye et al. 2023). This trend supports the notion that SARS-CoV-2 is evolving toward reduced miR-1307-3p binding. Despite the lack of knowledge regarding this interaction, we suspect that if miR-1307-3p is in fact a crucial player in the antiviral immune response, the Omicron variant of SARS-CoV-2 greatly benefits through an overall reduced miR-1307-3p binding. Further, if miR-1307-3p is involved in the “cytokine storm,” this variant also benefits through lower disease severity, which is seen in patients diagnosed with the Omicron variant. Sequence alignments of the SARS-CoV-2 variants to other coronaviruses, specifically human “cold-like” viruses, revealed that the Omicron deletion within the s2m is similar to those of the common cold (Arisan et al. 2022). This further supports that the removal of this region, and the decreased (or lack of) ability to bind miR-1307-3p may be a similar evolutionary step taken by SARS-CoV-2 to prevent severe, and potentially deadly, symptoms in the host by improving its evasion of the host immune response.

However, it is to be acknowledged that this data alone cannot prove directly if this s2m G15U mutant does indeed confer a selective advantage to SARS-CoV-2, or is a result of a hitchhiking effect, where the mutations within the S protein, which conferred an advantage to the Delta variant, happened on a genome that already had the s2m G15U mutation. Nonetheless, the findings that the Omicron variant, which is not a direct descendant of the Delta variant, acquired the 26 nt deletion that removed the miR-1307-3p binding sites suggest that the abolishment of the miR-1307-3p interactions with the genomic and subgenomic viral RNA is advantageous to the virus. Further research is needed to elucidate the details of the role of miR-1307-3p in the host immune response to SARS-CoV-2 infections and viral genome dimerization.

In summary, we showed here that the s2m G15U mutation is present in 99.32% of the SARS-CoV-2 Delta variant, and that while this mutation does not affect the s2m dimerization properties in the presence of the N protein, it does reduce its ability to interact with host miR-1307-3p, potentially having implications on the host response to SARS-CoV-2 infection and viral evolution to evade this response.