MMR-T1.1.1 Criterion B Details

This notebook provides step-by-step details on Criterion B calculations.

Setup

Import Python modules

import os

import ee
import geopandas as gpd
from lonboard import Map, PolygonLayer, BitmapTileLayer

from rle_python_gee.ee_auth import initialize_ee
from rle_python_gee.ee_rle import load_yaml, make_eoo, area_km2, get_aoo_grid_projection

Initialize Earth Engine

initialize_ee(project=os.environ['GOOGLE_CLOUD_PROJECT'])

Analysis

# Define the ecosystem that this notebook is analyzing
ecosystem_code = 'MMR-T1.1.1'
# Load the country config
country_config_path = os.environ['PIXI_PROJECT_ROOT'] + '/config/country_config.yaml'
country_config = load_yaml(country_config_path)

# Extract the GEE project path from the country config
gee_project_path = country_config['gee_project_path']

# Extract the class info for the ecosystem
class_info = [x for x in country_config['classified_image']['classes'] if x['code'] == ecosystem_code][0]
print(f'{class_info = }')

ecosystem_image = {
    'asset_id': f"{gee_project_path}/{ecosystem_code}/{class_info['id']}",
    'pixel_value': class_info['id']
}
print(f'{ecosystem_image = }')
class_info = {'id': 52, 'name': 'Tanintharyi island rainforest', 'code': 'MMR-T1.1.1'}
ecosystem_image = {'asset_id': 'projects/goog-rle-assessments/assets/MMR-T1.1.1/52', 'pixel_value': 52}
classified_image_asset_id = f"{country_config['gee_project_path']}/{country_config['classified_image']['asset_id']}"
print(f'{classified_image_asset_id = }')

class_img = (
    ee.Image(classified_image_asset_id)
      .eq(ecosystem_image['pixel_value'])
      .selfMask()
)
print(f'class_img: {class_img.getInfo()}')
classified_image_asset_id = 'projects/goog-rle-assessments/assets/mm_ecosys_v7b'
class_img: {'type': 'Image', 'bands': [{'id': 'b1', 'data_type': {'type': 'PixelType', 'precision': 'int', 'min': 0, 'max': 1}, 'dimensions': [14299, 23576], 'crs': 'EPSG:4326', 'crs_transform': [0.0008084837557075694, 0, 89.60991403135986, 0, -0.0008084837557075694, 28.548369897789982]}], 'properties': {'system:footprint': {'type': 'LinearRing', 'coordinates': [[94.66768735253133, 28.548775583027798], [93.22272474215634, 28.548775580809586], [91.41652148658464, 28.548775539452144], [89.60939625454633, 28.54877430628462], [89.60949715757825, 9.485667642440385], [91.41652148658464, 9.486595817478428], [92.86148408075269, 9.486595814990817], [94.30644671824773, 9.486595784089035], [95.39016867380931, 9.48715124036664], [96.47389058147114, 9.486595770655939], [97.91885317497325, 9.486595836813008], [99.36381579356966, 9.486595799973498], [101.1708401374696, 9.48566765739743], [101.17094098108427, 28.548774348455886], [99.00257513499045, 28.548775586792935], [97.55761254104706, 28.54877552452014], [95.75140925342585, 28.548775587887942], [94.66768735253133, 28.548775583027798]]}}}

Extent of occurrence (EOO) (subcriterion B1)

Set the scale (in meters) for reducing the image pixels to polygons. Use the image’s nominal scale unless is is less than 50 meters per pixel.

reduction_scale = max(class_img.projection().nominalScale().getInfo(), 50)
reduction_scale
90

Convert the classified image to vectors.

ecosystem_polygons = class_img.updateMask(1).reduceToVectors(
    scale=reduction_scale,
    geometry=class_img.geometry(),
    geometryType='polygon',
    maxPixels=1e12,
    bestEffort=False
)

Create a convex hull that encompasses the ecosystem polygons.

# convexHull() is called twice as a workaround for a bug
# (https://issuetracker.google.com/issues/465490917)
hull = ecosystem_polygons.geometry().convexHull(maxError=1).convexHull(maxError=1)
# Note:
#   - area() without the projection argument calculates the geodesic area.
#   - area(proj=) calculates the planar area. the area is calculated in projected units (10 km x 10 km grid cells)!
# This is different than the javascript version, which calculates the area without reprojection
aoo_area_km2 = hull.area(1, proj=get_aoo_grid_projection()).getInfo() * 1e2
print(f'The area of the hull is {aoo_area_km2:.2f} km²')
The area of the hull is 50141.45 km²
aoo_area_km2_geodesic = hull.area().getInfo() * 1e-6
print(f'The geodesic area of the hull is {aoo_area_km2_geodesic:.2f} km²')
The geodesic area of the hull is 50336.64 km²

Display EOO Layers

Ecosystem Tiles

Define a tile layer for the displaying the ecosystem image on a Lonboard map. The tiles are served by the Earth Engine API.

tile_url = class_img.getMapId(
    vis_params={
        'palette': ['blue'],
        'opacity': 0.5
    }
)['tile_fetcher'].url_format

tile_layer = BitmapTileLayer(
    data=tile_url,
    tile_size=256,
    max_requests=-1,
    min_zoom=0,
    max_zoom=19,
)

Ecosystem Polygons Layer

ecosystem_polygons_layer = PolygonLayer.from_geopandas(
    ecosystem_polygons_gdf,
    get_fill_color=[255, 0, 0, 127],
    stroked=True,
    get_line_width=2,
    get_line_color=[0, 0, 0, 150],
)

Ecosystem Hull Layer

type(hull_geom_projected)
shapely.geometry.polygon.Polygon
hull_gdf = gpd.GeoDataFrame.from_features(
    ee.FeatureCollection(hull).getInfo(),
    crs='EPSG:4326'
)

hull_layer = PolygonLayer.from_geopandas(
    hull_gdf,
    get_fill_color=[0, 0, 255, 63],
    stroked=True,
    get_line_width=200,
    get_line_color=[0, 0, 0, 150],
)

Display the map.

m = Map(
    layers=[
        ecosystem_polygons_layer,
        tile_layer,
        hull_layer,
    ],
)

m

Area of Occupancy (AOO) (subcriterion B2)

The protocol for this adjustment includes the following steps:

  • Intersect AOO grid with the ecosystem’s distribution map.
  • Calculate extent of the ecosystem type in each grid cell (‘area’) and sum these areas to obtain the total ecosystem area (‘total area’).
  • Arrange grid cells in ascending order based on their area (smaller first).
  • Calculate accumulated sum of area per cell (‘cumulative area’).
  • Calculate ‘cumulative proportion’ by dividing ‘cumulative area’ by ‘total area’ (cumulative proportion takes values between 0 and 1).
  • Calculate AOO by counting the number of cells with a ‘cumulative proportion’ greater than 0.01 (i.e. exclude cells that in combination account for up to 1% of the total mapped extent of the ecosystem type).

Intersect AOO grid with the ecosystem’s distribution map

Load the AOO grid projection

aoo_grid_proj = get_aoo_grid_projection()

aoo_grid_proj.getInfo()
{'type': 'Projection',
 'wkt': 'PROJCS["World_Cylindrical_Equal_Area", \n  GEOGCS["WGS 84", \n    DATUM["WGS_1984", \n      SPHEROID["WGS 84", 6378137.0, 298.257223563, AUTHORITY["EPSG","7030"]], \n      AUTHORITY["EPSG","6326"]], \n    PRIMEM["Greenwich", 0.0], \n    UNIT["degree", 0.017453292519943295], \n    AXIS["Longitude", EAST], \n    AXIS["Latitude", NORTH]], \n  PROJECTION["Cylindrical_Equal_Area"], \n  PARAMETER["central_meridian", 0.0], \n  PARAMETER["standard_parallel_1", 0.0], \n  PARAMETER["false_easting", 0.0], \n  PARAMETER["false_northing", 0.0], \n  UNIT["m", 1.0], \n  AXIS["Easting", EAST], \n  AXIS["Northing", NORTH], \n  AUTHORITY["ESRI","54034"]]',
 'transform': [10000, 0, 0, 0, 10000, 0]}

Extract the grid scale parameters

aoo_x_scale, _, _, _, aoo_y_scale, _ = aoo_grid_proj.getInfo()['transform']
print(f'{aoo_x_scale = } meters')
print(f'{aoo_y_scale = } meters')
aoo_x_scale = 10000 meters
aoo_y_scale = 10000 meters
  • Create an Earth Engine feature collection of AOO grid cells that intersect with the ecosystem, and calculate the fractional coverage of the ecosystem within the grid cell.
fractional_coverage_fc = class_img.unmask().reduceRegions(
  collection=class_img.geometry().coveringGrid(aoo_grid_proj),
  reducer=ee.Reducer.mean(),
).filter(ee.Filter.gt('mean', 0))

# Convert the Earth Engine feature collection to a GeoPandas GeoDataFrame.
fractional_coverage_gdf = ee.data.computeFeatures({
    "expression": fractional_coverage_fc,
    "fileFormat": "GEOPANDAS_GEODATAFRAME",
})
fractional_coverage_gdf.rename(columns={"mean": "coverage"}, inplace=True)

# Set the CRS (Earth Engine data is in EPSG:4326)
fractional_coverage_gdf = fractional_coverage_gdf.set_crs('EPSG:4326')
aoo_grid_cell_count = len(fractional_coverage_gdf)
aoo_grid_cell_count
206

Calculate grid cell area and the total ecosystem area

  • Calculate extent of the ecosystem type in each grid cell (‘area’) and sum these areas to obtain the total ecosystem area (‘total area’).
fractional_coverage_gdf['area'] = fractional_coverage_gdf['coverage'] * aoo_x_scale * aoo_y_scale

fractional_coverage_gdf.sort_values(by="area")[0:4]
geometry coverage area
175 POLYGON ((98.09603 12.57845, 98.18586 12.57845... 0.000035 3481.367444
125 POLYGON ((98.09603 11.93102, 98.18586 11.93102... 0.000079 7875.324046
111 POLYGON ((98.09603 11.74632, 98.18586 11.74632... 0.000079 7880.464175
74 POLYGON ((98.54519 11.1008, 98.63502 11.1008, ... 0.000079 7898.091697
total_area_km2 = fractional_coverage_gdf['area'].sum() / 1e6
print(f'Total ecosystem area: {total_area_km2:.0f} km²')
Total ecosystem area: 1937 km²

Calculate cumulative area in ordered cells

  • Arrange grid cells in ascending order based on their area (smaller first).
  • Calculate accumulated sum of area per cell (‘cumulative area’).
fractional_coverage_gdf = fractional_coverage_gdf.sort_values(by="area")
fractional_coverage_gdf["cumulative_area"] = fractional_coverage_gdf["area"].cumsum()

fractional_coverage_gdf.sort_values(by="area").head()
geometry coverage area cumulative_area
175 POLYGON ((98.09603 12.57845, 98.18586 12.57845... 0.000035 3481.367444 3481.367444
125 POLYGON ((98.09603 11.93102, 98.18586 11.93102... 0.000079 7875.324046 11356.691491
111 POLYGON ((98.09603 11.74632, 98.18586 11.74632... 0.000079 7880.464175 19237.155665
74 POLYGON ((98.54519 11.1008, 98.63502 11.1008, ... 0.000079 7898.091697 27135.247363
160 POLYGON ((97.91637 12.39331, 98.0062 12.39331,... 0.000123 12301.457861 39436.705224

Calculate the cumulative proportion

  • Calculate ‘cumulative proportion’ by dividing ‘cumulative area’ by ‘total area’ (cumulative proportion takes values between 0 and 1).
fractional_coverage_gdf["cumulative_proportion"] = fractional_coverage_gdf["cumulative_area"] / 1e6 / total_area_km2

fractional_coverage_gdf
geometry coverage area cumulative_area cumulative_proportion
175 POLYGON ((98.09603 12.57845, 98.18586 12.57845... 0.000035 3.481367e+03 3.481367e+03 0.000002
125 POLYGON ((98.09603 11.93102, 98.18586 11.93102... 0.000079 7.875324e+03 1.135669e+04 0.000006
111 POLYGON ((98.09603 11.74632, 98.18586 11.74632... 0.000079 7.880464e+03 1.923716e+04 0.000010
74 POLYGON ((98.54519 11.1008, 98.63502 11.1008, ... 0.000079 7.898092e+03 2.713525e+04 0.000014
160 POLYGON ((97.91637 12.39331, 98.0062 12.39331,... 0.000123 1.230146e+04 3.943671e+04 0.000020
... ... ... ... ... ...
108 POLYGON ((98.45536 11.65401, 98.54519 11.65401... 0.549146 5.491457e+07 1.689946e+09 0.872610
102 POLYGON ((98.45536 11.56174, 98.54519 11.56174... 0.561263 5.612625e+07 1.746072e+09 0.901591
164 POLYGON ((98.36552 12.39331, 98.45536 12.39331... 0.586257 5.862569e+07 1.804698e+09 0.931863
99 POLYGON ((98.18586 11.56174, 98.27569 11.56174... 0.607082 6.070816e+07 1.865406e+09 0.963210
93 POLYGON ((98.18586 11.46949, 98.27569 11.46949... 0.712497 7.124972e+07 1.936656e+09 1.000000

206 rows × 5 columns

Calculate AOO

  • Calculate AOO by counting the number of cells with a ‘cumulative proportion’ greater than 0.01 (i.e. exclude cells that in combination account for up to 1% of the total mapped extent of the ecosystem type).
aoo_grid_cells = fractional_coverage_gdf[fractional_coverage_gdf["cumulative_proportion"] > 0.01]
print(f"AOO (number of cells with cumulative proportion > 0.01): {len(aoo_grid_cells)}")
AOO (number of cells with cumulative proportion > 0.01): 139
aoo_grid_cells_dropped = fractional_coverage_gdf[fractional_coverage_gdf["cumulative_proportion"] <= 0.01]

Display the layers

aoo_grid_cells_layer = PolygonLayer.from_geopandas(
    aoo_grid_cells,
    get_fill_color=[0, 255, 0, 63],
)

aoo_grid_cells_layer_dropped = PolygonLayer.from_geopandas(
    aoo_grid_cells_dropped,
    get_fill_color=[255, 0, 0, 63],
)

m = Map(
    layers=[
        tile_layer,
        aoo_grid_cells_layer,
        aoo_grid_cells_layer_dropped
    ],
)

m

Criterion B Summary

from IPython.display import Markdown, display
display(
    Markdown(f'AOO and EOO were measured as '
      f'{aoo_grid_cell_count} 10 x 10 km grid cells '
      f'and {aoo_area_km2:.0f} km², respectively. '
      f'See [Criterion B Details](mmr-t1-1-1-crit-b-test) for more information. '
      f'There is no evidence that suggests this ecosystem should meet the criteria to '
      f'be listed as Near Threatened. The ecosystem is assessed as Least Concern under '
      f'Criterion B1 and B2. '
      f'**Least Concern.**'
    )
)

AOO and EOO were measured as 206 10 x 10 km grid cells and 50141 km², respectively. See Criterion B Details for more information. There is no evidence that suggests this ecosystem should meet the criteria to be listed as Near Threatened. The ecosystem is assessed as Least Concern under Criterion B1 and B2. Least Concern.