Zonal statistics

Uncomment the following line to install geemap if needed.

In [ ]:

# !pip install geemap

In [ ]:

import ee
import geemap
import os

Create an interactive map¶

In [ ]:

Map = geemap.Map()
Map

Add Earth Engine data¶

In [ ]:

# Add Earth Engine dataset
dem = ee.Image('USGS/SRTMGL1_003')

# Set visualization parameters.
dem_vis = {
  'min': 0,
  'max': 4000,
  'palette': ['006633', 'E5FFCC', '662A00', 'D8D8D8', 'F5F5F5']}

# Add Earth Engine DEM to map
Map.addLayer(dem, dem_vis, 'SRTM DEM')

# Add Landsat data to map
landsat = ee.Image('LE7_TOA_5YEAR/1999_2003')

landsat_vis = {
    'bands': ['B4', 'B3', 'B2'], 
    'gamma': 1.4
}
Map.addLayer(landsat, landsat_vis, "Landsat", False)

states = ee.FeatureCollection("TIGER/2018/States")
Map.addLayer(states, {}, 'US States')

Compute zonal statistics for one image¶

In [ ]:

out_dir = os.path.expanduser('~/Downloads')
if not os.path.exists(out_dir):
    os.makedirs(out_dir)

In [ ]:

out_dem_stats = os.path.join(out_dir, 'dem_stats.csv')  

# Allowed output formats: csv, shp, json, kml, kmz
# Allowed statistics type: MEAN, MAXIMUM, MINIMUM, MEDIAN, STD, MIN_MAX, VARIANCE, SUM
geemap.zonal_statistics(dem, states, out_dem_stats, statistics_type='MEAN', scale=1000)

In [ ]:

out_landsat_stats = os.path.join(out_dir, 'landsat_stats.csv')  
geemap.zonal_statistics(landsat, states, out_landsat_stats, statistics_type='SUM', scale=1000)

Compute zonal statistics for time-series images¶

In [ ]:

Map = geemap.Map()

collection = ee.ImageCollection('MODIS/MCD43A4_006_NDVI') \
                  .filter(ee.Filter.date('2018-04-01', '2018-05-01')) \
                  .select("NDVI")\

vis_params = {
  'min': 0.0,
  'max': 1.0,
  'palette': [
    'FFFFFF', 'CE7E45', 'DF923D', 'F1B555', 'FCD163', '99B718', '74A901',
    '66A000', '529400', '3E8601', '207401', '056201', '004C00', '023B01',
    '012E01', '011D01', '011301'
  ],
}

first_image = collection.first()

Map.addLayer(first_image, vis_params, "First image")
Map.setCenter(-7.03125, 31.0529339857, 2)
Map

In [ ]:

modis = collection.toBands()
Map.addLayer(modis, {}, "MODIS Time series", False)

In [ ]:

out_landsat_stats = os.path.join(out_dir, 'ndvi.csv')  
geemap.zonal_statistics(modis, states, out_landsat_stats, statistics_type='MEAN', scale=1000)

---

Last update: 2021-03-17