Skimage.graph: knight's move neighbourhood

[zoran-cuckovic]

Description:

Currently, the skimage.graph module operates on immediate adjacency neighbourhood, i.e. the graph connects pixels that touch each other. However, for cost surface calculation, this may introduce artefacts (image below). A possible solution is to extend the neighbourhood to the distance of two pixels, a.k.a. the "knight's move":

Octagonal shaped cost surface artefacts (lines correspond to isochrones):

The images are from https://grass.osgeo.org/grass83/manuals/r.cost.html (GRASS).

Implementing this feature would improve the usability of the module for spatial analysis, which has a very large community of potential users (ecology, geography, urban planning, etc.).

[lagru]

Thanks for the suggestion. At a first glance this sounds good but I'm unclear on the details.

• To implement this, would we need to modify skimage.graph.RAG such that the graph will contain these extra edges?
• Do you have an idea how an intuitive API would look like?

[zoran-cuckovic]

Hello,
Yes, perhaps the most intuitive solution would be to add the third option to the connectivity parameter (connectivity = 3) ?
I've checked the scipy.ndimage.generate_binary_structure function, which is mentioned in the manual for skimage.graph.RAG, and which has a connectivity parameter. It seems that only the immediate neighbourhood is taken into account; specifying 3 or more for a 2D matrix yields only a 3x3 window. This implies that the skimage.graph connectivity would behave differently than the ndimage connectivity. I don't know whether this may introduce confusion or not (?).
Thanks for your interest @lagru !

[lagru]

connectivity specifying the neighborhood as a squared distance from the footprint's center is pretty established. So deviating from that definition, so that connectivity=3 means use "knight's move neighborhood" would likely make for a very confusing API.

That said, we already have a convention for functions to specify a more complex arbitrary footprint: a footprint parameter taking a simple binary array which is True for connected pixels. That would make it possible to define a footprint like you suggest. We could even consider an helper function to generate such a footprint. :)

Do you have a concrete small example using scikit-image, where we could evaluate the effects of such a change? That would probably be very useful.

[zoran-cuckovic]

Hello,
Here is a snippet which should demonstrate the problem when using flat (or near flat) values:

import numpy as np
from skimage import graph
import matplotlib.pyplot as plt

fig, axs = plt.subplots(1, 3, figsize=(50, 3))

# Create a flat image, all ones
img = np.ones(2500).reshape(50,50)

# Introduce some random noise (half of the image), 
# to make one part more natural
img[25 :, : ] += np.random.rand(25, 50)

# Show the raw image,
axs[0].imshow(img,cmap='gray_r')

# Create a graph with diagonal correction (MCP_Geopetric)
lg = graph.MCP_Geometric (img )
lcd, traceback = lg.find_costs(starts=[(25,25)])

# Show accumulated cost surface
axs[1].imshow(lcd,cmap='gray_r')

# Round values to reveal the artefacts in the UPPER PART
# Lower part of the image, having random values, looks more natural
axs[2].imshow(np.round(lcd/5), cmap='gray_r') 

plt.show()