GLCM Texture Features

This example illustrates texture classification using texture classification using grey level co-occurrence matrices (GLCMs). A GLCM is a histogram of co-occurring greyscale values at a given offset over an image.

In this example, samples of two different textures are extracted from an image: grassy areas and sky areas. For each patch, a GLCM with a horizontal offset of 5 is computed. Next, two features of the GLCM matrices are computed: dissimilarity and correlation. These are plotted to illustrate that the classes form clusters in feature space.

In a typical classification problem, the final step (not included in this example) would be to train a classifier, such as logistic regression, to label image patches from new images.

../_images/plot_glcm_1.png

import matplotlib.pyplot as plt

from skimage.feature import greycomatrix, greycoprops
from skimage import data


PATCH_SIZE = 21

# open the camera image
image = data.camera()

# select some patches from grassy areas of the image
grass_locations = [(474, 291), (440, 433), (466, 18), (462, 236)]
grass_patches = []
for loc in grass_locations:
    grass_patches.append(image[loc[0]:loc[0] + PATCH_SIZE,
                               loc[1]:loc[1] + PATCH_SIZE])

# select some patches from sky areas of the image
sky_locations = [(54, 48), (21, 233), (90, 380), (195, 330)]
sky_patches = []
for loc in sky_locations:
    sky_patches.append(image[loc[0]:loc[0] + PATCH_SIZE,
                             loc[1]:loc[1] + PATCH_SIZE])

# compute some GLCM properties each patch
xs = []
ys = []
for i, patch in enumerate(grass_patches + sky_patches):
    glcm = greycomatrix(patch, [5], [0], 256, symmetric=True, normed=True)
    xs.append(greycoprops(glcm, 'dissimilarity')[0, 0])
    ys.append(greycoprops(glcm, 'correlation')[0, 0])

# create the figure
fig = plt.figure(figsize=(8, 8))

# display original image with locations of patches
ax = fig.add_subplot(3, 2, 1)
ax.imshow(image, cmap=plt.cm.gray, interpolation='nearest',
          vmin=0, vmax=255)
for (y, x) in grass_locations:
    ax.plot(x + PATCH_SIZE / 2, y + PATCH_SIZE / 2, 'gs')
for (y, x) in sky_locations:
    ax.plot(x + PATCH_SIZE / 2, y + PATCH_SIZE / 2, 'bs')
ax.set_xlabel('Original Image')
ax.set_xticks([])
ax.set_yticks([])
ax.axis('image')

# for each patch, plot (dissimilarity, correlation)
ax = fig.add_subplot(3, 2, 2)
ax.plot(xs[:len(grass_patches)], ys[:len(grass_patches)], 'go',
        label='Grass')
ax.plot(xs[len(grass_patches):], ys[len(grass_patches):], 'bo',
        label='Sky')
ax.set_xlabel('GLCM Dissimilarity')
ax.set_ylabel('GLVM Correlation')
ax.legend()

# display the image patches
for i, patch in enumerate(grass_patches):
    ax = fig.add_subplot(3, len(grass_patches), len(grass_patches)*1 + i + 1)
    ax.imshow(patch, cmap=plt.cm.gray, interpolation='nearest',
              vmin=0, vmax=255)
    ax.set_xlabel('Grass %d' % (i + 1))

for i, patch in enumerate(sky_patches):
    ax = fig.add_subplot(3, len(sky_patches), len(sky_patches)*2 + i + 1)
    ax.imshow(patch, cmap=plt.cm.gray, interpolation='nearest',
              vmin=0, vmax=255)
    ax.set_xlabel('Sky %d' % (i + 1))


# display the patches and plot
fig.suptitle('Grey level co-occurrence matrix features', fontsize=14)
plt.show()

STDOUT


        

STDERR


        

Python source code: download (generated using skimage 0.11dev)

IPython Notebook: download (generated using skimage 0.11dev)

aW1wb3J0IG1hdHBsb3RsaWIucHlwbG90IGFzIHBsdAoKZnJvbSBza2ltYWdlLmZlYXR1cmUgaW1wb3J0IGdyZXljb21hdHJpeCwgZ3JleWNvcHJvcHMKZnJvbSBza2ltYWdlIGltcG9ydCBkYXRhCgoKUEFUQ0hfU0laRSA9IDIxCgojIG9wZW4gdGhlIGNhbWVyYSBpbWFnZQppbWFnZSA9IGRhdGEuY2FtZXJhKCkKCiMgc2VsZWN0IHNvbWUgcGF0Y2hlcyBmcm9tIGdyYXNzeSBhcmVhcyBvZiB0aGUgaW1hZ2UKZ3Jhc3NfbG9jYXRpb25zID0gWyg0NzQsIDI5MSksICg0NDAsIDQzMyksICg0NjYsIDE4KSwgKDQ2MiwgMjM2KV0KZ3Jhc3NfcGF0Y2hlcyA9IFtdCmZvciBsb2MgaW4gZ3Jhc3NfbG9jYXRpb25zOgogICAgZ3Jhc3NfcGF0Y2hlcy5hcHBlbmQoaW1hZ2VbbG9jWzBdOmxvY1swXSArIFBBVENIX1NJWkUsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBsb2NbMV06bG9jWzFdICsgUEFUQ0hfU0laRV0pCgojIHNlbGVjdCBzb21lIHBhdGNoZXMgZnJvbSBza3kgYXJlYXMgb2YgdGhlIGltYWdlCnNreV9sb2NhdGlvbnMgPSBbKDU0LCA0OCksICgyMSwgMjMzKSwgKDkwLCAzODApLCAoMTk1LCAzMzApXQpza3lfcGF0Y2hlcyA9IFtdCmZvciBsb2MgaW4gc2t5X2xvY2F0aW9uczoKICAgIHNreV9wYXRjaGVzLmFwcGVuZChpbWFnZVtsb2NbMF06bG9jWzBdICsgUEFUQ0hfU0laRSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICBsb2NbMV06bG9jWzFdICsgUEFUQ0hfU0laRV0pCgojIGNvbXB1dGUgc29tZSBHTENNIHByb3BlcnRpZXMgZWFjaCBwYXRjaAp4cyA9IFtdCnlzID0gW10KZm9yIGksIHBhdGNoIGluIGVudW1lcmF0ZShncmFzc19wYXRjaGVzICsgc2t5X3BhdGNoZXMpOgogICAgZ2xjbSA9IGdyZXljb21hdHJpeChwYXRjaCwgWzVdLCBbMF0sIDI1Niwgc3ltbWV0cmljPVRydWUsIG5vcm1lZD1UcnVlKQogICAgeHMuYXBwZW5kKGdyZXljb3Byb3BzKGdsY20sICdkaXNzaW1pbGFyaXR5JylbMCwgMF0pCiAgICB5cy5hcHBlbmQoZ3JleWNvcHJvcHMoZ2xjbSwgJ2NvcnJlbGF0aW9uJylbMCwgMF0pCgojIGNyZWF0ZSB0aGUgZmlndXJlCmZpZyA9IHBsdC5maWd1cmUoZmlnc2l6ZT0oOCwgOCkpCgojIGRpc3BsYXkgb3JpZ2luYWwgaW1hZ2Ugd2l0aCBsb2NhdGlvbnMgb2YgcGF0Y2hlcwpheCA9IGZpZy5hZGRfc3VicGxvdCgzLCAyLCAxKQpheC5pbXNob3coaW1hZ2UsIGNtYXA9cGx0LmNtLmdyYXksIGludGVycG9sYXRpb249J25lYXJlc3QnLAogICAgICAgICAgdm1pbj0wLCB2bWF4PTI1NSkKZm9yICh5LCB4KSBpbiBncmFzc19sb2NhdGlvbnM6CiAgICBheC5wbG90KHggKyBQQVRDSF9TSVpFIC8gMiwgeSArIFBBVENIX1NJWkUgLyAyLCAnZ3MnKQpmb3IgKHksIHgpIGluIHNreV9sb2NhdGlvbnM6CiAgICBheC5wbG90KHggKyBQQVRDSF9TSVpFIC8gMiwgeSArIFBBVENIX1NJWkUgLyAyLCAnYnMnKQpheC5zZXRfeGxhYmVsKCdPcmlnaW5hbCBJbWFnZScpCmF4LnNldF94dGlja3MoW10pCmF4LnNldF95dGlja3MoW10pCmF4LmF4aXMoJ2ltYWdlJykKCiMgZm9yIGVhY2ggcGF0Y2gsIHBsb3QgKGRpc3NpbWlsYXJpdHksIGNvcnJlbGF0aW9uKQpheCA9IGZpZy5hZGRfc3VicGxvdCgzLCAyLCAyKQpheC5wbG90KHhzWzpsZW4oZ3Jhc3NfcGF0Y2hlcyldLCB5c1s6bGVuKGdyYXNzX3BhdGNoZXMpXSwgJ2dvJywKICAgICAgICBsYWJlbD0nR3Jhc3MnKQpheC5wbG90KHhzW2xlbihncmFzc19wYXRjaGVzKTpdLCB5c1tsZW4oZ3Jhc3NfcGF0Y2hlcyk6XSwgJ2JvJywKICAgICAgICBsYWJlbD0nU2t5JykKYXguc2V0X3hsYWJlbCgnR0xDTSBEaXNzaW1pbGFyaXR5JykKYXguc2V0X3lsYWJlbCgnR0xWTSBDb3JyZWxhdGlvbicpCmF4LmxlZ2VuZCgpCgojIGRpc3BsYXkgdGhlIGltYWdlIHBhdGNoZXMKZm9yIGksIHBhdGNoIGluIGVudW1lcmF0ZShncmFzc19wYXRjaGVzKToKICAgIGF4ID0gZmlnLmFkZF9zdWJwbG90KDMsIGxlbihncmFzc19wYXRjaGVzKSwgbGVuKGdyYXNzX3BhdGNoZXMpKjEgKyBpICsgMSkKICAgIGF4Lmltc2hvdyhwYXRjaCwgY21hcD1wbHQuY20uZ3JheSwgaW50ZXJwb2xhdGlvbj0nbmVhcmVzdCcsCiAgICAgICAgICAgICAgdm1pbj0wLCB2bWF4PTI1NSkKICAgIGF4LnNldF94bGFiZWwoJ0dyYXNzICVkJyAlIChpICsgMSkpCgpmb3IgaSwgcGF0Y2ggaW4gZW51bWVyYXRlKHNreV9wYXRjaGVzKToKICAgIGF4ID0gZmlnLmFkZF9zdWJwbG90KDMsIGxlbihza3lfcGF0Y2hlcyksIGxlbihza3lfcGF0Y2hlcykqMiArIGkgKyAxKQogICAgYXguaW1zaG93KHBhdGNoLCBjbWFwPXBsdC5jbS5ncmF5LCBpbnRlcnBvbGF0aW9uPSduZWFyZXN0JywKICAgICAgICAgICAgICB2bWluPTAsIHZtYXg9MjU1KQogICAgYXguc2V0X3hsYWJlbCgnU2t5ICVkJyAlIChpICsgMSkpCgoKIyBkaXNwbGF5IHRoZSBwYXRjaGVzIGFuZCBwbG90CmZpZy5zdXB0aXRsZSgnR3JleSBsZXZlbCBjby1vY2N1cnJlbmNlIG1hdHJpeCBmZWF0dXJlcycsIGZvbnRzaXplPTE0KQpwbHQuc2hvdygp