参考链接：图像纹理——灰度共生矩阵
知乎——提取图像的颜色、纹理特征(传统算法)

灰度共生矩阵（Gray-level Co-occurrence Matrix，GLCM）原理

灰度共生矩阵可反映灰度 值 和 空间 分布情况。
共生矩阵 $P$ 的描述方法：
规定一个方向（如水平，垂直，对角线）和距离（一个像素、两个像素），矩阵中 $P_{ij}$ 的值由灰度为 $i$ 和 $j$ 的像素对在该方向和距离上出现的次数除以N得到（归一化），N为对P有贡献的像素对总数。
共生矩阵的大小为 $L\ast L$ ，$L$ 是图像灰度级数目。

例子

在图像中任意一点$(x,y)$及偏离它的一点$(x+a,y+b)$ (其中 $a，b$ 为整数) 构成像素点对。

a为原灰度图像，灰度级数为16。为了表示方便，将灰度范围减小到4级，即为下图中样式

若取 a=1，b=0，即统计水平方向上相距为1像素的像素对出现的个数，这里的坐标系为

图b中黄色部分为图e中像素值对（1,0）的个数值


这是图e中（1,2）位置所对应的像素对

同理，图f为左下角方向一个位置相应灰度值对应的矩阵；图g为右侧方向两个像素对应的矩阵。

对比可以看出，图e中(0,1),(1,2),(2,3)和(3,0)均有较高的出现频数。图b表明，图像中存在明显的左上右下方向的纹理。

a，b取值较小对应于变化缓慢（粗而规则）的纹理图像，其灰度共生矩阵对角线上的数值较大。纹理的变化越快（细而规则），则对角线上的数值越小，而对角线两侧的值增大。

灰度共生矩阵特征

由于灰度共生矩阵的数据量较大，一般不直接作为区分纹理的特征，而是基于它构建的一些统计量作为纹理分类特征。Haralick曾提出了14种基于灰度共生矩阵计算出来的统计量：即：能量、熵、对比度、均匀性、相关性、方差、和平均、和方差、和熵、差方差、差平均、差熵、相关信息测度以及最大相关系数。

对比度 Contrast

度量矩阵的值是如何分布和图像中局部变化的多少，反应了图像的清晰度和纹理的沟纹深浅。纹理的沟纹越深，反差越大，效果越清晰；反之，对比值小，则沟纹浅，效果模糊。
反映某像素值及其领域像素值亮度的对比情况，图像亮度值变化快，换句话说纹理较深，它的对比度就越大，也就是它的纹理越清晰。

能量 Energy

能量变换反映了图像灰度分布均匀程度和纹理粗细度。若灰度共生矩阵的元素值相近，则能量较小，表示纹理细致；若其中一些值大，而其它值小，则能量值较大。能量值大表明一种较均一和规则变化的纹理模式。

熵 Entropy

图像包含信息量的随机性度量。当共生矩阵中所有值均相等或者像素值表现出最大的随机性时，熵最大；因此熵值表明了图像灰度分布的复杂程度，熵值越大，图像越复杂。

同质性/逆差距 Homogeneity

逆方差反映了图像纹理局部变化的大小，若图像纹理的不同区域间较均匀，变化缓慢，逆方差会较大，反之较小。与对比度或相异性相反，同质性的权重随着元素值与对角线的距离而减小，其减小方式是指数形式的。

相关性 Correlation

用来度量图像的灰度级在行或列方向上的相似程度，因此值的大小反应了局部灰度相关性，值越大，相关性也越大。

MATLAB代码

glcm = graycomatrix(I)

​glcms = graycomatrix(I,param1,val1,param2,val2,...)

​[glcms,SI] = graycomatrix(...)

关注的像素与其邻点之间的距离通过Offsets来进行调整，比如[0 1]代表是水平方向，[-1 1]代表是右上角45度方向，[-1 0]代表是竖直方向，即90度方向，而[-1 -1]则代表是左上角，即135度方向

这里的坐标系为水平向右为Y轴正方向，垂直向下为X轴正方向。

图像处理库scikits-image已经支持计算灰度共生矩阵和提取GLCM的纹理属性contrast、dissimilarity、homogeneity、ASM、energy、correlation

首先了解一下灰度共生矩阵是什么，下面介绍摘自百度百科。https://baike.baidu.com/item/%E7%81%B0%E5%BA%A6%E5%85%B1%E7%94%9F%E7%9F%A9%E9%98%B5

灰度共生矩阵，指的是一种通过研究灰度的空间相关特性来描述纹理的常用方法。1973年Haralick等人提出了用灰度共生矩阵来描述纹理特征。
由于纹理是由灰度分布在空间位置上反复出现而形成的，因而在图像空间中相隔某距离的两像素之间会存在一定的灰度关系，即图像中灰度的空间相关特性。
取图像(N×N)中任意一点 （x，y）及偏离它的另一点 （x+a，y+b），设该点对的灰度值为 （g1，g2）。令点（x，y） 在整个画面上移动，则会得到各种 （g1，g2）值，设灰度值的级数为 k，则（g1，g2） 的组合共有 k 的平方种。对于整个画面，统计出每一种 （g1，g2）值出现的次数，然后排列成一个方阵，再用（g1，g2） 出现的总次数将它们归一化为出现的概率P（g1，g2） ，这样的方阵称为灰度共生矩阵。距离差分值（a，b） 取不同的数值组合，可以得到不同情况下的联合概率矩阵。（a，b） 取值要根据纹理周期分布的特性来选择，对于较细的纹理，选取（1，0）、（1，1）、（2，0）等小的差分值。
当 a=1，b=0时，像素对是水平的，即0度扫描；当a=0，b=1 时，像素对是垂直的，即90度扫描；当 a=1，b=1时，像素对是右对角线的，即45度扫描；当 a=-1，b=1时，像素对是左对角线，即135度扫描。
这样，两个象素灰度级同时发生的概率，就将 （x，y）的空间坐标转化为“灰度对” （g1，g2）的描述，形成了灰度共生矩阵。

灰度共生矩阵的特征
对比度（contrast）

直接反映了某个像素值及其领域像素值的亮度的对比情况。如果偏离对角线的元素有较大值，即图像亮度值变化很快，则CON会有较大取值，这也符合对比度的定义。其中 。反映了图像的清晰度和纹理沟纹深浅的程度。纹理沟纹越深，其对比度越大，视觉效果越清晰；反之，对比度小，则沟纹浅，效果模糊。灰度差即对比度大的象素对越多，这个值越大。灰度公生矩阵中远离对角线的元素值越大，CON越大。

ASM 能量（angular second moment)

也即每个矩阵元素的平方和。
如果灰度共生矩阵中的值集中在某一块（比如对连续灰度值图像，值集中在对角线；对结构化的图像，值集中在偏离对角线的位置），则ASM有较大值，若G中的值分布较均匀（如噪声严重的图像），则ASM有较小的值。
能量是灰度共生矩阵元素值的平方和，所以也称能量，反映了图像灰度分布均匀程度和纹理粗细度。如果共生矩阵的所有值均相等，则ASM值小；相反，如果其中一些值大而其它值小，则ASM值大。当共生矩阵中元素集中分布时，此时ASM值大。ASM值大表明一种较均一和规则变化的纹理模式。

自相关（correlation）

其中

自相关反应了图像纹理的一致性。如果图像中有水平方向纹理，则水平方向矩阵的COR大于其余矩阵的COR值。它度量空间灰度共生矩阵元素在行或列方向上的相似程度，因此，相关值大小反映了图像中局部灰度相关性。当矩阵元素值均匀相等时，相关值就大;相反，如果矩阵像元值相差很大则相关值小。

百度百科对应的只有上述三个特征的公式，下面的另外三个特征的公式来自http://tonysyu.github.io/scikit-image/api/skimage.feature.html

相异性（Dissimilarity）

这里的P对应上面的G，levels对应k

计算对比度时，权重随矩阵元素与对角线的距离以指数方式增长，如果改为线性增长，则得到相异性。

同质性/逆差距（Homogeneity）

测量图像的局部均匀性，非均匀图像的值较低，均匀图像的值较高。与对比度或相异性相反，同质性的权重随着元素值与对角线的距离而减小，其减小方式是指数形式的。

能量(Energy)

是灰度共生矩阵各元素值的平方和，是对图像纹理的灰度变化稳定程度的度量，反应了图像灰度分布均匀程度和纹理粗细度。能量值大表明当前纹理是一种规则变化较为稳定的纹理。。

好了，下面开始开始实现上面的内容。通过skimage计算灰度共生矩阵并提取上述的六个纹理特征。

scikits-image库提供了两个模块，skimage.feature.greycomatrix(image, …[, …])计算灰度共生矩阵，skimage.feature.greycoprops(P[, prop])计算GLCM的纹理特征。

skimage.feature.greycomatrix(image, distances, angles, levels=256, symmetric=False, normed=False)

参数：
image：array_like of uint8

整数型输入图像。该图像将强制转换为uint8，因此最大值必须小于256。

distances：array_like

像素对距离偏移量列表。

angles：array_like

弧度的像素对角度列表。

级别：整数，可选

输入图像应在[0，levels-1]中包含整数，其中level表示所计数的灰度级数（对于8位图像，通常为256）。最大值为256。

对称：布尔型，可选

如果为True，则输出矩阵P [：，：，d，theta]是对称的。这是通过忽略值对的顺序来实现的，因此在给定偏移量遇到（i，j）时，将同时累积（i，j）和（j，i）。默认值为False。

标准化（归一化）：bool，可选

如果为True，则通过除以给定偏移量的累积共现总数，将每个矩阵P [：，：，d，theta]归一化。结果矩阵的元素总和为1。默认值为False。

举个例子演示一下
from skimage.feature import greycomatrix, greycoprops
import numpy as np
image = np.array([[0, 0, 1, 1],
[0, 0, 1, 1],
[0, 2, 2, 2],
[2, 2, 3, 3]], dtype=np.uint8)
#构建了一个44的数组
result = greycomatrix(image, [1], [0, np.pi/2], levels=4)#调用函数
#第一个参数是图像，第二个为距离1，第三个是扫描方向，这里加入了两个方向，0为向右扫描，np.pi/2,90°扫描即垂直扫描；最后一个levels=4即灰度级数。
运行这行代码，会返回一个四维的result结果。即灰度共生直方图，表示形式为p[i,j,d,theta],是灰度级j与灰度级i之间距离为d,角度为θ的情况下出现灰度级j的次数。

生成的result是一个4x4x1x2的数组，因为原始输入的image是一个44的数组，距离d这里只计算了一个，角度计算了两个，所以是4x4x1x2

读取一下result数组
result[:, :, 0, 0] #距离d为1时，扫描角度为0时的result
array([[2, 2, 1, 0],
[0, 2, 0, 0],
[0, 0, 3, 1],
[0, 0, 0, 1]], dtype=uint32)
result[:, :, 0, 1] #距离d为1时，扫描角度为90°时的result
array([[3, 0, 2, 0],
[0, 2, 2, 0],
[0, 0, 1, 2],
[0, 0, 0, 0]], dtype=uint32)
这样一看难以理解，画了个图以便大家理解

图中右边展示的是0°向右扫描，下边展示的是90°垂直扫描。生成的灰度共生矩阵里面的数是出现次数的统计。

尝试距离为2时的结果
result1 = greycomatrix(image,[2],[0, np.pi/2],levels=4)

result1[:, :, 0, 0]

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

result1[:, :, 0, 1]

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

就画了下0°扫描的（0，1），垂直扫描的示意图就不画了，图已经有点乱了。

result2 = greycomatrix(image,[1,2,3,4],[0, np.pi/2,np.pi/4,np.pi*3/4],levels=4)# 距离为1，2，3，4,扫描方向为0°水平扫描，90°垂直扫描，45°右对角线扫描，135°左对角线扫描。

#读取距离4，左对角线扫描的结果
result2[:, :, 3, 3]
array([[0, 0, 0, 0],
[0, 0, 1, 0],
[0, 0, 0, 0],
[0, 0, 0, 0]], dtype=uint32)

计算6个纹理特征

skimage.feature.greycoprops(P[, prop])

Parameters :
P : ndarray

Input array. P is the grey-level co-occurrence histogram for which to compute the specified property. The value P[i,j,d,theta] is the number of times that grey-level j occurs at a distance d and at an angle theta from grey-level i.

prop : {‘contrast’, ‘dissimilarity’, ‘homogeneity’, ‘energy’, ‘correlation’, ‘ASM’}, optional

The property of the GLCM to compute. The default is ‘contrast’.

Returns :
results : 2-D ndarray

2-dimensional array. results[d, a] is the property ‘prop’ for the d’th distance and the a’th angle.

参数：
P：ndarray

输入数组。P是用于计算指定属性的灰度共生直方图。值 P [i，j，d，θ]是灰度级j在距灰度级i的距离为d且角度为θ的角度出现的次数。

prop：{“对比度”，“相异性”，“同质性”，“能量”，“相关性”，“ ASM”}，可选

要计算的GLCM的属性。默认值为“对比度”。

返回值：
结果：2-D ndarray
二维数组。result [d，a]是第d个距离和第a个角度的属性’prop’。

image = np.array([[0, 0, 1, 1],
[0, 0, 1, 1],
[0, 2, 2, 2],
[2, 2, 3, 3]], dtype=np.uint8)
g = greycomatrix(image, [1, 2], [0, np.pi/2], levels=4,
normed=True, symmetric=True)#这里计算灰度共生矩阵，输出的矩阵P [：，：，d，theta]是对称的，并对结果进行了归一化
contrast = greycoprops(g, ‘contrast’)
contrast
array([[0.58333333, 1. ],
[1.25 , 2.75 ]])

dissimilarity = greycoprops(g, ‘dissimilarity’)
dissimilarity
array([[0.41666667, 0.66666667],
[1. , 1.5 ]])
#greycoprops()不能同时计算这六个特征，只能挨个计算。

官方文档见：http://tonysyu.github.io/scikit-image/api/skimage.feature.html

这两天因为任务刚接触这个，做了个整理，不足之处还请谅解。有错误之处还请指正，谢谢。