add binding for DockQ (https://github.com/bjornwallner/DockQ)

modules/bindings/doc/bindings.rst

@@ -24,3 +24,4 @@ So far, the binding module includes:
   naccess
   lga
   cadscore
+  dockq

modules/bindings/doc/dockq.rst

+:mod:`~ost.bindings.dockq` - Evaluate protein-protein interfaces
+================================================================
+
+.. module:: ost.bindings.dockq
+  :synopsis: Evaluate protein-protein interfaces
+
+
+.. autofunction:: ost.bindings.dockq.DockQ
+
+.. autoclass:: ost.bindings.dockq.DockQResult
+   :members:
+   :member-order: bysource

modules/bindings/pymod/CMakeLists.txt

@@ -16,6 +16,7 @@ cadscore.py
 kclust.py
 ialign.py
 mmseqs2.py
+dockq.py
 )
 
 set(OST_BINDINGS_PYMOD_SOURCES

modules/bindings/pymod/dockq.py

+import sys
+import os
+import subprocess
+import tempfile
+import shutil
+
+from ost import io
+from ost import mol
+
+def _Setup(mdl, ref, mdl_ch1, mdl_ch2, ref_ch1, ref_ch2):
+    """ Performs parameter checks and dumps files for DockQ
+
+    In case of dimeric interfaces the respective chains are selected from
+    mdl/trg, renamed to A/B and dumped to disk.
+
+    In case of interfaces with more chains involved, we simply select the
+    specified chains and do no renaming before dumping to disk.
+    """
+    if isinstance(mdl_ch1, str):
+        mdl_ch1 = [mdl_ch1]
+    if isinstance(mdl_ch2, str):
+        mdl_ch2 = [mdl_ch2]
+    if isinstance(ref_ch1, str):
+        ref_ch1 = [ref_ch1]
+    if isinstance(ref_ch2, str):
+        ref_ch2 = [ref_ch2]
+
+    if len(mdl_ch1) == 0:
+        raise RuntimeError("mdl_ch1 is empty")
+    if len(mdl_ch2) == 0:
+        raise RuntimeError("mdl_ch2 is empty")
+
+    if len(mdl_ch1) != len(ref_ch1):
+        raise RuntimeError("mdl_ch1/ref_ch1 inconsistent in size")
+    if len(mdl_ch2) != len(ref_ch2):
+        raise RuntimeError("mdl_ch2/ref_ch2 inconsistent in size")
+
+    for cname in mdl_ch1:
+        ch = mdl.FindChain(cname)
+        if not ch.IsValid():
+            raise RuntimeError(f"Chain {cname} specified in mdl_ch1 not "
+                               f"present in mdl")
+
+    for cname in mdl_ch2:
+        ch = mdl.FindChain(cname)
+        if not ch.IsValid():
+            raise RuntimeError(f"Chain {cname} specified in mdl_ch2 not "
+                               f"present in mdl")
+
+    for cname in ref_ch1:
+        ch = ref.FindChain(cname)
+        if not ch.IsValid():
+            raise RuntimeError(f"Chain {cname} specified in ref_ch1 not "
+                               f"present in ref")
+
+    for cname in ref_ch2:
+        ch = ref.FindChain(cname)
+        if not ch.IsValid():
+            raise RuntimeError(f"Chain {cname} specified in ref_ch2 not "
+                               f"present in ref")
+
+    mdl_to_dump = mdl.CreateFullView()
+    ref_to_dump = ref.CreateFullView()
+
+    if len(mdl_ch1) == 1 and len(mdl_ch2) == 1:
+        # Dimer processing of mdl => Create new entity only containing 
+        # the two specified chains and rename them to A, B
+        mdl_to_dump = mol.CreateEntityFromView(mdl_to_dump, True)
+        tmp = mol.CreateEntity()
+        ed = tmp.EditXCS()
+        ch1 = mdl_to_dump.FindChain(mdl_ch1[0])
+        ed.InsertChain("A", ch1, deep=True)
+        ch2 = mdl_to_dump.FindChain(mdl_ch2[0])
+        ed.InsertChain("B", ch2, deep=True)
+        mdl_ch1 = ["A"]
+        mdl_ch2 = ["B"]
+        mdl_to_dump = tmp
+
+        # Same for ref
+        ref_to_dump = mol.CreateEntityFromView(ref_to_dump, True)
+        tmp = mol.CreateEntity()
+        ed = tmp.EditXCS()
+        ch1 = ref_to_dump.FindChain(ref_ch1[0])
+        ed.InsertChain("A", ch1, deep=True)
+        ch2 = ref_to_dump.FindChain(ref_ch2[0])
+        ed.InsertChain("B", ch2, deep=True)
+        ref_ch1 = ["A"]
+        ref_ch2 = ["B"]
+        ref_to_dump = tmp
+    else:
+        # Interface with more chains...
+        raise NotImplementedError("DockQ computations beyond two interacting "
+                                  "chains has not been properly tested...")
+        #mdl_chain_names = mdl_ch1 + mdl_ch2
+        #ref_chain_names = ref_ch1 + ref_ch2
+        #mdl_to_dump = mdl_to_dump.Select(f"cname={','.join(mdl_chain_names)}")
+        #ref_to_dump = ref_to_dump.Select(f"cname={','.join(ref_chain_names)}")
+
+    # first write structures to string, only create a tmpdir and the actual
+    # files if this succeeds
+    mdl_str = io.EntityToPDBStr(mdl_to_dump)
+    ref_str = io.EntityToPDBStr(ref_to_dump)
+
+    tmp_dir = tempfile.mkdtemp()
+    with open(os.path.join(tmp_dir, "mdl.pdb"), 'w') as fh:
+        fh.write(mdl_str)
+    with open(os.path.join(tmp_dir, "ref.pdb"), 'w') as fh:
+        fh.write(ref_str)
+
+    return (tmp_dir, mdl_ch1, mdl_ch2, ref_ch1, ref_ch2)
+
+class DockQResult:
+    """ DockQ result object
+    """
+    def __init__(self, Fnat, Fnonnat, native_contacts, model_contacts, iRMS,
+                 LRMS, DockQ):
+        self._Fnat = Fnat
+        self._Fnonnat = Fnonnat
+        self._native_contacts = native_contacts
+        self._model_contacts = model_contacts
+        self._iRMS = iRMS
+        self._LRMS = LRMS
+        self._DockQ = DockQ
+
+    @property
+    def Fnat(self):
+        """ DockQ - Fnat output
+
+        :type: :class:`float`
+        """
+        return self._Fnat
+    
+    @property
+    def Fnonnat(self):
+        """ DockQ - Fnonnat output
+
+        :type: :class:`float`
+        """
+        return self._Fnonnat
+
+    @property
+    def native_contacts(self):
+        """ DockQ - number native contacts
+
+        :type: :class:`int`
+        """
+        return self._native_contacts
+
+    @property
+    def model_contacts(self):
+        """ DockQ - number model contacts
+
+        :type: :class:`int`
+        """
+        return self._model_contacts
+
+    @property
+    def iRMS(self):
+        """ DockQ - iRMS output
+
+        :type: :class:`float`
+        """
+        return self._iRMS
+
+    @property
+    def LRMS(self):
+        """ DockQ - LMRS output 
+
+        :type: :class:`float`
+        """
+        return self._LRMS
+
+    @property
+    def DockQ(self):
+        """ DockQ - DockQ output
+
+        :type: :class:`float`
+        """
+        return self._DockQ
+
+    def JSONSummary(self):
+        """ Returns JSON serializable summary
+        """
+        return {"Fnat": self.Fnat,
+                "Fnonnat": self.Fnonnat,
+                "native_contacts": self.native_contacts,
+                "model_contacts": self.model_contacts,
+                "iRMS": self.iRMS,
+                "LRMS": self.LRMS,
+                "DockQ": self.DockQ}
+
+    @staticmethod
+    def FromDockQOutput(output):
+        """ Static constructor from raw DockQ output
+
+        :param output: Raw output from DockQ executable
+        :type output: :class:`str`
+        :returns: Object of type :class:`DockQResult`
+        """
+        Fnat = None
+        Fnonnat = None
+        native_contacts = None
+        model_contacts = None
+        iRMS = None
+        LRMS = None
+        DockQ = None
+
+        for line in output.splitlines():
+            if line.startswith('*'):
+                continue
+            if line.startswith("Fnat"):
+                Fnat = float(line.split()[1])
+                native_contacts = int(line.split()[5])
+            elif line.startswith("Fnonnat"):
+                Fnonnat = float(line.split()[1])
+                model_contacts = int(line.split()[5])
+            elif line.startswith("iRMS"):
+                iRMS = float(line.split()[1])
+            elif line.startswith("LRMS"):
+                LRMS = float(line.split()[1])
+            elif line.startswith("DockQ"):
+                DockQ = float(line.split()[1])
+
+        return DockQResult(Fnat, Fnonnat, native_contacts, model_contacts,
+                           iRMS, LRMS, DockQ)
+
+
+def DockQ(dockq_exec, mdl, ref, mdl_ch1, mdl_ch2, ref_ch1,
+          ref_ch2):
+    """ Computes DockQ for specified interface 
+
+    DockQ is available from https://github.com/bjornwallner/DockQ - 
+    For this binding to work, DockQ must be properly installed and its
+    dependencies must be available (numpy, Biopython).
+
+    :param dockq_exec: Path to DockQ.py script from DockQ repository
+    :type dockq_exec: :class:`str`
+    :param mdl: Model structure
+    :type mdl: :class:`ost.mol.EntityView`/:class:`ost.mol.EntityHandle`
+    :param ref: Reference structure, i.e. native structure
+    :type ref: :class:`ost.mol.EntityView`/:class:`ost.mol.EntityHandle`
+    :param mdl_ch1: Specifies chain(s) in model constituting first part of
+                    interface
+    :type mdl_ch1: :class:`str`/:class:`list` of :class:`str`
+    :param mdl_ch2: Specifies chain(s) in model constituting second part of
+                    interface
+    :type mdl_ch2: :class:`str`/:class:`list` of :class:`str`
+    :param ref_ch1: ref equivalent of mdl_ch1
+    :type ref_ch1: :class:`str`/:class:`list` of :class:`str`
+    :param ref_ch2: ref equivalent of mdl_ch2
+    :type ref_ch2: :class:`str`/:class:`list` of :class:`str`
+    :returns: Result object of type :class:`DockQResult`
+    """
+    if not os.path.exists(dockq_exec):
+        raise RuntimeError(f"DockQ executable ({dockq_exec}) does not exist")
+
+    tmp_dir, mdl_ch1, mdl_ch2, ref_ch1, ref_ch2 = \
+    _Setup(mdl, ref, mdl_ch1, mdl_ch2, ref_ch1, ref_ch2)
+
+    cmd = [sys.executable, dockq_exec, os.path.join(tmp_dir, "mdl.pdb"),
+           os.path.join(tmp_dir, "ref.pdb")]
+
+    # add mdl/ref chains
+    cmd.append("-model_chain1")
+    cmd += mdl_ch1
+    cmd.append("-model_chain2")
+    cmd += mdl_ch2
+    cmd.append("-native_chain1")
+    cmd += ref_ch1
+    cmd.append("-native_chain2")
+    cmd += ref_ch2
+
+    proc = subprocess.run(cmd, stdout=subprocess.PIPE, stderr=subprocess.PIPE)
+
+    shutil.rmtree(tmp_dir) # cleanup, no matter if DockQ executed successfully 
+
+    if proc.returncode != 0:
+        raise RuntimeError("DockQ run failed - returncode: " + \
+                           str(proc.return_code))
+
+    if proc.stderr.decode() != "":
+        raise RuntimeError("DockQ run failed - stderr: " + proc.stderr.decode())
+
+    return DockQResult.FromDockQOutput(proc.stdout.decode())

modules/bindings/tests/CMakeLists.txt

@@ -10,6 +10,7 @@ set(OST_BINDINGS_UNIT_TESTS
   test_ialign.py
   test_lga.py
   test_hbplus.py
+  test_dockq.py
 )
 
 ost_unittest(MODULE bindings 

modules/bindings/tests/test_dockq.py

+''' Unit tests for the DockQ wrapper
+'''
+import sys
+import unittest
+
+import ost
+from ost import settings
+from ost import io
+from ost.bindings import dockq
+
+
+class TestDockQBinding(unittest.TestCase):
+
+    def testDimerExample(self):
+        # runs Dimer example from Bjoerns github repo and checks result
+        mdl = io.LoadPDB("testfiles/dockq_model.pdb")
+        ref = io.LoadPDB("testfiles/dockq_native.pdb")
+        result = dockq.DockQ(settings.Locate("DockQ.py"), mdl, ref,
+                             "A", "B", "A", "B")
+        self.assertEqual(result.Fnat, 0.533)
+        self.assertEqual(result.native_contacts, 60)
+        self.assertEqual(result.Fnonnat, 0.238)
+        self.assertEqual(result.iRMS, 1.232)
+        self.assertEqual(result.LRMS, 1.516)
+        self.assertEqual(result.DockQ, 0.700)
+
+    def testMultichainExample(self):
+        # multichain means one or both interface partner (ligand/receptor)
+        # consist of multiple chains. DockQ provides such functionality and the
+        # interface of the binding can deal with it. However, I just had no time
+        # for proper testing. The binding therefore raises a NotImplementedError
+        # for such cases.
+        mdl = io.LoadPDB("testfiles/dockq_model.pdb")
+        ref = io.LoadPDB("testfiles/dockq_native.pdb")
+        with self.assertRaises(NotImplementedError):
+            result = dockq.DockQ(settings.Locate("DockQ.py"), mdl, ref,
+                                 ["A", "B"], "B", ["A", "B"], "B")
+
+
+if __name__ == "__main__":
+
+    try:
+        settings.Locate("DockQ.py")
+    except:
+        print("Could not find DockQ.py, could not test binding...")
+        sys.exit(0)
+
+    from ost import testutils
+    testutils.RunTests()