The results of FEP are inconsistent with the results of the tutorial

[Zinc-X]

Dear Prof.

We finished the FEP calculations of FXR (12 to 74, 76, 84, 85, 88) without REST2, but got different results from you(tutorial01/tutorial/ddg_to_center.csv). We also noticed that the calculated results of 5 ligands are different from the experimental values in the article. Our calculation results are shown as follows：

  | Our ddg | Tutorial | Exp. Dg FXR_12 | 0 | 0 | -9.93271 FXR_74 | 0.070849 | 0.109407 | -8.48765 FXR_76 | 0.185977 | -0.23893 | -6.00989 FXR_84 | 0.009766 | -0.16105 | -7.33357 FXR_85 | 0.079347 | 0.511437 | -8.95884 FXR_88 | 0.272636 | 0.759893 | -8.60273

For FEP calculations, all parameters and procedures we use are based on those given in tutorial [OPLS/AA-m, Ligpargen, gromacs, plumed, without REST2/hrex] Our .ini configuration is shown below

# This is the configuration section for prepare_dual_topology.py (which is the only section in this file, but the [ section ] is mandatory)
[prepare_dual_topology]

# Read ligands and topologies from this folder
input_ligands = lig_data

# Read the macromolecule structure from this file
structure = receptor_data/5q17_processed.pdb

# This is the force field directory, it will be copied to each perturbation dir and used to prepare the MD systems
extradirs = oplsaam.ff

# Options controlling the core-constrained superimposition
# First select the use of it instead of reading all ligand poses
pose_loader = superimpose

#lambda set up
lambda_input = lambdas12

# Use this pose as the reference for the superimposition
poses_reference_pose_superimpose = receptor_data/9mv.pdb

# Name of the output. This will be a self-extracting bash file
perturbations_dir = tutorial-t1

# Sets the path to GROMACS executable in the run node, uncomment and modify if needed.
# gmx_bin_run = /usr/local/bin/gromacs

output_hidden_temp_dir=False

# Options controlling the output, see manual for more info. Uncomment as needed
# FEP legs are to be submitted to a slurm scheduler
#output_scripttype = slurm
# Run these commands at the beginnig of the jog (useful to load modules, importing libs)
# output_runbefore = module load python3; module load cuda
# Fine-tune job resources
#output_resources = all_cpus:48; all_gpus:4; all_time: 24
# Use a python file instead of a binary during the collect step
# output_collecttype = python

So we would like to ask you if you can give some guidance on how to calculate the simulation accurately. And qe would be deeply grateful if you could give more detailed configuration files of this case.

[luancarvalhomartins]

First, you are comparing experimental (absolute) ΔG to calculated (relative) ΔΔG. I redid your table to compare relative to relative free energies (below). Second, due to stochastic nature of MD, each run will sample a different part of the phase space, therefore, yielding a (hopefully) slightly different ΔΔG. So, comparing different runs must be done in this context. Third, comparing to experimental results is even more complex, but your FEP results pretty match ours, which suggests the method to the be consistent.

Lig	Our ΔΔG	Tutorial	Exp. ΔG	Exp. ΔΔG
FXR_12	0.0	0.0	-9.9	0.0
FXR_74	0.1	0.1	-8.5	1.4
FXR_76	0.2	-0.2	-6.0	3.9
FXR_84	0.0	-0.2	-7.3	2.6
FXR_85	0.1	0.5	-9.0	1.0
FXR_88	0.3	0.8	-8.6	1.3