Kevin J. Naidoo

Carbohydrate solution simulations: Producing a force field with experimentally consistent primary alcohol rotational frequencies and populations

We present a CHARMM Carbohydrate Solution Force Field (CSFF) suitable for nanosecond molecular dynamics computer simulations. The force field was derived from a recently published sugar parameter set. 1 Dihedral angle parameters for the primary alcohol as well as the secondary hydroxyl groups were adjusted. Free energy profiles of the hydroxymethyl group for two monosaccharides (β‐D‐glucose and β‐D‐galactose) were calculated using the new parameter set and compared with similar force fields. Equilibrium rotamer populations obtained from the CSFF are in excellent agreement with NMR data (glucose gg:gt:tg ≈ 66:33:1 and galactose gg:gt:tg ≈ 4:75:21). In addition, the primary alcohol rotational frequency is on the nanosecond time scale, which conforms to experimental observations. Equilibrium population distributions of the primary alcohol conformers for glucose and galactose are reached within 10 nanoseconds of molecular dynamics simulations. In addition, gas phase vibrational frequencies computed for β‐D‐glucose using this force field compare well with experimental frequencies. Carbohydrate parameter sets that produce both conformational energies and rotational frequencies for the pyranose primary alcohol group that are consistent with experimental observations should allow for increased accuracy in modeling the flexibility of biologically important (1‐6)‐linked saccharides in solution.

Haemozoin (β-haematin) biomineralization occurs by self-assembly near the lipid/water interface

Several blood‐feeding organisms, including the malaria parasite detoxify haem released from host haemoglobin by conversion to the insoluble crystalline ferriprotoporphyrin IX dimer known as haemozoin. To date the mechanism of haemozoin formation has remained unknown, although lipids or proteins have been suggested to catalyse its formation. We have found that β‐haematin (synthetic haemozoin) forms rapidly under physiologically realistic conditions near octanol/water, pentanol/water and lipid/water interfaces. Molecular dynamics simulations show that a precursor of the haemozoin dimer forms spontaneously in the absence of the competing hydrogen bonds of water, demonstrating that this substance probably self‐assembles near a lipid/water interface in vivo.

Ring puckering: a metric for evaluating the accuracy of AM1, PM3, PM3CARB-1, and SCC-DFTB carbohydrate QM/MM simulations.

The puckered conformations of furanose and pyranose carbohydrate rings are central to analyzing the action of enzymes on carbohydrates. Enzyme reaction mechanisms are generally inaccessible to experiments and so have become the focus of QM(semiempirical)/MM simulations. We show that the complete free energy of puckering is required to evaluate the accuracy of semiempirical methods used to study reactions involving carbohydrates. Interestingly, we find that reducing the free energy space to lower dimensions results in near meaningless minimum energy pathways. We analyze the furanose and pyranose free energy pucker surfaces and volumes using AM1, PM3, PM3CARB-1, and SCC-DFTB. A comparison with DFT optimized structures and a HF free energy surface reveals that SCC-DFTB provides the best semiempirical description of five- and six-membered carbohydrate ring deformation.

Water structure about the dimer and hexamer repeat units of amylose from molecular dynamics computer simulations

We have analyzed a set of molecular dynamics (MD) trajectories of maltose in vacuum and water for solute imposed structuring on the solvent. To do this, we used a novel technique to calculate water probability densities to locate the areas in which the solvent is most populated in the maltose solution. We found that only the layer of water within the first maltose hydration shell has a probability density 50% and greater than that of bulk water. On investigating this water layer using Voronoi polyhedra (VP) analysis it was seen that only the waters adjacent to the hydrophobic (CH and CH2) groups are more structured than bulk water. We found that in a maltose solution of approximately 1.0 g/cm3 the solute does not disrupt the structure of the surrounding water beyond the first hydration shell. Next we performed a 700‐ps MD simulation of a maltohexaose strand in a box of 4096 SPC/E waters. The water probability density calculations and the VP analysis of the maltohexaose solution show that the larger amylose repeat unit decreases the solvent configurational entropy of the water beyond the first hydration shell. Analysis of this trajectory reveals that the helical conformation of the maltohexaose strand is preserved via bridging intermolecular water hydrogen bonds, indicating that a single amylose helical turn in water is preserved by hydrophilic and not hydrophobic interactions. Using VP analysis we present a method to accurately determine the number of water molecules in the first hydration shell of dissolved solutes. In the case of maltose, there are 40 water molecules in this shell, while for maltohexaose the number is 98.

Implementation of an adaptive umbrella sampling method for the calculation of multidimensional potential of mean force of chemical reactions in solution

We describe the implementation of an adaptive umbrella sampling method, making use of the weighted histogram analysis method, for computing multidimensional potential of mean force for chemical reaction in solution. The approach is illustrated by investigating the effect of aqueous solution on the free energy surface for the proton transfer reaction of [H3N—H—NH3]+ using a combined quantum mechanical and molecular mechanical AM1/TIP3P potential.

Stereoelectronic and Solvation Effects Determine Hydroxymethyl Conformational Preferences in Monosaccharides

Although the conformational preferences in glucose and galactose have been studied since the early 1970s, only recently have the glucose and galactose hydroxymethyl populations been resolved by combining (3)J(HH) and (2)J(HH) NMR coupling data using a modified Karplus equation. A preference for gauche conformations is observed in monosaccharides, but the reasons for this are not understood. We calculated the free energy of rotation profiles for glucose and galactose primary alcohols using a semiempirical description of the monosaccharides in QM/MM simulations. From this we observed excellent agreement between our simulated population distributions for glucose gg/gt/tg = 35:57:3 and galactose gg/gt/tg = 4:86:7 with those measured from NMR. A stereoelectronic analysis of the minimum energy conformations using natural bond orbitals provides a clear description of the stabilizing contribution to the gauche conformers stemming from the C-H bonding and the C-O antibonding orbital interactions, specifically sigma(C6-H) --> sigma*(C5-O5) and sigma(C5-H) --> sigma*(C6-O6). Analysis of the solution trajectories reveals that persistent intramolecular hydrogen bonds and intermolecular bridging hydrogen bonds formed by water molecules between the ring oxygen and the hydroxymethyl group further stabilizes the gt conformation making it the preferred rotamer in both hydrated glucose and galactose. The hydroxymethyl quantum mechanics/molecular mechanics molecular dynamics trajectories and derived rotational free energies for these monosaccharides in water solutions explain that the experimental observations are due to a combination of competing stereoelectronic (gauche), electronic (intramolecular hydrogen bonding), and electrostatic (solvent-saccharide hydrogen bonding) factors.

Molecular details from computational reaction dynamics for the cellobiohydrolase I glycosylation reaction.

Glycosylation of cellobiose hydrolase I (CBHI), is a key step in the processing and degradation of cellulose. Here the pathways and barriers of the reaction are explored using the free energy from adaptive reaction coordinate forces (FEARCF) reaction dynamics method coupled with SCC-DFTB/MM. In many respects CBHI follows the expected general GH7 family mechanism that involves the Glu-X-Asp-X-X-Glu motif. However, critical electronic and conformational details, previously not known, were discovered through our computations. The central feature that ensures the success of the glycosylation reaction are the Glu212 nucleophiles hydrogen bond to the hydroxyl on C2, of the glucose in the -1 position of the cellulosic strand. This Glu212 function restricts the C2 hydroxyl in such a way as to favor the formation of the (4)E ring pucker of the -1 position glucose. A frontier molecular orbital analysis of the structures along the reaction surface proves the existence of an oxocarbenium ion, which has both transition state and intermediate character. The transition state structure is able to descend down the glycosylation pathway through the critical combination of Asp214 (HOMO), ring oxygen (LUMO), and Glu212 (HOMO), anomeric carbon (LUMO) interactions. Using the fully converged FEARCF SCC-DFTB/MM reaction surface, we find a barrier of 17.48 kcal/mol separating bound cellulose chain from the glycosylated CBHI. Taking recrossing into account gives k(cat) = 0.415 s(-1) for cellobiohydrolase glycosylation.

Acceleration of the GAMESS‐UK electronic structure package on graphical processing units

The approach used to calculate the two‐electron integral by many electronic structure packages including generalized atomic and molecular electronic structure system‐UK has been designed for CPU‐based compute units. We redesigned the two‐electron compute algorithm for acceleration on a graphical processing unit (GPU). We report the acceleration strategy and illustrate it on the (ss|ss) type integrals. This strategy is general for Fortran‐based codes and uses the Accelerator compiler from Portland Group International and GPU‐based accelerators from Nvidia. The evaluation of (ss|ss) type integrals within calculations using Hartree Fock ab initio methods and density functional theory are accelerated by single and quad GPU hardware systems by factors of 43 and 153, respectively. The overall speedup for a single self consistent field cycle is at least a factor of eight times faster on a single GPU compared with that of a single CPU.

Pyranose Ring Transition State Is Derived from Cellobiohydrolase I Induced Conformational Stability and Glycosidic Bond Polarization

Understanding carbohydrate ring pucker is critical to rational design in materials and pharmaceuticals. Recently we have generalized our adaptive reaction coordinate force biasing method to perform calculations on multidimensional reaction coordinates. We termed this the Free Energies from Adaptive Reaction Coordinate Forces (FEARCF) method. Using FEARCF in SCC-DFTB QM/MM non-Boltzmann simulations, we are able to calculate multidimensional ring pucker free energies of conformation. Here we apply this to the six-membered glucopyranose ring located in an eight-membered β 1-4 linked octaose oligosaccharide (cellooctaose). The cellooctaose was built following the conformation of the saccharides bound to cellobiohydrolase I (CBHI) of Trichoderma reesei as reported in the 7CEL crystal structure obtained from the PDB. We calculate the free energy of ring puckering of the glucopyranose ring at the -1 position in vacuum, in water, and bound to the protein. We find that the protein induces (4)E and (4)H(3) conformations that are much more stable than the usually preferred (4)C(1) conformer. Furthermore, for the (4)H(3) conformation in the catalytic binding domain, there is significant electronic rearrangement that drives the structure toward the transition state of the glycosylation reaction.

Free Energies from Adaptive Reaction Coordinate Forces (FEARCF): an application to ring puckering

Previously an adaptive reaction coordinate force biasing method based on probability distributions has been used to calculate the free energy of conformation, configuration and chemical reactions. This method has recently been generalised to perform calculations on multidimensional reaction coordinates. This paper presents details of this method, termed Free Energies from Adaptive Reaction Coordinate Forces (FEARCF). The efficiency of sampling of this method is demonstrated by applying it to the problem of sampling the many characteristic pucker conformations of a β-D-glucose pyranose using a semi-empirical PM3 Hamiltonian. The sampling ratio of the global minimum conformer (4C1) to the highest energy conformer (a planar hexapyranose ring) is 1.7:1. Pucker free energy surfaces such as the one presented here can be a useful tool in the analysis of enzymatic reactions involving molecular ring puckers in the Michaelis complex.