From f118694e33fe369b1d7ed330866f18a38af9e02f Mon Sep 17 00:00:00 2001
From: Mirko Birbaumer <mirko.birbaumer@hslu.ch>
Date: Tue, 11 Mar 2025 21:57:08 +0000
Subject: [PATCH] =?UTF-8?q?analytische=20Rechnung=20f=C3=BCr=20BF=20hinzug?=
 =?UTF-8?q?ef=C3=BCgt?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit

---
 notebooks/Model Comparison/MC_2_1.ipynb | 195 +++++++++++++-----------
 1 file changed, 108 insertions(+), 87 deletions(-)

diff --git a/notebooks/Model Comparison/MC_2_1.ipynb b/notebooks/Model Comparison/MC_2_1.ipynb
index d66ef90..2275e72 100644
--- a/notebooks/Model Comparison/MC_2_1.ipynb	
+++ b/notebooks/Model Comparison/MC_2_1.ipynb	
@@ -22,11 +22,19 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 13,
+   "execution_count": 3,
    "metadata": {
     "tags": []
    },
-   "outputs": [],
+   "outputs": [
+    {
+     "name": "stderr",
+     "output_type": "stream",
+     "text": [
+      "WARNING (pytensor.tensor.blas): Using NumPy C-API based implementation for BLAS functions.\n"
+     ]
+    }
+   ],
    "source": [
     "%matplotlib inline\n",
     "import pymc as pm\n",
@@ -42,53 +50,39 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 14,
+   "execution_count": 10,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from scipy.special import betaln\n",
+    "def beta_binom(prior, y):\n",
+    "    \"\"\"\n",
+    "    Calculate the marginal probability, analytically, for a beta-binomial model.\n",
+    "    prior : tuple\n",
+    "      alpha and beta parameters for the beta prior\n",
+    "    y : array\n",
+    "      array with \"1\" and \"0\" corresponding to success and failure respectively\n",
+    "    \"\"\"\n",
+    "    alpha, beta = prior\n",
+    "    h = np.sum(y)\n",
+    "    n = len(y)\n",
+    "    p_y = np.exp(betaln(alpha + h, beta + n - h) - betaln(alpha, beta))\n",
+    "\n",
+    "    return p_y"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 12,
    "metadata": {
     "tags": []
    },
    "outputs": [
     {
-     "name": "stderr",
-     "output_type": "stream",
-     "text": [
-      "Initializing SMC sampler...\n",
-      "Sampling 2 chains in 2 jobs\n"
-     ]
-    },
-    {
-     "data": {
-      "text/html": [
-       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
-      ],
-      "text/plain": []
-     },
-     "metadata": {},
-     "output_type": "display_data"
-    },
-    {
-     "name": "stderr",
-     "output_type": "stream",
-     "text": [
-      "We recommend running at least 4 chains for robust computation of convergence diagnostics\n",
-      "Initializing SMC sampler...\n",
-      "Sampling 2 chains in 2 jobs\n"
-     ]
-    },
-    {
-     "data": {
-      "text/html": [
-       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
-      ],
-      "text/plain": []
-     },
-     "metadata": {},
-     "output_type": "display_data"
-    },
-    {
-     "name": "stderr",
+     "name": "stdout",
      "output_type": "stream",
      "text": [
-      "We recommend running at least 4 chains for robust computation of convergence diagnostics\n"
+      "11\n"
      ]
     }
    ],
@@ -106,55 +100,77 @@
     "alpha_2 = 8\n",
     "beta_2 = 4\n",
     "\n",
-    "with pm.Model() as model_BF_0:\n",
-    "   θ = pm.Beta('θ', alpha_1, beta_1)\n",
-    "   y = pm.Bernoulli('y', θ, observed=y_d)\n",
-    "   trace_smc_0 = pm.sample_smc(chains=2)\n",
-    "\n",
-    "with pm.Model() as model_BF_1:\n",
-    "   θ = pm.Beta('θ', alpha_2, beta_2)\n",
-    "   y = pm.Bernoulli('y', θ, observed=y_d)\n",
-    "   trace_smc_1 = pm.sample_smc(chains=2)"
+    "BF= beta_binom((alpha_1, beta_1), y_d)/ beta_binom((alpha_2, beta_2), y_d)\n",
+    "print(round(BF))"
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 12,
-   "metadata": {
-    "tags": []
-   },
+   "execution_count": 17,
+   "metadata": {},
    "outputs": [
     {
      "name": "stdout",
      "output_type": "stream",
      "text": [
-      "11.16\n"
+      "Collecting preliz\n",
+      "  Downloading preliz-0.11.0-py3-none-any.whl.metadata (6.1 kB)\n",
+      "Requirement already satisfied: matplotlib>=3.5 in /opt/conda/lib/python3.10/site-packages (from preliz) (3.9.1.post1)\n",
+      "Collecting numba>=0.59 (from preliz)\n",
+      "  Downloading numba-0.61.0-cp310-cp310-manylinux2014_x86_64.manylinux_2_17_x86_64.whl.metadata (2.8 kB)\n",
+      "Requirement already satisfied: numpy>=1.22 in /opt/conda/lib/python3.10/site-packages (from preliz) (1.26.4)\n",
+      "Collecting scipy<1.13,>=1.9.1 (from preliz)\n",
+      "  Downloading scipy-1.12.0-cp310-cp310-manylinux_2_17_x86_64.manylinux2014_x86_64.whl.metadata (60 kB)\n",
+      "Requirement already satisfied: contourpy>=1.0.1 in /opt/conda/lib/python3.10/site-packages (from matplotlib>=3.5->preliz) (1.2.1)\n",
+      "Requirement already satisfied: cycler>=0.10 in /opt/conda/lib/python3.10/site-packages (from matplotlib>=3.5->preliz) (0.12.1)\n",
+      "Requirement already satisfied: fonttools>=4.22.0 in /opt/conda/lib/python3.10/site-packages (from matplotlib>=3.5->preliz) (4.56.0)\n",
+      "Requirement already satisfied: kiwisolver>=1.3.1 in /opt/conda/lib/python3.10/site-packages (from matplotlib>=3.5->preliz) (1.4.8)\n",
+      "Requirement already satisfied: packaging>=20.0 in /opt/conda/lib/python3.10/site-packages (from matplotlib>=3.5->preliz) (24.2)\n",
+      "Requirement already satisfied: pillow>=8 in /opt/conda/lib/python3.10/site-packages (from matplotlib>=3.5->preliz) (8.3.2)\n",
+      "Requirement already satisfied: pyparsing>=2.3.1 in /opt/conda/lib/python3.10/site-packages (from matplotlib>=3.5->preliz) (3.2.1)\n",
+      "Requirement already satisfied: python-dateutil>=2.7 in /opt/conda/lib/python3.10/site-packages (from matplotlib>=3.5->preliz) (2.9.0.post0)\n",
+      "Collecting llvmlite<0.45,>=0.44.0dev0 (from numba>=0.59->preliz)\n",
+      "  Downloading llvmlite-0.44.0-cp310-cp310-manylinux_2_17_x86_64.manylinux2014_x86_64.whl.metadata (4.8 kB)\n",
+      "Requirement already satisfied: six>=1.5 in /opt/conda/lib/python3.10/site-packages (from python-dateutil>=2.7->matplotlib>=3.5->preliz) (1.17.0)\n",
+      "Downloading preliz-0.11.0-py3-none-any.whl (514 kB)\n",
+      "Downloading numba-0.61.0-cp310-cp310-manylinux2014_x86_64.manylinux_2_17_x86_64.whl (3.8 MB)\n",
+      "\u001b[2K   \u001b[90m━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━\u001b[0m \u001b[32m3.8/3.8 MB\u001b[0m \u001b[31m31.7 MB/s\u001b[0m eta \u001b[36m0:00:00\u001b[0m\n",
+      "Downloading scipy-1.12.0-cp310-cp310-manylinux_2_17_x86_64.manylinux2014_x86_64.whl (38.4 MB)\n",
+      "\u001b[2K   \u001b[90m━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━\u001b[0m \u001b[32m38.4/38.4 MB\u001b[0m \u001b[31m66.8 MB/s\u001b[0m eta \u001b[36m0:00:00\u001b[0m:00:01\u001b[0m00:01\u001b[0m\n",
+      "Downloading llvmlite-0.44.0-cp310-cp310-manylinux_2_17_x86_64.manylinux2014_x86_64.whl (42.4 MB)\n",
+      "\u001b[2K   \u001b[90m━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━\u001b[0m \u001b[32m42.4/42.4 MB\u001b[0m \u001b[31m51.5 MB/s\u001b[0m eta \u001b[36m0:00:00\u001b[0m:00:01\u001b[0m00:01\u001b[0m\n",
+      "Installing collected packages: scipy, llvmlite, numba, preliz\n",
+      "  Attempting uninstall: scipy\n",
+      "    Found existing installation: scipy 1.13.1\n",
+      "    Uninstalling scipy-1.13.1:\n",
+      "      Successfully uninstalled scipy-1.13.1\n",
+      "Successfully installed llvmlite-0.44.0 numba-0.61.0 preliz-0.11.0 scipy-1.12.0\n"
      ]
+    },
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA/sAAAESCAYAAACirMVZAAAAP3RFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMS5wb3N0MSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8kixA/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAB4J0lEQVR4nO3dd3gU5drH8e/upncChITee0epIiiiKCLYG1UBCx5FX/sR8VjAjhz0WLChKAgqWFBQEBBEeu+d0EsgCenJ7rx/rEYmm0ACSSa7+X2uK1fI/czs3rtMNnPPPMVmGIaBiIiIiIiIiPgMu9UJiIiIiIiIiEjxUrEvIiIiIiIi4mNU7IuIiIiIiIj4GBX7IiIiIiIiIj5Gxb6IiIiIiIiIj1GxLyIiIiIiIuJjVOyLiIiIiIiI+Bi/wmzkcrk4dOgQ4eHh2Gy2ks5JRERERETKKMMwOH36NFWrVsVu171DkbKqUMX+oUOHqFGjRknnIiIiIiIiXmL//v1Ur17d6jREpACFKvbDw8MB9y90REREiSYkIiIiIiJlV3JyMjVq1MitEUSkbCpUsf931/2IiAgV+yIiIiIiouG9ImWcBtmIiIiIiIiI+BgV+yIiIiIiIiI+RsW+iIiIiIiIiI9RsS8iIiIiIiLiY1Tsi4iIiIiIiPiYQs3GLyJln8tlkJyRzcnULE6lZZOR7cTpMnAaBhgQHOAgLNCP0EA/ooL9iQrx1yy6IlI6DAPST0FG0j9fzmx3mw2w2SEwAoIrQFAUBEeB3WFhwiIiIt5Pxb6Il8nIdrLpUDIbDiSy63gqexNS2XMilcNJGThdRqEfJ9jfQVxUEFUjg6lbOZRGseE0jg2nYZVwwoP8S/AViIjPykqFo5vhyHo4uglO7YXEeEjaDzkZhX8cux9E1oAKtSG6DsQ0hbjWENsc/INLKHkRERHfYjMM45zVQXJyMpGRkSQlJREREVEaeYnIX9KznCzdncDC7cdZue8kWw+fJqcIRX1R2WzQqEo47etEc3HtaDrUjSYmPKjEnk9EvFjKMdj3B+z9w/392Bag5D6fsDmgSlOofSnUuwxqdYaA0JJ7PhHJl2oDEe+gYl+kDEpIyeTnjUeYs+kIy/acJCvHZWk+LatHcnnjGC5vHEOLapHq/i9SXhkGHN0IW3+CrT+67+BbyRHgLvib9oMm10FoRWvzESknVBuIeAcV+yJlRGaOk9kbj/DN6oP8sfNEkbrk5xXs7yA00IHdZsNht2ED0rKdpGTkXHCvgOoVgrmuVVX6talGwyrhF/RYIuIljm2FdV/CxhmQFH9hj2X3g6BI8AtyXzwAcOX8NY4/8/wf1+aAut2h9R3QpA/4BV5YniJSINUGIt5Bxb6IxeIT0vhyeTzTV+4nITWrUPvYbNAwJpymVSOoUymU2pVCqV0xhMrhgVQICSDIP/+JrQzDIDPHxYmUTA4nZXAoMZ34hDS2HT3NtiOn2X0itUgXGZrERXBHh5pc36YaYYGaAkTEp6Sfgg1fw9ov4dDqwu9ns0PFBhDbAio3hgq1IKomRFaH4Gj3mPuCegdlp7ufN/kwnNrj/jqxAw6vg+PbKPQQgZBK0KY/tBvsHvMvIsVKtYGId1CxL2KRjQeT+N+Cnfy88Qjn+i102G20rRnFpQ0qc3GdaJpXiyyR4joj28nGg0ks33uSFXtOsmLvKVIyc865X2iAgxvaVqd/x1o0itXdfhGvdnwbLHsP1k2F7LRzbx8QBjU7Qq0u7q+4liUziV5mirvo37sYds+HAyvcPQLOygZN+0LX/3PnJSLFQrWBiHdQsS9SylbHn2LCvB3M33b8rNsF+zvo0SSGa1rEcUmDSkRYMEN+Vo6LlXtPMm/rMeZtOcrehHOf+F/asDIjutejfZ1oje0X8RaG4S6g/3wHds499/ZRNaFRb2h8DdTsBA4LVvDISIZdv8Gmb2H7L5CTfvbtG1wJlz4ONS4unfxEfJhqAxHvoGJfpJTsPp7CK7O3MmfT0bNu16luRW5rX4OeTasQElB2usYbhsGmQ8l8t/Yg3687xNHks4+tbVerAvd3r8fljWNU9IuUVYYBO+fBgrFwcOXZtw2rAi1vgZa3QpXmBXfFt0JmCmz7GVZPgr2Lzr5tkz5wxX+gYr3SyU3EB6k2EPEOKvZFSlhCSibj5+3gy2XxBU6OFxbox80XVefODrWoHxNWyhkWndNlsHjnCSYv3ce8LUc52zD/FtUiebxXI7o2qFx6CYrI2RmG+w7+grFwcFXB29n9oHFvaN0f6l0OjrJzAbJAx7fByo9hzReQdTr/bex+cNFd0O0JCK1UuvmJ+ADVBiLeQcW+SAlxuQymrdzP2J+3kpSene820aEB3NWlNgM61SYy2IJusMXgYGI6U5fH8+Wy+LNOMNilfkWe6NWYltWjSi85EfF0aA38Mursd8CDK0C7IXDxUIisVnq5Faf0RFgxEf78H6SfzH+boCi48gX3xQy7vTSzE/Fqqg1EvIOKfZESsP3oaf49YwMr9p7Kt71CiD8jLqvPHR1qlqmu+hciI9vJtJX7eX/hbg4mFjx29tqWcTx9TROqRpXABF4iUrDE/fDbC7D+q4K3qVAbuox0d9UPCCmtzEpWViqs/AQWvVFw0V+zE1z7FsQ0LtXURLyVagMR76BiX6QY5ThdvDN/FxN+25Fvl/0APzt3danD/ZfVs2TCvdKQ7XTx/dpDvLNgJ7uPp+a7TbC/gwcur8/QrnUI9Mt/mUARKSbZGbD4TfhjPORk5L9NhTpw6WPuMflWTLZXGtITYfE4WPouOPOZc8Tu557Ar+v/ecdwBRELqTYQ8Q4q9kWKye7jKTw8bR3r9ifm235tyzieuqYJ1crJHe0cp4uvVx3grbk7OJKcf4FRu2IIz/ZpyuWNq5RydiLlxM558NOjcHJ3/u2RNaD7U+47+eWlwE3cD78+657FPz/VLoLr34dK9Us3LxEvotpAxDuo2Be5QIZhMHlZPC/N2kxGtsujvUZ0MC/2a0G3huVzgrqMbCeTluzlfwt2FTh3wVXNqvBC3+bERASVcnYiPur0EZj9VMEFbWAEdH0EOtwH/uX0927HXJj1CCTu82zzD3GP5b/o7rK16oBIGaHaQMQ7qNgXuQDJGdk8Nn1dvsvpOew2hl9alwcvb0BwgLqqn0rN4o1ft/HFsnjy+9SJCPJj1LVNualddS3VJ3K+DANWfeq+c52Z7Nluc7hnoe/+pGahB8hKg99fhT/+C4bTs73JddD3HQjSuY/ImVQbiHgHFfsi52nzoWTu/2IVexPSPNrqVApl3K2taV0jqvQTK+M2Hkzi2e82sjo+Md/2SxtWZsz1zalewUcmBxMpLUkH4ft/wa55+bdXbw/XjoPY5qWblzc4sAq+HQYnd3m2RdeDWz+HKs1KPy+RMkq1gYh3ULEvch6mrdzPqJkbyczx7LY/oGMtnrqmsc/Msl8SXC6Db9ccZMxPWziZz3J9oQEORvdpxs0X6S6/yDkZBqybAj8/CZlJnu1BUdDzeWgzQMvLnU1WqrtHxIoPPdv8guHaN6H1HaWfl0gZpNpAxDuo2BcpghynixdnbeHTJXs92qJC/Bl3S2suaxxT+ol5qYSUTP7zw2a+X3co3/arm8cy5voWVAgNKOXMRLxEynH44UHY9lP+7a1uhytfVJf9otj+C8wYDun5LJ3a4T646iWwa2iWlG+qDUS8g4p9kUJKzsjmgS/X8Pv24x5trapH8s6dbdX1/DzN3XyUf8/cwNFkz+WwqkQE8sbNrbmkgYoVEZNd8+Hb4ZB6zLMtPA6umwANepZ+Xr4gMR6mDYJDqz3b6veEmz7WOH4p11QbiHgHFfsihRCfkMZdk1aw81iKR9vATrX4d+8mWi/+AiVnZPPSj1v4auX+fNuHda3DY1c1JsBP3ZClnHNmw/wx7jXjyedPeMvb4OqXIbhCqafmU3IyYc7T+Xfrr9wE7pgKFWqXeloiZYFqAxHvoGJf5Bw2HEhi8CfLScgztjzAYWfsDS24sV11izLzTbM3HuGpb9dzKs1zmb7WNaJ45862VIsKtiAzkTLg1D745m44sMKzLbQy9BkPjXuXfl6+bO0U91AJZ575RUIqwZ3ToVpba/ISsZBqAxHvoGJf5CwW7zjBPZ+vJDXLvCRTxdAAPhjYjna1oi3KzLcdTc7g0enrWLTjhEeb5kaQcmvrLJhxX/6T8DW4Evq9q7H5JWXfn/DVnZCWYI4HhMGtk6HeZdbkJWIR1QYi3kH9YUUK8MO6Qwz5dLlHod84NpzvHuiiQr8EVYkIYtKQ9jx7bVOPbvuJadkM+XQFr87eSo7TczUEEZ/jcsK852HqHZ6Fvt0frhoDt3+lQr8k1eoEw35zd98/U1YKfHEzbPzGmrxERETOQnf2RfIxaclenvthE3l/O7rUr8h7/dsRHuRvTWLl0KZDSdz/xWr2JaR5tHWsG807d7SlYligBZmJlIK0k+5u+7t+82yrUAdu/gSqtin9vMqrjGSYNgB2L8jTYIOrX4UOw63ISqTUqTYQ8Q66sy+Sx3sLdzH6e89C/9qWcXw8+GIV+qWsWdVIfvjXJVzdPNajbenuk1z39h9sOpRPt2YRb3doLbzfLf9Cv/lNcM/vKvRLW1AE3DENmt2Qp8GAnx+DP8ZbkpaIiEh+VOyLnOHt33bw8s9bPeKDOtXiv7e10Yz7FokI8ud/d7bl2Wub4me3mdoOJqZz47tL+GHdIYuyEykB66bCR1dCUrw5bveDq1+DGz/U0m9W8QuEGz+C9vncxf/1WVj0RunnJCIikg8V+yJ/eWvudl7/ZbtH/P96NuS565phz1NkSumy2WzcdUkdpt3bidiIIFNbRraLf01Zwyuzt+J0nXNkkkjZ5XLB3P/AjHvAmWluC6sCg2e5u4rb9HlkKbvd3W3/sn97ts17Hha8Uvo5iYiI5KFiX8o9wzB485dtvDV3h0fbqGub8q8eDbDpxLrMaFuzAt//qwvtanmuIf7ugl0MnbSClMwcCzITuUCZKe7x4Ivf9Gyr0dHdbb9mx9LPS/Jns0G3x+HKlzzbFoyB317CYzyYiIhIKVKxL+Xe27/t5L+/7fSI/+e6Ztx9SR0LMpJziQkP4sthHbi9fQ2PtvnbjnPTu0s4nJRuQWYi5ynpAHzSC7b+6NnWfjgM+gHCPeetkDKg8wPQK587+b+/Cr+/Xvr5iIiI/EXFvpRrHy/ewxu/enbdf7FfcwZ1rl36CUmhBfo5GHtDS17s19xjHP/WI6fp984fbDyoifvECxxYBRMvhyMbzHG7H1w7Dq55DfwCrMlNCqfjvXBNPoX9/Bdh6buln4+IiAgq9qUcm7ZiP8//uNkjPvaGFvTvWMuCjOR89O9Yiy+GdqBCiHmVhKPJmdzy/p/M23LUosxECmHLD/DpNZCS5zgNioL+38JFd1mSlpyH9sPg2rc847OfhNWfl3o6IiIiKvalXPph3SGe+Ha9R/z5vs24vX1NCzKSC9GhbkVm3N+FupVCTfG0LCfDPlvJp3/ssSgzkbNY9gF8NQByMszxivVh6Dyo282avOT8XTTEvVpCXt//CzZ+U/r5iIhIuaZiX8qdBduO8fBXaz3mTXq8VyMGdqptSU5y4WpXCuXb+zvTvk60Ke4y4LkfNvOfHzbh0kz9Uha4XPDraPe67OQ5Jut0g6FzoVJ9S1KTYtBhOPR4Nk/QgG+Hw465lqQkIiLlk4p9KVc2Hkzi/i9Wk5On6BtxWT3u766Ta28XFRLA53e354Y21TzaPvljLw99tZasHJcFmYn8JScLZt4Lf7zl2dZ2IPT/BoI9V5oQL9P1/+CSR8wxVw5MGwiH1liTk4iIlDsq9qXc2H8yjcGfrCAty2mKD+5cm0evbGRRVlLcAv0cvHFLK0Ze0cCj7Yd1h7jrUy3NJxbJSIYvb4b1X3m2dX8a+vwXHP6ebeKdejwL7e8xx7JT4Yub4aSGFomISMlTsS/lwqnULAZ9spwTKZmmeN/WVXn22qbYbLYC9hRvZLPZGHlFQ8bd2gp/h/n/dvHOE9wxcanHsSBSok4fgU+ugd0LzHGbA66bAN2fcK/bLr7DZoNeL0OLm83x1OMw+UZITbAmLxERKTdU7IvPy8h2MvSzlew+nmqKd6pbkVdvaondrhNsX3V9m+p8NOhiQgIcpvj6A0nc/N6f7D+ZZlFmUq6c2gsfXwVH8yyt5x8Ct091d98X32S3Q993oM6l5vjJXTDlVsjSZ5CIiJQcFfvi01wug5FT17Jq3ylTvHFsOO8PbEegn6OAPcVXXNqwMlOGdSQ61LxO+Z4Tqdzw7hI2H0q2KDMpF45tgY+uchf8ZwqtDINnQcMrLUlLSpFfINw6Gao0N8cPrIAZw90TNoqIiJQAFfvi0177ZRuzNx0xxeIig/hkyMVEBGlsbHnRqkYUX9/biWpRwab48dOZ3PrBnx4Xg0SKxcFV8MnVkGL+DCK6Ltz9C1Rra01eUvqCIuHOryGiujm+5QdYMMaanERExOep2BefNWPNAd5dsMsUCw/y49Mh7YmLDC5gL/FVdSuH8e39nWkcG26Kn87IYcBHy1iy64RFmYlP2rMIJl0H6XkuJFVpAXfNcRf8Ur5ExLlXWwiKNMd/fw3WT7cmJxER8Wkq9sUnrY4/xRPfmMfH+tltvN+/HY3yFHtSflSJCOKrezrRvna0KZ6W5WTIJyuYv/WYRZmJT9n2s3sCtqwUc7x6exj8A4TFWJOXWC+mMdzymXtixjN9NwL2r7AmJxER8Vkq9sXnHEpMZ/hnqzzWU3++b3M6169kUVZSVkQG+/PZ3e25rFFlUzwzx8Xwz1fy04bDFmUmPmH9dJh6JzjzrPZQ73IYOBOCK1iSlpQhdbvDNa+aY85MmHoHJO63JCUREfFNKvbFp6Rl5TDss5Uey6oN7lybOzrUtCgrKWuC/B28P+AirmkRa4pnOw0e+HI136w6YFFm4tVWTYJvh4HhNMebXOeedT8g1Jq8pOy5eCi0H26OpR6DKbdDZkr++4iIiBSRin3xGYZh8Nj09WzKM7v6JfUr8UzvJhZlJWVVgJ+d/97WhhvaVjPFXQb83/R1fL50n0WZiVda8SH88CBgmOOt+8NNn7hnZBc501Vj3T0+znR0A3z/LzCM/PcREREpAhX74jPeXbiLWXm6YNepFMo7d7TFz6FDXTz5Oey8flMr+nf07PUxauZGPvh9Vz57ieSx7H2Y9X+e8Y73w3UTwOFX+jlJ2efwc18IqtjAHN/0Lfz5jjU5iYiIT1EFJD7hj50neH3ONlMsPMiPDwddRGSIltiTgtntNl7o25zhl3rOjj7mp60eKzqImPz5Dvz8uGf80sfgqjFg159ZOYvgKLjjK88Z+n991r2ig4iIyAXQWYh4vcNJ6Tw4ZQ2uM3o92mzw9h1tqVc5zLrExGvYbDaeuroxD1/R0KPtldlb+d+CnRZkJWXe4rdgztOe8e5Pw+XPuD+IRM6lYj24YaI5Zjhh+mBIOmhJSiIi4htU7ItXy8pxcf8Xq0lIzTLF/69nQ7o1rFzAXiKebDYbD13RgH9f4zm/w6uzt/HOfBX8cobfX4O5oz3jl4+C7k+Ufj7i3RpeBd2fMsfSTsC0gZCTmf8+IiIi56BiX7zai7M2syY+0RS7okkM93evb01C4vWGXVo33wkdX5ujgl9wT5y24GX47UXPtp7Pw6WPln5O4hsufRwaXGWOHVwJP+vikYiInB8V++K1Zqw5wGd/mmdMrxkdwhu3tMZuV/dZOX9Du9Zl1LVNPeKvzdnG27/tsCAjKRMMA+a/BAvGerZdNQa6PFT6OYnvsNvhhvehQh1zfNUnsH6aNTmJiIhXU7EvXmn70dM89e0GUyzQz857/dsRGawJ+eTC3X1JHZ7Np+B//Zft/HeeCv5yacFYd/f9vK5+DTqNKP18xPcEV4BbJ4NfsDn+w0g4oc8dEREpGhX74nXSs5yM+GI1GdkuU3zM9S1oWjXCoqzEF911SR2e6+NZ8L/563bGz9WJd7my8DVY+IpnvPeb0GF46ecjviu2OVz3X3MsO9U9YV92uiUpiYiId1KxL17nPz9sYsexFFPszg41ubFddYsyEl82uEsd/nNdM4/4uLkq+MuNxeNgfj5j9PuMh4vvLv18xPe1vAXaDjTHjm6E2U/lv72IiEg+VOyLV/l+3SGmrthvijWNi8h3fLVIcRnUuTbP982/4H93wS4LMpJS8+c7MPc5z/i1b0G7waWcjJQrvV6BmDx/21Z9Ahu/sSYfERHxOir2xWvsS0jl6Tzj9EMCHLx9RxuC/B0WZSXlxcBOtXkhn4L/ldlb+XjxHgsykhK37H2Y87Rn/JrX4aIhpZ+PlC8BIXDzp+AfYo5//xAk6CKjiIicm4p98QpZOS7+NWUNKZk5pviL/ZpTt3KYRVlJeTOggIL/+R83M3npvnz2EK+14iP4+XHP+FVjof2w0s9HyqfKjdzzQpwp6/Rf4/czLElJRES8h4p98Qqvzt7K+gNJptiNbatzQ1uN05fSNaBTbZ7p3cQj/szMjUxbuT+fPcTrrP4MZj3iGe/5PHS6v/TzkfKt9e3Q+k5z7Mh6mDvamnxERMRrqNiXMu+3rUf5ME836bqVQ/MdQy1SGoZ2rcvjvRp5xJ/4Zj3frT1oQUZSbNZ+Cd8/6Bm/fBR0eaj08xEBuOY1qJTnM2fZe7BzrjX5iIiIV1CxL2XasdMZPDp9vSkW4Gfn7dvbEhroZ1FWInB/9/o81KOBKWYY8Mi0dfy84bBFWckFWT8dZt4PGOZ496fg0kctSUkEgIBQ9/h9vyBzfOb9kHrCkpRERKTsU7EvZZZhGDzx9XpOpmaZ4qN6N6Fp1QiLshL5x8grGnBvt3qmmNNl8K8pa5i7+ahFWcl52ToLZtyDR6Hf9VHo9oQlKYmYVGkKV+ZZAjLlKHz/L/eVRhERkTxU7EuZNXlZPPO3HTfFrmxahf4da1mUkYiZzWbjiV6NuKtLHVM8x2Vw/xerWbj9eAF7Spmye4F7wjPDaY53fhAufwZsNiuyEvF08VBocKU5tu0nWPWpJemIiEjZpmJfyqRdx1N4adZmU6xyeCAv39gSm068pQyx2WyMurYJ/TvWNMWznC6Gf7aSJbvUxbZM278CptwBTnMPIjre756QT583UpbYbND3HQipZI7PfgpO7LAmJxERKbNU7EuZk+108fBXa8nIdpnir93UkujQAIuyEimYzWbj+euac3M78+oQmTkuhk5ayer4UxZlJmd1ZCN8cSNkp5rjbQbAVWNU6EvZFBbjLvjPlJMO3wyFnKz89xERkXJJxb6UOf+dt8Njmb2BnWrRvVGMRRmJnJvdbuPlG1vSt3VVUzwty8ngj5ez5XCyRZlJvhJ2wefXQ4b5s4Zm10Of8Sr0pWxr1AsuutscO7wWFoyxJB0RESmbVOxLmbJq30nemb/TFKtXOZSnrvZc11ykrHHYbbxxcyuubh5riidn5DDgo2XsPp5iUWZiknQAPusLqcfM8fo94foPwO6wJi+RorjyRajU0Bxb/BbsW2JJOiIiUvao2JcyIyUzh4e/WofrjEmF/ew23rq1DcEBOvkW7+DnsDP+tjZ0a1jZFD+RkkX/D5dxMDHdoswEgJTj8Fk/SNpvjtfsDLd8Bn4aKiReIiAEbvwQ7P5nBA33cnxZqQXuJiIi5YeKfSkzXpq1mfiTaabYwz0b0qJ6pEUZiZyfAD877/VvR/va0ab4oaQM+n+4jOOnMy3KrJxLT4TJ10NCnonM4lrDHV+5iycRbxLXyr1ixJlO7YFfR1uTj4iIlCkq9qVMWLj9OFOWm++0XVSrgsca5iLeIjjAwUeDL6JFNfPFqj0nUhnw0TKS0rItyqycykqFL2+BIxvM8UqNoP+3EBRhTV4iF6rzv6BGB3NsxUT3kpIiIlKuqdgXyyVnZPPkN+tNsdAAB2/e0hqHXZNkifcKD/Jn0l3taRATZopvPXKawZ8uJzUzx6LMypmcTPiqP+xfZo5H1YKBMyG0oiVpiRQLuwP6vQt+web4dw94TkApIiLliop9sdyLP27mcFKGKfZ07ybUrKguteL9okMDmDy0AzWjzcfzmvhEhn22koxsp0WZlRPOHPjmbtj1mzkeFgsDv4OIqvnvJ+JNKtaDnv8xx5L2w5ynrclHRETKBBX7Yqn5W48xbeUBU+yS+pW4o31NizISKX5VIoL4YmgHqkQEmuJLdiXwwJeryXa6LMrMx7lc8MODsOUHczy4gvuOfnQdS9ISKREXD4PaXc2xNZNh+xxr8hEREcup2BfLJKVl8+S35u77YYF+vHxjC2xa41p8TI3oEL4Y2oHoUPNs73O3HOPR6etwnbkMhVw4w4A5T8HaL8zxgHDo/w3EaDlP8TF2O/R9BwLMw4b4/l+QdtKanERExFIq9sUyz/+4maPJ5lnJ/927CdUrqPu++Kb6MeF8dld7wgP9TPHv1h7ime82Yhgq+IvN/DGw7D1zzC8I7pgK1dpZk5NISatQC64aY46lHIWfHrMmHxERsZSKfbHE3M1H+Wa1uft+1waVuO3iGhZlJFI6mleL5OMhFxPkb/74/XJZPC//vFUFf3FYMgF+f9Ucs/vBLZ9B7UusyUmktLQdCPV7mmMbv4ats6zJR0RELKNiX0pdYloWT88wL38VHujHKze2VPd9KRcurh3N+wMuwt9hPt7f/30378zfaVFWPmLVJPglz7rj2OCGD6DhVZakJFKqbDa47r8QZF72kx8fgfRES1ISERFrqNiXUvefHzZz7LS5+/6oa5tSNSq4gD1EfE+3hpWZcHsb8q4u+fov2/n0jz3WJOXtNn4DPzzkGe/zFjS/sdTTEbFMRFXo9Yo5lnIknwthIiLiy1TsS6n6betRZqw5aIp1b1SZmy+qblFGItbp1TyOV29q5RF/7ofNfL3qQD57SIG2/wLfDgfyDIPo+QK0G2xFRiLWanUb1L/CHFvzOeyab00+IiJS6lTsS6lJyczhmRkbTbHwID/G3qDZ96X8uqlddf5zXTOP+ONfr2P2xsMWZOSF9i6GaQPAlWOOX/oYdHnQmpxErGazwbVvec7O/8ODkJliSUoiIlK6VOxLqXl19lYOJWWYYqN6NyUuUt33pXwb1Lk2j17Z0BRzGfDglLUs2nHcoqy8xMHV8OVtkGP+bKH9PXDZv63JSaSsiKoBPf9jjiXGw28vWpOPiIiUKr9zbyJy4VbuPcnnS/eZYl3qV1T3fZG/jLisPqczcnj/9925sSyni+GfrWLy0Pa0qxVtYXZl1LGtMPlGyDptjre6A3q97L6zKVLetbsLNn4L+/74J7bsPWh2PdTsYF1eUi44nU6ys7OtTkPEZ/j7++NwOAq9vc0oxDpPycnJREZGkpSURERExAUlKOVPZo6Ta8YvYtfx1NxYkL+dOSMvpVbFUAszEylbDMPg6RkbmbI83hQPD/Jj6vCONKsaWcCe5dDJPfBxL/ekY2dqfC3cPAkcupYtkithF7zb2dwDplJDuGcR+AdZl5d4rXPVBoZhcOTIERITE0s/OREfFxUVRWxsbKGGQetsSErcO7/tNBX6AI/0bKhCXyQPm83Gi/2ak5KZww/rDuXGT2fkMPCj5Uy/txN1K4ed5RHKieTD8Flfz0K/7mVw08cq9EXyqlgPLnsafn32n9iJ7fD7a9BjlHV5ic/6u9CPiYkhJCREczOJFAPDMEhLS+PYsWMAxMXFnXMfnRFJidp6JJn/LdhlirWoFsldXepYlJFI2eaw23jzllakZubw29ZjufGE1Cz6f7iM6fd1plp5XqYyNQE+7weJ5mFB1OgAt30BfoGWpCVS5nUc4e7Of3jtP7HF46BpX4hraVla4nucTmduoV+xYkWr0xHxKcHB7nPAY8eOERMTc84u/ZqgT0qM02XwxDcbyHH9M1LEz27jlRtb4ufQoSdSEH+Hnf/d2ZYOdczj9A8lZdD/w2UcP51pUWYWy0iGL26E41vN8dgWcMc0CFBvIZECOfyg7ztgP+M+j+GE70aAU2Oqpfj8PUY/JCTE4kxEfNPfv1uFmQ9DFZeUmE+X7GXd/kRT7J5udWlaVfM+iJxLkL+DDwddRKvq5nH6e06kMuCjZSSllbOT8+x0mHIbHFpjjlesD/1nQHCUJWmJeJXY5tD1/8yxI+vhz7etyUd8mrrui5SMovxuqdiXErH/ZBqvz9lmitWtFMq/Lm9gUUYi3ic8yJ9Ph7SnQYx5nP7WI6cZ8ulyUjNzCtjTx+RkwbSB5tnEASJrwMDvIKyyNXmJeKOu/weVG5tjC15xT3opIiI+RcW+FDv3jOIbSM92muJjb2hBkH/hl4oQEagQGsDkoR2oGW3uDrk6PpF7Pl9FRp7fM5/jcsKMe2DHL+Z4aGUYMBMitXynSJH4Bbq783PGnaGcdJj1CJx7gSYREfEiKval2H2z+iCLdpwwxe7sUJMOdTVJi8j5qBIRxBdDO1Alwjz53OKdJ3hwyhpynC6LMithhgE/joRN35rjQZHuQr9SfSuyEvF+1S+C9sPMsV2/wYavrclHRERKhIp9KVYJKZm8OGuzKRYbEcSTVzcuYA8RKYwa0SFMvrsDFUL8TfFfNh/l8a/X43L52B05w4BfnoHVn5nj/qFw5zfusccicv4uHwXhVc2x2U9C2klr8hERDwkJCcTExLB3795C79O9e3dGjhxZpOc5n33KCm/OPa/bbruNN954o1gfU8W+FKuXftpCYp6Jw17s15zwIP8C9hCRwmpQJZxJd7UnLNC8auq3aw7y3A+bMHypC+7vr3tOGuYIcC+vV+Nia3IS8SVBEXDNq+ZY2gn4dZQ1+YiUAd27d8dms+V++fv7ExcXxy233MKePUWb12Lw4MHYbDa6d+9+3vm89NJL9O3bl9q1a5/3Y3iDl19+GZvNZlnRXloXDJxOJ6NGjaJOnToEBwdTr149Xnjhhdzzt2eeeYaXXnqJpKSkYntOv3NvIlI4S3ad4NvVB02x3i3iuKJpFYsyEvE9LatH8dGgixj48XIyc/7pvv/Zn/uICPLn0asaWZhdMVn6Hsx/0RyzOeCmT6DeZdbkJOKLmvSBxtfC1h//ia2ZDK1uh9qXWJeX+BSXy+BUWpZlz18hJAC7vWgrAwQEBNCmTRvS0tLYsGED06dPZ/PmzWzcuLGEsvSUlpbGRx99xJw5c0rtOa2wYsUK3n//fVq2bGl1KiXulVde4d1332XSpEk0a9aMlStXMmTIECIjI3nwwQdp3rw59erVY/LkyYwYMaJYnlPFvhSLzBwnz8wwfwCGB/oxuk9TizIS8V0d6lbkvf7tGPbZSnLO6L7/9vydhAf5cU+3ehZmd4HWTIbZT3jG+/0Pmlxb+vmI+LqrX4XdCyHr9D+xH0bCvYvBP8iytMR3nErLot2Lcy17/lXPXEHFsMBzb3iGuLg4li5dCsDQoUP56KOP2LRpEwkJCVSsWBGXy8WECROYOHEiO3fuJDg4mCuuuIJXX32VOnXqULt2bfbt2wfAwoULc5dKmz9/Pt27d2fcuHFMmjSJ+Ph4Tp8+TVRUFF27duXll1+mYcOGAPz0008EBgbSsWNHU26zZ8/mxRdfZOPGjTgcDjp16sT48eOpVy//v/3du3eneXP30LfPP/8cf39/7rvvPp5//nnTEm4ul4vHH3+cDz/8kICAAO69916ee+65837ewkhJSeHOO+9k4sSJvPjii+feoQA5OTk88MAD+b4+l8vFK6+8wgcffMCRI0do2LAho0aN4qabbgLcPTAWLlzIwoULGT9+PAB79uxh69atxf56lyxZQt++fenduzcAtWvXZsqUKSxfvjx3mz59+jB16tRiK/bVjV+KxXsLdrP7RKop9nivRsRE6ERBpCRc1jiGcbe2Ju9Sq2N/3sqXy+KtSepCbZoB3//LM37N69DqttLPR6Q8iKwGPfJ03U/YAYvftCYfkTIkLS2NgwfdvVYrV65MREQEAA888AAjR45k06ZN1K9fH4fDwddff03nzp05duwYbdq0oVKlSgCEh4fToUMHOnTokLv/woUL2blzJ7GxsTRu3JhTp04xY8YMevToQUZGBgCLFi2iXbt2HjmlpqbyyCOPsHLlSubNm4fdbuf666/H5Sp4st5Jkybh5+fH8uXLGT9+PG+++SYffvihxzahoaEsW7aMV199leeff55ff/210M87ZswYwsLCzvoVH28+PxkxYgS9e/fmiiuuKNL/S1Fe39ixY/nss89477332LRpEw8//DD9+/dn4cKFAIwfP55OnToxbNgwDh8+zOHDh6lRo0ah3ueivubOnTszb948tm/fDsC6detYvHgxV199de427du3Z/ny5WRmZl7Qe/I33dmXC7b7eArvzN9pirWqEcUdHWpZlJFI+dCnVVVSM3N48tsNpvi/Z24gNNBB39bVLMrsPGz/Bb4ZBkaek5Ueoz1nDReR4nXxUFj/FRxc9U9s0ZvQ/Eao7ANDg0SKaN++faa73gEBAUyePBl/f3/27NnDe++9B7iLzIEDB5KSkkKTJk04cOAAEyZMYMaMGQwePJhJkybRtm1bFixYYHr8MWPGMH36dPz93XNazZ07l549e3LgwAH++OMPevTowb59+6haNc8kmsCNN95o+vnjjz+mcuXKbN68OfcOfl41atRg3Lhx2Gw2GjVqxIYNGxg3bhzDhv3z97Vly5aMHj0agAYNGvD2228zb948evbsWajnvffee7nlllvO+r6e+XqmTp3K6tWrWbFixVn3KYyCXt/AgQMZM2YMc+fOpVOnTgDUrVuXxYsX8/7779OtWzciIyMJCAggJCSE2NjY3McszPtc1Nf85JNPkpycTOPGjXE4HDidTl566SXuvPNO0/ZZWVkcOXKEWrUuvJZSsS8XxDAMnpm5kawzlv5y2G2Mub45jiKOjxKRorutfU1SMnN4cdaW3JhhwP9NW0dYoB89mnjBnBl7F8O0AeAyT+7JJY9A10esyUmkPLE7oM9/4f1LwXC6Y65sd3f+wbPAro6gUr78PWY/JyeHzZs3k56ezpAhQ1i2bBkrV67MnVBt0KBBDBo0yLTv393/z2bfvn0MHz6c9evXk5KSYppg99ChQwCkp6cTFOTZQ3bHjh08++yzLFu2jBMnTuTeaY6Pjy+w2O/YsaPp4kWnTp144403cDqdOBwOAI8x83FxcRw7dqzQzxsdHU10dPQ5XzvA/v37eeihh/j111/zfY1FVdDr27lzJ2lpabkXLP6WlZVFmzZtzvqYhXmfi/KaAaZNm8YXX3zBl19+SbNmzVi7di0jR46katWqucdRcHAw4O5VUhxU7MsFmbn2IEt2JZhiQzrXplnVSIsyEil/hnatS3J6Nv/97Z8eNjkug/u/WM2nQ9rTqV5FC7M7hwOr4MtbISfDHL94GPR41pqcRMqj2ObQ+QH4Y/w/sfglsOYzaDfYsrTE+1UICWDVMxfWTftCn7+ozhyzv3nzZpo1a8ahQ4d47733aNWqVe52rVu3JjDQPB/Aue7G7t69m379+pGVlUV4eDjt2rUjJyeHtWvXAu4Z2wEqVarEqVOnPPbv06cPtWrVYuLEiVStWhWXy0Xz5s3JyrqwSRD/7mXwt7/Huxf2eceMGcOYMWPO+hybN2+mZs2arFq1imPHjtG2bdvcNqfTye+//87bb79NZmZm7kWIC5GSkgLArFmzqFbN3Nsx7/9bXoV5n4vymgEee+wxnnzySW67zT00sUWLFuzbt4+xY8fmFvsnT7qXP61cuXIRXmnBVOzLeUtMy+LFH7eYYlUjg3i4Z0OLMhIpvx7u2ZDkjBw+XbI3N5aZ42LopBV8OawjrWpEWZZbgY5ugsk3QFaKOd7qdvekYXknJBCRktXtSdg0ExL3/RP79VloeDWEe0EvISmT7HZbkSfIK6syMjJo164dNpsNwzAYPHgwDz30EODu7bp48WIiI903vEJCQgD3WPczrVmzJrdgnDNnDp06dWLq1Kncfvvtpu3atGnD5MmTTbGEhAS2bdvGxIkT6dq1KwCLFy8+Z97Lli0z/bx06VIaNGhQ6IK6MM9blC7tPXr0YMMG8xDEIUOG0LhxY5544okiF/oFvb6mTZsSGBhIfHw83bp1K3D/gICA3IssUPj3uajd+NPS0rDn6SnlcDhMF1U2btxI9erVc+d8uFAq9uW8vTJ7Kwmp5quI/+nbnNBAHVYipc1ms/HstU1Jyczh61UHcuOpWU4GfbKcr4Z3olFsuIUZ5pGwCz7rBxmJ5njja+G6t9VtWMQKASFw7Tj3Rbi/ZSTBnKfhpo+sy0uklB0+fJiOHTvmduMHsNvt9OnTh7p16zJs2DA++OADRo4cyfjx4wkLC2Pfvn0kJyfzySef0LJlSxo3bgzAypUradGiBaGhocyfP59mzZrljtfu1asXNWvW5MiRIx45XHXVVTz11FOcOnWKChUqAFChQgUqVqzIBx98QFxcHPHx8Tz55JPnfD3x8fE88sgj3HPPPaxevZoJEybwxhtvFPr9KMzzFqVLe3h4uMeQg9DQUCpWrFjgUISzKej1hYeH8+ijj/Lwww/jcrm45JJLSEpK4o8//iAiIiL3bnrt2rVZtmwZe/fuJSwsjOjo6EK9z0Xtxt+nTx9eeuklatasSbNmzVizZg1vvvkmd911V+42ixYt4sorryzye1AQnU3JeVmx9yRTlu83xa5sWoWeTXXlX8QqdruNl29oQa9msaZ4Ylo2Az5axr6E1AL2LGVJB9yFfuoxc7ze5XDTx+DQBUMRy9TvAS1uNsc2fg0751mTj4gFsrKyWLZsGatWrcLPz49OnTrx1Vdf5d4dfvfddxk3bhwtWrTg0KFD7Nu3j9q1a/PII4/QvXt3AO666y5uvPFGIiMj2bhxI8uWLcPpdNK4cWM+/vhj6tSpQ1ZWFpUqVWLKlCkeObRo0YK2bdsybdq03Jjdbmfq1KmsWrWK5s2b8/DDD/Paa6+d8/UMHDiQ9PR02rdvz4gRI3jooYcYPnx4od+P833eC/Xpp5+axuIX5Gyv74UXXmDUqFGMHTuWJk2a0KtXL2bNmkWdOnVy93/00UdxOBw0bdqUypUrEx8fXyKvd8KECdx0003cf//9NGnShEcffZR77rmHF154AXD3HJk5c6Zp4sQLZTPOnBGiAMnJyURGRpKUlJS7ZISUX1k5Lq6dsIjtR//pehsa4ODXR7pRNSrYwsxEBCAzx8nQSStZtOOEKV69QjBf39uZ2EgLl8RMOQ6f9IIE8woe1OwE/b+BgFBr8hKRf6Qch7cvMve8qVAb7l8K/vo7L2evDTIyMtizZw916tQplsnXyrNZs2bx2GOPsXHjRo/u34XVvXt3WrduzVtvvVW8yZWC0aNHs3DhQo/VDHzVu+++y4wZM/jll1/Oul1Rfsd0Z1+K7MPFu02FPrjHC6vQFykbAv0cvD+gHe1qVTDFD5xKp/9HyziZemGT+Jy39FPw+fWehX5cK7jjKxX6ImVFWGXo+bw5dmov/F7yd/JE5B+9e/dm+PDhHDx40OpULPHzzz/z6quvWp1GqfH392fChAnF+pgq9qVI9p9M47/zdphizapGMLhzbWsSEpF8hQT48fHgi2kaZ77jsvNYCoM+Xk5yRnYBe5aQzBT44mY4ap6Qh0qNoP+3EKQVPETKlDYDoEZHc+yP8XBsS/7bi0iJGDlyJDVq1LA6DUssX76c9u3bW51GqRk6dCiNGjUq1sdUsS+FZhgGz8zcSEb2PzNG2mww5voW+Dl0KImUNZHB/nx2d3vqVjLfMd9wMImhn64kPctZwJ7FLDsDpt4OB1aY41G1YOBMCC2eGWdFpBjZ7dDnLbCfMYeGKwd+GAlnzBwtImXbggULvLILvxQPVWhSaLM2HGbh9uOm2MCOtcrmkl4iAkClsEAmD+1AtTzDbJbvPcl9X6wiK6eET9pzsmDaANjzuzkeHgcDv4OIqvnvJyLWi2kCXR4yx/YvhTWfW5OPiIgUiYp9KZTkjGz+88NmUywmPJD/u6p4u5qISPGrGhXM5KEdqJRnneMF247z8FdryXGWUMHvzIavh8COPBPNBEfDgJkQXSff3USkDLn0MffkfGf6dRSkHMt3cxERKTtU7EuhvD5nG8dPZ5piz13XjIggf4syEpGiqFMplM/vbk9EkHlZu1kbDvPo9HU4XedcmKVoXE74djhs/dEcD4yAAd9CTOPifT4RKRn+wdD7TXMsIwnm/NuafEREpNBU7Ms5rd2fyOdL95lilzWqzNXNYwvYQ0TKoiZxEXx6V3tCAhym+My1h3jym/W4iqvgd7nguxGw6Vtz3D/Uvbxe1TbF8zwiUjrq94DmN5ljG6bBrt+syUdERApFxb6cVY7TxdPfbsA4owYI8rfzfN/m2Gw26xITkfPStmYFJg68iEA/88f/9FUHeOa7jRjGBRb8hgE/joR1U8xxv2C4cxrUKD+z6or4lKvGeK6a8eMjkJ1uTT4iInJOKvblrD5dspfNh5NNsYd6NKRGdIhFGYnIhepSvxLvD2hHQJ5VNL5cFs9/fth8/gW/YcDPT8DqSea4IxBunwK1LznPjEXEcuFV4Ir/mGOn9sCiN6zJR0REzknFvhToUGI6b/663RRrVCWcoV01qZaIt+veKIb/3dkWP7u5h86nS/Yy5qctRS/4DcM9adfy981xuz/c+jnUu+wCMxYRy7UdBDU6mGOL34JjWy1JR0REzk7FvhToue83kZZnHe4xNzTH36HDRsQXXNG0Cm/f0QZHnoJ/4qI9vP7LtqIV/PPHwJIJ5pjNATd/Ag2vKoZsRcRydjtcOw7sZ0z06cp2D91xlfAyniIiUmSq2iRfv24+yi+bj5pit7evSbta0RZlJCIloVfzOMbd2po89T7vzN/Ff+ftLNyD/P4a/P6qOWazw40ToUmf4klURMqGKs2g87/Msfg/Ye1ka/IREZECqdgXD2lZOTz3/SZTrGJoAE/20lJZIr7oulZVef3mVuSdc3Pc3O38b8E5Cv4lE+C3F/MEbdDvXWh+Y7HmKSJlxKWPQ1Qtc+yXUZBy3Jp8REQkXyr2xcP4uTs4mGieXfeZa5sQGeJvUUYiUtJuaFudl29o4RF/dfY2Ply0O/+d/vwf/PKMZ7zPW9DqtuJNUETKjoAQ6P2mOZaRmP/ngYict4SEBGJiYti7d2+h9+nevTsjR44s0vOczz5lhTfnntdtt93GG28U76SnKvbFZOuRZD5cvMcU61yvIv1aV7MoIxEpLbdeXJMX+jX3iL84awuTluw1B//8H8x5yvNBrn4N2g0ukfxEpAxpcAU0u8EcWz8Vdi+wJB2R4tC9e3dsNlvul7+/P3Fxcdxyyy3s2bPn3A9whsGDB2Oz2ejevft55/PSSy/Rt29fateufd6PUVY5nU5GjRpFnTp1CA4Opl69erzwwgsXvgTwebDigsHLL7+MzWYzPe8zzzzDSy+9RFJSUrE9j9+5N5HywuUyePrbDThd//ySBTjsvNCvOba8/XtFxCcN6FiL7BwXz/+42RQf/f0mbDYY2Kk2/PkOzHnac+eeL0CH4aWTqIhYr9dY2DkPMs84Mf3xEbhvCfgHWZeXlC0uF6SftO75g6Pdk0sWQUBAAG3atCEtLY0NGzYwffp0Nm/ezMaNG0soSU9paWl89NFHzJkzp9SeszS98sorvPvuu0yaNIlmzZqxcuVKhgwZQmRkJA8++KDV6ZWoFStW8P7779OyZUtTvHnz5tSrV4/JkyczYsSIYnku3dmXXF+t3M/q+ERT7N7u9ahXOcyahETEEnddUoenrvaco+PZ7zax/Mvn8y/0r3gOuvj2H2cRySM8Fq4YbY6d3AWL38x/eymf0k/Ca/Ws+zqPCw1xcXEsXbqU9evXc/fddwOwadMmEhISAHC5XIwfP57mzZsTFBREhQoVuPnmm3Pv/teuXZtJkyYBsHDhwtyeAgsWLABg3LhxtG7dmujoaPz9/alcuTI33HAD27f/s+T1Tz/9RGBgIB07djTlNnv2bC655BKioqKoWLEi1157Lbt27SrwtXTv3p0HHniABx54gMjISCpVqsSoUaM87qC7XC4ef/xxoqOjiY2N5bnnnrug5z2XJUuW0LdvX3r37k3t2rW56aabuPLKK1m+fHmRHysnJ6fA1+dyuRg7dmxuD4JWrVrx9ddf5+47ePBgFi5cyPjx43P/n/bu3Vvsr/dvKSkp3HnnnUycOJEKFSp4tPfp04epU6de8PP8TcW+AHAiJZOXfzavk1urYgj3d69nUUYiYqV7utXj0SsbmmJ3O2bRfns+Y8mu+A9c8nApZSYiZUq7IVD9YnNs0ZtwfHv+24t4kbS0NA4ePAhA5cqViYiIAOCBBx5g5MiRbNq0ifr16+NwOPj666/p3Lkzx44do02bNlSqVAmA8PBwOnToQIcOHXL3X7hwITt37iQ2NpbGjRtz6tQpZsyYQY8ePcjIyABg0aJFtGvXziOn1NRUHnnkEVauXMm8efOw2+1cf/31uM6y/OWkSZPw8/Nj+fLljB8/njfffJMPP/zQY5vQ0FCWLVvGq6++yvPPP8+vv/5a6OcdM2YMYWFhZ/2Kj4/PfbzOnTszb9683Asc69atY/HixVx99dVF+086x+sbO3Ysn332Ge+99x6bNm3i4Ycfpn///ixcuBCA8ePH06lTJ4YNG8bhw4c5fPgwNWrUKNT7XNTXDDBixAh69+7NFVdcke9rad++PcuXLyczM7PI70N+1I1fABgzawtJ6dmm2Iv9mhPk77AoIxGx2gOXN8BlwJu/buduxyxG+X/hudEV/4FLRpZ6biJSRtjtcO1b8P6lYDjdMVc2/DgSBs/CY5kPES+wb98+0xDWgIAAJk+ejL+/P3v27OG9994D3EXmwIEDSUlJoUmTJhw4cIAJEyYwY8YMBg8ezKRJk2jbtm3uHf2/jRkzhunTp+Pv7578eu7cufTs2ZMDBw7wxx9/0KNHD/bt20fVqlU9crvxRvNKNx9//DGVK1dm8+bNNG/uOe8OQI0aNRg3bhw2m41GjRqxYcMGxo0bx7Bhw3K3admyJaNHu3vqNGjQgLfffpt58+bRs2fPQj3vvffeyy233HLW9/XM1/Pkk0+SnJxM48aNcTgcOJ1OXnrpJe68886zPkZRXt/AgQMZM2YMc+fOpVOnTgDUrVuXxYsX8/7779OtWzciIyMJCAggJCSE2NjY3McszPtc1Nc8depUVq9ezYoVK866fVZWFkeOHKFWrVoFbldYKvaFJbtO8O2ag6bYda2q0rVBZYsyEpGy4sEeDWh7cDKX7PYs9JfVH0kHFfoiEtscOj8Af4z/J7bvD1j7BbTpb11eIufp7zH7OTk5bN68mfT0dIYMGcKyZctYuXJlbhfxQYMGMWjQINO+S5cuPefj79u3j+HDh7N+/XpSUlJMXeoPHToEQHp6OkFBnnNf7Nixg2effZZly5Zx4sSJ3DvN8fHxBRb7HTt2NF286NSpE2+88QZOpxOHw31jL+/48bi4OI4dO1bo542OjiY6Ovqcr/1v06ZN44svvuDLL7+kWbNmrF27lpEjR1K1alWP9/RcCnp9O3fuJC0tLfeCxd+ysrJo06bNWR+zMO9zUV7z/v37eeihh/j111/z/X/9W3BwMODuVVIcVOyXc5k5Tp6ZaZ5sJDzIj2eubWJRRiJSpiyZwCW73/IIj8m+nQ82tuephbu4p5uG+4iUe92egE0zIPGMLqu/PAMNe0FoJevyEusFR8NjFz7W+YKev4j+HrMPsHnzZpo1a8ahQ4d47733aNWqVe52rVu3JjAw0LTvue7G7t69m379+pGVlUV4eDjt2rUjJyeHtWvXAu5Z6gEqVarEqVOnPPbv06cPtWrVYuLEiVStWhWXy0Xz5s3Jysoq8us809+9DP5ms9lMXdbP9bxjxoxhzJgxZ32OzZs3U7NmTQAee+wxnnzySW67zb1Ub4sWLdi3bx9jx44tcrFfkJSUFABmzZpFtWrmlcXy/r/lVZj3uSivedWqVRw7doy2bdvmtjmdTn7//XfefvttMjMzcTgcnDzpnmOicuXiuemqYr+ce3/hbnYfTzXFHr+qETHhmkVXpNxbMiHfdbNfyr6Dic5rARj781ZcBtyn+T1EyreAUOj9Jnxx0z+x9FPwyyi4/l3r8hLr2e0+c8EnIyODdu3aYbPZMAyDwYMH89BDDwFgGAaLFy8mMjISgJCQEMA91v1Ma9asyS0Y58yZQ6dOnZg6dSq33367abs2bdowefJkUywhIYFt27YxceJEunbtCsDixYvPmfeyZctMPy9dupQGDRrk3tU/l8I8b1G7tKelpWHPs0qCw+E469wDBSno9TVt2pTAwEDi4+Pp1q1bgfsHBATkXmSBwr/PRXnNPXr0YMOGDaa2IUOG0LhxY5544onc/4uNGzdSvXr13DkfLpSK/XJs74lU3p6/0xRrVT2SOzpc+PgQEfFyi96Aec97hP+sN5KJm9qbYq/M3orLMBhxWf3Syk5EyqIGPaFpP9g885/Yui+h9e1Q51KrshIpssOHD9OxY8fcbvwAdrudPn36ULduXYYNG8YHH3zAyJEjGT9+PGFhYezbt4/k5GQ++eQTWrZsSePG7lVtVq5cSYsWLQgNDWX+/Pk0a9Ysd4x6r169qFmzJkeOHPHI4aqrruKpp57i1KlTubO2V6hQgYoVK/LBBx8QFxdHfHw8Tz755DlfT3x8PI888gj33HMPq1evZsKECbzxRj4T7hagMM9b1G78ffr04aWXXqJmzZo0a9aMNWvW8Oabb3LXXXcV+jH+VtDrCw8P59FHH+Xhhx/G5XJxySWXkJSUxB9//EFERERuD4LatWuzbNky9u7dS1hYGNHR0YV6n4vymsPDwz2GWYSGhlKxYkVTfNGiRVx55ZVFfg8Kotn4yynDMBj13Uaycv65ema3wUvXt8Bh12Q6IuWWYcBvL+Zb6HPli3Qa8B+e6e05zOe1Odv477wdHkv5iEg50+tlCIwwx358GLIzrMlH5DxkZWWxbNkyVq1ahZ+fH506deKrr77KvTv87rvvMm7cOFq0aMGhQ4fYt28ftWvX5pFHHqF79+4A3HXXXdx4441ERkayceNGli1bhtPppHHjxnz88cfUqVOHrKwsKlWqxJQpUzxyaNGiBW3btmXatGm5MbvdztSpU1m1ahXNmzfn4Ycf5rXXXjvn6xk4cCDp6em0b9+eESNG8NBDDzF8+PBCvx/n+7xnM2HCBG666Sbuv/9+mjRpwqOPPso999zDCy+8kLvNp59+ahqLX5Czvb4XXniBUaNGMXbsWJo0aUKvXr2YNWsWderUyd3/0UcfxeFw0LRpUypXrkx8fHyxv97CyMjIYObMmaaJEy+UzSjEmVlycjKRkZEkJSXlLhkh3u37dYd4cMoaU+yuLnV4tk9TizISEcsZhrvb/p9ve7Zd+SJ0/lfujx8t3sMLP2722OzebvV4olejQv1xFhEftXwi/PSoOdbtSbjsKWvykWJ3ttogIyODPXv2UKdOnbNORCbnNmvWLB577DE2btzo0eW9sLp3707r1q156623ije5UjB69GgWLlzosZqBr3r33XeZMWMGv/zyy1m3K8rvmO7sl0NJ6dkeJ+mxEUE8kmdNbREpR1wu98l5foX+1a+aCn2Auy+pw+h8Lg6+t3AXz32/CZdLd/hFyq2L7oJqF5lji9+EEzusyUfES/Xu3Zvhw4dz8ODBc2/sg37++WdeffVVq9MoNf7+/kyYMKFYH1Nj9suh1+ds4/jpTFNsdJ+mhAXqcBApl1xO+P5BWDs5T4MN+rwF7Qbnu9uQLnXws9sY9d0mU3zSn/tIzXLyyo0tNSxIpDyyO9yfHe93A+OvSa+cWe7u/IN+APX8ESm0kSNHWp2CZZYvX251CqVq6NChxf6Yqu7KmbX7E5m8bJ8pdnnjGHo1j7UoIxGxlDMbZtwDG78xx2126PcetLr1rLsP6FSbIH8HT3yznjNv5n+96gDp2U7eurU1/g51IhMpd2JbQKf73at6/G3vIlg3BVrfYV1eIuVMeekCL/nTGVg5kuN08e8ZGzhzloYgfzv/ua6ZxteKlEc5mTB9sGehb/eDmz45Z6H/t5svqsF/b2+DX567+LPWH+a+yavIyHYWsKeI+LTuT0FkDXNszr8hNcGafEREyhkV++XIp0v2sulQsin2YI8G1IgOsSgjEbFMVhpMvQO2/miOOwLg1i+gWb8iPdy1Lavy/oB2BPiZ/6zM3XKMoZNWkpaVc4EJi4jXCQiFa143x9JPwq/PWpOPiEg5o2K/nNh/Mo03ftluijWsEsawrnUtykhELJN+Cj7vBzvnmuN+wXDHV9Co13k9bI8mVfh40MUE+ztM8cU7TzDwo+UkZ2SfZ8Ii4rUa9YIm15ljayfDnkXW5CMiUo6o2C8HDMPgmZkbSc/Tlfal61toLK1IeXP6CHzSG/YvM8cDwqD/N1Dv8gt6+EsaVOKzu9sTnmfCz5X7TnH7B0s9JgcVkXLg6lcgINwc+/Fh91Ai8VmFWN1bRM5DUX63VOmVA9+vO8TC7cdNsTs71OTi2tEWZSQilji5Bz6+Co6ZZ88nKAoGfge1uxTL01xcO5ovhnUgKsTfFN90KJmb31vC/pNpxfI8IuIlIqpCjzxd9xN2wOK3LElHSpa/v/uzPy1Nn/UiJeHv362/f9fOxmYU4tJAcnIykZGRJCUlERERceEZSqlJTMuixxsLSUjNyo3FhAfy6yPdiAw+9wEiIj7iyEaYfAOkHDXHw+NgwAyIaVLsT7ntyGnu/HAZJ1LMd+8qhwfy2V3taRKnvyci5YbLCR9eAYdW/xNzBMB9f0Kl+tblJeflXLXB4cOHSUxMJCYmhpCQEE0ELVIMDMMgLS2NY8eOERUVRVxc3Dn3UbHv4x6bvo7pqw6YYu/e2ZarW5z74BARHxG/FL68BTKSzPHoujBgJlSoVWJPvfdEKgM+Xsb+k+mmeHiQHx8Nupj2ddTDSKTcOLwOPugOhuufWJ1LYeD3oGLQq5yrNjAMgyNHjpCYmFj6yYn4uKioKGJjYwt1EU3Fvg9bsvMEd3xoHpfbs2kVPhjQTldYRcqLHb/CVwMgx1xsE9sC+n8LYTElnsKx5AwGfrycrUdOm+KBfnbeuaMtVzStUuI5iEgZMeff8Ofb5tj170Or26zJR85LYWsDp9NJdrYmZxUpLv7+/jgcjnNv+BcV+z4qI9tJr7d+Z2/CP+OlwgL9+PWRS4mLDLYwMxEpNeunwcz7wJVn2buaneGOqRAUWWqpJKVnM2zSSpbvPWmKO+w2Xr6hBTdfVKOAPUXEp2SmwDsdIPmMXochFeGBlRCinj7eQrWBiHfQBH0+asJvO0yFPsDjvRqp0BcpDwwDFo+Db4d5FvoNrnLPul+KhT5AZLA/n93dniuamHsSOF0Gj329nvcX7irVfETEIoFh0Pt1cywtAX4dZU0+IiI+TMW+D9p6JJn3F+42xdrUjOLODiU3LldEygiXE356FOY+59nW8la47QsICCn1tACC/B28178dN7Wr7tE29uetvPDjZlwuLdUk4vMaXQ2NrzXH1kyGvX9Yk4+IiI9Sse9jnC6DJ7/ZQM4ZJ8x+dhtjb2iBw65x+iI+LTsdpg2EFR96tnW4D/q9Bw5rV+Hwc9h57aaW3NOtrkfbR4v3MOLL1WRkOy3ITERK1dWvQkCYOfbjSMjJzHdzEREpOhX7Pmby0n2s3Z9oit3brR6NYzWeSsSnpSbApOtg64+ebVe+BL3Ggr1sfOTbbDaeuroJT1/T2KPt541HuPPDZZw8Y7lQEfFBkdXg8jxd909shz/+a00+IiI+qGyc+UmxOJSYzmtztplidSqF8sDlWr9WxKed3AMf9YQDy81xRwDc9DF0fqBMLms1/NJ6vHlLK/zy9Dpate8UN767hH0JqRZlJiKlov0wiGttjv3+GiRoDg8RkeKgYt9HGIbB0zM2kJJpnozrpeubE+Rf+OUZRMTLHFzlLvRP5jk5DoqEATOg+Y3W5FVIN7StzqS72hMe6GeK7zmRyg3/W8Ka+FMWZSYiJc7ugD7jwXbG6agzE2Y94p5oVERELoiKfR/xzeqDLNh23BS7uV11OterZFFGIlLiNs2ET3pDqvl3n4jqcNccqH2JJWkVVZf6lZh+XyfiIoNM8YTULG6fuJTZG49YlJmIlLiqraHDvebY7gXupUNFROSCqNj3AceSM3j+h02mWEx4IM/0bmpRRiJSogzD3dV1+iDISTe3VWkOQ+dCTBNrcjtPjWMjmHF/F5rEmecXych2cd8Xq3h3wS4M3ekT8U2XPQ0R1cyx2U9CyjFr8hER8REq9r2cYRj8e+ZGkjPydt9vQWSItbNui0gJyMmEGffAby96ttXtDkN+hoi4Uk+rOMRGBjHtno50bWDukWQY8Mrsrfzf9HVk5mimfhGfExgO17xmjqWfVHd+EZELpGLfy/2w/jC/bj5qil3Xqio9m1axKCMRKTGpJ9wz7q//yrOt7SC482sI8u6VN8KD/Pl48MXcclF1j7ZvVx/kzonLOJGipblEfE7j3tDkOnNsyw+w6Vtr8hER8QEq9r1YQkomz31v7r5fMTSA565rZlFGIlJijm2BiZfD/qV5Gmxw1Rj3JFcO3+jN4++w88qNLXny6sYeiwis3HeKvm//wdYjydYkJyIlp/ebEBxtjs16FFKO57+9iIiclYp9L/bs95s81qJ+vm9zokMDLMpIRErEjrnw0ZWQuM8c9w+F26dApxFlcmm9C2Gz2bi3Wz3e79+OkADziiIHE9O58X9LmLflaAF7i4hXCqsMvV83x9SdX0TkvKnY91KzNx5m1vrDplivZrFc0yLWooxEpNgZBix6E764CTLz3MmOqA53z4FGV1uTWym5slksX9/bmap5ZupPzXIy9LOVfPC7Ju4T8SnNbsinO//3sGmGNfmIiHgxFfteKDEti2dmmrvvR4X483y/Zth87O6eSLmVmQLTB8O8/wB5itlq7WDYbxDbworMSl3TqhHMfKALbWpGmeKGAWN+2sqDU9eSlpWT/84i4l1sNuj9hmd3/p/UnV9EpKhU7HuhZ7/b5DFB1eg+TYkJDypgDxHxKid3w0c9YfNMz7Zm18PgWRBevibhjAkPYsqwjvRrXdWj7Yd1h7jhf0vYl5BqQWYiUuzCYjy786clwE//Z00+IiJeSsW+l/lh3SG+X3fIFOvROIZ+rasVsIeIeJUdc+GD7nBsc54GG/R4Fm76BPyDrcjMckH+Dsbd2ppHr2zo0bb1yGn6TFjMgm1al1vEJzS7AZr0Mcc2f6fu/CIiRaBi34scTc5g1HcbTbGIID9eur6Fuu+LeLszx+dnJJnbgiLhzunQ9f98biK+orLZbDxweQM+GnQR4YF+prbkjByGfLqCd+bv1Dh+EW9nsxUwO///qTu/iEghqdj3EoZh8MQ360lMyzbFX+jXnNhIdd8X8WoZSTBtYP7j8ys3gWHzoUFPS1Irq3o0qcL3/7qEBjFhprhhwGtztnHv5FWczsguYG8R8QphMXDNa+ZYWgL8OFKz84uIFIKKfS8xZfl+FmwzX8nu3TKO61p5jl8VES9yeB28380923ReTfvC0LlQsV7p5+UF6lQKZeaILvmuQjJn01H6vvMHW48k57OniHiN5jd6duff+iOsm2JNPiIiXkTFvheIT0jjxVnm8buVwwN5sW9zdd8X8VaGASs/gQ97wqk9eRr/Gp9/8yQIDMt3d3ELDfTjnTva8kSvxtjzfBzuPp5K37f/4KsV8erWL+Kt/u7OH1LRHP/pcTi115KURES8hYr9Ms7pMvi/6WtJy3Ka4q/e2JIKoQEWZSUiFyQzBb4d7u6K6jSvrEFQlMbnF5HNZuO+7vWYdFd7okL8TW2ZOS6e+GYDD3+1ltRMLc8n4pXCYqDPeHMs6zTMuBdczvz3ERERFftl3YeLdrNi7ylT7Pb2NbmscYxFGYnIBTm2BSZeBhumebZVuwjuXaTx+eepa4PK/PDAJTSvFuHRNnPtIfpMWMyWw+rWL+KVmvSB1v3Nsfg/Ycl/rclHRMQLqNgvw7YcTuaNX7abYjWjQ3imdxOLMhKR82YYsPpz+OAyOLHds73j/TDkZ4iqWfq5+ZAa0SF8c19nBneu7dG2+0Qq/d75gy+XqVu/iFe6+mWIqmWO/faSe+4TERHxoGK/jErPcvKvKWvIcrpyYzYbvHFLK0LzLDclImVc2kmYPgi+fwBy0s1tgRFwy+fQayz4aWhOcQj0c/Dcdc149862HsvzZea4eHrGBh6cupZkzdYv4l0Cw+H698F2xumrK9s9LCo7w7q8RETKKBX7ZdQLszaz81iKKTa8a10urh1dwB4iUibtWQTvXQKbv/Nsi2sF9yyEpteVfl7lwNUt4pj1YFdaVo/0aPth3SGufmsRy/ectCAzETlvtTpBl5Hm2PGt8OsoS9IRESnLVOyXQbM3HubLZfGmWLOqETxyZUOLMhKRIsvJgrnPwaQ+kHzQs/2iu+GuXyC6bqmnVp7UrBjC9Hs7MaRLbY+2g4np3PbBn7w2ZytZOS7PnUWkbOr+FMS2NMeWfwBbZ1mTj4hIGaViv4w5lJjOE99sMMWC/R389/Y2BPo5LMpKRIokYRd8fCUsHgfkGRseXAFu/QKufRP8gyxJr7wJ9HMwuk8z3uvfjoggc7d+lwHvzN/Fje8uYdfxlAIeQUTKFL8AuGEi+AWb49+NgKR8Lq6KiJRTKvbLEKfLYORXa0lKN48j/c91zahXWWtti5R5Lhcsn+jutn9ojWd7nW5w35/Q5NrSz03o1TyWn0deSoc6nsOhNhxM4tr/Lmby0n2avE/EG8Q0ds91cqb0U/DNUHBqmU0REVCxX6a8M3+nx/jR3i3juPmi6hZlJCKFdmoffHYd/PQoZKeZ2+z+0PMFGDATIuIsSU/cqkUF8+WwjjzRqzH+DpupLT3byTMzNzLgo+UcOJVWwCOISJnRbjA07WeOxS+B31+zIhsRkTJHxX4Z8eeuBN6aa16Oq1pUMGOub4HNZitgLxGxnGHAyo/h3c6wd5Fne8UGMHQudHkQ7PrILQscdhv3da/HjPu7UK9yqEf74p0n6PXWIi3RJ1LW2WzQZzxE5lmy9PdXYe9ia3ISESlDdOZZBhw7ncG/pqzBdcY5pd0G429rTWSwv3WJicjZJe6Hz/vBjw9DVj7jvS+6yz3bftXWpZ2ZFELzapH8+K+uDOhYy6MtJTOHp2ds0F1+kbIuOApu+ghsZ8xrZLjgm2GQesKytEREygIV+xbLcbp4cMoaTqRkmuIPX9GQi7TMnkjZ5HLBqk/hf51g9wLP9sga7i77146DAM87x1J2BAc4eKFfcz4dcjFxkZ4TJi7eeYKrxv3OF8s0ll+kzKrRHi7/tzl2+pB7/L7LaU1OIiJlgIp9i42bu52lu83j9C9tWJkRl9W3KCMROavj2+DT3vDDQ5B12rO97SC4bwnUu6z0c5Pz1r1RDHMevpRbL6rh0Zaa5eTfMzZy6wdL2XlMM/aLlEldHoa63c2x3fNhwdh8NxcRKQ9sRiFuVSQnJxMZGUlSUhIRERGlkVe58NvWo9z16UpTLC4yiFkPdiU6NMCirEQkX9kZsPhNWPQmuLI92yOqwXX/hfpXlH5uUqwWbDvGU99u4HBShkebv8PGfd3rc3/3egT5azlUkTIl5Ri81xVSjpjjd0yDhldZk5OPUm0g4h10Z98iB06l8fBX60wxP7uNt+9oq0JfpKzZu9i9nN7CV/Iv9Nv0h/v/VKHvI852lz/bafDfeTu4evwiluzUeGCRMiUsBm6ZBHY/c/zbYXByjzU5iYhYSMW+BTKyndz/xWqS0s1Fw9PXNKFdrQoWZSUiHlKOwcz73d32E3Z4tkfVgv7fQN93ICiy9POTEhMR5M8rN7Vk0l3tqRYV7NG+50Qqd3y4jEe+WktCnjlXRMRCNTu6lzo9U0YSTBsI2enW5CQiYhEV+6XMMAye/nYD6w8kmeJXN49lSJfa1iQlImbObFj6LkxoB2u/8Gy3OaDLQ3D/Ut3N93HdGlbm10cu5Z5udXHYPZdB/XbNQS57fQGTluwlx+myIEMR8dDxPmh2vTl2ZD3MetS9XKqISDmhMful7KPFe3jhx82mWJ1KoXz3QBcigrTMnojl9vwOPz0Ox7fk316tnXtd59gWpZuXWG7L4WSe+nYDa/cn5tveqEo4o/s0pXP9SqWbmIh4yjwNEy+HE9vN8V4vuy8GyAVRbSDiHVTsl6I/dp5g4MfLcbr+ectDAxzMHNGFBlXCLcxMREjcD7+Ogk0z8m8PCIMez8LFQ8GuidnKK5fL4Mvl8bwyeyunM3Ly3eaaFrE8fU0TqlcIKeXsRMTk2FZ3wZ+d+k/MZoc7p6tX1gVSbSDiHVTsl5J9Can0e+cPTqWZx+m/P6AdVzWLtSgrESEjGRaPg6X/gxzP2dcBaHkb9PwPhOt3VdyOJWfw4qwtfL/uUL7tgX52hl9al3u61SMs0C/fbUSkFGyaCdMHmWOBkTB0LlRuaElKvkC1gYh3ULFfChLTsrjh3SXsPp5qij/UowEP99QfGhFLOLNh1aew4GVIK2BW9dgWcM3r7gmfRPKxfM9JRn+/iS2Hk/NtrxQWwEM9GnBb+5r4OzRNjoglFrwCC8aYY9H1YNg8CNbEyOdDtYGId1CxX8KyclwM+GgZy/acNMV7Nq3C+/3bYc9nwicRKUGGAdt+hl+fzX+GfYCgKOgxCtoNUZd9OSeny2Dqinhen7PNo/fW3+pWCuXxXo25qlkVbDZ97ouUKsOAr4d4DtOq2x3u/AYc6n1TVKoNRLyDiv0SZBgG/zdtHd+uOWiKN44NZ/q9nQjXhHwipWvPIvjtRdi/NP92mwMuGgLdn4bQiqWbm3i9xLQsxv26nc+X7sNVwF/WdrUq8ESvxrSvE126yYmUd1lp8EkvOLzOHL/oLuj9JugiXJGoNhDxDir2S9D4uTsYN9c8C2xMeCAzR3Shaj7rNotICdm/An57AfYsLHibRr3d4/IrNSi9vMQn7Th6mldmb2XulmMFbnNJ/Uo83LMB7Wqp6BcpNUkHYeJlkHLUHO/xLHT9P2ty8lKqDUS8g4r9EjJzzUFGfrXWFAv2dzD93k40rxZpTVIi5c3h9TD/Jdg+u+BtqraBK1+C2l1KLy8pF5buTmDsT1tYdyCpwG0ubViZh69oQJuaGjcsUioOrIRPrgFnpjl+/fvQ6jZrcvJCqg1EvIOK/RKwcPtxhk5aQbbzn7fWboMPBlzEFU2rWJiZSDlxcBX8/gZsm1XwNlG13Hdzmt0Adk2cJiXDMAxmbTjMq7O3EX8yrcDtLm8cw8NXNKRFdV0MFilxm2bA9CHAGafAdj/3knz1LrcsLW+i2kDEO6jYL2ar9p2i/4fLSM92muLP9WnK4C51LMpKpBwwDNi7GBa9AbvnF7xdeFXo9hi07g9+AaWXn5RrWTkupiyP5535Ozl2OrPA7S5rVJn7utfn4toVNJGfSEn6838w5ylzLCAMBn0P1dpZk5MXUW0g4h1U7BejbUdOc8v7f5KUbp6NeXDn2jx3XTOLshLxcYYBO36FRa/D/mUFbxdSyT0m86Ih4K85M8QaGdlOvlgWz7sLdnEipeCiv12tCtzXrR6XN47Rqi0iJWXOv+HPt82x4Aow+Ceo0tSanLyEagMR76Biv5jEJ6Rx03tLPO7YXNeqKm/d2lonayLFLTsDNkyHpe/CsU0FbxcUBZ3/BR3uhcCwUktP5GzSs5xMXrqPdxfu4mRqVoHbNYgJ495u9biudVX8HRpuIlKsXC745m7Y9K05HhoDd82GivWsycsLqDYQ8Q4q9ovB/pNp3PbBUg4mppvi3RtVZuLAi3SCJlKcUo7Big9hxUeQdqLg7cKqQKcH3HfyA8NLLz+RIkjNzOGzP/fxwe+7OJWWXeB2VSICGdCxFre3r0nFsMBSzFDEx+VkwpTbYdc8czyiurvgj6phTV5lnGoDEe+gYv8CHTiVxq3vexb6F9WqwOd3dyA4wGFRZiI+5shG9138DdPAWfCdUCJrQJeHoM0A8A8qvfxELkBaVg5frdjPxN93cygpo8DtAvzs9G1VlcFdatOsqibzEykWWWkw+UaIX2KOR9WCQT9AhVrW5FWGqTYQ8Q4q9i/AwcR0bn3/Tw6cMhf6TeMimDK8I5HB/hZlJuIjstNh83ew6lOI//Ps21ZsAJc8DC1vAYd+98Q7ZTtdfL/2EO//vovtR1POum37OtH071iLq5pVIdBPF5ZFLkhGEky6Dg6vNccja7gL/mhNsnwm1QYi3kHF/nnadTyFAR8u87gD0zg2nCnDOlIhVLN8i5y349vcBf7aLyEj8ezb1ukGnUZA/Z5aQk98hstlMH/bMd5buIsVe0+dddsKIf7c0LY6t7evQf0YDVkROW9pJ+GTa+D4FnM8opq74NcY/lyqDUS8g4r987DxYBKDPl5OQp5JlRrHhvPlsI5Eq9AXKbrsdNjyA6z8xLMrZV6OAGhxC3S8D2Kbl05+IhbZcCCJT5fs5Yd1h8hyus667UW1KnBb+5r0bhGnYWQi5yPlOHzW13Pi17BYGPAtVNHqSqDaQMRbqNgvomW7Exg6aSWnM3NM8UZVwvlyWAdNnCRSFC4X7FsM675yd9fPOn327UMrw0V3w8V3Q1hM6eQoUkYcP53JlOXxfL50H8dPF7xsH0B4kB/Xtoyjb+tqtK8drRVhRIoiNQE+7wtHNpjjgRFw25dQp6s1eZUhqg1EvIOK/SL4ZdMR/jVlDZk55jsrrWpE8engi9V1X6Swjm2F9VNh/XRIPnDu7et2h3aDoVFv8NPvmZRvWTkuft54mCnL41m6++Q5t68aGUSf1lW5vk01GseW37/hIkWSdhI+v95zDL8jAG74AJpdb0laZYVqAxHvoGK/EAzDYOKi3Yz9eSt5363O9SrywcCLCAv0syY5EW+RsMt9937zTDi87tzbh1SCNv2h7UCNkxQpwJ4TqUxdEc83qw5wIuUsq1T8pXFsOH1bV+PalnHUiA4phQxFvFh6Inx5K+xfmqfBBr1eho73WpFVmVDeawMRb6Fi/xyynS5GzdzI1BX7PdqubFqF/97ehiB/jYsUydfxbX8V+N/B0Y2F2MEGdbtB20HQ+FrdxRcppKwcF/O2HGXKiv0s2nHc48J0fppVjaBXs1h6NY+lfkwYNpu6+ot4yE6Hb4bC1h892y662130l8O/VeW5NhDxJir2z+JESib/+nINf+5O8Gi7uV11xt7QAj+HZv8WyeVywsHVsGOOe7K941sLt19MM2h1K7S4GSKqlmyOIj7uYGI63609yMw1B8+5fN/f6lYK5armsfRqFkvL6pEq/EXO5HLCT4/Byo8822p1gVs+g9BKpZ+XhcprbSDibVTsF2DF3pM88OVqjiZ7ToL0eK9G3Netnk6GRADST8HOebDjF9g5F9I8L47lKywWWt4MLW/TjPoiJcAwDLYcPs13aw/y/bpDHM6zVGxBYsID6dawMt0aVaZr/cpEhviXcKYiXsAwYNHr8NuLnm2RNeCmT6DGxaWfl0XKY20g4o1U7OdhGAYfLtrDy7O34nSZ35ogfzvjbmnN1S3iLMpOpAxwOd0zFO/6zV3g718GxtmXA8sVGgNNroWmfaF2V7BrCIxIaXC5DJbtOcl3aw8yZ9MRTqVlF2o/h91GmxpRdG9Ume6NYmgaF6GZ/aV8Wz8dvn8AcvJcPLP7weWjoPODYPf9Xp/lqTYQ8WYq9s9wNDmDJ75Zz4Jtxz3aqkQEMnHgRbSsHlX6iYlYyTDcY+/3/A57FsLexZCRWPj9w+OgyXXuAr9mRxX4IhbLcbpYsfcUczYdYfbGIxxJLtwdf4Do0AA61ImmU72KdKpbUWP9pXw6tAam3gnJBz3b6vWA69/z+eVhy0ttIOLtVOzjvpv//bpDPPvdJpLSPe92dG1QibdubU3FsEALshMpZS4nHNvivmO/b4m7yE89VrTHqNQIGvR0F/nVLy4XdzlEvJHLZbD+YFJu4b/nRGqR9q8UFkCHuhXpWLciHepEU79ymO78S/lw+ihMHwzxSzzbgqPh6lehxU3goxfDfL02EPEV5b7YP5SYzvM/bGb2piMebTYbPHh5Ax7s0QCHTl7EV2UkwYGVsH+5u8A/sBKyThftMfyCoM6l0OBKd5FfoXaJpCoiJccwDPYmpLFg2zEWbDvO0t0JZOYUcojOX8KD/GhdI4q2NSvQpmYUbWpU0Jh/8V3OHFj4Cvz+GpDP6XTDXnDtOJ+ceNaXawMRX1Jui/3MHCcfLd7DhHk7Sc92erTHhAfy2s2t6NawsgXZiZSQjGT3ePvDa91r3R9e5+6in99JyrlE14N6l0GDq6D2JRCgNbtFfElGtpM/dyewcNtxFmw7xt6EtPN6nPoxYbSqHkXzahE0qxpJk7hwwoN0AUB8yO6F8O0wSDnq2RYYAZc+Bh3uAT/f6SHqi7WBiC8qd8W+YRj8svkor/y8ld0FdFfs27oq/7muGVEh5W/dVPERhgGnj8DxLX8V9+vg0Fo4uev8HzOiOtTt5r6DX7srRFYrtnRFpOzbfzKNpbsT+HN3Akt3JXCokLP756d2xRCaVYukWVX3BYCmcRFUCgvQ+H/xXinH4af/g83f5d9eoQ5c+QI0vtYnuvb7Um0g4svKTbFvGAYLth9n3K/bWX8gKd9tKoT489L1LbhGs+2Lt8gt6re6v45t+effGfkf54VWoTbU6AA1O7mL/Ap1fOIERUQunGEY7D+Zzp+7T7B090mW7b6w4h/cf4Prx4RRPyacBjFh1I8Jo0GVMGIjgnQRQLzH5u9h1v8VPNdN9Yvh0sfdQ968+Lj2hdpApDzw+WI/2+li9sYjfLR4D2v3J+a7jc0Gt7evyWNXNqJCqO7mSxljGJByDE7tgZO7Pb8utKgHcARAXGuo0d5d4NdoD+GxF/64IlJuHEnKYE38KVbHn2J1fCIbDiaRVcQx//kJC/SjXkwYtSuGUCs6hJoVQ6n1178rhwfqQoCUPWkn4ZdRsHZywdvEtYau/weNrgGHX6mlVly8uTYQKU98ttg/mpzB16sO8Pmf+866rFCrGlG80LeZltQT6zhz4PRh9xI+SQf++Uo+CIn73QV9dtFmyD4rux/ENIG4Vu6TjbhWENsS/IOK7zlEpNzLzHGy+VAya/cnsvFgMpsOJbHzWAo5rvOYI6QAwf4OakaHULNiCNWigqkaFURsZDBVI4OIiwomJjwQf4dWAxGLHFwFs59yT35bkIhq0HYQtB3gVRP5eWNtIFIe+VSxn5SezZyNR/hu3UGW7ErgbK+sbqVQHrqiAX1aVtUyQVIycjIh9bh7wp6Uv76nHvvn36cPu4v604fBuPC7X/nyD4HKjcyFfZVmPjVJkIh4j4xsJzuOprDpUBKbDrkvAGw5fDrfiXKLg90GlcMDiYsMJi4yiNjIICqFBVIpLOCv74FUCnf/HOjnKJEcpJwzDNj0Lcx9DhLjC97O5nBPetu0LzTqDaEVSy3F8+EttYFIeefVxX6O08WWw6f5fcdxFm47zqr4UzjPccegRnQwD/VoSL/WVfHT1X4pDJcTMk9DRiKkn/L8SssbO+nudp+RWHo5/l3UV27i/h7TBCo3hsgaWuNeRMo0p8vg4Kl0dhw7zY5jKew4msLO4ynsPHqa1KySuQiQn/AgPyr/dQEgKsSfyGD/3O+RIQFEBZtjUcEBhAf56YaBFI4zG9ZPg0VvnHuyXJsDanWGut3dk+JWbQOOsrWCRVmtDUTEzGuKfcMwOHAqnbX7E1m3P5F1B9zdAgt7N6BdrQoM6VKbq5rFqkufrzEMcGZBdjrkZBT+e1aaez35zNOQmQJZKe7vmcln/Ps05KRb/QrdHAHuSfKi6/719de/K9aDyJoq6kXEpxiGweGkDHYcS2HP8RT2nUwjPiGNvQmp7D+VXizzARSHkAAHoYF+hAX6ERroIDTAj9BAv79i//wcFuhHcICDIH8HQf52Av3c34P8HQT6ub8H+TkI9Lfnfg/0s2tOAl/jzIFNM2Dp/+DQ6sLtExDmLvirtnb30qvaxn0+YOHf/bJQG4jIuRWp2D9y5Mh5/UL7+/vj55f/5CMZGRkUIgX2JqRx9dtLi/S8oQEOrmxamTsurk7zqoXP2+FwEBCQ/0R9mZmZuFznd4JxtsfNysrC6Ty/Oxh2u53AwPy7ZV/I49psNoKC8h/HfSGP67/gBfxy0txXuV3ZZ3zPwZWTieHMAmcONtcZ7c4c979d2dhcOX/FssGVA85MbCXVDb40+QW5x+5FVsOIqE5OaBWM8GoYFWq7v8LiwF64bqYF/c65XC4yMzPPO8Wz/S6np5//RZHi+IzIz9l+50rqcfUZUbKPm5OTQ3Z29nk9LkBwcHCJPG5gYCD2fE68S+pxS+p32Vs+I1wugyPJGexLSGPX0ST2n0zjyOlMjiRlcDg5k6PJmWSWkYsBFyrQz130+9lt+Dns+Dvs+Dts+Dvsf/1sw89uw2EDh9321892/P6K+zls+Nvd24UG+vFYz/oez6HPiH+U6mfEoTW4VnyMbeN0bNlF+/3IuHc5RmRNyz4jkpOTiY2NVbEvUsYVafrPrVu3EhYWVuQniYqKol69evm2bdq0qVCPYRgGof42UrPPfnIe4LBzacPKtIrKpF1cIIEOIPEAmxMLn29QUBDNmjXLt23z5s3nfSLv7+9Py5Yt823bunXref8RsdvttGnTJt+2HTt2kJFxfssh2Ww22rZtm2/bnj17SElJOa/HbbX6c8hOzrfNF+9NGzYH2YEVyA6MJuev79mBFcgOiia6dgvCqjaGyOoQUjF3GZ7EU6fYvXu3+wFSgdRkIP/3LD/Vq1enSpUqHvHMzEw2b9583q8lJiaGGjVq5Nt2IY9bHJ8R+QkLC6NRo0b5tm3evPm8i319RriV1GcEQLt27fKN79u3j8TExGJ/3MOHD3PsWAFLZRVC06ZN8y0SEhISOHDgwHk/bt26dalQoYJHPCkp6Z/PiPPg7Z8RdruNqlHBVI0KJiBxL03igDiAICAIwzA4nWVwIs1JQrqThDQXJ9KdnEp3kpTpIjHDRVKmi6QMF87imy+wRGTmuIrtwkVYgI3e1bI84vqMOPfjlshnRNU2HO/4DIdibyHy6FIqHP6dyKPLsLvOXqTn+Eew6cBpOLjZss+IC/m/E5HS4zVrfdhsNupH+7PuqOcfqbqVQrm0YWW6NaxMh7rRhAT4sXr16vM+kZeSZRTy7nSZZPeD4Aq5X5n2EE47A3AGhJPjH44zIIIc/7C/inp3Ye/0Dwdb/pcxwuvWhXz+SIuIyPmz2WxEBNqICLRTt0LBY51DQ0OJrVmXEymZHD+dRUJqJidOZ5KUnsP2fQdIzXJxOssgNctFSpaLlCyDlCwXOV56euGnIQFlkssvhFPVLudUtcux56QTfnwVEQlrCD+xhuDTez22T41qmHtzQETkbLym2AdoEO3P7lPZ1I/2p0G0Py2rRXBNx+ZEh+bf7VXKJsNWSoedI9C9nJxfMDk2P7Lxw2UPwOUIxLAH4PQLwuUXgvOvL/e/g3H6heLyC/4r7v53vSatCIqMgcBw0x/YxKNHL+iunYiIWMdmsxEVEkBUSAD1Y8xtq1efzvemgWEYZDoNUrMMUrNdpOcYpGcb7u85LjKyDbJxEBZVkdTMHFIynaRm5pCalUN6lpOTySlk5bjIckGW0yDLaZDtNMguhVEHmrKo7HP5BZMUdwlJcZcA4Jd5ktBTWwhN3EZI0nZCEneQFtnQ4ixFxFt4VbF/U5MwbmsWljtZTVhYqAp9L5RQ/UriKoS4Z5Z1+IP97+9+nEg8TXpmNobdD8PmwGX3x7D5Ydgdf33/6+vvf9scGI4AXPZAXI5AWrZt7x777hdkmrjm8P79F9T9zoioBkH5j+MTEZHyw2azEeRnI8gPKpJ/TzX3UJ/G+batWbMm36E+LsMg2wlZLsN0ESDL6b4o4DQMsDmoUas22U6DHJeLHKdBttNFttNg3/79ZGU7yTHcKxw4XZBjGOS43D/nuCDEX3eDvU1OYDRJsV1Iiu3iDhgGNiPH2qRExGt4zQR9+dHkW26aWOcfmnzLzdsn3zqTPiPc9Bnhps8IN31G/EOfEW76jHDTZ4SbJugTEfCipfdERERERMR6qg1EvINGb4mIiIiIiIj4GBX7IiIiIiIiIj5Gxb6IiIiIiIiIj1GxLyIiIiIiIuJjVOyLiIiIiIiI+Jj817HJ4+8J+5OTk0s0GRERERERKdv+rgnOd3lKESkdhSr2T58+DUCNGjVKNBkREREREfEOp0+fJjIy0uo0RKQANqMQl+RcLheHDh0iPDwcm812Xk+UnJxMjRo12L9/v9bjlAum40mKk44nKW46pqQ46XiS4lQcx5NhGJw+fZqqVatit2tUsEhZVag7+3a7nerVqxfLE0ZEROgPlRQbHU9SnHQ8SXHTMSXFSceTFKcLPZ50R1+k7NOlOBEREREREREfo2JfRERERERExMeUWrEfGBjI6NGjCQwMLK2nFB+m40mKk44nKW46pqQ46XiS4qTjSaT8KNQEfSIiIiIiIiLiPdSNX0RERERERMTHqNgXERERERER8TEq9kVERERERER8jIp9ERERERERER+jYl9ERERERETExxRrsf/OO+9Qu3ZtgoKC6NChA8uXLz/r9tOnT6dx48YEBQXRokULfvrpp+JMR7xcUY6niRMn0rVrVypUqECFChW44oorznn8SflS1M+nv02dOhWbzUa/fv1KNkHxKkU9nhITExkxYgRxcXEEBgbSsGFD/c0Tk6IeU2+99RaNGjUiODiYGjVq8PDDD5ORkVFK2UpZ9fvvv9OnTx+qVq2KzWZj5syZ59xnwYIFtG3blsDAQOrXr8+nn35a4nmKSOkotmL/q6++4pFHHmH06NGsXr2aVq1acdVVV3Hs2LF8t1+yZAm33347d999N2vWrKFfv37069ePjRs3FldK4sWKejwtWLCA22+/nfnz5/Pnn39So0YNrrzySg4ePFjKmUtZVNTj6W979+7l0UcfpWvXrqWUqXiDoh5PWVlZ9OzZk7179/L111+zbds2Jk6cSLVq1Uo5cymrinpMffnllzz55JOMHj2aLVu28NFHH/HVV1/x9NNPl3LmUtakpqbSqlUr3nnnnUJtv2fPHnr37s1ll13G2rVrGTlyJEOHDmXOnDklnKmIlAqjmLRv394YMWJE7s9Op9OoWrWqMXbs2Hy3v+WWW4zevXubYh06dDDuueee4kpJvFhRj6e8cnJyjPDwcGPSpEkllaJ4kfM5nnJycozOnTsbH374oTFo0CCjb9++pZCpeIOiHk/vvvuuUbduXSMrK6u0UhQvU9RjasSIEcbll19uij3yyCNGly5dSjRP8S6AMWPGjLNu8/jjjxvNmjUzxW699VbjqquuKsHMRKS0FMud/aysLFatWsUVV1yRG7Pb7VxxxRX8+eef+e7z559/mrYHuOqqqwrcXsqP8zme8kpLSyM7O5vo6OiSSlO8xPkeT88//zwxMTHcfffdpZGmeInzOZ6+//57OnXqxIgRI6hSpQrNmzdnzJgxOJ3O0kpbyrDzOaY6d+7MqlWrcrv67969m59++olrrrmmVHIW36HzcRHf5lccD3LixAmcTidVqlQxxatUqcLWrVvz3efIkSP5bn/kyJHiSEm82PkcT3k98cQTVK1a1eMPmJQ/53M8LV68mI8++oi1a9eWQobiTc7neNq9eze//fYbd955Jz/99BM7d+7k/vvvJzs7m9GjR5dG2lKGnc8xdccdd3DixAkuueQSDMMgJyeHe++9V934pcgKOh9PTk4mPT2d4OBgizITkeKg2fjF57z88stMnTqVGTNmEBQUZHU64mVOnz7NgAEDmDhxIpUqVbI6HfEBLpeLmJgYPvjgA9q1a8ett97Kv//9b9577z2rUxMvtWDBAsaMGcP//vc/Vq9ezbfffsusWbN44YUXrE5NRETKkGK5s1+pUiUcDgdHjx41xY8ePUpsbGy++8TGxhZpeyk/zud4+tvrr7/Oyy+/zNy5c2nZsmVJpileoqjH065du9i7dy99+vTJjblcLgD8/PzYtm0b9erVK9mkpcw6n8+nuLg4/P39cTgcubEmTZpw5MgRsrKyCAgIKNGcpWw7n2Nq1KhRDBgwgKFDhwLQokULUlNTGT58OP/+97+x23UvRwqnoPPxiIgI3dUX8QHF8tcgICCAdu3aMW/evNyYy+Vi3rx5dOrUKd99OnXqZNoe4Ndffy1weyk/zud4Anj11Vd54YUXmD17NhdddFFppCpeoKjHU+PGjdmwYQNr167N/bruuutyZyquUaNGaaYvZcz5fD516dKFnTt35l40Ati+fTtxcXEq9OW8jqm0tDSPgv7vi0mGYZRcsuJzdD4u4uOKa6a/qVOnGoGBgcann35qbN682Rg+fLgRFRVlHDlyxDAMwxgwYIDx5JNP5m7/xx9/GH5+fsbrr79ubNmyxRg9erTh7+9vbNiwobhSEi9W1OPp5ZdfNgICAoyvv/7aOHz4cO7X6dOnrXoJUoYU9XjKS7Pxy5mKejzFx8cb4eHhxgMPPGBs27bN+PHHH42YmBjjxRdftOolSBlT1GNq9OjRRnh4uDFlyhRj9+7dxi+//GLUq1fPuOWWW6x6CVJGnD592lizZo2xZs0aAzDefPNNY82aNca+ffsMwzCMJ5980hgwYEDu9rt37zZCQkKMxx57zNiyZYvxzjvvGA6Hw5g9e7ZVL0FEilGxFfuGYRgTJkwwatasaQQEBBjt27c3li5dmtvWrVs3Y9CgQabtp02bZjRs2NAICAgwmjVrZsyaNas40xEvV5TjqVatWgbg8TV69OjST1zKpKJ+Pp1Jxb7kVdTjacmSJUaHDh2MwMBAo27dusZLL71k5OTklHLWUpYV5ZjKzs42nnvuOaNevXpGUFCQUaNGDeP+++83Tp06VfqJS5kyf/78fM+H/j5+Bg0aZHTr1s1jn9atWxsBAQFG3bp1jU8++aTU8xaRkmEzDPX3EhEREREREfElmsFFRERERERExMeo2BcRERERERHxMSr2RURERERERHyMin0RERERERERH6NiX0RERERERMTHqNgXERERERER8TEq9kVERERERER8jIp9ERERERERER+jYl9ERERERETEx6jYFxEREREREfExKvZFREREREREfMz/A5cjxd2bmk+BAAAAAElFTkSuQmCC",
+      "text/plain": [
+       "<Figure size 1000x300 with 1 Axes>"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
     }
    ],
    "source": [
-    "BF_smc = np.exp(trace_smc_0.sample_stats[\"log_marginal_likelihood\"].mean() - trace_smc_1.sample_stats[\"log_marginal_likelihood\"].mean())\n",
-    "print(np.round(BF_smc.item(),2))"
-   ]
-  },
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "Der Bayes-Faktor ist etwa 11 (es wurden 2 Simulationen (chains)  erzeugt). Der BF is somit grösser als 3, woraus wir schliessen, dass die erste Prior-Verteilung besser zu den Daten passt, was hier keine \\\"Uberraschung ist. \n",
-    "\n",
-    "\n",
-    "- Wir verwenden `pm.sample_smc()` anstatt `pm.sample()`.\n",
-    "- Per default werden 8 Simulationen durchgeführt. Mit der Option `chains=...` kann dies geändert werden.\n",
-    "- `pymc` kennt nur die log-Marginal-Verteilung, da diese bessere numerische Eigenschaften hat.\n"
+    "!pip install preliz\n",
+    "import preliz as pz\n",
+    "priors = ((alpha_1, beta_1),(alpha_2, beta_2))\n",
+    "with plt.rc_context({\"lines.linewidth\": 3}):\n",
+    "    for a, b in priors:\n",
+    "        ax = pz.Beta(a, b).plot_pdf()"
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 15,
-   "metadata": {
-    "tags": []
-   },
+   "execution_count": 13,
+   "metadata": {},
    "outputs": [
     {
      "name": "stderr",
@@ -202,27 +218,20 @@
     }
    ],
    "source": [
-    "trials = 30\n",
-    "head = 14 \n",
-    "\n",
-    "data = np.zeros(trials)\n",
-    "data[np.arange(head)]  = 1\n",
-    "y_d = data\n",
-    "\n",
     "with pm.Model() as model_BF_0:\n",
-    "   θ = pm.Beta('θ', 10, 10)\n",
+    "   θ = pm.Beta('θ', alpha_1, beta_1)\n",
     "   y = pm.Bernoulli('y', θ, observed=y_d)\n",
-    "   trace_BF_0 = pm.sample_smc(chains=2)\n",
+    "   trace_smc_0 = pm.sample_smc(chains=2)\n",
     "\n",
     "with pm.Model() as model_BF_1:\n",
-    "   θ = pm.Beta('θ', 200, 200)\n",
+    "   θ = pm.Beta('θ', alpha_2, beta_2)\n",
     "   y = pm.Bernoulli('y', θ, observed=y_d)\n",
-    "   trace_BF_1 = pm.sample_smc(chains=2)"
+    "   trace_smc_1 = pm.sample_smc(chains=2)"
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 16,
+   "execution_count": 14,
    "metadata": {
     "tags": []
    },
@@ -231,7 +240,7 @@
      "name": "stdout",
      "output_type": "stream",
      "text": [
-      "11.21\n"
+      "11.24\n"
      ]
     }
    ],
@@ -239,6 +248,18 @@
     "BF_smc = np.exp(trace_smc_0.sample_stats[\"log_marginal_likelihood\"].mean() - trace_smc_1.sample_stats[\"log_marginal_likelihood\"].mean())\n",
     "print(np.round(BF_smc.item(),2))"
    ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Der Bayes-Faktor ist etwa 11 (es wurden 2 Simulationen (chains)  erzeugt). Der BF is somit grösser als 3, woraus wir schliessen, dass die erste Prior-Verteilung besser zu den Daten passt, was hier keine \\\"Uberraschung ist. \n",
+    "\n",
+    "\n",
+    "- Wir verwenden `pm.sample_smc()` anstatt `pm.sample()`.\n",
+    "- Per default werden 8 Simulationen durchgeführt. Mit der Option `chains=...` kann dies geändert werden.\n",
+    "- `pymc` kennt nur die log-Marginal-Verteilung, da diese bessere numerische Eigenschaften hat.\n"
+   ]
   }
  ],
  "metadata": {
-- 
GitLab