A question recently came up on the mailing list about how to use custom fingerprints with the RDKit PostgreSQL cartridge without having to edit the cartridge code itself. Since the answer isn't trivial but may be useful to others, I'm doing a blog post with the answer.
We'll start the usual way, with a bunch of imports. In the interests of being maximally explicit and having this whole notebook be normal Python code, I'm handling the database connection with the usual psycopg2 connector to PostgreSQL.
from rdkit import Chem
from rdkit import DataStructs
from rdkit import rdBase
from rdkit.Chem.Draw import IPythonConsole
from rdkit.Chem import Draw
import psycopg2,time
print(time.asctime(),"rdkit=",rdBase.rdkitVersion)
Read in a set of molecule data that we've used before:
with open('../data/chembl16_2010-2012.smi') as inf:
data = [x.strip().split() for x in inf]
len(data)
Establish a connection to the database we'll use and insert that data:
conn = psycopg2.connect(database='rdkit_blog_demo')
curs = conn.cursor()
curs.execute('create extension if not exists rdkit')
curs.execute('drop table if exists raw_data')
curs.execute('create table raw_data (smiles text,molregno int)')
curs.executemany('insert into raw_data values (%s,%s)',data)
conn.commit()
Create a molecule table with the first 10K rows (we don't need all the data for the purposes of this post):
curs.execute('drop table if exists mols')
curs.execute('select molregno,mol_from_smiles(smiles::cstring) as m into mols from raw_data limit 10000')
conn.commit()
Now generate our custom fingerprints for those molecules. We're using one of the Atom Pair/Topological Torsion variants that the RDKit provides here:
from rdkit.Chem.AtomPairs import Sheridan
from rdkit.Chem import rdMolDescriptors
fps = []
# grab the molecules, we're pulling them out in their pickled form:
curs.execute('select molregno,mol_send(m) from mols limit 10000')
for molregno,pkl in curs.fetchall():
if pkl is None: continue
# construct a molecule
m = Chem.Mol(pkl.tobytes())
# now do our fingerprint. We're using Topological Torsions with Sheridan's binding properties
# to define atom types
fp = Sheridan.GetBTFingerprint(m,fpfn=rdMolDescriptors.GetHashedTopologicalTorsionFingerprintAsBitVect)
fps.append((molregno,fp))
Now insert the fingerprints into the database, we do this by sending a byte string to the database using DataStructs.BitVectToBinaryText() on the python side and then converting that into a bit vector in the database using the function bfp_from_binary_text():
curs.execute('drop table if exists fps')
curs.execute('create table fps (molregno int, btfp bfp)')
curs.executemany('insert into fps values (%s,bfp_from_binary_text(%s))',
[(x,DataStructs.BitVectToBinaryText(y)) for x,y in fps])
conn.commit()
If this were a larger database I'd also create an index on the fingerprint column in order to speed similarity searches up a bit, but since this one's small we'll skip that.
And here's how you do a query:
fp = fps[-1][-1]
curs.execute('set rdkit.tanimoto_threshold=0.6')
curs.execute('select molregno,m from mols join fps using (molregno) where btfp%%bfp_from_binary_text(%s)',
(DataStructs.BitVectToBinaryText(fp),))
res = curs.fetchall()
len(res)
It's always good to test that we're getting the right answers, so let's verify that by repeating the same process in Python (just to be sure we did everything right!):
len([x for (x,y) in fps if DataStructs.TanimotoSimilarity(y,fp)>=0.6])
Pull back some similarity values:
curs.execute('select molregno,tanimoto_sml(btfp,bfp_from_binary_text(%s)) from fps limit 10',
(DataStructs.BitVectToBinaryText(fp),))
curs.fetchall()
And again, just a test, ensure that we're getting the same thing we'd see in Python:
[(x,DataStructs.TanimotoSimilarity(y,fp)) for (x,y) in fps][:10]
That's it. Not particularly complicated, but still a useful thing to know how to do. Following this approach, any bit vector fingerprint (even one's generated outside of the RDKit) can be inserted into PostgreSQL tables and then searched using the RDKit cartridge.
