6
$\begingroup$

I am trying to write a tcl script to measure percentage helicity at each frame of a VMD .dcd trajectory. I found the following discussion on the VMD mailing list

to calculate the secondary structure for a selection at a given timestep, you should go to that timestep, and the use the command mol reanalyze top which runs that timestep through stride (among other things). Once this is done, you can access the secondary structure of an atomselect object with $sel get structure

I have tried to implement it as follows,

set outfile [open ./percent_helix.dat w]
set lookup {H G I}
set frame_num [molinfo top get numframes]
set full [atomselect top "name CA"]
set len [llength [$full get resid]]

for {set i 0} {$i < $frame_num} {incr i} {
    set helix 0
    $full frame i
    $full update
    mol reanalyze top
    set struc_string [$full get structure]
    foreach letter $lookup {
        set temp [expr {[llength [split $struc_string $letter]] - 1}]
        incr helix $temp
    }
    set percent [expr {double($helix) / double($len) * 100}]
    puts $outfile "$i\t$percent"
}
$full delete
close $outfile

However, the output i am the same percentage value at each timestep. I have figured that the reason behind that is the command mol reanalyze top is not changing the secondary structure values at each frame. What should I correct in my script?

Improve this question
$\endgroup$

1 Answer 1

Reset to default
6
$\begingroup$

I have found the solution to the problem. One needs to use the animate command to change the frame in the top molecule in VMD. The corrected script which produces the desired output is as follows:

set outfile [open ./percent_helix.dat w]
set lookup {H G I}
set frame_num [molinfo top get numframes]
set full [atomselect top "name CA"]
set len [llength [$full get resid]]
$full delete

for {set i 0} {$i < $frame_num} {incr i} {
    animate goto $i
    set sel [atomselect top "name CA"]
    mol ssrecalc top
    set struc_string [$sel get structure]
    set helix 0
    foreach letter $lookup {
        set temp [expr {[llength [split $struc_string $letter]] - 1}]
        incr helix $temp
    }
    set percent [expr {double($helix) / double($len) * 100}]
    puts $outfile "$i\t$percent"
    $sel delete
}
close $outfile
Improve this answer
$\endgroup$
7
  • $\begingroup$ +1. Thanks for coming back and sharing the solution so that others in the future might benefit from it! I'm sure they'll appreciate it! $\endgroup$
    – Nike Dattani
    Jun 19, 2021 at 3:38
  • $\begingroup$ @NikeDattani Thanks. I will mark this answer as accepted if no one else answers in a few days. A better way to do this may pop up in a new answer. $\endgroup$
    – Roni Saiba
    Jun 19, 2021 at 10:00
  • $\begingroup$ No point in marking your own answer as accepted: you only gain the +2 for marking other people's answers as accepted, and the accepted answer only goes to the top if it's by someone else. Accepting your own answer will just discourage others from trying to answer. $\endgroup$
    – Nike Dattani
    Jun 19, 2021 at 14:26
  • $\begingroup$ One more thing you can do is calculate the secondary structure at each frame from timeline plugin of vmd and use other coding platform to look up for letter H (for alpha helix). $\endgroup$
    – Vasista
    Jun 22, 2021 at 7:14
  • $\begingroup$ @AdupaVasista Hi. I think you mailed me on the vmd-l list as well. The problem with timeline tool is that I don't have any way of exporting it in a useful manner. $\endgroup$
    – Roni Saiba
    Jun 22, 2021 at 9:45

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .