VMD-L Mailing List

From: Peter Freddolino (petefred_at_ks.uiuc.edu)
Date: Thu Jan 11 2007 - 22:33:55 CST

Could you please elaborate exactly what selections are being fed to the
measure sasa commands in each of these cases? Certainly, if you mean
that the first selection being given was only the first five residues,
you'll get different results because portions of those residues will be
in contact with (and buried by) other regions, and you should use the
whole protein to get reasonable results.
Peter

alexandra.marques_at_fc.up.pt wrote:
>
> Hi again,
>
> I used the pdb file with and without hydrogens and yes, the servers
> donīt use
> hydrogens.
> I think I found the problem. If I use only selected residues of the
> protein to
> do the calculation, the sasa values are different and much higher than
> if I use
> the whole protein residues. For instance, for the first 5 aa:
>
> Sasa (using the whole protein) Sasa (using only the 5 aa in the
> selection)
> 0 124.216300964 148.595306396
> 1 2.31623363495 127.548957825
> 2 174.497711182 204.40007019
> 3 55.33360672 87.0589981079
> 4 49.2410736084 174.806564331
>
> And the sasa values computed using the whole protein are similar to
> the ones
> reported by the servers. Also, now I can obtain zero values for the
> buried
> residues.
> Do you know the reason for this difference? Do I have to use the whole
> protein
> residues to get correct sasa values?
>
> Thanks
> Alexandra
>
>
> Quoting Peter Freddolino <petefred_at_ks.uiuc.edu>:
>
>> Hi Alexandra,
>> as a first check, could you verify that you're using input that includes
>> explicit hydrogens? The servers you mentioned ignore them (I believe),
>> but VMD doesn't. As described in the VMD user's manual
>> (http://www.ks.uiuc.edu/Research/vmd/current/ug/node124.html), the sasa
>> is calculated by extending the radius of each atom in the target
>> selection by the solvent radius and calculating the area of the
>> resultant spheres which is not overlapping another atom. By this
>> definition even buried residues are likely to have some amount of
>> exposed surface (since they do not fit together perfectly); the question
>> of whether or not a given residue is actually in a solvent-accessible
>> crevice or not is not addressed here (this is part of the purpose of the
>> -restrict option).
>>
>> Peter
>>
>> alexandra.marques_at_fc.up.pt wrote:
>>> Hi
>>>
>>> I am using a script (_sasa_res.tcl) that I found in the mailing list
>>> archives to
>>> calculate the sasa by residue. However the values obtained are much
>>> higher than
>>> expected and buried residues do not have zero values. For instance,
>>> I compared
>>> the values obtained for a pdb file with those reported by the PQS or
>>> EBI PISA
>>> servers and they are very different.
>>> I am sending the script in attachment. Could someone please take a
>>> look to it
>>> and/or suggest me another way to calculate sasa by residue? What is
>>> the method
>>> used by vmd to calculate SASA?
>>>
>>> Thanks in advance
>>> Alexandra
>>>
>>> -------------------------------------------------------------
>>> A FCUP utiliza o sistema de webmail Horde/IMP (www.horde.org)
>>>
>>> Visite: http://www.fc.up.pt/
>>>
>>> ------------------------------------------------------------------------
>>>
>>>
>>> # version 1.2
>>> # calculate SASA for every residue of a selection
>>> # usage: getAllResSASA "selection" probe_radius <startframe <endframe>>
>>>
>>> proc getAllResSASA { selection radius args } {
>>> puts "SASA value calculator for $selection residues, probe radius
>>> $radius"
>>>
>>> # get number of frames
>>> set numframes [molinfo top get numframes]
>>> puts "total number of frames in trajectory: $numframes"
>>> if {[llength $args] == 0} {
>>> set startframe 1
>>> set endframe $numframes
>>> }
>>> if {[llength $args] == 1} {
>>> set startframe [lindex $args 0]
>>> set endframe $numframes
>>> set stepframe 1
>>> if {$startframe <= 0 || $startframe > $numframes} {
>>> puts "illegal value of startframe, changing to first frame"
>>> set startframe 1
>>> }
>>> }
>>> if {[llength $args] >= 2} {
>>> set startframe [lindex $args 0]
>>> set endframe [lindex $args 1]
>>> if {$startframe <= 0 || $startframe > $numframes} {
>>> puts "illegal value of startframe, changing to first frame"
>>> set startframe 1
>>> }
>>> if {$endframe < $startframe || $endframe > $numframes} {
>>> puts "illegal value of endframe, changing to last frame"
>>> set endframe $numframes
>>> }
>>> }
>>>
>>> set totframes [expr ($endframe - $startframe + 1)]
>>> puts "analysis will be performed on $totframes frame(s)
>>> ($startframe to $endframe)"
>>>
>>> # get list of residues
>>> set allsel [atomselect top $selection]
>>> set residlist [lsort -unique -integer [ $allsel get residue ]]
>>>
>>> # resid map for every atom
>>> set allResid [$allsel get residue]
>>> # resname map for every atom
>>> set allResname [$allsel get resname]
>>> # create resid->resname map
>>> foreach resID $allResid resNAME $allResname {
>>> set mapResidResname($resID) $resNAME
>>> }
>>>
>>> # set sasa for every residue to 0.0
>>> foreach r $residlist {
>>> set resSASA($r) 0.0
>>> # puts "setting 0.0 for residue $r"
>>> }
>>>
>>> # now subtract 1 from all frame indexes - numbering starts with 0
>>> set startframe [expr $startframe - 1]
>>> set endframe [expr $endframe - 1]
>>>
>>> puts "go and get coffee..."
>>> for {set i $startframe; set d 1} {$i <= $endframe} {incr i; incr d} {
>>> # show activity - one dot for every frame
>>> puts -nonewline "."
>>> if { [expr $d % 50] == 0 } {
>>> puts " "
>>> }
>>> flush stdout
>>>
>>> $allsel frame $i
>>> $allsel update
>>>
>>> foreach r $residlist {
>>> set sel [atomselect top "residue $r" frame $i]
>>> set rsasa [measure sasa $radius $allsel -restrict $sel]
>>> # set user value for this frame
>>> $sel set user $rsasa
>>> $sel delete
>>> set valuesasa 0.0
>>> foreach {tmp valuesasa} [split [array get resSASA $r] ] break
>>> # puts "residue $r, sasa: $rsasa, old: $valuesasa"
>>> # if sasa is above threshold, store it
>>> #if {$rsasa >= $sasa_limit} {
>>> set resSASA($r) [ expr { $valuesasa + $rsasa/$totframes} ]
>>> # }
>>> }
>>> }
>>>
>>> set fw [open "res_sasa.dat" w]
>>> foreach r $residlist {
>>> foreach {tmp resName} [split [array get mapResidResname $r] ] break
>>> foreach {tmp valuesasa} [split [array get resSASA $r] ] break
>>> puts $fw "$r $resName $valuesasa"
>>> }
>>> close $fw
>>> puts "done"
>>>
>>> # uncomment following code to change coloring method to "user based"
>>> #mol modcolor 0 [molinfo top] User
>>> #mol colupdate 0 [molinfo top] 1
>>> #mol scaleminmax [molinfo top] 0 auto
>>> }
>>>
>>
>
>
>
> -------------------------------------------------------------
> A FCUP utiliza o sistema de webmail Horde/IMP (www.horde.org)
>
> Visite: http://www.fc.up.pt/
>
>