GROMACS version:
GROMACS modification: Yes/No
Here post your question: I am trying to calculate free energy profiles of AquaporinZ and its mutants through umbrella sampling. But the pmf curve I got is far more different than the curves I found in previous literature. On left is the reference image and on right (with peaks and energy curve) is the one I got after umbrella sampling and using wham analysis.
Even the values I got on x-axis are different than the actual pore dimensions of Aquaporins ( which are between -20 to +20 angstrom). Is there any way to solve this issue? Thanking in advance for the responses.
100 per window? Also the x you get is just your reaction coordinate, what it is specifically depends on what you are biasing. If your CV is along the axis of the aquaporin then -20A/+20A is compatible with your reported reaction coordinate length of 3.x nm ca. But this is something that depends on how you set up the system, that is, how you generated the starting windows, if you do not populate correctly all the values of the reaction coordinate that you want to sample then you won’t explore all the space you are interested into.
This is a brief look of the whole system, cyan is the Aquaporin protein, embedded in a POPC lipid bilayer, submerged in a water box. I do not use a single water molecule and do not perform steered MD simulation. I directly use step6.6_equilibration.gro file (as shown in figure) as an input for umbrella sampling and use multiple water molecules (located inside the protein channel) and mention their indices in index.ndx file. Reaction coordinates are defined one by one on the basis of these water molecules for US simulations.
Each NPT step run for 100ps=0.1ns & US mdrun for 5 ns with 600kJ mol^-1 nm^-2 , by using following commands.
The resulted .tpr files and pullf files were used to generate histograms by using following command and the resulted graphs were attached in the query above.
So you run 5ns or sampling per window. Is this the same as the paper you are referring to? It seems very very short. Secondly, what pull groups are you specifying in the input file? The reaction coordinate depends on that. You are keeping the centre of mass of the molecule of water restrained with respect to…? The centre of mass of the lipid bilayer? Of the protein?
In some reaction coordinates, water molecule move in a correct direction(which is z-axis), but in some windows, water move on x-axis and thus sometimes move out of the pbc box.
How do I keep it at the pore centre or make it move in only on z-axis?
I do not know much about the sampling windows mentioned in previous literatures. But if 5ns are not enough then I would like to try for 10ns.
I specify the water molecule as pull group1. And a residue of protein as pull group2, located at the opening of channel passage.
The drift is due to the poor choice of spring constant. if i remeber, i used in between 500-2000 (unit) spring constnt for umbrella sampling with aquaporin.
I don’t fully understand your set up, but for sure one thing you should do is check previous literature and see how much sampling time do they need. 5ns seems short, but I may be mistaken. You have the reference, as you posted the plots.
You want to sample the PMF along the channel, right? Basically, along the trans-bilayer axis z, that you also report in your previous image. Then your reaction coordinate ξ is z, which is probably centered at the center of mass of the trans-bilayer protein channel. The first thing that seems out of order is then the value of your reported histograms - why from 0 to +3? Where are you setting the zero of the reaction coordinate? Are you sampling just half channel? If you want to sample along z, then the restraint must be applied along that axis, which means that you should see the water molecule oscillate around a certain z value (supposing you are dividing the systems in windows, sampling them individually, and then combining stuff together with something like WHAM, which seems to be the case looking at your posts). The molecule can (and, in principle, should) move freely along the xy plane, as you are not biasing those directions. What is your mdp set up?
Issue is not with the window or sampling time- water molecule is passing through outside the channel cavity. @EishaKhilji If you mail me the files, i can have a look.
Very well. What if my script isn’t in GROMACS but with OpenMM instead? Am I still allowed to make a post?
This is the only forum I was able to find that addresses the issue of unsymmetrical PMF plots. Best
