Mapping the Chemical Space of the RNA Cleavage and Its

Subscriber access provided by READING UNIV

Article

Mapping the Chemical Space of the RNA Cleavage and its Implications for Ribozyme Catalysis Vojtech Mlynsky, Petra Kührová, Petr Jurecka, Jiri Sponer, Michal Otyepka, and Pavel Banáš J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.7b09129 • Publication Date (Web): 08 Nov 2017 Downloaded from http://pubs.acs.org on November 9, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry B is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Mapping the Chemical Space of the RNA Cleavage and its Implications for Ribozyme Catalysis Vojtěch Mlýnský,1,2,# Petra Kührová,1,# Petr Jurečka,1 Jiří Šponer,1,3 Michal Otyepka,1* and Pavel Banáš,1,3* 1

Regional Centre of Advanced Technologies and Materials, Department of Physical Chemistry, Faculty of Science, Palacký University, tř. 17 listopadu 12, 771 46, Olomouc, Czech Republic 2

3

Scuola Internazionale Superiore di Studi Avanzati (SISSA), via Bonomea 265, 34136 Trieste, Italy

Institute of Biophysics of the Czech Academy of Sciences, Kralovopolská 135, 612 65 Brno, Czech Republic *

corresponding authors: [email protected], [email protected] #

both authors contributed equally

1 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract Ribozymes utilize diverse catalytic strategies. We report systematic quantum chemical calculations mapping the catalytic space of RNA cleavage by comparing all chemically feasible reaction mechanisms of RNA self-cleavage, using appropriate model systems including those chemical groups that may directly participate in ribozyme catalysis. We calculated the kinetics of uncatalyzed cleavage reactions proceeding via both monoanionic and dianionic pathways, and explicitly probed effects of various groups acting as general bases (GBs) and/or general acids (GAs), or electrostatic transition state stabilizers. In total, we explored 115 different mechanisms. The dianionic scenarios are generally preferred to monoanionic scenarios, although they may compete with one-another under some conditions. Direct GA catalysis seems to exert the dominant catalytic effect, while GB catalysis and electrostatic stabilization are less efficient. Our results indirectly suggest that the dominant part of the catalytic effect might be explained by the shift of the reaction mechanism from the mechanism of uncatalyzed cleavage to the mechanism occurring in ribozymes. This would contrast typical protein enzymes, primarily achieving catalysis by overall electrostatic effects in their catalytic center.

Introduction Ribozymes are nucleic acid counterparts of protein enzymes.1-2 While many features of RNA enzymatic activity are well understood,3-7 the general principles of RNA catalysis are still debated.8-13 Small self-cleaving ribozymes are particularly suitable model systems for studying the catalytic mechanisms of ribozymes because of their small size and the growing body of relevant structural and mechanistic data.14-19 They catalyze phosphodiester bond cleavage (and ligation) by an SN2-type reaction that is initiated by the nucleophilic attack of the O2´ oxygen on the adjacent scissile phosphate. The reaction proceeds via a pentacoordinate phosphorane transition state, generating a 2´,3´-cyclic phosphate and a 5´-OH terminated RNA fragment.9, 20 An enzyme’s catalytic effect can be defined as the difference between the reaction rates (or difference between the corresponding free-energy barriers, ∆G‡) of the catalyzed reaction and the uncatalyzed (reference) reaction in water.7, 21 However, interpretation of the overall catalytic effect is complicated by the fact that the reaction mechanisms (reaction pathways) in enzymes might differ from the mechanism of the corresponding reference uncatalyzed reaction in water. Therefore, Warshel et al. suggested to divide the catalytic effect into two general components.7 First, there are explicit contributions due to a change in the reaction mechanism. This part of the catalytic effect reflects the difference in the reaction barrier between the mechanism of the uncatalyzed reference reaction and the actual mechanism utilized by enzymes that directly involve active site groups; both considered in water environment. This explicit component thus corresponds to direct participation of (external) groups such as nucleobases and/or hydrated Mg2+ ions in the chemical reaction. The second effect stems from additional implicit influence of the enzyme environment,7, 21-24 which is associated with electrostatic contributions of the enzyme environment including electrostatic pre-organization of the active site and reduction of its polarity. It reflects the role of the surroundings of the active site not directly participating in the reaction pathway that provides transition state (TS) stabilization or ground state destabilization. Thus, to decipher a catalytic effect and in particular to separate the role of the explicit and implicit components we need to characterize, in water environment, the reference reaction and all 2 ACS Paragon Plus Environment

Page 2 of 31

Page 3 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


possible reaction pathways that include directly the proximal chemical groups. In protein enzymes, the dominant part of the catalytic effect is typically due to the implicit effects.7, 21, 25 As we show here, RNA catalysis appears to be a different story. Taking into account both above-defined components of the catalytic enhancement, four potentially relevant strategies for ribozyme catalysis have been suggested:20, 26-27 (i) preorganization of the active site in the in-line attack conformation, where the angle between the O2´ atom and the P-O5´ unit of the scissile phosphate approaches 180°; (ii) activation of the 2´-OH nucleophile by a general base (GB); (iii) TS stabilization by neutralizing the negative charge developed on the nonbridging oxygens (nbO) of the scissile phosphate during the cleavage reaction; and (iv) protonation of the leaving O5´ group by a general acid (GA). Uncatalyzed RNA backbone cleavage has been studied computationally and experimentally using small models in the past. Lönnenberg and co-workers suggested that the rate-limiting step under GA/GB catalysis is the exocyclic cleavage step, which occurs in concert with proton transfer to the departing alkoxy group.28-30 Perreault and Anslyn concluded that there are two possible and interchangeable mechanisms of RNA cleavage (known as the dianionic and monoanionic reaction paths) that differ with respect to the protonation state of the phosphorane intermediate (IN).31 Quantum mechanical (QM) calculations can provide unique atomistic insights into biocatalysis,32-35 enabling rationalization of catalytic effects. Numerous computational studies on the internal transesterification reaction have been conducted, but most of them used the simplest possible model of phosphorane systems, and only a few included the potentially crucial effects of the ribose ring.36-41 Charged nucleobases and hydrated Mg2+ ions have been identified as potential GB and GA in the self-cleavage mechanism of known ribozymes.15, 36, 42-52 In particular, in case of hydrated Mg2+ ions, the hydroxyl ion and water molecules coordinated to such Mg2+ ion can directly participate during the catalysis as GB and GA, respectively, as shown for the phosphodiester cleavage of small RNA models,37, 53 as well as RNA backbone cleavage by proteins54-55 and ribozymes.40, 44, 56-60 Further, the ions might also play a role of Lewis acid stabilizing negative charge developed on the scissile phosphate during cleavage reaction.44, 58-60 In the effort to probe explicit effects of the active site groups on the phosphodiester cleavage reaction and to map all possible reaction mechanisms of the sugar-phosphate backbone cleavage, we used high-level QM calculations to follow self-cleavage of a wide range of RNA cleavage models. We studied the reference (minimal unassisted) scenarios as well as variety of reaction paths containing additional (charged) nucleobases and hydrated Mg2+ ions, following reaction mechanisms anticipated in ribozymes. In this way we explicitly probed the direct component of the ribozyme’s catalytic effect as described above. The effects of external groups were divided according to the direct participation in proton transfers and/or their indirect stabilization of TS via the electrostatic potential during the cleavage reaction. We also probed the performance of the implicit solvent methods for description of the free-energy balance between various protonation states of the nucleobases and we suggested a careful adjustment of the implicit solvent model settings. In addition, for some reaction mechanisms we explicitly compared the description of the reaction profiles obtained by implicit solvation with those using solvation by explicit water molecules. Taking together, our finding supports the notion that the GA plays the dominant role in the catalytic effect whereas the presence of GB catalysis leads to subtle enhancements of the cleavage rates. Among all analyzed reaction paths, the direct involvement in proton transfers of the external groups in both GA and GB positions showed the largest rate enhancements, which were almost comparable to the catalytic effects of naturally occurring ribozymes. This might suggest that in contrast to protein enzymes a dominant part of the ribozyme’s catalytic effect can 3 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

be explained solely by shift of the reaction mechanism due to explicit involvement of GA/GB groups.

Methods The reference self-cleavage reaction was investigated as the cleavage of the 3´-(1´-amino-4´methylribose)-5´-methylphosphodiester. That minimal model containing 27 atoms was later extended by particular combinations of methyl-capped nucleobases and hexa-coordinated Mg2+ ions, namely the N9-methylguanine (Gua), N1-deprotonated N9-methylguanine (Gua−), N9methyladenine (Ade), N1-protonated N9-methyladenine (AdeH+), N9-methylcytosine (Cyt), N3protonated N9-methylcytosine (CytH+), and hexa-coordinated Mg2+ ions with both inner-shell ([Mg(H2O)5]2+ or [Mg(H2O)4OH]+) and outer-shell ([Mg(H2O)6]2+ or [Mg(H2O)5OH]+) coordination to the nbO atom of the scissile phosphate (Table 1). For each particular reaction mechanism we modeled geometry of TS or IN state that was subsequently used as an initial structure for following the reaction coordinate toward R and P states. All stationary points along the reaction path were subsequently optimized (to minima or first order saddle point transition state) and their free energies were calculated using Hessian matrix and harmonic oscillator and rigid rotor approximation (see below). The initial geometries used to probe the minimal monoanionic model (here labeled as unassisted monoanionic scenario V, Figure 1) were taken from previous works, where QM calculations of reference monoanionic reaction on small model systems were performed.39-40 The initial geometry of the TS of dianionic reaction path (here labeled as unassisted dianionic scenario I) was prepared by removing the hydrogen of the 2´-OH group from the TS geometry of monoanionic path and was subsequently optimized to transition state. Starting TS or IN geometries of reaction scenarios with extended models were prepared from IN or TS structures of uncatalyzed scenarios by placing external groups in positions allowing H-bond interactions with O2´, O5´ and nbO of the scissile phosphate. In particular, we modeled positions of suspected GB and GA to form H-bonds with O2´ and O5´, respectively and simultaneously forming the second H-bond contact with one of the nbOs of the scissile phosphate. Again, all these newly created models were fully relaxed and subsequently a detailed search for R, P, TS and IN states was performed by scanning the reaction coordinate. The exocyclic cleavage step corresponding to the final part of the reaction path (i.e., from the initial TS/IN toward P state) was scanned by elongating of the P…O5´ distance. Similarly, the nucleophile attack was scanned in reverse direction (i.e., from the initial TS/IN to R state) by elongating of the O2´…P distance. Finally, the rearrangement between two IN states of the monoanionic paths (scenarios V and VI) corresponding to the rotation of the phosphorane hydroxyl group toward the leaving O5´ oxygen was followed by scanning the O(pro-RP)-PO(pro-SP)-H and O(pro-SP)-P-O(pro-RP)-H dihedrals in endo and exo paths, respectively. The scans were performed in 0.1 Å step (or 10° in case of scanning the hydroxyl reorientation) and all remaining degrees of freedom including the ribose pucker were fully relaxed at each point. Note that none of the reaction models used in this study to probe the cleavage mechanisms are specific for any particular ribozyme and are entirely independent of any ribozyme structure (see section Classification of reactive states from catalytic scenarios of Results and Discussion for more details). All structures were optimized in implicit solvent model (conductor-like polarizable continuum model CPCM) by hybrid DFT functional developed for kinetics (modified Perdew-Wang 1parameter for kinetics MPW1K)61-62 with 6-311+G(2df,2p) triple-zeta basis set. For all localized 4 ACS Paragon Plus Environment

Page 4 of 31

Page 5 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


minima and TS states we subsequently calculated free energies using harmonic oscillator and rigid rotor approximation. Hessian matrix and corresponding frequencies were calculated at the same level for each optimized structure to estimate free energies of the corresponding states at 300 K and 1 atm. The differences between the gas-phase MPW1K/6-311+G(2df,2p) energies and the CPCM (εr = 78.4)/MPW1K/6-311+G(2df,2p) energies were used to estimate the solvation contributions to the energy profile of the reference reaction. The MPW1K functional was used because our previous studies have shown that the MPW1K/6-31+G(d,p) (gas-phase) error against CBS(T) method (MP263 energies extrapolated for the complete basis set (CBS) and corrected for the higher order correlation effects by CCSD(T) in a small basis set64-65) for the unassisted monoanionic scenario V (Figure 1) is ~1 kcal/mol.40-41 Thus, we selected the hybrid DFT MPW1K functional also for the present study. For the sake of completeness, we calculated the reference CBS(T) energies for the other unassisted path, i.e., the dianionic scenario I (see Accuracy of the MPW1K functional in the Supporting Information for details), in order to evaluate MPW1K gas-phase error also for this particular mechanism. For the purpose of this work, we have modified standard default setting for building up the cavity in the CPCM solvent model of Gaussian 09.66 We found that none of the standard settings of CPCM was able to correctly reproduce the relative balances of the (de)protonations of different ionizable groups due to unbalanced representation of the solvation term (see Supporting Information). Therefore, the standard UFF radii model (sphere around each solute atom multiplied by the electrostatic scaling factor of 1.1)66 was replaced by Bondi radii model66-67 scaled by a factor of 0.9 (see Tests and Adjustments of the Computational Protocol in Supporting Information for details). Beside implicit solvent calculations, we probed the uncatalyzed mechanisms also in explicit solvent. For this purpose, the structures of minimal models for dianionic and monoanionic mechanisms (of reaction scenarios I and V) were solvated by a shell of 53 explicit water molecules ensuring a single continuous layer of water molecules surrounding the solute atoms. Initial positions of explicit water molecules were generated by the Leap module of AMBER1668 using three R states geometries (one from dianionic scenario I and two from monoanionic scenarios V). These solvated systems were further optimized in Gaussian09 with the reparameterized implicit solvent (CPCM/Bondi radii scaled by factor 0.9), which modeled solvation outside the water droplet. Subsequently, we replaced geometry of the solute atoms by corresponding structures of TS, IN and P states from the particular mechanism and each state was further re-optimized. We note that the calculation of Hessian matrix required for TS optimization and calculations of the free-energy corrections was not feasible due to excessive computer demands, as the explicitly hydrated systems contained 188 and 189 QM atoms for dianionic scenario I and monoanionic scenarios V, respectively. Thus, the TS states were localized by constrained optimizations, in which we fixed the O2´-P and P-O5´ distances at values obtained in the corresponding implicit solvent calculations. In addition, the free-energy corrections for all states along the reaction coordinate from a particular reaction scenario were extrapolated from the corresponding implicit solvent calculation. The major protonation states at physiological conditions (pH ~7, T = 300 K) are ribose moiety with 2´-OH, canonical forms of nucleobases and Mg2+ hexahydrate. Since we also considered the model with deprotonated 2´-O− ribose moiety and additional extended models having charged methyl-capped nucleobases and/or OH− groups coordinated to Mg2+ ions, the calculated ∆G energies of the R states (and thus all energies along corresponding reaction paths) had to be corrected for the presence of only a minor equilibrium population of the assumed catalytically potent ionization form of the catalytic species. This has been done using empirical corrections for 5 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

the (de)protonation of particular ribose/nucleobase/complex Mg2+ ion moieties which are summarized in Table S2 in Supporting Information. We note that the pKa value of the 2´-OH group from phosphorylated ribose moieties is reported ambiguously in literature. Different approaches estimated pKa values between 12.3 and 14.9.69-76 We opted for the pKa of 13.7, which is in the middle of the above mentioned range and was measured at T=298 K and the physiological ionic concentration (Table S2 in Supporting Information).69

Results and Discussion The purpose of accurate QM calculation on small models is to evaluate the role of the shift of the chemical reaction Here we used advanced QM methods to compute free-energy activation barriers (∆G‡300K) related to kinetic rate constants (see Supporting Information for details) for all possible mechanisms of intramolecular phosphodiester cleavage of the RNA model system in water. In order to calculate free-energy barriers, we localized and minimized all geometries corresponding to R, TS, IN and P states along the reaction coordinates and subsequently we computed their free energies using harmonic oscillator and rigid rotor approximation. By performing the calculations on small model systems using accurate QM methods, we were able to investigate the explicit effects of the active site groups directly participating on the reaction, i.e., (charged) nucleobases and hydrated Mg2+ ions. In other words, we fully mapped part of the ribozyme’s catalytic effect, namely, the effect of the shift of the reaction mechanism from the mechanism of uncatalyzed cleavage to the mechanisms anticipated in ribozymes, which explicitly includes involvement of the active site groups. We would like to point out that probing the contributions standing beyond the local effects, such as electrostatic contribution of the active site preorganization (except of the effect of active site groups actively participating on the mechanism such as GA/GB groups), conformational strain, specific desolvation etc., is beyond the scope of this study. In order to include these nonlocal effects, a description of full ribozyme by QM/MM techniques would be required, as commonly done in the past in order to describe a reaction mechanism of various ribozymes. The QM/MM approach typically limits the study to probing only a few suggested reaction mechanisms and the results may depend on the reliability of the assumed reactant and product structures; their preparation in the context of the complete ribozyme is a challenging issue, due to the ruggedness of the conformational space of the ribozymes.77-79 In addition, some of the QM/MM studies may require description of the catalytic pocket by computationally less expensive semi-empirical methods in order to guarantee sufficient sampling. Unfortunately, the semi-empirical methods may compromise the quality of the QM description.80-82 We suggest that the present accurate QM computations on small models complement the experimental studies and QM/MM investigations of complete ribozymes. They aim to separate the catalytic effect caused by change of the reaction pathway by the ribozyme from the other contributions of the ribozyme environment. Studied mechanisms In the present computations, we probed unassisted mechanisms as well as the mechanisms, in which the nucleobases and hydrated Mg2+ ions either participated directly in proton transfers by functioning as a GA or GB, or assisted the reaction indirectly via electrostatic stabilization of the developing negative charge on the scissile phosphate during the cleavage reaction (Table 1). In total, we studied seven reaction scenarios (Figure 1). Four scenarios (I-IV) represented dianionic 6 ACS Paragon Plus Environment

Page 6 of 31

Page 7 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


mechanisms in which cleavage proceeded via a double-negatively charged deprotonated phosphorane. The remaining three scenarios (V-VII) corresponded to monoanionic mechanisms involving a negatively charged singly-protonated phosphorane (see Figure 1 and Table 1). Reaction scenarios II to IV, VI, and VII encompass several mechanisms with different moieties serving as GB or GA. Detailed definitions of the reference state of the reactants and free-energy adjustments (∆Gcorr) for the (de)protonation of specific ribose/nucleobase/Mg2+ ion moieties are summarized in Supporting Information. In total, we explored 115 different pathways. In the main text we discuss only overall rate-determining barrier of each mechanisms; for complete data including free energies of all states along the reaction coordinates for all studied mechanisms please see Tables S6-S8 in Supporting Information.

Figure 1: All the possible reaction scenarios for the internal transesterification. Reaction scenarios shown on yellow and white backgrounds correspond to dianionic and monoanionic mechanisms, respectively. Rare protonation states of phosphates and phosphoranes at pH ~7 (shown on a grey background) are assumed not to contribute to the catalysis due to enough high free-energy penalty associated with their rare protonation state that is not compensated by enhanced reactivity, and thus they were not considered in this study. In total, 115 distinct pathways were calculated. Table 1: Reaction scenarios probed in this study. The numbering of each scenario matches that used in Figure 1. Each scenario is characterized in terms of the chemical groups that are directly involved in proton transfer as GAs and/or GBs, or as electrostatically stabilizers at position of general base (es-GB) or general acid (es-GA). The calculated pathway depict idealized chemical



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 31

reactions of the used model systems immersed into water and are thus not related to specific ribozymes. Scenario I

II

Mechanism (number of paths) dianionic (1)

dianionic (4)

GB

GA

es-GB

es-GA

-

-

-

-

-

dianionic (4)

Gua−, Cyt, [Mg(H2O)4OH]+, [Mg(H2O)5OH]+

IV

dianionic (16)

Gua−, Cyt, [Mg(H2O)4OH]+, [Mg(H2O)5OH]+

V

monoanionic (2)

III

VI

VII

monoanionic (48)

monoanionic (40)

AdeH+, CytH+, [Mg(H2O)5]2+, [Mg(H2O)6]2+ AdeH+, CytH+, [Mg(H2O)5]2+, [Mg(H2O)6]2+

nbO

nbO

nbO

nbO

AdeH+, CytH+, [Mg(H2O)5]2+, [Mg(H2O)6]2+

nbO

-

Gua−, Cyt, [Mg(H2O)4OH]+, [Mg(H2O)5OH]+ Gua−, Cyt, [Mg(H2O)4OH]+, [Mg(H2O)5OH]+ -

AdeH+, CytH+, [Mg(H2O)5]2+, [Mg(H2O)6]2+ AdeH+, CytH+, [Mg(H2O)5]2+, [Mg(H2O)6]2+ -

Gua, CytH+, [Mg(H2O)5]2+, [Mg(H2O)6]2+


Gua, CytH+, [Mg(H2O)5]2+, [Mg(H2O)6]2+


Classification of reactive states in catalytic scenarios We characterized the conformational variability of R states of the sugar-phosphate backbone with an attempt to identify the variability of catalytically potent backbone geometries. We used inhouse script to generate all possible conformations of the sugar-phosphate backbone using the minimal model of the sugar-phosphate backbone as used in the monoanionic scenario V. The potentially reactive geometries require high values (>170°) of in-line attack angle (IAA, O2´-PO5´) and simultaneously avoiding steric clashes.20 The IAA is fully determined by two backbone dihedral angles, ζ (C3´-O3´-P-O5´) and ε (C4´-C3´-O3´-P), and the sugar-pucker of the ribose ring. Therefore, we systematically scanned all these parameters in order to generate all possible conformers of the studied model. Subsequently, we searched all these geometries for suspected catalytically potent conformation satisfying high values of IAA and simultaneously avoiding intramolecular clashes. This structure generation algorithm kept all bond distances and bond angles constant; the geometry of the ribose ring was generated in such a way that only the bond distances were assumed to be constant while the four remaining degrees of freedom of the ribose ring were systematically scanned (see ref. 83 for the definition of the ribose ring configurations). Using this rigid geometry search, we identified four maxima of the IAA in the ε/ζ/pucker space with IAA of ~175-180° (see Figure S3 in Supporting Information). The geometries corresponding 8 ACS Paragon Plus Environment

Page 9 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


to the obtained maxima of IAA can be fully characterized by their ε/ζ/pucker values. Two of them contained C3´-endo pucker of the ribose and differed by ε/ζ dihedrals (ε/ζ of either ~140°/~140° or ~50°/~-140°). The other two contained C2´-endo pucker and their ε/ζ dihedrals equaled to either ~-170°/~140° or ~90°/~-140°. Both pairs of IAA maxima corresponding to the same pucker were separated in the ε/ζ space by structures having only slightly smaller IAA than the corresponding maximal values but still keeping IAA values above ~170° (see Figure S3 in Supporting Information). Note that these conformations were identified solely based on geometrical criteria, so they represent distinct islands in the backbone dihedral conformational space characterized by high IAA angle values (surrounded by conformations with lower IAA value), but do not necessarily correspond to four well-defined energy minima. In addition, the shape of IAA function in the ε/ζ/pucker space and in particular the number of IAA maxima and their localization in such space might be biased by the assumption of rigid bonds and bond angles used in the geometry generation algorithm. Therefore, we used these geometries as initial structures for subsequent series of constrained geometry optimizations, in which we constrained ε/ζ values, while the other degrees of freedom including the ribose pucker were fully relaxed. We thus obtained a QM energy and IAA as a function of ε/ζ dihedrals. We performed such ε/ζ scans for both C2´-endo and C3´-endo puckers. These geometry minimizations were performed at the BLYP/6-31G(d,p) level with the reparameterized CPCM solvent model. In contrast to the rigid bond-angles geometry generation algorithm, the QM-optimized structures revealed only two maxima of IAA in the ε/ζ/pucker space, one for each ribose pucker (with ε/ζ equaling to ~160°/~160° and ~100°/~180° for C2´-endo and C3´-endo pucker, respectively) with relatively broad distribution of high-IAA structures in the ε/ζ space (Figure 2). None of these IAA maxima corresponded to an energy minimum, but both of them were structurally and energetically accessible from the neighboring energy minima (Figure 2). In addition, as the barrier between C2´-endo and C3´-endo ribose pucker is ~4kcal/mol,84 we suggest that the transitions between catalytically potent geometries are energetically accessible and a final reactive geometry is mostly dictated by the local structural environment, i.e., presence of the external chemical groups around the scissile phosphate.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 31

Figure 2: In-line attack angle (A and C panels, in deg.) and potential energy calculated on BLYP/6-31G(d,p)/CPCM level of theory (B and D panels, in kcal/mol) as a function of ε/ζ dihedral angles for structures with C2´-endo and C3´-endo puckering. The particular ε/ζ and pucker of the reactive conformations of the selected ribozymes41, 59, 85 and unassisted reaction paths probed in this study (scenarios I and V) are depicted by black dots in the right panels. The R state structures of the unassisted reaction paths are depicted on panel F. The energy maxima highlighted by stars on panels B and D correspond to structures with steric clashes of O2´ oxygen with one of non-bridging phosphate oxygens. We observed that R states of all unassisted reaction scenarios (dianionic scenario I and monoanionic scenarios V) possessed C3´-endo puckering. Catalyzed dianionic paths involving either only GB or only GA catalysis (scenarios II and III, respectively) have R state geometries belonging to C3´-endo (75% of the studied pathways) and C2´-endo (25% of the studied pathways) puckering, while mechanisms involving both GB and GA (scenario IV) revealed that all R states are unanimously belonging to the C2´-endo region. Finally, the catalyzed monoanionic paths (scenario VI) displayed also R states of both puckering (~73% and ~26% for C2´-endo and C3´-endo, respectively). We note that R states from scenarios VI and VII are identical, i.e., sharing the same first part of the cleavage reaction (see Figure 1). In summary, the different mechanisms investigated in our study featured several different catalytically relevant conformations of the RNA sugar-phosphate backbone differing by ribose puckering and ε/ζ values, depending on their structural context. While unassisted reaction mechanisms preferred C3´-endo pucker that is energetically preferred by the RNA backbone, a spontaneous repuckering 10 ACS Paragon Plus Environment

Page 11 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


due to structural context was observed in several mechanisms involving GA and/or GB groups. For the sake of completeness, we compared the observed structural variability of the sugarphosphate backbone in the studied mechanisms with the reactive conformations of small selfcleaving ribozymes identified by recent QM/MM studies.40-41, 59, 85 The comparison revealed that ribozymes utilized reactive conformations both with C2´-endo (glmS and HDV ribozymes) and C3´-endo (hairpin ribozyme) puckering (see Figure 2). Dianionic Scenarios The uncatalyzed scenario I corresponds to the nucleophilic attack of the activated 2´-O− group already deprotonated by solvent on the negatively charged phosphodiester. The reaction proceeds through a stable IN state, with the first TS corresponding to the nucleophilic attack (TSNA) and the rate-limiting TS corresponding to the exocyclic cleavage event (TSEC, ∆G‡300K = 27.1 kcal/mol, Table S6 in the SI). The inclusion of an external group acting as a GA (scenario II) reduced the barrier of TSEC to such extent that made the passage through TSNA the rate-limiting step with corresponding ∆G‡300K barrier equaling to 19.9–22.7 kcal/mol (the greatest and smallest reductions were achieved with [Mg(H2O)6]2+ and AdeH+ as the GA, respectively). In scenario III, a GB deprotonated the 2´-OH nucleophile, making the nucleophilic attack barrier-less. Consequently, the reaction proceeded through a single TS (TSEC). Most external GBs had no effect on the overall ∆G‡300K barrier (associated with TSEC), but Gua- acting as a GB reduced this ∆G‡300K barrier to 24.5 kcal/mol most likely due to its simultaneous electrostatic TS-stabilizing effect. Mechanisms involving both GB and GA, corresponding to scenario IV, reduced the ∆G‡300K barrier in 6 out of 16 cases. The rate-limiting step in this scenario was typically associated with passage through TSEC (see Figure 3, Tables 2, 5, and S6, and the detailed description of all dianionic scenarios in the Supporting Information). Table 2: Summary of all calculated rate-limiting ∆G‡300K barriers (in kcal/mol) for dianionic scenarios. Top left value (in bold) corresponds to unassisted scenario I, other values in the first row and column correspond to paths of scenarios II and III, respectively. The remaining values belong to paths from scenario IV. See Table S6 in Supporting Information for energies of all states along the reaction paths. -

AdeH+ as GA

CytH+ as GA

27.1

22.7

21.8

19.3 a

19.9

Gua− as GB

24.5

29.3

30.7

25.6

22.2

Cyt as GB

26.3

33.8

30.7

25.8

27.8

[Mg(H2O)4OH]+ as GB

27.7

30.7

31.7

28.3

24.6

[Mg(H2O)5OH]+ as GB

27.5

33.9

32.1

27.0

24.7

GB/GA -

a

[Mg(H2O)5]2+ as GA

[Mg(H2O)6]2+ as GA

the [Mg(H2O)5]2+ ion acted as both GB and GA during the reaction, thus this path corresponds to scenario IV.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 31

Figure 3: Structures of R, rate-limitting TS, and P states along most feasible pathways (with the lowest ∆G‡300K free energy barrier) in four dianionic scenarios I to IV (see Figure 1). Free energies in kcal/mol of each state with respect to ground state of reactants, i.e., including free energy correction due to rare protonation form, are depicted in brackets. Monoanionic Scenarios The uncatalyzed scenario V involves activation of the 2´-OH group by proton transfer to the pro-RP (endo path) or pro-SP (exo path) nbO of the scissile phosphate, reorientation of the resulting OH group of the phosphorane IN, and finally proton transfer to the leaving O5´ group (see the SI for details). The reaction coordinate features two stable phosphorane IN states (IN1, IN2) that are together with R and P states separated by three TS states (TSNA, TSreor, and TSEC). Passage through TSEC was the rate-limiting step, and the overall ∆G‡300K barriers were 34.8 kcal/mol and 33.6 kcal/mol for the endo and exo paths, respectively (Table S7 in the SI). Henceforth only the lower ∆G‡300K barriers for the endo and exo paths will be discussed (see Supporting Information for complete data). Compared to this reference scenario, the inclusion of external groups acting as es-GBs and/or es-GAs (scenario VI) reduced the ∆G‡300K barriers in 8 12 ACS Paragon Plus Environment

Page 13 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


out of 23 paths (Table 3 and Supporting Information Table S7). The most favorable path was [Mg(H2O)6]2+/[Mg(H2O)5]2+, where [Mg(H2O)6]2+ acted as an es-GB together with [Mg(H2O)5]2+ as an es-GA (∆G‡300K 30.1 kcal/mol). Passage through TSEC remained rate-limiting in this case (Table S7 in the Supporting Information). Table 3: Summary of all calculated rate-limiting ∆G‡300K barriers (in kcal/mol) for monoanionic scenarios without direct participation of external chemical group on proton transfers. The external groups, if present, provide only the electrostatic stabilization effect, i.e., act as es-GB and/or esGA. All barriers represent the minimal value of two competing endo and exo paths (see Table S7 in Supporting Information for full data, i.e., both routes and energies of all states along the reaction paths). The top left value (in bold) corresponds to path of the unassisted scenario V, while other values belong to paths from scenario VI. es-GB / es-GA

-

AdeH+ as es-GA

CytH+ as es-GA

[Mg(H2O)5]2

[Mg(H2O)6]2

as es-GA

as es-GA

+

+

-

33.6

35.8

35.2

30.3

30.5

Gua as es-GB

33.0

36.1

35.2

30.8

31.7

CytH+ as es-GB

38.6

40.9

39.7

34.4

35.6

[Mg(H2O)5]2+ as es-GB

38.8

38.6

37.5

-a

33.4


36.7

38.6

36.3

30.1

33.2

a

we were not able to obtain reaction path with two [Mg(H2O)5]2+ ions (both with inner-shell arrangement).

In scenario VII, the 2´-OH group was activated by an nbO, and the leaving O5´ group was protonated by an external GA. In total, 7 out of 19 paths provided lower barriers than scenario V, with the most favorable path (∆G‡300K = 30.1 kcal/mol) involving Gua and [Mg(H2O)5]2+ as the es-GB and GA, respectively (Table 4 and Supporting Table S8 in SI). Direct protonation by AdeH+ made passage through TSNA rate-limiting in three paths, but in all other cases passage through TSEC remained rate-limiting (see Figure 4, Table S8 and Supporting Information). Table 4: Summary of calculated rate-limiting ∆G‡300K barriers (kcal/mol) for scenario VII with direct protonation by GA in the second reaction step. All values represent the minimal barrier of two competing endo and exo paths (see Table S8 in Supporting Information for full data, i.e., both routes and energies of all states along the reaction paths). Each column corresponds to the direct GA role (also with the simultaneous electrostatic stabilization, i.e, the indirect es-GA role) of specific nucleobase or hydrated Mg2+ ion. Rows represent possible involvement of external groups acting as es-GB. GA / es-GB / es-GA

AdeH+ as GA / es-GA

CytH+ as GA / es-GA

[Mg(H2O)5]2+ as GA / es-GA

[Mg(H2O)6]2+ as GA / es-GA

-

33.5

33.5

33.0

31.1

Gua as es-GB

33.8

33.1

30.1

32.0



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

CytH+ as es-GB

38.7

36.7

37.2

35.9


36.2

35.2

-a

35.5


36.8

36.3

33.7

34.0

a

Page 14 of 31

we were not able to obtain reaction path with two [Mg(H2O)5]2+ ions (both with inner-shell arrangement).

Figure 4: Structures of R, rate-limitting TS, and P states along most feasible pathways (with the lowest ∆G‡300K free energy barrier) in three monoanionic scenarios V to VII (see Figure 1). Free energy in kcal/mol of each state with respect to ground state of reactants, i.e., including free energy correction due to rare protonation form is depicted in brackets. Table 5: Overall rate-limiting ∆G‡300K barriers summarizing data from Tables 2-4, where only the path with the lowest ∆G‡300K barrier is shown for a particular scenario (see Tables S6, S7 and S8 in the Supporting Information for complete data including all free energies along reaction paths of all studied mechanisms). Scenario

∆G‡300K

GB

GA

es-GB

es-GA

I

27.1

-

-

-

-

II

19.9

-

[Mg(H2O)6]2+

-

[Mg(H2O)6]2+


Page 15 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


III

24.5

Gua−

-

Gua−

-

IV

22.2

Gua−

[Mg(H2O)6]2+

Gua−

[Mg(H2O)6]2+

V

33.6

nbO

nbO

-

-

VI

30.1

nbO

nbO

VII

30.1

nbO

[Mg(H2O)5]2+

[Mg(H 2+ O) 2 6] Gua

[Mg(H2O)5]2+ [Mg(H2O)5]2+

Unassisted sugar-phosphate backbone cleavage Experimental measurements of uncatalyzed self-cleavage rates for small RNA phosphate models86 or chimeric DNA/RNA oligonucleotides69 in water have yielded estimated ∆G‡300K barriers of ~30 kcal/mol. The barriers we calculated for the uncatalyzed dianionic scenario I and monoanionic scenario V of the RNA self-cleavage mechanisms were 27.1 and 33.6 kcal/mol, respectively, both of which are close to the experimental value. These values are slightly lower than those obtained in previous computational studies, which were ~30.0 and ~33.9 for the dianionic and monoanionic unassisted pathways (Table S5 in the Supporting Information).37-41 While the dianionic scenario I is clearly preferred over the monoanionic scenario V under the given conditions (300 K and pH 7, Table 5), the difference between their barriers is comparable to the typical catalytic effect of ribozymes (6-8 kcal/mol, see below), therefore, the monoanionic mechanism may still be an important pathway for RNA cleavage that might compete with the dianionic mechanism in some circumstances. However, the dianionic mechanism is undoubtedly expected to be the dominant reaction channel under more basic conditions because the ∆Gcorr for the rare deprotonation of the 2´-OH nucleophile (and thus also the overall reaction barrier) will be lower in such conditions compared to physiological pH 7 (Table S2 in the SI). The calculated ∆G‡300K barriers for the uncatalyzed scenarios can also be used to estimate the catalytic effect of the RNA enzymes. Kinetic experiments on various naturally occurring small self-cleaving ribozymes have yielded estimated reaction barriers of 19.4 to 20.9 kcal/mol (Table 6),87-90 indicating that ribozymes reduce the cleavage barrier by at least ~6-8 kcal/mol. Ribozymes are thus efficient biocatalysts that are at least somewhat comparable to protein enzymes, which reduce reaction barriers by 6 to 34 kcal/mol.7, 91 Effect of exogenous chemical groups on the cleavage reaction A GB alone only reduced the ∆G‡300K barrier by a small amount (~0.6 kcal/mol) compared to scenario I. In keeping with previous reports, nucleobases were more effective GBs than hydrated Mg2+ ions56 (Table 2). The activation of the 2´-OH nucleophile by a GB is thus beneficial but not essential for RNA cleavage. Conversely, a GA with no other exogenous catalyst reduced the ∆G‡300K barrier by at least 4.4 kcal/mol compared to scenario I (see Tables 2 and 5), showing that GAs can strongly facilitate the RNA cleavage reaction and may contribute significantly to ribozymes’ catalytic effects. Because both GBs and GAs were beneficial for catalysis, one might expect a cooperative participation of both GB and GA (scenario IV) to yield significantly reduced ∆G‡300K barriers. However, combined GB + GA mechanisms provided lower ∆G‡300K barriers than the reference scenario I only in 6 out of 16 paths (Table 2). A similar trend was observed for the indirect cooperative participation of a GB and a GA via electrostatic stabilization of transition state (scenarios VI and VII, Tables 3 and 4). It thus appears that some combinations of nucleobases and 15 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 31

hydrated Mg2+ ions are even unfavorable for catalysis, especially if they only participate indirectly via electrostatic interactions. On the other hand, there were some reaction paths in which the direct participation of both GB and GA significantly reduced the ∆G‡300K barrier, mimicking the catalytic effect of naturally occurring ribozymes (see below). Most of the calculated ∆G‡300K barriers from scenarios IV, VI and VII were dominantly determined by the pKa correction applied to specific functional groups (∆Gcorr, Table S2 in Supporting Information). The accuracy of our description may therefore have been limited by the use of unperturbed rather than perturbed pKa values. Moreover, the ∆Gcorr values were estimated without accounting for potential cooperative effects whereby the (de)protonation of one group in the reaction complex might affect the (de)protonation propensity of another group. Consequently, there may be some residual inaccuracies in the estimated ∆G‡300K barriers despite the careful adjustment of the computational protocol. Still, the approximate nature of implicit solvent models is one of the main sources of uncertainty in our study (see section “Probing the accuracy of implicit solvent model by explicit solvent calculations” below). Implications for the self-cleavage of ribozymes As mentioned above, the aim of this study was to evaluate the part of the catalytic effect associated with the shift of the reaction mechanism and explicit involvement of the active site groups in the chemical reaction step. Although the implicit effects of the enzyme environment also contribute to the catalytic effects, they are beyond the scope of the present study, since the model systems were chosen to separate the explicit involvement of the active site groups from the rest of the catalytic effects. This also delineates limits of the comparison of the obtained results with the overall catalytic effects of real ribozyme systems. In addition, besides the discussed explicit and implicit catalytical effects affecting the chemical reaction step, the apparent reaction rate of ribozymes’ self-cleavage reaction as observed in experiments might be also affected by the kinetics of folding and structural rearrangements of the ribozymes required for formation of the proper reactive conformation. This effect is, e.g., responsible for different kinetics of the selfcleavage in minimal and full-length hammerhead ribozyme.92 The full-length hammerhead ribozyme is almost 1000-fold more active than the minimal sequence hammerhead ribozyme. The mechanism of full-length hammerhead ribozyme utilize deprotonated guanine and hydrated Mg2+ ion as GB and GA, respectively,93-94 while in the minimal hammerhead ribozyme both GA and GB roles are played by hydrated Mg2+ ions.89, 95 Our data suggest that the reactions (in water environment) following these two mechanisms are comparable in their free energy barriers, so although these two forms of hammerhead ribozymes differ in the reaction mechanism of selfcleavage, the difference in their activities cannot be attributed to the shift of reaction mechanism (see Table 2). In turn, as suggested elsewhere, the difference can be dominantly explained by coupling of the chemical reaction step with slow conformational rearrangement required for the formation of the reactive conformation in the minimal hammerhead ribozyme.96 Structural and biochemical studies have indicated that nucleobases and hydrated Mg2+ ions participate in the autocatalytic mechanism of small self-cleaving ribozymes (Table 6).15, 42-50, 52, 9799 Most small self-cleaving ribozymes (e.g., hairpin, hammerhead and glmS ribozymes) contain an active site Gua moiety in close proximity to the 2´-OH nucleophile, suggesting a role as a GB.8, 49, 52, 100-110 We found that direct involvement of the GB in the activation of 2´-OH nucleophile brought only marginal catalytical effect. In particular, the reaction pathways where a deprotonated Gua− acted as a GB had only slightly lower barriers than paths where Gua played only an indirect role in electrostatic stabilization of the TS (es-GB). Note that the electrostatically stabilizing role of Gua is implicit even in scenarios where it acts directly as a GB. Thus, although 16 ACS Paragon Plus Environment

Page 17 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


the activation of the 2´-OH nucleophile by deprotonated Gua− appears to be a slightly preferred mechanism, its contribution to the overall catalytic effect seems relatively modest compared to the effect of the same guanine in the electrostatic TS stabilization. For pathways involving deprotonated Gua−, the TS is stabilized by the withdrawal of electron density from the nbO via hydrogen bonding with the exocyclic N2 amino group. A similar electron flow occurs in the HDV ribozyme, where a hydrated Mg2+ ion acts as both a GB and a Lewis acid via direct innershell contact with the scissile phosphate.50, 77, 90, 111 Our calculations suggest that the dominant catalytic effect stems from the presence of a GA. Hairpin, neurospora Varkud satellite (VS), and twister ribozymes all contain an active site Ade residue in close proximity to the O5´ leaving group, suggesting a role as a GA. The direct involvement of protonated AdeH+ as a GA is strongly favored over the alternative pathway where it acts as an es-GA and contributes only electrostatic stabilization (see scenarios IV, VI, and VII). In HDV and hammerhead ribozymes, the catalytic GA is a protonated CytH+ or the 2´-OH (or a water molecule) coordinated to the Mg2+ ion, respectively.95, 112-114 Notably, the dianionic scenario IV, in which Gua− and [Mg(H2O)6]2+ act as a GB and GA, respectively, exhibited one of the lowest ∆G‡300K barriers of all the investigated reaction mechanisms (Table 2). Table 6: Nucleobases and hydrated Mg2+ ions suggested to act as general acid-base catalysts in small self-cleaving ribozymes, with the corresponding experimentally observed reaction barriers (∆G‡300K, kcal/mol). ribozyme

∆G‡300K

GB/es-GB

GA/es-GA

hairpin115

20.9

Gua−

AdeH+

hammerhead116

19.9

Gua−

Mg2+…2´-OH/Mg2+…H2O

HDV117

19.9

Mg2+…OH−

CytH+

glmS118

19.8

Gua−

GlcN6P

VS119

20.1

Gua−

AdeH+

twister52

19.4

Gua−

(N3-protonated) AdeH+

Probing the accuracy of implicit solvent model by explicit solvent calculations As noted above, the accuracy of present calculations is mainly limited by the approximate nature of the implicit solvent model. The overall accuracy of QM calculations is typically given by trade-off between level of the QM method, description of the solvation (implicit vs. explicit model), completeness of the studied model, and level of sampling. In principle, an extension of the studied model by an additional layer of explicit water molecules may increase completeness of the solvation description. However, this approach is not a panacea, as it increases conformational complexity of the studied system and consequently also the computational cost. Further, the explicit water molecules may lead to occurrence of multiple local minima on the potential energy surface, a common problem when immersing nucleic acids in a cluster of waters.120 This may add noise to the energy analyses. The free-energy contribution of the solvent, which is inherently included in the implicit solvent models, has to be calculated by proper sampling techniques. However, this might seriously limit the choice of the available QM methods to some less expensive semi-empirical techniques, which performance in turn might not be sufficient for comparison of various reaction mechanisms.80 Nevertheless, upon the suggestion 17 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 31

of the Reviewers, we have decided to perform a series of QM calculations with explicit water molecules, in order to estimate the accuracy of the used implicit solvent model. We compared free-energy profiles of both uncatalyzed (unassisted) reaction scenarios (i.e., scenarios I and V) calculated in implicit as well as explicit solvent models. We added a layer of explicit water molecules around the models used in the scenarios I and V, so that these extended models contained sugar-phosphate backbone immersed in a water droplet (~190 QM atoms in total) surrounded by implicit solvent as in other calculations (see the Methods section). Note that the added shell of explicit water molecules did not participate in any proton transfer. In all cases in which we assumed direct participation of water molecules, i.e., those scenarios where hydrated Mg2+ ion acts as general base or general acid, the complete first-shell solvation sphere of the Mg2+ ion has always been described explicitly even in cases denoted as implicit solvent calculations (see Methods). We observed that the explicit solvent description increased energies of both TS and IN states, and consequently the ∆G‡300K barriers, of the unassisted dianionic scenario I by ~3.9 kcal/mol compared to the implicit solvent calculations. The rate-limiting step remained connected with the exocyclic cleavage step (TSEC, ∆G‡300K = 30.8 kcal/mol, Table S6 in the Supporting Information). Similar effect was observed even in monoanionic scenarios V, in which energies of TS and IN states were increased by ~4 kcal/mol and ~3 kcal/mol, respectively. The rate-limiting steps again remained connected with the exocyclic cleavage step (TSEC, Table S7 in Supporting Information). Unfortunately, we were not able to perform these explicit solvent calculations for more complex systems including exogenous groups acting as general acid-base catalysts due to their enormous computational demands. Nevertheless, based on the above described comparison, we suggest that our implicit solvent model is accurate enough to distinguish correct reaction mechanism and ratelimiting step. Although it might have a tendency to underestimate free-energy barriers by ~3-4 kcal/mol, this underestimation seems to be systematic and thus the differences between reaction barriers obtained by implicit solvent calculations are suitable for meaningful comparison of various studied reaction mechanisms.

Conclusions Our calculations revealed that several phosphodiester cleavage reaction paths involving exogenous chemical groups (especially those featuring GAs that protonate the leaving O5´ group) yielded cleavage rates comparable to those of naturally occurring small self-cleaving ribozymes (0.2-2.8 min-1). We thus suggest that the dominant part of the ribozyme catalytic effect stems from shifting the reaction mechanism including the preorganization of the active site with appropriately positioned GA and GB moieties. In addition we found that the GA makes the biggest contribution to the ribozyme’s catalytic effect. In other words, the above-mentioned effect of the shifting the reaction mechanism by involvement of GA and GB moieties is sufficient to explain dominant part of the catalytic effect of ribozymes. In contrast, the catalytic effect in protein enzymes is dominantly achieved by implicit influence of the enzyme environment, in particular by the electrostatic preorganization of the active site. The implicit effect of protein enzyme environment may contribute to the catalytical effect by a few tens of kcal/mol (up to ~30 kcal/mol),7 which is much larger than even the overall catalytic effect of ribozymes (6-8 kcal/mol). Therefore, the ribozymes appear to exploit the implicit catalytic effects of enzyme environment less efficiently than protein enzymes.7 Thus we propose that the catalytic effect of ribozymes may be largely explainable by explicit involvement of active site groups directly participating in the self-cleavage mechanism. In other words, RNA and protein catalysis could 18 ACS Paragon Plus Environment

Page 19 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


feature different propensities to use explicit and implicit catalytic effects. We, however, in no case suggest that implicit effects, folding and other factors that are not included in our present computations do not contribute to the ribozyme catalysis. The difference between protein enzymes and ribozymes may be related to the RNA’s polyanionic nature, which is compensated by a diffuse shell of counter-ions. The corresponding strong electrostatic field generated by polyanionic RNA and its counter-ions is primarily dictated by overall folds of RNA molecules, i.e., by the disposition of the phosphate groups. In contrast to proteins, the RNA electrostatic field cannot be easily fine-tuned by the enzyme sequence. This hinders tailoring of a specific electric field in ribozymes’ active sites to support catalyzed reactions. We thus suggest that the limited capability of RNA enzymes to fine-tune the longrange electrostatic catalytic effects (and the need to rely only on direct involvement of active site chemical groups) may be one of the reasons why evolution ultimately preferred protein enzymes to ribozymes in most contexts.

Acknowledgements This work was supported by grant P208/12/1878 (J.S., M.O., P.K.) from the Czech Science Foundation. Further funding was provided by project LO1305 of the Ministry of Education, Youth and Sports of the Czech Republic (M.O., P.B., P.K.). JS acknowledges support from Praemium Academiae.

Supporting Information Supporting Information is available free of charge on the ACS Publications website and contains these sections: accuracy of MPW1K functional, reference state of the reactions, tests and adjustments of the computational protocol, detailed description of dianionic and monoanionic scenarios, detailed description how endogenous chemical groups affects the cleavage reaction, supporting tables, figures and Cartesian coordinates of all reaction states.

References (1) Grabowski, P. J.; Zaug, A. J.; Cech, T. R. The Intervening Sequence of the RibosomalRNA Precursor Is Converted to a Circular RNA in Isolated-Nuclei of Tetrahymena. Cell 1981, 23, 467-476. (2) Guerriertakada, C.; Gardiner, K.; Marsh, T.; Pace, N.; Altman, S. The RNA Moiety of Ribonuclease-P Is the Catalytic Subunit of the Enzyme. Cell 1983, 35, 849-857. (3) Miller, B. G.; Wolfenden, R. Catalytic proficiency: The unusual case of OMP decarboxylase. Annu. Rev. Biochem. 2002, 71, 847-885. (4) Kraut, D. A.; Carroll, K. S.; Herschlag, D. Challenges in enzyme mechanism and energetics. Annu. Rev. Biochem. 2003, 72, 517-571. (5) Benkovic, S. J.; Hammes-Schiffer, S. A perspective on enzyme catalysis. Science 2003, 301, 1196-1202. (6) Garcia-Viloca, M.; Gao, J.; Karplus, M.; Truhlar, D. G. How enzymes work: Analysis by modern rate theory and computer simulations. Science 2004, 303, 186-195. (7) Warshel, A.; Sharma, P. K.; Kato, M.; Xiang, Y.; Liu, H. B.; Olsson, M. H. M. Electrostatic basis for enzyme catalysis. Chem. Rev. 2006, 106, 3210-3235. 19 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 31

(8) Fedor, M. J. Comparative Enzymology and Structural Biology of RNA Self-Cleavage. Annu. Rev. Biophys. 2009, 38, 271-299. (9) Lilley, D. M. J. The origins of RNA catalysis in ribozymes. Trends Biochem. Sci. 2003, 28, 495-501. (10) Doudna, J. A.; Lorsch, J. R. Ribozyme catalysis: not different, just worse. Nat. Struct. Mol. Biol. 2005, 12, 395-402. (11) Fedor, M. J.; Williamson, J. R. The catalytic diversity of RNAS. Nat. Rev. Mol. Cell Bio. 2005, 6, 399-412. (12) Herschlag, D. Learning from ribozymes. RNA 2015, 21, 527-528. (13) Eckstein, F.; Lilley, D. M. J. Catalytic RNA. Springer-Verlag Berlin Heidelberg: 1997. (14) Ferre-D'Amare, A. R.; Scott, W. G. Small Self-cleaving Ribozymes. Cold Spring Harbor Perspectives in Biology 2010, 2. (15) Cochrane, J. C.; Strobel, S. A. Catalytic strategies of self-cleaving ribozymes. Accounts. Chem. Res. 2008, 41, 1027-1035. (16) Ward, W. L.; Plakos, K.; DeRose, V. J. Nucleic Acid Catalysis: Metals, Nucleobases, and Other Cofactors. Chem. Rev. 2014, 114, 4318-4342. (17) Jimenez, R. M.; Polanco, J. A.; Luptak, A. Chemistry and Biology of Self-Cleaving Ribozymes. Trends Biochem. Sci. 2015, 40, 648-661. (18) Müller, S.; Appel, B.; Balke, D.; Hieronymus, R.; Nübel, C. Thirty-five years of research into ribozymes and nucleic acid catalysis: where do we stand today? [version1; referees: 2 approved]. F1000Research 2016, 5(F1000 Faculty Rev). (19) Lee, K. Y.; Lee, B. J. Structural and Biochemical Properties of Novel Self-Cleaving Ribozymes. Molecules 2017, 22, 678. (20) Soukup, G. A.; Breaker, R. R. Relationship between internucleotide linkage geometry and the stability of RNA. RNA 1999, 5, 1308-1325. (21) Warshel, A. Computer Modeling of Chemical Reactions in Enzymes and Solutions. John Wiley & Sons, Inc.: New York, NY, 1997. (22) Gao, J.; Ma, S.; Major, D. T.; Nam, K.; Pu, J.; Truhlar, D. G. Mechanisms and free energies of enzymatic reactions. Chem. Rev. 2006, 106, 3188-3209. (23) Cramer, C. J. Essentials of Computational Chemistry: Theories and Models, 2nd Edition. John Wiley & Sons, Inc.: New York: NY, 2004. (24) Field, M. Simulating enzyme reactions: challenges and perspectives. J. Comput. Chem. 2002, 23, 48-58. (25) Olsson, M. H. M.; Parson, W. W.; Warshel, A. Dynamical contributions to enzyme catalysis: Critical tests of a popular hypothesis. Chem. Rev. 2006, 106, 1737-1756. (26) Emilsson, G. M.; Nakamura, S.; Roth, A.; Breaker, R. R. Ribozyme speed limits. RNA 2003, 9, 907-918. (27) Breaker, R. R.; Emilsson, G. M.; Lazarev, D.; Nakamura, S.; Puskarz, I. J.; Roth, A.; Sudarsan, N. A common speed limit for RNA-cleaving ribozymes and deoxyribozymes. RNA 2003, 9, 949-957. (28) Oivanen, M.; Lonnberg, H. Kinetics and Mechanisms for Reactions of Adenosine 2´- and 3´-Monophosphates in Aqueous Acid - Competition between Phosphate Migration, Dephosphorylation, and Depurination. J. Org. Chem. 1989, 54, 2556-2560. (29) Kosonen, M.; Yousefi-Salakdeh, E.; Stromberg, R.; Lonnberg, H. pH- and bufferindependent cleavage and mutual isomerization of uridine 2´- and 3´-alkyl phosphodiesters: implications for the buffer catalyzed cleavage of RNA. Journal of the Chemical Society-Perkin Transactions 2 1998, 1589-1595. 20 ACS Paragon Plus Environment

Page 21 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(30) Mikkola, S.; Kosonen, M.; Lonnberg, H. Cleavage and isomerization of RNA phosphodiester bonds: Nucleoside phosphotriesters and chimeric ribo/2´-O-methylribo oligonucleotides as tools for mechanistic studies. Curr. Org. Chem. 2002, 6, 523-538. (31) Perreault, D. M.; Anslyn, E. V. Unifying the current data on the mechanism of cleavage Transesterification of RNA. Angew. Chem. Int. Edit. 1997, 36, 432-450. (32) Wong, K. Y.; Gu, H.; Zhang, S. M.; Piccirilli, J. A.; Harris, M. E.; York, D. M. Characterization of the Reaction Path and Transition States for RNA Transphosphorylation Models from Theory and Experiment. Angew Chem Int Ed 2012, 51, 647-651. (33) Sponer, J.; Sponer, J. E.; Mladek, A.; Banas, P.; Jurecka, P.; Otyepka, M. How to understand quantum chemical computations on DNA and RNA systems? A practical guide for non-specialists. Methods 2013, 64, 3-11. (34) Duarte, F.; Aqvist, J.; Wiliams, N. H.; Kamerlin, S. C. L. Resolving Apparent Conflicts between Theoretical and Experimental Models of Phosphate Monoester Hydrolysis. J. Am. Chem. Soc. 2015, 137, 1081-1093. (35) Huang, M.; York, D. M. Linear free energy relationships in RNA transesterification: theoretical models to aid experimental interpretations. Phys. Chem. Chem. Phys. 2014, 16, 1584615855. (36) Nam, K. H.; Gaot, J. L.; York, D. M. Quantum mechanical/molecular mechanical simulation study of the mechanism of hairpin ribozyme catalysis. J. Am. Chem. Soc. 2008, 130, 4680-4691. (37) Boero, M.; Terakura, K.; Tateno, M. Catalytic role of metal ion in the selection of competing reaction paths: A first principles molecular dynamics study of the enzymatic reaction in ribozyme. J. Am. Chem. Soc. 2002, 124, 8949-8957. (38) Lopez, X.; Dejaegere, A.; Leclerc, F.; York, D. M.; Karplus, M. Nucleophilic attack on phosphate diesters: A density functional study of in-line reactivity in dianionic, monoanionic, and neutral systems. J. Phys. Chem. B 2006, 110, 11525-11539. (39) Liu, H. N.; Robinet, J. J.; Ananvoranich, S.; Gauld, J. W. Density functional theory investigation on the mechanism of the hepatitis delta virus ribozyme. J. Phys. Chem. B 2007, 111, 439-445. (40) Banas, P.; Rulisek, L.; Hanosova, V.; Svozil, D.; Walter, N. G.; Sponer, J.; Otyepka, M. General base catalysis for cleavage by the active-site cytosine of the hepatitis delta virus ribozyme: QM/MM calculations establish chemical feasibility. J. Phys. Chem. B 2008, 112, 11177-11187. (41) Mlynsky, V.; Banas, P.; Walter, N. G.; Sponer, J.; Otyepka, M. QM/MM Studies of Hairpin Ribozyme Self-Cleavage Suggest the Feasibility of Multiple Competing Reaction Mechanisms. J. Phys. Chem. B 2011, 115, 13911-13924. (42) Rupert, P. B.; Ferre-D'Amare, A. R. Crystal structure of a hairpin ribozyme-inhibitor complex with implications for catalysis. Nature 2001, 410, 780-786. (43) Bevilacqua, P. C. Mechanistic considerations for general acid-base catalysis by RNA: Revisiting the mechanism of the hairpin ribozyme. Biochemistry-Us 2003, 42, 2259-2265. (44) Lee, T. S.; Lopez, C. S.; Giambasu, G. M.; Martick, M.; Scott, W. G.; York, D. M. Role of Mg2+ in hammerhead ribozyme catalysis from molecular simulation. J. Am. Chem. Soc. 2008, 130, 3053-3064. (45) Guo, M.; Spitale, R. C.; Volpini, R.; Krucinska, J.; Cristalli, G.; Carey, P. R.; Wedekind, J. E. Direct Raman Measurement of an Elevated Base pK(a) in the Active Site of a Small Ribozyme in a Precatalytic Conformation. J. Am. Chem. Soc. 2009, 131, 12908-12909.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 31

(46) Wilcox, J. L.; Ahluwalia, A. K.; Bevilacqua, P. C. Charged Nucleobases and Their Potential for RNA Catalysis. Accounts. Chem. Res. 2011, 44, 1270-1279. (47) Viladoms, J.; Scott, L. G.; Fedor, M. J. An Active-Site Guanine Participates in glmS Ribozyme Catalysis in Its Protonated State. J. Am. Chem. Soc. 2011, 133, 18388-18396. (48) Wilson, T. J.; Lilley, D. M. J. Do the hairpin and VS ribozymes share a common catalytic mechanism based on general acid-base catalysis? A critical assessment of available experimental data. RNA 2011, 17, 213-221. (49) Suslov, N. B.; DasGupta, S.; Huang, H.; Fuller, J. R.; Lilley, D. M. J.; Rice, P. A.; Piccirilli, J. A. Crystal structure of the Varkud satellite ribozyme. Nat. Chem. Biol. 2015, 11, 840846. (50) Chen, J. H.; Yajima, R.; Chadalavada, D. M.; Chase, E.; Bevilacqua, P. C.; Golden, B. L. A 1.9 angstrom Crystal Structure of the HDV Ribozyme Precleavage Suggests both Lewis Acid and General Acid Mechanisms Contribute to Phosphodiester Cleavage. Biochemistry-Us 2010, 49, 6508-6518. (51) Wilson, T. J.; Liu, Y. J.; Domnick, C.; Kath-Schorr, S.; Liiley, D. M. J. The Novel Chemical Mechanism of the Twister Ribozyme. J. Am. Chem. Soc. 2016, 138, 6151-6162. (52) Liu, Y.; Wilson, T. J.; McPhee, S. A.; Lilley, D. M. Crystal structure and mechanistic investigation of the twister ribozyme. Nat. Chem. Biol. 2014, 10, 739-744. (53) Boero, M.; Tateno, M.; Terakura, K.; Oshiyama, A. Double-Metal-Ion/Single-Metal-Ion Mechanisms of the Cleavage Reaction of Ribozymes: First-Principles Molecular Dynamics Simulations of a Fully Hydrated Model System. J. Chem. Theory Comput. 2005, 1, 925-934. (54) De Vivo, M.; Dal Peraro, M.; Klein, M. L. Phosphodiester cleavage in ribonuclease H occurs via an associative two-metal-aided catalytic mechanism. J. Am. Chem. Soc. 2008, 130, 10955-10962. (55) Rosta, E.; Nowotny, M.; Yang, W.; Hummer, G. Catalytic mechanism of RNA backbone cleavage by ribonuclease H from quantum mechanics/molecular mechanics simulations. J. Am. Chem. Soc. 2011, 133, 8934-8941. (56) Chval, Z.; Chvalova, D.; Leclerc, F. Modeling the RNA 2´ OH Activation: Possible Roles of Metal Ion and Nucleobase as Catalysts in Self-Cleaving Ribozymes. J. Phys. Chem. B 2011, 115, 10943-10956. (57) Casalino, L.; Palermo, G.; Rothlisberger, U.; Magistrato, A. Who Activates the Nucleophile in Ribozyme Catalysis? An Answer from the Splicing Mechanism of Group II Introns. J. Am. Chem. Soc. 2016, 138, 10374-10377. (58) Ganguly, A.; Thaplyal, P.; Rosta, E.; Bevilacqua, P. C.; Hammes-Schiffer, S. Quantum mechanical/molecular mechanical free energy simulations of the self-cleavage reaction in the hepatitis delta virus ribozyme. J. Am. Chem. Soc. 2014, 136, 1483-1496. (59) Mlynsky, V.; Walter, N. G.; Sponer, J.; Otyepka, M.; Banas, P. The role of an active site Mg(2+) in HDV ribozyme self-cleavage: insights from QM/MM calculations. Phys. Chem. Chem. Phys. 2015, 17, 670-679. (60) Wong, K. Y.; Lee, T. S.; York, D. M. Active participation of Mg ion in the reaction coordinate of RNA self-cleavage catalyzed by the hammerhead ribozyme. J. Chem. Theory Comput. 2011, 7, 1-3. (61) Lynch, B. J.; Fast, P. L.; Harris, M.; Truhlar, D. G. Adiabatic connection for kinetics. J. Phys. Chem. A 2000, 104, 4811-4815. (62) Lynch, B. J.; Truhlar, D. G. How well can hybrid density functional methods predict transition state geometries and barrier heights? J. Phys. Chem. A 2001, 105, 2936-2941.


Page 23 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(63) Møller, C.; Plesset, M. S. Note on an approximation treatment for many-electron systems. Phys. Rev. 1934, 46, 618. (64) Jurecka, P.; Hobza, P. On the convergence of the (Delta E-CCSD(T)-Delta E-MP2) term for complexes with multiple H-bonds. Chem. Phys. Lett. 2002, 365, 89-94. (65) Jurecka, P.; Hobza, P. True stabilization energies for the optimal planar hydrogen-bonded and stacked structures of guanine···cytosine, adenine···thymine, and their 9-and 1-methyl derivatives: Complete basis set calculations at the MP2 and CCSD(T) levels and comparison with experiment. J. Am. Chem. Soc. 2003, 125, 15608-15613. (66) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A., et al. Gaussian 09. Gaussian, Inc.: Wallingford, CT, USA, 2009. (67) Bondi, A. Van der Waals Volumes and Radii. J. Phys. Chem.-Us 1964, 68, 441-451. (68) Case, D. A.; Cerutti, D. S.; Cheatham, T. E.; Darden, T. A.; Duke, R. E.; Giese, T. J.; Gohlke, H.; Goetz, A. W.; Green, D.; Homeyer, N., et al. AMBER 2017. University of California: San Francisco, 2017. (69) Li, Y. F.; Breaker, R. R. Kinetics of RNA degradation by specific base catalysis of transesterification involving the 2´-hydroxyl group. J. Am. Chem. Soc. 1999, 121, 5364-5372. (70) Lyne, P. D.; Karplus, M. Determination of the pK(a) of the 2´-hydroxyl group of a phosphorylated ribose: Implications for the mechanism of hammerhead ribozyme catalysis. J. Am. Chem. Soc. 2000, 122, 166-167. (71) Velikyan, I.; Acharya, S.; Trifonova, A.; Foldesi, A.; Chattopadhyaya, J. The pK(a)'s of 2´-hydroxyl group in nucleosides and nucleotides. J. Am. Chem. Soc. 2001, 123, 2893-2894. (72) Izatt, R. M.; Hansen, L. D.; Rytting, J. H.; Christensen, J. J. Proton Ionization from Adenosine. J. Am. Chem. Soc. 1965, 87, 2760-2761. (73) Acharya, S.; Foldesi, A.; Chattopadhyaya, J. The pK(a) of the internucleotidic 2´hydroxyl group in diribonucleoside (3 '-> 5 ') monophosphates. J. Org. Chem. 2003, 68, 19061910. (74) Usher, D. A.; Richardson, D. I.; Oakenfull, D. G. Models of Ribonuclease Action .2. Specific Acid, Specific Base, and Neutral Pathways for Hydrolysis of a Nucleotide Diester Analog. J. Am. Chem. Soc. 1970, 92, 4699-4712. (75) Jarvinen, P.; Oivanen, M.; Lonnberg, H. Interconversion and Phosphoester Hydrolysis of 2´,5´-Dinucleoside and 3´,5´-Dinucleoside Monophosphates - Kinetics and Mechanisms. J. Org. Chem. 1991, 56, 5396-5401. (76) Davies, J. E.; Doltsinis, N. L.; Kirby, A. J.; Roussev, C. D.; Sprik, M. Estimating pK(a) values for pentaoxyphosphoranes. J. Am. Chem. Soc. 2002, 124, 6594-6599. (77) Sripathi, K. N.; Tay, W. W.; Banas, P.; Otyepka, M.; Sponer, J.; Walter, N. G. Disparate HDV ribozyme crystal structures represent intermediates on a rugged free-energy landscape. RNA 2014, 20, 1112-1128. (78) Radak, B. K.; Lee, T. S.; Harris, M. E.; York, D. M. Assessment of metal-assisted nucleophile activation in the hepatitis delta virus ribozyme from molecular simulation and 3DRISM. RNA 2015, 21, 1566-1577. (79) Panteva, M. T.; Dissanayake, T.; Chen, H. Y.; Radak, B. K.; Kuechler, E. R.; Giambasu, G. M.; Lee, T. S.; York, D. M. Multiscale Methods for Computational RNA Enzymology. Computational Methods for Understanding Riboswitches 2015, 553, 335-374. (80) Mlynsky, V.; Banas, P.; Sponer, J.; van der Kamp, M. W.; Mulholland, A. J.; Otyepka, M. Comparison of ab Initio, DFT, and Semiempirical QM/MM Approaches for Description of Catalytic Mechanism of Hairpin Ribozyme. J. Chem. Theory Comput. 2014, 10, 1608-1622. 23 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 31

(81) Huang, M.; Giese, T. J.; York, D. M. Nucleic acid reactivity: Challenges for nextgeneration semiempirical quantum models. J. Comput. Chem. 2015, 36, 1370-1389. (82) Mlynsky, V.; Bussi, G. Understanding in-line probing experiments by modeling cleavage of nonreactive RNA nucleotides. RNA 2017, 23, 712-720. (83) Marzec, C. J.; Day, L. A. An exact description of five-membered ring configurations. I. Parameterization via an amplitude S, an angle gamma, the pseudorotation amplitude q and phase angle P, and the bond lengths. J. Biomol. Struct. Dyn. 1993, 10, 1091-1123. (84) Huang, M.; Giese, T. J.; Lee, T. S.; York, D. M. Improvement of DNA and RNA Sugar Pucker Profiles from Semiempirical Quantum Methods. J. Chem. Theory Comput. 2014, 10, 1538-1545. (85) Dubecky, M.; Walter, N. G.; Sponer, J.; Otyepka, M.; Banas, P. Chemical feasibility of the general acid/base mechanism of glmS ribozyme self-cleavage. Biopolymers 2015, 103, 550562. (86) Breslow, R.; Huang, D. L.; Anslyn, E. On the mechanism of action of ribonucleases: dinucleotide cleavage catalyzed by imidazole and Zn2+. Proc. Natl. Acad. Sci. USA 1989, 86, 1746-1750. (87) Fedor, M. J. Structure and function of the hairpin ribozyme. J. Mol. Biol. 2000, 297, 269291. (88) Nakano, S.; Chadalavada, D. M.; Bevilacqua, P. C. General acid-base catalysis in the mechanism of a hepatitis delta virus ribozyme. Science 2000, 287, 1493-1497. (89) Torres, R. A.; Himo, F.; Bruice, T. C.; Noodleman, L.; Lovell, T. Theoretical examination of Mg2+-mediated hydrolysis of a phosphodiester linkage as proposed for the hammerhead ribozyme. J. Am. Chem. Soc. 2003, 125, 9861-9867. (90) Das, S. R.; Piccirilli, J. A. General acid catalysis by the hepatitis delta virus ribozyme. Nat. Chem. Biol. 2005, 1, 45-52. (91) Radzicka, A.; Wolfenden, R. A Proficient Enzyme. Science 1995, 267, 90-93. (92) Scott, W. G.; Horan, L. H.; Martick, M. The hammerhead ribozyme: structure, catalysis, and gene regulation. Prog. Mol. Biol. Transl. Sci. 2013, 120, 1-23. (93) Lee, T. S.; Lopez, C. S.; Martick, M.; Scott, W. G.; York, D. M. Insight into the role of Mg in hammerhead ribozyme catalysis from X-ray crystallography and molecular dynamics simulation. J. Chem. Theory Comput. 2007, 3, 325-327. (94) Mir, A.; Chen, J.; Robinson, K.; Lendy, E.; Goodman, J.; Neau, D.; Golden, B. L. Two Divalent Metal Ions and Conformational Changes Play Roles in the Hammerhead Ribozyme Cleavage Reaction. Biochemistry-Us 2015, 54, 6369-6381. (95) Leclerc, F.; Karplus, M. Two-metal-ion mechanism for hammerhead-ribozyme catalysis. J. Phys. Chem. B 2006, 110, 3395-3409. (96) Radhakrishnan, R. Coupling of fast and slow modes in the reaction pathway of the minimal hammerhead ribozyme cleavage. Biophys. J. 2007, 93, 2391-2399. (97) Nam, K.; Gao, J. L.; York, D. M. Electrostatic interactions in the hairpin ribozyme account for the majority of the rate acceleration without chemical participation by nucleobases. RNA 2008, 14, 1501-1507. (98) Gebetsberger, J.; Micura, R. Unwinding the twister ribozyme: from structure to mechanism. Wiley Interdiscip Rev RNA 2016. (99) Roth, A.; Weinberg, Z.; Chen, A. G. Y.; Kim, P. B.; Ames, T. D.; Breaker, R. R. A widespread self-cleaving ribozyme class is revealed by bioinformatics. Nat. Chem. Biol. 2014, 10, 56-60.


Page 25 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(100) Mir, A.; Golden, B. L. Two Active Site Divalent Ions in the Crystal Structure of the Hammerhead Ribozyme Bound to a Transition State Analogue. Biochemistry-Us 2016, 55, 633636. (101) Lau, M. W. L.; Ferre-D'Amare, A. R. An in vitro evolved glmS ribozyme has the wildtype fold but loses coenzyme dependence. Nat. Chem. Biol. 2013, 9, 805-810. (102) Gong, B.; Klein, D. J.; Ferre-D'Amare, A. R.; Carey, P. R. The glmS ribozyme tunes the catalytically critical pK(a) of its coenzyme glucosamine-6-phosphate. J. Am. Chem. Soc. 2011, 133, 14188-14191. (103) Eiler, D.; Wang, J.; Steitz, T. A. Structural basis for the fast self-cleavage reaction catalyzed by the twister ribozyme. Proc. Natl. Acad. Sci. USA 2014, 111, 13028-13033. (104) Gaines, C. S.; York, D. M. Ribozyme Catalysis with a Twist: Active State of the Twister Ribozyme in Solution Predicted from Molecular Simulation. J. Am. Chem. Soc. 2016, 138, 30583065. (105) Banas, P.; Walter, N. G.; Sponer, J.; Otyepka, M. Protonation States of the Key Active Site Residues and Structural Dynamics of the glmS Riboswitch As Revealed by Molecular Dynamics. J. Phys. Chem. B 2010, 114, 8701-8712. (106) Viladoms, J.; Fedor, M. J. The glmS Ribozyme Cofactor is a General Acid-Base Catalyst. J. Am. Chem. Soc. 2012, 134, 19043-19049. (107) Mir, A.; Chen, J.; Robinson, K.; Lendy, E.; Goodman, J.; Neau, D.; Golden, B. L. Two Divalent Metal Ions and Conformational Changes Play Roles in the Hammerhead Ribozyme Cleavage Reaction. Biochemistry-Us 2015, 54, 6369-6381. (108) Weinberg, Z.; Kim, P. B.; Chen, T. H.; Li, S.; Harris, K. A.; Lunse, C. E.; Breaker, R. R. New classes of self-cleaving ribozymes revealed by comparative genomics analysis. Nat. Chem. Biol. 2015, 11, 606-610. (109) Nguyen, L. A.; Wang, J.; Steitz, T. A. Crystal structure of Pistol, a class of self-cleaving ribozyme. Proc. Natl. Acad. Sci. USA 2017, 114, 1021-1026. (110) Lilley, D. M. J. How RNA acts as a nuclease: some mechanistic comparisons in the nucleolytic ribozymes. Biochem. SocTrans. 2017, 45, 683-691. (111) Kapral, G. J.; Jain, S.; Noeske, J.; Doudna, J. A.; Richardson, D. C.; Richardson, J. S. New tools provide a second look at HDV ribozyme structure, dynamics and cleavage. Nucleic Acids Res. 2014, 42, 12833-12846. (112) Lott, W. B.; Pontius, B. W.; von Hippel, P. H. A two-metal ion mechanism operates in the hammerhead ribozyme-mediated cleavage of an RNA substrate. Proc. Natl. Acad. Sci. USA 1998, 95, 542-547. (113) Blount, K. F.; Uhlenbeck, O. C. The structure-function dilemma of the hammerhead ribozyme. Annu. Rev. Biophys. Biomol. Struct. 2005, 34, 415-440. (114) Leclerc, F. Hammerhead Ribozymes: True Metal or Nucleobase Catalysis? Where Is the Catalytic Power from? Molecules 2010, 15, 5389-5407. (115) Young, K. J.; Gill, F.; Grasby, J. A. Metal ions play a passive role in the hairpin ribozyme catalysed reaction. Nucleic Acids Res. 1997, 25, 3760-3766. (116) Hertel, K. J.; Herschlag, D.; Uhlenbeck, O. C. A kinetic and thermodynamic framework for the hammerhead ribozyme reaction. Biochemistry-Us 1994, 33, 3374-3385. (117) Shih, I. H.; Been, M. D. Kinetic scheme for intermolecular RNA cleavage by a ribozyme derived from hepatitis delta virus RNA. Biochemistry-Us 2000, 39, 9055-9066. (118) McCarthy, T. J.; Plog, M. A.; Floy, S. A.; Jansen, J. A.; Soukup, J. K.; Soukup, G. A. Ligand requirements for glmS ribozyme self-cleavage. Chem. Biol. 2005, 12, 1221-1226.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 31

(119) Wilson, T. J.; McLeod, A. C.; Lilley, D. M. J. A guanine nucleobase important for catalysis by the VS ribozyme. Embo J. 2007, 26, 2489-2500. (120) Gkionis, K.; Kruse, H.; Sponer, J. Derivation of Reliable Geometries in QM Calculations of DNA Structures: Explicit Solvent QM/MM and Restrained Implicit Solvent QM Optimizations of G-Quadruplexes. J. Chem. Theory Comput. 2016, 12, 2000-2016.


Page 27 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


TOC Graphic 82x44mm (300 x 300 DPI)

ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 1: All the possible reaction scenarios for the internal transesterification. Reaction scenarios shown on yellow and white backgrounds correspond to dianionic and monoanionic mechanisms, respectively. Rare protonation states of phosphates and phosphoranes at pH ~7 (shown on a grey background) are assumed not to contribute to the catalysis due to enough high free-energy penalty associated with their rare protonation state that is not compensated by enhanced reactivity, and thus they were not considered in this study. In total, 115 distinct pathways were calculated. 209x296mm (300 x 300 DPI)


Page 28 of 31

Page 29 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 2: In-line attack angle (A and C panels, in deg.) and potential energy calculated on BLYP/631G(d,p)/CPCM level of theory (B and D panels, in kcal/mol) as a function of ε/ζ dihedral angles for structures with C2’-endo and C3’-endo puckering. The particular ε/ζ and pucker of the reactive conformations of the selected ribozymes41, 59, 85 and unassisted reaction paths probed in this study (scenarios I and V) are depicted by black dots in the right panels. The R state structures of the unassisted reaction paths are depicted on panel F. The energy maxima highlighted by stars on panels B and D correspond to structures with steric clashes of O2’ oxygen with one of non-bridging phosphate oxygens. 89x120mm (300 x 300 DPI)



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 3: Structures of R, rate-limitting TS, and P states along most feasible pathways (with the lowest ∆G‡300K free energy barrier) in four dianionic scenarios I to IV (see Figure 1). Free energies in kcal/mol of each state with respect to ground state of reactants, i.e., including free energy correction due to rare protonation form, are depicted in brackets. 192x170mm (300 x 300 DPI)


Page 30 of 31

Page 31 of 31

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 4: Structures of R, rate-limitting TS, and P states along most feasible pathways (with the lowest ∆G‡300K free energy barrier) in three monoanionic scenarios V to VII (see Figure 1). Free energy in kcal/mol of each state with respect to ground state of reactants, i.e., including free energy correction due to rare protonation form is depicted in brackets. 192x135mm (300 x 300 DPI)


Mapping the Chemical Space of the RNA Cleavage and Its

Recommend Documents