1_Train_Random_Forest_ML.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": []
   },
   "source": [
    "![Digital Earth Pacific](dep.png)\n",
    "\n",
    "### Digital Earth Pacific Notebook 1 Train Random Forest Machine Learning (ML) Model\n",
    "\n",
    "The objective of this notebook is to train the machine learning model that will allow us to classify all of Fiji within the land cover classes defined through the training data explored in notebook 0. \n",
    "\n",
    "Through the notebook you will be working through the following steps: \n",
    "\n",
    "1. **Setting up your area of interest**  \n",
    "2. **Setting up your time of interest** \n",
    "3. **Training the machine learning model using the training data explored in notebook 0**\n",
    "4. **Exploring the training data in an interactive map**\n",
    "5. **Plotting the training data**\n",
    "6. **Saving the machine learning model in a small file that can be used later**  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "import geopandas as gpd\n",
    "import joblib\n",
    "import numpy as np\n",
    "import odc.geo.xr  # noqa: F401\n",
    "import pandas as pd\n",
    "import xarray as xr\n",
    "from odc.stac import load\n",
    "from pystac_client import Client\n",
    "from shapely import geometry\n",
    "from sklearn.ensemble import RandomForestClassifier\n",
    "from dask.distributed import Client as dask_client\n",
    "from depal_ck import get_image_values\n",
    "import warnings\n",
    "warnings.filterwarnings('ignore')\n",
    "import depal_ck as dep"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "%reload_ext autoreload\n",
    "%autoreload 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": []
   },
   "source": [
    "## Load study area\n",
    "\n",
    "Load data and set up your array to use for prediction"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "# Configure some things up front\n",
    "chunks = dict(x=100, y=100)\n",
    "year = \"2023\"\n",
    "country_code = \"ck\"\n",
    "training_file = \"CK/LULC_Initial.geojson\"\n",
    "aoi = dep.get_country_admin_boundary(\"Cook Islands\", \"Island Council\", \"Rarotonga\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "bbox = dep.get_bbox(aoi)\n",
    "bbox_geometry = geometry.box(*bbox)\n",
    "gdf = gpd.GeoDataFrame({'geometry': [bbox_geometry]}, crs='EPSG:4326')\n",
    "gdf.explore()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "catalog = \"https://stac.staging.digitalearthpacific.org\"\n",
    "client = Client.open(catalog)\n",
    "\n",
    "# Search for Sentinel-2 GeoMAD data\n",
    "items = client.search(\n",
    "    collections=[\"dep_s2_geomad\"],\n",
    "    bbox=bbox,\n",
    "    datetime=year\n",
    ").items()\n",
    "\n",
    "# Load the data\n",
    "data = load(items, chunks=chunks, bbox=bbox, resolution=10).squeeze(\"time\")\n",
    "\n",
    "#coastal clip\n",
    "#data = dep.do_coastal_clip(aoi, data, buffer=0)\n",
    "data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "training_data = gpd.read_file(training_file, bbox=bbox_geometry)\n",
    "training_data.sort_values(by=['lulc'], inplace=True)\n",
    "training_data.to_crs(4326, inplace=True)\n",
    "training_data.explore(\n",
    "    column=\"lulc\",\n",
    "    tiles=\"https://server.arcgisonline.com/ArcGIS/rest/services/World_Imagery/MapServer/tile/{z}/{y}/{x}\",\n",
    "    attr=\"Satellite\",\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# Count the values in a specific column (e.g., \"Class\" column)\n",
    "class_counts = training_data['lulc'].value_counts()\n",
    "\n",
    "# Display the counts\n",
    "class_counts.plot.bar()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# Subset the training data to a smaller sample size, if required\n",
    "#subset = training_data.sample(1000)\n",
    "subset = training_data\n",
    "\n",
    "with dask_client(\n",
    "    n_workers=8, threads_per_worker=8, memory_limit=\"8GB\"\n",
    "):\n",
    "    variables = get_image_values(subset, data)\n",
    "\n",
    "variables.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Sort the variables by name, so they're in a consistent order\n",
    "sorted_variables = variables.sort_index(axis=1)\n",
    "\n",
    "# Join the new variables to the original points and drop non-required columns\n",
    "training_array = pd.concat([training_data[\"class_id\"], sorted_variables], axis=1)\n",
    "training_array = training_array.drop(columns=['time','x','y','spatial_ref'])\n",
    "\n",
    "# Drop rows where there are any NaNs\n",
    "training_array = training_array.dropna()\n",
    "\n",
    "# Explore our data\n",
    "training_array.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "classifier = RandomForestClassifier(\n",
    "    n_estimators=100,\n",
    "    max_depth=10,\n",
    "    min_samples_leaf=10,\n",
    "    n_jobs=-1,\n",
    "    random_state=42,\n",
    ")\n",
    "\n",
    "training_data = np.array(training_array)[:, 1:]\n",
    "classes = np.array(training_array)[:, 0]\n",
    "\n",
    "model = classifier.fit(training_data, classes)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Print feature importances against column headings\n",
    "fields_importances = sorted(\n",
    "    zip(training_array.columns[1:], classifier.feature_importances_),\n",
    "    key=lambda x: x[1],\n",
    "    reverse=True,\n",
    ")\n",
    "\n",
    "for i in fields_importances:\n",
    "    # Format as a table to 2 decinal places\n",
    "    print(f\"{i[0]:<11}| {i[1]:.2f}\")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Export the model for use in the prediction notebook"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "joblib.dump(model, country_code + \"_lulc.model\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "print(\"Model Trained.\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.14"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}