专栏名称: 生信杂谈
生物信息学;生物信息;计算机辅助药物设计;测序分析;Python;R;机器学习;论文写作;网站制作;LOL;dota2。
目录
相关文章推荐
51好读  ›  专栏  ›  生信杂谈

Pymol加载对称单元

生信杂谈  · 公众号  ·  · 2018-05-05 20:13

正文

有时候有一些蛋白是对称单元组成的,如铁蛋白,但是晶体结构pdb坐标只有一部分,另外一部分写在信息里,我们如果想看全部或者需要全部的pdb坐标,可以使用如下方法实现:

Assembly

该功能需要PyMOL 1.8以上

设置对称单元,如下:

  1. set assembly, ""

  2. fetch 3bw1, asu, async=0

  3. set assembly, 1

  4. fetch 3bw1, assembly1, async=0

  5. set grid_mode

BiologicalUnit/Quat

该方法可以重构pdb的另外一边

python代码

  1. '''

  2. (c) 2010-2011 Thomas Holder, MPI for Developmental Biology

  3. Module for reading REMARK records from PDB files and in particular

  4. generate quaterny structure from REMARK 350.

  5. '''

  6. import sys, os

  7. from pymol import cmd, stored

  8. local_mirror_divided = '/mnt/bio/db/pdb.divided'

  9. def pdbremarks(filename):

  10.    '''

  11.    Read REMARK lines from PDB file. Return dictionary with remarkNum as key

  12.    and list of lines as value.

  13.    '''

  14.    remarks = dict()

  15.    if not isinstance(filename, basestring):

  16.        f = filename

  17.    elif filename[-3:] == '.gz':

  18.        import gzip

  19.        f = gzip.open(filename)

  20.    else:

  21.        f = open(filename)

  22.    for line in f:

  23.        recname = line[0:6]

  24.        if recname == 'REMARK':

  25.            num = int(line[7:10])

  26.            lstring = line[11:]

  27.            remarks.setdefault(num, []).append(lstring)

  28.    return remarks

  29. def quat350(rem350):

  30.    '''

  31.    Get transformation matrices for biomolecule 1 from REMARK 350.

  32.    '''

  33.    biomt = dict()

  34.    chains = tuple()

  35.    seenbiomolecule = False

  36.    for line in rem350:

  37.        if line.startswith('BIOMOLECULE:'):

  38.            if seenbiomolecule:

  39.                break

  40.            seenbiomolecule = True

  41.        elif line.startswith('APPLY THE FOLLOWING TO CHAINS:'):

  42.            chains = tuple(chain.strip() for chain in line[30:].split(','))

  43.        elif line.startswith('                   AND CHAINS:'):

  44.            chains += tuple(chain.strip() for chain in line [30:].split(','))

  45.        elif line.startswith('  BIOMT'):

  46.            row = int(line[7])

  47.            num = int(line[8:12])

  48.            vec = line[12:].split()

  49.            vec = map(float, vec)

  50.            biomt.setdefault(chains, dict()).setdefault(num, []).extend(vec)

  51.    return biomt

  52. def quat(name=None, filename=None, prefix=None , quiet=0):

  53.    '''

  54. DESCRIPTION

  55.    Read REMARK 350 from `filename` and create biological unit

  56.    (quaternary structure)

  57. USAGE

  58.    quat [name [, filename [, prefix]]]

  59. ARGUMENTS

  60.    name = string: name of object and basename of PDB file, if

  61.    filename is not given {default: first loaded object}

  62.    filename = string: file path {default: .pdb}

  63.    prefix = string: prefix for new objects {default: }

  64. EXAMPLE

  65.    fetch 1rmv, type=pdb

  66.    quat 1rmv

  67.    '''

  68.    quiet = int (quiet)

  69.    if name is None:

  70.        name = cmd.get_object_list()[0]

  71.    if prefix is None:

  72.        prefix = name

  73.    if filename is None:

  74.        candidates = [

  75.            '%s.pdb' % (name),

  76.            '%s/%s.pdb' % (cmd.get('fetch_path'), name),

  77.            '%s/%s/pdb%s.ent.gz' % (local_mirror_divided, name[1:3], name),

  78.        ]

  79.        for filename in candidates:

  80.            if os.path.exists(filename):

  81.                break

  82.        else:

  83.            print 'please provide filename'

  84.            return

  85.        if not quiet:

  86.            print 'loading from %s' % (filename)

  87.    remarks = pdbremarks(filename)

  88.    







请到「今天看啥」查看全文