COGMOD-HWI/HW03/hw3.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Problem 4\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "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": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Sample Mean (ybar): 4.84 hours\n",
      "Posterior Mean (µ_n): 4.89\n",
      "Posterior Standard Deviation (σ_n): 0.15\n"
     ]
    }
   ],
   "source": [
    "np.random.seed(42)\n",
    "n = 100 #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 = 7\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": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA04AAAIjCAYAAAA0vUuxAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjEsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvc2/+5QAAAAlwSFlzAAAPYQAAD2EBqD+naQAApGdJREFUeJzs3Xd4m9Xd//G3JNvylPdIYidx7BCyyGRkkAGBkLAClNk2hF0aKPwoUNI+ZfOk8DBLKYUWSFmFll12yGAkYWVBQgjZ03bivYes+/eHfMtWZDuxLVu29Xldly7Z59ySjmxw/PE553sshmEYiIiIiIiISIusgR6AiIiIiIhId6fgJCIiIiIichgKTiIiIiIiIoeh4CQiIiIiInIYCk4iIiIiIiKHoeAkIiIiIiJyGApOIiIiIiIih6HgJCIiIiIichgKTiIiIiIiIoeh4CQi0oWmTZvGtGnTPJ/v3LkTi8XCokWLAjamnmLgwIHMmzcv0MPoNHfeeScWiyXQwxARkRYoOImItGLbtm1cc801DBo0iPDwcBwOB5MmTeKxxx6jqqoq0MPrMjt37uSyyy4jKyuL8PBw0tLSmDJlCnfccYfXdX/96197dAh8//33ufPOO/36nOXl5dxxxx2MGDGCqKgoEhMTGT16NDfccAP79+/362t1tnnz5hEdHd1iv8Vi4brrruvCEYmIdJ2QQA9ARKS7eu+99zj//POx2+3MnTuXESNGUFtbyxdffMEtt9zCxo0befrppwM9zE63detWjj32WCIiIrj88ssZOHAgOTk5rFmzhvvvv5+77rrLc+1f//pXkpKSeuzM0Pvvv88TTzzht/BUV1fHlClT+PHHH7n00ku5/vrrKS8vZ+PGjbz88succ8459O3b1y+vJSIinUvBSUSkGTt27OCiiy5iwIABLF26lD59+nj65s+fz9atW3nvvfcCOMKu88gjj1BeXs66desYMGCAV9+BAwcCNKqe4a233mLt2rW89NJLXHLJJV591dXV1NbWBmhkvUtFRQVRUVGBHoaI9HJaqici0owHHniA8vJynnnmGa/QZMrOzuaGG27wfO50OrnnnnvIysrCbrczcOBAfv/731NTU9Ou1//xxx/52c9+RkJCAuHh4YwfP5533nnH57rvvvuOqVOnEhERQXp6Ovfeey/PPfccFouFnTt3el37wQcfcOKJJxIVFUVMTAynn346GzduPOxYtm3bRnp6uk9oAkhJSfF8PHDgQDZu3Minn36KxWLBYrF49nO1tH9n0aJFPmM1DIN7772X9PR0IiMjmT59eovjLC4u5sYbbyQjIwO73U52djb3338/LpfLc425j+zBBx/k6aef9nyPjj32WL755hvPdfPmzeOJJ54A8Iy/6ZhfeeUVxo0bR0xMDA6Hg5EjR/LYY48d9msHMGnSJJ8+c+nn4bz44ouMGzeOiIgIEhISuOiii9izZ4/PdV999RWnnXYasbGxREZGMnXqVFasWOF1jfl9+PHHH7ngggtwOBwkJiZyww03UF1dfdixtMeBAwe44oorSE1NJTw8nFGjRvHPf/7T65rly5djsVhYvny5V3tzewDN5YLbtm1j9uzZxMTE8POf/xyALVu2cN5555GWlkZ4eDjp6elcdNFFlJSUdMp7E5HgohknEZFm/Pe//2XQoEFMnDjxiK6/8sor+ec//8nPfvYzfvvb3/LVV1+xcOFCNm3axJtvvtmm1964cSOTJk2iX79+3HbbbURFRfHvf/+bOXPm8Prrr3POOecAsG/fPqZPn47FYmHBggVERUXxj3/8A7vd7vOcL7zwApdeeikzZ87k/vvvp7KykieffJLJkyezdu1aBg4c2OJ4BgwYwCeffMLSpUs56aSTWrzu0Ucf5frrryc6Opo//OEPAKSmprbpvQPcfvvt3HvvvcyePZvZs2ezZs0aTj31VJ/ZmcrKSqZOncq+ffu45ppr6N+/PytXrmTBggXk5OTw6KOPel3/8ssvU1ZWxjXXXIPFYuGBBx7g3HPPZfv27YSGhnLNNdewf/9+Fi9ezAsvvOD12MWLF3PxxRdz8sknc//99wOwadMmVqxY4RWgD2WGzeeff57/+Z//aXPxh/vuu48//vGPXHDBBVx55ZUcPHiQxx9/nClTprB27Vri4uIAWLp0KbNmzWLcuHHccccdWK1WnnvuOU466SQ+//xzjjvuOK/nveCCCxg4cCALFy7kyy+/5M9//jNFRUU8//zzRzSu/Pz8I7quqqqKadOmsXXrVq677joyMzP5z3/+w7x58yguLm71a9cap9PJzJkzmTx5Mg8++CCRkZHU1tYyc+ZMampquP7660lLS2Pfvn28++67FBcXExsb267XEhHxMERExEtJSYkBGGefffYRXb9u3ToDMK688kqv9ptvvtkAjKVLl3rapk6dakydOtXz+Y4dOwzAeO655zxtJ598sjFy5Eijurra0+ZyuYyJEycagwcP9rRdf/31hsViMdauXetpKygoMBISEgzA2LFjh2EYhlFWVmbExcUZV111ldf4cnNzjdjYWJ/2Q23YsMGIiIgwAGP06NHGDTfcYLz11ltGRUWFz7XDhw/3en+mO+64w2jun5znnnvOa6wHDhwwwsLCjNNPP91wuVye637/+98bgHHppZd62u655x4jKirK+Omnn7ye87bbbjNsNpuxe/duwzAav8aJiYlGYWGh57q3337bAIz//ve/nrb58+c3O84bbrjBcDgchtPpbP6L1ILKykpjyJAhBmAMGDDAmDdvnvHMM88YeXl5Ptce+jXauXOnYbPZjPvuu8/ruu+//94ICQnxtLtcLmPw4MHGzJkzvb5mlZWVRmZmpnHKKaf4vMZZZ53l9Zy//vWvDcBYv359q+/n0ksvNYBWb/Pnz/dc/+ijjxqA8eKLL3raamtrjQkTJhjR0dFGaWmpYRiGsWzZMgMwli1b5vV6zf3/YY7htttu87p27dq1BmD85z//afU9iIi0l5bqiYgcorS0FICYmJgjuv79998H4KabbvJq/+1vfwvQpr1QhYWFLF26lAsuuICysjLy8/PJz8+noKCAmTNnsmXLFvbt2wfAhx9+yIQJExg9erTn8QkJCZ5lS6bFixdTXFzMxRdf7Hm+/Px8bDYbxx9/PMuWLWt1TMOHD2fdunX84he/YOfOnTz22GPMmTOH1NRU/v73vx/xezsSn3zyCbW1tVx//fVeszM33nijz7X/+c9/OPHEE4mPj/d6XzNmzKC+vp7PPvvM6/oLL7yQ+Ph4z+cnnngiANu3bz/suOLi4qioqGDx4sVtej8RERF89dVX3HLLLYB7aeIVV1xBnz59uP7661tdyvnGG2/gcrm44IILvN5fWloagwcP9nzf1q1bx5YtW7jkkksoKCjwXFdRUcHJJ5/MZ5995rV0Edz79Jq6/vrrgcb/llsTHh7O4sWLm70d6v333yctLY2LL77Y0xYaGspvfvMbysvL+fTTTw/7ei259tprvT43Z5Q++ugjKisr2/28IiIt0VI9EZFDmPtOysrKjuj6Xbt2YbVayc7O9mpPS0sjLi6OXbt2HfFrb926FcMw+OMf/8gf//jHZq85cOAA/fr1Y9euXUyYMMGn/9BxbNmyBaDFZXZHss/mqKOO4oUXXqC+vp4ffviBd999lwceeICrr76azMxMZsyYcdjnOBLm12rw4MFe7cnJyV6hB9zv67vvviM5ObnZ5zq0cEX//v29Pjefr6io6LDj+vWvf82///1vZs2aRb9+/Tj11FO54IILOO200w772NjYWB544AEeeOABdu3axZIlS3jwwQf5y1/+QmxsLPfee2+zj9uyZQuGYfh8LUyhoaGe6wAuvfTSFsdQUlLi9fU79DmzsrKwWq0+++KaY7PZjvj7vWvXLgYPHozV6v132qFDh3r62yMkJIT09HSvtszMTG666SYefvhhXnrpJU488UTOOussfvGLX2iZnoj4hYKTiMghHA4Hffv2ZcOGDW16nD8OLzVnBm6++WZmzpzZ7DWHBqMjfc4XXniBtLQ0n/6QkCP/p8BmszFy5EhGjhzJhAk
      "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",
    "# Prior distribution\n",
    "prior = norm.pdf(x, mu0, sigma0)\n",
    "\n",
    "# Likelihood distribution (centered at sample mean)\n",
    "likelihood = norm.pdf(x, ybar, sigma / np.sqrt(n))\n",
    "\n",
    "# Posterior distribution\n",
    "posterior = norm.pdf(x, mu_n, sigma_n)\n",
    "\n",
    "# Plot\n",
    "plt.figure(figsize=(10, 6))\n",
    "plt.plot(x, prior, label=\"Prior\", linestyle=\"--\")\n",
    "plt.plot(x, likelihood, label=\"Likelihood\", linestyle=\"-.\")\n",
    "plt.plot(x, posterior, label=\"Posterior\", linestyle=\"-\")\n",
    "plt.title(\"College Students Sleep Hours\")\n",
    "plt.xlabel(\"Sleep Time (hours)\")\n",
    "plt.ylabel(\"Density\")\n",
    "plt.legend()\n",
    "plt.show()\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "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    7.0               4.844      4.892\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
}