New Methodology for Computer-Aided Modelling of Biomolecular Structure and Dynamics 2. Local Deformations and Cycles

Abstract

A new methodology for the conformational modelling of biomolecular systems (1) is extended to local deformations of chain molecules and to flexible molecular rings. It is shown that these two cases may be reduced to considering an equivalent molecular model with a regular tree-like topology. A simple procedure is developed to analyze any flexible rings (the five- and six-membered suguar rings of carbohydrates and nucleic acids, in particular) and local deformation regions by energy minimization. Dynamic equations are also derived for such molecular systems. As a result, a unified approach is proposed for the efficient energy minimization and simulation of dynamic behavior of multimolecular systems having any set of variable internal coordinates, local deformation regions and cycles. Advantages and domains of applicability of the approach are discussed.     

Dynamics of the aromatic amino acid residues in the globular conformation of the basic pancreatic trypsin inhibitor (BPTI). II. Semi-empirical energy calculations.

The molecular conformation of the basic pancreatic trypsin inhibitor (BPTI) is known in considerable detail from both X-ray studies in single crystals and NMR studies in solution. The NMR experiments showed that the aromatic rings of the phenylalanyl and tyrosyl residues can undergo rapid rotational motions about the C beta--Cv bond. The present paper describes a model investigation of the mechanistic aspects of these intramolecular rotational motions. From calculations of the conformational energies for molecular species derived from the X-ray structure by rotations of individual aromatic rings, it was apparent that the rotational motions of the aromatics could only be understood in a flexible structure. Flexibility was simulated by allowing the protein to relax to an energetically favorable conformation for each of the different rotation states of the aromatic rings. It was then of particular interest to investigate how the perturbations caused by different rotation states of the aromatic rings were propagated in the protein structure. It was found that the rotation axes C beta--Cv were only slightly affected (delta X1 approximately less than 20 degrees. The most sizeable perturbations are caused by through space interactions with nearby atoms, which move away from the ring center and thus release the steric hindrance opposing the rotational motions. The values for the energy barriers obtained from the energy minimization are of the same order of magnitude as those measured by NMR.     

The combinatorial distance geometry method for the calculation of molecular conformation. I. A new approach to an old problem

A new approach to the long-standing local minimum problem of molecular energy minimization is proposed. The approach relies upon a field of computer mathematics known as combinatorial optimization, together with methods of conformational analysis derived from distance geometry. The advantages over the usual numerical techniques of optimization are, first, that the algorithms derived are globally convergent, and second, that the mathematical problems involved are well-posed and suitable for study within the modern theory of computational complexity. In this paper we introduce the approach, and describe a computer program based on it.     

Supercoiled DNA energetics and dynamics by computer simulation.

A new formulation is presented for investigating supercoiled DNA configurations by deterministic techniques. Thus far, the computational difficulties involved in applying deterministic methods to supercoiled DNA studies have generally limited computer simulations to stochastic approaches. While stochastic methods, such as simulated annealing and Metropolis-Monte Carlo sampling, are successful at generating a large number of configurations and estimating thermodynamic properties of topoisomer ensembles, deterministic methods offer an accurate characterization of the minima and a systematic following of their dynamics. To make this feasible, we model circular duplex DNA compactly by a B-spline ribbon-like model in terms of a small number of control vertices. We associate an elastic deformation energy composed of bending and twisting integrals and represent intrachain contact by a 6-12 Lennard Jones potential. The latter is parameterized to yield an energy minimum at the observed DNA-helix diameter inclusive of a hydration shell. A penalty term to ensure fixed contour length is also included. First and second partial derivatives of the energy function have been derived by using various mathematical simplifications. First derivatives are essential for Newton-type minimization as well as molecular dynamics, and partial second-derivative information can significantly accelerate minimization convergence through preconditioning. Here we apply a new large-scale truncated-Newton algorithm for minimization and a Langevin/implicit-Euler scheme for molecular dynamics. Our truncated-Newton method exploits the separability of potential energy functions into terms of differing complexity. It relies on a preconditioned conjugate gradient method that is efficient for large-scale problems to solve approximately for the search direction at every step. Our dynamics algorithm is numerically stable over large time steps. It also introduces a frequency-discriminating mechanism so that vibrational modes with frequencies greater than a chosen cutoff frequency are essentially frozen by the method. With these tools, we rapidly identify corresponding circular and interwound energy minima for small DNA rings for a series of imposed linking-number differences. These structures are consistent with available electron microscopy data. The energetic exchange of stability between the circle and the figure-8, in very good agreement with analytical results, is also detailed. Molecular dynamics trajectories at 100 femtosecond time steps then reveal the rapid folding of the unstable circular state into supercoiled forms. Significant bending and twisting motions of the interwound structures are also observed. Such information may be useful for understanding transition states along the folding pathway and the role of enzymes that regulate supercoiling.(ABSTRACT TRUNCATED AT 400 WORDS)     

An exact expression to calculate the derivatives of position-dependent observables in molecular simulations with flexible constraints

In this work, we introduce an algorithm to compute the derivatives of physical observables along the constrained subspace when flexible constraints are imposed on the system (i.e., constraints in which the constrained coordinates are fixed to configuration-dependent values). The presented scheme is exact, it does not contain any tunable parameter, and it only requires the calculation and inversion of a sub-block of the Hessian matrix of second derivatives of the function through which the constraints are defined. We also present a practical application to the case in which the sought observables are the Euclidean coordinates of complex molecular systems, and the function whose minimization defines the flexible constraints is the potential energy. Finally, and in order to validate the method, which, as far as we are aware, is the first of its kind in the literature, we compare it to the natural and straightforward finite-differences approach in a toy system and in three molecules of biological relevance: methanol, N-methyl-acetamide and a tri-glycine peptide.     

A general approach to the optimization of the conformation of ring molecules with an application to valinomycin

A general and efficient methodology is presented which allows molecules containing one or many rings of any size to be manipulated within energy minimization procedures. Variables describing the conformation of the molecules concerned are limited to dihedral and ring valence angles and the ring closure conditions are treated as equality constraints. An application is made to the ion transporter valinomycin and its complexes with K+ and Na+ which illustrates the possibilities of the approach and leads to results which allow a better understanding of the conformational mechanics of this important ionophore.     

Calculation of the conformation of glycosphingolipids. 2. GM1- and GM2-gangliosides

The conformations of GM1- and GM2-gangliosides have been predicted by energy minimization techniques including an orbital force field approach. The global energy minimum conformers for these two gangliosides show marked differences, particularly in the relative orientation of the sugar rings. The predicted structures are compared with those postulated from NMR spectroscopy in relation to the formation of a cation-binding site. The minimum energy conformer of GM2-ganglioside is able to form this binding site whereas this conformer of GM1-ganglioside is not. The nature of the specific interaction of gangliosides with activator proteins is discussed.     

New general approach for determining the solution structure of a ligand bound weakly to a receptor: structure of a fibrinogen Aalpha-like peptide bound to thrombin (S195A) obtained using NOE distance constraints and an ECEPP/3 flexible docking program.

A new approach incorporating flexible docking simulations and NMR data is presented for calculating the bound conformation of a ligand that interacts weakly with an enzyme. This approach consists of sampling directly the conformation of a flexible ligand inside a receptor active site containing surrounding flexible loops. To make this sampling efficient, a ligand-growing procedure has been adopted. Optimization of the ECEPP/3-plus-NOE constraint function is carried out by using a collective variable Monte Carlo minimization technique. Numerous energy minimizations are made possible for such a large system by using a Bezier splines energy grid technique. This new flexible docking approach was applied to determine the structure of a fibrinogen Aalpha-like peptide (7DFLAEGGGVRGPRV20) bound to an active site mutant of thrombin [thrombin(S195A)]. Structure calculations of the bound ligand, using 2D-transferred NOESY distance constraints in the DIANA program, showed that the N-terminal portion of the peptide (D7-R16) involves a chain reversal, whereas the C-terminal portion (G17-V20) adopts a fold that exists in several different orientations. In addition, the ECEPP/3 flexible docking package was used to assess the conformational variability of the ligand and surrounding 60D-insertion loop of thrombin. Amino acid residues (17-20) of the peptide interact with a region of the enzyme that exhibits broad specificity, with a preferred direction between the 60D-insertion loop and Pro37 of thrombin.     

The importance of being exhaustive. Optimization of bridging structural water molecules and water networks in models of biological systems

Algorithms and protocols are described for the optimization for H-bonding of isolated singular H2O molecules and entire networks of H2O molecules. Unlike other approaches that are prone to being trapped in local energy minima, these methods rely on exhaustive searches of orientation space for the H2O molecules. The results are scored with the HINT hydropathic interaction model, but the algorithms should be general for any energy-scoring computation. Two examples are provided: 1) the tightly-bound H2O molecule 301 of HIV-1 protease is shown to be more reasonably oriented in terms of forming H-bonds with this method than with a molecular mechanics energy minimization method; and 2) the H2O network surrounding carbonmonoxymyoglobin is constructed and analyzed for a 1.80-A neutron-diffraction structure. The H-atom positions calculated with this method show a somewhat better agreement with the experimental results than do the H-atom positions calculated with molecular mechanics, and both are considerably better than random.     

Molecular dynamics and free energy analysis of neuraminidase-ligand interactions

We report molecular dynamics calculations of neuraminidase in complex with an inhibitor, 4-amino-2-deoxy-2,3-didehydro-N-acetylneuraminic acid (N-DANA), with subsequent free energy analysis of binding by using a combined molecular mechanics/continuum solvent model approach. A dynamical model of the complex containing an ionized Glu119 amino acid residue is found to be consistent with experimental data. Computational analysis indicates a major van der Waals component to the inhibitor-neuraminidase binding free energy. Based on the N-DANA/neuraminidase molecular dynamics trajectory, a perturbation methodology was used to predict the binding affinity of related neuraminidase inhibitors by using a force field/Poisson-Boltzmann potential. This approach, incorporating conformational search/local minimization schemes with distance-dependent dielectric or generalized Born solvent models, correctly identifies the most potent neuraminidase inhibitor. Mutation of the key ligand four-substituent to a hydrogen atom indicates no favorable binding free energy contribution of a hydroxyl group; conversely, cationic substituents form favorable electrostatic interactions with neuraminidase. Prospects for further development of the method as an analysis and rational design tool are discussed.     

Elastic energy storage in unmineralized and mineralized extracellular matrices (ECMs): a comparison between molecular modeling and experimental measurements

In order to facilitate locomotion and limb movement many animals store energy elastically in their tendons. In the turkey, much of the force generated by the gastrocnemius muscle is stored as elastic energy during tendon deformation and not within the muscle. As limbs move, the tendons are strained causing the collagen fibers in the extracellular matrices to be strained. During growth, avian tendons mineralize in the portions distal to the muscle and show increased tensile strength, modulus, and energy stored per unit strain as a result. In this study the energy stored in unmineralized and mineralized collagen fibers was measured and compared to the amount of energy stored in molecular models. Elastic energy storage values calculated using the molecular model were slightly higher than those obtained from collagen fibers, but display the same increases in slope as the fiber data. We hypothesize that these increases in slope are due to a change from the stretching of flexible regions of the collagen molecule to the stretching of less flexible regions. The elastic modulus obtained from the unmineralized molecular model correlates well with elastic moduli of unmineralized collagen from other studies. This study demonstrates the potential importance of molecular modeling in the design of new biomaterials.     

A novel approach to decoy set generation: designing a physical energy function having local minima with native structure characteristics

We suggest a new approach to the generation of candidate structures (decoys) for ab initio prediction of protein structures. Our method is based on random sampling of conformation space and subsequent local energy minimization. At the core of this approach lies the design of a novel type of energy function. This energy function has local minima with native structure characteristics and wide basins of attraction. The current work presents our motivation for deriving such an energy function and also tests the derived energy function.Our approach is novel in that it takes advantage of the inherently rough energy landscape of proteins, which is generally considered a major obstacle for protein structure prediction. When local minima have wide basins of attraction, the protein's conformation space can be greatly reduced by the convergence of large regions of the space into single points, namely the local minima corresponding to these funnels. We have implemented this concept by an iterative process. The potential is first used to generate decoy sets and then we study these sets of decoys to guide further development of the potential. A key feature of our potential is the use of cooperative multi-body interactions that mimic the role of the entropic and solvent contributions to the free energy.The validity and value of our approach is demonstrated by applying it to 14 diverse, small proteins. We show that, for these proteins, the size of conformation space is considerably reduced by the new energy function. In fact, the reduction is so substantial as to allow efficient conformational sampling. As a result we are able to find a significant number of near-native conformations in random searches performed with limited computational resources.     

Determining minimum energy conformations of polypeptides by dynamic programming

A combinatorial optimization approach is used for solving the multiple-minima problem when determining the low-energy conformations of short polypeptides. Each residue is represented by a finite number of discrete states corresponding to single residue local minima of the energy function. These precomputed values constitute a search table and define the conformational space for discrete minimization by a generalized dynamic programming algorithm that significantly limits the number of intermediate conformations to be generated during the search. Since dynamic programming involves stagewise decisions, it results in buildup-type procedures implemented in two different forms. The first procedure predicts a number of conformations by a completely discrete search and these are subsequently refined by local minimization. The second involves limited continuous local minimization within the combinatorial algorithm, generally restricted to two dihedral angles in a buildup step. Both procedures are tested on 17 short peptides previously studied by other global minimization methods but involving the same potential energy function. The discrete method is extremely fast, but proves to be successful only in 14 of the 17 test problems. The version with limited local minimization finds, however, conformations in all the 17 examples that are close to the ones previously presented in the literature or have lower energies. In addition, results are almost independent of the cutoff energy, the most important parameter governing the search. Although the limited local minimization increases the number of energy evaluations, the method still offers substantial advantages in speed.     

Thermitase, a thermostable subtilisin: comparison of predicted and experimental structures and the molecular cause of thermostability

The Subtilisin family of proteases has four members of known sequence and structure: subtilisin Carlsberg, Subtilisin novo, proteinase K, and thermitase. Using thermitase as a test case, we ask two questions. How good are methods for model building a three-dimensional structure of a protein based on sequence homology to a known structure? And what are the molecular causes of thermostability? First, we compare predicted models of thermitase, refined by energy minimization and varied by molecular dynamics, with the preliminary crystal structure. The predictions work best in the conserve structural core and less well in seven loop regions involving insertions and deletions relative to Subtilisin. Here, variation of loop regions by molecular dynamics simulation in vacuo followed by energy minimization does not improve the prediction since we find no correlation between in vacuo energy and correctness of structure when comparing local energy minima. Second, in order to identify the molecular case of thermostability we confront hypotheses erived by calculation of the details of interatomic interactions with inactivation experiments. As a result, we can exclude salt bridges and hydrophobic interactions as main cause of thermostability. Based on a combination of theoretical and experimental evidence, the unusually tight binding of calcium by thermitase emerges as the most likely single influence responsible for its increased thermostability.     

A General Approach to the Optimization of the Conformation of Ring Molecules with an Application to Valinomycin

Abstract

A general and efficient methodology is presented which allows molecules containing one or many rings of any size to be manipulated within energy minimization procedures. Variables describing the conformation of the molecules concerned are limited to dihedral and ring valence angles and the ring closure conditions are treated as equality constraints. An application is made to the ion transporter valinomycin and its complexes with K⁺ and Na⁺ which illustrates the possibilities of the approach and leads to results which allow a better understanding of the conformational mechanics of this important ionophore.     

Activated Conformations of the ras-gene-Encoded p21 Protein. 1. An Energy-Refined Structure for the Normal p21 Protein Complexed with GDP

Abstract

A complete three-dimensional structure for the ras-gene-encoded p21 protein with Gly 12 and Gin 61, bound to GDP, has been constructed in four stages using the available α-carbon coordinates as deposited in the Brookhaven National Laboratories Protein Data Bank. No all-atom structure has been made available despite the fact that the first crystallographic structure for the p21 protein was reported almost four years ago. In the p21 protein, if amino acid substitutions are made at any one of a number of different positions in the amino acid sequence, the protein becomes permanently activated and causes malignant transformation of normal cells or, in some cell lines, differentiation and maturation. For example, all amino acids except Gly and Pro at position 12 result in an oncogenic protein; all amino acids except Gin, Glu and Pro at position 61 likewise cause malignant transformation of cells. We have constructed our all-atom structure of the non-oncogenic protein from the x-ray structure in order to determine how oncogenic amino acid substitutions affect the three-dimensional structure of this protein. In Stage 1 we generated a poly-alanine backbone (except at Gly and Pro residues) through the α-carbon structure, requiring the individual Ala, Pro or Gly residues to conform to standard amino acid geometry and to form trans-planar peptide bonds. Since no a-carbon coordinates for residues 60–65 have been determined these residues were modeled by generating them in the extended conformation and then subjecting them to molecular dynamics using the computer application DISCOVER and energy minimization using DISCOVER and the ECEPP (Empirical Conformational Energies for Peptides Program). In Stage 2, the positions of residues that are homologous to corresponding residues of bacterial elongation factor Tu (EF-Tu) to which p21 bears an overall 40% sequence homology, were determined from their corresponding positions in a high-resolution structure of EF-Tu. Non-homologous loops were taken from the structure generated in Stage 1 and were placed between the appropriate homologous segments so as to connect them. In Stage 3, all bad contacts that occurred in this resulting structure were removed, and the coordinates of the α-carbon atoms were forced to superimpose as closely as possible on the corresponding atoms of the reference (x-ray) structure. Then the side chain positions of residues of the nonhomologous loop regions were modeled using a combination of molecular dynamics and energy minimization using DISCOVER and ECEPP respectively. All of the residues of the structure were then allowed to move under restrained energy minimization where the restraints were gradually removed. In Stage 4, the nucleotide GOP was added to the model and further energy minimization was carried out. The energy of the protein-GOP complex was minimized by allowing the atoms of GOP to move with the protein held fixed and then by allowing both the nucleotide and the residues of the protein to move together. The reconstructed model agrees with the published features of the p21 protein-GOP complex including the hydrogen bonding scheme, the distribution of backbone dihedral angles, the residues contacting the nucleotide, and the orientation of loops with respect to one another in the protein. The structure also agrees with one that was predicted previously (Chen, J.M. et al., J. Biomol. Struct. Dynamics 6, 850–875 (1989)). In our molecular dynamics-energy minimization procedures, we also have been able to place all residues except Ala 66, which occurs in a poorly-defined region crystallographically, in local single residue minima, including residues reported to be in high energy regions in the x-ray structure. The constructed model can explain observed physical phenomena such as autophosphorylation by GTP on Thr 59 in proteins containing Thr in place of Ala 59.     

Conformational characteristics of mixed sugar puckered deoxytetranucleoside triphosphate d-GpCpGpC from energy minimization studies

As a continuation of our theoretical studies on nucleic acid subunit systems, in this article we consider the case of the tetranucleoside d-GpCpGpC, the minimally ideal representative unit for analyzing the relative stabilities of different forms of homo- and mixed helical conformation of polynucleotides. The four sugar rings are kept so as to generate B-genus, B+A genus and Z-genus conformations. Twenty five helical conformational states which resulted from judicious mixing of A-, B-, C-, W-, and Z-, states locally are subjected to energy minimization permitting the 19 dihedral angles to vary simultaneously. Conformational states corresponding to regular helical forms and mixed helical forms, when analyzed provide valuable information as to the local conformational flexibility and transitions available to polynucleotides.     

The Study of Possible A and B Conformations of Alternating DNA Using a New Program for Conformational Analysis of Duplexes (CONAN)

Abstract

A new program, CONAN has been designed for CONformational ANalysis of oligonucleotide duplexes with natural and modified bases. It allows to model both regular DNA fragments with different types of symmetry and irregular ones including bends, junctions, mismatched pairs and base lesions. Computations and minimization of the energy are performed in a space of internal structural variables chosen to build start structure easier and conveniently analyze the results obtained. These internal structural variables determine mutual base-base and base-sugar arrangement and sugar puckering. The analytical closure procedure is applied both to sugar rings and to backbone fragments between adjacent sugars. For more effective energy minimization, analytical gradient is calculated. The CONAN was applied to the search for low-energy conformations of poly(dA-dT)·poly(dA-dT) and poly(dG-dC)·poly(dG-dC) duplexes. Extended regions of low-energy A and B conformations are revealed and characterized. These regions contain structures with different relative values of helical twist, τ, for pur-pyr and pyr-pur steps, namely, conformations with τ(pur-pyr)>τ(pyr-pur) and with τ(pur-pyr)<τ(pyr-pur). Two types of sugar puckering were found for B-form low-energy conformations, the first type with all C2′-endo sugar residues and the second one—;with C2′-endo purines and O1′-endo pyrimidines. The calculated conformations are compared with X-ray diffraction data for crystals and fibers and NMR data for solution.     

Nucleic acid base-pairing and N-methylacetamide self-association in chloroform: affinity and conformation.

A recently developed computational method, 'mining minima', is used to examine the hydrogen-bonding interactions of nucleic acid base-pairs and of the N-methylacetamide homodimer in chloroform. The mining minima algorithm aggressively samples molecular conformations, identifies the most important local minima, and computes their contributions to the overall free energy of the system. Here, the CHARMM 98 parameter set is used for the potential energy and the generalized Born/surface area solvent model is used to account for the influence of the solvent. Good agreement with experiment is obtained for the non-covalent binding affinities of a series of complexes. The computational approach used here is applicable to a range of molecular systems.