0

我正在尝试浏览一个 .pdb 文件，计算蛋白质复合物链 A 和 B 上不同残基的 α 碳原子之间的距离，然后将距离连同链标识符和残基编号一起存储在字典中。

例如，如果在链 A 上的残基 100 上发现第一个 α 碳（“CA”），并且它所结合的在链 BI 上的残基 123 上，我希望我的字典看起来像 d={(A, 100) :[B, 123, distance_between_atoms]}

from Bio.PDB.PDBParser import PDBParser
parser=PDBParser()
struct = parser.get_structure("1trk", "1trk.pdb")

def getAlphaCarbons(chain):
    vec = []
    for residue in chain:
        for atom in residue:
            if atom.get_name() == "CA":
                vec = vec + [atom.get_vector()]
    return vec

def dist(a,b):
    return (a-b).norm()


chainA = struct[0]['A']
chainB = struct[0]['B']

vecA = getAlphaCarbons(chainA)
vecB = getAlphaCarbons(chainB)

t={}
model=struct[0]

for model in struct:
    for chain in model:
        for residue in chain:
            for a in vecA:
               for b in vecB:
                if dist(a,b)<=8:
                    t={(chain,residue):[(a, b, dist(a, b))]}

     break
print t

它已经运行程序很长时间了，我不得不中止运行（我在某处做了一个无限循环吗？？）

我也试图这样做：

t = {i:[((a, b, dist(a,b)) for a in vecA) for b in vecB if dist(a, b) <= 8] for i in chainA}
print t

但它以以下格式打印有关残留物的信息：

<Residue PHE het=  resseq=591 icode= >: []

它没有打印与距离相关的任何内容。

非常感谢，我希望一切都清楚。

4

1 回答 1

0

强烈建议在计算距离时使用 C 库。我将mdtraj用于这类事情，它比 BioPython 中的所有 for 循环都快得多。

要获得所有对 alpha-Carbons：

import mdtraj as md
def get_CA_pairs(self,pdbfile):
  traj = md.load_pdb(pdbfile)
  topology = traj.topology 
  CA_index = ([atom.index for atom in topology.atoms if (atom.name == 'CA')])
  pairs=list(itertools.combinations(CA_index,2))
return pairs

然后，为了快速计算距离：

  def get_distances(self,pdbfile,pairs):
  #returns list of resid1, resid2,distances between CA-CA
   traj = md.load_pdb(pdbfile)
   pairs=self.get_CA_pairs(pdbfile)
   dist=md.compute_distances(traj,pairs)
#make dictionary you desire. 
   dict=dict(zip(CA, pairs))
  return dict

这包括所有 α-碳。mdtraj 中也应该有一个链标识符来从每个链中选择 CA。

于 2017-07-22T16:19:16.213 回答