Custom Search
-- Jason M. Swails > On Sep 25, 2024, at 2:59 PM, Adrian Roitberg via AMBER <amber.ambermd.org> wrote: > > Hi, > > > This is what I sent Nick. It is more complicated than just the charge smearing or the intrinsic radii. > > > ------------------- > > Hi Nikolay. > > My group wrote the CpH code for Amber, modifying the original from Case and others. > > Cysteine is a BIG problem, and we do not trust the results. The problem seems technical, so let me explain. > > In the Amber implementaion, when we change protonation state, we just change the charges in the atom, nothing else. We keep everything else in the force field the SAME for protonated and unprotonated forms, and we "hide" all other differences inside the so-called reference energy, which is fitted to get the correct pKa for the reference compound, which is the simplest peptide . > > As you might have seen, it works great for almost everything, except Cysteine. > > The issue is simple. Cysteinate (S-), is MUCH larger than Cysteine (SH). (3.92 A vs 3.56 A) > > See https://doi.org/10.1063/1.5038010 for how they fixed this in CHARMM. > > Amber uses the SAME vdw radius for S- and S in SH, which makes waters to be too close to the S in deprotonated cysteines. Makes sense ? > > So, the solution is to have a different set of radii for S- and S in SH, and to change vdw radii AND charges once one changes protonation states (or attempts to). > > On top of that, on the regular CpH version in Amber, the actual montecarlo step to decide if the protonation change is accepted is done in GB, not in explicit solvent, which would require also changing the intrinsic GB radii for S- and S in SH. > > We more or less know what needs to be done to at least test this, and we know it is important, but at this point, my group does not have the people to work on this problem. We will try find time at some point, but honestly, it does not look good. > > If you are interested, we can discuss things, and you can take over, but it will require going knee deep into the amber code. > > Let me know if you want to chat. > > Adrian > > ------------------------------- > >> On 9/25/24 2:30 PM, Jason Swails via AMBER wrote: >> [External Email] >> >> Cysteines are not widely titrated to my knowledge, and it wasn't until I >> added support for them that the cpin generator would create input files for >> them. >> >> That said, there are three tunable "sets" of parameters for cysteine in >> constant pH simulations: >> >> 1. The force field parameters for the atoms. These are the same as the >> non-titratable parameters, as is true for every other titratable residue in >> Amber. >> 2. The partial charge sets for the protonated and deprotonated states. >> There is a prescribed recipe for assigning these used by every other >> titratable residue, and it is described in Mongan's original paper from >> 2004. In short, the backbone charges are the same as the "protonated" >> state, the sidechain charges are the same as the standard (non-titratable) >> version of that residue for each given protonation state, and the >> "difference" is dumped on the side-chain carbon atom in the beta position >> to give the correct formal charge for the entire residue. >> 3. The reference energy. >> >> However, the reference energy is fully determined by the first two (and >> *must* be equal to the value required to yield the correct titration curve >> of the reference compound). There are additional constraints on how much >> flexibility you have in modifying either of the first two parameter sets. >> First, the code only changes charges and only evaluates energy terms >> involving the charges when computing the energy difference used in the >> Monte Carlo evaluation, meaning you can't assign different bonded >> parameters (like torsion parameters) to different charge states of a single >> residue. Second, you *can't* change the charge of any atom that is within 2 >> bonds (3 atoms away) from one of the atoms that bonds with an adjacent >> residue (the carboxyl carbon or amide nitrogen of the backbone). The reason >> for this is also explained in Mongan's paper, and comes down to the fact >> that the scaled 1-4 nonbonded interactions are effectively part of the 1-4 >> bonded interaction, and if these are free to change between charge states >> then this bonded term starts to depend on the actual residue sequence. As >> there are three different sets of backbone charges in the default Amber >> amino acid parameters (one set for all +1 amino acids, another for all >> neutral amino acids, and a third for all -1 amino acids), the reference >> energies for titratable residues whose backbone atoms change charges would >> become sequence dependent (and the complexity absolutely *explodes* when >> you think about neighboring titratable residues). >> >> I tried this during my PhD, and it failed spectacularly for all the reasons >> I described. The only other strategy I can even think of that would be >> reasonable would be to "smear" the excess charge over the entire >> side-chain, rather than dumping it all on the carbon-beta (and, of course, >> you would need to recalculate the reference energy). You could also try >> other GB models, but you'd need to compute the free energy for these >> choices, too. In my experience, the impact of these choices on the computed >> pKa is minimal (less than a pK unit), which is an absolutely miniscule >> deviation in free energy. It will certainly not resolve the failure of a >> residue to titrate at all at a certain pH. >> >> The last knob you can turn would be the intrinsic solvation radius, which >> was a problem I identified with the original implementation (the "dummy" >> hydrogens led to a systematic overestimation of the carboxylate oxygen >> atoms on aspartate and glutamate residues). Correcting this systematic >> overestimation, (reducing the intrinsic radius by 0.1 Angstroms), shifted >> the pKa of buried carboxylates by about 2 pK units (it makes sense >> conceptually). However, the correction for these carboxylate residues was >> "obvious" (the distribution of effective radii for carboxylate oxygens was >> almost exactly shifted by 0.1 Angstroms between the AS4 and ASP residues, >> so the adjustment yielded largely the correct distributions for the >> deprotonated forms). >> >> This is a long-winded way to say that "parameter improvement" is a large >> and complex project that I've found to require a very intimate knowledge of >> the models involved (GB, CpHMD, charge fitting, TI/FEP, and probably >> others), and a good grasp on the underlying implementation. >> >> An important piece of knowledge is "when does CYS titrate?" -- keep >> increasing the pH until you start to see transitions. Maybe you need to go >> up to pH=20. Maybe higher. Eventually you will see CYS -> CYM transitions. >> You really should try to narrow down the *actual* pKa (which you can only >> get when you have a pH or two where you have protonation fractions that are >> between 0 and 1). Without knowing that, you won't even know the scale of >> the problem you're facing, and you'll basically be stabbing blindly at >> solutions. >> >> Good luck, >> Jason >> >> On Wed, Sep 25, 2024 at 8:57 AM Nikolay Kuzmich < >> nikolay.kuzmich.weizmann.ac.il> wrote: >> >>> Hi Jason, thanks a lot! >>> >>> I have started a 80ns-CpH run at pH=9.5. >>> >>> And yes, the CYS residues are in the cpin file. >>> >>> But since Adrian said that the cysteines are sort of underdeveloped in the >>> Amber CpH MD implementation, >>> >>> I'll do that rather from curiosity and I'll try to improve the >>> CYM parameters for CpH MD. >>> >>> Kind regards >>> >>> Nick >>> ------------------------------ >>> *From:* Jason Swails <jason.swails.gmail.com> >>> *Sent:* Tuesday, September 24, 2024 7:52:14 PM >>> *To:* Nikolay Kuzmich; AMBER Mailing List >>> *Subject:* Re: [AMBER] CpH MD in explicit solvent: no switch for CYS >>> although it was in the titratable residues list >>> >>> *Caution: External Sender. Do not click on links or open attachments >>> unless you recognize the sender.* >>> >>> In explicit solvent, protonation state changes are attempted for every >>> residue (in random order) every ntcnstph steps, which means that >>> protonation state changes were attempted for every CYS residue. If the >>> protonation state never changed, then every attempt was rejected. This >>> means that the energy penalty was far too large to be overcome by the pH >>> adjustment (and that in this simulation, the CYS residues are highly >>> unlikely to be found in the deprotonated form). >>> >>> A few things to check: >>> >>> * Make sure that CYS residues are actually specified in the cpin file >>> * Run some short sample simulations at much higher pH values (keep >>> increasing solvph until you start to see CYS transitions). >>> >>> You can also generally get a good idea of why protonation changes are >>> constantly rejected by looking at the local environment of the CYS residues >>> in the trajectory. >>> >>> Good luck, >>> Jason >>> >>> On Mon, Sep 23, 2024 at 4:13 AM Nikolay Kuzmich via AMBER < >>> amber.ambermd.org> wrote: >>> >>>> Dear Amber users and developers, >>>> >>>> just wanted to ask you about one detail: >>>> >>>> I carried out constant pH MD in explicit solvent, ASP, GLU, HIS and CYS >>>> were included in the list of titratable residues (the former 3 were renamed >>>> in AS4, GL4 and HIP, accordingly) and pH was set to 7.35. >>>> >>>> >>>> cpinutil.py -igb 2 -resnames HIP AS4 GL4 CYS -p >>>> cpH_prot_RC_tip3p18oct_hmr.prmtop -op cpH_prot_RC_tip3p18oct_hmr_cpi.prmtop >>>> -o cpH_prot_RC_tip3p18oct_hmr.cpin >>>> >>>> >>>> the mdin file was: >>>> >>>> &cntrl >>>> imin = 0, nstlim = 125000000, dt = 0.004, >>>> irest = 1, ntx = 5, ig = -1, >>>> temp0 = 300.0, icnstph=2, >>>> ntcnstph=100, ntrelax=200, solvph=7.35, >>>> ntc = 2, ntf = 2, >>>> ntwx = 12500, ntwr = 12500, ntpr = 2500, >>>> cut = 9.0, iwrap = 1, saltcon=0.15, >>>> ntt =3, gamma_ln=2.0, ntb = 1, ntp = 0, >>>> nscm = 1000, >>>> ioutfm=1, ntxo=2, >>>> / >>>> I applied the HMR technique as well. >>>> >>>> >>>> Having performed analysis after 500ns, I found that the cysteines >>>> underwent no transition and remained all the time in SH state; then I >>>> carried out the simulation 80ns long at pH=8.5 and again no transitions of >>>> CYS into anion were observed. From literature it is known that one CYS can >>>> have unusually low side-chain pKa. >>>> >>>> I guess some additional work-around for cysteines should have been done... >>>> >>>> Can you please prompt me which one exactly? >>>> >>>> Kind regards >>>> >>>> Nick >>>> >>>> Computational chemistry and >>>> molecular modeling >>>> >>>> The Nancy and Stephen Grand Israel National Center for Personalized >>>> Medicine >>>> Weizmann Institute of Science >>>> >>>> M +972 52 7439671 >>>> 234 Herzl Street, P.o. Box 26, >>>> Rehovot, 7610001 Israel >>>> >>>> >>>> [1694421514730] >>>> _______________________________________________ >>>> AMBER mailing list >>>> AMBER.ambermd.org >>>> http://lists.ambermd.org/mailman/listinfo/amber >>>> >>> >>> -- >>> Jason M. Swails >>> >> >> -- >> Jason M. Swails >> _______________________________________________ >> AMBER mailing list >> AMBER.ambermd.org >> http://lists.ambermd.org/mailman/listinfo/amber > > -- > Dr. Adrian E. Roitberg > Frank E. Harris Professor > Department of Chemistry > University of Florida > roitberg.ufl.edu > 352-392-6972 > > > _______________________________________________ > AMBER mailing list > AMBER.ambermd.org > http://lists.ambermd.org/mailman/listinfo/amber _______________________________________________ AMBER mailing list AMBER.ambermd.org http://lists.ambermd.org/mailman/listinfo/amberReceived on Wed Sep 25 2024 - 15:00:01 PDT