wies.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "\n",
       "<style>\n",
       ".output_wrapper button.btn.btn-default,\n",
       ".output_wrapper .ui-dialog-titlebar {\n",
       "  display: none;\n",
       "}\n",
       "</style>\n"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# %matplotlib notebook\n",
    "\n",
    "from IPython.display import HTML\n",
    "\n",
    "HTML('''\n",
    "<style>\n",
    ".output_wrapper button.btn.btn-default,\n",
    ".output_wrapper .ui-dialog-titlebar {\n",
    "  display: none;\n",
    "}\n",
    "</style>\n",
    "''')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Introduction\n",
    "\n",
    "In this activity, we'll be using satellite imagery to automatically determine how the water level of a particular lake changes over time. We can then use this to predict what will happen to the lake in the future! \n",
    "\n",
    "In order to do this, we need to develop a system that will detect water in an image."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Colour Selection\n",
    "\n",
    "Our satellite images are colour images. All colour images consist of three components: red, green and blue. The first step is to select which colours (or ratio of colours) we want to work with.\n",
    "\n",
    "Click on the buttons below to select which colours you want. The graph on the left is the original image, while the one on the right is the colour component you selected (note that the individual colour components render in greyscale!) \n",
    "\n",
    "Think about the three colours and how they relate to land and water. We want an image that will separate the water from the land as best as possible. But how do we know which one is best?\n",
    "\n",
    "The histogram graph below shows the magnitude frequency of the new image. **The best choices will result in a graph with two peaks (a bimodal graph), which correspond to land and water respectively!**   \n",
    "\n",
    "Once you've made your choice, move on to the next step."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from __future__ import print_function\n",
    "from ipywidgets import interact, interactive, fixed, interact_manual\n",
    "import ipywidgets as widgets\n",
    "from skimage.color import label2rgb\n",
    "\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy import ndimage as ndi\n",
    "from skimage import data\n",
    "from skimage import io\n",
    "from skimage import color\n",
    "from IPython.display import HTML\n",
    "\n",
    "global image_rgb\n",
    "global image\n",
    "global image_cleaned\n",
    "global threshold_value\n",
    "global choice\n",
    "\n",
    "choice = \"All Bands\"\n",
    "threshold_value = 0\n",
    "image_rgb = io.imread('1984.png')\n",
    "image = image_rgb\n",
    "training_years = range(1984, 2008)\n",
    "testing_years = range(2008, 2030) \n",
    "\n",
    "import warnings\n",
    "warnings.filterwarnings(\"ignore\")\n",
    " \n",
    "def _select_bands(image_rgb, choice):\n",
    "    if choice == \"All Bands\":\n",
    "        image = color.rgb2grey(image_rgb) * 255\n",
    "    elif choice == \"Red Band\":\n",
    "        image = image_rgb[:, :, 0] * 255\n",
    "    elif choice == \"Green Band\":\n",
    "        image = image_rgb[:, :, 1] * 255\n",
    "    elif choice == \"Blue Band\":\n",
    "        image = image_rgb[:, :, 2] * 255\n",
    "    elif choice == \"Red-Green Ratio\":\n",
    "        image = image_rgb[:, :, 0] / image_rgb[:, :, 1]\n",
    "        image = image / image.max() * 255   \n",
    "    elif choice == \"Red-Blue Ratio\":\n",
    "        image = image_rgb[:, :, 0] / image_rgb[:, :, 2]\n",
    "        image = image / image.max() * 255   \n",
    "    elif choice == \"Green-Red Ratio\":\n",
    "        image = image_rgb[:, :, 1] / image_rgb[:, :, 0]\n",
    "        image = image / image.max() * 255   \n",
    "    elif choice == \"Green-Blue Ratio\":\n",
    "        image = image_rgb[:, :, 1] / image_rgb[:, :, 2]\n",
    "        image = image / image.max() * 255   \n",
    "    elif choice == \"Blue-Red Ratio\":\n",
    "        image = image_rgb[:, :, 2] / image_rgb[:, :, 0]\n",
    "        image = image / image.max() * 255   \n",
    "    elif choice == \"Blue-Green Ratio\":\n",
    "        image = image_rgb[:, :, 2] / image_rgb[:, :, 1]\n",
    "        image = image / image.max() * 255   \n",
    "    else:\n",
    "        raise ValueError\n",
    "    return image\n",
    "    \n",
    "    \n",
    "\n",
    "\n",
    "\n",
    "def segment_image(image_rgb, spectrum, threshold):\n",
    "    \"\"\"\n",
    "    Return the segmented image and area\n",
    "    \"\"\"\n",
    "    im = _select_bands(image_rgb, spectrum)\n",
    "    clean_im = morphology.remove_small_objects((im < threshold), 50)\n",
    "    segmentation = ndi.binary_fill_holes(clean_im)\n",
    "    labeled_image, _ = ndi.label(segmentation)\n",
    "    image_label_overlay = label2rgb(labeled_image, image=clean_im)\n",
    "    plt.imshow(im_rgb, cmap=plt.cm.gray, interpolation='nearest')\n",
    "    plt.contour(segmentation, [0.5], linewidths=1.2, colors='r')\n",
    "    return plt, image_label_overlay[:,:,1].sum()\n",
    "\n",
    "\n",
    "def _calculate_area(clean_im):\n",
    "    segmentation = ndi.binary_fill_holes(clean_im)\n",
    "    labeled_image, _ = ndi.label(segmentation)\n",
    "    image_label_overlay = label2rgb(labeled_image, image=clean_im)\n",
    "    return image_label_overlay[:,:,1].sum()\n",
    "\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-02-19T14:54:52.965362Z",
     "start_time": "2018-02-19T14:54:52.912496Z"
    },
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "9e9dd9f64c0741a299aa74134a3a79f1",
       "version_major": 2,
       "version_minor": 0
      },
      "text/html": [
       "<p>Failed to display Jupyter Widget of type <code>interactive</code>.</p>\n",
       "<p>\n",
       "  If you're reading this message in the Jupyter Notebook or JupyterLab Notebook, it may mean\n",
       "  that the widgets JavaScript is still loading. If this message persists, it\n",
       "  likely means that the widgets JavaScript library is either not installed or\n",
       "  not enabled. See the <a href=\"https://ipywidgets.readthedocs.io/en/stable/user_install.html\">Jupyter\n",
       "  Widgets Documentation</a> for setup instructions.\n",
       "</p>\n",
       "<p>\n",
       "  If you're reading this message in another frontend (for example, a static\n",
       "  rendering on GitHub or <a href=\"https://nbviewer.jupyter.org/\">NBViewer</a>),\n",
       "  it may mean that your frontend doesn't currently support widgets.\n",
       "</p>\n"
      ],
      "text/plain": [
       "interactive(children=(ToggleButtons(description='Colours:', options=('All Bands', 'Red Band', 'Green Band', 'Blue Band', 'Red-Green Ratio', 'Red-Blue Ratio', 'Green-Red Ratio', 'Green-Blue Ratio', 'Blue-Red Ratio', 'Blue-Green Ratio'), value='All Bands'), Output()), _dom_classes=('widget-interact',))"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "def select_bands(c):\n",
    "    \n",
    "    global image\n",
    "    global image_rgb\n",
    "    global choice\n",
    "    choice = c\n",
    "    image = _select_bands(image_rgb, c)\n",
    "    \n",
    "\n",
    "    fig = plt.figure(figsize=(12, 10))\n",
    "    ax1 = fig.add_subplot(2,2,1)\n",
    "    ax1.set_title(\"Original Image\")\n",
    "    ax1.imshow(image_rgb, cmap=plt.cm.gray, interpolation='nearest')\n",
    "    ax2 = fig.add_subplot(2,2,2)\n",
    "    ax2.set_title(\"New Image\")\n",
    "    ax2.imshow(image, cmap=plt.cm.gray, interpolation='nearest')\n",
    "    ax3 = fig.add_subplot(2,2,3)\n",
    "    hist = np.histogram(image, bins=np.arange(0, 256))\n",
    "    ax3.set_title(\"Histogram\")\n",
    "   \n",
    "    plt.plot(hist[1][:-1], hist[0], lw=2)\n",
    "    \n",
    "    \n",
    "    #plt.imshow(image_rgb, cmap=plt.cm.gray, interpolation='nearest')\n",
    "    #plt.show()\n",
    "    #plt.imshow()\n",
    "    #plt.show()\n",
    "    ##hist = np.histogram(image, bins=np.arange(0, 256))\n",
    "    #plt.plot(hist[1][:-1], hist[0], lw=2)\n",
    "    #plt.show()\n",
    "\n",
    "\n",
    "buttons = widgets.ToggleButtons(\n",
    "    options=['All Bands', 'Red Band', 'Green Band', 'Blue Band', 'Red-Green Ratio', 'Red-Blue Ratio', 'Green-Red Ratio', 'Green-Blue Ratio','Blue-Red Ratio', 'Blue-Green Ratio'],\n",
    "    description='Colours:',\n",
    "    disabled=False,\n",
    "    selected=2,\n",
    "    button_style='',\n",
    ")\n",
    "\n",
    "interact(select_bands, c=buttons);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Note:\n",
    "\n",
    "Why didn't selecting blue work at all? Lakes are blue, right?! Well, that's because the land, which is a bit greyish, also contains blue. Additionally, plantlife is green, which is quite similar to blue. Red, on the other hand, has absolutely nothing at all to do with water! And thus it's a good choice of colour, since the **absence** of red is a good indication of water!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Thresholding\n",
    "\n",
    "We now pick which sections of the image are part of the lake, and which are part of the land. Drag the slider to select a value. All pixels with values greater than the threshold will be considered as the lake. Everything else is the land.\n",
    "\n",
    "**Use the histogram above to select the threshold by picking a value in between the two peaks!**\n",
    "\n",
    "Try pick the value so the image looks like the shape of the lake (remember that the lake has two islands!). Once you're satisfied, move to the next step. If the thresholding doesn't result in a lake-shaped image, you may have to pick different colours."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-02-19T14:56:22.213479Z",
     "start_time": "2018-02-19T14:56:21.823732Z"
    },
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "50ae1f5b0e7c46ecae7dc1a84783f950",
       "version_major": 2,
       "version_minor": 0
      },
      "text/html": [
       "<p>Failed to display Jupyter Widget of type <code>interactive</code>.</p>\n",
       "<p>\n",
       "  If you're reading this message in the Jupyter Notebook or JupyterLab Notebook, it may mean\n",
       "  that the widgets JavaScript is still loading. If this message persists, it\n",
       "  likely means that the widgets JavaScript library is either not installed or\n",
       "  not enabled. See the <a href=\"https://ipywidgets.readthedocs.io/en/stable/user_install.html\">Jupyter\n",
       "  Widgets Documentation</a> for setup instructions.\n",
       "</p>\n",
       "<p>\n",
       "  If you're reading this message in another frontend (for example, a static\n",
       "  rendering on GitHub or <a href=\"https://nbviewer.jupyter.org/\">NBViewer</a>),\n",
       "  it may mean that your frontend doesn't currently support widgets.\n",
       "</p>\n"
      ],
      "text/plain": [
       "interactive(children=(IntSlider(value=0, continuous_update=False, description='Threshold:', max=255), Output()), _dom_classes=('widget-interact',))"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "<function __main__.threshold>"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "def threshold(val):\n",
    "    global threshold_value\n",
    "    threshold_value = val\n",
    "    \n",
    "    fig = plt.figure(figsize=(9, 9))\n",
    "    ax = fig.add_subplot(1,1,1)\n",
    "    ax.set_title(\"Image after Thresholding\")\n",
    "    ax.imshow(image < val, cmap=plt.cm.gray, interpolation='nearest')\n",
    "    \n",
    "    #plt.imshow(image > val, cmap=plt.cm.gray, interpolation='nearest')\n",
    "    \n",
    "    \n",
    "    \n",
    "    #plt.show()\n",
    "    #plt.tight_layout()\n",
    "\n",
    "threshold_slider = widgets.IntSlider(\n",
    "    value=0,\n",
    "    min=0,\n",
    "    max=255,\n",
    "    step=1,\n",
    "    description='Threshold:',\n",
    "    disabled=False,\n",
    "    continuous_update=False,\n",
    "    orientation='horizontal',\n",
    "    readout=True,\n",
    "    readout_format='d'\n",
    ")\n",
    "\n",
    "interact(threshold, val=threshold_slider)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Remove small objects\n",
    "\n",
    "There may be small specks in the above image. To remove these artefacts, click on the button below which removes shapes smaller than a certain size. This will clean up the image a bit."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-02-19T14:57:01.024278Z",
     "start_time": "2018-02-19T14:57:00.863971Z"
    },
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "c70bfb9352884f0e933281fbe54dc403",
       "version_major": 2,
       "version_minor": 0
      },
      "text/html": [
       "<p>Failed to display Jupyter Widget of type <code>Button</code>.</p>\n",
       "<p>\n",
       "  If you're reading this message in the Jupyter Notebook or JupyterLab Notebook, it may mean\n",
       "  that the widgets JavaScript is still loading. If this message persists, it\n",
       "  likely means that the widgets JavaScript library is either not installed or\n",
       "  not enabled. See the <a href=\"https://ipywidgets.readthedocs.io/en/stable/user_install.html\">Jupyter\n",
       "  Widgets Documentation</a> for setup instructions.\n",
       "</p>\n",
       "<p>\n",
       "  If you're reading this message in another frontend (for example, a static\n",
       "  rendering on GitHub or <a href=\"https://nbviewer.jupyter.org/\">NBViewer</a>),\n",
       "  it may mean that your frontend doesn't currently support widgets.\n",
       "</p>\n"
      ],
      "text/plain": [
       "Button(description='Click to clean image', style=ButtonStyle())"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from skimage import morphology\n",
    "\n",
    "\n",
    "def remove_small_objects(button):\n",
    "    from IPython import display as dsp\n",
    "    dsp.clear_output(wait=True)\n",
    "    display(button)\n",
    "    global image\n",
    "    global threshold_value\n",
    "    global image_cleaned\n",
    "    image_cleaned = morphology.remove_small_objects((image < threshold_value), 50)\n",
    "    fig = plt.figure(figsize=(9, 9))\n",
    "    ax = fig.add_subplot(1,1,1)\n",
    "    ax.set_title(\"Cleaned Image\")\n",
    "    \n",
    "    # fig, ax = plt.subplots(figsize=(12, 8))\n",
    "    ax.imshow(image_cleaned, cmap=plt.\n",
    "              cm.gray, interpolation='nearest')\n",
    "    ax.axis('off')\n",
    "    \n",
    "    #ax.set_adjustable('box-forced')\n",
    "\n",
    "\n",
    "button1 = widgets.Button(description=\"Click to clean image\")\n",
    "display(button1)\n",
    "button1.on_click(remove_small_objects)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Segmentation\n",
    "\n",
    "Given the above image, we now find different connected clusters in the image, which we then label with different colours. Looking at the original images, we can expect three clusters: the land, the lake and the islands.\n",
    "\n",
    "Click on the button below to segment the image into these different classes.\n",
    "\n",
    "The left image shows the different classes. If everything has gone to plan, the right image should provide a nice outline of the lake.\n",
    "\n",
    "If you're happy with the outline, move to the next step.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-02-19T14:57:04.267463Z",
     "start_time": "2018-02-19T14:57:03.904013Z"
    },
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "728fb764ded143a8b5e95fecf302231b",
       "version_major": 2,
       "version_minor": 0
      },
      "text/html": [
       "<p>Failed to display Jupyter Widget of type <code>Button</code>.</p>\n",
       "<p>\n",
       "  If you're reading this message in the Jupyter Notebook or JupyterLab Notebook, it may mean\n",
       "  that the widgets JavaScript is still loading. If this message persists, it\n",
       "  likely means that the widgets JavaScript library is either not installed or\n",
       "  not enabled. See the <a href=\"https://ipywidgets.readthedocs.io/en/stable/user_install.html\">Jupyter\n",
       "  Widgets Documentation</a> for setup instructions.\n",
       "</p>\n",
       "<p>\n",
       "  If you're reading this message in another frontend (for example, a static\n",
       "  rendering on GitHub or <a href=\"https://nbviewer.jupyter.org/\">NBViewer</a>),\n",
       "  it may mean that your frontend doesn't currently support widgets.\n",
       "</p>\n"
      ],
      "text/plain": [
       "Button(description='Click to segment', style=ButtonStyle())"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "def segment(button):\n",
    "    from IPython import display as dsp\n",
    "    dsp.clear_output(wait=True)\n",
    "    display(button)\n",
    "    global image\n",
    "    global image_rgb\n",
    "    global image_cleaned\n",
    "    segmentation = ndi.binary_fill_holes(image_cleaned)\n",
    "    labeled_image, _ = ndi.label(segmentation)\n",
    "    \n",
    "    image_label_overlay = label2rgb(labeled_image, image=image_cleaned)\n",
    "    # print(image_label_overlay[:,:,1].sum())  \n",
    "    \n",
    "    fig = plt.figure(figsize=(12, 10))\n",
    "    ax1 = fig.add_subplot(1,2,2)\n",
    "    \n",
    "    ax1.imshow(image_rgb, cmap=plt.cm.gray, interpolation='nearest')\n",
    "    ax1.contour(segmentation, [0.5], linewidths=1.2, colors='r')\n",
    "    ax1.set_title(\"Contour Overlay\")\n",
    "    \n",
    "    ax2 = fig.add_subplot(1,2,1)\n",
    "    ax2.imshow(image_label_overlay, interpolation='nearest')\n",
    "    ax2.set_title(\"Segmentation\")\n",
    "    #plt.show()\n",
    "    #plt.imshow(image_label_overlay, interpolation='nearest')\n",
    "    plt.tight_layout()\n",
    "    \n",
    "    \n",
    "button2 = widgets.Button(description=\"Click to segment\")\n",
    "display(button2)\n",
    "button2.on_click(segment)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Apply to Multiple Images\n",
    "\n",
    "We have satellite images from a number of consecutive years. So let's apply your colour choice and threshold to these images too. Click on the button to generate an animation of these images. Be careful! What worked for the first image may not work for the rest of them! If you're happy with the result, move to the final step. Otherwise, go back and change your colour or threshold choice. \n",
    "\n",
    "**If you make a change, you must execute all the steps again, including removing small objects and the segmentation**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "39799356b0f549ddbd733d740527f453",
       "version_major": 2,
       "version_minor": 0
      },
      "text/html": [
       "<p>Failed to display Jupyter Widget of type <code>Button</code>.</p>\n",
       "<p>\n",
       "  If you're reading this message in the Jupyter Notebook or JupyterLab Notebook, it may mean\n",
       "  that the widgets JavaScript is still loading. If this message persists, it\n",
       "  likely means that the widgets JavaScript library is either not installed or\n",
       "  not enabled. See the <a href=\"https://ipywidgets.readthedocs.io/en/stable/user_install.html\">Jupyter\n",
       "  Widgets Documentation</a> for setup instructions.\n",
       "</p>\n",
       "<p>\n",
       "  If you're reading this message in another frontend (for example, a static\n",
       "  rendering on GitHub or <a href=\"https://nbviewer.jupyter.org/\">NBViewer</a>),\n",
       "  it may mean that your frontend doesn't currently support widgets.\n",
       "</p>\n"
      ],
      "text/plain": [
       "Button(description='Build Animation', style=ButtonStyle())"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from matplotlib import animation  \n",
    "import matplotlib.image as mgimg\n",
    "import imageio\n",
    "from IPython.display import Image\n",
    "\n",
    "def segment_image(year):\n",
    "    \n",
    "    global choice\n",
    "    global threshold_value\n",
    "\n",
    "    file = '{}.png'.format(year)\n",
    "    im_rgb = io.imread(file)\n",
    "    im = _select_bands(im_rgb, choice)\n",
    "    \n",
    "    clean_im = morphology.remove_small_objects((im < threshold_value), 50)\n",
    "    segmentation = ndi.binary_fill_holes(clean_im)\n",
    "    labeled_image, _ = ndi.label(segmentation)\n",
    "    image_label_overlay = label2rgb(labeled_image, image=clean_im)\n",
    "    \n",
    "    plt.clf()\n",
    "    plt.title(str(year))\n",
    "    plt.imshow(im_rgb, cmap=plt.cm.gray, interpolation='nearest')\n",
    "    plt.contour(segmentation, [0.5], linewidths=1.2, colors='r')\n",
    "    fname = \"temp_{}.png\".format(year)\n",
    "    plt.savefig(fname)   \n",
    "    return fname\n",
    "\n",
    "\n",
    "def build_animation(button):\n",
    "    \n",
    "    from IPython import display as dsp\n",
    "    dsp.clear_output(wait=True)\n",
    "    display(button)\n",
    "    \n",
    "    \n",
    "    progress = widgets.IntProgress(\n",
    "    value=0,\n",
    "    min=0,\n",
    "    max=24 + 1,\n",
    "    step=1,\n",
    "    description='Building:',\n",
    "    bar_style='', # 'success', 'info', 'warning', 'danger' or ''\n",
    "    orientation='horizontal'\n",
    "    )\n",
    "    display(progress)\n",
    "    \n",
    "    \n",
    "    #from IPython import display as dsp\n",
    "    #dsp.clear_output(wait=True)\n",
    "    #display(button)\n",
    "        \n",
    "    # segment_image(1984)\n",
    "    #from IPython import display as dsp\n",
    "    #dsp.clear_output(wait=True)\n",
    "    #interact(segment_image, year=year_slider)\n",
    "    files = []\n",
    "    for i in range(1984, 2008):\n",
    "        files.append(segment_image(i))\n",
    "        progress.value += 1\n",
    "    \n",
    "    imageio.mimsave('temp.gif', [imageio.imread(f) for f in files], 'GIF', duration=0.5)\n",
    "    progress.value += 1\n",
    "    display(Image(filename='temp.gif', width=750, height=500))\n",
    "    #ani = animation.ArtistAnimation(fig, ims, interval=1000, blit=True, repeat_delay=1000)\n",
    "    #plt.show()\n",
    "    plt.clf()\n",
    "\n",
    "button3 = widgets.Button(description=\"Build Animation\")\n",
    "display(button3)\n",
    "\n",
    "button3.on_click(build_animation)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Model the Future\n",
    "\n",
    "Once you've settled on your choices, we can now use your model to build a prediction for what will happen to the lake.\n",
    "\n",
    "We do this by first calculating the area of the outlines in the above animation. We then try to fit a quadratic function through these points as best we can. \n",
    "\n",
    "The generated graph shows the calculated areas in blue, and the quadratic function in yellow. We can use the yellow line to estimate the trend in the amount of water in the lake and use it to make predictions.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "84cedf58f00242579e4225d7f0178beb",
       "version_major": 2,
       "version_minor": 0
      },
      "text/html": [
       "<p>Failed to display Jupyter Widget of type <code>Button</code>.</p>\n",
       "<p>\n",
       "  If you're reading this message in the Jupyter Notebook or JupyterLab Notebook, it may mean\n",
       "  that the widgets JavaScript is still loading. If this message persists, it\n",
       "  likely means that the widgets JavaScript library is either not installed or\n",
       "  not enabled. See the <a href=\"https://ipywidgets.readthedocs.io/en/stable/user_install.html\">Jupyter\n",
       "  Widgets Documentation</a> for setup instructions.\n",
       "</p>\n",
       "<p>\n",
       "  If you're reading this message in another frontend (for example, a static\n",
       "  rendering on GitHub or <a href=\"https://nbviewer.jupyter.org/\">NBViewer</a>),\n",
       "  it may mean that your frontend doesn't currently support widgets.\n",
       "</p>\n"
      ],
      "text/plain": [
       "Button(description='Build predictive model', style=ButtonStyle())"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "\n",
    "def calculate_area(im, choice, threshold_value, progress=None):\n",
    "    im_rgb = io.imread(im)\n",
    "    im = _select_bands(im_rgb, choice)\n",
    "   \n",
    "    \n",
    "    clean_im = morphology.remove_small_objects((im < threshold_value), 50)\n",
    "    segmentation = ndi.binary_fill_holes(clean_im)\n",
    "    labeled_image, _ = ndi.label(segmentation)\n",
    "    image_label_overlay = label2rgb(labeled_image, image=clean_im)\n",
    "    if progress is not None:\n",
    "        progress.value += 1\n",
    "    return image_label_overlay[:,:,1].sum()                      \n",
    "    \n",
    "def predict(button):\n",
    "    global choice\n",
    "    global threshold_value\n",
    "    \n",
    "    from IPython import display as dsp\n",
    "    dsp.clear_output(wait=True)\n",
    "    display(button)\n",
    "    \n",
    "    \n",
    "    progress = widgets.IntProgress(\n",
    "    value=0,\n",
    "    min=0,\n",
    "    max=24,\n",
    "    step=1,\n",
    "    description='Building:',\n",
    "    bar_style='', # 'success', 'info', 'warning', 'danger' or ''\n",
    "    orientation='horizontal'\n",
    "    )\n",
    "    display(progress)\n",
    "    \n",
    "    \n",
    "    x = np.array([i for i in range(1984, 2008)])\n",
    "    vals = [calculate_area('{}.png'.format(i), choice, threshold_value, progress) for i in x]\n",
    "\n",
    "    z = np.poly1d(np.polyfit(x, np.array(vals), 2))\n",
    "    #plt.clf()\n",
    "    fig = plt.figure(figsize=(8, 6))\n",
    "    plt.plot(vals)\n",
    "    plt.plot(z(x))\n",
    "    plt.xlabel('Year')\n",
    "    plt.ylabel('Total Area')\n",
    "    plt.title('Lake Area Changes Over Time')\n",
    "    plt.show()\n",
    "    \n",
    "    \n",
    "button4 = widgets.Button(description=\"Build predictive model\")\n",
    "display(button4)\n",
    "button4.on_click(predict)\n",
    "\n",
    "#image_label_overlay[:,:,1].sum()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Visualise the Future\n",
    "\n",
    "Now that we can predict the lake size in the future, we can visualise this. Click on the button below to see what happens to the lake if it continues to follow this trend. We make this prediction by *eroding* the lake at the specified rate."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-02-19T14:57:11.518376Z",
     "start_time": "2018-02-19T14:57:11.511091Z"
    },
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "3ee4baaa97244cc8a7299734853c9463",
       "version_major": 2,
       "version_minor": 0
      },
      "text/html": [
       "<p>Failed to display Jupyter Widget of type <code>Button</code>.</p>\n",
       "<p>\n",
       "  If you're reading this message in the Jupyter Notebook or JupyterLab Notebook, it may mean\n",
       "  that the widgets JavaScript is still loading. If this message persists, it\n",
       "  likely means that the widgets JavaScript library is either not installed or\n",
       "  not enabled. See the <a href=\"https://ipywidgets.readthedocs.io/en/stable/user_install.html\">Jupyter\n",
       "  Widgets Documentation</a> for setup instructions.\n",
       "</p>\n",
       "<p>\n",
       "  If you're reading this message in another frontend (for example, a static\n",
       "  rendering on GitHub or <a href=\"https://nbviewer.jupyter.org/\">NBViewer</a>),\n",
       "  it may mean that your frontend doesn't currently support widgets.\n",
       "</p>\n"
      ],
      "text/plain": [
       "Button(description='Visualise predictions', style=ButtonStyle())"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "\n",
    "def predict_image(future, image_rgb, spectrum, threshold, area):\n",
    "    \"\"\"\n",
    "    Return the segmented image and area\n",
    "    \"\"\"\n",
    "    im = _select_bands(image_rgb, spectrum)\n",
    "    clean_im = morphology.remove_small_objects((im < threshold), 50)\n",
    "        \n",
    "    #plt.imshow(clean_im, cmap=plt.cm.gray, interpolation='nearest')\n",
    "    \n",
    "    \n",
    "    #erode/dilate cleaned image until correct area is achieved\n",
    "    current_area = _calculate_area(clean_im)\n",
    "    current_im = clean_im\n",
    "    while True:\n",
    "        \n",
    "        diff = current_area - area\n",
    "        \n",
    "        if diff > 0:\n",
    "            #erode\n",
    "            new_clean_im = ndi.binary_erosion(current_im)\n",
    "        else:\n",
    "            #dilate\n",
    "            new_clean_im = ndi.binary_dilation(current_im)\n",
    "            \n",
    "        \n",
    "        new_clean_im = morphology.remove_small_objects(new_clean_im, 50)\n",
    "        \n",
    "        #plt.imshow(new_clean_im, cmap=plt.cm.gray, interpolation='nearest') \n",
    "        #plt.show()\n",
    "        \n",
    "        new_area = _calculate_area(new_clean_im)\n",
    "        new_diff = new_area - area\n",
    "        if np.sign(diff) != np.sign(new_diff):\n",
    "            # overshoot. Pick the better of the two\n",
    "            if np.abs(new_diff) < np.abs(diff):\n",
    "                clean_im = new_clean_im\n",
    "            else:\n",
    "                clean_im = current_im\n",
    "            break\n",
    "        else:\n",
    "            # keep going\n",
    "            current_im = new_clean_im\n",
    "            current_area = new_area\n",
    "    \n",
    "    segmentation = ndi.binary_fill_holes(clean_im)\n",
    "    labeled_image, _ = ndi.label(segmentation)\n",
    "    image_label_overlay = label2rgb(labeled_image, image=clean_im)\n",
    "    \n",
    "    plt.clf()\n",
    "    #fig = plt.figure(figsize=(8, 6))\n",
    "    temp = \"s\" if future != 1 else \"\"\n",
    "    plt.title(\"{} year{} in future\".format(future, temp))\n",
    "    plt.imshow(image_rgb, cmap=plt.cm.gray, interpolation='nearest')\n",
    "    plt.contour(segmentation, [0.5], linewidths=1.2, colors='r')\n",
    "    fname = \"_temp_{}.png\".format(future)\n",
    "    plt.savefig(fname)   \n",
    "    return fname, image_label_overlay[:,:,1].sum()               \n",
    "    \n",
    "def predict(button):\n",
    "    global choice\n",
    "    global threshold_value\n",
    "    \n",
    "    from IPython import display as dsp\n",
    "    dsp.clear_output(wait=True)\n",
    "    display(button)   \n",
    "    \n",
    "    \n",
    "    progress2 = widgets.IntProgress(\n",
    "    value=0,\n",
    "    min=0,\n",
    "    max=2030-2008,\n",
    "    step=1,\n",
    "    description='Animating:',\n",
    "    bar_style='', # 'success', 'info', 'warning', 'danger' or ''\n",
    "    orientation='horizontal'\n",
    "    )\n",
    "    display(progress2)\n",
    "    \n",
    "    x = np.array([i for i in range(1984, 2008)])\n",
    "    vals = [calculate_area('{}.png'.format(i), choice, threshold_value) for i in x]\n",
    "    z = np.poly1d(np.polyfit(x, np.array(vals), 2))\n",
    "    \n",
    "    image_rgb = io.imread('2007.png')\n",
    "    files = []\n",
    "    \n",
    "\n",
    "    \n",
    "    for i, year in enumerate(testing_years):\n",
    "        \n",
    "        area = z(year)\n",
    "        if area < 0:\n",
    "            break\n",
    "        x, area = predict_image(i, image_rgb, choice, threshold_value, area)\n",
    "        files.append(x)\n",
    "        progress2.value += 1\n",
    "    progress2.value = 2030-2008\n",
    "    plt.clf()\n",
    "    imageio.mimsave('_temp.gif', [imageio.imread(f) for f in files], 'GIF', duration=0.5)\n",
    "    display(Image(filename='_temp.gif', width=750, height=500))    \n",
    "    \n",
    "button5 = widgets.Button(description=\"Visualise predictions\")\n",
    "display(button5)\n",
    "button5.on_click(predict)\n",
    "\n",
    "#image_label_overlay[:,:,1].sum()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Congrats!\n",
    "\n",
    "You've just used machine learning to building a predictive model on satellite imagery :)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "\n",
       "<script>\n",
       "code_show=true; \n",
       "function code_toggle() {\n",
       " if (code_show){\n",
       " $('div.input').hide();\n",
       " } else {\n",
       " $('div.input').show();\n",
       " }\n",
       " code_show = !code_show\n",
       "} \n",
       "$( document ).ready(code_toggle);\n",
       "</script>\n",
       "The raw code for this IPython notebook is by default hidden for easier reading.\n",
       "You can toggle the code on and off by clicking <a href=\"javascript:code_toggle()\">here</a>."
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "HTML('''\n",
    "<script>\n",
    "code_show=true; \n",
    "function code_toggle() {\n",
    " if (code_show){\n",
    " $('div.input').hide();\n",
    " } else {\n",
    " $('div.input').show();\n",
    " }\n",
    " code_show = !code_show\n",
    "} \n",
    "$( document ).ready(code_toggle);\n",
    "</script>\n",
    "The raw code for this IPython notebook is by default hidden for easier reading.\n",
    "You can toggle the code on and off by clicking <a href=\"javascript:code_toggle()\">here</a>.''')"
   ]
  }
 ],
 "metadata": {
  "hide_input": false,
  "kernelspec": {
   "display_name": "Python 3",
   "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.6.4"
  },
  "latex_envs": {
   "LaTeX_envs_menu_present": true,
   "autocomplete": true,
   "bibliofile": "biblio.bib",
   "cite_by": "apalike",
   "current_citInitial": 1,
   "eqLabelWithNumbers": true,
   "eqNumInitial": 1,
   "hotkeys": {
    "equation": "Ctrl-E",
    "itemize": "Ctrl-I"
   },
   "labels_anchors": false,
   "latex_user_defs": false,
   "report_style_numbering": false,
   "user_envs_cfg": false
  },
  "toc": {
   "colors": {
    "hover_highlight": "#DAA520",
    "navigate_num": "#000000",
    "navigate_text": "#333333",
    "running_highlight": "#FF0000",
    "selected_highlight": "#FFD700",
    "sidebar_border": "#EEEEEE",
    "wrapper_background": "#FFFFFF"
   },
   "moveMenuLeft": true,
   "nav_menu": {
    "height": "12px",
    "width": "252px"
   },
   "navigate_menu": true,
   "number_sections": true,
   "sideBar": true,
   "threshold": 4,
   "toc_cell": false,
   "toc_section_display": "block",
   "toc_window_display": false,
   "widenNotebook": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}