05a numpy and interact.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Ordinary python lists"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "[3,4,5,\"Bob\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mylist = [2,5,7,\"Alice\"]\n",
    "mylist"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mylist[3]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mylist[4]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Building lists\n",
    "\n",
    "With a `for` loop:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "evennumbers = []\n",
    "for i in range(200):\n",
    "    if i%2 == 0:\n",
    "        evennumbers.append(i)\n",
    "len(evennumbers)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With a `while` loop"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "evennumbers = [0]\n",
    "while len(evennumbers) < 100:\n",
    "    evennumbers.append(evennumbers[-1]+2)\n",
    "len(evennumbers)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Another `for` loop"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "evennumbers = []\n",
    "for i in range(100):\n",
    "    evennumbers.append(2*i)\n",
    "\n",
    "    \n",
    "# We can slice lists\n",
    "evennumbers[3:7]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With a *list comprehension*:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "evennumbers = [2*i for i in range(100)]\n",
    "evennumbers[-4:-1]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With a different list comprehension:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "evennumbers = [i for i in range(200) if i%2 == 0]\n",
    "evennumbers[50:-48]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sum(evennumbers) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Nested clauses"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "v = 0\n",
    "for e in evennumbers:\n",
    "    #print(\"This will happen every time\")\n",
    "    if e > 50:\n",
    "        #print(\"This will happen fifty times\")\n",
    "        if e%7 == 0:\n",
    "            #print(\"This will happen only for even multiples of 7 bigger than 50\")\n",
    "            v += 1\n",
    "v"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## About numpy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A `numpy` `array` is just a dressed-up list:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "evennumbers = np.arange(0,200,2)\n",
    "evennumbers"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With some extra functionality:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "evennumbers.mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "evennumbers[evennumbers < 10]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### and pandas"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`pandas` just dresses `numpy` in some more functionality:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "a = np.array([1,2,3])\n",
    "b = np.array([\"bob\",\"gene\",\"tina\"])\n",
    "df = pd.DataFrame(data = {\n",
    "    'number': a,\n",
    "    'name': b\n",
    "})\n",
    "df"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Back to numpy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's generate a normally distributed population:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "population = np.random.normal(loc=50, scale=3, size=10000)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`loc` is $\\mu$, `scale` is $\\sigma$, `size` is $N$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Take a sample from this population:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sample = np.random.choice(population, size=100, replace=False)\n",
    "\n",
    "print(\"Population mean = {}\".format(population.mean().round(2)))\n",
    "print(\"Sample mean = {}\".format(sample.mean().round(2)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import seaborn as sns\n",
    "import matplotlib.pyplot as plt\n",
    "import warnings\n",
    "warnings.simplefilter(action='ignore', category=FutureWarning)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, axs = plt.subplots()\n",
    "sns.distplot(population, ax = axs)\n",
    "axs.axvline(population.mean())\n",
    "sns.distplot(sample, ax = axs)\n",
    "axs.axvline(sample.mean())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The interact widget"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ipywidgets import interact"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def update(n):\n",
    "    print(n)\n",
    "\n",
    "interact(update,n=(1,100))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The continuous updating is annoying so..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ipywidgets import IntSlider"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def update(n):\n",
    "    print(n)\n",
    "\n",
    "interact(update,n=IntSlider(min=1, max=100, step=1, continuous_update=False))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "\n",
    "def update(n):\n",
    "    fig, axs = plt.subplots()\n",
    "    sns.distplot(population, ax = axs)\n",
    "    axs.set_ylim(0,0.2)\n",
    "    axs.set_xlim(30,70)\n",
    "    axs.axvline(population.mean())\n",
    "    sample = np.random.choice(population, size=n, replace=False)\n",
    "    sns.distplot(sample, ax = axs)\n",
    "    axs.axvline(sample.mean())\n",
    "\n",
    "interact(update,n=IntSlider(value=10, min=2, max=1000, step=1, continuous_update=False))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So the sample mean $\\bar{x}$ is a good (unbiased) estimator for the population mean $\\mu$.\n",
    "\n",
    "The same is **not** true of the standard deviation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "def update(n):\n",
    "    # make a set of axes\n",
    "    fig, axs = plt.subplots()\n",
    "    # set the axes\n",
    "    axs.set_ylim(0,0.2)\n",
    "    axs.set_xlim(30,70)\n",
    "    # plot the population\n",
    "    sns.distplot(population, ax = axs)\n",
    "    # find the mean and sd for the population\n",
    "    mu = population.mean()\n",
    "    sigma = population.std()\n",
    "    # color one sd from the mean\n",
    "    axs.axvspan(mu-sigma, mu+sigma, facecolor=\"lightsteelblue\", alpha=0.4)\n",
    "    # draw a sample\n",
    "    sample = np.random.choice(population, size=n, replace=False)\n",
    "    # plot the sample\n",
    "    sns.distplot(sample, ax = axs)\n",
    "    # find the sample mean and the 'wrong' sd\n",
    "    xbar = sample.mean()\n",
    "    s = sample.std()\n",
    "    axs.axvspan(xbar-s,xbar+s, facecolor=\"wheat\", alpha=0.2)\n",
    "    print(\"Standard deviation of sample with /n = {}\".format(s))\n",
    "    print(\"Population standard deviation = {}\".format(sigma))\n",
    "\n",
    "\n",
    "interact(update,n=IntSlider(value=10, min=2, max=1000, step=1, continuous_update=False))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**This is not quite as compelling a visualisation as we wanted.**\n",
    "\n",
    "We want to show that using $\\dfrac{\\sum(x-\\bar{x})^2}{n}$ tends to *underestimate* the standard deviation of the population, which is why we use $n-1$ instead.\n",
    "\n",
    "But of course, at these sample sizes, the difference between dividing by $n$ and dividing by $n-1$ is not going to be visible on a graph."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Instead, we could draw a lot of samples and plot the distribution of the uncorrected and corrected standard deviations compared to the true value."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "hide_input": false
   },
   "outputs": [],
   "source": [
    "# Our population\n",
    "population = np.random.normal(loc=50, scale=5, size=10000)\n",
    "# Make 1000 samples, each of size 5\n",
    "samples = [np.random.choice(population, size=5, replace=False) for _ in range(1000)]\n",
    "# Make an array of the means of each sample\n",
    "sample_means = np.array([sample.mean() for sample in samples])\n",
    "# And their (uncorrected) standard deviations\n",
    "sample_stds =np.array([sample.std() for sample in samples])\n",
    "# Amd their (corrected) standar deviations\n",
    "sample_stds_corrected = np.array([sample.std(ddof=1) for sample in samples])\n",
    "# We want three sets of axes\n",
    "fig, axs = plt.subplots(1,3)\n",
    "fig.set_figwidth(18)\n",
    "# plot the distribution of means around the true mean\n",
    "axs[0].set_title(\"Mean\")\n",
    "axs[0].get_yaxis().set_visible(False)\n",
    "sns.distplot(sample_means, ax=axs[0], color=\"goldenrod\")\n",
    "axs[0].axvline(population.mean(), color=\"steelblue\")\n",
    "axs[0].axvline(sample_means.mean(), color=\"orange\")\n",
    "# and the distribution of uncorrected standard deviations around the true standard deviation\n",
    "axs[1].set_title(\"Uncorrected standard deviation\")\n",
    "axs[1].get_yaxis().set_visible(False)\n",
    "sns.distplot(sample_stds, ax=axs[1], color=\"goldenrod\")\n",
    "axs[1].axvline(population.std(), color=\"steelblue\")\n",
    "axs[1].axvline(sample_stds.mean(), color=\"orange\")\n",
    "# and the distribution of corrected standard deviations around the true standard deviation\n",
    "axs[2].set_title(\"Corrected standard deviation\")\n",
    "axs[2].get_yaxis().set_visible(False)\n",
    "sns.distplot(sample_stds_corrected, ax=axs[2], color=\"goldenrod\")\n",
    "axs[2].axvline(population.std(), color=\"steelblue\")\n",
    "axs[2].axvline(sample_stds_corrected.mean(), color=\"orange\")"
   ]
  }
 ],
 "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}