{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Delayed Nonlinear Negative Feedback Oscillator\n",
    "\n",
    "Implement the 3-state model of Brian Goodwin (B.C. Goodwin, Oscillatory behavior in enzymatic control processes, Adv. Enzym. Regul. 3 (1965) 425e437).\n",
    "\n",
    "The parameters are\n",
    "- alpha_1: maximal rate of X's transcription\n",
    "- K : activation threshold for Z on X\n",
    "- n : Hill coefficient\n",
    "- gamma_1: rate constant for degradation of X\n",
    "- alpha_2: rate constant for X causing production of Y\n",
    "- gamma_2: rate constant for degradation of Y\n",
    "- alpha_3: rate constant for Y causing production of Z\n",
    "- gamma_3: rate constant for degradation of Z\n",
    "\n",
    "The model is\n",
    "\\begin{eqnarray*}\n",
    "\\frac{dX}{dt} &=& \\alpha_1 \\frac{K^n}{K^n+Z^n} - \\gamma_1 X\\\\\n",
    "\\frac{dY}{dt} &=& \\alpha_2 X - \\gamma_2 Y\\\\\n",
    "\\frac{dZ}{dt} &=& \\alpha_3 Y - \\gamma_3 Z \\\\\n",
    "\\end{eqnarray*}\n",
    "\n",
    "(Equations as presented by D Gonze & P Ruoff, The Goodwin Oscillator and its Legacy, Acta Biotheoretica, 2020)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import scipy.integrate\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def goodwin_oscillator( t, statevars, params ):\n",
    "    X = statevars[0];\n",
    "    Y = statevars[1];\n",
    "    Z = statevars[2];\n",
    "    # set parameters here\n",
    "    alpha_1 = params[0]\n",
    "    K = params[1]\n",
    "    n = params[2]\n",
    "    gamma_1 = params[3]\n",
    "    alpha_2 = params[4]\n",
    "    gamma_2 = params[5]\n",
    "    alpha_3 = params[6]\n",
    "    gamma_3 = params[7]\n",
    "    \n",
    "    dydt = np.zeros( (3,) )\n",
    "    # Write code here\n",
    "    \n",
    "    return dydt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Simulate the model\n",
    "\n",
    "Simulate the model with parameters from Figure 2 (with n modified) from B. Ananthasubramaniam, C. Schmal and H. Herzel, Amplitude Effects Allow Short Jet Lags and Large Seasonal Phase Shifts in Minimal Clock Models..., Journal of Molecular Biology.\n",
    "\n",
    "- alpha_1 = 1\n",
    "- K = 0.61\n",
    "- n = 11\n",
    "- gamma_1 = 0.15\n",
    "- alpha_2 = 1\n",
    "- gamma_2 = 0.15\n",
    "- alpha_3 = 1\n",
    "- gamma_3 = 0.25\n",
    "\n",
    "We will assume the model reaches its limit cycle within 200 hours of simulation. Simulate the model first for 200 hours. Then start a new simulation from the end time point of the first simulation and simulate for 100 hours."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Simulate here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Plot the simulation\n",
    "Plot both simulations in separate subplots. The first should show the transient (the path it takes to reach the oscillation's limit cycle). The second should show the oscillation on its limit cycle)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "#plot here\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualize the transient in X-Y-Z space\n",
    "\n",
    "Examine the transient as it approaches the limit cycle. Visualize the limit cycle in X-Y-Z space."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Write code here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### What are the effects of n on the oscillations?\n",
    "Plot the simulation for n = 5,7,9,11,13"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Write code here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Examine the activation/repression over time\n",
    "\n",
    "Choose a parameter set of interest (e.g. the original one), and determine how K^n/(K^n+Z^n) changes over time"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Write code here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### What are the effects of K on the oscillations?\n",
    "Plot the simulation for several values of K ranging from 0.1 to 1. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Write code here"
   ]
  }
 ],
 "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.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}