【数据分析】基于GEE的2000-2023年逐年归一化差值水分指数（NDMI）获取-以成都市为例

<think>好的，用户想用GEE计算并导出2002-2022年桂林市的RSEI数据，并实现可视化。首先，我需要回忆一下RSEI的组成。RSEI包括四个指标：NDVI（植被指数）、LST（地表温度）、Wet（湿度指数）和NDBSI（建筑指数），通常通过主成分分析（PCA）来合成。 首先，用户需要定义研究区域，也就是桂林市的边界。根据引用[2]中的方法，可能需要使用GEE的行政区划数据来筛选桂林市的范围。比如，引用[2]提到通过筛选中国省市边界数据构建研究区域，这里应该类似，需要找到桂林市的FeatureCollection。 接下来是数据选择部分。RSEI通常使用Landsat数据，因为其空间分辨率较高。但引用[3]中提到使用MODIS数据，而用户的问题可能希望更高分辨率的结果，所以可能需要混合使用不同Landsat卫星的数据，比如Landsat 5、7、8，像引用[2]中那样分时间段处理。例如，2002-2011年用Landsat 5，2012-2013用Landsat 7，2014-2022用Landsat 8。需要注意不同卫星的传感器差异，可能需要做辐射校正或波段调整。 然后，影像的过滤条件，比如云覆盖率，引用[2]中用了20%的云量过滤，这里可能也需要类似的处理。此外，时间范围设定为每年，可能需要逐年处理，比如每年生成一个合成影像，如中值合成。 接下来是计算各个指数。NDVI是近红外和红波段的差值比和值之和。LST需要根据热红外波段计算，可能需要大气校正，但GEE可能有现成的产品。Wet指数可能是指湿度分量，通常用缨帽变换的湿度分量，或者使用特定的波段组合。NDBSI是建筑指数，可能由SI（裸土指数）和IBI（建筑指数）组合而成，具体公式可能需要查阅相关文献。 然后，各指数归一化处理，将四个指标标准化到0-1范围，以便进行PCA分析。主成分分析在GEE中可以用ee.Array.covariance等方法实现，或者使用内置的PCA模块。通常RSEI使用第一主成分，但需要确认是否需要反转，因为有时第一主成分可能代表的是负面生态状况。 计算完每年的RSEI后，需要逐年导出到Google Drive或Google Cloud Storage。导出时需要注意设置区域、尺度、分辨率等参数，确保数据正确。 可视化部分，需要在GEE中使用地图图层添加每年的RSEI结果，设定颜色梯度，比如从红色到绿色表示生态指数从低到高。同时，可能需要生成时间序列图表，显示RSEI的年际变化。 可能的挑战包括不同Landsat卫星数据的波段名称和位置不同，需要统一处理。例如，Landsat 5的波段名称与Landsat 8不同，需要调整代码以适应不同卫星。此外，云遮蔽处理可能需要额外的步骤，如使用QA波段或现有的云掩膜算法。 用户可能还关心如何处理缺失年份的数据，比如某一年云覆盖过多导致无法获得有效影像，可能需要使用插值或其他方法填补。此外，RSEI的计算是否有现成的GEE脚本可供参考，引用[3]提到了一个博客可能有相关代码，可能需要参考该资源。 最后，确保代码高效运行，因为处理长时间序列的数据可能需要较大的计算资源，合理使用GEE的批量处理功能和导出任务，避免超时或内存不足的问题。</think>### 分步解决方案：GEE计算桂林市RSEI遥感生态指数 #### **一、数据准备** 1. **定义研究区域** ```javascript var roi = ee.FeatureCollection("users/your_account/Guilin_boundary"); // 或使用GEE内置行政区划数据筛选桂林市 var china = ee.FeatureCollection("USDOS/LSIB_SIMPLE/2017"); var guilin = china.filter(ee.Filter.eq('name', 'Guilin')); ``` 2. **数据选择与过滤** - 整合Landsat 5/7/8数据，按年份分段处理： ```javascript var landsat5 = ee.ImageCollection("LANDSAT/LT05/C02/T1_L2") .filterDate('2002-01-01', '2012-01-01') .filterBounds(roi) .filter(ee.Filter.lte('CLOUD_COVER', 20)); var landsat7 = ee.ImageCollection("LANDSAT/LE07/C02/T1_L2") .filterDate('2012-01-01', '2014-01-01') .filterBounds(roi) .filter(ee.Filter.lte('CLOUD_COVER', 20)); var landsat8 = ee.ImageCollection("LANDSAT/LC08/C02/T1_L2") .filterDate('2014-01-01', '2023-01-01') .filterBounds(roi) .filter(ee.Filter.lte('CLOUD_COVER', 20)); ``` #### **二、计算RSEI组成指数** 1. **NDVI（归一化植被指数）** $$ \text{NDVI} = \frac{\text{NIR} - \text{Red}}{\text{NIR} + \text{Red}} $$ ```javascript var addNDVI = function(image) { var nir = image.select('SR_B4'); // Landsat 5/7/8波段需对应调整 var red = image.select('SR_B3'); return image.addBands(nir.subtract(red).divide(nir.add(red)).rename('NDVI')); }; ``` 2. **LST（地表温度）** - 使用Landsat热红外波段（需辐射定标）： $$ \text{LST} = \frac{\text{TIRS}}{10} $$ [^3] ```javascript var addLST = function(image) { var lst = image.select('ST_B6').multiply(0.1); return image.addBands(lst.rename('LST')); }; ``` 3. **Wet（湿度指数）** - 缨帽变换湿度分量公式： $$ \text{Wet} = 0.1511B_1 + 0.1972B_2 + 0.3283B_3 + 0.3408B_4 - 0.7117B_5 - 0.4559B_7 $$ 4. **NDBSI（建筑指数）** $$ \text{NDBSI} = \frac{\text{SI} + \text{IBI}}{2} $$ - SI（裸土指数）： $$ \text{SI} = \frac{(B_5 + B_3) - (B_4 + B_1)}{(B_5 + B_3) + (B_4 + B_1)} $$ - IBI（建筑指数）： $$ \text{IBI} = \frac{2B_5/(B_5 + B_4) - [B_2/(B_2 + B_5) + B_4/(B_4 + B_3)]}{2B_5/(B_5 + B_4) + [B_2/(B_2 + B_5) + B_4/(B_4 + B_3)]} $$ #### **三、主成分分析与RSEI合成** 1. **归一化处理** ```javascript var normalize = function(image) { var minMax = image.reduceRegion({ reducer: ee.Reducer.minMax(), geometry: roi, scale: 30, maxPixels: 1e13 }); return image.subtract(minMax.get('min')).divide(minMax.get('max').subtract(minMax.get('min'))); }; ``` 2. **PCA计算** ```javascript var calculateRSEI = function(year) { var annualCollection = mergedCollection.filter(ee.Filter.calendarRange(year, year, 'year')); var composite = annualCollection.median(); var pca = ee.Array(composite.select(['NDVI','LST','Wet','NDBSI']).toArray()).pca(); var rsei = ee.Image(pca.project([0])).rename('RSEI'); return rsei; }; ``` #### **四、数据导出与可视化** 1. **批量导出** ```javascript for (var y = 2002; y <= 2022; y++) { var rseiImg = calculateRSEI(y); Export.image.toDrive({ image: rseiImg, description: 'RSEI_Guilin_' + y, scale: 30, region: roi, maxPixels: 1e13 }); } ``` 2. **可视化** ```javascript var visParams = {min: -1, max: 1, palette: ['red', 'yellow', 'green']}; Map.addLayer(rseiImg.clip(roi), visParams, 'RSEI ' + year); ```