Skip to content

Instantly share code, notes, and snippets.

Show Gist options
  • Save davipatti/69ae07fd202e6422d7cf2cd35ebd7f52 to your computer and use it in GitHub Desktop.
Save davipatti/69ae07fd202e6422d7cf2cd35ebd7f52 to your computer and use it in GitHub Desktop.
[Bayesian prior distributions don't need to look like what you think the posterior distribution will look like]
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Bayesian prior distributions don't need to look like what you think the posterior distribution will look like"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"import pymc as pm\n",
"import matplotlib.pyplot as plt\n",
"\n",
"np.random.seed(42)\n",
"plt.style.use(\"seaborn-v0_8-whitegrid\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Simulate some data from a Normal(mu=5, sigma=1) distribution:"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAiMAAAGsCAYAAAAPJKchAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAytUlEQVR4nO3df3zN9f//8fs6sx+sTRhhFROiZpZlkr1DCm/e7zIpKT+KJr/2rSQbhcJHkn75FSGJoiHlR4SPVJ83Uxjrh7GNENZWGbPZ4WzfP1x23k7bmY0zz/24XS+XXTjP83w9X4/zep5zdt/r9Tqv45aXl5cnAAAAQ64zXQAAAKjcCCMAAMAowggAADCKMAIAAIwijAAAAKMIIwAAwCjCCAAAMIowAgAAjCKMAJfgGoAAcO0RRlDhrFq1Sk2bNtWxY8eKvYzVatX//M//aM2aNaVYWdHi4uLUtGlTDR06tND7r+RxmeCszn379qlz586yWq2GKisoNjZWTZs2LfDz6quvunQ90dHR6tixo0vHvNYyMzPVsWNHrVq1yqF9xYoVioyMNFQVKgrCCCDp999/14cffqgLFy6YLkVbtmzRF198YboMl8rJydHo0aM1atQoeXh4mC7H7pdfflHDhg21fPlyh5+nnnrKpesZOnSoZs6c6dIxr6WMjAwNHjxYv/32W4H7evbsqbS0NK1YscJAZago3E0XAMCRr6+vJk+erLZt26pWrVqmy3GJjz/+WO7u7urUqZPpUhz88ssvCgoKUsuWLUt1PTfffHOpjl+atmzZosmTJ+vs2bOF3u/m5qbBgwfr1VdfVffu3eXl5XWNK0RFwJ4RlGu5ubmaPXu22rdvr+DgYA0dOlQZGRkF+m3evFl9+vRRSEiI7rjjDnXp0kVLly6VJB07dkz33XefJCkmJsZhd3psbKwiIiLUsmVLtWjRQg8++KC+/PLLUn1Mzz33nLKysjRhwoTL9j18+LCioqJ0zz33qGXLlurbt6927dplv//YsWNq2rSpPvjgA3Xp0kXBwcFauXKlZsyYoS5dumjTpk3q3r27goKC9OCDD2rPnj2Kj49Xr1691KJFC3Xv3l3bt293WGdR27IwVqtVH3zwgbp3716grg0bNmjo0KFq2bKl2rZtq9mzZyszM1NjxoxRq1at1LZtW02bNs1+Lk/+oay4uDiHdfTt21d9+/Z1GNvZT36/vLw8JSYmqlmzZpfdzpfz448/qn///mrVqpVCQkI0YMAAxcfH2+//+2Ga8+fP64033tA//vEPtWjRQgMHDtTq1asdDm9FR0dr4MCBWr58uTp16qQWLVqod+/eOnTokLZu3ap//etfCg4OVq9evfTLL7841HO55210dHSR2yh/+54+fVrDhw/XXXfdpfnz5zt9/B06dFBOTo5Wrlx51dsSlRN7RlCuTZs2TYsXL9aQIUMUHBysL7/8UtOnT3fo8/XXX2vYsGHq16+fRowYoXPnzunjjz/Wq6++qjvuuEPNmjXTzJkzNXz4cA0ZMkQPPPCAJGnp0qWaNGmSRowYoVatWikjI0Pvv/++XnjhBYWEhOjGG28slcfUqFEjjRgxQtOnT9fatWsdfolfKikpSY888ogaNGigl156SVWqVNHixYvVv39/LVy4UK1bt7b3nTFjhsaOHSsfHx8FBwcrNjZWJ0+e1GuvvabnnntOVatW1cSJExUVFaUqVaromWeeUd26de33f/311/Ly8rrstgwODi5QZ1xcnFJTU+3b9VIvvfSSnnjiCfXt21crVqzQO++8oy+++EJt27bVzJkztXHjRs2fP1933HGHunbtWqztV7t2bS1fvtzp/T4+PpKkI0eO6OzZs0pISFDnzp117NgxBQQEaMiQIXrooYeKtS7p4rkUgwYNUps2bTRjxgxZrVbNmTNHAwcO1Ndff63rr7++wDLjxo3T2rVrNWLECDVr1kxr167Vyy+/XKDfnj179Pvvvys6Olo5OTmaMGGCIiMj5ebmpqioKHl7e2v8+PF64YUXtG7dOknFe94OHTpUvXv3dvqYbr31VkmSl5eX1q1bp8DAwCLPVfL09FSHDh20Zs0aPf7448XedkA+wgjKrdOnT+ujjz7Sk08+qeHDh0uSwsPD9fvvv+vbb7+190tKSlKPHj00duxYe1tISIjCwsIUFxen4OBg+1/HN998s5o3by5JOnr0qAYOHOhwQmn9+vUVERGhXbt2qVu3bgVqOn/+vN566y198cUXyszMVPPmzdW1a1d17txZkjRv3jw9+eSTql+/fpGPbeDAgdq0aZMmTpyoNm3aFHq4ZubMmfLw8NDixYvtv2Dbt2+v7t276/XXX3c4ht+1a1f17NnTYfns7GyNHz9e//jHP+zbafr06Zo8ebIefvhhSVJWVpaioqJ06NAhNWvWrFjb8u927NghX19fNWzYsMB94eHhevbZZyVJjRs31tq1a1WzZk2NGzdOktSmTRutWbNGu3fvLnYY8fDwKNZhl/y9CceOHVN0dLTc3d21evVqjR49WlarVY888kix1peUlKS//vpL/fr105133ilJCgwM1PLly3X27NkCYeTIkSP67LPPNHr0aD355JOSLm6H9PR0fffddw59z549q7fffluNGjWSJO3cuVPLli3TokWLdPfdd0uSfv31V02dOlWnT5+Wr69vsZ63N998c7EOHXl4eCgwMLBY2yEoKEjr169XZmam/fkIFBdhBOVWfHy8zp8/rw4dOji0d+3a1SGMDBo0SNLFN/ZDhw7pyJEjSkhIkKQiP9kRHR0t6WLoSUlJ0a+//mrffe1sudTUVP3888+aMGGCqlSpoq1bt2rGjBmaNGmSpIu/XKtXr37Zx2axWDRlyhT16NFDr7zyimbMmFGgz86dO9WhQweHN353d3d169ZNs2bNcjjG7+xQRP4vT0n2wHNpoMiv9fTp05KubFsePXrUafgKCQkpsP4WLVrY29zc3OTn56czZ84UurwzRZ2I7ObmJovForvuukvvvfeewsLCVLVqVUkXQ8Gff/6pd999V7169ZKbm9tl19W4cWPVqFFDzzzzjLp06aLw8HDdc889GjVqVKH94+LilJeXpy5duji0d+/evUAY8fPzswcR6fJz5OvrW6znbW5urnJzc50+JovFUqzHfqn69evLZrPp5MmT9j0rQHERRlBu5Z8bcsMNNzi0+/v7O9z+888/NX78eG3evFlubm665ZZbFBoaKqno64ocOXJE48aN0/bt21WlShUFBgbqtttuK3K5OnXqaOHChbruuounY917770aO3asjhw5Ih8fH9WpU6fYj+/WW2/V8OHD9eabb9p3wV8qIyOj0D0mtWrVUl5enjIzM+1t+b9s/66wv2C9vb2d1nQl2zIzM9PpmIWt31mtxXXpOUCFad26tT766CPVrFmzQJCVLs7Zf/7zH6Wnpxd4LhWmWrVqWrp0qebMmaMvv/xSy5cvl5eXlx588EG99NJLBT499Oeff0qSatas6dD+99tS4dtHKnobFed5O2bMGH322WdOx1i8eLHCwsKc3l9UTSUNjoBEGEE5lh9C/vjjD4ddyadOnXLo98ILLyglJUWLFi1SSEiIPDw8lJ2drU8//dTp2Lm5uYqMjFSVKlW0YsUKNWvWTO7u7kpKStLnn3/udLkqVaoU2nbpX7clMWjQIH311VeaOHGiBg4c6HCfn5+f0tPTCyyTlpYm6eL2+f33369ovc5cybZ0ZR35f63//a/6s2fPqlq1apIunjNS1MdM8/v98MMPOnr0qHr06OFwf05OjiwWi/z8/IpdV2BgoKZNmyabzaZ9+/bp888/1yeffKKbb77ZvjcpX34gTU9PV7169ezt+SHlahT3eTt8+PAiz+0o7JDa5Tj74wAoDsIIyq2QkBB5eXlpw4YNuuuuu+ztW7dudei3a9cuPfroow5/6X3zzTeS/vtLzWKxOCzz119/6dChQxozZoyCgoKcLlfaLBaLXnvtNfXo0UNz5851uO+uu+7S1q1bHY7R22w2rVu3TkFBQaVyPY/ibMu/q1evnrZt26a8vLwS7/r/u/zHefLkSXtbRkaGkpOT7Yd3PDw8HObMmR07dmjGjBlq2bKl/Zdvbm6uNm7caA9axbFhwwZNmDBBa9askb+/v0JCQhQSEqJ169bp+PHjBfq3atVKFotFmzZtUv/+/e3tX331VbHWV5TiPm8DAgIUEBBw1eu7VGpqqiwWS4n2/gH5CCMot6pVq6ahQ4fq7bfflre3t9q0aaNt27YVCCMtWrTQmjVrdPvtt+vGG2/U7t27NW/ePLm5uSk7O1uS7CcZbt++XY0aNVJwcLDq16+vpUuX6sYbb5Svr6++/fZbLV68WJLsy10LjRs31rBhw/T22287tA8fPlzffPON+vXrZ/9reMmSJTp69GiRH8O8GsXZln93zz33aN68eTpw4ICaNm16Vetv2rSp6tatq1mzZsnHx0dubm6aO3dukYeWnOndu7eWLVumZ555Rv/v//0/eXt76+OPP9aBAwccPqqclJQkq9VqP7H57+68807l5uZq2LBhioyMVLVq1fTll1/qzJkzhX6C6KabblLPnj315ptv6vz587rtttu0adMm+/M2/xDflahZs6ax5+2uXbsUGhp6RXMBcJ0RlGuDBw/WmDFjtGHDBg0ZMkSJiYkaPXq0Q5/XXntNwcHBmjhxooYNG6YtW7bolVdeUbt27fTDDz9IuvgX95NPPqnNmzfr6aef1vnz5zV79mzVqVNH0dHRevbZZ7V3717NmTNHgYGB9uWulaefflq33367Q1vjxo318ccfq2bNmoqJidGoUaOUl5enxYsXq23btqVSR3G25d+FhoaqZs2a2rZt21Wv32Kx6N1331WtWrX0/PPPa/LkyerWrVuhv/Qvp1atWlq6dKmaNm2qSZMm6dlnn1V2drYWLVrkcILoK6+8Yv+0VmFq166t+fPn6/rrr9fYsWM1ePBg/fTTT5oxY4batGlT6DIvv/yyevfurYULF2ro0KE6efKkhgwZIunqz5kx8bzNyclRXFxcgZNygeJyy+ObwQCUsoULF+qTTz7RV199ddWHaq41q9WqiIgIrV271iXjnTp1St98843Cw8Mdzq+YOnWqVq1aVeCCbuXB6tWr9cYbb2jz5s1cgRVXhD0jAEpdnz59lJubqw0bNpgupcTmz59f4k+WFMXb21uTJ0/Wc889p61btyouLk5z587VkiVL7FeHLU9yc3O1cOFCDR8+nCCCK8aeEQDXxO7duxUdHa21a9eWqS/Lu5zExEQ1atRI7u6uO8Xul19+0dtvv634+HhlZ2fr5ptvVu/evfX444+Xuz1HsbGx2rBhgxYsWGC6FJRjhBEAAGAUh2kAAIBRhBEAAGAUYQQAABhVLi56duHCBWVkZMjT0/OqLggEAACundzcXOXk5MjPz6/Ik8DLRRjJyMjQ4cOHTZcBAACuQIMGDQr9Msh85SKMeHp6Srr4YFx9qWGbzaYDBw6oSZMmBb6fBNce81G2MB9lC/NRtjAfl5edna3Dhw/bf487Uy7CSP6hGW9v76u+VPLf2Ww2SRcvwcyTyTzmo2xhPsoW5qNsYT6K73KnWHACBgAAMIowAgAAjCKMAAAAowgjAADAKMIIAAAwijACAACMIowAAACjCCMAAMAowggAADCKMAIAAIwijAAAAKMIIwAAwCjCCAAAMIowAgAAjHI3XQCAsq1n7EkpdoPpMq7a4de6mS4BgBPsGQEAAEYRRgAAgFGEEQAAYBRhBAAAGHXFJ7BarVZFRETo5ZdfVlhYmKKjo/XZZ58V6BcWFqbFixcXaM/IyFDr1q0d2qpXr664uLgrLQkAAJRDVxRGcnJyNHLkSB08eNDeNnbsWI0cOdJ++7ffflPfvn3Vr1+/QsdISkpS9erVtXbtWnvbddexowYAgMqmxGEkKSlJI0eOVF5enkP79ddfr+uvv95+Ozo6Wl26dFGnTp0KHSclJUUNGzaUv79/SUsAAAAVSIl3RezcuVNhYWFavny50z7bt2/X999/r+eff95pn6SkJDVo0KCkqwcAABVMifeM9OnT57J95s2bpx49eqhu3bpO+yQnJ+vChQt6+OGHlZqaqtDQUMXExKh27dpOl7HZbLLZbCUtuUj547l6XFwZ5qNsqUjzUBEeC6+PsoX5uLzibhuXX4H16NGj2rFjh8aOHVtkv5SUFNWoUUMxMTHKy8vTW2+9pWeeeUaxsbGyWCyFLnPgwAFXl2uXkJBQamOj5JgPuFp8fLzpElyG10fZwnxcPZeHkY0bN6pZs2a69dZbi+y3bt06ubm5ycvLS5L07rvvql27dtq7d6/uvPPOQpdp0qSJqlat6tJ6bTabEhISFBQU5DQE4dphPsoWm80mxZ40XYZLtGzZ0nQJV43XR9nCfFxeVlZWsXYkuDyMfPvtt7rvvvsu28/b29vhds2aNVW9enWlpqY6XcZisZTahJfm2Cg55gOuVpGeT7w+yhbmw7nibheXfpY2Ly9PCQkJTvds5MvMzNRdd92lHTt22NtSU1P1119/KTAw0JUlAQCAMs6lYeS3337T2bNnCz1Ec+7cOaWlpUmSfHx81KpVK02ZMkX79u3TTz/9pOeee07h4eFq2rSpK0sCAABlnEvDyB9//CFJ8vPzK3Df+vXr1a5dO/vtqVOnqnnz5oqMjFTfvn1Vv359vfHGG64sBwAAlANXdc5IYmKiw+3g4OACbfkiIiIUERFhv+3n56cpU6ZczeoBAEAFwPXXAQCAUYQRAABgFGEEAAAYRRgBAABGEUYAAIBRhBEAAGAUYQQAABhFGAEAAEYRRgAAgFGEEQAAYBRhBAAAGEUYAQAARhFGAACAUYQRAABgFGEEAAAYRRgBAABGEUYAAIBRhBEAAGAUYQQAABhFGAEAAEYRRgAAgFGEEQAAYBRhBAAAGEUYAQAARhFGAACAUYQRAABgFGEEAAAYRRgBAABGEUYAAIBRhBEAAGAUYQQAABhFGAEAAEYRRgAAgFGEEQAAYBRhBAAAGEUYAQAARhFGAACAUVccRqxWq7p37664uDh726RJk9S0aVOHnyVLljgdY9GiRQoPD1dISIjGjBmj7OzsKy0HAACUU+5XslBOTo5GjhypgwcPOrQnJydr5MiR6tGjh73Nx8en0DE2btyomTNnatq0aapZs6ZiYmI0bdo0jRs37kpKAgAA5VSJ94wkJSXpkUce0ZEjRwrcl5ycrObNm8vf39/+4+3tXeg4ixcvVv/+/dWhQwe1aNFCr7zyilauXMneEQAAKpkSh5GdO3cqLCxMy5cvd2jPzMxUamqqGjRocNkxbDabEhISFBoaam9r2bKlzp8/r/3795e0JAAAUI6V+DBNnz59Cm1PTk6Wm5ub3nvvPX3zzTeqXr26nnzySYdDNvlOnz6tnJwc1a5d+7+FuLurevXqOnnypNN122w22Wy2kpZcpPzxXD0urgzzUbZUpHmoCI+F10fZwnxcXnG3zRWdM1KYlJQUubm5KTAwUE888YS+//57vfzyy/Lx8dH999/v0PfcuXOSJA8PD4d2Dw8PWa1Wp+s4cOCAq8otICEhodTGRskxH3C1RmM3mC7BJVb2upHXRxnDfFw9l4WRhx56SB06dFD16tUlSbfddpsOHz6sTz75pEAY8fT0lKQCwcNqtTo9x0SSmjRpoqpVq7qqZEn/PWQUFBQki8Xi0rFRcsxH2WKz2aRY53srYQavj7KB96vLy8rKKtaOBJeFETc3N3sQyRcYGKgdO3YU6Fu9enV5enoqPT1djRo1kiRduHBBp06dkr+/v9N1WCyWUpvw0hwbJcd8AM7x+ihbmA/nirtdXHbRs3feeUcDBgxwaNu/f78CAwMLrvS66xQUFKRdu3bZ2+Lj4+Xu7q7bbrvNVSUBAIBywGVhpEOHDvr++++1YMECHTlyRB9//LFWr16tp556StLF80TS0tLs/fv06aMFCxZo8+bN2rdvnyZMmKBHHnmkyMM0AACg4nHZYZoWLVronXfe0bvvvqt33nlH9evX1/Tp0xUSEiJJWr9+vWJiYpSYmChJ6tatm3777TeNGzdOVqtVDzzwgEaNGuWqcgAAQDlxVWEkP1jk69Spkzp16lRo34iICEVERDi0RUZGKjIy8mpKAAAA5RxflAcAAIwijAAAAKMIIwAAwCjCCAAAMIowAgAAjCKMAAAAowgjAADAKMIIAAAwijACAACMIowAAACjCCMAAMAowggAADCKMAIAAIwijAAAAKMIIwAAwCjCCAAAMIowAgAAjCKMAAAAowgjAADAKMIIAAAwijACAACMIowAAACjCCMAAMAowggAADCKMAIAAIwijAAAAKMIIwAAwCjCCAAAMIowAgAAjCKMAAAAowgjAADAKMIIAAAwijACAACMIowAAACjCCMAAMAowggAADCKMAIAAIy64jBitVrVvXt3xcXF2dvi4+PVu3dvhYSEqHPnzoqNjS1yjNDQUDVt2tTh5+zZs1daEgAAKIfcr2ShnJwcjRw5UgcPHrS3paWl6emnn9Zjjz2m1157TT/99JNiYmLk7++v9u3bFxgjNTVVZ86c0ebNm+Xl5WVvr1q16pWUBAAAyqkSh5GkpCSNHDlSeXl5Du2bN29WrVq19Pzzz0uSGjRooLi4OK1Zs6bQMJKcnCx/f3/ddNNNV1Y5AACoEEocRnbu3KmwsDA999xzatmypb09PDxczZo1K9A/MzOz0HGSkpLUsGHDkq4eAABUMCUOI3369Cm0PSAgQAEBAfbbf/zxh9atW6cRI0YU2j85OVnZ2dnq27evDh06pGbNmmnMmDFFBhSbzSabzVbSkouUP56rx8WVYT7KFuahbGJeygbery6vuNvmis4ZuZxz585pxIgRqlWrlh599NFC+6SkpCgjI0PPP/+8fHx89P7772vAgAFat26dfHx8Cl3mwIEDpVGuJCkhIaHUxkbJMR+Ac7w+yhbm4+q5PIycPXtWQ4cO1eHDh/Xxxx/L29u70H4LFizQ+fPnVa1aNUnSG2+8oXvvvVdbt27Vv/71r0KXadKkictPcLXZbEpISFBQUJAsFotLx0bJMR9li81mk2JPmi4Df8Pro2zg/erysrKyirUjwaVhJDMzU4MGDdKRI0f04YcfqkGDBk77enh4yMPDw37b09NTAQEBSk1NdbqMxWIptQkvzbFRcswH4Byvj7KF+XCuuNvFZRc9y83N1fDhw3Xs2DF99NFHaty4sdO+eXl56tSpk1atWmVvy8rK0q+//qrAwEBXlQQAAMoBl+0ZWbFiheLi4jRnzhz5+voqLS1NklSlShVVr15dVqtVGRkZqlGjhiwWi9q3b68ZM2aofv36qlGjht555x3deOONuvfee11VEgAAKAdcFkY2btyo3NxcDR482KG9devW+uijj7Rnzx7169dPW7ZsUUBAgEaNGiV3d3eNHDlSmZmZatOmjebNm8euLgAAKpmrCiOJiYn2/y9YsKDIvmFhYQ79PT09FR0drejo6KspAQAAlHN8UR4AADCKMAIAAIwqlYueAZAaRK8zXQIAlAvsGQEAAEYRRgAAgFGEEQAAYBRhBAAAGEUYAQAARhFGAACAUYQRAABgFGEEAAAYRRgBAABGEUYAAIBRhBEAAGAUYQQAABhFGAEAAEYRRgAAgFGEEQAAYBRhBAAAGEUYAQAARhFGAACAUYQRAABgFGEEAAAYRRgBAABGEUYAAIBRhBEAAGAUYQQAABhFGAEAAEYRRgAAgFGEEQAAYBRhBAAAGEUYAQAARhFGAACAUYQRAABgFGEEAAAYRRgBAABGEUYAAIBRVxxGrFarunfvrri4OHvb0aNHNWDAALVs2VL//Oc/9d133xU5xtq1a9WpUycFBwdr2LBh+vPPP6+0HAAAUE5dURjJycnR888/r4MHD9rb8vLyNGzYMNWqVUsrV67Ugw8+qOHDh+v48eOFjrFv3z6NHTtWw4cP1/Lly3X69GnFxMRc2aMAAADllntJF0hKStLIkSOVl5fn0L5jxw4dPXpUy5YtU9WqVdWoUSNt375dK1eu1IgRIwqMs2TJEnXt2lUPPfSQJOn1119Xhw4ddPToUd10001X9mgAAEC5U+I9Izt37lRYWJiWL1/u0L537141b95cVatWtbe1atVK8fHxhY6zd+9ehYaG2m/XrVtX9erV0969e0taEgAAKMdKvGekT58+hbanpaWpdu3aDm01a9bUyZMnC+3/+++/l6i/JNlsNtlsthJWXLT88Vw9Lq4M8wFcHq+PsoH3q8sr7rYpcRhxJjs7Wx4eHg5tHh4eslqthfY/d+5cifpL0oEDB66+UCcSEhJKbWyUHPMBOMfro2xhPq6ey8KIp6enTp065dBmtVrl5eXltP/fg4fVapW3t7fTdTRp0sThMJAr2Gw2JSQkKCgoSBaLxaVjo+Qq1HzEbjBdASqoCvH6qAAq1PtVKcnKyirWjgSXhZE6deooKSnJoS09Pb3AoZhL+6enpxfo7+/v73QdFoul1Ca8NMdGyTEfgHO8PsoW5sO54m4Xl130LDg4WD/99JPOnTtnb9u1a5eCg4Od9t+1a5f99okTJ3TixAmn/QEAQMXksjDSunVr1a1bVzExMTp48KDmzZunffv26eGHH5Z08RBMWlqa/WSWxx57TJ9//rliY2O1f/9+vfjii2rfvj0f6wUAoJJxWRixWCyaPXu20tLSFBERoS+++EKzZs1SvXr1JEl79uxRu3btdOLECUlSSEiIXn31Vc2aNUuPPfaY/Pz8NGXKFFeVAwAAyomrOmckMTHR4fYtt9yiJUuWFNo3LCysQP+IiAhFRERcTQkAAKCc44vyAACAUYQRAABgFGEEAAAYRRgBAABGEUYAAIBRhBEAAGAUYQQAABhFGAEAAEYRRgAAgFGEEQAAYBRhBAAAGEUYAQAARhFGAACAUYQRAABgFGEEAAAYRRgBAABGEUYAAIBRhBEAAGAUYQQAABhFGAEAAEYRRgAAgFGEEQAAYBRhBAAAGEUYAQAARhFGAACAUYQRAABgFGEEAAAYRRgBAABGEUYAAIBRhBEAAGAUYQQAABhFGAEAAEYRRgAAgFGEEQAAYBRhBAAAGEUYAQAARhFGAACAUe6uHGzVqlWKiYkp0O7m5qb9+/cXaP/3v/+txMREh7Y1a9aoSZMmriwLAACUYS4NI//85z8VHh5uv33hwgX1799f7du3L9DXZrPp8OHDWrJkiRo0aGBvv+GGG1xZEgAAKONcGka8vLzk5eVlvz137lzl5eXphRdeKND32LFjOn/+vFq0aCFPT09XlgEAAMoRl4aRS506dUrvv/++Jk2aJA8PjwL3JyUlqW7duiUKIjabTTabzZVl2sdz9bi4MswHcHm8PsoG3q8ur7jbptTCyCeffKLatWurS5cuhd6fnJysKlWqaPDgwfrxxx/VsGFDvfjii2rRooXTMQ8cOFBa5SohIaHUxkbJMR+Ac7w+yhbm4+qVShjJy8tTbGysBg0a5LTPoUOHlJGRoV69eikqKkqffvqp+vfvr/Xr16tu3bqFLtOkSRNVrVrVpbXabDYlJCQoKChIFovFpWOj5CrUfMRuMF0BKqgK8fqoACrU+1UpycrKKtaOhFIJIwkJCUpNTVW3bt2c9pk4caLOnTsnHx8fSdKECRO0e/duff7553rmmWcKXcZisZTahJfm2Cg55gNwjtdH2cJ8OFfc7VIq1xn59ttvFRoaKj8/P6d93N3d7UFEuvjx38DAQKWmppZGSQAAoIwqlTCyb98+3XnnnUX26du3r2bOnGm/nZubq8TERAUGBpZGSQAAoIwqlTBy8OBB3XrrrQ5tNptNaWlpslqtkqSOHTtq0aJF2rJli1JSUvTqq6/qzJkz6tGjR2mUBAAAyqhSOWckPT1dvr6+Dm0nTpzQfffdp8WLFyssLEwDBgxQTk6OJk2apPT0dAUHB+uDDz5wOHQDAAAqvlIJI/v27SvQFhAQ4HDpdzc3Nz3zzDNOT1YFAACVA1+UBwAAjCKMAAAAowgjAADAKMIIAAAwijACAACMIowAAACjCCMAAMAowggAADCKMAIAAIwijAAAAKMIIwAAwCjCCAAAMIowAgAAjCKMAAAAowgjAADAKMIIAAAwijACAACMIowAAACjCCMAAMAowggAADCKMAIAAIwijAAAAKMIIwAAwCjCCAAAMIowAgAAjCKMAAAAowgjAADAKMIIAAAwijACAACMIowAAACjCCMAAMAowggAADCKMAIAAIwijAAAAKMIIwAAwCjCCAAAMIowAgAAjHJpGNm0aZOaNm3q8BMVFVVo3//85z/q3r27goOD1a9fPx09etSVpQAAgHLC3ZWDJSUlqUOHDpo4caK9zdPTs0C/48ePa9iwYRoxYoTCw8M1a9YsDR06VF988YXc3NxcWRIAACjjXBpGkpOT1aRJE/n7+xfZLzY2VnfccYeeeuopSdKUKVN0zz33aOfOnQoLC3NlSQAAoIxzeRhp27btZfvt3btXoaGh9tve3t66/fbbFR8fX2QYsdlsstlsLqn10jEv/RdmMR/A5fH6KBt4v7q84m4bl4WRvLw8HTp0SN99953mzp0rm82mLl26KCoqSh4eHg5909LSVLt2bYe2mjVr6uTJk0Wu48CBA64qt4CEhIRSGxslx3wAhesZe1KKLfq9sjxY2etG0yW4DO9XV89lYeT48ePKzs6Wh4eH3n77bR07dkyTJk3SuXPn9NJLLzn0ze93KQ8PD1mt1iLX0aRJE1WtWtVVJUu6mNoSEhIUFBQki8Xi0rFRchVqPmI3mK4AKLNatmxpuoSrVqHer0pJVlZWsXYkuCyM1K9fX3FxcfLz85Obm5uaNWum3NxcjRo1SjExMQ4T5enpWSB4WK1W+fr6FrkOi8VSahNemmOj5JgPoGKrSK9v3q+cK+52celHe6tXr+7waZhGjRopJydHGRkZDv3q1Kmj9PR0h7b09PTLnvgKAAAqHpeFkW+//VZhYWHKzs62t/3yyy+qXr26atSo4dA3ODhYu3btst/Ozs7Wzz//rODgYFeVAwAAygmXhZGQkBB5enrqpZdeUkpKirZt26bXX39dgwYNks1mU1pamv3QTM+ePbV7927NmzdPBw8eVExMjAICAvhYLwAAlZDLwoiPj48WLFigP//8Uz179tTYsWP16KOPatCgQTpx4oTatWunPXv2SJICAgI0Y8YMrVy5Ug8//LBOnTqlWbNmccEzAAAqIZdeZ6Rx48b64IMPCrQHBAQoMTHRoe3ee+/Vvffe68rVAwCAcogvygMAAEYRRgAAgFGEEQAAYBRhBAAAGEUYAQAARhFGAACAUYQRAABgFGEEAAAYRRgBAABGEUYAAIBRhBEAAGAUYQQAABhFGAEAAEYRRgAAgFGEEQAAYBRhBAAAGEUYAQAARhFGAACAUYQRAABgFGEEAAAYRRgBAABGEUYAAIBRhBEAAGAUYQQAABhFGAEAAEYRRgAAgFGEEQAAYBRhBAAAGEUYAQAARhFGAACAUYQRAABgFGEEAAAYRRgBAABGEUYAAIBRhBEAAGAUYQQAABhFGAEAAEa5NIykpqYqKipKrVu3Vnh4uKZMmaKcnJxC+w4ZMkRNmzZ1+Nm6dasrywEAAOWAu6sGysvLU1RUlHx9fbV06VJlZGRozJgxuu666zR69OgC/ZOTkzVt2jTdfffd9jY/Pz9XlQMAAMoJl4WRlJQUxcfH6//+7/9Uq1YtSVJUVJSmTp1aIIxYrVYdO3ZMQUFB8vf3d1UJAACgHHJZGPH399f8+fPtQSRfZmZmgb4pKSlyc3PTTTfdVKJ12Gw22Wy2q6qzsDEv/RdmMR9A5VARXuO8X11ecbeNy8KIr6+vwsPD7bdzc3O1ZMkStWnTpkDflJQU+fj46MUXX9TOnTt14403asSIEbr33nuLXMeBAwdcVW4BCQkJpTY2So75ACq2+Ph40yW4DO9XV89lYeTvpk2bpp9//lkrVqwocF9KSorOnTundu3aKTIyUps2bdKQIUO0fPlyBQUFOR2zSZMmqlq1qkvrtNlsSkhIUFBQkCwWi0vHRslVqPmI3WC6AqDMatmypekSrlqFer8qJVlZWcXakVAqYWTatGn68MMP9dZbb6lJkyYF7h86dKj69u1rP2H1tttu008//aRPP/20yDBisVhKbcJLc2yUHPMBVGwV6fXN+5Vzxd0uLr/OyMSJE/XBBx9o2rRp6ty5c+Erve66Ap+cCQwMVGpqqqvLAQAAZZxLw8jMmTO1bNkyvfnmm+rWrZvTftHR0YqJiXFo279/vwIDA11ZDgAAKAdcFkaSk5M1e/ZsPf3002rVqpXS0tLsP5KUlpamc+fOSZI6duyoNWvWaPXq1fr11181c+ZM7dq1S0888YSrygEAAOWEy84Z2bJli2w2m+bMmaM5c+Y43JeYmKh27dppypQpioiI0AMPPKDx48drzpw5On78uBo3bqz58+crICDAVeUAAIBywmVhJDIyUpGRkU7vT0xMdLjdq1cv9erVy1WrBwAA5RRflAcAAIwijAAAAKMIIwAAwCjCCAAAMIowAgAAjCKMAAAAowgjAADAKMIIAAAwijACAACMIowAAACjCCMAAMAowggAADCKMAIAAIwijAAAAKMIIwAAwCjCCAAAMIowAgAAjCKMAAAAowgjAADAKMIIAAAwijACAACMIowAAACjCCMAAMAowggAADCKMAIAAIwijAAAAKMIIwAAwCjCCAAAMIowAgAAjCKMAAAAowgjAADAKMIIAAAwijACAACMIowAAACjCCMAAMAowggAADCKMAIAAIxyaRjJycnRmDFjFBoaqnbt2mnhwoVO+/7888/q1auXgoOD1bNnT/3444+uLAUAAJQTLg0jr7/+un788Ud9+OGHGj9+vGbOnKkNGzYU6JeVlaXIyEiFhoZq1apVCgkJ0eDBg5WVleXKcgAAQDngsjCSlZWl2NhYjR07Vrfffrvuv/9+DRo0SEuXLi3Qd/369fL09NSLL76oRo0aaezYsapWrVqhwQUAAFRs7q4aaP/+/bpw4YJCQkLsba1atdJ7772n3NxcXXfdf3PP3r171apVK7m5uUmS3NzcdOeddyo+Pl4REREFxs7NzZUknT17VjabzVUlO4ydmZnpUCPMqEjz0bC6y15eQIVz5swZ0yVctYr0flVazp07J+m/28oZl71bpqWl6YYbbpCHh4e9rVatWsrJydGpU6dUo0YNh7633nqrw/I1a9bUwYMHCx07JydHknTkyBFXlVtAUlJSqY2NkqsI8/HG/bVMlwCUWQcOHDBdgstUhPer0paTkyMfHx+n97ssjGRnZzsEEUn221artVh9/94vn5+fnxo0aCBPT0/SJwAA5URubq5ycnLk5+dXZD+XhRFPT88CYSL/tpeXV7H6/r2fvUh3d9WsWdNVpQIAgGukqD0i+Vy2m6FOnTr666+/dOHCBXtbWlqavLy85OvrW6Bvenq6Q1t6erpq167tqnIAAEA54bIw0qxZM7m7uys+Pt7etmvXLgUFBRU4tBIcHKw9e/YoLy9PkpSXl6fdu3crODjYVeUAAIBywmVhxNvbWw899JAmTJigffv2afPmzVq4cKH69esn6eJekvyzart06aLTp09r8uTJSkpK0uTJk5Wdna2uXbu6qhwAAFBOuPRs0JiYGN1+++3q37+/XnnlFY0YMUIPPPCAJKldu3Zav369pIvHj+bOnatdu3YpIiJCe/fu1bx581S1alVXllOkklwtFqUvNTVVUVFRat26tcLDwzVlyhT7p6hgVmRkpKKjo02XUelZrVa98soruuuuu9S2bVu9+eab9r3LuPZOnDihwYMH684771THjh21aNEi0yWVay69EIK3t7emTp2qqVOnFrgvMTHR4XaLFi302WefuXL1JXLp1WKPHz+u0aNHq169eurSpYuxmiqrvLw8RUVFydfXV0uXLlVGRobGjBmj6667TqNHjzZdXqW2bt06bdu2TT169DBdSqU3adIkxcXFacGCBTp79qyee+451atXT7179zZdWqX07LPPql69elq1apWSkpL0wgsvqH79+rr//vtNl1YuVcrPyZbkarEofSkpKYqPj9eUKVPUuHFjhYaGKioqSmvXrjVdWqV26tQpvf766woKCjJdSqV36tQprVy5UhMnTlSLFi10991366mnntLevXtNl1YpZWRkKD4+XkOGDFGDBg3UqVMnhYeHa/v27aZLK7cqZRhxdrXYvXv3XvYqcXA9f39/zZ8/X7VqOV4kLDMz01BFkKSpU6fqwQcfLHCBQlx7u3btko+Pj1q3bm1vi4yM1JQpUwxWVXl5eXnJ29tbq1at0vnz55WSkqLdu3erWbNmpksrtyplGLnc1WJxbfn6+io8PNx+Ozc3V0uWLFGbNm0MVlW5bd++XT/88IOGDh1quhRIOnr0qOrXr6/Vq1erS5cuuu+++zRr1iz+eDLE09NT48aN0/LlyxUcHKyuXbvqH//4h3r16mW6tHKrUn55RkmuFotrb9q0afr555+1YsUK06VUSjk5ORo/frzGjRvn9EKEuLaysrL066+/atmyZZoyZYrS0tI0btw4eXt766mnnjJdXqWUnJysDh066Mknn9TBgwc1ceJE3X333fr3v/9turRyqVKGkZJcLRbX1rRp0/Thhx/qrbfeUpMmTUyXUynNnDlTd9xxh8PeKpjl7u6uzMxMTZ8+XfXr15ckHT9+XJ988glhxIDt27drxYoV2rZtm7y8vBQUFKTU1FTNmTOHMHKFKmUYufRqse7uFzeBs6vF4tqZOHGiPvnkE02bNk2dO3c2XU6ltW7dOqWnp9vPqcoP6hs3btSePXtMllZp+fv7y9PT0x5EJKlhw4Y6ceKEwaoqrx9//FG33HKLwx+vzZs313vvvWewqvKtUoaRS68WGxoaKsn51WJxbcycOVPLli3Tm2++ycerDfvoo48cvtbhjTfekCS98MILpkqq9IKDg5WTk6NDhw6pYcOGki5+Cu3ScIJrp3bt2vr1119ltVrth/hTUlIUEBBguLLyq1L+5r3c1WJxbSUnJ2v27Nl6+umn1apVK6Wlpdl/cO3Vr19ft9xyi/2nWrVqqlatmm655RbTpVVagYGBat++vWJiYrR//359++23mjdvnh577DHTpVVKHTt2VJUqVfTSSy/p0KFD+t///V+999576tu3r+nSyi23vEp6Cb/s7GxNmDBBX331lXx8fDRw4EANGDDAdFmV0rx58zR9+vRC7/v7xfJw7eVfffW1114zXEnldubMGU2cOFGbNm2St7e3+vTpo2HDhsnNzc10aZVS/leZ7Nu3TzVq1NDjjz+u/v37Mx9XqNKGEQAAUDZUysM0AACg7CCMAAAAowgjAADAKMIIAAAwijACAACMIowAAACjCCMAAMAowggAADCKMAIAAIwijAAAAKMIIwAAwKj/D6HklMyXJ8rPAAAAAElFTkSuQmCC",
"text/plain": [
"<Figure size 640x480 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"data = np.random.randn(50) + 5\n",
"plt.hist(data, bins=np.arange(10))\n",
"plt.title(\"data ~ Normal(mu=5, sigma=1)\")\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Model this data using an _Exponential_ distribution for it's mean. (This may seem\n",
"slightly counter intuitive.)"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"Sampling: [likelihood, mu, sigma]\n"
]
}
],
"source": [
"with pm.Model() as model:\n",
" mu = pm.Exponential(\"mu\", lam=1)\n",
" sigma = pm.Exponential(\"sigma\", lam=1)\n",
" pm.Normal(\"likelihood\", mu=mu, sigma=sigma, observed=data)\n",
" idata = pm.sample_prior_predictive()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"If we look at the prior distribution of mu, it looks like Exponential(1):"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAhwAAAGsCAYAAACW3H6UAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAtsElEQVR4nO3dfVxUdd7/8Tc3AYOkohil9dNLTaNpVMJyK3qkljdtmZaPrtwe3pSZtRrs1c2GYKtgNwq2diOlqdHqxlaLmt1obnlV1tVaJhsCraBpWl6m4ZWYCgwx8PvDB1OnAeUoXweY1/Px8LHNd875zud8wJ2333PmTFBdXV2dAAAADAr2dwEAAKDtI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAaDV4D6FQOtF4EDAmTBhgvr27Wv5c8kll2jw4MHKyMjQ4cOHT7j/3r171bdvX61evfqM1Dt06FCfen/55/777z8jdfjTjz/+qIcfflhbtmzxjk2YMEETJkywNc/QoUM1Y8YMn/F7771XeXl5pzzvmfT++++rb9++PuPjx4/XunXr/FAR0DSh/i4A8IeLL75Ys2fP9j7+6aef9OWXX2rBggXatm2bXnnlFQUFBTW47znnnKPXXntN/+///b8zVa6uueYaTZs2rcHnoqOjz1gd/rJt2za98cYbGjt2rHfslz+/07F69WodOHDAMndL9dlnn+nBBx9s8Lm0tDTdddddGjRokDp37nyGKwNOjsCBgBQVFaUBAwZYxi677DIdO3ZMzz77rLZu3erzfL2wsLBGnzOlU6dOZ/w1W7revXuf9hxVVVV68sknNXv2bAUHt9wF36NHj2rp0qVaunSpzj77bFVUVPhsc/HFF6tfv35atGiRHnnkET9UCZxYy/0bBvjBJZdcIknat2+fpOPL6w899JCSk5M1YMAA3XnnnQ2eUtm9e7eSk5N11VVXacCAAZowYYLy8/O9z9fv89JLL2nkyJHq37+/Vq1a1ay1r1ixwqeuTz/9VBdddJGee+45SdKMGTM0YcIErVy5UkOGDFF8fLwmTZqkkpISy1xNPZ533nlHycnJio+P1+WXX65HHnnE580wLy9PN9xwg/e01cKFC+XxeLzPz5gxQ3fccYdWrVqlESNG6JJLLtHo0aP10UcfSTr+r/qJEydKkiZOnOg93fHrUx8//PCDMjIyNGTIEF1yySW6/PLLNX36dO3du7fRnq1atUput1tDhgxpdJumzDthwgTNmjVLzz//vK6++mr1799fd999tw4ePKhVq1Zp2LBhio+P1x133GHZ70SnyoYOHerdbuXKlfr73/+uWbNmafz48Y3WOmrUKK1cuVI//PBDo9sA/sIKB/ALX3/9tSTpggsu8I698847uummm7Ro0SLV1tb67PPVV1/pP//zP9WjRw898sgjOuuss7RixQpNmjRJOTk5uvzyy73bLly4UDNnzlRUVJT69+/f5Lrq6upUU1PT4HOhocf/Gk+YMEHvvvuuMjMzNXjwYIWFhSktLU0DBgzQvffe691+27Zt2rVrlx544AF16NBBzz77rPf8/znnnGPreGbPnq2xY8fq+eefV2FhoZ566ilFR0d7l/1feOEFPfXUUxo/frxSU1O1bds2LVy4UN99952eeOIJ7zzFxcX6/vvvlZycrKioKD3zzDNKSkrSRx99JKfTqVmzZmnOnDmaNWuWBg0a1GB/7rnnHh0+fFgPPfSQYmJiVFpaqqefflqzZ8/Wiy++2GDv3nzzTW+vGut7U+d9++235XQ69fjjj2v//v2aM2eOxo8fr/DwcKWkpKiystJ7HEuWLJEkvfbaaw2+riRLTUOHDtVtt90mh8OhhQsXNrrP0KFD5fF49N577+m2225rdDvAHwgcCEi/fgM/fPiwNm/erEWLFik+Pt670iFJZ511ljIyMrxvAL/+F3N2drbCwsK0YsUKRUVFSZIGDx6sG2+8UVlZWVq5cqV32+uvv/6UrhVYs2aN1qxZ0+BzK1eulMvlUlBQkObOnaubbrpJ8+fPV0hIiMrLy7V8+XKFhIR4tz9y5IgWL16sgQMHSpL69eun6667TitWrNBDDz1k63iuueYapaSkSJKuuOIKffLJJ/rwww/14IMP6siRI3r++ed12223eZf4ExMT1bFjRz3yyCO68847deGFF3prWr16tfe6mMjISI0fP16ffvqpRowY4T190rt37wZPpXz//fdyOBxKSUnxHtegQYP0zTffNPqmfvToURUVFen6669vtO925q2pqVF2drY6dOggSXr33Xf18ccfa8OGDd4AW1BQoDfeeMO7T1NPkzX1eqHIyEj16tVLmzZtInCgxSFwICB9/vnncjqdlrHg4GBdeeWVmjNnjuWC0Z49ezb6L2BJ2rx5s4YMGeJ9c5aOrzrccMMNeu6553Ts2DHveFxc3CnVO2TIEE2fPr3B5375BnzBBRfooYce0qOPPqq6ujrNnTvXslojSeeff773zVM6fhFsfHy8Pv/8c9vH8+s3zHPPPVf/+7//K0n64osvVFVVpaFDh1rCXf2pgk8++cQbODp16mR5Uz333HMlSZWVlSfpzHGxsbFasWKF6urqtHfvXu3Zs0e7du3Sv/71L1VXVze4z3fffSePx6Pzzz+/Webt1auXN2xIUkxMjKKjoy3979ixo44cOeJ93NiqlSQFBQVZgmJTdevW7YSnkQB/IXAgIDmdTmVkZEg6/n/s4eHhOu+88yxvsvXatWt3wrkOHz6smJgYn/GYmBjV1dXp6NGj3rHIyMhTqrdjx45yuVxN2va3v/2t5s2bJ0m66qqrfJ6PjY31GevcubO+/PJLSfaOx+FwWLYJDg723iujvLxckjR16tQG6/z+++8bnac+8DV0Cqsxb775phYsWKDvvvtOHTt2VFxcnCIiIhrdvv6N/2Q/k6bO29Dvzsnm/nXo/aVu3brp/fffP+H+DXE4HJZQA7QUBA4EpHbt2jX5DfxkOnTooIMHD/qMl5WVSTr+sdVfvrma9thjj6ldu3YKCwvTrFmz9MILL1ieP3TokM8+Bw8e9H6UsrmOp3379pKkJ598Uj169PB5vqFQc6q2bNmilJQUTZgwQXfddZc3VGVlZVkudv2l+o8T//jjj806rx2/PD31aydaVTuRH3/8MSA+Ko3Wh8ABnKbLLrtMH3zwgY4ePer9V67H49HatWvlcrlO+Y3jVLz77rt6++23lZWVpXbt2mn69OlatWqV5bqR3bt3a+fOnerVq5ck6cCBA/riiy+8KxHNdTz9+/fXWWedpQMHDmjUqFHe8W3btikrK0vTpk3Teeed16S5TnZq4YsvvlBtba2SkpJ09tlne2v+5z//Ken4SsmvP/YaGxurkJAQ7d+/v1nntaO5Qu8v7d+/33uqCmhJCBzAabrvvvv00UcfaeLEiZo6darOOussvfzyy/r222+1bNmyE+5bUFDgc/1CQ3744QcVFBQ0+FxISIhcLpd++OEHpaenKzExUaNHj5YkXXfddZo7d66uuuoq73URdXV1uvfee3X//fcrJCTEe6Fj/UdMT+d4fik6OlpTpkzRM888o6NHj2rQoEE6cOCAnnnmGQUFBemiiy5q8lz1b/YffvihOnTo4LNvv379JElz5szR2LFjdfjwYeXm5no/7ltRUeFzyiMyMlKXXnqp8vPzdccddzT4uqcyrz8dOXJEO3bs0OTJk/1dCuCDwAGcpgsvvFB/+9vftGDBAqWmpiooKEj9+vXTihUrLBdnNuS2227TzTff7L3mojEbN27Uxo0bG3zu7LPP1pYtW5SRkaHKykrvtSmSNGvWLP32t7/VzJkzvR/h7Nq1qyZPnqwnnnhClZWVuvLKK7Vo0SJ17NjxtI/n1/7rv/5LXbp00d/+9jctW7ZMHTp00BVXXKEHHnjAGyKa4sILL9SNN96o3Nxcffzxx3r77bctzw8aNEizZs3SSy+9pPXr1ysmJkaDBg1Sdna2pk+frvz8fF1zzTU+844YMUILFy6U2+1WeHi4z/OnOq+/fPzxxzrrrLM0ePBgf5cC+Aiq49uQAL/ZtGmT3nnnHc2ZM+eMvN6MGTO0efPmU7oYsS2qrKzUddddpz/+8Y8aM2aMv8s5bZMmTVKfPn00c+ZMf5cC+OBOo4Cf1NbWatmyZQ1+kgRnhsPhUFJSkl588UXL3U9bo6KiIpWUlDT6qSDA3wgcgJ8EBwfrgQce0IgRI/xdSkAbN26czj33XO+3xbZWc+fO1Z/+9Cd16dLF36UADeKUCgAAMI4VDgAAYByBAwAAGEfgAAAAxrWI+3DU1NTo8OHDCg8PP6279gEAgDOntrZWbrdbHTp0UGjoiSNFiwgchw8f1u7du/1dBgAAOAU9evTwfh9TY1pE4Ki/w1+PHj18vjXydHk8Hm3fvl19+vQ5pa96bmvohxX98EVPrOiHFf3wFcg9qays1O7duxu8U++vtYjAUX8axeFwnPLXdzem/mY+kZGRAfeL0BD6YUU/fNETK/phRT980RM16XIILpgAAADGETgAAIBxBA4AAGAcgQMAABhH4AAAAMYROAAAgHEEDgAAYByBAwAAGEfgAAAAxhE4AACAcQQOAABgHIEDAAAYR+AAAADGETgAAIBxLeLr6U0bm7dfylvv7zJs2T3vBn+XAABAs2GFAwAAGEfgAAAAxhE4AACAcQQOAABgHIEDAAAYR+AAAADGETgAAIBxBA4AAGAcgQMAABhH4AAAAMYROAAAgHEEDgAAYByBAwAAGEfgAAAAxhE4AACAcQQOAABgHIEDAAAYR+AAAADGETgAAIBxBA4AAGAcgQMAABhH4AAAAMYROAAAgHEEDgAAYJztwOF2u5WWlqaBAwcqMTFROTk5jW5bWlqq3/3ud+rXr59GjRqlTz/99LSKBQAArZPtwJGVlaXi4mItX75cs2fPVnZ2ttavX++z3ZEjRzR58mT17t1bb731loYNG6b77rtP//d//9cshQMAgNbDVuCoqKhQXl6eZs6cKafTqWHDhmnKlCnKzc312fb1119XZGSk0tPT1b17dyUnJ6t79+4qLi5utuIBAEDrEGpn45KSEtXU1Cg+Pt47lpCQoMWLF6u2tlbBwT/nl82bN+vaa69VSEiId2zVqlUnnN/j8cjj8dgp6aSae74zxVTd9fO21r40N/rhi55Y0Q8r+uErkHti55htBY6ysjJFR0crLCzMOxYTEyO3263y8nJ16tTJO/7tt9+qX79++tOf/qT3339f3bp1U0pKihISEhqdf/v27XbKadMKCgqMzl9UVGR0/taGfviiJ1b0w4p++KInJ2YrcFRWVlrChiTv4+rqast4RUWFlixZookTJ2rp0qVau3at7rrrLr3zzjs677zzGpy/T58+ioyMtFPSSXk8Hilvf7POeSYMGDDAyLwej0dFRUVyuVyW1adART980RMr+mFFP3wFck8qKiqavFhgK3CEh4f7BIv6xxEREZbxkJAQxcXFKTk5WZJ08cUX65NPPtEbb7yhe++9t8H5Q0JCAu6H1RjTfaDXVvTDFz2xoh9W9MNXIPbEzvHaumg0NjZWhw4dUk1NjXesrKxMERERat++vWXbLl26qGfPnpaxHj166LvvvrPzkgAAoA2wFTji4uIUGhpqub4gPz9fLpfLcsGodPyUQGlpqWVs165d6tat26lXCwAAWiVbgcPhcGjMmDFKT09XYWGhNmzYoJycHE2cOFHS8dWOqqoqSdK4ceNUWlqqhQsXas+ePXrmmWf07bffavTo0c1/FAAAoEWzfeOv1NRUOZ1OTZo0SRkZGUpKStLw4cMlSYmJiVq3bp0kqVu3blq2bJk++OAD3Xjjjfrggw+0ZMkSxcbGNu8RAACAFs/WRaPS8VWOzMxMZWZm+jz361MoCQkJWr169alXBwAA2gS+vA0AABhH4AAAAMYROAAAgHEEDgAAYByBAwAAGEfgAAAAxhE4AACAcQQOAABgHIEDAAAYR+AAAADGETgAAIBxBA4AAGAcgQMAABhH4AAAAMYROAAAgHEEDgAAYByBAwAAGEfgAAAAxhE4AACAcQQOAABgHIEDAAAYR+AAAADGETgAAIBxBA4AAGAcgQMAABhH4AAAAMYROAAAgHEEDgAAYByBAwAAGEfgAAAAxhE4AACAcQQOAABgHIEDAAAYR+AAAADGETgAAIBxBA4AAGAcgQMAABhH4AAAAMYROAAAgHEEDgAAYByBAwAAGEfgAAAAxtkOHG63W2lpaRo4cKASExOVk5PT6La///3v1bdvX8ufDz744LQKBgAArU+o3R2ysrJUXFys5cuXa9++fUpJSVHXrl01cuRIn2137typ+fPn64orrvCOdejQ4fQqBgAArY6twFFRUaG8vDwtXbpUTqdTTqdTO3bsUG5urk/gqK6u1t69e+VyudSlS5dmLRoAALQutgJHSUmJampqFB8f7x1LSEjQ4sWLVVtbq+Dgn8/Q7Nq1S0FBQbrggguaPL/H45HH47FTUpPmbI1M1V0/b2vtS3OjH77oiRX9sKIfvgK5J3aO2VbgKCsrU3R0tMLCwrxjMTExcrvdKi8vV6dOnbzju3btUlRUlB5++GFt3rxZ5557rpKSknTNNdc0Ov/27dvtlNOmFRQUGJ2/qKjI6PytDf3wRU+s6IcV/fBFT07MVuCorKy0hA1J3sfV1dWW8V27dqmqqkqJiYmaOnWq3nvvPf3+97/Xa6+9JpfL1eD8ffr0UWRkpJ2STsrj8Uh5+5t1zjNhwIABRub1eDwqKiqSy+VSSEiIkddoTeiHL3piRT+s6IevQO5JRUVFkxcLbAWO8PBwn2BR/zgiIsIyPm3aNE2YMMF7kehFF12kL7/8Un//+98bDRwhISEB98NqjOk+0Gsr+uGLnljRDyv64SsQe2LneG19LDY2NlaHDh1STU2Nd6ysrEwRERFq3769deLgYJ9PpPTs2VMHDhyw85IAAKANsBU44uLiFBoaarm+ID8/Xy6Xy3LBqCTNmDFDqamplrGSkhL17Nnz1KsFAACtkq3A4XA4NGbMGKWnp6uwsFAbNmxQTk6OJk6cKOn4akdVVZUkaejQoXrrrbe0Zs0a7dmzR9nZ2crPz9f48eOb/ygAAECLZvtOo6mpqXI6nZo0aZIyMjKUlJSk4cOHS5ISExO1bt06SdLw4cM1e/ZsLVq0SDfeeKPef/99LVu2TOeff37zHgEAAGjxbN9p1OFwKDMzU5mZmT7PlZaWWh7feuutuvXWW0+9OgAA0Cbw5W0AAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAOAIHAAAwjsABAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAOAIHAAAwjsABAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAOAIHAAAwjsABAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDjbgcPtdistLU0DBw5UYmKicnJyTrrP3r17FR8fr88+++yUigQAAK1bqN0dsrKyVFxcrOXLl2vfvn1KSUlR165dNXLkyEb3SU9PV0VFxWkVCgAAWi9bgaOiokJ5eXlaunSpnE6nnE6nduzYodzc3EYDx5tvvqljx441S7GBpMeMtWZfIG99s0+5e94NzT4nAKBtsHVKpaSkRDU1NYqPj/eOJSQkaOvWraqtrfXZ/tChQ5o/f77mzJlz+pUCAIBWy9YKR1lZmaKjoxUWFuYdi4mJkdvtVnl5uTp16mTZft68ebr55pt14YUXNml+j8cjj8djp6QmzYkzozX2ur7m1li7KfTEin5Y0Q9fgdwTO8dsK3BUVlZawoYk7+Pq6mrL+D//+U/l5+fr7bffbvL827dvt1MOWpiCggJ/l3DKioqK/F1Ci0NPrOiHFf3wRU9OzFbgCA8P9wkW9Y8jIiK8Y1VVVZo1a5Zmz55tGT+ZPn36KDIy0k5JJ+XxeKS8/c06Jxo2YMAAf5dgm8fjUVFRkVwul0JCQvxdTotAT6zohxX98BXIPamoqGjyYoGtwBEbG6tDhw6ppqZGoaHHdy0rK1NERITat2/v3a6wsFDffvutkpOTLfvffffdGjNmTKPXdISEhATcD6stac0/O373fNETK/phRT98BWJP7ByvrcARFxen0NBQFRQUaODAgZKk/Px8uVwuBQf/fP1pv3799O6771r2HT58uB577DFdddVVdl4SAAC0AbYCh8Ph0JgxY5Senq4nnnhC33//vXJycjR37lxJx1c7zj77bEVERKh79+4++8fGxqpz587NUzkAAGg1bN9pNDU1VU6nU5MmTVJGRoaSkpI0fPhwSVJiYqLWrVvX7EUCAIDWzfadRh0OhzIzM5WZmenzXGlpaaP7neg5AADQtvHlbQAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAOAIHAAAwjsABAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAOAIHAAAwjsABAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAOAIHAAAwjsABAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAOAIHAAAwznbgcLvdSktL08CBA5WYmKicnJxGt33zzTc1YsQI9evXT+PGjVNhYeFpFQsAAFon24EjKytLxcXFWr58uWbPnq3s7GytX7/eZ7stW7Zo5syZmjZtmtauXav4+HjdfffdOnbsWLMUDgAAWg9bgaOiokJ5eXmaOXOmnE6nhg0bpilTpig3N9dn27KyMk2bNk2jR4/WBRdcoOnTp6u8vFw7d+5stuIBAEDrEGpn45KSEtXU1Cg+Pt47lpCQoMWLF6u2tlbBwT/nl+uvv97731VVVfrLX/6izp07q1evXs1QNgAAaE1sBY6ysjJFR0crLCzMOxYTEyO3263y8nJ16tTJZ59NmzZp8uTJqqur05NPPql27do1Or/H45HH47FT0kk193xoXGvsdX3NrbF2U+iJFf2woh++Arkndo7ZVuCorKy0hA1J3sfV1dUN7nPhhRdq9erV+uCDDzRjxgydf/75GjBgQIPbbt++3U45aGEKCgr8XcIpKyoq8ncJLQ49saIfVvTDFz05MVuBIzw83CdY1D+OiIhocJ+YmBjFxMQoLi5OW7du1auvvtpo4OjTp48iIyPtlHRSHo9HytvfrHOiYY39XFsyj8ejoqIiuVwuhYSE+LucFoGeWNEPK/rhK5B7UlFR0eTFAluBIzY2VocOHVJNTY1CQ4/vWlZWpoiICLVv396ybWFhoUJCQuR0Or1jvXr1OuFFoyEhIQH3w2pLWvPPjt89X/TEin5Y0Q9fgdgTO8dr61MqcXFxCg0NtSyd5+fny+VyWS4YlaSVK1dqwYIFlrEvv/xSPXv2tPOSAACgDbAVOBwOh8aMGaP09HQVFhZqw4YNysnJ0cSJEyUdX+2oqqqSJN1222369NNPtXz5cu3evVvPPvusCgsLdccddzT7QQAAgJbN9o2/UlNT5XQ6NWnSJGVkZCgpKUnDhw+XJCUmJmrdunWSJKfTqezsbK1cuVI33XSTNm7cqBdffFGxsbHNewQAAKDFs3UNh3R8lSMzM1OZmZk+z5WWlloeDxkyREOGDDn16gAAQJvAl7cBAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAOAIHAAAwjsABAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAOAIHAAAwjsABAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAuFB/F4C2o8eMtf4uwbadj4/0dwkAEBBY4QAAAMYROAAAgHEEDgAAYByBAwAAGEfgAAAAxhE4AACAcQQOAABgnO3A4Xa7lZaWpoEDByoxMVE5OTmNbvvhhx9q9OjRio+P16hRo/Tf//3fp1UsAABonWwHjqysLBUXF2v58uWaPXu2srOztX79ep/tSkpKdN9992ns2LFas2aNxo0bpz/84Q8qKSlplsIBAEDrYetOoxUVFcrLy9PSpUvldDrldDq1Y8cO5ebmauRI6x0b3377bf3mN7/RxIkTJUndu3fX+++/r3feeUcXXXRR8x0BAABo8WwFjpKSEtXU1Cg+Pt47lpCQoMWLF6u2tlbBwT8vmNx888366aeffOY4cuTIaZQLAABaI1uBo6ysTNHR0QoLC/OOxcTEyO12q7y8XJ06dfKO9+rVy7Lvjh07tGnTJo0bN67R+T0ejzwej52STqq550PbUv/7we/Jz+iJFf2woh++Arkndo7ZVuCorKy0hA1J3sfV1dWN7vfDDz8oKSlJl156qa699tpGt9u+fbudcoDTVlRUZPlf/IyeWNEPK/rhi56cmK3AER4e7hMs6h9HREQ0uM/Bgwd15513qq6uTs8++6zltMuv9enTR5GRkXZKOimPxyPl7W/WOdF2uFwuFRUVyeVyKSQkxN/ltAgej4ee/AL9sKIfvgK5JxUVFU1eLLAVOGJjY3Xo0CHV1NQoNPT4rmVlZYqIiFD79u19tj9w4ID3otEVK1ZYTrk0JCQkJOB+WPCv+t83fvd80RMr+mFFP3wFYk/sHK+tj8XGxcUpNDRUBQUF3rH8/Hy5XC6flYuKigpNmTJFwcHBevnllxUbG2vnpQAAQBtiK3A4HA6NGTNG6enpKiws1IYNG5STk+NdxSgrK1NVVZUk6YUXXtA333yjzMxM73NlZWV8SgUAgABk65SKJKWmpio9PV2TJk1SVFSUkpKSNHz4cElSYmKi5s6dq1tuuUX/+Mc/VFVVpVtvvdWy/80336x58+Y1T/UAAKBVsB04HA6HMjMzvSsXv1RaWur974buPgoAAAITX94GAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAOAIHAAAwjsABAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAOAIHAAAwjsABAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAOAIHAAAwjsABAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwL9XcBgD/1mrn++H/krfdvITbsnneDv0sAANtsr3C43W6lpaVp4MCBSkxMVE5Ozkn32bJli6699tpTKhAAALR+tlc4srKyVFxcrOXLl2vfvn1KSUlR165dNXLkyAa3Ly0t1R/+8AeFh4efdrEAAKB1srXCUVFRoby8PM2cOVNOp1PDhg3TlClTlJub2+D2r776qsaNG6fOnTs3S7EAAKB1srXCUVJSopqaGsXHx3vHEhIStHjxYtXW1io42JpfPvroI2VmZuro0aPKzs4+6fwej0cej8dOSU2aE2hLTP9O18/P353j6IcV/fAVyD2xc8y2AkdZWZmio6MVFhbmHYuJiZHb7VZ5ebk6depk2f7555+XJK1evbpJ82/fvt1OOUBAKigoOCOvU1RUdEZep7WgH1b0wxc9OTFbgaOystISNiR5H1dXV592MX369FFkZORpz/NLHo9HytvfrHMC/jRgwACj83s8HhUVFcnlcikkJMToa7UG9MOKfvgK5J5UVFQ0ebHAVuAIDw/3CRb1jyMiIuxM1aCQkJCA+2EBdp2pvyP8fbSiH1b0w1cg9sTO8dq6aDQ2NlaHDh1STU2Nd6ysrEwRERFq3769nakAAEAAsRU44uLiFBoaajmHnJ+fL5fL5XPBKAAAQD1bKcHhcGjMmDFKT09XYWGhNmzYoJycHE2cOFHS8dWOqqoqI4UCAIDWy/ayRGpqqpxOpyZNmqSMjAwlJSVp+PDhkqTExEStW7eu2YsEAACtm+07jTocDmVmZiozM9PnudLS0gb3ueWWW3TLLbfYrw4AALQJXHgBAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMI7AAQAAjCNwAAAA4wgcAADAOAIHAAAwjsABAACMI3AAAADjCBwAAMA4AgcAADCOwAEAAIwjcAAAAOMIHAAAwDgCBwAAMC7U3wUAsKfHjLVn5oXy1jfbVLvn3dBscwFonVjhAAAAxhE4AACAcQQOAABgHIEDAAAYR+AAAADGETgAAIBxBA4AAGAc9+EAYNwZu3dIM+LeIUDzYoUDAAAYR+AAAADGETgAAIBxBA4AAGAcgQMAABhH4AAAAMYROAAAgHEEDgAAYByBAwAAGMedRgGgAT53R81b759CbODuqGjJWOEAAADGETgAAIBxnFIBgDbijH5JXjOeYmqNp4Ia7HULP+3m7z7bXuFwu91KS0vTwIEDlZiYqJycnEa3/fe//61bb71V/fv319ixY1VcXHxaxQIAgNbJ9gpHVlaWiouLtXz5cu3bt08pKSnq2rWrRo4cadmuoqJCU6dO1ahRozRv3jy98soruueee/Tee+8pMjKy2Q4AANC6ndGVGfiNrRWOiooK5eXlaebMmXI6nRo2bJimTJmi3Nxcn23XrVun8PBwPfzww+rVq5dmzpypdu3aaf36lr3kBAAAmp+tFY6SkhLV1NQoPj7eO5aQkKDFixertrZWwcE/55etW7cqISFBQUFBkqSgoCBdeumlKigo0C233GKZt7a2VpJ07NgxeTyeUz6YhtTW1uo/OnKpCgAgsB05cqTZ56yqqpL08/v4idh6Jy4rK1N0dLTCwsK8YzExMXK73SovL1enTp0s2/bu3duyf+fOnbVjxw6fed1utyTpm2++sVNOkz05LMbIvAAAtBbbt283Nrfb7VZUVNQJt7EVOCorKy1hQ5L3cXV1dZO2/fV2ktShQwf16NFD4eHhllUSAADQctXW1srtdqtDhw4n3dZW4AgPD/cJDPWPIyIimrTtr7eTpNDQUHXu3NlOKQAAoAU42cpGPVvLCbGxsTp06JBqamq8Y2VlZYqIiFD79u19tj148KBl7ODBgzrnnHPsvCQAAGgDbAWOuLg4hYaGqqCgwDuWn58vl8vlcyqkf//++uKLL1RXVydJqqur07/+9S/179//9KsGAACtiq3A4XA4NGbMGKWnp6uwsFAbNmxQTk6OJk6cKOn4akf9FasjR47Ujz/+qMcff1xfffWVHn/8cVVWVur6669v/qMAAAAtmu0rNFNTU+V0OjVp0iRlZGQoKSlJw4cPlyQlJiZq3bp1ko6f03nhhReUn5+vW265RVu3btWSJUvO2E2/7NwRNdBUV1frxhtv1GeffebvUvzqwIEDSk5O1uWXX66rr75ac+fO9X5iKhDt2bNHd911l+Lj4zV48GAtW7bM3yW1GFOnTtWMGTP8XYbfvffee+rbt6/lT3Jysr/L8pvq6mplZGTosssu05VXXqkFCxZ4V/Xhy/YNKhwOhzIzM5WZmenzXGlpqeVxv3799Prrr596daehqXdEDTRut1sPPvhggx9PDiR1dXVKTk5W+/btlZubq8OHDystLU3BwcFKSUnxd3lnXG1traZOnSqXy6XXX39de/bs0QMPPKDY2FiNGjXK3+X51dq1a7Vx40bdfPPN/i7F77766isNGTJEjz76qHcsPDzcjxX512OPPabPPvtML774oo4dO6b7779fXbt21bhx4/xdWovUJu+IVX9H1KVLl8rpdMrpdGrHjh3Kzc0N6MDx1Vdf6cEHHySBS9q1a5cKCgr0ySefKCbm+H1akpOTlZmZGZCB4+DBg4qLi1N6erqioqLUo0cPXXHFFcrPzw/owFFeXq6srCy5XC5/l9Ii7Ny5U3369FGXLl38XYrflZeXa9WqVXrppZfUr18/SdLkyZO1detWAkcj2uRNLxq7I+rWrVubdDe0tmrz5s0aNGiQXnvtNX+X4nddunTRsmXLvGGj3tGjR/1UkX+dc845evrppxUVFaW6ujrl5+fr888/1+WXX+7v0vwqMzNTo0eP9rmJYaDauXOnevTo4e8yWoT8/HxFRUVZ/o5MnTpVc+fO9WNVLVubDBwnuyNqoLr99tuVlpYmh8Ph71L8rn379rr66qu9j2tra/Xyyy/rN7/5jR+rahmGDh2q22+/XfHx8RoxYoS/y/GbTZs2acuWLZo2bZq/S2kR6urq9PXXX+t//ud/NGLECF133XV68sknG7yZYyD49ttv1a1bN61Zs0YjR47Utddeq+eeey6g/1F7Mm0ycNi5IyogSfPnz9e///1v3X///f4uxe+effZZLV68WNu2bQvYf6253W7Nnj1bs2bNavBmhYFo37593v9vffrpp5WSkqK33npLWVlZ/i7NLyoqKrRnzx69+uqrmjt3rlJSUvTXv/5Vf/nLX/xdWovVJq/hsHNHVGD+/Plavny5nnrqKfXp08ff5fhd/fUKbrdbDz30kB5++GGfAN/WZWdn65JLLrGsggW6bt266bPPPlOHDh0UFBSkuLg41dbW6o9//KNSU1MVEhLi7xLPqNDQUB09elR//vOf1a1bN0nHQ9krr7yiyZMn+7m6lqlNBo5f3hE1NPT4ITZ2R1QEtkcffVSvvPKK5s+fH9CnDw4ePKiCggJdd9113rHevXvrp59+0tGjRy1fzBgI1q5dq4MHD3qvA6v/B8s//vEPffHFF/4sza86duxoedyrVy+53W4dPnw44H5HunTpovDwcG/YkKT/+I//0HfffefHqlq2NnlKxc4dURG4srOz9eqrr2rBggW64YYb/F2OX+3du1f33XefDhw44B0rLi5Wp06dAu6NRJL++te/6q233tKaNWu0Zs0aDR06VEOHDtWaNWv8XZrffPzxxxo0aJAqKyu9Y9u2bVPHjh0D8nekf//+crvd+vrrr71ju3btsgQQWLXJd9+T3REV2Llzp55//nndfffdSkhIUFlZmfdPIHK5XHI6nUpLS9NXX32ljRs3av78+br33nv9XZpfdOvWTd27d/f+adeundq1a6fu3bv7uzS/iY+PV3h4uB555BHt2rVLGzduVFZWlqZMmeLv0vyiZ8+eGjx4sFJTU1VSUqKPP/5YS5Ys0e9+9zt/l9ZiBdW10ZsyVFZWKj09Xe+++66ioqJ011136Y477vB3WS1G3759tWLFCg0aNMjfpfjFkiVL9Oc//7nB5359A7tAceDAAT366KPatGmTHA6Hxo8fr3vuuUdBQUH+Ls3v6u8yOm/ePD9X4l87duzQE088oYKCArVr107jxo3T9OnTA/Z35MiRI3r00Uf13nvvyeFw6Pbbbw/ofpxMmw0cAACg5WiTp1QAAEDLQuAAAADGETgAAIBxBA4AAGAcgQMAABhH4AAAAMYROAAAgHEEDgAAYByBAwAAGEfgAAAAxhE4AACAcf8f+D8RMNgVttcAAAAASUVORK5CYII=",
"text/plain": [
"<Figure size 640x480 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"plt.hist(idata.prior[\"mu\"].values[0], density=True)\n",
"plt.title(\"Prior, Exponential(lam=1)\")\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Sample from the posterior..."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"Auto-assigning NUTS sampler...\n",
"Initializing NUTS using jitter+adapt_diag...\n",
"Multiprocess sampling (4 chains in 4 jobs)\n",
"NUTS: [mu, sigma]\n"
]
},
{
"data": {
"text/html": [
"\n",
"<style>\n",
" /* Turns off some styling */\n",
" progress {\n",
" /* gets rid of default border in Firefox and Opera. */\n",
" border: none;\n",
" /* Needs to be in here for Safari polyfill so background images work as expected. */\n",
" background-size: auto;\n",
" }\n",
" progress:not([value]), progress:not([value])::-webkit-progress-bar {\n",
" background: repeating-linear-gradient(45deg, #7e7e7e, #7e7e7e 10px, #5c5c5c 10px, #5c5c5c 20px);\n",
" }\n",
" .progress-bar-interrupted, .progress-bar-interrupted::-webkit-progress-bar {\n",
" background: #F44336;\n",
" }\n",
"</style>\n"
],
"text/plain": [
"<IPython.core.display.HTML object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"\n",
" <div>\n",
" <progress value='8000' class='' max='8000' style='width:300px; height:20px; vertical-align: middle;'></progress>\n",
" 100.00% [8000/8000 00:00&lt;00:00 Sampling 4 chains, 0 divergences]\n",
" </div>\n",
" "
],
"text/plain": [
"<IPython.core.display.HTML object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"Sampling 4 chains for 1_000 tune and 1_000 draw iterations (4_000 + 4_000 draws total) took 1 seconds.\n"
]
}
],
"source": [
"with model:\n",
" idata.extend(pm.sample())"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Now compare the prior distribution (Exponential(1)) to the posterior, which looks like\n",
"Normal(5, 1):"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Text(0.5, 1.0, 'Posterior')"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAhwAAAGsCAYAAACW3H6UAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABA3ElEQVR4nO3deXRU9f3/8VcWshFBCAFlEQwQljgESASRIIoVQqFWwdiKIHULDUq+Iv5YLYuURoiCSyBsYmmlEDY5RVIqCl+0fpFoICFAIZsiHJQOpxDEJBMzye8PDoPjkJAJc2eS4fk4Z07P3PnM577n3dG8vHPv5/pUV1dXCwAAwEC+ni4AAAB4PwIHAAAwHIEDAAAYjsABAAAMR+AAAACGI3AAAADDETgAAIDhCBwAAMBwBA4AAGA4AgfQwI0bN07dunWze9xxxx269957NW/ePJWUlNT6/lOnTqlbt27aunWroXWuWbNGL730kt0+a3usX7/ekDoOHTqkYcOGqaKiwuG13bt3q1u3bg7bx44dq8zMTLtt+/bt069//Wv9+OOPhtQJ3Gj8PV0AgGvr2bOn5syZY3v+448/6siRI1q8eLH+/e9/a/369fLx8bnqe1u3bq2MjAzddttthtVXVFSkFStW6O9//7vd9qSkJN17771XfU+HDh1cXofFYtG0adP0//7f/1NAQIDda/v379eUKVOu+r6ZM2fq6aefVv/+/RUWFiZJGjBggNq1a6dly5bpf/7nf1xeK3CjIXAAjUBoaKh69+5tt+3OO+/UDz/8oLfeeku5ubkOr18WEBBQ42uukpqaqpEjR6pNmzZ222+77TbD9/1Tf/vb3+Tv769f/OIXtm0XL17UqlWrtGrVKt10000qLS11eF/Pnj3Vq1cvpaen6+WXX7ZtT0pK0pgxY/TYY4+pdevWbvkMgLfiJxWgEbvjjjskSadPn5Z06eeXl156ScnJyerdu7eefPLJq/6k8vXXXys5OVkDBw5U7969NW7cOGVnZ9tev/yed999V/Hx8YqOjtaWLVuuWkN+fr7+93//VyNHjqzXZ3j++edlMplUXFxs2/b222+rR48eysrKkiQNGTJES5Ys0Z/+9Cfdeeed6t+/v6ZOnarz58/b3lNRUaF3333XoY7Nmzdr48aNmj17tsaOHVtjHb/61a+0efNm/fe//7VtM5lMatu2rd599916fTYAVxA4gEbsq6++kmT/88Q//vEPNW3aVOnp6XrmmWcc3lNYWKhRo0bp1KlTevnll/Xaa6/Jx8dH48ePt/2Bv+ztt9/Ws88+q0WLFmngwIFXrWH79u0KDw+/6pGMqqoqVVZWOjysVqttzNy5cxUSEmL7yejw4cNavny5nnrqKfXr18827m9/+5sOHDiglJQUTZkyRXv37tWECRN0+YbX+/fv15kzZzR06FC7GoYMGaLdu3frt7/9bW2t1JAhQ2S1WrVr1y677fHx8frggw9qfS+Aa+MnFaARqK6uVmVlpe15SUmJsrKylJ6erj59+tiOdEhSkyZNNG/ePNs5DKdOnbKbKy0tTQEBAfrLX/6i0NBQSdK9996rkSNHatGiRdq8ebNt7PDhwzV69Ohaa/v8889lMpmueg7JrFmzNGvWLIftISEhOnjwoCSpVatWmjNnjiZPnqxNmzZp7dq1ioyMdDhvwtfXV++++65uuukmSVLLli313HPP6dNPP9U999yjzz//XM2aNdPtt99u9766nrsSEhKizp07a9++ffrNb35j224ymbR8+XIVFRWpc+fOdZoLgCMCB9AIfPHFF4qKirLb5uvrq7vvvluvvPKK3R/7iIgIhxMmfyorK0v33XefLWxIkr+/v0aMGKGlS5fqhx9+sG3v0aPHNWs7efKk+vTpc9XXnn/++aueNOrn52f3/Je//KV27typ2bNnKyAgQFu3bnX4DEOGDLGFjcvP/f399cUXX+iee+7RyZMn1a5du2vWW5t27do5BLT27dtLuhTcCBxA/RE4gEYgKipK8+bNkyT5+PgoMDBQt956q11ouKxp06a1zlVSUqJWrVo5bG/VqpWqq6t18eJF27aQkJBr1nbx4kUFBwdf9bV27drJZDJdcw5Jevjhh/XPf/5TnTp1cjhKIcnhhFRfX1+1aNHCdllwbXXUVXBwsL7//nuHbZIctgNwDoEDaASaNm1a5z/c19K8eXOdPXvWYbvZbJYktWjRQv/5z3/qPN/NN9983X+My8rKlJKSosjISOXn52vNmjUO55+cO3fO7rnVatW5c+fUsmXLetV9NRcuXFCLFi3stl0OND/fDsA5nDQK3GDuvPNO7dmzx+5IhtVq1Y4dO2QymWr9OeZq2rVrp2+//fa6anr99df13Xff6e2339bYsWP11ltvqaioyG7MJ598YreY18cff6zKykoNGDBAktS2bVt99913tpNI6+O7775z+FnmzJkztvkB1B+BA7jBPP/887JYLHriiSe0c+dOffzxx3rmmWd08uRJvfjii07PN3DgQB08ePCqf+i/+eYb5eTkXPVx+QqbrKwsvffee3ruuefUqVMnvfDCC2rZsqWmT59udzXLt99+q6SkJO3du1cbNmzQyy+/rEGDBql///62Or7//nvl5+fXqy/ff/+9CgoKNGjQILvt2dnZat++/VV/5gFQd/ykAtxgunbtqr/97W9avHixZsyYIR8fH/Xq1Ut/+ctfFBsb6/R8Q4cO1dKlS3Xo0CFFR0fbvZaenq709PSrvu/+++/Xa6+9phkzZigyMlJPP/20pEs/H82ePVtJSUlavXq1JkyYIEkaMWKEmjVrphdeeEEhISF6+OGHNXnyZNt8sbGxCgsL0969e6+6fPm1fPrpp2rSpInDSa6ffvqp4uPjnZ4PgD2f6us5/ggAkn7/+9+rRYsWSklJMWT+IUOGqF+/fnr11VdrHbdmzRqtX79eH374YY1Lvddk/PjxioyMtLuM98svv9RTTz2ljz76iJVGgevETyoArtvkyZP14Ycf2lY89ZQxY8aoqqpKO3fudOp9eXl5OnbsmBITE+22r169WuPHjydsAC5A4ABw3bp166YJEybotdde82gdQUFBSk1N1ZIlS656t9iapKSk6A9/+IPCw8Nt2/bt26fTp09r0qRJRpQK3HD4SQUAABiOIxwAAMBwBA4AAGA4AgcAADBcg1iHo7KyUiUlJQoMDJSvLxkIAIDGoKqqShaLRc2bN5e/f+2RokEEjpKSEn399deeLgMAANRDp06dFBYWVuuYBhE4AgMDJV0q+Hrv9vhzVqtV+fn5ioyMdLglNlyHPrsHfXYP+uw+9No9jOpzWVmZvv76a9vf8do0iMBx+WeU4ODgOt0O2xmX78UQEhLCl9lA9Nk96LN70Gf3odfuYXSf63I6BCdMAAAAwxE4AACA4QgcAADAcAQOAABgOAIHAAAwHIEDAAAYjsABAAAM53TgOHPmjJKTk9WvXz8NGjRIKSkpslgsVx2blJSkbt262T327Nlz3UUDAIDGxamFv6qrq5WcnKxmzZpp3bp1Kikp0cyZM+Xr66tp06Y5jC8qKlJqaqoGDBhg29a8efPrrxoAADQqTgWO4uJi5eTk6LPPPlOrVq0kScnJyVq4cKFD4KioqNCpU6dkMpkUHh7uuooBAECj41TgCA8P1+rVq21h47KLFy86jC0uLpaPj486dOhQ5/mtVqtt+VVXuTyfq+eFPfrsHvTZPeiz+9Br9zCqz87M51NdXV1d3x1VVVVpzJgxatGihdLT0+1ey8zM1Lx58zRw4EBlZWXplltu0aRJkzR48GCHeUpLS/Xvf/+7vmUAAAAP6tGjxzXvhXZdN29LTU3V0aNHtXnzZofXiouLVV5erri4OCUmJmrXrl1KSkpSRkaGTCbTVeeLjIw05OZteXl5MplM3BjIQPTZPeize9Bn96HX7mFUn0tLS5Wfn1+nsfUOHKmpqVq7dq2WLFmiyMhIh9cnTpyocePG2U4S7d69u44cOaKNGzfWGDj8/PwM+8IZOTeuoM/uQZ/dgz67D712D1f32Zm56rUOx/z58/Xuu+8qNTVVw4YNu/rEvr4OV6RERETozJkz9dklAABoxJwOHGlpadqwYYMWL16sESNG1Dhu+vTpmjFjht22Y8eOKSIiwvkqAQBAo+ZU4CgqKtKyZcv07LPPKiYmRmaz2faQJLPZrPLycknSkCFDtH37dm3btk0nTpxQWlqasrOzNXbsWNd/CgAA0KA5dQ7Hxx9/LKvVqvT0dIerUo4fP664uDilpKRo1KhRGjp0qObMmaP09HSdPn1aXbt21erVq9W+fXuXfgAAANDwORU4EhMTlZiYWOPrx48ft3uekJCghISE+lUGAAC8BjdvAwAAhiNwAAAAwxE4AACA4QgcAADAcAQOAABgOAIHAAAwHIEDAAAYjsABAAAMR+AAAACGI3AAAADDETgAAIDhCBwAAMBwBA4AAGA4AgcAADAcgQMAABiOwAEAAAxH4AAAAIYjcAAAAMM5FTjOnDmj5ORk9evXT4MGDVJKSoosFstVxx49elQJCQmKjo7W6NGjdfjwYZcUDAAAGp86B47q6molJyerrKxM69at05IlS7Rnzx698cYbDmNLS0uVmJio2NhYbd26VX369NGECRNUWlrqytoBAEAjUefAUVxcrJycHKWkpKhr166KjY1VcnKyPvjgA4exmZmZCgwM1NSpU9W5c2fNmjVLTZs21c6dO11aPAAAaBzqHDjCw8O1evVqtWrVym77xYsXHcbm5uYqJiZGPj4+kiQfHx/17dtXOTk511ctAABolPzrOrBZs2YaNGiQ7XlVVZXee+893XXXXQ5jzWazunTpYrctLCxMBQUFte7DarXKarXWtaQ6uTyfq+eFPfrsHvTZPeiz+9Br9zCqz87MV+fA8XOpqak6evSoNm/e7PBaWVmZAgIC7LYFBASooqKi1jnz8/PrW8415eXlGTY3rqDP7kGf3YM+uw+9dg9P9rlegSM1NVVr167VkiVLFBkZ6fB6YGCgQ7ioqKhQUFBQrfNGRkYqJCSkPiXVyGq1Ki8vTyaTSX5+fi6dG1fQZ/egz+5Bn92HXruHUX0uLS2t88ECpwPH/PnztX79eqWmpmrYsGFXHdOmTRudPXvWbtvZs2fVunXrWuf28/Mz7Atn5Ny4gj67B312D/rsPvTaPVzdZ2fmcmodjrS0NG3YsEGLFy/WiBEjahwXHR2tgwcPqrq6WtKlS2oPHDig6OhoZ3YHAAC8RJ0DR1FRkZYtW6Znn31WMTExMpvNtod06UTR8vJySVJ8fLwuXLigBQsWqLCwUAsWLFBZWZmGDx9uzKcAAAANWp0Dx8cffyyr1ar09HTFxcXZPSQpLi5OmZmZkqTQ0FCtWLFC2dnZGjVqlHJzc7Vy5UqXn58BAAAahzqfw5GYmKjExMQaXz9+/Ljd8169eun999+vf2UAAMBrcPM2AABgOAIHAAAwHIEDAAAYjsABAAAMR+AAAACGI3AAAADDETgAAIDhCBwAAMBwBA4AAGA4AgcAADAcgQMAABiOwAEAAAxH4AAAAIYjcAAAAMMROAAAgOEIHAAAwHAEDgAAYDgCBwAAMByBAwAAGK7egaOiokIjR47U/v37axyTlJSkbt262T327NlT310CAIBGyr8+b7JYLJoyZYoKCgpqHVdUVKTU1FQNGDDAtq158+b12SUAAGjEnA4chYWFmjJliqqrq2sdV1FRoVOnTslkMik8PLzeBQIAgMbP6Z9UsrKy1L9/f2VkZNQ6rri4WD4+PurQoUO9iwMAAN7B6SMcY8aMqdO44uJihYaGaurUqcrKytItt9yiSZMmafDgwTW+x2q1ymq1OltSrS7P5+p5YY8+uwd9dg/67D702j2M6rMz89XrHI66KC4uVnl5ueLi4pSYmKhdu3YpKSlJGRkZMplMV31Pfn6+UeUoLy/PsLlxBX12D/rsHvTZfei1e3iyz4YFjokTJ2rcuHG2k0S7d++uI0eOaOPGjTUGjsjISIWEhLi0DqvVqry8PJlMJvn5+bl0blxBn92DPrsHfXYfeu0eRvW5tLS0zgcLDAscvr6+DlekREREqLCwsMb3+Pn5GfaFM3JuXEGf3YM+uwd9dh967R6u7rMzcxm28Nf06dM1Y8YMu23Hjh1TRESEUbsEAAANlEsDh9lsVnl5uSRpyJAh2r59u7Zt26YTJ04oLS1N2dnZGjt2rCt3CQAAGgGXBo64uDhlZmZKkoYOHao5c+YoPT1dI0eO1O7du7V69Wq1b9/elbsEAACNwHWdw3H8+PFanyckJCghIeF6dgEAALwAN28DAACGI3AAAADDETgAAIDhCBwAAMBwBA4AAGA4AgcAADAcgQMAABiOwAEAAAxH4AAAAIYjcAAAAMMROAAAgOEIHAAAwHAEDgAAYDgCBwAAMByBAwAAGI7AAQAADEfgAAAAhiNwAAAAwxE4AACA4eodOCoqKjRy5Ejt37+/xjFHjx5VQkKCoqOjNXr0aB0+fLi+uwMAAI1YvQKHxWLRiy++qIKCghrHlJaWKjExUbGxsdq6dav69OmjCRMmqLS0tN7FAgCAxsnpwFFYWKhHH31U33zzTa3jMjMzFRgYqKlTp6pz586aNWuWmjZtqp07d9a7WAAA0Dj5O/uGrKws9e/fX5MnT1bv3r1rHJebm6uYmBj5+PhIknx8fNS3b1/l5ORo1KhRV32P1WqV1Wp1tqRaXZ7P1fPCHn12D/rsHvTZfei1exjVZ2fmczpwjBkzpk7jzGazunTpYrctLCys1p9h8vPznS2nzvLy8gybG1fQZ/egz+5Bn92HXruHJ/vsdOCoq7KyMgUEBNhtCwgIUEVFRY3viYyMVEhIiEvrsFqtysvLk8lkkp+fn0vnxhX02T3os3vQZ/eh1+5hVJ9LS0vrfLDAsMARGBjoEC4qKioUFBRU43v8/PwM+8IZOTeuoM/uQZ/dgz67D712D1f32Zm5DFuHo02bNjp79qzdtrNnz6p169ZG7RIAADRQhgWO6OhoHTx4UNXV1ZKk6upqHThwQNHR0UbtEgAANFAuDRxms1nl5eWSpPj4eF24cEELFixQYWGhFixYoLKyMg0fPtyVuwQAAI2ASwNHXFycMjMzJUmhoaFasWKFsrOzNWrUKOXm5mrlypUuPykUAAA0fNd10ujx48drfd6rVy+9//7717MLAADgBbh5GwAAMByBAwAAGI7AAQAADEfgAAAAhiNwAAAAwxE4AACA4QgcAADAcAQOAABgOAIHAAAwHIEDAAAY7rqWNm8sRm/6Ttq0s9YxX786wk3VAABw47khAkdddJq+45pjCCUAANQPP6kAAADDETgAAIDhCBwAAMBwBA4AAGA4AgcAADAcV6k4gStZAACoH6ePcFgsFs2cOVOxsbGKi4vTmjVrahyblJSkbt262T327NlzXQUDAIDGx+kjHIsWLdLhw4e1du1anT59WtOmTVPbtm0VHx/vMLaoqEipqakaMGCAbVvz5s2vr2IAANDoOBU4SktLtWnTJq1atUpRUVGKiopSQUGB1q1b5xA4KioqdOrUKZlMJoWHh7u0aAAA0Lg49ZPKsWPHVFlZqT59+ti2xcTEKDc3V1VVVXZji4uL5ePjow4dOrimUgAA0Gg5dYTDbDarRYsWCggIsG1r1aqVLBaLzp8/r5YtW9q2FxcXKzQ0VFOnTlVWVpZuueUWTZo0SYMHD65xfqvVKqvVWo+PUTNXz+eq/XWeVfu9XSSpaIHjz1QN1eXP7e5+32jos3vQZ/eh1+5hVJ+dmc+pwFFWVmYXNiTZnldUVNhtLy4uVnl5ueLi4pSYmKhdu3YpKSlJGRkZMplMV50/Pz/fmXIapJycnAY5l7vk5eV5uoQbAn12D/rsPvTaPTzZZ6cCR2BgoEOwuPw8KCjIbvvEiRM1btw420mi3bt315EjR7Rx48YaA0dkZKRCQkKcKemarFartOk7l85Zm9Eu3Ffv3r1dNpfRrFar8vLyZDKZ5Ofn5+lyvBZ9dg/67D702j2M6nNpaWmdDxY4FTjatGmjc+fOqbKyUv7+l95qNpsVFBSkZs2a2Y319fV1uCIlIiJChYWFNc7v5+fHF+4nGmMv+P/QPeize9Bn96HX7uHqPjszl1Mnjfbo0UP+/v52h/qzs7NlMpnk62s/1fTp0zVjxgy7bceOHVNERIQzuwQAAF7AqcARHByshx56SHPnztWhQ4f00Ucfac2aNXriiSckXTraUV5eLkkaMmSItm/frm3btunEiRNKS0tTdna2xo4d6/pPAQAAGjSnVxqdMWOGoqKiNH78eM2bN0+TJk3S0KFDJUlxcXHKzMyUJA0dOlRz5sxRenq6Ro4cqd27d2v16tVq3769az8BAABo8JxeaTQ4OFgLFy7UwoULHV47fvy43fOEhAQlJCTUvzoAAOAVuFssAAAwHIEDAAAYjsABAAAMR+AAAACGI3AAAADDOX2VCtyn0/Qd1xzz9asj3DYPAAD1xREOAABgOAIHAAAwHIEDAAAYjnM4Grm6nJ/hynk41wMAUB8c4QAAAIYjcAAAAMMROAAAgOEIHAAAwHCcNAqnXPPk0k07ObEUAOCAIxwAAMBwHOGAy7nqUl2OlACA9+AIBwAAMByBAwAAGM7pwGGxWDRz5kzFxsYqLi5Oa9asqXHs0aNHlZCQoOjoaI0ePVqHDx++rmIBAEDj5PQ5HIsWLdLhw4e1du1anT59WtOmTVPbtm0VHx9vN660tFSJiYn61a9+pVdffVXr16/XhAkTtGvXLoWEhLjsA+DG5sol2esyF+eVAED9OBU4SktLtWnTJq1atUpRUVGKiopSQUGB1q1b5xA4MjMzFRgYqKlTp8rHx0ezZs3SJ598op07d2rUqFEu/RDwTq46+dSVc7kqlHDvGgA3GqcCx7Fjx1RZWak+ffrYtsXExGj58uWqqqqSr++VX2hyc3MVExMjHx8fSZKPj4/69u2rnJwcAge8mluC0qadTs3TEIOLO8NbQ/z8wI3GqcBhNpvVokULBQQE2La1atVKFotF58+fV8uWLe3GdunSxe79YWFhKigocJi3qqpKkvTDDz/IarU69QGupaqqSrffzNW/uLHd9+o/3bq/vz939zXH1OWfy++//95t80jSg0v/z37Dro8cxtTlsznMU0912Vdjd/nf/xcvXrT7j1a4llF9Li8vt5u/Nk79JS4rK7MLG5JszysqKuo09ufjpEsnokrSN99840w5dfbaA60MmRfA1eXn519zTF3+uXTnPO6uqS7qWrc3KCws9HQJNwSj+myxWBQaGlrrGKcCR2BgoENguPw8KCioTmN/Pk6Smjdvrk6dOikwMJCECwBAI1FVVSWLxaLmzZtfc6xTgaNNmzY6d+6cKisr5e9/6a1ms1lBQUFq1qyZw9izZ8/abTt79qxat27tWIS/v8LCwpwpBQAANADXOrJxmVOHE3r06CF/f3/l5OTYtmVnZ8tkMjkcmYiOjtbBgwdVXV0tSaqurtaBAwcUHR3tzC4BAIAXcCpwBAcH66GHHtLcuXN16NAhffTRR1qzZo2eeOIJSZeOdlw+gSQ+Pl4XLlzQggULVFhYqAULFqisrEzDhw93/acAAAANmtMnTMyYMUNRUVEaP3685s2bp0mTJmno0KGSpLi4OGVmZkq6dIhlxYoVys7O1qhRo5Sbm6uVK1e6ddEvZ1ZFRf2dOXNGycnJ6tevnwYNGqSUlBTbicBwvcTERE2fPt3TZXitiooKzZs3T3feeafuvvtuLV682HakFq7z7bffasKECerbt6+GDBmiP//5z54uyetUVFRo5MiR2r9/v23byZMn9bvf/U69e/fWL3/5S/3rX/9yWz1OXy8aHByshQsXauHChQ6vHT9+3O55r1699P7779e/uutU11VRUX/V1dVKTk5Ws2bNtG7dOpWUlGjmzJny9fXVtGnTPF2e19mxY4f27t2rhx9+2NOleK0//vGP2r9/v9555x398MMPmjx5stq2bavf/va3ni7Nq7zwwgtq27attm7dqsLCQr300ktq166dHnjgAU+X5hUsFoumTJlitxRFdXW1nnvuOUVGRmrLli366KOP9PzzzyszM1Nt27Y1vCavvSTk8qqos2bNUlRUlB544AE988wzWrdunadL8yrFxcXKyclRSkqKunbtqtjYWCUnJ+uDDz7wdGle5/z581q0aJFMJpOnS/Fa58+f15YtWzR//nz16tVLAwYM0FNPPaXc3FxPl+ZVSkpKlJOTo6SkJHXq1Em/+MUvNGjQIO3bt8/TpXmFwsJCPfroow5LTXz++ec6efKkXnnlFXXu3FkTJkxQ7969tWXLFrfU5bWBo6ZVUXNzc+u0QAnqJjw8XKtXr1arVvbrDly8eNFDFXmvhQsX6te//rXDgnpwnezsbIWGhqpfv362bYmJiUpJSfFgVd4nKChIwcHB2rp1q3788UcVFxfrwIED6tGjh6dL8wpZWVnq37+/MjIy7Lbn5uaqZ8+edqc2xMTE2F0IYiSvDRzXWhUVrtGsWTMNGjTI9ryqqkrvvfee7rrrLg9W5X327dunL7/8UhMnTvR0KV7t5MmTateunbZt26b4+Hjdf//9Wrp0Kf+R4mKBgYGaPXu2MjIyFB0dreHDh+uee+5RQkKCp0vzCmPGjNHMmTMVHBxst91sNjssTREWFqbvvvvOLXV57ZrfzqyKCtdJTU3V0aNHtXnzZk+X4jUsFovmzJmj2bNnX3XhPLhOaWmpTpw4oQ0bNiglJUVms1mzZ89WcHCwnnrqKU+X51WKiop033336cknn1RBQYHmz5+vAQMG6MEHH/R0aV7LmRXAjeC1gcOZVVHhGqmpqVq7dq2WLFmiyMhIT5fjNdLS0nTHHXfYHUmCMfz9/XXx4kW9/vrrateunSTp9OnTWr9+PYHDhfbt26fNmzdr7969CgoKkslk0pkzZ5Senk7gMFBgYKDDEf6aVgA3gtcGDmdWRcX1mz9/vtavX6/U1FQNGzbM0+V4lR07dujs2bO285EuB+d//vOfOnjwoCdL8zrh4eEKDAy0hQ1Juv322/Xtt996sCrvc/jwYXXs2NHuD13Pnj21fPlyD1bl/dq0aeNwL5WaVgA3gtcGjp+uihobGyup5lVRcX3S0tK0YcMGLV68mEuODfDXv/5VlZWVtuevvfaaJOmll17yVEl1Ul1dLR8fH0+X4ZTo6GhZLBZ99dVXuv322yVduhLrpwEE169169Y6ceKEKioqbIf4i4uL1b59ew9X5t2io6O1cuVKlZeX28Jedna2YmJi3LJ/r/3Le61VUeEaRUVFWrZsmZ599lnFxMTIbDbbHnCNdu3aqWPHjrZH06ZN1bRpU3Xs2LHOc4wbN07dunWze9xxxx269957NW/ePJWUlLis3oqKCv3pT3/S9u3bXTLfuHHjNG7cOJfMdS0RERG69957NWPGDB07dkyffvqpVq5cqccee8wt+79RDBkyRE2aNNHLL7+sr776Srt379by5cvd9v/zjapfv3669dZbNWPGDBUUFGjlypU6dOiQHnnkEbfs32uPcEiXVkWdO3euxo8fr9DQULtVUeEaH3/8saxWq9LT05Wenm732s8XgoNn9ezZU3PmzLE9//HHH3XkyBEtXrxY//73v7V+/XqXHJH4z3/+o7Vr17rsUtKf1uwOr732mubPn6/HHntMwcHBevzxx/lD6GI33XST/vznP2vBggV65JFH1LJlSyUlJek3v/mNp0vzan5+flq2bJlmzZqlUaNGqWPHjlq6dKlbFv2SJJ9q1uwFvN7lP5h//etfHV5bunSp3nrrLWVkZKh3797Xva9Tp07p/vvvV0pKikaNGnXd8wHwDl77kwqAurnjjjskXboaQ5IyMzM1atQo9enTRwMHDtTs2bPtfnIpLy/X3Llzdc899+iOO+5QfHy83nnnHUlXwoZ06QjjkCFDbO/78ssvNXbsWEVHR6tfv36aNm2a/vvf/9pe37p1q3r27KlNmzZp4MCB6tevnwoLCx1+UrFYLFq6dKni4+NlMpk0dOhQrVy50m6tjHHjxumll15ScnKyevfurSeffNKAzgFwhlf/pALg2r766itJUocOHbRs2TK99dZbGjNmjCZPnqyTJ0/qzTffVE5OjjZu3KigoCD96U9/0r/+9S9NmzZNrVq10ieffKJFixbp5ptv1q9+9SulpaXp+eefV1JSku0nzC+++EJPPvmk7rrrLr3xxhsqKSnRm2++qSeeeEKbN2+2ncBmtVq1Zs0aLViwQOfOnVPnzp3taq2urtbvf/975eTk6Pnnn1f37t21f/9+vfHGGzp58qTmz59vG/uPf/xDDz74oNLT01m4C2gACBzADaK6utruapeSkhJlZWUpPT1dffr00W233ab09HQ9+uijmj17tm1cZGSkHn/8cW3ZskWPP/64srKyNHDgQI0YMUKS1L9/f4WEhCgsLEwBAQG25alvu+029ezZU5L0+uuv6/bbb9eKFSvk5+cn6dIZ8yNGjLDNe9nvf/973XvvvVf9DJ988on+7//+T4sXL7btf+DAgQoKCrIFmK5du0qSmjRponnz5jksdATAMwgcwA3iiy++UFRUlN02X19f3X333XrllVeUk5Nju531T8XGxqpdu3bKysrS448/rv79+2vDhg367rvvNHjwYA0ePFjPPfdcjfstKytTbm6unn76abvQ06FDB3Xu3FmfffaZXeCo7X4aWVlZ8vf3d7j8+sEHH9Sbb76prKwsW+CIiIggbAANCIEDuEFERUVp3rx5kiQfHx8FBgbq1ltvVWhoqKRL1+NLcrgR3+Vt33//vSRp1qxZuuWWW/T3v/9d8+fP1/z589WnTx/NnTtX3bt3d3jvhQsXVFVVpVWrVmnVqlUOrwcGBto9/+mNpX6upKRELVq0sB0luSw8PFySbDVKUtOmTWucB4D7ETiAG0TTpk1rvbV98+bNJV1aeTAiIsLuNbPZrA4dOki6dO+FpKQkJSUl6fTp09qzZ4+WLVumKVOmaMeOHVfdr4+Pj373u9/Zfgb5qZ/fYKo2zZs317lz52S1Wu1Cx3/+8x9JUosWLeo8FwD34ioVAJIunVMREBCgDz74wG77l19+qdOnT6tv374qLy/XsGHDtGbNGklS27Zt9fjjj2vEiBG2q1x+fvQhNDRUPXv2VHFxsUwmk+3RtWtXvf3229q/f3+da+zXr58qKyu1c+dOu+1///vfJcltKyYCcB5HOABIkm6++WYlJiZq6dKlatKkie677z6dOnVKb775prp06aKHH35YQUFBioqKUlpampo0aaJu3brpq6++0vvvv2+7h85NN90k6dINujp37qzo6Gi9+OKLSkxM1JQpU/Tggw/arkbJzc3VxIkT61zjPffco/79++vll1/WmTNn1L17d2VlZWnVqlV6+OGH1aVLF0N6A+D6ETgA2EyaNEmtWrXSe++9p4yMDN18882Kj4/XCy+8YDu34pVXXtEbb7yhNWvWyGw2KywsTI888oj+53/+R9KlIxpPPvmkMjIytHfvXn322WeKi4vTO++8o7S0NCUnJ6tJkyaKiorSu+++69RiYz4+PlqxYoXeeust/fnPf9Z///tftW/fXi+++CJrbQANHCuNAgAAw3EOBwAAMByBAwAAGI7AAQAADEfgAAAAhiNwAAAAwzWIy2IrKytVUlKiwMBA+fqSgQAAaAyqqqpksVjUvHlz+fvXHikaROAoKSnR119/7ekyAABAPXTq1ElhYWG1jmkQgePyzZs6derk1H0V6sJqtSo/P1+RkZEOSy7Ddeize9Bn96DP7kOv3cOoPpeVlenrr792uAnj1TSIwHH5Z5Tg4OBa7xRZH1arVdKlO1DyZTYOfXYP+uwe9Nl96LV7GN3nupwOwQkTAADAcAQOAABgOAIHAAAwHIEDAAAYjsABAAAM1yCuUgHQOHSavqNO475+dYTBlQBobDjCAQAADOd04Dhz5oySk5PVr18/DRo0SCkpKbJYLFcdm5SUpG7dutk99uzZc91FAwCAxsWpn1Sqq6uVnJysZs2aad26dSopKdHMmTPl6+uradOmOYwvKipSamqqBgwYYNvWvHnz668aAAA0Kk4FjuLiYuXk5Oizzz5Tq1atJEnJyclauHChQ+CoqKjQqVOnZDKZFB4e7rqKAQBAo+PUTyrh4eFavXq1LWxcdvHiRYexxcXF8vHxUYcOHa6vQgAA0Og5dYSjWbNmGjRokO15VVWV3nvvPd11110OY4uLixUaGqqpU6cqKytLt9xyiyZNmqTBgwfXOL/VarWt9+4ql+dz9bywR5/do7H0uaHXdy2Npc/egF67h1F9dma+67osNjU1VUePHtXmzZsdXisuLlZ5ebni4uKUmJioXbt2KSkpSRkZGTKZTFedLz8//3rKqVVeXp5hc+MK+uweDb3POTk5ni7BJRp6n70JvXYPT/a53oEjNTVVa9eu1ZIlSxQZGenw+sSJEzVu3DjbSaLdu3fXkSNHtHHjxhoDR2RkpCF3i83Ly5PJZOJOhAaiz+7h8T5v2lmnYb179za2DoN5vM83EHrtHkb1ubS0tM4HC+oVOObPn6/169crNTVVw4YNu+oYX19fhytSIiIiVFhYWOO8fn5+hn3hjJwbV9Bn92jofW7ItTmjoffZm9Br93B1n52Zy+l1ONLS0rRhwwYtXrxYI0bUvJrg9OnTNWPGDLttx44dU0REhLO7BAAAjZxTgaOoqEjLli3Ts88+q5iYGJnNZttDksxms8rLyyVJQ4YM0fbt27Vt2zadOHFCaWlpys7O1tixY13/KQAAQIPm1E8qH3/8saxWq9LT05Wenm732vHjxxUXF6eUlBSNGjVKQ4cO1Zw5c5Senq7Tp0+ra9euWr16tdq3b+/SDwAAABo+pwJHYmKiEhMTa3z9+PHjds8TEhKUkJBQv8oAAIDX4OZtAADAcAQOAABgOAIHAAAwHIEDAAAYjsABAAAMR+AAAACGI3AAAADDETgAAIDhCBwAAMBwBA4AAGA4AgcAADAcgQMAABiOwAEAAAxH4AAAAIYjcAAAAMMROAAAgOEIHAAAwHAEDgAAYDgCBwAAMJxTgePMmTNKTk5Wv379NGjQIKWkpMhisVx17NGjR5WQkKDo6GiNHj1ahw8fdknBAACg8alz4KiurlZycrLKysq0bt06LVmyRHv27NEbb7zhMLa0tFSJiYmKjY3V1q1b1adPH02YMEGlpaWurB0AADQSdQ4cxcXFysnJUUpKirp27arY2FglJyfrgw8+cBibmZmpwMBATZ06VZ07d9asWbPUtGlT7dy506XFAwCAxqHOgSM8PFyrV69Wq1at7LZfvHjRYWxubq5iYmLk4+MjSfLx8VHfvn2Vk5NzfdUCAIBGyb+uA5s1a6ZBgwbZnldVVem9997TXXfd5TDWbDarS5cudtvCwsJUUFBQ6z6sVqusVmtdS6qTy/O5el7Yo8/u0Vj63NDru5bG0mdvQK/dw6g+OzNfnQPHz6Wmpuro0aPavHmzw2tlZWUKCAiw2xYQEKCKiopa58zPz69vOdeUl5dn2Ny4gj67R0Pvs7cczWzoffYm9No9PNnnegWO1NRUrV27VkuWLFFkZKTD64GBgQ7hoqKiQkFBQbXOGxkZqZCQkPqUVCOr1aq8vDyZTCb5+fm5dG5cQZ/dw+N93lS387B69+5tbB0G83ifbyD02j2M6nNpaWmdDxY4HTjmz5+v9evXKzU1VcOGDbvqmDZt2ujs2bN2286ePavWrVvXOrefn59hXzgj58YV9Nk9GnqfG3JtzmjoffYm9No9XN1nZ+Zyah2OtLQ0bdiwQYsXL9aIESNqHBcdHa2DBw+qurpa0qVLag8cOKDo6GhndgcAALxEnQNHUVGRli1bpmeffVYxMTEym822h3TpRNHy8nJJUnx8vC5cuKAFCxaosLBQCxYsUFlZmYYPH27MpwAAAA1anQPHxx9/LKvVqvT0dMXFxdk9JCkuLk6ZmZmSpNDQUK1YsULZ2dkaNWqUcnNztXLlSpefnwEAABqHOp/DkZiYqMTExBpfP378uN3zXr166f33369/ZQAAwGtw8zYAAGA4AgcAADAcgQMAABiOwAEAAAxH4AAAAIYjcAAAAMMROAAAgOEIHAAAwHAEDgAAYDgCBwAAMByBAwAAGI7AAQAADEfgAAAAhiNwAAAAwxE4AACA4QgcAADAcAQOAABgOAIHAAAwHIEDAAAYrt6Bo6KiQiNHjtT+/ftrHJOUlKRu3brZPfbs2VPfXQIAgEbKvz5vslgsmjJligoKCmodV1RUpNTUVA0YMMC2rXnz5vXZJQAAaMScDhyFhYWaMmWKqqurax1XUVGhU6dOyWQyKTw8vN4FAgCAxs/pwJGVlaX+/ftr8uTJ6t27d43jiouL5ePjow4dOtR5bqvVKqvV6mxJ15zzp/8LY9Bn92gsfW7o9V1LY+mzN6DX7mFUn52Zz+nAMWbMmDqNKy4uVmhoqKZOnaqsrCzdcsstmjRpkgYPHlzje/Lz850tp87y8vIMmxtX0Gf3aOh9zsnJ8XQJLtHQ++xN6LV7eLLP9TqHoy6Ki4tVXl6uuLg4JSYmateuXUpKSlJGRoZMJtNV3xMZGamQkBCX1mG1WpWXlyeTySQ/Pz+Xzo0r6LN7eLzPm3bWaVhtRz8bA4/3+QZCr93DqD6XlpbW+WCBYYFj4sSJGjdunO0k0e7du+vIkSPauHFjjYHDz8/PsC+ckXPjCvrsHg29zw25Nmc09D57E3rtHq7uszNzGRY4fH19Ha5IiYiIUGFhoVG7BHAdOk3f4ekSAHgxwxb+mj59umbMmGG37dixY4qIiDBqlwAAoIFyaeAwm80qLy+XJA0ZMkTbt2/Xtm3bdOLECaWlpSk7O1tjx4515S4BAEAj4NLAERcXp8zMTEnS0KFDNWfOHKWnp2vkyJHavXu3Vq9erfbt27tylwAAoBG4rnM4jh8/XuvzhIQEJSQkXM8uAACAF+DmbQAAwHAEDgAAYDgCBwAAMByBAwAAGI7AAQAADEfgAAAAhiNwAAAAwxE4AACA4QgcAADAcAQOAABgOAIHAAAwHIEDAAAYjsABAAAMR+AAAACGI3AAAADDETgAAIDhCBwAAMBwBA4AAGC4egeOiooKjRw5Uvv3769xzNGjR5WQkKDo6GiNHj1ahw8fru/uAABAI1avwGGxWPTiiy+qoKCgxjGlpaVKTExUbGystm7dqj59+mjChAkqLS2td7EAAKBxcjpwFBYW6tFHH9U333xT67jMzEwFBgZq6tSp6ty5s2bNmqWmTZtq586d9S4WAAA0Tk4HjqysLPXv318ZGRm1jsvNzVVMTIx8fHwkST4+Purbt69ycnLqVSgAAGi8/J19w5gxY+o0zmw2q0uXLnbbwsLCav0Zxmq1ymq1OltSrS7P5+p5YY8+u0dj6XNDr+9aGkufvQG9dg+j+uzMfE4HjroqKytTQECA3baAgABVVFTU+J78/HyjylFeXp5hc+MK+uweDb3P3nIks6H32ZvQa/fwZJ8NCxyBgYEO4aKiokJBQUE1vicyMlIhISEurcNqtSovL08mk0l+fn4unRtX0Gf3MLTPm1x3flXv3r1dNpcn8H12H3rtHkb1ubS0tM4HCwwLHG3atNHZs2fttp09e1atW7eu8T1+fn6GfeGMnBtX0Gf3aOh9bsi1OaOh99mb0Gv3cHWfnZnLsIW/oqOjdfDgQVVXV0uSqqurdeDAAUVHRxu1SwAA0EC5NHCYzWaVl5dLkuLj43XhwgUtWLBAhYWFWrBggcrKyjR8+HBX7hIAADQCLg0ccXFxyszMlCSFhoZqxYoVys7O1qhRo5Sbm6uVK1e6/BwNAADQ8F3XORzHjx+v9XmvXr30/vvvX88uAACAF+DmbQAAwHAEDgAAYDgCBwAAMByBAwAAGI7AAQAADEfgAAAAhiNwAAAAwxE4AACA4QgcAADAcAQOAABgOAIHAAAwHIEDAAAYjsABAAAMd113iwWAq+k0fcc1x3z96gg3VAKgoeAIBwAAMByBAwAAGI7AAQAADEfgAAAAhiNwAAAAwzkdOCwWi2bOnKnY2FjFxcVpzZo1NY5NSkpSt27d7B579uy5roIBAEDj4/RlsYsWLdLhw4e1du1anT59WtOmTVPbtm0VHx/vMLaoqEipqakaMGCAbVvz5s2vr2IAANDoOBU4SktLtWnTJq1atUpRUVGKiopSQUGB1q1b5xA4KioqdOrUKZlMJoWHh7u0aAAA0Lg49ZPKsWPHVFlZqT59+ti2xcTEKDc3V1VVVXZji4uL5ePjow4dOrimUgAA0Gg5dYTDbDarRYsWCggIsG1r1aqVLBaLzp8/r5YtW9q2FxcXKzQ0VFOnTlVWVpZuueUWTZo0SYMHD65xfqvVKqvVWo+PUbPL87l6Xtijz+7hTX1uyJ/Bm/rc0NFr9zCqz87M51TgKCsrswsbkmzPKyoq7LYXFxervLxccXFxSkxM1K5du5SUlKSMjAyZTKarzp+fn+9MOU7Jy8szbG5cQZ/dwxv6nJOT4+kSrskb+txY0Gv38GSfnQocgYGBDsHi8vOgoCC77RMnTtS4ceNsJ4l2795dR44c0caNG2sMHJGRkQoJCXGmpGuyWq3Ky8uTyWSSn5+fS+fGFfTZPQzt86adrp3vGnr37u3W/TmD77P70Gv3MKrPpaWldT5Y4FTgaNOmjc6dO6fKykr5+196q9lsVlBQkJo1a2Y31tfX1+GKlIiICBUWFtY4v5+fn2FfOCPnxhX02T28oc+NoX5v6HNjQa/dw9V9dmYup04a7dGjh/z9/e0OhWZnZ8tkMsnX136q6dOna8aMGXbbjh07poiICGd2CQAAvIBTgSM4OFgPPfSQ5s6dq0OHDumjjz7SmjVr9MQTT0i6dLSjvLxckjRkyBBt375d27Zt04kTJ5SWlqbs7GyNHTvW9Z8CAAA0aE6vNDpjxgxFRUVp/PjxmjdvniZNmqShQ4dKkuLi4pSZmSlJGjp0qObMmaP09HSNHDlSu3fv1urVq9W+fXvXfgIAANDgOb3SaHBwsBYuXKiFCxc6vHb8+HG75wkJCUpISKh/dQAAwCtw8zYAAGA4AgcAADAcgQMAABiOwAEAAAxH4AAAAIYjcAAAAMMROAAAgOEIHAAAwHAEDgAAYDgCBwAAMByBAwAAGI7AAQAADEfgAAAAhiNwAAAAwxE4AACA4fw9XQAA43WavsPTJQC4wXGEAwAAGI7AAQAADOd04LBYLJo5c6ZiY2MVFxenNWvW1Dj26NGjSkhIUHR0tEaPHq3Dhw9fV7EAAKBxcjpwLFq0SIcPH9batWs1Z84cpaWlaefOnQ7jSktLlZiYqNjYWG3dulV9+vTRhAkTVFpa6pLCAQBA4+FU4CgtLdWmTZs0a9YsRUVF6YEHHtAzzzyjdevWOYzNzMxUYGCgpk6dqs6dO2vWrFlq2rTpVcMJAADwbk5dpXLs2DFVVlaqT58+tm0xMTFavny5qqqq5Ot7Jb/k5uYqJiZGPj4+kiQfHx/17dtXOTk5GjVqlIvKB9BY1eXKma9fHeGGSgC4g1OBw2w2q0WLFgoICLBta9WqlSwWi86fP6+WLVvaje3SpYvd+8PCwlRQUOAwb1VVlSTphx9+kNVqdeoDXMvluS9evGgXiOBa9Nk96tvn229unFfAf//99x7ZL99n96HX7mFUn8vLy+3mr41T/xYqKyuzCxuSbM8rKirqNPbn46RLJ6JK0jfffONMOU4pLCw0bG5cQZ/dw9k+v/ZAK4MqMVZ+fr5H98/32X3otXsY1WeLxaLQ0NBaxzgVOAIDAx0Cw+XnQUFBdRr783GS1Lx5c3Xq1EmBgYEkXAAAGomqqipZLBY1b978mmOdChxt2rTRuXPnVFlZKX//S281m80KCgpSs2bNHMaePXvWbtvZs2fVunVrxyL8/RUWFuZMKQAAoAG41pGNy5w6nNCjRw/5+/srJyfHti07O1smk8nhyER0dLQOHjyo6upqSVJ1dbUOHDig6OhoZ3YJAAC8gFOBIzg4WA899JDmzp2rQ4cO6aOPPtKaNWv0xBNPSLp0tOPyCSTx8fG6cOGCFixYoMLCQi1YsEBlZWUaPny46z8FAABo0HyqLx+CqKOysjLNnTtXH374oUJDQ/X000/rd7/7nSSpW7duSklJsV32eujQIc2ZM0dFRUXq1q2b5s2bp549e7r8QwAAgIbN6cDRmFgsFs2bN08ffvihgoKC9NRTT+mpp57ydFle58yZM1qwYIE+//xzBQYG6pe//KVefPFFBQYGero0r5SYmKiWLVvq1Vdf9XQpXqmiokIpKSn64IMP1KRJEz3yyCOaPHmybU0huMa3336ruXPn6osvvtDNN9+sJ554wvYfr3CNiooKjRo1Sn/4wx/Uv39/SdLJkyf1hz/8QTk5OWrbtq1mzpypuLg4t9TTOC/Or6OfLsN++vRpTZs2TW3btlV8fLynS/Ma1dXVSk5OVrNmzbRu3TqVlJRo5syZ8vX11bRp0zxdntfZsWOH9u7dq4cfftjTpXitP/7xj9q/f7/eeecd/fDDD5o8ebLatm2r3/72t54uzau88MILatu2rbZu3arCwkK99NJLateunR544AFPl+YVLBaLpkyZYrf2VXV1tZ577jlFRkZqy5Yt+uijj/T8888rMzNTbdu2Nbwmr70G1Zll2FF/xcXFysnJUUpKirp27arY2FglJyfrgw8+8HRpXuf8+fNatGiRTCaTp0vxWufPn9eWLVs0f/589erVSwMGDNBTTz2l3NxcT5fmVUpKSpSTk6OkpCR16tRJv/jFLzRo0CDt27fP06V5hcLCQj366KMOa1t9/vnnOnnypF555RV17txZEyZMUO/evbVlyxa31OW1gaOmZdhzc3PrtCIa6iY8PFyrV69Wq1b2C0tdvHjRQxV5r4ULF+rXv/61wwq+cJ3s7GyFhoaqX79+tm2JiYlKSUnxYFXeJygoSMHBwdq6dat+/PFHFRcX68CBA+rRo4enS/MKWVlZ6t+/vzIyMuy25+bmqmfPngoJCbFti4mJsbvy1EheGziutQw7XKNZs2YaNGiQ7XlVVZXee+893XXXXR6syvvs27dPX375pSZOnOjpUrzayZMn1a5dO23btk3x8fG6//77tXTpUv4jxcUCAwM1e/ZsZWRkKDo6WsOHD9c999yjhIQET5fmFcaMGaOZM2cqODjYbrvZbHZYCyssLEzfffedW+ry2nM4nFmGHa6Tmpqqo0ePavPmzZ4uxWtYLBbNmTNHs2fPvupKvXCd0tJSnThxQhs2bFBKSorMZrNmz56t4OBgTjh3saKiIt1333168sknVVBQoPnz52vAgAF68MEHPV2a13LmliNG8NrA4cwy7HCN1NRUrV27VkuWLFFkZKSny/EaaWlpuuOOO+yOJMEY/v7+unjxol5//XW1a9dOknT69GmtX7+ewOFC+/bt0+bNm7V3714FBQXJZDLpzJkzSk9PJ3AYKDAw0OEIf023HDGC1wYOZ5Zhx/WbP3++1q9fr9TUVA0bNszT5XiVHTt26OzZs7bzkS4H53/+8586ePCgJ0vzOuHh4QoMDLSFDUm6/fbb9e2333qwKu9z+PBhdezY0e4PXc+ePbV8+XIPVuX92rRp43DztppuOWIErw0cP12GPTY2VlLNy7Dj+qSlpWnDhg1avHgxlxwb4K9//asqKyttz1977TVJ0ksvveSpkrxWdHS0LBaLvvrqK91+++2SLl2J9dMAguvXunVrnThxQhUVFbZD/MXFxWrfvr2HK/Nu0dHRWrlypcrLy21hLzs7WzExMW7Zv9f+5b3WMuxwjaKiIi1btkzPPvusYmJiZDabbQ+4Rrt27dSxY0fbo2nTpmratKk6duzo6dK8TkREhO69917NmDFDx44d06effqqVK1fqscce83RpXmXIkCFq0qSJXn75ZX311VfavXu3li9frnHjxnm6NK/Wr18/3XrrrZoxY4YKCgq0cuVKHTp0SI888ohb9u/VK43Wtgw7XGPlypV6/fXXr/ra8ePH3VzNjWH69OmSxEqjBvn+++81f/587dq1S8HBwRozZoyee+45Vhp1scv32Dp06JBatmypxx9/XOPHj6fPLtatWzf95S9/sa00euLECc2aNUu5ubnq2LGjZs6cqbvvvtsttXh14AAAAA2D1/6kAgAAGg4CBwAAMByBAwAAGI7AAQAADEfgAAAAhiNwAAAAwxE4AACA4QgcAADAcAQOAABgOAIHAAAwHIEDAAAY7v8DyE4u+nh2I+oAAAAASUVORK5CYII=",
"text/plain": [
"<Figure size 640x480 with 2 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"kwds = dict(bins=np.linspace(0, 10), density=True)\n",
"\n",
"top = plt.subplot(211)\n",
"plt.hist(idata.prior[\"mu\"].values[0], **kwds)\n",
"plt.title(\"Prior (Exp(1))\")\n",
"\n",
"plt.subplot(212, sharex=top, sharey=top)\n",
"plt.hist(idata.posterior[\"mu\"].values.ravel(), **kwds)\n",
"plt.title(\"Posterior\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"So, Bayesian prior distributions don't _have_ to look like what you think the posterior\n",
"distribution will look like.\n",
"\n",
"Sometimes there may be reasons to prefer a certain prior distribution; for instance if\n",
"you know your data cannot be negative, then an Exponential may be better than a Normal."
]
}
],
"metadata": {
"kernelspec": {
"display_name": "pymc-env",
"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.10.4"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment