Skip to content

Instantly share code, notes, and snippets.

@akelleh
Created July 21, 2018 19:26
Show Gist options
  • Save akelleh/7d559d948ba80a1c5f9469be9d9ae117 to your computer and use it in GitHub Desktop.
Save akelleh/7d559d948ba80a1c5f9469be9d9ae117 to your computer and use it in GitHub Desktop.
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"%matplotlib inline"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"import pandas as pd\n",
"import numpy as np"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"First, lets do an example to show that the graph D --> Y --> S still has selection bias! There is an implicit error term on each node. We can say $U_D$ is the error term on $D$, $U_Y$ is the error term on $Y$, and $U_S$ is the error term on $S$ There is a collider of $U_Y$ and $D$ at $Y$, and sampling conditions on a descendant of the collider, $S$! We can see this by estimating the sample conditional mean, $E[Y|D=d, S=1]$, and comparing with the population quantity, $E[Y|D=d]$. We'll see they're different! "
]
},
{
"cell_type": "code",
"execution_count": 28,
"metadata": {},
"outputs": [],
"source": [
"N = 1000\n",
"\n",
"d = np.random.normal(size=N)\n",
"\n",
"epsilon = np.random.normal(size=N)\n",
"y = d + epsilon\n",
"\n",
"p_s = 1. / (1. + np.exp(-5.*y))\n",
"s = np.random.binomial(1., p=p_s)\n",
"\n",
"X = pd.DataFrame({'D': d, 'Y': y, 'S': s})\n",
"X_sampled = X[X['S'] == 1]"
]
},
{
"cell_type": "code",
"execution_count": 29,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>D</th>\n",
" <th>S</th>\n",
" <th>Y</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>0</th>\n",
" <td>0.644067</td>\n",
" <td>1</td>\n",
" <td>0.764798</td>\n",
" </tr>\n",
" <tr>\n",
" <th>1</th>\n",
" <td>-0.773107</td>\n",
" <td>1</td>\n",
" <td>1.145900</td>\n",
" </tr>\n",
" <tr>\n",
" <th>2</th>\n",
" <td>-1.705162</td>\n",
" <td>0</td>\n",
" <td>-5.040828</td>\n",
" </tr>\n",
" <tr>\n",
" <th>3</th>\n",
" <td>-0.610889</td>\n",
" <td>0</td>\n",
" <td>-1.376712</td>\n",
" </tr>\n",
" <tr>\n",
" <th>4</th>\n",
" <td>-0.621720</td>\n",
" <td>0</td>\n",
" <td>0.282031</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" D S Y\n",
"0 0.644067 1 0.764798\n",
"1 -0.773107 1 1.145900\n",
"2 -1.705162 0 -5.040828\n",
"3 -0.610889 0 -1.376712\n",
"4 -0.621720 0 0.282031"
]
},
"execution_count": 29,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"X.head()"
]
},
{
"cell_type": "code",
"execution_count": 30,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>D</th>\n",
" <th>S</th>\n",
" <th>Y</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>0</th>\n",
" <td>0.644067</td>\n",
" <td>1</td>\n",
" <td>0.764798</td>\n",
" </tr>\n",
" <tr>\n",
" <th>1</th>\n",
" <td>-0.773107</td>\n",
" <td>1</td>\n",
" <td>1.145900</td>\n",
" </tr>\n",
" <tr>\n",
" <th>6</th>\n",
" <td>-0.332474</td>\n",
" <td>1</td>\n",
" <td>0.052216</td>\n",
" </tr>\n",
" <tr>\n",
" <th>7</th>\n",
" <td>-0.337974</td>\n",
" <td>1</td>\n",
" <td>1.451546</td>\n",
" </tr>\n",
" <tr>\n",
" <th>8</th>\n",
" <td>0.037773</td>\n",
" <td>1</td>\n",
" <td>0.056914</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" D S Y\n",
"0 0.644067 1 0.764798\n",
"1 -0.773107 1 1.145900\n",
"6 -0.332474 1 0.052216\n",
"7 -0.337974 1 1.451546\n",
"8 0.037773 1 0.056914"
]
},
"execution_count": 30,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"X_sampled.head()"
]
},
{
"cell_type": "code",
"execution_count": 31,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"<matplotlib.axes._subplots.AxesSubplot at 0x7f114fab24e0>"
]
},
"execution_count": 31,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAEKCAYAAAAGvn7fAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJztvXtwW+d5Jv5+uBEEeCd4kUjxIoqUrFscW5Rtydc4F9mxkzrrbNqmrd124/nNtJk2bbaXZLrt7EzbX2azbXfS3enIbTfpL+mktROtHUeWb3EcO44lSrJkWVeKoiiBVwAkQBLEhQDO749H736HIECCJHgB+T4zGhDAwXe+c5w833ue732fVxmGQQKBQCBYP7Cs9gQEAoFAkF8IsQsEAsE6gxC7QCAQrDMIsQsEAsE6gxC7QCAQrDMIsQsEAsE6gxC7QCAQrDMIsQsEAsE6gxC7QCAQrDPYVuOkHo/HaGlpWY1TCwQCQcHi1KlTfsMwauY7blWIvaWlhU6ePLkapxYIBIKChVKqL5fjRIoRCASCdQYhdoFAIFhnEGIXCASCdYZV0dgFAoFgtTA9PU1er5ei0ehqTyUrnE4nNTY2kt1uX9TvhdgFAsGGgtfrpdLSUmppaSGl1GpPZxYMw6BAIEBer5daW1sXNYZIMQKBYEMhGo1SdXX1miR1IiKlFFVXVy/piUKIXSAQbDisVVJnLHV+IsUIBBsUXi9RVxeRz0dUU0PU2UnU2LjasxLkAxKxCwQbEF4v0QsvEE1NEdXV4fWFF/C5YHlhGAbde++99PLLL//fz5577jk6dOhQ3s4hEbtAsAHR1UVUUUFUVob3/NrVJVF7OvL9ZKOUon/4h3+gz3/+8/TQQw9RIpGgr33ta3Ts2LG8zVmIXSDYgPD5EKmbUVJCNDy8OvNZq+Anm4oK3K/JSbz/7GeXRu67d++mxx9/nL7xjW9QOBym3/iN36C2tra8zVuIXSDYgKipAUlxpE6E9zXz2kttLCznk82f//mf0x133EEOhyPv3llC7ALBBkRnJyJPIkTqk5NEwSDRAw+s7rzWGpbzycbtdtMXvvAFKikpoaKioqUPaIJsngoEGxCNjZATXC6QlMu1dHlhPYKfbMzI55ONxWIhiyX/NCwRu0CwQdHYKEQ+Hwr1yUYidoFAIMiCQn2yyVvErpSyEtFJIuo3DOOxfI0rEAgEq4nlfLL5i7/4i2UZN58R++8R0cU8jicQCASCRSAvxK6UaiSiTxPRP+ZjPIFAIBAsHvmK2P+OiP6IiFJ5Gk8gEAiWDYZhrPYU5sRS57dkjV0p9RgRjRiGcUop9eAcxz1DRM8QETU1NS31tAKB4BbEzGthcDqdFAgE1qx1L/uxO53ORY+hlroyKKX+moh+nYgSROQkojIi+qFhGL+W7Tf79u0z8l1pJRBsRJhL3s3peGsxc2OtLECF3EFJKXXKMIx98/1+ycSedtIHieir82XFCLELBPnBkSNwZjRbA4yPIy3viScWP26+SbiQFqC1jFyJXfLYBYIChs8HojSjpASfLxbLYelr9lyxWPBaUYHPBflHXondMIyfSg67QLByWI6S9+Ug4eVYgATZIZYCAkEBYzlK3nM1vlqIXCNukisLkWIEggLGcpS85/IUsFC5prMTC874OFEqhddgEJ8L8g+J2AWCAkc+S969XqJAgOiNN0DYe/YQOZ2znwIW6lPOC1BXFxagmhqMJxunywMhdoGgwLHUDBb+/eXLRNeuEe3eTfTxjxOdOweCf+ih2U8Bi/EpFzfJlYMQu0BQwFhq6zbz74NBIpuN6MIForvuInr4YZ06mT7WcmrmayXfvZAhGrtAsIbh9SJX/fBhvKZr2EvNYDH/fmKCqKqKyO0m6u7G99kyV5ZLM1+OVMuNCCF2gWCNIheSW2oaoc9HFI0S/eIXkGE++IAoHicKhfB9tih8uXzKJd89PxApRiBYo8hlgzKTJNLXR9Tfjyh/PilDKaK33yaqriZqbyf68EOis2exacpReLbUyeXQzJezx+hGgkTsAsEaRS7ReLok0ttL9POfEzU05CZlKEXEriKlpUTbthFZrVpbX+mS/+XuMbpRIMQuEKxRpJOc30/05ptEp09rvT1dEunvJzp4kKi1NTcpI5Uiuv9+oqIiLBAeD9Gv/RrRfffBa2alNy0l3z0/EClGIFijMFeVDgwQHT1KFIsRffSjRDdv4jOOqJmADx9eWNXo6dPIU//IR/B5dzfRu++C4HnhWElIvnt+IMQuEKxRMMkdO0b0b/9GlEgQVVaC1EMhov37ZxcE5ZKGaE5x7Owk+tnPsGg4HETl5Uh5bGhYWNpkPiH57kuHSDECwRpGYyM0cLudqLkZBG2xEA0OoqAoPfulsxM6+xtvEL38Ml57e2dKGeZN2dpaogcfJIpEiK5fR5QcjRINDRElk5KNUqiQiF0gWOM4fRrSiFL4V1QEsr94ERWi6UhvCpT+Pj3zxONB/nooRNTSAmkmGiU6fx6Enw1SSLR2IcQuEKxxGAZRfT0kGCJIJvE40fT07E3Fri6Q8969+rPx8flTJEdHsQHrcuG9ywVS7+vDRm06eedS8SrEv3oQKUYgWOPYtw9E3tQE/XtsDBH144/PJsrFpEiOjyPF0eNBeqRh4DUUAjlnKpCar5BIKkhXFxKxCwRpWGykuVwR6qFD0L59PhB7YyPGf+ih2dF0LpunmTJPHnsMEfrQEBaO8nJIMjU1mQuk5iskWqj7oyC/EGIXCExYrKmW10v0ne+A8OJxyCUXLhA99dTSiayxEeOYF43GRqITJ2bPs7mZ6K23kEFTV0e0aRM8YGprZ1eimufF171rl27Ycf060Y4dsBsIhUD2bW1E4fD8C4hUkK4uhNgFAhOOHUMu9/Q0iKy9XUsMcxH0sWNEV65AzqisBLG9/jrRpUtEn/vc0qP3dCI+cmR2RBwIEL34IlIVz5xB1ozFArJuadGEzQsV0czFYv9+EPzwMH5ntRL98If4rqUFOfRvv40CqPk6N5mJ3+/HPR0eXr38+I0G0dgFglvwepEeaLGANGMxouPHoWfPZ6p18iT8VlwuENqNG5AyRkaWR1/OpKUPDiKyHhwk2roVUo3LBYKPx7UWnkwS/dVfEf3n/4zo3mrFHE+cAGE/+iiu3eEgKi5G9H/lCiJ/w0CWzXwmYKzj9/YSvfeetgTm/HjR2pcXErELBLfQ1QX5QCmQIGeInDs3fw9Rs+fKwABI3TBAkEvVlzNp95mkkOFhELjbreduteJfdzeiZb8faYzd3US33455d3XBf928+VlRgWvYvRsLRTCIJ4JHH8XCQDR3IRET/7e+hacfNhnzeGZn6QjyDyF2geAWfD64GjK5FRcja2RkZH6vkjvugPmWxQIN2m4H8d52G76fT182kzfnq6dSGG9oCN4vZi19/35E2Dz25CQiYo6yGTYbiJ1teLu78b6oCMdZLPrzu+7Sc6yrgxQVi0FnT6VA7k6nXjTmQ2MjNPkDB/R5crkXgqVDpBiB4BZqakBcd90F4hsbQ9T90EPzR5ePPAJnxGQSJBiLYePyzjvx/VwOhebUQKsVC8Tbb+PvDz8kunp1ppRSUYHfpEshX/oSyHh0FHOYmsLxTicWmlQKxyaTINxoFOcvLgbx8xz5aaC9HYvU1BQyZuz2hRtyiVvj6kAidoHgFnhDsKIC5M4bgo88Mv9vGxuJnn4aUfemTUQ9PZAxqqrm9zU3pwaePw/ZgghjTE9jM5alFCId8WaTQp59FvJJXR0ie86K4c3LhgZY9B4/juNTKU3aPEezl8y5c/jtww8j9XIhEsp8m6yC5YEyWBhcQezbt884efLkip9XIJgP+cpFzyatpFdvdnUR/eAHeL99OzZhKyvxpBAMIgKPRvEEcOgQxmav9CeeWPg1mNM5o1GQ9sgInkoeeWR5qkalAjV/UEqdMgxj33zHScQuWLdYDKHky1mQx/jXf0UKYlERtOrt27G5yhp5RQXRli2QQo4fhybO/iycbvnTn4LsR0ZmRs/Z0gazXQPfD87aqaxE5LzcRCtujSsPidgF6xLmyNQsAayUDa3XS/TtbyN9sqgI/yYmIIN0dqIhRkMDNkb7+/FaWwuJhMv6d+/Gd729kGnGx7GJumcPdPOFXM9C7sdcxw4OEj3/vJ7/k09KE4yVhETsgg2NuUra+XU5pYGuLqQW2mw4t1KIxC9dQuQ9PIzPi4oQmVdX4/PJSaQUjo0Rvf8+yP4zn0GkHo0ip/z0afymvn5m2mBX10zSve8+aPQ+H/T6hobcSvxffhma/tgY5uN2Y6/gf/0vjFNdDd+aYJDom98k+upXhdzXGoTYBesS2UraL12CFLJQy4BsyCb3sLVAeTlek0kQfTKJTJNAABE4V6kqhWi8vR3ZLUeOzCTisTH8JhZDxJ5uq9vVRfRf/yuIXCk4Qf74x0T/8T+iO9Lx4yDi0tLZm7Dp1/Pmm9Dwh4aQiTM5iQXolVdwfby5y6/PPy/EvtYg6Y6CdYlsaXajo3O7Ei4EczkY1tQgp7yqCiTMBMpEWVmJz4NBRMSpFLxlON87vbKUyTWZ1MVTViuuh4jon/8ZOr3TiYXk8mXo6P/7fyNdsq5OFyqZx0xPO+zqwlPCyAjSIMvK8Do8jDmGwzOPr6jAE4JgbUGIXbAuka0pcmXl/La2uYIli1deIfqXfyF69VW8f/llnN/jAck2NiKyDgYRcVssMytSYzGQe0UF5ko0c2Hy+0Hg166BcEMhLCLJJK6HCIuCw4HmGxydGwZ+d+QISD2ZRBQ+MgLt/8c/xtjm8n4u0gqF8Hv+FwziicLvn3kPgkE8WQjWFoTYBesS2bxMtm/PT8EMSxaTk9oPZmgI7998E8c8/TQMs+x2RO719UT33gsfF8OAvJJK4e9IBK8DA3Bh9PvhrsheKw4HForychh8BQLQ269fB3FzU4zBQSwcNhsWjEgE83vtNaLNmxH1v/465vfww7gvZu8WLtLaswfzCYfxuncv0T33YHEIBPRrIIANVMHagmjsgnWLbGl2CymYyZaP3tMDiWJkBERYVAQiHRkBcXd1Ic/8S1/COM8+i4rS0lIc/+GHIGbDwGsspu1+Wfvn9neJBLT3YBDf9fcj3720FJWtN29iXrxQFBVpXd/pxHejo4jq77xzpnbP4E1ULijq6IBez5F+Rwf+/pM/QVXsjRsY5zd/U/T1tQghdsGGQqYmEw88kD3vm9P+rFYQmmEQ3X8/IupYDGPU18+ULPbsmS3tGAayVHp6oK3X1+vv2tshfcRiWCxGRyHpXL0KV8W9e/G5zUZ06hQ+Ly5GURH7vtfWItqfmkJapVIgdasVTxJ1dbDv7e3FU4sZ5k1U8/2JRDCXykrk2vPG8GOP5fe/iSD/EGIXbDjkWjAzV6l/XR1IvLpaSxZWKzJQMhllKYWURfZ5b24GoY+O6gh9716Q8vHjiLxDIYx7/DiidpcLJGwYWByuXYM/jduNce12EPDICCJ2tgooKUF7vf5+EH8uHZakoKiwsWRiV0ptIaJ/IaI6IjKI6LBhGP9jqeMKBKsNc8pkKDSz1H/fPmjfJSX4nCULlkzM0o7XC7IdG8O/Y8fwvqQE6Yi3346o+oMPQMh79+I49kQPhfDK0gxvtPp8+pwTE5hrcTFeP/hAb9QeOADTr0AA5wwGMa/FereIRcDaRz4i9gQR/aFhGKeVUqVEdEop9ZphGBfyMLZAQETLSybZxjZ7npeXzyz193ggbfT3g9jHxrBB6nKBgI8e1WN1daED0cgIMmhsNhwbi6FD0eXL2qQrlcImKEsqjY0g51AI8g/nqadSWASKiiDXKIVOTefPg9QdDsyloQH6OGfRtLbmLkVlu1eLaR0oWFksmdgNwxgkosFbf08opS4SUQMRCbEL8oLlJJO5xjY7E7a1IdL2+yF3vPEGSPHLX55Z+fnssyDQmhrIJQMDIOjbboOkUleHhcIwMJbfj88/9jFo5hcuaAIvLcUiMjwM4nY4oM1Ho5rUw2G8/8xndJ47zyccxr9Ll3A+mw0LEFFmA7FcFk9pUl0YyGu6o1KqhYg+SkTH8zmuYGPDTCYLLSryepEOePgwXtNbss01tjll0ueDXr11q86DN9sseb0gdZsNtr3T0yDpZFKX5o+OQkYhAuHG4xjL6URE7/NBNtmxAxo8b6ROT4PUlcJnnJJYWgo/maefJvrVX4U0dP06jmPt327Hk0JdHa4rW2u6uYqtzMjUkm+xdQCC5UPeNk+VUiVE9AMi+n3DMMYzfP8MET1DRNTU1JSv0wo2ABbb8T6XSH++sXkj8cgRROrmTUdzi7euLmxw1tSAWHnzdGAAskswiIi6vx8LyNQU5JxweGYru0AA5DsygnF8Puj327Zh/oODkFO2b4c01N6uPdKtVlxrIqHz3i9eREQ/Oop9gelpbLr+5V8SffKTOirPNRLP1JJPGmesPeSF2JVSdgKpf88wjB9mOsYwjMNEdJgI7o75OK9gY2CxZJILWeUyttcLbZybXHPvTvMCcPkyyPvmTRyzaRM+P3MGGndLC7T4wUEQOcsiSkGmmZoC8ff04HcPPYSova+P6Ec/QuRdW4tFyDBwvvFxLBSvv44o3+fDdblc0NnZbsDlAqG/8w6Kp9rbUazEUflnP5v74imNMwoD+ciKUUT0T0R00TCMv1n6lASCmVgsmeRCVvONzbp5T48uUHr9dUgdZWUg7K4uRMElJYjAJyeJzp7F36WliJTPnUOk/dGPQlPv6wOR796NMS5c0Fkwjz6KY4nweu+9qBzduxeLxuXL0M1378Ym6sQEcuwrK1H4VFIC0o7HEcUXF+PvqSl8NzSkrQGqqmCBoBRIn9Mx29sxl/TFcyF1AILVQz409oNE9OtE9DGl1Jlb/x7Nw7gCARFltweYj0xy6bc519hm3bypCdF4by8iZa8X0bLbjWMaGvDb5mZ85vNBAnnsMUTa09NYDNxuoq9/HVa3TU0g9L4+/L6pCWTKOjwjEkGEXVGB7Jj+fpwnGARp8yLDG6fDw1hUkkkQt1J42ojH8ZRw8ybm4/fju5de0ja9RUU431tv4VqlqrQwkY+smHeISOVhLgJBViymaCbXSD/b2Gbd3OdDZN3TA+Jzu4l27tQe6ZEI+qR2d0M2mZwEoXd0YCxOl+S0RXZZrK7GGN3dGD8chp7/0EPQ20MhROc7dsCrhQgpjR4PIvctWxCJB4MY9xOf0Lq8ywViHxxEVD89jd8XF+M8585hsUgmQfoPPoh5hEL4vL5+7iYcku64diGVp4J1i6XKBizlRKMgz/JykGIige9HRyGfbN6Mzc6DB7XX+Rtv6HH8fkTZbObV06OfGnbtwnlu3sQ5du6EnPLcc5Baysogt4yMYByPB4tBMIjFKhTShUpVVdDfGxq0WVhfnzYDq6rCnF0uzCMeRzOPnTvxt8ej559KZd6clnTHwoAQu2DdgvOyr1wBofFmItHcPUE5j1spbIKeOKG7Gw0PQ67YtEkXEbW1IUIeH9dPBjU1INneXmjqN2+CkF0uRMVTU4job97EGGVl+N3p0yBilwueMAcOED3+ODT7s2cRye/diz6q/AQxPQ0i3rULr9u3I5rfuRNjDA3NtCyORvGvpgbnYTvhY8fm1teJ9NPA+fNYVMrLcf3pPu2C1YUQu2BdgiWDZBIbm1YrCMzl0s2kvV5N4o2Nurk0SwwjI5r0yssRdVutujCInRi9Xrg4er36yeCppyCB/NmfIUPFbkfK4q5dmEMoBFLmTU6nEwtDMqnPFw5jMWhvx3t+8ti7Fy6Nb7+NdMapKRD6e+9hPtev4/hf+iV9P/7P/8H1eTy4ptpanDMcxpjbt8/U12tq8NnhwzOLlSwWop/9DONUVGCB+NnPYHAmWDsQYhcUDBZiK8CSwfnziFJdLu2ZXl+PDc8HHtAk/rd/q423DANk2t6u+4fW1uK4iQlIIhcuoFhp2zbtfNjZqef43e+C0IeHQYJOJ+Z9/DjI8coVkP2OHSBir1fb/yaTuIZIBJH500/ju099SlsUDAzoOQ4OgmBtNsxxYAAbsSzdECHyv3IFqZU2G845NoY57N+PzBvW161WXPftt+snENbRDQNPMkQzXw1JYF5TEGIXFAQWumnHksHp0yBStxvySTSqo2W7XTeSvnQJkXNbG35//jyILxwGsY6OYrMxHofOPj4O/ToeB4mmzzEYBMGOjoKoXS6cb2hIOzVGIpBMWlrwN3up2+04pqICvz97FgTe0THzHrz5pu6oVFKCeXHDDosFRM3E7nQiQ6enB08EO3dib+AXvwDpm/X1n/8c42bS0Q0DTx3vvYfN3epqorvvFmJfaxBiFxQEFrppx5IBl+KPj4OY2Z2R/cljMbS04w1RpxOvXPXpcmER6e7GbycmtHTicIDcf/AD/GZ0VDexGBwEOZaWgni5+lMpveHa3g6in5jAcT4f9PLaWp3ayF4wvAmcTGp9u6dHd0siwjkNQxM8m4pxNtAXv6jvGT/1uFz6mhlz5f8rhfO3tCD6j0Tw/uDBJf3nFeQZQuyCgsBCbQVYMqivhwTR3w8CtdkQEUcikEr27gX5cZMLLhLidnV33AFS5Nxw9kKPxTCW1QrivHEDmrfXi89TKf2kYLHgyWBiQvcpbW4GebtcWBwsFpB7WxsINx7XLek+9SksXt/7ni6EqqjQlr6Tk/hnteIfz8vhmJkNRDRbyiLCU0AggCcMnw/Xku76wRvCgYCOzs2vShKe1xSk56mgIJBLsZEZ3LGI5QXDAKFaLNgULC4GKQ4M4G8iyBN2uybwzZt1xDs6CpKcmABZszd7LAai/+lPQdCsN09MzGym7ffjNZVCVB+Pa7/05mZo7Z/+NOYRCOC4ykro+EzAo6M6P91iQVWqxYLzcI56JKJTG++4g+iZZ7STY7rJ17e/jUwYrxe2Bb29eFq47z7sIfT2zmwE3tmJ9/ffD8mI/W/uv1834RasDUjELigIzFdslClV0elEUU8opDVunw+RscuFKN7vR5TMtris33/4oe47Ojmpi44sFt1UI5nE98PD+MwwQJgs6/CiQIRji4vx+0AApG23g9y3bcM8Dx1Cef/p0xhr3z5t8EWE3wSDOAe3ymM9v6wMhM4+NS0tOqL2eom+9S1ca10dJCAipEL6/bhP7e2Yb0eHznjp78d35vz/mhqcn4uliED86R2jBKsLIXbBmkJ67nllJSLszs7sxUZeL9F3vgPSDgSw6ckl/HffDSJnT5XbbtMOiw0NIHWPB+c4dw6kGgyCFJ1ORK7l5fB4cbt1FE6EMW02Lc9YLFggmNiJdEPo6mqc3+nUG547doBQp6Y0CXPza74PR49qj5rr13Xj62gUn3FhUlOT9pfp70dxUkkJxjhxAvemvh6/PX4c866sxJPD9DTGiET0hmtzM+b6zDMz//uICVhhQIhdsGaQKffc6wXZPPccinMeeWT2ZumxY1gIiopAmuyLkkoh8q6qwnehEMiROxLV1iK6TSRA8tu3Q5I4eRIEGAyC7Hp7ddZMIqGjc8stIdPhADkmElpbJ9LaM4/jdmvJIhhE1Wd9PWwAnE4sTrW1WNB6elB56nbrJtp79mDxmZzExu/584j6d+7EJuzJk3oPwe3GNT37LMZhUufI+vx5/K68HO8jEcyB2+Zlk7nEBKwwIMQuWDNIzz1PJJDpUVKCyPr8eZBreorjyZOIiG/ehERRVASyDQZRtRmJIGq/dg3RbDyO6PYrX8HCMTWls2yuXsX30SgicyZsll44Ane5dHu6WAy/TaVmpv0phblz0RG7LhYVYUyLBdfHC9GVK4isOavlxz/GHBIJkO7kJLJPbtyAZcGOHbi+7m6Mffo0nla2bIHXemsrUiUHB0HAr7yCufBGbiAAfxgiRPGRCL5jTT1bFC7Nrtc+hNgFawac+cImVFeugJATCUSgY2MzuxsxeMNyakqnKXKWSDQK4t+5E2MnEvgumSR6/nlE8zt26LHGx3VlKRO1xYLjuYMRpywqhWNtNhB9IKBL6y0W/JueBmlzJyXOXDEM/OvrIzp1CgtLdTXGvXkTEXs4jOidnyqIcE86O6GNP/QQzlFVBXJ3uyGjfOELetO4pgYL2tgY5JXRUdwfiwVRvMOBRYE9ajhTR6LwwoYQu2DNgDNfLBZIDteugUhra/H3tWsgn+JikA5ni9xxB4pq4nFIEpyuWF2NXPa+PhTisJVtfT2I9MIFEP/p0zgHb5ROTWHc0lKQNVeCsqbt8YCcOT++uBjyh92OvznKZ22c5ZdUSm+42u342+8n+vd/x7Xv2YOxe3r0gqIU5lRaiuPdbtyb0tLZ/uk8byJcL+fPDwzgXng8+J3HAyIvLtZ2xVu2oBl2NjJfzmbigvxD0h0FawadndDABwd1NsnQELTo117T+eCGQfTNb+qep488AmLmoqCxMRBWayuI+eJFEC0X8vT3a0fF+nr87v33QfS1tSDtVEq3uLPZtLzDTwFKgSQ5X9zpxBOB1aqj6/TIXil8b7frSJ51++FhLE5DQ4igLRZII+yzPjoKoj57lujdd3EPzP7pR4/i3l28qDeSHQ781unEufn4u+5C9B4I5PbfJdd+qIK1AyF2wZpBYyOItbER8kIqpUnQ4UDEPT4O6aS6GlIKo6wMpHPbbbqRs8MB0mTdPRYDyVksIPeSEmSTsJxhtYKs77xTpz2ytFJVpas0i4owp2QSPiubN2Oc5mYtzfB5lEKUzc0zkkmt2Suli4nsdiwY77yjNXWrVY8xMaEXu9JSHLt3L44LBLQvzLZtmMPly5jjAw9go9VuR+rkPffguL4+PAHlQtZLaSYuWB2IFCNYUzAMaMfHj4MoEwmiH/5QR71ut45ib9zAb7q6EJ1PTenMD96knJhAWl8yicg2HNYbmvX1IEinEzIIR7VVVYjeKyo0iY6PY8HxeHR6Y0UFznfjhs5tJ8L5DUNH7yzNsD7PGTvFxTg3Z7IQITLnzJpt20DasRjOtXUrzltfj2sKBEDUv/gFFrN4HGPs26erY0+dgnTDlbcHDuCcH36o/d6J5rZoWGwzccHqQYhdsKbAvTfPnAGJbt4MIkskdIchImRtNDTgbyae9nYsCEQg9aEhvbFqGBiHUw9+yx4YAAAgAElEQVSLikDA4TB+F43id+XlIPXGRhzDEgibax04gM+vXoV8Y7Np6YRbzRHpjVyWZXgz1WLR3iyGATKORHQmDPcoNQycg7Nh6uqwqG3eDEJXCvMiwmtRkU5djEaxmHR1aXOwlhZE6H4/zMDa2rBwmpGNrBfbTFywehBiFywJi9lUy/Ybr1dr5KWlII/z5xGNcqOM8nIQWyBA9Ju/ifGYeDwekPS770I3r6tDYVFXF4iOW8SlUjjO6cQmotuNTValUE7/1ls6ArfbEa1brTj/pUuQMXgD1+nUVr/9/fg7Pe2RCIQei+E7txskOj6O47hRhlK4Bm5nx77pra24bpZoSkp09SrLVempi2wHwBvJu3ZhEeLNX48nd7KWoqTCgxC7YNFYTP/Lri4UzSQSuu3cwIAuemlpQVTKEoLTiejynnuIXn8d55ycJHr4YaQBEmni6emBRs0547ffDiLeuRM53EQYu7IS5z9wAJuNXImaSuFpYWwMJM5VpSUlWucfGABJnjqF8/OG7siIzrrJhKIi3Ws0mQShFhfrFnSpFAi9vh7HT0/r6tCaGmz01tWBsNkorLYWv921C+c3py6eP69NyNhxMpXCtfl8RI8+qjdZ43HdMempp2bPXYqSCg9C7IJFY6FWul4vSN1mAzlEo5A9du7UEXxdHQjpU5/CJuaVKyC1228HEcdiIKKBARDTU0/hXPv3E33ta/iupkZHvS4XZJ2PfARPAUQYY3AQzTXuvBPnu3oVJN7YCFKdnASpcr/Tigqij30MxFpWhs8uXgSpV1UhiucNUU6FZGMujv55oeBCp+lpjMtGWzYbiJOlHy5k4qeD0VEsDJ/5zOwKXH4K4tTF3/gNkDsXUxHhWswt79KfKubyVF+OoiRJoVw+CLELFo2Fbqqxn3htrU4lJALJstlUXx+0ce6nuWkTSDsQwOceDwjX50MK5KVLcBc8fRqkt3kz/pWXg3wjEejKDQ1a0w4GIclMT+Ozt97C53V18FXhyJdllVQKpPrqqzjH+fNYJGIxLAB+P+bHhJ5OkJzPzla8ySTGicUwx2hUN6w2R/xFRTgPR+9lZUS/8ivaT8aMdOLlNn1Xr2qZJxCAyRd3YWptxYLHGB9fuabUi3naE+QOIXbBorHQTTWOzCIRTepOJ4h9/378H/r734emzq6FPT1EX/0qMmOCQZAoEQjR6cRC8OqryO9mn3W/H6mAZWX4zY4duvrT5YI8YRggzTNnMF5pKcg2EMD3XEiUSunmz/x+ehrET4Qnh4kJRNZVVZiX2dSLwRupnNpos2GunM7JG7w2m14cnE4sYlyklUggUp8PHAlzwZTfj2u67z7tFnn0KAjV79ct8UpLcd9XAgt92hMsDELsgkVjoZtqZvmFCKQzNobPAgHo4CUlINjeXp0l8td/jd/U1elmzZz2NzyMMSIRnVsejyNSra5GlN/RgQ1FNuHilnYeD0iZG2yEwxgjGASpcs65maQdDowzPY0ovroapNvfj3lYLDOP5wImh0NLKzwvu13bFXAnJG64wZq/y4WFr6oKkTqT3lwb0BwJ79iBz4LB2ZEwPx1duKD95XmD2utdfnKVFMrlhRQoCRYN3lTjsnSXa+5H6c5OEB03tODMjy1btFe5UojSKyt1/0/WsoeGcDzr1X19kA/sdpB8PK790i9dwuJw332wHKirw+IxOqrT/2IxTa7BIIiUCZ17kLJOToRFIBQCgXOKIdsKcGUrF0ZZLBi7qQnWwVztybYEXOg0OoqxOXedrQSmp3WXpnvvJfov/0VbKMxVCWqOhEdHIRudOgU/dnPxUWcnctk5/ZKvd/fulSk8WmjjFMHCIBG7YElYyKaaObvC6QTZBQIg2nhca8IWCyont2zB7xwORLk2G0iVZYtUCoRvt+s+pJEIiNPhgA3vBx8Qffe7uiPRjh2QaX70I03OgQBea2pAfiUlOr88HezoyNG8zabfs9zidmNOFRV4Krj7bqRTtrRgcYhEMks2XGDETyq1tXj/5JM6Gu/qwpON0wl9nCtBiWZuQPv9+imFn2zMGnZjI7JrgkH8Ky8HqVdVrUzULCmUywshdkHekEuWQ/pCcPgwZJWuLpD3xYsg2YEBXebf0IAI9s47tYvhpUsg2ZISnC8aBRFyn894HJ93d+v2dKEQpIef/ETLMjYb/m5u1tkr3IeUdfV0sHQyPa2fGMwt9Ww2XQkaCsH6oKUFufNWKzZoOXedo3TW1dlfxm7H9ezcifu6aRNa2fX2wp89FiN66SXk6dfV6dz5O+7Qzbfdbl2FW1c32xmzo2OmZTERrmclomZJoVxeCLEL8oLFZjnU1CArhXOsrVaQIcsjBw4ggu3rQzS8dSvGvu023Vc0HMZ5rVYQptuNcS9cwNMASyvcLSgaxYLAG5klJSDKBx8EIZ85gyrNbOl/nNs+Oal9X/gpgv1pIhH8bbMheu7oAGm3teG33/++JnJ2f0wmdVFRRwcIlxuI+P3YIB4Y0N2awmE8CWzbpgutRkYw7tCQjtT7+vB3KjVzc3S1o2bxdV8+CLEL8oLFZjl0diIS9fu102E0CoLlXHClQPDsysht79i7ncGbjqWlIDJ2dUyltGcLEyln5oyPI5rmzdjmZvi3Z5JhGFar3gg1N+Fgfdww8HlzMxaf6mqQNcPp1FkvxcX4LpUCGVutiMJ5z+DmTdxLfoLga+eq1WQSZNzVpVMXN2/GwtbTg++am3VnJvPmqETN6xdC7IK8YLFZDo2NIK7Ll3V6X00NIk+3GyRdWambUPzkJ7qpxdQUFoTpaUTed9+NBSAQAKE1NoIQrVaci/PJGckkjk0kUL158iQqT19/Xbe3ywY+p1J45bRI3oDduhXVsqdOgaTDYZB8fz82LdmvnSUdToFk07GeHszbZkMk7XbjetmT3W7HeMmkXpwaG7Egvf8+0R/8AYrBqqtx/9I3R5m8c4mapZCo8CDELsgLlmIUZRiIMsvKEL3G4yC0zZuJvvEN3az63XdBfCzNcNaMxTKzwxAbZT30ENH3vqc1cN5wZT2aI+2pKZDWBx/M7Jc6F1wujOnz4T2fo6hIy0mXL+M9l/afOweNfHoa2TLRKK4lGsV41dX4bmxMR/OMffswHmfqWCz4DXu/1NfrrKK6OtyzTJujqRQ2X3MlaSkkKkxIuqMgL+jsBIHwRiX3zeQUPSKQxJEj2DA9ckSn35WWat+TixfxWlurLQA4WrTbsZHq8+luR+Pj2qr33Dmcc2oKRU333gsrAvY9ZzLnTkaJBMa0WqGpx+NYEFgSygSHQ7frCwRAeOzWWFqKDd94HPOKRECGxcXY8AyH9Sbvli0g282bUQHKVans+T45iX+jo5j/oUMo4iouxjXGYlgIiPCbTZtw3eEwLIh9Puj0u3fjt+3teHr4p3/Cfbdac2uYIV7shQkhdkFeMF9O+1y51y0tKDriopqKCrxvacFv2aiqvFynFrLXOeeus1TR20v03nuQbIiIfvd3ER2zNwtvWBLpHqZsNcB6+Vzg83MxFHchYl/1oiIQLi9K3HP1/HlcKzfgGBvDouL1ws+djcGmpqCrT0zgb7cbi933vgdy/qVfAsFzU5EHH8Q9SyZ1dyS2Z+DFlu9JdzeO2bwZxByPz0/SPp+usmVwJpJg7UKkGEHeMJdem2lzNRAg+qu/guY8MgIiY1JsbNQRaU0NIuWqKmjPnDtOpCs7h4ZAgpOT+O7aNUS13POU+5KyTzpXg5r90tlTZS4zLCJN7kR6kamq0s01uBn34CBImWUV9pMJBrUdMG+eVleDzJXSx0cieF9VpfcM3G5U4qYvmBUVMzNbOjpwzycmkM/ucGCsPXu0j053N+Sat9/OLs2IF3thIi/ErpQ6RET/g4isRPSPhmH8v/kYV7B+4POB1M6fB/EphdS9wUGQVVsbyDkchnzQ1IQ0yJMnQfo3b4JEecOQG18QgchsNkgcbHVbUoIo9epVEJHDAS06mcT5uZvR6KjOaInFcrsWthlwOnVkzro3O0GWlmLOnJrIMhD7uhPhMzNhVlTgnpSV4VpGRjCnjg4sAvE4yPgv/5Lok5/UJJye2dLRgVz5igpE9VeuYCGpqdGdmpxOLAr9/bhf2fTz1U6JFCwOSyZ2pZSViP4nEX2CiLxE1KWUetEwjAtLHVuwfmCxIHr2eEA4586BpBob8X56GuTOG5snToB03G5E3IaBv8fHcWxZmfZ24Y3O4WGdRRKL6Y1KTm+MRBAZFxdjgSgt1SX+2TT1TDBb73Krvp07QbzsW8Ot6Ox2fMZz5ObW3DCbG1aPjOBe8FPFzZs6z76vD2NdvYoFj49JryRlHDky8+mIO0cR6Q5UhoEnCI8ncwWrOWtGUiILD/mI2PcT0VXDMK4RESmlvk9EnyUiIXbB/0W6th2JgOQCARBcIADSsNmQrufzaVtbrtAMh/VGZTKJDUO24GWNeWoKhHjliu6Byu3mEglExJxGGIvhnH4/Iuts5G6x6CYZDE6dNAx8x9WwExOIbJn8Obpn10ZuasHmX4YBouYIOhzWtsHT01pGKi/HsdeuIQrn6P1P/1Q3D9m+HRF2euppezueXmIxpIOeOweSrqyEl47Ho4/NlKIqhUSFh3wQewMR3TS99xLRXXkYV1CgyJT3bBggkZ4eRLA2m05tZB+YgQF8xvnaiYQmU87i4E1GbobBpMkyDUfsk5OaeHnjNJ2MlcK5o1Gt2WcCZ9E4HFpb53NGItrwixeYyUndeDsSwbH8lMCLVFHRTG8YznZxu0GioRDuH2vjRPg7FtO6eSIB4i8txW9dLn0Pzbq4x4M8/f5+zO+BB7Qn+9TUzGsV/Xx9YMU2T5VSzxDRM0RETU1NK3VaQZ4xX7FKtrxnhwNR6T334LhXXkGUzpG706nL3W+/HU00zJo3Ezg3ymAdfXhYHzc1pQ3B+AmBNyEZnLbHG6deL8jUTO4c0af7uZhtdjnKTiSwUI2NIarm33DWTGWlrnDlAiR2jGSL4nAYmjpn+HAOv92Oa7t2DQtBUREKn/r6tG9ORQUWp6kp7FHs2oXrZemFdXGrlejLX54deYt+vj6Rj3THfiLaYnrfeOuzGTAM47BhGPsMw9hXIyFBQWKulEVGtrxnpWbmuUcikAjuvFO3u9u7F4vFgQP4LUe4nNZotWpbXiZlszxCpEmdSDsxskVucbHuUzo+rvPE2c6AI3rbrXCHI2WrFaTNDTjKyxEFcyOLJ57AIsXXODmpm0xHItg7OHgQBJ9KgdydTlwfZ9Rw5Wp5OY7nCD4QwFyJMK7Ph8WM6wSam/FdcTGi/JISjJWLnfJCbZcFhYN8ROxdRNSulGolEPovE9Gv5mFcwRpDLn4w2awFvF4Q2DvvgHiLirBp19qqu/gMDeE7rxfn6e/X1aJEIHEmW/aI4TRFIm0CxuCNUSZ6TjlkAjfLPVzeb26Uwa+cScMpkkRYjCoqcA3RKObKUTU/NSgF8v7gAyxilZV4v2kTFpSpKZDp9DSundMdy8v1YhmPY0FRCmMPDmIzdHAQ5+b5RCL4HUspueriop+vTyyZ2A3DSCilfpeIXiGkO/6zYRjnlzwzwZpDLn4w5rxnJuyeHhDZ44+jCnJyEhkrExN45abLdjtkmB/8ABFpUZGWUVgr57J5MzkTzbQVYDBZM/ErpYt+uFkHuz/G45BDKiuhU3M6Jcs2rHUXFekm1Gyx+/77iJw5dXByEtcWjeIcU1O4B3V1GKe5GQuFzYZjeHOXM224QxS37+P8fKUw/s6duNb2dtw79qBvbhYpRQDkRWM3DOMoER3Nx1iCtYu5ilVYe798GZpwYyOIzmqF/lxfDxvd0lLIGK2tILyrV3W7uPZ2PSaRttNlwmbdnHuHckERb16yXMPHcOokl+qzvMFRMJEew2YDGbMzI5/L6dRSzuSkllfGx7FwhUIg1aefRqHP6KieP19XRQXOl0ggXfHsWfy+tFTPr7ZWa+PcTaq4GNczOYl7wYvFe+8R/Yf/AFuCSATnrKzEezHoEhBJ5amAQMjPPw8ibmhAxx6zxwsjW7FKR4feML3tNpDoj34E6aWtDaTe1AQS6u4GsfPGYVub1tSJiH7xC5AgSyI2mzbuIppJ8izTsOZONJOozQ6NTPJMkOam0RwNc5TMm5hFRfjt1JTuahQK6fOy/4xSWKAeeAD30VzZGo9rfX9yEk8DmzbpVn1cPBUM4ri2NrhU9vbis/5+7e9usWD+ySQ2UH/nd6DvCwTpEK+YDY6uLqJvfhNk19SE129+M7N/SLbNNtbEecO0tRXRI1sCsHMidxQiAjH19BCdPg2rXL8fn4dCuiNQKDS7mTSDydPpnKmlm8meSEfrlZXafZHPMzmpC4rYctdigVbNiwUTOPu5c4u94mJIH7ffjkYXH3yAe3fvvTpirq7W7fqmp3V3J4dDLyBEugFIWRmqbjdvxhwOHtQOlDwvux2ZRcXF+THiymbMJihsSMS+wfH88yAgJmF+ff75mc2T50pxPHp0tvaeSBC98YaulOzvR1RbWwsSHRmBFe3kJHKyz54levRR3SvUZgPxmzc/05FIaMkiGgVpszUA561z1D81heO52Id92PlJgF/ZnyUW0xuTJSUYm3PkN2/GU4jLhe/NBVd2O4y4WEYZHtZ59SwXBYN4mikt1d42ZWUwPevowO++9CXc9+5u3Hs2GKuvx99VVUs34hJL3vULIfYNjv5+kK8ZFRVwHCTK7f/86do7b5qyZ/jUFCJjlkaGh0G4TPSdndDmX3yR6LHHMN4HH4BY5+pkxBIHt5NzOnV3Iy7d58Il1uJZs7fbdWEPb8DyhizLOhUVIN1Nm0DmLS1Er76qz8NR99gYeo1+8pO6h2hrK44tK9O2xN3dSO/kRcTt1kVM+/bhv4NZJ+fXr3xFNyGZnMRexSc+sfRCosV2vRKsfQixb3A0NIDcOFInwvuyMjyav/qqTk3M5ieSrr2fPQsC3bQJC8PoqG4YzVa8N24gxS8aBRmWloJkOjow7rFj2rAqE9jV0WLRlZ4TEzrK5hZ4nN9NBHLkvHiWcriAibV5TqdkjX//fpA2EeQjTlWcmMB8i4rwBPDpT+N63ngDC+CePfjsww91URGbf/EGr92O+3nXXUR//MeZr9PrBYm/+65umOHxYEH+7d/O+T9zRiy265Vg7UOIfYPjySehqROBqIJBGFDt2KEtbi0WyCV33aU3Ps3/5083ivL7tXzR1oaIWClkgpSX4zdjY5rwOZ0vkYDHCxEI/uJFnd6Y7pPOxGve5GTwRiU3lzanRjKh84YsLwB8DKddsidNdTV+5/UitfAzn8FTyFtv4Vp37oSE1NeH6x4awp7Bt7+NheyLX8QxL72ExYBTJS0WbDSHwzjHkSOZpS6fDzp+UxMi/lBIF30tNaoWS971CyH2DY7OTnQbev55RNENDXoDkAkkFgP5ckZLpv/zmwtdenpAklzR6XKBkCYmdErjiy+CUNmHXSlEttev419Pj/ZZIZpZqMRkXlWl34+OgpTZqyWVQjSaSOjG06y52+04xu8HgfOTAfvXJJNYHG6/Hdd94waO3bYN8y0vR7qhw6ELjHp6cE1cWVtcDFI+dozo3/8d45SXY9Hs78d9nJrCPRsawvGZpC4mX49Hm3WNj2t9f6n/7cVSYH1CsmIE1NmJ3qLf/S5e2YeFCEQcDoOs2BIgveWdGV4vUvr6+mZ2OEomoR9XVYEQKytBiH4/CDwUQlXqq68iah8a0r/L1NUoGtV+4lYrJBKu0HQ6QYrshc4bm2xdEA6D7Fkr5wbUTufM9nSjo3i/d6+uFLVaMc7x4xjf58N8f/ITXI/NNjMN8sYN/btAAJITV8D290Oaam3N3noul5aDi4VYCqxfSMQumAXzI7rHA3J/801dVPPkk5n/z88brdXVIPWrV0FC27ZBq7bbtaf37t3YUOzpAdFOToKAuXiIo1JuWWfumsRpjbwJyn0+6+v174JBNOmYmNDz44ifJZqqKiw2qRQWkqkpLQnFYtpXvaNDR//d3ThPIED0r/+KTc/x8Zmkb7eD9GMxjFdVBQuApiadt+5247djY5lbz7HUtdx+6GIpsD4hxC6YBfMjejQK/+7qahTDOJ1ogrFp02xC4CyLpiY0Tq6p0S3ifD6iz38eenRjI6Lcnh5E7kzG7N0yOKj7atpsIEpzRajTCbLkVMF4HHO8dg3RcGkpxpia0m6JZlteHouvIZUCabLNQEmJ9nkhQrHQpz8NUh8fh/bvdmNMjwc+7Gx/a7HgnrElgdWK4yYnsYg0NODvZBKFWV1donML8g8hdsEsmKPEd97RZles8QYCRN/6FjZGzZt9nGVx/jw2BgcGEPHabIj6336b6MwZoocfxkJRWQmpgqsziXShEKcJsjFXIqGzWoqK8H0ioatGuXI0FNIavMWS2WedHSIHB3WzjqkpkDYXD7GHezIJ0u/rw324dEnvG4TDyAByOKCfj4xol0Z+muCiIrsdTwVcPPWpT+F69u2bbbHb24ungsOHMdehIcg1kmsuyBVC7IKM4Ed0JmvOQff7QdzT04g4zUTDEk4ohL/HxpAZwpubqRRkjQ8/xCIRi+FY3kQ12wWMjYG83W6Qb3k5fscVqkTaotfhQCTMksjIiM5/5w1VIr15arXqTc+xMUTinD3DLfR4HkT4/fnzIHIuYOK+qaOj2rtm1y6cu7dXa/Xt7dDRb9zApnBLC8770kt4fewxoo99DDLW8LBepIqLQfRvvgnib2jInm4qEKRDiH2DYb4q0nSkp8R1d4PEqqtnEw1LOGyUFQzqiJslEXOGSTg8M4uF0xE5ymZ3x3PnsEBwNWk8rqUWIhBrfz8WCpsN43KjDHPVKmfElJToJhmGgWthz3PW4c3gbk3BoCZsLnwqKcH9sFpxb7lptMuFe8OplEzUly9jfo2N2GO4eBHjP/00PjtyRP/W78cTQiKBheCxxzKnmwoE6ZCsmA2EXBplmI89cgRa+FtvIQo1a9Gctkik9XCWcHbvBvEWFUFvj8dB0JOTOKfHg3ECAcgx3NuU5Q8iEKjfr9MWy8q03zrn1rMlAFsAmDNo3G6QLBcxMVwuXeRzzz24rps3tetjJnBPVl4YOEJnPd1igXTCnYo6OvTi1t0NUr56Fedl35rWVjyJVFdjDt/6FqSXV17BmH4/Mm84pXJ8HO/9ftHgBfNDIvYNhFxLyM02Ajt2INr88ENNyg0NMxsgm4mmsZHoP/0n+K4fO4ZKzFRKSzGhECo5r1wBSQ8OzvRRZ32bJQmHAzo0R6mbNoEYR0ZAdvw0EIvp33K/0WQS0TFX1zLpc4XqjRtaF2fTsHSYq1F5kTD70CSTiLyZ5LnZRjiM87CfTVWVXgTKyrD/wL1bu7sx1oEDeJL42c+0D3trK55YSkpA8GfPYlGVXHPBXJCIfQPB58ucWpduJpXe3q61FUSyfTv6Zlqt8+dVM8H/wR8g4g6HdREOl9fb7ToqT0c8rqP1S5eIfv5z3Te0qgpE6HZrd0ez5MJSidOpPWQ2bwYhsoRUWak3arnRRjqY0Hl87q3KBU1sKMZPEmVluE6fT0f1ySTmc+edevOTSMtI169jnPp6zOsjH8H5Ll7UxmMNDXjy4RRK2TgVzAeJ2DcQci0hn8tDZK686kz6vdeLVMGqKkTIyaTuedrfn5nUibRuPTICMkulQIahkLYn4DRIsyEXE7zZ0XF0VD9pEIFM2UKYbQRY3zcvEGw3wFIOSzG1tchUUQopjJs2YfG4cgWLHMtALMeUlGABaG6eKUFxD9OtW7W05fEQ3XcfqlWHh/Hf4WMfw+ecFiqkLpgPQuwbCLmWkGdyazx7FoT0/vu6o88dd+jN12wukKzfT0yAoKurtWUuZ6JkAmvoTJREunH18DDGikbxe49H565z9Dw4qDcyuVl0OKw3cbnJBUstrNWbq1y5GQdLQnY7IutIRHeC4g5J99+PMVha4gXC48H7qSn8fvdu3bQjGoVtwW23zZS2nE7k/HNf1ZIS/WQkEowgFwixbyDkWsWYXqD0s5/pzc3r10F4bW2QR4aHkdHR1aXTAvv7dfHQjRuQF7ZsgbTwzjsgOvZwsdtnWwawLQBnkhAhsuVo3+/X3uvcN5SzVaxW3diC/WqSSUTIo6OI0sfHtV866/NM6lYrxqyowDUXFeE6Ght12mdfH2SpTZuQjx4MEj3yCJqGVFbiHCUlSKUsL9da/PAwcvgPHdL3nBfE8fGZi+1nP4vvl6viVLC+IcRe4Fho+mIuJeTmBeDtt0FWpaUgabb3DYWg+/r9OO7KFVR+GgbIf2JCdwu6dEl7sUxNgZyZxMwSCG+EsrOiwwHZw27H+bjKlFvfEc3MUXe5NEkTaWmGvW7Yi503NFnnZpMwblzNjUcqKhDtDw3NHGfTJlx3fz+eajo6cA/uuAOLXVMTiH9kBBujW7aAlDP9t5lvsRUiFywGQuwFjOXsgJNeoPTqqyA/jqDDYUTFly+D5G7e1NYBbE3LfUZDIRBeIqErQxnsn04EAud+p52dIMg33tBROBM/R9r8WleH6Jh913kj1GLR2nwqhc/tdmjdY2N6M7S4WMs03KiDfdqvX4fPTSikK099Phzb2IjFbnoax37xi9q2eHQUC83DD+sc9fnutUCQLwixFzBWogMO6+3l5bM9Vz78UKcAXrmiDbSIQH5cmNPait9xmiFnjJjb0pWUgNRLSlC89LnP4buTJ3Wkzd7sDgfIlJtt1NeDeFl+YXlnelovIsmkbmzB6ZAlJTrd0e3G52wvEI3qFMXf+i0sMMEgFgTeRA0EMHZFBYi8sVHLUvM9QS30SUsgWAiE2AsYK9EBh/X2+npoyzduaIJnZ8LhYXzmdmtfdG4GHYthnuy5breD7IuKdLVnTQ1K7X0+fBYIwESMCLJMMqkrPKNR3VSaJZdgEFKRYeC3HNmzhzv70LBlwNWrmAOP7XbDDuDKFW0rPDKCeey6/4QAABqLSURBVH3iEyDh7dtxDVevojUdZ+NMT0OC4oybXKJv6TUqWG4IsRcwsmWvxGKoGp0vCswWNaZ/vn8/PuN+niMjOE9LC7oDvfaaJnFz7jenCobDugCImzp7PPh+ehpkPTgIQq2rQ/R94QLGaG0FcdbW6swazlRRCoRaUYHvamq0oyNH9EVFSCW8dg3nYU2eq12Vwu/4NzYbrqW8nOiXfxlEPzwMV8rvfAdVotPT2me+vV1n3eQK6TUqWG4IsRcwMmWvKIU8aLYLyBYFZosa9++HLa/58xMnMM4TT+C3R45gQ5DNtCoqEL3zZiRngRCBxCsqECFbLFqD93pBxjabtu7l8UZG8HdxMTYoq6t1RWllpW6usXkz5phKgeDjcXR/OnsW5MzReCQC8uQMGt58LSnRzbAHBjBPttx1uXBNDofO82fL4FAIY7N7ZCqFeeUKftLipt/mdncCQT4gxF7AyJS9YrbXJcoeBWaLGp9/Hh2DskWTXi82UsNhRLLNzfj83DmQXUUFCJDJkf3R6+pwfEUFSNXv1w0nbDYcm0ppT3VuejE2BrJuaQGJHjqEOY6O6oId3tB9/338dts23TC6thYSEi8go6M62ueOQX19WCwGBnCO8nL87kc/Ivr4x4meegrX39pK9NGPwno4HkfUHwzimszNwOdDTQ3OeeECrp8XxkAA91eidsFSIcRe4Mhmr0s0t96eTZ/v74dnSfrnw8Mgt2efRYMMl0s3mQgEQNLsqFhUhGiaN1oTCXze2IjPuFK0tRW/i0RA9GNjOp2QUxDjcWSmVFXpRtgNDSDSGzcwzuXL+C1n1UQiWGT270e3pz/7M+21/pGPYO7RKBaO/fsxbiSC7zlLxzBwTG0t5n30KO7XxAQibH664NdshVaZ0NlJ9PLLuF+cfZNKoXhJ5BhBPiDEXiCYL4tioR3nsx3PHX7SP7dYQOo2G3T1CxdA6JyDzl4n3DSDDbHq60G2LLuMj4OkYzHo6ps364yVaFQvGCMjGLu+HlG834/I/8wZjHfhAsi8o0N7sFutIOI9e3CuWAyyTE0NyH90FJ/V1YFMt2zRHjfd3dhAZYKemtJNN/h+9fVBY9+8WRM8V5MuRGNvbISNAGfZlJdjDNbzBYKlQoi9AOD1YuPO5wNxOhwgtqee0uS+0I7z2Y5/8klo6umfOxzayVApENG774IoudNQLIbjeGPR6QQZ9/VB5nC7ZzbEINIVpdzAgl0X+bp6ekB2mzdDSuHNy0QCYw8MaIsDhwPHcw6530/03HPauqC4GOebmIDc9JWv6PO8/DJ+x+Zg4TBkJl4YOcrmbJriYix+TU14ytm2bWH/TTnLxryAjo+LHa8gPxB3xwLAsWNIxbNaQTxWK94fO6aPWWjH+WzHd3Zm/twwQDrcXaisDO9LSxG1Op14z/o5k20qpStIDQOLRDQKYtuxA4vA6Ch+v2ULro+rSrmpRksLZBuvF5+VlmIR2bsX82PnRZ8PTwHd3ch/v3ED5OlyaTvdrVvxxPHpT+tovbGR6Etf0v1WHQ4cY7XOPKatDQQeDOK62tu1TJXubjkfOjsxznwumQLBYiARewHg5Eloyi4X3rtceH/yJKxxGQutYMx2fKbPa2pAyBcu4H1xMSLtUAhzKS+HrEAEsquu1n1KrVZE69PTIHOvF3LL1JQu3d++HQsCk393N8jVbgexs0Oiy4Wo3+/Hv0BAuzOGw/g9p1qy8dfICBYMLmbiRcaMzk5o+3PJXR0deL9vn85mSaVQXbpQXTxX3x6BYDEQYi8ApNvJEmkCWwl0dSFXvatLb5KaG0iUl+NvjuBdLhCuzQZiZqMtdmlk3xW3W5N/PA6SLy0l+uADHO9yQQ5h/drl0tkybW3wcCHC76NR/G2z4XibTRco8Sbu6Ch+09iYWfKYb2Fk+aqiguiuu7RMdejQ4u6rWAkIlgtC7AUANpeyWHQWxdgY0cGDy3M+80bt2BgcGbdsga7+3nsg0fvuQ+T63nuIlgcHQXgdHdioNAy8DwZ1Y2tu8VZailePR+eOj48jMr9+HcfdcQciaKW09S9Xto6NIUp+5RXt+cKaOxuNMbgoic3GgkGcN5PkMd8G9UKjbLENEKwWhNgLAI88MtNcyuGA1vvIIzOPyweRpBcuvfgiyHDrVmjWpaVa43/oIaK774ZnzMGD2Mj0+RB5790L+eTnP0dka7Ph82gUG6FsErZ5M8blYiQubCopwUZobS3km9FRPAXs3w8pprcXsspHPwrC5tz6RAILgLmHqc0G2YSfIOrr57/uTGX+6fe3sRHvjx6dfb/FNkCwmlBG+jP+CmDfvn3GyZMnV/y8hYz5SNtMJOm+3gshkiNHZmZr/Pf/jr9TKZApe7aMjWkrWp8PGrl5boODSI9kcy67HZLNl76ko2W+pitXQNxnziDSb2kBoY+MIKOGCA0pqqrwRBAIYCEYHkbUnkxiwejr09E569/sYeNy4QmjvR3pkOn3Jv26iXQB1BNPzL6/fX1YtA4ehFyUfr/nG08gWAyUUqcMw9g333FLitiVUv+NiB4nojgR9RDRbxqGEVzKmILMmE+PXaz/SPqCceUKImRGdbXegHS7ER1zZx+3GxWnu3ZlHveBB2YT29mz+O7yZeSn796N801OolFFfb0unOKy+54evCoFe+A33sD5e3t1tgsRiL+4WPvGu1x40hgd1QVJ27dnvjfzGaql39+hIZxjaAgZO+ljroRBm0CQDUuVYl4joj81DCOhlPoGEf0pEf3x0qclWCjmIpK5zL7S5YKeHm21S4Qq1OeeA5nX10MDTybxeSqFjkkjI9C26+pAoAMDyBW/7TY9F78fjo3HjyPKZcOtCxd0mzku0KmtxXg+H76vqIDeHolofbuvD5GyUjjn1q1YJGIx6PNVVVhIBgch8djtmB/bLaST7HwFXun3l/3Z2QM+fcyFFowJBPnEkvLYDcN41TCMWw7c9B4RiXq4SmAiMYMrRl94AbJAXZ02B2Oy5yjUYsHr7t3QzDm/ur4eG6VbtoDEmpuhq/NGaDSqiTceBxEnk5BqeD5+P6Lsd94B8V67hnOzTNLdjeO2bkVaIvuiBwIYu7FRe7uEQlgIuP+pw4ExT5zA8deu4X0kAuln0yYsVLEYfuP363tjJtlMeeW9vRjz8GEseCwLEWHsYFDbHKSPKXnqgtVEPjdPf4uI/i2P4wkWgGyVpA4HyDEeR7TMZfAvvww9Oz3Kb27WRT2c+fF7v4fv0jX8H/8YpHXxIojT6cR3NhsWguAtUY6LhSIRbPpaLPj7wgWkDfJxmzaBSHftwjgvvojIv6oK34+PQ4qJRvH91BS+Z1/4mhosCr29+M5qnenr4nYji2fXLnxmrspNz3hhozDuuxqNQlPne1RfD7LfsQP3IL3SV/LUBauJeTdPlVKvE1GGPAL6umEYL9w65utEtI+IPmdkGVAp9QwRPUNE1NTUdGefOfwR5AWZJJejR0FiXV0gNm4DNzyMrBaXK/cNPvP4ShG99BIIdHgYi4VhQOqYmCD67d9GfndXF9Hf/Z3egLXbdUR+8ybRgw9qr5RgUHu/+3wgTrcb0k4qhcWhvx/n2bsXv3c6IeuEw5h3SwsWh0gEr9PTmNMDDyBjZmgI9+bLX56bZDNtfvb24vxtbTor5swZ7A0YBjZnzY2qBYJ8I2+bp4ZhfHyeEz1NRI8R0cPZSP3WOIeJ6DARsmLmO69g4chWMfrWWzozhAgRM3u+cLSci78Mj+/1En3rWzie/Wu4tN/n09a9fPwrr+D7ZBKbpkQgd5cLMkxpKYj7ySexGJmzZvgp4ac/BVlXVmLTctMmXBPLO21tiJ77+kDyTU14X1mJ3yUSRPfcg/MMD89Pvpn2LJqbsZA884ye34kTsBXm+ycpjYK1gKVmxRwioj8iogcMw5ia73jB8iNTrvXwMIgwldIGV52dINpc5QIel7NZQiFsjnLxEDd0Hh/Hxuvp09pbfN8++MV7PEg3vHQJYygFG90DB0CYJ05gnnx+s5zhciHVcft2fNfdraP/zk6QLhHGLSrSTpWRiNb+zW6M8yGXzU/phCRYq1iqxv73RFRERK8p1Le/ZxjG/7PkWQkWhHTS3b1b51afOIEinoEBvdm3eze0d5dr4T06g0Fo6DzWtm26QYXViii5pQWLCEevhw5h4fD5QK6JBKLhXbtAhl1d0NorKmaTonl+ZmnE48EiMjGhG3MUF2MOU1N4ImlvJ/rJT5DJU1GBax4bgxwzX0OLXNwyJaVRsFaxJGI3DGOBZqWCfINJN5kEifv9kCM+9SkQGxEiWO7Uk42k5iqAMkemvJmZSOiUw0QCv6+rQ2bL1BSI2uHAb594QncheuUVROmjoyBfbgzS3Y3fZCPFbER7//2IyoeGQNpNTRijqko3prZadY773r16XnMRey6bn5LSKFirEEuBNY75Kk67ukDqFy6AdOvqILW88gr05aoqvJ+LpDLls3/nO5Atrl1D1svu3ZBUuNtQTQ20bY8HUTjr1iy3cLNqJur0Tk/Hj+vGGsXFOj0yGylmI1oiom9/G4uJYWg3yO3b9abur/86FhGGeV5zIdvTzFxPSHPtUQgEKwUh9jWMXPxG2IOcI/LpaUS0Y2O6KxBr7dki1HStOB5HtofPh2KfsjJUmIZCiMhPnQK5ezy6QxJ3RTIjE1FzlNveDnIn0tky85Fipmtgj3YzKip0dgpnt8w3r1x9dsz/TW67DQvThx/ifnR0SEqjYG1AGm2sYWQqIGIdmlFTg+jT6cTGYyyGyL2sDPJELkUxPp9uCE2EBSEa1T1E2Z3xxAnklrvdOJ/djrxw3khlcn7vPaQGZjo3F+44HHoDd3gYC9Biskm6uhChP/wwTNEefhjv+R7lUijEZJ2piCvT+cz/TVpbQeYdHZCchNQFawESsa9h5LI5xy3bxsYQNW/ZAu3bZgMhZyLL9OhUqZlacSiECJRL5rk36ZUrkHXCYWj4HR2oKI3HQXDT00Svvw5SP36c6Fd+ZfY1mSWVcFgbiS3GhbKri+gHP8Bvt2/PbBeQi1a+kOwW2TAVFAKE2Ncwctmc47Zuzz4LAq6pgeZrtWYn9XR5Z2QEUez0NEj6xg2Q7tatGJNNtkpLQebNzSj+IdLWuH4/5JtYTDeYvngR+vfTT2fPdGGCzmR9mw3ma2Crg+PHsfnq8WS+R3ONuRCylg1TQSFAiH2FsBiv9FwbVOfS1o2RKTpNJiGf2O3agiAaReZKMAgdORzWm7G82UkEIh8exuat14vfc3l/dTUIP1sGymI9y83X0NEBUrdYsJnpcCx8A3MhZL3QpuECwWpAiH0FsFgCm0tGyLRQ5OLznR6d+v2o3rRYUJnJBUw7d0Kjn5jA+9tuQwTP3ullZZBcOKodGADJT08jO8Xlwt+JBM6ZCYst8DFfg8eDSP3KFVgM3HXXwjcwF0LW4gEjKAQIsa8AllKhmC0TZLHdedKj0+5uEHJNDcidbQciETS3+NrXZpp/2WzIAqmogG/KJz8JUh0ehixjsegMnevXsRhkkykWIoGYF7L338eiUVSE9Mv2dshP+/cvronFQslaepUK1jqE2FcA+dxwY58Wnw+phpwzTpTbQpEenQ4NgSDZQZEIGS+DgyDKdNLbsoXoc5/D54cP47ouX9Yt8vx+kG4qhTnu3589KydXCcS8kFmtOMeNG+iEFI3CR6ajA0VQi4WQtWA9QYh9BZCvDTcmOL8fpB6L6U1DblIxHzJFp83NIHIuyx8bQ2TOhJyN9Pi6yssxlz17sGEaCODz22+fvXFqRq4SiPmJ5/x5VJdWVuI8RUX4u7Z27mIiaSgt2EgQYl8B5GvDjQmurg5EyrKJuRApF6RnpbzwAshxYAAEb7Ui0ybX6L++HoRrtUJ6efzxzFk5mUh2/36i55+HrNPQAIfHuVIMQyHM1elEtM4NL8ymY+bzSUNpwUaEEPsKIF8bbkxw5qrNoiLIKQ0Ni8vMMM/N6UT0n2tUa/5tJIJN1cpKyDVzNds22xawt/qBA9q0zOzwSDTziYctDSYm9L10OmeajmXyuCES90XBxoEQ+wohHxouExxngnR3a3JbShSaq8NjJkkj1+vKRLKcLfORj+jP+FjzmOYnnrY2WABfvw5nSaLZpmO5NqgWCNYrhNgLCGaCq6qC/NLQkDupL1ZvzoekkYlkubjJjEzEm16tevAgcu8dDjyx7N4923SMSIqJBBsXQuwFhKVIOgslZ/Mi0NODBWQpkkYmknU4Zh+XjXjTnww8ntmt6zI1qJZiIsFGhBB7gSFbXvt8kfhC9Ob0ReD4cRBiaWlmP5ZckIlka2qgsY+PL5x4cyFtKSYSbFQIsRc4co3EF6I3py8CdXUgze5uTewLlTQykSznnS+GeHMlbclPF2xECLEXOHKNxBeiN6cvAu3t8JIZGoKOvVhJIxvJZiPe+Z5EhLQFgswQP/YCR7qXOhHep/uz5OJLzuBFgOHx6Dz54WHkzy93LvhCPNIFAsFMSMRe4Mg1El+I3pxJv7Zaib785ZWLkOd6EuFXqSYVCDJDiL3AsVBnwnzq19mQjzL+bHsCFy+iQlaqSQWC7BBiL3As1No3V/JLXwS8XvQPXUhP0KUQb7YnkbExVLZKNalAkB1C7OsA+bb2TcdCxspXGX+2J5Gqqsx7ClJNKhBoyObpOkUujbCXY6xcN3PnAz+JuFwzN2w7OmZu7BJJNalAkA6J2Ncp8umTslo9QbPtCUg1qUAwN4TYCwy56ub5JNh89QTNx6aqVJMKBPNDGYax4ifdt2+fcfLkyRU/b6HDrHWbSTOT1p3rsbmQ7ULOm21MooWNIRAIZkMpdcowjH3zHifEXjg4cmS28dX4OPTnTL0+5yPthS4US4m2Fzp3gUAwG7kSu0gxBYSF6ubz5a0vJINlqeX74o0uEKwchNgLCPn2F19JsjXP3e+HodjQED5Pb2knEAiWBkl3LCAsxO8lF6R7whAtX+ogz723F4ZiwSCR3Q6fd/GAEQjyCyH2AkK23O7FRrtLXSi4GvXwYbzORc489/5+okQCEtDddxO1ti4+v14gEGSGSDEFhnxa1a5kRyY+X1sbGldbTCGFaO0CQX4hxL7BsdiFYrHWAdKHVCBYfuRFilFK/aFSylBKefIxnmDtY7HWAfneJxAIBLOxZGJXSm0hok8S0Y2lT0dQKFjsxmu+9wkEAsFs5EOK+Vsi+iMieiEPYwkKBAvxgU+HtLQTCJYXS4rYlVKfJaJ+wzDO5mk+ggKBRN4CwdrFvBG7Uup1IqrP8NXXiehrBBlmXiilniGiZ4iImpqaFjBFwVqFOfJmy4GjR6VdnUCw2li0V4xSag8RvUFEU7c+aiSiASLabxjG0Fy/Fa+Y9YWFmoQJBILFYdm9YgzDOEdEtaYTXieifYZh+Bc7pqAwka+uSQKBID+QylPBkpGvrkkCgSA/yBuxG4bRItH6xsRKes4IBIL5IRG7YMmQoiOBYG1BiF2wZEjqo0CwtiBeMYK8QIqOBIK1A4nYBQKBYJ1BiF0gEAjWGYTYBQKBYJ1BiF0gEAjWGYTYBQKBYJ1BiF0gEAjWGYTYBQKBYJ1BiF0gEAjWGYTYBQKBYJ1BiF0gEAjWGYTYBQKBYJ1BiF0gEAjWGcQETDAvuJ+pzyf9TAWCQoBE7II5wf1Mp6aI6urw+sIL+FwgEKxNLLqZ9ZJOqpSPiPpW/MQz4SEi6fgEzHEvqiqILBaiVFJ/ZrGipcZocEVmt7KQ/11oyL3QWCv3otkwjHl7k60Ksa8FKKVO5tLteyNA7oWG3AsNuRcahXYvRIoRCASCdQYhdoFAIFhn2MjEfni1J7CGIPdCQ+6FhtwLjYK6FxtWYxcIBIL1io0csQsEAsG6hBA7ESml/lApZSilPKs9l9WCUuq/KaUuKaU+UEodUUpVrPacVhpKqUNKqctKqatKqT9Z7fmsFpRSW5RSbyqlLiilziulfm+157TaUEpZlVLvK6VeWu255IINT+xKqS1E9EkiurHac1llvEZEuw3D2EtEV4joT1d5PisKpZSViP4nET1CRDuJ6FeUUjtXd1arhgQR/aFhGDuJ6G4i+p0NfC8Yv0dEF1d7ErliwxM7Ef0tEf0REW3ozQbDMF41DCNx6+17RLTRTAP2E9FVwzCuGYYRJ6LvE9FnV3lOqwLDMAYNwzh96+8JAqE1rO6sVg9KqUYi+jQR/eNqzyVXbGhiV0p9loj6DcM4u9pzWWP4LSJ6ebUnscJoIKKbpvde2sBkxlBKtRDRR4no+OrOZFXxd4TgL7XaE8kV694ETCn1OhHVZ/jq60T0NYIMsyEw170wDOOFW8d8nfAo/r2VnJtg7UEpVUJEPyCi3zcMY3y157MaUEo9RkQjhmGcUko9uNrzyRXrntgNw/h4ps+VUnuIqJWIziqliCA9nFZK7TcMY2gFp7hiyHYvGEqpp4noMSJ62Nh4ebD9RLTF9L7x1mcbEkopO4HUv2cYxg9Xez6riINE9Bml1KNE5CSiMqXUdw3D+LVVnteckDz2W1BKXSeifYZhrAWjnxWHUuoQEf0NET1gGIZvteez0lBK2Qibxg8TCL2LiH7VMIzzqzqxVYBCpPMdIho1DOP3V3s+awW3IvavGobx2GrPZT5saI1dMAN/T0SlRPSaUuqMUuofVntCK4lbG8e/S0SvEDYL/30jkvotHCSiXyeij93638KZWxGroEAgEbtAIBCsM0jELhAIBOsMQuwCgUCwziDELhAIBOsMQuwCgUCwziDELhAIBOsM675ASSCYD0qpJBGdIyI7oer2X4jobw3DKJgScoHADCF2gYAoYhjG7URESqlaIvpXIiojoj9f1VkJBIuE5LELNjyUUpOGYZSY3m8lVJ56NqC1gmAdQDR2gSANhmFcIyIrEdWu9lwEgsVAiF0gEAjWGYTYBYI03JJikkQ0stpzEQgWAyF2gcAEpVQNEf0DEf296OuCQoVsngo2PDKkO/5/RPQ3ku4oKFQIsQsEAsE6g0gxAoFAsM4gxC4QCATrDELsAoFAsM4gxC4QCATrDELsAoFAsM4gxC4QCATrDELsAoFAsM4gxC4QCATrDP8/Qp8wYUYlfs0AAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"X.plot(x='D', y='Y', style='bo', alpha=0.3, xlim=(-5,5), ylim=(-5, 5))"
]
},
{
"cell_type": "code",
"execution_count": 32,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"<matplotlib.axes._subplots.AxesSubplot at 0x7f114fad8ac8>"
]
},
"execution_count": 32,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAEKCAYAAAAGvn7fAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJztnXtwnFd9939ntSutVndpdbElW3Ic2Y5v5CIlkBASEwq2CUnToaWltAkdyD9vmbZDStswLZ3OFF7mZUo7tB3G4W2Bl3RoubhckjgkISSBEEe2c3F8j2LLXt21ut+l1Xn/+PrXs5J1WWlXWunR9zPjWe3u85zn7BP4nt/zO7+LsdYKIYQQ7+BL9wQIIYSkFgo7IYR4DAo7IYR4DAo7IYR4DAo7IYR4DAo7IYR4DAo7IYR4DAo7IYR4DAo7IYR4DH86LhoOh21NTU06Lk0IIWuW48ePd1lrSxc6Li3CXlNTI8eOHUvHpQkhZM1ijGlK5Di6YgghxGNQ2AkhxGNQ2AkhxGOkxcdOCCHpYmJiQiKRiIyOjqZ7KnMSDAalqqpKAoHAks6nsBNC1hWRSETy8vKkpqZGjDHpns41WGslGo1KJBKRLVu2LGkMumIIIeuK0dFRKSkpWZWiLiJijJGSkpKknigo7ISQdcdqFXUl2fnRFUPIOiUSEWloEOnsFCktFamvF6mqSvesSCqgxU7IOiQSEfnRj0SGh0XKy/H6ox/hc7K8WGvlve99rzz11FP/89n3vvc92b9/f8quQYudkHVIQ4NIYaFIfj7e62tDA632maT6ycYYI1//+tflt3/7t2Xfvn0yOTkpjz76qBw5ciRlc6awE7IO6eyEpR5Pbq5Ie3t65rNa0SebwkLcr8FBvL///uTEfffu3fKRj3xEvvzlL8vQ0JD84R/+oWzdujVl86awE7IOKS2FSKmlLoL3pQuWl1pfLOeTzRe+8AW5+eabJTMzM+W1syjshKxD6utheYrAUh8cFOntFbnrrvTOa7WxnE82OTk58rGPfUxyc3MlKysr+QHj4OYpIeuQqiq4E0IhiFQolLx7wYvok008qXyy8fl84vOlXoZpsROyTqmqopAvxFp9sqHFTgghc7BWn2xSZrEbYzJE5JiINFtr703VuIQQkk6W88nmb//2b5dl3FRa7H8iImdSOB4hhJAlkBJhN8ZUiciHReQbqRiPEELI0kmVxf6PIvI5EZlK0XiEELJsWGvTPYV5SXZ+SfvYjTH3ikiHtfa4MebueY57WEQeFhHZvHlzspclhFyFxbwWRzAYlGg0umpL92o99mAwuOQxTLIrgzHmSyLyByIyKSJBEckXkR9aaz8x1zl1dXU21ZlWhKxH4lPe48PxVmPkxmpZgNZyByVjzHFrbd1C5yct7DMuereIPLJQVAyFnZDUcPgwKjPGlwbo70dY3gMPLH3cVIvwWlqAVjOJCjvj2AlZw3R2Qijjyc3F50tlOUr6xtdc8fnwWliIz0nqSamwW2t/wRh2QlaO5Uh5Xw4RXo4FiMwNSwoQsoZZjpT3RAtfLcZdw2qSKwtdMYSsYZYj5T2Rp4DFumvq67Hg9PeLTE3htbcXn5PUQ4udkDVOKlPeIxGRaFTkuecg2Hv2iASD1z4FLLZOuS5ADQ1YgEpLMR43TpcHCjsha5xkI1j0/HPnRN55R2T3bpEPfEDk5EkI/L591z4FLKVOOatJrhwUdkLWMMm2bos/v7dXxO8XOX1a5LbbRO65x4VOzhxrOX3mqyXefS1DHzshq5hIBLHqhw7hdaYPO9kIlvjzBwZEiotFcnJELlzA93NFriyXz3w5Qi3XIxR2QlYpiYhcsmGEnZ0io6Miv/413DBvvikyPi7S14fv57LCl6tOOePdUwNdMYSsUhLZoJzNJdLUJNLcDCt/IVeGMSIvvSRSUiJSWyvy1lsib7yBTVO1wucKnVwOn/ly9hhdT9BiJ2SVkog1PtMlcvGiyK9+JVJZmZgrwxgRrSqSlydy/fUiGRnOt77SKf/L3WN0vUBhJ2SVMlPkurpEnn9e5MQJ52+f6RJpbha54w6RLVsSc2VMTYm8730iWVlYIMJhkU98QuTOO1FrZqU3LRnvnhroiiFklRKfVdrSIvLkkyJjYyI33SRy5Qo+U4taBfjQocVljZ44gTj1d70Ln1+4IPLyyxB4XThWEsa7pwYKOyGrFBW5I0dE/vM/RSYnRYqKIOp9fSK33nptQlAiYYjxIY719SIvvohFIzNTpKAAIY+VlYsLm0wljHdPHrpiCFnFVFXBBx4IiFRXQ6B9PpHWViQUzYx+qa+Hn/2550SeegqvFy9Od2XEb8qWlYncfbfIyIjIpUuwkkdHRdraRGIxRqOsVWixE7LKOXECrhFj8C8rC2J/5gwyRGcysynQzPczI0/CYcSv9/WJ1NTANTM6KnLqFAR/LphItHqhsBOyyrFWpKICLhgRuEzGx0UmJq7dVGxogDjv3es+6+9fOESyuxsbsKEQ3odCEPWmJmzUzhTvRDJeKfzpg64YQlY5dXUQ8s2b4f/u6YFF/ZGPXCuUSwmR7O9HiGM4jPBIa/Ha1wdxni1BaqFEImaQphda7ITMYKmW5nJZqPv3w/fd2Qlhr6rC+Pv2XWtNJ7J5Olvkyb33wkJva8PCUVAAl0xp6ewJUgslEi22+iNJLRR2QuJYalGtSETkW9+C4I2Pw11y+rTIgw8mL2RVVRgnftGoqhJ59dVr51ldLfLCC4igKS8X2bABNWDKyq7NRI2fl/7uXbtcw45Ll0R27EC5gb4+iP3WrSJDQwsvIMwgTS8UdkLiOHIEsdwTExCy2lrnYphPoI8cETl/Hu6MoiII27PPipw9K/Jbv5W89T5TiA8fvtYijkZFfvxjhCq+/jqiZnw+iHVNjRNsXahEpi8Wt94KgW9vx3kZGSI//CG+q6lBDP1LLyEBaqHOTfHC39WFe9renr74+PUGfeyEXCUSQXigzwfRHBsTOXoU/uyFimodO4Z6K6EQBO3yZbgyOjqWx788my+9tRWWdWuryHXXwVUTCkHgx8edLzwWE/niF0X+/M9h3WdkYI6vvgrBPngQvz0zUyQ7G9b/+fOw/K1FlM1CRcDUj3/xosgrr7iSwBofT1/78kKLnZCrNDTAfWAMRFAjRE6eXLiHaHzNlZYWiLq1EMhk/cuz+e5nc4W0t0PAc3Lc3DMy8O/CBVjLXV0IY7xwQeTGGzHvhgbUX4/f/CwsxG/YvRsLRW8vnggOHsTCIDJ/IpEK/9e+hqcfLTIWDl8bpUNSD4WdkKt0dqKqoYpbdjaiRjo6Fq5VcvPNKL7l88EHHQhAeG+4Ad8v5F+OF2+NV5+awnhtbaj9Eu9Lv/VWWNg69uAgLGK1shW/H8KuZXgvXMD7rCwc5/O5z2+7zc2xvByuqLEx+NmnpiDuwaBbNBaiqgo++dtvd9dJ5F6Q5KErhpCrlJZCuG67DcLX0wOre9++ha3LAwdQGTEWgwiOjWHj8pZb8P18FQrjQwMzMrBAvPQS/n7rLZG3357uSiksxDkzXSGf/jTEuLsbcxgexvHBIBaaqSkcG4tBcEdHcf3sbAi/zlGfBmprsUgNDyNiJhBYfEEuVmtMD7TYCbmKbggWFkLcdUPwwIGFz62qEnnoIVjdGzaINDbCjVFcvHBd8/jQwFOn4LYQwRgTE9iMVVeKiLN453KFPPYY3Cfl5bDsNSpGNy8rK1Gi9+hRHD815URb5xhfS+bkSZx7zz0IvVyMC2WhTVayPBirjsEVpK6uzh47dmzFr0vIQqQqFn0u18rM7M2GBpEf/ADvt2/HJmxREZ4UenthgY+O4glg/36MrbXSH3hg8b8hPpxzdBSi3dGBp5IDB5Yna5QZqKnDGHPcWlu30HG02IlnWYqgpKqyoI7xH/+BEMSsLPiqt2/H5qr6yAsLRTZtgivk6FH4xLU+i4Zb/uIXEPuOjunW81xhg3P9Br0fGrVTVATLebmFltUaVx5a7MSTxFum8S6AlSpDG4mIfPObCJ/MysK/gQG4Qerr0RCjshIbo83NeC0rg4tE0/p378Z3Fy/CTdPfj03UPXvgN1/M71nM/Zjv2NZWke9/383/ox9lE4yVhBY7WdfMl9Kur8vpGmhoQGih349rGwNL/OxZWN7t7fg8KwuWeUkJPh8cREhhT4/Ia69B7O+7D5b66Chiyk+cwDkVFdPDBhsapovunXfCR9/ZCX99ZWViKf5PPQWffk8P5pOTg72Cf/1XjFNSgro1vb0iX/mKyCOPUNxXGxR24knmSmk/exaukMWWDJiLudw9WlqgoACvsRiEPhZDpEk0Cgtcs1SNgTVeW4volsOHpwtxTw/OGRuDxT6zrG5Dg8jf/R2E3BhUgnziCZHf+R10Rzp6FEKcl3ftJuzM3/P88/Dht7UhEmdwEAvQ00/j9+nmrr5+//sU9tUGwx2JJ5krzK67e/6qhIthvgqGpaWIKS8uhgirgKpQFhXh895eWMRTU6gto/HeMzNLVVxjMZc8lZGB3yMi8m//Bj99MIiF5Nw5+NH//d8RLlle7hKV4secGXbY0ICnhI4OhEHm5+O1vR1zHBqafnxhIZ4QyOqCwk48yVxNkYuKFi5rmyjqsnj6aZFvf1vkZz/D+6eewvXDYYhsVRUs695eWNw+3/SM1LExiHthIeYqMn1h6uqCgL/zDgS3rw+LSCyG3yOCRSEzE8031Dq3FucdPgxRj8VghXd0wPf/xBMYOz69X5O0+vpwvv7r7cUTRVfX9HvQ24snC7K6oLATTzJXLZPt21OTMKMui8FBVw+mrQ3vn38exzz0EApmBQKw3CsqRN77XtRxsRbulakp/D0ygteWFlRh7OpCdUWttZKZiYWioAAFvqJR+NsvXYJwa1OM1lYsHH4/FoyREczvmWdENm6E1f/ss5jfPffgvsTXbtEkrT17MJ+hIbzu3SvynvdgcYhG3Ws0ig1Usrqgj514lrnC7BaTMDNXPHpjI1wUHR0QwqwsCGlHB4S7oQFx5p/+NMZ57DFklObl4fi33oIwW4vXsTFX7ld9/9r+bnISvvfeXnzX3Ix497w8ZLZeuYJ56UKRleX8+sEgvuvuhlV/yy3TffeKbqJqQtG2bfDXq6W/bRv+/su/RFbs5csY55OfpH99NUJhJ+uK2ZpM3HXX3HHfGvaXkQFBs1bkfe+DRT02hjEqKqa7LPbsuda1Yy2iVBob4VuvqHDf1dbC9TE2hsWiuxsunbffRlXFvXvxud8vcvw4Ps/ORlKR1n0vK4O1PzyMsEpjIOoZGXiSKC9H+d6LF/HUEk/8Jmr8/RkZwVyKihBrrxvD996b2v8mJPVQ2Mm6I9GEmflS/cvLIeIlJc5lkZGBCJTZCmUZg5BFrfNeXQ1B7+52FvrevRDlo0dheff1YdyjR2G1h0IQYWuxOLzzDurT5ORg3EAAAtzRAYtdSwXk5qK9XnMzhD+RDktMKFrbJC3sxphNIvJtESkXESsih6y1/5TsuISkm/iQyb6+6an+dXXwfefm4nN1WajLJN61E4lAbHt68O/IEbzPzUU44o03wqp+800I8t69OE5rovf14VVdM7rR2tnprjkwgLlmZ+P1zTfdRu3tt6PoVzSKa/b2Yl5Lrd3CEgGrn1RY7JMi8llr7QljTJ6IHDfGPGOtPZ2CsQkRkeUVk7nGjq95XlAwPdU/HIZro7kZwt7Tgw3SUAgC/OSTbqyGBnQg6uhABI3fj2PHxtCh6Nw5V6RragqboOpSqaqCOPf1wf2jcepTU1gEsrLgrjEGnZpOnYKoZ2ZiLpWV8I9rFM2WLYm7oua6V0tpHUhWlqSF3VrbKiKtV/8eMMacEZFKEaGwk5SwnGIy39jxlQm3boWl3dUFd8dzz0EUP/OZ6Zmfjz0GAS0thbukpQUCfcMNcKmUl2OhsBZjdXXh8/e/Hz7z06edgOflYRFpb4dwZ2bCNz866kR9aAjv77vPxbnrfIaG8O/sWVzP78cCJDJ7AbFEFk82qV4bpDTc0RhTIyI3icjRVI5L1jfxYrLYpKJIBOGAhw7hdWZLtvnGjg+Z7OyEv/q661wcfHyZpUgEou73o2zvxAREOhZzqfnd3XCjiEBwx8cxVjAIi76zE26THTvgg9eN1IkJiLox+ExDEvPyUE/moYdEPv5xuIYuXcJx6vsPBPCkUF6O3zVXa7r5kq3ima0l31LzAMjykbLNU2NMroj8QET+1FrbP8v3D4vIwyIimzdvTtVlyTpgqR3vE7H0FxpbNxIPH4alHr/pGN/iraEBG5ylpRBW3TxtaYHbpbcXFnVzMxaQ4WG4c4aGpreyi0Yhvh0dGKezE/7766/H/Ftb4U7Zvh2uodpaVyM9IwO/dXLSxb2fOQOLvrsb+wITE9h0/fu/F/ngB51VnqglPltLPjbOWH2kRNiNMQGBqD9urf3hbMdYaw+JyCERVHdMxXXJ+mCpYpKIWCUydiQC37g2udbenfELwLlzEO8rV3DMhg34/PXX4eOuqYEvvrUVQq5uEWPgphkehvA3NuK8fftgtTc1ifzkJ7C8y8qwCFmL6/X3Y6F49llY+Z2d+F2hEPzsWm4gFIKg//KXSJ6qrUWyklrl99+f+OLJxhlrg1RExRgR+b8icsZa+w/JT4mQ6SxVTBIRq4XGVr95Y6NLUHr2Wbg68vMh2A0NsIJzc2GBDw6KvPEG/s7Lg6V88iQs7Ztugk+9qQlCvns3xjh92kXBHDyIY0Xw+t73InN0714sGufOwW++ezc2UQcGEGNfVITEp9xciPb4OKz47Gz8PTyM79raXGmA4mKUQDAGoq/hmLW1mMvMxXMxeQAkfaTCx36HiPyBiLzfGPP61X8HUzAuISIyd3mAhcQkkX6b840d7zffvBnW+MWLsJQjEVjLOTk4prIS51ZX47POTrhA7r0XlvbEBBaDnByRz38epW43b4agNzXh/M2bIabqh1dGRmBhFxYiOqa5Gdfp7YVo6yKjG6ft7VhUYjEItzF42hgfx1PClSuYT1cXvvvpT12Z3qwsXO+FF/BbmVW6NklFVMwvRcSkYC6EzMlSkmYStfTnGjveb97ZCcu6sRHCl5MjsnOnq5E+MoI+qRcuwG0yOAhB37YNY2m4pIYtapXFkhKMceECxh8agj9/3z742/v6YJ3v2IFaLSIIaQyHYblv2gRLvLcX4/7Gbzi/fCgEYW9thVU/MYHzs7NxnZMnsVjEYhD9u+/GPPr68HlFxfxNOBjuuHph5inxLMm6DdSVMzoK8SwogChOTuL77m64TzZuxGbnHXe4WufPPefG6eqCla3FvBob3VPDrl24zpUruMbOnXCnfO97cLXk58Pd0tGBccJhLAa9vVis+vpcolJxMfzvlZWuWFhTkysGVlyMOYdCmMf4OJp57NyJv8NhN/+pqdk3pxnuuDagsBPPonHZ589D0HQzUWT+nqAax20MNkFffdV1N2pvh7tiwwaXRLR1Kyzk/n73ZFBaCpG9eBE+9StXIMihEKzi4WFY9FeuYIz8fJx34gSEOBRCTZjbbxf5yEfgs3/jDVjye/eij6o+QUxMQIh37cLr9u2w5nfuxBhtbdNLFo+O4l9pKa6j5YSPHJnfvy7ingZOncKiUlCA3z+zTjtJLxR24knUZRCLYWMzIwMCFgq5ZtKRiBPxqirXXFpdDB0dTvQKCmB1Z2S4xCCtxBiJoIpjJOKeDB58EC6Qv/5rRKgEAghZ3LULc+jrgyjrJmcwiIUhFnPXGxrCYlBbi/f65LF3L6o0vvQSwhmHhyHor7yC+Vy6hON/8zfd/fjv/8bvC4fxm8rKcM2hIYy5fft0/3ppKT47dGh6spLPJ/LiixinsBALxIsvosAZWT1Q2MmaYTFlBdRlcOoUrNRQyNVMr6jAhudddzkR/+pXXeEtayGmtbWuf2hZGY4bGIBL5PRpJCtdf72rfFhf7+b4ne9A0NvbIYLBIOZ99CjE8fx5iP2OHRDiSMSV/43F8BtGRmCZP/QQvvvQh1yJgpYWN8fWVgis3485trRgI1ZdNyKw/M+fR2il349r9vRgDrfeisgb9a9nZOB333ijewJRP7q1eJIRmf5qGcC8qqCwkzXBYjft1GVw4gSENCcH7pPRUWctBwKukfTZs7Cct27F+adOQfiGhiCs3d3YbBwfh5+9vx/+6/FxiOjMOfb2QmC7uyHUoRCu19bmKjWOjMBlUlODv7WWeiCAYwoLcf4bb0DAt22bfg+ef951VMrNxby0YYfPB6FWYQ8GEaHT2Igngp07sTfw619D9OP967/6FcadzY9uLZ46XnkFm7slJSLvfjeFfbVBYSdrgsVu2qnLQFPx+/shzFqdUeuTj42hpZ1uiAaDeNWsz1AIi8iFCzh3YMC5TjIzIe4/+AHO6e52TSxaWyGOeXkQXs3+NMZtuNbWQugHBnBcZyf85WVlLrRRa8HoJnAs5vzbjY2uW5IIrmmtE3gtKqbRQL//++6e6VNPKOR+szJf/L8xuH5NDaz/kRG8v+OOpP7zkhRDYSdrgsWWFVCXQUUFXBDNzRBQvx8W8cgIXCV790L8tMmFJglpu7qbb4Yoamy41kIfG8NYGRkQzsuX4fOORPD51JR7UvD58GQwMOD6lFZXQ7xDISwOPh/EfetWCO74uGtJ96EPYfF6/HGXCFVY6Er6Dg7iX0YG/um8MjOnRwOJXOvKEsFTQDSKJ4zOTvyWmVU/dEM4GnXWefyrYcDzqoI9T8maIJFko3i0Y5G6F6yFoPp82BTMzoYotrTgbxG4JwIBJ+AbNzqLt7sbIjkwALHW2uxjYxD6X/wCAq3+5oGB6c20u7rwOjUFq3583NVLr66Gr/3DH8Y8olEcV1QEP74KcHe3i0/3+ZCV6vPhOhqjPjLiQhtvvlnk4YddJceZRb6++U1EwkQiKFtw8SKeFu68E3sIFy9ObwReX4/373sfXEZa/+Z973NNuMnqgBY7WRMslGw0W6hiMIiknr4+5+Pu7IRlHArBiu/qgpWsZXHVf//WW67v6OCgSzry+VxTjVgM37e34zNrIZjq1tFFQQTHZmfj/GgUoh0IQNyvvx7z3L8f6f0nTmCsujpX4EsE5/T24hraKk/9+fn5EHStU1NT4yzqSETka1/Dby0vhwtIBKGQXV24T7W1mO+2bS7ipbkZ38XH/5eW4vqaLCUC4Z/ZMYqkFwo7WVXMjD0vKoKFXV8/d7JRJCLyrW9BtKNRbHpqCv+73w0h15oqN9zgKixWVkLUw2Fc4+RJiGpvL0QxGITlWlCAGi85Oc4KF8GYfr9zz/h8WCBU2EVcQ+iSElw/GHQbnjt2QFCHh50Ia/NrvQ9PPulq1Fy65Bpfj47iM01M2rzZ1ZdpbkZyUm4uxnj1Vdybigqce/Qo5l1UhCeHiQmMMTLiNlyrqzHXhx+e/t+HRcDWBhR2smqYLfY8EoHYfO97SM45cODazdIjR7AQZGVBNLUuytQULO/iYnzX1wdx1I5EZWWwbicnIfLbt8MlcewYBLC3F2J38aKLmpmcdNa576ojMzMT4jg56XzrIs73rOPk5DiXRW8vsj4rKlAGIBjE4lRWhgWtsRGZpzk5ron2nj1YfAYHsfF76hSs/p07sQl77JjbQ8jJwW967DGMo6KulvWpUzivoADvR0YwB22bN5ebi0XA1gYUdrJqmBl7PjmJSI/cXFjWp05BXGeGOB47Bov4yhW4KLKyILa9vcjaHBmB1f7OO7Bmx8dh3f7Zn2HhGB52UTZvv43vR0dhmatgq+tFLfBQyLWnGxvDuVNT08P+jMHcNelIqy5mZWFMnw+/Txei8+dhWWtUyxNPYA6TkxDdwUFEn1y+jJIFO3bg9124gLFPnMDTyqZNqLW+ZQtCJVtbIcBPP4256EZuNIr6MCKw4kdG8J361OeywtnsevVDYSerBo180SJU589DkCcnYYH29EzvbqTohuXwsAtT1CiR0VEI/86dGHtyEt/FYiLf/z6s+R073Fj9/S6zVIXa58Px2sFIQxaNwbF+P4Q+GnWp9T4f/k1MQLS1k5JGrliLf01NIsePY2EpKcG4V67AYh8agvWuTxUiuCf19fCN79uHaxQXQ9xzcuBG+djH3KZxaSkWtJ4euFe6u3F/fD5Y8ZmZWBS0Ro1G6tAKX9tQ2MmqQSNffD64HN55B0JaVoa/33kH4pOdDdHRaJGbb0ZSzfg4XBIarlhSglj2piYk4mgp24oKCOnp0xD+EydwDd0oHR7GuHl5EGvNBFWfdjgMcdb4+OxsuD8CAfytVr76xtX9MjXlNlwDAfzd1SXyX/+F375nD8ZubHQLijGYU14ejs/Jwb3Jy7u2frrOW0VdBJE9x4/jXoTDOC8chpBnZ7tyxZs2oRn2XGK+nM3ESephuCNZNdTXwwfe2uqiSdra4It+5hkXD26tyFe+4nqeHjgAYdakoJ4eCNaWLRDmM2cgtJrI09zsKipWVOC8116D0JeVQbSnplyLO7/fuXf0KcAYiKTGiweDeCLIyHDW9UzL3hh8Hwg4S1799u3tWJza2mBB+3xwjWid9e5uPMm88YbIyy/jHsTXT3/ySdy7ixfhpunowKKTkYG6Mvn57vjbboP1Ho0m9t8l0X6oZPVAYSerhqoqCGtVFdwLU1NOBDMzYXH398N1UlICV4qSnw/RueEG18g5MxOiqX73sTEIsM8Hcc/NRTSJujMyMiDWt9ziwh7VtVJc7LI0s7Iwp1gMdVY2bsQ41dXONaPXMQZWtjbPiMWcz94Yl0wUCGDB+OUvnU89I8ONMTDgFru8PBy7dy+Oi0ZdXZj77sMxzz0HAb7/fgj5nj0InXzPe3BcUxOegBIR62SaiZP0QFcMWVVYC9/x0aMQyslJkR/+0Fm9OTnOir18Gec0NMA6Hx52kR+6STkwgLC+WAyW7dCQ29CsqIBABoMQPmNg1RYXw3ovLHQi2t+PBSccduGNhYW43uXLLrZdBNe31lnv6ppR/7xG7GRn49oaySICy1wja66/HqI9NoZrXXcdrltRgd8UjUKof/1rLGbj41gY77kH8x0ZQVz8iy/Cmt+2DXVhgkG4tLTeu8j8JRqW2kycpA8KO1lVaO/N11+HiG7cCCGbnHQdhkQQtVHRadj7AAAcCElEQVRZib9VeGprsSCIQNTb2tzGqrUYR0MPs7IgwENDOG90FOcVFEDUq6pwjLpAtLjW7bfj87ffhvvG73euE201J+I2ctUto5upPp+rzWItxHhkxEXCaI9Sa3ENjYYpL8eitnEjBN0YzEsEr1lZLnRRBL/nJz/BAqDVJTUh6d57EWVUXT393s8l1kttJk7SB4WdJMVSNtXmOicScT7yvDyIx6lTsEa1UUZBAYQtGhX55CcxngpPOAyRfvll+M3Ly5FY1NAAq1lbxE1N4bhgEJuIOTmwao1BOv0LLzgLPBBwvurubsSzNzW5Ddxg0JX6bWtzJQTUeld8PljeU1O4Xm4ujrXWNcowBr9B29lp3fQtW/C71UWTm+uyV9VdFR+6KIIN1ljMNffIzXVPOlrJMVGxZlLS2oPCTpbMUvpfNjQgaWZy0rWda2lxSS81NS6S4+RJCNnWrXA5PPssrjk4CHfDhg0YU4WnsRE+ao0Zv/FGCPHOnYjhFsHYRUW4/u23Y2NVM1GnpvC00NMDEdes0txc5+dvaUFy0PHjuL5u6Pb3u7BGLcQVn4GaleV6jcZiENTsbNeCbmoKgl5RgeMnJlx2aGkpNnrLy7E5qoXCyspw7q5dWBAzM10Zg44OnKt1cERwL3t6sKAePOiydcfHXcekBx+89r8Zk5LWHhR2smQWW0o3EoGo+/0Qh9FRuD127nQWfHk5xPJDH8Im5vnzELUbb4QQj41BiFpaIEwPPohr3XqryKOP4rvSUmf1hkJw67zrXXgKEMEYra1ornHLLbje229DxKuqIKqDgxBG7XdaWCjy/vdDWPPz8dmZMxD14mLM3RhXRtfnc6UCNA5eFwpNdJqYwLhaaMvvh3Cq60cTmfTpoLsbC8N9912bgatPQSq8+/bhaWdkxGWbjo5Ob3k3s4b6fDXVlyMpiSGUyweFnSyZxW6qaT3xsjIXSigCkdViU01NcGloP80NGyDa0Sg+D4chuJ2dCIE8exbVBU+cgOht3Ih/BQUQ1pER+JUrK51Pu7cXLpmJCXz2wgv4vLwcdVU0aaevz1Vk7O5G3fahIQhmXh4WiMFBjN/dDZdIVpbbOFX/vsazayneWAzjjI1hjqOjrmF1vAsnKwvXUes9P1/k937P1ZOJZ6bwapu+t992bp5oFBuo2oVpyxYseEp//8o1pV7K0x5JHAo7WTKL3VRTyyzeigwGIey33or/Q3/3u/ALa9XCxkaRRx5BZExvL0RUBIIYDGIh+NnPEN+tdda7uhAKmJ+Pc3bscG6SUAg+bmshmq+/jvHy8iC20Si+10SiqSnX/FnfT0xA+EXw5DAwgGvfcAMEq6sLoh2LuboxupGqoY1+P+aq3ZIUTVzSImTqTikrw3EHDiT236aqCp2f4qtF3nmnqxb55JMQ1K4u1xIvLw/3fSVY7NMeWRwUdrJkFrupFu9+EYFg9fTgs2gUfvDcXAjsxYsuSuRLX8I55eWuWbOG/bW3Y4yRERdbPj4OS7WkBFb+tm2IltEiXNrSLhyGKGuDjaEhjNHbC1HVmPN4F0VmJsaZmICVXlIC0e3rg2tGS/Oqrz0ry52nrhWdl2arqrvG78e/jRtdSYNQCAtfcTEsdRW9RNwYVVWzW/f636KpCfdV68vrBnUksvziyhDK5YUJSmTJ6KaapqWHQvM/StfXQ+y0oYVGfmza5ATRGFjpRUWu/6f6stvacLz6q5ua4D4IBCDy4+OuXvrZs1gc7rwTJQfKy7F4dHfjWjU1EPNAwDWrmJpygq49SLWBhQgWrv5+LCQaYqhlBfr7ca52YvL7XYXFd7/bZXtqWYJg0DX70AUkPmtVo2j8fjSa/pu/cSUUEs0EjUREDh8WOXQIr/Hf19cjll3DL/X37t69MolHi22cQhYHLXaSFIvZVIuPrggGIXbRKARufNz5hH0+NHnetAnnZWbCyvX7Iaoalz41BcEPBFwf0pERCGVmJsrwvvmmyHe+4zoS7dgBN81PfuLEORrFa2kpxC8318WXxxO/8Tk56coUaE1zbbFXXAwB1wzN/fsR837wILJl29rwxDEy4jZcjXGx9P39mH9ZGd5/9KMuHLShAU82wSD845oJKjLdjbGQD7uqCtE1vb34V1AAUS8uXhmrmSGUywuFnaSMRN0D8Z8dOgRBa2iAeJ85A5FtaXFp/pWVEM5bbnFVDM+ehZDm5uJ6o6MQQu3zOT6Ozy9ccLHlfX1wPfz8584to5Z1dbWLXdc+pLrZObPtm7pOJibcE4O21NNSvpoJOjKCevF79mCz8v77IWjt7ThfSxBo7Lvf77owXX89nm4iEWwif/ObeAp56SUsID/9KeL0y8tdJUq954n4sLdtm16yWAS/ZyWsZoZQLi8UdpISlhrlUFqKqJSpKVfWtq/PuUduvx3i2NQEsbvuOox9ww2ur+jQEK6bkQGxzMnBuKdP42lAXStqWY+OYkHQjczcXAjl3XdDkF9/HXVUZiYZKRrbPjjo6r7oU4TWpxkZcU8aTU0iH/84zt22DYlVX/oSfrNGz4i4p43cXPzO/HzXQKSrCxvELS2uW9PQEBKrdAHIyXH3XH3Y822OpttqZl335YPCTlLCUqMc6uthiXZ1uUqHo6MQRGMgxMZA4LUqo7a909rtSl4eBC8vDz53reo4NeVqtmgZXY3M6e+Hv103Y6urUb99phtGQxa1nrpuhMY34VBL3lp8Xl2NuWzc6MoIKFpjPhSCGyQnBz7v8XG4WAoK8PeVK7iXTU2Yo/52zVqNxWB1X74s8ru/i8WkoSGxzVFazd6Fwk5SwlKjHKqqIFznzrnwvtJSxHTn5EAYi4pcE4qf/9x1L9KGGKOjsN4rK7EARKOwPKuqIG4ZGbiWirMSi+HYyUlkbx47hszTZ5915yhqkWvp34kJWPvG4FXDInUD9rrrkC3b1YUN3N5eXOvcORTt0jIEfj+OUTdOYeH0bFG/H5Z0Tg4EXGuyBwL4zRonX1GBKB/NZD14EKGOWmly5uaoinciVjMTidYeFHaSEpIpFGUtrNr8fFic4+Mu7O/LX3bNql9+GUJvDKz37m4IcEGBC4U8f94Vytq3T+Txx6eHFWrCkVq7avF2dmKjVdvUacGumf51jUsPhTBmZyfe6zWysjAndats2wY3TGsrsm5PncKisGGDy0wdGYEob9sGYVd/fSiEDc6MDJG6OiwKGqnj87kaMFVVriCa3vO5NkenprD5mqhIM5FobUJhJykhEX/tXJZfXh42PpuacF5urnNjiLhzAgEImMa4t7S4+ug33ghL1e+HCD3yCK536hTcEX6/SxpS61ut7YwM1zxba7koauFrNmlennPD9PfjWqOjeJ+Xh7lrQ4tdu1zq/1NPQYh7e13pYW28rfH0paU4r6gI19OErAMHEFlz4QI2TqNRjFVSgrlmZUHEZ/Yqjd8c7epCfZsTJ7ApXV/vwiQXqu3DRKK1B4WdpISF/LXzWX41NQgHLCrC8YODSEK66Sacq4WqtNiWWtPagWjDBvzLzMS5Kly6cDz6KMaPT+/Xv7OzIcrj4y7rM5745KRYzIm5Wvz9/dNdUFlZsIx37cKYra2IinnsMXfc5CTEPxDAYqaWfn09nhq0/Z829hgfx3mf+QyqUr7wAu5JOIzQzeJiF80Tf891sY1GscBduoTrbNyI/0633TZ7D9l4mEi0NqGwk5Qxn792NssvGhX54hchZB0dsCBLSmD5VlXhbxEsEpmZELDGRnymm6Aa793Whhj1kRG4aESwmDz9tOsZmpkJ37ZugM5MCNIolbmiYUQgsmNjuGZJCYTbGMTcqz98bAzCfeIE2t2VlblG2iL4nYEA5qzWe0WF+/3BIJ5YtJepbog+8IDIpz6FfwuhT0cDA8i6zczEPPfscXV0LlyApf/SS3O7ZliLfW2SEmE3xuwXkX8SkQwR+Ya19n+nYlziHTo7IZqnTkFcjYErpbUV4rZ1K4RuaAgW74034pzDh+E37+hwvnfd8MzOhmD39Lj471gMlr8ILOVXXoF7whi4ccrLIfyhEASqu9sJuW6IzibsKvzGQCTz8zHexo2uf6n6vnfsgPCdPOlcSDU1sMa1QcbkJH5rZSXmu3MnfsOlS9h43b/fXVs3RBPdxIx/OrrhBty/QADnaKemYBDHNTfj+nP5z9MdEkmWRtLCbozJEJF/EZHfEJGIiDQYY35srT2d7NjEO/h8iLkOhyE4J0/CYqyqwvuJCYj71JQryPXmmxD78XGI+PAwPs/Lc6GO2u5Oo2S0afOXv4wEHmsxfkeHE0RN5c/Kwj/1kWsW6Uw03n18HMdPTWGs3FzMORAQ+eAHXYMOLXCmm6KDgzjP58O18/Mx5/jN2ZYWuJMGBhC6eOTIdIvd50t8E3Pm05F2jhJxHaisxQIZDs+fwcqQyLVJKiz2W0XkbWvtOyIixpjvisj9IkJhJ/+D+rRF8DoyApdFNApxi0YhGn4/xPzSJYivJgH5/a5MgBbG0mJawSCs/mAQTwR33OGeBNQi1Y1LrQ6pIYPhsEvx15h4bTId38M0HIaYa2SNLgJ9fRB1jcePzyAdG4NLZWICvvRg0FnMWkJ4bAyCOj6OMEjt8DQygjmeP4+9hrKyxDcxZ/rFa2vx5DI2hnDQkych0kVFCMUMh92xs/nPmUi09kiFsFeKyJW49xERuS0F4xIPoWVjGxud60RDG7UOTEsLBGvvXljzWtZ2YMCl2efkwK0xMOAaR+TnY/wnnsB3589DvMbHsWCo0OXkuGtrBEpFBYS7sdE1lR4bw3jafHpiwvm/dZOyshLXDgRcKV1jkFA0MeEWg74+Vwu+tRXRN9dd5wR0bMw13dDs0D17XK/VgYHpnZziSbRHaTiMzdzmZsz/rrtcTfbh4enn0n/uDVZs89QY87CIPCwisnnz5pW6LFkllJZCRN7zHrx/+mmIrlruwSAs0k2bEP3xqU/Bcu3qcmn//f0QsqEhFyY5OAjh3rQJ7hYtJaCNrIuLXX10a2Fd5+biO81irahwwhoIYNxwGNdTURcR2bwZ9V7y8lw998JCzOUb34AAl5Q4AdbiYrt3I9QwKwst/crLsbFaVeWKmtXXu8gcjcnv78dCcPgw3DSBAMoH1NY6F0+iPUozMnBfZ1re9J97k1QIe7OIbIp7X3X1s2lYaw+JyCERkbq6unmacBEvMlNsRkYgUBodYgws9c2bIT7hMGqjZGa6/qHWOis7MxOCm5uL74aH8fm5cxAo7Z60YQPGrKhwrodoFIJZVAQXiXYyKi3FolFR4Zp2aG/UcBgbmjMLZqk//dgxXKeyEk8ew8MYPyMDYl9djaeIigoXEz86ioUuFsPr0aMu67S/H79FXUTbtmHOWo1y1y6MMZsIJ+oXp//cu6RC2BtEpNYYs0Ug6L8rIh9PwbjEQ8wUkXAYIrhliytU1dbmapns3YtQQa2/kpfnQhY3bnSdkIJBCOLwMOKyH3/clfKdnIRVXlAAwd63z7WHy86Gy2Z0FCI6MACrf3TURcv4/Vg4Sksx1sWLmO9s1q3GxufnO/FvbcUGsLbey8pymaebN8M1cv48znvuOYyrC5Y22I5Gce+qq/Fbo1Ec09w8uwUef78TEWj6z71J0sJurZ00xvyxiDwtCHf8N2vtqaRnRjxHvIhoSN7FixA7jTyprMQm5MDA9LR7rb0+Po5zd+7E36OjELyaGrg1NPNTI1G01ntODhaP7m40pT51CsfcfDNE+bXXsNjs2IFa7fn5OLewEIvIzp141aYiM63bm2/GQqSRLyMjcAlp44+iIkTQaAu/m25yPnatO79jh9tU7ezEHEIhfC6CDVSNwGlvpyCTuUmJj91a+6SIPJmKscjaZLGFotSC/9rXXLs49R0/9xyO+ehH4YuPxSC+WtgrM9NFn+za5cIEX3wR4qfFtXp6XAckrZDY14dCXxUVWFTOnMHY1dUQ+oEBbF729eF627ejDrw2oHjggdl/z4ED+F4bW+vG8E03OXdNeTmeRC5fdovYwYMuKqW/HwtCSYlbsKqrMZYIviso4AYnWRhmnpKkWWqhKC1UdfvtrriWiIslr62FpXvhAjYfQyEIoQg+a293LomnnoIVHm8BHzvmmllkZeHzyUmRH/8Y7eb8fiwMExMIr2xrwwZsdTUiV971rumiO5+YavPo+MUtGp09kuW661yDjPjfrf1eH3jA7UnEYvhdWlSsupobnGRhKOwkaZIpFDVbyrpaqCIQVo01LytzQhtforaqClb6Aw/gmjk5cIdcuQILes8e12BCa5qfOwdLf3IS7hFNdIrFYM0XFrong0SjRWb6q7u6rnXP9PQgzl6jWuZK1Y/fk9AyCUVFWJxYNpcsBIWdJE0yhaJmC80Lh12TaP1My/jGEy+EGk55222w5nt63KIRDDrf9eAgfNatrRDx8+ed6yQQgABrga25/OmJMpt75vrrXdz7QqGG3NgkS4XCTpImmUJRs4XcPfQQvmtoQG/T7m74nd96C5+rPzxeCHWBKCyEuA8OYiF49VVnpfv9sNK3b4eYVlbie58P4p+bi99QVATXzF/8RXL3ZTb3TLy1zVBDslxQ2EnSzGZ1X7wIl8ahQ3Nvps7ccD148NpjWlrwWW4uNj/feguW+fbt1wphZqbIL38J67yuTuTBBxFB8thjcLmUl0PYMzKwMfvqq9gUzcmBuI+OYjHQ0MVUMJ/VTYucLBcUdpI0M61urXmene2EfuZmaiIbrjN991u2wHIPhaZHp8SPtX+/s+ZFsKBs2DC71bxhAzYmX3sNn9fWulZ1d965cvePkFRDYScpId76PHwY4jvfZmoiG66J+u4XGmsuy1ibcHzrW66Zh7aziy+bS8hag8JOUk4igpzIMYn67pPZvK2qgsuGzZqJl6Cwk5STiCAnckyiTR6S7fJDXzfxGr6FDyFkcdTXQ4D7+xEHrk2W6+sXd4z67jXsMBSaPekpkbEIWU8Ym6rt/0VQV1dnjx07tuLXJStHIiUGFluGINnrEbLWMcYct9bWLXQcXTFkWUjEvZFKF0iqxuICQbwAXTGEXEXDJoeHsRk7PIz3kUi6Z0bI4qCwE3KV+LBJbe5cWIjPCVlLUNgJuUpn5+zVGDs70zMfQpYKhZ2Qq2jYZDysfU7WIhR2Qq7CsEniFSjshFwl0bh5QlY7DHckJA5moRIvQIudEEI8BoWdEEI8BoWdEEI8BoWdEEI8BoWdEEI8BoWdEEI8BoWdEEI8BoWdEEI8BoWdEEI8BoWdEEI8BoWdEEI8BoWdEEI8BoWdEEI8BoWdEEI8BoWdEEI8BoWdEEI8RlLCboz5P8aYs8aYN40xh40xhamaGCGEkKWRrMX+jIjsttbuFZHzIvJXyU+JEEJIMiQl7Nban1lrJ6++fUVE2FSMEELSTCp97H8kIk+lcDxCCCFLYMFm1saYZ0WkYpavPm+t/dHVYz4vIpMi8vg84zwsIg+LiGzevHlJkyWEELIwCwq7tfYD831vjHlIRO4VkXustXaecQ6JyCERkbq6ujmPI4QQkhwLCvt8GGP2i8jnROQua+1waqZECCEkGZL1sf+ziOSJyDPGmNeNMV9PwZwIIYQkQVIWu7X2+lRNhBBCSGpg5ikhhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHgMCjshhHiMlAi7MeazxhhrjAmnYjxCCCFLJ2lhN8ZsEpEPisjl5KdDCCEkWVJhsX9VRD4nIjYFYxFCCEmSpITdGHO/iDRba99I0XwIIYQkiX+hA4wxz4pIxSxffV5EHhW4YRbEGPOwiDwsIrJ58+ZFTJEQQshiMNYuzYNijNkjIs+JyPDVj6pEpEVEbrXWts13bl1dnT127NiSrksIIesVY8xxa23dQsctaLHPhbX2pIiUxV3wkojUWWu7ljomIYSQ5GEcOyGEeIwlW+wzsdbWpGosQgghS4cWOyGEeAwKOyGEeAwKOyGEeAwKOyGEeAwKOyGEeAwKOyGEeAwKOyGEeAwKOyGEeAwKOyGEeAwKOyGEeAwKOyGEeAwKOyGEeAwKOyGEeAwKOyGEeAwKOyGEeAwKOyGEeAwKOyGEeIwlN7NO6qLGdIpI04pfeDphEWF/VsB74eC9cPBeOFbLvai21pYudFBahH01YIw5lki37/UA74WD98LBe+FYa/eCrhhCCPEYFHZCCPEY61nYD6V7AqsI3gsH74WD98Kxpu7FuvWxE0KIV1nPFjshhHgSCruIGGM+a4yxxphwuueSLowx/8cYc9YY86Yx5rAxpjDdc1ppjDH7jTHnjDFvG2P+Mt3zSRfGmE3GmOeNMaeNMaeMMX+S7jmlG2NMhjHmNWPMT9M9l0RY98JujNkkIh8UkcvpnkuaeUZEdltr94rIeRH5qzTPZ0UxxmSIyL+IyAER2Skiv2eM2ZneWaWNSRH5rLV2p4i8W0T+1zq+F8qfiMiZdE8iUda9sIvIV0XkcyKyrjcbrLU/s9ZOXn37iohUpXM+aeBWEXnbWvuOtXZcRL4rIveneU5pwVrbaq09cfXvAYGgVaZ3VunDGFMlIh8WkW+key6Jsq6F3Rhzv4g0W2vfSPdcVhl/JCJPpXsSK0yliFyJex+RdSxmijGmRkRuEpGj6Z1JWvlHgfE3le6JJIo/3RNYbowxz4pIxSxffV5EHhW4YdYF890La+2Prh7zecGj+OMrOTey+jDG5IrID0TkT621/emeTzowxtwrIh3W2uPGmLvTPZ9E8bywW2s/MNvnxpg9IrJFRN4wxojA9XDCGHOrtbZtBae4Ysx1LxRjzEMicq+I3GPXXxxss4hsintfdfWzdYkxJiAQ9cettT9M93zSyB0icp8x5qCIBEUk3xjzHWvtJ9I8r3lhHPtVjDGXRKTOWrsaCv2sOMaY/SLyDyJyl7W2M93zWWmMMX7BpvE9AkFvEJGPW2tPpXViacDA0vmWiHRba/803fNZLVy12B+x1t6b7rksxLr2sZNp/LOI5InIM8aY140xX0/3hFaSqxvHfywiTws2C/9rPYr6Ve4QkT8Qkfdf/d/C61ctVrJGoMVOCCEegxY7IYR4DAo7IYR4DAo7IYR4DAo7IYR4DAo7IYR4DM8nKBGyEMaYmIicFJGAIOv22yLyVWvtmkkhJyQeCjshIiPW2htFRIwxZSLyHyKSLyJfSOusCFkijGMn6x5jzKC1Njfu/XWCzNPwOiytQDwAfeyEzMBa+46IZIhIWbrnQshSoLATQojHoLATMoOrrpiYiHSkey6ELAUKOyFxGGNKReTrIvLP9K+TtQo3T8m6Z5Zwx/8nIv/AcEeyVqGwE0KIx6ArhhBCPAaFnRBCPAaFnRBCPAaFnRBCPAaFnRBCPAaFnRBCPAaFnRBCPAaFnRBCPMb/BywbnyBl5XzeAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"X_sampled.plot(x='D', y='Y', style='bo', alpha=0.3, xlim=(-5,5), ylim=(-5, 5))"
]
},
{
"cell_type": "code",
"execution_count": 33,
"metadata": {},
"outputs": [],
"source": [
"from statsmodels.nonparametric.kernel_regression import KernelReg\n",
"import matplotlib.pyplot as pp"
]
},
{
"cell_type": "code",
"execution_count": 34,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Text(0,0.5,'Y')"
]
},
"execution_count": 34,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAEWCAYAAABmE+CbAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAGXBJREFUeJzt3X20XHV97/H39xwSuDwVIcEAITm29VZyfUBuVmJ8KKlggoqCtaxrDSrUa6DVCrdIAlItVQs0bYG7roI5leJFaC1eLKLFRUCb5QMDGJ5qgbYiDQUkGAEFAQkJ3/vH3uMZTs7MmZyTM3vmzPu11qycmb1n5nv2TOZzvr/f3nsiM5EkaaDqAiRJ3cFAkCQBBoIkqWQgSJIAA0GSVDIQJEmAgdD3IuKEiPhOw/WfR8Svtlj/rohYOsU1DUVERsQuU/k8nRSFSyPi8Yi4pep66iJiaUQ8OIn7fzYiPrYza2ryPJOqU+0xELpURLw7IjaUH9APR8TXI+L1U/28mblnZt5X1vD5iPjUqOX/LTPXT3Ud09DrgTcBczNzUdXFTMToPx4AMvPkzPxkVTWNZaw61R4DoQtFxB8BFwLnAC8G5gEXAcdUWZcmZT6wMTOfqroQqanM9NJFF+BXgJ8Dx7VYZ1eKwPhRebkQ2LVcthR4EDgN+DHwMHBiw333A64BngBuAT4JfKdheQK/DqwEngO2lPV8tVy+EThyJ9TxVuD2so4HgLMblg2VdezS5PffCJwO/DPwFHAJRXB+HXgSuAF4UcP6rwFuBH4K3AksbVh2InBPeb/7gJMalrX8Hcao68By2z4G3At8oLz9/cAvgG3ltvzTMe57AvBd4NPAz4B/BY4Y77HLZWcD/w/4+/L3uA141ejXtOH654FPNf6ODcvOAH5YPs7dwDvK2w8Z9Tv8dPRjldc/UNb3WFnvgaPqOBn4QflafAaIJtvyv5SP/XhZx+mTrLPp+81Lw3avugAvo14QOArYSpMPw3KdTwA3AfsDs8sPu0+Wy5aW9/8EMAN4C/A05Qck8EXgSmAP4OXAQ4wRCOXPL/jPXt62kZFAmEwdS4FXUHSprwQeAY4tlw0xfiDcRBECB1F8WN8GvBrYDfgm8CflugcBj5bPP0AxbPMoMLtc/lbg14AADi9rPKyd32GMur5F0cntBhwKbAbeWC47oXE7j3HfE8rn+l/lc/0PimDYt43HPpsivH+nvO9HgP8AZox+TUe/rmwfCMdRhM9AWcNTwAHNfodRj/VG4CfAYRR/LPwf4Fuj3ltfA/ah6Ho3A0c12R7nAd8G9gUOBv5lknUupcn7zUvDdqq6AC+jXhBYAWwaZ50fAm9puL6cYjii/sZ/hoYPU4oPzNcAg+UHx8salp3DxANhQnU0+Z0uBC4ofx5i/EBY0XD9KuDihut/CFxd/rwa+MKo+18HvK/JY18NnLKjv0P5obUN2KvhtnOBz5c/b/chNer+J1B0WdFw2y3Ae9p47LOBmxqWDVB0M28Y/ZqOfl0ZFQhj1HUHcEyz32HUY10CrGlYtmf5fhtqqOP1DcuvBM5o8rz30RAWFB3rhOts9X7zMnJxDqH7PArMGmcPmwOB+xuu31/e9svHyMytDdefpvjPORvYhaJlbrzvRE20DiJicUT8U0RsjoifUQwlzNqB536k4ednxri+Z/nzfOC4iPhp/UIxwXtAWcebI+KmiHisXPaWUXU0/R1GORB4LDOfbLjtfooOpV0PZflp1XD/A9t87F++ppn5PMVQV+Nr0ZaIeG9E3NGwrV5O+6/LC94PmflzivdzY52bGn5uti3rj9X0fbqjde6E91tfMBC6Tw14Fji2xTo/ovigq5tX3jaezRTDEgePum8z2WLZZOoA+FuKMeaDM/NXgM9SDNvsbA9QdAj7NFz2yMzzImJXiu7iL4EXZ+Y+wLUTrONHwL4RsVfDbfMohuTadVBEND53fXu289i/fE0jYgCYy8hr8TSwe8O6c8Z68oiYD/w18CFgv3J7/Asj22OH3g8RsQfFnNWObIO6h2nyPp1gnZ16v/U0A6HLZObPgI8Dn4mIYyNi94iYUf4lu6Zc7e+AP46I2RExq1z/8jYeexvwZeDs8nEXAO9rcZdHgKbHJEy0jtJeFH/1/iIiFgHvbvN+O+py4G0RsTwiBiNit3Kf9rnATIqx7s3A1oh4M7BsIk+SmQ9QzKGcWz7HKykmk9vdHlDMxXy4fL2Po5ggvbbNx/7vEfHbZWd5KsUfFTeVy+4A3l3+/kdRzJWMZQ+KD9PNABFxIsVf3nWPAHMjYmaT+/8dcGJEHFqG7TnAzZm5cQe2Qd2VwJkR8aLytfrDSdbZqfdbTzMQulBm/hXwR8AfU7zpH6D4a+jqcpVPARso9rL5PsWE6qe2f6QxfYiiTd9EMf57aYt1LwEWlG351WMsn0wdfwB8IiKepAiSK9u83w4pP0yPAT7KyLY8HRgoh2A+XD734xQfEtdM4ul+l2L+40fAP1BMbN+wA/e/GXgpxcTsnwG/k5mPtvnYX6GYXH2cYt7htzPzuXLZKcDbKPbsWcHI++gFMvNu4K8outRHKCZhv9uwyjeBu4BNEfGTMe5/A/Axiq7rYYrJ+ne1/du/0J9SDBP9B7AO+MIk6+zI+63XxQuHLCVVISJOAP5nZu7wwYcRcTbFpPHxO7su9Rc7BEkSYCBIkkoOGUmSADsESVKpp04vPGvWrBwaGqq6DEnqKbfeeutPMnP2eOv1VCAMDQ2xYcOGqsuQpJ4SEW2dkcAhI0kSYCBIkkoGgiQJMBAkSSUDQZIEGAiSpJKBIEkCDARJUslAkCQBBoIkqWQgSJIAA0GSVDIQJEmAgSBJKhkIkiTAQJAklQwESRLQBYEQEYMRcXtEfK3qWiSpn1UeCMApwD1VFyFJ/a7SQIiIucBbgc9VWYckqfoO4UJgFfB8sxUiYmVEbIiIDZs3b+5cZZLUZyoLhIg4GvhxZt7aar3MHM7MhZm5cPbs2R2qTpL6T5UdwuuAt0fERuCLwBsj4vIK65GkvlZZIGTmmZk5NzOHgHcB38zM46uqR5L6XdVzCJKkLrFL1QUAZOZ6YH3FZUhSX7NDkCQBBoIkqWQgSJIAA0GSVDIQJEmAgSBJKhkIkiTAQJAklQwESRJgIEiSSgaCJAkwECRJJQNBkgQYCJKmkVoNzj23+Fc7ritOfy1Jk1WrwRFHwJYtMHMmfOMbsGRJ1VX1FjsESdPC+vVFGGzbVvy7fn3VFfUeA0HStLB0adEZDA4W/y5dWnVFvcchI0k9p1YrOoClS0eGhZYsKYaJRt+u9hkIknpKq7mCJUsMgslwyEhST3GuYOrYIUjqeo1DRPW5gnqH4FzBzmMgSOpqYw0ROVcwNQwESV2rVoOzz4Znn4Xnnx8ZIjrzTINgKhgIkrpOrQaXXQaXXgrPPVeEwcCAQ0RTzUCQ1FXqQ0S/+AVkFrcNDMCRRxbdgp3B1HEvI0ldpb4XUT0MImDXXQ2DTjAQJHWV0Uccn3SS5yXqFIeMJHUVjziujoEgqeOGh+Gqq+Cd74SVK7df7hHH1TAQJHXM6tXFnkObNxfX160r/h0rFNR5ziFI6ojly2HNmpEwqLvqqmrq0fYMBElTqlaDd7xjpBsY7Z3v7Gw9as4hI0lTZngYPvhB2Lp1+2UvehGcd57DRd3EQJC009VqcMYZ8O1vjxxP0GhgAP7xH5047jaVBUJEHAxcBrwYSGA4M/93VfVI2jmWLx97eGhgAH7jN4rLqlWGQTeqskPYCpyWmbdFxF7ArRFxfWbeXWFNkiaoVoP3vQ9+8IPtl0XAxRc7PNTtKptUzsyHM/O28ucngXuAg6qqR9LEDA/DgQfCa187dhgAnH66YdALumIOISKGgFcDN1dbiaQd0Wx4qG6//eCccwyDXlF5IETEnsBVwKmZ+cQYy1cCKwHmzZvX4eokjaVWg+OOg4cear7OsmVw3XWdq0mTV+lxCBExgyIMrsjML4+1TmYOZ+bCzFw4e/bszhYoaTvDw8XwULMwOPRQuPFGw6AXVRYIERHAJcA9mXl+VXVIak99ruCkk8Zevvvuxd5Dt9/uHkS9qsoho9cB7wG+HxF3lLd9NDOvrbAmSWMYb65g0SK42RnAnldZIGTmd4Co6vkljW94uNhD6IntZvdGHHKIYTBdVD6pLKk7jdcVzJxZTCxffnnnatLUMhAkvUCtBm96Ezz1VPN1VqwwCKYjA0HSLw0Nwf33N18+dy5ceaWTxtOVp7+WxOrVxeklWoXBihXwwAOGwXRmhyD1ucWL4ZZbmi/fc89iLsEgmP4MBKlPHX88XHFF63Xmz4eNGztSjrqAQ0ZSH1q8uHUYRMDatYZBv7FDkPrI8DB8+MPw7LPN1znkELjbk9D3JTsEqQ/UasUeQied1DwMZs4sugLDoH/ZIUjT3HgHmIGnnlDBDkGapoaHYXCwdRjMn19857FhILBDkKal8XYlBbsCbc8OQZpGVq8uuoJWYTBnjl2BxmaHIE0Te+8NTz7Zeh2PK1ArdghSjxseLo4bGC8MVq0yDNSaHYLUwxYsgHvuab2OcwVqlx2C1IMWLCi6glZhsM8+xXcbGwZql4Eg9Zi99x6/K1i1Ch5/3BPSacc4ZCT1iNWrYc2a1uvMmAFbtnSmHk0/BoLUA/bYA55+uvU6noNIk+WQkdTFFi8u5grGC4NMw0CTZyBIXahWK4JgvKONZ8wowkDaGRwykrpMO7uSgkGgnc9AkLrILrvAtm2t15kzBx5+uDP1qL84ZCR1gaGhYoioVRjstVfRFRgGmioGglSxgQG4//7W66xdC0880Zl61L8cMpIqUqvBG97Qei5g333h0Uc7V5P6mx2CVIHddoPXvrb1ENGKFYaBOssOQeqg/faDxx4bfz33IFIV7BCkDqifonq8MKh/paVUBTsEaYoNDY0/aewX16gb2CFIU6TeFYwXBitWGAbqDnYI0hSYOROee671Oh5gpm5jhyDtRPVzEI0XBsuWGQbqPnYI0k7Szimq/b4CdbNKO4SIOCoi/i0i7o2IM6qsRZqoelcwXhjceKNhoO7WtEOIiGuBP8jMjVPxxBExCHwGeBPwIPC9iLgmMz2ru3pGO3MF4K6k6g2tOoRLgXURcVZEzJiC514E3JuZ92XmFuCLwDFT8DzSlGhnrmDFCsNAvaNph5CZX4qIrwMfAzZExBeA5xuWnz/J5z4IeKDh+oPA4tErRcRKYCXAvHnzJvmU0uQdcABs2tR6nQh4/vnW60jdZrw5hC3AU8CuwF6jLh2RmcOZuTAzF86ePbtTTyttZ/Xq4oN+vDBYtcowUG9qNYdwFHA+cA1wWGaOM2W2wx4CDm64Pre8Teo67ZyDaK+9PEW1elurDuEs4LjMPGMKwgDge8BLI+IlETETeBdF+Ehdo/7FNeOFgd9XoOmg1RzCG6byiTNza0R8CLgOGAT+JjPvmsrnlHZExPjreLSxppNKj0PIzGsz879m5q9l5p9VWYtUd/zx44fB7rv7dZaafjx1hVSq1Yp5gCuuaL3e2rXw1FOdqUnqJE9dIdH+F9fceCMsWTL19UhVMBDU95wrkAoOGalv1Y8rGI9zBeoXBoL6zvAwDA7CmjWt11u2zNNOqL84ZKS+svfe8OSTrdcZHIStWztTj9RN7BDUNyLGD4NVqwwD9S87BE17y5fDunXjr+fwkPqdHYKmrfoBZuOFwSGHGAYS2CFommpnrmD33T3ATGpkh6Bppb4r6XhhsGKFYSCNZoegaaOdrsAvrpGas0NQz1u8uP09iAwDqTk7BPW0wcHxP+RnzIAtWzpTj9TL7BDUk+p7ELUKg4GBoiswDKT22CGo5yxYAPfc03qdZcvguus6U480XdghqGfUarD//q3DYM6c4vsKDANpx9khqCe00xWsWAGXX96ZeqTpyEBQ19ttN3j22ebLI+C73/WLa6TJcshIXev444uJ4VZhMH9+MbFsGEiTZyCo66xeXewqesUVzc8xNH9+sWzjxo6WJk1rDhmpqzhXIFXHDkFdYXgY9tijdRjUuwLDQJoadgiq3AEHwKZNrddZuxZWruxMPVK/skNQZepnJm0VBi99Kdx4o2EgdYIdgiqxeDHcckvz5bvvDhdcYBBInWQgqKOGh+Hss+Hhh5uvM3++ew9JVTAQ1DGeg0jqbs4haMqtXg277NI6DJYtK/YgMgyk6tghaErttx889ljz5a96FVx8sUcaS93ADkFTon5m0mZhMGtWsffQHXcYBlK3sEPQTrd8Oaxb13y5k8ZSd7JD0E61eHHzMIgovsHMMJC6kx2CJm14GK66CmbPHvvYggh42cvg7rs7X5uk9lUSCBHxF8DbgC3AD4ETM/OnVdSiyRlveMjdSKXeUdWQ0fXAyzPzlcC/A2dWVIcmYbwwWLTIMJB6SSWBkJnrMnNrefUmYG4VdWjiVq8eOwxWrCi6grVr4eabO1+XpInrhjmE3wP+vuoi1J5aDS67rPjAH23ZMk9NLfWyKQuEiLgBmDPGorMy8yvlOmcBW4ErWjzOSmAlwLx586agUrWjVoM1a+CrXy2+snL0N5k5VyD1vikLhMw8stXyiDgBOBo4IrPZFyVCZg4DwwALFy5sup6mRj0IrrmmCIJGEcXlIx+BP//zauqTtPNUtZfRUcAq4PDMfLqKGjS+Wg2WLoUtW7ZfNjgIH/gAvPe9HmksTRdVzSF8GtgVuD4iAG7KzJMrqkVjqNWK01Q3C4OLLvK7CqTpppJAyMxfr+J51Z5aDY44Ap599oW3DwzA299eHG1sVyBNP92wl5G6QH3vobotW4o5g4EBWLgQDjvM4SFpujMQxPAw/P7vj0waz5hRfH8BwMyZcOGFBoHUDwyEPlarwRlnwLe+9cLbt26Fk06CefOKSWXDQOoPBkKfqtXg8MPhuee2XzYw4PCQ1I88/XUfqu9B1CwMLrrIMJD6kR1CHxkehksugdtvh23btl/+m78J551nGEj9ykDoA/Wjja+++oW3DwwU31Ow557w/vd7XIHU7wyEaa5Wg9/6re2PKQDYdVf43OfsCCQVnEOYxmo1OPXUscPg2GPhG98wDCSNsEOYpupHGz/zzAtvP+gg+PjHHR6StD0DYZpav3778xDNnAlf+pJdgaSxOWQ0TS1dWgTA4GDx78knFyFhGEhqxg5hmlqypJgjWL/eo40ltcdAmMaWLDEIJLXPISNJEmAgSJJKBoIkCTAQekKtBueeW/wrSVPFSeUuVz/AbMuWYvdRjy6WNFXsELpc/QCzbduKf9evr7oiSdOVgdDlRh9gtnRp1RVJmq4cMupyHmAmqVMMhB7gAWaSOsEhI0kSYCBIkkoGgiQJMBAkSSUDQZIEGAiSpJKBIEkCDARJUslAkCQBBoIkqWQgSJIAA0GSVKo0ECLitIjIiJhVZR2SpAoDISIOBpYB/1lVDZKkEVV2CBcAq4CssAZJUqmSQIiIY4CHMvPONtZdGREbImLD5s2bO1CdJPWnKfuCnIi4AZgzxqKzgI9SDBeNKzOHgWGAhQsX2k1I0hSZskDIzCPHuj0iXgG8BLgzIgDmArdFxKLM3DRV9UiSWuv4V2hm5veB/evXI2IjsDAzf9LpWiRJIzwOQZIEVNAhjJaZQ1XXIEmyQ5AklQwESRJgIEiSSgaCJAkwECRJJQNBkgQYCJKkkoEgSQIMBElSyUCQJAEGgiSpZCBIkgADQZJUMhAkSYCBIEkqGQiSJAAis3e+tz4iNgP3V1zGLMCv+yy4LUa4LUa4LUZ0y7aYn5mzx1uppwKhG0TEhsxcWHUd3cBtMcJtMcJtMaLXtoVDRpIkwECQJJUMhB03XHUBXcRtMcJtMcJtMaKntoVzCJIkwA5BklQyECRJgIEwKRFxWkRkRMyqupaqRMRfRMS/RsQ/R8Q/RMQ+VdfUaRFxVET8W0TcGxFnVF1PVSLi4Ij4p4i4OyLuiohTqq6pahExGBG3R8TXqq6lHQbCBEXEwcAy4D+rrqVi1wMvz8xXAv8OnFlxPR0VEYPAZ4A3AwuA342IBdVWVZmtwGmZuQB4DfDBPt4WdacA91RdRLsMhIm7AFgF9PWsfGauy8yt5dWbgLlV1lOBRcC9mXlfZm4BvggcU3FNlcjMhzPztvLnJyk+CA+qtqrqRMRc4K3A56qupV0GwgRExDHAQ5l5Z9W1dJnfA75edREddhDwQMP1B+njD8G6iBgCXg3cXG0llbqQ4o/G56supF27VF1At4qIG4A5Yyw6C/goxXBRX2i1LTLzK+U6Z1EMGVzRydrUfSJiT+Aq4NTMfKLqeqoQEUcDP87MWyNiadX1tMtAaCIzjxzr9oh4BfAS4M6IgGKI5LaIWJSZmzpYYsc02xZ1EXECcDRwRPbfgS0PAQc3XJ9b3taXImIGRRhckZlfrrqeCr0OeHtEvAXYDdg7Ii7PzOMrrqslD0ybpIjYCCzMzG44o2HHRcRRwPnA4Zm5uep6Oi0idqGYTD+CIgi+B7w7M++qtLAKRPEX0v8FHsvMU6uup1uUHcJHMvPoqmsZj3MImqxPA3sB10fEHRHx2aoL6qRyQv1DwHUUk6hX9mMYlF4HvAd4Y/leuKP8C1k9wg5BkgTYIUiSSgaCJAkwECRJJQNBkgQYCJKkkgemSRMUEduA7wMzKI7Svgy4IDN75lQFUiMDQZq4ZzLzUICI2B/4W2Bv4E8qrUqaII9DkCYoIn6emXs2XP9ViiOVZ/XhKTw0DTiHIO0kmXkfMAjsX3Ut0kQYCJIkwECQdppyyGgb8OOqa5EmwkCQdoKImA18Fvi08wfqVU4qSxM0xm6nXwDOd7dT9SoDQZIEOGQkSSoZCJIkwECQJJUMBEkSYCBIkkoGgiQJMBAkSaX/D08FPl6ZEmyfAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"model = KernelReg(X['Y'], X['D'], 'c')\n",
"result = model.fit(X['D'])[0]\n",
"\n",
"pp.plot(X['D'], result, 'b.'); pp.xlim(-5,5); pp.ylim(-5,5); pp.title(\"Conditional mean of population data\"); pp.xlabel('D'); pp.ylabel(\"Y\")"
]
},
{
"cell_type": "code",
"execution_count": 35,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Text(0,0.5,'Y')"
]
},
"execution_count": 35,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAEWCAYAAABmE+CbAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAF0VJREFUeJzt3X2UXHWd5/H3N51AiAQCJIgQQnRwd2HExTnZID6RIwxBZAaU0RXNKKNj8Dgo7KAwyOyAjjvgw0DYgyi1sijCjONMFF1HlvDU7vFMAwYEPTzouBgGMGjkGSRAwnf/+P16u9L0c1J9q9Pv1zl1uqruvVXfulV9P/X7/e69FZmJJEkzmi5AktQdDARJEmAgSJIqA0GSBBgIkqTKQJAkAQaChhERJ0bED9puPxURrxhh/jsjYlmHa1ocERkRMzv5PJMpissi4tGIuKXpeoayNet9vMtGxFci4tPjr1LbgoEwxUTEuyNibd1Ar4+IqyPiDZ1+3szcOTPvrTW86J82M383M3s7Xcd26A3A7wMLM3Np08VMJRHRGxF/2nQd2xMDYQqJiD8HVgF/A7wUWARcDBzbZF3aKvsB6zLz6aYLkchML1PgAuwKPAW8Y4R5dqQExi/rZRWwY522DHgAOA34NbAe+JO2ZfcAvgM8AdwC/DXwg7bpCewPrASeB56r9fyvOn0dcMQ2qOOtwI9qHfcD57RNW1zrmDnM618HfBz4MfA0cCklOK8GngSuA3Zrm/+1wL8AjwF3AMvapv0JcHdd7l7gpLZpI76GIerau67bR4CfAx+s938A2Ahsruvyk0Msuz/wfeBx4DfAP7RNu7CuoyeAW4E3tk07B/hH4Ir6Gn4C/DvgzFrz/cCRbfP3AufW9/4J4NvA7kOtd8pn8dL6uh8EPg301Gk9wOdrrfcCfzbKe/Ya4LZa4z8AXwc+XaftBnwX2AA8Wq8vrNP+W11vG+u6u2i0deJlDNuZpgvwMsY3Co4CNg33j1Xn+RRwE7AnsKBu7P66TltWl/8UMAs4GvgtdQNZ/xG/AbwEeFX9R39RINTrX+n/p22bvo6BQNiaOpYBB1Far68GfgUcV6dtsWEa4vWvq8/7UmAfyobvtrrRmQ3cAJxd590HeLg+/wxKt83DwII6/a3A7wABHFZr/L2xvIYh6vo/lJbcbODguoF7c512Yvt6HmLZvwfOqjXOBt7QNm0FJchnUsLpIWB2nXYOZWO5vE6/HPhFfaxZwAeBX7Q9Vm99z19VPwOrgSuGWu/At4BL6nx7UkLkpDrtQ8A9wL7A7sCNw71nwA7AfcB/qTX9EeXLRn8g7AEcD8wB5lIC7qpBNf/poMccdp14GcN2pukCvIzxjYL3AA+NMs//BY5uu72c0h3RvxF7pv0fk7LBfC3lW93zwH9om/Y3TDwQJlTHMK9pFXBBvb7FhmmIedcB72m7vRr4Ytvtj/RvUIAzgK8NWv4a4H3DPPZVwCnjfQ11w7gZmNt237nAV+r1Exk5EC4HWtRvxqO8/48C/7FePwe4tm3aH1C+Sfd/k59b1+W8ersXOK9t/gMprcCe9vVOCdtngZ3a5j0BuLFevwH4UNu0I4d7z4A3UVqQ0Xbfvwz+bLVNOxh4tO12L4MCYaR14mX0i2MIU8fDwPxR9tbYm/KNq9999b7//xiZuant9m+BnSnf4mdSmtrty07UROsgIg6JiBsjYkNEPE75xjl/HM/9q7brzwxxe+d6fT/gHRHxWP+FMsD7slrHWyLipoh4pE47elAdw76GQfYGHsnMJ9vuu4/SQhmL0ymtlFvqnlzv758QER+LiLsj4vFa466Dahz82n+TmZvbbjOo5sHv/yxevO73q/evb1tvl1BaClBe71g/R3sDD2bdcg+ePyLmRMQlEXFfRDxBaWnNi4ie4R5wDOtEIzAQpo4+yjez40aY55eUf9h+i+p9o9lA6QLZd9CywxntFLkTrQPg7yj97ftm5q7AlygbxG3tfkoLYV7b5SWZeV5E7EhpXXweeGlmzgO+N8E6fgnsHhFz2+5bROmeGVVmPpSZH8zMvYGTgIsjYv+IeCMlLN5J6aqaRxln2Jp1Nfj9f54yFtDufsrncH7betslM3+3Tl8/xOMMZz2wT0S019w+/2nAvwcOycxdKC0KGHiNW3wOO7ROphUDYYrIzMeBvwK+EBHH1W9Ps+o32c/W2f4e+MuIWBAR8+v8V4zhsTcD3wTOqY97IPC+ERb5FTDsMQkTraOaS/lGvTEilgLvHuNy43UF8AcRsTwieiJidkQsi4iFlL7tHalBGRFvoXR9jFtm3k/pBjm3PserKYPJY1ofEfGOWhOU7o8EXqCsp021xpkR8VfALhOpsc2KiDgwIuZQxkf+qa1F0f961gNrgL+NiF0iYkZE/E5EHFZn+Qbw0YhYGBG7AX8xwvP11dfw0fpZfjvQvuvtXEpL5rGI2B04e9Dygz+HnVgn04qBMIVk5t8Cfw78JeVDfz9wMqV/G8reHmspe9n8hDKgOtaDfE6mdB88RBkjuGyEeS8FDqxdBlcNMX1r6vgw8KmIeJISJN8Y43LjUjfUxwKfYGBdfhyYUbt3Plqf+1FKKH1nK57uBEo//C8pA7JnZ+Z1Y1z2PwE3R8RTtYZTshwPcg3wv4GfUbpZNrJlV81EfI3y3j9EGcD+6DDzvZcSmndR1s8/UbvagP9Ra7uD8r5/c7gny8zngLdTxlEeAf7zoPlXATtRWik3UV5vuwuBP6oH9f13OrNOppXYsvtO0nQUEb2UvYq+3HQtao4tBEkSYCBIkiq7jCRJgC0ESVI1pU4jPH/+/Fy8eHHTZUjSlHLrrbf+JjMXjDbflAqExYsXs3bt2qbLkKQpJSLGdOYBu4wkSYCBIEmqDARJEmAgSJIqA0GSBBgIkqTKQJAkAQaCJKkyECRJgIEgSaoMBEkSYCBIkioDQZIEGAiSpMpAkCQBBoIkqTIQJElAFwRCRPRExI8i4rtN1yJJ01njgQCcAtzddBGSNN01GggRsRB4K/DlJuuQJDXfQlgFnA68MNwMEbEyItZGxNoNGzZMXmWSNM00FggRcQzw68y8daT5MrOVmUsyc8mCBQsmqTpJmn6abCG8HvjDiFgHfB14c0Rc0WA9kjStNRYImXlmZi7MzMXAu4AbMnNFU/VI0nTX9BiCJKlLzGy6AIDM7AV6Gy5DkqY1WwiSJMBAkCRVBoIkCTAQJEmVgSBJAgwESVJlIEiSAANBklQZCJIkwECQJFUGgiQJMBAkSZWBIEkCDARJUmUgSJIAA0GSVBkIkrYbfX1w7rnlr8avK34xTZK2Vl8fHH44PPcc7LADXH89HHpo01VNLbYQJG0XentLGGzeXP729jZd0dRjIEiacobqGlq2rLQMenrK32XLmqpu6rLLSNKUMlzX0KGHluu9vSUM7C4aPwNBUtfr6xvY0A/VNdS/8e8PBk2MgSCpqw1uEaxaVf7237ZraNsxECR1rb4+OOccePZZeOGFEgIPP2zXUKcYCJK6Tl8fXH45XHYZPP98CYMZMwZaBHYNdYaBIKmr9HcRbdwImeW+GTPgiCNKa8Eg6Bx3O5XUVfoHjfvDIAJ23NEwmAwGgqSuMvh4gpNO8qjjyWKXkaRJ12rB6tVw/PGwcuWW0zyeoDkGgqRJ0WrB2WfDI4+ULiGANWvK36FCwSCYfHYZSeqovj5429tK189DDw2EQb/Vq5upSy9mC0FSR7RacOGFcM89ZbfR4Rx//OTVpJEZCJK2qRUryrf+jRtHnm+vveCTn3xxd5GaYyBI2ib6+uB974N//dehp8+YAfPmlaOOjzsOrrhicuvT6BobQ4iIfSPixoi4KyLujIhTmqpF0sS1WrD33vC6140cBl/8YjntxFNPGQbdqskWwibgtMy8LSLmArdGxLWZeVeDNUkag74++Oxn4YYb4IknRp73TW+C885zr6GpoLFAyMz1wPp6/cmIuBvYBzAQpC7VHwRXXTX6vK98JXz1qwbBVNIVYwgRsRh4DXBzs5VIGkpfH3z4w3D77aPPe/DBcPHFBsFU1HggRMTOwGrg1Mx8UeMzIlYCKwEWLVo0ydVJWrECrrxy9PnmzIGTT4bPfKbzNakzGj0wLSJmUcLgysz85lDzZGYrM5dk5pIFCxZMboHSNNXXB4cdVs4lNFoY7LADnH46PP20YTDVNdZCiIgALgXuzszzm6pD0oDxdA3Nnl0OKnOPoe1Hky2E1wN/DLw5Im6vl6MbrEeats44o5xi+nWvGz0Mdt8dLrkEnnnGMNjeNLmX0Q+AaOr5JcHy5XDttQO/PTCSXXaBz33OI4u3Z57cTpqG+vpgt93K2UZHC4OenjJG8PjjhsH2rvG9jCRNnlYLPvYxePLJsc1/5JFwzTWdrUndw0CQtnN9ffDOd8IDD4xt/hkz4IQTHB+YjuwykrZTrRbssUcZKB4tDCJg8eIyWLx5s2EwXdlCkLYzrRZ8/OOjn2Oo33veYwCoMBCk7cR49hiKgP3391xD2pKBIE1xBx4Id9899vltEWg4jiFIU9CKFeVAsoixhUFE2WMo0zDQ8AwEaQo54wyYNaucX2jwj9UPZe7cMlD8wgvuPqrR2WUkTQHj7RaaMwcuuMADyTQ+thCkLrZ8+di7hWCgRfD004aBxs8WgtRlli+H668vxwOM1Yc+BO99r3sMaesYCFKXWLwY7rtvfMvsuCNs3NiRcjQN2WUkNajVgl13Ld1CYw2DiLLraKZhoG3LFoLUgIm0Bg44AO66qyPlSIAtBGnStB87MJ4w2Guv0howDNRpBoLUYWecUUJgrMcOAOy3XwmBTFi/vrP1Sf3sMpI6YLzHDfSbM6fsMio1wRaCtI20WuXXxcZz3ACULqFLLimtAcNATbKFIG2lQw6BW24Z3zKzZsFFF3nwmLqLgSBNwPLlcOON8Pzz41/Wn6VUt7LLSBqjVgv23rt0Ca1ZM74wOOCAgUFiw0DdykCQRtBqwUteUkLgpJPGt8fPnDkDYwPuMqqpwC4jaQgTGRcA2Gkn+MhH4DOf2fY1SZ1mIEjV8uXw/e+Xk8pt2jS+ZT2KWNsDu4w0bfX1wdveNtAltGYNPPvs2MJgxoyBXyCzS0jbC1sImlZaLVi1Cp55ppw+Yiw/SN/u4IPh4os9zbS2TwaCtnutFpx5Jjz66PgDoN/SpXDzzdu2LqnbGAjaLq1YAatXl98SHuv5g/rNmwfvepc/OKPpx0DQdqO/O+gXv5jY7wTMnQuf/7xHD2v6MhA0pbVacM458PDD428JLF1aWgPHH28ISGAgaIrp64PeXrjzTvjnf4bHHhvf8jvtBAcdBB/4gCEgDWYgaEro7w766U/LuMB4HHBAOc3E29/uAWPSSAwEdZ3+VsCyZWVQt9Uqp40Yj5e9rBxtfPrpDgxLY2UgqCv0h8Aee8Cpp5bxgB12gOuvL3sLjWb2bHjTm0qI9AeJpPFpNBAi4ijgQqAH+HJmntdkPZpc/SHw2GNwwQXllBE9PeVv/+6ivb1l0HfNmoHlZsyAI44oRwe/4hVw3nkGgLQtDBsIEfE94MOZua4TTxwRPcAXgN8HHgB+GBHfyUxPAjAN9PXB4YeXU0UMHhOYMaOcSmKHHbb8tn/ppeX003YDSZ0xUgvhMmBNRHwV+GxmTuCnQEa0FPh5Zt4LEBFfB44FDIRpoLe3tAAGh0FPT/klsYcf3jIMVq50ryCp04YNhMz8x4i4GvivwNqI+BrwQtv087fyufcB7m+7/QBwyOCZImIlsBJg0aJFW/mU6hbLlpUWQH8LIWIgDNzwS80YbQzhOeBpYEdgLm2BMFkyswW0AJYsWTLBM9Go2xx6aBkw7h9IHtwikDT5RhpDOAo4H/gO8HuZ+dtt/NwPAvu23V5Y79M0ceihBoDUTUZqIZwFvCMz7+zQc/8QeGVEvJwSBO8C3t2h55IkjWKkMYQ3dvKJM3NTRJwMXEPZ7fR/djB8JEmjaPQ4hMz8HvC9JmuQJBX+hKYkCTAQJEmVgSBJAgwESVJlIEiSAANBklQZCJIkwECQJFUGgiQJMBAkSZWBIEkCDARJUmUgSJIAA0GSVBkIkiTAQJAkVQaCJAkwECRJlYEgSQIMBElSZSBIkgADQZJUGQiSJMBAkCRVBoIkCTAQJEmVgSBJAgwESVJlIEiSAANBklQZCJIkwECQJFUGgiQJMBAkSVUjgRARn4uIeyLixxHxrYiY10QdkqQBTbUQrgVelZmvBn4GnNlQHZKkqpFAyMw1mbmp3rwJWNhEHZKkAd0whvB+4Oqmi5Ck6W5mpx44Iq4D9hpi0lmZ+e06z1nAJuDKER5nJbASYNGiRR2oVJIEHQyEzDxipOkRcSJwDHB4ZuYIj9MCWgBLliwZdj5J0tbpWCCMJCKOAk4HDsvM3zZRgyRpS02NIVwEzAWujYjbI+JLDdUhSaoaaSFk5v5NPK8kaXjdsJeRJKkLGAiSJMBAkCRVBoIkCTAQJEmVgSBJAgwESVJlIEiSAANBklQZCJIkwECQJFUGgiQJMBAkSZWBIEkCDARJUmUgSJIAA0GSVBkIkiTAQJAkVQaCJAkwECRJlYEgSQIMBElSZSBIkgADQZJUGQiSJMBAkCRVBoIkCTAQJEmVgSBJAgwESVJlIEiSAANBklQZCJIkwECQJFWNBkJEnBYRGRHzm6xDktRgIETEvsCRwL81VYMkaUCTLYQLgNOBbLAGSVLVSCBExLHAg5l5xxjmXRkRayNi7YYNGyahOkmanmZ26oEj4jpgryEmnQV8gtJdNKrMbAEtgCVLltiakKQO6VggZOYRQ90fEQcBLwfuiAiAhcBtEbE0Mx/qVD2SpJF1LBCGk5k/Afbsvx0R64Almfmbya5FkjTA4xAkSUADLYTBMnNx0zVIkmwhSJIqA0GSBBgIkqTKQJAkAQaCJKkyECRJgIEgSaoMBEkSYCBIkioDQZIEGAiSpMpAkCQBBoIkqTIQJEmAgSBJqgwESRIAkTl1frc+IjYA9zVcxnzAn/ssXBcDXBcDXBcDumVd7JeZC0abaUoFQjeIiLWZuaTpOrqB62KA62KA62LAVFsXdhlJkgADQZJUGQjj12q6gC7iuhjguhjguhgwpdaFYwiSJMAWgiSpMhAkSYCBsFUi4rSIyIiY33QtTYmIz0XEPRHx44j4VkTMa7qmyRYRR0XETyPi5xHxF03X05SI2DciboyIuyLizog4pemamhYRPRHxo4j4btO1jIWBMEERsS9wJPBvTdfSsGuBV2Xmq4GfAWc2XM+kioge4AvAW4ADgRMi4sBmq2rMJuC0zDwQeC3wZ9N4XfQ7Bbi76SLGykCYuAuA04FpPSqfmWsyc1O9eROwsMl6GrAU+Hlm3puZzwFfB45tuKZGZOb6zLytXn+SsiHcp9mqmhMRC4G3Al9uupaxMhAmICKOBR7MzDuarqXLvB+4uukiJtk+wP1ttx9gGm8E+0XEYuA1wM3NVtKoVZQvjS80XchYzWy6gG4VEdcBew0x6SzgE5TuomlhpHWRmd+u85xF6TK4cjJrU/eJiJ2B1cCpmflE0/U0ISKOAX6dmbdGxLKm6xkrA2EYmXnEUPdHxEHAy4E7IgJKF8ltEbE0Mx+axBInzXDrol9EnAgcAxye0+/AlgeBfdtuL6z3TUsRMYsSBldm5jebrqdBrwf+MCKOBmYDu0TEFZm5ouG6RuSBaVspItYBSzKzG85oOOki4ijgfOCwzNzQdD2TLSJmUgbTD6cEwQ+Bd2fmnY0W1oAo35C+CjySmac2XU+3qC2Ej2XmMU3XMhrHELS1LgLmAtdGxO0R8aWmC5pMdUD9ZOAayiDqN6ZjGFSvB/4YeHP9LNxevyFrirCFIEkCbCFIkioDQZIEGAiSpMpAkCQBBoIkqfLANGmCImIz8BNgFuUo7cuBCzJzypyqQGpnIEgT90xmHgwQEXsCfwfsApzdaFXSBHkcgjRBEfFUZu7cdvsVlCOV50/DU3hoO+AYgrSNZOa9QA+wZ9O1SBNhIEiSAANB2mZql9Fm4NdN1yJNhIEgbQMRsQD4EnCR4weaqhxUliZoiN1Ovwac726nmqoMBEkSYJeRJKkyECRJgIEgSaoMBEkSYCBIkioDQZIEGAiSpOr/AcRkVOakYD/dAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"model = KernelReg(X_sampled['Y'], X_sampled['D'], 'c')\n",
"result = model.fit(X_sampled['D'])[0]\n",
"\n",
"pp.plot(X_sampled['D'], result, 'b.'); pp.xlim(-5,5); pp.ylim(-5,5); pp.title(\"Conditional mean of sampled data\"); pp.xlabel('D'); pp.ylabel(\"Y\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## The sampled conditional mean is biased due to conditioning on a descendant of a collider involving an error term!! "
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.5.2"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment