From Hydrogen Peroxide-Responsive Boronated Nucleosides Towards Antisense Therapeutics – A Computational Mechanistic Study

We used a combination of MD simulations and DFT calculations to reveal the precise chemical mechanism underlying the conversion of boronated nucleosides to natural nucleosides in the presence of hydrogen peroxide, which was recently experimentally demonstrated by Morihiro and Obika et al. (Chem. Sci. 2018, 9, 1112). Our results show that this process is initiated by the H2O2 deprotonation to a base concerted with the nucleophilic attack of the resulting OOH– anion onto the boron atom as the rate-limiting step of the overall transformation. This liberates a free base, followed by the 1,2-rearrangement to the C–OOH– adduct. Lastly, breaking of the O–O bond within the peroxide moiety cleaves the boron–carbon bond, giving boronic acid ester and the matching ketone as the final products. The obtained reaction profiles successfully interpret a much higher conversion rate of the thymine derivative over its guanine analogue, and rationalize why t-Buhydroperoxide is hindering the conversion, thus placing both aspects in firm agreement with experiments. The offered insight represents a promising tool for the future synthetic approaches of stimuli-responsive biomolecules, especially chemically caged prodrug-type nucleic acid therapeutics, bearing significant importance due to their application potential in diagnostics and therapy of various genetic disorders.


INTRODUCTION
TIMULI-responsive biomolecules have attracted growing attentions because of the possibility to control their functions and processes in biological environments, which opens the door for their potential application in the fields of chemical and synthetic biology. As a typical example, photo-caged systems, which feature a photo-removable chemical group bound to a target biomolecule, have been explored for decades. [1] These are based on temporarily masking the activity of an investigated biomolecule, which can be restored by an external stimulus, such as photo-irradiation, to remove the introduced chemical groups (decaging). Moreover, recent progress and deeper understanding of bio-orthogonal chemical reactions have extended the scope of the stimuli towards other means of the bond cleavage reactions, [2] including Diels-Alder reaction-triggered decaging to release active proteins, [3] a phosphine-mediated Staudinger reduction capable of activating proteins, [4] or a hydrogen peroxide (H2O2) mediated decaging of proteins with phenylboronic acid derivatives. [5] Considering these seminal works on proteins and their progress, chemically caged nucleic acids are not yet sufficiently explored.
Deoxyribonucleic acid (DNA) is widely regarded as the cornerstone of biological heredity and a conveyor of genetic information. [6] The highly sensitive detection of distinct genetic sequences plays a pivotal role in gene therapy, medical diagnostics, food safety, biodefense applications, and environmental monitoring. [7] Highly sensitive DNA probes may thus be used as an effective tool for the early stage detection of genetic disorders, such as cancer, [8] where both the sheer volume of the number of reported cases and the huge death toll indicate the importance of the development of efficient strategies for the early detection of these diseases to safeguard human health and further the advancement of allied research areas. [9] With this in mind, it comes as a no surprise that several antisense oligonucleotides (ASOs) have been under clinical development for the last few decades due to their promising therapeutic and diagnostic applications. [10] ASOs hybridize with complementary mRNA to regulate specific gene expression either via RNase H-induced strand scission, translational arrest caused by the steric blockade of the factors involved in translation, splicing and metabolism, or via splice switching. [10,11] Several antisense drugs have been approved by pharmaceutical regulatory agencies in recent years. [10,12] Although ASOs offer great potential to be a simple approach for the rational design of effective therapeutics, their poor nuclease stability, binding affinity, physicochemical and pharmacokinetic properties are often cited as barriers to antisense drug development. [13] In order to improve in vivo pharmacokinetics, regulate the activity in a desirable manner and reduce the off-target effects of ASOs, using a desired stimuli-responsive moieties could be a valuable tool.
Small molecules can trigger the selective decaging to restore active biomolecules even in living cells. [2b,14] In addition, the chemical environment inside cells, which varies depending on the physiological condition, can be employed to alter the properties of nucleic acids under a specific biomarker. It is known that while normal cells show very low reactive oxygen species (ROS) levels, cancer cells and inflamed cells exhibit increased levels of intracellular ROS including hydrogen peroxide, hydroxyl radical and superoxide anion. [15] H2O2 -a prominent ROS with a crucial role in carcinogenesis and cancer development [16] -is present in significant concentrations in vivo, and is thus an ideal candidate as a biomarker of tumors. The H2O2 levels in cancer cells are found to be in the range of 1 -5 mM, significantly higher than those in normal cells (≈ 0.5 μM -1 mM). [17] These biological features of cancer and inflamed cells provide an opportunity for the selective prodrug design. [18] In synthetic organic chemistry, H2O2 is a good uncharged nucleophile and can be used as a two-electron electrophilic oxidant. [18c] Moreover, it is well established that alkyl-or aryl boronic acids and their esters can be easily dissociated by H2O2 to the products which do not show any cytotoxicity, suggesting that the reaction between boronate moiety and H2O2 is chemospecific, bioorthogonal and biocompatible, [18] on top of the fact that many boroncontaining compounds are already in the clinical use. [19] Along these lines, Morihiro and Obika et al. [18b] very recently prepared a series of boronated nucleosides ( Figure   1A) and introduced O4-(4-boronobenzyl)-modified thymine into ASOs for the H2O2-triggered gene silencing in mammalian cells. Namely, they masked the Watson-Crick face of the nucleobase with pinacol borane that sterically hindered duplex formation between the ASO and target mRNA. In the presence of H2O2, the boronate group was swiftly removed, resulting in gene silencing through hybridization of the ASO to the target mRNA. Although the authors proposed that phenylboronic acid pinacol ester is first oxidized to phenol derivative by H2O2 and subsequent 1,6-elimination releases the active molecule ( Figure 1B), the mechanism of the activation is not fully clarified. In this work, we used a combination of molecular dynamics (MD) simulations and density functional theory (DFT) calculations to provide deeper insights into the H2O2-triggered ASO activation and reveal that its precise chemical mechanism follows a different reaction sequence. Specifically, the reaction is initiated by the H2O2 deprotonation to a base concerted with the nucleophilic attack of the resulting HOOonto the boron atom, being the rate-limiting process of the overall transformation. This liberates a free base, while the rest of the system undergoes an internal rearrangement from the boron-bound HOOto a neighboring carbon-bound HOOderivative. The following breaking of the O-O bond within the peroxide cleaves the carbon-boron bond, yielding boronic acid pinacol ester and the corresponding ketone as the final products. The obtained insight represents a promising tool for future synthetic approaches of prodrug-type nucleic acid therapeutics for selective disease treatment.

Molecular Dynamics Simulations
Molecular dynamics simulations were performed in order to sample the conformational space of the reacting complexes and elucidate representative structures for the subsequent density functional theory analysis. For all reactants, geometry optimization and RESP charge calculations were performed with the Gaussian 16 program [23] at the HF/6-31G(d) level to be consistent with the employed GAFF force field. The formed reaction complexes were solvated in a truncated octahedral box of TIP3P water molecules spanning a 10 Å thick buffer and submitted to geometry optimization in AMBER16 program, [24] employing periodic boundary conditions in all directions. Optimized systems were gradually heated from 0 to 300 K and equilibrated during 30 ps using NVT conditions. This was followed by productive and unconstrained MD simulations of 300 ns employing a time step of 2 fs at a constant pressure (1 atm) and temperature (300 K), the latter held constant using Langevin thermostat with a collision frequency of 1 ps -1 . The long-range electrostatic interactions were calculated employing the Particle Mesh Ewald method, [25] while the nonbonded interactions were truncated at 10.0 Å. Following MD simulations, the obtained structures were clustered according to the distance between the reactive oxygen atom on ROS and the boron atom on arylboronatemodified nucleosides on every tenth structure from the last 100 ns of simulations, and representative snapshots submitted to mechanistic DFT calculations.

Density Functional Theory Calculations
In order to minimize errors associated with the initial selection of starting geometries from MD trajectories, we tried several conformations of each reactive complex and proceeded with the mechanistic calculations using the most stable complexes. As a good compromise between accuracy and feasibility, all geometries were optimized in the gas-phase by the very efficient M06-2X/6-31+G(d) method with thermal Gibbs free energy corrections extracted from the corresponding frequency calculations without the scaling factors. To account for the solvent effects, the final gasphase energies were corrected for the Gibbs solvation energies employing the SMD polarizable continuum model with all parameters corresponding to pure water, giving rise to (SMD)/M06-2X/6-31+G(d)//M06-2X/6-31+G(d) model employed here. The choice of this computational setup was prompted by its success in reproducing kinetic and thermodynamic parameters of various organic, organometallic and enzymatic reactions, [26] being particularly accurate for relative trends among similar reactants, which is the focus here. It is clear that such a methodology could be further improved by considering explicit solvation through any of the established QM/MM or EVB techniques, [27] yet this is beyond the scope of current work, especially since experimental kinetic and thermodynamic parameters for the investigated processes are not known. All of the transition state structures were verified to have the appropriate imaginary frequency from which the corresponding reactants and products were determined using the Intrinsic Reaction Coordinate (IRC) procedure. [28] Atomic charges were obtained by natural bond orbital (NBO) analyses [29] as the single-point calculations at the (SMD)/M06-2X/6-31+G(d) level. All of the calculations were performed using the Gaussian 16 software. [20]

RESULTS AND DISCUSSION
Structures of the investigated systems are depicted in Figure 1A and these involve boronated thymine (1) and guanine (2), as well as hydrogen peroxide (H2O2) and tbutyl-hydroperoxide (TBHP) as reactive oxygen species. For the sake of simplicity, both nucleosides were truncated at the N-glycosidic position, where the sugar moiety was replaced by the N-methyl group. This is a reasonable approximation, since this position is very distant from the reactive boron-carbon bond in investigated nucleosides, and this modification should not have a significant impact on the obtained results. In this context, although our calculations were done on modified nitrogen bases, the obtained insight and conclusions are closely applicable to actual nucleosides as well, which is a terminology that we will keep throughout the text. The choice of these nucleosides was prompted by the fact that Morihiro and Obika et al. demonstrated that, out of all four boronated systems, thymine derivative is most efficiently converted to the corresponding natural nucleoside, while guanine undergoes the least efficient conversion, both with H 2O2.
[18b] In addition, thymine is a pyrimidine, while guanine belongs to a purine class of nucleobases, which justifies their separate analysis. Lastly, both nucleosides were shown to be practically inactive towards TBHP, which deemed it worthwhile to investigate and rationalize such a dramatic influence of the t-Bu substituent effect on the reactivity. Approach of H2O2 towards nucleosides generates reactive complexes in which one of the O-H groups forms a hydrogen bond with a base fragment, while the other is oriented in a way to stack above the aromatic phenyl ring (SP1, Figure 2). In 1, this hydrogen bond is established with the carbonyl oxygen on thymine, with the corresponding O•••O distance of 2.785 Å (SP1(1)). On the other hand, in 2, the latter is formed with the imino nitrogen in the fivemembered ring and the matching O•••N distance is slightly longer at 2.889 Å (SP1 (2)). In both cases, this interaction  turns out to be crucial for the reactivity, as it enables the deprotonation of H2O2 to a base and the subsequent nucleophilic reactivity of the formed perhydroxyl HOOanion. The latter is accomplished onto the electrophilic site in nucleosides, which is a region around the B-C bond ( Figure  S1), which is exactly what needs to be cleaved to afford the final products. [18b] Without this deprotonation, the unionized H2O2 appears not to have enough nucleophilicity to perform the reaction, as all of our attempts to locate transition states (or even stable products) for its attack on the B atom or its neighboring aromatic C atom failed. Nevertheless, this is fully in line with the fact that, in organic chemistry, a typical H2O2-mediated oxidation is performed in the alkaline media, [20] where OHperforms the H2O2 deprotonation (pKa = 11.54) [21] to give a nucleophilic HOOanion.
With this in mind, we first investigated the approach of H2O2 towards the B atom in both 1 and 2 ( Figure 3). In the transition structure TSSP1-SP2 there is a simultaneous formation of the new B-O bond accompanied with the peroxide deprotonation to a base ( Figure 2). This reaction is feasible for both systems, yet the intrinsic barrier for 2 is by 0.9 kcal mol -1 lower. This is easily explained knowing that guanine is a much stronger base than thymine, [22] therefore deprotonation to guanine is easier. However, guanine is a larger and more polar system, thus its solvation energy is higher, which makes it is less prone to oxidation. To firmly support that, our calculations for 1 and 2 gave the watersolution solvation energies of -15.4 and -20.9 kcal mol -1 , respectively. As a consequence, it is by 2.4 kcal mol -1 less favorable to form a reactive complex with 2 than it is with 1, which increases the total activation free energy for 2 to 33.2 kcal mol -1 , being 1.5 kcal mol -1 higher than that for 1 ( Figure 3). This reaction gives a tetrahedral B-OOHadduct SP2, and is slightly more favorable for guanine 2, although, for both nucleosides, the reaction is significantly endergonic.
Following SP2, we investigated several mechanistic possibilities. A direct cleavage of the B-C bond, as proposed by Morihiro and Obika et al. to potentially give a phenol derivative, [18b] does not offer good products. Namely, in this case, the mentioned C atom picks up hydrogen from the OOH group to form a phenyl moiety, which terminates the reaction. Analogously, the cleavage of the O-O bond within the peroxide and the formation of the OH adduct on the aromatic carbon atom is linked with a further barrier of 28.8 kcal mol -1 for 1, which turns out to be one of the less feasible pathways, and later does not offer a viable route towards the final products. Interestingly, the formation of SP2 weakens the O-C bond connecting a base and a benzyl moiety, as evidenced in its elongation from 1.441 and 1.438 Å in SP1 for 1 and 2, respectively, to 1.510 and 1.469 Å in SP2, in the same order. This allows a very efficient release of a base, having a small barrier of 7.9 kcal mol -1 for 1 and 17.3 kcal mol -1 for 2, to give a non-covalent complex between a base and the rest of the system (SP3), in which a base spontaneously tautomerizes to a more stable physiological form, as in 4 and 5. It is worth noting that a significant difference in the height of these barriers in TSSP2-SP3 comes as a result of different extents of the O-C bond elongation, and is lower for a much longer O-C bond in SP2(1). A free base, then, leaves the complex, being favorable for both 1 (by 5.2 kcal mol -1 ) and 2 (by 6.4 kcal mol -1 ), which allows it to undergo specific antisense pairing with the target DNA, thus forming basis for its potential therapeutic effect.
At this point on a reaction pathway, several observations are in place. Release of a free base leaves an intermediate SP4, which is high in energy (Figure 3), leading to a conclusion that it is not the end point of the whole process. Also, further steps, which should enable the thermodynamic feasibility of the overall conversion, no longer depend on the initial substrate, as SP4 will be the same for all analogously boronated nucleosides. Yet, our calculations clearly show that the formation of SP4 is significantly more favorable for thymine 1 than it is for the guanidine analogue 2. Namely, the overall barrier leading to SP4 is ΔG ‡ = 33.9 kcal mol -1 for 1, being 2.5 kcal mol -1 lower than for 2. Analogously, the reaction Gibbs free-energy until this point is ΔGR = 22.7 kcal mol -1 for 1, which is 5.3 kcal mol -1 more favorable than for 2 ( Figure 3). Therefore, both kinetic and thermodynamic aspects strongly favor the transformation of 1 over 2. This is found in excellent agreement with experiments by Morihiro and Obika et al., [18b] and nicely explains why thymine analogue 1 is much favorably converted to a free thymine 4 than 2 is to 5. The underlying reasons for that are several: (1) guanidine is larger and is better solvated in the solution, making it less prone to forming a reactive complex with H2O2 and undergoing the subsequent oxidation, and (2) its size allows the guanidine fragment to more effectively distribute an excess positive charge formed upon its protonation in SP2, which leaves less positive charge around the bond connecting guanidine with the rest of the nucleoside, making its release less effective. The later is evidenced in the NBO charges at the corresponding oxygen atom in SP1 being -0.55 |e| for both 1 and 2, respectively, which are changed to -0.58 and -0.60 |e| in SP2, in the same order.
SP4 intermediate features OOHbound to the boron atom. From that, we considered several routes, but the only feasible pathway involved the cleavage of the B-OOHbond and the 1,2-rearrangement to the C-OOHadduct SP5. This is facilitated by the fact that the latter is by 5.1 kcal mol -1 more stable, while the barrier to achieve it is moderate at 8.7 kcal mol -1 . In the last step there is a cleavage of the O-O bond in the peroxide fragment accompanied by a simultaneous boron-carbon cleavage and the B-O bond forming reaction (TSSP5-SP6, Figure 2). The barrier for this process is 23.5 kcal mol -1 , and the reaction is highly exergonic to give SP6 involving a complex of boronic acid pinacol ester 6 and the matching ketone 7. Dissociation of both products is thermodynamically further favorable, giving the overall reaction free-energies of ΔGR = -91.1 kcal mol -1 for 1 and ΔGR = -85.8 kcal mol -1 for 2. It is worth noting that, on top of its demonstrated kinetic feasibility, the calculated thermodynamic aspects also support the conversion of 1 as more favorable, in line with experiment. [18b] As an alternative, following the reactive complex SP1, the reaction could potentially proceed through the nucleophilic attack of H2O2 onto the aromatic carbon atom, since this site is also significantly electrophilic ( Figure S1). This is an interesting pathway, since this process spontaneously liberates a base and immediately positions the resulting OOHanion onto the aromatic carbon atom (SP2', Figure 4). Consequently, this sequence of events reduces the number of steps for the already demonstrated feasible pathway connecting SP1 to SP5 (Figure 3), which initially required three steps, into a single reaction. However the barrier for this process is very high, being ΔG ‡ = 56.0 kcal mol -1 for 1 and ΔG ‡ = 53.0 kcal mol -1 for 2, which is much higher than those calculated for the previous pathway. All of this renders the mechanism, which is initiated through the addition of H2O2 onto the aromatic carbon atom, as very implausible.
Unlike H2O2, which is a symmetrical molecule, its tert-butyl analogue TBHP is an unsymmetrical system, with the oxygen atom bearing the t-Bu group being considerably more nucleophilic than its O-H counterpart ( Figure S1). Nevertheless, its initial reactive complex is identical to that of H2O2 in a way that a stable hydrogen bonding with a base is formed (SP1, Figure 2 However, although our attempts to model the analogous approach of the O(t-Bu) atom in TBHP on either boron or carbon atom in 1 and 2 gave stable products, such a reaction proceeds by cleaving the O-O bond in TBHP followed by the departure of the reactive hydroxyl anion OH -. The latter spontaneously binds to different positions on each nucleoside, thus preventing further conversion, and not leading to final products. As an illustration, the attack of the O(t-Bu) atom on the aromatic carbon in 1 is associated with a barrier of ΔG ‡ = 89.4 kcal mol -1 and the produced OHforms a chemical bond with the neighboring ortho-carbon on the phenyl ring, thus is clearly not feasible and unproductive. On the other hand, a scenario analogous to that depicted in Figure 3 would involve the attack of the less nucleophilic O(H) moiety on the boron atom. Still, during this approach, TBHP rotates in a way to lose the hydrogen bonding with a base for the hydrogen bonding with one of the O-atom on the boronic ester (SP1'(2), Figure 2), thus giving the unproductive orientation. This is a stable structure, which is facilitated by the favorable C-H•••••π interactions that t-Bu hydrogen atoms form with aromatic moieties on both the phenyl group and a base. Still, in this orientation, the attack of the O(H) from TBHP onto the Batom is feasible and the barrier for that is reasonable at ΔG ‡ = 34.4 kcal mol -1 , with ΔGR = 13.2 kcal mol -1 . Nevertheless, the product of that reaction is HOObound to the B-atom with a proton residing on the O(B) moiety. This prevents a favorable release of a base, and all other mechanistic options failed to give the right products. Alternatively, forcing TBHP to approach the electrophilic B-atom in an orientation that maintains the hydrogen bonding to a base, allows a mechanism analogous to that in Figure 3. However, the barrier for the transition SP1→SP2 with 1 is ΔG ‡ = 41.2 kcal mol -1 , being 9.5 kcal mol -1 higher than in a reaction with H2O2. This is found in excellent agreement with experiments, [18b] which revealed that the selectivity of 1 towards H2O2 is around 16-times higher than towards TBHP. All of this lends strong credence to our results and supports the proposed mechanism. Moreover, this kinetic aspect is also combined with a better thermodynamic feasibility of the whole transformation, which for H2O2 is ΔGR = -91.1 kcal mol -1 (Figure 3), while for TBHP it is ΔGR = -83.4 kcal mol -1 . Lastly, although Morihiro and Obika et al., [18b] did not investigate the reactivity of TBHP with boronated guanidine 2, our calculations show that this process would be even less favorable than that for the mentioned reaction between H2O2 and 2. Specifically, the first transition state TSSP1-SP2 would be 44.0 kcal mol -1 above individual reactants, 2.8 kcal mol -1 higher than with H2O2, while the reaction Gibbs free-energy would be the least favorable of all examined here at ΔGR = -78.1 kcal mol -1 . Along these lines, considering the mechanism presented in Figure 3, we can reasonably speculate that very low conversion rates of 1 with other reactive oxygen species examined experimentally, [18b] such as OCl • , HO • , tBuO • , O2 • and NO, are likely attributed to (i) the non-existence of the acidic hydrogen in these systems, which would protonate the base and allows its effective release, and also to (ii) a well-demonstrated non-selective high reactivity of these open-shell systems, [15] which is certainly not limited to the most electrophilic regions around the B-C bond in 1-2. Any process other than that would give non-productive conversions of such and similar boronated nucleosides, and would disfavor the appropriate antisense pairing with target DNA molecules, thus is of no interest here.

CONCLUSIONS
In this work we used a combination of molecular dynamics simulations and density functional theory calculations to reveal the precise molecular mechanism governing the hydrogen peroxide-triggered conversion of boronated nucleosides into natural nucleosides, which was very recently experimentally demonstrated by Morihiro and Obika et al. [18b] By investigating several mechanistic possibilities, we showed that this process involves several steps and is initiated by the nucleophilic attack of H2O2 onto the boron atom of the substrate, where the resulting B-O bond formation is concerted with the hydrogen peroxide deprotonation to a base fragment. This is followed by the release of a free base, and the 1,2-rearrangement of the formed B-OOHadduct into the C-OOHadduct. Lastly, the latter undergoes the O-O bond cleavage within the peroxide moiety, which gives boronic acid pinacol ester and the corresponding ketone as the final products. In this way our results contradict those by Morihiro and Obika et al., [18b] which suggested that initial nucleosides are first oxidized to phenol derivatives and that a subsequent 1,6-elimination releases the active molecules.
The obtained kinetic and thermodynamic parameters are consistent in confirming that the H2O2-mediated conversion of the thymine derivative 1 is more favorable than that of the guanine analogue 2, as its both activation and reaction free energies are by ΔΔG ‡ = 2.5 kcal mol -1 and ΔΔGR = 5.3 kcal mol -1 , respectively, more favorable than those calculated for 2. In this way, our results are found in excellent agreement with those reported experimentally, [18b] which revealed the same trend among nucleosides. In addition, we showed that the reaction with tBuderivative TBHP is hindered by the unproductive hydrogen bonding of its O-H moiety with the oxygen atom of the pinacol fragment, which does not allow the protonation of the base moiety and its efficient release. Nevertheless, the overall conversion of 1 with TBHP is still possible, yet the calculated barrier is 9.5 kcal mol -1 higher than in the reaction with H2O2. This also nicely ties in with Morihiro and Obika et al., [18b] which revealed that the selectivity of 1 towards H2O2 is around 16-times higher than towards TBHP. Taken all together, all of this lends credence to our results and strongly supports the proposed mechanism.
Given the growing attention and significant progress in the development, deeper understanding and widespread application of stimuli-responsive proteins, we hope that presented results will aid in directing the attention of researchers towards reaching the same level of knowledge and sophistication with chemically caged nucleic acids. These may likely find their use in the early stage detection of genetic disorders, such as cancer, which is becoming of utmost importance and urgency from both the sheer volume of the number of reported cases of genetic disorders and the huge death toll.