Merge 4b65ba0 into bd7616b

CHANGELOG.rst

@@ -17,6 +17,8 @@ Added
 - ``Sequence.seq_encoded`` to allow turning a sequence into an encoded list
 - ``Sequence.encoded_matrix`` to give the entire alignment as encoded matrix
 - ``SequenceFile.filter_gapped`` to filter sequences with a certain threshold of gaps
+- ``SequenceFile.to_string`` and ``ContactMap.to_string`` methods
+- ``ContactMapMatrixFigure`` added to illustrate prediction signal of entire ``ContactMap``
 Changed
 ~~~~~~~
 - Changed API interface for ``conkit.plot`` in accordance to necessary changes for above

conkit/command_line/conkit_plot.py

@@ -144,6 +144,25 @@ def add_contact_map_chord_args(subparsers):
     subparser.set_defaults(which='contact_map_chord')
 
 
+def add_contact_map_matrix_args(subparsers):
+    description = u"""
+This command will plot a contact map matrix using the provided contacts
+alongside any additional information.
+
+"""
+    subparser = subparsers.add_parser(
+        'cmat',
+        help="Plot a contact map matrix",
+        description=description,
+        formatter_class=argparse.RawDescriptionHelpFormatter)
+    subparser.add_argument('-e', dest='otherfile', default=None, help='a second contact map to plot for comparison')
+    subparser.add_argument('-ef', dest='otherformat', default=None, help='the format of the second contact map')
+    _add_default_args(subparser)
+    _add_sequence_default_args(subparser)
+    _add_contact_default_args(subparser)
+    subparser.set_defaults(which='contact_map_matrix')
+
+
 def add_contact_density_args(subparsers):
     description = u"""
 This command will plot a contact density plot using the provided
@@ -331,6 +350,24 @@ def altloc_remove(cmap):
         figure = conkit.plot.ContactMapChordFigure(con_sliced, use_conf=args.confidence, legend=True)
         figure_aspect_ratio = 1.0
 
+    elif args.which == 'contact_map_matrix':
+
+        seq = conkit.io.read(args.seqfile, args.seqformat)[0]
+        con = conkit.io.read(args.confile, args.conformat)[0]
+
+        con.sequence = seq
+        con.assign_sequence_register()
+
+        if args.otherfile:
+            other = conkit.io.read(args.otherfile, args.otherformat)[0]
+            other.sequence = seq
+            other.assign_sequence_register()
+        else:
+            other = con
+
+        figure = conkit.plot.ContactMapMatrixFigure(con, other=other, legend=True)
+        figure_aspect_ratio = 1.0
+
     elif args.which == 'contact_density':
         logger.info('Min sequence separation for contacting residues: %d', args.dtn)
         logger.info('Contact list cutoff factor: %f * L', args.dfactor)

conkit/core/contactmap.py

@@ -40,6 +40,7 @@
 
 import collections
 import numpy as np
+import os
 import sys
 
 if sys.version_info.major < 3:
@@ -467,12 +468,10 @@ def calculate_jaccard_index(self, other):
            [doi: 10.1093/bib/bbw106].
 
         """
-        intersection = np.sum([1 for contact in self if contact.id in other])
-        union = len(self) \
-            + np.sum([1 for contact in other if contact.id not in self])
-        # If self and other are both empty, we define J(x,y) = 1
+        union = len(self) + np.sum([1 for contact in other if contact.id not in self])
         if union == 0:
             return 1.0
+        intersection = np.sum([1 for contact in self if contact.id in other])
         return float(intersection) / union
 
     def calculate_kernel_density(self, *args, **kwargs):
@@ -845,6 +844,11 @@ def sort(self, kword, reverse=False, inplace=False):
         contact_map._sort(kword, reverse)
         return contact_map
 
+    def to_string(self):
+        """Return the :obj:`ContactMap <conkit.core.contactmap.ContactMap>` as :obj:`str`"""
+        content = ["%d\t%d\t%.5f" % (c.res1_seq, c.res2_seq, c.raw_score) for c in self]
+        return os.linesep.join(content)
+
     @staticmethod
     def _adjust(contact_map, keymap):
         """Adjust res_altseq entries to insertions and deletions"""

conkit/core/sequencefile.py

@@ -39,6 +39,7 @@
 __version__ = "1.0"
 
 import numpy as np
+import os
 import sys
 
 if sys.version_info.major < 3:
@@ -448,6 +449,11 @@ def sort(self, kword, reverse=False, inplace=False):
         sequence_file._sort(kword, reverse)
         return sequence_file
 
+    def to_string(self):
+        """Return the :obj:`SequenceFile <conkit.core.sequencefile.SequenceFile>` as :obj:`str`"""
+        content = [s.seq for s in self]
+        return os.linesep.join(content)
+
     def trim(self, start, end, inplace=False):
         """Trim the :obj:`SequenceFile <conkit.core.sequencefile.SequenceFile>`
 

conkit/plot/__init__.py

@@ -50,6 +50,12 @@ def ContactMapChordFigure(*args, **kwargs):
     return ContactMapChordFigure(*args, **kwargs)
 
 
+def ContactMapMatrixFigure(*args, **kwargs):
+    """:obj:`ContactMapMatrixFigure <conkit.plot.contactmatrix.ContactMapMatrixFigure>` instance"""
+    from conkit.plot.contactmapmatrix import ContactMapMatrixFigure
+    return ContactMapMatrixFigure(*args, **kwargs)
+
+
 def ContactDensityFigure(*args, **kwargs):
     """:obj:`ContactDensityFigure <conkit.plot.contactdensity.ContactDensityFigure>` instance"""
     from conkit.plot.contactdensity import ContactDensityFigure

conkit/plot/contactmap.py

@@ -204,13 +204,10 @@ def draw(self):
         else:
             self_radius = other_radius = 0.48
 
-        self._patch_scatter(
-            reference_data[:, 0], reference_data[:, 1], facecolor=reference_colors, radius=0.5, linewidth=0)
-        self._patch_scatter(
-            reference_data[:, 1], reference_data[:, 0], facecolor=reference_colors, radius=0.5, linewidth=0)
-        self._patch_scatter(self_data[:, 1], self_data[:, 0], facecolor=self_colors, radius=self_radius, linewidth=0)
-        self._patch_scatter(
-            other_data[:, 0], other_data[:, 1], facecolor=other_colors, radius=other_radius, linewidth=0)
+        self._patch_scatter(reference_data[:, 0], reference_data[:, 1], symbol="o", facecolor=reference_colors, radius=0.5, linewidth=0)
+        self._patch_scatter(reference_data[:, 1], reference_data[:, 0], symbol="o", facecolor=reference_colors, radius=0.5, linewidth=0)
+        self._patch_scatter(self_data[:, 1], self_data[:, 0], symbol="o", facecolor=self_colors, radius=self_radius, linewidth=0)
+        self._patch_scatter(other_data[:, 0], other_data[:, 1], symbol="o", facecolor=other_colors, radius=other_radius, linewidth=0)
 
         if self.lim:
             min_max_data = np.arange(self.lim[0], self.lim[1] + 1)

conkit/plot/contactmapchord.py

@@ -173,8 +173,8 @@ def draw(self):
             xy_highlight.append(xy.tolist())
 
         radius = ContactMapChordFigure.get_radius_around_circle(coords)
-        self._patch_scatter(coords[:, 0], coords[:, 1], facecolor=colors, linewidth=0.0, radius=radius)
-        self._patch_scatter(*zip(*xy_highlight), facecolor="none", edgecolor="#000000", radius=radius)
+        self._patch_scatter(coords[:, 0], coords[:, 1], symbol="o", facecolor=colors, linewidth=0.0, radius=radius)
+        self._patch_scatter(*zip(*xy_highlight), symbol="o", facecolor="none", edgecolor="#000000", radius=radius)
 
         arrow_x, arrow_y = (npoints + npoints / 5, 0)
         self.ax.arrow(arrow_x, arrow_y, 0, npoints / 10, head_width=1.5, color="#000000")

conkit/plot/contactmapmatrix.py

@@ -0,0 +1,192 @@
+# BSD 3-Clause License
+#
+# Copyright (c) 2016-18, University of Liverpool
+# All rights reserved.
+#
+# Redistribution and use in source and binary forms, with or without
+# modification, are permitted provided that the following conditions are met:
+#
+# * Redistributions of source code must retain the above copyright notice, this
+#   list of conditions and the following disclaimer.
+#
+# * Redistributions in binary form must reproduce the above copyright notice,
+#   this list of conditions and the following disclaimer in the documentation
+#   and/or other materials provided with the distribution.
+#
+# * Neither the name of the copyright holder nor the names of its
+#   contributors may be used to endorse or promote products derived from
+#   this software without specific prior written permission.
+#
+# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
+# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
+# DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
+# FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+# DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+# SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
+# CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
+# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
+# OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+"""A module to produce a contact map plot"""
+
+from __future__ import division
+from __future__ import print_function
+
+__author__ = "Felix Simkovic"
+__date__ = "10 Jan 2018"
+__version__ = "1.0"
+
+import matplotlib.collections as mcoll
+import matplotlib.pyplot as plt
+import numpy as np
+
+from conkit.core._struct import _Gap
+from conkit.misc import normalize
+from conkit.plot.figure import Figure
+from conkit.plot.tools import ColorDefinitions, _isinstance
+
+
+class ContactMapMatrixFigure(Figure):
+    """A Figure object specifically for a Contact Map Matrix
+
+    This figure will illustrate the contacts in a contact
+    map matrix. This plot is a very common representation of contacts.
+    With this figure, you can illustrate either your contact
+    map by itself, compared against a second contact map, and/or
+    matched against contacts extracted from a contact map.
+
+    Attributes
+    ----------
+    hierarchy : :obj:`ContactMap <conkit.core.contactmap.ContactMap>`
+       The default contact map hierarchy
+    other : :obj:`ContactMap <conkit.core.contactmap.ContactMap>`
+       The second contact map hierarchy
+    altloc : bool
+       Use the res_altloc positions [default: False]
+
+    Examples
+    --------
+    >>> import conkit
+    >>> cmap = conkit.io.read('toxd/toxd.mat', 'ccmpred').top_map
+    >>> conkit.plot.ContactMapMatrixFigure(cmap)
+
+    """
+
+    def __init__(self, hierarchy, other=None, altloc=False, lim=None, **kwargs):
+        """A new contact map plot
+
+        Parameters
+        ----------
+        hierarchy : :obj:`ContactMap <conkit.core.contactmap.ContactMap>`
+           The default contact map hierarchy
+        other : :obj:`ContactMap <conkit.core.contactmap.ContactMap>`, optional
+           The second contact map hierarchy
+        altloc : bool, optional
+           Use the res_altloc positions [default: False]
+        lim : tuple, list, optional
+           The [min, max] residue numbers to show
+        **kwargs
+           General :obj:`Figure <conkit.plot.figure.Figure>` keyword arguments
+
+        """
+        super(ContactMapMatrixFigure, self).__init__(**kwargs)
+
+        self._hierarchy = None
+        self._other = None
+        self._lim = None
+
+        self.altloc = altloc
+
+        self.hierarchy = hierarchy
+        if other:
+            self.other = other
+        if lim:
+            self.lim = lim
+
+        self.draw()
+
+    def __repr__(self):
+        return "{0}(file_name=\"{1}\")".format(self.__class__.__name__, self.file_name)
+
+    @property
+    def hierarchy(self):
+        return self._hierarchy
+
+    @hierarchy.setter
+    def hierarchy(self, hierarchy):
+        if hierarchy and _isinstance(hierarchy, "ContactMap"):
+            self._hierarchy = hierarchy
+        else:
+            raise TypeError("Invalid hierarchy type: %s" % hierarchy.__class__.__name__)
+
+    @property
+    def other(self):
+        return self._other
+
+    @other.setter
+    def other(self, hierarchy):
+        if hierarchy and _isinstance(hierarchy, "ContactMap"):
+            self._other = hierarchy
+        else:
+            raise TypeError("Invalid hierarchy type: %s" % hierarchy.__class__.__name__)
+
+    @property
+    def lim(self):
+        return self._lim
+
+    @lim.setter
+    def lim(self, lim):
+        if (isinstance(lim, list) or isinstance(lim, tuple)) and len(lim) == 2:
+            self._lim = lim
+        elif (isinstance(lim, list) or isinstance(lim, tuple)):
+            raise ValueError("A list with 2 entries is required!")
+        else:
+            raise TypeError("A list with [min, max] limits is required!")
+
+    def draw(self):
+        _hierarchy = self._hierarchy.rescale()
+
+        self_data = np.array([c for c in _hierarchy.as_list() if all(ci != _Gap.IDENTIFIER for ci in c)])
+        self_colors = ContactMapMatrixFigure._determine_color(_hierarchy)
+        self_rawsc = np.array([c.raw_score for c in _hierarchy if all(ci != _Gap.IDENTIFIER for ci in [c.res1_seq, c.res2_seq])])
+
+        if self._other:
+            _other = self._other.rescale()
+            other_data = np.array([c for c in _other.as_list() if any(ci != _Gap.IDENTIFIER for ci in c)])
+            other_colors = ContactMapMatrixFigure._determine_color(_other)
+            other_rawsc = np.array([c.raw_score for c in _other if all(ci != _Gap.IDENTIFIER for ci in [c.res1_seq, c.res2_seq])])
+        else:
+            other_data = self_data
+            other_colors = self_colors
+            other_rawsc = self_rawsc
+
+        self._patch_scatter(self_data[:, 1], self_data[:, 0], symbol="s", facecolor=self_colors, radius=1.0, linewidth=0)
+        self._patch_scatter(other_data[:, 0], other_data[:, 1], symbol="s", facecolor=other_colors, radius=1.0, linewidth=0)
+
+        if self.lim:
+            min_max_data = np.arange(self.lim[0], self.lim[1] + 1)
+            self.ax.set_xlim(self.lim[0] - 0.5, self.lim[1] + 0.5)
+            self.ax.set_ylim(self.lim[0] - 0.5, self.lim[1] + 0.5)
+        else:
+            min_max_data = np.append(self_data[:, 0], self_data[:, 1])
+            min_max_data = np.append(min_max_data, other_data[:, 0])
+            min_max_data = np.append(min_max_data, other_data[:, 1])
+
+        self.ax.set_xlim(min_max_data.min(), min_max_data.max() + 1.)
+        self.ax.set_ylim(min_max_data.min(), min_max_data.max() + 1.)
+        gap = int(10 * (min_max_data.max() - min_max_data.min()) / 100)
+        tick_range = np.arange(min_max_data.min(), min_max_data.max(), gap, dtype=np.int64)
+
+        self.ax.set_xticks(tick_range + 0.5)
+        self.ax.set_xticklabels(tick_range)
+        self.ax.set_yticks(tick_range + 0.5)
+        self.ax.set_yticklabels(tick_range)
+
+        self.ax.set_xlabel('Residue number')
+        self.ax.set_ylabel('Residue number')
+
+    @staticmethod
+    def _determine_color(h):
+        """Determine the color of the contacts in order"""
+        greys = plt.get_cmap("Greys")
+        return [greys(contact.raw_score) for contact in h]

conkit/plot/figure.py

@@ -66,7 +66,7 @@ def __init__(self, ax=None, legend=True):
     def __repr__(self):
         return self.__class__.__name__
 
-    def _patch_scatter(self, x, y, facecolor="#ffffff", edgecolor="#000000", radius=0.5, linewidth=1.0):
+    def _patch_scatter(self, x, y, symbol="o", facecolor="#ffffff", edgecolor="#000000", radius=0.5, linewidth=1.0):
         """Draw scatter points as :obj:`Circles <matplotlib.pyplot.Circle>` to control width for discrete data"""
         if len(x) != len(y):
             raise ValueError("Unequal x and y data provided")
@@ -99,13 +99,22 @@ def _patch_scatter(self, x, y, facecolor="#ffffff", edgecolor="#000000", radius=
                 raise ValueError("Unequal x/y data and radii provided")
             r = radius
 
+        if symbol not in ["o", "s"]:
+            raise ValueError("Symbol needs to be circle (\"o\") or square (\"s\")")
+
         # Credits to https://stackoverflow.com/a/48174228/3046533
-        circles = [
-            plt.Circle((xi, yi), facecolor=fci, edgecolor=eci, radius=ri, linewidth=lwi)
-            for xi, yi, fci, eci, ri, lwi in zip(x, y, fc, ec, r, lw)
-        ]
-        if len(circles) > 0:
-            patch_collection = mcoll.PatchCollection(circles, match_original=True)
+        if symbol == "o":
+            patches = [
+                plt.Circle((xi, yi), facecolor=fci, edgecolor=eci, radius=ri, linewidth=lwi)
+                for xi, yi, fci, eci, ri, lwi in zip(x, y, fc, ec, r, lw)
+            ]
+        elif symbol == "s":
+            patches = [
+                plt.Rectangle((xi, yi), facecolor=fci, edgecolor=eci, height=ri, width=ri, linewidth=lwi)
+                for xi, yi, fci, eci, ri, lwi in zip(x, y, fc, ec, r, lw)
+            ]
+        if len(patches) > 0:
+            patch_collection = mcoll.PatchCollection(patches, match_original=True)
             self.ax.add_collection(patch_collection)
 
     def savefig(self, filename, dpi=300, overwrite=False):