Dear GMX users, I am trying to calculate the PMF of a drug inside a membrane protein which forms a channel. The reaction coordinate is the z-axis. In order to do so, I would like to use Umbrella sampling of the drug in many positions along the z-axis, followed by WHAM. That should give me the PMF for each position along the pore, inside-out, right?
I am having trouble defining the parameters for the simulations. As I understand, I should generate many starting structures, in which the position of the drug differ by, say, 0.5 Angstrom apart on the z axis. So let's say I have 30 starting structures (protein + drug + lipid bilayer + water). Then Umbrella sampling should be performed on each of those structures, so I should also have 30 files of pull.ppa. The difference between them is the position to which the pulled group (drug) is restrained to. Does that sound reasonable so far? This is how my pull.ppa looks like: -- verbose = no runtype = umbrella ngroups = 1 group_1 = AMAN ;reference_group = weights_1 = reference_weights = reftype = com reflag = 0 pulldim = Y Y N ; DYNAMIC REFERENCE GROUP OPTIONS r = 0 rc = 0 update = 1 ; UMBRELLA SAMPLING OPTIONS K1 = 1000 Pos1 = 33.482 31.225 33.173 -- K1 is the force constant acting on the pulled group in order to restrain it, right? how do I determine the best value for it? For Pos1 - is that OK to input the COM (center of mass) of the ligand, as it is in each starting structure? Since this is the position I would like to restrain it to. I did not input a reference group, since I'm not sure if the reference should be the protein or if it'll be better to work in absolute coordinates. Pulldim = Y Y N, since for each starting structure, the point is to see if it moves along x,y due to the forces acting on it by the protein. As I understand, since the ligand's z position is different in each starting structure, we can analyze it with WHAM and calculate the PMF. Does that sound right, or am I completely wrong? The result of a simulation with the aforementioned parameters is that the ligand escapes the protein *into* the lipid environment, which doesn't make sense. When I use K1=10, the drug does not escape the pore. So I'm confused here. I thought K1=1000 was supposed to put more force on the drug to stay in its original position. Any help and ideas regarding why this is wrong? Thanks in advance, Hadas. _______________________________________________ gmx-users mailing list gmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to [EMAIL PROTECTED] Can't post? Read http://www.gromacs.org/mailing_lists/users.php