Python：采用非矩形区域的 GLCM答案

【问题标题】：Python: taking the GLCM of a non-rectangular regionPython：采用非矩形区域的 GLCM
【发布时间】：2017-04-03 20:09:06
【问题描述】：

我一直在使用 skimage 的 SLIC 实现来分割超像素中的图像。我想使用 GLCM 从这些超像素中提取额外的特征来解决分类问题。这些超像素不是矩形的。在 MATLAB 中，您可以将像素设置为 NaN，算法将忽略它们 (link)。我可以使用它在超像素周围制作边界框，然后将未使用的像素设置为 NaN。

但是，skimage 中的 greycomatrix 函数与 MATLAB 实现的工作方式并不完全相同。将像素设置为 NaN 时，该函数无法通过断言检查所有值是否大于 0。

是否有可以使用非矩形 ROI 的 Python 实现？

【问题讨论】：

标签： python image-processing scikit-image mahotas glcm

【解决方案1】：

虽然mahotas 也是一个出色的计算机视觉库，但没有必要停止使用skimage 来做到这一点。

正如@Tonechas 指出的那样，需要将这些 NaN 值设置为整数，因为 np.nan 的类型为 float，而 greycomatrix 函数需要一个数组整数。

最简单的选择是将那些NaN 设置为零，但是，如果您的像素中已经有零值并且不想混合它们，您可以选择任何其他常量。之后，您所要做的就是从 GLCM 中过滤掉选择的值（再一次，通常为零）。

要了解这意味着什么，让我们看看skimage 告诉我们greycomatrix function 的输出：

4-D 数组

[...] 值 P[i,j,d,theta] 是灰度 j 在距离 d 和与灰度 i 成角度 theta 处出现的次数。如果 normed 为 False，则输出类型为 uint32，否则为 float64。维度是：级别 x 级别 x 距离数 x 角度数。

换句话说，数组的前两个维度定义了一个矩阵，它告诉我们两个不同的值有多少次相距一定距离。请注意，GLCM 确实不保持输入数组的形状。这些行和列告诉我们这些值是如何关联的。

知道了这一点，就很容易过滤掉 ROI 之外的值（假设我们已将这些 NaN 设置为零）：

glcm = greycomatrix(img, [1], [0])  # Calculate the GLCM "one pixel to the right"
filt_glcm = glcm[1:, 1:, :, :]           # Filter out the first row and column

现在您可以轻松计算过滤后的 GLCM 的 Haralick 属性。例如：

greycoprops(filt_glcm, prop='contrast')

【讨论】：

【解决方案2】：

问题是您必须将整数数组传递给greycomatrix，但np.nan 的类型为float（详情请查看this thread）。因此，您无法将 ROI 之外的像素编码为 NaN。

处理非矩形 ROI 的近似解决方法是将 ROI 外的像素设置为 0 并使用 mahotas 库中的函数 haralick。该函数返回从四个不同的 GLCM 中提取的 13 个 Haralick 特征，对应于四个二维方向和距离参数的特定值。

来自文档：

ignore_zeros 可用于让函数忽略任何零值像素（作为背景）。

总而言之，您需要屏蔽掉那些落在 ROI 之外的像素，并在对haralick 的调用中将ignore_zeros 设置为True。

演示

首先，让我们生成一些模拟数据：

In [213]: import numpy as np

In [214]: shape = (3, 4)

In [215]: levels = 8

In [216]: np.random.seed(2017)

In [217]: x = np.random.randint(0, levels, size=shape)

In [218]: x
Out[218]: 
array([[3, 1, 6, 5],
       [2, 0, 2, 2],
       [3, 7, 7, 7]])

然后我们必须从图像中删除所有零，因为在这种方法中，零强度级别是为 ROI 之外的像素保留的。值得指出的是，将强度 0 和 1 合并为单个强度 1 会导致结果不准确。

In [219]: x[x == 0] = 1

In [220]: x
Out[220]: 
array([[3, 1, 6, 5],
       [2, 1, 2, 2],
       [3, 7, 7, 7]])

下一步包括为 ROI 之外的像素（在此玩具示例中为图像的四个角）定义掩码并将这些像素设置为 0。

In [221]: non_roi = np.zeros(shape=shape, dtype=np.bool)

In [222]: non_roi[np.ix_([0, -1], [0, -1])] = True

In [223]: non_roi
Out[223]: 
array([[ True, False, False,  True],
       [False, False, False, False],
       [ True, False, False,  True]], dtype=bool)

In [224]: x[non_roi] = 0

In [225]: x
Out[225]: 
array([[0, 1, 6, 0],
       [2, 1, 2, 2],
       [0, 7, 7, 0]])

我们现在可以从非矩形 ROI 的 GLCM 中提取特征：

In [226]: import mahotas.features.texture as mht

In [227]: features = mht.haralick(x, ignore_zeros=True)

In [228]: features.size
Out[228]: 52

In [229]: features.ravel()
Out[229]: array([ 0.18      ,  5.4       ,  0.5254833 , ...,  0.81127812,
       -0.68810414,  0.96300727])

检查共现矩阵的外观可能很有用。例如，“右侧像素”GLCM 将是：

In [230]: mht.cooccurence(x, 0)
Out[230]: 
array([[0, 1, 0, ..., 0, 1, 2],
       [1, 0, 2, ..., 0, 1, 0],
       [0, 2, 2, ..., 0, 0, 0],
       ..., 
       [0, 0, 0, ..., 0, 0, 0],
       [1, 1, 0, ..., 0, 0, 0],
       [2, 0, 0, ..., 0, 0, 2]])

【讨论】：