空间域二阶统计纹理

纹理是影像的局部特征，其度量方法有：基于空间一阶和二阶统计特征、基于傅里叶功率谱、基于分形、基于小波变换等。这里主要介绍空间域的二阶统计纹理，资料来源于Mryka Hall-Beyer的主页，包括以下内容：

• 定义和创建灰度共生矩阵（GLCM）
• 基于GLCM的纹理运算
• 其他相关问题

1. GLCM的定义及创建

1.1 定义

灰度共生矩阵（Grey Level Co-occurrence Matrix, GLCM）是像元灰度值的不同组合（combination）出现的频率的列联表，也叫Grey Tone Spatial Dependency Matrix(灰度色调空间依赖矩阵)。

测试影像及其像元值如下：

1.2 普通GLCM的结构

1. 计算窗口
   一般为\((2k-1)\times(2k-1)\)的窗口，纹理特征是该窗口内影像的局部特征。本文的测试影像是$ 4 \times 4 $的一小幅影像，直接使用它计算纹理而不采用任何窗口。

2. 位移(Offset)
   位移是一个向量\(\vec{d}=(\Delta x, \Delta y)\)，不同的位移决定了不同方向的纹理特征，理论上可以计算任意方向的纹理特征。

3. 两个像元间的空间关系
   GLCM纹理考虑了两个像元间的关系，一个为参考像元（reference pixel），另一个为邻域像元（neighbour pixel），两个像元的灰度值按照这种空间关系会形成一个组合。
   计算某个组合的频数时，从窗口左上角到右下角，每个像元依次作为参考像元。若参考像元的位置为\((i,j)\)，则处于\((i+Δx,j+Δy)\)位置的像元就是其邻域像元。比如\(\vec{d}=(1,0)\)时，邻域像元为参考像元右侧的第一个像元，此时纹理的方向是向东。

4. 普通GLCM
   假设测试影像灰度值的范围为\([0,3]\)，\(\vec{d}=(1,0)\)，则所有可能的组合如下表：

	0	1	2	3
0	(0,0)	(0,1)	(0,2)	(0,3)
1	(1,0)	(1,1)	(1,2)	(1,3)
2	(2,0)	(2,1)	(2,2)	(2,3)
3	(3,0)	(3,1)	(3,2)	(3,3)

将每个组合出现的频数作为对应位置的矩阵元素，得到的矩阵就是测试影像的一个GLCM：

	0	1	2	3
0	2	2	1	0
1	0	2	0	0
2	0	0	3	1
3	0	0	0	1

可以记为：

\[ G=\begin{bmatrix} 2&2&1&0\\ 0&2&0&0\\ 0&0&3&1\\ 0&0&0&1 \end{bmatrix}\]

1.3 对称GLCM

以上GLCM并不对称，而纹理运算要求GLCM是对称矩阵，因为非对称矩阵按照行和列计算得到的纹理值不同。

GLCM对称意味着两个灰度值顺序相反的组合有相同的频数，比如(2,3)和(3,2)。也就是说将它们视为同一组合，计算过程中每个像元组合都按两次计数则能获得对称矩阵：一次正向、一次反向。以\(\vec{d}=(1,0)\)为例，正向是向东，反向则向西，对应\(\vec{d}=(0,1)\)，此时其实得到的是水平方向的GLCM。

这样以来，计算纹理的方向就会在0°-180°范围内。然而由于窗口大小有限，只需要选择有代表性的几个角度即可例如下图所示的四个方向：

获得对称阵的一个简单方法是加上自身的转置矩阵。测试影像对应的对称GLCM为：

\[ V=G+G^T= \begin{bmatrix} 2&2&1&0\\ 0&2&0&0\\ 0&0&3&1\\ 0&0&0&1 \end{bmatrix}+\begin{bmatrix} 2&0&0&0\\ 2&2&0&0\\ 1&0&3&0\\ 0&0&1&1 \end{bmatrix}=\begin{bmatrix} 4 & 2 & 1 & 0\\ 2 & 4 & 0 & 0\\ 1 & 0 & 6 & 1\\ 0 & 0 & 1 & 2 \end{bmatrix} \]

1.4 归一化GLCM

\(V\)的元素是不同组合出现的次数，也就是频数，\(V\)是一个列联表。为了得到灰度值组合出现的频率（或概率），可以进一步进行归一化，得到归一化GLCM。归一化的过程就是每个组合的频数除以总频数，若归一化GLCM记为P，则：

\[P_{i,j} = \frac{V_{i,j}} {\sum\limits_{i,j=0}^{N-1} V_{i,j}} \]

其中\(V_{i,j}\)是矩阵元素，\(N\)是行数或列数（两者相等）。符号\(\sum\limits_{i,j=0}^{N-1}\)是\(\sum\limits_{i=0}^{N-1}\sum\limits_{j=0}^{N-1}\)的简写。

归一化得到的是联合概率分布：

\[P = \begin{bmatrix} .166 & .083 & .042 & 0\\ .083 & .166 & 0 & 0\\ .042 & 0 & .25 & .042\\ 0 & 0 & .042 & .083 \end{bmatrix} \]

如无特别说明，下文中GLCM均指归一化的GLCM。

1.5 GLCM的属性

1. GLCM是方阵

2. 行数和列数与影像的量化等级相同
   对于量化为\(k\) bit的影像而言，其量化等级为\(2^k\)，灰度值的取值范围是为\([0,2^k-1]\)，对应的GLCM大小为\(2^k \times 2^k\)。

3. GLCM对称

   • 主对角线上是值相等的组合（0-0, 1-1, 2-2...），主对角线元素和越大，说明影像值越均一，差异越小，纹理越不明显
   • 平行于主对角线的直线（副对角线）上，是与参考像元有差值的组合，越远插值越大，最近的两条线上差值为1，次近的两条线差值为2，依此类推

---

2. 基于GLCM的纹理运算

大多数纹理运算（texture calculation）是归一化GLCM元素值的加权平均。加权平均的主要目的是强调归一化GLCM中不同值的相对重要性。

按照不同纹理特征所表达的纹理信息，将其分为三组：

• 反差（Contrast）
• 有序性（Orderliness）
• 描述性统计量（Descriptive Statistics）

2.1 反差

与反差相关的特征使用的是和GLCM对角线距离相关的权重，也就是说，权重是距离的函数。

为了强调窗口内的局部反差，权重应该使反差大的影像的运算结果得到较大的数值。对角线元素没有反差，离对角线越远反差越大。因此，权重应该随距离而增加。

• 对比度（Contrast）

\[Con = \sum\limits_{i,j=0}^{N-1} P_{i,j}(i-j)^2 \]

以前文为例，对应的权重矩阵为：

\[ \begin{bmatrix} 0&1&4&9\\ 1&0&1&4\\ 4&1&0&1\\ 9&4&1&0 \end{bmatrix}\]

因此

\[Con = 0\cdot.166 + 1\cdot.083+\cdots+0\cdot.083=0.586 \]

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

\[Dis = \sum\limits_{i,j=0}^{N-1} P_{i,j}|i-j| \]

P对应的权重矩阵为：

\[ \begin{bmatrix} 0&1&2&3\\ 1&0&1&2\\ 2&1&0&1\\ 3&2&1&0 \end{bmatrix}\]

且$ Con = 0.418 $

• 同质性（Homogeneity）
  与对比度或相异性相反，同质性的权重随着元素值与对角线的距离而减小，其减小方式是指数形式的。

\[Hom = \sum\limits_{i,j=0}^{N-1} \frac{P_{i,j}}{1+(i-j)^2} \]

P对应的权重矩阵为：

\[ \begin{bmatrix} 0&0.5&0.2&0.1\\ 0.5&1&0.5&0.2\\ 0.2&0.5&1&0.5\\ 0.1&0.2&0.5&0 \end{bmatrix}\]

且\(Hom=0.807\)

2.2 有序性

有序性（Orderliness）是指窗口内像元值的规律程度或有序程度。比如以下两个窗口A, B，其水平方向有相同的对比度（都为1），因为每个像元都比左侧像元大1。但是A, B的有序性却很不相同，A很明显比B更有规律。

\[ A=\begin{bmatrix} 1&2&3&4\\ 1&2&3&4\\ 1&2&3&4\\ 1&2&3&4 \end{bmatrix} \; \; \; \; \; B=\begin{bmatrix} 3&4&5&6\\ 1&2&3&4\\ 2&3&4&5\\ 4&5&6&7 \end{bmatrix} \]

在GLCM中，元素的值是各种组合的概率，用GLCM元素本身来构造权重就可以得到度量有序性的纹理特征。

• 角度二阶矩（Angular Second Moment, ASM）和能量（Energy）
  这两个特征直接以\(P_{i,j}\)作为权重，W=P，有序性越高，ASM和能量越大。能量有时也被称为均匀性（Uniformity）。

\[ ASM = \sum\limits_{i,j=0}^{N-1} P_{i,j}^2 \; \; \; \; \; Energy=\sqrt{ASM}\]

对于测试影像，\(ASM=0.145,Energy=0.381\)。

• 熵（Entropy）
  以\(P_{i,j}\)的自然对数的相反数作为权重。有序性越低，熵越大。

\[Ent= \sum\limits_{i,j=0}^{N-1} P_{i,j}(-lnP_{i,j}) \]

当元素为0时，由于ln0没有定义，这时的熵被定义为0。上文示例的熵为2.095。

度量有序性的三个特征都有一定的物理意义。角度二阶矩是旋转加速度(额，其实，咱家也不是很理解...)的一种度量，能量则是一个常见的物理量，而熵在热动力学中代表的是永久损失掉无法再利用的热量。熵可以不严格的理解为混沌或无序，可以认为和能量是相反的一对概念。

2.3 描述性统计量

描述性统计量包括GLCM均值、GLCM方差（或GLCM标准差）和GLCM相关系数，这三个特征的计算过程和一般的均值、方差、相关系数的计算方式并没有本质区别，只是这里使用的是GLCM矩阵的元素值而不是原始影像的灰度值，也就是说，这三个特征度量的是灰度值组合出现的期望，标准差和相关系数。

如果使用原始影像的灰度值，得到的是一阶统计纹理。

• GLCM均值（GLCM Mean）

\[\mu_i=\sum\limits_{i,j=0}^{N-1} iP_{i,j} \; \; \; \; \; \mu_j=\sum\limits_{i,j=0}^{N-1} jP_{i,j} \]

**左侧的权重为参考像元的灰度值**，得到的是基于参考像元的GLCM均值，而**右侧的权重为邻域像元的灰度值**，得到的是基于邻域像元的GLCM均值。由于GLCM的对称性，这两个结果是相等的。

GLCM均值是像元的灰度值（i或j）与该像元的所有组合出现的频率的加权平均，其实质是离散型随机变量的期望，因此得到的是GLCM均值。

*文中，求和从0到N-1而不是1到N，这对于GLCM均值、方差和相关系数的计算很方便，行列号就是参考像元和邻域像元的灰度值。*

• GLCM方差（GLCM Variance）

\[\sigma_i^2=\sum\limits_{i,j=0}^{N-1} P_{i,j}(i-\mu_i)^2 \; \; \; \; \; \sigma_j^2=\sum\limits_{i,j=0}^{N-1} P_{i,j}(j-\mu_j)^2 \]

GLCM标准差（GLCM Standard Deviation）为：

\[\sigma_i=\sqrt{\sigma_i^2} \; \; \; \; \; \sigma_j=\sqrt{\sigma_j^2} \]

与一般方差或标准差一样，GLCM方差或标准差度量的是GLCM元素值以其均值为中心的分散程度，类似于熵。

• GLCM相关系数（GLCM Correlation）

\[r=\sum\limits_{i,j=0}^{N-1} P_{i,j} \left[ \frac{(i-\mu_i)(j-\mu_j)}{\sqrt{\sigma_i^2\cdot\sigma_j^2}} \right] \]

GLCM相关系数度量的是像元灰度值与其邻域像元灰度值的线性相关性。

GLCM相关性的计算和其他纹理特征的计算有很大的不同，因此，它**独立与其他纹理特征**，提供不同的信息。将它和其他特征组合使用往往比较好。它的值有更直观的意义：0表示无关，1表示完全相关。

*当窗口内的灰度值相等时，GLCM标准差为0，相关系数分母为0。这时，不同的软件有不同的处理方式，定义为1（完全相关）更符合逻辑。*

---

3. 其他相关问题

3.1 数据压缩（重采样）

• 减少计算量
  前文说过，GLCM的大小取决于影像的量化等级，量化等级越高矩阵越大。由于GLCM创建时计算复杂度很高，一般将影像压缩为8bit。

• 改善统计结果
  另外，如果直接使用原始灰度值，可能得到很稀疏的GLCM，因为很多灰度值组合不会出现在影像上。压缩后的GLCM稀疏性更低，有利于改善统计有效性（statistical validity）。

3.2 如何建立纹理影像

纹理运算得到窗口内某一纹理特征的值，最终将该值赋予窗口中心位置的像元，得到一景纹理影像。

影像边缘
对于\(N\times N\)的窗口，在影像边缘上将有宽度为(N-1)/2的条带不能处于窗口中心，通常的处理方式是用最近的纹理值来填充。

3.3 GLDV

GLDV(Grey Level Difference Vector，灰度差值向量)是GLCM每条对角线上元素的和，分别是差值为0, 1, 2...的像元组合的频数或频率。测试影像的GLDV如下表：

差值	频数	频率
0	16	.666
1	6	.250
2	2	.083
3	0	0

3.4 多通道影像处理

• 每个通道单独计算纹理特征
• 使用主成份变换或其他方法进行组合
• 任意三个纹理特征可以投影到色彩空间查看
• 可以获取多个纹理影像的主成份，但结果的解释十分困难

3.5 参数选择

纹理运算涉及多个参数的选择：

• 窗口大小
• 位移方向
• 位移大小
• 影像通道
• 纹理特征

下面是一些经验法则：

• 窗口大小应该小于物体大小，同时应该足够大以包含对象的特征变化
• 通过目视检查或用主成份来减少通道数量，有些对象在不同的通道上表现不同
• 目视检查确定重要的方向，如果没有明显的方向性，那么选择所有方向求平均
• 位移几乎总是为1
• 纹理特征间存在相关性，其中只有4到5个是完全独立的（见下文）
• 最后，剩余的参数组合可以使用特征选择方法，或者直接比较其精度

3.6 纹理特征间的相关性

由于纹理特征的计算方法只是权重不同，因此有些特征有很大的相关性。比如对比度和相异性，本质上包含了相同的信息。而GLCM方差和度量反差的特征也有很高的相关性。

下面是特征间的相关系数，假设具有相同的窗口，不同的影像可能不尽相同，但大体趋势是可以肯定的：

• 同质性 vs 对比度，\(r=-0.80\)
• 同质性 vs 相异性，\(r=-0.95\)
• GLCM方差 vs 对比度，\(r=0.89\)
• GLCM方差 vs 相异性，\(r=0.91\)
• GLCM方差 vs 同质性，\(r=0.83\)
• 熵 vs ASM，\(r=-0.87\)

而GLCM均值和GLCM相关系数的无关性要高一点，在测试影像中与其他特征的相关系数分别小于0.1和0.5.

因此，实际可以选择一个反差特征，一个有序性特征，两个或三个描述性特征。