Inactivation and reactivation of ribonuclease A studied by computer simulation

The year 2011 marked the half-centenary of the publication of what came to be known as the Anfinsen postulate, that the tertiary structure of a folded protein is prescribed fully by the sequence of its constituent amino acid residues. This postulate has become established as a credo, and, indeed, no contradictions seem to have been found to date. However, the experiments that led to this postulate were conducted on only a single protein, bovine ribonuclease A (RNAse). We conduct molecular dynamics (MD) simulations on this protein with the aim of mimicking this experiment as well as making the methodology available for use with basically any protein. There have been many attempts to model denaturation and refolding processes of globular proteins in silico using MD, but only a few examples where disulphide-bond containing proteins were studied. We took the view that if the reductive deactivation and oxidative reactivation processes of RNAse could be modelled in silico, this would provide valuable insights into the workings of the classical Anfinsen experiment.     

Gastric flow and mixing studied using computer simulation

The fed human stomach displays regular peristaltic contraction waves that originate in the proximal antrum and propagate to the pylorus. High-resolution concurrent manometry and magnetic resonance imaging (MRI) studies of the stomach suggest a primary function of antral contraction wave (ACW) activity unrelated to gastric emptying. Detailed evaluation is difficult, however, in vivo. Here we analyse the role of ACW activity on intragastric fluid motions, pressure, and mixing with computer simulation. A two-dimensional computer model of the stomach was developed with the 'lattice-Boltzmann' numerical method from the laws of physics, and stomach geometry modelled from MRI. Time changes in gastric volume were specified to match global physiological rates of nutrient liquid emptying. The simulations predicted two basic fluid motions: retrograde 'jets' through ACWs, and circulatory flow between ACWs, both of which contribute to mixing. A well-defined 'zone of mixing', confined to the antrum, was created by the ACWs, with mixing motions enhanced by multiple and narrower ACWs. The simulations also predicted contraction-induced peristaltic pressure waves in the distal antrum consistent with manometric measurements, but with a much lower pressure amplitude than manometric data, indicating that manometric pressure amplitudes reflect direct contact of the catheter with the gastric wall. We conclude that the ACWs are central to gastric mixing, and may also play an indirect role in gastric emptying through local alterations in common cavity pressure.     

The effect of amino acid substitution on protein-folding and -unfolding transition studied by computer simulation

Protein-folding and -unfolding transitions were studied by the method of computer simulation. The protein was modeled as a two-dimensional lattice polymer. Various energy terms were assumed to be operative between units composing the polymer. But hydrophobic interactions were neglected explicitly. Both thermodynamic and kinetic quantities were obtained from the simulation, and from their temperature dependence in the transition zone characteristics of the conformational transition of proteins were discussed. Two amino acid substituted models, differing in the location of substitution, were studied and compared with the original in order to clarify the effect of substitution on conformational transition of proteins. The following conclusions were reached in this study: (1) The relaxation time of the slow mode, which reflects the overall folding and unfolding processes, shows a peak near the transition temperature, while that of the fast mode is almost independent of temperature. The peak of the slow mode occurs at a slightly lower temperature than the transition temperature. (2) The dependence of the logarithm of the rate constants on the inverse of temperature (Arrhenius plot) is linear. Therefore, the plot of the free energy of activation vs temperature is linear. (3) The values of kinetic parameters obtained suggest that in the activated state the intramolecular interactions are half broken, while the state is close to the native state on the entropy axis. (4) The amino acid substitution, which is modeled as having slightly unfavorable short-range interactions, causes the substituted ones to be slightly unstable. Moreover, it causes the folding transition to slow. From the analysis of the way slowing down is observed in the two substituted models, we conclude that a structure, designed to model a β-sheet, is formed before it gets assembled with other structures, which are designed to model α-helices. The process of assembly occurs nearly at the activated state of the folding and unfolding transition. (5) It is suggested from this study that the maximum of folding rate constant in the Arrhenius plot that has been observed experimentally in real proteins is likely due to hydrophobic interactions.     

Peripheral binding mode and penetration depth of cobra cardiotoxin on phospholipid membranes as studied by a combined FTIR and computer simulation approach

Cobra cardiotoxin, a cytotoxic beta-sheet basic polypeptide, is known to cause membrane leakage in many cells including human erythrocytes. Herein, we demonstrate that the major cobra cardiotoxin from Naja atra, CTX A3, can cause leakage of vesicle contents in phosphatidylglycerol (PG) and phosphatidylserine containing, but not in pure phosphatidylcholine (PC), membrane bilayers. By the combined polarized attenuated total reflection infrared spectroscopy and computer simulation studies, CTX A3 is shown to peripherally bind to both zwitterionic and anionic monolayers in a similar edgewise manner with a tilted angle of approximately 48 +/- 20 degrees between the beta-sheet plane of the CTX molecule and the normal of the membrane surface. The average surface area expansion induced by CTX A3 binding to the PG monolayer, however, is two times larger than that of the PC monolayer as determined by the Langmuir minitrough method. Interaction energy considerations of CTX A3 on neutral and negatively charged membrane surfaces suggests that the electrostatic interaction between anionic lipid and cationic CTXs plays a role in modulating the penetration depth of CTX molecules on the initial peripheral binding mode and reveals a pathway leading to the formation of an inserted mode in negatively charged membrane bilayers.     

Dynamics of folding and unfolding transition in a globular protein studied by time correlation functions from computer simulation

A method of calculating time correlation functions from records of computer simulated equilibrium conformational fluctuations in a globular protein is discussed. Use of the calculated time correlation function for discussions of dynamics of folding and unfolding transition in the two-dimensional lattice model of proteins. The time correlation functions can be approximated in general by a sum of two simple exponential terms. The relaxation time of the slower mode does not depend on the nature of the physical quantity with respect to which the time correlation function is calculated. This time characterizes the overall folding and unfolding transition. The relaxation time of the faster mode depends on the nature of the physical quantity and characterizes conformational fluctuations within each of the native and denatured states. The mechanism of a previously observed phenomenon of the acceleration of the folding and unfolding transition by short-range interactions is discussed.     

Formation of the chick primitive streak as studied in computer simulations

We have used a computer simulation system to examine formation of the chick primitive streak and to test the proposal (Wei and Mikawa Development 127 (2000) 87) that oriented cell division could account for primitive streak elongation. We find that this proposal is inadequate to explain elongation of the streak. In contrast, a correctly patterned model streak can be generated if two putative mechanisms are operative. First, a subpopulation of precursor cells that is known to contribute to the streak is assigned a specific, but simple, movement pattern. Second, additional cells within the epiblast are allowed to incorporate into the streak based on near-neighbor relations. In this model, the streak is cast as a steady-state system with continuous recruitment of neighboring epiblast cells, egress of cells into deeper layers and an internal pattern of cell movement. The model accurately portrays elongation and maintenance of a robust streak, changes in the composition of the streak and defects in the streak after experimental manipulation.     

Superhelical DNA studied by solution scattering and computer models

We present here recent results on the structure of superhelical DNA and its changes with salt concentration between 0.01 and 1.5 M NaCl. Scattering curves of two different superhelical DNAs were determined by static light scattering. The measured radii of gyration do not change significantly with salt concentration. Small-angle neutron scattering, together with calculations from a Monte Carlo model, allows to determine the superhelix diameter. Measured and simulated scattering curves agreed almost quantitatively. Experimentally we find that the diameter decreases from 16.0±0.9 nm at 10 mM to 9.0±0.7 nm at 100 mM NaCl. The superhelix diameter from the simulated conformations decreased from 18.0±1.5 nm at 10 mM to 9.4±1.5 nm at 100 mM NaCl. At higher salt concentrations up to 1.5 M NaCl, the diameter stays constant at 9 nm.This revised version was published online in October 2005 with corrections to the Cover Date.     

Thyroid homeostasis in animals with iodine deficiency

The content of total iodine, its hormonal and nonhormonal fractions as well the level of protein-bound iodine in blood and basic tissue targets in representatives of 4 classes of animals: Esox lucius L., Rana esculenta, Streptopelia decaocto Priv., Lepus europaeus Pall. inhabiting the mountain regions with iodine deficiency in environment and in the lowlands of Transcarpathia with higher iodine provision have been investigated. A considerable decrease of general and hormone iodine level in the animal tissues of the mountain area accompanied by the suppression of the thyroid function has been stated. The utilization of thyroid hormones under the iodine deficiency condition is increased in the majority of cases and the level of protein-bound iodine is lowered that testifies to the transition of animal organism in the iodine-deficient areas to the lower level of thyroid homeostasis.     

Cold-active enzymes studied by comparative molecular dynamics simulation

Enzymes from cold-adapted species are significantly more active at low temperatures, even those close to zero Celsius, but the rationale of this adaptation is complex and relatively poorly understood. It is commonly stated that there is a relationship between the flexibility of an enzyme and its catalytic activity at low temperature. This paper gives the results of a study using molecular dynamics simulations performed for five pairs of enzymes, each pair comprising a cold-active enzyme plus its mesophilic or thermophilic counterpart. The enzyme pairs included α-amylase, citrate synthase, malate dehydrogenase, alkaline protease and xylanase. Numerous sites with elevated flexibility were observed in all enzymes; however, differences in flexibilities were not striking. Nevertheless, amino acid residues common in both enzymes of a pair (not present in insertions of a structure alignment) are generally more flexible in the cold-active enzymes. The further application of principle component analysis to the protein dynamics revealed that there are differences in the rate and/or extent of opening and closing of the active sites. The results indicate that protein dynamics play an important role in catalytic processes where structural rearrangements, such as those required for active site access by substrate, are involved. They also support the notion that cold adaptation may have evolved by selective changes in regions of enzyme structure rather than in global change to the whole protein. Figure Collective motions in C^α atoms of the active site of cold-active xylanase Electronic supplementary material Supplementary material is available in the online version of this article at and is accessible for authorized users.     

Thermostability of protein studied by molecular dynamics simulation

The thermostability of protein thermostable cathechol 2,3-dixoygenase (TC23O) has been studied by the parallel molecular dynamics simulations. By analysis of the exponent beta, which is related to the scattering spectrum and constant-pressure heat capacity Cp, we reveal the respective contribution of a specific residue 228 proline; a specific salt bridge, Lys188N-Glu291OE1; four ions; and a different water environment to the thermostability of TC23O. The dynamic transition temperature of the mutants, Pro228Ser and Glu291Gly of the TC23O, was decreased about 10 degrees C and 19 degrees C respectively. The displacement of the four ions had no significant effect on the thermostability of TC23O. Water affects the thermostability by influencing the changes of accessible conformation to a certain extent. All these results agree with the known experimental results.     

Stability of polypeptide conformational states as determined by computer simulation of the free energy

A method is developed to extract the entropy of polypeptides and proteins from samples of conformations. It is based on techniques suggested previously by Meirovitch, and has the advantage that it can be applied not only to states in which the molecule undergoes harmonic or quasiharmonic conformational fluctuations, but also to the random coil, as well as to mixtures of these extreme states. In order to confine the search to a region of conformational space corresponding to a stable state, the transition probabilities are determined not by “looking to the future,” as in the previous method [H. Meirovitch and H. A. Scheraga (1986) J. Chem. Phys. 84 , 6369–6375], but by analyzing the previous steps in the generation of the chain. The method is applied to a model of decaglycine with rigid geometry, using the potential energy function ECEPP (Empirical Conformational Energy Program for Peptides). The model is simulated with the Metropolis Monte Carlo method to generate samples of conformations in the α-helical and hairpin regions, respectively, at T = 100 K. For the α-helix, the four dihedral angles of the N- and C-terminal residues are found to undergo full rotational variation. The results show that the α-helix is a more stable structure than the hairpin. Both its Helmholtz free energy F and energy E are lower than those of the hairpin by ΔF ～ 0.4 and ΔE ～ 0.3 kcal/mole/residue, respectively. It should be noted that the contribution of the entropy ΔS to ΔF is significant (TΔS ～ 0.1 kcal/mole/residue). Also, the entropy of the α-helix is found to be larger than that of the hairpin. This is a result of the extra entropy arising from the rotational freedom about the four terminal single bonds of the α-helix.     

Conformational equilibria of valine studied by dynamics simulation.

The conformational probability distribution of a valine residue in the valine dipeptide and of the valine side chain in an alpha-helix, as well as the change in helix stability for replacing alanine with valine, has been calculated by molecular dynamics simulations of explicitly hydrated systems: dipeptide, tetrapeptide and 10-, 14- and 18-residue oligoalanine helices. All computed free-energy differences are means from at least eight separate slow-growth simulations, four in each direction and are reported with their root-mean-square deviations. Different values for the change in free energy of folding (delta delta G degrees) have been calculated with the use of forcefields having an all-atom and a central-atom representation of methyl groups, etc. The value obtained with the all-atom forcefield agrees well with new experimental values (3 kJ/mol = 0.7 kcal/mol). Furthermore, the most stable valine side-chain rotamer in the helix is different for these two representations. The most stable rotamer for the all atom conformation is the same one that predominates for valines in alpha-helices in proteins of known conformation. The lower conformational freedom of the valine side chain in the helix contributes 1 kJ/mol to the difference in stability computed with the all-atom potential; unfavorable interactions of the side chain with helix, even in the most stable conformation, further increase delta delta G degrees.     

Lipid bilayer polypeptide interactions studied by molecular dynamics simulation

A model membrane with a polypeptide alpha-helix inserted has been simulated by molecular dynamics at a temperature well above the gel/liquid crystalline phase transition temperature. Order parameters of the lipids and other equilibrium and dynamic quantities have been calculated. Three systems, polyglycine constrained into an alphahelical configuration, glycophorin with similarly conformationally constrained backbone and finally glycophorin free to change its backbone conformation, have been studied. In all cases there was an ordering of the chains close to the helix. This effect was, however, much smaller for glycophorin with its rather bulky side chains than for polyglycine. The dynamics of the lipids were affected by the neighbouring helix, not drastically however. Lateral diffusion and reorientational time correlations of lipids close to the helix were slower than for the bulk ones, but not more than two or three times. Thus, we did not find any evidence of bound or frozen boundary lipids.     

Acetylcholine receptor pore permeability studied by molecular dynamics simulation

A dynamic model of the closed-state pore of an acetylcholine receptor (five M2 α-helices stabilized with a (CH₂)₁₀₅ ring) is used to examine the migration of uncharged and charged probe particles equivalent to a hexahydrated sodium ion (van der Waals diameter 7.27 Å) propelled by varied external force along the channel axis. Ion movement through the pore is hindered by steric constraints and electrostatic interactions. The van der Waals gate is formed by helix residues 13′ (A-Val255, B-Val261, C-Val269, D-Val255, and E-Ile264), whereas the negatively charged residues in the upper part of the channel are important for ion selectivity.     

Formation of the branching pattern of blood vessels in the wall of the avian yolk sac studied by a computer simulation

The present study was performed to provide data to support the notion previously believed but not proved experimentally or theoretically, that blood vessels are formed by the selection of capillaries in the network. In an attempt to understand the mechanism of formation of blood vessel branching structures, the transformation of a capillary network to a branching system in the wall of quail yolk sac was successively recorded by a series of photographs, and a computer simulation was carried out for the process of in vivo vascularization based on the photographs. The simulation demonstrated that a positive feedback system participated in the formation of a branching structure. That is, vessels which had been much used were enlarged, whereas less used vessels were reduced in their size and finally extinguished. The enlarged vessels became major components of the branching system. As the body of an embryo grew, it was observed that polygonal capillary networks enlarged, which led each polygon of the network to divide into a few finer polygons. Then, some of the capillary vessels were again selected and formed a branching system. This process repeated during the body growth, indicating that the vascular system developed adaptively to the body growth. A region where the growth was fast, received much blood flow and produced finer networks of capillaries. Thus, it was experimentally demonstrated for the first time that capillaries in the network are successively selected by a positive feedback mechanism and form blood vessels.     

Dioxygen affinity in heme proteins investigated by computer simulation

We present an investigation of the molecular basis of the modulation of oxygen affinity in heme proteins using computer simulation. QM-MM calculations are applied to explore distal and proximal effects on O(2) binding to the heme, while classical molecular dynamics simulations are employed to investigate ligand migration across the polypeptide to the active site. Trends in binding energies and in the kinetic constants are illustrated through a number of selected examples highlighting the virtues and the limitations of the applied methodologies. These examples cover a wide range of O(2)-affinities, and include: the truncated-N and truncated-O hemoglobins from Mycobacterium tuberculosis, the mammalian muscular O(2) storage protein: myoglobin, the hemoglobin from the parasitic nematode Ascaris lumbricoides, the oxygen transporter in the root of leguminous plants: leghemoglobin, the Cerebratulus lacteus nerve tissue hemoglobin, and the Alcaligenes xyloxidans cytochrome c'.     

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)