Hi Justin, As you advised I reduced number of my windows and I obtined histogram:
http://speedy.sh/2b9dT/histo.JPG Which looks really good. The corresponding profile: http://speedy.sh/y8Ssz/profile.JPG I do not understand it. Does my deltaG of binding correspond to everything above 0 kcal/mol which is 6 kcal/mol? Thank you in advance, Steven On Thu, Mar 8, 2012 at 2:42 PM, Steven Neumann <s.neuman...@gmail.com>wrote: > > > On Thu, Mar 8, 2012 at 2:18 PM, Justin A. Lemkul <jalem...@vt.edu> wrote: > >> >> >> Steven Neumann wrote: >> >>> Dear Gmx Users, Dear Justin, >>> I pulled my ligand away from my protein. Ligand was attached to lower >>> part of my protein, I pulled in Z coordinate it using: >>> >>> ; Run parameters >>> >>> integrator = md ; leap-frog integrator >>> >>> nsteps = 5000000 ; 2 * 5000000 = 10 ns >>> >>> dt = 0.002 ; 2 fs >>> >>> tinit = 0 >>> >>> nstcomm = 10 >>> >>> ; Output control >>> >>> nstxout = 50000 ; save coordinates every 100 ps >>> >>> nstvout = 50000 ; save velocities every >>> >>> nstfout = 5000 >>> >>> nstxtcout = 5000 ; every 10 ps >>> >>> nstenergy = 5000 >>> >>> ; Bond parameters >>> >>> continuation = yes ; first dynamics run >>> >>> constraint_algorithm = lincs ; holonomic constraints >>> >>> constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained >>> >>> ; Neighborsearching >>> >>> ns_type = grid ; search neighboring grid cells >>> >>> nstlist = 5 ; 10 fs >>> >>> rlist = 0.9 ; short-range neighborlist cutoff (in nm) >>> >>> rcoulomb = 0.9 ; short-range electrostatic cutoff (in nm) >>> >>> rvdw = 0.9 ; short-range van der Waals cutoff (in nm) >>> >>> ewald_rtol = 1e-5 ; relative strenght of the Ewald-shifted potential >>> rcoulomb >>> >>> ; Electrostatics >>> >>> coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics >>> >>> pme_order = 4 ; cubic interpolation >>> >>> fourierspacing = 0.12 ; grid spacing for FFT >>> >>> fourier_nx = 0 >>> >>> fourier_ny = 0 >>> >>> fourier_nz = 0 >>> >>> optimize_fft = yes >>> >>> ; Temperature coupling is on >>> >>> tcoupl = V-rescale ; modified Berendsen thermostat >>> >>> tc_grps = Protein_LIG Water_and_ions ; two coupling groups - more >>> accurate >>> >>> tau_t = 0.1 0.1 ; time constant, in ps >>> >>> ref_t = 298 298 ; reference temperature, one for each group, in K >>> >>> ; Pressure coupling is on >>> >>> pcoupl = Parrinello-Rahman ; pressure coupling is on for NPT >>> >>> pcoupltype = isotropic ; uniform scaling of box vectors >>> >>> tau_p = 1.0 ; time constant, in ps >>> >>> ref_p = 1.0 ; reference pressure, in bar >>> >>> compressibility = 4.5e-5 ; isothermal compressibility of water, bar^-1 >>> >>> ; Periodic boundary conditions >>> >>> pbc = xyz ; 3-D PBC >>> >>> ; Dispersion correction >>> >>> DispCorr = EnerPres ; account for cut-off vdW scheme >>> >>> ; Velocity generation >>> >>> gen_vel = no ; assign velocities from Maxwell distribution >>> >>> ; These options remove COM motion of the system >>> >>> ; Pull code >>> >>> pull = umbrella >>> >>> pull_geometry = distance >>> >>> pull_dim = N N Y >>> >>> pull_start = yes >>> >>> pull_ngroups = 1 >>> >>> pull_group0 = Protein >>> >>> pull_group1 = LIG182 >>> >>> pull_init1 = 0 >>> >>> pull_rate1 = 0.0 >>> >>> pull_k1 = 200 ; kJ mol^-1 nm^-2 >>> >>> pull_nstxout = 1000 ; every 2 ps >>> >>> pull_nstfout = 1000 ; every 2 ps >>> >>> Following Justin's tutorial I used perl script to extract coordinate for >>> each window. >>> >>> 0 2.4595039 >>> >>> 1 2.4745028 >>> >>> ... >>> >>> 500 8.74 >>> >>> My ligand at the begining was at such distance as it was in the lower >>> part of the protein. Then I used 0.1 nm spacing at the begining (till 4 nm) >>> and 0.2 nm later on. >>> >>> And following equilibration in each window I run umbrella sampling for >>> 10ns in app 49 windows: >>> >>> Run parameters >>> >>> integrator = md ; leap-frog integrator >>> >>> nsteps = 5000000 ; 2 * 5000000 = 10 ns >>> >>> dt = 0.002 ; 2 fs >>> >>> tinit = 0 >>> >>> nstcomm = 10 >>> >>> ; Output control >>> >>> nstxout = 50000 ; save coordinates every 100 ps >>> >>> nstvout = 50000 ; save velocities every >>> >>> nstfout = 5000 >>> >>> nstxtcout = 5000 ; every 10 ps >>> >>> nstenergy = 5000 >>> >>> ; Bond parameters >>> >>> continuation = yes ; first dynamics run >>> >>> constraint_algorithm = lincs ; holonomic constraints >>> >>> constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained >>> >>> ; Neighborsearching >>> >>> ns_type = grid ; search neighboring grid cells >>> >>> nstlist = 5 ; 10 fs >>> >>> rlist = 0.9 ; short-range neighborlist cutoff (in nm) >>> >>> rcoulomb = 0.9 ; short-range electrostatic cutoff (in nm) >>> >>> rvdw = 0.9 ; short-range van der Waals cutoff (in nm) >>> >>> ewald_rtol = 1e-5 ; relative strenght of the Ewald-shifted potential >>> rcoulomb >>> >>> ; Electrostatics >>> >>> coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics >>> >>> pme_order = 4 ; cubic interpolation >>> >>> fourierspacing = 0.12 ; grid spacing for FFT >>> >>> fourier_nx = 0 >>> >>> fourier_ny = 0 >>> >>> fourier_nz = 0 >>> >>> optimize_fft = yes >>> >>> ; Temperature coupling is on >>> >>> tcoupl = V-rescale ; modified Berendsen thermostat >>> >>> tc_grps = Protein_LIG Water_and_ions ; two coupling groups - more >>> accurate >>> >>> tau_t = 0.1 0.1 ; time constant, in ps >>> >>> ref_t = 298 298 ; reference temperature, one for each group, in K >>> >>> ; Pressure coupling is on >>> >>> pcoupl = Parrinello-Rahman ; pressure coupling is on for NPT >>> >>> pcoupltype = isotropic ; uniform scaling of box vectors >>> >>> tau_p = 1.0 ; time constant, in ps >>> >>> ref_p = 1.0 ; reference pressure, in bar >>> >>> compressibility = 4.5e-5 ; isothermal compressibility of water, bar^-1 >>> >>> ; Periodic boundary conditions >>> >>> pbc = xyz ; 3-D PBC >>> >>> ; Dispersion correction >>> >>> DispCorr = EnerPres ; account for cut-off vdW scheme >>> >>> ; Velocity generation >>> >>> gen_vel = no ; assign velocities from Maxwell distribution >>> >>> ; These options remove COM motion of the system >>> >>> ; Pull code >>> >>> pull = umbrella >>> >>> pull_geometry = distance >>> >>> pull_dim = N N Y >>> >>> pull_start = yes >>> >>> pull_ngroups = 1 >>> >>> pull_group0 = Protein >>> >>> pull_group1 = LIG182 >>> >>> pull_init1 = 0 >>> >>> pull_rate1 = 0.0 >>> >>> pull_k1 = 200 ; kJ mol^-1 nm^-2 >>> >>> pull_nstxout = 1000 ; every 2 ps >>> >>> pull_nstfout = 1000 ; every 2 ps >>> >>> >>> My PMF profile: >>> >>> http://speedy.sh/zerqZ/**profile.JPG<http://speedy.sh/zerqZ/profile.JPG> >>> >>> My histogram: >>> http://speedy.sh/PyhnN/Histo.**JPG<http://speedy.sh/PyhnN/Histo.JPG> >>> >>> Why g_wham takes into account distances below 2.45 nm as the 1st >>> structure was at 2.45. If I get rid of the distances below 2.45 (those >>> weird values PMF values) I obtain beautiful profile: >>> >>> http://speedy.sh/TUXGC/**profile1.JPG<http://speedy.sh/TUXGC/profile1.JPG> >>> >>> Please, explain! >>> >>> >> The way you're thinking about distance is not consistent. Again, this is >> a hazard of trying to map my tutorial onto your problem. You say you have >> a ligand bound to the "lower part" of your protein, and then you're pulling >> in the z-direction. The COM distance (as measured by g_dist and extracted >> using my script) is not equivalent to the distance along the reaction >> coordinate, if that reaction coordinate is only one dimension. In the >> tutorial, it was. Here, it is not, hence the massive sampling defects that >> you're observing and considerable redundancy in many of your windows. >> >> Check the output of grompp for the actual restraint distances that mdrun >> will interpret. They are printed to the screen. >> >> -Justin >> > > Thank you Justin!!! > > Steven > > > > >> >> -- >> ==============================**========== >> >> Justin A. Lemkul >> Ph.D. Candidate >> ICTAS Doctoral Scholar >> MILES-IGERT Trainee >> Department of Biochemistry >> Virginia Tech >> Blacksburg, VA >> jalemkul[at]vt.edu | (540) 231-9080 >> http://www.bevanlab.biochem.**vt.edu/Pages/Personal/justin<http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin> >> >> ==============================**========== >> -- >> gmx-users mailing list gmx-users@gromacs.org >> http://lists.gromacs.org/**mailman/listinfo/gmx-users<http://lists.gromacs.org/mailman/listinfo/gmx-users> >> Please search the archive at http://www.gromacs.org/** >> Support/Mailing_Lists/Search<http://www.gromacs.org/Support/Mailing_Lists/Search>before >> posting! >> Please don't post (un)subscribe requests to the list. Use the www >> interface or send it to gmx-users-requ...@gromacs.org. >> Can't post? Read >> http://www.gromacs.org/**Support/Mailing_Lists<http://www.gromacs.org/Support/Mailing_Lists> >> > >
-- gmx-users mailing list gmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/Support/Mailing_Lists