KNMI LOESS TG

{
 "metadata": {
  "name": "",
  "signature": "sha256:e9fa647959fcfc66d16a54b6fbff46073f16597ea6cc222e30de541053408583"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import pandas as pd\n",
      "import matplotlib.pyplot as plt\n",
      "import statsmodels.api as sm\n",
      "%matplotlib inline"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 1
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "url = r'http://www.knmi.com/klimatologie/onderzoeksgegevens/homogeen_260/tg_hom_mnd260.txt'"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 2
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# load data from the website\n",
      "df = pd.read_csv(url, skiprows=28, delim_whitespace=True, date_parser=lambda x: pd.datetime.strptime(x, '%Y%m%d'), \n",
      "                 header=None, names=['STN', 'Date', 'TG'], parse_dates=['Date'], index_col='Date')\n",
      "\n",
      "# calculate the yearly mean\n",
      "# not weighted for the length of the month!\n",
      "dfyear = df.groupby(df.index.year).mean()\n",
      "\n",
      "# loose 2014 since its not a full year\n",
      "dfyear = dfyear[dfyear.index < 2014]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 3
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# Lowess fit, after Cleveland, W.S. (1979)\n",
      "lowess = sm.nonparametric.lowess\n",
      "z = lowess(dfyear.TG, dfyear.index, frac=0.5)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 4
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# plot the results\n",
      "fig, ax = plt.subplots(figsize=(10,5), subplot_kw={'xlabel': 'Year', \n",
      "                                                   'ylabel': 'Temperature ($^\\circ$C)',\n",
      "                                                   'yticks': range(7,12),\n",
      "                                                   'ylim': (7,11.5)})\n",
      "\n",
      "ax.step(dfyear.index, dfyear.TG, 'k-', where='mid')\n",
      "ax.grid(axis='y')\n",
      "\n",
      "ax.plot(dfyear.index, z, color='#0080FF', lw=2)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 5,
       "text": [
        "[<matplotlib.lines.Line2D at 0x13465f60>,\n",
        " <matplotlib.lines.Line2D at 0x132e9438>]"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAmoAAAFCCAYAAABBxfG5AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xl8VNXdx/HvZGUJIQlCQLYgsiaBsIhGKsUqoM+jKFit\niFqkdSmVVmsVtbXFahVbtGJbl7qV4o6KPm5oFQER2SOLCAImhn2RAAmQ/Tx/3EwgISGTyZ25d2Y+\n79drXsm9meU3ZyaZb8695xyPMcYIAAAArhPldAEAAACoG0ENAADApQhqAAAALkVQAwAAcCmCGgAA\ngEsR1AAAAFwqxukCfJGVlaXVq1c7XQYAAECD+vfvry+//NKW+wqJHrXVq1fLGMMliJc//vGPjtcQ\naRfanDaPhAttTptHwsXOzqWQCGoAAACRiKAGAADgUgQ11Gn48OFOlxBxaPPgo82DjzYPPto8tHmM\nMa5f69Pj8SgEygQAALA1t9CjBgAA4FIENQAAAJciqAEAALgUQQ0AAMClCGoAAAAuRVADAABwKYIa\nAACASxHUAAAAXIqgBgAA4FIENQAAAJciqAEAALgUQQ0AAMClCGoAAAAuRVADAABwKYIaAACASxHU\nAAAAXIqgBgAA4FIENQAAAJciqAEAALgUQQ0AAMClghbUJk6cqNTUVGVmZlbvmz17ttLT0xUdHa1V\nq1YFqxQAAGCzlJQUeTyek15SUlKcLjPkBC2oXXfddZo7d26NfZmZmZozZ46GDRsWrDIAAEAAFBQU\nyBhz0ktBQYHTZYacmGA90DnnnKO8vLwa+3r37h2shwcAAAg5nKMGAADgUgQ1AAAAlwraoc+mmjBh\ngtLS0iRJSUlJysrK0vDhwyVJ8+fPlyS22WabbbbZZtvF2959bqnHrm3v97VP8bKDxxhjbL/XeuTl\n5eniiy/W2rVra+w/99xzNX36dA0aNKjO23k8HgWxTAAA0Ei+fFZHyue5nc8zaEFt3LhxWrBggfbt\n26fU1FTde++9SklJ0eTJk7Vv3z61bt1aAwYM0AcffHBikRHywgIAEKoIaseEZFBrikh5YQEACFUE\ntWPsfJ5RttwLAAAAbEdQAwAAcCmCGgAAgEsR1AAAAFyKoAYAAOBSBDUAAACXIqgBAAC4FEENAADA\npQhqAAAALkVQAwCgHikpKfJ4PDUuKSkpTpd1glCp01/h/vxOhiWkAACoR12fP278TKqrppSUFBUU\nFFRvJycna//+/UGtwZ/r+Ho7N74OXnbWFmPLvQAAAFepHco8Ho9DlaApOPQJAADgUgQ1AAAAlyKo\nAQCAsHCwWHp4sZR3wOlK7MM5agAAIKTlHZBmLJGeyZGKSqXthdIjo5yuyh4ENQAAEJKWbJMe+UJ6\n42upsmqQ5blp0qjuTlZlL4IaAABV6prSAu5SUSm9vVHSxEXKftbaFxMljc+Ubj1LGtDB0fJsxzxq\nAABUCeRcYIHkhroDPY9aUYnR819Kjy6RtlRl6aRm0o2DpMlDpI6J/lQdGMyjBgAAmqx2D6IU+Ilx\nG2tnoaQf3a/Of5MKiq193ZKk3Jcma+tnf1dCnKPlBRw9agAAVHFDz5Q//K3bzhn/7W679XutEZwv\nrJVKK6x9Z3WSbsuWxvSWYqLd9zp40aMGAADCjjHS/Dxp+hfS+5usfR5J+vpNff7XsTq7s4PFOYQe\nNQAAqtCj5kyPWnml9MZ66a+LpZU7rX3NY6TrsqRbzpJ6nsJanwAAAEF1uFR6Nkf625Jjk9S2bSHd\nPESadIZ0Sgtn63MDghoAAAiq3UXSP5ZJ/1x+bIBAjxTr/LNr+0vNY52tz00IagAAoFHqGy3aoDY9\ndOM70szVUknVAIHsTtLtZ0uje0nRLGx5As5RAwCgCueoBaYNlmyzzj97c32l5ImSR1Ywu/1saWiX\nhm9v57l0wcA5agAAwNUqjfTeN1ZA+yy/amdFmX5+RrxuO1vqfYqj5YUMetQAAKgSyB61QE4u66Ye\ntZJy6aW1VkD7ep+1r3W8NTjgwUs7yBTubODZ+F+7W9CjBgBAiCkoKKgzbISLg8XSv1ZKjy6VdhRa\n+zolWutvXj9QahUvPRlbcsJzdttKCG5DUAMAoInqWsw9UsLHjkJpxhLpyZXSoRJrX0Y76Y6zpSsz\npNjoY9etq03CKawGAkENAIAmqt1bFgnh4+u90vTF0qw1UlmltW94mjVA4MLTpQhogqAgqAEAAJ8t\n3irpyrfU93Fr2yPpx32tgDako5OVhScGEwAAUMWuk/IDvVyTLzX5ch1f91Ua6d1vpL98Ln2+1doX\nHy1NyJJ+e7Z0eoo7ancLBhMAAICAq2sEZ1Iz6cBH9+u7t3+v1ISmP0ZycnKNQ8WRdH6fL+hRAwCg\nCj1q1r6DxUZPrThxBOdvzpJ+PlBKbBa4z+Vgt10g0KMGAEAYcFtv0s5CSec/qM5/a3gEJ4KDHjUA\nAKoEu0fNn8f393Ynq2njPmsE53/WSKVVa3AO6ypNGVr3CM5Afi7To1YTPWoAAESqTmdpzKvS2xsk\nI2sEp9a/riUP/1hndnK4NkiiRw0AgGqR0KPmXYPzL4ulRVVrcHpHcN6WLfU8xdmF6elRq4keNQAA\nQkh9a4Y2KDpOz+VYhzi9Izh1tEC/G5msyUNkywhO2I8eNQAAqoRCj1pjPxMPFktPrZSmzNkhtTpV\nktQ5UfpNtnTrD1vJlBQGvIbGoEet1n0R1AAAkai+nqmGRl2GSlDbfkiasVR6coVUWGrty2xnrSDg\nHcEZ7EO0viCo1bovghoAIBLZeT6Ym4Lauj3W4c2X1h5bg/PcNGn5o5er6MvXa1y3rmBKUGs6ghoA\nAE0UTkHNGGnBd9YKAu9vsvZFeaw1OH+bLZ3RiDU46+pprC2Q870R1GpiMAEAACGqolJ682sroC3f\nYe1rHiNNHCDdepbU3Y81OENl+Sa3TRYcKAQ1AABCzJEy6d9fSg9/IX1b1fnVprl08xDpl2dIbVs6\nW18w1HXINhwR1AAACBF7DksaPlVd/iZ9f9Ta1z3ZGsE5IUtqEetoeQgAzlEDAESkUDlHLSUlRQWe\nNtLZt0n9fyrFNpckDelojeAc01uKjmr003AtN45EbSzOUQMAIAIs3ioVjHxWnj5j5P3Yv7inFdB+\n0OXENTgRfuhRAwBEJLf2qJVXGP3fRmuAwBfbrP1x0dI1/awlnvq0bXTJIYUetZroUQMAwAWOlEka\nfKP6/FPaVHWefHIzqeCj+/Xd279Xe5Z4ikj0qAEAIpJbetT2Hpb+udy67Dti7UtLkn5zlnTdAKlV\nfGR9BtKjVlPQTj+cOHGiUlNTlZmZWb1v//79GjFihHr27KmRI0fqwIEDwSoHqFdKSoo8Hk+NS0qK\nH5MRAcBJbNwn3fSu1OVR6d4FVSFt+zK99mNp02Rp8plSQpzTVcJpQQtq1113nebOnVtj37Rp0zRi\nxAh98803Ou+88zRt2rRglQPUq6CgQMaYGpeGZukGAF8YI332nXTJK1Kff1qLpReXSxf1lBZMkPT0\nmbo8XYoJo1GcaJqgHvrMy8vTxRdfrLVr10qSevfurQULFig1NVW7du3S8OHDtWHDhhOLdFF3JsJf\nqC1VAsA/wTz0WV4pxWZeriG3ztay7da++Gjp2v7WHGi9T/H9vsMdhz5rcnQwwe7du5WamipJSk1N\n1e7du50sB0AYqmvdwnBdagbuU1QqPZcjPbpE0hVWSGvTXJp0hrWCQCoDBNAA14z69J4LBAB28h7K\nPh5/axBwCe119yfSkyukguKqfd9v0uPX9NBPWUEAjeBoUPMe8mzfvr127typdu3a1XvdCRMmKC0t\nTZKUlJSkrKwsDR8+XJI0f/58SWKbbVu2vftq//z4n7mpXrZPvu3dx+sZmtutWrVSUVGRvBISEvTO\nO+84Xp9X7e3n58zXa19JuiVPDy6SlDtfGe2kP103XGPTe6vP2E+07HPer/62r1vqq6ve+fPnKy8v\nT3Zz9By1O+64Q23atNGUKVM0bdo0HThwoM4BBW467ozwxzki4YXXM7T5Mu2FXfft7+2MkaK6/0gX\n3jdPH2z27qzUZX2jdFu2lN3Z98fj/co5aifcV2ODWnFxsTwej+Lj4xv1QOPGjdOCBQu0b98+paam\n6k9/+pMuueQSXXHFFcrPz1daWppee+01JSUlnVikixof4Y8/lOGF1zN0+HI+oZuCWlmFNHu9NH2x\nlLPL+lnzGGniAOmf15wu8/3mOm/X2Joi7f1KUKt1Xw0FtcrKSr311lt6+eWXtXjxYlVWVsoYo+jo\naGVnZ2v8+PG69NJLA3rOh5saH+GPP5ThhdczdAT7w9jvoBbfSo/ML9SjS6X8g1U7i3brvotT9YvB\nUpsW9oaNSHu/EtRq3VdDQW3YsGE655xzNHr0aGVlZVX3pJWUlCgnJ0f/93//p0WLFmnhwoW2FFRn\nkS5qfIQ//lCGF17P0OH2oLazUHpsqTTt4wKpebIkqVcba/3NG85qJlNWXH1dgpr/CGq17quhoFZS\nUtLgYU5frtMUbmp8hD/+UIYXXs/Q4dagtmGfdXhz1hqptMLa94Mu0m+zpYt7SVEeexdlj/T3K0Gt\npqiGrpCfn69FixadsH/RokXasmWLJAU0pAFwh9pLa7GsFsLdkm3SmFelvv+Uns2RyiqksX0kPXOW\nPrtOuqS3FdKAQGowqN1yyy1KTEw8YX9iYqJuueWWgBQFwH1qL63FslrhLxLXvTVG+mCT9MN/S9nP\nSm9tkOKipRsGShtvlt64QtK2pU6XiQjS4Dxqu3fvVr9+/U7Y369fP+Xm5gakKACA8yJpsuCKSun1\n9dKDi6TVVYvktI63VhD41ZlSe1YQgEMaDGoHDhyo92fFxcX1/gyA81g+CTi50gpp1mrpoc+lTVW/\nFh0SpFvPkm4cLCVyZg8c1mBQGzx4sP71r3/phhtuqLH/6aef1qBBgwJWGICmi6QeEaAxSsolnfEL\nnf6YtPWQta9bkjRlqPTTLKmZaxZYRKRrcNTnrl27NGbMGMXFxVUHs5UrV6qkpERz5sxRhw4dAl+k\ni0ZyRJpI7JEJp1FXdj6XQM4Q76/a78+63pvh9HoGW7DbLhgj+4rLpWdWSdMWSdsLrX3pbaW7z5Gu\nSJdiGjxz27d2YdSn/xj1Weu+fFmZwBijTz/9VOvWrZPH41F6erp+9KMf2VKAL9zU+JEmEv9ohNNz\nDveg5u+Ho9O1h8o/QOEU1EorpGdXSfd/Ju2oCmjavUazb+6nsX0aN3qToBZYBLVa9+VLUHOamxrf\nX6Hyh7m2SPyjEU7PmaDm275g/36GynssHIJaRaX00lrpj/Ol3KpTrvunSn/4oXRZepRMZaUtdRLU\n7ENQq4mj8EHCuUKAe/H7GX6MsabW+P2n0vq91r7ep0j3n2vNhebxVF0JcDmfgtqRI0fUokWLQNcC\nAHC55OTkGiHWbUcGjJHe/cbqQfMulN61tXTvcOnqflK0D+egAW7SYFD7/PPPFR0drYqKCg0dOjQY\nNQEAHFDX4Iza6hqs4QbGSHM3S3+YL63YYe3rkGANErh+oBTP8SOEqAbfumVlZRo6dKg++eSTYNQD\nwAXqO2cL4a2uQ8BuZ4z00RZp6gJrySdJSm0p3fkD6cZBUvNYZ+sDmqrBoNanTx/NmzdP6enpwagH\ngAuE4gc2Iosx0sffWoc4v6gKaKe0sOZBm3SG1IKAhjDRYFBLTU1VampqMGoBAKBBH38r3btAWpRv\nbbdpLt1RFdAS4pytDbBbg0HNGNPgOQi+XAcAAH8ZI324RdLPPteIWda+lObS7WdLNw8hoOHEgS7h\nclSgwaA2fPhwXXTRRbrkkkvUs2fPGj/buHGj3nrrLb333ntauHBhwIoEAEQm7yjO+xZKy3dI6ny2\nTmkh3ZYt/fIMqRVrcaKKm0Yf26nBoPbRRx/pxRdf1C9/+UutW7dOrVq1kjFGRUVFysjI0Pjx4/Xx\nxx8Ho1YALheqEzvDfSoqpdnrpQcXSWt2W/vatZT2zPmtcj+dTg8aIkajViaoqKjQvn37JEmnnHKK\noqOjA1bY8dw027C/QnW26VCtuynC6TkHcoZvNyzXZNfKBE7XHejH85U/Ndi6VI7Ho5Jyo1mrpWmf\nS5ur8v2praTfZks3DpZaxtn7eHas0lHXPlYm8F84tIFjKxNER0czsCDI6KFAuHHDe7quSVsjjdum\nYDlUIuns29RtxrG1OE9LtkZx/rQ/86AhcrHWZ5CE6n9XTj++E8LpObuxRy3Y64/a2WtiFze8x+x6\nvKbez85CacZS6YkVVWFNUnpb6a4fSD/JkGJqrSRgdw8ePWruEw5twFqfAICQtma39Lcl1oLppRVV\nO/Pm6727h+vC06vW4gQgn1c9q6ys1KxZs/SnP/1JkpSfn69ly5YFrLBIlJKSIo/HU+MS7odkaj/n\nlJQUp0sCECCVRnp/k3T+f6T+T0r//lIqq7AWSV/yM0n/Plf/04OQBhzP50OfN910k6KiojRv3jxt\n2LBB+/fv18iRI7VixYpA1xhyXZ518ber3Nf7CpRAd0H72wa12XneUzh0u3tx6JNDn42pIVD3c6hE\nmvml9I/l0jffW/taxko/GyD96kype4rv98Whz/AXDm3gyKHPpUuXKicnRwMGDJBkfTCWlZXZUgTQ\nVHUteWTnJMx1nXzeUAh0w0nzgJM27JP+sUyauVoqKrX2dUqUfjVEun6QlNTM2fqAUOBzUIuLi1NF\nRUX19t69exUV5fORUyCk1Q5XvoTAQIdHwI1KK6S3N0hPrZQ+yT22/4ddpclDpEt6nzhAIBzV/ufO\nuw9oLJ+D2uTJkzVmzBjt2bNHd999t15//XXdf//9gawNwEnU7rGjtw5OyjsgPb1SejZH2n3Y2tc8\nRrq6n7XEU78Im9mJ30XYxaegZozRsGHDNGjQIH3yySeSpLffflt9+vQJaHEA6le7x47eOgRbcbmk\njJ9oxCzpk28l77sxva1002ArpHF4E2ganwYTGGOUmZmpdevWBaOmE4TaSYR1YTCBb/cfyJN7m3Jf\ndl0nkIMe7FwpgMEEDCaojzFSzi7p+RzpxbVSQbG1Pz5a+nFfK6AN7ezfyM1wGkxgZ03h8BnYGOHQ\nBkEfTODxeDRo0CAtW7ZMQ4YMseWBgUhUVyCjJwyhIO+A9OIa6YW11iCBajtW6J8/H6xxGVJyc8fK\nA8KWz9Nz9OrVS5s3b1bXrl3VsmVL68Yej9asWRPQAr2PE0pJui70qPl2/274DzeQPU5NuV5Dt6NH\njR61ptZQ285C6Y2vpVe/khblH9vftoU0LkOaOEDK6hDcHi561MJfOLSBI9NzfPjhh7Y8IADAPnZP\nA7P9kBXOZq+XPs8/dt5Z8xhpTB/p6kzp/NOk2OgmFg7AJz4HtbS0tACWAQDwR1Ongak00qqd0rvf\nSO98Y33vFR8tXXC6dHlfaXQvqVW8XVWHNqbeQDD5HNTuvffeE/Z5PB794Q9/sLWgcFDff7gA4Aa7\niiRlXqWJb0tzN0s7i479rHnMsXD2vz2lRMLZCZh6A8Hkc1Br2bJl9X8QR48e1bvvvqu+ffsGrLBQ\nVtd/uADglN1F0udbpYXfWZPQrtsj6bIX9fyX1s87J0oX9bQu56ZJzWMdLBZADT4PJqitpKREI0eO\n1IIFC+yu6QThehIhgwnCazCBvxPQOj2YwN9znBhM4L9A1llWIcV1GqB/vZujxdus88w21XopW8RK\nR9Z/oOmTLtSI7lJmO/+m02hK3f7eV6DX9XVaOJxI31Th0Aa2/l74G9T279+vIUOGaPPmzbYUcjLh\n+gIR1EI3qLnhw8KuoBboNieo+Xbf/jzewWIpqefZemrOYq3eJa3YKa3eJZVU1Lxey1jprE7W/GY/\n6iZld5biYwI3SjHY9+XG19hf4RBSmioc2sCRUZ+ZmZnV31dWVmrPnj2cn4aIxTkqCJayCin/oNUr\ntul7afN+6Zv90ld7pK2HJP1ssW58t9aNvv9G437YU0M6Sud0kfq3j4z1NYFw5HNQe/fdd6vTYUxM\njFJTUxUby4kMgFswEs296josbox0oNhaF3NXkaT0KzRjibS90Apg+Qety45Ca2RmXZrFSMVbV+ma\nUQOV2U4afKo0sIOU1LyXXnosdHofANTP50OfU6ZM0UMPPdTgvkAI1y5PDn2G7qFPOzl9uIdDn77d\nd11KK6RDJVJhiXSwxPr+YLH1/cFia3mlgqPSI48/qzFX/Uz7j0r7jliX749K5ZU+1CGpY6LUI8W6\nnJ4i9Wgj9W0rdU+WYqIDN6mxLzj0aa9wOOzXVOHQBo6cozZgwADl5OTU2JeZmam1a9faUsjJhOsL\nRFAjqDXl8QhqjRsYUV4pFZVKh0ulw2XHvi+q5/LAw4/puht/pcJSK4gVlar6e+/X2ueCNVarOCk1\nQWqfIC36YLZunnC5OrSSurS2RmJ2aW2FtLiTTC7rb3sS1Nz52RIOIaWpwqENgnqO2hNPPKHHH39c\nW7ZsqXGeWmFhoYYOHWpLEUCw+DsyE/Uzxjo0V14pKa6lDhRb35dVSGWVkpK7aeM+6/vSCkkdh2hR\nvvV9WUXVvj5j9fLaqn2Vx/Z7b1N7u7RCKimXNOY/umK2FZiKy6WMh/eruFw6Wi4dLZOOVF3i7quq\npTHO/FX19BX1ifZY84y1ire+to6XWjer+hpvrX2Z3Ey641c/1+uznlFKc+mUFtalTYuaAcwz8Qr9\n/bXQ+SACEBwN9qgdPHhQBQUFuvPOO/XQQw9VJ8RWrVqpTZs2wSkyTJM0PWrh1bvjL7t7A4yxen8K\niq1zoAqOWl+9vUCTbr1D99z3Fx0utULMk8/O1I+v+qmOlFnhxxt6SsprBifv1/JK63L4aLFi4pr5\ndPjODaI8UkKcNQKy5XFfW8VZ+70/axVv7fv9Hb/Ws0/MUMJx1/H+zPu1WYzk8WFai0D2cNGjRo9a\nuAmHNnDk0KdkTeS6adMmFRcXV+8bNmyYLYWcTLi+QAQ1gpqvj1dppD2HrXUYdxZZJ59ff+s9mnTH\nfdpz5Nh5T95LsMNTlMcaVVh69LASW7VUTJQUG2WtB7ntu2/Vo/tpiou2tr9csVRDs89UbLS1RFFs\ntPTuW2/oJ5dfpthoKS5K1deNrfreux1XtS8+xvr++uuu0asvzVJctDWjfvNYKzzFR1tzhbWMs762\niLVu50uo8gp2ACGoEdQkd/xNclo4tIEjQe3pp5/WY489pm3btikrK0tLlixRdna25s2bZ0shJxOu\nLxBBjaDmfbyScqP8g1JugZR3oOpyUNp60BoBuP1Q4w7dtYyVkppZl+Tm1mG4xKpeoH/9Y7qm3vXb\n6l6lSTdM0Ksv/FvNY6pCTlXQ8X71Bqbjv0Z7pMSE5iotPqroKCuoeZ+Lk+eo2Ymg5v/jB/u+3Pj+\n8Zcb/iY5LRzawJGglpGRoeXLlys7O1tffvmlNmzYoLvuuktz5syxpZCTCdcXKJhBzdcJWn05hyvY\nQc2X2kMhqBWVWnNgeS9b9ktbCqRPv/xOnqSuaugR2zSXOiVKHVpJHRKk5//+gGb8+W61aym1bSG1\nbVl17lNzK2TVx65w7Ot7iqDmfFCzc4JmXx4/2PflxvePv8IhpDRVOLSBI0Ft8ODBWrFiRXVvWrNm\nzdS3b1+tX7/elkJOJlgvkF0nmrsxqPkbZOx6fF/vK9CBK9BB7fgw5p2cdFPV9vELX9cW5bFG+XVL\nltKSpLTW1tcura1w1inxxPUX7fpwCvTvF0HN+aAWTu3khtsFUjiElKYKhzaws16fJ7zt3LmzCgoK\ndOmll2rEiBFKTk5WWlqaLUXYzd/AVXsxdU9jTmjxQ+0JSpmc1P0qKq0JSPMOSLkHrEOVmwuO9Y7t\nOVz/beOirXmvTk85dumeLF1w5ukq3rNZsSeZggEIdYy4BvzjU4+aMUZbt25Vly5dJEnz58/XoUOH\ndMEFFyguLi7wRTYymdrVU+XW/9L8OaxBj5pvtysqlVp17K2Plmyonh3+u6oZ4r87YH092bli8dFW\nr1jtyUl7pFi9YtF1LOPjdC8CPWqBv29/f6/8PRXBjT1qgTwc68b3j7/CoTepqcKhDYJ+6NMYo8zM\nTK1bt86WB20sgtrJ7z+QQcaXx2+KYAU1Y6yZ4Nue1l/vLlyt7YXWCfp/evQZjbr859peKG07ZE1j\n0ZDUlscOUXZLsnrFulf1jnVMPHZivb9tEOzbufH96u9928kNQS2Qvx92CeQpG3bXYNftAikcQkpT\nhUMbBP3Qp8fj0aBBg7Rs2TINGTLElgcG7FJpJCWkavl2K2x5Lxr7gn74b+v77YeqZpH/xWpd9PJx\nNx74c3245dhmfLRUsmeThg/qoU6JUtfW1qVL1aVrkjXVAwAAweDzYIJevXpp8+bN6tq1q1q2bGnd\n2OPRmjVrAlqg93HoUav//sO5R62iUopp3V45m3ZpR6E1XUX+QemBf87SsNHXKP+g71NXtI6XDm5d\np5HZGerYyjoUed/tN+qdl55Sx1ZS59bWiMmoqOD+50aPmjt/1+hR8w09avYKh96kpgqHNgj6oU9J\nysvLq/PB02wYUDBjxgw988wzMsbo+uuv169//euaRRLUTnr/oRrUDhZLSV0z9P6idVavV6F07yNP\n68KfXK/thdakrvuOVPWYNeTIPmWddoo6JVqjJzslSr+bfK0+mfOf6lDWMi6wocFfBDV3/q4FO6gF\nexoauxDU7OWGv0lOC4c2cCSoVVZW6sUXX1Rubq7+8Ic/KD8/X7t27WryodB169Zp3LhxWr58uWJj\nY3XBBRfoySefVPfu3Y8VSVA76f27MahVGml3kTU68ruDVV8PHDsxP/+gtaSRTw7vUUa3duqQcOwQ\n5B9vmaCP3/i3Orf2hrDAnNsWDAQ1d/6uBTuo+Xs/BDV3vn/85Ya/SU4LhzZwJKjddNNNioqK0rx5\n87Rhwwbt379fI0eO1IoVK5pUwOuvv665c+fqmWeekSTdf//9io+P1+23336sSILaSe/fiaBWWWm0\n/6g1RcW3BdY0FbkHjs2q/91Ba73Ik2keIx3dtUHnD+mtTolSx1bSn6f8Qm+/8IQ6trImdm3bQoqL\nCe4hIYLqIyMqAAAWqElEQVSavQhqBDW31GDX7QLJDX+TnBYObWBnvT7Po7Z06VLl5ORowIABkqxu\n+rKysiYXkJGRod/97nfav3+/mjVrpvfee48BCy5RXG71gn1bcCyQfVsg6aYcJT0kHSo5+e1PaWGd\niJ+WVHVSflLNE/NTmktRUX303+PezH9e8aRG93oisE8MCFG151707gMQvnwOanFxcaqoqKje3rt3\nr6Ki6pgUqpF69+6tKVOmaOTIkWrZsqUGDBhgy/3CNzsLawYxXfq8hj1vbW8/pLqXNWqfpUMl1rqR\n3ZKlblXTVFR/n2wFslbxNW9W3zk4AHzDBLGIBPxDUpPPQW3y5MkaM2aM9uzZo7vvvluvv/667r//\nfluKmDhxoiZOnChJuvvuu6sn1j3ehAkTqgcuJCUlKSsrS8OHD5dkTcArqXrbu6++7drX924ff926\nthu6fWMfz99tX5+fJB0ulV57b77UOVszllhhbNni+dpZKOl3h3XqI5Jyq67fbbiUNUGfLbC2o08b\nri6tpeRd89U+QTrnh8PVLUm6cngvvfXCUxo9arg8nuPqO/PY4++ro37vyg92tefxP/Nl2y2vn12P\nZ9fzC/Tz9fXx/H0+gdr2pZ66/vlISEjQO++84/jzq/14gXp8775AvV8D/fvgtm3vvnB9fr5sv/nm\nm66qx5dt7/fegZd28vkcNUn6+uuvNW/ePBljdN5556lPnz62FLFnzx61a9dO+fn5GjVqlJYuXarE\nxMRjRXKOWr33X1YhxbU9TR+v+LZ6SaNvD0ivzF2idr3OOumSRpI1HUW3ZOm0ZOm0JGnalOv18eyn\n1S3ZGj1Z17JGgWyXQJ9rxjlqnKNm533bee6nnThHLbC3CyQ3vH/QdI6co3b06FG9//77WrRokTwe\nj8rKytStWzc1a9asyUX8+Mc/1vfff6/Y2Fg9/vjjNUJapKuotBbzzi04tr6kLnlO58609m09JOnX\n3+r8WbVu2MkKaXHR1jli3yz5QJPGXXgslCVLA9Jaa1/xwRo3m7bqGZ132tPBenoAAOAkfO5Ru/zy\ny5WYmKirr75axhi99NJLOnjwoGbPnh3oGpvco+brmnJO9Kgdv8j38VNZeC8NrS3pkWQObtU5mZ1r\nnC824ZJztHXtZzq1lbWkUbCn5/D3dvSoBf52wV4cmx41599TTreTW2qw63aB5Ib3D5rOkek5+vbt\nq/Xr1ze4LxCaGtR8vU4gglpxec0Fvb/zXqpC2bZDUnkDs+qntrRGTHpD2LQ7b9CHr/5L3ZKs0ZPN\nYoM/PQdBzT7h9CFTF4Ka8++pQD6+v/8I28mN7x9/ueH9g6Zz5NDnwIED9cUXXyg7O1uStGTJEg0a\nNMiWIkKVMdZC394JXL1hTJe/piFPW9/vbuAcMUlqn3BsCgvvV++i311an7i25LSVT2tk938F5DkB\nTqhvlBejHN2vrteo9msJwH8+96j17t1b33zzjTp37iyPx6P8/Hz16tVLMTEx8ngCu+anYz1qsc20\naXdx9fqS3h6x/OMuRxuY1DUmStWLe3dpXXM+sa7eHjGf43I9dQawx8mXx/e3bl+vQ4+avbcLNn97\ndwL5XvQFPWr2PX6gawin3yE3vH/QdI70qM2dO9eWB3QNj0e7ilQjhGnUI7rsNev7rQcl/b5YPf5+\n8rtpHX9sAldvELvzF+P0+fsvq0trqUOCFM20cAAQVME+FxMIlEZNz+EUf3rUDhWb6gC29VCtrwel\nLXtLpJj4k99RRZm6psSqc+uaM+p3aW1NXdGltdS6jkGvgf7vhx41etSacrtgo0fN+feU048f6BpC\n5XfdF6FSJ07OkR615cuX64EHHlBeXp7Ky8urCwnkIc/6lFVI2wuP9YRtrRXIdGeBEqc1cCcx8WrT\nXNUhrHOi9I8/365Xnvpr9b5OSc2UV1nRwB0hlDEDNgDAzXzuUevZs6emT5+ujIyMGks8pVWtFhBI\nxyfTtzZIY1+tZ2mj4zSLscJX16RjvV/er11aS71PbSFTeqTex6lr2596A4EetcA+v2Bz6/vMLvSo\nOd+j5fTjB7oGetTgNo70qLVt21ajR4+25UGbIrWl5PFIpyZUha/WUpfEqq9Vl0Gnn6IjRft00oFH\nZUeDVjMAAIA/fO5R++ijj/Tqq6/q/PPPV1xcnHVjj0djx44NaIHex/GWWVEpVRopLsZKYf7O1+Pv\nf2CNrTcQ6FGjR60ptws2etSc79Fy+vEDXQM9anAbR3rUZs6cqY0bN6q8vLzGoc9gBLXjRUdJ0VJ1\nA9Q+vwihjXPGAAA4xuegtmLFCm3YsIFghIBi+Hz4qR2+Cd6BR5sD4cPnoHb22Wdr/fr1Sk9PD2Q9\nAMIM4Tv4aPMT0VuPUOVzUPviiy+UlZWlbt26KT7emn/Mqek5AABoDMIrQlWjVybgpEbncVgDAIDI\n4PPiRl26dNFnn32mmTNnKi0tTVFRUdqzZ08gawsob9g5/hIqgWf//v0yxlRf+E8RAIDw5HNQmzRp\nkr744gu99NJLkqSEhARNmjQpYIUFWu2wQ+AB3KX2P1MpKSlOlwQAQefzoc+lS5cqJydHAwYMkGQt\neFtWVhawwgBEtrrmRwSASONzj1pcXJwqKo6te7l3794a86kBAADAXj4nrcmTJ2vMmDHas2eP7r77\nbg0dOlR33XVXIGsDAACIaA0uIVVWVqbY2FhJ0tdff61PPvlEknTeeeepT58+ga9QJx9pGsjlNkJ5\naZ9gLyGVkpKigoKC6m1/l/bylS/Pr3ZN9dXltFB+nwWTG5Yl8rcmXivnl5AKFeH0XCKZrZ93DQW1\ngQMHatWqVbY8mL8Iao0X7KAW7HX1nF7L0E6h/D4LJjd8+PtbE68VQc1X4fRcIpmdr1mDhz55cwAA\nADijwVGfe/fu1SOPPFJnYPN4PPrNb34TkMIAAAAiXYNBraKiQoWFhcGoxXXqWxvObec4AXAP1pQE\nYKcGz1EbMGCAcnJyglVPnZw6R62xtQTj8X3FOWqho3btvg6CCOXn7A9+10OHG16rUH1twum5RDI7\nXzOfJ7wFEBwFBQV1/qEGAESeBgcTfPzxx8GoAwAAALU0GNTatGkTjDoAAABQC2tAAQAAuBRBrZG8\nI7qOv6SkpDhdFgAACEMMJmikuqbm4ERvAHBOOE2JEk7PBfYgqCEg+GMDIFjCaW7LcHousAdBDQHB\nHxsgcvGPGmAfghoAwFb8owbYh8EEAAAALkVQAwAAcCmCGgAAgEsR1AAAAFyKoAYAAOBSBDUAAACX\nYnoOG9SeM4j5ggAAgB0IajZgziAAABAIBLUwRS8fAAChj6AWpujlAwAg9DGYAAAAwKVCvkeNxX8B\nAEC4CvmgxiE+AAAQrkI+qCEyMVgCABAJCGpoNDccbo60nlQ3tDkAIPg8xhjjdBEN8Xg8CoEyQ05d\n7UpbB1/tNuc1qFuw36++3DevFYC62Pm3gVGfAAAALkVQAwAAcCmCGgAAgEu5Iqg9+OCDSk9PV2Zm\npq666iqVlJQ4XRIAAIDjHA9qeXl5evrpp7Vq1SqtXbtWFRUVeuWVV5wuCwAAwHGOB7XExETFxsbq\nyJEjKi8v15EjR9SxY0enywIQ4bxTohx/SUlJcbosABHG8aCWkpKi2267TV26dNGpp56qpKQknX/+\n+U6XBcBl6gpOgZxLbv/+/TLG1LgUFBQE7PEAoC6OB7UtW7bo0UcfVV5ennbs2KGioiK9+OKLTpcV\nEYL9wQc0RV3BKdImPgYQeRxfmWDFihU6++yz1aZNG0nS2LFjtXjxYo0fP77G9SZMmKC0tDRJUlJS\nkrKysjR8+HBJ0vz58yWJ7UZuez/k3FJPpG5799W37XR9bLPNNttsn3zb+31eXp7s5vjKBKtXr9b4\n8eO1fPlyNWvWTBMmTNCQIUP0y1/+svo6zP6NcJaSklLjkFpycjI9RS7FKhIAfGHn3wbHe9T69++v\na6+9VoMHD1ZUVJQGDhyoG264wemygKAhlAEA6uN4j5ov+K8VgBvQowbAF6z1CQAAEAEIagAAAC5F\nUAMAAHApghoAAIBLEdQAAABciqAGAADgUgQ1AAAAlyKoAQAAuBRBDQAAwKUIagAAAC7l+FqfABAq\nkpOT5fF4amwDQCCx1icAAICNWOsTAAAgAhDUAAAAXIqgBgAA4FIENQAAAJciqAEAALgUQQ0AAMCl\nCGoAAAAuRVADAABwKYIaAACASxHUAAAAXIqgBgAA4FIENQAAAJciqAEAALgUQQ0AAMClCGoAAAAu\nRVADAABwKYIaAACASxHUAAAAXIqgBgAA4FIENQAAAJciqAEAALgUQQ0AAMClCGoAAAAuRVADAABw\nKYIaAACASxHUAAAAXIqgBgAA4FIENQAAAJciqAEAALgUQQ0AAMClCGoAAAAuRVADAABwKYIaAACA\nSxHUAAAAXIqgBgAA4FIENQAAAJciqAEAALgUQQ0AAMClCGoAAAAuRVADAABwKYIaAACASzke1DZu\n3KgBAwZUX1q3bq3HHnvM6bIAAAAc53hQ69Wrl3JycpSTk6OVK1eqRYsWGjNmjNNlRbz58+c7XULE\noc2DjzYPPto8+Gjz0OZ4UDvexx9/rO7du6tz585OlxLx+MUOPto8+Gjz4KPNg482D22uCmqvvPKK\nrrrqKqfLAAAAcAXXBLXS0lK98847uvzyy50uBQAAwBU8xhjjdBGS9Pbbb+uJJ57Q3LlzT/jZ6aef\nri1btjhQFQAAQON0795dmzdvtuW+Ymy5Fxu8/PLLGjduXJ0/s+vJAgAAhBJX9KgdPnxYXbt2VW5u\nrlq1auV0OQAAAK7giqAGAACAEzkymGDixIlKTU1VZmZm9b7Vq1crOztb/fr10+jRo1VYWFj9swcf\nfFA9evRQ79699dFHH1XvX7lypTIzM9WjRw/9+te/DupzCDWNafP//ve/Gjx4sPr166fBgwfr008/\nrb4Nbe67xr7PJSk/P18JCQl6+OGHq/fR5r5rbJuvWbNG2dnZysjIUL9+/VRaWiqJNm+MxrR5cXGx\nxo0bp379+qlv376aNm1a9W1oc99t3bpV5557rtLT05WRkVE9Sfz+/fs1YsQI9ezZUyNHjtSBAweq\nb8PnaNM0ts1t/Rw1Dli4cKFZtWqVycjIqN43ePBgs3DhQmOMMc8995y55557jDHGfPXVV6Z///6m\ntLTU5Obmmu7du5vKykpjjDFnnHGGWbp0qTHGmAsvvNB88MEHQX4moaMxbZ6Tk2N27txpjDFm3bp1\npmPHjtW3oc1915g297rsssvMFVdcYaZPn169jzb3XWPavKyszPTr18+sWbPGGGPM/v37TUVFhTGG\nNm+MxrT5888/b6688kpjjDFHjhwxaWlp5rvvvjPG0OaNsXPnTpOTk2OMMaawsND07NnTrF+/3tx+\n++3moYceMsYYM23aNDNlyhRjDJ+jdmhsm9v5OepIUDPGmNzc3Bq/2K1bt67+Pj8/3/Tt29cYY8wD\nDzxgpk2bVv2zUaNGmS+++MLs2LHD9O7du3r/yy+/bG688cYgVB66fG3z41VWVpqUlBRTWlpKm/uh\nMW0+Z84cc/vtt5upU6dWBzXavPF8bfP33nvPXH311SfcnjZvPF/bfO7cuebiiy825eXlZu/evaZn\nz56moKCANm+iSy65xPz3v/81vXr1Mrt27TLGWMGiV69exhg+RwOhoTY/XlM/R10zj1p6errefvtt\nSdLs2bO1detWSdKOHTvUqVOn6ut16tRJ27dvP2F/x44dtX379uAWHeLqa/PjvfHGGxo0aJBiY2O1\nfft22ryJ6mvzoqIi/eUvf9HUqVNrXJ82b7r62vybb76Rx+PRBRdcoEGDBumvf/2rJNrcDvW1+ahR\no5SYmKgOHTooLS1Nt99+u5KSkmjzJsjLy1NOTo7OPPNM7d69W6mpqZKk1NRU7d69WxKfo3bzpc2P\n19TPUdcEteeee06PP/64Bg8erKKiIsXFxTldUthrqM2/+uor3XnnnXrqqaccqjD81NfmU6dO1a23\n3qoWLVrIML7HVvW1eXl5uRYtWqSXXnpJixYt0pw5czRv3jx5PB6HKw599bX5Cy+8oKNHj2rnzp3K\nzc3V9OnTlZub63C1oauoqEiXXXaZZsyYccKMCR6Ph/dyADS2ze34HHXNPGq9evXShx9+KMn6T/e9\n996TZKXN43t6tm3bpk6dOqljx47atm1bjf0dO3YMbtEhrr42l6z2HDt2rGbNmqVu3bpJEm1ug9pt\n/v7770uSli1bpjfeeEN33HGHDhw4oKioKDVv3lxjx46lzZuovvd5586dNWzYMKWkpEiS/ud//ker\nVq3S1VdfTZs3UX3v88WLF2vMmDGKjo5W27ZtNXToUK1cuVI/+MEPaPNGKisr02WXXaZrrrlGl156\nqSSrR2fXrl1q3769du7cqXbt2knic9QujWlzyb7PUdf0qO3du1eSVFlZqfvvv1+/+MUvJEmjR4/W\nK6+8otLSUuXm5mrTpk0aMmSI2rdvr8TERC1dulTGGM2aNau64eCb+tr8wIED+t///V899NBDys7O\nrr5+hw4daPMmqt3mN910kyRp4cKFys3NVW5urm655Rb97ne/06RJk3if26C+9/moUaO0du1aHT16\nVOXl5VqwYIHS09NpcxvU9z7v3bu35s2bJ8maP3PJkiXq3bs3bd5Ixhj97Gc/U9++fXXLLbdU7x89\nerRmzpwpSZo5c2Z1G/I52nSNbXNbP0ftOKmusa688krToUMHExsbazp16mSeffZZM2PGDNOzZ0/T\ns2dPc9ddd9W4/p///GfTvXt306tXLzN37tzq/StWrDAZGRmme/fuZvLkycF+GiGlMW1+3333mZYt\nW5qsrKzqy969e40xtHljNPZ97jV16lTz8MMPV2/T5r5rbJu/8MILJj093WRkZFSP1jKGNm+MxrR5\ncXGxGT9+vMnIyDB9+/atMbqZNvfdZ599Zjwej+nfv3/13+gPPvjAfP/99+a8884zPXr0MCNGjDAF\nBQXVt+FztGka2+Z2fo4y4S0AAIBLuebQJwAAAGoiqAEAALgUQQ0AAMClCGoAAAAuRVADAABwKYIa\nAACASxHUAIQdY4zOOecczZ07t3rf7NmzdeGFFzpYFQA0HvOoAQhLX331lS6//HLl5OSorKxMAwcO\n1Icffli9lEtjlJeXKybGNSvuAYggBDUAYWvKlClq0aKFDh8+rISEBH333Xdat26dysrKNHXqVI0e\nPVp5eXm69tprdfjwYUnSP/7xD2VnZ2v+/Pm65557lJKSog0bNmjjxo0OPxsAkYigBiBsHTlyRAMH\nDlRcXJwuuugipaena/z48Tpw4IDOPPNM5eTkyOPxKCoqSvHx8dq0aZOuuuoqLV++XPPnz9dFF12k\nr776Sl27dnX6qQCIUPTlAwhbLVq00E9+8hMlJCTotdde0zvvvKPp06dLkkpKSrR161a1b99eN998\ns1avXq3o6Ght2rSp+vZDhgwhpAFwFEENQFiLiopSVFSUjDF688031aNHjxo/nzp1qjp06KBZs2ap\noqJCzZo1q/5Zy5Ytg10uANTAqE8AEWHUqFF67LHHqrdzcnIkSYcOHVL79u0lSf/5z39UUVHhSH0A\nUBeCGoCw5/F4dM8996isrEz9+vVTRkaG/vjHP0qSJk2apJkzZyorK0sbN25UQkJCjdsBgJMYTAAA\nAOBS9KgBAAC4FEENAADApQhqAAAALkVQAwAAcCmCGgAAgEsR1AAAAFyKoAYAAOBSBDUAAACX+n/T\n76pMACeloAAAAABJRU5ErkJggg==\n",
       "text": [
        "<matplotlib.figure.Figure at 0x3f8aa58>"
       ]
      }
     ],
     "prompt_number": 5
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": []
    }
   ],
   "metadata": {}
  }
 ]
}