Package rdkit :: Package Chem :: Package Pharm2D :: Module Generate
[hide private]
[frames] | no frames]

Source Code for Module rdkit.Chem.Pharm2D.Generate

  1  # $Id$ 
  2  # 
  3  # Copyright (C) 2002-2008 greg Landrum and Rational Discovery LLC 
  4  # 
  5  #   @@ All Rights Reserved @@ 
  6  #  This file is part of the RDKit. 
  7  #  The contents are covered by the terms of the BSD license 
  8  #  which is included in the file license.txt, found at the root 
  9  #  of the RDKit source tree. 
 10  # 
 11  """ generation of 2D pharmacophores 
 12   
 13  **Notes** 
 14   
 15    - The terminology for this gets a bit rocky, so here's a glossary of what 
 16      terms used here mean: 
 17   
 18        1) *N-point pharmacophore* a combination of N features along with 
 19           distances betwen them. 
 20   
 21        2) *N-point proto-pharmacophore*: a combination of N feature 
 22           definitions without distances.  Each N-point 
 23           proto-pharmacophore defines a manifold of potential N-point 
 24           pharmacophores. 
 25   
 26        3) *N-point scaffold*: a collection of the distances defining 
 27           an N-point pharmacophore without feature identities. 
 28   
 29    See Docs/Chem/Pharm2D.triangles.jpg for an illustration of the way 
 30    pharmacophores are broken into triangles and labelled. 
 31   
 32    See Docs/Chem/Pharm2D.signatures.jpg for an illustration of bit 
 33    numbering 
 34   
 35  """ 
 36  from __future__ import print_function 
 37  from rdkit.Chem.Pharm2D import Utils,SigFactory 
 38  from rdkit.RDLogger import logger 
 39  logger = logger() 
 40   
 41  _verbose = 0 
 42   
43 -def _ShortestPathsMatch(match,featureSet,sig,dMat,sigFactory):
44 """ Internal use only 45 46 """ 47 if _verbose: 48 print('match:',match) 49 nPts = len(match) 50 distsToCheck = Utils.nPointDistDict[nPts] 51 nDists = len(distsToCheck) 52 dist = [0]*nDists 53 bins = sigFactory.GetBins() 54 minD,maxD = bins[0][0],bins[-1][1] 55 56 for i in range(nDists): 57 pt0,pt1 = distsToCheck[i] 58 minSeen=maxD 59 for idx1 in match[pt0]: 60 for idx2 in match[pt1]: 61 minSeen=min(minSeen, dMat[idx1,idx2]) 62 if minSeen==0 or minSeen<minD: return 63 # FIX: this won't be an int if we're using the bond order. 64 d = int(minSeen) 65 # do a quick distance filter 66 if d == 0 or d < minD or d >= maxD: 67 return 68 dist[i] = d 69 70 idx = sigFactory.GetBitIdx(featureSet,dist,sortIndices=False) 71 if _verbose: 72 print('\t',dist,minD,maxD,idx) 73 74 if sigFactory.useCounts: 75 sig[idx] = sig[idx]+1 76 else: 77 sig.SetBit(idx)
78 79
80 -def Gen2DFingerprint(mol,sigFactory,perms=None,dMat=None):
81 """ generates a 2D fingerprint for a molecule using the 82 parameters in _sig_ 83 84 **Arguments** 85 86 - mol: the molecule for which the signature should be generated 87 88 - sigFactory : the SigFactory object with signature parameters 89 NOTE: no preprocessing is carried out for _sigFactory_. 90 It *must* be pre-initialized. 91 92 - perms: (optional) a sequence of permutation indices limiting which 93 pharmacophore combinations are allowed 94 95 - dMat: (optional) the distance matrix to be used 96 97 """ 98 if not isinstance(sigFactory,SigFactory.SigFactory): 99 raise ValueError('bad factory') 100 featFamilies=sigFactory.GetFeatFamilies() 101 if _verbose: 102 print('* feat famillies:',featFamilies) 103 nFeats = len(featFamilies) 104 minCount = sigFactory.minPointCount 105 maxCount = sigFactory.maxPointCount 106 if maxCount>3: 107 logger.warning(' Pharmacophores with more than 3 points are not currently supported.\nSetting maxCount to 3.') 108 maxCount=3 109 110 # generate the molecule's distance matrix, if required 111 if dMat is None: 112 from rdkit import Chem 113 useBO = sigFactory.includeBondOrder 114 dMat = Chem.GetDistanceMatrix(mol,useBO) 115 116 # generate the permutations, if required 117 if perms is None: 118 perms = [] 119 for count in range(minCount,maxCount+1): 120 perms += Utils.GetIndexCombinations(nFeats,count) 121 122 # generate the matches: 123 featMatches = sigFactory.GetMolFeats(mol) 124 if _verbose: 125 print(' featMatches:',featMatches) 126 127 sig = sigFactory.GetSignature() 128 for perm in perms: 129 # the permutation is a combination of feature indices 130 # defining the feature set for a proto-pharmacophore 131 featClasses=[0] 132 for i in range(1,len(perm)): 133 if perm[i]==perm[i-1]: 134 featClasses.append(featClasses[-1]) 135 else: 136 featClasses.append(featClasses[-1]+1) 137 138 # Get a set of matches at each index of 139 # the proto-pharmacophore. 140 matchPerms = [featMatches[x] for x in perm] 141 if _verbose: 142 print('\n->Perm: %s'%(str(perm))) 143 print(' matchPerms: %s'%(str(matchPerms))) 144 145 # Get all unique combinations of those possible matches: 146 matchesToMap=Utils.GetUniqueCombinations(matchPerms,featClasses) 147 for i,entry in enumerate(matchesToMap): 148 entry = [x[1] for x in entry] 149 matchesToMap[i]=entry 150 if _verbose: 151 print(' mtM:',matchesToMap) 152 153 for match in matchesToMap: 154 if sigFactory.shortestPathsOnly: 155 _ShortestPathsMatch(match,perm,sig,dMat,sigFactory) 156 return sig
157