CoCalc -- GEE_Workshop

GitHub Repository: giswqs/geemap
Path: blob/master/docs/workshops/GEE_Workshop_2025.ipynb
²³¹³ views

Kernel: Python 3 (ipykernel)

Geospatial Data Analysis and Visualization with Earth Engine

Notebook: https://geemap.org/workshops/GEE_Workshop_2025
Earth Engine: https://earthengine.google.com
Geemap: https://geemap.org

Introduction

This notebook contains the materials for the workshop Geospatial Data Analysis and Visualization with Earth Engine at the 第八届地球空间大数据与云计算研讨会.

Agenda

This workshop covers the following topics:

Colab setup
Data visualization in 3D
Creating timelapse animations
Exporting Earth Engine data
Charts

Prerequisites

To use geemap and the Earth Engine Python API, you must register for an Earth Engine account and follow the instructions here to create a Cloud Project. Earth Engine is free for noncommercial and research use. To test whether you can use authenticate the Earth Engine Python API, please run this notebook on Google Colab.

Technical requirements

Install packages

conda create -n gee python=3.12
conda activate gee
conda install -c conda-forge mamba
mamba install -c conda-forge geemap leafmap maplibre

In [ ]:

# %pip install "geemap[workshop]" leafmap maplibre

Import libraries

In [ ]:

import ee
import geemap
import leafmap.maplibregl as leafmap

In [ ]:

ee.Authenticate()

In [ ]:

ee.Initialize(project="YOUR-PROJECT-ID")

Colab setup

Uncomment the following line to get the Earth Engine authorization token. Please treat your token with care and don't share it with anyone. Copy the token to the clipboard.

In [ ]:

# geemap.get_ee_token()

Open your Google Colab notebook and click on the secrets tab.
Create a new secret with the name EARTHENGINE_TOKEN.
Paste the content from the clipboard into the Value input box of the created secret.
Toggle the button on the left to allow notebook access to the secret.

Data Visualization in 3D

Globe projection

In [ ]:

m = leafmap.Map(center=[-100, 40], zoom=3, style="liberty")
m.add_globe_control()
m

Overture buildings

In [ ]:

m = leafmap.Map(center=[-100, 40], zoom=3, style="positron", projection="globe")
m.add_overture_3d_buildings()
m

Add basemaps

In [ ]:

m = leafmap.Map(center=[-100, 40], zoom=3, style="positron", projection="globe")
m.add_basemap("Esri.WorldImagery")
m.add_overture_3d_buildings()
m.add_layer_control()
m

Visualize Earth Engine data

To use the 3d-terrain basemap, you will need an API key from MapTiler. Once you have the API key, you can uncomment the following code block and replace YOUR_API_KEY with your actual API key. Then, run the code block code to set the API key as an environment variable.

In [ ]:

# import os
# os.environ["MAPTILER_KEY"] = "YOUR_API_KEY"

In [ ]:

m = leafmap.Map(style="3d-terrain", projection="globe")
m.add_ee_layer(asset_id="ESA/WorldCover/v200", opacity=0.5)
m.add_legend(builtin_legend="ESA_WorldCover", title="ESA Landcover")
m.add_overture_3d_buildings()
m.add_layer_control()
m

In [ ]:

m.layer_interact()

Land cover data.

In [ ]:

m = leafmap.Map(style="3d-terrain", projection="globe")
dataset = ee.ImageCollection("ESA/WorldCover/v200").first()
vis_params = {"bands": ["Map"]}
m.add_ee_layer(dataset, vis_params, name="ESA Worldcover", opacity=0.5)
m.add_legend(builtin_legend="ESA_WorldCover", title="ESA Landcover")
m.add_layer_control()
m

Night time light data.

In [ ]:

m = leafmap.Map(style="darkmatter", projection="globe")
dataset = ee.ImageCollection("NOAA/VIIRS/DNB/ANNUAL_V22").filter(
    ee.Filter.date("2022-01-01", "2023-01-01")
)
nighttime = dataset.select("maximum")
nighttimeVis = {"min": 0.0, "max": 60.0}
m.add_ee_layer(nighttime, nighttimeVis, name="Nighttime")

countries = ee.FeatureCollection("USDOS/LSIB_SIMPLE/2017")
styleParams = {
    "fillColor": "00000000",
    "color": "ff0000",
    "width": 1.0,
}
countries = countries.style(**styleParams)
m.add_ee_layer(countries, {}, name="Country Boundaries")
m.add_layer_control()
m

Creating timelapse animations

Landsat

In [ ]:

m = geemap.Map(center=[22.95459, 113.94078], zoom=10)
m

Pan and zoom the map to an area of interest. Use the drawing tools to draw a rectangle on the map. If no rectangle is drawn, the default rectangle shown below will be used.

In [ ]:

m.user_roi_coords()

In [ ]:

roi = m.user_roi
if roi is None:
    roi = ee.Geometry.BBox(113.5059, 22.6805, 114.2939, 23.1195)
    m.add_layer(roi)
    m.center_object(roi)

In [ ]:

timelapse = geemap.landsat_timelapse(
    roi,
    out_gif="Dongguan.gif",
    start_year=1986,
    end_year=2025,
    start_date="01-01",
    end_date="12-31",
    bands=["SWIR1", "NIR", "Red"],
    frames_per_second=5,
    title="Landsat Timelapse",
    progress_bar_color="blue",
    mp4=True,
)
geemap.show_image(timelapse)

In [ ]:

m = geemap.Map(center=[22.95459, 113.94078], zoom=10)
m.add_gui("timelapse")
m

In [ ]:

m = geemap.Map()
roi = ee.Geometry.BBox(113.6278, 22.639, 113.8983, 22.8024)
m.add_layer(roi)
m.center_object(roi)
m

In [ ]:

timelapse = geemap.landsat_timelapse(
    roi,
    out_gif="Donguan2.gif",
    start_year=1986,
    end_year=2025,
    start_date="01-01",
    end_date="12-31",
    bands=["SWIR1", "NIR", "Red"],
    frames_per_second=5,
    title="Dongguan, China",
    font_color="red",
)
geemap.show_image(timelapse)

In [ ]:

m = geemap.Map()
roi = ee.Geometry.BBox(113.8252, 22.1988, 114.0851, 22.3497)
m.add_layer(roi)
m.center_object(roi)
m

In [ ]:

timelapse = geemap.landsat_timelapse(
    roi,
    out_gif="hong_kong.gif",
    start_year=1990,
    end_year=2022,
    start_date="01-01",
    end_date="12-31",
    bands=["SWIR1", "NIR", "Red"],
    frames_per_second=3,
    title="Hong Kong",
)
geemap.show_image(timelapse)

In [ ]:

m = geemap.Map()
roi = ee.Geometry.BBox(-115.5541, 35.8044, -113.9035, 36.5581)
m.add_layer(roi)
m.center_object(roi)
m

In [ ]:

timelapse = geemap.landsat_timelapse(
    roi,
    out_gif="las_vegas.gif",
    start_year=1984,
    end_year=2025,
    bands=["NIR", "Red", "Green"],
    frames_per_second=5,
    title="Las Vegas, NV",
    font_color="blue",
)
geemap.show_image(timelapse)

MODIS

MODIS vegetation indices

In [ ]:

m = geemap.Map()
m

In [ ]:

roi = m.user_roi
if roi is None:
    roi = ee.Geometry.BBox(-18.6983, -36.1630, 52.2293, 38.1446)
    m.add_layer(roi)
    m.center_object(roi)

In [ ]:

timelapse = geemap.modis_ndvi_timelapse(
    roi,
    out_gif="ndvi.gif",
    data="Terra",
    band="NDVI",
    start_date="2000-01-01",
    end_date="2025-12-31",
    frames_per_second=3,
    title="MODIS NDVI Timelapse",
    overlay_data="countries",
)
geemap.show_image(timelapse)

MODIS temperature

In [ ]:

m = geemap.Map()
m

In [ ]:

roi = m.user_roi
if roi is None:
    roi = ee.Geometry.BBox(-171.21, -57.13, 177.53, 79.99)
    m.add_layer(roi)
    m.center_object(roi)

In [ ]:

timelapse = geemap.modis_ocean_color_timelapse(
    satellite="Aqua",
    start_date="2018-01-01",
    end_date="2020-12-31",
    roi=roi,
    frequency="month",
    out_gif="temperature.gif",
    overlay_data="continents",
    overlay_color="yellow",
    overlay_opacity=0.5,
)
geemap.show_image(timelapse)

GOES

In [ ]:

roi = ee.Geometry.BBox(167.1898, -28.5757, 202.6258, -12.4411)
start_date = "2022-01-15T03:00:00"
end_date = "2022-01-15T07:00:00"
data = "GOES-17"
scan = "full_disk"

In [ ]:

timelapse = geemap.goes_timelapse(
    roi, "goes.gif", start_date, end_date, data, scan, framesPerSecond=5
)
geemap.show_image(timelapse)

In [ ]:

roi = ee.Geometry.BBox(-159.5954, 24.5178, -114.2438, 60.4088)
start_date = "2021-10-24T14:00:00"
end_date = "2021-10-25T01:00:00"
data = "GOES-17"
scan = "full_disk"

In [ ]:

timelapse = geemap.goes_timelapse(
    roi, "hurricane.gif", start_date, end_date, data, scan, framesPerSecond=5
)
geemap.show_image(timelapse)

In [ ]:

roi = ee.Geometry.BBox(-121.0034, 36.8488, -117.9052, 39.0490)
start_date = "2020-09-05T15:00:00"
end_date = "2020-09-06T02:00:00"
data = "GOES-17"
scan = "full_disk"

In [ ]:

timelapse = geemap.goes_fire_timelapse(
    roi, "fire.gif", start_date, end_date, data, scan, framesPerSecond=5
)
geemap.show_image(timelapse)

NAIP

In [ ]:

m = geemap.Map(center=[40, -100], zoom=4)
m

In [ ]:

roi = m.user_roi
if roi is None:
    roi = ee.Geometry.BBox(-99.1019, 47.1274, -99.0334, 47.1562)
    m.add_layer(roi)
    m.center_object(roi)

In [ ]:

timelapse = geemap.naip_timelapse(
    roi,
    out_gif="naip.gif",
    bands=["N", "R", "G"],
    frames_per_second=3,
    title="NAIP Timelapse",
)
geemap.show_image(timelapse)

Sentinel-1

In [ ]:

m = geemap.Map()
m

In [ ]:

roi = m.user_roi
if roi is None:
    roi = ee.Geometry.BBox(117.1132, 3.5227, 117.2214, 3.5843)
    m.add_layer(roi)
    m.center_object(roi)

In [ ]:

timelapse = geemap.sentinel1_timelapse(
    roi,
    out_gif="sentinel1.gif",
    start_year=2019,
    end_year=2019,
    start_date="04-01",
    end_date="08-01",
    bands=["VV"],
    frequency="day",
    vis_params={"min": -30, "max": 0},
    palette="Greys",
    frames_per_second=3,
    title="Sentinel-1 Timelapse",
    add_colorbar=True,
    colorbar_bg_color="gray",
)
geemap.show_image(timelapse)

Exporting Earth Engine data

Exporting images

Add a Landsat image to the map.

In [ ]:

m = geemap.Map()

image = ee.Image("LANDSAT/LC08/C02/T1_TOA/LC08_044034_20140318").select(
    ["B5", "B4", "B3"]
)

vis_params = {"min": 0, "max": 0.5, "gamma": [0.95, 1.1, 1]}

m.center_object(image)
m.add_layer(image, vis_params, "Landsat")
m

Add a rectangle to the map.

In [ ]:

region = ee.Geometry.BBox(-122.5955, 37.5339, -122.0982, 37.8252)
fc = ee.FeatureCollection(region)
style = {"color": "ffff00ff", "fillColor": "00000000"}
m.add_layer(fc.style(**style), {}, "ROI")

To local drive

In [ ]:

geemap.ee_export_image(image, filename="landsat.tif", scale=30, region=region)

Check image projection.

In [ ]:

projection = image.select(0).projection().getInfo()
projection

In [ ]:

crs = projection["crs"]
crs_transform = projection["transform"]

Specify region, crs, and crs_transform.

In [ ]:

geemap.ee_export_image(
    image,
    filename="landsat_crs.tif",
    crs=crs,
    crs_transform=crs_transform,
    region=region,
)

To Google Drive

In [ ]:

geemap.ee_export_image_to_drive(
    image, description="landsat", folder="export", region=region, scale=30
)

In [ ]:

geemap.download_ee_image(image, "landsat.tif", scale=90)

Exporting image collections

In [ ]:

point = ee.Geometry.Point(-99.2222, 46.7816)
collection = (
    ee.ImageCollection("USDA/NAIP/DOQQ")
    .filterBounds(point)
    .filterDate("2008-01-01", "2018-01-01")
    .filter(ee.Filter.listContains("system:band_names", "N"))
)

In [ ]:

collection.aggregate_array("system:index")

To local drive

In [ ]:

geemap.ee_export_image_collection(collection, out_dir=".", scale=10)

To Google Drive

In [ ]:

geemap.ee_export_image_collection_to_drive(collection, folder="export", scale=10)

Exporting feature collections

In [ ]:

m = geemap.Map()
states = ee.FeatureCollection("TIGER/2018/States")
fc = states.filter(ee.Filter.eq("NAME", "Alaska"))
m.add_layer(fc, {}, "Alaska")
m.center_object(fc, 4)
m

To local drive

In [ ]:

geemap.ee_to_shp(fc, filename="Alaska.shp")

In [ ]:

geemap.ee_export_vector(fc, filename="Alaska.shp")

In [ ]:

geemap.ee_to_geojson(fc, filename="Alaska.geojson")

In [ ]:

geemap.ee_to_csv(fc, filename="Alaska.csv")

In [ ]:

gdf = geemap.ee_to_gdf(fc)
gdf

In [ ]:

df = geemap.ee_to_df(fc)
df

To Google Drive

In [ ]:

geemap.ee_export_vector_to_drive(
    fc, description="Alaska", fileFormat="SHP", folder="export"
)

Charts

Feature and FeatureCollection

Import libraries

In [ ]:

import calendar
from geemap import chart

feature_by_feature

Features are plotted along the x-axis, labeled by values of a selected property. Series are represented by adjacent columns defined by a list of property names whose values are plotted along the y-axis.

In [ ]:

ecoregions = ee.FeatureCollection("projects/google/charts_feature_example")
features = ecoregions.select("[0-9][0-9]_tmean|label")

In [ ]:

geemap.ee_to_df(features)

In [ ]:

x_property = "label"
y_properties = [str(x).zfill(2) + "_tmean" for x in range(1, 13)]

labels = calendar.month_abbr[1:]  # a list of month labels, e.g. ['Jan', 'Feb', ...]

colors = [
    "#604791",
    "#1d6b99",
    "#39a8a7",
    "#0f8755",
    "#76b349",
    "#f0af07",
    "#e37d05",
    "#cf513e",
    "#96356f",
    "#724173",
    "#9c4f97",
    "#696969",
]
title = "Average Monthly Temperature by Ecoregion"
x_label = "Ecoregion"
y_label = "Temperature"

In [ ]:

fig = chart.feature_by_feature(
    features,
    x_property,
    y_properties,
    colors=colors,
    labels=labels,
    title=title,
    x_label=x_label,
    y_label=y_label,
)
fig

feature.by_property

In [ ]:

ecoregions = ee.FeatureCollection("projects/google/charts_feature_example")
features = ecoregions.select("[0-9][0-9]_ppt|label")

In [ ]:

geemap.ee_to_df(features)

In [ ]:

keys = [str(x).zfill(2) + "_ppt" for x in range(1, 13)]
values = calendar.month_abbr[1:]  # a list of month labels, e.g. ['Jan', 'Feb', ...]

In [ ]:

x_properties = dict(zip(keys, values))
series_property = "label"
title = "Average Ecoregion Precipitation by Month"
colors = ["#f0af07", "#0f8755", "#76b349"]

In [ ]:

fig = chart.feature_by_property(
    features,
    x_properties,
    series_property,
    title=title,
    colors=colors,
    x_label="Month",
    y_label="Precipitation (mm)",
    legend_location="top-left",
)
fig

feature_groups

In [ ]:

ecoregions = ee.FeatureCollection("projects/google/charts_feature_example")
features = ecoregions.select("[0-9][0-9]_ppt|label")

In [ ]:

features = ee.FeatureCollection("projects/google/charts_feature_example")
x_property = "label"
y_property = "01_tmean"
series_property = "warm"
title = "Average January Temperature by Ecoregion"
colors = ["#cf513e", "#1d6b99"]
labels = ["Warm", "Cold"]

In [ ]:

chart.feature_groups(
    features,
    x_property,
    y_property,
    series_property,
    title=title,
    colors=colors,
    x_label="Ecoregion",
    y_label="January Temperature (°C)",
    legend_location="top-right",
    labels=labels,
)

feature_histogram

In [ ]:

source = ee.ImageCollection("OREGONSTATE/PRISM/Norm91m").toBands()
region = ee.Geometry.Rectangle(-123.41, 40.43, -116.38, 45.14)
features = source.sample(region, 5000)

In [ ]:

geemap.ee_to_df(features.limit(5).select(["07_ppt"]))

In [ ]:

property = "07_ppt"
title = "July Precipitation Distribution for NW USA"

In [ ]:

fig = chart.feature_histogram(
    features,
    property,
    max_buckets=None,
    title=title,
    x_label="Precipitation (mm)",
    y_label="Pixel Count",
    colors=["#1d6b99"],
)
fig

Image charts

image_by_region

In [ ]:

ecoregions = ee.FeatureCollection("projects/google/charts_feature_example")
image = (
    ee.ImageCollection("OREGONSTATE/PRISM/Norm91m").toBands().select("[0-9][0-9]_tmean")
)

In [ ]:

labels = calendar.month_abbr[1:]  # a list of month labels, e.g. ['Jan', 'Feb', ...]
title = "Average Monthly Temperature by Ecoregion"

In [ ]:

fig = chart.image_by_region(
    image,
    ecoregions,
    reducer="mean",
    scale=500,
    x_property="label",
    title=title,
    x_label="Ecoregion",
    y_label="Temperature",
    labels=labels,
)
fig

image_regions

In [ ]:

ecoregions = ee.FeatureCollection("projects/google/charts_feature_example")
image = (
    ee.ImageCollection("OREGONSTATE/PRISM/Norm91m").toBands().select("[0-9][0-9]_ppt")
)

In [ ]:

keys = [str(x).zfill(2) + "_ppt" for x in range(1, 13)]
values = calendar.month_abbr[1:]  # a list of month labels, e.g. ['Jan', 'Feb', ...]

In [ ]:

x_properties = dict(zip(keys, values))
title = "Average Ecoregion Precipitation by Month"
colors = ["#f0af07", "#0f8755", "#76b349"]

In [ ]:

fig = chart.image_regions(
    image,
    ecoregions,
    reducer="mean",
    scale=500,
    series_property="label",
    x_labels=x_properties,
    title=title,
    colors=colors,
    x_label="Month",
    y_label="Precipitation (mm)",
    legend_location="top-left",
)

image_by_class

In [ ]:

ecoregions = ee.FeatureCollection("projects/google/charts_feature_example")

image = (
    ee.ImageCollection("MODIS/061/MOD09A1")
    .filter(ee.Filter.date("2018-06-01", "2018-09-01"))
    .select("sur_refl_b0[0-7]")
    .mean()
    .select([2, 3, 0, 1, 4, 5, 6])
)

wavelengths = [469, 555, 655, 858, 1240, 1640, 2130]

In [ ]:

fig = chart.image_by_class(
    image,
    class_band="label",
    region=ecoregions,
    reducer="MEAN",
    scale=500,
    x_labels=wavelengths,
    title="Ecoregion Spectral Signatures",
    x_label="Wavelength (nm)",
    y_label="Reflectance (x1e4)",
    colors=["#f0af07", "#0f8755", "#76b349"],
    legend_location="top-left",
    interpolation="basis",
)
fig

image_histogram

In [ ]:

image = (
    ee.ImageCollection("MODIS/061/MOD09A1")
    .filter(ee.Filter.date("2018-06-01", "2018-09-01"))
    .select(["sur_refl_b01", "sur_refl_b02", "sur_refl_b06"])
    .mean()
)

region = ee.Geometry.Rectangle([-112.60, 40.60, -111.18, 41.22])

In [ ]:

fig = chart.image_histogram(
    image,
    region,
    scale=500,
    max_buckets=200,
    min_bucket_width=1.0,
    max_raw=1000,
    max_pixels=int(1e6),
    title="MODIS SR Reflectance Histogram",
    labels=["Red", "NIR", "SWIR"],
    colors=["#cf513e", "#1d6b99", "#f0af07"],
)
fig

ImageCollection charts

image_series

In [ ]:

# Define the forest feature collection.
forest = ee.FeatureCollection("projects/google/charts_feature_example").filter(
    ee.Filter.eq("label", "Forest")
)

# Load MODIS vegetation indices data and subset a decade of images.
veg_indices = (
    ee.ImageCollection("MODIS/061/MOD13A1")
    .filter(ee.Filter.date("2010-01-01", "2020-01-01"))
    .select(["NDVI", "EVI"])
)

In [ ]:

title = "Average Vegetation Index Value by Date for Forest"
x_label = "Year"
y_label = "Vegetation index (x1e4)"
colors = ["#e37d05", "#1d6b99"]

In [ ]:

fig = chart.image_series(
    veg_indices,
    region=forest,
    reducer=ee.Reducer.mean(),
    scale=500,
    x_property="system:time_start",
    chart_type="LineChart",
    x_cols="date",
    y_cols=["NDVI", "EVI"],
    colors=colors,
    title=title,
    x_label=x_label,
    y_label=y_label,
    legend_location="right",
)
fig

image_series_by_region

In [ ]:

# Import the example feature collection.
ecoregions = ee.FeatureCollection("projects/google/charts_feature_example")

# Load MODIS vegetation indices data and subset a decade of images.
veg_indices = (
    ee.ImageCollection("MODIS/061/MOD13A1")
    .filter(ee.Filter.date("2010-01-01", "2020-01-01"))
    .select(["NDVI"])
)

In [ ]:

title = "Average NDVI Value by Date"
x_label = "Date"
y_label = "NDVI (x1e4)"
x_cols = "index"
y_cols = ["Desert", "Forest", "Grassland"]
colors = ["#f0af07", "#0f8755", "#76b349"]

In [ ]:

fig = chart.image_series_by_region(
    veg_indices,
    regions=ecoregions,
    reducer=ee.Reducer.mean(),
    band="NDVI",
    scale=500,
    x_property="system:time_start",
    series_property="label",
    chart_type="LineChart",
    x_cols=x_cols,
    y_cols=y_cols,
    title=title,
    x_label=x_label,
    y_label=y_label,
    colors=colors,
    stroke_width=3,
    legend_location="bottom-left",
)
fig

image_doy_series

In [ ]:

# Import the example feature collection and subset the grassland feature.
grassland = ee.FeatureCollection("projects/google/charts_feature_example").filter(
    ee.Filter.eq("label", "Grassland")
)

# Load MODIS vegetation indices data and subset a decade of images.
veg_indices = (
    ee.ImageCollection("MODIS/061/MOD13A1")
    .filter(ee.Filter.date("2010-01-01", "2020-01-01"))
    .select(["NDVI", "EVI"])
)

In [ ]:

title = "Average Vegetation Index Value by Day of Year for Grassland"
x_label = "Day of Year"
y_label = "Vegetation Index (x1e4)"
colors = ["#f0af07", "#0f8755"]

In [ ]:

fig = chart.image_doy_series(
    image_collection=veg_indices,
    region=grassland,
    scale=500,
    chart_type="LineChart",
    title=title,
    x_label=x_label,
    y_label=y_label,
    colors=colors,
    stroke_width=5,
)
fig

image_doy_series_by_year

In [ ]:

# Import the example feature collection and subset the grassland feature.
grassland = ee.FeatureCollection("projects/google/charts_feature_example").filter(
    ee.Filter.eq("label", "Grassland")
)

# Load MODIS vegetation indices data and subset years 2012 and 2019.
veg_indices = (
    ee.ImageCollection("MODIS/061/MOD13A1")
    .filter(
        ee.Filter.Or(
            ee.Filter.date("2012-01-01", "2013-01-01"),
            ee.Filter.date("2019-01-01", "2020-01-01"),
        )
    )
    .select(["NDVI", "EVI"])
)

In [ ]:

title = "Average Vegetation Index Value by Day of Year for Grassland"
x_label = "Day of Year"
y_label = "Vegetation Index (x1e4)"
colors = ["#e37d05", "#1d6b99"]

In [ ]:

fig = chart.doy_series_by_year(
    veg_indices,
    band_name="NDVI",
    region=grassland,
    scale=500,
    chart_type="LineChart",
    colors=colors,
    title=title,
    x_label=x_label,
    y_label=y_label,
    stroke_width=5,
)
fig

image_doy_series_by_region

In [ ]:

# Import the example feature collection and subset the grassland feature.
ecoregions = ee.FeatureCollection("projects/google/charts_feature_example")

# Load MODIS vegetation indices data and subset a decade of images.
veg_indices = (
    ee.ImageCollection("MODIS/061/MOD13A1")
    .filter(ee.Filter.date("2010-01-01", "2020-01-01"))
    .select(["NDVI"])
)

In [ ]:

title = "Average Vegetation Index Value by Day of Year for Grassland"
x_label = "Day of Year"
y_label = "Vegetation Index (x1e4)"
colors = ["#f0af07", "#0f8755", "#76b349"]

In [ ]:

fig = chart.image_doy_series_by_region(
    veg_indices,
    "NDVI",
    ecoregions,
    region_reducer="mean",
    scale=500,
    year_reducer=ee.Reducer.mean(),
    start_day=1,
    end_day=366,
    series_property="label",
    stroke_width=5,
    chart_type="LineChart",
    title=title,
    x_label=x_label,
    y_label=y_label,
    colors=colors,
    legend_location="right",
)
fig

Array and list charts

Scatter plot

In [ ]:

# Import the example feature collection and subset the forest feature.
forest = ee.FeatureCollection("projects/google/charts_feature_example").filter(
    ee.Filter.eq("label", "Forest")
)

# Define a MODIS surface reflectance composite.
modisSr = (
    ee.ImageCollection("MODIS/061/MOD09A1")
    .filter(ee.Filter.date("2018-06-01", "2018-09-01"))
    .select("sur_refl_b0[0-7]")
    .mean()
)

# Reduce MODIS reflectance bands by forest region; get a dictionary with
# band names as keys, pixel values as lists.
pixel_vals = modisSr.reduceRegion(
    **{"reducer": ee.Reducer.toList(), "geometry": forest.geometry(), "scale": 2000}
)

# Convert NIR and SWIR value lists to an array to be plotted along the y-axis.
y_values = pixel_vals.toArray(["sur_refl_b02", "sur_refl_b06"])


# Get the red band value list; to be plotted along the x-axis.
x_values = ee.List(pixel_vals.get("sur_refl_b01"))

In [ ]:

title = "Relationship Among Spectral Bands for Forest Pixels"
colors = ["rgba(29,107,153,0.4)", "rgba(207,81,62,0.4)"]

In [ ]:

fig = chart.array_values(
    y_values,
    axis=1,
    x_labels=x_values,
    series_names=["NIR", "SWIR"],
    chart_type="ScatterChart",
    colors=colors,
    title=title,
    x_label="Red reflectance (x1e4)",
    y_label="NIR & SWIR reflectance (x1e4)",
    default_size=15,
    xlim=(0, 800),
)
fig

In [ ]:

x = ee.List(pixel_vals.get("sur_refl_b01"))
y = ee.List(pixel_vals.get("sur_refl_b06"))

In [ ]:

fig = chart.array_values(
    y,
    x_labels=x,
    series_names=["SWIR"],
    chart_type="ScatterChart",
    colors=["rgba(207,81,62,0.4)"],
    title=title,
    x_label="Red reflectance (x1e4)",
    y_label="SWIR reflectance (x1e4)",
    default_size=15,
    xlim=(0, 800),
)
fig

Transect line plot

In [ ]:

# Define a line across the Olympic Peninsula, USA.
transect = ee.Geometry.LineString([[-122.8, 47.8], [-124.5, 47.8]])

# Define a pixel coordinate image.
lat_lon_img = ee.Image.pixelLonLat()

# Import a digital surface model and add latitude and longitude bands.
elev_img = ee.Image("USGS/SRTMGL1_003").select("elevation").addBands(lat_lon_img)

# Reduce elevation and coordinate bands by transect line; get a dictionary with
# band names as keys, pixel values as lists.
elev_transect = elev_img.reduceRegion(
    reducer=ee.Reducer.toList(),
    geometry=transect,
    scale=1000,
)

# Get longitude and elevation value lists from the reduction dictionary.
lon = ee.List(elev_transect.get("longitude"))
elev = ee.List(elev_transect.get("elevation"))

# Sort the longitude and elevation values by ascending longitude.
lon_sort = lon.sort(lon)
elev_sort = elev.sort(lon)

In [ ]:

fig = chart.array_values(
    elev_sort,
    x_labels=lon_sort,
    series_names=["Elevation"],
    chart_type="AreaChart",
    colors=["#1d6b99"],
    title="Elevation Profile Across Longitude",
    x_label="Longitude",
    y_label="Elevation (m)",
    stroke_width=5,
    fill="bottom",
    fill_opacities=[0.4],
    ylim=(0, 2500),
)
fig

Metadata scatter plot

In [ ]:

# Import a Landsat 8 collection and filter to a single path/row.
col = ee.ImageCollection("LANDSAT/LC08/C02/T1_L2").filter(
    ee.Filter.expression("WRS_PATH ==  45 && WRS_ROW == 30")
)

# Reduce image properties to a series of lists; one for each selected property.
propVals = col.reduceColumns(
    reducer=ee.Reducer.toList().repeat(2),
    selectors=["CLOUD_COVER", "GEOMETRIC_RMSE_MODEL"],
).get("list")

# Get selected image property value lists; to be plotted along x and y axes.
x = ee.List(ee.List(propVals).get(0))
y = ee.List(ee.List(propVals).get(1))

In [ ]:

colors = [geemap.hex_to_rgba("#96356f", 0.4)]
print(colors)

In [ ]:

fig = chart.array_values(
    y,
    x_labels=x,
    series_names=["RMSE"],
    chart_type="ScatterChart",
    colors=colors,
    title="Landsat 8 Image Collection Metadata (045030)",
    x_label="Cloud cover (%)",
    y_label="Geometric RMSE (m)",
    default_size=15,
)
fig

Mapped function scatter & line plot

In [ ]:

import math

start = -2 * math.pi
end = 2 * math.pi
points = ee.List.sequence(start, end, None, 50)


def sin_func(val):
    return ee.Number(val).sin()


values = points.map(sin_func)

In [ ]:

fig = chart.array_values(
    values,
    points,
    chart_type="LineChart",
    colors=["#39a8a7"],
    title="Sine Function",
    x_label="radians",
    y_label="sin(x)",
    marker="circle",
)
fig

Data table charts

In [ ]:

import pandas as pd

Manual DataTable chart

In [ ]:

data = {
    "State": ["CA", "NY", "IL", "MI", "OR"],
    "Population": [37253956, 19378102, 12830632, 9883640, 3831074],
}

df = pd.DataFrame(data)
df

In [ ]:

fig = chart.Chart(
    df,
    x_cols=["State"],
    y_cols=["Population"],
    chart_type="ColumnChart",
    colors=["#1d6b99"],
    title="State Population (US census, 2010)",
    x_label="State",
    y_label="Population",
)
fig

Computed DataTable chart

In [ ]:

# Import the example feature collection and subset the forest feature.
forest = ee.FeatureCollection("projects/google/charts_feature_example").filter(
    ee.Filter.eq("label", "Forest")
)

# Load MODIS vegetation indices data and subset a decade of images.
veg_indices = (
    ee.ImageCollection("MODIS/061/MOD13A1")
    .filter(ee.Filter.date("2010-01-01", "2020-01-01"))
    .select(["NDVI", "EVI"])
)

# Build a feature collection where each feature has a property that represents
# a DataFrame row.


def aggregate(img):
    # Reduce the image to the mean of pixels intersecting the forest ecoregion.
    stat = img.reduceRegion(
        **{"reducer": ee.Reducer.mean(), "geometry": forest, "scale": 500}
    )

    # Extract the reduction results along with the image date.
    date = geemap.image_date(img)
    evi = stat.get("EVI")
    ndvi = stat.get("NDVI")

    # Make a list of observation attributes to define a row in the DataTable.
    row = ee.List([date, evi, ndvi])

    # Return the row as a property of an ee.Feature.
    return ee.Feature(None, {"row": row})


reduction_table = veg_indices.map(aggregate)

# Aggregate the 'row' property from all features in the new feature collection
# to make a server-side 2-D list (DataTable).
data_table_server = reduction_table.aggregate_array("row")

# Define column names and properties for the DataTable. The order should
# correspond to the order in the construction of the 'row' property above.
column_header = ee.List([["Date", "EVI", "NDVI"]])

# Concatenate the column header to the table.
data_table_server = column_header.cat(data_table_server)

In [ ]:

data_table = chart.DataTable(data_table_server, date_column="Date")
data_table.head()

In [ ]:

fig = chart.Chart(
    data_table,
    chart_type="LineChart",
    x_cols="Date",
    y_cols=["EVI", "NDVI"],
    colors=["#e37d05", "#1d6b99"],
    title="Average Vegetation Index Value by Date for Forest",
    x_label="Date",
    y_label="Vegetation index (x1e4)",
    stroke_width=3,
    legend_location="right",
)
fig

Interval chart

In [ ]:

# Define a point to extract an NDVI time series for.
geometry = ee.Geometry.Point([-121.679, 36.479])

# Define a band of interest (NDVI), import the MODIS vegetation index dataset,
# and select the band.
band = "NDVI"
ndvi_col = ee.ImageCollection("MODIS/061/MOD13Q1").select(band)

# Map over the collection to add a day of year (doy) property to each image.


def set_doy(img):
    doy = ee.Date(img.get("system:time_start")).getRelative("day", "year")
    # Add 8 to day of year number so that the doy label represents the middle of
    # the 16-day MODIS NDVI composite.
    return img.set("doy", ee.Number(doy).add(8))


ndvi_col = ndvi_col.map(set_doy)

# Join all coincident day of year observations into a set of image collections.
distinct_doy = ndvi_col.filterDate("2013-01-01", "2014-01-01")
filter = ee.Filter.equals(**{"leftField": "doy", "rightField": "doy"})
join = ee.Join.saveAll("doy_matches")
join_col = ee.ImageCollection(join.apply(distinct_doy, ndvi_col, filter))

# Calculate the absolute range, interquartile range, and median for the set
# of images composing each coincident doy observation group. The result is
# an image collection with an image representative per unique doy observation
# with bands that describe the 0, 25, 50, 75, 100 percentiles for the set of
# coincident doy images.


def cal_percentiles(img):
    doyCol = ee.ImageCollection.fromImages(img.get("doy_matches"))

    return doyCol.reduce(
        ee.Reducer.percentile([0, 25, 50, 75, 100], ["p0", "p25", "p50", "p75", "p100"])
    ).set({"doy": img.get("doy")})


comp = ee.ImageCollection(join_col.map(cal_percentiles))

# Extract the inter-annual NDVI doy percentile statistics for the
# point of interest per unique doy representative. The result is
# is a feature collection where each feature is a doy representative that
# contains a property (row) describing the respective inter-annual NDVI
# variance, formatted as a list of values.


def order_percentiles(img):
    stats = ee.Dictionary(
        img.reduceRegion(
            **{"reducer": ee.Reducer.first(), "geometry": geometry, "scale": 250}
        )
    )

    # Order the percentile reduction elements according to how you want columns
    # in the DataTable arranged (x-axis values need to be first).
    row = ee.List(
        [
            img.get("doy"),
            stats.get(band + "_p50"),
            stats.get(band + "_p0"),
            stats.get(band + "_p25"),
            stats.get(band + "_p75"),
            stats.get(band + "_p100"),
        ]
    )

    # Return the row as a property of an ee.Feature.
    return ee.Feature(None, {"row": row})


reduction_table = comp.map(order_percentiles)

# Aggregate the 'row' properties to make a server-side 2-D array (DataTable).
data_table_server = reduction_table.aggregate_array("row")

# Define column names and properties for the DataTable. The order should
# correspond to the order in the construction of the 'row' property above.
column_header = ee.List([["DOY", "median", "p0", "p25", "p75", "p100"]])

# Concatenate the column header to the table.
data_table_server = column_header.cat(data_table_server)

In [ ]:

df = chart.DataTable(data_table_server)
df.head()

In [ ]:

fig = chart.Chart(
    df,
    chart_type="IntervalChart",
    x_cols="DOY",
    y_cols=["p0", "p25", "median", "p75", "p100"],
    title="Annual NDVI Time Series with Inter-Annual Variance",
    x_label="Day of Year",
    y_label="Vegetation index (x1e4)",
    stroke_width=1,
    fill="between",
    fill_colors=["#b6d1c6", "#83b191", "#83b191", "#b6d1c6"],
    fill_opacities=[0.6] * 4,
    labels=["p0", "p25", "median", "p75", "p100"],
    display_legend=True,
    legend_location="top-right",
    ylim=(0, 10000),
)
fig