保姆级教程：在Google Earth Engine (GEE) 上复现遥感生态指数RSEI（Landsat 8数据）

从零到一GEE平台遥感生态指数(RSEI)全流程实战解析遥感生态指数RSEI作为评估区域生态环境质量的综合指标近年来在学术研究和实际应用中备受关注。对于刚接触Google Earth EngineGEE平台的研究者而言如何利用Landsat 8数据快速实现RSEI计算往往面临代码理解困难、参数设置模糊、报错排查无门等痛点。本文将采用理论代码调试三维度教学法带您完整走通从数据预处理到结果可视化的全链路流程。1. 环境准备与数据加载1.1 研究区定义技巧在GEE中定义研究区时建议先通过交互式地图工具绘制大致范围再手动调整坐标值。以下代码展示了如何定义青岛西海岸新区的研究范围var roi ee.Geometry.Polygon( [[[120.121, 35.975], [120.121, 35.886], [120.257, 35.886], [120.257, 35.975]]], null, false ); Map.centerObject(roi, 10); // 第二个参数控制缩放级别常见问题排查坐标顺序必须为逆时针方向最后一个点需与第一个点闭合若出现Geometry is invalid错误检查坐标格式是否正确1.2 数据加载与云掩膜Landsat 8 TOA数据需特别注意云量筛选和去云处理。推荐使用质量评估波段(BQA)进行精确去云function removeCloud(image) { var qa image.select(BQA); var cloudMask qa.bitwiseAnd(1 4).eq(0); // 第4位表示云 var cloudShadowMask qa.bitwiseAnd(1 8).eq(0); // 第8位表示云阴影 return image.updateMask(cloudMask.and(cloudShadowMask)); } var L8 ee.ImageCollection(LANDSAT/LC08/C01/T1_TOA) .filterBounds(roi) .filterDate(2018-01-01, 2019-12-31) .filterMetadata(CLOUD_COVER, less_than, 30) // 云量阈值建议设为30% .map(removeCloud);注意不同Landsat数据集的BQA位掩码可能不同需查阅官方文档确认2. 生态指标计算详解2.1 四大核心指标实现RSEI包含绿度(NDVI)、湿度(Wet)、热度(LST)、干度(NDBSI)四个分量指标湿度分量计算基于缨帽变换系数var Wet img.expression( B*0.1509 G*0.1973 R*0.3279 NIR*0.3406 SWIR1*(-0.7112) SWIR2*(-0.4572), { B: img.select(B2), // Blue G: img.select(B3), // Green R: img.select(B4), // Red NIR: img.select(B5), // Near Infrared SWIR1: img.select(B6), // Shortwave Infrared 1 SWIR2: img.select(B7) // Shortwave Infrared 2 } );干度指数(NDBSI)复合计算var ibi img.expression( (2*SWIR1/(SWIR1NIR) - (NIR/(NIRRED) GREEN/(GREENSWIR1))) / (2*SWIR1/(SWIR1NIR) (NIR/(NIRRED) GREEN/(GREENSWIR1))), { SWIR1: img.select(B6), NIR: img.select(B5), RED: img.select(B4), GREEN: img.select(B3) } ); var si img.expression( ((SWIR1RED) - (NIRBLUE)) / ((SWIR1RED) (NIRBLUE)), { SWIR1: img.select(B6), NIR: img.select(B5), RED: img.select(B4), BLUE: img.select(B2) } ); var ndbsi ibi.add(si).divide(2);2.2 温度数据融合策略由于Landsat 8地表温度产品精度限制推荐融合MODIS温度数据var lst ee.ImageCollection(MODIS/006/MOD11A1) .filterDate(2018-01-01, 2019-12-31) .select([LST_Day_1km, LST_Night_1km]) .mean() .reproject(EPSG:4326, null, 1000); // 重投影保证分辨率一致 var lst_mean lst.expression( (Day Night) / 2, { Day: lst.select(LST_Day_1km), Night: lst.select(LST_Night_1km) } );3. 数据标准化与PCA分析3.1 归一化处理关键参数归一化是消除量纲影响的关键步骤需注意scale参数与maxPixels设置var img_normalize function(img) { var minMax img.reduceRegion({ reducer: ee.Reducer.minMax(), geometry: roi, scale: 30, // 应与原始数据分辨率匹配 maxPixels: 1e13 // 大区域需增大此值 }); return ee.ImageCollection.fromImages( img.bandNames().map(function(name) { var band img.select(name); return band.unitScale( ee.Number(minMax.get(name.cat(_min))), ee.Number(minMax.get(name.cat(_max))) ); }) ).toBands().rename(img.bandNames()); };3.2 主成分分析实现细节PCA是RSEI的核心算法需理解协方差矩阵计算过程var pca function(img) { // 均值中心化 var meanDict img.reduceRegion({ reducer: ee.Reducer.mean(), geometry: roi, scale: 30, maxPixels: 1e13 }); var centered img.subtract(ee.Image.constant(meanDict.values(img.bandNames()))); // 协方差矩阵计算 var arrays centered.toArray(); var covar arrays.reduceRegion({ reducer: ee.Reducer.centeredCovariance(), geometry: roi, scale: 30, maxPixels: 1e13 }); // 特征分解 var covarArray ee.Array(covar.get(array)); var eigens covarArray.eigen(); // 主成分计算 var pcImage ee.Image(eigens.slice(1, 1)) .matrixMultiply(arrays.toArray(1)) .arrayProject([0]) .arrayFlatten([[PC1, PC2, PC3, PC4]]); return pcImage; };提示第一主成分(PC1)通常包含最大方差信息是构建RSEI的基础4. RSEI计算与可视化优化4.1 指数构建与标准化最终RSEI计算需注意方向一致性生态质量与PC1通常负相关var RSEI PCA_img.select(PC1).multiply(-1).add(1).rename(RSEI); var RSEI_normalized img_normalize(RSEI); // 再次归一化到[0,1]范围4.2 可视化参数调优推荐使用渐变色带突出生态质量差异var visParams { min: 0, max: 1, palette: [ d7191c, fdae61, ffffbf, a6d96a, 1a9641 // 红-黄-绿渐变 ] }; Map.addLayer(RSEI_normalized, visParams, RSEI);进阶技巧使用image.updateMask()去除水体等干扰区域结合image.sample()提取典型区域值进行验证通过Export.image.toDrive()导出GeoTIFF结果5. 全流程代码封装与调试将各模块封装为独立函数便于复用和调试function calculateRSEI(roi, startDate, endDate) { // 数据加载与预处理 var L8 ee.ImageCollection(LANDSAT/LC08/C01/T1_TOA) .filterBounds(roi) .filterDate(startDate, endDate) .map(removeCloud); // 生态指标计算 var indicators calculateIndicators(L8, roi); // 标准化与PCA var normalized indicators.map(img_normalize); var pcaResults normalized.map(pca); // RSEI计算 var rsei pcaResults.map(function(img) { return img.select(PC1).multiply(-1).add(1).rename(RSEI); }); return rsei.mean(); // 返回时间序列均值 } // 调用示例 var finalRSEI calculateRSEI(roi, 2018-01-01, 2018-12-31);调试建议使用print()输出中间结果分步执行代码确认各阶段输出对异常值检查波段计算表达式缩小研究区范围测试运行效率实际项目中遇到过最棘手的问题是PCA计算超时最终通过优化scale参数和适当缩小研究区范围解决。建议初次运行时先使用小区域测试确认无误后再扩展到大范围计算。