Nano-modeling of insulin-like growth factor 1 (IGF-1) by computational methods

The empirical force fields have great difficulty in simulating folding of insulin-like growth factor 1 (IGF1). In an effort to understand the conformational preferences that may be attributed to stereoelectronic effects, a number of computational studies are carried out. Monte Carlo, molecular dynamics and Langevin simulation methods by MM+, amber and optimized potential for liquid simulations (OPLS) force fields of calculations have been performed on IGF-1 as growth factor. The parameters of minimized structure of IGF-1, calculated potential energy for important dihedral angles and the effect of temperature on geometry of optimized structure have been calculated. In this work, we have used different temperatures at gas and water media and we have seen that in simulation approaches, scaling up the interaction energy has a similar effect to lowering temperature. This study has demonstrated that the simple model including an approximate average solvent effect can simulate the qualitative feature of the IGF-1. The key research was to find dynamics of biomolecular structure and an appropriate effective stabilized energy.


INTRODUCTION
Insulin-like growth factor (IGF-1), also known as somatomedin C, mediates the growth promoting activity of growth hormone.IGF-1 is autocrine regulator of cell proliferation, paracrine growth and survival factor for mammalian embryo development (Emmitte et al., 2009).Recent NMR studies have revealed that IGF-I has three -helical regions surrounding a hydrophobic core (Laajok et al., 2000).The over expression or auto activation of the insulin-like growth factor-1 receptor (IGF-1R) tyrosine kinase has been associated with various cancers.Insulinlike growth factor (IGF-1) is an anti-apoptosis factor of multiple cell types and the anti-apoptotic effects are mediated through mitochondrial and cytochrome-c pathway (Li et al., 2003).Development of faster computers that are within the reach of the widest scientific community as well as efficient computational *Corresponding author.E-mail: m_monajjemi@yahoo.com.methods allows investigating systems between 50 to 100 atoms in the frame of quantum mechanics and up to 50,000 atoms with molecular dynamics.Since the models become increasingly realistic, direct comparison with experimental data becomes possible (Na´ray-Szabo´ and Berenteb, 2003).In addition to hit identification, docking techniques are increasingly used to support lead optimazation efforts (Kitchen et al., 2004).Recently, constant temperature molecular simulations of peptide folding have been reported using implicit solvent models and explicit solvent models (Monajjemi et al., 2006;Dzubiell et al., 2006).Recently, however, several computer simulations have demonstrated a strong coupling between hydrophobicity, solute-solvent dispersion attractions and electrostatics.For example, simulations of explicit water between plate like solutes revealed that hydrophobic attraction and dewetting phenomena are strongly sensitive to the nature of solute-solvent dispersion interactions (Monajjem et al., 2006).
The competing effects of the solvent such as the van der Waals (VDW) attraction and hydrogen bonding between the protein and solvent reduce the strength of the interactions and consequently reduce the energy barrier related to the multiple minima problem (Ozkan et al., 2004).In solution, the intramolecular VDW interactions of a protein molecule are balanced by the intermolecular VDW interactions with solvent molecules.The possible difference between the protein intramolecular VDW attraction and that with water may be included in the hydrophobic interaction energy (Ozkan et al., 2004).Kurochkina and Lee have shown that the pair wise sum of the buried surface area is linearly related to the true buried area (Sung, 1999).Since the specific interactions between the residues and solvent play an important role in the stability of the native structure, it is useful to carry out such simulations at atomistic detail.This comes with the problem of timescale of folding/ unfolding that is several orders of magnitude larger than those currently attainable by MD simulations (Ozkan et al., 2004).Water plays a crucial role for the stability, dynamics and function of proteins.For this reason molecular dynamics (MD), Monte Carlo (MC) and Langevin dynamic (LD) simulations must account for the effects that this solvent has both on protein structure and on protein dynamics (Hamaneh and Buck, 2007).The aim of the present work was to describe and characterize the molecular structure vibrational properties IGF-1 crystalline-structure.In this work, the structures of a coordination compound modeling the IGF-1 computationally.Thus, it is worthwhile to collect information on their structures by the means of computational chemistry as well.

THEORETICAL BACKGROUND AND COMPUTATIONAL METHODS
The crystal structures of proteins were from the Brookhaven protein data bank.The structure of protein IGF-1 was selected from the protein data bank (PDB code 1B9G).These studies provided insights into the steric, electrostatic, hydrophobic and hydrogen bonding properties and other structural features influencing the IGF-1.In vacuum the system was simulated using Monte Carlo, molecular dynamic and Langevin dynamics with 100 ps step and without any constraints.Temperature was kept constant at 300 K.In water, simulations, the system was placed in a box (3 x 3 x 3 nm) containing one molecule of solute and 884 TIP3P water molecules (Figure 1).The system was simulated using Newtonian dynamics with 100 ps step and no constraints applied to the solute.

MC simulation
Monte Carlo simulations are based on pair wise additive potentials of the form (Monajjemi et al., 2008): Where rij is the distance between atoms i and j, Aij and Bij are coefficients associated with the particular atom pair and qi and qj are the partial charges associated with each of the atomic sites.
Each distinct atom i in the system is assigned a set of parameters Aii, Bii and qi.The coefficients Aij and Bij can then be obtained from the mixing rules Aij = (AiiAjj) 1/2 and Bij = (BiiBjj) 1/2 (Monajjemi et al., 2008).

MD simulation
In concepts and algorithms of classical MD simulations the atoms of a biopolymer move according to the Newtonian equations of motion (Berendsen, 1990): ..
Where m is the mass of atom , r is its position, and Etotal is the total potential energy that depends on all atomic positions and, thereby, couples the motion of atoms.For an all-atom MD simulation, one assumes that every atom experiences a force specified by a model force field accounting for the interaction of that atom with the rest of the system.
The electrostatic potential energy is represented as a pairwise summation of Coulombic interactions as described in Equation 4( Phillips et al., 2005): In Equation 4, N is the number of atoms in molecules A and B respectively, and q the charge on each atom.

VDW interaction
The van der Waals potential energy for the general treatment of non-bonded interactions is often modeled by a Lennard-Jones 12 to 6 function as shown in Equation 5 (Kitchen et al., 2004): In Equation 5, is the well depth of the potential and is the collision diameter of the respective atoms i and j (Kitchen et al., 2004).We can consider an effective Hamiltonian operator constructed for molecule in a given geometry of it and the solvent: Where H0 is the Hamiltonian in gas phase (the unperturbed Hamiltonian), Velec is the perturbation from the permanent charge distribution of water, represented as a set of point-charges, Vind is the perturbation from the induced dipoles in the solvent and Vnon-elec is a non-electrostatic perturbation which models the effect of the anti-symmetry between the solute and solvent (Na´ray-Szabo´ and Berenteb, 2003).

Langevin dynamics
Using Langevin dynamics, you can model solvent effects and study the dynamical behavior of a molecular system in a liquid environment.These simulations can be much faster than molecular dynamics.These simulations can be used to study the same kinds of problems as molecular dynamics: time dependent properties of solvated systems at non-zero temperatures.Because of the implicit treatment of the solvent, this method is particularly well-suited for studying large molecules in solution.Langevin dynamics simulates the effect of molecular collisions and the resulting dissipation of energy that occur in real solvents, without explicitly including solvent molecules.This is accomplished by adding a random force and a frictional force to each atom at each time step.Mathematically, this is expressed by the Langevin equation of motion (Berendsen, 1990): Here, is the friction coefficient of the solvent in units of ps -1 and Ri is the random force imparted to the solute atoms by the solvent.
The friction coefficient is related to the diffusion constant D of the solvent by Einstein's relation: The random force is calculated as a random number, taken from a Gaussian distribution, with an average value of zero and no correlation with the atom's velocity.Molecular mechanics (MM) force fields rely on the combination of Coulomb and Lennard-Jones interactions to describe all nonbonded interactions (Ponder and Case, 2003).Even though the functional form of the potential energy is quite simple, it depends on a large number of empirical parameters which must be obtained from ab initio electronic structure calculations on small molecules and/or experimental data.
Because each new term in the MM potential function requires additional empirical parameters, it is quite appealing to keep the functional form of the potential function as simple as possible.While most widely used current force fields such as amber, OPLS do not employ explicit hydrogen bonding terms, this was not always the case (Monajjemi and Chahkandi, 2005;Hagler and Lifson, 1974;Cornell et al., 1995;Monajjemi et al., 2005;Weiner et al., 1984).

RESULTS AND DISCUSSION
The complex was solvated by added water molecules.The systems were first energy minimized steps with the conjugate gradient algorithm.Then, the positionrestrained MC, MD and LD simulation were run 100 ps afterwards, 1 ps simulations were carried out at a time step of 100 ps (Figure 1).Several simulations were carried out as listed in Table 1.MC, MD and LD simulations of the IGF-1 were performed with the Hyper-  HyperChem, 2001).The geometries, and the interaction energies, bonds, angles, stretch-bends, electrostatic and the VDW interactions were carried out in solution and gas phase (Table 1 and Figure 2).In solution, the intramolecular VDW interactions of a protein molecule are balanced by the intermolecular VDW interactions with solvent molecules.Thus, when solvent molecules are not explicitly included, the intramolecular VDW interactions must be adjusted accordingly.The longer-range attractive VDW interactions provide a nearly uniform background potential (Chandler et al., 1983) and therefore can serve as the reference for the VDW energy calculation (McCammon et al., 1980).The following text describes methods for generating and evaluating representative molecular conformations, particularly for peptides and small proteins, based on 'molecular mechanics' energy functions.On the other hand, 'molecular mechanics' describes molecules as atoms linked with springs (harmonic bond stretches and bond angle wagging), each atom having finite volume and relatively sharp boundaries ("6 to 12" hard spheres potentials), with sinusoidal torsional energies.The force field for a typical protein can be given as a sum of the various components including bond stretching and bending, torsional potentials and non-bonded interactions.In this paper, we have used Monte Carlo methods to study IGF-1 in the bulk and in confined environments.Results are presented in Table 1 effects on the specific media of the structure.The potential energy was for the growth factor (IGF-1) with water during the MC simulation is shown in Figure 3. Molecular dynamics simulations were carried out on the two systems, gas and solvent IGF-1 molecule.All simulations were performed at 300 K.Each solvent system was immersed in a periodic water box and the structures of water molecules were maintained.A 100 ps time step was used in all the simulations.The potential energy, represented through the MD "force field," is the most crucial part of the simulation since it must faithfully represent the interaction between atoms yet be cast in the form of a simple mathematical function that can be calculated quickly.The system was well equilibrated and 500 ps in the range of the MD equilibration were selected for further processing analysis.After equilibration, the MD simulation was very stable and in order to compare the difference between the relation coefficients (R2 = 0.8173 in gas and R2 = 0.7558 in water), we have shown in Figure 3 respectively.The theoretically possible stable conformers of free molecule were searched by means of a molecular dynamics calculation performed in a temperature interval from 0 to 500 K; for example, the iterative calculation with time step of "100 ps", carried out by utilizing the software "Chem3D", the experimental x-ray geometrical data reported for IGF-1 in crystalline structure were used as input geometrical data (Hartung et al., 1992).At the next step, the appropriate ones carefully selected from the structures obtained throughout this calculation were optimized using MM+, amber and OPLS force field parameters included into the same software.In this paper, a comprehensive conformational search on free molecule was carried out.The obtained results have demonstrated that the free molecule has a very flexible macro-cyclic structure.On the basis of the theoretical results obtained for the determined most stable, the dependencies of the geometrical and force constants parameters of the free molecule to its conformational structure were discussed.Furthermore, we have used MD and LD methods to study protein in the bulk and in confined environments.The structures obtained throughout this calculation were optimized using MM+, amber and OPLS force field parameters.Also, all these  approaches were included discrete particles moving in a defined energy landscape according to Langevin dynamics (LD).These results have shown that the force field of AMBER has convenient relation in all of simulation methods and various media (Figure 3).
From the simulations we have shown that the kinetic temperature of the system is properly bounded around the prescribed equilibrium temperature.The length of each simulation was 100 ps.We have measured the relative drift of molecular temperature denoted by T in percent with respect to mean temperature, T in Kelvin.In the simulation of the small water system with temperature of 295, 297, 299, 301, 303 and 305 K.The instantaneous kinetic temperature is given by Equation 8 (Español and Warren, 1995): Where k B is the Boltzmann constant, N f is the degrees of freedom (N f = 3N − 3 for a system of N particles with fixed total momentum), m i is the atom weight for atom i, v i is the velocity of atom i.We block-average over many instantaneous values to get an accurate estimate of the temperature.These methods, which rely upon uniform sampling of energy space, can yield thermodynamic data over the entire temperature range of interest and have been shown to overcome large free energy barriers.
We have reported findings for six different temperatures of various sizes and topologies.Results are presented in Table 2 that were indicated potential energies of IGF-1 in various temperatures.For certain confining environments, individual proteins do exhibit power-law dependence, but the relationship is different for each molecule.In other cases, the increase in stability upon confinement interestingly demonstrates nonmonotonic behavior.Several molecular dynamics simulations could be performed over a wide range of temperature and the data could be combined using a weighted histogram approach (Weiner et al., 1984); however, the statistical error associated with the tails of the sampled distributions is usually large and can propagate when data from simulations at different temperatures are merged.Potential energies for the three force fields of MM+, amber and OPLS at Monte Carlo simulation were compared in Figure 4.The average energies are in good agreement within the simulation accuracy.As expected, amber demonstrates much smoother energy profiles than the other two simulation methods due to higher-order energy conservation in themodified Hamiltonian (Figure 4a).The magnitudes of energy fluctuations in both MM+ and OPLS approaches are significantly smaller than the other (Figure 4b and c).The sampling results of step-size of MD and LD methods are presented in Figures 5 and 6 respectively.Observed data are almost identical for both choices of the MD simulation length, which suggests that the MD simulation have affected much more the acceptance rate at least for this particular model than MC and LD approaches.This potential does not have any terms describing angular dependencies of hydrogen bonds and is similar to the 10 to 12 hydrogen bonding potential originally proposed by McGuire et al. (1972).They found that hydrogen bonding energies were represented adequately by a sum of Lennard-Jones and electrostatic interactions plus the 10 to 12 hydrogen bonding term with empirical constants adjusted according to the hydrogen bond type.
It was because the functional form of such a hydrogen bonding term was very close to the Lennard-Jones component of the force field, the second-generation amber force field omitted it altogether (Cornell et al., 1995), relying instead on the combination of Lennard-Jones and Coulomb interactions to model hydrogen bonded complexes, thus the data of this force field in three simulation methods have shown the changes of potential energy via time at various temperatures more better than MM+ and OPLS force fields (Figures 4a,5a and 6a).Similarly, the widely used OPLS force field does not contain an explicit hydrogen bonding term: the emphasis of OPLS parameterization is on reproducing thermodynamic properties of organic liquids such as enthalpies of vaporization, densities and free energies of hydration (Monajjemi et al., 2005;Jorgensen and Tirado-Rives, 1988) (Figures 4b, 5b and 6b).Because each new term in the MM+ potential function requires additional empirical parameters, it is quite appealing to keep the functional form of the potential function as simple as possible (Figures 4c, 5c and 6c).The effect of confinement on the thermodynamic properties of several statement proteins was investigated by performing simulations over a large range of temperatures.We have computed the transition temperature for the IGF-1 molecule.The results are summarized in Table 1 for 300 K in gas or solvent and in Table 2 for 295, 297, 299, 301, 303 and 305 K temperatures.Figures have shown the function of the reduced temperature.Low reduced temperatures promote complex structure stability, whereas high reduced temperatures oppose it.The major part of this difference is due to the interaction of IGF-1 with solvent molecules correspond to various simulation methods and force fields.A difficult task in computational study of stabilized structure is to find a proper energy function that can lead to a unique structure.
Our simulations showed that the simple energy function modified to include solvent effect has a parameter range that can simulate indicated structure at constant temperature of 300 K.

Conclusion
In this work we have used molecular dynamic models to explore the stability of IGF-1 by comparing theoretical methods of simulation.A highly selective on effect of temperature and environment was discovered in chemical structure and it has been investigated the standard constant temperature at MC, MD and LD simulations.We have employed the molecular dynamics simulation method as the main tool to study conformational dynamics of biomolecules.One of the force field designed for treating macromolecules can be simplified by not considering explicitly -the so-called united atom approach is 'amber'.It was appeared that solvent effects influence the calculated potential energy surface, by lowering potential energy barriers on angle.
This means that the parameterizations that have been developed for small molecules with considerable effort can be carried over into macromolecular calculations with little or no change.Also, we have applied MM+ and OPLS force fields parameters for IGF-1 model in gas and for water environment.Also, the possible difference between IGF-1 intramolecular VDW attraction and that with water has been included in the hydrophobic interaction energy.The shortrange repulsion represents the exclusive volume of each atom and needs to be calculated explicitly.The measurement of the potential of solvation under similar conditions of temperature in solution along with investigation of energetic and structural aspects of solution have been used to gain insight into the molecular level interaction with IGF-1.
Solute-solvent pair interaction of potential energies was shown that the greater stability of solvent observed over all states investigated in this study is related to the MD/amber approach.

Figure 2 .
Figure 2. Geometry optimized variables of bond length (B), bond angle (A) and dihedral angle (D) in gas and water media at 300 K.

Table 1 .
Calculation various variables in 300 K temperature for IGF-1 at MM+, amber and OPLS.