VMDで二次構造の表示をステップごとに更新する

#はじめに

タンパク質のフォールディングやアンフォールディングでは時間経過と共にαヘリックスやβシートなどの二次構造が変化しますが、VMDはデフォルトではこれを描画してくれません。
これを描画する方法をご紹介します。スクリプトは参考リンクからそのまま使っています。

ちなみに、VMDについては以下の記事も書いたので参考にしてください。

#デフォルトの表示

Webから適当なペプチド(PDB : 1e0m)を持ってきてMD計算を実行します。初期状態では以下のようにβシートを形成しています。

フォールディングは大変なので、700Kの超高温で1ns計算して構造をバラバラにしました。なお、計算にはGROMACSを利用しました。以下がデフォルトのVMDの表示です。

すべての状態がランダムコイルで表示されています。これは最初に読み込むgroファイルに依存するということです。フォールドした状態のgroファイルを読み込むと以下のようになります。

すべての状態でβシートの表示になってしまいました。これだと構造の変化がわかりにくいです。

#二次構造の表示を更新するスクリプト

参考リンクに書いてある方法を紹介します。

##方法1

以下のスクリプトをコンソールに打ち込んで実行します。自動で実行したければ、.vmdrcなど起動時に読み込むファイルに書き込んでもOKです。

proc structure_trace {name index op} {
      vmd_calculate_structure $index
}

trace variable vmd_frame w structure_trace

これを使うと表示が毎フレームごとに変化します。

問題解決です。

##方法2

方法1でも十分な気もしますが、より高速なのが方法2です。

• 一度計算した二次構造をキャッシュに保存しておくことで何度も同じフレームの二次構造を計算することを防ぎ、
• 複数分子を読み込んだ場合でも、特定の分子だけ二次構造を計算させることができる

のがメリットだそうです。
何はともあれ、方法1と同様にしてコンソールに打ち込むか、.vmdrcに記述して実行します。

# Cache secondary structure information for a given molecule

# reset the secondary structure data cache
proc reset_sscache {{molid top}} {
    global sscache_data
    if {! [string compare $molid top]} {
      set molid [molinfo top]
    }
    if [info exists sscache_data($molid)] {
        unset sscache_data
    }
}

# start the cache for a given molecule
proc start_sscache {{molid top}} {
    if {! [string compare $molid top]} {
      set molid [molinfo top]
    }
    global vmd_frame
    # set a trace to detect when an animation frame changes
    trace variable vmd_frame($molid) w sscache
    return
}

# remove the trace (need one stop for every start)
proc stop_sscache {{molid top}} {
    if {! [string compare $molid top]} {
      set molid [molinfo top]
    }
    global vmd_frame
    trace vdelete vmd_frame($molid) w sscache
    return
}


# when the frame changes, trace calls this function
proc sscache {name index op} {
    # name == vmd_frame
    # index == molecule id of the newly changed frame
    # op == w
    
    global sscache_data

    # get the protein CA atoms
    set sel [atomselect $index "protein name CA"]

    ## get the new frame number
    # Tcl doesn't yet have it, but VMD does ...
    set frame [molinfo $index get frame]

    # see if the ss data exists in the cache
    if [info exists sscache_data($index,$frame)] {
        $sel set structure $sscache_data($index,$frame)
        return
    }

    # doesn't exist, so (re)calculate it
    vmd_calculate_structure $index
    # save the data for next time
    set sscache_data($index,$frame) [$sel get structure]

    return
}

二次構造を計算するmolid（デフォルトでは0）を指定して実行します。

start_sscache 0

このコマンドごと.vmdrcに書いてしまうこともありだと思います。なお、表示は変わりませんので、動画は割愛します。

毎フレーム更新をやめるには下記コマンドを実行します。

stop_sscache 0

ちなみに私は今のところ更新をやめる場面が思いつきませんが、少し処理を軽くしたければ、実行してもよいかもしれません。

以上、長くなりましたが皆様の参考になれば幸いです。

#参考