ndvi2gif — Land Cover Classification Examples#

This notebook demonstrates supervised and unsupervised land-cover classification workflows using:

  • NdviSeasonality to build multi-temporal feature stacks, and

  • LandCoverClassifier to train/apply classifiers and evaluate accuracy.

Requirements: earthengine-api, geemap, geopandas (if loading local shapefiles), matplotlib, seaborn.

import warnings
warnings.filterwarnings("ignore", category=UserWarning, module="geemap.conversion")
import sys
print(f"Python executable: {sys.executable}")
print(f"Python version: {sys.version}")

# Para ver el entorno conda específico
import os
print(f"CONDA_DEFAULT_ENV: {os.environ.get('CONDA_DEFAULT_ENV', 'No conda env')}")
Python executable: /home/diego/miniconda3/envs/ndvi2gif_clean/bin/python
Python version: 3.9.23 | packaged by conda-forge | (main, Jun  4 2025, 17:57:12) 
[GCC 13.3.0]
CONDA_DEFAULT_ENV: ndvi2gif_clean
import ee
import geemap
import matplotlib.pyplot as plt
from ndvi2gif import NdviSeasonality, TimeSeriesAnalyzer, SpatialTrendAnalyzer, LandCoverClassifier
#from ndvi2gif import NdviSeasonality
#from ndvi2gif import NdviSeasonality
#from ndvi2gif.timeseries import TimeSeriesAnalyzer, SpatialTrendAnalyzer
print("✓ Módulos importados correctamente")

# Authenticate & initialize Earth Engine
# ee.Authenticate() # Just first time 
ee.Initialize()
✓ Módulos importados correctamente
# Create interactive map
Map = geemap.Map()
Map
roi = Map.draw_last_feature
# --- 1) Configure NdviSeasonality
# Example ROI (Lebrija, Spain). Replace with your AOI.
roi = roi #ee.Geometry.Rectangle([-4.80, 37.80, -4.60, 37.95])

proc = NdviSeasonality(
    roi=roi,
    start_year=2020,
    end_year=2023,      # exclusive -> generates 2020–2022
    periods=4,          # seasonal
    sat='Landsat',           # Sentinel-2
    index='ndvi',
    key=['median']
)
There we go again...
Using MODIS Terra + Aqua LST (maximum coverage)
Using all Sentinel-1 orbits (ascending + descending).
Applying S1 ARD preprocessing:
  - Speckle filter: REFINED_LEE
  - Terrain correction: True
  - Terrain model: VOLUME
Landsat collection includes thermal bands: ST_B10 (L8/9) and ST_B6 (L4/5/7)
# --- 2) Build multi-temporal feature stack
clf = LandCoverClassifier(proc)

feature_img = clf.create_feature_stack(
    indices=['ndvi', 'mndwi', 'lst'],   # adjust as needed
    include_statistics=True,   # adds mean/std/max/min per index
    normalize=True             # scales all bands to [0,1]
)
feature_img.bandNames().size().getInfo()
LandCoverClassifier initialized for Landsat
Period: 2020-2022, 4 periods/year
Creating feature stack...
  Processing ndvi...
Year 2020: Successfully processed 4 periods using ndvi index
Year 2021: Successfully processed 4 periods using ndvi index
Year 2022: Successfully processed 4 periods using ndvi index
  Processing mndwi...
Year 2020: Successfully processed 4 periods using mndwi index
Year 2021: Successfully processed 4 periods using mndwi index
Year 2022: Successfully processed 4 periods using mndwi index
  Processing lst...
Year 2020: Successfully processed 4 periods using lst index
Year 2021: Successfully processed 4 periods using lst index
Year 2022: Successfully processed 4 periods using lst index
  Adding temporal statistics...
  Normalizing features...
Feature stack ready: 48 bands
48

Supervised Classification#

Choose one of the following paths to provide training labels:

  • A) Polygons (each polygon has a class field) — points will be sampled per class.

  • B) Points already labeled.

  • C) Stratified sampling from a public land-cover product (quick-start).

# --- 3A) Training from polygons (recommended for balanced classes)
# Provide your polygons file (shapefile or GeoJSON) with a 'class' attribute.
# Example:
# clf.add_training_data(
#     training_polygons="/path/to/polygons.shp",  # or .geojson
#     class_property="class",
#     points_per_class=200
# )

# --- 3B) Training from labeled points
# clf.add_training_data(
#     training_points="/path/to/points.geojson",
#     class_property="class"
# )

# --- 3C) Quick-start: stratified sample from ESA WorldCover (optional)
# worldcover = ee.Image("ESA/WorldCover/v100/2020").select('Map')
# # Map ESA codes to your class IDs (example mapping below; adjust!)
# reclass = worldcover.remap([10,20,30,40,50,60,70,80,90,95],
#                            [1, 2, 3, 4, 5, 6, 7, 8, 9, 10])
# samples = reclass.stratifiedSample(
#     numPoints=1000, classBand='remapped', region=roi, scale=10, geometries=True
# ).rename(['class'])  # ensure property name is 'class'
# clf.add_training_data(training_points=samples, class_property='class')
# --- 4) Train and classify (Random Forest by default)
# Make sure you ran one of 3A/3B/3C above to populate training data.
classified = clf.classify_supervised(
    algorithm='random_forest',
    params={'numberOfTrees': 200, 'bagFraction': 0.6}  # tweak as needed
)
classified
# --- 5) Visualize
m = geemap.Map(center=[37.88, -4.7], zoom=12)
m.addLayer(feature_img, {}, 'Feature stack (NDVI/EVI)')
m.addLayer(classified.randomVisualizer(), {}, 'Classification (RF)')
m.addLayer(roi, {'color': 'yellow'}, 'ROI', False)
m
# --- 6) Accuracy report (computed if validation split exists)
clf.accuracy_results
# --- 7) Plot confusion matrix
# Provide class labels in the same order as your class IDs in training data.
labels = ['Water', 'Urban', 'Crop', 'Forest', 'Bare']  # example
ax = clf.plot_confusion_matrix(labels)
# --- 8) Feature importance (Random Forest only)
try:
    fi = clf.get_feature_importance()
    import pandas as pd
    pd.Series(fi).sort_values(ascending=False).head(20).plot(kind='bar', figsize=(10,3), title='Top Feature Importance');
except ValueError as e:
    print(e)
# --- 9) Export to Google Drive
# task = clf.export_results(
#     description='ndvi2gif_landcover_RF',
#     scale=10,
#     region=roi
# )
# print('Export task started:', task.id)

Unsupervised Classification (Clustering)#

# --- 10) K-means clustering
clusters = clf.classify_unsupervised(
    algorithm='kmeans',
    n_clusters=8,
    max_iterations=30
)
m2 = geemap.Map(center=[37.88, -4.7], zoom=12)
m2.addLayer(clusters.randomVisualizer(), {}, 'K-means (k=8)')
m2.centerObject(roi, zoom=11)
display(m2)
Performing kmeans clustering with 8 clusters...
Clustering complete: 8 clusters
# Export Classification
proc.get_export_single(
    image=clusters.toInt16(),            
    name="landcover_clusters.tif",       
    crs="EPSG:32629",                    
    scale=30                             
)
Generating URL ...
Downloading data from https://earthengine.googleapis.com/v1/projects/16456038103/thumbnails/7ede0554391fa974f51a0b6f159f971e-915f0c5ad372a9501520e291531f897d:getPixels
Please wait ...
Data downloaded to /home/diego/git/Ndvi2Gif/examples_notebooks/landcover_clusters.tif
Image have been exported
# --- 11) Cascade kmeans clustering
clusters = clf.classify_unsupervised(
    algorithm='cascade_kmeans',
    n_clusters=8,
    max_iterations=30
)
m3 = geemap.Map(center=[37.88, -4.7], zoom=12)
m3.addLayer(clusters.randomVisualizer(), {}, 'K-means (k=8)')
m3
Performing cascade_kmeans clustering with 8 clusters...
Clustering complete: 8 clusters
# --- 11) Lda clustering
clusters = clf.classify_unsupervised(
    algorithm='lda',
    n_clusters=8,
    max_iterations=30
)
m4 = geemap.Map(center=[37.88, -4.7], zoom=12)
m4.addLayer(clusters.randomVisualizer(), {}, 'K-means (k=8)')
m4
Performing lda clustering with 8 clusters...
Clustering complete: 8 clusters

Exportar clasificación no supervisada a Drive (Landsat → 30 m)#

clusters = clf.classified_image.rename('class').toInt16()
task = proc.export_to_drive(
    image=clusters,
    description="landcover_kmeans_2020_2022",
    folder="ndvi2gif",
    scale=30,             # 30 m para L5/L7/L8/L9
    crs="EPSG:4326"
)
print("Drive task:", task.id)

Exportar supervisada a un Asset (con política de piramidado para clases)#

classified = clf.classified_image.rename('class').toInt16()
task = proc.export_to_asset(
    image=classified,
    asset_id="users/tu_usuario/ndvi2gif/landcover_rf_2020_2022",
    pyramiding_policy={"class": "mode"},  # importante para clases
    scale=30,
    overwrite=True
)
print("Asset task:", task.id)

Exportar un composite NDVI (S2 → 10 m)#

ndvi_max_2021 = proc.get_year_composite(2021).select('ndvi_max')
proc.export_to_drive(
    image=ndvi_max_2021,
    description="ndvi_max_2021_s2",
    folder="ndvi2gif",
    scale=10,
    crs="EPSG:4326"
)