ch13-code-listing.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Chapter 13: Statistics"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Robert Johansson\n",
    "\n",
    "Source code listings for [Numerical Python - A Practical Techniques Approach for Industry](http://www.apress.com/9781484205549) (ISBN 978-1-484205-54-9).\n",
    "\n",
    "The source code listings can be downloaded from http://www.apress.com/9781484205549"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Imports"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from scipy import stats"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from scipy import optimize"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import random"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import seaborn as sns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "sns.set(style=\"whitegrid\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Descriptive statistics"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "x = np.array([3.5, 1.1, 3.2, 2.8, 6.7, 4.4, 0.9, 2.2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "np.mean(x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "np.median(x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "x.min(), x.max()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "x.var()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "x.std()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "x.var(ddof=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "x.std(ddof=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Random numbers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "random.seed(123456789)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "random.random()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "random.randint(0, 10)  # 0 and 10 inclusive"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "np.random.seed(123456789)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "np.random.rand()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "np.random.randn()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "np.random.rand(5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "np.random.randn(2, 4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "np.random.randint(10, size=10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "np.random.randint(low=10, high=20, size=(2, 10))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "fig, axes = plt.subplots(1, 3, figsize=(12, 3))\n",
    "\n",
    "axes[0].hist(np.random.rand(10000))\n",
    "axes[0].set_title(\"rand\")\n",
    "axes[1].hist(np.random.randn(10000))\n",
    "axes[1].set_title(\"randn\")\n",
    "axes[2].hist(np.random.randint(low=1, high=10, size=10000), bins=9, align='left')\n",
    "axes[2].set_title(\"randint(low=1, high=10)\")\n",
    "\n",
    "fig.tight_layout()\n",
    "fig.savefig(\"ch13-random-hist.pdf\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#random.sample(range(10), 5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "np.random.choice(10, 5, replace=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "np.random.seed(123456789)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "np.random.rand()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "np.random.seed(123456789); np.random.rand()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "np.random.seed(123456789); np.random.rand()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "prng = np.random.RandomState(123456789)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "prng.rand(2, 4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "prng.chisquare(1, size=(2, 2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "prng.standard_t(1, size=(2, 3))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "prng.f(5, 2, size=(2, 4))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "prng.binomial(10, 0.5, size=10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "prng.poisson(5, size=10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Probability distributions and random variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "np.random.seed(123456789)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X = stats.norm(1, 0.5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X.mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X.median()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X.std()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "X.var()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "[X.moment(n) for n in range(5)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X.stats()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X.pdf([0, 1, 2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X.cdf([0, 1, 2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X.rvs(10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "stats.norm(1, 0.5).stats()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "stats.norm.stats(loc=2, scale=0.5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X.interval(0.95)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X.interval(0.99)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def plot_rv_distribution(X, axes=None):\n",
    "    \"\"\"Plot the PDF, CDF, SF and PPF of a given random variable\"\"\"\n",
    "    if axes is None:\n",
    "        fig, axes = plt.subplots(1, 3, figsize=(12, 3))\n",
    "    \n",
    "    x_min_999, x_max_999 = X.interval(0.999)\n",
    "    x999 = np.linspace(x_min_999, x_max_999, 1000)\n",
    "\n",
    "    x_min_95, x_max_95 = X.interval(0.95)\n",
    "    x95 = np.linspace(x_min_95, x_max_95, 1000)\n",
    "\n",
    "    if hasattr(X.dist, 'pdf'):\n",
    "        axes[0].plot(x999, X.pdf(x999), label=\"PDF\")\n",
    "        axes[0].fill_between(x95, X.pdf(x95), alpha=0.25)\n",
    "    else:\n",
    "        x999_int = np.unique(x999.astype(int))\n",
    "        axes[0].bar(x999_int, X.pmf(x999_int), label=\"PMF\")\n",
    "    axes[1].plot(x999, X.cdf(x999), label=\"CDF\")\n",
    "    axes[1].plot(x999, X.sf(x999), label=\"SF\")\n",
    "    axes[2].plot(x999, X.ppf(x999), label=\"PPF\")\n",
    "    \n",
    "    for ax in axes:\n",
    "        ax.legend()\n",
    "    \n",
    "    return axes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "fig, axes = plt.subplots(3, 3, figsize=(12, 9))\n",
    "\n",
    "X = stats.norm()\n",
    "plot_rv_distribution(X, axes=axes[0, :])\n",
    "axes[0, 0].set_ylabel(\"Normal dist.\")\n",
    "X = stats.f(2, 50)\n",
    "plot_rv_distribution(X, axes=axes[1, :])\n",
    "axes[1, 0].set_ylabel(\"F dist.\")\n",
    "X = stats.poisson(5)\n",
    "plot_rv_distribution(X, axes=axes[2, :])\n",
    "axes[2, 0].set_ylabel(\"Poisson dist.\")\n",
    "\n",
    "fig.tight_layout()\n",
    "fig.savefig(\"ch13-distributions.pdf\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def plot_dist_samples(X, X_samples, title=None, ax=None):\n",
    "    \"\"\" Plot the PDF and histogram of samples of a continuous random variable\"\"\"\n",
    "    if ax is None:\n",
    "        fig, ax = plt.subplots(1, 1, figsize=(8, 4))\n",
    "\n",
    "    x_lim = X.interval(.99)\n",
    "    x = np.linspace(*x_lim, num=100)\n",
    "\n",
    "    ax.plot(x, X.pdf(x), label=\"PDF\", lw=3)    \n",
    "    ax.hist(X_samples, label=\"samples\", normed=1, bins=75)\n",
    "    ax.set_xlim(*x_lim)\n",
    "    ax.legend()\n",
    "    \n",
    "    if title:\n",
    "        ax.set_title(title)\n",
    "    return ax"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "fig, axes = plt.subplots(1, 3, figsize=(12, 3))\n",
    "X = stats.t(7.0)\n",
    "plot_dist_samples(X, X.rvs(2000), \"Student's t dist.\", ax=axes[0])\n",
    "X = stats.chi2(5.0)\n",
    "plot_dist_samples(X, X.rvs(2000), r\"$\\chi^2$ dist.\", ax=axes[1])\n",
    "X = stats.expon(0.5)\n",
    "plot_dist_samples(X, X.rvs(2000), \"exponential dist.\", ax=axes[2])\n",
    "fig.tight_layout()\n",
    "fig.savefig(\"ch13-dist-sample.pdf\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X = stats.chi2(df=5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X_samples = X.rvs(500)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "df, loc, scale = stats.chi2.fit(X_samples)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "df, loc, scale"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Y = stats.chi2(df=df, loc=loc, scale=scale)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(8, 3))\n",
    "\n",
    "x_lim = X.interval(.99)\n",
    "x = np.linspace(*x_lim, num=100)\n",
    "\n",
    "ax.plot(x, X.pdf(x), label=\"original\")\n",
    "ax.plot(x, Y.pdf(x), label=\"recreated\")\n",
    "ax.legend()\n",
    "\n",
    "fig.tight_layout()\n",
    "fig.savefig(\"ch13-max-likelihood-fit.pdf\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "fig, axes = plt.subplots(1, 2, figsize=(12, 4))\n",
    "\n",
    "x_lim = X.interval(.99)\n",
    "x = np.linspace(*x_lim, num=100)\n",
    "\n",
    "axes[0].plot(x, X.pdf(x), label=\"original\")\n",
    "axes[0].plot(x, Y.pdf(x), label=\"recreated\")\n",
    "axes[0].legend()\n",
    "\n",
    "axes[1].plot(x, X.pdf(x) - Y.pdf(x), label=\"error\")\n",
    "axes[1].legend()\n",
    "\n",
    "fig.tight_layout()\n",
    "fig.savefig(\"ch13-max-likelihood-fit.pdf\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Hypothesis testing"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "np.random.seed(123456789)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "mu, sigma = 1.0, 0.5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "X = stats.norm(mu-0.2, sigma)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "n = 100"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "X_samples = X.rvs(n)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "z = (X_samples.mean() - mu)/(sigma/np.sqrt(n))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "z"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "t = (X_samples.mean() - mu)/(X_samples.std(ddof=1)/np.sqrt(n))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "t"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "stats.norm().ppf(0.025)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "2 * stats.norm().cdf(-abs(z))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "2 * stats.t(df=(n-1)).cdf(-abs(t))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "t, p = stats.ttest_1samp(X_samples, mu)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "t"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "p"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(8, 3))\n",
    "\n",
    "sns.distplot(X_samples, ax=ax)\n",
    "x = np.linspace(*X.interval(0.999), num=100)\n",
    "ax.plot(x, stats.norm(loc=mu, scale=sigma).pdf(x))\n",
    "\n",
    "fig.tight_layout()\n",
    "fig.savefig(\"ch13-hypothesis-test-dist-sample-mean.pdf\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "n = 50"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "mu1, mu2 = np.random.rand(2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "X1 = stats.norm(mu1, sigma)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "X1_sample = X1.rvs(n)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "X2 = stats.norm(mu2, sigma)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "X2_sample = X2.rvs(n)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "t, p = stats.ttest_ind(X1_sample, X2_sample)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "t"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "p"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "mu1, mu2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sns.distplot(X1_sample)\n",
    "sns.distplot(X2_sample)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Nonparameteric methods"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "np.random.seed(0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X = stats.chi2(df=5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X_samples = X.rvs(100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "kde = stats.kde.gaussian_kde(X_samples)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "kde_low_bw = stats.kde.gaussian_kde(X_samples, bw_method=0.25)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "x = np.linspace(0, 20, 100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "fig, axes = plt.subplots(1, 3, figsize=(12, 3))\n",
    "\n",
    "axes[0].hist(X_samples, normed=True, alpha=0.5, bins=25)\n",
    "axes[1].plot(x, kde(x), label=\"KDE\")\n",
    "axes[1].plot(x, kde_low_bw(x), label=\"KDE (low bw)\")\n",
    "axes[1].plot(x, X.pdf(x), label=\"True PDF\")\n",
    "axes[1].legend()\n",
    "sns.distplot(X_samples, bins=25, ax=axes[2])\n",
    "\n",
    "fig.tight_layout()\n",
    "fig.savefig(\"ch13-hist-kde.pdf\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "kde.resample(10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def _kde_cdf(x):\n",
    "    return kde.integrate_box_1d(-np.inf, x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "kde_cdf = np.vectorize(_kde_cdf)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(8, 3))\n",
    "\n",
    "sns.distplot(X_samples, bins=25, ax=ax)\n",
    "x = np.linspace(0, 20, 100)\n",
    "ax.plot(x, kde_cdf(x))\n",
    "\n",
    "fig.tight_layout()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def _kde_ppf(q):\n",
    "    return optimize.fsolve(lambda x, q: kde_cdf(x) - q, kde.dataset.mean(), args=(q,))[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "kde_ppf = np.vectorize(_kde_ppf)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "kde_ppf([0.05, 0.95])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X.ppf([0.05, 0.95])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Versions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%reload_ext version_information"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%version_information numpy, scipy, matplotlib, seaborn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "lines = \"\"\"alpha             -- Alpha\n",
    "       anglit            -- Anglit\n",
    "       arcsine           -- Arcsine\n",
    "       argus             -- Argus\n",
    "       beta              -- Beta\n",
    "       betaprime         -- Beta Prime\n",
    "       bradford          -- Bradford\n",
    "       burr              -- Burr (Type III)\n",
    "       burr12            -- Burr (Type XII)\n",
    "       cauchy            -- Cauchy\n",
    "       chi               -- Chi\n",
    "       chi2              -- Chi-squared\n",
    "       cosine            -- Cosine\n",
    "       crystalball       -- Crystalball\n",
    "       dgamma            -- Double Gamma\n",
    "       dweibull          -- Double Weibull\n",
    "       erlang            -- Erlang\n",
    "       expon             -- Exponential\n",
    "       exponnorm         -- Exponentially Modified Normal\n",
    "       exponweib         -- Exponentiated Weibull\n",
    "       exponpow          -- Exponential Power\n",
    "       f                 -- F (Snecdor F)\n",
    "       fatiguelife       -- Fatigue Life (Birnbaum-Saunders)\n",
    "       fisk              -- Fisk\n",
    "       foldcauchy        -- Folded Cauchy\n",
    "       foldnorm          -- Folded Normal\n",
    "       frechet_r         -- Deprecated. Alias for weibull_min\n",
    "       frechet_l         -- Deprecated. Alias for weibull_max\n",
    "       genlogistic       -- Generalized Logistic\n",
    "       gennorm           -- Generalized normal\n",
    "       genpareto         -- Generalized Pareto\n",
    "       genexpon          -- Generalized Exponential\n",
    "       genextreme        -- Generalized Extreme Value\n",
    "       gausshyper        -- Gauss Hypergeometric\n",
    "       gamma             -- Gamma\n",
    "       gengamma          -- Generalized gamma\n",
    "       genhalflogistic   -- Generalized Half Logistic\n",
    "       gilbrat           -- Gilbrat\n",
    "       gompertz          -- Gompertz (Truncated Gumbel)\n",
    "       gumbel_r          -- Right Sided Gumbel, Log-Weibull, Fisher-Tippett, Extreme Value Type I\n",
    "       gumbel_l          -- Left Sided Gumbel, etc.\n",
    "       halfcauchy        -- Half Cauchy\n",
    "       halflogistic      -- Half Logistic\n",
    "       halfnorm          -- Half Normal\n",
    "       halfgennorm       -- Generalized Half Normal\n",
    "       hypsecant         -- Hyperbolic Secant\n",
    "       invgamma          -- Inverse Gamma\n",
    "       invgauss          -- Inverse Gaussian\n",
    "       invweibull        -- Inverse Weibull\n",
    "       johnsonsb         -- Johnson SB\n",
    "       johnsonsu         -- Johnson SU\n",
    "       kappa4            -- Kappa 4 parameter\n",
    "       kappa3            -- Kappa 3 parameter\n",
    "       ksone             -- Kolmogorov-Smirnov one-sided (no stats)\n",
    "       kstwobign         -- Kolmogorov-Smirnov two-sided test for Large N (no stats)\n",
    "       laplace           -- Laplace\n",
    "       levy              -- Levy\n",
    "       levy_l\n",
    "       levy_stable\n",
    "       logistic          -- Logistic\n",
    "       loggamma          -- Log-Gamma\n",
    "       loglaplace        -- Log-Laplace (Log Double Exponential)\n",
    "       lognorm           -- Log-Normal\n",
    "       lomax             -- Lomax (Pareto of the second kind)\n",
    "       maxwell           -- Maxwell\n",
    "       mielke            -- Mielke's Beta-Kappa\n",
    "       moyal             -- Moyal\n",
    "       nakagami          -- Nakagami\n",
    "       ncx2              -- Non-central chi-squared\n",
    "       ncf               -- Non-central F\n",
    "       nct               -- Non-central Student's T\n",
    "       norm              -- Normal (Gaussian)\n",
    "       norminvgauss      -- Normal Inverse Gaussian\n",
    "       pareto            -- Pareto\n",
    "       pearson3          -- Pearson type III\n",
    "       powerlaw          -- Power-function\n",
    "       powerlognorm      -- Power log normal\n",
    "       powernorm         -- Power normal\n",
    "       rdist             -- R-distribution\n",
    "       reciprocal        -- Reciprocal\n",
    "       rayleigh          -- Rayleigh\n",
    "       rice              -- Rice\n",
    "       recipinvgauss     -- Reciprocal Inverse Gaussian\n",
    "       semicircular      -- Semicircular\n",
    "       skewnorm          -- Skew normal\n",
    "       t                 -- Student's T\n",
    "       trapz              -- Trapezoidal\n",
    "       triang            -- Triangular\n",
    "       truncexpon        -- Truncated Exponential\n",
    "       truncnorm         -- Truncated Normal\n",
    "       tukeylambda       -- Tukey-Lambda\n",
    "       uniform           -- Uniform\n",
    "       vonmises          -- Von-Mises (Circular)\n",
    "       vonmises_line     -- Von-Mises (Line)\n",
    "       wald              -- Wald\n",
    "       weibull_min       -- Minimum Weibull (see Frechet)\n",
    "       weibull_max       -- Maximum Weibull (see Frechet)\n",
    "       wrapcauchy        -- Wrapped Cauchy\"\"\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "len(lines.split(\"\\n\"))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "help(stats)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3.6 (NPM ed2)",
   "language": "python",
   "name": "py3.6-npm-e2"
  },
  "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
}