https://colab.research.google.com/assets/colab-badge.svg

Remote Sensing Time Series Anomaly#

Setup#

Ecoscope#

[ ]:
ECOSCOPE_RAW = "https://raw.githubusercontent.com/wildlife-dynamics/ecoscope/master"

# !pip install ecoscope
[ ]:
import os
import sys

import geopandas as gpd
import matplotlib.pyplot as plt
import pandas as pd

import ecoscope

ecoscope.init()

Google Drive Setup#

[ ]:
output_dir = "Ecoscope-Outputs"

if "google.colab" in sys.modules:
    from google.colab import drive

    drive.mount("/content/drive/", force_remount=True)
    output_dir = os.path.join("/content/drive/MyDrive/", output_dir)

os.makedirs(output_dir, exist_ok=True)

Earth Engine#

[ ]:
import ee

try:
    EE_ACCOUNT = os.getenv("EE_ACCOUNT")
    EE_PRIVATE_KEY_DATA = os.getenv("EE_PRIVATE_KEY_DATA")
    if EE_ACCOUNT and EE_PRIVATE_KEY_DATA:
        ee.Initialize(credentials=ee.ServiceAccountCredentials(EE_ACCOUNT, key_data=EE_PRIVATE_KEY_DATA))
    else:
        ee.Initialize()

except ee.EEException:
    ee.Authenticate()
    ee.Initialize()

Load Region Data#

[ ]:
ecoscope.io.download_file(
    f"{ECOSCOPE_RAW}/tests/sample_data/vector/maec_4zones_UTM36S.gpkg",
    os.path.join(output_dir, "maec_4zones_UTM36S.gpkg"),
)

aoi_file = os.path.join(output_dir, "maec_4zones_UTM36S.gpkg")
regions_gdf = gpd.GeoDataFrame.from_file(aoi_file).to_crs(4326)
regions_gdf.set_index("ZONE", drop=True, inplace=True)
[ ]:
regions_gdf.explore()

Calculate Albedo Anomaly#

[ ]:
result = ecoscope.io.eetools.calculate_anomaly(
    gdf=regions_gdf,
    img_coll=ee.ImageCollection("MODIS/006/MCD43C3").select("Albedo_BSA_vis"),
    historical_start="2000-01-01",
    start="2010-01-01",
    end="2022-01-01",
    scale=5000,
)
result["img_date"] = pd.to_datetime(result["img_date"])

Optionally Clean Data#

[ ]:
result = result.dropna()
result = result.reset_index()
[ ]:
result

Plot#

[ ]:
fig, ax = plt.subplots(figsize=(25, 10))
ax.grid(True)
ax.set_ylabel("Albedo Anomaly")
ax.set_xlabel("Time")
ax.set_title("Albedo anomaly 2010-2022 based on mean reference from 2000-2010")
result.groupby("ZONE").apply(lambda df: ax.plot(df["img_date"], df["Albedo_BSA_vis"], label=df.name))
ax.legend()

Export to CSV#

[ ]:
result.to_csv(os.path.join(output_dir, "mara_albedo_anomaly.csv"), header=True)