optimizing_designs.ipynb

{
 "metadata": {
  "name": "",
  "signature": "sha256:733a2b40aaaf8752d00addf0ef59507fb53370b4f7d58e004e507b30750127ad"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "# Optimizing experimental designs for power"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from __future__ import print_function\n",
      "from __future__ import division"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 1
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import numpy as np\n",
      "import numpy.linalg as npl\n",
      "import matplotlib.pyplot as plt"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 2
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%matplotlib inline"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 3
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Remember that the students-t statistic is:\n",
      "\n",
      "$$\n",
      "t = \\frac{c^T \\hat\\beta}{\\sqrt{\\mathrm{var}(c^T \\hat\\beta)}}\n",
      "$$\n",
      "\n",
      "which works out to:\n",
      "\n",
      "$$\n",
      "t = \\frac{c^T \\hat\\beta}{\\sqrt{\\hat{\\sigma}^2 c^T (X^T X)^+ c}}\n",
      "$$\n",
      "\n",
      "where $\\hat{\\sigma}^2$ is our estimate of variance in the residuals, and $(X^T X)^+$ is the [pseudo-inverse](https://en.wikipedia.org/wiki/Moore%E2%80%93Penrose_pseudoinverse) of $X^T X$."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from nipy.modalities.fmri.design import block_design, natural_spline, stack2designs\n",
      "from nipy.modalities.fmri import hrf"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 4
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "times = np.arange(0, 240, 2)  # TR onset times\n",
      "T = len(times)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 5
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def random_onsets(times, p_thresh = 0.5):\n",
      "    is_onset = np.random.uniform(size=times.shape) > p_thresh\n",
      "    return times[is_onset]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 6
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "isis = np.diff(random_onsets(times))\n",
      "print(isis)\n",
      "_ = plt.hist(isis)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[ 6 10  4  2  4 10  2  6  2  6  6  2  2 10  4  6  4  6  2  2  2  2  2  6  2\n",
        "  2  2  2  8  8  2  2  2  2  2  4  2  6  4  2  2  2  4  2  2  6  2  2  2  6\n",
        "  2  4  2  2  6  2  4  2  4  2  2  4  4  2]\n"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAXIAAAD9CAYAAAChtfywAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAGEJJREFUeJzt3X9QVPe9//H3Qck1BBVJYSFCCpOGyC/ZNUQyNtRVXDQ3\nkS/W1K+k6o5iptPetDGmucFOneAkY9bmmmhs7h0nJelqcrX5ppGSH1IxuAa10WvYvRJ/hDaFEcyy\n/kBUBAeB8/3DkFECLLuyLJ/k+Zj5zOyu53zOa3bH1+x89hyOpuu6AADUFRLsAACAm0ORA4DiKHIA\nUBxFDgCKo8gBQHEUOQAoblBF3tXVNcpkMjnnzp37nohIc3NzpMViqUhKSqrNzc3d1dLSEhHYmACA\n/gyqyDdu3PhESkrKMU3TdBERm81WZLFYKmpra5NycnI+stlsRYGNCQDoj9cib2xsjPvwww//dfny\n5X/QdV0TESkrK8uzWq12ERGr1WovLS3ND3RQAEDfvBb5k08++fKLL774dEhISHfPax6Px2AwGDwi\nIgaDwePxeAyBDAkA6N/ogf7x/ffffzg6Ovq0yWRyOhwOc1/baJqm9yy59PVvQ5ARAL5zelZABmPA\nb+QHDhyYVlZWlpeYmFhXUFCwrbKycubixYu3GgwGT1NTU4yIiNvtjo2Ojj49QJgRP5599tmgZ/g2\nZCQnOUf6UCWnrwYs8rVr1/6moaEhvq6uLnH79u0LZ86cWbl169bFeXl5ZXa73SoiYrfbrfn5+aU+\nHxkAMCR8Oo+8Z6mkqKjIVlFRYUlKSqqtrKycWVRUZAtMPACANwOukV9v+vTpe6dPn75XRCQyMrJ5\n9+7dswIXa3iZzeZgR/BKhYwi5Bxq5BxaquT0lebPesygJ9c0PZDzA8C3kaZpog/Vj50AgJGPIgcA\nxVHkAKA4ihwAFEeRA4DiKHIAUBxFDgCKo8gBQHEUOQAojiIHAMVR5ACgOIocABRHkQOA4ihyAFAc\nRQ4Aihv0jSX8tXz544E+xICys+8Xq3VRUDMAQCAF/MYSIpsCNr93n4vJdFSqqyuDmAEAfOPrjSUC\n/o1cJJjfyCtF5GgQjw8Aged1jfzKlStjsrKyDhqNRldKSsqxVatWvSAiUlxcXBwXF9doMpmcJpPJ\nWV5ePifwcQEAvXn9Rj5mzJgre/bsmREWFtbW2dk5+oEHHti3b9++BzRN01euXPnSypUrXxqOoACA\nvg3qrJWwsLA2EZGOjo5burq6Rk2YMOG8iPi0hgMACIxBFXl3d3eI0Wh0GQwGz4wZM/akpqYeFRHZ\ntGnTLzMyMv63sLCwpKWlJSKwUQEAffHprJULFy6Mnz179l9tNltRSkrKsaioqDMiIqtXr37O7XbH\nlpSUFN4wuabpIs9e94r5qzFcKsVkep6zVgCMaA6HQxwOx9fP16xZ49OKh8+nHz733HOrb7311vZf\n//rX/9HzWn19fcLcuXPfq6mpSb9hck3TRQJ3eqN3FDkA9fh6+qHXpZWzZ89+r2fZpL29/daKigqL\nyWRyNjU1xfRss2PHjnnp6ek1/kUGANwMr2etuN3uWKvVau/u7g7p7u4OWbx48dacnJyPlixZssXl\nchk1TdMTExPrNm/e/LPhCAwAuNEwXNnJ0goA+GLIl1YAACMbRQ4AiqPIAUBxFDkAKI4iBwDFUeQA\noDiKHAAUR5EDgOIocgBQHEUOAIqjyAFAcRQ5ACiOIgcAxVHkAKA4ihwAFEeRA4DiKHIAUBxFDgCK\no8gBQHEDFvmVK1fGZGVlHTQaja6UlJRjq1atekFEpLm5OdJisVQkJSXV5ubm7mppaYkYnrgAgN4G\nLPIxY8Zc2bNnzwyXy2U8cuTI5D179szYt2/fAzabrchisVTU1tYm5eTkfGSz2YqGKzAA4EZel1bC\nwsLaREQ6Ojpu6erqGjVhwoTzZWVleVar1S4iYrVa7aWlpfmBDgoA6Ntobxt0d3eHTJkypfqLL764\n6+c///l/paamHvV4PAaDweARETEYDB6Px2Pof4bi6x6bvxoAgB4Oh0McDoff+3st8pCQkG6Xy2W8\ncOHC+NmzZ/91z549M67/d03TdE3T9P5nKPY7HAB8F5jNZjGbzV8/X7NmjU/7D/qslfHjx1946KGH\nPvj000/vNRgMnqamphgREbfbHRsdHX3ap6MCAIbMgEV+9uzZ7/WckdLe3n5rRUWFxWQyOfPy8srs\ndrtVRMRut1vz8/NLhyMsAOCbBlxacbvdsVar1d7d3R3S3d0dsnjx4q05OTkfmUwm54IFC94uKSkp\nTEhIqH/77bcXDFdgAMCNNF0fYHn7ZifXNF0kcPN7Vykm0/NSXV0ZxAwA4BtN00TXdW2w23NlJwAo\njiIHAMVR5ACgOIocABRHkQOA4ihyAFAcRQ4AiqPIAUBxFDkAKI4iBwDFUeQAoDiKHAAUR5EDgOIo\ncgBQHEUOAIqjyAFAcRQ5ACiOIgcAxVHkAKA4r0Xe0NAQP2PGjD2pqalH09LSPnvllVd+JSJSXFxc\nHBcX12gymZwmk8lZXl4+J/BxAQC9jfa2QWho6NWXX375SaPR6GptbQ2/9957P7VYLBWapukrV658\naeXKlS8NR1AAQN+8FnlMTExTTExMk4hIeHh4a3Jy8vFTp05NFBGf7vIMAAgMr0V+vfr6+gSn02m6\n//77P9m/f/8PN23a9MstW7YsyczMPLx+/fqnIiIiWr65V/F1j81fDQBAD4fDIQ6Hw+/9NV3XB7Vh\na2truNlsdvz2t799Pj8/v/T06dPRUVFRZ0REVq9e/Zzb7Y4tKSkpvGFyTdNFBjd/YFSKyfS8VFdX\nBjEDAPhG0zSfVjwGddbK1atXQ+fPn//nRYsWvZmfn18qIhIdHX1a0zRd0zR9+fLlfzh06NBUf0MD\nAPzntch1XdcKCwtLUlJSjq1YsWJDz+tutzu25/GOHTvmpaen1wQqJACgf17XyPfv3//DN998c9Hk\nyZOPmEwmp4jI2rVrf7Nt27YCl8tl1DRNT0xMrNu8efPPAh8XANDboNfI/ZqcNXIA8FlA1sgBACMX\nRQ4AiqPIAUBxFDkAKI4iBwDFUeQAoDiKHAAUR5EDgOIocgBQHEUOAIqjyAFAcRQ5ACiOIgcAxVHk\nAKA4ihwAFEeRA4DiKHIAUBxFDgCKo8gBQHFei7yhoSF+xowZe1JTU4+mpaV99sorr/xKRKS5uTnS\nYrFUJCUl1ebm5u5qaWmJCHxcAEBvXos8NDT06ssvv/zk0aNHUz/55JP7X3311X87fvx4ss1mK7JY\nLBW1tbVJOTk5H9lstqLhCAwAuJHXIo+JiWkyGo0uEZHw8PDW5OTk46dOnZpYVlaWZ7Va7SIiVqvV\nXlpamh/osACAbxrty8b19fUJTqfTlJWVddDj8RgMBoNHRMRgMHg8Ho+h772Kr3ts/moAAHo4HA5x\nOBx+7z/oIm9tbQ2fP3/+nzdu3PjE2LFjL13/b5qm6Zqm6X3vWex3OAD4LjCbzWI2m79+vmbNGp/2\nH9RZK1evXg2dP3/+nxcvXrw1Pz+/VOTat/CmpqYYERG32x0bHR192qcjAwCGhNci13VdKywsLElJ\nSTm2YsWKDT2v5+XlldntdquIiN1ut/YUPABgeGm63s+KyFf27dv3wI9+9KOPJ0+efKRn+eSFF15Y\nNXXq1EMLFix4++TJk3cmJCTUv/322wsiIiJabphc03SRgecPrEoxmZ6X6urKIGYAAN9omia6rmuD\n3t5bkd9kGIocAHzka5FzZScAKI4iBwDFUeQAoDiKHAAUR5EDgOIocgBQHEUOAIqjyAFAcRQ5ACiO\nIgcAxVHkAKA4ihwAFEeRA4DiKHIAUBxFDgCKo8gBQHEUOQAojiIHAMVR5ACgOK9FvmzZstcNBoMn\nPT29pue14uLi4ri4uEaTyeQ0mUzO8vLyOYGNCQDoj9ciX7p06Ru9i1rTNH3lypUvOZ1Ok9PpNM2Z\nM6c8cBEBAAPxWuTZ2dlVEyZMON/7dV/u8AwACJzR/u64adOmX27ZsmVJZmbm4fXr1z8VERHR0veW\nxdc9Nn81AAA9HA6HOBwOv/fXdF33ulF9fX3C3Llz36upqUkXETl9+nR0VFTUGRGR1atXP+d2u2NL\nSkoKvzG5puki3ucPnEoxmZ6X6urKIGYAAN9omubTqodfZ61ER0ef1jRN1zRNX758+R8OHTo01Z95\nAAA3z68id7vdsT2Pd+zYMe/6M1oAAMPL6xp5QUHBtr17904/e/bs9+Lj4xvWrFnzrMPhMLtcLqOm\naXpiYmLd5s2bfzYcYQEA3zSoNXK/J2eNHAB8Nixr5ACAkYMiBwDFUeQAoDiKHAAUR5EDgOIocgBQ\nHEUOAIqjyAFAcRQ5ACiOIgcAxVHkAKA4v28sAd+MGxcply5940ZLw2rs2Aly8WJzUDMAGHoU+TC5\nVuLB/ANiIpcucXc+4NuIpRUAUBxFDgCKo8gBQHGskeM7aST8+CzCD9AYGhQ5vpNGwo/P13LwAzRu\nHksrAKA4r0W+bNmy1w0Ggyc9Pb2m57Xm5uZIi8VSkZSUVJubm7urpaUlIrAxAQD98VrkS5cufaO8\nvHzO9a/ZbLYii8VSUVtbm5STk/ORzWYrClxEAMBAvBZ5dnZ21YQJE274VaisrCzParXaRUSsVqu9\ntLQ0P1ABAQAD8+vHTo/HYzAYDB4REYPB4PF4PIb+ty6+7rH5qwEA6OFwOMThcPi9/02ftaJpmq5p\n2gA//xff7CEA4FvNbDaL2Wz++vmaNWt82t+vs1YMBoOnqakpRkTE7XbHRkdHn/ZnHgDAzfOryPPy\n8srsdrtVRMRut1vz8/NLhzYWAGCwvBZ5QUHBtmnTph34/PPP74mPj2944403lhYVFdkqKiosSUlJ\ntZWVlTOLiopswxEWAPBNmq4H7uq2a2vnwbx6rlJMpueluroyiBmu0TRNgn8loSaB/LxVMjI+DxE+\nE/RF0zTRdX3Ql/1yZScAKI4iBwDFUeQAoDiKHAAUR5EDgOIocgBQHEUOAIqjyAFAcRQ5ACiOIgcA\nxVHkAKA4ihwAFEeRA4DiKHIAUBxFDgCKo8gBQHEUOQAojiIHAMWNvpmdExIS6seNG3dx1KhRXaGh\noVcPHTo0daiCAQAG56aKXNM03eFwmCMjI5uHKhAAwDc3vbTiyw1CAQBD76aKXNM0fdasWbszMzMP\nv/baa48NVSgAwODd1NLK/v37fxgbG+s+c+ZMlMViqZg0adKJ7Ozsqhu3Kr7usfmrAQDo4XA4xOFw\n+L2/puv6kARZs2bNs+Hh4a1PPfXU+q8n1zRdZGjm90+lmEzPS3V1ZRAzXKNpmgT3vRAR0WSoPm/V\njYzPQ4TPBH3RNM2nZWu/l1ba2trCLl26NFZE5PLly7ft2rUrNz09vcbf+QAA/vF7acXj8RjmzZu3\nQ0Sks7Nz9E9/+tO3cnNzdw1dNADAYPhd5ImJiXUul8s4lGEAAL7jyk4AUBxFDgCKo8gBQHEUOQAo\njiIHAMVR5ACguJu6RB8Avk3GjYuUS5fOBzuGzyhyAPjKtRIfCX8ywbc/KsvSCgAojiIHAMVR5ACg\nOIocABTHj53Ad9xIOVNj7NgJcvEit//1B0UOfMeNlDM1Ll3i9r/+YmkFABRHkQOA4ihyAFAcRQ4A\niqPIAUBxN1Xk5eXlcyZNmnTi7rvv/vu6deueGapQw83hcAQ7wiA4gh1gUNR4L9WhzvvpCHaAQVHn\n/fSN30Xe1dU16vHHH/99eXn5nGPHjqVs27at4Pjx48lDGW64qPHhOoIdYFDUeC/Voc776Qh2gEFR\n5/30jd9FfujQoak/+MEP/pGQkFAfGhp6deHChdv/8pe//J+hDAcA8M7vC4JOnTo1MT4+vqHneVxc\nXOPBgwezem83fvxcfw9x0zo7z0po6G1BOz4ADAtd1/0a77zzzvzly5e/1vN869atix5//PFN128j\n1y4XYzAYDIaPw5c+9vsb+cSJE081NDTE9zxvaGiIj4uLa7x+G13XueYWAALM7zXyzMzMw3//+9/v\nrq+vT+jo6LjlT3/60//Ny8srG8pwAADv/P5GPnr06M7f//73j8+ePfuvXV1dowoLC0uSk5OPD2U4\nAMAg+LtG3t84efJkvNls3pOSknI0NTX1s40bN/5qqI8xFKO9vX3M1KlTD2ZkZLiSk5OPFRUVvRDs\nTAONzs7OUUaj0fnwww+/F+ws/Y3vf//79enp6UeMRqPzvvvuOxTsPP2N8+fPR8yfP/+dSZMmHU9O\nTj72t7/97f5gZ+o9Tpw4cY/RaHT2jHHjxl0Yif+X1q5duyolJeVoWlpaTUFBwX9fuXLlX4Kdqa+x\nYcOGJ9LS0mpSU1M/27BhwxPBztMzli5d+np0dLQnLS2tpue1c+fORc6aNavi7rvvrrVYLLvOnz8f\n4W2eIQ/mdrtjnE6nUdd1uXTpUnhSUtLnx44dSw72G9bXuHz5cpiu63L16tXRWVlZn1RVVT0Q7Ez9\njfXr16989NFH35o7d25ZsLP0NxISEurOnTsXGewc3saSJUvsJSUly3o++5aWlvHBzjTQ6OrqComJ\niXGfPHkyPthZrh91dXUJiYmJ/+wp7wULFvzpj3/8ozXYuXqPmpqatLS0tJr29vYxnZ2do2bNmlXx\nj3/8465g59J1XT7++OPs6upq0/VF/vTTT/9u3bp1/67ruthstmeeeeYZm7d5hvwS/ZiYmCaj0egS\nEQkPD29NTk4+/uWXX94x1McZCmFhYW0iIh0dHbd0dXWNioyMHJF/1b6xsTHuww8//Nfly5f/QR/h\nPyCP9HwXLlwYX1VVlb1s2bLXRa4tEY4fP/5CsHMNZPfu3bPuuuuuL64/3XckGDdu3MXQ0NCrbW1t\nYZ2dnaPb2trCJk6ceCrYuXo7ceLEpKysrINjxoy5MmrUqK7p06fvfffdd38c7FwiItnZ2VUTJky4\n4a4eZWVleVar1S4iYrVa7aWlpfne5gno31qpr69PcDqdpqysrIOBPI6/uru7Q4xGo8tgMHhmzJix\nJyUl5ViwM/XlySeffPnFF198OiQkpDvYWQaiaZo+a9as3ZmZmYdfe+21x4Kdpy91dXWJUVFRZ5Yu\nXfrGlClTqh977LHX2trawoKdayDbt29f+Oijj/53sHP0FhkZ2fzUU0+tv/POO0/ecccdX0ZERLTM\nmjVrd7Bz9ZaWlvZZVVVVdnNzc2RbW1vYBx988FBjY2NcsHP1x+PxGAwGg0dExGAweDwej8HbPgEr\n8tbW1vBHHnnknY0bNz4RHh7eGqjj3IyQkJBul8tlbGxsjPv4449/5HA4zMHO1Nv777//cHR09GmT\nyeQc6d929+/f/0On02nauXPng6+++uq/VVVVZQc7U2+dnZ2jq6urp/ziF7/4z+rq6im33XbbZZvN\nVhTsXP3p6Oi45b333pv7k5/85P8FO0tvX3zxxV0bNmxYUV9fn/Dll1/e0draGv7WW2/9NNi5eps0\nadKJZ555Zl1ubu6uBx98cKfJZHKO9C9FPTRN0zVN071tF5Aiv3r1auj8+fP/vGjRojfz8/NLA3GM\noTR+/PgLDz300AeHDx/ODHaW3g4cODCtrKwsLzExsa6goGBbZWXlzCVLlmwJdq6+xMbGukVEoqKi\nzsybN2/HoUOHpgY7U29xcXGNcXFxjffdd9//iIg88sgj71RXV08Jdq7+7Ny588F7773306ioqDPB\nztLb4cOHM6dNm3bg9ttvPzd69OjOH//4x+8eOHBgWrBz9WXZsmWvHz58OHPv3r3TIyIiWu65557P\ng52pPwaDwdPU1BQjIuJ2u2Ojo6NPe9tnyItc13WtsLCwJCUl5diKFSs2DPX8Q+Xs2bPfa2lpiRAR\naW9vv7WiosJiMpmcwc7V29q1a3/T0NAQX1dXl7h9+/aFM2fOrNyyZcuSYOfqra2tLezSpUtjRUQu\nX758265du3LT09Nrgp2rt5iYmKb4+PiG2traJJFr68+pqalHg52rP9u2bSsoKCjYFuwcfZk0adKJ\nTz755P729vZbdV3Xdu/ePWukLk+ePn06WkTk5MmTd+7YsWPeSFyq6pGXl1dmt9utIiJ2u906qC/D\nQ/0rbFVV1QOapnVnZGS4ek6d2rlz55xg/zrcexw5ciTdZDJVZ2RkuNLT04/87ne/ezrYmbwNh8Mx\nfaSetfLPf/4zMSMjw5WRkeFKTU39bO3atauCnam/4XK5MjIzM/9n8uTJ/ztv3rx3R+pZK62trbfd\nfvvtZy9evDg22Fn6G+vWrfv3ntMPlyxZYu/o6AgNdqa+RnZ29scpKSlHMzIyXJWVlTOCnadnLFy4\ncFtsbOyXoaGhHXFxcQ2vv/760nPnzkXm5OTs9uX0Q03XvS6/AABGMO4QBACKo8gBQHEUOQAojiIH\nAMVR5ACgOIocABT3/wHnTKUsDdRYXwAAAABJRU5ErkJggg==\n",
       "text": [
        "<matplotlib.figure.Figure at 0x512ad10>"
       ]
      }
     ],
     "prompt_number": 7
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "onsets = random_onsets(times)\n",
      "n_onsets = len(onsets)\n",
      "offsets = onsets + 1.\n",
      "block_spec = np.ones(n_onsets, dtype=([('start', float), ('end', float), ('type', int)]))\n",
      "block_spec['start'] = onsets\n",
      "block_spec['end'] = offsets\n",
      "print(block_spec)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[(4.0, 5.0, 1) (6.0, 7.0, 1) (8.0, 9.0, 1) (10.0, 11.0, 1) (14.0, 15.0, 1)\n",
        " (16.0, 17.0, 1) (18.0, 19.0, 1) (20.0, 21.0, 1) (22.0, 23.0, 1)\n",
        " (24.0, 25.0, 1) (26.0, 27.0, 1) (28.0, 29.0, 1) (30.0, 31.0, 1)\n",
        " (32.0, 33.0, 1) (40.0, 41.0, 1) (48.0, 49.0, 1) (50.0, 51.0, 1)\n",
        " (60.0, 61.0, 1) (62.0, 63.0, 1) (64.0, 65.0, 1) (66.0, 67.0, 1)\n",
        " (68.0, 69.0, 1) (74.0, 75.0, 1) (76.0, 77.0, 1) (78.0, 79.0, 1)\n",
        " (82.0, 83.0, 1) (86.0, 87.0, 1) (90.0, 91.0, 1) (102.0, 103.0, 1)\n",
        " (106.0, 107.0, 1) (114.0, 115.0, 1) (118.0, 119.0, 1) (120.0, 121.0, 1)\n",
        " (122.0, 123.0, 1) (124.0, 125.0, 1) (126.0, 127.0, 1) (130.0, 131.0, 1)\n",
        " (132.0, 133.0, 1) (148.0, 149.0, 1) (158.0, 159.0, 1) (160.0, 161.0, 1)\n",
        " (166.0, 167.0, 1) (168.0, 169.0, 1) (172.0, 173.0, 1) (182.0, 183.0, 1)\n",
        " (184.0, 185.0, 1) (190.0, 191.0, 1) (192.0, 193.0, 1) (196.0, 197.0, 1)\n",
        " (204.0, 205.0, 1) (208.0, 209.0, 1) (210.0, 211.0, 1) (214.0, 215.0, 1)\n",
        " (216.0, 217.0, 1) (218.0, 219.0, 1) (222.0, 223.0, 1) (226.0, 227.0, 1)\n",
        " (230.0, 231.0, 1) (234.0, 235.0, 1) (236.0, 237.0, 1) (238.0, 239.0, 1)]\n"
       ]
      }
     ],
     "prompt_number": 8
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def random_design(times, duration=1):\n",
      "    onsets = random_onsets(times)\n",
      "    offsets = onsets + duration\n",
      "    n_onsets = len(onsets)\n",
      "    block_spec = np.ones(n_onsets, dtype=([('start', float), ('end', float), ('type', int)]))\n",
      "    block_spec['start'] = onsets\n",
      "    block_spec['end'] = offsets\n",
      "    return block_design(block_spec, times, hrfs=[hrf.spm])"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 9
    },
    {
     "cell_type": "heading",
     "level": 3,
     "metadata": {},
     "source": [
      "An example of output of random_design"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "X, contrasts = random_design(times)\n",
      "plt.plot(times, X)\n",
      "print(contrasts)\n",
      "print(X.shape)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "{'constant_0': array(1.0000000000000007)}\n",
        "(120,)\n"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD9CAYAAAC1DKAUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXtcVFeW71fx8v3CB5oqlERQwAeiqFFjmm5jlHQH7Y7T\nktvxOmmTMJlxOrm3u5P0TM+N5nZntB/Tkxv7k6E7j04nHXXy6GCmFR8xdEii4gM1CkZAKhYIGEQI\niAoU5/6xZqcOxXk/6uyqWt/Px49U1eHU5tSp/dvrt/Ze2yUIAhAEQRDRR4zTDSAIgiCcgQSAIAgi\nSiEBIAiCiFJIAAiCIKIUEgCCIIgohQSAIAgiSjEtACUlJSvT09PPpaWlVW/duvVJqWNKS0tzs7Oz\nK2bOnHkmNze31Ox7EgRBEOZxmVkH4Pf7Y6dPn/7ZgQMH7nK73Q3z588/un379vszMjKq2DFtbW2j\nlyxZ8vHevXtXeDye+paWlnHjxo1rsaT1BEEQhGFMRQDl5eULUlNTa1JSUrzx8fE9BQUFO4qLi1eJ\nj3njjTf+x3333fe2x+OpBwCgzp8gCIIPTAlAQ0ODOzk52cceezye+oaGBrf4mOrq6rTW1tbEr3/9\n6x/k5OQce+2119aZeU+CIAjCGuLM/LLL5VL1j3p6euJPnDgx9/3331/W1dU1dNGiRYduv/32w2lp\nadV6z0UQBEH0RxAEl9HfNRUBuN3uBp/Pl8we+3y+ZGb1MJKTk3133333viFDhlwfO3bslTvvvPPD\nU6dOZUmdTxCEqP/39NNPO94GXv7RtaBrQddC+Z9ZTAlATk7Oserq6jSv15vS3d2dsHPnzrX5+fm7\nxMesWrWq+KOPPrrD7/fHdnV1DT1y5MjCzMzMSnPNJgiCIMxiygKKi4vr3bZt28YVK1bs9fv9sRs2\nbHgpIyOjqqioqBAAoLCwsCg9Pf3cypUrS2bPnn06Jiam7+GHH/49CQBBEITzmJoGaiUul0vgpS1O\nUlpaCrm5uU43gwvoWgSgaxGArkUAl8sFgokcAAkAQRBEmGJWAKgUBEEQRJRCAkAQBBGlkAAQBEFE\nKSQABEEQUQoJAEEQRJRCAkAQBBGlkAAQBEFEKSQABEEQUQoJAEEQRJRCAkAQBBGlkAAQBEFEKSQA\nBEEQUQoJAEEQRJRCAkAQBBGlkAAQBEFEKSQABEEQUQoJAEEQRJRCAkA4Rl2d0y0giOiGBIBwhDNn\nABYscLoVBBHdkAAQjnDgAMDVqwC0DTRBOAcJAOEI778P4PcDdHWZP9cHHwDs22f+PAQRbZAAECGn\ntxegrAxg2DCA9nbj5/F6Ae67D+Cb3wT47W8tax5BRA2mBaCkpGRlenr6ubS0tOqtW7c+Gfx6aWlp\n7qhRo9qzs7MrsrOzK372s5/91Ox72kVvL8D+/WRL2M3RowApKQAej3EBuHEDYO5cgOxsgDffBPjy\nS0ubSBBRQZyZX/b7/bEbN27cduDAgbvcbnfD/Pnzj+bn5+/KyMioEh/3ta997a+7du3KN9dU+ykt\nBVixAuBrXwP43e8A0tKcblFkcuAAwLJlAB99ZFwALlwAGD8e4Kc/BTh2jASAIIxgKgIoLy9fkJqa\nWpOSkuKNj4/vKSgo2FFcXLwq+DhBEFxm3idUlJUB/PCHAKtWASxaBPDee063KDJ5/30UgFGjjAtA\nTQ1Aair+PHIkCQBBGMFUBNDQ0OBOTk72sccej6f+yJEjC8XHuFwu4ZNPPlmclZV1yu12N/zqV7/6\nUWZmZqXU+TZt2vTVz7m5uZCbm2umebopKwP40Y8A7rkHE5QHDwLce29ImxDxdHXhiP3OOwFeeYUE\ngCD0UFpaCqWlpZadz5QAuFwuVbd87ty5J3w+X/LQoUO79uzZk7d69ep3z58/P03qWLEAhJrubvSm\nFy/Gx6NHA1RKyhRhho8+ApgzB2D4cPMRQEYG/kwCQEQLwQPjzZs3mzqfKQvI7XY3+Hy+ZPbY5/Ml\nezyeevExI0aM6Bg6dGgXAEBeXt6enp6e+NbW1kQz72sHJ04ATJ2KHT8AdSp2wewfABQAo9e4tjYQ\nAQwZAtDTgyJOEKGipQWj2XDGlADk5OQcq66uTvN6vSnd3d0JO3fuXJufn79LfExzc3MSywGUl5cv\nEATBlZiY2Grmfe2grAxg6dLAYxIAezh0CO0fAOtyAC4Xfl4dHda0kSC08IMfAPzzPzvdCnOYsoDi\n4uJ6t23btnHFihV7/X5/7IYNG17KyMioKioqKgQAKCwsLHrrrbfWvPDCC4/GxcX1Dh06tGvHjh0F\n1jTdWsrKAB54IPB4xAgSADtobga45Rb8edQonM2jl+5ugIYGgClTAs8xwR471pp2EoQSR44AvPUW\n2pnhjEvgZNK7y+USnGpLXx9OKTxzBmDSJHzuzBmAtWsBzp51pEkRy7hxmFuZMAGTwH/9K8Af/qDv\nHNXVACtXog3EmD0b4LXXALKyLG0uQQxAEADuuAMdg7ffxvvRKVwul6lZlhG7EvhnPwN48UVtx1ZW\nAowZE+j8AcgCsoO+PoC2NrzWAMYtoJoazNeIoc9LHkEAeOcdp1sRObz9NsC1awCPP471rMIZUxYQ\nzxw9CjB9urZjy8oCvjSDPGXraW/H8g/x8fjYjAAw/59Bn5c8zc1YMqO7O3DtCWMIAsBTTwEUFaHd\n2NaGz7nCYqXTQCI2AqitBfjiC23HBieAATAH0NFBZSGs5MqV/h691QJAEYA0LM/S1uZsOyKBlhYc\n9S9bhmI6eHB4DzwiUgD6+vCmb2nRdvzRowC3397/udhYnF547Zr17YtWWlsBEkUTgI0KgHgKKIME\nQB4mAOFuV/CA14t1rBhjxoT3dY1IAWhsBLh+XVsE4PcDXLwIcOutA1+jTsVaKAJwBpYsb+Vu8nX4\nQQIQBtTW4iwTLQLQ2Iid0uDBA1+jTsVapARA7/X1+/FLeNtt/Z+nz0oeFgGQAJgnWAASE0kAuKOm\nBou5aRGA4A9UDK0FsJbW1v4CMGgQJs9u3NB+Dp8PxT1YsEkA5KmtBXC7SQCsgCKAMKC2FmvF37gB\ncPOm8rFKAkCdirVcudI/BwCg3waSsn8A6LNS4sIFgJwcEgArIAEIA1iScNw49SiABCB0BFtAACQA\ndtPVhR3/rFnh3VHZzXe/i9V/1ZASgHAW1ogUALZQaPx49ZlAJAChI9gCAsBrTAJgH3V1eH+PHRve\nHZXdVFerr+gVBOwvxCVIKALgDEHoLwBmI4BwnuPLG1ZYQBcuDEwAA5AAyMGuV2IiCYASTU34T4mW\nFsw9jRwZeI4EgDPYTT52rDUCQJ2KdVhhAV2+DJCUNPB5+qykuXABB0MkAPL4/XhfNTYqHyfVV5AA\ncEZtLd7wLpe6APj9OKtk8mTp16lTsRYpC0ivAEiJCAB9VnLU1gYigHDuqOykpQUXj5IARADiVaJq\nAqC0BgCAOhWrscICIgHQB7OAwj1ZaSdNTThg1CIAwQtGtQrAv/4rwO7dhptoGxEnAOJKkWoCoGT/\nANA6ACvp6cGyGqNG9X9ejwAIAn7ZgkUEALeYvHYNR3JEABYRkwUkT1MTQFqavRFASQnAp58abaF9\nRJwAiCMAtWmgagJAo0rruHoVvywxQXecntXA7e0AQ4dKV7SMicFKo52d5tsaKfT1BUatrKOi4oYD\naWrCjV2am5UHEFL9hVZrraqKTwGOOAEIjgCUpoGSAIQOKfsHQF8EIGf/MOjz6k9jI17fYcMAEhLC\nv3KlXTQ14dTOESPwHpNDqr8YPTpQElqOlhYciCqd2ykiTgD05ABIAEKHVAIYgATATpj9w6BEsDRN\nTQATJ+KGUHI2kNQaAACMRocMURbWqir8nyIAm7l2DTsTtrMXCQA/yHXeJAD2EbxmghLB0jQ2qgtA\nSwt29CNGDHxNLQ9QVYX9DEUANlNbi34n85kTE7Fz8fulj9ciABQyW4PS7B0SAHsIFgBKBEujJQJQ\n6ivUhLWqCvcQ5vHaR5QANDZi1UNGbCx6dFLK6/cD1NfLrwEAoA7FSoI3g2FQBGAfUhYQj52Q01gh\nAGoRAAlACGhpwZk/YuRsILU1AABYrrivT72iKKGOFRZQSwsJgB4uXeo/ICIBkCYUArBkCX4HeJuF\nFXECMH58/+fkpoKq2T8AuDiE7Q1sNX6/vjr44Y5VOYBggRdDAtCf4OsV7qtW7eD6dfw3Zow9AtDZ\nif1PRgb2J9evW9Fq6zAtACUlJSvT09PPpaWlVW/duvVJueOOHj06Py4urvedd975jtn3lEMuApCa\nCqpFAADs61TeeAOgsND68/KKnAU0ZAhAby9Ad7f6OcgC0kdwxEQRwECam7G2lMulLgBydrGSAJw7\nh4vMYmPxs+AtEWxKAPx+f+zGjRu3lZSUrKysrMzcvn37/VVVVRlSxz355JNbV65cWSIIgsvMeyrx\nxRfaLSCpKV1S2NWp1NQAfPaZ9eflFbnO2+XSHgWQAGhHEAZeLxKAgTD7B0BZAC5elO8vlKbXVlXh\n6J8dx9v1NyUA5eXlC1JTU2tSUlK88fHxPQUFBTuKi4tXBR/3/PPP/+OaNWveGj9+vIZNGo2jJwfg\ndATg82EbogWlzlvramASAO10dgLExWGExeCxA3IaKQGQ8umVBEApAhALAI8RQJyZX25oaHAnJyf7\n2GOPx1N/5MiRhcHHFBcXrzp48OA3jh49Ot/lcsmmQTZt2vTVz7m5uZCbm6urPXICILXRg9eLuwCp\nYacANDfjjk1Dh1p/ft6Qs4AAKAKwA6nvAuUABiIWgOHD0ar58sv+Nava21EUgutYMdQE4P778Wcr\nBLi0tBRKS0vNnUSEKQFQ6swZjz/++L9v2bLlKZfLJQiC4FKygMQCYASpJPD48QCffDLwWD0WkB1J\nYJ8PR2iffx4YIUQyahEACYC1SCXMKQIYiFgAAAJRgLiz//xz9P9dMj2XmgBkZuLPVuzKFjww3rx5\ns6nzmRIAt9vd4PP5ktljn8+X7PF46sXHHD9+fF5BQcEOAICWlpZxe/bsyYuPj+/Jz8/fZea9pdCa\nA+jrU18DwLCjUxEEfP+cHNyyL9IF4Pp1vOZykY4WAbhxA5PFw4bJH0PVWwNITZklARhIUxNAVlbg\nMROA9PTAcxcvKvcVcgLQ3Y0DzbQ0fJyYyJ8FZCoHkJOTc6y6ujrN6/WmdHd3J+zcuXNtcMd+4cKF\n2+rq6m6tq6u7dc2aNW+98MILj9rR+ff1SdebkZoG2tiIH5rYH5XDDgFoa8NQc/bs6MgDMPtHbgSl\nZTUwG/3LnYOdx8rPKi+Pvy+sVigC0IZcBCBGiwBIXdeaGvy9QYPwMY/X35QAxMXF9W7btm3jihUr\n9mZmZlauXbt2Z0ZGRlVRUVFhUVFRSCc5trfj6DC4VLDUNFCt9g+APaNKnw8gORnLVtTVWXtuHlGz\nbrREAGrnALBWAAQB4P33AY4cseZ8oUYqAhg2DPdloIWNAawSAKkI4Px5gGnTAo8jLgkMAJCXl7cn\nLy9vj/i5wsLCIqljX3nllQfNvp8cUkkvAHyupQW/0Gz0+Pnn2mYAAWCnYnUnzQQgJQXg+HFrz80j\n4SgA165hZ3nsGMA991hzzlAiFQG4XIHOStzpRTNaBUDpHhg9Gu/fvr7++10E7yAWcREAT3zxxcAE\nMACGX0OG9O9gtE4BBbDHAvL5ADye6IkA5EpBM6wSABatWbHcno3ojh41fy4nkCubwWMn5BSCgAKQ\nlBR4Tk4AlBwDuZLQn3/e//esSAJbTcQIgFwEAABwyy3Y6TL0WEB2CQCLAKIhByC3GQxDyzXWIgAJ\nCfhltGK5/dWrOLI7epS/+i1akCubQQIQoL0dB4jiyQlSAsBmASkhZQMFC0DEJYF5QkkAsrP7Wy16\nLSCrBaC+HgVgwgRcBxDpJafZdpByDBuGlosSWgQAwLrPq7UVYNYsDOsbGsyfL9RQBKAO2wdATLAA\n9PTgep1bblE+l9RqYIoAQoiSAMyfj14ugwcLKDkZPdloiAI6OqQ30mCwDd2VCLUAMNGaPz88bSC5\n7wMtBgsQ7P8DoABcuhSI+i5dQotIah9qMVoigDFj+KsIGlECIJUDAMD59uxL3NenntUXY8dCMCYA\nANGRB1ATAC0RgFopaIaVApCYiPeOePAQLpAFpE6w/w+AnfTIkYHqAVr7imAB6OzE6H7ChMBzgwej\nkKjd66EkYgRAahEYY+5cgDNncGHG5cvYGSktKBJjdQTAFoGJBSDaIwBeLaBwjQAEgSwgLcjNlLrr\nLoADB/CxVgFwu/FYhtzqYd6mgkaMAChZQMOG4dZ4Z87os38ArF8H0NKCSSeWeEpJoQhg2DAcMSmh\nRwC07i+gBLOAWATAU9iuRlcXTkeUWnlNAhAguOYPw4gAzJwJ8OmngcfB9g+Dt+sfFQIAELCB9MwA\nAgj40319ppsIAIEpoAyKALRHAEqfL0NLPkELzAKaOBHbd+GC+XOGCiW7jDaGD9DejgOGYJYtA/jg\nA9y0Sa4jD2b2bG0CwFsiOGoEgIXyemYAAWDJhqFD1UeoWhHbPwAUAQBYmwQeOhRHwGZhFhBA+NlA\nSmJplUBGAu3t0hHAxIk4SDt+XF8EUFmJogEg38/wNhU0YgRAbiEYg80E0hsBAFibBxAngAEoAgBQ\njwD8fvyyKk0lFZ/LCgEQT10VTyIIB5QiAC3RVrQgJwAAARtIqwCMGIEJ39pafEwRQAjp6cGbWu7D\nBMAQ7fx5LM+qJwIAsFcAxoxBeymSp+apCQAbtcv57G1t+BnExqq/l1URALOAAAAWLAA4fNj8OUOF\nUgRAAhBALgcAgAKwf7+2RWAMsQ0kN9CkCMAG2ErTGIW/ZtAgrMv90Uf6BcDKjeGDBcDlivwoQE0A\nYmNxFa/cCl6t9g8ACoAVHZzYArrjDpxAILW3NI9QBKANuRwAAMCdd2IhQLZlqRZmzQI4fRp/piRw\nCFHz/xnz56OdoNcCsmpUCTBQAAAAUlMje39gNQEAUO6Y9AqA1RbQ4MGYGPzLX8yfNxQoRQBWCWQk\noGQBjRgBMG+e8kYwwcyahRHAzZv4GUitHiYLyAaU1gCIyclBBVbrjIKxWwAWLgQ4dMia8/NGdzda\nXKwmuhw8CUBf38Ccw6pVAMXF5s4bKtSmRJMAIEoCAACwfLl2+wcgYAH5fLguQMqyJAvIBpRWAYvJ\nzcVNPvRi1Zemrw+Xlrvd/Z9fskR628pIgI3+1UZRSrNTQi0AX36J54kTFUv/5jdxf4AbN8ydOxQ4\nYQE99lj47camJgAPPwzw059qP19aGtaNqqyUdxkoArABrRbQ1KkAr7+u//xWRQCHDmEbBg/u//y8\neZicjsSRmRb7B0C5Y2prw8qcWrDis5IqXjduHMCcOSgCvKOWBFZKuBvlpZcAdlm+z5999PXh1G65\nHAAA1gVavFj7OePiAKZPR6tQTgAoArABrQJgFKsE4Pe/B9iwYeDzgwdj+Fhebv49eMMKAWhv1y4A\nVkwDFc8AEpOfHx6dnFIEoJZwN0JfH352b71l3TnN0NamfkxnJ9bw1zKzTA+zZwP8139RBBBStOYA\njGJFp9LeDvDuuwDr1km/vmQJwMcfm3sPHtEjAHKL7dratM/EsEKsxTOAxKxahQJg1apwu1BbNW21\nDXTtGorKwYPO20DNzTjLT+0eULN/jDJrFtq8cjMN2UpsXkqLRIQAhCICMPuF2b4dk0ri6oBi9OYB\nKisBfvlLc20KBaGOAKz4rOT2L0hNxciA90VhapVTrRaAzk68LkuXOj9TqrUV75fdu5WPs1MAAOQj\ngIQEjPidFkpGxAiAliSwUawYVb74IsBDD8m/vmgR5gi0ji5PnQL4t3/jZyQhh1U5gFBGAEob2Kxc\niSNdXunqwntIqdqtVaulGZ2dmMRfswbgzTetO6/RtgAA7NypfJzSIjAzzJ6N/ytNNR81igTAUuyO\nAMx+YSoqsI3Ll8sfk5SEf0NlpbZztrdjPfOaGuPtCgVaBUBpFpDeCMCuHAAAQEYGrijnFTZjSmnW\nldURAPuMV63CJLlVdbOMtmX2bIB9+5TbYVcEMHEiwD/8g/L00REjnL1GYiJCAC5f5tsC+o//APj+\n95VXKgPoywOwEcRf/2q8XaEgHGcByeUAANAGYpuF8IiWqql2WEDDh6NoLl7srA3U0YGd75IlAO+9\nJ3+c0ipgM7hcANu2Ke8gZmVlAbOEvQBcvow3YPDiKisx06n813/hF+Lv/k792MWLtecB2tvxRo8k\nAeAlCaxkAaWl8S0AWnZOs0sAANAGevtt685tpC0jRgCsXatsA9kVAWhh+PAIEoCSkpKV6enp59LS\n0qq3bt36ZPDrxcXFq7Kysk5lZ2dXzJs37/jBgwe/YfY9xRw+jCtp1UbXZjBqAZ0/jyP/N9+UT/6K\n0RMBtLfjtMQPP9TfrlDiRBLY7Dx3JQvolluwk+HFww1Gix1qpwDMmeOsLcnut9Wrsaa/3OZAduUA\ntBAxEYDf74/duHHjtpKSkpWVlZWZ27dvv7+qqipDfMxdd9114NSpU1kVFRXZf/jDH/72kUce+Z25\nJvfn0CFMoNqJEQuoowPg298GeOYZ7e3LyMD9ArTM0W5vx9pGN2/yXUjOrAAIgr4IgM1zv3lTXzvF\nKFlALhffNlBtrXqtK7tyAADObzrPxGjUKFz5L1e+w8kIIGJyAOXl5QtSU1NrUlJSvPHx8T0FBQU7\niouLV4mPGTZs2Fe3Wmdn5/Bx48ZZWlMxVAKgNwLYsQNX/RYWav+dmBgsE9HQoH4su4HvvJNvG8is\nAFy/jn5qQoL29zSbs1GygAD4toGOH8eV5UrYGQE4LQDi++3uu7H6rxROCwAvEUCc+iHyNDQ0uJOT\nk33sscfjqT9y5MjC4OPefffd1T/5yU/+tbGxcdK+ffvuljvfpk2bvvo5NzcXcnNzFd+/pwc3eVk4\n4B2txYgAtLUBTJumvZIgw+PBKCA1Vfk4dgN/7WtoA61fr+99QoXZWUB6Rv8M9nlprR8UjJIFBMC3\nAJw4AbBli/IxVlcEFQvAqFH42O+3fpWt1rZMmoQ/jx8vv+rWriSwFszkAEpLS6G0tNSytpgSAJfL\npclpXb169burV69+t6ysbOm6dete++yzz6ZLHScWAC2cPo0r7rT6w0YxMmK6dk15LrYcTADUEEcA\nv/mN/vcJFR0d2r5ocklgPTOAGGYTwUoWEAAKAI9R15Ur2Ha1wYMdEQD7jGNiUPDb25VF1C46OnDg\nBaAcjYRrBBA8MN68ebOptpiygNxud4PP5/tq/o3P50v2eDyy3dfSpUvLent7465cuWJwbNafUNg/\nAMY6lFAIwMiRADNmYCepxTZyArMWkJ4EMMOMAPj92KEpdQ68RgAnTgBkZ6tPiLAzBwDgrA3EZgGp\ntcPpJHBE5ABycnKOVVdXp3m93pTu7u6EnTt3rs3Pz+9XLqu2tnaqIAguAIATJ07MBQAYO3as4Xp4\n4uReJAqA260vAoiJ4buOkFkBMGMBGYG9n1InyrMAzJ2rfpydOQAAZwWgo0NbPiJcIwCrMSUAcXFx\nvdu2bdu4YsWKvZmZmZVr167dmZGRUVVUVFRYVFRUCADw9ttv3zdr1qxPs7OzKx577LHnduzYUWD0\n/UpLMbHa1ISPP/lEX7lWoxiZBmpnBCAI/UcwyclYBItHnIgAzKzcVrN/AHDV9o0b2qpOhpITJ9QT\nwACRLwDsflPaftHJHABPAmAqBwAAkJeXtycvL2+P+LnCwsIi9vMTTzzxiyeeeOIXZt8HAMskdHUB\nfPe7AH/6E36IzO+zk/h4rK/S06O8wk/MtWs4EtWLFgHo6sJZMawto0bJz3d2GieTwEZQmwEE0H8q\n6Pz5xt7HDo4fB3j6afXj7LCAxAIwerSzFpCWhDQtBEPCaiXwhQsATzyBHUp+PsDtt9u7AIzhcunv\nVMxEAGp+fvDNy6sACMLA0aEcVucAjHZwajOAGLzZQG1tGBlPl5xe0R87IgBecgDiAUdMDI7ygyO1\n4Ag61PAUAYSVANTVBXb16uhA7ztU6LUVjApAUhLO5ujulj8mXASARSpxGuJMJrDB1VBDHQFosYAA\n+BOAkycBsrK0Tb20qxooY8wY5+yx4LYkJg4Uoxs3UBzU9qm2C56SwKYtoFBy4QLAbbfhDXbkCO7o\nEyr0jiqNCkBsLFYUbGyUX9EZLAAjR/IpAFrtHwD8Qg4ejAu/xNetrU1/nSe7LSAAFIADB4y9hx1o\nTQADRE8OgLUlOA/gpP0DQBGAIQQBBeDWW/Hx2LHGPHajhMoCAlDPA4RLBKBHAACkO6ZQTwPVIwBO\nRwANDYGyIcePOycAwTkApwRAEPDvUmuLkwlgAMoBGOLqVfTitXw57YA3ARDfwJEsAE5YQHpzADdu\nOLPPa0EBJqOLinBFvJYZQACRmwPo6kJbR2yDyQkARQBI2AgAs3/0llawCr1fmq4uigD0CoDUTKBQ\nTwO9ckWbAIwfj7NLsrOxfSpVS2zhwgWA55/HarMNDVhMUAtWCoDfj2tzxHasUwIgdb9J5QCcTAAD\n4ADl5k28dk4TNjkAsf3jBKGOAHw++delBIDH8sRGIoDg5FioI4AvvtBWutvlAnjnHXyv9HQsE93d\nra9onRlu3sTSz6tWAXznOxiBaJ2ibKUAsPtcPDBzSgCkZpzxmANwubCdaivOQ0HYRAB1dRgBOIWe\nTsXvxzUDgwcbey+1CCB4BBMpEYCVOQCjHdwXX2jfX3rZMlyJPmYMbs4TypyAz4erxpndoafuTkJC\nYF2LWYL9fwDrBOD6dX1CLnW/8WgBAfCTBwgbAWAWkFPoGTWxRWBG7Sq9FtDw4fhF4SGkFGNVDiCU\nSeDLl7VFAMFkZgKcPWvsPY3g9WIhRCOwdS1WRAHB/j+AdQLw7LMAP/6x9W1xOgkMwE8eIKwEIFws\nIDP2D4B+AWAVGHmzgcwKQE8PJlj1XksnBGDGDFypHio+/1x94xclrLKBpGwXZkkGr+nQy5kzAO++\nq/08UtEIjzkAAH7WAoSNAISTBWRWACZNwto+vb3Sr0uFsDyuBTArAOyLqjeSMioA167hVEIjn92M\nGeETAQAvmtirAAAgAElEQVRYJwBSnW5cHH4GZgck587hOY4f13a8XATAWw4AgCIAXfT2oudpZsRj\nFj0zS8wKQHw87usqV+BN6gbmMQ9gZBaQeFRkJAEMYFwA2OjfiHWXmUkRgBizNlBPDw76HnoIowAt\naM1H8CAAlAPQQX09fjGdWroNoM8zNSsAAMploSNVAII7JSMJYHYeIwKgJwEczPTpaFMqlfCwEl4i\nAKlRN4D5chC1tWiFrl0rv69vMHqSwJQDQMJCAJy2fwBCawEBKOcBpG7gSBQApyIAIwwejCUramqM\n/b5eIj0COHcOp9cuWIBrM7RcV6m2UA5AmbAQAKcTwAChtYAA1AVAKgKItCSw0QjA6AwXMxEAQOjy\nAD09WPnT4zF+DjtzAADmBaCqCgUgJgbg3nu1RQFS99uIETidVDzllQcLiCIAHTg9BRQg9BaQnAAI\nQvRYQEamgAI4EwEAhG4mUEMDVo3VuvBLCqsqgtoZAbCVzatXaxMAKTvK5Rq4PwEPAkA5AB1EowUk\nlwO4eRNv6uBFZpEqAGYsIEHQ93tmI4BQrQUw6/8DhCYHYIUFBADwjW8AnDqFn48SWqMRHgSAIgAd\n8GABhVoAxo/Hpf7ByN28kTANNHgWkFELKD4eRVLvSlcrIoBQCIBZ/x+A7xyAIPQXgMGD8dqeP6/e\nFqn7LTgPwEsSmHIAGqmrc14A9K4ENisAY8di8isYOQGgCKA/RmwgsxEAmwlkRYkFJXiKAOzIATQ2\n4oy/sWMDz8l9H7S2ha0F6OnBmVpmv59moQhAI9euYXJz0iRn2xHqCIAEwHgEwM6lVwDMRgBsJpDd\nNYEiPQIQj/4ZWgVAzY5qbsY1Nk5VFWaQAGiE3exOf2ChFoBx48JbAPx+/WUcwj0CAAhNIpinCEDO\ndjGzMbycAEhZolraIhaAigqAOXOMtctKKAmsEStudivQOw3U7G5lw4djqHrzZv/nw0UAOjsHlglW\nw8oIQO9UUEHACMCsAEyfDvDZZ+bOoQZPEYAdFlBV1cC9DeQGRFraIs4B6Nk6004oB6CRujo+BCDU\n00BdLumwVy6BxZsA6LV/AAbuB2B0GiiA/gigsxNLK5v93CZMUB+pmsHvx9lhkyebO4+V1UB5sYCU\nIgCWAzh+XPvOaXYSMRZQSUnJyvT09HNpaWnVW7dufTL49T/96U/fy8rKOjV79uzTS5Ys+fj06dOz\n9Zzf63U+AQygb2qhFQIAIC8AchEATwvBjNgpwTuChdICMuv/MxIT7d0e8tIlvC/MlkUJRQ5AXApC\nT6lyMzkANTHiKQIIewHw+/2xGzdu3FZSUrKysrIyc/v27fdXVVX1C95uu+22Cx9++OGdp0+fnv0v\n//Iv//eRRx75nZ734MUCio/HVYlaZnjYKQByy9h5iwAuXcJdsvQwZAjmDVj5X7MWkB4BsML/B7Bf\nAD7/3JrvQ6hqAQkC2pi33QZQUqJ+vo4OvOeDLS61HACrwSQljMwCam7Ge4KH/iQicgDl5eULUlNT\na1JSUrzx8fE9BQUFO4qLi1eJj1m0aNGhUaNGtQMALFy48Eh9fb2uBey8CACA9rDZiQiAjSjM1mC3\nCiMCEBODIsAiLTPztSM1AvB6ramKa3cOID4eO+POTlzFGxMD8PDD6gXiKiowkR4T1DOp5QDk2gEQ\niADY6N/pCSUAeP1v3HB+EydTewI3NDS4k5OTv9q91uPx1B85cmSh3PEvvfTShnvuuWe33OubNm36\n6ufc3FzIzc3lTgC6uvCGUsJuAcjMHHhsXBx2np2dzi9yATAmAACBjqm5Gb/QRssd6J0GalUEMHas\n/QLASwTQ06O89SnreF9+GeBnPwP4+GOAxx8H+MMf5M9ZUgJw990Dn1ezgJRyTsECwAMuV+Az0PN9\nLS0thdLSUsvaYUoAXC6X5sX2H3zwwddffvnl73/88cdL5I4RCwAAdmbXrlkzMrMCrZ2KExEAQMAG\n4kUAsrP1/96wYfhlfvhhgH/6J+Pv72QEoOZVm6GuDmCh7BBLO1YIwLVrKNJyI+oxYwBOnwY4ehTg\nz3/GDeyzsgB27QLIz5f+nT17AJ5/fuDz7LsgCNLvJ2dFsXa0tqIA/M3faPvbQgGL2vV8X9nAmLF5\n82ZTbTBlAbnd7gafz5fMHvt8vmSPxzOggs3p06dnP/zww7/ftWtX/pgxYzTPDWCjHR5CNgC+LSAA\nvvIARiOA4cMBfvEL9I1/+EPj7693lssXX1gjAGy0qbcOkVasWhVvhQAodboAeC1+8xus6T9kCH62\nzz8P8Mwz0sc3NuJ3/vbbB76WkICRhtxEByULiOUAeJkBxOAhD2BKAHJyco5VV1eneb3elO7u7oSd\nO3euzc/P3yU+5uLFi5O/853vvPP6668/kJqaqqtaOk/2D4D2UWVXFwmAGQvo9dfRJoiNNf7+RiIA\nKywg1lHZ9cW2UgDMVgNV6nQBUAAOHgTYsCHw3Lx5mMiWYu9egLvuQjtTCqU8gJIYDR2KuwpevQow\ndap8e0MNDzOBTFlAcXFxvdu2bdu4YsWKvX6/P3bDhg0vZWRkVBUVFRUCABQWFhY988wz/+fq1atj\nHn300RcAAOLj43vKy8sXaDl/OAqAIOAxZheCAehbBwDgnACcOYNfMPEKS6MC4HYDPPAAwLRp5to0\ndKg+K8aqCAAgkAi22orr7cXranYNAABen+vXcdJAcMJVK3JTQBljxgDMnt3fdx8/HkfxN24MzB3s\n2QOQlyd/PvZ9kKoMrCRGLhe2JSPD+N9qBzwsBjMlAAAAeXl5e/Ly8vaInyssLCxiP7/44osPvfji\niw8ZOTcvi8AYWsLm69dx9oOZ0SvDSARg1VqAnh6ckbFAg1T/8pf4ZS7670+9txfbbaRD/c//tOba\nDR2K+0hrxaoIACDwuVl97/p8uA9AQoL5c8XE4Gd2/brxaFVNAObNw1LOYgs3Jgbrel261L8j7+0F\n2L8fLSM5lKaCarGjeEkAM3iIADjSw4HwsgiMoSUCsMr/BxgoAHKbwTCsjAD27sWRuBqCgGG+uFQv\nK7glF8orYUXnD2BsHYDVEYDVWF0V12weQG2198aNAOvWDXze48FNbcQcOYKRjVLUqGQBqbUlMZE/\nAeAhB2A6ArCTcLSA7BSA1laMLuTsJSv3BCgrwxGn3KwLRm0ttktcAbOx0fnqrXo8bkGwbhooQPQI\ngFoEIIfUZkd79gCsXKn8e0pTQdXyEc8+y58AUASgAm8CoKVTsVIAEhNx4Qxb3FVbi0ksuQ7Zygjg\no4/Qp1Xz0Q8exC37rlwJdCZG/X8r0RMBfPklCqvcfHa9kAAoIxUB7NtnTgDULKDcXD6mR4vhIQfA\nrQCwRNG4cU63JICWqYVWCkBcHJ6LdepqeyNbJQDXrwOcPInvpeajHzyIMzemTgWo+e85XrwIgNbO\nzUr/H0Bb3RojWC0AZgvCWRUBCALupKY2RVMpB6AWAfAIRQAK8LYGACD0FhBA/87kwgXlaWxWCcCx\nY7gcPz1dWQAEAeCDDzDRl5YWsIF4EQCtEYCV/j9A9EQARiq+AgyMAC5dwvOoncvoNFBe4SEHwLUA\n8JQABgi9BQTQXwBqa81FAH19AHfeqV6mt6wM4I47cHcrJQGorMSbeMoUnLbJEsHhJgCNjTi7xiq0\nCsA99+jbPYw3AbAqAqiuxgGEGmZyADxCEYACvPn/AKG3gACsjQBOnsTOvbZW+T0/+kibABw8iKN/\ngPCOAOrr8W+1Ci0C0NmJM63OnNF2zuvXUbitvK5mOyCjnW5wBKBHAIxOA+URygEowKsAhHME8Je/\n4P9KnbrfD/DJJ8YEIFwjACcE4NgxjMjkVsUG4/XiNEkrFzKJd8oygtENeyZNAmhqClTCrK4GSE1V\n/z2KAKyHWwFoagKYONHpVvTHSQG4eRPn1yt1VCNHKi8E270bff2LF+WPOXsW7ZAJE5QFwO8H+Otf\nAb7+dXw8bRpfEYCeaaA+H45KrUJLEvjwYayPo1UArLZ/AMznKq5eVa+MK0VCAr735cv4WGsEoJQD\nMNoWJ6EcgALt7fx9oFo8U7sEwOvFDllpcZVSBNDSgp79976nPKpn9g+AsgCcP49fYibSEyeiTfHF\nFzgytHJWjRH0RgBWCoCWjvXQIYBvfUufAFgdEZutXGqm0xXnAbQKwNChGDVJfa4tLXzNGNTC4sUA\n773nbBu4FQAz2wHahZMRgNoUUICAAEhVoty7F0frqanKAsASwADYKV66JL3JTFNT/07T5cIvcVkZ\nRg9O11xJSMDyAlp2cLNaAFj5YbmKoIKAEcDatc5GAGp7FzzzDMD778u/fvUqiogRWB6grw+tTS0W\nkNw+2YKAf8fYscba4hQJCc6vTeBaAIxuB2gXTgoAWwSmREICbqAiFaXs3o2zTpKTlS2g8vJAOd7B\ng1FUmpsHHidVPz8tDW0hp+0fAOws1CwxAOyArLasBg0K7IYlhdeLkdwdd/BrAQkC1nY6fFj+962I\nAOrr8Rxa/XspAejowOttRY2kaINbAVCqeeMUWi0gKyqBMvREAABovQTPlPD7MQLIy1NP7AaP7JOT\nBy7bB5AXgNJSPgQAQNu6iMuXcaBh1SpghlIe4PBhFNkJE/B+0TIV0y4BkGvjmTMojHIC1ddnbpDG\nIgCt9g9DKg9w5Ur4jf55gVsB4DUCcCoHoCUCAMAEbvCIvbwcZ14kJ2PnfPmytDXS1YViIR6NeTzS\ngiFVO2faNNwBihcBGD1aXQCstn8YSqPrQ4dQAFwunNmjJQoItQVUUoIDDq9X+vWODkxiG92yk0UA\negVASljD0f/nBS4FoLcXy0DwNq1r5Ej1rL2TOQAAFICmpv7PlZQE6qzHxeExly4N/F3WqYtXX8tF\nDHIRAAA/AqAlArB6CihDSQAOHwZYtAh/njJFXQDa21Gwre7klNpYUgLwyCPybTM768ZoBCC1FoAi\nAONwKQBs0xOeykAAKCdZGXYIQEuLdgGYOHFgBPDZZ/03a5k8WX5UH9zJ6BEAtokLTwLQ1qZ8jNVT\nQBlynev16zjVllWmnDJFfpTNuHwZRdvq7wPbvjI4yd/ZiVHj3/4t5ouk7ncrBMCqCIAEwDhcCgCP\n9g8AhruDBinbQFYLwLBhaMsMHaptxoCUBRSc5JRLBLe0DLR1lAQg+NixY7FT4EUAnLSA5OyVigrc\nmYrlibREAGZm2ygRH4/3V3Ci/IMPcCOgpCSMwqUmAZgVALfbuhwAWUDG4VIA2tv5FAAAbJfSqNJq\nAWBT37TuZSplAUkJgFZfX+lYqQJq69YBzJypra12o9UCsisCkEqwlpcDLFwYeKxVAOxaEyMVqZSU\nBEozy7XPbJuGD0cBqq7Wt08vWUDWwqUA8LgGgBFqAQDAm1uL/QMwMAIQBBQA8QYtkydLRwB6BEDK\nAgIAeO45ezpUI2ixgEKdBPb5+n+WWgSgtTV0AiAI/TdnSUmRFwCzUYnbjZalnllzZAFZC7cCQBFA\nAD0RQHAOoL0dE7/iQll6IgC3G8/X2xt4rrsbfWJePyMGj7OAgoXT6QgguEOtqcGyIyyKk8tRWNEm\nj0ef/QMgPwuIBMAYXAoAj2sAGGq2gh0CkJSk/YsSbAFJLXKSSwJLeanx8fhcY2PgOZYsdnq1rxpq\nn5UgoA8dagEQi6zbjdezu1v+XHblAAAGtvPoUYAlSwIJZ7ssIAD82/UKwPjxeL3EXLlCOQCjcLkn\ncLhGAIJgjwD89rfal4wHW0BSAiCXBJbbF5dFDGy6pNUbqNiFmgC0tKAXPWSI9e8ttxAsOAKIi8PP\nx+eTj/KuXrV2vwIxwQJw8SJ2+oyUFNyuMZjWVvPCuWyZ/gV4t9wSKE/CBiBkARmHyzFcuCaBe3pw\n5GT1kvRx47Sfc/RoDOGvX8fHUgIwbhwu+gqezaQmAAw5/5831HIAdk0BBdBuAQGo20B25gCCZyuJ\nhZ61zS4L6HvfA7jvPn2/M2QI3uPiQQ5ZQMYxLQAlJSUr09PTz6WlpVVv3br1yeDXz507l75o0aJD\ngwcPvvHrX//6h1rOyXMSWGlUacfoXy8uV/8oQEoAXC7pPIAeAXC62qcW1HIAdvn/ANICIAjS11hN\nAOy2gMSRipQAfP75wLUATpZfDp7EQBGAcUwJgN/vj924ceO2kpKSlZWVlZnbt2+/v6qqKkN8zNix\nY688//zz//ijH/3oV1rPG64W0Jdf8rErkTgPIFfoTEoA5OZTp6YGNnwHCK8IwGkBEHecbW0442XQ\noP7HahGAUM0CChaA0aMBYmMHbhzDiwBcv452kNMDr3DFlACUl5cvSE1NrUlJSfHGx8f3FBQU7Cgu\nLl4lPmb8+PFf5OTkHIuPj9dQmBfhOQmsJgA8tFs8E0irAPT0YPulRprp6QDnzgUeh1MOQMkCslMA\nWHVKcUVQOeHk2QJi7Qu2gXgRADb6561qQLhgKgnc0NDgTk5O/qob8Xg89UeOHFmo9DtKbNq0CQBw\n79qlS3MBINdM82xBSQB4ES41CwhAOoxOTJSe2RMsAJcv61u84xRaLCC2paUdMHuFRYVKAvCnP8mf\nx+4IgFlAXV0oWMEWFVsLwMpXsDbZZUupMXlyQJCizf8vLS2F0tJSy85nSgBcLpdCVRz9MAHYtQsg\nN9fKM1uHkq3AkwBosYDEtd7l/H8A/P1r1wIdUbjkAAYPRnvg5s2BtguAfYXgGOxzYDt5GY0AQjUN\ntL4ep2YGDwKCI4C+PmcnakyeDPDhh/hztE0Bzc3NhVxR57h582ZT5zNlAbnd7gafz/fVV8jn8yV7\nPB6J6vH64DkJHA4RALOABAHn74tXATOCIwCleiouF0YBn32Gj8MlB+ByKQu2nRYQwMCOXe66TZ6M\n6xHYJuliurtRwOyqjCsWACn7B2Dg39HRgbkMpe1J7UTKAiKMYUoAcnJyjlVXV6d5vd6U7u7uhJ07\nd67Nz8/fJXWsIAiaXbpwTQKzKqZOwyygK1cwOSY1zz0tLdChAyhHAAD9baBwyQEAKOcBGhrsLVyn\nVQAGD8ZOTKpE99WreM/Z5XEnJgYqgsoJQHA5CKc3YBcLQLRZQFZjSsPj4uJ6t23btnHFihV7/X5/\n7IYNG17KyMioKioqKgQAKCwsLGpqapo4f/78o19++eXImJiYvueee+6xysrKzOHDh0tumCcImIzk\noSOVQslX5iUCYNaD0laHKSk48mPXWo8AhIsFBCD/ed28iaNrO/ecSEkZmDyfPl3+WK93YAdst9ce\nFxeoCKoUAYgtIKcFYMIEzFV0dVEEYBbTQVxeXt6evLy8PeLnCgsLi9jPEydObBLbRGp0duKIyOhO\nQ3ajNKLkRQCYBaQkADExAJmZWJt+0SJtAvDaa5gL6O3lY7qrFuQsINax2jl7ZMoUrKzJuHwZYOlS\n6WOZAAS/HorOltlAPh9AdvbA14MjGTtnJWlBvI7lypX+K5cJfXC3EpjnVcAAKE6CgDuWBcPLNFBm\nAaltdj5zJu79CiC9F4AYFgEw+ydcpt3JCXYoOlatFhCAfNXNUHS2bCqoXAQwbhxOE2bX0ekIACBg\nA5EFZA7uBIDnBDAAdnxyeQBeIoBRo9DeqK3VLgBSu4GJSU3FDurSpfDx/wHkLaBQCgBbDKYkAHaW\nXFCDTQWVEwCXC2DWLIBTp0LXJjWYAJAFZA4uBYDnCABA3lbgRQBcLuxoKiqUBWDGjP4CoBQBDBqE\nncOhQ+Hj/wPIf1atrfbPYx89Gq02topWLQKQEwC72ym2gOSmxc6Zg+tzQtUmNSZPRnElATAHdwLA\nSyeqBO8RAADmAdQEIDgCUOvY09Nx/nU4RQBOWkAAgSigtxfvD7mOU0kAQmEBeb3YRrn3ys7G+ylU\nbVJDHAFE0zoAq+FOAMIhAlASAF5mL7GZQEoCcMstaBVdvqyeAwBAASgrCy8BcNICAgh4+y0t2PnH\nxkofx/ZoCN6gPRQ5gMREtHeSk+VzO8ERAC8CQDkAc3AnALwngQHkOxWeIgBWP15JAFyuQBSg5YuU\nno5f/kixgEIVAXi96ovnhgzBjli88Q5A6CygkyeVV0XPnInrRrq7+RGA2lqclcbLdy4c4U4AeE8C\nA8jbCjwJwMSJ/f+XY+ZMgI8+kq5SGUx6Ov4fThGA2jRQu2EWkJbV01I2UKgsoOpqZQEYMgT3Mj57\nlg8BYJsaydWvIrTB3aULVwtIEHCJPE8W0Pjx6hvJzJwJ8MEH2kb1bBFTuAkADzkAowIQKguor0+9\nLlJ2NkYKPAjAkCF4z5L9Yw7uBICnUbQcUgJw7RqOoHlZwJaUpK3MwcyZ2mf2jB2Lx4WTBeR0DkAs\nAGrXzakIgEVCagLA8gBOLwRjTJ5MAmAW7gQgHCIAKVuBN+FasABg/Xr142bMwLIIWjv1114DyMoy\n17ZQwkMOwKwFZLdVxTpRLQJQUcFHBABAAmAF3AlAuCSBgyMA3gQgJQXgf/0v9ePGjcNoQetUuhUr\nrN/z2E6ULKBQ5AAmTMDosK5OXQCc2nhFbwTw5Zd8fEcnT6YpoGbhTgDCIQksJwC8+P96mTkzvGwd\nPYwahR2WU3vaulzYsR89qi0CEJeDYNsdSlVztZIxY7AonJoAjBuH13PYMOdKQYu54w6AhYa3nyIA\nLCgGZzXhYAGFQwSgh8WL7d0YxUkSEjAv09UV2DdWEELrY0+ZArBvn7YI4OJF7PTZCmK7C9YBYGd+\n9qy2An9z5gB8+qm97dHKmjVOtyD84U4AwqEjDYccgB6eecbpFtgL+7yYAFy/jh2s3SNrBqtWqSYA\nQ4bg4IIt4Aul1z5tmrbj5szBjXSIyIBLCygcIwBeKoESAwnOA4Q6ialVAAD6J4J5mW0jJicnvKYB\nE8pwJQA3b+K2eKEamRkl0iygSCd4KqgTAjBokDaLRSwAvMy2EXPvvQBvvOF0Kwir4EoAWCfKe635\n4cNxP4De3sBzJAD8EmzZhXpkPWWK9j0UggXA6aqbwcTE0NTLSIIrAQgH+wcAv8gjR/bvVEgA+EXK\nAgplx5qTA/CLX2g7dvFigHfeCX2imohOuBKAcFgDwAi2gUgA+MVpC2jIEICCAm3HfutbuPvWe+/x\naQERkQVXAtDaGj6dqJQAhOs6gEjHaQtIDzExAJs24b9QbFpDRDdcCcA77+DijnAguFOhWUD8EvxZ\n8eiti1m9Gi2gt9/mV6iIyIArAXjzTYC//3unW6ENsoDCh+DPindrxeUCePppXA/AczuJ8IcrAViz\nJnzmGJMAhA9SEQDvHeuqVQB33w1w661Ot4SIZEwLQElJycr09PRzaWlp1Vu3bn1S6pgf/OAH/y8t\nLa06KyvrVEVFRbbcubQUL+OF4MQiCQC/jB2L1TgZPOcAGC4XwN69WK2VIOzClAD4/f7YjRs3bisp\nKVlZWVmZuX379vurqqoyxMfs3r37npqamtTq6uq03/3ud488+uijL8idLyND7hX+EE8tFAQSAJ6Z\nNQv3vGUF4XjPARBEqDAlAOXl5QtSU1NrUlJSvPHx8T0FBQU7iouLV4mP2bVrV/769etfBQBYuHDh\nkba2ttHNzc1JZt6XB8aPR48WABeFxcSob6lIOIPHg51/QwM+DgcLiCBCgalicA0NDe7k5GQfe+zx\neOqPHDmyUO2Y+vp6T1JSUnPw+TZt2vTVz7m5uZCbm2umebaSnQ3wyiv4M00B5RuXC2DuXIATJ1AM\nwsECIggpSktLobS01LLzmRIAl8slqB8FIAhCv0Xwcr8nFgDeyc4GqKzE0T9NAeWfefMAjh/HWjZt\nbSQARHgSPDDevHmzqfOZsoDcbneDz+f7qpK8z+dL9ng89UrH1NfXe9xud4OZ9+WBIUMA0tNxizzy\n//ln3jyMADo70aoLp13NCMIuTAlATk7Oserq6jSv15vS3d2dsHPnzrX5+fm7xMfk5+fv+uMf//g/\nAQAOHz58++jRo9uk7J9wZOFCgPJyEoBwYO5cjADI/yeIAKYsoLi4uN5t27ZtXLFixV6/3x+7YcOG\nlzIyMqqKiooKAQAKCwuL7rnnnt27d+++JzU1tWbYsGHXXnnllQetabrzLFgAcOAA+sokAHwzZQqW\nG6+sJAEgCIZLCN4s1SFcLpfAS1u0cvYsLtv/yU8AysoCSWGCT5Yvxx2tjh4FsDCPRhCO4XK5BuRY\n9cDVSuBwIz0doLkZoK6OIoBwYN48jNgoAiAIhATABLGxgU6FpoHyz9y5ACdPkgAQBIMEwCQsEUwR\nAP/Mm4f/0ypggkBIAEyyYAFAXx8JQDhw2234OVEEQBAICYBJFizA/0kA+IetCCYBIAjE1DRQAsDt\nBpg0iQQgXPjpT3HjdYIgaBqoJbzxBtZuHzfO6ZYQBBFNmJ0GSgJAEAQRptA6AIIgCMIQJAAEQRBR\nCgkAQRBElEICQBAEEaWQABAEQUQpJAAEQRBRCgkAQRBElEICQBAEEaWQABAEQUQpJAAEQRBRCgkA\nQRBElEICQBAEEaWQABAEQUQpJAAEQRBRimEBaG1tTVy+fPn+adOmnb/77rv3tbW1jZY67vvf//7L\nSUlJzbNmzfrUeDOjh9LSUqebwA10LQLQtQhA18I6DAvAli1bnlq+fPn+8+fPT1u2bNn7W7ZseUrq\nuAcffPCVkpKSlcabGF3QzR2ArkUAuhYB6FpYh2EB2LVrV/769etfBQBYv379q+++++5qqeOWLl1a\nNmbMmKtG34cgCIKwB8MC0NzcnJSUlNQMAJCUlNTc3NycZF2zCIIgCLtR3BJy+fLl+5uamiYGP//z\nn//8n9evX//q1atXx7DnEhMTW1tbWxOlzuP1elPuvffe9z799NNZsg1xuWg/SIIgCJ2Y2RIyTunF\n/fv3L5d7LSkpqbmpqWnixIkTmxobGydNmDDhstFGAJj7IwiCIAj9GLaA8vPzd7366qvrAQBeffXV\n9atXr37XumYRBEEQdmNYAJ566qkt+/fvXz5t2rTzBw8e/MZTTz21BQDg0qVLt3zzm9/8Czvu/vvv\n3/wmpIYAAASgSURBVL548eJPzp8/Py05Odn3yiuvPGhFwwmCIAiTCILg6L89e/asnD59+rnU1NTq\nLVu2POl0e0L9b8qUKd5Zs2adnjNnTsX8+fPLBUGAK1euJN51113709LSzi9fvnzf1atXRzvdTjv+\nPfjggy9PmDCheebMmZ+y55T+9mefffYnqamp1dOnTz+3d+/eu51uv93X4umnn97kdrvr58yZUzFn\nzpyK3bt350XDtbh48WJybm7uB5mZmWdnzJhx5rnnnvtBtN4bctfCqnvD0T+ut7c3durUqTV1dXUp\n3d3d8VlZWScrKysznL7oofyXkpJSd+XKlUTxcz/+8Y9/sXXr1icEQYAtW7Y8+eSTT25xup12/Pvw\nww+XnjhxIlvc6cn97WfPns3Myso62d3dHV9XV5cyderUGr/fH+P032Dntdi0adPTv/71r/938LGR\nfi0aGxsnVlRUzBEEATo6OoZPmzbts8rKyoxovDfkroVV94ajpSDKy8sXpKam1qSkpHjj4+N7CgoK\ndhQXF69ysk1OIAQlwLWusQh3pNaIyP3txcXFq+6///7t8fHxPSkpKd7U1NSa8vLyBU602w7k1ssE\n3xsAkX8tJk6c2DRnzpyTAADDhw/vzMjIqGpoaHBH470hdy0ArLk3HBWAhoYGd3Jyso899ng89eyP\nixZcLpdw1113HcjJyTn2+9///mGA6F5jIfe3X7p06RaPx1PPjouWe+X555//x6ysrFMbNmx4iZVb\niaZr4fV6UyoqKrIXLlx4JNrvDXYtbr/99sMA1twbjgoAzf0H+Pjjj5dUVFRk79mzJ++3v/3tP5SV\nlS0Vv+5yuYRovU5qf3ukX5dHH330hbq6ultPnjw5Z9KkSY0//OEPfy13bCRei87OzuH33Xff2889\n99xjI0aM6BC/Fm33Rmdn5/A1a9a89dxzzz02fPjwTqvuDUcFwO12N/h8vmT22OfzJYvVKxqYNGlS\nIwDA+PHjv/j2t7/95/Ly8gVsjQUAgBVrLMIJub89+F6pr6/3uN3uBqfaGQomTJhwmXV0Dz300Iss\nlI+Ga9HT0xN/3333vb1u3brX2BTzaL032LV44IEHXmfXwqp7w1EByMnJOVZdXZ3m9XpTuru7E3bu\n3Lk2Pz9/l5NtCiVdXV1DOzo6RgAAXLt2bdi+ffvunjVr1qfRvMZC7m/Pz8/ftWPHjoLu7u6Eurq6\nW6urq9MWLFhQ7mxr7aWxsXES+/nPf/7zt1lF3Ui/FoIguDZs2PBSZmZm5eOPP/7v7PlovDfkroVl\n94bTWe7du3fnTZs27bOpU6fWPPvssz9xuj2h/HfhwoVbs7KyTmZlZZ2cMWPGGfb3X7lyJXHZsmUH\nIn0aaEFBwfZJkyZdio+P7/Z4PL6XX375QaW//ec///k/TZ06tWb69OnnSkpKVjjdfjuvxUsvvfT9\ndevW/XHWrFmnZ8+efWrVqlXvNjU1JUXDtSgrK7vD5XL1ZWVlnWTTHPfs2bMyGu8NqWuxe/fuPKvu\nDcVaQARBEETkQjuCEQRBRCkkAARBEFEKCQBBEESUQgJAEAQRpZAAEARBRCkkAARBEFHK/wdhxevi\nbl4IcwAAAABJRU5ErkJggg==\n",
       "text": [
        "<matplotlib.figure.Figure at 0x5278850>"
       ]
      }
     ],
     "prompt_number": 10
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "contrasts"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 11,
       "text": [
        "{'constant_0': array(1.0000000000000007)}"
       ]
      }
     ],
     "prompt_number": 11
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "spline = (natural_spline(times, order=3))\n",
      "print(spline.shape)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "(120, 4)\n"
       ]
      }
     ],
     "prompt_number": 23
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "best_design = None\n",
      "best_efficiency = 0\n",
      "n_iters = 1000\n",
      "efficiencies = np.zeros(n_iters)\n",
      "splines = natural_spline(times)\n",
      "for i in range(n_iters):\n",
      "    X_ev, contrasts_ev = random_design(times)\n",
      "    X, contrasts = stack2designs(X_ev, splines, contrasts_ev)\n",
      "    C = contrasts['constant_0'][:, None]\n",
      "    efficiency = 1. / C.T.dot(npl.pinv(X.T.dot(X)).dot(C))\n",
      "    efficiencies[i] = efficiency\n",
      "    if efficiency > best_efficiency:\n",
      "        best_efficiency = efficiency\n",
      "        best_design = X"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 12
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "stack2designs?"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 13
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "natural_spline?"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 14
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import scipy.stats as sst\n",
      "m = sst.scoreatpercentile(efficiencies, 95)\n",
      "_ = plt.hist(efficiencies[efficiencies>m], bins=100)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAD9CAYAAAC7iRw+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAEuBJREFUeJzt3X9QlPedwPHPUyGdtmiURBYOuCzyKwgIOybqpTq3OQX1\njyiMiZFpOEZIxnqJqdXJhbbTCndXJZ0kFo3XOikmtF5IzGUg3o2gknRba2qYXJfEqImmLhPEZZPw\nI4F0zII890ezgRoEw/5g5fN+zTwzy7PP8nyfr0/e4rP7EMM0TQEA6PG1yR4AACC0CD8AKEP4AUAZ\nwg8AyhB+AFCG8AOAMmOGv729PfHOO+/8bWZm5qmsrKy3d+3a9bCISEVFRUVCQsIFm83mtNlszqam\nphWhGS4AwF/GWJ/j7+zsjO3s7IzNzc1t7e/vj5o/f/7/NTQ0FBw4cGDt9OnT+7Zs2fJkCMcKAAiA\niLGejI2N7YyNje0UEYmKiurPyMg409HRES8iYpqmEYoBAgAC65qv8be1tVmdTqdt0aJFJ0REdu/e\nvSknJ+fNsrKymt7e3pnBGyIAIKBM0xx36evri5o/f/4b9fX1BaZpisfjiRkaGjKGhoaMH/3oR/9R\nWlpac+VrRMRkYWFhYfnqy7V02Z9l3A28Xm9kfn7+4Z07d24e7XmXy2XNyso6OVr4a2trzaioYlPE\nNEUOmPn5d5sabdu2bbKHEDaYi2HMxTDmYlgowj/mpR7TNI2ysrKauXPnnt68efPPfevdbnec73F9\nfX1hdnb2ybG+DwAgfIz55u7x48e/vX///vvmzZv3ls1mc4qIbN++/Yd1dXVFra2tuYZhmElJSa69\ne/duCM1wAQD+GjP8ixcv/sPQ0NCX/lWwcuXKxuANaeqx2+2TPYSwwVwMYy6GMRehxZ27IcBJPYy5\nGMZcDGMuQovwA4AyhB8AlCH8AKAM4QcAZQg/AChD+AFAGcIPAMoQfgBQhvADgDKEHwCUIfwAoAzh\nBwBlCD8AKEP4AUAZwg8AyhB+AFCG8AOAMoQfAJQh/ACgDOEHAGUIPwAoQ/gBQBnCDwDKEH4AUIbw\nA4AyhB8AlCH8AKAM4QcAZQg/AChD+AFAGcIPAMoQfgBQhvADgDKEHwCUIfwAoMyY4W9vb0+88847\nf5uZmXkqKyvr7V27dj0sItLd3R2dl5d3NC0t7Wx+fv6R3t7emaEZLgDAX2OGPzIycmDnzp3fP3Xq\nVOaJEycW7dmz58EzZ85kVFVVlefl5R09e/Zs2tKlS1+pqqoqD9WAAQD+GTP8sbGxnbm5ua0iIlFR\nUf0ZGRlnOjo64g8ePLiqpKSkVkSkpKSktqGhoSAUgwUA+C/iWjdsa2uzOp1O28KFC1/3eDwWi8Xi\nERGxWCwej8djGe019fX14vWeF5GKwIwWAKYYh8MhDocjpPu8pvD39/dHrVmz5qXq6urvTZ8+vW/k\nc4ZhmIZhmKO9rrCwUJqbm8XrrRCRF0XklN8DBoCpxG63i91u/+LrysrKoO9z3E/1DAwMRK5Zs+al\n4uLi3xQUFDSI/PWn/M7OzlgREbfbHRcTE/NBsAcKAAiMMcNvmqZRVlZWM3fu3NObN2/+uW/9qlWr\nDtbW1paIiNTW1pb4/kIAAIS/MS/1HD9+/Nv79++/b968eW/ZbDaniMiOHTt+UF5eXrV27doDNTU1\nZVarte3AgQNrQzNcAIC/xgz/4sWL/zA0NDTqvwqam5uXBWdIAIBg4s5dAFCG8AOAMoQfAJQh/ACg\nDOEHAGUIPwAoQ/gBQBnCDwDKEH4AUIbwA4AyhB8AlCH8AKAM4QcAZQg/AChD+AFAGcIPAMoQfgBQ\nhvADgDKEHwCUIfwAoAzhBwBlCD8AKEP4AUAZwg8AyhB+AFCG8AOAMoQfAJQh/ACgDOEHAGUIPwAo\nQ/gBQBnCDwDKEH4AUIbwA4AyhB8AlBk3/KWlpfssFosnOzv7pG9dRUVFRUJCwgWbzea02WzOpqam\nFcEdJgAgUMYN//r165+5MuyGYZhbtmx50ul02pxOp23FihVNwRsiACCQxg3/kiVLjs2aNavnyvWm\naRrBGRIAIJgmfI1/9+7dm3Jyct4sKyur6e3tnRnIQQEAgidiIi/auHHjL37yk5/8m4jIj3/843/f\nunXrEzU1NWVXbldfXy9e73kRqfBvlAAwRTkcDnE4HCHd54TCHxMT84Hv8f333/+ru+66639G266w\nsFCam5vF660QkRdF5NSEBgkAU5Xdbhe73f7F15WVlUHf54Qu9bjd7jjf4/r6+sKRn/gBAIS3cX/i\nLyoqqvvd7373jx999NHNiYmJ7ZWVldscDoe9tbU11zAMMykpybV3794NoRgsAMB/44a/rq6u6Mp1\npaWl+4IzHABAsHHnLgAoQ/gBQBnCDwDKEH4AUIbwA4AyhB8AlCH8AKAM4QcAZQg/AChD+AFAGcIP\nAMoQfgBQhvADgDKEHwCUIfwAoAzhBwBlCD8AKEP4AUAZwg8AyhB+AFCG8AOAMoQfAJQh/ACgDOEH\nAGUIPwAoQ/gBQBnCDwDKEH4AUIbwA4AyhB8AlCH8AKAM4QcAZQg/AChD+AFAGcIPAMoQfgBQZtzw\nl5aW7rNYLJ7s7OyTvnXd3d3ReXl5R9PS0s7m5+cf6e3tnRncYQIAAmXc8K9fv/6ZpqamFSPXVVVV\nlefl5R09e/Zs2tKlS1+pqqoqD94QAQCBNG74lyxZcmzWrFk9I9cdPHhwVUlJSa2ISElJSW1DQ0NB\nsAYIAAisCV3j93g8FovF4hERsVgsHo/HYwnssAAAwRLh7zcwDMM0DMMc7bn6+nrxes+LSIW/uwGA\nKcnhcIjD4QjpPicUfovF4uns7IyNjY3tdLvdcTExMR+Mtl1hYaE0NzeL11shIi+KyCk/hgoAU4/d\nbhe73f7F15WVlUHf54Qu9axatepgbW1tiYhIbW1tSUFBQUNghwUACJZxw19UVFR3xx13vPbuu++m\nJyYmtj/zzDPry8vLq44ePZqXlpZ29tVXX/2n8vLyqlAMFgDgv3Ev9dTV1RWNtr65uXlZ4IcDAAg2\n7twFAGUIPwAoQ/gBQBnCDwDKEH4AUIbwA4AyhB8AlCH8AKAM4QcAZQg/AChD+AFAGcIPAMoQfgBQ\nhvADgDKEHwCUIfwAoAzhBwBlCD8AKEP4AUAZwg8AyhB+AFAmpOF/5ZX/FcMwxDAMmTEjOpS7BgB8\nLiKUO7t8+ZKImCIi0tdnhHLXAIDPcakHAJQh/ACgDOEHAGUIPwAoQ/gBQBnCDwDKEH4AUIbwA4Ay\nhB8AlCH8AKAM4QcAZQg/AChD+AFAGb9+O6fVam2bMWPGJ9OmTbscGRk50NLSsiBQAwMABIdf4TcM\nw3Q4HPbo6OjuQA0IABBcfl/qMU2TX6wPANcRv3/iX7ZsWfO0adMub9iwYe8DDzzw9Mjn6+vrxes9\nLyIV/uwGAKYsh8MhDocjpPs0TNOc8IvdbndcXFyc+8MPP5ydl5d3dPfu3ZuWLFlyTOSvfynU1tbK\ngw82S3//r0XkRRFZK77/A5eIIf7sGwCmIsMwgn4lxa9LPXFxcW4RkdmzZ39YWFhYz5u7ABD+Jhz+\nv/zlL9/s6+ubLiLy6aeffuvIkSP52dnZJwM3NABAMEz4Gr/H47EUFhbWi4gMDg5GfOc73/mv/Pz8\nI4EbGgAgGCYc/qSkJFdra2tuIAcDAAg+7twFAGUIPwAoQ/gBQBnCDwDKEH4AUIbwA4AyhB8AlCH8\nAKAM4QcAZQg/AChD+AFAGcIPAMoQfgBQhvADgDKEHwCUIfwAoAzhBwBlCD8AKEP4AUAZwg8AyhB+\nAFBmEsMfIYZhiGEYMmNG9BdrZ8yIHnX91Vxt+5Hrr/V7AYAGEZO360ERMUVEpK/P+GJtX1/PqOuv\n5mrbj1x/rd8LADTgUg8AKEP4AUAZwg8AyhB+AFCG8AOAMoQfAJQh/ACgTJiEf/hmrqutD+QNWH97\nc9cNQb3J66vekKbF9TIvwR4nNxoGRiD/nK6Xc9Mfk3gD10jDN3OJGKOuD+QNWH97c5cRlH2Mti9u\nIht2vcxLsMfJjYaBEcg/p+vl3PRHmPzEDwAIFcIPAMoQfgBQZsLhb2pqWnHrrbe+k5qaeu6xxx57\nNJCDmnockz2AsOFwOCZ7CGGDuRjGXITWhMJ/+fLlaQ899NBTTU1NK06fPj23rq6u6MyZMxmBHtzU\n4ZjsAYQN/gMfxlwMYy5Ca0Lhb2lpWZCSkvKe1Wpti4yMHFi3bt3zL7/88upADw4AEHgTCn9HR0d8\nYmJiu+/rhISECx0dHfFf+uZf+5oMDr4qN954l3zjG1X+jBMAECAT+hy/YRjm+FuJFBcXi4jIpUsd\nI1894cdfvsFrvG3+dvurP3ct+/BH5edLsL7/9aWysvLzR8Ge90AJ9jjHOk/1GD4vJiqQf07Xy7k5\nMRMKf3x8fEd7e3ui7+v29vbEhISECyO3MU1z6s0WAEwBE7rUc9ttt71x7ty51La2NqvX673hhRde\nuHfVqlUHAz04AEDgTegn/oiIiMGnnnrqoeXLlx++fPnytLKyspqMjIwzgR4cACAITNMM+NLY2Lgi\nPT39nZSUlHNVVVWPBmMfk7HccsstbdnZ2W/l5uY6b7/99hbTNKWrqyt62bJlR1NTU8/m5eUd6enp\nmenbfvv27T9ISUk5l56e/s7hw4fzfevfeOON+VlZWSdTUlLOPfzww9W+9ZcuXfr62rVrX0hJSTm3\ncOHCE21tbbdM9jH7lvXr1++LiYnxZGVlnfStC9WxP/vssyWpqalnU1NTz9bW1v5zOM7Ftm3bKuLj\n4y/k5uY6c3NznYcOHVqpYS7ef//9RLvd/tu5c+eeyszMfLu6uvphrefG1eYiHM+NgB/84ODgtOTk\n5PdcLpfV6/VG5uTktJ4+fTpjsk/QQCxWq9XV1dUVPXLdI4888rPHHnvsX03TlKqqqkcfffTRKtM0\n5dSpU3NzcnJavV5vpMvlsiYnJ783NDRkmKYpt99+e8vrr7++wDRNWbly5aHGxsYVpmnKnj17/mXj\nxo3/aZqmPP/88/fee++9z0/2MfuW3//+90v+9Kc/2UbGLhTH3tXVFT1nzpw/9/T0zOzp6Znpexxu\nc1FRUbHtiSee2HLltlN9Ltxud6zT6cw1TVP6+vqi0tLS3j19+nSGxnPjanMRjudGwA/+tdde+4fl\ny5c3+b7esWNH+Y4dO8on8+QM1GK1Wl0fffTRTSPXpaenv9PZ2Wnx/cGnp6e/Y5p//Zt85L92li9f\n3vTHP/5x0cWLF+NuvfXWM771dXV16zZs2PBL3zYnTpxYaJqmDAwMRNx8880fTvYxj1xcLpd1ZOxC\ncezPPfdc0Xe/+91f+F6zYcOGX9bV1a0Lt7moqKjY9vjjj2+9cjsNczFyWb16dcPRo0eXaT43rpyL\ncDw3Av67eq71M/7XI8MwzGXLljXfdtttbzz99NMPiIh4PB6LxWLxiIhYLBaPx+OxiIhcvHjx70Z+\n0sk3D1euj4+P7/DNz8i5i4iIGLzxxhs/7u7uDttfCB7sY+/q6rrpat8rVMf4VezevXtTTk7Om2Vl\nZTW9vb0zRXTNRVtbm9XpdNoWLlz4uvZzwzcXixYtOiESfudGwMN/rZ/xvx4dP378206n09bY2Lhy\nz549Dx47dmzJyOcNwzCn8vGPRfOxi4hs3LjxFy6XK6m1tTU3Li7OvXXr1icme0yh1N/fH7VmzZqX\nqqurvzd9+vS+kc9pOzf6+/uj7r777v+urq7+XlRUVH84nhsBD/+1fMb/ehUXF+cWEZk9e/aHhYWF\n9S0tLQssFouns7MzVkTE7XbHxcTEfCDy5Xm4cOFCQkJCwoX4+PiOCxcuJFy53vea999//+9FRAYH\nByM+/vjjG6Ojo7tDeYxfRbCP/aabbuq6Xs6nmJiYD3yBu//++3/V0tKyQETHXAwMDESuWbPmpeLi\n4t8UFBQ0iOg9N3xzcd999+33zUVYnhuBvq41MDAQMWfOnD+7XC7rZ599dsNUeXP3008//eYnn3wy\n3TRN6e/v/9Ydd9xx/PDhw/mPPPLIz3zX6Xbs2FF+5ZtYn3322Q3nz59PmjNnzp99b9wsWLDg9RMn\nTiwcGhoyrnzjxnedrq6ubl04vblrml++rh2KY+/q6opOSko639PTM7O7u3uW73G4zcXFixfjfI+f\nfPLJ7xcVFT2nYS6GhoaM4uLiX2/evHnnyPUaz42rzUU4nhtBmYBDhw6tTEtLezc5Ofm97du3/2Ay\nT8xALefPn0/KyclpzcnJac3MzHzbd1xdXV3RS5cubR7tY2s//elPf5icnPxeenr6O01NTct9630f\n1UpOTn5v06ZNu3zrL1269PV77rnngO+jWi6XyzrZx+1b1q1bVxcXF3cxMjLSm5CQ0L5v3771oTr2\nffv2rU9JSTmXkpJy7tlnny0Jt7moqakpLS4u/nV2dvZb8+bNe3P16tUNvjc2p/pcHDt2bLFhGEM5\nOTmtvo8rNjY2rtB4bow2F4cOHVoZjueGYZpqLr0BAIT/AxcAqEP4AUAZwg8AyhB+AFCG8AOAMoQf\nAJT5f852YFLmXI+FAAAAAElFTkSuQmCC\n",
       "text": [
        "<matplotlib.figure.Figure at 0x5124850>"
       ]
      }
     ],
     "prompt_number": 15
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "efficiency"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 16,
       "text": [
        "array([[ 109.45761821]])"
       ]
      }
     ],
     "prompt_number": 16
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "best_efficiency"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 17,
       "text": [
        "array([[ 227914.54247173]])"
       ]
      }
     ],
     "prompt_number": 17
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "plt.plot(best_design[:,0])"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 18,
       "text": [
        "[<matplotlib.lines.Line2D at 0x22cbccd0>]"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD9CAYAAAC1DKAUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXtYVNe5/9/hoijeABUjM4oREPCCGIwxtxKNUWyCJqYN\npjWexJNw7LFpeprG5Jym1fNrrLRNT9PY5tDmfqnS5oZNI4nGEG2NoIIaAyrooAMKKgIqigPD/P54\nz3I2w+w9+36ZeT/PwwMz7Jm19sze67vey3qXzev1AkEQBBF+RBjdAYIgCMIYSAAIgiDCFBIAgiCI\nMIUEgCAIIkwhASAIgghTSAAIgiDCFMUCUFZWtiA9Pf1wampqXVFR0epAx5SXl+dmZ2dXT5ky5VBu\nbm650jYJgiAI5diUrAPweDyRkyZNOrJt27Y7k5KSmmbOnLln48aNSzMyMmrZMe3t7SNuueWWf37y\nySfz7XZ747lz50aOHDnynCq9JwiCIGSjyAKorKy8MSUlpT45ObkhOjq6u6CgYFNpaeki7jF//vOf\nH1yyZMl7dru9EQCABn+CIAhzoEgAmpqakhwOh4s9ttvtjU1NTUncY+rq6lLPnz8ff8cdd3yek5Oz\n96233lqmpE2CIAhCHaKUvNhmswX1H3V3d0dXVVXN+Oyzz+Zevnx58OzZs7+86aabdqemptZJfS+C\nIAiiL16v1yb3tYosgKSkpCaXy+Vgj10ul4O5ehgOh8N11113fTpo0KArCQkJrbfffvuOAwcOZAV6\nP6/XG5I/P/vZzwzvA50fnR+dX+j9KEWRAOTk5Oytq6tLbWhoSHa73QNKSkoeyM/P38w9ZtGiRaX/\n+Mc/bvV4PJGXL18eXFFRMSszM7NGWbcJgiAIpShyAUVFRfVs2LBh1fz58z/xeDyRK1aseCUjI6O2\nuLi4EACgsLCwOD09/fCCBQvKpk2bdjAiIqL30Ucf/RMJAEEQhPEoSgNVE5vN5jVLX9SmvLwccnNz\nje6GZtD5WRs6P+tis9nAqyAGQAJAEARhUZQKAJWCIAiCCFNIAAiCIMIUEgCCIIgwhQSAIAgiTCEB\nIAiCCFNIAAiCIMIUEgCCIIgwhQSAIAgiTCEBIAiCCFNIAAiCIMIUEgCCIIgwhQSAIAgiTCEBIAiC\nCFNIAAiCIMIUEgCCIIgwhQSAIAgiTCEBIAiCCFNIAAgCALxegJYWo3tBmJkLFwCam43uhbqQABAE\nAJSVAdx+u9G9IMzM888DTJoE8O67+rftdAJ0d6v/viQABAEAv/sdwNmzRveCMDN79wI8/jjA6tUA\n3/sewNWr+rW9cCHAjh3qvy8JABH2HDkCsG8fmvi9vUb3hjAjXi9eI48+ClBVhWLwwQf6tH3qFMDh\nwwCtreq/NwkAEfZs2ADw2GMAsbEoAgThz6lTAB4PgMMBMHw4QE4OwPnz+rT9+ef4W4v2FAtAWVnZ\ngvT09MOpqal1RUVFq/3/X15enjt8+PCO7Ozs6uzs7Oqf//znP1HaJkGoxYULAO+8A7ByJUBcHEBb\nm9E9IszIvn0AN9wAYLPh42HD9JssbN8OEB+vjQBEKXmxx+OJXLVq1YZt27bdmZSU1DRz5sw9+fn5\nmzMyMmq5x33jG9/4YvPmzfnKukoQ6vP66wB33QWQlOQTgAkTjO4VYTaqqlAAGHoLwOLFJrQAKisr\nb0xJSalPTk5uiI6O7i4oKNhUWlq6yP84r9drU9IOQWjFSy8BfP/7+DdZAAQfzAJgDB2qjwA4nQBd\nXQA332xCC6CpqSnJ4XC42GO73d5YUVExi3uMzWbz7tq16+asrKwDSUlJTb/+9a+fzMzMrAn0fmvW\nrLn2d25uLuTm5irpHkEI0t0NUF8PcNNN+JgEgOBj3z6AF1/0PR42DODiRe3b3b4dYM4cgIQEFIDy\n8nIoLy9X7f0VCYDNZvMGO2bGjBlVLpfLMXjw4MtbtmzJW7x48YdHjx5NC3QsVwAIQmvOnAEYORIg\nMhIfkwAQgTh9GlM+x4/3PaeXC4gJQFwcCoD/xHjt2rWK3l+RCygpKanJ5XI52GOXy+Ww2+2N3GOG\nDh16cfDgwZcBAPLy8rZ0d3dHnz9/Pl5JuwShBi0tAGPG+B6TABCBYP5/G8eRrYcAeL0+AdAqCKxI\nAHJycvbW1dWlNjQ0JLvd7gElJSUP5Ofnb+Ye09LSkshiAJWVlTd6vV5bfHy8TglUBMFPSwtAYqLv\nMQkAEQh//z+APgJw+DBATAwmJZgyCygqKqpnw4YNq+bPn/+Jx+OJXLFixSsZGRm1xcXFhQAAhYWF\nxe++++79L7300sqoqKiewYMHX960aVOBOl0nCGU0N/cXgJMnjesPYU727QP47nf7PqeHAGzfDnDH\nHfg3EwCvt68lohSb1xvUja8LNpvNa5a+EOHB+vU44y8qwsebNuHqzpISY/tFmAu7HcswXH+977nT\npwFmzMDfWrFqFUBaGpafAAAYNAhXAw8e7DvGZrMpyrKklcBE2BLIAiAXEMGlpQXg8uX+a0P0sADO\nn8ckBYYWbiASACJsoSAwEYyaGoApU/q7XQYPxsygnh7t2m5txUGfQQJAECpCFgARjNOncZW4PzYb\nLgbTci1Aayvm/zPi49W/PkkAiLCFLIDQwePRppJrc3Pfa4SL1m6g8+fJAiD8+OADTA8jlONvAYwY\nAdDRQSWhrcizzwL8/vfqv+/p08YJQCALgAQgzHn5ZV95WEI+bjfApUt9Z1hRUejb1WOJP6EuJ05g\nWQ+1aW4GuO66wP/TUgC6uwE6O7ENhhYCoGgdAKE/LS20c5UanDkDMGoUQITfFIi5gYYPN6ZfhDzO\nnsWiaWpjlAuorQ2vRe71SRYAAc3NJABq4O/+YVAcwJqcOwfQ1KT++xolAP7uHwCyAMKe3l6cuZIA\nKMc/AMwgAbAmZ8/iKlm1MSoG4B8ABiABCHva2tA3SAKgnFCwAD79FOMWc+YY3RNj8XrxnvB4MC8/\nSqVRze3GAZ67GItLKFgA5AKyEM3N6BM8c8bonlifULAA/vpXgL//3eheGE9nJ5b0HjkSv1e14IsT\nMULBAiABsBAtLQApKWQBqIGQBSDmJtu9W/0+ScXpxPMId86excHfbgdobAx+vFiE/P8AZAEQOtPc\njMvSW1spV10pSiyA1laA2bOxRoyRkAAg587hTD0pSd1AsJD/H4AEgNCZlhYAhwMgNhagvd3o3lgb\n/70AGGIEoLoafyu5GVtbAT78UP7rPR4sXU0CgBbAqFHaWAB8awAA9HcBDRmCqa5ut3rtkABYCGaS\njhpFbiClKAkC79+Pv5UIwJdfAvz0p/Jf39gIMGAACQCAzwWUlBTaLiCbTf0YFQmAhWCzVhIA5Shx\nATELoLVVfvutrQANDfJTF51OgOxsXLV89ar8foQCzAVkt6vrAjJSAAJZAADqu4FIACwEWQDq0NWF\nmSNxcf3/J9YCmDBBuQvo4kX5szmnEzcoGT2assK0cgGZLQYAQAIQ1pAAqMOZM2hJBdpaL5gAXL6M\ng+9ttym3AACwho0cnE4UoTFjyA3EdQGpbQEEiwFoVTfKfy8ABglAGEMuIHXg8/8DBBeAQ4cAJk3C\ngVepBQCAbiA5OJ0AycnYDzVz362IfxaQWiuCjXYBkQVAXMPj8V3oJADK4PP/AwQvCV1djb73+Hjl\nFkBiojIBIAsAYS6g2FiAmBh1BkivV3iiAIBZOZcuqV+CoqsLVzTHxvb/n9qbwpAAWITWVhycoqNJ\nAJQidGNHR+Pm23ym/f79ANOn4+xMqQVwww3yBaChAQUgMZEEgLmAANSLA1y4gKuLhwzhPyYyEgWn\ns1N5e1xYADiQi5IsgDCFa46OGkWBPyUIWQAAwm4gNS0AuQJw9SoOenY7WQAAPssYQD0BCOb/Z2jh\nBuILAAOQAIQt3FkrWQDK4FsExuATAI8HYwBZWcZaACdO4EAXGUkC0N2NbpgRI/CxWoHgYP5/hlYC\nECgADGBCASgrK1uQnp5+ODU1ta6oqGg133F79uyZGRUV1fP+++/fp7TNcIQ7ax09mgRACcF8u3wC\ncPQofgfDhqEAKLUAZsyQJwDM/w9AAsAGS1awTU0LwCgB4AsAA5hMADweT+SqVas2lJWVLaipqcnc\nuHHj0tra2oxAx61evbpowYIFZV6vN4BniwiGvwvo3Dlt6p+HA2fPoojywScA+/ej+wdA2Y14+TJ+\nd3Y7/pZa1oMEwAcLADPUsgCCrQFghLULqLKy8saUlJT65OTkhujo6O6CgoJNpaWli/yPe/HFF79/\n//33vztq1Ciat8qE67aIicEyAFpuSB3KtLcHXgTG4BOA6moMAAP4YgByRJjd4DYbpnJKtQJIAHz4\nC0CoxAD0cgEp2jqhqakpyeFwuNhju93eWFFRMcv/mNLS0kXbt2+fs2fPnpk2m433llmzZs21v3Nz\ncyE3N1dJ90KK5mb0PTNYHID2rpVOW5vPZxwIIQvgBz/Av2NiMGOos1M4UyQQ3BkeEwAmLGJwOgHu\nvRf/HjIERejSJen9CAXOneu7YYta9YCam3G9RzD0dgHV1JRDU1M5cIZKRSgSAKHBnPHEE0/8dv36\n9U/bbDav1+u1CbmA1qh1ViGIv9+aCUBKinF9sirt7fIE4MgRgMxM32MWB1BDAKTAFoEBoBXBrIBw\nvBYCWQChEAROTQ38v29+Mxe6u3Ph2WcxCWDt2rWK2lIkAElJSU0ul8vBHrtcLofdbu+jv/v27buh\noKBgEwDAuXPnRm7ZsiUvOjq6Oz8/f7OStsMN/9RFygSSh8eDs/Zhw/iPiYvrP4h0dfnKcTOYOT5+\nvLQ+cAVg/HhlLiAAEgCuBRAXh+WSlVpERsYAhCyAyEiAn/0MF4pFRipvS1EMICcnZ29dXV1qQ0ND\nstvtHlBSUvKA/8B+/Pjx651O5wSn0znh/vvvf/ell15aSYO/dPxnJCQA8rhwAWDoUP5t/gDQ0jp9\nuu9zx4/jrJu736zcTCAlFsDFixhE5lqD4RwH4K4BAECLSI1AsFljAAAAzz4LMHCgOm0pEoCoqKie\nDRs2rJo/f/4nmZmZNQ888EBJRkZGbXFxcWFxcXGhOl0kurvRbcGdFYSqANxxB0B9vXbvH8z9A4C+\n38OH+z5XX99/hi03IOcvAFIKwjU04Gu4q0TDWQD8XUAAaFU5nfLfs6cHv1e+zeC56J0FpDaKXEAA\nAHl5eVvy8vK2cJ8rLCwsDnTsa6+99rDS9sIRZuZyTb5Ro/S/6Z97DmD16r6zYLX56iuAf/xDO3dG\nsAAwAEBaGg60bjdmWwEEFgAlFoDdjn9LtQBOnPD5/xnhLAD+QWAAjNN8/TXAggXy3rOtDZMrxFzn\neruA1IZWAluAQCtX9bYA3G6An/wEF0NpRW8v3nx792rXhhgLYOBAgHHj+loiWlkA8fE44xS7FiDQ\ntZCYGL4VQQNZAJMnowDIRcoMXG0B8HqDu4DUhATAAgTKSNC7HhAbYGpqtGuDVeHUWgCE1gAwMjIA\namt9j+vq1LUA2ADD1gKIdQMFGvDC2QII9HlMmWJdAejsxPTimBj13lMIEgALwCcAeloAbIBRcmMF\ng5VI/uordTe+5iLGAgBAAeCKnVYWAIC0TKAzZ/qvYg5XAWCzZf/BevJk/O74SnoHQ8oMXG0B0HP2\nD0ACYAkCXeR61wNqbsbZqpYWQGsrul6Sk7UTGjExAAD0IzML4OpVzAry973LtQDOnev7fUqJA5AF\n4KOjA0t3+2fEjBiBA/PJk/Le10gLQM8AMAAJgCW4cKH/il9mAehVD+j0aSxeprUFkJAAMHMmwJ49\n2rQhxQJgAtDQgPn//kFBtSyAiRPFx1YCCQCLAcid8VqVQJ8FQ4kbyEgB4NsMXitIACzAhQv9Fy7F\nxuKMXO3NKPhobgbIzQU4dgzTUrWA3Xg5OdrFAcTGANLTcVDu7Q3s/weQZwH09GAuP7cPM2YAVFWJ\ne32gQW/gQFz0pOZOUVYgUAYQQ0kgWIoADB2K96daEzESAKIfgQQAQN84QHMzrj6127XL02f+z5wc\n4y2AoUOxLydOBPb/A8gTAJZiyE3pzc4GOHgQxSEYfLPecHQDCVkAegnAwIG4qPDqVXlt+XPxovAq\ndbUhAbAAbPWqP8OH61cRlC2NZznWWsDyn7OysO7OlSvqtyFWAADwXGtq+AWA1QyS4noJNLgMG4Yu\npmDxFa9XWABOnRLfj1BA6LucPBk375GD1Dx8Nd1AFy8Gvte1ggTAAvDNCmJj9XUBXXedL8NCC9jg\nOGgQrsY9cED9NsQGgQF8cYD6+sDFuaKj8TuQcvPzzS7FuL06O3G2GWiz8ECrl0OdS5f4B8vMTPw8\n5MRFpAZimRtIDS5e1LeqKwmABeBzAQ0ZgjeBHrBUVC0tAO6NN3OmNnEAKRYAEwC+GACA9ECwkADs\n2yf82jNn+F0eU6di+mw40dkZWAwB8H4ZOVJeSQipAhAbi/WZ1IAsAKIfRguA14suoMREfSwAAO3i\nAGKDwAAoAAcOYH15/xRQhtQ4gBILIFjWi1yXh1UJVvFTrhtIai7+4MHqWeJCVo0WkABYAKMFoKMD\na+LExqKrob5em0wgfwtAKwGQEgOoqsLqkqwmkD9qWQDTp+NgJfS5Bgt6HjoUXtuEBhMAOamgfIvL\nhCALgNAUPgHQKwbALY07eDAOiMeOqd8O98abMgVz210u4ddIobsb6/rzuQ38GTkS+yNUmE4tC2DI\nkOAL4IQEICEB30PNz8vsiLEApArA5cuYXj14sPjXqHkfUgyA6IPQdn96WQD+m2MoLbbFB9f0jo4G\nWLQI4K9/Ve/92ezfxrsnXX8yMvh3ZwJQzwIACO4GEhIAgPCLAwjFAADkXadyVuIOHkwWAKERnZ2Y\nFRNo9x+9BMC/FhFLj1STq1ex/g/34v/2twH+8hf12pDi/2fMmIHWCB9qWQAA4gTAvw4Ql3CLAwSz\nADIyMJ1YzPoKhhwBUNsCIAEgrsHn/gHQ1wWktQXAVkByZ+dz56KrSeqWiXxI8f8znn8eoFBgayM1\nLYAbbhDOBApmAZAA9CU2Fr8fKesjjLYAKAhM9EFIAPR0AXG3x5syBWD/fnXbCHTjRUcD3Huvem4g\nOQIQGSm8faSaFsD06SisfKtKhdJAAcLPBSRm31+pO67J2YyFYgCEZphBAPwtgGnT8EJVc+ER38D4\nwAMAJSXqtCFlEZhY1LQABg/GgDPfLD6YBSDH5WFlxAjA+PHSBMBoC4BcQEQfzCgAEREA99+vboCW\n78b7xjcws0WNrCM5FkAwpFgAYlIMs7L4Z/HBBCA2FjO0tNxT2UwECwIDSNtrAUBePX61LACvlwSA\n8MOMMQAAgG99S90ALd/AGBUFsGSJOm3JCQIHY9gw8SJ8+TKK56BB/Mdcfz3/6tVgAgAQXnEALVxA\nRloAXV3o9oyOVv5eYiEBMDl8heAAjIsBAADcfDO6PrjbJipBaOa1aBHAJ58ob0MLCyA2Vvx3ICbA\nd/31AMeP93/+8mWsaxNswAunOIBYF5BUC8CoGIDes38AEgDTI1QeVg8B6O7GgdO/7rrabiChG2/i\nRGmzOD60iAFI+Q7EDFh8AsBm/8HWMISLBdDTg2nDwfbOtVIMQO8AMAAJgOkxOgZw5gwO/oHWIXz7\n2+oJgFD2hd2OqXxKd7zSwgLQSwCCZQAxwkUAmP8/mCCOH48xJLHXDlkAEikrK1uQnp5+ODU1ta6o\nqGi1//9LS0sXZWVlHcjOzq6+4YYb9m3fvn2O0jbDCaNjAIHcP4zZs3FWrYYbSOjGi4nBvQ9aWpS1\noUUMYMAAHFzEbGLf2RlcAK67DvvpP6MU4/8HAEhLw71wtdhLwUyI+SwBcHY+dKj4a8dIC0DvNQAA\nCgXA4/FErlq1akNZWdmCmpqazI0bNy6tra3N4B5z5513bjtw4EBWdXV19uuvv/4vjz322B+VdTm8\nMNoCCBQAZqjpBgp2440bp7zOjRYWgM2G34MYIb50KXjWSkQEBi79A8FiBSA6Gme9gayIUEKMNcWQ\n4gYych2A5SyAysrKG1NSUuqTk5MboqOjuwsKCjaVlpYu4h4TGxt77aO5dOnSkJEjR55T0ma4ISQA\ngwZh5oDHo137QgIAAHD77erU7Q8mAA6HcgHQIgYAIF6IxQ5agdxAYgWA7/WhhhQBEJsJ5PFg5Vup\nVqKVYwBRSl7c1NSU5HA4rt2Wdru9saKiYpb/cR9++OHiZ5555henT5++7tNPP72L7/3WrFlz7e/c\n3FzIzc1V0r2QQCgLiO0O1dmp3T6i/oXg/Bk/Hl0OStFDALSwAADMJwATJ2pTrdVMSLUAxGQCtbfj\nvRYo3iWEnhZAeXk5lJeXK2/s/1AkADabTVT18cWLF3+4ePHiD3fu3HnbsmXL3jpy5MikQMdxBYBA\ngm0SrbUANDcDpKfz/19qlkUgvF5fLSA+HA7lQqOlAIgZAMQsXALgF4C0NHH9CRcBEFvWe/x4cXEq\nOf5/AHUtgGAC4D8xXrt2raI2FbmAkpKSmlwul4M9drlcDrvd3sh3/G233bazp6cnqrW1VcbHHJ4I\nuYAAtI8DBHMBJSRgAFTJnqgXLmCgl2/TFQDlFkBXF/4WWoQlF7NZAOHgAhIbBAYQ7wKSKwBqWQCW\nCwLn5OTsraurS21oaEh2u90DSkpKHsjPz9/MPebYsWMTvV6vDQCgqqpqBgBAQkKChPJZwfF4QveC\nN1oA6utxBsWHzabcDSQm8KZUALSa/QOIXwymRADEpoEChI8FoLYLSK4AxMTgBENpmrLlYgBRUVE9\nGzZsWDV//vxPPB5P5IoVK17JyMioLS4uLgQAKCwsLH7vvfeWvPnmmw9FR0d3Dxky5NKmTZsK1Om6\nj+3bsWTvsWPSNvuwAsEEQMtU0AsXUACmTxc+btw4nGEJ1c0XQsyNp1QAtAoAA0izAPhSarlMmIBZ\nQF6v73qWYgFMmIADnscj3Z9tFeRkAXE/z0DIFQBW3uPKFfFuqUBcvIjXuZ4oEgAAgLy8vC15eXlb\nuM8VFhYWs7+feuqpXz711FO/VNqOEHV1eMMcPYp71oYSRloAlZUA2dnCrhkA5XEAMTfe2LE4CHZ3\ny6uVoqUFIPY7EBsDGDIEf7hbcUoRgMGDfXXw9R5QGG++iVbhT36izftLiQEMH47XTGtr/xXtXOQK\nAIBvIqZUACzlAjILx47hF1xWZnRP1OXqVZzFDRzIf4yWAvDll1jzJxhKXUBibryoKIDERGmbe3DR\nYhEYQ+0YAEBfN1BXF8ZZpAT6jXYDVVQAvP66dpvUS/ksAcRNUpQIgBqBYMvFAMzCsWNYnXLLluDH\nWgmWASRktmopALt2iRMA5gKSi9gbT4kbyAwWgFwB2L0b9yWW4t40WgCcTmxf7a1DGVKCwADiBEDO\nIjCGGq5YsgBkcuwYxgD++U/1NmYwA8FSQAHUiwEcOgTwu9/5Hvf24sAze3bw16rhAhJTg11JKqhZ\nYgByBKCoCOCJJ6T1yehMoIYGgNxcgM2bgx0pDzkWQLBAsJy9ABhqWABUDE4GXi9e6NnZuIH3F18Y\n3SP1COb/B1DPAnjjDYBnnsGZMgDu9hUfj26XYOgRAwCwvgUgxUfMBvD9+7G887Jl0vpkpAXg9eJg\n+4MfAJSWatOGlBgAgLhUUDViAEogC0AGzc344Q8dCrBgQWi5gfQUgK1bcTept9/Gx19+KW72D+AL\n0IopiBaIc+eEg3MMJQKgZHYXDCm1gKRYAE4nwPr1AD/8oXAcKBBGCkBLC96TCxfiFpWnT6vfhhYx\ngJYW8YF2fygGYBDHjuHFDgCQl0cCIIczZ3DG9vvfAxQX4wxOrP8fAAO0110H0Mi7BFAYoYqjXJQU\nhBPbhhzUXgcAgAJQXQ3w2WcAjz0mvU9GuoAaGnDGPWAATso++kj9NqTGANLTAXbsAPjNb3xWrv/7\nOZ24r7IcyAIwCK4AZGXhTRYqe6KKEQA1LrzPPsO9d++8E2fxu3aJzwBiKHEDiR2clVgAWgqAFjGA\npCTMAlu5Ut6gMGoUfpeBBjutcTpxLQIAQH6+Nm4gqRbA5Mk4Ody7F/v26qt9/79vH+6mJtXSYii1\nANh+wBQDkAhXAGw2nHGESjqoUCE4hhoWwNatAPPm4edXWAiwbh3O5qUs7FKSChouAiAlBhAZiTn0\njz8ur082m3FuIK4A5OXhzFvtxYpSBQAAYNYsgD//GeC99wBefLHv/yoq8P9yUToR6+rC7zzYmhu1\nCSkBAACYOTN0dkQSkwWkVAC8Xp8AAAAsX44WwcyZ6NoRi9xU0M5OXNw1fHjwY0eNws9EzmYnRguA\nx4M3uZRaRD/9qbjYCB9GuYG4AjBiBEBODoCKBSwBQHoQmMutt+K40copSKNUAJRaAEa4fwBCUAAS\nEvp+sVZGjxjAkSM4W2SVJhMSAB58EFP4pCDXBcQGZjE57hER6BqRagVcvozuECOzgNjsP0LHO84o\nC4DFALj9kBsf4kOOBcAYMABFgCtKFRUAN94ovz9KLQAjAsAAJACmRo8YANf9wyguxpRQKch1AUmd\nmctxA0kRGTlIEQA94QpAS4s6+zaIgWsBAGhzT0oNAvszZw7WEAPA1eVXrvQdR6SihgWgt/8fwOIC\ncOECfnHcXPX4eFzRFwroYQFs24bBXy7R0dLcPwDKLQCxyMkEOnVKO/cPgLjvQMmMVS7XX48z24ce\nwhn5v/2b9m16PPj9cCvIqi0ALN1Yib987lx0dQL4Zv9KJghKJ2LkApLBsWN4kXO/uHCzAOQKQEcH\nZvp88QXeDEphM3OpJXH1tAC0gs3+hOreGCEAU6diYHHqVFwlv3+/9m2eOoVxi5gY33Nq35NK/P+M\nrCxcu9LUpNz/D2DdGIDiaqBGcuwYQEpK3+dCTQCCXRRyZh5z5uBFn5GBWSajR8vvI7cfQ4bgmgKh\nDWT8kTo4jx4t3a+ttQBERuKAd+UKDgSBMEIArrsO1xIAoDhdviytqqgcnM6+/n8AbQRA6WcZEYFx\nrs8/x3th9Wpl72fVGIDlBcDfbxcbi2bolSva7P6kJ1pYAFeu4Gyws1O6mycYLA4gVQCklPAeMUJ6\nbrvWAgCg58DeAAAgAElEQVTgWwzGJwBGxAC42Gw46z1woL/LT038/f8A5hQAALR8P/0U1wAoCQAD\nWNcCsLwLyF8AbLbQsQLEpIHGxGAaZXe3uPc8cQL96GoP/gDyUkGlDs5xcVjYTcs25BBMiI2wAPzJ\nytLeDaSHACgNADPmzAH4619xwqK0TIgaMQAKAkskkAAAhE4gWIwFYLOJr0UDEPgGVQu7HX2qUpA6\nOJvVArCCAEyfjhaAljQ0WMcCmDQJxwql/n8AtAAoCKwzfAIQKhaAGAEAkDb7OH4cA+daoMfsnCwA\n+ehlAfjHAOLjUbSV7pnLUCMIDICTp7w8LIOilNhYZS4gigFIxO3GSqDjxvX/XygIQG+v+EFDShxA\nSwsgPl5agNbtRpGTstpVjgCcOoUVS7Uk2HdgdAwAAOvh1NdjjSG5NW+CEej6iorCz6e9XZ2KrGqK\n6Z/+pM77qGEB2O3q9EUKlrUATpzADyyQLzsUBKCzE4PYYjb1NosASB2cm5sxq0fK6lipLqCrV/Hm\nklvnXSxWsABiYtBi1mqXLrcbF5wF2odYzXtSzc/SZlNngaBSC4BiABIRcmWEggCIdf8ASBMAM7mA\n5LhmYmNxoBG790BzMy4U1LoEgxUEAADjAFq5gVwutLS0npSpFQRWE1oIpjPBBMDqQWAxGUAMKRef\nmSwAOQJgs0mzAvTw/wMED8SbRQBYKqgW+NcA4qK2BWC0O80fMYsBhSABkIiQAMTHW98C6OhQ3wJo\na8PYglY7Y+khAFLb0UsAgm0KY4YYAIC2FsCJE31LQHAxqwtILSIjsYTK1avyXm/ZYnBlZWUL0tPT\nD6emptYVFRX1W0/3zjvvfCcrK+vAtGnTDt5yyy3/PHjw4DSlbQKEhgvoyhXc8i8hAeDrr/v+T0rg\nUqwAsM9Mq6JocXHSLK9QEgCruICYBSB3pioEW2MSiFAXAABli8EsGQPweDyRq1at2lBWVragpqYm\nc+PGjUtra2v7bKp2/fXXH9+xY8ftBw8enPbss8/+v8cee+yPyrqMWF0APvsMSzDv3Yv5yIcP9/3/\nyZOBg2mBECsAWrp/AHwDs9jBRe7gbFYXkBUEYPRoDAZrURk0nC0AAGVxAEu6gCorK29MSUmpT05O\nboiOju4uKCjYVFpauoh7zOzZs78cPnx4BwDArFmzKhobGxUnO3m91heA//1fLLn87ru4DN1/Be3J\nk/yzKX/EXnhaC8CgQWhdiN2whSwAY9BqQdjJk/oIgBmDwADKLQDLrQNoampKcjgc12oz2u32xoqK\nCt51da+88sqKhQsXfsz3/zVr1lz7Ozc3F3J5diU5fx59bnwbfFghCHz8OMBTT+HfgUoouFy4K5cY\npLiAJk+W1k+psMGZryYOFz0sAD3WAABYYx0AIyMDNwJSGz1dQGb5LLnItQC8XvExgPLycihXcXs1\nRQJgs9lEexI///zzO1599dVH/vnPf97CdwxXAIQIlsoYH+9zRWjl71YK9xzGjwfYubPv/6VYAEOG\nYJXHYDidAHffLa2fUmECkJQU/FiyAIxh9GiAc+fUfc/eXtz1i2IA0l939SqmKYvZ38B/Yrx27Vrp\nDXJQ5AJKSkpqcrlc1zzVLpfLYbfb+23+dvDgwWmPPvronzZv3pwfFxcncR1nf4IJQHQ0uiM6OpS2\npA3+2Th8FoDVXEAAPvENhseDosXdzEcsFANQRkKC+gLQ3IzCzN0HwL/NUBcAuRaAUQFgAIUCkJOT\ns7euri61oaEh2e12DygpKXkgPz9/M/eYkydPjrvvvvvef/vtt7+bkpJSr6y7CNsIRgi1LrglS3BR\nllx6egB27+77nH82jv92im43Do5iBy4xLqDeXhQZvjxttRA7Oz97Fo+NjtaujZ4evAbU2O8gGFZy\nAY0cqb4ACLl/ACgGIISUlG+1USQAUVFRPRs2bFg1f/78TzIzM2seeOCBkoyMjNri4uLC4uLiQgCA\n//7v//5pW1tb3MqVK1/Kzs6uvvHGGyuVdlrMalY1Lrj2doD331fmL925E+Dee/s+59//UaPwomaz\nh6YmHPzFlIEAECcAp07hwCnGN68EsYOzkpm52DbOnMHrQIvS1/4ILQRj5br5Zsd6M3Kk+kkSQgFg\nAIoBCHHmjDxLWA0U3xp5eXlb8vLytnCfKywsLGZ/v/zyy//68ssv/6vSdrgcPw5QUCB8jBqB4Npa\n/O10ig/I+nPwIJrH5875ip75C4DNhrOnkycxQOdyiU8BBRAnAHq4fwD0EQCxLiC93D8AwgvB2IzV\nLPEorSwAIQFQc6Mms7qA5FoALS3GCYAlVwLrZQFwBUAuBw/i76++8j0XqP/cOICUADCAuJmHljWA\nuJjJAtB6M3guQiJstgHLCBeQWhs1eb3mcqdxkWsBtLTo46YMhOUEwO3GwSPYDFmNchC1tZhCqFQA\nsrKCCwA3DhAOFsDx4/LjEWItgLo6gNRUeW1IReg7MNuAFReHn5/Ho957BnMBAYgTgKYmFG4+rlzB\nbBk93HpSkWsBGOkCspwAnDyJKYbBgodqWQALF8oXAI8HS+8uXeqzBADUtwDMJgBiXG9HjgCkp8tv\nQ4zIHD4svw2pDByIQedAW3OazQKIisKgo9R9FYQI5gICEHdPrl4N8D//w/9/swaAAZRZACQAIhHr\nyjCDANTX436jN9/sswB6ejBf2v9m4VoAcgTg4kXhY2pq9BkMpQzOUjaD5zJ8OGZmBdthSk8BENqa\n02wCAKB+IDiYCwgg+D3pdgN89BFawHxcvGgua4qL3E1hSAAkIEUAlASBr1xBc3TuXByQ5Wxnd/Ag\nwLRpAFOmYLG33l68uMeM6b/og2sBSHUBDR+Og0xPT+D/d3dj+1lZ0s9BKmIEwOPBVN60NHltREXh\nzRZM9PQUAAB+S8ysAqBWHKCjA6/tuDjh44IJwBdf4LUqVKeorS14O0Yhd1MYEgAJ6GUBHD2K7Qwb\nhheckF+SDyYAcXH443Ty91+JBRAZiSLAN/DW1OD76zFzEiMAJ05g0EtJSmqwdlpbcYXlmDHy25AK\nnwCYLQYAoK4AsNl/sCynYPfkBx8ALF8uLADnz2u/u5tcyALQAb0EoLYWIDMT/54wQZ4biAkAAMDU\nqegG4uu/3Y4i09aGsyCpsxyh862uBpgxQ9r7yUXMSmAl7h9GsEAwizHomXppNReQmgIQzP8PIHyN\n9vYClJYC/Pu/Y1A0UCwFAF+v1X4WSlFiAVAWkEiEdh3iojQLqLYWc/IB9BGAAQNwQVhFhbjZlD9C\nPt2qKv0EQExJaCUBYP92+NDb/QPAvxbAjAKgZjkIMRlArE2+a3TPHrS2J09Gq43P4jazBSAnCNzV\nhe5mo9xalhOA1lbfgioh1LAAlAhARweWO2CDfTABAMCB/x//kOb+YQjd0NXVANnZ0t9TDjExWNhK\nqCS0HhaAEQJgtRiAWkFgMQFgAOF78sMPfSvmHQ5+N5CZLQA5aaBnzuDs36hFgpYTgPZ2/jLQXIYP\nx0GIz5QMhlIBOHQIg7+snIMYARg/HgVASgCYwXdz9fbiFoB6CQBA8Nl5qFoA4RoDUMMC+OADnwCw\nVfGBCDULwEj/P4DFBKC3FzM/hg8PfqzNJn2LQkZPD6ZwslnqhAk4cEvhwAGf+wcAByOn0xdcDsS4\ncT4XkFT4ZnR1dfg/PU1MMYOzUgvASgJgVgvALDGAw4fxM7rhBnwsJAChZgGQAEjgwgVUWbFF0uS6\ngZxO9EOyLBU5FgDX/w+AC4UmTkQR43NhjR+PPkG5FkCgG1rPADBDaHDu6MCbXcx+AUIIuYCuXsUB\nZOJEZW1IJZwFQIkLaP9+XCsT8X+j0bhx/GsByAJQF0sJgFj3D0NuIJjr/gHAAbmlBQcWsfgLAAC6\ngYQ2ZWc3kdwYQKBz1TMAzBCyvI4cwdm/Up+nkMgcO4ZiKmaDDTUREgCzuYDUCgK73fg+YnZdi4vD\nCYB/CYpTp/pOCKwaAxCzINMfI8tAAIS4AMi1APwFICoK0zSlbKT99df9t1+cNk04hZWZ0XIsAD4X\nkJ4BYIbQ4KyG+wdA2AIwwv0DIBwDMKMFoEYQuL0dXbJirPKoKNz20P978y8MaNUYwLBh0uOORqaA\nAoS4AIwaJW6rRH/8BQBAmhvo6lX0BfrPVJYtA3j6af7XsZm/Wi4gr9c4C4BPANQIAAdro7bWOAEI\nNAM0owuIzcb5Vo+LpaNDXEyOMXo0znq5+O/bbNUYgJi443e+A7B9u+8xuYAk0N4uLZiZmIgfsFSa\nmvq7YaQIAOunv5vD4QCYNYv/dSNGYD60nHrpgaydkyfRDaLnalgAfSyAYG34C7ge8N38ZhSAyEi8\n3pQWhJMqAHY71sLi4l+2Oy4Ohcl/J77eXry3zCoAAMJeB7cb01137PA9RwIggbY2aRZAYmL/2YYY\nAlkaUgRASb2SnBx5rwtk0hsRAAbQxwIwowuIz+K8cAFdH2ZDjUCwVAFwOPoHeE+f7msB2GyBj2NJ\nIGYsBc2Ij+e3ACor0UW0b5/vORIACUh1Acm1ADo6jBMAubALj7sCl+1FoDd85SBYETg1avTziYzX\nq56VIRU+ATh7Fv9nNswiAIE27gnkBmptNa//nyFkAWzfDnDffeiWZZAASECqAIweLU8AWGCLi9kF\nYMAAdB11dPiea2jQZw8Af/gG54YG5UXgGHwWQHMzfg5GLK0PJABuN7qAzFjBUo2S6XJcQFwBYFVs\n/d8jkACcP29u9w9AcAF45BFM9W5uxvPu6DBW1CwnAFrHALzewBe12QUAoL8bSOwCHbXhE4ATJ9QT\npEGD0Cfc1dX3+dOnla8xkMuoUTij5lphZ8/i9xJhwjtNDQsg0GRJCH8LgGUA+cfLAq0FMHMAmMEn\nAJcvA+zdC3DbbeiWrarCayM+Xvy6Ji0w4WXJjx4uoK4uvBhjYvo+P3o0ZniIWeknNVahFv6ZQGYT\nADXNXZstsBXQ3GycSR0Tg5YYN3hpdJ63EGZwAfn7/7nHBbIAzO4C4osB7NqF7tihQ1EA9u0z3v0D\nYDEBkDqwxsfjoO12i38N3wUdESG8QMW/n0ZYANzZR28vZlvIWVSmFCEBUDMjKVA7archFX83ECv2\nZUbUEgAp9yQTAGYlBfL/A/DHAKxqAWzfDjBnDv59ww1oAZAASESqBRARgTeklEwgoTa4m7YIYQYX\nUEsLCpmclFKlsHRI/5LQal/wgQSguZkEQCxGWADsWGYl8VkAgVxAVrAAxAgAcwGFhACUlZUtSE9P\nP5yamlpXVFS02v//hw8fTp89e/aXMTExXc8///yPlLQlVQAApLuBhC7o8eN92zYKYaQFwG5oo9w/\nAOgKiYzs7y5T2z0TKI/dSBcQgC8OwDCzABgRBPZP8eSzANh6Ae5WrFa1ADo6sDrw7Nn4eOJEfK6m\nxuIC4PF4IletWrWhrKxsQU1NTebGjRuX1tbW9lmCk5CQ0Priiy9+/8knn/y1sq5KDwIDSBeAYBaA\nGAGQ00814F58RgoAAL97Rs0LPtB3Sy4g8RhhAQD0FQA+CyAmBu9D7vdrBQsgUAxg505cAMriijYb\nlmfZssXiAlBZWXljSkpKfXJyckN0dHR3QUHBptLS0kXcY0aNGnU2Jydnb3R0tMzK/D6MtgC4G7cL\nYQYXkNgKjVqhhwAEyiknF5B4zCAAfBYAQP84gFUtgIoKgFtv7fvcjBm4TsdoAVC0pq6pqSnJ4XBc\nuwXtdntjRUWFQLEDYdasWXPt79zcXMjNzb32uKcHXQpSl9TLEQCrxgC4F9/JkwBpafr3gcE+96lT\nfc+pPTt3OHB1JRczuIC41xsJQH+4awH4LAAALJxYX+8rn2IFC4Ddg16vL7X15EkAzlAGAL69D6Re\nG+Xl5VBeXq60m9dQJAA2m01g51fpcAXAn44OrLYnNZ86MbF/7REhhPKarRYDmDdP/z4w/Gdvvb3q\nD4Z2O8B77/V9zgwuoEOHfI/NLAAjRmCWXE+P/PIKci2AnTvxbyELYMqUvp+lFRaCDRqEY9Tly74S\n4I2N/Qs8shItUicr/hPjtWvXyu8sKHQBJSUlNblcrmun5nK5HHa7XcJwKx457h8AdV1AdjvOWIJV\nUDRqHYC/C8jIGIC/ALS14Q0xcKB6bTgcfcW9qwtvPCNX3VrJBRQRIX/XPIYSF1BnJ6Zo890r/gJg\nhVIQAP33IXG5cOzgkpaGx4nZR0FLFAlATk7O3rq6utSGhoZkt9s9oKSk5IH8/PzNgY71er2KtgCR\nG1hVMwg8YAAOsqdP87++uxsLPhlR/MtMQWB/AdAi5c0/BsBqqxu1wTZAXwHwelEAzFgHiKHEDdTV\nhefov2gyGEy4+VYBM6ZO9QmAx4Opo0ZMrKSSkOATVa8Xz9VfACIi0L1lpLUKoNAFFBUV1bNhw4ZV\n8+fP/8Tj8USuWLHilYyMjNri4uJCAIDCwsLi5ubmMTNnztxz4cKFYREREb0vvPDCD2pqajKHDBkS\nYOsMfpRYAFLWAQSb0TA3EF/NfuZCMmLpP3MBtbfjhWfkzTJuHMCmTb7HWgjA8OF9S3cY7f4B6CsA\nly5hOqzZdgPjwmIWmZnSX8s+d6mCy4Tbfx8AfyZMwHuXLeYcNszYsgli4U7Ezp9HqzdQ7NIM9aEU\nF1bNy8vbkpeXt4X7XGFhYTH7e8yYMc1cN5Fc5LpVpBaECyY0TAD8o/oMo/z/AFhkzWbDDVHGjTN2\nJqyHBWCz+QKKw4cbnwEE0FcAzOz+YSQl4UAsBznuHwAcDAcOxNm9kABERuK+Dl9/jYOq2f3/DK4A\nBPL/mwnLrASWawGMHImDstidj8RaAHwYKQAAeL5VVca6fwD6L/nXanDmuoGMzgACwNm+x4OxCKsI\nQFOTvNfKFQAAXwYXXwCYweIAVkgBZXBjAIH8/2Yi5AUgKgoHZLF+zmDtCG1Xx15vpAAkJJhDAGJj\n8YfNhrVa9s4NBJvBBWSz+awAKwjA2LHGCUBFRfAgKBMAK6SAMrgxALIAVELJwColEGx1C8AsAgDQ\nVyy1FAAzWQAAPgE4e9b8AmCUBWC348Y9oWgBcF1AZAGohFwLAEC6AIiJAfBhtACMHIk3TLgKgNEW\nAIC1LAAjXUAAwS2AqVMBvvrKehYAVwDIAlABJbn1YgWgtxczDoRSOFk5CP9Kl2r0Uw0SEjDeES4C\nwF1VagYXEID1BEDvIDCAb1AMZgGMHYsZQIcPW8cC4MYAAqWAmgnLCIAeFsDFi+i3Fko1GzYM1wPw\nLZ4x2gJgsyQj6wAxuAIQLkFgAGsJwNix+Llxq26KRepeAFzEWgA2G7qBduywlgXAxgeyAFRCDwEQ\n24aQG8hoARg5EiA6OvjMSg+YALAFUVoMhiwI7PWSC0gOAwfipCbQZvbBkLodJBeHA9sWc79NmYKl\nk61iAXDrAZEFoBJ6BIHFmrRmFoCEBLy5zLAHLdvUo60N1yhIXTEqhqFDUfBY7XgjVmD7YyUBAJAf\nB1DiArr+eoC33hK3VmXKFPxtJQugtRV/Bg0y90JAEwwT4jCbBcCXCmq0ADgcAJMmGdc+F2YBaL3z\nkd0OsGcPtmHk4jcGCUBwIiMBvvUtcceyirJWsQDi4nAsOXnS3LN/ABVWAuuFHkFgsRe00L4ARq8D\nuP12385DRjNmDPpCT57UVgAcDoC9e83h/gHwlVdobUWXnNmRuxZAiQBIYfJk/G0VCyAqClc7Hzpk\nbv8/gEUsgKtXscja4MHyXq+2BcDqlAfCaAvAZsMgtRmIiMDBRevB2eFAC8BMAnDkCA6Ocsss64kR\nFoAUEhIAli0zR1xLLAkJAAcOmN8CsIQAsGwDueY926c1WKaD2As6OxsXWwXCaAEwG+PG+dwzWsEs\nADNkAAHg9dbRYQ33D4D5BQAA4M035U8AjYAJAFkAKqDUrTJgAAYHg9U9F2sBJCdjKVz/stAeD6aS\n6nVTWAG9BKC93TwWwIgROPMnAQhf4uMB9u8nC0AV1FhcNW4cgNMpfIzYC9pmA5g5Ewc2/9cPHWqO\nDByzMG4cLjTSOggMYB4BsNnQ908CEL6wTCCyAFRASQYQY8oULCsrhJQLOienvwCQ+6c/bEGa1haA\n1m1IZdQoawmA1NXAXV0odFqk9oYCLGBNFoAKqCUA3O3llLYzcyb6nbmQAPRHTwEwiwUAYC0BSEjA\n8tVXroh/jZJFYOEACYCK6CUAciwAbk0gEoD+MAHQcnAePBgHXKP3V+UyerR1BMBmk54KSu4fYeLj\ncSww8yIwAIusA1Ajt15tC2DsWFzK3tCAW9cBkAAEggmA1oNhVZW5Zls//7l18tYBfHGAlBRxx5MA\nCMNW5JsdS1gAagWBL1zA9+JD6kXt7wYyehGYGRk6FD+jQYO0bcdMgz8AwMSJ1tjAnEEWgLpMngww\nf77RvQiOJQTg1Cnly8AjIvBLEbICpLqa/DOByAIIzA03GN0DIhhSM4FIAITJygL45S+N7kVwTC8A\nbW0AH30EcM89yt8rmBtI6kXtnwlEAkBYFRKA8MT0AlBcDHD33eoE+IQE4OpVXMglxVWRk4O+Z7bC\n2OjNYAhCLiQA4YmpBcDtBnjxRYD/+A913k9IANgFLaXcREICLvg5ehQfkwVAWBU5AkCTHeujWADK\nysoWpKenH05NTa0rKipaHeiYxx9//Hepqal1WVlZB6qrq7PFvndJCUB6OsD06Up7iUyZgvuLBtrO\nUe6MJicH4NlnAZ56Ct1BJACEFZG6GIwsgNBAkQB4PJ7IVatWbSgrK1tQU1OTuXHjxqW1tbUZ3GM+\n/vjjhfX19Sl1dXWpf/zjHx9buXLlS2Le2+sF+M1vAH70IyU97Mvo0ViH3L+GD4D8tQZPPgkwbRpa\nAz/8IcDNNyvvJ0HozdixeF94POKOJwEIDRStA6isrLwxJSWlPjk5uQEAoKCgYFNpaemijIyMWnbM\n5s2b85cvX/4GAMCsWbMq2tvbR7S0tCQmJib2K9D8ve/5/r5wAf3yCxYo6WFf2P6ihw71jynIvaBn\nzsQfgrAyMTFY5PDQIcxgCQatBA4NFAlAU1NTksPhcLHHdru9saKiYlawYxobG+2BBODEiTXX/k5N\nzYW//CVX9cJqU6fiRX7XXX2fV2O1MUFYmdtuA9i5U5wAkAVgDOXl5VBeXq7a+ykSAJvNFsCb3h+v\n19sntMr3ur//fY2S7ohiyhSA3bv7P08XNBHu3H47plyvWhX8WLIAjCE3Nxdyc3OvPV67dq2i91M0\nv05KSmpyuVzXFjy7XC6H3W5vFDqmsbHRnpSUJKP4rDqwQLA/ZAEQ4Q6zAAIlSXC5dAl3xEtL06df\nhHYoEoCcnJy9dXV1qQ0NDclut3tASUnJA/n5+Zu5x+Tn529+8803HwIA2L17900jRoxoD+T+0Yup\nUwHq6voHgvfswa0eCSJcSU7GJIljx4SP++wzgFmzAIYN06VbhIYocgFFRUX1bNiwYdX8+fM/8Xg8\nkStWrHglIyOjtri4uBAAoLCwsHjhwoUff/zxxwtTUlLqY2NjO1977bWH1em6PIYOxf1Ff/tbgKIi\nfO7YMYBt2wD++Ecje0YQxmKzoRtoxw7honAffQTwzW/q1y9CO2zeYPaeTthsNq9efTlxAmDGDBz4\nR4wAWLkS0zh//nNdmicI0/LSSwCVlQCvvRb4/729WHgvmEgQ+mCz2frFWKVg6pXAWjF+PJaX+P3v\nAZqbccHZ448b3SuCMB4WB+CjuhpdPzT4hwaW2A9AC55+GiA3F6ClBeDBB62zeQdBaElmJpY0OXUq\ncP2tv/+d3D+hRFhaAAAAGRkAt9wC8Ic/qLvamCCsTEQEwK238lsBH32E1jMRGoStAACgz7+oyLej\nF0EQ/G6g5mbMoLv1Vv37RGhDWAaBCYLgZ/9+gLw8gPJygEmTfM+/+ipAWRnAX/5iWNcIPygITBCE\nqkyfDrBuHcbIqqpwYdhbbwH8538CfPe7RveOUJOwDQITBMHPww9jivSCBViS/dIlgL/9jQofhhrk\nAiIIgpcvvsDiiYWFAFE0XTQdSl1AJAAEQRAWhWIABEEQhCxIAAiCIMIUEgCCIIgwhQSAIAgiTCEB\nIAiCCFNIAAiCIMIUEgCCIIgwhQSAIAgiTCEBIAiCCFNIAAiCIMIUEgCCIIgwhQSAIAgiTCEBIAiC\nCFNkC8D58+fj582btzUtLe3oXXfd9Wl7e/uIQMc98sgjryYmJrZMnTr1K/ndtDbl5eVGd0FT6Pys\nDZ1f+CJbANavX//0vHnzth49ejRt7ty5n61fv/7pQMc9/PDDr5WVlS2Q30XrE+oXIJ2ftaHzC19k\nC8DmzZvzly9f/gYAwPLly9/48MMPFwc67rbbbtsZFxfXJrcdgiAIQhtkC0BLS0tiYmJiCwBAYmJi\nS0tLS6J63SIIgiC0RnBHsHnz5m1tbm4e4//8c88991/Lly9/o62tLY49Fx8ff/78+fPxgd6noaEh\n+Z577vnbV199NZW3IzYbbQdGEAQhESU7ggnu8rl169Z5fP9LTExsaW5uHjNmzJjm06dPXzd69Ogz\ncjsBoOwkCIIgCOnIdgHl5+dvfuONN5YDALzxxhvLFy9e/KF63SIIgiC0RrYAPP300+u3bt06Ly0t\n7ej27dvnPP300+sBAE6dOjX2m9/85t/ZcUuXLt1488037zp69Giaw+Fwvfbaaw+r0XGCIAhCIV6v\n19CfLVu2LJg0adLhlJSUuvXr1682uj9Kf06ePOnIzc39PDMz8+vJkycfeuGFFx73er3Q2toaf+ed\nd25NTU09Om/evE/b2tpGGN1XuT89PT2R06dPr7777rv/Fmrn1tbWNmLJkiXvpqen12ZkZNTs3r17\nViid37p1657JzMz8esqUKV8tXbr0z11dXQOtfH4PP/zwq6NHj26ZMmXKV+w5ofNZt27dMykpKXWT\nJk06/Mknn9xldP/lnN+TTz75q/T09Npp06YduPfee99vb28fLvf8DD25np6eyIkTJ9Y7nc5kt9sd\nnUCCh0kAAATuSURBVJWVtb+mpibD6A9dyc/p06fHVFdXT/d6vXDx4sUhaWlpR2pqajJ+/OMf/7Ko\nqOgpr9cL69evX7169er1RvdV7s/zzz//Hw8++OA799xzz2av1wuhdG4PPfTQG6+88sojXq8Xuru7\no9rb24eHyvk5nc7kCRMmHO/q6hro9Xrh29/+dsnrr7++3Mrnt2PHjtuqqqqyuQMk3/l8/fXXmVlZ\nWfvdbne00+lMnjhxYr3H44kw+hyknt+nn346j/V79erV65Wcn6Ent2vXrtnz588vY49/8YtfPP2L\nX/ziaaM/dDV/Fi1a9OHWrVvvnDRp0uHm5uZErxdFYtKkSYeN7pucH5fLZZ87d+627du338EsgFA5\nt/b29uETJkw47v98qJxfa2trfFpa2pHz58/HdXd3R919991/+/TTT+dZ/fycTmcyd4DkO59169Y9\nw/UyzJ8/v+zLL7+8yej+Sz0/7s/7779/73e+85235Z6fobWAmpqakhwOh4s9ttvtjU1NTUlG9klN\nGhoakqurq7NnzZpVESrrJn74wx/+z69+9asfR0RE9LLnQuXcnE7nhFGjRp19+OGHX5sxY0bVo48+\n+qfOzs7YUDm/+Pj48z/60Y+eHzdu3MmxY8eeGjFiRPu8efO2hsr5MfjO59SpU2PtdnsjOy4UxptX\nX331kYULF34MIO/8DBWAUM79v3Tp0pAlS5a898ILL/xg6NChF7n/s9lsXiue+0cffXT36NGjz2Rn\nZ1d7edJ2rXpuAAA9PT1RVVVVM773ve/9oaqqakZsbGynf4kTK5/fsWPHJv72t799oqGhIfnUqVNj\nL126NOTtt9/+LvcYK59fIIKdj5XP9bnnnvuvAQMGuB988ME/8x0T7PwMFYCkpKQml8vlYI9dLpeD\nq2BWpbu7O3rJkiXvLVu27C2WHsvWTQAAqLFuwgh27dp18+bNm/MnTJjgXLp06cbt27fPWbZs2Vuh\ncG4AOGOy2+2NM2fO3AMAcP/9979bVVU1Y8yYMc2hcH579+7Nufnmm3clJCS0RkVF9dx3333vf/nl\nl7ND5fwYfNej/3jT2NhoT0pKajKqn0p4/fXX/+Xjjz9e+M4773yHPSfn/AwVgJycnL11dXWpDQ0N\nyW63e0BJSckD+fn5m43sk1K8Xq9txYoVr2RmZtY88cQTv2XPh8K6iXXr1v2ny+VyOJ3OCZs2bSqY\nM2fO9rfeemtZKJwbAMCYMWOaHQ6H6+jRo2kAANu2bbtz8uTJX99zzz1/C4XzS09PP7x79+6brly5\nMsjr9dq2bdt2Z2ZmZk2onB+D73rMz8/fvGnTpgK32z3A6XROqKurS73xxhsrje2tdMrKyhb86le/\n+nFpaemimJiYLva8rPMzOsDx8ccf56WlpR2ZOHFi/bp1654xuj9Kf3bu3HmrzWbrzcrK2j99+vTq\n6dOnV2/ZsmVBa2tr/Ny5c7dZMdUu0E95efk3WBZQKJ3b/v37s3JycvZwU+xC6fyKioqeYmmgDz30\n0BtutzvayudXUFCw8brrrjsVHR3tttvtrldfffVhofN57rnn/nPixIn1kyZNOlxWVjbf6P5LPb9X\nXnnlkZSUlLpx48adYOPLypUr/yD3/ARrAREEQRChC+0IRhAEEaaQABAEQYQpJAAEQRBhCgkAQRBE\nmEICQBAEEaaQABAEQYQp/x+C+dfnjxjwpgAAAABJRU5ErkJggg==\n",
       "text": [
        "<matplotlib.figure.Figure at 0x22af4110>"
       ]
      }
     ],
     "prompt_number": 18
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": []
    }
   ],
   "metadata": {}
  }
 ]
}