1. 为什么说传统采样方式已经“过时”了如果你做过大范围的遥感土地利用分类尤其是那种需要分析好几年、甚至十几年变化的研究我猜你一定对“选样本点”这个步骤又爱又恨。爱的是样本选得好分类精度就高论文数据就漂亮恨的是这个过程实在是太折磨人了。我最早做项目的时候就是对着高分辨率影像一个像素一个像素地手动勾画这块是林地那块是耕地那片是水体……眼睛都快看花了一天下来也画不了多少。更头疼的是当你需要做时序分析比如分析过去十年某个区域的森林变化你需要为每年的影像都准备一套样本。先不说工作量是指数级增长光是保证不同年份样本位置和类别解释的一致性就足以让人崩溃。去年你觉得这块地是“旱地”今年再看可能因为作物轮作或者影像质量你又觉得它像“草地”了。这种主观性带来的误差最后都会直接反映在你的分类结果里让整个研究的可靠性大打折扣。所以我们一直在寻找一种更聪明、更“懒”的办法。核心思路就是能不能让机器帮我们自动、客观地生成样本这个想法在Google Earth EngineGEE这个强大的云端平台成为主流后变得前所未有的可行。因为GEE里集成了越来越多全球性的、持续更新的权威土地覆盖数据产品比如我们今天要重点用到的Google Dynamic World。这些产品本身就是由顶尖团队利用海量数据和先进算法生产出来的我们可以把它们看作是“标准答案库”或者“超级老师”。我们的新思路就是“借力”这位老师让它来告诉我们在某个位置、某个时间点最可能的地类是什么然后我们基于这些信息自动化、随机化地“抄作业”——生成成千上万个训练和验证样本。这彻底把我们从繁重、主观的目视解译中解放了出来把精力更多地投入到模型构建和结果分析上。下面我就手把手带你走通这套高效、可复现的样本构建与分类全流程。2. 认识我们的“权威标签源”Google Dynamic World在开始动手写代码之前我们得先了解清楚我们借力的“权威数据”到底靠不靠谱。这里我首推Google Dynamic World (DW)。为什么是它而不是MODIS土地覆盖或者ESA WorldCover呢我对比使用过不少产品DW有几个让我觉得特别顺手的特点。首先它的时空分辨率非常友好。DW是基于Sentinel-2影像生成的这意味着它提供了近乎全球覆盖的10米分辨率数据而且更新频率很高。相比之下很多30米或500米分辨率的产品在分析地块破碎的农田、城市细节时就不够用了。对于做精细分类10米是一个很实用的尺度。其次它是近实时Near Real-Time的。DW数据产出很快这意味著你可以获取到最近几天、几周的地表覆盖情况。这对于监测快速变化比如火灾、洪水后的土地覆盖变化或者作物生长季的动态非常有价值。你不再需要苦等年度产品发布。最关键的是它的分类体系很实用。DW将每个像素的概率预测为9个类别水体、树木、草地、淹没植被、作物、灌木/灌丛、建设用地、裸地、冰雪/云。这个分类体系基本覆盖了全球主要的地表类型而且“淹没植被”这个类别对于湿地研究特别有用。它输出的不是一个简单的“硬分类”结果而是每个类别的概率0-1之间这给了我们很大的灵活性。我们可以直接使用概率最高的类别作为标签也可以设定一个概率阈值比如0.7来获取更高置信度的样本这对于提升样本质量很有帮助。在GEE中调用它非常简单它的数据资产ID是GOOGLE/DYNAMICWORLD/V1。你可以把它想象成一个巨大的、覆盖全球、持续更新的“地类字典”。我们接下来的所有自动化操作都基于从这个“字典”里查询信息。当然没有任何一个产品是完美的DW在复杂山区或者同物异谱/同谱异物现象严重的区域也可能有误判。但作为一个自动化、大批量生成样本的“起点”它的客观性、一致性和易用性已经远远超过了传统人工采样。我们的策略是用它生成基础样本再结合少量的本地验证数据进行校准这比完全从零开始要高效和稳健得多。3. 搭建自动化样本工厂从数据准备到样本生成理论说清楚了咱们直接上干货。我会把整个流程拆解成几个关键步骤并附上详细的代码和解释。你完全可以把这段代码复制到GEE的代码编辑器里替换掉你的研究区就能直接跑起来。3.1 第一步划定战场并获取干净的影像做任何遥感分析第一步都是确定你的研究区ROI和获取预处理好的影像。这里我们以Landsat 8地表反射率数据为例因为它时间序列长且经过了大气校正。// 1. 定义研究区 - 这里以中国某个区域为例你可以通过上传矢量文件或画图工具定义你的roi var roi ee.FeatureCollection(你的研究区资产路径); // 或者用 ee.Geometry.Rectangle([经度, 纬度, 经度, 纬度]) // 2. 定义Landsat 8影像去云和缩放函数 function maskL8sr(image) { // 利用QA_PIXEL波段进行云、云阴影、卷云掩膜 var qaMask image.select(QA_PIXEL).bitwiseAnd(parseInt(11111, 2)).eq(0); // 利用QA_RADSAT波段去除饱和像素 var saturationMask image.select(QA_RADSAT).eq(0); // 对光学波段和热红外波段分别进行缩放转换 var opticalBands image.select(SR_B.).multiply(0.0000275).add(-0.2); var thermalBands image.select(ST_B.*).multiply(0.00341802).add(149.0); // 替换原始波段并应用掩膜 return image.addBands(opticalBands, null, true) .addBands(thermalBands, null, true) .updateMask(qaMask) .updateMask(saturationMask); } // 3. 获取并处理指定时间和区域的影像合成一幅镶嵌图 var img ee.ImageCollection(LANDSAT/LC08/C02/T1_L2) .filterDate(2024-05-01, 2024-10-30) // 选择时间范围例如一个生长季 .filterBounds(roi) .filter(ee.Filter.lte(CLOUD_COVER, 10)) // 过滤掉云量高于10%的影像 .map(maskL8sr) // 应用去云函数 .median() // 使用中位数合成能有效抑制异常值如残余云、噪声 .clip(roi); // 4. 可视化一下看看影像质量 Map.centerObject(roi, 9); // 地图中心定位到研究区缩放级别9 var visualization { min: 0.0, max: 0.3, bands: [SR_B4, SR_B3, SR_B2] }; // 标准假彩色 Map.addLayer(img, visualization, 2024 Landsat 8 RGB);这里我用了.median()中位数合成而不是原文中的.mosaic()。在实际应用中中位数合成对于去除时序中的噪声比如没完全掩膜掉的薄云效果更好能得到更干净、更具代表性的光谱特征。这是个小技巧实测下来比简单镶嵌更稳。3.2 第二步请出“权威老师”获取土地覆盖标签现在我们请出核心角色——Dynamic World数据来为我们的影像提供“标准答案”。// 1. 获取与研究区、时间相匹配的Dynamic World数据 // 注意DW数据更新快建议选取与影像期相近的日期以减少因地表变化导致的标签错误 var dw ee.ImageCollection(GOOGLE/DYNAMICWORLD/V1) .filterDate(2024-07-01, 2024-09-30) // 选取夏季稳定时期的数据 .filterBounds(roi) .select(label) // 选择‘label’波段这是概率最大的类别硬分类 .mode() // 对时间序列取众数得到期间最稳定的地类 .clip(roi); print(Dynamic World 标签图像, dw); // 2. 可选但推荐如果你想要更高质量的样本可以使用概率阈值进行过滤 // 先获取最大概率波段和对应的标签 var dw_with_prob ee.ImageCollection(GOOGLE/DYNAMICWORLD/V1) .filterDate(2024-07-01, 2024-09-30) .filterBounds(roi) .max(); // 取时间序列内每个像素各类别的最大概率 var max_prob dw_with_prob.select(water,trees,grass,flooded_vegetation,crops, shrub_and_scrub,built,bare,snow_and_ice).max(); // 计算最大概率值 var dw_label dw_with_prob.select(label); // 对应的类别标签 // 只保留置信度高于0.7的像素作为有效标签区域 var dw_high_confidence dw_label.updateMask(max_prob.gt(0.7)); // 你可以选择使用 dw_high_confidence 代替下面的 dw这里我引入了两个改进点一是使用.mode()众数来合成时间序列的DW数据这能确保我们得到的是研究时段内最稳定的土地覆盖类型避免了单一时相可能的偶然误差。二是提供了一个概率阈值过滤的选项。DW的‘label’波段虽然是硬分类但其背后是有每个类别的概率值的。通过只保留最大概率大于0.7这个阈值你可以调整的像素我们能得到一份置信度更高的“黄金标准”标签用这个来生成样本理论上质量会更好。你可以根据自己研究区的情况决定是否启用这一步。3.3 第三步自动化随机采样构建样本库这是最核心的一步我们将影像特征光谱、指数等和土地覆盖标签结合起来进行分层随机采样。// 1. 将土地覆盖标签作为新增波段添加到影像中 var img_with_label img.addBands(dw.rename(landcover)); // 将DW标签波段命名为‘landcover’ // 2. 定义分类体系对应的数值DW的label默认是0-8 var classValues [0, 1, 2, 3, 4, 5, 6, 7, 8]; // 对应DW的9个类别 // 3. 执行分层随机采样 // 关键参数 // - numPoints: 总样本点数或每类点数取决于‘classPoints’参数 // - classBand: 包含类别信息的波段名 // - region: 采样区域 // - scale: 采样尺度米应与影像分辨率匹配或略大 // - geometries: 是否保留样本点的地理坐标 var sample_points img_with_label.stratifiedSample({ numPoints: 100, // 这里指**每个类别**采集100个点总样本数 100 * 类别数 classBand: landcover, region: roi, scale: 30, // 与Landsat分辨率一致 geometries: true // 设为true后面才能导出或可视化 }); print(生成的样本点集合, sample_points); print(样本数量, sample_points.size()); // 4. 将样本集随机分割为训练集80%和验证集20% var sample_with_random sample_points.randomColumn(); // 给每个样本添加一个0-1的随机数列 var training_sample sample_with_random.filter(ee.Filter.lte(random, 0.8)); var validation_sample sample_with_random.filter(ee.Filter.gt(random, 0.8)); print(训练样本数, training_sample.size()); print(验证样本数, validation_sample.size()); // 5. 可视化查看样本分布 Map.addLayer(training_sample, {color: FF0000}, 训练样本 (红)); Map.addLayer(validation_sample, {color: 0000FF}, 验证样本 (蓝));踩过的一个坑stratifiedSample里的numPoints参数。如果你直接写numPoints: 500GEE可能会尝试在所有类别间总共采集500个点但某些小面积类别可能根本采不到。更稳妥的做法是使用numPoints: 100并配合classPoints参数本例中未显式使用但默认行为是每类点数这能确保每个地类都有足够的样本对于样本不平衡问题比如研究区内水体很少特别重要。你可以通过调整这个数字来控制总样本量。4. 训练与评估让随机森林模型跑起来样本准备好了我们就可以训练分类器了。随机森林Random Forest在遥感分类里是出了名的“万金油”对参数不敏感抗过拟合能力强效果通常不错。// 1. 定义用于分类的特征波段即影像的所有反射率波段 var input_bands img.bandNames(); // 获取影像的所有波段名 // 2. 训练随机森林分类器 var classifier ee.Classifier.smileRandomForest({ numberOfTrees: 100, // 树的数量通常50-100足够更多可能提升有限但计算更慢 variablesPerSplit: null, // 每次分裂随机选择的变量数null表示sqrt(总特征数) minLeafPopulation: 1, // 叶节点最小样本数对于大样本可设为5或10以防过拟合 bagFraction: 0.5, // 每棵树训练时使用的样本比例 seed: 0 // 随机种子设为固定值可使结果可重复 }).train({ features: training_sample, // 训练样本 classProperty: landcover, // 样本中的类别属性字段 inputProperties: input_bands // 用于训练的特征 }); // 3. 在验证集上评估分类器性能 var validation_assessed validation_sample.classify(classifier); var confusion_matrix validation_assessed.errorMatrix(landcover, classification); // 计算混淆矩阵 // 4. 打印一系列精度指标 print(混淆矩阵行真实值列预测值:, confusion_matrix); print(总体分类精度 (OA):, confusion_matrix.accuracy()); print(Kappa系数:, confusion_matrix.kappa()); // 5. 计算并打印每类的生产者精度漏分误差和用户精度错分误差 var consumers_accuracy confusion_matrix.consumersAccuracy(); // 生产者精度 var producers_accuracy confusion_matrix.producersAccuracy(); // 用户精度 print(生产者精度 (PA - 各类别被正确分类的比例):, producers_accuracy); print(用户精度 (UA - 分类结果中各类别的正确比例):, consumers_accuracy); // 6. 将整个研究区的影像进行分类 var classified_image img.classify(classifier).clip(roi);这里我稍微调整了随机森林的参数。numberOfTrees设为100是一个比较稳健的选择。minLeafPopulation我习惯设为1但如果你发现模型在训练集上精度接近100%而在验证集上差很多过拟合可以尝试把它调到5或10。bagFraction是自助采样率0.5是一个常用值意味着每棵树只用50%的样本进行训练这有助于增加树之间的差异性提升模型泛化能力。精度评估部分除了总体精度和Kappa一定要看生产者精度PA和用户精度UA。PA低说明该类地物很多被分到了其他类漏分UA低说明分到该类的像元很多不是这一类错分。这能帮你精准定位分类模型在哪些地类上表现不佳。5. 结果可视化、导出与深度优化思路模型训练评估完我们就可以看看成果并把它保存下来。// 1. 定义可视化配色与Dynamic World官方配色一致 var dw_palette [ #419bdf, // 0: 水 #397d49, // 1: 树 #7a87c6, // 2: 草 #e49635, // 3: 淹没植被 #dfc35a, // 4: 作物 #c4281b, // 5: 灌木/灌丛 #a59b8f, // 6: 建设用地 #b39fe1, // 7: 裸地 #ffffff // 8: 雪/冰 ]; // 2. 将Dynamic World标签和我们的分类结果添加到地图进行对比 Map.addLayer(dw, {min:0, max:8, palette: dw_palette}, Dynamic World 参考); Map.addLayer(classified_image, {min:0, max:8, palette: dw_palette}, 随机森林分类结果); // 3. 导出分类结果到Google云盘 Export.image.toDrive({ image: classified_image, description: LandUse_Classification_2024_RF, folder: GEE_Exports, // 你云盘里的文件夹名 fileNamePrefix: RF_Result_2024, region: roi, scale: 30, crs: EPSG:4326, // 坐标系 maxPixels: 1e13 // 允许导出大量像素 }); // 4. 可选导出样本点用于其他软件或后续分析 Export.table.toDrive({ collection: training_sample, // 也可以导出验证样本或全部样本 description: Training_Samples_2024, folder: GEE_Exports, fileNamePrefix: Training_Points, fileFormat: CSV });到这一步一个完整的自动化分类流程就走通了。但作为实战经验分享我觉得还有几个可以深度优化的点能让你的结果更上一层楼。第一特征工程。我们目前只用了原始的光谱波段。实际上加入一些光谱指数会极大提升分类能力特别是对于区分植被类型。你可以在合成影像img之后计算NDVI归一化植被指数、NDWI归一化水体指数、NDBI归一化建筑指数等作为新增波段加入。在GEE里计算NDVI只需要一行代码var ndvi img.normalizedDifference([SR_B5, SR_B4]).rename(NDVI);然后把ndvi波段加到img里。我试过加入NDVI和NDWI后水体、植被大类的区分度明显提升。第二后处理。直接分类的结果往往会有很多“椒盐噪声”特别是一些散点。我们可以使用简单的众数滤波来平滑结果var smoothed classified_image.focal_mode({radius: 3, units: pixels});。这能去除很多孤立的小斑块让图面更干净。不过要注意滤波也会抹掉一些真实的细节半径要根据你的应用尺度来调整。第三处理类别不平衡。如果你的研究区里90%是森林只有1%是水体那么即使模型把所有像素都分成森林也能达到90%的总体精度但这显然不是我们想要的。我们之前在采样时用分层采样确保了每类有固定数量的样本这在一定程度上缓解了问题。但在训练时可以进一步使用classWeight参数给少数类别更高的权重或者使用smileRandomForest的classProperties参数进行更精细的控制。第四时序特征的威力。本文示例用了单一时相的影像。但GEE最大的优势就是处理时序数据。你可以计算某个时间段内比如一年的NDVI最大值、最小值、平均值、振幅等这些物候特征对于区分作物类型如水稻、小麦、玉米、常绿林和落叶林等有奇效。思路是先构建一个时间序列影像集合然后通过.reduce(ee.Reducer)来生成这些统计特征影像最后将它们一起作为分类特征。这相当于让模型不仅看“长相”还看“生长节奏”分类精度通常会有一个质的飞跃。这套“借力权威数据自动化构建样本”的思路我把它应用到了好几个省级尺度的土地覆盖制图项目中效率提升是肉眼可见的。以前一个熟练的研究生需要一两周才能完成的样本选取工作现在喝杯咖啡的功夫代码就跑完了而且样本的空间分布更均匀完全避免了人为选择样本时无意识的空间偏好。更重要的是这套流程是完全可复现的。你只要保存好代码明年用同样的流程跑一遍就能得到一套时间上可比的结果这对于做变化检测研究来说价值巨大。当然它也不是银弹在权威数据本身质量不高的区域或者地类非常复杂混淆的地区可能还需要引入少量的人工修正样本。但无论如何它已经为我们解决了80%以上的问题让我们能把宝贵的精力投入到那20%更需要人类智能的环节中去。