CompCB@Python

Biopython使用4:PDB结构

2019-01-19  本文已影响7人  爱折腾的大懒猪

Biopython 第十一章:走向3D:PDB模块

PDBParser和MMCIFParser

PDBParser 解析器

使用p=PDBParser(参数)来构建一个PDB解释器. 接收以下参数 (常用QUIET=True ) :

>>> from Bio.PDB.PDBParser import PDBParser
>>> p = PDBParser()
>>> s = p.get_structure('1d3z','1D3Z.pdb')

解释器对象有以下方法和属性:

读取后获得的Structure对象有多项属性. 其中header属性同PDBParser获得. id属性也是解析时指定.

Bio.PDB.parse_pdb_header(file) 也可以用来读取文件头. 包括head, deposition_date, release_date, structure_method, resolution, structure_reference, journal_reference, author and compound.

读取mmCIF文件

类似PDBParser, 参数支持structure_builder与QUIET. 解释器支持get_structure方法, 有line_counter属性. 要获取类似文件头信息, 要用MMCIF2Dict 方法. 获取的也是一个字典.

>>> from Bio.PDB.MMCIFParser import MMCIFParser
>>> parser = MMCIFParser()
>>> structure = parser.get_structure('1fat', '1fat.cif')
### 要解释mmCIF的额外信息, 需要MMCIF2Dict
>>> from Bio.PDB.MMCIF2Dict import MMCIF2Dict
>>> mmcif_dict = MMCIF2Dict('1fat.cif')

PDBIO类对象保存文件

该类主要用于输出PDB文件. from Bio.PDB import PDBIO来加载后, 直接构建一个io对象. 用set_structure来指定结构后, 再用save方法来保存.

如果想保存部分结构而非全部结构, 可以用Select类来实现. Select支持四种方法来进行选择: accept_model(model), accept_chain(chain), accept_residue(residue), accept_atom(atom). 默认情况下, 每种方法的返回值都为1(即被保存), 可以继承该类自行构建子类来构建自己的选择器, 从而实现强大的功能. 当经常需要选择性保存时, 这种做法比较方便.

>>> from Bio.PDB import PDBIO
### 简单的读写
>>> io = PDBIO()
>>> io.set_structure(s)
>>> io.save('out.pdb')

#### 构建一个选择器来选择GLY. 
>>> class GlySelect(Select):
...     def accept_residue(self, residue):
...         if residue.get_name()=='GLY':
...             return True
...         else:
...             return False
...
>>> io.save('gly_only.pdb', GlySelect())

Structure类对象

Biopython的结构采用SMCRA体系架构（Structure/Model/Chain/Residue/Atom，结构/模型/链/残基/原子）.

一般地，一个实体子类（即原子，残基，链，模型）能通过id作为键来从父类（分别为残基，链，模型，结构）中提取。可以使用len(object)来获取实体的子类的数量.

可以使用child_list = parent_entity.get_list()来获取子对象的列表 (顺序有讲究), 也可以用get_parent()来获取父类. 对于残基和原子, 其id是一个元组, 比较奇怪, 建议获取全列表后再用索引获取.

在SMCRA的所有层次水平，你还可以提取一个完整id 。完整id是包含所有从顶层对象（结构）到当前对象的id的一个元组。一个残基对象的完整id可以这么得到：

>>> full_id = residue.get_full_id()
>>> print full_id
("1abc", 0, "A", ("", 10, "A"))
>>> atom0=residue.get_list()[0]
>>> print atom0.get_full_id()
("1abc", 0, "A", ("", 10, "A"), ('N', ' '))

一些通用的方法:

以下属性方法基本都有, 但一般不用. 除非很清楚父子类间操作.

## transform(rot, tran) 用法
## rot是右乘旋转矩阵, 3x3 array
## tran是平移向量, 大小3的向量
>>> rotation = rotmat(pi, Vector(1, 0, 0))
>>> translation = array((0, 0, 1), 'f')
>>> entity.transform(rotation, translation)

Structure/Model/Chain

这三层次都比较高. Structure的ID靠指定字符串, Model的ID源自于文件位置, 从0开始的整数. 一般只有一个模型, 有些结构含有多个模型(例如NMR). Chain的ID一般是链的标识符, 单字符(一般是字母), 如'A'. 每个模型的链具有唯一ID. 有多个模型时会出现多个相同ID的链.

Residue

残基的ID是一个三元组: (het, resseq, icode):

残基存在disordered状态. 例如存在Thr 80 A, Ser 80 B, Asn 81, 这个80残基可能有突变, 有另外一种可能情况. 这个会在icode中反映出来.

残基类对象有这些特殊方法 :

Bio.PDB.is_aa(residue, standard=False) 方法可以检查氨基酸. residue可以使残基对象, 也可以是3字母字符串. standard设置True, 只检查20个标准氨基酸. 例如FME默认True, 设置后就False.

Atom

原子属性有很多, 除了之前很多SMCRA对象的共有方法外, 还有很多方法获取的是对应的属性.

原子的ID就是他的名称name. 要注意唯一性. 一般就是PDB文件中原子名称去除空格来创建. 而实际在PDB原子名称是可以有空格的.例如'CA '代表的是钙而非Cα (' CA '). 此时为了防止歧义, 会保留空格 (所以要慎防这种奇葩情况的出现, 尽管Ca不常见).

注意: 原子序号项和电荷没有存储!!!

Bio.PDB.Vector 模块及相关运算

除了上述方法, Bio.PDB里面也有一些和向量操作相关的方法:

####  获得CB的估算位置
## 获得原子坐标向量
>>> n = residue['N'].get_vector()
>>> c = residue['C'].get_vector()
>>> ca = residue['CA'].get_vector()
## 计算C-CA, N-CA的向量. 
>>> n = n - ca
>>> c = c - ca
## 将N-CA向量绕C-CA轴旋转-120度, 大致可以到CB的位置
# 计算旋转矩阵用于沿CA-C旋转-120度
>>> rot = rotaxis(-pi * 120.0/180.0, c)
# 将旋转矩阵用于N-CA, 从而获得旋转后N-CA向量
>>> cb_at_origin = n.left_multiply(rot)
# 根据旋转后向量, 加上CA的位置, 可以计算出CB的位置.
>>> cb = cb_at_origin+ca

紊乱Disordered残基和原子

紊乱残基和原子分别使用DisorderedResidue和DisorderedAtom来处理. 在列出链内残基时, 相应的会显示紊乱残基, 会选择其中一个残基显示. 但是, 当有A,B多个可选残基时, 他不一定选的是A的残基, 默认会选择最后的一个残基作为缺省!

最好获得链的残基的方法是使用get_unpacked_list()来获取残基, 此时显示的会是A链的残基且是Residue而不是DisorderedResidue.

一个PDB例子是1EN2. 该例子残基10, 14, 16, 80, 81均是紊乱残基, 而1, 90-93是het残基, 94-169是水.

对于紊乱残基, 会有一些特殊的方法:

from Bio.PDB import PDBParser
p=PDBParser()
dm=p.get_structure('1en2','pdb1en2.ent')
dc=dm[0]['A']
dr=dc.get_list()[9] # 或者 dc[10]
dr.disordered_get_list()
dr.disordered_select('SER')

遍历SMCRA

>>> p = PDBParser()
>>> structure = p.get_structure('X', 'pdb1fat.ent')
>>> for model in structure:
...     for chain in model:
...         for residue in chain:
...             for atom in residue:
...                 print atom
>>> atoms = chain.get_atoms()
>>> for atom in atoms:
...     print atom

另外, Bio.PDB.Selection里的一些方法可以帮助选择出相应的一些实例.

例如:

## 选择一个结构的所有残基, 原子的话相应 'A'
## A=atom, R=residue, C=chain, M=model, S=structure
>>> res_list = Selection.unfold_entities(structure, 'R')
## 这个方法也可以帮助选择一些原子对应的唯一的残基或链等. 
>>> chain_list = Selection.unfold_entities(atom_list, 'C')
## 有专用于选择相应实例的父类实例的
>>> residue_list = Selection.get_unique_parents(atom_list)
## 这个方法用来选择出unique的实例
>>> unique_atom=Selection.uniqueify(atom_list)

多肽对象 Polypeptide

这是一个多肽相关的类. 多肽类实际是列表list继承的类, 储存一系列的残基对象.

该类有两个Builder, 一个是CaPPBuilder, 一个是PPBuilder (也可以在Bio.PDB内获得), 前者根据Cα-Cα距离来建立多肽, 后者通过C-N距离来构建多肽.

# Using C-N
>>> ppb=PPBuilder()
>>> for pp in ppb.build_peptides(structure):
...     print pp.get_sequence()
...
# Using CA-CA
>>> ppb=CaPPBuilder()
>>> for pp in ppb.build_peptides(structure):
...     print pp.get_sequence()

RCSB 蛋白结构数据库

Biopython的PDB模块还能在PDB数据库进行下载. 通过PDBList对象来实现.

PDBList(self, server='ftp://ftp.wwpdb.org', pdb='/Users/hom', obsolete_pdb=None, verbose=True)

下载PDB文件: retrieve_pdb_file(pdb_code, obsolete=False, pdir=None, file_format=None, overwrite=False)

可以在命令行运行python PDBList.py 1fat 来下载. python PDBList.py all /data/pdb 可以下载整个PDB数据库到/data/pdb. 使用-d选项可以将下载所有文件到同一目录而非对应子目录.

>>> from Bio.PDB import PDBList
>>> pdbl = PDBList()
>>> pdbl.retrieve_pdb_file('1FAT', file_format='pdb')

以下是PDBList对象其他一些方法:

上一篇 下一篇

猜你喜欢

热点阅读