GEE：探索北运河流域水体面积 30 年来的变化 [逐年变化

---

前言

使用GEE分析流域内水体面积变化特征，以北运河流域为例

同区域分析参见：
GEE：北运河流域土地利用的重分类

分析区域如下：

---

一、分析步骤

主要分析步骤：

二、python代码

1、加载流域范围

python代码如下

# 添加显示底图，选取高德地图
Map = geemap.Map()
Map.add_basemap('HYBRID')

# 导入流域范围
basin_shp = '../world_basins/byh.shp'
basin_bj = geemap.shp_to_ee(basin_shp)
roi = basin_bj.geometry()
# 添加显示
Map.addLayer(ee.Image().paint(roi, 0, 2), {'palette':'red'}, 'roi')
Map.centerObject(basin_bj)
Map

结果如下：

2、影像截取与合成

• geemap.landsat_timeseries()

python代码如下：

# 导入影像数据，年度合成
images = geemap.landsat_timeseries(roi=roi,
                                    start_year=1984,
                                    end_year=2020,
                                    start_date='01-01',
                                    end_date='12-31')
# 查看影像波段
first = images.first()
first.bandNames().getInfo()
# ['Blue', 'Green', 'Red', 'NIR', 'SWIR1', 'SWIR2', 'pixel_qa', 'NBR']
# 添加显示
Map.addLayer(first, {'bands':['SWIR2', 'SWIR1', 'Green'], 'min':0, 'max':3000}, 'first image')

结果如下：

3、NDWI计算与水面提取

• normalizedDifference()

• ee.Image.pixelArea()
  - ee.Image.reduceRegion()
  

python代码如下：

def ndwi_cac(image):
    ndwi_image = image.normalizedDifference(['Green', 'SWIR1']).rename('ndwi')
    return ndwi_image.gt(0).selfMask()
    
# 获得水面
ndwi_water = images.map(ndwi_cac)
Map.addLayer(ndwi_water.first(), {'palette':'blue'}, 'first water')

# 计算水面面积，单位按照km2计
def cal_area(image):
    pixel_area = image.multiply(ee.Image.pixelArea()).divide(1e6)
    img_area = pixel_area.reduceRegion(**{
        'reducer':ee.Reducer.sum(),
        'geometry':roi,
        'scale':30,
        'maxPixels':1e13
    })
    return image.set({'water_area':img_area})

water_areas = ndwi_water.map(cal_area)
water_stats = water_areas.aggregate_array('water_area').getInfo()
water_stats

结果如下：
[{‘ndwi’: 35.52171285398048},
{‘ndwi’: 49.82769552462259},
{‘ndwi’: 95.52709095289183},
{‘ndwi’: 88.0812142082916},
{‘ndwi’: 93.93069886911051},
{‘ndwi’: 80.60827039904903},
{‘ndwi’: 136.78296227087307},
{‘ndwi’: 108.33794441999756},
{‘ndwi’: 171.22477916253243},
{‘ndwi’: 127.52933278774321},
{‘ndwi’: 86.3748384304267},
{‘ndwi’: 112.22821871601487},
{‘ndwi’: 101.4823960050295},
{‘ndwi’: 89.49672596974511},
{‘ndwi’: 295.2877416742835},
{‘ndwi’: 80.7449547841241},
{‘ndwi’: 120.11593294853824},
{‘ndwi’: 134.62668476349447},
{‘ndwi’: 149.84151247550966},
{‘ndwi’: 236.29430841558016},
{‘ndwi’: 158.22021181601244},
{‘ndwi’: 70.20089856972656},
{‘ndwi’: 97.93744088662496},
{‘ndwi’: 108.04200576387466},
{‘ndwi’: 107.22519927537127},
{‘ndwi’: 63.4901824073373},
{‘ndwi’: 85.06409974042677},
{‘ndwi’: 74.55227930805393},
{‘ndwi’: 221.75472506976308},
{‘ndwi’: 85.00424085790755},
{‘ndwi’: 60.99660440626263},
{‘ndwi’: 68.4423834487381},
{‘ndwi’: 67.66773335986407},
{‘ndwi’: 84.48483915371168},
{‘ndwi’: 76.99267855924732},
{‘ndwi’: 58.5695615691552},
{‘ndwi’: 87.93523549212745}]

4、统计分析

python代码如下：

# 时间序列变化柱状图显示
import matplotlib.pyplot as plt

x = list(range(1984, 2021))
y = [item.get('ndwi') for item in water_stats]

plt.bar(x, y, align='center', alpha=0.5)
plt.ylabel('Area hm2')
plt.title('surface water dynamics in BYH')
plt.show()

---

三、对比分析

拟采用JRC数据直接对水面变化进行分析，其python代码如下：
【有待分析。。。】

结果如下：
二者计算结果的差异比较分析：

1. 可以看出，

---

四、小结与讨论

1. 分析了北运河流域长时间序列水面变化特征
2. 有待探索的问题，可考虑不同分区的水面变化，分布不同行政区或者流域上下游区域的水面时空变化

---

参考：

• https://geemap.org/notebooks/tn_surface_water/
• https://developers.google.com/earth-engine/apidocs/ee-image-pixelarea?hl=en