{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "4d03d7c7",
   "metadata": {},
   "source": [
    "# ex1: practical exercises for [Example Sheet 1](https://www.cl.cam.ac.uk/teaching/current/DataSci/ex/ex1.pdf).\n",
    "\n",
    "The example sheet asks you to implement the three classes given below: `PoissonModel`, `PiecewiseLinearModel`, and `StepPeriodModel`. \n",
    "The class skeletons are given, and you should fill in the missing pieces.\n",
    "To test your answers on [Moodle](https://www.cl.cam.ac.uk/teaching/current/DataSci/redirect/ex1moodle.html),\n",
    "please upload either a Jupyter notebook called `ex1.ipynb` or a plain Python file called `ex1.py`. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b04bc7a9",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import scipy.optimize\n",
    "import sklearn.linear_model"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ada21ca5",
   "metadata": {},
   "source": [
    "**Poisson model.** Suppose we're given a dataset $[x_1,\\dots,x_n]$. We wish to fit the model that says each $x_i$ is an independent sample from \n",
    "the Poisson(λ) distribution. Estimate λ using [`scipy.optimize.fmin`](https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.fmin.html#scipy.optimize.fmin).\n",
    "\n",
    "_Note. If the tester reports that your answer is a little bit off, try increasing the precision that `scipy.optimize.fmin` is using by e.g. passing\n",
    "in the argument `xtol=0.00001`._\n",
    "\n",
    "```python\n",
    "class PoissonModel():\n",
    "    def __init__(self):\n",
    "        self.λ_ = np.nan\n",
    "    def fit(self, x):\n",
    "        # Input: x is a numpy vector of integers\n",
    "        # TODO: set self.λ_\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6986f935",
   "metadata": {},
   "source": [
    "**Piecewise linear response.** Suppose we're given a dataset of $(x_i,y_i)$ pairs. We wish to fit a model for $y$ as a function of $x$, made up of two straight lines. The function must be continuous, i.e. the two straight lines must meet at an inflection point. The $x$-coordinate of the inflection point is given.\n",
    "<br>\n",
    "<img src=\"res/ex1_piecewiselinear.png\" style=\"height:12em\"/>\n",
    "\n",
    "```python\n",
    "class PiecewiseLinearModel():\n",
    "    def fit(self, x, y, inflection_x):\n",
    "        # Input: x and y are numpy vectors of real numbers, inflection_x is a real number\n",
    "        # TODO: fit the model, and store its parameters\n",
    "    def predict(self, x):\n",
    "        # Input: x is a numpy vector of real numbers\n",
    "        # TODO: return a numpy vector of real numbers, with the predicted y values\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c5ea24a1",
   "metadata": {},
   "source": [
    "**Stepwise climate model.** We're given a time series consisting of $(t_i,\\text{temp}_i)$ pairs, recording average daily temperatures in &deg;C. Here time is measured in years, and readings are monthly. We wish to fit a model which describes the temperature as a sinusoid plus a step response function, \n",
    "$$\n",
    "\\text{temp} \\approx \\beta_1 \\sin(2\\pi t)+\\beta_2 \\cos(2\\pi t) + ???\n",
    "$$\n",
    "where the last term tells us the decadal average temperature. (Take decades to be represented by their start year, `np.floor(t/10)*10`.)\n",
    "<br>\n",
    "<img src=\"res/ex1_climatestep.png\" style=\"height:12em\"/>\n",
    "\n",
    "```python\n",
    "class StepPeriodicModel():\n",
    "    def fit(self, t, temp):\n",
    "        # Input: x and y are numpy vectors of real numbers\n",
    "        # TODO: fit the model, and store its parameters\n",
    "    def predict_step(self, t):\n",
    "        # Input: t is a numpy vector of real numbers\n",
    "        # TODO: return a numpy vector of real numbers, with the predicted decadal average temperatures\n",
    "        # (i.e. the ??? part of the equation above, not including the sinusoidal part).\n",
    "        # It should return np.nan for timepoints outside the range where we have data.\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "08a1794c",
   "metadata": {},
   "source": [
    "# TEST\n",
    "_The Moodle checker will look for a markdown cell with the contents `# TEST`, and ignore everything beneath it. Put your working code above this cell, and put any experiments and tests below._"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7157490a",
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c6f51e9a",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Test the Poisson model (where we know the true answer from algebra)\n",
    "\n",
    "x = [3,2,8,1,5,0,8]\n",
    "m = PoissonModel()\n",
    "m.fit(x)\n",
    "assert np.isclose(m.λ_, np.mean(x), rtol=1e-3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c2453089",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot the piecewise linear model\n",
    "\n",
    "df = pandas.read_csv('https://www.cl.cam.ac.uk/teaching/current/DataSci/ticks/res/piecewiselinear.csv')\n",
    "\n",
    "m = PiecewiseLinearModel()\n",
    "m.fit(df.x, df.y, 2.8)\n",
    "\n",
    "fig,ax = plt.subplots(figsize=(3,2))\n",
    "ax.scatter(df.x, df.y, alpha=.5)\n",
    "xnew = np.linspace(0,5,100)\n",
    "ax.plot(xnew, m.predict(xnew), color='black', linestyle='--')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "dd06c3cf",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot the stepwise climate model\n",
    "\n",
    "climate = pandas.read_csv('https://www.cl.cam.ac.uk/teaching/current/DataSci/ticks/res/climate_202610.csv')\n",
    "df = climate.loc[(climate.station=='Oxford') & (~pandas.isna(climate.temp))]\n",
    "\n",
    "m = StepPeriodicModel()\n",
    "m.fit(df.t, df.temp)\n",
    "\n",
    "# For Oxford, t=1855, the expected answer is 9.55417.\n",
    "# If your code produces an answer around 3.039, you've probably included the sinusoid in your prediction,\n",
    "# which is not what's wanted. You are meant to fit a model that includes the sinusoid,\n",
    "# but your predict_step function should only return the non-sinusoid part of the fitted model.\n",
    "# If your answer is off by a little bit, perhaps you've specified an unidentifiable model\n",
    "# (i.e. linearly dependent features), causing sklearn.linear_model to lose precision.\n",
    "test_years = [1855,1965]\n",
    "pred = m.predict_step(np.array(test_years))\n",
    "print(pred)\n",
    "assert np.allclose(pred, [9.55417, 9.855]), f\"Answers for {test_years} not what was expected\"\n",
    "\n",
    "fig,ax = plt.subplots(figsize=(5,2.5))\n",
    "tnew = np.arange(1850,2030,10)\n",
    "ax.step(tnew, m.predict_step(tnew), where='post')\n",
    "plt.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.10.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}