请问如何用脚本实现两个分子对齐align？

wei

各位好！我想用Tcl脚本实现一批分子的对齐，再计算其RMSD。我选择了第一个分子作为参考（index=0）,后续分子依次和这第一个分子对齐再做RMSD。RMSD功能已实现，但对齐未能实现，请各位指教。以下为部分代码。

set j 1
while {$j <=$num_mol} {
# select ref and compare molecule
set mol_comare [atomselect 1 "backbone"]
set mol_ref [atomselect 0 "backbone"]
# compute the transformation
set trans_mat [measure fit $mol_compare $mol_ref]
# do the alignment
$mol_com move $trans_mat
set mol_rmsd_all [measure rmsd $mol_ref $mol_com]
puts $myfile "$j    $mol_rmsd_all "
incr j 1
    }
  close $myfile

lyj714

有很多程序都可以实现，不一定要用tcl，因为如果分子大，叠合过程可能慢。比如crest有个-rmsd的选项就可以直接算rmsd(一般的rmsd都是指的先叠合再计算)。
但是你的描述也可能是指的是不同分子之间计算rmsd，如果你指的是这种，那就得考虑更多，比如选着相同数目的原子进行叠合等等。

wei

感谢回复！事实上，vmd的GUI界面也是可以完成对齐和RMSD操作的，但由于种种原因，我目前的应用场景需要用Tcl脚本+vmd实现。

lyj714

wei 发表于 2021-4-12 21:29
感谢回复！事实上，vmd的GUI界面也是可以完成对齐和RMSD操作的，但由于种种原因，我目前的应用场景需要用Tc ...

https://www.ks.uiuc.edu/Research/vmd/vmd-1.9.3/ug/node203.html
粗略看了下你写的那个，不是做了叠合了吗。参考分子不变，比较分子每次更新帧数即可。

ene

1. #----------------------------------------------------
   
2. #                     PARAMETERS
   
3. #----------------------------------------------------
   
4. set align_frame     0
   
5. set align_reference "protein"
   
6. set align_selection "protein"
   
7. #----------------------------------------------------
   
8. set nf [molinfo top get numframes]
   
9. set f0 [atomselect top "$align_reference" frame $align_frame]
   
10. set f  [atomselect top "$align_reference"]
    
11. set sel [atomselect top "$align_selection"]
    
12. for { set i 0 } { $i < $nf } { incr i } {
    
13.     $f frame $i
    
14.     $sel frame $i
    
15.     set mv [measure fit $f $f0]
    
16.     $sel move $mv
    
17. }
    
18. puts "All Done!"
    

复制代码

wei

谢谢各位！我测试一下。