choderalab / yank

An open, extensible Python framework for GPU-accelerated alchemical free energy calculations.
http://getyank.org
MIT License
177 stars 70 forks source link

Positive Enthalpy of binding #775

Open mmagithub opened 7 years ago

mmagithub commented 7 years ago

Hi, I am trying to carry out an absolute FE calculation with YANK for a small molecule binding to a receptor. I tried to increase the number of iteration and cycles/iteration to get well converged results. Here is the final output: ####################### Free energy of binding: -26.699 +- 0.831 kT ( -15.917 +- 0.496 kcal/mol) DeltaG complex : 162.771 +- 0.649 kT DeltaG restraint : 0.508 kT DeltaG solvent : 136.580 +- 0.519 kT Enthalpy of binding: 420.726 +- 424.983 kT ( 250.821 +- 253.359 kcal/mol) ########################

What looks strange for me is that the Enthalpy of binding is positive and has a very high error bar. Is it normal?, I am not sure. here is my YAML script ####################

options: minimize: yes verbose: yes output_dir: . number_of_iterations: 1000 nsteps_per_iteration: 1000 temperature: 300kelvin pressure: 1atmosphere timestep: 1.5*femtoseconds output_dir: ligand5_similar2ima resume_setup: yes resume_simulation: yes

systems: bms-135: phase1_path: [input/prod_npt-1.5ns_com.inpcrd,input/5-49_com.prmtop] phase2_path: [input/prod_npt-1.5ns_norec.inpcrd,input/5-49_norec.prmtop] ligand_dsl: resname MOL solvent: PME-solvent

solvents: PME-solvent: nonbonded_method: PME nonbonded_cutoff: 10angstrom clearance: 12angstroms

protocols: absolute-binding: complex: alchemical_path: lambda_restraints: [0.00, 0.05, 0.10, 0.15, 0.20, 0.25, 0.50, 0.75, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.000, 1.00, 1.000, 1.00, 1.00] lambda_electrostatics: [1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 0.95, 0.90, 0.85, 0.80, 0.75, 0.70, 0.65, 0.60, 0.55, 0.50, 0.45, 0.40, 0.35, 0.30, 0.25, 0.20, 0.15, 0.10, 0.05, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.000, 0.00, 0.000, 0.00, 0.00] lambda_sterics: [1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 0.95, 0.90, 0.85, 0.80, 0.75, 0.70, 0.65, 0.60, 0.55, 0.50, 0.45, 0.40, 0.35, 0.30, 0.25, 0.20, 0.15, 0.10, 0.075, 0.05, 0.025, 0.01, 0.00] solvent: alchemical_path: lambda_electrostatics: [1.00, 0.95, 0.90, 0.85, 0.80, 0.75, 0.70, 0.65, 0.60, 0.55, 0.50, 0.45, 0.40, 0.35, 0.30, 0.25, 0.20, 0.15, 0.10, 0.05, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00] lambda_sterics: [1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 0.95, 0.90, 0.85, 0.80, 0.75, 0.70, 0.65, 0.60, 0.55, 0.50, 0.45, 0.40, 0.35, 0.30, 0.25, 0.20, 0.15, 0.10, 0.00]

experiments: system: bms-135 protocol: absolute-binding restraint: type: Harmonic

#################### The ligand is neutral and the protein required 2 chlorine to neutralize. Before running the simulation, the system was prepared and relaxed in AMBER. any clue?

Thanks, Marawan

jchodera commented 7 years ago

This sounds like it might be an overlap issue. Can you generate the Jupyter analysis notebook and post this (or render it to a PDF)? This should help us give useful suggestions for how to improve overlap/convergence.