空间域二阶统计纹理

纹理是影像的局部特征，其度量方法有：基于空间一阶和二阶统计特征、基于傅里叶功率谱、基于分形、基于小波变换等。这里主要介绍空间域的二阶统计纹理，资料来源于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.

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