levels parameter in greycomatrix scikit-image python

Question

I'm moving my Matlab image processing algorithms to Python using scikit-image tools, and I'm calculating the gray level co-occurrence matrix (GLCM) using greycomatrix. I have a problem if the parameter levels is lesser than the maximum value of the intensity image (image.max()). For instance:

import numpy as np
from skimage.feature import greycomatrix
image = np.array([[0, 0, 1, 1],[0, 0, 1, 1],[0, 2, 2, 2],[2, 2, 3, 3]], dtype=np.uint8)
result = greycomatrix(image, distances = [1], angles = [0], levels = 4, symmetric=True)

The output is:

glcm = result[:,:,0,0]

array([[4, 2, 1, 0],
   [2, 4, 0, 0],
   [1, 0, 6, 1],
   [0, 0, 1, 2]], dtype=uint32)

which is correct, a 4x4 matrix. But if levels=3, I can't calculate the GLCM, and the error is:

result = greycomatrix(image, distances = [1], angles = [0], levels = 3, symmetric=True)

Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/usr/lib/python3.4/site-packages/skimage/feature/texture.py", line 97, in greycomatrix
assert image.max() < levels
AssertionError

And of course ... I get the error, but I should be able to calculate a GLCM (3x3 matrix) with levels lesser than image.max(). For instance, for:

result = greycomatrix(image, distances = [1], angles = [0], levels = 3, symmetric=True)

I should get the following GLCM (I can do it in Matlab):

4     3     0
3    10     1
0     1     2

When I work with huge images I reduce the the levels of the GLCM in order to reduce the calculation time. Is there any problem with the greycomatrix or I'm thinking wrong? Thanks in advance.

I think the only reason is that no one's asked for it before. The three options are: 1) discuss it on the scikit-image mailing list 2) file an issue against scikit-image on GitHub or 3) submit a pull request that fixes the problem. — Stefan van der Walt, Jul 22 '14 at 19:53

Tonechas · Answer 1 · 2017-06-05T19:09:01.983

Maybe I'm late to the party but hopefully my answer might be useful for someone else in the future...

According to Matlab documentation 'NumLevels' can be less than max(image(:)) because grayscale values are clipped:

Grayscale values less than or equal to low are scaled to 1. Grayscale values greater than or equal to high are scaled to 'NumLevels'.

You could use the same workaround in Python through NumPy's clip:

In [11]: import numpy as np

In [12]: from skimage.feature import greycomatrix

In [13]: image = np.array([[0, 0, 1, 1], 
    ...:                   [0, 0, 1, 1], 
    ...:                   [0, 2, 2, 2], 
    ...:                   [2, 2, 3, 3]], dtype=np.uint8)
    ...: 

In [14]: clipped = np.clip(image, a_min=0, a_max=2)

In [15]: clipped
Out[15]: 
array([[0, 0, 1, 1],
       [0, 0, 1, 1],
       [0, 2, 2, 2],
       [2, 2, 2, 2]], dtype=uint8)

In [16]: result = greycomatrix(clipped, 
    ...:                       distances=[1], 
    ...:                       angles=[0], 
    ...:                       levels=3, 
    ...:                       symmetric=True)
    ...: 

In [17]: result[:, :, 0, 0]
Out[17]: 
array([[ 4,  2,  1],
       [ 2,  4,  0],
       [ 1,  0, 10]], dtype=uint32)

If you wish to reduce the size of the gray level co-ocurrence matrix I would recommend requantizing the image rather than clipping the grayscale values. This approach can be easily implemented through NumPy's digitize:

In [59]: nlevels = 256

In [60]: nbins = 64

In [61]: arr = np.arange(nlevels).reshape((16, 16)).astype(np.uint8).T

In [62]: np.set_printoptions(threshold=300, linewidth=100)

In [63]: arr
Out[63]: 
array([[  0,  16,  32,  48,  64,  80,  96, 112, 128, 144, 160, 176, 192, 208, 224, 240],
       [  1,  17,  33,  49,  65,  81,  97, 113, 129, 145, 161, 177, 193, 209, 225, 241],
       [  2,  18,  34,  50,  66,  82,  98, 114, 130, 146, 162, 178, 194, 210, 226, 242],
       [  3,  19,  35,  51,  67,  83,  99, 115, 131, 147, 163, 179, 195, 211, 227, 243],
       [  4,  20,  36,  52,  68,  84, 100, 116, 132, 148, 164, 180, 196, 212, 228, 244],
       [  5,  21,  37,  53,  69,  85, 101, 117, 133, 149, 165, 181, 197, 213, 229, 245],
       [  6,  22,  38,  54,  70,  86, 102, 118, 134, 150, 166, 182, 198, 214, 230, 246],
       [  7,  23,  39,  55,  71,  87, 103, 119, 135, 151, 167, 183, 199, 215, 231, 247],
       [  8,  24,  40,  56,  72,  88, 104, 120, 136, 152, 168, 184, 200, 216, 232, 248],
       [  9,  25,  41,  57,  73,  89, 105, 121, 137, 153, 169, 185, 201, 217, 233, 249],
       [ 10,  26,  42,  58,  74,  90, 106, 122, 138, 154, 170, 186, 202, 218, 234, 250],
       [ 11,  27,  43,  59,  75,  91, 107, 123, 139, 155, 171, 187, 203, 219, 235, 251],
       [ 12,  28,  44,  60,  76,  92, 108, 124, 140, 156, 172, 188, 204, 220, 236, 252],
       [ 13,  29,  45,  61,  77,  93, 109, 125, 141, 157, 173, 189, 205, 221, 237, 253],
       [ 14,  30,  46,  62,  78,  94, 110, 126, 142, 158, 174, 190, 206, 222, 238, 254],
       [ 15,  31,  47,  63,  79,  95, 111, 127, 143, 159, 175, 191, 207, 223, 239, 255]], dtype=uint8)

In [64]: binned = np.uint8(np.digitize(arr, np.arange(0, nlevels, nbins))) - 1

In [65]: binned
Out[65]: 
array([[0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3],
       [0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3]], dtype=uint8)

In [66]: glcm = greycomatrix(binned, 
    ...:                     distances=[1], 
    ...:                     angles=[0], 
    ...:                     levels=binned.max()+1, 
    ...:                     symmetric=True)
    ...: 

In [67]: glcm[:, :, 0, 0]    
Out[67]: 
array([[96, 16,  0,  0],
       [16, 96, 16,  0],
       [ 0, 16, 96, 16],
       [ 0,  0, 16, 96]], dtype=uint32)

levels parameter in greycomatrix scikit-image python

1 Answers1