1 #!/bin/env python 2 # 3 # File: RDKitGenerateMolecularFrameworks.py 4 # Author: Manish Sud <msud@san.rr.com> 5 # 6 # Copyright (C) 2026 Manish Sud. All rights reserved. 7 # 8 # The functionality available in this script is implemented using RDKit, an 9 # open source toolkit for cheminformatics developed by Greg Landrum. 10 # 11 # This file is part of MayaChemTools. 12 # 13 # MayaChemTools is free software; you can redistribute it and/or modify it under 14 # the terms of the GNU Lesser General Public License as published by the Free 15 # Software Foundation; either version 3 of the License, or (at your option) any 16 # later version. 17 # 18 # MayaChemTools is distributed in the hope that it will be useful, but without 19 # any warranty; without even the implied warranty of merchantability of fitness 20 # for a particular purpose. See the GNU Lesser General Public License for more 21 # details. 22 # 23 # You should have received a copy of the GNU Lesser General Public License 24 # along with MayaChemTools; if not, see <http://www.gnu.org/licenses/> or 25 # write to the Free Software Foundation Inc., 59 Temple Place, Suite 330, 26 # Boston, MA, 02111-1307, USA. 27 # 28 29 from __future__ import print_function 30 31 import os 32 import sys 33 import time 34 import re 35 36 # RDKit imports... 37 try: 38 from rdkit import rdBase 39 from rdkit import Chem 40 from rdkit.Chem import AllChem 41 from rdkit.Chem.Scaffolds import MurckoScaffold 42 except ImportError as ErrMsg: 43 sys.stderr.write("\nFailed to import RDKit module/package: %s\n" % ErrMsg) 44 sys.stderr.write("Check/update your RDKit environment and try again.\n\n") 45 sys.exit(1) 46 47 # MayaChemTools imports... 48 sys.path.insert(0, os.path.join(os.path.dirname(sys.argv[0]), "..", "lib", "Python")) 49 try: 50 from docopt import docopt 51 import MiscUtil 52 import RDKitUtil 53 except ImportError as ErrMsg: 54 sys.stderr.write("\nFailed to import MayaChemTools module/package: %s\n" % ErrMsg) 55 sys.stderr.write("Check/update your MayaChemTools environment and try again.\n\n") 56 sys.exit(1) 57 58 ScriptName = os.path.basename(sys.argv[0]) 59 Options = {} 60 OptionsInfo = {} 61 62 63 def main(): 64 """Start execution of the script.""" 65 66 MiscUtil.PrintInfo( 67 "\n%s (RDKit v%s; MayaChemTools v%s; %s): Starting...\n" 68 % (ScriptName, rdBase.rdkitVersion, MiscUtil.GetMayaChemToolsVersion(), time.asctime()) 69 ) 70 71 (WallClockTime, ProcessorTime) = MiscUtil.GetWallClockAndProcessorTime() 72 73 # Retrieve command line arguments and options... 74 RetrieveOptions() 75 76 # Process and validate command line arguments and options... 77 ProcessOptions() 78 79 # Perform actions required by the script... 80 GenerateMolecularFrameworks() 81 82 MiscUtil.PrintInfo("\n%s: Done...\n" % ScriptName) 83 MiscUtil.PrintInfo("Total time: %s" % MiscUtil.GetFormattedElapsedTime(WallClockTime, ProcessorTime)) 84 85 86 def GenerateMolecularFrameworks(): 87 """Generate Bemis Murcko molecular framworks.""" 88 89 Infile = OptionsInfo["Infile"] 90 Outfile = OptionsInfo["Outfile"] 91 92 UseChirality = OptionsInfo["UseChirality"] 93 94 RemoveDuplicateFrameworks = OptionsInfo["RemoveDuplicateFrameworks"] 95 UseGraphFrameworks = OptionsInfo["UseGraphFrameworks"] 96 97 SortFrameworks = OptionsInfo["SortFrameworks"] 98 if SortFrameworks: 99 FrameworkMolIDs = [] 100 FrameworkMolIDToMolMap = {} 101 FrameworkMolIDToAtomCountMap = {} 102 103 DuplicateFrameworkMolIDs = [] 104 DuplicateFrameworkMolIDToMolMap = {} 105 DuplicateFrameworkMolIDToAtomCountMap = {} 106 107 DuplicatesOutfile = "" 108 if RemoveDuplicateFrameworks: 109 DuplicatesOutfile = OptionsInfo["DuplicatesOutfile"] 110 111 # Setup a molecule reader... 112 MiscUtil.PrintInfo("\nProcessing file %s..." % Infile) 113 Mols = RDKitUtil.ReadMolecules(Infile, **OptionsInfo["InfileParams"]) 114 115 # Set up a molecular framework writer... 116 Writer = RDKitUtil.MoleculesWriter(Outfile, **OptionsInfo["OutfileParams"]) 117 if Writer is None: 118 MiscUtil.PrintError("Failed to setup a writer for output fie %s " % Outfile) 119 120 # Set up a duplicate molecular framework writer... 121 if RemoveDuplicateFrameworks: 122 DuplicatesWriter = RDKitUtil.MoleculesWriter(DuplicatesOutfile, **OptionsInfo["OutfileParams"]) 123 if Writer is None: 124 MiscUtil.PrintError("Failed to setup a writer for duplicates output fie %s " % DuplicatesOutfile) 125 126 if RemoveDuplicateFrameworks: 127 MiscUtil.PrintInfo("Generating files: %s and %s..." % (Outfile, DuplicatesOutfile)) 128 else: 129 MiscUtil.PrintInfo("Generating file %s..." % Outfile) 130 131 # Process molecules... 132 MolCount = 0 133 ValidMolCount = 0 134 135 FrameworksCount = 0 136 UniqueFrameworksCount = 0 137 DuplicateFrameworksCount = 0 138 139 CanonicalSMILESMap = {} 140 141 Compute2DCoords = OptionsInfo["OutfileParams"]["Compute2DCoords"] 142 143 for Mol in Mols: 144 MolCount += 1 145 146 if Mol is None: 147 continue 148 149 if RDKitUtil.IsMolEmpty(Mol): 150 MolName = RDKitUtil.GetMolName(Mol, MolCount) 151 MiscUtil.PrintWarning("Ignoring empty molecule: %s" % MolName) 152 continue 153 154 ValidMolCount += 1 155 156 if UseGraphFrameworks: 157 FrameworksMol = MurckoScaffold.MakeScaffoldGeneric(Mol) 158 else: 159 FrameworksMol = MurckoScaffold.GetScaffoldForMol(Mol) 160 161 if Compute2DCoords: 162 AllChem.Compute2DCoords(FrameworksMol) 163 164 if SortFrameworks: 165 HeavyAtomCount = FrameworksMol.GetNumHeavyAtoms() 166 167 FrameworksCount += 1 168 169 if RemoveDuplicateFrameworks: 170 CanonicalSMILES = Chem.MolToSmiles(FrameworksMol, isomericSmiles=UseChirality, canonical=True) 171 if CanonicalSMILES in CanonicalSMILESMap: 172 DuplicateFrameworksCount += 1 173 if SortFrameworks: 174 # Track duplicate frameworks... 175 DuplicateFrameworkMolIDs.append(DuplicateFrameworksCount) 176 DuplicateFrameworkMolIDToMolMap[DuplicateFrameworksCount] = FrameworksMol 177 DuplicateFrameworkMolIDToAtomCountMap[DuplicateFrameworksCount] = HeavyAtomCount 178 else: 179 # Write it out... 180 DuplicatesWriter.write(FrameworksMol) 181 else: 182 UniqueFrameworksCount += 1 183 CanonicalSMILESMap[CanonicalSMILES] = CanonicalSMILES 184 if SortFrameworks: 185 # Track unique frameworks... 186 FrameworkMolIDs.append(UniqueFrameworksCount) 187 FrameworkMolIDToMolMap[UniqueFrameworksCount] = FrameworksMol 188 FrameworkMolIDToAtomCountMap[UniqueFrameworksCount] = HeavyAtomCount 189 else: 190 # Write it out... 191 Writer.write(FrameworksMol) 192 elif SortFrameworks: 193 # Track for sorting... 194 FrameworkMolIDs.append(FrameworksCount) 195 FrameworkMolIDToMolMap[FrameworksCount] = FrameworksMol 196 FrameworkMolIDToAtomCountMap[FrameworksCount] = HeavyAtomCount 197 else: 198 # Write it out... 199 Writer.write(FrameworksMol) 200 201 if SortFrameworks: 202 ReverseOrder = OptionsInfo["DescendingSortOrder"] 203 SortAndWriteFrameworks( 204 Writer, FrameworkMolIDs, FrameworkMolIDToMolMap, FrameworkMolIDToAtomCountMap, ReverseOrder 205 ) 206 if RemoveDuplicateFrameworks: 207 SortAndWriteFrameworks( 208 DuplicatesWriter, 209 DuplicateFrameworkMolIDs, 210 DuplicateFrameworkMolIDToMolMap, 211 DuplicateFrameworkMolIDToAtomCountMap, 212 ReverseOrder, 213 ) 214 215 Writer.close() 216 if RemoveDuplicateFrameworks: 217 DuplicatesWriter.close() 218 219 MiscUtil.PrintInfo("\nTotal number of molecular frameworks: %d" % FrameworksCount) 220 if RemoveDuplicateFrameworks: 221 MiscUtil.PrintInfo("Number of unique molecular frameworks: %d" % UniqueFrameworksCount) 222 MiscUtil.PrintInfo("Number of duplicate molecular frameworks: %d" % DuplicateFrameworksCount) 223 224 MiscUtil.PrintInfo("\nTotal number of molecules: %d" % MolCount) 225 MiscUtil.PrintInfo("Number of valid molecules: %d" % ValidMolCount) 226 MiscUtil.PrintInfo("Number of ignored molecules: %d" % (MolCount - ValidMolCount)) 227 228 229 def SortAndWriteFrameworks(MolWriter, MolIDs, MolIDToMolMap, MolIDToAtomCountMap, ReverseOrder): 230 """Sort frameworks and write them out.""" 231 SortedMolIDs = sorted(MolIDs, key=lambda MolID: MolIDToAtomCountMap[MolID], reverse=ReverseOrder) 232 for MolID in SortedMolIDs: 233 FrameworksMol = MolIDToMolMap[MolID] 234 MolWriter.write(FrameworksMol) 235 236 237 def ProcessOptions(): 238 """Process and validate command line arguments and options.""" 239 240 MiscUtil.PrintInfo("Processing options...") 241 242 # Validate options... 243 ValidateOptions() 244 245 OptionsInfo["Infile"] = Options["--infile"] 246 OptionsInfo["InfileParams"] = MiscUtil.ProcessOptionInfileParameters( 247 "--infileParams", Options["--infileParams"], Options["--infile"] 248 ) 249 250 OptionsInfo["Outfile"] = Options["--outfile"] 251 OptionsInfo["OutfileParams"] = MiscUtil.ProcessOptionOutfileParameters( 252 "--outfileParams", Options["--outfileParams"], Options["--infile"], Options["--outfile"] 253 ) 254 255 OptionsInfo["Overwrite"] = Options["--overwrite"] 256 257 OptionsInfo["Mode"] = Options["--mode"] 258 OptionsInfo["UseGraphFrameworks"] = False 259 if re.match("^GraphFrameworks$", OptionsInfo["Mode"], re.I): 260 OptionsInfo["UseGraphFrameworks"] = True 261 262 OptionsInfo["RemoveDuplicates"] = Options["--removeDuplicates"] 263 OptionsInfo["RemoveDuplicateFrameworks"] = False 264 if re.match("^Yes$", OptionsInfo["RemoveDuplicates"], re.I): 265 OptionsInfo["RemoveDuplicateFrameworks"] = True 266 267 # Setup outfile for writing out duplicates... 268 OptionsInfo["DuplicatesOutfile"] = "" 269 if OptionsInfo["RemoveDuplicateFrameworks"]: 270 FileDir, FileName, FileExt = MiscUtil.ParseFileName(OptionsInfo["Outfile"]) 271 OptionsInfo["DuplicatesOutfile"] = "%sDuplicates.%s" % (FileName, FileExt) 272 273 OptionsInfo["Sort"] = Options["--sort"] 274 OptionsInfo["SortFrameworks"] = False 275 if re.match("^Yes$", OptionsInfo["Sort"], re.I): 276 OptionsInfo["SortFrameworks"] = True 277 278 OptionsInfo["SortOrder"] = Options["--sortOrder"] 279 OptionsInfo["DescendingSortOrder"] = False 280 if re.match("^Descending$", OptionsInfo["SortOrder"], re.I): 281 OptionsInfo["DescendingSortOrder"] = True 282 283 OptionsInfo["UseChirality"] = False 284 if re.match("^yes$", Options["--useChirality"], re.I): 285 OptionsInfo["UseChirality"] = True 286 287 288 def RetrieveOptions(): 289 """Retrieve command line arguments and options.""" 290 291 # Get options... 292 global Options 293 Options = docopt(_docoptUsage_) 294 295 # Set current working directory to the specified directory... 296 WorkingDir = Options["--workingdir"] 297 if WorkingDir: 298 os.chdir(WorkingDir) 299 300 # Handle examples option... 301 if "--examples" in Options and Options["--examples"]: 302 MiscUtil.PrintInfo(MiscUtil.GetExamplesTextFromDocOptText(_docoptUsage_)) 303 sys.exit(0) 304 305 306 def ValidateOptions(): 307 """Validate option values.""" 308 309 MiscUtil.ValidateOptionFilePath("-i, --infile", Options["--infile"]) 310 MiscUtil.ValidateOptionFileExt("-i, --infile", Options["--infile"], "sdf sd smi txt csv tsv") 311 312 MiscUtil.ValidateOptionFileExt("-o, --outfile", Options["--outfile"], "sdf sd smi") 313 MiscUtil.ValidateOptionsOutputFileOverwrite( 314 "-o, --outfile", Options["--outfile"], "--overwrite", Options["--overwrite"] 315 ) 316 MiscUtil.ValidateOptionsDistinctFileNames( 317 "-i, --infile", Options["--infile"], "-o, --outfile", Options["--outfile"] 318 ) 319 320 MiscUtil.ValidateOptionTextValue("-m, --mode", Options["--mode"], "GraphFrameworks AtomicFrameworks") 321 322 MiscUtil.ValidateOptionTextValue("-r, --removeDuplicates", Options["--removeDuplicates"], "yes no") 323 MiscUtil.ValidateOptionTextValue("-s, --sort", Options["--sort"], "yes no") 324 MiscUtil.ValidateOptionTextValue("--sortOrder", Options["--sortOrder"], "ascending descending") 325 326 MiscUtil.ValidateOptionTextValue("--useChirality", Options["--useChirality"], "yes no") 327 328 329 # Setup a usage string for docopt... 330 _docoptUsage_ = """ 331 RDKitGenerateMolecularFrameworks.py - Generate Bemis Murcko molecular frameworks 332 333 Usage: 334 RDKitGenerateMolecularFrameworks.py [--infileParams <Name,Value,...>] 335 [--mode <GraphFrameworks or AtomicFrameworks> ] 336 [ --outfileParams <Name,Value,...> ] [--overwrite] [--removeDuplicates <yes or no>] 337 [--sort <yes or no>] [--sortOrder <ascending or descending>] 338 [--useChirality <yes or no>] [-w <dir>] -i <infile> -o <outfile> 339 RDKitGenerateMolecularFrameworks.py -h | --help | -e | --examples 340 341 Description: 342 Generate Bemis Murcko [ Ref 133 ] molecular frameworks for molecules. Two types of molecular 343 frameworks can be generated: Graph or atomic frameworks. The graph molecular framework 344 is a generic framework. The atom type, hybridization, and bond order is ignore during its 345 generation. All atoms are set to carbon atoms and all bonds are single bonds. The atom type, 346 hybridization, and bond order is preserved during generation of atomic molecular frameworks. 347 348 The supported input file formats are: SD (.sdf, .sd), SMILES (.smi, .csv, .tsv, .txt) 349 350 The supported output file formats are: SD (.sdf, .sd), SMILES (.smi) 351 352 Options: 353 -e, --examples 354 Print examples. 355 -h, --help 356 Print this help message. 357 -i, --infile <infile> 358 Input file name. 359 --infileParams <Name,Value,...> [default: auto] 360 A comma delimited list of parameter name and value pairs for reading 361 molecules from files. The supported parameter names for different file 362 formats, along with their default values, are shown below: 363 364 SD: removeHydrogens,yes,sanitize,yes,strictParsing,yes 365 SMILES: smilesColumn,1,smilesNameColumn,2,smilesDelimiter,space, 366 smilesTitleLine,auto,sanitize,yes 367 368 Possible values for smilesDelimiter: space, comma or tab. 369 -m, --mode <GraphFrameworks or AtomicFrameworks> [default: GraphFrameworks] 370 Type of molecular frameworks to generate for molecules. Possible values: GraphFrameworks 371 or AtomicFrameworks. The graph molecular framework is a generic framework. The atom type, 372 hybridization, and bond order is ignore during its generation. All atoms are set to carbon atoms 373 and all bonds are single bonds. The atom type, hybridization, and bond order is preserved 374 during the generation of atomic molecular frameworks. 375 -o, --outfile <outfile> 376 Output file name. 377 --outfileParams <Name,Value,...> [default: auto] 378 A comma delimited list of parameter name and value pairs for writing 379 molecules to files. The supported parameter names for different file 380 formats, along with their default values, are shown below: 381 382 SD: compute2DCoords,auto,forceV3000,no 383 SMILES: smilesDelimiter,space,smilesTitleLine,yes 384 385 Default value for compute2DCoords: yes for SMILES input file; no for all other 386 file types. 387 --overwrite 388 Overwrite existing files. 389 -r, --removeDuplicates <yes or no> [default: no] 390 Remove duplicate molecular frameworks. Possible values: yes or no. The duplicate 391 molecular franworks are identified using canonical SMILES. The removed frameworks 392 are written to a separate output file. 393 -s, --sort <yes or no> [default: no] 394 Sort molecular frameworks by heavy atom count. Possible values: yes or no. 395 --sortOrder <ascending or descending> [default: ascending] 396 Sorting order for molecular frameworks. Possible values: ascending or descending. 397 -u, --useChirality <yes or no> [default: yes] 398 Use stereochemistry for generation of canonical SMILES strings to identify 399 duplicate molecular frameworks. 400 -w, --workingdir <dir> 401 Location of working directory which defaults to the current directory. 402 403 Examples: 404 To generate graph molecular framworks for molecules and write out a SMILES file, 405 type: 406 407 % RDKitGenerateMolecularFrameworks.py -i Sample.smi -o SampleOut.smi 408 409 To generate graph molecular framworks, remove duplicate frameworks for molecules 410 and write out SD files for unique and duplicate frameworks, type: 411 412 % RDKitGenerateMolecularFrameworks.py -m GraphFrameworks -r yes 413 -i Sample.sdf -o SampleOut.sdf 414 415 To generate atomic molecular framworks, remove duplicate frameworks, sort 416 framworks by heavy atom count in ascending order, write out SMILES files for 417 unique and duplicate frameworks, type: 418 419 % RDKitGenerateMolecularFrameworks.py -m AtomicFrameworks -r yes 420 -s yes -i Sample.smi -o SampleOut.smi 421 422 To generate graph molecular framworks for molecules in a CSV SMILES file, 423 SMILES strings in column 1, name in olumn 2, emove duplicate frameworks, 424 sort framworks by heavy atom count in decending order and write out a SD 425 file, type: 426 427 % RDKitGenerateMolecularFrameworks.py -m AtomicFrameworks 428 --removeDuplicates yes -s yes --sortOrder descending --infileParams 429 "smilesDelimiter,comma,smilesTitleLine,yes,smilesColumn,1, 430 smilesNameColumn,2" --outfileParams "compute2DCoords,yes" 431 -i SampleSMILES.csv -o SampleOut.sdf 432 433 Author: 434 Manish Sud(msud@san.rr.com) 435 436 See also: 437 RDKitConvertFileFormat.py, RDKitDrawMolecules.py, RDKitSearchFunctionalGroups.py, 438 RDKitSearchSMARTS.py 439 440 Copyright: 441 Copyright (C) 2026 Manish Sud. All rights reserved. 442 443 The functionality available in this script is implemented using RDKit, an 444 open source toolkit for cheminformatics developed by Greg Landrum. 445 446 This file is part of MayaChemTools. 447 448 MayaChemTools is free software; you can redistribute it and/or modify it under 449 the terms of the GNU Lesser General Public License as published by the Free 450 Software Foundation; either version 3 of the License, or (at your option) any 451 later version. 452 453 """ 454 455 if __name__ == "__main__": 456 main()
