Alaska 2024 Part2
Geospatial Cloud Computing with the GEE Python API - Part 2
- Notebook: https://geemap.org/workshops/Alaska_2024_Part2
- Earth Engine: https://earthengine.google.com
- Geemap: https://geemap.org
Introduction¶
This notebook contains the materials for the first part of the workshop Geospatial Cloud Computing with the GEE Python API at the University of Alaska Fairbanks.
This workshop provides an introduction to cloud-based geospatial analysis using the Earth Engine Python API. Attendees will learn the basics of Earth Engine data types and how to visualize, analyze, and export Earth Engine data in a Jupyter environment with geemap. In addition, attendees will learn how to develop and deploy interactive Earth Engine web apps with Python. Through practical examples and hands-on exercises, attendees will enhance their learning experience. During each hands-on session, attendees will walk through Jupyter Notebook examples on Google Colab with the instructors. At the end of each session, they will complete a hands-on exercise to apply the knowledge they have learned.
Agenda¶
The workshop is divided into three parts. The second part will cover the following topics:
- Processing of vector data (shapefiles, json, conversion from one format to another)
- Processing of raster data: extract pixel value, raster calculator, zonal statistics etc.
- Working with local geospatial data in Geemap
- Accessing Cloud Optimized GeoTIFF
- Exporting EE Image and Feature data
- Creating timelapse animations using Landsat or Sentinel 2 for Alaska
- Time series analysis: Forest cover change for a test site in Alaska (e.g. Bonanza Creek LTER or Caribou-Poker Creeks Research Watershed)
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.11
conda activate gee
conda install -c conda-forge mamba
mamba install -c conda-forge pygis
# %pip install geemap pygis
Import libraries¶
import ee
import geemap
geemap.ee_initialize()
in_geojson = "https://github.com/gee-community/geemap/blob/master/examples/data/countries.geojson"
m = geemap.Map()
fc = geemap.geojson_to_ee(in_geojson)
m.add_layer(fc.style(**{"color": "ff0000", "fillColor": "00000000"}), {}, "Countries")
m
From Shapefile¶
url = "https://github.com/gee-community/geemap/blob/master/examples/data/countries.zip"
geemap.download_file(url, overwrite=True)
in_shp = "countries.shp"
fc = geemap.shp_to_ee(in_shp)
m = geemap.Map()
m.add_layer(fc, {}, "Countries")
m
From GeoDataFrame¶
import geopandas as gpd
gdf = gpd.read_file(in_shp)
fc = geemap.gdf_to_ee(gdf)
m = geemap.Map()
m.add_layer(fc, {}, "Countries")
m
To GeoJSON¶
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
geemap.ee_to_geojson(fc, filename="Alaska.geojson")
To Shapefile¶
geemap.ee_to_shp(fc, filename="Alaska.shp")
To GeoDataFrame¶
gdf = geemap.ee_to_gdf(fc)
gdf
gdf.explore()
To DataFrame¶
df = geemap.ee_to_df(fc)
df
To CSV¶
geemap.ee_to_csv(fc, filename="Alaska.csv")
m = geemap.Map(center=[40, -100], zoom=4)
dem = ee.Image("USGS/SRTMGL1_003")
landsat7 = ee.Image("LANDSAT/LE7_TOA_5YEAR/1999_2003")
vis_params = {
"min": 0,
"max": 4000,
"palette": ["006633", "E5FFCC", "662A00", "D8D8D8", "F5F5F5"],
}
m.add_layer(
landsat7,
{"bands": ["B4", "B3", "B2"], "min": 20, "max": 200, "gamma": 2},
"Landsat 7",
)
m.add_layer(dem, vis_params, "SRTM DEM", True, 1)
m
in_shp = "us_cities.shp"
url = "https://github.com/giswqs/data/raw/main/us/us_cities.zip"
geemap.download_file(url)
in_fc = geemap.shp_to_ee(in_shp)
m.add_layer(in_fc, {}, "Cities")
geemap.extract_values_to_points(in_fc, dem, out_fc="dem.shp")
geemap.shp_to_gdf("dem.shp")
geemap.extract_values_to_points(in_fc, landsat7, "landsat.csv")
geemap.csv_to_df("landsat.csv")
Extracting pixel values along a transect¶
m = geemap.Map(center=[40, -100], zoom=4)
m.add_basemap("TERRAIN")
image = ee.Image("USGS/SRTMGL1_003")
vis_params = {
"min": 0,
"max": 4000,
"palette": ["006633", "E5FFCC", "662A00", "D8D8D8", "F5F5F5"],
}
m.add_layer(image, vis_params, "SRTM DEM", True, 0.5)
m
line = m.user_roi
if line is None:
line = ee.Geometry.LineString(
[[-120.2232, 36.3148], [-118.9269, 36.7121], [-117.2022, 36.7562]]
)
m.add_layer(line, {}, "ROI")
m.centerObject(line)
reducer = "mean"
transect = geemap.extract_transect(
image, line, n_segments=100, reducer=reducer, to_pandas=True
)
transect
geemap.line_chart(
data=transect,
x="distance",
y="mean",
markers=True,
x_label="Distance (m)",
y_label="Elevation (m)",
height=400,
)
transect.to_csv("transect.csv")
m = geemap.Map(center=[40, -100], zoom=4)
# Add NASA SRTM
dem = ee.Image("USGS/SRTMGL1_003")
dem_vis = {
"min": 0,
"max": 4000,
"palette": ["006633", "E5FFCC", "662A00", "D8D8D8", "F5F5F5"],
}
m.add_layer(dem, dem_vis, "SRTM DEM")
# Add 5-year Landsat TOA composite
landsat = ee.Image("LANDSAT/LE7_TOA_5YEAR/1999_2003")
landsat_vis = {"bands": ["B4", "B3", "B2"], "gamma": 1.4}
m.add_layer(landsat, landsat_vis, "Landsat", False)
# Add US Census States
states = ee.FeatureCollection("TIGER/2018/States")
style = {"fillColor": "00000000"}
m.add_layer(states.style(**style), {}, "US States")
m
out_dem_stats = "dem_stats.csv"
geemap.zonal_stats(
dem, states, out_dem_stats, statistics_type="MEAN", scale=1000, return_fc=False
)
out_landsat_stats = "landsat_stats.csv"
geemap.zonal_stats(
landsat,
states,
out_landsat_stats,
statistics_type="MEAN",
scale=1000,
return_fc=False,
)
Zonal statistics by group¶
m = geemap.Map(center=[40, -100], zoom=4)
# Add NLCD data
dataset = ee.Image("USGS/NLCD_RELEASES/2019_REL/NLCD/2019")
landcover = dataset.select("landcover")
m.add_layer(landcover, {}, "NLCD 2019")
# Add US census states
states = ee.FeatureCollection("TIGER/2018/States")
style = {"fillColor": "00000000"}
m.add_layer(states.style(**style), {}, "US States")
# Add NLCD legend
m.add_legend(title="NLCD Land Cover", builtin_legend="NLCD")
m
nlcd_stats = "nlcd_stats.csv"
geemap.zonal_stats_by_group(
landcover,
states,
nlcd_stats,
statistics_type="SUM",
denominator=1e6,
decimal_places=2,
)
nlcd_stats = "nlcd_stats_pct.csv"
geemap.zonal_stats_by_group(
landcover,
states,
nlcd_stats,
statistics_type="PERCENTAGE",
denominator=1e6,
decimal_places=2,
)
Zonal statistics with two images¶
m = geemap.Map(center=[40, -100], zoom=4)
dem = ee.Image("USGS/3DEP/10m")
vis = {"min": 0, "max": 4000, "palette": "terrain"}
m.add_layer(dem, vis, "DEM")
m
landcover = ee.Image("USGS/NLCD_RELEASES/2019_REL/NLCD/2019").select("landcover")
m.add_layer(landcover, {}, "NLCD 2019")
m.add_legend(title="NLCD Land Cover Classification", builtin_legend="NLCD")
stats = geemap.image_stats_by_zone(dem, landcover, reducer="MEAN")
stats
stats.to_csv("mean.csv", index=False)
geemap.image_stats_by_zone(dem, landcover, out_csv="std.csv", reducer="STD")
Map algebra¶
m = geemap.Map()
# Load a 5-year Landsat 7 composite 1999-2003.
landsat_1999 = ee.Image("LANDSAT/LE7_TOA_5YEAR/1999_2003")
# Compute NDVI.
ndvi_1999 = (
landsat_1999.select("B4")
.subtract(landsat_1999.select("B3"))
.divide(landsat_1999.select("B4").add(landsat_1999.select("B3")))
)
vis = {"min": 0, "max": 1, "palette": "ndvi"}
m.add_layer(ndvi_1999, vis, "NDVI")
m.add_colorbar(vis, label="NDVI")
m
# Load a Landsat 8 image.
image = ee.Image("LANDSAT/LC08/C02/T1_TOA/LC08_044034_20140318")
# Compute the EVI using an expression.
evi = image.expression(
"2.5 * ((NIR - RED) / (NIR + 6 * RED - 7.5 * BLUE + 1))",
{
"NIR": image.select("B5"),
"RED": image.select("B4"),
"BLUE": image.select("B2"),
},
)
# Define a map centered on San Francisco Bay.
m = geemap.Map(center=[37.4675, -122.1363], zoom=9)
vis = {"min": 0, "max": 1, "palette": "ndvi"}
m.add_layer(evi, vis, "EVI")
m.add_colorbar(vis, label="EVI")
m
Exercise 1 - Zonal statistics¶
Find out which state has the highest mean temperature in the United States on June 28, 2023. Relevant Earth Engine assets:
url = "https://github.com/giswqs/data/raw/main/raster/srtm90.tif"
filename = "dem.tif"
geemap.download_file(url, filename)
Multi-band raster¶
m = geemap.Map()
m.add_raster(filename, cmap="terrain", layer_name="DEM")
vis_params = {"min": 0, "max": 4000, "palette": "terrain"}
m.add_colorbar(vis_params, label="Elevation (m)")
m
url = "https://github.com/giswqs/data/raw/main/raster/cog.tif"
filename = "cog.tif"
geemap.download_file(url, filename)
m = geemap.Map()
m.add_raster(filename, indexes=[4, 1, 2], layer_name="False color")
m
in_geojson = (
"https://github.com/opengeos/datasets/releases/download/vector/cables.geojson"
)
m = geemap.Map()
m.add_geojson(in_geojson, layer_name="Cable lines", info_mode="on_hover")
m
m = geemap.Map()
m.add_basemap("CartoDB.DarkMatter")
callback = lambda feat: {"color": feat["properties"]["color"], "weight": 2}
m.add_geojson(in_geojson, layer_name="Cable lines", style_callback=callback)
m
url = "https://github.com/opengeos/datasets/releases/download/world/countries.geojson"
m = geemap.Map()
m.add_geojson(
url, layer_name="Countries", fill_colors=["red", "yellow", "green", "orange"]
)
m
import random
m = geemap.Map()
def random_color(feature):
return {
"color": "black",
"weight": 3,
"fillColor": random.choice(["red", "yellow", "green", "orange"]),
}
m.add_geojson(url, layer_name="Countries", style_callback=random_color)
m
m = geemap.Map()
style = {
"stroke": True,
"color": "#0000ff",
"weight": 2,
"opacity": 1,
"fill": True,
"fillColor": "#0000ff",
"fillOpacity": 0.1,
}
hover_style = {"fillOpacity": 0.7}
m.add_geojson(url, layer_name="Countries", style=style, hover_style=hover_style)
m
Shapefile¶
url = "https://github.com/opengeos/datasets/releases/download/world/countries.zip"
geemap.download_file(url, overwrite=True)
m = geemap.Map()
in_shp = "countries.shp"
m.add_shp(in_shp, layer_name="Countries")
m
GeoDataFrame¶
import geopandas as gpd
m = geemap.Map(center=[40, -100], zoom=4)
gdf = gpd.read_file("countries.shp")
m.add_gdf(gdf, layer_name="Countries")
m
GeoPackage¶
m = geemap.Map()
data = "https://github.com/opengeos/datasets/releases/download/world/countries.gpkg"
m.add_vector(data, layer_name="Countries")
m
CSV to vector¶
data = "https://github.com/gee-community/geemap/blob/master/examples/data/us_cities.csv"
geemap.csv_to_df(data)
geemap.csv_to_geojson(
data, "cities.geojson", latitude="latitude", longitude="longitude"
)
geemap.csv_to_shp(data, "cities.shp", latitude="latitude", longitude="longitude")
geemap.csv_to_vector(data, "cities.gpkg", latitude="latitude", longitude="longitude")
gdf = geemap.csv_to_gdf(data, latitude="latitude", longitude="longitude")
gdf
cities = (
"https://github.com/gee-community/geemap/blob/master/examples/data/us_cities.csv"
)
m = geemap.Map(center=[40, -100], zoom=4)
m.add_points_from_xy(cities, x="longitude", y="latitude")
m
regions = "https://github.com/gee-community/geemap/blob/master/examples/data/us_regions.geojson"
m = geemap.Map(center=[40, -100], zoom=4)
m.add_geojson(regions, layer_name="US Regions")
m.add_points_from_xy(
cities,
x="longitude",
y="latitude",
layer_name="US Cities",
color_column="region",
icon_names=["gear", "map", "leaf", "globe"],
spin=True,
add_legend=True,
)
m
m = geemap.Map(center=[40, -100], zoom=4)
m.add_circle_markers_from_xy(
data,
x="longitude",
y="latitude",
radius=8,
color="blue",
fill_color="black",
fill_opacity=0.5,
)
m
m = geemap.Map()
url = "https://opendata.digitalglobe.com/events/california-fire-2020/pre-event/2018-02-16/pine-gulch-fire20/1030010076004E00.tif"
m.add_cog_layer(url, name="Fire (pre-event)")
m
geemap.cog_center(url)
geemap.cog_bands(url)
url2 = "https://opendata.digitalglobe.com/events/california-fire-2020/post-event/2020-08-14/pine-gulch-fire20/10300100AAC8DD00.tif"
m.add_cog_layer(url2, name="Fire (post-event)")
m = geemap.Map(center=[39.4568, -108.5107], zoom=12)
m.split_map(left_layer=url2, right_layer=url)
m
SpatioTemporal Asset Catalog (STAC)¶
url = "https://tinyurl.com/22vptbws"
geemap.stac_bounds(url)
geemap.stac_center(url)
geemap.stac_bands(url)
m = geemap.Map()
m.add_stac_layer(url, bands=["pan"], name="Panchromatic")
m.add_stac_layer(url, bands=["B3", "B2", "B1"], name="False color")
m
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.
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
geemap.ee_export_image(image, filename="landsat.tif", scale=30, region=region)
Check image projection.
projection = image.select(0).projection().getInfo()
projection
crs = projection["crs"]
crs_transform = projection["transform"]
Specify region, crs, and crs_transform.
geemap.ee_export_image(
image,
filename="landsat_crs.tif",
crs=crs,
crs_transform=crs_transform,
region=region,
)
To Google Drive
geemap.ee_export_image_to_drive(
image, description="landsat", folder="export", region=region, scale=30
)
geemap.download_ee_image(image, "landsat.tif", scale=90)
Exporting image collections¶
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"))
)
collection.aggregate_array("system:index")
To local drive
geemap.ee_export_image_collection(collection, out_dir=".", scale=10)
To Google Drive
geemap.ee_export_image_collection_to_drive(collection, folder="export", scale=10)
Exporting feature collections¶
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
geemap.ee_to_shp(fc, filename="Alaska.shp")
geemap.ee_export_vector(fc, filename="Alaska.shp")
geemap.ee_to_geojson(fc, filename="Alaska.geojson")
geemap.ee_to_csv(fc, filename="Alaska.csv")
gdf = geemap.ee_to_gdf(fc)
gdf
df = geemap.ee_to_df(fc)
df
To Google Drive
geemap.ee_export_vector_to_drive(
fc, description="Alaska", fileFormat="SHP", folder="export"
)
m = geemap.Map(center=[64.838721, -147.763366], zoom=11)
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.
roi = m.user_roi
if roi is None:
roi = ee.Geometry.BBox(-147.9701, 64.7733, -147.5849, 64.8717)
m.add_layer(roi)
m.center_object(roi)
timelapse = geemap.landsat_timelapse(
roi,
out_gif="Fairbanks.gif",
start_year=2000,
end_year=2023,
start_date="06-01",
end_date="09-01",
bands=["SWIR1", "NIR", "Red"],
frames_per_second=5,
title="Landsat Timelapse",
progress_bar_color="blue",
mp4=True,
)
geemap.show_image(timelapse)
m = geemap.Map(center=[64.838721, -147.763366], zoom=11)
m.add_gui("timelapse")
m
m = geemap.Map()
roi = ee.Geometry.BBox(-115.5541, 35.8044, -113.9035, 36.5581)
m.add_layer(roi)
m.center_object(roi)
m
timelapse = geemap.landsat_timelapse(
roi,
out_gif="las_vegas.gif",
start_year=1984,
end_year=2023,
bands=["NIR", "Red", "Green"],
frames_per_second=5,
title="Las Vegas, NV",
font_color="blue",
)
geemap.show_image(timelapse)
m = geemap.Map()
roi = ee.Geometry.BBox(113.8252, 22.1988, 114.0851, 22.3497)
m.add_layer(roi)
m.center_object(roi)
m
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)
Sentinel-2¶
m = geemap.Map(center=[64.838721, -147.763366], zoom=11)
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.
roi = m.user_roi
if roi is None:
roi = ee.Geometry.BBox(-147.9701, 64.7733, -147.5849, 64.8717)
m.add_layer(roi)
m.center_object(roi)
timelapse = geemap.sentinel2_timelapse(
roi,
out_gif="sentinel2.gif",
start_year=2017,
end_year=2023,
start_date="06-01",
end_date="09-01",
frequency="year",
bands=["SWIR1", "NIR", "Red"],
frames_per_second=3,
title="Sentinel-2 Timelapse",
)
geemap.show_image(timelapse)
MODIS¶
MODIS vegetation indices
m = geemap.Map()
m
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)
timelapse = geemap.modis_ndvi_timelapse(
roi,
out_gif="ndvi.gif",
data="Terra",
band="NDVI",
start_date="2000-01-01",
end_date="2022-12-31",
frames_per_second=3,
title="MODIS NDVI Timelapse",
overlay_data="countries",
)
geemap.show_image(timelapse)
MODIS temperature
m = geemap.Map()
m
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)
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¶
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"
timelapse = geemap.goes_timelapse(
roi, "goes.gif", start_date, end_date, data, scan, framesPerSecond=5
)
geemap.show_image(timelapse)
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"
timelapse = geemap.goes_timelapse(
roi, "hurricane.gif", start_date, end_date, data, scan, framesPerSecond=5
)
geemap.show_image(timelapse)
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"
timelapse = geemap.goes_fire_timelapse(
roi, "fire.gif", start_date, end_date, data, scan, framesPerSecond=5
)
geemap.show_image(timelapse)
Time series analysis¶
Visualizing forest cover¶
We will use the Hansen Global Forest Change v1.10 (2000-2022) dataset.
dataset = ee.Image("UMD/hansen/global_forest_change_2022_v1_10")
dataset.bandNames()
Select the imagery for 2000.
m = geemap.Map()
first_bands = ["first_b50", "first_b40", "first_b30"]
first_image = dataset.select(first_bands)
m.add_layer(first_image, {"bands": first_bands, "gamma": 1.5}, "Landsat 2000")
Select the imagery for 2022.
last_bands = ["last_b50", "last_b40", "last_b30"]
last_image = dataset.select(last_bands)
m.add_layer(last_image, {"bands": last_bands, "gamma": 1.5}, "Landsat 2022")
Select the tree cover imagery for 2000.
treecover = dataset.select(["treecover2000"])
treeCoverVisParam = {"min": 0, "max": 100, "palette": ["black", "green"]}
name = "Tree cover (%)"
m.add_layer(treecover, treeCoverVisParam, name)
m.add_colorbar(treeCoverVisParam, label=name, layer_name=name)
m.add("layer_manager")
m
Extract tree cover 2000 by using the threshold of 10%.
threshold = 10
treecover_bin = treecover.gte(threshold).selfMask()
treeVisParam = {"palette": ["green"]}
m.add_layer(treecover_bin, treeVisParam, "Tree cover bin")
Visualizing forest gain and loss¶
Visualize forest loss.
m = geemap.Map(center=[64.864983, -147.840441], zoom=4)
m.add_basemap("Esri.WorldImagery")
treeloss_year = dataset.select(["lossyear"])
treeLossVisParam = {"min": 0, "max": 22, "palette": ["yellow", "red"]}
layer_name = "Tree loss year"
m.add_layer(treeloss_year, treeLossVisParam, layer_name)
m.add_colorbar(treeLossVisParam, label=layer_name, layer_name=layer_name)
m.add("layer_manager")
m
Compare forest loss and gain.
m = geemap.Map(center=[64.864983, -147.840441], zoom=4)
m.add_basemap("Esri.WorldImagery")
treeloss = dataset.select(["loss"]).selfMask()
treegain = dataset.select(["gain"]).selfMask()
m.add_layer(treeloss, {"palette": "red"}, "Tree loss")
m.add_layer(treegain, {"palette": "yellow"}, "Tree gain")
m.add("layer_manager")
m
Calculating forest cover change¶
Compute zonal statistics to find out which county in Alaska has the largest forest area in 2000.
Add a county boundary layer to the map.
counties = ee.FeatureCollection("TIGER/2018/Counties").filter(
ee.Filter.eq("STATEFP", "02")
)
df = geemap.ee_to_df(counties)
df.head()
style = {"color": "0000FFFF", "fillColor": "00000000"}
m.add_layer(counties, {}, "Counties Vector", False)
m.add_layer(counties.style(**style), {}, "Counties Raster")
m
Compute zonal statistics by county.
geemap.zonal_stats(
treecover_bin,
counties,
"forest_cover.csv",
stat_type="SUM",
denominator=1e6,
scale=300,
)
Create a pie chart to visualize the forest area by county.
geemap.pie_chart(
"forest_cover.csv", names="NAME", values="sum", max_rows=20, height=400
)
Create a bar chart to visualize the forest area by county.
geemap.bar_chart(
"forest_cover.csv",
x="NAME",
y="sum",
max_rows=20,
x_label="County",
y_label="Forest area (km2)",
)
Calculate the forest loss area by county.
geemap.zonal_stats(
treeloss.gt(0),
counties,
"treeloss.csv",
stat_type="SUM",
denominator=1e6,
scale=300,
)
Create a bar chart to visualize the forest loss area by county.
geemap.pie_chart("treeloss.csv", names="NAME", values="sum", max_rows=20, height=600)
Create a bar chart to visualize the forest loss area by county.
geemap.bar_chart(
"treeloss.csv",
x="NAME",
y="sum",
max_rows=20,
x_label="County",
y_label="Forest loss area (km2)",
)
Exercise 3 - Analyzing forest cover gain and loss¶
Find out which US state has the largest forest gain and loss between 2000 and 2022. Create pie charts and bar charts to show the results. Relevant Earth Engine assets: