linear-models.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "from sklearn.linear_model import LogisticRegression"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Linear models\n",
    "\n",
    "## Resources\n",
    "\n",
    "Vincent Warmerdam: Winning with Simple, even Linear, Models PyData London 2018 - [youtube](https://www.youtube.com/watch?v=68ABAU_V8qI) - [slides](https://koaning.io/theme/notebooks/simple-models.pdf)\n",
    "\n",
    "## Linear models\n",
    "\n",
    "Good\n",
    "- simple\n",
    "- interpretable\n",
    "\n",
    "Bad\n",
    "- only linear relationships between feature & target\n",
    "- no relationship between multiple features (such as if this, then that etc)\n",
    "\n",
    "When to use linear models\n",
    "- linear relationship between features & target\n",
    "- when you have good features\n",
    "- interpretability is important\n",
    "\n",
    "Requirement of no-collinearity\n",
    "- collinearity = correlation between features\n",
    "- this cause instability of parameters -> can't interpret the parameters\n",
    "\n",
    "We can use PCA with `whiten=True` to remove this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.decomposition import PCA\n",
    "\n",
    "samples = 100\n",
    "data = np.random.uniform(0, 1, size=samples)\n",
    "\n",
    "dataset = np.vstack([\n",
    "    data, \n",
    "    data + np.random.uniform(0, 0.5, size=samples),\n",
    "    data - np.random.uniform(0, 0.5, size=samples)\n",
    "]).T\n",
    "\n",
    "assert dataset.shape[0] == samples"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The covariance matrix shows collinearity:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.corrcoef(dataset, rowvar=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can use PCA to remove this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tr = PCA(whiten=True)\n",
    "clean = tr.fit_transform(dataset)\n",
    "np.corrcoef(clean, rowvar=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Non-linear relationships with linear models\n",
    "\n",
    "Clever feature engineering can allow you to solve non-linear problems with a linear model.\n",
    "\n",
    "First create a dataset:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "num_pts = 100\n",
    "\n",
    "red_x = np.concatenate([np.random.normal(3, 0.5, num_pts), np.random.normal(-3, 0.5, num_pts)])\n",
    "red_y = np.concatenate([np.random.normal(3, 0.5, num_pts), np.random.normal(-3, 0.5, num_pts)])\n",
    "\n",
    "blue_x = np.concatenate([np.random.normal(0, 1, num_pts), np.random.normal(0, 1, num_pts)])\n",
    "blue_y = np.concatenate([np.random.normal(0, 1, num_pts), np.random.normal(0, 1, num_pts)])\n",
    "\n",
    "plt.scatter(red_x, red_y, color='red')\n",
    "plt.scatter(blue_x, blue_y, color='blue')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can show that a linear model works poorly:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "features = np.zeros((num_pts*4, 2))\n",
    "target = np.zeros((num_pts*4))\n",
    "\n",
    "for idx, (x, y) in enumerate(zip(red_x, red_y)):\n",
    "    features[idx, :] = x, y\n",
    "    target[idx] = 0\n",
    "    \n",
    "for idx, (x, y) in enumerate(zip(blue_x, blue_y), idx):\n",
    "    features[idx, :] = x, y\n",
    "    target[idx] = 1\n",
    "    \n",
    "mdl = LogisticRegression(solver='lbfgs')\n",
    "\n",
    "mdl.fit(features, target)\n",
    "\n",
    "mdl.score(features, target)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's add a new feature:\n",
    "- x1 * x2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "features = np.concatenate(\n",
    "    [features, \n",
    "     (features[:, 0] * features[:, 1]).reshape(-1, 1)], \n",
    "    axis=1\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now that we add this feature, our linear classifier does well:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mdl = LogisticRegression(solver='lbfgs')\n",
    "\n",
    "mdl.fit(features, target)\n",
    "\n",
    "mdl.score(features, target)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Radial basis functions\n",
    "\n",
    "[Wikipedia](https://en.wikipedia.org/wiki/Radial_basis_function)\n",
    "\n",
    "Class of functions where the value depends only on the distance between the input and a fixed point\n",
    "\n",
    "There are a number of different kernels - we will focus on the Gaussian RBF kernel:\n",
    "\n",
    "$$ r(x) = e^{-\\epsilon (x-x_{0})^2}$$\n",
    "\n",
    "where\n",
    "- $\\epsilon$ is a hyperparameter\n",
    "- $x_{0}$ is the center\n",
    "\n",
    "Let's try this on a real time series dataset - describes the minimum daily temperatures over 10 years (1981-1990) in the city Melbourne, Australia [from ML-mastery](https://machinelearningmastery.com/time-series-datasets-for-machine-learning/).\n",
    "\n",
    "We will do year ahead prediction"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = pd.read_csv(\n",
    "    'https://raw.githubusercontent.com/jbrownlee/Datasets/master/daily-min-temperatures.csv', \n",
    "    index_col=0, parse_dates=True\n",
    ")\n",
    "\n",
    "data = data.iloc[-6*365:, :]\n",
    "\n",
    "ax = data.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A typical feature would be to use the month.\n",
    "\n",
    "Let's train a simple linear model, with the month one-hot encoded:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y = data.loc[:, 'Temp']\n",
    "from sklearn.preprocessing import OneHotEncoder\n",
    "\n",
    "from sklearn.linear_model import Ridge\n",
    "enc = OneHotEncoder(sparse=False)\n",
    "x = enc.fit_transform(data.index.month.values.reshape(-1, 1))\n",
    "\n",
    "mdl = Ridge()\n",
    "mdl.fit(x, y)\n",
    "\n",
    "y_p = mdl.predict(x)\n",
    "\n",
    "plt.plot(y)\n",
    "plt.plot(y.index, y_p)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.metrics.pairwise import rbf_kernel"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "days = np.arange(1, 365).reshape(-1, 1)\n",
    "\n",
    "k = rbf_kernel(days, days, gamma=0.1)\n",
    "print(k.shape)\n",
    "plt.plot(k[128])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rbf_data = data.copy()\n",
    "\n",
    "days = rbf_data.index.dayofyear.values.reshape(-1, 1)\n",
    "k = rbf_kernel(days, days, gamma=0.1)\n",
    "plt.plot(k[5])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for day in range(1, 365):\n",
    "    rbf_data.loc[:, 'rbf{}'.format(day)] = k[day]\n",
    "    \n",
    "rbf_data.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y = rbf_data.loc[:, 'Temp']\n",
    "x = rbf_data.drop('Temp', axis=1)\n",
    "\n",
    "mdl = Ridge()\n",
    "mdl.fit(x, y)\n",
    "\n",
    "y_p_r = mdl.predict(x)\n",
    "\n",
    "f, a = plt.subplots(nrows=2, sharex=True, figsize=(20, 10))\n",
    "\n",
    "a[0].plot(y, label='true', color='lightblue')\n",
    "a[1].plot(y, label='true', color='lightblue')\n",
    "\n",
    "a[0].plot(y.index, y_p, label='one-hot', color='black')\n",
    "a[1].plot(y.index, y_p_r, label='rbf', color='black')\n",
    "\n",
    "for ax in a:\n",
    "    ax.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Feature selection with linear models\n",
    "\n",
    "As we have a linear model, we can check to see which features are significant (note we need to use `statsmodels` to get this):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import statsmodels.api as sm\n",
    "\n",
    "ols = sm.OLS(y, x)\n",
    "ols_result = ols.fit()\n",
    "#ols_result.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Conversion of coordinates\n",
    "\n",
    "class1 = np."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.datasets import make_circles\n",
    "\n",
    "x, y = make_circles(factor=0.5)\n",
    "\n",
    "plt.scatter(x[:, 0], x[:, 1], c=y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mdl = LogisticRegression(solver='lbfgs')\n",
    "mdl.fit(x, y)\n",
    "mdl.score(x, y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "https://stackoverflow.com/questions/20924085/python-conversion-between-coordinates"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def cartesian_to_polar(x, y):\n",
    "    rho = np.sqrt(x**2 + y**2).reshape(-1, 1)\n",
    "    phi = np.arctan2(y, x).reshape(-1, 1)\n",
    "    return np.hstack([rho, phi])\n",
    "\n",
    "xp = cartesian_to_polar(x[:, 0], x[:, 1])\n",
    "\n",
    "mdl = LogisticRegression(solver='lbfgs')\n",
    "mdl.fit(xp, y)\n",
    "mdl.score(xp, y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.scatter(xp[:, 0], xp[:, 1], c=y)"
   ]
  },
  {
   "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.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}