# Umbrella sampling - WHAM analysis **URL:** \ **Category:** User discussions\ **Tags:** analysis-tools, free-energy\ **Created:** [April 13, 2023, 3:18pm UTC](https://gromacs.bioexcel.eu/t/umbrella-sampling-wham-analysis/6179 "2023-04-13T15:18:40Z")\ **Posts on this page:** 1\ **Page:** 1
**Author:** ![yasser.almeida](https://avatars.discourse-cdn.com/v4/letter/y/76d3ee/32.png) [@yasser.almeida](https://gromacs.bioexcel.eu/u/yasser.almeida)\ **Post date:** [April 13, 2023, 3:18pm UTC](https://gromacs.bioexcel.eu/t/umbrella-sampling-wham-analysis/6179/1 "2023-04-13T15:18:40Z")
GROMACS version: 2023 GROMACS modification: No Dear all, I am doing a WHAM analysis of an umbrella sampling of an enzyme-inhibitor complex and I have some questions. I did SMD simulations to pull the inhibitor from the enzyme along the X-axis with `pull_coord1_k = 1000` kJ mol^-1 nm^-2, at different pulling velocities 0.1, 0.05, 0.02, 0.01, 0.005 Å/ps. I got frames along the pulling axis for US sampling during 10 ns. Then, I run the wham command as: `gmx_mpi wham -it tpr-files.dat -ix pullx-files.dat -o profile_pullx -hist -unit kCal -b 1000 -min 2.375790 -max 5 -nBootstrap 100` Below is a figure of the histograms and the PMF for each pulling velocity. ![US_histograms_PMFs](https://europe1.discourse-cdn.com/flex017/uploads/bioexcel1/original/2X/6/61bdc8fd260be0cfe52564cce93738900b1d890b.jpeg) At first glance, it seems that the pulling velocities of 0.1 and 0.005 Å/ps give the best histogram overlapping. The ∆G was calculated manually as the difference between the minimum and the maximum value of the PMF values (∆G = min(PMF)-max(PMF)). 1 - Does the wham command give the ∆G value directly or should it be calculated as before manually? 2 - Although the ∆G values are reasonable, these are 2x the expected value ≈ 14 kcal/mol. I am aware that this analysis can be difficult and tricky. Could you give some tips to improve this calculation? My mdp file for the US simulations is: ```auto title = Umbrella pulling simulation define = -DPOSRES_X ; Restraining the enzyme ; Run parameters integrator = md dt = 0.002 tinit = 0 nsteps = 5000000 ; 10 ns nstcomm = 10 ; Output parameters nstxout-compressed = 5000 ; every 10 ps nstenergy = 5000 ; Bond parameters constraint_algorithm = lincs constraints = all-bonds continuation = yes ; Single-range cutoff scheme cutoff-scheme = Verlet nstlist = 20 ns_type = grid rlist = 1.4 rcoulomb = 1.4 rvdw = 1.4 ; PME electrostatics parameters coulombtype = PME fourierspacing = 0.12 fourier_nx = 0 fourier_ny = 0 fourier_nz = 0 pme_order = 4 ewald_rtol = 1e-5 optimize_fft = yes ; Berendsen temperature coupling is on in two groups Tcoupl = Nose-Hoover tc_grps = Protein Non-Protein tau_t = 1.0 1.0 ref_t = 310 310 ; Pressure coupling is on Pcoupl = Parrinello-Rahman pcoupltype = isotropic tau_p = 1.0 compressibility = 4.5e-5 ref_p = 1.0 refcoord_scaling = com ; Generate velocities is off gen_vel = no ; Periodic boundary conditions are on in all directions pbc = xyz ; Long-range dispersion correction DispCorr = EnerPres ; Pull code pull = yes pull_ncoords = 1 ; only one reaction coordinate pull_ngroups = 2 ; two groups defining one reaction coordinate pull_group1_name = emzyme pull_group2_name = inhibitor pull_coord1_type = umbrella ; harmonic potential pull_coord1_geometry = distance ; simple distance increase pull_coord1_dim = Y N N pull_coord1_groups = 1 2 pull_coord1_start = yes ; define initial COM distance > 0 pull_coord1_rate = 0.0 ; restrain in place pull_coord1_k = 1000 ; kJ mol^-1 nm^-2 ``` I am considering repeating the calculations with a different pulling constant to improve the spacing of the initial configurations for the US. Thanks in advance.