diff cheminfolib.py @ 15:4242b4d68e9c draft default tip

planemo upload for repository https://github.com/bgruening/galaxytools/tree/master/chemicaltoolbox/openbabel commit d9c51279c061a1da948a2582d5b502ca7573adbf
author bgruening
date Thu, 15 Aug 2024 11:06:27 +0000
parents 1400d1977e7b
children
line wrap: on
line diff
--- a/cheminfolib.py	Tue Nov 10 20:33:21 2020 +0000
+++ b/cheminfolib.py	Thu Aug 15 11:06:27 2024 +0000
@@ -11,28 +11,32 @@
 import tempfile
 from multiprocessing import Pool
 
-
 try:
     from galaxy import eggs
-    eggs.require('psycopg2')
+
+    eggs.require("psycopg2")
 except ImportError:
     psycopg2 = None
-    print('psycopg2 is not available. It is currently used in the pgchem wrappers, that are not shipped with default CTB')
+    print(
+        "psycopg2 is not available. It is currently used in the pgchem wrappers, that are not shipped with default CTB"
+    )
 
 try:
     from openbabel import openbabel, pybel
+
     openbabel.obErrorLog.StopLogging()
 except ImportError:
     openbabel, pybel = None, None
-    print('OpenBabel could not be found. A few functions are not available without OpenBabel.')
+    print(
+        "OpenBabel could not be found. A few functions are not available without OpenBabel."
+    )
 
 
 def CountLines(path):
-    out = subprocess.Popen(['wc', '-l', path],
-                           stdout=subprocess.PIPE,
-                           stderr=subprocess.STDOUT
-                           ).communicate()[0]
-    return int(out.partition(b' ')[0])
+    out = subprocess.Popen(
+        ["wc", "-l", path], stdout=subprocess.PIPE, stderr=subprocess.STDOUT
+    ).communicate()[0]
+    return int(out.partition(b" ")[0])
 
 
 def grep(pattern, file_obj):
@@ -49,15 +53,15 @@
     for line_counter, line in enumerate(open(filepath)):
         if line_counter > 10000:
             break
-        if line.find('$$$$') != -1:
-            return 'sdf'
-        elif line.find('@<TRIPOS>MOLECULE') != -1:
-            return 'mol2'
-        elif line.find('ligand id') != -1:
-            return 'drf'
-        elif possible_inchi and re.findall('^InChI=', line):
-            return 'inchi'
-        elif re.findall(r'^M\s+END', line):
+        if line.find("$$$$") != -1:
+            return "sdf"
+        elif line.find("@<TRIPOS>MOLECULE") != -1:
+            return "mol2"
+        elif line.find("ligand id") != -1:
+            return "drf"
+        elif possible_inchi and re.findall("^InChI=", line):
+            return "inchi"
+        elif re.findall(r"^M\s+END", line):
             mol = True
         # first line is not an InChI, so it can't be an InChI file
         possible_inchi = False
@@ -65,99 +69,128 @@
     if mol:
         # END can occures before $$$$, so and SDF file will
         # be recognised as mol, if you not using this hack'
-        return 'mol'
-    return 'smi'
+        return "mol"
+    return "smi"
 
 
 def db_connect(args):
     try:
-        db_conn = psycopg2.connect("dbname=%s user=%s host=%s password=%s" % (args.dbname, args.dbuser, args.dbhost, args.dbpasswd))
+        db_conn = psycopg2.connect(
+            "dbname=%s user=%s host=%s password=%s"
+            % (args.dbname, args.dbuser, args.dbhost, args.dbpasswd)
+        )
         return db_conn
     except psycopg2.Error:
-        sys.exit('Unable to connect to the db')
+        sys.exit("Unable to connect to the db")
 
 
 ColumnNames = {
-    'can_smiles': 'Canonical SMILES',
-    'can': 'Canonical SMILES',
-    'inchi': 'InChI',
-    'inchi_key': 'InChI key',
-    'inchi_key_first': 'InChI key first',
-    'inchi_key_last': 'InChI key last',
-    'molwt': 'Molecular weight',
-    'hbd': 'Hydrogen-bond donors',
-    'donors': 'Hydrogen-bond donors',
-    'hba': 'Hydrogen-bond acceptors',
-    'acceptors': 'Hydrogen-bond acceptors',
-    'rotbonds': 'Rotatable bonds',
-    'logp': 'logP',
-    'psa': 'Polar surface area',
-    'mr': 'Molecular refractivity',
-    'atoms': 'Number of heavy atoms',
-    'rings': 'Number of rings',
-    'set_bits': 'FP2 bits',
-    'id': 'Internal identifier',
-    'tani': 'Tanimoto coefficient',
-    'spectrophore': 'Spectrophores(TM)',
-    'dist_spectrophore': 'Spectrophores(TM) distance to target',
-    'synonym': 'Entry id',
+    "can_smiles": "Canonical SMILES",
+    "can": "Canonical SMILES",
+    "inchi": "InChI",
+    "inchi_key": "InChI key",
+    "inchi_key_first": "InChI key first",
+    "inchi_key_last": "InChI key last",
+    "molwt": "Molecular weight",
+    "hbd": "Hydrogen-bond donors",
+    "donors": "Hydrogen-bond donors",
+    "hba": "Hydrogen-bond acceptors",
+    "acceptors": "Hydrogen-bond acceptors",
+    "rotbonds": "Rotatable bonds",
+    "logp": "logP",
+    "psa": "Polar surface area",
+    "mr": "Molecular refractivity",
+    "atoms": "Number of heavy atoms",
+    "rings": "Number of rings",
+    "set_bits": "FP2 bits",
+    "id": "Internal identifier",
+    "tani": "Tanimoto coefficient",
+    "spectrophore": "Spectrophores(TM)",
+    "dist_spectrophore": "Spectrophores(TM) distance to target",
+    "synonym": "Entry id",
 }
 
 OBDescriptor = {
-    'atoms': ["atoms", "Number of atoms"],
-    'hatoms': ["hatoms", "Number of heavy atoms"],  # self defined tag hatoms in plugindefines.txt
-    'can_smiles': ["cansmi", "Canonical SMILES"],
-    'can_smilesNS': ["cansmiNS", "Canonical SMILES without isotopes or stereo"],
+    "atoms": ["atoms", "Number of atoms"],
+    "hatoms": [
+        "hatoms",
+        "Number of heavy atoms",
+    ],  # self defined tag hatoms in plugindefines.txt
+    "can_smiles": ["cansmi", "Canonical SMILES"],
+    "can_smilesNS": ["cansmiNS", "Canonical SMILES without isotopes or stereo"],
     # ["abonds", "Number of aromatic bonds"],
     # ["bonds", "Number of bonds"],
     # ["dbonds", "Number of double bonds"],
     # ["formula", "Chemical formula"],
-    'hba': ["HBA1", "Number of Hydrogen Bond Acceptors 1 (JoelLib)"],
-    'hba2': ["HBA2", "Number of Hydrogen Bond Acceptors 2 (JoelLib)"],
-    'hbd': ["HBD", "Number of Hydrogen Bond Donors (JoelLib)"],
-    'inchi': ["InChI", "IUPAC InChI identifier"],
-    'inchi_key': ["InChIKey", "InChIKey"],
+    "hba": ["HBA1", "Number of Hydrogen Bond Acceptors 1 (JoelLib)"],
+    "hba2": ["HBA2", "Number of Hydrogen Bond Acceptors 2 (JoelLib)"],
+    "hbd": ["HBD", "Number of Hydrogen Bond Donors (JoelLib)"],
+    "inchi": ["InChI", "IUPAC InChI identifier"],
+    "inchi_key": ["InChIKey", "InChIKey"],
     # ["L5", "Lipinski Rule of Five"],
-    'logp': ["logP", "octanol/water partition coefficient"],
-    'mr': ["MR", "molar refractivity"],
-    'molwt': ["MW", "Molecular Weight filter"],
+    "logp": ["logP", "octanol/water partition coefficient"],
+    "mr": ["MR", "molar refractivity"],
+    "molwt": ["MW", "Molecular Weight filter"],
     # ["nF", "Number of Fluorine Atoms"],
     # ["s", "SMARTS filter"],
     # ["sbonds", "Number of single bonds"],
     # ["smarts", "SMARTS filter"],
     # ["tbonds", "Number of triple bonds"],
     # ["title", "For comparing a molecule's title"],
-    'psa': ["TPSA", "topological polar surface area"],
-    'rotbonds': ['ROTATABLE_BOND', 'rotatable bonds'],
+    "psa": ["TPSA", "topological polar surface area"],
+    "rotbonds": ["ROTATABLE_BOND", "rotatable bonds"],
 }
 
 
 def print_output(args, rows):
-    if args.oformat == 'table':
-        outfile = open(args.output, 'w')
-        requested_fields = (filter(lambda x: x not in ["[", "]", "'"], args.fetch)).split(', ')
+    if args.oformat == "table":
+        outfile = open(args.output, "w")
+        requested_fields = (
+            filter(lambda x: x not in ["[", "]", "'"], args.fetch)
+        ).split(", ")
         if args.header:
-            outfile.write('Identifier\t' + '\t'.join([ColumnNames[key] for key in requested_fields]) + '\n')
+            outfile.write(
+                "Identifier\t"
+                + "\t".join([ColumnNames[key] for key in requested_fields])
+                + "\n"
+            )
         for row in rows:
-            outfile.write(row['synonym'] + '\t' + '\t'.join([str(row[key]) for key in requested_fields]) + '\n')
+            outfile.write(
+                row["synonym"]
+                + "\t"
+                + "\t".join([str(row[key]) for key in requested_fields])
+                + "\n"
+            )
 
-    elif args.oformat in ['sdf', 'mol2']:
+    elif args.oformat in ["sdf", "mol2"]:
         outfile = pybel.Outputfile(args.oformat, args.output, overwrite=True)
         for row in rows:
             try:
-                mol = pybel.readstring('sdf', row['mol'])
-                if args.oformat == 'sdf':
-                    keys = filter(lambda x: x not in ["[", "]", "'"], args.fetch).split(', ')
-                    mol.data.update({ColumnNames['synonym']: row['synonym']})
-                    if 'inchi_key' in keys:
-                        keys = (', '.join(keys).replace("inchi_key", "inchi_key_first, inchi_key_last")).split(', ')
-                    [mol.data.update({ColumnNames[key]: row[key]}) for key in keys if key]
+                mol = pybel.readstring("sdf", row["mol"])
+                if args.oformat == "sdf":
+                    keys = filter(lambda x: x not in ["[", "]", "'"], args.fetch).split(
+                        ", "
+                    )
+                    mol.data.update({ColumnNames["synonym"]: row["synonym"]})
+                    if "inchi_key" in keys:
+                        keys = (
+                            ", ".join(keys).replace(
+                                "inchi_key", "inchi_key_first, inchi_key_last"
+                            )
+                        ).split(", ")
+                    [
+                        mol.data.update({ColumnNames[key]: row[key]})
+                        for key in keys
+                        if key
+                    ]
                 outfile.write(mol)
             except OSError:
                 pass
     else:
-        outfile = open(args.output, 'w')
-        outfile.write('\n'.join(['%s\t%s' % (row[args.oformat], row['synonym']) for row in rows]))
+        outfile = open(args.output, "w")
+        outfile.write(
+            "\n".join(["%s\t%s" % (row[args.oformat], row["synonym"]) for row in rows])
+        )
     outfile.close()
 
 
@@ -167,31 +200,37 @@
 
 def get_properties_ext(mol):
     HBD = pybel.Smarts("[!#6;!H0]")
-    HBA = pybel.Smarts(("[$([$([#8,#16]);!$(*=N~O);"
-                        "!$(*~N=O);X1,X2]),$([#7;v3;"
-                        "!$([nH]);!$(*(-a)-a)])]"
-                        ))
+    HBA = pybel.Smarts(
+        (
+            "[$([$([#8,#16]);!$(*=N~O);"
+            "!$(*~N=O);X1,X2]),$([#7;v3;"
+            "!$([nH]);!$(*(-a)-a)])]"
+        )
+    )
     calc_desc_dict = mol.calcdesc()
 
     try:
-        logp = calc_desc_dict['logP']
+        logp = calc_desc_dict["logP"]
     except KeyError:
-        logp = calc_desc_dict['LogP']
+        logp = calc_desc_dict["LogP"]
 
-    return {"molwt": mol.molwt,
-            "logp": logp,
-            "donors": len(HBD.findall(mol)),
-            "acceptors": len(HBA.findall(mol)),
-            "psa": calc_desc_dict['TPSA'],
-            "mr": calc_desc_dict['MR'],
-            "rotbonds": mol.OBMol.NumRotors(),
-            "can": mol.write("can").split()[0].strip(),  # tthis one works fine for both zinc and chembl (no ZINC code added after can descriptor string)
-            "inchi": mol.write("inchi").strip(),
-            "inchi_key": get_inchikey(mol).strip(),
-            "rings": len(mol.sssr),
-            "atoms": mol.OBMol.NumHvyAtoms(),
-            "spectrophore": OBspectrophore(mol),
-            }
+    return {
+        "molwt": mol.molwt,
+        "logp": logp,
+        "donors": len(HBD.findall(mol)),
+        "acceptors": len(HBA.findall(mol)),
+        "psa": calc_desc_dict["TPSA"],
+        "mr": calc_desc_dict["MR"],
+        "rotbonds": mol.OBMol.NumRotors(),
+        "can": mol.write("can")
+        .split()[0]
+        .strip(),  # tthis one works fine for both zinc and chembl (no ZINC code added after can descriptor string)
+        "inchi": mol.write("inchi").strip(),
+        "inchi_key": get_inchikey(mol).strip(),
+        "rings": len(mol.sssr),
+        "atoms": mol.OBMol.NumHvyAtoms(),
+        "spectrophore": OBspectrophore(mol),
+    }
 
 
 def get_inchikey(mol):
@@ -206,10 +245,12 @@
     spectrophore = pybel.ob.OBSpectrophore()
     # Parameters: rotation angle = 20, normalization for mean and sd, accuracy = 3.0 A and non-stereospecific cages.
     spectrophore.SetNormalization(spectrophore.NormalizationTowardsZeroMeanAndUnitStd)
-    return ', '.join(["%.3f" % value for value in spectrophore.GetSpectrophore(mol.OBMol)])
+    return ", ".join(
+        ["%.3f" % value for value in spectrophore.GetSpectrophore(mol.OBMol)]
+    )
 
 
-def split_library(lib_path, lib_format='sdf', package_size=None):
+def split_library(lib_path, lib_format="sdf", package_size=None):
     """
     Split a library of compounds. Usage: split_library(lib_path, lib_format, package_size)
     IT currently ONLY WORKS FOR SD-Files
@@ -217,18 +258,39 @@
     pack = 1
     mol_counter = 0
 
-    outfile = open('/%s/%s_pack_%i.%s' % ('/'.join(lib_path.split('/')[:-1]), lib_path.split('/')[-1].split('.')[0], pack, 'sdf'), 'w')
+    outfile = open(
+        "/%s/%s_pack_%i.%s"
+        % (
+            "/".join(lib_path.split("/")[:-1]),
+            lib_path.split("/")[-1].split(".")[0],
+            pack,
+            "sdf",
+        ),
+        "w",
+    )
 
-    for line in open(lib_path, 'r'):
+    for line in open(lib_path, "r"):
         outfile.write(line)
-        if line.strip() == '$$$$':
+        if line.strip() == "$$$$":
             mol_counter += 1
             if mol_counter % package_size == 0:
                 outfile.close()
                 pack += 1
-                outfile = open('/%s/%s_pack_%i.%s' % ('/'.join(lib_path.split('/')[:-1]), lib_path.split('/')[-1].split('.')[0], pack, 'sdf'), 'w')
+                outfile = open(
+                    "/%s/%s_pack_%i.%s"
+                    % (
+                        "/".join(lib_path.split("/")[:-1]),
+                        lib_path.split("/")[-1].split(".")[0],
+                        pack,
+                        "sdf",
+                    ),
+                    "w",
+                )
                 if mol_counter * 10 % package_size == 0:
-                    print('%i molecules parsed, starting pack nr. %i' % (mol_counter, pack - 1))
+                    print(
+                        "%i molecules parsed, starting pack nr. %i"
+                        % (mol_counter, pack - 1)
+                    )
     outfile.close()
 
     return True
@@ -242,7 +304,7 @@
     output_files = []
     tfile = tempfile.NamedTemporaryFile(delete=False)
 
-    smiles_handle = open(smiles_file, 'r')
+    smiles_handle = open(smiles_file, "r")
     for count, line in enumerate(smiles_handle):
         if count % structures_in_one_file == 0 and count != 0:
             tfile.close()
@@ -257,16 +319,19 @@
 
 def mp_run(input_path, regex, PROCESSES, function_to_call):
     paths = []
-    [paths.append(compound_file) for compound_file in glob.glob(str(input_path) + str(regex))]
+    [
+        paths.append(compound_file)
+        for compound_file in glob.glob(str(input_path) + str(regex))
+    ]
     paths.sort()
 
     pool = Pool(processes=PROCESSES)
-    print('Process initialized with', PROCESSES, 'processors')
+    print("Process initialized with", PROCESSES, "processors")
     result = pool.map_async(function_to_call, paths)
     result.get()
 
     return paths
 
 
-if __name__ == '__main__':
+if __name__ == "__main__":
     print(check_filetype(sys.argv[1]))