Biskit.Fold

1 ## 2 ## Biskit, a toolkit for the manipulation of macromolecular structures 3 ## Copyright (C) 2004-2006 Raik Gruenberg & Johan Leckner 4 ## 5 ## This program is free software; you can redistribute it and/or 6 ## modify it under the terms of the GNU General Public License as 7 ## published by the Free Software Foundation; either version 2 of the 8 ## License, or any later version. 9 ## 10 ## This program is distributed in the hope that it will be useful, 11 ## but WITHOUT ANY WARRANTY; without even the implied warranty of 12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 13 ## General Public License for more details. 14 ## 15 ## You find a copy of the GNU General Public License in the file 16 ## license.txt along with this program; if not, write to the Free 17 ## Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. 18 ## 19 ## 20 ## $Revision: 2.11 $ 21 ## last $Date: 2007/04/06 10:16:07 $ 22 ## last $Author: graik $ 23 """ 24 Get binding energy from fold_X 25 """ 26 27 from Biskit import Executor, TemplateError 28 from Biskit import BiskitError 29 30 import Biskit.tools as T 31 32 import re, string, tempfile 33 34 import Biskit.molUtils as molUtils 35 36

37 -class Fold_XError( BiskitError ):

38 pass

39 40

41 -class Fold_X( Executor ):

42 """ 43 Run fold_X with given PDBModel 44 ============================== 45 Returns a dictionary with energy terms. 46 47 Example usage 48 ------------- 49 >>> x = Fold_X( model, verbose=1 ) 50 >>> result = x.run() 51 52 Energy terms (kcal/mol) 53 ----------------------- 54 This is the order of the energy terms in the output file. There 55 are two additional values that are unlabeled in the output 56 (version 2.5.2). The first is the entropy cost of making a complex. 57 This is zero if not 'AnalyseComplex' is run. The final value are 58 the number of residues. 59 60 We'll try to keep the energy labels of Fold-X version 1 of possible. 61 :: 62 ene - total energy 63 bb_hb - Backbone Hbond 64 sc_hb - Sidechain Hbond 65 vw - Van der Waals 66 el - Electrostatics 67 sol_p - Solvation Polar 68 sol_h - Solvation Hydrophobic 69 vwcl - Van de Waals clashes 70 sc - entropy side chain 71 mc - entropy main chain 72 sloop - sloop_entropy 73 mloop - mloop_entropy 74 cis - cis_bond 75 tcl - torsional clash 76 bbcl - backbone clash 77 dip - helix dipole 78 wtbr - water bridge 79 disu - disulfide 80 el_kon - electrostatic kon 81 p_cov - partial covalent bonds 82 ( - complex entropy ) 83 ( - number of residues ) 84 85 Reference 86 --------- 87 U{http://foldx.embl.de} for more info 88 89 @note: September 11 - 2006: 90 Tested with Fold-X 2.5.2 91 """ 92

93 - def __init__(self, model, **kw ):

94 """ 95 @param model: reference PDBModel 96 @type model: PDBModel 97 98 @param kw: additional key=value parameters for Executor: 99 @type kw: key=value pairs 100 :: 101 debug - 0|1, keep all temporary files (default: 0) 102 verbose - 0|1, print progress messages to log (log != STDOUT) 103 node - str, host for calculation (None->local) NOT TESTED 104 (default: None) 105 nice - int, nice level (default: 0) 106 log - Biskit.LogFile, program log (None->STOUT) (default: None) 107 """ 108 self.temp_pdb = tempfile.mktemp('_foldx_.pdb') 109 self.temp_command = tempfile.mktemp('_foldx_.command') 110 self.temp_option = tempfile.mktemp('_foldx_.option') 111 self.temp_result = tempfile.mktemp('_foldx_.result') 112 self.temp_runlog = tempfile.mktemp('_foldx_.log') 113 self.temp_errlog = tempfile.mktemp('_foldx_.err') 114 115 Executor.__init__( self, 'fold_X', args='-manual %s %s %s'\ 116 %(self.temp_pdb, self.temp_option, 117 self.temp_command), **kw ) 118 119 self.model = model.clone() 120 121 ## fold-X-allowed atoms for each res in standard order 122 self.aminoAcidDict = molUtils.aaAtoms 123 for k in self.aminoAcidDict: 124 if 'HN' not in self.aminoAcidDict[ k ]: 125 self.aminoAcidDict[ k ] += ['HN']

126 127

128 - def __prepareModel( self, model, f_pdb_out ):

129 """ 130 Prepare model for fold_X and write to file:: 131 - remove all hydrogens but H 132 - rename H to HN 133 - rename H1 to HN (remove H2 and H3) 134 - ( terminal oxygens should be OXT ) 135 - sort atoms according to self.aminoAcidDict 136 137 @param model: reference PDBModel 138 @type model: PDBModel 139 @param f_pdb_out: remaned pdb file name 140 @type f_pdb_out: str 141 142 @raise Fold_XError: if preparation fail 143 """ 144 145 def __cmpAtoms( a1, a2 ): 146 """ 147 compare atoms. 148 149 @param a1: atom name 150 @type a1: str 151 @param a2: atom name 152 @type a2: str 153 """ 154 res = a1['residue_name'] 155 target = self.aminoAcidDict[ res ] 156 try: 157 return cmp(target.index( a1['name'] ), 158 target.index( a2['name'] )) 159 except ValueError, why: 160 s = "Unknown atom for %s %i: %s or %s" % \ 161 (res, a1['residue_number'], a1['name'], a2['name'] ) 162 raise Fold_XError( s )

163 164 ## make a copy 165 # model = model.take( range(model.lenAtoms()) ) 166 167 ## mask for all heavy atoms, H and H3 168 heavy_mask = model.maskHeavy() 169 # HN_mask = model.mask( lambda a: a['name'] == 'H' ) 170 # H3_mask = model.mask( lambda a: a['name'] == 'H3' ) 171 172 ## rename H and H3 -> HN 173 # atm_dic = model.getAtoms() 174 # for i in N.nonzero( HN_mask + H3_mask ): 175 # atm_dic[i]['name'] = 'HN' 176 177 ## remove all none backbone hydrogens 178 keep_mask = heavy_mask #+ HN_mask + H3_mask 179 model = model.compress( keep_mask ) 180 181 ## consecutive residue numbering 182 model.renumberResidues() 183 184 ## sort atoms 185 model = model.sort( model.argsort( __cmpAtoms ) ) 186 187 model.writePdb( f_pdb_out )

188 189

190 - def prepare( self ):

191 """ 192 Write a Fold-X compatible pdb file to disc. 193 194 @note: Overrides Executor method. 195 """ 196 self.__prepareModel( self.model, self.temp_pdb ) 197 198 f_com = open( self.temp_command, 'w') 199 f_com.writelines(['<TITLE>FOLDX_commandfile;\n', 200 '<Stability>%s;\n'%self.temp_result]) 201 f_com.close() 202 203 f_opt = open( self.temp_option, 'w') 204 f_opt.writelines(['<TITLE> FOLDX_optionfile;\n', 205 '<logfile_name> %s;'%self.temp_runlog, 206 '<errorfile_name> %s;'%self.temp_errlog]) 207 f_opt.close()

208 209

210 - def cleanup( self ):

211 """ 212 Remove files created for and by the calculation. 213 """ 214 Executor.cleanup( self ) 215 216 if not self.debug: 217 T.tryRemove( self.temp_pdb ) 218 T.tryRemove( self.temp_command ) 219 T.tryRemove( self.temp_option ) 220 T.tryRemove( self.temp_result ) 221 T.tryRemove( self.temp_runlog ) 222 T.tryRemove( self.temp_errlog ) 223 ## Fold-X writes a file called "runlog.txt" 224 ## to local directory. Try to remove it. 225 T.tryRemove( 'runlog.txt' ) 226 ## and even though the error log is supposed 227 ## to be written to self.temp_errlog, I get a 228 ## 'errorfile.txt' in the local directory. Remove. 229 T.tryRemove( 'errorfile.txt' )

230 231

232 - def parse_foldx( self):

233 """ 234 Extract energies from output. 235 236 @return: dictionary with energy terms as keys 237 @rtype: dict 238 """ 239 energy = {} 240 241 ## extract energy terms from output 242 l= open(self.temp_result).readlines() 243 244 if len(l) == 0: 245 raise Fold_XError, 'ERROR: fold_X Segmentation fault' 246 247 if re.match('^%s'%self.temp_pdb, l[-1]): 248 249 ## get the energies 250 E = string.split( l[-1])[1:] 251 252 ## energy labels can't be parsed from file, so 253 ## we'll have to use this list (and lets hope 254 ## they don't change the order of the items). 255 keys = [ 'ene', 'bb_hb', 'sc_hb', 'vw', 256 'el', 'sol_p', 'sol_h', 'vwcl', 257 'sc', 'mc', 'sloop', 'mloop', 258 'cis', 'tcl', 'bbcl', 'dip', 259 'wtbr', 'disu', 'el_kon', 'p_cov' ] 260 261 ## add energies to dictionary 262 for i in range( len( keys ) ): 263 energy[ keys[i] ] = float( E[i] ) 264 else: 265 print 'No fold_X energy to extract' 266 energy = 0 267 268 return energy

269 270

271 - def isFailed( self ):

272 """ 273 @note: Overrides Executor method 274 """ 275 return not self.error is None

276 277

278 - def finish( self ):

279 """ 280 @note: Overrides Executor method 281 """ 282 Executor.finish( self ) 283 self.result = self.parse_foldx()

284 285 ############# 286 ## TESTING 287 ############# 288 import Biskit.test as BT 289

290 -class Test(BT.BiskitTest):

291 """Fold_X test""" 292 293 TAGS = [ BT.EXE ] 294

295 - def test_Fold_X( self):

296 """Fold_X test""" 297 from Biskit import PDBModel 298 299 ## Loading PDB... 300 self.m = PDBModel( T.testRoot() + '/rec/1A2P.pdb' ) 301 self.m = self.m.compress( self.m.maskProtein() ) 302 303 ## Starting fold_X 304 self.x = Fold_X( self.m, debug=self.DEBUG, 305 verbose=(self.VERBOSITY>2) ) 306 307 ## Running 308 self.r = self.x.run() 309 310 if self.local: 311 print "Result: ", self.r 312 313 for k, v in self.r.items(): 314 self.assertAlmostEqual( v, self.EXPECTED[k], 6 )

315 316 317 EXPECTED = {'el': -13.766400000000001, 'wtbr': -4.8059700000000003, 'ene': -18.475000000000001, 'mc': 160.28800000000001, 'sloop': 0.0, 'dip': 0.00177626, 'sol_p': 167.71100000000001, 'disu': 0.0, 'tcl': 0.72696700000000003, 'cis': 0.0, 'p_cov': 0.0, 'sol_h': -134.613, 'bb_hb': -87.362499999999997, 'sc_hb': -48.350000000000001, 'vw': -116.67100000000001, 'sc': 58.089300000000001, 'el_kon': 0.0, 'mloop': 0.0, 'vwcl': 0.27728599999999998, 'bbcl': 0.37019200000000002}

318 319 320 321 if __name__ == '__main__': 322 323 BT.localTest() 324

Source Code for Module Biskit.Fold_X