synthetic_experiment2_with_A_as_feature.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Fair Dummies: Synthetic Example\n",
    "\n",
    "This notebook implements a the fair dummies framework for learning predictive models that approximately satisfy the equalized odds notion of fairness.\n",
    "\n",
    "Paper: \"Achieving Equalized Odds by Resampling Sensitive Attributes,\" Y. Romano, S. Bates, and E. J. Candès, 2020\n",
    "\n",
    "### Proposed approach\n",
    "\n",
    "__Core idea__: fit a regression function, minimizing\n",
    "\n",
    "$$ \\text{loss = prediction error + distance to equalized odds}$$\n",
    "\n",
    "\n",
    "__Input__: $  \\{(X_i,A_i,Y_i)\\}_{i=1}^n \\sim P_{XAY}$ training data\n",
    "    \n",
    "\n",
    "__Step 1__: sample dummy protected attributes\n",
    "$$\n",
    "\\tilde{A}_i \\sim P_{A|Y}(A \\mid Y=y_{i}) \\quad \\forall \\ i=1,2,\\dots,n \\nonumber\n",
    "$$\n",
    "\n",
    "$A \\in \\{0,1\\}$? generate $\\tilde{A}$ using a biased coin-flip, with\n",
    "$$\n",
    "P\\{A=1|Y=y\\} = \\frac{P\\{y \\mid A=1\\}P\\{A=1\\}}{P\\{y \\mid A=1\\}P\\{A=1\\} + P\\{y \\mid A=0\\}P\\{A=0\\}}\n",
    "$$\n",
    "    \n",
    "__Step 2__: fit a regression function on $\\{(X_i, A_i, Y_i)\\}_{i=1}^n$\n",
    "\n",
    "$$\n",
    "        \\hat{f}(\\cdot) \\,= \\, \\underset{f \\in \\mathcal{F}}{\\mathrm{arg min}} \\, \\frac{1-\\lambda}{n} \\sum_{i=1}^n (Y_i - f(X_i))^2  + \\lambda \\mathcal{D}\\left( [\\hat{\\mathbf{Y}}, \\mathbf{A}, \\mathbf{Y}] , [\\hat{\\mathbf{Y}}, \\tilde{\\mathbf{A}}, \\mathbf{Y}] \\right) \\nonumber\n",
    "$$\n",
    "\n",
    "where\n",
    "    \n",
    "$$\n",
    "\\hat{\\mathbf{Y}} = \\left[f(X_{1}), f(X_{2}), \\dots, f(X_{n})\\right]^T \\ ; \\ \\mathbf{A} = \\left[A_{1}, A_{2}, \\dots, A_{n}\\right]^T \\ ; \\ \\tilde{\\mathbf{A}} = \\left[\\tilde{A}_{1}, \\tilde{A}_{2}, \\dots, \\tilde{A}_{n}\\right]^T \\ ; \\ \\mathbf{Y} = \\left[Y_{1}, Y_{2}, \\dots, Y_{n}\\right]^T\n",
    "$$\n",
    "\n",
    "and $\\mathcal{D}\\left( \\mathbf{Z}_1, \\mathbf{Z}_2 \\right)$ tests whether $P_{Z_1} = \\ P_{Z_2}$ given $\\mathbf{Z}_1, \\mathbf{Z}_2$,  here it implemented as a classifier two-sample test"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.preprocessing import StandardScaler\n",
    "\n",
    "import matplotlib\n",
    "matplotlib.use('Agg')\n",
    "%matplotlib inline\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "np.warnings.filterwarnings('ignore')\n",
    "import random\n",
    "\n",
    "import torch\n",
    "import torch.nn as nn\n",
    "import pandas as pd\n",
    "\n",
    "import os\n",
    "import sys\n",
    "\n",
    "\n",
    "sys.path.append(os.path.abspath(os.path.join(os.getcwd() + '/others/third_party/cqr')))\n",
    "base_path = os.getcwd() + '/data/'\n",
    "    \n",
    "from fair_dummies import fair_dummies_learning\n",
    "from fair_dummies import utility_functions\n",
    "\n",
    "seed = 123\n",
    "\n",
    "random.seed(seed)\n",
    "np.random.seed(seed)\n",
    "torch.manual_seed(seed)\n",
    "\n",
    "# maximal number of testpoints to plot\n",
    "max_show = 100\n",
    "\n",
    "# save figures?\n",
    "save_figures = False\n",
    "\n",
    "# display results\n",
    "if save_figures:\n",
    "    SMALL_SIZE = 26\n",
    "else:\n",
    "    SMALL_SIZE = 18\n",
    "\n",
    "MEDIUM_SIZE = SMALL_SIZE\n",
    "BIGGER_SIZE = SMALL_SIZE\n",
    "\n",
    "plt.rc('font', size=SMALL_SIZE)          # controls default text sizes\n",
    "plt.rc('axes', titlesize=SMALL_SIZE)     # fontsize of the axes title\n",
    "plt.rc('axes', labelsize=MEDIUM_SIZE)    # fontsize of the x and y labels\n",
    "plt.rc('xtick', labelsize=SMALL_SIZE)    # fontsize of the tick labels\n",
    "plt.rc('ytick', labelsize=SMALL_SIZE)    # fontsize of the tick labels\n",
    "plt.rc('legend', fontsize=SMALL_SIZE)    # legend fontsize\n",
    "plt.rc('figure', titlesize=BIGGER_SIZE)  # fontsize of the figure title\n",
    "plt.rc('text', usetex=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_func(x,\n",
    "              y,\n",
    "              y_u=None,\n",
    "              y_l=None,\n",
    "              pred=None,\n",
    "              point_color=\"b.\",\n",
    "              shade_color=\"\",\n",
    "              method_name=\"\",\n",
    "              title=\"\",\n",
    "              filename=None,\n",
    "              save_figures=False,\n",
    "              show_legend=True):\n",
    "    \n",
    "    \"\"\" Scatter plot of (x,y) points along with the constructed prediction interval \n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    x : numpy array, corresponding to the feature of each of the n samples\n",
    "    y : numpy array, target response variable (length n)\n",
    "    pred : numpy array, the estimated prediction. It may be the conditional mean,\n",
    "           or low and high conditional quantiles.\n",
    "    shade_color : string, desired color of the prediciton interval\n",
    "    method_name : string, name of the method\n",
    "    title : string, the title of the figure\n",
    "    filename : sting, name of the file to save the figure\n",
    "    save_figures : boolean, save the figure (True) or not (False)\n",
    "    \n",
    "    \"\"\"\n",
    "    \n",
    "    inds = np.random.permutation(x.shape[0])[:max_show]\n",
    "    x_ = x[inds]\n",
    "    y_ = y[inds]\n",
    "    if y_u is not None:\n",
    "        y_u_ = y_u[inds]\n",
    "    if y_l is not None:\n",
    "        y_l_ = y_l[inds]\n",
    "    if pred is not None:\n",
    "        pred_ = pred[inds]\n",
    "\n",
    "    fig = plt.figure()\n",
    "    inds = np.argsort(np.squeeze(x_))\n",
    "    \n",
    "    if (y_u is not None) and (y_l is not None):\n",
    "        plt.fill_between(x_[inds].squeeze(),\n",
    "                         y_u_[inds],\n",
    "                         y_l_[inds],\n",
    "                         alpha=.2,\n",
    "                         color=shade_color,\n",
    "                         edgecolor=\"\",\n",
    "                         label = u'Prediction interval')\n",
    "    \n",
    "    if pred is not None:\n",
    "        if pred_.ndim == 2:\n",
    "            plt.plot(x_[inds,:], pred_[inds,0], 'k', lw=2, alpha=0.2,\n",
    "                     label=u'Predicted low and high quantiles')\n",
    "            plt.plot(x_[inds,:], pred_[inds,1], 'k', lw=2, alpha=0.2)\n",
    "        else:\n",
    "            plt.plot(x_[inds,:], pred_[inds], point_color, lw=2, alpha=0.8, markersize=15,\n",
    "                     fillstyle='none')\n",
    "    min_val = -8\n",
    "    max_val = 8 \n",
    "    \n",
    "    ident = [min_val, max_val]\n",
    "    plt.xlim(ident[0], ident[1])\n",
    "    plt.ylim(ident[0], ident[1])\n",
    "    \n",
    "    plt.plot(ident,ident, ls=\"--\", c=\".3\")\n",
    "    \n",
    "    plt.xlabel(r'$Y$')\n",
    "    plt.ylabel(r'$\\hat{Y}$')\n",
    "    \n",
    "    ax = plt.gca()\n",
    "    ax.set_aspect('equal', 'box')\n",
    "    \n",
    "    plt.title(title)\n",
    "\n",
    "\n",
    "    if save_figures and (filename is not None):\n",
    "        plt.savefig(filename, bbox_inches='tight', dpi=300)\n",
    "    \n",
    "    plt.show()\n",
    "    \n",
    "def plot_groups_pointwise(x_axis_0,x_axis_1,A,Y,Yhat,\n",
    "                          y_u=None,y_l=None,point_color=\"b.\",shade=\"\",\n",
    "                          filename_0=None,filename_1=None,save_figures=False):\n",
    "   \n",
    "    inds = A==0\n",
    "    \n",
    "    y_u0 = None\n",
    "    y_l0 = None\n",
    "    if (y_u is not None) and (y_l is not None):\n",
    "        y_u0=y_u[inds]\n",
    "        y_l0=y_l[inds]\n",
    "        \n",
    "    plot_func(x_axis_0.reshape(x_axis_0.shape[0],1),\n",
    "              Y[inds], y_u=y_u0, y_l=y_l0, pred=Yhat[inds],\n",
    "              point_color=point_color,shade_color=shade, method_name=\"\", title=r'$A=0$',\n",
    "              filename=filename_0, save_figures=save_figures,show_legend=True)\n",
    "\n",
    "    inds = A==1\n",
    "    \n",
    "    y_u1 = None\n",
    "    y_l1 = None\n",
    "    if (y_u is not None) and (y_l is not None):\n",
    "        y_u1=y_u[inds]\n",
    "        y_l1=y_l[inds]\n",
    "        \n",
    "    plot_func(x_axis_1.reshape(x_axis_1.shape[0],1),\n",
    "              Y[inds], y_u=y_u1, y_l=y_l1, pred=Yhat[inds],\n",
    "              point_color=point_color,shade_color=shade, method_name=\"\", title=r'$A=1$',\n",
    "              filename=filename_1, save_figures=save_figures,show_legend=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Generate data\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(500, 3)\n",
      "(500,)\n"
     ]
    }
   ],
   "source": [
    "def init_seed(seed = 0):\n",
    "    random.seed(seed)\n",
    "    np.random.seed(seed)\n",
    "    torch.manual_seed(seed)\n",
    "\n",
    "def synthetic_example(include_A=False, n = 6000):\n",
    "        \n",
    "    p0 = 0.9\n",
    "    p1 = 1 - p0\n",
    "\n",
    "    A = np.random.binomial(1,p1,n).T\n",
    "    X = np.random.randn(n,2)\n",
    "\n",
    "    beta0 = [2,1]\n",
    "    beta1 = [1,2]\n",
    "\n",
    "    Y = np.random.randn(n)\n",
    "    Y[A==0] = Y[A==0] + np.dot(X[A==0],beta0)\n",
    "    Y[A==1] = Y[A==1] + np.dot(X[A==1],beta1)\n",
    "\n",
    "    x_axis_0 = Y[A==0]\n",
    "    x_axis_1 = Y[A==1]\n",
    "    \n",
    "    if include_A:\n",
    "        X = np.concatenate((X,A[:,np.newaxis]),axis=1)\n",
    "    \n",
    "    return X, A, Y, x_axis_0, x_axis_1, beta0, beta1\n",
    "\n",
    "include_A = True\n",
    "\n",
    "init_seed(seed)\n",
    "X, A, Y, x_axis_0, x_axis_1, beta0, beta1 = synthetic_example(include_A=include_A, n = 500)\n",
    "X_cal, A_cal, Y_cal, x_axis_0, x_axis_1, beta0, beta1 = synthetic_example(include_A=include_A, n = 2000)\n",
    "X_test, A_test, Y_test, x_axis_0_test, x_axis_1_test, beta0, beta1 = synthetic_example(include_A=include_A, n = 2000)\n",
    "\n",
    "# compute input dimensions\n",
    "n = X.shape[0]\n",
    "in_shape = X.shape[1]\n",
    "    \n",
    "print(X.shape)\n",
    "print(Y.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Estimate $P_{A|Y}$ using kernel density estimator, then compute $P\\{A=1|Y=y\\}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYEAAAD+CAYAAADcWrmEAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAALCklEQVR4nO3dz20b2ZYH4HMGHQBb3dp6DIYgKwQ7A+uFIIUwIRh2Bq0QWs6gHYKs/SyegFn3g5qLtz9vwdKYVlNiiRJFk+f7AC2qyHvrAmXUr+4fXmdVBQA9/de2GwDA9ggBgMaEAEBjQgCgMSEA0NhP227AY/3666/1+vXrbTcDYKd8/fr1X1V1ePf8zoXA69ev4/LyctvNANgpmfl/y84bDgJoTAgANCYEABoTAgCNCQGAxoQAQGNCAKAxIQDQmBAAaGznfjEMP6qLi+1d++Rke9dmt+kJADQmBAAaEwIAjQkBgMZGTQxn5mlE3AyH06r6NKLMJCL+ERHvqurkzmdvI+IsIj5ExCwi3kfErKrOH9F2AJ5oZU/gNgCq6nNVfY6Iz5n524oyRxHxNubBMV3ylclw/uvw94sAAHh5Y3oCZ1X15vagqq6HN/l7VdVVRFwNYXDfd97c9xkAL+PBnsAwpLPsQT5bFQQA/PhW9QSmMR+zv+u+YZ7RhhCZDPUfjZlnAOB5rQqBg/g2IbxoFvMH+LquIuZDSxERmXmTmX9U1btlXx7mJU4jIl69evWEywKwaCtLRKvq+jYAhuOriDjOzKW9i6o6r6rjqjo+PDx8sXYC7LsxIXCw5NxTegH3uY75iiIAXsiqELiM5Q/8gxiGdB4rM6eZ+dc6ZQF4Xg/OCVTVLDOvM3NSVYsTxJOq+vKE635Ycm4aEU+pEyJiu7t5wq4ZMxz0MYZJ2Yj//yHYl4XjaWZeDMtJ7/rbUNIwF/DdiqPMfB8Rvy/OEwCweSt/LFZV55l5urCkc1pVZwtfmcZ8LP8ghof7MMH7PiLeRcRRZn6MiH/e/ir4ts6h/GQ4t1gnAC9g1N5BD23pMAwL/Xzn3HVEfBr+Hl0nAC/DLqIAjQkBgMaEAEBjQgCgMSEA0JgQAGhMCAA0JgQAGhMCAI0JAYDGhABAY0IAoDEhANCYEABoTAgANCYEABoTAgCNCQGAxoQAQGNCAKAxIQDQmBAAaEwIADQmBAAaEwIAjQkBgMaEAEBjQgCgMSEA0JgQAGhMCAA0JgQAGhMCAI0JAYDGftp2A4Cnu7jYznVPTrZzXZ6PngBAY0IAoDEhANCYEABoTAgANCYEABoTAgCNCQGAxoQAQGNCAKAxIQDQmBAAaEwIADQmBAAaEwIAjQkBgMaEAEBjQgCgsVH/vWRmnkbEzXA4rapPI8pMIuIfEfGuqv72n9CtUycAz2tlCNw+rKvq83A8zczfqursgTJHETGN+UN++hx1AvD8xvQEzqrqze1BVV1n5tuHClTVVURcDWHwLHUC8PwenBMYhnSWPchn6z60N1EnAOtZNTE8jYjZkvNLh3lG2kSdAKxhVQgcxLfJ20WziJisec1H15mZp5l5mZmXf/7555qXBeCunVgiWlXnVXVcVceHh4fbbg7A3hgzMXyw5Ny6vYBN1skP5OJi2y0AxljVE7iM5Q/ng4i4WvOam6gTgDU8GAJVNYuI62FFz6JJVX1Z54KbqBOA9YyZE/gYEae3B8Pa/y8Lx9PMvFjyUI9YPuyzsk4AXsbKOYGqOh9W57yN+TDO9M4ve6cR8TbmD/xZxDwYIuJ9RLyLiKPM/BgR/6yq85F1AvACsqq23YZHOT4+rsvLy203gxVMDLNJJ3/bjYxVMvNrVR3fPb8TS0QB2AwhANCYEABoTAgANCYEABoTAgCNCQGAxoQAQGNCAKAxIQDQmBAAaEwIADQmBAAaEwIAjQkBgMaEAEBjQgCgMSEA0JgQAGhMCAA0JgQAGhMCAI0JAYDGhABAY0IAoDEhANCYEABoTAgANCYEABoTAgCNCQGAxoQAQGNCAKAxIQDQmBAAaEwIADQmBAAaEwIAjQkBgMaEAEBjQgCgMSEA0JgQAGhMCAA0JgQAGhMCAI0JAYDGhABAY0IAoDEhANCYEABoTAgANPbTmC9l5mlE3AyH06r69JQymfk2Is4i4kNEzCLifUTMqur8EW0H4IlW9gRuH+ZV9bmqPkfE58z87YllJhExjYivw98vAgDg5Y3pCZxV1Zvbg6q6Ht7kn1Rm8XMAtuPBnkBmTiLiaMlHs/uCYJ0yAGzHqp7ANOZj9nfdDJ+tXWYIhMnw3aMx8wwAPK9VIXAQ3yZ3F81i/gBft8xVxHyYKCIiM28y84+qereswmGO4TQi4tWrVyuaDMBYW1kiWlXXtwEwHF9FxHFmLu1dVNV5VR1X1fHh4eGLtRNg340JgYMl5+7rBTylzHVEmDMAeEGrQuAylj+8D2IY0nlsmcycZuZf45sIwKY8GAJVNYuI62HFz6JJVX15QpkPS4pOI2JpnQBsxpjhoI8xTMpGRGTmUSw8rIc3+4s7D/17ywxzAd+tHsrM9xHx++I8AQCbt/LHYlV1npmnC0s6p1V1tvCVaczH8g9ieLivKnP7+XA4Gc4t1gnACxi1d9BDWzoMQzw/P6bMmM8B2Dy7iAI0JgQAGhMCAI0JAYDGhABAY0IAoDEhANCYEABoTAgANCYEABobtW0EwI/k4mJ71z452d61N0FPAKAxIQDQmBAAaEwIADQmBAAaszpoz21zFQXw49MTAGhMCAA0JgQAGhMCAI0JAYDGhABAY0IAoDEhANCYEABoTAgANCYEABoTAgCNCQGAxoQAQGNCAKAxIQDQmBAAaEwIADQmBAAaEwIAjQkBgMaEAEBjQgCgMSEA0JgQAGhMCAA09tO2GwCwSy4utnPdk5PN1KsnANCYEABoTAgANCYEABoTAgCNCQGAxiwRfQHbWlIGsIqeAEBjrXoC3sgBvjcqBDLzNCJuhsNpVX16apl16gTgea0MgduHdVV9Ho6nmflbVZ2tW2adOgF4fmN6AmdV9eb2oKquM/PtE8usUycAz+zBieHMnETE0ZKPZvc9tFeVWadOADZjVU9gGhGzJedvhs/WKXOzRp0AbMCqEDiIb5O3i2YRMVmzzKPrHOYQTofDf2fm/z7QZiJ+jYh/bbsRjOZ+7Y5dvlf/vezkTiwRrarziDjfdjt2RWZeVtXxttvBOO7X7tjHezXmx2IHS87d1wsYW2adOgF4ZqtC4DKWP5wPIuJqzTLr1AnABjwYAlU1i4jrYUXPoklVfVmnzDp18miGznaL+7U79u5ejRkO+hjfJmUjM48i4svC8TQzL+481B8sM+JznmCYQ2FHuF+7Yx/vVVbV6i/NV+dcx3wY57stHoa1/RcR8aaqrseUGfM5AJs3KgQA2E87sUQU9pFNFHfDMNpxFhEfYv57pvcRMduXoSH/nwBsweImisNGip8z87dtt4ulJjHfzeDr8PfLvgRAhJ7A3tv3t5gdZhPFHbJ4r/aNnsD+2+u3mF1kE0V+JHoCDezzW8yOWmdjRrZoCOdJzO/b0T7N3wgBeHnrbMzI9lxFzIfsIiIy8yYz/6iqd9tt1vMQAg3s81sMbNri75+G46vMPM7M6d3PdpEQ2H97/Razw2yiuNuuI+Jt7ME2EiaG91xVXS++rVTVVUQcZ6ax5+2xieKOGLbF+Wvb7dgkPYEdMQzp/M/Ir5+t6KbuzVvMLqqqWWZeZ+Zk2FDxlk0Uf0wflpybxp7sdyYEdsTwcHjUP7rhbf9rVf28mVbxBLebKH6KsInij2r4/cZ3K7ky831E/L4P8wER9g7aa0MIvF+yed9fcWfDP16eTRR3x3CvIoZhvH26V0Jgz2Xm6eKPw4a3mHdVdbbFZgE/CCHQwD6/xQBPIwQAGrNEFKAxIQDQmBAAaEwIADQmBAAaEwIAjQkBgMb+A7nZiVyx/AjuAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYEAAAD+CAYAAADcWrmEAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAWZklEQVR4nO3d6XJc1bnG8WclQEgqcGQZQyqpGKdTGSqppIitG0hh7sDnXIJ8B+QSKLgDfAnBX/PJrlQIkASQBFgIeZBb82hLao3W2O/5sHbHbbmlnnvtvdf/V9UF3a299dLY++m9RmdmAgDE6XuhCwAAhEMIAEDECAEAiBghAAARIwQAIGIvhC6gWa+99ppdunQpdBkAkCnDw8OPzezCydczFwKXLl3S0NBQ6DIAIFOcc9O1Xqc5CAAiRggAQMQIAQCIGCEAABEjBAAgYoQAAESMEACAiBECABCxzE0WA/KsXJamp6XFRenVV6Xf/U76Hl/V0EWEABBIuSytr/vH1JR/bXdX2t/3/z4/Lz16JP3pT1JfX7AykXOEANBjZtLSkvTNN9LW1tk/+/ixdOuW9Mc/Sr/5TW/qQ1wIAaCHNjakL7/03/6bceeO9MMfShcvdqcuxIsQAHrg4EAqFqWxMd8M1Iovv5ReeUU6d66ztSFuhADQRUdHvtlnZsb/ezvKZWl4WHr7bcm5ztQHMO4A6BIz6fPP/R1AuwFQsb7u+xOATmkoBJxzg865a8nj3QaP6UuO+6jGe1edcx855y475wrOuXedc4PNFg+k2f370sJCd84LdErdEEguzmtmdtPMbkq66Zz7sM4xlyVdlbQmqVDjR/qS14eTx3kzu9Fs8UBara5Ko6PdOffKirS52Z1zIz6N9AlcN7MrlSdmVnTOXT3rADMbkTSShMFpP3PltPeALFtclL74wjcHdcvUlB82CrTrzDsB51yfpFoX8lK9IABic3wsjY9Ln33mRwN109xcd0MG8ah3J1CQVKrx+mnNPA1LQqQvOf9lM/ugnfMBvWbm2/zX1vxM36Wl7l/8K3Z2/O89f743vw/5VS8E+uUv+CeV5C/grRqRfNOSJDnn1pxzt8zsnVo/nPRLDErSRWbLIKCtLWlvz7f5z801P+mrk6anCQG0L8g8gcrFv+r5iHNuwDlXOPle8v4NSTckaWBggJtg9JyZ9PXX0sRE6Eqemp2V3nqLBebQnkb++PTXeK0by1kV5UcUAalycCB9+mm6AkDydS0uhq4CWVcvBIZU+4Lfr6RJp1nJvICAN9FA48x8R29aJ2hVVh8FWnVmCJhZSVIxGSVUrc/Mbrfxe9+r8VpBUjvnBDpufNyv5JlWi4u+jwJoVSPNQe8r6ZSV/jsR7HbV80Iy+/e0O4ZnJG3+z4w4cs5dk/TXWv0BQCgLC37BtzQz8+sSAa2q2zFsZjeS5R8qQzoLZna96kcK8m35/Uou7s65gqRrkt6RdNk5976kh5VZwZVzJsf3Ja9VnxMIqlSS/v3v0FU0ZnJS+tWvWFQOrXGWsRknAwMDNjQ0FLoM5Nwnn6S3H6CWP/9Zeu210FUgzZxzw2Y2cPJ1BpcBJ6ysZCsApPSNXEJ2EAJAlXLZzwfImrk5aXs7dBXIIkIAqPLggd8CMmvMpHv3QleBLCIEgMT2dvpHA52lE7uXIT6EACD/TXpoyK8EmlVHR34pCaAZhAAgvxjbo0ehq2hfkZk2aBIhgOjt7UlffRW6is5YW8tmnwbCIQQQvfHxfLWlT0+HrgBZQgggaisr0sOHoavorJkZdh1D4wgBRGt3V/rXv/J3wXzyJN2L3iFdCAFEqVz2awMdHoaupDsmJ0NXgKwgBBCl8XHfiZpXc3P5DTh0FiGA6Ozt5X927fGxXwobqIcQQHTu3s32pLBGMUoIjSAEEJW9vfyNBjrN8jKLyqE+QgBRuXfPdwrHIpbAQ+sIAURjfz++i+LUVBxNX2gdIYBoTEzEd0E8OJDm50NXgTQjBBCFcjnexdWmpkJXgDQjBBCFxUXfKRyjlZV4/9tRHyGAKMQ8XNLMTx4DaiEEkHsHB/5OIGaTk/lbIwmdQQgg92Zm4hoWWkupxKJyqI0QQO7F3BRU7f790BUgjQgB5Nr2dr4XimvGwoJfPhuoRggg1xge+axYh8nidIQAcqtcZl39k4pF+kfwLEIAubWwwPj4k/b3pdnZ0FUgTQgB5NbMTOgK0mliInQFSBNCALnE3IDTra3RWY6nCAHk0uwsbd9n4W4AFYQAcolRMGebnaW/BB4hgNxZX/czZHG6ctlvswkQAsgdhoU2ZmKC7SdBCCBnjo8ZFdQoM2lsLHQVCI0QQK7Mz0uHh6GryI7ZWZaSiB0hgFyhQ7g5ZowUih0hgNzY2pIePQpdRfZMTkpHR6GrQCiEAHLjwYPQFWTTwQGfXcwIAeTCwQGjgtpx96705EnoKhACIYBcmJpihnA7jo6k0dHQVSAEQgCZR+dmZ0xPswVljAgBZN78vLSzE7qKfPj8c7/cNOJBCCDTmPDUWbu70qef+j4WxIEQQKZNT0ubm6GryJe1Nenjj5l0FwtCAJllJo2Ph64in0olaWTEf8bIN0IAmTU/zwJo3TQz45vaCIJ8IwSQSWbSvXuhq8i/8XHps8/oLM4zQgCZNDvLFom9srjoO4uPj0NXgm54oZEfcs4NSqr8lSuY2QcNHNMn6f8kvWNm/9uJcwKSH8EyPBy6irisrfnJZG+9FboSdFrdEKhcrM3sZvK84Jz70Myun3HMZUkF+Yt8oRPnBCq+/ZYFz0J48EDq65MuXQpdCTqpkTuB62Z2pfLEzIrOuatnHWBmI5JGkjDoyDkByQ8HnZ4OXUW8Rkak8+elV14JXQk65cw+gaRJp9aFvNTqRbsb50Q8mBgW1vGx9MkndBTnSb2O4YKkWlt212zmaVA3zokIbG5Kc3Ohq8DOjjQ0FLoKdEq9EOjX087baiVJfS3+zqbP6ZwbdM4NOeeGHrFrSLRY8z49Fhb8qCFkXyaGiJrZDTMbMLOBCxcuhC4HAezv++WikR6jo0wky4NGQqC/xmut3gV085zIsZkZ9gtIm40NaWUldBVoV70QGFLti3O/pJEWf2c3zomcm50NXQFqYdZ29p0ZAmZWklRMRvRU6zOz2638wm6cE/m2uyutroauArUsL7OKa9Y10hz0vqTBypNk7P/tqucF59xHNS7qUu1mn7rnBKpxF5Bu3A1kW90QMLMbSsbwO+euSbp6YmZvQdJVVV3wk2B4V9JfJF12zr2fzBJu9JzAf83MhK4AZ5meZme3LHOWse79gYEBG2KQcjR2d6W//S10FainUJCuXKn/cwjHOTdsZgMnX8/EEFHEa2kpdAVoxPQ0W1JmFSGAVCMEsuH4mGa7rCIEkFrHx4RAlkxOhq4ArSAEkForK2xkkiWlkrS+HroKNIsQQGpxF5A9LO2RPYQAUosQyB6W98geQgCptLMjbW+HrgLNOjggvLOGEEAqcSHJLnZ+yxZCAKm0vBy6ArRqfp6dx7KEEEDqlMuEQJaZMWcgSwgBpM7qqnR0FLoKtINRQtlBCCB1FhZCV4B2lUosMZ0VhABSxYwQyAuahLKBEECqbG0xNDQvZmbYgzgLCAGkCncB+bGzw45wWUAIIFXm50NXgE56+DB0BaiHEEBq7O1Ja2uhq0Anzc76jYGQXoQAUmNxMXQF6DQz6cGD0FXgLIQAUoOmoHx6+FA6PAxdBU5DCCAVDg+ZJZxXx8esJ5RmhABSYWGBJYjzrFhkuGhaEQJIBb4p5tvGBp3+aUUIILgnT2gKigHDRdOJEEBwLC8Qh9lZv+kM0oUQQFBmNAXFolxmddE0IgQQ1MaGfyAODx/SQZw2hACCoikoLtvb0spK6CpQjRBAMGbS3FzoKtBrzCBOF0IAwZRKfqVJxGVx0S8ZjnQgBBAMdwHxuns3dAWoIAQQBJuRx216ms2D0oIQQBBraywxHDMz6bvvQlcBiRBAIDQFYXqaPqE0IATQc2YsGw1vcjJ0BSAE0HOMCkLF9DSTx0IjBNBzbCaPit1d6fHj0FXEjRBAzy0tha4AacIosbAIAfTU/j7ryuNZc3N+9zGEQQigp9g3ACcdHNBEGBIhgJ6iKQi1sMR0OIQAesaMEEBtS0tMHgyFEEDPlEq+TwCohbuBMAgB9Ax3ATjL5CRzBkIgBNAzhADOsrsrPXoUuor4EALoiYMDaXU1dBVIu4cPQ1cQH0IAPbG8zK0+6pufl/b2QlcRF0IAPcE4cDTCzK8nhN4hBNB15TL9AWjc1BR3jb1ECKDr1tZ8nwDQiM1NlhbppRca+SHn3KCkyv+Wgpl90M4xzrmrkq5Lek9SSdI1SSUzu9FE7cgImoLQrMlJ6fz50FXEoe6dQOVibmY3zeympJvOuQ/bPKZPUkHScPI4TwDk1+Ji6AqQNTMz3D32SiN3AtfN7ErliZkVk2/ybR1T/T7ya2fH394DzTg+lopF6be/DV1J/p15J+Cc65N0ucZbpdOCoJVjkF80BaFV9+9LR0ehq8i/es1BBfk2+5PWkvdaPsY5d9U5dy3557uNFIvsoSkIrdrflyYmQleRf/Wag/r1tHO3Wkm+Xb/VY0Yk30wkSc65NefcLTN7p9YJkz6GQUm6ePFinZKRFkdHLAOA9ty7JxUK0ksvha4kv4IMETWzYiUAkucjkgacczXvLszshpkNmNnAhQsXelYn2rOy4ucIAK06OJC++y50FfnWSAj013jttLuAdo4pSqLPIEdoCkInPHwoPXkSuor8qhcCQ6p98e5X0qTT7DHOuYJzbr3xEpFFZoQAOqNclsbGQleRX2eGgJmVJBWTET/V+szsdhvHvFfj0IKkmudE9mxt8e0NnTM5KW1shK4inxppDnpfSaesJDnnLqvqYp18s//oxEX/1GOSvoBnRg85565J+mt1PwGyjQ3l0Wmjo6EryKe6k8XM7IZzbjAZ498nvwTE9aofKci35fcrubjXO6byfvK0L3mt+pzIOEIAnba46AcbvP566ErypaG1g85a0iFp4jnXzDGNvI/sKpcZGoruuHNHevttybnQleQHq4ii4x4/ZqYnumN9XZqdDV1FvhAC6DhGBaGbxsbYb6CTCAF0lBnrBaG7trfZfayTCAF01NaW/0sKdNPYmF9pFO0jBNBR3AWgF3Z3/UxitI8QQEfRH4BeuXuXAQidQAigY/b3/cggoBf29/2m9GgPIYCO4S4AvfbgASOF2kUIoGPm50NXgNhsb/Plo12EADri6IilIhAGu4+1hxBARywuMmQPYSwvs8JoOwgBdMTMTOgKELO7d0NXkF2EANq2v0+7LMKanZV2dkJXkU2EANo2Pc0IDYRlRt9AqwgBtMWMmZtIh8lJJo+1ghBAW5aXWSsI6XB4yOSxVhACaMuDB6ErAJ5i8ljzCAG0bHNTWloKXQXwFJPHmkcIoGX37oWuAHgefVTNIQTQEjb2QFotLTFctBmEAFry7be0vSK9isXQFWQHIYCmra2x2TfSbXKSZUwaRQigKWbSN9+ErgI42/6+NDcXuopsIATQlKUlNo5BNjBctDGEABpmJo2Ohq4CaMz6OsubN4IQQMNmZ1myF9ly5w53A/UQAmhIucxdALJnY4OlJOohBNCQqSlpdzd0FUDzRkf9ukKojRBAXWbS/fuhqwBas7/v57WgNkIAdS0sSFtboasAWjcx4ee34HmEAM5kJo2Ph64CaN+XXzKBrBZCAGdaWPBD7YCs29yUxsZCV5E+hABOdXzsh9gBeXHvHnMHTiIEcKr799k1DPnzn//w57oaIYCatrel774LXQXQeQcH0j//Ke3tha4kHQgBPMdM+uILP0EMyKOdHR8EBwehKwmPEMBz7tyRVldDVwF018aG9I9/cEdACOAZ4+NMDEM8Njakjz+WnjwJXUk4hAD+a3ycmZWIz+am9Pe/S6VS6ErCIAQgiQBA3HZ3fRAsLISupPcIARAAgPy8mM8+838XYhoUQQhErLJJDAEAPDU+Lt265SeVxbAXwQuhC0AYq6t+FBBbRQLP29z0Q0hffVX65S+lixell14KXVV3EAIR2dvzewQXiwwBBRqxuSl99ZX09dfShQvSG2/4f/7P/0gv5OTqmZP/jDiUy35J5/V135G1v+8v7Pv70o9+JP3gB9LLL0vf/76/jT089O/v7PihcDEPgwPaYSatrPhHxY9/LP30p9KlSz4UsooQSKnjY99Us7Hhh66VSj4AYuqwAtJse9vPqbl/X/rJT6Rf/1p6/XXJudCVNYcQSJG9Pf9NY2nJD1VjSzwgG5aW/KOvz4fBz38ufS8jw24IgYAq3/aXl/0j1skqQF6USn7drdFRHwa/+IX04ouhqzobIRDA6qq/hVxcZKcjII+ePJG++cYPN/3DH3y/QVrvDBoKAefcoKTKDp0FM/ug3WNaOWeWHR35Jp6JCUbmALE4OJCGh/2Xvt//XvrZz9IXBnVDoHKxNrObyfOCc+5DM7ve6jGtnDNLzPzonVLpacfu8rIPAgDx2drym9m89JKfc/Dmm9K5c+noRG7kTuC6mV2pPDGzonPuapvHtHLO4Mz8iIDNTd9pe3z8dJjmyQfNPABOOjjwrQETE36I6Ztv+pFFIQPhzBBwzvVJulzjrZJz7qqZ3W72GElDzZ6z18z8/6ytrWe/zW9scHEH0Bnb237j+7Exf4fw+us+EN54w8/76ZV6dwIFSbXGrKwl77VyzFoL52ybmW+OqUyuqv7WXuv5aePx03D7BiBfDg+l+Xn/kPxyFW+8IfX3+39/+WUfFN24/tQLgX497bytVpLU1+IxTZ8z6UMYlKSLFy+eXfEpnPNDtV58UXrllZZOAQC5k7J+6trM7IaZDZjZwIULF0KXAwC50UgI9Nd47bS7gEaPaeWcAIAOqxcCQ6p9ce6XNNLiMa2cEwDQBWeGgJmVJBWTET/V+k4bxVPvmFbOCQDojkaag95X0ikrSc65y5JuVz0vOOc+OnFRP/OYBt4HAPRA3cliZnbDOTeYjPHvk1/ioXpmb0HSVfnmnFIjxzRwTgBADzS0dpCZ3TjjvduSzjVzTCPvAwC6LxNDRAEA3UEIAEDEnJmFrqEpzrlHkqZD15Ezr0l6HLqICPG5hxHr5/6mmT032zZzIYDOc84NmdlA6Dpiw+ceBp/7s2gOAoCIEQIAEDFCAJLEcN0w+NzD4HOvQp8AAESMOwEAiBghAAARa2jZCACdk+yUV9ldr2BmH4SsJwbJOmXXJb0nv8bZNUkllq/hTgDoqUoAmNlNM7sp6aZz7sPQdUWgT36xy+HkcZ4A8LgTwDP4xtR1183sSuWJmRWTzxxdVv254ynuBHAS35i6JNlz43KNt0oEAULhTgDP4RtT1xSU7LlxwlryHrqoav+SkqTL9MV4hADQO/162iFcraTa+26jc0Yk3/wmSc65NefcLTN7J2xZ4RECeA7fmJA3lYt/1fMR59yAc65w8r3YEAI4iW9M3dVf4zXuAsIoym+NG3WfFx3DeIaZFau/GZnZiKQB5xxt1u0bUu0Lfr+S8EXnOecKzrn10HWkFXcCOZU06fylwR+/XueWmG9MHWBmJedc0TnXZ2bVHcR9yV7d6J73arxWkBT9504I5FRyUWnqD3jybX/YzM51pypIel/SoKQPJMk5d1lciLoqmYvxzKgs59w1SX+NvT9AYhVRVElC4NrJjuDkVvoKf2E6I5k1XFQyJ4OO995IPncpaZLjc/cIATzDOTdYPTks+cb0jpldD1gWgC4hBPAcvjEB8SAEACBiDBEFgIgRAgAQMUIAACJGCABAxAgBAIgYIQAAESMEACBi/w8rWFEhhB854AAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from sklearn.neighbors import KernelDensity\n",
    "\n",
    "def compute_density(Z,A,show_graphs=False,Z_test=[]):\n",
    "\n",
    "    bandwidth = np.sqrt(np.median(np.abs(Z)))\n",
    "\n",
    "    kde_0 = KernelDensity(kernel='linear', bandwidth=bandwidth).fit(Z[A==0][:, np.newaxis])\n",
    "    kde_1 = KernelDensity(kernel='linear', bandwidth=bandwidth).fit(Z[A==1][:, np.newaxis])\n",
    "    \n",
    "    if show_graphs:\n",
    "        plt.clf()\n",
    "        plt.hist(Z[A==0], fc='#AAAAFF', density=True)\n",
    "        plt.show()\n",
    "\n",
    "        X_plot = np.linspace(min(Z[A==0])-1, max(Z[A==0])+1, 1000)[:, np.newaxis]\n",
    "        log_dens = kde_0.score_samples(X_plot)\n",
    "\n",
    "        plt.clf()\n",
    "        plt.fill(X_plot[:, 0], np.exp(log_dens), fc='#AAAAFF')\n",
    "        plt.show()\n",
    "\n",
    "    log_dens_0 = np.exp(np.squeeze(kde_0.score_samples(Z[:, np.newaxis])))\n",
    "    log_dens_1 = np.exp(np.squeeze(kde_1.score_samples(Z[:, np.newaxis])))\n",
    "    p_0 = np.sum(A==0) / A.shape[0]\n",
    "    p_1 = 1 - p_0\n",
    "\n",
    "    # p(A=1|y) = p(y|A=1)p(A=1) / (p(y|A=1)p(A=1) + p(y|A=0)p(A=0))\n",
    "    p_success = (log_dens_1*p_1) / (log_dens_1*p_1 + log_dens_0*p_0)\n",
    "\n",
    "    p_success_test = []\n",
    "    if len(Z_test) > 0:\n",
    "        log_dens_0_test = np.exp(np.squeeze(kde_0.score_samples(Y_test[:, np.newaxis])))\n",
    "        log_dens_1_test = np.exp(np.squeeze(kde_1.score_samples(Y_test[:, np.newaxis])))\n",
    "        p_success_test = (log_dens_1_test*p_1) / (log_dens_1_test*p_1 + log_dens_0_test*p_0)\n",
    "    \n",
    "    return p_success, p_success_test\n",
    "\n",
    "p_success, p_success_test = compute_density(Y,A,True,Y_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Implementation of a fairness-unaware baseline algorithm"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# step size\n",
    "lr = 0.01\n",
    "\n",
    "# inner epochs to fit adversary\n",
    "dis_steps = 1\n",
    "\n",
    "# inner epochs to fit loss\n",
    "loss_steps = 1\n",
    "\n",
    "batch_size = 6000\n",
    "\n",
    "# utility loss\n",
    "cost_pred = torch.nn.MSELoss()\n",
    "in_shape = X.shape[1]\n",
    "out_shape = 1\n",
    "\n",
    "model_type = \"deep_model\"\n",
    "\n",
    "# equalized odds penalty\n",
    "mu = 0\n",
    "second_moment_scaling = 0\n",
    "\n",
    "# total number of epochs\n",
    "epochs = 500\n",
    "\n",
    "base_reg = fair_dummies_learning.EquiRegLearner(lr=lr,\n",
    "                                                pretrain_pred_epochs=0,\n",
    "                                                pretrain_dis_epochs=0,\n",
    "                                                epochs=epochs,\n",
    "                                                loss_steps=loss_steps,\n",
    "                                                dis_steps=dis_steps,\n",
    "                                                cost_pred=cost_pred,\n",
    "                                                in_shape=in_shape,\n",
    "                                                batch_size=batch_size,\n",
    "                                                model_type=model_type,\n",
    "                                                lambda_vec=mu,\n",
    "                                                second_moment_scaling=second_moment_scaling,\n",
    "                                                out_shape=out_shape)\n",
    "\n",
    "init_seed(seed)\n",
    "\n",
    "input_data_train = np.concatenate((A[:,np.newaxis],X),1)\n",
    "base_reg.fit(input_data_train, Y)\n",
    "\n",
    "input_data_cal = np.concatenate((A_cal[:,np.newaxis],X_cal),1)\n",
    "Yhat_out_cal = base_reg.predict(input_data_cal)\n",
    "\n",
    "input_data_test = np.concatenate((A_test[:,np.newaxis],X_test),1)\n",
    "Yhat_out_test = Yhat_test = base_reg.predict(input_data_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Baseline (All): Test Error = 1.0248804225307333\n",
      "Baseline (A=0): Test Error = 1.0075439989539974\n",
      "Baseline (A=1): Test Error = 1.1756597527296029\n"
     ]
    }
   ],
   "source": [
    "# bhat, res, rank, s = np.linalg.lstsq(X, Y)\n",
    "\n",
    "# Yhat_out_cal = np.dot(X_cal,bhat)\n",
    "# Yhat_out_test = np.dot(X_test,bhat)\n",
    "\n",
    "print(\"Baseline (All): Test Error = \" + str(np.sqrt(np.mean((Yhat_out_test - Y_test)**2))))\n",
    "print(\"Baseline (A=0): Test Error = \" + str(np.sqrt(np.mean((Yhat_out_test[A_test==0] - Y_test[A_test==0])**2))))\n",
    "print(\"Baseline (A=1): Test Error = \" + str(np.sqrt(np.mean((Yhat_out_test[A_test==1] - Y_test[A_test==1])**2))))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARkAAAElCAYAAADQqJbBAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO2deXiU5dn2z5uwD0sSQHaQDAFUIsqmFRRlEcWlgiJCodbIEiqfaKWAIGiCgkO17ftaUfA1ClisSAyHay3U2oXyUAJFBtRSSRDFBCIhiEESIPf3xzUP88xkJpkss5+/48gR59nmTsycXPe1Kq01CCEkWDQK9wIIIbENRYYQElQoMoSQoEKRIYQEFYoMISSoUGQIIUGFIkOCjlJqi1IqMdzrIOGhcbgXQGIbpdRdAEYDSAZQGsL3nQmgxPUyRWu9MlTvTTyhJUOCzSTX95RQvaEpMFrrTVrrTQA2KaVWh+r9iScUGRI0lFLzASxwvQzldmmWS1wAAFrrfIg1RcIARYYEBaVUCoBS1wccCJEl4/L9DPRxqlQpRaEJAxQZEixmaa3XuP67FIA9RO+bAt++nxKEcMtG3NDxSxocl8Ww03IoH+L4rem+gP0mWutZfk4lw+3wtVKK0G7ZiAuKDAkGE71EICArohrhIFEMRYY0KK7ITrJSymE5HOptii+riVZMmKDIkAbD5XS1a60neh0/DsDh+64GJw++BSUZwO4QrYFYoMiQhmSm1nqBj+MBJeEppd4M8H1K/G2ttNalSql8pVSi1tr6vola660BPp80IBQZ0iC4Mnvz/ZzOd12TYglpV8HbAqoHDgAzAax0ve9AABSYMMEQNqkXSqkUlwXyJoAq1oXLR2NaN6tDkaviCp2XKqVGm2UNdCqHD8Uev4SQYEJLhhASVCgyhJCgQpEhhAQVigwhJKjERQi7ffv2+uKLLw73MgiJScrKylBaWoqioqJvtdYdvM/HhchcfPHFyMvLC/cyCIk5DMOAw+GA3W7Hs88++6Wva7hdIoTUCavAZGZm+r2OIkMIqTWffPKJh8DYbDa/18bFdokQ0rCkpqZi7NixmDZtWrUCA9CSIYTUgn379uHMmTNo2bIlMjIyahQYgCJDCAkQwzCwZMkSvPbaa7W6jyJDCKkRq5N38uTJtbqXIkMIqRbvKFIgWyQrFBlCiF8qKiqwZs2aOgsMwOgSIaQamjZtimXLliExMbFOAgPQkiGE+MAwDKxbtw5aa3Tt2rXOAgPQkiGEeGH1wVRUVKBZs2b1eh4tGULIBbydvPUVGICWDCHRg9MJ5OYChYVA587A+PFAWlqDPb6+USR/0JIhJBrIzgYWLgR69QLuu0++L1woxxuI8+fPo2/fvg0qMAAtGUIiH6cTyMkB1q8Hkl3DMYcOBW65BZg6FRgypF4WTWlpKRITEzFs2DD86Ec/QqNGDWt70JIhJNLJzQXuucctMCbJycDkyXK+jhiGgRkzZmDPnj0A0OACA1BkCIl8CguBvn19n+vTBygqqtNjTR9Mz549kZqaWo8FVg9FhpBIp3Nn4D//8X3uwAGgU6daPzJYTl5fUGQIiXTGjwf+8AegpMTzeEkJ8Prrcr4WHDp0KGQCA9DxS0jkk5YG3HmnOHknT5Yt0oEDIjB33VVrp2/Pnj1x//3344Ybbgi6wABxMqZ28ODBmo3ESdRj5skUFckWqZZ5Mjt37kTnzp3RrVu3oCxPKbVLaz3Y+zgtGUKihbS0OoeqTR/MkCFDsGjRogZeWPXQJ0NIjGN18s6dOzfk70+RISSGCWUUyR8UGUJiFK013nvvvbAKDECfDCExidYaSiksWrQIlZWVYRMYgJYMITGHOVXg9OnTaNGiRVgFBqDIEBJTmD6YM2fOIFLSUygyhMQIkeDk9QVFhpAYYOfOnREpMAAdv4TEBF26dMGQvn0xt2lT2ObNC0rnvLpCS4aQKKagoEAmCnz4IRY5nbClpgatc15doSVDSJRi+mB+NnYsfvzOO0HrnFdfaMkQEoVYnbyjv/kmaJ3zGoKoEhml1Gil1JtKqYFKqRSl1Hyl1Mxwr4uQUFIlilRcHJTOeQ1FVIkMgEQAKQB2ub7aaa3XhHdJhISOkydP4tfPPAN78+bIPHUKtmefBRo1avDOeQ1J1PlktNaDwr0GQsJF25wcLP7hB/QeNw62/v1FXLZsAf7+d/HBWLdMZuc8hyN8C0YUigwh8YhhGDh78CCuzcnBAF8O3hEjgB//GJg5s96d8xqaqBMZpdRoyLapFMBArfXKMC+JkKBi+mD6tGiBYZMmoZEvB+/8+YBhAAUFwPbtskVyOMIuMED0icxuANBa5wOAUqpEKbVFaz3G+0KXQ3gmAPTo0SOkiySkobA6eZeeOoVG/fr5vrBPHxGXpUtDu8AAiCrHr9Y63xQY1+vdAAYrpVJ8XLtGaz1Yaz24Q4cOIV0nIX5xOoGsLGD2bPnudPq9tEoUqWvXiHbw+iOqRMYP+QBGh3sRhNRILedZf/HFF561SA08GiVURM20Ape1sktrneR1fBeA1dWFsjmtgIQdp1MExeq0BUQgpk718J+Ul5ejWbNm0FqjoqICzZo1c1+fnQ1s2lR1NMrVV8v5wsKw1S35m1YQbZbMCh/HUgBsDfVCCKkVAc6zNgwDGRkZ+Prrr6GU8hQYAEhPF0EqKADWrpXvV18N7NgRsIUUaqLG8au1zldKlVqPKaXuArDR6qchJCIpLARuusn3OZfT1uqDSUpK8n0t4DkaxZeFFEF1S0AUiQwgzlxLGUGi69isMC6JkMAw51kPHVr13IEDMJo1q1s/mEAsJIpM7WAZAYlKxo8Xi8NHVu6+11+Ho2lT2Hv3rn3DqQAspHATdSJDSFRSzTzr1DvuwK0JCbjnnntq39GuBgspEsLaFBlCQkV6uvhIcnOB7dvxScuW6J2ZCduQIbi/rs+sxkKKhLolgCJDSGhxOW1NJ++NeXmYPWSIOHBzc2sfgq7GQoqEuiWAIkNIyLFGkX76059KqDknRxy4N90k25+FC0U80tNrfqCXhRRJdUsARYaQkFKlVCA/XwSmviFoa1g7woi2ZDxCopazZ88iOzvbM0wdYJJeNENLhpCGIACfSpMmTbBs2TK0atXKHUWKghB0faElQ0h9qaHw0TAMrF69GpWVlejYsaNnmNoMQfsiQkLQ9YWWDCH1wems1qdi2Gxw/OEPsNvtqKioQPPmzT3vj4IQdH2hyBBSH6rxqRjXXQfHhg2wp6YiMzOzqsAAURGCri8UGUJ8EWjeih+fimEYcBgG7E2a1FwqEOEh6PpCkSEE8BSVwkLZrsyYUXPeip+0/iZNmuCSdu2wuGvXwEoFIjgEXV/o+CXE6ri99lrg0CGgWTPg7FkRj2nTxOeyaVPVdple3eqOHz8OABjUqxeeKiyEbeLEEP8wkQctGRLfeDtus7KARx6pmgznr3WCxadiXHcdfrVjBxZefTWG/PWvUDHiU6kvtGRIfOPtuC0slJGvvpLh/I18TU+Hce+9cGzfjl4JCbi0rEx8KoGUBMQBtGRIfOPtuLX6WLyT4cy8FS+nsNG3r4SpXVGkWrdriHFoyZD4xjsZzupjsSbDmXkrgEfi3dft2sHx2muwt2lDgfEDLRkS33gnw6WlAVddJb6Y0lLghhuA5cuBf/zD3bDbknjXbehQzO7SBcNeekmKHemDqULUjESpDxyJQqrFOmbkP/8B3n0XOHECUAo4dw44dUr8NKdPA5WVQP/++Fe3bmg/dixSxo6VZ6xfL5MDInCCY6jwNxKFlgwh6elA27birP3sM9kitW8PPPSQiM7mzeKHSUgAbrwRxrFjcJw5gyv/+18sPXJE7o+RYsZgQJ8MIdnZ8tW9u4StbTbZDhUUAP/+N/Dii4DdDjRtCiMlBY6kJNi1xiMPPeTOnYmRYsZgwO0SiV+cThGQd9+VrVJ+PpCYKAl5t9wivpn77pMtU2UljI8+gqNtW+kHc/nlsH39tTiA9+8H9u6NqVKAusDtEiFWzJaXzZuLkNjtwFtvyZZo+nSxZHr0AHbvBr78EhrA1jZtYC8rQ+bmzbB98olcW1kpUafHH/ccuFaXfr0xCrdLJP6wZvkmJADffAP885/ukoK335brzp0DPv4YlZ06QR0+jPmffYaszz+H7cYbge+/B/LygHfeAX70I+Drr4HZs4E77gAeeCBiR8aGA4oMiT/MLN/Nm4E9e4BGjUQQrr9eHMC/+hXw2GPAoUMw2rfHgvJynKqoQNPUVLTcuFG2Vc2bAxs2iC/mv/8FGjcGiosl1H3kiPhzaqp7ihMoMiT+KCwUCyYnRyyMr74CevcWQVi3TgTkuedgNGoER2oq9NmzaNS4sYjHa68BX34p4ewVK8SH07EjsG0boDUwZw4wbx7wu99JHRQQU/166wJ9MiT+6NzZPYJk+HDPplEHDwJ9+sD49ls4+veHvVkzZB46BNuIESJMBw4APXsC330HXH450L8/8P77EoVavBgYN87d9uE3v3H7Y+I4xE1LhsQf48eLP8UMOaenS2To00+B1auRB8Bx2WWwN2qEzAULpF3DwIHAZZdJUl7btpIRXFkJfP45MHiwWCvWEoWBAyUkblovcRzipiVD4o+0NODKK2VbM2mSOHsPHgT27QPatkWPAQNwzbff4udffgnbokWSN/PPfwKZmeLknTdP/C+tWkno+rnn5LnWEoUDByRiVVQUU/166wJFhsQny5bJNunll4GLLwaSkvBFx45IadcOF+3YgV9u2ADs3Am89JLUK91xh+TUAMBttwFdu8r2SCng5Ek5bvaWmTQJOHYMuOYa4OhR2YrFcW8ZigyJX5o3By66CLj3XhhbtsBx8cX4SUUF7lJKzqenA7t2if/lz38WJ6/NJtGkpk3Fotm6VXJkANkiNWkCnDkjIvT++8CttwIZGXErMABFhsQStUmCy80FfvlL4IorYGRnw3HiBOxNm+LmJUvEz2J2wPv0U3H49usH3H+/OIwNAygvByZOlOPXXy9O3u7dZYt0880Syn788eobV8VJ0h5FhsQGvobW//znQLt28gG2foidTkmi694dxu7dcBQUwN6qFTLHjIHNdOhu3y7XHTki257SUuD3v5fnZ2SIFbN4MdCihYjMmDGyDq0lia+mEgNf6/XXrDzKociQ6MfpBF59VcTgn/+URDhABKCgQPrA7Nwp/pfkZKB1a6BHD3x/ySX47aefwn7+PDJTU2HLyRHfixkJys0F7r4bWLVKokrvvedu03n55SIS3btLpvD589LsKhCRqGEg3IW+wjECRYZEHuY2Yt8+ifycPi2CMXSob//GkiXiNzFDzFu3Ar/9rWTt7t0rVsOSJcDYsWKBXHYZMGAAWu3YgceXLkWPrl1hmzlTxMh02k6eDKxZA5SVyeSCkhIRmT59xD/z7rvig5k1C/jiC+kjE6hIVDMQzmez8iiHeTIksjDHkxw+LMWJR4+K/6NLF8mq/clPPOuAnE5JhPvVryRjd+hQ2a5kZgJvvin5MMOHy7miImDKFBgDB+JPO3YAV12FSzIzYfvLX+S+Tz4Rp22jRsALL4hjeMIE4NJLReRycoC1a0XERo8GbrzRs7l4oJm9ZrNyX/hrVh7FUGRI5GBuIx59VHwUKSkiIDt2yAyk3/1OokGvvOKuA8rNlWQ46wezsBAYNEgskOJiuTYrC9i3D0b37nB89RW2XnQRzmstvpOCArnmq6+kUVWXLvK+b7whYjd0qFgsFRVSBNmzp/hqxo+vmmQXiEh49xW2EoNJexQZEjmY24iPPhIxmTZNrAPTQvjoI+Dee+WcaS0UFoofxDJgDZ07y1bH6ZRq6YcfBnr1gvHFFzI61mbD4z/7GRKOHpVtydKl4ou59VbpDWNuZcx+vxs2uLdkDz8s6ygvB379axGvyy5z/wyBiITXQLgLmEl748c33O80AqDIkMjB3EYUFspr65bCtBD69JHXprXQubM4Xc36o/XrJYdl40Y53ru3bJFSU+Ho0QP2kyeROXIkbIcPV51EMH6851YmO1usqFmzREicTuBvfxOB+fZb8RmNGOHurBeoSFgGwmH9endz8hhN2qPjl0QO5jaic2dxvlpnTJsWwoED8toUCDOVf/1699D6d9+VPi6HDokl86c/4chf/ypRJAC2p5+WkoD58+W+1193f7jNNZg+GDMCtGgR8OCDUoXduLHkw5w9KyLTsaPk3LzyirSMCEQk0tPd692+XX6eGO2sR0uGRA7mNmLkSBGZ9evFOjAthJEjxfF67JjbWrBaBXv2SASpdWvJzH3kEfzQrx8wYwbu/PxzrJg4EbaRIyVa1auXOJILCjynPZpr2LDBMwJUUgJ88IFEp1askNyYVavk/s2bxWK6/PK65bjEeAtcWjIkcjAFY/ly2bJ8/DFwxRXiwL3kEuk4p7U4Z63/4ntbBVoDt94KY+RIPP/f/yLrpZfQyzDQxDAkWjVhAvD889Wv4Ykn5Lk7doj19PrrQLduElXSWt4nLc29jh07RAADhcl4hIQJq2AMGiSWwunT4oMZPtx/HZD1A+90wli4EI6dO2G323HRsGESbi4pEYsnI6PmNezaJT6YY8fcW5ncXHdUyNu5W5uoEJPxCAkzVsGoA0ZZGRyNG8N+7hwyL70Utn373NZIII5VMzyelydRLWtN0SOPSNnBb37jvr62rRziLBmPIkMij3oUDn722WdwOBywp6Yic9Ik2D74wN2gyhSBrCz/z7ZuY5KSpMbpj3+U8LYZ+VJK/D+nT9dOvEwKC2WL5IsY7KBHkSGRRT19Fb1798aECRMwYcIE2Gw22XoE+mxf25hJk8QJvH69CM2GDXK8PlEhM4JlRs6sxGAyHrTWMf81aNAgTaKAvXu1HjdO6+PHPY8fP671zTfLeT/s3r1bnzx5sn7PzszUet063/evWyfnG4J6/JyRDIA87ePzR0uGhJbqtkJ19FUYhgGHw4Hrr78ec+fO9f2+3s+2rqNFC+l6V1kZmm2MNew+ebI8uy7briiBeTIkdGRlyYfrr3+VbNzGjT0Hn5mjSrKypEYoK8vthPVTE2QKjN1ux/Tp0/2/t3cm78KF7gFsaWmSwFdYGLqaIrN5eUGBhL6983ViCFoyJDRkZUmBY2amdPvPyQE+/FDS9V95RXwnhYXAypUSwfH2mTRpUuVDbhWYzEmTYHv2Wf8OXX+ZvIAcnzxZkvNeeklCyVZrKliNwOsZRYsWlI7xbEMAGDx4sM7Lywv3MuIXp1OE4uGHJbHOdL526iQZtGZP3LNnRQT+8AfPTNu775Z6oVWrLnwoz58/jwcffBAtWrSQhlNvvy3P7NtXRMO7gZTTKYJ1+eXSumHaNPfzp04VAdmzR9ZWUeF7GxODVkZDopTapbUe7H2clgwJPrm50kGubVtpYWm1IubNA7ZskfnT8+eL8Hj7KoqLZWtj+Vc/ISEBWVlZaJ6fD9sTT9Sc2FZdJq/pBzl9Wto8zJ4dFzVFoaJGkVFKXay1PhSCtZBYpbBQCgpNC8a6FenTR1pnJidLQeOiRZ4lAkoBqanAZ58BWVkw+vbFjqIizJkzB+3atZPygECdxf4yec3zpt8lTrYxoSIQS8YBYFKwF0JimM6dxcn74YdSwGhNhquokGZRAwfKgDXA/SE381q6dweGDYPRtCkcr70Ge3IyysvL0bJly9ontmVkyLbpySdD43chAUWXJiqlBgR9JSR2GT9enKpt28r26JtvxFrJz5e2mZdeKiLz1VfuRk5mYtz//i9w+DCMUaPg2LFDnLzFxWhpClJtu8zFWS+XSKBGx69SaguAPABvaK33hGRV1a9nJgCzpViK1nplTffQ8RsBPPCAWCZNm4p/xWzOfeedkvZfXi5jRQxDtjl//7v02j18GMawYXDk5YnAZGbC9tZbEvJdutTt0LX6ZABPh64v4TDzZIqKRIhidOZRKKmV41cp1UZr/R0AaK3HuI7dqZRK0Vq/Fdyl+scUGK31JtfrFKXUaq31rHCtiQRAdrY4d2+9VXwh+/fL9un22yWE/cEHIjymaOTmiuUzbBjgcKCVUrjshx/w6KOPSqmAdRtU18Q2+l1CRhWRUUqNAqABfGQ9rrXOUUpdqZSap7V+JlQL9GKW1nqQZU35SqnRYVoLCQRz2zNsmERthg6VYxs2yFeTJhJV2rZNrjc//EeP4lhhIS5atQr9O3fGsrvvhrLZ5BrvbVAcdZmLRjxERik1HVJ/4HNbpLX+t1IqXyn1NIDlprUTCpRSiQAG+jhVqpQarbXeGqq1kFpgpvMXFLiLAtPSpLvcpZfK8cpKT9HIzobx2WdYabPhF2lpGN66NdSjj4rF0qePOI6vvFK+m9uc2lgmcTIeNlLwdvwWAEip7gat9Umt9UIAM305hJVSFzfY6jxJAVDq43gJalgzCSNmOr+vDv19+kjY2tp82+mEkZMDR9u2SElOxpXvvivH586V2dK33y59defNky2WtSwhELxLCuryDFIrPERGa/1nAP9WSk2o7ial1EgAQwHsVkp9qJSaZxGcYPlHkuF2+FopBZDoY40zlVJ5Sqm84uLiIC2J1IgZ/fEV1XnlFRlMb/pOnE4Yv/gFHI0awd68OTLnzIHNrO95+WV34+4ePeTeggKZ0bRpk7vGqTqsrRzMQXDTpsnrQJ9Bak2VELbWusDq3FVKrXB9H6mUekEpdRzAVgCjAbwE4CSADIg4nQcwPyQrrwGt9Rqt9WCt9eAOHTqEeznxi9WCsRYFrl4t415XrJDj2dkoWrgQK5s3hz05GZljxkgm786d8oziYpno+Je/iNPYtEJWrJAs3ZqmNgKBVXmTBieQZLxZrqhOEsRq2AjgTZfVcwGlVFsAYwA83eCrdJPs41gVK4ZEEL6iP40bS4i6WzeJNAFATg46rV+PBzMzMWTAANjS02Xg2ogR0qGuvFx8N9deK07iTp1EnG65RaJWp07VvJY460gXKQSSjJcI4M8Axmitk7XWGd4CA1zw1WwCsKmhF+kiD74FJRnA7iC9J2kIrBbMwoVSRnDffXKsVy/seOopfN6+PZCcjOunT4ctJ0csnyNHpPXDNddIfg0gkwys25vkZCl6DGRLHGfjYSOFQCyZNVrrGtq7e7C6roupDq11qSuylai1tjqAExlZigKso0MM48KWxaishKNbN6QdPYosp9PT8mnRQuZc79wJfPGFzFLq2tVze9O1K7B3r5Qe1IQ5CC5UrRwIgMBrlwJGa11Qx7UEggPATAArAUApNRDiHyLRgJdP5EI/GJsNCzp1chczmnkv6eky8L51a5mVVF7u3nZVVEjzK8OQUHiPHjW/f5x1pIsUahSZIItGrdBar3FFjUZDtk4pzPaNIiw+kSoNpx56CDh3zjOHpXlzyRA2w8sLF0qV9kcficBoLa+XL5cQdyAwcS/kRF0/Ga31mnCvgdQRS5f+bdu2uWuRbDaZEPnOO8Du3eLw7dhRok9FRRJRevBB93TJ224DfvhBHL7Ll9feCmFJQUiJOpEhUcz48Ti/cCESbrkFc+fORUVFhbRrKCkR8enWTYokzV4vOTkyqnbZMuD774FRo0R0fv1rufbcOVohUQBFhtSPQFL0XdcY+fl4PSEBT9xxB5JmzEBjq0+kXTtg5kx3W0yTtDQJT2/dCnz9tYjPW29RWKIIigypO/6GpV11lZwvLJSvkhIY48bBUVoKe3Iymp44Ifd16eL2iaxa5Z4m4M2oUSIwq1aF7mcjDQZFJhYJRQGgv6HxpaVSYzRwoDhm9++H0bcvHNu3y+jYzEzYzCjR7NnudcXbVMU4gnOXYo1gFgA6ne6ZSHPnSvsGa76J0ykjRbp0kS54djt233QTHG3awH7mjESRbDbfafy+CigBdw6LWUBJog6KTCwRzAJAb/FSCnjvPbd4OZ3A/fdLXkvjxsAnnwD//jd6XXstrhs2DJlHj8L26qvu53kPa2NbzJiF26VYoo5jXmskJwd4+mmxXAoKgCuukJqiDh2kGvq116QJ+LffSuJccjI+P3MGvSsrkfTcc3j4wQclLP3+++6tnJnnYmb5AsxhiVFoycQS1lGs3vgZ81oj2dnSTmHkSNkmmdsvAHjhBXnPTz8FFi8W0WjeHEZZGR7t2xcbExLkPZ98Uq4pLpZ7O3SQPJfRo6tu5dLSpA3nqlXynQIT9dCSiSUaynnqdMoA+o8/liLFzp2B666T55qD0378Y5n8+MMPYinl5QG9esE4fRqODh1gP3sWP77iCqk/OnhQrKFWrSSB7p13gOnTxXIxG36bQ9hIzEFLJpZoCOdpdjbwk59IO4WLLgJuuAE4cQJYsMBtcSQny7m0NEmKa9kS2LZNokj9+8N+6hQyExJgu+oq8c8kJsrkAaWAffs8B8uzl0vMQ0smlqhvAaDTKT6WTp1ErBYvFidvRYX4VNascVscWkunupIS4L77cLq8HM/t3w+7Usj88kvYGjWSrVJZmWTmdusG3HuviI33OtjLJaahJRNrWHu3rF0r362WQ3Xk5oqFMm2aWBjm9mv4cBGZQ4dklvSOHTKIbedOuXbbNrQcMwaZP/yAzKQk2IYOBc6fBz7/XIocW7UCxo2T9/DlF2IeTExDSybWqCkRr7rzhYXy3XQeW/uvzJghIvPVV9I68+hRoGtXGNdcg2Pl5bh9+XL07tJFfC9lZSIaI0fKPVpLweNbb1XtTMdeLjEPLZlYorpEPKdTrJFbb5UclmuvrZqo17mzfDe7x1m3X6+8ItZJ9+5ijTzxBIwJE+D4n//B344cwbnJk0VQyspkltLAgcCgQRLe/vhj4G9/k3PbtjEPJs6gJRMr+Evzv+UWyWlp2lRE4v/+T0Ti978XATE/6EOGiOWyZYscM7vHpaeLz+SnP5VI0oQJQEYGjLIyON55B/YePZDZrh0ab98uPXifeEKsmE2bRMROn5akveJimRa5bBnzYOIMikys4C8Rz+yT27KlVDmPGSPHb7nFLS5mdGfpUrGAfvMbiSqZ0ajcXPGrPPYYkJ7u2XDK7AfjdIog5eVJst7990uTcFNMzCmR7OUSd1BkYgV/nfhzc2Wb9OabIjZZWW5/zPDhcn7sWHd0x8y6ffFFydAF5LqMjAviUFJSgt69e4ez6hoAAA8dSURBVOOJJ54QgTGrsXv0kGt69ZJjd94pwgWIdUTnblxCkYkV/CXiFRaK07WyEli5EnjkEXdbhhdflG1Ur16eApCWBjz/fJW3KCsrg81mw7hx4zB27FgkJCR4btOOHAEeeECS9C66SGYitW0rVhGdu3ELHb+xgr9EvMREGWzfsqVk3+7f756+eMMNUtC4Zk2NiXqGYWD69Ok4cOAAAIjAAJ7btJ07pdXDBx8APXsC/fuL6Fx/PZ27cQwtmVjBXyLeP/4BnDwpjaQ++UQEYPx4aWf54osiPu3aVSsAhmHA8fTTsLdoga7PPitjSMzQt7lNMy0asxQhN1cS77p0kYfs2iUixAH3cQdFJpbwVcW8apUcf/ddabrdpo3kq+Tlia/lyBGgWTO/jzQMA44VK2CvrJTRsf37uzvg3Xmne5tWUOC2aJKTRUTWrxfhOXRIIloLFnjeG0iCIIl6KDKxhq/oTevWUvFcXCwWx4ABEmpOSxOHrmsL5M2BAwfEgqmsROYLL8DWrZucMEPjU6dKFCk7W3wwVsdzSYlsw5QSX9DmzZ4FliyKjBsoMvFA+/ZisaxeXXVyYjXTF+12OyZ16oTbhg51C4xJcrJYQs8/Lxm9b78tW7B58zybg995p2QHWx3L9e1vQ6IKOn7jgf79Za6Rr65zl14qSXIWdu3ahePHjyMhIQH3KCVbJG+ys4GNG2UrlJwMpKbKluyRR9z1Up07i7j4qgCva38bEnVQZOKB8ePFQbtokWfh5KJFwDffeAiAYRh48skn8arZKtPXkHqnE/jtb6W6+sorxdfywANASoo4lxs3lkzfY8fEsvEVWWJRZNzA7VKsUF3hoxl5Wr5ctiljx8qH3Gv6ojWTNyMjQ+71NaT+xRcl76ZTJ7dvx/S1XH21+F/MGUnFxZIMaIVFkXEFRSYW8Df/yBrBqaF/rs9SAcB3aPyDDyR7eMoUTwslOVmOvf++e0ZSdjYH3Mc5FJlop7rCSO8Ijp+6ofPnz2Pjxo1VBcbEW6DKy4Hbbw8sBM3m4HEPRSbaqW5CwfDhMh+pb99qk+ASEhLw+OOPo3HjxlUFxsQqUEpJsWNJiTvxrrBQsotzctxFmL7uJXEHHb/Rjr8JBdnZ0mJBKb9D3swt0tmzZ9G2bVv/AuNNRob4ZIYPl854lZVSH/X22xIx6tixgX44EgtQZKIdf9GfnBzZLo0Y4XPImykwxcXFOHv2bO3eMy1NfCpFRZJj43TKV3KyjD8xjPoNkiMxBbdL0Y6v6E9urnv0iDWC40qCM7Kz4Th06IIPpmXLlnV773nzJIxdVCS+FnM7lpTERDtyAYpMtJOWJsWPV18tVoXdLrkqCQky28jrg74DgOPgQdj79PHt5A2UwkLZhvma8cTpA8QCRSZU1NTgu673ZWdLy8wePSSxrqBAHLKjRvmM/iSXlODyNm0wvz4CAzTcIDkS89AnEwqqa/Bdn/vMzNsWLWSm0auvApmZUkLwpz9JmwcX33zzDVBSgtQPP0TmggX1ExigYQbJkfhAax3zX4MGDdJhY+9erceN0/r4cc/jx49rffPNcr6u9/3851oPGOD7mm7dtO7RQ+t16/T2tWv1+Ntv11vHjdP65Zcb7md7+WVZy7p1WhuGfL/55oZ9DxI1AMjTPj5/3C4Fm+ryWKqrRA7kvn/9S1L2fV1z//3Apk0wdu8WH4zNhqt/8QtJjGsomGhHAoAiE2z8NfgGqneQ1vU+C0ZSkkSR6uvkrQ4m2pEaoE8m2PjKYzHxdpA6nTJNYPZsuWfr1urvGzpUihF9+EWK334bK5OS/JcKEBIiaMkEG195LEDVSmTvIsetW8Wp26oV8OCD/u/7+9+Bu+8Wx69ZgLh2LTqcO4dfTpmCy2+7jQJDwgpFJtj4a/BtrUT2V+TYurVMXPz+ewlJ+6pgfughmT6weTMAwKisRIvz5zHgoYfwo3vuCdMPTYgbJU7h2Gbw4ME6Ly8vvIsw8128s2MBafhUVCR9cr1zYZYvF6umX7+q93k92ygogKOkBJd27Ignu3aFKiridAASMpRSu7TWg72P05IJFf4cpNnZ0rbyvvuAceOq9oIZNUoaQJn9Wfw82ygrk34wiYlYdOQI1LXXAjffzOkAJOxQZMKJuU2aMkXKAXx18w8ge/ZCw6kuXZB56BBsgfSWISREMLoUTsxcmClTPLNnzVyY3/8+oOzZvLw8iSK1awdbTbk1hIQYikw4MXvBWJ3D5jSBgwdrbFN5bs8eICsLs/ftwzKbDbbPP/fdWwbgdAASNigy4cSaQ5OeLmFpc5qA0yntGvz4UYynnsL/W7IE33bsiIT0dLTo3RvYs0cGqvmCRYskTNAnE068c2hM53BJiVg15sQAL4w33oDDMGC329Fi6lTAZhPfS69ewIQJUhg5fLj7Bk4HIGGEIhNOAsmh8cIwDDg2bIC9XTtkPvWUZ6Ld8OHuBt9LlnA6AIkIKDLhphZFhnv27JEoUuPGyExP953JO2MGcOqUbLtYtEgiAIpMJBBgkaHdbseoUaNw37FjsH35pe+LDhyQsbNLlzbwIgmpG3T8RgGffvopKioq0Lp1a8yZMwe2iRPZMIpEDRSZCMcwDCxevBgbNmxwH/QV8l6/Xl7T90IiDG6XIhjr6NiJEyd6nmTDKBIlUGQiFL+zqa2wYRSJArhdikDOnDmDF154gQ2nSEwQNZaMUmo0gFkAVgAoBXAXgFKttZ8U1+ilefPmyMrKQvv27SkwJOqJJksmEUAKgF2ur3axJjCGYeCNN94AAPTs2ZMCQ2KCaBIZaK0Haa2V1jpJa70g3OtpSEwfzM6dO2s/m5qQCCZqtksRTV2nQ7rwdvI2adIkiIslJLRElSWjlBqtlLrL9X1+DdfOVErlKaXyiouLg7eouk6HdBFQFImQKCaaLJndAKC1zgcApVSJUmqL1nqMr4td/po1gPT4DcqK/DUAr0UnutOnTyM1NRWPP/44BYbEJFEjMqa4WF7vVkoNVkqleJ8LGXWdDgngu+++Q5s2bTBy5EiMGDECCQkJIVgwIaEn5Nsl11ZnS4BfKTU8Lh/A6FCs2ydmZztfVNOJzjAMzJgxA/v27QMACgyJaUJuyWittwLwMxrRNy6x2aW1TgrOquqI2dlu6NCq5/x0orP6YHr16hWCRRISXqLJ8bvCx7EU1FKwGpTx42tVDU0nL4lHosIno7XOV0qVWo8ppe4CsDFs/higVp3tDh48SIEhcUlUTZBUSs10/WciAGitVwZyX9AnSFY3HdJFZWUlNm/ejLFjx1JgSEzib4JkVIlMXQnnmNq8vDx069YNnTgpgMQ4/kQmmnwyUYdhGHjqqaewdu3acC+FkLBBkQkSVifvnDlzwr0cQsIGRSYIMIpEiBuKTANjOngpMIQIURHCjha01mjUqBEee+wxKKUoMISAlkyDYRgGli1bhvLycrRq1YoCQ4gLikwDYPpgvvvuO5w7dy7cyyEkoqDI1BM6eQmpHopMPfjXv/5FgSGkBuj4rQcdO3bEoEGD8PDDD1NgCPEDLZk6cPjwYWit0bNnTzz22GMUGEKqgSJTSwzDwNy5c/HHP/4x3EshJCqgyNQCq5P3uuuuC/dyCIkKKDIBwigSIXWDIhMAJ06cwDPPPEOBIaQOMLoUAElJSXj00UfRr18/CgwhtYQiUw2GYUAphauuugqDBg0K93IIiUq4XfKD6YPJzc1FPHQPJCRYUGR8YHXyLlmyBEqpcC+JkKiFIuMFo0iENCwUGS/2799PgSGkAaHj18XZs2fRpEkTpKeno6KiAs2aNQv3kgiJCWjJQLZIs2fPRmFhIZRSFBhCGpC4FxnTB5OYmIg2bdqEezmExBxxLTJ08hISfOJWZPbu3UuBISQExK3jt3fv3rjpppswdepUCgwhQSTuLBmn04kzZ86gZcuWmDVrFgWGkCATVyJjGAaWLl2K9evXh3sphMQNcSMyVifvlClTwr0cQuIGFQ/Ff5dcconu168fnbyEBBGl1C6t9WDv43FhyRQXF1NgCAkTcWHJKKWKAXwZ7nXUQHsA34Z7EREEfx+eRMPvo6fWuoP3wbgQmWhAKZXny9SMV/j78CSafx9xsV0ihIQPigwhJKhQZCKHNeFeQITB34cnUfv7oE+GEBJUaMkQQoIKRYYQElTitgqbRCZKqZkASlwvU7TWK8O5nnCilBoNYBaAFQBKAdwFoFRrHVX+GVoyEYhSarRS6k2l1EClVIpSar7rwxfTmAKjtd6ktd4EYJNSanW41xVGEgGkANjl+moXbQID0JKJVKx/XKUA1sTJv+iztNYXRnVqrfNd/5rHLdbfR7RCkYlQYuGPqzYopRIBDPRxqlQpNVprvTXUayINA0WGRAopEKvNmxLXubjEZcklQn43A6PRoqXIRCix8MdVS5LhdvhaKYX8HuKR3YBsGwFAKVWilNqitR4T3mXVDopMZBITf1ykfpj//y2vdyulBiulUrzPRTKMLkUgWut86x+R1no3gMFKqVjfNiT7OBavVow/8gFElTOclkwQcW15FgR4+awa/nUy/7iiLoQZIHnwLSjJcFl28YTrH5RdWuukcK+lvlBkgogrIlKrqEgs/XHVBq11qVIqXymVqLW2OoAT4ziytMLHsRTU8m8q3HC7FJnExB9XHXAAuJB0qJQaiNj/mX3ismo9om1KqbsAbIwmfwzAKuyIRCk105rZ6frjGqO1nhXGZYUEV9ZvPlwJiXEQVasWS6Z3IgBE4++DIhOhxMIfFyEARYYQEmTokyGEBBWKDCEkqFBkCCFBhSJDCAkqFBkSdlyNuVYrpbRS6oRSyuHnui2ua96MgxKLmIHRJRIxKKXeBDDaX7azS3zecNVykSiBlgyJJFYDSPTVDc+VkEiBiUIoMiRicNUolQKYaD3uKi8ABSY6ociQSGMjgLvNF662nKNdjcVJFEKRIZHGm/DcMj3Kkoroho5fEnEopU5Aqq93QiY1+Or9S6IEWjIkEtkIGWS2lQIT/VBkSCSSDBEYOnpjAIoMiUQGAtgS7kWQhoEiQyIKVzQpHroAxg0UGRJpjAaYExNLUGRIROCqX3JA+vxCKeVgfVJswBA2ISSo0JIhhAQVigwhJKhQZAghQYUiQwgJKhQZQkhQocgQQoIKRYYQElQoMoSQoEKRIYQEFYoMISSo/H+CX4qm6X2D2gAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARkAAAElCAYAAADQqJbBAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO2de3hU1bnG3wWEAEMx3ISAXJIh3EMpEJSKqIgVUatcFMGkFjxyOdpSURFQ1MRHcdBae06rJX2MQqhFBKFaq6eirfXCRhIKBFRQCCAaIhgCQgoJYZ0/vtnMzmSSzCSz5/r+nocnZs++rImZN996v299S2mtQQghdtEs3AMghMQ2FBlCiK1QZAghtkKRIYTYCkWGEGIrFBlCiK1QZIjtKKXeUUolhXscJDxQZIitKKWmABgHoEOIn5uklJqllHo1lM8ltaHIELuZ6v6aGqoHKqWGQYStLJTPJb6hyBDbUEotAPCA+9uQTZe01lu11msB7AvVM0ndUGSILSilUgGUa63NDzojijiFIkPsYrbWOtf93+UAnOEcDAkfLcI9ABJ7KKXGAdhiObQPfhi/Sqnl/j5Daz27EUMjYYAiQ+zgZi8R8MuApXDEJhQZElSUUrMAdFBKuSyH6cfEMRQZEjTcBXdOrfXNXse/A+DyfRWJdSgyJJjM0lo/4ON4uT8XB1A4V8apVfRAkSFBwV3ZW1ddyj73OamWlHYtvCOgIBDSKmPiG6awSZNQSqW6I5BXAdSKLtwejRndLHdnnkIxJrMQcJhSyuUeBwkDij1+CSF2wkiGEGIrFBlCiK1QZAghtkKRIYTYSlyksDt16qR79+4d7mEQEpOcOnUK5eXlOHz48FGtdWfv1+NCZHr37o2CgoJwD4OQmMMwDLhcLjidTvz6178+4OscTpcIIY3CKjDZ2dl1nkeRIYQEzPbt22sIjMPhqPPcuJguEUKCS1paGq655hpkZWXVKzAAIxlCSADs3LkTp0+fRps2bTBnzpwGBQagyBBC/MQwDCxZsgSrVq0K6DqKDCGkQawm77Rp0wK6liJDCKkX7yySP1MkKxQZQkidVFZWIjc3t9ECAzC7RAiph5YtW+Kxxx5DUlJSowQGYCRDCPGBYRhYuXIltNbo3r17owUGYCRDCPHC6sFUVlYiMTGxSfdjJEMIOY+3ydtUgQEoMoQQN03NItUFp0uEEABAdXU1+vXrhyVLlgRNYACKDCFxT3l5OZKSknDppZdi1KhRaNYsuBMcTpcIiWMMw8Cdd96Jbdu2AUDQBQagyBASt5geTK9evZCWlmbbcygyhMQhdpm8vqDIEBJn7N+/P2QCA9D4JSTu6NWrF+644w5ceeWVtgsMwEiGkLhhy5YtOHToEJRSuP7660MiMAAjGULiAtODycjIwOLFixt3k6IiYP16oKQESE4GJk4E0tMbvIyRDCExjtXknTdvXuNukpcHLFwIpKQAM2bI14UL5XgDMJIhJIYJShapqAhYtw7Izwc6dJBjI0cC110HZGYCGRn1RjSMZAiJUbTWePPNN5ueRVq/Hrj1Vo/AmHToAEybJq/XAyMZQmIQrTWUUli8eDHOnTvXNJO3pAQYP973a337Aps21Xs5IxlCYgxzV4GKigq0bt266Vmk5GRg927fr+3ZA3TtWu/lFBlCYgjTgzl9+jS01sG56cSJwOrVQFlZzeNlZcCf/yyv1wOnS4TECLYtFUhPByZPFpN32jSZIu3ZIwIzZUqDaWyKDCExwJYtW+xdKjBzpmSR1q8XD6ZrV8Dl8qtOhiJDSAzQrVs3ZGRkYN68efZV8qan+yUq3tCTISSKKS4uPr+jwOLFi0O2VCAQKDKERCmGYWD+/Pl4/fXXwz2UeqHIEBKFWE3ecePGhXs49RJVIqOUGqeUelUpNUwplaqUWqCUmhXucRESSkLZcCoYRJXIAEgCkAqg0P2vo9Y6N7xDIiR0HD9+HM8880zUCAwQhdklrfXwcI+BkHBxwQUX4MEHH0SfPn2iQmCAKBQZQuIRwzBQVVWFyy67DD/84Q/DPZyAiDqRUUqNg0ybygEM01ovC/OQCLEV04Pp27cvLr30Ulu2LbGTaBOZrQCgtd4HAEqpMqXUO1rrq71PdBvCswCgZ8+eIR0kIcHCavI+/PDDUScwQJQZv1rrfabAuL/fCmCEUirVx7m5WusRWusRnTt3Duk4CQkG0ZZFqouoEpk62AcgsgsFCGkEX375ZdQLDBBF0yV3tFKotW4f7rEQYidnzpxBYmIibrvtNtx8881ITEwM95CaRLRFMkt9HEsFsDHUAyHEDgzDwJw5c85vXRLtAgNEkci4vZhy6zGl1BQAa6w+DSHRiunBdOzYEe3bx07AHjXTJUDMXMsygiT3sdlhHBIhQSFWTF5fRJXIACI04R4DIcFk586dMSswQBSKDCGxRlpaGq6//nrceuutMScwQBR5MoTEGtu3b8epU6eQmJiIO+64IyYFBqDIEBIWDMPAo48+ipUrV4Z7KLZDkSEkxFhN3p/97GfhHo7tUGQICSGxnEWqCxq/hNhFUZFsIVJSAiQno+qGG5CXlxdXAgNQZAixh7w8YN062ah+/Hhg924kPPQQHhs/Hm1nzowbgQE4XSIk+BQVicDk5wNZWTDOncPyigqcW7ECXd56C4598VWgzkiGkGBSVATMmwcoBfzudzD69YNr9Wo4nU5Utm2LVqNHy+v9+slG9hMnNmrDtGiCkQwhwSIvD1i4UAQmMxNGy5ZwrVoFZ7t2yM7ORquXXwbefFNenzEDSEmR8/Pywj1yW6HIEBIMrFOkyy+HUVoK1+bNYvIeOQLH22/L65mZwOWXAyNHAllZcv7atXJ9jEKRISQYrF8vJm+HDsDEiUj44AMM6NMH2Y8/Dse0acDvfw9cfz3wxhsyRTLp0AGYNk2uj1EoMoQEg5ISoF8/fPfdd0B6OoZPnozHd+2C47XXgMpKYM8eiVqmTKntwfTtCxw+HJ5xhwAavyQ+8KpZCbrhmpwM47338JTLhYULFyJj5kyojAx55vvvAwkJEsnMnFn72j17gK5dG/dcu99XEGAkQ2If05BNSbHNcDX69YPro4+Q0qMHBg4cKAfT04G77wYSE4GnnwY++ggoK6t5YVkZ8Oc/15xC+UsI3lcwYCRDYhurIduhgxwbORK47joxYTMymvyX3zAMSVN36IDsPXtkitS3r0Qof/6zTJEmTwaOH5dnTptW+/VAxxCC9xUsKDIktrEaslashmsTPoyHDh2quRZp3z6556ZNMgVyuTz3nzlTPvx1vR5B7yuYUGRIbFNSImX9vujbVz7sjcHthVxUUoK57drh0lOn4LjvvoZ9kfT04Hz47XpfNkBPhsQ2ycnA7t2+X2us4ZqXh08WLsS+9u2Bjh3xk+3b4di0CejYMXS+iB3vyyYoMiS2mTgRWL06eIZrURGMdeuwtGVLrDpwAPj3v4H33pN/W7cCQ4eGpsAu2O/LRigyJLZJTxfTNTNTPvx/+pN8P3Qo0LJlwLcz8vLgatECTqcT97Zp4/FFrF5IKArsvN/X5s3yNTOzcUayjVBkSOwzc6YYrOvWAYsXy7Enn5QPaQBTG8Mw4Nq7F85u3cTkPXJEFjqaWIvqQlFgZ76v4mJgxQr56nL5rsUJIzR+SfxQVSXTG2tGxs+Ur9YaGzduhNPhQHZGhvSDMX2RkSPlJKsXEipfJFhGso0wkiHxgT8p3zo4d+4clFJYsGABcubPh2PdOvE+rL6I1QuJQF8knDCSIfFBI1O+hmFgXX4+Hu7YET84cgQtk5OBiy/2FNX96EfAlVfKyTfcAGzbBjzwQMT5IuGEIkPiA++pjZU6pjaGYcC1dCmc586h2fDhsvZo926JXi65RDyQsjJg9Gi5oKzM44sAQE5ORK8pChUUGRIfTJwoJu9119WcMplTG1MY3JwXmIoKZLdpA8eqVUCnTsDgwcCiRcATT9Rdreujvy8WLhSjOcJM2VBAkSHxgTXl28DaoYKCAo/AVFXBcfSovF5UBLRrByxdCnTr5rt0P4rWFIUKigyJH/xcO9SzogI/rq7Gf7dpIwLzxhsiGGVlIhSLFwMLFgBffVV7ShRFa4pCBUWGxBf1pHy//PJLpKam4sJ//hP3T50KvPiinGsKhikU770nhXwFBSI61ilRRUXUrCkKFUxhEwLxYO6//3689tpr57vcoVMnmf5YS/f79pXlA9u3A5ddJn16rf16Dx0CNm70/ZAIW1MUKhjJkNgjwG5x1q1jr732WuDTTyUyGTxYPBirj/Pii8AHHwBDhgDDhtV+Xq9ewPLlwJw5fhnM8QAjGRJbBNgt7rzAdOsmWaT77gNKS4E//hEYO1aEY/FiSU0vXw68+y7wwx9K9fDEibWfN3o0cPIkcMUVEb+mKFQwkiGxQ4CZnZMnT+LZZ5+VfZH274eje3egulqmNbt3i1E8ZoyYvEOGAJ9/LtuZvP++fAVqP2/3buDOO6XVpmE0vTlVDMBIhsQOAS4daNu2LR7JzET20aNwjB8PvP028N13wAUXiDiVlgL79gE9eojw9OgB9OwpiytLSoCXX675PHNKdNttwKxZQJcuwHPPAQ8/HLcCAzCSIbGEn0sHDMPAiRMn8JOf/AQDtm4FLr1UpkddugA33SSm7+7dMkXatUuMXG+ROH4cePRRiXY2b65dc1NREZeZJF9QZEjsUN/SgXffBT7/HMbPfw7XsWNI69EDV111FZqXlABHjgDNmgGvvFJzmtWlCzB9OvCHP8jmbFZmzgQKC2WK9u23tadEcZpJ8gWnSyR2sK6KLiqSQrm5c4EJE4Cnn4YxbBhc5eVwtm+PR4qL0XzFChEmwwDatwcefFCuMTvalZaKaGzZ4vt5Y8fKuf/5T83jXIVdA0YyJHYwlw5cfjnQvLlMfQDgww9lX6QvvoAzLU0aTp05I2ZwSoqIydChkh0yi+quuUY8mv79JWLxxlyfNGaMpLwTEyVtPWAA8M03cZtJ8gVFhsQWGRkSlfz4xzKNcYvG1z17wvnKK8ieOlUaTjkckm5+6SVZUf2vf4mYDBsmq60feUSmRG+95VllbeKdxTLrZCorJfNkdt0jACgyJNIJdBvW9eslhZyVBQD4z+zZaD1uHCaPHImfVlcj4a23RIgAYP9+yRg9+6yIwtNPA717A06nREOvvgq0aSMRivczrFkl61KF/Hwxiyky56EnQyKXxmzDai4JKCqCcc89mHXoEIqffBIoKkJC//41++7u3SuCAgADBwKnTsn1+/YBR49Khuiee2qLmvkMX4Sit2+UwUiGRCbWKcnXX3uimSFDpLS/rpYJyclAbi6MkhLZVeDCC3Hhtm1St1JZKVmknBwxbQ8elJqYhQslMrnlFnlmQYGsW5o0yXf/l0Y0wIpnGMmQyMSckmzYUDOaGTgQ+P57YMkS39cNGgTjgw/gSkgQk/d3v4Nj6FD58H/6qYhTUZEYuxUVwP/9H3DXXTK9mj5dROa996SNw9ixvp8RRXseRQKMZEhkUlIC9OkjH2bvZQJduojvUlRUK5r5bPNmuPr1g7O6Gtnt2sFx++3AP/8JtGol1731FpCQINOgFSskG/TznwM//anc0yyq++lP6/ZWAmiARSgyJFJJTva0sPReJlBaCqSlAfPmiTdiMYT7nDiBSZdeikmnTknLzLNn5UN/9dXAmjVy7RdfyD03bRJP5tlnJUr6/ntg0CApqquoEBGqCz8bYBGKDIlUJk4Enn9eUs0ff+wRku7dgWXLpA6mY8fztS3/XrgQzuuvR7vmzZH1zjsiJNOmiYl7330SAW3aJKbs11+LyDidYv6OHi3ri4qL5Ssg0VND3koU7HkUCdCTIZHJli3AuXMy1enSRf47K0tqWpQCZs+WNPPIkTDS0pDTsiVeXLNGqnd37pQIaNAg2aIkN1f8ku3bPS0yAUlXf/WVvGbNCtFbCSoUGRJ5mJmldetECE6flsK6a68FzpwRv+SNN4CJE2v0g/mv8nKZFg0YII2lAODee2WaNGGCFNp16SJiUlYGfPih+CqZmZKxOns27nu/2IHSWod7DLYzYsQIXVBQEO5hxBaBFskFQk6OZJOysqQmZu1aj8E6c6Z4J48+CmPgwPMd7bLbtIGjWTPPYsUWLUQ0Dh+WVdJHj0r089JL4uNUV8szPvtMpkmlpbJR2zXXeN6Lne8xBlFKFWqtR3gfpydDAsfufYWsLRu8DdbWrYG2bVF9++3I/+UvRWCys6WjXZcuIjDmHktmVmrzZuCZZ0ScSktFiBISgH/8Q0SjbVvpE/PZZyJk6em+3+N//7f4QMnJFJ0AoMiQwGjKvkL+RgbexW5Wg/WJJ4CXXkLz48eRk5ODVq1ayVqkpCS596pVtVPMe/dKVPPXv0oP3sGD5dyhQ6VJVWamGL579wKjRokgvftuzfe4c6cIXHGxFO1VV8f1hm2B0KAno5Tqbf8wSNTQ2I3rA1kiUE+xm/HRR/jtpZeiOjMTHf/2Nzh27hQx+PBDQGvJPgHywXe5pADvxRclcmnfXozkw4eB22+XyMgw5F9RkWSbJk0C/vd/a75HU1hXr5ZM1ZdfenYnWLvW0xqC+MQf4zf+2quTumnMuh1r9DN0qLRQ+PhjzxIB7w+pNRKxNOM2srLgatECX3XsiDM5ORJVrFghX597DvjVr2pes22bZJQeeUQqhe+8U4SjtBS4+GJ5lrc4ZmRI1zvre7QKq/U9NiSsBIB/06WblVJPaK232z4aEvk0Zt2OdYmAt8/x1luyRGDDhprXeHkxRmKiLBXo0wfZ2dlo43AAI7w8xvT0ugvk5s71+DxduogIzZgh31s3XduyRdYzWd+j1SPyfo9xumFbIPgjMu8CuFUppbTW2+weUEMopWYBMOPoVK31snCOJ+4IcON6FBVJurldO6lJycvz9GfxtUTA9G127pSMUKdOMLp3h6u4+LzAOBwOz719eTwN+Tx33SWic//9MkUyhWPvXuC11+T9rV7teY/mtX361H6PXBDZID5FRinVTmt9AgC01le7j01WSqVqrV8L5QC9xjULQJnWeq37+1Sl1HKt9exwjSnuCGTdjpmh6dlTRGjsWNms3mqWlpZKRLJ+vUQR69bJh/rrr8+LTttWrTCouhqL0tI8AhNohssqjpMniwk8apTUzxQWSge8UaOk2dXJk7Ks4JJLZF2T0wk89ZRMzaZP97zHON6wLSC01jX+AbgKwFjv4+7XfgTgPl+vheIfgEIfx/Y2dN3w4cM1CTI7dmidna313LnydceO2q9PmKD1d9/Jf/fsqfXf/y7fX3utHDP/e9UqrW+5Rc7/4IPz15WWlp4/59y//uW5znpvK9Z7++KFF+T1lSu1Ngytf/ELrTt00LpbN62vvFLey4QJ8vrmzVo//rjWaWlaX3GF1jfeqPXo0Z5rV66Ue73wgj0/3ygEQIH28fmrUYynlPov94l1TouUUhcAWATgCe2OdkKBUioJwDGttfI6XgjgAa11HRsQsxgvLFgL6gDpt1tcLNmZvXs9DbinTJHMz4svikdSXAykpMBIS8OyZcswf/58jC4uPn8cxcVyP+u9reTn11yD5I05xTp82FNTY07TrLU1JmVlErWZ0YqvawkA/4vxigGkAqhTZLTWxwEsVErdp5R6R3sZwkqp3lrr/UEYszepAMp9HC9zv0YiCe89kB57TLyQXbtk8eJXXwEvvCAp58xMaRLVrx/w8ccwUlLOLxX40d//Lv7MoUOy6vrwYUlV+7G/kk/q8mz8Sc3H+SZtjaVGCltr/S6AfyulJtV3kVJqLICRALYqpf7PLTg/dL9slz/SAR7D10o5gCQfY5yllCpQShUcOXLEpiGROmnWTP76z50rUQ0g/sZHH0lHupMnZXO0G2+UaGbwYGD3bhitWsG1dq1n69i0NEl19+ghvsg333iMWF801ohlS03bqFUno7Uu1hZzVym11P11rFLqeaXUdwA2AhgH4I8AjgOYAxGnagALQjLyBtBa52qtR2itR3Tu3Dncw4kvcnLEWH3/fembe/SomK3Llknz7uJioFs3OVe5Z78TJ+Lw6tVYVlwMZ1UVsktK4MjPF6P2k0+kyVSnTlKhO2hQ8DvT2SFcBIB/KezZ7qxOe0jUsAbAq+6o5zxur+ZqAE8GfZQeOvg4ViuKIWEkJwf43e+A7GzJEK1YIcVt6elSHJecLJW306dLFsj0PDIy0HXyZPxy3TpkfPstHMePA4sXAzt2SCHdE0/INQkJno51wexMF2hqnviNPyKTBGAtgOXewmLF7dWsVUrVMn6CRAF8C0oHAFtteiYJhKIiWTuUnS3TpKIi2c9o2DBZ/dy7N3DFFeLPWNY5bb78clyQn4/+y5bhiowMEZ9evYCtWyWVvH271LdkZMj3f/2rPK9tWzGMO3f2dLRrrGfClpq24Y/I5Gqt5zR82nmWN3Yw9aG1LldK7VNKJWmtrQZwUn2ZJdIEGrPnUY8ewPDhnu/NPZCKioAf/EBMW4uZapw6BdemTUhPSEAOIPe/4QbxbUpKxIzt189TB/P99+KPXH+9mL8bN4qAffutZxyNFQS21LSFoK9d0loXN3Is/uACMMv8Rik1DOIPkWDT2D2PnE6Pt2E1U62VtQDQty+M4mLJIiUl4YEBAzz3GTQIeP11yUYVF0u0UlwsU5nNm2UalZUlWaePPhK/5swZ6SHT0BgbIj1dskjPPcdsUpBoMJKxWTQCQmud684ajYNMnVI1q32DT2PbOSQnywfdLMm3lvJfeCFw4MD5bUaMf/wDrrIyOFNSkL1nDxz33uu5z65dUoH785/LquiRI8X8dblkynXwYO0xtm0rQmR2tquv5QQJKVHXT0ZrnRvuMcQ8/tSM+PoAm+bp+PHyQR89GvjDH2TpQG6utMZctAhISsJH330HZ/v2yC4shGPGjJr327lT+rW89prsgVRQII3De/aUaZG564B1jGZ9TENjJCGHPX5JbRpbM2Kap+Ym9fv3S3bmkUeAyy4DBg9G9bffAlu2YN6ZM8hJTITD132OHpV7mbsIjBgh9xo4UKZt586J+WtNOXtNxVjXEjlQZEhtmlIzkpEhQrB+vaSsx44FVq4EunaF8Ze/YP6QITi2fj1aFBaizWuvAX/5S+3GT506yfdlZTWnRffcI0sS+vaVdPYbb3jOs9bHsK4looi66RIJEvVljhpbM2JdGX377SJUq1cDx4/DaNkSrtRUOHv1QsshQzzX+JreDB4srSEyM6XlZUYG8Oab8uz0dOl2N3KkrH1avFimVmaamXUtEQdFJh7xbpOwcaNMc7p3B668UkQm0JqResxiIysLrubN4ezevWY/GBPTT7H2ktm2TQTk97+Xc5o1k++feELaRezaJSnxDz4AFiyQjFR+PutaIhCKTLzhLQZ5eZ40cH6+Jw08ebJEA/7WjNRhFm/dvx+uFi3gPHsW2RkZtQUGEAH75ht5ril8ubkiHj16yDNTUkRgpkyRsU2eLMcMQ5qEr1jBupYIhSITb1jFoKE08AUXyDXWvbnqmmZ5r7p2k5KSgjGDB2PWF1/AsW6dtHzwnoLl5soapg0bakZBH34o3ktuLvCzn9UUEOu0iKIS0VBk4g2rGHhHH9Y0cLdukm5essTTee622+S8+++vvRfR119LevnRR4H0dHz++efo06cP2rdvj3v69wcSE4GLLvKktvfvFxP34EGJVmbMqJ0yHz0aePxxEcIdO2RKVFHRuHJ/btQWNphdijWKimSRotliwXsnAGvmyDtVbWZliopko7PRo6WytnVrmZYcOyYtGtq0kUijqsqzF9HkySIa8+fDePxxLFq0CGvWrKmZ+Zk5U1pavvSSnOt0irjs3Vt3NqtvXxE8l6vm7gQul//7HTWmepkEDUYysYQ/fW+tmSNrRa51+mFmerp189yzVSvgjjuk/mTRIjFeN2+W7NGbb8oH/5e/hPHKK3AZBpwdO+LG06dr7itdVCTXGEbNqOX77yUdPXVq7ejCFL66mk01RFM2oyNBgZFMrGD9MGVlyQfJ1wZk1tXGLVpIRe7zz9cUg+JiEYJPP5UV1enpUnE7YYIU1bVtK71hPvsMmDNHMjy7dsne1K1bw9msGbI/+wyOQ4dqRhx1VRJPny6ezJ/+VPN4U/rDmDR2MzoSNBjJxAqBLAWwrjZu2RL4zW9knZCZBv7b3+R4VZVML5xOKfG/7z7gxAkRnBYtxGf5+mugtBQVpaX436eekm1LbrgBjldeqd1ntw5zGOnpsqo6P1/GEGibhfr8lrqeCXDPpBDBSCZWCHQpgLna+B//kAjITAMbhhi0AwZIQdyECRIRPfigvFZVBTgcwJgxsvK5d2+goABtbrsN2f/5D7KnToXj4EHfFbf1VRL36ydCE6jv0pDfwo53YYeRTKwQyM6Ovv7ym1FHTo5s91pVJUVviYmSwn79dZkmHTggTaK2b5dK3gMH8O2Pf4yf5uejz5gxnkyQr4pbfyqJA/FH/PFb2PEu7DCSiRXq2aS+hq/R0F9+MyKaORN48klZBb18uaSaf/UraZ355ZdAdTWMjAy42rfHvzp2xNl+/aSn74oVdU9x6tjjuoYfFAj+ThGD+UwSMIxkopG6PIiGlgL485ffGhFNniz9edeulTYL//43cPo00Lo1jNmz4SoogLNbN2SXlKBFaSmQlCQi5D3F8R7vHXdIdqqp3ef89VvY8S6sUGSijYbS1PV9mPz5y+89vTDv+eyzIk7Nmsm2JYYBZ+fOyJ40CY45c8SwLSyUdgwNjTcvT8Zb1wZs/hLIFLGxKXDSdHxtKxlr/2Jmm9rGbs9qMmeObL/qC8OQbVq1rr2d68qVWnfvrnWrVlq3a6ffHDlS33fllfpkq1Zat2yp9ZgxWo8aJV/T0jxbtzZ1vA1h9/1JQKCObWrpyUQTTa358DfTMnNmzQrbN98EKipwassWIC8PE44dw5MAHK1aee67bJlknfLyZHq1bp3s+Hj6tGyRYq08DlaNCv2WqIAiE000dZdDf81hK1oDmzfDuOIK/Nfjj2NPejqwcCGaf/WVVOhOmCDZJHOF9OjRnnVPSskH3lcZf7C613kLYqBLDojt0JOJJgLxIHzhbQ43ayZFdlu2yN5IJlYfpfseDSsAAA5LSURBVE8fGG+/DVdlJZxnzqD7b34j5yQmSpvM8nL5+sorHnPZXPfUu7cU7WVl1S7jD2aNCv2WiIYiE00MGiQRwnvvSXRgZpUCqfkwjdwlSyRbNGKEpKqrq2VFdUKCLFicNg0oLoaxeTNcgwfDWVmJ7E8+gWPXLtl4raJCWjFcdJF0srOay+a6J28T2Zwide/OGpU4giITSprSbsCMLsaMkTVFiYmybmjAAGn4FKgHUVUlImP6O3l5sqJ6xw4ZV7t22PPHP8I1ZAicbdsie8MGOK6+Wt5DeTkwa5ZcU1gI3HKL577FxRJt3XVX7cipslJqaQyDnkkcQZEJFf6skK4L7/oWU6zMD+2TT8p9/MXbQDbvv3o1cO21Mn0qKIBz8GBMrajADd26wdG/vxjAqakSCa1bJ2KZlCRj+NOfZBeBDz8ELr/c9zqp998Xj4c1KnEFRSYUNLXdgLcoWD2I/HzPBvT+4l3EZt5/wwagpASFF16I3j16oOOQIbh1wwap8L36amlO9cUXMl07c0b6wpw+LVOv3/5W1h49+aQIallZzfF27y4RDAUm7qDIhILGbpZmEuyVxN4GckkJ0KcPsHo1jBtvhGv/fox2OHCv0ymp6VtukWlRYqKsYdq2zdOoKj8fGDdODF6zuO74cW5cT85DkQkFO3cChw4BH39c24vxRyRMUWjduran4ytL05D3423IJicD69bJWqTCQjjbtsWc/fuBX/8auOYaiUKKiiRl3b078MADEpGYhvPo0bKK24Rl/MQCRcZu8vI8m5zNmFHbi/EnlTtxovTX7dJFGmqbns78+bJdq7XZkz/ej7ch63TCyM+XNHWHDsh+7jk49u0TMfmf/5Hs0UUXyUZsW7fKeLZtk9fHjwfefrt2pohpZeKGImMnpheTlydtE/r0qenF9O3buFSudfcA750E/PV+LNFGdUkJ1qSkwHn6NLIvvxyOnTtF/ABJmW/fLoV1vXtLenvTJmkaPnq0CAynQaQeKDJ2Ynoxo0fXXiHds6d80BcubPgDun697BAwdGjNKcgzz0hEYXo6gXo/7mijOYBHDAMtFiyAY88e2eDenOJ07y7jfuEFT7Hd+vVSrXv2LKdBpEEoMnZiNWy9fYpz54Af/ci/8nfzPr6mIBUVHk8nQIPYMAy8//77mD9/Pi645BJpwbl2rUcIzSmRNVLhNIgECEXGTryzON6p5+Lixt3HitXT8T7PagB/+20N78cwDLhcLjidTlRVVSEhIYGGLbEHX0uzY+1f2Fo9BKsVgb/3sZ73wgvy3ytXav33v2s9ZIjWl12m9Qsv6E2bNumbbrpJ33vvvfrkyZNNf5+E6LpbPTCSsZuEBPFSMjKASZNkmhRozYg/Xe+s5914o+xltGyZZJ9+8xtpu3DTTdiclQXX66/LrgLZ2b73piYkiChtzU7EKCNGjNAFBQWhfag1ldy8ufx3QYH4MI891viNykzTtWvXutc+3XWXnNOlS63zvvjtb7FqyxYsWL6cAkOCilKqUGs9wvs4Ixk78JVKnj5ditcyMxt/X+tK55IST9Mnb6E5d04MW4uH880336Bbt25Iu+QSZO/eLQ2mCAkBbFplB3btWujvns5eHfAMw8Ddd9+Nd999l3sNkZBDkbGDpnaw84W/29ACNTrgmVmk1NRUXNKvX9O3fSUkQDhdsoOmdrDzRSCFdm4D2MjKgqtFCziTk5E9cCAcs2axOpeEHEYydtCYXroNEWB0dOSGG7CsZUs427RBNgDHoUPsfUvCAiMZO/A35RwIAUZHnTt3xv0LFmDIkCHMIpGwwkjGLoLdRd/P6MgwDGzfvh0AMGrUKAoMCTuMZOwkmOt8/IiOTJN34MCBGDJkCJRSwXk2IU2AIhNN1LO2yLoWafHixRQYEjFQZKINH9GRVWC4VIBEGvRkYoCCggIKDIlYGMlEMWfPnkWLFi0wd+5cVFZWonXr1uEeEiG1YCRjJ0VFQE4OMHeufLVW5TYRwzDwi1/8AkePHkXz5s0pMCRiocjYhb/rjBqB6cE4HA6KC4l4OF2yg6Zu5lYPNHlJtMFIxg5sWoW9bds2CgyJOhjJ2EF9Db2bNQP++te6N16rB6fTiauuugozZsygwJCogZGMHXj1czlPXh7w1FNAjx4B+TSffvopKisr8YMf/AB33303BYZEFRQZO/C1zqioSI516gQ8+mj9/WAsGIaBBx98EC+//HJoxk5IkKHI2IF1nVF+PrB5swhLaam04bROj+rxaawm78033xy68RMSROjJ2IX3OqODB6Xv7vTptc+tY+M1mrwkFqDI2Il1nVFODlBd7fs8r34wp0+fxvPPP0+BITFB1EyXlFLjlFKvKqWGKaVSlVILlFKzwj0uvwmgW16rVq2Qk5NDgSExQTRFMkkAUgEUAigHkKu1XhbeIQWAn/1gDhw4gKlTp6JXr16ea63bzQaY9iYk3ESTyEBrPTzcY2gSfvaDmTRpkuxNDdTcJG78eEmNL1wogsV+vSQKiCqRiQn86AdzXmBsXJ5ASKiIGk8GOO/LTHF/XdDAubOUUgVKqYIjR46EaogBU28Wya5N4ggJIdEUyWwFAK31PgBQSpUppd7RWl/t62StdS6AXED2wg7KCGzwRioqKpCWloZHHnmktslb3/IEH2lvQiKRqIlktNb7TIFxf78VwAilVGpIBhDk1g0nTpwAAIwdOxZLly71nUWqa3kCwO1mSdQQcpFxT3Xe8fNfQwKyD8A42wcdyBaxfmAYBu68807s3LkTANC8eXPfJ9qxSRwhISbk0yWt9UYAGwO5xi02hVrr9vaMqgEC2SK2AaweTEpKSv0n27FJHCEhJpo8maU+jqUiQMFqFEHyRhq1VKCetDch0UBUiIzWep9Sqtx6TCk1BcAaq09jGwFuEeuLvXv3Nn4tUjA3iSMkxESFyACSLbIsI0hyH5sdkodPnCgm73XX1Zwymd6Iy9XgLVJSUpCVlYVrrrmGSwVIXBE1IgOcT0uHniZ4IwUFBbjooovQtWtXTJo0KYSDJiQyiCqRCSuN8EZMD+aSSy7BAw88EMLBEhI5UGQCIQBvxGry3n333TYPjJDIJWqK8aIJNpwixANFJsicO3cOGzZsoMAQ4obTpSCitUazZs3w0EMPQSlFgSEEjGSChmEYeOyxx3DmzBm0bduWAkOIG4pMEDA9mBMnTuDs2bPhHg4hEQVFponQ5CWkfigyTeCTTz6hwBDSADR+m0CXLl0wfPhw3HPPPRQYQuqAkUwjOHjwILTW6NWrFx566CEKDCH1QJEJEMMwMG/ePLz99tvhHgohUQFFJgCsJu+YMWPCPRxCogKKjJ8wi0RI46DI+MGxY8fw9NNPU2AIaQTMLvlB+/btsWjRIvTv358CQ0iAUGTqwTAMKKVw8cUXY/jw6N4hl5BwwelSHZgezPr166F1cPaGIyQeocj4wGryLlmyBEqpcA+JkKiFIuMFs0iEBBeKjBe7du2iwBASRGj8uqmqqkJCQgJmzpyJyspKJCYmhntIhMQEjGQgU6S5c+eipKQESikKDCFBJO5FxvRgkpKS0K5du3APh5CYI65FhiYvIfYTtyKzY8cOCgwhISBujd8+ffpg/PjxyMzMpMAQYiNxF8kUFRXh9OnTaNOmDWbPnk2BIcRm4kpkDMPAww8/jPz8/HAPhZC4IW5ExmryTp8+PdzDISRuUPGw+G/AgAG6f//+NHkJsRGlVKHWeoT38biIZI4cOUKBISRMxEUko5Q6AuBAuMfRAJ0AHA33ICII/jxqEg0/j15a687eB+NCZKIBpVSBr1AzXuHPoybR/POIi+kSISR8UGQIIbZCkYkccsM9gAiDP4+aRO3Pg54MIcRWGMkQQmyFIkMIsZW4XYVNIhOl1CwAZe5vU7XWy8I5nnCilBoHYDaApQDKAUwBUK61jip/hpFMBKKUGqeUelUpNUwplaqUWuD+8MU0psBorddqrdcCWKuUWh7ucYWRJACpAArd/zpGm8AAjGQiFesvVzmA3Dj5iz5ba31+q06t9T73X/O4xfrziFYoMhFKLPxyBYJSKgnAMB8vlSulxmmtN4Z6TCQ4UGRIpJAKidq8KXO/Fpe4I7kkyM9mWDRGtBSZCCUWfrkCpAM8hq+VcsjPIR7ZCsi0EQCUUmVKqXe01leHd1iBQZGJTGLil4s0DfP/v+X7rUqpEUqpVO/XIhlmlyIQrfU+6y+R1norgBFKqVifNnTwcSxeo5i62AcgqsxwRjI24p7yPODn6bMb+Otk/nJFXQrTTwrgW1A6wB3ZxRPuPyiFWuv24R5LU6HI2Ig7IxJQViSWfrkCQWtdrpTap5RK0lpbDeCkOM4sLfVxLBUB/k6FG06XIpOY+OVqBC4A54sOlVLDEPvv2SfuqLZGtk0pNQXAmmjyYwCuwo5IlFKzrJWd7l+uq7XWs8M4rJDgrvrdB3dBYhxk1erFUumdBADR+POgyEQosfDLRQhAkSGE2Aw9GUKIrVBkCCG2QpEhhNgKRYYQYisUGRJ23I25liultFLqmFLKVcd577jPeTUOlljEDMwukYhBKfUqgHF1VTu7xecV91ouEiUwkiGRxHIASb664bkLEikwUQhFhkQM7jVK5QButh53Ly8ABSY6ociQSGMNgFvMb9xtOce5G4uTKIQiQyKNV1FzyrSISyqiGxq/JOJQSh2DrL7eAtmpwVfvXxIlMJIhkcgayEZmGykw0Q9FhkQiHSACQ6M3BqDIkEhkGIB3wj0IEhwoMiSicGeT4qELYNxAkSGRxjiANTGxBEWGRATu9UsuSJ9fKKVcXJ8UGzCFTQixFUYyhBBbocgQQmyFIkMIsRWKDCHEVigyhBBbocgQQmyFIkMIsRWKDCHEVigyhBBbocgQQmzl/wEVY0dl0G7RqQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Init Loss = 4.8537\n",
      "Final Loss = 0.8792094\n",
      "Fair dummies test (regression score), p-value: 0.000999000999000999\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY0AAAElCAYAAAAY6L5SAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAYn0lEQVR4nO3dPXAj533H8d9fJxUejSc4KCluVFjCKT3Fo1tnYuHcO0PqyjQWqEpqEjKXGanMiRenjohLipQnUq4zJuRMnDIkfak0GYdw3Pt4cJwy8j/F8yy5xC2AZxdvC/L7mcHcYRe7+yxf9sfnZZ81dxcAACleWXYBAACrg9AAACQjNAAAyQgNAEAyQgMAkIzQAAAkIzSwUGbWMTMf8TpbYhlemNmRmbUWUYY6MrPN+HWY9DpZdlmxPK8uuwC4se5L6ktq5JYNUjaMF61zd78/ZRluu/sg7nNd0kNJJ2a25e69Kfc9FzM89yI9Se8NLTuR1JW0P4fjYQURGliW4+yCXcFMLmD547v7qaQtM9uXdCDp9iyOMQdzu3jHr8dpfpmZSdJZ/PpMxcwaktrufjjtvrA8NE9h5bh71927c9r3tiSZ2d489j+teZ77AmwoBDJWGKEBvOwLSZvLLgRQR4QGlqVpZo38K1thZm0zO8l1UF/5q9/MDszsIPd+38z24nZHZrYzZdlOJF10iA8fLy5zM9vMvd+Prx0zO4vl3onndpR19Me+k2m2KSrLXq6Tei+3fCcuy/bTmfLrckWZ45rZkaSj3NfOzaw9y/JgMQgNLMuZpBf5V27kUkvSrqS7krYk7Uy44DUVagYHChemaZtv+pJUciRVU1JWxvuxDHuSvor/3lPo6D+YcpsrYoBsKnRgb0najBfpdtzHewpfx21J5yXOZ6yyx40d91uS5O4WX7UcbIDx6AjHstwe1RE+1GbfN7NDhQvOuDBoSbo3iw7buC+5e7/kdgN3fyxJZvZI0o5Ch38vLtvXyx3ZVbZRXNdSuHBfnLeZ7SqMAnukEDj9+HUuey4jLeu4qAdCA7UUm34eKFzA1xWGg45zOqPAkMJfyFVGdh1n/3H3QRx5lL+noegv/SrbZLLmnZO4XV4vbvvCzHqSjrJwmoFlHRc1QPMUaifei/BQ0lOFGkbKEM1Z/kXbUegML6soaCY1CVXZ5sr2ueae/Gvg7lnz3kDS3nBfyJSWdVwsGTUN1Ers9F2XdDdrHir4a3aex99TuOFwd2hVvqO+oXroSWqY2fqoWla8J+Iw1txmdfFe1nFRA4QG6iarMezGC/i6Qvv5XDtNY1htK9Qy7hf0t7RiW/65pCfzLEsqd++bWVfSgZltK3zt1hWa9Z5K+m78d6DLO/CXedzsj4DNuE50hq8emqdQK/Fi/Vjh4n2icNE5VLU+hrFyQz9dYcRSU6GGM3whexTXnSn81fxI4QI48zKVFW9GzDrLz3TZrHeq0B/0VVy+oTh6aVnHjbWSU4WvYS1vnsRktuhnhMeq/fsKf80l/RDH4ZZZO2+LjjUAWI6FNk/FJoCsip80Bj4LjGy+GjNrmdl+Nt0DAGBxFl7TkC7C44m730v47Mnw58zsLI7QAAAsUK37NGJT1nrBqgFTEADA4tU6NBSasEaNY7+xD8sBgGWp+5DbpopvdBro6sN7LsQ+kI4kfetb37r31ltvVT/6734n/f73kz/3yivSt79d/TgAUCNff/31b9z9j4rW1T00SovzFnUlaWNjw4+PjydsMWZff/I9ffPhn0/83K3P/0n2rz+vfBwAqBMz+/WodXVvnpJCbWNYXe7IBYAbpe6hcazigGhq6LGUAID5q3VoZFMrF8z102D6AQBYvGWFRlGTU3bj3sFQSOzp8kE12T0eBAYALMGi7wjPHt5yX9J6nJDuLPfQnZbCXP1NXU5o1s09EayhMI0Id4MDwBIsNDTiVNeP46tofU/S7YLl0z6+EwAwA7Xu0wAA1AuhAQBIRmgAAJIRGgCAZIQGACAZoQEASEZoAACSERoAgGSEBgAgGaEBAEhGaAAAkhEaAIBkhAYAIBmhAQBIRmgAAJIRGgCAZIQGACAZoQEASEZoAACSERoAgGSEBgAgGaEBAEhGaAAAkhEaAIBkhAYAIBmhAQBIRmgAAJIRGgCAZIQGACAZoQEASPbqMg5qZh1J5/Fty90fJ26TaUjquvtgHuUDABRbeGhkgeHuh/F9y8z23X17zDY7GgoJM9uXNHIbAMDsLaN5ajsLDEly976k9oRtvltQqxiYWWPmpQMAjLTQ0IgX+fWCVQMzGxccrYL1DZqnAGCxFl3TaEkqutCfx3Wj7Eo6MrM96aKJa3/2xQMAjLPo0GjqsgM8b6DQuV3I3XuS7knqmJlL6rv76XyKCAAYZSmjp8oys5akDUlvS3qoUOvYdvduwWc7kjqSdOfOHT179qz6gdtt+WuvTy5fuy1NcxwAWBHLCI1mwbJJHdq7udFVu2b2VNJXZtaPtZALMUi6krSxseFra2uVC+off6Rv3nlz4udu9XqyTz6tfBwAWBWLbp46VnFANCUVNjfFDvCj/LLYNLUl6f6sCwgAGG2hoRFHO/ULhso2hmsMCY4lPZ9NyQAAKZZxn8aeYp+DJJnZuqRe7n3LzA6yYIlh8qBgPx3FZigAwGIsvE/D3btm1onNTg2FaUTyd3a3FG72a+pyeO4HcbjtWXzfkHTIfRoAsFhLGT1VNOopt64n6fbQsoHCvRoAgCVillsAQDJCAwCQjNAAACQjNAAAyQgNAEAyQgMAkIzQAAAkIzQAAMkIDQBAMkIDAJCM0AAAJCM0AADJCA0AQDJCAwCQjNAAACQjNAAAyQgNAEAyQgMAkIzQAAAkIzQAAMkIDQBAMkIDAJCM0AAAJCM0AADJSoWGmX1/XgUBANRf2ZpGz8yem9nfm9mfzqVEAIDaKhsatyU9lHRX0lcxQB6Z2drsiwYAqJtSoeHuv3X3rrv/QFJT0l8pBMiBmf0yBsh35lFQAMDyVe4Id/eBpOeSziVZfN2T9Csz+2cze2sWBQQA1Efp0DCztdin8VzSgUKT1ba7v5Orgfxa0k9nW1QAwLK9WubDZvZfkt6W9AuFpqkv3P23+c+4+8DM9iR9MLNSAgBqoVRoSOpJ+tzdn0343HNJ90etNLOOQrOWJLXc/XHKwc1sR9Ig29bdD1O2AwDMRtnQOJH0QNKV0DCzn0p66u7/KIUOc0lfFe0gC4zsgm9mLTPbd/ftcQc2syNJW7EvRWb2wsx62XsAwPyV7dPYkfTvBcu7kj5M3Md2vobg7n1J7XEbxBrGwVBA3CMwAGCxyobGXUmnBctPJK1P2tjMGiM+NzCzccHxUKFp7EIMGwDAApUNjb6kdwuWr0v6VcL2LYU+iWHncd1LYtA04v83zaxtZjtxOQBggcr2aXQl/YOZDdz9X6SL+ai6kh4lbN/UZQd43kAxGApsZOtz/SDHCsN9X+psj30mHUm6c+eOnj2b1Gc/Rrstf+31iR+zdlua5jgAsCJKhYa7PzazNxSmEPHcqifu/uPZFu2KhkItJyvHwMyaZrbu7leay9y9qxBi2tjY8LW16jOc+Mcf6Zt33pz4uVu9nuyTTysfBwBWRdmahtx918z+RuGvfJd06u4pTVOZZsGycU1N/Xjc4Watc4UO9KI+FgDAHJQODeliSG2VeySOVRwQTY24+Lt738xG7Y/RUwCwQKVDI84p1VLBxd/dfzJu29is1DezxlDNoeHuvZEbSqdm1hoaMdVSCCEAwIKUnUbkLyV9pjA54TCXdCthN3sKHdWP4z7XlRtOa2at+JkPcsGyG1/buW36w/0ZAID5KlvTeCjpbyU9Gp5zKpW7d82sE+/LaChMI5K/G7yl0FfRVGx+cveemTXiTX6S9Ia7j5ymBAAwH2VD41xh7qlKgZGJI5xGrespzJw7vJx5pgBgycre3LcvaXMeBQEA1F/ZmsZdSR0ze6jQCX1l9JK7P5hVwQAA9VNlyG2+mWjkWFgAwPVT9o7w1JlsAQDXUJXHvf7IzI7N7P9yy56a2Q9nWzQAQN2UCg0z+0zhfolHQ9seSPrrGZYLAFBDZfs0OpI23f1nQxMWJj1PAwCw2ko3Tyk8/1u62gneUtrzNAAAK6xsaDyRdGhmawrThsjM3pX0eXwBAK6xsqOnduPcUKeSZGbPFaYC6c75eRoAgBqo8jyNrRgc78VFvZLP0wAArKiqz9PoK/ckPQDAzVB2avQ/G7N64O4/m7I8AIAaK1vTGDXTrCvUPP54uuIAAOqs1Ogpd39l+KXw3IsvJW3NpYQAgNqocp/GFfHpeh8oPG0PAHCNTR0aUVPhaXsAgGusbEf404LFDYU7wgcF6wAA10jZmoYVvH4r6ReipgEA117ZO8Lfn1dBAAD1N6s+DQDADVC2T+P7ZT7PzX4AcL2Uvbmvp3Aj3/CzwUctu1WxXACAGirbPLUl6StJd4du8PtS4eFM+Rv/CAwAuGbKhsZnknYKZrXdlfRwNkUCANRV2dB4Q9K9guUt8bhXALj2yobGgaSumf0wWxA7x79Q6O8AAFxjZe/T2DazpqQvzczjYlN4kh8TFgLANVf1yX1v67I56pQn9wHAzVA6NMzsR5I+lLTm7q/GZU8lPXX3n8y4fACAGinVp2FmnymMlHqkq/dlHKjE6Ckz65jZZnztlClD3P6g7DYAgOmVrWl0FO7H+JnZlXv5TpQ4esrMOpLO3f0wvm+Z2b67byduvy5ps1yxAQCzULp5StLz+G8+NVqSUvs1tt39Ytiuu/fNrMwMua0SnwUAzFDZIbdPJB2a2ZrCNCEys3clfR5fY5lZQ8U1kkFKcJjZZlZDAQAsXtkht7tm1lIYYisze67wEKauu/84YRejHtZ0rgk1iHjcfpnyAgBmq+qQ25ak9+KiXokht02FgBg2UAifcdapZQDAcpWdGv0vFGoVfS3wr/7YdJV0x3nsaO9I0p07d/Ts2bPqB2635a+9PvmY7bY0zXEAYEWUrWn8QKEv4++mOGazYNnIWkbsB5G7Jz2D3N27krqStLGx4Wtra1XKGPb18Uf65p03J37uVq8n++TTyscBgFVRNjT2JH1hZj13/48KxztWcUA0FftJCnSki6G2F+L9HYMYEgCABSgbGvclvZB0amanGmqicvcH4zZ294GZ9c2sMVRzaLh7YfOTuz8eXmZme0XLAQDzVXbIbUOhRvClwn0ZNvRKsadYe5AuahC93PuWmR1kzVIAgPooO+T2w2kP6O7dOI1IWyGEWkN3g7cktRWarK70Y8RttuL/9yUdjKqhAABmb2RomNkvJd139/8eWv7W8LKyxvVDxBC4PWZdT1LSlCMAgNka1zx1V0Od1mb2B5LO5loiAEBtle3TkNL7LgAA10yV0PDJHwEAXEeTQiPphjoAwM0wafTU22aWnyuqKUlm9m0VNFO5+//MsGwAgJqZFBpFw1lNo2sgt6YrDgCgzsaFBsNaAQBXjAwNd3+yyIIAAOqvyugpAMANRWgAAJIRGgCAZIQGACAZoQEASEZoAACSERoAgGSEBgAgGaEBAEhGaAAAkhEaAIBkhAYAIBmhAQBIRmgAAJIRGgCAZIQGACAZoQEASEZoAACSERoAgGSEBgAgGaEBAEhGaAAAkhEaAIBkry7joGbWkXQe37bc/XHiNpJ0V1JD0q67D+ZURABAgYWHRhYY7n4Y37fMbN/dt8dt4+7d3PtNSScKAQIAWJBlNE9tZ4EhSe7el9Qe9WEza2koHOL2zRgeAIAFWWhomFlD0nrBqoGZjQwOSZ2CZeeSmjMpGAAgyaKbp1qSivohzuO6l8SayO0R+zqeXdEAAJMsOjSauuwAzxsodG4nif0iPXc/HbGuI0l37tzRs2fPKhZVUrstf+31yeVpt6VpjgMAK2Ipo6emEfs4tt39XtH62GHelaSNjQ1fW1urfCz/+CN9886bEz93q9eTffJp5eMAwKpYRkd4UT9Eci1D0p6k92ZUFgBACYsOjWMVB0RT0ktNTcPMbE/cnwEAS7PQ0IgX+34cRZXXcPfeuG1jX8V+7BjPlo0bcQUAmLFlNE/tKTeE1szWJfVy71tmdpAPlhgOx1lgmFmDwACAxVt4R7i7d82sEy/6DYVpRPJ3g7cUbvZrKty/0ZJ0JElmNry7oqG4AIA5WcroqfyUIAXresqFQaxdvJQWAIDFY5ZbAEAyQgMAkIzQAAAkIzQAAMkIDQBAMkIDAJCM0AAAJCM0AADJCA0AQDJCAwCQjNAAACQjNAAAyQgNAEAyQgMAkIzQAAAkIzQAAMkIDQBAMkIDAJCM0AAAJCM0AADJCA0AQDJCAwCQjNAAACQjNAAAyQgNAEAyQgMAkIzQAAAkIzQAAMkIDQBAMkIDAJDs1WUc1Mw6ks7j25a7P57HNgCA2Vp4aGQXf3c/jO9bZrbv7tuz3AYAMHvLqGlsu/u97I27982sPYdtAAAzttA+DTNrSFovWDUYFQJVtgEAzMeiO8JbkgYFy8/julltAwCYg0WHRlOXndl5A0mNGW4DAJiDpYyemqfYad6Jb//XzP5zqh3+/N+Gl/yhpN8UHHiqw9Rc8TlfbzftnG/a+Uqc8zjfGbViGaHRLFg2qcaQvI27dyV1yxYqlZkdu/vGvPZfR5zz9XfTzlfinKtadPPUsYov9k1JpzPcBgAwBwsNDXcfSOrHEVF5DXfvzWobAMB8LGMakT1d9jnIzNYl9XLvW2Z2MBQSY7dZsLk1fdUY53z93bTzlTjnSszdZ1GQcgcNndV9hWanK1OCxHsvDiTdc/d+yjYAgMVYSmgAAFYTs9wCAJJdu/s0pnETZ9+d4pwl6a5Cc+FuHLCwEqb9npnZgbtvzb5k81P1nM1sR+FG2nNJyiYNXQVT/mxL4We7uwo/27EP+H1J91N/Niv/Hrg7r9BE15G0mXvfkrQ/623q9Kp6zkPvNyWdLftc5nnOQ9uvh1+b5Z/LvM9Z0pHCKMXs/Yv8+zq/Kv5s7wyf3yr8Psefyc34Opnnz4S7Exq5L9pLX+xJF8Mq29TpVbb88Qdrr2D5i/wPYJ1f037P4i+mL/s85n3O8QI6/AdCa9nnMudzPihYtrdCQbleIjQq/x7Qp6GbOfvuFOXvFCw7V/Fd+7Uy7ffMzDZ9hZpnpKnO+aGGhrV7bjRjnU1xzq2C9Q1fgeapMqb9PSA0gps4+27p8rt7391vj9jX8QzLNi+Vv2dm1lIY8r1qSp9zvKg04v83zaxtZjsFN9jWVdXv866kIzPbky7a/PdnX7ylm+raRWgEN3H23ZmUP/5i9dx9FaZ0meac11fkHIdVOeeNbL27H3qYeaGrcP/UKqj0fY7neU9Sx8xcUn9Fv+eTTPW7T2igsvjX97a73192WeYpVtlv2pQ1DeVqVrGJphlnY7iW4s/zhqS3JT1WqHUUNcfeaITGpbnOvltT05Z/T9J7MyrLopQ656xJZsXbtct+n/tS4TmfS6p9f11U5Wd719277j5w912FWsfeKvRRVlD5d5/7NIKbOPvuVOWP7b4rdX+Gqp1zR7qY7+xCdv+Ch6n466z0Obt730Y/H2YVvt+lzzkGw1F+mbufmtmWpPu6XjXNqX73CQ2Fv6jMrG9mwyMlxs6+W3abOpmm/FkHoV+dG6xd9/Ou+H1+6YYnM9srWl5HU3yfT82sNTRiaiUGPMz4d/NYxSONVta0Xx+apy6t+uy7VZQ+5/gX2XF2MTGzxopV36t8n1ddlXPeja/8NqvUMVzqnOPF8kHBfjpandlwC4e9z/raxYSFOTdx9t0y5xw7Cs9G7Or2qjRVVfk+59Zt6fJCclD32lWm4s/2pi6HYL4R2/lXRtlzjhfVh7r8GW9IOqz7/Snx93JToRmtrdCJf5Y1nc762kVoAACS0TwFAEhGaAAAkhEaAIBkhAYAIBmhAQBIRmgAicxsP05kN+4zB2Z2kri/AzNblUkAAUmEBlDGgXRx/8Iom7qe02kDkggNIFm8kW+g4juH82HyxcIKBSwYoQGU01WoTRR5oPBskZW4Mx6ogtAAynkqjWyiomkK1x6hAZQQJ+zrS9rOL89CJP8M8fiY1BMzczN7kT1GdJSijvG47ebQsr24v4n7BGaNqdGB8vYVHs6Tn1p6W9Lh0OdaCjPF9uP/j8zsYiK5KmKorCs8/KopaX/afQJlUNMAysvC4f3csraGmqbiU+B67t6PneiHCrPkVpKbzXTL3U/jPnc1VOsB5omaBlBSnCb+VCEAunGK6UHRNOmxaemBQk1j2uetZM8tORnzZD1grggNoJp9XdYstlQwzDZ3k98jhcdo7mn6Z8gP3P32lPsAKqN5CqjmC+miJpE95OZCfBLaukJTUpkH+eSfkjgcMD1JjeHnlQOLRGgAFcQO8J6kJypumspCYjc+bnNTo+/vyGvFzzfivvPH7Cs+MTCOzGqZ2SZTkWCRCA2gugOFmsFLI5diqDxWeDTsicKjOA8V7igf5ZHCiKizuO9HCuFzsY27b+uyaexM4fGkT6c/FSANj3sFACSjpgEASEZoAACSERoAgGSEBgAgGaEBAEhGaAAAkhEaAIBkhAYAIBmhAQBIRmgAAJL9P1bRZ+u69scCAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plot_groups_pointwise(x_axis_0_test,\n",
    "                      x_axis_1_test,\n",
    "                      A_test,\n",
    "                      Y_test,\n",
    "                      Yhat_out_test,\n",
    "                      y_u=None,\n",
    "                      y_l=None,\n",
    "                      point_color='r.',\n",
    "                      shade=None,\n",
    "                      filename_0=\"supp_synth_use_A_baseline_0.png\",\n",
    "                      filename_1=\"supp_synth_use_A_baseline_1.png\",\n",
    "                      save_figures=save_figures)\n",
    "\n",
    "p_val_base = utility_functions.fair_dummies_test_regression(Yhat_out_cal,\n",
    "                                            A_cal,\n",
    "                                            Y_cal,\n",
    "                                            Yhat_out_test,\n",
    "                                            A_test,\n",
    "                                            Y_test,\n",
    "                                            num_reps = 1,\n",
    "                                            num_p_val_rep=1000,\n",
    "                                            reg_func_name=\"Net\",\n",
    "                                            lr = 0.01,\n",
    "                                            return_vec=True)\n",
    "\n",
    "\n",
    "plt.clf()\n",
    "x, bins, p=plt.hist(x=p_val_base, bins=30, color='salmon', alpha=0.7, range=(0, 1), density=True, ec='red', align='left', label=u'Baseline')\n",
    "plt.grid(axis='y', alpha=0.75)\n",
    "plt.xlabel('Value')\n",
    "plt.ylabel('Frequency')\n",
    "plt.title('Fair Dummies Test')\n",
    "for item in p:\n",
    "    item.set_height(item.get_height()/sum(x))\n",
    "plt.ylim(0, 1)\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# step size\n",
    "lr = 0.01\n",
    "\n",
    "# inner epochs to fit adversary\n",
    "dis_steps = 40\n",
    "\n",
    "# inner epochs to fit loss\n",
    "loss_steps = 40\n",
    "\n",
    "batch_size = 6000\n",
    "\n",
    "# utility loss\n",
    "cost_pred = torch.nn.MSELoss()\n",
    "in_shape = X.shape[1]\n",
    "out_shape = 1\n",
    "\n",
    "model_type = \"deep_model\"\n",
    "\n",
    "# equalized odds penalty\n",
    "mu = 0.99\n",
    "second_moment_scaling = 0\n",
    "\n",
    "# total number of epochs\n",
    "epochs = 500"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "fair_reg = fair_dummies_learning.EquiRegLearner(lr=lr,\n",
    "                                                pretrain_pred_epochs=0,\n",
    "                                                pretrain_dis_epochs=0,\n",
    "                                                epochs=epochs,\n",
    "                                                loss_steps=loss_steps,\n",
    "                                                dis_steps=dis_steps,\n",
    "                                                cost_pred=cost_pred,\n",
    "                                                in_shape=in_shape,\n",
    "                                                batch_size=batch_size,\n",
    "                                                model_type=model_type,\n",
    "                                                lambda_vec=mu,\n",
    "                                                second_moment_scaling=second_moment_scaling,\n",
    "                                                out_shape=out_shape)\n",
    "\n",
    "init_seed(seed)\n",
    "\n",
    "input_data_train = np.concatenate((A[:,np.newaxis],X),1)\n",
    "fair_reg.fit(input_data_train, Y)\n",
    "\n",
    "input_data_cal = np.concatenate((A_cal[:,np.newaxis],X_cal),1)\n",
    "Yhat_out_cal = fair_reg.predict(input_data_cal)\n",
    "\n",
    "input_data_test = np.concatenate((A_test[:,np.newaxis],X_test),1)\n",
    "Yhat_out_test = Yhat_test = fair_reg.predict(input_data_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Fair Dummies (All): Test Error = 1.0606429853029742\n",
      "Fair Dummies (A=0): Test Error = 1.061409485353968\n",
      "Fair Dummies (A=1): Test Error = 1.0533977557298904\n"
     ]
    }
   ],
   "source": [
    "print(\"Fair Dummies (All): Test Error = \" + str(np.sqrt(np.mean((Yhat_out_test - Y_test)**2))))\n",
    "print(\"Fair Dummies (A=0): Test Error = \" + str(np.sqrt(np.mean((Yhat_out_test[A_test==0] - Y_test[A_test==0])**2))))\n",
    "print(\"Fair Dummies (A=1): Test Error = \" + str(np.sqrt(np.mean((Yhat_out_test[A_test==1] - Y_test[A_test==1])**2))))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARkAAAElCAYAAADQqJbBAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO2deXiU5bnG7zcBQhgKYZMEZEkmi0WjCCQu2FYwiqfipWGpoMRqrCxFj4oWAQUNqDgctdqjVLEGNS4oYanLwWrUcqzth2HTIApCQCoGBULAw07ynj+e+ZhvJjPJZJn9/l0XVzLfMvMa89153mdVWmsQQkigiAv1Aggh0Q1FhhASUCgyhJCAQpEhhAQUigwhJKBQZAghAYUiQwKOUuoDpVRSqNdBQkObUC+ARDdKqTEA8gB0BVATxM+dCKDa+TJNa70gWJ9N3KElQwLNdc6vacH6QFNgtNalWutSAKVKqeeC9fnEHYoMCRhKqekA7nW+DOZ2aZJTXAAAWutKiDVFQgBFhgQEpVQagBrnAw4EyZJx+n4GeTlVo5Si0IQAigwJFJO01ouc39cAsAfpc9Pg3fdTjSBu2YgLOn5Jq+O0GMothyohjt/G7vPbb6K1nuTjVFe4HL5WahDcLRtxQpEhgWCshwj4ZUU0IBwkgqHIkFbFGdnpqpRyWA4He5vizWqiFRMiKDKk1XA6Xe1a67Eex/cDcHi/q9VZC++C0hXA+iCtgVigyJDWZKLW+l4vx/1KwlNKLfXzc6p9ba201jVKqUqlVJLW2vq5SVrrMj/fn7QiFBnSKjgzeyt9nK50XpNmCWnXw9MCagEOABMBLHB+7iAAFJgQwRA2aRFKqTSnBbIUQD3rwumjMa2b54KRq+IMndcopfLMsgY6lUOHYo9fQkggoSVDCAkoFBlCSEChyBBCAgpFhhASUGIihN29e3fdv3//UC+DkKjk8OHDqKmpwZ49e/ZprXt4no8Jkenfvz/Wrl0b6mUQEnUYhgGHwwG73Y7HH3/8W2/XcLtECGkWVoEpKiryeR1FhhDSZD7//HM3gbHZbD6vjYntEiGkdcnIyMCIESNQUFDQoMAAtGQIIU1g06ZNOHbsGDp06IDJkyc3KjAARYYQ4ieGYWD27Nl45ZVXmnQfRYYQ0ihWJ+/48eObdC9FhhDSIJ5RJH+2SFYoMoQQn5w4cQKLFi1qtsAAjC4RQhqgXbt2mDdvHpKSkpolMAAtGUKIFwzDwMsvvwytNXr37t1sgQFoyRBCPLD6YE6cOIGEhIQWvR8tGULIaTydvC0VGIAiQwhx0tIoki+4XSKEAABqa2uRlZWF2bNnt5rAABQZQmKempoaJCUlYejQobjooosQF9e6GxxulwiJYQzDwK233oqNGzcCQKsLDECRISRmMX0w/fr1Q0ZGRsA+hyJDSAwSKCevNygyhMQYO3fuDJrAAHT8EhJz9OvXD7fccguGDRsWcIEBaMkQEjOUl5fju+++g1IKI0eODIrAALRkCIkYKiqAFSuAqiogJQXIzweys/271/TB5OTkYNasWYFdqAe0ZAiJAIqLgRkzgNRU4Oab5euMGXK8MaxO3jvuuCPwi/WAlgwhYU5FBbBsGVBSAnTtKsdyc4GrrgImTABycnxbNMGMIvmClgwhYc6KFcC4cS6BMenaFRg/Xs57Q2uNd999N6QCA9CSISTsqaoCrrzS+7nMTOBf/6p/XGsNpRRmzZqFurq6kAkMQEuGkLAnJQXYssX7ua1bgeRk92PmVIEjR44gMTExpAIDUGQICXvy84ElS4Dqavfj1dXA66/LeRPTB3Ps2DForYO7UB9wu0RImJOdDYweLU7e8eNli7R1qwjMmDEup284OHm9QZEhJAIoLJQo0ooV4oNJTgYcDpfAlJeXh6XAABQZQiKG7GzfoepevXohJycHd9xxR1gJDECfDCERzY4dO05PFJg1a1bYCQxAkSEkYjEMA9OmTcNbb70V6qU0CEWGkAjE6uTNy8sL9XIaJKJERimVp5RaqpQapJRKU0pNV0pNDPW6CAkm4RpF8kVEiQyAJABpANY5/3XTWi8K7ZIICR4HDx7EE088ETECA0RgdElrPTjUayAkVHTu3Bn33Xcf0tPT6wlMS1pBBJJIs2QIiUkMw8Ann3wCADjvvPPqCUxLWkEEmoizZJRSeZBtUw2AQVrrBSFeEiEBxfTBZGZmYujQofXGlrSkFUQwiDSRWQ8AWutKAFBKVSulPtBaX+55odMhPBEA+vbtG9RFEtJaWJ28c+bMQVxcXL1t0Q8/NN4KIpQiE1HbJa11pSkwztfrAQxRSqV5uXaR1nqI1npIjx49grpOQloDb1Ekb9uid97xXaWdmQns2RPcdXsSaZaMNyoB5AFglIlEFdu2bXMTGF/bos2bgbffBq67rr7F4q0VRLCJGEvGmRdzINTrICTQHD9+HABwww034OGHHz7t5PXVIe/66wGlgFdfdT/urRVEKIgYkXEy38uxNABlwV4IIYHAMAxMnjz59OiShISE0+eqqoCsrPr3ZGcDI0eKoJSUAGvWyNcJE9xbQYSKiNkuaa0rlVI11mNKqTEA3rT6aQiJVKw+mC5dutQ7b3bIy82tf29WlgjNjh3eW0GEkoiyZLTWi5RSE53/pgNI01pPCvW6CGkp/pQKNNQhb5HTI1lVJQITLol4QARZMiYsIyCRjLesXKU2+VWL5KtD3oIF4pO58EKxaLZskQjU6NHS7CrUqHDpAxpIhgwZoteuXRvqZZAYp7hYokPjxrnEYMkS4JprjiM+/hWMGzfOr1okU6j27BFxqagAVq50dwhXV4sYBXPLpJRap7Ue4nk84iwZQiIRb+HnhITPMWxYOiZOtMHhuAX+1jpaO+TNnSsWTLgm4gER5pMhJFLxDD8bhoEHH3wQb7/9coMD2hrDV8QJCI9EPIAiQ0hQsIqB1cl74403tkgMmjqTKRRQZAgJAqYYeIsitUQMmjKTKVTQJ0NIEMjPB2bMOImUlGI3gTHFwOFo3vv6O5MplFBkCAkCIgZtUVo6Dxdf3BGbNtlaTQwam8kUaigyhAQYwzDw+eef49Zbb0VOTk+sWAGsXdu6YtDQTKZQQ5EhJIBYfTAnTpxAdnb7sBWDQEGRIaQV8JbJe/iwu5O3ffv2Tbo/WsSIIkNIC7Fm8l55pZnWb6BNGwcyMhqfKuD9/vApC2gpFBlCWoCvRlJ9+rTF44//HNddd1+DAhPu/XlbA+bJENICPDN59+/fDwC49NLBGDv2Yaxa1XCtgK9GVNaygEiHIkNIC/DM5J04cSLKy8sBAFlZqtFM3kgoC2gpFBlCWoBnJm9qaioGDBgAwL+0/kgoC2gp9MkQ4gV/oz2SyWugXTtHszJ55X7xwXi2amhJJnA4QZEhxIOmRHu6dPkObdo4cOqUHQMGFDU5kzcSygJaCptWEQKX5bJpE7BxowjNJZe4zjfUBOr9999Hly5DsWqVDXv2NK/9pbURVbi1z/QXX02rKDIk5rFaLp98AsTFAf/+d33LpaREGnXn5wOLF3+Ggwe7o1+/tIgUhEDgS2To+CUxjTVPpaAAiI8XYSkpAUpL5bxJZibw8cfig9m+fT46d34lrAbbhysUGRLTeOapmNEeb3kqH34IHDwoTt6MDDvmzbsbBQXeBYm4oOOXxDRVVeLcNTGjPampsnX69FM5Pnw4sHSpgb59608VCKd+uuEILRkS03jmqWRnA23bAiNGAO++C9TWioUyapRGly5l6NXLey1StCTOBQJaMiSmsVouH30kg+s3bxan78cfA507A59/XofCwjgsWTIdcXGncM89HerlzkRL4lwgoCVDYhYzbLx9u1gu//gHsGsXMGgQ8N57QF4e8NRTBoYPvxcvv/wT9u9vB8PogPPPh5vDN5z66YYjFBkSkxQXi0i0aQN06wbcdRewbh1w4gTw858Dy5cDP/5oYP58B3bt0khIiMMvfgH06QPcfz+wejVwxx3AU08B11wTPYlzgYDbJRJzWMPWTz8NTJ4s4evdu4GyMqB/f6BNG+kHc+SIHeefX4TcXBs++UQcvMOHAzffDPz0E5CeLtskaz5NNDegag60ZEjMYQ1bW6ugU1NlGuPSpWvhcDjQpo0ddXVF+OknG44elQS9/HzJBJ4zBzj7bGD6dMCaz2paSKmpIkTMo6ElQ2IQa9g6Lk5KBc44Q77/5hugU6e+OHnyYnz77e/Ro4cNH30kmb7WMoPMTJkMsHWrzKOeOhX4+9/FGho9Ghg4UKyXaGtA1RxoyZCYwwxbFxcDX3whjt9Ro4CsrG3Yv78OW7acgYED/4Bt22wwDGDfPtkSjR3rskq2bgU6dQIWLAA++EDyac44Axg2THw7N9zgsl6iqQFVc6AlQ2KO/Hzg978HEhOBv/4VWLkSePJJ8cHY7Tdg584xeOEFoF07oHt34KWXRFRKS4Grrwb+/Gfg0CHAZgOOHxcfzpIlwH33yRYpPR34zW+AxYtd1otp+cQitGRIzJGdLRGlH34Ann0W+PvfDcTFOXDokB1Hj/4HBg+WKNLTT4twPPKIJNp16AD8138BX34pUaiLLpL3KigQa8VakvDb34plY1ovsZxHQ5EhMUlKCnDBBUBpqYEDBxxo186OYcOK0K6dDfv2Adu2SVlBTo44g198EThwQHJnfvtboK5OkvYAl+PYOpc6M1OO7dnDPBpul0hMEhcHGMb/ITPzSfTpI6UCb7xhQ2Ii8PXXIkJxcbKtOngQMAxXEWVJiWQCL1smVdtbtoiD19qAqm9fERdAXsdyHg1FhsQsWnfEtGkP4Jxz+qKy0oZly8TfMmqUJOT9+9/A4MHAqlXA7NnA3r0Smdq2Ta75j/8Qn05Jiat9ZmGhWDE33ggcPSrXTZ4cuwIDsGkViWJ8JcUZhoGFCw+hb98rsH69RH7MZlVmWUF1NTBkCPCHP0jBZFKSbJGOHROH8KWXAl99BezfL+fj4lzboRUrJHfmzjsbHs4WbUl7vppW0ZIhUYm12116unz/wgvAwIEGAAdOncrA7t2X4ayz4mEYEoIeOlRyZjZulLyYnBygY0exULZtk4hR+/bAkSNSVnDeebK1GjJE/DWlpeIcvuSSxq2XaJ8aaYUiQ8KWpv6l99and+tWccaOGwdccYWBd95x4Kef7OjS5QHU1MTjnHOAd94BMjKAX/wC6N0buPdeEZsVK4BJk4A//QlISxMn7rhx4ujdulXODxoE9OsnfpclS/wTiViYGmmF0SUSljQ1PX/uXHnAV68GKivFypg5E3jySXmYMzIMrFrlwIkTdmRnF+H4cRtGjJDWDiNHiuXy0ENirYwZI+/59tvAP/8J1NRIwt7IkeJvKSuTyNKQIbJVOngQTeqQFwtTI61QZEjYsWwZ8OijkmeyY4ckzTX0EM+dKzktd90FzJ8vkZ/vvwd+9jNJ+d+9G9i9ezcSE+3IzS1Cba0N48fLvQ6HOHS7dxe/yvHjIiyTJ8t5M+/lZz8D/vY3ScwrK5NQ9lVXAeeeK/cD/otELEyNtEKRIWFFcbFYIMOHA1OmuFsw3h7iigrglVeAoiK5PjcX+NWvZNuxYQOQkHAUhYXAl1+OxqlT8xEfb0Nysuthzs6WYsdHHpGSgDvvBP73f0UEcnNlm5abKyHp//s/+Yx+/US4hg+XsoQePVzr8UckYmFqpBWKDAkbTF/F0KHi18jNrW/BeD7EK1ZIdu7gwa5j+fkSiu7UyUD37hPRp88OpKYC69e3xV/+Iuc9H2bz9ZdfSri6uBh45hngppvkWGIi8O23sq6yMtkmPfIIMGCAVGN7vk9DWJP2rERr0h5FhoQNpq8iNdX9L73VgvF8iKuqALu9vmWQmGggM9OBw4d7wmY7A4BYHIcOSf6L9WG2PtyeW5nCQtm6ff890KuXfB0wQK677TZ57e19GsKatFdSAqxZI1+jNWmPIkPCBvMB9/aXPjMT2Lmz/kOckuIqUDSvLy42kJzsgM1mx8aNRSgrs2HlSuD22yUaNHasCMaRI/Ufbm9bmdGjZctWUyMlBjfeKD6a0aOl/umLL5ouEoWF4g/asUP8PDt2yOtoC18DFBkSRpgPuLe/9IsXSw9ez4c4P19yXK68Uq7/4x+/wvbtDiQk2LF6dRHuvtuG99+XaNOnn4pYXXKJ+Fi8Pdy+tjLXXiutHQ4cAGbNEmF74AH5OmuWbPOaKxLRng/LPBkSNpiTA666Sh7WnBzZIn38sWTkPvqoiI8VU5BKS0U8Nm9Ox48/jsKuXaMwZYoNc+a4rgNEtM48E6ePe2IVuPHjRZS2bhUL6vLLRfA2bHAPP5tzsptCLCXjsayAhBXFxSIY5gP+4YdixZx5pkR/fCXkLV++ARs22LF/fycoJU7ilSvri8E110jYua6u4QQ/M7Fvzx7xAeXny+vUVHFGe2LOyfYlXp7vPWOGezKeub4JE8QiikS/DMsKSERQWCh5Ls88Iyn7x49Ll7kJE3z/tTcMAyUlDlx66aWYN+8OACJWntbIggWSN3PhheL7ach6yM6u/6AvXOg+bdJKU5pS+ZOMF4ki4wuKDAkrzG1EXp4IzIQJkva/aZMIgWfqvWEYcDhkdOzvfve70+9j3W79618iLl26uFs3vlL5fZUzmD6j3Nz6625KfovnaFwr0dhBjyJDgkpD9Ui+RpVcd527EJh/7Q8fdgmMt9GxVmtk7lyxYEyBsa4jMVE65D3zTH1fSVmZWDq9e0vouqzM1dbBxAxdOxz+/QxaS6wiBUaXSNCw1hfV1spgNWs90ooVkoj39NPAG2+Is7eion6mb2YmUFVVi5KSEp8C44k1/8WzLio7W6yluXNdIldQINbTp59KucLx4yI0NTXS5qEl+S2xloxHS4YEBbO+qKjI1VXub3+TbFmz4fbHH0vv3MmTRWzi4lw+k7PPdm0jtm4FUlLiMXXqXLRv3x42m63Rim3TekhMrF8BvWWLiFhJiYhO1671K6U7dhTH7t//Ls5jw5D1JCc33VHbUAQrJpPxlFL9A78MEs1Y64sSEoBXX5W8k7/8RWqAfvxRRr7u3i0PXkGBJM79+9/SZqG0VKJMyclAWZmBN954CtdcU4tu3brBZrP5VbFtWg+vvebudDWthxtukNyZnTvluKdz1ixn6NoVmDgR6NlTHMFz5jRPFGIpGc8fS8YB4LpAL4REL2Z9UefO4vu4+GKpdE5JEX/L118Db70lmbTvvCPHzL/2t98u1djPPgv85jcGnnrKgW7d7MjIOA6gg9+9Wcz3e/BBeZDXrKlvPdjt0tIBqO+ctfpKWss56y2CFY3445MZq5Q6L+ArIVGLWV/0+ONSO3T22S6L4913gfJyEYhjx9wzfc8+WxLsPvwQ6NjRwDffONCnjx0LFxahQ4cOAJrWm6WwUOYmVVR4tx769xfrqbravbzA01cSjc7ZQOKPJfMhgHFKKaW13hjoBTWGUmoiANNllqa1XhDK9ZDGSUmRKYybN8tWadQoOZ6bK7U/DzwgM4y2bwcWLRLr49lngf/+b0m5HzDAQEKCAxkZ9Z28TQ0HT54sW6mHHqofIfrHP0TgJkwQcXv2WWnv8PbbLmunqZEk4kNklFKdtNaHAEBrfbnz2GilVJrWenkwF+ixrokAqrXWpc7XaUqp57TWk0K1JtI4+flioWRmil/m8GH5fv16SZC78krpLmdaEeXl4hu5/XaxGO67ryO6dTsbGRkz60WRmhoObszpWljoCm+3awf88Y/S7uHss8W6ilbnbCCpV1aglLoMgNZaf1TvYqXOB3CZ1vqxIK3P8/PXaa0HexzbrrW2N3QfywpCz7BhIiDHjslQtd27ZTJAQYHkr9x6K3DLLTJX+qefRHwqK3/Eq6+egW7dgKQkjU8/Vbj2Wpk1/cMP4nwdOVK2U01N0fdWNtCS64jvsgI3kVFK/Q7A2oa2RUqpzgBmAnjEtHaCgVIqCcABrbXyOL4OwL1a6zJf91JkQs/UqWJxbN8uW6cePaQ5d+fO8iAPGCDZuFOnygPdqZOB/fsX4PDhabjppkuQlSWCtG2bbGVuukmct8uXizXSqZNvy4QEB18i4+n43QEgraE30lof1FrPADDRm0M4gCHvNAA1Xo5Xo5E1k9BSXCw9V/btE3+M2cvlwAHJi9m3T/q8VFQAn30G1NYa2L/fgVOn0vDSS+dj4EDgiSdkC5WVJcWNQ4aIz+Rf/xJRycuLjXBwJOJtu5QK4PyGfC9KqeEAJgMYDaAMwAcAPtBaf66Umq+1ntnqC1UqD8BznlsjpdRSAOWeDmCn/2YiAPTt23fwt99+29pLIn5grTh++mlXQp7dLoPTOnWSJlLPPy+1Rd26GaitdUApO3btKsLAgTb8+KOroXj//mIN9e8vlg8gztytW4GP6m3wSTDx15KB1nqHVWCUUvOdX4crpf6slNoPEZY8AM8DOAgRnA1KqVoA0wP039AktNaLtNZDtNZDelg7PZOg4hlivvNO6e6/cqUMXTv3XPHXxMcDQ4bsgdYLcPSoHd98U4Q+fWz461+lsXhdnVg9Dz4o/poNG1xTC3JyJKGPhCf+hLAnOa2CLpDtypsAlmqtP7Re5PTVXA7g0VZfpYuuXo4lBfDzSAuxhpirqiQ/xowErVkj25sVKyQDePXqZNhs/4k+fXLw9ddSKtC2rWQD790LnHOOOF03bpTtktkSobxcEvZIeOJPMl4SJFfmcq11V631ZE+BAU77akoBlLb2Ip2shXdB6QpgfYA+k7QQa1KbZ/9cM8RcWbkGu3Z9jd27gZdfvhTFxTY8+aT4bFJTZXtUWCizp9etE6fuqFHiIN6+XZy/U6eG5D+P+IE/IrNIa/0bb8Lig+dasiBfaK1rAFQ6o0xWkhqKLJHQYq04tn5vJrVlZRmorn4U27e/hptucm2rRo+Waufvv5cRsgcPAmedJSUJhw9LMWV5uSTxXXVV/bacJHzwR2SalNuotd7RzLX4gwNOZy4AKKUGQfxDJEyxJr9t3Aicf774YIYNA+x2A6+95sCpU3Z89dW9yMlxv/eeeyTjNiVFap3y8mRCZF2dtIJQSuYrLV4cmv824h8R1+PX6R+qhGyd/CorYJ5M6LEmtSkFHD1qYN8+GVsybVoRHn7Yhh07RFisuS49e0pmcHo6MH166+TANNYWgjQPv5LxohWKTPBo7AE2z3/55eOIi6vCtGlFyMkRJ+/UqbIdOnTIlV3bu7er+bfWLc+6tXa+M/v8LlkSnVMCgg1FhiLTYhoTkMYeYDlfi3Hj4pGefgpffXUCy5Z1cDtvnVTQ2pm70TolIFzgtALSInzNCbrgAjm/aZP4XIqLpVTAFKNzzxWfiXTDM9C//+sYOfJBdOnSBRdd1AbXXuvq++LZ/Ls5XecaItamBIQLFJkopLV9DsuWyWC1oUMl83bgQKkjqqkB5s2TBLu6OhGSa66Re84/X7J5f/xRqqtnzjSQleVAdbUd27a1O+3k9XzAA9nIKdamBIQLbCQeZfjTirIxKiqkJ++UKZIkN20aMHy4vDbfb+5c4L33pETgs88knDxkiAjNiRMiLAMGSKlAfLyBzEwHOna04/LLi/Dggza39ZitLQONtznXJmxEFTgoMlFERQXw4otiCfzznyICAweKD6K01JWG3xBWkfrFL6QosbYW6N5dMnULCuT9SkrEsunUSVL8r78e+OYbeb1ggVhPq1cDlZXrcc45Dvz0kx179xahsNBWbz3BesBjbUpAuECRiSJmz67f3nLGDKkT8mxF6Q1rv1yzrcLdd8v3b7/tEoWuXV1Nt5cvFwvm+uuBtWuBESOkovrQIUn179AhFbW1v8QZZxRh7Vobqqvdt0jBfMCtOTstGWlCmgZ9MlFCRYVYFH/5iwyGB9wbat9wQ+NbEk/HqOnDGDgQWLVKWmc+6qxMs9ulAfg334hFs3Ch5L+8/rpk5/7449dYsyYdBw50wRVX3IXu3WXrZHakO3FCLB3DCO4DHmjnMqkPRSZKWLFCLApPITGthuXLgfMaaQfv6Rg1fRgFBdKBzmzunZkpdUVr1sj2qHt3sZy++UYe3L/9zcDZZzuQkjIWvXpdjxEjgHvvlfU5HLLW1asl7yUUD3isTAkIF7hdihKqqmQr4M3n0LOnbGUa25J4OkatPoysLBGaHTuA554TgUlOlkxcw5BI0/DhQPfuIjC1tXbs2nUNduwQP0/v3sDJk/K+t90m85eeeooPeyxASyZKSEkRB623JtmPPSYh5cYe6Px88eGYs55NH8Z110koevp0CVUbhvhdRo+WrdTs2dKjd9AgA+ee60BNjR0bNhThzDNtOH5c8mcuukhE6MEHgaNH6QOJJSgyUYIpECUl7j6Hzp1lSzNvXuP5M946+bdtK82/U1NlLrTpw1i4UKyb7GxxLL/++hG8+up/49AhO777rgjJyTYMHy7X9+snI0j+8Q9xCL/wAgUmlqDIRAmeAjFihCst/+abJY/FW8auZ82O6Ri1zj3KzZUWl95mS+fmmtnAHZCQUIT4+BSkptrw3XfigzHfOy5OxOiXv6TAxBr0yUQRvuYr5+S4h6at+S7e8mesc4+eeUZGlviaLb14sYHly99CSQnQs2c6xo+3YdMm4Mkn3d/bX78QiT5YIBkDzJ0r252CgvpbpjZtgFOnZHA80LQiwocfNmAYDrRrZ8fddz+Kxx9vg23bxF/zww+SS1NXJ0l9jz0mazCbf5Pow+9G4iT6qKoS/4m3koN335Uucyb+zpY2DANr1zoQF2dHly5FeOWVNkhOFh/Q44+LGNXVSebx5s0uvxCJPSgyMUBKCvDKK5JId8YZrvEiBQVineze7drWmILkDbPGyDAMOBwO2O12nHdeEX71KxsWLpSt1U03ubrg9eol7/X55yJq9MXEJhSZGOGFF6QhlLXI8U9/knKBm25yWSj+FBFWV1cjPT0dRUVFGDvW5pabY/qFNm+WFg/JyTJ+9rvv5LPnzvWvhopEDxSZKKeiQhLnzj9frAsz4zYuTnJWLrwQuOwyV6Zww0WEh5GfD/z617/Go48+CpvN5rUeaONGsV4eeAAYPFi2aS2pCieRDUPYUY7pY1m2TGYXAVJj1K6dzKNev14efLMK2luuzNatwJIlBhITn0JCQhbtIiAAAA4rSURBVBGATMTHx7s5kVNTJUnPWg8E1HciW+upcnK4hYoFKDJRTlWVpP7v3y9Jcc8/73rg339ftjf790tynYlnEWFCgoG2bR1ITbWjd+/eAMQSefFFESqtRbj27pWJAoC835Yt8pqd6GIbbpeinJQUsWImThSLxrqtWb5cwtfdutV/2LOzJax9440Gdu50oHdvOzp0KMI999gwdao4kRMSpKnVjBnydd8+CVW3aSNbI6UkeuVtaxSsRlUk9FBkopz8fEmCS052T9Z77jngww+ldWavXt7v3bp1KxwOBzp1smPnziJkZNhw881ioezZA1x9tSu5b+BA4MwzgT59JFqVmwv86lciat4S/tiJLnbgdilK8FWXlJ0tTt977nHNNDL9JzNnSm3SiRPe39Nut+PSS6/D++9fjZIS2+ltz8GDEil67z0Z0padLZ9dUCBjY//nf+Q6s57q6qvdt0ZmoypHk8YGkkiFlkwU0Fhf33nzJEnuyy/dyw2uvdZ7V7p169Zh//79iI+Px8GD4zBunK2eXyUx0T05z1t+jelELimR/jHsRBebUGQiHM+Wmd7qkrKzJRfmiy+k5cKIERJm9vawG4aBhx56CC+++CIA7+KRmyvlAT17uvwqKSnAunUiOtZxs4WF0odGa3eB4yC12IHbpQjH31lC/rSdtGbyTp48GYB7tbXJ5MnAJ5/IFIOsLLFQ2rSRvjI9esh5k+pqafHABlWxC0UmwmnKLKGG2k5aBaaoqAg2mw1A/UZW5vvceqsk23XsKFu0Y8fEWqmtFSvpyBH3CZAUmNiFIhPheLM0TPyN4NTW1uLNN9+sJzCA7+S8VaukN8zmzRJRstuB/v3l+LJlErFik24CsNVDxNNa850PHjyINm3auAmM5+esWCE+mORkaSheXMy50sQFZ2FHKb4sDX+2KYZhYPXq1Zg2bRo6d+7c6OdY32vuXM6VJv5BkYkCcnIksvOnP0mWbU6OKwdl7lzvPX2tPpiTJ0+ibdu2TfpMzpUm/kKRiXCKi6WlwhlnSCc6QELVTzwhEwa89fQdMMDdyduhQ4cmf25r+IJIbECRiWAqKqSXbnKyJNZlZcmDv3Ch5MisWgVccolca1Y/FxSswVtvOZCeXt/J2xS8RZ0AZvOS+lBkIphnn5W+MEuWuLdS2LxZ6odef90lMoBcc8UVXVFefi6KiqY3W2CAlvmCSGxBkYlgPvtMLBhP52tNjWT1lpe7jn3//ffo1asXLrwwA1u2FKEF+nIazpUm/kCRCRKNDVZrzftSUoANG1yvDcPAggULMHXqVHz//WWt6i/hXGnSGBSZICDDzxofrOYpKICk7Pu6z6wh+s//dLdm+vUDHn5YCiXvusvAjh3ig8nKuvD0JAFCggWT8QKMv8lyViHKygLKysSpe//9IiLe7gOAG26QyNJvfyt+keefF+GpqwOysgz06OHAyZN2pKQUYc8eG8aMYXEiCQxMxgsR/s4xMiupzeveew8oKpJpAmbPFs/75syRplOLF4uwHDwIbNok7TYnTNiLDz9cgPbt7UhIKMInn9jw6KNiBRESTCgyAcafpDVvQmTe17Fj/exZa7Kb1fn6zjvA0KEyhSA7uwcGD/4Dzj33XNhsNpSUSD8ZigwJNhSZAONP0tqmTTKX6J//dDl3zfu8Zc96JruZzteqKuDiiw3U1SUCOA8XXXTR6WuYhUtCBZtWBZiG5xjJ9xs3Sr6LtasdIPetX+8uKOZ93gbXt29voLTUgTfffBOevjZm4ZJQQUsmwDSUtHbhhRI9Ki4G5s8XX4p1LtF550nPljvvlOs8k92s0aj27Q1UVjpQV2fHlCmzoJQ6vQZm4ZJQwuhSkPBslZCfL69TU6VdZnGxlAKYQrR4MfDRRyI2gPt9ntEoQCyYujo7unYtwt69Nq9ZuIwqkUDiK7pEkQkhU6bIFsn011iF6NQpiRa98Ub9+zzD4k8//TR27tyJO+8swsSJNtxyizh5PYWJkEDCEHYY4ukUtmbPlpRI021vmNGoTp1OAWiDKVOm4MSJE0hMTMT48SIwc+YE5T+BkEah4zeENOYU9ubcBcQHAxi4/fbbsW/fPsTHxyMxMREAJzOS8IOWTAhpbiWzGUXKyLCfFhcTRpFIuEGRCTFNrWQ2DKlFqquz48473fvBMIpEwhGKTBjgbyXzxo0b4XBIsWNGhjh52cuFhDsUmQjCbrfjsssuw8033wybzYaLL2YvFxL+MIQdAWzevBnp6elo165dqJdCiE98hbAZXQpzDMPAfffdh9deey3USyGkWXC71Ao0t+tdY1jHlowdO7blb0hICKAl00KKiyX7NjXVvcCxuLhl7+trNjUhkQYtmRZQUVG/2ZS1wDEnp3kWzbFjx/DnP/+ZAkOigogRGaVUHoBJAOYDqAEwBkCN1npRqNbkT9e75ohM+/btMXfuXHTv3p0CQyKeSNouJQFIA7DO+a9bKAUGEB9MVpb3c81J7zcMA284KyL79etHgSFRQSSJDLTWg7XWSmvdRWt9b6jXYxY4eqOp6f2mD6a8vBwnT55snQUSEgZElMiEG80tcPTE08nbtm3b1l8sISEiYnwywGm/TBLEJzNIa72ggWsnApgIAH3NSfStTGuMamUUiUQ7kSQy6wFAa10JAEqpaqXUB1rry71d7PTXLAIk4zdQi2rpqNYjR44gIyMDDzzwAAWGRCURXVaglDoAYLApPL4Ix7KCQ4cOoVOnTgCA2tpaxMfHh3hFhLSMsCkrUErlKaU+8PNfWiNvVwkgLxjrbk0Mw8Ctt96KTZs2AQAFhkQ1Qd8uaa3LAJQ15R6n2KzTWncJzKqCh9UHk5qaGurlEBJwIim6NN/LsTQ0UbBCCZ28JBaJCJFx+lxqrMeUUmMAvNmYPyZc2L59OwWGxCQRE13SWi9yhqUBCWNDaz0phEtqEqmpqSgoKMCIESMoMCSmiBiRAU6HpSOKtWvX4swzz0RycjJGjRoV6uUQEnQiYrsUqRiGgYcffhgvvfRSqJdCSMigyAQIq5P3tttuC/VyCAkZFJkAwCgSIS4oMq1MXV0dVq5cSYEhxElEOX7DHa014uLicP/990MpRYEhBLRkWg3DMDBv3jwcP34cHTt2pMAQ4oQi0wqYPphDhw7h1KlToV4OIWEFRaaF0MlLSMNQZFrAZ599RoEhpBHo+G0BPXv2xODBg3HXXXdRYAjxAS2ZZrBr1y5ordGvXz/cf//9FBhCGoAi00QMw8Add9yB9957L9RLISQioMg0AauT95e//GWol0NIRECR8RNGkQhpHhQZPzhw4AAee+wxCgwhzYDRJT/o0qULZs6cibPOOosCQ0gTocg0gGEYUErhggsuwODBg0O9HEIiEm6XfGD6YFasWIFInk1FSKihyHjB6uSdPXs2lFKhXhIhEQtFxgNGkQhpXSgyHnz55ZcUGEJaETp+nZw8eRJt27ZFYWEhTpw4gYSEhFAviZCogJYMZIs0ZcoUVFVVQSlFgSGkFYl5kTF9MElJSejUqVOol0NI1BHTIkMnLyGBJ2ZF5osvvqDAEBIEYtbxm56ejiuvvBITJkygwBASQGLOkqmoqMCxY8fQoUMHTJo0iQJDSICJKZExDANz5sxBSUlJqJdCSMwQMyJjdfJef/31oV4OITGDioXiv5///Of6rLPOopOXkACilFqntR7ieTwmLJm9e/dSYAgJETFhySil9gL4NtTraITuAPaFehFhBH8e7kTCz6Of1rqH58GYEJlIQCm11pupGavw5+FOJP88YmK7RAgJHRQZQkhAociED4tCvYAwgz8PdyL250GfDCEkoNCSIYQEFIoMISSgxGwVNglPlFITAVQ7X6ZprReEcj2hRCmVB2ASgPkAagCMAVCjtY4o/wwtmTBEKZWnlFqqlBqklEpTSk13PnxRjSkwWutSrXUpgFKl1HOhXlcISQKQBmCd81+3SBMYgJZMuGL95aoBsChG/qJP0lqfHtWpta50/jWPWaw/j0iFIhOmRMMvV1NQSiUBGOTlVI1SKk9rXRbsNZHWgSJDwoU0iNXmSbXzXEzitOSSID+bQZFo0VJkwpRo+OVqIl3hcvhaqYH8HGKR9YBsGwFAKVWtlPpAa315aJfVNCgy4UlU/HKRlmH+/7e8Xq+UGqKUSvM8F84wuhSGaK0rrb9EWuv1AIYopaJ929DVy7FYtWJ8UQkgopzhtGQCiHPLc6+fl09q5K+T+csVcSFMP1kL74LSFU7LLpZw/kFZp7XuEuq1tBSKTABxRkSaFBWJpl+upqC1rlFKVSqlkrTWVgdwUgxHluZ7OZaGJv5OhRpul8KTqPjlagYOAKeTDpVSgxD9/81ecVq1btE2pdQYAG9Gkj8GYBV2WKKUmmjN7HT+cl2utZ4UwmUFBWfWbyWcCYkxEFVrEEumdxIAROLPgyITpkTDLxchAEWGEBJg6JMhhAQUigwhJKBQZAghAYUiQwgJKBQZEnKcjbmeU0pppdQBpZTDx3UfOK9ZGgMlFlEDo0skbFBKLQWQ5yvb2Sk+bzhruUiEQEuGhBPPAUjy1g3PmZBIgYlAKDIkbHDWKNUAGGs97iwvAAUmMqHIkHDjTQC/MV8423LmORuLkwiEIkPCjaVw3zLNZElFZEPHLwk7lFIHINXX5ZBJDd56/5IIgZYMCUfehAwyK6PARD4UGRKOdIUIDB29UQBFhoQjgwB8EOpFkNaBIkPCCmc0KRa6AMYMFBkSbuQBzImJJigyJCxw1i85IH1+oZRysD4pOmAImxASUGjJEEICCkWGEBJQKDKEkIBCkSGEBBSKDCEkoFBkCCEBhSJDCAkoFBlCSEChyBBCAgpFhhASUP4fwdPWXSRafTsAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARkAAAElCAYAAADQqJbBAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO2de3gTZfbHv69QoA1CablVEKSBgkAFC2VVEBetK0h3kZty1QWVizfWGwK7ICmKxp+4q+uVdctiARG5eQUFUVfUQUoFqq5UbMUFW4GWIizSAn1/f5xMM00maZJmkkxyPs/Tp+0kM/O2dL6c23uOkFKCYRjGKM4L9wIYholuWGQYhjEUFhmGYQyFRYZhGENhkWEYxlBYZBiGMRQWGcZwhBBbhBCJ4V4HEx5YZBhDEUKMAZAFICnE900UQkwTQrweyvsy7rDIMEZzk+NzaqhuKITIAAlbRSjvy+jDIsMYhhBiNoCHHN+GzF2SUhZIKdcCKA7VPRnPsMgwhiCESAVQKaVUH3S2KGIUFhnGKKZLKZc6vq4EYA3nYpjw0TjcC2CiDyFEFoCdmkPF8CHwK4R4ydd7SCmnB7A0JgywyDBGMNZFBHwKwLJwRCcsMkxQEUJMA5AkhLBrDnM8JoZhkWGChqPgziqlHOtyvByAXf8sJtphkWGCyTQp5UM6xyt9OdmPwrkKdq3MA4sMExQclb2e6lKKHe9J1aS03XC1gIJASKuMGX04hc00CCFEqsMCeR2Am3XhiNGo1s1LjsxTKNakFgJmCCHsjnUwYUBwj1+GYYyELRmGYQyFRYZhGENhkWEYxlBYZBiGMZSYSGG3bt1aXnTRReFeBsNEJf/73/9QWVmJsrKyo1LKNq6vx4TIXHTRRcjPzw/3Mhgm6lAUBXa7HVarFUuWLDmg9x52lxiGCQitwNhsNo/vY5FhGMZv9uzZU0dgLBaLx/fGhLvEMExw6datG6677jpMnjzZq8AAbMkwDOMHX331FU6fPo2EhATMmDGjXoEBWGQYhvERRVEwf/58rFixwq/zWGQYhqkXbZB3/Pjxfp3LIsMwjFdcs0i+uEhaWGQYhvFIdXU1li5dGrDAAJxdYhjGC02aNMGiRYuQmJgYkMAAbMkwDKODoih45ZVXIKVEhw4dAhYYgC0ZhmFc0MZgqqur0bRp0wZdjy0ZhmFqcQ3yNlRgALZkGIZxUF8WqbAQ2LABKC0FUlKAkSOB9PT6r8uWDMMwAIBz586he/fuugKTmwvMmQN06QJMmUKf58yh4/XBlgzDxDiVlZVITEzEwIEDcfnll+O88+raHoWFwLp1QF4ekOQYMjNgADB8ODBpEpCZ6d2iYUuGYWIYRVFw++23Y/fu3QDgJjAAuUjjxjkFRiUpCRg/nl73BosMw8Qoagymc+fO6Natm8f3lZYC3bvrv5aWBpSVeb8PiwzDxCD+bBVISQH27dN/ragIaN/e+71YZBgmxvjhhx/82os0ciSwejVQUVH3eEUF8Oqr9Lo3OPDLMDFG586dceutt2LIkCE+VfKmpwOjR1OQd/x4cpGKikhgxoypP40dE2Nq+/fvL7mROBPr7Ny5EykpKejYsWNA56t1MmVl5CK51skIIXZJKfu7nseWDMPEAGoMJjMzE/PmzQvoGunpvhXfucIxGYaJcrRB3lmzZoX8/iwyDBPFNLThVDBgkWGYKEVKiXfeeSesAgNwTIZhohIpJYQQmDdvHmpqasImMABbMgwTdahTBU6dOoX4+PiwCgzAIsMwUYUagzl9+jQipTyFRYZhooRICPLqwSLDMFHAzp07I1JgAA78MkxUcMEFFyAzMxOzZs2KKIEB2JJhGFNTUlJSO1Fg3rx5EScwAIsMw5gWRVFw33334c033wz3UrzCIsMwJkQb5M3Kygr3crxiKpERQmQJIV4XQmQIIVKFELOFENPCvS6GCSWRmkXyhKlEBkAigFQAuxwfyVLKpeFdEsOEjuPHj+Opp54yjcAAJswuSSn7hXsNDBMuWrZsiT//+c/o2rWrKQQGMKHIMEwsoigKzpw5gyuvvBJ9+vQJ93L8wnQiI4TIArlNlQAypJRPhHlJDGMoagwmLS0NAwcO1B1bEsmYTWQKAEBKWQwAQogKIcQWKeW1rm90BISnAUCnTp1CukiGCRbaIO+CBQtMJzCAyQK/UspiVWAc3xcA6C+ESNV571IpZX8pZf82bdqEdJ0MEwzMlkXyhKlExgPFACK7UIBhAmD//v2mFxjARO6Sw1rZJaVsFe61MIyRVFVVoWnTppg4cSLGjh2Lpk2bhntJDcJslsxjOsdSAWwN9UIYxggURcGMGTNw8OBBCCFMLzCAiUTGEYup1B4TQowBsEYbp2EYs6LGYJKTk9GqVfQY7KZxlwAK5mq2ESQ6jk0P45IYJihES5BXD1OJDEBCE+41MEww+eqrr6JWYAATigzDRBvdunVDdnY2xo0bF3UCA7DIMExQUedFl5YCat1cTQ2QkuI+O3rPnj21e5BuvfXW8Cw4BJgm8MswkU5uLjBnDtClC5CcDHz6KX0kJ9OxOXPoPQDFYBYuXIhXXnklvIsOAWzJMEwQKCwE1q0D8vKAQ4eA1auBbdvotUmTgJtuAoYPp68tFgWrV1MM5uabbw7vwkMAWzIMEwQ2bADGjQOSkup+nZQEjB9Px5KSgMGDFaxaFb1BXj1YZBgmCJSWAt27u38NAGlpQFkZcObMGezdm4uqKitOnLBhyRILCgvDs95QwiLDMEEgJQXYt8/9awAoKgLatwfy8uKwb98iJCTYcOutFrc4TbTCIsMwQWDkSIrDVFTU/bqiAli9WsHJky9h9eoatGrVDjabBQMGAJMnUwxn7VpEtUXDgV+GCQLp6cDo0RTYHT8euPRSYMgQID5eQdu2duzbZ8Xhw9W4555mddLY2piN9rivaFPmemnySIBFhokJQvEwTp0KZGbSfSoqgH79FBw9akeTJlb8+KMNs2c3w4QJ7uelpQGff+7//XJzKaM1bhwwdCi5aHPmkNhNndrwnydYsLvERD3a+pUpU9xrVoJJejqwYAFw880Kjh2zIy3NiuXLbcjOtuDcOf1z1JiNP2hT5pMnI6LdLxYZJqoJ18MYFxeHiy++uDZNrY3TaKmoAF59lSwrf9CmybVo3a9Igd0lJqrx5WEMpttUXl6O5ORk9OvXDxkZGRBCAHCP2aSlkQXz6qvAmDH+r6G0lFwkPQJ1v4yCLRkmqnGtWdGi1q8EC0VRMG3aNOzcuRMAagVGZepUwG4HSkqA5cvps90eWPzENU2uJRD3y0jYkmGiGvVhHDDA/bVgPozafjA9e/b0+L709OBYTiNHUlxp+PC6VprqftntDb9HsGCRYaIaox5GbbaqqoqySHFxViQk2FBcbDE8jRxs98tIWGSYqMaIh1GbOv7ll4M4eNCO6morevWyoVs3S8jSyNqU+eefk1Vmt0eWwAAsMkwMEMyH0X23dUfcdttMZGYOxLRpFkyY4NxtnZlp/AMfLPfLSISUMtxrMJz+/fvL/Pz8cC+DiQACKcrTnrNvH5CVBWRlfYH161ujV69UTJ5M78vLo2DuggXA4sXA1q0UdI7UStxgI4TYJaXs73qcs0tMzOBvUV5hIXDDDUB2NrBnD3DllYAQwPvvK3j00cfw/fcrdHdb5+YC77xD7zW6+M8MsLvExARaN0cNAA8Y4Nm1yc0Fli0DTpwAXn6ZxGPlSqB5cwWAHTU1VrRufX+dzFVREQnLunV0zSNH6DVv94kF2JJhYgJ/KmRVQRo0CLj/forhlJQAiYkKamrsqKqy4tprbTjvPEud3davvkrnZ2cDb71Vt4o3EitxQwWLDBMT+FOUpwpSZaVz0+FFF0mkpGxFkyZW7N1rwzPPWPDtt87d1kOGABkZwLffkrWkl7kKdvGfWWB3iYkJ/CnKU0v2FQV4+21g69YatG59HqqrZ0NRzqJlywRs306uVEUFWTwAfS0lWTJ66etIq8QNFWzJMDGBPxsUtSX78fEKHn30IZw4cQJNmjTBgQMJaNsWiI8HUlNJNGpqgHbtgJkzgaefpgkFwdoIGQ2wJcPEBP4U5alVwomJCtq0sWP/fitef/08dOwIPPkk0KIFcOGFVHMzfLh7LxezVOKGCq6TYWIKtealrIysEE/1K48+qkBRqOFUdrYNb79tQX4+xWAmTwYWLSJxefxx5zVLSoDt2+lYWppv94kmPNXJsCXDxBS+VMjm5+cjP9+ONm2s+OwzG1q2tKBPH2DhQjp37lyKvUyc6N6drmlTen3OHCrKY1hkGMaNTp064YorrsAdd9yB116zYO1ayh6dOkWZo5UrAXUmm2vtjZRAdTU1xIrFmhg9WGSYqCPQfr6bN+/Hjh2pKCtri5SUB1FcrL/v6fe/p3S4Xu1NURFw0UUkSsFuiGVWWGSYqCLQ5tpqDOaSSyZiypQxbudpxaKwkF5r27Zudzpt+4hTpyKrO104YZFhogZ/tw6ovPYaCUxcnBVbtgzD++/TeXPn0kZH1/PUTNVjj1EMRkr3DFJeXmzWxOjBIsNEDYH081UUBStX2nHqlBUVFTaMGEGzqTduBD75hMRG77ypU4GWLUmIqqvJRVLbR0Rid7pwwiLDRA3+Ntc+efIkliz5G06etOLUKRs++MBSK1D33APceCOQn0+VvZ6u2bMntXTIzASsVmD37tiuidGDK36ZqMFbc+0PPqB9RTNnAjk55Fo1b94cXbo8jMOHbejb11LHAkpKAm65BYiLo93UrqhtI0aPproYAJg3j9y1QJuDRytsyTBRg7af76FDzgzTgQO0D+mBB6jh1LZtCubO/QWjRv0OJ05cjObNgXffdfZ+UbNR7doB//0vcNVVde+jF/uZMIHcpEmTQv9zRzosMkzUoAZkr7oKaNSIGk4BVIXbsaO6z0jBjh12pKV1w9q11+D06UbYvx9o3hzYtIm2C6xfD/TvD+zcCXToAPTqVfc+oZ7lZHbYXWKiisxMoFUris0cPuxMYW/fDqxbp+Dxx2lsyaJFD2PQoEbYswc4d45EJicH6NQJKC8H3niDskYJCe6bGkM5yykaYJFhoooNG4Dbb6eS//btgZ9/Bn74AdiyRUHjxnbEx1trR8cWFFAnu5wc+vzcc2SBjB8PnDxJLteUKe5WiZkGq0UC7C4xEY2/1bulpZRSXr2aXJqBA4HzzgPy8g6heXMrWrYkgQGAggK61t13A7/9LfDiixSbAYDevYFff9UP4JppsFokwJYME7GoGZzGjcml+fhjirnk5Hg+57zzqNFUXh7Qty/QqNGv2LQJ6Np1NAoKHkOjRiQwFRXkFvXoQeelp5Ml88UX9DF8OLlKemjbRuTlATt20OdJkzh1rQdbMkxEomZwhg4FNm8mq+S224Bdu4CHH6b3eNvlvHIlsHmzgoSE5zBmTA7eeKMLDhyIw5o15Op89x3FYTZsINdK7W4HkABt2FD3mCtmGawWCbDIMBHJhg3k6mze7L5NAAD++ld916mmht7z7LMKevSwIznZimPH2qKigqpyDxyggrnkZMpCbdpEWagbbqBYTlERsHw5BX1nzPC+RjMMVosEWGSYiKS0lFwkvVRxRgalmvVSxSkpwJdfksAkJNBWgd27LXjzTeDf/waeeQb44x+B2bPJvbn9dhKat96iAHHLlkBVFfCnP9H1cnL8383N1IVjMkxEkpICfP+9fqq4qIhK+PVSxRkZ/8HZs3a0bWvFyy/b0Lu3BfPnU/n/smUkFL/8QsI1aBDw0UckWMnJdL0+fYDnn6dr+TMIjvEMWzJMRDJyJLBiBcVgtBMG1AzOoEHA2bPu5113XVcsXToKO3eOwvr1FuzeTS7UpElUkNeuHQWH1SmPzZpRJfC775IIdexI19Hbza0KzubNlH1iy8Y3uMcvEzL8TUfn5ADPPgvYbOQiqe0Uhg4FXn8duOQSEpCUFKB37y+RlWVFixYtUFgI3HkncMUVwLZtZKksXAisWkVismABicyllwKffUbW0uHDQGIira1LF+Cyy1A74xpw9qm58EK655VXUpq8vj41sQTPwmbCir9zqAESod696X3jx1MV7qBBwD/+ARw/TkIwZQrQpImCZctyMG/esloha9SItgdcdhm5Qdu308evvwJ/+xtll/76V7r+qFHA/v107IILKIWtddO0e5WmTKGU+uTJ9P3atfQ64xl2lxjDCaSZlGo5TJlCqet162gvUUUFbRvYuJGupSi0F8lqtWLbttuwdy8Fc4cOpRYM//oXuUeqoJSWkgtWU0MFeP/+N8Vl/vQnyjBlZ5MQaQfBafcqaSt6ea+Sb7AlwwREYSG5M9rWCZ7wZw61em1VlCZPph3O69YBX35Ju6KHDnUKjN1OAjNxIs2mvuIKOmfAAGq9oCg0hG3GDOCrr4DBg0mkrFY6PmwYZZQyM+mal1xCRXjaQXDqXiW9AW28V6l+WGQYv/HX9fF3Q6E3UerUiVLN586dQ15eHqxW2ou0aZMFN9xAmSPXc8aPpyrg+fNpjQMHUvzlkUeohea0aXTPigpg7176ebQVvefO0Xl6Fb28V6l+2F1i/CIQ18efOdSA9w53Viulths1aoScnBw0a9YMFosFpaWUOdJ74NPSaKOkKnS9e9MUSHXKY3U1bVlQFEp1d+pUt6L3+HGyonJz3SuDea9S/dRryQghLjJ+GYxZ8Nf1AfybQw143+XcooWCM2eextGj55CcnFy72TExke6tN2u6qIgESL3myJEkZPPm0dTHVauownfePOCnn5zXSE+nTNRrr5Gltngx71UKBF8sGTuAm4xeCGMO/O2jC/g3hxrwvMt561YF+/bZ0aGDFTffXIVBgxLwww9k2Xz3HdW8dOhQ91qqkN15J1kiw4c717N4Mc1Q+vVXCvguXuxZNHivUuD4IjJjhRCLpZR7DF8NE/H46/qoqA/piy8Cf/87WQ4DBtAxV/RE6cMPFXz6qR1JSVY8/7wNS5Yk4F//oroVqxW45hqyln77W+DBB92FbPRocnvUa/bqRbGgp56iAryzZ+sXDd6rFBi+iMwHAMYJIYSUcrfRC6oPIcQ0AKrhnSqlfCKc64k1GtJLZedOCtrefTfFR7wNXlNHjjz3HHDihIL27Wk29d//bkNxsQU7dlAMRbuGGTOAESPouJ61oWeNrF/PwmE0uiIjhGghpfwFAKSU1zqOjRZCpEop14dygS7rmgagQkq51vF9qhDiJSnl9HCtKdbw1/VR8TdgrK2TadasOd57rxdKS+fitdcsOHjQc1xo2jSKs3hqA8HWSOhxExkhxDUAJIBt2uNSynVCiEuFEA9IKZ8M1QJdmC6l7KdZU7EQIitMa4lZ/IlPqBW4b71FWZtDh+A2ekRb0FZYSC7V228Do0cfRt++bZGe3htjxizCsWMCkyZRHxh/40K+EugcbcYzdbJLQojbAJRLKbfpvVlK+SWAfwghHhdCtAjFAjVrSwSQofNSJQtN6FEzL88/T5/1HkRtPU2nTvQevXoatVZGfX9ZGTBkiIL9+2dg4sTtuOEG4KuvRK0gHT1qTI/dQLY+MPXjmsIuAZDq7QQp5XEp5RwA04QQfVxfNzDlnQqgUud4BepZMxN6XKt2+/ShAK3efp+iImrkrb7/l18UlJfbkZSUisWLL0VJiTM7lJZGA9cWLSIh0FYbe0qJB7LeAQN4f1KwqCMyUsoPAHwphBjl7SQhxNUABgAoEEK8J4R4QCM4RsVHkuAM+GqpBJCos8ZpQoh8IUT+Eb0RgIxhFBYCs2YBp0/TLurCQmetDFC3nkYVBoDiLFu2KGjWzI7ERCuee86G7GwLHniAdlSvXUtuWUkJbQ/Yt4/qWmbMoP1NDalbCaT+h/ENt5iMlLIEZNEAAIQQj0kp5zqEZSyAGwG0Aj3c/wA9/DMAPCGEUPtGzDV64fUhpVwKYClArR7CvJyYQQ3YCkEPfePGwB13UFOoU6eovULbttTEe88earEwZQplntq2LcPzzz+B88+njnZVVRZYLM44S0YGbXR87z2qvFXjJ2rF7uOPU1A6EAKp/2F8w5cU9nRHVkcVljUAXndYPbUIIVoCuBbA40FfpZMknWNuVgwTHrQux7PPksCcOQPEx5P18ZvfkOVRVAQ4CnUhBH1OSQEOH26PxMR7cOutmdi/31Kbwfr+e9oWoAaP1dJ+baYoLw/4+uvARSbQ+h+mfnzZIJkIqpW5VkqZJKWc4SowQG2sZi2AtcFepIN86AtKEoACg+7J+IHW5ejVC/jLX2iyQEYGNe3+6CNqGpWZSfGZhQupZcO6dTvQps23WL0auOCC3+LAAQumTiXX6JtvqGNd+/bU6+XKK/Xv3dDd0P5ufWB8xxdLZqmUsp6+7XV4KdDFeENKWSmEKBZCJEoptQHgRCnlViPuyfiH6nKoLlPr1rQrurLS2Ujq7rupdUNcHB373e9osuOnn6Zj9OgcLFtGjb3btQPy8ykOk5pKGxrLyqjYTovqMn38MVURFxYGFpMJtP6HqZ96228KIbo44jRhx+G2JapVvkKIDFDtjNdgM7ffDIxA2mU2bgx8+im5L3/+M1ke+/cDL7xAMZncXHqYd+wAXnxRQWWlHe3aWQHY8OKLFhQWUkuGTz4hV6pLF7r2jz/S3qSWLckiSkpyill2Nt0vO5vu3ZCWmOrPXFZG1hPXyfiOp/ab9VoykSIwAAVzHVmjLJDrlMrVvsagPsDjxpF14m0LgMrIkfT6vfeSCKSkUC+Wu+4iERg0yBk3+fBDBRUVdnTrZkXPnjYcPEhBmvR0Sk+PHk3ntmhBrtXIkcCHHwL/+Q9ZM8OGUSPwSZMoVnPbbbSuigrPLSd8gSuCgw83EmfcKCwkQdFuAQCcD7C3jYRDhtDcopkzqe3l//0fuU2DB1NWqE0b4OKLgU8+WYL27Usxc6YNS5ZY6lzzzjvJItm2zf3+Q4bQ+YcPk6Vz1VXu1kZenvetBYwxBGzJMLGHLzUj3kSmcWN6yMvKaGpjeTm5PadOAZs2nUNBQSNcffUsNGlSjfHjE2rbL6h88QX129W7/8iRFDzu149S33rZIE45RxbcfpNxw992mVpGjiQr5K67aMvBxo30+aOPgPJyBddddx+uv/4YEhIao2PHBKxfT+lpfytqvTW24pRzZMEiw7jRkAdYm6VRu8jt3g0UFCjo29eOtLQ4LFnSpHbP06BB7hW1AwaQOOmlkzdsoHgLp5zNA4sM40ZDHuDCQuDgQdopvWwZNYUqKFBw4YXU0c5ms9W2zFRxtY5mzKCRJTfd5BSqxYtJXA4edL7PVcy4JWZkwjEZxo1Aa0b0MlKrVxeguNgOi8WKzEx3gQHcraP0dJqDtGwZWTQ//EAtIlq2pIBy9+7OTJfd7l/LCW7hEHo4uxSj+PLQ+VMz4ikjVVx8DPfc8y9MnDgNK1ZYvGasgLpr6tWLMkxvv01iN3Gi8/6+ZLpUtOKnduTjEbPBh7NLTC3+1sBo/x/yJE6uGalvv/0WXbt2RWpqK9x007347jvv1tHOne5rys2lyuBHHqk7lxrwfXpjICNcmODCIhNj+PrQ6QnRxIn0/gcfdBcn7S5mdbLj2LFjMWHChNqU8oIF+h31AHcrSF1T377AjTfq/yy+pKobko5nggOLTIyh99Cp1klVFfWBufNOdyGKjydRqKmhBz89nYRA7SLXtCkVyN14o4LVq2l07AjHRiNtzEWvojYnx7MQZGbSWiZMcP9ZfElVcwuH8MMiE2XUF2txfei0FkvHjsDKlcDcuVShq33oN2wgl+W//yUh6t6drlVRQeNIjh4FDh1SsGIFjS1Rs0i+TDHwJgSjRtHQtYoK/6cjANzCIRLgFHYU4UuPWm0NjGvLybg4KtMfOJBaLGgL5EpL6bx33qFy/iuvpKxP06a0baB161O44IK/Q0or9uyxYeNGi88pZW91OTU11Ogq0FQ119NEAFLKqP/o16+fjHb27pXy+uulLC+ve7y8XMphw+h11/fZbFK+8or7+2w2KadPp8/qOZmZUiYnS5mdLeUddzjPLS+X8pJL6Nh3330nv/jipBw1it5vsznv29C1q+uaOdP366r88590nVdekVJR6POwYXScCR4A8qXO88fuUpTga4BTWwNTVUUuUl6eew3MjBnU1lJ1p1q0oDqV4mKaxPjNN+TiFBQoOP/8wwD+gK5duwIAZs8Gli/3vEFRz6XzpS7HW/rcm4vII2bDC4tMlOBPgFN96GbNohjMVVfVfejS02mn85YtQEEB8MQTzthHo0YkBDYbcOedNNmxRQsrjhy5Huqfk7dYh7f0ua+Fdb5eT5uO5xYO4YNFJkrwN8CZng48/TQ9kHfd5R5U/ekn4PLLqUfvxo1UK3PrrTSgfv58ID5eQbt2dpw7Z8XRozZ8/XVj5ObS7mlPAVlf0uf+tGfgGhhzwIHfKCGQAKfeZkY1qHrZZdTAW0oSKLXHy/HjQGmpgt697aiqsqK83IayMgvOP5/6+Y4Y4R6QLSykNPXUqdTd7tChuusIdOwIjzExBywyUYI3wfCWhVEbdpeUUBylpIQEZscOmgxwySWUpcrNpVEmEycC1dUVSEjoivPPt6FJEwsyMsiCOHMG6NChrpvi7xRJf2hISwomdLC7FEUEGuDUxiu0e5AOHaKv77lHdUH+h7g4CyyW66Eo1yE7uxFsNue5CxdSQykVV3empITE5p573N2ZQGpWuAbGHLDIRBkNDXBqXZCkJKd1NHiwgqZNn0ZVlQ0//5yGhx9uVO/mQld3ZuRIEq3hw+tmvHwtrHNFe71ACvWY0MAiE4N4S/m6ZqmmTgUsFgWrVtlx9qwVx493wIkTZClox4+oDaXUwWt619K6dAMG0Pl66XNf4TEm5oBFJsaoL+Xr6oIoCu1FSky04osvbGjd2oLycrrOhg2U/i4vpxjOmTPA1Vc776Xnzqgu3cKFtEWhpKRhNStcAxP5sMjEEL6kfFUXpEsX4M03i/Dtt3Y0aWLF55/bUFFhwX33UUPwNWsozb1sGVkQrVsD/fuT+Bw/Ttf6+Wfgn/+kMSjaXjAdOgC//kqvBUMMuAYmsmGRiSFcYyRatyk+HnjxReC55yiLNGoUMHKkFT163IRNm36P4mILbryRNivWWrsAAA6uSURBVCsCtJ/ooYfITbJYnIJRUUHWTatWwO230+e33gI2b6bha927szsTa7DImJBAW0lqYySubtO775JV0q4dcPjwLixffhF27kxGWdk4VFTQBkp1qkB6OlknDz9MtTTvvuu8/6FDVBV8xRXORlM33QSsWuWc8sjuTGzBImMyApnsqKLGSOLj3d2mffvUjI+CTp3s+OabQViw4H4A1LWuc2egWzdnm4f8fJoU6cqGDVT1e/iw81h6OvDYY0DPnhSDYYGJLbgYz0SsWwc8/ji5MyUlJBaTJ5NYrF1b/+witSp41aq6bpOa8m3XTkHHjnacOmXF8eMzaq83YAC5Q2+8QW0epkyhwrrFi+l4ZqbzHqWl9FmvRoUL5GITFhmTkJtLzaSuvpo69mt7xfhaRq+mfFeuJDdHWxXcqpWCbdvsqKmx4vhxG7p1s2DOHNoOUF5OrtL33wP/+x/w8ss0mqS4mJpVaTNKiYm0Dr1tDFwgF5uwu2QC1KzQwIHkEg0Y4J4V8rWV5NSpwK5ddM3Dh+mhnzLlHF55ZQ2sViuSkmxYscKC4cOBykpg0SLaWnDNNTQ5oKCArJEuXaj9w/HjZNF07Egisn077Xfq0KHufblALnZhkTEBalaopKRu3YnWgunSxXcrYcYMsoIeeYSukZPTCCNGPIyhQxtj/HgLOnak+AlAQvbxx8D551Nc5vrr6X49ewLPPAN89hllmZ56ikaYPP88xXC4QI5RYZExAWpWqG9f9zL6tDTgww8BRfHdSlDdpsmTFXTu/DFKS+9DTU1LDBtGMZcHH6T4y4EDwKZN1NwqNRW45RaKAfXr5wzg7t5NVk6vXs42DenpXCDHOGGRMQFqVmjyZPcy+mXLyEWZM8e/h7hnTwVvvmnHsWNW/PTTGZw4EYdWrah3TFISWSizZ5Nr9fjjlLJWLSXVNVNdoIkTqQ2EFi6QY1RYZEyAdiOgtoz+ww+BTz4hERg92vn++upo1LlIXbvSVIHi4gSMHk0padVCSkkhgXnrLWD6dAoQP/AAfXz/PQWC1TYSNTUc0GU8w9klE+DaK+bUKYrBlJVRxkkrMPVNLNixYwfsdpqLpI4tSU+nQG1enrMXTePGVGw3dCi5Q5dcQj1+d+0i60l1gdROeNz1n/EEz8I2EfXNpvY0j1o7N7pZs++wYsUKzJ49GxaLpfY9OTkkLGfPOq8PUKynUyeyVmpqgDffBP7wB9oyoA3o8kxpxtMsbBaZKCInhywX17nRAPDCCz/hyJELvE4Q0BOo7dvJIsrIoOBur14Un/EkdEzs4klkOCYTQgLdc+QrniYWKIqC9957Ai1b3gngGt1zvfVmmTu3rqWidc8Ypj44JhMCCgspdpGdDezZQ9MX9aY7NhS9SYxqkDcxMRWdO1/m9Xy9fr92O7tCTMNgS8ZgcnMpUHriBJXjl5VRWf/o0c6S/mCN7nBtR6kKTOfOVhQV2XD//ZZ6r8GpZybYsMgYiLodYNAgqpC99lo6rt0OoO1168v1vLlbWpdnxIgj2Lz5CbRsSQIzZoyFxYMJC+wuGYi6HaCysu7oDu12AF93JteXmlZRXZ6ff26DxMQH0amTDXa7hV0eJmywJWMgaiDWdc8R4Kya9WVnsj+TEhVFQXx8PBYs6APgckN+LobxB7ZkDEQNxOpNdywqol3MvhSy+TopUY3BrFmzBrFQmsCYAxYZA1HFpUOHuhW7778PPPkk7ffxZWeyL5MSVYGxWq2YN28ehBDB/4EYJgDYXTIQ19qTiROB9eupdeWll1KvFl+CsfVNSmzaVHHbKgAYX5fDML7AlozBaGtPPv0U6NMHePtt2u3s6wOv524Bzl3QrVrluwmMr4FihjEatmQMxNWSmDkzMEvCczXuWYwZ0xi33DIT1dXViI+Pr72vr4FihjEatmQMItiWhGs1bkGBgvbt78Yf/nAUjRo1qhUYwPdAMcOEArZkDMAoS0KtxtUGebXiouJpDxPgey9ghgkWbMkYgJGWhFZgtDEYLXp7mFR4YgATalhkDMCXlHMg7N69u16BAeoPFHODKSaUsLtkAPWlnAO1JKxWK6655hpMmTLFo8AA3ts28MQAJtRw0yoD8KVDnT8P+jfffIOuXbuiSZMmfq/DWyc9hgkm3LQqhATTklBjMCNGjMAf//hHv9fBosKEGxYZg9BOFQh09pA2yDt27FgAXMXLmA8WGQNpiCWhl0XKzaXU+LhxlKLet4/cstGjuXsdE7mwyEQgp0+fxgsvvFBHYLiKlzErphEZIUQWgOkAHgNQCWAMgEop5dKwLsxPfHF3mjVrhpycHLRu3bo2i+RL7Q2LDBOJmKlOJhFAKoBdjo9kswlMfVsNFEXBa6+9BgDo3LlznTS1UbU3DGM0prFkAEBK2S/cawiU+twdi0XB6tUUgxk1ahTi4uLqnG9U7Q3DGI2pRMZsaF2jffuArCx9d2fwYAWrVtnRrRvFYFwFBnCfRKCiVvHa7Qb/MAwTIGZylyCEyBJCjHF8nl3Pe6cJIfKFEPlHjhwJ1RJrcXWNhADeecd9F7aiKFAUO+LivG8VcJ2HvWOHc6QKV/EykYyZLJkCAJBSFgOAEKJCCLFFSnmt3psd8ZqlAFX8hmyV0HeNrroKaNMGWLu2bibo1KlTaNWqGy688GGvWwWA4NTeMEyoMY3IqOKi+b5ACNFfCJHq+lq40csEqe7O739Pr3fu/AtatGiBvn2vxpNPXoX77mvk8XpcgMeYmZCLjCMV/ZCPb59ej4AUA8iCw2KJFPT6uajuzssvA82bK9i796+44or52Lq1N8aMaeRRNLgAjzE7IRcZKeVWAFv9OUcIkQpgl5SylTGr8g1fLQpPmaCpU4FvvlGwb58dgBWVlV28ujtcgMdEA2YK/D6mcywVfgpWoPjTTtNTP5etW0lgOnWyYvlyG3JyvI+O5TaaTDRgipiMlLJYCFGpPSaEGANgTSjiMf5aFHq7sHfs+B5bttiRlGTFE094ziJp4TaaTDRgCpEBKFskhJjm+DbRcWx6KO4dSEm/ayaoXbsuGDZsMm655TqfBAbgAjwmOjCNyAC1aemQE6hFkZ4OVFXlo2PHjmjfvj2AUX7dlwvwmGjATDGZsBFoY25FUfDoo49i+fLlAd2XC/CYaIBFxgcCacyt7Qdz1113BXxv13lLJSX0PaevGbNgKncpXPjbTtOXsSX+3p+tFsassMj4iK8l/TU1Ndi4cWPQBIZhzA5PKwgiUkoIIXDy5EkIIVhgmJjC07QCjskECUVRsGjRIlRVVaF58+YsMAzjgEUmCKgxmF9++QVnz54N93IYJqJgkWkgwQ7yMky0wSLTAL744gsWGIapB84uNYB27dqhX79+uPfee1lgGMYDbMkEwI8//ggpJTp37oy//OUvLDAM4wUWGT9RFAWzZs3C5s2bw70UhjEFLDJ+oA3yDh48ONzLYRhTwCLjI5xFYpjAYJHxgWPHjuHJJ59kgWGYAODskg+0atUKc+fORY8ePVhgGMZPWGS8oCgKhBD4zW9+g379TDshl2HCCrtLHlBjMBs2bEAsbCJlGKNgkdFBG+SdP38+hBDhXhLDmBYWGRc4i8QwwYVFxoWvv/6aBYZhgggHfh2cOXMGcXFxmDp1Kqqrq9G0adNwL4lhogK2ZEAu0syZM1FaWgohBAsMwwSRmBcZNQaTmJiIFi1ahHs5DBN1xLTIcJCXYYwnZkVm7969LDAMEwJiNvDbtWtXDB06FJMmTWKBYRgDiTlLprCwEKdPn0ZCQgKmT5/OAsMwBhNTIqMoChYsWIC8vLxwL4VhYoaYERltkHfChAnhXg7DxAwxMUHy4osvlj169OAgL8MYSExPkDxy5AgLDMOEiZiwZIQQRwAcCPc66qE1gKPhXkQEwb+Pupjh99FZStnG9WBMiIwZEELk65masQr/Pupi5t9HTLhLDMOEDxYZhmEMhUUmclga7gVEGPz7qItpfx8ck2EYxlDYkmEYxlBYZBiGMZSY3YXNRCZCiGkAKhzfpkopnwjnesKJECILwHQAjwGoBDAGQKWU0lTxGbZkIhAhRJYQ4nUhRIYQIlUIMdvx8EU1qsBIKddKKdcCWCuEeCnc6wojiQBSAexyfCSbTWAAtmQiFe0fVyWApTHyP/p0KWXtqE4pZbHjf/OYRfv7MCssMhFKNPxx+YMQIhFAhs5LlUKILCnl1lCviQkOLDJMpJAKstpcqXC8FpM4LLlE0O8mw4wWLYtMhBINf1x+kgRnwFdLJej3EIsUAOQ2AoAQokIIsUVKeW14l+UfLDKRSVT8cTENQ/3313xfIIToL4RIdX0tkuHsUgQipSzW/hFJKQsA9BdCRLvbkKRzLFatGE8UAzBVMJwtGQNxuDwP+fj26fX876T+cZkuhekj+dAXlCQ4LLtYwvEfyi4pZatwr6WhsMgYiCMj4ldWJJr+uPxBSlkphCgWQiRKKbUB4MQYziw9pnMsFX7+TYUbdpcik6j44woAO4DaokMhRAai/2fWxWHV1sm2CSHGAFhjpngMwLuwIxIhxDRtZafjj+taKeX0MC4rJDiqfovhKEiMgayaVzSV3okAYMbfB4tMhBINf1wMA7DIMAxjMByTYRjGUFhkGIYxFBYZhmEMhUWGYRhDYZFhwo6jMddLQggphDgmhLB7eN8Wx3tej4EtFlEDZ5eYiEEI8TqALE/Vzg7xec2xl4sxCWzJMJHESwAS9brhOQoSWWBMCIsMEzE49ihVAhirPe7YXgAWGHPCIsNEGmsA3Kh+42jLmeVoLM6YEBYZJtJ4HXVdprm8pcLccOCXiTiEEMdAu693giY16PX+ZUwCWzJMJLIGNMhsKwuM+WGRYSKRJJDAcKA3CmCRYSKRDABbwr0IJjiwyDARhSObFAtdAGMGFhkm0sgCuCYmmmCRYSICx/4lO6jPL4QQdt6fFB1wCpthGENhS4ZhGENhkWEYxlBYZBiGMRQWGYZhDIVFhmEYQ2GRYRjGUFhkGIYxFBYZhmEMhUWGYRhDYZFhGMZQ/h+ilh9Ewxb1ugAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Init Loss = 4.1642585\n",
      "Final Loss = 0.7604558\n",
      "Fair dummies test (regression score), p-value: 0.4015984015984016\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlkAAAElCAYAAADX8PdgAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3dbZCV5Z3n8d/fRtjB6uZ02w4SsTCHxNJuLEg3UxYOaJk0E5pKWE1BspAaZ7HWxp4qyIsQWFdipSawBsK+GH0x0GZhZWqhpIm6CQU4tGNFULGG7sFKaGuNnEjUiY4t3cBCRaTz3xf3feBwOP1wHu7z0Hw/Vafg3I/Xfbrp/nFd//u6zd0FAACAwrqu1A0AAAAYiwhZAAAAESBkAQAARICQBQAAEAFCFgAAQAQIWQAAABEgZKGozKzNzHyI14kStqHfzA6aWbwYbShHZrY4/BxGenWXuq0AUAnGlboBuGbNl5SQFEtZNjCaHcNf8qfcfX6ebah194HwmE2SHpPUbWZL3L0rz2NHooDXnkmXpK+lLeuW1CFpawTnA4AxjZCFUjmaDDg5KMgv/NTzu3uPpCVmtlVSp6TaQpwjApGFnfDz6EldZmaSdCL8fPJiZjFJLe6+J99jAUAlYLgQFcfdO9y9I6Jjr5AkM9sYxfHzFeW1F8FsBQEWAK4JhCzgarslLS51IwAAlY2QhVKpM7NY6iu5wsxazKw7pSD9il4lM+s0s86U91vNbGO430EzW5Nn27olXSqATz9fuMzNbHHK+63ha42ZnQjbvSa8toPJwv6w9iuffTK1ZWNKUfrGlOVrwmXJ47Tl+blcIZvzmtlBSQdTPjs3s5ZCtgcAyg0hC6VyQlJ/6ivlzr64pLWSpktaImnNCAGhTkHPU6eCX+T5DqclJCnLOw3rJCXbOD9sw0ZJL4d/Niso7O/Mc58rhIFrsYKC9SWSFoehpiU8xtcUfI4rJJ3K4nqGle15w0L9JZLk7ha+yvLmAgAoFArfUSq1QxW+p9UcJcxsj4Jf0MOFp7ik5kIUaIfHkrsnstxvwN03SZKZPSlpjYIC/65w2VZdXbieyz4K18UVBJ1L121maxXcJfmkgoCWCD/nbK9lSKU6LwBUGkIWylI4FPcdBYGnScH0AsPpKVDAkoIemFzufDya/Iu7D4R35qXOKZWpJymXfZKSw23d4X6pusJ9+82sS9LBZJgrgFKdFwAqCsOFKDvhXFCPSXpOQQ/WaG75L2SPSZuC4vdsZQpmIw3R5bLPFfunDL+lvgbcPTncOiBpY3otV55KdV4AqBj0ZKGshEXeTZKmJ4frMvSWRHn+jQomSF2btiq1MD+m8tAlKWZmTUP14oVzUu0JewYLFXZKdV4AqCiELJSbZI/U2jDwNCmo/4m0SDoMdysU9GLNz1AvFg9rkU5JeibKtoyWuyfMrENSp5mtUPDZNSkYZn1O0l+Efw7o8gz7pTxvMjQvDteJ4ncAYxnDhSgrYbjZpCDsdCv4Jb1HudVIDStlKgFXcEdfnYIetPRf/E+G604o6JV5UkFgKHibshVOnposjj+hy8OsPQrq2V4Ol89WeHdfqc4b9nr1KPgMy3KyVwAoJHP34p4wGGr5toLeglH90A9v30/WqcQppAUAAOWuqMOF4ZBMcshlVHMQJQNW8nlnZhY3s63Jx58AAACUo6L3ZEmXwtYz7t48im2707czsxPhHUwAAABlqaxrssKhxaYMqwZ4JAcAAChnZR2yFAwpDjWPUDaPPAEAACiqcp/CoU6ZJ2YcUMq8RanCGq42SfqzP/uz5ttuuy33s589K/3pTyNvd911UnV17udByZw9Kw2O4kssSVV8mXGNePvtt/vc/aZStwOodOUesrIWPveuQ5Jmz57tR48eHWGPYY51370afPRvRtyuasuzsl+9mvN5UDpz5l3U/UtHN8H6K7vq9MahMfdPBriKmZ0sdRuAsaDchwuloDcrXbnMuA0AAJBRuYeso8ocqOoUTGoIAABQlso6ZIWzfycyPCsuxuM4AABAOStVyMo0BJicaLQzLVRtVFjIHm7TpIifYwcAAJCvYs/4HlfwsN/5kprCBwCfCIvVpWBahhYFISz5ANkOM2sL58WKKXisDrO9AwCAslbUkOXuCQUP/8347MFwCLA2w/KODJsDAACUrbKuyQIAAKhUTPoDACiZ7u7u26qqqtquu+66Vne/aiQDKFdm1v+nP/1p/+DgYEdzc/N7mbYhZAEASqK7u/u266+//vnJkyfHYrHY2fHjx/eZWambBYzI3XXhwoXrBwYG/tPHH3+8oLu7+1uZghbDhQCAkqiqqmqbPHlybPLkyacmTJjwOQELlcLMNGHChM8nT558avLkybGqqqq2TNsRsgAAJXHddde1xmKxs6VuB5CPWCx29rrrrmvNtI6QBQAoCXevHT9+/OelbgeQj/Hjx38+VD0hIQsAUDIMEaLSDfc9TMgCAACIACELAAAgAoQsAAAitHnz5nozax7q1dfXV1XqNmbS2toab21tjee7TTkxs+YXX3yxejTbFuLamCcLAFB2Fn5zcPofPvLxpW5H0pSb7cK+X1adyOcYL7zwwju33377Z6dOnbrid299ff1gfq27UmNj452TJk26+Prrr/+2kMdF9ghZAICy84ePfPz93/20bO48fOV/35h34Js7d+75MFBdKECThrR8+fJPUt/39fVV/fKXv6xZvnx5f5TnxdUIWQAAjCGrV6/uS31/+PDhiQ8//HB8+fLl3aVq07WKmiwAAIAIELIAACiCf//3f6/q6+u74pW+zbJly6bV1NTMuvXWW2ds3769Nr1QO1Mxtpk1b9++vTbTNvfcc8+XH3zwwduT26Ue78UXX6xubGy808yaa2pqZrW3t9+S3p7Tp09Xtbe333LrrbfOGGqbdO3t7bfU1NTMGs32y5Ytm7Zs2bJp69atm5w8x7p16yb39fVV3XPPPV82s+Zbb711xuHDhyem7pe6fWtrazzTZ9nX11fV2toaN7PmxsbGOzdv3lyfb3uzRcgCAKAIGhsb77rppptmpb5SA1Nra2t8165d9U888cQHu3btSjzzzDMZQ0E2Xn/99d9u27YtIUnu3u3u3Q888MBZSXr33XcnbNiw4YPjx4//eseOHSe2bNlyc3oQeeONN2qmT5/+2UsvvfTOqlWr/pBpm1Stra3xvXv31u7bt++dHTt2nNi7d2/tcNv39/dX7dq1q16SXnrppXe++93vfrJhw4ap99133+1r1qz56NChQ2/X1NQMLl269NLntGzZsmnPPvvsTU8//fTJI0eO9NbW1g7G4/G70oPWokWL4q+99lrNtm3bEp2dnSd6enpuyLe92aImCwCAIvjkk0+ODXUnYW9v7/gDBw7Ubtu2LZEsUP/Zz352srGx8a6o2pNau9XQ0HBhwYIF/c8//3xt6vI5c+acSb5fv379x4lE4j/8/d///c3pdV+p13Do0KG3586de16S+vv7P9y8eXPG7ZOqq6sH169f/7Ek/fjHP/5oy5YtN8+cOfN8MgwuX778kx/84AfTpKC+bNeuXfXHjx//dUNDwwVJ2rlz58nGxsaJq1atmrpz586Tyba88cYbNS+88MI7yePs3LnzZDLQ5dPebBCyAAAosTfffPMGSUq9AzAZIqK0ffv22t27d9f+/ve/n9Db2ztxzpw5Z4bbfv78+WdSg0qqffv21UjSvHnz7symDTNmzDiX/HsyhDY1NV1aduONN14KpkeOHJlYXV09mP7Z3HvvvWdeffXVmuT75OeZDFiFbG82CFkAAFyDGhsb75Sk1atXf3T33Xef+/73vz/19OnTw06M+umnnw67vrq6evDMmTPHsmnHpEmTrurdSw1WqQYGBoY8/5kzZ7Ke1DWX9maDmiwAAErs7rvvPicFxegjbZsahHKdLf7w4cMTe3t7J3Z2dp5Yvnx5/2h7zV5++eWahoaG85nWLVy48MzZs2er0ovUC2nBggVnM51j7969tfPmzbvUa5Xp80z/rIrRXnqyAAAogsOHD0/MVJN1xx13fNbQ0HChoaHh/EMPPTR9x44dJ26//fbPvv/970/NdJz3339/Qm9v7/g///M/H/zrv/7raSOd98tf/vJnUjA0WFtbe1EKJkaVpPXr109Zt27dH958880bDhw4UJs+XPjGG2/UbN68uX7hwoVndu7cWZusG8t0noaGhgtLly7tW7p0afzpp58+efvtt3/25ptv3rB79+7a/fv3Z9wnW3Pnzj2/YMGC/uQ56uvrBzds2HDz6dOnxz311FMfpLZlzpw5Z1auXDlN0sn6+vrBNWvWXHHnYDHaS8gCAKAIklMppEsWu//qV796Z9GiRfEHH3zw9oaGhvOrV6/+6MCBA7Wp2z7++OMfLVy48PbGxsa75syZc2bTpk0f/uY3v5mYDE+ZzJ0793xDQ8P5hx9+ON7Q0HB+w4YNH9TX1w8++uijH23ZsuXmvXv31n7jG9/oX7BgwRUzwtfW1g4uWLCgv6en54a/+7u/mzpp0qSLP/3pT08ON3P8zp07T65bt27yypUrp33wwQcTkteR7Wc1nP379yfa29tvWbly5bTTp0+PmzFjxrkjR470pgfYX/ziF4lFixbFH3rooekzZsw498gjj/S9//77E4rZXnP3Qh2r7MyePduPHj2a8/5+370afPRvRtyuasuzsl+9mvN5UDpz5l3U/UtPjWrbV3bV6Y1D/L8EY5+Zdbv77KjP89Zbb703c+bMjHdxfeUvLt5ZZo/Vuf5f/2Xc28U+r5k1p94hh/L01ltv1c+cOfO29OX8xgAAlJ0pN9uFQjwvsFCm3GyR3+mHsYeQBQAoO/t+WXWi1G0A8sXdhQAAABGgJwsAgDLl7t2lbgNyR08WAABABAhZAAAAESBkAQAARKAkNVlm1iYpOTlR3N03jXKfpJikDncfiKJ9AAAA+Sp6yEoGLHffE76Pm9lWd18xzD5rlBaqzGyrpCH3AQAAKKVSDBeuSAYsSXL3hKSWEfb5iwy9VgNmFit46wAAAAqgqCErDEVNGVYNmNlwQSueYX2M4UIAAFCuit2TFZeUKRidCtcNZa2kg2a2Ubo05Li18M0DAAAojGKHrDpdLnhPNaCgmD0jd++S1CypzcxcUsLde6JpIgAAQP4qYsZ3M4tLmi3pi5IeU9CrtcLdOzJs2yapTZKmTJmiY8eO5X7ilhb59TeM3L6WFimf86Bkvj7f9YXYH0e57Yc6dswibhEAYKwoRciqy7BspAL2tSl3H641s+ckvWxmibCX65IweHVI0uzZs33WrFk5N9S/t0qDX7plxO2qurpkP3wi5/OgdNpXXtT99Zk6V6/2ysE6/eiJivh/CVDx/D8umq6PPxpf6nZcMvnmC/Z/fpHTQ6s3b95c/4Mf/GBa6rKpU6d+Nm/evLNPPfXUB/X19YOFaWRmra2tcUnav39/IsrzlJqZNb/wwgvvPPDAA2dH2rZYn0mxf2McVeZAVScp4/BfWPB+MHWZu/eY2RJJ8yV1ZdoPAFDBPv5o/GD7f/681M1IqvqH/5V34Es+h7Cvr6/q8OHDEx9//PGp8Xj8rkQi8euog1aqxsbGOydNmnTx9ddf/22xznmtKmpNVng3YCLD1Aux9B6pUTgq6dPCtAwAgOKor68ffOCBB84eP3787UmTJl1ctWrV1GKef/ny5Z9861vf6i/mOa9VpRj72KigZmqTJJlZk1J6o8L6q42SHnH3AXfvMrNOSXvSjtOmcFgQAIBKNGPGjPPvvfdeUYdFV69e3VfM813Lih6y3L3DzNrCYcCYgsfqpM7cHlcwOWmdLk/38Eg4fUNyPDwmaQ/zZAEAKlVfX1/Va6+9VrNq1ao/lLotiEZJHhDt7h3u3uXue9KfWxgurw1ngk8uG3D3teF+He6+KXU9AADlrq+vr6qvr6+qt7d3/Pbt22u/8pWv3PmXf/mXZ9avX/9xcpsXX3yxurGx8U4za66pqZnV3t5+xd1X69atm1xTUzPLzJpvvfXWGZs3b65Prmtvb7+lpqZmVqb9UrW2tsaThd+StGzZsmnLli2btm7dusm33nrrjJqamlnr1q2bnL7faI+faqhj9/X1Vd1zzz1fTl7H4cOHJ6ZfZ3L71tbWeF9fX1X6Z9na2ho3s+bGxsY7Uz+HfNtcSCUJWQAAXGtuuummWTfddNOsxsbGux5++OH4vHnzzv7jP/7jydRt3n333QkbNmz44Pjx47/esWPHiS1bttycDBAvvvhi9YYNG6bu27fvnePHj//66aefPnnjjTcOSkFw2rt3b+2+ffve2bFjx4m9e/fWDhU80vX391ft2rWrXpJeeumld7773e9+smHDhqm9vb2XhjFzPf5Qx77vvvtuX7NmzUeHDh16u6amZnDp0qVXhL5nn332pqeffvrkkSNHemtrawfj8fhdqUFr0aJF8ddee61m27Ztic7OzhM9PT1XzbeUz2dSKNyPDgBAESTvLpSk3t7e8evXr58Sj8fv2rFjx4nktAOp9VINDQ0XFixY0P/888/Xrl69uq+/v39cdXX14B133PFZfX39YENDw4XksQ4cOFB76NCht+fOnXtekvr7+z/cvHnzzaOtv6qurh5M9qj9+Mc//mjLli03v/POOxMaGhou5Hv8TMeeOXPm+eQ1L1++/JPkFBeHDx+euGvXrvrjx4//Onl9O3fuPNnY2Dhx1apVU3fu3Hmyt7d3/BtvvFGTOl3Dzp07TybDXKE+k0IgZAEAUGQNDQ0Xdu7ceXLZsmV66KGHpp85c+bSjNbbt2+v3b17d+3vf//7Cb29vRPnzJlzRpK++c1vnnniiScu3nTTTbPmzJlz5qtf/eqZ9evXf7xv374aSZo3b96dubZnxowZ55J/T04n0d/fP06S8j1+pmM3NTVdWpbsjZOkI0eOTKyurr4UIJPuvffeM6+++mqNJL355ps3SNJw82EV4jMpBIYLAQAokb/927/95OzZs1XJmqSwvujmb3/72/2dnZ0nFixYcGmqhfr6+sH333//N9u2bUtMmjRpcMOGDVOTtVXV1dWD7t6d/hptOyZNmjTsPF35HD/TsVODVaqBgYGqTMsl6cyZM0OuyyTfz6QQCFkAAJTIgQMHqiVp7ty55w8fPjyxt7d3Ymdn54nly5f3p/fmJC1fvrx///79iW3btiUOHDhQu3DhwjOpQa3Qoj5+qgULFpzNdK69e/fWzps376wk3X333eekoEYtuT69ML6YbR4OIQsAgCLr7e0d397efsuGDRumPv744x9I0h133PGZJK1fv35K8g7EAwcO1Cb32b59e217e/stYRgbf/DgwZqpU6d+1tDQcGHp0qV9S5cujb/44ovVyX1T7yDMR9THTzV37tzzCxYs6E+e6/DhwxNbW1vjp0+fHvfUU099kGzPnDlzzqxcuXJacptFixZd0ZZitnk41GQBAFAEZtac/Ht1dfXgjBkzzqUWb9fX1w8++uijH23ZsuXmvXv31n7jG9/oTx0uvPvuu8/t3r27duHChbefPXu2qqGh4fyuXbsSUlD4vW7duskrV66c9sEHH0xoaGg4v3r16o8K1faoj59q//79ifb29ltWrlw57fTp0+NmzJhx7siRI72pjx76xS9+kVi0aFH8oYcemj5jxoxzjzzySN/7778/oVRtHoq5ezHPV1SzZ8/2o0eP5ry/33evBh/9mxG3q9ryrOxXr+Z8HpTOnHkXdf/SUT4geled3jjE/0sw9plZt7vPjvo8b7311nszZ87MeKfXWHpANMa+t956q37mzJm3pS/nNwYAoOwQaDAWUJMFAAAQAUIWAABABAhZAAAAESBkAQAARICQBQAAEAFCFgCgVHwsTyOEa0P4PZzxG5mQBQAoCTN779y5cyV97AmQr3Pnzk00s99lWkfIAgCUxOeff77pvffeG9fX1xe7cOHCOHq1UCncXRcuXBjX19cXe++998Z9/vnnmzJtx2SkAICSaG5u3t/d3f3uhx9++Ni//du/zXT3OvGff1SGP5lZv7u/evHixSebm5t/m2kjQhYAoGTCX04Pl7odQBT4HwMAAEAECFkAAAARIGQBAABEgJAFAAAQAUIWAABABLIKWWb21agaAgAAMJZk25PVZWafmtk/mNn9kbQIAABgDMg2ZNVKekzSdEkvh4HrSTObVfimAQAAVK6sQpa7n3b3Dnf/K0l1kv6rgsDVaWa/DQPXtCgaCgAAUElyLnx39wFJn0o6JcnCV7Ok35nZATO7rRANBAAAqERZhywzmxXWZH0qqVPBEOIKd/9SSg/XSUn/VNimAgAAVI6snl1oZu9K+qKkf1UwVLjb3U+nbuPuA2a2UdIjBWslAABAhcn2AdFdkra4+7ERtvtU0vyhVppZm4JhRkmKu/um0ZzczNZIGkju6+57RrMfAABAsWUbsrolfUfSFSHLzP5J0nPu/j+loEBe0suZDpAMWMmAZGZxM9vq7iuGO7GZHZS0JKwFk5n1m1lX8j0AAEA5ybYma42kf8mwvEPSo6M8xorUHih3T0hqGW6HsAerMy1QNROwAABAuco2ZE2X1JNhebekppF2NrPYENsNmNlwQesxBUOVl4ThDAAAoCxlG7ISkr6SYXmTpN+NYv+4gpqqdKfCdVcJg1ks/PtiM2sxszXhcgAAgLKUbU1Wh6SfmdmAu78iXXqeYYekJ0exf50uF7ynGlAYpDKYnVyfUsd1VMH0EVcV14c1X22SNGXKFB07NlKN/jBaWuTX3zDiZtbSIuVzHpTM1+e7vhD74yi3/VDHjlnELQIAjBVZhSx332RmNyp4pI6nrHrG3TcXtmlXiCnoRUu2Y8DM6sysyd2vGL509w4FoU+zZ8/2WbNyf+KPf2+VBr90y4jbVXV1yX74RM7nQem0r7yo++sz5f6rvXKwTj96Itv/lwAArlVZ/8Zw97Vm9t8V9CK5pB53H81QYVJdhmXDDf0lwvOmDzOeUlAwn6lGDAAAoKRy+m95OEVDLnNUHVXmQFWnIcKSuyfMhhyi4e5CAABQlrIOWeEzCePKEJbc/fnh9g2H+RJmFkvrmYq5e9eQO0o9ZhZPu6MwriC0AQAAlJ1sH6vzA0k/UfAw6HQuqWoUh9mooDB9U3jMJqVMz2Bm8XCbR1KC2NrwtSJln0R6PRYAAEC5yLYn6zFJP5X0ZPozC0fL3TvMrC2cFyum4LE6qbO9xxXUWtUpHA509y4zi4WTkkrSje4+5GN7AAAASi3bkHVKwbMLcwpYSeEdgEOt65JUm2E5zykEAAAVI9vJSLdKWhxFQwAAAMaSbHuypktqM7PHFBSdX3F3n7t/p1ANAwAAqGS5TOGQOmzH9NcAAAAZZDvj+6NRNQQAAGAsybYmS2b2X8zsqJldTFn2nJk9WNimAQAAVK6sQpaZ/UTBfFVPpu3bKem/FbBdAAAAFS3bmqw2SYvd/Z/THhDdLampcM0CAACobFkPF0r6NPwzteg9Limbh0QDAACMadmGrGck7TGzWQoeoyMz+4qkLeELAAAAyv7uwrXhswV7JMnMPlXwaJwOd98cQfsAAAAqUtbzZLn7kjBofS1c1OXuDBUCAACkyGUyUrl7QlKiwG0BAAAYM7IKWWb2rWFWD7j7P+fZHgAAgDEh256sPUMsdwU9W1/OrzkAAABjQ1Z3F7r7dekvSXWSfi5pSSQtBAAAqEC5zJN1BXcfkPSIpI35NwcAAGBsyDtkheoktRToWAAAABUv28L35zIsjimY8X2gIC0CAAAYA7LtybIMr9OS/lX0ZAEAAFyS7Yzv346qIQAAAGNJoWqyAAAAkCLbmqyvZrM9k5MCAIBrVbaTkXYpmHjU0pYPtawqx3YBAABUtGyHC5dIelnS9LQJSX8uaXHaRKUELAAAcM3KNmT9RNIad/9d2vK1kh4rTJMAAAAqX7Yh60ZJzRmWxyU15d8cAACAsSHbkNUpqcPMHkwuCIvhdyuo1wIAAICynydrhZnVSfq5mXm42CT1iAdEAwAAXJLt3YVy9yVm9kVdHh7syVCjBQAAcE3LOmSZ2X+R9KikWe4+Llz2nKTn3P35ArcPAACgImVVk2VmP1FwJ+GTunJerE5lcXehmbWZ2eLwtSabNoT7d2a7DwAAQDFl25PVpmA+rH82u2Lu0W6N8u5CM2uTdMrd94Tv42a21d1XjHL/JkmLs2s2AABAcWU9XCjp0/DP1JQVlzTauqwV7n5pGgh3T5hZSxbnj2exLQAAQElkO4XDM5L2mNksBY/NkZl9RdKW8DUsM4spc4/XwGiClpktTvaAAQAAlLNsp3BYa2ZxBVM2yMw+lRST1OHum0dxiLikgQzLT2mEHqrwvIls2gsAAFAquU7hEJf0tXBRVxZTONQpCFTpBhSEteE00YsFAAAqRVYhy8xWK+i1SqiIvUrhUOKoZpQPC+vbJGnKlCk6duxY7iduaZFff8PI52xpkfI5D0rm6/NdX4j9cZTbfqhjx2zkDQEAUPY9WX+loBbrf+RxzroMy4bsxQrruOTumYYZr+LuHZI6JGn27Nk+a9asXNoYHOt7qzT4pVtG3K6qq0v2wydyPg9Kp33lRd1fn6lz9WqvHKzTj57I5V4RAMC1KNvfGBsl7TazLnd/K4fzHVXmQFWnsM4rgzbp0tQNl4Tzaw2EoQoAAKCsZBuy5kvql9RjZj1KGzJ09+8Mt7O7D5hZwsxiaT1TMXfPOBzo7pvSl5nZxkzLAQAAykW2UzjEFPQ4/VzBvFiW9hqNjQp7p6RLPVRdKe/jZtaZHCYEAACoRNlO4fBovid0947wsTotCkJbPG2297ikFgVDiFfUYYX7LAn/vlVS51A9YAAAAKU0ZMgys99Kmu/u76Utvy19WbaGq6MKQ1PtMOu6JI3qETwAAAClMtxw4XSlFamb2SRJJyJtEQAAwBiQbU2WNPraKwAAgGtWLiHLC94KAACAMWakkDWqCUABAABwpZHuLvyimaVOh10nSWZWrQzDhu5+poBtAwAAqFgjhaxM0yOYhu7hqsqvOQAAAGPDcCGLaRIAAAByNGTIcvdnitkQAACAsSSXuwsBAAAwAkIWAABABAhZAAAAESBkAQAARICQBQAAEAFCFgAAQICNyVYAAArASURBVAQIWQAAABEgZAEAAESAkAUAABABQhYAAEAECFkAAAARIGQBAABEgJAFAAAQAUIWAABABAhZAAAAESBkAQAARICQBQAAEAFCFgAAQAQIWQAAABEgZAEAAESAkAUAABABQhYAAEAECFkAAAARGFeKk5pZm6RT4du4u28a5T6SNF1STNJadx+IqIkAAAB5KXrISgYsd98Tvo+b2VZ3XzHcPu7ekfJ+saRuBYELAACg7JRiuHBFMmBJkrsnJLUMtbGZxZUWpsL968KwBQAAUHaKGrLMLCapKcOqATMbMmhJasuw7JSkuoI0DAAAoMCKPVwYl5SpjupUuO4qYU9X7RDHOlq4pgEAABROsUNWnS4XvKcaUFDMPiphXVeXu/cMsa5NkqZMmaJjx47l2FRJLS3y628YuT0tLVI+50HJfH2+6wuxP45y2w917JhF3CIAwFhRkrsL8xHWaK1w9+ZM68MC+Q5Jmj17ts+aNSvnc/n3VmnwS7eMuF1VV5fsh0/kfB6UTvvKi7q/PlPuv9orB+v0oycq7p8MAKBESlH4nqmOatS9WJI2SvpagdoCAAAQiWKHrKPKHKjqJF019JfOzDaK+bEAAEAFKGrICsNRIrzLMFXM3buG2zestdoaFsInlw13RyIAAEDJlGK4cKNSpmQwsyZJXSnv42bWmRrEwjB1NBmwzCxGwAIAAOWs6FW87t5hZm1hSIopeKxO6mzvcQWTk9YpmD8rLumgJJlddWdXpqkdAAAASq4kt0qlPiInw7oupYSnsPeK++YBAEBFKcVwIQAAwJhHyAIAAIgAIQsAACAChCwAAIAIELIAAAAiQMgCAACIACELAAAgAoQsAACACBCyAAAAIkDIAgAAiAAhCwAAIAKELAAAgAgQsgAAACJAyAIAAIgAIQsAACAChCwAAIAIELIAAAAiQMgCAACIACELAAAgAoQsAACACBCyAAAAIkDIAgAAiAAhCwAAIAKELAAAgAgQsgAAACJAyAIAAIgAIQsAACAChCwAAIAIELIAAAAiQMgCAACIwLhSnNTM2iSdCt/G3X1TFPsAAACUStFDVjIsufue8H3czLa6+4pC7gMAAFBKpejJWuHuzck37p4ws5YI9gEAACiZotZkmVlMUlOGVQNDhaZc9gEAACi1Yhe+xyUNZFh+KlxXqH0AAABKqtghq06Xi9dTDUiKFXAfAACAkirJ3YVRCovk28K3/8/M/m9eB3z1UPqSekl9GU6c12nKXOZrHiOOHM64OOM1j+Ev85j+Gg+Bax7atKgbAlwLShGy6jIsG6lHatT7uHuHpI5sGzVaZnbU3WdHdfxyxDWPfdfa9UpcM4DoFXu48Kgyh6M6ST0F3AcAAKCkihqy3H1AUiK8YzBVzN27CrUPAABAqZXisTobdblmSmbWJKkr5X3czDrTQtWw+xRZZEORZYxrHvuuteuVuGYAETN3L/5Jg+L0hIJhwCsekRPOfdUpqdndE6PZBwAAoNyUJGQBAACMdaUYLgQAABjzxtw8WflIPog6fDuqIclc9ikneVyzJE1XMHy7NrxBoSLk+zUzs053X1L4lkUn12s2szUKJv49JUnJh7RXgjy/t6Xge7ujEr63wxrWb0uaP9rvzUr/2QVUAkJWKPkDJ/lLJCzA3+ruKwq5TznJ9ZrDuciS7xdL6lYQuMpevl+z8KaLxVG2sdByvWYzOyhpSTJkmFm/mXVVSOjI5Xt7jdJClZltlVTW/57D78m4snjUWKX/7AIqBTVZITPrdvfmtGUn3H3I8JDLPuUk2/abWVzSCndfm7a8X9IjldDLke/XLAyVne5eMXO/5/i9vUbSQFqgjqfejFLOcrzmq3oozWyjpCcrJFg2SXom/bqH2Laif3YBlYKaLF3qam/KsGogvNuxIPuUkzza35Zh2SllnpW/rOT7NTOzxZUQJFPlcc2PKW2alAoKWLleczzD+lglBKxsVPrPLqCSELICcQV1J+mG637PZZ9yknX73T3h7rVDHOtoAdsWlZy/ZmEvXkWEjDRZX3P4SzgW/n2xmbWY2ZoMEwKXq1y/zmslHQx7r5JDalsL37ySq/SfXUDFIGQF6nS5ADTVgIZ+rmIu+5STgrQ//EXU5e6V8IijfK65qUKuMV0u1zw7ud7d94RPVuhQMH9dJcjp6xxeZ7OkNjNzSYkK/ZqPpNJ/dgEVg5CFnKXUaM0vdVuiFA6hXGuPcIoppecuHDKrC+t+xqTw+3m2pC9K2qSgVyvT8DgAjAoh67JMNUUj/a8ul33KSb7t3yjpawVqS7Fkdc3JIbIKr8vJ9uuckDJe8ylJlVKzk8v39lp373D3gfDmjmZJG8donVKl/+wCKgJTOASOKvMPmDpJQw0X5LJPOcmr/WHdSkXNj6XcrrlNunTn1iWZ7r4rU1lfs7snzIa8ebISvt5ZX3MYpA6mLnP3HjNbImm+xlZPZqX/7AIqBiFLwf/YzSxhZul3EsXCOo2C7FNO8ml/siA47dmSLeV+3Tl+na+aoNHMNlbKxI15fJ17MkzZUBE3OBT43+ZRZb4Tr2JV+s8uoJIwXHjZRqVMTxD2XHSlvI+bWWfaHVbD7lMBsr7m8H/8R5O/fM0sVmHDKbl8nStdLte8Nnyl7lNJheBZXXMYLr6T4ThtCor+K0HGaVTG6M8uoCIwGWmKsIcmoaAr/YrHTIRBolNSc1oPzpD7VIJsrjksDD4xxKFqK2XoMJevc8q6Jbr8i7ezUv7nn+P39mJdvqX/xvRJaMtdttcchpDHdPl7PCZpT7nPDxb+u1ysYFizRUHR/onkUPZY/dkFVAJCFgAAQAQYLgQAAIgAIQsAACAChCwAAIAIELIAAAAiQMgCAACIACELKDPhnEaV8jBmAMAQCFkAAAARIGQBAABEgJAFAAAQAUIWAABABAhZQAZmttHMtprZGjM7YWb9ZrZ1FPttNbMTactiZubJB2mbWYuZdYfL+s1s4wjHvKoQPtx3cYY294/mmACA6BGygMziCh4ELQUPhe6Q1GZma0bYr1NSPHxob9K3JQ2kPEw6LmmtpOnhsdeED+vNWRjCFkv6WnjMxfkeEwCQH0IWMLQud9/k7j3uvlbSHkkrhtshDFIDadutkLQ7ZZsOd+9y90S4/R4FwSgnYaBbLGlJ2NYuBSFu2LYCAKJFyAJG76CCXqiRdCgIPTKzmKQmSVcMNZrZ4nAYsDu5bR5awj+TQ5CuoEetKc/jAgDyQMgCchDWaXnKK7VmaquCIcMmXR4q7EnZt1vSY5KeU9CDtacATRpwd0t/FeC4AIAcjSt1A4AKskRSjyS5+/ShNnL3hJklFPQwzVfQsyVJCoNXk6Tp7p4Il43m3LGUY8TS1nVJiplZU2qYAwCUFj1ZwNBazKzNzOLh3Xotkp4c5b5bJX0n3Cd1qDAR/rk2PO5ijW64MB5uH5P0TOqKMKx1SOoM71yMJ4cjR9lWAEAECFnA0HokNUtK1k2tcPfRDu11KOixSiR7rCTJ3QckbVJw52K3gp6uPQqK5YfypKQ6SScU1Fo9qSCsXdrH3VcoCHNbw+2Sw5EAgBIxdy91G4CyE/YCxdx9fqnbAgCoTPRkAQAARICQBQAAEAFCFgAAQASoyQIAAIgAPVkAAAARIGQBAABEgJAFAAAQAUIWAABABAhZAAAAESBkAQAAROD/A7VSBqxq1rVSAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plot_groups_pointwise(x_axis_0_test,\n",
    "                      x_axis_1_test,\n",
    "                      A_test,\n",
    "                      Y_test,\n",
    "                      Yhat_out_test,\n",
    "                      y_u=None,\n",
    "                      y_l=None,\n",
    "                      point_color='b.',\n",
    "                      shade=None,\n",
    "                      filename_0=\"supp_synth_use_A_equi_0.png\",\n",
    "                      filename_1=\"supp_synth_use_A_equi_1.png\",\n",
    "                      save_figures=save_figures)\n",
    "\n",
    "\n",
    "p_val_equi = utility_functions.fair_dummies_test_regression(Yhat_out_cal,\n",
    "                                                            A_cal,\n",
    "                                                            Y_cal,\n",
    "                                                            Yhat_out_test,\n",
    "                                                            A_test,\n",
    "                                                            Y_test,\n",
    "                                                            num_reps = 1,\n",
    "                                                            num_p_val_rep=1000,\n",
    "                                                            reg_func_name=\"Net\",\n",
    "                                                            lr = 0.01,\n",
    "                                                            return_vec=True)\n",
    "\n",
    "\n",
    "plt.clf()\n",
    "x_equi, bins_equi, p_equi=plt.hist(x=p_val_equi, bins=30, color='royalblue', alpha=0.7, range=(0, 1), density=True, ec='blue', align='left', label=u'Equitable model')\n",
    "x_base, bins_base, p_base=plt.hist(x=p_val_base, bins=30, color='salmon', alpha=0.7, range=(0, 1), density=True, ec='red', align='left', label=u'Baseline model')\n",
    "\n",
    "plt.grid(axis='y', alpha=0.75)\n",
    "plt.xlabel('p-value')\n",
    "plt.ylabel('Frequency')\n",
    "plt.title('Fair Dummies Test')\n",
    "for item in p_equi:\n",
    "    item.set_height(item.get_height()/sum(x_equi))\n",
    "for item in p_base:\n",
    "    item.set_height(item.get_height()/sum(x_base))\n",
    "\n",
    "plt.ylim(0, 1)\n",
    "plt.legend(bbox_to_anchor=(1.05, 1.0), loc='upper left') \n",
    "plt.show()\n"
   ]
  },
  {
   "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.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}