【遥感学习整理 05】水质火点森林农业与矿物蚀变遥感监测

水质、火点、森林、农业与矿物蚀变遥感监测

📝 摘要： 本文以湖泊叶绿素反演、MODIS秸秆焚烧监测、林冠健康变化、农业耕作变化和Crosta矿物蚀变提取五类案例为主线，重点解释各类方法为什么成立、数据为什么要这样处理、阈值和模型为何不能直接照搬，以及如何把ENVI操作转化为可验证、可复现的遥感分析流程。

关键词： 水质反演、MODIS火点、GLT、林冠变化、分类后比较、Crosta、PCA、遥感质量控制

01 从五个案例理解应用型遥感的共同结构

1.1 五个专题不是五套互不相干的软件操作

湖泊水质、秸秆焚烧、林冠状态、农业耕作和矿物蚀变分别属于水环境、大气环境、森林生态、农业生产和地质调查。若只按照菜单顺序学习，它们会表现为五套彼此割裂的ENVI操作；若从遥感信息形成机制出发，它们实际上都在完成同一件事：把传感器记录的辐射差异转换为具有业务含义的空间对象或连续参数。

水质专题将水体光谱差异转换为叶绿素a浓度；火点专题将热异常像元与农用地叠加，形成秸秆焚烧候选点；林冠专题将多时相植被指数差异转换为健康下降区域；农业专题把两个时期的分类结果转换为耕作类型转移；蚀变专题则将波段间的吸收与反射差异压缩到特定主成分中，再识别异常区域。

共同流程并不意味着可以用同一套阈值处理所有任务。不同专题的差异主要来自目标信号的形成机制：

专题	被观测对象	主要信号	主要干扰	输出
湖泊水质	水体组分	可见光—近红外反射变化	大气、邻近效应、风浪、悬浮物	连续浓度
秸秆焚烧	高温火点	中红外和热红外异常	工业热源、云边、林火	点状候选目标
林冠状态	森林冠层	多时相NDVI变化	物候、阴影、干旱、配准误差	状态变化等级
农业耕作	耕地利用状态	光谱、纹理与物候	裸地混淆、收获时间差异	类别转移
矿物蚀变	蚀变矿物组合	波段吸收和反射差异	植被、水体、岩性背景、大气	异常区

1.2 应用型遥感的真正难点在“解释链条”

影像中出现一个高值或低值，并不自动具有业务意义。一个高温像元可能是秸秆焚烧，也可能是工业烟囱；一片NDVI下降区可能是病虫害，也可能是收获、干旱或阴影；一个主成分异常区可能对应蚀变，也可能只是裸岩、盐碱地或地形照度差异。

因此，应用型遥感必须建立完整解释链条：

1. 目标为什么会改变某些波段或特征；
2. 传感器是否能够观测这种变化；
3. 预处理是否消除了主要非目标因素；
4. 模型输出与业务对象之间是否存在明确映射；
5. 是否有独立数据证明这种映射可靠。

::: important 判断遥感方法是否成立
一个方法不仅要“能算出结果”，还要说明结果为何与目标有关、在哪些条件下有关、在哪些条件下会失效。
:::

1.3 五类案例对应五种典型推理方式

这五个专题几乎覆盖了传统应用遥感的主要推理模式：

• 水质反演属于从光谱特征到连续参数的经验回归；
• 火点监测属于产品编码解释与空间条件筛选；
• 林冠监测属于同类地物内部状态的连续变化检测；
• 农业监测属于离散类别的分类后比较；
• 蚀变提取属于有目标约束的统计异常检测。

理解这些推理方式，比记住某个版本ENVI中的工具路径更重要。软件版本会变化，但模型假设、误差来源和验证逻辑不会因为界面变化而消失。

---

02 湖泊水质遥感：从水体光谱到叶绿素a浓度

2.1 水体遥感为什么比陆地遥感更依赖预处理

陆地植被和裸土在可见光与近红外中通常具有较强反射，而清洁水体的离水辐射较弱。传感器接收到的水体像元中，大气路径辐射、邻近陆地散射和水面镜面反射可能占有很高比例。也就是说，影像中“看见的水”并不完全是水体自身发出的信号。

传感器接收到的水体辐射可概念化表示为：

$$L_{sensor}=L_{path}+T_{up}{L}_{water}+L_{adj}+L_{glint}$$

其中，$L_{path}$为大气路径辐射，$T_{up}{L}_{water}$为水体离水辐射经大气传输后的部分，$L_{adj}$为邻近陆地贡献，$L_{glint}$为太阳耀斑或镜面反射。

当目标信号本身较弱时，任何一个附加项处理不当都会显著改变波段比值。因此，水质反演不能把“大气校正”理解为一般性的图像美化步骤，它决定了模型输入是否仍然代表水体光学性质。

2.2 叶绿素a与水色变化之间的关系

叶绿素a是浮游植物光合作用色素，也是评价水体富营养化和藻类生物量的重要指标。叶绿素a在蓝光和红光区域具有较强吸收，在绿光附近形成相对反射峰；当藻类浓度升高时，红光吸收、绿光反射和红边附近的谱形都可能发生变化。

但在太湖等内陆复杂水体中，叶绿素a并不是唯一影响水色的组分。悬浮泥沙、有色可溶性有机物、浅水底质和藻类漂移都可能同时改变反射率。这类水体常被称为光学复杂水体，同一叶绿素浓度在不同悬浮物条件下可能产生不同的光谱响应。

这意味着经验模型并不是“叶绿素浓度的普适物理定律”，而是特定水体、特定季节和特定传感器条件下的统计映射。

2.3 环境小卫星案例的预处理顺序

案例使用HJ-1B CCD影像，并以太湖叶绿素a为目标参数。原始流程包括数据读取、辐射定标、几何校正、大气校正、太湖区域裁剪、星地样本匹配、模型建立和结果验证。

先裁剪再配准和大气校正，主要是为了减少计算量。但裁剪范围不能仅保留水面本身。大气校正和自动配准往往需要一定的地物多样性和空间上下文，裁剪过小可能削弱稳定特征点数量或场景统计代表性。

2.4 为什么几何精度会直接影响反演模型

水质采样点通常由GPS坐标表示，而卫星像元代表一个具有面积的空间单元。若影像存在半个到一个像元的位置误差，采样点可能落入岸边混合像元；若水体近岸分布着水草、码头、船舶或高浑浊河口，提取到的光谱将与实测水样不再对应。

因此，星地匹配不应只是“把点叠加到影像上并读取一个像元”。更稳健的方法包括：

• 检查采样点是否位于纯水像元；
• 以采样点为中心提取$3\times 3$或$5\times 5$窗口；
• 剔除窗口内离群像元；
• 使用中位数或均值代表卫星观测；
• 对岸边样点单独标记；
• 检查卫星过境与采样时间差。

窗口不能无限增大。窗口越大，越容易混入空间异质性；窗口过小，则容易受定位误差和单像元噪声影响。

2.5 波段比值模型如何建立

案例采用近红外与红光比值作为自变量：

$$R_{NR}=\frac{B_{NIR}}{B_{RED}}$$

叶绿素a经验模型写为：

$$Chla=a{R}_{NR}+b$$

案例给出的示例系数为：

$$Chla=0.0424R_{NR}+0.0164$$

这个公式的意义是：在案例样本范围内，近红外与红光比值增大时，叶绿素a浓度呈线性增加趋势。系数$a$控制变化速度，截距$b$表示比值为零时模型外推的理论值。由于比值为零通常不在真实样本范围内，截距不应被解释为具有明确水体物理意义。

波段比值的优势在于能够削弱共同照度变化，但它也可能放大红光分母较小时的噪声。因此，建模前应检查：

• 分母是否接近零；
• 是否存在负反射率；
• 比值是否出现极端值；
• 高浓度样本是否对回归线形成过强杠杆；
• 残差是否随浓度增大而增大。

2.6 为什么不能只看决定系数

散点图中的$R^2$只能描述样本内线性拟合程度，并不能证明模型具备空间外推能力。一个样本量很小、浓度范围很窄的模型，也可能获得较高$R^2$。

水质模型至少应报告：

$$RMSE=\sqrt{\frac{1}{n}\sum _{i=1}^n({\hat{y}}_i-y_i{)}^2}$$

$$MAE=\frac{1}{n}\sum _{i=1}^n|{\hat{y}}_i-y_i|$$

$$Bias=\frac{1}{n}\sum _{i=1}^n({\hat{y}}_i-y_i)$$

其中，RMSE对大误差更敏感，MAE反映平均绝对偏差，Bias用于识别系统性高估或低估。验证时应使用未参与建模的独立样本，不能把回归拟合样本再次作为“验证点”。

::: warning 经验模型的迁移限制
案例模型只能代表特定影像、特定采样期和特定水体条件。更换季节、传感器或湖泊后，应重新进行星地匹配与参数拟合。
:::

2.7 从浓度图到水环境解释

叶绿素a空间分布图生成后，还需要把连续数值放回湖泊水动力和污染输入背景中解释。高值可能与湖湾封闭、河口营养盐输入、风场驱动藻类聚集或局部水华有关；低值可能出现在开阔水面，也可能是薄云、耀斑掩膜或大气校正异常造成的假低值。

合理的后续分析包括：

• 按湖区或管理单元统计均值与高值面积；
• 对比河口、湖湾和湖心；
• 结合风向、风速和水位解释藻类漂移；
• 对高值异常区回看原始波段与真彩色影像；
• 与同期监测站数据进行空间验证。

2.8 水质反演中的波段选择不能只看相关系数

在经验模型建立阶段，研究者常会计算不同波段或波段组合与叶绿素a浓度的相关系数，并选择相关性最高的特征。但相关性高并不等于模型稳定。某个比值可能只是碰巧与当前样本中的悬浮物梯度一致，当水体成分结构发生变化后，关系就会失效。

波段选择至少要同时考虑三方面：

第一，目标参数的生物光学机理。叶绿素a的吸收和散射特征决定哪些波段具有解释基础。若某个波段组合在光谱机理上无法解释，即使样本内相关性很高，也可能只是偶然关系。

第二，传感器波段宽度和响应函数。环境小卫星CCD波段并不是某个单一波长的观测，而是传感器响应函数覆盖范围内的加权积分。地面光谱仪上的窄波段特征，未必能在卫星宽波段中完整保留。

第三，其他水体组分的干扰。悬浮物、有色可溶性有机物和水深可能与叶绿素浓度共同变化，造成虚假相关。若条件允许，应同步记录浊度、悬浮物和透明度，分析模型是否真正响应叶绿素。

因此，特征筛选不能只输出一个最高相关系数，还应比较不同模型在独立验证集中的误差、残差结构和跨区域稳定性。

2.9 水体邻近效应为何在太湖边缘尤其明显

邻近效应是指传感器观测某个水体像元时，周围高反射陆地经过大气散射后对该像元产生额外贡献。太湖岸线附近存在城市、农田、堤岸和水生植被，这些目标在近红外中的反射通常远高于水体，因此近岸水体像元容易被“抬亮”。

邻近效应可能造成两类问题：

• 近红外反射率异常升高，使$B_{NIR}/B_{RED}$比值偏大；
• 高浓度区域看似沿岸分布，但实际部分来自岸边散射影响。

处理时可以设置一定宽度的岸线缓冲区，将最靠近陆地的若干像元排除；也可以比较近岸与湖心样本的残差，判断模型是否存在明显空间偏差。对于分辨率较低的数据，岸线缓冲区通常比简单水体掩膜更重要。

2.10 回归模型应检查残差，而不是只输出拟合直线

残差定义为：

$$e_i=y_i-{\hat{y}}_i$$

若残差在低浓度和高浓度区随机分布，说明线性模型可能基本适用；若浓度越高，残差幅度越大，则可能存在异方差；若残差呈弯曲趋势，则说明线性关系不足，可能需要对数、指数或分段模型。

需要重点检查：

• 是否存在少量样本决定整条回归线；
• 高浓度水华样本是否超出模型线性范围；
• 残差是否与采样位置、近岸距离或悬浮物相关；
• 训练集与验证集浓度范围是否一致；
• 模型是否出现负浓度预测。

若模型预测出现负值，不应只把负值截断为0，还应判断模型是否在低浓度区发生外推。

---

03 MODIS秸秆焚烧监测：从热异常到农田火点

3.1 火点遥感识别的物理基础

燃烧目标的温度远高于周围地表。根据普朗克辐射定律，高温目标在中红外和热红外波段会产生强烈辐射。即使火焰只占一个像元很小比例，也可能显著抬高该像元的中红外亮温。

一个混合火点像元可概念化表示为：

$$L_{\lambda }=f{L}_{\lambda ,fire}+(1-f)L_{\lambda ,bg}$$

其中，$f$为火点面积比例，$L_{\lambda ,fire}$为火焰或高温目标辐射，$L_{\lambda ,bg}$为背景辐射。由于高温辐射增长很快，较小的$f$也可能形成可检测异常。

但“热异常”不是“秸秆焚烧”的同义词。工业烟囱、油气燃烧、森林火灾、垃圾焚烧甚至部分高温裸地都可能产生类似信号。秸秆焚烧识别需要增加土地利用和空间背景约束。

3.2 MOD14、MOD03与MCD12Q1分别解决什么问题

案例使用三个MODIS相关产品：

产品	提供的信息	在流程中的作用
MOD14	热异常与火点编码	确定可能存在燃烧的像元
MOD03	经纬度和观测几何	对MOD14进行地理定位
MCD12Q1	土地覆盖分类	判断火点是否位于农用地

这三个产品解决的是不同问题。MOD14回答“哪里有热异常”，MOD03回答“这些异常在地球上的什么位置”，MCD12Q1回答“这些位置属于什么地表类型”。

3.3 为什么需要GLT

MODIS扫描数据不是简单的规则地图网格。原始数据中的每个像元具有对应经纬度，像元形状和位置会随扫描角度变化。GLT，即Geometry Lookup Table，通过经纬度数组建立原始行列号到目标地图网格的查找关系。

使用GLT时必须保证MOD14和MOD03来自同一过境时刻。若时间不匹配，经纬度数组与火点数据的扫描结构不同，结果会出现整体错位或局部扭曲。

3.4 产品编码不是一般影像灰度

案例根据产品说明，将DN值7、8、9解释为火点等级。这里的DN不是普通辐射亮度，而是分类编码。使用这类产品时，必须先阅读数据集属性，明确：

• 每个编码代表什么；
• 哪些值表示无数据；
• 哪些值表示云；
• 哪些值表示低、中、高置信火点；
• 是否需要同时读取质量控制层。

若把分类编码当成连续物理量进行插值，会生成不存在的中间类别。因此，地理校正和重投影时应优先使用最近邻方法。

3.5 农用地提取与空间约束

MCD12Q1是全球土地覆盖产品，案例根据分类说明提取农用地编码。经过地理定位、投影变换和镶嵌后，农用地掩膜与火点栅格进行空间相交。

设火点掩膜为$F(x,y)$，农用地掩膜为$A(x,y)$，则候选秸秆焚烧点为：

$$S(x,y)=F(x,y)\times A(x,y)$$

只有当两个掩膜均为1时，结果才为1。

这一逻辑看似简单，但空间分辨率差异会引入不确定性。火点像元和土地覆盖像元可能大小不同，一个火点像元可能覆盖农田、村庄和道路。若直接按像元中心判断，容易遗漏边界火点或误纳固定热源。

3.6 从“相交”到“业务候选点”

更完整的秸秆焚烧筛选应增加：

• 距工业企业和电厂的距离；
• 距村庄和垃圾处理设施的距离；
• 周边农田比例；
• 连续多时相是否重复出现；
• 是否处于收获季；
• 是否存在烟羽；
• 是否与高分辨率影像或现场记录一致。

可以构建一个规则分数：

$$Score=w_1{C}_{fire}+w_2{P}_{crop}-w_3{D}_{industry}-w_4{R}_{repeat}$$

其中，$C_{fire}$代表火点置信度，$P_{crop}$代表周边农田比例，$D_{industry}$代表工业热源影响，$R_{repeat}$代表固定热源重复出现特征。该表达式仅表示分析思路，具体权重需要业务数据校准。

3.7 火点结果的时间含义

MODIS火点产品反映卫星过境瞬间或短时间窗口内的热异常。一次未检测到火点，并不代表当天没有焚烧，因为燃烧可能发生在两次过境之间，也可能被云遮挡。

因此，火点统计应注明：

• 过境时间；
• 白天或夜间；
• 云量；
• 观测次数；
• 置信等级；
• 重复火点去重规则。

::: important 火点产品的正确定位
MOD14适合大范围快速筛查，但不能替代全天候连续监控。结果更适合作为执法巡查和现场核查的候选线索。
:::

3.8 MODIS火点的空间分辨率并不等于最小火场面积

火点产品的像元尺寸通常远大于实际燃烧区域，但高温目标能够通过亚像元热异常被检测。材料中提到可识别较小火点，这反映的是高温目标对中红外辐射的强烈贡献，而不是卫星能够直接描绘几十平方米火场的几何边界。

因此，火点产品中的一个像元应解释为“该像元范围内存在热异常证据”，而不是“整个像元都在燃烧”。将火点栅格转为矢量后，通常应保留像元中心点或置信范围，不宜把整个像元多边形直接当作实际燃烧面积。

3.9 固定热源识别可以利用时间重复性

工业烟囱、电厂和油气燃烧设施常在相近位置重复出现热异常，而秸秆焚烧通常具有季节性和短时性。可以基于多年或多日火点记录计算重复出现次数：

$$R_i=\sum _{t=1}^T{I}_{i,t}$$

其中，$I_{i,t}$表示位置$i$在时刻$t$是否出现火点。若某位置在非收获季也频繁出现，且空间位置高度稳定，应优先怀疑固定热源。

这种时间重复性分析能够显著减少仅靠单景农用地叠加造成的误判。

3.10 农用地产品本身也存在分类误差

MCD12Q1是全球土地覆盖产品，其分类体系和空间分辨率决定了它只能提供区域尺度土地覆盖背景。城乡交错区、小块农田和设施农业可能被混合到其他类别中。

因此，秸秆焚烧候选点筛选至少包含两层不确定性：

1. 火点产品是否正确识别热异常；
2. 土地覆盖产品是否正确识别农用地。

较严格的业务流程可进一步叠加更高分辨率耕地数据、永久基本农田边界或当年作物分布图，以提高空间约束精度。

3.11 火点数量不能直接等同于焚烧事件数量

相邻像元可能属于同一场燃烧，同一火场也可能在连续过境中被重复记录。统计前应明确事件聚合规则，例如：

• 空间距离小于一定阈值；
• 时间间隔小于一定小时数；
• 位于同一农田或行政村；
• 火点置信等级相近。

只有完成时空聚类后，才能从“火点像元数”进一步估算“焚烧事件数”。二者属于不同统计口径。

---

04 林冠状态动态监测：识别“同一地类内部的改变”

4.1 林冠变化属于状态改变，而不只是地类转化

土地覆盖变化通常分为两类。第一类是转化，例如森林变为耕地；第二类是改变，例如森林仍然是森林，但冠层密度、叶量、健康或物种组成发生变化。

林冠监测关注的正是第二类问题。虫害、干旱、林火、风倒和衰退都可能使树冠由绿色转为红褐色，或使林隙增加，却未必立即导致“森林”类别消失。

这类变化比土地覆盖转化更难，因为目标差异通常小于森林与非森林之间的差异，对多时相一致性要求更高。

4.2 为什么先提取林区

若直接在整景影像上计算NDVI差值，农田收获、城市施工、水位变化和裸土变化都会产生大幅差异。案例先通过面向对象方法提取林区，再把林区作为掩膜应用于后续变化检测。

面向对象林区提取利用的不只是NDVI，还包括：

• 对象平均光谱；
• 纹理；
• 面积；
• 紧致度；
• 邻接关系；
• 冠层空间连续性。

分割尺度和合并尺度控制对象边界。尺度过小会把单个树冠拆碎，尺度过大则可能将林地与道路、阴影或草地合并。案例中的参数只能视为特定影像上的试验结果，不能直接套用到其他分辨率和季节。

4.3 两期影像必须满足哪些条件

林冠差异只有在以下条件尽量一致时才可解释：

• 成像季节接近；
• 植被物候接近；
• 太阳高度角差异可控；
• 大气影响得到处理；
• 空间配准精确；
• 影像分辨率一致；
• 林区掩膜定义一致。

案例采用QUAC快速大气校正，适合减少明显大气差异。但QUAC依赖场景统计，若两期场景组成差异很大，校正后仍可能存在系统偏差。更稳健的做法是检查稳定地物在两期影像中的反射率是否接近。

4.4 NDVI差值的计算和解释

两期NDVI差值为：

$$\Delta NDVI=NDV{I}_{t2}-NDV{I}_{t1}$$

若$\Delta NDVI$显著为负，说明后期植被绿度或冠层活力低于前期；若接近零，说明变化不明显；若显著为正，可能表示恢复、生长增强或物候差异。

关键问题是“显著”如何定义。不能仅凭颜色表主观判断，而应分析差值分布。若稳定林区差值近似服从某种集中分布，可以使用均值和标准差确定初始阈值：

$$T_{low}=\mu -k\sigma$$

其中，$\mu$为稳定林区差值均值，$\sigma$为标准差，$k$控制异常程度。实际阈值还需结合已知虫害区或人工判读样本调整。

4.5 阈值法为何容易失效

固定阈值可能受到以下因素影响：

• 阴坡与阳坡照度差异；
• 不同树种的季节变化；
• 局部云影；
• 高山阴影；
• 配准误差；
• 林缘混合像元；
• 正常采伐与病虫害混淆。

因此，较好的策略是将差值阈值与其他条件组合，例如：

• 仅在连续林区内部判断；
• 排除阴影和云；
• 结合红光、近红外和红边变化；
• 对持续多个时相下降区域提高置信度；
• 使用林业小班或调查数据验证。

4.6 两期差值向时间序列扩展

两期差值只能告诉我们起点和终点不同，却不能回答变化发生在何时、是否持续、是否恢复。若拥有多年数据，可分析林冠指数时间序列：

$$V{I}_t=Tren{d}_t+Seaso{n}_t+Disturbanc{e}_t+{\epsilon }_t$$

其中，$Tren{d}_t$表示长期趋势，$Seaso{n}_t$表示季节性，$Disturbanc{e}_t$表示扰动，${\epsilon }_t$表示噪声。通过识别突变点，可区分正常季节波动与病虫害、火灾等持续性异常。

4.7 林冠健康结果的业务解释

“健康下降”应被视为遥感异常类别，而不是病因诊断。进一步判定病虫害需要结合：

• 虫情调查；
• 枯死木分布；
• 气象干旱；
• 林火记录；
• 采伐许可；
• 风灾记录；
• 无人机或高分影像。

遥感能够快速确定“哪里值得去看”，但不能仅凭NDVI差值确定生物学病因。

4.8 面向对象林区提取为何有助于降低变化噪声

高分辨率影像中，单棵树冠内部可能同时存在受光面、阴影面和枝叶间隙。如果按像元阈值提取森林，同一林冠容易被分成多个零散区域。面向对象方法先把相邻像元组合成具有空间意义的对象，再计算对象平均光谱和纹理，可降低局部阴影对林区识别的干扰。

对象尺度还决定了变化监测的最小单元。若对象接近林分或树冠尺度，变化结果更容易与林业小班和调查单元对应；若对象过大，则局部虫害会被平均掉；若对象过小，则结果又退化为像元级噪声。

4.9 NDVI饱和会削弱高覆盖森林的变化敏感性

在高叶面积指数和高冠层覆盖条件下，NDVI可能进入饱和区。此时，林冠叶量发生一定下降，NDVI变化仍可能很小。也就是说，NDVI对中低覆盖植被较敏感，但对非常茂密森林的结构变化未必足够敏感。

可补充使用：

• 红边指数；
• 短波红外水分指数；
• 可见光颜色变化；
• 纹理；
• 冠层阴影比例；
• SAR后向散射；
• LiDAR冠层高度。

山松甲虫侵害等过程往往伴随叶片颜色、水分和冠层结构变化，多特征组合通常比单一NDVI差值更稳定。

4.10 变化阈值应区分森林类型和地形条件

同一阈值应用于所有林区，默认不同树种、坡向和海拔具有相同的正常变化范围。现实中，常绿林与落叶林、阴坡与阳坡、低海拔与高海拔的物候和光照差异明显。

可以按林型或地形分区分别估计稳定样本的$\mu$和$\sigma$，再设置分区阈值。这样做的本质是让“正常变化基线”与局地生态条件匹配。

4.11 林冠变化验证应采用空间独立样本

若验证样本来自阈值调参时使用的同一片受灾区，结果会偏乐观。应在不同空间区域选择：

• 已知健康林；
• 已知受灾林；
• 轻度变化区；
• 阴影和林缘混淆区；
• 采伐或火灾区域。

验证指标不仅要报告总体精度，还应报告“健康下降”类别的用户精度和生产者精度，分别回答误报和漏报问题。

---

05 农业耕作变化监测：分类后比较与误差传播

5.1 农业变化为什么强烈依赖物候

农田在播种、出苗、生长、成熟、收获、留茬和休耕阶段呈现完全不同的光谱状态。两个年份的影像即使相差几天，也可能因为农时不同产生明显差异。

因此，农业变化监测首先要比较的是“相同物候阶段”，其次才是“相同日历日期”。例如，一年气温偏高导致提前收获，使用固定日期影像可能把正常农时差异误判为耕作类型变化。

5.2 案例分类体系的业务含义

案例将地表划分为：

• 绿植耕地；
• 留茬耕地；
• 休耕地或裸地；
• 水体；
• 居民用地。

这一体系不是一般土地覆盖分类，而是针对耕作状态设计。绿植耕地反映正在生长的作物，留茬耕地反映收获后的残茬，休耕地或裸地反映无明显植被覆盖的耕地状态。

类别体系必须在两个时相保持一致。若前期把裸地归为休耕地，后期却把相同地表归为建设用地，转移矩阵将混入分类规则差异。

5.3 训练样本不应只选择“典型区域”

只选颜色纯正、边界清晰的样本，会导致模型在真实复杂区域中失效。农业样本应覆盖：

• 不同作物类型；
• 不同长势；
• 不同土壤背景；
• 不同灌溉状态；
• 阴影；
• 田块边缘；
• 留茬密度差异；
• 裸地与建设用地混淆区。

训练集和验证集还应空间分离。若同一田块相邻像元同时进入训练和验证，精度会被高估，因为这些像元具有高度相似的光谱和纹理。

5.4 分类后比较的核心逻辑

分类后比较先对每个时相独立分类，再组合前后类别：

$$T(x,y)=\left(C_{t1}(x,y),C_{t2}(x,y)\right)$$

例如，“绿植耕地→居民用地”表示可能的建设占用，“休耕地→绿植耕地”表示恢复耕作，“绿植耕地→留茬耕地”可能只是收获阶段变化。

5.5 转移矩阵如何阅读

转移矩阵的行表示前期类别，列表示后期类别。对角线表示未变化，非对角线表示转移。

前期\后期	绿植耕地	留茬耕地	休耕地	居民用地
绿植耕地	持续生长	收获或季节差异	休耕或退化	建设占用候选
留茬耕地	复种或再生长	保持	休耕	建设占用候选
休耕地	恢复耕作	恢复后收获	保持	建设占用候选
居民用地	需重点复核	需重点复核	拆除或误分	保持

转移矩阵不能只看净变化。例如，耕地总面积不变，可能同时存在大量转入和转出。净变化会掩盖空间重组过程。

5.6 分类误差如何传递到变化图

若两个时相各自存在误分，变化图会把这些误差组合起来。变化类别的正确率不是简单等于两个总体精度的乘积，但可以明确：分类后比较的误差通常大于单期分类误差。

假设某像元前期真实为绿植耕地，后期真实为留茬耕地。只要任一时相被误分，最终变化类型就可能错误。尤其当两个时期采用不同样本、不同分类器或不同后处理规则时，非真实变化会显著增加。

因此，应分别验证：

1. 前期分类；
2. 后期分类；
3. 变化类别。

5.7 后处理为何不能过度

多数滤波和小斑块删除能够减少椒盐噪声，但农业田块中真实存在狭长地块、田埂、小型设施和碎片化耕作。过度平滑会：

• 合并不同农田；
• 删除小地块；
• 改变边界；
• 低估变化面积。

后处理参数应与最小业务制图单元一致，而不是以“图看起来整齐”为目标。

5.8 从两期分类向农业时序监测扩展

农业生产具有完整生长曲线。多时相数据可以提取：

• 返青时间；
• 植被指数峰值；
• 生长季长度；
• 收获下降时间；
• 复种次数；
• 灌溉脉冲；
• 休耕持续时间。

相比两期分类，时间序列更能区分“作物类型变化”和“同一作物不同生育阶段”。

5.9 不同耕作状态的光谱差异来自哪里

绿植耕地通常在红光中吸收较强、近红外中反射较高；留茬耕地同时包含作物残体和裸土，光谱处于植被与土壤之间；休耕地或裸地主要由土壤亮度、含水量和粗糙度控制；居民用地则表现出建筑材料、道路和阴影的混合特征。

难点在于这些类别并非稳定光谱端元。降雨后的裸地可能明显变暗，浅色屋顶可能与干燥裸土接近，稀疏绿植耕地可能与留茬地混淆。因此，农业分类应尽量加入纹理、田块形状和时间信息。

5.10 分类体系必须服务于变化问题

若研究目标是识别建设占用，分类体系应确保耕地与居民用地可稳定区分；若目标是识别休耕，则需要细分绿植、留茬和裸地状态；若目标是作物类型变化，则需要引入多时相物候。

类别越多并不一定越好。类别过细会增加样本需求和类间混淆，最终可能使变化矩阵充满不可靠的小类别转移。合理做法是先根据业务问题定义最小必要类别，再评估影像是否具备区分这些类别的能力。

5.11 变化面积应考虑分类误差

直接统计变化图像元数量会把分类误差当作真实面积。较严谨的做法是利用误差矩阵进行面积校正，或至少为变化面积报告置信区间。

例如，某类“耕地转居民用地”在样本中用户精度较低，说明图中标出的很多变化可能是假变化。此时，未经校正的像元面积不能直接作为正式统计结果。

5.12 地块级分析通常比像元级结果更接近农业管理

农业管理对象通常是田块。将像元分类结果叠加地块边界后，可以按多数比例或面积比例确定地块状态。

设地块$p$中类别$c$的像元面积为$A_{p,c}$，地块总面积为$A_p$，则类别比例为：

$$P_{p,c}=\frac{A_{p,c}}{A_p}$$

当某类别比例超过设定阈值时，可将地块判为该状态。地块级汇总能够减少零散噪声，也便于与承包地、作物调查和农业统计数据关联。

---

06 Crosta矿物蚀变提取：从波段选择到异常解释

6.1 围岩蚀变为何能够成为遥感找矿线索

成矿热液与围岩发生反应后，会形成新的矿物组合。部分羟基矿物、黏土矿物、碳酸盐矿物和铁氧化物在可见光、近红外和短波红外范围具有特征吸收。

遥感不能直接看到矿床本身，但可以识别与矿化活动相关的蚀变信息，从而缩小野外调查范围。其业务定位是“圈定找矿靶区”，而不是直接证明地下存在矿体。

6.2 Crosta方法不是普通PCA

普通PCA通常把所有波段同时输入，目标是最大化总体方差；Crosta方法则先根据目标矿物的光谱机理选择波段，再对目标波段子集进行PCA，以提高目标信息在某个主成分中的集中程度。

案例选择TM1、TM4、TM5和TM7，针对羟基和碳酸根相关蚀变进行分析。其核心不是“第几个主成分固定对应蚀变”，而是读取特征向量，判断波段系数之间的符号和大小关系。

6.3 为什么必须先去除植被和水体

植被在红光和近红外之间具有强烈反差，水体在近红外和短波红外中具有低反射。若不进行掩膜，这些大幅度变化会主导协方差结构，使PCA优先解释植被覆盖和水陆差异，而不是较弱的矿物蚀变差异。

因此，预处理通常包括：

• 辐射定标；
• FLAASH大气校正；
• 工程区裁剪；
• 植被掩膜；
• 水体掩膜；
• 云和阴影排除。

背景掩膜不是为了“让结果更漂亮”，而是为了改变统计总体，使PCA的方差结构更接近地质目标。

6.4 PCA的数学含义

设输入波段向量为$\mathbf{x}$，均值为$\mathbf{\mu }$，协方差矩阵为$\mathbf{\Sigma }$。PCA求解：

$$\mathbf{\Sigma }{\mathbf{a}}_k={\lambda }_k{\mathbf{a}}_k$$

其中，${\lambda }_k$为第$k$个特征值，${\mathbf{a}}_k$为对应特征向量。主成分为：

$$P{C}_k={\mathbf{a}}_k^T(\mathbf{x}-\mathbf{\mu })$$

特征值越大，说明该主成分解释的总体方差越多。但蚀变信号可能只占较小区域，因此未必出现在第一主成分。目标成分常位于较后的主成分中，只要其波段系数组合符合目标吸收特征。

6.5 如何利用特征向量判断目标成分

案例给出的判断原则是：TM5系数与TM7、TM4系数符号相反，TM1通常与TM5符号相同。这种相反符号意味着主成分在强化一个波段的同时抑制另一个波段，从而突出目标波段间的反差。

假设某主成分为：

$$PC=a_{1}TM1+a_{4}TM4+a_{5}TM5+a_{7}TM7$$

若$a_5$与$a_7$符号相反，且绝对值较大，说明该成分强调TM5与TM7之间的差异。是否对应蚀变，还需结合目标矿物光谱机理、异常区亮暗和已知地质信息。

需要注意：特征向量整体乘以$-1$仍然表示同一个主成分。因此，“亮异常”与“暗异常”可能因软件输出符号而反转，解释时应看相对符号而不是固定亮暗。

6.6 标准差异常分级

案例采用标准差倍数划分异常等级。设目标主成分均值为$\mu$，标准差为$\sigma$，高值异常阈值为：

$$T_k=\mu +k\sigma$$

常见分级形式为：

• 一级异常：大于$\mu +3\sigma$；
• 二级异常：$\mu +2.5\sigma$至$\mu +3\sigma$；
• 三级异常：$\mu +2\sigma$至$\mu +2.5\sigma$。

这种方法假设背景值集中分布，异常位于分布尾部。若研究区岩性差异很大或主成分分布明显偏斜，单一标准差阈值可能把正常岩性边界识别为异常。

6.7 为什么异常必须地质验证

遥感异常可能由多种因素造成：

• 蚀变矿物；
• 裸岩差异；
• 盐碱地；
• 土壤含水量；
• 地形阴影；
• 建筑材料；
• 大气残差；
• 传感器噪声。

验证可结合：

• 已知矿点；
• 地质图；
• 岩性图；
• 构造线；
• 地球化学异常；
• 野外光谱；
• 岩石样品；
• 高光谱数据。

当遥感异常与构造、岩性接触带、地球化学异常和已知矿点同时吻合时，其找矿意义才显著提高。

::: warning 遥感蚀变异常的正确表述
应表述为“蚀变信息异常区”或“找矿候选靶区”，不能直接表述为“矿体分布区”。
:::

6.8 波段组合的选择应回到矿物吸收机理

Crosta方法的波段组合并非任意。目标矿物在某些波段表现为相对高反射，在另一些波段表现为吸收低值。选择波段时，应确保组合能够同时覆盖目标吸收与背景反射差异。

对于宽波段Landsat TM，单个波段往往包含多个矿物特征，因此提取结果通常是“蚀变矿物组合异常”，而不是某一种纯矿物的精确识别。若需要区分高岭石、绢云母、绿泥石等具体矿物，通常需要更高光谱分辨率数据和光谱匹配方法。

6.9 协方差PCA与相关系数PCA的差别

若不同波段数值范围差异很大，使用协方差矩阵PCA时，方差较大的波段会占据更高权重；使用相关系数矩阵相当于先对各波段标准化，使每个波段具有相同方差。

两种方法可能得到不同特征向量。Crosta分析时应记录PCA采用的统计方式，不能只记录输入波段和主成分序号。

6.10 主成分异常的正负反转如何处理

特征向量$\mathbf{a}$和$-\mathbf{a}$都是同一特征值的合法解，因此软件可能输出正负相反的主成分。若目标信息表现为暗异常，可以将主成分乘以$-1$后再进行统一分级：

$$P{C}_{rev}=-PC$$

这种反转只改变显示方向，不改变空间异常结构。真正的判别依据是目标波段系数之间的相对符号关系。

6.11 标准差分级前应检查分布形态

均值加标准差方法隐含了背景分布较集中、异常位于尾部的假设。若主成分直方图具有明显偏态、多峰或长尾，$\mu +3\sigma$可能过严或过松。

可以结合：

• 分位数；
• 自然断点；
• 局部背景统计；
• 已知矿点响应；
• ROC分析；
• 分区阈值。

在岩性复杂地区，按岩性单元或地貌单元分别估计背景值，通常比全区统一阈值更合理。

6.12 从像元异常到找矿靶区需要空间综合

零散异常像元通常不具有直接找矿意义。后处理可包括：

• 连通域分析；
• 最小面积筛选；
• 异常等级合并；
• 沿构造线缓冲分析；
• 与岩体接触带叠加；
• 与地球化学异常叠加；
• 形成靶区优先级。

靶区评价可概念化表示为：

$$Priority=w_1{A}_{alteration}+w_2{D}_{fault}+w_3{G}_{geochem}+w_4{K}_{known}$$

其中，各项分别表示蚀变异常、断裂邻近性、地球化学异常和已知矿点支持。该式表达综合判断思路，实际权重需由地质专家和样本数据确定。

---

07 不同方法的误差结构为何不能混为一谈

7.1 连续反演、离散分类和异常检测的验证对象不同

水质反演输出连续浓度，应使用RMSE、MAE、Bias和残差分析；农业分类和火点筛选输出离散类别，应使用混淆矩阵、用户精度、生产者精度和F1；Crosta方法输出统计异常，其验证更接近目标发现问题，应关注已知蚀变带命中率、背景误报率和空间合理性。

若只用一类指标评价所有任务，会掩盖方法真正的问题。例如，水质模型的$R^2$不能替代浓度误差，火点数量与现场记录的相关性也不能替代单点正确率。

7.2 用户精度与生产者精度对应不同业务风险

对于秸秆焚烧监测，“用户精度低”意味着候选点中误报较多，会增加巡查成本；“生产者精度低”意味着真实焚烧被漏掉，会削弱监管覆盖。

对于矿物蚀变靶区，误报意味着增加野外验证成本，漏报则可能错过潜在找矿线索。不同业务对两类错误的容忍度不同，因此阈值不能仅以总体精度最大为目标。

7.3 阈值本质上是在误报和漏报之间取舍

无论是NDVI变化阈值、火点置信等级还是标准差异常倍数，阈值越严格，通常误报越少、漏报越多；阈值越宽松，候选区增多，但现场核查负担也增加。

因此，阈值应根据任务场景选择：

• 科研分析强调误差可解释和独立验证；
• 应急监测可能优先降低漏报；
• 执法巡查需要控制误报成本；
• 找矿初筛可接受较宽松候选区；
• 正式成果统计需要更严格的精度保证。

---

08 ENVI实施过程中最容易被忽略的技术细节

7.1 数据读取必须保留产品语义

HJ-1B、MODIS和Landsat数据都包含专用元数据。使用传感器专用读取器的意义在于自动获取：

• 波段顺序；
• 定标系数；
• 波长；
• 投影；
• 缩放因子；
• NoData；
• 产品编码含义。

把HDF中的分类层直接当普通灰度影像打开，或忽略角度层缩放因子，都会使后续分析失去物理意义。

7.2 连续数据和离散数据不能使用同一种重采样方法

反射率、NDVI和温度属于连续变量，可使用双线性或三次卷积；土地覆盖、火点等级和分类结果属于离散变量，应优先使用最近邻。

若对火点等级使用双线性插值，原本的7、8、9会产生7.4、8.2等无定义数值，产品编码将被破坏。

7.3 掩膜应用前后要检查NoData

掩膜后的背景可能被赋值为0，也可能被赋值为NoData。若目标参数本身允许0值，使用0作为背景会混淆“真实零值”和“无数据”。

更合理的做法是：

• 使用明确NoData值；
• 在统计时排除NoData；
• 输出前检查直方图；
• 在颜色分级中单独设置背景。

7.4 波段运算应明确数据类型

比值、回归和指数计算应输出浮点型。若输入为整数且软件按整数运算，可能发生截断。

例如，计算近红外与红光比值时应使用浮点转换：

float(B_NIR) / B_RED

同理，叶绿素反演、NDVI差值和PCA异常指数都应检查输出精度和缩放因子。

7.5 模型参数和阈值应保存为工程元数据

一个可复现的专题工程至少应记录：

内容	示例
输入数据	文件名、产品级别、日期
软件环境	ENVI版本、扩展模块
预处理	定标、大气、配准、重采样
模型	公式、系数、样本数
阈值	火点编码、NDVI阈值、标准差倍数
掩膜	水体、林区、农用地、植被
验证	样本来源、指标、结果
输出	坐标系、分辨率、NoData、数据类型

---

09 五类专题的质量控制与复核路径

8.1 水质反演检查链

• 影像与实测点是否同日或接近卫星过境时间？
• 波段是否完成辐射定标和大气校正？
• 水岸混合像元是否剔除？
• 模型是否保留独立验证样本？
• 是否检查高浓度和低浓度区的残差？
• 输出浓度是否存在负值或异常极值？

8.2 火点监测检查链

• MOD14与MOD03是否来自同一过境时刻？
• GLT地理定位是否与海岸线、道路或行政边界吻合？
• 火点编码是否依据对应产品说明解释？
• 土地覆盖数据是否完成正确投影和镶嵌？
• 是否排除工业固定热源和重复异常点？
• 是否注明过境时间、云量和火点置信度？

8.3 林冠变化检查链

• 两期影像物候是否一致？
• 林区掩膜是否在两期采用一致规则？
• 配准误差是否小于目标变化尺度？
• 稳定林区差值是否接近零？
• 阈值是否由样本或统计分布支持？
• 健康下降是否有林业调查或高分影像验证？

8.4 农业变化检查链

• 两期分类体系和编码是否一致？
• 训练与验证样本是否空间独立？
• 影像日期是否对应相近农时？
• 两期分类是否分别验证？
• 变化图是否独立抽样验证？
• 后处理是否保留小型真实田块？

8.5 蚀变异常检查链

• 植被、水体、云和阴影是否掩膜？
• 输入波段是否与目标矿物机理一致？
• 是否查看特征向量而不是只看主成分序号？
• 是否考虑主成分正负号整体反转？
• 异常阈值是否结合区域背景校准？
• 是否与地质、构造、矿点和地球化学数据叠加验证？

---

10 常见失败情形及其原因

9.1 水质图出现大面积负值

可能原因包括大气校正过度、模型外推、波段顺序错误、缩放因子遗漏或水体掩膜不完整。负值不应简单全部改为0后结束，而应先确认错误来自模型还是数据。

9.2 火点整体偏移

通常与MOD14和MOD03时间不匹配、经纬度数组读取错误、GLT输出投影设置错误或数据行列不一致有关。应先用海岸线和大型地理目标检查地理定位，再进行火点提取。

9.3 林冠变化沿道路和林缘呈现细线

这类“边缘变化”通常是配准误差造成的。道路、林缘和树冠边界的空间梯度很大，即使小于一个像元的错位，也会产生正负成对差值。

9.4 农业变化面积异常巨大

可能是两期影像物候差异、分类规则不一致、训练样本偏差或后处理尺度不同造成的。应先比较两期类别直方图，再检查转移矩阵中是否出现不合理大规模转移。

9.5 蚀变异常与阴坡高度重合

说明主成分可能仍受地形照度影响。可以检查DEM、坡向和阴影掩膜，或进一步进行地形校正。若异常只出现在同一坡向而与地质构造无关，其找矿意义通常有限。

---

11 术语表

术语	说明
叶绿素a	表征浮游植物和水体营养状态的重要色素参数
星地匹配	将地面采样与卫星像元按时间和空间对应
波段比值	利用两个波段之比构建遥感特征
MOD14	MODIS热异常和火点相关产品
MOD03	MODIS地理定位和观测几何产品
MCD12Q1	MODIS年度土地覆盖产品
GLT	Geometry Lookup Table，几何查找表
林冠状态	冠层叶量、密度、健康和林隙等状态
改变	同一土地覆盖类型内部结构或功能变化
转化	不同土地覆盖类型之间的变化
分类后比较	分别分类两个时相后组合类别进行变化分析
转移矩阵	记录前后类别转换数量的矩阵
Crosta方法	基于目标波段组合和PCA的蚀变异常提取方法
特征向量	PCA中各输入波段对主成分的权重
标准差异常	根据均值和标准差倍数提取分布尾部异常
羟基蚀变	与含羟基矿物相关的遥感蚀变信息
地表反射率	尽量消除大气影响后的地表反射特性
NoData	表示无有效观测值的像元编码
RMSE	对大误差较敏感的均方根误差
Bias	反映系统性高估或低估的平均偏差

---