{ "nbformat": 4, "nbformat_minor": 0, "metadata": { "colab": { "name": "Earth_Engine_TensorFlow_logistic_regression.ipynb", "private_outputs": true, "provenance": [], "collapsed_sections": [], "toc_visible": true }, "kernelspec": { "name": "python3", "display_name": "Python 3" } }, "cells": [ { "cell_type": "code", "metadata": { "id": "fSIfBsgi8dNK" }, "source": [ "#@title Copyright 2021 Google LLC. { display-mode: \"form\" }\n", "# Licensed under the Apache License, Version 2.0 (the \"License\");\n", "# you may not use this file except in compliance with the License.\n", "# You may obtain a copy of the License at\n", "#\n", "# https://www.apache.org/licenses/LICENSE-2.0\n", "#\n", "# Unless required by applicable law or agreed to in writing, software\n", "# distributed under the License is distributed on an \"AS IS\" BASIS,\n", "# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.\n", "# See the License for the specific language governing permissions and\n", "# limitations under the License." ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "aV1xZ1CPi3Nw" }, "source": [ "
\n", "\n", " Run in Google Colab\n", "\n", " View source on GitHub
" ] }, { "cell_type": "markdown", "metadata": { "id": "A5_vH_K4PM6X" }, "source": [ "# Introduction\n", "\n", "## Logistic regression\n", "Logistic regression is a classical machine learning method to estimate the probability of an event occurring (sometimes called the \"risk\"). Specifically, the probability is modeled as a sigmoid function of a linear combination of inputs. This can be implemented as a very simple neural network with a single trainable layer.\n", "\n", "Here, the event being modeled is deforestation in 2016. If a pixel is labeled as deforesetation in 2016 according to the [Hansen Global Forest Change dataset](https://developers.google.com/earth-engine/datasets/catalog/UMD_hansen_global_forest_change_2018_v1_6), the event occurred with probability 1. The probability is zero otherwise. The input variables (i.e. the predictors of this event) are the pixel values of two Landsat 8 surface reflectance median composites, from 2015 and 2017, assumed to represent before and after conditions.\n", "\n", "The model will be hosted on [Google AI Platform](https://cloud.google.com/ai-platform) and used in Earth Engine for interactive prediction from an `ee.Model.fromAIPlatformPredictor`. See [this example notebook](http://colab.research.google.com/github/google/earthengine-api/blob/master/python/examples/ipynb/Earth_Engine_TensorFlow_AI_Platform.ipynb) for background on hosted models.\n", "\n", "**Running this demo may incur charges to your Google Cloud Account!**" ] }, { "cell_type": "markdown", "metadata": { "id": "KiTyR3FNlv-O" }, "source": [ "# Setup software libraries\n", "\n", "Import software libraries and/or authenticate as necessary." ] }, { "cell_type": "markdown", "metadata": { "id": "HsyDopq-yy2b" }, "source": [ "## Authenticate to Colab and Cloud\n", "\n", "*This should be the same account you use to login to Earth Engine*." ] }, { "cell_type": "code", "metadata": { "id": "sYyTIPLsvMWl", "cellView": "code" }, "source": [ "from google.colab import auth\n", "auth.authenticate_user()" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "Ejxa1MQjEGv9" }, "source": [ "## Authenticate to Earth Engine\n", "\n", "*This should be the same account you used to login to Cloud previously*." ] }, { "cell_type": "code", "metadata": { "id": "HzwiVqbcmJIX", "cellView": "code" }, "source": [ "import ee\n", "ee.Authenticate()\n", "ee.Initialize()" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "xhZiXmnSyy2l" }, "source": [ "## Test the TensorFlow installation" ] }, { "cell_type": "code", "metadata": { "cellView": "code", "id": "WjOh_CJeyy2m" }, "source": [ "import tensorflow as tf\n", "print(tf.__version__)" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "w6hPSVmYyy2p" }, "source": [ "## Test the Folium installation" ] }, { "cell_type": "code", "metadata": { "id": "frodQp2syy2q" }, "source": [ "import folium\n", "print(folium.__version__)" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "DrXLkJC2QJdP" }, "source": [ "# Define variables" ] }, { "cell_type": "code", "metadata": { "id": "GHTOc5YLQZ5B" }, "source": [ "# REPLACE WITH YOUR CLOUD PROJECT!\n", "PROJECT = 'your-project'\n", "\n", "# Output bucket for trained models. You must be able to write into this bucket.\n", "OUTPUT_BUCKET = 'your-bucket'\n", "\n", "# Cloud Storage bucket with training and testing datasets.\n", "DATA_BUCKET = 'ee-docs-demos'\n", "\n", "# This is a good region for hosting AI models.\n", "REGION = 'us-central1'\n", "\n", "# Training and testing dataset file names in the Cloud Storage bucket.\n", "TRAIN_FILE_PREFIX = 'logistic_demo_training'\n", "TEST_FILE_PREFIX = 'logistic_demo_testing'\n", "file_extension = '.tfrecord.gz'\n", "TRAIN_FILE_PATH = 'gs://' + DATA_BUCKET + '/' + TRAIN_FILE_PREFIX + file_extension\n", "TEST_FILE_PATH = 'gs://' + DATA_BUCKET + '/' + TEST_FILE_PREFIX + file_extension\n", "\n", "# The labels, consecutive integer indices starting from zero, are stored in\n", "# this property, set on each point.\n", "LABEL = 'loss16'\n", "# Number of label values, i.e. number of classes in the classification.\n", "N_CLASSES = 3\n", "\n", "# Use Landsat 8 surface reflectance data for predictors.\n", "L8SR = ee.ImageCollection('LANDSAT/LC08/C01/T1_SR')\n", "# Use these bands for prediction.\n", "OPTICAL_BANDS = ['B1', 'B2', 'B3', 'B4', 'B5', 'B6', 'B7']\n", "THERMAL_BANDS = ['B10', 'B11']\n", "BEFORE_BANDS = OPTICAL_BANDS + THERMAL_BANDS\n", "AFTER_BANDS = [str(s) + '_1' for s in BEFORE_BANDS]\n", "BANDS = BEFORE_BANDS + AFTER_BANDS\n", "\n", "# Forest loss in 2016 is what we want to predict.\n", "IMAGE = ee.Image('UMD/hansen/global_forest_change_2018_v1_6')\n", "LOSS16 = IMAGE.select('lossyear').eq(16).rename(LABEL)\n", "\n", "# Study area. Mostly Brazil.\n", "GEOMETRY = ee.Geometry.Polygon(\n", " [[[-71.96531166607349, 0.24565390557980268],\n", " [-71.96531166607349, -17.07400853625319],\n", " [-40.32468666607349, -17.07400853625319],\n", " [-40.32468666607349, 0.24565390557980268]]], None, False)\n", "\n", "# These names are used to specify properties in the export of training/testing\n", "# data and to define the mapping between names and data when reading from\n", "# the TFRecord file into a tf.data.Dataset.\n", "FEATURE_NAMES = list(BANDS)\n", "FEATURE_NAMES.append(LABEL)\n", "\n", "# List of fixed-length features, all of which are float32.\n", "columns = [\n", " tf.io.FixedLenFeature(shape=[1], dtype=tf.float32) for k in FEATURE_NAMES\n", "]\n", "\n", "# Dictionary with feature names as keys, fixed-length features as values.\n", "FEATURES_DICT = dict(zip(FEATURE_NAMES, columns))\n", "\n", "# Where to save the trained model.\n", "MODEL_DIR = 'gs://' + OUTPUT_BUCKET + '/logistic_demo_model'\n", "# Where to save the EEified model.\n", "EEIFIED_DIR = 'gs://' + OUTPUT_BUCKET + '/logistic_demo_eeified'\n", "\n", "# Name of the AI Platform model to be hosted.\n", "MODEL_NAME = 'logistic_demo_model'\n", "# Version of the AI Platform model to be hosted.\n", "VERSION_NAME = 'v0'" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "YmBGDcKYO7S2" }, "source": [ "# Generate training data\n", "\n", "This is a multi-step process. First, export the image that contains the prediction bands. When that export completes (several hours in this example), it can be reloaded and sampled to generate training and testing datasets. The second step is to export the traning and testing tables to TFRecord files in Cloud Storage (also several hours)." ] }, { "cell_type": "code", "metadata": { "id": "Lt9OMetJkDe6" }, "source": [ "# Cloud masking function.\n", "def maskL8sr(image):\n", " cloudShadowBitMask = ee.Number(2).pow(3).int()\n", " cloudsBitMask = ee.Number(2).pow(5).int()\n", " qa = image.select('pixel_qa')\n", " mask1 = qa.bitwiseAnd(cloudShadowBitMask).eq(0).And(\n", " qa.bitwiseAnd(cloudsBitMask).eq(0))\n", " mask2 = image.mask().reduce('min')\n", " mask3 = image.select(OPTICAL_BANDS).gt(0).And(\n", " image.select(OPTICAL_BANDS).lt(10000)).reduce('min')\n", " mask = mask1.And(mask2).And(mask3)\n", " return image.select(OPTICAL_BANDS).divide(10000).addBands(\n", " image.select(THERMAL_BANDS).divide(10).clamp(273.15, 373.15)\n", " .subtract(273.15).divide(100)).updateMask(mask)\n", "\n", "# Make \"before\" and \"after\" composites.\n", "composite1 = L8SR.filterDate(\n", " '2015-01-01', '2016-01-01').map(maskL8sr).median()\n", "composite2 = L8SR.filterDate(\n", " '2016-12-31', '2017-12-31').map(maskL8sr).median()\n", "\n", "stack = composite1.addBands(composite2).float()\n", "\n", "export_image = 'projects/google/logistic_demo_image'\n", "\n", "image_task = ee.batch.Export.image.toAsset(\n", " image = stack, \n", " description = 'logistic_demo_image', \n", " assetId = export_image, \n", " region = GEOMETRY,\n", " scale = 30,\n", " maxPixels = 1e10\n", ")" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "3NYcEjex2Yjw" }, "source": [ "First, export the image stack that contains the predictors." ] }, { "cell_type": "code", "metadata": { "id": "FTNhKtTGv5Jn" }, "source": [ "image_task.start()" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "jHH59m6q-OeZ" }, "source": [ "Wait until the image export is completed, then sample the exported image." ] }, { "cell_type": "code", "metadata": { "id": "5eKAILt-cC3i" }, "source": [ "sample = ee.Image(export_image).addBands(LOSS16).stratifiedSample(\n", " numPoints = 10000,\n", " classBand = LABEL,\n", " region = GEOMETRY,\n", " scale = 30,\n", " tileScale = 8\n", ")\n", "\n", "randomized = sample.randomColumn()\n", "training = randomized.filter(ee.Filter.lt('random', 0.7))\n", "testing = randomized.filter(ee.Filter.gte('random', 0.7))\n", "\n", "train_task = ee.batch.Export.table.toCloudStorage(\n", " collection = training,\n", " description = TRAIN_FILE_PREFIX,\n", " bucket = OUTPUT_BUCKET,\n", " fileFormat = 'TFRecord'\n", ")\n", "\n", "test_task = ee.batch.Export.table.toCloudStorage(\n", " collection = testing,\n", " description = TEST_FILE_PREFIX,\n", " bucket = OUTPUT_BUCKET,\n", " fileFormat = 'TFRecord'\n", ")" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "sFQ7vwpU2gER" }, "source": [ "Export the training and testing tables. This also takes a few hours." ] }, { "cell_type": "code", "metadata": { "id": "EjAyMkEFt1W8" }, "source": [ "train_task.start()\n", "test_task.start()" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "QNfaUPbcjuCO" }, "source": [ "# Parse the exported datasets\n", "\n", "Now we need to make a parsing function for the data in the TFRecord files. The data comes in flattened 2D arrays per record and we want to use the first part of the array for input to the model and the last element of the array as the class label. The parsing function reads data from a serialized `Example` proto (i.e. [`example.proto`](https://github.com/tensorflow/tensorflow/blob/r1.12/tensorflow/core/example/example.proto)) into a dictionary in which the keys are the feature names and the values are the tensors storing the value of the features for that example. ([Learn more about parsing `Example` protocol buffer messages](https://www.tensorflow.org/programmers_guide/datasets#parsing_tfexample_protocol_buffer_messages))." ] }, { "cell_type": "code", "metadata": { "id": "x2Q0g3fBj2kD", "cellView": "code" }, "source": [ "def parse_tfrecord(example_proto):\n", " \"\"\"The parsing function.\n", "\n", " Read a serialized example into the structure defined by FEATURES_DICT.\n", "\n", " Args:\n", " example_proto: a serialized Example.\n", "\n", " Returns:\n", " A tuple of the predictors dictionary and the label, cast to an `int32`.\n", " \"\"\"\n", " parsed_features = tf.io.parse_single_example(example_proto, FEATURES_DICT)\n", " labels = parsed_features.pop(LABEL)\n", " return parsed_features, tf.cast(labels, tf.int32)\n", "\n", "\n", "def to_tuple(inputs, label):\n", " \"\"\" Convert inputs to a tuple.\n", "\n", " Note that the inputs must be a tuple of tensors in the right shape.\n", "\n", " Args:\n", " dict: a dictionary of tensors keyed by input name.\n", " label: a tensor storing the response variable.\n", "\n", " Returns:\n", " A tuple of tensors: (predictors, label).\n", " \"\"\"\n", " # Values in the tensor are ordered by the list of predictors.\n", " predictors = [inputs.get(k) for k in BANDS]\n", " return (tf.expand_dims(tf.transpose(predictors), 1),\n", " tf.expand_dims(tf.expand_dims(label, 1), 1)) \n" ], "execution_count": null, "outputs": [] }, { "cell_type": "code", "metadata": { "id": "T3PKyDQW8Vpx", "cellView": "code" }, "source": [ "# Load datasets from the files.\n", "train_dataset = tf.data.TFRecordDataset(TRAIN_FILE_PATH, compression_type='GZIP')\n", "test_dataset = tf.data.TFRecordDataset(TEST_FILE_PATH, compression_type='GZIP')\n", "\n", "# Compute the size of the shuffle buffer. We can get away with this\n", "# because it's a small dataset, but watch out with larger datasets.\n", "train_size = 0\n", "for _ in iter(train_dataset):\n", " train_size+=1\n", "\n", "batch_size = 8\n", "\n", "# Map the functions over the datasets to parse and convert to tuples.\n", "train_dataset = train_dataset.map(parse_tfrecord, num_parallel_calls=4)\n", "train_dataset = train_dataset.map(to_tuple, num_parallel_calls=4)\n", "train_dataset = train_dataset.shuffle(train_size).batch(batch_size)\n", "\n", "test_dataset = test_dataset.map(parse_tfrecord, num_parallel_calls=4)\n", "test_dataset = test_dataset.map(to_tuple, num_parallel_calls=4)\n", "test_dataset = test_dataset.batch(batch_size)\n", "\n", "# Print the first parsed record to check.\n", "from pprint import pprint\n", "pprint(iter(train_dataset).next())" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "Nb8EyNT4Xnhb" }, "source": [ "Note that each record of the parsed dataset contains a tuple. The first element of the tuple is a dictionary with bands for keys and the numeric value of the bands for values. The second element of the tuple is the class label, which in this case is an indicator variable that is one if deforestation happened, zero otherwise." ] }, { "cell_type": "markdown", "metadata": { "id": "t9pWa54oG-xl" }, "source": [ "# Create the Keras model\n", "\n", "This model is intended to represent traditional logistic regression, the parameters of which are estimated through maximum likelihood. Specifically, the probability of an event is represented as the sigmoid of a linear function of the predictors. Training or fitting the model consists of finding the parameters of the linear function that maximize the likelihood function. This is implemented in Keras by defining a model with a single trainable layer, a sigmoid activation on the output, and a crossentropy loss function. Note that the only trainable layer is convolutional, with a 1x1 kernel, so that Earth Engine can apply the model in each pixel. To fit the model, a Stochastic Gradient Descent (SGD) optimizer is used. This differs somewhat from traditional fitting of logistic regression models in that stocahsticity is introduced by using mini-batches to estimate the gradient." ] }, { "cell_type": "code", "metadata": { "id": "OCZq3VNpG--G", "cellView": "code" }, "source": [ "from tensorflow import keras\n", "\n", "# Define the layers in the model.\n", "model = tf.keras.models.Sequential([\n", " tf.keras.layers.Input((1, 1, len(BANDS))),\n", " tf.keras.layers.Conv2D(1, (1,1), activation='sigmoid')\n", "])\n", "\n", "# Compile the model with the specified loss function.\n", "model.compile(optimizer=tf.keras.optimizers.SGD(momentum=0.9),\n", " loss='binary_crossentropy',\n", " metrics=['accuracy'])\n", "\n", "# Fit the model to the training data.\n", "model.fit(x=train_dataset, \n", " epochs=20,\n", " validation_data=test_dataset)" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "shbr6cSXShRg" }, "source": [ "## Save the trained model\n", "\n", "Save the trained model to `tf.saved_model` format in your cloud storage bucket." ] }, { "cell_type": "code", "metadata": { "id": "Sgg7MTXfS1PK" }, "source": [ "model.save(MODEL_DIR, save_format='tf')" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "keijPyVQTIAq" }, "source": [ "# EEification\n", "\n", "The first part of the code is just to get (and SET) input and output names. Keep the input name of 'array', which is how you'll pass data into the model (as an array image)." ] }, { "cell_type": "code", "metadata": { "id": "w49O7n5oTS4w" }, "source": [ "from tensorflow.python.tools import saved_model_utils\n", "\n", "meta_graph_def = saved_model_utils.get_meta_graph_def(MODEL_DIR, 'serve')\n", "inputs = meta_graph_def.signature_def['serving_default'].inputs\n", "outputs = meta_graph_def.signature_def['serving_default'].outputs\n", "\n", "# Just get the first thing(s) from the serving signature def. i.e. this\n", "# model only has a single input and a single output.\n", "input_name = None\n", "for k,v in inputs.items():\n", " input_name = v.name\n", " break\n", "\n", "output_name = None\n", "for k,v in outputs.items():\n", " output_name = v.name\n", " break\n", "\n", "# Make a dictionary that maps Earth Engine outputs and inputs to \n", "# AI Platform inputs and outputs, respectively.\n", "import json\n", "input_dict = \"'\" + json.dumps({input_name: \"array\"}) + \"'\"\n", "output_dict = \"'\" + json.dumps({output_name: \"output\"}) + \"'\"\n", "print(input_dict)\n", "print(output_dict)" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "AX2icXa1UdFF" }, "source": [ "## Run the EEifier\n", "\n", "Use the command line to set your Cloud project and then run the eeifier." ] }, { "cell_type": "code", "metadata": { "id": "IYmH_wCOUhIv" }, "source": [ "!earthengine set_project {PROJECT}\n", "!earthengine model prepare --source_dir {MODEL_DIR} --dest_dir {EEIFIED_DIR} --input {input_dict} --output {output_dict}" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "0_uTqQAaVTIK" }, "source": [ "# Deploy and host the EEified model on AI Platform\n", "\n", "**If you change anything about the model, you'll need to re-EEify it and create a new version!**" ] }, { "cell_type": "code", "metadata": { "id": "8RZRRzcfVu5T" }, "source": [ "!gcloud ai-platform models create {MODEL_NAME} \\\n", " --project {PROJECT} \\\n", " --region {REGION}\n", "\n", "!gcloud ai-platform versions create {VERSION_NAME} \\\n", " --project {PROJECT} \\\n", " --region {REGION} \\\n", " --model {MODEL_NAME} \\\n", " --origin {EEIFIED_DIR} \\\n", " --framework \"TENSORFLOW\" \\\n", " --runtime-version=2.3 \\\n", " --python-version=3.7" ], "execution_count": null, "outputs": [] }, { "cell_type": "markdown", "metadata": { "id": "5aTGza-rWIjp" }, "source": [ "# Connect to the hosted model from Earth Engine\n", "\n", "Now that the model is hosted on AI Platform, point Earth Engine to it and make predictions. These predictions can be thresholded for a rudimentary deforestation detector. Visualize the after imagery, the reference data and the predictions." ] }, { "cell_type": "code", "metadata": { "id": "WtWJMvzAo279" }, "source": [ "# Turn into an array image for input to the model.\n", "array_image = stack.select(BANDS).float().toArray()\n", "\n", "# Point to the model hosted on AI Platform. If you specified a region other\n", "# than the default (us-central1) at model creation, specify it here.\n", "model = ee.Model.fromAiPlatformPredictor(\n", " projectName=PROJECT,\n", " modelName=MODEL_NAME,\n", " version=VERSION_NAME,\n", " # Can be anything, but don't make it too big.\n", " inputTileSize=[8, 8],\n", " # Keep this the same as your training data.\n", " proj=ee.Projection('EPSG:4326').atScale(30),\n", " fixInputProj=True,\n", " # Note the names here need to match what you specified in the\n", " # output dictionary you passed to the EEifier.\n", " outputBands={'output': {\n", " 'type': ee.PixelType.float(),\n", " 'dimensions': 1\n", " }\n", " },\n", ")\n", "\n", "# Output probability.\n", "predictions = model.predictImage(array_image).arrayGet([0])\n", "\n", "# Back-of-the-envelope decision rule.\n", "predicted = predictions.gt(0.7).selfMask()\n", "\n", "# Training data for comparison.\n", "reference = LOSS16.selfMask()\n", "\n", "# Get map IDs for display in folium.\n", "probability_vis = {'min': 0, 'max': 1}\n", "probability_mapid = predictions.getMapId(probability_vis)\n", "\n", "predicted_vis = {'palette': 'red'}\n", "predicted_mapid = predicted.getMapId(predicted_vis)\n", "\n", "reference_vis = {'palette': 'orange'}\n", "reference_mapid = reference.getMapId(reference_vis)\n", "\n", "image_vis = {'bands': ['B4', 'B3', 'B2'], 'min': 0, 'max': 0.3}\n", "image_mapid = composite2.getMapId(image_vis)\n", "\n", "# Visualize the input imagery and the predictions.\n", "map = folium.Map(location=[-9.1, -62.3], zoom_start=11)\n", "folium.TileLayer(\n", " tiles=image_mapid['tile_fetcher'].url_format,\n", " attr='Map Data © Google Earth Engine',\n", " overlay=True,\n", " name='image',\n", ").add_to(map)\n", "folium.TileLayer(\n", " tiles=probability_mapid['tile_fetcher'].url_format,\n", " attr='Map Data © Google Earth Engine',\n", " overlay=True,\n", " name='probability',\n", ").add_to(map)\n", "folium.TileLayer(\n", " tiles=predicted_mapid['tile_fetcher'].url_format,\n", " attr='Map Data © Google Earth Engine',\n", " overlay=True,\n", " name='predicted',\n", ").add_to(map)\n", "folium.TileLayer(\n", " tiles=reference_mapid['tile_fetcher'].url_format,\n", " attr='Map Data © Google Earth Engine',\n", " overlay=True,\n", " name='reference',\n", ").add_to(map)\n", "map.add_child(folium.LayerControl())\n", "map" ], "execution_count": null, "outputs": [] } ] }