Python实战：基于skimage的灰度共生矩阵(GLCM)纹理特征分析与应用
1. 灰度共生矩阵(GLCM)基础入门
第一次接触灰度共生矩阵这个概念时，我也被这个拗口的名字吓到了。但实际理解后发现，它就是一种描述图像纹理特征的数学工具。想象一下，当你观察一块布料时，会注意到它的纹理走向、粗细程度、规则性等特征，GLCM就是用数学方法来量化这些视觉感受。
GLCM的核心思想是统计图像中特定距离和角度上的像素对出现的频率。举个例子，假设我们有一张简单的4x4像素图像：
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)这张图像有4个灰度级(0-3)。如果我们统计水平方向(0度)相邻像素(距离=1)的组合出现次数，就会得到一个4x4的矩阵，这就是最基本的GLCM。
在skimage中计算这个矩阵非常简单：
from skimage.feature import greycomatrix glcm = greycomatrix(image, distances=[1], angles=[0], levels=4)得到的glcm是一个4维数组，理解起来可能有点抽象。我建议新手可以这样可视化：
print(glcm[:, :, 0, 0]) # 距离=1，角度=0时的GLCM输出结果会显示每种灰度组合出现的次数。比如(0,0)位置的值表示灰度值0和0水平相邻出现的次数。
2. skimage中的GLCM实战技巧
在实际项目中，我发现skimage的greycomatrix函数有几个关键参数需要特别注意：
levels参数：这个参数新手最容易出错。它应该设置为图像中的最大灰度级数+1。对于8位图像通常是256，但如果你的图像已经归一化到0-15，就应该设为16。
距离和角度选择：根据纹理特性选择很重要。我做过一个卫星图像分析项目，发现对于细密的城市纹理，距离1-3效果最好；而对于广阔的农田，可能需要更大的距离值。
一个实用的多角度计算示例：
import numpy as np from skimage import io, color from skimage.feature import greycomatrix # 加载图像并转为灰度 image = io.imread('texture.jpg') gray = color.rgb2gray(image) gray = (gray * 255).astype(np.uint8) # 转为0-255的uint8 # 计算多角度GLCM distances = [1, 2, 3] angles = [0, np.pi/4, np.pi/2, 3*np.pi/4] glcm = greycomatrix(gray, distances=distances, angles=angles, levels=256)这里有个小技巧：如果图像太大，可以先裁剪感兴趣区域(ROI)进行计算，能显著提升效率。
3. 六大纹理特征详解与提取
GLCM本身只是统计矩阵，真正有用的是从中提取的纹理特征。skimage提供了六个最常用的特征：
对比度(Contrast)：衡量图像的清晰度。在医学影像分析中，我发现肿瘤区域的对比度通常比正常组织高。
相异性(Dissimilarity)：类似于对比度，但对微小变化更敏感。在工业质检中特别有用。
同质性(Homogeneity)：反映局部均匀性。在卫星图像分类时，水体通常有很高的同质性。
提取这些特征的代码很简单：
from skimage.feature import greycoprops # 计算所有特征 features = {} props = ['contrast', 'dissimilarity', 'homogeneity', 'energy', 'correlation', 'ASM'] for prop in props: features[prop] = greycoprops(glcm, prop)但要注意，greycoprops一次只能计算一个特征。在我的项目中，通常会封装一个函数来批量计算：
def extract_glcm_features(image, distances, angles, levels=256): glcm = greycomatrix(image, distances, angles, levels) features = {} for prop in ['contrast', 'dissimilarity', 'homogeneity', 'energy', 'correlation', 'ASM']: features[prop] = greycoprops(glcm, prop) return features4. 实际应用案例与参数调优
在遥感图像分类项目中，我发现参数选择直接影响分类效果。经过多次实验，总结出以下经验：
- 医学影像：距离1-3，四个角度(0,45,90,135)，levels=256
- 工业检测：距离1-5，重点关注0度和90度
- 卫星图像：距离3-7，所有角度
一个完整的遥感图像分类示例：
import numpy as np from skimage import io, color, exposure from sklearn.ensemble import RandomForestClassifier # 1. 数据准备 image = io.imread('satellite.tif') gray = color.rgb2gray(image) gray = exposure.rescale_intensity(gray, out_range=(0, 255)).astype(np.uint8) # 2. 提取纹理特征 def extract_features(roi): glcm = greycomatrix(roi, distances=[3], angles=[0, np.pi/4, np.pi/2, 3*np.pi/4], levels=256) features = [] for prop in ['contrast', 'dissimilarity', 'homogeneity', 'energy']: features.append(greycoprops(glcm, prop).ravel()) return np.concatenate(features) # 3. 滑动窗口提取特征 window_size = 30 features = [] for i in range(0, gray.shape[0]-window_size, window_size//2): for j in range(0, gray.shape[1]-window_size, window_size//2): roi = gray[i:i+window_size, j:j+window_size] features.append(extract_features(roi)) # 4. 分类器训练(假设已有标签) clf = RandomForestClassifier() clf.fit(features, labels)这个流程在多个项目中都取得了不错的效果，关键是要根据具体问题调整窗口大小和特征组合。
5. 常见问题与性能优化
在实际使用中，我踩过不少坑，这里分享几个常见问题的解决方法：
内存不足问题：计算大图像的GLCM时，经常会遇到内存错误。解决方法：
- 降低灰度级数(levels)，比如从256降到64
- 分块处理图像
- 使用sparse矩阵
特征选择问题：不是所有六个特征都有用。我通常先用随机森林评估特征重要性，然后只保留重要的几个。
一个实用的特征重要性评估代码：
from sklearn.ensemble import RandomForestClassifier import matplotlib.pyplot as plt # 假设X是特征矩阵，y是标签 clf = RandomForestClassifier() clf.fit(X, y) # 绘制特征重要性 importances = clf.feature_importances_ features = ['Contrast_0', 'Contrast_45', ..., 'ASM_135'] # 根据实际特征命名 plt.barh(range(len(importances)), importances) plt.yticks(range(len(importances)), features) plt.show()参数敏感性问题：GLCM对参数很敏感。建议先用网格搜索确定最佳参数组合：
from sklearn.model_selection import GridSearchCV param_grid = { 'distances': [[1], [1,2], [1,3]], 'angles': [[0], [0, np.pi/2], [0, np.pi/4, np.pi/2, 3*np.pi/4]] } best_score = 0 best_params = {} for distances in param_grid['distances']: for angles in param_grid['angles']: features = extract_glcm_features(train_images, distances, angles) score = cross_val_score(clf, features, labels, cv=5).mean() if score > best_score: best_score = score best_params = {'distances': distances, 'angles': angles}6. 进阶技巧与可视化方法
为了让纹理分析结果更直观，我经常使用以下可视化技巧：
热力图可视化：将纹理特征值映射回图像位置
def visualize_feature_map(image, feature_name, window_size=30): gray = color.rgb2gray(image) gray = (gray * 255).astype(np.uint8) feature_map = np.zeros_like(gray, dtype=np.float32) for i in range(0, gray.shape[0]-window_size, window_size//2): for j in range(0, gray.shape[1]-window_size, window_size//2): roi = gray[i:i+window_size, j:j+window_size] glcm = greycomatrix(roi, distances=[3], angles=[0], levels=256) feature_val = greycoprops(glcm, feature_name)[0,0] feature_map[i:i+window_size, j:j+window_size] = feature_val plt.imshow(feature_map, cmap='hot') plt.colorbar() plt.title(f'{feature_name} Feature Map') plt.show()多特征联合分析：有时候单个特征区分力不足，我会组合多个特征：
# 计算能量和同质性的比值 glcm = greycomatrix(image, distances=[1], angles=[0], levels=256) energy = greycoprops(glcm, 'energy')[0,0] homogeneity = greycoprops(glcm, 'homogeneity')[0,0] ratio = energy / homogeneity这个比值在区分某些相似纹理时特别有效，比如在木材缺陷检测项目中，它能更好地区分裂纹和木节。
7. 与其他纹理分析方法的对比
虽然GLCM很强大，但它不是唯一的纹理分析方法。在我的项目中，经常会结合其他方法：
- LBP(局部二值模式)：计算更快，但对噪声更敏感
- Gabor滤波器：能捕捉多尺度纹理特征，但参数更多
- 深度学习：效果最好，但需要大量标注数据
一个实用的混合特征提取方案：
from skimage.feature import local_binary_pattern def extract_hybrid_features(image): # GLCM特征 glcm = greycomatrix(image, distances=[1,2], angles=[0], levels=256) glcm_features = [greycoprops(glcm, prop)[0,0] for prop in ['contrast', 'homogeneity']] # LBP特征 lbp = local_binary_pattern(image, P=8, R=1) lbp_hist, _ = np.histogram(lbp, bins=256) lbp_hist = lbp_hist / lbp_hist.sum() # 归一化 return np.concatenate([glcm_features, lbp_hist])这种混合方法在多个项目中都取得了比单一方法更好的效果，特别是在数据量不大的情况下。