Issue in Probability distribution overlap for D2E transformation

pmx
1

Dear All,

I am performing amino acid transformation for my protein-peptide interaction studies.
When I do D2E mutation in both peptide and complex arms with 51 lambdas (5 ns each), there is a poor Probability distribution overlap in the terminals. So we increased lambda points to 55 (5 ns each) for terminal convergence. (Attached both Old and New Overlap Image along with BAR data and MDP files), Even after increasing, still no proper terminal convergence.

Old Parameters used:

; Alchemical Transformations
free-energy = yes ; use free energy code
init_lambda_state = 0 ; start from state A
;lambda_states = 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50
fep_lambdas = 0.0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2 0.22 0.24 0.26 0.28 0.3 0.32 0.34 0.36 0.38 0.4 0.42 0.44 0.46 0.48 0.5 0.52 0.54 0.56 0.58 0.6 0.62 0.64 0.66 0.68 0.7 0.72 0.74 0.76 0.78 0.8 0.82 0.84 0.86 0.88 0.9 0.92 0.94 0.96 0.98 1.0
calc_lambda_neighbors = 1 ; only immediate neighboring windows
couple-intramol = yes
;
; soft-core function
sc-function = gapsys ;
sc-gapsys-scale-linpoint-lj = 0.85 ;
sc-gapsys-scale-linpoint-q = 0.3 ;
sc-gapsys-sigma-lj = 0.3 ;
;
nstdhdl = 10 ; write to dhdl at each step
dhdl-derivatives = yes
dhdl-print-energy = total
separate-dhdl-file = yes

Kindly do the needful. Thank You

New_md.mdp (3.5 KB)

New_Overlap
New_Overlap1695×890 207 KB

New_Res.txt (23.5 KB)
Old_Overlap
Old_Overlap1815×890 215 KB

Old_Res.txt (21.8 KB)

2

Could it be a plotting issue? Why is there only one distribution in the first and last lambda window? From the text file, I don’t see a particularly large uncertainty for the terminal window dG estimates. Also, in one of the plots, lambda values go up to 1.08

3

Thank you for your response, and my sincere apologies for the delayed reply.

I have attached the entire BAR data for both new and old approaches, for reference.
(https://drive.google.com/file/d/1BqCvgyxzVoDXf--SdQlzk6Bc5RV2XBct/view?usp=sharing)

Narration: In my complex, the aspartic acid of the peptide interacts with the arginine of the Protein via a salt bridge. Once we mutate aspartic acid to glutamic acid (increasing the spacer length by one carbon), the interaction is lost. I believe that the lack of overlap is due to this irreversible loss of salt-bridge interaction.

Coming to your response,

  1. No, Sir, the issue is not with the Python script (attached). Here, I plotted the distribution for each of fifty-one lambda.xvg, based on the following formula:
    dH(i-1)λ - dH(i)λ : orange
    dH(i)λ - dH(i+1)λ : green
  2. For the terminal xvgs, I have either the next or previous dH λ data. That is the reason for the single peak.
  3. Since I used the same Python script for both 51 and 54 lambdas, the lambdas go to 1.08, kindly read as 1.0 (Final lambda)

My request here is, if there are any ways to have a smoother transition for D2E.

Another doubt:
What is the significance of the “rotations” in mutres.mtp? Since both wild-type and mutant amino acids have proper descriptions in the ffbonded.itp file.

Thank you for your time and consideration.

4

I don’t see lack of overlap in your figure where lambda goes properly from 0 to 1.0. Which lambda window do you find problematic there?

5

Sir, in the case of fifty-five lambdas, the lambdas 2, 3, 4, 51, 52, 53 and in the case of fifty-one, 2 and 49. Are they proper, as it is?

Thank you.

6

In the 51 lambda case, it seems that the overlap is fine. Strange that for more lambdas the overlap decreases. Do you use HREX? It should help with convergence

7

I will definitely consider HREX in my future, Sir.

Any word about “rotations”, Sir? The reason is, how is the end state taken care of when the rotations are focused only on the initial state? For example, in the F2L mutation

[ rotations ]
CA-CB HB1 HB2 CG CD1 HD1 CD2 HD2 CE1 HE1 CE2 HE2 CZ HZ DCG HV1 DCD1 HV2 HV3 HV4 DCD2 HV5 HV6 HV7
CB-CG CD1 HD1 CD2 HD2 CE1 HE1 CE2 HE2 CZ HZ

I expect,

CB-DCG …
DCG-DCD1 …
DCG-DCD2 …

Though I assume “rotations” help to ensure smooth alchemical transitions and avoid structural artifacts. I am unable to understand the relations. I have performed with and without rotations, and the results are similar; with both end states are fine for a tripeptide (in water). Kindly enlighten me.

Thank You.

8

The rotations in .mtp are there to position dummy amino acid in the initial configuration. They do not affect the alchemical transitions.

9

Thank you so much for your response and time, Sir.