Merge branch 'master' of git://github.com/scikit-image/scikit-image

doc/examples/segmentation/plot_rag_boundary.py

@@ -0,0 +1,27 @@
+"""
+==========================
+Region Boundary based RAGs
+==========================
+
+This example demonstrates construction of region boundary based RAGs with the
+`rag_boundary` function.
+"""
+from skimage.future import graph
+from skimage import data, segmentation, color, filters, io
+from skimage.util.colormap import viridis
+
+
+img = data.coffee()
+gimg = color.rgb2gray(img)
+
+labels = segmentation.slic(img, compactness=30, n_segments=400)
+edges = filters.sobel(gimg)
+edges_rgb = color.gray2rgb(edges)
+
+g = graph.rag_boundary(labels, edges)
+
+out = graph.draw_rag(labels, g, edges_rgb, node_color="#999999",
+                     colormap=viridis)
+
+io.imshow(out)
+io.show()

skimage/future/graph/__init__.py

@@ -1,5 +1,5 @@
 from .graph_cut import cut_threshold, cut_normalized
-from .rag import rag_mean_color, RAG, draw_rag
+from .rag import rag_mean_color, RAG, draw_rag, rag_boundary
 from .graph_merge import merge_hierarchical
 ncut = cut_normalized
 
@@ -9,4 +9,5 @@ __all__ = ['rag_mean_color',
            'ncut',
            'draw_rag',
            'merge_hierarchical',
+           'rag_boundary',
            'RAG']

skimage/future/graph/rag.py

@@ -2,6 +2,7 @@ import networkx as nx
 import numpy as np
 from numpy.lib.stride_tricks import as_strided
 from scipy import ndimage as ndi
+from scipy import sparse
 import math
 from ... import draw, measure, segmentation, util, color
 try:
@@ -11,6 +12,43 @@ except ImportError:
     pass
 
 
+def _edge_generator_from_csr(csr_matrix):
+    """Yield weighted edge triples for use by NetworkX from a CSR matrix.
+
+    This function is a straight rewrite of
+    `networkx.convert_matrix._csr_gen_triples`. Since that is a private
+    function, it is safer to include our own here.
+
+    Parameters
+    ----------
+    csr_matrix : scipy.sparse.csr_matrix
+        The input matrix. An edge (i, j, w) will be yielded if there is a
+        data value for coordinates (i, j) in the matrix, even if that value
+        is 0.
+
+    Yields
+    ------
+    i, j, w : (int, int, float) tuples
+        Each value `w` in the matrix along with its coordinates (i, j).
+
+    Examples
+    --------
+
+    >>> dense = np.eye(2, dtype=np.float)
+    >>> csr = sparse.csr_matrix(dense)
+    >>> edges = _edge_generator_from_csr(csr)
+    >>> list(edges)
+    [(0, 0, 1.0), (1, 1, 1.0)]
+    """
+    nrows = csr_matrix.shape[0]
+    values = csr_matrix.data
+    indptr = csr_matrix.indptr
+    col_indices = csr_matrix.indices
+    for i in range(nrows):
+        for j in range(indptr[i], indptr[i + 1]):
+            yield i, col_indices[j], values[j]
+
+
 def min_weight(graph, src, dst, n):
     """Callback to handle merging nodes by choosing minimum weight.
 
@@ -311,6 +349,70 @@ def rag_mean_color(image, labels, connectivity=2, mode='distance',
     return graph
 
 
+def rag_boundary(labels, edge_map, connectivity=2):
+    """ Comouter RAG based on region boundaries
+
+    Given an image's initial segmentation and its edge map this method
+    constructs the corresponding Region Adjacency Graph (RAG). Each node in the
+    RAG represents a set of pixels within the image with the same label in
+    `labels`. The weight between two adjacent regions is the average value
+    in `edge_map` along their boundary.
+
+    labels : ndarray
+        The labelled image.
+    edge_map : ndarray
+        This should have the same shape as that of `labels`. For all pixels
+        along the boundary between 2 adjacent regions, the average value of the
+        corresponding pixels in `edge_map` is the edge weight between them.
+    connectivity : int, optional
+        Pixels with a squared distance less than `connectivity` from each other
+        are considered adjacent. It can range from 1 to `labels.ndim`. Its
+        behavior is the same as `connectivity` parameter in
+        `scipy.ndimage.filters.generate_binary_structure`.
+
+    Examples
+    --------
+    >>> from skimage import data, segmentation, filters, color
+    >>> from skimage.future import graph
+    >>> img = data.chelsea()
+    >>> labels = segmentation.slic(img)
+    >>> edge_map = filters.sobel(color.rgb2gray(img))
+    >>> rag = graph.rag_boundary(labels, edge_map)
+
+    """
+
+    conn = ndi.generate_binary_structure(labels.ndim, connectivity)
+    eroded = ndi.grey_erosion(labels, footprint=conn)
+    dilated = ndi.grey_dilation(labels, footprint=conn)
+    boundaries0 = (eroded != labels)
+    boundaries1 = (dilated != labels)
+    labels_small = np.concatenate((eroded[boundaries0], labels[boundaries1]))
+    labels_large = np.concatenate((labels[boundaries0], dilated[boundaries1]))
+    n = np.max(labels_large) + 1
+
+    # use a dummy broadcast array as data for RAG
+    ones = as_strided(np.ones((1,), dtype=np.float), shape=labels_small.shape,
+                      strides=(0,))
+    count_matrix = sparse.coo_matrix((ones, (labels_small, labels_large)),
+                                     dtype=np.int_, shape=(n, n)).tocsr()
+    data = np.concatenate((edge_map[boundaries0], edge_map[boundaries1]))
+
+    data_coo = sparse.coo_matrix((data, (labels_small, labels_large)))
+    graph_matrix = data_coo.tocsr()
+    graph_matrix.data /= count_matrix.data
+
+    rag = RAG()
+    rag.add_weighted_edges_from(_edge_generator_from_csr(graph_matrix),
+                                weight='weight')
+    rag.add_weighted_edges_from(_edge_generator_from_csr(count_matrix),
+                                weight='count')
+
+    for n in rag.nodes():
+        rag.node[n].update({'labels': [n]})
+
+    return rag
+
+
 def draw_rag(labels, rag, img, border_color=None, node_color='#ffff00',
              edge_color='#00ff00', colormap=None, thresh=np.inf,
              desaturate=False, in_place=True):

skimage/future/graph/tests/test_rag.py

@@ -196,3 +196,23 @@ def test_generic_rag_3d():
     assert h.has_edge(0, 1) and h.has_edge(0, 3) and not h.has_edge(0, 7)
     k = graph.RAG(labels, connectivity=3)
     assert k.has_edge(0, 1) and k.has_edge(1, 2) and k.has_edge(2, 5)
+
+
+def test_rag_boundary():
+    labels = np.zeros((16, 16), dtype='uint8')
+    edge_map = np.zeros_like(labels, dtype=float)
+
+    edge_map[8, :] = 0.5
+    edge_map[:, 8] = 1.0
+
+    labels[:8, :8] = 1
+    labels[:8, 8:] = 2
+    labels[8:, :8] = 3
+    labels[8:, 8:] = 4
+
+    g = graph.rag_boundary(labels, edge_map, connectivity=1)
+    assert set(g.nodes()) == set([1, 2, 3, 4])
+    assert set(g.edges()) == set([(1, 2), (1, 3), (2, 4), (3, 4)])
+    assert g[1][3]['weight'] == 0.25
+    assert g[2][4]['weight'] == 0.34375
+    assert g[1][3]['count'] == 16