Python_Notebook_Tutorial.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introduction\n",
    "\n",
    "This document produces a short, minimalist worked-example of an Analysis Notebook documenting the process from data entry to production of publication-quality figures.\n",
    "\n",
    "## Hypothesis and design\n",
    "A good analysis follows any *a priori* hypotheses and follows from an appropriate experimental design. In this case, Darwin's hypothesis was that self-pollination leads to inbreeding depression. Accordingly, the design compares the fitness (using plant height as a surrogate) of pairs of self- and cross-pollinated plants. \n",
    "\n",
    "Adding SMPracticals darwin dataset from file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "df = pd.read_csv(\"SMPractical_darwindataset.csv\")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exploratory graphics"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "plt.scatter(df[\"type\"],df[\"height\"])\n",
    "plt.margins(1)\n",
    "plt.grid()\n",
    "plt.xlabel(\"type\",fontsize=12)\n",
    "plt.ylabel(\"height\",fontsize=12)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Analysis\n",
    "Linear model analysis:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                            OLS Regression Results                            \n",
      "==============================================================================\n",
      "Dep. Variable:                 height   R-squared:                       0.175\n",
      "Model:                            OLS   Adj. R-squared:                  0.146\n",
      "Method:                 Least Squares   F-statistic:                     5.940\n",
      "Date:                Thu, 23 May 2019   Prob (F-statistic):             0.0214\n",
      "Time:                        12:27:51   Log-Likelihood:                -73.889\n",
      "No. Observations:                  30   AIC:                             151.8\n",
      "Df Residuals:                      28   BIC:                             154.6\n",
      "Df Model:                           1                                         \n",
      "Covariance Type:            nonrobust                                         \n",
      "================================================================================\n",
      "                   coef    std err          t      P>|t|      [0.025      0.975]\n",
      "--------------------------------------------------------------------------------\n",
      "Intercept       20.1917      0.759     26.596      0.000      18.637      21.747\n",
      "type[T.Self]    -2.6167      1.074     -2.437      0.021      -4.816      -0.417\n",
      "==============================================================================\n",
      "Omnibus:                       16.447   Durbin-Watson:                   2.106\n",
      "Prob(Omnibus):                  0.000   Jarque-Bera (JB):               18.083\n",
      "Skew:                          -1.568   Prob(JB):                     0.000118\n",
      "Kurtosis:                       5.152   Cond. No.                         2.62\n",
      "==============================================================================\n",
      "\n",
      "Warnings:\n",
      "[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.\n"
     ]
    }
   ],
   "source": [
    "from statsmodels.formula.api import ols\n",
    "linmodel = ols(formula=\"height ~ type\", data=df)\n",
    "result = linmodel.fit()\n",
    "print(result.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Analysis of variance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "            df      sum_sq    mean_sq         F    PR(>F)\n",
      "type       1.0   51.352083  51.352083  5.939518  0.021414\n",
      "Residual  28.0  242.083333   8.645833       NaN       NaN\n"
     ]
    }
   ],
   "source": [
    "from statsmodels.stats import anova\n",
    "table = anova.anova_lm(result)\n",
    "print(table)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Assumptions\n",
    "Diagnostic plot of residuals vs. fitted values:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYcAAAENCAYAAADkNanAAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi40LCBodHRwOi8vbWF0cGxvdGxpYi5vcmcv7US4rQAAGV9JREFUeJzt3X2UXHWd5/H31wAaVyAwCUoa2gAqHo8ZDbaA4uiirCjiEnHXEcXxgTWjs7sjOoMSYWfXh1kcmXF1j7vjRPGsszCoI5hRh3Mi+DRnHIN0DBh5CILy1KCE0cjDZAXDd/+o26STW9V9O13V996u9+ucOqn+1a1b35vqqk/f+/vd343MRJKkqR5XdwGSpOYxHCRJJYaDJKnEcJAklRgOkqQSw0GSVGI4SJJKDAdJUonhIEkq2afuAvbW0qVLc8WKFXWXIUmtsmnTpvsyc9lMy7U2HFasWMH4+HjdZUhSq0TE7VWW87CSJKnEcJAklTQqHCJiUURsjoiv1V2LJA2zRoUD8C7gxrqLkKRh15hwiIjDgFcBn6m7Fkkado0JB+DjwHuBR+suRJKGXSPCISJOBe7NzE0zLLcmIsYjYnzbtm3zVJ0kDZ9GhANwAvBvI+I24PPASyPi4j0Xysx1mTmWmWPLls14DockaS81Ihwyc21mHpaZK4DXA9/MzDNrLkuShlYjwkGS1CyNmz4jM78NfLvmMiRpqLnnIEkqMRwkSSWGgySppHF9DpK00K3fPMGFG7Zy9/YdLF+ymHNOPprVq0bqLms3hoMkzaP1mydYe/kWdjyyE4CJ7TtYe/kWgEYFhIeVJGkeXbhh62PBMGnHIzu5cMPWmirqznCQpHl09/Yds2qvi+EgSfNo+ZLFs2qvi+EgSfPoxGd2nxeuV3tdDAdJmkffuqn7jNK92utiOEjSPGpLn4NDWStow5hkSe2wfMliJroEgX0OLTM5Jnli+w6SXWOS12+eqLs0SS1kn8MC0ZYxyZLawT6HBaItxwcltUNbvlMMhxm0ZUyypHZoy3eK4TCDc04+msX7LtqtbfG+izjn5KNrqkhSm7XlO8VwmMHqVSMcM3rgbm3HjB7oaCVJe2X1qhEuOH0lI0sWE8DIksVccPrKxn2nRGbWXcNeGRsby/Hx8YG/zvnrt3DxxjtK7WceP8qHV68c+OtLUj9FxKbMHJtpOfccZnDp1XfOql2SFgLDYQY7e+xZ9WqXpIXAcJjBoohZtUvSQmA4zOCM4w6fVbskLQSNCIeIODwivhURN0TE9RHxrrprmvTh1Ss58/jRx/YUFkXYGS1pwWvEaKWIOBQ4NDN/EBH7A5uA1Zl5Q6/nzNdoJUnqt/PXb+HSq+9kZyaLIjjjuMPn7Q/OVo1Wysx7MvMHxf0HgBuBZg36laQ+mBwePzmoZWcmF2+8g/PXb6m5st01IhymiogVwCrg6i6PrYmI8YgY37atWZNUSVIVbRke36hwiIgnAZcBZ2fm/Xs+npnrMnMsM8eWLWvW9LaSVEVbhsc3JhwiYl86wXBJZl5edz2SNMwaEQ4REcBFwI2Z+bG665GkYdeIcABOAN4EvDQiri1up9RdlCT120iPqbl7tdelEeGQmf+YmZGZv52Zzy1uV9RdlyT1m1N2S5JKVq8a4bXPG9ntxNrXPm+kcVN2Gw6SNI/Wb57gsk0Tu53ncNmmCdZvnqi5st0ZDpI0jy7csJUdj+zcrW3HIzu5cMPWmirqznCQpHk0sX3HrNrrYjhI0jxqy2UADAdJmkdtOUN6n7oLaIM6Z1CUtLCMLFnc9RCS5zm0TFtmUJTUDp7nsEC0ZQZFSe2wetUIF5y+kpEliwk6ewwXnL6ycec5eFhpBm05PiipPVavat5Jb3syHGawKKJrEDRtZIGk9mhDP6aHlWZw/JEHzapdkqbTln5Mw2EGt/1z9xNTerVL0nTa0o9pOMzg7h5nLfZql6TptKUf03CYwfIeY497tUvSdDxDeoFoy5hkSe1wxnGHz6q9LobDDFavGuGwg56wW9thBz2h8cPQJDXTh1ev5MzjR3e7nsOZx482brRSZMOOc1U1NjaW4+PjA3+dN376e3z31l+U2k846mAuefsLBv76ktRPEbEpM8dmWs49hxl0C4bp2iVpITAcJEklniEtSfNs/eYJLtywlbu372D5ksWcc/LRjevHdM9hBiccdfCs2iVpOus3T7D28i1MbN9B0rkC3NrLt3gN6V4i4hURsTUibomIc+uuZ9Ilb39BKQjsjJa0t9pyDelGHFaKiEXA/wL+DXAXcE1EfCUzb6i3sg6DQFK/tGXWhabsORwL3JKZP8nMh4HPA6fVXJMk9V1bZl1oSjiMAFNnnbqraJOkBaUtsy404rBSVRGxBlgDMDo6Om+v24aRBZLaYfK7o+nfKU0Jhwlg6sQihxVtu8nMdcA66JwhPR+FTY4smOxAmhxZADTuzZTUDm24EtxeHVaKiCMjYkUf67gGeHpEHBER+wGvB77Sx/XvtbaMLJCkfqoUDhFxaUS8sLj/VuB64PqIOKsfRWTmb4D/BGwAbgS+mJnX92Pdc9WWkQWS1E9VDyu9DHhzcf89wEnAdmA9cFE/CsnMK4Ar+rGuflq+ZDETXYKgaSMLJKmfqh5W2i8zH46IEeDgzPxu8Zf9kwdYWyO0ZWSBJPVT1T2HayNiLfBU4O8BiqC4f1CFNUVbRhZIUj9VDYezgA8BjwDnFG0vAC4ZRFFN04aRBZLUT5XCITNvBd6wR9uXgC8NoihJUr2qjlaKiHh7RHwjIn5YtL04Il432PIkSXWo2iH9QTqHlj4NTJ6afBfwvkEUJUmqV9VweAtwamZ+Hpg8M/mnwJGDKEqSVK+q4bAIeLC4PxkOT5rSJklaQKqOVroC+FhEvBs6fRB0Ri99dVCFSdJC1YbJPKuGw3uAzwG/Avals8fwdeD3BlRXo7ThjZTUDm2ZzLPSYaXMvD8zX0PnJLjjgaMy8zWZ+cBAq2uAtlzvVVI7tGUyz6pDWR8XEY8DtgGbgHuntC1obXkjJbVDWybzrPrl/hs6Z0d3uy1obXkjJbXDQrtM6BF0hq1O3k6g0xm9ZkB1NUZb3khJ7dCWyTyr9jncvsdtI50pvBf8SXAnPnPZrNolaTqrV41wwekrGVmymABGlizmgtNXNqozGuZ2mdADgAX/Dfmtm7bNql2SZtKGyTwrhUNE/F92nfwG8ETgxcDFgyiqSexzkDSMqu453LLHzw8Bn8rMq/pcT+N4JThJw6jqlN0fGHQhTXXOyUfvdsIKNLPzSJL6qWc4RMTbqqwgMz/bv3KaxyvBSRpG0+05vKnC8xNY0OEA7eg8kqR+6hkOmXnifBYiSWqOWQ9lLWZkjcmfM/PRuRQQERcCrwYeBm4F3pqZ2+eyTklqsjZM5ll1bqWRiPhyRPwz5ak05upK4NmZ+dvAzcDaPqxTkhqpLZN5Vp0+41N0/rJ/GZ3puo8BvgK8Y64FZObXM/M3xY8bgcPmus5+W795ghM+8k2OOPfvOeEj32zcmyipPdoymWfVw0ovBEYz86GIyMy8LiLOAv6JznWl++VtwBf6uL45a8vc65LaoS0n1lbdc9hJ53ASwPaIWEbnRLhK344RcVVE/KjL7bQpy5xXvMYl06xnTUSMR8T4tm3zM31FW1JeUju0ZTLPquFwNXBKcX8Dnb/uLwfGqzw5M0/KzGd3uf0dQES8BTgVeGNm5jTrWZeZY5k5tmzZ/Ezr1JaUl9QObZmVtephpTexK0jOBv4I2B/4+FwLiIhXAO8FXpKZ/zLX9fWb02dI6qe2nFhbNRz2ycz7ADJzB/DhPtbwSeDxwJWdUbJszMw5d3T3y4nPXMbFG+/o2i5Je6MNJ9ZWDYc7IuLbwN8AX87Mh/pVQGY+rV/rGgSn7JY0jKr2OYwCX6MzdPVnEXFpRLw6IuZyPYhWsM9B0jCqeiW4+zLzf2fmi4BnA9cBfwrcM8jimqAtIwskqZ+q7jlMdQjwZGApsOCnuWjLyAJJ6qeqV4J7FnBGcVsMfBFYnZnfH2BtjdCWkQWS1E9V+wy+C1wG/D7wrblOttc2bRhZIEn9VDUcnpyZDw+0EklSY1TtkH4sGCJiy+DKkSQ1wd50SD+171VIkhplb8IhZl5EktRmexMOr+x7FZKkRunZIR0RvYLjn6Y+NmwjlyRpGEw3Wuk3QM/ps6dYNPMikqQ2mS4cjphy/1XAvwMuAG6n0yn9PjrnPkiSFpie4ZCZt0/ej4j3AGOZOTldxs0RMU7nYj9/OdgSJUnzrWqH9IHAE/doe2LRLklaYKqeIf054KqI+DhwJ3A48IdFuyRpgakaDu8FbgF+F1hOZ6ruTwKfHlBdkqQaVQqHYrjqp4qbJGmBq9TnEB1vj4hvRMQPi7YXR8TrBlueJKkOVTukPwicRecw0mjRdhed4aySpAWmaji8BTg1Mz/PrhPjfgocOYiiJEn1qhoOi4AHi/uT4fCkKW2SpAWkajhcAXwsIh4PnT4I4EPAVwdVmCSpPlXD4T3AocCv6Jz49iCdKTTO7VchEfFHEZERsbRf65Qk7Z2qQ1nvB14TEYfQCYU7M/Nn/SoiIg4HXg7c0a91SpL2XtWhrJsBMvPezLxmMhiK+ZX64X/QOdGuyiywkqQBq3pY6Wl7NhT9DnMerRQRpwETmXndXNclSeqPaQ8rRcRfF3f3m3J/0grg+iovEhFXAU/p8tB5wPvpHFKqsp41wBqA0dHRGZaWJO2tmfocbu1xP4HvAn9b5UUy86Ru7RGxks51I67r7IhwGPCDiDi2W59GZq4D1gGMjY15CEqSBmTacMjMDwBExM2Zeemej0dEt72ByjJzC3DIlPXdRue6EffNZb2SpLmp2ufwVz3ab+hXIZKk5qg6ZXeUGiIOAB7tZzGZuaKf65Mk7Z2ZOqTvpNO/sDgi9jwH4beA0qEmSVL7zbTncCadvYYrgDdNaU/g55m5dVCFSZLqM1OH9HcAImJpZv7L/JQkSapbz3CIiPMy80+LH88thpqWZOafDKIwSVJ9pttz+AAwGQ5HAQ8PvhxJUhNMFw5TDyO9OjMPGHQxTbV+8wQXbtjK3dt3sHzJYs45+WhWrxqpuyxJGpjpwuHWiPgLOlNk7BMRb6XLkNbM/OygimuC9ZsnWHv5FnY8shOAie07WHv5FgADQtKCNV04/C6dmVLPAPYDfq/LMgks6HC4cMPWx4Jh0o5HdnLhhq2Gg6QFq2c4ZObNwH8AiIhvZObL5q2qBrl7+45ZtUvSQlBp+oxhDQaA5UsWz6pdkhaCqnMrDa1zTj6axfsu2q1t8b6LOOfko2uqSJIGr+rcSkNrsl/B0UqShonhUMHqVSOGgaSh4mElSVKJ4SBJKjEcJEklhoMkqcRwkCSVGA6SpBLDQZJUYjhIkkoMB0lSSSPCISL+c0TcFBHXR8RH665HkoZd7dNnRMSJwGnAczLz1xFxSN01SdKwa8KewzuBj2TmrwEy896a65GkodeEcHgG8DsRcXVEfCcinl93QZI07OblsFJEXAU8pctD5xU1HAwcDzwf+GJEHJmZ2WU9a4A1AKOjo4MrWJKG3LyEQ2ae1OuxiHgncHkRBt+PiEeBpcC2LutZB6wDGBsbK4WHJKk/mnBYaT1wIkBEPAPYD7iv1ookacjVPloJ+Czw2Yj4EfAw8OZuh5QkSfOn9nDIzIeBM+uuQ5K0SxMOK0mSGsZwkCSVGA6SpBLDQZJUYjhIkkoMB0lSieEgSSoxHCRJJYaDJKnEcJAklRgOkqQSw0GSVGI4SJJKDAdJUonhIEkqMRwkSSWGgySpxHCQJJUYDpKkEsNBklRiOEiSSgwHSVKJ4SBJKqk9HCLiuRGxMSKujYjxiDi27pokadjVHg7AR4EPZOZzgT8pfpYk1agJ4ZDAAcX9A4G7a6xFkgTsU3cBwNnAhoj4czph9cJeC0bEGmANwOjo6PxUJ0l9dv76LVx69Z3szGRRBGccdzgfXr2y7rJ2My/hEBFXAU/p8tB5wMuAd2fmZRHxOuAi4KRu68nMdcA6gLGxsRxQuZI0MOev38LFG+947OedmY/93KSAiMx6v2Mj4lfAkszMiAjgV5l5wEzPGxsby/Hx8cEXKEl9dNTaK9jZ5Xt3UQS3XnDKwF8/IjZl5thMyzWhz+Fu4CXF/ZcCP66xFkkaqG7BMF17XZrQ5/B24BMRsQ/w/yj6FCRpIVoU0XPPoUlq33PIzH/MzOdl5nMy87jM3FR3TZI0KGccd/is2uvShD0HSRoak53OTR+tVHuH9N6yQ1qSZq9NHdKSpIYxHCRJJYaDJKnEcJAklRgOkqQSw0GSVGI4SJJKDAdJUonhIEkqMRwkSSWGgySpxHCQJJUYDpKkEsNBklTS2im7I2IbcHufV7sUuK/P66yb29QOblM7LIRtempmLptpodaGwyBExHiVec7bxG1qB7epHRbiNvXiYSVJUonhIEkqMRx2t67uAgbAbWoHt6kdFuI2dWWfgySpxD0HSVLJUIRDRHw2Iu6NiB9NaftCRFxb3G6LiGt7PPe2iNhSLDc+f1VPr8c2PTciNk7WGhHH9njumyPix8XtzfNX9fTmuE07p7yfX5m/qqfXY5ueExHfK36vvhoRB/R47isiYmtE3BIR585f1dOb4zY17vMUEYdHxLci4oaIuD4i3lW0HxwRVxafkysj4qAez2/k52nOMnPB34AXA8cAP+rx+F8Af9LjsduApXVvQ5VtAr4OvLK4fwrw7S7POxj4SfHvQcX9g+renrlsU/HYg3XXP4ttugZ4SXH/bcCHujxvEXArcCSwH3Ad8Ky6t2cu21Q81rjPE3AocExxf3/gZuBZwEeBc4v2c4E/6/Lcxn6e5nobij2HzPwH4BfdHouIAF4HXDqvRc1Rj21KYPIvtgOBu7s89WTgysz8RWb+ErgSeMXACp2FOWxTY/XYpmcA/1DcvxJ4bZenHgvckpk/ycyHgc8Dpw2s0FmYwzY1Umbek5k/KO4/ANwIjND5//5csdjngNVdnt7Yz9NcDUU4zOB3gJ9n5o97PJ7A1yNiU0Ssmce69sbZwIURcSfw58DaLsuMAHdO+fmuoq2pqmwTwBOKw04bI6Lbh7hJrmfXF/2/Bw7vskzb3qcq2wQN/zxFxApgFXA18OTMvKd46GfAk7s8pW3vU2WGA5zB9HsNL8rMY4BXAv8xIl48P2XtlXcC787Mw4F3AxfVXE8/VN2mp2bnzNU3AB+PiKPmq8C98DbgDyJiE53DGA/XXE8/VN2mxn6eIuJJwGXA2Zl5/9THsnMMaaiGdg51OETEPsDpwBd6LZOZE8W/9wJfprO731RvBi4v7v8t3WudYPe/6g4r2pqqyjZNfZ9+Anybzl9/jZSZN2XmyzPzeXT+MLm1y2Ktep8qblNjP08RsS+dYLgkMyd/334eEYcWjx8K3Nvlqa16n2ZjqMMBOAm4KTPv6vZgRPyriNh/8j7wcuBH3ZZtiLuBlxT3Xwp0O1S2AXh5RBxUjL54edHWVDNuU7Etjy/uLwVOAG6YtwpnKSIOKf59HHA+8Kkui10DPD0ijoiI/YDXA40ZhbWnKtvU1M9T0e94EXBjZn5sykNfofPHCcW/f9fl6W37PFVXd4/4fNzo/CVzD/AInWOCZxXt/wd4xx7LLgeuKO4fSWeUyHV0jqmeV/e2TLdNwIuATUW9VwPPK5YdAz4z5blvA24pbm+te1vmuk3AC4EtxTJbJt/fJtx6bNO76IyIuRn4CLtORn3sd6/4+ZRimVtb8Ls34zY19fNU/I4l8EPg2uJ2CvBbwDfo/EFyFXDwnr97xc+N/DzN9eYZ0pKkkmE/rCRJ6sJwkCSVGA6SpBLDQZJUYjhIkkoMB2kWImI0Ih6MiEU9Hv9vEXFxn14rI+Jp/ViXNFv71F2A1CaZeQfwpLrrkAbNPQcNpWLqFEk9GA4aGsWFZt4XET8EHioOEV0WEdsi4qcR8YdTlj22mOX1/oj4eUR8rGhfURzu2af4+YiI+E5EPBARVwJLp6zjX0fEXV1qOGnKa3wvIrZHxD0R8cliqoxutZ9SXIzmgYiYiIg/7v//kLSL4aBhcwbwKjoXZ/kynakcRoCXAWdHxMnFcp8APpGZBwBHAV/ssb6/oTO9x1LgQ+yai6eKnXRmml0KvKCo4Q96LHsR8PuZuT/wbOCbs3gdadYMBw2b/5mZd9L5gl2WmR/MzIezM5vrp+lMcAedeYOeFhFLM/PBzNy454oiYhR4PvBfMvPX2bkIzlerFpKZmzJzY2b+JjNvA/6KXZMM7ukR4FkRcUBm/jKLi9NIg2I4aNhMXpjlqcDy4pDO9ojYDryfXRd0OYvO1c1uiohrIuLULutaDvwyMx+a0nZ71UIi4hkR8bWI+FlE3A/8d6YcltrDa+lMBnd7cRjrBVVfR9obhoOGzeRMk3cCP83MJVNu+2fmKQCZ+ePMPAM4BPgz4EvFNNNT3QMctEf76JT7DwFPnPyhGP66bMrjfwncBDy9OHz1fiC6Fp15TWaeVtSznt6HuaS+MBw0rL4PPFB0UC+OiEUR8eyIeD5ARJwZEcsy81Fge/GcR6euIDNvB8aBD0TEfhHxIuDVUxa5mc7lS19VXEzmfODxUx7fH7gfeDAinknnqnclxbrfGBEHZuYjxXMe7bas1C+Gg4ZSZu4ETgWeC/wUuA/4DHBgscgrgOsj4kE6ndOvz8wdXVb1BuA44BfAfwX+espr/IpOB/Nn6Fwd7CE61z+Y9MfF8x+g09/R84qEwJuA24rDT+8A3jiLzZVmzes5SJJK3HOQJJUYDpKkEsNBklRiOEiSSgwHSVKJ4SBJKjEcJEklhoMkqcRwkCSV/H9CMxdBY1OnUAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "from statsmodels.api import graphics\n",
    "plt.scatter(result.fittedvalues,result.resid)\n",
    "plt.margins(0.1)\n",
    "plt.xlabel(\"residuals\",fontsize=12)\n",
    "plt.ylabel(\"fitted-values\",fontsize=12)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAEKCAYAAAASByJ7AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi40LCBodHRwOi8vbWF0cGxvdGxpYi5vcmcv7US4rQAAIABJREFUeJzt3Xl8VPX1//HXAUWkigtutZJEkLpQQQUVl7ovdcNat69GrXXBpeKOG/251biBWhGrDYoVjWtxQYt1rdbWIgZFAbeCJghiDVoWRRDI+f3xuZEQkpmbZGbuTPJ+Ph7zmJl7Z+49GXHO3M9yPubuiIiIdEg6ABERyQ9KCCIiAighiIhIRAlBREQAJQQREYkoIYiICKCEICIiESUEEREBlBBERCSyWtIBNMcGG2zgJSUlSYchIlJQJk2aNNfdN0z3uoJKCCUlJVRWViYdhohIQTGz6jivU5ORiIgASggiIhJRQhAREUAJQUREIkoIIiICKCGIiEhECUFERAAlBBGRwldZCZdcAq1cElkJQUSkELnDa6/BgQfCjjvCqFFQHWv+WZOUEERECok7PPcc/PznsNdeMHky3HRTSAatLO1TUKUrRETardpaePJJuP56ePtt6N4d7rgDTj0V1lwzI6dQQhARyWdLl8LDD8MNN8CHH0KvXjB6NJSWQqdOGT2VEoKISD5avBj+/OfQHFRVBX36wCOPwFFHQceOWTmlEoKISD755hsoL4fhw2HOHNh5ZxgxAg49FMyyemolBBGRfDBvXugTuP12+Oor2GcfePBB2HvvrCeCOkoIIiJJ+vJLuO02uPNOWLgwXAkMHQoDBuQ8FCUEEZEkfPYZDBsW5g8sWQJHHw1XXAF9+yYWkhKCiEguTZ8ON94IY8aEOQUnngiXXQY//WnSkSkhiIjkxJQpYejoo4/C6qvDoEEwZAgUFycd2Q+UEEREsmniRCgrg3HjYK214OKL4YILYJNNko5sFUoIIiKZVldnqKwMXnoJ1lsPrr4aBg+G9ddPOromKSGIiGSKO4wfH8pLvPFGuAoYNgzOOAPWXjvp6NJSQhARaa3ly+GJJ0IimDwZiorCMNJTToHOnZOOLrbEq52aWUcze8fMnk06FhGRZlm6FO6/H3r3hmOOgUWL4L77wkiis88uqGQA+XGFcB7wAdA16UBERGJZvDgUmLv55lB2um9feOwx+NWvslZnKBcSvUIws82AQ4B7koxDRCSWb74JNYY23xx++1vYdFN49ll4550wsayAkwEkf4XwB+ASIP97W0Sk/fr6axg5MtQZ+vpr2HdfeOihsEBNjuoM5UJiVwhmdijwpbtPSvO6QWZWaWaVNTU1OYpORAT473/h0kvD5LGrroLdd4cJE8JQ0hwWncuVJJuMdgMGmlkV8Aiwj5k92PBF7l7u7v3dvf+GG26Y6xhFpD2aOTPMGSgpCU1Ehx4K774LTz8dylG3UYklBHe/3N03c/cS4P+AV9z9hKTiERHh44/DUNGePeHuu+H448MqZQ8/HBaoaeOS7kMQEUnee++FOQSPPx6WpTzrrFBioqgo6chyKi8Sgru/CryacBgi0t68+WYoL/HMM6HO0JAhoc7QxhsnHVki8iIhiIjkjDv8/e/hiuDll0NtoWuvhXPOCTWH2jElBBFpH9zhr38NVwQTJoQ6Q8OHhzpDa62VdHR5IfHSFSIiWbV8eZhFvP32cNhhYeH6P/4RPv0ULrqo4JNBRUUYDNWhQ7ivqGj5sXSFICJt09KlYZH6G28Mo4e22irUHTruuLBATRtQURHW2Vm0KDyvrg7PAUpLm388XSGISNvy3Xeh0ugWW4QhpF26hNFDU6fCSSe1mWQAMHToimRQZ9GisL0llBBEpG1YuDAUm9t889BBvNlmoc/g7bfhqKNyVmcok0046cyc2bzt6ajJSEQK29dfw4gR4fa//8H++4efyHvskfPSEpluwkmnqCico7HtLaErBBEpTF98AZdcEuoMXXNNSAATJ8ILL8CeeyZSZyjTTTjplJWFFrH6unQJ21tCCUFECkt1dSg9XVICt9wCAweGmcZPPQU77phoaJluwkmntBTKy0NONAv35eUtvxpRQhCRwvDRR3DyyaGzeNQoOPHEsK2iArbdtlVt95lq92+qqSabFTBKS6GqCmprw31rmqaUEEQkv02eHJan3HrrMJ/g7LNhxoyQFLbYAljRdl9dHeaf1bXdx/lib817G8p0E06umbsnHUNs/fv398rKyqTDEJFc+Pe/wzfpX/8Ka68dmokuuAA22miVl5aUNN65WlwcfjWn0pr3NqaiIvQZzJwZrgzKyrLTodwcZjbJ3funfZ0SgojkDfdQX6isDF59Fbp1g/PPD8NI1123ybd16BDe2pBZaEpJpTXvLRRxE4KajEQkebW1YfGZAQPCsNGPPgodxlVV8LvfpUwG0Lq2+yTa/fOVEoKIJGf58rD4TN++8MtfQk1NWJjmk0/gwgtj1xlqTdt9obf7Z5ISgojkTN1onjXsey7pdi8LfrJVWJVs+XIYMybUHDrjDCrGdm7WqJ/WDL/M9NDNQqY+BBHJiYoKOO/0RRz/3T0MYRjdmcVk254Fg4eyx21HhG9/Vp3tC+EXe3v9ks4E9SGISP5YsIBZ59zI+9+VMILzqKaYgxjP9j6Jk54+8odkALmf7SsrKCGISPZ89RVceSUUF3PpvMt5h+3Zk1f5Oa/zNw4CbJVZvLme7SsrqLidiGTe55/DrbeGDuJvv4UjjmDgvy/nmS9WLS3RcDRPpgu2SXyJXSGYWXcz+7uZvW9m08zsvKRiEZEMqaqCs84KJahvuy2MHJo6FZ54gmOH7xhrNI9G/SQnySajZcBF7r4NMAD4rZltk2A8ItJSH34Iv/51KCUxenSoOfTxx2HFst69gfijeTTqJzmJJQR3n+Pub0ePFwIfAD9JKh4RaVzKwm/vvANHH41vsw3fPfA4f1g+mJ03/ISKPf4EPXuucqy4hdgyWbBN4suLPgQzKwG2B95MNhIRqa+pBV82/PhfHPBWGTz3HN+v2ZXbOl7B8GXnMZcNYXZ2F4WR7El8lJGZrQWMBc539wWN7B9kZpVmVllTU5P7AEXasHRln1ceAursx4v8ddFeHHDt7vDWW1BWRr8NZnLZsutCMohomGhhSnRimpmtDjwLPO/ut6Z7vSamiWROnAlgHToAXsthPMNQytiJt5jFT7iFi7ntm9PhRz9qF8XhCl3eT0wzMwPuBT6IkwxEJLPSTgBbtoxzuj3Eu/TlaX5JN77idMrpyQyeLD4ffvQjQMXh2pIkm4x2A04E9jGzydHt4ATjEWlXmproNaf6e7jnHthqK0bMLaWj1VLKg2zJR9zD6azWZY2VhoBqmGjbkeQoo3+6u7l7H3ffLrqNTyoekfam4S/4NVnEYEbwaceecPrpoeT0E0/wzv1T+FdxKbW2WqNDQDVMtO1IvFNZRLInVadx3S/7rsznMm6gilBnqGOvHvD886HT+IgjKD2xQ9ohoBom2jbkxbBTEcm8poaMQvjCLj1wLn32v53icXfQ1efzaudfMHnI0DCCSNolXSGIFLimrgKa6jS+/dLPw+IzxcVsO66Mrr/aDyor2eu755QM2rm0Vwhm1hOY5e5LzGwvoA8wxt3nZTs4EUkt1VVAw07jEj7lUm7iN7PvgxHLw8I0l10G26hijARxrhDGAsvNbAugHOgOPJTVqEQkllRDR+s6jbfmfe7nJP5DL37DfTy+1imhztCYMUoGspI4CaHW3ZcBRwB3uPsQ4MfZDUtE4ki1dsBdp03iyY5HMpWfcSRjuZ3z2Kbzp9jdd0GPHrkNVApCnISw1MyOA35NmFUMsHr2QhJpf1KNBkq1r7HJX7vzOq+scRAH/b/+HLTGK4xcZyglVHNH8S1ce8+mGgEkTYozyug3wJlAmbt/amabAw9kNyyR9iNVPwCkHilUVla33zmAFxhKGXvwOos7bQhX3cAaZ5/NuV27cm5u/yQpULFqGZnZmkCRu3+U/ZCaplpG0haVlDS+Qlhxcbhval9VFVBby2sXPs16d5XR5/tJfN5xMz4/fgj97z5t1enD0m5lrJaRmR0GTAb+Fj3fzszGtT5EEYHU/QBN7ZtdvSxcWvTpw563/4o+3efBqFFsumgG/cecq2QgLRKnD+FqYCdgHoC7TwbUIyXtVkvb+5val6o4XMN9nVjC6ZQzfbUt4YQTQq2Ihx4KK5addhp06pSpP1PaoTh9CEvdfX4oTvoDFbWVdqml7f2p9q3oB1jx2vrF4QYNAhZ9y+mM4mKGsxmzmVu0Y1jE/rDDohrVIq0XJyFMM7PjgY5m1gs4F3gju2GJ5Kd0JaNbsq+qasWxZ84MVwVlZVE9oHnz6H3InXQf+we61c5lwhp78tH597HvDfuFqwORDErbqWxmXYChwAGAAc8Dv3f3xdkPb2XqVJakpVoMBlq2r9FFZGpq4A9/gJEjYcECOPhguOIK2G23Fscu7VfcTuW0VwjuvoiQELQgnrR7RUWNj/qpa+tv6b4fzJoFw4eH+tGLF8ORR4ZEsP32rYpbJI4mGx/N7BkzG9fULZdBirREuvWCW9IBnGoxmJbuA2DGjNBZ0KNHuCo45hiYNg0ef1zJQHLH3Ru9AXumujX1vmze+vXr5yJxPPige5cu7qGhJty6dAnb0+2P897iYnezcF+3vUX7pk51P/549w4d3NdYw/3MM90//TT7H5C0K0Clx/iOjTUxLV+oD0HiSjXZq6qqlZPBMqGyMlwePPVUWJv4zDPhoovgxyoTJpnX6j4EM3vM3Y8xsynAKlnD3fu0MkaRrEk12SvO/ubui+0f/wiJ4IUXwhKVV14J554L3bpl4OAirZNqAPN50f2hwGGN3EQS15LJXun2p3tvs7nDc8/Bz38Oe+4JkyfDjTeGy5BrrlEykLzRZEJw9znRw7Pdvbr+DTg7N+GJNK1uklh1dfjOrZvsVVGRvhO3VR3AcdXWwtix0L9/GDZaXQ133AGffgqXXgpduzbzgCJZlq6TAXi7kW3vxemgiHHsXwAfAdOBy9K9Xp3KUl9x8codv3W34uKwP1UHb7r96d6b0tKl7mPGuG+9dQioVy/3e+91X7KkxX+rSGvQ2k5lMzuLcCXQA5hRb9fawL/c/YTWJCIz6wh8DOwPzALeAo5z9/ebeo86laW+VJPEGp3slW2LF8P998NNN4WrgG23DXMIjj4aOnZMICCRIBMT0x4CngNuAC6rt32hu3/dyvggFMyb7u6fAJjZI8DhQJMJQaS+dJPEcubbb+FPfwoTyubMgZ13DrOMDztM5SWkoKTqQ5jv7lXufhzhF/xSwmijtcwsE//L/QT4rN7zWdG2lZjZIDOrNLPKmpqaDJxWCkmqyWMZa+tvqXnz4LrrwnjUiy6CrbaCl16Cf/8bBg5UMpCCk7Z0hZmdQyiB/V9WVDl1ICfDTt29HCiH0GSUi3NKfkhVWbS0lB+Wgmy0KFw2ffkl3HYb3HknLFwIhx4amoZ22SXLJxbJrjjF7aYDO7v7Vxk9sdkuwNXufmD0/HIAd7+hqfeoD6F9STe5LOdmzYJhw2DUqNBfcPTRIRH07ZtAMCLxZay4HaFZZ37rQ1rFW0CvaI3m2cD/Acdn4TxSoFoyeSwrpk8P8wbGjAm92CecAJddBltumeNARLIrTkL4BHjVzP4KLKnb6O63tubE7r4sao56HugIjHb3aa05prQtiXcaT50K118Pjz4Kq68e2quGDFlR30KkjYmTEGZGt07RLWPcfTwwPpPHlLYj3UpiWTNxYkgETz8Na60VOowvvBA22STLJxZJVpz1EK7JRSAiDeW009gdXnstnOCll2C99eDqq2HwYFh//SycUCT/xBlltCFwCdAb6Fy33d33yWJcIsDKo4myoq7OUFkZvPEGbLwx3HxzqD669tpZPLFI/omzOncF8CGwOXANUEXoEBYpXMuXh8VndtgBDjkkjCAaOTLMMB4yRMlA2qU4CaGbu98LLHX319z9FEBXB1KYli4N5SV69w6rkn37LYweHUYS/fa3sOaaSUcokpg4ncpLo/s5ZnYI8DmgRlUpLIsXw333heagqqowd+DRR8OaxaozJALESwjXmdk6wEXAHUBX4IKsRiWSKd98A3ffDbfcAl98AQMGhBLUhxyi0hIiDcQZZfRs9HA+sHd2wxHJkP/9L3zx3347fP017LsvPPQQ7LWXEoFIE9L2IZjZfWY2uuEtF8FJ4UtVnC7O/mb773/D4jNFRXDVVbD77jBhQhhKuvfeSgYiKcRpMnq23uPOwBGEfgSRlNIVp0u3v1lmzgx1hu65B77/PnQYX3459NHS3yJxpS1ut8obzDoA/3T3XbMTUtNU3K6wpCtOl5HidR9/HOoMPfBAeH7SSaHOUK9eLYpZpC3KZHG7hnoBG7XgfdLOpCtO16ride+9F8pLPP44dOoEZ50FF1+cwOo4Im1HnD6EhWa2oO4eeAa4NPuhSb5L1/7f1Hdz3fZ0+xs1YUJYiaxvXxg/Pkwiq6qCESOUDERaKW1CcPe13b1rvfufuvvYXAQn+auu/b+6OlR/qGv/b86KZrFXPHOHV14JI4V22SWUmLj22nDSG28M5SZEpPXcvckbsCZwGnBrdDse6JTqPdm89evXzyU/FBe7h2/qlW/FxSu/7sEHwzazcP/gg83YX1vr/swz7gMGhINvson78OHuCxdm8S8TaXuASo/xHdtkp7KZbQuMA/4BTIo29wO2AfYHLnb332UzWTWkTuX80aFDSAENmUFt7arbm2X5cvjLX0IfwXvvhfaoSy6B3/wGOndO+3YRWVkmOpVHAIPc/cUGB94PmApoMZt2LCuL13z/PTz4YGgG+s9/wqL1998Pxx0XFqgRkaxK1Yfw44bJAMDdXyLUNzoia1FJ3ovd/h/Hd9+FSqNbbAGnnhoWpXn88bBi2UknKRmI5EiqhNDBzNZouNHMOhMqny5q5D3STpSWQnl5mDNgFu7Ly5s5oWzBArjpptAkNHhwuLwYPx4mTYKjjlLROZEcS9VkNAYYa2a/dfdqADMrITQlPZD90CTftXjxmq++CsNER4yAefNg//3hiitgzz1VWkIkQU0mBHe/zszOAV43s7rGgW+B4e5+R06ik7Zlzhy49Va4666wDsHhh4dEsNNOSUcmIqSZh+DuI929iLBa2ubuXpyJZGBmw8zsQzN7z8yeNLN1W3tMyZyMF5yrrg6Lz2y+eUgIAweG0UNPPaVkIJJH4qyYhrsvdPeFGTzvi8DP3L0P8DFweQaPLa0QZ8JZbB99BCefHDqLR42CE08M2x56CLbdNtOhi0grxUoImebuL7j7sujpBGCzJOJoT+L+6h86dEX10TqLFoXtsU2eHKqNbr01PPYYnH02zJgRksIWW7TwLxCRbGtJcbtMOwV4NOkg2rLmlJluVcG5N94I407Hj4euXUPV0fPPh41UC1GkEMQpbtfFzP6fmY2Knvcys0NjvO8lM5vayO3weq8ZCiwDmmyQMLNBZlZpZpU1NTXx/ipZSXN+9Te74Jz7isVndtsN3nwTfv/7kHWuv17JQKSAxLlCuI9QumKX6Pls4HFWXjhnFe6+X6r9ZnYycCiwrzdVPyMcpxwoh1C6Ika80kBzfvWXla18NQFNTDirrYVnnw07Jk6ETTcNHcaDBsGPfpSx2EUkd+L0IfR095sJs5OJJqS1arC4mf0CuAQYqAluLdOckUDN+dWfdsLZ8uXwyCOw3XZh2GhNTVjE/pNP4IILlAxEClichPC9ma0JOICZ9QSWtPK8I4G1gRfNbLKZ3d3K47UrzR0J1NwyE6WlYYmB2tpwX1pKqDN0772hvtBxx8GyZWGVso8/hjPOgDVWmdQuIgUmTpPRVcDfgO5mVgHsBpzcmpO6u4aatEKqPoHGZg7XbRs6NDQTFRWFZBBrlvF334V1iocNg88+gx12gLFj4Ze/DJcnItJmxFpT2cy6AQMITUUT3H1utgNrjMpfB1ktPV1nwYIwo/jWW+HLL2H33UNGOfBAlZcQKTCtLn9tZjs02DQnui8ysyJ3f7s1AUrLZaX0dJ2GdYYOPDCUl9hjjwwcXETyWaomo1tS7HNgnwzHIjHFHgnUHHPmwC23hA7ib7+FI44IiaB/2h8VItJGpCput3cuA5H4WtUn0FBVFdx8M4weDUuXhg7jyy+H3r0zGbKIFIC0ncrR+gdnA7sTrgxeB+5298VZjk1SaHHp6Toffgg33BCGJnXoEGoOXXop9OyZqRBFpMDEGWU0BlgI1FU5PZ6wHsLR2QpKsuidd8IM4rFjw/rEgwfDRRfBZionJdLexUkIP3P3beo9/7uZvZ+tgCRL/vWvkAjq6gxdfnmoM7ThhklHJiJ5Is5A8rfNbEDdEzPbGdDYz0LgDi++CHvtFYaNTpwI110XhiiVlSkZiMhK4lwh9APeMLO6yjdFwEdmNgXwaE0DySe1tTBuXLgieOutUGfottvg9NNVWkJEmhQnIfwi61FIZixbFtYfuP56mDYNevQIhYhOOkmlJUQkrbQJwd2rzWw9oHv912tiWh5ZsgTGjIGbbgoL0WyzDTz4IBx7LKyWD0teiEghiDPs9PeE2kUziArcoYlp+WHRorAK2bBhMHs29OsHTzwRqpCqzpCINFOcn4/HEEpgf5/tYCSm+fPhj38M/QI1NaGsxOjRsP/+qjMkIi0WJyFMBdYFvsxyLJLO3Lnwhz/AyJEhKfziF2G68u67Jx2ZiLQBcRLCDcA7ZjaVeusguPvArEUlK5s9O9QZ+tOfQjnqX/0q1BnaoWH9QRGRlouTEO4HbgKmAJkqrixxfPJJ6Cj+85/DSmXHHx8mlG29ddKRiUgbFCchLHL3EVmPRFZ4//1QZ+jhh6FjRzjlFBgyJAwjFRHJkjgJ4XUzuwEYx8pNRhp2mmmTJoU5BE88EepZn3deqDO06aZJRyYi7UCchLB9dD+g3jYNO82k118PpSSefx7WWQd+97uQDDbYIOnIRKQdiTMxTesiZIM7vPBCSASvvx7qCt1wA5x9dig+JyKSY7GmsZrZIUBvoHPdNne/NltBtWm1tfDUU6FpaNKkUHb69tvhtNNCM5GISELSTmc1s7uBY4HBgBHWQSjOxMnN7CIzczNr+20jy5aFchLbbgtHHhnmEdxzTyg1ce65SgYikrg49Q12dfeTgP+5+zXALsBPW3tiM+sOHADMTPfagrZkSSgwt+WWcOKJoaTEQw/BBx/AqadCp05JRygiAsRLCN9F94vMbFNgKfDjDJz7NuASVtRHalu+/TaUlujRA844A7p1C01F774b1i1W0TkRyTNxvpWeNbN1gWHA24Qv8FGtOamZHQ7Mdvd3LU3tHTMbBAwCKCoqas1pc2PePLjzzlBiYu7csDjN/ffDvvuqzpCI5DVzj/8D3czWADq7+/wYr30J2KSRXUOBK4AD3H2+mVUB/d19brpj9u/f3ysr83SxtpqacEVw552wYAEcfHCoM7TrrklHJiLtnJlNcvf+6V7X5BWCme0IfObuX0TPTwKOBKrN7Gp3/zrVgd19vyaOuy2wOVB3dbAZYZnOnerOVVBmzYLhw0M/weLFcNRRobzE9tunf6+ISB5J1YfwJ+B7ADPbA7gRGAPMB8pbekJ3n+LuG7l7ibuXALOAHQouGcyYAYMGhT6CkSPhmGNCyYnHHlMyEJGClKoPoWO9q4BjgXJ3HwuMNbPJ2Q8tT02btqLO0Oqrh3WKhwyBkpKkIxMRaZWUCcHMVnP3ZcC+RB27Md7XLNFVQv6rrAyzip96KixUf+GF4fbjTAy4EhFJXqov9oeB18xsLmHo6esAZrYFodmoffjHP0IieOEFWHdduPLKMJGsW7ekIxMRyagmE4K7l5nZy4Q5By/4iuFIHQizltsud/jb30J5iX/+EzbaKKxLcOaZqjMkIm1WyqYfd5/QyLaPsxdOwmpr4cknQyJ4+23o3h3uuCPMKF5zzaSjExHJqjgzldu+pUthzBjo3TsMG124EO69F6ZPh3POyVoyqKgIfdEdOoT7ioqsnEZEJJb2XT9h8eKwPOVNN0FVFfTpA488EpJCx45ZPXVFRRi1umhReF5dHZ4DlJZm9dQiIo1qn1cI33wDt94a5hCcdRZsvDGMGweTJ8Oxx2Y9GUCYxFyXDOosWhS2i4gkoX1dIcybF/oEbr8dvvoK9tkHHngg3Oe4ztDMJmq8NrVdRCTb2kdC+PLLFXWGFi6EQw8NP8UHDEj/3iwpKgrNRI1tFxFJQttuMvrsszBnoLg49BMcdFBoFnrmmUSTAYSpDQ3XxOnSJWwXEUlC20wI06eHJSl79oS77grrD3z4ITz6KPTtm3R0QOg4Li8Pucos3JeXq0NZRJLTtpqMpkwJdYYefTTUGRo0KNQZKs7Iip8ZV1qqBCAi+aNtJISJE0Nby7hxsNZacPHFcMEFsEljyzGIiEhjCj8hlJeHJSrXWw+uvhoGD4b11086KhGRglP4CeHww2H+/FBnaO21k45GRKRgFX5C2Hjj0E8gIiKt0jZHGYmISLMpIYiICKCEICIiESUEEREBlBBERCSSWEIws8Fm9qGZTTOzm5OKQ0REgkSGnZrZ3sDhQF93X2JmGyURh4iIrJDUFcJZwI3uvgTA3b9MKA4REYkklRB+CvzczN40s9fMbMeE4hARkUjWmozM7CWgsepyQ6Pzrg8MAHYEHjOzHu7ujRxnEDAIoEirx4iIZE3WEoK779fUPjM7C3giSgATzawW2ACoaeQ45UA5QP/+/VdJGCIikhlJNRk9BewNYGY/BToBcxOKRURESK643WhgtJlNBb4Hft1Yc5GIiOROIgnB3b8HTkji3CIi0jjNVBYREUAJQUREIkoIIiICKCGIiEhECUFERAAlBBERiSghiIgIoIQgIiIRJQQREQGUEEREJKKEICIigBKCiIhElBBERARQQhARkYgSgoiIAEoIIiISUUIQERFACUFERCJKCCIiAighiIhIJJGEYGbbmdkEM5tsZpVmtlMScYiIyApJXSHcDFzj7tsBV0bPRUQkQUklBAe6Ro/XAT5PKA4REYmsltB5zweeN7PhhKS0a0JxiIhIJGtXCGb2kplO3nveAAAJNUlEQVRNbeR2OHAWcIG7dwcuAO5NcZxBUT9DZU1NzUr7KiqgpAQ6dAj3FRXZ+muSPaeISC6Yu+f+pGbzgXXd3c3MgPnu3jXd+/r37++VlZVA+CIeNAgWLVqxv0sXKC+H0tLsxJ3EOUVEWsvMJrl7/3SvS6oP4XNgz+jxPsB/mnuAoUNX/mKG8Hzo0FbHllfnFBHJlaT6EE4Hbjez1YDFwKDmHmDmzOZtz4QkzikikiuJJAR3/yfQrzXHKCqC6urGt2dLEucUEcmVgp2pXFYW2u/r69IlbG9L5xQRyZWCTQilpaEzt7gYzMJ9tjt3kziniEiuJDLKqKXqjzISEZF48n2UkYiI5BklBBERAZQQREQkooQgIiKAEoKIiEQKapSRmdUAjUwNS9QGwNykg4ipkGKFwoq3kGKFwoq3kGKF/Iy32N03TPeigkoI+cjMKuMM58oHhRQrFFa8hRQrFFa8hRQrFF689anJSEREACUEERGJKCG0XnnSATRDIcUKhRVvIcUKhRVvIcUKhRfvD9SHICIigK4QREQkooTQTGZ2tJlNM7NaM2tyJIGZVZnZFDObbGaJVORrRqy/MLOPzGy6mV2WyxgbxLG+mb1oZv+J7tdr4nXLo891spmNy3GMKT8rM1vDzB6N9r9pZiW5jK+ReNLFe7KZ1dT7PE9LIs4oltFm9qWZTW1iv5nZiOhvec/Mdsh1jPViSRfrXmY2v97nemWuY2wRd9etGTdga2BL4FWgf4rXVQEb5HusQEdgBtAD6AS8C2yTULw3A5dFjy8Dbmridd8kFF/azwo4G7g7evx/wKMJ/vePE+/JwMikYmwQyx7ADsDUJvYfDDwHGDAAeDOPY90LeDbpz7S5N10hNJO7f+DuHyUdRxwxY90JmO7un7j798AjwOHZj65RhwP3R4/vB36ZUBxNifNZ1f8b/gLsa2aWwxjry6f/tmm5+z+Ar1O85HBgjAcTgHXN7Me5iW5lMWItSEoI2ePAC2Y2ycyavWZ0Dv0E+Kze81nRtiRs7O5zosdfABs38brOZlZpZhPMLJdJI85n9cNr3H0ZMB/olpPoVhX3v+2RURPMX8yse25Ca5F8+rcaxy5m9q6ZPWdmvZMOJo5E1lTOd2b2ErBJI7uGuvvTMQ+zu7vPNrONgBfN7MPoV0VGZSjWnEkVb/0n7u5m1tQQuOLos+0BvGJmU9x9RqZjbSeeAR529yVmdgbh6mafhGNqC94m/Dv9xswOBp4CeiUcU1pKCI1w9/0ycIzZ0f2XZvYk4fI94wkhA7HOBur/Ktws2pYVqeI1s/+a2Y/dfU7UFPBlE8eo+2w/MbNXge0JbeXZFuezqnvNLDNbDVgH+CoHsTUmbbzuXj+2ewj9OPkqp/9WW8PdF9R7PN7M/mhmG7h7vtU4WomajLLAzH5kZmvXPQYOABodjZAH3gJ6mdnmZtaJ0BGa05E79YwDfh09/jWwyhWOma1nZmtEjzcAdgPez1F8cT6r+n/DUcArHvUyJiBtvA3a4AcCH+QwvuYaB5wUjTYaAMyv18SYV8xsk7q+IzPbifBdm9QPg/iS7tUutBtwBKHtcgnwX+D5aPumwPjocQ/CiI53gWmE5pu8jDV6fjDwMeFXdiKxRnF0A14G/gO8BKwfbe8P3BM93hWYEn22U4BTcxzjKp8VcC0wMHrcGXgcmA5MBHok/O81Xbw3RP9G3wX+DmyVYKwPA3OApdG/21OBM4Ezo/0G3Bn9LVNIMcovD2I9p97nOgHYNcl/B3FvmqksIiKAmoxERCSihCAiIoASgoiIRJQQREQEUEIQEZGIEoLknJl1q1cF8gszmx09nmdmuZpTUBfLdtFM0rrnA1ta8TWqcLtBI9vXMbMxUZXOGWZW0VQl19ZI9beY2dVmdnGmzyltixKC5Jy7f+Xu27n7dsDdwG3R4+2A2kyfL5ox3JTtCGP162Ib5+43ZjiEe4FP3H0Ld+9JmKPw5wyfA3Lzt0gbpoQg+aajmY2ysI7DC2a2JoCZ9TSzv0XFAl83s62i7SVm9kpUnO1lMyuKtv/ZzO42szeBm6PZ46PNbKKZvWNmh0ezd68Fjo2uUI61sD7AyOgYG5vZk1GBsnfNbNdo+1NRHNPSFS40sy2AfsDv622+FuhrZltGdfOfrff6kWZ2cvT4SjN7y8ymmll5vZmvr5rZTdHf8rGZ/Tzd39IgpqY+y6Ojc71rZhkvsyL5TwlB8k0v4E537w3MA46MtpcDg929H3Ax8Mdo+x3A/e7eB6gARtQ71maEGaIXEornveLuOwF7A8OA1YErCWsWbOfujzaIZQTwmrv3JdS+nxZtPyWKoz9wrpmlqma6DTDZ3ZfXbYgev0NYryKVke6+o7v/DFgTOLTevtWiv+V84CoP5a1T/S31NfVZXgkcGP29A9PEJm2QittJvvnU3SdHjycBJWa2FqFkxeO2YmmBNaL7XYBfRY8fYOXibI/X+yI+ABhYrx29M1CUJpZ9gJPghy/x+dH2c83siOhxd0ISy0admr3N7BKgC7A+ISE9E+17IrqfBJTEPWCaz/JfwJ/N7LF6x5d2RAlB8s2Seo+XE34ZdwDmRf0MzfFtvccGHOkNFgwys52bc0Az2wvYD9jF3RdZqLbaOcVb3ge2M7MO7l4bHaMD0JdQIrmIla/UO0ev6Uz45d7f3T8zs6sbnKfuc1pO8/4/bvKzdPczo8/jEGCSmfXzlauhShunJiPJex5KCX9qZkfDD2vr9o12v0Go4glQCrzexGGeBwbXa4ffPtq+EFi7ife8DJwVvb6jma1DKGf9vygZbEVYyjFV7NMJzUO/q7f5d8DL7j4TqAa2sbAW87rAvtFr6r7850a/6o9KdZ4Yf0tdPE1+lmbW093fdPcrgRpWLjUt7YASghSKUuBUM6urIFu3FORg4Ddm9h5wInBeE+//PaHP4D0zm8aKTt6/E76QJ5vZsQ3ecx6h2WYKoWlmG+BvwGpm9gFwI6GSZTqnEMpQzzCzGkISORPA3T8DHiOUR3+MkDxw93nAqGj784RS1umk+lvqa+qzHGZmUywsHP8GoVKntCOqdiqSQ2a2JfBX4Fx3H590PCL1KSGIiAigJiMREYkoIYiICKCEICIiESUEEREBlBBERCSihCAiIoASgoiIRP4/ozfYVKRyhRUAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "from statsmodels.api import graphics\n",
    "fig = graphics.qqplot(result.resid,line=\"s\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Figures for publication"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAENCAYAAADgwHn9AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi40LCBodHRwOi8vbWF0cGxvdGxpYi5vcmcv7US4rQAAIABJREFUeJzt3XuYXVV9//H3h5BLITEIsUNCKAm0seQumYDBUgYpiKBQLIooLVfDrVSUPCI3raSAlYj8UPuTFEjA8hMoF6Fgq0FzIDwEkEkD5gJBcCITMoEkgJmETG7f3x97n+NhmMvJZM51Pq/nOc/Ze519+R42Od9Za+29liICMzMzgN3KHYCZmVUOJwUzM8txUjAzsxwnBTMzy3FSMDOzHCcFMzPLcVIwM7McJwUzM8txUjAzs5zdyx3Azho2bFiMGjWq3GFYD23cuJE999yz3GFYD/jaVbfGxsa1EfGh7raruqQwatQonnvuuXKHYT2UyWRoaGgodxjWA7521U3SykK2c/ORmZnlOCmYmVmOk4KZmeVUXZ+CmVl3tm7dSnNzM5s3by53KCU3aNAgRo4cSf/+/Xu0v5OCmdWc5uZmhgwZwqhRo5BU7nBKJiJYt24dzc3NjB49ukfHcPORmdWczZs3s88++/SphAAgiX322WeXakhOCmZWk/paQsja1e/tpGBmZjlOCmZmRdLS0sLnP/95DjroIKZMmcLxxx/PihUryh1Wl5wUzMyKICI4+eSTaWho4JVXXqGxsZHrr7+eNWvW5LbZtm1bGSPsWEmSgqT9Jc2XtEzSUklfbvf5pZJC0rBSxGNmVmzz58+nf//+nH/++bmySZMmsX37do444ghOPPFExo4dC8CNN97I+PHjGT9+PDfddBOQjDV1wgknMGnSJMaPH88999wDwNe//nXGjh3LxIkTmTFjRq/HXapbUrcBl0bEIklDgEZJ8yJimaT9gWOB35coFjOz91u4EDIZaGiAadN2+XBLlixhypQpHX62aNEilixZwujRo2lsbGTOnDk888wzRASHHXYYRx55JK+++iojRozg0UcfBeCdd95h3bp1PPjgg7z44otI4u23397lONsrSU0hIlZHxKJ0eQOwHNgv/fh7wNeAKEUsZmbvs3AhHH00XH118r5wYVFPd+ihh+aeI3jyySc5+eST2XPPPRk8eDCf+cxnWLBgARMmTGDevHlcdtllLFiwgKFDhzJ06FAGDRrEOeecwwMPPMAee+zR67GVvE9B0ijgI8Azkk4CVkXE86WOw8wsJ5OBLVtg+/bkPZPZ5UOOGzeOxsbGDj8rZAjyMWPGsGjRIiZMmMBVV13FNddcw+67786zzz7LKaecwiOPPMJxxx23y3G2V9InmiUNBu4HLiFpUrqCpOmou/2mA9MB6urqyPTCBbPyaG1t9fWrUtV07YYOHcqGDRsK3n63qVPZY8CAJCEMGMCmqVPZsRP7d2Tq1Kls2rSJm2++mbPOOgtImpQee+wxtm3blovvkEMO4YILLuCiiy4iIrj//vuZPXs2K1as4IMf/CAnnXQSAwYM4M4772T16tW8++67HHHEEUycOJGJEyd2+D03b97c82sVESV5Af2BnwNfTdcnAG8ATelrG0m/wr5dHWfKlClh1Wv+/PnlDsF6qJqu3bJly3Z+p6eeirjuuuS9l6xatSo++9nPxoEHHhhjx46N448/PmbPnh0nnHDCe7b77ne/G+PGjYtx48bF9773vYiI+J//+Z+YMGFCTJo0Kerr6+PXv/51vP766zF16tSYMGFCjB8/PubOndvheTv6/sBzUcBvtZJti0vJI3Z3AOsj4pJOtmkC6iNibVfHqq+vD0+yU708UUv1qqZrt3z5cg4++OByh1E2HX1/SY0RUd/dvqXqU/gY8PfAxyUtTl/Hl+jcZmZWoJL0KUTEk0CXA3JExKhSxGJmZp3zE81mZpbjpGCls3AhtLQU/R5wM+s5JwUrjezDQatWleThIDPrGScFK43sw0HQaw8HmVnvc1Kw0mhogAEDkuUBA5J1sxp37bXXMm7cOCZOnMjkyZN55plnOt32zDPP5L777gNgwYIFjBs3jsmTJ/Puu++WKlzAczRbqUybBjfdBFLy3gsDjplVsoULF/LII4+waNEiBg4cyNq1a9mSrS1346677uLyyy/n9NNPL3KU7+eagpXGwoVwySXwhz8k7+5TsAqzuq2Ng55+mpa2tt453urVDBs2jIEDBwIwbNgwRowYQWNjI0ceeSRTpkzhE5/4BKtXr37Pfrfeeiv33nsvV199NV/84hd7JZad4aRgpeE+BatwM5uaaNq8mZkrV/bK8Y499lhee+01xowZw4UXXsjjjz/O1q1bufjii7nvvvtobGzk7LPP5sorr3zPfueeey4nnngiN9xwA3fddVevxLIz3HxkpdHQAP36Jcv9+rlPwSrK6rY25qxZww5gTksLVx9wAPumf+H31ODBg2lsbGTBggXMnz+fU089lauuuoolS5ZwzDHHALB9+3aGDx/eC9+g9zgpWOlI7303qxAzm5rYkY4Dtz2CmStX8sMxY3b5uP369aOhoYGGhgYmTJjAD3/4Q8aNG8fCCm4+dfORlUYmA9n5aLdtc/ORVYxsLWFLmhS2RDCnpWWX+xZeeuklXn755dz64sWLOfjgg3nzzTdzSWHr1q0sXbp0l87T25wUrDR8S6pVqPxaQla2trArWltbOeOMM3LzKS9btoxrrrmG++67j8suu4xJkyYxefJknnrqqV06T29z85GVxrRp8MtfwsqVybtvSbUK8fC6dblaQtaWCB5au3aXmpCmTJnS4Q/+sGHDeOKJJ95XPnfu3A6XS81JwUpn2jRoa3NCsIrSfPjh5Q6horj5yMzMcpwUzKwmlWJWyUq0q9/bScHMas6gQYNYt25dn0sMEcG6desYNGhQj4/hPgUzqzkjR46kubmZN998s9yhlNygQYMYOXJkj/d3UjCzmtO/f39Gjx5d7jCqkpuPzMwsx0nBzMxynBTMzCynJElB0v6S5ktaJmmppC+n5TdIelHSC5IelLRXKeKxMpk9G15+OXk3s4pUqprCNuDSiBgLfBS4SNJYYB4wPiImAiuAy0sUj5Xa7Nlw3nnJJDvnnefEYFahSpIUImJ1RCxKlzcAy4H9IuIXEZEOncnTQM/vo7LKdv/9Xa+bWUUoeZ+CpFHAR4D2M1ifDfx3qeOxEpk8uet1M6sIKuUTf5IGA48D10bEA3nlVwL1wGeig4AkTQemA9TV1U25++67SxSx9ZqWFli1itaRIxnc3Az77Qf77lvuqGwntLa2Mnjw4HKHYT101FFHNUZEfXfblSwpSOoPPAL8PCJuzCs/EzgPODoiNnV3nPr6+njuueeKFqcVycKFcPTRZGbOpOHqqz18dhXKZDI0eB6MqiWpoKRQkieaJQm4DVjeLiEcB3wNOLKQhGBVzPMpmFWFUg1z8THg74HfSFqcll0B3AwMBOYleYOnI+L8EsVkpeb5FMwqXkmSQkQ8CXQ0W/vPSnF+qxALFyZ9CwsXOjGYVSg/0WylkfYpsGpV8p5OXG5mlcVJwUojk0majiB5z2TKGY2ZdcJJwUpjn31gx45keceOZN3MKo6TgpXGunWwW/q/2267JetmVnGcFKw0Ghpg4MBkeeDAZN3MKo6TgpVG9jmF/fbzcwpmFczTcVrp+DkFs4rnmoKZmeU4KVjpeJIds4rnpGCl4Ul2zKqCk4KVhifZMasKTgpWGn/3d12vm1lF8N1HVhrTpyfvEtxyyx/XzayidJkUJO0OnAicAEwC9gLeBp4nmTrzp3lzLJt1bfr0ZMwjP7hmVrE6bT6SdD7wKsmsaK8A1wLnp++vAF8CXk23MzOzGtBVTeHPgUMjoqWDzx4ErpM0HLi0KJGZmVnJdZoUImJGdztHxGqg2+3MzKw6FHT3kaSxkurS5SGSviXpm5L2KG54ZlYJVre1sWTjRlqyc2JYzSr0ltSfkHQyA9wA/DXwUeCWYgRlNSp/Ok6rKjObmmiLYObKleUOxYqs0KQwKiJekiTgM8BngVOATxQtMqstno6zaq1ua2POmjUAzGlpcW2hxhWaFDZLGgIcCvw+ItYCbcCgokVmtSWTgS1bkuUtWzwdZxWZ2dTEjggAtru2UPMKTQr/D/gVcAcwNy07BPhdEWKyWtTQAAMGJMsDBvhZhSqRrSVsSZPClgjXFmpcQUkhIr4CXAlcEBE/SIt3AF8pZH9J+0uaL2mZpKWSvpyW7y1pnqSX0/cP9uRLWBXwJDtVKb+WkOXaQm0reOyjiPgF8FtJH03Xn4uIXxW4+zbg0ogYS9JBfZGkscDXgV9GxF8Av0zXrVZNmwb77uuEUEUeXrcuV0vI2hLBQ2vXlikiK7aCxj6S9GckdyBNBgIYLOkU4LiIOLe7/dPnGVanyxskLQf2A04CGtLN7gAywGU79xXMrFiaDz8cgPHPPgutrYzbYw+WHHpomaOyYiq0pnAL8CgwBNials0DjtnZE0oaBXwEeAaoSxMGQAtQt7PHM7PiWrxhA0s3bQJg6aZNvLBhQ5kjsmJStKsadriRtA74UETskLQ+IvZOy9+OiL262T3/OIOBx4FrI+KB9vtLeisi3tevIGk6MB2grq5uyt13313oKa3CtLa2Mnjw4HKHYTth2aZNvLtjByO3b6e5Xz/+ZLfdGLuHn1utNkcddVRjRNR3t12hQ2evIRkLaUW2IO0T+H2hAUnqD9wP3BURD2SPK2l4RKxOx1F6o6N9I2I2MBugvr4+GnznStXKZDL4+lWPxRs2cFFjIwCzWluZkSb056dMYeKQIeUMzYqk0OajWcAjks4Cdpd0GnAP8K+F7Jw+9HYbsDwibsz76GHgjHT5DOChAuMxsxI4ffnyDsu/0Em5Vb+CagoRcXvahHQe8BrJD/jVEfHTAs/zMeDvgd9IWpyWXQF8G7hX0jnASuBzOxO8mRXX8rQvodByq34Fz7wWEQ/Rw7/kI+JJQJ18fHRPjmlmxbc9r6kvk8kQbvqreQUnBUnHktyS+p5ewoj4Rm8HZWZm5VHocwo/IGnamQ/k1xu7v3XJzMyqRqE1hS8AkyLitWIGY2Zm5VXo3UdrgbeLGYiZmZVfpzUFSQfmrX4XuEvS9STPLORExKtFis3MzEqsq+aj35L0GeTfNfSpdtsE0K+3gzIzs/LoNClERMEjqJqZWW0o6Idf0n7t5zqQ9EFJI4oTltUkz9FsVvEKrQ38FBjZrmwk8GDvhmM1y3M0m1WFQpPCmIj4TX5Buv6XvR+S1STP0WxWFQpNCm9K+vP8gnR9Xe+HZDXJczSbVYVCk8LtwP2SPiVprKRPA/cBtxYvNKspnqPZrCoU+kTzt0lmXJsF7E8yUuqtwI1d7WT2HtOmQVubE4JZBSt06OwdwA3py8z6mMUbNrCotZW9N2zw5Do1bmdGSf0wMIn3j5J6e28HZWaV5fTlyzmLZHKdJYceWu5wrIgKHSX1CuAbwPO8f5RUJwWzGrZ4wwaWppPqLN20iRdcW6hphdYULgEOjYgXihmMmVWe9lNyurZQ2wq9++hd4MViBmJmlSe/lpCVrS1YbSo0KVwNfF/ScEm75b+KGZyZlVf7WkLWFzopt+pXaPPR3PT93Lwy4VFSzWra8na1hO7KrfoV+pf+6PR1YN4ru25WGA+IV3W2NzRwwfDhDFAygv4AiQtHjGC7n0ivWQUlhYhY2dmr2AFajfCAeFVpdVsbc9asYUsk07FviWBOSwstbW1ljsyKpcvmI0nXdHeAiPhGd9tIup1kgp43ImJ8WjYZ+BEwCNgGXBgRzxYStFWhjgbE85PNFW9mUxM70oSQtT2CmStX8sMxY8oUlRVTdzWF/Qt4FWIucFy7su8A34qIySTPQHynwGNZNfKAeFXp4XXrcrWErC0RPLR2bZkismLrsqYQEWf1xkki4glJo9oXAx9Il4cCr/fGuaxCZQfEW7nSA+JVkebDDwdg7wULAPhgv36sP+KIcoZkRdZpUpBUFxFrujtAodt14BLg55JmkdRYDu/BMayaeEC8qjRv/Xre2r4dgLe2b+dX69fz8b33LnNUVixd1RR+Jelx4MfAM+mgeACkzyccCvwD8NfA+B6c+wLgKxFxv6TPAbcBf9PRhpKmA9MB6urqyHiClqrV2trq61dllra2MgsYuX07s1pbeeGpp9ht8OBu97PqpGjXXpj7QBpA8kN8Psntp68CG4Ah6fpvgVuA2yJiS7cnSpqPHsnraH4H2CsiQpKAdyLiA10cAoD6+vp47rnnuv9mVpEymQwN7k+oGvPWr+fYF5LRbWa1tjIjTQa/nDjRtYUqI6kxIuq7267TmkL6Q/8D4AeS9gcmAHsBbwEvRMSqXYzxdeBIIAN8HHh5F49nZr3s1KVLOyw/ZelS9y3UqELnU3iNZGKdHpH0E6ABGCapGfgm8CXg/0jaHdhM2jxkZpUj25dQaLlVv4LnU9gVEXFaJx9NKcX5zaxnIq+pL5PJvGfdapMHtDMzsxwnBTMzyykoKUia0Un5V3s3HDOrRKvb2liycaPHPOoDCq0pdDa+0VW9FYiZVa6ZTU20pWMeWW3rbkC8j6eL/SQdRTKHQtaBJM8tmFkNy46U+i/AnJYWrj7gAPYdOLDcYVmRdHf30W3p+yDg9rzyAFqAi4sRlJlVjvyRUj1Cau3rbkC80QCS7oyIfyhNSGZWKTqbT8G1hdpV6CQ7uYTgOZrN+o6u5lOw2lTo3UeHSFooaSOwNX1tS9/NrEZ5PoW+p9Anmu8A/gs4G/CM3WZ9RHY+BfATzX1FoUnhAODK6GxIVTMzqwmF9gk8CBxbzEDMzKz8upp57cckt54CDAQelPQkya2oOb4rycysdnTVfPTbduvLihmImZmVX1eT7HyrlIGYmVn5FdTRnDfcRXttQHNE+KZlsxqWHRDvL9va/NBajSv07qPbgBHp8jpgn3T5DWBfSS8An48IT6lpVoNmNjVxkIe46BMKvfvoNuBmYK+IGEEyV/NNwI/S5V8D/1aUCM2srLJDXUAyIJ6Hz65thSaFLwOXR8S7AOn7VcAlEbERuBSoL06IZlZOHQ2IZ7Wr0KSwEZjarmwKf3y6eUevRWRmFaOzAfFcW6hdhfYpfAP4haSHgdeAkcCn+ePQ2UcD9/V+eGZWTl0NiOe+hdpU6CipdwKHAS8CQ4EVwLS0nIh4JCK+VLQozawsPCBe31Pw0NcRsSwiZkbEBRFxTUQU/DCbpNslvSFpSbvyiyW9KGmppO/sTOBmVnzNhx/OBcOHM0DJpIsDJC4cMeI9A+VZbelqmIvZETE9Xc4f8uI9ChzmYi7wA+DOvOMfBZwETIqINkl/uhNxm1kJeJKdvqerPoXf5S23H/Jip0TEE5JGtSu+APh2RLSl27yxK+cws97nPoW+p6thLq7PWy7GkBdjgCMkXQtsBmZExK+LcB4z66Gu+hScFGpToXcfIekY4PPAn0bEpyXVAx+IiF/twrn3Bj5KcrvrvZIO7GjOBknTgekAdXV1ZDKZHp7Syq21tdXXr4r8R95yKzA/u7Jli69jjSp07KOLSR5guxU4JS1+l+Qp5572ODUDD6RJ4FlJO4BhwJvtN4yI2cBsgPr6+mjw7E9VK5PJ4OtXnXzt+oZC7z66BPibiPg2f3xQ7UXgw7tw7p8CRwFIGgMMAHyfm5lZGRXafDSE5KE1+ONdSP2BLYXsLOknQAMwTFIz8E3gduD29DbVLcAZnu7TzKy8Ck0KTwBfB67NK/sn8poYuxIRp3Xy0ekFnt/MzEqg0Oaji4GTJTUBQyS9BHwO+GqxAjOzypGdT8FjHtW+Qoe5WE1yh9CpwBeAM4BDI6Klyx3NrCbMbGqizSOk9gk7M8xFRMQzEfGfEfF0RHhkVLM+wPMp9C1d9ilIWkAnw1tkRcRf92pEZlZROppPwQ+u1a7uOppvzVsW8EPgwuKFY2aVxGMf9T1dJoWIuCN/XdKN7cvMrHZ57KO+p+A+BTPrezyfQt/jpGBmnfJ8Cn1Pdx3NH2+/fToPgrIFuzAgnplVOPcp9D3ddTTf1m59HcnwFFkBHNirEZlZxXCfQt/TXUfz6FIFYmaVx/Mp9D0Fz6dgZn1Pft9BJpMhPHR2zXNHs5mZ5TgpmJlZjpOCmZnlOCmYmVmOk4KZmeU4KZiZWY6TgpmZ5TgpmJlZjpOCmZnlOCmYmVmOk4KZmeWUJClIul3SG5KWdPDZpZJC0rBSxGJmZp0rVU1hLnBc+0JJ+wPHAr8vURxmZtaFkiSFiHgCWN/BR98DvkYyL4OZmZVZ2YbOlnQSsCoinpfU3bbTgekAdXV1ZDKZ4gdoRdHa2urrV6V87fqGsiQFSXsAV5A0HXUrImYDswHq6+ujwWO6V61MJoOvX3XytesbynX30UHAaOB5SU3ASGCRpH3LFI+ZmVGmmkJE/Ab40+x6mhjqI2JtOeIxM7NEqW5J/QmwEPiwpGZJ55TivGZmtnNKUlOIiNO6+XxUKeIwM7Ou+YlmMzPLcVIwM7McJwUzM8txUjAzsxwnBTMzy3FSMDOzHCcFMzPLcVIwM7McJwUzM8txUjAzsxwnBTMzy3FSMLNurW5rY8nGjbS0tZU7FCsyJwUz69bMpibaIpi5cmW5Q7Eic1Iwsy6tbmtjzpo1AMxpaXFtocY5KZhZl2Y2NbEjAoDtri3UPCcFM+tUtpawJU0KWyJcW6hxTgpm1qn8WkKWawu1zUnBzDr18Lp1uVpC1pYIHlrr6dRrVUmm4zSz6tR8+OG55UwmQzQ0lC8YKwnXFMzMLMdJwczMcpwUzMwspyRJQdLtkt6QtCSv7AZJL0p6QdKDkvYqRSxmZta5UtUU5gLHtSubB4yPiInACuDyEsVi5TJ7Nrz8cvJuZhWpJEkhIp4A1rcr+0VEbEtXnwZGliIWK5PZs+G88+APf0jenRjMKlKl3JJ6NnBPZx9Kmg5MB6irqyOTyZQoLOs1EsyaRevIkWRmzUrWfR2rSmtrq//t9QFlTwqSrgS2AXd1tk1EzAZmA9TX10eD75WuPitWwIwZZGbNomHGDLjlFvB1rCqZTAb/26t9ZU0Kks4EPgUcHdHusUmrLdOnJ+9SkhCy62ZWUcqWFCQdB3wNODIiNpUrDiuh6dOTJiP/tWlWsUp1S+pPgIXAhyU1SzoH+AEwBJgnabGkH5UiFjMz61xJagoRcVoHxbeV4txmZlY4P9FsZmY5TgpmZpbjpGBmZjlOCmZmluOkYGZmOU4KZmaW46RgZmY5TgpmZpbjpGBmZjlOCmZmlqNqG5xU0pvAynLHYT02DFhb7iCsR3ztqtsBEfGh7jaquqRg1U3ScxFRX+44bOf52vUNbj4yM7McJwUzM8txUrBSm13uAKzHfO36APcpmJlZjmsKZmaW46RgPSZpX0l3S3pFUqOkn0kaU+64rGckXSlpqaQX0ilyD+ti27mSTkmXj0j3WyzpT0oXsRVDSabjtNojScCDwB0R8fm0bBJQB6xI13ePiG3li9IKJWka8CngkIhokzQMGFDg7l8Ero+I/yhagFYyrilYTx0FbI2IH2ULIuJ5oJ+kBZIeBpYBSPqqpCXp65K0bE9Jj0p6Pi0/NS3/tqRl6V+rs8rwvfqq4cDaiGgDiIi1EfG6pCmSHk9rgj+XNDx/J0nnAp8DZkq6qwxxWy9zTcF6ajzQ2MlnhwDjI+J3kqYAZwGHAQKekfQ4cCDwekScACBpqKR9gJOBv4yIkLRX0b+FZf0C+IakFcBjwD3AU8D3gZMi4s00cV8LnJ3dKSJulfRXwCMRcV8Z4rZe5pqCFcOzEfG7dPmvgAcjYmNEtAIPAEcAvwGOkfSvko6IiHeAd4DNwG2SPgNsKkfwfVF6baYA04E3SZLCeSTJf56kxcBVwMiyBWkl4ZqC9dRS4JROPtvY3c4RsULSIcDxwL9I+mVEXCPpUODo9Nj/CHy8twK2rkXEdiADZCT9BrgIWBoR08oamJWUawrWU78CBkqani2QNJGkFpBvAfC3kvaQtCdJ89ACSSOATWnn5A3AIZIGA0Mj4mfAV4BJpfgiBpI+LOkv8oomA8uBD6Wd0EjqL2lcWQK0knFNwXokbfM/GbhJ0mUkzT5NwE/bbbdI0lzg2bTo1oj4X0mfAG6QtAPYClwADAEekjSIpP/hqyX5MgYwGPh+2o+zDfgtSVPSbOBmSUNJfi9uIqklWo3yE81mZpbj5iMzM8txUjAzsxwnBTMzy3FSMDOzHCcFMzPLcVKwmiDpnyX9R7o8SlJI2j1d/29JZxTpvFdIurUYxzYrBycFqyiSmiS9K6lV0pp0iObBu3LMiPhkRNzRC7E1SGpud+zrIuLcXT12u/NckX7/VkmbJW3PW/czAlZUTgpWiT4dEYNJBtarJxlzp89IE83g9L/B+cDC7HpE+IliKyonBatYEbEK+G+SQdmQNELSw5LWS/qtpC8VchxJmXSIZySdKelJSbMkvSXpd5I+mbftWZKWS9og6VVJ56Xle6axjMj7q31EfrNVut2J6YQzb6fnPTjvsyZJM9Jhwd+RdE/69PZOkXSLpH9tV/YzSReny82SLku/x1uSbpM0sF2Mz6cxPilp/M7GYLXLScEqlqT9SQbM+9+06G6gGRhBMmDedZJ6MmDeYcBLwDDgOySjsir97A2SyWY+QDLk9/ckHRIRG4FPkgz3nf2r/fV28Y4BfgJcAnwI+BnwX5LyJ6v5HHAcMBqYCJzZg/jvAE7LxiypDmhIz531ReAY4C+AccDl6bZTgX8HzgX2AW4nGVqk0Al1rMY5KVgl+qmkt4EngcdJfvz3Bz4GXBYRmyNiMXAr8A89OP7KiPj3dFTQO0gmmKkDiIhHI+KVSDxOMs9A+0H+OnMq8GhEzIuIrcAs4E+Aw/O2uTkiXo+I9cB/kQw8t1Mi4imSsaYa0qLTgMciYm278zSnZdel20AyntG/RcSvI2J7RNyelk/d2TisNjkpWCX624jYKyIOiIgLI+JdktrB+ojYkLfdSmC/Hhy/JbsQEdk5GwYDSPqkpKfTJqq3SWoqwwo87og0puyxdwCvtYuxJW95U/a8PXAncHq6fDrw43afv5a3vDKNDeAA4LK06eirBEGIAAABkElEQVTt9DsOp2f/Ha0GOSlYtXgd2FvSkLyyPwNW9dYJ0nb3+0n+wq+LiL1ImoCyTUvdjR75OsmPbvZ4AvbvzRjz/Bg4WdJHgINIah359s9b/rM0NkiSxbfSpJt97RER9xYhRqtCTgpWFSLiNZLpIa+XNCidu+EcoDcnix8ADCSZeWxb2gF9bN7na4B90mGkO3IvcIKkoyX1By4F2tK4e1VErASeJ2n++s+I2Nxuk3+UtJ+SKU4vJ5lJDZL+hIskTVVisKRPpx3pZk4KVlVOA0aR/NX7IPDNiHistw6eNk39E8mP+1vAF4CH8z5/kaQz99W06WVEu/1fImnK+T6wFvg0ye21W3orxnbuACbw/qYj0jgfA14h6VS/Lo3xaZK5K/4vyXdcwR+bocw8n4JZtUrvvLoNODDy/iGnD9idHhGZcsVm1cs1BbMqlN5C+mXg38N/2VkvclIwqzKSJpA0/ewN3FzmcKzGuPnIzMxyXFMwM7McJwUzM8txUjAzsxwnBTMzy3FSMDOzHCcFMzPL+f9ypq2v0gKW9QAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "plt.scatter(df[\"type\"][df[\"type\"]==\"Cross\"],df[\"height\"][df[\"type\"]==\"Cross\"],marker=\".\",color=\"r\",label=\"Cross\")\n",
    "plt.scatter(df[\"type\"][df[\"type\"]==\"Self\"],df[\"height\"][df[\"type\"]==\"Self\"],marker=\"^\",color=\"c\",label=\"Self\")\n",
    "plt.margins(1)\n",
    "plt.grid()\n",
    "plt.ylim(11,25)\n",
    "plt.xlabel(\"Pollination Type\",fontsize=12)\n",
    "plt.ylabel(\"Height (Inches)\",fontsize=12)\n",
    "plt.legend()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Summary\n",
    "The mean heights of the cross- and self-pollinated plants were 20.2 (95% CI: 18.6, 21.7) and 17.6 (16.0, 19.1) inches respectively, with a difference in height of -2.6 inches. A 95% CI for the difference in height did not contain zero.\n",
    "## Software versions:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "pandas: 0.24.2\n",
      "matplotlib: 2.2.4\n",
      "numpy: 1.16.3\n",
      "statsmodels: 0.9.0\n"
     ]
    }
   ],
   "source": [
    "print(\"pandas: \"+str(pd.__version__))\n",
    "import matplotlib\n",
    "print(\"matplotlib: \"+str(matplotlib.__version__))\n",
    "import numpy as np\n",
    "print(\"numpy: \"+str(np.__version__))\n",
    "import statsmodels\n",
    "print(\"statsmodels: \"+str(statsmodels.__version__))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.15rc1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}