From: Brian Radak (brian.radak_at_gmail.com)
Date: Thu Feb 14 2019 - 12:11:47 CST
Hi Hristina,
There are a number of ways to set up the calculation you want. In my
opinion, the main difficulty is keeping the two alchemical regions in a
similar space near the lambda endpoints. One way to achieve this is to add
extraBonds between between atoms in opposing alchemical regions - this
essentially mimics what "single topology" methods do.
Also, because you are looking at ions, if you are comparing to old CHARMM
calculations, make sure you are not comparing absolute solvation free
energies with PBCs against SSBPs - these will not match.
HTH,
Brian
On Wed, Feb 13, 2019, 7:59 PM Hristina Zhekova <hzhekova_at_ucalgary.ca wrote:
> Hello everyone,
>
>
> I have the following issue. A protein I'm interested in prefers to bind
> some anions over other anions. I want to evaluate the preference of the
> binding site for the different anions. The anions may have different
> charges and numbers of atoms, e.g. Cl- to HCO3- or HCO3- to CO3^2-. I tried
> to follow the NAMD tutorial for alchemical perturbations (the part for
> mutation of Tyr to Ala) with dual topology for these two alchemical
> transformations in a water box (no protein, only ion in about 3000 water
> molecules) to see if I'm getting sensible numbers. I have free energies of
> hydration for all 3 anions from previous calculations that were done with
> CHARMM, so I was comparing the DeltaDeltaG(hydration) values with the
> DeltaDeltaG(alchemical perturbation) values. Presumably, they should be the
> same or very similar (from thermodynamic cycle considerations). However,
> I'm getting significant difference (~10 kcal/mol), which means that I'm
> doing something wrong. I've used the same force fields parameters in all
> calculations (all parameters are available in the current CHARMM toppar
> files). The input files were modified from the tutorial files or Tyr to Ala
> alchemical perturbation. I haven't changed anything in the FEP block. These
> are the dual topologies I used:
>
>
> RESI ZERO -1.00
> GROUP
> ATOM CA CG2O6 0.69 ! H3 O1
> ATOM H3A HGP1 0.43 ! | /
> ATOM O1A OG2D2 -0.76 ! O3--C (-)
> ATOM O2A OG2D2 -0.76 ! Cl \
> ATOM O3A OG311 -0.60 ! O2
> GROUP
> ATOM CLB CLA -1.00
>
> BOND CA O1A CA O2A CA O3A O3A H3A
> IMPR CA O1A O2A O3A
> ! seed = O1A CA O2A
> IC O2A O1 *CA O3A 0.0000 0.00 180.00 0.00 0.0000
> IC O1A CA O3A H3A 0.0000 0.00 0.00 0.00 0.0000
>
>
>
> RESI ZERO -2.00
> GROUP
> ATOM CA CG2O6 0.69 ! O1
> ATOM H3A HGP1 0.43 ! H3 /
> ATOM O1A OG2D2 -0.76 ! O3--C (-)
> ATOM O2A OG2D2 -0.76 ! \
> ATOM O3A OG311 -0.60 ! O2
> GROUP
> ATOM C1B CG2O6 1.42
> ATOM O1B OG2D2 -1.14
> ATOM O2B OG2D2 -1.14
> ATOM O3B OG2D2 -1.14
>
> BOND CA O1A CA O2A CA O3A O3A H3A
> BOND C1B O1B C1B O2B C1B O3B
> IMPR CA O1A O2A O3A
> ! seed = O1A CA O2A
> ! required for out-of-plane vibrations
> IMPR C1B O1B O2B O3B
> IC O2A O1A *CA O3A 0.0000 0.00 180.00 0.00 0.0000
> IC O1A CA O3A H3A 0.0000 0.00 0.00 0.00 0.0000
> IC O2B O1B *C1B O3B 0.0000 0.00 180.00 0.00 0.0000
> IC O3B O2B *C1B O1B 0.0000 0.00 180.00 0.00 0.0000
> !redundant definition needed to enable seeding.
>
>
> The FEP part in the input (forward-off.namd In-aqua) looks like this:
>
>
> source ../tools/fep.tcl
>
> alch on
> alchType FEP
> alchFile solvate.fep
> alchCol B
> alchOutFile forward-off.fepout
> alchOutFreq 10
>
> alchVdwLambdaEnd 1.0
> alchElecLambdaStart 0.5
> alchVdWShiftCoeff 4.0
> alchDecouple off
>
> alchEquilSteps 100
> set numSteps 500
>
> runFEP 0.0 1.0 0.05 $numSteps
>
>
> Increasing the numbers of steps for equilibration and forward/backward FEP
> calculation didn't change the results.
>
>
> While I was checking the NAMD mailing list for answers, I also noticed
> that people use residue patches (PRES) for alchemical perturbations with
> change in charge (protonation/deprotonation and pKa calculations). I didn't
> try this, but it makes me wonder if that would be another possibility for
> such calculations, where we have charged species.
>
>
> I will appreciate any comments and suggestions. Please let me know if you
> need more information.
>
>
> Thank you!
>
>
> Hristina
>
This archive was generated by hypermail 2.1.6 : Thu Dec 31 2020 - 23:17:10 CST