CohensCongressExample.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#Cohen's Congress example"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##What this is about"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In 1994 Jacob Cohen argued that the reliance on p-values in psychology is a bad thing ([Cohen, J. The earth is round (p<.05). American Psychologist, 1994, 49, 997-1003](http://dx.doi.org/10.1037/0003-066X.49.12.997)). In his argument he used the following to exemplify that the logic of null hypothesis significance testing is wrong:\n",
    "\n",
    "> If a person is an American, then he is probably not a member of Congress.\n",
    "> This person is a member of Congress.\n",
    "\n",
    "> Therefore, he is probably not an American.\n",
    "\n",
    "Cohen's point was that the conclusion in this example, although it follows the [modus tollens](https://en.wikipedia.org/wiki/Modus_tollens) in logic, is not sensible. At the same time, the structure of the example is exactly that of null hypothesis significance testing:\n",
    "\n",
    "> If H0 is true, then this result (statistical significance) would probably not occur.\n",
    "> This result has occurred.\n",
    "\n",
    "> Then H0 is probably not true and therefore formally invalid.\n",
    "\n",
    "In the discussion associated with the article Baril and Cannon criticised this example as \"a particularly artificial and irrelevant illustration\".\n",
    "\n",
    "Here, we'll explore why Baril and Cannon came to this conclusion by framing Cohen's example properly with probabilities. Actually, [Edwin Thompson Jaynes](https://en.wikipedia.org/wiki/Edwin_Thompson_Jaynes) (Jaynes, E. T. Probability Theory: The Logic of Science. Cambridge University Press, 2003, ch. 1, [@Google books](http://books.google.de/books?id=tTN4HuUNXjgC)) provided a nice introduction to the extension of logic with probabilities and already pointed out that the modus tollens does not apply, when the implication (if A then B) is uncertain."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##Formalising the Congress example"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can formalise the example using binary random variables. Especially we define\n",
    "\n",
    "$A \\in \\{0, 1\\} :=$ whether a person is a citizen of the USA (American)\n",
    "\n",
    "$B \\in \\{0, 1\\} :=$ whether a person is a member of the Congress of the USA\n",
    "\n",
    "These random variables can take truth values 0 and 1 such that when, e.g., $A = 0$ a person is not an American and when $A = 1$ a person is an American. We can now make intuitive guesses about the probabilities which implement Cohen's example."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###Likelihood"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To state that \"If a person is an American, then he is probably not a member of Congress.\" we can set the probability of $K=0$ given $A=1$ to a very high value (this should actually just be the number of members of Congress divided by the number of Americans, but we choose 0.99 out of convenience):\n",
    "\n",
    "$$P(K=0|A=1) = 0.99.$$\n",
    "\n",
    "Additionally, we need the probability of not being in Congress under the premise that a person is not American which should be 1, but we allow for some, very unlikely mistakes and set it to:\n",
    "\n",
    "$$P(K=0|A=0) = 0.9999999$$\n",
    "\n",
    "These two probabilities fully specify the likelihood $P(K|A)$, because $P(K=1|A) = 1 - P(K=0|A)$ irrespective of the value of $A$. Now store these probabilities in Python variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "pnotKA = 0.99\n",
    "pnotKnotA = 0.9999999"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###Prior"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To specify the full probabilistic model we need to set the prior probability of a person being American. Of course, this depends on the group of people you consider. We take the whole world and use estimates of 300 million Americans vs. 6 billion people on earth:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.05\n"
     ]
    }
   ],
   "source": [
    "pA = 300.0 / 6000.0\n",
    "print pA"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###Posterior"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In his example, Cohen is interested in the plausibility of the statement that a person is not an American given that he is in Congress. We can now compute this posterior probability of $A$ given $K$ by applying Bayes' formula:\n",
    "\n",
    "$$P(A|K) = \\frac{P(K|A)P(A)}{P(K)}.$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###Marginal likelihood"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We see that the posterior depends on the probability $P(K)$ which is the so-called *marginal likelihood* and can be computed from the generative model $P(K|A)P(A)$ by *marginalising* variable $A$. *Marginalising* means that you remove all dependencies on $A$ from the model by summing all probabilities of $K$ under all possible hypotheses about $A$ weighted by the prior probabilities that these hypotheses are true. Here, there are two possible hypotheses: $A=0$ and $A=1$ leading to the following equation\n",
    "\n",
    "$$P(K) = P(K|A=1)P(A=1) + P(K|A=0)p(A=0).$$\n",
    "\n",
    "As we have above defined the probabilities $P(K=0|A)$ in Python, we now can easily compute $P(K=0)$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.999499905\n"
     ]
    }
   ],
   "source": [
    "pnotK = pnotKnotA * (1-pA) + pnotKA * pA\n",
    "print pnotK"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###Final probability"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, we can calculate the posterior probability of not being American under the information that a person is a member of Congress which represents the plausibility of the last sentence in Cohen's Congress example (notice that I plug the appropriate truth values for variables $A$ and $B$, i.e., 0 or 1, into the Bayes' formula above):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.000189963906758\n"
     ]
    }
   ],
   "source": [
    "pKnotA = 1 - pnotKnotA\n",
    "pnotAK = pKnotA*(1-pA) / (1-pnotK)\n",
    "print pnotAK"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "##Conclusion"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We see that the plausibility is very low, as it should be when we simply stick with our intuition. It is clear that we cannot apply modus tollens in this situation, but probability theory, which extends binary logic, leads to sensible inferences, just as Jaynes had already pointed out. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "##Null hypothesis significance testing"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Cohen's main point with this example was that people have incorrect intuitions about p-values and null hypothesis significance testing. Especially, Cohen wanted to point out that it is not a sensible conclusion to reject $H_0$ based on the low probability of the data under the assumption that $H_0$ is true. The computed probability `pnotAK` above is exactly the probability with which that conclusion would be true. In Cohen's example this is extremely low, but does this example actually reflect the typical situation in null hypothesis significance testing? Baril and Cannon argued against this.\n",
    "\n",
    "So let's replace the particulars of the Congress example with typical values from null hypothesis significance testing. First, let's get clear about what $A$ and $K$ mean here\n",
    "\n",
    "$A = 1$ := $H_0$ is true\n",
    "\n",
    "$K = 1$ := test is significant\n",
    "\n",
    "Test statistics are usually defined by two important quantities: the significance level of a test, $\\alpha$, and the power of a test, $1-\\beta$ where $\\alpha$ is the false positive rate, i.e., the probability of rejecting $H_0$ although it is true, and $\\beta$ is the false negative rate, i.e., the probability of accepting $H_0$ although it is false. $H_0$ is rejected when the test is significant, i.e., $K=1$. $\\alpha$ is often set to 0.01 such that"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.99\n"
     ]
    }
   ],
   "source": [
    "pKA = 0.01\n",
    "pnotKA = 1 - pKA\n",
    "print pnotKA"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Whereas $\\alpha$ is usually fixed to a given value, it is hard to get a good estimate of the underlying power of a test, because it depends on size of the effect under study which is unknown. [Button et al. (2013)](http://dx.doi.org/10.1038/nrn3475) have done a meta-analysis of studies in neuroscience and have found median estimated powers of 8-31%. Let's be positive and take the upper limit of these values such that $\\beta=0.69$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "pnotKnotA = 0.69"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, we need an estimate of the prior probability with which $H_0$ is true. Obviously this is not known and it depends on the skill of the researcher to generate sensible hypotheses. I would guess that across all possible scientific hypotheses the by far largest number of considered hypotheses are false. Because $H_0$ is often taken to be the opposite of the actual hypothesis under consideration, my guess corresponds to $p(A=1) \\gg 0.9$. But let's remain positive and assume a researcher with great skill whose hypotheses have a 50/50 chance of being true:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "pA = 0.5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now I can compute the posterior probability of correctly rejecting the null hypothesis by using the same equations as before:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.96875\n"
     ]
    }
   ],
   "source": [
    "pnotK = pnotKnotA * (1-pA) + pnotKA * pA\n",
    "pKnotA = 1 - pnotKnotA\n",
    "pnotAK = pKnotA*(1-pA) / (1-pnotK)\n",
    "print pnotAK"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Say what? Null hypothesis significance testing works!? The posterior of correctly rejecting $H_0$ is a little bit smaller than the original $\\alpha$, but it's still reasonably close! That was exactly the point Baril and Cannon made in their reply to Cohen. The examples are only extreme cases. In typical situations, null hypothesis significance testing is not so bad. Well, that's comforting. Also, there must be a reason why it stuck around for so long and is loved by so many. But, wait, what if the prior probability of hypotheses being true is not 0.5? That was the most uncertain quantity in our calculations. Here is some code which plots the probability of correctly rejecting $H_0$ based on different values of $P(A=1)$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Populating the interactive namespace from numpy and matplotlib\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEVCAYAAADzUNLBAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XmcXXV9//HXZ7Y7+8ZMJmSYkBDWGAIJWxDQQRECRWyR\nimhFsVWqxfqzLaL+tOanVktt1YfVCqJY2qpYFWWRpYAdkZ1AgAAJ2cnKZJvMvs/n98c5kwxjkjl3\n5txlZt7Px+M+zj3nfuecT85j5n5yvqu5OyIiImPJyXQAIiIyOShhiIhIJEoYIiISiRKGiIhEooQh\nIiKRKGGIiEgkGUkYZnarmTWb2crDlPm2ma01sxfMbFE64xMRkT+UqSeMHwFLD/WhmV0CHOvuxwEf\nBb6XrsBEROTgMpIw3P33QMthilwG3BaWfQqoNLO6dMQmIiIHl61tGPXAlhH7W4GjMhSLiIiQvQkD\nwEbtaw4TEZEMyst0AIewDWgYsX9UeOwNzExJRERkHNx99H/Kx5StCeMu4DrgdjNbAuxz9+aDFcy2\nyROHhpy+gSH6BofoHRgK3o949fYPBp/1jzg2MBh+FvzMyDIjP+vpD973hMd7+oOyvQNDrH/o36l5\n8/vpHRiif3CIRF4OifxcCvNzSOQF28KCXArzh185+98XFQT7Rfm5FBYE+0X7j4f7I9+POJafa5gl\n/XuXUsuWLWPZsmWZDiMr6F4coHtxwHj/ZjOSMMzsp8BbgRoz2wJ8EcgHcPeb3f1eM7vEzNYBncA1\nmYhzPHJyLPhiJjet1102+AjLll0EBEmrdyBIMMHrwPve/jce7+4L3/cN0t0/yJ6OPrr7Bvcf7+4b\npLd/kK7w856+8H34cqA4TB7BNm//fkki2A+2uZQk8ihO5FJSEG4TeRSHnxcngu1wmUReeu+fiIwt\nIwnD3a+KUOa6dMQyFeXk2P4ngFTrCxNTV+8gXX0DdIcJpat3INiGiaWzd4Cu3kF2tfWyua+LjnC/\nq3eAjnDbGZbr7BnADEoSeZQk8igtzKM0kUdJuC0tzKU0kR8cD19lI7Z72nvZsa+bssJ8ShK5WfcE\nJDJZZWuVlCSpsbExI9ctyMuhIC+H8qL82M7pHlTrdfQO0NkzGGx7B+joCbfh+46eAZpbe1jXPEB7\nT//+Y6+11/En33qM9u4BegeGKEnkUl6UT1lhPmVFeZQX5VNRlE95+H7/fnGwrSjOpzJ8X5iGpJtK\nmfq9yEa6FxNn2dYGkAwz88kcv6TewOAQHT0DtPcM0NbdT1t3/4j3A7R29dMaHm/tOrDd1xVsLQcq\ni/OpLC7Yn0gqiwuoKjmwrSopoLqkgKrSYFtRlE9Ojp5qJHuZ2bgavZUwRA7B3enpH6K1q5+Wrr4w\nkfTR0tnPvs4+Wrr6aensY29H3/7jezv76OgZoKIon+rSAqpLCzhi/yvBEaUF1JQlgvdlBdSWJSgr\nzFO1maSVEoZIlhgYHKKlq5+9HX3s7ehlT0df+Oplb0cfu9t72d0e7O9u76NvcIjasgS15QlqyxLM\nCLe15YXUVSSoC7fVpQly9eQiMVDCEJmkuvuCzgC7O3rZ1dbLzrYedrb1hq8emluDbVt3P0eUJphZ\nWcjMivBVWcisyqJgW1VEXXmCvNxsHo8r2UAJQ2SK6xsYChNID6+39vL6vm52tPawo6WH11u72d7S\nw56OXmrLEsyqKqK+qoj66iLqq4o5qrqIhiOKqa8qoiBPCWW6U8IQEfoHh2hu7WHr3m62tYSvvV1s\n3dvNlj1dNLf2UlNWwOyaYmYfUczsI0o4uqaYOTUlHF1bTFlhfL3dJHspYYjImAYGh9je0sPmPV1s\n3tPFa7s72byni427gm1JIo+5tSXMrS3hmBklzJtRyry6Uhqqi1TVNYUoYYjIhLg7za29bNzVwYad\nnWzY2cn6nR2s39nBrvZe5tSUcGxdKcfPLOO4maWccGQZR9eUqCF+ElLCEJGU6e4bZMPODta83s7a\n1zt4dUc7r+5oZ29nH8fVlXLirHLm15dz0qxyTqovU9VWllPCEJG0a+/pZ82OdlZtb2fVtjZe2d7G\nmh3t1FUUsuCoChY0lLOwoZKTGyooSWhiiWyhhCEiWWFgcIj1Ozt5aUsrK7e28uLmfaze0U5DdRGn\nHl3F4jmVLJ5TxbwZpRoRnyFKGCKStfoGhnh1RzsrXmvhuU0tPLexhbaeAU6bU8UZx1Rz5rxqTm6o\nIF8N62mhhCEik0pzaw/LN7bwzIa9PL1+L1v2dLF4bhVnH3sE5xxfw5vqy/UEkiJKGCIyqbV09vH0\n+r08vnY3j63Zzd7OPs45voa3nFBL40m11JYXZjrEKUMJQ0SmlO0t3Ty2Zje/W72LR9fsZnZ1MY3z\na3nHgpksOKpcEzZOgBKGiExZ/YNDPLuxhd++spMHX2qmr3+Qd5w8k6ULZ3LGMdUaC5IkJQwRmRbc\nnfU7O3jgxWbue2EHO9t6ufiUmbxz0SxOm1ulJ48IlDBEZFrasLOD3zy/g7uf207vwBB/fFo9f3J6\nPXNqSzIdWtZSwhCRac3deXlrG3cs38pdz21nXl0pV509m4sXziSRP7mX2o2bEoaISKhvYIiHX27m\np09s5qWtrVx+xlFcfc4cZtcUZzq0rKCEISJyEFv2dPFfj73Gz5/ewulzq7nmLXNYcuwR07qtQwlD\nROQwunoH+NXybdz6yEZKC/P4+NuP5R0L6qbl4EAlDBGRCIaGnP956XW+99B6uvoG+cSFx3LpqbOm\nVeJQwhARSYK78+ia3XzzvjV09Q3yNxcfzzsW1E2LqiolDBGRcXB3Hn55J/9876sUF+TyhT+ez6I5\nVZkOK6WUMEREJmBoyPnV8m18/d7VnH1sDTdceiIzK6fm/FVKGCIiMejsHeDfHlrHTx7fzMcumMeH\n3zJ3yq1nroQhIhKj13Z38n//eyWt3f187T0LWdBQkemQYqOEISISM3fnl89s5R/vXs2VSxr45EXH\nU5A3+Z82lDBERFJkV3svn7n9RXa19fLND5zKvBmlmQ5pQpQwRERSyN35r8de45v3r+HTl57IlWc1\nTNouuEoYIiJpsPb1dq677TlOPbqSL717waSc2FAJQ0QkTTp6Brj+py+wY18P37tmMUdWFmU6pKQo\nYYiIpJG7c9PD6/n332/ipmtOm1SD/ZQwREQy4LcvN/Pp21/kn967kLe9qS7T4UQy3oSRkf5hZrbU\nzFab2Vozu+Egn9eY2f1m9ryZvWRmH8pAmCIiY3rbm+q45c9P54afreQXT2/JdDgplfYnDDPLBV4F\nLgC2Ac8AV7n7qhFllgEJd/+smdWE5evcfWDUufSEISJZYX1zBx+8+WmuPu9oPnr+vEyHc1iT6Qnj\nTGCdu29y937gduBdo8rsAMrD9+XAntHJQkQkm8yrK+UXf/1mfvr4Zn7QtCHT4aREJhJGPTDyuW1r\neGykW4A3mdl24AXgk2mKTURk3GZWFvLjjy/hP36/idt+vynT4cQuLwPXjFKH9DngeXdvNLN5wINm\ndoq7t48uuGzZsv3vGxsbaWxsjCtOEZGkzaoq4scfX8J7v/sk+bk5vO/NszMdEk1NTTQ1NU34PJlo\nw1gCLHP3peH+Z4Ehd79xRJl7gX9w98fC/YeBG9x9+ahzqQ1DRLLSpl2dvPe7T/Cldy/gwpNnZjqc\nN5hMbRjLgePMbI6ZFQBXAneNKrOaoFEcM6sDTgCmZqWgiExJc2pL+P6HT+dz/72Sl7e2ZjqcWKQ9\nYYSN19cBDwCvAD9z91Vmdq2ZXRsW+ypwupm9ADwEfNrd96Y7VhGRiVg4u5IvX7GAj/xwOa/v68l0\nOBOmgXsiIin2bw+t474XdvCz686mOJGJpuM30khvEZEs5e787U9eAOAb7z81w9GkuQ3DzArNLDGe\nnxURmW7MjK9csYAXN+/jV8u3ZjqccYuUMMwsx8wuN7Ofm9k2YCPwmpltM7NfmNmf2GSdGF5EJA2K\nE3l8++rFfOXXq9i0qzPT4YxL1CeMJuA04J+BY9z9SHefCRwTHjsD+F1KIhQRmSLm15fziQuP5ZP/\nuYK+gaFMh5O0SG0YZpZw996Jlomb2jBEZLJxdz7yw+XMm1HKZy87KSMxpLQNYzgRmNmNoz8bPpbu\nZCEiMhmZGf/03lO4Y/k2Vm6ZXOMzkm30vvAgxy6JIxARkemiurSAz1x6Ip//+UoGhyZPLUnURu+P\nmdlK4AQzWznitQl4MaURiohMQZefUU9RQS4/fvy1TIcSWdQ2jAqgCvga8BnACCYRbM/kCGy1YYjI\nZLauuZ0r//VJ7rv+PGZUFKbtuqluw2h1903Ar4CW8P3VwA/MbHGyFxURETi2roz3nt3AV+5cNXbh\nLJBsG8YX3L3NzM4F3g7cCtwUf1giItPDde84judfa+HJdXsyHcqYkk0Yg+H2UuAWd78HyI83JBGR\n6aOoIJdPXXw8/3Lvq2R7FXuyCWObmX2fYEry35hZ4TjOISIiI1y2uJ7Wrn4eWb0r06EcVrJf9u8B\n7gcudPd9BA3h18celYjINJKbY+FTxpqsfspINmF0AyXAVeF+PrAv1ohERKahi06eyZA7D6x8PdOh\nHFKyCePfgCXA+8L9jvCYiIhMQE6O8beXnMA371uTtYP5kk0YZ7n7XxE8aRCOwVCjt4hIDBpPqqW0\nMI+7V2zPdCgHlWzC6DOz3OEdM6sFJt+UiyIiWcjM+OsLj+P7v92QlW0ZySaMfyUYvDfDzL4KPEYw\n+ltERGJw3gm1dPcP8syGlkyH8geSXqLVzE4iGLQH8LC7Z2yIoqYGEZGp6EePbGT5hha++6HUTKSR\nliVazexGd1/l7t8JX6sONuW5iIiM3xVnHsVja3ezY193pkN5A01vLiKSZcoK87ls8Sx+8vjmTIfy\nBpreXEQkC33w3Dnc/uQWegcGxy6cJnkRy/0EuA/4R+AGgunNIZjePPtnzBIRmWTm1ZVy4pFl/GbF\nDi4/46hMhwOMo9E7m6jRW0SmsodeauZfH1zLnZ86N9bzpqXRW0RE0uf8+TPY2drLuub2TIcCKGGI\niGSt3BzjklOO5J4VOzIdCjCObrVRjomISDwuXXwkd6/YnhUjv9WtVkQki506u5K+gSFWbW/LdCjq\nVisiks3MjHcumsXdz2W+WipSLykzqyBYLCmrutWql5SITAevbGvj2luX88jnz8cs6c5Nf2C8vaQi\njcNw91agFXhvshcQEZGJOWlWGQV5OTy/eR+Ljq7KWBxRB+4BEK7h/W5gzoifdXf/UsxxiYhIaLha\n6p7ndmQ0YSTb6H0ncBnQT7DaXgfQGXdQIiLyRpcuOpLfvLA9o6vxJfWEAdS7+0UpiURERA7p2Loy\nqksSPLNhL0uOPSIjMST7hPG4mS2c6EXNbKmZrTaztWZ2wyHKNJrZCjN7ycyaJnpNEZHJ7oIFM/jd\nql0Zu36yCeM84FkzWzOia21S3WrDJV6/AywF5gNXhYsyjSxTCXwXeKe7LwCuSDJOEZEp59zja3l0\nTeYSRrJVUhcf5FiyFWpnAuvcfROAmd0OvAsYuXLf+4BfuvtWAHffneQ1RESmnEVzKnltdxd7O/qo\nLi1I+/WTesIIv+RbgRnA7PB1dJLXrAe2jNjfGh4b6Tig2sz+18yWm9kHkryGiMiUk5+bwxnHVPP4\n2sz8HzrZbrUfAf4aaABWAEuAJ4C3JXGaKE8k+cBigrXDi4EnzOxJd1+bTLwiIlPNuSfU8PtXd3Pp\nollpv3ayVVKfBM4AnnD3883sROBrSZ5jG0HCGdZA8JQx0hZgt7t3A91m9ghwCvAHCWPZsmX73zc2\nNtLY2JhkOCIik8e5x9fwg6aNuHvkUd9NTU00NTVN+NpJLaBkZsvd/XQzex5Y4u49ZvaKu89P4hx5\nwKsETw/bgaeBq9x91YgyJxI0jF8EJICngCvd/ZVR59LUICIyrbg7Z/+/h/nJx5dwzIzScZ0jpVOD\njLDVzKqAXwMPmlkLsCmZE7j7gJldBzwA5AI/dPdVZnZt+PnN7r7azO4nmNhwCLhldLIQEZmOzCzs\nLbV73Alj3Nce7//QzawRKAfud/e+OINKIgY9YYjItPPrZ7dx3ws7uPnDp4/r59OyRKuZ3RY+YeDu\nTcCjwE3JXlRERMbvnONreHLdHgYGh9J63WQH7p3i7i3DO+6+l6A3k4iIpEltWYJZVUW8uKU1rddN\nNmGYmVWP2KkmaIcQEZE0Ovf4Gh59Nb3jMZJNGP9CMCbiy2b2FYIxGF+PPywRETmcc46vSfsAvqR6\nSbn7f5jZswQD9Ry43N1fTklkIiJySKfMruTlrW0MDTk5ORNfhS+KZLvVEiYIJQkRkQyqKimgqiSf\njbs6mVeXnu61kaqkzOzRcNt+kFdbakMUEZGDWTi7khc270vb9aK2YWwKt19w97JRr/IUxSYiIodx\nckMFK9PYUypqwlhsZrOAD5tZ9ehXKgMUEZGDW5jmhBG1DeMm4GHgGODZg3w+N7aIREQkkgUNFaza\n3sbA4BB5ucl2ek1epCu4+7fd/STgR+4+d/QrxTGKiMhBlBXmM7OykLXNHWm5XrILKP1lqgIREZHk\npbMdI2ovqcfCbYd6SYmIZI9TGipZuSU9PaUitWG4+znhNr1z6YqIyGGd3FDBr5/dlpZrpb6VRERE\nUmZ+fTlrXm+nd2Aw5ddKdk3vQuDdwJwRP+vu/qWY4xIRkQiKE3kcXVPCmh0dnNxQkdJrJfuEcSdw\nGdAPdISvzriDEhGR6BY2VPBiGkZ8JzuXVL27X5SSSEREZFxObqjgxS2tvD/F10n2CeNxM1uYkkhE\nRGRcFs6uTEvX2mSfMM4DrjGzjUBveMzdXUlERCRDTpxVxsZdHXT3DVJUkLo17ZJNGEsBI1gLQ0RE\nskAiL5eja0rYuKuD+fWpa/hONmFcQZAshlfrcKDVzJ519+djjUxERCKbW1vChp2dKU0YybZhnAb8\nJTALqAeuBS4GbjGzG2KOTUREIppbW8LGXanttJpswmgAFrv737r73xAkkBnAW4EPxRybiIhElI0J\noxboG7HfD9S5exfQE1tUIiKSlHQkjGTbMH4MPGVmvyZox3gn8BMzKwFeiTs4ERGJZu6MEjbs7MDd\nMbOxf2AckkoY7v5lM7sfOIegwftad18efpzqMSMiInII1SUFALR09lNdWpCSayT7hIG7PwM8k4JY\nRERknMyMubWlbNzVQXVpalbO1noYIiJTRKrbMbQehojIFDF3RmoTRlK9pMzsT82sLHz/BTO7w8wW\npyY0ERFJxjHh4L1USbZb7d+7e7uZnQu8HbgVuCn+sEREJFmprpJKNmEML+l0KXCLu98D5McbkoiI\njMec2hJe293J0FBqpvtLNmFsM7PvA1cCvwlX4NMyryIiWaAkkUdFcT47WlMzjjrZL/v3AA8AF7r7\nPqAKuD72qEREZFzm1JSwMUXtGEklDHfvdPdfuvvacH+Hu/9PSiITEZGkBe0YHSk5d0aqk8xsqZmt\nNrO1h5vl1szOMLMBM7s8nfGJiExWqexam/aEYWa5wHcIFmOaD1xlZicdotyNwP0cWH9DREQOIxjt\nnZqEEWngnpkdTcRV9tx98xhFzgTWufum8Ny3A+8CVo0q9wngF8AZUa4rIiLBWIyMJgzgNqIvy3r+\nGJ/XA1tG7G8FzhpZwMzqCZLI2wgShpaEFRGJoOGIYna09NA3MERBXryVSFGnBmmM8ZpRvvy/BXzG\n3d2CeXpVJSUiEkFBXg5HVhWyZU8X8+rinc0p6dlqY7CNYOW+YQ0ETxkjnQbcHs7pXgNcbGb97n7X\n6JMtW7Zs//vGxkYaGxtjDldEZHKZW1vChl2d+xNGU1MTTU1NEz6vuae3tsfM8oBXCaYW2Q48DVzl\n7qPbMIbL/wi4293vOMhnnu74RUSy3Zd+9TIzKwv56PnzDvq5meHuSdfcpP0Jw90HzOw6ggGAucAP\n3X2VmV0bfn5zumMSEZlKZlUVsb2lO/bzZqJKCne/D7hv1LGDJgp3vyYtQYmITBF15YWs2LQv9vNG\nThhmdh5Br6WZBJMQ7gKe0EhvEZHsMqMiwc62+OeTijoO43MEs9KuADoJqpLKgbeb2dvc/TOxRyYi\nIuNSV15IcwomIIz6hPHSwXooAb8wsyviDEhERCZmRnmC5rZe3J2wt2ksoiaMU8zsVOA5oIugSqoE\nWAjUEozIFhGRLFCcyCORl0NrVz+VJQWxnTfqwL0vm9kFwJuBGQRzUDUDjwK/jS0aERGJRV1FIc1t\nvelPGADu/hDwEICZVbr7PjOr0kAIEZHsU1eeoLm1hxOOLIvtnOOdaOSD4fbquAIREZH41JYXxt5T\nSsuriohMQXUVCXa29cZ6TiUMEZEpKBVda5UwRESmoBkVhXrCEBGRsQ03esdJCUNEZAoa7lYbp/Em\njAdHbUVEJIvMKE+wq62HoaH4Rj6MK2G4+ysjtyIikl0S+bmUJPJo6eqL7ZyqkhIRmaLqKgppbo2v\nWkrTm4uITFEzyoNpzufXl8dyPk1vLiIyRQVPGPH1lNL05iIiU9TwNOdx0fTmIiJTVF1FIWtfb4/t\nfJGfMIA24Bw0vbmIyKRQV17Io6/uju18URPGV4GfAc8AT7v7ruEPzKwCaI0tIhERicWMigTNMc5Y\nGzVh3ANsBq4AfmBm3cDTBAlkHvDx2CISEZFY1JUXsjMD3Wo/7+69wK1m9iHgTmARsBg4K7ZoREQk\nNrXlCfZ09DI45OTmTHxt76gD994z4r27e4u7/9bd/xn4/ISjEBGR2OXn5lBelM/ejnieMqI+Ydxo\nZucCjxE0eo/0YiyRiIhI7IZHe9eWF074XFETxjeA5QTVT2eY2fMEI71XAA3AVROOREREYheMxehh\nARUTPlekhBFWPQE0DR8zs5kECUQN3iIiWWpGeXzzSUWeS2o0d38duNPMdo1ZWEREMiJY2zuerrVj\nNnqbWcLMag71ubs/PqLs7FiiEhGRWMS5tveYCSPsTrvEzK4ys6KDlTGzKjP7KHB0LFGJiEgs4lzb\nO2qV1O+ATwHXm1kDwVxS+eG2C9gK3OLuGvEtIpJF6srjq5KKmjC+TjD9RwNQD1zi7p2xRCAiIilT\nU5Zgd3s8q+5FTRgr3f27AGZ2JHAlcGssEYiISMqUF+XT1t0fy7mijvTeXwHm7jsIZq4VEZEsV5LI\npXdgiP7BoQmfK+oTxmdGrIexAvDhD8yszt2bJxyJiIjEzswoL8qjvXuA6tKCCZ0rasK4jWBm2iXA\nu4FFZnY9wVQhtcDVE4pCRERSZrhaKi0Jw92/HL69f/iYmc0jGOn9kfFc2MyWAt8iWB/8B+5+46jP\n3w98GjCgHfiYu2veKhGRJJUXxtOOMZGR3uuB9Wa2NdmfNbNc4DvABcA24Bkzu8vdV40otgF4i7u3\nhsnl+wRPOCIikoTyojxaY0gYURu9D8ndHxnHj50JrHP3Te7eD9wOvGvUeZ8YMa7jKeCoiUUqIjI9\nlRXl09Y9MOHzTDhhjFM9sGXE/tbw2KH8OXBvSiMSEZmi4upaO+4qqQnysYsEzOx84MPAOQf7fNmy\nZfvfNzY20tjYOMHQRESmlt3rnuPH9z3Hqw9UT+g85h75uzs2ZrYEWObuS8P9zwJDB2n4XgjcASx1\n93UHOY9nIn4RkcnkOw+upbtvkOv/6EQg6Grr7kmv2ZqpKqnlwHFmNsfMCghGjt81skA48+0dwJ8d\nLFmIiEg05TG1YWSkSsrdB8zsOuABgm61P3T3VWZ2bfj5zcDfA1XA98wMoN/dz8xEvCIik1l5UT7t\nk7gNA3e/D7hv1LGbR7z/C+Av0h2XiMhUU16YF0ujd6aqpEREJE3Ki/KzYxyGiIhkt7jaMJQwRESm\nuLjGYShhiIhMcRXFShgiIhJBYX4Og0NO78DghM6jhCEiMsUFa2Lk0z7BdgwlDBGRaaAshnYMJQwR\nkWkgjinOlTBERKaBYBElVUmJiMgY4pgeRAlDRGQaKC+a+PQgShgiItNAMBZDVVIiIjKGOEZ7K2GI\niEwDShgiIhKJutWKiEgkccxYq4QhIjINqEpKREQiUbdaERGJpEJVUiIiEoWqpEREJJJEfi4Avf3j\nXxNDCUNEZJqY6FOGEoaIyDQx0bEYShgiItPERMdiKGGIiEwTqpISEZFIJjoWQwlDRGSamOhYDCUM\nEZFpQlVSIiISiRKGiIhEUl6UR2uXEoaIiIxB3WpFRCSSMlVJiYhIFOVFebT1KGGIiMgY1K1WREQi\nKS/Kp32yVUmZ2VIzW21ma83shkOU+Xb4+QtmtijdMYqITDVlk62XlJnlAt8BlgLzgavM7KRRZS4B\njnX344CPAt9Ld5yTTVNTU6ZDyBq6FwfoXhygewGJvFzycm3cP5+JJ4wzgXXuvsnd+4HbgXeNKnMZ\ncBuAuz8FVJpZXXrDnFz0x3CA7sUBuhcH6F4Eyovyx/2zmUgY9cCWEftbw2NjlTkqxXGJiEx5ky1h\neMRyo5+bov6ciIgcQnlR3rh/1tzT+z1sZkuAZe6+NNz/LDDk7jeOKHMT0OTut4f7q4G3unvzqHMp\niYiIjIO7J92YMf5UM37LgePMbA6wHbgSuGpUmbuA64DbwwSzb3SygPH9g0VEZHzSnjDcfcDMrgMe\nAHKBH7r7KjO7Nvz8Zne/18wuMbN1QCdwTbrjFBGRN0p7lZSIiExOk2Kktwb6HTDWvTCz94f34EUz\ne8zMFmYiznSI8nsRljvDzAbM7PJ0xpcuEf8+Gs1shZm9ZGZNaQ4xbSL8fdSY2f1m9nx4Lz6UgTDT\nwsxuNbNmM1t5mDLJfW+6e1a/CKqt1gFzgHzgeeCkUWUuAe4N358FPJnpuDN4L84GKsL3S6fzvRhR\n7rfAPcC7Mx13hn4nKoGXgaPC/ZpMx53Be7EM+NrwfQD2AHmZjj1F9+M8YBGw8hCfJ/29ORmeMDTQ\n74Ax74W7P+HureHuU0zd8StRfi8APgH8AtiVzuDSKMp9eB/wS3ffCuDuu9McY7pEuRc7gPLwfTmw\nx93HPxtfFnP33wMthymS9PfmZEgYGuh3QJR7MdKfA/emNKLMGfNemFk9wRfG8NQyU7HBLsrvxHFA\ntZn9r5mJzfd+AAAENklEQVQtN7MPpC269IpyL24B3mRm24EXgE+mKbZslPT3Zia61SZLA/0OiPxv\nMrPzgQ8D56QunIyKci++BXzG3d3MjD/8HZkKotyHfGAx8HagGHjCzJ5097UpjSz9otyLzwHPu3uj\nmc0DHjSzU9y9PcWxZaukvjcnQ8LYBjSM2G8gyISHK3NUeGyqiXIvCBu6bwGWuvvhHkknsyj34jSC\nsTwQ1FdfbGb97n5XekJMiyj3YQuw2927gW4zewQ4BZhqCSPKvXgz8A8A7r7ezDYCJxCMD5tukv7e\nnAxVUvsH+plZAcFAv9F/8HcBV8P+keQHHeg3BYx5L8xsNnAH8Gfuvi4DMabLmPfC3Y9x97nuPpeg\nHeNjUyxZQLS/jzuBc80s18yKCRo4X0lznOkQ5V6sBi4ACOvrTwA2pDXK7JH092bWP2G4BvrtF+Ve\nAH8PVAHfC/9n3e/uZ2Yq5lSJeC+mvIh/H6vN7H7gRWAIuMXdp1zCiPg78VXgR2b2AsF/mD/t7nsz\nFnQKmdlPgbcCNWa2BfgiQfXkuL83NXBPREQimQxVUiIikgWUMEREJBIlDBERiUQJQ0REIlHCEBGR\nSJQwREQkEiUMERGJRAlDREQiUcIQicDMEhm4ZmG6rylyOEoYIiEz+4KZrTOzvzCzvzKzm8yszMwu\nBcpGlPs7M9tuZh8ws4SZfd3MVpnZH0/w+jlm9o0Rh44ys3dM5Jwiccr6uaRE0ugZoM7dfwBgZjcC\nnwLWjlp0aDlwt7v/Z1huFdDp7r8e74XNrIpgLp+3Dh9z93Vm9kdm9ri7d4733CJx0ROGyAFLgOdG\n7A9P/Tw6EZwFPDZi/y3A0xO5sLu3uPs3gLZRH90DvH8i5xaJixKGyAFnACsAzOxIghXJGsJ1JEaX\nqzOzD5rZB4GLmGDCOBR3Xw8sSMW5RZKlKimRAxYBx4QrsZUDlwD/cpByx7r7FbC/KumLY62TbWbz\ngUO1R9zm7vsO8+O5Y0YukgZKGCJAmCQ2u/svRx3PH7U/E9g14tAiwqcLMzuLoHF83uj1OML1J8a7\nBoV6S0lWUMIQCYxulxg2eJByK0bsLyZoLAf4U3f/OzObb2az3X3zcKExnjD+Y4yldIcOH7pIeihh\nyLRnZm8BPgJsMbO6UctUdo0o92bg48BeM6sHFgJXAU+ZWQVQHBbtAOqA/QljrCcMMysBPgqcZGb/\nh2BVvE4Llk3siOGfKTJhShgy7bn7I8D5h/h4q5lVhb2YHido4B62DbhveMfMWsO3lUBSa8qH3Wa/\nGb5GOgV4MplziaSKekmJHN4twJ9GLHuvmZ0PDI2sjpqgC4Cfx3QukQnRmt4iYzCz84DXYkwCUa97\nMpDr7s+n87oih6KEISIikahKSkREIlHCEBGRSJQwREQkEiUMERGJRAlDREQiUcIQEZFIlDBERCQS\nJQwREYnk/wO1bPD0RKuTqwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7fc51aca9490>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# import numpy and matplotlib and allow for inline plots\n",
    "%pylab --no-import-all inline\n",
    "\n",
    "# make a range of prior probabilities\n",
    "pA = np.linspace(0.0, 1.0, 100)\n",
    "\n",
    "# compute posterior with previous assumptions\n",
    "pnotK = pnotKnotA * (1-pA) + pnotKA * pA\n",
    "pKnotA = 1 - pnotKnotA\n",
    "pnotAK = pKnotA*(1-pA) / (1-pnotK)\n",
    "\n",
    "# plot results\n",
    "ax = plt.axes()\n",
    "ax.plot(pA, pnotAK, color=plt.cm.Blues(0.75));\n",
    "ax.set_xlabel('$P(H_0=1)$');\n",
    "ax.set_ylabel('$P(H0=0$ | significant test)');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The plot shows that the entire validity of null hypothesis significance testing hinges on the prior probability with which the null hypothesis is true, i.e., the researchers ability to identify sensible hypotheses. For less likely scientific hypotheses ($P(H_0=1) \\rightarrow 1$) it is also less likely that they are true, despite a significant test. For example, if only 1 in 10 of the questions researchers pose about their data turn out to be true on average, a significant test only suggests a probability that the corresponding hypothesis is true of:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "p(H0=1)=0.90 then p(H0=0|significant test)=0.78\n"
     ]
    }
   ],
   "source": [
    "print \"p(H0=1)=%4.2f then p(H0=0|significant test)=%4.2f\" % (pA[89], pnotAK[89])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**How many of your hypotheses have turned out to be true in the past?** Remember, I don't mean how often your tests turned out to be significant.\n",
    "\n",
    "There was another uncertain quantity in our computations: the typical power of tests. [Button et al. (2013)](http://dx.doi.org/10.1038/nrn3475) estimated a very low median power of studies in neuroscience of 8-31% and argued that this was a bad thing, because it leads to irreproducable results. Baril and Cannon, on the other hand, argued the opposite, saying that the typical low power of studies prevents us from finding trivial, i.e., tiny effects with null hypothesis significance testing. They gave an estimate of power for medium sized effects in psychology of 0.57 (see [Rossi, 1990](http://dx.doi.org/10.1037/0022-006X.58.5.646)). With ever increasing sample sizes, e.g., from online surveys, it is also possible to achieve very high powers of, let's say, 0.9. What happens to the probability of correctly rejecting a null hypothesis then? To check this we just replot the curve above for all the values mentioned here:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEVCAYAAADzUNLBAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3Xd4VFXewPHvmZo26Y2EEHqvIgKiEhUVEXVXsaBrf3ft\nuvraViy4rorv6uruui4qNlZdXcvaC4oCFqQoSO8kpPee6fN7/5hhCCRAJpkklPN5nvvcdu6dk4jz\ny+lKRNA0TdO0gzF0dwY0TdO0w4MOGJqmaVqb6IChaZqmtYkOGJqmaVqb6IChaZqmtYkOGJqmaVqb\ndEvAUEq9pJQqVUqtPUCavymltiqlflFKjenK/GmapmktdVcJ42Vg6v5uKqWmAf1FZADwO+CfXZUx\nTdM0rXXdEjBE5Fug+gBJzgFeDaRdBsQrpdK6Im+apmla6w7VNoxMIL/ZeQHQs5vyommapnHoBgwA\ntc+5nsNE0zStG5m6OwP7UQhkNTvvGbi2F6WUDiKapmntICL7/lF+UIdqwPgQuAl4Uyk1AagRkdLW\nErZl8kSPx0tFTSNllXWUVdVTVllPaVUdZZWB40r/cUlFLdX1TSTFRZOaGEt6SixpSbGkJ8eSlhxL\nenIcPVLi6BE4joq0tPgsn09weXy4vD6cHp//uNnm9Hj9e7dvTzr33vcczc733PPhcO951unx+t/v\n9r9j25cvk37SZbi9glcEq8mA1WwiwmzAYjZiNSmsZiMWkxGLSWExGTGbFBaTAbPRgMm4e68wGhUm\nowGjQWE0KAxKYTAojEqhDGBAgVIo5f/9+4TAJvhEEAGvz7/3EdjLnrQGBSaDwmQwYAp8hslA8Nio\ndn82mNTuPIBRGTAY8OdH4c+PCpwDKnA+769z+N3v70G1KKQC4i+qikhg77/sfw8o9rwzeE3tcy3w\nWbvv7b3fkz+lFAbDnnuGfdIaAvcNzc5V4H64zJ49m9mzZ4ftfYcz/bvYo73/xrolYCil/g1MBpKV\nUvnAg4AZQESeE5FPlVLTlFLbgEbgqo58nslkJD3Z/8V/MB6Pl/LqBkoq6iirrKOkoo6Silq27yrn\n+5+3B89LKuqwmk2kp8T6g0hKXDCgZKT695mp8fRIicNi7vxf82zXImbPngL4g9buAOPf/MdOd8tr\nDpf/2B489tLY5MW++3rgmt0V2Jofu7wIEGkxEmUxBvdRFpP/2Bo4t5oC141EW01EWgxEmk1EWAxY\nLSasJgMWkwGL2YDFZMRkVCilcHsFj2/35sPT7NzpFTxeLx6v4PYJbq8Pj09we4XiBic/FtYGz3ff\ncwWOd7/XZFCYjQqL0YDVpIgw+oOo1WjAYgwEU0Pg2GjAbPQHObNRYVIKk1E1C3KBwAYYAkFWQSBY\nSnDfPHD6RPD59lzfna550Nk3qOwJ4Luv+QPrnqAe+PzAsdPjo8np9T9n8D+vCG9Q0o4e3RIwRGRm\nG9Lc1BV52ZfJZAwGgAMREWrq7RSX11JSXktxeS1F5bXsyC/n+5+3Bc9LK+qIt0WSkRpPRmocGanx\nZKbGk5kW2FITyEyLxxYdEbafwWBQRAa+wDubKxCYmpxemlwe7C4vTa7WzxtdHsrrnTS5vDQ6PP69\n00Oj0xO83+jwX1MKoq2mwGbEFmEm2mokJsJ/LSbCv9kiTMRYTcRGmIiJNmOLMLM6IYIZw1ODz7T2\n5egT2SugNA8mLq+/1Lb72O31B2C3V2j0eHE5A6U+757S5L7nu9Objcpf2jMZsBgNRASOrabAsXnP\ntQiTvyQYYTZgNe55xmxUmI0GFOwVaLwiiE/wBkp0Lq8veN3n8wehmkY3W0sa8Yrg9fn/mNj9b8QY\nCCC7g92ejb3Og0GxWYA0Gf3BSju6qMN5PQyllBzq+ff5fJRVNVBUVkNRWQ2FgX1BSQ1F5TUUldVS\nUFKNyWSgZ1oCPdP8ASQrPYGe6f7z3ccRVvN+P2fRokXk5OR03Q/WiUT8X7wNTg+NTn9waXB6aHB4\nAgHGS73DTUOz6/UO/77B4SFv/QrMGcOpd3hwenxEW43ERpqxRZiwRZiJjTQRF2UOHsdGmomLMhMX\naSY2sI8P7CM6EHR9IsFg4/T4cASqJR2Bc2fwOFDiC5w73LvvebEHrtvdXgxKEWk2EGk27rWPMhv3\nHFv8x1FmA9EWI6t+/I7TTj2FCLMh+AXvE8Hr828+H8Fjb7Pr3sB1j0/wemXPcWDzeAV2Vy0GAonJ\n6A8mZqNhz3Xj3ufmQOmxOxxJ/490lFKqXW0YOmAcAnaXVgpKqiksrSG/pIqC0moKSnYf11BYWkO8\nLZKe6Qn06pFIVnoCWT0S6RXYsjMSSYiN0lUN+/B4fcGAUmd3U2f3UO/w7+vsbmrtbuqa3NQ5PNQ2\n7TmvaXJT2+RGKfzBI8pMfJSF+Cgz8dH+44TAPjHGf5wYbSEh2kJcpBmDIbz/HSQQfOyBIGIPVCPu\n2XtpcvlocgdKdIF7jYESnsvrI8ocqBa0GIkJ7KOtRmIsJmKs/msxVv+xzeqvKjxQfnyBkk2wqjBQ\n9efx7rnm9gaue6VFVaDJqDCb/FV+ZqPB36ZmNPivmfxBSP977hw6YBzhfD4fJRV1FJTWsKu4ivyS\nKnYVVZFfXMWu4mryiisRgeyMRLJ7JNIrI4nsjER6ZyYH9knE26K6+8c4rIgIDreP2iY31U0uagOB\npKbRRXVgX9XoorrRTU2Tf1/V4KLB6SEu0kxijIXEaAtJtsA+xkKSzUpys32yzYotwtTpX4xen9Dk\n9vpLbK49W73TQ6PLS4PTS4PLQ0PguN7hrxa0WU3+LcJIbKDqLzbCRFyzLcZqwtjGACkiweCxuwrQ\n7fXh9gSq/5rtvSJYjIZAJw1/O5c10Na1u1rPZNQBpT10wNCoqW8ir6iKvKJK8oqqyC2sZFdRJblF\nleQWVmIyGuidmUTvzGT69EymT2YSfXom07dnMlnpiZjNnd/mcTTweH1UN/mDR1WDP6hU1jupbHBR\n2eiiot5JZb2Ligb/3uX1kWKz+rfYwGazkhrr39LiIkiLjSAxxtLmL+aOEvFXpdU7/UGlzuGhzumh\n3uGl1uGh1uEvldU5/AEnxmoiPtJEXISZhCgTCZFm4iPNJEaaSIzyVweG2ubhC7TLuDwS7Cnocu+p\nznO6fSgVaPMJtAVFWAxEmI1EmA3dWv11qNMBQzsgEaGyppGdhRXkFlaSW1DBzsJKdhRUsCO/nNLK\nejLT4unbM5l+WSn0Cez79/IfH6j9ROsYu8tLeb2Tinon5XVOynbv6xyU1vr3ZXVO6uxukmKspMdF\nkBZnpUd8JOlxEfRIiCA9LoKMhEjSYq2YjF07HtfrE+ocHmocHmrsbqrtbmqaPFQHjqua3NjdPuID\nwSMpykJStJnkaAspgX1UO9qKJNBxwd/7zx9A7Lt7ALp8COJv27EY/T3zAj33LDqQ6IChdYzT5Sav\nqIod+RXsKChne74/kGzPLyevqIrUJBsDeqXSv1cK/XqlMCA7lYHZaWRnJOmSSRdxeXyU1zkpqXVQ\nWuuguNZBSY2D4ho7JbUOiqodVDY4SbZZyUiIJDMhkp6Je/Y9E6PITIjEcoC2iU7L++5SV5ObyiY3\nFY0uKhv9+/JGNyaDIjXGQkqMhdQYC2kxFtJsVlJiLJjaWapye33YXT4crkD7TqAruE8g2mrca4sw\nG46qIKIDhtZpPB4vu4qr2LarnG27ytiaV8b2/HK25pVRXF5Lrx6JDMxOZUDvNAZmpzGwTxqD+6SR\nFB/T3Vk/6ri9PkprHRRW2Sms9m8FVXYKqu0UVDZRWuskyWYhOymKrKQoeiVH0Ts5muzkKLKTo7BF\ndH1JUkSod3opa3BR3uCitMFJab2L0non1XYPSVFmesRa6RFrJSPWSmZcBEnR5nZ363V7fHvachxe\nGpwevD78Df8RRmyR/u7aXVX91x10wNC6hdPlZkd+BZtzS9maW8qWvDI255ayeWcJFpMpGDwG9Uln\ncJ90BvdNJys94aj6a+5Q4vH6KKpxsKuyifzKJvIqmsiraCS3ooldlU1EW430SYmmT0o0fVNi6JcW\nTd/UGLISI7u8qgv8AbCswUVxnZPiOidFga3J5SUj1krP+Ah6xUeQlRBBWoy13V/yrkAQqbf7e9Q1\nOr1EWozERZqIizJhiwy9DeZQpgOGdkgREUoq6tiSW8qmnSVs2lES3Dc0ORnUJ52h/dIZ0q8HQ/v2\nYGj/DHqmxetA0o1EhNJaJzvLG9lR3sCOskZ2lDWyvayB8jon2clRDEi3MSAthgHpMQzqYSM7Obpb\n/hJvcnkpqHVQUONgV42D/BoHtQ4PWfER9E6IpHdiJH2TIrFZ2zc22eeTQFdsD7VN/gGotkgTCdFm\nEqLN3VKtF046YGiHjeq6JjbuKGbTjhLWbytm4/ZiNu4oxu5wM7R/D4b168GwARkM75/JsAEZJMTq\n7sDdze7ysqOsga0lDWwpqWdLSQObi+upanTRPy2GIRk2hmbEMjgzliEZtm6p2mpyecmrtrOzyk5u\nlX8fF2GiX3IU/ZOjGJQSTWxE+wKIx+vzN+Q3+hv0I8wGkmPMJNksh2Xw0AFDO+xVVDewYXsR67YW\nsX6bf79hezHxsVGMGJDB8AGZjBiYyahBPemXlYKxG6pItL3VO9xsKW5gY1GdfyusZ0tJPWlxEQzr\nGcuInnGMzIpjeFYc0e38a7+9fCIU1jrZXtHE1opGtlY0ER9hZlBqFEPTYuifHIW5Hf+GfCLU2T1U\n1rupanQTZTGQGmslMcZ82LR76IChHZF8Ph95RVWs2VLI2i0FrN1axJrNBZRX1TOsfwajB2cxanBP\nRg3qybD+Gbr77yHA4/WxvayRdQW1rMuvZU1+LZuK68lKjGR0djxjshM4pnc8/VJjwj4i/kB8IuTX\nONhU1siG0gaKap30T45ieHoMIzNs7aq+8olQ0+imrM5Fg8NLss1MWpy1S+Zx6wgdMLSjSm29nbVb\nClm9OZ9fNhXwy6YCtuWXMaBXKqOHZDFmSBZjhvRi5MBMIiNaTkOvdS2Xx8fm4npW5VWzKreGn3Or\nqbN7GNsngXF9ExjXN5ERWXHt+ou/vRpdXjaVNbCmqIGNZQ1kxFoZnRnLMZmx7aq6crp9lNU5Katz\nERNhJCM+AlvkobmChA4Y2lHP4XSzblsRqzfm8/OGXazauIvNuaUM6JXKMUN7cezwbMYOy2ZYvww9\nduQQUFrr4Ked1azYUcXyHVXkV9o5pnc8E/onMWlgMsMyY7usBOL2+thc1sjPhfWsK6knOyGScVlx\njMqwHXBOrdb4fEJZvYviaicWk6JXUuQhFzh0wNC0VjicbtZuLeSn9Xn8tH4XK9flkV9SxchBPRk3\nvDfjRvTmuBG9dVffQ0B1o4vl26tYuq2S77dUUNXoYtKAZE4cnEzO4FRSYq1dkg+Xx8faknqW76ol\nt8rO2J5xHN8nnp5xoS1BICJU1LvJr7ITbTXRKynikKmq0gFD09qott7Oqo27WLE2l+WBzWBQHDey\nDxNG9WHCyL6MGZKl20O6WVG1ne+3VrBkUwXfba4gKymSk4ekMmV4GsN7xnZJgK9ucrM0r4aleTUk\nRJo5pX8iIzNsIY3J8PmE4lonxdVOUmItZCVGdGnbTWt0wNC0dhIR8ooq+XHNTpavyWXZLzvYtLOU\nEQMzmTi6L8eP7sfE0X1JTtAj17uL2+vj553VfL2hjC/Xl+J0+zh9eBpnjExnXN/ETu+d5PUJa4rr\n+XprFQ0uDyf3T2RidnxIbS5uj4+dFXaanF76pUZ1azWVDhiaFkaNdicr1uaydPUOfli9neVrc8lI\nieeEY/oxaWx/TjimPz3TEro7m0clEWF7WSNfrC3h819KKKtzcuaodKaP6cHY3p1ftbijsokvt1RS\nUOvgjEHJTMiOD2m+q8oGF7nldpJtFrKSIrplBLkOGJrWiTweL2u3FvHdz9v4/udtfP/zdmwxEZw4\ntj8njR3ASeMGkpWuA0h32FHWwKe/lPDRz0U4PT7OHZvBr8dm0jslulM/N7fKzqebyimtdzF9aArH\nhlBN5vb62F7ahNcnDOwR3aW9w0AHDE3rUj6fj007S/n2p60sXrGVb1duJdYWQc64geQcN4iccQNJ\nSbR1dzaPKiLC+sI63ltRyEeriuiXGsPFE7M4c2Q61k7sFbe9ool315ZiNihmjEojKz6yzfktqHJQ\nXu9iUHo00e0chd4eOmBoWjfy+Xxs2F7MohVbWLRsM9/+vI3eGUmcPH4QUyYOYdKYfno8SBdyeXws\nXF/Kmz/ms66glvOO7cllJ2TTK6lzppnxifBjXg0fbyhnVIaNc4elEtHGIFXZ4GJnmZ1+aVEkRHdN\nRwsdMDTtEOJ2e1m5Po9vlm1m4bJNrNlcwHEjenPa8UM5fdJQhvRN1914u0h+ZROv/7CLt5fnc2yf\nRK48sTcT+id2yu+/yeXlv+tK2VrexKXH9GBAG6vF6h0eNhc30j81ivguCBo6YGjaIayuwc6iFVv4\n6oeNLPhhIx6vlzMmDeWME4Zx8nGDsEWH1sdfC12T08N/fyrk5SW5xFhNXH9qP04bntYpXVzXFdfz\n5uoSRmf6SxttaaPoyqChA4amHSZEhK15ZXz+3Xq++G49y9fmctyI3px50nDOOmkEfXomd3cWj2g+\nn7BgXSn/XLgdu8vLTaf1Z/roHmEPHI0uL/9eVUyN3c01x/UkIergQaDe7mFzSSP906KIb0P69tIB\nQ9MOU/WNDr5etplPl6zl82/XkxQfzfSckUzPGcmxw3phMOhZeTuDiPDdlgqe/nwrTS4vt585kCnD\nUsNaVSUiLNxaxTfbKrliXCYD21BFVWf3sLWkkWE9Y9rcDhIqHTA07Qjg8/lYsS6PTxav5aNv1lBT\n38TZOSM599RRnDR2oJ4DqxOICF9vKOOJT7cQZTFy37lDGNM7vF2kN5U1Mn9lIWcOTuHEvgd/d0mN\nfxLD4T07Z0ZfHTA07Qi0Na+UD79Zw/sLV7Mjv5xpJ43gvNPGcOqEwVjMh9aEdoc7n0/470+FPPHJ\nZiYOSOKuswaTHh++tqWKRhf/+H4XE7PjOW1g0gFLMiLCttImDErRNzUy7A30OmBo2hEuv6SaD75e\nzXtfrmLTjhKmnTScC84YyynjB+uSRxg1Oj38c+F23vhhF9ef2o+rTuodtvXMa+1u/vF9PsPSYzhn\nWMoBA4HXJ6wrqCc9zkpaXHgnXtQBQ9OOIoVlNbz/1WreWfAz2/LKOPfUUVw4dSwnHNNft3mESV5F\nI7PeXked3c2jF45geM+4sLy30enh2aX59IqP4MJRB+5ebXd5WV/QwPCs8LZn6IChaUepvKJK3v7i\nZ97+fCWVtY1cOPVYZk4bx4iBmd2dtcOeiPDuikIe/3gTF47P4tYzBoRlDW+H28vfv9vFsPQYpg1J\nOWDawioHDU4vg3qEb6oTHTA0TWP9tiLe/HQFb362ksTYKC49ezwXTxtHqp6mpEPK653c89Yayuuc\nPP2b0fRN7fjMxXUOD08tyeW0gUkcf4BGdp9P+CW/nj4pkWHraqsDhqZpQT6fjyUrt/LaR8v4eNFa\nThzbnyt+NZEzJg3T7R3tJCK89sMunv58C3edNZgLx/fscGN0Wb2Tv36bxyXHZDAsff9BqLrRTV6F\nnZG9QluLY390wNA0rVUNTU7e+/JnXnl/KTsLKrh0+niu+vXx9Ot14KoQrXVbS+q5ef4qRvWK54/n\nD+vwxIY7q+w8vzSfm0/oRcZ+VvUTETYVNxIfZaJHGHpu6YChadpBbd5Zwsv//YE3Pl7O8AGZXH3+\nJM45eaTuohuiBoeHu95cQ3GNnWevPIYebZyhdn+W7arhqy2V3HlyHyz76ZG1uwF8ZC9bh9tRdMDQ\nNK3NnC43H3z9C/Pe+Y6teWVc9evjufr8SXpRqBCICHO/3sGr3+byzyuP6dBgPxHh1ZVFRJoNXDS6\nx37T7SxvwmhQ9ErqWIDSAUPTtHbZuL2Y597+lv98tpKc4wZxw8zJTBrTT8+m20Zfbyjj7jfX8PjF\nIzllaGq732N3e3n86538ekQaozJa76Rgd3nZUNjAmOzYDo0AP6wChlJqKvA0YATmicjj+9xPBl4D\n0gET8ISIvNLKe3TA0LQwqW908NpHy3j234uIibJy48wcLpg6Fqula9ZoOJytyq3m2pd/4q6zBjPj\nuJ7tfs/OyiZeWFbAnSf3ISGy9d/7xqIGkmMspMS2f32VwyZgKKWMwGZgClAIrABmisjGZmlmA1YR\n+UMgeGwG0kTEs8+7dMDQtDDz+Xx8+cNG/v76N6zfVsQNF0/mmhknkBjXuUueHu62lzZw5fMruPyE\nbH57ct92v+ezTeXsqnZw7cSsVu9XN7opqHIwIqv9XaXbGzC6Y0joccA2EckVETfwJnDuPmmKgdjA\ncSxQuW+w0DStcxgMBs44YRgf//MmPvzHjWzJK2P4OQ9xx5/fIb+kuruzd8jqlxbD2zdP5I2lu5i3\naEe73zNlQBKl9U42lja0ej8+yoTHKzQ4uv4rsTsCRiaQ3+y8IHCtuReAYUqpIuAX4NYuypumac2M\nGJjJC3+8jBVv34vJaGD8RY/xuwf/xeadJd2dtUNSenwEr18/nn99l8f8b3Pb9Q6z0cCvhqfx3tpS\nvL6WNShKKdLiLJTUOjuY29B1R1+6ttQh3QusFpEcpVQ/4Eul1CgRqd834ezZs4PHOTk55OTkhCuf\nmqYFZKbGM+f287jrmjN47q0lnHbN00weN5B7fjuVYf0zujt7h5SMhEheu348M5/9EbPJwMyJvUJ+\nx4geMSzeUcV3O6uZ3C+xxf2UWAur8+pxe31tWs1v0aJFLFq0KOR87Ks72jAmALNFZGrg/A+Ar3nD\nt1LqU+AREfk+cL4QuFtEVu7zLt2GoWndoKHJyfNvf8vf/rWQCaP6ct910xg+QM9d1VxueSMz//Ej\nD50/jNNHpIf8fFGtg79/t4v7TutHtKXl4MDtpU1EWAxkJoQ+kO9wasNYCQxQSvVWSlmAi4AP90mz\nCX+jOEqpNGAQ0P5KQU3TwiomysrtV0xhw0cPMXF0X6Zf/wyX3f2SrqpqpndKNM9dM5ZZb69jfUFt\nyM9nxEUwOtPGZ5vKW72fGmuhot7V0WyGpMsDRqDx+ibgC2AD8JaIbFRKXauUujaQ7FHgWKXUL8BX\nwF0iUtXVedU07cCiIi3cetmprPtwNiMH9eS0a57mtw/8i7wi/b8rwMgs//Qhv3vpJ0pqHCE/P21I\nCit21VJtd7e4FxNhxO0RnG5fOLLaJnrgnqZpYVNbb+fp+Qt5/u0lXDp9PHddcwbJCR2f2fVw9+xX\n2/h8TQlv3jiBKGtoTcfvrCnBYjRwzrCWgwK3ljQSF2kiNcQFlg6nKilN045QcbZIHrxxOj+9cx8u\nt4cx5/2Jv7z6FQ5ny7+QjybXn9qPAek27n93fcjPTu6byA+5NTg9LUsS8VFmapq6rnttuwKGUipC\nKRXeNQM1TTtipCfH8vQfLuLrV25j6eodjD7vT7z9xU8crTUCSikePn8Ya3bV8P7KwpCeTYmx0C8p\nkmW7alrci4syUWt34+ui32ubqqSUUgbgV8BM4Hj8gUYBXmAp8DrwflfXD+kqKU07PCxZuYW7n3yP\nqEgrf7n7AkYNav/0GYezjYV1XDZ3Oe/cMpHeKW0fOb+took3fi7mvtP6tlgPY21+PdnJkcRGtr2q\nq7OrpBYBY4EngL4i0kNE0oG+gWvjgMWhfrimaUeHk44dyHev3cWl04/j3Bv/wS2PvEllTesjmY9k\nQzJjufn0/tz62mpcrVQx7U+/pEgizAY2lLT8ncVHmahp6poqv7aWMKwicsBhhW1JE266hKFph5/q\nuiYe/ucnvLvgZ/548zlcds54DIajpzlVRPjdSz/RLzWGe84e3ObnVubXsjSvhptPyN7rep3dQ265\nfzW+turUEsbuQKCUenzfe7uvdXWw0DTt8JQQG8Vf7r6AD/5xA/Pe/Y4p1zzN+m1F3Z2tLqOU4vGL\nRvLeygLW5rd9fMaYzFjK6l0U1O7dPdcWYcTp8YVUYmmvUMP66a1cmxaOjGiadnQZPTiLRa/8LzOn\njWPqb//GQ//4GKfr6OhNlRhj4e7pg7nvnXWtzhfVGqNBMT47jpW79g4ySiniuqhaqk0BQyl1vVJq\nLTBIKbW22ZYLrOnUHGqadsQyGg389oITWfbWPWzYUcz4i+fww+rt3Z2tLnHesZlEWYy8/kNem58Z\nkxnLqqL6Fr3N4qNM1HZB99q2ljDeAM4GPgCmB46nA8eIyKWdlDdN044SGanxvPXkb3nwhun85q6X\nuOPP79Bk79ppL7qaUoqHZwznb19so6yubaPAM2KtGJUif59R4/FRZmqbPJ3ebbmtbRi1IpIL/Beo\nDhxfDsxTSh3TednTNO1o8uspY1j59r1UVDcwYeYclq3Z2d1Z6lT902K4eGIWf/pg48ET4w8yozNt\nrCrce+Jui8mAyaiwuzq3HSPUNoz7RaROKXUCcCrwEjA3/NnSNO1olRgXzSuPXsnsG8/motuf54G/\nf4jLfeSun3bjlP78klfDsm2VbUo/JjOW1UV1LUoTURYjTS5vZ2QxKNSAsTs304EXRORjQC/4q2la\n2J132hiWvfUH1m4p5OQr/8K2vLLuzlKniLQYuW3qQJ74bEubqpR6BuaNKtxnAaVIiwH7IRYwCpVS\nz+OfkvwTpVREO96haZrWJmlJsbz3t+v4zdnjOfnKv/Dq+0uPyOlFzj4mg7omN0s2VRw0rVKK0Rn+\nUkZzUdZDr4RxIfA5cLqI1AAJwJ1hz5WmaVqAUorrL57M5y/cwt9f/5qr75tPQ9ORNezLaFD8fuoA\nnvxsc5sC4u52jOZpoyxGmpyHVhuGHYjGP6cU+KujWs6IpWmaFmbD+mewZP6dRFhMTLr0/1i3NbRJ\n/A51Z4xIRwS+WFt60LS94iNwe30U1+0JnBFmA26vr83jOtoj1IDxLDABuCRw3hC4pmma1umiIi38\n88FLuevq0znzd3/ntQ9/7O4shY3BoLj9zIE8/fmWg37p+3tL+cdkNL8WYTZ2ajtGqAFjvIjciL+k\nQWAVPN3k7nzJAAAgAElEQVTorWlal7r07PEsmHcr//fSAn4/5z9HTC+qnCEpxESY+GjVwadKGdnD\nxvp9JiOMsho6tR0j1IDhUkoFVyNXSqUAXbc+oKZpWsCQfj349l93UlBcxZnX/p2SirqDP3SIU0px\n82kDeOGbHQdty8hOiKC03rnXwkqRFmOnjsUINWD8Hf/gvVSl1KPA98BjYc+VpmlaG8TZIvnPU7/j\n1PGDOeE3/8fPG3Z1d5Y67MRBydjdXlbsrD5gOrPRQGZcBLnV9uA1f8P3IVLCEJHXgLvxB4ki4FwR\n+U9nZEzTNK0tDAYD9157Jk/cOYNzb3yW/361qruz1CEGg+LySb2Z/23uQdP2TYpkZ2VT8LyzB++F\nFDCUUo+LyEYReSawbWxtynNN07Su9qtTR/Phszdy5xPv8vi8zw/r8RrnH5fJ91sqKa6xHzBd38Qo\ntlfuSWMxKXw+we3tnGopPb25pmlHjDFDslgy/w4+/PoXrnvoddzuzh3I1llsEWbOOSaDfy/NP2C6\nPkmR5FbZg2t6K6U6tR1DT2+uadoRJSM1ni/m/Z6yynouuO25w3aQ3+UnZPPmj7twevYf9GxWE7ER\nJoqajcfozBHfoU5v/iF7pjc/GxirpzfXNO1QExNl5e2nfkePlDjO+J+nKa08/HpQ9UuLYXCPWD5d\nXXLAdH2TItnRrB0j0mLE3kkN3yFNby4iF4tIXuA4V0TaNr2ipmlaFzOZjDz7wCVMPXE4p171FHlF\nVd2dpZBdcWI2rx6k8btvYhQ7Kpv3lOq8sRh64kBN045YSinuv/4srrv4JKZc8xRbcg8+7cahJGdI\nKmV1TraVNuw3zb4lDH9PKV+nNPrrgKFp2hHvpktO5v7rpnHGb//Kms0F3Z2dNjMaFGeOSufjA4z8\nTo2x4PYK1YE1vc0mA0qB29vNAaO1LrS6W62maYeDy8+dyBN3zeDsG/7BT+vbvo52dzt7TAYfry7e\nb4lBKUWfpEh2VO1TyuiEdgzdrVbTtKPG+acdwzP3zeS8W+YeNkFjVK84XB4fG4vq95umq9oxdLda\nTTtC5ebmcvLJJxMdHc2QIUNYuHDhAdPffffdJCcnk5yczD333LPXvfXr15OTk0N8fDxZWVn86U9/\n6sysd6qzTx55WAUNpRTTR/c44ISE+7ZjWEwGXJ7uq5LS3Wo1rZOJSFgbKmfOnMnYsWOpqqrikUce\nYcaMGVRUtL6i23PPPccHH3zAmjVrWLNmDR999BHPPfdc8P5ll13GiSeeSHV1NYsXL+bZZ5/lo48+\nClteu1rzoLFq44EHxx0Kzh6TwScHqJbKio+grMGFIzBQ0Ww0dM5o793/SA/HzZ99TWu/7Oxseeyx\nx2To0KGSkJAgV111lTgcjuD9559/Xvr37y+JiYlyzjnnSFFRkYiIPPDAA3LzzTeLiIjL5ZKoqCi5\n8847RUSkqalJrFarVFdXi4jI0qVLZeLEiRIfHy+jRo2SRYsWBd8/efJkmTVrlhx//PESGRkp27dv\nD8vPtXnzZrFardLQ0BC8dtJJJ8ncuXNbTT9x4kR54YUXgucvvfSSTJgwIXhutVpl48aNwfMLLrhA\n5syZE5a8dqf3F66W3lP+IBu3F3d3Vg7I5/PJqY8tklW51ftNM2fhdtlZ2SQiItUNLtlQUL/ftIHv\nzpC/c0Nt9I5QSl2qlJqllHowsD0Q/jCmaV3njTfeYMGCBWzfvp0tW7YEq1u+/vpr7r33Xt5++22K\ni4vJzs7m4osvBiAnJ4dFixYBsGLFCnr06MGSJUsAWLp0KUOGDCE+Pp7CwkKmT5/OAw88QHV1NU88\n8QTnn38+lZV7hjC99tprzJs3j4aGBnr16tUif9OnTychIaHV7Zxzzmn1Z1q/fj19+/YlOjo6eG3U\nqFGsX7++1fQbNmxg1KhRwfORI0fulfb000/n1VdfxePxsGnTJpYuXcqUKVPa8us9pJ17yigeufVX\nTL/hGXILD76edndpS7VUaoyV8gYXAGaT6pQShinE9B/gX5L1J8AR9txoR6XIMTeF5T32Vc+E/IxS\niptuuonMzEwAZs2axc0338zDDz/M66+/zjXXXMPo0aMBeOyxx0hISGDXrl1MmDCBrVu3UlVVxbff\nfss111zDs88+S2NjI4sXL2by5MmAPxhMmzaNqVOnAjBlyhSOPfZYPvnkEy6//HKUUlx55ZUMGTIE\n8M+8uq+PP/445J+roaGBuLi4va7FxsZSWNj6sqb7po+NjaWhYU/f/6eeeorTTjuNJ554Aq/Xy4MP\nPsjYsWNDzteh6JLpx1HbYOes65/hqxdvo0dK3MEf6gbTx2Rw2T+Xce85QzAaVIv7qTYLpQ3+KUL8\nVVLhb8MINWBkisgZYc+FdlRrzxd9OGVlZQWPe/XqRVGR/6+44uJijj322OC96OhokpKSKCwspFev\nXhx77LEsXryYJUuWMGvWLFavXs3333/PkiVLuOWWWwDIy8vj7bff3qu+3+PxcMopp7T6+eESExND\nXd3e02HU1NQQGxvbpvS1tbXExMQA0NTUxCmnnMIjjzzCJZdcQnFxMTNmzCA1NZXrr78+7HnvDtdf\nPJm6Bjvn3PgPvnrxNuJskd2dpRb6p8WQGGNh5Y4qxvdPanE/LcbCmmJ/TyqTUeHx+quRlGoZXNor\n1G61PyilRnb0Q5VSU5VSm5RSW5VSd+8nTY5SapVSap1SalFHP1PT9mfXrl17He8ubWRkZJCbmxu8\n19jYSGVlZfD+5MmTWbhwIatWrWLcuHFMnjyZzz//nOXLl3PSSScB/gB02WWXUV1dHdzq6+u56667\ngu892P/QZ555JjabrdXtrLPOavWZYcOGsWPHjr1KCb/88gvDhg3bb/rVq1fvlXb48OGAv3qrvr6e\n3/zmNxgMBjIzM7nooov49NNPD5jvw81d15zBpGP6M/OOeYfskq9ThqWxaFN5q/dSYyyUBaqkDEph\nNPiDRliF0uABbATcwBZgbWBbE+I7jMA2oDf+9cBXA0P2SRMPrAd6Bs6T9/Oug7cWadoBZGdny8iR\nI6WgoEAqKytl0qRJMmvWLBER+eqrryQlJUVWr14tDodDbrnlFjnxxBODzy5YsEBsNptMmTJFRETW\nr18vNptNhg8fHkyTn58v6enp8sUXX4jH4xG73S7ffPONFBQUiIhITk6OzJs3r1N+tgkTJsgdd9wh\ndrtd3n33XYmPj5eKiopW086dO1eGDBkihYWFUlBQIEOHDpXnnntORESqqqokOjpa3njjDfF6vVJc\nXCwTJkwI/p6OJB6PV2b8/jm5etar4vP5ujs7LSzbXilnP/ltq/eaXB65/YON4g3ke3VurTQ6PK2m\npZ2N3qEGjN6tbNkhvmMi8Hmz83uAe/ZJcwPwxza8q42/Zk1rXe/evWXOnDkydOhQiY+PlyuvvFLs\ndnvw/ty5c6Vfv36SmJgoZ599thQWFgbv1dfXi9lslj/+8Y8i4u/JkpqaKjfccMNen7Fs2TKZPHmy\nJCYmSkpKikyfPl3y8/NFxB8wXnzxxU752XJzcyUnJ0ciIyNl8ODBsnDhwuC9JUuWSExMzF7p77rr\nLklMTJTExES5++6797r36aefypgxYyQ2NlbS09Pld7/73V6/pyNJY5NTTrrsz/LgMx92d1ZacHm8\nMuIPX0hlvbPV+/d+ukUqG10iIrK+oF5qAsf7am/AUP5n204plQAMACKalVKWhPD8DOAMEflt4Pw3\nwHgRublZmqfwlz6GATbgryLyr1beJaHmX9Oa69OnDy+++OJebQqaVl5VT86VT3LP/0zlsnMmdHd2\n9vI/81byq7EZTB+T0eLe377N4/RBSQxOjWFrSSMJ0WaSbZYW6ZRSiEjIjRuhdqv9LbAEWAA8BHwB\nzA7xM9vyDW8GjsE/7cgZwP1KqQEhfo6maVq7pCTaePfp65j19PssXb2ju7OzlxMGJvPdlta7ADdv\nxzAbw9+1NtReUrcC44ClInKyUmow8FiI7ygEmncLyQL2nT4yH6gQETtgV0otAUYBW/d92ezZs4PH\nOTk55OTkhJgdTdO0lgb3Tef5hy7j0rteZPH8O8hKT+juLAEwaWAS8xbvaLUHVGqMhbL63QHDgDsw\nPciiRYuC44Y6IqQqKaXUShE5Vim1GpggIg6l1AYRGRrCO0zAZuBUoAhYDswUkY3N0gwGnsFfurAC\ny4CLRGTDPu/SVVKapnWqp+Z/xX8+W8nCl24nKrJl9U5XExGOf+hrXr9hPH1TY/a6t66knsXbq7lx\nUi/Kap3UO7z0S4tq8Y4uqZICCgJtGO8DXyqlPgRyQ3mBiHiAm/BXZ20A3hKRjUqpa5VS1wbSbAI+\nxz+x4TLghX2DhaZpWlf4/WWnMqx/Btf/8XUOhT9QlVKcMKj1aqnUGGuwSspkCv98UiE3egcfVCoH\niMXf48kVzkyFkAddwtA0rdPZHS5yrniSq847nusumtzd2eH9nwr5/JcS5l6992h7r0+446PNPD59\nIE63j9xyOyOybC2e76pG71cDJQxEZBHwHTA31A/VNE07nERGWHjjz//Do899xvK1ud2dHSYNTObH\n7ZV49ilBGA2KpCgzFQ2uTmn0DrVKapSIVO8+EZEq/L2ZNE3Tjmj9eqXwzP0z+c1dL1JZs/81trtC\nis1KRnwka/JrW9xLjbFQ2uAKNnqHsxYm1IChlFKJzU4S8Y/c1jRNO+Kdc/IoZpwxlqtnzcfn64T1\nJkJwwqBkvtvcSjuGzd+11mhQKAXhLGSEGjCeBJYqpR5WSv0JWAr8OXzZ0TRNO7T98aazqW2w8/fX\nv+nWfBw/IIml2ypbXE+Lse7pWhvmhu+QAoaIzAfOA8qAEuC8wDVN0w4xoSzR+tRTT9GvXz9iY2NJ\nS0vjqquuor5+zxrS999/PyNGjMBsNvPQQw91RfYPWSaTkVcevYInX/6S1Zu6b7W+Ub3iWV9Qh8+3\nd5WTf/De7mnOVVinOQ+1hIGIrBeRv4vIMyLS+mosmqaFTCS89c2hLNF67rnnsnLlSurq6ti0aRO7\ndu3ikUceCd4fMGAAf/7znznrrLPCOl324ap3ZjKP33E+V977Ck32bukkSkK0hYRoMzvLG/e67l8X\nw4WIhH2p1jYFDKXUd4F9fStb3cGe17RDVe/evZkzZw7Dhg0jMTGRq6++GqfTGbz/wgsvMGDAAJKS\nkjj33HMpLi4G4MEHHwyueeF2u4mOjg5OWW6324mIiKCmpgaAH3/8keOPP56EhARGjx7N4sWLg+/P\nycnhvvvuY9KkSURHR7Nz586w/Fxbtmxh1apVPPTQQ1itVs477zxGjhzJu+++22r6vn37kpDgH8ns\n8/kwGAz06NEjeP/yyy9n6tSp2Gy2Q2IswqFg5rRxHDO0F3c/2frvtCuM7BXfouE7xuJvVm5webut\nhJEb2N8vIrZ9ttZXZNG0w4ReonXP7yEuLo6UlBRSUlK49dZbQ/tFHoWevudCvvpxE58sXtstnz+i\nZxxr9wkYSinSAlOEmI0Ktyd8JYy2ziV1jFIqA7haKdWizSLQvVbT2qXv7eFZiGfHX6aF/IxeonWP\nSy65hEsuuYRt27ZxwQUX8NRTT3HbbbeF/NlHk9iYSF546DIu/8PLTBzdl8S46IM/FEYjsuL4cl1p\ni+spgRHffeIjaXJ5w/Z5bQ0Yc4GFQF/863nvq0/YcqQdddrzRR9OeonWvfXv35977rmHOXPm6IDR\nBieM7c+vTh3NnX9+lxf/dHmXfvbwrFg2FtXh8fowGff8sZEcbaayycXA5Cjc9i6ukhKRv4nIEOBl\nEemz7xa23GhaN9BLtLbkdruJimo5aV1b8ns0eviWc1j6y44ur5qyRZhJj49ga+neAwljI0zUO72Y\njYYWo8E7ItRutdeF7ZM17RAgIjz77LMUFhYGexNddNFFgL+X0csvv8wvv/yC0+nk3nvvZcKECcF2\nhsmTJzN//nyGDRuG2WwmJyeHefPm0bdvX5KSkgD4zW9+w0cffcSCBQvwer04HA4WLVq0V9XQwRqR\nP/vsM+rr61vdPvnkk1afGThwIKNHj+ahhx7C4XDw3nvvsW7dOs4///xW08+bN4/ycv9a0Rs2bGDO\nnDl7pfV4PDgcDrxeL263G4fD0e0D1w4l0ZFW5j54Cbc8+hbVdU1d+tmttWPEWIzUOzxhb/Ru67Kq\n3wf2DUD9Pltde5b6C8eGXqJV6yC9RKvfVVddJWlpaRITEyMDBw6Uxx9/fK81ra+44gpRSu21vfrq\nq52S78PZ7x97S377wPwu/cyXFu+Q+95eu9e17RWN8sQ3O8Xt8crybdUtnqGrlmg9lOjZarWO0ku0\nauFU3+hgzPl/4uU/XcGJx3bNIqErd1bx8Psb+eC2ScFrZQ0unv1+Fw+e3o/l22sZ1zcOg2FPVWJX\nrYehaZqm7YctOoIn7pzBLY++hcvt6ZLPHJoRy9aSepyePb2hbFYj9U4PSinMpvBVS4U6vXmEUupS\npdQspdSDge2BsORE0zTtCHDuKaPo0zOJv87f/1Qs4RRlNZGdHM2W4j0N3xEmAz4Bp8cX1tHeoZYw\nPgDOAdz42zMagMYDPqFph7CdO3fq6igtrJRS/OXuC/jra1+TW9j6VCzhNiIrjjX5NXvlYXcpI5wN\n320dh7FbpoicEZZP1jRNO0L1zkzm1stO5bY5b/Pe367r9K7II7Ja9pSyWU00BLrWdlcJ4wel1Miw\nfLKmadoR7NbLTmFnYQWff9f5c7SO3E/ACJYwPN3QhgGcCPyklNqilFob2NaEJSeapmlHEIvZxOO3\nn8fdT76H2x2+6TlaMzjDxs7yRuyulg3f4VyqNdSAMRUYAJwGTA9sZ4clJ5qmaUeYM04YRu/MJJ57\ne0mnfo7VZKRXUtReU53bdo/2Nhm6rQ1jBiDA7go5AWqVUj+JyOqw5EjTNO0IMuf285j6278yc9o4\nkuJjOu1z+qREs7O8kaGZ/vnCYqxGqprcYW30DrWEMRa4DsgAMoFrgTOBF5RSd4clR5qmaUeQof16\ncN5pY3jkuc869XN2B4zdYq0m6h1eTN3Y6J0FHCMi/ysit+MPIKnAZODKsORI07SwCGWJ1tmzZ2M2\nm4OTGsbGxgYnXty1a1eLSQ8NBgNPPfVUF/0kh7/7rjuLtz//iU07SjrtM/qk7h0wYgKN3iaDwuvr\nnhJGCtB8PUI3kCYiTYAjLDnStKOU7JkjLSxCWaJVKcXMmTODkxrW1dXRu3dvwD/jbvMJD9euXYvB\nYNjvRIZaS8kJMdx+1RQeeObDTvuMPinR7CxrXsLwN3obFIQpXoQcMF4HlgVGeM8GfgDeUEpFAxvC\nkyVN6zp6iVa/UILVq6++yuTJk1tdHVDbv+suPImf1u9ixbrcTnl/n5RodpQ3BP87xgTGYSiDQrqj\nhCEiDwO/A2qBauBaEXlIRBpF5NKw5EjTupheotVfwvjoo49ISkpi+PDhzJ07t9V0IsL8+fO54oor\nDvJb1fYVGWHhD7+dyuxnPjp44nZIjLYAUN3oBvyN3k1uLyISthJGqL2kEJEVwIrwfLymwc3/3RiW\n9/z910NCfkYv0ep34YUXcu2115KWlsaPP/7I+eefT3x8fDBA7vbdd99RVlbGjBkzQs6TBlecO5Gn\n5n/FouWbyTluUFjfrZQKNnwnxlgwKEWU2UiTy0tgPvMOjzhvU8BQSn0vIpOUUg34u9I2JyJy8HUf\nNW0/2vNFH056iVaCAQtg4sSJ3HrrrbzzzjstAsarr77KjBkz9rsan3ZgZrOR+68/iwef+YhFrw4M\n+5QhuwPG2D4JQGC0t8sbbMcwdvDj2rpE66TAPkZEbPtsOlhohzW9RGvb2O123nnnHV0d1UEXnjGW\nJoeLT5esC/u79+1aa4sw0hCY5twXjg4Voay2BFwA2ALH9wPv4e9mq1fc0w5L2dnZMnLkSCkoKJDK\nykqZNGmSzJo1S0REvvrqK0lJSZHVq1eLw+GQW265RU488cTgswsWLBCbzSZTpkwREZH169eLzWaT\n4cOHB9Pk5+dLenq6fPHFF+LxeMRut8s333wjBQUFIuJfcW/evHmd8rNNmDBB7rjjDrHb7fLuu+9K\nfHy8VFRUtJr2/fffl6qqKvH5fLJs2TLJyMiQ+fP3Xjnu9ddflz59+nRKXo82H33zi4y74NG9VjUM\ny3t/LpTrX/4peP7S8gJZnlcjK3fUiNPtDV6nnSvuhdpL6gERqVdKnQCcCrwEtN46pmmHAaUUl1xy\nCaeffjr9+vVjwIAB3HfffQCceuqpPPzww5x//vlkZGSwc+dO3nzzzeCzEydOxOFwBEsTQ4YMITIy\nMngO0LNnTz744AMeffRRUlNT6dWrF08++eRePZI6aybTN998k5UrV5KYmMisWbN49913g2uNf/vt\nt9hstmDat956iwEDBhAbG8sVV1zBH/7wBy677LK93jd//vwW17T2OWvyCAwGFfZSRsvBe7u71oan\nhBHSEq1KqdUiMlopNQdYKyKvK6VWiciYDuekHfQSrVpH6SVate7y7pc/89d/fc3iV/83bH80NDo9\njHvgK9Y9dgYGg+KLzRU43D56x0YwIC2aKKsR6LolWguVUs8DFwGfKKUi2vEOTdO0o96vThlNbb2d\nRcu3hO2d0VYTcZFmimv946hjA6O9w1XCCPXL/kLgC+B0EakBEoA7O5wLTdO0o4zRaODOq07j8Re/\nCOt7ezerlooJru0N4aiMCXXgXqOIvCsiWwPnxSKyoOPZ0LTuoZdo1brTRWeOI7ewkh9/2RG2dzZv\nx/AvouTtthJGWCilpiqlNimlth5ollul1DillEcpdV5X5k/TNK0rmM1Gbr9iCv8XxlJG8zmlbFZ/\nt9pwzSfV5QFDKWUEnsG/GNNQYKZSqsXIrUC6x4HP2bP+hqZp2hHl8nMnsHpTAWu3tD4KP1StlTCU\noutKGEqpbKVUr7ZsbXjdccA2EckVETfwJnBuK+luBt4Bytv802iaph1mIqxmrr/oJJ5545uwvK9P\nSjS5gYBhMRkwGhQi4SlhtHUuqVdpOSXI/px8kPuZQH6z8wJgfPMESqlM/EHkFGBcCJ+taZp22Ln6\n/EkMP+chHr7lXFITbQd/4ACykqIornHg8viwmAzYrEa8PgnLjLVtChgiktPhT2r2ujakeRq4R0RE\n+Tso6yopTdOOWEnxMZx/+jG88Pa3zLp2WofeZTEZ6BEfQX5lE/3SYrBZTbh94ZmxNuTZasOgEP/K\nfbtl4S9lNDcWeDMwmCUZOFMp5RaRFquPzJ49O3ick5NDTk5OmLOraZrW+W6cmcOZ1/6NO646DavF\n3KF37e5a6w8YRn78YQlbV/2ILbKDX/ntmU+kIxv+ILUd6A1YgNXAkAOkfxk4bz/32j7JiqYdZXbu\n3Ck5OTkSFRUlgwcPlq+++mq/aaurq+Xyyy+X1NRUSU1NldmzZ7f7XVr7nX3DMzL/g6Udfs9D762X\n57/eLiIi//65SL7bWiX5lfbgfbpoLqkOExEPcBP+AYAbgLdEZKNS6lql1LVdnR9NO1TInj+EwiKU\nJVpvu+02HA4HeXl5LF++nH/961+88sor7XqX1n43XXIyz7z+TYf/HWQkRFBc4x/tbbOacHp94fm3\n1Z4oc6hs6BKG1kHZ2dny2GOPydChQyUhIUGuuuoqcTgcwfvPP/+89O/fXxITE+Wcc86RoqIiERF5\n4IEH5OabbxYREZfLJVFRUXLnnXeKiEhTU5NYrVaprq4WEZGlS5fKxIkTJT4+XkaNGiWLFi0Kvn/y\n5Mkya9YsOf744yUyMlK2b98elp9r8+bNYrVapaGhIXjtpJNOkrlz57aaPjk5WVasWBE8f/TRR4Mz\n84b6Lq39fD6fjP71w7J4xeYOvefDnwvlxlf8s9Yu2lYpn64rk9zypuB9OruEoZQ6MbCW9z+VUs8E\njk/veMjStO6ll2j1k2Z/gfp8PtatW9fud2nto5TihpmT+ce/Fx888QGkxUZQWudfmz42woTd4wvL\nOIy2rrh3L2AGVgGNgBGIBU5VSp0iIvd0OCfaUevHbTVhec+E/vEhP6OXaPWbOnUqjz/+OK+88gol\nJSW89NJL2O32dr1L65iLp43jgb9/SHF5LT1S4g7+QCtSY62UBSYgjLEasbu9+Hwdz1tbm8zXSSs9\nlIB3lFJ6cV+tQ9rzRR9OeolW+Nvf/sbNN9/MgAEDSE5O5pJLLuHf//53u96ldYwtOoLzThvDax/+\nyJ3XnNGud6TGWimtcyIi2CwmmlzeLp1LapRS6gGl1HSl1ClKqclKqWlKqXuAiR3OhaZ1I71EKyQk\nJPDaa69RXFzM2rVr8Xq9jB8/vl3v0jru6vMm8dJ/f8DXzmJBlNWE1WSgtslNhNkQaPQOQ8ba2tgB\nTAEewD8P1LPAg/hX3VPtaTwJx4Zu9NY6SC/R6rd9+3apqKgQj8cjn376qSQnJ8uGDRva9S6t43w+\nn0y4+DH58ocNB0+8H6fNWSybiuqk0emRJ77eIRsL64P3aGejd3u/qOMD+4T2PB+uTQcMraN69+4t\nc+bMkaFDh0p8fLxceeWVYrfv6a8+d+5c6devnyQmJsrZZ58thYWFwXv19fViNpvlj3/8o4j4/ydP\nTU2VG264Ya/PWLZsmUyePFkSExMlJSVFpk+fLvn5+SLiDxgvvvhip/xsubm5kpOTI5GRkTJ48GBZ\nuHBh8N6SJUskJiYmeP6f//xHMjIyJCoqSsaMGSMLFixo87u0zvH8f5bIxf/7Qrufv/TZH2XJpjJx\ne33yyIJtsr6g4wEjpCVad1NK3Soif92973Axp530Eq1aR+klWrVDVV2DnUHTHmDVe/eTnhx6e9Ht\nr6/m+AHJzDiuJw8v2MavhqQxIss/T1VXLdGqaZqmdYHYmEh+PWUM//rwx3Y9nxobQVmdv6eUQSl8\nYZhMSgcMTdO0Q9TV50/i5XY2fqfFWSkLjMUwGBTeLuwlpWlHJL1Eq3YoGzu0F5FWM0t/2Rnys2mx\nEXtFxLMAABBJSURBVJQGxmIYD9cV9zRN07S2UUpx8bRx/PuT5SE/m9qshGE0HMZremuapmltc/G0\ncfz3q1U4Xe6QnvNPD+IvYZiNhrCMw2hvwPhyn72maZrWCbLSExg+IJPPvg1t7q7UWCvldU58PsFo\noPsChohsaL7XNE3TOs8lZx3Hvz9dEdIzVrORaKuJ6iYX5sAcZR0dhqCrpDRN0w5xvzp1NIuXb6Gq\ntjGk5/yTEDoxmwIBo4P5aPN6fUqpE4FTgHTAC5QDS0VkQQfzoGmaph1AnC2SKccP4b0vV/E/M05o\n83Npcf52DKtRIfjX9TaEPFxvjzaVMALTm5+Cf3rzd4APgfX4pzef0/6P1zRN09pi5rRxvBFib6nU\nWCultU4sxkAJo4N9a/X05pqmaYeB048fynWzX2dnQQV9eia36Zm0OP9o72ybFZGOj8XQ05trmqYd\nBsxmI+edNoZ3vvipzc+kxvrHYliMCp/Q4bEYbQ0Y64DvgWOA84CLgHHACuCODuVA0zRNa5NfTxnN\n+wtXtzl9WmxEIGAYgm0YHdHWKqlHgbfwB4jlIlK++4ZSKg6o7Vg2NE3TtIM54Zj+5JdUk1dUSXZG\n0kHTp8X5pwexmAx4XdJlJYyPgV3ADGCNUmqHUupNpdT/Ao91KAeapmlam5hMRqbnjOS/X7WtlLGn\nSsqAV6TDg/faGjDuE5GXROQq4A/AWOB5QAHjO5YFTdM0ra1CqZZKibVS2eDEZFB4fdLhKc7bGjAu\nbHYsIlItIl+LyBPAfR3KgaZpmtZmk8cNZGteKYVlNQdNazYaiI0043R7/QGji0oYjyulnlNKXQ6k\n7nNvTceyoGmaprWVxWzizP9v7/6DrCrvO46/P7t7WS7KRdC4UWARFQmrAUH5pUGhWkXTNmmbTGuM\nTWusTls7bTpJTTrTSCdtOjaJ7bR2TMb+GNuZxk4apxrHH+NMJDYqGJVfCigrGECIo6KiILB777d/\nnAt7ufvr7F72XPbyec0w3HPOs+d895l773ef5znPeZZcwAMpWxltE8bywYEuukvZjWHcCXwfOANY\nKGmtpMcl/R3w7ZoiMDOzIfn1K+em75Ya38re/UnCyGQMIyK+HRErI+KOiPhMRFwI3EByq+2k2kIw\nM7OhuGLRx9jwyuu88fbeQcu2TWjl3f3ZtjB6iYhfRMQDwF/VFIGZmQ3J2NYcV13awY+eGHxEoK0w\nlj3vH6SrVBr5MQxJrZL6nYceEU9XlG2vLRwzM0vj01dcyAM/Hrxb6vQJrez54BBdxQxaGBFxEFgk\n6TpJ+b7KSJoo6WZgWk3RmJlZKlcunsXq9a/xwf6DA5ZLZnsf4FCx9ttq0870/gnwJeArkqaSPN48\nV/5/P7ATuCciPOPbzCwD408ay8XnT+OJ1S/zq8tm91vu9AnJynvdpRLFUm3XTJswvkXy+I+pwGTg\n2ogY2koeZmZ2TF39iQ4ee+qlARPGaSe38tb7hwDRXaotY6Qd9N4QEbdFxOeA3yZ5+KCZmdXRNUsu\n4NGfvjTg0quFfI69H3YhBcWMZnof6SSLiN3A4PdymZnZiJox7XRax7Sw4ZXX+y1zUmszB7tLSKK7\nmE3C+KqkuyTdKGkuFUvDSmqrKQIzMxsWSVx96fk89tOXBixTyLdABMWM5mHcS/LE2nbgr4F/krRK\n0ndIxjfMzKwOrllyPo8MkDAACmNzBNQ86J12pvc3IuLRiFgREZ+MiDOB64HnSQbCh0zSckmbJW2R\ndFsfx6+XtE7SeklPSep/VMfM7AS15KIZvLhlF3ve6/8+pPH5liRh1HGm96sR8V/A7UP9WUnNwF3A\ncqADuE7SrKpiW4HLImI28A2Sx6mbmVmFsa05LrvoXB5/elO/ZQr5HKUMH2/er4h4chg/tgDojIjX\nIqILuA/4VNV5n6mY17EamFJbpGZmjWl5+W6p/hTyOYpZPXxwBEwGdlRs7yzv688XgYdHNCIzs1Hq\nqks7ePzpjRT7GaQo5FvoLtb+LKm0E/eOtdRhS1oG3Ahc2tfxFStWHHm9dOlSli5dWmNoZmajS/sZ\nk2g7tcALm7Yz/4Kzeh1/q3MNzz/0PC9OzPPYpD6f8JSKBprwMVIkLQJWRMTy8vbXgFJE3FFVbjZw\nP7A8Ijr7OE/UI34zs+PNn93xA874yAS+cuNVvY798+OdvNNd4pfO/QiXzJiIJCJCQ71GvbqkngNm\nSDpL0hiSmeMPVhYoP/n2fuDzfSULMzPrsWzBTH7ys1f6PDY+38KBQ8War1GXhBER3cCtwGPARuC/\nI2KTpFsk3VIu9nVgInC3pDWSnq1HrGZmo8GSi89l9fptHDzU1etYIZ/jw2OQMOo1hkFEPAI8UrXv\nexWvbwJuyjouM7PR6JTx45g5/aM8u/41llw846hjhXyOfQe7GHIfVJV6dUmZmdkxtnT+eazso1uq\nkG9h34GkhVHLuK8ThplZg1i28DyeePblXvuTFkY3ATXNxXDCMDNrEIvnnMP6l3f2WoWvMDbHvgPd\nlIKalml1wjAzaxDj8mOYO6udp144+sbSCeNyvH+gi1JETZP3nDDMzBrIsgUze41jjM010dVdohjh\nFoaZmSUuX3Ber/kYksjnmimW3MIwM7Oy+RdMo3P7m70ed54f00R3KYgaMoYThplZAxmTa2HxnLN5\n8rktR+0/KdfiFoaZmR3tsvkzeiWMk1tb6C55DMPMzCosvvBsVq/betS+Qr6F7giK7pIyM7PD5s1q\nZ/O2N46aj5GsugddRScMMzMrG9ua4+PnTeb5l35+ZF8hn6MUwaF+FllKwwnDzKwBLb7wbJ5Z29Mt\nVcgng95uYZiZ2VEWzalOGDm6S0G3E4aZmVVaNGc6z27YRqmUdEEV8jmKpRJdJXdJmZlZhbZTC5x6\nysls2voLIOmS6iq6hWFmZn1YNGc6q8q31xbyObqKJd9Wa2ZmvVUOfBfyOQ6VH0A4XE4YZmYNavGF\n5/BMuYUxIZ/jQFeRGu6qdcIwM2tUH5vexp539/HG23sp5Fs40OUuKTMz60NTUxMLZ09n1bpttOaa\nkzUxnDDMzKwvyTjGqwCUokS3E4aZmfXl6Al8ouSEYWZmfbno/Gm8uGUXXV1FBF4Pw8zM+nbyuFba\nz5jExq27aVIQvq3WzMz6M6+jnTUbt9PS3IRqOI8ThplZg5vbMZU1m3aQaxbS8FOGE4aZWYOb19HO\n8xu309rcRA35wgnDzKzRzZk5hU2v7qa1pYkmtzDMzKw/J+VbmXbmJPbv/5BmtzDMzGwg8zra2bNn\nr1sYZmY2sLmz2tm5802am5wwzMxsAPM62tn62q6azuGEYWZ2Apgzcwpbt+32wwfNzGxg4/JjmHza\nyaNvASVJyyVtlrRF0m39lPnH8vF1kuZmHaOZWaOZN/PM0dXCkNQM3AUsBzqA6yTNqipzLXBuRMwA\nbgbuzjrO0WblypX1DuG44bro4bro4bqAizqmjrrHmy8AOiPitYjoAu4DPlVV5teAewEiYjVwiqS2\nbMMcXfxh6OG66OG66OG6SAa+R1vCmAzsqNjeWd43WJkpIxyXmVlDm33eFIo1LOpdj4SRNr1V3yxc\nw1PczcxsXH5MTS0M1fJs9GFdUFoErIiI5eXtrwGliLijosx3gZURcV95ezNweUS8UXUuJxEzs2GI\niCHP4GsZiUAG8RwwQ9JZwC7gt4Drqso8CNwK3FdOMO9WJwsY3i9sZmbDk3nCiIhuSbcCjwHNwL9G\nxCZJt5SPfy8iHpZ0raROYB/we1nHaWZmR8u8S8rMzEanUTHT2xP9egxWF5KuL9fBeklPSZpdjziz\nkOZ9US43X1K3pN/IMr6spPx8LJW0RtKLklZmHGJmUnw+TpP0qKS15br43TqEmQlJ/ybpDUkbBigz\ntO/NiDiu/5F0W3UCZwE5YC0wq6rMtcDD5dcLgVX1jruOdbEYmFB+vfxErouKcj8GHgJ+s95x1+k9\ncQrwEjClvH1aveOuY12sAP72cD0AbwMt9Y59hOpjCTAX2NDP8SF/b46GFoYn+vUYtC4i4pmIeK+8\nuZrGnb+S5n0B8MfA/wBvZhlchtLUw+eAH0bEToCIeCvjGLOSpi52A4Xy6wLwdkR0ZxhjZiLi/4B3\nBigy5O/N0ZAwPNGvR5q6qPRF4OERjah+Bq0LSZNJvjAOP1qmEQfs0rwnZgCTJD0h6TlJN2QWXbbS\n1MU9wPmSdgHrgD/JKLbj0ZC/N+txW+1QeaJfj9S/k6RlwI3ApSMXTl2lqYt/AL4aESFJ9H6PNII0\n9ZAD5gFXAOOAZyStiogtIxpZ9tLUxV8AayNiqaRzgMclzYmI90c4tuPVkL43R0PCeB2YWrE9lSQT\nDlRmSnlfo0lTF5QHuu8BlkfEQE3S0SxNXVxEMpcHkv7qayR1RcSD2YSYiTT1sAN4KyI+BD6U9CQw\nB2i0hJGmLi4B/gYgIl6VtA2YSTI/7EQz5O/N0dAldWSin6QxJBP9qj/wDwK/A0dmkvc50a8BDFoX\nktqB+4HPR0RnHWLMyqB1ERFnR8T0iJhOMo7xBw2WLCDd5+MB4BOSmiWNIxng3JhxnFlIUxebgSsB\nyv31M4GtmUZ5/Bjy9+Zx38IIT/Q7Ik1dAF8HJgJ3l/+y7oqIBfWKeaSkrIuGl/LzsVnSo8B6oATc\nExENlzBSvie+Cfy7pHUkfzD/eUTsqVvQI0jS94HLgdMk7QBuJ+meHPb3pifumZlZKqOhS8rMzI4D\nThhmZpaKE4aZmaXihGFmZqk4YZiZWSpOGGZmlooThpmZpeKEYWZmqThhmKUgqbUO1xyb9TXNBuKE\nYVYm6S8ldUq6SdIfSfqupPGSfgUYX1Huy5J2SbpBUqukb0naJOnTNV6/SdKdFbumSPrlWs5pdiwd\n98+SMsvQz4C2iPgXAEl3AF8CtlQtOvQc8KOI+M9yuU3Avoj43+FeWNJEkmf5XH54X0R0SvqkpKcj\nYt9wz212rLiFYdZjEfBCxfbhRz9XJ4KFwFMV25cBz9Zy4Yh4JyLuBPZWHXoIuL6Wc5sdK04YZj3m\nA2sAJJ1BsiLZ1PI6EtXl2iR9QdIXgKupMWH0JyJeBS4YiXObDZW7pMx6zAXOLq/EVgCuBb7TR7lz\nI+IzcKQr6fbB1smW1AH0Nx5xb0S8O8CPNw8auVkGnDDMgHKS2B4RP6zan6va/ijwZsWuuZRbF5IW\nkgyOn1O9Hkd5/YnhrkHhu6XsuOCEYZaoHpc4rNhHuTUV2/NIBssBPhsRX5bUIak9IrYfLjRIC+M/\nBllKtzRw6GbZcMKwE56ky4DfB3ZIaqtapnJ/RblLgD8E9kiaDMwGrgNWS5oAjCsX/QBoA44kjMFa\nGJJOAm4GZkn6U5JV8fYpWTbxg2Pwa5rVzAnDTngR8SSwrJ/DOyVNLN/F9DTJAPdhrwOPHN6Q9F75\n5SnAkNaUL982+/flf5XmAKuGci6zkeK7pMwGdg/w2ZRlH5a0DChVdkfV6ErgB8foXGY18ZreZoOQ\ntAT4+TFMAmmv+3GgOSLWZnlds/44YZiZWSrukjIzs1ScMMzMLBUnDDMzS8UJw8zMUnHCMDOzVJww\nzMwsFScMMzNLxQnDzMxS+X+AzWa82XgbgwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7fc535544b50>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# betas corresponding to different values of power\n",
    "pnotKnotA = 1 - np.array([[0.08, 0.31, 0.57, 0.9]]).transpose()\n",
    "nump = pnotKnotA.shape[0]\n",
    "\n",
    "# compute posterior with previous assumptions\n",
    "pnotK = pnotKnotA * (1-pA) + pnotKA * pA\n",
    "pKnotA = 1 - pnotKnotA\n",
    "pnotAK = pKnotA*(1-pA) / (1-pnotK)\n",
    "\n",
    "# plot results\n",
    "ax = plt.axes()\n",
    "ax.set_color_cycle( plt.cm.Blues_r(np.linspace(0, 0.75, nump)) )\n",
    "ax.plot(pA, pnotAK.transpose());\n",
    "\n",
    "# text labels\n",
    "ax.set_xlabel('$P(H_0=1)$');\n",
    "ax.set_ylabel('$P(H0=0$ | significant test)');\n",
    "ax.legend( [\"power = %4.2f\" % x for x in 1-pnotKnotA.flatten()], \n",
    "           loc='center', frameon=False );"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The posterior probability of correctly rejecting the null hypothesis increases with power. This is what Button et al. (2013) mean when they say that low power leads to low reproducability of results. Baril and Cannon's (and Cohen's) point is that with high power you really need to take care to come up with sensible hypotheses, because you'll correctly reject your null hypothesis for tiny effects that may not mean anything, or may mean something very different from what you are actually interested in, e.g., some unforeseen, tiny, but reliable property of your particular sample of data."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###Lessons learnt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In sum, Cohen's Congress example has taught us that we can trust probability theory in making sensible inferences in uncertain situations. Especially, probability theory helps us to understand the circumstances under which our decisions about scientific findings hold and when they don't. We have seen that decisions based on null hypothesis significance testing can be reliable, but all depends on the prior probability with which a considered hypothesis is true."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}