How to set-up three-step transformation path in FEP calculation in GROMACS mdp file

User discussions
1

GROMACS version: 2022.5
GROMACS modification: No

Hi, Everyone!

I am new to FEP calculation and recently perform some FEP calculations of protein single mutation using GROMACS + pmx.

In my benchmarks, the calculation details are listed below:

  1. Mutation
    ala-scaning mutation(residue mutated to Ala) and normal mutation(residue mutated to other non-Ala residues)

  2. Transformation Path
    “one-step” and “two-step” transformation path for. And “one-step” transformation path means ele and vdw change simultaneously with 12 or 21 lambda states, while “two-step” transformation path means ele term changes firstly, then vdw and 21 lambda states are used.

  3. Alchemical Ion is added and fixed at box corner to keep charge neutralization for charge-change mutation

In my results, my potocol had a normal performance in ala-scaning mutation with “two-step” transformation path. But in many normal mutations cases FEP calculation with “two-step” tranformation path, the simulations in some lambda windows are easy to fail due to atom collapse, although Beutler soft-core potential was used(I am not sure whether Gapsys soft-core will alleviate this issue). The reasone is in some windows, vdw is too weak, while ele is stong to attract water or other atoms. Using “one-step” transformation path, FEP calculation can be done, but the performance is worse than ala-scaning mutation. Analysis of overlap between lambda windows with pymbar shows small overlap(< 0.1, even 0) in “one-step” transformation path.

So My Question is:
I think whether the performation and the success FEP calculation can be improved with “three-step” transformation, that is turn off wild type ele firstly, then vdw change from wild type to mutant, and turn on mutant’s ele in the last. This maybe need to set coul_lambdas, vdw_lambdas in mdp independently for wild type and mutant residues, like this
coul_lambda_A = 1 0.8 0.6 0.4 0.2 0 0 0 0 0 0 0 0 0 0 0
coul_lambda_B = 0 0 0 0 0 0 0 0 0 0 0 0.2 0.4 0.6 0.8 1
vdw_lambda_A = 1 1 1 1 1 1 0.8 0.6 0.4 0.2 0 0 0 0 0 0
vdw_lambda_B = 0 0 0 0 0 0 0.2 0.4 0.6 0.8 1 1 1 1 1 1

Is it possible or Is there some ways to realize “three-step” transformation path in GROMACS?

Or Do you have suggestions/advices/experience for me to for normal mutation FEP calculation setup?

Looking forward to your reply and Really Appreciated!

FYI, my FEP setup in mdp is listed below:

;#### fep setting

free_energy = yes
init_lambda_state = 0
delta_lambda = 0
couple-lambda0 = vdw-q
couple-lambda1 = vdw-q
couple-intramol = yes

; fep analysis related
calc_lambda_neighbors = 1
dhdl-print-energy = potential
dhdl-derivatives = yes
nstdhdl = 2000

; lambda_index 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

coul_lambdas = 0.00 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00
vdw_lambdas = 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 1.00
bonded_lambdas = 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50 0.55 0.60 0.65 0.70 0.75 0.80 0.85 0.90 0.95 1.00
mass_lambdas = 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50 0.55 0.60 0.65 0.70 0.75 0.80 0.85 0.90 0.95 1.00
restraint_lambdas = 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 0.00
temperature_lambdas = 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 0.00

; soft-potential
sc-alpha = 0.5
sc-coul = no
sc-power = 1
sc-sigma = 0.3