Science topic

Molecular Dynamics Simulation - Science topic

Explore the latest questions and answers in Molecular Dynamics Simulation, and find Molecular Dynamics Simulation experts.
Questions related to Molecular Dynamics Simulation
Question
6 answers
How can I do MMPBSA analysis on Gromacs trajectory using the MMPBSA tool of amber software packge?
Relevant answer
Answer
Not sure it is possible, since these two codes use different output formats for their trajectories. An alternative, I think that you can convert the gromacs trajectory files (xtc or trr) into AMBER trajectory files with VMD. You will also need to convert your GROMACS topology file into the AMBER file format.
HTH
Question
11 answers
I would like to calculate auto-correlation function for tilt angle for the lipids in a bilayer from molecular dynamic simulation. I have the angles for all lipids from each frames for all subsequent frames but I need a code to calculate the auto-correlation function. Where can I find or get the code apart from writing my own?
Relevant answer
Answer
Autocorrelation function is implemented in many (commercial/open) software.
Depending on the MD code you have used, there may be a tool already included like 'g_analyze' which comes with GROMACS: http://manual.gromacs.org/online/g_analyze.html
Question
9 answers
I have simulated the protein in water and 2) protein with ligand
now I want to check the g energy in gromacs for both system individually and then combine to one plot. How could I do this?
Relevant answer
Answer
With grompp, you are not performing a simulation run, you just generate the input file for the simulation. So proceed by running e.g.
mdrun -s md_0_1.tpr -deffnm md_out
and you will get a "md_out.edr" file for your simulation run, which you can then analyze with g_energy.
Question
2 answers
How and with what program can I calculate the MD simulation of the proteins modified by polyGlu-tails of different lengths. Will suitable the Gromaks to solve such problem? Thanks for the help!
Relevant answer
Answer
Hi,
In addition to Hiqmet's response. It is also depend (a little) of the force field you want to use to model your poly Glu tails. Did you read the literature to see if somebody has do similar simulations?
HTH
Question
1 answer
Is there a tutorial for this and for using viparr.py?
Relevant answer
Answer
It's in the "Desmond 2.2 User Manual, Chapter 8: Using Alternate Force Field Parameters and Constraints". You need to run the viprr.py script in command line using cms files as input.
Question
2 answers
I am doing MD simulation of ag-ab complexes using GROMACS and want to know the eletrostatic potential of a surface atom in a coordinate. Please tell me if there is any method to calculate the same that would be helpful to me.
Relevant answer
Answer
Your question is not clear, but I think you mean how to calculate the electrostatic potential at a certain coordinate R originating from certain atom r_i. If this is the case, you must take the vector difference (R - r_i), then calculate the length of that vector, call that d_i. The electrostatic potential is (1/(4 Pi e0))*(q_i/d_i). Where e0 is the dielectric constant, and q_i is the point charge assigned to atom i by GROMACS. Then you must add over all atoms you want (sum over index i).
Question
3 answers
I tried to do MD refinement for a protein-membrane Bilayer complex by using Gromacs. In manual they mentioned only distance restraints between two atoms with in the molecule but here it is a complex with two molecules. Anyone could provide references to this topic is really appreciation.
Relevant answer
Answer
One method is to rotate the bilayer so it is in the x-y plane . A harmonic restraint can then be applied in the z direction. From the Gromacs 4.6.1 manual, section 4.3.1
"Using three different force constants the position restraints can be turned on or off in each spatial dimension; this means that atoms can be harmonically restrained to a plane or a line. Position restraints are applied to a special fixed list of atoms. Such a list is usually generated by the pdb2gmx program"
Question
2 answers
I'm getting two different types of electrostatic potential maps in presence, and in the absence of water molecules.
Relevant answer
Answer
APBS calculate electrostatic energy by using the implicit solvent model. So, you do not need to include water molecules explicitly. APBS uses dielectric constant of water to calculate the electrostatic potential. Even in absence of any water molecules, obtained electrostatic potentials of solute is assumed to be similar to that of the solute present in the bulk water. If water molecule will be present explicitly, certainly it will affect the electrostatic potential, because APBS will consider water molecules as a solute similar to the protein/ligand.
Question
1 answer
Can anyone advise how to modify the water model and which force field is used?
Relevant answer
Answer
Hello
For simulating DPPC there wont be any inbuilt forcefield available in gromacs.
You have to either manually add these lipid parameters to an existing forcefield or can download some of the user submitted forcefields containing lipid parameters submitted in gromacs website.
And as far as water model is concerned this link might help you
Good luck !!
Question
2 answers
#T PM6 Geom=ModRedundant ADMP
...
1 F
2 F
3 F
...
the modredundant option was read but the atoms still move.
Relevant answer
Answer
Does the distance between two frozen atoms remain constant throughout the trajectory? If so, the atoms are actually frozen and just translated during the trajectory.
If this is not the case, there are a few ways to freeze atoms in g09 using cartesian coordinates. Although I have no experience in ADMP calculations in g09, you may try to see if any of these options work for your purpose:
1) try changing Geom=ModRedundant with Geom=AddRedundant: the specification of the frozen atoms is the same. From the g09 Manual (http://www.gaussian.com/g_tech/g_ur/k_geom.htm):
"AddRedundant is synonymous with ModRedundant. This option may be used for job types other than optimizations."
2) try changing Geom=ModRedundant with Geom=ReadOptimize or Geom=ReadFreeze: you can find more on this keyword in the Manual (http://www.gaussian.com/g_tech/g_ur/k_geom.htm). In my experience, Geom=ReadFreeze works with ONIOM calculations, even though it's said to be deprecated in the Manual, while Geom=ReadOptimize doesn't work (at least in G09 rev A.02).
3) try using the freeze-code -1 in the molecule specification after the element label (http://www.gaussian.com/g_tech/g_ur/m_molspec.htm).
Hope this helps,
Daniele
Question
6 answers
I have simulated the protein-ligand system by using the Gromacs force-field 43a1 after obtaining the docking results for semi-flexible and flexible docking.
My ligand is a steroid, and they only work related to this published have used amber package, which is not free.
Can I continue my analysis by using this Gromacs 43a1 force-field?
Relevant answer
Answer
The AMBER force field can be used with GROMACS (and it is included in the standard distribution), so you do not need to use the AMBER package if all you need is the same force field. It is usually considered among the best force fields for protein molecular dynamics.
The GROMOS force fields are generally considered (from what I've heard from many, many users: I didn't compare them myself) to be not as good as AMBER, but they can be useful nonetheless.
What is critical however is the *parametrization* of your ligand. It has to be consistent with the force field you use and it has to be of good quality. Therefore your choice of force field will depend on how you can parametrize your ligand. Do you already have a known, good topology for the ligand? If yes, use the force field compatible with that topology. If not, look around to create an AMBER or GROMOS compatible ligand topology (I am unaware of tools able to create OPLS-friendly topologies but they may well exist).
Question
40 answers
Can anyone suggest a web based server or free software which is comparable to GastroPlus (GastroPlusTM, a physiologically-based simulation tool that predicts the absorption, pharmacokinetics, pharmacodynamics, and drug-drug interactions for drugs administered through intravenous, oral, ocular, or pulmonary routes)?
Relevant answer
Answer
Can you please send me the link in sbhandari157388@troy.edu
Sagar.
Question
17 answers
I have four proteins A, B, C and D. A is a transmembrane receptor and its the substrate that gets phosphorylated in its cytosolic domain. B is the kinase, therefore, Mg++ and ATP dependent. B also depends on C and D for its activity. There might as well be a domain swapping between B and C. Is it at all possible to build a reliable model of the ternary complex of BCD? I only have the sequence data. Direct structural information of these proteins is not available.
Relevant answer
Answer
This would be a tall order. Like Felix said, you're best hope would be to find a model of a homologous multi-protein complex that you could use as a template for building a predicted model. If such a model does not exist, you could NOT perform de novo predictions of how a multi-protein complex assembles, even in the case where a reasonable model for the structure of each subunit exists.
A fruitful approach for modeling multi-protein complexes that cannot be studied directly by X-ray crystallography has been to perform docking studies in cryo-EM images of multiprotein complexes. The minimum criteria you need are high-resolution cyro-EM, and structural models for each subunit.
Question
51 answers
I've been told to use Xeno View to build my polymer and assign the Force Field parameters. The problem is that I want to use my own set of parameters and it seems that the program can only read parameters from a default list of frc file. Can you suggest a software that lets you build a polymer from a monomer, assign the atom types, add water and save the result in a topology file containing the parameters I need?
Relevant answer
Answer
Some ideas that I am aware of ...
Commercial -- These tools allow you to either build or import individual molecules, assign force-fields, build complex molecular systems (initial condition) from the individuals, and then automate the running of simulations.
Materials Studio: http://accelrys.com/
Of these, I have experience with the Scienomics software. It is very user friendly, has excellent molecular construction tools including a polymer builder and a crystal builder, great force-field support, and interfaces with a lot of different open source MD and Monte-Carlo programs, LAMMPS, GROMACS, etc.
They all probably have academic licenses.
Free
I have not personally used this software, but understand it has some, but not all of the capabilities of the commercial entities above. Has its own MD engine and great for GPU is my understanding.
Question
2 answers
I am just new to gromacs. Can anyone please help me in correcting these charges for the itp file, i.e. how to adjust them? These have been generated by prodrg, I have read the article for correcting it but couldn't make it. At further steps the error is generated while running the energy pressure and temperature.
Relevant answer
Answer
You can use SWISSPARAM for generating topology for your ligands. But it uses charmm 27 forcefield. Its mentioned on the site how to use this topology with gromacs. Hope this helps you.
Energy minimztion probelm in Gromacs.
Question
14 answers
When I run the energy minimization step for Groamcs it gives the following error: Double precision normally gives you higher accuracy. You might need to increase your constraint accuracy, or turn off constraints alltogether (set constraints = none in mdp file) Steepest Descents converged to machine precision in 43 steps, but did not reach the requested Fmax < 1000. Potential Energy = -1.4074464e+06 Maximum force = 8.7960381e+05 on atom 2201 Norm of force = 3.9092959e+03 writing lowest energy coordinates. How can I correct this?
Relevant answer
Answer
Hello. First, you have to try to understand why that is happening. 1) Are you sure that the protein is well built (i.e. no bonded atoms inside aromatic rings, clashes, etc? 2) What are you minimizing? Protein, DNA? Small molecules? Please carefully check the structure. Open the failed minimized GRO file and check for distortions or any information that can help you to evaluate the initial structure. Besides steepest descent method, you can try minimizing the structure by conjugated gradient method (CG), a more thorough method. Simply change integrator to cg (inside minim.mdp): integrator = cg Best regards, Ricardo Ferreira
Question
1 answer
I am planning to make a molecular dynamics simulation by embedding my protein into nuclear membrane, to study interactions in active site. But I couldn't find any solution to put lipids different for nuclear membrane rather than cellular membrane. If possible I am interested in building system with Maestro to use Desmond.
Relevant answer
Answer
Hi,
I think the differences between cellular and nuclear membrane is only for their lipid compositions. So in the simulations point of view, you can use the same approach to embedded the MP in the membranes.
To construct your models of membrane and insert your membrane protein in it, you could use CHARMM-GUI (http://www.charmm-gui.org/)
AFAIK Desmond supports the CHARMM force field
HTH
Plz suggest some best workstation for modeling studies ?
Question
5 answers
Hello all, I wanted to purchase workstation for modeling study. If any one suggest some software names/ workstation suitable package for docking and molecular modeling in one. Thanks
Relevant answer
Answer
Hi look at these workstation from Dell.... They may be costly but service is very good....... May be you can also try IBM workstations http://www.dell.com/in/business/p/precision-t7600/pd http://www.dell.com/in/business/p/precision-t5600/pd http://www.dell.com/in/business/p/precision-t3600/pd
Gromacs run problem
Question
65 answers
I want to do 3ns simulations using gromacs on an i3 laptop. My problem is that I have to face the load shedding problem. After every 2 hours, the light goes for 1 hour, sometimes the light goes for more than an hour and a half and the laptop battery diminishes. I have heard that the graph will be effected if I shut the lid when only 10% battery is left. Is there any solution for this problem? Can I split my simulations into different parts so that they are not effected by the light problem? If so, how can I do this? Or is there any other way to solve my problem?
Relevant answer
Answer
yes, in Gromacs mdrun will save your run every 15 mins. please look into documentation for mdrun "-cpt" property of mdrun will help you to restart your simulation, when power goes off or PC crashes. just append -cpt to your mdrun command more info: http://www.gromacs.org/Documentation/How-tos/Doing_Restarts all the best
Question
3 answers
I want to buy a workstation for bio molecules simulation for minimum of 100ns. Can any one suggest the best configuration?
xenon E5 1650, 3.02 Ghz,12 Mb cache, 12 GB of Sata RAM, DDR3, NVIDA Quardro 4000 2Gb
Relevant answer
Answer
If you want to use GPU for simulation, then buy Nvidia Tesla. Cards like GTX or Quadro is nice for testing and debugging purpose. Simulations may run continuously for many days and weeks, which may affect the hardware. GTX or Quadro are not made for heavy continuous use, while simulation is long and computationally expensive process.
Question
8 answers
I would like to build a polymer membrane and do MD studies. Is there any free software available to build a polymer membrane and do MD simulation? In the future I want to study the interaction of polymer membranes with biomolecules. What is the procedure to build a model in gromacs?
Relevant answer
Answer
Hi Javier,
Thanks for answer. But VMD can only be used for bio membrane. Do you think for polymer membrane also U can use it. Thanks
Can anyone suggest how to perform RMSIP analysis over the gromacs trajectories ?
Question
8 answers
I want to perform the RMSIP analysis over gromacs trajectories, which I came across in one of my reference papers. Is there any specific tools or scripts for performing this analysis
Relevant answer
Answer
You can calculate inner product of covariance matrix by g_anaeig tool using -inpr option. You can also calculate overlap in covariance matrix using same tool. This article may be helpful. http://pubs.acs.org/doi/full/10.1021/jp062548j
Question
2 answers
I would like to setup unrestrained MD simulations on gromacs for protein-protein complex . The protein complex was calculated by HADDOCK by using NMR Chemical shift restraints?
Could you please suggest, or give the reference to, unrestrained MD simulations?
Relevant answer
Answer
Hello.
I think you are a little bit confused. All MD production runs are unrestrained, i.e. the atoms spatial position is allowed to move.
However, in all MD simulations, the first steps (1-10 ns) must be restrained in order to allow a correct temperature (NVT) and pressure (NPT) equilibration. Only then, the restrains can be removed and the system is allowed to equilibrate on a MD production run.
Mount your GROMACS system, do an energy minimization EM) followed by a short NVT to equilibrate pressure (10-100 ps) and a restrained NPT (2-10 ns) to allow pressure equilibration. Only after this you can remove restrains and perform an unrestricted MD run.
Best regards,
Ricardo Ferreira
Is it possible to perform a molecular dynamics simulation (GROMACS, AMBER ...) at a pH other than 7?
Question
16 answers
I want to do a molecular dynamics simulation with my protein. In vitro experiments were performed at pH 8. Is it possible to alter the parameters of the molecular dynamics properties with explicit solvent at pH=8? Or do I have to alter the protonation of my protein in a way that it mimics a pH of 8?
Relevant answer
Answer
My experimental colleagues have no (little) difficulty in adjusting the pH of their solutions, as pH is a macroscopic parameter, like temperature and pressure. At the microscopic level of MD simulations, the macroscopic pH will be reflected in your choice of protonation states for titratable residues. It is generally not possible to provide much in the way of computational insights into pH-dependence of reactions. The pH is an ensemble average value and you are typically simulating systems with a copy number of one. If your concern is about a particular residue in the active site, say a histidine, which is much less likely to be protonated at pH 8 than at pH 6.5, then you will need to run simulations with it singly protonated and others with it doubly protonated. Even at pH 8, where the probability of protonation is reduced by an order of magnitude from pH 7, remember that there are probably 10^15 copies of your protein in a small vial and millions of those histidines will remain protonated even at high pH.
Question
3 answers
I'm studying stability of certain enzymes with their mesophilic and thermophilic homologues using MD Simulation methods. Using the information obtained from the study I'm trying to design point mutations which would help in increasing stability. Is it possible to suggest some appropriate methods of studying such mutants quantitatively. e.g. in terms of melting temperatures or something like that?
Relevant answer
Answer
A similar discussion is in the link below
You may also take a look at this paper for some ideas on how structural factors affect thermal stability
DNA simulations with groamcs
Question
20 answers
I want to simulate DNA in water. All gromacs force fields do not support DNA. So I chose the amber force fields from Gromacs package. I want to ask, are all the remaining steps the same for DNA simulations like the lyzome example? I have done up to this step pdb2gmx -ff amber99sb -f DNAAmber.pdb -o DNA2.pdb -p DNA.top -water spce -ignh editconf -bt octahedron -f DNA2.pdb -o Complex.pdb -d 1.0 genbox -cp Complex.pdb -cs spc216.gro -o Complex_b4ion.pdb -p DNA.top grompp -f em.mdp -c Complex_b4ion.pdb -p DNA.top -o Complex_b4ion.tpr genion -s Complex_b4ion.tpr -o Complex_b4em.pdb –nname Na –nn 8 –g trp_ion.log For the rest of the steps, can I continue from: http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/lysozyme/05_EM.html I am stuck at this step
Relevant answer
Answer
Yes, parmbsc0 is still not included in new version. So, you have to add it yourself. You can use following command to add 0.150 M NaCl : genion -s Complex_b4ion.tpr -o Complex_b4em.pdb -conc 0.150 -neutral It will add sufficient number of Na+ and Cl- ions such that ion concentration will be 0.150 M and whole molecular system will be neutral.
Question
9 answers
I am working with an enzyme which has Zinc in its catalytic site.
I want to further continue my study using molecular dynamics simulation, as the zinc is involved in the catalysis process. I tried the simulation with dummy atoms approach but the mass gets distributed between dummy atom (DZ) and zinc, and has an effect on catalysis.
So if anyone has any idea about zinc parameter please let me know.
Relevant answer
Answer
There are parameters for Zn in the CHARMM, AMBER, and GROMOS force fields. There is no zinc in the OPLS force field files that come with GROMACS, so I don't know about it. If you are really interested in the reaction, you will have to do QM/MM and will then not need the parameters since the zinc will be in the QM region.
Question
10 answers
RMSD calculation using g_rms in Gromacs. My plan is to compare between deviation of the protein with respect to crystal structure. However I am a bit confused.
Gromacs ask for:
1) group for atleast square fit
2) group for RMSD calculation.
If I want to plot backbone vs carbon alpha group, which one should I pick so that the results can be compared?
Square fit | RMSD calculation
c-alpha --------c-alpha
c-alpha --------backbone
or
c-alpha--------c-alpha
backbone-----backbone
Relevant answer
Answer
g_rms indicates: Each structure from a trajectory (-f) is compared to a reference structure. The reference structure is taken from the structure file (-s).
So you need to put the crystal structure under the flag -s and the structure that you want to compare with the crystal structure under the flag -f.
Then, when you select the groups, the first that the tool asks for is the group for least square fit (so, the structures under the -f flag) and the second is the group for RMSD calculation (structure under the -s flag).
For better information you can see section 8.9 of Gromacs Manual.
The g_rms is a graphic of RMSD vs time. You should also see root mean square fluctuation (g_rmsf), perhaps this tool is what you are looking for: http://manual.gromacs.org/online/g_rmsf.html
Question
6 answers
Using NAMD, I want to simulate aquaporin channel in POPC membrane. To induce the water flow, pressure difference across the membrane is needed. I am stuck with this point. Could you please suggest how to do it?
Relevant answer
Answer
Actually using PBC You cant just use higher pressure on one side of membrane. Recently I heared a conference talk, where guy did it by putting two membranes in system - obtained a kind of pseudo-liposome where water pressure between two membranes was higher than outside.
Question
2 answers
I want to simulate a protein using PDB file extracted from Protein Data Bank. After inspection, I noted that the protein atoms in the PDB has ANISOU atoms. In my understanding, this type of atoms has two coordinates as they are flexible (correct me if I am wrong).
My question is, can I run a molecular dynamics on this PDB file? Do I need to do any pretreatment?
Relevant answer
Answer
Hi
you have to actually remove the ANISOU atoms data from the PDB file before going for simulation.
You can use the following python script to remove them
Question
50 answers
I did the PCA of one protein after MD simulation and now want to visualize the movement of protein in VMD using PorcupinePlot.tcl., but I am not able to.. If any one knows how to do this please share.
Thanks
Jitendra
Relevant answer
Answer
Hi Jitendra
See this crude video that i made to demo the steps involved to you
The video is crude. u need to pause and read,
the step or also outlined below
BTW: modevector in pymol gives better quality images then Porcupine.tcl
There are more then one way of doing this. I am describing one method here.
Requirement:
-------------
You should have downloaded a Porcupine.tcl from
How to use Porcupine.tcl and generate a porcupine plot
------------------------------------------------------
Step 1: Get extreme1.pdb using the following command
"g_anaeig -s protein_CAlpha.pdb -f combinedItajectory.xtc -extr extreme1.pdb -first 1 -last 1 -nframes 50"
Step2: Using a text editor, save the first frame in "extreme1.pdb" to "first_frame.pdb"
Similarly, do for last frame and save it to "last_frame.pdb"
Step3: Open VMD, and Load first_frame.pdb in New Molecule.
Similarly load, last frame as a newmolecule.
Step4:
Open Tcl/Tk Console
Step5: type "play Porcupine.tcl" and enter
Step6: Then in VMD Main, click on first_frame and type in console : atomselect top "name CA"
Step7:
Then in VMD Main, click on last_frame and type in console:
atomselect top "name CA"
Step8: Then type
porcupineplot::Diff atomselect1 atomselect0
Output looks like this:
-----------------------
>Main< (ras-raf) 17 % atomselect top "name CA"
atomselect0
>Main< (ras-raf) 20 % atomselect top "name CA"
atomselect1
>Main< (ras-raf) 21 % porcupineplot::Diff atomselect1 atomselect0
running through 242 residues
drawing arrows for distances between 0.14418759278366697 and 5.825031451230942 (colors from scheme BGR)
drawing color scale bar
done
view the plot in VMD display window.
Question
8 answers
I want to simulate Piroxicam (NSAID drug) molecules with GROMACS. I have the GROMACS topology file generated by PRODRG server. As mentioned here http://pubs.acs.org/doi/abs/10.1021/ci100335w the topology generated by PRODRG server contains numerous inaccuracies. I want to validate this topology. How can I do this?
Relevant answer
Answer
I think the easiest way is try to construct a gromacs topology file using pdb2gmx , using your PRODRG generated file as input and selecting one of the gromos forcefield options ( as the paper that you cite is based on this ff). If works, will mean that the topology generated with PRODRG is valid for run a simulation with gromos.
But without doing this, and only reading the abstract of this paper, I think that you are not going to be able, at elast not in a good way, but in the paper are proposing several points as “best practices” for parametrization of small molecules under the GROMOS force fields, try to take advantage of them.
Question
3 answers
Please consider the adsorption of some organic molecules at the surface of semiconductor.
Can the dipole dipole interaction between the adsorbent molecules affect the surface coverage or adsorption energy of the system?
Any reference or suggestion will be appreciated.
Thanks in advance.
Relevant answer
Answer
Dear Brian and Manoranjan
thanks for your good references
Question
2 answers
Is there any correlation between the value of Delta G and the stability of the protein complex?
Relevant answer
Answer
The higher, in absolute value, of delta G, the better the energy of the protein-ligand complex. However, it doesn't correlate with the activity presenting different ligands. A priori, what should correlate, is if you have a ligand A, and the ligand B, the difference in whose delta G (AAG a-b = | dG(a) - dG(b)| ) should correlate better with the activity. No values a priori can be told, as depending of the size of the binding site and the size of your ligand, more contacts, hence, better or worse dG
If the complex protein-ligand is stable, then, molecular dynamics is required. What you need to do is monitorize the RMSD of your ligand, but it won't give so much information.
Essentially, you need to combine dG and RMSD information both together. You could have very good binders that present a good energetic profile (dG) but are less stable (RMSD high) and viceversa. At the end, both descriptors could give you a hint about if the ligand is stable and energetically well stabilized by the protein
Question
11 answers
I am a new user of GOMACS. My system contains Polymer and Water. I used gromos53a6 force field and spc water model and it ran nicely. But now I want to use TIP3P water model with the same force field (gromos53a6), in that case, I am getting a error message - "Fatal error: Atomtype OWT3 not found".
In addition, I added OWT3 and HW at the bottom of the .atp file and also corrected nonbonded.itp file of the corresponding force field. But I couldn't get rid of such kind of error.
Would you please give me suggestions to fix it?
Thanks in advance ...
Relevant answer
Answer
Juganta,
you are right. Foolish me to mix atomtype with atomname... dah!
Nonetheless, I've tried and gave me the same error. This means that GROMOS forcefield in GROMACS isn't capable of using TIP3P because, in ffnonbonded.itp, there isn't any OWT3 atomtype.
You can add it, of course, but is never a good idea to mix/change/alter forcefields if you are planning to publish your work later on.
In GROMACS, the only FF that have TIP3P are GMX, GROMOS43A1 and OPLS. This information is given by /opt/gromacs-4.6.3/share/gromacs/top/tip3p.itp:
; This file is for backward compatibility only.
; Please directly include the tip3p.itp file from your force field directory.
#ifdef _FF_GROMACS
#include "gmx.ff/tip3p.itp"
#endif
#ifdef _FF_GROMOS96
#include "gromos43a1.ff/tip3p.itp"
#endif
#ifdef _FF_OPLS
#include "oplsaa.ff/tip3p.itp"
#endif
Hope it helps this time!
Ricardo
Question
4 answers
I need the complete time series of a protein dynamics (hopefully also provided by someone wanting to collaborate). The main aim of our group is to extract topological descriptors of the corresponding protein contact network to be compared with energetic variables of the MD.
Relevant answer
In which kind of protein are you interested in? Globular protein, small peptide, membrane protein?
Question
1 answer
Which of the potential energy surfaces is the best interatomic potential to describe graphene when using LAMMPS?
Relevant answer
Answer
Dear Mohsen
The adaptive intermolecular reactive bond order (AIREBO) potential is able to model covalent C-C bond breaking with associated changes in hybridization of atomic orbitals. In addition, the AIREBO potential reproduces the bonding characteristics (such as binding energy and bond length) and forces associated with the rotation about dihedral angles for carbon double bonds. However, If you are trying mechanical study of carbon family, the potential file should be slightly modified or it predict large strength before fracture.
Question
3 answers
The average structure calculated during analysis would be the average of all the atomic positions/co-ordinates during the length of simulation. But, I wanted to know how are the bond angles treated for computation of the average structure? Since, we use the structure for further analysis too. Are the bond lengths distorted or are they retained while the average is being computed? Also, what are the co-ordinates used for this computation?
Relevant answer
Answer
@Tarak: Thank you for the explanation. So, if my force field parameters are correct, then the average structure computed via g_rmsf fucntion of gromacs should give me a structure which obeys the standard chemical conformations, i.e. the bond angles and bond lengths? Also, could you please guide so as to how to check if the average structure is correct as obtained?
Question
5 answers
The MD simulations using gromacs 96 force field was done on a compound obtained after semi-flexible docking and flexible docking. The protein was also simulated in water. Comparing the docking results with the original protein structure the total energy was high, so should I conclude that after docking the structure was unstable?
Relevant answer
Answer
Hello.
It's rather difficult to say. It could be yes, but it could also mean that the binding pose and the sidechain orientation, although good for the ligand, it's not ideal for the system. It will depend on the system, because it is to be expected an energetic barrier, but it cannot be calculated like that. Only with metadynamics or free energy profiles.
Best regards,
Ricardo
Question
3 answers
Can anyone comment on the reliability of the parameters generated from SWISSPARAM server?
I read few of them, saying that the parameters are good enough to start with.
Thanks
Relevant answer
Answer
Like all theoretical models, they are good enough to start with, but if you need accuracy then you need to measure it yourself. If you are unsure about relative reliability of various computational methods, it's common to start with the most popular (so even when it's wrong it can be compared to everyone else who used the same software), and then get some experimental data to refine your estimate. Just remember, all models are wrong, but some are useful.
Question
34 answers
I am using Gaussian for molecular optimization and energy calculation but now I want to learn how to do the same for dimers. In the case of dimers, orientation and different types of interaction play a major role in optimization of system, so how to deal with all these parameters?
Relevant answer
Answer
Dear Nikhil,
don't be to affraid of studying these systems. Thanks to DFT-D3 it is much easier nowadays then some years ago. In the beginnig I would suggest you prepare some input structures with different intermolecular orientations (angle, slide, distance), this will help you to see if you end up in the same minimum and to compare the energies. Regarding solvation, why not starting in vacuo to get a feeling for these calculations and later include (maybe different) solvent models? HF will not do a good job for vdw / dispersive interactions, here you need at least MP2 in a reasonable basis set including CP, but forget about that with G..
With DFT you can use BLYP with Grimmes Dispersion correction and Density Fitting to make it faster and stepwise increase the basis set size. Counterpoise correction is not recommended for DFT-D, use the QZVP basis set for the final single points, this will give you essentially no BSSE. If you want to give a try to a different software you can try Turbomole, where each system will take some days on a single core...
Best wishes,
Malte
Question
4 answers
Fatal error:
There is a dangling bond at at least one of the terminal ends. Select a proper terminal entry.
Relevant answer
Answer
I don't know what software you are using but it sounds like there is an atom with an unfilled valence (i.e. not enough bonds attached to it, but also hasn't been assigned a negative charge). Or it might be that it is expecting a specific functinality on the atom but sees a hydrogen instead (i.e. it expects methyl group but there is a hydrogen instead.
It all depends on the specific software you are using and the test it is performing that your molecule is failing, but certainly at first glance it appears as though part of you molecule is incorrect in some way.
Question
5 answers
I have seen in a lot of papers the probability distribution of a particular eigenvector and comparison with Gaussian PC for the same set of analyses. They have projected the histogram of the fluctuations. I was trying to get the same for my set, I do have the rmsf projection along individual eigenvectors, say 1-8 but I was little confused with the histogram of the fluctuations as to how to obtain those.
Also, for obtaining the Gaussian PC from a list of eigenvectors, is there a straightforward measure (which I have obviously missed reading) to find which one is Gaussian in my set, or does one have to go through hit and trial along with the eigenvalue spectrum? Say, choose PC 100 and plot the projection and obtain the histogram and then see.
Relevant answer
Answer
Dear Ankita,
sorry for the delay in answering your question. I had a very busy schedule in the past weeks. I took a look to the paper you mentioned and now I have a better understanding of what you need.
Firstly, about producing histograms from the fluctuations around some projection over principal coordinates you have, of course, a number of possible solutions. Usually I produce histograms with Octave (a free version of Matlab). Say that you have an ASCII text file ("data.txt") containing 2 rows of numbers. The first being time, the second being the rmsf of the projection along some eigenvalue. These are the basic commands to plot and save an histogram (in what follows ">>>" means the Octave prompt). I will assume that you launch Octave from the same directory where you have the "data.txt" file.
>>> a = load("./data.txt");
>>> [h, x] = hist(a(:,2),50,1.0);
>>> plot(x,h);
>>> M= [x' h'];
>>> save hist.dat M;
The command "hist" creates an histogram using data contained in the second row of the matrix a [a(:,2)], with 50 bins and normalizing to 1. I suggest you to refer to Octave user guide about the hist command.
The plot command just plots for you the histogram.
The following commands are used to save your histogram into a file named "hist.dat", which is a 2-column ASCII file with the first row being the value of the rmsf associated to the n-th bin (n = 50 in the example), while the second row contains the value of the histogram in that bin.
About the "Gaussianity" of the distributions, I think that there is no other way than trial-and-error. There are different ways to establish if the distribution is Guassian or not. One possibility is to fit the histogram with a Gaussian as both Haiguang and I suggested before. Here I suggest you a solution using Octave (not gnuplot as I wrote in the previous comment).
This is the idea: if your variable has a Gaussian distribution
p(x) = N exp[-(x-x0)^2/2s^2]
where N is the normalization constant, x0 the average value of x and s^2 the standard deviation, then the natural logarithm of p(x) is a second order polynomial
ln(p(x)) = A + Bx + Cx^2
where A, B and C trivially depend upon N, x0 and s^2. What you can do, after having produced the histogram (and so the x and h arrays) with the previous Octave commands is to do a polynomial fit or order 2 of ln(p), recalculate the value of the distribution with the fitting parameters and obtain a chi-square value. The closer the chi-square is to 0, the closer your distribution is to a Gaussian one. Here are the Octave commands:
>>> z = polyfit(x, log(h + 1.0e-15), 2);
>>> y = poly ], which allows one to treat different portions of a molecular system with different ab initio many-body methods. With CIM you choose the portion that is chemically active to be treated with the highest level ab initio many-body method, while the rest that is not chemically active, and will not play a significant role in determining the correlation energy, is treated with a lower level method. A good example of how this can be used for real systems is P.M. Kozlowski, M. Kumar, P. Piecuch, W. Li, N.P. Bauman, J.A. Hansen, P. Lodowski, and M. Jaworska, \The Cobalt-Methyl Bond Dissociation in Methyl-cobalamin: New Benchmark Analysis Based on Density Functional Theory and Completely Renormalized Coupled-Cluster Calculations," J. Chem. Theory Comput. 8, 1870-1894 (2012). There are some nice illustrations that demonstrate how the chemically active portion is chosen, and most importantly, it is chosen not by the user, thus making it easier to use. I know this is different from ONIOM, but the thinking for choosing the higher-level and lower-level portions are the same. Which portions are chemically important and which are mainly spectators, if you will. Also, if you're interested, CIM will be freely available in GAMESS most likely in the next version that is released.
Question
5 answers
Why is it used in the negative and not the positive of gradient?
Relevant answer
Answer
All of these answers are correct.
There is a very simple way of articulating this convention. Objects roll down-hill, the negative sign just says that: Force is in the direction of DECREASING potential energy; this convention makes potential energy diagrams intuitive!
I hope Luis finds this slightly helpful.
Question
11 answers
The protein is of 543 amino acids, its MD simulation is done using Gromacs, the problem is occurring when the duration is set for 1000 PS and contemporary steps. The problem occurs in 457 step showing "molecules too much far in the trajectory". Could any one suggest how to overcome this issue?
Relevant answer
Answer
try looking into the Gromacs users mailing list (http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List) - it offers solution to many simulation issues, and advice both from the users and developers
Question
7 answers
I don't know how to find the 0001 and 1012 hematite surfaces. Where I find the surfaces? I need to draw them with Avogadro or some draw tool?
Relevant answer
Answer
Actually I think you can do it, why not. The only problem is to find a suitable force field to reproduce hematite accurately. It is a pretty common iron oxide, so you will probably find something in the literature.
Regarding the surfaces, yes, you should create them with a visualization tool or with a drawing tool. The second one is easier of course. I recommend you the free software GDIS, or the commercial code Material Studio if you have access to it.
Question
15 answers
Gromacs, namd.. Etc
Relevant answer
Answer
In terms of how the force fields compare w.r.t structure and dynamics, rather than speed factors, there have been a few recent articles that compared a number of the more popular force fields and compared them against experiment.
Here's a couple:
doi: 10.1021/ct300323g
doi: 10.1021/ct2007814
Question
40 answers
I have installed gromacs on ubuntu 12.10, and when I run the command
/usr/share/gromacs/tutor/gmxdemo
It says it's a directory
But when I write the following comand by sudo
/usr/share/gromacs/tutor/gmxdemo/demo, it says command not found....
And by cd result is permsion is denied,
How can I run the demo for gromacs correctly?
Relevant answer
Answer
Hello.
Enter the directory and run the command as follows:
$ ./demo
Question
2 answers
How can I measure the formation of water mediated hydrogen bonds formed between a ligand and a protein complex during simulation. I'm using gromacs 4.5v. Is there any specific way by which we can calculate the number of hydrogen bonds?
Relevant answer
Answer
Water mediated hydrogen bonds are very difficult to monitor in a molecular dynamics simulation. 1) They are much less stable (around 4 kcal.mol⁻¹) than protein-protein hydrogen-bonds (~8 kcal.mol⁻¹), and 2) as they are more labile, they break more easily and are rapidly replaced by another water molecule. You have to select a particular residue and all water molecules to find out how many HB are established, but you will not find which specific water molecules (and how many) are involved). If your water molecules are relatively stable and 'trapped' between ligand and protein, maybe you can check them. If not, is more accurate to chech ligand-protein interactions (with LigPlot software) that focus only in water molecules...
For further informations on g_hbond, check: "Van der Spoel, D., Van Maaren, P. J., Larsson, P., & Tîmneanu, N. (2006). Thermodynamics of hydrogen bonding in hydrophilic and hydrophobic media. J Phys Chem B, 110(9), 4393–8. doi:10.1021/jp0572535"
Question
9 answers
I began to run a MD simulation using the GROMACS program. In the first step, pdb2gmx, it was asked to select the force filed. How do I select the right force field?
Relevant answer
Answer
Many force fields can be used to study protein-ligand interactions, for example CHARMM with CGenff, AMBER with GAFF, OPLS, etc.). The choice of the force field for your simulations *should not* depend on the MD code. It is particularly true if the program you want to use accepts different potentials (such as GROMACS, NAMD). To choose the better force field for your simulations, you will need to answer some questions and also read the literature !!! For example:
- Does my ligand or my protein have been already studied with MD simulations. if yes, with which force field ? If not, can I used (say GAFF and AMBER parameters to model accurately my ligand and my protein, respectively?) If yes, good !! , Try these force fields in your simulations and see if you obtain accurate results with them. If not, try different combination of force fields . If all the parameters available in literature are not able to model accurately your protein-ligand interactions, you will need to parameterize your own force field. It is a long and complicated task !!!
HTH
Question
53 answers
Please tell me tools for molecular dynamics simulation studies. Software for Windows and Linux platform, free and commercially available.
Relevant answer
Answer
Hello.
Check GROMACS (www.gromacs.org). In my profile, you will find an article about membrane protein simulation (DOI: 10.1021/ct300083m), more specifically ABC transporters (P-glycoprotein) on a POPC membrane. This paper have all the steps/methods you will need.
Any help on membrane building or any other issue, please ask.
Best regards.
Question
2 answers
I want to simulate dobule stranded DNA using Gromacs. I read that I should use an amber forcefield for DNA simulation and I have crossed a website explaining some parameters that need to be changed in the Gromacs package to become compatible with this amber force filed. Is it possible to use the pdb2gmx command on PDF file for DNA (taken from pdb data bank) and choose amber forcefield for that? Are those changes necessary? And if pdb2gmx does not give error, how can I verify the correctness of the procedure?
Relevant answer
Answer
Make sure that the double strands do not come apart. Some DNA simulations will require constraints to keep the two strands together.
DSSP instalation error.
Question
3 answers
I want to do secondary structure analysis using DSSP. I am getting an error during the execution. error command: Failed to execute command: /usr/local/bin/dssp -na dd02brnL dd1YGM3g > /dev/null 2> /dev/null. . Please help me out.
Relevant answer
Answer
Use -v, that gives you a verbose output. Try to put the extension (dd02brnL.pdb). You have dssp installed, right?
Question
5 answers
Materials studio 3.4
Relevant answer
Answer
Hi, I don't think it is possible, since MS 6.1 supports only COMPASS, CVFF, PCFF, Dreiding, and Universal, and thus probably your MS version is too old for using coarse grains force fields.
HTH
Question
6 answers
I will be glad to learn from anyone who has had experience (or knowledge) of how to estimate configurational entropy from MD trajectories. In particular, I will like to obtain an estimate from cluster analyses obtained from trajectories from different force fields.
Relevant answer
Answer
Dear Olujide,
I have recently written a tool to calculate angular spread of single residues from MD trajectories In fact the output is a "flexibility spectrum", where the flexibility of each residue is used to characterize the given trajectory. The tool can also distinguish local fluctuations from conformational changes. Let me know if you want to try it.
Question
1 answer
I have a question about Material Studio. Actually I'm a beginner to deal with this software. I would like to use forcite calculation to do my simulation. I have a problem with this calculation. I would like to use COMPASS as a force field but one of my compounds (Boron) can't be found in this force field. Then, I tried to change the force field using Universal Force Field and it works. My questions are:
1. What is the different between using COMPASS and Universal Force Field for MD simulation?
2. If I want to use COMPASS Force Field, it means I have to add the compound (Boron) to the database. How to add the database of this compound? Is there any equation that we have to add?
3. Do you have any suggestions on how to do this simulation? Should I change to another method (GULP or the others)?
Relevant answer
Answer
A force field is a set of force constants obtained from experimental data.
Usually this set is calculated from the vibrational spectra of a set of chosen molecules so the resulting force field is limited to the same family of molecules.
UFF has ben developped in a different way: see the article of Rappe et al 1992.
Question
6 answers
is there any server which uses amber force field to make gro file. acpype commnad can convert amber topology to gromacs. If theres another direct way to deal with this
Relevant answer
Answer
U must first install ambertools (its fortunately free!). Then make your .prmtop and .crd files using tleap or xleap (included in Ambertools). then use this file to convert .prmtop and .crd files to .top and .gro files
How to retrieve the structures for the energy minima from the generated Gibbs free energy landscapes in Gromacs?
Question
20 answers
I have retrieved the Gibbs free energy landscapes for my principal components and located the energy minima. Can anyone suggest how to retrieve the structures from the specific energy clusters on the free energy landscape? I have attached the paper which I am using as reference.
Relevant answer
Answer
Thanks Rajendra !!
Question
23 answers
I want to simulate a modeled protein which is around 7 angstrom from native and I want to see whether it goes towards the native with time. What could be the parameters so that the protein does not change much?
Relevant answer
Answer
As you probably known folding of small protein with classical explicit MD simulations (without use enhanced sampling approach, such as T-REMD, metadynmics, etc..) take a (very) long time and need a lot computational resources.
For me, It normal that RMSF of C and N terminal are highers, since generally these parts are more flexible than other part of the protein. To conclude I will suggest you to do more refinements of your modeled protein and carry out other MD.
HTH
Based on what criteria should we select different temperature levels while simulating a point mutation to study the folding pattern?
Question
27 answers
I want to perform MD simulations for two point mutations of my protein at two or three different temperature levels and compare the folding at different temperatures. I want to know whether there is any criteria for selecting different temperatures. Any suggestions?
Relevant answer
Answer
You could predict the temperature stability change of your mutated protein and choose temperatures for simulation below and above the predicted melting temperature. Tools that predict the melting temperature change according to mutations are: Auto-Mute: http://proteins.gmu.edu/automute/AUTO-MUTE.html Prethermut: http://www.mobioinfor.cn/prethermut/
Question
3 answers
How can i measure the hydrophobic interaction between the ligand and protein using gromacs
Relevant answer
Answer
Hello.
Between the ligand and protein, unfortunately gromacs don't have any tool. But, you can download LigPLOT (or LigPLOT+ with GUI), a simple program that quantifies the number of non-bonded interactions and hydrogen-bonds between any ligand and the protein.
Question
2 answers
In my study I want to simulate the solvent effect on Diels Alder reaction by MD method using different solvents.
Relevant answer
Answer
I don't have expertise on SYBYL, but in GROMACS you cannot break chemical bonds (as seen in a chemical reaction). The only way you can simulate a similar effect is a multi-step quantum calculation with GAUSSIAN or employing a QM/MM simulation, with explicit solvent that you can change. But, in QM, solvent effects can be difficult to simulate.
Question
3 answers
I am new to MD simulation for Protein - protein complex , sorry for asking some basic questions !
1. What Kind of information it will provide if I do MD simulation for protein complex after docking ?
2. How to analyses results if i use MD softwares like AMBER , NAMD or Gromacs ?
3.How to compare these simulation result with Experimental result ( if any ) ?
Relevant answer
Answer
Dear Mr. Penumutchu, hear are my answer for your questions. Hope these answer satisfies you.. :)
1. If you conduct a MD simulation for protein complex after docking, you may predict its stability based on the curve of simulation time versus potential energy of the complex. This complex stability calculation thus very important to obtain the best ligands or lead compounds.
2. To answer this question maybe you have to see the software manual. Usually, to analyze the MD simulation results you have to draw curve between potential energy of the complex versus simulation time or RMSD/RMSF (root mean square deviation/fluctuation) versus simulation time.
3. You could compare the RMSD/RMSF curve between simulation result and experimental result from protein data bank (PDB) and then you could conclude is your simulation result valid or not.
Question
9 answers
I am doing molecular dynamic simulation of one compound to see the binding mode in its protein target, but the absolute configuration of one chiral center of this compound changed during the dynamic simulation. How can I fix its absolute configuration to be unchanged?
Relevant answer
Answer
the chiral molecule is not corrected parameterized. Therefore, you must be modified parameters. For example, you can use swissparam.ch, such as describe by Ashar.
Question
8 answers
Protein Simulation is giving too many LINCS errors and the program terminates in 14 ps. I am using parameters as in sample mdp files.
Relevant answer
Answer
GO through tutorial it will helpful for you
Question
1 answer
I want to do a simulation of this protein (metallprotease) with triterpenoid saponin inhibitors with antihemorragic properties from Pentaclethra macroloba.
Relevant answer
Answer
it seems that there is no structure available for your protein in pdb database, but if you make a blast on its sequence, you will find this .
which could be a good start for homology modelling due to the 60 % of sequence identity from your protein : blast scoring : 497 evalue : 1E-140
I join attached files , 3DSL , and homology construction of your prot ( made with Chimera and Modeller V9.11)
I hope they help you
Question
2 answers
I have the data for a gene network obtained from HumanNet but I am unable to interpret the scores been given from HumanNet. I would like to know on what basis these scores have been calculated.
Also I am unable to find out the exact cutoff value associated with the HumanNet score.
Relevant answer
Answer
HumanNet scores were designed to indicate pairs of genes that are more or less likely to share the same GO functional annotation (GO "biological process"). For an example of construction of one of the individual features: for CX, or co-expression, lots of microarray experiments in different conditions, calculating the pearson correlation between all pairs of genes for their expression over different conditions, and using regression to get a score from those correlation values related to GO functional annotation. The individual features were integrated based on their power in predicting GO functional annotation, resulting in an "integrated" score, the last column in the dataset.
There's a better explanation of the derivation methods in the WormNet (2008) paper: http://www.nature.com/ng/journal/v40/n2/abs/ng.2007.70.html
>>Also I am unable to find out the exact cutoff value associated with the HumanNet score.
Not sure what you mean by this. Glad to try to answer if you could clarify.
[Context: I'm not an author of HumanNet, but it's from my lab (http://marcottelab.org), and I've used it a lot.]
Question
20 answers
I am trying a MD simulation on fructose (as a practice) using VMD and NAMD. I am having difficulties in generating the *.PSF file. I am always getting an error saying some atoms can’t be parameterized. Can you please help me to fix this issue?
I am using toppar CHARMM36 all-atom carbohydrate force field files.
Relevant answer
Answer
Hi, I personally think glycam force field is better for carbohydrates. and u can get some help in Amber mailing list
Can anyone suggest any tools or softwares which would aid in visualizing the eigenvectors obtained from PCA of gromacs trajectories.
Question
12 answers
I want to visualize the eigen vectors individually and compare the relative motion. Initially i tried the dynamite server , but it is down from past few months. so, if anyone is using a similar or better tool please let me know
Relevant answer
Answer
You can get trajectory for the given eigenvector through g_analyze using "-extr " option and frame number in this trajectory is control by "-nframes" options. Then, you can use VMD to visualize the trajectory, which will give you the motion along the given eigenvector.
Question
18 answers
I have calculated modulus of elasticity but I am confused to choose which value?
Do I have to choose the maximum value or do I have to get the value via trend line (plot stress vs strain)?
Relevant answer
Answer
I think taking an average would be a better idea.try to get the best curve fitting. It may exclude some points which are out of the path.
Question
2 answers
Could you please help me out how to do MD simulation using Gromacs for a protein which has Threonine residue being phosphorylated (TPO). Ending up with Error "Residue TPO not found in Residue topology File."
Relevant answer
Answer
For proteins containing the phosphorylated residues, you need to edit the forcefield so that the new residues can be recognized by the forcefield or u can download either one of this forcefield (ffG43a1p.tar.gz/gromos43a1p-4.5.1.tgz)
and use it which already has the related info about the phospho-(Ser,Tyr,Thr)
Best regards
bhavani
What will be the best docking tool for a DNA drug or DNA metal docking
Question
40 answers
.
Relevant answer
Answer
If the metal has any kind of activity in the active site, as far as I know, any docking program can calculate that. In contrast, in the case of the metal has only "structural effects" (such as sustain the ligand which participates in the reaction), you can do it doing some tricks as change the name to the metal from X to Zn (in the case of GOLD, that works) or create a new atom type. In one of the works where I collaborate, my colleges use a "huge proton" with 6 valences... it's chemically meaningful but can model the structural effect of the metal.
Getting error in autodock execution for autogrid running.
Question
32 answers
When I am running autogrid run, I am getting error that I cannot open or find Grid parameter file, I am using Windows 7.
Relevant answer
Answer
Thanks @Rajnikant for ur kind response, I have some more doubts, I am facing problem when performing docking of metal with protein, can any body suggest me how to make metal pdb file, when i am running autogrid of complex, i am getting error as unknown ligand atom type x, add parameters for it to the parameter library first. what should i do now, I am using Win 7,I have one more doubt like, do we need cygwin for autodock 4.2, or not, as I have cygwin and while execution i am following the path from directly the scrippt research institute instead of from cygwin folder. I have saved my input files in user in one folder, do I need to keep it in C drive itself the one which I am giving to program pathname for autogrid. exe selection at time of execution.
Question
12 answers
I have found that it is the spatial fourier transform of gs(r,t). But gs(r,t) is calculated at a particular time. So what time should I select? And then for the fourier transform which k value should I use? Since I am not familiar with fourier transformation can you suggest me which algorithm should I use?
Relevant answer
Answer
Either you can do the Fourier transform of this quantity or you can directly calculate intermediate scattering function by the formula attached here. The <> averaging is the ensemble average which is generally being done from the time trajectory given that your system is ergodic. Basically, one can theoretically derive those expression from the spatial fourier transform of pair-pair correlation function, gs(r,t).
If you want to follow the direct path (calculating from gs(r,t)), then you can choose the time interval 0,5, 10, 15 or any regular timestep interval upto a large value - lets say 100000 timesteps. Please remember you have to have the trajectory snapshots which are consistent with these interval.
Thanks.
Question
2 answers
Do you know any tutorial for MD on DNA-Ligand Complex? Can I use Gromacs for DNA ?
Relevant answer
Answer
Hi,
1- See http://ambermd.org/tutorials/basic/tutorial1/. It is for AMBER, but with some changes easely transferable for simulations with GROMACS.
By googling, I have also found this one http://wwwuser.gwdg.de/~ggroenh/SaoCarlos2008/html/build.html#top. They are old, but they provide good starting points.
2- Probably yes with the adequate force field (for instance CHARMM, GROMOS and AMBER (bsc0)) All these force fields are supported by GROMACS.
HTH
Question
5 answers
I want to have comprehensive information on the force fields present in gromacs, in the PDB2GMX command.
Any article on comparative analysis of ff will be welcomed.
Also, info. on which ff to use on a particular system and its explanation would be useful.
Relevant answer
Answer
Hello,
GROMACS contains a lot of force fields for many types of biomolecular simulations (i.e. protein, detergents, lipids, ionic solutions, etc.) so it is difficult to give you a definitive response. Be more precise in your question. For protein force fields you can see the reference attached:
In case you second question, there is no clear/definitive response since it depends a lot of what do you want to do. A good start is to (always) read the literature about previous MD works of similar system.
HTH
Question
35 answers
I want to know what type of analysis we can do from gromacs results, specifically the graphs?
And what is the best time step after which we should observe simulations?
Relevant answer
Answer
Hello,
I would suggest u to go throgh this link, where u can have an idea about the analysis of simulation
Question
7 answers
I am trying to run a quench molecular dynamics on a molecule. This molecule consists of Cyclopentadienyl rings with Fe together with a (CH)2 chain. I am using Materials Studio 6.1 (Accelrys) version program package.
The problem is, when I click on the run bottom then it is showing "the dynamics can't integrate the equation of motion for a particle of zero (or small) mass". But if I cut the Cyclopentadienyl rings together with Fe then there is no problem. That means there is no problem with the rest of the (CH)2 chain.
Does anyone know where the problem is?
Relevant answer
Answer
The problem that you are facing is that you have a cyclopentadienyl ring which has a dummy atom representing the Pi system. This is not recognised by Forcite or Discover as a valid particle for molecular dynamics because it has no mass associated with it and hence the program gives an error.
I think that it is pretty hard to study MD of organometallics using classical simulations techniques. I did a quick search for "organometallics" on
This turned up many hits but publications were nearly all using DMol3. I guess you could potentially delete the Pi system and create bonds between the cyclopentadienyl fragment and the metal atom but I am not convinced this would give good results (picking a forcefield might be challenging!)
Question
5 answers
I have read a couple of tutorials and am trying to incorporate the basic analysis tools available in gromacs like g_rmsd, g_rmsf, g_rms, g_covar, etc. But as I have recently started, I am confused about how to proceed in an appropriate direction. I am interested in studying the allosteric movements in my protein. Any help for me to channelize properly would be really kind.
Relevant answer
Answer
Hi,
Sorry but your question is not very precise. A good start is to read the literature about similar system as yours Probably you have ligand(s) in your system (correct ?), so you could examine the influence of the ligand on the protein dynamic against the apo protein. For instance, by examining
- the stability of the protein vs. time : g_rmsd
- the fluctuation of the protein domains w/o ligand : g_rmsf
- the protein domain motions w/o ligand : g_covar
- Interaction between protein residues (in the pocket ) - ligand : g_hbonds and g_saltbr
- etc...
HTH
Question
10 answers
I want to use Gromacs for energy minimization and molecular dynamics simulation.
I see different simulation times (like; 200ps, 1000ps, 5ns) in papers. Does anyone have suggestions about the best time for simulation, and the factors that affect the decision?
Relevant answer
Answer
Hi, a simulation lenght depends from what information you want to get from your MD.
Generally a MD simulation should be carried out untill the property you are looking for reaches convergence. For example, if you are interested in generating a set of snapshots for MMPBSA analysis, the property you might look at can be the protein RMSD, so you can run your simulation untill the RMSD vs time plot reaches a plateau. If you start from a well defined crystal structure, generally you'll reach convergence in a relatively short time. Otherwise, if you have to test different ligands docked in the same protein structure, simulation time will depend on how the protein need to "adapt" to the different ligands.
Differently, if you need or want to observe large protein motion following a particular "perturbation" (i.e. activation following the binding of an agonists to an inactive state of the protein), simulation time might be decidedly longer.
You might want to give a look to http://dx.doi.org/10.1016/j.jmgm.2013.02.002
Question
9 answers
By different time scale such as 5000ps, 10000ps etc. for structural analysis by different server.
Relevant answer
Answer
Hi,
You don't need a webserver to do that, you can use the g_angle, g_hbonds tools of GROMACS with the appropriate snapshots and index file of your MD trajectory, since you did your MD with GROMACS. Go to the gromacs website and start the helps of these tools.
Additionnally, i don't think you did 5000 ns or 10000 ns of MD ;)
HTH
Question
11 answers
I want to get the psi and phi angles on a particular residue and plot them on ramachandran plot. I went through some tutorials and was informed to use ndx file for that. But can someone give me an actual command to do so.
Relevant answer
Answer
Hi
I think this is what u are looking for
To extract the angles for a single residue, e.g. ILE-32, type:
------grep ILE-32 ramachandran.xvg | awk '{print $1 " " $2'} > rama-ILE-32.xvg
How to retrieve a phosphorylated form of a protein insilico to perform molecular dynamics studies.
Question
64 answers
In my protein, tyrosine kinase (TYR) is phosphorylated. I have the structure of the normal form. Is there a server, tool, or any other method to replace the TYR with phospho-TYR? I want to use this phosphorylated form for my MD studies.
Relevant answer
Answer
Dear Bhavaniprasad, I am guessing GROMACS does not recognize your phosphotyrosine modification because its forcefield parameters are not defined. You can download Phosphotyrosine forcefield parameters for GROMACS here: http://www.gromacs.org/@api/deki/files/12/=ffG43a1p.tar.gz Its a modified GROMOS96 43a1 forcefield that includes phosphoserine, phosphothreonine, and phosphotyrosine parameters Briefly below are the steps you need to take. 1. untar the downloaded file in your current directory where you setup your simulations. 2. You can use a variety of tools like the ones suggested by others to add the phosphotyrosine residue into your PDB file 3. To let GROMACS recognize the residue correctly, the important thing is to make sure the phosphotyrosine residue in your input PDB file has its residue name as PTR and the atoms are named as listed below and strictly in that order. N H CA CB CG CD1 HD1 CD2 HD2 CE1 HE1 CE2 HE2 CZ OH P O1P O2P O3P C O NOTE: If the phosphotyrosine residue you inserted using VMD or any other tool has a different atom naming convention, you can manually edit the PDB file and change them to the ones listed above. NOTE: Explore the downloaded forcefield directory to find PTR residue listing in the ffG43a1p.rtp file. 4. using this modified PDB file as your input run pdb2gmx -ff ffG43a1p <including the other usual command line options) You should now be able to interactively choose the forcefield ffG43a1p to include phosphotyrosine in your simulation box. Hope this helps Venky
Question
4 answers
I will simulate PbTiO3 structure, how should I use the distance of atomic range
for the 3D model?
Relevant answer
Answer
Begin with a fairly sparse k-point grid, or even Gamma, and increase the energy cut-off in steps of 20% from a low value until the total energy of the model system stops decreasing on each step by, say 0.1 meV. Then, increase the grid density monotonically until the total energy stops changing by the same amount as chosen previously for the convergence criterion, and that's it.
Question
2 answers
I am working on RNA simulations and having a problem. When I start simulation from full stretched sequence, it never goes to native structure (a T-loop). It is a very small RNA (20 nt), and it is never showing any base-pairing during simulation. I am using amber99sb.
Can anybody suggest to me something which could be helpful for my work?
Relevant answer
Answer
Hello,
I think the problem is your force field. The AMBER99SB force is not optimized for DNA and RNA simulations . It is mainly (but not only) used for proteins. For simulations of RNA/DNA, I would suggest you to use the recent bsc0χOL3 force field (J Phys Chem B 116(33):9899-916 (2012), PMID 22809319)
For more information, you can also read for instance the following papers
Effect of Guanine to Inosine Substitution on Stability of Canonical DNA and RNA Duplexes: Molecular Dynamics Thermodynamics Integration Study
Miroslav Krepl, Michal Otyepka, Pavel Banáš, and Jiří Šponer
The Journal of Physical Chemistry B
2013 117 (6), 1872-1879 DOI: 10.1021/jp311180u
The Nuclear Magnetic Resonance of CCCC RNA Reveals a Right-Handed Helix, and Revised Parameters for AMBER Force Field Torsions Improve Structural Predictions from Molecular Dynamics
Jason D. Tubbs, David E. Condon, Scott D. Kennedy, Melanie Hauser, Philip C. Bevilacqua, and Douglas H. Turner
Biochemistry
2013 52 (6), 996-1010
DOI: 10.1021/bi3010347
Performance of Molecular Mechanics Force Fields for RNA Simulations: Stability of UUCG and GNRA Hairpins J. Chem. Theory Comput., 2010, 6 (12), pp 3836–3849
DOI: 10.1021/ct100481h
HTH
Question
2 answers
I made a model for my protein using Modeller with template having 46% identity. Now I want to use Gromacs for energy minimization and molecular dynamics simulation.
I see different simulation time (like; 200 ps, 1000ps, 5 ns,...) in papers. Does anyone have suggestions about the best time of simulation, and factors that affect the decision?
Relevant answer
Answer
There is no easy answer for this question. The amount of time required for a simulation depends a lot upon the properties that you want to study, and how long the system will need to provide you with sufficiently high-quality data that your answers will be reliable. A time appropriate for one system may be completely insufficient or overkill for another simulation of a different material or under different state conditions.
There are a few basic rules of thumb that you can follow:
(1) You need to an equilibrate a system until there is no discernible trace of how you "built" the system in the data that you're collecting.
(2) Your "production runs" should be long enough to reduce the uncertainty in your your statistical averaging to an "acceptable" level. (You can specify what that level of uncertainty is. The uncertainty in the averaging can be predicted using classical tools, or can be estimated with methods such as those developed by Flyvbjerg and Petersen.)
Question
11 answers
Looking to buy the best computational and molecular docking studies for my QSAR studies for my molecule.
Relevant answer
Answer
In my opinion Schrodinger is the best option if you can afford it. However, if your budget is limited than ICM (Molsoft) is one of the best options. It costs only around 3000$ and can achieve better docking results than Glide. Other good options are Molegro's MVD for docking and MDM for QSAR. MVD is also CUDA ported. I personally tested it and the CUDA version achieve really good results and 17x speed up!
Good Luck
Question
148 answers
I am new to simulations and gromacs. I have docked the protein-ligand complex and now want to simulate it using Gromacs. I am using the best docked conformation of protein-ligand structure. Should I remove hydrogens from it, as hydrogens are causing problems, and ignore the missing residues and generate the topology of protein ( and generate the ligand topology by PRODRUG server)? Is this the correct method?
Relevant answer
Answer
Hello.
Instead of ProDRG, try using ATB: Automated Topology Builder and Repository (http://compbio.biosci.uq.edu.au/atb/). You can submit your own molecule or browse through their database.
But, if you want to go deep, calculate partial charges by ab-initio (Gaussian).
Best regards.
Question
28 answers
I am new to gromacs. When I run my protein file to generate topology for it, by following the tutorial example mentioned here (http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/complex/02_topology.html).
I am getting the following error for my protein pdb file:
Fatal error:
Atom HA in residue SER 705 was not found in RTP entry SER with 10 atoms
while sorting atoms.
Can anybody please help me in removing this error. What should I do for this SER HA atom?
Relevant answer
Answer
Hello.
That is no missing residue. Simply, in GROMACS, atom names must match. For instance, if the acidic proton in serine is named H, but in the database is named HA, they don't match and give that specific error.
Open the aminoacid database (aminoacids.rtp) of the correspondent FF and check the atom (probably HA). Then, go to your protein and see if the correspondent proton is also HA:
- If not, change it to HA;
- If missing, you have a protein protonation problem. Re-protonate the protein.
Best regards.
Question
8 answers
I have been running a membrane protein simulation. I followed the tutorial on: http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/membrane_protein/index.html
But I faced a problem in the shrinking phase. I have done up to 14 cycles of the shrinking phase using inflate gro. I got only 2.14513126335007 nm^2 area per lipid. The protein was in the membrane up to 13 cycles but after the 14th cycle it came away from the lipid bilayer. I was unable to find the error. Please help me in this regard.
Relevant answer
Answer
Hello.
Since GROMACS 4.5.4, you have an alternative method to insert protein in membranes: g_membed (Wolf, MG et al. J Comput Chem, 31(11), 2169–74. doi:10.1002/jcc.21507), and works a lot better.
Check the GROMACS manual, you'll find instructions. If you need further assistance, I can send to you more detailed instructions and the .mdp file.
Best Regards.
Question
4 answers
Can it say anything about changes in domains if the protein ligand interacts with only one domain?
Relevant answer
Answer
Hello.
It is possible, but will depend on the flexibility of your receptor. With MD, the whole structure will adapt towards the perturbation made by each specific ligand. However, it can be difficult to calculate binding energies. One way is through Gibbs Free Energy: if you simulate ligand and protein alone, and then the complex ligand-protein, you can calculate the difference between energies and extrapolate the binding energy. Another thing that you can search is for evidences of 'induced-fit' mechanism between ligand and receptor.
Best regards.
Question
8 answers
MD file generated need to be analyzed. What could be the best program and where it is available with all the details required?
Relevant answer
Answer
Hello.
First of all, which features do you want to analyse? You need to be more specific on this issue.
Best regards.
Question
13 answers
I have a protein (crystal structure) complexed with two molecules of the same ligand (inhibitor). I am trying to build a pharmacophore model for searching potential hits. But the problem is that two molecules are bound in same binding pocket but in different conformation (image attached). Both the molecules are interacting with hydrophobic interaction. The binding pocket consists of exclusively hydrophobic residues. Could some one please suggest me how can I develop a pharmacophore model by using these? Shall I straightforwardly consider any one of the ligand molecules? I was thinking of merging two similar pharmacophore features to make one point. But the problem is that the binding conformation of both the inhibitor molecules is different.
Also, if I want to do docking of several inhibitors in the same site before going for pharmacophore modelling (for building common pharmacophore model), how should I do this?
Please help me in solving this. Thanking you in anticipation.
Regards.
Relevant answer
Answer
Given your ligand and the pocket size, there may not be a simple answer or a straightforward way use existing pharmacophore generating software. You have to decide if you want to find a ligand that is of similar size, where two of them will bind next to each other. Or do you want to find a single ligand that spans the whole binding site region.
An initial place to start is to try to find similar ligand based on simple molecular descriptors (e.g. via ShaEP, Feature trees, MOE’s MACCS). Carefully selected descriptors might allow you to find a few similar compounds that would bind and stack the way the native ligand does.
However, I would recommend starting even more basic, and form some hypotheses that you can computationally test. A better understanding of the system may help you develop a suitable approach for searching for new ligands. By investigating the system, you will develop some chemical intuition concerning it. Some basic questions I would have for you system are
1) Could the ligands adopt different binding poses under physiological conditions?
2) What are the ligand-protein binding strengths of each ligand alone and when together?
3) What is the ligand-ligand interaction energy?
4) How dynamic is the binding pocket?
5) What are the more dynamically stable ligand-protein interactions?
6) How important are all of those sulfur atoms (if that is what the yellow atoms are), and what changes occur if they are replaced with oxygen atoms?
7) What forces are important for binding? How large is the molecule’s dipole moment?
8) How conformationally labile is the molecule?
Insight into many of these questions can be gain through good MD simulations. Such an investigation would be time consuming, but potentially rewarding.