Some form of fingerprint is often used to make substructure searching more efficient. The idea is to use a fingerprinting algorithm with the property:
FP(query) & FP(mol) = FP(query)
if query is a substructure of mol. In words: if query
is a substructure of mol
then every bit set in the fingerprint of query
should also be set in mol
.
A bunch of different approaches to this have been developed, I'm not going to review them all here. :-)
Andrew Dalke has done some writing on the topic: http://dalkescientific.com/writings/diary/archive/2012/06/11/optimizing_substructure_keys.html
One of the best-known approaches is the Daylight algorithm: http://www.daylight.com/dayhtml/doc/theory/theory.finger.html
As Andrew mentions in the post I link to above, one of the big problems here is the lack of a reasonable query dataset for benchmarking. The approach I've taken is to use collections of small molecules from ZINC as well as some fragments of PubChem molecules. The database I query is a set of ChEMBL molecules from an earlier post.
Here's the atom-number count distributions of the queries and molecules:
It's not perfect, but it's a start.
Using the RDKit pattern fingerprinter (the default used in the postgresql cartridge, more on that in another post), a screenout accuracy of around 60% is achievable with these three query sets. This means that 60% of molecules passing the fingerprint screen actually have a substructure match.
To put that in perspective: with the leads query set and the pattern fingerprint, 4650 substructure matches are found in a total of 7601 searches, that's 61% accuracy for the fingerprint. If the fingerprint had not been used to pre-screen, 25000000 (500*50000) searches would have needed to be done, so we've reduced the number of substructure searches by 99.97%. By this logic, even a comparatively poor fingerprint with an accuracy of 5% would help enormously: reducing the number of searches by 99.7%. This will become important when actually using the fingerprints in a database index, but that's a different post.
On the technical side: there's some info below about using IPython's ability to run code in parallel.
from rdkit import Chem
from rdkit.Chem import rdMolDescriptors
from rdkit.Avalon import pyAvalonTools
from rdkit.Chem import Draw
from rdkit.Chem.Draw import IPythonConsole
from rdkit import rdBase
from rdkit import DataStructs
import cPickle,random,gzip,time
from __future__ import print_function
print(rdBase.rdkitVersion)
2014.03.1pre
Here we will use the 50K molecules that make up the set of 25K reference pairs generated in an earlier post: http://rdkit.blogspot.ch/2013/10/building-similarity-comparison-set-goal.html
As a quick reminder: these are pairs of molecules taken from ChEMBL with MW<600 and a count-based MFP0 similarity of at least 0.7 to each other.
ind = [x.split() for x in gzip.open('../data/chembl16_25K.pairs.txt.gz')]
mols = []
for i,row in enumerate(ind):
m1 = Chem.MolFromSmiles(row[1])
mols.append(m1)
m2 = Chem.MolFromSmiles(row[3])
mols.append(m2)
We'll use three sets:
These sets were discussed in this thread on the mailing list: http://www.mail-archive.com/rdkit-discuss@lists.sourceforge.net/msg02066.html and this presentation: http://www.hinxton.wellcome.ac.uk/advancedcourses/MIOSS%20Greg%20Landrum.pdf
frags = [Chem.MolFromSmiles(x.split()[0]) for x in file('../data/zinc.frags.500.q.smi')]
leads = [Chem.MolFromSmiles(x.split()[0]) for x in file('../data/zinc.leads.500.q.smi')]
pieces = [Chem.MolFromSmiles(x) for x in file('../data/fragqueries.q.txt')]
Look at the size distributions:
figure(figsize=(16,4))
subplot('141')
hist([x.GetNumAtoms() for x in pieces],bins=20)
title('Pieces')
xlabel('NumAtoms')
subplot('142')
hist([x.GetNumAtoms() for x in frags],bins=20)
title('Fragments')
xlabel('NumAtoms')
subplot('143')
hist([x.GetNumAtoms() for x in leads],bins=20)
title('Leads')
_=xlabel('NumAtoms')
subplot('144')
hist([x.GetNumAtoms() for x in mols],bins=20)
title('Mols')
_=xlabel('NumAtoms')
I'm going to be doing a fair amount of embarassingly parallel computation, so I'll take advantage of IPython's support for parallel computing. For this to work, you need to have a cluster of workers set up. This is easy in the Notebook (there's a tab in the dashboard), and it's doable from the command line (http://ipython.org/ipython-doc/stable/parallel/parallel_process.html)
It took a bit of time to figure out how to do this, but once over that hurdle, this thing is super easy and useful.
from IPython.parallel import Client
rc = Client()
dview = rc[:]
with dview.sync_imports():
from rdkit import Chem
from rdkit import DataStructs
from rdkit.Avalon import pyAvalonTools
importing Chem from rdkit on engine(s) importing DataStructs from rdkit on engine(s) importing pyAvalonTools from rdkit.Avalon on engine(s)
Verify that it works:
t1=time.time()
_ = [Chem.RDKFingerprint(m) for m in mols[:10000]]
t2=time.time()
print('Serial: %.2f'%(t2-t1))
Serial: 15.39
fn = lambda x:Chem.RDKFingerprint(x)
t1=time.time()
_ = dview.map_sync(fn,mols[:10000])
t2=time.time()
print('Parallel: %.2f'%(t2-t1))
Parallel: 5.24
The above runtimes were collected with a cluster of 4 workers, so it looks like the parallization is working.
def subTest(mol,queries,qfps,fpf,verify,silent):
nFound=0
nScreened=0
nTot=0
nFailed = 0
molfp = fpf(mol,False)
for j,queryfp in enumerate(qfps):
nTot += 1
if DataStructs.AllProbeBitsMatch(queryfp,molfp):
nScreened += 1
if mol.HasSubstructMatch(queries[j]):
nFound += 1
elif verify:
if mol.HasSubstructMatch(queries[j]):
nFailed += 1
if not silent:
print('Failure: %s %s'%(Chem.MolToSmiles(mol,True),Chem.MolToSmiles(queries[j],True)))
return nFound,nScreened,nTot,nFailed
def testScreenout(mols,queries,fpf,dview,verify=False,silent=False):
nFound=0
nScreened=0
nTot=0
nFailed = 0
#if not silent: print("Building query fingerprints")
qfps = [fpf(m,True) for m in queries]
#if not silent: print("Running Queries")
t = lambda x,qfps=qfps,queries=queries,fpf=fpf,verify=verify,silent=silent,subTest=subTest:subTest(x,queries,qfps,fpf,verify,silent)
res=dview.map_async(t,mols)
for entry in res:
nFound+=entry[0]
nScreened+=entry[1]
nTot+=entry[2]
nFailed+=entry[3]
if not silent:
accuracy = float(nFound)/nScreened
print("Found %d matches in %d searches with %d failures. Accuracy: %.3f"%(nFound,nScreened,nFailed,accuracy))
return nTot,nScreened,nFound,nFailed
methods = {
'Avalon-2K':lambda x,y: pyAvalonTools.GetAvalonFP(x,nBits=2048,isQuery=y,bitFlags=pyAvalonTools.avalonSSSBits),
'Pattern-1K':lambda x,y: Chem.PatternFingerprint(x,fpSize=1024),
'Pattern-2K':lambda x,y: Chem.PatternFingerprint(x,fpSize=2048),
'Pattern-4K':lambda x,y: Chem.PatternFingerprint(x,fpSize=4096),
'Layered': lambda x,y:Chem.LayeredFingerprint(x,layerFlags=Chem.LayeredFingerprint_substructLayers),
}
leadResults={}
fragResults={}
pieceResults={}
for method,func in methods.iteritems():
print("----------------------------")
print("Doing %s"%method)
print("Leads")
leadResults[method]=testScreenout(mols,leads,func,dview)
print("Frags")
fragResults[method]=testScreenout(mols,frags,func,dview)
print("Pieces")
pieceResults[method]=testScreenout(mols,pieces,func,dview)
---------------------------- Doing Avalon-2K Leads Found 886 matches in 2044 searches with 0 failures. Accuracy: 0.433 Frags Found 4479 matches in 13656 searches with 0 failures. Accuracy: 0.328 Pieces Found 1926728 matches in 4141259 searches with 0 failures. Accuracy: 0.465 ---------------------------- Doing Layered Leads Found 1274 matches in 6797 searches with 0 failures. Accuracy: 0.187 Frags Found 4659 matches in 44258 searches with 0 failures. Accuracy: 0.105 Pieces Found 1935315 matches in 4642461 searches with 0 failures. Accuracy: 0.417 ---------------------------- Doing Pattern-2K Leads Found 1274 matches in 1781 searches with 0 failures. Accuracy: 0.715 Frags Found 4659 matches in 7892 searches with 0 failures. Accuracy: 0.590 Pieces Found 1935315 matches in 3381436 searches with 0 failures. Accuracy: 0.572 ---------------------------- Doing Pattern-1K Leads Found 1274 matches in 7527 searches with 0 failures. Accuracy: 0.169 Frags Found 4659 matches in 16372 searches with 0 failures. Accuracy: 0.285 Pieces Found 1935315 matches in 3874208 searches with 0 failures. Accuracy: 0.500 ---------------------------- Doing Pattern-4K Leads Found 1274 matches in 1748 searches with 0 failures. Accuracy: 0.729 Frags Found 4659 matches in 7742 searches with 0 failures. Accuracy: 0.602 Pieces Found 1935315 matches in 3048867 searches with 0 failures. Accuracy: 0.635
mns = sorted(methods.keys())
for mn in mns:
print(mn,'Pieces','%.3f'%(float(pieceResults[mn][2])/pieceResults[mn][1]),'Fragments','%.3f'%(float(fragResults[mn][2])/fragResults[mn][1]),
'Leads','%.3f'%(float(leadResults[mn][2])/leadResults[mn][1]))
Avalon-2K Pieces 0.465 Fragments 0.328 Leads 0.433 Layered Pieces 0.417 Fragments 0.105 Leads 0.187 Pattern-1K Pieces 0.500 Fragments 0.285 Leads 0.169 Pattern-2K Pieces 0.572 Fragments 0.590 Leads 0.715 Pattern-4K Pieces 0.635 Fragments 0.602 Leads 0.729
The results from the Avalon fingerprint shown above are not directly comparable to the others since structures are being filtered out that shouldn't be.
To demonstrate this, we need a serial form of the screenout test:
def testScreenoutSerial(mols,queries,fpf,verify=False,silent=False):
nFound=0
nScreened=0
nTot=0
nFailed = 0
#if not silent: print("Building query fingerprints")
qfps = [fpf(m,True) for m in queries]
#if not silent: print("Running Queries")
t = lambda x,qfps=qfps,queries=queries,fpf=fpf,verify=verify,silent=silent,subTest=subTest:subTest(x,queries,qfps,fpf,verify,silent)
res=map(t,mols)
for entry in res:
nFound+=entry[0]
nScreened+=entry[1]
nTot+=entry[2]
nFailed+=entry[3]
if not silent:
accuracy = float(nFound)/nScreened
print("Found %d matches in %d searches with %d failures. Accuracy: %.3f"%(nFound,nScreened,nFailed,accuracy))
return nTot,nScreened,nFound,nFailed
And now we can get some examples:
testScreenoutSerial(mols[:500],leads,methods['Avalon-2K'],verify=True)
Failure: Cc1ccc(C)c(S(=O)(=O)c2nnn3c4ccsc4c(Nc4ccc(C)c(C)c4)nc23)c1 c1cn[nH]n1 Failure: COc1ccc(OCC#Cc2cn([C@H](C)C[C@@H]3CC[C@H]([C@@H](C)C(=O)N(C)Cc4ccccc4)O3)nn2)cc1 c1cn[nH]n1 Failure: Cc1c(C(=O)Nc2ccccc2)nnn1Cc1ccccc1O c1cn[nH]n1 Failure: c1c(CCc2ccccc2)nnn1C[C@H]1CC[C@@H]([C@@H]2CC[C@H](Cn3cc(CCc4ccccc4)nn3)O2)O1 c1cn[nH]n1 Failure: O=C(NCc1ccc(F)cc1)[C@@H](CCN1[C@@H]2CC[C@@H]1CC(n1nnc3cccnc31)C2)c1ccccc1 c1cn[nH]n1 Found 10 matches in 22 searches with 5 failures. Accuracy: 0.455
(250000, 22, 10, 5)
Look at a specific example:
mfp = pyAvalonTools.GetAvalonFP('Cc1ccc(C)c(S(=O)(=O)c2nnn3c4ccsc4c(Nc4ccc(C)c(C)c4)nc23)c1',True,2048,isQuery=False)
qfp = pyAvalonTools.GetAvalonFP('c1cn[nH]n1',True,2048,isQuery=True)
DataStructs.AllProbeBitsMatch(qfp,mfp)
False
Chem.MolFromSmiles('Cc1ccc(C)c(S(=O)(=O)c2nnn3c4ccsc4c(Nc4ccc(C)c(C)c4)nc23)c1')
The normal explanation for this would be differences in the aromaticity model. In this case, however, that turns out not be exactly it.
Here's the molecule we're querying with:
Chem.MolFromSmiles('c1cn[nH]n1')
qfp2 = pyAvalonTools.GetAvalonFP('c1cnn[nH]1',True,2048,isQuery=True)
qfp==qfp2
False
The deciding factor is the tautomer. The Avalon fingerprinter considers this significant, and that difference explains the different screenout.
If we use the fingerprint from the second tautomer we see the expected screenout behavior:
DataStructs.AllProbeBitsMatch(qfp2,mfp)
True
This ends up mattering because the RDKit's substructure matcher ignores the H completely:
Chem.MolFromSmiles('c1cnn[nH]1').HasSubstructMatch(Chem.MolFromSmiles('c1cn[nH]n1'))
True
Chem.MolFromSmiles('c1cnnn1C').HasSubstructMatch(Chem.MolFromSmiles('c1cn[nH]n1'))
True
Unless, of course, you construct the query from SMARTS:
Chem.MolFromSmiles('c1cnnn1C').HasSubstructMatch(Chem.MolFromSmarts('c1cn[nH]n1'))
False
fn = lambda x,v:Chem.RDKFingerprint(x,minPath=1,maxPath=5,nBitsPerHash=1,useHs=False)
testScreenoutSerial(mols,leads,fn,verify=True)
Failure: COc1ccc(-n2c(SC(C)=O)nc3sc4c(c3c2=O)CCCC4)cc1 CSC(C)=O Failure: CCCCCC[C@H]1SC(=O)c2ccccc2[C@@H]1C(=O)O CSC(C)=O Failure: COc1ccccc1SC(=O)c1cccc(C=O)n1 CSC(C)=O Failure: COc1ccc(-n2c(SC(C)=O)nc3sc4c(c3c2=O)CCCC4)cc1 CSC(C)=O Failure: CC(C)(C)OC(=O)n1c(-c2ccc3c(c2)CC(NS(=O)(=O)c2ccccc2)C3)cc2cc(O)ccc21 Cc1cc2cc(O)ccc2[nH]1 Failure: Cc1c2cc(O)ccc2n(Cc2ccc(OC[C@H](C)N3CCCC3)cc2)c1-c1ccc(O)cc1 Cc1cc2cc(O)ccc2[nH]1 Failure: CC1(C)Oc2ccc(-c3nc(-c4cccc5c(CCC(=O)O)c[nH]c54)no3)cc2O1 Cc1noc(-c2ccc(O)cc2)n1 Failure: COc1ccc(CC(=O)Nc2ccc(-c3noc(-c4cc(OC)c(OC)c(OC)c4)n3)cc2)cc1 Cc1noc(-c2ccc(O)cc2)n1 Failure: CO[C@H]1C(=O)N2C(C(=O)C(C)(C)C)=C(C)[C@@H](S(=O)(=O)c3ccccc3)S(=O)(=O)C12 CSCSC Failure: Cc1c2c(nc3ccc(OCCN4CCCCC4)cc13)-c1cc3cc(O)ccc3n1C2 Cc1cc2cc(O)ccc2[nH]1 Failure: COc1ccc2[nH]c(-c3ccc4ccccc4c3)c(-c3cc(OC)c(OC)c(OC)c3)c2c1 Cc1cc2cc(O)ccc2[nH]1 Failure: O=C(O)C1CCCN1C(=O)CC(SC(=O)c1ccccc1)C(=O)c1ccccc1 CSC(C)=O Failure: CC(C)(C)OC(=O)n1c(-c2ccc3c(c2)CC(NS(=O)(=O)c2ccccc2)C3)cc2cc(O)ccc21 Cc1cc2cc(O)ccc2[nH]1 Failure: O=Cc1cccc(C(=O)Sc2ccc(Cl)cc2)n1 CSC(C)=O Failure: O=C1NCc2cccc(-c3cc4cc(OCCN5CCCCC5)ccc4[nH]3)c21 Cc1cc2cc(O)ccc2[nH]1 Failure: COc1ccc(C#CCC2(S(=O)(=O)c3ccc(C)cc3)SC(=O)NC2=O)cc1 CSCSC Failure: CCCC(=O)Sc1nc2sc3c(c2c(=O)n1-c1ccc(OC)cc1)CCCC3 CSC(C)=O Failure: COc1ccc(-n2c(SC(C)=O)nc3sc4c(c3c2=O)CCCC4)cc1 CSC(C)=O Failure: COc1ccc(-c2nc(-c3ccc(NC(=O)c4ccc(Br)o4)cc3)no2)cc1OC Cc1noc(-c2ccc(O)cc2)n1 Failure: COc1cc2c(cc1OC)C(=O)N1CSCC1C(=O)S2 CSC(C)=O Failure: Oc1nc2ccccc2cc1-c1cc2cc(OCCN3CCCCC3)ccc2[nH]1 Cc1cc2cc(O)ccc2[nH]1 Failure: CC(C)(C)OC(=O)[C@@H]1CCCN1C(=O)[C@H]1CCC[C@@H]1SC(=O)c1ccccc1 CSC(C)=O Failure: O=C1Sc2ccccc2C(=O)N2CCSC12 CSC(C)=O Failure: CC(C)(C)OC(=O)n1c(-c2ccc3c(c2)CC(NS(=O)(=O)c2ccccc2)C3)cc2cc(O)ccc21 Cc1cc2cc(O)ccc2[nH]1 Failure: O=c1n(Cc2ccccc2)nc2n1-c1ccccc1OC2 CCn1nc2n(c1=O)-c1ccccc1OC2 Failure: COc1ccc(-n2c(SC(C)=O)nc3sc4c(c3c2=O)CCCC4)cc1 CSC(C)=O Failure: CCCOC(=O)c1c(CCC)c(C(=O)SCC)c(CC)nc1-c1cccc(F)c1 CSC(C)=O Failure: O=C1Sc2sc(=S)sc2SC1c1ccccc1 CSC(C)=O Failure: COc1ccc(-n2c(SC(C)=O)nc3sc4c(c3c2=O)CCCC4)cc1 CSC(C)=O Failure: COc1ccc2c(c1)c(CCN(C)C)c1n2S(=O)(=O)c2ccccc2-1 Cc1cc2cc(O)ccc2[nH]1 Failure: CC(C)(C)OC(=O)n1c(-c2ccc3c(c2)CC(NS(=O)(=O)c2ccccc2)C3)cc2cc(O)ccc21 Cc1cc2cc(O)ccc2[nH]1 Failure: CCCC(=O)Sc1nc2sc3c(c2c(=O)n1-c1ccc(OC)cc1)CCCC3 CSC(C)=O Failure: COc1ccc(-c2cc3cc(OC)ccc3[nH]2)cc1 Cc1cc2cc(O)ccc2[nH]1 Failure: O=C(O)[C@@H](Cc1ccc(-c2ccccc2)cc1)SC(=O)c1ccccc1 CSC(C)=O Found 1240 matches in 1399 searches with 34 failures. Accuracy: 0.886
(25000000, 1399, 1240, 34)
The failures happen here because the RDKit fingerprinter uses information about aromaticity in the calculation of atom invariants and while hashing bonds. That extra information is a big part of why the accuracy is so high.
Turning off all bond order information (finer grained control is unfortunately not possible) and using atomic number as the atom invariants gives no errors, but a much lower screenout rate:
def func(mol,v):
invs = [x.GetAtomicNum() for x in mol.GetAtoms()]
return Chem.RDKFingerprint(mol,minPath=1,maxPath=5,nBitsPerHash=1,useHs=False,useBondOrder=False,atomInvariants=invs)
testScreenoutSerial(mols,leads,func,verify=True)
Found 1274 matches in 125651 searches with 0 failures. Accuracy: 0.010
(25000000, 125651, 1274, 0)
Using a shorter maxPath length also hurts accuracy:
fn = lambda x,v:Chem.RDKFingerprint(x,minPath=1,maxPath=4,nBitsPerHash=1,useHs=False)
testScreenout(mols,leads,fn,dview)
Found 1240 matches in 2331 searches with 0 failures. Accuracy: 0.532
(25000000, 2331, 1240, 0)
It is most likely luck that there are no failures here.
Andrew Dalke's chemfp package (http://chemfp.com/) includes a pattern-based fingerprinter which is based on the PubChem/CACTVS substructure keys:
from chemfp import rdkit_patterns
dalke_fp=rdkit_patterns.SubstructRDKitFingerprinter_v1.make_fingerprinter()
fn = lambda x,v,fp=dalke_fp:DataStructs.CreateFromBinaryText(fp(x))
testScreenoutSerial(mols[:500],leads,fn,verify=True)
Failure: CC(C)(C)c1cc2c(cc1SC1=C(O)OC(CCc3ccccc3)(c3ccccc3)CC1=O)CCCC2 C1=COCCC1 Failure: CC(C)(C)c1cc2c(cc1SC1=C(O)OC(CCc3ccccc3)(c3ccccc3)CC1=O)CCCC2 O=C1C=CO[C@H](c2ccccc2)C1 Failure: Cc1ccc(C)c(S(=O)(=O)c2nnn3c4ccsc4c(Nc4ccc(C)c(C)c4)nc23)c1 c1cn[nH]n1 Failure: CCOC(=O)Cc1c(C)n(CC2CCCCC2)c2ccc(OC)cc12 Cc1cc2cc(O)ccc2[nH]1 Failure: COc1ccc(OCC#Cc2cn([C@H](C)C[C@@H]3CC[C@H]([C@@H](C)C(=O)N(C)Cc4ccccc4)O3)nn2)cc1 c1cn[nH]n1 Failure: COC1=C(OC)C(=O)C2=C(C[C@@H]3[C@@]4(C)C[C@H]4C[C@@]3(C)C2)C1=O COC1=CCC=CC1 Failure: COc1ccc2[nH]c3c(c2c1)CC(NC(=O)C1CC1)CC3 Cc1cc2cc(O)ccc2[nH]1 Failure: Cc1c(C(=O)Nc2ccccc2)nnn1Cc1ccccc1O c1cn[nH]n1 Failure: COc1c(=O)[nH]c2ccccc2c1-c1ccccc1 Oc1cnc2ccccc2c1 Failure: Cc1ccc2c(nc(-c3ccc(NC(=O)C(F)(F)F)cc3)c(O)c2C(=O)O)c1C Oc1cnc2ccccc2c1 Failure: c1c(CCc2ccccc2)nnn1C[C@H]1CC[C@@H]([C@@H]2CC[C@H](Cn3cc(CCc4ccccc4)nn3)O2)O1 c1cn[nH]n1 Failure: O=C(NCc1ccc(F)cc1)[C@@H](CCN1[C@@H]2CC[C@@H]1CC(n1nnc3cccnc31)C2)c1ccccc1 c1cn[nH]n1 Found 3 matches in 108 searches with 12 failures. Accuracy: 0.028
(250000, 108, 3, 12)
It's not particularly effective as a screening fingerprint and, as seen above, creates some holes.
Another key-based fingerprint that isn't particularly effective and results in some misses:
fn = lambda x,v:rdMolDescriptors.GetMACCSKeysFingerprint(x)
testScreenoutSerial(mols[:500],leads,fn,verify=True)
Failure: CC(C)(C)c1cc2c(cc1SC1=C(O)OC(CCc3ccccc3)(c3ccccc3)CC1=O)CCCC2 C1=COCCC1 Failure: Cc1ccc(C)c(S(=O)(=O)c2nnn3c4ccsc4c(Nc4ccc(C)c(C)c4)nc23)c1 c1cn[nH]n1 Failure: CCOC(=O)C1=COC(C)(CCc2ccccc2)CC1=O C1=COCCC1 Failure: CCOC(=O)Cc1c(C)n(CC2CCCCC2)c2ccc(OC)cc12 Cc1cc2cc(O)ccc2[nH]1 Failure: COc1ccc(OCC#Cc2cn([C@H](C)C[C@@H]3CC[C@H]([C@@H](C)C(=O)N(C)Cc4ccccc4)O3)nn2)cc1 c1cn[nH]n1 Failure: COC1=C(OC)C(=O)C2=C(C[C@@H]3[C@@]4(C)C[C@H]4C[C@@]3(C)C2)C1=O COC1=CCC=CC1 Failure: COc1ccc2[nH]c3c(c2c1)CC(NC(=O)C1CC1)CC3 Cc1cc2cc(O)ccc2[nH]1 Failure: Cc1c(C(=O)Nc2ccccc2)nnn1Cc1ccccc1O c1cn[nH]n1 Failure: COc1c(=O)[nH]c2ccccc2c1-c1ccccc1 Oc1cnc2ccccc2c1 Failure: CC(=O)NCCCC[C@H](NC(C)=O)C(=O)NCC(=O)S[C@H](C)C(=O)O CSC(C)=O Failure: c1c(CCc2ccccc2)nnn1C[C@H]1CC[C@@H]([C@@H]2CC[C@H](Cn3cc(CCc4ccccc4)nn3)O2)O1 c1cn[nH]n1 Failure: O=C(NCc1ccc(F)cc1)[C@@H](CCN1[C@@H]2CC[C@@H]1CC(n1nnc3cccnc31)C2)c1ccccc1 c1cn[nH]n1 Found 3 matches in 134 searches with 12 failures. Accuracy: 0.022
(250000, 134, 3, 12)
A first pass at improving is to combine some selected patterns with an algorithmic fingerprinter like the pattern fingerprinter. Here's a thread describing an experiment:
http://www.mail-archive.com/rdkit-discuss@lists.sourceforge.net/msg02078.html
And a reproduction of that with this test harness:
qs="""0 2 times O largest 55458
1 2 times Ccc largest 29602
2 1 times CCN largest 16829
3 1 times cnc largest 11439
4 1 times cN largest 8998
5 1 times C=O largest 7358
6 1 times CCC largest 6250
7 1 times S largest 4760
8 1 times c1ccccc1 largest 4524
9 2 times N largest 2854
10 1 times C=C largest 2162
11 1 times nn largest 1840
12 2 times CO largest 1248
13 1 times Ccn largest 964
14 1 times CCCCC largest 857
15 1 times cc(c)c largest 653
16 3 times O largest 653
17 1 times O largest 466
18 2 times CNC largest 464
19 1 times s largest 457
20 1 times CC(C)C largest 335
21 1 times o largest 334
22 1 times cncnc largest 334
23 1 times C=N largest 321
24 2 times CC=O largest 238
25 4 times Ccc largest 238
26 1 times Cl largest 230
27 4 times O largest 149
28 2 times ccncc largest 149
29 6 times CCCCCC largest 76
30 2 times c1ccccc1 largest 76
31 1 times F largest 75
32 3 times CCOC largest 44
33 3 times N largest 44
34 1 times c(cn)n largest 44
35 1 times N largest 41
36 9 times C largest 41
37 1 times CC=C(C)C largest 33
38 1 times c1ccncc1 largest 26
39 1 times CC(C)N largest 26
40 1 times CC largest 26
41 4 times CCC(C)O largest 25
42 2 times ccc(cc)n largest 21
43 6 times C largest 21
44 1 times C1CCCC1 largest 18
45 1 times C largest 18
46 5 times O largest 18
47 2 times Ccn largest 14
48 1 times CNCN largest 13
49 3 times cncn largest 13
50 1 times CSC largest 13
51 3 times CC=O largest 11
52 1 times CCNCCCN largest 11
53 1 times CccC largest 11
54 3 times ccccc(c)c largest 10"""
def _initPatts():
ssPatts=[]
for q in qs.split('\n'):
q = q.split(' ')
count = int(q[1])
patt = Chem.MolFromSmiles(q[3],sanitize=False)
patt.UpdatePropertyCache(strict=False) # <- github #149, without this we cannot pickle the queries
ssPatts.append((patt,count))
return ssPatts
_ssPatts=None
def GetCombinedSubstructFP(m,ssPatts=None,fpSize=1024,verbose=False):
if ssPatts is None:
ssPatts = _initPatts()
sz = len(ssPatts)
lfp=Chem.PatternFingerprint(m,fpSize=fpSize)
res = DataStructs.ExplicitBitVect(fpSize+sz)
obls = [x+sz for x in lfp.GetOnBits()]
res.SetBitsFromList(obls)
for i,(p,count) in enumerate(ssPatts):
matches = m.GetSubstructMatches(p,uniquify=True)
if len(matches)>=count:
res.SetBit(i)
return res
ssPatts = _initPatts()
fn = lambda x,v,ssPatts=ssPatts:GetCombinedSubstructFP(x,ssPatts=ssPatts,fpSize=2048)
testScreenoutSerial(mols[:500],leads,fn,verify=True)
Found 15 matches in 23 searches with 0 failures. Accuracy: 0.652
(250000, 23, 15, 0)
leadResults={}
fragResults={}
pieceResults={}
print("----------------------------")
fn = lambda x,v,ssPatts=ssPatts,fpf=GetCombinedSubstructFP:fpf(x,ssPatts=ssPatts,fpSize=1024)
method='1024'
print("Doing %s"%method)
print("Pieces")
pieceResults[method]=testScreenout(mols,pieces,fn,dview)
print("Frags")
fragResults[method]=testScreenout(mols,frags,fn,dview)
print("Leads")
leadResults[method]=testScreenout(mols,leads,fn,dview)
print("----------------------------")
fn = lambda x,v,ssPatts=ssPatts,fpf=GetCombinedSubstructFP:fpf(x,ssPatts=ssPatts,fpSize=2048)
method='2048'
print("Doing %s"%method)
print("Pieces")
pieceResults[method]=testScreenout(mols,pieces,fn,dview)
print("Frags")
fragResults[method]=testScreenout(mols,frags,fn,dview)
print("Leads")
leadResults[method]=testScreenout(mols,leads,fn,dview)
---------------------------- Doing 1024 Pieces Found 1935315 matches in 3023990 searches with 0 failures. Accuracy: 0.640 Frags Found 4659 matches in 13765 searches with 0 failures. Accuracy: 0.338 Leads Found 1274 matches in 6857 searches with 0 failures. Accuracy: 0.186 ---------------------------- Doing 2048 Pieces Found 1935315 matches in 2810633 searches with 0 failures. Accuracy: 0.689 Frags Found 4659 matches in 7669 searches with 0 failures. Accuracy: 0.608 Leads Found 1274 matches in 1761 searches with 0 failures. Accuracy: 0.723
Comparing to the previous pattern fp results:
Pattern-1K Pieces 0.500 Fragments 0.285 Leads 0.169
Combine-1K Pieces 0.640 Fragments 0.338 Leads 0.186
Pattern-2K Pieces 0.572 Fragments 0.590 Leads 0.715
Combine-2K Pieces 0.689 Fragments 0.608 Leads 0.723
Pattern-4K Pieces 0.635 Fragments 0.602 Leads 0.729
These help a bit, particularly with the pieces, but not a whole lot.