Using internal and collective variables in Monte Carlo simulations of nucleic acid structures: chain breakage/closure algorithm and associated Jacobians

This article describes a method for solving the geometric closure problem for simplified models of nucleic acid structures by using the constant bond lengths approximation. The resulting chain breakage/closure equations, formulated in the space of variable torsion and bond angles, are easy to solve, and have only two solutions. The analytical simplicity is in contrast with the high complexity of the closure problem in the torsion angle space with at most 16 solutions, which has been dealt with by several authors and was solved analytically by Wu and Deem (J. Chem. Phys. 1999, 111, 6625). The discussion on the choice of variables and associated Jacobians is focussed on the question of how conformational equilibration is affected in Monte Carlo simulations of molecular systems. In addition to the closure of the phosphate backbone, it is necessary to also solve the closure problem for the five-membered flexible furanose sugar ring. Explicit closure equations and the resulting Jacobians are given both for the complete four-variable model of the furanose ring and simulations in the phase-amplitude space of the five-membered ring, which are based on the approximate two-variable model of furanose introduced by Gabb et al. (J. Comput. Chem. 1995, 16, 667). The suggested closure algorithm can be combined with collective variables defined by translations and rotations of the monomeric nucleotide units. In comparison with simple internal coordinate moves, the resulting concerted moves describe local structural changes that have high acceptance rates and enable fast conformational equilibration. Appropriate molecular models are put forward for prospective Monte Carlo simulations of nucleic acids, but can be easily adapted to other biomolecular systems, such as proteins and lipid structures in biological membranes.     

Explicit factorization of external coordinates in constrained statistical mechanics models

If a macromolecule is described by curvilinear coordinates or rigid constraints are imposed, the equilibrium probability density that must be sampled in Monte Carlo simulations includes the determinants of different mass-metric tensors. In this work, the authors explicitly write the determinant of the mass-metric tensor G and of the reduced mass-metric tensor g, for any molecule, general internal coordinates and arbitrary constraints, as a product of two functions; one depending only on the external coordinates that describe the overall translation and rotation of the system, and the other only on the internal coordinates. This work extends previous results in the literature, proving with full generality that one may integrate out the external coordinates and perform Monte Carlo simulations in the internal conformational space of macromolecules.     

Potential energy constrained molecular dynamics simulations

A method for carrying out molecular dynamics simulations in which the potential energy U of the molecular system is constrained at its initial value is developed and thoroughly tested. The constraint is not introduced within the framework of the Lagrange multipliers technique, rather it is fulfilled in a natural way by carrying out the simulations in terms of suitable sets of delocalized coordinates. Such coordinates are defined by an appropriate tuning of the Baker, Kessi, and Delley internal delocalized nonredundant coordinates technique [J. Chem. Phys. 105, 192 (1996)]. The proposed method requires multiple evaluations of energy and gradients in each step of the molecular dynamics simulation, so that constant U simulations suffer some overhead compared to ordinary simulations. But the particular formulation of the delocalized coordinates and of the equations of motion greatly simplifies all the various steps required by the Baker's technique, thus allowing for the efficient implementation of the method itself. The technique is reliable and allows for very high accuracy in the potential energy conservation during the whole simulation. Moreover, it proved to be free of drift troubles which can occur when standard constraint methods are straightforwardly implemented without the application of appropriate correcting techniques.     

Rotational symmetry of the molecular potential energy in the Cartesian coordinates

We consider the molecular Born-Oppenheimer potential energy as a function of atomic Cartesian coordinates and discuss the non-stationary Hessian properties arising due to rotational symmetry. A connection with the extended Hessian theory is included. New applications of Cartesian representation for examining and correcting raw numerical Hessian data and a simple formulation of harmonic vibrational analysis of partially optimized systems are proposed. Exemplary calculations for the porphyrin molecule with an internal proton transfer are reported. We also develop the normal transformation method to incorporate the rotational symmetry into the approximate analytical potentials, which are parametrized in the Cartesian coordinates. The transformation converts the coordinates from the space fixed frame to the frame which translates and rotates with the molecule and is determined by the Eckart conditions. New simple analytical formulas for the first and second derivatives of the transformed potential are derived. This fast method can be used to calculate the potential and its derivatives in the simulations of chemical reaction dynamics in the space fixed Cartesian frame without the need to constrain the molecular rotation or to define the local non-redundant internal coordinates.     

Molecular structure of free canonical 2′-deoxyribonucleosides: a density functional study

The molecular structure of isolated canonical 2′-deoxyrinobucleosides was calculated using the density functional theory. It was demonstrated that the geometry of the base unit (BU) is almost unchanged compared to free nucleobases. Only slight out-of-plane deformation of the pyrimidine ring in deoxy-cytidine is observed. The conformation of the furanose ring strongly depends on the nature and orientation of the nucleobase. All nucleosides possess different conformations of this ring. Significant influence of the steric repulsion between the nucleobase and the sugar unit (SU) on puckering of the furanose ring and variation of the C–O and glycosyl bond lengths was demonstrated. The C(3′)-endo conformer of the furanose ring is more stable at the anti-orientation BU with respect to SU. An opposite trend is observed for the syn-orientation which is additionally stabilized by an intramolecular hydrogen bond with participation of the C(5′)OH group.     

Graphical approach for defining natural internal coordinates

A simple, customizable connectivity scheme is rigorously defined in which pairs of atoms are classified into three categories. The tools of graph theory are used to analyze the molecular graph and to efficiently find rings and ring assemblies through a combination of pruning and homeomorphic reduction. The definition of natural internal coordinates is extended in a nonredundant fashion for the various cases of weakly interacting components and for fused ring systems. The ring system coordinates were tested and found to be superior to Z-matrix coordinates. © 1998 John Wiley & Sons, Inc. J Comput Chem 19: 504–511, 1998     

Solid-state nuclear magnetic resonance spectroscopy studies of furanose ring dynamics in the DNA HhaI binding site

The dynamics of the furanose rings in the GCGC moiety of the DNA oligomer [d(G 1A 2T 3A 4 G 5 C 6 G 7 C 8T 9A 10T 11C 12)] 2 are studied by using deuterium solid-state NMR (SSNMR). SSNMR spectra obtained from DNAs selectively deuterated on the furanose rings of nucleotides within the 5'-GCGC-3' moiety indicated that all of these positions are structurally flexible. The furanose ring within the deoxycytidine that is the methylation target displays the largest-amplitude structural changes according to the observed deuterium NMR line shapes, whereas the furanose rings of nucleotides more remote from the methylation site have less-mobile furanose rings (i.e., with puckering amplitudes < 0.3 A). Previous work has shown that methylation reduces the amplitude of motion in the phosphodiester backbone of the same DNA, and our observations indicate that methylation perturbs backbone dynamics through the furanose ring. These NMR data indicate that the 5'-GCGC-3' is dynamic, with the largest-amplitude motions occurring nearest the methylation site. The inherent flexibility of this moiety in DNA makes the molecule more amenable to the large-amplitude structural rearrangements that must occur when the DNA binds to the HhaI methyltransferase.     

5‐Phenyluridine trihydrate

In the title compound, C₁₅H₁₆N₂O₆·~3H₂O, the substituted uracil ring is oriented in the anti position relative to the ribose ring, and the phenyl and uracil rings are oriented in a noncoplanar fashion. The furanose ring adopts a conformation close to 3′‐endo, in contrast to the furanose conformation seen in the crystal structure of the synthetic precursor 5‐bromouridine, which is close to 2′‐endo. The molecule is involved in an extensive hydrogen‐bonding network with several water molecules, some of which are disordered.     

Geometry optimization of periodic systems using internal coordinates

An algorithm is proposed for the structural optimization of periodic systems in internal (chemical) coordinates. Internal coordinates may include in addition to the usual bond lengths, bond angles, out-of-plane and dihedral angles, various "lattice internal coordinates" such as cell edge lengths, cell angles, cell volume, etc. The coordinate transformations between Cartesian (or fractional) and internal coordinates are performed by a generalized Wilson B-matrix, which in contrast to the previous formulation by Kudin et al. [J. Chem. Phys. 114, 2919 (2001)] includes the explicit dependence of the lattice parameters on the positions of all unit cell atoms. The performance of the method, including constrained optimizations, is demonstrated on several examples, such as layered and microporous materials (gibbsite and chabazite) as well as the urea molecular crystal. The calculations used energies and forces from the ab initio density functional theory plane wave method in the projector-augmented wave formalism.     

Internal‐to‐Cartesian back transformation of molecular geometry steps using high‐order geometric derivatives

In geometry optimizations and molecular dynamics calculations, it is often necessary to transform a geometry step that has been determined in internal coordinates to Cartesian coordinates. A new method for performing such transformations, the high‐order path‐expansion (HOPE) method, is here presented. The new method treats the nonlinear relation between internal and Cartesian coordinates by means of automatic differentiation. The method is reliable, applicable to any system of internal coordinates, and computationally more efficient than the traditional method of iterative back transformations. As a bonus, the HOPE method determines not just the Cartesian step vector but also a continuous step path expressed in the form of a polynomial, which is useful for determining reaction coordinates, for integrating trajectories, and for visualization. © 2013 Wiley Periodicals, Inc.     

Quasi-Hamiltonian equations of motion for internal coordinate molecular dynamics of polymers

Conventional molecular dynamics simulations of macromolecules require long computational times because the most interesting motions are very slow compared to the fast oscillations of bond lengths and bond angles that limit the integration time step. Simulation of dynamics in the space of internal coordinates, that is, with bond lengths, bond angles, and torsions as independent variables, gives a theoretical possibility of eliminating all uninteresting fast degrees of freedom from the system. This article presents a new method for internal coordinate molecular dynamics simulations of macromolecules. Equations of motion are derived that are applicable to branched chain molecules with any number of internal degrees of freedom. Equations use the canonical variables and they are much simpler than existing analogs. In the numerical tests the internal coordinate dynamics are compared with the traditional Cartesian coordinate molecular dynamics in simulations of a 56 residue globular protein. For the first time it was possible to compare the two alternative methods on identical molecular models in conventional quality tests. It is shown that the traditional and internal coordinate dynamics require the same time step size for the same accuracy and that in the standard geometry approximation of amino acids, that is, with fixed bond lengths, bond angles, and rigid aromatic groups, the characteristic step size is 4 fs, which is 2 times higher than with fixed bond lengths only. The step size can be increased up to 11 fs when rotation of hydrogen atoms is suppressed. © 1997 by John Wiley & Sons, Inc. J Comput Chem 18 : 1354–1364, 1997     

A comparative study of molecular dynamics in Cartesian and in internal coordinates: dynamical instability in the latter caused by nonlinearity of the equations of motion

The stability of a general molecular dynamics (MD) integration scheme is examined for simulations in generalized (internal plus external) coordinates (GCs). An analytic expression is derived for the local error in energy during each integration time step. This shows that the explicit dependence of the mass-matrix on GCs, which makes the system's Lagrange equations of motion nonlinear, causes MD simulations in GCs to be less stable than those in Cartesian coordinates (CCs). In terms of CCs, the corresponding mass-matrix depends only on atomic masses and thus atomistic motion is subject to the linear Newton equations, which makes the system more stable. Also investigated are two MD methods in GCs that utilize nonzero elements of the vibrational spectroscopic B-matrices. One updates positions and velocities in GCs that are iteratively adjusted so as to conform to the velocity Verlet equivalent in GCs. The other updates positions in GCs and velocities in CCs that are adjusted to satisfy the internal constraints of the new constrained WIGGLE MD scheme. The proposed methods are applied to an isolated n-octane molecule and their performances are compared with those of several CCMD schemes. The simulation results are found to be consistent with the analytic stability analysis. Finally, a method is presented for computing nonzero elements of B-matrices for external rotations without imposing the Casimir-Eckart conditions.     

New way of describing static and dynamic deformations of the Jahn-Teller type in ring molecules

A new method is presented to describe deformations of an N-membered planar ring (N-ring) molecule in terms of deformation vectors that can be expressed by a set of 2N-3 deformation amplitudes and phase angles. The deformation coordinates are directly derived from the normal vibrational modes of the N-ring and referenced to a regular polygon (N-gon) of unit length. They extend the conceptual approach of the Cremer-Pople puckering coordinates (J. Am. Chem. Soc. 1975, 97, 1354) to the planar ring and make it possible to calculate, e.g., a planar ring of special deformation on a Jahn-Teller surface. It is demonstrated that the 2N-3 deformation parameters are perfectly suited to describe the pseudorotation of a bond through the ring as it is found in cyclic Jahn-Teller systems. In general, an N-membered planar ring can undergo N-2 different bond pseudorotations provided the energetics of such a process is feasible. The Jahn-Teller distortions observed in ring compounds correspond either directly to the basic pseudorotation modes or to linear combinations of them. Any deformed ring molecule can be characterized in terms of the new ring deformation coordinates, which help to identify specific electronic effects. The usefulness of the ring deformation coordinates is demonstrated by calculating the Jahn-Teller surfaces for bond pseudorotation in the case of the cyclopropyl radical cation and cyclobutadiene as well as the ring deformation surfaces of disulfur dinitride and its dianion employing multireference averaged quadratic coupled cluster (MR-AQCC) theory, equation-of-motion coupled cluster theory in form of EOMIP-CCSD, and single determinant coupled cluster theory in form of CCSD(T).     

PMR investigation into the structure of some n-substituted 1-amino-1-deoxy-D-fructoses (amadori rearrangement products) : Evidence for a preferential conformation in solution

The PMR spectra at 220 MHz of some Amadori rearrangement products deduced from D-glucose with p-toluidine (1), N-methylphenylamine (2), di-butylamine (4), piperdine (5), and morpholine (6) have been studied in detail.Compounds 1-6 appear to exist in solution predominantly as an equilibrium mixture of the furanose and pyranose ring. The pyranose ring occurs exclusively in the β(D)-²C₅-conformation (corresponds to Reeves 1C-conformation). The furanose ring probably exists as a mixture of both the β- and α-anomer, in which the β-anomer is favoured.     

The phase behavior of polyethylene ring chains

The equilibrium properties of an isolated polyethylene ring chain are studied by using molecular dynamics (MD) simulations. The results of an 80-bond linear chain are also presented, which are in agreement with previous studies of square-well chains and Lennard-Jones (LJ) homopolymers. Mainly, we focus on the collapse of polyethylene ring chains. At high temperatures, a fully oblate structure is observed for the ring chains with different chain lengths. For such an oblate structure, a shape factor of delta(*)=0.25 and a rodlike scaling relation between the radius of gyration and chain lengths could be deduced easily in theory, and the same results are obtained by our MD simulations. Such an oblate structure can be obtained by Monte Carlo simulation only for sufficient stiff ring chains. When the temperature decreases, an internal energy barrier is observed. This induces a strong peak in the heat capacity, denoting a gas-liquid-like transition. This energy barrier comes mainly from the local monomer-monomer interactions, i.e., the bond-stretching, the bond-bending, and the torsion potentials. A low temperature peak is also observed in the same heat capacity curve, representing a liquid-solid-like transition. These numerical simulation results support a two-stage collapse of polyethylene ring chains; however, the nature should be different from the square-well and LJ ring chains.     

Synthesis and conformational studies of peptidomimetics containing furanoid sugar amino acids and a sugar diacid

Furanoid sugar amino acids (1) were synthesized and used as dipeptide isosteres to induce interesting turn structures in small linear peptides. They belong to a new variety of designed hybrid structures that carry both amino and carboxyl groups on rigid furanose sugar rings. Four such molecules, 6-amino-2,5-anhydro-6-deoxy-D-gluconic acid (3, Gaa) and its mannonic (4, Maa), idonic (5, Iaa), and a 3,4-dideoxyidonic (6, ddIaa) congeners were synthesized. The synthesis followed a novel reaction path in which an intramolecular 5-exo S(N)2 opening of the hexose-derived terminal aziridine ring in 2 by the gamma-benzyloxy oxygen with concomitant debenzylation occurred during pyridinium dichromate oxidation of the primary delta-hydroxyl group to carboxyl function, leading to the formation of furanoid sugar amino acid frameworks in a single step. Incorporation of these furanoid sugar amino acids into Leu-enkephalin replacing its Gly-Gly portion gave analogues 8-11. Detailed structural analysis of these molecules by circular dichroism (CD) and various NMR techniques in combination with constrained molecular dynamics (MD) simulations revealed that two of these analogues, 8a and 10a, have folded conformations composed of an unusual nine-membered pseudo beta-turn-like structure with a strong intramolecular H-bond between LeuNH --> sugarC3-OH. This, in turn, brings the two aromatic rings of Tyr and Phe in close proximity, a prerequisite for biological activities of opioid peptides. The analgesic activities of 8a,b determined by mouse hot-plate and tail-clip methods were similar to that of Leu-enkephalin methyl ester. The syn disposition of the beta-hydroxycarboxyl motif on the sugar rings appears to be the driving force to nucleate the observed turn structures in some of these molecules (8 and 10). Repetition of the motif on both sides of a furanose ring resulted in a novel molecular design of sugar diacid, 2,5-anhydro-D-idaric acid (7, Idac). Bidirectional elongation of the diacid moieties of 7 with identical peptide strands led to the formation of a C2-symmetric reverse-turn mimetic 12 which displayed a very ordered structure consisting of identical intramolecular H-bonds at two ends between LeuNH --> sugar-OH, the same as in 8 and 10.