Introduction to DEA Mangroves 065bef4b605d47be8be0e681f30dc96a

Background

Mangroves are unique, valuable and vulnerable woody plant species that inhabit intertidal regions around much of the Australian coastline.

They provide a diverse array of ecosystem services such as:

  • coastal protection

  • carbon storage

  • nursery grounds and habitat for a huge variety of avian, coastal and marine animal species.

However, mangrove ecosystems are impacted upon by both natural and anthropogenic drivers of change such as sea-level rise, coastal land reclaimation and severe tropical cyclone damage. In Australia, mangroves are protected by law and consequently, changes in their extent and canopy density are driven predominantly by natural drivers.

Modelling mangrove canopy density offers an effective mechanism to assess how mangroves are responding to these external influences and to monitor their recovery to ensure that they can continue to thrive and support our vital coastal ecosystems.

What this product offers

The DEA Mangroves data product maps the annual canopy cover density of Australian mangroves within a fixed extent around the entire continental coastline. The extent represents a union of Global Mangrove Watch layers for multiple years, produced by the Japanese Aerospace Exploration Agency.

Within this extent, mangroves are identified by leveraging a relationship between the 10th percentile green vegetation component of the DEA Fractional Cover data product with Light Detection And Ranging (LiDAR)-derived Planimetric Canopy Cover% (PCC). More detail on the method can be found here.

Three cover classes are identified within the product which are defined as:

  • Closed Forest - pixels with more than 80 % mangrove canopy cover

  • Open Forest - pixels with between 50 % and 80 % canopy cover

  • Woodland - pixels with between 20 % and 50 % canopy cover

Publications

  • Lymburner, L., Bunting, P., Lucas, R., Scarth, P., Alam, I., Phillips, C., Ticehurst, C., & Held, A., (2020). Mapping the multi-decadal mangrove dynamics of the Australian coastline. Remote Sensing of Environment, 238, 111185. Available at: https://doi.org/10.1016/j.rse.2019.05.004

Description

This notebook introduces the DEA Mangroves data product and steps through how to:

  1. View the product name and associated measurements in the DEA database

  2. Load the dataset and view the data classes within it

  3. Plot a single timestep image

  4. Create and view an animation of the whole timeseries

  5. Plot change over time by graphing the timeseries of each class

  6. Identify hotspot change areas within each class

Note: Visit the DEA Mangroves product documentation for detailed technical information including methods, quality, and data access.


Getting started

To run this analysis, run all the cells in the notebook, starting with the “Load packages” cell.

Load packages

Import Python packages that are used for the analysis.

[1]:
%matplotlib inline

import datacube
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd

from datacube.utils import masking
from IPython.core.display import Video
from matplotlib.colors import LinearSegmentedColormap

import sys

sys.path.insert(1, "../Tools/")
from dea_tools.plotting import display_map, xr_animation

Connect to the datacube

Connect to the datacube so we can access DEA data. The app parameter is a unique name for the analysis which is based on the notebook file name.

[2]:
dc = datacube.Datacube(app="DEA_Mangroves")

Available products and measurements

List products available in Digital Earth Australia

[3]:
dc_products = dc.list_products()
dc_products.loc[["ga_ls_mangrove_cover_cyear_3"]]
[3]:
name description license default_crs default_resolution
name
ga_ls_mangrove_cover_cyear_3 ga_ls_mangrove_cover_cyear_3 Geoscience Australia Landsat Mangrove Cover Ca... CC-BY-4.0 None None

List measurements

View the list of measurements associated with the DEA Mangroves product. Note the single measurement unit, the nodata value and the flags_definitions.

[4]:
dc_measurements = dc.list_measurements()
dc_measurements.loc[["ga_ls_mangrove_cover_cyear_3"]]
[4]:
name dtype units nodata aliases flags_definition
product measurement
ga_ls_mangrove_cover_cyear_3 canopy_cover_class canopy_cover_class uint8 1 255 NaN {'woodland': {'bits': [0, 1, 2, 3, 4, 5, 6, 7]...

Loading data

Select and view your study area

If running the notebook for the first time, keep the default settings below. This will demonstrate how the analysis works and provide meaningful results. Zoom around the displayed map to understand the context of the analysis area. To select a new area, click on the map to reveal the Latitude (y) and Longitude (x) for diagonally opposite corners and place these values into the query

Replace the y and x coordinates to try the following locations:

  • Bowling Green Bay, Qld

    • “y”: (-19.50688750115376, -19.27501266742088)

    • “x”: (147.05183029174404, 147.47617721557216)

  • Pellew Islands, NT

    • “y”: (-15.6786, -16.0075)

    • “x”: (136.5360, 137.0682)

[5]:
# Set up a region to load data
query = {
    "y": (-15.6806, -15.8193),
    "x": (136.5202, 136.6885),
    "time": ("1988", "2022"),
}
display_map(x=query["x"], y=query["y"])
[5]:
Make this Notebook Trusted to load map: File -> Trust Notebook

Load and view DEA Mangroves

[6]:
# Load data from the DEA datacube catalogue
mangroves = dc.load(product="ga_ls_mangrove_cover_cyear_3", **query)
mangroves
[6]:
<xarray.Dataset>
Dimensions:             (time: 35, y: 531, x: 624)
Coordinates:
  * time                (time) datetime64[ns] 1988-07-01T23:59:59.999999 ... ...
  * y                   (y) float64 -1.675e+06 -1.675e+06 ... -1.691e+06
  * x                   (x) float64 4.871e+05 4.872e+05 ... 5.058e+05 5.058e+05
    spatial_ref         int32 3577
Data variables:
    canopy_cover_class  (time, y, x) uint8 255 255 255 255 ... 255 255 255 255
Attributes:
    crs:           EPSG:3577
    grid_mapping:  spatial_ref

View the DEA Mangroves class values and definitions

You’ll see that four classes are identified in this dataset. The notobserved class was separated from nodata pixels in this workflow to identify locations where mangroves have been observed but there is poor observation density. This class is usually insignificant in size and has not been included in the remainder of this notebook analysis.

[7]:
# Extract the flags information from the dataset
flags = masking.describe_variable_flags(mangroves)
flags["bits"] = flags["bits"].astype(str)
flags = flags.sort_values(by="bits")

# Append the class descriptions to each class
descriptors = {
    "closed_forest": "> 80 % canopy cover",
    "open_forest": "50 - 80 % canopy cover",
    "woodland": "20 - 50 % canopy cover",
    "notobserved": "Fewer than 3 clear observations",
}
flags = flags.rename(columns={"description": "class"})
flags["description"] = pd.Series(data=descriptors, index=flags.index)

# View the values in the dataset that are associated with each class
flags
[7]:
bits values class description
woodland [0, 1, 2, 3, 4, 5, 6, 7] {'1': True} Woodland 20 - 50 % canopy cover
notobserved [0, 1, 2, 3, 4, 5, 6, 7] {'0': True} Mangroves not observed Fewer than 3 clear observations
open_forest [0, 1, 2, 3, 4, 5, 6, 7] {'2': True} Open Forest 50 - 80 % canopy cover
closed_forest [0, 1, 2, 3, 4, 5, 6, 7] {'3': True} Closed Forest > 80 % canopy cover

Plotting data

Plot a single timestep

[8]:
# Firstly, mask out the nodata values
mangroves["masked"] = mangroves.canopy_cover_class.where(mangroves.canopy_cover_class != 255)

# Plot the most recent timestep
mangroves["masked"].isel(time=-1).plot(cmap="Greens")
[8]:
<matplotlib.collections.QuadMesh at 0x7f25cde04100>
../../../_images/notebooks_DEA_products_DEA_Mangroves_22_1.png

View all timesteps as an animation

[9]:
# Produce a time series animation of mangrove canopy cover
xr_animation(
    ds=mangroves,
    bands=["masked"],
    output_path="DEA_Mangroves.mp4",
    annotation_kwargs={"fontsize": 20, "color": "white"},
    show_date="%Y",
    imshow_kwargs={"cmap": "Greens"},
    show_colorbar=True,
    interval=1000,
    width_pixels=800,
)

# Plot animation
plt.close()

# View and interact with the animation
Video("DEA_Mangroves.mp4", embed=True)
Exporting animation to DEA_Mangroves.mp4
[9]:

Plot the area (number of pixels) occupied by each class at every timestep in the area of interest

[10]:
# Set up the data for plotting: count the number of pixels per class in the loaded location
mangroves["closed_forest"]= mangroves.canopy_cover_class.where(mangroves.canopy_cover_class == 3)
mangroves["open_forest"] = mangroves.canopy_cover_class.where(mangroves.canopy_cover_class == 2)
mangroves["woodland"] = mangroves.canopy_cover_class.where(mangroves.canopy_cover_class == 1)

y1 = mangroves["woodland"].count(dim=["y", "x"])
y2 = mangroves["open_forest"].count(dim=["y", "x"])
y3 = mangroves["closed_forest"].count(dim=["y", "x"])

# Simplify the date labels for the x-axis
x = np.arange(int(query["time"][0]), int(query["time"][-1]) + 1, 1)

# Prepare the figures
fig, axs = plt.subplots(1, 2, figsize=(12, 5))
fig.suptitle("Regional mangrove canopy change over time")

# Plot the single class summaries
axs[0].plot(x, y1, color="#aedea7", label="Woodland")
axs[0].plot(x, y2, color="#37a055", label="Open Forest")
axs[0].plot(x, y3, color="#00441b", label="Closed Forest")
axs[0].set_title("Canopy class change")
axs[0].set(ylabel="Pixel count", xlabel="Year")

# Stack the classes to plot a snapshot of the region at each time step
axs[1].bar(x, y1, color="#aedea7", label="Woodland")
axs[1].bar(x, y2, color="#37a055", label="Open Forest", bottom=y1)
axs[1].bar(x, y3, color="#00441b", label="Closed Forest", bottom=(y1 + y2))
axs[1].legend(bbox_to_anchor=(1.0, 0.5), loc="center left", frameon=False)
axs[1].set_title("Class distribution")
axs[1].set(xlabel="Year")

plt.tight_layout()
plt.show()
../../../_images/notebooks_DEA_products_DEA_Mangroves_26_0.png

Example applications

Detect change hotspots in each class

Identify change hotspots in the analysis area by assessing the frequency that each class was identified during the length of the timeseries.

Pale yellow coloured areas signify pixels that have a high frequency of being identified as the class of interest over the length of the timeseries and were consistently modelled as that class. Orange to red areas are pixels that were less frequently identified as the class of interest during the timeseries and thus, experienced change between classes during the timeseries.

[11]:
# Create a boolean mask of the `Closed Forest` class
mangroves["closed_forest"] = masking.make_mask(mangroves.canopy_cover_class, closed_forest=True)
# Calculate the frequency of each pixel being identified as `Closed Forest` during the timeseries
mg_cf = mangroves["closed_forest"].mean(dim=["time"])
# Mask out values of 0 to ensure a clean visualisation
mg_cf.where(mg_cf > 0).plot(cmap='YlOrRd_r')
# Title the plot
plt.title("Closed forest frequency (%)")
[11]:
Text(0.5, 1.0, 'Closed forest frequency (%)')
../../../_images/notebooks_DEA_products_DEA_Mangroves_28_1.png
[12]:
# Create a boolean mask of the `Open Forest' class
mangroves["open_forest"] = masking.make_mask(mangroves.canopy_cover_class, open_forest=True)
# Calculate the frequency of each pixel being identified as `Open Forest` during the timeseries
mg_of = mangroves["open_forest"].mean(dim=["time"])
# Mask out values of 0 to ensure a clean visualisation
mg_of.where(mg_of > 0).plot(cmap='YlOrRd_r')
# Title the plot
plt.title("Open forest frequency (%)")
[12]:
Text(0.5, 1.0, 'Open forest frequency (%)')
../../../_images/notebooks_DEA_products_DEA_Mangroves_29_1.png
[13]:
# Create a boolean mask of the `Woodland' class
mangroves["woodland"] = masking.make_mask(mangroves.canopy_cover_class, woodland=True)
# Calculate the frequency of each pixel being identified as `Woodland` during the timeseries
mg_wl = mangroves["woodland"].mean(dim=["time"])
# Mask out values of 0 to ensure a clean visualisation
mg_wl.where(mg_wl > 0).plot(cmap='YlOrRd_r')
# Title the plot
plt.title("Woodland frequency (%)")
[13]:
Text(0.5, 1.0, 'Woodland frequency (%)')
../../../_images/notebooks_DEA_products_DEA_Mangroves_30_1.png

Additional information

License: The code in this notebook is licensed under the Apache License, Version 2.0. Digital Earth Australia data is licensed under the Creative Commons by Attribution 4.0 license.

Contact: If you need assistance, please post a question on the Open Data Cube Discord chat or on the GIS Stack Exchange using the open-data-cube tag (you can view previously asked questions here). If you would like to report an issue with this notebook, you can file one on GitHub.

Last modified: December 2023

Compatible datacube version: 1.8.13

[14]:
print(datacube.__version__)
1.8.13

Tags

Tags: sandbox compatible, landsat, fractional_cover, DEA_Mangroves, time series