season-2/cadcad-hack-21-profiling-models.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Table of Contents\n",
    "\n",
    "<ul>\n",
    "    <li><a href='#Dependencies'>Dependencies</a></li>\n",
    "    </li>\n",
    "    <li><a href='#Modelling'>Modelling</a>\n",
    "        <ol style='margin-top: 0em;'>\n",
    "            <li><a href='#1.-State-Variables'>State Variables</a></li>\n",
    "            <li><a href='#2.-System-Parameters'>System Parameters</a></li>\n",
    "            <li><a href='#3.-Policy-Functions'>Policy Functions</a></li>\n",
    "            <li><a href='#4.-State-Update-Functions'>State Update Functions</a></li>\n",
    "            <li><a href='#5.-Partial-State-Update-Blocks'>Partial State Update Blocks</a></li>\n",
    "        </ol>\n",
    "    </li>\n",
    "    <li><a href='#Simulation'>Simulation</a>\n",
    "        <ul style='margin-top: 0em; start=\"6\"'>\n",
    "            <li><a href='#6,-7,-and-8:-Configuration,-Execution-and-Output-Preparation'>6, 7, and 8: Configuration, Execution and Output Preparation</a></li>\n",
    "            <li><a href='#9.-Analysis'>Analysis</a></li>\n",
    "        </ul>\n",
    "    </li>\n",
    "</ul>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Dependencies"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# No need for cadCAD standard dependencies\n",
    "# we're going to use the cadCAD_tools library\n",
    "# Install it through `pip install cadCAD-tools`\n",
    "\n",
    "from cadCAD_tools import profile_run"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Additional dependencies\n",
    "\n",
    "# For visualizing how much each substep is impacting the simulation\n",
    "from cadCAD_tools.profiling.visualizations import visualize_substep_impact\n",
    "\n",
    "# For visualizing how much elapsed time takes for each timestep\n",
    "from cadCAD_tools.profiling.visualizations import visualize_elapsed_time_per_ts\n",
    "\n",
    "# For analytics\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "# For visualization\n",
    "import plotly.express as px"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Modelling"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. State Variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "initial_state = {\n",
    "    'predator_population': 15, # Initial number of predators\n",
    "    'prey_population': 100, # Initial number of preys\n",
    "}\n",
    "initial_state"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. System Parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "system_params = {\n",
    "    # Parameters describing the interaction between the two populations:\n",
    "    \n",
    "    # A parameter used to calculate the rate of predator birth\n",
    "    \"predator_birth_parameter\": [0.01], \n",
    "    # A parameter used to calculate the rate of predator death\n",
    "    \"predator_death_parameter\": [1.0],\n",
    "    # A parameter used to calculate the rate of prey birth\n",
    "    \"prey_birth_parameter\": [0.6, 1.0],\n",
    "    # A parameter used to calculate the rate of prey death\n",
    "    \"prey_death_parameter\": [0.03],\n",
    "\n",
    "    # Parameters used for random Numpy variable tuning\n",
    "    # These parameters scale the random variable used to create the random birth / death rate\n",
    "    \n",
    "    \"random_predator_birth\": [0.0002],\n",
    "    \"random_predator_death\": [0.005],\n",
    "    \"random_prey_birth\": [0.1],\n",
    "    \"random_prey_death\": [0.1],\n",
    "\n",
    "    # Parameter used as conversion factor between 1 unit of time and 1 timestep\n",
    "    # 10 timesteps == 1 unit of time, i.e. every 10 cadCAD model timesteps, 1 unit of actual model time passes\n",
    "    \n",
    "    \"dt\": [0.1],\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Policy Functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def p_predator_births(params, substep, state_history, previous_state):\n",
    "    '''Predator Births Policy Function\n",
    "    The predator birth rate (rate of predators born per unit of time) is a product of\n",
    "    the prey population and the predator birth parameter plus a random variable.\n",
    "    \n",
    "    i.e. the larger the prey population, the higher the predator birth rate\n",
    "    '''\n",
    "    # Parameters   \n",
    "    dt = params['dt']\n",
    "    predator_birth_parameter = params['predator_birth_parameter']\n",
    "    random_predator_birth = params['random_predator_birth']\n",
    "\n",
    "    # State Variables\n",
    "    predator_population = previous_state['predator_population']\n",
    "    prey_population = previous_state['prey_population']\n",
    "\n",
    "    # Calculate the predator birth rate\n",
    "    birth_rate = prey_population * (predator_birth_parameter + np.random.random() * random_predator_birth)\n",
    "    # Calculate change in predator population\n",
    "    births = birth_rate * predator_population * dt\n",
    "\n",
    "    return {'add_to_predator_population': births}\n",
    "\n",
    "\n",
    "def p_predator_deaths(params, substep, state_history, previous_state):\n",
    "    '''Predator Deaths Policy Function\n",
    "    The predator death rate (rate of predators that die per unit of time) is a function of\n",
    "    the predator death parameter plus a random variable.\n",
    "    \n",
    "    i.e. the larger the predator death parameter, the higher the predator death rate\n",
    "    '''\n",
    "    # Parameters\n",
    "    dt = params['dt']\n",
    "    predator_death_parameter = params['predator_death_parameter']\n",
    "    random_predator_death = params['random_predator_death']\n",
    "\n",
    "    # State Variables\n",
    "    predator_population = previous_state['predator_population']\n",
    "\n",
    "    # Calculate the predator death rate\n",
    "    death_rate = predator_death_parameter + np.random.random() * random_predator_death\n",
    "    # Calculate change in predator population\n",
    "    deaths = death_rate * predator_population * dt\n",
    "\n",
    "    return {'add_to_predator_population': -1.0 * deaths}\n",
    "\n",
    "\n",
    "def p_prey_births(params, substep, state_history, previous_state):\n",
    "    '''Prey Births Policy Function\n",
    "    The prey birth rate (rate of preys born per unit of time) is a function of\n",
    "    the prey birth parameter plus a random variable.\n",
    "    \n",
    "    i.e. the larger the prey birth parameter, the higher the prey birth rate\n",
    "    '''\n",
    "    # Parameters\n",
    "    dt = params['dt']\n",
    "    prey_birth_parameter = params['prey_birth_parameter']\n",
    "    random_prey_birth = params['random_prey_birth']\n",
    "\n",
    "    # State Variables\n",
    "    prey_population = previous_state['prey_population']\n",
    "\n",
    "    # Calculate the prey birth rate\n",
    "    birth_rate = prey_birth_parameter + np.random.random() * random_prey_birth\n",
    "    # Calculate change in prey population\n",
    "    births = birth_rate * prey_population * dt\n",
    "\n",
    "    return {'add_to_prey_population': births}\n",
    "\n",
    "\n",
    "def p_prey_deaths(params, substep, state_history, previous_state):\n",
    "    '''Prey Deaths Policy Function\n",
    "    The prey death rate (rate of preys that die per unit of time) is a product of\n",
    "    the predator population and the prey death parameter plus a random variable.\n",
    "    \n",
    "    i.e. the larger the predator population, the higher the prey death rate\n",
    "    '''\n",
    "    # Parameters\n",
    "    dt = params['dt']\n",
    "    prey_death_parameter = params['prey_death_parameter']\n",
    "    random_prey_death = params['random_prey_death']\n",
    "\n",
    "    # State Variables\n",
    "    prey_population = previous_state['prey_population']\n",
    "    predator_population = previous_state['predator_population']\n",
    "\n",
    "    # Calculate the prey death rate\n",
    "    death_rate = predator_population * (prey_death_parameter + np.random.random() * random_prey_death)\n",
    "    # Calculate change in prey population\n",
    "    deaths = death_rate * prey_population * dt\n",
    "\n",
    "    return {'add_to_prey_population': -1.0 * deaths}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. State Update Functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def s_predator_population(params, substep, state_history, previous_state, policy_input):\n",
    "    '''Predator Population State Update Function\n",
    "    Take the Policy Input `add_to_predator_population`\n",
    "    (the net predator births and deaths)\n",
    "    and add to the `predator_population` State Variable.\n",
    "    '''\n",
    "    # Policy Inputs\n",
    "    add_to_predator_population = policy_input['add_to_predator_population']\n",
    "\n",
    "    # State Variables\n",
    "    predator_population = previous_state['predator_population']\n",
    "\n",
    "    # Calculate updated predator population\n",
    "    updated_predator_population = predator_population + add_to_predator_population\n",
    "\n",
    "    return 'predator_population', updated_predator_population\n",
    "\n",
    "\n",
    "def s_prey_population(params, substep, state_history, previous_state, policy_input):\n",
    "    '''Prey Population State Update Function\n",
    "    Take the Policy Input `add_to_prey_population`\n",
    "    (the net prey births and deaths)\n",
    "    and add to the `prey_population` State Variable.\n",
    "    '''\n",
    "    # Policy Inputs\n",
    "    add_to_prey_population = policy_input['add_to_prey_population']\n",
    "\n",
    "    # State Variables\n",
    "    prey_population = previous_state['prey_population']\n",
    "\n",
    "    # Calculate updated prey population\n",
    "    updated_prey_population = prey_population + add_to_prey_population\n",
    "\n",
    "    return 'prey_population', updated_prey_population"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. Partial State Update Blocks"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "partial_state_update_blocks = [\n",
    "    {   \n",
    "        'label': \"Lotka-Volterra equations\",\n",
    "        # Configure the model Policy Functions\n",
    "        'policies': {\n",
    "            # Calculate the predator birth rate and number of births\n",
    "            'predator_births': p_predator_births,\n",
    "            # Calculate the predator death rate and number of deaths\n",
    "            'predator_deaths': p_predator_deaths,\n",
    "            # Calculate the prey birth rate and number of births\n",
    "            'prey_births': p_prey_births,\n",
    "            # Calculate the prey death rate and number of deaths\n",
    "            'prey_deaths': p_prey_deaths,\n",
    "        },\n",
    "        # Configure the model State Update Functions\n",
    "        'variables': {\n",
    "            # Update the predator population\n",
    "            'predator_population': s_prey_population,\n",
    "            # Update the prey population\n",
    "            'prey_population': s_predator_population\n",
    "        }\n",
    "    },\n",
    "    {\n",
    "        'label': 'Do nothing',\n",
    "        'policies': {},\n",
    "        'variables': {}\n",
    "    }\n",
    "]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Simulation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6, 7, and 8: Configuration, Execution and Output Preparation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df = profile_run(initial_state,\n",
    "                 system_params,\n",
    "                 partial_state_update_blocks,\n",
    "                 20, # 20 timesteps\n",
    "                 5, # 5 MC runs\n",
    "                 use_label=True, # Assign PSUB labels to the results df\n",
    "                 assign_params=True # Have the parameters on the results df\n",
    "                )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 9. Analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# See how much each substep is slowing the sim down.\n",
    "# Pass relative=True to see the fraction\n",
    "# relative=False will give the absolute time\n",
    "visualize_substep_impact(df, relative=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# See how much time has elapsed for each timestep\n",
    "# Pass relative=True to see how much time it took for just one timestep\n",
    "# relative=False will give the total elapsed time\n",
    "visualize_elapsed_time_per_ts(df, relative=False)"
   ]
  }
 ],
 "metadata": {
  "jupytext": {
   "formats": "ipynb,py:percent"
  },
  "kernelspec": {
   "display_name": "Python (cadCAD Hacks)",
   "language": "python",
   "name": "python-cadcad-hacks"
  },
  "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.8.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}