手把手教你用Google Earth Engine分析全球森林变化（附Python代码）

基于Google Earth Engine的全球森林动态监测实战指南当我们需要分析跨越数十年的全球森林覆盖变化时传统遥感处理方法往往面临数据量大、计算资源有限、处理流程繁琐等挑战。Google Earth EngineGEE的出现彻底改变了这一局面它将PB级的卫星影像与强大的云计算能力相结合使得全球尺度的环境监测变得前所未有的高效。1. GEE平台核心优势与森林监测应用场景GEE作为遥感云计算平台的代表主要具备三大核心优势海量数据即时访问集成超过40年的Landsat系列数据、Sentinel-1/2卫星数据以及MODIS等关键数据集所有数据已经过预处理和校正并行计算能力借助Google云基础设施可在几分钟内完成传统工作站需要数周计算的全球尺度分析开放API生态同时支持JavaScript和Python接口便于集成到科研工作流中在森林监测领域GEE特别适合以下应用场景应用方向典型指标数据需求森林覆盖变化年度变化率、砍伐热点Landsat/Sentinel-2森林健康评估NDVI/EVI指数、叶面积指数MODIS/Landsat碳储量估算生物量、树高模型GEDI/LiDAR 光学数据火灾影响评估烧毁严重度指数(NBR)Sentinel-2/Landsat提示GEE的数据目录每天新增约4000张影像建议定期检查是否有新数据集发布2. 环境配置与数据准备2.1 Python API初始化配置首先需要安装Earth Engine Python API包并完成认证pip install earthengine-api然后进行账户认证和初始化import ee # 触发认证流程 ee.Authenticate() # 初始化API ee.Initialize()2.2 关键数据集加载GEE提供了多种森林相关数据集以下是常用数据的加载方式# 加载 Hansen 全球森林变化数据集 hansen ee.Image(UMD/hansen/global_forest_change_2022_v1_10) # 加载Landsat 8地表反射率数据 landsat ee.ImageCollection(LANDSAT/LC08/C02/T1_L2) # 加载Sentinel-2 MSI数据 sentinel2 ee.ImageCollection(COPERNICUS/S2_SR)对于特定区域的数据筛选可以使用以下方法# 定义研究区域以亚马逊雨林为例 amazon ee.Geometry.Rectangle([-74, -12, -44, 4]) # 筛选2020-2023年的Landsat数据 landsat_filtered landsat.filterBounds(amazon) \ .filterDate(2020-01-01, 2023-12-31) \ .filter(ee.Filter.lt(CLOUD_COVER, 20))3. 森林变化检测核心技术实现3.1 基于NDVI的植被动态监测NDVI归一化植被指数是监测植被健康的经典指标def calculate_ndvi(image): # Landsat 8的NDVI计算 ndvi image.normalizedDifference([SR_B5,SR_B4]).rename(NDVI) return image.addBands(ndvi) # 对筛选后的影像集合计算NDVI landsat_ndvi landsat_filtered.map(calculate_ndvi) # 生成NDVI时间序列动画 ndvi_params { min: 0, max: 1, palette: [red, yellow, green] } Map.addLayer(landsat_ndvi.select(NDVI), ndvi_params, NDVI Time Series)3.2 森林损失与增益分析利用Hansen数据集可以直接提取森林变化信息# 提取2000-2022年森林损失数据 loss hansen.select(loss) loss_area loss.gt(0).multiply(ee.Image.pixelArea()) # 计算研究区域内总损失面积(平方公里) total_loss loss_area.reduceRegion( reduceree.Reducer.sum(), geometryamazon, scale30, maxPixels1e13 ).get(loss) print(Total forest loss (km²):, ee.Number(total_loss).divide(1e6).getInfo())3.3 变化检测高级方法对于更精确的变化检测可以使用光谱时间序列分析方法# 构建年度合成NDVI影像集 def create_year_composite(year): start ee.Date.fromYMD(year, 1, 1) end start.advance(1, year) composite landsat.filterDate(start, end) \ .map(calculate_ndvi) \ .median() return composite.set(year, year) years ee.List.sequence(2000, 2022) annual_composites ee.ImageCollection(years.map(create_year_composite)) # 计算NDVI趋势 trend annual_composites.select(NDVI) \ .toArray() \ .reduce(ee.Reducer.linearFit())4. 结果可视化与导出4.1 交互式地图可视化GEE支持多种可视化方式展示分析结果# 定义可视化参数 loss_vis {palette: [red], min: 0, max: 1} gain_vis {palette: [green], min: 0, max: 1} # 创建交互地图 Map geemap.Map() Map.centerObject(amazon, 6) Map.addLayer(loss, loss_vis, Forest Loss) Map.addLayer(hansen.select(gain), gain_vis, Forest Gain) Map.addLayer(trend.select(scale), {min: -0.05, max: 0.05, palette: [red, white, green]}, NDVI Trend) Map4.2 数据导出到本地对于需要进一步分析的结果可以导出到Google Drive# 导出NDVI趋势结果 export_task ee.batch.Export.image.toDrive( imagetrend.select(scale), descriptionNDVI_Trend_Export, scale30, regionamazon, fileFormatGeoTIFF, maxPixels1e13 ) export_task.start()5. 实际应用中的优化策略在长期使用GEE进行森林监测的过程中有几个关键优化点值得注意数据质量控制云覆盖是光学遥感的主要干扰因素推荐使用以下云掩膜方法def mask_clouds(image): # Sentinel-2云掩膜示例 qa image.select(QA60) cloud_mask qa.bitwiseAnd(1 10).Or(qa.bitwiseAnd(1 11)) return image.updateMask(cloud_mask.Not())计算效率优化对小区域分析使用reduceRegion()而非reduceRegions()对大型计算任务使用Export而非getInfo()合理设置scale和maxPixels参数时间序列分析技巧使用median()而非mean()进行时间合成考虑物候变化使用相同季节影像进行比较对缺失数据使用线性插值填补在亚马逊雨林监测项目中采用上述方法后原本需要数周完成的2000-2022年森林变化分析现在只需15分钟即可获得初步结果且可以直接生成交互式可视化图表供决策参考。