您好, 欢迎来到 !    登录 | 注册 | | 设为首页 | 收藏本站

如何计算具有多个模型/构象的蛋白质的平均结构

5b51 2022/1/14 8:21:36 python 字数 5153 阅读 505 来源 www.jb51.cc/python

我有一个PDB文件'1abz'(https://files.rcsb.org/view/1ABZ.pdb),包含一个蛋白质结构的坐标,有23种不同的模型(编号为模型1-23).请忽略标题备注,有趣的信息从第276行开始,标题为“模型1”.我想计算蛋白质的平均结构.蛋白质的PDB文件包含多个构象/模型,我想计算每个残基的单个原子的平均坐标,这样我最终得到一个构

概述

我有一个PDB文件’1abz'(https://files.rcsb.org/view/1ABZ.pdb),包含一个蛋白质结构的坐标,有23种不同的模型(编号为模型1-23).请忽略标题备注,有趣的信息从第276行开始,标题为“模型1”.

我想计算蛋白质的平均结构.蛋白质的PDB文件包含多个构象/模型,我想计算每个残基的单个原子的平均坐标,这样我最终得到一个构象/模型.

我无法弄清楚如何使用Biopython来做到这一点,所以我尝试使用Pandas计算平均坐标.我想我已经设法计算平均值,但现在的问题是我有一个不再是PDB格式的csv文件,所以我无法将此文件加载到PyMol中.

我的问题是,如何将我的csv文件转换为PDB格式.更好的是,如何在不损害原始pdb文件格式的情况下获取Biopython或Python中的平均坐标?

这是我用来计算熊猫平均坐标的代码.

#I first converted my pdb file to a csv file

import pandas as pd
import re

pdbfile = '1abz.pdb'
df = pd.DataFrame(columns=['Model','Residue','Seq','Atom','x','y','z']) #make dataframe object
i = 0 #counter

b = re.compile("MODEL\s+(\d+)")
regex1 = "([A-Z]+)\s+(\d+)\s+([^\s]+)\s+([A-Z]+)[+-]?\s+([A-Z]|)"
regex2 = "\s+(\d+)\s+([+-]?\d+\.\d+\s+[+-]?\d+\.\d+\s+[+-]?\d+\.\d+)"
reg = re.compile(regex1+regex2)

with open(pdbfile) as f:
    columns = ('label','ident','atomName','residue','chain','sequence','z','occ','temp','element')
    data = []
    for line in f:
        n = b.match(line)
        if n:
            modelNum = n.group(1)

        m = reg.match(line)
        if m:
            d = dict(zip(columns,line.split()))
            d['model'] = int(modelNum)
            data.append(d)

df = pd.DataFrame(data)
df.to_csv(pdbfile[:-3]+'csv',header=True,sep='\t',mode='w')

#Then I calculated the average coordinates

df = pd.read_csv('1abz.csv',delim_whitespace = True,usecols = [0,5,7,8,10,11,12]) 
df1 = df.groupby(['atomName','sequence'],as_index=False)['x','z'].mean()
df1.to_csv('avg_coord.csv',mode='w')

首先使用所有模型解析pdb文件

import Bio.PDB
import numpy as np

parser = Bio.PDB.PDBParser(QUIET=True)  # Don't show me warnings
structure = parser.get_structure('1abz','1abz.pdb')  # id of pdb file and location

好的,现在我们有了文件内容,假设你的所有模型都有相同的原子,那么得到一个列表,其中包含每个原子的唯一标识符(例如:链残基pos原子名):

atoms = [a.parent.parent.id + '-' + str(a.parent.id[1]) + '-' +  a.name for a in structure[0].get_atoms() if a.parent.id[0] == ' ']  # obtained from model '0'

请注意,我忽略了a.parent.id [0] ==”的杂质残基.现在让我们得到每个原子的平均值:

atom_avgs = {}
for atom in atoms:
    atom_avgs[atom] = []
    for model in structure:
        atom_ = atom.split('-')
        coor = model[atom_[0]][int(atom_[1])][atom_[2]].coord
        atom_avgs[atom].append(coor)
    atom_avgs[atom] = sum(atom_avgs[atom]) / len(atom_avgs[atom])  # average

现在让我们使用结构的一个模型创建新的pdb:

ns = Bio.PDB.StructureBuilder.Structure('id=1baz')  #  new structure
ns.add(structure[0])  # add model 0
for atom in ns[0].get_atoms():
    chain = atom.parent.parent
    res = atom.parent
    if res.id[0] != ' ':
        chain.detach_child(res)  # detach hetres
    else:
        coor = atom_avgs[chain.id + '-' + str(res.id[1]) + '-' + atom.name]
        atom.coord = coor

现在让我们写一下pdb

io = Bio.PDB.PDBIO()
io.set_structure(ns)
io.save('new_1abz.pdb')

总结

以上是编程之家为你收集整理的如何计算具有多个模型/构象的蛋白质的平均结构全部内容,希望文章能够帮你解决如何计算具有多个模型/构象的蛋白质的平均结构所遇到的程序开发问题。


如果您也喜欢它,动动您的小指点个赞吧

除非注明,文章均由 laddyq.com 整理发布,欢迎转载。

转载请注明:
链接:http://laddyq.com
来源:laddyq.com
著作权归作者所有。商业转载请联系作者获得授权,非商业转载请注明出处。


联系我
置顶