From 99ea3e555ae6eddd5fc0c0fa0cdbb64e75bb4dca Mon Sep 17 00:00:00 2001 From: DarioMarzella Date: Thu, 15 Jul 2021 20:28:24 +0200 Subject: [PATCH] Add docstrings to deeprank/features/Edesolv.py. Add Edesolv to docs/deeprank.features.rst --- deeprank/features/Edesolv.py | 134 +++++++++++++++++++---------------- docs/deeprank.features.rst | 10 +++ 2 files changed, 84 insertions(+), 60 deletions(-) diff --git a/deeprank/features/Edesolv.py b/deeprank/features/Edesolv.py index 4a429649..d854ce20 100755 --- a/deeprank/features/Edesolv.py +++ b/deeprank/features/Edesolv.py @@ -3,7 +3,6 @@ from pdb2sql import pdb2sql, interface from deeprank.features import FeatureClass from copy import deepcopy -from os import system from io import StringIO import numpy as np @@ -16,11 +15,20 @@ class Edesolv(FeatureClass): -#class Edesolv(): # init the class - def __init__(self,pdb_data): - #super().__init__('Atomic') + def __init__(self, pdb_data): + """Compute the energy of desolvation feature. + + Biopython 1.79 or later is required for this feature. + + Args: + pdb_data (list(byte) or str): pdb data + + Example: + >>> edesolv = Edesolv('1AK4.pdb') + >>> edesolv.get_feature() + """ self.pdb = pdb_data print('PDB DATA: ') print(pdb_data) @@ -29,26 +37,27 @@ def __init__(self,pdb_data): self.feature_data = {} self.feature_data_xyz = {} + # the feature extractor - def get_feature(self, verbose=False): - #def __compute_feature__(self, verbose=False): + """ Computes the actual Edesolv feature""" + esolcpx = 0.0 esolfree = 0.0 + # create a sql database pdb_db = pdb2sql(self.pdb) self.db = interface(pdb_db) + + #Create a temporary pdb-like string to calculate complex and free Edesolv temp_pdb = StringIO(('\n').join(pdb_db.sql2pdb())) temp_pdb.seek(0) - #Export PDB to parse it in BioPython - #pdb_db.exportpdb(fname = '_temp_.pdb') # Parse structure in Bio.PDB p = PDBParser(QUIET=1) - #struct = p.get_structure('temp', '_temp_.pdb')#(self.pdb.split('/')[-1].split('.')[0], self.pdb) - #system('rm _temp_.pdb') struct = p.get_structure('temp', temp_pdb) + # Identify chains in the Bio.PDB object chains = [x.id for x in struct[0] if x.id != ' '] chainA = chains[0] chainB = chains[1] @@ -68,6 +77,7 @@ def get_feature(self, verbose=False): # get the contact atoms indA,indB = list(self.db.get_contact_atoms(chain1=chainA, chain2=chainB).values()) contact = indA + indB + # extract the atom keys and xyz of the contact CA atoms keys = self.db.get('serial,chainID,resName,resSeq,name',rowID=contact) xyz = self.db.get('x,y,z',rowID=contact) @@ -85,7 +95,6 @@ def get_feature(self, verbose=False): # create the dictionary of human readable and xyz-val data - hread, xyzval = {},{} self.edesolv_data = {} self.edesolv_data_xyz = {} @@ -93,98 +102,85 @@ def get_feature(self, verbose=False): brk_flag = False disulf_bond = False atom = Atom(key) - #if verbose: - print('key: ', key) - #print('Values for atom n. ', atom.serial) - #print('PDB2SQL values:') - #print('Name: ', atom.name, 'Resn: ', atom.resn, 'Pos: ', atom.position) # atom.position specifies the position in the residue, if Back Bone or Side Chain + # Calculate complex sasa try: atom.cpx_sasa = struct[0][atom.chainID][atom.resid][atom.name].sasa except KeyError: #Handling alternative residues error print('Alternative residue found at:', key) brk_flag = True - #raise Exception('KeyError') pass + + # Calculate unbound (translated) sasa try: atom.free_sasa = free_struct[0][atom.chainID][atom.resid][atom.name].sasa except KeyError: #Handling alternative residues error print('Alternative residue found at:', key) brk_flag = True - #raise Exception('KeyError') pass if not brk_flag: - #if verbose: - #print('Complex SASA: ', atom.cpx_sasa) - #print('Free SASA: ', atom.free_sasa) - + # check if the atom is a sulfide involved in a disulfide bond for cys in disulfides_cys: if cys == (atom.chainID, atom.resid): disulf_bond = True - assign_solv_param(atom, disulf_bond) - #if verbose: - #print('Solv: ', atom.solv) + # Assign solvation parameters + Assign_solv_param(atom, disulf_bond) + #Calculate atom complex edesolv, free edesolv and final edesolv (cpx - free) atom.esolcpx = atom.cpx_sasa * atom.solv atom.esolfree = atom.free_sasa * atom.solv atom.edesolv = atom.esolcpx - atom.esolfree esolcpx += atom.esolcpx esolfree += atom.esolfree - #if verbose: - #print('atom_esolcpx: ', atom.esolcpx) - #print('atom_esolfree: ', atom.esolfree) - #print('atom_Edesolv: ', atom.edesolv) - #print('\n') - #raise Exception('OK.') - - # human readable - # { (chainID,resName,resSeq,name) : [val] } - hread[tuple(key)] = [atom.edesolv] #Put Edesolv here! # xyz-val # { (0|1,x,y,z) : [val] } chain = [{chainA:0,chainB:1}[key[1]]] - #k = tuple(chain + xyz) k = tuple(chain + coords) - #print(k, type(k)) - #xyzval[k] = [atom.edesolv] + #Human readable self.edesolv_data[(key[1], key[3], key[2], key[4])] = [atom.edesolv] + #XYZ data self.edesolv_data_xyz[k] = [atom.edesolv] - - #self.feature_data['Edesolv'] = hread - #self.feature_data_xyz['Edesolv'] = xyzval - print('edesolv_data: ', self.edesolv_data) - print('edesolv_data_xyz: ', self.edesolv_data_xyz) self.feature_data['Edesolv'] = self.edesolv_data self.feature_data_xyz['Edesolv'] = self.edesolv_data_xyz - #Edesolv = esolcpx - esolfree - #return Edesolv - ################################################################## ################################################################## -def __compute_feature__(pdb_data, featgrp, featgrp_raw, chainA, chainB): +def __compute_feature__(pdb_data, featgrp, featgrp_raw): + """Main function called in deeprank for the feature calculations. - edesolv_feat = Edesolv(pdb_data) - edesolv_feat.get_feature() + Args: + pdb_data (list(bytes)): pdb information + featgrp (str): name of the group where to save xyz-val data + featgrp_raw (str): name of the group where to save human readable data + """ + + edesolv_feat = Edesolv(pdb_data) + edesolv_feat.get_feature() - # export in the hdf5 file - edesolv_feat.export_dataxyz_hdf5(featgrp) - edesolv_feat.export_data_hdf5(featgrp_raw) + # export in the hdf5 file + edesolv_feat.export_dataxyz_hdf5(featgrp) + edesolv_feat.export_data_hdf5(featgrp_raw) - # close - edesolv_feat.db._close() + # close + edesolv_feat.db._close() class Atom(): - def __init__(self,atom_specs): # [498, 'A', 'HIS', 54, 'CB'] + def __init__(self,atom_specs): + """Atom class necessary to calculate atomic-level surface area, + solvation parameters and desolvation energy. + + Args: + atom_specs(""" + serial = atom_specs[0] chainID = atom_specs[1] resn = atom_specs[2] @@ -202,7 +198,18 @@ def __init__(self,atom_specs): # [498, 'A', 'HIS', 54, 'CB'] self.position = 'SC' -def assign_solv_param(atom, disulf_bond=False): +def Assign_solv_param(atom, disulf_bond=False): + """Assigns solvation parameter to an Atom object + + Args: + atom (Atom): Edesolv.Atom object. + disulf_bond (bool, optional): Disulfide bond flag. If true, the atom is a sulfide involved in a disulfie bond. Defaults to False. + + Returns: + None. + + """ + arofac = 6.26 alifac = 1.27 polfac = 2.30 @@ -236,7 +243,7 @@ def assign_solv_param(atom, disulf_bond=False): atom.solv = 0.0000 ''' - # Corarse-Grained model + # Corarse-Grained model from HADDOCK code for eventual future implementation if atom.position == "BB" and atom.resn=="ALA": atom.solv = -0.0107 elif atom.position == "BB" and atom.resn=="GLY": @@ -367,8 +374,14 @@ def get_disulfide_bonds(pdb_data): return disulfide_cys def get_dihedral(p): - """Praxeolitic formula - 1 sqrt, 1 cross product""" + """Get dihedral angle from four points (atoms) 3D coordinates. + + Args: + p (list): list of four points 3D coordinates. Each element must be a list of floats. + + Returns: + dihedral (float): Dihedral angle. + """ p0 = p[0] p1 = p[1] p2 = p[2] @@ -394,4 +407,5 @@ def get_dihedral(p): # v and w may not be normalized but that's fine since tan is y/x x = np.dot(v, w) y = np.dot(np.cross(b1, v), w) - return np.degrees(np.arctan2(y, x)) + dihedral = np.degrees(np.arctan2(y, x)) + return dihedral diff --git a/docs/deeprank.features.rst b/docs/deeprank.features.rst index 28f53f3d..ef67db3e 100644 --- a/docs/deeprank.features.rst +++ b/docs/deeprank.features.rst @@ -35,6 +35,16 @@ Buried Surface Area .. autofunction:: deeprank.features.BSA.__compute_feature__ +Energy of desolvation +------------------- + +.. automodule:: deeprank.features.Edesolv + :members: + :undoc-members: + :private-members: + +.. autofunction:: deeprank.features.Edesolv.__compute_feature__ + FullPSSM --------