GEE KMeans非监督分类完整链路与结果解译

发布时间:2026/10/27 22:53:24
GEE KMeans非监督分类完整链路与结果解译

跑通 KMeans 的代码只要十分钟,看懂它输出的那张彩色图却可能要花上一整天——这是我带过的新人最常卡住的地方。地图上摊开六种颜色的斑块,每个斑块下面挂着一个从 0 到 5 的数字,然后呢?哪个是水体,哪个是建筑,哪个是刚收割完的耕地?GEE 的非监督分类(Unsupervised Classification)恰恰是把"分"这件事做得又快又漂亮,却把"认"这件事完整地甩回给了使用者。它不像监督分类那样给你一套现成的类别标签和精度报告,它只负责按光谱相似度把像素抱成团,剩下的解释工作全靠你自己。我见过太多人在这里掉头去用监督分类,其实只是没搞明白聚类结果该怎么读。

这篇内容讲的就是 GEE 平台上非监督分类的完整链路:从影像怎么预处理成聚类器能吃的格式,到ee.Clusterer.wekaKMeans的参数怎么调,再到结果里那一堆数字标签怎么还原成真实地物。适合已经能跑通 GEE 基础代码、做过简单指数计算、但对cluster()和train()这条线还没摸透的人。全文基于 Landsat 8/9 Collection 2 地表反射率产品演示,同样的思路换成 Sentinel-2 完全通用。里面会写到我实际踩过的几个坑,包括量纲不统一导致聚类结果完全失效、采样数量设了五千实际只拿到八百这种,都是文档里不会强调但一定会遇到的东西。

1. 先搞清楚 GEE 里的聚类到底在做什么

1.1 聚类不是分类,这个区别决定了后面所有操作

很多人把这两个词混着用,结果在参数设置上反复走弯路。分类(Classification)的前提是你事先知道有哪些类别,比如"水体、植被、建筑、裸地"这四类,然后拿标注好的样本去教模型认识它们,模型学的是"什么样的光谱组合对应哪一类"。聚类(Clustering)完全反过来,你事先不知道有几类、也不知道是什么类,算法只做一件事:把光谱特征相近的像素归到一组,让组内差异尽可能小、组间差异尽可能大。

这个差别带来的直接后果是:聚类输出的标签是纯粹的编号,0、1、2除了表示"它们不是同一组"之外没有任何含义。编号为0的簇在你的研究区可能是水体,换一个研究区可能就变成了裸地,甚至同一景影像换一次运行顺序都可能重排。所以任何试图把cluster == 2直接当成"植被"写进后续代码的做法,都是在给自己埋雷。

我在实际项目里更愿意把非监督分类当成一个探索工具而不是制图工具。它的强项是:你对研究区一无所知的时候,快速看出这里大概有几种光谱上可区分的地物;或者你手上没有可靠的训练样本,想先摸清楚各类地物的光谱分布范围,为后面做监督分类准备样本。把它当成最终成果直接出图,那是对它的误用。

1.2 GEE 把聚类器包装成了一个"可以 train 的对象"

GEE 里的聚类器都在ee.Clusterer这个命名空间下面,本质是把 Weka 的几个经典聚类算法搬到了云端。常用的几个我列一下各自的适用场景,方便你选:

聚类器核心机制是否需要指定类别数适用场景
wekaKMeansK-均值,迭代找质心需要,nClusters最常用,地物光谱差异明显时效果好
wekaCascadeKMeans级联 K-均值,自动搜索 k不需要,给范围即可不确定该分几类时先跑它探路
wekaXMeans扩展 K-均值给上下界想比 KMeans 更灵活一点,但计算量更大
wekaLVQ学习向量量化需要有先验的初始中心时可以用
wekaCobweb概念聚类不需要输出类别数不确定,实用性一般

它们的调用方式高度统一,都是先.train()再.cluster():

var clusterer = ee.Clusterer.wekaKMeans({ nClusters: 6 }).train(trainingFeatures); var result = image.cluster(clusterer);

.train()吃的是一个ee.FeatureCollection,返回的是一个训练好的Clusterer对象;.cluster()吃一张Image,返回一张单波段、整型、波段名为cluster的影像。这两步之间没有任何魔法,train()做的是在多维特征空间里找质心,.cluster()做的是把每个像素往最近的质心上贴。理解了这个,你就知道为什么训练样本的代表性比数量更重要——质心的位置完全由样本决定,样本偏了,质心就偏了,整张图就废了。

1.3 为什么结果是一张"整数地图"而不是概率图

cluster()输出的像素值是簇的编号,数据类型是整数。这一点和分类器的classify()不同,后者可以输出概率波段(把ee.Classifier的setOutputMode('PROBABILITY')打开)。聚类没有概率的概念,一个像素要么属于簇 3,要么属于簇 5,没有中间状态。

这带来一个隐性限制:你没法用"置信度"来筛掉那些光谱上模棱两可的像素。比如城市边缘的混合像元,它在光谱上可能介于建筑簇和裸地簇之间,KMeans 会按距离硬性把它分给更近的那个,误差就此产生。想缓解这个问题,只能靠增加簇的数量、加入纹理或指数特征来拉开距离,或者干脆放弃聚类转投监督分类。

提示:cluster()返回的影像没有 mask 之外的信息,原始影像被掩膜掉的区域(云、云影、条带)在结果里同样是无效值,出图时会显示为透明,这是正常的,不用去"修复"它。

2. 送进聚类器之前,影像必须先变成一张干净的表

2.1 波段选择:少而准,远胜全而杂

GEE 的聚类器对波段数量是敏感的。你塞进去 12 个波段和塞进去 5 个波段,跑出来的结果可能天差地别,而且不是"多了更好"。原因在于 KMeans 用的是欧氏距离,每个波段在距离计算里的权重是一样的——波段越多,那些信息量低、噪声大的波段(比如大气校正后仍然残留的短波红外噪声)对距离的贡献就越大,反而把真正有区分度的可见光波段淹没了。

我的习惯是先做一轮波段筛选。对于 Landsat 8/9 的地表反射率产品,SR_B2到SR_B7这六个多光谱波段是基础盘,其中SR_B2(蓝)、SR_B3(绿)、SR_B4(红)、SR_B5(近红外)是区分植被、水体、建筑的核心,SR_B6和SR_B7(两个短波红外)对裸地、建筑、云有额外区分力。SR_B1(海岸气溶胶)我一般直接扔掉,它在大陆内部的信噪比太差。

var bands = ['SR_B2', 'SR_B3', 'SR_B4', 'SR_B5', 'SR_B6', 'SR_B7'];

除了原始波段,我强烈建议加上一两个指数波段。NDVI对植被的区分力比任何单波段都强,加上它之后,植被簇往往会从一个"浑圆的堆"分裂成"茂密植被"和"稀疏植被"两簇,这对后续解译帮助极大。

var ndvi = img.normalizedDifference(['SR_B5', 'SR_B4']).rename('NDVI'); var stack = img.select(bands).addBands(ndvi);

NDBI(归一化建筑指数)在城镇化区域也很好用,公式是(SWIR1 - NIR) / (SWIR1 + NIR),对应SR_B6和SR_B5。但要注意,加上指数之后波段之间的量纲差距会拉大,这就引出下一个必须处理的问题。

2.2 量纲统一:被忽略最多、后果最严重的一步

这是我看过的非监督分类失败案例里,出现频率最高的一个原因,而且它失败得很隐蔽——代码不报错,图也能出,只是结果毫无意义。

举个具体的例子。假设你不仅加了 NDVI,还加了一个高程波段(SRTM,取值范围 0 到 8000 米左右)。NDVI 的范围是 -1 到 1,反射率范围是 0 到 1,而高程的范围是 0 到 8000。在 KMeans 计算欧氏距离的时候,一个像素和质心在高程上差 100 米,产生的距离平方是 10000;而 NDVI 差 0.5,产生的距离平方才 0.25。高程直接以压倒性优势主导了整个聚类过程,结果就是你得到的不是地物分类图,而是一张按海拔分层的等高线图。

解决办法是把所有波段缩放到同一个量级。GEE 里最方便的是unitScale():

var dem = ee.Image('USGS/SRTMGL1_003').select('elevation').unitScale(0, 3000);

unitScale(min, max)会把值线性拉伸到 0 到 1 区间,超出范围的值被截断。这样高程和反射率就站在了同一起跑线上。

对于反射率本身,也要注意 Fractional 的缩放因子。Landsat Collection 2 Level-2 的地表反射率产品,官方存储的是经过放大 10000 再偏移的整数,要还原成真实反射率必须做一次线性变换:

var scaled = img.select(bands) .multiply(0.0000275) .add(-0.2);

这个系数不是随便来的,是产品说明里给的scale factor和offset。如果你忘了这一步,波段值域是 0 到 40000 左右,虽然各波段量纲一致、聚类照样能跑,但后面和 NDVI(-1 到 1)混在一起的时候,NDVI 就被彻底压制了。同样,Sentinel-2 的S2_SR_HARMONIZED要除以 10000。

2.3 掩膜与合成:别让云把样本污染了

聚类对异常值不敏感(毕竟它只找质心),但云和云影会在光谱上形成极端值,这些极端值会自己聚成一簇,白白浪费你的nClusters名额。所以预处理阶段的去云是必须的。

Landsat Collection 2 的QA_PIXEL是位掩码,第 3 位是云、第 4 位是云影、第 1 位是膨胀云。用bitwiseAnd把它们挑出来:

function maskClouds(img) { var qa = img.select('QA_PIXEL'); var clear = qa.bitwiseAnd(1 << 1).eq(0) // 膨胀云 .and(qa.bitwiseAnd(1 << 3).eq(0)) // 云 .and(qa.bitwiseAnd(1 << 4).eq(0)); // 云影 return img.updateMask(clear); }

时间合成我一般用median()而不是mean(),中位数对残余云的鲁棒性明显更好,一个时间窗口里偶尔混进一两景带薄云的影像,中位数几乎不受影响,均值就会被拉偏。窗口长度选在目标地物生长季内,比如做植被分类就取 6 到 9 月。

var composite = col.map(maskClouds).median();

这里有个容易忽略的细节:median()合成之后,某些像素在所有影像里都被掩膜掉了,最终结果里这些像素是无效的。这会直接影响后面采样——sample()会跳过无效像素,你的实际样本数会远小于设定值。这一点我们在下一节细说。

3. sample() 这一步,决定了聚类结果的上限

3.1 numPixels 不是你想的那样,采样前先算笔账

很多人写numPixels: 5000,就觉得拿到了五千个样本。实际上numPixels是"最多返回这么多",实际数量受三个因素制约:研究区内的有效像素总数、scale参数对应的分辨率、以及掩膜情况。

先算个账。假设你的研究区是 0.4 度 × 0.3 度的一块区域,纬度 40 度左右,经度方向约 34 公里,纬度方向约 33 公里,总面积大约 1100 平方公里。在scale: 30的情况下,有效像素数是 1100 / (0.03 × 0.03) ≈ 122 万个。理论上采 5000 个完全没问题。

但如果你的研究区只有 10 公里见方,那就是 100 平方公里,约 11 万个像素,采 5000 个也还行。真正出问题的是研究区再小一点,或者云掩膜覆盖率很高(比如南方的雨季,一景影像 60% 被云挡掉),那numPixels: 5000可能只能拿到几百个。

所以我的做法是先探一下实际能拿到多少:

var train = stack.sample({ region: roi, scale: 30, numPixels: 5000, seed: 42, geometries: false, tileScale: 4 }); print('实际样本数:', train.size());

如果这个数字和你设的差得远,说明你的研究区里有效像素不够,这时候要么扩大研究区,要么降低scale(用更粗的分辨率,比如 60 米),要么接受现实——样本少一点也不是不能跑,只是质心的估计会不那么稳。

3.2 scale 和 tileScale,一个控制精度一个控制内存

scale决定采样时按多大分辨率去取像素。这里有个反直觉的点:采样用的 scale 不必等于影像的原生分辨率。如果你只是想要代表性的光谱样本,用 60 米甚至 100 米采样反而更好——它会自动做一次聚合,样本更能代表"一片区域"而不是"一个像素",对噪声的鲁棒性更强,而且计算量小得多。

而tileScale是另一个维度。它控制的是 GEE 在计算时把影像切成多大的块。默认tileScale: 1,遇到内存报错时把它调到 4 甚至 8,能显著降低单块内存占用。代价是计算时间变长。这个参数在sample()、cluster()、reduceRegion()里都可以用,凡是涉及大区域统计的地方都建议显式写上。

我在做全省尺度的聚类时,scale通常设在 60 到 100 之间,tileScale设 4,sample数量设 8000 到 10000。这个组合在大多数区域都能在几分钟内跑完,出来的质心也足够稳定。

3.3 采样分布不均匀,是个结构性缺陷

sample()的采样策略是在整个region内随机均匀采样。听起来很公平,但对地物分类来说恰恰不公平——因为不同地物的面积占比本来就差得远。一个典型的研究区里,耕地可能占 60%,建成区占 15%,水体占 3%,随机采样得到的样本里耕地的样本数是水体的 20 倍。

KMeans 的质心计算受样本分布影响。耕地样本多,它的质心估计就准;水体样本只有寥寥几十个,质心就飘。更糟的情况是某些稀有地物在整个研究区只有零星分布,随机采样可能一个点都抽不到,这个地物就完全没有对应的簇。

缓解办法有两个。一是分层采样:把研究区按大致的亚区分块,每个块单独采样再merge()。比如你想确保山区和城区都有足够样本,就分别框两个矩形区域各采 3000 个点。二是提高总采样量,用数量换稀有地物的覆盖概率,但这会带来内存压力。我通常用第一种,可控性更强:

var trainA = stack.sample({region: roiUrban, scale: 60, numPixels: 3000, seed: 1, geometries: false}); var trainB = stack.sample({region: roiMountain, scale: 60, numPixels: 3000, seed: 2, geometries: false}); var train = trainA.merge(trainB);

注意geometries: false这个参数。设为false时返回的FeatureCollection不带几何信息,体积小很多,而train()本来就只用属性值,几何信息纯属浪费内存。这个细节能省掉不少"内存超限"的麻烦。

4. 训练与出图:wekaKMeans 的完整调用链

4.1 参数逐个拆开说

ee.Clusterer.wekaKMeans的参数不多,但每一个都影响结果:

参数类型说明我的常用值
nClustersInteger簇的数量5 到 10,按研究区复杂度定
maxIterationsInteger最大迭代次数50,够用,再多收益很小
distanceFunctionString距离度量,Euclidean 或 ManhattanEuclidean
seedInteger随机种子任意固定值,比如 1
dropLastClusterBoolean是否丢弃不完整的簇默认即可

nClusters是最需要斟酌的。设小了,不同地物被强行合并——水体和深色建筑经常被并成一簇;设大了,同一地物被拆成好几簇,比如"茂密林地"和"稀疏林地"分开,你解译的时候还得手动合并回去。经验值是:研究区里明显可区分的地物有几种,就设几加二。比如一个典型的城乡混合区,水体、林地、耕地、建成区、裸地这五类,设 7 是比较合理的,多出来的两个名额留给光谱变异(比如深水浅水、新建筑旧建筑)。

maxIterations设 50 基本不会成为瓶颈。KMeans 在遥感影像这种数据上通常二三十次迭代就收敛了,设太大只是浪费时间。

4.2 从样本到训练好的聚类器

完整的一步是这样的:

var clusterer = ee.Clusterer.wekaKMeans({ nClusters: 7, maxIterations: 50, distanceFunction: 'Euclidean', seed: 1 }).train(train, 4); // 第二个参数是并行数

train()的第二个可选参数是parallelScale,控制训练时的并行度。样本量大的时候(上万条)把它设成 4 甚至 8 能明显提速。

这里有一个隐蔽的坑:train()要求FeatureCollection里每个 feature 的所有属性都是数值型。如果你的样本是从sampleRegions()带过来的,里面可能混进了字符串属性(比如地类名'water'),训练时就会报Property 'xx' is missing或者直接抛类型错误。两种处理方式:采样时用geometries: false并且只保留影像波段对应的属性,或者在训练前.select()一下。

var cleanTrain = train.select(['SR_B2','SR_B3','SR_B4','SR_B5','SR_B6','SR_B7','NDVI']);

train()内部做的事,用大白话说就是:把所有样本点扔进一个七维空间(六个反射率加一个 NDVI),反复迭代找七个"重心",让所有点到自己最近重心的距离平方和最小。找到之后,这七个重心坐标就固定下来了,成为后续分类的依据。

4.3 应用聚类器与固定配色

训练好之后,应用到整张影像是一行的事:

var clustered = stack.cluster(clusterer, 'cluster');

第二个参数是指定输出波段名,默认就是cluster,写不写都行。

出图时很多人图省事用randomVisualizer():

Map.addLayer(clustered.randomVisualizer(), {}, 'clusters');

这个函数确实方便,它给每个簇随机分配颜色。但它有两个问题:一是颜色每次刷新都可能变,你今天出的图和昨天出的图对不上;二是随机的颜色可能让两个光谱相近的簇颜色也相近,肉眼根本分不出来。做正式成果的时候,一定要用固定调色板:

var palette = ['#1f78b4', '#33a02c', '#e31a1c', '#ff7f00', '#6a3d9a', '#b15928', '#a6cee3']; Map.addLayer(clustered, {min: 0, max: 6, palette: palette}, 'clusters');

min: 0, max: 6对应七个簇的编号范围。这里要注意,如果你用了remap改了编号,min/max也要跟着改。

导出到 Drive 的时候,记住cluster()的输出虽然是整数类型,但导出时显式转一下toInt8()更稳妥,文件体积也能小一半:

Export.image.toDrive({ image: clustered.toInt8(), description: 'kmeans_unsupervised', region: roi, scale: 30, maxPixels: 1e13, fileNamePrefix: 'kmeans_result' });

5. 结果怎么读:从一堆编号到真正的地物

5.1 用统计特征给每个簇"挂牌"

这是整篇文章我认为最有价值的一步。拿到七色地图之后,不要凭眼睛猜,用数据说话:算出每个簇在各个波段上的均值,然后对着光谱特征判断它是什么。

var withCluster = stack.addBands(clustered); var stats = withCluster.reduceRegion({ reducer: ee.Reducer.mean().group({ groupField: stack.bandNames().length(), groupName: 'cluster' }), geometry: roi, scale: 60, maxPixels: 1e13, tileScale: 4 }); print(stats);

groupField这里填的是簇波段在合并影像里的索引位置,stack有 7 个波段(6 个反射率加 NDVI),所以簇波段的索引就是 7。这个索引一定要算准,填错了分组就完全乱套。

拿到统计表之后,判断逻辑大致是这样的:

簇编号NIR 均值NDVI 均值SWIR1 均值推断地物
0低< 0低水体
1高> 0.6低茂密植被
2中0.3 到 0.5中稀疏植被/农田
3中< 0.15高建成区
4中< 0.1很高裸土/沙地

水体的判定最干脆:近红外波段反射率极低(水体强吸收),NDVI 为负值。这个组合几乎是判别水体的充分条件,任何簇只要 NDVI 全为负,基本可以直接认定是水体。

建成区和裸地的区分要小心,两者在不同波段的差异很微妙。一般来说建成区在两个短波红外上都有较高反射,而裸土在 SWIR1 上比 SWIR2 更高。如果用肉眼看图能分辨,那结合空间形态判断更靠谱——建成区呈规则的块状或线状,裸土边界更自然。

提示:不要一次把所有簇都定性完。先挑出你最确定的(水体、茂密植被),剩下的反复对照高分底图或者同时期的真彩色合成图,一个一个确认。这个过程急不得,做一次能管很久。

5.2 混淆矩阵:非监督也有正规的精度评价方法

非监督分类没有训练样本,但它可以有验证样本。你手上有一些已知地类的点(哪怕只有二三十个),就能做一个混淆矩阵。

前提是先把簇编号和真实类别对应起来。假设经过上一节的统计判断,你确认簇 1 是植被、簇 3 是建成区、簇 0 是水体,那就用remap转成统一的类别编码:

var classed = clustered.remap( [0, 1, 3], // 原始簇编号 [1, 2, 3] // 目标类别编码:1 水体 2 植被 3 建成区 ).rename('class');

然后拿验证点去采样并计算:

var pred = classed.sampleRegions({ collection: validationPoints, properties: ['class'], scale: 30, tileScale: 4 }); var cm = pred.errorMatrix('class', 'class_1'); print('总体精度:', cm.accuracy()); print('Kappa:', cm.kappa()); print('混淆矩阵:', cm);

后面的class_1是sampleRegions()给影像波段自动加的后缀,你可以在print(pred.first())里确认一下实际字段名。这一步我踩过坑,字段名写错会直接报Property not found。

非监督分类的精度通常不会很高,总体精度能到 70% 到 80% 就算不错了,主要误差来源是混合像元和类别内部的巨大变异(比如同样是建成区,老城区和新开发区的光谱能差出天际)。所以如果项目对精度有硬要求,最终还是要走监督分类。

5.3 类别的合并与后处理

有时候聚类跑出来的七个簇里,有两三个其实是一类地物的不同变异(比如深水、浅水、浑浊水)。这时候需要在编号层面合并:

var merged = clustered.remap( [0, 4, 5], // 三个水体相关的簇 [1, 1, 1] // 全部映射为类别 1 );

remap()的行为要留意:没有出现在from列表里的原始值,默认会变成0(无效值),除非你把defaultValue参数设成别的。所以在做多类合并的时候,要么一次性把七个簇全列出来,要么显式设置defaultValue。

var merged = clustered.remap([0,4,5], [1,1,1], null, 'cluster');

合并完成后,通常还要做一次去除小斑块的处理。GEE 里可以用connectedPixelCount()找出孤立的小连通区再掩膜掉:

var minSize = 8; var count = merged.connectedPixelCount(minSize, true); var cleaned = merged.updateMask(count.gte(minSize));

这个操作对视觉效果的提升非常明显,尤其是建成区边缘那些零散的分类噪声,清掉之后整张图立刻"干净"了。

6. 必然遇到的报错与性能问题

6.1 User memory limit exceeded:九成是采样或统计的问题

这个报错在 GEE 里几乎是家常便饭,尤其是在sample()和reduceRegion()上。它的本质是单个计算任务的内存超过了配额(大约 250MB 量级),注意是单个任务而不是总量,所以解决办法往往不是"少算点",而是"拆开算"。

三个最有效的调整,按优先级排:

第一个,把tileScale从 1 提到 4。这个改动对结果没有任何影响,只是把计算切成更小的块,是所有内存报错的第一反应。

第二个,把scale调粗。scale: 30换成scale: 60,参与计算的像素数直接降到四分之一。对聚类来说,用 60 米分辨率训练质心完全够用,因为你要的是光谱特征的代表性,不是像素级的精度。

第三个,缩小region。如果你的研究区是一个地级市那么大的范围,先把它裁成几个子区域分别跑,再合并结果。这个方案最麻烦但最有效。

还有一种情况是reduceRegion()里的maxPixels设得太小。默认值是 1e7(一千万像素),对于 30 米分辨率的区域很容易超。养成习惯,涉及大区域统计时显式写maxPixels: 1e13。

6.2 Computation timed out:导出比交互式计算更靠谱

在 Code Editor 里点print()或者Map.addLayer()触发的计算,有一个大约五分钟的超时限制。超过就报Computation timed out。这不是代码写错了,纯粹是任务太大。

遇到这个,第一步是把交互式计算改成导出任务:

Export.image.toDrive({ image: clustered, ... });

导出任务的时限宽松得多(几十分钟级别),而且是异步的,提交完可以干别的。先在 Code Editor 里用小区域试通参数,确认结果合理了,再提交大区域的导出任务,这个工作流能省掉大量反复等待的时间。

第二步是检查有没有"不必要的全图计算"。典型的例子是reduceRegion()里没有指定scale,GEE 会试图用影像的原生分辨率去算,像素数瞬间爆炸。任何时候都要显式指定统计的scale,哪怕它和影像分辨率一样。

6.3 每次跑出来结果不一样,是种子的问题还是别的

理论上seed固定了,结果就应该是确定的。但实际使用中,我确实遇到过改了无关参数之后簇编号整体重排的情况。原因在于样本集变了——你改了numPixels、scale、或者影像的时间范围,采样点就变了,质心跟着变,编号重排是必然的。

还有一种更微妙的情况:即使所有参数完全一致,在非常大的区域上重复运行,偶尔也会出现个别簇的编号交换。这跟 GEE 后端的分布式计算有关,不是你能完全控制的。

所以工程上的正确做法是:永远不要在代码里硬编码"簇 3 是水体"这种假设。要么把解译结果手动 remap 一次并存成常量,要么把整条链路写成一个可复现的流程,每次运行后重新确认一遍标签对应关系。我自己的习惯是在代码里留一段注释,写上"本区域 nClusters=7 时,簇 0 水体 / 簇 1 林地 / ……(2024 年 6 月确认)",换数据就要重新确认。

7. 一套可以直接改的模板和参数速查

7.1 完整可运行的最小示例

下面这段代码把前面所有环节串起来了,改一改区域坐标就能跑。用的是 Landsat 8 Collection 2 地表反射率,7 个簇,研究区是一块城乡混合区。

// ===== 1. 研究区 ===== var roi = ee.Geometry.Rectangle([116.20, 39.70, 116.60, 40.00]); // ===== 2. 影像准备 ===== var col = ee.ImageCollection('LANDSAT/LC08/C02/T1_L2') .filterBounds(roi) .filterDate('2023-06-01', '2023-09-30') .filter(ee.Filter.lt('CLOUD_COVER', 15)); function maskClouds(img) { var qa = img.select('QA_PIXEL'); var clear = qa.bitwiseAnd(1 << 1).eq(0) .and(qa.bitwiseAnd(1 << 3).eq(0)) .and(qa.bitwiseAnd(1 << 4).eq(0)); return img.updateMask(clear); } var bands = ['SR_B2', 'SR_B3', 'SR_B4', 'SR_B5', 'SR_B6', 'SR_B7']; var composite = col.map(maskClouds).median() .select(bands) .multiply(0.0000275) .add(-0.2); var ndvi = composite.normalizedDifference(['SR_B5', 'SR_B4']).rename('NDVI'); var stack = composite.addBands(ndvi); // ===== 3. 采样 ===== var train = stack.sample({ region: roi, scale: 60, numPixels: 8000, seed: 42, geometries: false, tileScale: 4 }); print('实际样本数', train.size()); // ===== 4. 训练聚类器 ===== var clusterer = ee.Clusterer.wekaKMeans({ nClusters: 7, maxIterations: 50, distanceFunction: 'Euclidean', seed: 1 }).train(train, 4); // ===== 5. 应用 ===== var clustered = stack.cluster(clusterer, 'cluster'); // ===== 6. 出图 ===== var palette = ['#1f78b4', '#33a02c', '#e31a1c', '#ff7f00', '#6a3d9a', '#b15928', '#a6cee3']; Map.centerObject(roi, 11); Map.addLayer(composite, {bands: ['SR_B4','SR_B3','SR_B2'], min: 0, max: 0.3}, '真彩色'); Map.addLayer(clustered, {min: 0, max: 6, palette: palette}, '聚类结果'); // ===== 7. 各簇光谱统计 ===== var withCluster = stack.addBands(clustered); var stats = withCluster.reduceRegion({ reducer: ee.Reducer.mean().group({ groupField: stack.bandNames().length(), groupName: 'cluster' }), geometry: roi, scale: 60, maxPixels: 1e13, tileScale: 4 }); print('各簇波段均值', stats); // ===== 8. 导出 ===== Export.image.toDrive({ image: clustered.toInt8(), description: 'kmeans_unsupervised', region: roi, scale: 30, maxPixels: 1e13, fileNamePrefix: 'kmeans_result' });

这段代码我在几个不同的气候区跑过,nClusters: 7是个比较稳的起点。南方的植被覆盖区可以加到 8 到 9,西北的荒漠区 5 到 6 就够了。

7.2 参数速查与调整思路

现象最可能的原因调整动作
结果几乎全是同一色波段量纲差异过大检查multiply/add缩放,加unitScale
水体和建筑混在一起nClusters偏小加到 8 到 10 再试
同一地物被拆成多簇nClusters偏大减到 5 到 6,或后期 remap 合并
样本数远小于设定值研究区有效像素不足扩大区域或调粗scale
报内存超限单块计算量太大tileScale: 4,scale调粗,缩小region
出图颜色每次不同用了randomVisualizer()换固定palette
结果呈明显条带状时间窗口内影像太少放宽日期范围或降低云量阈值

7.3 什么时候该放弃非监督,转投监督分类

跑完这一整套流程,你大概会有两个感受:一是非监督确实快,不用准备样本;二是它的精度上限真的有限,而且解译成本高。我的判断标准很直接——如果这个分类成果需要写进报告、需要用面积统计、需要和其他数据做定量比对,那非监督分类不够用,必须做监督分类。

非监督的真正价值在于给监督分类铺路。它能帮你快速回答"这个区域大概有几种可区分的地物"这个问题,还能帮你找到每一类的光谱中心位置。你可以拿聚类结果当成初始样本,在每个簇内部随机选几个点,人工核对后作为训练样本,这样选出来的样本比盲目框选代表性强得多,能显著缩短监督分类的样本准备时间。

另外一个我自己常用的组合是:先用wekaCascadeKMeans扫一遍,让它自己在 3 到 12 的范围内搜索最优簇数,看看它给出的结果。如果它建议 6 个簇,那就说明你这个区域的光谱可分性大概就是六类左右,用普通的 KMeans 时设 6 或者 7 都不会错太远。这个探路的过程成本很低,但能省掉不少反复试参数的功夫。

提示:wekaCascadeKMeans的参数是minClusters和maxClusters,不要传nClusters,传了会被忽略或者报错。这个细节文档里写得不算显眼,我第一次用的时候卡了好一会儿。

后续如果要做监督分类,这条链路上的大半准备代码可以直接复用——影像预处理、采样、导出,唯一要换的就是把Clusterer换成Classifier(比如ee.Classifier.smileRandomForest),再加上样本标注这一步。从这个角度看,把非监督分类的流程走一遍,本身就是在给后面更复杂的工作打地基。