简介：面向遥感、地信专业学生和研究者的最大似然法监督分类Matlab实践资源，以八波段遥感影像中的建筑物、道路、植被、水四类地物为对象，覆盖训练样本读取、分类器构建、像素分类及基于真实类别的精度评价全过程，旨在帮助理解监督分类中参数估计与精度验证的核心环节。包内共九个文件，包含Matlab源码脚本（.m）、训练样本与待分类像素的Excel数据表（.xls/.xlsx）以及结果说明与输出文档（.docx/.txt），压缩包整体约270KB，结构清晰便于按需查阅。目前已有5678人学习浏览。资源中已附有计算得出的总体精度、用户精度、制图精度与Kappa系数结果文件，使用者只需在代码中调整读取路径，即可在Matlab中直接运行并复现完整实验过程，非常适合作为遥感图像处理课程设计、实验报告或监督分类方法学习的配套材料。 算是个老话题了，但直到今天，每当有做遥感应用的小伙伴拿着影像来找我问“怎么分地类”，我脑子里第一时间跳出来的，往往还是最大似然法（Maximum Likelihood Classification, MLC）。不是因为它最炫，而是它在绝大多数场景下，是性价比和可解释性最均衡的一个选择。这篇用一篇完整项目的形式，把“最大似然法 + 监督分类”从原理到实操、从踩坑到优化，一次性讲透。
1. 整体设计与技术原理剖析
1.1 为什么在深度学习时代还要用最大似然法
先说个扎心的现实：很多刚入门的朋友，一上来就想着用U-Net、Random Forest、SVM，觉得“传统方法过时了”。但你真到生产环境下跑一遍就会发现，深度学习模型对样本量的要求、对硬件的要求、对数据预处理的要求，绝不是拿一块研究区随便跑跑就能落地的。反倒是最大似然法，在样本数量有限、类别光谱可分性尚可、计算资源紧张的情况下，经常能交出稳定且可解释的结果。
最大似然法的核心假设是：每一类地物的光谱特征（在各个波段上的DN值或反射率）服从多维正态分布。然后，对于影像中的每一个像元，我们计算它属于每一个类别的概率（似然值），最后把它判给概率最大的那一类。这个思路，本质上就是贝叶斯决策理论里“后验概率最大化”的简化版。
这个思路有个很直观的生活类比：你正走在一条街道上，闻到一股味道，你需要判断它是来自咖啡店还是面包房。如果你平时去咖啡店的频率更高、咖啡味在你记忆中的浓度分布更稳定，而面包房的糖香味波动又大，那你大概率会判定“这是咖啡味”。最大似然法做的事，就是把这种基于经验和统计的判断，用多维正态分布的概率密度函数精确算出来。
还有个常被忽略但非常关键的优势：最大似然法输出的不只是分类结果，它还会生成一张“规则/似然度”图层。你可以用这个图层判断哪些像元是“分得不自信”的，从而辅助人工修正或标注新样本。这个特性对一个生产型项目来说极其珍贵。
1.2 监督分类的整体逻辑流程
监督分类的基本逻辑可以归纳为四个字：先教后判。“教”就是人工选择训练样本（ROI，Region of Interest），告诉算法“哪块区域是水体”“哪块区域是农田”；“判”就是用算法计算整个影像上每个像元与这些样本统计特征的相似程度，然后完成归类。
完整的流程一般长这样：
数据获取与预处理。包括辐射定标、大气校正、几何校正，必要时要进行镶嵌和裁剪。
波段分析与特征选择。不是所有波段都参与分类效果就最好，噪声波段的加入反而会拉低精度。
训练样本采集。这是整个监督分类里最花时间、最影响精度的环节，没有之一。
统计特征参数计算。最大似然法需要计算每一类的均值向量、协方差矩阵等统计量。
分类执行。对每个像元计算概率，完成归类。
精度评价。用混淆矩阵、Kappa系数来量化分类效果。
后处理与制图。包括小图斑去除、类别合并、平滑滤波等。
需要特别强调的是：很多人觉得第3步样本采集“就是画几个多边形而已”，这是最大的误解。样本的质量和数量直接决定了协方差矩阵估得准不准，而协方差矩阵又是最大似然法的“命根子”——它一旦病态，后面所有概率计算都会失真。
1.3 最大似然法与其它监督分类方法的对比
为了让你更清楚它处在什么位置，我用一张表把主流的几种监督分类算法放在一起做个对比：
| 算法 | 核心思想 | 优点 | 缺点 | 适合场景 |
|---|---|---|---|---|
| 最大似然法 | 概率统计，多维正态分布 | 原理清晰、稳定、可解释性强，输出概率辅助信息 | 计算量大，对样本统计量敏感，对非正态数据适应性差 | 样本质量高、地类光谱差异明显、生产型项目 |
| 最小距离法 | 距离度量，欧氏/马氏距离 | 计算快、原理简单 | 对光谱方差和协方差结构不敏感，精度一般偏低 | 快速预分类、类别少且光谱紧凑 |
| 支持向量机 | 寻找最大间隔超平面 | 小样本下泛化能力较强，非线性能力好 | 参数调优复杂，核函数选择有门槛，结果可解释性弱 | 高维数据、样本不均衡、光谱重叠严重 |
| 随机森林 | 集成学习，多棵决策树投票 | 抗过拟合、能处理高维和非线性 | 模型较大，解释性弱，调参需要经验 | 特征维度高、类别复杂的中大型项目 |
| 深度学习（CNN等） | 多层神经网络自动提取特征 | 精度天花板高，能利用空间上下文信息 | 需要大量样本和算力，训练时间长，部署复杂 | 样本充足、研究型项目、复杂地物场景 |
从表里能看出来，最大似然法的核心优势区间是“样本有限但质量高、类别可分性尚可、需要稳定复现”的场景。换句话说，它更像是监督分类里的“稳健老手”，赢在综合稳定性而非单项上限。
2. 样本选择与预处理的核心细节
2.1 训练样本的“质”与“量”，如何平衡
关于最大似然法的样本选择，业内有一个经典的经验法则：每个类别的训练样本数至少要达到波段数的10到30倍。什么意思？如果你用的是Landsat 8的6个反射率波段（不算热红外），那每类至少要有60到180个像元。这个底线的本质是：协方差矩阵的估计需要足够的自由度。样本太少，协方差矩阵会奇异或近似奇异，导致计算时出现负数概率或者无法求逆，整个算法直接崩掉或用诡异方式输出。
但“量”不是唯一标准，“质”同样关键。这里说的质，包括三方面：
纯净性：样本尽量落在“纯粹”的地物像元上，避免混合像元。比如提取农田样本时，不要包含田埂、水渠、零星的裸土和树木，这些在影像上往往表现为“脏”的边界像元。
代表性：样本要覆盖该类地物在不同地理位置、不同地形条件、不同生长期的光谱差异。比如同一块水稻田，早期和晚期光谱完全不一样，如果你只在某一个时相上选样本，分类时就会漏掉其他时期的“水稻”。
独立性：训练样本和验证样本必须严格分开。有人一套样本既做训练又做验证，得出的精度高得吓人，但一拿到真实场景去用就翻车，这就是典型的“自欺欺人”。
我在实际项目中，一般会把训练样本和验证样本控制在7:3到8:2的比例，而且验证样本尽量用完全不同的图斑或区域来选。如果是小范围研究区，也可以用“留一法”去做交叉验证，结果更客观。
2.2 预处理中那些容易拉低精度的隐形坑
预处理的目的是让影像的光谱信息尽量真实、稳定地反映地物特征，但很多细节容易被忽视：
辐射定标的单位坑：有些影像的DN值范围是0到255，有些是0到65535（16位），还有些已经转成了表观反射率。最大似然法本身对数据的绝对尺度没要求，因为它用的距离是马氏距离，自带协方差归一化。但你如果混用了不同尺度来源的数据，比如一部分从USGS下载的Landsat 8 L1级产品，一部分是从别的平台拿到的表面反射率产品，那统计特征就乱了。
大气校正的必要性：做单幅影像分类时，有时候不做大气校正，直接用表观反射率也能分出个大概。但你只要想把这个分类器移植到另一景影像上，就会发现不同时相的大气和光照条件差异让光谱统计特征完全对不上。所以我一直强调，即使只在单景上分类，也建议至少做一次Quick Atmospheric Correction（QUAC）或FLAASH，让光谱特征更有物理意义。
波段选择不是越多越好：以Landsat 8为例，很多人觉得把所有波段都喂进去“总多点信息”，但实际效果往往适得其反。短波红外波段对土壤和植被的区别能力很强，但近红外和红光波段在植被覆盖度高的区域高度相关，同时放进模型会引入共线性问题。这会直接影响协方差矩阵的条件数，让矩阵更接近奇异，分类结果就更容易不稳定。具体能不能全波段用，建议先跑一次主成分分析或相关性矩阵，把高度相关的波段对找出来再做取舍。
地形阴影处理：在山区，阴坡和阳坡的同一地物光谱差异可能比不同地物之间的差异还大。这种情况下，不做地形校正，最大似然法会把很多阴坡的森林分成水体或阴影类。最简单的应急办法是把阴影单独作为一个类别参与分类，更彻底的办法是引入地形校正（如SCS+C校正）或地形阴影指数作为辅助波段。
2.3 训练样本的可视化检查与光谱可分性评估
在正式跑分类之前，强烈建议先做一步“预检”，也就是评估样本之间的光谱可分性。ENVI里自带一个工具叫“计算ROI可分离性”，输出的是Jeffries-Matusita距离和转换散度，值都在0到2.0之间。
经验判断标准很直接：
大于1.9：可分性优秀，放心跑。
1.7到1.9：可分性尚可，分类结果会出现零星错分，需要后续后处理。
1.4到1.7：可分性较差，结果会明显出现混分，建议重新选样本或合并类别。
小于1.4：这俩类别基本不可分，强行分类就是折磨自己。
我自己的操作习惯是：做完ROI后先看均值曲线，再看可分性报告，两者结合判断。均值曲线能直观反映“这一类在哪个波段上有辨识度”，可分性报告则用数值告诉你能不能分得开。
比如有一次做城市土地利用分类，我发现“草地”和“林地”的可分性只有1.35，查均值曲线后看到两类在近红外波段虽然有差异，但方差太大——因为公园草地和郊区草地的养护状态完全不一样。后来我把“草地”拆成“城市绿地”和“自然草地”两类，质量马上上来了。
3. 实操过程与核心环节实现
3.1 基于ENVI的完整分类实操
ENVI对最大似然法的支持非常成熟，我一般这么操作：
第一步：数据准备
把预处理后的影像加载进ENVI，建议将多波段文件按标准格式（如GeoTIFF）保存。如果影像比较大（例如超过2GB），优先考虑用ENVI的DAT格式，读写效率更高。打开影像后先做个快速拉伸显示，用真彩色合成的方式初步观察地物分布。
第二步：ROI绘制
在ROI Tool里新建类，类名建议用英文或拼音（中文容易在后续批处理和二次开发中出编码问题），比如Water、Forest、Cropland、BuiltUp、BareLand。绘制时用多边形工具沿着地物边界小心勾画，尽量避开混合区域。每个类别至少画3到5个分散的图斑，覆盖研究区内的不同位置。
第三步：计算可分离性
用“ROI可分离性”工具，勾选所有类别执行计算。此时要重点看两两组合值。处理方法是：低于1.4的类别要么合并、要么重新选样本、要么增加波段信息。在这个环节多花半小时，后面能省半天。
第四步：执行最大似然法
在工具箱中搜索“Maximum Likelihood Classification”，打开界面后选择输入影像和ROI文件。这里有个参数叫“似然度阈值”，默认是None，意思是所有像元必须归到某一类。但生产项目中我建议设定一个阈值（比如0.05或0.01），让那些对任何类别都没有足够把握的像元单独分到一个“未分类”类。这样后续检查和处理起来会从容很多。
第五步：精度评价
生成分类结果后，用“Confusion Matrix Using Ground Truth ROIs”做精度评定，选择之前留出的验证样本。重点关注Overall Accuracy和Kappa Coefficient，前者反映整体正确率，后者排除了随机一致性的干扰，更能体现分类器的真实水平。
第六步：后处理
用Majority/Minority Analysis去除小图斑，用形态学开闭运算平滑边界。注意核大小的选择，3×3到5×5通常够用，过大的核会抹掉细碎但重要的地物边界。
3.2 用Python从零实现最大似然分类
如果你不想完全依赖商业软件，或者希望在特定项目里把分类器嵌入自动化流程，用Python自己写一个最大似然分类器其实并不复杂。核心步骤如下：
import numpy as np from scipy.stats import multivariate_normal from osgeo import gdal # 1. 读取影像和样本像元 def read_image(path): ds = gdal.Open(path) data = ds.ReadAsArray() # shape: (bands, rows, cols) rows, cols = data.shape[1], data.shape[2] data_2d = data.reshape(data.shape[0], -1).T # shape: (pixels, bands) return data_2d, rows, cols # 2. 计算每类的均值向量和协方差矩阵 def fit_class_distribution(samples, labels): classes = np.unique(labels) class_params = {} for c in classes: class_samples = samples[labels == c] mean_vec = np.mean(class_samples, axis=0) cov_mat = np.cov(class_samples, rowvar=False) # 为了防止奇异矩阵，加一个小的正则项 cov_mat += np.eye(cov_mat.shape[0]) * 1e-6 class_params[c] = (mean_vec, cov_mat) return class_params # 3. 最大似然分类 def ml_classify(image_2d, class_params, threshold=0.01): all_classes = list(class_params.keys()) # 用零向量初始化一个足够小的概率值 prob_matrix = np.zeros((image_2d.shape[0], len(all_classes))) for idx, c in enumerate(all_classes): mean_vec, cov_mat = class_params[c] rv = multivariate_normal(mean=mean_vec, cov=cov_mat, allow_singular=True) prob_matrix[:, idx] = rv.pdf(image_2d) # 找出每个像元的最大概率 max_prob = np.max(prob_matrix, axis=1) max_class_idx = np.argmax(prob_matrix, axis=1) # 低于阈值的归为未分类 max_class_idx[max_prob < threshold] = -1 return max_class_idx有几个细节值得展开说一下：
np.cov默认是除以n-1，这是无偏估计，没问题。但如果某个类别的样本非常少（比如少于波段数），协方差矩阵就会变成奇异矩阵。这时加正则项（上面代码中的1e-6）只是应急方案，最根本的办法还是回去补充样本。multivariate_normal.pdf用的是概率密度值而不是概率值，所以得到的“概率”可能大于1。这不影响分类决策，因为我们是比较相对大小，不是看绝对数值。但如果要用阈值来定“未分类”，阈值本身需要根据实际概率密度值的数量级来调整。你的影像如果是0到255的DN值范围，概率密度值可能极小；如果是表观反射率（0到1），密度值又会很大。实践上我一般先跑一次，输出max_prob的直方图，再决定合理的阈值区间。allow_singular=True这个参数要慎用，它会把近似奇异的矩阵强行求逆，结果可能不准确。生产环境里宁可报错提示样本不足，也比输出一个看似正常实则不可靠的结果好。
3.3 关键参数的选择与调优建议
最大似然法看似“无参”，实际上还是有几个决策点直接影响结果质量：
似然度阈值：这个阈值控制了分类器对“低置信度像元”的处理方式。阈值设得越小，越多的像元会被“勉强”分到某一类，整体的未分类像元越少；设得越大，未分类像元越多，但每一类的可信度更高。我的经验是：初始设为0.01，看结果直方图后微调，目标是在“分类完整性”和“准确率”之间找到平衡点。
协方差矩阵的处理方式：有些软件允许你选择“使用完整协方差矩阵”或“只使用对角线（即各波段独立）”。绝大多数情况下应该用完整协方差矩阵，因为它捕捉了波段之间的相关性信息。只有当样本量不足时，才退而求其次，用对角线版本。
先验概率：最大似然法理论上允许设置各类别的先验概率。比如你已知研究区里水体占了20%，农田占了40%，那么按面积比例设置先验概率会让分类结果更符合实际情况。但使用先验概率要谨慎——它会把分类结果向先验“拉扯”，如果先验本身不准，反而会降低精度。初学阶段建议先用等先验（Equal Priors），跑出基线结果再说。
4. 常见问题与排查技巧实录
4.1 分类结果出现大量“椒盐”噪声，怎么办
现象：结果图上散布着密密麻麻的孤立像元，和周围类别不一致。这在最大似然法中非常常见，原因主要是像元级分类没有考虑空间上下文，加上地物光谱本身的随机波动。
排查步骤和解决办法：
先用3×3或5×5窗口的Majority Filter做一次众数滤波，能消掉大部分孤立点。注意不要重复多次大核滤波，否则边界会变得非常圆滑，反而失真。
如果滤波后仍然噪声严重，回头检查样本的光谱方差。常见问题是某个类别的训练样本里混入了异质成分，导致协方差矩阵过大，特征空间里“四处开花”。
如果以上两步都不奏效，考虑是否某些类别本身在光谱上就纠缠不清。此时最好的办法是增加辅助数据，比如引入NDVI波段、DEM或者纹理特征，增大维数来提升可分性。
4.2 某些类别整体错分，尤其是“阴影”和“水体”
这是几乎所有遥感分类项目都会撞上的经典问题，因为阴影和水体在可见光波段的特征非常相似——都表现为低反射率、在近红外波段也偏低。
我的排查方法分两步：
第一步，查看均值曲线。在ENVI里把“水”和“阴影”的训练样本均值曲线叠加显示，如果两个曲线几乎重合，那基本没有通过算法挽救的余地。
第二步，增加辅助信息。常见做法有三种：
引入DEM和太阳入射角，计算地形阴影指数，把它作为一个新波段参与分类。
增加纹理特征，比如基于灰度共生矩阵（GLCM）的方差和对比度波段。阴影区域往往纹理平滑，水体在大型湖泊中也平滑，但水体的空间延展性要强得多，纹理上的差异能帮上忙。
如果数据源允许，结合短波红外（SWIR）波段。虽然水体在SWIR上是强烈吸收的，阴影在SWIR上的反射率会受大气散射影响，两者有一定细微差异。这个差异不够大，但在部分场景下能起到关键的一票。
4.3 分类精度验证不理想，但目视效果“看着还行”
这种场景也很常见：输出的图上地物轮廓分明、颜色鲜艳，看起来完美无瑕，一算混淆矩阵却发现Overall Accuracy只有70%多。
这时候通常问题出在验证样本上，而不是分类器本身。常见坑：
验证ROI和训练ROI距离太近，甚至重叠。由于地物光谱在空间上有自相关性，邻近像元的光谱高度相似，这会导致验证结果虚高；反之，如果验证样本选在与训练样本完全不同的地区，由于地域光谱差异，结果又会虚低。合理的做法是验证样本在空间上独立于训练样本。
验证样本类别定义模糊。比如“裸地”和“滩涂”验证样本边界不清，评估时互相计为误差，混淆矩阵就很难看。这时候重新审视类别体系的定义，做到互斥、完备，比闷头调参数更有效。
样本量不足导致精度估计方差过大。验证样本太少（比如每类只有20个像元），精度估算的置信区间会非常宽，结果的偶然性也大。建议每类至少50到100个验证像元。
4.4 两张不同时相的影像用同一套样本，效果天差地别
有人做过一个项目：把某地区5月份的Landsat影像分成两景，一景做训练和分类，另一景想复用同样的样本，结果第二景的分类精度惨不忍睹。原因很简单：两个月之间的植被物候差异让“农田”“林地”的光谱特征发生了显著偏移。
应对思路有几种：
如果用多时相数据做分类，尽量选择物候差异最小的时相，或在每个时相重新选择样本。
如果硬要用一套样本，必须做归一化处理，比如将DN值转化为地表反射率，并确保两景影像的大气校正参数一致性。
稳妥的方案是采用“迁移学习”的思路：先做光谱归一化（如直方图匹配），再复用样本。当然这不属于最大似然法的标准能力范围，需要额外的数据预处理工作。
4.5 分类结果里出现“未分类”空洞，怎么补
设置了似然度阈值之后，会有一部分像元因为概率低于阈值而变成“未分类”。这在生产制图中是个麻烦事。我的处理方式分两步：
先看未分类像元的空间分布。如果它们集中分布在混合像元较多的过渡地带或阴影区，说明这是正常的“光谱模糊区”，可以用邻域类别众数填充。
如果未分类像元分散得毫无规律，那可能是影像上存在云、阴影、条带或传感器噪声等异常区域。这些像元的光谱特征完全超出正常地物的分布范围，不应强行归类。常规做法是先定义Mask剔除异常区域，再做分类，最后在制图时单独标注，而不是强行填充。
5. 从单幅分类到批量生产的工程化建议
5.1 批处理时如何保证结果一致性
如果你面对的不是一幅图，而是全省甚至全国的影像分类任务，最大的挑战不是算法本身，而是“一致性”。同一套分类规则，在前一景表现良好，到后一景可能就崩了，这种“崩”不是随机误差，而是数据源差异导致的系统性偏移。
解决办法是在工程层面建一套标准化流程：
明确统一的预处理参数（大气校正模型、云掩膜阈值、地形校正策略）。
建立完好的样本库管理机制，允许不同区域微调样本，但保持类别体系的一致。
为每一景影像自动生成分类后的质量报告，包括总体精度、Kappa系数、各类别的制图精度和用户精度，并用阈值自动化标记“疑似问题影像”，让人工只检查标记出来的部分。
这套流程做下来，单景分类的时间可能没有缩短，但整体交付效率和质量的稳定性会提升非常明显。
5.2 最大似然法的sample库管理与扩展
在实际项目里，我们把ROI样本当成项目资产来管理，而不是散落在各个工程文件里。具体做法：
将样本文件统一存储为Shapefile或GeoJSON，几何和属性独立，方便跨平台使用。
每个样本属性字段记录：类别编码、类别名称、采集日期、影像时相、数据源标识、采集团队等信息。这样当发现某一批样本影响了精度时，可以快速追踪问题源头。
样本版本化更新。每次分类项目结束后，把经人工修正的“高置信度像元”增量加入样本库，让样本库越用越强大。这里要注意，增量加入前必须做好光谱一致性检查，防止把错误像元引入。
5.3 结合其他算法的融合策略
用最大似然法作为“主体框架”、用SVM或RF作为“补充校验”，是我在业务中经常采用的策略。具体操作是：先用最大似然法完成主体分类，然后在“光谱模糊区”（比如概率小于某个阈值的像元）或特定类别对（比如易混类别之间）上，引入SVM或RF重新判断。这样既保留了最大似然法的稳定性和可解释性，又借助更灵活的算法降低了混分率。
还可以用最大似然法的概率输出做软分类，把每个像元的类别归属概率作为中间数据，再做亚像元尺度的混合像元分解或者生态参数反演。这属于把最大似然法从“分类器”升级成“估算器”的思路，在植被覆盖度估算、不透水面比例估算等场景中有实际价值。
6. 实操中的一些个人心得
用最大似然法做了这些年项目，有几个体会比较深。
第一，这个算法对“干净数据”的依赖远大于对“复杂算法”的依赖。很多分类效果不好的根源，早就在数据预处理和样本设计阶段埋下了。与其纠结换什么高级算法，不如先把样本和预处理打磨扎实，往往能收到立竿见影的效果。
第二，最大似然法的概率输出是一个被严重低估的信息。多数人只用了最终的分类结果，却忽略了“概率图”里藏着的价值。比如在一幅土地利用分类图中，概率低的区域往往对应混合像元较多的地方，这些地方恰恰是野外核查最需要优先布点的位置。用这个思路来指导野外采样，效率能翻倍。
第三，没有万能的分类器，只有合适的流程。最大似然法在光谱分离度好的场景下，精度不输很多机器学习的算法，而且你完全可以讲清楚每一个决策的根据——这在撰写技术报告和对接业务人员时，真的非常有用。业务方通常更想听到“这个像元因为光谱特征接近水体，所以分为水体”，而不是“模型根据256维特征综合打分得出结果”。
最后想说的是，无论你手头最终选用了哪种分类器，都建议把最大似然法当作理解监督分类的一个“锚点”。把它的原理吃透了，再去看SVM、随机森林、深度学习等方法，你会更容易看懂这些算法在解决哪些问题、牺牲了哪些东西。基础打牢了，上层建筑才稳。
本文还有配套的精品资源，点击获取
