{"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"},{"sourceId":8064546,"sourceType":"datasetVersion","datasetId":4757656}],"dockerImageVersionId":30673,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"Hengck23 found some more data in the literature for sEH and asked if someone could advise on how to label the DNA attachment point as Dy so that it matches the Kaggle competition data. See discussion: https://www.kaggle.com/competitions/leash-BELKA/discussion/491908\n\nThis wasn't completely straightforward, so I decided it was easiest to create this notebook to show one way to approach this problem.\n\nOne important note: the structures in this dataset are quite different from those in the Kaggle dataset. This is because the Kaggle dataset was created from a library that enumerated a triazine core (a 6-membered aromatic ring with alternating carbons and nitrogens in the ring) while this dataset was enumerated from different cores.\n\nEDIT: Fixed the link.\nThis notebook uses a dataset from here: https://doi.org/10.26434/chemrxiv-2023-pq197","metadata":{}},{"cell_type":"code","source":"import numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nimport os","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2024-04-08T18:51:19.143278Z","iopub.execute_input":"2024-04-08T18:51:19.143681Z","iopub.status.idle":"2024-04-08T18:51:19.149972Z","shell.execute_reply.started":"2024-04-08T18:51:19.143652Z","shell.execute_reply":"2024-04-08T18:51:19.148434Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#This is the cheminformatics package that will do most of the heavy lifting for us\n!pip install rdkit","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:51:19.152685Z","iopub.execute_input":"2024-04-08T18:51:19.153409Z","iopub.status.idle":"2024-04-08T18:51:31.590026Z","shell.execute_reply.started":"2024-04-08T18:51:19.153370Z","shell.execute_reply":"2024-04-08T18:51:31.588885Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from rdkit import Chem\nfrom rdkit.Chem import Draw, AllChem\nfrom rdkit import RDLogger","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:51:31.591420Z","iopub.execute_input":"2024-04-08T18:51:31.591820Z","iopub.status.idle":"2024-04-08T18:51:31.599704Z","shell.execute_reply.started":"2024-04-08T18:51:31.591783Z","shell.execute_reply":"2024-04-08T18:51:31.598420Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import contextlib\nimport joblib\nfrom tqdm import tqdm\nfrom joblib import Parallel, delayed\n\n#Code from this nice Stack Overflow discussion: https://stackoverflow.com/questions/37804279/how-can-we-use-tqdm-in-a-parallel-execution-with-joblib\n@contextlib.contextmanager\ndef tqdm_joblib(tqdm_object):\n    \"\"\"Context manager to patch joblib to report into tqdm progress bar given as argument\"\"\"\n    class TqdmBatchCompletionCallback(joblib.parallel.BatchCompletionCallBack):\n        def __call__(self, *args, **kwargs):\n            tqdm_object.update(n=self.batch_size)\n            return super().__call__(*args, **kwargs)\n\n    old_batch_callback = joblib.parallel.BatchCompletionCallBack\n    joblib.parallel.BatchCompletionCallBack = TqdmBatchCompletionCallback\n    try:\n        yield tqdm_object\n    finally:\n        joblib.parallel.BatchCompletionCallBack = old_batch_callback\n        tqdm_object.close()","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:51:31.603158Z","iopub.execute_input":"2024-04-08T18:51:31.603593Z","iopub.status.idle":"2024-04-08T18:51:31.613661Z","shell.execute_reply.started":"2024-04-08T18:51:31.603559Z","shell.execute_reply":"2024-04-08T18:51:31.612368Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data = pd.read_csv('/kaggle/input/additional-data/total_compounds.csv')\n#data = data[0:100000]","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:51:31.615271Z","iopub.execute_input":"2024-04-08T18:51:31.616135Z","iopub.status.idle":"2024-04-08T18:51:46.659370Z","shell.execute_reply.started":"2024-04-08T18:51:31.616070Z","shell.execute_reply":"2024-04-08T18:51:46.658190Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data.head()","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:51:46.660679Z","iopub.execute_input":"2024-04-08T18:51:46.661067Z","iopub.status.idle":"2024-04-08T18:51:46.681230Z","shell.execute_reply.started":"2024-04-08T18:51:46.661035Z","shell.execute_reply":"2024-04-08T18:51:46.678932Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def replace_DNA_attachment_point(structure_SMILES: str, bb_SMILES: str) -> str:\n    \"\"\"\n    This function takes the final structure and the first building block as SMILES strings, then does some substructure matching and replacement to substitute the N-Me where\n    the DNA is attached for [Dy].\n    \n    This function works when the first building block in the synthesis is an acid that gets conjugated to the DNA linker or an amine that gets conjugated to the DNA linker.\n    \n    It is assumed that if there is an acid present in the first building block, that this is the site of attachment to the DNA linker.\n    \n    Args:\n        structure_SMILES: The SMILES string for the final molecule.\n    \n        bb_SMILES: The SMILES string for the first building block (i.e. the one that should have the DNA attachment point).\n    \n    Returns: \n        \n        A SMILES string representing the final molecule where the DNA attachment point is labeled with Dy.\n\n    \"\"\"\n    #print(structure_SMILES)\n    #print(bb_SMILES)\n    #Get the molecule object\n    mol = Chem.MolFromSmiles(structure_SMILES)\n\n    #Get the building block molecule object\n    bb_mol = Chem.MolFromSmiles(bb_SMILES)\n    \n    #Define a molecule object from a SMARTS pattern matching a carboxylic acid\n    acid_pattern = Chem.MolFromSmarts('C(=O)[O;h1]')\n    \n    #Just in case RDKit had trouble reading any of the SMILES\n    try:\n        if bb_mol.HasSubstructMatch(acid_pattern):\n            #Enumerate the methyl group onto the building block\n            new_bb = AllChem.ReplaceSubstructs(bb_mol, acid_pattern, Chem.MolFromSmarts('C(=O)[N]-[C;h3]'))[0]           \n            \n            #Sanity checkup\n            Chem.SanitizeMol(new_bb)\n            \n            #Convert to string\n            bb_string = Chem.MolToSmiles(new_bb)\n            \n            #Define a recursive SMARTS that checks for an N-Methyl group coming off of the new_building block.\n            string = f\"[N]-[$([C;h3]);$({bb_string})]\"\n            \n            #Get a molecule object from the SMARTS\n            bb_smarts = Chem.MolFromSmarts(string)\n            \n            #Replace the identified N-Me with N-Dy\n            final_mol = AllChem.ReplaceSubstructs(mol, Chem.MolFromSmarts(string), Chem.MolFromSmarts('[N]-[Dy]'))[0]\n            \n            Chem.SanitizeMol(final_mol)\n\n            #Return the string\n            return Chem.MolToSmiles(final_mol)\n            \n            \n        else:\n            #Enumerate the methyl group onto the building block\n            new_bb = AllChem.ReplaceSubstructs(bb_mol, Chem.MolFromSmarts('[N;h2]'), Chem.MolFromSmarts('[N]-[C;h3]'))[0]\n        \n            #Sanity checkup\n            Chem.SanitizeMol(new_bb)\n        \n            #Convert to string\n            bb_string = Chem.MolToSmiles(new_bb)\n        \n            #Define a recursive SMARTS that checks for an N-Methyl group coming off of the new_building block.\n            string = f\"[N;h0]-[$([C;h3]);$({bb_string})]\"\n        \n            #Get a molecule object from the SMARTS\n            bb_smarts = Chem.MolFromSmarts(string)\n        \n            #Replace the identified N-Me with N-Dy\n            final_mol = AllChem.ReplaceSubstructs(mol, Chem.MolFromSmarts(string), Chem.MolFromSmarts('[N]-[Dy]'))[0]\n            \n            Chem.SanitizeMol(final_mol)\n\n            #Return the string\n            return Chem.MolToSmiles(final_mol)\n    \n    except:\n        \n        return None","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:51:46.683887Z","iopub.execute_input":"2024-04-08T18:51:46.684569Z","iopub.status.idle":"2024-04-08T18:51:46.697781Z","shell.execute_reply.started":"2024-04-08T18:51:46.684531Z","shell.execute_reply":"2024-04-08T18:51:46.696319Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Get the structures of the molecules and first building blocks as lists, then zip them\nstructures = data.structure.tolist()\n\nbbs = data.bb1.tolist()\n\nzipped_structures = list(zip(structures, bbs))","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:51:46.699763Z","iopub.execute_input":"2024-04-08T18:51:46.700199Z","iopub.status.idle":"2024-04-08T18:51:46.743279Z","shell.execute_reply.started":"2024-04-08T18:51:46.700154Z","shell.execute_reply":"2024-04-08T18:51:46.741861Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Let's see how many CPUs are available to us\nprint(f\"This notebook has {os.cpu_count()} CPUs available\")","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:51:46.744549Z","iopub.execute_input":"2024-04-08T18:51:46.745155Z","iopub.status.idle":"2024-04-08T18:51:46.764377Z","shell.execute_reply.started":"2024-04-08T18:51:46.745100Z","shell.execute_reply":"2024-04-08T18:51:46.763188Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Paralellizing our DNA attachment point replacement code across all of the available CPUs\nwith tqdm_joblib(tqdm(desc=\"Labeling the DNA attachment point\", total=len(zipped_structures))) as progress_bar:\n    values = Parallel(n_jobs=-1)(delayed(replace_DNA_attachment_point)(x,y) for x,y in zipped_structures)","metadata":{"_kg_hide-output":true,"execution":{"iopub.status.busy":"2024-04-08T18:51:46.767618Z","iopub.execute_input":"2024-04-08T18:51:46.768702Z","iopub.status.idle":"2024-04-08T18:52:37.401423Z","shell.execute_reply.started":"2024-04-08T18:51:46.768640Z","shell.execute_reply":"2024-04-08T18:52:37.400438Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Update the dataframe with the new values\ndata['new_structure'] = values","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:52:37.402389Z","iopub.execute_input":"2024-04-08T18:52:37.402699Z","iopub.status.idle":"2024-04-08T18:52:37.417755Z","shell.execute_reply.started":"2024-04-08T18:52:37.402672Z","shell.execute_reply":"2024-04-08T18:52:37.416444Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Check out the updated dataframe\ndata.head()","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:52:37.419932Z","iopub.execute_input":"2024-04-08T18:52:37.420544Z","iopub.status.idle":"2024-04-08T18:52:37.442033Z","shell.execute_reply.started":"2024-04-08T18:52:37.420508Z","shell.execute_reply":"2024-04-08T18:52:37.440819Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"###Visualize some of the original data\nDraw.MolsToGridImage([Chem.MolFromSmiles(x) for x in data['structure'][0:4]], molsPerRow=4, subImgSize=(400,300))","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:52:37.443485Z","iopub.execute_input":"2024-04-08T18:52:37.444656Z","iopub.status.idle":"2024-04-08T18:52:37.503437Z","shell.execute_reply.started":"2024-04-08T18:52:37.444612Z","shell.execute_reply":"2024-04-08T18:52:37.502173Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"###Visualize some of the newly labeled data\nDraw.MolsToGridImage([Chem.MolFromSmiles(x) for x in data['new_structure'][0:4]], molsPerRow=4, subImgSize=(400,300))","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:52:37.505421Z","iopub.execute_input":"2024-04-08T18:52:37.506144Z","iopub.status.idle":"2024-04-08T18:52:37.568799Z","shell.execute_reply.started":"2024-04-08T18:52:37.506100Z","shell.execute_reply":"2024-04-08T18:52:37.567575Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#check to see if we have any NaN\nprint(f\"Number of structures that did not convert: {data['new_structure'].isna().sum()}\")","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:52:37.570451Z","iopub.execute_input":"2024-04-08T18:52:37.571097Z","iopub.status.idle":"2024-04-08T18:52:37.586714Z","shell.execute_reply.started":"2024-04-08T18:52:37.571055Z","shell.execute_reply":"2024-04-08T18:52:37.585056Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Write the file to .csv format\ndata.to_csv('DNA_Labeled_Data.csv', index=False)","metadata":{"execution":{"iopub.status.busy":"2024-04-08T18:52:37.588001Z","iopub.execute_input":"2024-04-08T18:52:37.588396Z","iopub.status.idle":"2024-04-08T18:52:38.540713Z","shell.execute_reply.started":"2024-04-08T18:52:37.588355Z","shell.execute_reply":"2024-04-08T18:52:38.539783Z"},"trusted":true},"execution_count":null,"outputs":[]}]}