[gmx-users] Re: Wierd results from Umbrella sampling

[Thomas Schlesier]

The picture of pullx.xvg is not really helpful. The reason is that in 
pullx.xvg the second column (after the time) is the position of the 
reference group (which naturally moves; but it's not the thing in which 
wie are interested) and the third column is the distance from the 
reference group to the pulled group (this is the one which interestes you).
Another thing you could try is to pull in all 3 dimensions, if you pull 
only in the z-direction your pulled group can move freely in the x-y 
plane and wouldn't be affected by the umbrella potential. But i think 
for the histogram only the direction in which one pulls are accounted 
for, so the wide histograms are somewhat strange.
For the problems with the constraints i have absolute no idea, never had 
this error message. Even don't know if the problems you get from there 
follow you into the pull-mode simulation. Does the constraint simulation 
run without the pull-code, or do you get there the same problem?
greetings
thomas



> I increased the potential to 2350 instead of 1000 in my md_umbrella.mdp (the other parameters kept unchanged), then the movement of my protein was less than the previous one, which we can see from pullx.xvg (attached). But it still not reach to the designed region. So I am going to increase the potential higher again.
>
> It still does not work if I switched on the LINCS. The error is that the number of constraint is more than the number of freedom.
>
> I am running another trail with high potential (let's say 3000, ok?) and dt=2 ft, which value is definitely not recommended by Martini developers.
>
> With Regards,
>
> Jiangfeng.
>
>
>
> Message: 3
> Date: Tue, 22 May 2012 18:53:54 +0200
> From: Thomas Schlesier<schlesi at uni-mainz.de>
> Subject: [gmx-users] Re: Wierd results from Umbrella sampling,  (Justin
>          A. Lemkul)
> To:<gmx-users at gromacs.org>
> Message-ID:<4FBBC4A2.9020506 at uni-mainz.de>
> Content-Type: text/plain; charset="ISO-8859-1"; format=flowed
>
> I never worked with the MARTINI (or other coarse-grained) force field,
> but this in the umbrella.mdp
>
> title       = Umbrella pulling simulation
> integrator  = md
> dt          = 0.019
>
> looks suspicious. The dt is about an order of magnitude greater than one
> uses in normal (bond-)constrainted md-simulations. I know that in
> coarse-grained simulations the potentials are smoother than in atomistic
> force fields and one therefore could use a higher timestep, but i don't
> know if you can go so high.
> Anyway you should look for the umbrella potential. If the force constant
> and the timestep are both too high you could get problems:
> Assume a particle in a harmonic well. If the force constant is high and
> the timestep too high, you wont sample the harmonic well, but 'jump'
> each time from one side to the other (high force + great timestep ->
> long movement).
> So i would test it first with a timestep of 0.002 (same as you used for
> pulling sim).
>
>
>
>
>
> ------------------------------
>
> Message: 4
> Date: Tue, 22 May 2012 19:15:41 +0200
> From: "Justin A. Lemkul"<jalemkul at vt.edu>
> Subject: Re: [gmx-users] Re: Wierd results from Umbrella sampling,
>          (Justin A. Lemkul)
> To: Discussion list for GROMACS users<gmx-users at gromacs.org>
> Message-ID:<4FBBC9BD.8070702 at vt.edu>
> Content-Type: text/plain; charset=ISO-8859-1; format=flowed
>
>
>
> On 5/22/12 6:53 PM, Thomas Schlesier wrote:
>> I never worked with the MARTINI (or other coarse-grained) force field, but this
>> in the umbrella.mdp
>>
>> title = Umbrella pulling simulation
>> integrator = md
>> dt = 0.019
>>
>> looks suspicious. The dt is about an order of magnitude greater than one uses in
>> normal (bond-)constrainted md-simulations. I know that in coarse-grained
>> simulations the potentials are smoother than in atomistic force fields and one
>> therefore could use a higher timestep, but i don't know if you can go so high.
>> Anyway you should look for the umbrella potential. If the force constant and the
>> timestep are both too high you could get problems:
>> Assume a particle in a harmonic well. If the force constant is high and the
>> timestep too high, you wont sample the harmonic well, but 'jump' each time from
>> one side to the other (high force + great timestep ->  long movement).
>> So i would test it first with a timestep of 0.002 (same as you used for pulling
>> sim).
>>
>
> Testing this would be good, though the MARTINI developers routinely report
> timesteps of 20-40 fs as "normal" use.  I have never obtained stable
> trajectories above 20 fs, but I also do not do much coarse graining.  It could
> very well be that the pull code is not stable with such an integration step.
>
> -Justin
>
>>
>> Am 22.05.2012 18:37, schrieb gmx-users-request at gromacs.org:
>>> Thank you for your reply, Thomas. Under your explanation, I am well understood
>>> of the terms: pull_k and pull_rate.
>>> Here I upload both md_umbrella.mdp and md_pulling.mdp.
>>> I have to mention that it is coarse gained system with Martini force field.
>>>
>>> At the same time, i am going to run a simulation without pull code but with
>>> LINCS constraint, and also run another system with a huge pull_k but without
>>> LINCS. Hope I could get some interesting information.
>>>
>>> With Regards,
>>>
>>> Jiangfeng.
>>>
>>> Message: 3
>>> Date: Tue, 22 May 2012 16:36:47 +0200
>>> From: Thomas Schlesier<schlesi at uni-mainz.de>
>>> Subject: [gmx-users] Re: Wierd results from Umbrella, sampling (Justin
>>> A. Lemkul)
>>> To:<gmx-users at gromacs.org>
>>> Message-ID:<4FBBA47F.5010503 at uni-mainz.de>
>>> Content-Type: text/plain; charset="ISO-8859-1"; format=flowed
>>>
>>> Think it would be best to show the .mdp file, else we can only guess
>>> what the parameters are.
>>>  From the histogram it looks like that the force constant of the
>>> restraining potential is too low, since the histograms are really wide,
>>> but pull_k1=1000 is a 'normal' value.
>>> On thing which confueses me is you said that fluctuations from g_dist
>>> are about 0.4nm, but in the histogram i looks like the distance
>>> fluctuates about 1nm or even more. for this be sure, that you compare
>>> the right distances ->  if you pull only in z-direction, the only look
>>> into the z-direction from the output from g_dist. Since the x, and
>>> y-directions would be unaffected by the umbrella potential.
>>>
>>> for the error message with the constraints:
>>> what happens if you run the system with constraints but without the
>>> pull-code?
>>>
>>> for pull_k1>0 and pull_rate1=0:
>>> if you're pulling with an umbrella potential pull_k1 defines the force
>>> constant of the potential (hook'sches law). Let's assume you put a
>>> spring to your molecule
>>> with pull_rate1=0, it means that the other end of the spring doesn't
>>> move, and the spring restrains the position of the molecule (molecule
>>> cn't diffuse away). in umbrella sampling you want to restrian of
>>> molecule (or part of it) relative to another one / part ->  so pull_rate1=0
>>> with pull_rate1>0, it's like you pull the spring away for the molecule,
>>> spring gets strechted ->  wants to relax ->  molecule follows spring (like
>>> that what happens in afm pulling)
>>>
>>>
>>> Am 22.05.2012 15:40, schriebgmx-users-request at gromacs.org:
>>>>> Hi Justin, As for the maximum of 0.4 nm fluctuation of my pulled
>>>>> protein, I used g_dist to calculate the distance between my pulled
>>>>> protein and stable part in a window, where the distance is treated as
>>>>> the fluctuation of the protein along z-axis. Well, from pullx.xvg, the
>>>>> position moved a lot (3nm for instance.) As for the windows simulation,
>>>>> I didn't apply constraint but only the internal constraint in the itp
>>>>> file. I still don't understand why it have to do constraint? why not
>>>>> give a fully freedom to run the simulation? In the md_umbrella.mdp, I
>>>>> set pull_k1=1000; pull_rate1= 0.0, but apparently, I am confused with
>>>>> pull_k1>0 combined with pull_rate1=0. In my mdp, i set "none" to LINCS,
>>>>> because if I use "all-bonds", an error of "1099 constraints but degrees
>>>>> of freedom is only1074" occurs. Actually, there is no any window with a
>>>>> designed distance. Here I attach the histo.xvg, profile.xvg. Thank you
>>>>> with regards, Jiangfeng. On 5/22/12 9:36 AM, Du Jiangfeng (BIOCH) wrote:
>>>>>>>>> Dear Justin,
>>>>>>>>>
>>>>>>>>> Based on your questions to my simulation, I posted here yesterday
>>>>>> hopefully
>>>>>>>>> it was the correct way to reply in this forum.
>>>>>>>>>
>>>>> You've still replied to the entire digest message (which I've cut out); please
>>>>> make sure to keep replies free of superfluous posts in the future. The archive
>>>>> is already pretty hopeless, but let's not make it worse:)
>>>>>
>>>>>>>>> In this morning I got a list of new windows of umbrella sampling, the
>>>>>> overlap
>>>>>>>>> is sufficient enough, but I saw another problem: In the histogram figure,
>>>>>>>>> the base of peak covers the distance of 2 nm instead of 0.2 nm., that's
>>>>>>>>> horrible! However, when I checked back to the simulation results of each
>>>>>>>>> window, the fluctuation of my pulled protein is only 0.4nm in maximum.
>>>>>> So the
>>>>>>>>> base of peak shouldn't cover such long distance, right?
>>>>>>>>>
>>>>> If the peaks aren't corresponding to the desired restraint distances, then
>>>> there
>>>>> are several potential problems:
>>>>>
>>>>> 1. Your restraints aren't set up the way you think they are (check grompp
>>>> output
>>>>> and .tpr file contents to be sure)
>>>>> 2. Your restraints are ineffectual (in which case you may need to revisit the
>>>>> force constant)
>>>>>
>>>>> I can't determine from your description what's going on. What do you mean by a
>>>>> maximum of 0.4 nm fluctuation? In what quantity? What do the contents of
>>>>> pullx.xvg show you for the problematic window, and for that matter, the
>>>> others?
>>>>> Are any of them producing the desired restraint distances?
>>>>>
>>>>> Again I will ask you to share an image of the histogram and PMF profile; these
>>>>> would be very helpful to see.
>>>>>
>>>>> -Justin
>>>>>
>>>>> --
>>>>> ========================================