1
2
3
4
5
6 from __future__ import generators
7
8 import os
9 import tempfile
10 from Bio.PDB import *
11 from PDBExceptions import PDBException
12 from AbstractPropertyMap import AbstractResiduePropertyMap
13 import re
14
15
16 __doc__="""
17 Use the DSSP program to calculate secondary structure and accessibility.
18 You need to have a working version of DSSP (and a license, free for
19 academic use) in order to use this. For DSSP, see U{http://www.cmbi.kun.nl/gv/dssp/}.
20
21 The DSSP codes for secondary structure used here are:
22
23 - H Alpha helix (4-12)
24 - B Isolated beta-bridge residue
25 - E Strand
26 - G 3-10 helix
27 - I pi helix
28 - T Turn
29 - S Bend
30 - - None
31 """
32
33
34 _dssp_cys=re.compile('[a-z]')
35
36
37
38
39 MAX_ACC={}
40 MAX_ACC["ALA"]=106.0
41 MAX_ACC["CYS"]=135.0
42 MAX_ACC["ASP"]=163.0
43 MAX_ACC["GLU"]=194.0
44 MAX_ACC["PHE"]=197.0
45 MAX_ACC["GLY"]=84.0
46 MAX_ACC["HIS"]=184.0
47 MAX_ACC["ILE"]=169.0
48 MAX_ACC["LYS"]=205.0
49 MAX_ACC["LEU"]=164.0
50 MAX_ACC["MET"]=188.0
51 MAX_ACC["ASN"]=157.0
52 MAX_ACC["PRO"]=136.0
53 MAX_ACC["GLN"]=198.0
54 MAX_ACC["ARG"]=248.0
55 MAX_ACC["SER"]=130.0
56 MAX_ACC["THR"]=142.0
57 MAX_ACC["VAL"]=142.0
58 MAX_ACC["TRP"]=227.0
59 MAX_ACC["TYR"]=222.0
60
61
63 """
64 Secondary structure symbol to index.
65 H=0
66 E=1
67 C=2
68 """
69 if ss=='H':
70 return 0
71 if ss=='E':
72 return 1
73 if ss=='C':
74 return 2
75 assert(0)
76
78 """
79 Create a DSSP dictionary from a PDB file.
80
81 Example:
82 >>> dssp_dict=dssp_dict_from_pdb_file("1fat.pdb")
83 >>> aa, ss, acc=dssp_dict[('A', 1)]
84
85 @param in_file: pdb file
86 @type in_file: string
87
88 @param DSSP: DSSP executable (argument to os.system)
89 @type DSSP: string
90
91 @return: a dictionary that maps (chainid, resid) to
92 amino acid type, secondary structure code and
93 accessibility.
94 @rtype: {}
95 """
96 out_file=tempfile.mktemp()
97 os.system(DSSP+" %s > %s" % (in_file, out_file))
98 dict, keys=make_dssp_dict(out_file)
99
100
101 return dict, keys
102
104 """
105 Return a DSSP dictionary that maps (chainid, resid) to
106 aa, ss and accessibility, from a DSSP file.
107
108 @param filename: the DSSP output file
109 @type filename: string
110 """
111 dssp={}
112 fp=open(filename, "r")
113 start=0
114 keys=[]
115 for l in fp.readlines():
116 sl=l.split()
117 if sl[1]=="RESIDUE":
118
119 start=1
120 continue
121 if not start:
122 continue
123 if l[9]==" ":
124
125 continue
126 resseq=int(l[5:10])
127 icode=l[10]
128 chainid=l[11]
129 aa=l[13]
130 ss=l[16]
131 if ss==" ":
132 ss="-"
133 acc=int(l[34:38])
134 res_id=(" ", resseq, icode)
135 dssp[(chainid, res_id)]=(aa, ss, acc)
136 keys.append((chainid, res_id))
137 fp.close()
138 return dssp, keys
139
140
141 -class DSSP(AbstractResiduePropertyMap):
142 """
143 Run DSSP on a pdb file, and provide a handle to the
144 DSSP secondary structure and accessibility.
145
146 Note that DSSP can only handle one model.
147
148 Example:
149 >>> p=PDBParser()
150 >>> structure=parser.get_structure("1fat.pdb")
151 >>> model=structure[0]
152 >>> dssp=DSSP(model, "1fat.pdb")
153 >>> # print dssp data for a residue
154 >>> secondary_structure, accessibility=dssp[(chain_id, res_id)]
155 """
156 - def __init__(self, model, pdb_file, dssp="dssp"):
157 """
158 @param model: the first model of the structure
159 @type model: L{Model}
160
161 @param pdb_file: a PDB file
162 @type pdb_file: string
163
164 @param dssp: the dssp executable (ie. the argument to os.system)
165 @type dssp: string
166 """
167
168 dssp_dict, dssp_keys=dssp_dict_from_pdb_file(pdb_file, dssp)
169 dssp_map={}
170 dssp_list=[]
171
172
173
174 for key in dssp_keys:
175 chain_id, res_id=key
176 chain=model[chain_id]
177 res=chain[res_id]
178 aa, ss, acc=dssp_dict[key]
179 res.xtra["SS_DSSP"]=ss
180 res.xtra["EXP_DSSP_ASA"]=acc
181
182 resname=res.get_resname()
183 rel_acc=acc/MAX_ACC[resname]
184 if rel_acc>1.0:
185 rel_acc=1.0
186 res.xtra["EXP_DSSP_RASA"]=rel_acc
187
188
189 resname=to_one_letter_code[resname]
190 if resname=="C":
191
192
193 if _dssp_cys.match(aa):
194 aa='C'
195 if not (resname==aa):
196 raise PDBException, "Structure/DSSP mismatch at "+str(res)
197 dssp_map[key]=((res, ss, acc, rel_acc))
198 dssp_list.append((res, ss, acc, rel_acc))
199 AbstractResiduePropertyMap.__init__(self, dssp_map, dssp_keys, dssp_list)
200
201
202 if __name__=="__main__":
203
204 import sys
205
206 p=PDBParser()
207 s=p.get_structure('X', sys.argv[1])
208
209 model=s[0]
210
211 d=DSSP(model, sys.argv[1])
212
213 for r in d:
214 print r
215
216 print d.keys()
217
218 print len(d)
219
220 print d.has_key(('A', 1))
221
222 print d[('A', 1)]
223
224 print s[0]['A'][1].xtra
225