{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.13","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":67356,"databundleVersionId":8006601,"sourceType":"competition"}],"dockerImageVersionId":30715,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Mining Fingerprint Collisions\n\nWe know that molecular fingerprints suffer from hash collisions where different molecular substructures are mapped to the same bit value. We can't determine ahead of time which features will have collisions, but we can mine them empirically from existing data.\n\nThis notebook shows how to extract SMARTS string associated with fingerprint features to analyze collisions.","metadata":{}},{"cell_type":"code","source":"!pip install rdkit","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:41:33.239417Z","iopub.execute_input":"2024-06-05T04:41:33.239814Z","iopub.status.idle":"2024-06-05T04:41:53.186998Z","shell.execute_reply.started":"2024-06-05T04:41:33.239782Z","shell.execute_reply":"2024-06-05T04:41:53.185608Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import os\n\nimport numpy as np\nimport matplotlib.pyplot as plt\n\nfrom rdkit import Chem\nfrom rdkit.Chem import rdMolDescriptors\nfrom rdkit.Chem import DataStructs\nfrom rdkit.Chem import rdFingerprintGenerator\n\nfrom collections import defaultdict\n\nimport datasets\nimport pandas as pd","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:41:56.807503Z","iopub.execute_input":"2024-06-05T04:41:56.807970Z","iopub.status.idle":"2024-06-05T04:41:59.130996Z","shell.execute_reply.started":"2024-06-05T04:41:56.807933Z","shell.execute_reply":"2024-06-05T04:41:59.129755Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We'll start by grabbing some molecules. For this analysis, we only need SMILES strings.\n\nTo keep the notebook quick, we'll only look at 2000 molecules from the test dataset. For a more in dept analysis, you can run this on the full dataset.","metadata":{}},{"cell_type":"code","source":"df = pd.read_csv('/kaggle/input/leash-BELKA/test.csv')\nsmiles = df.molecule_smiles.unique()\n\n# take 2k smiles\ndataset = datasets.Dataset.from_dict({'smiles' : smiles[:2000]})","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:00.816529Z","iopub.execute_input":"2024-06-05T04:42:00.817815Z","iopub.status.idle":"2024-06-05T04:42:08.374558Z","shell.execute_reply.started":"2024-06-05T04:42:00.817767Z","shell.execute_reply":"2024-06-05T04:42:08.373304Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"To start, lets look at a single fingerprint.\n\nCreating the fingerprint with `AdditionalOutput` allows us to inspect the features that map to each fingerprint bit","metadata":{}},{"cell_type":"code","source":"def get_fp_with_ao(mol, radius=3, fpSize=2048):\n    fpg = rdFingerprintGenerator.GetMorganGenerator(radius=radius,fpSize=fpSize)\n    \n    ao = rdFingerprintGenerator.AdditionalOutput()\n    ao.AllocateAtomCounts()\n    ao.AllocateAtomToBits()\n    ao.AllocateBitInfoMap()\n    \n    fp = fpg.GetFingerprint(mol, additionalOutput=ao)\n    return fp, ao","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:08.376462Z","iopub.execute_input":"2024-06-05T04:42:08.376857Z","iopub.status.idle":"2024-06-05T04:42:08.383933Z","shell.execute_reply.started":"2024-06-05T04:42:08.376823Z","shell.execute_reply":"2024-06-05T04:42:08.382591Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Lets look at a random molecule.\n\nNote that here we are removing chirality from the molecule (`useChiral=0`). The analysis in this notebook will not look at chirality. This is to keep things consistent with our fingerprint function. By default, Morgan fingerprints do not look at chirality. You can change this by passing `includeChirality=True` to the fingerprint generator.\n\nSince our fingerprint function doesn't look at chirality, we expect enantiomers to map to the same bit value. Removing chirality prevents us from counting enantiomers as a bit collision. If you do want to treat enantiomers as a bit collision, you can run this notebook using chiral fingerprints.","metadata":{}},{"cell_type":"code","source":"mol = Chem.MolFromSmiles(Chem.CanonSmiles(dataset[0]['smiles'], useChiral=0))\nmol","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:12.412105Z","iopub.execute_input":"2024-06-05T04:42:12.412525Z","iopub.status.idle":"2024-06-05T04:42:12.441493Z","shell.execute_reply.started":"2024-06-05T04:42:12.412489Z","shell.execute_reply":"2024-06-05T04:42:12.440168Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fp, ao = get_fp_with_ao(mol, 3, 2048)","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:13.524438Z","iopub.execute_input":"2024-06-05T04:42:13.524877Z","iopub.status.idle":"2024-06-05T04:42:13.531948Z","shell.execute_reply.started":"2024-06-05T04:42:13.524841Z","shell.execute_reply":"2024-06-05T04:42:13.530218Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We can inspect the bit info map of the fingerprint additional output. This is a dict of bit values with `(atom_index, radius)` values denoting the central atom and radius of the feature that mapped to that specific bit.","metadata":{}},{"cell_type":"code","source":"bit_map = ao.GetBitInfoMap()\nbit_map","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:19.062157Z","iopub.execute_input":"2024-06-05T04:42:19.062579Z","iopub.status.idle":"2024-06-05T04:42:19.082964Z","shell.execute_reply.started":"2024-06-05T04:42:19.062544Z","shell.execute_reply":"2024-06-05T04:42:19.081746Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def get_bit_atoms_and_bonds(mol, atom_idx, radius):\n    \n    env = Chem.FindAtomEnvironmentOfRadiusN(mol, radius, atom_idx)\n\n    atoms = set((atom_idx, ))\n    bonds = set(env)\n    \n    for bond_idx in env:\n        bond = mol.GetBondWithIdx(bond_idx)\n        atoms.update([bond.GetBeginAtomIdx(), \n                      bond.GetEndAtomIdx()])\n        \n    for atom_idx in atoms:\n        atom = mol.GetAtomWithIdx(atom_idx)\n        for bond in atom.GetBonds():\n            bonds.add(bond.GetIdx())\n            \n    return atoms, bonds","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:20.284207Z","iopub.execute_input":"2024-06-05T04:42:20.284735Z","iopub.status.idle":"2024-06-05T04:42:20.296564Z","shell.execute_reply.started":"2024-06-05T04:42:20.284690Z","shell.execute_reply":"2024-06-05T04:42:20.295166Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We can extract the substructure defined by the `(atom_index, radius)` tuple and visualize it on the molecule","metadata":{}},{"cell_type":"code","source":"bit_idx = 1433\natom_idx, radius = bit_map[bit_idx][0]\nbit_atoms, bit_bonds = get_bit_atoms_and_bonds(mol, atom_idx, radius)\n\nChem.Draw.MolToImage(mol, \n                     highlightAtoms=list(bit_atoms), \n                     highlightBonds=list(bit_bonds))","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:21.455820Z","iopub.execute_input":"2024-06-05T04:42:21.456224Z","iopub.status.idle":"2024-06-05T04:42:21.500545Z","shell.execute_reply.started":"2024-06-05T04:42:21.456191Z","shell.execute_reply":"2024-06-05T04:42:21.499217Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We can also use this information to extract the substructure as a SMARTS string","metadata":{}},{"cell_type":"code","source":"def get_morgan_smarts(mol, atom_idx, radius):\n    mol_copy = Chem.Mol(mol)\n    \n    atoms, bonds = get_bit_atoms_and_bonds(mol_copy, atom_idx, radius)\n            \n    feature_neighbors = []\n    for atom_idx in atoms:\n        neighbors = mol_copy.GetAtomWithIdx(atom_idx).GetNeighbors()\n\n        for neighbor in neighbors:\n            neighbor_idx = neighbor.GetIdx()\n            if neighbor_idx not in atoms:\n                feature_neighbors.append(neighbor_idx)\n\n    for atom_idx in feature_neighbors:\n        mol_copy.GetAtomWithIdx(atom_idx).SetAtomicNum(0)\n\n    submol = Chem.PathToSubmol(mol_copy, list(bonds))\n    bit_smarts = Chem.MolToSmarts(submol).replace('[#0]', '*')\n    return bit_smarts","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:23.142125Z","iopub.execute_input":"2024-06-05T04:42:23.142549Z","iopub.status.idle":"2024-06-05T04:42:23.151549Z","shell.execute_reply.started":"2024-06-05T04:42:23.142513Z","shell.execute_reply":"2024-06-05T04:42:23.150256Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Here we get the SMARTS string for the feature visualized above. We can see it matches the highlighted region.","metadata":{}},{"cell_type":"code","source":"bit_smarts = get_morgan_smarts(mol, atom_idx, radius)\nprint(bit_smarts)","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:24.041995Z","iopub.execute_input":"2024-06-05T04:42:24.042375Z","iopub.status.idle":"2024-06-05T04:42:24.049366Z","shell.execute_reply.started":"2024-06-05T04:42:24.042346Z","shell.execute_reply":"2024-06-05T04:42:24.047955Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"Chem.Draw.MolsToGridImage([mol, Chem.MolFromSmarts(bit_smarts)], \n                          highlightAtomLists=[list(bit_atoms), []],\n                          highlightBondLists=[list(bit_bonds), []],\n                          legends=['Input Molecule', f'Bit {bit_idx} SMARTS'],\n                          molsPerRow=2,\n                          subImgSize=(300,300)\n                         )","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:24.827000Z","iopub.execute_input":"2024-06-05T04:42:24.827430Z","iopub.status.idle":"2024-06-05T04:42:24.871589Z","shell.execute_reply.started":"2024-06-05T04:42:24.827394Z","shell.execute_reply":"2024-06-05T04:42:24.869742Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We can use this smarts string to run a substructure query against the full molecule. In this case, we see two matches. This matches the information in the bit map showing that `(14, 3)` and `(25, 3)` both map to bit 1433. In this case, we don't have a collision - we have the same feature in two parts of the molecule mapping to the same fingerprint bit.","metadata":{}},{"cell_type":"code","source":"print(mol.GetSubstructMatches(Chem.MolFromSmarts(bit_smarts)))\nprint(bit_map[bit_idx])","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:25.955417Z","iopub.execute_input":"2024-06-05T04:42:25.955987Z","iopub.status.idle":"2024-06-05T04:42:25.962738Z","shell.execute_reply.started":"2024-06-05T04:42:25.955948Z","shell.execute_reply":"2024-06-05T04:42:25.961505Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def plot_bit_features(mol, bit_info, bit_idx, **draw_kwargs):\n    bit_tuples = bit_info[bit_idx]\n\n    plot_mols = []\n    atom_highlights = []\n    bond_highlights = []\n\n    for (atom_idx, radius) in bit_tuples:\n        bit_atoms, bit_bonds = get_bit_atoms_and_bonds(mol, atom_idx, radius)\n        atom_highlights.append(list(bit_atoms))\n        bond_highlights.append(list(bit_bonds))\n        plot_mols.append(mol)\n\n    img = Chem.Draw.MolsToGridImage(plot_mols, \n                              highlightAtomLists=atom_highlights, \n                              highlightBondLists=bond_highlights,\n                              **draw_kwargs\n                             )\n    return img","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:26.540517Z","iopub.execute_input":"2024-06-05T04:42:26.540928Z","iopub.status.idle":"2024-06-05T04:42:26.549353Z","shell.execute_reply.started":"2024-06-05T04:42:26.540896Z","shell.execute_reply":"2024-06-05T04:42:26.548020Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We can also plot all substructures associated with a given bit","metadata":{}},{"cell_type":"code","source":"plot_bit_features(mol, bit_map, bit_idx, molsPerRow=2, subImgSize=(300,300))","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:27.630664Z","iopub.execute_input":"2024-06-05T04:42:27.631061Z","iopub.status.idle":"2024-06-05T04:42:27.667984Z","shell.execute_reply.started":"2024-06-05T04:42:27.631030Z","shell.execute_reply":"2024-06-05T04:42:27.666502Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Now lets extract all the smarts associated with every bit and see if there are any collisions","metadata":{}},{"cell_type":"code","source":"def get_fp_smarts(mol, radius, fp_size, canonicalize=False):\n    fp, ao = get_fp_with_ao(mol, radius, fp_size)\n    \n    output = {}\n    \n    for bit_idx, fp_pairs in ao.GetBitInfoMap().items():\n        output[bit_idx] = set()\n        \n        for (atom_idx, radius) in fp_pairs:\n            bit_smarts = get_morgan_smarts(mol, atom_idx, radius)\n            if canonicalize:\n                bit_smarts = Chem.MolToSmiles(\n                                Chem.MolFromSmarts(bit_smarts))\n            output[bit_idx].update([bit_smarts])\n            \n    return output","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:29.132988Z","iopub.execute_input":"2024-06-05T04:42:29.133379Z","iopub.status.idle":"2024-06-05T04:42:29.141804Z","shell.execute_reply.started":"2024-06-05T04:42:29.133348Z","shell.execute_reply":"2024-06-05T04:42:29.140178Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"bit_smarts = get_fp_smarts(mol, 3, 2048)","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:29.614751Z","iopub.execute_input":"2024-06-05T04:42:29.615695Z","iopub.status.idle":"2024-06-05T04:42:29.683277Z","shell.execute_reply.started":"2024-06-05T04:42:29.615637Z","shell.execute_reply":"2024-06-05T04:42:29.682082Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Here it looks like we have a good number of collided bits. But if we look closely, most of these are near-identical substructures","metadata":{}},{"cell_type":"code","source":"collided_bits = [k for (k,v) in bit_smarts.items() if len(v)>1]\ncollided_bits","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:30.924323Z","iopub.execute_input":"2024-06-05T04:42:30.924808Z","iopub.status.idle":"2024-06-05T04:42:30.934193Z","shell.execute_reply.started":"2024-06-05T04:42:30.924770Z","shell.execute_reply":"2024-06-05T04:42:30.932816Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"For example, look at bit 1855. Here we have what is essentially the same structure. These structures map to different smarts because of how RDKit stores bond information.","metadata":{}},{"cell_type":"code","source":"bit_idx = 1855\nsmarts = [get_morgan_smarts(mol, atom_idx, radius) \n          for (atom_idx, radius) in bit_map[bit_idx]]\n\nplot_bit_features(mol, bit_map, bit_idx, molsPerRow=3, subImgSize=(300,300), legends=smarts)","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:32.500661Z","iopub.execute_input":"2024-06-05T04:42:32.501089Z","iopub.status.idle":"2024-06-05T04:42:32.546706Z","shell.execute_reply.started":"2024-06-05T04:42:32.501055Z","shell.execute_reply":"2024-06-05T04:42:32.545421Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"If we look at the smarts themselves, we see they each pick up a slightly different aromatic bond orientation that gives them different smarts strings. You could probably fix this if you got really into the weeds on bond assignment.","metadata":{}},{"cell_type":"code","source":"Chem.Draw.MolsToGridImage([Chem.MolFromSmarts(i) for i in smarts], \n                          legends=smarts,\n                          molsPerRow=3,\n                          subImgSize=(300,300)\n                         )","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:34.301335Z","iopub.execute_input":"2024-06-05T04:42:34.302086Z","iopub.status.idle":"2024-06-05T04:42:34.332341Z","shell.execute_reply.started":"2024-06-05T04:42:34.302048Z","shell.execute_reply":"2024-06-05T04:42:34.331198Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We can fix this by converting the SMARTS to SMILES strings, which includes canonicalization. The results are not strictly SMARTS patterns, but they help resolve this duplication issue.\n\nAfter SMILES conversion, we no longer have any collided bits within the molecule\n\nNote that this canonicalization approach works for the subset used in this notebook, but you may have canonicalization issues on a larger corpus. If this occurs, you will need to find a different way of resolving the same structure mapping to different SMARTS.","metadata":{}},{"cell_type":"code","source":"bit_smarts = get_fp_smarts(mol, 3, 2048, canonicalize=True)\n\ncollided_bits = [k for (k,v) in bit_smarts.items() if len(v)>1]\ncollided_bits","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:36.070079Z","iopub.execute_input":"2024-06-05T04:42:36.070480Z","iopub.status.idle":"2024-06-05T04:42:36.162965Z","shell.execute_reply.started":"2024-06-05T04:42:36.070450Z","shell.execute_reply":"2024-06-05T04:42:36.161770Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"To find more real fingerprint collisions, we'll need to mine them from a corpus of molecules.\n\nHere we'll extract fingerprints from all molecules and look over the full set for collisions. We'll do this at three fingerprint sizes - 2048, 32000 and 102400 - to see how the collision rate evolves with fingerprint size.","metadata":{}},{"cell_type":"code","source":"def process_fp_vals(fp_dict):\n    fp_ints = []\n    fp_smarts = []\n    for k,v in fp_dict.items():\n        fp_ints.append(k)\n        fp_smarts.append(v)\n        \n    return fp_ints, fp_smarts\n\n\n\ndef compute_fp_vals(row):\n    smile = Chem.CanonSmiles(row['smiles'], useChiral=0)\n    mol = Chem.MolFromSmiles(smile)\n    \n    fp_smarts_small = get_fp_smarts(mol, 3, 2048, canonicalize=True)\n    fp_smarts_med = get_fp_smarts(mol, 3, 32000, canonicalize=True)\n    fp_smarts_large = get_fp_smarts(mol, 3, 102400, canonicalize=True)\n    \n    output = {}\n    names = ['small', 'med', 'large']\n    fp_dicts = [fp_smarts_small, fp_smarts_med, fp_smarts_large]\n    for i in range(len(fp_dicts)):\n        fp_ints, fp_smarts = process_fp_vals(fp_dicts[i])\n        output[f'{names[i]}_ints'] = fp_ints\n        output[f'{names[i]}_smarts'] = fp_smarts\n    \n    return output","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:37.483928Z","iopub.execute_input":"2024-06-05T04:42:37.484354Z","iopub.status.idle":"2024-06-05T04:42:37.498018Z","shell.execute_reply.started":"2024-06-05T04:42:37.484311Z","shell.execute_reply":"2024-06-05T04:42:37.496826Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"dataset = dataset.map(compute_fp_vals, num_proc=os.cpu_count())","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:42:38.997357Z","iopub.execute_input":"2024-06-05T04:42:38.997826Z","iopub.status.idle":"2024-06-05T04:46:20.369128Z","shell.execute_reply.started":"2024-06-05T04:42:38.997790Z","shell.execute_reply":"2024-06-05T04:46:20.367639Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Now we compile the results down to a single dict for each fingerprint size","metadata":{}},{"cell_type":"code","source":"def compile_fp_vals(batch):\n    names = ['small', 'med', 'large']\n    compiled_dicts = [defaultdict(set) for i in range(len(names))]\n    output = {}\n    \n    for i, name in enumerate(names):\n        compiled = compiled_dicts[i]\n        \n        fp_ints = batch[f'{name}_ints']\n        fp_smarts = batch[f'{name}_smarts']\n        \n        for j in range(len(fp_ints)):\n            fpi = fp_ints[j]\n            fps = fp_smarts[j]\n            \n            for k in range(len(fpi)):\n                compiled[fpi[k]].update(fps[k])\n                \n        compiled = [(k,v) for k,v in compiled.items()]\n        compiled = sorted(compiled, key=lambda x: x[0])\n        output[f'{name}_ints'] = [[i[0] for i in compiled]]\n        output[f'{name}_smarts'] = [[i[1] for i in compiled]]\n                \n    return output","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:20.371694Z","iopub.execute_input":"2024-06-05T04:46:20.372108Z","iopub.status.idle":"2024-06-05T04:46:20.384107Z","shell.execute_reply.started":"2024-06-05T04:46:20.372063Z","shell.execute_reply":"2024-06-05T04:46:20.382997Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"compiled_data = dataset.map(\n                    compile_fp_vals, \n                    batched=True, \n                    batch_size=500, \n                    remove_columns=dataset.column_names, num_proc=os.cpu_count())\n\ncompiled_data = compiled_data.map(compile_fp_vals, \n                                  batched=True, \n                                  batch_size=len(compiled_data)+1, \n                                  remove_columns=compiled_data.column_names)","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:20.387712Z","iopub.execute_input":"2024-06-05T04:46:20.388168Z","iopub.status.idle":"2024-06-05T04:46:24.427719Z","shell.execute_reply.started":"2024-06-05T04:46:20.388137Z","shell.execute_reply":"2024-06-05T04:46:24.426221Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"small_collision_dict = {k:v for (k,v) in \n                    zip(compiled_data[0]['small_ints'], compiled_data[0]['small_smarts'])}\n                        \nmed_collision_dict = {k:v for (k,v) in \n                    zip(compiled_data[0]['med_ints'], compiled_data[0]['med_smarts'])}\n\nlarge_collision_dict = {k:v for (k,v) in \n                    zip(compiled_data[0]['large_ints'], compiled_data[0]['large_smarts'])}","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:25.755723Z","iopub.execute_input":"2024-06-05T04:46:25.756858Z","iopub.status.idle":"2024-06-05T04:46:26.690518Z","shell.execute_reply.started":"2024-06-05T04:46:25.756803Z","shell.execute_reply":"2024-06-05T04:46:26.689255Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Now we can look at the distribution of how many collisions occur per bit.\n\nWe can see that collisions are common at 2048 bits, but very rare at 102400 bits","metadata":{}},{"cell_type":"code","source":"small_counts = np.array([len(i) for i in small_collision_dict.values()])\nmed_counts = np.array([len(i) for i in med_collision_dict.values()])\nlarge_counts = np.array([len(i) for i in large_collision_dict.values()])","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:28.498310Z","iopub.execute_input":"2024-06-05T04:46:28.498743Z","iopub.status.idle":"2024-06-05T04:46:28.510429Z","shell.execute_reply.started":"2024-06-05T04:46:28.498708Z","shell.execute_reply":"2024-06-05T04:46:28.508795Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"counts = [small_counts, med_counts, large_counts]\nlabels = ['2048 Bits', '32000 Bits', '102400 Bits']\n\nfor i in range(len(counts)):\n    print(f'{labels[i]} mean collision count: {counts[i].mean():.3f}')","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:29.202473Z","iopub.execute_input":"2024-06-05T04:46:29.202902Z","iopub.status.idle":"2024-06-05T04:46:29.209851Z","shell.execute_reply.started":"2024-06-05T04:46:29.202869Z","shell.execute_reply":"2024-06-05T04:46:29.208723Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.hist(small_counts, density=True, alpha=0.5, bins=20, label='2048 Bits')\nplt.hist(med_counts, density=True, alpha=0.5, bins=20, label='32000 Bits')\nplt.hist(large_counts, density=True, alpha=0.5, bins=20, label='102400 Bits')\nplt.legend()\nplt.xlabel('Collisions Per Bit')","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:30.683864Z","iopub.execute_input":"2024-06-05T04:46:30.684241Z","iopub.status.idle":"2024-06-05T04:46:31.253663Z","shell.execute_reply.started":"2024-06-05T04:46:30.684213Z","shell.execute_reply":"2024-06-05T04:46:31.252514Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def plot_bit_collision(collision_dict, bit_idx, **draw_kwargs):\n    bit_smarts = collision_dict[bit_idx]\n    bit_mols = [Chem.MolFromSmarts(i) for i in bit_smarts]\n    img = Chem.Draw.MolsToGridImage(bit_mols, legends=bit_smarts, **draw_kwargs)\n    return img","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:33.464262Z","iopub.execute_input":"2024-06-05T04:46:33.464704Z","iopub.status.idle":"2024-06-05T04:46:33.471470Z","shell.execute_reply.started":"2024-06-05T04:46:33.464671Z","shell.execute_reply":"2024-06-05T04:46:33.470041Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Here we can inspect a collided bit to see the structures that map to that bit. Looking at bit 53, we can see very different structures mapping to the same bit","metadata":{}},{"cell_type":"code","source":"plot_bit_collision(small_collision_dict, 53, molsPerRow=4, subImgSize=(300,300))","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:35.433623Z","iopub.execute_input":"2024-06-05T04:46:35.434102Z","iopub.status.idle":"2024-06-05T04:46:35.535062Z","shell.execute_reply.started":"2024-06-05T04:46:35.434068Z","shell.execute_reply":"2024-06-05T04:46:35.533614Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Looking at a collided bit in the large fingerprint, we can see there are still issues with near-identical structures being considered different. In this case, due to the presence/absence of hydrogens in the structure. This can be fixed with additional post-processing.\n\nThat said, we still have a genuine collision","metadata":{}},{"cell_type":"code","source":"plot_bit_collision(large_collision_dict, 78699, molsPerRow=4, subImgSize=(300,300))","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:37.296061Z","iopub.execute_input":"2024-06-05T04:46:37.296478Z","iopub.status.idle":"2024-06-05T04:46:37.388340Z","shell.execute_reply.started":"2024-06-05T04:46:37.296444Z","shell.execute_reply":"2024-06-05T04:46:37.387037Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"To look at this more systematically, we can calculate the tanimoto similarity for all SMARTS in a collided bit and see how different the largest discrepancy is ","metadata":{}},{"cell_type":"code","source":"def measure_bit_discrepancy(collision_dict, bit_idx, radius, fpSize):\n    bit_smarts = collision_dict[bit_idx]\n    bit_mols = [Chem.MolFromSmarts(i) for i in bit_smarts]\n    for mol in bit_mols:\n        mol.UpdatePropertyCache()\n        Chem.GetSymmSSSR(mol)\n    \n    fpg = rdFingerprintGenerator.GetMorganGenerator(radius=radius,fpSize=fpSize)\n    fps = [fpg.GetFingerprint(i) for i in bit_mols]\n    \n    tanimotos = np.array([DataStructs.BulkTanimotoSimilarity(i, fps) for i in fps])\n\n    return tanimotos.min()","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:40.185252Z","iopub.execute_input":"2024-06-05T04:46:40.186491Z","iopub.status.idle":"2024-06-05T04:46:40.194524Z","shell.execute_reply.started":"2024-06-05T04:46:40.186445Z","shell.execute_reply":"2024-06-05T04:46:40.193143Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"For bit 53 of the small fingerprint, we see a minimum tanimoto similarity of 0.02, indicating two very different molecules are colliding in that bit","metadata":{}},{"cell_type":"code","source":"measure_bit_discrepancy(small_collision_dict, 53, 3, 2048)","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:41.977551Z","iopub.execute_input":"2024-06-05T04:46:41.978024Z","iopub.status.idle":"2024-06-05T04:46:41.988627Z","shell.execute_reply.started":"2024-06-05T04:46:41.977988Z","shell.execute_reply":"2024-06-05T04:46:41.987346Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We also see a similarly low value for bit 78699 of the large fingerprint","metadata":{}},{"cell_type":"code","source":"measure_bit_discrepancy(large_collision_dict, 78699, 3, 2048)","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:42.862368Z","iopub.execute_input":"2024-06-05T04:46:42.862795Z","iopub.status.idle":"2024-06-05T04:46:42.875439Z","shell.execute_reply.started":"2024-06-05T04:46:42.862754Z","shell.execute_reply":"2024-06-05T04:46:42.874216Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"On the other hand, bit 95844 of the large fingerprint has a minimum similarity of 0.525, indicating most things within that bit are broadly similar","metadata":{}},{"cell_type":"code","source":"measure_bit_discrepancy(large_collision_dict, 95844, 3, 2048)","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:44.098280Z","iopub.execute_input":"2024-06-05T04:46:44.098744Z","iopub.status.idle":"2024-06-05T04:46:44.108959Z","shell.execute_reply.started":"2024-06-05T04:46:44.098707Z","shell.execute_reply":"2024-06-05T04:46:44.107508Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We can confirm this visually","metadata":{}},{"cell_type":"code","source":"plot_bit_collision(large_collision_dict, 95844, molsPerRow=4, subImgSize=(300,300))","metadata":{"execution":{"iopub.status.busy":"2024-06-05T04:46:44.842758Z","iopub.execute_input":"2024-06-05T04:46:44.843174Z","iopub.status.idle":"2024-06-05T04:46:44.922220Z","shell.execute_reply.started":"2024-06-05T04:46:44.843144Z","shell.execute_reply":"2024-06-05T04:46:44.920788Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Hopefully this notebook helps you visualize and analyze fingerprint collisions.","metadata":{}},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}