{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Binomial and negative binomial distributions" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Today's post is prompted by [this question from Reddit](https://www.reddit.com/r/statistics/comments/4ijzkm/number_of_selections_with_replacement_before/):\n", "\n", ">How do I calculate the distribution of the number of selections (with replacement) \n", "I need to make before obtaining `k`? For example, let's say I am picking marbles from \n", "a bag with replacement. There is a 10% chance of green and 90% of black. I want `k=5` green \n", "marbles. What is the distribution number of times I need to take a marble before getting 5? \n", "\n", ">I believe this is a geometric distribution. I see how to calculate the cumulative \n", ">probability given `n` picks, but I would like to generalize it so that for any value of `k` \n", ">(number of marbles I want), I can tell you the mean, 10% and 90% probability for the \n", ">number of times I need to pick from it.\n", "\n", ">Another way of saying this is, how many times do I need to pull on a slot machine \n", ">before it pays out given that each pull is independent?\n", "\n", "Note: I've changed the notation in the question to be consistent with convention.\n" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false }, "outputs": [], "source": [ "from __future__ import print_function, division\n", "\n", "import thinkplot\n", "from thinkstats2 import Pmf, Cdf\n", "\n", "from scipy import stats\n", "from scipy import special\n", "\n", "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Solution\n", "\n", "There are two ways to solve this problem. One is to relate the desired distribution to the binomial distribution. \n", "\n", "If the probability of success on every trial is `p`, the probability of getting the `k`th success on the `n`th trial is\n", "\n", " PMF(n; k, p) = BinomialPMF(k-1; n-1, p) p\n", "\n", "That is, the probability of getting `k-1` successes in `n-1` trials, times the probability of getting the `k`th success on the `n`th trial.\n", "\n", "Here's a function that computes it:" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": false }, "outputs": [], "source": [ "def MakePmfUsingBinom(k, p, high=100):\n", " pmf = Pmf()\n", " for n in range(1, high):\n", " pmf[n] = stats.binom.pmf(k-1, n-1, p) * p\n", " return pmf" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "And here's an example using the parameters in the question." ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYcAAAEACAYAAABYq7oeAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xl8VfWd//HXJxthJ6wBEsK+CAgIAopLLIqorVhtrV1G\nbZ1f7bRM/dXpPLSd+Y3aR2da7W/s1J91qq1tta1LtaNi3XCLdWMX2QIJO2EJgbCEACEk398f9+bk\nZL9JbnJu7n0/H488vN+T77n3c68h75zzPd/vMeccIiIifklBFyAiIrFH4SAiIg0oHEREpAGFg4iI\nNKBwEBGRBhQOIiLSQEThYGYLzWyzmRWY2V1N9HnIzArNbK2ZTQ9vyzKzd8xso5mtN7Pv+vrfY2ZF\nZrYm/LUwOm9JRETaK6WlDmaWBDwMzAf2ASvN7CXn3GZfn6uAMc65cWY2B/gVMBc4C9zpnFtrZr2A\n1Wa21Lfvg865B6P8nkREpJ0iOXKYDRQ653Y55yqBZ4BF9fosAp4EcM4tB/qa2RDn3AHn3Nrw9hNA\nPjDct5+19w2IiEj0RRIOw4E9vnYRdX/BN9Znb/0+ZjYSmA4s921eHD4N9Rsz6xthzSIi0sE6ZUA6\nfErpeeCO8BEEwCPAaOfcdOAAoNNLIiIxosUxB0JHASN87azwtvp9shvrY2YphILhD865l2o6OOdK\nfP1/Dbzc2IubmRZ/EhFpA+dcm0/dR3LksBIYa2Y5ZpYG3AQsqddnCXAzgJnNBY4654rD3/stsMk5\n9wv/DmaW6WteD2xoqgDnnL6i9HXPPfcEXkO8fOmz1OcZy1/t1eKRg3OuyswWA0sJhcnjzrl8M7s9\n9G33mHPuVTO72sy2AuXArQBmNg/4KrDezD4BHPBD59zrwAPhS16rgZ3A7e1+NyIiEhWRnFYi/Mt8\nQr1tj9ZrL25kvw+B5Cae8+bIyxQRkc6kGdIJJjc3N+gS4oY+y+jS5xlbLBrnpjqSmblYr1FEJNaY\nGa6DB6RFRCTBKBxERKQBhYOIiDSgcBARkQYUDiIi0oDCQUREGlA4iIhIAwoHERFpQOEgIiINKBxE\nRKQBhYOIiDSgcBARkQYUDiIi0oDCQUREGlA4iIhIAwoHERFpQOEgIiINKBxERKQBhYOIiDSgcBAR\nkQYUDiIi0oDCQUREGlA4iIhIAwoHERFpQOEgIiINKBxERKQBhYOIiDSgcBARkQYUDiIi0oDCQURE\nGlA4iIhIAylBFyAt23vwKGs27qas/DSWZIzJHsS544eT3i016NJEJE4pHGLYzr2HePTP71Ows7jB\n91JTkrn6kil8YcFMenRPC6A6EYln5pwLuoZmmZmL9Ro7wotvr+VPLy+nuoX3ntGnB9//+gImjs7s\npMpEpCswM5xz1tb9IxpzMLOFZrbZzArM7K4m+jxkZoVmttbMpoe3ZZnZO2a20czWm9l3ff0zzGyp\nmW0xszfMrG9b30Q8cc7x9Ksr+cOSZV4wJCcnMXvqSL501SyuvWwaWUMyvP5Hjp/k3x5ewt9WFQRV\nsojEoRaPHMwsCSgA5gP7gJXATc65zb4+VwGLnXPXmNkc4BfOublmlglkOufWmlkvYDWwyDm32czu\nBw475x4IB06Gc+7uRl4/oY4cXn53Hb9/8SOvPWFUJt/5Si7DB/fztjnn+GD1Vh7/nw8pKz8NgAH/\n+5bLuei8sZ1dsojEoM44cpgNFDrndjnnKoFngEX1+iwCngRwzi0H+prZEOfcAefc2vD2E0A+MNy3\nzxPhx08A17X1TcSLLTsO8OSSZV57xqRs7v3OZ+sEA4T+p188axz3/9P1ZGeGjiIc8NAf3+HTLUWd\nWbKIxKlIwmE4sMfXLqL2F3xTffbW72NmI4HpQM1vv8HOuWIA59wBYHCkRcej8lMVPPjEW1RXVwMw\nLmcwd922kLTUpq8ZGDKgDz/6x2u900xVVdU8+Ps3OVha1ik1i0j86pSrlcKnlJ4H7nDOlTfRrclz\nR/fee6/3ODc3l9zc3GiWFxOefW0Vh46cAKBHehp33noFqanJLe7Xp1d3/s8/XMMPfv4CpcfKOXGy\ngv/726X8+x3XRbS/iMSHvLw88vLyovZ8kYw5zAXudc4tDLfvBpxz7n5fn18B7zrnng23NwOXOueK\nzSwF+CvwmnPuF7598oHccJ/M8P6TGnn9uB9z2LWvlO8/8Jw3AP29my/nopmtGzvYsuMA//rQEu/I\n4/rLZ/DVz82Jeq0i0jV0xpjDSmCsmeWYWRpwE7CkXp8lwM3hguYCR2tOGQG/BTb5g8G3z63hx7cA\nL7W+/Pjwuxc+9IJhyrhhzDtvTKufY8KoTG6+dq7XfvHttWzddTBqNYpIYmkxHJxzVcBiYCmwEXjG\nOZdvZreb2TfDfV4FdpjZVuBR4B8AzGwe8FXgM2b2iZmtMbOF4ae+H7jCzLYQuhLqp1F+b13Cpm37\nWV+wF4AkM75x/UWYtS3sP5s7lcljhwFQ7RwPP53H2bNVUatVRBKHJsEF7EeP/NW7wuiyORNY/JXL\n2vV8+0uOcef9z3Gm8iwAt1x3AddeNq3ddYpI19Ipk+CkYxTuKvaCwYAbrjiv3c85dFBfbrr6fK/9\n7GurOHL8ZLufV0QSi8IhQC++/an3+KKZ4xg6KDqTxK+5ZIo3N+J0RSV/+uvyqDyviCQOhUNASo+V\ns2LdDq/9+cunR+25U1KSue0LF3ntvOVb2LWvNGrPLyLxT+EQkLc+zveuUDpnzFByhg2I6vNPm5DF\njEnZQGgCyVM6ehCRVlA4BKCqqpo3P8r32lfOm9whr/O1z82hZjRq1cZdbN5+oENeR0Tij8IhAGvy\nd1N6LDRRvE+v7sydNqpDXmfk8IFcNHOc1/7z66s65HVEJP4oHALw3spC7/H8ORNISem4ZS5uXDjT\nO3r4dEuRJsaJSEQUDp2s/FQFqzbs9NqXnD++Q19v2OB+XOhbxvsvb67p0NcTkfigcOhkyz/dQWV4\n1vLI4QMZMbR/h7+mf/7EivU7deWSiLRI4dDJ3vPdse2SWeOa6Rk9OcP6c/6UkV77f97S0YOINE/h\n0ImOlp1kY+E+IDQj+uJWrrzaHl9YUHv08OHqrewvOdZpry0iXY/CoROtXL/Tu2nFpDFD6d+3Z6e9\n9ticwUybkAWE5j288NYnnfbaItL1KBw60XLfjOjZUzvm8tXm3OA7enhvVSHHyk51eg0i0jUoHDpJ\n+akK1oWX5gaY00FzG5pzzpihjMkeBMDZs1Us/WhTp9cgIl2DwqGTfLJpD1VVobu0jRw+kMH9e3d6\nDWbGZ3Oneu03Ptio+z2ISKMUDp1k+fraU0pzzh0ZWB0XTh9DRp8eABw5fpKP124PrBYRiV0Kh05Q\nVVXN2vw9Xnv21JGB1ZKSksyCeed47b++tz6wWkQkdikcOkHhroOcPH0GgP59e0Z9BdbWunLeZJKT\nQ//rt+4+SMHO4hb2EJFEo3DoBGs27fYeT5+Y3eZ7REdL397dudi3IJ+OHkSkPoVDJ1iTXxsO550z\nIsBKan320tqB6Y/XbvdWiRURAYVDhzty/CQ7ig4BkJSU5E1EC9qorIGcM2YoANXV1by9bHPAFYlI\nLFE4dLBPN9cORE8cNYQe3dMCrKYu/02G3vxoE9XV1QFWIyKxROHQwdaH11ICmDYxO8BKGppz7ih6\n90wH4PDRctb4rqgSkcSmcOhAzjnWFxR57XPHDw+wmoZSU5OZP3ei1176gWZMi0iIwqED7Ss5xuGj\noYHe7ulp3tIVseTyCyZ5j9ds2kVJaVmA1YhIrFA4dKD1W2rXUpo8Zqg3tyCWDB3Ul3PH167W+pYG\npkUEhUOH8p9Smhpjp5T8/DOm3/44X+stiYjCoaM459iwtXYweur42LiEtTHnT8mhX+/a9ZZWbtgV\ncEUiEjSFQwfZc+AIJ05WANC7ZzojhmYEXFHTUlLqDky/vSw/wGpEJBYoHDrI5u0HvMeTRmcGvmRG\nS+ZfUBsOa/P3cOjIiQCrEZGgKRw6SP72/d7jiaOHBlhJZIYM6MOUccOA0MB03sqCYAsSkUApHDpI\n/SOHruDyubWXtb6zbDPOuWZ6i0g8Uzh0gNJj5RwMzxdITUlmdNbAgCuKzJxpo+iRHlreo/jwcTb4\nZneLSGJROHSAfN9Rw7icwaSkJAdYTeTSUlO4ZFbtUt5ajE8kcSkcOsBm33jDpC4w3uDnnzG97NPt\nlJ+qCLAaEQlKROFgZgvNbLOZFZjZXU30ecjMCs1srZnN8G1/3MyKzWxdvf73mFmRma0Jfy1s31uJ\nHZt31N5ZbWIXGW+oMSprICOHh06DVZ6t4v1VWwOuSESC0GI4mFkS8DBwJTAZ+LKZTazX5ypgjHNu\nHHA78N++b/8uvG9jHnTOnRf+er0tbyDWnK6oZGf4/g0GTBg1JNiC2mD+3Ane47eX69SSSCKK5Mhh\nNlDonNvlnKsEngEW1euzCHgSwDm3HOhrZkPC7Q+AI008d2xf/N8GBTuLqQ5f5ZM9tD89u3cLuKLW\nu2TWeG+cZPueEnbuPRRwRSLS2SIJh+GAf6H/ovC25vrsbaRPYxaHT0P9xsz6RtA/5uXXuYS1a403\n1OjVoxtzzh3ltTUwLZJ4ghyQfgQY7ZybDhwAHgywlqjxz2+YOLrrnVKq4V9O472VhZypPBtgNSLS\n2VIi6LMXGOFrZ4W31e+T3UKfOpxzJb7mr4GXm+p77733eo9zc3PJzc1t7qkDU1VVzZad/sHornnk\nAKEbEw3K6E3JkTLKT1WwYv1OLjpvbNBliUgT8vLyyMvLi9rzRRIOK4GxZpYD7AduAr5cr88S4DvA\ns2Y2FzjqnCv2fd+oN75gZpnOuZo/s68HNjRVgD8cYtmufYepOFMJwIB+PRmU0SvgitrOzPjM3Ak8\n+9oqIDRjWuEgErvq/+F83333tev5Wjyt5JyrAhYDS4GNwDPOuXwzu93Mvhnu8yqww8y2Ao8C367Z\n38yeAj4CxpvZbjP7evhbD5jZOjNbC1wKfK9d7yQG+I8aJoyK/cX2WnLZ7Aleoq/bUqS7xIkkkEiO\nHAhfZjqh3rZH67UXN7HvV5rYfnOENXYZW3fXnikbn9N1xxtqDOrfm6njs1hXUOQtxvfFK2cGXZaI\ndALNkI6irbsOeo/Hjoi9+0W3hX9g+t3lW7QYn0iCUDhEyanTZ9hbHJrOYYRmGseD2eeOrLMY38at\nWoxPJBEoHKJk254Sav6mzh7an/RuqYHWEy1pqSlcPLN2Mb53lm8JsBoR6SwKhyjxjzeMHTE4wEqi\nz39q6aNPtnHy1JkAqxGRzqBwiBJ/OIzLia9wGJ09kBFD+wOhxfg+/ESL8YnEO4VDlMTjYHQNM+Mz\nc2qPHrSchkj8UzhEwbGyU5QcCc0BSElJ9v7KjieXzBpHcnLox6Vw10H2HGhqLUURiQcKhyjYurv2\nqGHU8AFd5s5vrdG3d3fOn5zjtd/VUt4icU3hEAXxPN7gd5lvYDpvZQFnz1YFWI2IdCSFQxRsi+Mr\nlfxmTMwmo08PIHQqbU3+nhb2EJGuSuHQTs45Cn2nlcbG8ZFDcnISueeP99o6tSQSvxQO7VRy5ATH\nT5wCoHt6GsMGxcU9i5rkP7W0auNujpadDLAaEekoCod28g9Gj8ke2OVXYm3J8MH9mDAqE4Dq6mre\nW1kYcEUi0hEUDu3kH28YF8fjDX7z59Yu0PvOss1ajE8kDikc2qnOkUOChMOF08fQLS20dlRR8ZE6\nn4GIxAeFQzs459hRdNhrj86Oj5VYW9I9PY0Lpo/22poxLRJ/FA7tcOjICcpPVQDQIz2Nwf17B1xR\n5/EvxvfBmm3e7VFFJD4oHNphx97ao4ZRWfE/GO03aXQmmQP7AKF7WSz7dEfAFYlINCkc2mFH0SHv\n8ajhiXFKqYaZcZlvMb53NOdBJK4oHNqhTjhkDQiwkmDknj+emmOlDYX7KD58PNB6RCR6FA7tsLPe\naaVEMzCjF9MnZXtt3SVOJH4oHNqorPx0nWW6hw/uF3BFwfCfWspbsUVzHkTihMKhjfxHDSOG9o/L\nZbojMXvKSHr16AaErt5aV7A34IpEJBoUDm20Y2/teMPIYYk33lAjNTWZS2aN89pvfpQfYDUiEi0K\nhzbyD0YnyuS3plx+wSTv8Yr1OzhWdirAakQkGhQObVRnjkOCXcZaX86wAd5Njqqqqnl3hQamRbo6\nhUMbnKk8y97wPZQNGDk8cU8r1Vhw4Tne4zc/2qSBaZEuTuHQBrv3lVId/uU3dFBf0rulBlxR8Oad\nN4Ye6WkAHDh0nA2F+wKuSETaQ+HQBv7B6JwEP6VUo1taKpf67hL3xoebAqxGRNpL4dAGdVZiTcDJ\nb0254kINTIvEC4VDG9S5jFXjDZ6cYQMYP3IIEBqY1npLIl2XwqGVqqur60yAS/TLWOvzD0y/9XG+\nBqZFuiiFQyvtKznGmcqzAGT06UG/3j0Crii2XDhjdJ2B6fWaMS3SJSkcWmmnb7xBp5Qa6paWSu7s\n2oHppZoxLdIlKRxayT/eMDprUICVxK7LL6g9tbR83Q6Olp0MsBoRaQuFQyv5l83I0ZFDo3KG9WfC\nqEwgNEbzzjLNmBbpaiIKBzNbaGabzazAzO5qos9DZlZoZmvNbIZv++NmVmxm6+r1zzCzpWa2xcze\nMLO+7XsrHc85V2fZDF3G2rQFvsta3/hwI1VV1QFWIyKt1WI4mFkS8DBwJTAZ+LKZTazX5ypgjHNu\nHHA78N++b/8uvG99dwNvOecmAO8AP2jTO+hEpcfKOX4idO1+erdU7x7K0tCFM8bQu2c6EFrKe+WG\nncEWJCKtEsmRw2yg0Dm3yzlXCTwDLKrXZxHwJIBzbjnQ18yGhNsfAEcaed5FwBPhx08A17W+/M7l\nP2oYOXwAZtZM78SWlppS57LWV/+2IcBqRKS1IgmH4cAeX7sovK25Pnsb6VPfYOdcMYBz7gAwOIJa\nAlVnmW6dUmrRgnnnkBQO0I1b97Fr3+EW9hCRWBFLA9IxP1tqZ5FmRrfGwIxezJk22mvr6EGk60iJ\noM9eYISvnRXeVr9Pdgt96is2syHOuWIzywQONtXx3nvv9R7n5uaSm5vbctUdoO5gtC5jjcTVl0zh\n47XbAHhvZQFf+9wcbyxCRKInLy+PvLy8qD2ftbS8gZklA1uA+cB+YAXwZedcvq/P1cB3nHPXmNlc\n4L+cc3N93x8JvOycm+rbdj9Q6py7P3wFVIZz7u5GXt/FwhIM5acquPnu3wGQlJTEUw/cRmpqYt43\nujWcc/zTA897p5RuXnQBiz4zLeCqROKfmeGca/PAaIunlZxzVcBiYCmwEXjGOZdvZreb2TfDfV4F\ndpjZVuBR4Nu+Ap8CPgLGm9luM/t6+Fv3A1eYWU3w/LStb6Iz+NdTys7MUDBEyMy45tIpXvu1v22g\nulqXtYrEukhOK+Gcex2YUG/bo/Xai5vY9ytNbC8FLo+szOD5B6NHaTC6VS6eOY4nX1rGiZMVlBwp\nY9XG3cyeOjLoskSkGbE0IB3T6t4zWoPRrZGWmsIVF9ROinvlvXXN9BaRWKBwiJCOHNrnyosmU3Py\nc0PhPnbvLw20HhFpnsIhAmfPVlFUXDuPT5extt6g/r2Zc+4or/3Ke+sDrEZEWqJwiMCeA0e8tYGG\nDOhDz+7dAq6oa7r6Uu9iNfJWFmi1VpEYpnCIwPaiEu+xxhva7pwxQxmTHZofcvZsFa+9vzHgikSk\nKQqHCPgvYx2p8YY2MzMWzZ/utV9/fwOnKyoDrEhEmqJwiMB2DUZHzdxzRzG4f28ATpys4J3lmwOu\nSEQao3BogXOuzpGDTiu1T3JyEp/NPddrv/zuOt3rQSQGKRxacODQce/UR59e3enft2fAFXV98+dO\npFeP0KD+wdIylq3bEXBFIlKfwqEFdU4p6R4OUZHeLZWFF0322i+9vZZYWD9LRGopHFqwy39KSeMN\nUXPVJVNISQmtT7VtTwnrClpaxFdEOpPCoQV1L2NVOERLv949uHxu7d1mn39jdYDViEh9CocW1L2M\nVYPR0XTd/OkkJYV+BDdt28/GrfsCrkhEaigcmnG07CRHjodm8aalpjBsUN+AK4ovg/r3Jvf88V77\n+TfWBFiNiPgpHJqxfU/d24LW/JUr0XP9FTO8BfnWFRRRsLM40HpEJES/7Zrhv1JptAajO8TQQX25\naOY4r/2XpTp6EIkFCodm7NhTOxg9Olvh0FFuWHCe93jVxl11lkcXkWAoHJpR98hhUICVxLfszAzm\n+pbzfvqVlQFWIyKgcGjSiZMVHCwtA0JLPmRnZgRcUXy78apZ3tjD6k27NPYgEjCFQxP8pzZGDO3v\nTdiSjpEzbAAXnjfWaz/1yooAqxERhUMTNBjd+b501SySwsuTrC/Yy3rNmhYJjMKhCf6Z0Rpv6BzD\nB/cjd/YEr/3UKyu05pJIQBQOTdjhm+OgK5U6zxcXziQ5OfRjWbCzmNWbdgdckUhiUjg04nRFJfsO\nHgXAgJxh/YMtKIEM7t+bBRee47X/+PJyqqt1vweRzqZwaMTOvYepOZmRlZlBt7TUQOtJNDcsOM/7\nzPfsL+Wd5VsCrkgk8SgcGlFnJVYNRne6jD49uG7+NK/99Csrda9pkU6mcGjEjqLalVg1GB2May+b\nRkafHkBoAcQX3l4bcEUiiUXh0Ig6l7FqMDoQ6d1S+co1s732S2+v5fDREwFWJJJYFA71VFZWsXt/\nqdceOVz3cAhK7uzx5AwLff6VZ6t4SstqiHQahUM9u/eXelfHZA7sQ8/u3QKuKHElJSVx63UXeO33\nVmzRshoinUThUE/dwWiNNwTt3AlZzJqcA4ADHnvufV3aKtIJFA71FO466D3Wshmx4Rs3zCM1vLbV\njqJDLP0wP+CKROKfwqEefziMHzk4wEqkxpABfbj+ihle+09/Xc6xslMBViQS/xQOPqcrKtkTHow2\nYEy2TivFiuvmT2fIgD4AnDx9hj++vDzgikTim8LBZ9uektqZ0UP70z09LdB6pFZaagq33TDPa7+z\nfDP52/YHWJFIfFM4+GzdXTsYPW6ETinFmpmTc5g9daTXfuTpPM5Ung2uIJE4FlE4mNlCM9tsZgVm\ndlcTfR4ys0IzW2tm01va18zuMbMiM1sT/lrY/rfTPv7LJMflKBxi0W03XER6t9C6S/tKjvHn11YF\nXJFIfGoxHMwsCXgYuBKYDHzZzCbW63MVMMY5Nw64HfhVhPs+6Jw7L/z1ejTeUHts3V07GK1wiE0D\nM3pxy6LauQ8vvr2Wrb6LCEQkOiI5cpgNFDrndjnnKoFngEX1+iwCngRwzi0H+prZkAj2NWLEkeMn\nOXQktDxDakqy7hkdw664cBJTxg0DQnMffvl0HmfPVgVblEiciSQchgN7fO2i8LZI+rS07+Lwaajf\nmFnfiKvuAHXmN2QP0j2jY5iZ8a0vXUpaagoQmtX+3NI1AVclEl9SOuh5IzkieAT4kXPOmdmPgQeB\n2xrreO+993qPc3Nzyc3NjUKJdflPTYzXKaWYN3RQX75yzWx+/+JHAPzljdXMmJjNxNGZAVcmEoy8\nvDzy8vKi9nyRhMNeYISvnRXeVr9PdiN90pra1zlX4tv+a+Dlpgrwh0NH8R85jFU4dAnXXDqFFet3\nsGnbfhzwX0++zX/e9QWthyUJqf4fzvfdd1+7ni+S00orgbFmlmNmacBNwJJ6fZYANwOY2VzgqHOu\nuLl9zcz/J971wIZ2vZN2cM5pMLoLSkpK4o6/m0+P8HyUkiNlPPbc+wFXJRIfWgwH51wVsBhYCmwE\nnnHO5ZvZ7Wb2zXCfV4EdZrYVeBT4dnP7hp/6ATNbZ2ZrgUuB70X3rUVuX8kxTp4+A0DvnukM7t87\nqFKklQZm9OJbN13qtT9YvZX3VhYEWJFIfIhozCF8memEetserddeHOm+4e03R15mxyr0zW8YnzME\ns5i5iEoiMG/GGD7J38274XtNP/rn9xkzYhBZQ3TFmUhbaYY09ccbtJ5SV3Tb9fPIHBhae6niTCU/\ne3wpp8JHgyLSegoHoMAfDlo2o0vqnp7G97++wFvau6j4CL98+j2ccy3sKSKNSfhwOHX6DDv2hC6c\nMmDCqCHBFiRtNiprIN/60iVe++O123g5b12AFYl0XQkfDgW7DnorsWYP7a/LILu43NkTuHLeZK/9\nh5eW8emWogArEumaEj4cNvmWfT5nzNAAK5Fo+cb1F3qXI1c7x89+u5Td4ft0iEhkEj4cNm+vDYdJ\noxUO8SAlJZl//sYC+vftCYROHf7Ho69xtOxkwJWJdB0JHQ5nz1axZUftZaxaeiF+DOjXix9+8yq6\npYWW9y45UsZPHnudijOVAVcm0jUkdDhsLzpEZXg1z0EZvRmY0SvgiiSaRmUN5M5bL/cW+tq6+yAP\n/v4treAqEoGEDgf/eMOkMTpqiEezJudw2xcu8tqrNu7ioT+9S3V1dYBVicS+hA6HDYW16wfW3B9A\n4s9VF0/h8/O9mxPy4Zqt/OrZv2kOhEgzEjYczp6tYuPW2iOHqeOzAqxGOtpXPzeHhRfVXuL69rLN\n/PZ/PlRAiDQhYcOhYNdB7+b0Qwb00WJ7cc7M+PsvXETu7Nplvl792wYee+59BYRIIxI2HNYX6JRS\nojEzvn3TpcydNtrbtvTDTTz0x3eoqtIYhIifwgE4V6eUEkZychJ33nI5F88c523726pCHvz9m1RW\n6iomkRoJGQ6nKyop2FU7v2HKeB05JJLk5CTu+LvPcMWFk7xty9bt4N5HXub4iVMBViYSOxIyHNYV\n7PVOI4wY2p9+vXsEXJF0NjPj9hsv4XO553rbNm8/wA9+/gJ7Dx4NsDKR2JCQ4bB64y7v8cxzRjTT\nU+KZmXHLdRdw86ILvIlyBw4d5wcPvlDntKNIIkq4cHDOsWbTbq89c3JOgNVI0MyMRZ+Zxve/UXsv\niPJTFdz3y5f5y5trdCWTJKyEC4edew9TeqwcgJ7duzF+pO7fIDB32mh+/N1F9O3dHQAHPPXXFfzk\nsdc5cbIi2OJEApBw4bDad9QwfVI2yckJ9xFIE8bmDOZn37+BCaNql1JZvWkX33/geTZu3RdgZSKd\nL+F+M66fs9akAAAKSUlEQVRcv9N7PGuyxhukrgH9evGjxZ+rM1BdcqSMe/7fEp548WNv4qRIvEuo\ncDhYWsbW3aH7RSclJTFjksJBGkpJSebWz1/I97++wLszoAOWvPsp//yzv7B5+4FgCxTpBAkVDh+u\n2eo9nj4xi9490wOsRmLdBdNH8/O7v8j0idnetqLiI/zLL17kkafzKCs/HWB1Ih0rscLhk23e4wun\njwmwEukqBvTrxb9+62q++cWLSUtN8ba/vWwz//jvz7D0w01aekPiksX6pXpm5qJR4/6SYyz+8dNA\naIbs7/79Fu+UgUgkDpaW8du/fMjKDTvrbM8aksHXrp3DrMk5mFnjO4t0MjPDOdfmH8iECYdnXlvJ\nc6+vBuD8KSO5+38tbPdzSmJasX4nj//lAw4dOVFn+4RRmdy4cCbTJmQpJCRwCocIVFVV8637/uTN\nb7jz1iuYN0OnlaTtKs5U8nLeel546xNOV9S9L/W4nMHcsOA8HUlIoBQOEVixfif3/+Z1APr06s6v\n7/saKeHZsCLtcazsFM8vXc0bjYw9DBvUl6sumcJlsyfQPT0toAolUSkcIvDjX73CJ/l7ALj+8hl8\n9XNzolGaiKektIwX3lrL28s3c/Zs3aW/u6encdF5Y7hs9gTGjxyiownpFAqHFuw5cITv/eRZHGDA\nL//tKwwZ0Cdq9Yn4lR4r5+V31/HWx/mcPH2mwfeHDerLpbMncMmscbr7oHQohUMLHnziLW9+w8xz\ncvjh7VdFqzSRJp2uqCRvRQGvvLeOfSXHGu0zLmcw508dyeypo8ga0k9HFBJVCodm7N5fyp0//TM1\ne//0zs8zLkcL7Unncc6Rv/0A7y7fwkdrtzUYvK6RObAPs6eOYtrELCaOyiS9W2onVyrxRuHQBOcc\n9//mDe+adB01SNAqzlSyYt1O3l2xhXVbimjqpzo5OYmxIwYzddwwpowbzricwQoLaTWFQxM+/GQb\nD/7+Ta99/53XMzZncDRLE2mzsvLTrN64i5Xrd/LJ5iIqzjR+RAGhsbKszAzG5gxmbPZgxo4YxIhh\n/evM2BapT+HQiGNlp7jjJ896a99cceEkvvWlSzuiPJF2O1N5lvUFe/kkfw8bCvey58CRFvcxYOig\nvowY2p+sof0ZMbQ/wwb1ZfCA3pr5L0AnhYOZLQT+i9BaTI875+5vpM9DwFVAOXCrc25tc/uaWQbw\nLJAD7ARudM41GLlrbTicOn2Gex5+mW17SgAY0K8nP7/7Rv2DkS7jWNkpNmzdx4bCveRvP0DR/tIm\nT0E1plePbmQO7MuQgX3IHNCHIQN7k9GnJwP69aR/35706tFNg98JoMPDwcySgAJgPrAPWAnc5Jzb\n7OtzFbDYOXeNmc0BfuGcm9vcvmZ2P3DYOfeAmd0FZDjn7m7k9SMOh9Jj5Tz4+7fI374/tC/wL9+6\nhhmTspvfMYHk5eWRm5sbdBlxobM+y9MVlWzbU8LW3SVs3X2Q7XtKKD50vFWB4ZecnET/Pj3p368n\n/Xp3p1ePbvTpmU6vnun07tmN3j27h9vd6N0jne7pqaSmJHd4oOhnM7raGw6RnLScDRQ653aFX/AZ\nYBGw2ddnEfAkgHNuuZn1NbMhwKhm9l0E1JzreQLIAxqEQyROnKwgb8UWnl+6ps4yyt+88RIFQz36\nBxg9nfVZpndLZfLYYUweO8zbVnGmkr3FR9m9v5Td+0spOnCU4sPHOXD4eINJePVVVVVTcqSMkiNl\nEdeQZEZ6t1S6p6fSvVtanf+md0slPS2VbmkppKYkk5qaTFpqCqkpSeH/+rclk5aSTEpKEinJySQl\nGUlJSSQnGa+89gZTps0iOTmJlOQkkpKM5KTQ4+TkJB3tdLJIwmE4sMfXLiIUGC31Gd7CvkOcc8UA\nzrkDZtbkaPF/PPoa1a6aqipHtaumutpRVe2orq6mrPx0g7+iDPjatXNZMO+cCN6eSNfTLS2V0dmD\nGJ09qM525xyHj5ZTfPg4xYeOc+DQcUqOlHHk+ElKj5ZTevwkpxqZnNeSauc4efpMeGJfeZTeRV2b\nlq1n++k/Nvl9A5KSk0hOCoVFkhlmob+QLfw4yffYzDDC25PC26nbv/Zx7f7g699EIPk3h5614fa6\n/VvZh4ad6rxmC883a/LIxl+kFTrqcoe2RHyTR8mrN+2K+Eky+vTgjr+bz9Txw9tQgkjXZmYMzOjF\nwIxedY40/E5XVFJ6rJzSY+UcO3GaE+WnKTtZQdmJ0xwvP8WJkxWUlZ+mrPw0J05WcKqiMibuWeEI\nHfVUVVVD0xd3CZA5sG/7n8Q51+wXMBd43de+G7irXp9fAV/ytTcDQ5rbF8gndPQAkAnkN/H6Tl/6\n0pe+9NX6r5Z+vzf3FcmRw0pgrJnlAPuBm4Av1+uzBPgO8KyZzQWOOueKzexQM/suAW4F7gduAV5q\n7MXbM6AiIiJt02I4OOeqzGwxsJTay1Hzzez20LfdY865V83sajPbSuiE5Neb2zf81PcDfzazbwC7\ngBuj/u5ERKRNYn4SnIiIdL6koAtoipktNLPNZlYQngchrWRmO83sUzP7xMxWhLdlmNlSM9tiZm+Y\nWRRGruKTmT1uZsVmts63rcnPz8x+YGaFZpZvZguCqTp2NfF53mNmRWa2Jvy10Pc9fZ5NMLMsM3vH\nzDaa2Xoz+254e9R+PmMyHMKT5x4GrgQmA182s4nBVtUlVQO5zrkZzrmaS4jvBt5yzk0A3gF+EFh1\nse93hH4G/Rr9/MzsHEKnRicRWingEdOF+fU19nkCPOicOy/89TqAmU1Cn2dzzgJ3OucmAxcA3wn/\njozaz2dMhgO+iXfOuUqgZvKctI7R8P/xIkKTDgn/97pOragLcc59ANRf6Kipz+9a4Bnn3Fnn3E6g\nkIbzgRJaE58nNH7p+yL0eTbJOXegZoki59wJQld/ZhHFn89YDYemJtVJ6zjgTTNbaWZ/H95WZ/Ih\noKVqW2dwE59f/Z/ZvehnNlKLzWytmf3GdxpEn2eEzGwkMB1YRtP/vlv9ecZqOEh0zHPOnQdcTeiw\n82JCgeGnKxLaR59f+zwCjHbOTQcOAP8ZcD1dipn1Ap4H7ggfQUTt33eshsNeYISvnRXeJq3gnNsf\n/m8J8CKhw8ji8LpXmFkmcDC4Crukpj6/vYB/IS/9zEbAOVfiW1nz19Se6tDn2QIzSyEUDH9wztXM\nE4vaz2eshoM38c7M0ghNnlsScE1dipn1CP9VgZn1BBYA66mdfAjNTD4Uj1H3nHhTn98S4CYzSzOz\nUcBYYEVnFdmF1Pk8w7/AalwPbAg/1ufZst8Cm5xzv/Bti9rPZ0zeSqqFyXMSmSHAC2bmCP1//pNz\nbqmZrUKTDyNiZk8BucAAM9sN3AP8FHiu/ufnnNtkZn8GNhFa+efbbbqFYRxr4vO8zMymE7qybidw\nO+jzbImZzQO+Cqw3s08InT76IU1MLm7L56lJcCIi0kCsnlYSEZEAKRxERKQBhYOIiDSgcBARkQYU\nDiIi0oDCQUREGlA4iIhIAwoHERFp4P8DJk4Ndblud50AAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "pmf = MakePmfUsingBinom(5, 0.1, 200)\n", "thinkplot.Pdf(pmf)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can solve the same problem using the negative binomial distribution, but it requires some translation from the parameters of the problem to the conventional parameters of the binomial distribution.\n", "\n", "The negative binomial PMF is the probability of getting `r` non-terminal events before the `k`th terminal event. (I am using \"terminal event\" instead of \"success\" and \"non-terminal\" event instead of \"failure\" because in the context of the negative binomial distribution, the use of \"success\" and \"failure\" is often reversed.)\n", "\n", "If `n` is the total number of events, `n = k + r`, so\n", "\n", " r = n - k\n", "\n", "If the probability of a terminal event on every trial is `p`, the probability of getting the `k`th terminal event on the `n`th trial is\n", "\n", " PMF(n; k, p) = NegativeBinomialPMF(n-k; k, p) p\n", "\n", "That is, the probability of `n-k` non-terminal events on the way to getting the `k`th terminal event.\n", "\n", "Here's a function that computes it:" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false }, "outputs": [], "source": [ "def MakePmfUsingNbinom(k, p, high=100):\n", " pmf = Pmf()\n", " for n in range(1, high):\n", " r = n-k\n", " pmf[n] = stats.nbinom.pmf(r, k, p)\n", " return pmf" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Here's the same example:" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYcAAAEACAYAAABYq7oeAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xl8VfWd//HXJxthJ6wBEsK+CAgIAopLLIqorVhtrV1G\nbZ1f7bRM/dXpPLSd+Y3aR2da7W/s1J91qq1tta1LtaNi3XCLdWMX2QIJO2EJgbCEACEk398f9+bk\nZL9JbnJu7n0/H488vN+T77n3c68h75zzPd/vMeccIiIifklBFyAiIrFH4SAiIg0oHEREpAGFg4iI\nNKBwEBGRBhQOIiLSQEThYGYLzWyzmRWY2V1N9HnIzArNbK2ZTQ9vyzKzd8xso5mtN7Pv+vrfY2ZF\nZrYm/LUwOm9JRETaK6WlDmaWBDwMzAf2ASvN7CXn3GZfn6uAMc65cWY2B/gVMBc4C9zpnFtrZr2A\n1Wa21Lfvg865B6P8nkREpJ0iOXKYDRQ653Y55yqBZ4BF9fosAp4EcM4tB/qa2RDn3AHn3Nrw9hNA\nPjDct5+19w2IiEj0RRIOw4E9vnYRdX/BN9Znb/0+ZjYSmA4s921eHD4N9Rsz6xthzSIi0sE6ZUA6\nfErpeeCO8BEEwCPAaOfcdOAAoNNLIiIxosUxB0JHASN87azwtvp9shvrY2YphILhD865l2o6OOdK\nfP1/Dbzc2IubmRZ/EhFpA+dcm0/dR3LksBIYa2Y5ZpYG3AQsqddnCXAzgJnNBY4654rD3/stsMk5\n9wv/DmaW6WteD2xoqgDnnL6i9HXPPfcEXkO8fOmz1OcZy1/t1eKRg3OuyswWA0sJhcnjzrl8M7s9\n9G33mHPuVTO72sy2AuXArQBmNg/4KrDezD4BHPBD59zrwAPhS16rgZ3A7e1+NyIiEhWRnFYi/Mt8\nQr1tj9ZrL25kvw+B5Cae8+bIyxQRkc6kGdIJJjc3N+gS4oY+y+jS5xlbLBrnpjqSmblYr1FEJNaY\nGa6DB6RFRCTBKBxERKQBhYOIiDSgcBARkQYUDiIi0oDCQUREGlA4iIhIAwoHERFpQOEgIiINKBxE\nRKQBhYOIiDSgcBARkQYUDiIi0oDCQUREGlA4iIhIAwoHERFpQOEgIiINKBxERKQBhYOIiDSgcBAR\nkQYUDiIi0oDCQUREGlA4iIhIAwoHERFpQOEgIiINKBxERKQBhYOIiDSgcBARkQYUDiIi0oDCQURE\nGlA4iIhIAylBFyAt23vwKGs27qas/DSWZIzJHsS544eT3i016NJEJE4pHGLYzr2HePTP71Ows7jB\n91JTkrn6kil8YcFMenRPC6A6EYln5pwLuoZmmZmL9Ro7wotvr+VPLy+nuoX3ntGnB9//+gImjs7s\npMpEpCswM5xz1tb9IxpzMLOFZrbZzArM7K4m+jxkZoVmttbMpoe3ZZnZO2a20czWm9l3ff0zzGyp\nmW0xszfMrG9b30Q8cc7x9Ksr+cOSZV4wJCcnMXvqSL501SyuvWwaWUMyvP5Hjp/k3x5ewt9WFQRV\nsojEoRaPHMwsCSgA5gP7gJXATc65zb4+VwGLnXPXmNkc4BfOublmlglkOufWmlkvYDWwyDm32czu\nBw475x4IB06Gc+7uRl4/oY4cXn53Hb9/8SOvPWFUJt/5Si7DB/fztjnn+GD1Vh7/nw8pKz8NgAH/\n+5bLuei8sZ1dsojEoM44cpgNFDrndjnnKoFngEX1+iwCngRwzi0H+prZEOfcAefc2vD2E0A+MNy3\nzxPhx08A17X1TcSLLTsO8OSSZV57xqRs7v3OZ+sEA4T+p188axz3/9P1ZGeGjiIc8NAf3+HTLUWd\nWbKIxKlIwmE4sMfXLqL2F3xTffbW72NmI4HpQM1vv8HOuWIA59wBYHCkRcej8lMVPPjEW1RXVwMw\nLmcwd922kLTUpq8ZGDKgDz/6x2u900xVVdU8+Ps3OVha1ik1i0j86pSrlcKnlJ4H7nDOlTfRrclz\nR/fee6/3ODc3l9zc3GiWFxOefW0Vh46cAKBHehp33noFqanJLe7Xp1d3/s8/XMMPfv4CpcfKOXGy\ngv/726X8+x3XRbS/iMSHvLw88vLyovZ8kYw5zAXudc4tDLfvBpxz7n5fn18B7zrnng23NwOXOueK\nzSwF+CvwmnPuF7598oHccJ/M8P6TGnn9uB9z2LWvlO8/8Jw3AP29my/nopmtGzvYsuMA//rQEu/I\n4/rLZ/DVz82Jeq0i0jV0xpjDSmCsmeWYWRpwE7CkXp8lwM3hguYCR2tOGQG/BTb5g8G3z63hx7cA\nL7W+/Pjwuxc+9IJhyrhhzDtvTKufY8KoTG6+dq7XfvHttWzddTBqNYpIYmkxHJxzVcBiYCmwEXjG\nOZdvZreb2TfDfV4FdpjZVuBR4B8AzGwe8FXgM2b2iZmtMbOF4ae+H7jCzLYQuhLqp1F+b13Cpm37\nWV+wF4AkM75x/UWYtS3sP5s7lcljhwFQ7RwPP53H2bNVUatVRBKHJsEF7EeP/NW7wuiyORNY/JXL\n2vV8+0uOcef9z3Gm8iwAt1x3AddeNq3ddYpI19Ipk+CkYxTuKvaCwYAbrjiv3c85dFBfbrr6fK/9\n7GurOHL8ZLufV0QSi8IhQC++/an3+KKZ4xg6KDqTxK+5ZIo3N+J0RSV/+uvyqDyviCQOhUNASo+V\ns2LdDq/9+cunR+25U1KSue0LF3ntvOVb2LWvNGrPLyLxT+EQkLc+zveuUDpnzFByhg2I6vNPm5DF\njEnZQGgCyVM6ehCRVlA4BKCqqpo3P8r32lfOm9whr/O1z82hZjRq1cZdbN5+oENeR0Tij8IhAGvy\nd1N6LDRRvE+v7sydNqpDXmfk8IFcNHOc1/7z66s65HVEJP4oHALw3spC7/H8ORNISem4ZS5uXDjT\nO3r4dEuRJsaJSEQUDp2s/FQFqzbs9NqXnD++Q19v2OB+XOhbxvsvb67p0NcTkfigcOhkyz/dQWV4\n1vLI4QMZMbR/h7+mf/7EivU7deWSiLRI4dDJ3vPdse2SWeOa6Rk9OcP6c/6UkV77f97S0YOINE/h\n0ImOlp1kY+E+IDQj+uJWrrzaHl9YUHv08OHqrewvOdZpry0iXY/CoROtXL/Tu2nFpDFD6d+3Z6e9\n9ticwUybkAWE5j288NYnnfbaItL1KBw60XLfjOjZUzvm8tXm3OA7enhvVSHHyk51eg0i0jUoHDpJ\n+akK1oWX5gaY00FzG5pzzpihjMkeBMDZs1Us/WhTp9cgIl2DwqGTfLJpD1VVobu0jRw+kMH9e3d6\nDWbGZ3Oneu03Ptio+z2ISKMUDp1k+fraU0pzzh0ZWB0XTh9DRp8eABw5fpKP124PrBYRiV0Kh05Q\nVVXN2vw9Xnv21JGB1ZKSksyCeed47b++tz6wWkQkdikcOkHhroOcPH0GgP59e0Z9BdbWunLeZJKT\nQ//rt+4+SMHO4hb2EJFEo3DoBGs27fYeT5+Y3eZ7REdL397dudi3IJ+OHkSkPoVDJ1iTXxsO550z\nIsBKan320tqB6Y/XbvdWiRURAYVDhzty/CQ7ig4BkJSU5E1EC9qorIGcM2YoANXV1by9bHPAFYlI\nLFE4dLBPN9cORE8cNYQe3dMCrKYu/02G3vxoE9XV1QFWIyKxROHQwdaH11ICmDYxO8BKGppz7ih6\n90wH4PDRctb4rqgSkcSmcOhAzjnWFxR57XPHDw+wmoZSU5OZP3ei1176gWZMi0iIwqED7Ss5xuGj\noYHe7ulp3tIVseTyCyZ5j9ds2kVJaVmA1YhIrFA4dKD1W2rXUpo8Zqg3tyCWDB3Ul3PH167W+pYG\npkUEhUOH8p9Smhpjp5T8/DOm3/44X+stiYjCoaM459iwtXYweur42LiEtTHnT8mhX+/a9ZZWbtgV\ncEUiEjSFQwfZc+AIJ05WANC7ZzojhmYEXFHTUlLqDky/vSw/wGpEJBYoHDrI5u0HvMeTRmcGvmRG\nS+ZfUBsOa/P3cOjIiQCrEZGgKRw6SP72/d7jiaOHBlhJZIYM6MOUccOA0MB03sqCYAsSkUApHDpI\n/SOHruDyubWXtb6zbDPOuWZ6i0g8Uzh0gNJj5RwMzxdITUlmdNbAgCuKzJxpo+iRHlreo/jwcTb4\nZneLSGJROHSAfN9Rw7icwaSkJAdYTeTSUlO4ZFbtUt5ajE8kcSkcOsBm33jDpC4w3uDnnzG97NPt\nlJ+qCLAaEQlKROFgZgvNbLOZFZjZXU30ecjMCs1srZnN8G1/3MyKzWxdvf73mFmRma0Jfy1s31uJ\nHZt31N5ZbWIXGW+oMSprICOHh06DVZ6t4v1VWwOuSESC0GI4mFkS8DBwJTAZ+LKZTazX5ypgjHNu\nHHA78N++b/8uvG9jHnTOnRf+er0tbyDWnK6oZGf4/g0GTBg1JNiC2mD+3Ane47eX69SSSCKK5Mhh\nNlDonNvlnKsEngEW1euzCHgSwDm3HOhrZkPC7Q+AI008d2xf/N8GBTuLqQ5f5ZM9tD89u3cLuKLW\nu2TWeG+cZPueEnbuPRRwRSLS2SIJh+GAf6H/ovC25vrsbaRPYxaHT0P9xsz6RtA/5uXXuYS1a403\n1OjVoxtzzh3ltTUwLZJ4ghyQfgQY7ZybDhwAHgywlqjxz2+YOLrrnVKq4V9O472VhZypPBtgNSLS\n2VIi6LMXGOFrZ4W31e+T3UKfOpxzJb7mr4GXm+p77733eo9zc3PJzc1t7qkDU1VVzZad/sHornnk\nAKEbEw3K6E3JkTLKT1WwYv1OLjpvbNBliUgT8vLyyMvLi9rzRRIOK4GxZpYD7AduAr5cr88S4DvA\ns2Y2FzjqnCv2fd+oN75gZpnOuZo/s68HNjRVgD8cYtmufYepOFMJwIB+PRmU0SvgitrOzPjM3Ak8\n+9oqIDRjWuEgErvq/+F83333tev5Wjyt5JyrAhYDS4GNwDPOuXwzu93Mvhnu8yqww8y2Ao8C367Z\n38yeAj4CxpvZbjP7evhbD5jZOjNbC1wKfK9d7yQG+I8aJoyK/cX2WnLZ7Aleoq/bUqS7xIkkkEiO\nHAhfZjqh3rZH67UXN7HvV5rYfnOENXYZW3fXnikbn9N1xxtqDOrfm6njs1hXUOQtxvfFK2cGXZaI\ndALNkI6irbsOeo/Hjoi9+0W3hX9g+t3lW7QYn0iCUDhEyanTZ9hbHJrOYYRmGseD2eeOrLMY38at\nWoxPJBEoHKJk254Sav6mzh7an/RuqYHWEy1pqSlcPLN2Mb53lm8JsBoR6SwKhyjxjzeMHTE4wEqi\nz39q6aNPtnHy1JkAqxGRzqBwiBJ/OIzLia9wGJ09kBFD+wOhxfg+/ESL8YnEO4VDlMTjYHQNM+Mz\nc2qPHrSchkj8UzhEwbGyU5QcCc0BSElJ9v7KjieXzBpHcnLox6Vw10H2HGhqLUURiQcKhyjYurv2\nqGHU8AFd5s5vrdG3d3fOn5zjtd/VUt4icU3hEAXxPN7gd5lvYDpvZQFnz1YFWI2IdCSFQxRsi+Mr\nlfxmTMwmo08PIHQqbU3+nhb2EJGuSuHQTs45Cn2nlcbG8ZFDcnISueeP99o6tSQSvxQO7VRy5ATH\nT5wCoHt6GsMGxcU9i5rkP7W0auNujpadDLAaEekoCod28g9Gj8ke2OVXYm3J8MH9mDAqE4Dq6mre\nW1kYcEUi0hEUDu3kH28YF8fjDX7z59Yu0PvOss1ajE8kDikc2qnOkUOChMOF08fQLS20dlRR8ZE6\nn4GIxAeFQzs459hRdNhrj86Oj5VYW9I9PY0Lpo/22poxLRJ/FA7tcOjICcpPVQDQIz2Nwf17B1xR\n5/EvxvfBmm3e7VFFJD4oHNphx97ao4ZRWfE/GO03aXQmmQP7AKF7WSz7dEfAFYlINCkc2mFH0SHv\n8ajhiXFKqYaZcZlvMb53NOdBJK4oHNqhTjhkDQiwkmDknj+emmOlDYX7KD58PNB6RCR6FA7tsLPe\naaVEMzCjF9MnZXtt3SVOJH4oHNqorPx0nWW6hw/uF3BFwfCfWspbsUVzHkTihMKhjfxHDSOG9o/L\nZbojMXvKSHr16AaErt5aV7A34IpEJBoUDm20Y2/teMPIYYk33lAjNTWZS2aN89pvfpQfYDUiEi0K\nhzbyD0YnyuS3plx+wSTv8Yr1OzhWdirAakQkGhQObVRnjkOCXcZaX86wAd5Njqqqqnl3hQamRbo6\nhUMbnKk8y97wPZQNGDk8cU8r1Vhw4Tne4zc/2qSBaZEuTuHQBrv3lVId/uU3dFBf0rulBlxR8Oad\nN4Ye6WkAHDh0nA2F+wKuSETaQ+HQBv7B6JwEP6VUo1taKpf67hL3xoebAqxGRNpL4dAGdVZiTcDJ\nb0254kINTIvEC4VDG9S5jFXjDZ6cYQMYP3IIEBqY1npLIl2XwqGVqqur60yAS/TLWOvzD0y/9XG+\nBqZFuiiFQyvtKznGmcqzAGT06UG/3j0Crii2XDhjdJ2B6fWaMS3SJSkcWmmnb7xBp5Qa6paWSu7s\n2oHppZoxLdIlKRxayT/eMDprUICVxK7LL6g9tbR83Q6Olp0MsBoRaQuFQyv5l83I0ZFDo3KG9WfC\nqEwgNEbzzjLNmBbpaiIKBzNbaGabzazAzO5qos9DZlZoZmvNbIZv++NmVmxm6+r1zzCzpWa2xcze\nMLO+7XsrHc85V2fZDF3G2rQFvsta3/hwI1VV1QFWIyKt1WI4mFkS8DBwJTAZ+LKZTazX5ypgjHNu\nHHA78N++b/8uvG99dwNvOecmAO8AP2jTO+hEpcfKOX4idO1+erdU7x7K0tCFM8bQu2c6EFrKe+WG\nncEWJCKtEsmRw2yg0Dm3yzlXCTwDLKrXZxHwJIBzbjnQ18yGhNsfAEcaed5FwBPhx08A17W+/M7l\nP2oYOXwAZtZM78SWlppS57LWV/+2IcBqRKS1IgmH4cAeX7sovK25Pnsb6VPfYOdcMYBz7gAwOIJa\nAlVnmW6dUmrRgnnnkBQO0I1b97Fr3+EW9hCRWBFLA9IxP1tqZ5FmRrfGwIxezJk22mvr6EGk60iJ\noM9eYISvnRXeVr9Pdgt96is2syHOuWIzywQONtXx3nvv9R7n5uaSm5vbctUdoO5gtC5jjcTVl0zh\n47XbAHhvZQFf+9wcbyxCRKInLy+PvLy8qD2ftbS8gZklA1uA+cB+YAXwZedcvq/P1cB3nHPXmNlc\n4L+cc3N93x8JvOycm+rbdj9Q6py7P3wFVIZz7u5GXt/FwhIM5acquPnu3wGQlJTEUw/cRmpqYt43\nujWcc/zTA897p5RuXnQBiz4zLeCqROKfmeGca/PAaIunlZxzVcBiYCmwEXjGOZdvZreb2TfDfV4F\ndpjZVuBR4Nu+Ap8CPgLGm9luM/t6+Fv3A1eYWU3w/LStb6Iz+NdTys7MUDBEyMy45tIpXvu1v22g\nulqXtYrEukhOK+Gcex2YUG/bo/Xai5vY9ytNbC8FLo+szOD5B6NHaTC6VS6eOY4nX1rGiZMVlBwp\nY9XG3cyeOjLoskSkGbE0IB3T6t4zWoPRrZGWmsIVF9ROinvlvXXN9BaRWKBwiJCOHNrnyosmU3Py\nc0PhPnbvLw20HhFpnsIhAmfPVlFUXDuPT5extt6g/r2Zc+4or/3Ke+sDrEZEWqJwiMCeA0e8tYGG\nDOhDz+7dAq6oa7r6Uu9iNfJWFmi1VpEYpnCIwPaiEu+xxhva7pwxQxmTHZofcvZsFa+9vzHgikSk\nKQqHCPgvYx2p8YY2MzMWzZ/utV9/fwOnKyoDrEhEmqJwiMB2DUZHzdxzRzG4f28ATpys4J3lmwOu\nSEQao3BogXOuzpGDTiu1T3JyEp/NPddrv/zuOt3rQSQGKRxacODQce/UR59e3enft2fAFXV98+dO\npFeP0KD+wdIylq3bEXBFIlKfwqEFdU4p6R4OUZHeLZWFF0322i+9vZZYWD9LRGopHFqwy39KSeMN\nUXPVJVNISQmtT7VtTwnrClpaxFdEOpPCoQV1L2NVOERLv949uHxu7d1mn39jdYDViEh9CocW1L2M\nVYPR0XTd/OkkJYV+BDdt28/GrfsCrkhEaigcmnG07CRHjodm8aalpjBsUN+AK4ovg/r3Jvf88V77\n+TfWBFiNiPgpHJqxfU/d24LW/JUr0XP9FTO8BfnWFRRRsLM40HpEJES/7Zrhv1JptAajO8TQQX25\naOY4r/2XpTp6EIkFCodm7NhTOxg9Olvh0FFuWHCe93jVxl11lkcXkWAoHJpR98hhUICVxLfszAzm\n+pbzfvqVlQFWIyKgcGjSiZMVHCwtA0JLPmRnZgRcUXy78apZ3tjD6k27NPYgEjCFQxP8pzZGDO3v\nTdiSjpEzbAAXnjfWaz/1yooAqxERhUMTNBjd+b501SySwsuTrC/Yy3rNmhYJjMKhCf6Z0Rpv6BzD\nB/cjd/YEr/3UKyu05pJIQBQOTdjhm+OgK5U6zxcXziQ5OfRjWbCzmNWbdgdckUhiUjg04nRFJfsO\nHgXAgJxh/YMtKIEM7t+bBRee47X/+PJyqqt1vweRzqZwaMTOvYepOZmRlZlBt7TUQOtJNDcsOM/7\nzPfsL+Wd5VsCrkgk8SgcGlFnJVYNRne6jD49uG7+NK/99Csrda9pkU6mcGjEjqLalVg1GB2May+b\nRkafHkBoAcQX3l4bcEUiiUXh0Ig6l7FqMDoQ6d1S+co1s732S2+v5fDREwFWJJJYFA71VFZWsXt/\nqdceOVz3cAhK7uzx5AwLff6VZ6t4SstqiHQahUM9u/eXelfHZA7sQ8/u3QKuKHElJSVx63UXeO33\nVmzRshoinUThUE/dwWiNNwTt3AlZzJqcA4ADHnvufV3aKtIJFA71FO466D3Wshmx4Rs3zCM1vLbV\njqJDLP0wP+CKROKfwqEefziMHzk4wEqkxpABfbj+ihle+09/Xc6xslMBViQS/xQOPqcrKtkTHow2\nYEy2TivFiuvmT2fIgD4AnDx9hj++vDzgikTim8LBZ9uektqZ0UP70z09LdB6pFZaagq33TDPa7+z\nfDP52/YHWJFIfFM4+GzdXTsYPW6ETinFmpmTc5g9daTXfuTpPM5Ung2uIJE4FlE4mNlCM9tsZgVm\ndlcTfR4ys0IzW2tm01va18zuMbMiM1sT/lrY/rfTPv7LJMflKBxi0W03XER6t9C6S/tKjvHn11YF\nXJFIfGoxHMwsCXgYuBKYDHzZzCbW63MVMMY5Nw64HfhVhPs+6Jw7L/z1ejTeUHts3V07GK1wiE0D\nM3pxy6LauQ8vvr2Wrb6LCEQkOiI5cpgNFDrndjnnKoFngEX1+iwCngRwzi0H+prZkAj2NWLEkeMn\nOXQktDxDakqy7hkdw664cBJTxg0DQnMffvl0HmfPVgVblEiciSQchgN7fO2i8LZI+rS07+Lwaajf\nmFnfiKvuAHXmN2QP0j2jY5iZ8a0vXUpaagoQmtX+3NI1AVclEl9SOuh5IzkieAT4kXPOmdmPgQeB\n2xrreO+993qPc3Nzyc3NjUKJdflPTYzXKaWYN3RQX75yzWx+/+JHAPzljdXMmJjNxNGZAVcmEoy8\nvDzy8vKi9nyRhMNeYISvnRXeVr9PdiN90pra1zlX4tv+a+Dlpgrwh0NH8R85jFU4dAnXXDqFFet3\nsGnbfhzwX0++zX/e9QWthyUJqf4fzvfdd1+7ni+S00orgbFmlmNmacBNwJJ6fZYANwOY2VzgqHOu\nuLl9zcz/J971wIZ2vZN2cM5pMLoLSkpK4o6/m0+P8HyUkiNlPPbc+wFXJRIfWgwH51wVsBhYCmwE\nnnHO5ZvZ7Wb2zXCfV4EdZrYVeBT4dnP7hp/6ATNbZ2ZrgUuB70X3rUVuX8kxTp4+A0DvnukM7t87\nqFKklQZm9OJbN13qtT9YvZX3VhYEWJFIfIhozCF8memEetserddeHOm+4e03R15mxyr0zW8YnzME\ns5i5iEoiMG/GGD7J38274XtNP/rn9xkzYhBZQ3TFmUhbaYY09ccbtJ5SV3Tb9fPIHBhae6niTCU/\ne3wpp8JHgyLSegoHoMAfDlo2o0vqnp7G97++wFvau6j4CL98+j2ccy3sKSKNSfhwOHX6DDv2hC6c\nMmDCqCHBFiRtNiprIN/60iVe++O123g5b12AFYl0XQkfDgW7DnorsWYP7a/LILu43NkTuHLeZK/9\nh5eW8emWogArEumaEj4cNvmWfT5nzNAAK5Fo+cb1F3qXI1c7x89+u5Td4ft0iEhkEj4cNm+vDYdJ\noxUO8SAlJZl//sYC+vftCYROHf7Ho69xtOxkwJWJdB0JHQ5nz1axZUftZaxaeiF+DOjXix9+8yq6\npYWW9y45UsZPHnudijOVAVcm0jUkdDhsLzpEZXg1z0EZvRmY0SvgiiSaRmUN5M5bL/cW+tq6+yAP\n/v4treAqEoGEDgf/eMOkMTpqiEezJudw2xcu8tqrNu7ioT+9S3V1dYBVicS+hA6HDYW16wfW3B9A\n4s9VF0/h8/O9mxPy4Zqt/OrZv2kOhEgzEjYczp6tYuPW2iOHqeOzAqxGOtpXPzeHhRfVXuL69rLN\n/PZ/PlRAiDQhYcOhYNdB7+b0Qwb00WJ7cc7M+PsvXETu7Nplvl792wYee+59BYRIIxI2HNYX6JRS\nojEzvn3TpcydNtrbtvTDTTz0x3eoqtIYhIifwgE4V6eUEkZychJ33nI5F88c523726pCHvz9m1RW\n6iomkRoJGQ6nKyop2FU7v2HKeB05JJLk5CTu+LvPcMWFk7xty9bt4N5HXub4iVMBViYSOxIyHNYV\n7PVOI4wY2p9+vXsEXJF0NjPj9hsv4XO553rbNm8/wA9+/gJ7Dx4NsDKR2JCQ4bB64y7v8cxzRjTT\nU+KZmXHLdRdw86ILvIlyBw4d5wcPvlDntKNIIkq4cHDOsWbTbq89c3JOgNVI0MyMRZ+Zxve/UXsv\niPJTFdz3y5f5y5trdCWTJKyEC4edew9TeqwcgJ7duzF+pO7fIDB32mh+/N1F9O3dHQAHPPXXFfzk\nsdc5cbIi2OJEApBw4bDad9QwfVI2yckJ9xFIE8bmDOZn37+BCaNql1JZvWkX33/geTZu3RdgZSKd\nL+F+M66fs9akAAAKSUlEQVRcv9N7PGuyxhukrgH9evGjxZ+rM1BdcqSMe/7fEp548WNv4qRIvEuo\ncDhYWsbW3aH7RSclJTFjksJBGkpJSebWz1/I97++wLszoAOWvPsp//yzv7B5+4FgCxTpBAkVDh+u\n2eo9nj4xi9490wOsRmLdBdNH8/O7v8j0idnetqLiI/zLL17kkafzKCs/HWB1Ih0rscLhk23e4wun\njwmwEukqBvTrxb9+62q++cWLSUtN8ba/vWwz//jvz7D0w01aekPiksX6pXpm5qJR4/6SYyz+8dNA\naIbs7/79Fu+UgUgkDpaW8du/fMjKDTvrbM8aksHXrp3DrMk5mFnjO4t0MjPDOdfmH8iECYdnXlvJ\nc6+vBuD8KSO5+38tbPdzSmJasX4nj//lAw4dOVFn+4RRmdy4cCbTJmQpJCRwCocIVFVV8637/uTN\nb7jz1iuYN0OnlaTtKs5U8nLeel546xNOV9S9L/W4nMHcsOA8HUlIoBQOEVixfif3/+Z1APr06s6v\n7/saKeHZsCLtcazsFM8vXc0bjYw9DBvUl6sumcJlsyfQPT0toAolUSkcIvDjX73CJ/l7ALj+8hl8\n9XNzolGaiKektIwX3lrL28s3c/Zs3aW/u6encdF5Y7hs9gTGjxyiownpFAqHFuw5cITv/eRZHGDA\nL//tKwwZ0Cdq9Yn4lR4r5+V31/HWx/mcPH2mwfeHDerLpbMncMmscbr7oHQohUMLHnziLW9+w8xz\ncvjh7VdFqzSRJp2uqCRvRQGvvLeOfSXHGu0zLmcw508dyeypo8ga0k9HFBJVCodm7N5fyp0//TM1\ne//0zs8zLkcL7Unncc6Rv/0A7y7fwkdrtzUYvK6RObAPs6eOYtrELCaOyiS9W2onVyrxRuHQBOcc\n9//mDe+adB01SNAqzlSyYt1O3l2xhXVbimjqpzo5OYmxIwYzddwwpowbzricwQoLaTWFQxM+/GQb\nD/7+Ta99/53XMzZncDRLE2mzsvLTrN64i5Xrd/LJ5iIqzjR+RAGhsbKszAzG5gxmbPZgxo4YxIhh\n/evM2BapT+HQiGNlp7jjJ896a99cceEkvvWlSzuiPJF2O1N5lvUFe/kkfw8bCvey58CRFvcxYOig\nvowY2p+sof0ZMbQ/wwb1ZfCA3pr5L0AnhYOZLQT+i9BaTI875+5vpM9DwFVAOXCrc25tc/uaWQbw\nLJAD7ARudM41GLlrbTicOn2Gex5+mW17SgAY0K8nP7/7Rv2DkS7jWNkpNmzdx4bCveRvP0DR/tIm\nT0E1plePbmQO7MuQgX3IHNCHIQN7k9GnJwP69aR/35706tFNg98JoMPDwcySgAJgPrAPWAnc5Jzb\n7OtzFbDYOXeNmc0BfuGcm9vcvmZ2P3DYOfeAmd0FZDjn7m7k9SMOh9Jj5Tz4+7fI374/tC/wL9+6\nhhmTspvfMYHk5eWRm5sbdBlxobM+y9MVlWzbU8LW3SVs3X2Q7XtKKD50vFWB4ZecnET/Pj3p368n\n/Xp3p1ePbvTpmU6vnun07tmN3j27h9vd6N0jne7pqaSmJHd4oOhnM7raGw6RnLScDRQ653aFX/AZ\nYBGw2ddnEfAkgHNuuZn1NbMhwKhm9l0E1JzreQLIAxqEQyROnKwgb8UWnl+6ps4yyt+88RIFQz36\nBxg9nfVZpndLZfLYYUweO8zbVnGmkr3FR9m9v5Td+0spOnCU4sPHOXD4eINJePVVVVVTcqSMkiNl\nEdeQZEZ6t1S6p6fSvVtanf+md0slPS2VbmkppKYkk5qaTFpqCqkpSeH/+rclk5aSTEpKEinJySQl\nGUlJSSQnGa+89gZTps0iOTmJlOQkkpKM5KTQ4+TkJB3tdLJIwmE4sMfXLiIUGC31Gd7CvkOcc8UA\nzrkDZtbkaPF/PPoa1a6aqipHtaumutpRVe2orq6mrPx0g7+iDPjatXNZMO+cCN6eSNfTLS2V0dmD\nGJ09qM525xyHj5ZTfPg4xYeOc+DQcUqOlHHk+ElKj5ZTevwkpxqZnNeSauc4efpMeGJfeZTeRV2b\nlq1n++k/Nvl9A5KSk0hOCoVFkhlmob+QLfw4yffYzDDC25PC26nbv/Zx7f7g699EIPk3h5614fa6\n/VvZh4ad6rxmC883a/LIxl+kFTrqcoe2RHyTR8mrN+2K+Eky+vTgjr+bz9Txw9tQgkjXZmYMzOjF\nwIxedY40/E5XVFJ6rJzSY+UcO3GaE+WnKTtZQdmJ0xwvP8WJkxWUlZ+mrPw0J05WcKqiMibuWeEI\nHfVUVVVD0xd3CZA5sG/7n8Q51+wXMBd43de+G7irXp9fAV/ytTcDQ5rbF8gndPQAkAnkN/H6Tl/6\n0pe+9NX6r5Z+vzf3FcmRw0pgrJnlAPuBm4Av1+uzBPgO8KyZzQWOOueKzexQM/suAW4F7gduAV5q\n7MXbM6AiIiJt02I4OOeqzGwxsJTay1Hzzez20LfdY865V83sajPbSuiE5Neb2zf81PcDfzazbwC7\ngBuj/u5ERKRNYn4SnIiIdL6koAtoipktNLPNZlYQngchrWRmO83sUzP7xMxWhLdlmNlSM9tiZm+Y\nWRRGruKTmT1uZsVmts63rcnPz8x+YGaFZpZvZguCqTp2NfF53mNmRWa2Jvy10Pc9fZ5NMLMsM3vH\nzDaa2Xoz+254e9R+PmMyHMKT5x4GrgQmA182s4nBVtUlVQO5zrkZzrmaS4jvBt5yzk0A3gF+EFh1\nse93hH4G/Rr9/MzsHEKnRicRWingEdOF+fU19nkCPOicOy/89TqAmU1Cn2dzzgJ3OucmAxcA3wn/\njozaz2dMhgO+iXfOuUqgZvKctI7R8P/xIkKTDgn/97pOragLcc59ANRf6Kipz+9a4Bnn3Fnn3E6g\nkIbzgRJaE58nNH7p+yL0eTbJOXegZoki59wJQld/ZhHFn89YDYemJtVJ6zjgTTNbaWZ/H95WZ/Ih\noKVqW2dwE59f/Z/ZvehnNlKLzWytmf3GdxpEn2eEzGwkMB1YRtP/vlv9ecZqOEh0zHPOnQdcTeiw\n82JCgeGnKxLaR59f+zwCjHbOTQcOAP8ZcD1dipn1Ap4H7ggfQUTt33eshsNeYISvnRXeJq3gnNsf\n/m8J8CKhw8ji8LpXmFkmcDC4Crukpj6/vYB/IS/9zEbAOVfiW1nz19Se6tDn2QIzSyEUDH9wztXM\nE4vaz2eshoM38c7M0ghNnlsScE1dipn1CP9VgZn1BBYA66mdfAjNTD4Uj1H3nHhTn98S4CYzSzOz\nUcBYYEVnFdmF1Pk8w7/AalwPbAg/1ufZst8Cm5xzv/Bti9rPZ0zeSqqFyXMSmSHAC2bmCP1//pNz\nbqmZrUKTDyNiZk8BucAAM9sN3AP8FHiu/ufnnNtkZn8GNhFa+efbbbqFYRxr4vO8zMymE7qybidw\nO+jzbImZzQO+Cqw3s08InT76IU1MLm7L56lJcCIi0kCsnlYSEZEAKRxERKQBhYOIiDSgcBARkQYU\nDiIi0oDCQUREGlA4iIhIAwoHERFp4P8DJk4Ndblud50AAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "pmf2 = MakePmfUsingNbinom(5, 0.1, 200)\n", "thinkplot.Pdf(pmf2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "And confirmation that the results are the same within floating point error." ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "8.6736173798840355e-17" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "diffs = [abs(pmf[n] - pmf2[n]) for n in pmf]\n", "max(diffs)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Using the PMF, we can compute the mean and standard deviation:" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "(49.998064403376738, 21.207570382894403)" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "pmf.Mean(), pmf.Std()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "To compute percentiles, we can convert to a CDF (which computes the cumulative sum of the PMF)" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAEACAYAAAC9Gb03AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAEctJREFUeJzt3X+s3fVdx/HXq4OSOCaZI+tiO5jCFtziIAvWLpt6BgZu\nidldRiJQ4wRdUpRqExNWNDGcP5YIMSYD62w6r3Nsknau225dRuhid2Iw/CgRWjZ6oUjo2gId002z\nJdOuvv3jfFvOPT2/7/ec7/f7+T4fyQ3nx/ee8+Gbc1/33ff38/lcR4QAAOlaVfQAAADTRdADQOII\negBIHEEPAIkj6AEgcQQ9ACRuaNDbXrB90vahAcfcb/uI7adtX5XvEAEAKzFKRf9ZSdf3e9L2RkmX\nRcQ7JW2WtCOnsQEAcjA06CPiEUnfH3DIvKQHsmMfl3SR7TX5DA8AsFJ59OjXSjrWcf9E9hgAoAS4\nGAsAiTsvh9c4IentHffXZY+dwzYb6wDABCLCk37vqEHv7KuXvZLukLTb9gZJP4iIk/1eiE3U8tNs\nNtVsNoseRjImPZ+L+w9q90NP6n/+91T+g6qoZx/7mt694TeKHkYyvnz/76/o+4cGve0HJTUkvcX2\ndyTdLWm1pIiInRHxdds32H5B0o8k3baiEQElVJYwv2D1+bpp49Wav+bKQscxTLP5qprN24seRjI8\n7aCPiE0jHLNlRaMASmSaoV6VoEZa8ujRoyCNRqPoIVReZ6i/dvy/9czWyZaBEODL8dksF8+yZ247\n6NGjaJNW7IQ5imJ7JhdjgcqaJNgJdaSEoEeSxgl3Qh2pI+iRDMId6I2gR+WNEvAEO+qMoEclEe7A\n6Ah6VMqwgCfcgXMR9KgEAh6YHEGPUhsU8IQ7MBqCHqW1uP+gHlh89JzHCXhgPAQ9SqdfFU/AA5Mh\n6FEag9o0H5t/PwEPTIigRynQpgGmh6BHoWjTANNH0KMQtGmA2SHoUQiqeGB2CHrM3OL+g8tCnoAH\npougx8z0atdcsPp8PfgXv1fgqID0rSp6AKiPXu2amzZeXdBogPqgosfU9avkadcAs0HQY6p6zY+n\nXQPMFq0bTNXuh55cdv9MJQ9gdqjoMRW92jXMjweKQdAjV4NWuhLyQDFo3SBXgxZCASgGFT1yw0Io\noJwIeuSm88IrM2uA8iDosWK9+vK0aoDyoEePFeu1GIp2DVAeVPSY2KAVrwDKg6DHxNigDKgGWjeY\nSL8ZNgDKh4oeE2GGDVAdBD3GwgwboHpo3WAszLABqmekoLc9Z3vJ9vO2t/V4/qdt77X9tO1nbN+a\n+0hRqMX9B7XpzgX68kAFDW3d2F4labukayW9LOmA7cWIWOo47A5J346ID9u+WNJztr8QET+Zyqgx\nc8ywAaprlIp+vaQjEXE0Ik5J2iVpvuuYkPSm7PabJP0HIZ8WKnmguka5GLtW0rGO+8fVDv9O2yXt\ntf2ypAsl3ZTP8FAGi/sPLrtPJQ9US16zbq6X9FREXGP7MknfsP3eiPhh94HNZvPs7UajoUajkdMQ\nMC3dUykBTFer1VKr1crt9RwRgw+wN0hqRsRcdv8uSRER93Yc8zVJfx4R/5rd/2dJ2yLiya7XimHv\nh/Lgr0QB5WBbEeFJv3+UHv0BSZfbvtT2akk3S9rbdcxRSb+eDWiNpHdJenHSQaEcmEoJpGFo6yYi\nTtveImmf2r8YFiLisO3N7adjp6RPSvp724eyb/tERPzn1EaNmeACLJCGoa2bXN+M1k0l9GrZ7Lnv\n9gJHBNTbLFo3qJleLRsA1UXQYxl2pQTSw6ZmWIZdKYH0UNFjGXalBNJDRQ9Jr1+A7cRUSiANVPSQ\nxAVYIGUEPSQxZx5IGa2bmuvVsuECLJAWKvqao2UDpI+grzlaNkD6aN3UGPvMA/VARV9j7DMP1ANB\nX2MsjgLqgdZNDbE4CqgXKvoaYqYNUC8EfQ0x0waoF1o3NcLiKKCeqOhrhJYNUE8EfY3QsgHqidZN\nTdGyAeqDir4mulfBAqgPgr4mWAUL1BdBXxOsggXqix594lgFC4CKPnFMqQRA0CeOKZUAaN3UCFMq\ngXoi6BPVqzcPoJ5o3SSK3jyAMwj6RNGbB3AGrZsE8bdgAXSiok8Qq2ABdCLoE8QqWACdCPrEsQoW\nAD36hDClEkAvI1X0tudsL9l+3va2Psc0bD9l+1u2v5nvMDEKplQC6GVoRW97laTtkq6V9LKkA7YX\nI2Kp45iLJP21pOsi4oTti6c1YPTHlEoAvYzSulkv6UhEHJUk27skzUta6jhmk6Q9EXFCkiLie3kP\nFONhSiWAM0YJ+rWSjnXcP652+Hd6l6Tzs5bNhZLuj4jP5zNEDENvHsAgeV2MPU/S+yRdI+mNkh61\n/WhEvJDT62MAevMABhkl6E9IuqTj/rrssU7HJX0vIn4s6ce2/0XSlZLOCfpms3n2dqPRUKPRGG/E\nOAe9eSAtrVZLrVYrt9dzRAw+wH6DpOfUvhj7iqQnJN0SEYc7jrlC0l9JmpN0gaTHJd0UEc92vVYM\nez+M78atO87e3nPf7QWOBMA02FZEeNLvH1rRR8Rp21sk7VN7OuZCRBy2vbn9dOyMiCXbD0s6JOm0\npJ3dIY/p6N7XBgC6Da3oc30zKvrcbbpz4Wzr5oLV5zPbBkjQSit6tkCoOPa1ATAMQZ8Q9rUB0At7\n3VQUc+cBjIqKvqKYOw9gVAR9RTF3HsCoaN0kgJk2AAYh6CuG3jyAcdG6qRh68wDGRdBXDL15AOOi\ndVNh9OYBjIKKvkLY1wbAJAj6Cum8CEtvHsCoCPoKYV8bAJMg6CuKfW0AjIqLsRXA3HkAK0FFXwHM\nnQewEgR9BTB3HsBK0LqpGObOAxgXFT0AJI6KvsS4CAsgD1T0JcZFWAB5IOhLjIuwAPJA66YiuAgL\nYFJU9CXFBmYA8kLQlxQbmAHIC0FfUmxgBiAvBH0FsIEZgJUg6AEgccy6KRkWSQHIGxV9ybBICkDe\nCPqSYZEUgLzRuikxFkkByANBXxL05gFMC62bkqA3D2BaCPqSoDcPYFpo3ZQQvXkAeRqporc9Z3vJ\n9vO2tw047pdsn7L90fyGCABYiaFBb3uVpO2Srpf0Hkm32L6iz3H3SHo470Gmjp0qAUzTKBX9eklH\nIuJoRJyStEvSfI/j/lDSlyR9N8fx1QI7VQKYplGCfq2kYx33j2ePnWX7ZyV9JCL+RpLzG149sFMl\ngGnKa9bNpyR19u4J+wmxUyWAvI0y6+aEpEs67q/LHut0taRdti3pYkkbbZ+KiL3dL9ZsNs/ebjQa\najQaYw45HSySAtBLq9VSq9XK7fUcEYMPsN8g6TlJ10p6RdITkm6JiMN9jv+spH+KiC/3eC6GvV+d\nbLpz4Zz580ytBNDNtiJi4k7J0Io+Ik7b3iJpn9qtnoWIOGx7c/vp2Nn9LZMOpm5YJAVgFoZW9Lm+\nGRX9Mjdu3XH29p77bi9wJADKbKUVPVsgAEDi2AKhAFyEBTBLVPQFYKdKALNE0BeAi7AAZonWTcGY\nTglg2qjoASBxVPQzxEVYAEWgop8hLsICKAJBP0NchAVQBFo3BeEiLIBZoaKfEf6KFICiEPQzwl+R\nAlAUgn5G+CtSAIpC0BeAvyIFYJYIegBIHLNupoxFUgCKRkU/ZSySAlA0gn7KWCQFoGi0bmaIRVIA\nikDQTwm9eQBlQetmSujNAygLgn5K6M0DKAtaNzNAbx5AkajoASBxBP0UsFMlgDIh6KeAnSoBlAlB\nPwXsVAmgTAj6KWOnSgBFY9ZNjlgkBaCMqOhzxCIpAGVE0OeIRVIAyojWzZSwSApAWRD0OaA3D6DM\naN3kgN48gDIj6HNAbx5AmdG6yRm9eQBlM1JFb3vO9pLt521v6/H8JtsHs69HbP9i/kMFAExiaNDb\nXiVpu6TrJb1H0i22r+g67EVJvxoRV0r6pKTP5D3QsmIDMwBlN0pFv17SkYg4GhGnJO2SNN95QEQ8\nFhH/ld19TNLafIdZXmxgBqDsRgn6tZKOddw/rsFB/nFJD61kUFXCBmYAyi7Xi7G2PyTpNkkf7HdM\ns9k8e7vRaKjRaOQ5hEKxgRmAPLRaLbVardxezxEx+AB7g6RmRMxl9++SFBFxb9dx75W0R9JcRPx7\nn9eKYe9XFWcWSXVW9Hvuu73AEQFIlW1FhCf9/lFaNwckXW77UturJd0saW/XIC5RO+R/u1/Ip4ZF\nUgCqYmjrJiJO294iaZ/avxgWIuKw7c3tp2OnpD+T9DOSPm3bkk5FxPppDrxoLJICUBVDWze5vllC\nrZsbt+44e5uWDYBpWmnrhpWxY2IDMwBVw143Y6I3D6BqCPox0ZsHUDW0bsbQvd0BG5gBqAIq+jGw\n3QGAKiLox8B2BwCqiKCfENsdAKgKevQjYEolgCqjoh8BUyoBVBlBPwKmVAKoMlo3Y2JKJYCqIegH\noDcPIAW0bgagNw8gBQT9APTmAaSA1k0fbHcAIBVU9H2w3QGAVBD0fbDdAYBU0Lrp0mumDdsdAKgy\nKvouzLQBkBqCvgszbQCkhtZNplfLhpk2AFJARZ+hZQMgVQS92tU8LRsAqaJ1o3PnzNOyAZASKnox\nZx5A2mpd0TNnHkAd1Lqi5wIsgDqoZUV/ppLnAiyAOqhl0PcKeS7AAkhV7Vo3TKUEUDe1q+iZSgmg\nbmoT9L368lTyAOqgNq2bXn15plICqIPkK3pm2ACou6SDfnH/QT2w+Oiyx+jLA6ibkYLe9pykT6nd\n6lmIiHt7HHO/pI2SfiTp1oh4Os+BjqNXFS9RyQOop6FBb3uVpO2SrpX0sqQDthcjYqnjmI2SLouI\nd9r+ZUk7JG2Y0pj76hfwkvSx+fcn15NvtVpqNBpFDyMZnM/8cC7LZZSLseslHYmIoxFxStIuSfNd\nx8xLekCSIuJxSRfZXpPrSPtY3H9Qm+5c0I1bd+iBxUd7VvEphrzU/mFCfjif+eFclssorZu1ko51\n3D+udvgPOuZE9tjJFY2uy6CKvduZNk2KAQ8A45j5xdgbt+6Y6usT8ACwnCNi8AH2BknNiJjL7t8l\nKTovyNreIembEbE7u78k6dci4mTXaw1+MwBATxHhSb93lIr+gKTLbV8q6RVJN0u6peuYvZLukLQ7\n+8Xwg+6QX+lAAQCTGRr0EXHa9hZJ+/T69MrDtje3n46dEfF12zfYfkHt6ZW3TXfYAIBRDW3dAACq\nbWZ73dies71k+3nb22b1vqmw/ZLtg7afsv1E9tibbe+z/Zzth21fVPQ4y8r2gu2Ttg91PNb3/Nn+\nE9tHbB+2fV0xoy6vPufzbtvHbf9b9jXX8Rznsw/b62zvt/1t28/Y/qPs8fw+nxEx9S+1f6G8IOlS\nSedLelrSFbN471S+JL0o6c1dj90r6RPZ7W2S7il6nGX9kvRBSVdJOjTs/El6t6Sn1G5tviP77Lro\n/4cyffU5n3dL+uMex/4C53PguXybpKuy2xdKek7SFXl+PmdV0Y+y6AqDWef+C2xe0uey25+T9JGZ\njqhCIuIRSd/verjf+fuwpF0R8ZOIeEnSEZ27dqTW+pxPqf057TYvzmdfEfFqZFvGRMQPJR2WtE45\nfj5nFfS9Fl2tndF7pyIkfcP2Adsfzx5bE9nspoh4VdJbCxtdNb21z/nrtwAQw22x/bTtv+1oNXA+\nR2T7HWr/S+kx9f/5Hvt81mY/+gR8ICLeJ+kGSXfY/hW1w78TV9ZXhvO3Mp+W9PMRcZWkVyX9ZcHj\nqRTbF0r6kqStWWWf28/3rIL+hKRLOu6vyx7DiCLiley/r0n6qtr/VDt5Zk8h22+T9N3iRlhJ/c7f\nCUlv7ziOz+sIIuK1yJrIkj6j19sJnM8hbJ+ndsh/PiIWs4dz+3zOKujPLrqyvVrtRVd7Z/TelWf7\np7Lf9rL9RknXSXpG7XN4a3bY70ha7PkCOMNa3kPud/72SrrZ9mrbPyfpcklPzGqQFbLsfGZhdMZH\nJX0ru835HO7vJD0bEfd1PJbb53Mme91En0VXs3jvRKyR9JVsC4nzJP1DROyz/aSkL9r+XUlHJf1m\nkYMsM9sPSmpIeovt76g9Q+QeSf/Yff4i4lnbX5T0rKRTkv6go1KF+p7PD9m+StL/SXpJ0maJ8zmM\n7Q9I+i1Jz9h+Su0WzZ+qPevmnJ/vSc4nC6YAIHFcjAWAxBH0AJA4gh4AEkfQA0DiCHoASBxBDwCJ\nI+gBIHEEPQAk7v8BQAWBA2MLe/oAAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "cdf = Cdf(pmf)\n", "scale = thinkplot.Cdf(cdf)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "And here are the 10th and 90th percentiles." ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "(26, 78)" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "cdf.Percentile(10), cdf.Percentile(90)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": true }, "source": [ "Copyright 2016 Allen Downey\n", "\n", "MIT License: http://opensource.org/licenses/MIT" ] }, { "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.11" } }, "nbformat": 4, "nbformat_minor": 0 }