04-WOFOST_playground.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "toc": true
   },
   "source": [
    "<h1>Table of Contents<span class=\"tocSkip\"></span></h1>\n",
    "<div class=\"toc\"><ul class=\"toc-item\"><li><span><a href=\"#Introduction\" data-toc-modified-id=\"Introduction-1\"><span class=\"toc-item-num\">1&nbsp;&nbsp;</span>Introduction</a></span></li><li><span><a href=\"#A-simulation-with-WOFOST\" data-toc-modified-id=\"A-simulation-with-WOFOST-2\"><span class=\"toc-item-num\">2&nbsp;&nbsp;</span>A simulation with WOFOST</a></span><ul class=\"toc-item\"><li><span><a href=\"#Your-first-WOFOST-run\" data-toc-modified-id=\"Your-first-WOFOST-run-2.1\"><span class=\"toc-item-num\">2.1&nbsp;&nbsp;</span>Your first WOFOST run</a></span></li></ul></li><li><span><a href=\"#Some-important-WOFOST-parameters\" data-toc-modified-id=\"Some-important-WOFOST-parameters-3\"><span class=\"toc-item-num\">3&nbsp;&nbsp;</span>Some important WOFOST parameters</a></span></li><li><span><a href=\"#Remarks\" data-toc-modified-id=\"Remarks-4\"><span class=\"toc-item-num\">4&nbsp;&nbsp;</span>Remarks</a></span></li></ul></div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div style=\"float:right\">\n",
    "<table>\n",
    "<tr>\n",
    "    <td> <img src=\"figs/nceo_logo.png\" alt=\"Drawing\" style=\"width:200px;height:40px;\"/> </td>\n",
    "    <td> <img src=\"figs/gssti_logo.png\" alt=\"Drawing\" style=\"width:200px;height:40px;\"/> </td>\n",
    "    <td> <img src=\"figs/multiply_logo.png\" alt=\"Drawing\" style=\"width:40px;height:40px;\"/> </td>\n",
    "</tr>\n",
    "    </table>\n",
    "</div>\n",
    "\n",
    "# The WOFOST crop model\n",
    "\n",
    "**Author:** Jose Gómez-Dans (NCEO & UCL)  `j.gomez-dans@ucl.ac.uk`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "This notebook was built with:\n",
      "python version: 3.6.7 | packaged by conda-forge | (default, Jul  2 2019, 02:18:42) \n",
      "[GCC 7.3.0] \n",
      "PCSE version: 5.4.2\n"
     ]
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "%load_ext autoreload\n",
    "%autoreload 2\n",
    "import sys, os\n",
    "import copy\n",
    "import datetime as dt\n",
    "import pcse\n",
    "from pcse.models import Wofost71_WLP_FD\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib\n",
    "import matplotlib.pyplot as plt\n",
    "print(\"This notebook was built with:\")\n",
    "print(\"python version: %s \" % sys.version)\n",
    "print(\"PCSE version: %s\" %  pcse.__version__)\n",
    "import datetime\n",
    "\n",
    "from pandas.plotting import register_matplotlib_converters\n",
    "register_matplotlib_converters()\n",
    "\n",
    "from wofost_utils import *"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introduction\n",
    "\n",
    "The WOFOST crop model has been around for a long time. It is widely used in both research and operational applications, and is fairly easy to use. It also has been used to simulate a large number of crops, and there are a number of extensions that deal with other crops.\n",
    "\n",
    "While WOFOST lacks some of the more complex processes that more sophisticated models provide, it's more than adequate to familiarise oneself with a crop model.\n",
    "\n",
    "The following diagram shows what WOFOST does\n",
    "\n",
    "![schematic of the WOFOST model](figs/wofost_schematic.png)\n",
    "\n",
    "Some of those boxes will already be familiar to you. We can see at the top that WOFOST is driven by radiation intercepted by the canopy. The canopy interception itself is not driven by a satellite fAPAR product as we saw earlier, but rather uses the predictions of LAI from the previous time step. Once photosynthesis is calculated and both growth and maintenance respiration are lost, partitioning of the assimilated C to the different plant organs takes place. These organs are:\n",
    "\n",
    "1. Leaves\n",
    "2. Stems\n",
    "3. Roots\n",
    "4. Storage organs (e.g. grain!)\n",
    "\n",
    "A number of processes deal with the age and eventual death of leaves, stems and roots.\n",
    "\n",
    "WOFOST implements a large number of processes to model a crop. Although not shown here, it does have a soil moisture module, and many of the processes are controlled by the development stage of the crop (its *phenology*). The phenology is simulated often by thermal time (there are other options, using e.g. daylength and stuff).\n",
    "\n",
    "WOFOST doesn't implement many processes that are to do with management, pests, complex soil chemistry, etc. However, some management and additional information is needed to run the model, such as selecting the crop and variety, providing meteo data, configuring things around... This is all done for you (lucky you!). Since our aim here is to build towards linking the model with EO data, we'll consider some parameters that are important to diagnose what's happening. These will include LAI, soil moisture, total storage organ weight, and the development stage. \n",
    "\n",
    "\n",
    "## A simulation with WOFOST\n",
    "\n",
    "There are two ways of using the model: potential and water limited. the former ignores the effect of water stress, whereas the latter does some soil water calculations and includes the effect of soil moisture in crop development. \n",
    "\n",
    "A simulation starts by the user selecting a soil type, a crop variety, some meteorological inputs, and some other parameters such as either the sowing or emergence date, the amount of water in the soil at the start of the simulation, etc. Once this is done, the model takes over, calculating soil hydraulics and any limitations to growth due to lack of soil water, and calculating radiation interception, photosynthesis, removing respiration losses and eventually partioning the assimilate. Defined by two crop variety parameters, the model calculates the thermal time between emergence and flowering, and between flowering and harvest, this is the development stage (DVS), which goes from 0 (emergence) to 1 (flowering) to 2 (harvest). Many processes are impacted by DVS.\n",
    "\n",
    "### Your first WOFOST run\n",
    "\n",
    "The next grid will run WOFOST and plot some graphs of some state variables with respect to time as a function of sowing and harvest date. The code will run the simulation changing the sowing date from sowing to harvest date, every `n_days`. You can select different regions and/or years."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "c92a5e2c89f4421e971bf06efd4ae41c",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "interactive(children=(DatePicker(value=datetime.date(2011, 7, 1), description='Earliest possible sowing date')…"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "change_sowing_slider();"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can probably see some relationships emerging between some parameters. You can also see that changing the sowing date results in massively different trajectories for LAI and biomass allocation to the different organs, but that ultimately, the variation in yield seems to be quite small. A reason for this is that once the crop reaches maturity (DVS=2), there's not assimilation going on. Also, as the leaves age and senesce, there's isn't that much C going into the system, which produces a saturation pattern.\n",
    "\n",
    "The example above is a gross simplification. We have given you the ability to change the meteorological inputs, but we still use the same crop parameters, same soil properties and same starting values for everything. Changing these settings would results in different trajectories. Let's have a play with some important parameters\n",
    "\n",
    "## Some important WOFOST parameters\n",
    "\n",
    "WOFOST has got tons of parameters, and it would be unwise to just list all of them. We will play around with just a handful to see the effect they have on crop behaviour. Here's a list of them.\n",
    "\n",
    "* `WAV`: Initial amount of available water in total rootable zone (moisture content above wilting point) (cm)\n",
    "* `RDMSOL`: Maximum rooting depth allowed by the soil\n",
    "* `SPAN`: life span of leaves growing at an average temperature of 35 °C\n",
    "* `TDWI`: Initial total dry weight of crop\n",
    "* `RGRLAI`: Relative growth in LAI per day\n",
    "* `TSUMEM`: Length (in effective degree days) from sowing to emergence \n",
    "* `TSUM1`: Length (in effective degree days) of the vegetative period\n",
    "* `TSUM2`: Length (in effective degree days) of the reproductive period\n",
    "* `CVO`: Conversion efficiency of C into storage organs\n",
    "* `CVL`: Conversion efficienc of C into leaves\n",
    "\n",
    "`TSUM1` and `TSUM2` are critical in controlling phenology. `SPAN` and `RGRLAI` tend to control the LAI dynamics, whereas `TDWI` and `TSUMEM` control the crop emergence and that first phase of the crop development. A good first starting point is ensuring that your LAI curve looks sensible before you try to tweak other parameters, in this case by ensuring that there's enough water in the soil. Once you get a sensible LAI curve, see the effet that the other parameters have on the different organs etc."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "e1d16e2e54a44c4b9396db6433c272c1",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "interactive(children=(FloatSlider(value=40.0, description='span', max=50.0, min=20.0), FloatSlider(value=20.0,…"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "results = {}\n",
    "wofost_parameter_sweep();\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It is often a bit tricky to see the differences in so many parameters. Each run dumps the state above into a CSV file, and the following code will plot the different files parameter by parameter so that you can spot the changes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "2d7a1cfb6e414bc898e77b2f044df4f4",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "interactive(children=(Dropdown(description='state', options=('DVS', 'LAI', 'TAGP', 'TWSO', 'TWLV', 'TWST', 'TW…"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "WOFOST_STATE = ['DVS', 'LAI', 'TAGP', 'TWSO', 'TWLV', 'TWST',\n",
    "                'TWRT', 'TRA', 'RD', 'SM']\n",
    "\n",
    "@widgets.interact(state=widgets.Dropdown(options=WOFOST_STATE))\n",
    "def plot_results(state):\n",
    "    plt.figure(figsize=(16,8))\n",
    "    results = {f.name:pd.read_csv(f) for f in Path().cwd().glob(\"*-???.csv\")}\n",
    "    for k, v in results.items():\n",
    "        \n",
    "        plt.plot_date(v.index, v[state], '-', label=k[:-4].replace(\"_\", \" \"))\n",
    "    plt.legend(loc=\"best\")\n",
    "    plt.ylabel(state)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Remarks\n",
    "\n",
    "The above examples show that many different model trajectories are possible. Some of them will be broadly within your expectations of how the crop behaves, some others will just be crazy. Which begs an important question... **What are the *correct* parameter values? And how do we select them?** \n",
    "\n",
    "Before we address that important thought, observe too that we have a lot of parameters (>100), and from your tinkering, you're probably aware that many parameter combinations may result in the same model behaviour. The fancy word for that is **equifinality**.\n",
    "\n",
    "So in order to find the best parameters, we tend to use field measurements to **calibrate** the model: e.g. we fit the model to the data using the measurements. This is an non-linear fitting exercise, which is hard to do with such a large parameter space. Using prior distributions for the parameters, and only considering the most sensitive parameters to our output fields of interest is the way to go."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a rel=\"license\" href=\"http://creativecommons.org/licenses/by-nc/4.0/\"><img alt=\"Creative Commons Licence\" style=\"border-width:0\" src=\"https://i.creativecommons.org/l/by-nc/4.0/88x31.png\" /></a><br />This work is licensed under a <a rel=\"license\" href=\"http://creativecommons.org/licenses/by-nc/4.0/\">Creative Commons Attribution-NonCommercial 4.0 International License</a>."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "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.7"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": true,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": true,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}