{"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":8462195,"sourceType":"datasetVersion","datasetId":5044537}],"dockerImageVersionId":30698,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"**paper**  \n- 'High-Quality Conformer Generation with CONFORGE: Algorithm and Performance Assessment'- T Seidel , acs 2023  \n- https://pubs.acs.org/doi/10.1021/acs.jcim.3c00563\n- https://github.com/molinfo-vienna/CDPKit\n\nThis one looks good, fast and accurate.","metadata":{}},{"cell_type":"code","source":"!pip install rdkit\n!pip install py3Dmol\n!pip install cdpkit\n","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2024-05-21T02:36:25.661792Z","iopub.execute_input":"2024-05-21T02:36:25.662155Z","iopub.status.idle":"2024-05-21T02:37:09.560811Z","shell.execute_reply.started":"2024-05-21T02:36:25.662126Z","shell.execute_reply":"2024-05-21T02:37:09.559617Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import sys\nimport CDPL.Chem as Chem\nimport CDPL.ConfGen as ConfGen\nimport py3Dmol\n\n# dictionary mapping status codes to human readable strings\nstatus_to_str = {\n    ConfGen.ReturnCode.UNINITIALIZED: 'uninitialized',\n    ConfGen.ReturnCode.TIMEOUT: 'max. processing time exceeded',\n    ConfGen.ReturnCode.ABORTED: 'aborted',\n    ConfGen.ReturnCode.FORCEFIELD_SETUP_FAILED: 'force field setup failed',\n    ConfGen.ReturnCode.FORCEFIELD_MINIMIZATION_FAILED: 'force field structure refinement failed',\n    ConfGen.ReturnCode.FRAGMENT_LIBRARY_NOT_SET: 'fragment library not available',\n    ConfGen.ReturnCode.FRAGMENT_CONF_GEN_FAILED: 'fragment conformer generation failed',\n    ConfGen.ReturnCode.FRAGMENT_CONF_GEN_TIMEOUT: 'fragment conformer generation timeout',\n    ConfGen.ReturnCode.FRAGMENT_ALREADY_PROCESSED: 'fragment already processed',\n    ConfGen.ReturnCode.TORSION_DRIVING_FAILED: 'torsion driving failed',\n    ConfGen.ReturnCode.CONF_GEN_FAILED: 'conformer generation failed'\n}\n\ndef gen3DStructure(mol: Chem.Molecule, struct_gen: ConfGen.StructureGenerator) -> int:\n    # prepare the molecule for 3D structure generation\n    ConfGen.prepareForConformerGeneration(mol)\n\n    # generate the 3D structure\n    status = struct_gen.generate(mol)\n\n    # if successful, store the generated conformer ensemble as\n    # per atom 3D coordinates arrays (= the way conformers are represented in CDPKit)\n\n    if status == ConfGen.ReturnCode.SUCCESS:\n        struct_gen.setCoordinates(mol)\n\n    # return status code\n    return status\n\nprint('import ok!!!')","metadata":{"execution":{"iopub.status.busy":"2024-05-21T02:37:48.206425Z","iopub.execute_input":"2024-05-21T02:37:48.206839Z","iopub.status.idle":"2024-05-21T02:37:48.316534Z","shell.execute_reply.started":"2024-05-21T02:37:48.206807Z","shell.execute_reply":"2024-05-21T02:37:48.315485Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"kaggle_smiles = [\n    'C#CCOc1ccc(CNc2nc(NCC3CCCN3c3cccnn3)nc(N[C@@H](CC#C)CC(=O)NC)n2)cc1',\n    'C#CCOc1ccc(CNc2nc(NCc3cccc(Br)n3)nc(N[C@@H](CC#C)CC(=O)NC)n2)cc1',\n    'C#CCOc1ccc(CNc2nc(Nc3cc(C(C)(C)C)[nH]n3)nc(N[C@@H](CC#C)CC(=O)NC)n2)cc1',\n    'C#CCOc1ccc(CNc2nc(Nc3nc(Cl)nc4cc(OC)c(OC)cc34)nc(N[C@@H](CC#C)CC(=O)NC)n2)cc1',\n    'C#CCOc1cccc(CNc2nc(NCc3nc4ccccc4s3)nc(N[C@@H](CC#C)CC(=O)NC)n2)c1',\n    'C#CC[C@@H](CC(=O)NC)Nc1nc(NCC23CCC(C)(CO2)C3)nc(Nc2cccc(S(=O)(=O)NC(C)(C)C)c2)n1',\n    'C#CC[C@@H](CC(=O)NC)Nc1nc(NCc2cccnc2N(C)C)nc(Nc2cnc(Cl)nc2OC)n1',\n    'C#CC[C@@](C)(Nc1nc(NCc2ccccc2N2CCOCC2)nc(NCc2ncccc2N(C)C)n1)C(=O)NC',\n    'CNC(=O)c1cc(OC)c(Nc2nc(NCCC3CC3)nc(NCc3nnc4c(=O)[nH]ccn34)n2)cc1N',\n    'CNC(=O)[C@H](CCCN=[N+]=[N-])Nc1nc(Nc2noc3ccc(F)cc23)nc(Nc2noc3ccc(F)cc23)n1'\n]\n\n# export only a single 3D structure (in case of multi-conf. input molecules)\nwriter = Chem.MolecularGraphWriter('xxx.sdf')\nChem.setMultiConfExportParameter(writer, False)\n\n# create and initialize an instance of the class ConfGen.StructureGenerator which will\n# perform the actual 3D structure generation work\nstruct_gen = ConfGen.StructureGenerator()\nstruct_gen.settings.timeout = 3600 * 1000  # apply the -t argument\n\n# create an instance of the default implementation of the Chem.Molecule interface\nmol = Chem.BasicMolecule()\n\nfor t,s in enumerate(kaggle_smiles):\n\n    # compose a simple molecule identifier\n    mol_id = f'xxx{t}'\n    mol = Chem.parseSMILES(s)\n\n    #print('- Generating 3D structure of molecule %s...' % mol_id) \n    try:\n\n        # generate 3D structure of the read molecule \n        status = gen3DStructure(mol, struct_gen)\n\n        # check for severe error reported by status code \n        if status != ConfGen.ReturnCode.SUCCESS: \n                print('Error: 3D structure generation for molecule %s failed: %s' % (\n                mol_id, status_to_str[status]))\n\n        else:\n\n            # enforce the output of 3D coordinates in case of MDL file formats\n            Chem.setMDLDimensionality(mol, 3)\n\n            # output the generated 3D structure\n            if not writer.write(mol):\n                sys.exit('Error: writing 3D structure of molecule %s failed' % mol_id)\n\n\n    except Exception as e:\n        sys.exit('Error: 3D structure generation or output for molecule %s failed: %s' % (mol_id, str(e)))\n\nwriter.close()\nprint('conformer ok!')\n","metadata":{"execution":{"iopub.status.busy":"2024-05-21T02:37:51.730945Z","iopub.execute_input":"2024-05-21T02:37:51.731296Z","iopub.status.idle":"2024-05-21T02:38:00.256753Z","shell.execute_reply.started":"2024-05-21T02:37:51.731269Z","shell.execute_reply":"2024-05-21T02:38:00.255770Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#visualization\nfrom rdkit import Chem as rdChem\nkaggle_mol = rdChem.SDMolSupplier('xxx.sdf')\n\n#print one for debug \n#print(kaggle_smiles[7])\n#print(kaggle_mol[7])\n\np = py3Dmol.view(width=1200, height=400, viewergrid=(1,3))\n#\nfor j,t in enumerate([7,0,5]):\n    print(f'view {j}:',kaggle_smiles[t])\n    p.removeAllModels(viewer=(0,j))\n    p.addModel(rdChem.MolToMolBlock(kaggle_mol[t], confId=0), 'sdf', viewer=(0,j))\n    p.setStyle({'stick':{}}, viewer=(0,j))\np.zoomTo()\np.show()\n     \n","metadata":{"execution":{"iopub.status.busy":"2024-05-21T02:38:39.536698Z","iopub.execute_input":"2024-05-21T02:38:39.537083Z","iopub.status.idle":"2024-05-21T02:38:39.745225Z","shell.execute_reply.started":"2024-05-21T02:38:39.537053Z","shell.execute_reply":"2024-05-21T02:38:39.744129Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## example application\n[paper] 'Spherical Message Passing for 3D Graph Networks' - Yi Liu, Arvix 2022\n\nhttps://arxiv.org/abs/2102.05013  \nhttps://github.com/divelab/DIG/  ","metadata":{}},{"cell_type":"code","source":"#!pip install dive-into-graphs","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from dig.threedgraph.method import SphereNet\nfrom torch_geometric.data import Data, Batch\nfrom rdkit import Chem\nfrom rdkit.Chem.rdchem import BondType\nimport numpy as np\nimport torch\n\nmodel = SphereNet(energy_and_force=False, cutoff=5.0, num_layers=4,\n                  hidden_channels=128, out_channels=1, int_emb_size=64,\n                  basis_emb_size_dist=8, basis_emb_size_angle=8, basis_emb_size_torsion=8, out_emb_channels=256,\n                  num_spherical=3, num_radial=6, envelope_exponent=5,\n                  num_before_skip=1, num_after_skip=2, num_output_layers=3)\n\n\n\nsdf_file = '<your sdf file from conformer generation>.sdf'\nmols = Chem.SDMolSupplier(sdf_file, removeHs=False, sanitize=False)\natom_type_list, position_list, con_mat_list = [], [], []\n\natomic_num_to_type={\n    1:0, 6:1, 7:2, 8:3, 9:4, 'unknown':8 #<todo> \n}\nbond_to_type = {\n    BondType.SINGLE: 1, BondType.DOUBLE: 2, BondType.TRIPLE: 3\n}\n\n#https://github.com/divelab/DIG/blob/21476b079c9226f38915dcd082b5c2ee0cddaac8/dig/threedgraph/dataset/PygQM93D.py#L11\ndata_list = []\nfor idx in range(100):#len(mols)\n\tprint('\\r',idx,end='',flush=True)\n\tmol = mols[idx]\n\tnum_atoms = mol.GetNumAtoms()\n\tposition = mols.GetItemText(idx).split('\\n')[4:4 + num_atoms]\n\tposition = np.array([[float(x) for x in line.split()[:3]] for line in position], dtype=np.float32)\n\tatom_type = np.array([atomic_num_to_type.get(atom.GetAtomicNum(),8) for atom in mol.GetAtoms()])\n\n    #bond edge not used: edge creation is done in SphereNet forward()\n    \n# \tcon_mat = np.zeros([num_atoms, num_atoms], dtype=int)\n# \tfor bond in mol.GetBonds():\n# \t\tstart, end = bond.GetBeginAtomIdx(), bond.GetEndAtomIdx()\n# \t\tbond_type = bond_to_type[bond.GetBondType()]\n# \t\tcon_mat[start, end] = bond_type\n# \t\tcon_mat[end, start] = bond_type\n\n# \tif atom_type[0] != 1:\n# \t\tcarbon_idxs = np.nonzero(atom_type == 1)\n# \t\tperm = np.arange(len(atom_type))\n# \t\tif len(carbon_idxs[0]) > 0:\n# \t\t\tfirst_carbon = carbon_idxs[0][0]\n# \t\t\tperm[0] = first_carbon\n# \t\t\tperm[first_carbon] = 0\n# \t\tatom_type, position = atom_type[perm], position[perm]\n# \t\tcon_mat = con_mat[perm][:, perm]\n \n\tdata = Data(\n\t\tpos=torch.tensor(position),\n\t\tz=torch.tensor(atom_type),\n\t)\n\tdata_list.append(data)\n\nbatch = Batch.from_data_list(data_list)\nbatch = batch.cuda()\nmodel = model.cuda()\nmodel = model.eval()\n\nwith torch.no_grad():\n\ty = model(batch)\nprint(y.shape)","metadata":{},"execution_count":null,"outputs":[]}]}