Skip to content

Instantly share code, notes, and snippets.

@statkclee
Created November 30, 2015 05:45
Show Gist options
  • Save statkclee/f81a55631769c906548d to your computer and use it in GitHub Desktop.
Save statkclee/f81a55631769c906548d to your computer and use it in GitHub Desktop.
Think Stats 2, 9장 연습문제
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"통계적 사고 (2판) 연습문제 ([thinkstats2.com](thinkstats2.com), [think-stat.xwmooc.org](http://think-stat.xwmooc.org))<br>\n",
"Allen Downey / 이광춘(xwMOOC)"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"%matplotlib inline\n",
"from __future__ import print_function\n",
"\n",
"import thinkstats2\n",
"import thinkplot\n",
"\n",
"import math\n",
"import random\n",
"import numpy as np\n",
"\n",
"from scipy import stats\n",
"from estimation import RMSE, MeanError"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 연습문제 8.1\n",
"\n",
"이번 장에서, $\\mu$를 추정하는데 xbar 와 중위수를 사용했고, xbar가 MSE 하한을 산출함을 알아냈다. 또한, $\\sigma$를 추정하는데 $S^2$ 와 $S_{n-1}^2$ 을 사용했고, $S^2$ 은 편향되었고, $S_{n-1}^2$은 불편향임을 알아냈다.\n",
"유사한 실험을 실행해서, xbar와 중위수가 $\\mu$의 편향된 추정값임을 알아내라. 또한, $S^2$ 혹은 $S_{n-1}^2$ 가 MSE 하한을 산출하는지 검사하라."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Experiment 1\n",
"mean error xbar -0.000501853563251\n",
"mean error median -5.92353811989e-05\n",
"Experiment 2\n",
"RMSE biased 0.515294279494\n",
"RMSE unbiased 0.577482542019\n"
]
}
],
"source": [
"\"\"\"Mean error for xbar and median as estimators of population mean.\n",
"\n",
" n: sample size\n",
" m: number of iterations\n",
"\"\"\"\n",
"n = 7\n",
"m=100000\n",
"\n",
"mu = 0\n",
"sigma = 1\n",
"\n",
"means = []\n",
"medians = []\n",
"for _ in range(m):\n",
" xs = [random.gauss(mu, sigma) for i in range(n)]\n",
" xbar = np.mean(xs)\n",
" median = np.median(xs)\n",
" means.append(xbar)\n",
" medians.append(median)\n",
"\n",
"print('Experiment 1')\n",
"print('mean error xbar', MeanError(means, mu))\n",
"print('mean error median', MeanError(medians, mu))\n",
"\n",
"\"\"\"RMSE for biased and unbiased estimators of population variance.\n",
"\n",
"n: sample size\n",
"m: number of iterations\n",
"\"\"\"\n",
"mu = 0\n",
"sigma = 1\n",
"\n",
"estimates1 = []\n",
"estimates2 = []\n",
"for _ in range(m):\n",
" xs = [random.gauss(mu, sigma) for i in range(n)]\n",
" biased = np.var(xs)\n",
" unbiased = np.var(xs, ddof=1)\n",
" estimates1.append(biased)\n",
" estimates2.append(unbiased)\n",
"\n",
"print('Experiment 2')\n",
"print('RMSE biased', RMSE(estimates1, sigma**2))\n",
"print('RMSE unbiased', RMSE(estimates2, sigma**2))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"1. xbar 와 중위수가 m이 증가함에따라 더 낮은 평균오차를 산출한다. 실험으로부터 식별할 수 있는한 어떤 것도 분명히 편향되지 않았다.\n",
"1. 분산에 대한 불편 추정량이 불편 추정량보다 더 낮은 약 10% 낮은 RMSE를 산출한다. 그리고 m이 증가함에 따라 차이는 유지된다."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 연습문제 8.2 \n",
"모수 $\\lambda=2$를 갖는 지수분포에서 표본 $n=10$개를 추출한다고 가정하자. 실험을 1000번 모의시험하고 추정값 lamhat의 표본 분포를 도식화한다. 추정값의 표준오차와 90% 신뢰구간을 계산하라.\n",
"\n",
"다른 $n$ 값을 갖는 실험을 반복하고, $n$ 값과 표준오차를 도식화한다."
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"standard error 0.784058129255\n",
"confidence interval (1.2783641952567839, 3.5793120939847975)\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYEAAAEZCAYAAABxbJkKAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XmYFNXZBfBzGDYHBVQUWcUFEVwRJSioYzCKC+75FBUF\n10QxromKRkfjrolGccENiBskrmAEAoFRUIOiAqKDgogw7KAssjO83x9VjF093TPD0DW3qvv8noeH\nutV3qk/P0m/XreXSzCAiIrmplusAIiLijoqAiEgOUxEQEclhKgIiIjlMRUBEJIepCIiI5DAVAcka\nJAtJvuQvtya5miRr4Hn7kJyQ0F5Nsk2Gtn0ryef85TYkt5DMyN9tTX6PJLpUBGS7kexG8iOSK0gu\nJzmR5OEOopRd9GJmc81sJ3NwIYz/vHMq6kOygOS8KmzrfjO7PBO5SM4h+euEbTv7Hkl01HYdQOKN\nZEMA7wK4EsA/AdQDcDSADS7iOHjO0JDMM7PSDG7SkGXfI9l+2hOQ7bUfADOzYeZZb2ZjzOxLACC5\nD8lxJJeRXEryZZKNtn6x/+n0JpLT/KGJF0g2JTmS5EqSY0g29vtuHQ65nOR8kgtI3pgqVPLQCcki\nknf7eymrSI4muWtC/4tI/uDnvN3P1T3NtnclOdzPNwnAPkmPbyG5t798Msmv/OcsIXkDyXwAIwE0\n91/zKpLN/OGs10m+RHIlgD6JQ1wJLk31+kkOJvmXhHbZ3oa/jdYARvjPeVOK71Fz/3UtJzmT5GUJ\n2yok+U+SQ/y800l2Sv0rIXGiIiDb6xsApf4bUA+SO6focy+AZgDaA2gFoDDhMQNwFoDuANoBOBXe\nG+QtAHaH9zv6h6TtFQDYF8AJAG5O92adQi8Affzt1gVwEwCQ7ADgSf/xZgAaAWiOhOGlJE8CWAtg\nDwCXAOhbQd8XAFxhZg0BHABgvJmtBdADwAJ/OKahmS30+58G4F9m1gjAK2m2W4DUr9/S5TCz3gDm\nAjjVf85HUnQb6vdpBuAcAPeRPC7h8Z4AXoP3/RkOYECa1ywxoiIg28XMVgPoBu/N5zkAS0i+Q3J3\n//HvzOy/ZrbJzJYBeBTAsUmbecLMlprZAgATAHxsZlPNbAOAtwB0TOp/l5mtM7PpAAbBe/OuNCqA\nQWY2y8zWwxu6OtR/7BwAw83sIzPbBOAOpHkzJZkHr2jd4Wf4CsAQpB9m2QjgAJINzWylmX2xdVNp\n+n9kZsMBwM+Zql9Fr79awz0kWwE4CsDNZrbRzKYCeB7ARQndJpjZKP8YwssADqnOc0m0qAjIdjOz\nGWbW18xaATgQ3qfoxwDAH9oZ6g+FrATwEoBdkzaxOGF5XVJ7PYAdk/onHlCd6z9fVSxKep6t220O\noCTh9awDsDzNNnaDdywtOUM6ZwM4GcAcf0iqSyUZSyp5HCmeu6qvvyLNAfxoZmuStt0ioZ34c1kL\noH6mzlQSd/QDlIwys2/gfTI+0F91H4BSAAf6Qxy9UfnvXWWfZlsnLc+vRtRECwC0LHtycgeUL1Rb\nLQWwOUWGlMxsspmdAa94vA1vDwRIvaeRajgnVb90r38NgPyEx/aowra2WgBgF5KJBbc1qlaUJMZU\nBGS7kGznH+xs4bdbwRue+NjvsiO8N6dVfp8/ZuBpbye5A8kD4I3xD6tq3DTr3wDQk+SRJOvCO2aR\nsq9/ts6bAAr9DB0AXJzyycg6JC8g2cj/utXwCiLgfare1T+7qqJ8qdale/1TAJxMcmeSewC4Lunr\nFiPpIHbC65oH4CMA95OsR/JgeMc7Xk7VX7KHioBsr9UAfgVgEsmf4b35TwOw9ayVuwAcBmAlgBHw\n3nArOy/dkpaT+78PYBaAsQAeNrOxafpW9Km6rK8/rn8NvAOjC/zXtATpT3PtB6+4LQLwov8v3fNe\nCOB7fyjsCgAX+M85A95B1tkkfyTZLM1rTfWa0r3+lwBMBTAHwCj/9SR+7f3wCshPJG9IkbUXgDb+\n9+BNeMc9xqXJkfy1ElMM8zoRki8COAXAEjM7KE2fxwGcBG+MsU/CgTORAHpX4c4GUNvMtoT4PDsC\n+AnAvmb2Q1jPIxIFYe8JDIJ3KlxKJE+G94fWFt6npKdDziOSEsmeJPNJNgDwCIBpKgCSC0ItAmY2\nAd4nqnROg3cQEWY2CUBjkk3DzCSxF9au62nwDrDOhzdufl5IzyMSKa5vG9ECwdPdSuCdpbE4dXfJ\nZf79ePJC2vblADJyjx6ROInCgeHksx90sElEpIa43hOYD+82Alu1RIpzvkmqMIiIVIOZVXjdjesi\nMBze6XZD/SspV5hZyqGgONzttrCwEIWFhZX2Ky4uLltu3759iIlSq2pO18LMmamfQXUymhmWLl+F\n9Rs2YeZ3C1CrVi18OWMuNm8uRe3aeXj7vU+wT5umqFXL+9udXlzRBclVM/fbCWi939HBlemmEfD/\n1g7qsGf67ZUswyEHtkHzpqluFZWwKQDr1m1A545t0z5dosEvPIU+l15VeUfHqppzv32ao/keu9RA\notRYhW96qEWA5Gvw7hPTxL+b4Z0A6gCAmQ00s/f8uyzOgndBUd8w84jUpNLSLfj32M8wcdIMfDip\nGA0a1AdJ/Pzzukq/9utvKp1qIK3junlnY69YuQYtmu2CvfZsiteHLsI5550CANgxvx5O6n4Y8vKi\nMBocVDRmdxR0PbDyjo7FJWdVhFoEzKzSG3uZWb8wM4iEpbR0C9au24BlP64CADzw+FuY9f0irFu/\nMeUb/Zo16zPyvB3atcK3sxei6xHt0O1X7dFu3+ZosktDNGqYn/ZrZkwZg/PO6JqR55fs4no4KKsU\nFBS4jlAlyrltVv+8DiULvPvJ/VCyFI88NRzr1nkXE69cPhcfTn9gu7afn18PnTu2xaIlK3DwAXti\n910bYtPmUhzUfs+yT+u1ahHt9mmOOnWq9ycble9lZZSz5oV6xXCmkMyqGfBcHxOQin8G74z6FM+/\nPBZr1m7Ahg0bM/J8Jxx3KNrv1xJdO++PvFq1sHOjBqhXr05Gti2SDsnIHxgWiZRPp8zCQ0+8tc1f\nt+suO6G0dAtWrFyDay4/Be3btsC+e+2BHerXRa1a0Rt7F9lKRUBy3vxFP+KhJ97C3PnLsGTpyrT9\n9th9ZzRqmI/Zcxdjv72b4/Lex+OQDm1Qt67+jCS+9NsrOWnp8tV4cvB4zPx+MXaoXz9lnzatd8dN\nV52O/fZphgb5qfuIxJ2KgOSUz6fNxj1/ex1z5i2qtO8rTyffjl8k+6gISM64/IanKz3//qarT8dh\nB++NPVvuVkOpRNxSEZCs90PJUpx/5aNpH7/t+nNwXLcDsUP9ujWYSiQaVAQk66UqAPXq1sYd1/dE\nm1ZNdJqu5DQVAclKZoZxE77E2yM/KffYX27phWZN9KsvAqgISJaZPmMurrzxmbSPfzD8HuTl1Qpc\nLCaSy1QEJGt0PaV/hY+327dFJG+aJuKSioDE3qIlK3B234fSPt7l8HboeUKnrLnro0gmqQhI7KUr\nAGNevxP5O9Sr4TQi8aIiILG1cPFPOOeSh8utH/DA5eh40F4OEonEj4qAxM7sOYvR++q/p3xs4rv3\nVmk2JRHxqAhIbJgZHh34Lt4Y8XHKxwc93k8FQGQbqQhIbPz+T8/iy69/SPnYiFf6Y5fGO9ZwIpH4\nUxGQWJg4qThlAZgw4h7dr19kO6gISORNnvodbr77pcC6wj+di98ce4ijRCLZQx+hJNJmfb8Q1/Z/\nIbCuc6e2KgAiGaIiIJFlZri43xOBdY0aNsCjd/d1lEgk+6gISGR1O/W2cuvee638OhGpPhUBiaRf\nn1VYbt2H/76v5oOIZDkVAYmc8ROnY8OGjYF1rz17g6M0ItlNRUAiZdmPq3D7/a8G1j3w5wvRukUT\nR4lEspuKgERGaekWnN77gcC6gw9og6O7dHCUSCT7qQhIZDz/ythy65568HIHSURyh4qARMKs7xfi\nH8OKAus+GH6P7gUkEjIVAYmE5OsBrurbQ7OAidQA/ZWJcx99+k25deeffbSDJCK5R0VAnDIz/LFw\nSGCdhoFEao6KgDh11c3PBdoXn3echoFEapD+2sSpaV/NCbSv6P0bN0FEcpSKgDjT9ZT+gfb9t1/o\nKIlI7lIRECcefOKtcuuOOVIXhYnUtFCLAMkeJGeQnEny5hSPNyE5iuQUktNJ9gkzj0TDlOnfY/io\nTwPrXnryWkdpRHJbaDOLkcwDMADA8QDmA/iU5HAzK07o1g/AF2Z2K8kmAL4h+bKZbQ4rl7h150PD\nMPb9qYF1l5zfHXu3aeookUhuC3NPoDOAWWY2x8w2ARgK4PSkPgsBNPSXGwJYrgKQvf732bflCgAA\nXHpBdwdpRAQItwi0ADAvoV3ir0v0HIADSC4AMBWAxgSy1Pr1G3HjHYMD644+sgM+GH6Pm0AiAiDc\nieatCn36A5hiZgUk9wEwhuQhZrY6uWNhYWHZckFBAQoKCjKVU0JWWroF3c8uLLf+AZ0NJJJRRUVF\nKCoq2qavCbMIzAfQKqHdCt7eQKKjANwLAGb2HcnvAbQDMDl5Y4lFQOLlpruGlFunWcJEMi/5A/Jd\nd91V6deEORw0GUBbkm1I1gVwLoDhSX1mwDtwDJJN4RWA2SFmEgc++WxmoD3uzcp/MUWkZoS2J2Bm\nm0n2AzAaQB6AF8ysmOSV/uMDAdwHYBDJqfAK0p/M7MewMknN27QpeJz/wTt6o169Oo7SiEiyMIeD\nYGYjAYxMWjcwYXkZgJ5hZhC3vp+7JNDu2nl/R0lEJBVdMSyh6vuHAYG27g4qEi0qAhKKTZs2l7s3\nUPNmuzpKIyLpqAhIxpWWbkHBGXeUW//YPX0dpBGRiqgISMZ1P6ew3LqnHroCLfbYpebDiEiFQj0w\nLLnn2ZfGYNPG4BlBrzx9Hdq03t1RIhGpiIqAZEzyMQAAGPb8jWipYwEikaXhIMmIcRO+LLfu0gu6\nqwCIRJz2BCQj/vzAa4F2v0tPQq+zjnaURkSqSnsCst2WLFsZaJ96wuEqACIxoSIg223I0PGB9q3X\nnuUoiYhsKxUB2S4bN27G2yM/KWvn1c5zmEZEtpWKgGyX484MXhR2w+90KyiROFERkGob+8G0cuvO\nOKmzgyQiUl0qAlItZoY7HxwaWDfilfLXCYhItKkISLV0O/W2YLtLe+zSeEdHaUSkulQEZJt9/8Pi\ncus0X7BIPKkIyDa746FhgfaIV/prngCRmFIRkG22du2GQFvDQCLxpSIg22zRkp/Klh/4s4aBROJM\nRUC2yew5weMBHdq1cpRERDJBRUC2Se+r/x5o77rzTo6SiEgmqAhIlY0ePyXQbt1yN0dJRCRTVASk\nyu5+5J+B9pAnrnGUREQyRUVAqmT5T6sD7TNP6YK6dTUdhUjcqQhIlRR/Oz/Qvumq0xwlEZFMUhGQ\nKhn02jjXEUQkBCoCUiUzZpaULeu0UJHsoSIglfpuzqJA+5rLTnaUREQyTUVAKmRmuOjqxwPrDu6w\np6M0IpJpKgJSocFJ8weLSHZREZAKPf/y2ED7g+H3OEoiImFQEZC0ShYuD7R/1+dE5OXpV0Ykm+gv\nWtI697K/BtoXnnOMoyQiEhYVAUkp+QphAJo4RiQLqQhISvc++kag/d5rtztKIiJhCrUIkOxBcgbJ\nmSRvTtOngOQXJKeTLAozj1TdpM++DbQbNcx3lEREwhTaHcBI5gEYAOB4APMBfEpyuJkVJ/RpDOBJ\nACeaWQnJJmHlkao7/3ePBdqP3NXHTRARCV2YewKdAcwyszlmtgnAUACnJ/U5H8AbZlYCAGa2LMQ8\nUgUbNmzCD/OWBNYdefh+jtKISNjCLAItAMxLaJf46xK1BbALyfEkJ5PsHWIeqYKeve8PtB+79xJH\nSUSkJoR5Q3irQp86AA4D0B1APoCPSf7PzGaGmEvS+LL4B6xZsz6w7ohD93WURkRqQphFYD6AxNtN\ntoK3N5BoHoBlZrYOwDqSHwA4BEC5IlBYWFi2XFBQgIKCggzHlacGjQ60X33mekdJRKQ6ioqKUFRU\ntE1fQ7OqfGDfdiRrA/gG3qf8BQA+AdAr6cDw/vAOHp8IoB6ASQDONbOvk7ZlYeV0obi47FuA9u3b\nO0wS1PWU/mXL+fn1MOZfdzpME66o/gxEMokkzKzCC3xC2xMws80k+wEYDSAPwAtmVkzySv/xgWY2\ng+QoANMAbAHwXHIBkJqxeXNpoP3gn3V4RiQXhDpJrJmNBDAyad3ApPYjAB4JM4dU7qnBwaEg3S5a\nJDfoimEBAAx7a2KgXbt2nqMkIlKTVAQE74z6NNC+77YLHCURkZqmIiB46Im3Au1jjuzgKImI1DQV\nAQko6Hqg7hYqkkNUBHLczNkLA+0/Xp18Zw8RyWYqAjmuzzVPBNqNGzVwlEREXFAREBHJYSoCOWzj\nxs2B9uAnrnGURERcURHIYa++OSHQbrt3M0dJRMQVFYEc9txLY1xHEBHH0hYBkoMTli+ukTRSY9at\n3xhoH3vUAY6SiIhLFe0JHJKwfF3YQaRmvfufzwJtXSUskps0HJSjHhs4omy5dp1Q7yMoIhFW0V9/\nS5KPAyCAFgnLAGBm9ofQ00koVq5aG2if+ptOjpKIiGsVFYE/wpsikgA+S3ose2Z4yUH/Gv5RoH3j\n73s6SiIirqUtAmY2uAZzSA0a9Nq4QLtWLY0KiuSqCv/6SfYh+TnJtf6/yTpTKN6Sp+k8/thD0vQU\nkVyQdk/Af7O/FsANAL6ANyzUEcDD/py//6iZiJJJk6d+F2jf8oczHSURkSioaE/gKgBnmdl4M1th\nZj+Z2TgAZwO4umbiSaa9/d4ngfYO9es6SiIiUVBREdjJzL5PXmlmcwDsFFoiCVXRh9PLlo88op3D\nJCISBRUVgfXVfExioucJh7uOICKOVXSKaHuSX6Z5bJ8wwki45i/6MdDudIh+jCK5rqIicDCApgBK\nkta3ArCwfHeJujfe/V+gvWOD+o6SiEhUVDQc9BiAlWY2J/EfgJUAHq2RdJJRK1auKVvOz6/nMImI\nREVFRaCpmZUbDjKzaQD2Ci+ShGX0uC/Kli/+vwJ3QUQkMioqAo0reEzjCDGzeXNpoN2ubQtHSUQk\nSioqApNJXpG8kuTlKH8vIYm4Dz7+OtDudPDejpKISJRUdGD4OgBvkbwAv7zpdwJQD4AuM42Z5IPC\nul+QiAAV30BuEcmjABwH4EB4dw59179qWGJmyvRfrvvreaKuDxART4WziZh3t7Fx/j+JqZKFywPt\nbr/q4CiJiESNxgSynJnh3Mv+GljXpVNbR2lEJGpUBLLc2yM/Kbeudu08B0lEJIpUBLLcI0++E2hP\nGHGPoyQiEkUqAjnk1BMO11lBIhKgd4QstmjJikD7qr49HCURkahSEchiTw8eHWg3apjvKImIRFWo\nRYBkD5IzSM4keXMF/Y4guZnkWWHmyTVj359atly3bh2HSUQkqkIrAiTzAAwA0ANABwC9SLZP0+9B\nAKPgzWMsGbBu/cZA++ZrznCURESiLMw9gc4AZvm3oN4EYCiA01P0uwbA6wCWhpgl5/z1qeGBdo9f\nd3SURESiLMwi0ALAvIR2ib+uDMkW8ArD0/4qCzFPTvno029cRxCRGAizCFTlDf0xALf4t6cgNByU\nMStX/TKBTJ/zjnOYRESirMJ7B22n+fCmotyqFcpPVdkJwFCSANAEwEkkN5nZ8KR+KCwsLFsuKChA\nQUFBhuNmr26/KncoRkSyUFFREYqKirbpa+h9CM88krUBfAOgO4AFAD4B0MvMitP0HwRghJm9meIx\nCyunC8XFv3wL2rfP/Bt00YfTcdt9r5a1x715F+rV09lBicL+GYhEAUmYWYUjLKHtCZjZZpL9AIwG\nkAfgBTMrJnml//jAsJ471yUWAAAqACKSVpjDQTCzkQBGJq1L+eZvZn3DzJIr1iedGnqizgoSkQro\niuEs8+6Y4Myft19/tqMkIhIHKgJZ5tFnRgTaumGciFRE7xBZTKeGikhlVASyyA/zghddX3J+d0dJ\nRCQuVASySOEjwwLtvDz9eEWkYnqXyBJmhm9nLShr16tX12EaEYkLFYEsMfaDaYH2gPsvdZREROJE\nRSBLPD0oOIFMh3at0vQUEfmFikCWWLz0l6kkzzyli8MkIhInKgJZYPlPqwPt3552pKMkIhI3KgJZ\n4PNpswPtPVvu5iiJiMSNikAW+N/kb11HEJGYUhHIAqPGfVG2fMD+rR0mEZG4URHIMqee0Ml1BBGJ\nERWBmEs+KHxiwaGOkohIHKkIxFzyhPKaQEZEtoWKQMz9O2H+gMaNGjhMIiJxpCIQc013a1y23P2Y\ngx0mEZE4UhGIubHvTy1b7txxX4dJRCSOVASySMOd8l1HEJGYURGIsS1btgTaLZrt4iiJiMSVikCM\njZs4PdDedeedHCURkbhSEYix5DkERES2lYpAjE34+Ouy5SOPaOcwiYjElYpATP244udA+4yTOjtK\nIiJxpiIQU6+P+DjQ7tp5f0dJRCTOVARiasnSlYE2SUdJRCTOVARiauR/Py9b1lCQiFSXikAMrVm7\nPtDeb98WjpKISNypCMTQ04P/E2j31BwCIlJNKgIxNHz0p4F2rVr6MYpI9ejdI4ZKN5eWLV/422Md\nJhGRuFMRiJnxSbeKOPvULo6SiEg2UBGImdvvfzXQ3r1JI0dJRCQbqAjEiJkF2r8++iBHSUQkW6gI\nxEjiVJIAcNNVpztKIiLZIvQiQLIHyRkkZ5K8OcXjF5CcSnIayQ9Jao7ENB4c8Hag3aihJpERke0T\nahEgmQdgAIAeADoA6EWyfVK32QCOMbODAfwFwLNhZoqzLaW/TCKju4aKSCaEvSfQGcAsM5tjZpsA\nDAUQGMMws4/NbOuNcCYBaBlypli686FhgfYNv+vpKImIZJOwi0ALAPMS2iX+unQuBfBeqIliKnFC\neQBovoemkhSR7Vc75O1b5V08JI8DcAmArqkeLywsLFsuKChAQUHBdkaLj02bNgfa/S49yVESEYmy\noqIiFBUVbdPXhF0E5gNoldBuBW9vIMA/GPwcgB5m9lOqDSUWgVwz5as5gfZ5Z3ZzE0REIi35A/Jd\nd91V6deEPRw0GUBbkm1I1gVwLoDhiR1ItgbwJoALzWxWyHli6eEn3wm0NXeAiGRKqHsCZraZZD8A\nowHkAXjBzIpJXuk/PhDAHQB2BvC0/+a2ycx0g/wE8xcsL1tuu09zh0lEJNuEPRwEMxsJYGTSuoEJ\ny5cBuCzsHNnigrOPdh1BRLKIrhiOuHXrNwbaR+n6ABHJIBWBiJs8JXiYpEF+fUdJRCQbqQhE3C1/\nedl1BBHJYioCETZ3/rJAu2WLJo6SiEi2UhGIsBvvGBxov/jYVW6CiEjWUhGIsAWLfgy0dTxARDJN\nRSCivv1uQaD9xP06i1ZEMk9FIKKuvGlgoH3YwXs7SiIi2UxFIILMDBs3biprN2igYSARCYeKQAQV\nPvzPQHvw4/0cJRGRbKciEEGaO0BEaoqKQMRs3BicO+DVgdc7SiIiuUBFIGKOO/OOQHvPlrs5SiIi\nuUBFIELMqjwRm4hIRqgIRMirb0wItJ/72+8dJRGRXKEiECFPDRoVaHdo1ypNTxGRzFARiIj5SbeI\nuOzC4x0lEZFcoiIQEUUTpwfafc47zlESEcklKgIR8fwr/y1bbrJrQ00mLyI1QkUgAlatXhu4TcQJ\nBYc6TCMiuURFIAJOOu+eQPvCc45xlEREco2KgGOz5ywut65Rw3wHSUQkF6kIONb76r8H2m8M+pOj\nJCKSi1QEHFqybFW5dXvs3thBEhHJVSoCjmwu3YLrC4cF1k18915HaUQkV6kIOPLggJHl1um0UBGp\naSoCjnw9MziH8Pvv/MVREhHJZSoCDlzQ77lA++/3XYratfMcpRGRXKYiUIPMDKf1vr/c+k6aRF5E\nHKntOkCuWL9+I7qfXVhu/eAnrtGxABFxRnsCNSRVATjjxI5ou3ezmg8jIuLTnkDIli5fhTMueqDc\n+kvO64bu3do7SCQi8gsVgRBNnDQDN9/9j3LrB/2tL+rW1bdeRNzTcFCIUhWA1569QQVARCJDRSAE\npaVb0PWU/uXWF719N1q3aOIgkYhIaqEWAZI9SM4gOZPkzWn6PO4/PpVkxzDz1IQfSpbimNNuL7f+\ng+H3oE4d7QGISLSEVgRI5gEYAKAHgA4AepFsn9TnZAD7mllbAFcAeDqsPGExMxR/W4LHn38PB3a5\nAOdf+Wi5Pg/deRHy8qKz01VUVOQ6QpXEIWccMgLKmWlxyVkVYb4zdQYwy8zmmNkmAEMBnJ7U5zQA\nQwDAzCYBaEyyaYiZttvyn1Zj7vxlGD1+Crqe0h/dTr0Nl13/FIa9NRErl88t1//tf9yCrp33d5A0\nvbj8AschZxwyAsqZaXHJWRVhjk+0ADAvoV0C4FdV6NMSQPmZViLgvCsfxbySpVXqW7duHfz71f7I\n36FeyKlERKovzCJgVeyXfLlsVb+uxlVWAJru1gi/73Miep11dKSGf0RE0qFZOO+5JLsAKDSzHn77\nVgBbzOzBhD7PACgys6F+ewaAY81scdK2IlsYRESizMwqvC9NmHsCkwG0JdkGwAIA5wLoldRnOIB+\nAIb6RWNFcgEAKn8RIiJSPaEVATPbTLIfgNEA8gC8YGbFJK/0Hx9oZu+RPJnkLABrAPQNK4+IiJQX\n2nCQiIhEX6SPXlblYrMoIPkiycUkv3SdJR2SrUiOJ/kVyekk/+A6Uyok65OcRHIKya9Jlp+AIUJI\n5pH8guQI11nSITmH5DQ/5yeu86RDsjHJ10kW+z/7Lq4zJSPZzv8+bv23Mop/SyRv9f/WvyT5Ksm0\npylGdk/Av9jsGwDHA5gP4FMAvcys2GmwFEgeDeBnAP8ws4Nc50mF5B4A9jCzKSR3BPAZgDMi+v3M\nN7O1JGsDmAjgJjOb6DpXKiRvANAJwE5mdprrPKmQ/B5AJzP70XWWipAcAuB9M3vR/9k3MLOVrnOl\nQ7IWvPemzmY2r7L+NcU/DjsOQHsz20ByGID3zGxIqv5R3hOoysVmkWBmEwD85DpHRcxskZlN8Zd/\nBlAMoLkUnr91AAAE10lEQVTbVKmZ2Vp/sS6840mRfPMi2RLAyQCeR/lTnaMm0vlINgJwtJm9CHjH\nFKNcAHzHA/guSgXAtwrAJgD5fjHNh1esUopyEUh1IVkLR1myiv9JoSOASW6TpEayFskp8C4aHG9m\nX7vOlMajAP4IYIvrIJUwAGNJTiZ5ueswaewFYCnJQSQ/J/kcyXzXoSpxHoBXXYdI5u/x/RXAXHhn\nZq4ws7Hp+ke5CERznCrm/KGg1wFc6+8RRI6ZbTGzQ+FdPX4MyQLHkcoheSqAJWb2BSL+KRtAVzPr\nCOAkAFf7w5dRUxvAYQCeMrPD4J0teIvbSOmRrAugJ4B/uc6SjOQ+AK4D0Abe3v6OJC9I1z/KRWA+\ngFYJ7Vbw9gakmkjWAfAGgJfN7G3XeSrjDwf8G8DhrrOkcBSA0/zx9tcA/Jpk+QkkIsDMFvr/LwXw\nFryh1qgpAVBiZp/67dfhFYWoOgnAZ/73NGoOB/CRmS03s80A3oT3+5pSlItA2cVmftU9F97FZVIN\n9GazfwHA12b2mOs86ZBsQrKxv7wDgN8A+MJtqvLMrL+ZtTKzveANC4wzs4tc50pGMp/kTv5yAwAn\nAIjcWWxmtgjAPJL7+auOB/CVw0iV6QWv+EfRDABdSO7g/90fDyDtkGpkb3Cf7mIzx7FSIvkagGMB\n7EpyHoA7zGyQ41jJugK4EMA0klvfVG81s1EOM6XSDMAQ/8yLWgBeMrP/Os5UFVEdvmwK4C3vvQC1\nAbxiZv9xGymtawC84n/o+w4RvXjUL6bHA4jk8RUzm+rvlU6Gd7zqcwDPpusf2VNERUQkfFEeDhIR\nkZCpCIiI5DAVARGRHKYiICKSw1QERERymIqAiEgOUxEQSYHkxSSbJbSfI9k+A9vdk2TyDHsizqgI\niKTWBwl3WTWzyzN0seJeAM7PwHZEMkJFQHIKyQv9SWu+IPmMf8fSwf7kG9NIXkfybHj3X3nFv6Nl\nfZJFJA/zt/EzyYf8yXnGkOxC8n2S35Hs6fdpQ/IDkp/5/470IzwA4Gj/+a/1n/9hkp+QnEryCjff\nGclVumJYcoY/nPMggDPNrJTkkwCWwLvL5gl+n4ZmtorkeAA3mtnn/vqyNsktAE4ys9Ek3wSwI7wb\nih0AYIiZdfTve7TFn9SjLYBXzewIksfCmyRna7G4AsBuZnavP/vTRAC/NbM5NfedkVwW2XsHiYSg\nO7xZwCb799KpD+/eVHuTfBzeHUsT76uT7hbRG81stL/8JYD1flGZDu/2vYA3Ic4AkocAKAXQNs02\nTwBwEMlz/HZDAPsCmLPNr06kGlQEJNcMMbP+iStI9gfQA8DvAPwfgEv9h9LtJm9KWN4CYCPgzYPg\nz+QEANcDWGhmvf2pUtdXkKmfmY3Ztpchkhk6JiC55L8AziG5GwCQ3IVkawC1zexNAH+GN+MaAKyG\n96m8uhoCWOQvXwTvTrhbt7tTQr/RAK7aWjxI7heDGbUki2hPQHKGmRWTvB3Af/xbVW8EcCOAR/02\n8MtsVoMBPENyLcpPyJG8h2Aplp8C8AbJiwCMArB1FrepAEr96TMHAXgc3hDS5/6935cAOLPaL1Jk\nG+nAsIhIDtNwkIhIDlMREBHJYSoCIiI5TEVARCSHqQiIiOQwFQERkRymIiAiksNUBEREctj/A3ny\nwzyHx4XiAAAAAElFTkSuQmCC\n",
"text/plain": [
"<matplotlib.figure.Figure at 0xf21afd0>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"<matplotlib.figure.Figure at 0x7adeb50>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"\"\"\"Sampling distribution of L as an estimator of exponential parameter.\n",
"\n",
"lam: parameter of an exponential distribution\n",
"n: sample size\n",
"m: number of iterations\n",
"\"\"\"\n",
"lam=2\n",
"n=10\n",
"m=1000\n",
"\n",
"def VertLine(x, y=1):\n",
" thinkplot.Plot([x, x], [0, y], color='0.8', linewidth=3)\n",
"\n",
"estimates = []\n",
"for j in range(m):\n",
" xs = np.random.exponential(1.0/lam, n)\n",
" lamhat = 1.0 / np.mean(xs)\n",
" estimates.append(lamhat)\n",
"\n",
"stderr = RMSE(estimates, lam)\n",
"print('standard error', stderr)\n",
"\n",
"cdf = thinkstats2.Cdf(estimates)\n",
"ci = cdf.Percentile(5), cdf.Percentile(95)\n",
"print('confidence interval', ci)\n",
"VertLine(ci[0])\n",
"VertLine(ci[1])\n",
"\n",
"# plot the CDF\n",
"thinkplot.Cdf(cdf)\n",
"thinkplot.Show(xlabel='estimate',\n",
" ylabel='CDF',\n",
" title='Sampling distribution')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"1. 표본크기 10인경우\n",
" - 표준오차 : 0.896717911545\n",
" - 신뢰구간: (1.2901330772324622, 3.8692334892427911)\n",
"1. 표본크기가 증가함에 따라, 표준오차와 CI 구간범위가 줄어든다:\n",
" - 10 0.90 (1.3, 3.9)\n",
" - 100 0.21 (1.7, 2.4)\n",
" - 1000 0.06 (1.9, 2.1)\n",
"상기 신뢰구간 세개 모두 실제 값 2를 포함하고 있다."
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": false
},
"source": [
"## 연습문제 8.3 \n",
"하키와 축구같은 스포츠 게임에서 득점 사이 시간은 대체로 지수를 따른다. 그래서 게임에서 한 팀이 득점한 골을 관측함으로써 득점을 추정할 수 있다. 이 추정 과정은 득점 사이 시간을 표집하는 것과 약간 다르다. 그래서 작동방법을 살펴보다.\n",
"\n",
"게임당 골로 득점률 ${\\tt lam}$을 인자로 받고, 전체 시간이 1 게임 경과할 때까지 득점사이 시간을 생성함으로서 게임을 모의시험하고 나서, 득점한 점수를 반환하는 함수를 작성하라.\n",
"\n",
"많은 게임을 모의시험하고, ${\\tt lam}$ 추정값을 저장하고 나서 평균 오차와 RMSE를 계산하는 또다른 함수를 작성하라.\n",
"\n",
"추정값을 이와 같은 방식으로 만드는 것이 편향됐을까? 추정값에 대한 표본분포와 90% 신뢰구간을 도식화하시오. 표준오차는 얼마인가? ${\\tt lam}$ 값을 크게하면, 표집오차에 무슨일이 생길까?"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Experiment 4\n",
"rmse L 1.41495194265\n",
"mean error L -0.000607\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXsAAAEACAYAAABS29YJAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAER9JREFUeJzt3X+MZWV9x/H3p0NpC9SqiYUKYyC6KhipkmZLSmuvFs2G\nVNekfyDRaEUof3SVNqYFbSK7adJQ+9tupKug0Ra7TSyYbcJmgbY3oUqAVQpUd5UN2bq7/BBBqZiQ\nLuHbP+7Z7d3Zmbl3dmbn3uV5v5LNnPOc5zn3O7Mzn3nmuefcm6pCkvTi9hOTLkCSdPwZ9pLUAMNe\nkhpg2EtSAwx7SWqAYS9JDRgZ9knWJdmd5OEk18xzfH2SB5Lcn+TrSd427lhJ0urIYtfZJ5kBvg1c\nDBwA7gMuq6pdQ31Oraofd9tvBG6tqteMM1aStDpGzezXAnuqam9VHQS2AuuHOxwK+s5pwPfHHStJ\nWh2jwv5MYN/Q/v6u7QhJ3p1kF7Ad+MhSxkqSjr9RYT/WaylU1Veq6lzgncDfJ8myK5MkrZiTRhw/\nAMwO7c8ymKHPq6ruSnIS8PKu38ixSXxxHkk6BlU19sR61Mx+J7AmydlJTgYuBbYNd0jy6kMz+SQX\ndAU8Nc7YoYKn/t9111038Rqs0zpP1Bqtc+X/LdWiM/uqej7JBmAHMAPcVFW7klzVHd8C/Bbw/iQH\ngWeB9yw2dskVSpKWbdQyDlW1ncETr8NtW4a2Pwl8ctyxkqTV5x20Y+r1epMuYSzWubJOhDpPhBrB\nOidt0ZuqVqWApCZdgySdaJJQK/gErSTpRcCwl6QGGPaS1ADDXpIaYNhLUgMMe0lqgGEvSQ0w7CWp\nAYa9JDXAsJekBhj2ktQAw16SGmDYS1IDDHtJaoBhL0kNGPlOVTo2G669cUXOs/n6K1b0fMPnlNQO\nZ/aS1ADDXpIaYNhLUgMMe0lqgGEvSQ0w7CWpAYa9JDXAsJekBhj2ktQAw16SGjAy7JOsS7I7ycNJ\nrpnn+HuTPJDkwSRfTXL+0LG9Xfv9Se5d6eIlSeNZ9LVxkswAm4GLgQPAfUm2VdWuoW6PAG+pqmeS\nrAM+A1zYHSugV1VPr3zpkqRxjZrZrwX2VNXeqjoIbAXWD3eoqrur6plu9x7grDnnyIpUKkk6ZqPC\n/kxg39D+/q5tIR8CbhvaL+DOJDuTXHlsJUqSlmvUSxzXuCdK8lbgcuCioeaLquqxJK8A7kiyu6ru\nOoY6JUnLMCrsDwCzQ/uzDGb3R+ielP0ssK6qfnCovaoe6z4+meRWBstCR4X9xo0bD2/3ej16vd7Y\nn4AktaDf79Pv9495/Kiw3wmsSXI28ChwKXDZcIckrwJuAd5XVXuG2k8BZqrqR0lOBd4BbJrvQYbD\nXpJ0tLkT4U2b5o3TBS0a9lX1fJINwA5gBripqnYluao7vgX4BPAy4IYkAAerai1wBnBL13YScHNV\n3b6k6iRJK2Lk2xJW1XZg+5y2LUPbVwBHvc9dVT0CvGkFapQkLZN30EpSAwx7SWqAYS9JDTDsJakB\nhr0kNcCwl6QGGPaS1ADDXpIaYNhLUgMMe0lqgGEvSQ0w7CWpAYa9JDXAsJekBhj2ktQAw16SGmDY\nS1IDDHtJaoBhL0kNMOwlqQGGvSQ1wLCXpAYY9pLUAMNekhpg2EtSAwx7SWqAYS9JDTDsJakBhr0k\nNWBk2CdZl2R3koeTXDPP8fcmeSDJg0m+muT8ccdKklbHomGfZAbYDKwDzgMuS3LunG6PAG+pqvOB\nPwY+s4SxkqRVMGpmvxbYU1V7q+ogsBVYP9yhqu6uqme63XuAs8YdK0laHaPC/kxg39D+/q5tIR8C\nbjvGsZKk4+SkEcdr3BMleStwOXDRUsdu3Ljx8Hav16PX6407VJKa0O/36ff7xzx+VNgfAGaH9mcZ\nzNCP0D0p+1lgXVX9YClj4ciwlyQdbe5EeNOmTUsaP2oZZyewJsnZSU4GLgW2DXdI8irgFuB9VbVn\nKWMlSatj0Zl9VT2fZAOwA5gBbqqqXUmu6o5vAT4BvAy4IQnAwapau9DY4/i5SJIWMGoZh6raDmyf\n07ZlaPsK4Ipxx0qSVp930EpSAwx7SWqAYS9JDTDsJakBhr0kNcCwl6QGGPaS1ADDXpIaYNhLUgMM\ne0lqwMiXS9CL14Zrb1yxc22+ft5XzJA0JZzZS1IDDHtJaoBhL0kNMOwlqQGGvSQ1wLCXpAYY9pLU\nAMNekhpg2EtSAwx7SWqAYS9JDTDsJakBhr0kNcCwl6QGGPaS1ADDXpIaYNhLUgNGhn2SdUl2J3k4\nyTXzHH99kruTPJfko3OO7U3yYJL7k9y7koVLksa36NsSJpkBNgMXAweA+5Jsq6pdQ92eAj4MvHue\nUxTQq6qnV6heSdIxGDWzXwvsqaq9VXUQ2AqsH+5QVU9W1U7g4ALnyPLLlCQtx6iwPxPYN7S/v2sb\nVwF3JtmZ5MqlFidJWhmLLuMwCOvluKiqHkvyCuCOJLur6q65nTZu3Hh4u9fr0ev1lvmwkvTi0u/3\n6ff7xzx+VNgfAGaH9mcZzO7HUlWPdR+fTHIrg2WhRcNeknS0uRPhTZs2LWn8qGWcncCaJGcnORm4\nFNi2QN8j1uaTnJLkZ7vtU4F3AA8tqTpJ0opYdGZfVc8n2QDsAGaAm6pqV5KruuNbkpwB3Ae8BHgh\nydXAecDPA7ckOfQ4N1fV7cfvU5EkLWTUMg5VtR3YPqdty9D24xy51HPIs8CbllugJGn5vINWkhpg\n2EtSAwx7SWqAYS9JDTDsJakBhr0kNWDkpZct2HDtjStyns3XX7Ei55GklebMXpIaYNhLUgMMe0lq\ngGEvSQ0w7CWpAYa9JDXAsJekBhj2ktQAw16SGmDYS1IDDHtJaoBhL0kNMOwlqQGGvSQ1wLCXpAYY\n9pLUAMNekhpg2EtSAwx7SWqAYS9JDTDsJakBI8M+yboku5M8nOSaeY6/PsndSZ5L8tGljJUkrY5F\nwz7JDLAZWAecB1yW5Nw53Z4CPgz8+TGMlSStglEz+7XAnqraW1UHga3A+uEOVfVkVe0EDi51rCRp\ndYwK+zOBfUP7+7u2cSxnrCRpBZ004ngt49xjj924cePh7V6vR6/XW8bDStKLT7/fp9/vH/P4UWF/\nAJgd2p9lMEMfx9hjh8NeknS0uRPhTZs2LWn8qGWcncCaJGcnORm4FNi2QN8sY6wk6ThadGZfVc8n\n2QDsAGaAm6pqV5KruuNbkpwB3Ae8BHghydXAeVX17Hxjj+cnI0ma36hlHKpqO7B9TtuWoe3HOXK5\nZtGxkqTV5x20ktQAw16SGjByGUdaig3X3rhi59p8/RUrdi6pdc7sJakBhr0kNcCwl6QGGPaS1ADD\nXpIaYNhLUgMMe0lqgGEvSQ0w7CWpAYa9JDXAsJekBhj2ktQAw16SGmDYS1IDDHtJaoBhL0kNMOwl\nqQGGvSQ1wLCXpAYY9pLUAMNekhpg2EtSAwx7SWqAYS9JDTDsJakBI8M+yboku5M8nOSaBfp8qjv+\nQJI3D7XvTfJgkvuT3LuShUuSxnfSYgeTzACbgYuBA8B9SbZV1a6hPpcAr6mqNUl+GbgBuLA7XECv\nqp4+LtVLksYyama/FthTVXur6iCwFVg/p8+7gC8AVNU9wEuTnD50PCtVrCTp2IwK+zOBfUP7+7u2\ncfsUcGeSnUmuXE6hkqRjt+gyDoOwHsdCs/dfrapHk7wCuCPJ7qq6a/zyJEkrYVTYHwBmh/ZnGczc\nF+tzVtdGVT3afXwyya0MloWOCvuNGzce3u71evR6vbGKl6RW9Pt9+v3+MY8fFfY7gTVJzgYeBS4F\nLpvTZxuwAdia5ELgh1X1RJJTgJmq+lGSU4F3AJvme5DhsJckHW3uRHjTpnnjdEGLhn1VPZ9kA7AD\nmAFuqqpdSa7qjm+pqtuSXJJkD/Bj4IPd8DOAW5Icepybq+r2JVUnSVoRo2b2VNV2YPucti1z9jfM\nM+4R4E3LLVCStHzeQStJDTDsJakBhr0kNcCwl6QGGPaS1ADDXpIaMPLSS2nSNlx744qda/P1V6zY\nuaQTiTN7SWqAYS9JDTDsJakBhr0kNcCwl6QGGPaS1ADDXpIaYNhLUgMMe0lqgGEvSQ0w7CWpAYa9\nJDXAsJekBhj2ktQAw16SGmDYS1IDDHtJaoDvVKUm+e5Xao0ze0lqgGEvSQ0w7CWpASPDPsm6JLuT\nPJzkmgX6fKo7/kCSNy9lrCTp+Fs07JPMAJuBdcB5wGVJzp3T5xLgNVW1Bvgd4IZxx55IDnz325Mu\nYSzWubJOhDr7/f6kSxiLdU7WqKtx1gJ7qmovQJKtwHpg11CfdwFfAKiqe5K8NMkZwDljjD1hHPju\ntznzVa+bdBkjWefKWkqdk7rCp9/v0+v1VuyxjxfrnKxRyzhnAvuG9vd3beP0eeUYYyVJq2DUzL7G\nPE+WW4iko43z18K9//ENvv/c6H7eD9C2VC2c50kuBDZW1bpu/2PAC1X1p0N9/g7oV9XWbn838OsM\nlnEWHdu1j/sLRZI0pKrGnmiPmtnvBNYkORt4FLgUuGxOn23ABmBr98vhh1X1RJKnxhi7pGIlScdm\n0bCvqueTbAB2ADPATVW1K8lV3fEtVXVbkkuS7AF+DHxwsbHH85ORJM1v0WUcSdKLw1TcQZvkz5Ls\n6m7KuiXJz026pmEnws1hSWaT/HuSbyb5ryQfmXRNC0kyk+T+JP8y6VoW0l1C/OXu+/Jb3RLl1Eny\nse7//KEkX0ryU5OuCSDJ55I8keShobaXJ7kjyXeS3J7kpZOssatpvjqnLo/mq3Po2EeTvJDk5Yud\nYyrCHrgdeENV/SLwHeBjE67nsBPo5rCDwO9X1RuAC4HfndI6Aa4GvsX4V3tNwt8At1XVucD5TOH9\nId3zYVcCF1TVGxksl75nkjUN+TyDn5lh1wJ3VNVrgX/t9idtvjqnMY/mq5Mks8Dbgf8edYKpCPuq\nuqOqXuh27wHOmmQ9cxy+sayqDgKHbg6bKlX1eFX9Z7f9LINweuVkqzpakrOAS4AbmdJLdruZ3K9V\n1edg8PxTVT0z4bLm8z8MfsmfkuQk4BTgwGRLGqiqu4AfzGk+fANm9/Hdq1rUPOarcxrzaIGvJ8Bf\nAn84zjmmIuznuBy4bdJFDBnnxrKp0s343szgG3Xa/BXwB8ALozpO0DnAk0k+n+QbST6b5JRJFzVX\nVT0N/AXwXQZXvP2wqu6cbFWLOr2qnui2nwBOn2QxY5q2PDosyXpgf1U9OE7/VQv7bq3uoXn+vXOo\nzx8B/1tVX1qtusYwzUsNR0lyGvBl4Opuhj81kvwm8L2qup8pndV3TgIuAD5dVRcwuMpsGpYcjpDk\n1cDvAWcz+CvutCTvnWhRY6rBlSFT/bM1pXkEQDf5+Dhw3XDzYmNW7Z2qqurtix1P8tsM/rz/jVUp\naHwHgNmh/VkGs/upk+QngX8G/qGqvjLpeubxK8C7uhfP+2ngJUm+WFXvn3Bdc+1nMGO6r9v/MlMY\n9sAvAV+rqqcAktzC4Gt880SrWtgTSc6oqseT/ALwvUkXtJApzqNDXs3gl/wDSWCw1PT1JGurat6v\n61Qs4yRZx+BP+/VV9dyk65nj8I1lSU5mcHPYtgnXdJQM/sdvAr5VVX896XrmU1Ufr6rZqjqHwROJ\n/zaFQU9VPQ7sS/Laruli4JsTLGkhu4ELk/xM9/9/MYMnvqfVNuAD3fYHgGmckEx7HgFQVQ9V1elV\ndU7387SfwRP1C/4CnYqwB/4WOA24o7sk79OTLuiQqnqewR3COxj8IP3TlN4cdhHwPuCt3dfw/u6b\ndppN85/xHwZuTvIAg6tx/mTC9Rylqh4AvshgQnJo3fYzk6vo/yX5R+BrwOuS7EvyQeB64O1JvgO8\nrdufqHnqvJwpzKOhOl879PUcNvJnyZuqJKkB0zKzlyQdR4a9JDXAsJekBhj2ktQAw16SGmDYS1ID\nDHtJaoBhL0kN+D9EoxgPDDPaIgAAAABJRU5ErkJggg==\n",
"text/plain": [
"<matplotlib.figure.Figure at 0xebbc870>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"<matplotlib.figure.Figure at 0xebbc910>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"def SimulateGame(lam):\n",
" \"\"\"Simulates a game and returns the estimated goal-scoring rate.\n",
"\n",
" lam: actual goal scoring rate in goals per game\n",
" \"\"\"\n",
" goals = 0\n",
" t = 0\n",
" while True:\n",
" time_between_goals = random.expovariate(lam)\n",
" t += time_between_goals\n",
" if t > 1:\n",
" break\n",
" goals += 1\n",
"\n",
" # estimated goal-scoring rate is the actual number of goals scored\n",
" L = goals\n",
" return L\n",
"\n",
"def Estimate4(lam=2, m=1000000):\n",
"\n",
" estimates = []\n",
" for i in range(m):\n",
" L = SimulateGame(lam)\n",
" estimates.append(L)\n",
"\n",
" print('Experiment 4')\n",
" print('rmse L', RMSE(estimates, lam))\n",
" print('mean error L', MeanError(estimates, lam))\n",
" \n",
" pmf = thinkstats2.Pmf(estimates)\n",
"\n",
" thinkplot.Hist(pmf)\n",
" thinkplot.Show()\n",
"\n",
"Estimate4() "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"1. 이런 방식으로 lambda 추정에 대한 RMSE는 1.4.\n",
"1. 평균오차는 작고, m이 커짐에 따라 줄어든다. 그래서, 추정량은 편향되지 않은 것으로 보인다.\n",
"주의: 득점 사이 시간이 지수분포라면, 게임에서 득점된 골의 분포는 포아송이다. [https://en.wikipedia.org/wiki/Poisson_distribution](https://en.wikipedia.org/wiki/Poisson_distribution) 참조한다."
]
},
{
"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
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment