Effect of arsenate substitution on phosphate repository of cell: a computational study

The structural analogy with phosphate derives arsenate into various metabolic processes associated with phosphate inside the organisms. But it is difficult to evaluate the effect of arsenate substitution on the stability of individual biological phosphate species, which span from a simpler monoester form like pyrophosphate to a more complex phosphodiester variant like DNA. In this study, we have classified the physiological phosphate esters into three different classes on the basis of their structural differences. This classification has helped us to present a concise theoretical study on the kinetic stability of phosphate analogue species of arsenate against hydrolysis. All the calculations have been carried out using QM/MM methods of our Own N-layer Integrated molecular Orbital molecular Mechanics (ONIOM). For quantum mechanical region, we have used M06-2X density functional with 6-31+G(2d,2p) basis set and for molecular mechanics we have used the AMBER force field. The calculated rate constants for hydrolysis show that none of the phosphate analogue species of arsenate has a reasonable stability against hydrolysis.


Introduction
Maintaining negative charge at physiological pH leads phosphate (PO 2À 4 or P i ) ester to dominate the living world as it prevents chemically stored energy and genetic material from escaping the cell [1]. (V) is also capable of retaining negative charge in its phosphate analogue arsenate (AsO 2À 4 or As i ) over a range of physiological pH conditions [2] and it is already demonstrated in bacteria that at a high concentration ratio of As i /P i the phosphate transporters allow the transport of As i on the basis of structural similarity [3]. This structural analogy further leads to the integration of As i into various metabolic processes associated with P i inside the organisms [ concern related to the P i analogue species of As i (PASA) is their rapid rate of hydrolysis [1]. However, the scientific world is really surprised by the discovery of arsenate substituted DNA backbone in place of phosphate (As i -DNA) in a bacterium strain GFAJ-1 of the Halomonadaceae family, which was isolated from Mono Lake, California with a very high As i /P i ratio (greater than 103) [4]. The theoretical investigations on the possibility of As i -DNA also revealed that both arsenodiester and phosphodiester linkages have similar geometric and conformational properties [5,6]. Although these theoretical studies support geometrical similarities between As i -DNA and P i -DNA, neither of them claims the stability of As i -DNA against hydrolysis. Whereas the models described by Fekry et al. and Mlá dek et al. to estimate the kinetic stability of phosphodiester and arsenodiester linkages in DNA lack the fundamental base stacking interactions, also called p 2 p stacking [6,7]. In order to evaluate the stability of a molecule like DNA, a more realistic model is required because the work of Wang et al. reveals that the stacking between neighbouring bases in DNA single strands raises the activation energy barrier for hydrolysis of DNA [8]. Along with this, it is found that the discovery of As i -DNA in GFAJ-1 captured the attention of the scientific community from other PASA. However, it is not clear how As i could be directly incorporated into P i -DNA because neither As i nor P i are the substrates for the DNA synthesizing enzymes and neither of these are directly incorporated into the diester backbone of the DNA. There are various metabolic steps through which monosaccharides are processed using P i biomolecules and metabolized into nucleoside triphosphate (NTP's). NTP's are the monomeric units of nucleic acid polymers (DNA and RNA), and P i is able to form the phosphodiester backbone of the DNA through the nucleophilic substitution reaction between deoxy-NTP molecules. Similarly, before reaching DNA, As i must have a reasonable stability in its dNTP analogues, whereas biosynthesis of As i analogue of dNTP again depends upon the availability of some other PASA like Ribose-5-As i , Glucose-1-As i , Pyro-As i , etc. [9]. Physiologically, these are dianionic monoesters of P i and highly sensitive to As i substitution [10]. Hence, to investigate the effect of As i substitution on the P i repository of the cell a comparative study at different biomolecular level is needed, that is not merely confined to analysing the impact of As i substitution on the stability of a particular P i biomolecule. However, due to a vast number of P i species inside organisms, it is difficult to compare the rate of hydrolysis for each and every PASA against their respective P i analogue. In this study, through a unique classification of P i biomolecules, we examined the consequence of As i substitution on the kinetic stability of P i species spanning from a simpler monoester form like pyrophosphate to a more complex phosphodiester variant like DNA.

Hydrolysis pathways
Pathways of the reaction for the hydrolysis of As i and P i are well documented in the literature [11 -13]. The associative pathway is suggested to be dominant for both the esters [11] and follows the S N 2 mechanism deformed reactants TS water As/P-ester Figure 1. The initial step of associative pathway for the hydrolysis of (As/P) esters. Inset showing major geometrical distortion (in red) observed in ester geometry. (E ¼ As/P, L ¼ leaving group, O p ¼ oxygen atom of ester capturing proton from water, rsos.royalsocietypublishing.org R. Soc. open sci. 5: 181565 with first-order kinetics. The initial step of the pathway as shown in figure 1 involves the attack of water as a nucleophile on the central atom (As/P) of the ester. It is found that the attack of water has the highest barrier in terms of activation energy and hence determines the rate of hydrolysis of the ester [6,11,12]. In addition to the nature of reactant and attacking nucleophile, the chemistry of the leaving group is also the key factor for the kinetics of nucleophilic substitution reactions. However, the complete mechanism of the reaction for the nucleophilic attack on different types of phosphate ester was investigated in detail by Lopez et al. who concluded that the rate determining step is always that for the attack of nucleophile rather than the departure of the leaving group [12]. In other words, the activation energy barrier for nucleophilic substitution is mainly influenced by attacking nucleophile and steric congestion around a central atom [14,15]. We employed the activation strain model [14] to relate the height of the activation energy barrier with the geometrical deformation and rigidity of reactants (figure 1). Our study aims at comparing the kinetic stability of As and P esters against water, hence we only modelled the reactants and transition state (TS) structures which belong to the rate determining step. In an attempt to understand the TS, it was found that after reaching the TS for monoanionic monoesters and diesters of P form a reaction intermediate, having a pentacoordinated centre with trigonal bipyramidal geometry [12]. However, direct product is observed in the case of dianionic monoester of P [12]. We also observed that the TS is followed by a pentacovalent intermediate, during the hydrolysis of all classes of PASA, irrespective of the type of ester and anionic form. The reaction via pentacovalent intermediate is further proceeded by internal proton transfer that leads to the breaking of the As/P-O bond [6,11,12]. The structures of reactant and TS as shown in figure 1 for all three classes were modelled (available with the electronic supplementary material). We carried out intrinsic reaction coordinate calculations (movie files available with the electronic supplementary material) to confirm coordinates of the reactant and intermediate for the initial step as shown in figure 1.

Classification
In order to determine the effect of As i substitution on the P i repository of the cell, we first need to understand the physiological functions of P i species. P i has a variety of roles inside a cell, and is incorporated into biomolecules mainly through kinases. This process is also known as phosphorylation and the reverse of this is de-phosphorylation, removal of P i by hydrolases. This cycle of addition and removal of P i governs biological phenomena such as activation or deactivation of proteins, synthesis of DNA and removing falsely paired nucleotide bases from it, breaking down hexoses into trioses, synthesis and utilization of ATP, etc. However, all of these happen at specific sites in a well-controlled way. Phosphorylation/de-phosphorylation on proteins usually takes place on the amino acids that have a hydroxyl group on their side chains, e.g. serine, threonine and tyrosine. Similarly, on hexoses the hydroxyl group at terminal carbon atoms is the site of phosphorylation/de-phosphorylation. Whereas in ADP, the hydroxyl group attached to the terminal beta phosphorus (Pb) atom is used to synthesise ATP [9]. These specific sites are conserved throughout the organisms, and the rate of phosphorylation/de-phosphorylation at these sites is well synchronized with each biological process.
We classified physiological P i species based on the type of ester and conserved sites through which P i anchors with biomolecules (table 1). The same sites also behaved as leaving groups during dephosphorylation or hydrolysis of P i ester. The modelled compound from each class is shown in figure 2. Our classification broadly covers most of the P i biomolecules, Class 1A representing the energy currency of the cell, like ATP and ADP, while examples of Class 1B members are sugar or glycerol moiety containing biomolecules like Glucose-1-P i , Fructose-1,6-Bi-P i , Phosphoglycerate, etc. Whereas the molecules capable of storing genetic information belong to Class 2, like DNA and RNA. This classification not only enabled us to determine the relative stability of PASAs against hydrolysis with their respective P i analogue, it also  helped us to understand the consequence of As i substitution at two different levels of biological esters (monoester and diester) and the role of conserved sites for phosphorylation/de-phosphorylation in the kinetics of hydrolysis.

Methods
For quantum chemical calculations, we need a suitable method to describe the thermochemical kinetics accurately for hydrolysis of inorganic systems like pyro-P i /As i and for bigger biomolecular systems like DNA that have non-covalent interactions. Zhao and Truhlars' work recommended Minnesota density functionals for thermochemistry, thermochemical kinetics and non-covalent interactions [16]. The presence of a polar solvent (1 ¼ 78.4, water) was mimicked with the polarizable continuum model using the integral equation formalism variant [17]. We carried out computations using M06-2X Minnesota density functional with a 6-31þG(2d,2p) basis set for the hydrolysis of pyroarsenate, and calculated the rate of the hydrolysis at 298 K using equation (2.1) where c8 ¼ 1 and DG y represents the free energy of activation. The calculated rate constant value of 0.08 s 21 is in good agreement with the available experimental value of 0.05 s 21 for the non-enzymatic hydrolysis of pyroarsenate. [18] 2.1. Modelling of P i -DNA and As i -DNA The crystal structure of B-DNA dodecamer on PDB (id: 1bna.pdb) was used to construct a single stranded P i -DNA model system consisting of d(TpTpApA) nucleotide bases. For modelling the As i -DNA d(TAsTAsApA) structure, we replaced the respective two phosphorus atoms of P i -DNA d(TpTpApA) with arsenic. It is difficult to deal with the whole system quantum mechanically (QM) due to the size complexity of P i -DNA and As i -DNA. Hence, QM was only applied on the portion which is involved in the reaction and the remaining part of the structure was dealt with using molecular mechanics (MM). We used two layers (our Own N-layer Integrated molecular Orbital molecular Mechanics) ONIOM [19], integrated into a Gaussian 09 [20] package, for QM/MM calculations. To avoid the problem of cancellation in ONIOM, a non-polar bond between two sp 3 hybridized C atoms, which is more than four bonded atoms away from the centre of the reaction, was The MM region was treated with AMBER [21] and frozen to maintain the structural constraints for diester backbone. For the QM region we chose M06-2X with 6-31þG(2d,2p) basis set which provides a good description of base stacking interaction for single-stranded DNA [22]. Our modelled P i -DNA with M06-2X 6-31þG(2d,2p)/AMBER optimization of the geometry shows a reasonable geometrical similarity with the identical portion of B-DNA crystal structure with a RMSD value of 1.54 Å as shown in figure 4.

Results and discussion
The oxygen atoms bonded to central atom (As/P) of the esters are labelled (figure 2) to discuss and compare the parameters of the geometry between modelled reactants and TS compounds. Comparison of the main geometrical data between computed arsenate and phosphate models belong to Class 1A, 1B and 2 is presented in tables 2 -4, respectively. Although the tabular data shows differences among the selected bonded parameters of the reactant and TS structures, these parameters were rarely used to relate the two different ester species. However, using the activation strain model we represented these geometrical changes in terms of the following energy terms (DE # dist ¼ distortion energy, DE # int ¼ interaction energy, DE # act ¼ activation energy), were basically derived from single point energy calculations [14]. The calculated values for distortion, interaction and activation energies for each species of three classes were shown in figure 5. Throughout the classification, irrespective of the anionic form, geometry of P i esters were more deformed in TS compared to the As i esters, which is   figure 5. The higher deformity in TS for P i esters might be due to the smaller size of the central P atom, compared to the bigger size of the As atom which allows more space for the incoming water nucleophile. Similarly, we also found that the geometry of dianionic P i esters was more deformed than its monoionic form because two negatively charged oxygen atoms on the smaller central P atom possess a higher steric barrier for incoming nucleophile. However, the effect of higher deformity in the geometry of dianionic P i monoester was not reflected in terms of the activation energy barrier due to the higher contribution from the interaction energy. In fact, it  is the interaction energy which lowers the activation energy barrier for the hydrolysis of all the dianionic As i as well as P i monoesters. Interestingly, in the case of diester of class 2, the activation energy barrier achieved a comparable height equal to the interaction energy and this signifies the importance of diester linkages inside a living system. Compared to As i -DNA a lot of energy is required to deform the structure of  Figure 5. Graph of distortion, interaction and activation energies for the rate determining step of reaction between As i /P i esters and water. (green: distortion energy, blue: distortion energy, red: activation energy, black: interaction energy, in kcal mol 21 ). P i -DNA during hydrolysis and this causes the biology to rely upon the smaller P atom for more rigid and stable geometry compared to the larger As atom. Apart from the comparison of geometrical deformity, the kinetics of hydrolysis for the esters is also discussed further in detail. We calculated the rate constant for the hydrolysis of different classes of As/P-esters as shown in table 5. Both Class 1A and 1B dianionic ester of P hydrolyzed much faster than the monoanionic form as reported earlier [12]. We observed the same trend in the case of As esters of Class 1A and Class 1B. Higher stability of monoanionic ester over dianionic ester has biological significance. For instance, a stable molecule, e.g. DNA, with a monoanionic phosphodiester backbone, preserves the genetic information. On the other hand, comparatively facile hydrolyzable Class 1A and 1B esters are available physiologically in dianionic form to provide energetic requirements and assist in several other anabolic processes for an organism [1]. In comparison to the P i counterpart of Class 1A and 1B, we found that none of the PASA have reasonable stability against hydrolysis. We calculated the rate for the hydrolysis of As i -DNA to be 1.23 Â 10 26 s 21 which is slower than other PASA like pyroarsente and ribose-1-arsenate, which have first-order rate constant of 0.08 s 21 and 0.07 s 21 , respectively, in dianionic form and 0.0003 s 21 and 0.0007 s 21 , respectively, in monoanionic form. Higher kinetic stability of As i -DNA against hydrolysis in comparison to other monoanionic esters of As i belonging to Class 1A and 1B, was expected because of steric hindrance and base stacking interactions, [8] but it is still a long way from competing with the stability of P i -DNA against hydrolysis which has a firstorder rate constant value of 5.74 Â 10 220 s 21 . Another aspect of selecting P i over As i is that P i is 10 6 times more stable in monoanionic diester form than in dianionic monoester variants of Class 1B, whereas this difference was found to be narrower in the case of As i . This difference might favour the integration of dianionic monoesters of P i into a relatively more stable biomolecule like DNA. We also found that there is no significant difference among the rate constants of monoanionic pyro-As i and ribose-1-As i , monoionic pyro-P i and ribose-1-P i , similar to the dianionic species of both the esters, which signifies that the conserved sites through which these two esters bonded covalently have nothing to do with the kinetics of the hydrolysis. However, these conserved sites have well-defined regulatory roles in protein -protein interactions [23].
We reported that the different classes of PASA have much less stability against hydrolysis as compared to their P i analogue species. Now, the rapid hydrolysis of PASA results in free As i which is again going to compete with P i . In order to avoid the structural confusion between As i and P i , the living organisms require metabolic machineries for the methylation of As i that results into organoarsenical compounds which have one or more As-C bonds instead of As-O bonds. Methylation not only resolves the structural ambiguity between As i and P i , it might also prevent conversion of As i into a highly reactive form which is known as Arsenite(III). On the other hand, organoarsenicals are highly stable against oxidation as well as against hydrolysis [24], as a result a number of marine organisms accept arsenolipids as a component of their cell membrane which are functionally, but not structurally, similar to phospholipids [25][26][27]. Arsenolipids primarily consist of As-C bonds which indicates the significant biological stability of As-C bond-based compounds as compared to kinetically unstable As-O bonds containing PASA.

Conclusion
Our study concluded that all classes of PASA have a higher hydrolysis rate compared to their respective P i -analogues. In biological systems, where the rate of hydrolysis of P i biomolecules is synchronously coupled with another biological process, any alteration in the kinetics of hydrolysis could have a severe effect on the metabolic activity of a cell. Although As i -DNA has a higher resistance against hydrolysis compared to PASA of Class 1A and 1B, it is still a kinetically unstable molecule, and organisms would face severe consequences by allowing such a highly unstable molecule to substitute P i -DNA for its genetic material. The rapid rate of the hydrolysis of dianionic PASAs raised another concern regarding their potential to reach and integrate into DNA. However, the mechanism behind the stability of organoarsenicals needs to be investigated for a further understanding of the metabolic process involving methylation of As i .