diff --git a/module3/exo3/Around Simpson's Paradox.ipynb b/module3/exo3/Around Simpson's Paradox.ipynb index a1ab52efb216e2990948010361d7788244d59ba6..5ea54bfdcf0b5527f129ef916d42b4ddb1d2744e 100644 --- a/module3/exo3/Around Simpson's Paradox.ipynb +++ b/module3/exo3/Around Simpson's Paradox.ipynb @@ -2,7 +2,7 @@ "cells": [ { "cell_type": "code", - "execution_count": 49, + "execution_count": 1, "metadata": {}, "outputs": [], "source": [ @@ -10,7 +10,8 @@ "import pandas as pd\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", - "import seaborn as sns" + "import seaborn as sns\n", + "from sklearn.linear_model import LogisticRegression" ] }, { @@ -22,7 +23,7 @@ }, { "cell_type": "code", - "execution_count": 50, + "execution_count": 2, "metadata": {}, "outputs": [ { @@ -102,7 +103,7 @@ "4 Yes Alive 81.4" ] }, - "execution_count": 50, + "execution_count": 2, "metadata": {}, "output_type": "execute_result" } @@ -123,7 +124,7 @@ }, { "cell_type": "code", - "execution_count": 51, + "execution_count": 3, "metadata": { "scrolled": false }, @@ -151,7 +152,7 @@ }, { "cell_type": "code", - "execution_count": 52, + "execution_count": 4, "metadata": { "scrolled": false }, @@ -194,7 +195,7 @@ }, { "cell_type": "code", - "execution_count": 53, + "execution_count": 5, "metadata": { "scrolled": true }, @@ -237,7 +238,7 @@ }, { "cell_type": "code", - "execution_count": 54, + "execution_count": 6, "metadata": {}, "outputs": [ { @@ -264,7 +265,7 @@ " Smoker\n", " Status\n", " Age\n", - " Binned age\n", + " Age group\n", " \n", " \n", " \n", @@ -308,7 +309,7 @@ "" ], "text/plain": [ - " Smoker Status Age Binned age\n", + " Smoker Status Age Age group\n", "0 Yes Alive 21.0 18-34 years\n", "1 Yes Alive 19.3 18-34 years\n", "2 No Dead 57.5 55-64 years\n", @@ -316,7 +317,7 @@ "4 Yes Alive 81.4 Over 65 years" ] }, - "execution_count": 54, + "execution_count": 6, "metadata": {}, "output_type": "execute_result" } @@ -334,7 +335,7 @@ " else:\n", " return \"Over 65 years\"\n", "\n", - "data['Binned age'] = data['Age'].apply(bin_age)\n", + "data['Age group'] = data['Age'].apply(bin_age)\n", "data.head()" ] }, @@ -347,7 +348,7 @@ }, { "cell_type": "code", - "execution_count": 55, + "execution_count": 7, "metadata": {}, "outputs": [ { @@ -362,44 +363,111 @@ "print(\"Number of rows with missing values:\", data.isnull().any(axis=1).sum())" ] }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Let's visualize the number on women alive and dead after twenty years, according to their smoking habits and age. Different colors correspond to different couples of smoking habits and status." + ] + }, { "cell_type": "code", - "execution_count": 56, + "execution_count": 8, "metadata": {}, "outputs": [ { "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAA3EAAAKZCAYAAAAWFKIXAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzs3XmYFdWZ+PHvC4gQQUFtFUUFd5RFtAUyLoiiomJQRwQmLkQTE7fEJKODSyIa/WmcxKjJZNQYhTiK4DaJiInJiEmIiqLgCnEBBZQgEFFQkYDn90cV7aXp5Ta0tpX+fp6nn7636tSpt07VrXvfe07VjZQSkiRJkqRiaNHUAUiSJEmSymcSJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEid9SiJiTERc0UTrjoi4LSLeiYgn12P51yNiYP74ooi4pfGjXGedIyNiyqdUd9X2NDcRcUVELI6IvzV1LPr0RESKiF3yxzdGxPeaOqZyRMToiPif9Vz2Mzk3NVRd55uIeDEiDv4U1rne58+IODgi5tcx/3N5PDXn87oEJnFqRvIT/sKI2KRk2lcj4tEmDOvTcgBwGNA5pdRnQypKKf2/lNJXGycsiIiXI2K3xqpPtYuI7YHvAnumlLaJiC75h/1WjbiO3SLi1xGxKCL+HhG/i4jdq5X5dkT8LSLejYhbI2Ljxlq/1pVS+kZK6QeNXe+GHj/1JQsN1djnps9CSmmvlNKjTR1HQ5QeT429DyWtP5M4NTetgG81dRANFREtG7jIjsDrKaX3P4141ldE7Ay0SCm93NSxNBM7AktSSm83RmW1fHjvAPwG2B3YGngS+HXJMkcAo4BDgS7ATsBljRHP+mrMJPazVuTYP0u2k6R/diZxam7+E/j3iOhQfUZN3zJHxKMR8dX88ciI+EtE/CQilkbE7Ij4l3z6vIh4OyJOrVbtlhHx+4hYFhF/jIgdS+reI5/394j4a0ScWDJvTET8d0RMioj3gQE1xLttRPwmX/7ViPhaPv104BbgixGxPCLW+cAcETtHxCMRsSQfandHTW2Sl60a7hQRv42Ic6rNfzYijq9vm3JHA5Pyslvk8b+XD/ncuVq9dbXP0RExPV92XkSMrrbsyRHxRr59F9e0XXm5rvm+bJE/vyUi3i6Z/z8RcV5d7V3SRnfn5ZdFxPN5D9WF+XExLyIOLym/WUT8MiIWRMSbkQ15bJnPGxkRUyLiR5ENh50TEUfWsQ2jIuK1fL0vRcRx+fSBwO+BbfPjYAzwp3yxpfm0L+ZlT4uImfn6flftOE0RcXZEvAK8Un39KaUnU0q/TCn9PaX0D+AnwO4RsUVe5FTglymlF1NK7wA/AEbWsi3/FRE/rjbtgWr74N7Iev3mRMQ3S8r1iYjH8/25ICJ+FhGta9uOyPwk3z/vRsRzEdG9tnauFtP1+T59LyKejogDS+a1jGyY35p98nRkPaJExF4lx/TCiLgon75xRFwXEW/lf9dF3lsZec9HRPxHZENib8unn59v51sRcVq1+KqGcpcs/918WxdExFdKym6Rt/F7EfFUfizWNixvneMnIlpExCWRvd7ejohfRcRmNbTZJsBDfHI8Lo+IbfPZrfPllkU23LCyZLm69nnpuWnN+fv0iJgLPFJDDFtGxMT8GPl7RPw5Pnntv5636XMR8X5kr8+tI+KhPK4/RETHkrq+lMe6NLL3iW41NVhk57E5ETG8ZD1rhqqPjogJdWz7PpGd55ZFdn4ZH/UM0Y9azhsR8ZXIXuPLInvv+noNy14U2fvB6xHx5ZLpY/LjosZ9GNlrb1p+DC2MiGtria1j3v6L8hgnRkTnkvmPRsQPInufXRYRD0fEliXzyzqv52Xre484paSu71XbLy3ik/PqknwfbV7X+qQmkVLyz79m8Qe8DgwE7gOuyKd9FXg0f9wFSECrkmUeBb6aPx4JrAK+ArQErgDmAv8FbAwcDiwD2uXlx+TPD8rnXw9MyedtAszL62oF7AMsBvYqWfZdYH+yL1va1LA9fwR+DrQB9gYWAYeWxDqljrbYhWy45cZABdmHs+uqt1X+eDTwP/njU4C/lJTbE1ia11PnNuXlfwsckT++C5iQL9cdeLMB7XMw0CNvm57AQuDYkpiWl7T7tfl+G1hLW8wF9s0f/xWYDXQrmde7jPYeDawAjsjj/RUwB7gY2Aj4GjCnZJ3/C9yUb+dWZL1XXy/Zd//Il2kJnAm8BUQt8Q8Fts3bYhjwPtCppJ3ml5TtwrrH+LHAq0C3PPZLgMdK5ieyZHBzoG0Zr7NjgQUlz58FhpU83zKvc4salu2Tb2uLkrIfkPXwtQCeBr4PtCbr0ZvNJ8fTvkC/fBu6ADOB82rbjnxfPU3Wkxj59ncq81xyErBFvq7vAn8jf40C5wPPk/VMBtArL9seWJCXb5M/75svcznwRH4sVACPAT8o2YergB+SHc9tgUFkx3z3/Bi6M9++XUrOH1dUW/5ysmPxqLxNO5a8Du8CvkD22plHLeeOWo6f0/LjZyegHdn59fZalj+YkuOx2mvnKLLj/SrgiXxefft8NJ+cm9bE9qu8TdY5VvO6b8zbYSPgQPLXFdk57wmyY2074G3gGaB33u6PAJfmZXcje50dltdzQd4Grau91+xDdg4ZXMe5tbZtbw28QTZyZCPgeGDlmv1aw7aNpI7zBtkXaDuTHZP982Ngn2rHyLX5tvbPt2/3Wo6n6vvwceDk/HE7oF8tMW4B/CvZsdYeuBv435L5jwKv5e3bNn9+9Xqe1w+m/veIA/J2/lHedmv2y3n5sdA5X9dNwLhyzg3++fdZ/jV5AP7591n9lbyxdidLkCpoeBL3Ssm8Hnn5rUumLQH2zh+PAe4qmdcOWA1sT/Zh+8/V4ruJTz4kjAF+Vce2bJ/X1b5k2lXAmJJYa03iaqjvWGB69bbKH4/mkw9K7cne3HfMn18J3Jo/rm+bvpC3TxuyDxn/APYoKfv/+CSJq7OuGuK/DvhJ/vj71dp9E7IPP7W92d8OfAfYhiyJuwb4BtCVLEFtUUZ7jwZ+XzLvGLIPCS1L2i2RJQxbAx9R8iETGAFMLtl3r5bM+0K+7DZl7ssZwJD88cHUn8Q9BJxe8rwF2Qe8Nfs4AYeUue7OZMn4iJJprwGDSp5vlNfZpZY6ZgKH5Y/PASblj/sCc6uVvRC4rZZ6zgPuL3m+1nYAhwAvkyV+Lcp9rdSyrneAXvnjv65p/2plRlDyGqs27zXgqJLnR5ANh16zD1dS8kUOcCv5h9v8+W7UncR9WG2fv51v95rX4e4l866gYUnc/wFnlTzfPa+zVQ3Lr3U8lrx2/lDyfE/gw3L2OTUncTvVsZ8uJxvqu0sN814Hvlzy/F7gv0uen0uecADfAyZUe828CRxcUtdlwHxgQA3rGVjGth+U1xkl86dQdxJX9nmD7Iukb5Xsl1XAJiXzJwDfq+V4qr4P/5Rv75YNfN3sDbxT8vxR4JKS52cBv80fN+i8XsO6qr9HjCuZ94XSusjOQYeWzO9ELce0f/415Z/DKdXspJReACaSXafTUAtLHn+Y11d9WruS5/NK1rsc+DtZr8mOQN98KM7SiFgKfJkskVhn2RpsC/w9pbSsZNobZN8g1ysitoqIuyIbyvce8D9kvR51ytf3IDA8nzQcuCN/XN82HUrWw7OCLIFuVW0b3yh5XGddEdE3Iibnw3LeJUu61sS/LWu3+/tkyWNt/kj2weQgsg8jj5J9E92fLJH8mPLau/pxsDiltLrkOWTHxo5kicyCkm27iawXZo2qO0mmlD4oWXYd+bCgGSV1daeMfVliR+D6kuX/TvZtfem21XUsromjAngY+HlKaVzJrOXApiXP1zwubctSY8l6usj/314S57bVjomLyJLiNTdYmRjZDVTeI/tSoHo7lB4XjwA/I+tJXxgRN0fEppQhsqGJMyMbhrkU2KxkXduTJWXV1TYdsuOr9Ph/I5+2xqL8dVNavrbXTk2WpJRWlTz/gOx4qul1WO++rqam2FuR75cyld459QOgTWTD2uvc57WoK/7/JOsxezgfUlj9PaD6a7i2c/ta25yfI+ax9mvmG2Tnu8l1xAO1b/u2wJsppVQyv759U+t5IyKOjIgnIhtGupSs96/09fFOWvsa6urHYF1OJ/siYVZkQ3IH11QoIr4QETflwxjfIzvfdoi1r/mu3h6lbV72eb2B7xEfVKtrR+D+kmNuJtmXeA05pqVPnUmcmqtLyYadlL7prnkD+0LJtNKkan1sv+ZBRLQjG8r1FtkbyB9TSh1K/tqllM4sWTZRu7eAzSOifcm0Hci+uS3HVXn9PVNKm5J9WI4ylx0HjIjseqq2wJoPKfVt01FkCSBkQxFXUdI+efyUWdedZDfT2D6ltBnZEKk18S9g7Xb/Atkwntr8kWxY1cH54ylkw1j7589hw9u71DyynrgtS7Zt05TSXg2tKLJr135B1mO1RUqpA/ACte/Lmo6peWRDOUvbum1K6bF6liuNoyNZAveblNKV1Wa/SDakcI1ewMKUUm0fwP4HGBIRvciGOP5vSZxzqsXZPqV0VD7/v4FZwK75MX0R67bDWtuRUrohpbQvsBfZh9Dz69rOfFsPBP4DOJFsSGIHsp79NeuaR7XrO+uZDtnxtWPJ8x3yaTXGTbVjnLVfOw2x5nXYuWTa9rWUrSkOqDn2VaydANW1fF3q2+flxpjNSGlZSum7KaWdyHrLvxMRhzYwJqi2zRERZO1Wej74BrBDRPxkPeqHbB9vl9e9Rl37plaRXV95L9mwwa3zY3YSa78+OkbJnZtZ9xhcY532TSm9klIaQfZF1A+Be6rVtcZ3yXpq++av0YPWhFjGZjT0vF7fe0TptXhtq9U1Dziy2nHXJqW0Pud76VNjEqdmKaX0KjAe+GbJtEVkb8InRXZzgtOo/UNXuY6KiAMiu8HCD4CpKaV5ZD2Bu+UXam+U/+0XtVwcX0P888ium7kqItpERE+yb0PvqHvJKu3JekiWRsR2lPHhtcQksg8wlwPj82+hKWObjsyXJe+hug8YnX87uyfZDTAos672ZD1jKyKiD/BvJcveAwwuaffLqeNcl1J6hexb9pOAP6WU3iP7APqv5ElcI7R36foWkCU8P46ITfOL6HeOiP4NrYtsSFEi+zBOZDesqOvmHIuAj8muLVrjRuDCiNgrr2OziBhabgB579XvyK6VrKl3+1fA6RGxZ57sXUI2PKtGKaX5wFNkPXD3ppTW9GI+CbwX2Q0+2uav0e4RsV8+vz3wHrA8IvYguyaorrj3y7+t34jsC5wVZN+2r7m5zOu1LNqeLElZBLSKiO+zdk/jLcAPImLXyPSM7CYvE4FtIuK8yG5k0j4i+ubLjAMuiYiKyG7k8H2yZLY2E4CReZt+gexLqQar4XW4B9l1r7Wp6fgZB3w7spsEtSPrAR1fredvjYXAFlHDjU9qUd8+b5CIGBwRu+SJ0Xtk+3t1PYvVZAJwdEQcmh8/3yX7Yqb0i49lZNcuHhQRV6/HOh7PYzsnIlpFxBCya0bXR2uya7sWAasiu+HJ4TWUuywiWudfVAwmu2atunX2YUScFBEV+XvB0nxyTe3anuxcuzSyG4U05Lht0Hmd+t8jjonsxmStyYaCliaSNwJX5l+Skb8uhzQgVukzYRKn5uxysg/Bpb5GltAsIft2/rHqCzXQnWRvVH8nu/HCl6FqWOLhZMMR3yIbQrLmxgXlGkF2HchbwP1k14v9vsxlLyO76P5dst6x+8pdaUrpo7z8QLLtWzO91m2K7K5/y1NKc0uqOodsqMzfyD7U31ZOXXmRs4DLI2IZ2QfeCSXLvgicnce2gOx6pfp+1+iPZEPO5pY8D2B6SZkNae/qTiH7YPVSHt89ZNddNEhK6SXgx2Qf+BaSXaf5lzrKf0B2HeNf8qFC/VJK95O17V2RDXF6gSzhLtdxwH7AV+KTO9Ytj4gd8nX+luw6w8lkQ7TeoP4Pb2PzbVkzlHJNwnEM2XU0c8hudHML2VBGgH8n+6C2jKx3cnw969g0L/dOHtMSsp4KyL7xr60df0d2HeHL+XIrWHuY27Vkx+PDZInCL8muf1xGdiOMY8iO51f45K6zVwDTgOfIboryTD6tRimlh8iu8XmEbHjgOndibIBzyNrwb2TtPY4sIalpvescP2TX591ONjRuDll7nFvL8rPy+mfny9c5XK+Mfd5QuwJ/IPsC63Gyob+PNrSSlNJfyb70+Wke0zHAMSmlldXKLSXb50dGRIN+ty+v63iyL4uW5uubSC37pp66lpF9YTmB7Hj/N7JeqlJ/y+e9Rfbl1Dfy/VW9rpr24SDgxYhYTnYDr+HVhv+ucR3Z6I3FZDcO+W0DtqGh5/X63iPOJbuhzwKyc8bbfNK215O1z8P58k+QXZ8pfa6suWuRJH1qIuICsuGDFzR1LPr8i4iDyHqiupT09H6W63+Y7KYPMz/rdTe1iPgh2c0wTq23sD5TETEVuDGldFu9hVW2vAd5KdlQ7DlNHY9ULnviJH0WXqekp02qTT487VvALU2RwAGklA5vLglcZL9j1jMf+tmHrOfn/qaOSxAR/SNim3w45alkt8ovu/dKtYuIY/IhxJuQ9cA/T/Y+JRVGq/qLSNKGSSlNqL+Umrv8msdpZL8t95V6iqtxtCcbHrct2ZCyH5Pdhl9Nb3eyYYDtyO5sekJ+Ta023BCyYcBBds4ZnhyapoJxOKUkSZIkFYjDKSVJkiSpQEziJEmSJKlAPhfXxG255ZapS5cuTR2GJEmSJDWJp59+enFKqaKcsp+LJK5Lly5MmzatqcOQJEmSpCYREW+UW9bhlJIkSZJUICZxkiRJklQgJnGSJEmSVCCfi2viJEmSJDWNf/zjH8yfP58VK1Y0dSjNQps2bejcuTMbbbTRetdhEidJkiQ1Y/Pnz6d9+/Z06dKFiGjqcP6ppZRYsmQJ8+fPp2vXrutdj8MpJUmSpGZsxYoVbLHFFiZwn4GIYIstttjgXk+TOEmSJKmZM4H77DRGW5vESZIkSSrblVdeyV577UXPnj3Ze++9mTp16gbV9+ijjzJ48OBGiq72+K677jo++OCDepcvt1xT8po4SZIkSWV5/PHHmThxIs888wwbb7wxixcvZuXKlU0Wz6pVq2jV6pOUpq74rrvuOk466SS+8IUv1FlnueWakj1xkiRJksqyYMECttxySzbeeGMAttxyS7bddlu6dOnCRRddxBe/+EUqKyt55plnOOKII9h555258cYbgeymHueffz7du3enR48ejB8/fp36n3rqKXr37s3s2bN5//33Oe2009hvv/3o3bs3v/71rwEYM2YMQ4cO5ZhjjuHwww8vK74bbriBt956iwEDBjBgwAAAzjzzTCorK9lrr7249NJLAWos165du6r677nnHkaOHAnA3XffTffu3enVqxcHHXRQYzVxeVJKTf637777JkmSJEmfvZdeeqnsssuWLUu9evVKu+66azrzzDPTo48+mlJKaccdd0w///nPU0opnXfeealHjx7pvffeS2+//XaqqKhIKaV0zz33pIEDB6ZVq1alv/3tb2n77bdPb731Vpo8eXI6+uij01/+8pe0zz77pDfeeCOllNKFF16Ybr/99pRSSu+8807adddd0/Lly9Ntt92Wtttuu7RkyZKy41sT46JFi6qer1l+1apVqX///unZZ5+tsdwmm2xS9fjuu+9Op556akoppe7du6f58+dXxdcQNbU5MC2VmT/ZEydJkiSpLO3atePpp5/m5ptvpqKigmHDhjFmzBgAvvSlLwHQo0cP+vbtS/v27amoqKBNmzYsXbqUKVOmMGLECFq2bMnWW29N//79eeqppwCYOXMmZ5xxBg888AA77LADAA8//DBXX301e++9NwcffDArVqxg7ty5ABx22GFsvvnmDYqvugkTJrDPPvvQu3dvXnzxRV566aUGtcX+++/PyJEj+cUvfsHq1asbtOyG8po4SZIkSWVr2bIlBx98MAcffDA9evRg7NixAFVDGFu0aFH1eM3zVatWkXU21axTp06sWLGC6dOns+222wLZiMF7772X3Xfffa2yU6dOZZNNNmlQfGuGQK4xZ84cfvSjH/HUU0/RsWNHRo4cWett/0vvJlla5sYbb2Tq1Kk8+OCD7L333syYMYMtttii1rgakz1xkiRJksry17/+lVdeeaXq+YwZM9hxxx3LWvaggw5i/PjxrF69mkWLFvGnP/2JPn36ANChQwcefPBBLrroIh599FEAjjjiCH76059WJX/Tp0+vte499tij3vjat2/PsmXLAHjvvffYZJNN2GyzzVi4cCEPPfRQ1TKl5QC23nprZs6cyccff8z9999fNf21116jb9++XH755Wy55ZbMmzevrHZoDPbESZIkSSrL8uXLOffcc1m6dCmtWrVil1124eabb2bixIn1Lnvcccfx+OOP06tXLyKCa665hm222YZZs2YBWbL0wAMPcOSRR3Lrrbfyve99j/POO4+ePXuSUqJLly41rmfx4sVViV5t8QGcccYZHHnkkXTq1InJkyfTu3dv9tprL3baaSf233//qvqql7v66qsZPHgw22+/Pd27d2f58uUAnH/++bzyyiuklDj00EPp1avXBrdvuaKubs3PSmVlZZo2bVpThyFJkiQ1OzNnzqRbt25NHcZ6mzhxIrNnz+ab3/xmU4dStpraPCKeTilVlrO8PXGSJEmSCqsxfyi8KLwmTpIkSZIKxCROkiRJkgrEJE6SJEmSCsQkTpIkSZIKxCROkiRJkgrEJE6SJEmSCsSfGJAkSZJUpcuoBxu1vtevPrreMh9++CGDBg3ikUceYd68eXTt2pUbbriBc889F4BzzjmHyspKRo4cud5xDBkyhLfffpvHH3+8atro0aNp164d//7v/873v/99DjroIAYOHNjguocPH84PfvADdt111/WOryFM4iR9Ls3co7g/OtpcdJs1s6lDkCT9k7j11ls5/vjjadmyJQBbbbUV119/PV//+tdp3br1Bte/dOlSnnnmGdq1a8ecOXPo2rXrOmUuv/zy9a7/zDPP5JprruEXv/jFhoRZNodTSpIkSWpSd9xxB0OGDKl6XlFRwaGHHsrYsWPXKTtjxgz69etHz549Oe6443jnnXfqrf/ee+/lmGOOYfjw4dx11101lhk5ciT33HMPDz30ECeeeGLV9EcffZRjjjkGgIcffpgvfvGL7LPPPgwdOpTly5cDcOCBB/KHP/yBVatWNWi715dJnCRJkqQms3LlSmbPnk2XLl3Wmj5q1Ch+/OMfs3r16rWmn3LKKfzwhz/kueeeo0ePHlx22WX1rmPcuHGMGDGCESNGMG7cuDrLHnbYYTzxxBO8//77AIwfP55hw4axePFirrjiCv7whz/wzDPPUFlZybXXXgtAixYt2GWXXXj22WcbsOXrzyROkiRJUpNZvHgxHTp0WGd6165d6dOnD3feeWfVtHfffZelS5fSv39/AE499VT+9Kc/1Vn/woULefXVVznggAPYbbfdaNWqFS+88EKt5Vu1asWgQYN44IEHWLVqFQ8++CBDhgzhiSee4KWXXmL//fdn7733ZuzYsbzxxhtVy2211Va89dZbDd389eI1cZIkSZKaTNu2bVmxYkWN8y666CJOOOEEDjrooPWuf/z48bzzzjtV18G999573HXXXVxxxRW1LjNs2DD+67/+i80335z99tuP9u3bk1LisMMOq7Unb8WKFbRt23a942wIk7hPiTdl+HzzhgySJEmfDx07dmT16tWsWLGCNm3arDVvjz32YM8992TixIn06dOHzTbbjI4dO/LnP/+ZAw88kNtvv72qV+5nP/sZkN3JstS4ceP47W9/yxe/+EUA5syZw2GHHVZnEnfwwQdz+umn84tf/IJhw4YB0K9fP84++2xeffVVdtllFz744APmz5/PbrvtBsDLL7/MXnvt1TiNUg+TOEmSJElVyvlJgMZ2+OGHM2XKlBpv73/xxRfTu3fvqudjx47lG9/4Bh988AE77bQTt912GwCzZs1i//33X2vZ119/nblz59KvX7+qaV27dmXTTTdl6tSptcbTsmVLBg8ezJgxY6purlJRUcGYMWMYMWIEH330EQBXXHEFu+22GwsXLqRt27Z06tRp/RuhASKl9JmsqC6VlZVp2rRpTR1Go7In7vPNnrjPP19Dn3++jiTpn8PMmTPp1q1p33enT5/Otddey+23377edQwePJj77ruvUX6SoKF+8pOfsOmmm3L66aeXVb6mNo+Ip1NKleUsb0+cJEmSpCbVu3dvBgwYwOrVq6t+K66hJk6c2MhRla9Dhw6cfPLJn9n6TOIkSZIkNbnTTjutqUNYb1/5ylc+0/X5EwOSJEmSVCAmcZIkSZJUICZxkiRJklQgXhP3KTnxQpv28+z5pg5AkiRJWk9mGpIkSZI+MXqzRq7v3XqLfPjhhwwaNIhHHnmEefPm0bVrV2644QbOPfdcIPsB78rKSkaOHNng1Y8ZM4bzzz+fzp07s3z5cnbaaScuvfRS/uVf/qXBddWkS5cuTJs2jU033ZSBAwfyyCOP0KrVp5tmOZxSkiRJUpO69dZbOf7446t+XmCrrbbi+uuvZ+XKlY1S/7Bhw5g+fTqvvPIKo0aN4vjjj2fmzMb9vdPWrVtz6KGHMn78+EattyYmcZIkSZKa1B133MGQIUOqnldUVHDooYcyduzYdcrOmDGDfv360bNnT4477jjeeeedBq1rwIABnHHGGdx8880AvPbaawwaNIh9992XAw88kFmzZgHwwAMP0LdvX3r37s3AgQNZuHAhAEuWLOHwww+nd+/efP3rXyelVFX3scceyx133NHg7W8okzhJkiRJTWblypXMnj2bLl26rDV91KhR/PjHP2b16tVrTT/llFP44Q9/yHPPPUePHj247LLLGrzOffbZpypZO+OMM/jpT3/K008/zY9+9CPOOussAA444ACeeOIJpk+fzvDhw7nmmmsAuOyyyzjggAOYPn06X/rSl5g7d25Vvd27d+epp55qcDwN5TVxkiRJkprM4sWL6dChwzrTu3btSp8+fbjzzjuVIQQOAAAgAElEQVSrpr377rssXbqU/v37A3DqqacydOjQBq9zTe/Z8uXLeeyxx9aq46OPPgJg/vz5DBs2jAULFrBy5Uq6du0KwJ/+9Cfuu+8+AI4++mg6duxYtWzLli1p3bo1y5Yto3379g2Oq1z2xEmSJElqMm3btmXFihU1zrvooov44Q9/yMcff9yo65w+fTrdunXj448/pkOHDsyYMaPqb821cueeey7nnHMOzz//PDfddNNaMUZErXV/9NFHtGnTplHjrc4kTpIkSVKT6dixI6tXr64xkdtjjz3Yc889mThxIgCbbbYZHTt25M9//jMAt99+e1Wv3M9+9jN+9rOf1bu+P/7xj9x888187WtfY9NNN6Vr167cfffdQNZD9+yzzwJZr992220HsNa1eQcddFDVdW8PPfTQWtfkLVmyhIqKCjbaaKMGt0NDlD2cMiJaAtOAN1NKgyNic2A80AV4HTgxpfROXvZC4HRgNfDNlNLvGjluSZIkSZ+GMn4SoLEdfvjhTJkyhYEDB64z7+KLL6Z3795Vz8eOHcs3vvENPvjgA3baaSduu+02AGbNmsX+++9fY/3jx49nypQpfPDBB3Tt2pV7772Xbt26AdlNVc4880yuuOIK/vGPfzB8+HB69erF6NGjGTp0KNtttx39+vVjzpw5AFx66aWMGDGCffbZh/79+7PDDjtUrWfy5MkcddRRjdYutYnSu6nUWTDiO0AlsGmexF0D/D2ldHVEjAI6ppT+IyL2BMYBfYBtgT8Au6WUVtdWd2VlZZo2bdqGbsvnSo+xPZo6BNXh+VP9ue/Pu5l7dGvqEFSPbrMa99bMkqSmMXPmzKqEpqlMnz6da6+9lttvv3296xg8eDD33XcfrVu3bsTIGub444/nqquuYvfdd6+zXE1tHhFPp5Qqy1lPWcMpI6IzcDRwS8nkIcCafsWxwLEl0+9KKX2UUpoDvEqW0EmSJEnSOnr37s2AAQPWuRNlQ0ycOLFJE7iVK1dy7LHH1pvANYZyr4m7DrgAKL2icOuU0gKA/P9W+fTtgHkl5ebn09YSEWdExLSImLZo0aIGBy5JkiTpn8dpp51W9WPfRdS6dWtOOeWUz2Rd9SZxETEYeDul9HSZddZ0q5Z1xmymlG5OKVWmlCorKirKrFqSJEmSmrdybmyyP/CliDgKaANsGhH/AyyMiE4ppQUR0Ql4Oy8/H9i+ZPnOwFuNGbQkSZIkNVf19sSllC5MKXVOKXUBhgOPpJROAn4DnJoXOxX4df74N8DwiNg4IroCuwJPNnrkkiRJktQMlf0TAzW4GpgQEacDc4GhACmlFyNiAvASsAo4u647U0qSJEmSytegJC6l9CjwaP54CXBoLeWuBK7cwNgkSZIkfcYa+6eyyvlppw8//JBBgwYxadIk9t13X+6++2569MjiuOaaa5g9ezY33nhjo8bVvXt3evfuvdbPGpx00kmccMIJHHvssXzlK19h1KhR63W3yUMPPZT77ruPzTbbrDFDrlLu3SklSZIk6VNx6623cvzxx7PJJptw3XXXcdZZZ5FS4s033+Smm27iqquuatT1Pffcc7Rq1YpHHnmEDz/8sMYyt91223r/XMC//du/NXrSWcokTpIkSVKTuuOOOxgyZAgAgwYNolOnTvzqV7/i29/+NqNHj6Zjx44AXH311fTp04eePXty+eWXA7Bs2TKOPPJIevXqRffu3bnnnnvqXd+4ceM45ZRTOOSQQ5g4cWKNZQ444ABmzJjBT3/6Uy666KKq6bfccgvf/va3ARg7dix9+vRh77335qyzzuLjj7NfZBsyZAh33nnn+jdIPUziJEmSJDWZlStXMnv2bLp06VI17brrruPiiy9m0aJFnHzyyQBMmjSJuXPnMnXqVGbMmMFjjz3GY489xqRJk+jSpQvPPvssL7zwAocddli965wwYQLDhg1jxIgRjBs3rs6yQ4cOXSsxHD9+PMOGDeOFF17g/vvv57HHHmPGjBmsWrWKu+66C4Att9ySZcuWsXTp0vVokfptyI1NJEmSJGmDLF68mA4dOqw1bdttt+WQQw5h8ODBVdMefvhhHnroIXr37g3A8uXLefnll+nbty+jRo1i1KhRHHPMMey///51ru/xxx+nc+fObLfddmy11VZ87Wtf49133631+rVtttmGzp07M23aNHbYYQfmzJlD3759uf7663nqqaeorKwEsuv6tt/+k19aq6ioYMGCBetsW2MwiZMkSZLUZNq2bcuKFSvWmd6iRQtatPhk4GBKiUsuuYTTTz99nbLTpk1j0qRJnH/++QwePHit4Y/VjRs3jhdeeKGq5++9997j/vvvZ+TIkbUuM2zYMCZMmECXLl3413/9VyKClBKnnXYaP/jBD2pcZsWKFbRt27bWOjeEwyklSZIkNZmOHTuyevXqGhO5UkcccQS//OUvef/99wGYP38+ixcv5s0336Rdu3acfPLJfOc73+GZZ54B4IILLuCBBx5Yq47Vq1dz77338tJLL/H666/z+uuvc99999U7pPKEE07gvvvu46677mLYsGEADBw4kAkTJrB48WIAlixZwty5cwH4+OOPWbx48Vo9c43JnjhJkiRJVcr5SYDGdvjhhzNlyhQGDhxYa5mjjjqKWbNm0a9fPwDat2/PnXfeyUsvvcSoUaNo0aIFrVu3rror5HPPPcfQoUPXqmPy5Ml07dqVrbfeumragAEDOOmkk1i4cGGt695iiy3YZZddeO2119hnn30A6NGjB5deeikDBw7k448/ZqONNuLGG29khx124Mknn+SAAw6gZcuW690mdYmU0qdScUNUVlamadOmNXUYjaqxf19DjaspTk5qmJl7dGvqEFSPbrNmNnUIkqRGMHPmTLp1a9r33enTp3Pttdeu9ZttGyKlxKBBg/jd737XKPU11Nlnn82JJ55I//79a5xfU5tHxNMppcpy6nc4pSRJkqQm1bt3bwYMGMDq1asbpb6IaLIEDrLtqS2BawwOp5QkSZLU5E477bSmDqHRfPWrX/1U67cnTpIkSZIKxCROkiRJkgrEJE6SJEmSCsQkTpIkSZIKxBubSJIkSarS2D/zU85P0nz44YcMGjSISZMmse+++3L33XfTo0f2k13XXHMNs2fPrvr9tw11ySWXcNttt1FRUcH7779Pz549ufLKK9ljjz02uO5Vq1ax5ZZbsnTpUv72t79x+umn8+CDDzZC1GuzJ06SJElSk7r11ls5/vjj2WSTTbjuuus466yzSCnx5ptvctNNN3HVVVc16vrOP/98ZsyYwSuvvMIJJ5zAgAEDWLJkSaOuY5tttmHzzTdn6tSpjVovmMRJkiRJamJ33HEHQ4YMAWDQoEF06tSJX/3qV3z7299m9OjRdOzYEYCrr76aPn360LNnTy6//HIAli1bxpFHHkmvXr3o3r0799xzT4PWPWLECAYMGMBdd90FwFNPPUX//v3Zd999OfLII1m4cCEAN954I/vttx+9evVi6NChfPjhhwC89tpr9O3bl/3224/Ro0evVfexxx7LHXfcsd7tUhuTOEmSJElNZuXKlcyePZsuXbpUTbvuuuu4+OKLWbRoESeffDIAkyZNYu7cuUydOpUZM2bw2GOP8dhjjzFp0iS6dOnCs88+ywsvvMBhhx3W4Bj22WcfZs2axUcffcS3vvUt7r33Xp5++mlOOukkvve97wEwdOhQnnrqKZ599ll23nlnxowZA8C5557Lt771LZ566ikqKirWqreyspI///nP69cwdfCaOEmfSyde6Onp8+75pg5AkvRPYfHixXTo0GGtadtuuy2HHHIIgwcPrpr28MMP89BDD9G7d28Ali9fzssvv0zfvn0ZNWoUo0aN4phjjmH//fdvcAwpJQBmzpzJiy++yMCBAwFYvXo1nTt3BuC5557j+9//PkuXLmXZsmVVsT3++OM88MADAJx88slceumlVfVutdVWvPXWWw2Opz5+SpIkSZLUZNq2bcuKFSvWmd6iRQtatPhk4GBKiUsuuYTTTz99nbLTpk1j0qRJnH/++QwePJiLLrqoQTFMnz6dAw44gJQSPXv2rLH37JRTTuGhhx6ie/fu3HLLLTzxxBMARAQRUWO9K1asoG3btg2KpRwOp5QkSZLUZDp27Mjq1atrTORKHXHEEfzyl7/k/fffB2D+/PksXryYN998k3bt2nHyySfzne98h2eeeQaACy64oKqHrC4TJkxg8uTJDBs2jD333JM333yTJ598EsiGer744osAvP/++2yzzTb84x//4M4776xavl+/fkyYMAFgnevfXn75Zbp3715mS5TPnjhJkiRJVcr5SYDGdvjhhzNlypSqYYw1Oeqoo5g1axb9+vUDoH379tx555289NJLjBo1ihYtWtC6deuqnyJ47rnnGDp0aI11/ed//idjxozh/fffp0ePHkyePJktttgCgHvuuYdvfvObLFu2jFWrVvHd736Xvfbai8svv5w+ffqwww470L1796qk84YbbuDLX/4y1157Lccdd9xa65k8eTJHH330BrdPdbFm/GdTqqysTNOmTWvqMBpVj7E9mjoE1eH5U72a5/PO19Dnn68jSfrnMHPmTLp1a9zfhmuo6dOnc+2113L77bc3Sn0pJQYNGsTvfve7RqlvfWM48MADefDBB9lss83WmldTm0fE0ymlynLqdjilJEmSpCbVu3dvBgwYwOrVqxulvoho0gQO4O233+aCCy5YJ4FrDA6nlCRJktTkTjvttKYOoVFtvfXWfOlLX/pU6rYnTpIkSWrmPg+XWDUXjdHWJnGSJElSM9amTRuWLFliIvcZSCmxZMkS2rRps0H1OJxSkiRJasY6d+7M/PnzWbRoUVOH0iy0adOm6gfE15dJnCRJktSMbbTRRnTt2rWpw1ADOJxSkiRJkgrEJE6SJEmSCsQkTpIkSZIKxCROkiRJkgrEJE6SJEmSCsQkTpIkSZIKxCROkiRJkgrEJE6SJEmSCsQkTpIkSZIKxCROkiRJkgrEJE6SJEmSCsQkTpIkSZIKxCROkiRJkgrEJE6SJEmSCsQkTpIkSZIKxCROkiRJkgrEJE6SJEmSCsQkTpIkSZIKxCROkiRJkgrEJE6SJEmSCsQkTpIkSZIKxCROkiRJkgrEJE6SJEmSCqTeJC4i2kTEkxHxbES8GBGX5dNHR8SbETEj/zuqZJkLI+LViPhrRBzxaW6AJEmSJDUnrcoo8xFwSEppeURsBEyJiIfyeT9JKf2otHBE7AkMB/YCtgX+EBG7pZRWN2bgkiRJktQc1dsTlzLL86cb5X+pjkWGAHellD5KKc0BXgX6bHCkkiRJkqTyromLiJYRMQN4G/h9SmlqPuuciHguIm6NiI75tO2AeSWLz8+nSZIkSZI2UFlJXEppdUppb6Az0CciugP/DewM7A0sAH6cF4+aqqg+ISLOiIhpETFt0aJF6xW8JEmSJDU3Dbo7ZUppKfAoMCiltDBP7j4GfsEnQybnA9uXLNYZeKuGum5OKVWmlCorKirWK3hJkiRJam7KuTtlRUR0yB+3BQYCsyKiU0mx44AX8se/AYZHxMYR0RXYFXiyccOWJEmSpOapnLtTdgLGRkRLsqRvQkppYkTcHhF7kw2VfB34OkBK6cWImAC8BKwCzvbOlJIkSZLUOOpN4lJKzwG9a5h+ch3LXAlcuWGhSZIkSZKqa9A1cZIkSZKkpmUSJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFUm8SFxFtIuLJiHg2Il6MiMvy6ZtHxO8j4pX8f8eSZS6MiFcj4q8RccSnuQGSJEmS1JyU0xP3EXBISqkXsDcwKCL6AaOA/0sp7Qr8X/6ciNgTGA7sBQwCfh4RLT+N4CVJkiSpuak3iUuZ5fnTjfK/BAwBxubTxwLH5o+HAHellD5KKc0BXgX6NGrUkiRJktRMlXVNXES0jIgZwNvA71NKU4GtU0oLAPL/W+XFtwPmlSw+P58mSZIkSdpAZSVxKaXVKaW9gc5An4joXkfxqKmKdQpFnBER0yJi2qJFi8qLVpIkSZKauQbdnTKltBR4lOxat4UR0Qkg//92Xmw+sH3JYp2Bt2qo6+aUUmVKqbKiomI9QpckSZKk5qecu1NWRESH/HFbYCAwC/gNcGpe7FTg1/nj3wDDI2LjiOgK7Ao82diBS5IkSVJz1KqMMp2AsfkdJlsAE1JKEyPicWBCRJwOzAWGAqSUXoyICcBLwCrg7JTS6k8nfEmSJElqXupN4lJKzwG9a5i+BDi0lmWuBK7c4OgkSZIkSWtp0DVxkiRJkqSmZRInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQVSbxIXEdtHxOSImBkRL0bEt/LpoyPizYiYkf8dVbLMhRHxakT8NSKO+DQ3QJIkSZKak1ZllFkFfDel9ExEtAeejojf5/N+klL6UWnhiNgTGA7sBWwL/CEidksprW7MwCVJkiSpOaq3Jy6ltCCl9Ez+eBkwE9iujkWGAHellD5KKc0BXgX6NEawkiRJktTcNeiauIjoAvQGpuaTzomI5yLi1ojomE/bDphXsth8akj6IuKMiJgWEdMWLVrU4MAlSZIkqTkqO4mLiHbAvcB5KaX3gP8Gdgb2BhYAP15TtIbF0zoTUro5pVSZUqqsqKhocOCSJEmS1ByVlcRFxEZkCdwdKaX7AFJKC1NKq1NKHwO/4JMhk/OB7UsW7wy81XghS5IkSVLzVc7dKQP4JTAzpXRtyfROJcWOA17IH/8GGB4RG0dEV2BX4MnGC1mSJEmSmq9y7k65P3Ay8HxEzMinXQSMiIi9yYZKvg58HSCl9GJETABeIruz5dnemVKSJEmSGke9SVxKaQo1X+c2qY5lrgSu3IC4JEmSJEk1aNDdKSVJkiRJTcskTpIkSZIKxCROkiRJkgrEJE6SJEmSCqScu1NKkiRJzcrMPbo1dQiqR7dZM5s6hCZjT5wkSZIkFYhJnCRJkiQViEmcJEmSJBWISZwkSZIkFYhJnCRJkiQViEmcJEmSJBWISZwkSZIkFYi/EydJ0j8hf+Pq8685/8aVpA1jT5wkSZIkFYhJnCRJkiQViEmcJEmSJBWISZwkSZIkFYhJnCRJkiQViEmcJEmSJBWISZwkSZIkFYhJnCRJkiQViEmcJEmSJBWISZwkSZIkFUirpg7gn9Xzc+Y2dQiSJEmS/gnZEydJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQXiTwx8SrqsuLOpQ1AdXm/qACRJkqT1ZE+cJEmSJBWISZwkSZIkFYhJnCRJkiQViEmcJEmSJBWISZwkSZIkFYhJnCRJkiQViEmcJEmSJBWISZwkSZIkFYhJnCRJkiQViEmcJEmSJBWISZwkSZIkFYhJnCRJkiQViEmcJEmSJBWISZwkSZIkFYhJnCRJkiQViEmcJEmSJBWISZwkSZIkFYhJnCRJkiQViEmcJEmSJBWISZwkSZIkFYhJnCRJkiQVSL1JXERsHxGTI2JmRLwYEd/Kp28eEb+PiFfy/x1LlrkwIl6NiL9GxBGf5gZIkiRJUnNSTk/cKuC7KaVuQD/g7IjYExgF/F9KaVfg//Ln5POGA3sBg4CfR0TLTyN4SZIkSWpu6k3iUkoLUkrP5I+XATOB7YAhwNi82Fjg2PzxEOCulNJHKaU5wKtAn8YOXJIkSZKaowZdExcRXYDewFRg65TSAsgSPWCrvNh2wLySxebn0yRJkiRJG6hVuQUjoh1wL3BeSum9iKi1aA3TUg31nQGcAbDDDjuUG4YkSZL0qTvxwrI/JquJPN/UATShsnriImIjsgTujpTSffnkhRHRKZ/fCXg7nz4f2L5k8c7AW9XrTCndnFKqTClVVlRUrG/8kiRJktSslHN3ygB+CcxMKV1bMus3wKn541OBX5dMHx4RG0dEV2BX4MnGC1mSJEmSmq9y+on3B04Gno+IGfm0i4CrgQkRcTowFxgKkFJ6MSImAC+R3dny7JTS6kaPXJIkSZKaoXqTuJTSFGq+zg3g0FqWuRK4cgPikiRJkiTVoEF3p5QkSZIkNS2TOEmSJEkqEJM4SZIkSSoQkzhJkiRJKhB/xVCSJEmq5vk5c5s6BKlW9sRJkiRJUoHYEyfpc8lvQCVJkmpmT5wkSZIkFYhJnCRJkiQViEmcJEmSJBWISZwkSZIkFYhJnCRJkiQViEmcJEmSJBWIPzEgSdI/oRMv9C3+8+75pg5AUmHZEydJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEidJkiRJBWISJ0mSJEkFYhInSZIkSQViEifp/7d359GWnWWZwJ+XFARCIAwZiBCSGIN2EDuQIoDQDRINiAMBGRIbRKZ0t9BLGl0ubOlmkLgio4CiMsfGEEFFkGYxSINIY0zCICTMEAgxQMIsMwlv/3FOya3LrVRR+9zs+qp+v7Vq3X32mZ5krTep5+69vw0AwECUOAAAgIEocQAAAAPZMncAgI0c9c2z547ATnxi7gAAsI9yJA4AAGAgShwAAMBAlDgAAICBKHEAAAADUeIAAAAGosQBAAAMRIkDAAAYiBIHAAAwECUOAABgIEocAADAQJQ4AACAgShxAAAAA1HiAAAABqLEAQAADESJAwAAGIgSBwAAMBAlDgAAYCBKHAAAwECUOAAAgIEocQAAAANR4gAAAAaixAEAAAxEiQMAABjITktcVb24qi6vqgvX7HtiVf1LVb1n+edea5777ar6aFV9qKrusVnBAQAA9kW7ciTupUnuucH+Z3X38cs/r0uSqjouyalJbr18z/Oqar9VhQUAANjX7bTEdffbknxhFz/v3knO6e5vdffFST6a5MQJ+QAAAFhjyjVxj66q9y5Pt7zxct/Nk3xqzWsuXe4DAABgBXa3xP1xkmOSHJ/k00mesdxfG7y2N/qAqjq9qi6oqguuuOKK3YwBAACwb9mtEtfdn+3uq7r7u0lekO+dMnlpkiPWvPQWSS7bwWc8v7u3dvfWQw45ZHdiAAAA7HN2q8RV1eFrHt4nybaVK1+T5NSq2r+qjk5ybJLzpkUEAABgmy07e0FVvTzJ3ZIcXFWXJnlCkrtV1fFZnCr5iST/OUm6+6KqekWS9ye5MsmjuvuqzYkOAACw79lpievu0zbY/aKref0ZSc6YEgoAAICNTVmdEgAAgGuYEgcAADAQJQ4AAGAgShwAAMBAlDgAAICBKHEAAAADUeIAAAAGosQBAAAMRIkDAAAYiBIHAAAwECUOAABgIEocAADAQJQ4AACAgShxAAAAA1HiAAAABqLEAQAADESJAwAAGIgSBwAAMBAlDgAAYCBKHAAAwECUOAAAgIEocQAAAANR4gAAAAaixAEAAAxEiQMAABiIEgcAADAQJQ4AAGAgShwAAMBAlDgAAICBKHEAAAAD2TJ3AABg9d538SVzRwBgkzgSBwAAMBAlDgAAYCBKHAAAwECUOAAAgIEocQAAAANR4gAAAAaixAEAAAxEiQMAABiIEgcAADAQJQ4AAGAgShwAAMBAlDgAAICBKHEAAAADUeIAAAAGosQBAAAMRIkDAAAYiBIHAAAwECUOAABgIEocAADAQJQ4AACAgShxAAAAA1HiAAAABrLTEldVL66qy6vqwjX7blJVb6qqjyx/3njNc79dVR+tqg9V1T02KzgAAMC+aFeOxL00yT3X7Xtckjd397FJ3rx8nKo6LsmpSW69fM/zqmq/laUFAADYx+20xHX325J8Yd3ueyc5a7l9VpJT1uw/p7u/1d0XJ/lokhNXlBUAAGCft7vXxB3W3Z9OkuXPQ5f7b57kU2ted+ly3/epqtOr6t4lCnwAABTcSURBVIKquuCKK67YzRgAAAD7llUvbFIb7OuNXtjdz+/urd299ZBDDllxDAAAgL3T7pa4z1bV4Umy/Hn5cv+lSY5Y87pbJLls9+MBAACw1u6WuNckechy+yFJXr1m/6lVtX9VHZ3k2CTnTYsIAADANlt29oKqenmSuyU5uKouTfKEJGcmeUVVPTzJJUnunyTdfVFVvSLJ+5NcmeRR3X3VJmUHAADY5+y0xHX3aTt46qQdvP6MJGdMCQUAAMDGVr2wCQAAAJtIiQMAABiIEgcAADAQJQ4AAGAgShwAAMBAlDgAAICBKHEAAAADUeIAAAAGosQBAAAMRIkDAAAYiBIHAAAwECUOAABgIEocAADAQJQ4AACAgShxAAAAA1HiAAAABqLEAQAADESJAwAAGMiWuQMAAKt31DfPnjsCO/GJuQMAw3IkDgAAYCBKHAAAwECUOAAAgIEocQAAAANR4gAAAAaixAEAAAxEiQMAABiIEgcAADAQJQ4AAGAgShwAAMBAlDgAAICBKHEAAAADUeIAAAAGosQBAAAMRIkDAAAYiBIHAAAwECUOAABgIEocAADAQJQ4AACAgShxAAAAA1HiAAAABqLEAQAADESJAwAAGIgSBwAAMBAlDgAAYCBKHAAAwECUOAAAgIEocQAAAANR4gAAAAaixAEAAAxEiQMAABiIEgcAADAQJQ4AAGAgW6a8uao+keRfk1yV5Mru3lpVN0nyF0mOSvKJJA/o7i9OiwkAAECymiNxP9Xdx3f31uXjxyV5c3cfm+TNy8cAAACswGacTnnvJGctt89KcsomfAcAAMA+aWqJ6yRvrKp3VtXpy32Hdfenk2T589CJ3wEAAMDSpGvikty5uy+rqkOTvKmqPrirb1yWvtOT5Ja3vOXEGAAAAPuGSUfiuvuy5c/Lk7wqyYlJPltVhyfJ8uflO3jv87t7a3dvPeSQQ6bEAAAA2GfsdomrqutX1Q22bSc5OcmFSV6T5CHLlz0kyaunhgQAAGBhyumUhyV5VVVt+5yzu/v1VXV+kldU1cOTXJLk/tNjAgAAkEwocd398ST/foP9n09y0pRQAAAAbGwzbjEAAADAJlHiAAAABqLEAQAADESJAwAAGIgSBwAAMBAlDgAAYCBKHAAAwECUOAAAgIEocQAAAANR4gAAAAaixAEAAAxEiQMAABiIEgcAADAQJQ4AAGAgShwAAMBAlDgAAICBKHEAAAADUeIAAAAGosQBAAAMRIkDAAAYiBIHAAAwECUOAABgIFvmDgAAAHuao7559twR2IlPzB1gRo7EAQAADESJAwAAGIgSBwAAMBAlDgAAYCBKHAAAwECUOAAAgIEocQAAAANR4gAAAAaixAEAAAxEiQMAABiIEgcAADAQJQ4AAGAgShwAAMBAlDgAAICBKHEAAAADUeIAAAAGosQBAAAMRIkDAAAYiBIHAAAwECUOAABgIEocAADAQJQ4AACAgShxAAAAA1HiAAAABqLEAQAADESJAwAAGIgSBwAAMBAlDgAAYCBKHAAAwECUOAAAgIFsWomrqntW1Yeq6qNV9bjN+h4AAIB9yaaUuKraL8kfJfnZJMclOa2qjtuM7wIAANiXbNaRuBOTfLS7P97d305yTpJ7b9J3AQAA7DOqu1f/oVX3S3LP7n7E8vGDk9yhux+95jWnJzl9+fBHk3xo5UFYpYOTfG7uEDAwMwTTmCGYzhzt2Y7s7kN25YVbNilAbbBvu7bY3c9P8vxN+n5WrKou6O6tc+eAUZkhmMYMwXTmaO+xWadTXprkiDWPb5Hksk36LgAAgH3GZpW485McW1VHV9V1kpya5DWb9F0AAAD7jE05nbK7r6yqRyd5Q5L9kry4uy/ajO/iGuPUV5jGDME0ZgimM0d7iU1Z2AQAAIDNsWk3+wYAAGD1lDgAAICBKHEAAAADUeLYqaq6cVX9xNw5YFRmCKYxQ7B7qurXq+qGtfCiqnpXVZ08dy6mU+LYUFW9dTn0N0nyz0leUlXPnDsXjMIMwTRmCFbiYd39lSQnJzkkyUOTnDlvJFZBiWNHDloO/X2TvKS7T0jy0zNngpGYIZjGDMF0tfx5ryzm6J/X7GNgShw7sqWqDk/ygCSvnTsMDMgMwTRmCKZ7Z1W9MYsS94aqukGS786ciRXYlJt9s1d4UhY3a397d59fVT+c5CMzZ4KRmCGYxgzBBFVVSf5XFqdRfry7v15VN83ilEoGp8TxfapqvyRHdPe/XUTe3R9P8kvzpYJxmCGYxgzBdN3dVfU3y1ORt+37fJLPzxiLFXE6Jd+nu69K8otz54BRmSGYxgzBypxbVbefOwSrV909dwb2QFV1RpKDkvxFkq9t29/d75otFAzEDME0Zgimq6r3J7lVkk9mMUeVxUE6t+wYnBLHhqrqLRvs7u6++zUeBgZkhmAaMwTTVdWRG+3v7k9e01lYLSUOAAD2YlV1aJLrbnvc3ZfMGIcVsLAJO1RVP5fk1tl+6J88XyIYixmCacwQTFNVv5jkGUl+KMnlSY5M8oEs5oqBWdiEDVXVnyR5YJL/lsX50/fPYvCBXWCGYBozBCvxu0numOTD3X10kpOS/L95I7EKShw78pPd/StJvtjdT0pypyRHzJwJRmKGYBozBNN9Z3lbgWtV1bW6+y1Jjp87FNM5nZId+cby59er6oeyuKfI0TPmgdGYIZjGDMF0X6qqA5P8Q5I/r6rLk1w5cyZWQIljR15bVTdK8rQk70rSSV44byQYihmCacwQTHfvLH4h8pgk/ymL23a4rnQvYHVKdqqq9k9y3e7+8txZYERmCKYxQ7D7lrcZOLa7/66qDkiyX3f/69y5mMY1cWyoqg6oqv9ZVS/o7m8lObSqfn7uXDAKMwTTmCGYrqoemeQvk/zpctfNk/zNfIlYFSWOHXlJkm9lcSF5klya5CnzxYHhmCGYxgzBdI9KcuckX0mS7v5IkkNnTcRKKHHsyDHd/dQk30mS7v5GFks8A7vGDME0Zgim+1Z3f3vbg6raksX1pQxOiWNHvl1V18ty0KvqmCx+IwrsGjME05ghmO7vq+p/JLleVf1Mklcm+duZM7ECFjZhQ1V1cpLfSXJckjdmcSj+V7v7rXPmglGYIZjGDMF0VXWtJA9PcnIWR7LfkOSFrQAMT4ljh6rqpknumMXQn9vdn5s5EgzFDME0ZgimWS4G9Lru/u7cWVgtJY4NVdVfJnlxktcbfPjBmSGYxgzBdFX1siwWB/qrJC/p7g/MHIkVcU0cO/InWdwU8iNVdWZV/djcgWAwZgimMUMwUXc/KMltk3wsyUuq6h+r6vSqusHM0ZjIkTiuVlUdlOS0LK5L+FSSFyR5WXd/Z9ZgMAgzBNOYIZiuqg5O8qAkj0nygSQ/kuQ53f3cWYOx25Q4dmh5LcKDkjw4yWVJ/jzJXZLcprvvNmM0GIIZgmnMEExTVb+Q5GFJjknyv5Oc1d2XV9UBST7Q3UfOGpDdpsSxoar66yQ/lsXAv7S7P73muQu6e+ts4WAAZgimMUMwXVX9WRarUb5tg+dO6u43zxCLFVDi2FBV3b27/+/cOWBUZgimMUMAO6bEAQAADMTqlAAAAANR4gAAAAaixLFTVfVrc2eAkZkhmMYMwQ+mqu65ZvugqnpRVb23qs6uqsPmzMZqbJk7AHuWqnrs+l1Jfruqrpsk3f3Maz4VjMMMwTRmCFbi95K8frn9jCSfTvILSe6b5E+TnDJTLlZEiWO9JyV5XZKLsvgfZ5Lsl+QGsyWCsZghmMYMwWpt7e7jl9vPqqqHzJqGlbA6JdupqlsmeWaSjyV5Und/vao+3t0/PHM0GIIZgmnMEExXVZdmMUeV5FFJjunlX/qr6r3d/RNz5mM618Sxne6+pLvvl+QdSd5UVfebOxOMxAzBNGYIVuIFWRy9PjDJWUkOTpKqulmS98yYixVxJI4dqqoDsjit5Q7d/R/nzgOjMUMwjRkC2JgSBwDs0arqpt39+blzwCiq6g5JPtDdX6mq6yV5XJLbJXl/kt/r7i/PGpDJnE7JdqrqZlX1x1X1R1V106p64nJJ2ldU1eFz54M93bplnW9UVS+0rDPsuqo6s6q2nfq1tao+nuTcqvpkVd115ngwihcn+fpy+9lJDkry+8t9L5krFKujxLHeS7P4Lc2nkrwlyTeS/HySf0jyJ/PFgmH83prtpyf5TBbLOp+fxbLOwNX7ue7+3HL7aUke2N3HJvmZLJZKB3buWt195XJ7a3c/prvf3t1PSmKRoL2AEsd6h3X3c7v7zCQ36u7fX15k/twkR84dDgaztbsf392f7O5nJTlq7kAwgGtX1bZbIF2vu89Pku7+cJL954sFQ7mwqh663P7nqtqaJFV1qyTfmS8Wq+I+cay3ttj/2dU8B2zs0OXNiivJDauq+nsXH5sh2Lk/SvK6qjozyeur6g+S/HWSk2JVPdhVj0jy7Kp6fJLPJfnHqvpUFmdaPWLWZKyEEsd6r66qA7v7q939+G07q+pHknx4xlwwim3LOiffW9b5Css6w67p7udW1fuS/Nckt8ri7yo/muRVSZ4yZzYYxXLhkl+tqhtkcfrkliSXdvdn503GqlidEgAAYCBO7WGnquq1c2eAkZkhmMYMAWxPiWNX3HzuADA4MwTTmCGANZQ4dsW75w4AgzNDMI0Zgh9QVe1XVX83dw42h2viAABgL1RVr0ny4OVCJ+xFHIljO1W1tareUlUvq6ojqupNVfWlqjq/qm47dz7Y01XVPddsH1RVL6qq91bV2VV12JzZYARmCFbqm0net5yj52z7M3copnMkju1U1XlJnpDkRkmemuS/d/dfVtVJSZ7S3XeaNSDs4arqXd19u+X2C5N8JovbDtw3yV27+5Q588GezgzB6lTVQzba391nXdNZWC0lju1U1bu7+7bL7Uu6+5YbPQdsbN1fQN/T3ceveW67x8D3M0OwWlV1vSS37O4PzZ2F1XGzb9b7ZlWdnOSgJF1Vp3T331TVXZNcNXM2GMGhVfXYJJXkhlVV/b3fljmFHXbODMGKVNUvJHl6kuskObqqjk/y5O7+xXmTMZX/GLLef0nyG0keluQeSX6qqr6U5HlJfn3OYDCIFyS5QZIDk5yV5OAkqaqbJXnPjLlgFGYIVueJSU5M8qUk6e73JDl6zkCshtMpAQBgL1RV/9Tdd1h3ucx7u/sn5s7GNI7Escuq6qFzZ4ARVNWJVXX75fZxVfXYqrrX3LlgRFV1l+UMnTx3FhjQhVX1y0n2q6pjq+q5Sd4xdyimcySOXbZ+oRPg+1XVE5L8bBbXHL8pyR2SvDXJTyd5Q3efMV862PNV1XndfeJy+5FJHpXkVUlOTvK33X3mnPlgJFV1QJLfyWJ+kuQNWaw2/s35UrEKShzbqar37uipJLfq7v2vyTwwmqp6X5Ljk+yfxdLot+juryxXB/snp7DA1Vt32tf5Se7V3VdU1fWTnNvdt5k3IYyjqm7b3e+eOwerZ3VK1jssiwVNvrhuf8Xhd9gVV3b3VUm+XlUf6+6vJEl3f6OqvjtzNhjBtarqxllc8lHdfUWSdPfXqurKeaPBcJ5ZVYcneWWSc7r7orkDsRpKHOu9NsmBy9WLtlNVb73m48Bwvl1VB3T315OcsG1nVR2URImDnTsoyTuz+OVhV9XNuvszVXXgch+wi7r7p5Yruz4gyfOr6oZJ/qK7nzJzNCZyOiXAClXV/t39rQ32H5zk8O5+3wyxYHjLa3sO6+6L584CI6qq2yT5rSQP7O7rzJ2HaZQ4AADYC1XVv0vywCT3T/K5JOck+avuvnzWYEymxAEAwF6oqs5N8vIkr+zuy+bOw+oocQAAsBdarox8TJJO8jG3Fth7uNk3AADsRapqS1U9NcklSc5K8rIkn6qqp1bVtedNxyoocQAAsHd5WpKbJPnh7j5hee/FY5LcKMnTZ03GSjidEgAA9iJV9ZEkt+p1f9Gvqv2SfLC7j50nGaviSBwAAOxden2BW+68Kovr4xicEgcAAHuX91fVr6zfWVUPSvLBGfKwYk6nBACAvUhV3TzJXyf5RpJ3ZnH07fZJrpfkPt39LzPGYwWUOAAA2AtV1d2T3DpJJbmou988cyRWRIkDAAAYiGviAAAABqLEAQAADESJA2AIVXWfquqq+rG5swDAnJQ4AEZxWpK3Jzn1mvzS5c1xAWCPocQBsMerqgOT3DnJw7OmxFXVtarqeVV1UVW9tqpeV1X3Wz53QlX9fVW9s6reUFWHb/C5x1TVuVV1flU9uaq+utx/t6p6S1WdneR9y32PraoLl38es9x3VFVduObzfrOqnrjcfmtV/UFVvWP5nhM37V8QAPsUJQ6AEZyS5PXd/eEkX6iq2y333zfJUUluk+QRSe6UJFV17STPTXK/7j4hyYuTnLHB5z47ybO7+/ZJLlv33IlJfqe7j6uqE5I8NMkdktwxySOr6ra7kPv63f2TSX5tmQEAJlPiABjBaUnOWW6fs3ycJHdJ8sru/m53fybJW5b7fzTJjyd5U1W9J8njk9xig8+9U5JXLrfPXvfced198ZrveVV3f627v5rFTXT/wy7kfnmSdPfbktywqm60C+8BgKu1Ze4AAHB1quqmSe6e5MerqpPsl6Sr6reyuIHthm/L4sa2d5rw1V9b93kbuTLb/0L0uuueX38zVjdnBWAyR+IA2NPdL8mfdfeR3X1Udx+R5OIsjo69PckvLa+NOyzJ3Zbv+VCSQ6rq306vrKpbb/DZ5yb5peX21S2Y8rYkp1TVAVV1/ST3SfIPST6b5NCqumlV7Z/k59e974HL779Lki9395d/oH9yANiAI3EA7OlOS3Lmun1/leSXkzwqyUlJLkzy4ST/lEVZ+vZygZPnVNVBWfz/7g+SXLTucx6T5GVV9RtJ/k+SDUtWd7+rql6a5Lzlrhd297uTpKqevPzei5N8cN1bv1hV70hywyQP+0H+oQFgR6rbmR0AjKuqDuzury5PuzwvyZ2X18ftynsPSPKN7u6qOjXJad197xXlemuS3+zuC1bxeQCwjSNxAIzutcsFQ66T5Hd3tcAtnZDkD6uqknwpjpYBMABH4gAAAAZiYRMAAICBKHEAAAADUeIAAAAGosQBAAAMRIkDAAAYiBIHAAAwkP8PK3JKkUPrKpcAAAAASUVORK5CYII=\n", "text/plain": [ - "array([[221],\n", - " [ 6],\n", - " [182],\n", - " [ 7],\n", - " [172],\n", - " [ 19],\n", - " [190],\n", - " [ 39],\n", - " [ 81],\n", - " [ 40],\n", - " [ 64],\n", - " [ 51],\n", - " [ 28],\n", - " [165],\n", - " [ 7],\n", - " [ 42]])" + "
" ] }, - "execution_count": 56, - "metadata": {}, - "output_type": "execute_result" + "metadata": { + "needs_background": "light" + }, + "output_type": "display_data" + } + ], + "source": [ + "_, ax = plt.subplots(figsize=(15, 10))\n", + "data.pivot_table(index=['Age group'], columns=['Smoker', 'Status'], aggfunc='size').plot(kind='bar', stacked=True, ax=ax)\n", + "ax.set_title(\"Number of alive/dead women after 20 years, according to their smoking habits and age\")\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "One can see that most women that were over 65 in 1972 are dead twenty years after and that, at that a great majority of polled older women were non-smokers. For the other age groups results look similar but are difficult to interpret.\n", + "\n", + "Let's therefore compute the mortality rates for smokers and non-smokers in different age groups." + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Age group: 18-34 years\n", + "\tMortality rate for smokers:\t3.70% ± 2.69%\n", + "\tMortality rate for non smokers:\t2.64% ± 2.09%\n", + "\n", + "Age group: 35-54 years\n", + "\tMortality rate for smokers:\t17.03% ± 4.87%\n", + "\tMortality rate for non smokers:\t9.95% ± 4.24%\n", + "\n", + "Age group: 55-64 years\n", + "\tMortality rate for smokers:\t44.35% ± 9.08%\n", + "\tMortality rate for non smokers:\t33.06% ± 8.38%\n", + "\n", + "Age group: Over 65 years\n", + "\tMortality rate for smokers:\t85.71% ± 9.80%\n", + "\tMortality rate for non smokers:\t85.49% ± 4.97%\n", + "\n" + ] } ], "source": [ - "np.array(data.groupby(['Binned age', 'Smoker', 'Status']).count())" + "def ci_per_age(age_group):\n", + " pop = data[data['Age group'] == group]\n", + " print(f\"Age group: {group}\")\n", + " \n", + " num_smokers = sum(pop['Smoker'] == \"Yes\")\n", + " num_dead_smokers = sum(np.logical_and(pop['Smoker'] == \"Yes\", pop['Status'] == \"Dead\"))\n", + " rate_smokers = num_dead_smokers / num_smokers\n", + " ci_smokers = z * (rate_smokers * (1 - rate_smokers) / num_smokers) ** 0.5\n", + " print(f\"\\tMortality rate for smokers:\\t{rate_smokers:.2%} \" + u\"\\u00B1\" + f\" {ci_smokers:.2%}\")\n", + "\n", + " num_non_smokers = len(pop) - num_smokers\n", + " num_dead_non_smokers = sum(np.logical_and(pop['Smoker'] == \"No\", pop['Status'] == \"Dead\"))\n", + " rate_non_smokers = num_dead_non_smokers / num_non_smokers\n", + " ci_non_smokers = z * (rate_non_smokers * (1 - rate_non_smokers) / num_non_smokers) ** 0.5\n", + " print(f\"\\tMortality rate for non smokers:\\t{rate_non_smokers:.2%} \" + u\"\\u00B1\" + f\" {ci_non_smokers:.2%}\")\n", + " \n", + "for group in sorted(data['Age group'].unique()):\n", + " ci_per_age(group)\n", + " print()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now, the mortality rate is considerably higher for smokers, especially for people between 35 and 65 years old. This might seem like a contradiction, as before the rate was higher for non-smokers. However, from the previous bar chart it is clear that the percentage of polled smokers/non-smokers is different in different age groups, in particular for older women as mentioned. In addition, the fact that most polled women over 65 are non-smokers can be an argument in favor of the hypothesis that smoking is dangerous for health, but this can't be proven through statistics.\n", + "\n", + "---\n", + "\n", + "The age groups are fixed a-priori. In order to have more flexible results and reduce the introduced bias it is possible to try to perform a logistic regression, studying the probability of death in the two groups according to the age." ] }, { "cell_type": "code", - "execution_count": 57, + "execution_count": 10, "metadata": {}, "outputs": [ { @@ -423,133 +491,138 @@ " \n", " \n", " \n", - " \n", - " \n", - " Age\n", - " \n", - " \n", - " Binned age\n", " Smoker\n", " Status\n", - " \n", + " Age\n", + " Age group\n", + " Death\n", " \n", " \n", " \n", " \n", - " 18-34 years\n", - " No\n", - " Alive\n", - " 221\n", - " \n", - " \n", - " Dead\n", - " 6\n", - " \n", - " \n", - " Yes\n", - " Alive\n", - " 182\n", - " \n", - " \n", - " Dead\n", - " 7\n", - " \n", - " \n", - " 35-54 years\n", - " No\n", - " Alive\n", - " 172\n", - " \n", - " \n", - " Dead\n", - " 19\n", - " \n", - " \n", - " Yes\n", - " Alive\n", - " 190\n", - " \n", - " \n", - " Dead\n", - " 39\n", - " \n", - " \n", - " 55-64 years\n", - " No\n", - " Alive\n", - " 81\n", - " \n", - " \n", - " Dead\n", - " 40\n", - " \n", - " \n", - " Yes\n", - " Alive\n", - " 64\n", - " \n", - " \n", - " Dead\n", - " 51\n", + " 0\n", + " Yes\n", + " Alive\n", + " 21.0\n", + " 18-34 years\n", + " 0\n", " \n", " \n", - " Over 65 years\n", - " No\n", - " Alive\n", - " 28\n", + " 1\n", + " Yes\n", + " Alive\n", + " 19.3\n", + " 18-34 years\n", + " 0\n", " \n", " \n", - " Dead\n", - " 165\n", + " 2\n", + " No\n", + " Dead\n", + " 57.5\n", + " 55-64 years\n", + " 1\n", " \n", " \n", - " Yes\n", - " Alive\n", - " 7\n", + " 3\n", + " No\n", + " Alive\n", + " 47.1\n", + " 35-54 years\n", + " 0\n", " \n", " \n", - " Dead\n", - " 42\n", + " 4\n", + " Yes\n", + " Alive\n", + " 81.4\n", + " Over 65 years\n", + " 0\n", " \n", " \n", "\n", "" ], "text/plain": [ - " Age\n", - "Binned age Smoker Status \n", - "18-34 years No Alive 221\n", - " Dead 6\n", - " Yes Alive 182\n", - " Dead 7\n", - "35-54 years No Alive 172\n", - " Dead 19\n", - " Yes Alive 190\n", - " Dead 39\n", - "55-64 years No Alive 81\n", - " Dead 40\n", - " Yes Alive 64\n", - " Dead 51\n", - "Over 65 years No Alive 28\n", - " Dead 165\n", - " Yes Alive 7\n", - " Dead 42" + " Smoker Status Age Age group Death\n", + "0 Yes Alive 21.0 18-34 years 0\n", + "1 Yes Alive 19.3 18-34 years 0\n", + "2 No Dead 57.5 55-64 years 1\n", + "3 No Alive 47.1 35-54 years 0\n", + "4 Yes Alive 81.4 Over 65 years 0" ] }, - "execution_count": 57, + "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ - "data.groupby(['Binned age', 'Smoker', 'Status']).count()" + "death = []\n", + "for status in data['Status']:\n", + " if status == \"Alive\":\n", + " death.append(0)\n", + " elif status == \"Dead\":\n", + " death.append(1)\n", + "assert len(death) == len(data['Age'])\n", + "data['Death'] = death\n", + "data.head()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Barchart with one bar for each age, divided in four parts(smoker or not, dead or alive)" + "It will be done using statsmodels: sklearn does not allow to derive prediction confidence intervals and uses regularization by default." + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": {}, + "outputs": [], + "source": [ + "smokers = data[data['Smoker'] == \"Yes\"]\n", + "non_smokers = data[data['Smoker'] == \"No\"]\n", + "log_smokers = LogisticRegression().fit(smokers['Age'].values.reshape([-1, 1]), smokers['Death'])\n", + "log_non_smokers = LogisticRegression().fit(non_smokers['Age'].values.reshape([-1, 1]), non_smokers['Death'])" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "At this point, we can plot the predicted probabilities of death as a function of the age." + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlYAAAF1CAYAAAAqdaQaAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzs3Xd0VNXax/HvTiGBAKGDEHqT3kJTmkqzgKD0DipivdjvtRe8Yq+vF1Q6UgQRaQp2UJCO9F5DTYAEkhBI2e8fZ4KTmApJJiS/z1pZmZmzZ59nTptn9tlnH2OtRURERESunpenAxARERHJK5RYiYiIiGQRJVYiIiIiWUSJlYiIiEgWUWIlIiIikkWUWImIiIhkESVW+YQxpq0xZpen40iJMaaDMSYkh+ZVxRhjjTE+V/h+a4ypkcq0gcaYZSmVNcaMM8a8cGVRZzrGB4wxJ40xkcaYkhkof9AY0zGbYsm2uuXKGWOGGWN+d3seaYyp5smYPCWt4092HjeNMb8aY+69wvdONsaMSWN6vl2fuYESq1zO9cV0wbWjJP59koH3JUkArLUrrLW1synGNHfy/MJa+6W1tnMq00ZZa1+D7E0kjTG+wHtAZ2ttYWvt6eyYTyrz1nZwjXJtK/uzul5jzMvGmOlZXW9Oyc7jZnZyX5/aL3PeFf1qlxzXzVr7o6eDuBYYY3ystXGejsODygL+wDZPB+JJ2g7+ZozxttbGezoOkfxCLVbXMGNMDWPMb8aYCGNMmDFmtuv15a4if7lauPombyVxtYQ9ZYzZbIyJMsZMMMaUNcZ8Z4w5b4z50RhT3K38HGPMCde8lhtj6rleHwkMBJ52zWuh6/XyxpivjTGhxpgDxphH3eoq6PoVddYYsx1ons7ntMaYR40x+12f821jjJdr2jBjzB/GmPeNMWeAl40xXsaY540xh4wxp4wxU40xgcmqHWGMOWaMOW6MecJtXi2MMauMMeGuaZ8YYwoke+9tacTyOylI/NVojAkAvgPKu7VAljfGRLuftjPGNHMtO98U6vIzxnzgiv+Y67GfMaYWkHjaItwY83MqsQx2LZvTxpjnkk3zMsb82xizzzX9K2NMCbfpmdoOXBq7trMIY8xsY4x/KnElrsuPXWV3GmNucZse6NpOjxtjjrqWp3ey917eDlKoP811a4zpbIzZ5Zr3p8bZt+51mz7CGLPDtd0uNcZUTulzpDDf6saYn13LM8wY86Uxppjb9IrGmHmu9X3auLVIG2Puc83zvDFmuzGmqev1OsY5lRRujNlmjOnu9p7Jxpj/GWOWGGOigJuMMSWNMQuMMeeMMWuA6slidD9tPdkY83/GmMWu+a42xlR3K5vmcnIr1xV4Fujr2ib+cr1e3hXLGWPMXmPMfWksu9tcn/u8a50/6Xq9gzEmxBjztHH28ePGmB6u8rtddT/rVk+K+0wq83zUNc8gk/Jx88nUtmdXPMdd87jXpNF1wKWya7s9b4xZZowp5VZXivuam1LGmB9c7/3NfXtMnK9J/fj8jGt5nnety1uQrGOt1V8u/gMOAh1TmTYTeA4nQfYH2rhNs0ANt+cdgJBk9f6J08JRATgFbACaAH7Az8BLbuVHAEVc0z4ANrlNmwyMcXvuBawHXgQKANWA/UAX1/SxwAqgBFAR2OoeWwqf0wK/uMpXAnYD97qmDQPigEdwWmALumLd65pvYWAeMM1VvoqrvplAANAACE1cxkAzoJWrrirADmB0JmL5PaV14L6Mkq8L12tLgAfcnr8PfJzK8njVte7KAKWBlcBryT6fTyrvrQtEAu1c6/I91/JL/PyjXXUHuaaPB2ZeyXbgtp2tAcq7ltkOYFQqsSWuy8cAX6AvEAGUcE2f74onwPXZ1wD3p7YdpFB/qusWKAWcA+5yTf8XEOu2bnvgbFN1XNOfB1ZmcB+uAXRyLbPSwHLgA9c0b+Av1/oOwG0/BnoDR3F+eBhXPZVdy2YvTtJSALgZOA/UdlsPEcCN/H1smAV85ZpHfVe9aW2rZ4AWrs/6JTArI8sphc/+MjA92Wu/AZ+64mqMs//dksr7jwNtXY+LA03d9qE4nGOML3Cfq54ZONtnPSAGqJaBfaYDrv0ReAHnOFg6jeNmitsz0BU44Zp3IWAayY7DyT7br8A+oBbOcetXYGwm9rXz/L0ff5jO+nQ/PtcGjgDl3Y4Z1TOyLesvY38eD0B/6awgZ0eOBMLd/u5zTZsKfAYEpfC+jCRWA92efw38z+35I8D8VGIq5qo/0PU8+Y7bEjic7D3/ASa5Hu8HurpNG0n6iZV7+QeBn1yPh6Uwr5+AB92e18Y5+Cd+oVrgerfpbwETUpn3aOCbTMRypYlVX+AP12NvnAN0i1Ri2gfc5va8C3DQ9Tjx86WWWL2I60vS9TwAuMTfidUO3L7kgOsSl11mtwO37WxQsmU9LpXYhgHHAOP22hpgMM4PgIu4JUxAf+CX1LaDDOxbl9ctMARY5TbN4Hz5JCZW3wH3uE33AqKBypmZp+u9PYCNrsetcRKClJbvUuBfKbze1rV9eLm9NhN42W09THWb5u1ah+7b/H/T2Va/cJt2G7AzI8sphVhfxi2xwvkhFQ8UcXvtDWByKu8/DNwPFE32egfgAuDtel7E9RlaupVZD/TIwD7TASfRfA/4Hdf2nNK+ShrbMzAReMNtWg3ST6yed3v+IPB9KmVT2tfc9+PCruVaMZX1OSZZXKeAjoBvZrdf/aX/p1OB14Ye1tpibn+fu15/GufAtsZ1OmBEJus96fb4QgrPC4PTR8MYM9Y4p4fO4RxcwPn1mpLKOKe6whP/cH5dl3VNL49zME50KAOxJi9fPpVpifW713kIJ6kq6/ZaivUZY2oZYxa5muDP4XwBJf+cacVypb4F6hrnSp5OQIS1dk0qZVP6fBmNIcmyt9ZGAe4d3CsD37ittx04B+yyV7AdJDrh9jga13aViqPWdfR3SfxsiS01x91iG4/TApEo+XaQRDrrNvlysYD7BQaVgQ/d5n0GZ9+rkNY8XfMtY4yZ5Tr1cg6Y7jbfisAhm3J/sIo4CUFy5YEj1toEt9cOJYvFfVmUxtn+M7PPpbbO0ltO6SkPnLHWnk8WS2rL8W6cxO6Q63RXa7dpp+3ffccuuP6neAwj/X2mGM4PvDestRHpfIYMLRvS2R7TqiuD+5r7eojE2SbTPQ5Ya/fi/Kh4GTjl2jaz4hgmLkqsrmHW2hPW2vusteVxftV9ms75/Cs1ALgT5xdOIE6rCDhfLOD8OnJ3BDiQLBksYq29zTX9OM6XRqJKGYghefljbs+Tz/8Yzhehe/k4kh50U6vvf8BOoKa1tihOQmhIKq1YMiJ5vFhrY3BO1QzEaaGZlsb7U/p8GY0hybI3xhQC3IdkOALcmmzd+Vtrj5L57eBKVDDGuC/vxM92BKfFqpRbXEWtte79TtKbf1rr9jjO6U8AXDEEub33CM5pR/flUtBauzIDn+kNV2wNXfMd5DbfI0Alk/LwH0dI1hfK5RhQ0bj69rlUwml1SeS+LEJxtv/M7nMpSW85JZfSvlnCGFMkWSxHSYG1dq219k6cBHo+zj5yJdLbZ84CdwCTjDE3XuE8kiwbki7vzEpvX0tSvzGmMM6pyZSOAykdb2ZYa9vgLBMLvHkVsUoySqyuYcaY3saYxB35LM4OkvgL7iROH6OsUATnS+00Tt+B/yabnnxea4Bzrg6SBV2/vuobYxI7qX8F/McYU9wV/yMZiOEpV/mKOP06ZqdRdibwmDGmquuA819gdrJWgReMMYVcHUKHu9VXBKcPSaQx5nrggauMJSUngZLmnx3qp+Kc0uqO06qRmpnA88aY0q7Ori+mU97dXOAOY0wb43TcfpWkx4FxwOuJHWFd87jTNS2z28GVKAM8aozxNcb0xunTtMRaexxYBrxrjClqnE721Y0x7TNRd1rrdjHQwDgdoH2Ah4BybtPH4WyziZ31A13x4Xr+qzHm5TTmG4lzQUEF4Cm3aWtwvpDHGmMCjDH+bl/sXwBPGudCBuPqjFwZWA1E4XRI9jXGdAC64fSj+gdXq848nAs7Chlj6gJD01pQaUhvOSV3EqiSmARaa4/g9G96w/VZGwL34PTjSsIYU8A4Y8MFWmtjcdbdlV7dmO4+Y639FeeHzTfGmJZXMI+vgOHGubCgkGseVyq9fQ2ci2gS9+PXgNWu5Ztckv3SGFPbGHOzcTrvx+C07Omq0SykxOrasNAkHcfqG9frzYHVxphIYAFOf4wDrmkvA1Ncpy76XOX8p+I0nR8FtuN0AnU3Aec0VrgxZr7rQN4Np2PqASAM50siMZF4xVXfAZwvy7RaZxJ9i9NnYhPOwX1CGmUnuupc7ppHDP9M3n7D6QD8E/COtTZxYM8ncX4tngc+J+WkKTOx/IO1difOgX6/a5mVd73+B5AAbLDWHkyjijHAOmAzsAWns22Gxqmx1m7D+TKcgfOFfpakp3I+xNmWlhljzuOs68QvmUxtBxmJJwWrgZo428zrQC/791hcQ3A6a293xT0Xpw9YRqW6bq21YTidxd/C+TKri7OML7qmf4Pzq36W69TMVuBWt7orAn+kMt9XgKY4HcoX4yQ5ifNN3Fdq4PQnCsHpb4e1do5rGcxwxTwfpyP/JZzk+1ac5fQpMMS1XaXmYZzTTCdw+txMSqNsqtJbTimY4/p/2hizwfW4P04LzDHgG5yLZH5I5f2DgYOuZT4Kp7XvSmRon3HFMRxYYIxplpkZWGu/Az7CubhlL7DKNSm1ZZOW9PY1cLaLl3BOATbDSQpTkny/9MO5gCgMZ3sog9N6K1nEJO3OIJL7GGMszumbvZ6OJbsZZ4iEGdbaLzwdS04zxgzD6QTdJhfE4oWT5Ay01v6STtkgYI61tnVa5fKizCyn/MYYUwcnAfdLpQ+d5FFqsRLJJVynSpuS+VOLkgWMMV2MMcVcp0gS+1+l1FKQhLU2JD8lVVe6nPIDY0xP1ynM4jgtnAuVVOU/SqxEcgFjzBTgR5xxlc6nV16yRWucq/DCcE7P9bDWXkj7LfmSllPq7se5WGAfTr+llPpoSh6nU4EiIiIiWUQtViIiIiJZRImViIiISBZJaVC6JIwxE3EGTjtlra2fwnSDc4n2bTgjxw6z1m5IXi65UqVK2SpVqmQ6YBEREZGctn79+jBrben0yqWbWOGMefIJzrgaKbkVZ9yZmjjj3fyPv8e9SVWVKlVYt25dBmYvIiIi4lnGmIzcfi39U4HW2uU4A5Cl5k6cG35aa+2fQDFjTGYG7RMRERHJE7Kij1UFkt5sMoRUbqhpjBlpjFlnjFkXGhqaBbMWERERyT2yIrFKfoNaSOVmqNbaz6y1wdba4NKl0z1NKSIiInJNyUgfq/SEkPQu3kGkfIftdMXGxhISEkJMTEwWhCVZwd/fn6CgIHx9fT0dioiISK6XFYnVAuBhY8wsnE7rEa470WdaSEgIRYoUoUqVKjgXG4onWWs5ffo0ISEhVK1a1dPhiIiI5HoZGW5hJtABKGWMCcG5m7YvgLV2HLAEZ6iFvTjDLQy/0mBiYmKUVOUixhhKliyJ+sOJiIhkTLqJlbW2fzrTLfBQVgWkpCp30foQERHJOI28noLXX3+devXq0bBhQxo3bszq1auvqr5ff/2VO+64I4uiExERkdwqK/pY5SmrVq1i0aJFbNiwAT8/P8LCwrh06ZLH4omLi8PHR6tJRETkWqAWq2SOHz9OqVKl8PPzA6BUqVKUL1+eKlWq8Oyzz9K6dWuCg4PZsGEDXbp0oXr16owbNw5wOns/9dRT1K9fnwYNGjB79ux/1L927VqaNGnC/v37iYqKYsSIETRv3pwmTZrw7bffAjB58mR69+5Nt27d6Ny5M8ePH6ddu3Y0btyY+vXrs2LFipxbICIiIpJhubYp5JWF29h+7FyW1lm3fFFe6lYvzTKdO3fm1VdfpVatWnTs2JG+ffvSvn17ACpWrMiqVat47LHHGDZsGH/88QcxMTHUq1ePUaNGMW/ePDZt2sRff/1FWFgYzZs3p127dpfrXrlyJY888gjffvstlSpV4tlnn+Xmm29m4sSJhIeH06JFCzp27Ag4LWebN2+mRIkSvPvuu3Tp0oXnnnuO+Ph4oqOjs3S5iIiISNbItYmVpxQuXJj169ezYsUKfvnlF/r27cvYsWMB6N69OwANGjQgMjKSIkWKUKRIEfz9/QkPD+f333+nf//+eHt7U7ZsWdq3b8/atWspWrQoO3bsYOTIkSxbtozy5csDsGzZMhYsWMA777wDOFdFHj58GIBOnTpRokQJAJo3b86IESOIjY2lR48eNG7cOKcXi4iISO4TfQZO7QD/QChX39PRALk4sUqvZSk7eXt706FDBzp06ECDBg2YMmUKwOXTg15eXpcfJz6Pi4vDuUAyZddddx0xMTFs3LjxcmJlreXrr7+mdu3aScquXr2agICAy8/btWvH8uXLWbx4MYMHD+app55iyJAhWfZ5RUREcrXYGAjdASe3wYmtcHIrhO6EKNdwQE2HQvePPBujS65NrDxl165deHl5UbNmTQA2bdpE5cqV2bJlS7rvbdeuHePHj2fo0KGcOXOG5cuX8/bbb7Nz506KFSvGhAkT6Ny5MwEBAXTo0IEuXbrw8ccf8/HHH2OMYePGjTRp0uQf9R46dIgKFSpw3333ERUVxYYNG5RYiYhI3mMtnD/hJE4nt/6dRIXtARvvlPEtBGXqQK0uULoOlL4+17RWgRKrf4iMjOSRRx4hPDwcHx8fatSowWeffcaiRYvSfW/Pnj1ZtWoVjRo1whjDW2+9Rbly5di5cycAZcuWZeHChdx6661MnDiRF154gdGjR9OwYUOstVSpUiXF+fz666+8/fbb+Pr6UrhwYaZOnZrln1tERCTHRZ6Coxvg6Hrn7/hfEB329/TAilC2PtTp5vwv1wCKVwEvb4+FnB6T1umr7BQcHGzXrVuX5LUdO3ZQp04dj8QjqdN6ERGRq3YxEo5v+juJOroRIpx+xRgvKFMXyjeGsg2cFqiy9aBgcc/G7MYYs95aG5xeObVYiYiISNYLPwyHVsKhPyBkndMnyiY404pVhqBgaDkSKjSD6xpBgYC067tGKLESERGRq2MtnN7rJFGHVjp/EUecaX6BULE51OnuJFEVmkJAKc/Gm42UWImIiEjmWOtcoXfoj7+TqcQr9AJKQ+Ub4IZHnP9l6ubqPlFZTYmViIiIpO/ccdj/C+z7Bfb/ClGnnNcDK0H1W5wkqvKNULI6GOPRUD1JiZWIiIj8U2wMHPwd9v3kJFOhO5zXA0pDtQ5Q7Sao2haKVfJYiBHRsaw/fIYSAX40rljMY3G4U2IlIiIijnPHYc8y2L3UaZ2KjQYff6jUGhr3h+o3Q5l64OWZWw2fOh/D6v1nWH3gNGsPnGXXyfMA3NW0Ao0r5o67kiixSsYYw+OPP867774LwDvvvENkZCQvv/yyZwNLRYcOHXjnnXcIDk73ClAREZGkrHXGjtq1BHZ/7zwGZ/yoxgOdQTirtAHfgh4JLyzyIqv3n2HV/jD+3H+GvaciASjs50PTysW5o+F1BFcpkWtaq0CJ1T/4+fkxb948/vOf/1CqVN69agEgPj4eb+/806FQRESAhAQIWQs7Fjh/4YedcaSCWsAtL0Gtrs7I5h7oJ3Um6hKr959m1f7T/Ln/NLtPOolUQAFvmlctQa9mQbSuVpJ65Yvi4+2ZVrP0KLFKxsfHh5EjR/L+++/z+uuvJ5l26NAhRowYQWhoKKVLl2bSpElUqlSJYcOGUbRoUdatW8eJEyd466236NWr1z/qnjNnDq+88gre3t4EBgayfPlyJk+ezPz584mPj2fr1q088cQTXLp0iWnTpuHn58eSJUsoUaIEmzZtYtSoUURHR1O9enUmTpxI8eJ/D5yWkJDA8OHDqVixImPGjGHZsmW89NJLXLx4kerVqzNp0iQKFy5MlSpVGDFiBMuWLePhhx/m1KlTjBs3Dh8fH+rWrcusWbOyfRmLiEgOS4h3+kvtWAA7FkHkCfDyheo3QbunofZtEFAyx8OKjU9gw6Gz/LY7lN92h7Lt2DkAChXwJrhKCXo0qUDraiWpXyEQ31yaSCWXexOr7/4NJ9K/P1+mlGsAt45Nt9hDDz1Ew4YNefrpp5O8/vDDDzNkyBCGDh3KxIkTefTRR5k/fz4Ax48f5/fff2fnzp107949xcTq1VdfZenSpVSoUIHw8PDLr2/dupWNGzcSExNDjRo1ePPNN9m4cSOPPfYYU6dOZfTo0QwZMoSPP/6Y9u3b8+KLL/LKK6/wwQcfABAXF8fAgQOpX78+zz33HGFhYYwZM4Yff/yRgIAA3nzzTd577z1efPFFAPz9/fn9998BKF++PAcOHMDPzy9JTCIico2z1rldzJY5sPVr5yo+n4JQsyPUuRNqdQb/wBwP68iZaJbvCeW3XaGs3HeayItx+HgZmlYuzpOda9G6eikaBl07iVRyuTex8qCiRYsyZMgQPvroIwoW/Pu88qpVq5g3bx4AgwcPTpJ49ejRAy8vL+rWrcvJkydTrPfGG29k2LBh9OnTh7vuuuvy6zfddBNFihShSJEiBAYG0q1bNwAaNGjA5s2biYiIIDw8nPbt2wMwdOhQevfuffn9999/P3369OG5554D4M8//2T79u3ceOONAFy6dInWrVtfLt+3b9/Ljxs2bMjAgQPp0aMHPXr0uLIFJiIiuUfobtg610mozuwH7wJOX6n6vaBmpxwf4TwmNp7VB87w265Qftt9in2hUQBUKFaQ7o3L075WaW6oXpIi/r45Gld2yb2JVQZalrLT6NGjadq0KcOHD0+1jHE7/+zn53f5ceL9F5977jkWL14MwKZNmxg3bhyrV69m8eLFNG7cmE2bNv3jvV5eXpefe3l5ERcXl26sN9xwA7/88gtPPPEE/v7+WGvp1KkTM2fOTLF8QMDfO9XixYtZvnw5CxYs4LXXXmPbtm34+OTezUJERFIQfQa2zIVN010d0A1UbQdtHnduYFwwZzt3n426xM87T7Fs+wmW7w7jQmw8fj5etKxWkgEtK9O+Vmmqlw5I8j2aV+gbNBUlSpSgT58+TJgwgREjRgBOAjNr1iwGDx7Ml19+SZs2bdKs4/XXX0/ST2vfvn20bNmSli1bsnDhQo4cOZKhWAIDAylevDgrVqygbdu2TJs27XLrFcA999zD8uXL6d27N9988w2tWrXioYceYu/evdSoUYPo6GhCQkKoVatWknoTEhI4cuQIN910E23atGHGjBlERkZSrFjuubpCRERSkRDvDNS5cTrsXATxl6BcQ+jyX6h3FxS9LkfDOXImmh+2n2TZ9hOsPXiW+ARLuaL+9GoWxC11ytCqWkn8ffP+BVNKrNLwxBNP8Mknn1x+/tFHHzFixAjefvvty53XM+Opp55iz549WGu55ZZbaNSo0eVWq/RMmTLlcuf1atWq/WPejz/+OBEREZeTvsmTJ9O/f38uXrwIwJgxY/6RWMXHxzNo0CAiIiKw1vLYY48pqRIRye3OHoJNX8KmGc79+AoWh+ARzvAI1zXMsTCstWw7do5l20/yw/aT7DjudDyvXbYID7SvTud6ZWlQITBPtkqlxSSetsppwcHBdt26dUle27FjB3Xq1PFIPJI6rRcREQ9LSIC9P8Laz2HPD85r1W+GJoOcK/p8/XMkDGst24+fY+Ffx1m0+RghZy9gDDSvXIJOdcvSqW5ZqpTK2T5cOcUYs95am+6gkWqxEhERya2izzin+tZNgLMHIaAMtHsKmg6BYhVzLIy9pyJZ+NcxFm4+xv7QKHy8DG1qluLRm2tyc50ylCrsl34l+YQSKxERkdzm+F+w+jPn6r64GKh0A9zyIlzfDXwK5EgIR85Es3DzMRb+dZwdx89hDLSqWpJ721Sja/1ylAjImTiuNUqsREREcoPE030rP4KDK8A3ABr1h+b3Qrn6ORJCWORFFmxyWqY2HnbGNmxaqRgvdavL7Q2uo0zRnDnleC3LdYmVtTbfdXTLzTzVB09EJN+IjYEtX8HKTyBsFxStAJ1eg2ZDc2QAz0txCfy88yRz1x/l112niEuw1L2uKP++9Xpub3AdFUsUyvYY8pJclVj5+/tz+vRpSpYsqeQqF7DWcvr0afz99QtFRCTLRZ9x+k6t/swZFb1cA7jrc6jXE7yzd7DMxCv65q4P4dtNRzkbHUuZIn7c07YqvZoGUbNskWydf16WqxKroKAgQkJCCA0N9XQo4uLv709QUJCnwxARyTuiwmDVJ7Dmc7gUCTU6wg2PQNX22X7j44joWL7ZGMKstUfYeeI8BXy86Fy3LL2aBdGmRqlce2Pja0muSqx8fX2pWrWqp8MQERHJeudPOv2n1k2E2AtOy1S7J6FsvWydrbWW9YfOMmPNYRZvPs7FuAQaBQUypkd9ujUsT2ChvHErmdwiVyVWIiIiec65Y/DHh7B+sjM6eoPe0PZJKF0r3bdejfDoS8zbcJSZaw6z51Qkhf186B0cRL/mlahfIedvvpxfKLESERHJDpGhsOIdp4XKJkDDftD2cShZPVtnu+P4OaasPMj8TUeJiU2gUcVivHV3Q25veB0Bfvraz25awiIiIlkpJsK5wm/V/zljUDUe4JzyK14l22YZG5/Asm0nmbLyIGsOnsHf14ueTSowqFVl6pVX61ROUmIlIiKSFWIvOB3Sf38PLpx1+lDd9DyUqpFtszwdeZGZaw4z/c/DnDgXQ8USBXnutjr0Ca6ovlMeosRKRETkasTHOTdF/nUsnD8G1W+BW16A8k2ybZb7QiP5YsUBvt4QwqW4BNrWLMXrPevToXYZvL00XJEnKbESERG5Uvt+gaXPwqntENQc7voMqrbNlllZa1l94AxfrNjPjztO4efjxd1Ng7inTVVqlCmcLfOUzFNiJSIiklmn98HS52D3d1CsMvSZCnW6Z8s4VHHxCSzZeoIvVuxnc0gEJQIKMLpjTQa3qkxJ3fw411FiJSIiklEXwmH527B6PPj4Q8dXoOUo8M36O1RcjIvn6/VHGffbPg6fiaZaqQBe71mfu5sG4e/rneXzk6yhxEpERCQ9CfGwfhL8/LrTMb3pELjUqeZkAAAgAElEQVT5eShcJstnFX0pjhmrD/P5iv2cPHeRRkGBPHd7MzrVKYuX+k/lekqsRERE0nJ0Ayx6DI5vgiptoct/4bqGWT6biAuxTF15kIl/HOBsdCytqpXg3d6NubGG7p97LVFiJSIikpIL4fDzGFj7BRQuC70mQr27srwfVcSFWCb8foBJvx/g/MU4br6+DA/dVJ1mlUtk6XwkZyixEhERcWctbJnjdE6PDoOW98NNz4J/1g60eT4mlkl/HOTzFfs5HxNH13rleOSWGhrQ8xqnxEpERCRR2B5Y/DgcWA4VmsHAOVC+cZbOIupiHJNXOglVeHQsneqWZXTHmkqo8gglViIiIvFxsPJDZ5BP34Jw+3vQbBh4Zd3VdzGx8UxZeZDxy/dzJuoSt1xfhtEda9EgSAlVXqLESkRE8rcTW+Hbh5zO6XXvhNveydKr/eLiE5i7PoQPftzDiXMxtKtVmsc71aJxxWJZNg/JPZRYiYhI/hR3ybmv3/J3oGAx6D0F6vXIsuqttSzbfpK3vt/JvtAomlQqxof9GtOyWsksm4fkPkqsREQk/zm2Eb59GE5uhQZ9oOtYCMi6hGf1/tO8+f1ONhwOp3rpAMYNakaXemU1bEI+oMRKRETyj/hY+O1NWPGec7qv/yyofWuWVb/3VCRvLNnBTztPUbaoH2PvakCvZkH4eHtl2Twkd1NiJSIi+UPYHph3n9Na1XigM9Bnwazp53Q26hIf/rSH6X8ewt/Xm6e71mb4DVUpWEC3nslvlFiJiEjeZi2smwBLn3fu6ddnGtTtniVVx8YnMG3VIT78aQ/nY2Lp16ISj3eqRSndHDnfUmIlIiJ5V+Qppy/VnqVQ/Wa481Moet1VV2ut5acdp/jvkh3sD4uibc1SPHd7Ha4vVzQLgpZrmRIrERHJm3YugQWPwKVIuPUtaH4feF19X6e9pyJ5ZeE2VuwJo1rpACYOC+am2mXUMV2ADCZWxpiuwIeAN/CFtXZssumBwHSgkqvOd6y1k7I4VhERkfTFxsCy55x7/JVrAHd9DmXqXHW1URfj+OjnPUz8/QD+vt68eEddBreujK86poubdBMrY4w38H9AJyAEWGuMWWCt3e5W7CFgu7W2mzGmNLDLGPOltfZStkQtIiKSktP7YM4wOLEZWj8Mt7wIPlfX38lay6LNx3l98Q5OnIuhV7Mgnul6PaWLqB+V/FNGWqxaAHuttfsBjDGzgDsB98TKAkWM0w5aGDgDxGVxrCIiIqnb+jUs+JdzG5osGkZhz8nzvLRgGyv3naZe+aL838AmNKtcIguClbwqI4lVBeCI2/MQoGWyMp8AC4BjQBGgr7U2IUsiFBERSUtsDCx91rnyL6g59JoExSpeVZUXLsXz4U97+GLFfgoV8Oa1O+sxoGVlvL3Uj0rSlpHEKqWtyCZ73gXYBNwMVAd+MMassNaeS1KRMSOBkQCVKlXKfLQiIiLuTu+DOUPhxBa44VHn1J+371VV+euuUzw/fyshZy/Qu1kQ/771ekpq+ATJoIwkViGAe+ofhNMy5W44MNZaa4G9xpgDwPXAGvdC1trPgM8AgoODkydnIiIiGbd9Acx/ELx9YMBXUKvLVVV36nwMry3awcK/jlGtdACzRraile7rJ5mUkcRqLVDTGFMVOAr0AwYkK3MYuAVYYYwpC9QG9mdloCIiIgAkxMPPY5wbKFdoBn2mQmDQlVeXYJm19ghjv9tBTGwCozvW5IEO1fHz0ajpknnpJlbW2jhjzMPAUpzhFiZaa7cZY0a5po8DXgMmG2O24Jw6fMZaG5aNcYuISH4UfQa+vhf2/QRNh8Jtb1/VVX97T0Xyn3mbWXvwLK2qleD1ng2oXrpwFgYs+U2GxrGy1i4BliR7bZzb42NA56wNTURExM2JrTB7IEQchTs+gODhV1xVXHwCn684wPs/7qagrzdv9WpI72ZBGuRTrppGXhcRkdxvy1xnFHX/QBi+BCq2uOKqdp44x9NzN7M5JIKu9crxao96lCnin4XBSn6mxEpERHKvhHj44UVY9QlUag29p0CRsldUVWx8Ap/+so9PftlDUX9f/m9AU25rUE6tVJKllFiJiEjuFHMO5o6AvT9A83uhyxvgU+CKqtp6NIKn5m5mx/FzdG9Unpe61dUQCpItlFiJiEjuc/YgzOgHYbvhjvcheMQVVZPYSvXxz3soHlCAzwY3o3O9clkbq4gbJVYiIpK7HF4NswZAQiwM+hqq33RF1ewLjeTx2Zv4KySCOxuX55Xu9ShW6MpavEQySomViIjkHpu/gm8fcsalGvAVlKqZ6SoSEixTVx3kje92UrCAN58MaMIdDctnfawiKVBiJSIinpeQAL+8DivegSptnUE/C2X+ZsfHwi/w1Ny/+GPvaW6qXZo3725ImaK64k9yjhIrERHxrNgY+OZ+2D4fmgyG29+7ok7q8zce5YVvtxKfYHnjrgb0a15RV/xJjlNiJSIinhN9xulPdXgVdHoNbngEMpkMnYuJ5cX5W5m/6RjBlYvzbp9GVC4ZkE0Bi6RNiZWIiHhG+GGY3gvOHoBeE6H+3ZmuYv2hs4yevZFj4TE83qkWD91UA28vtVKJ5yixEhGRnHdii5NUxV6AQfOgattMvT0+wfLpL3v54Kc9XBfoz1f3t6ZZ5eLZFKxIximxEhGRnLXvF5g9GPyLwojvoWzdTL39aPgFHpu1iTUHz3Bn4/K81qM+Rf19sylYkcxRYiUiIjnnr9nw7YNQqjYMnAOBFTL19u+2HOeZrzeTYOH9vo3o2SQomwIVuTJKrEREJGf88aFz378qbaHfl84NlTMoJjae/y7ZwdRVh2hUsRgf9WusDuqSKymxEhGR7GUt/PiSk1jVuwt6jgOfjN+n72BYFA/N2MC2Y+e4r21VnupyPQV8vLIxYJErp8RKRESyT0I8LHoMNkyB4HvgtnfAK+NJ0aLNx/j311vw9jJ8MSSYjnXLZmOwIldPiZWIiGSPuEsw7z5n4M+2T8LNz2d4jKqY2HjGLN7O9D8P06RSMT4Z0JQKxQpmc8AiV0+JlYiIZL1LUc6Vf/t+gs5jnIE/M+hgWBQPfrmB7cfPcX+7ajzZpTa+3jr1J9cGJVYiIpK1LpyFGX0hZC10/wSaDs7wW3/cfpLHvtqEt5dh4rBgbr5ep/7k2qLESkREsk7kKZjWE8J2Q+/JUPfODL0tPsHywY+7+fjnvTSoEMj/BjUlqHih7I1VJBsosRIRkaxx7hhM6eb8HzAbqt+cobedjbrEo7M2smJPGH2Cg3j1zvr4+3pnc7Ai2UOJlYiIXL3wI05SFRUGg7+BSq0y9LYtIRGMmr6e0PMXeeOuBvRvUSmbAxXJXkqsRETk6pw5AFO6Q0wEDJkPQcEZettXa4/w/LdbKRVQgDmjWtOoYrFsDlQk+ymxEhGRK3d6n9NSFRsNQ7+F8k3SfcvFuHheXrCdmWsOc2ONknzUrwklC2d8wFCR3EyJlYiIXJnQXU5LVUIsDF0I5Rqk+5aj4Rd4cPp6/gqJ4IEO1XmiUy18NJSC5CFKrEREJPNOboep3QEDwxZDmTrpvmXl3jAenrmRS3EJjBvUjK71y2V/nCI5TImViIhkzoktTkuVj5/TUlWqZprFrbVMWXmQ1xbvoGqpAMYPbkb10oVzKFiRnKXESkREMu7EViep8i0IwxZBiWppFr8Ul8BLC7Yyc80ROtYpw/t9G1PE3zeHghXJeUqsREQkYxJP//n4ZyipOh15kQe+3MCaA2d4oEN1nuxcG2+vjN0rUORapcRKRETSd2qnk1R5+WYoqdpx/Bz3TV3HqfMX+aBvY3o0qZBDgYp4lhIrERFJW+huZ0gF4+UkVSWrp1l86bYTPDZ7E4X9fJhzv8ankvxFiZWIiKQubK+TVGFh6KI0O6pba/nk5728+8NuGgUF8tmQYMoW9c+5WEVyASVWIiKSstP7YModkBDntFSVrp1q0QuX4nlq7l8s2nycOxuX5827G+p+f5IvKbESEZF/OnvQaamKu+gkVWmMU3U84gIjp65n67EInul6PaPaV8MYdVKX/EmJlYiIJHXuGEy9Ey5FOeNUla2XatEtIRGMmLKW6ItxfD44mI51y+ZgoCK5jxIrERH5W1SYk1RFhcGQBXBdw1SLLtt2gn/N2kSJgAJMf/BGapcrkoOBiuROSqxERMRxIRym9YDwwzDoawhqlmIxay0T/zjImMXbaVghkM+HBlOmiDqpi4ASKxERAbgYCV/2csar6j8LqrRJsVhcfAKvLtrO1FWH6FKvLB/0bULBAuqkLpJIiZWISH4XewFm9oOjG6DPFKjZMcVikRfjeGTGBn7ZFcrIdtX4d9fr8dJI6iJJKLESEcnP4i7BV0Ph4O/QczzU6ZZiseMRFxgxeR27T55nTI/6DGpVOYcDFbk2KLESEcmvEuJh3n2wZync8T406ptisa1HI7hnylqiLsYzcVhz2tcqncOBilw7lFiJiORH1sLiJ2D7fOj0GgSPSLHYTztO8sjMjRQr6MvcB1pzfbmiORyoyLVFiZWISH70y+uwfhK0eQxufDTFIlNWHuSVhduoVz6QCUODKaPb04ikS4mViEh+8+c4WP42NBkMt7z0j8kJCZax3+/ks+X76VS3LB/2a0yhAvq6EMkI7SkiIvnJ5jnw/TNw/R1wxweQ7NYzF+PieWrOZhb8dYwhrSvzUrd6eOvKP5EMU2IlIpJf7PkB5o+CKm3h7gngnfQr4FxMLKOmrWflvtO655/IFVJiJSKSHxxZA7MHQ5m60G8G+CbtL3XyXAxDJ65h76lI3uvTiLuaBnkoUJFrmxIrEZG87tQO+LI3FL3OuVWNf9Ir+/aeOs/QiWsJj77EpOHNaVtTwymIXCklViIieVn4YZh2F/j4w+BvoHCZJJPXHjzDvVPWUcDHi9n3t6Z+hUAPBSqSNyixEhHJq6LCYFpPiI2C4d9B8SpJJn+/9TiPztpEULGCTBnRgoolCnkmTpE8RImViEhedPG8c1PliBAY8i2UrZdk8pSVB3l54TaaVCzGF0ObUyKggIcCFclblFiJiOQ1cRdh1kA4vhn6z4RKrS5Pstby1tJd/O/XfXSsU5aP+zehYAFvDwYrkrcosRIRyUsSEmD+A3DgN+gxDmp1uTwpLj6Bf8/bwtz1IQxoWYlXu9fDx9vLg8GK5D1KrERE8pIfX4KtX0PHl6Fx/8svx8TG89CXG/hp5ylGd6zJv26pqTGqRLJBhn6qGGO6GmN2GWP2GmP+nUqZDsaYTcaYbcaY37I2TBERSdfq8bDyI2h+L9w4+vLLERdiGTxhNT/vOsVrPeozumMtJVUi2STdFitjjDfwf0AnIARYa4xZYK3d7lamGPAp0NVae9gYUybl2kREJFvsWAjfPQO1b4db37p8q5pT52IYMnEN+0Ij+bh/E+5oWN7DgYrkbRk5FdgC2Gut3Q9gjJkF3AlsdyszAJhnrT0MYK09ldWBiohIKg6vhq/vhaBguPsL8HI6ox8Mi2LwxNWcjrzExGEa+FMkJ2TkVGAF4Ijb8xDXa+5qAcWNMb8aY9YbY4akVJExZqQxZp0xZl1oaOiVRSwiIn8L2wMz+0LRCtB/NhRwxqLadiyCXuNWERkTx8z7WimpEskhGUmsUjoRb5M99wGaAbcDXYAXjDG1/vEmaz+z1gZba4NLl9ZOLiJyVc6fhOl3gZePc6uagJIA/Ln/NP3G/0kBb8OcUTfQqGIxDwcqkn9k5FRgCFDR7XkQcCyFMmHW2iggyhizHGgE7M6SKEVEJKmLkTCjtzO6+rBFUKIqAEu3neCRmRupVKIQU0e0oHyxgh4OVCR/yUiL1VqgpjGmqjGmANAPWJCszLdAW2OMjzGmENAS2JG1oYqICADxsTBnKJzYCr2nQIVmAHy19ggPTF9P3euKMuf+1kqqRDwg3RYra22cMeZhYCngDUy01m4zxoxyTR9nrd1hjPke2AwkAF9Ya7dmZ+AiIvmStbBoNOz9Ebp9BLU6AzD+t3288d1O2tYsxbhBzQjw0zCFIp6QoT3PWrsEWJLstXHJnr8NvJ11oYmIyD/8OhY2Tof2z0CzoVhreXvpLj79dR93NLyO9/o0poCPRlMX8RT9pBERuVZsnA6/jYXGA6HDf0hIsLy8cBtTVx2if4uKjOnRAG8vDfwp4klKrERErgX7f4OF/4JqN0G3D4lLsDw9dzPzNh5lZLtq/OfW6zWaukguoMRKRCS3C90FswdDyZrQZwoxCV48OnMDy7af5MnOtXjophpKqkRyCSVWIiK5WWQofNkbfArAgNlEmQBGTlnLH3tP83K3ugy7saqnIxQRN0qsRERyq9gYmDUAIk/CsCVE+JVn+ITVbDoSzru9G3F3syBPRygiySixEhHJjRISYP4DELIG+kwlNLA+gz9bxf7QKD4d2Iyu9ct5OkIRSYESKxGR3OiXMbBtHnR8haPlOzNo/CpORMQwYViw7vsnkospsRIRyW02TocV70LToeyvdQ+D/reS8xfjmHZPC4KrlPB0dCKSBiVWIiK5iduwCtuavsjQz/7EWpg1shX1ygd6OjoRSYeG5xURyS1Cd8NXg6FkDTa2/oB+X6zH19uLr0a1VlIlco1QYiUikhtEhcGXvcC7AGtvGM+AqTspGVCAOaNaU710YU9HJyIZpFOBIiKeFhsDM/tD5ElWtZ3K0LnHqVY6gKn3tKBMEX9PRycimaDESkTEk9yGVfiz2XsMWhpPw6BAJg9rQWAhX09HJyKZpMRKRMSTfnkdts1jfa3R9PujHDdUL8HnQ4IJ8NPhWeRapD1XRMRTNn4JK95hW7ke3L25OZ3qluXj/k3w9/X2dGQicoWUWImIeMKB5diF/+JQYAvuPHgXPRpX4O3ejfD11jVFItcyJVYiIjktdDd29iBCCwTR/eR99G5Rjdd71MfLy3g6MhG5SkqsRERyUlQYdkZvImO9uCt6NH3a1Oe52+tgjJIqkbxAiZWISE6JjSFhZn/iwo8xOOZ57r75BkZ3rKmkSiQPUWIlIpITEhKI/+YBvEPW8K9L/+K2W+9gZLvqno5KRLKYEisRkRxw6ccxFNg+j7Gx/bih+z0MblXZ0yGJSDZQYiUiks2i10yj0Mp3mR3fgZo9n+fu4IqeDklEsokSKxGRbHRuxy8UXDKaPxLqE9jrY7o2UlIlkpdpwBQRkWwSdnArzB7EIVuWhN5T6NqokqdDEpFsphYrEZFscOzYERKm3I2xXpzrOYO2DWp4OiQRyQFqsRIRyWIHjocR+nkvStnThHWbTNPGjT0dkojkECVWIiJZaOfxCHZ9NpRGdiehHT+idvAtng5JRHKQEisRkSzy15Fwfhv/GF3t74S1+g8V2wzwdEgiksOUWImIZIHV+08z+4u3uJ+viaw7gFJdnvF0SCLiAeq8LiJylX7bHcqEaVOY4D2eixXbUvjuj0C3qRHJl5RYiYhchaXbTvDejMXMLfA+png1/AZMB29fT4clIh6ixEpE5ArN33iUMXNWsLDg2wT4+eM1eA4ULObpsETEg5RYiYhcgRmrD/PK/PV8W+RDyiWcxfRfBMWreDosEfEwJVYiIpn0xYr9vL54GzNLTOb66O3QezJUbO7psEQkF1BiJSKSQdZaPv55L+/9sJtPy39PqzO/QseXoV5PD0cmIrmFEisRkQyw1jL2u52MX76fsdU2c9ux6dB0CNw42tOhiUguosRKRCQdCQmWFxdsZfqfh3mxXhh9D7wD1W6C29/TsAoikoQSKxGRNMTFJ/D015uZt+Eoz7YwDN/1PKZkTegzRcMqiMg/KLESEUnFpbgE/jVrI99tPcHz7Utyz677MD7+MPAr8A/0dHgikgspsRIRSUFMbDyjpq/n112hvNS1KsP3PAKRoTB8MRSr5OnwRCSXUmIlIpLM+ZhY7pu6jtUHzvBGz3r0P/gCHF0PfadDhWaeDk9EcjElViIibs5EXWLYpDVsO3aOD/o25s7Q8bBjAXT5L9S5w9PhiUgup8RKRMTlREQMgyas5siZaD4b3IxbopbAHx9C83uh1YOeDk9ErgFKrEREgINhUQyasJrw6FimjGhBq4SNsPgJqNEJur6pYRVEJEOUWIlIvrfj+DkGT1hDfEICM+5rSUPfozBhGJSpC70ngbcOlSKSMV6eDkBExJPWHzpL3/Gr8PEyfHV/axoWvQBf9gG/wjBgNvgV8XSIInIN0c8wEcm3VuwJZeTU9ZQp6sf0e1pSsbCFSX3hwlkY8R0EVvB0iCJyjVFiJSL50vdbj/PozE1UKx3A1HtaUCbAF2YPghNboN9MuK6Rp0MUkWuQEisRyXfmrDvCM19vpnHFYkwa1oLAgj7w3dOwawnc+jbU7urpEEXkGqXESkTylQm/H+C1RdtpU6MU4wc3I8DPB/74CNZ8Bq0fhpYjPR2iiFzDlFiJSL5greX9H/fw0U976FqvHB/2b4yfjzdsmQs/vAD1ekKn1zwdpohc45RYiUiel5BgeXXRdiavPEivZkGMvasBPt5ecGAFzH8AKt8IPcaBly6UFpGro8RKRPK0uPgEnv56M/M2HGXEjVV5/vY6eHkZOLUDZg2E4lWh35fg6+/pUEUkD8jQzzNjTFdjzC5jzF5jzL/TKNfcGBNvjOmVdSGKiFyZmNh4Rk3fwLwNR3msYy1euMOVVJ07BtN7gW9BGDQXChb3dKgikkek22JljPEG/g/oBIQAa40xC6y121Mo9yawNDsCFRHJjIjoWO6dupZ1h87ySvd6DL2hijMh5hx82RtiwmH4d1CskkfjFJG8JSMtVi2Avdba/dbaS8As4M4Uyj0CfA2cysL4REQy7eS5GPqMX8WmI+F81K/J30lV3CX4ajCE7oQ+U+G6hh6NU0TynowkVhWAI27PQ1yvXWaMqQD0BMalVZExZqQxZp0xZl1oaGhmYxURSdf+0Eju+nQlIWejmTSsBd0alXcmWAsLH4X9v0L3j6HGLR6NU0TypowkVind0t0me/4B8Iy1Nj6tiqy1n1lrg621waVLl85ojCIiGbI5JJxe41YRExvPzJGtaFOz1N8Tf3kd/poJNz0PjQd4LkgRydMyclVgCFDR7XkQcCxZmWBgljEGoBRwmzEmzlo7P0uiFBFJx+97wrh/2jqKFSrAtHtaUK104b8nrpsEy9+GpkOh3ZOeC1JE8ryMJFZrgZrGmKrAUaAfkOTnnrW2auJjY8xkYJGSKhHJKYs2H+Ox2ZuoXrowU0a0oGxRt6ETdn0Pix+Hmp3h9vfApNQILyKSNdJNrKy1ccaYh3Gu9vMGJlprtxljRrmmp9mvSkQkO01ddZCXFmwjuHJxvhjanMCCvn9PPLQK5gx1bqjcaxJ4a+g+EcleGTrKWGuXAEuSvZZiQmWtHXb1YYmIpM39FjUd65TlkwFN8Pf1/rvAyW0wsy8EVoSBc8GvcOqViYhkEf18E5FrTlx8Ai98u42Zaw7TJziI//Z03aIm0dlDMO0u8A2AwfMgoFTqlYmIZCElViJyTYm+FMcjMzby085TPNihOk91qY1x7zcVGQrTekJcDIz4XgOAikiOUmIlIteMsMiL3DN5LVuORvBaj/oMblU5aYGL5+HLXs4ta4Z8C2XqeCZQEcm3lFiJyDXhQFgUQyeu4dT5GMYPDqZT3bJJC8RddG6qfGIL9J8FlVp6JlARydeUWIlIrrfh8FnumbwWYwwz72tFk0rJbpqcEA/z7oMDv0HP8VCrs2cCFZF8T4mViORqy7ad4JGZGykX6M+U4S2oUiogaQFrYcmTsP1b6PJfaNTPM4GKiKDESkRysWmrnDGqGgQVY8LQYEoV9vtnoV/HwrqJ0OYxaP1QjscoIuJOiZWI5DoJCZa3lu5i3G/76FinDB/1b0KhAikcrtZ8Dr+NhSaD4JaXcj5QEZFklFiJSK5yMS6ep+du5ttNxxjYshKvdK+XdIyqRFu/hiVPQe3b4I4PdasaEckVlFiJSK4RER3LqOnrWbX/NE93rc0D7asnHaMq0e6lMG8kVGoNvSbqVjUikmvoaCQiucLBsChGTFnLkTPRvN+3ET2bBKVc8MBymD0YyjWAAbPBt2DOBioikgYlViLicWsOnOH+aeuwwPR7WtKyWsmUCx5ZCzP6QYlqMGge+BfN0ThFRNKjxEpEPOqbjSE8M3cLQcULMnFY838Op5DoxBb48m4oUhaGzIdCJXI2UBGRDFBiJSIeYa3l/R9289HPe2ldrST/G9SUYoUKpFw4bI9z/78ChZ1b1RQpl7PBiohkkBIrEclxMbHOlX8L/jpGn+AgxvRoQAGfFK78Azh7CKbe6TweskA3VRaRXE2JlYjkqLDIi4ycuo4Nh8N5puv1jGpfLeUr/wDOn3CSqkuRMGwJlKqRs8GKiGSSEisRyTF7Tp5n+OS1hEVe5H8Dm3Jrg+tSLxx5CqZ0g6hQ5/Rfufo5F6iIyBVSYiUiOWLFnlAenL4B/wLezB7ZmkYVi6VeOCoMpnSHiBAY9DUEBedcoCIiV0GJlYhkK2stU1Ye5LXFO6hZpjAThjWnQrE0xp6KPuOc/jt7EAZ+BZVvyLFYRUSulhIrEck2F+PieXH+NmavO0LHOmX5oF9jCvulcdi5cNZJqk7vhf6zoGq7nAtWRCQLKLESkWxx6nwMo6atZ8PhcB69uQajO9bCyyuN+/ldCHeGVAjdyf+3d9/xWZX3/8dfVwIZJCGDDDIhYe89RIsCDhQV66ggrg791mqtrR1f7bDt79uHdtcObXEVK0oRaUErUkerCJUlIBtCgCRk7z3v6/fHSTBFkEDukdx5Px+PPO517nM+PC6Jb865zudi8UswZI73ihURcRMFKxFxu49yK7j7+e1U1jfzxJLJXPVpk9QBGqrghRugYA8sWg5DL/VOoSIibqZgJSJu9bcdufzvK7uJDQ/mlXtmMTrpLMvONFbD8hshfyd87nkYfoV3ChUR8QAFK9+bOAQAABvhSURBVBFxi1aX5advHGDpe1nMSI/hiSWTGRAe/OlfaqiEF26EE9vhpudg5ALvFCsi4iEKViLSZZV1zdy/YgfvHirm9gsG8f2rR9M38Ayd1NvVV8AL10P+LrjpzzD6Wq/UKiLiSQpWItIlmUXV3PX8dnLL63j0+nEsnt6JJWfqyuAv10HhPvjcX2DkVZ4vVETECxSsROS8vbmvkK//dSchfQN48a6ZTBscc/Yv1ZbA89dBySFY9CIMv9zzhYqIeImClYics1aX5VdvHuQP/zrCuORI/njblE9v+tmuphievxbKspyWCkPneb5YEREvUrASkXNSVtvE/S/t4P3MEhZPT+WRa8YQ0jfw7F+sLmhbpiYHblkJGRd7vlgRES9TsBKRTtuZU8FXXthOSW0TP71hHDdP68R8KnDW/Ht+IVTlw5JVMPhCzxYqIuIjClYiclbWWl7aksMP1+4lLiKYV748i3EpkZ37ckmmM1G9oRJu+xukzfBssSIiPqRgJSKfqr6plR+s2cPL23OZPTyOx2+eSHRYUOe+XLDbWabGWrjzNUic4NliRUR8TMFKRM4os6iGe5d/yMHCar7att5f4Ket99dR9mZ48SYIioDb/w6xwzxbrIhIN6BgJSKntWbnCR5avZuQvoEs+8J0Lh4e1/kvH3kHViyBiES4fQ1EpXquUBGRbkTBSkT+S0NzKz96dS8vbclh2uBofrt4EomRnWil0G7fWnjlixA7Am5bDeHxnitWRKSbUbASkZOyimv4yvIPOVBQzT2XDOHBy4bT52xL03S0YzmsvQ+Sp8KSlRAa7bliRUS6IQUrEQGcS38Pr95NUJ8Anvv8NOaMOIczTdbC+7+Gt38EGXNg0XIICvNcsSIi3ZSClUgvV9vYwiNr97Jqey5TBkXzu8WTSOpMF/V2rlZY9x3Y+hSMuwkWPgF9OnnXoIiIn1GwEunFdudWcv+KHRwrreWrc4dy/7xh9D2XS3/NDbD6S7D/VZh1P1z6Iwg4h++LiPgZBSuRXsjlsjz9fhY/X3+Q2PBgXrprJjMzBpzbTurL4aVbIPs/MP8xmHmPZ4oVEelBFKxEepmi6gYeXLmLDYdLuGJMAj+9YTxR/c7x0l1FDiy/0VlM+cZnYez1nilWRKSHUbAS6UXeOVDIt17+iNqmFn7y2bHcMj0NYzrZ8LNd4V544UZoqoVbV0P6ZzxTrIhID6RgJdIL1Da28JPX9/Pi5mxGDoxgxeKZDEuIOPcdHX4LXr4TgsPhC+sgYYzbaxUR6ckUrET83Pbj5Xxj5U6yy+r4n9kZfOPy4QT3CTz3HW1eCm98xwlTi/8KkcnuL1ZEpIdTsBLxU00tLn779mGe+HcmiZGhrLhrJjPOdYI6QGsLrH8ItiyFEVfB9U85Z6xEROQTFKxE/NDhwmq+vnIne05UceOUFB65ZjQRIX3PfUcNVbDqC5D5JlxwH1z2Ywg4j7NdIiK9hIKViB9pdVme23iUn60/SHhwH/546xTmjx14fjuryIYXb4aSQ3DN4zDlTrfWKiLijxSsRPzEkeIavvXyLj7MruDSUfE8ev144iKCz29nOVthxWJoaYJbX4GMS9xZqoiI31KwEunhWl2Wpzdk8cs3DxHaN5Df3DyRhROTzr2NQrvty+AfD0JkCty5EuKGu7dgERE/pmAl0oMdLqzmm6s+YldOBZePTuD/PjuW+IiQ89tZS5Nz19+2Z2HIXLjhGegX496CRUT8nIKVSA/U0upi6YYsfvPmYcKCA/nt4klcMz7x/M9SVRfAytshZzNc+ADM+4EmqYuInAcFK5EeZldOBQ+t3s2+/CquHDuQHy8ce/5zqcCZT7XyNmiohBuf0/I0IiJdoGAl0kPUNLbwy38eZNmmY8SGB/PkkslcOS6xazv98HlnPlVEInzxTRg41j3Fioj0Up0KVsaY+cDjQCDwtLX2sVM+XwJ8p+1lDXCPtXaXOwsV6c3e2lfID9bsIb+qgSUz0vj2/JH0P5++VO2a6+H1b8GOv0DGHGchZc2nEhHpsrMGK2NMIPAH4DIgF9hqjFlrrd3XYbOjwMXW2nJjzJXAUmCGJwoW6U2Kqhr44at7eX13AcMTwll1ywVMGdTFAFSSCS/fAYV74DMPwiUPQ6BOXouIuENnfptOBzKttVkAxpgVwELgZLCy1m7qsP0HQIo7ixTpbVpaXbzwwXF++eYhGltcfPPy4dw9ewhBfQK6tuPdq+DVr0FgECx5BYZd6p6CRUQE6FywSgZyOrzO5dPPRn0RWNeVokR6s+3Hy/je3/eyP7+Ki4bG8v+uG0t6bFjXdtrcAOsfhm3PQOoM59JfpP79IyLibp0JVqe7f9uedkNj5uAEq4vO8PndwN0AaWlpnSxRpHcoqWnksXUHWLU9l8TIEJ5YMpkrxw48/xYK7cqy4OU7IX8XzLrfaaUQ2IX5WSIickadCVa5QGqH1ylA3qkbGWPGA08DV1prS0+3I2vtUpz5V0ydOvW04Uykt2l1WZZvPs4v1h+kvrmVey4Zwn1zhhIW7IZ5Tx+9DP/4BpgAWLwCRlzZ9X2KiMgZdeY391ZgmDEmHTgBLAJu6biBMSYNWA3cZq095PYqRfzU5qxSfvzaPvbmOZf9fnjtGIbGh3d9xw2VThuF3S87l/6ufwqiB3V9vyIi8qnOGqystS3GmPuA9TjtFp611u41xny57fM/Aj8ABgBPtF22aLHWTvVc2SI9W3ZpHY+u28+6PQUkRYbw+1smsWBcFzqnd3R8E6z+H6g6AXO+Cxd9Q3f9iYh4ibHWN1fkpk6dardt2+aTY4v4SnVDM0/8+wjPbDhKYIDhK5cM4a7ZGYT0dcPyMa3N8O9H4f1fQ9QguOFpSNG/b0RE3MEYs70zJ430z1gRL2h1WVZtz+Hn6w9RUtPI9ZOT+fYVIxkYeZ4LJp+q9Ai88iXI+xAm3QrzH4PgCPfsW0REOk3BSsSDrLW8d7iEx9YdYH9+FVMGRfPMHVOZkBrlngO4XLBlKbz1Q+gTDJ97HkYvdM++RUTknClYiXjIR7kVPLbuAJuOlJIaE8pvF0/imvFumkcFzlmqNfdB9iYYdjlc8zj0T3LPvkVE5LwoWIm42fHSWn6+/iCvfZRPTFgQj1wzmiUzBnW9a3o7Vyts/hO8/WOng/p1T8KExeCuwCYiIudNwUrETUpqGvnd24dZvjmbvoEBfHXuUO6enUFEVxZL/sRBMmHNvZDzAQyfD1f/WmepRES6EQUrkS4qr23iT+9lsWzTMZpaXdw8LZUH5g0jvr+bJqaDc8ffB0/Cv34CfULgs0th/Od0lkpEpJtRsBI5T5X1zTyzIYtnNx6jtqmFayck8bV5w8iIc0ODz45ytsJrD0DhHhixAK7+FUQMdO8xRETELRSsRM5RTWMLz71/lKc2ZFHV0MJV4wbywKXDGZ7g5vYG9eXw1o9g+5+dy303L4eRC3SWSkSkG1OwEumkqoZm/vKf4zy9IYvyumYuHZXA1y8bxpikSPceyFpnKZr1D0NdGVxwL1zyEAS7+UyYiIi4nYKVyFmU1Tbx3Maj/HnTMaobWrhkRBxfv3S4+3pRdVRy2Fnj7+i7kDwFbl0NiePdfxwREfEIBSuRMyiqauCpDVks35xNXVMr88cM5L65Qxmb7OYzVAD1FfDuz2DLn6BvGCz4JUz5PAS4YakbERHxGgUrkVNkl9bx1IYs/roth5ZWF9dOSOIrc4a6fw4VOD2pPlwG7/yfc9lv8u0w9/sQHuf+Y4mIiMcpWIm02ZVTwdL3sli3J5/AAMMNk1P48sVDGBwb5pkDHn0P3njIudtv0IUw/1FInOCZY4mIiFcoWEmv5nJZ/nWwiD+9l8WWo2VEhPTh7tlDuHPWYPctkHyq0iPw1iOw/1WITIObljnr++luPxGRHk/BSnqlhuZW1uw8wVMbjpJZVENSZAjfWzCKRdPTCA/20F+L6kJ496fOpb/AIJjzPZh1H/QN9czxRETE6xSspFc5UVHPCx8cZ8WWbMrrmhmd2J/HF03kqnGJ9A1001p+p2qohI2/hQ+egNYmmHwHXPxtNfkUEfFDClbi96y1bD5axrJNx1i/twCAy0cP5I5Zg5mZEYPx1CW45gbY+jRs+IXT7HPsDTDnuzBgiGeOJyIiPqdgJX6rtrGFtbvyWLbpGAcKqonq15e7Zw/h1plppET389yBWxphxwuw4VdQlQtD5sK8RyBpoueOKSIi3YKClfidvXmVvLg5mzU786hpbGHkwAgeu34cCycmExrkwb5QzQ2w4y/w/q+h6gSkTIPrnoCMiz13TBER6VYUrMQv1DW18NqufJZvyWZXTgXBfQJYMD6RJTPSmJwW7bnLfeAEqg+Xwfu/geo8SJ0JC38PGXN0p5+ISC+jYCU9lrWWXbmVvLwth7U786hubGFofDg/uHo0109OJqpfkGcLaKyBD5+HjY9DTQGkzYLPPgnpFytQiYj0UgpW0uMUVTXwtx0nWLU9l8NFNQT3CeCqcYksnp7GtMEePjsFUFMEm//kTExvqIBBF8ENT8HgzyhQiYj0cgpW0iM0trTy9v4iVm3P5d1DxbS6LFMGRfPo9eNYMD6R/iF9PV9EyWHY9DvYtcJpmzByAVz4NUid7vlji4hIj6BgJd1Wq8uyOauUNTvzeH1PPtUNLST0D+bu2RncOCWFIXHhni/CWsj+D2z6PRx83WnsOfEWuOA+iB3q+eOLiEiPomAl3Yq1lr15Vfx9xwle/SiPwqpGwoICuWLMQBZOSuaiobEEBnjhcltTHex+GbY8BYW7ITQaZn8Lpt+tBZJFROSMFKzE56y17Muv4o09Bfxjdz5ZxbX0DTRcPDye7y1I4tJRCZ5tk9BRWRZsfcZpm9BQCfFj4JrHYdxNEOShxZhFRMRvKFiJT1hr2ZlTwRt7Cli3p4DssjoCDExPj+FLF2Vw1biBnr+rr11rCxx52wlUh/8JJgBGX+ucnUq7QBPSRUSk0xSsxGtaXZbtx8tZtyef9XsKyKtsoE+AYdbQWO65ZAiXjU4gNjzYewWVHXU6pO9cDtX5EBbvrOE35U7on+S9OkRExG8oWIlHNTS38kFWKW/tL+SNPYWU1DQS1CeA2cPiePDyEVw6KoHIfl64o69dcz3sf81p6Hlsg3N2auhlcNXPYfh8CPRiLSIi4ncUrMTt8ivreedAEf86UMTGzFLqm1sJ7RvInJFxzB+byNyR8YQHe/E/PZcLcj6Aj1bC3tXO3KmoQTD3ezDhFohM9l4tIiLi1xSspMtaXZYd2eW8c6CIdw4UcaCgGoCU6FBumprCnJHxXJAxgJC+XpqA3q5ovxOmdq+Cymzo28/pPTXpNqeZZ0CAd+sRERG/p2Al56WoqoGNR0r498Fi3j1UTEVdM4EBhqmDonnoypHMHRnP0Phwz3dBP1VFNuz9O+xeCQW7wQTCkDkw7/sw4ioI9kLvKxER6bUUrKRTqhqa2ZxVxsbMEjZmlnC4qAaAmLAg5o6IZ+6oeD4zLI7IUB/MUSo9AvvWwP61kLfDeS95Klz5MxjzWQiP935NIiLSKylYyWk1trSy/Xg5mzJLeT+zhI9yK3BZCOkbwPT0Adw4JYULh8YyOrE/Ad5o2NmRtVB80AlS+9Y6DTwBkqfAZT+GUddCTLp3axIREUHBStrUN7WyI6ecrUfL2XqsjG3Hy2hodhEYYBifEsm9c4Yya0gskwdFEdzHy3OlAFqaIHsTHHwDDq2D8mPO+6kz4YpHYdQ1EJXq/bpEREQ6ULDqpcprm9oCVDlbjpax50QlLS6LMTAiIYJF09K4cGgsMzJivLPA8enUlkLmm3BwHRx5BxqrIDAY0mfDrK86c6bUb0pERLoRBatewFrLsdI6duaUs/VYOVuPlp2cIxUUGMD4lEjump3BtMHRTEmL8W5fqY5amiB3ixOijrwDeTsBC+EJMOY6p89UxiVaWkZERLotBSs/VFLTyK6cCnblVLAzt5JdORVU1jcDEBHch8mDorluUjLTBscwPiXS+20Q2rXPlTr6rhOkjm6A5lrnTr7U6TDnYRg6DxInqTWCiIj0CApWPVxtYwv78qvYlVPBjrYwlVteD0CAgeEJEVw5diATU6OYkBrF8IQIAr092bydywVFe+HYRjj+PhzfBHWlzmfR6TBhkROkBl8EIZG+qVFERKQLFKx6kOLqRvblV7E3r5K9eVXsy6viWGkt1jqfJ0eFMjE1itsvGMSElCjGpUTSL8iHQ9zaAgUfwfGNTpjK3uR0PQeITINhl8OgC2HwhRCT4bs6RURE3ETBqhtqdVmOl9ZyoKCafXkfB6mi6saT26REhzImqT/XTUxmTFJ/xqdGEh8R4sOqgdoSyN0GJ7bDiW2QsxWanC7sxAyB0QudIDVoFkSl+bZWERERD1Cw8iGXy5JbXs+hwmoOFVVzqKCaQ4U1ZBbX0NTiAiAwwDA0LpyLhsYyOqk/o5P6MyYx0ncTzNs110P+LidEtYepiuPOZyYA4sfA+M85Z6PSZkH/RN/WKyIi4gUKVl7Q1OIiu6yOrOIaskpqySyq4VBhNYcLa6hvbj25XVJkCMMSIrhoWCzD4sMZMTCC4QkRvptc3q6xBor2OUvEFOx2QlTRPnC1OJ9HpjrNOad9CVKmQuIE3bknIiK9koKVm1hrKa5pJKu4tu3HCVFHS2rJLquj1WVPbhsXEez0ipqeyvAEJzwNSwj3Xb+odtZCVZ4TngrbQlTBHijLAtrqD4mEpElw4QNOmEqeAhEJPi1bRESku1CwOgeNLa3klteTXVZHblkd2Sd/6skpq6OmseXktsF9AkiPDWNUYgQLxiWSERdGRlw46bFhvllPryNroaYISg467Q6KD0LxASjcA/XlH28XnQ4Dxzp36yWMdZ5HpoK3F1YWERHpIRSsOmh1WYqrGzlRUUdOWX2H4FRHTlkdBVUNJ+/AAyc8pcX0Iy2mHzPSYxg8oB8ZceFkxIWRFBnq/TX0TuVqhYpsKDnkhKeSg1B8yHlsvzsPICgC4kY4a+wNHOf8xI+GkP6+q11ERKQH6jXBylpLaW0T+RUN5FXWk19RT35lA3mVDSefF1Y10NLhkh3AwP4hpMX0Y9aQWNJi+pEaE3oyTMVFBGN8ffampckJT2VZUH7UeSzLgrKjznp6ruaPtw2Lg7iRMPZGJ0jFDnceIxJ1FkpERMQN/DZYZZfW8Zu3D5Ff0UB+pROcGtvutGsXFBhAYlQIiZEhzEiPaXseSnJUKKkx/UiJDvX9xPHWFqjOh8ocqMiBymyozHVCU9lR533b4c8VFA4x6RA/CkYucPpDtYeofjE++2OIiIj0Bn4brFqt5YMjpSRGhTI2OZLLxwwkMTKEpKhQkiJDSYwKYUBYkG/POLlaobYYqgvafvKg8kSHEJXjTCa3rf/9vbA4pw9U6nRn/lN0uhOgYtKdz3T2SURExCf8Nlilx4ax6aF5vjm4ywX1Zc6ZpuqCMz/WFH0yNJlA6J8MUalOM82oVIhMcSaNR6U5z/uG+ubPJSIiIp/Kb4OVW7lczt1ydSVOd/Ha4rbnpR2et/3UlUBd2ScDE0C/Ac58poiBkDDm4+ftj+EDITwBAjUsIiIiPVHv+T+4tdBcB/UV0FDR4bH8NO+d8lhX+t/zmDoKiXQuv/WLhQFDnMtzYXHOT3to6p/oBKY+wd79M4uIiIhX+W+wKj4Ea77y3yGp4x1yn2CckBQaBSFRzmNkivO834C2sBT7yeeBPu5JJSIiIt1Gp4KVMWY+8DgQCDxtrX3slM9N2+dXAXXAndbaD91c67npE+wsq9I/+b/DUkgUhEZ/8r3g/hAQ4NOSRUREpGc7a7AyxgQCfwAuA3KBrcaYtdbafR02uxIY1vYzA3iy7dF3ogfB7Wt8WoKIiIj0Lp05RTMdyLTWZllrm4AVwMJTtlkIPG8dHwBRxphEN9cqIiIi0q11JlglAzkdXue2vXeu24iIiIj4tc4Eq9N1m7TnsQ3GmLuNMduMMduKi4s7U5+IiIhIj9GZYJULpHZ4nQLkncc2WGuXWmunWmunxsXFnWutIiIiIt1aZ4LVVmCYMSbdGBMELALWnrLNWuB245gJVFpr891cq4iIiEi3dta7Aq21LcaY+4D1OO0WnrXW7jXGfLnt8z8Cr+O0WsjEabfwec+VLCIiItI9daqPlbX2dZzw1PG9P3Z4boF73VuaiIiISM+ijpgiIiIibqJgJSIiIuImClYiIiIibqJgJSIiIuImClYiIiIibqJgJSIiIuImxumU4IMDG1MMHPfJwT0vFijxdRHiMRpf/6Wx9V8aW//mjfEdZK0967IxPgtW/swYs81aO9XXdYhnaHz9l8bWf2ls/Vt3Gl9dChQRERFxEwUrERERETdRsPKMpb4uQDxK4+u/NLb+S2Pr37rN+GqOlYiIiIib6IyViIiIiJsoWHWRMSbVGPMvY8x+Y8xeY8zX2t6PMca8aYw53PYY7eta5fwYYwKNMTuMMa+1vdbY+gFjTJQxZpUx5kDb398LNLb+wxjz9bbfyXuMMS8ZY0I0vj2TMeZZY0yRMWZPh/fOOJbGmIeMMZnGmIPGmCu8Xa+CVde1AA9aa0cBM4F7jTGjgf8F3rbWDgPebnstPdPXgP0dXmts/cPjwBvW2pHABJwx1tj6AWNMMnA/MNVaOxYIBBah8e2p/gzMP+W9045l2/9/FwFj2r7zhDEm0HulKlh1mbU231r7YdvzapxfzsnAQmBZ22bLgOt8U6F0hTEmBVgAPN3hbY1tD2eM6Q/MBp4BsNY2WWsr0Nj6kz5AqDGmD9APyEPj2yNZa98Dyk55+0xjuRBYYa1ttNYeBTKB6V4ptI2ClRsZYwYDk4DNQIK1Nh+c8AXE+64y6YLfAN8GXB3e09j2fBlAMfBc22Xep40xYWhs/YK19gTwCyAbyAcqrbX/ROPrT840lslAToftctve8xoFKzcxxoQDrwAPWGurfF2PdJ0x5mqgyFq73de1iNv1ASYDT1prJwG16LKQ32ibb7MQSAeSgDBjzK2+rUq8xJzmPa+2P1CwcgNjTF+cULXcWru67e1CY0xi2+eJQJGv6pPzdiFwrTHmGLACmGuMeQGNrT/IBXKttZvbXq/CCVoaW/9wKXDUWltsrW0GVgOz0Pj6kzONZS6Q2mG7FJzLwF6jYNVFxhiDM09jv7X2Vx0+Wgvc0fb8DmCNt2uTrrHWPmStTbHWDsaZDPmOtfZWNLY9nrW2AMgxxoxoe2sesA+Nrb/IBmYaY/q1/Y6ehzP/VePrP840lmuBRcaYYGNMOjAM2OLNwtQgtIuMMRcBG4DdfDwP52GceVYrgTScv+Q3WWtPnXwnPYQx5hLgm9baq40xA9DY9njGmIk4NyUEAVnA53H+samx9QPGmB8BN+Pcub0D+BIQjsa3xzHGvARcAsQChcAjwN85w1gaY74LfAFn7B+w1q7zar0KViIiIiLuoUuBIiIiIm6iYCUiIiLiJgpWIiIiIm6iYCUiIiLiJgpWIiIiIm6iYCUiIiLiJgpWIiIiIm6iYCUiIiLiJv8foblOLb3QKzYAAAAASUVORK5CYII=\n", + "text/plain": [ + "
" + ] + }, + "metadata": { + "needs_background": "light" + }, + "output_type": "display_data" + } + ], + "source": [ + "plt.figure(figsize=(10,6))\n", + "data_pred = np.linspace(start=15, stop=100, num=150)\n", + "smoker_proba = log_smokers.predict_proba(data_pred.reshape([-1, 1]))\n", + "plt.plot(data_pred, smoker_proba[:, 1], label=\"Smokers\")\n", + "non_smoker_proba = log_non_smokers.predict_proba(data_pred.reshape([-1, 1]))\n", + "plt.plot(data_pred, non_smoker_proba[:, 1], label=\"Non-smokers\")\n", + "plt.title(\"Estimated probability of death per age, according to smoking habits\")\n", + "plt.legend()\n", + "plt.show()" ] } ],