Hello,
looking at the update list on
http://ambermd.org/bugfixesat.html I don't
think this patch has ever been released for AT 13. Is that true? Is the
behavior fixed in AT 14?
Cheers,
Jan-Philip
On 03.12.2013 21:21, wmsmith.uci.edu wrote:
> The patch you mentioned was constructed, but its changes were not yet
> pushed to the branch as we are awaiting approval. I can easily amend the
> patch so that pb_read will no longer force radiopt = 2 when phiout=1 and
> will instead display a warning if it is set to some other value.
> Be aware that if you use radiopt other than 2, it will be up to you to
> find and use the appropriate atomic radii in your visualization software
> since the molecular surface will no longer correspond to the one defined
> by the radii in your topology file... If you use the levelset description
> (current default I believe), you can always print out the level set as a
> dx file. If not, you will have to construct the surface yourself.
> -Wes
>> Jan-Philip,
>>
>> Wes has submitted the patch. I think it's pending approval …
>>
>> Also as for radiopt upon the printout of phi map, we can certainly do it
>> as
>> you suggested … The original argument for hardwiring the radii for
>> visualization is this. It is not straightforward to set radii in VMD/PyMol
>> and the phi map is supposed to use values on the outside of the molecular
>> surface, i.e. solvated phi values since the purpose is to get some
>> intuition on the effect of a molecular probe approaching the surface from
>> the outside. Thus we want the best possible match between the radii used
>> in
>> the visualization software and PBSA. Furthermore. PBSA reduces the radii
>> by
>> 0.1 A to make sure the projected phi values on the molecular surfaces are
>> solvated values.
>>
>> However, this is not necessary when the newer DX volumetric data output
>> from Amber is used for surface rendering. The plan is that we'll still set
>> radiopt == 2 as the default condition and remind users to reset radii in
>> the visualization software to those used by Amber if radiopt != 2.
>>
>> Ray
>>
>> --
>> Ray Luo, Ph.D.
>> Professor,
>> Biochemistry, Molecular Biophysics, and
>> Biomedical Engineering
>> University of California, Irvine, CA 92697-3900
>>
>>
>> On Tue, Dec 3, 2013 at 4:54 AM, Jan-Philip Gehrcke
>> <jgehrcke.googlemail.com>wrote:
>>
>>> Hello Wes and Ray,
>>>
>>> thanks for your immediate feedback. I have seen that you have a patch
>>> since last week -- any news about when this will be made public?
>>>
>>> Also, I would like to discuss again the fact that radiopt is internally
>>> overridden to 2 when phiout == 1. In understand that there were some
>>> arguments for taking this decision. However, in my opinion, the current
>>> behavior must be changed for two reasons:
>>>
>>> - the current behavior contradicts documentation. The docs say that via
>>> radiopt = 1 PBSA uses the radii from the prmtop file. When no error
>>> appears as of a missing atom type, the user has no chance to realize
>>> that PBSA magically decided to use its hard-coded radii rather than the
>>> ones from the prmtop.
>>>
>>> - one could fix the docs or the code -- I vote for code. There should be
>>> a way to make PBSA take radii from the prmtop file *without being forced
>>> to manually change code*. And I don't see why this statement should be
>>> invalid depending on an output format.
>>>
>>>
>>> Thanks a lot for looking into this,
>>>
>>> Jan-Philip
>>>
>>>
>>> On 11/26/2013 07:49 PM, Ray Luo, Ph.D. wrote:
>>>> Jan-Philip,
>>>>
>>>> Okay, it turns out someone overwrote that change in the AmberTools13
>>>> release. Wes will need to issue a patch to put it back.
>>>>
>>>> Thanks a lot for reporting!
>>>>
>>>> Ray
>>>>
>>>> --
>>>> Ray Luo, Ph.D.
>>>> Professor,
>>>> Biochemistry, Molecular Biophysics, and
>>>> Biomedical Engineering
>>>> University of California, Irvine, CA 92697-3900
>>>>
>>>>
>>>> On Mon, Nov 25, 2013 at 10:52 AM, Jan-Philip Gehrcke <
>>>> jgehrcke.googlemail.com> wrote:
>>>>
>>>>> On 11/25/2013 07:13 PM, Ray Luo, Ph.D. wrote:
>>>>>> Jan-Philip,
>>>>>>
>>>>>> Indeed, when you ask for potential printout, pbsa would set radii as
>>> used
>>>>>> in VMD for consistency. That's why it overwrites your option in the
>>> input
>>>>>> file. If you insist in using radiopt=0, you can certainly comment
>>> the
>>>>> line
>>>>>> "radiopt=2" output and do it again.
>>>>>>
>>>>>> Ray
>>>>>
>>>>> Hey Ray,
>>>>>
>>>>> thanks for the confirmation. Two things left to be discussed:
>>>>>
>>>>> 1) I outcommented radiopt=2. In another attempt I left it as is and
>>>>> simply added the missing atom types to phi_aaradii(). In both cases,
>>>>> pbsa runs without error and produces the same result. So far so good.
>>>>> However, I end up with electrostatic potential values that seem to be
>>>>> too large by about two orders of magnitude.
>>>>>
>>>>> I have used pbsa earlier (1-2 years ago, from AmberTools12), with
>>>>> basically the same input (radiopt=0 is the only difference here)
>>> where I
>>>>> got values between 1 and 10 kcal/mol for most of the space regions
>>> near
>>>>> the protein surface. Now, I get something between 500 and 2000
>>> (whatever
>>>>> units). So, the units seem to have changed.
>>>>>
>>>>> Back in 2012, we also had an issue with the estatic potential units
>>> in
>>>>> the DX output file: http://archive.ambermd.org/201207/0466.html. At
>>>>> these times, Wes took care of this and released a bug fix making sure
>>>>> that the unit would be kcal/(mol*e). Did this accidentally change
>>> again?
>>>>>
>>>>> I have further proof that something changed. Found the
>>>>> AmberTools12-based pbsa binary that I used earlier. I used both, the
>>> old
>>>>> binaries and the new binaries with *the same input*. The first few
>>>>> values of the DX file corresponding to the old binary:
>>>>>
>>>>> 1.3067049144E+00 1.3232708113E+00 1.3399177244E+00
>>>>> 1.3566168815E+00 1.3733524282E+00 1.3900526596E+00
>>>>>
>>>>> Using the new AmberTools13-based pbsa binary:
>>>>>
>>>>> 5.4732938513E+02 5.5427901231E+02 5.6121241210E+02
>>>>> 5.6816177213E+02 5.7512228560E+02 5.8210218482E+02
>>>>>
>>>>> The factor is about 418.
>>>>>
>>>>>
>>>>> 2) I would suggest to change the current code behavior, because
>>>>> otherwise what happens is entirely intransparent to the user (we
>>> define
>>>>> radiopt=0 and magically get something undocumented). Would it be
>>>>> possible to use pymol/vmd-radii only in case radiopt was not
>>> explicitly
>>>>> defined in the input file (as kind of fallback)? Then, one should
>>> also
>>>>> add a corresponding warning/information to the output file, that
>>> these
>>>>> radii will now be applied. A missing-radius-error right after such an
>>>>> info line would make sense.
>>>>>
>>>>>
>>>>> Thanks for your support,
>>>>>
>>>>> Jan-Philip
>>>>>
>>>>>
>>>>>>
>>>>>> --
>>>>>> Ray Luo, Ph.D.
>>>>>> Professor,
>>>>>> Biochemistry, Molecular Biophysics, and
>>>>>> Biomedical Engineering
>>>>>> University of California, Irvine, CA 92697-3900
>>>>>>
>>>>>>
>>>>>> On Mon, Nov 25, 2013 at 10:04 AM, Jan-Philip Gehrcke <
>>>>>> jgehrcke.googlemail.com> wrote:
>>>>>>
>>>>>>> On 11/25/2013 06:38 PM, Jan-Philip Gehrcke wrote:
>>>>>>>
>>>>>>>> I used 'set default pbradii bondi' in leap when creating that
>>> topology
>>>>>>>> file. I am using AmberTools 13 with the latest patches for pbsa.
>>> The
>>>>>>>> whole pbsa input file is this:
>>>>>>>>
>>>>>>>> &cntrl
>>>>>>>> inp=0,
>>>>>>>> /
>>>>>>>> &pb
>>>>>>>> npbopt=0,
>>>>>>>> radiopt=0,
>>>>>>>> npbverb=0,
>>>>>>>> space=1,
>>>>>>>> fillratio=1.0,
>>>>>>>> sasopt=2,
>>>>>>>> phiout=1,
>>>>>>>> phiform=2,
>>>>>>>> npbverb=1,
>>>>>>>> /
>>>>>>>>
>>>>>>>> So, what is the issue here? Why is the given radius not used?
>>>>>>>>
>>>>>>>
>>>>>>> In pb_read.F90 there is this code:
>>>>>>>
>>>>>>> ! set phimap output options when requested
>>>>>>>
>>>>>>> if ( phiout == 1 ) then
>>>>>>> outphi = .true.
>>>>>>> radiopt = 2
>>>>>>> donpsa = .false.
>>>>>>> npopt = 0
>>>>>>> end if
>>>>>>>
>>>>>>>
>>>>>>> Since I am using phiout=1, radiopt becomes set to the
>>> (undocumented)
>>>>>>> value 2, which in turn makes pb_init() end up in this section:
>>>>>>>
>>>>>>> else if ( radiopt == 2 ) then
>>>>>>> call phi_aaradi( natom, isymbl, radi )
>>>>>>>
>>>>>>> So, if this quick analysis is right, the question would be why we
>>> set
>>>>>>> radiopt = 2 when phiout == 1.
>>>>>>>
>>>>>>>
>>>>>>> _______________________________________________
>>>>>>> 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/amber
>>>>>>
>>>>>
>>>>>
>>>>> _______________________________________________
>>>>> 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/amber
>>>>
>>>
>>>
>>> _______________________________________________
>>> 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/amber
>>
>>
>
>
>
> _______________________________________________
> 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/amber
Received on Sat May 24 2014 - 11:00:02 PDT
