PDB 形式のファイルの ATOM レコードまたは SEQRES レコードにはアミノ酸配列の情報が記載されている。このページでは、PDB ファイルから FASTA ファイルを作成する Python スクリプトの例を示す。

ディクショナリを利用した方法

ファイルを開いて ATOM 行を探し、C_α を見つけたら、そのアミノ酸の 3 文字コードを取得し、1 文字に変換する。最後に、それらを出力する。

import sys
argvs = sys.argv

pdb_file = '1ALK.pdb'
fa = {}

amino_code = {
	'ALA':'A', 'ARG':'R', 'ASN':'N', 'ASP':'D',
	'CYS':'C', 'GLN':'Q', 'GLU':'E', 'GLY':'G',
	'ILE':'I', 'LEU':'L', 'LYS':'K', 'MET':'M',
	'PHE':'F', 'PRO':'P', 'SER':'S', 'THR':'T',
	'TRP':'W', 'TYR':'Y', 'VAL':'V', 'HIS':'H',
	'ASX':'B', 'GLX':'Z', 'UNK':'K'
}

with open(pdb_file) as fh:
	for buff in fh:
		if (buff[0:4] != 'ATOM'):
			continue
		
		chain_name = buff[21:22]
		res_number = int(buff[22:26])
		amino_acid = buff[17:20]
		
		if not (chain_name in fa):
			fa[chain_name] = []

		aa = 'X'
		if (amino_acid in amino_code):
			aa = amino_code[amino_acid]
		
		if (len(fa[chain_name]) != res_number):
			fa[chain_name] += ['X'] * (res_number - len(fa[chain_name]))
		
		fa[chain_name][res_number - 1] = aa

for k, v in sorted(fa.items()):
	print ('>' + k)
	print (''.join(v))

PDB から 1ALK.pdb をダウンロードし、上のスクリプトを実行させると、FASTA 形式が画面に出力される。

BioPython を利用した方法

BioPython にある PDBParser などのモジュールを利用する場合、以下のように FASTA 形式のファイルを作成できる。

from Bio.PDB import *

amino_code = {
    'ALA':'A', 'ARG':'R', 'ASN':'N', 'ASP':'D',
    'CYS':'C', 'GLN':'Q', 'GLU':'E', 'GLY':'G',
    'ILE':'I', 'LEU':'L', 'LYS':'K', 'MET':'M',
    'PHE':'F', 'PRO':'P', 'SER':'S', 'THR':'T',
    'TRP':'W', 'TYR':'Y', 'VAL':'V', 'HIS':'H',
    'ASX':'B', 'GLX':'Z', 'UNK':'K'
}

p = PDBParser()
structure = p.get_structure('X', '1ALK.pdb')

for model in structure:
    for chain in model:
        chain_id = chain.get_id()
        chain_seq = ''
        for residue in chain:
            aaa = residue.get_resname()
            a = ''
            if (aaa in amino_code):
                a = amino_code[aaa]
            chain_seq += a
        print('>' + chain_id)
        print(chain_seq)