[gmx-users] Simulation parameter problem about protein unfolding

[chris.neale at utoronto.ca]

Generally, forcefields are not parameterized for temperatures other  
than 298K, so simulations are not expected to reproduce the expected  
properties (like boiling water and the correct temperature  
denaturation of proteins).

There's almost certainly other issues here (including the fact that  
I'm entirely sure that you can get a lot more than 24 ns of simulation  
on a 54 aa protein; and 26 atom of pressure seems pretty arbitrary)  
but it will come down to this eventually.

Just because you found a paper in which they get a denatured state  
does not imply that they got the correct denatured state.

Chris.

-- original message --

Hi All,

I met a problem when I try to unfold a protein using Gromacs, It  
seemed the protein cannot be totally unfolded!

The simulated system has one Engrailed Homeodomain (En) protein (a  
three helix bundle protein with 54 residues, 629 atoms), total 4848  
spce waters, and 7 Cl- used to neutralize the system in a 5.752(nm)^3  
water BOX. I use the NTP ensemble with T=498K and P=26atm. The system  
has 1nm thick water in each side of the En protein, and the density of  
the system has been adjusted to0.829 g/cm3.

The simulation lasted 24ns. The helixes disappeared at about 4ns. And  
after that some beta sheet formed in the N terminal of the protein.  
However, the protein kept in a compact state during the 24ns  
simulation. The radius of gyration of the protein fluctuated around  
1.1nm during the simulation.

I've also noticed similar simulation done by others. For example,  
David Becka and Valerie Daggett reviewed their En protein unfolding in  
paper "Methods for molecular dynamics simulations of protein  
folding/unfolding in solution. Methods 34 (2004) 112¨C120". The  
simulations were performed with ENCAD and ilmm, and used an 0.8nm  
cutoff range. And the ensemble they used is NVE as I know. A stretched  
unfolding state occurred at about 5ns in their 60ns simulation in  
498K, with little helix structure.

I was wondering whether the difference is caused by using different MD  
software and force field, or by some wrong parameters in my .mdp file.  
I've also conducted another 18ns simulation, and the result is almost  
the same. I list he mdp file below. Any comment is appreciated!