Mercurial > repos > bgruening > rdconf

comparison rdkit_descriptors.py @ 0:5c501eb8d56c draft default tip

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
"planemo upload for repository https://github.com/bgruening/galaxytools/tree/master/chemicaltoolbox/rdkit commit c1d813d3f0fec60ea6efe8a11e59d98bfdc1636f"
author bgruening
date Sat, 04 Dec 2021 16:39:31 +0000
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:5c501eb8d56c
1 #!/usr/bin/env python
2
3 import argparse
4 import inspect
5 import sys
6
7 from rdkit import Chem
8 from rdkit.Chem import Descriptors
9
10
11 def get_supplier(infile, format="smiles"):
12 """
13 Returns a generator over a SMILES or InChI file. Every element is of RDKit
14 molecule and has its original string as _Name property.
15 """
16 with open(infile) as handle:
17 for line in handle:
18 line = line.strip()
19 if format == "smiles":
20 mol = Chem.MolFromSmiles(line, sanitize=True)
21 elif format == "inchi":
22 mol = Chem.inchi.MolFromInchi(
23 line,
24 sanitize=True,
25 removeHs=True,
26 logLevel=None,
27 treatWarningAsError=False,
28 )
29 if mol is None:
30 yield False
31 else:
32 mol.SetProp("_Name", line.split("\t")[0])
33 yield mol
34
35
36 def get_rdkit_descriptor_functions():
37 """
38 Returns all descriptor functions under the Chem.Descriptors Module as tuple of (name, function)
39 """
40 ret = [
41 (name, f)
42 for name, f in inspect.getmembers(Descriptors)
43 if inspect.isfunction(f) and not name.startswith("_")
44 ]
45 # some which are not in the official Descriptors module we need to add manually
46 ret.extend([("FormalCharge", Chem.GetFormalCharge), ("SSSR", Chem.GetSSSR)])
47 ret.sort()
48 return ret
49
50
51 def descriptors(mol, functions):
52 """
53 Calculates the descriptors of a given molecule.
54 """
55 for name, function in functions:
56 yield (name, function(mol))
57
58
59 if __name__ == "__main__":
60 parser = argparse.ArgumentParser()
61 parser.add_argument("-i", "--infile", required=True, help="Path to the input file.")
62 parser.add_argument("--iformat", help="Specify the input file format.")
63
64 parser.add_argument(
65 "-o",
66 "--outfile",
67 type=argparse.FileType("w+"),
68 default=sys.stdout,
69 help="path to the result file, default is stdout",
70 )
71
72 parser.add_argument(
73 "-s",
74 "--select",
75 default=None,
76 help="select a subset of comma-separated descriptors to use",
77 )
78
79 parser.add_argument(
80 "--header",
81 dest="header",
82 action="store_true",
83 default=False,
84 help="Write header line.",
85 )
86
87 args = parser.parse_args()
88
89 if args.iformat == "sdf":
90 supplier = Chem.SDMolSupplier(args.infile)
91 elif args.iformat == "smi":
92 supplier = get_supplier(args.infile, format="smiles")
93 elif args.iformat == "inchi":
94 supplier = get_supplier(args.infile, format="inchi")
95 elif args.iformat == "pdb":
96 supplier = [Chem.MolFromPDBFile(args.infile)]
97 elif args.iformat == "mol2":
98 supplier = [Chem.MolFromMol2File(args.infile)]
99
100 functions = get_rdkit_descriptor_functions()
101 if args.select and args.select != "None":
102 selected = args.select.split(",")
103 functions = [(name, f) for name, f in functions if name in selected]
104
105 if args.header:
106 args.outfile.write(
107 "%s\n" % "\t".join(["MoleculeID"] + [name for name, f in functions])
108 )
109
110 for mol in supplier:
111 if not mol:
112 continue
113 descs = descriptors(mol, functions)
114 try:
115 molecule_id = mol.GetProp("_Name")
116 except KeyError:
117 molecule_id = Chem.MolToSmiles(mol)
118 args.outfile.write(
119 "%s\n"
120 % "\t".join([molecule_id] + [str(round(res, 6)) for name, res in descs])
121 )