计算原子坐标之间的距离

Question

我有一个文本文件,如下所示

ATOM    920  CA  GLN A 203      39.292 -13.354  17.416  1.00 55.76           C 
ATOM    929  CA  HIS A 204      38.546 -15.963  14.792  1.00 29.53           C
ATOM    939  CA  ASN A 205      39.443 -17.018  11.206  1.00 54.49           C  
ATOM    947  CA  GLU A 206      41.454 -13.901  10.155  1.00 26.32           C
ATOM    956  CA  VAL A 207      43.664 -14.041  13.279  1.00 40.65           C 
.
.
.

ATOM    963  CA  GLU A 208      45.403 -17.443  13.188  1.00 40.25           C  

我想计算两个α碳原子之间的距离,即计算第一和第二原子之间的距离,然后计算第二和第三原子之间的距离等等......两个原子之间的距离可以表示为:distance = sqrt((x1-x2)^2+(y1-y2)^2+(z1-z2)^2) .

列7,8和9分别代表x,y和z坐标.我需要打印距离和相应的残差对(第4列),如下所示.(距离值不是真实的)

GLN-HIS   4.5
HIS-ASN   3.2
ASN-GLU   2.5

如何用perl或python进行计算？

Answer 1

如果您的数据是由空格分隔的,那么简单split就可以完成这项任务.缓冲线以顺序地将它们相互比较.

use strict;
use warnings;

my @line;
while (<>) {
    push @line, $_;            # add line to buffer
    next if @line < 2;         # skip unless buffer is full
    print proc(@line), "\n";   # process and print 
    shift @line;               # remove used line 
}

sub proc {
    my @a = split ' ', shift;   # line 1
    my @b = split ' ', shift;   # line 2
    my $x = ($a[6]-$b[6]);      # calculate the diffs
    my $y = ($a[7]-$b[7]);
    my $z = ($a[8]-$b[8]);
    my $dist = sprintf "%.1f",                # format the number
                   sqrt($x**2+$y**2+$z**2);   # do the calculation
    return "$a[3]-$b[3]\t$dist"; # return the string for printing
}

输出(带样本数据):

GLN-HIS 3.8
HIS-ASN 3.8
ASN-GLU 3.9
GLU-VAL 3.8

如果您的数据以制表符分隔,则可以拆分/\t/而不是' '.

Answer 2

不要在空格上拆分

这里给出的其他答案是一个有缺陷的假设 - 坐标将以空格分隔.根据PDB规范ATOM,这不是必需的情况:PDB记录值由列索引指定,并且可以相互流动.例如,您的第一条ATOM记录为:

ATOM    920  CA  GLN A 203      39.292 -13.354  17.416  1.00 55.76           C 

但这也完全有效:

ATOM    920  CA  GLN A 203      39.292-13.3540  17.416  1.00 55.76           C 

更好的方法

由于列指定的索引以及PDB文件中可能出现的其他问题的数量,您不应编写自己的解析器.PDB格式很乱,并且有很多特殊情况和格式错误的文件需要处理.相反,使用已经为您编写的解析器.

我喜欢Biopython的PDB.PDBParser.它将为您解析Python对象的结构,并提供方便的功能.如果您更喜欢Perl,请查看BioPerl.

PDB.Residue对象允许按名称键入访问Atoms,PDB.Atom对象使-运算符超载以返回两个Atom之间的距离.我们可以用它来编写简洁明了的代码:

码

from Bio import PDB
parser = PDB.PDBParser()

# Parse the structure into a PDB.Structure object
pdb_code = "1exm"
pdb_path = "pdb1exm.ent"
struct = parser.get_structure(pdb_code, pdb_path)

# Grab the first two residues from the structure
residues = struct.get_residues()
res_one = residues.next()
res_two = residues.next()

try:
    alpha_dist = res_one['CA'] - res_two['CA']
except KeyError:
    print "Alpha carbon missing, computing distance impossible!"