尧图网页设计 
简介一份基于Python实现的遥感影像滑坡场景分类代码与项目文档面向毕业设计、课程设计及项目开发场景适合想快速掌握SVM分类流程并扩展应用的学习者。源码已经过严格测试可放心参考并在此基础上继续延伸。资源压缩包共21个文件、约1.86MB包括4个py源码、训练与测试txt数据、pkl特征/模型文件、lda模型相关文件及说明md目录按功能模块划分便于对照学习。整体流程清晰先从影像提取光谱统计特征与GLCM纹理特征再用K-Means聚类生成视觉单词接着通过LDA模型将词袋抽象为主题分布并保存为libsvm格式最后采用Libsvm完成SVM分类。项目文档对各环节原理与运行方式进行了说明便于理解特征工程—词袋模型—主题模型—分类器这条完整链路。目前已有56人学习下载用作课程报告或论文实验的参考基线都很合适。1. 遥感影像滑坡场景分类不是语义分割而是特征SVM的经典管线第一次拿遥感影像做滑坡场景分类时我的第一反应是上深度学习语义分割但实际一跑才发现场景分类要的是「这张图块是不是滑坡」不需要逐像元分割用语义分割反而把问题做复杂了。这个基于Python实现的项目走的是另一条完整可复现的路线从影像里提取光谱特征均值、标准差和GLCM纹理特征用K-Means把低层特征聚类成视觉单词再用LDA主题模型对视觉词袋做抽象把每张影像变成主题分布最后交给LibSVM分类器完成滑坡场景识别。整套流程纯CPU就能跑起来适合毕业设计、课程设计做闭环演示也适合想在遥感场景分类方向上做快速验证的从业者。资源里源码、中间模型文件和项目文档都齐顺着代码能一路推进到精度调优。2. 特征提取光谱统计量与GLCM纹理参数怎么定2.1 光谱特征为什么取均值和标准差特征提取是整条管线的底座也是最容易被跳过思考的一步。项目对每张影像块提取光谱特征取的是三个波段的平均值和标准差组合起来就是一个六维向量。为什么是这两个统计量而不是直方图或者更复杂的颜色描述子因为场景分类的输入是固定大小的图块均值能反映图块的整体亮度分布标准差能反映波段内的波动程度。滑坡体的裸土区域在可见光波段上通常亮度偏高、标准差偏大而同面积的草坡、林地方向均值低、波动小这两个统计量已经能把正负样本拉开。直方图虽然信息更充足但维度会膨胀到几十甚至上百对后续K-Means和SVM都不友好。从像素级到向量级的转换常见做法是先读入影像块再按通道统计。下面的代码是整个特征提取流程的第一段import numpy as np from PIL import Image def spectral_features(patch): 输入三波段影像块返回光谱均值与标准差特征6维 arr np.asarray(patch).astype(np.float32) means arr.mean(axis(0, 1)) # 每个波段一个均值 stds arr.std(axis(0, 1)) # 每个波段一个标准差 return np.concatenate([means, stds])这里arr的shape是(H, W, 3)axis(0,1)表示把所有像素在高度和宽度方向上做统计保留波段维。返回的向量顺序是三个均值接三个标准差调用方必须固定这个顺序否则后面生成libsvm格式时特征列会整体错位这是个最容易犯的低级错误。如果输入是单波段灰度图这段代码同样能跑只是输出的均值、标准差各一个维度变成2。光谱特征这一段虽然简单但在整个项目里承担的任务是「颜色分布的粗颗粒度描述」。后续的GLCM纹理特征负责描述空间关系两者互补之后一个图块的特征向量才具备区分滑坡和非滑坡的能力。2.2 GLCM纹理特征窗口、灰度级、方向与距离纹理特征这块项目选的是GLCM灰度共生矩阵的多个统计值。GLCM的核心思想很直接统计某个灰度值在指定方向和距离上与另一个灰度值相邻出现的频率。如果图块纹理粗糙像素对规律性强共生矩阵会呈现明显的主对角分布如果纹理细碎矩阵分布就均匀扩散。参数是GLCM落地时最需要较真的地方。我一般会把灰度级压缩到16或32距离取1方向取0度、45度、90度、135度四个方向然后对四个方向的统计量取均值。灰度级如果取256共生矩阵会变成256×256绝大多数格子都是零统计结果全是稀疏值方差大且不稳定距离取1能刻画相邻像素的最短关系取大了反而把纹理结构打成碎块。窗口大小通常用7×7或9×9窗口太小对边缘敏感窗口太大容易把两类地物的边界抹平。四个方向取均值是为了削弱影像旋转对纹理统计的影响这也是场景分类任务里很实用的做法。import numpy as np from skimage.feature import graycomatrix, graycoprops def glcm_features(gray_patch, levels16, distance1): 灰度图块输出GLCM四个方向统计值的均值特征 img (gray_patch / 255.0 * (levels - 1)).astype(np.uint8) angles [0, np.pi / 4, np.pi / 2, 3 * np.pi / 4] glcm graycomatrix(img, [distance], angles, levelslevels, symmetricTrue, normedTrue) stats [] for name in [contrast, dissimilarity, homogeneity, energy, correlation, asm]: prop graycoprops(glcm, name) # 输出 shape: (1,4) stats.extend(prop.mean(axis1)) # 四个方向取平均 return np.array(stats, dtypenp.float32)代码里两个关键点symmetricTrue让共生矩阵对称化这样0度和90度这类互为反向的组合会被合并矩阵更稳定normedTrue把频率归一化成概率避免图块尺寸不同导致数值量级不一致。graycoprops输出的每个统计量在四个方向上各有一个值对axis1取均值后返回的向量长度就是统计量个数乘4方向归一后的数量。这里有个容易踩的坑skimage的graycomatrix要求输入是uint8类型而且levels参数不能小于输入数组的最大值加1。如果之前做了归一化到[0,1]的浮点操作直接丢进去会报错或者得到全零矩阵。所以我在代码里先乘回[0,255]区间再转uint8这样levels16时数据能正确落在0到15的灰度级上。特征向量最终是6维光谱加6个GLCM统计量总计12维进入下一层K-Means时的特征空间并不算大这也是经典特征工程的典型做法。2.3 样本与标签生成滑窗取样和训练测试划分特征函数准备好之后下一步是把原始影像变成一批带标签的样本。项目里的generate_samples_and_labels.py干的就是这件事。场景分类一般用滑窗从大影像上切图块每个图块对应一个标签。滑坡场景是典型的小样本正类所以采样策略要特别设计正样本图块尽量全部保留负样本按正样本数量的一定比例采样否则SVM很容易被大量负样本带偏。图块尺寸常见做法是64×64或者128×128步长取等于窗口尺寸时样本之间没有重叠取小于窗口尺寸时有重叠相当于做了滑动增强。我建议先把图块坐标和标签存成中间文件再做训练测试划分这样后面调参时不用重新切图。import cv2 import numpy as np def extract_patches_with_label(image, patch_size64, stride64): 按滑窗切图块返回坐标列表标签由外部掩膜对应坐标决定 h, w image.shape[:2] patches [] for y in range(0, h - patch_size 1, stride): for x in range(0, w - patch_size 1, stride): patches.append((y, x, image[y:y patch_size, x:x patch_size])) return patches def make_train_test(patches, labels, test_ratio0.2): 按顺序切分训练测试注意先打乱索引 idx np.arange(len(patches)) rng np.random.RandomState(42) rng.shuffle(idx) split int(len(idx) * (1 - test_ratio)) train_idx, test_idx idx[:split], idx[split:] return train_idx, test_idx这里在切分前用固定RandomState打乱索引是保证实验结果可复现的关键。实际业务中更严格的划分是按影像文件划分而不是按图块划分同一张大图切出的图块共享光照和地物背景如果不按影像划分训练测试之间会有严重的信息泄漏精度虚高。这个细节在第五章会用专门一节展开讲。样本坐标和图块内容生成后要对每一块做特征提取。常见做法是对每个图块分别调用光谱特征和GLCM特征然后横向拼接成一个特征行所有图块的特征行纵向堆叠就形成了后续K-Means聚类的输入矩阵。这一步如果直接用Python循环跑速度会比较慢但项目样本量在几千的量级时完全够用。数据准备做到这里管线才走完第一层接下来进入视觉词袋阶段。3. 特征抽象K-Means视觉词袋与LDA主题分布为什么中间要多一层3.1 为什么低层特征不能直接喂SVM很多人拿到12维特征后会想既然维度不高为什么不直接训练SVM分类器还要绕一圈做K-Means和LDA答案是「一个图块的特征」和「一张影像的语义」之间还差着一个层级。每个64×64图块的特征描述的是局部颜色和纹理而一张遥感影像往往有一万个图块把这些图块特征直接平均成一条向量会丢失掉「哪些纹理模式在空间上频繁共现」这个重要信号。视觉词袋Bag of Visual Words的思路就是把影像类比成文档每个图块的特征向量经过K-Means聚类得到一个簇编号作为视觉单词一张影像内所有视觉单词的出现频次组成词袋直方图。这样影像被描述成「哪些视觉模式出现了多少次」而不是一堆没名字的连续特征。LDA主题模型再在这个词袋上做一层抽象输出多个主题的分布概率得到的就是项目描述里说的「高级特征」。这一层抽象的价值在分类精度上体现得很明显。直接用12维物理特征训练SVM碰上不同光照、不同季节的影像特征分布偏移很厉害而经过视觉词袋和主题分布两级抽象后分类器面对的是「语义模式比例」而非「像素统计量」泛化能力明显更强。3.2 K-Means聚类参数簇数量、初始化与标准化K-Means这一步把大量图块特征聚成若干个簇每个簇中心就是一个视觉单词。簇数量是关键参数常见做法是取200到500之间数量太小时不同地物被强行合并区分能力弱数量太大时每个簇里样本太少词袋直方图变得特别稀疏反而放大噪声。我一般先用256跑通流程再看后续分类精度的变化趋势决定是否需要调大。K-Means对特征尺度非常敏感。光谱特征取值在0到255之间GLCM统计值基本在0到1附近如果直接送进K-Means距离计算会被光谱维度主导等于纹理特征白提了。所以训练聚类器之前需要先用StandardScaler对特征矩阵做标准化。from sklearn.preprocessing import StandardScaler from sklearn.cluster import KMeans import numpy as np scaler StandardScaler() X_scaled scaler.fit_transform(feature_matrix) kmeans KMeans(n_clusters256, initk-means, n_init10, max_iter300, random_state42) kmeans.fit(X_scaled) word_ids kmeans.predict(X_scaled) # 每个图块一个视觉单词ID参数里n_init10表示用10个不同初始中心跑K-Means取最优结果这一步能有效避免聚类结果陷入局部最优random_state42固定了初始状态保证重复运行得到相同的单词划分。后续生成词袋时训练集和测试集必须用同一个kmeans对象做predict绝不能对测试集重新拟合聚类器否则视觉单词的编号语义就对不上了。词袋的构建在图像检索里也叫直方图量化。每张影像里统计每个簇编号出现的次数得到一个长度为n_clusters的整数向量from collections import Counter def build_bow(word_ids, n_clusters256): counter Counter(word_ids) bow np.zeros(n_clusters, dtypenp.float32) for wid, cnt in counter.items(): bow[wid] cnt return bow这个向量绝大部分位置是零属于高维稀疏表示。LDA模型接受的输入正是这种离散计数型矩阵所以词袋结果可以直接喂给下一层不需要额外做归一化。3.3 LDA主题模型从词袋到主题分布再写入libsvm格式LDALatent Dirichlet Allocation在这里是主题模型不是线性判别分析。这是一个特别容易混淆的坑sklearn里LinearDiscriminantAnalysis和LatentDirichletAllocation是两个完全不同的类前者是监督降维后者是无监督主题建模。项目文档里写的「利用LDA模型对视觉词袋进行主题分析」指的一定是后者。LDA把每张影像看成一篇文档视觉单词看成词主题数就是影像中可能存在的「隐含场景模式」。主题数太少会把不同地表类型揉在一起太多会产生大量重复主题。我做这个项目的经验是先从50个主题开始观察训练集和验证集精度的交叉验证曲线再决定加还是减。from sklearn.decomposition import LatentDirichletAllocation lda LatentDirichletAllocation( n_components50, learning_methodonline, max_iter30, random_state0 ) lda.fit(bow_matrix) topic_features lda.transform(bow_matrix) # 每张影像一个50维主题分布learning_methodonline适合样本量较大的情况训练速度快且迭代更稳定max_iter控制迭代轮数轮数太少主题分布还没收敛轮数太多耗时线性增长。topic_features每一行元素之和接近1含义是该影像在50个主题上的分布概率这个连续向量就是SVM要吃的最终特征。由于最后要训练LibSVM分类器需要把主题分布和类别标签写成libsvm文本格式。libsvm格式的约定是每行第一列为标签后面是“特征编号:特征值”编号从1开始连续递增def write_libsvm(topic_features, labels, output_path): with open(output_path, w) as f: for lab, vec in zip(labels, topic_features): items [{}:{:.6f}.format(i 1, v) for i, v in enumerate(vec)] line {} {}\n.format(int(lab), .join(items)) f.write(line)这段代码里标签用int强制转换防止numpy的float标签写入文件后变成“1.0”让LibSVM解析报错特征值保留6位小数就足够写太长文件会膨胀对精度没有实际帮助。项目里training.txt、test.txt以及training2.txt、test2.txt就是不同特征组合或样本划分下生成的libsvm文件classify.py读取这些文件完成后续训练。资源里的visual_lda.model、text_clt.pkl、spec_clt.pkl这些文件对应的正是这一层的中间产物。visual_lda.model是训练好的LDA模型各个clt.pkl是不同特征分支训练出的K-Means聚类器。把它们落盘的意义在于测试新影像时可以直接加载同一套模型保证单词编号和主题空间完全对齐。4. 分类器落地LibSVM格式、C/gamma搜索与模型持久化4.1 LibSVM文本格式与数据读取整个管线到了最后一步才进入分类环节但LibSVM数据格式的读取和转换反而是很多新手卡壳的地方。LibSVM工具包是SVM训练最常用的开源实现项目里生成的libsvm样本格式每行的结构是「标签 索引:值 索引:值...」索引从1开始而不是Python习惯的0。这个格式设计是为了支持稀疏存储对LDA输出的50维主题分布来说这个向量是稠密的不需要稀疏存储但仍然要按索引编号写入。读取libsvm格式的方式有两种。如果直接用LibSVM自带的命令行工具txt文件直接喂给svm-train就能训练如果整个管线都在Python里跑更常见的做法是用sklearn的svm.SVC它底层同样调用了LibSVM实现训练效果一致还能和前面的K-Means、LDA无缝衔接。分类器要读取libsvm文本先把文本转成稠密矩阵import numpy as np def load_libsvm(path, dim50): X_list, y_list [], [] with open(path, r) as f: for line in f: parts line.strip().split() y_list.append(int(float(parts[0]))) vec np.zeros(dim, dtypenp.float32) for item in parts[1:]: idx, val item.split(:) vec[int(idx) - 1] float(val) X_list.append(vec) return np.array(X_list), np.array(y_list)解析时每行的第一个元素要先用float再转int兼容标签写成“1”和“1.0”两种情况索引列要减1因为文件里是从1开始编号而numpy数组下标从0开始。如果dim设置得比文件里最大的索引还要小运行时会报索引越界反之则会在向量末尾留出一堆零这两种错误都不容易一眼发现。4.2 核函数选择与C/gamma参数搜索SVM分类器的核函数这个项目场景下直接选RBF径向基核就行。RBF核能把低维空间里线性不可分的样本映射到高维空间滑坡场景的主题分布特征虽然已经是高层抽象但正负样本的边界仍然是非线性的线性核在这个任务上精度会差一截。RBF核有两个核心参数C是误分类惩罚系数gamma控制RBF核的宽度。C太大容易过拟合C太小又欠拟合gamma太大会导致每个训练样本只影响自己附近很小的区域gamma太小则决策边界过于平滑。C和gamma怎么定我一般直接网格搜索。范围先按数量级大步长扫C取[1, 10, 100]gamma取[0.001, 0.01, 0.1]找到最优值附近再加密。如果正负样本不平衡网格搜索的scoring不能只盯着accuracy看滑坡场景里负样本占比高模型全部预测成负样本也可能有70%以上的准确率必须把scoring换成f1。from sklearn.svm import SVC from sklearn.model_selection import GridSearchCV param_grid { C: [1, 10, 100], gamma: [0.001, 0.01, 0.1], } svm SVC(kernelrbf, probabilityTrue, random_state42) grid GridSearchCV(svm, param_grid, cv5, scoringf1, n_jobs-1) grid.fit(X_train, y_train) print(best params:, grid.best_params_) print(best f1:, grid.best_score_)这里的probabilityTrue让LibSVM额外计算出概率形式的预测结果后面做混淆矩阵和预测图时直接使用predict_proba输出代价是训练时间会增加但对这个项目的数据量来说完全能接受。交叉验证cv5把训练数据切成5份轮流用其中4份训练、1份验证避免用同一批数据既调参又评估造成的过拟合。4.3 模型持久化与预测流程网格搜索结束后最优模型要连同前面的scaler、kmeans、lda一起保存下来不然测试阶段重新加载时模型之间会产生不一致。这是个特别容易忽略的命题LDA的主题编号、K-Means的簇编号、SVM的特征索引都是一一对应的任何一个模型重新训练过整条管线就断掉了。常见做法是把所有对象打包成一个字典用pickle或joblib存成一个文件。import pickle pipeline { scaler: scaler, kmeans: kmeans, lda: lda, svm: grid.best_estimator_, } with open(landslide_pipeline.pkl, wb) as f: pickle.dump(pipeline, f) # 预测阶段 with open(landslide_pipeline.pkl, rb) as f: loaded pickle.load(f) X_test_scaled loaded[scaler].transform(test_features) bow_test build_bow_from_kmeans(loaded[kmeans], test_features) topic_test loaded[lda].transform(bow_test) y_pred loaded[svm].predict(topic_test)我见过不少人在这个环节只保存SVM的模型文件然后测试时重新K-Means、重新LDA拟合结果预测准确率直接掉到随机猜测水平。要记住从特征提取到SVM每个环节都有训练和转换两个阶段训练时fit转换时transform全程保持同一套对象。5. 避坑排查特征到分类的六个常见翻车点与解决方法5.1 现象准确率卡在70%上下不去先描述症状训练集准确率能做到85%以上但验证集一直卡在70%左右无论怎么调C和gamma都没用。这个现象在特征差异大的项目里非常典型。原因基本是特征没有归一化。光谱特征取值0到255GLCM统计值在0到1附近SVM优化时大数值特征主导了距离计算纹理特征的贡献被压没了。解决方法是像前面写的那样在K-Means聚类之前对特征矩阵做StandardScaler标准化并且测试阶段用同一个scaler做transform不能在测试集上重新fit。做完这一步我的验证集准确率一般能直接拉升5到10个百分点。5.2 现象K-Means每次跑出的视觉单词都不一样还有一类问题是不固定随机种子造成的。K-Means的初始化是随机的LDA的求解也带随机性SVM训练同样有随机因素。如果这三处都不设置random_state每次运行的结果都会有差异有时精度甚至能差出好几个百分点。我要先说明这类现象不一定是代码写的bug而是算法本身的随机性。解决方式很简单就是在这三个模型的构造函数里都显式指定random_stateK-Means和SVM用random_state42LDA用random_state0。固定种子后如果两次运行结果仍有差异那就要怀疑数据切分时是否也用了随机。5.3 现象样本文件里的特征行数对不上训练时发现training.txt的行数和标签列表长度不一致或者某一行解析出的特征个数比预设维度多。原因是生成libsvm文件时特征向量的顺序没有和标签顺序对齐。最常见的是在读取图块时用了两个不同的目录遍历顺序一个目录按文件名排序另一个目录按文件修改时间排序两边顺序不一致标签和特征就错位了。解决方法是生成样本阶段就把「文件名、坐标、标签、特征」放在同一个DataFrame或列表里所有操作基于同一个索引不分别维护两个列表。5.4 现象训练集准确率极高测试集却远低于预期如果训练集准确率超过95%但测试集只有60%左右第一反应应该是过拟合。但在遥感场景里还有另一个常见原因数据泄漏。滑窗切图块的模式中同一张遥感影像切出的图块背景相似光照条件一致。如果直接随机划分训练集和测试集同一个影像切出的图块会同时出现在两边模型学到了「这张影像的整体风格」而不是「滑坡区域的视觉模式」分数自然虚高。解决方式是划分数据集时按影像ID进行分组整张影像的所有图块要么全进训练集要么全进测试集不要跨组混合。5.5 现象LDA训练时报负值错误有时LatentDirichletAllocation.fit会报“Negative values in data passed to LatentDirichletAllocation”或者训练完成后主题分布全是nan。这个现象的原因通常是输入给LDA的矩阵不是词袋计数而是被强行归一化成了小数或者某个特征列出现了负数。LDA的输入必须是矩阵元素非负的整数计数词袋直方图天然满足这个条件。但如果你在词袋之前做了StandardScaler特征值就会变成负数和浮点数这时候丢给LDA必然报错。解决方法是把标准化只用于K-Means聚类步骤词袋构建完成后不再做任何缩放。如果仍然报错用np.isnan检查一下词袋矩阵看看是不是某些图块在K-Means预测时产生了空簇。5.6 现象GLCM统计量返回nanGLCM特征提取时某个图块的所有GLCM统计值返回nan导致整条训练管线崩溃。原因是输入graycomatrix的数据类型和levels不对。skimage要求输入是uint8数组如果传入float64或者数据最大值超过了levels减1共生矩阵会越界或者稀疏异常graycoprops就会输出nan。解决方式是在进入GLCM前把图像数据先归一化到[0,1]再乘上(levels-1)然后转成uint8。另外levels不要设成25616或32足够灰度级越多共生矩阵越稀疏统计值越不稳定。这个坑在调试时最迷惑的地方在于它报错的位置距离真正出错的代码很远你可能会先怀疑分类器、再怀疑LDA最后才发现是特征提取环节的某个图块出了问题。我的排查习惯是从管线末端往上游逐层打印shape和合法性检查哪一层出现nan问题就在哪一层或它的下一层输入。6. 进阶验证混淆矩阵与预测类别图发现模型悄悄学错特征6.1 混淆矩阵逐类看错误只报告准确率在场景分类里是远远不够的尤其是滑坡这类正样本稀缺的任务准确率会被大量负样本稀释。我做完分类后一定会先看混淆矩阵逐类检查模型到底错在什么地方。from sklearn.metrics import confusion_matrix, classification_report import seaborn as sns import matplotlib.pyplot as plt y_pred grid.predict(X_test) cm confusion_matrix(y_test, y_pred) sns.heatmap(cm, annotTrue, fmtd, xticklabels[non-landslide, landslide], yticklabels[non-landslide, landslide]) plt.show()看混淆矩阵不是看对角线有多漂亮而是看非对角线的分布。滑坡被判成非滑坡意味着漏检这对灾害场景是致命的非滑坡被判成滑坡意味着虚警说明模型把某些地物纹理误当成了滑坡体。如果虚警集中在某一块下一步就要回到特征提取层看看这一块地物到底是什么样子而不是盲目调SVM参数。6.2 把预测结果映射回原始影像混淆矩阵只能告诉你错多少不能告诉你在哪里出错。遥感影像这个场景的特殊性在于错误往往是空间相关的模型可能把某条河流的冲积扇整体误判成滑坡。要看到这个现象最直接的做法是把每个图块的预测概率映射回原始影像坐标生成一张预测类别图叠加在原始影像上观察。这个映射图的价值在于它能直观暴露模型的学习偏好。我印象很深的一次测试混淆矩阵显示虚警率偏高把概率图叠加回原始影像后才发现模型把影像中一块裸土河漫滩判成了滑坡。河漫滩和滑坡体在光谱上确实非常接近差异主要在纹理的灰度共现结构上。看到类别图之后我返回去把GLCM窗口从7×7调整到9×9增加了correlation统计量虚警才压下去。如果不做这一层可视化验证这个问题会一直藏在指标背后。从那以后我每做完一版模型都会强制走一遍「混淆矩阵预测类别图」双验证一次都不跳过。这个习惯让我避开了至少三次看似精度达标、实则学习偏置的翻车状况。希望帮到你。本文还有配套的精品资源点击获取 遥感影像滑坡场景分类：从特征提取到SVM的完整实现与避坑指南
简介一份基于Python实现的遥感影像滑坡场景分类代码与项目文档面向毕业设计、课程设计及项目开发场景适合想快速掌握SVM分类流程并扩展应用的学习者。源码已经过严格测试可放心参考并在此基础上继续延伸。资源压缩包共21个文件、约1.86MB包括4个py源码、训练与测试txt数据、pkl特征/模型文件、lda模型相关文件及说明md目录按功能模块划分便于对照学习。整体流程清晰先从影像提取光谱统计特征与GLCM纹理特征再用K-Means聚类生成视觉单词接着通过LDA模型将词袋抽象为主题分布并保存为libsvm格式最后采用Libsvm完成SVM分类。项目文档对各环节原理与运行方式进行了说明便于理解特征工程—词袋模型—主题模型—分类器这条完整链路。目前已有56人学习下载用作课程报告或论文实验的参考基线都很合适。1. 遥感影像滑坡场景分类不是语义分割而是特征SVM的经典管线第一次拿遥感影像做滑坡场景分类时我的第一反应是上深度学习语义分割但实际一跑才发现场景分类要的是「这张图块是不是滑坡」不需要逐像元分割用语义分割反而把问题做复杂了。这个基于Python实现的项目走的是另一条完整可复现的路线从影像里提取光谱特征均值、标准差和GLCM纹理特征用K-Means把低层特征聚类成视觉单词再用LDA主题模型对视觉词袋做抽象把每张影像变成主题分布最后交给LibSVM分类器完成滑坡场景识别。整套流程纯CPU就能跑起来适合毕业设计、课程设计做闭环演示也适合想在遥感场景分类方向上做快速验证的从业者。资源里源码、中间模型文件和项目文档都齐顺着代码能一路推进到精度调优。2. 特征提取光谱统计量与GLCM纹理参数怎么定2.1 光谱特征为什么取均值和标准差特征提取是整条管线的底座也是最容易被跳过思考的一步。项目对每张影像块提取光谱特征取的是三个波段的平均值和标准差组合起来就是一个六维向量。为什么是这两个统计量而不是直方图或者更复杂的颜色描述子因为场景分类的输入是固定大小的图块均值能反映图块的整体亮度分布标准差能反映波段内的波动程度。滑坡体的裸土区域在可见光波段上通常亮度偏高、标准差偏大而同面积的草坡、林地方向均值低、波动小这两个统计量已经能把正负样本拉开。直方图虽然信息更充足但维度会膨胀到几十甚至上百对后续K-Means和SVM都不友好。从像素级到向量级的转换常见做法是先读入影像块再按通道统计。下面的代码是整个特征提取流程的第一段import numpy as np from PIL import Image def spectral_features(patch): 输入三波段影像块返回光谱均值与标准差特征6维 arr np.asarray(patch).astype(np.float32) means arr.mean(axis(0, 1)) # 每个波段一个均值 stds arr.std(axis(0, 1)) # 每个波段一个标准差 return np.concatenate([means, stds])这里arr的shape是(H, W, 3)axis(0,1)表示把所有像素在高度和宽度方向上做统计保留波段维。返回的向量顺序是三个均值接三个标准差调用方必须固定这个顺序否则后面生成libsvm格式时特征列会整体错位这是个最容易犯的低级错误。如果输入是单波段灰度图这段代码同样能跑只是输出的均值、标准差各一个维度变成2。光谱特征这一段虽然简单但在整个项目里承担的任务是「颜色分布的粗颗粒度描述」。后续的GLCM纹理特征负责描述空间关系两者互补之后一个图块的特征向量才具备区分滑坡和非滑坡的能力。2.2 GLCM纹理特征窗口、灰度级、方向与距离纹理特征这块项目选的是GLCM灰度共生矩阵的多个统计值。GLCM的核心思想很直接统计某个灰度值在指定方向和距离上与另一个灰度值相邻出现的频率。如果图块纹理粗糙像素对规律性强共生矩阵会呈现明显的主对角分布如果纹理细碎矩阵分布就均匀扩散。参数是GLCM落地时最需要较真的地方。我一般会把灰度级压缩到16或32距离取1方向取0度、45度、90度、135度四个方向然后对四个方向的统计量取均值。灰度级如果取256共生矩阵会变成256×256绝大多数格子都是零统计结果全是稀疏值方差大且不稳定距离取1能刻画相邻像素的最短关系取大了反而把纹理结构打成碎块。窗口大小通常用7×7或9×9窗口太小对边缘敏感窗口太大容易把两类地物的边界抹平。四个方向取均值是为了削弱影像旋转对纹理统计的影响这也是场景分类任务里很实用的做法。import numpy as np from skimage.feature import graycomatrix, graycoprops def glcm_features(gray_patch, levels16, distance1): 灰度图块输出GLCM四个方向统计值的均值特征 img (gray_patch / 255.0 * (levels - 1)).astype(np.uint8) angles [0, np.pi / 4, np.pi / 2, 3 * np.pi / 4] glcm graycomatrix(img, [distance], angles, levelslevels, symmetricTrue, normedTrue) stats [] for name in [contrast, dissimilarity, homogeneity, energy, correlation, asm]: prop graycoprops(glcm, name) # 输出 shape: (1,4) stats.extend(prop.mean(axis1)) # 四个方向取平均 return np.array(stats, dtypenp.float32)代码里两个关键点symmetricTrue让共生矩阵对称化这样0度和90度这类互为反向的组合会被合并矩阵更稳定normedTrue把频率归一化成概率避免图块尺寸不同导致数值量级不一致。graycoprops输出的每个统计量在四个方向上各有一个值对axis1取均值后返回的向量长度就是统计量个数乘4方向归一后的数量。这里有个容易踩的坑skimage的graycomatrix要求输入是uint8类型而且levels参数不能小于输入数组的最大值加1。如果之前做了归一化到[0,1]的浮点操作直接丢进去会报错或者得到全零矩阵。所以我在代码里先乘回[0,255]区间再转uint8这样levels16时数据能正确落在0到15的灰度级上。特征向量最终是6维光谱加6个GLCM统计量总计12维进入下一层K-Means时的特征空间并不算大这也是经典特征工程的典型做法。2.3 样本与标签生成滑窗取样和训练测试划分特征函数准备好之后下一步是把原始影像变成一批带标签的样本。项目里的generate_samples_and_labels.py干的就是这件事。场景分类一般用滑窗从大影像上切图块每个图块对应一个标签。滑坡场景是典型的小样本正类所以采样策略要特别设计正样本图块尽量全部保留负样本按正样本数量的一定比例采样否则SVM很容易被大量负样本带偏。图块尺寸常见做法是64×64或者128×128步长取等于窗口尺寸时样本之间没有重叠取小于窗口尺寸时有重叠相当于做了滑动增强。我建议先把图块坐标和标签存成中间文件再做训练测试划分这样后面调参时不用重新切图。import cv2 import numpy as np def extract_patches_with_label(image, patch_size64, stride64): 按滑窗切图块返回坐标列表标签由外部掩膜对应坐标决定 h, w image.shape[:2] patches [] for y in range(0, h - patch_size 1, stride): for x in range(0, w - patch_size 1, stride): patches.append((y, x, image[y:y patch_size, x:x patch_size])) return patches def make_train_test(patches, labels, test_ratio0.2): 按顺序切分训练测试注意先打乱索引 idx np.arange(len(patches)) rng np.random.RandomState(42) rng.shuffle(idx) split int(len(idx) * (1 - test_ratio)) train_idx, test_idx idx[:split], idx[split:] return train_idx, test_idx这里在切分前用固定RandomState打乱索引是保证实验结果可复现的关键。实际业务中更严格的划分是按影像文件划分而不是按图块划分同一张大图切出的图块共享光照和地物背景如果不按影像划分训练测试之间会有严重的信息泄漏精度虚高。这个细节在第五章会用专门一节展开讲。样本坐标和图块内容生成后要对每一块做特征提取。常见做法是对每个图块分别调用光谱特征和GLCM特征然后横向拼接成一个特征行所有图块的特征行纵向堆叠就形成了后续K-Means聚类的输入矩阵。这一步如果直接用Python循环跑速度会比较慢但项目样本量在几千的量级时完全够用。数据准备做到这里管线才走完第一层接下来进入视觉词袋阶段。3. 特征抽象K-Means视觉词袋与LDA主题分布为什么中间要多一层3.1 为什么低层特征不能直接喂SVM很多人拿到12维特征后会想既然维度不高为什么不直接训练SVM分类器还要绕一圈做K-Means和LDA答案是「一个图块的特征」和「一张影像的语义」之间还差着一个层级。每个64×64图块的特征描述的是局部颜色和纹理而一张遥感影像往往有一万个图块把这些图块特征直接平均成一条向量会丢失掉「哪些纹理模式在空间上频繁共现」这个重要信号。视觉词袋Bag of Visual Words的思路就是把影像类比成文档每个图块的特征向量经过K-Means聚类得到一个簇编号作为视觉单词一张影像内所有视觉单词的出现频次组成词袋直方图。这样影像被描述成「哪些视觉模式出现了多少次」而不是一堆没名字的连续特征。LDA主题模型再在这个词袋上做一层抽象输出多个主题的分布概率得到的就是项目描述里说的「高级特征」。这一层抽象的价值在分类精度上体现得很明显。直接用12维物理特征训练SVM碰上不同光照、不同季节的影像特征分布偏移很厉害而经过视觉词袋和主题分布两级抽象后分类器面对的是「语义模式比例」而非「像素统计量」泛化能力明显更强。3.2 K-Means聚类参数簇数量、初始化与标准化K-Means这一步把大量图块特征聚成若干个簇每个簇中心就是一个视觉单词。簇数量是关键参数常见做法是取200到500之间数量太小时不同地物被强行合并区分能力弱数量太大时每个簇里样本太少词袋直方图变得特别稀疏反而放大噪声。我一般先用256跑通流程再看后续分类精度的变化趋势决定是否需要调大。K-Means对特征尺度非常敏感。光谱特征取值在0到255之间GLCM统计值基本在0到1附近如果直接送进K-Means距离计算会被光谱维度主导等于纹理特征白提了。所以训练聚类器之前需要先用StandardScaler对特征矩阵做标准化。from sklearn.preprocessing import StandardScaler from sklearn.cluster import KMeans import numpy as np scaler StandardScaler() X_scaled scaler.fit_transform(feature_matrix) kmeans KMeans(n_clusters256, initk-means, n_init10, max_iter300, random_state42) kmeans.fit(X_scaled) word_ids kmeans.predict(X_scaled) # 每个图块一个视觉单词ID参数里n_init10表示用10个不同初始中心跑K-Means取最优结果这一步能有效避免聚类结果陷入局部最优random_state42固定了初始状态保证重复运行得到相同的单词划分。后续生成词袋时训练集和测试集必须用同一个kmeans对象做predict绝不能对测试集重新拟合聚类器否则视觉单词的编号语义就对不上了。词袋的构建在图像检索里也叫直方图量化。每张影像里统计每个簇编号出现的次数得到一个长度为n_clusters的整数向量from collections import Counter def build_bow(word_ids, n_clusters256): counter Counter(word_ids) bow np.zeros(n_clusters, dtypenp.float32) for wid, cnt in counter.items(): bow[wid] cnt return bow这个向量绝大部分位置是零属于高维稀疏表示。LDA模型接受的输入正是这种离散计数型矩阵所以词袋结果可以直接喂给下一层不需要额外做归一化。3.3 LDA主题模型从词袋到主题分布再写入libsvm格式LDALatent Dirichlet Allocation在这里是主题模型不是线性判别分析。这是一个特别容易混淆的坑sklearn里LinearDiscriminantAnalysis和LatentDirichletAllocation是两个完全不同的类前者是监督降维后者是无监督主题建模。项目文档里写的「利用LDA模型对视觉词袋进行主题分析」指的一定是后者。LDA把每张影像看成一篇文档视觉单词看成词主题数就是影像中可能存在的「隐含场景模式」。主题数太少会把不同地表类型揉在一起太多会产生大量重复主题。我做这个项目的经验是先从50个主题开始观察训练集和验证集精度的交叉验证曲线再决定加还是减。from sklearn.decomposition import LatentDirichletAllocation lda LatentDirichletAllocation( n_components50, learning_methodonline, max_iter30, random_state0 ) lda.fit(bow_matrix) topic_features lda.transform(bow_matrix) # 每张影像一个50维主题分布learning_methodonline适合样本量较大的情况训练速度快且迭代更稳定max_iter控制迭代轮数轮数太少主题分布还没收敛轮数太多耗时线性增长。topic_features每一行元素之和接近1含义是该影像在50个主题上的分布概率这个连续向量就是SVM要吃的最终特征。由于最后要训练LibSVM分类器需要把主题分布和类别标签写成libsvm文本格式。libsvm格式的约定是每行第一列为标签后面是“特征编号:特征值”编号从1开始连续递增def write_libsvm(topic_features, labels, output_path): with open(output_path, w) as f: for lab, vec in zip(labels, topic_features): items [{}:{:.6f}.format(i 1, v) for i, v in enumerate(vec)] line {} {}\n.format(int(lab), .join(items)) f.write(line)这段代码里标签用int强制转换防止numpy的float标签写入文件后变成“1.0”让LibSVM解析报错特征值保留6位小数就足够写太长文件会膨胀对精度没有实际帮助。项目里training.txt、test.txt以及training2.txt、test2.txt就是不同特征组合或样本划分下生成的libsvm文件classify.py读取这些文件完成后续训练。资源里的visual_lda.model、text_clt.pkl、spec_clt.pkl这些文件对应的正是这一层的中间产物。visual_lda.model是训练好的LDA模型各个clt.pkl是不同特征分支训练出的K-Means聚类器。把它们落盘的意义在于测试新影像时可以直接加载同一套模型保证单词编号和主题空间完全对齐。4. 分类器落地LibSVM格式、C/gamma搜索与模型持久化4.1 LibSVM文本格式与数据读取整个管线到了最后一步才进入分类环节但LibSVM数据格式的读取和转换反而是很多新手卡壳的地方。LibSVM工具包是SVM训练最常用的开源实现项目里生成的libsvm样本格式每行的结构是「标签 索引:值 索引:值...」索引从1开始而不是Python习惯的0。这个格式设计是为了支持稀疏存储对LDA输出的50维主题分布来说这个向量是稠密的不需要稀疏存储但仍然要按索引编号写入。读取libsvm格式的方式有两种。如果直接用LibSVM自带的命令行工具txt文件直接喂给svm-train就能训练如果整个管线都在Python里跑更常见的做法是用sklearn的svm.SVC它底层同样调用了LibSVM实现训练效果一致还能和前面的K-Means、LDA无缝衔接。分类器要读取libsvm文本先把文本转成稠密矩阵import numpy as np def load_libsvm(path, dim50): X_list, y_list [], [] with open(path, r) as f: for line in f: parts line.strip().split() y_list.append(int(float(parts[0]))) vec np.zeros(dim, dtypenp.float32) for item in parts[1:]: idx, val item.split(:) vec[int(idx) - 1] float(val) X_list.append(vec) return np.array(X_list), np.array(y_list)解析时每行的第一个元素要先用float再转int兼容标签写成“1”和“1.0”两种情况索引列要减1因为文件里是从1开始编号而numpy数组下标从0开始。如果dim设置得比文件里最大的索引还要小运行时会报索引越界反之则会在向量末尾留出一堆零这两种错误都不容易一眼发现。4.2 核函数选择与C/gamma参数搜索SVM分类器的核函数这个项目场景下直接选RBF径向基核就行。RBF核能把低维空间里线性不可分的样本映射到高维空间滑坡场景的主题分布特征虽然已经是高层抽象但正负样本的边界仍然是非线性的线性核在这个任务上精度会差一截。RBF核有两个核心参数C是误分类惩罚系数gamma控制RBF核的宽度。C太大容易过拟合C太小又欠拟合gamma太大会导致每个训练样本只影响自己附近很小的区域gamma太小则决策边界过于平滑。C和gamma怎么定我一般直接网格搜索。范围先按数量级大步长扫C取[1, 10, 100]gamma取[0.001, 0.01, 0.1]找到最优值附近再加密。如果正负样本不平衡网格搜索的scoring不能只盯着accuracy看滑坡场景里负样本占比高模型全部预测成负样本也可能有70%以上的准确率必须把scoring换成f1。from sklearn.svm import SVC from sklearn.model_selection import GridSearchCV param_grid { C: [1, 10, 100], gamma: [0.001, 0.01, 0.1], } svm SVC(kernelrbf, probabilityTrue, random_state42) grid GridSearchCV(svm, param_grid, cv5, scoringf1, n_jobs-1) grid.fit(X_train, y_train) print(best params:, grid.best_params_) print(best f1:, grid.best_score_)这里的probabilityTrue让LibSVM额外计算出概率形式的预测结果后面做混淆矩阵和预测图时直接使用predict_proba输出代价是训练时间会增加但对这个项目的数据量来说完全能接受。交叉验证cv5把训练数据切成5份轮流用其中4份训练、1份验证避免用同一批数据既调参又评估造成的过拟合。4.3 模型持久化与预测流程网格搜索结束后最优模型要连同前面的scaler、kmeans、lda一起保存下来不然测试阶段重新加载时模型之间会产生不一致。这是个特别容易忽略的命题LDA的主题编号、K-Means的簇编号、SVM的特征索引都是一一对应的任何一个模型重新训练过整条管线就断掉了。常见做法是把所有对象打包成一个字典用pickle或joblib存成一个文件。import pickle pipeline { scaler: scaler, kmeans: kmeans, lda: lda, svm: grid.best_estimator_, } with open(landslide_pipeline.pkl, wb) as f: pickle.dump(pipeline, f) # 预测阶段 with open(landslide_pipeline.pkl, rb) as f: loaded pickle.load(f) X_test_scaled loaded[scaler].transform(test_features) bow_test build_bow_from_kmeans(loaded[kmeans], test_features) topic_test loaded[lda].transform(bow_test) y_pred loaded[svm].predict(topic_test)我见过不少人在这个环节只保存SVM的模型文件然后测试时重新K-Means、重新LDA拟合结果预测准确率直接掉到随机猜测水平。要记住从特征提取到SVM每个环节都有训练和转换两个阶段训练时fit转换时transform全程保持同一套对象。5. 避坑排查特征到分类的六个常见翻车点与解决方法5.1 现象准确率卡在70%上下不去先描述症状训练集准确率能做到85%以上但验证集一直卡在70%左右无论怎么调C和gamma都没用。这个现象在特征差异大的项目里非常典型。原因基本是特征没有归一化。光谱特征取值0到255GLCM统计值在0到1附近SVM优化时大数值特征主导了距离计算纹理特征的贡献被压没了。解决方法是像前面写的那样在K-Means聚类之前对特征矩阵做StandardScaler标准化并且测试阶段用同一个scaler做transform不能在测试集上重新fit。做完这一步我的验证集准确率一般能直接拉升5到10个百分点。5.2 现象K-Means每次跑出的视觉单词都不一样还有一类问题是不固定随机种子造成的。K-Means的初始化是随机的LDA的求解也带随机性SVM训练同样有随机因素。如果这三处都不设置random_state每次运行的结果都会有差异有时精度甚至能差出好几个百分点。我要先说明这类现象不一定是代码写的bug而是算法本身的随机性。解决方式很简单就是在这三个模型的构造函数里都显式指定random_stateK-Means和SVM用random_state42LDA用random_state0。固定种子后如果两次运行结果仍有差异那就要怀疑数据切分时是否也用了随机。5.3 现象样本文件里的特征行数对不上训练时发现training.txt的行数和标签列表长度不一致或者某一行解析出的特征个数比预设维度多。原因是生成libsvm文件时特征向量的顺序没有和标签顺序对齐。最常见的是在读取图块时用了两个不同的目录遍历顺序一个目录按文件名排序另一个目录按文件修改时间排序两边顺序不一致标签和特征就错位了。解决方法是生成样本阶段就把「文件名、坐标、标签、特征」放在同一个DataFrame或列表里所有操作基于同一个索引不分别维护两个列表。5.4 现象训练集准确率极高测试集却远低于预期如果训练集准确率超过95%但测试集只有60%左右第一反应应该是过拟合。但在遥感场景里还有另一个常见原因数据泄漏。滑窗切图块的模式中同一张遥感影像切出的图块背景相似光照条件一致。如果直接随机划分训练集和测试集同一个影像切出的图块会同时出现在两边模型学到了「这张影像的整体风格」而不是「滑坡区域的视觉模式」分数自然虚高。解决方式是划分数据集时按影像ID进行分组整张影像的所有图块要么全进训练集要么全进测试集不要跨组混合。5.5 现象LDA训练时报负值错误有时LatentDirichletAllocation.fit会报“Negative values in data passed to LatentDirichletAllocation”或者训练完成后主题分布全是nan。这个现象的原因通常是输入给LDA的矩阵不是词袋计数而是被强行归一化成了小数或者某个特征列出现了负数。LDA的输入必须是矩阵元素非负的整数计数词袋直方图天然满足这个条件。但如果你在词袋之前做了StandardScaler特征值就会变成负数和浮点数这时候丢给LDA必然报错。解决方法是把标准化只用于K-Means聚类步骤词袋构建完成后不再做任何缩放。如果仍然报错用np.isnan检查一下词袋矩阵看看是不是某些图块在K-Means预测时产生了空簇。5.6 现象GLCM统计量返回nanGLCM特征提取时某个图块的所有GLCM统计值返回nan导致整条训练管线崩溃。原因是输入graycomatrix的数据类型和levels不对。skimage要求输入是uint8数组如果传入float64或者数据最大值超过了levels减1共生矩阵会越界或者稀疏异常graycoprops就会输出nan。解决方式是在进入GLCM前把图像数据先归一化到[0,1]再乘上(levels-1)然后转成uint8。另外levels不要设成25616或32足够灰度级越多共生矩阵越稀疏统计值越不稳定。这个坑在调试时最迷惑的地方在于它报错的位置距离真正出错的代码很远你可能会先怀疑分类器、再怀疑LDA最后才发现是特征提取环节的某个图块出了问题。我的排查习惯是从管线末端往上游逐层打印shape和合法性检查哪一层出现nan问题就在哪一层或它的下一层输入。6. 进阶验证混淆矩阵与预测类别图发现模型悄悄学错特征6.1 混淆矩阵逐类看错误只报告准确率在场景分类里是远远不够的尤其是滑坡这类正样本稀缺的任务准确率会被大量负样本稀释。我做完分类后一定会先看混淆矩阵逐类检查模型到底错在什么地方。from sklearn.metrics import confusion_matrix, classification_report import seaborn as sns import matplotlib.pyplot as plt y_pred grid.predict(X_test) cm confusion_matrix(y_test, y_pred) sns.heatmap(cm, annotTrue, fmtd, xticklabels[non-landslide, landslide], yticklabels[non-landslide, landslide]) plt.show()看混淆矩阵不是看对角线有多漂亮而是看非对角线的分布。滑坡被判成非滑坡意味着漏检这对灾害场景是致命的非滑坡被判成滑坡意味着虚警说明模型把某些地物纹理误当成了滑坡体。如果虚警集中在某一块下一步就要回到特征提取层看看这一块地物到底是什么样子而不是盲目调SVM参数。6.2 把预测结果映射回原始影像混淆矩阵只能告诉你错多少不能告诉你在哪里出错。遥感影像这个场景的特殊性在于错误往往是空间相关的模型可能把某条河流的冲积扇整体误判成滑坡。要看到这个现象最直接的做法是把每个图块的预测概率映射回原始影像坐标生成一张预测类别图叠加在原始影像上观察。这个映射图的价值在于它能直观暴露模型的学习偏好。我印象很深的一次测试混淆矩阵显示虚警率偏高把概率图叠加回原始影像后才发现模型把影像中一块裸土河漫滩判成了滑坡。河漫滩和滑坡体在光谱上确实非常接近差异主要在纹理的灰度共现结构上。看到类别图之后我返回去把GLCM窗口从7×7调整到9×9增加了correlation统计量虚警才压下去。如果不做这一层可视化验证这个问题会一直藏在指标背后。从那以后我每做完一版模型都会强制走一遍「混淆矩阵预测类别图」双验证一次都不跳过。这个习惯让我避开了至少三次看似精度达标、实则学习偏置的翻车状况。希望帮到你。本文还有配套的精品资源点击获取