{ "cells": [ { "cell_type": "markdown", "metadata": { "id": "p72wQoOWuWnU" }, "source": [ "# Agricultural drought\n", "\n", "####Environmental indicator in the STRATA platform https://unepstrata.org/\n", "\n", "The agricultural drought indicator assesses vegetation conditions using NDVI (Normalised Difference Vegetation Index) data. The default threshold is NDVI ≤25th percentile, thereby flagging if data is in the lowest 25% of data for the location historically for the time of year considered.\n" ] }, { "cell_type": "markdown", "metadata": { "id": "uXHJ011_UX2P" }, "source": [ "### Package imports\n", "In order for the code to run, Earth Engine needs to be imported and authenticated.\n", "\n", "For map visualisations, geemap and folium are used." ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "id": "XhLIe0asOyAd", "colab": { "base_uri": "https://localhost:8080/", "height": 211 }, "outputId": "5c972dcc-9093-4156-9a4b-ffe2b69a52a6" }, "outputs": [ { "output_type": "stream", "name": "stdout", "text": [ "geemap package not installed. Installing ...\n", "Installation complete\n", "To authorize access needed by Earth Engine, open the following URL in a web browser and follow the instructions. If the web browser does not start automatically, please manually browse the URL below.\n", "\n", " https://accounts.google.com/o/oauth2/auth?client_id=517222506229-vsmmajv00ul0bs7p89v5m89qs8eb9359.apps.googleusercontent.com&scope=https%3A%2F%2Fwww.googleapis.com%2Fauth%2Fearthengine+https%3A%2F%2Fwww.googleapis.com%2Fauth%2Fdevstorage.full_control&redirect_uri=urn%3Aietf%3Awg%3Aoauth%3A2.0%3Aoob&response_type=code&code_challenge=gLN4Ov2dl9N1qcdwKreArBPpEjvaigdhhM8DJrAjJaY&code_challenge_method=S256\n", "\n", "The authorization workflow will generate a code, which you should paste in the box below. \n", "Enter verification code: 4/1AX4XfWiEKJR2irRE9JGNKhXEOrNmev2GnNLriiK7mZSWOBEcarfDkNegcc0\n", "\n", "Successfully saved authorization token.\n" ] } ], "source": [ "# Installs geemap package\n", "import subprocess\n", "\n", "try:\n", " import geemap\n", "except ImportError:\n", " print('geemap package not installed. Installing ...')\n", " subprocess.check_call([\"python\", '-m', 'pip',\n", " 'install', 'geemap'])\n", "\n", "print('Installation complete')\n", "\n", "# Checks whether this notebook is running on Google Colab\n", "try:\n", " import google.colab\n", " import geemap.eefolium as emap\n", "except:\n", " import geemap as emap\n", "\n", "# Authenticates and initializes Earth Engine\n", "import ee\n", "ee.Authenticate(auth_mode='paste')\n", "Map = emap.Map()\n", "\n", "#Using folium for visualisation\n", "import folium \n", "\n", "import pprint\n", "pp = pprint.PrettyPrinter(indent=4)\n", "\n", "from IPython.display import Image\n" ] }, { "cell_type": "markdown", "metadata": { "id": "D1nFckQOVgtC" }, "source": [ "### Data visualisation\n", "The following code allows Earth Engine data to be displayed on a map by specifying the following: \n", "* The layer to visualise (ee_object)\n", "* Visualisation parameters (vis_params) \n", "* A name for the layer (name)\n" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "id": "wYL0ZHB3PmIT" }, "outputs": [], "source": [ "def add_ee_layer(self, ee_object, vis_params, name):\n", " try: \n", " # display ee.Image()\n", " if isinstance(ee_object, ee.image.Image): \n", " map_id_dict = ee.Image(ee_object).getMapId(vis_params)\n", " folium.raster_layers.TileLayer(\n", " tiles = map_id_dict['tile_fetcher'].url_format,\n", " attr = 'Map Data © \\\n", " Google Earth Engine',\n", " name = name,\n", " overlay = True,\n", " control = True\n", " ).add_to(self)\n", " # display ee.ImageCollection()\n", " elif isinstance(ee_object,\n", " ee.imagecollection.ImageCollection): \n", " ee_object_new = ee_object.mosaic()\n", " map_id_dict = ee.Image(ee_object_new) \\\n", " .getMapId(vis_params)\n", " folium.raster_layers.TileLayer(\n", " tiles = map_id_dict['tile_fetcher'].url_format,\n", " attr = 'Map Data © \\\n", " Google Earth Engine',\n", " name = name,\n", " overlay = True,\n", " control = True\n", " ).add_to(self)\n", " except Exception as e:\n", " print(\"Could not display {}\".format(name), e )\n", " \n", "# Add EE drawing method to folium.\n", "folium.Map.add_ee_layer = add_ee_layer" ] }, { "cell_type": "markdown", "metadata": { "id": "oow_5O3eXVYp" }, "source": [ "### Bitmap generation function \n", "Function to convert an image into a bitmap using supplied threshold parameters" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "id": "cVmS-TebxyTY" }, "outputs": [], "source": [ "def generateBitmapWithComparison(image, comparison,\n", " threshold, weight):\n", " if (comparison == 'lt'):\n", " bit_img = image.lt(threshold)\n", " else:\n", " bit_img = image\n", " \n", " weighted_image = image.expression('WGT * BTW', {\n", " 'WGT': weight,\n", " 'BTW': bit_img\n", " })\n", " return weighted_image" ] }, { "cell_type": "markdown", "metadata": { "id": "PanBw2aHWii1" }, "source": [ "### User input variables\n", "\n", "\n", "* Selected to and from dates \n", "* Selected area of interest\n", "\n" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "id": "ozeIAXaiWd5R" }, "outputs": [], "source": [ "#Current analysis period\n", "start_date = ee.Date('2020-10-15')\n", "end_date = ee.Date('2021-02-02')\n", "\n", "# Area of interest - Somalia\n", "country_code = 226" ] }, { "cell_type": "markdown", "metadata": { "id": "6zyBEvSpsRKq" }, "source": [ "###Datasets\n", "\n", "\n", "* [GAUL country boundaries](https://developers.google.com/earth-engine/datasets/catalog/FAO_GAUL_2015_level0) \n", "* [MODIS Terra NDVI](https://developers.google.com/earth-engine/datasets/catalog/MODIS_006_MOD13A1)\n", "\n" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "id": "gqMWQRXBWxXa" }, "outputs": [], "source": [ "#Global Administrative Unit Layers filtered to Somalia \n", "somalia_boundary = ee.FeatureCollection( \\\n", " \"FAO/GAUL/2015/level0\") \\\n", " .filter(ee.Filter.eq('ADM0_CODE', country_code))\n", "\n", "#MODIS Terra vegetation indices (NDVI selected) \n", "NDVI = ee.ImageCollection('MODIS/006/MOD13A1').select('NDVI') " ] }, { "cell_type": "markdown", "metadata": { "id": "XAmN3t1SY-yd" }, "source": [ "# Start of indicator code " ] }, { "cell_type": "markdown", "metadata": { "id": "qf7rCr8pstVr" }, "source": [ "### NDVI for user selected analysis period" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "id": "gOxMbwgxoBHy" }, "outputs": [], "source": [ "#Determine days of year for analysis period\n", "start_day = start_date.format('DD')\n", "end_day = end_date.format('DD')\n", "\n", "start_day_of_month = start_date.format('dd')\n", "end_day_of_month = end_date.format('dd')\n", "\n", "start_month = start_date.format('MM')\n", "end_month = end_date.format('MM')\n", "\n", "diff_days = end_date.difference(start_date, 'days')\n", "\n", "#Find median NDVI composite across specified\n", "# analysis period and clip to boundary\n", "analysis_period_NDVI = NDVI.filter \\\n", " (ee.Filter.date(start_date, end_date)).median() \n", "analysis_period_NDVI = analysis_period_NDVI \\\n", " .clipToCollection(somalia_boundary)\n" ] }, { "cell_type": "markdown", "metadata": { "id": "u3PAht8_pzlQ" }, "source": [ "NDVI for user selected period visualisation " ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "colab": { "base_uri": "https://localhost:8080/", "height": 421 }, "id": "NrNWp91vqFN0", "outputId": "4883a8c4-e618-4b98-8f1b-f1f3cc1bcb69" }, "outputs": [ { "output_type": "execute_result", "data": { "text/plain": [ "" ], "text/html": [ "" ] }, "metadata": {}, "execution_count": 7 } ], "source": [ "palette = [\n", " '#FFFFFF', '#CE7E45', '#DF923D', '#F1B555', \n", " '#FCD163', '#99B718', '#74A901',\n", " '#66A000', '#529400', '#3E8601', '#207401', \n", " '#056201', '#004C00', '#023B01',\n", " '#012E01', '#011D01', '#011301'\n", " ]\n", "Image(url = analysis_period_NDVI.getThumbURL(\n", " {'min': 0, 'max': 9000, \n", " 'region': somalia_boundary.geometry(), 'dimensions': 400,\n", " 'palette':palette }))" ] }, { "cell_type": "markdown", "metadata": { "id": "_JEXs88WuEA1" }, "source": [ "###Historic NDVI\n", "This code aims to find NDVI images in the selected time period per year between 2003 and the selected year in order to calculate a historic baseline NDVI. " ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "colab": { "base_uri": "https://localhost:8080/" }, "id": "9UNU7o3UnMgt", "outputId": "598356ea-50e9-4bf9-8ea6-88ab7cd38701" }, "outputs": [ { "output_type": "stream", "name": "stdout", "text": [ "start_year: 2020\n", "Years of historic data:\n", "[ 2003,\n", " 2004,\n", " 2005,\n", " 2006,\n", " 2007,\n", " 2008,\n", " 2009,\n", " 2010,\n", " 2011,\n", " 2012,\n", " 2013,\n", " 2014,\n", " 2015,\n", " 2016,\n", " 2017,\n", " 2018,\n", " 2019]\n" ] } ], "source": [ "#Select analysis period over which to look for\n", "# historic NDVI images\n", "#The same time period each year is selected \n", "historic_NDVI_images = NDVI.filter(\n", " ee.Filter.calendarRange(\n", " ee.Number.parse(start_day),\n", " ee.Number.parse(end_day),\"day_of_year\"))\n", "\n", "#Determine start year \n", "start_year = start_date.format(\"Y\")\n", "print(\"start_year:\", start_year.getInfo())\n", "\n", "#List years from 2003 to year before analysis\n", "# period start date\n", "years = ee.List.sequence(\n", " 2003, ee.Number.parse(start_year)\n", " .subtract(1))\n", "print(\"Years of historic data:\")\n", "pp.pprint(years.getInfo())\n", "\n", "#For each year, find median NDVI \n", "def get_yearly_NDVI(year):\n", " start = ee.Date.fromYMD(\n", " year, ee.Number.parse(start_month), \n", " ee.Number.parse(start_day_of_month))\n", " #If start date and end date straddle calendar year \n", " #we can capture same period for each year\n", " end = start.advance(diff_days , 'day')\n", " annual = historic_NDVI_images.filter(\n", " ee.Filter.date(start, end)).median()\n", " return annual.clipToCollection(somalia_boundary)\n", " \n", "#Create image collection for each year \n", "historic_NDVI_images_year_aggregation = \\\n", " ee.ImageCollection.fromImages(years.map(get_yearly_NDVI))\n", "\n", "#Create historic threshold image using percentile reducer \n", "historic_threshold = \\\n", " historic_NDVI_images_year_aggregation.reduce(\n", " ee.Reducer.percentile([25]))" ] }, { "cell_type": "markdown", "metadata": { "id": "f3TVXFjMkqnj" }, "source": [ "Historic NDVI aggregation visualisation " ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "colab": { "base_uri": "https://localhost:8080/", "height": 421 }, "id": "xBvSfAyPklUn", "outputId": "4810e459-bc80-4307-d3a7-91af945cb7f8" }, "outputs": [ { "output_type": "execute_result", "data": { "text/plain": [ "" ], "text/html": [ "" ] }, "metadata": {}, "execution_count": 9 } ], "source": [ "Image(url = historic_NDVI_images_year_aggregation.median()\n", " .getThumbURL(\n", " {'min': 0, \n", " 'max': 9000, \n", " 'region': somalia_boundary.geometry(), \n", " 'dimensions': 400,\n", " 'palette':palette }))" ] }, { "cell_type": "markdown", "metadata": { "id": "quFQXdZVkdHP" }, "source": [ "Historic threshold visualisation" ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "colab": { "base_uri": "https://localhost:8080/", "height": 421 }, "id": "7IIh9YVrqYK8", "outputId": "bff691bc-c743-41cc-8d44-5d598af76b67" }, "outputs": [ { "output_type": "execute_result", "data": { "text/plain": [ "" ], "text/html": [ "" ] }, "metadata": {}, "execution_count": 10 } ], "source": [ "Image(url = historic_threshold.getThumbURL(\n", " {'min': 0, \n", " 'max': 9000,\n", " 'region': somalia_boundary.geometry(), \n", " 'dimensions': 400,\n", " 'palette':palette }))" ] }, { "cell_type": "markdown", "metadata": { "id": "xfSt9xKcvQtL" }, "source": [ "###Produce bitmap\n", "Change is defined as lower than 75% historical average\n", "\n", "Generate bitmap using bitmap function " ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "id": "Cyvx2VeaoXZb" }, "outputs": [], "source": [ "bitmap = generateBitmapWithComparison(\n", " analysis_period_NDVI,'lt',historic_threshold, 1) \\\n", " .clip(somalia_boundary)" ] }, { "cell_type": "markdown", "metadata": { "id": "LOf5KjnWw9ub" }, "source": [ "###Visualisation parameters\n", "\n" ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "id": "tjMukBCSrqVk" }, "outputs": [], "source": [ "bitmap_params = {\n", " 'min': 0,\n", " 'max': 1,\n", " 'palette': ['#000000', '#FFFFFF']\n", "}" ] }, { "cell_type": "markdown", "metadata": { "id": "jHXSba0ZUf5u" }, "source": [ "###Add layers to the map" ] }, { "cell_type": "code", "execution_count": 13, "metadata": { "colab": { "base_uri": "https://localhost:8080/", "height": 918 }, "id": "s-kw2n5-sgHZ", "outputId": "afb1a7cf-2f48-42a4-db98-36c3ae7a964b" }, "outputs": [ { "output_type": "execute_result", "data": { "text/plain": [ "" ], "text/html": [ "
Make this Notebook Trusted to load map: File -> Trust Notebook
" ] }, "metadata": {}, "execution_count": 13 } ], "source": [ "#Add Google satellite basemap\n", "basemap = { \n", " 'Google Satellite': folium.TileLayer(\n", " tiles = \\\n", " 'https://mt1.google.com/vt/lyrs=s&x={x}&y={y}&z={z}',\n", " attr = 'Google',\n", " name = 'Google Satellite',\n", " overlay = True,\n", " control = True\n", " )\n", "}\n", "\n", "Map = folium.Map(\n", " location=[5.1, 46.1], \n", " width=600,\n", " height=600,\n", " zoom_start=6)\n", "basemap['Google Satellite'].add_to(Map)\n", "\n", "#NDVI for analysis period \n", "Map.add_ee_layer(analysis_period_NDVI,\n", " {'min':0, \n", " 'max':9000,\n", " 'palette':palette}, \n", " \"NDVI for analysis period\")\n", "\n", "#Historic NDVI aggregation \n", "Map.add_ee_layer(historic_NDVI_images_year_aggregation\n", " .median(),\n", " {'min':0,\n", " 'max':9000, \n", " 'palette':palette},\n", " \"Historic NDVI aggregation\")\n", "\n", "#Historic threshold \n", "Map.add_ee_layer(historic_threshold,\n", " {'min':0,\n", " 'max':9000,\n", " 'palette':palette},\n", " \"Historic threshold\")\n", "\n", "#Add final bitmap to map \n", "Map.add_ee_layer(bitmap,bitmap_params,\n", " \"Agricultural drought bitmap\")\n", "\n", "Map.add_child(folium.LayerControl())\n", "\n", "Map" ] } ], "metadata": { "colab": { "collapsed_sections": [], "name": "STRATA Agricultural Drought Indicator Code Example.ipynb", "provenance": [] }, "kernelspec": { "display_name": "Python 3", "name": "python3" }, "language_info": { "name": "python" } }, "nbformat": 4, "nbformat_minor": 0 }