Mapping Vegetation Recovery with A Custom Dataset

In this example, we’ll use sankee to visualize vegetation recovery following the 2014 Happy Camp Complex fire in California, USA. Rather than using a premade classified dataset like NLCD or MODIS Land Cover, we’ll use our own custom data derived from Landsat imagery.

Setup

import sankee
import ee

ee.Authenticate()
ee.Initialize()

Building a Custom Dataset

sankee contains built-in classified datasets like NLCD and LCMS, but what if we want to map changes in our own data? To do so, we’ll build a time series of classified Landsat images showing vegetation health.

First, we’ll load the fire perimeter polygon from Monitoring Trends in Burn Severity (MTBS) data. This will be our area of interest.

fire = (ee.FeatureCollection("USFS/GTAC/MTBS/burned_area_boundaries/v1")
    .filterMetadata("Event_ID", "equals", "CA4179612337420140814")
)

Now we can load Landsat 8 imagery from immediately after the fire and five years post-fire to see how vegetation regrew.

immediate = ee.Image("LANDSAT/LC08/C01/T1_TOA/LC08_045031_20150718").clip(fire)
recovery = ee.Image("LANDSAT/LC08/C01/T1_TOA/LC08_046031_20200807").clip(fire)

Finally, we need to turn continuous reflectance data into categorical classes. We’ll do that by calculating NDVI and binning it into three classes–unhealthy, moderate, and healthy.

def calculate_ndvi(img):
    """Calculate NDVI and apply some arbitrary thresholds to classify reflectance into three
    classes with the following pixel values:

    1: Unhealthy (NDVI <= 0.3)
    2: Moderate (NDVI > 0.3 and NDVI <= 0.5)
    3: Healthy (NDVI > 0.5)

    """
    ndvi = img.normalizedDifference(["B5", "B4"])
    thresholds = [-1, 0.3, 0.5]

    return ndvi.gt(thresholds).reduce(ee.Reducer.sum()).rename("health")

immediate_health = calculate_ndvi(immediate)
recovery_health = calculate_ndvi(recovery)

Parameters

Now that we have our classified data, we need to define some parameters that sankey will use to create the diagram. Essentially, we need to tell it how the pixel values we defined above correspond to class labels and colors, as well as what band we used to store our classified band.

labels = {
    1: "Unhealthy",
    2: "Moderate",
    3: "Healthy",
}

palette = {
    1: "#e5f5f9",
    2: "#99d8c9",
    3: "#2ca25f",
}

band = "health"

We also need to build a list of images to use and (optional) corresponding labels to add to the plot.

img_list = [immediate_health, recovery_health]
label_list = ["Immediate", "Recovery"]

Finally, we’re ready to generate our Sankey plot, passing in all the parameters we defined above.

sankee.sankify(
    img_list,
    region=fire,
    band=band,
    labels=labels,
    palette=palette,
    label_list=label_list,
    scale=30
)

As you can see above, much of the unhealthy and moderate areas recovered to higher NDVI values within five years of fire. We could also include more intermediate images, add more classes, use supervised classification, etc.