Freezing previously-mutated residues in a staged/sequential multi-mutation pmx workflow

pmx
1

Hi all,

I’m running a sequential (staged) multi-mutation alchemical calculation with pmx + GROMACS: I mutate one residue per stage (e.g. T2V), equilibrate both end states, then use the state-B equilibrated structure as the starting point for the next stage’s mutation, and so on through 16 stages.

Problem: starting from stage 2, GROMACS/pmx appears to treat all previously mutated residues (still named T2V in the structure/topology) as active dual-topology (A/B) mutations again, not just the current stage’s target. I confirmed this by checking [ atoms ] in each stage’s topol.top — every prior T2V residue still carries distinct typeB/chargeB entries, not just the new one. As supporting evidence, the peak dH/dλ magnitude for forward window 1 scales roughly linearly with stage number (~2,000 kJ/mol·λ⁻¹ at stage 1 up to ~22,000 at stage 16), consistent with an increasing number of residues being simultaneously coupled to λ rather than just one per stage.

I believe this happens because gmx pdb2gmx -ff charmm36m-mut (and pmx gentop) attach B-state parameters to any residue named T2V, with no way to distinguish “this stage’s new target” from “an already-completed prior mutation.”

Question: What’s the correct/recommended way to “freeze” a previously mutated residue as a plain, single-topology amino acid (e.g. VAL) before starting the next stage’s mutation, so it’s no longer treated as dual-topology? I want to avoid simply renaming T2V→VAL in the coordinate file, since the state-B structure still contains the extra (now non-interacting) dummy atoms unique to the original residue (e.g. THR’s OG1/HG1), which wouldn’t match VAL’s .rtp entry.

Is there a pmx utility/workflow for this, or is the expected approach to manually strip the dummy atoms and regenerate the topology for that residue via the plain (non-mut) force field before proceeding to the next stage? Section 4.3 of Aldeghi et al. (Methods Mol Biol 2019) discusses reusing equilibrium trajectories across a mutation scan, but I haven’t found guidance for this specific chained/sequential case.

Thanks,
Sankalp Jadhav

2

you’re right. Every time you run gmx pdb2gmx -ff charmm36m-mut (and pmx gentop) a new hybrid structure + topology is added, and so subsequent stages in your workflow will define more and more morphs.

I am guessing that you are doing the mutations sequentially rather than all at once to give the system more time to equilibrate after each individual mutation and therefore get a more converged result? I can see the reasoning behind that, but please also keep in mind that the overall error propagates over all 16 individual steps, which might become substantial as well.

So if that were the rationale, I would consider as alternative to do all mutations at once. This is a massive perturbation, so you may need long non-equilibrium transitions (of perhaps 10 ns or more) to get sufficient overlap between forward and backward runs, but if it converges you do not need to worry about reverting mutated residues back to normal residues in the structures and topologies, or about error propagation across the 16 stages.

Bert