skan.csr: CSR Graph representation of skeletons#

class skan.csr.JunctionModes(value)[source]#

Modes for cleaning up junctions in skeletons.

NONE:

Junctions are left as is. This is equivalent to unique_junctions=False in skan < 0.10.

Centroid:

Junctions are consolidated into the centroid of the contributing nodes. This is equivalent to unique_junctions=True in skan < 0.10.

MST:

Junctions are replaced with the minimum spanning tree. This is the new default in skan==0.10.0 and will become the only mode in skan>=0.11.

class skan.csr.Skeleton(skeleton_image, *, spacing=1, source_image=None, _buffer_size_offset=None, keep_images=True, junction_mode=<JunctionModes.MST: 'mst'>, unique_junctions=None)[source]#

Object to group together all the properties of a skeleton.

In the text below, we use the following notation:

  • N: the number of points in the pixel skeleton,

  • ndim: the dimensionality of the skeleton

  • P: the number of paths in the skeleton (also the number of links in the junction graph).

  • J: the number of junction nodes

  • Sd: the sum of the degrees of all the junction nodes

  • [Nt], [Np], Nr, Nc: the dimensions of the source image

Parameters

skeleton_image (array) – The input skeleton (1-pixel/voxel thick skeleton, all other values 0).

Other Parameters
  • spacing (float or array of float, shape (ndim,)) – The scale of the pixel spacing along each axis.

  • source_image (array of float, same shape as skeleton_image) – The image that skeleton_image represents / summarizes / was generated from. This is used to produce visualizations as well as statistical properties of paths.

  • keep_images (bool) – Whether or not to keep the original input images. These can be useful for visualization, but they may take up a lot of memory.

  • unique_junctions (bool, optional) – If True, adjacent junction nodes get collapsed into a single conceptual node, with position at the centroid of all the connected initial nodes.

graph#

The skeleton pixel graph, where each node is a non-zero pixel in the input image, and each edge connects adjacent pixels. The graph is represented as an adjacency matrix in SciPy sparse matrix format. For more information see the scipy.sparse documentation as well as scipy.sparse.csgraph. Note: pixel numbering starts at 1, so the shape of this matrix is (N + 1, N + 1) instead of (N, N).

Type

scipy.sparse.csr_matrix, shape (N + 1, N + 1)

nbgraph#

A thin Numba wrapper around the csr_matrix format, this provides faster graph methods. For example, it is much faster to get a list of neighbors, or test for the presence of a specific edge.

Type

NBGraph

coordinates#

skeleton_pixel_id i -> coordinates[i] The image coordinates of each pixel in the skeleton. Some values in this matrix are non-sensical — you should only access them from node ids.

Type

array, shape (N, ndim)

paths#

A csr_matrix where element [i, j] is on if node j is in path i. This includes path endpoints. The number of nonzero elements is N - J + Sd.

Type

scipy.sparse.csr_matrix, shape (P, N + 1)

n_paths#

The number of paths, P. This is redundant information given n_paths, but it is used often enough that it is worth keeping around.

Type

int

distances#

The distance of each path. Note: not initialized until path_lengths() is called on the skeleton; use path_lengths() instead

Type

array of float, shape (P,)

skeleton_image#

The input skeleton image. Only present if keep_images is True. Set to False to preserve memory.

Type

array or None

source_image#

The image from which the skeleton was derived. Only present if keep_images is True. This is useful for visualization.

Type

array or None

path(index)[source]#

Return the pixel indices of path number index.

Parameters

index (int) – The desired path.

Returns

path – The indices of the pixels belonging to the path, including endpoints.

Return type

array of int

path_coordinates(index)[source]#

Return the image coordinates of the pixels in the path.

Parameters

index (int) – The desired path.

Returns

path_coords – The (image) coordinates of points on the path, including endpoints.

Return type

array of float

path_label_image()[source]#

Image like self.skeleton_image with path_ids as values.

Returns

label_image – Image of the same shape as self.skeleton_image where each pixel has the value of its branch id + 1.

Return type

array of ints

path_lengths()[source]#

Return the length of each path on the skeleton.

Returns

lengths – The length of all the paths in the skeleton.

Return type

array of float

path_means()[source]#

Compute the mean pixel value along each path.

Returns

means – The average pixel value along each path in the skeleton.

Return type

array of float

path_stdev()[source]#

Compute the standard deviation of values along each path.

Returns

stdevs – The standard deviation of pixel values along each path.

Return type

array of float

path_with_data(index)[source]#

Return pixel indices and corresponding pixel values on a path.

Parameters

index (int) – The desired path.

Returns

  • path (array of int) – The indices of pixels on the path, including endpoints.

  • data (array of float) – The values of pixels on the path.

paths_list()[source]#

List all the paths in the skeleton, including endpoints.

Returns

paths – The list containing all the paths in the skeleton.

Return type

list of array of int

skan.csr.branch_statistics(graph, pixel_values=None, *, buffer_size_offset=0)[source]#

Compute the length and type of each branch in a skeleton graph.

Parameters
  • graph (sparse.csr_matrix, shape (N, N)) – A skeleton graph.

  • pixel_values (array of float, shape (N,)) – A value for each pixel in the graph. Used to compute total intensity statistics along each branch.

  • buffer_size_offset (int, optional) – The buffer size is given by the sum of the degrees of non-path nodes. This is usually 2x the amount needed, allowing room for extra cycles of path-only nodes. However, if the image consists only of such cycles, the buffer size will be 0, resulting in an error. Until a more sophisticated, expandable-buffer solution is implemented, you can manually set a bigger buffer size using this parameter.

Returns

branches – An array containing branch endpoint IDs, length, and branch type. The types are: - tip-tip (0) - tip-junction (1) - junction-junction (2) - path-path (3) (This can only be a standalone cycle) Optionally, the last column contains the average pixel value along each branch (not including the endpoints).

Return type

array of float, shape (N, {4, 5})

skan.csr.compute_centroids(image)[source]#

Find the centroids of all nonzero connected blobs in image.

Parameters

image (ndarray) – The input image.

Returns

label_image – The input image, with each connected region containing a different integer label.

Return type

ndarray of int

Examples

>>> image = np.array([[1, 0, 1, 0, 0, 1, 1],
...                   [1, 0, 0, 1, 0, 0, 0]])
>>> labels, centroids = compute_centroids(image)
>>> print(labels)
[[1 0 2 0 0 3 3]
 [1 0 0 2 0 0 0]]
>>> centroids
array([[0.5, 0. ],
       [0.5, 2.5],
       [0. , 5.5]])
skan.csr.make_degree_image(skeleton_image)[source]#

Create a array showing the degree of connectivity of each pixel.

Parameters

skeleton_image (array) – An input image in which every nonzero pixel is considered part of the skeleton, and links between pixels are determined by a full n-dimensional neighborhood.

Returns

degree_image – An image containing the degree of connectivity of each pixel in the skeleton to neighboring pixels.

Return type

array of int, same shape as skeleton_image

skan.csr.skeleton_to_csgraph(skel, *, spacing=1, value_is_height=False, junction_mode=<JunctionModes.MST: 'mst'>, unique_junctions=None)[source]#

Convert a skeleton image of thin lines to a graph of neighbor pixels.

Parameters
  • skel (array) – An input image in which every nonzero pixel is considered part of the skeleton, and links between pixels are determined by a full n-dimensional neighborhood.

  • spacing (float, or array-like of float, shape (skel.ndim,)) – A value indicating the distance between adjacent pixels. This can either be a single value if the data has the same resolution along all axes, or it can be an array of the same shape as skel to indicate spacing along each axis.

Other Parameters
  • value_is_height (bool, optional) – If True, the pixel value at each point of the skeleton will be considered to be a height measurement, and this height will be incorporated into skeleton branch lengths. Used for analysis of atomic force microscopy (AFM) images.

  • junction_mode (JunctionModes.{NONE,MST,CENTROID}) – If NONE, junction pixels are not collapsed. If MST, junction pixels are replaced by their minimum spanning tree, resulting in a single junction pixel. If CENTROID, junction pixels are collapsed to their centroid.

  • unique_junctions (bool, optional) – DEPRECATED: Use junction_mode=JunctionModes.Centroid to get behavior equivalent to unique_junctions=True. If True, adjacent junction nodes get collapsed into a single conceptual node, with position at the centroid of all the connected initial nodes.

Returns

  • graph (sparse.csr_matrix) – A graph of shape (Nnz + 1, Nnz + 1), where Nnz is the number of nonzero pixels in skel. The value graph[i, j] is the distance between adjacent pixels i and j. In a 2D image, that would be 1 for immediately adjacent pixels and sqrt(2) for diagonally adjacent ones.

  • pixel_coordinates (array of float) – An array of shape (Nnz + 1, skel.ndim), mapping indices in graph to pixel coordinates in degree_image or skel. Array entry (0,:) contains currently always zeros to index the pixels, which start at 1, directly to the coordinates.

skan.csr.submatrix(M, idxs)[source]#

Return the outer-index product submatrix, M[idxs, :][:, idxs].

Parameters
  • M (scipy.sparse.spmatrix) – Input (square) matrix

  • idxs (array of int) – The indices to subset. No index in idxs should exceed the number of rows of M.

Returns

Msub – The subsetted matrix.

Return type

scipy.sparse.spmatrix

Examples

>>> Md = np.arange(16).reshape((4, 4))
>>> M = sparse.csr_matrix(Md)
>>> print(submatrix(M, [0, 2]).toarray())
[[ 0  2]
 [ 8 10]]
skan.csr.summarise(image, *, spacing=1, using_height=False)[source]#

Compute statistics for every disjoint skeleton in image.

Note: this function is deprecated. Prefer summarize().

Parameters
  • image (array, shape (M, N, .., P)) – N-dimensional array, where nonzero entries correspond to an object’s single-pixel-wide skeleton. If the image is of type ‘float’, the values are taken to be the height at that pixel, which is used to compute the skeleton distances.

  • spacing (float, or array-like of float, shape (skel.ndim,)) – A value indicating the distance between adjacent pixels. This can either be a single value if the data has the same resolution along all axes, or it can be an array of the same shape as skel to indicate spacing along each axis.

  • using_height (bool, optional) – If True, the pixel value at each point of the skeleton will be considered to be a height measurement, and this height will be incorporated into skeleton branch lengths, endpoint coordinates, and euclidean distances. Used for analysis of atomic force microscopy (AFM) images.

Returns

df – A data frame summarising the statistics of the skeletons in image.

Return type

pandas DataFrame

skan.csr.summarize(skel: skan.csr.Skeleton, *, find_main_branch=False)[source]#

Compute statistics for every skeleton and branch in skel.

Parameters
  • skel (skan.csr.Skeleton) – A Skeleton object.

  • find_main_branch (bool, optional) – Whether to compute main branches. A main branch is defined as the longest shortest path within a skeleton. This step is very expensive as it involves computing the shortest paths between all pairs of branch endpoints, so it is off by default.

Returns

summary – A summary of the branches including branch length, mean branch value, branch euclidean distance, etc.

Return type

pandas.DataFrame