Applicability of the nearest-neighbour model for pseudoknot RNAs

Abstract

The validity of the nearest-neighbour (NN) model was verified in an RNA pseudoknot (PK) structure. The thermodynamic parameters of the second hairpin stem (S2) region, which separates the PK from a hairpin structure, were monitored using CD and UV melting. Different PKs with identical NN base pairs in the S2 region exhibited similar thermodynamic parameters, highlighting the validity of the NN model in this RNA tertiary structure motif.

---

The function of living cells largely depends on the structural stability of nucleic acids because of their central role in biological processes, such as replication, transcription, and translation.¹ Thus, predicting structural stability remains one of the most challenging and critical aspects of nucleic acids. Generally, the stability of nucleic acids is defined as the free energy change (ΔG°) associated with the folding and unfolding process, which is governed by the changes in enthalpy (ΔH°) and entropy (ΔS°) of the corresponding structure.²

Thermodynamic parameters (ΔH°, ΔS°, and ΔG°) that are based on the nearest-neighbour (NN) model, proposed by Tinoco et al.,³ are often used in calculating the thermodynamic stabilities of nucleic acid duplexes consisting of any Watson–Crick base pair composition. This model works on the simple principle that the summation of the free energy changes for the formation of all adjacent NN base pairs along with that for helix initiation results in the thermodynamic stability of the entire duplex. Although the parameters were originally determined in a buffer containing 1 M NaCl,⁴ several research groups, including ours, have quantitatively measured the NN parameters for DNA/DNA,⁵ RNA/RNA,⁶ and RNA/DNA⁷ duplexes at 100 mM NaCl under physiological conditions. The NN model has also been applied to other simple structural motifs such as hairpins,⁸ mismatch base pairs,⁹ and a few triplex structures.¹⁰ The NN parameters for DNA/DNA^5b,11 and RNA/RNA duplexes¹² were measured under molecular crowding conditions, highlighting the applicability of the NN model in different media with varying physicochemical properties. However, the NN model has not been applied to predict complex and non-canonical structures because of the tertiary interactions involved in these motifs.

Being intrinsically single-stranded, RNAs form various functional motifs with different tertiary interactions and play key roles in biological processes.¹³ These tertiary interactions are built on the preformed secondary structures consisting of continuous Watson–Crick base pairs. Thus, predicting the thermodynamic stability of secondary structure units involved in tertiary structure formation is key to perceiving the physicochemical treatment between simple secondary structures and the complex tertiary structures originating from them.

Pseudoknots (PKs) are one of the prevalent tertiary structures formed on RNAs (Fig. 1) and contain a minimum of two helical stems (S1 and S2) connected through a single-stranded loop (L3).¹⁴ Almost all naturally occurring functional RNAs exhibit a huge abundance of PK structures.^14b For instance, nearly 66% of known riboswitches comprise PKs, which play a vital role in their gene regulatory functions.¹⁵ The reliance of PK functions on their structural stabilities indicates the importance of predicting the stability of PKs.^14,16 The structural stabilities of PKs have been explored by various computational models¹⁷ and experimental techniques such as NMR spectroscopy,¹⁸ UV spectroscopy,¹⁹ laser temperature-jump spectroscopy,²⁰ and single-molecule atomic force microscopy.²¹ However, most studies have investigated only a particular PK structure, highlighting the need for a general approach for predicting the PK structural stability of interest.

Fig. 1 Schematic representation of the PK formation and secondary structures of PKs analysed in this study. In the secondary structures of PK-A and PK-B, “X–Y” (in the orange square) indicates any composition of the Watson–Crick base pair in the S2 region. HP-A and HP-B are control hairpins of PK-A and PK-B, respectively. A14 in the PK-A series and A16 in the PK-B series, which are numbered in red, are the L2 nucleotides.

Among the PK structures, simple H-type PKs make up the majority of the natural PKs.^14a,22 Generally, a H-type PK structure undergoes a stepwise formation from unfolded single-stranded RNA to a hairpin structure with the S1 stem followed by the formation of the PK, which consists of both S1 and S2 stems (Fig. 1).²³ Thus, the formation of the S2 stem region separates a conventional hairpin from the complex PK, and the structural stability of the S2 stem region is expected to dictate PK formation. Hence, we hypothesised that the NN model could be used to predict the thermodynamic parameters of the S2 stem region. However, before its implementation for PKs, the validity of the NN model needs to be verified as there are very few examples of an NN model being applied to higher-order RNA structures.¹⁰

In this study, the thermodynamic parameters for the formation of S2 stem regions in two different PKs (PK-A and PK-B, Fig. 1) were analysed under physiological buffer conditions containing 10 mM Na₂HPO₄ (pH 7), 100 mM NaCl, and 1 mM Na₂EDTA through CD and UV melting. PK-A and PK-B differ in the length of the S1 stem region with 5 and 7 G–C base pairs, respectively. All the designed PKs were derived from a naturally occurring PK in mouse mammary tumour virus (MMTV), i.e. MMTV VPK (Fig. S1 and Table S1, ESI†).²⁴ In our study, a uracil nucleobase at position 13 of MMTV VPK was replaced by cytosine to form a G–C Watson–Crick base pair at the end of the S2 stem. This eliminated any probable G-bulge formation by forming an A–U base pair with the adjacent adenine nucleobase in L3 (Fig. S1, ESI†). The dangling uridine nucleobase at the 3′ end of MMTV VPK was also removed to prevent non-specific interaction with L1 (Fig. S1, ESI†). A set of two PKs, PK-A-1a and PK-A-1b, were designed by placing compositions with identical NN base pairs in the S2 stem region. Similarly, seven sets of identical NN base pairs were placed at the S2 stem to design PK-A-1 to PK-A-7 (Table S1, ESI†). In addition to the PK-A series, four sets of identical NN base pairs were placed at the S2 region of PK-B to design PK-B-1 to PK-B-4 (Table S1, ESI†). Importantly, PK-A-1a and PK-B-1a have the same NN base pairs in their S2 region, which allowed us to explore the effect of S1 on the stability of S2. We also designed hairpin structures containing only the S1 regions of PK-A and PK-B (i.e., HP-A and HP-B, respectively, Fig. 1) for control experiments.

CD spectra of PK (PK-A, PK-B) and HP (HP-A, HP-B) (20 μM) at 0 °C, which showed ellipticity patterns with a positive peak at 265 nm and two negative peaks at 237 and 210 nm (Fig. 2A and Fig. S2, ESI†), corresponding to the signature of an RNA duplex.²⁵ As the PK structure is a combination of two or more hairpin structures, CD spectra of PKs are usually similar to the RNA duplex structure, as observed previously.²⁶ Differential spectra of PK from those of the corresponding HP showed a positive peak at 270 nm and negative peaks at 248 and 210 nm (Fig. S3, ESI†). As the structural difference of PK and HP is the formation of the S2 region, the differential spectra were expected to show the spectra of the S2 region in the PKs.

Fig. 2 (A) CD spectra of PK-A-5a (red) and HP-A (green) at 0 °C in a buffer containing 10 mM Na₂HPO₄ (pH 7), 100 mM NaCl, and 1 mM Na₂EDTA. (B) CD melting profile of PK-A-5a after subtraction of the HP-A signal. In (A) and (B), oligos of 20 μM were analysed.

Previous studies have indicated that the two stem regions, S1 and S2, show different melting temperatures.²³ For MMTV VPK, the S2 region is less stable than the S1 region and shows a lower melting temperature.^19a After the melting of the S2 region, the PK can be considered as a hairpin structure with a long 3′ tail in which the structure and subsequent melting behaviour of the remaining S1 region should be similar to that of the control HP (Scheme S1, ESI†). Here, it is important to note that the HP stability mainly depends on the loop length and nucleobase compositions close to the stem, which is the same for the HP structure formed after S2 dissociation in the designed PKs and control HP.²⁷ Thus, we envisioned that the subtraction of the melting profiles of the control HP from those of PK enables the extraction of the melting profile of the S2 region. The melting profiles of the PKs and HPs were collected at 210 nm as the differential spectra at 0 and 90 °C showed a maximum difference at 210 nm (Fig. S3, ESI†). Representative CD melting profiles of PK-A-5a, HP-A, PK-B-2a, and HP-B are shown in Fig. S4 (ESI†). Following this hypothesis, the CD melting profiles of PK showed two transitions, although the second transition, which corresponds to the melting of the S1 region, did not reach saturation in the experimental range. The melting transition of HP did not reach saturation as well as the second transition of PK. Thus, the first melting transition of PK around 50 °C corresponds to the melting of the S2 region. Fig. 2B shows the typical CD melting profile of the S2 region of PK-A-5a after subtracting the signal obtained from HP-A. The melting profiles of the PKs after subtracting the corresponding HP signal showed similar two-state transitions in both the PK-A and PK-B series, enabling calculation of the thermodynamic parameters using eqn (1), in the ESI,† after normalising the signals between 0–80 °C (Fig. S5, ESI†).

Table 1 shows the thermodynamic parameters ΔH°, TΔS°, and ΔG° at 37 °C for the S2 regions of all PKs from the PK-A and PK-B series and the parental MMTV VPK. The melting temperature (T_m) was calculated using the ΔH° and ΔS° values. The calculated parameters of the MMTV VPK were consistent with those from the previous studies, where the parameters were evaluated at 50 and 200 mM salt concentrations with different buffer reagents,^18a,19a,20 highlighting the accuracy of our analysis. Importantly, PKs in the same set, named “a” and “b”, having identical NN base pairs in the S2 region, exhibited similar thermodynamic parameters in both PK-A and PK-B (Table 1). For example, differences in between PK-A-5a and PK-A-5b and between PK-B-2a and PK-B-2b were 0.4 and 0.0 kcal mol⁻¹, respectively. The average between the PKs with the same NN base pairs in the PK-A series was 0.17 ± 0.13 kcal mol⁻¹, and that in the PK-B series was 0.03 ± 0.05 kcal mol⁻¹. In contrast, the maximum differences in in different sets in the PK-A (between PK-A-1b and PK-A-7b) and PK-B (between PK-B-4a or PK-B-4b and PK-B-2a or PK-B-2b) series were 2.5 and 0.7 kcal mol⁻¹, respectively, which are significantly larger than the average between the PKs with the same NN base pairs. The observation that the set of S2 stem regions with identical NN base pairs showed similar thermodynamic parameters, whereas the differences in parameters between different sets were relatively large, suggested the validity of the NN model in the PK structure.

Table 1 Thermodynamic parameters for the formation of S2 in the various PK-As and PK-Bs were analysed using CD melting at 210 nm

RNA sequences	ΔH° ^a (kcal mol^−1)	TΔS° ^a (kcal mol^−1)	^a (kcal mol^−1)	T m (°C)
a Values and errors are the average ± SD of at least five or more experiments. b Averaged values have a slight deviation from a value calculated from ΔH° and TΔS° values in the left columns due to the rounding off process of each replicated data.
MMTV VPK	−32.9 ± 6.6	−31.6 ± 6.3	−1.4 ± 0.3^b	50.4 ± 1.3
PK-A-1a	−44.7 ± 5.1	−45.1 ± 5.2	0.4 ± 0.2	34.5 ± 1.5
PK-A-1b	−37.9 ± 4.9	−38.6 ± 5.0	0.6 ± 0.2^b	31.9 ± 1.3
PK-A-2a	−39.0 ± 3.6	−38.3 ± 3.6	−0.7 ± 0.1	42.4 ± 1.0
PK-A-2b	−40.0 ± 4.1	−39.5 ± 4.1	−0.5 ± 0.1	40.6 ± 0.9
PK-A-3a	−42.2 ± 2.4	−42.3 ± 2.4	0.1 ± 0.0	36.3 ± 0.1
PK-A-3b	−35.8 ± 5.2	−35.9 ± 5.2	0.1 ± 0.2	36.1 ± 1.8
PK-A-4a	−38.2 ± 2.9	−38.5 ± 2.9	0.3 ± 0.1	34.6 ± 1.0
PK-A-4b	−42.2 ± 3.2	−42.4 ± 3.2	0.2 ± 0.1	35.9 ± 0.4
PK-A-5a	−38.2 ± 3.5	−36.8 ± 3.5	−1.4 ± 0.1	48.6 ± 0.9
PK-A-5b	−30.2 ± 2.7	−29.2 ± 2.7	−1.0 ± 0.1	47.4 ± 1.8
PK-A-6a	−36.3 ± 0.9	−35.5 ± 0.9	−0.8 ± 0.1	44.2 ± 1.1
PK-A-6b	−38.5 ± 1.6	−37.8 ± 1.5	−0.7 ± 0.1	43.0 ± 1.0
PK-A-7a	−29.4 ± 5.4	−27.7 ± 5.2	−1.7 ± 0.3	56.7 ± 2.5
PK-A-7b	−29.3 ± 6.6	−27.4 ± 6.3	−1.9 ± 0.3	58.4 ± 2.5
PK-B-1a	−43.3 ± 5.4	−43.2 ± 5.4	−0.1 ± 0.2	37.8 ± 1.6
PK-B-1b	−43.2 ± 5.7	−43.1 ± 5.6	−0.1 ± 0.1	37.8 ± 1.0
PK-B-2a	−43.9 ± 4.4	−43.3 ± 4.3	−0.6 ± 0.2	41.6 ± 1.0
PK-B-2b	−38.6 ± 6.2	−38.0 ± 6.0	−0.6 ± 0.3	42.0 ± 2.2
PK-B-3a	−35.6 ± 4.8	−35.6 ± 4.7	0.0 ± 0.1	37.4 ± 0.8
PK-B-3b	−30.5 ± 6.1	−30.4 ± 6.2	−0.1 ± 0.2	38.1 ± 2.7
PK-B-4a	−40.1 ± 7.7	−40.2 ± 8.1	0.1 ± 0.6	37.0 ± 4.1
PK-B-4b	−41.2 ± 4.9	−41.3 ± 4.8	0.1 ± 0.2	36.3 ± 1.4

---

UV melting profiles were also measured at 260 nm in the same buffer as that used in the CD melting analyses with an oligonucleotide concentration of 2 μM. The UV melting profiles of the PKs consisted of two melting transitions where the second transition was not saturated (Fig. S6, ESI†). Unsaturated UV melting was also observed for HP-A and HP-B (Fig. S6, ESI†). Thus, the corresponding HP was subtracted from the PK to calculate the thermodynamic parameters as well as the CD melting (Fig. S7, ESI†). The thermodynamic parameters calculated from the UV melting profiles are tabulated in Table S2 (ESI†). A set of PKs having identical NN base pairs in the S2 stem region resulted in similar thermodynamic parameters for both the PK-A and PK-B series. Similar to the results in Table 1, the average values of PKs with identical NN base pairs were 0.37 ± 0.30 and 0.35 ± 0.40 kcal mol⁻¹ for the PK-A and PK-B series, respectively, which are both significantly smaller than the maximum difference between PKs with different NN base pairs. The maximum difference for the PK-A series was 3.5 kcal mol⁻¹ between PK-A-7b and PK-A-1b, and that for the PK-B series was 1.3 kcal mol⁻¹ between PK-B-2a and PK-B-1a. Recently, we have reported a similar magnitude of average free energy difference of 0.2 ± 0.2 kcal mol⁻¹ for DNA duplexes having identical NN base pairs in molecular crowding conditions.¹¹ It is important to note that the concentration of the oligonucleotide used for UV melting analyses (2.0 μM) is 10-fold less than that used for CD melting analyses. The calculated thermodynamic parameters for all the PKs through CD (20.0 μM) and UV (2 μM) melting were closely related with an average difference of 0.27 ± 0.24 kcal mol⁻¹, indicating the intramolecular nature of the calculated parameters (Table 1 and Table S2, ESI†). These results also highlight the validity of the NN model in the PK structure.

The validity of the NN model was verified in both the PK-A and PK-B series, which had different S1 stem lengths. Both PK-A and PK-B, with identical NN base pairs of the S2 stem, exhibited similar thermodynamic parameters. The average between PKs with the same S2 and different S1 (e.g., PK-A-1a/PK-B-1a), as obtained from CD melting, was 0.25 ± 0.23 kcal mol⁻¹. The similar between PKs with the same S2 stem suggests that the stability of the S2 stem is less affected by the length of the S1 stem when the length of the L3 region is sufficient to form the S2 stem.

The of the S2 stem region, when accounting for the formation of an RNA duplex, was predicted by the established NN model under 100 mM NaCl conditions (Table S3 and Scheme S2, ESI†).^4f,28 The stability of the S2 stem in PK was relatively lower than that of the same sequence forming a canonical duplex. Based on the tertiary structure,^18b it is considered that the L1 region consisting of two nucleotides is too short to keep the continuous helicity of S1 and S2 through the coaxial stacking, which generally provides additional stabilisation energy to the stacked stems. In addition, due to the presence of L2 nucleotide, which functions as a hinge between the S1 and S2 stems, the two stems bend resulting in reduced stability of the S2 stem (Fig. S8, ESI†). In contrast to the overall reduced stability of the S2 stem, the relative stability between NN base pairs changed following the predicted stability based on the NN model. The values in different PKs of the PK-A series were plotted with the predicted values of the corresponding S2 duplexes (Fig. S9, ESI†). The linear correlation (y = 0.99x + 5.97), in which the slope value is close to 1, suggests that the internal physicochemical properties of the S2 stem region are similar to those of the canonical RNA duplex, highlighting the central role of NN base pairs in determining the S2 stem stability in PK. In contrast, the intercept value of the linear correlation, shown in Fig. S9 (ESI†), suggests that the extra factors destabilised the S2 stems by approximately 6 kcal mol⁻¹ as compared to the canonical duplexes. The destabilisation is attributed to the different bulk interactions of the S2 stem. The nucleobases in the L1 region are placed at the major groove side of the S2 stem, orienting away from the groove,^18b suggesting the potential of L1 in reducing hydration around the S2 stem. Reduced hydration drastically destabilises the duplex, as previously observed with modified nucleobases at the major groove side.²⁹ Reduced hydration caused by only two modified nucleobases destabilised the DNA duplex by 4–5 kcal mol⁻¹ at 20 °C. In addition to the bulk interaction, the overall S2 stem stability in PK is also affected by some additional factors such as bending of the S2 stem, inefficient coaxial stacking etc. from the complex tertiary structure. The influences of extra factors such as length and sequence in the loop region can be systematically analysed by varying the factors in different PK structures, in which their tertiary structures are available.

In summary, the applicability of the NN model in complex PK structures was investigated by monitoring S2 formation, which separates a PK from a conventional hairpin structure. Identical NN base pairs in the S2 stem of different PKs exhibited similar thermodynamic parameters, highlighting the validity of the NN model. Predicting the stability of PK stems has crucial implications in various biological processes. For instance, programmed ribosomal frameshifting (PRF) caused by a PK structure produces different proteins at a specific ratio.³⁰ The PRFs related to PK stability are receiving more attention for their importance in the production of viral proteins due to the coronavirus that caused the COVID-19 pandemic, SARS-CoV-2.³¹ Thus, predicting the thermodynamic stability of PK remains a key aspect to investigate the PRF processes.³² The thermodynamic parameters determined in this study will help demonstrate the relationship between PK stability and its contribution to biological processes.

This work was supported by Grants-in-Aid for Scientific Research from the Ministry of Education, Culture, Sports, Science and Technology (MEXT) and Japan Society for the Promotion of Science (JSPS) (KAKENHI Grant No. 17H06351, 18KK0164, 19H00928, and 21H05108), The Hirao Taro Foundation of KONAN GAKUEN for Academic Research, and Chubei Itoh Foundation.

Conflicts of interest

There are no conflicts of interest.

Notes and references

1. (a) P. Belmont, et al. , Chem. Soc. Rev., 2001, 30, 70; (b) P. H. von Hippel and E. Delagoutte, Cell, 2001, 104, 177; (c) N. Sugimoto, Chemistry and Biology of Non-canonical Nucleic Acids, Wiley-VCH, Weinheim, Germany, 2021.
2. (a) S. Takahashi and N. Sugimoto, Chem. Soc. Rev., 2020, 49, 8439; (b) S. Takahashi and N. Sugimoto, Acc. Chem. Res., 2021, 54, 2110.
3. (a) I. Tinoco, et al. , Nature, 1971, 230, 362; (b) I. Tinoco, et al. , Nature, New Biol., 1973, 246, 40.
4. (a) S. M. Freier, et al. , Proc. Natl. Acad. Sci. U. S. A., 1986, 83, 9373; (b) N. Sugimoto, et al. , Biochemistry, 1995, 34, 11211; (c) J. SantaLucia, et al. , Biochemistry, 1996, 35, 3555; (d) N. Sugimoto, et al. , Nucleic Acids Res., 1996, 24, 4501; (e) J. SantaLucia, Jr., Proc. Natl. Acad. Sci. U. S. A., 1998, 95, 1460; (f) T. Xia, et al. , Biochemistry, 1998, 37, 14719.
5. (a) J. M. Huguet, et al. , Proc. Natl. Acad. Sci. U. S. A., 2010, 107, 15431; (b) S. Ghosh, et al. , Proc. Natl. Acad. Sci. U. S. A., 2020, 117, 14194.
6. I. Ferreira, et al. , Chem. Phys., 2019, 521, 69.
7. (a) D. Banerjee, et al. , Nucleic Acids Res., 2020, 48, 12042; (b) D. Banerjee, et al. , Nucleic Acids Res., 2021, 49, 10796.
8. D. H. Mathews, et al. , J. Mol. Biol., 1999, 288, 911.
9. (a) N. Peyret, et al. , Biochemistry, 1999, 38, 3468; (b) L. He, et al. , Biochemistry, 1991, 30, 11124.
10. R. W. Roberts and D. M. Crothers, Proc. Natl. Acad. Sci. U. S. A., 1996, 93, 4320.
11. S. Ghosh, et al. , Nucleic Acids Res., 2019, 47, 3284.
12. M. S. Adams and B. M. Znosko, Nucleic Acids Res., 2019, 47, 3658.
13. D. D. Licatalosi and R. B. Darnell, Nat. Rev. Genet., 2010, 11, 75.
14. (a) A. Peselis and A. Serganov, Wiley Interdiscip. Rev.: RNA, 2014, 5, 803; (b) I. Brierley, et al. , Nat. Rev. Microbiol., 2007, 5, 598; (c) D. W. Staple and S. E. Butcher, PLoS Biol., 2005, 3, e213.
15. (a) A. Roth and R. R. Breaker, Annu. Rev. Biochem., 2009, 78, 305; (b) C. E. Hajdin, et al. , Proc. Natl. Acad. Sci. U. S. A., 2013, 110, 5498.
16. T. M. Hansen, et al. , Proc. Natl. Acad. Sci. U. S. A., 2007, 104, 5830.
17. (a) M. Kuwahara, et al. , Nucleic Acids Res., 2006, 34, 5383; (b) S. Cao and S.-J. Chen, Nucleic Acids Res., 2006, 34, 2634.
18. (a) R. L. Gonzalez and I. Tinoco, J. Mol. Biol., 1999, 289, 1267; (b) L. X. Shen and J. I. Tinoco, J. Mol. Biol., 1995, 247, 963.
19. (a) C. A. Theimer and D. P. Giedroc, RNA, 2000, 6, 409; (b) C. Reiling, et al. , J. Phys. Chem. B, 2015, 119, 1939.
20. R. Narayanan, et al. , J. Am. Chem. Soc., 2011, 133, 18767.
21. G. Chen, et al. , RNA, 2007, 13, 2175.
22. (a) C. W. A. Pleij, Curr. Opin. Struct. Biol., 1994, 4, 337; (b) F. H. D. van Batenburg, et al. , Nucleic Acids Res., 2001, 29, 194.
23. D. P. Giedroc, et al. , J. Mol. Biol., 2000, 298, 167.
24. X. Chen, et al. , EMBO J., 1995, 14, 842.
25. M. Sarkar, et al. , Biochemistry, 1996, 35, 4678.
26. K. H. Johnson and D. M. Gray, J. Biomol. Struct. Dyn., 1992, 9, 733.
27. M. J. Serra, et al. , Nucleic Acids Res., 1993, 21, 3845.
28. S.-i. Nakano, et al. , Nucleic Acids Res., 1999, 27, 2957.
29. M. Ganguly, et al. , J. Am. Chem. Soc., 2009, 131, 12068.
30. T. Endoh and N. Sugimoto, Chem. Rec., 2017, 17, 817.
31. J. A. Kelly, et al. , J. Biol. Chem., 2020, 295, 10741.
32. (a) K. Neupane, et al. , Nat. Commun., 2021, 12, 4749; (b) S. I. Omar, et al. , PLoS Comput. Biol., 2021, 17, e1008603.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: 10.1039/d1cc07094k

---