<jupyter_start><jupyter_text>MEGADOCK POST-PROCESSING Utility Class and functions Table of Contents  * [The megadock format](chapter-1)  * [Chapter 2](chapter-2)  * [Chapter 3](chapter-3)   The megadock format A typical megadock results file is derived from the [ZDOCK format](http://zdock.umassmed.edu/zdock_output_format.php).First line feature the number grid cells and the size of the cubic cell. Line 3 indicates the initial position of the **receptor** barycenter within the original system of coordinates. Line 4 indicates similar information for the **ligand**.From line 12 to 16, docking poses information are provided. In order to reconstruct a docking pose both **receptor** and **ligand** molecule barycenters have to be set to the origin.The first field is the Euler angles triplet to apply rotation to the **ligand** molecule to reconstruct its orientation in the pose. The second triplet is the translation vector in mesh units to move the **ligand** molecule to its position in  the pose. 1 126     1.20 2         0.000000        0.000000        0.000000 3 1a2k_u1.pdb     27.372499       8.099500        82.649002 4 1a2k_u2.pdb     59.376999       0.796000        75.828003 5  6 Elec score ratio:   1.00 7 ACE  score ratio:   1.00 8  9 *** docking results ***1011 Rank   Angle               Trans         total       rPSC     rPSC(+)  rPSC(-)  ELEC       RecACE     LCARMSD NearNat12     1: (-2.72, 1.87, 0.62) ( 19, 15, 12)     4376.78     2977     7522    -4545     315.14    1084.64  56.351       013     2: (-2.72, 1.74, 0.66) ( 19, 14, 13)     4288.76     3066     7566    -4500     140.73    1082.03  56.310       014     3: (-2.93, 1.98,-3.13) ( 27,  8,  7)     4288.50     3394     7579    -4185     -27.59     922.10  64.465       015     4: ( 1.68, 1.19, 2.71) (121,104,106)     4270.28     3384     5994    -2610      94.39     791.89  52.556       016     5: ( 2.41, 1.63, 1.84) (  7, 96,117)     4269.77     3563     4688    -1125     337.08     369.69  61.792       0So for each pose, the ccmap libray has to: 1. move ligand and receptor to center2. rotate ligand3. translate ligand4. compute interaction matrix  The megadock container     1* Each pose keeps a dictionary of its amino-acids contacts 10000 poses are treated in 4mn, with a memory cost ~ 9286M    2* Releasing the GIL w/in C extension does reduce total computation cost to 3mn40s. Bottleneck is single structure computation in C<jupyter_code>import sys, os, re, pickle, json
#sys.path.append("/Users/jmartin/dev/pyproteinsExt/src")
#sys.path.append("/Users/jmartin/dev/pyproteins/src")

sys.path.append("/Users/guillaumelaunay/work/DVL/pyproteinsExt/src")
sys.path.append("/Users/guillaumelaunay/work/DVL/pyproteins/src")

%load_ext autoreload
import pyproteinsExt.structure.coordinates as PDB
import pyproteinsExt.structure.operations as PDBop
import ccmap

from multiprocessing import Pool

reL1 = ur'^([\d]+)[\s]+([\d\.]+)[\s]*$'
reL2 = ur'^[\s]*([\.\d-]+)[\s]+([\d\.-]+)[\s]+([\.\d-]+)[\s]*$'
reL3 = ur'^[\s]*([\S]+)[\s]+([\.\d-]+)[\s]+([\d\.-]+)[\s]+([\.\d-]+)[\s]*$'
reZPOSE = ur'^[\s]*[\d]+:[\s]+\([\s]*([\.\d]+),[\s]*([\.\d]+),[\s]*([\.\d]+)\)[\s]+\([\s]*([\.\d]+),[\s]*([\.\d]+),[\s]*([\.\d]+)\)$'

reZPOSE = ur'^[\s]*([\d-]+):[\s]+\([\s]*([\.\d-]+),[\s]*([\.\d-]+),[\s]*([\.\d-]+)\)[\s]+\([\s]*([\d-]+),[\s]*([\d-]+),[\s]*([\d-]+)\)' 

parserPDB = PDB.Parser()


def unpickle(fileName):
    with open(fileName, 'r') as f:
        megadockObj = pickle.load(f)
    return megadockObj

def parse(fileName, maxPose = 0):
    megadockObj = Megadock()#ncpu and pSize arguments can be provided
    dockBool = False
    with open(fileName, 'r') as f:
        for line in f:
            if '*** docking results ***' in line:
                dockBool = True
            
            m = re.match(reL1, line)
            if m :
                megadockObj.nCells = int(m.groups()[0])
                megadockObj.step = float(m.groups()[1])
                continue
            m = re.match(reL2, line)
            if m :
                megadockObj.eulerREC = ( float(m.groups()[0]), float(m.groups()[1]),float(m.groups()[2]) )
                continue
            m = re.match(reL3, line)
            if m :
                if not megadockObj.fileREC:
                    megadockObj.fileREC = m.groups()[0]
                    megadockObj.baryREC = (float(m.groups()[1]), float(m.groups()[2]), float(m.groups()[3]))
                    continue
                megadockObj.fileLIG = m.groups()[0]
                megadockObj.baryLIG = (float(m.groups()[1]), float(m.groups()[2]), float(m.groups()[3]))
                continue
            m = re.match(reZPOSE, line)
            if m and dockBool:
                euler = (float(m.groups()[1]), float(m.groups()[2]), float(m.groups()[3]))
                tr = [int(m.groups()[4]), int(m.groups()[5]), int(m.groups()[6])]
                tr=tuple([ t - megadockObj.nCells if t > megadockObj.nCells / 2 else t for t in tr ])
                
                megadockObj.push( int(m.groups()[0]), euler, tr)
                if len(megadockObj.zList) == maxPose:
                    return megadockObj
        
        return megadockObj


# Helper function for multiprocessing
def mpCcmap(datum):
    zPack, dist = datum
    for z in zPack:
            z.ccmap(dist=dist)
    return zPack
# Megadock docking run results class
# Takes a flat result file or a serialized, once ccmap have been computed
class Megadock(object):
    
    def dump(self):
        return pickle.dumps(self)
    
    
    def __init__(self, **kwargs):
        self.step = None
        self.nCells = None
        self.eulerREC = None
        self.fileREC = None
        self.fileLIG = None
        self.baryREC = None
        self.baryLIG = None
        self.zList = []
        self.pdbObjLigand = None
        self.pdbObjReceptor = None
        

    #Post processing / statistics operations on set of poses
    def contactSummary(self):
        pass # access individual pose through self.zList
    
    
    def ccmap(self, **kwargs):
        start = kwargs['start'] if 'start' in kwargs else 0
        stop  = kwargs['stop']  if 'stop' in kwargs else len(self.zList)
        dist  = kwargs['dist']  if 'dist' in kwargs else 4.5
        
        self.threadPacketSize =  kwargs['pSize'] if 'pSize' in kwargs else 200
        self.ncpu = kwargs['ncpu'] if 'ncpu' in kwargs else None
        
        elem = []
        indices= slice(start, stop, self.threadPacketSize).indices(len(self.zList))
        for i in range(*indices):
            j = i + self.threadPacketSize if i + self.threadPacketSize <= stop else stop
            elem.append( (self.zList[slice(i, j)], dist) )
        
        print("Created " + str(len(elem)) + " data packets (" +  str(self.threadPacketSize) + " zObjects each) for process pool")
        pool = Pool(self.ncpu) if self.ncpu else Pool() #note the default will use the optimal number of workers
        mpData = pool.map(mpCcmap, elem)
        pool.close()
        self._unpackProcessPool(mpData)
    
    def _unpackProcessPool(self, poolResults):
         #print self.zList[0]._ccmap
        #print hope[0][0]._ccmap
        print ('unpacking')
        i=0
        for pack in poolResults:
            for zPoseClone in pack: 
                if zPoseClone.id != self.zList[i].id:
                    raise ValueError('mismatching id ' + str(self.zList[i].id) + ' '  + str(zPoseClone.id) ) 
                self.zList[i]._ccmap = zPoseClone._ccmap
                self.zList[i].dictorizedReceptor = zPoseClone.dictorizedReceptor
                self.zList[i].dictorizedLigand = zPoseClone.dictorizedLigand
                i += 1
        
    
    def push(self, *args):
        self.zList.append( zPose(self, *args) )
            
    def setReceptor(self, pdbFile): # Optional chain arguments ?
        self.pdbObjReceptor = parserPDB.load(file=pdbFile)
    def setLigand(self, pdbFile):
        self.pdbObjLigand = parserPDB.load(file=pdbFile)

    def __str__(self):
        return str({ 'step' : self.step, 'nCells' : self.nCells , 'EulerREC' : self.eulerREC,
                       'fileREC' : self.fileREC, 'baryREC' : self.baryREC,
                       'fileLIG' : self.fileLIG, 'baryLIG' : self.baryLIG,
                        'zList'  : self.zList
                   })
    
    def __getitem__(self, key):      
        return self.zList[key]
    
    def __iter__(self):
        for zPose in self.zList:
            yield zPose
    
class zPose(object):
    def __init__(self, belongsTo, id, euler, tr):
        self.id = id
        self.euler = euler
        self.translate = tuple([ (-1) * t * belongsTo.step for t in tr ])
        self._ccmap = None
        self.belongsTo = belongsTo
        self.ligOffset = tuple( [ (-1) * o for o in belongsTo.baryLIG] ) 
        self.recOffset = tuple( [ (-1) * o for o in belongsTo.baryREC] ) 
        self.dictorizedReceptor = None
        self.dictorizedLigand = None
        
    def __str__(self):
        return str(self.id) + ') ' + str(self.euler) + ' ' + str(self.translate)
    def __repr__(self):
        return str(self.id) + ') ' + str(self.euler) + ' ' + str(self.translate)
    
    def ccmap(self, dist=4.5):
        #print self._ccmap
        if not self._ccmap:
            pdbObjRec = self.belongsTo.pdbObjReceptor
            pdbObjLig = self.belongsTo.pdbObjLigand
            #tmp_ccmap = ccmap.zmap(dist, self.euler)
            self.dictorizedReceptor = pdbObjRec.atomDictorize
            self.dictorizedLigand = pdbObjLig.atomDictorize
            #print dist
            ccmapAsString = ccmap.zmap( (self.dictorizedReceptor, self.dictorizedLigand), dist, 
                                   self.euler, self.translate, self.recOffset,  self.ligOffset)
            self._ccmap = json.loads(ccmapAsString)
        return self._ccmap

    def dump(self):
        pdbObjRec = self.belongsTo.pdbObjReceptor.clone()
        pdbObjLig = self.belongsTo.pdbObjLigand.clone()

        if self.dictorizedReceptor:
            pdbObjRec.setCoordinateFromDictorize(self.dictorizedReceptor)
        if self.dictorizedLigand:
            pdbObjLig.setCoordinateFromDictorize(self.dictorizedLigand)

        return str(pdbObjRec) + str(pdbObjLig)
                
                
#zDockRes = parse('/Users/jmartin/CHOKO/2017_10_16_Notebook_ccmap_megadock/1a2k_u.detail')
#zDockRes.setReceptor('/Users/jmartin/CHOKO/2017_10_16_Notebook_ccmap_megadock/1a2k_u1_rename.pdb')
#zDockRes.setLigand('/Users/jmartin/CHOKO/2017_10_16_Notebook_ccmap_megadock/1a2k_u2_rename.pdb')

zDockRes = parse('/Users/guillaumelaunay/work/DVL/pythonExtend/ccmap/data/1a2k_megadock.all')
zDockRes.setReceptor('/Users/guillaumelaunay/work/DVL/pythonExtend/ccmap/data/1a2k_u1.pdb')
zDockRes.setLigand('/Users/guillaumelaunay/work/DVL/pythonExtend/ccmap/data/1a2k_u2.pdb')


zDockRes.ccmap(stop=10000,dist=5, pSize=200)<jupyter_output>Created 50 data packets (200 zObjects each) for process pool
unpacking<jupyter_text>printing ccmap as lists of contacts<jupyter_code>with open('start.temp',"w") as f:
    print 'start' 
    
zDockRes.ccmap(stop=10000,dist=5)

class Item(object):
    def __init__(self, datum):
        self.chainID = datum['chainID']
        self.resID = datum['resID']
        
    def __hash__(self):
        return hash(self.chainID + self.resID)

    def __str__(self):
        return self.resID + self.chainID

    def __repr__(self):
        return self.__str__()
    
for i in range(10000):
    mapRegistry = zDockRes[i].ccmap()
    with open('complex.'+str(i+1)+'_contact.txt', "w") as f:
        for cclist in mapRegistry['data']:
            rootRes = Item(cclist['root'])
            for p in cclist['partners']:
                partnerRes = Item(p)
                f.write(str(rootRes)+','+str(partnerRes)+'\n')
with open('finish.temp',"w") as f:
    print 'finish'<jupyter_output>start<jupyter_text>Printing PDB coordinates of specific poses<jupyter_code>zDockRes.ccmap(stop=100,dist=5, pSize=25)
#Dumping PDB files
# pdbCoordinateString = zDockRes[1].dump()
# print pdbCoordinateString
for pose in zDockRes[:100]:
    print pose
    pdbFile = '/Users/guillaumelaunay/work/DVL/pythonExtend/ccmap/data/1a2k_pose_' + str(pose.id) + '.pdb'
    pdbString = pose.dump()
    with open(pdbFile, "w") as f:
        f.write(pdbString)<jupyter_output>Created 1 data packets (200 zObjects each) for process pool
unpacking
1) (-2.72, 1.87, 0.62) (-22.8, -18.0, -14.399999999999999)
2) (-2.72, 1.74, 0.66) (-22.8, -16.8, -15.6)
3) (-2.93, 1.98, -3.13) (-32.4, -9.6, -8.4)
4) (1.68, 1.19, 2.71) (6.0, 26.4, 24.0)
5) (2.41, 1.63, 1.84) (-8.4, 36.0, 10.799999999999999)
6) (1.57, 1.28, 2.81) (7.199999999999999, 24.0, 25.2)
7) (0.1, 1.82, -1.56) (1.2, 12.0, 31.2)
8) (1.99, 2.53, 2.74) (22.8, -15.6, -18.0)
9) (-1.68, 1.31, 1.63) (31.2, 6.0, 3.5999999999999996)
10) (2.41, 1.51, 1.88) (-9.6, 34.8, 12.0)
11) (-2.41, 2.95, 1.13) (2.4, 12.0, 30.0)
12) (3.04, 0.89, -1.87) (-9.6, -25.2, -19.2)
13) (-0.94, 0.05, -2.41) (2.4, 9.6, 31.2)
14) (0.0, 0.05, -2.41) (10.799999999999999, 16.8, 30.0)
15) (1.47, 1.28, 2.81) (9.6, 25.2, 25.2)
16) (2.3, 1.63, 1.84) (-7.199999999999999, 36.0, 10.799999999999999)
17) (-1.88, 2.04, -2.65) (-33.6, -12.0, -7.199999999999999)
18) (0.31, 1.3, 3.04) (34.8, -4.8, 0.0)
19) (-2.09, 2.18, 2.6) (2.4, -33.6, -18.0)
20) (-0.1, 0.75[...]<jupyter_text>Serializing // Deserializing<jupyter_code>zDockRes.ccmap(stop=10)
serializedFile = '/Users/jmartin/CHOKO/2017_10_16_Notebook_ccmap_megadock/1a2k_megadock_10.pkl'
with open(serializedFile, 'w') as f:
    f.write(zDockRes.dump())
zDockResUpkl = unpickle(serializedFile)
zDockResUpkl.pdbObjReceptor.atomDictorize
print(zDockResUpkl[0]._ccmap)<jupyter_output>None
{'fraction_non_alphanumeric': 0.1496653084982538, 'fraction_numerical': 0.0725407450523865, 'mean_word_length': 2.3963429701013097, 'pattern_counts': {'":': 0, '<': 14, '<?xml version=': 0, '>': 14, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 8, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '4', 'dir': 'jupyter-structured-clean-dedup', 'id': '570448', 'license': 'MIT', 'n_tokens_mistral': 5686, 'n_tokens_neox': 5062, 'n_words': 1118, 'path': 'notebook/--= MEGADOCK POST PROCESSING =--.ipynb', 'repo_name': 'MMSB-MOBI/ccmap', 'size': '13744'}