Add nix-shell and prob things

4 files changed, 643 insertions, 0 deletions

diff --git a/sem6/prob/miniexam/Untitled.ipynb b/sem6/prob/miniexam/Untitled.ipynb
new file mode 100644
index 0000000..3d81069
--- /dev/null
+++ b/sem6/prob/miniexam/Untitled.ipynb
@@ -0,0 +1,346 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": 41,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import numpy as np\n",
+    "import matplotlib.pyplot as plt\n",
+    "from scipy import stats"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Problem 1\n",
+    "\n",
+    "Man kan sige at distribution af et enkelt kast af den 8 kantede terning er $1/8$ for hver værdi i $[1, 8]$.\n",
+    "Ved at summe dem sammen kommer de til at ligne en gaussian random variable.\n",
+    "\n",
+    "Man ved også de to limits for distributionen, nemlig $[50 \\cdot 1, 50 \\cdot 8]$, og hver imellem må mean ligge.\n",
+    "Ud fra intervallet og mean kan man regne varience."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 10,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Man kan først finde størrelsen af intervallet, eller var * 2\n",
+    "var2 = 50 * 8 - 50 * 1\n",
+    "sigma2 = var2/2\n",
+    "mean = 50 * 1 + sigma2\n",
+    "\n",
+    "# Ved coin flip får man egentlig det samme distribution"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Part A\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 8,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "<BarContainer object of 500 artists>"
+      ]
+     },
+     "execution_count": 8,
+     "metadata": {},
+     "output_type": "execute_result"
+    },
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAD4CAYAAAD1jb0+AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAANXklEQVR4nO3db4hl9X3H8fenu5s/RCFrd5RFpWODhEpoVhmsYAk2JulqSzWQQITKPrBsHkRQGiibBNrkWVoa0yclsKmShRqDoKJoaLNsDRII2lmz6i6rXZNuU+OyM1aC5kla9dsH9wwZxpm9d++fmfndeb/gcs753XPnfL8X9jNnz/2dO6kqJEnt+a2NLkCSNBwDXJIaZYBLUqMMcElqlAEuSY3avp4H27VrV83Ozq7nISWpeUePHn2tqmZWjq9rgM/OzjI/P7+eh5Sk5iX5r9XGvYQiSY0ywCWpUQa4JDXKAJekRhngktQoA1ySGmWAS1KjDHBJapQBLkmNMsA1tWYPPLHRJUgTZYBLUqMMcElqVN8AT/K+JM8keS7JiSRf68YvSnI4yaluuXPy5UqSlgxyBv5r4ONV9VFgD7A3yXXAAeBIVV0JHOm2JUnrpG+AV8+vus0d3aOAW4BD3fgh4NZJFChJWt1A18CTbEtyDFgADlfV08AlVXUGoFtevMZr9yeZTzK/uLg4prKld3PWibaagQK8qt6uqj3AZcC1ST4y6AGq6mBVzVXV3MzMu/6ghCRpSOc1C6Wqfgn8ENgLnE2yG6BbLoy7OEnS2gaZhTKT5IPd+vuBTwAvAo8B+7rd9gGPTqhGSdIqBvmbmLuBQ0m20Qv8B6vq8SQ/Bh5Mcgfwc+CzE6xTkrRC3wCvqueBq1cZ/x/gxkkUJUnqzzsxJalRBria5JRByQCXpGYZ4JLUKANckhplgEtSowxwSWqUAS5JjTLAJalRBrgkNcoAl6RGGeCS1CgDXJIaZYBLUqMMcElqlAEuSY0ywCWpUQa4JDXKAJekRhngktQoA1ySGmWAa+r59zM1rQxwSWqUAS5JjTLAJalRfQM8yeVJnkxyMsmJJHd1419N8oskx7rHzZMvV5K0ZPsA+7wFfLGqnk1yIXA0yeHuuW9W1d9PrjxJ0lr6BnhVnQHOdOtvJjkJXDrpwiRJ53Ze18CTzAJXA093Q3cmeT7JfUl2rvGa/Unmk8wvLi6OVq00JKcSahoNHOBJLgAeAu6uqjeAbwEfAvbQO0P/xmqvq6qDVTVXVXMzMzOjVyxJAgYM8CQ76IX3/VX1MEBVna2qt6vqHeDbwLWTK1OStNIgs1AC3AucrKp7lo3vXrbbp4Hj4y9PkrSWQWahXA/cDryQ5Fg39mXgtiR7gAJOA5+fQH2SpDUMMgvlR0BWeer74y9HkjQo78SUpEYZ4JLUKANckhplgEtSowxwSWqUAS5JjTLAJalRBriacT5fSLW0r19ipWlmgEtSowxwSWqUAS5JjTLAJalRBrgkNcoAl6RGGeCS1CgDXJIaZYBLUqMMcElqlAEuSY0ywCWpUQa4JDXKAFfTZg884TcOassywCWpUQa4JDWqb4AnuTzJk0lOJjmR5K5u/KIkh5Oc6pY7J1+uJGnJIGfgbwFfrKrfA64DvpDkKuAAcKSqrgSOdNuSpHXSN8Cr6kxVPdutvwmcBC4FbgEOdbsdAm6dUI2SpFWc1zXwJLPA1cDTwCVVdQZ6IQ9cvMZr9ieZTzK/uLg4Yrna6pxxIv3GwAGe5ALgIeDuqnpj0NdV1cGqmququZmZmWFqlCStYqAAT7KDXnjfX1UPd8Nnk+zunt8NLEymREnSagaZhRLgXuBkVd2z7KnHgH3d+j7g0fGXJ0lay/YB9rkeuB14IcmxbuzLwNeBB5PcAfwc+OxEKpQkrapvgFfVj4Cs8fSN4y1HkjQo78SUpEYZ4GrCOKYPOgVR08YAl6RGGeCS1CgDXJIaZYBLUqMMcElqlAEuSY0a5E5MqRlOFdRW4hm4JDXKAJekRhngktQoA1ySGmWAS1KjDHBJapQBLkmNMsAlqVEGuCQ1ygCXpEYZ4JLUKANckhplgGtLWfqyK7/0StPAAJekRhngktSovgGe5L4kC0mOLxv7apJfJDnWPW6ebJmSpJUGOQP/DrB3lfFvVtWe7vH98ZYlSeqnb4BX1VPA6+tQiyTpPIxyDfzOJM93l1h2jq0iSdJAhg3wbwEfAvYAZ4BvrLVjkv1J5pPMLy4uDnk4SdJKQwV4VZ2tqrer6h3g28C159j3YFXNVdXczMzMsHVKklYYKsCT7F62+Wng+Fr7SpImY3u/HZI8ANwA7EryCvA3wA1J9gAFnAY+P7kSJUmr6RvgVXXbKsP3TqAWSdJ58E5MSWqUAa7m+EVUUo8BLkmNMsAlqVEGuCQ1ygCXpEYZ4JLUKANckhplgGvTm/S0QaclqlUGuCQ1ygCXpEYZ4JLUKANckhplgEtSowxwSWqUAa5NzSl+0toMcElqlAEuSY0ywCWpUQa4JDXKAJekRhngktQoA1ySGmWAS1KjDHBJalTfAE9yX5KFJMeXjV2U5HCSU91y52TLlCStNMgZ+HeAvSvGDgBHqupK4Ei3LUlaR30DvKqeAl5fMXwLcKhbPwTcOt6yJEn9DHsN/JKqOgPQLS9ea8ck+5PMJ5lfXFwc8nCSpJUm/iFmVR2sqrmqmpuZmZn04SRpyxg2wM8m2Q3QLRfGV5IkaRDDBvhjwL5ufR/w6HjKkSQNapBphA8APwY+nOSVJHcAXwc+meQU8MluW5K0jrb326GqblvjqRvHXIsk6Tx4J6YkNcoA16Yze+AJ/xamNAADXJIaZYBLUqMMcElqlAEuSY0ywCWpUQa4hDNf1CYDXJIaZYBLUqMMcElqlAEuSY0ywCWpUQa4JDXKANeGcMqeNDoDXJIaZYBLUqMMcElqlAEuSY0ywCWpUQa4JDXKANeG2sjphE5lVOsMcElqlAEuSY3aPsqLk5wG3gTeBt6qqrlxFCVJ6m+kAO/8UVW9NoafI0k6D15CkaRGjRrgBfwgydEk+1fbIcn+JPNJ5hcXF0c8nKbVZpkRslnqkAYxaoBfX1XXADcBX0jysZU7VNXBqpqrqrmZmZkRDydJWjJSgFfVq91yAXgEuHYcRUmS+hs6wJN8IMmFS+vAp4Dj4ypMknRuo8xCuQR4JMnSz/luVf3LWKqSJPU1dIBX1c+Aj46xFknSeXAaoSQ1ygDXhlk5ZW8zTeHbTLVIazHAJalRBrgkNcoAl6RGGeCS1CgDXJIaZYBLUqMMcGkNTiXUZmeAS1KjDHBJapQBLkmNMsAlqVEGuCQ1ygDXutvsszs285dsScsZ4JLUKANckhplgEtSowxwSWqUAS5JjTLAJalRBrjGbvm0u6X1lctpMo09qQ0GuCQ1ygCXpEYZ4JLUqJECPMneJC8leTnJgXEVJUnqb+gAT7IN+EfgJuAq4LYkV42rMEnSuY1yBn4t8HJV/ayq/hf4HnDLeMqSJPWTqhruhclngL1V9Rfd9u3AH1TVnSv22w/s7zY/DLw0ZK27gNeGfG2r7HlrsOetYZSef6eqZlYObh+hmKwy9q7fBlV1EDg4wnF6B0vmq2pu1J/TEnveGux5a5hEz6NcQnkFuHzZ9mXAq6OVI0ka1CgB/u/AlUmuSPIe4HPAY+MpS5LUz9CXUKrqrSR3Av8KbAPuq6oTY6vs3Ua+DNMge94a7HlrGHvPQ3+IKUnaWN6JKUmNMsAlqVFNBPi03rKf5L4kC0mOLxu7KMnhJKe65c5lz32pew9eSvLHG1P18JJcnuTJJCeTnEhyVzc+zT2/L8kzSZ7rev5aNz61PS9Jsi3JT5I83m1Pdc9JTid5IcmxJPPd2GR7rqpN/aD3AelPgd8F3gM8B1y10XWNqbePAdcAx5eN/R1woFs/APxtt35V1/t7gSu692TbRvdwnv3uBq7p1i8E/qPra5p7DnBBt74DeBq4bpp7Xtb7XwLfBR7vtqe6Z+A0sGvF2ER7buEMfGpv2a+qp4DXVwzfAhzq1g8Bty4b/15V/bqq/hN4md5704yqOlNVz3brbwIngUuZ7p6rqn7Vbe7oHsUU9wyQ5DLgT4B/WjY81T2vYaI9txDglwL/vWz7lW5sWl1SVWegF3jAxd34VL0PSWaBq+mdkU51z92lhGPAAnC4qqa+Z+AfgL8C3lk2Nu09F/CDJEe7rxCBCfc8yq3062WgW/a3gKl5H5JcADwE3F1VbySrtdbbdZWx5nquqreBPUk+CDyS5CPn2L35npP8KbBQVUeT3DDIS1YZa6rnzvVV9WqSi4HDSV48x75j6bmFM/Ctdsv+2SS7AbrlQjc+Fe9Dkh30wvv+qnq4G57qnpdU1S+BHwJ7me6erwf+LMlpepc8P57kn5nunqmqV7vlAvAIvUsiE+25hQDfarfsPwbs69b3AY8uG/9ckvcmuQK4EnhmA+obWnqn2vcCJ6vqnmVPTXPPM92ZN0neD3wCeJEp7rmqvlRVl1XVLL1/r/9WVX/OFPec5ANJLlxaBz4FHGfSPW/0J7cDfrp7M70ZCz8FvrLR9YyxrweAM8D/0fuNfAfw28AR4FS3vGjZ/l/p3oOXgJs2uv4h+v1Dev9NfB441j1unvKefx/4SdfzceCvu/Gp7XlF/zfwm1koU9szvVlyz3WPE0s5NemevZVekhrVwiUUSdIqDHBJapQBLkmNMsAlqVEGuCQ1ygCXpEYZ4JLUqP8HvlS2LghG5EEAAAAASUVORK5CYII=\n",
+      "text/plain": [
+       "<Figure size 432x288 with 1 Axes>"
+      ]
+     },
+     "metadata": {
+      "needs_background": "light"
+     },
+     "output_type": "display_data"
+    }
+   ],
+   "source": [
+    "slots = np.zeros(500)\n",
+    "x = np.arange(500)\n",
+    "for exp in range(1000):\n",
+    "    s = 0\n",
+    "    for i in range(50):\n",
+    "        r = np.random.randint(1, 9)\n",
+    "        #r = [1, 8][np.random.randint(0, 2)]\n",
+    "        s += r\n",
+    "    slots[s] += 1\n",
+    "\n",
+    "plt.bar(x, slots)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Problem 2\n",
+    "\n",
+    "The distribution function is given by:\n",
+    "$$\n",
+    "f(x)_\\theta = \\frac x {\\theta^2} e ^ {- \\frac x \\theta }\n",
+    "$$\n",
+    "\n",
+    "Then we can find the likelyhood function $f(x_1,...,x_n | \\theta)$.\n",
+    "\n",
+    "$$\n",
+    "f(x_1,...,x_n | \\theta) = f(x_1)_\\theta \\cdot ... \\cdot f(x_n)_\\theta = \\frac {\\prod_i x_i} {\\theta^{2n}} \\cdot \\exp \\left( -\\frac {\\sum_i} \\theta \\right)\n",
+    "$$\n",
+    "\n",
+    "Her kan man tage log på hver side\n",
+    "\n",
+    "$$\n",
+    "\\log f(x_1,...,x_n | \\theta) = \\log \\left( \\frac {\\prod_i x_i} {\\theta^{2n}} \\right) - \\frac 1 \\theta \\sum_i x_i\n",
+    "$$\n",
+    "\n",
+    "Og diff i forhold til $\\theta$.\n",
+    "\n",
+    "$$\n",
+    "\\frac d {d\\theta} \\log f(x_1,...,x_n | \\theta) = -\\frac {2n} \\theta + \\frac 1 {\\theta^2} \\sum_i x_i\n",
+    "$$\n",
+    "\n",
+    "Og man kan løse for $0$.\n",
+    "\n",
+    "$$\n",
+    "\\theta = \\frac {\\sum_i x_i} {2n}\n",
+    "$$\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 34,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "Estimated theta = 1.2944444444444445\n"
+     ]
+    }
+   ],
+   "source": [
+    "# Part B\n",
+    "# Her indsættes samples\n",
+    "\n",
+    "samples = np.array([3.2, 1.4, 6.5, 2.2, 1.8, 2.6, 3.9, 0.5, 1.2])\n",
+    "N = len(samples)\n",
+    "\n",
+    "sum_samples = np.sum(samples)\n",
+    "\n",
+    "theta_est = sum_samples / (2 * N)\n",
+    "print(f\"Estimated theta = {theta_est}\")"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 39,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "argmax: 1.2962962962962963\n"
+     ]
+    },
+    {
+     "data": {
+      "text/plain": [
+       "[<matplotlib.lines.Line2D at 0x7faf4d186640>]"
+      ]
+     },
+     "execution_count": 39,
+     "metadata": {},
+     "output_type": "execute_result"
+    },
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEDCAYAAAAlRP8qAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAAj/ElEQVR4nO3deXScd33v8fd3NFqsxZK12JYtyWtsx3glirOTACELJKSl0JPALYUbbnALHLhcLnB7Tttzy+npdtsLtKQhDbRwaUlDSblprklCIAtJSLCdeInj2JZ3LbY2S9a+ze/+MTPORJbssT0zz6LP6xxFmplnnueryP74p9/zW8w5h4iIBF/E6wJERCQzFOgiIiGhQBcRCQkFuohISCjQRURCQoEuIhISnga6mX3XzNrN7PUMnOvdZrYj5WPYzH4jA2WKiASCeTkO3czeBfQD33fOrcngeSuBJqDOOTeYqfOKiPiZpy1059zzQHfqc2a2zMyeMLPtZvZLM1t1Eaf+MPBThbmIzCR+7EN/EPicc+4K4EvA/RdxjruBH2a0KhERn4t6XUAqMysFrgV+ZGbJpwsTr30I+JMp3tbinLs15Ry1wFrgyexWKyLiL74KdOK/MfQ45zZMfsE59yjwaBrn+G3g351zYxmuTUTE13zV5eKcOw0cNrOPAFjc+gs8zT2ou0VEZiCvhy3+EPgVsNLMms3sXuBjwL1mthPYA9x1AedbDNQDz2WhXBERX/N02KKIiGSOr7pcRETk4nl2U7S6utotXrzYq8uLiATS9u3bO51zNVO95lmgL168mG3btnl1eRGRQDKzo9O9pi4XEZGQUKCLiISEAl1EJCQU6CIiIaFAFxEJCQW6iEhIKNBFRELCb6stSpZ09I3w+K5WDPjAugXUlBV6XZKIZJgCfQZ4vaWXjz30Cr1D8RWF//fTB/jO7zbSuLjS48pEJJPU5RJy3QOjfOIft1JaGOXJL7yLJ7/wLipLCrj3e9s42jXgdXkikkHnDXQz+66ZtZvZ69O8bmb2TTNrMrNdZvbOzJcpF+uvntxHz+AoD/1uIyvnl7Fyfhnf/8+biMUcf/Dvu9FqmyLhkU4L/Z+A287x+u3AZYmP+4C/v/SyJBOOdQ3y8NZjfPyaxVxeO/vM8/WVxXz5tpW82NTFU2+c9LBCEcmk8wa6c+55oPsch9wFfN/FvQxUJPb1FI/900tHyDPj0zcuPeu1ezY1sKS6hK8/fYBYTK10kTDIRB/6QuB4yuPmxHNnMbP7zGybmW3r6OjIwKVlOsNjE/xo+3E+sK6WebOLzno9mhfhM+9ezt6207x4sNODCkUk0zIR6DbFc1M2+ZxzDzrnGp1zjTU1Uy7nKxnyywOd9A2P86F31k17zJ3ra6ksKeAHL0+7GqeIBEgmAr2Z+D6eSXVAawbOK5dgy+42KorzuXZZ1bTHFEbz+EhjHU/vbaetdyiH1YlINmQi0B8DPp4Y7XI10Ouca8vAeeUiDY9N8LM3TnLr6vnk5537R/zRTQ1MxBw/eU3/BosEXTrDFn8I/ApYaWbNZnavmW02s82JQ7YAh4Am4B+A389atZKWbUdO0T8yzq1r5p332EVVJWxsqOA/dirQRYLuvDNFnXP3nOd1B3wmYxXJJXuhqZP8POOqJdN3t6S6c90C/uTxN2hq72f53NIsVyci2aKZoiH0YlMnG+vnUFKY3soOd6yrxQwe36VWukiQKdBDpmdwlNdbe7lueXXa75k7u4grGubw9F5NMhIJMgV6yLx8qAvn4Lrl6XW3JL3n8rm83nKaE73DWapMRLJNgR4yrx7roSAaYV1dxQW97+bL4zdQf/FmexaqEpFcUKCHzGvHTrF2YTkF0Qv70V42t5S6ObP4xZvqdhEJKgV6iIxNxNjV3MuG+ooLfq+ZcfPl83ihqZPhsYnMFyciWadAD5E32/oYGY+xsaHiot5/44oahsdibDtyKrOFiUhOKNBD5LXj8SDe2DDnot6/aUkl0YhpsS6RgFKgh8jO471UlxayoPzs1RXTUVIYZWNDBS82KdBFgkiBHiJ7206zesFszKZaADM91y6rZndLL72DYxmsTERyQYEeEmMTMZra+7m8tuySznPd8mqcg18d6spQZSKSKwr0kDjY0c/oRIzVKVvNXYwN9RXMys/jJfWjiwSOAj0k9radBnjb3qEXoyAa4aqllbygfnSRwFGgh8Tetj4KohGWVpdc8rmuXlrFoY4BOvpGMlCZiOSKAj0k9radZsW8UqLn2dAiHVcurgRg+9Fz7Q0uIn6jQA+JvW2nuXz+pXW3JK1dWE5hNMKvD2uCkUiQKNBDoKNvhM7+UVZdYv95UkE0wob6CraphS4SKAr0EGhq7wdgxbzM7Ta0aUkle1pPMzAynrFzikh2KdBD4GBHPNCX1WQu0BsXVzIRc7x6TN0uIkGhQA+BpvZ+igvyqL3IKf9TeWdDBRGDrVqoSyQwFOghcLCjn2U1pZc05X+ysqJ8Vi+YzdbD6kcXCQoFeggc6hhgWc2ljz+frHFRJa8dP8XYRCzj5xaRzFOgB9zg6DgtPUMZ7T9Palw8h+GxGG+29WX83CKSeQr0gDvUMQDAsrmZD/TkuurJddZFxN8U6AGXHOGyPAuBvqC8iLllhbx2rCfj5xaRzFOgB9zB9n4iBouqijN+bjNjY0MFr2nookggKNAD7mDnAPWVxRRG87Jy/o0NczjSNUj3wGhWzi8imaNAD7ijXQMsqsr8CJekjfUVAOxQP7qI7ynQA8w5x9GuQRZVZr67JWltXTl5EVM/ukgAKNADrGdwjL7h8az0nycVF0RZNb9MgS4SAGkFupndZmb7zKzJzL46xevlZvYfZrbTzPaY2SczX6pMdrR7ECCrXS4AGxsq2HG8h4mYy+p1ROTSnDfQzSwP+BZwO7AauMfMVk867DPAG8659cBNwF+bWUGGa5VJjnbFx6Bns4UOsLF+Dv0j42eGSIqIP6XTQt8ENDnnDjnnRoGHgbsmHeOAMosvJlIKdANadzXLjnbFW+gNWexDh3gLHdDwRRGfSyfQFwLHUx43J55L9XfA5UArsBv4vHPurAVAzOw+M9tmZts6OjousmRJOto1yPzZRRTlZ2fIYtKS6hLKZ+WrH13E59IJ9KmW8JvcmXorsANYAGwA/s7Mzto+xzn3oHOu0TnXWFNTc4GlymRHuwZoyHJ3C6ROMOrJ+rVE5OKlE+jNQH3K4zriLfFUnwQedXFNwGFgVWZKlOkc7R5kcQ4CHeL96Pvb++gbHsvJ9UTkwqUT6FuBy8xsSeJG593AY5OOOQa8F8DM5gErgUOZLFTebnB0nI6+kayPcEna0FCBc7C7uTcn1xORC3feQHfOjQOfBZ4E9gKPOOf2mNlmM9ucOOxrwLVmthv4OfAV51xntooWONadmxuiSRvqKgB47XhPTq4nIhcums5BzrktwJZJzz2Q8nUrcEtmS5NzOdKZHIOem0AvL85naXUJOxToIr6lmaIBdTw5qagyN10uABvq4xOMnNMEIxE/UqAHVEvPEGWFUcqL83N2zQ0NFXT0jdDaO5yza4pI+hToAdXaM8SCilk5veaG5MqLGr4o4ksK9IBq7R2itqIop9dcNX82BdGIltIV8SkFekC19gznvIVeEI3wjgWzdWNUxKcU6AE0NDpB98AoC3Mc6BDvdtnd0svYxFkrO4iIxxToAdTWOwRAbXluu1wgHujDYzH2nejL+bVF5NwU6AHU2hMfZZLrLheILwEAsLO5J+fXFpFzU6AHUGtPvIXuRZdLfeUsKksKNNJFxIcU6AHU2juEGcybnfsuFzM7M8FIRPxFgR5ArT1D1JQWUhD15se3ob6Cpo5+rbwo4jMK9ADyYshiqg318ZUXd2nlRRFfUaAHUGvPkCf950nrkzNG1e0i4isK9IBxzsVniXowZDGpfFY+S2tKtIORiM8o0APm1OAYw2MxT7tcIL4+ulZeFPEXBXrAJIcseh7oDRV09o/QkqhHRLynQA+YFg/HoKfaoH50Ed9RoAdMWyLQc73S4mRnVl5UP7qIbyjQA6a1d5iCaISqkgJP6yiIRlizYLaWABDxEQV6wLQkhiyamdelsKF+jlZeFPERBXrAtPV4O2Qx1YYGrbwo4icK9IDxepZoqo26MSriKwr0ABmbiHGyzz+BXjdnFlUlBQp0EZ9QoAfIid5hnIOFHo9wSdLKiyL+okAPkLbe+MYWteX+aKFDfDz6wY5+TmvlRRHPKdADxC+zRFNtaEisvHhcKy+KeE2BHiAtZwLdH10uAOvqKgDYcfyUt4WIiAI9SFp7hqgozqe4IOp1KWckV15UP7qI9xToAdLWO8wCH/WfJyVvjGrlRRFvpRXoZnabme0zsyYz++o0x9xkZjvMbI+ZPZfZMgXiLXQ/9Z8nbayvoLN/VCsvinjsvIFuZnnAt4DbgdXAPWa2etIxFcD9wAedc+8APpL5UiU+7d8//edJG+rnAJpgJOK1dFrom4Am59wh59wo8DBw16RjPgo86pw7BuCca89smdI3PEbf8Di1Pmyhr6oto1ArL4p4Lp1AXwgcT3ncnHgu1Qpgjpk9a2bbzezjU53IzO4zs21mtq2jo+PiKp6hkmPQ/djlkp8XYc3CcrXQRTyWTqBPtazf5LtfUeAK4APArcAfmtmKs97k3IPOuUbnXGNNTc0FFzuTvbWxhf+6XCB+Y1QrL4p4K51AbwbqUx7XAa1THPOEc27AOdcJPA+sz0yJAm9NKvLTLNFUG+orGBnXyosiXkon0LcCl5nZEjMrAO4GHpt0zP8FbjCzqJkVA1cBezNb6szW1jNMXsSYW1bodSlTSm5J95q6XUQ8c95Ad86NA58FniQe0o845/aY2WYz25w4Zi/wBLAL+DXwkHPu9eyVPfO09gwxf3YR0Tx/Th2omzOL6tIC3RgV8VBaUw6dc1uALZOee2DS478C/ipzpUmqlp4hX035n+ytlRe1BICIV/zZ3JOztPYO+bb/PGl9XQUHOwboHdLKiyJeUKAHQCzmONHrn40tprOhoQKAXdo4WsQTCvQA6OwfYWzC+XbIYlJy5cWdujEq4gkFegC0+HzIYlL5rHyWaeVFEc8o0APAz7NEJ9tQP0crL4p4RIEeAK1nZokGINAb4isvNp/SyosiuaZAD4CWniFKCvKYPcs/G1tMZ2NigpG6XURyT4EeAK09Q9RWzMJsqmV1/GXl/DJm5eex/ajGo4vkmgI9ANoCMGQxKT8vwsaGCrYe6fa6FJEZR4EeAK0+3dhiOlcurmRv22n6hjXBSCSXFOg+Nzw2QWf/qC/3Ep3OpiWVxBzqdhHJMQW6zyWHLPpxp6LpbGyoIBoxfn1Y3S4iuaRA97m2xJBFPy/MNVlxQZQ1C8vVjy6SYwp0n2sJ0Bj0VJuWVLLzeC/DYxNelyIyYyjQfa61J97lMr88OC10iN8YHZ2Isau51+tSRGYMBbrPtfUOUV1aSGE0z+tSLkjjojkA6nYRySEFus+1BGzIYtKckgJWzCvVjVGRHFKg+1xrz1BgJhVNduXiSl49eoqJmBbqEskFBbqPOedo7Rn2/bK509m0pJK+kXH2tp32uhSRGUGB7mO9Q2MMjU0EashiqisXVwKo20UkRxToPhbUIYtJCypm0VBZzMuHurwuRWRGUKD7WHLIYpBmiU527bIqXj7UpX50kRxQoPtYawBniU52zbIqTg+Ps6dV49FFsk2B7mOtvUMU5EWoLin0upSLdu2yagBeOqhuF5FsU6D7WGvPMLUVRUQi/t/YYjo1ZYWsmFfKi02dXpciEnoKdB9r7RmiNmBT/qdy7bJqth7pZnQ85nUpIqGmQPextgBPKkp17bIqhsdi2mdUJMsU6D41PhHjxOnhwA5ZTHXV0ioihrpdRLJMge5TJ/tGiDkCO0s0VfmsfNYsLOdXujEqklUKdJ9qOZWYVDQn+IEO8eGLrx0/xeDouNeliIRWWoFuZreZ2T4zazKzr57juCvNbMLMPpy5Emem1oDPEp3sumXVjE04XtEyACJZc95AN7M84FvA7cBq4B4zWz3NcX8BPJnpImeilhBMKkq1aUklhdEIz+/v8LoUkdBKp4W+CWhyzh1yzo0CDwN3TXHc54AfA+0ZrG/GaukZorKkgOKCqNelZERRfh5XL63iuX0KdJFsSSfQFwLHUx43J547w8wWAr8JPHCuE5nZfWa2zcy2dXToL/a5tJwaCk3rPOmmlTUc6hzgWNeg16WIhFI6gT7VNMXJKy19HfiKc+6cOwI75x50zjU65xpramrSLHFmiu9UFI7+86QbV8R/5s/t1y9xItmQTqA3A/Upj+uA1knHNAIPm9kR4MPA/Wb2G5kocCaKb2wRjklFqZZUl9BQWcyz6nYRyYp0An0rcJmZLTGzAuBu4LHUA5xzS5xzi51zi4F/A37fOfeTTBc7U/QMjjE4OhG6FrqZceOKGl462MXI+Dl/mRORi3DeQHfOjQOfJT56ZS/wiHNuj5ltNrPN2S5wJkqOcKkLyRj0VDetrGFobIKth095XYpI6KQ1hMI5twXYMum5KW+AOuc+cellzWxvDVkMX6Bfs6yKgrwIz+1v5/rLqr0uRyRUNFPUh8I2qShVcUGUTUsq1Y8ukgUKdB9qOTVEUX6EypICr0vJiptW1nCgvZ/j3Rq+KJJJCnQfau2Nj3AxC+7GFufyvtXzAHjqjZMeVyISLgp0H2o5Fb4x6KkWVZWwan4ZT+054XUpIqGiQPehlp5wrIN+LresnsfWI910D4x6XYpIaCjQfWZ4bILO/pHQB/r7Vs8n5uDne9XtIpIpCnSfaesdBsI5ZDHVmoWzqS0v4mfqRxfJGAW6z4RtY4vpmBm3rJ7H8wc6GBrVrFGRTFCg+0xLT3woX9i7XABuecd8hsdiPH9AY9JFMkGB7jMtPcOYwfzycC2dO5VNSyqpKM7np7vbvC5FJBQU6D7T3D3IgvJZ5OeF/0eTnxfh9jXz+dkbJ9XtIpIB4U+NgDnWPUh9Zfi7W5LuXLeAgdEJntmnNdJFLpUC3WeOdQ/SUFnsdRk5c9XSKqpLC/mPnZOX2BeRC6VA95Gh0Qna+0ZmVKDnRYw71tXyizfb6Rse87ockUBToPtI86n4CJf6GRToAHeur2VkPMbTmmQkckkU6D5yrHtmBvrG+jksrJjFYzvU7SJyKRToPpIM9JnU5QIQiRh3rK/l+QOddPSNeF2OSGAp0H3kePcQxQV5VIV0HfRz+cgVdUzEHD95rcXrUkQCS4HuI8kRLmFdB/1cls8tY2NDBT/afhznnNfliASSAt1HjncPzrj+81QfuaKe/Sf72dXc63UpIoGkQPcJ59yMG4M+2R3raynKj/DItuNelyISSAp0n+jsH2VobGJGB/rsonxuX1PLYztbGR7TUgAiF0qB7hMzdYTLZB+5oo6+4XF++roW7BK5UAp0nzjWPQAwo9ZxmcrVS6tYWlPC91466nUpIoGjQPeJw52DRGzmTSqaLBIxPn71InYc72Hn8R6vyxEJFAW6TxzuHKBuTjGF0TyvS/Hcb11RR0lBHt/71RGvSxEJFAW6Txzq6GdJdYnXZfhCWVE+v3VFHY/vbKOrXzNHRdKlQPcB5xyHOwdYWqNAT/r4NYsYnYjx8FYNYRRJlwLdB9r7RhgcnWCpWuhnLJ9bxvXLq/neS0cYGdcQRpF0pBXoZnabme0zsyYz++oUr3/MzHYlPl4ys/WZLzW8DnXER7gsqS71uBJ/+fSNS2nvG+HH27W+i0g6zhvoZpYHfAu4HVgN3GNmqycddhi40Tm3Dvga8GCmCw2zw52JQFeXy9tcv7yatQvL+fbzBxmfiHldjojvpdNC3wQ0OecOOedGgYeBu1IPcM695Jw7lXj4MlCX2TLD7VBHP0X5EWpnF3ldiq+YGZ959zKOdg2y5fUTXpcj4nvpBPpCIPXOVHPiuencC/x0qhfM7D4z22Zm2zo6OtKvMuQOdw6wuKqESGTmrbJ4Presns+ymhLuf6ZJqzCKnEc6gT5Vykz5N8vM3k080L8y1evOuQedc43Oucaampr0qwy5Q50DGrI4jUjE2HzjMt480cdTb2iLOpFzSSfQm4H6lMd1wFl7hZnZOuAh4C7nXFdmygu/4bEJjnYNcNm8Mq9L8a3f3LiQpdUl/K8n9zERUytdZDrpBPpW4DIzW2JmBcDdwGOpB5hZA/Ao8DvOuf2ZLzO8Dnb0E3OwUoE+rWhehC/dupID7f08+mqz1+WI+NZ5A905Nw58FngS2As84pzbY2abzWxz4rA/AqqA+81sh5lty1rFIbP/ZB8AK+dryOK53L5mPmsXlvP1pw9oaV2RaaQ1Dt05t8U5t8I5t8w596eJ5x5wzj2Q+PpTzrk5zrkNiY/GbBYdJvtO9FOQF2FRlfrQz8XM+Mptq2jpGeIHL2slRpGpaKaox/af7GNpTQn5efpRnM/1l1Vzw2XVfOPnB+jo0xovIpMpRTy270QfK+er/zxdf3znOxgem+AvnnjT61JEfEeB7qG+4TFaeoZYoRuiaVs+t5R7r1/Kv21vZvvRbq/LEfEVBbqHDrT3AyjQL9Dn3rOc2vIi/vAne7QkgEgKBbqH9p1IjHBRoF+QksIof3THat5oO80Dzx30uhwR31Cge+j1ll7KiqIzfh/Ri3H72lruWFfLN35+gDdaT3tdjogvKNA9tLull7ULyzHTGi4X42t3raF8VgFffGQHo+PqehFRoHtkdDzGm219rF1Y7nUpgTWnpIA/+9Ba3jzRx18/tc/rckQ8p0D3yP6TfYxOxFijQL8k71s9j49e1cC3nz/EU3u0xK7MbAp0j+xu6QVQCz0D/uiO1axdWM5/+9FOjnYNeF2OiGcU6B7Z1dxDWVGURVXFXpcSeEX5edz/sXcSMWPzD15lYGTc65JEPKFA98j2o6d4Z8Mc3RDNkPrKYr5x9wb2n+zjs//yqsany4ykQPfAqYFR9p/s58rFc7wuJVRuWjmXr921hmf2dfDHj+3RDkcy40S9LmAm2n40vv3qlYsrPa4kfD56VQPNpwa5/9mDVJUU8MVbVnpdkkjOKNA9sPVoN/l5xvr6Cq9LCaUv3bKSrv5RvvmLJiIR4ws3r/C6JJGcUKB7YOvhbtYuLKcoP8/rUkIpEjH+7ENriTnH158+QMzBf735Mt2vkNBToOdY79AYO5t7+b0bl3ldSqhFIsZf/NY6AL6ZWD/9a3e9g6jWnZcQU6Dn2AsHOpmIOW5aWeN1KaEXiRh/+eF11JQVcv+zBzl5ephv3rOR0kL9sZdwUnMlx57d187soigb1H+eE2bGl29bxZ/+5hqe3dfOB//uhTOrXIqEjQI9h5xzPLe/gxtW1OhX/xz72FWL+OdPXc3poXHu+tYL/Nv2Zg1rlNBRquTQq8d6aO8b4b2r5npdyox0zbIqtnz+ejbUV/ClH+3k0/9nO+2nh70uSyRjFOg59NiOFgqjEd63ep7XpcxYc8uK+OdPXc0fvH8Vz+3v4Oa/eY5/3XqMiZha6xJ8CvQcGZ+I8fiuNm6+fB5lRflelzOj5UWM+961jJ9+/gZWzZ/NV368mzv/9gVeaur0ujSRS6JAz5HnD3TQNTDKnesXeF2KJCytKeVfP30137xnI71DY3z0oVf4ne+8wsuHutS/LoGk8Vs58o8vHmFuWSHvUf+5r5gZH1y/gFtWz+N7Lx3hH355iLsffJnGRXP41A1Lee/lc8nXDWwJCAV6Drx54jS/PNDJl25ZQUFU4eBHRfl5fPrGZfzutYt5ZNtxvv3cITb/YDtzywr57cZ6fruxngYtdSw+Z179atnY2Oi2bdvmybVz7VPf28Yrh7v45ZffTUVxgdflSBrGJ2I8t7+Df3nlGM/sayfm4puR3L52PrevqWVJdYnXJcoMZWbbnXONU72mFnqWPbe/g6f3nuS/37pSYR4g0bwI7718Hu+9fB6tPUM8vquVLbtP8JdP7OMvn9jH0poSrltWzXXLq7hmaTXlxbrRLd5TCz2LuvpHuPNvX6C4MMrjn7tei3GFQEvPEE++foJfHujglcPdDI5OYAYr5paxrq6c9fUVrK+rYOX8MnWvSVacq4WuQM+SvuExPvmPW9nd0suPNl/DuroKr0uSDBsdj7GzuYcXmzrZcbyHXc29dA+MAvGhkYsqi1k2t5Tlc0tZXlPK4uoSFlbMYm5ZIZGIVn6Ui3PJXS5mdhvwDSAPeMg59+eTXrfE6+8HBoFPOOdevaSqA+yN1tN88ZEdNLX38817NirMQ6ogGuHKxZVnNipxztF8aoidzT282dZHU3s/TR39PPNmO+MpE5fy84z55UUsrJjFgvJZVJUWUFVaSGVJAVUl8a+rSgqoKM6npCCq8Je0nTfQzSwP+BbwPqAZ2Gpmjznn3kg57HbgssTHVcDfJz7PCMNjE7T0DLHjWA9P7DnB03tPUlVSyHc+cSU3rtCqijOFmVFfWUx9ZTF3rHvr+bGJGEe7BjnePUhLzxAtPUO0Jj5eOdxNZ/8II+NT74FqBqUFUUqLopQWvvV5dlE+xQV5FOXnUZQfoTCaR2E0QmF+hKL8xNeJ54ry88jPi5AXMaJ5Rl7EyI+8/XE0kvwceetxXvxzxAwz4p/hzGOtL+8/6bTQNwFNzrlDAGb2MHAXkBrodwHfd/H+m5fNrMLMap1zbZku+Nl97Xzt8filXeI/ybaPcw4HOAcu8axz8Q/Odwypx7kzz531nsRBqdfuT9llvrq0gN+7cRn/5YalzCnRTVCB/LxIvNtlbumUrzvnGBydoHtglK6BUbr6R+gaGKVncJT+kQn6h8fpHxmjf2ScvuH4R1vvMAMj44yMxxgem2B4bAIvVi+YHPQYRAwMi39Ohj/x5YyN+HPxXzqSx8SPn3zetz0+67rn/8fkrHOcdc5zXzOd6571ljTrvvvKej51w9KzL3iJ0gn0hcDxlMfNnN36nuqYhcDbAt3M7gPuA2hoaLjQWgEoK8pn1fzZZ/5PJf+AJP/HJf/wpD5H4g+Mve09b/1AzZI/zPMcw9t/oMnXK0vyWThnFivnzWbV/DL9iiwXxMwoKYxSUhilvvLix7qPT8QYHo8xMjZxJuhHxmOMjMcYHY8xEXOMx2KMxxwTEy7+OfFc/LN76/NE/LjkcxD/hyeWaCDFzjRy3NsexxItodiZ5+MNp2RjKPVxLNEqik365cTx9n+ZJt/mm+rfrbOPOfebJp9jqnuJZx9zYec464wpT1SXFp51vUxIJ9CnSqfJtaZzDM65B4EHIX5TNI1rn+WKRXO4YtGci3mrSKhF8yKU5kW0gccMls64qmagPuVxHdB6EceIiEgWpRPoW4HLzGyJmRUAdwOPTTrmMeDjFnc10JuN/nMREZneeX83c86Nm9lngSeJD1v8rnNuj5ltTrz+ALCF+JDFJuLDFj+ZvZJFRGQqaXW2Oee2EA/t1OceSPnaAZ/JbGkiInIhNDdZRCQkFOgiIiGhQBcRCQkFuohISHi22qKZdQBHL/Lt1cBM29FX3/PMoO95ZriU73mRc27KRaI8C/RLYWbbpls+Mqz0Pc8M+p5nhmx9z+pyEREJCQW6iEhIBDXQH/S6AA/oe54Z9D3PDFn5ngPZhy4iImcLagtdREQmUaCLiIRE4ALdzG4zs31m1mRmX/W6nmwzs++aWbuZve51LbliZvVm9oyZ7TWzPWb2ea9ryiYzKzKzX5vZzsT3+z+9rilXzCzPzF4zs8e9riUXzOyIme02sx1mti3j5w9SH3piw+r9pGxYDdwzacPqUDGzdwH9xPdsXeN1PblgZrVArXPuVTMrA7YDvxHWn7PF9zUscc71m1k+8ALweefcyx6XlnVm9kWgEZjtnLvD63qyzcyOAI3OuaxMpApaC/3MhtXOuVEguWF1aDnnnge6va4jl5xzbc65VxNf9wF7ie9RG0ourj/xMD/xEZyW1kUyszrgA8BDXtcSFkEL9Ok2o5aQMrPFwEbgFY9LyapE18MOoB34mXMu1N9vwteBLwOx8xwXJg54ysy2m9l9mT550AI9rc2oJRzMrBT4MfAF59xpr+vJJufchHNuA/H9eDeZWai718zsDqDdObfd61py7Drn3DuB24HPJLpUMyZoga7NqGeIRF/yj4F/ds496nU9ueKc6wGeBW7ztpKsuw74YKJP+WHgPWb2A29Lyj7nXGviczvw78S7kTMmaIGezobVEnCJm4TfAfY65/7G63qyzcxqzKwi8fUs4GbgTU+LyjLn3P9wztU55xYT/3v8C+fcf/K4rKwys5LETX7MrAS4Bcjo6LVABbpzbhxIbli9F3jEObfH26qyy8x+CPwKWGlmzWZ2r9c15cB1wO8Qb7XtSHy83+uisqgWeMbMdhFvtPzMOTcjhvHNMPOAF8xsJ/Br4P85557I5AUCNWxRRESmF6gWuoiITE+BLiISEgp0EZGQUKCLiISEAl1EJCQU6CIiIaFAFxEJif8P0+PrXvkm6SQAAAAASUVORK5CYII=\n",
+      "text/plain": [
+       "<Figure size 432x288 with 1 Axes>"
+      ]
+     },
+     "metadata": {
+      "needs_background": "light"
+     },
+     "output_type": "display_data"
+    }
+   ],
+   "source": [
+    "# Prøver lige at løse den programmerbart også\n",
+    "thetas = np.linspace(0, 5, 1000)[1:]\n",
+    "res = np.empty(thetas.shape)\n",
+    "\n",
+    "for i, theta in enumerate(thetas):\n",
+    "    f = lambda x: ((x / (theta**2)) * np.exp(-x / theta))\n",
+    "    score = 1\n",
+    "    for sample in samples:\n",
+    "        score *= f(sample)\n",
+    "    \n",
+    "    res[i] = score\n",
+    "\n",
+    "maxtheta = thetas[np.argmax(res)]\n",
+    "print(f\"argmax: {maxtheta}\")\n",
+    "plt.plot(thetas, res)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Problem 3\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 60,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "radial = np.array([5.2, 5.7, 7.6, 8.0, 7.7, 5.5, 6.7, 7.0, 8.4, 5.9])\n",
+    "\n",
+    "belted = np.array([5.1, 5.9, 7.2, 7.9, 7.8, 5.4, 6.7, 6.8, 7.9, 5.7])\n",
+    "belted_mean = np.mean(belted)\n",
+    "N = len(belted)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Part A\n",
+    "\n",
+    "Here we don't know the varience, and will therefore have to use the t distribution method.\n",
+    "The interval is given by:\n",
+    "$$\n",
+    "\\bar X \\pm t_{\\alpha/2,n-1} \\frac S {\\sqrt{n}}\n",
+    "$$\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 59,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "Then the interval becomes [6.174311863003706, 7.105688136996291]\n"
+     ]
+    }
+   ],
+   "source": [
+    "alpha = 0.1\n",
+    "# Lets start by calculating the t value\n",
+    "tval = stats.t.ppf(1 - alpha, N-1)\n",
+    "\n",
+    "S = np.sqrt(np.sum((belted - belted_mean)**2) / (N - 1))\n",
+    "\n",
+    "diff = tval * S / np.sqrt(N)\n",
+    "print(f\"Then the interval becomes [{belted_mean - diff}, {belted_mean + diff}]\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Part B\n",
+    "\n",
+    "We can take the difference of the two and say that this should be less than or equal to $0$.\n",
+    "Thus the $H_0$ is $(\\mu_{belted} - \\mu_{radial}) \\leq 0$.\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 79,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "[ 0.1 -0.2  0.4  0.1 -0.1  0.1  0.   0.2  0.5  0.2] 0.12999999999999998\n",
+      "ts: 1.9475670608117976\n",
+      "p-value: 0.04164574347828587\n",
+      "H_0 is accepted, thus radial does not make it better\n"
+     ]
+    }
+   ],
+   "source": [
+    "mu_0 = 0\n",
+    "alpha = 0.05\n",
+    "diff = radial - belted\n",
+    "diff_mean = np.mean(diff)\n",
+    "print(diff, diff_mean)\n",
+    "\n",
+    "S = np.sqrt(np.sum((diff - diff_mean)**2) / (N - 1))\n",
+    "\n",
+    "ts = np.sqrt(N) * (diff_mean - mu_0) / S\n",
+    "print(f\"ts: {ts}\")\n",
+    "p_value = 1 - stats.t.cdf(ts, N - 1)\n",
+    "print(f\"p-value: {p_value}\")\n",
+    "\n",
+    "if p_value > alpha:\n",
+    "    print(\"H_0 is rejected thus radial tires have better economy\")\n",
+    "else:\n",
+    "    print(\"H_0 is accepted, thus radial does not make it better\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Part C\n",
+    "\n",
+    "Drivers can drive radically different in the way they brake, shift gears and apply power.\n",
+    "The driver is therefore as important a parameter as the car driven.\n",
+    "\n",
+    "When creating such tests it is very important to only change the variable in question, in this case the tires.\n",
+    "By keeping all other variables the same, one can say that a change in performance likely comes from the other tires."
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3",
+   "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.8.8"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}
diff --git a/sem6/prob/shell.nix b/sem6/prob/shell.nix
new file mode 100644
index 0000000..11c21ce
--- /dev/null
+++ b/sem6/prob/shell.nix
@@ -0,0 +1 @@
+import ../../shells/pymath.nix
diff --git a/sem6/prob/stat6/Opgaver.ipynb b/sem6/prob/stat6/Opgaver.ipynb
new file mode 100644
index 0000000..81a2213
--- /dev/null
+++ b/sem6/prob/stat6/Opgaver.ipynb
@@ -0,0 +1,280 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": 53,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import numpy as np\n",
+    "import matplotlib.pyplot as plt\n",
+    "from IPython.display import display, Math\n",
+    "from scipy import stats"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Problem 1\n",
+    "\n",
+    "> The following data indicate the relationship between x, the specific gravity of\n",
+    "a wood sample, and Y , its maximum crushing strength in compression parallel to\n",
+    "the grain.\n",
+    "> 1. Plot a scatter diagram. Does the linear relationship seem reasonable\n",
+    "> 2. Estimate the regression coefficients\n",
+    "> 3. Predict the maximum crushing strength of a wood sample whose specific gravity is $0.43$\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 25,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "x_i = np.array([0.41,  0.46,  0.44,  0.47,  0.42,  0.39,  0.41,  0.44,  0.43,  0.44])\n",
+    "y_i = np.array([1.850, 2.620, 2.340, 2.690, 2.160, 1.760, 2.500, 2.750, 2.730, 3.120])\n",
+    "x_mean = np.mean(x_i)\n",
+    "y_mean = np.mean(y_i)\n",
+    "\n",
+    "t = np.linspace(np.min(x_i), np.max(x_i))\n",
+    "\n",
+    "n = len(x_i)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 7,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "<matplotlib.collections.PathCollection at 0x7fa6ad4ab2e0>"
+      ]
+     },
+     "execution_count": 7,
+     "metadata": {},
+     "output_type": "execute_result"
+    },
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAATd0lEQVR4nO3df2zc933f8edrMoMwkVMVMRtYsj31R6ZkS5vIYOMM3rJs2SrbzRqlC7CfDuYt8Ap4hY0Zmif/sa4ItqLQZgSF0QhGnDpB3QRdzWhtmkbzMCeOkdkF9SOmbVqB67SJKG+ml6p2EmKRlPf+uKNN0yTvjjreHb98PgBCd5/v544vkdRLX36+3/teqgpJ0ub3l4YdQJLUHxa6JDWEhS5JDWGhS1JDWOiS1BCXDOsTX3bZZbV79+5hfXpJ2pSOHTv2QlVNrLRtaIW+e/dupqenh/XpJWlTSvJnq21zyUWSGsJCl6SGsNAlqSEsdElqCAtdkhpiaGe5SOrekRNzHDp6ijNnF9i5Y5wD+/awf++uYcfSiLHQpRF35MQcB6dmWDh3AYC5swscnJoBsNT1Ki65SCPu0NFTL5f5ooVzFzh09NSQEmlUWejSiDtzdqGncW1dFro04nbuGO9pXFuXhS6NuAP79jA+tu1VY+Nj2ziwb8+QEmlUeVBUGnGLBz49y0WdWOjSJrB/7y4LXB255CJJDWGhS1JDWOiS1BAWuiQ1hIUuSQ1hoUtSQ1joktQQFrokNYSFLkkNYaFLUkN0LPQkr0/yx0m+nuTJJL+6wpwk+Y0kzyR5PMnVGxNXkrSabq7l8v+Av1NV300yBjyS5I+q6tElc64H3tr+uAb4RPtPSdKAdNxDr5bvtu+OtT9q2bQPAp9pz30U2JHk8v5GlSStpas19CTbkpwEngcerKrHlk3ZBXx7yf3T7bHlz3Nzkukk0/Pz8+uMLElaSVeFXlUXqupdwBXAu5O8Y9mUrPSwFZ7nnqqarKrJiYmJnsNKklbX01kuVXUW+DJw3bJNp4Erl9y/AjhzMcEkSb3p5iyXiSQ72rfHgb8LPL1s2u8DH2mf7fIe4C+q6rl+h5Ukra6bs1wuBz6dZBut/wB+t6q+kOSXAKrqMPBF4AbgGeD7wE0blFeStIqOhV5VjwN7Vxg/vOR2Abf0N5okqRe+UlSSGsJCl6SGsNAlqSEsdElqCAtdkhrCQpekhrDQJakhLHRJaggLXZIawkKXpIaw0CWpISx0SWoIC12SGsJCl6SGsNAlqSEsdElqCAtdkhrCQpekhrDQJakhLHRJaoiOhZ7kyiQPJZlN8mSSW1eY8yNJ/iDJ19tzbtqYuJKk1VzSxZzzwO1VdTzJpcCxJA9W1VNL5twCPFVVfz/JBHAqyf1V9YONCC1Jeq2Oe+hV9VxVHW/ffgmYBXYtnwZcmiTAduA7tP4jkCQNSE9r6El2A3uBx5Ztuht4O3AGmAFuraofrvD4m5NMJ5men59fX2JJ0oq6LvQk24EHgNuq6sVlm/cBJ4GdwLuAu5O8aflzVNU9VTVZVZMTExPrDi1Jeq2uCj3JGK0yv7+qplaYchMwVS3PAN8E3ta/mJKkTro5yyXAvcBsVd21yrRvAe9vz38LsAd4tl8hJUmddXOWy7XAjcBMkpPtsTuBqwCq6jDwMeC+JDNAgDuq6oX+x5W2piMn5jh09BRnzi6wc8c4B/btYf/e5ecmaKvrWOhV9Qitkl5rzhng5/oVStIrjpyY4+DUDAvnLgAwd3aBg1MzAJa6XsVXikoj7tDRUy+X+aKFcxc4dPTUkBJpVHWz5CJtGaO4tHHm7EJP49q63EOX2haXNubOLlC8srRx5MTcUHPt3DHe07i2LgtdahvVpY0D+/YwPrbtVWPjY9s4sG/PkBJpVLnkIrWN6tLG4pLPqC0FafRY6FLbzh3jzK1Q3qOwtLF/7y4LXB255CK1ubShzc49dKnNpQ1ttI0+i8pCl5ZwaUMbZRAvEHPJRZIGYBBnUVnokjQAgziLykKXpAEYxAvELHRJGoBBnEXlQVFJGoBBnEVloUvSgGz0WVQuuUhSQ1joktQQLrlIapxRvK79IFjokhplK79ln0sukhplVK9rPwgdCz3JlUkeSjKb5Mkkt64y731JTrbnfKX/USWps1G9rv0gdLPkch64vaqOJ7kUOJbkwap6anFCkh3AbwLXVdW3kvzYxsSVpLWN8nXtN1rHPfSqeq6qjrdvvwTMAssXov4JMFVV32rPe77fQSWpG1v5uvY9raEn2Q3sBR5btumvAD+a5MtJjiX5yCqPvznJdJLp+fn5dQWWpLXs37uLX/vFn2bXjnEC7Noxzq/94k83/oAoQKqqu4nJduArwH+sqqll2+4GJoH3A+PA/wJ+vqq+sdrzTU5O1vT09HpzS9KWlORYVU2utK2r0xaTjAEPAPcvL/O208ALVfU94HtJHgbeCaxa6JKk/urmLJcA9wKzVXXXKtP+G/A3k1yS5A3ANbTW2iVJA9LNHvq1wI3ATJKT7bE7gasAqupwVc0m+RLwOPBD4JNV9cQG5JUkraJjoVfVI0C6mHcIONSPUJKk3vlKUUlqCAtdkhrCQpekhrDQJakhLHRJagivh95wW/VC/9JWZKE32Fa+0L+0Fbnk0mBb+UL/0lZkoTfYVr7Qv7QVWegNttoF/bfChf6lrchCb7CtfKF/aSvyoGiDLR749CwXaWuw0Btu/95dFri0RbjkIkkNYaFLUkNY6JLUEBa6JDWEhS5JDWGhS1JDWOiS1BAdCz3JlUkeSjKb5Mkkt64x92eTXEjy4f7GlCR10s0Li84Dt1fV8SSXAseSPFhVTy2dlGQb8OvA0Q3IKUnqoOMeelU9V1XH27dfAmaBlV56+MvAA8DzfU0oSepKT2voSXYDe4HHlo3vAj4EHO7w+JuTTCeZnp+f7zGqJGktXRd6ku209sBvq6oXl23+OHBHVV14zQOXqKp7qmqyqiYnJiZ6DitJWl1XF+dKMkarzO+vqqkVpkwCn0sCcBlwQ5LzVXWkX0ElSWvrWOhptfS9wGxV3bXSnKr68SXz7wO+YJlL0mB1s4d+LXAjMJPkZHvsTuAqgKpac91ckjQYHQu9qh4B0u0TVtU/v5hAkqT18Q0upE3gyIk533lKHVno0og7cmKOg1MzLJxrnUQ2d3aBg1MzAJa6XsVruUgj7tDRUy+X+aKFcxc4dPTUkBJpVFno0og7c3ahp3FtXRa6NOJ27hjvaVxbl4UujbgD+/YwPrbtVWPjY9s4sG/PkBJpVHlQVBpxiwc+PctFnVjo0iawf+8uC1wdueQiSQ1hoUtSQ1joktQQFrokNYSFLkkNYaFLUkNY6JLUEBa6JDWEhS5JDWGhS1JDWOiS1BAWuiQ1RMdCT3JlkoeSzCZ5MsmtK8z5p0keb398Lck7NyauJGk13Vxt8Txwe1UdT3IpcCzJg1X11JI53wT+VlX9eZLrgXuAazYgryRpFR0LvaqeA55r334pySywC3hqyZyvLXnIo8AVfc4pSeqgpzX0JLuBvcBja0z7l8AfrfL4m5NMJ5men5/v5VNLkjroutCTbAceAG6rqhdXmfO3aRX6HSttr6p7qmqyqiYnJibWk1eStIqu3rEoyRitMr+/qqZWmfMzwCeB66vq//YvoiSpGx0LPUmAe4HZqrprlTlXAVPAjVX1jf5GVBMdOTHne2RKfdbNHvq1wI3ATJKT7bE7gasAquow8O+BNwO/2ep/zlfVZN/TqhGOnJjj4NQMC+cuADB3doGDUzMAlrp0Ebo5y+URIB3mfBT4aL9CqdkOHT31cpkvWjh3gUNHT1no0kXwlaIauDNnF3oal9QdC10Dt3PHeE/jkrpjoWvgDuzbw/jYtleNjY9t48C+PUNKJDVDV6ctSv20uE7uWS5Sf1noGor9e3dZ4FKfueQiSQ1hoUtSQ1joktQQFrokNYSFLkkNYaFLUkNY6JLUEBa6JDWEhS5JDWGhS1JDWOiS1BAWuiQ1hIUuSQ1hoUtSQ1joktQQHQs9yZVJHkoym+TJJLeuMCdJfiPJM0keT3L1xsSVJK2mmze4OA/cXlXHk1wKHEvyYFU9tWTO9cBb2x/XAJ9o/ylJGpCOe+hV9VxVHW/ffgmYBZa/1cwHgc9Uy6PAjiSX9z2tJGlVPa2hJ9kN7AUeW7ZpF/DtJfdP89rSJ8nNSaaTTM/Pz/cYVZK0lq4LPcl24AHgtqp6cfnmFR5SrxmouqeqJqtqcmJiorekkqQ1dVXoScZolfn9VTW1wpTTwJVL7l8BnLn4eJKkbnVzlkuAe4HZqrprlWm/D3ykfbbLe4C/qKrn+phTktRBN2e5XAvcCMwkOdkeuxO4CqCqDgNfBG4AngG+D9zU96SSpDV1LPSqeoSV18iXzingln6FkiT1zleKSlJDWOiS1BAWuiQ1hIUuSQ1hoUtSQ1joktQQFrokNYSFLkkNYaFLUkNY6JLUEBa6JDWEhS5JDWGhS1JDdHP5XG1iR07McejoKc6cXWDnjnEO7NvD/r2veXdASQ1goTfYkRNzHJyaYeHcBQDmzi5wcGoGwFKXGsgllwY7dPTUy2W+aOHcBQ4dPTWkRJI2koXeYGfOLvQ0Lmlzs9AbbOeO8Z7GJW1uFnqDHdi3h/Gxba8aGx/bxoF9e4aUSNJG8qBogy0e+PQsF2lr6FjoST4FfAB4vqrescL2HwF+G7iq/Xz/uap+q99BtT779+6ywKUtopsll/uA69bYfgvwVFW9E3gf8F+SvO7io0mSetGx0KvqYeA7a00BLk0SYHt77vn+xJMkdasfB0XvBt4OnAFmgFur6ocrTUxyc5LpJNPz8/N9+NSSpEX9KPR9wElgJ/Au4O4kb1ppYlXdU1WTVTU5MTHRh08tSVrUj0K/CZiqlmeAbwJv68PzSpJ60I/TFr8FvB/4apK3AHuAZ/vwvK/hhaYkaXXdnLb4WVpnr1yW5DTwK8AYQFUdBj4G3JdkBghwR1W90O+gXmhKktbWsdCr6h932H4G+Lm+JVrFWheastAlaRO99N8LTUnS2jZNoXuhKUla26YpdC80JUlr2zQX5/JCU5K0tk1T6OCFpiRpLZtmyUWStDYLXZIawkKXpIaw0CWpISx0SWqIVNVwPnEyD/zZOh9+GdD368X0wajmgtHNZq7emKs3Tcz1l6tqxeuPD63QL0aS6aqaHHaO5UY1F4xuNnP1xly92Wq5XHKRpIaw0CWpITZrod8z7ACrGNVcMLrZzNUbc/VmS+XalGvokqTX2qx76JKkZSx0SWqIkSv0JNclOZXkmST/boXtH0zyeJKTSaaT/I0l225N8kSSJ5PcNshcS+b9bJILST7c62OHkOtTSZ5P8kQ/M11MriRXJnkoyWz7+3jriOR6fZI/TvL1dq5fHYVcS8a3JTmR5AujkivJnyaZWfy3OkK5diT5vSRPt3/O/vqwcyXZ0/46LX68uK4Oq6qR+QC2AX8C/ATwOuDrwF9dNmc7r6z9/wzwdPv2O4AngDfQuizw/wDeOqhcS+b9T+CLwId7eeygc7XH3wtcDTwx6O/jGl+vy4Gr27cvBb4xCl8vWm+Avr19ewx4DHjPsHMt2fZvgN8BvjAK38f2+J8Cl/XzZ6tPuT4NfLR9+3XAjlHItWz7/6b1AqKeMozaHvq7gWeq6tmq+gHwOeCDSydU1Xer/bcG3ggs3n478GhVfb+qzgNfAT40qFxtvww8ADy/jscOOhdV9TDwnT5l6Uuuqnquqo63b78EzAL9ugj+xeSqqvpu++5Y+6NfZxRc1PcxyRXAzwOf7FOevuTaQOvOleRNtHZk7gWoqh9U1dlh51rm/cCfVFXPr6QftULfBXx7yf3TrPCPOcmHkjwN/CHwL9rDTwDvTfLmJG8AbgCuHFSuJLto/QdyuNfHDinXRupLriS7gb209oaHnqu9rHGS1j/EB6tqJHIBHwf+LfDDPuXpV64C/nuSY0luHpFcPwHMA7/VXqL6ZJI3jkCupf4R8Nn1BBi1Qs8KY6/ZC6qqz1fV24D9wMfaY7PArwMPAl+i9evO+QHm+jhwR1VdWMdj1+ticm2ki86VZDutvZjbqurFUchVVReq6l3AFcC7k7xj2LmSfAB4vqqO9SlLX3K1XVtVVwPXA7ckee8I5LqE1jLjJ6pqL/A9oF/Htfrxc/864BeA/7qeAKP2FnSnefVe9RXAmdUmV9XDSX4yyWVV9UJV3Uv7V6kk/6n9fIPKNQl8Lgm0LrxzQ5LzXT524Lmq6kifMvQ9V5IxWmV+f1VNjUquxQlVdTbJl4HraP1mOLRcwDXALyS5AXg98KYkv11V/2yYuarqSFWdAaiq55N8ntaSxMPDzAU8Cpxe8tvV79G/Qu/Hz9f1wPGq+j/rStCPgwH9+qD1H8yzwI/zykGFv7Zszk/xykHRq4G5Jfd/rP3nVcDTwI8OKtey+ffxysG0nh47qFxLxnbT/4OiF/P1CvAZ4OPD+PlaI9cE7YNnwDjwVeADw861bPx99Peg6MV8vd4IXLrk9teA64adq33/q8Ce9u3/ABwahVztsc8BN603w0jtoVfV+ST/GjhK60jvp6rqySS/1N5+GPgHwEeSnAMWgH9Y7a8E8ECSNwPngFuq6s8HmKunxw47F0CSz9IqgcuSnAZ+pVq/5Qwz17XAjcBMe70a4M6q+uKQc10OfDrJNlpLlb9bVX05RfBiv48b5SJzvQX4fHtP9BLgd6rqSyOQC1oHJe9vL288C9w0Crnax/7+HvCv1pvBl/5LUkOM2kFRSdI6WeiS1BAWuiQ1hIUuSQ1hoUtSQ1joktQQFrokNcT/B2/O3py1bbiSAAAAAElFTkSuQmCC\n",
+      "text/plain": [
+       "<Figure size 432x288 with 1 Axes>"
+      ]
+     },
+     "metadata": {
+      "needs_background": "light"
+     },
+     "output_type": "display_data"
+    }
+   ],
+   "source": [
+    "# Part A\n",
+    "plt.scatter(x_i, y_i)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 28,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "A: -2.8259168241966486, B: 12.245746691871574\n"
+     ]
+    },
+    {
+     "data": {
+      "text/plain": [
+       "[<matplotlib.lines.Line2D at 0x7fa6ad289250>]"
+      ]
+     },
+     "execution_count": 28,
+     "metadata": {},
+     "output_type": "execute_result"
+    },
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAAkWElEQVR4nO3deXzV9Zn28c9NFggkEPYlEGKCgoBAaFwRg7UVd3Gp1bZ2tLWI085Mn+nYjs48M9On01lKp1NnOpNIwSJTrF2k1LF1qNM2IRQQQUAsi3qSsCRhJ4QlIcu5nz9yaCMGyHKS38nJ9X698uLkd74nuchy5eR37vONuTsiItLz9Qk6gIiIRIcKXUQkTqjQRUTihApdRCROqNBFROJEYlDveNiwYZ6VlRXUuxcR6ZE2bdp02N2Ht3ZdYIWelZXFxo0bg3r3IiI9kpntPt91OuUiIhInVOgiInFChS4iEidU6CIicUKFLiISJwKbchGRtlu5uYKFq3ZRWV3LmPQUnpw7kXm5GUHHkhijQheJcSs3V/DUim3UNjQBUFFdy1MrtgGo1OV9dMpFJMYtXLXr92V+Vm1DEwtX7QookcQqFbpIjKusrm3Xcem9VOgiMW5Mekq7jkvvpUIXiXFPzp1ISlLC+46lJCXw5NyJASWSWKUHRUVi3NkHPjXlIhejQhfpAeblZqjA5aJ0ykVEJE6o0EVE4oQKXUQkTqjQRUTihApdRCROqNBFROKECl1EJE6o0EVE4oQKXUQkTqjQRUTixEUL3cz6mdkGM9tqZr8zs6+2ssbM7N/M7D0ze8vMZnZNXBEROZ+27OVyBviwu580syRgjZm96u7rW6y5Fbg08nI1UBD5V0REuslF76F7s5ORV5MiL37OsruBZZG164F0Mxsd3agiInIhbTqHbmYJZrYFOAi85u6vn7MkA9jb4vV9kWPnvp35ZrbRzDYeOnSog5FFRKQ1bSp0d29y9xnAWOAqM5t6zhJr7WatvJ1F7p7n7nnDhw9vd1gRETm/dk25uHs1UATccs5V+4BxLV4fC1R2JpiIiLRPW6ZchptZeuRyCvARYOc5y14GPh2ZdrkGOO7uVdEOKyIi59eWKZfRwPNmlkDzD4AfufsrZrYAwN0LgV8AtwHvAaeBR7sor4iInMdFC93d3wJyWzle2OKyA5+PbjQREWkPPVNURCROqNBFROKECl1EpBu5O03hD0x1R4UKXUSkG4TDzqrf7efegrW88PruLnkfbZlyERGRDqpvDPOzLRUUFocIHTpF5pD+pPdP7pL3pUIXEekCp8408uIbe1lcUkrV8TouHz2Qf3sol9umjiIxoWtOjqjQRUSi6OipepauLWfZunKqTzdwTfYQ/vHeK8i/bDhmre2SEj0qdBGRKNh37DSLS8p48Y091DWE+ejkkTwxJ4eZmYO7LYMKXUSkE3bur+HZ4lJe3lqJAfNyM1iQn82EEWndnkWFLiLSAW+UH6WwKMSvdh6kf3ICj1yXxWevv4Qx6SmBZVKhi4i0UTjs/GbXQQqKQmzcfYzB/ZP4Px+5jE9fO57BA7pmcqU9VOgiIhfR0BTmlbcqKSwqZdeBE2Skp/B3d07mgSvH0T85dmo0dpKIiMSY2vomfvjGHr5bUkZFdS2XjUzlWw9M587pY0jqotHDzlChi4ico/p0PcvW7Wbp2nKOnqrnQ+MH89W7pvDhSSPo06drRw87Q4UuIhJRdbyWJSVlvLBhD6frm/jwpBE8MSeHK7OGBB2tTVToItLrvXfwJM8Wh1i5pYKww13Tx/B4fjaTRg0MOlq7qNBFpNfavOcYhcUhfrn9AH0T+/CJqzJ5bHY244b0Dzpah6jQRaRXcXdWv3uYwqIQ60qPMCgliT+5cQJ/dF0WQ1P7Bh2vU1ToItIrNDaFefXt/RQUhdheVcOogf3469sv58GrMkntGx9VGB//CxGR86hraOInm/axaHUpe46eJnv4AL5x3zTm5WaQnBh7o4edcdFCN7NxwDJgFBAGFrn7M+esGQR8H8iMvM1vuvv3oh9XRKRtauoa+P763Ty3ppzDJ88wfVw6T992OTdPHhnTo4ed0ZZ76I3Al9z9TTNLAzaZ2Wvuvr3Fms8D2939TjMbDuwys+XuXt8VoUVEzudgTR1LflvG8vV7OHmmkRsuG84T+Tlckz2ky7evDdpFC93dq4CqyOUTZrYDyABaFroDadb80UoFjtL8g0BEpFuUHz7Fs6tLeWnTPhrDYW67YjQL8nOYmjEo6Gjdpl3n0M0sC8gFXj/nqu8ALwOVQBrwcXcPt3L7+cB8gMzMzA7EFRF5v7crjlNQHOLVbVUkJvThY3ljmX9DNuOHDgg6Wrdrc6GbWSrwEvBFd6855+q5wBbgw0AO8JqZlZy7zt0XAYsA8vLyuubPXotI3HN31oWOUFAcouTdw6T1TeTx/BwenZXFiLR+QccLTJsK3cySaC7z5e6+opUljwL/5O4OvGdmZcAkYEPUkopIr9cUdl7b3jx6uHXfcYan9eUvb53EJ67OZGC/pKDjBa4tUy4GLAF2uPu3zrNsD3ATUGJmI4GJQGnUUopIr3amsYmVmyt4triU0sOnyBran3+45wrunZlBv6SEoOPFjLbcQ58FPAxsM7MtkWNP0zyiiLsXAl8DlprZNsCAr7j74ejHFemdVm6uYOGqXVRW1zImPYUn505kXm5G0LG63Mkzjbzw+m6WrCnjQM0ZpowZyH98Yia3TB1FQpyOHnZGW6Zc1tBc0hdaUwncHK1QIvIHKzdX8NSKbdQ2NAFQUV3LUyu2AcRtqR8+eYalvy1n2bpyauoauS5nKN/82HSunzAs7kcPO0PPFBWJcQtX7fp9mZ9V29DEwlW74q7Q9x49zaLVpfxo417qm8LcMmUUC/JzmD4uPehoPYIKXaSFWDy1UVld267jPdGOqhoKi0O88lYVfQzuyc3g8fwccoanBh2tR1Ghi0TE6qmNMekpVLRS3kH+dflocHfeKD9GQdF7/GbXIQYkJ/CZWVl89vpsRg3qvaOHnaFCF4mI1VMbT86d+L4fNAApSQk8OXdiYJk6Ixx2frXzIIXFITbtPsbQAcn8xc2X8fA1WQzqr9HDzlChi0TE6qmNsz9MYu1UUHs1NIV5eUslhcUh3j14krGDU/jqXVN4IG8cKckaPYwGFbpIRCyf2piXm9HjCvys0/WNvLhhL4tLSqk8XsekUWk88+AMbr9iNIkJ8bV9bdBU6CIR8XZqI2jHTtXz/Lpylq4tp/p0A1ddMoSv33MFcyYO1+hhF1Ghi0TEy6mNoFVU17K4pJQXN+yltqGJj1w+gifm5PCh8UOCjha4rp6iUqGLtNCTT20E7d0DJygsLuVnWyoAuGvGGBbk53DZyLSAk8WG7piiUqGLSKds2n2MgqIQ/7vjAClJCTx87Xgem51NRgw89hBLumOKSoUuIu3m7hTtOkRBcYgNZUdJ75/En950KY9cl8WQAclBx4tJ3TFFpUIXkTZrbArz821VFBSF2Ln/BGMG9eNv7pjMx68cx4C+qpML6Y4pKn0GROSi6hqa+PHGvSwqKWXv0VouHZHKNz82nbumjyE5UaOHbdEdU1QqdBE5r+OnG/iv9eV877flHDlVz8zMdP7mjincNGkEfbR9bbt0xxSVCl1EPuBATR1L1pSxfP1uTtU3cePE4TwxZwJXZg3WDHkndPUUlQpdRH4vdOgki4pL+enmCprcuXPaaB7Pz+Hy0QODjiZtoEIXEbburaagKMSq7ftJTujDg1eN43Ozsxk3pH/Q0aQdVOgivZS7U/LuYQqLQ6wNHWFgv0S+cOME/ui6LIal9g06XqfE4r723UGFLtLLNIWdV99uHj38XWUNIwf25a9uu5yHrs4kNQ5GD2N1X/vu0PM/eyLSJnUNTbz05j4WrS5l95HTZA8bwD/fdwXzcjPomxg/29fG6r723eGihW5m44BlwCggDCxy92daWTcH+DaQBBx29/xoBhWRjqmpa2D5+j0sWVPG4ZNnmD52EE99aiYfnTyKhDgcPYzVfe27Q1vuoTcCX3L3N80sDdhkZq+5+/azC8wsHfhP4BZ332NmI7omroi01cETdTy3ppzl63dz4kwjsy8dxhP5M7g2Z2hcjx7G8r72Xe2ihe7uVUBV5PIJM9sBZADbWyz7BLDC3fdE1h3sgqwi0gblh0+xqKSUn2zaR2NTmFuvGM0T+TlMzRgUdLRu0Zv3tW/XOXQzywJygdfPueoyIMnMioA04Bl3X9bK7ecD8wEyMzM7EFdEzuftiuMUFof4xbYqEvv04b4PjWX+DdlcMmxA0NG6VW/e197cvW0LzVKBYuDr7r7inOu+A+QBNwEpwDrgdnd/53xvLy8vzzdu3NjR3CJC8+jhutIjFBSFKHn3MKl9E/nUNeP5zKwsRgzsF3Q86QJmtsnd81q7rk330M0sCXgJWH5umUfso/mB0FPAKTNbDUwHzlvoItJx4bDzy+37KSguZeveaoal9uXLt0zkk1ePZ1BKUtDxJCBtmXIxYAmww92/dZ5lPwO+Y2aJQDJwNfCvUUspIgDUN4ZZubmCwtUhSg+dInNIf75+z1TumzmWfknxM3ooHdOWe+izgIeBbWa2JXLsaSATwN0L3X2Hmf0P8BbNo42L3f3tLsgr0iudPNPIixv2sLikjP01dUwePZB/fyiXW6eOIjFB29dKs7ZMuawBLjrj5O4LgYXRCCUizY6cPMPSteU8v7acmrpGrs0eyjfun8bsS4fF9eihdIyeKSoSg/YePc3iklJ+uHEvZxrD3Dx5JAvyc8jNHBx0NIlhKnSRGLJzfw3PFpfy8tZK+hjcPSODBfnZTBiRFnQ06QFU6CIx4I3yoxQUhfj1zoP0T07g0euy+OzsSxg9KP6f3SjRo0IXCUg47Px650EKi0Ns3H2MIQOS+fOPXsanrx1Pev/koONJD6RCF+lmDU1h/ntrJYXFId45cJKM9BS+etcUHsgbR0qyRg+l41Toca63bvQfi07XN/LDN/ayuKSMiupaJo5M418/Pp07po0hSaOHEgUq9DjWmzf6jyXHTtWzbN1ulq4t49jpBq7MGszX5k3hxokjNHooUaVCj2O9eaP/WFBZXcuSNWX8YMMeTtc3cdOkETwxJ4e8rCFBR5M4pUKPY715o/8gvXfwBIXFpazcXIEDd08fw+P5OUwcpdFD6Voq9DjWmzf6D8Kbe45RWBTil9sP0C+pD5+6ZjyPzb6EsYP7Bx1NegkVehzrzRv9dxd3p/idQxQUhXi97CiDUpL40w9P4JFZlzBkgEYPpXup0ONYb97ov6s1NoX5+bYqCotL2VFVw+hB/fi/d0zmwSvHMaCvvq0kGPrKi3PzcjNU4FFU19DEjzftY9HqEHuP1pIzfAAL75/G3TMySE7U6KEES4Uu0gbHaxv4/vrdfO+3ZRw+Wc+Mcen81W2TuXnySPr00eihxAYVusgFHKip47k1ZSx/fQ8nzzSSf9lwnpiTw9WXDNEMucQcFbpIK0oPnWTR6lJWvFlBYzjM7dPGsCA/myljBgUdTeS8VOgiLby1r5rC4hCvvr2fpIQ+PHDlWD43O5vxQwcEHU3kolTo0uu5O2tDRygoCrHmvcOk9Uvkj+fk8Mh1lzA8rW/Q8UTaTIUuvVZT2Fn1u/0UFIXYVnGcEWl9efq2STx0VSZp/ZKCjifSbip06XXONDax4s0KFq0upezwKS4ZNoB/vPcK7p2ZQd9EbV8rPddFC93MxgHLgFFAGFjk7s+cZ+2VwHrg4+7+k2gGFemsE3UNvPD6HpasKePgiTNckTGIgk/O5OYpo0jQ6KHEgbbcQ28EvuTub5pZGrDJzF5z9+0tF5lZAvDPwKouyCnSYYdOnOF7vy3jv9bv5kRdI9dPGMa/fnwG1+UM1eihxJWLFrq7VwFVkcsnzGwHkAFsP2fpnwAvAVdGO6RIR+w5cppFJSF+tHEfDU1hbp06igX5OUwbmx50NJEu0a5z6GaWBeQCr59zPAO4B/gwFyh0M5sPzAfIzMxsZ1SRttleWUNhcYhX3qoksU8f7vtQBp+bnU328NSgo4l0qTYXupml0nwP/IvuXnPO1d8GvuLuTRf6FdbdFwGLAPLy8rzdaUXOw915vewoBUUhit85RGrfRD43O5vPXH8JIwf2CzqeSLdoU6GbWRLNZb7c3Ve0siQPeDFS5sOA28ys0d1XRiuoSGvCYee1HQcoKAqxZW81w1KTeXLuRD51zXgGpWj0UHqXtky5GLAE2OHu32ptjbtf0mL9UuAVlbl0pfrGMCu3VPBscYjQoVNkDunP38+byv0fGku/JI0eSu/Ulnvos4CHgW1mtiVy7GkgE8DdC7smmsgHnTrTyA82NI8eVh2v4/LRA/m3h3K5beooEhO0fa30bm2ZclkDtHm2y90f6UwgkdYcOXmG59eW8/y63RyvbeCa7CH8033TuOHSYRo9FInQM0Ulpu07dprFJWW8+MYe6hrC3Dx5JAvm5DAzc3DQ0brVys0V+stTclEqdIlJu/afoLA4xMtbKzHgntwMHs/PZsKItKCjdbuVmyve97dhK6preWrFNgCVuryPCl1iyhvlRyksCvGrnQfpn5zAI9dl8dnrL2FMekrQ0QKzcNWu9/2hb4DahiYWrtqlQpf3UaFL4MJh5ze7DlJQFGLj7mMM7p/En3/0Mj597XjS+ycHHS9wldW17TouvZcKXQLT0BTmlbcqKSwqZdeBE2Skp/B3d07mgSvH0T9ZX5pnjUlPoaKV8u7Nv7VI6/RdI92utr6JH76xh++WlFFRXctlI1P51gPTuXP6GJI0evgBT86d+L5z6AApSQk8OXdigKkkFqnQpdtUn65n2brdLF1bztFT9eSNH8z/u3sKN04cQR9tX3teZ8+Ta8pFLkaFLl2u6ngtS0rKeGHDHk7XN/HhSSN4Yk4OV2YNCTpajzEvN0MFLhelQpcu897BkzxbHGLllgrCDndNH8Pj+dlMGjUw6GgicUmFLlG3ec8xCotD/HL7Afom9uETV2Xy2Oxsxg3pH3Q0kbimQpeocHdWv3uYgqL3WF96lEEpSXzhxgk8cl0WQ1P7Bh1PpFdQoUunNDaF+cXb+yksCrG9qoZRA/vx17dfzoNXZZLaV19eIt1J33HSIXUNTfxk0z4WrS5lz9HTZA8fwDfun8a8GRkkJ2r0UCQIKnRpl5q6Br6/fjfPrSnn8MkzTB87iKdvm8nNk0dp9FAkYCp0aZODNXUs+W0Zy9fv4eSZRm64bDgL8rO5Nnuotq8ViREqdLmg8sOneHZ1KS9t2kdjOMzt08bw+A3ZTM0YFHQ0ETmHCl1atW3fcQqLQ7z6dhWJCX34WN5Y5t+QzfihA4KOJiLnoUKX33N31oWOUFAcouTdw6T1TWRBfg6PzMpiRFq/oOOJyEWo0IWmsPPa9v0UFIXYuu84w9P68pe3TuITV2cysF9S0PFEpI1U6L3YmcYmVm6u4NniUkoPnyJraH/+8d4ruCc3g35JCUHHE5F2umihm9k4YBkwCggDi9z9mXPWfBL4SuTVk8AT7r41ylklSk6eaeQHr+9h8ZpSDtScYWrGQP7jEzO5ZeooEjR6KNJjteUeeiPwJXd/08zSgE1m9pq7b2+xpgzId/djZnYrsAi4ugvySiccPnmGpb8tZ9m6cmrqGpk1YSj/8rEZzJqg0UOReHDRQnf3KqAqcvmEme0AMoDtLdasbXGT9cDYKOeUTth79DTfLSnlh2/spb4pzC1TRrEgP4fp49KDjiYiUdSuc+hmlgXkAq9fYNlngVfPc/v5wHyAzMzM9rxr6YAdVTUUFod45a0q+hjcmzuW+fnZ5AxPDTqaiHSBNhe6maUCLwFfdPea86y5keZCv7616919Ec2nY8jLy/N2p5WLcnc2lB2loDhE0a5DDEhO4DOzsvjs9dmMGqTRQ5F41qZCN7Mkmst8ubuvOM+aacBi4FZ3PxK9iNIW4bDzq50HKSh6jzf3VDN0QDJ/cfNlPHxNFoP6a/RQpDdoy5SLAUuAHe7+rfOsyQRWAA+7+zvRjSgX0tAU5mdbKnm2OMS7B08ydnAKX7t7Ch/LGxfTo4crN1fob2SKRFlb7qHPAh4GtpnZlsixp4FMAHcvBP4GGAr8Z2RaotHd86KeVn7vdH0jL27Yy+KSUiqP1zFpVBrPPDiD268YTWJCbG9fu3Jzxfv+in1FdS1PrdgGoFIX6YS2TLmsAS440+bujwGPRSuUnN/RU/U8v7ac59eVU326gasuGcLX772COZcN7zGjhwtX7fp9mZ9V29DEwlW7VOginaBnivYQFdW1LC4p5cUNe6ltaOIjl4/kiTnZfGj8kKCjtVtldW27jotI26jQY9w7B05QWBzi5S2VANw1YwxP5Odw6ci0gJN13Jj0FCpaKe8x6SkBpBGJHyr0GLVp9zEKikL8744DpCQl8PC143lsdjYZcVB6T86d+L5z6AApSQk8OXdigKlEej4Vegxxd4p2HaKgOMSGsqOk90/iix+5lD+6NovBA5KDjhc1Z8+Ta8pFJLpU6DGgsSnMz7dVUVAUYuf+E4we1I+/uWMyD141jv7J8fkpmpeboQIXibL4bIseora+iR9v2sui1aXsO1bLpSNS+ebHpnPX9DEkJ8b26KGIxB4VegCOn25g2bpylq4t58ipemZmpvO3d07hpkkj6KPta0Wkg1To3Wj/8TqWrCnlhdf3cKq+iTkTh/PHcyZwZdbgHjNDLiKxS4XeDUKHTrKouJQVm/cRdrhz2mgez8/h8tEDg44mInFEhd6Ftu6tpqAoxKrt+0lO6MNDV2XyudnZjBvSP+hoIhKHVOhR5u6UvHuYwuIQa0NHGNgvkc/PmcAjs7IYlto36HgiEsdU6FHSFHZefbt59PB3lTWMHNiXv7rtch66OpPUvvowi0jXU9N0Ul1DEy+9uY/vri6l/MhpsocN4Bv3TePu3DH0TYzd7WtFJP6o0Duopq6B5ev3sGRNGYdPnmH62EEUfHImN08ZRYJGD0UkACr0djp4oo7n1pSzfP1uTpxpZPalw3gifwbX5gzV6KGIBEqF3kblh0+xqKSUn2zaR2NTmFuvGM0T+TlMzRgUdDQREUCFflFvVxynoDjEq9uqSOzTh/vzxjJ/djZZwwYEHU1E5H1U6K1wd9aVHqGgKETJu4dJ65vI/Bty+MysLEYM7Bd0PBGRVqnQWwiHnV9u309BcSlb91YzLLUvX75lIp+6ZjwD+yUFHU9E5IJU6EB9Y5iVmysoXB2i9NApMof05+/nTeX+D42lX5JGD0WkZ7hooZvZOGAZMAoIA4vc/Zlz1hjwDHAbcBp4xN3fjH7c6Dp5ppEXN+xhcUkZ+2vqmDx6IP/+UC63Th1FYoK2rxWRnqUt99AbgS+5+5tmlgZsMrPX3H17izW3ApdGXq4GCiL/xqQjJ8+wdG05z68tp6aukWuzh/KN+6cx+9JhGj0UkR7rooXu7lVAVeTyCTPbAWQALQv9bmCZuzuw3szSzWx05LYxY+/R0ywuKeWHG/dypjHMzZNHsiA/h9zMwUFHExHptHadQzezLCAXeP2cqzKAvS1e3xc59r5CN7P5wHyAzMzMdkbtuJ37a3i2uJSXt1bSx+Ce3Azm35DDhBGp3ZZBRKSrtbnQzSwVeAn4orvXnHt1KzfxDxxwXwQsAsjLy/vA9dH2RvlRCopC/HrnQfonJ/DodVl8dvYljB6U0tXvWkSk27Wp0M0sieYyX+7uK1pZsg8Y1+L1sUBl5+O1Xzjs/HrnQQqLQ2zcfYwhA5L50kcv4+Frx5PePzmISCIi3aItUy4GLAF2uPu3zrPsZeALZvYizQ+GHu/u8+cNTWH+e2slhcUh3jlwkoz0FL561xQeyBtHSrJGD0Uk/rXlHvos4GFgm5ltiRx7GsgEcPdC4Bc0jyy+R/PY4qNRT3oep+sb+eEbe1lcUkZFdS0TR6bx7Y/P4PZpo0nS6KGI9CJtmXJZQ+vnyFuuceDz0QrVFsdO1bNs3W6Wri3j2OkGrswazNfmTeHGiSM0eigivVKPe6Zo1fFavru6jB9s2ENtQxMfuXwEC/JzyMsaEnQ0EZFA9bhC37KnmufXlXP39DE8np/DxFFpQUcSEYkJPa7Qb54yitVfvpGMdI0eioi01OMeNUzoYypzEZFW9LhCFxGR1qnQRUTihApdRCROqNBFROKECl1EJE6o0EVE4oQKXUQkTqjQRUTihApdRCRO9Lin/kv7rNxcwcJVu6isrmVMegpPzp3IvNyMoGOJSBdQocexlZsreGrFNmobmgCoqK7lqRXbAFTqInFIp1zi2MJVu35f5mfVNjSxcNWugBKJSFdSocexyuradh0XkZ5NhR7HxpxnV8rzHReRnk2FHseenDuRlKT3/4HslKQEnpw7MaBEItKV9KBoHDv7wKemXER6h4sWupk9B9wBHHT3qa1cPwj4PpAZeXvfdPfvRTuodMy83AwVuEgv0ZZTLkuBWy5w/eeB7e4+HZgD/IuZJXc+moiItMdFC93dVwNHL7QESDMzA1IjaxujE09ERNoqGg+Kfge4HKgEtgF/5u7h1haa2Xwz22hmGw8dOhSFdy0iImdFo9DnAluAMcAM4DtmNrC1he6+yN3z3D1v+PDhUXjXIiJyVjQK/VFghTd7DygDJkXh7YqISDtEY2xxD3ATUGJmI4GJQGkU3u4HaKMpEZHza8vY4g9onl4ZZmb7gL8FkgDcvRD4GrDUzLYBBnzF3Q9HO6g2mhIRubCLFrq7P3SR6yuBm6OW6DwutNGUCl1EpAc99V8bTYmIXFiPKXRtNCUicmE9ptC10ZSIyIX1mM25tNGUiMiF9ZhCB200JSJyIT3mlIuIiFyYCl1EJE6o0EVE4oQKXUQkTqjQRUTihLl7MO/Y7BCwu4M3HwZEfb+YKIjVXBC72ZSrfZSrfeIx13h3b3X/8cAKvTPMbKO75wWd41yxmgtiN5tytY9ytU9vy6VTLiIicUKFLiISJ3pqoS8KOsB5xGouiN1sytU+ytU+vSpXjzyHLiIiH9RT76GLiMg5VOgiInEi5grdzG4xs11m9p6Z/WUr199tZm+Z2RYz22hm17e47s/M7G0z+52ZfbE7c7VYd6WZNZnZ/e29bQC5njOzg2b2djQzdSaXmY0zs9+Y2Y7I5/HPYiRXPzPbYGZbI7m+Ggu5WhxPMLPNZvZKrOQys3Iz23b2ezWGcqWb2U/MbGfk6+zaoHOZ2cTIx+nsS02HOszdY+YFSABCQDaQDGwFJp+zJpU/nPufBuyMXJ4KvA30p3lb4P8FLu2uXC3W/Rr4BXB/e27b3bkix28AZgJvd/fn8QIfr9HAzMjlNOCdWPh40fwH0FMjl5OA14Frgs7V4ro/B14AXomFz2PkeDkwLJpfW1HK9TzwWORyMpAeC7nOuX4/zU8galeGWLuHfhXwnruXuns98CJwd8sF7n7SI/9rYABw9vLlwHp3P+3ujUAxcE935Yr4E+Al4GAHbtvduXD31cDRKGWJSi53r3L3NyOXTwA7gGhtgt+ZXO7uJyOvJkVeojVR0KnPo5mNBW4HFkcpT1RydaEO5zKzgTTfkVkC4O717l4ddK5z3ASE3L3dz6SPtULPAPa2eH0frXwzm9k9ZrYT+Dnwmcjht4EbzGyomfUHbgPGdVcuM8ug+QdIYXtvG1CurhSVXGaWBeTSfG848FyR0xpbaP5GfM3dYyIX8G3gy0A4SnmilcuBX5rZJjObHyO5soFDwPcip6gWm9mAGMjV0oPADzoSINYK3Vo59oF7Qe7+U3efBMwDvhY5tgP4Z+A14H9o/nWnsRtzfRv4irs3deC2HdWZXF2p07nMLJXmezFfdPeaWMjl7k3uPgMYC1xlZlODzmVmdwAH3X1TlLJEJVfELHefCdwKfN7MboiBXIk0n2YscPdc4BQQrce1ovF1nwzcBfy4IwFi7U/Q7eP996rHApXnW+zuq80sx8yGufthd19C5FcpM/uHyNvrrlx5wItmBs0b79xmZo1tvG2353L3lVHKEPVcZpZEc5kvd/cVsZLr7AJ3rzazIuAWmn8zDCwXcDVwl5ndBvQDBprZ9939U0HmcveV7l4J4O4HzeynNJ+SWB1kLmA9sK/Fb1c/IXqFHo2vr1uBN939QIcSROPBgGi90PwDphS4hD88qDDlnDUT+MODojOBihavj4j8mwnsBAZ3V65z1i/lDw+mteu23ZWrxbEsov+gaGc+XgYsA74dxNfXBXINJ/LgGZAClAB3BJ3rnONziO6Dop35eA0A0lpcXgvcEnSuyOslwMTI5b8DFsZCrsixF4FHO5ohpu6hu3ujmX0BWEXzI73PufvvzGxB5PpC4D7g02bWANQCH/fIRwJ4ycyGAg3A5939WDfmatdtg84FYGY/oLkEhpnZPuBvvfm3nCBzzQIeBrZFzlcDPO3uvwg412jgeTNLoPlU5Y/cPSojgp39PHaVTuYaCfw0ck80EXjB3f8nBnJB84OSyyOnN0qBR2MhV+Sxv48Cj3c0g576LyISJ2LtQVEREekgFbqISJxQoYuIxAkVuohInFChi4jECRW6iEicUKGLiMSJ/w9Xw43AN/vRfQAAAABJRU5ErkJggg==\n",
+      "text/plain": [
+       "<Figure size 432x288 with 1 Axes>"
+      ]
+     },
+     "metadata": {
+      "needs_background": "light"
+     },
+     "output_type": "display_data"
+    }
+   ],
+   "source": [
+    "# Part B\n",
+    "\n",
+    "S_xY = sum((x * y for (x, y) in zip(x_i, y_i))) - n * x_mean * y_mean\n",
+    "S_xx = sum((x * x for x in x_i)) - n * x_mean * x_mean\n",
+    "\n",
+    "B = S_xY / S_xx\n",
+    "A = y_mean - B * x_mean\n",
+    "print(f\"A: {A}, B: {B}\")\n",
+    "\n",
+    "plt.scatter(x_i, y_i)\n",
+    "plt.plot(t, B * t + A)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 33,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "No 'a' percentage is given, so we just fint the regression value at 0.43\n",
+      "The predicted maximum strength is 2.4397542533081285\n",
+      "Wait this makes no sense, 0.43 is in the sample set.b\n"
+     ]
+    }
+   ],
+   "source": [
+    "# Part C\n",
+    "where = 0.43\n",
+    "print(f\"No 'a' percentage is given, so we just fint the regression value at {where}\")\n",
+    "\n",
+    "print(f\"The predicted maximum strength is {B * where + A}\")\n",
+    "\n",
+    "# TODO TODO GET HELP\n",
+    "print(\"Wait this makes no sense, 0.43 is in the sample set.b\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Problem 2\n",
+    "\n",
+    "This seemes like a continuation.\n",
+    "\n",
+    "> In Problem 4.\n",
+    "> 1. Estimate the variance of an individual response.\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 41,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "The estimated variance is 0.10566006616257498\n"
+     ]
+    }
+   ],
+   "source": [
+    "S_YY = sum((y*y for y in y_i)) - n * y_mean * y_mean\n",
+    "SS_R = (S_xx * S_YY - (S_xY)**2) / S_xx\n",
+    "sigma2 = SS_R / (n - 2)\n",
+    "print(f\"The estimated variance is {sigma2}\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Problem 3\n",
+    "\n",
+    "> The following data set presents the heights of 12 male law school classmates whose\n",
+    "law school examination scores were roughly equal. It also gives their annual salaries\n",
+    "5 years after graduation. Each of them went into corporate law. The height is in\n",
+    "inches and the salary in units of $1,000.\n",
+    "> 1. Do the above data establish the hypothesis that a lawyer’s salary is related to\n",
+    "his height? Use the 5 percent level of significance.\n",
+    "> 2. What was the null hypothesis in part (a)?\n",
+    "\n",
+    "\n",
+    "## Part A\n",
+    "\n",
+    "We let the null hypothesis be that there is no relation between the salary and height, thus that $\\beta = 0$."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 45,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "x_i = np.array([64, 65, 66, 67, 69, 70, 72, 72, 74, 74, 75, 76])\n",
+    "y_i = np.array([91, 94, 88, 103, 77, 96, 105, 88, 122, 102, 90, 114])\n",
+    "\n",
+    "alpha = 0.05\n",
+    "\n",
+    "n = len(x_i)\n",
+    "\n",
+    "x_mean = np.mean(x_i)\n",
+    "y_mean = np.mean(y_i)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 56,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "B: 1.4571428571428193\n",
+      "TS: 1.7483281067751004\n",
+      "p_value: 0.11098065306186689\n",
+      "h_0 is accepted and thus there is no connection, which kind of makes sense.\n"
+     ]
+    }
+   ],
+   "source": [
+    "# Part A\n",
+    "S_xx = sum((x * x for x in x_i)) - n * x_mean**2\n",
+    "S_YY = sum((y * y for y in y_i)) - n * y_mean**2\n",
+    "S_xY = sum((x * y for (x, y) in zip(x_i, y_i))) - n * x_mean * y_mean\n",
+    "SS_R = (S_xx * S_YY - (S_xY)**2) / S_xx\n",
+    "\n",
+    "print(f\"B: {(B := S_xY / S_xx)}\")\n",
+    "TS = np.sqrt(((n - 2) * S_xx) / SS_R) * np.abs(B)\n",
+    "print(f\"TS: {TS}\")\n",
+    "\n",
+    "p_value = 2 * (1 - stats.t.cdf(TS, n-2))\n",
+    "print(f\"p_value: {p_value}\")\n",
+    "if p_value < alpha:\n",
+    "    print(\"h_0 is rejected and there is a connection between salary and height. Why tho\")\n",
+    "else:\n",
+    "    print(\"h_0 is accepted and thus there is no connection, which kind of makes sense.\") "
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3",
+   "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.9.2"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}
diff --git a/shells/pymath.nix b/shells/pymath.nix
new file mode 100644
index 0000000..5fa8b1e
--- /dev/null
+++ b/shells/pymath.nix
@@ -0,0 +1,16 @@
+{ pkgs ? import <nixpkgs> {}, pythonPackages ? pkgs.python38Packages }:
+
+pkgs.mkShell {
+    buildInputs = with pythonPackages; [
+        jupyterlab
+        jupyter-c-kernel
+        numpy
+        scipy
+        matplotlib
+    ];
+
+    shellHook = ''
+        echo Lets do some python
+        alias start="jupyter lab"
+        '';
+}