mirror of
https://github.com/ION606/COGMOD-HWI.git
synced 2026-05-14 22:16:57 +00:00
150 lines
29 KiB
Plaintext
150 lines
29 KiB
Plaintext
{
|
||
"cells": [
|
||
{
|
||
"cell_type": "markdown",
|
||
"metadata": {},
|
||
"source": [
|
||
"# Problem 4\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 12,
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"import numpy as np\n",
|
||
"import matplotlib.pyplot as plt\n",
|
||
"from scipy.stats import norm\n",
|
||
"\n",
|
||
"import matplotlib.pyplot as plt\n",
|
||
"from scipy.stats import norm\n",
|
||
"import pandas as pd\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 17,
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"Sample Mean (ybar): 4.66 hours\n",
|
||
"Posterior Mean (µ_n): 4.70\n",
|
||
"Posterior Standard Deviation (σ_n): 0.21\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"np.random.seed(42)\n",
|
||
"n = 50 #sample size\n",
|
||
"mu_true = 5 #avg sleep time\n",
|
||
"sigma = 1.5 #standard deviation\n",
|
||
"data = np.random.normal(mu_true, sigma, n)\n",
|
||
"ybar = np.mean(data) \n",
|
||
"print(f\"Sample Mean (ybar): {ybar:.2f} hours\")\n",
|
||
"\n",
|
||
"# Prior parameters (initial belief about sleep time)\n",
|
||
"mu0 = 5.5\n",
|
||
"sigma0 = 1.0\n",
|
||
"\n",
|
||
"# Posterior parameters\n",
|
||
"sigma_n = np.sqrt(1 / (1 / sigma0**2 + n / sigma**2))\n",
|
||
"mu_n = (mu0 / sigma0**2 + n * ybar / sigma**2) / (1 / sigma0**2 + n / sigma**2)\n",
|
||
"\n",
|
||
"print(f\"Posterior Mean (µ_n): {mu_n:.2f}\")\n",
|
||
"print(f\"Posterior Standard Deviation (σ_n): {sigma_n:.2f}\")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 18,
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAzoAAAIjCAYAAADLH25TAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjEsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvc2/+5QAAAAlwSFlzAAAPYQAAD2EBqD+naQAATBhJREFUeJzt3XlcFuX+//H3DbLJ6gIiiohg7opbhuZWGq5pWprZwb1OkkuWJZ7cUiNL006p6am0XH62uhzTzH3LOm64Zy64VO4KCCoozO8Pv951CygoeMPwej4e9+PBzFwz85l7LHl7zXWNxTAMQwAAAABgIg72LgAAAAAAchtBBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAKoWbNmqlZs2bW5WPHjslisWj27Nl2q6mgKF++vHr27GnvMvLM6NGjZbFY7F0GANw3gg4AU5g9e7YsFou2bduW6fZmzZqpevXqD7iq3HPkyBG9+OKLqlChglxdXeXl5aVGjRrpgw8+0NWrV+1d3gNz7Ngx9erVSyEhIXJ1dZW/v7+aNGmiUaNG2bSbNm1agQ5ty5Yt0+jRo3P1mElJSRo1apSqV68ud3d3lShRQmFhYRo0aJD+/PPPXD0XAOQHRexdAADgzr7//ns988wzcnFxUWRkpKpXr67U1FRt2rRJQ4cO1b59+zRz5kx7l5nnDh8+rPr168vNzU29e/dW+fLlderUKe3YsUMTJkzQmDFjrG2nTZumkiVLFtiel2XLlmnq1Km5FnauX7+uJk2a6Ndff1WPHj00YMAAJSUlad++fZo/f76eeuopBQQE5Mq5ACC/IOgAQB5KTk6Wu7v7Pe8fFxenZ599VkFBQVqzZo1Kly5t3RYVFaXDhw/r+++/z41S873JkycrKSlJsbGxCgoKstl29uxZO1VVMCxatEg7d+7UvHnz9Nxzz9lsu3btmlJTU+1UGQDkHR5dA1Bo3bhxQ2PHjlVISIhcXFxUvnx5DR8+XCkpKTbtLBZLpv+yfvtYjVuPz61fv179+/eXn5+fypYtK0m6fPmyBg8erPLly8vFxUV+fn5q2bKlduzYccca3333XSUlJenTTz+1CTm3hIaGatCgQTm+puz69ddf9fTTT6t48eJydXVVvXr1tGTJkgztdu/eraZNm8rNzU1ly5bVuHHjNGvWLFksFh07dsym7fLly9W4cWO5u7vL09NTbdu21b59++5ay5EjR1S2bNkMIUeS/Pz8rD+XL19e+/bt0/r162WxWGSxWKzjkbIaf3Lr3v29VsMwNG7cOJUtW1ZFixZV8+bNs6wzPj5egwcPVmBgoFxcXBQaGqoJEyYoPT3d2ubWOKiJEydq5syZ1ntUv359bd261dquZ8+emjp1qiRZ6/97zQsWLFDdunXl6ekpLy8v1ahRQx988MFdvztJatSoUYZttx6FvJu5c+eqbt26cnNzU/HixfXss8/q5MmTGdr98ssvatWqlby9vVW0aFE1bdpUmzdvtmlz6z78+uuv6tKli7y8vFSiRAkNGjRI165du2stAJAd9OgAMJWEhASdP38+w/rr169nWNe3b199/vnnevrpp/Xqq6/ql19+UUxMjA4cOKCFCxfecw39+/eXr6+vRo4cqeTkZEnSP//5T33zzTd6+eWXVbVqVV24cEGbNm3SgQMHVKdOnSyP9d///lcVKlRQw4YNs3Xu3Lymffv2qVGjRipTpoyGDRsmd3d3ffXVV+rYsaO+/fZbPfXUU5KkP/74Q82bN5fFYlF0dLTc3d31ySefyMXFJcMx58yZox49eigiIkITJkzQlStXNH36dD366KPauXOnypcvn2U9QUFBWrVqldasWaPHHnssy3ZTpkzRgAED5OHhoX/961+SpFKlSuXo2iVp5MiRGjdunNq0aaM2bdpox44deuKJJzL0fly5ckVNmzbVH3/8oRdffFHlypXTTz/9pOjoaJ06dUpTpkyxaT9//nxdvnxZL774oiwWi95991116tRJR48elZOTk1588UX9+eefWrlypebMmWOz78qVK9WtWzc9/vjjmjBhgiTpwIED2rx5s03gvd2tcPjFF1/ozTffzPFkA+PHj9eIESPUpUsX9e3bV+fOndOHH36oJk2aaOfOnfLx8ZEkrVmzRq1bt1bdunU1atQoOTg4aNasWXrssce0ceNGPfzwwzbH7dKli8qXL6+YmBj9/PPP+ve//61Lly7piy++yFF9AJApAwBMYNasWYakO36qVatmbR8bG2tIMvr27WtznNdee82QZKxZs8a6TpIxatSoDOcMCgoyevTokaGGRx991Lhx44ZNW29vbyMqKipH15SQkGBIMjp06JCt9jm5pqZNmxpNmza1LsfFxRmSjFmzZlnXPf7440aNGjWMa9euWdelp6cbDRs2NCpWrGhdN2DAAMNisRg7d+60rrtw4YJRvHhxQ5IRFxdnGIZhXL582fDx8TH69etnU9/p06cNb2/vDOtvt3fvXsPNzc2QZISFhRmDBg0yFi1aZCQnJ2doW61aNZvru2XUqFFGZn/13bp3t2o9e/as4ezsbLRt29ZIT0+3ths+fLghyea+jx071nB3dzd+++03m2MOGzbMcHR0NE6cOGEYxl/fcYkSJYyLFy9a2y1evNiQZPz3v/+1rouKisq0zkGDBhleXl4Z/nzdzZUrV4xKlSoZkoygoCCjZ8+exqeffmqcOXMmQ9vbv6Njx44Zjo6Oxvjx423a7dmzxyhSpIh1fXp6ulGxYkUjIiLC5ju7cuWKERwcbLRs2TLDOZ588kmbY/bv39+QZOzatStH1wcAmeHRNQCmMnXqVK1cuTLDp2bNmjbtli1bJkkaMmSIzfpXX31Vku5r3Eu/fv3k6Ohos87Hx0e//PJLjma3SkxMlCR5enpmq31uXtPFixe1Zs0adenSRZcvX9b58+d1/vx5XbhwQRERETp06JD++OMPSdIPP/yg8PBwhYWFWfcvXry4unfvbnPMlStXKj4+Xt26dbMe7/z583J0dFSDBg20du3aO9ZUrVo1xcbG6vnnn9exY8f0wQcfqGPHjipVqpT+85//ZPvasmPVqlVKTU3VgAEDbHo/Bg8enKHt119/rcaNG6tYsWI219WiRQulpaVpw4YNNu27du2qYsWKWZcbN24sSTp69Ohd6/Lx8VFycrJWrlyZo+txc3PTL7/8oqFDh0q6+ahenz59VLp0aQ0YMOCOjzZ+9913Sk9PV5cuXWyuz9/fXxUrVrTet9jYWB06dEjPPfecLly4YG2XnJysxx9/XBs2bLB5lE+6Oc7s7wYMGCDprz/LAHA/eHQNgKk8/PDDqlevXob1t34JveX48eNycHBQaGioTTt/f3/5+Pjo+PHj91xDcHBwhnXvvvuuevToocDAQNWtW1dt2rRRZGSkKlSokOVxbo2buHz5crbOm5vXdPjwYRmGoREjRmjEiBGZtjl79qzKlCmj48ePKzw8PMP22+s4dOiQJGX52Fl2xok89NBDmjNnjtLS0rR//34tXbpU7777rl544QUFBwerRYsWdz1Gdtz6ripWrGiz3tfX1yakSDeva/fu3fL19c30WLdPlFCuXDmb5VvHu3Tp0l3r6t+/v7766iu1bt1aZcqU0RNPPKEuXbqoVatWd93X29tb7777rt59910dP35cq1ev1sSJE/XRRx/J29tb48aNy3S/Q4cOyTCMDN/FLU5OTtZ2ktSjR48sa0hISLD5/m4/ZkhIiBwcHDKM6wKAe0HQAVCo3c+LEdPS0jJd7+bmlmFdly5d1LhxYy1cuFA//vij3nvvPU2YMEHfffedWrdunelxvLy8FBAQoL179+aortx42eOtf3l/7bXXFBERkWmb24NMdo85Z84c+fv7Z9hepEj2/0pydHRUjRo1VKNGDYWHh6t58+aaN2/eXYNOVt9NVvcyO9LT09WyZUu9/vrrmW5/6KGHbJZv7+27xTCMu57Lz89PsbGxWrFihZYvX67ly5dr1qxZioyM1Oeff57tmoOCgtS7d2899dRTqlChgubNm5dl0ElPT5fFYtHy5cszrd3Dw8PaTpLee+89m969zNpmhReVAshNBB0AhVJQUJDS09N16NAhValSxbr+zJkzio+Pt5nZq1ixYoqPj7fZPzU1VadOncrROUuXLq3+/furf//+Onv2rOrUqaPx48dnGXQkqV27dpo5c6a2bNmSaa/JvV7T3dzqaXJycrpreAgKCtLhw4czrL99XUhIiKSbv6znVs+LJGsP3t/vR1a/MN/qTYiPj7cOoJeUobfr1nd16NAhm163c+fOZeh5CQkJUVJSUq5e051+4Xd2dlb79u3Vvn17paenq3///poxY4ZGjBiR4/BZrFgxhYSE3DFMh4SEyDAMBQcHZwhtt7eTbgb07H4Xhw4dsukBPXz4sNLT0+84KQUAZBdjdAAUSm3atJGkDDNivf/++5Kktm3bWteFhIRkGGcxc+bMbPcCpKWlKSEhwWadn5+fAgIC7jrt8+uvvy53d3f17dtXZ86cybD9yJEj1qmFc3JNd+Pn56dmzZppxowZmQa6c+fOWX+OiIjQli1bFBsba1138eJFzZs3z2afiIgIeXl56e233850Fry/HzMzGzduzHS/W+M5KlWqZF3n7u6eIZxKf/0y/vf7mZycnKE3pEWLFnJyctKHH35o09Ny+3cr3eyt27Jli1asWJFhW3x8vG7cuHHH68rMrXcv3X4NFy5csFl2cHCwjj+705+lXbt2ZTob4fHjx7V//36b7+52nTp1kqOjo8aMGZOh18kwDGtNdevWVUhIiCZOnKikpKQMx8ns/t6aRvuWDz/8UJLuGP4BILvo0QFQKNWqVUs9evTQzJkzFR8fr6ZNm+p///ufPv/8c3Xs2FHNmze3tu3bt6/++c9/qnPnzmrZsqV27dqlFStWqGTJktk61+XLl1W2bFk9/fTTqlWrljw8PLRq1Spt3bpVkyZNuuO+ISEhmj9/vrp27aoqVaooMjJS1atXV2pqqn766Sd9/fXX1nf55OSasmPq1Kl69NFHVaNGDfXr108VKlTQmTNntGXLFv3+++/atWuXpJthbO7cuWrZsqUGDBhgnV66XLlyunjxorV3wsvLS9OnT9c//vEP1alTR88++6x8fX114sQJff/992rUqJE++uijLOuZMGGCtm/frk6dOll/ud+xY4e++OILFS9e3GaigLp162r69OkaN26cQkND5efnp8cee0xPPPGEypUrpz59+mjo0KFydHTUZ599Zq3jFl9fX7322muKiYlRu3bt1KZNG+3cuVPLly/PcN+HDh2qJUuWqF27durZs6fq1q2r5ORk7dmzR998842OHTuW7T8rf69fkgYOHKiIiAg5Ojrq2WefVd++fXXx4kU99thjKlu2rI4fP64PP/xQYWFhNr14t1u5cqVGjRqlJ598Uo888og8PDx09OhRffbZZ0pJScn0PVG3hISEaNy4cYqOjtaxY8fUsWNHeXp6Ki4uTgsXLtQLL7yg1157TQ4ODvrkk0/UunVrVatWTb169VKZMmX0xx9/aO3atfLy8tJ///tfm2PHxcXpySefVKtWrbRlyxbNnTtXzz33nGrVqpWj7wsAMmXHGd8AINfcmh5469atmW5v2rSpzfTShmEY169fN8aMGWMEBwcbTk5ORmBgoBEdHW0znbJhGEZaWprxxhtvGCVLljSKFi1qREREGIcPH85yeunba0hJSTGGDh1q1KpVy/D09DTc3d2NWrVqGdOmTcv29f32229Gv379jPLlyxvOzs6Gp6en0ahRI+PDDz+0qTe715Sd6aUNwzCOHDliREZGGv7+/oaTk5NRpkwZo127dsY333xj027nzp1G48aNDRcXF6Ns2bJGTEyM8e9//9uQZJw+fdqm7dq1a42IiAjD29vbcHV1NUJCQoyePXsa27Ztu+N3sHnzZiMqKsqoXr264e3tbTg5ORnlypUzevbsaRw5csSm7enTp422bdsanp6ehiSba92+fbvRoEEDw9nZ2ShXrpzx/vvvZ5he2jBu3vcxY8YYpUuXNtzc3IxmzZoZe/fuzXDfDePm1NnR0dFGaGio4ezsbJQsWdJo2LChMXHiRCM1NdXmO37vvfcyXJtum8L8xo0bxoABAwxfX1/DYrFYp3v+5ptvjCeeeMLw8/Oz1v/iiy8ap06duuN3d/ToUWPkyJHGI488Yvj5+RlFihQxfH19jbZt29pMO24YWU/B/e233xqPPvqo4e7ubri7uxuVK1c2oqKijIMHD9q027lzp9GpUyejRIkShouLixEUFGR06dLFWL16dYZz7N+/33j66acNT09Po1ixYsbLL79sXL169Y7XAgDZZTGMbIx+BAAghwYPHqwZM2YoKSkpywH4KJxGjx6tMWPG6Ny5cznu7QKA7GKMDgDgvl29etVm+cKFC5ozZ44effRRQg4AwC4YowMAuG/h4eFq1qyZqlSpojNnzujTTz9VYmJilu/gAQAgrxF0AAD3rU2bNvrmm280c+ZMWSwW1alTR59++qmaNGli79IAAIUUY3QAAAAAmA5jdAAAAACYDkEHAAAAgOkUiDE66enp+vPPP+Xp6Wl98RwAAACAwscwDF2+fFkBAQFycMi636ZABJ0///xTgYGB9i4DAAAAQD5x8uRJlS1bNsvtBSLoeHp6Srp5MV5eXnauBgAAAIC9JCYmKjAw0JoRslIggs6tx9W8vLwIOgAAAADuOqSFyQgAAAAAmA5BBwAAAIDpEHQAAAAAmE6BGKMDAAAAZCYtLU3Xr1+3dxnIRY6OjipSpMh9v1aGoAMAAIACKSkpSb///rsMw7B3KchlRYsWVenSpeXs7HzPxyDoAAAAoMBJS0vT77//rqJFi8rX15eXypuEYRhKTU3VuXPnFBcXp4oVK97xpaB3QtABAABAgXP9+nUZhiFfX1+5ubnZuxzkIjc3Nzk5Oen48eNKTU2Vq6vrPR2HyQgAAABQYNGTY0732otjc4xcqAMAAAAA8hWCDgAAAADTYYwOAAAATGPyyt8e6PleaflQnp+jfPnyGjx4sAYPHpzn5zITenQAAACAB6Rnz56yWCyyWCxydnZWaGio3nrrLd24cSPLfbZu3aoXXnjhAVZpDvToAAAAAA9Qq1atNGvWLKWkpGjZsmWKioqSk5OToqOjbdqlpqbK2dlZvr6+93W+W8cpbOjRAQAAAB4gFxcX+fv7KygoSC+99JJatGihJUuWqGfPnurYsaPGjx+vgIAAVapUSdLNR9emTJli3f/EiRPq0KGDPDw85OXlpS5duujMmTPW7aNHj1ZYWJg++eQTBQcH3/P0zAUdPToAAACAHbm5uenChQuSpNWrV8vLy0srV67MtG16ero15Kxfv143btxQVFSUunbtqnXr1lnbHT58WN9++62+++47OTo6PojLyHcIOgAAAIAdGIah1atXa8WKFRowYIDOnTsnd3d3ffLJJ1k+arZ69Wrt2bNHcXFxCgwMlCR98cUXqlatmrZu3ar69etLuvm42hdffHHfj70VZDy6BgAAADxAS5culYeHh1xdXdW6dWt17dpVo0ePliTVqFHjjuNpDhw4oMDAQGvIkaSqVavKx8dHBw4csK4LCgoq1CFHokcHAAAAeKCaN2+u6dOny9nZWQEBASpS5K9fyd3d3XPlHLl1nIKMoAMAAAA8QO7u7goNDb2nfatUqaKTJ0/q5MmT1l6d/fv3Kz4+XlWrVs3NMgs8gg4Ac1kbY7vcPDrzdgAAFEAtWrRQjRo11L17d02ZMkU3btxQ//791bRpU9WrV8/e5eUrBB0AAACYxistH7J3CXnKYrFo8eLFGjBggJo0aSIHBwe1atVKH374ob1Ly3cshmEY9i7ibhITE+Xt7a2EhAR5eXnZuxwA+Rk9OgBQKFy7dk1xcXGF+j0xZnan+5vdbMCsawAAAABMh6ADAAAAwHQIOgAAAABMh6ADAAAAwHQIOgAAAABMh6ADAAAAwHQIOgAAAABMh6ADAAAAwHQIOgAAAABMp4i9CwAAAAByzdqYB3u+5tG5diiLxaKFCxeqY8eOOnbsmIKDg7Vz506FhYXl+vHWrVun5s2b69KlS/Lx8cm1a8iO0aNHa9GiRYqNjc3T8xB0AAAAgAekZ8+eio+P16JFizJsO3XqlIoVK5Yn5w0MDNSpU6dUsmTJPDl+fkTQAQAAAPIBf3//PDu2o6Njnh4/P2KMDgAAAJAPWCyWTHt6JCktLU29e/dW5cqVdeLECUnS4sWLVadOHbm6uqpChQoaM2aMbty4ken+x44dk8ViyfC42Pbt21WvXj0VLVpUDRs21MGDB222T58+XSEhIXJ2dlalSpU0Z84cm+0nTpxQhw4d5OHhIS8vL3Xp0kVnzpyxafPOO++oVKlS8vT0VJ8+fXTt2rUcfCv3jqADAAAA5GMpKSl65plnFBsbq40bN6pcuXLauHGjIiMjNWjQIO3fv18zZszQ7NmzNX78+Bwd+1//+pcmTZqkbdu2qUiRIurdu7d128KFCzVo0CC9+uqr2rt3r1588UX16tVLa9eulSSlp6erQ4cOunjxotavX6+VK1fq6NGj6tq1q/UYX331lUaPHq23335b27ZtU+nSpTVt2rTc+WLugkfXAAAAgHwqKSlJbdu2VUpKitauXStvb29J0pgxYzRs2DD16NFDklShQgWNHTtWr7/+ukaNGpXt448fP15NmzaVJA0bNkxt27bVtWvX5OrqqokTJ6pnz57q37+/JGnIkCH6+eefNXHiRDVv3lyrV6/Wnj17FBcXp8DAQEnSF198oWrVqmnr1q2qX7++pkyZoj59+qhPnz6SpHHjxmnVqlUPpFeHHh0AAAAgn+rWrZuSk5P1448/WkOOJO3atUtvvfWWPDw8rJ9+/frp1KlTunLlSraPX7NmTevPpUuXliSdPXtWknTgwAE1atTIpn2jRo104MAB6/bAwEBryJGkqlWrysfHx6ZNgwYNbI4RHh6e7fruBz06AAAAQD7Vpk0bzZ07V1u2bNFjjz1mXZ+UlKQxY8aoU6dOGfZxdXXN9vGdnJysP1ssFkk3H0kzA3p0AAAAgHzqpZde0jvvvKMnn3xS69evt66vU6eODh48qNDQ0AwfB4fc+RW/SpUq2rx5s826zZs3q2rVqtbtJ0+e1MmTJ63b9+/fr/j4eJs2v/zyi80xfv7551yp727o0QEAAAAeoISEhAyzn5UoUSLL9gMGDFBaWpratWun5cuX69FHH9XIkSPVrl07lStXTk8//bQcHBy0a9cu7d27V+PGjcuVOocOHaouXbqodu3aatGihf773//qu+++06pVqyRJLVq0UI0aNdS9e3dNmTJFN27cUP/+/dW0aVPVq1dPkjRo0CD17NlT9erVU6NGjTRv3jzt27dPFSpUyJUa74SgAwAAAPNoHm3vCu5q3bp1ql27ts26W4P1szJ48GClp6erTZs2+uGHHxQREaGlS5fqrbfe0oQJE+Tk5KTKlSurb9++uVZnx44d9cEHH2jixIkaNGiQgoODNWvWLDVr1kzSzUfdFi9erAEDBqhJkyZycHBQq1at9OGHH1qP0bVrVx05ckSvv/66rl27ps6dO+ull17SihUrcq3OrFgMwzDy/Cz3KTExUd7e3kpISJCXl5e9ywGQn62NsV0uAH/hAQBy7tq1a4qLi1NwcHCOxqSgYLjT/c1uNmCMDgAAAADTIegAAAAAMB2CDgAAAADTYTICAPg/02KnZVjXP6y/HSoBAAD3ix4dAAAAAKaTo6ATExOj+vXry9PTU35+furYsaMOHjx41/2+/vprVa5cWa6urqpRo4aWLVt2zwUDAAAAwN3kKOisX79eUVFR+vnnn7Vy5Updv35dTzzxhJKTk7Pc56efflK3bt3Up08f7dy5Ux07dlTHjh21d+/e+y4eAAAAADKTozE6P/zwg83y7Nmz5efnp+3bt6tJkyaZ7vPBBx+oVatWGjp0qCRp7NixWrlypT766CN9/PHH91g2AAAAAGTtviYjSEhIkCQVL148yzZbtmzRkCFDbNZFRERo0aJFWe6TkpKilJQU63JiYuL9lAmgMOMFogAAFEr3HHTS09M1ePBgNWrUSNWrV8+y3enTp1WqVCmbdaVKldLp06ez3CcmJkZjxoy519IAAABQSGU2g2ZeMuvsnKNHj9aiRYsUGxtr71Lu2T3PuhYVFaW9e/dqwYIFuVmPJCk6OloJCQnWz8mTJ3P9HAAAAMCD1rNnT1ksFlksFjk7Oys0NFRvvfWWbty4cV/HXbdunSwWi+Lj43Olztdee02rV6/OlWPZyz316Lz88staunSpNmzYoLJly96xrb+/v86cOWOz7syZM/L3989yHxcXF7m4uNxLaQAAAEC+1qpVK82aNUspKSlatmyZoqKi5OTkpOho+z9ebRiG0tLS5OHhIQ8Pj/s61vXr1+Xk5JRLleVcjnp0DMPQyy+/rIULF2rNmjUKDg6+6z7h4eEZ0uDKlSsVHh6es0oBAAAAE3BxcZG/v7+CgoL00ksvqUWLFlqyZIkuXbqkyMhIFStWTEWLFlXr1q116NAh637Hjx9X+/btVaxYMbm7u6tatWpatmyZjh07pubNm0uSihUrJovFop49e0q6OdwkJiZGwcHBcnNzU61atfTNN99Yj3mrJ2j58uWqW7euXFxctGnTJo0ePVphYWHWdunp6XrrrbdUtmxZubi4KCwszGaismPHjslisejLL79U06ZN5erqqnnz5uXtF3kXOerRiYqK0vz587V48WJ5enpax9l4e3vLzc1NkhQZGakyZcooJubmAOBBgwapadOmmjRpktq2basFCxZo27ZtmjlzZi5fCgAAAFDwuLm56cKFC+rZs6cOHTqkJUuWyMvLS2+88YbatGmj/fv3y8nJSVFRUUpNTdWGDRvk7u6u/fv3y8PDQ4GBgfr222/VuXNnHTx4UF5eXtbfzWNiYjR37lx9/PHHqlixojZs2KDnn39evr6+atq0qbWGYcOGaeLEiapQoYKKFSumdevW2dT4wQcfaNKkSZoxY4Zq166tzz77TE8++aT27dunihUr2hxn0qRJql27tlxdXR/I95eVHAWd6dOnS5KaNWtms37WrFnW1HjixAk5OPzVUdSwYUPNnz9fb775poYPH66KFStq0aJFd5zAAAAAADA7wzC0evVqrVixQq1bt9aiRYu0efNmNWzYUJI0b948BQYGatGiRXrmmWd04sQJde7cWTVq1JAkVahQwXqsW7Mg+/n5ycfHR9LNmYzffvttrVq1yvo0VYUKFbRp0ybNmDHDJui89dZbatmyZZa1Tpw4UW+88YaeffZZSdKECRO0du1aTZkyRVOnTrW2Gzx4sDp16pQL3879y1HQMQzjrm1uT3+S9Mwzz+iZZ57JyakAAAAAU1q6dKk8PDx0/fp1paen67nnnlOnTp20dOlSNWjQwNquRIkSqlSpkg4cOCBJGjhwoF566SX9+OOPatGihTp37qyaNWtmeZ7Dhw/rypUrGQJMamqqateubbOuXr16WR4nMTFRf/75pxo1amSzvlGjRtq1a1e2j/Og3dd7dADA7DKbptSsU4kCAB6M5s2ba/r06XJ2dlZAQICKFCmiJUuW3HW/vn37KiIiQt9//71+/PFHxcTEaNKkSRowYECm7ZOSkiRJ33//vcqUKWOz7faJv9zd3e/xamzl1nFywz1PLw0AAAAg59zd3RUaGqpy5cqpSJGb/Q5VqlTRjRs39Msvv1jbXbhwQQcPHlTVqlWt6wIDA/XPf/5T3333nV599VX95z//kSQ5OztLktLS0qxtq1atKhcXF504cUKhoaE2n8DAwGzX6+XlpYCAAG3evNlm/ebNm21qy2/o0QEAAADsrGLFiurQoYP69eunGTNmyNPTU8OGDVOZMmXUoUMHSTfHv7Ru3VoPPfSQLl26pLVr16pKlSqSpKCgIFksFi1dulRt2rSRm5ubPD099dprr+mVV15Renq6Hn30USUkJGjz5s3y8vJSjx49sl3f0KFDNWrUKIWEhCgsLEyzZs1SbGys3WdWuxOCDgAAAEyjID9ePGvWLA0aNEjt2rVTamqqmjRpomXLllnfRZOWlqaoqCj9/vvv8vLyUqtWrTR58mRJUpkyZTRmzBgNGzZMvXr1UmRkpGbPnq2xY8fK19dXMTExOnr0qHx8fFSnTh0NHz48R7UNHDhQCQkJevXVV3X27FlVrVpVS5YssZlxLb+xGNmZYcDOEhMT5e3trYSEBHl5edm7HAD52dqYO29vnvXL2DIbj5OZgvyXKACYxbVr1xQXF6fg4GC7T2OM3Hen+5vdbMAYHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDrMugbAdKbF786wrr9P1m+OBgAUXAVgXi3cg9y4r/ToAAAAoMBxdHSUJKWmptq5EuSFK1euSJJ1au17QY8OAAAACpwiRYqoaNGiOnfunJycnOTgwL/fm4FhGLpy5YrOnj0rHx8fa6C9FwQdAAAAFDgWi0WlS5dWXFycjh8/bu9ykMt8fHzk7+9/X8cg6AAo0DK85DOT8TkAAHNydnZWxYoVeXzNZJycnO6rJ+cWgg4AAAAKLAcHB7m6utq7DORDPMwIAAAAwHQIOgAAAABMh6ADAAAAwHQIOgAAAABMh6ADAAAAwHQIOgAAAABMh6ADAAAAwHQIOgAAAABMh6ADAAAAwHQIOgAAAABMh6ADAAAAwHQIOgAAAABMh6ADAAAAwHQIOgAAAABMh6ADAAAAwHQIOgAAAABMh6ADAAAAwHQIOgAAAABMh6ADAAAAwHSK2LsAAHig1sZkXNc8+sHXAQAA8hQ9OgAAAABMh6ADAAAAwHQIOgAAAABMhzE6AAqFafG7bZb7+9T8a+HWuJ1bbYIbP6CqAABAXqFHBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpEHQAAAAAmA5BBwAAAIDpFLF3AQBgD9Pid9u7BAAAkIfo0QEAAABgOgQdAAAAAKZD0AEAAABgOozRAYDbxW3MuC648YOvAwAA3DN6dAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOnwwlAABcfaGNvl5tH2qQMAAOR79OgAAAAAMB2CDgAAAADTIegAAAAAMB2CDgAAAADTIegAAAAAMJ0cB50NGzaoffv2CggIkMVi0aJFi+7Yft26dbJYLBk+p0+fvteaAQAAAOCOchx0kpOTVatWLU2dOjVH+x08eFCnTp2yfvz8/HJ6agAAAADIlhy/R6d169Zq3bp1jk/k5+cnHx+fHO8HAAAAADn1wMbohIWFqXTp0mrZsqU2b958x7YpKSlKTEy0+QAAAABAduV50CldurQ+/vhjffvtt/r2228VGBioZs2aaceOHVnuExMTI29vb+snMDAwr8sEAAAAYCI5fnQtpypVqqRKlSpZlxs2bKgjR45o8uTJmjNnTqb7REdHa8iQIdblxMREwg4AAACAbMvzoJOZhx9+WJs2bcpyu4uLi1xcXB5gRQAAAADMxC7v0YmNjVXp0qXtcWoAAAAAhUCOe3SSkpJ0+PBh63JcXJxiY2NVvHhxlStXTtHR0frjjz/0xRdfSJKmTJmi4OBgVatWTdeuXdMnn3yiNWvW6Mcff8y9qwAAAACAv8lx0Nm2bZuaN29uXb41lqZHjx6aPXu2Tp06pRMnTli3p6am6tVXX9Uff/yhokWLqmbNmlq1apXNMQAAAAAgN+U46DRr1kyGYWS5ffbs2TbLr7/+ul5//fUcFwYAAAAA98ouY3QAAAAAIC8RdAAAAACYDkEHAAAAgOkQdAAAAACYjl1eGAoABdm02Gk2y/3D+tupEgAAkBV6dAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYThF7FwAA92xtjBS/295VAACAfIgeHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDoEHQAAAACmU8TeBQBAgRS38a+fLyVIzaPtVwsAAMiAHh0AAAAApkPQAQAAAGA6BB0AAAAApsMYHQD51rTYaTbL/e1UBwAAKHjo0QEAAABgOgQdAAAAAKZD0AEAAABgOgQdAAAAAKbDZAQA8q+/v5RTknxq2qcOAABQ4NCjAwAAAMB0CDoAAAAATIegAwAAAMB0CDoAAAAATIegAwAAAMB0CDoAAAAATIegAwAAAMB0CDoAAAAATIegAwAAAMB0CDoAAAAATIegAwAAAMB0CDoAAAAATIegAwAAAMB0iti7AAAoEOI22rsCAACQA/ToAAAAADAdgg4AAAAA0yHoAAAAADAdgg4AAAAA0yHoAAAAADAdgg4AAAAA02F6aQAFxrT43fYuAQAAFBD06AAAAAAwHYIOAAAAANMh6AAAAAAwHYIOAAAAANMh6AAAAAAwHYIOAAAAANNhemkA+cPaGHtXAAAATIQeHQAAAACmQ9ABAAAAYDoEHQAAAACmwxgdAMgNmY0xah794OsAAACS6NEBAAAAYEIEHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDoEHQAAAACmk+Ogs2HDBrVv314BAQGyWCxatGjRXfdZt26d6tSpIxcXF4WGhmr27Nn3UCoAAAAAZE+Og05ycrJq1aqlqVOnZqt9XFyc2rZtq+bNmys2NlaDBw9W3759tWLFihwXCwAAAADZUSSnO7Ru3VqtW7fOdvuPP/5YwcHBmjRpkiSpSpUq2rRpkyZPnqyIiIicnh4AAAAA7irPx+hs2bJFLVq0sFkXERGhLVu2ZLlPSkqKEhMTbT4AAAAAkF15HnROnz6tUqVK2awrVaqUEhMTdfXq1Uz3iYmJkbe3t/UTGBiY12UCAAAAMJF8OetadHS0EhISrJ+TJ0/auyQAAAAABUiOx+jklL+/v86cOWOz7syZM/Ly8pKbm1um+7i4uMjFxSWvSwMAAABgUnneoxMeHq7Vq1fbrFu5cqXCw8Pz+tQAAAAACqkcB52kpCTFxsYqNjZW0s3po2NjY3XixAlJNx87i4yMtLb/5z//qaNHj+r111/Xr7/+qmnTpumrr77SK6+8kjtXAAAAAAC3yXHQ2bZtm2rXrq3atWtLkoYMGaLatWtr5MiRkqRTp05ZQ48kBQcH6/vvv9fKlStVq1YtTZo0SZ988glTSwMAAADIMzkeo9OsWTMZhpHl9tmzZ2e6z86dO3N6KgAAAAC4J/ly1jUAAAAAuB8EHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDoEHQAAAACmQ9ABAAAAYDoEHQAAAACmU8TeBQCAJE2L323vEgAAgInQowMAAADAdAg6AAAAAEyHoAMAAADAdAg6AAAAAEyHoAMAAADAdAg6AAAAAEyHoAMAAADAdHiPDgDcp8zeAdTfp6YdKgEAALfQowMAAADAdAg6AAAAAEyHoAMAAADAdAg6AAAAAEyHoAMAAADAdAg6AAAAAEyHoAMAAADAdAg6AAAAAEyHoAMAAADAdAg6AAAAAEyHoAMAAADAdAg6AAAAAEyHoAMAAADAdAg6AAAAAEyHoAMAAADAdAg6AAAAAEyHoAMAAADAdAg6AAAAAEyHoAMAAADAdAg6AAAAAEyHoAMAAADAdAg6AAAAAEyHoAMAAADAdAg6AAAAAEyHoAMAAADAdIrYuwAAMKNp8bul2Gl3bdc/rP8DqAYAgMKHHh0AAAAApkPQAQAAAGA6BB0AAAAApkPQAQAAAGA6BB0AAAAApkPQAQAAAGA6BB0AAAAApkPQAQAAAGA6BB0AAAAApkPQAQAAAGA6BB0AAAAApkPQAQAAAGA6BB0AAAAAplPE3gUAKKTWxti7AgAAYGL06AAAAAAwHYIOAAAAANPh0TUAD9y02GlS/G57lwEAAEyMHh0AAAAApkPQAQAAAGA6BB0AAAAApkPQAQAAAGA6BB0AAAAApkPQAQAAAGA6TC8NAA9K3Ebb5eDG9qkDAIBCgB4dAAAAAKZD0AEAAABgOgQdAAAAAKZD0AEAAABgOgQdAAAAAKZD0AEAAABgOgQdAAAAAKZD0AEAAABgOgQdAAAAAKZD0AEAAABgOgQdAAAAAKZD0AEAAABgOgQdAAAAAKZTxN4FAIBpxW20dwUAABRa9OgAAAAAMB2CDgAAAADTIegAAAAAMB3G6AB4MNbG/PVz/G771QEAAAqFe+rRmTp1qsqXLy9XV1c1aNBA//vf/7JsO3v2bFksFpuPq6vrPRcMAAAAAHeT46Dz5ZdfasiQIRo1apR27NihWrVqKSIiQmfPns1yHy8vL506dcr6OX78+H0VDQAAAAB3kuOg8/7776tfv37q1auXqlatqo8//lhFixbVZ599luU+FotF/v7+1k+pUqXuq2gAAAAAuJMcBZ3U1FRt375dLVq0+OsADg5q0aKFtmzZkuV+SUlJCgoKUmBgoDp06KB9+/bd8TwpKSlKTEy0+QAAAABAduVoMoLz588rLS0tQ49MqVKl9Ouvv2a6T6VKlfTZZ5+pZs2aSkhI0MSJE9WwYUPt27dPZcuWzXSfmJgYjRkzJielAchP/j7xAAAAgB3k+fTS4eHhioyMVFhYmJo2barvvvtOvr6+mjFjRpb7REdHKyEhwfo5efJkXpcJAAAAwERy1KNTsmRJOTo66syZMzbrz5w5I39//2wdw8nJSbVr19bhw4ezbOPi4iIXF5eclAYAAAAAVjnq0XF2dlbdunW1evVq67r09HStXr1a4eHh2TpGWlqa9uzZo9KlS+esUgAAAADIphy/MHTIkCHq0aOH6tWrp4cfflhTpkxRcnKyevXqJUmKjIxUmTJlFBNz8xn9t956S4888ohCQ0MVHx+v9957T8ePH1ffvn1z90oAAAAA4P/kOOh07dpV586d08iRI3X69GmFhYXphx9+sE5QcOLECTk4/NVRdOnSJfXr10+nT59WsWLFVLduXf3000+qWrVq7l0FAAAAAPyNxTAMw95F3E1iYqK8vb2VkJAgLy8ve5cD4G7uMuvatPjdD6iQfC64sfqH9bd3FQAAFCjZzQZ5PusaAAAAADxoBB0AAAAApkPQAQAAAGA6BB0AAAAApkPQAQAAAGA6OZ5eGgCQe6bFTrNZZhY2AAByBz06AAAAAEyHoAMAAADAdHh0DQCQqyav/C3Pjv1Ky4fy7NgAAHOhRwcAAACA6RB0AAAAAJgOQQcAAACA6RB0AAAAAJgOkxEAwD04GX/1vo/x+5ELGdalnLs5kJ9B9wAA3B+CDgAUMnk5KxoAAPkFj64BAAAAMB2CDgAAAADTIegAAAAAMB2CDgAAAADTIegAAAAAMB1mXQMAFBh5PWMc03oDgHnQowMAAADAdOjRAQDg/+RljxG9RQDwYNGjAwAAAMB0CDoAAAAATIegAwAAAMB0CDoAAAAATIegAwAAAMB0mHUNgGmdjL9q7xIAAICd0KMDAAAAwHQIOgAAAABMh0fXACAfyssXVwIAUBjQowMAAADAdAg6AAAAAEyHoAMAAADAdBijAwB2UjZxe4Z1OzJpV8era94XAwCAydCjAwAAAMB0CDoAAAAATIegAwAAAMB0GKMD4P6tjbF3BaZx+7id373q2qkSAAAKNnp0AAAAAJgOQQcAAACA6RB0AAAAAJgOQQcAAACA6TAZAYB7Ni122s0f4nfbtxAAAIDbEHQAAHgAJq/8LU+P/0rLh/L0+ABQ0PDoGgAAAADTIegAAAAAMB2CDgAAAADTYYwOALs5GX/V3iUAAACTokcHAAAAgOkQdAAAAACYDkEHAAAAgOkwRgdAzqyN+etnXhQKAADyKYIOAAAmkJcvJOVlpAAKIh5dAwAAAGA6BB0AAAAApkPQAQAAAGA6BB0AAAAApkPQAQAAAGA6BB0AAAAApkPQAQAAAGA6vEcHwB3d/m6OR05csP580uHqgy6n0CmbuF2PxCdYl38u94IdqwEAoOCgRwcAAACA6dCjAxRwefk2dAAAgIKKoAMAAO4or/9B5ZWWD+Xp8QEUTjy6BgAAAMB06NEBCpFHTsy0WWZgOwAAMCt6dAAAAACYDkEHAAAAgOkQdAAAAACYDmN0AACAXeXlrG7M6AYUXgQdADaYsKDg4Z4BAJARj64BAAAAMB16dAAAgGnxslOg8KJHBwAAAIDp0KMDAAXI7eNxAABA5ujRAQAAAGA69OgAyGCJw2Hrz78nfilJquPV1V7lFHp/vx938mR6aB5XAgBAwUHQAR6AvB4MCwAAAFs8ugYAAADAdAg6AAAAAEyHR9cAZMuO/xurczab40UAAADsiaADAABwj/JyDCYvIwXuD4+uAQAAADAdenQAMSsazO/2F43+XO4FO1UCAMCDQdABYB1/IzEGpyC79b6d3/92P2+5dV951w4AoLDg0TUAAAAApkPQAQAAAGA6PLoGAABQCDFjHMyOoIMCgwkDct+tAep3GpdTNnH7gyoHD0Bm43jqeHW96347Mhn3k539AACwl3t6dG3q1KkqX768XF1d1aBBA/3vf/+7Y/uvv/5alStXlqurq2rUqKFly5bdU7EAAAAAkB057tH58ssvNWTIEH388cdq0KCBpkyZooiICB08eFB+fn4Z2v/000/q1q2bYmJi1K5dO82fP18dO3bUjh07VL169Vy5CGQfvSIAABQM/J0N3J8c9+i8//776tevn3r16qWqVavq448/VtGiRfXZZ59l2v6DDz5Qq1atNHToUFWpUkVjx45VnTp19NFHH9138QAAAACQmRz16KSmpmr79u2Kjo62rnNwcFCLFi20ZcuWTPfZsmWLhgwZYrMuIiJCixYtyvI8KSkpSklJsS4nJCRIkhITE3NSbp6Zuibv3jMS9VjevuPiWnJSnh4f+Vvy1ZRM16c4XH/AlSAvpRa5lmFdytXrWba55nj3/y+kXsl4zOzsB6Bwilm0w94l3LO8/l0sLxXk31Fz4lYmMAzjju1yFHTOnz+vtLQ0lSpVymZ9qVKl9Ouvv2a6z+nTpzNtf/r06SzPExMTozFjxmRYHxgYmJNyC6Th9i4AgAlkZxzkX22+0tv3dJZ73Q8A8jN+F8tcfvxeLl++LG9v7yy358tZ16Kjo216gdLT03Xx4kWVKFFCFosl030SExMVGBiokydPysvL60GViixwP/IP7kX+wb3IX7gf+Qf3Iv/gXuQf3IusGYahy5cvKyAg4I7tchR0SpYsKUdHR505c8Zm/ZkzZ+Tv75/pPv7+/jlqL0kuLi5ycXGxWefj45OtGr28vPjDkI9wP/IP7kX+wb3IX7gf+Qf3Iv/gXuQf3IvM3akn55YcTUbg7OysunXravXq1dZ16enpWr16tcLDwzPdJzw83Ka9JK1cuTLL9gAAAABwv3L86NqQIUPUo0cP1atXTw8//LCmTJmi5ORk9erVS5IUGRmpMmXKKCYmRpI0aNAgNW3aVJMmTVLbtm21YMECbdu2TTNnzszdKwEAAACA/5PjoNO1a1edO3dOI0eO1OnTpxUWFqYffvjBOuHAiRMn5ODwV0dRw4YNNX/+fL355psaPny4KlasqEWLFuX6O3RcXFw0atSoDI+8wT64H/kH9yL/4F7kL9yP/IN7kX9wL/IP7sX9sxh3m5cNAAAAAAqYHL8wFAAAAADyO4IOAAAAANMh6AAAAAAwHYIOAAAAANMp8EEnJiZG9evXl6enp/z8/NSxY0cdPHjQ3mUVStOnT1fNmjWtL7YKDw/X8uXL7V0WJL3zzjuyWCwaPHiwvUsplEaPHi2LxWLzqVy5sr3LKrT++OMPPf/88ypRooTc3NxUo0YNbdu2zd5lFUrly5fP8N+GxWJRVFSUvUsrdNLS0jRixAgFBwfLzc1NISEhGjt2rJizyj4uX76swYMHKygoSG5ubmrYsKG2bt1q77IKnBxPL53frF+/XlFRUapfv75u3Lih4cOH64knntD+/fvl7u5u7/IKlbJly+qdd95RxYoVZRiGPv/8c3Xo0EE7d+5UtWrV7F1eobV161bNmDFDNWvWtHcphVq1atW0atUq63KRIgX+f78F0qVLl9SoUSM1b95cy5cvl6+vrw4dOqRixYrZu7RCaevWrUpLS7Mu7927Vy1bttQzzzxjx6oKpwkTJmj69On6/PPPVa1aNW3btk29evWSt7e3Bg4caO/yCp2+fftq7969mjNnjgICAjR37ly1aNFC+/fvV5kyZexdXoFhuumlz507Jz8/P61fv15NmjSxdzmFXvHixfXee++pT58+9i6lUEpKSlKdOnU0bdo0jRs3TmFhYZoyZYq9yyp0Ro8erUWLFik2NtbepRR6w4YN0+bNm7Vx40Z7l4JMDB48WEuXLtWhQ4dksVjsXU6h0q5dO5UqVUqffvqpdV3nzp3l5uamuXPn2rGywufq1avy9PTU4sWL1bZtW+v6unXrqnXr1ho3bpwdqytYCvyja7dLSEiQdPMXbNhPWlqaFixYoOTkZIWHh9u7nEIrKipKbdu2VYsWLexdSqF36NAhBQQEqEKFCurevbtOnDhh75IKpSVLlqhevXp65pln5Ofnp9q1a+s///mPvcuCpNTUVM2dO1e9e/cm5NhBw4YNtXr1av3222+SpF27dmnTpk1q3bq1nSsrfG7cuKG0tDS5urrarHdzc9OmTZvsVFXBZKpnJ9LT0zV48GA1atRI1atXt3c5hdKePXsUHh6ua9euycPDQwsXLlTVqlXtXVahtGDBAu3YsYNnevOBBg0aaPbs2apUqZJOnTqlMWPGqHHjxtq7d688PT3tXV6hcvToUU2fPl1DhgzR8OHDtXXrVg0cOFDOzs7q0aOHvcsr1BYtWqT4+Hj17NnT3qUUSsOGDVNiYqIqV64sR0dHpaWlafz48erevbu9Syt0PD09FR4errFjx6pKlSoqVaqU/t//+3/asmWLQkND7V1egWKqoBMVFaW9e/eSdu2oUqVKio2NVUJCgr755hv16NFD69evJ+w8YCdPntSgQYO0cuXKDP8ihAfv7/8iWrNmTTVo0EBBQUH66quveKzzAUtPT1e9evX09ttvS5Jq166tvXv36uOPPybo2Nmnn36q1q1bKyAgwN6lFEpfffWV5s2bp/nz56tatWqKjY3V4MGDFRAQwH8bdjBnzhz17t1bZcqUkaOjo+rUqaNu3bpp+/bt9i6tQDFN0Hn55Ze1dOlSbdiwQWXLlrV3OYWWs7Oz9V8b6tatq61bt+qDDz7QjBkz7FxZ4bJ9+3adPXtWderUsa5LS0vThg0b9NFHHyklJUWOjo52rLBw8/Hx0UMPPaTDhw/bu5RCp3Tp0hn+4aVKlSr69ttv7VQRJOn48eNatWqVvvvuO3uXUmgNHTpUw4YN07PPPitJqlGjho4fP66YmBiCjh2EhIRo/fr1Sk5OVmJiokqXLq2uXbuqQoUK9i6tQCnwY3QMw9DLL7+shQsXas2aNQoODrZ3Sfib9PR0paSk2LuMQufxxx/Xnj17FBsba/3Uq1dP3bt3V2xsLCHHzpKSknTkyBGVLl3a3qUUOo0aNcrwCoLffvtNQUFBdqoIkjRr1iz5+fnZDLzGg3XlyhU5ONj+Wujo6Kj09HQ7VQRJcnd3V+nSpXXp0iWtWLFCHTp0sHdJBUqB79GJiorS/PnztXjxYnl6eur06dOSJG9vb7m5udm5usIlOjparVu3Vrly5XT58mXNnz9f69at04oVK+xdWqHj6emZYZyau7u7SpQowfg1O3jttdfUvn17BQUF6c8//9SoUaPk6Oiobt262bu0QueVV15Rw4YN9fbbb6tLly763//+p5kzZ2rmzJn2Lq3QSk9P16xZs9SjRw+mXbej9u3ba/z48SpXrpyqVaumnTt36v3331fv3r3tXVqhtGLFChmGoUqVKunw4cMaOnSoKleurF69etm7tAKlwP8fZfr06ZKkZs2a2ayfNWsWAxofsLNnzyoyMlKnTp2St7e3atasqRUrVqhly5b2Lg2wq99//13dunXThQsX5Ovrq0cffVQ///yzfH197V1aoVO/fn0tXLhQ0dHReuuttxQcHKwpU6Yw4NqOVq1apRMnTvALtZ19+OGHGjFihPr376+zZ88qICBAL774okaOHGnv0gqlhIQERUdH6/fff1fx4sXVuXNnjR8/Xk5OTvYurUAx3Xt0AAAAAKDAj9EBAAAAgNsRdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAEzGYrFo0aJF9i7jrpo1a6bBgwfb7fxNmjTR/Pnzrcv59Xs7f/68/Pz89Pvvv9u7FAAoUAg6AFCAnDt3Ti+99JLKlSsnFxcX+fv7KyIiQps3b7Z3aVbHjh2TxWK542f27Nn67rvvNHbsWLvUuGTJEp05c0bPPvusXc6fEyVLllRkZKRGjRpl71IAoEApYu8CAADZ17lzZ6Wmpurzzz9XhQoVdObMGa1evVoXLlywd2lWgYGBOnXqlHV54sSJ+uGHH7Rq1SrrOm9vb7m5udmjPEnSv//9b/Xq1UsODvb/977r16/Lycnpjm169eqlunXr6r333lPx4sUfUGUAULDZ///wAIBsiY+P18aNGzVhwgQ1b95cQUFBevjhhxUdHa0nn3wyy/1OnjypLl26yMfHR8WLF1eHDh107NgxmzaffPKJqlSpIldXV1WuXFnTpk2zbrvVQ7NgwQI1bNhQrq6uql69utavX5/p+RwdHeXv72/9eHh4qEiRIjbr3NzcMjy6Vr58eY0bN06RkZHy8PBQUFCQlixZonPnzqlDhw7y8PBQzZo1tW3bNpvzbdq0SY0bN5abm5sCAwM1cOBAJScnZ/l9nDt3TmvWrFH79u0zbDt//ryeeuopFS1aVBUrVtSSJUtstq9fv14PP/ywXFxcVLp0aQ0bNkw3btywuYYpU6bY7BMWFqbRo0dbly0Wi6ZPn64nn3xS7u7uGj9+vC5duqTu3bvL19dXbm5uqlixombNmmXdp1q1agoICNDChQuzvC4AgC2CDgAUEB4eHvLw8NCiRYuUkpKSrX2uX7+uiIgIeXp6auPGjdq8ebM8PDzUqlUrpaamSpLmzZunkSNHavz48Tpw4IDefvttjRgxQp9//rnNsYYOHapXX31VO3fuVHh4uNq3b5/rPUmTJ09Wo0aNtHPnTrVt21b/+Mc/FBkZqeeff147duxQSEiIIiMjZRiGJOnIkSNq1aqVOnfurN27d+vLL7/Upk2b9PLLL2d5jk2bNqlo0aKqUqVKhm1jxoxRly5dtHv3brVp00bdu3fXxYsXJUl//PGH2rRpo/r162vXrl2aPn26Pv30U40bNy7H1zl69Gg99dRT2rNnj3r37q0RI0Zo//79Wr58uQ4cOKDp06erZMmSNvs8/PDD2rhxY47PBQCFlgEAKDC++eYbo1ixYoarq6vRsGFDIzo62ti1a5dNG0nGwoULDcMwjDlz5hiVKlUy0tPTrdtTUlIMNzc3Y8WKFYZhGEZISIgxf/58m2OMHTvWCA8PNwzDMOLi4gxJxjvvvGPdfv36daNs2bLGhAkT7lrzqFGjjFq1amVY37RpU2PQoEHW5aCgIOP555+3Lp86dcqQZIwYMcK6bsuWLYYk49SpU4ZhGEafPn2MF154wea4GzduNBwcHIyrV69mWs/kyZONChUqZFgvyXjzzTety0lJSYYkY/ny5YZhGMbw4cMzfJdTp041PDw8jLS0NOs1TJ482ea4tWrVMkaNGmVznsGDB9u0ad++vdGrV69M673llVdeMZo1a3bHNgCAv9CjAwAFSOfOnfXnn39qyZIlatWqldatW6c6depo9uzZmbbftWuXDh8+LE9PT2uPUPHixXXt2jUdOXJEycnJOnLkiPr06WPd7uHhoXHjxunIkSM2xwoPD7f+XKRIEdWrV08HDhzI1eurWbOm9edSpUpJkmrUqJFh3dmzZ63XN3v2bJvaIyIilJ6erri4uEzPcfXqVbm6ut71/O7u7vLy8rKe68CBAwoPD5fFYrG2adSokZKSknI8I1q9evVsll966SUtWLBAYWFhev311/XTTz9l2MfNzU1XrlzJ0XkAoDBjMgIAKGBcXV3VsmVLtWzZUiNGjFDfvn01atQo9ezZM0PbpKQk1a1bV/PmzcuwzdfXV0lJSZKk//znP2rQoIHNdkdHxzyp/07+Pij/VqDIbF16erqkm9f34osvauDAgRmOVa5cuUzPUbJkSV26dOmu5791vlvnyg4HBwfrY3W3XL9+PUM7d3d3m+XWrVvr+PHjWrZsmVauXKnHH39cUVFRmjhxorXNxYsX5evrm+1aAKCwo0cHAAq4qlWrZjn4vk6dOjp06JD8/PwUGhpq8/H29lapUqUUEBCgo0ePZtgeHBxsc6yff/7Z+vONGze0ffv2TMe5PEh16tTR/v37M9QeGhoqZ2fnTPepXbu2Tp8+nWXYyUqVKlW0ZcsWmyCzefNmeXp6qmzZspJuhse/zziXmJiYZc/S7Xx9fdWjRw/NnTtXU6ZM0cyZM2227927V7Vr185RzQBQmBF0AKCAuHDhgh577DHNnTtXu3fvVlxcnL7++mu9++676tChQ6b7dO/eXSVLllSHDh20ceNGxcXFad26dRo4cKD1casxY8YoJiZG//73v/Xbb79pz549mjVrlt5//32bY02dOlULFy7Ur7/+qqioKF26dEm9e/fO8+u+kzfeeEM//fSTXn75ZcXGxurQoUNavHjxHScjqF27tkqWLJnjdw/1799fJ0+e1IABA/Trr79q8eLFGjVqlIYMGWKdpvqxxx7TnDlztHHjRu3Zs0c9evTIVs/YyJEjtXjxYh0+fFj79u3T0qVLbULklStXtH37dj3xxBM5qhkACjMeXQOAAsLDw0MNGjTQ5MmTdeTIEV2/fl2BgYHq16+fhg8fnuk+RYsW1YYNG/TGG2+oU6dOunz5ssqUKaPHH39cXl5ekqS+ffuqaNGieu+99zR06FC5u7urRo0aNlM/S9I777yjd955R7GxsQoNDdWSJUsyzAz2oNWsWVPr16/Xv/71LzVu3FiGYSgkJERdu3bNch9HR0f16tVL8+bNU7t27bJ9rjJlymjZsmUaOnSoatWqpeLFi6tPnz568803rW2io6MVFxendu3aydvbW2PHjs1Wj46zs7Oio6N17Ngxubm5qXHjxlqwYIF1++LFi1WuXDk1btw42/UCQGFnMW5/mBgAgL85duyYgoODtXPnToWFhdm7nFxx+vRpVatWTTt27FBQUJC9y7mrRx55RAMHDtRzzz1n71IAoMDg0TUAQKHj7++vTz/9VCdOnLB3KXd1/vx5derUSd26dbN3KQBQoNCjAwC4IzP26AAAzI+gAwAAAMB0eHQNAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYDkEHAAAAgOkQdAAAAACYzv8HOa+5RcFMObkAAAAASUVORK5CYII=",
|
||
"text/plain": [
|
||
"<Figure size 1000x600 with 1 Axes>"
|
||
]
|
||
},
|
||
"metadata": {},
|
||
"output_type": "display_data"
|
||
}
|
||
],
|
||
"source": [
|
||
"# Define x-axis range\n",
|
||
"x = np.linspace(3, 8, 1000)\n",
|
||
"\n",
|
||
"# Generate random samples for histograms\n",
|
||
"prior_samples = norm.rvs(loc=mu0, scale=sigma0, size=1000)\n",
|
||
"likelihood_samples = norm.rvs(loc=ybar, scale=sigma / np.sqrt(n), size=1000)\n",
|
||
"posterior_samples = norm.rvs(loc=mu_n, scale=sigma_n, size=1000)\n",
|
||
"\n",
|
||
"# Plot histograms\n",
|
||
"plt.figure(figsize=(10, 6))\n",
|
||
"plt.hist(prior_samples, bins=30, density=True, alpha=0.5, label=\"Prior\")\n",
|
||
"plt.hist(likelihood_samples, bins=30, density=True, alpha=0.5, label=\"Likelihood\")\n",
|
||
"plt.hist(posterior_samples, bins=30, density=True, alpha=0.5, label=\"Posterior\")\n",
|
||
"plt.title(\"Hours College Students Sleep\")\n",
|
||
"plt.xlabel(\"Sleep Time (hours)\")\n",
|
||
"\n",
|
||
"plt.legend()\n",
|
||
"plt.show()"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 15,
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
" Prior Data (Sample Mean) Posterior\n",
|
||
"0 Precision 1.0 44.444 45.444\n",
|
||
"1 SD 1.0 0.150 0.148\n",
|
||
"2 Mean 5.5 4.844 4.859\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"# Create a DataFrame\n",
|
||
"df = pd.DataFrame({\n",
|
||
" \"\": [\"Precision\", \"SD\", \"Mean\"],\n",
|
||
" \"Prior\": [1 / sigma0**2, sigma0, mu0],\n",
|
||
" \"Data (Sample Mean)\": [n / sigma**2, sigma / np.sqrt(n), ybar],\n",
|
||
" \"Posterior\": [1 / sigma_n**2, sigma_n, mu_n]\n",
|
||
"})\n",
|
||
"\n",
|
||
"# Display the DataFrame\n",
|
||
"print(df.round(3))"
|
||
]
|
||
}
|
||
],
|
||
"metadata": {
|
||
"kernelspec": {
|
||
"display_name": "venv",
|
||
"language": "python",
|
||
"name": "python3"
|
||
},
|
||
"language_info": {
|
||
"codemirror_mode": {
|
||
"name": "ipython",
|
||
"version": 3
|
||
},
|
||
"file_extension": ".py",
|
||
"mimetype": "text/x-python",
|
||
"name": "python",
|
||
"nbconvert_exporter": "python",
|
||
"pygments_lexer": "ipython3",
|
||
"version": "3.12.3"
|
||
}
|
||
},
|
||
"nbformat": 4,
|
||
"nbformat_minor": 2
|
||
}
|