Hi Fanyu,
maybe running on CPU is sufficient then, pmemd.cuda is really more
optimized for PBC condensed phase systems, small molecules in vacuum may
not even see much speedup there.
Id say your equation is missing the "periodic box to single molecule in
vacuum" correction factor, but even with that in place, my guess is that it
would not work very well, the MMPBSA solvation free energy model is quite
different from what TI will give you and adding them up like that adds
considerable inaccuracy to your result.
If you want the ddG of your peptide mutation in solvent, you are way more
likely to get a useful result by explicitly simulating that change in
solvent (or figuring out why that failed before), than by calculating a
vacuum transformation and correcting that somehow.
Kind Regards,
Thomas
On Mon, Jun 5, 2023 at 10:35 AM Fanyu Zhao <fz2113.nyu.edu> wrote:
> Hi Thomas,
>
> Thanks for the suggestion!
>
> I just tried running igb=0 with pmemd on CPU and it worked fine. So maybe
> it's something wrong with CUDA?
>
> I think you are right that I don't have a fully closed thermodynamic
> cycle. I came up with an equation, like ddG(wt to mutant, in solvent, from
> TI)=ddG(wt to mutant, gas phase, from TI)+ddG(difference of mutant
> solvation energy and wild type solvation energy, from MM/GBSA). But it
> seems like not accessible.
>
> I'll try the alchemical ion approach again.
>
> Thanks,
> Fanyu
>
> On Mon, Jun 5, 2023 at 4:02 PM Steinbrecher, Thomas <
> thomas.steinbrecher.roche.com> wrote:
>
>> Hi Fanyu,
>>
>> sorry to hear igb=0 didnt work, but have you checked running the same job
>> using pmemd or sander on CPU? If it segfaults there as well, it likely is
>> something about the input you are using. If you want to make sure the
>> charge change is the source of the problem, could you try the same
>> calculation with a protonated Asp, just to check your input (plus one H)
>> works in that case?
>>
>> I still dont understand what you need the vacuum simulation for, if this
>> is about peptide binding, then you would normally simulate the peptide
>> transformation in solvent and bound to whatever it is binding to. A vacuum
>> TI transformation is only needed if you care about the change in solvation
>> free energy (from ddG(TI,solv) - ddG(TI,vac) ) and since that is probably
>> not something experimentally accessible for a peptide, what would you use
>> the result for?
>>
>> If you really care about this transformation in vacuum, it may be hard to
>> get a meaningful answer. Changing a charge in vacuum is not the same as
>> changing it in a periodic boundary environment, due to the charge-charge
>> interactions between periodic images, so your thermodynamic cycle doesnt
>> quite close. I'm not 100% sure what the current best practice is to handle
>> such cases, but the alchemical ion approach is intended to get rid of such
>> box size effects. Since the alchemical ion approach does not work in
>> vacuum, Id stear clear of total system charge changes for such a case.
>>
>> Kind Regards,
>>
>> Thomas
>>
>> On Mon, Jun 5, 2023 at 9:42 AM Fanyu Zhao <fz2113.nyu.edu> wrote:
>>
>>> Hi Thomas,
>>>
>>> Thanks for your reply!
>>>
>>> Yes, I tried using igb=0, and the simulation just stopped at reading
>>> restart file. The error message gives a segmentation fault so I think there
>>> might be something wrong with the input file I'm using.
>>>
>>> Actually, I'm trying to calculate the relative binding free energy of a
>>> peptide mutation involving charge-changing (e.g. mutating Asp to Ala, with
>>> the charge changing from -1 to 0). So far I've tried using co-alchemical
>>> ions to balance the charge between the two systems, or running the
>>> simulation without any correction of the changed charge. Both failed.
>>>
>>> So I'm thinking if I could firstly calculate the energy change of
>>> mutation on gas phase, and secondly calculate the solvation energy
>>> difference of the wild-type and mutant, and then sum the two terms up.
>>>
>>> I cannot figure out any other methods in TI to calculate the mutation
>>> involving charge-changing :(. It would be appreciated if you could give me
>>> any suggestions on this.
>>>
>>> Best,
>>> Fanyu
>>>
>>> On Mon, Jun 5, 2023 at 3:29 PM Steinbrecher, Thomas <
>>> thomas.steinbrecher.roche.com> wrote:
>>>
>>>> Hi Fanyu,
>>>>
>>>> have you tried running with igb = 0 ? It should be the same effect as
>>>> igb=6 based on what I see in the manual, but maybe there is a check in
>>>> pmemd about "igb > 0" that is missing that igb = 6 would also correspond to
>>>> vacuum?
>>>>
>>>> As a side note, based on the masks, you seem to have a 10+ residue
>>>> peptide in your simulation, are you sure a vacuum simulation of this
>>>> species is meaningful?
>>>>
>>>> Kind Regards,
>>>>
>>>> Thomas
>>>>
>>>> On Mon, Jun 5, 2023 at 4:26 AM Fanyu Zhao via AMBER <amber.ambermd.org>
>>>> wrote:
>>>>
>>>>> Dear amber users,
>>>>>
>>>>> I'm trying to run TI in the gas phase. But it always gives me an error
>>>>> message of "The CUDA code currently cannot handle TI with implicit
>>>>> solvent
>>>>> models! Please use the serial code instead".
>>>>>
>>>>> I was wondering if gas phase TI can only be performed with CPU and if
>>>>> I'm
>>>>> correct with the process of creating a pdb file in gas phase.
>>>>>
>>>>> Here's my input for a minimization.
>>>>>
>>>>> Minimization with position restraint
>>>>> &cntrl
>>>>> imin = 1, ntmin = 2, maxcyc = 10000,
>>>>> ntpr = 1000, ntwe = 1000,
>>>>> ntr = 1, restraint_wt = 100.0,
>>>>> restraintmask = '!.H= ',
>>>>> ntb=0,igb=6, cut = 9999.0, ioutfm = 1,
>>>>> icfe = 1, ifsc = 1, logdvdl = 0,
>>>>> timask1=':10',
>>>>> timask2=':3',
>>>>> scmask1=':10 & !.CA,C,N,O',
>>>>> scmask2=':3 & !.CA,C,N,O',
>>>>> clambda=0.5,
>>>>> &end
>>>>> /
>>>>>
>>>>>
>>>>> And here's the tleap file I used to generate a gas phase pdb file.
>>>>>
>>>>> source leaprc.protein.ff14SB
>>>>>
>>>>> com=loadpdb 1tvb_aftermd.pdb
>>>>> unbound=loadpdb peptide.pdb
>>>>>
>>>>> savepdb com wild_md_bound.pdb
>>>>> saveamberparm com wild_md_bound.prmtop wild_md_bound.inpcrd
>>>>>
>>>>> savepdb unbound wild_md_unbound.pdb
>>>>> saveamberparm unbound wild_md_unbound.prmtop wild_md_unbound.inpcrd
>>>>> quit
>>>>>
>>>>>
>>>>> Thank you so much!
>>>>>
>>>>> Best,
>>>>> Fanyu
>>>>> _______________________________________________
>>>>> AMBER mailing list
>>>>> AMBER.ambermd.org
>>>>> http://lists.ambermd.org/mailman/listinfo/amber
>>>>> <https://urldefense.proofpoint.com/v2/url?u=http-3A__lists.ambermd.org_mailman_listinfo_amber&d=DwMFaQ&c=slrrB7dE8n7gBJbeO0g-IQ&r=TA6Yv1GqGslhLaifgggoUQ&m=4WjHVK2ypm7vH0Okye8p2Ny0JAsKRcIz9NXV9bVSYGgs24GTPqJvmd0_PYIvgSFY&s=jj3P4Wq5mxqWq6hYdbE91qLTi7Stm6DKikZPbRtmspw&e=>
>>>>>
>>>>
>>>>
>>>> --
>>>> *Dr. Thomas Steinbrecher*
>>>> Principal Scientist CADD
>>>>
>>>> Roche Pharma Research and Early Development
>>>> Roche Innovation Center Basel
>>>> F. Hoffmann-La Roche Ltd
>>>> Bldg. 092/3.92
>>>> Grenzacherstrasse 124
>>>> 4070 Basel
>>>> Switzerland
>>>>
>>>> Phone +41 61 682 1319
>>>> mailto: thomas.steinbrecher.roche.com
>>>>
>>>
>>
>> --
>> *Dr. Thomas Steinbrecher*
>> Principal Scientist CADD
>>
>> Roche Pharma Research and Early Development
>> Roche Innovation Center Basel
>> F. Hoffmann-La Roche Ltd
>> Bldg. 092/3.92
>> Grenzacherstrasse 124
>> 4070 Basel
>> Switzerland
>>
>> Phone +41 61 682 1319
>> mailto: thomas.steinbrecher.roche.com
>>
>
--
*Dr. Thomas Steinbrecher*
Principal Scientist CADD
Roche Pharma Research and Early Development
Roche Innovation Center Basel
F. Hoffmann-La Roche Ltd
Bldg. 092/3.92
Grenzacherstrasse 124
4070 Basel
Switzerland
Phone +41 61 682 1319
mailto: thomas.steinbrecher.roche.com
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Mon Jun 05 2023 - 02:30:02 PDT