{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Python: Difference-in-Differences\n",
    "\n",
    "In this example, we illustrate how the [DoubleML](https://docs.doubleml.org/stable/index.html) package can be used to estimate the average treatment effect on the treated (ATT) under the conditional parallel trend assumption. The estimation is based on [Chang (2020)](https://doi.org/10.1093/ectj/utaa001), [Sant'Anna and Zhao (2020)](https://doi.org/10.1016/j.jeconom.2020.06.003) and [Zimmert et al. (2018)](https://arxiv.org/abs/1809.01643).\n",
    "\n",
    "In this example, we will adopt the notation of [Sant'Anna and Zhao (2020)](https://doi.org/10.1016/j.jeconom.2020.06.003).\n",
    "\n",
    "In the whole example our treatment and time variable $t\\in\\{0,1\\}$ will be binary. \n",
    "Let $D_i\\in\\{0,1\\}$ denote the treatment status of unit $i$ at time $t=1$ (at time $t=0$ all units are not treated) and let $Y_{it}$ be the outcome of interest of unit $i$ at time $t$.\n",
    "Using the potential outcome notation, we can write $Y_{it}(d)$ for the potential outcome of unit $i$ at time $t$ and treatment status $d$. Further, let $X_i$ denote a vector of pre-treatment covariates.\n",
    "In these difference-in-differences settings [Abadie (2005)](https://doi.org/10.1111/0034-6527.00321) showed that the ATTE\n",
    "\n",
    "$$\\theta = \\mathbb{E}[Y_{i1}(1)- Y_{i1}(0)|D_i=1]$$\n",
    "\n",
    "is identified when panel data are available or under stationarity assumptions for repeated cross-sections. Further, the basic assumptions are \n",
    "\n",
    " - **Parallel Trends:** We have $\\mathbb{E}[Y_{i1}(0) - Y_{i0}(0)|X_i, D_i=1] = \\mathbb{E}[Y_{i1}(0) - Y_{i0}(0)|X_i, D_i=0]\\quad a.s.$\n",
    "\n",
    "- **Overlap:** For some $\\epsilon > 0$, $P(D_i=1) > \\epsilon$ and $P(D_i=1|X_i) \\le 1-\\epsilon$ a.s.\n",
    "\n",
    "For a detailed explanation of the assumptions see e.g. [Sant'Anna and Zhao (2020)](https://doi.org/10.1016/j.jeconom.2020.06.003) or [Zimmert et al. (2018)](https://arxiv.org/abs/1809.01643).\n"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Panel Data (Repeated Outcomes)\n",
    "\n",
    "At first, we will consider two-period panel data, where we observe i.i.d. data $W_i = (Y_{i0}, Y_{i1}, D_i, X_i)$.\n",
    "\n",
    "### Data\n",
    "\n",
    "We will use the implemented data generating process `make_did_SZ2020` to generate data according to the simulation in [Sant'Anna and Zhao (2020)](https://doi.org/10.1016/j.jeconom.2020.06.003) (Section 4.1). \n",
    "\n",
    "In this example, we will use `dgp_tpye=4`, which corresponds to the misspecified settings in [Sant'Anna and Zhao (2020)](https://doi.org/10.1016/j.jeconom.2020.06.003) (other data generating processes are also available via the `dgp_type` parameter). In all settings the true ATTE is zero."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To specify a corresponding `DoubleMLData` object, we have to specify a single outcome `y`. For panel data, the outcome consists of the difference of \n",
    "\n",
    "$$\\Delta Y_i = Y_{i1}- Y_{i0}.$$\n",
    "\n",
    "This difference will then be defined as outcome in our `DoubleMLData` object. The data generating process `make_did_SZ2020` already specifies the outcome `y` accordingly."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "================== DoubleMLData Object ==================\n",
      "\n",
      "------------------ Data summary      ------------------\n",
      "Outcome variable: y\n",
      "Treatment variable(s): ['d']\n",
      "Covariates: ['X1', 'X2', 'X3', 'X4']\n",
      "Instrument variable(s): None\n",
      "No. Observations: 1000\n",
      "\n",
      "------------------ DataFrame info    ------------------\n",
      "<class 'pandas.core.frame.DataFrame'>\n",
      "RangeIndex: 1000 entries, 0 to 999\n",
      "Columns: 6 entries, X1 to d\n",
      "dtypes: float64(6)\n",
      "memory usage: 47.0 KB\n",
      "\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "from doubleml.datasets import make_did_SZ2020\n",
    "from doubleml import DoubleMLData\n",
    "\n",
    "np.random.seed(42)\n",
    "n_obs = 1000\n",
    "x, y, d = make_did_SZ2020(n_obs=n_obs, dgp_type=4, cross_sectional_data=False, return_type='array')\n",
    "dml_data = DoubleMLData.from_arrays(x=x, y=y, d=d)\n",
    "print(dml_data)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### ATTE Estimation\n",
    "\n",
    "To estimate the ATTE with panel data, we will use the `DoubleMLDID` class. \n",
    "\n",
    "As for all `DoubleML` classes, we have to specify learners, which have to be initialized first.\n",
    "Here, we will just rely on a tree based method. \n",
    "\n",
    "The learner `ml_g` is used to fit conditional expectations of the outcome $\\mathbb{E}[\\Delta Y_i|D_i=0, X_i]$, whereas the learner `ml_m` will be used to estimate the propensity score $P(D_i=1|X_i)$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from lightgbm import LGBMClassifier, LGBMRegressor\n",
    "\n",
    "n_estimators = 30\n",
    "ml_g = LGBMRegressor(n_estimators=n_estimators)\n",
    "ml_m = LGBMClassifier(n_estimators=n_estimators)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `DoubleMLDID` class can be used as any other `DoubleML` class. \n",
    "\n",
    "The score is set to `score='observational'`, since the we generated data where the treatment probability depends on the pretreatment covariates. Further, we will use `in_sample_normalization=True`, since normalization generally improved the results in our simulations (both `score='observational'` and `in_sample_normalization=True` are default values).\n",
    "\n",
    "After initialization, we have to call the `fit()` method to estimate the nuisance elements."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "================== DoubleMLDID Object ==================\n",
      "\n",
      "------------------ Data summary      ------------------\n",
      "Outcome variable: y\n",
      "Treatment variable(s): ['d']\n",
      "Covariates: ['X1', 'X2', 'X3', 'X4']\n",
      "Instrument variable(s): None\n",
      "No. Observations: 1000\n",
      "\n",
      "------------------ Score & algorithm ------------------\n",
      "Score function: observational\n",
      "DML algorithm: dml2\n",
      "\n",
      "------------------ Machine learner   ------------------\n",
      "Learner ml_g: LGBMRegressor(n_estimators=30)\n",
      "Learner ml_m: LGBMClassifier(n_estimators=30)\n",
      "Out-of-sample Performance:\n",
      "Learner ml_g0 RMSE: [[11.43627032]]\n",
      "Learner ml_g1 RMSE: [[11.66989604]]\n",
      "Learner ml_m RMSE: [[0.48873663]]\n",
      "\n",
      "------------------ Resampling        ------------------\n",
      "No. folds: 5\n",
      "No. repeated sample splits: 1\n",
      "Apply cross-fitting: True\n",
      "\n",
      "------------------ Fit summary       ------------------\n",
      "       coef   std err         t     P>|t|     2.5 %    97.5 %\n",
      "d  1.386988  1.827375  0.759006  0.447849 -2.194601  4.968577\n"
     ]
    }
   ],
   "source": [
    "from doubleml import DoubleMLDID\n",
    "dml_did = DoubleMLDID(dml_data, \n",
    "                      ml_g=ml_g, \n",
    "                      ml_m=ml_m,\n",
    "                      score='observational',\n",
    "                      in_sample_normalization=True,\n",
    "                      n_folds=5)\n",
    "\n",
    "dml_did.fit()\n",
    "print(dml_did)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As usual, confidence intervals at different levels can be obtained via"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "      5.0 %    95.0 %\n",
      "d -1.618776  4.392752\n"
     ]
    }
   ],
   "source": [
    "print(dml_did.confint(level=0.90))"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Coverage Simulation\n",
    "\n",
    "Here, we add a small coverage simulation to highlight the difference to the linear implementation of [Sant'Anna and Zhao (2020)](https://doi.org/10.1016/j.jeconom.2020.06.003). We generate multiple datasets, estimate the ATTE and collect the results (this may take some time). "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Iteration: 0/200\n",
      "Iteration: 20/200\n",
      "Iteration: 40/200\n",
      "Iteration: 60/200\n",
      "Iteration: 80/200\n",
      "Iteration: 100/200\n",
      "Iteration: 120/200\n",
      "Iteration: 140/200\n",
      "Iteration: 160/200\n",
      "Iteration: 180/200\n"
     ]
    }
   ],
   "source": [
    "n_rep = 200\n",
    "ATTE = 0.0\n",
    "\n",
    "ATTE_estimates = np.full((n_rep), np.nan)\n",
    "coverage = np.full((n_rep), np.nan)\n",
    "ci_length = np.full((n_rep), np.nan)\n",
    "\n",
    "np.random.seed(42)\n",
    "for i_rep in range(n_rep):\n",
    "    if (i_rep % int(n_rep/10)) == 0:\n",
    "        print(f'Iteration: {i_rep}/{n_rep}')\n",
    "    dml_data = make_did_SZ2020(n_obs=n_obs, dgp_type=4, cross_sectional_data=False)\n",
    "\n",
    "    dml_did = DoubleMLDID(dml_data, ml_g=ml_g, ml_m=ml_m, n_folds=5)\n",
    "    dml_did.fit()\n",
    "\n",
    "    ATTE_estimates[i_rep] = dml_did.coef.squeeze()\n",
    "    confint = dml_did.confint(level=0.95)\n",
    "    coverage [i_rep] = (confint['2.5 %'].iloc[0] <= ATTE) & (ATTE <= confint['97.5 %'].iloc[0])\n",
    "    ci_length[i_rep] = confint['97.5 %'].iloc[0] - confint['2.5 %'].iloc[0]"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us take a look at the corresponding coverage and the length of the confidence intervals."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Coverage: 0.925\n",
      "Average CI length: 5.32236455588136\n"
     ]
    }
   ],
   "source": [
    "print(f'Coverage: {coverage.mean()}')\n",
    "print(f'Average CI length: {ci_length.mean()}')"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here, we can observe that the coverage is still valid, since we did not rely on linear learners, so the setting is not misspecified in this example. \n",
    "\n",
    "If we know the conditional expectation is correctly specified (linear form), we can use this to obtain smaller confidence intervals but in many applications, we may want to safeguard against misspecification and use flexible models such as random forest or boosting.\n",
    "\n",
    "The distribution of the estimates takes the following form"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGdCAYAAAD60sxaAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAABa0klEQVR4nO3deXhU5d0+8HtmktmyTPaZhEx2IIQlQZaIomKNBlfQ2hf8+RZJLbZur74RrbQKVu0bUKRURbC0CLiBVYut1QhNCYqELWEnhC17MpONmck+ycz5/REZjSSQZZIzmbk/1zVXm5MzZ75jdHLnPN/neSSCIAggIiIi8iBSsQsgIiIiGm4MQERERORxGICIiIjI4zAAERERkcdhACIiIiKPwwBEREREHocBiIiIiDwOAxARERF5HC+xC3BFdrsdVVVV8PPzg0QiEbscIiIi6gNBENDY2IiIiAhIpZe/x8MA1IOqqiro9XqxyyAiIqIBKC8vR2Rk5GXPYQDqgZ+fH4Cuf4D+/v4iV0NERER9YbFYoNfrHb/HL4cBqAcXh738/f0ZgIiIiEaYvrSvsAmaiIiIPA4DEBEREXkcBiAiIiLyOOwBIiIij2Wz2dDR0SF2GdRHMpkMXl5eTlmihgGIiIg8UlNTEyoqKiAIgtilUD+o1WqEh4dDLpcP6joMQERE5HFsNhsqKiqgVqsRGhrKRW9HAEEQYLVaUVtbi+LiYowePfqKix1eDgMQERF5nI6ODgiCgNDQUKhUKrHLoT5SqVTw9vZGaWkprFYrlErlgK/FJmgiIvJYvPMz8gzmrk+36zjlKkREREQjCIfAiIiIvlNpasWFZuuwvV6gjxyjAjgEJwaXCEBr1qzBq6++CoPBgOTkZLzxxhuYPn36FZ+3ZcsW3HfffZgzZw62bdvmOC4IApYtW4b169fDZDLh2muvxdq1azF69OghfBdERDSSVZpacdNruWjrsA/bayq9pch5apaoIWjjxo148sknYTKZRKtBDKIHoK1btyIzMxPr1q1DamoqVq9ejfT0dBQVFSEsLKzX55WUlGDx4sW47rrrLvneK6+8gtdffx2bNm1CbGwsnn/+eaSnp+PkyZODapgiIiL3daHZirYOOx69MWFYAkmlqRVrdp7FhWZrv15v4cKF2LRp0yXH09PTkZ2dfdnnxsTE4Mknn8STTz7pODZv3jzcdtttfX79gXK1oCV6AFq1ahUWLVqEjIwMAMC6devwr3/9Cxs2bMCzzz7b43NsNhvuv/9+/P73v8c333zT7R+mIAhYvXo1nnvuOcyZMwcAsHnzZmi1Wmzbtg3z588f8vdEREQj16gAFWJDfMQu47Jmz56Nd955p9sxhUIxoGupVCqPnAknahO01WpFfn4+0tLSHMekUinS0tKQl5fX6/NefPFFhIWF4cEHH7zke8XFxTAYDN2uqdFokJqa2us129vbYbFYuj2IiIhclUKhgE6n6/YIDAyEIAh44YUXEBUVBYVCgYiICPzP//wPAGDWrFkoLS3F//7v/0IikThmwG3cuBEBAQGOa7/wwgtISUnBhg0bEBUVBV9fXzzyyCOw2Wx45ZVXoNPpEBYWhj/84Q/dalq1ahUmTpwIHx8f6PV6PPLII2hqagIA5ObmIiMjA2az2fHaL7zwAoCu38GLFy/GqFGj4OPjg9TUVOTm5g75P0NR7wDV1dXBZrNBq9V2O67VanHq1Kken7N792789a9/xeHDh3v8vsFgcFzjx9e8+L0fy8rKwu9///t+Vk9EnmDf+Xqs/vcZnDJ0/WEU6qfA/6aNwewJOk6hJpfzySef4I9//CO2bNmC8ePHw2Aw4MiRIwCATz/9FMnJyXjooYewaNGiy17n3Llz+PLLL5GdnY1z587h3nvvxfnz5zFmzBjs2rULe/bswS9+8QukpaUhNTUVQNcNjNdffx2xsbE4f/48HnnkETzzzDN46623cM0112D16tVYunQpioqKAAC+vr4AgMceewwnT57Eli1bEBERgb///e+YPXs2jh07NqS9u6IPgfVHY2Mjfv7zn2P9+vUICQlx2nWXLFmCzMxMx9cWiwV6vd5p1yeikcfS1oEnPjyEnUW1iA3xwc1JOkgkQJHBgoffL8C0mEC8df8UhPoNbNiBaDA+//xzR4C46Le//S2USiV0Oh3S0tLg7e2NqKgox6SioKAgyGQy+Pn5QafTXfb6drsdGzZsgJ+fH5KSknDjjTeiqKgIX3zxBaRSKcaOHYsVK1Zg586djgD0w76imJgYvPzyy/j1r3+Nt956C3K5HBqNBhKJpNtrl5WV4Z133kFZWRkiIiIAAIsXL0Z2djbeeecd/N///Z8z/nH1SNQAFBISAplMBqPR2O240Wjs8Ydz7tw5lJSU4M4773Qcs9u7uvW9vLxQVFTkeJ7RaER4eHi3a6akpPRYh0KhGPDYKRG5H3NLB37+1304V9eEJ9NGY3pM0A/u9ozCkXIT3v76HB7ceABbfzUDKrlM1HrJ89x4441Yu3Ztt2NBQUFobm7G6tWrERcXh9mzZ+O2227DnXfeCS+v/v26j4mJgZ+fn+NrrVYLmUzWbRFCrVaLmpoax9f//ve/kZWVhVOnTsFisaCzsxNtbW1oaWmBWq3u8XWOHTsGm82GMWPGdDve3t6O4ODgftXcX6L2AMnlckyZMgU5OTmOY3a7HTk5OZgxY8Yl5ycmJuLYsWM4fPiw43HXXXfhxhtvxOHDh6HX6xEbGwudTtftmhaLBfv27evxmkREP2Ru6cD89Xkorm/Gc7cnITU2+JKhrmR9AJ5OT8TpmkY8seUQbHZupknDy8fHBwkJCd0eQUFB0Ov1KCoqwltvvQWVSoVHHnkE119/fb93vPf29u72tUQi6fHYxZsQJSUluOOOOzBp0iR88sknyM/Px5o1awB09fv2pqmpCTKZDPn5+d1+txcWFuJPf/pTv2ruL9GHwDIzM/HAAw9g6tSpmD59OlavXo3m5mbHrLAFCxZg1KhRyMrKglKpxIQJE7o9/2Lj1g+PP/nkk3j55ZcxevRoxzT4iIgIzJ07d7jeFhGNQIIg4DefHEFZQwuW3TEe+qCe/2oFgNgQHzx+42i8tqMI63adw6M3JgxjpUS9U6lUuPPOO3HnnXfi0Ucfddw8uOqqqyCXy2Gz2Zz+mvn5+bDb7Xjttdccd4k++uijbuf09NqTJ0+GzWZDTU1Nj8vaDCXRA9C8efNQW1uLpUuXwmAwICUlBdnZ2Y4m5rKysn7v+/HMM8+gubkZDz30EEwmE2bOnIns7GyuAUREl/X+vjJknzDif9PGXDb8XHRVdCBmTwjHG/85g3unRELrz88Yd1BpanX512lvb79kYo+Xlxc+//xz2Gw2pKamQq1W47333oNKpUJ0dDSArqGtr7/+GvPnz4dCoXBaP21CQgI6Ojrwxhtv4M4778S3336LdevWdTsnJiYGTU1NyMnJQXJyMtRqNcaMGYP7778fCxYswGuvvYbJkyejtrYWOTk5mDRpEm6//Xan1NcT0QMQ0NUB/thjj/X4vStNhdu4ceMlxyQSCV588UW8+OKLTqiOiDzBKYMFL/7zJG5O0mJ6bFCfn3fP5FH45kwtXttehFfuTR7CCmmoBfrIofSWYs3Os8P2mkpvKQJ95P1+XnZ2drc+VwAYO3Ysli9fjuXLlyMzMxM2mw0TJ07EP//5T0c/zYsvvohf/epXiI+PR3t7OwTBOcO3ycnJWLVqFVasWIElS5bg+uuvR1ZWFhYsWOA455prrsGvf/1rzJs3D/X19Vi2bBleeOEFvPPOO3j55Zfx1FNPobKyEiEhIbj66qtxxx13OKW23kgEZ717N2KxWKDRaGA2m+Hv7y92OUQ0xOx2Afes3YO6pnb8Ye5EyL36d9d5+wkDNu4pwef/MxPjIzRDVCU5U1tbG4qLixEbG9ttdIB7gbm+3n52QP9+f7vEHSAiIjF9eqgSh8tNeP6OpH6HHwD4ybgwbD9pxCvZRdj0iyvvY0iua1SAioHEQ4g6C4yISGyWtg5kfVmIGXHBSAof2B1fL6kUc1IisOt0Lc7XNjm5QiIaCgxAROTR3sg5g6a2TtyfGjWo66TGBsNf6YX39pY5qTIiGkoMQETksarNrdi4pwR3Jkcg2Hdwi6HKvaSYNTYMfztYjhZrp5MqJKKhwgBERB5rzc6zUHrJcNuE8Cuf3Ac3JYahqb0T/zhc5ZTr0dDjPKCRx1k/MwYgIvJIFRdasGV/Oe6YFO60rSzC/JWYHBWATXkl/MXq4mSyrp/55VYpJtfU0tIC4NLVqvuLs8CIyCOt2XkWPgov3DL+8ptC9tdNiVq8ur0Ip41NGKvzu/ITSBReXl5Qq9Wora2Ft7d3vxfcpeEnCAJaWlpQU1ODgIAAR4gdKAYgIvI4laZWfHSwAvOn6aH0du5GppMiNfCRy/CvY9UMQC5MIpEgPDwcxcXFKC0tFbsc6oeAgIAr7mbfFwxARORxNuwuhspbhrRxWqdf20smxVXRgfjX0Spk3jzmyk8g0cjlcowePZrDYCOIt7f3oO/8XMQAREQexdzSgQ/3l+GWJK3T7/5cdHVs8HfDYI0Yo+VdIFcmlUq5T6SH4qAnEXmU9/aVosNmR7qTe39+aGKkBmq5DJ8frR6y1yCiwWEAIiKP0d5pwzvfFuO60aEIUPd/A8q+8pZJMSWqaxiMiFwTAxAReYx/HK5CXZMVt090zro/l5MaF4xztc04bWwc8tciov5jACIijyAIAjbuKUGKPgARw7DZ5cRRGii8pNh5qmbIX4uI+o8BiIg8QkGZCSeqLLglyfkzv3oi95IiKdwfuUW1w/J6RNQ/DEBE5BE27SlGuEaJZH3AsL3mpEgNDpY2cG8wIhfEAEREbq/G0oYvjhmQNk4LqUQybK87KTIAHTYBe8/XD9trElHfMAARkdv7cH85vGQS3DAmdFhfN1yjRJifAl+frhvW1yWiK2MAIiK31mmz48P9Zbg2PgQ+iuFd+1UikWDiKA1yi9gITeRqGICIyK3lFtXCYGnDTUOw7UVfJEcGoKS+BeUNLaK8PhH1jAGIiNzaB/vLEBfqg9gQH1Fef/wof0glwK7TnA1G5EoYgIjIbVWZWpFbVIOfJIaJVoNa7oXRYX7IO8dGaCJXwgBERG5ry4FyKLxkuCYuRNQ6EsP9sK+4HoIgiFoHEX2PAYiI3JLNLuCjA+W4Jj4YKvnQ7PreV4k6f9Q1WVFSzz4gIlfBAEREbumbM13NzzeKOPx10RitL6QSYH8xh8GIXAUDEBG5pY8OliMqSIU4kZqff0gt90JsiA/2FTeIXQoRfYcBiIjcTkOzFdtPGHHDmDBIhnHl58sZq/XDvvMMQESuggGIiNzOZ4crIQCYmSBu8/MPJYb7o9LUikpTq9ilEBEYgIjIzQiCgC0HyjElKhD+Km+xy3EYq/MDABzgMBiRS2AAIiK3cqLKgiJD47Dv+3Ul/kpvRAaq2AdE5CIYgIjIrXxaUAmNyhvJ+gCxS7lEos6PM8GIXAQDEBG5jQ6bHZ8drsS18cGQSV2j+fmHEsL8cK62GebWDrFLIfJ4LhGA1qxZg5iYGCiVSqSmpmL//v29nvvpp59i6tSpCAgIgI+PD1JSUvDuu+92O2fhwoWQSCTdHrNnzx7qt0FEIvvmTC3qm624zsWGvy4aHeYLADhaYRK3ECISPwBt3boVmZmZWLZsGQoKCpCcnIz09HTU1NT0eH5QUBB+97vfIS8vD0ePHkVGRgYyMjLw1VdfdTtv9uzZqK6udjw+/PDD4Xg7RCSiT/IrERWkRnSQWuxSeqTTKOGjkOFwmUnsUog8nugBaNWqVVi0aBEyMjKQlJSEdevWQa1WY8OGDT2eP2vWLNx9990YN24c4uPj8cQTT2DSpEnYvXt3t/MUCgV0Op3jERgYOBxvh4hEYm7twI6TRlw3OsRl1v75MalEgoRQXxwquyB2KUQeT9QAZLVakZ+fj7S0NMcxqVSKtLQ05OXlXfH5giAgJycHRUVFuP7667t9Lzc3F2FhYRg7diwefvhh1Nf33njY3t4Oi8XS7UFEI8sXx6rRabfjWhda+6cn8WG+OFRu4saoRCITNQDV1dXBZrNBq9V2O67VamEwGHp9ntlshq+vL+RyOW6//Xa88cYbuPnmmx3fnz17NjZv3oycnBysWLECu3btwq233gqbzdbj9bKysqDRaBwPvV7vnDdIRMNm26FKTIjQIFAtF7uUyxod5osLLR0oa+DGqERi8hK7gIHw8/PD4cOH0dTUhJycHGRmZiIuLg6zZs0CAMyfP99x7sSJEzFp0iTEx8cjNzcXN9100yXXW7JkCTIzMx1fWywWhiCiEaTK1Ir9xQ341Q3xYpdyRfGhXY3Qh8tNiA4Wf58yIk8lagAKCQmBTCaD0WjsdtxoNEKn0/X6PKlUioSEBABASkoKCgsLkZWV5QhAPxYXF4eQkBCcPXu2xwCkUCigUCgG/kaISFT/PFIFb5kU02Jcv9fPT+mNcI0Sh8pMmJMySuxyiDyWqENgcrkcU6ZMQU5OjuOY3W5HTk4OZsyY0efr2O12tLe39/r9iooK1NfXIzw8fFD1EpFr2na4EldFB0AtHxk3teNDfVHARmgiUYk+CywzMxPr16/Hpk2bUFhYiIcffhjNzc3IyMgAACxYsABLlixxnJ+VlYUdO3bg/PnzKCwsxGuvvYZ3330X//3f/w0AaGpqwtNPP429e/eipKQEOTk5mDNnDhISEpCeni7KeySioXPG2IjC6kZcE+/azc8/lBDmi5NVFrR19NyXSERDT/Q/l+bNm4fa2losXboUBoMBKSkpyM7OdjRGl5WVQSr9Pqc1NzfjkUceQUVFBVQqFRITE/Hee+9h3rx5AACZTIajR49i06ZNMJlMiIiIwC233IKXXnqJw1xEbmjb4Ur4KGRIccGtL3oTH+qLTruAIkOjS27ZQeQJJALnYl7CYrFAo9HAbDbD399f7HKIqBeCIOD6V3YiIcwXD13v+g3QF1k77cjYuB8vzpmA/746WuxyiNxGf35/iz4ERkQ0UMcqzSi/0Iqr44LFLqVf5F5S6IPUOF5pFrsUIo/FAEREI9bnR6uhUXljfIRG7FL6LSbYB0crGICIxMIAREQjkiAI+PxIFabFBLrkzu9XEhfig9PGRrR3shGaSAwMQEQ0IhWUmVBlbsOMETb8dVFsiI+jEZqIhh8DEBGNSJ8frUKg2huJupE5USE62AdSSVcfExENPwYgIhpx7HYBXxyrxvTYYEhH4PAX0NUIHRnIRmgisTAAEdGIc6jcBKOlHVfHBoldyqDEhrARmkgsDEBENOJkH69GgNobY7R+YpcyKLEhPigysBGaSAwMQEQ0ogiCgC+OGTA1OnDEDn9ddLER+rShSexSiDwOAxARjSjHKy2oNLVieuzInP31Q9HBakglwPEqDoMRDTcGICIaUb48Xg0/hRfGhY/s4S8AUHjJEBGgQmG1RexSiDwOAxARjRhdw1/VuCo6EF5S9/j4igpS40QVAxDRcHOPTxAi8ghnappQUt+C6SN89tcPRQepcaraArud+1ITDScGICIaMbafMEDlLcPEUSNv76/eRAf7oNlqQ/mFFrFLIfIoDEBENGJsP2HEpEgNvGXu89EVHawGAJzkMBjRsHKfTxEicmtGSxuOVpoxJTpQ7FKcKkAtR4Dam43QRMOMAYiIRoQdJ42QSoDJevcKQEBXH9AJBiCiYcUAREQjwvYTBowL94ev0kvsUpwuOtiHM8GIhhkDEBG5vMa2Duw5V+92w18XRQerYTC3wdRiFbsUIo/BAERELm/X6Vp02gVMddcAFOQDADjJYTCiYcMAREQu7z+FNYgKUiPUTyl2KUMiXKOEwkvKmWBEw4gBiIhcms0uYGdRDSZHBYhdypCRSiXQB6p4B4hoGDEAEZFLO1xuwoWWDrec/fVD+iA1TlU3il0GkcdgACIil7bzVA38lF4YHeYrdilDKipIjbM1Tei02cUuhcgjMAARkUvLOWXEpMgASKUSsUsZUlFBalhtdpTUN4tdCpFHYAAiIpdlMLehsLoRk/UBYpcy5PRBXVtiFHIYjGhYMAARkcvaWVQDqQRIjgwQu5Qh56f0RpCPHEUGBiCi4cAAREQua+epGozR+rnl6s890QeqUGjgTDCi4cAAREQuydppx7dn6zzi7s9FnAlGNHwYgIjIJRWUXUCz1YZkD+j/uSgqSI1KUyssbR1il0Lk9hiAiMglfX26FhqVN6KD1WKXMmyivmuEPs0+IKIhxwBERC4pt6gWE0ZpIJW49/T3H4oIUEEmleAUAxDRkGMAIiKXU9vYjpPVFiRHasQuZVh5y6SICFDiFBuhiYacSwSgNWvWICYmBkqlEqmpqdi/f3+v53766aeYOnUqAgIC4OPjg5SUFLz77rvdzhEEAUuXLkV4eDhUKhXS0tJw5syZoX4bROQku8/WAgAmeVAD9EX6QDZCEw0H0QPQ1q1bkZmZiWXLlqGgoADJyclIT09HTU1Nj+cHBQXhd7/7HfLy8nD06FFkZGQgIyMDX331leOcV155Ba+//jrWrVuHffv2wcfHB+np6Whraxuut0VEg7CrqBaxIT7QqLzFLmXY6YPUKDI0QhAEsUshcmuiB6BVq1Zh0aJFyMjIQFJSEtatWwe1Wo0NGzb0eP6sWbNw9913Y9y4cYiPj8cTTzyBSZMmYffu3QC67v6sXr0azz33HObMmYNJkyZh8+bNqKqqwrZt24bxnRHRQNjtAr4+U4eJozxr+OuiqEA1Gts7UW3mH2xEQ0nUAGS1WpGfn4+0tDTHMalUirS0NOTl5V3x+YIgICcnB0VFRbj++usBAMXFxTAYDN2uqdFokJqa2us129vbYbFYuj2ISBynDI1oaLZ6XP/PRfogFQCgyMhhMKKhJGoAqqurg81mg1ar7XZcq9XCYDD0+jyz2QxfX1/I5XLcfvvteOONN3DzzTcDgON5/blmVlYWNBqN46HX6wfztohoEL49WweFlxSjtX5ilyKKYF8FVN4yToUnGmKiD4ENhJ+fHw4fPowDBw7gD3/4AzIzM5Gbmzvg6y1ZsgRms9nxKC8vd16xRNQvu8/WIVHnB2/ZiPx4GjSpRILIQBX3BCMaYqJusBMSEgKZTAaj0djtuNFohE6n6/V5UqkUCQkJAICUlBQUFhYiKysLs2bNcjzPaDQiPDy82zVTUlJ6vJ5CoYBCoRjkuyGiwWrvtGFfcT3umRwpdimiigxUcy0goiEm6p9YcrkcU6ZMQU5OjuOY3W5HTk4OZsyY0efr2O12tLe3AwBiY2Oh0+m6XdNisWDfvn39uiYRDb9DZSa0ddgxwUMboC/SB6lwtqYJNjtnghENFdG3WM7MzMQDDzyAqVOnYvr06Vi9ejWam5uRkZEBAFiwYAFGjRqFrKwsAF39OlOnTkV8fDza29vxxRdf4N1338XatWsBABKJBE8++SRefvlljB49GrGxsXj++ecRERGBuXPnivU2iagPvj1bBz+ll0dtf9ETfaAaVpsdJfXNiA/1FbscIrckegCaN28eamtrsXTpUhgMBqSkpCA7O9vRxFxWVgap9PsbVc3NzXjkkUdQUVEBlUqFxMREvPfee5g3b57jnGeeeQbNzc146KGHYDKZMHPmTGRnZ0OpVA77+yOivvvmTB3GR/h71PYXPdH/YE8wBiCioSERuNrWJSwWCzQaDcxmM/z9/cUuh8gjWNo6kPL77fjFzFjclKi98hPc3K/fy0fGtTF4Mm2M2KUQjRj9+f3tmdMsiMjl7DvfALsATIjw7P6fizgTjGhoMQARkUvIO1ePUD8FtP4cqgY4E4xoqDEAEZFL2HOuDknhHHK+SB+kQml9M9o6bGKXQuSWGICISHQXmq04ZWjEOAYgB32gGnYBOFvTJHYpRG6JAYiIRLevuB4AMD6CAeiiyMCuPcHO1HAYjGgoMAARkejyztVD669AiC9XZL9ILfdCiK8cRQbeASIaCgxARCS6Pefq2f/Tg8hANU5zV3iiIcEARESiqmtqx5maJiRx+vslOBWeaOgwABGRqPadbwAA3gHqQWSgGpWmVjS3d4pdCpHbYQAiIlHtPV+PcI0SQT5ysUtxOd83QrMPiMjZGICISFR7z9cjUce7Pz0ZFaCCBF17ghGRczEAEZFoGpqtOFPThHHhfmKX4pKU3jJo/ZVshCYaAgxARCSa/cVd/T9cALF3owJVKGIAInI6BiAiEs2+4nqE+XH9n8vRB6p4B4hoCDAAEZFo9p1vQKKOw1+XExmohtHSDnNLh9ilELkVBiAiEoW5tQOF1RYkcvjrsi7OBDvNLTGInIoBiIhEcbCkAQK4/s+VRASoIJWACyISORkDEBGJYn9xA4J95AjzY//P5XjLpAjXqHCGfUBETsUARESiyDtfj7E6P0gkErFLcXmRgSoUGbkYIpEzMQAR0bBrsXbiZJWFDdB9FMmZYEROxwBERMPuUJkJnXaBK0D3kT5QjYZmK+qa2sUuhchtMAAR0bDbX9wAX4UXRn03w4kuLzJQDQC8C0TkRAxARDTs9hc3YKzWD1L2//SJVqOAl0yCM+wDInIaBiAiGlYdNjsOlV/AWPb/9JmXVIoIDbfEIHImBiAiGlYnqixo67CzAbqfIgNVXAuIyIkYgIhoWB0oboDCS4rYEB+xSxlRIgPVOG1shCAIYpdC5BYYgIhoWO0vaUBCmC+8ZPz46Q99oAqNbZ0wWjgTjMgZ+AlERMNGEAQcKG5g/88AcCYYkXMxABHRsDlX2wRTawfGahmA+ivMTwGFl5QBiMhJGICIaNgcKLkAqQQYHcYA1F9SqQSjAtgITeQsDEBENGwOlDQgJsQHKrlM7FJGpFGBnApP5CwMQEQ0bA4UN2AMh78GTB+oxhljE+x2zgQjGiwGICIaFkZLG8ovtCKRAWjA9EEqtHbYUGlqFbsUohHPJQLQmjVrEBMTA6VSidTUVOzfv7/Xc9evX4/rrrsOgYGBCAwMRFpa2iXnL1y4EBKJpNtj9uzZQ/02iOgyDpQ0AADGcAbYgOm/mwnGPiCiwRM9AG3duhWZmZlYtmwZCgoKkJycjPT0dNTU1PR4fm5uLu677z7s3LkTeXl50Ov1uOWWW1BZWdntvNmzZ6O6utrx+PDDD4fj7RBRLw6WXIDOX4lAtVzsUkasIB851HIZ+4CInED0ALRq1SosWrQIGRkZSEpKwrp166BWq7Fhw4Yez3///ffxyCOPICUlBYmJifjLX/4Cu92OnJycbucpFArodDrHIzAwcDjeDhH1Yn9xA8ZofcUuY0STSCTQf7ciNBENjqgByGq1Ij8/H2lpaY5jUqkUaWlpyMvL69M1Wlpa0NHRgaCgoG7Hc3NzERYWhrFjx+Lhhx9GfX19r9dob2+HxWLp9iAi52ls68ApgwWJOn+xSxnxIgNVOMUhMKJBEzUA1dXVwWazQavVdjuu1WphMBj6dI3f/OY3iIiI6BaiZs+ejc2bNyMnJwcrVqzArl27cOutt8Jms/V4jaysLGg0GsdDr9cP/E0R0SUOlZlgF9j/4wyRgWqcr21Ch80udilEI5qX2AUMxvLly7Flyxbk5uZCqVQ6js+fP9/x/ydOnIhJkyYhPj4eubm5uOmmmy65zpIlS5CZmen42mKxMAQROdHBkgb4K70QoVFe+WS6LH2QCh02AaX1zUjggpJEAybqHaCQkBDIZDIYjcZux41GI3Q63WWfu3LlSixfvhzbt2/HpEmTLntuXFwcQkJCcPbs2R6/r1Ao4O/v3+1BRM6zv+QCxur8IJFIxC5lxNM79gRrErkSopFN1AAkl8sxZcqUbg3MFxuaZ8yY0evzXnnlFbz00kvIzs7G1KlTr/g6FRUVqK+vR3h4uFPqJqK+s3bacbjsAhdAdBJ/lTcCVN6cCk80SKLPAsvMzMT69euxadMmFBYW4uGHH0ZzczMyMjIAAAsWLMCSJUsc569YsQLPP/88NmzYgJiYGBgMBhgMBjQ1df011NTUhKeffhp79+5FSUkJcnJyMGfOHCQkJCA9PV2U90jkyU5UmdHWaUci+3+cJjJQxZlgRIMkeg/QvHnzUFtbi6VLl8JgMCAlJQXZ2dmOxuiysjJIpd/ntLVr18JqteLee+/tdp1ly5bhhRdegEwmw9GjR7Fp0yaYTCZERETglltuwUsvvQSFQjGs742Iutb/UXhJERPsI3YpbiMyUI3Cas5WJRoMiSAI3FTmRywWCzQaDcxmM/uBiAbpoc0HUWlqxXO3J4lditvIOWXEht3FOPnibCi9ubEs0UX9+f0t+hAYEbkvQRBwoKQBY9n/41RRgWrYBeBsDRuhiQaKAYiIhsz5umZcaOnAWPb/OFXkdzPBuCAi0cAxABHRkDlQ3ACpBEgI4xYYzqSSy6D1V6DIwD4gooFiACKiIXOg5AKig32glos+38LtRAaqeQeIaBAYgIhoyOwrrmf/zxDRMwARDcqAAtD58+edXQcRuRmDuQ0VF1qRGM4ANBSiglSobWzHhWar2KUQjUgDCkAJCQm48cYb8d5776Gtrc3ZNRGRG9hf0gAAvAM0RNgITTQ4AwpABQUFmDRpEjIzM6HT6fCrX/0K+/fvd3ZtRDSC7S+uR0SAEgFqudiluKXwACW8pBKuCE00QAMKQCkpKfjTn/6EqqoqbNiwAdXV1Zg5cyYmTJiAVatWoba21tl1EtEIs7+Y6/8MJS+pFKMCVbwDRDRAg2qC9vLywj333IO//e1vWLFiBc6ePYvFixdDr9djwYIFqK6udladRDSCmFqsOG1sQqKOK6kPpa6ZYJwKTzQQgwpABw8exCOPPILw8HCsWrUKixcvxrlz57Bjxw5UVVVhzpw5zqqTiEaQAyUXAIAboA6xqEAVThsawR2NiPpvQItzrFq1Cu+88w6Kiopw2223YfPmzbjtttscm5bGxsZi48aNiImJcWatRDRCHChpQLCPHKF+3IB4KOmD1Gi22lBxoRX6ILXY5RCNKAMKQGvXrsUvfvELLFy4EOHh4T2eExYWhr/+9a+DKo6IRqZ9xfUYo/ODRCIRuxS3FvVd6CkyNDIAEfXTgALQjh07EBUV5bjjc5EgCCgvL0dUVBTkcjkeeOABpxRJRCNHU3snjldY8MA1MWKX4vaCfOTwVXihsNqCtCSt2OUQjSgD6gGKj49HXV3dJccbGhoQGxs76KKIaOTKL70AmyAgKYIN0ENNIpEgKogrQhMNxIACUG8Nd01NTVAqlYMqiIhGtr3n6xGg9kaEhp8FwyEqSI2T1ZwJRtRf/RoCy8zMBND1V8fSpUuhVn8/5myz2bBv3z6kpKQ4tUAiGln2nq9HIvt/hk1UkBrbTxrQarVBJZeJXQ7RiNGvAHTo0CEAXXeAjh07Brn8+xVe5XI5kpOTsXjxYudWSEQjRou1E8cqzFgwI1rsUjxGVLAadgE4U9OISZEBYpdDNGL0KwDt3LkTAJCRkYE//elP8PfnGD8RfS+/9AI67QKSwjVil+IxIgNVkEqAwmoLAxBRPwxoFtg777zj7DqIyA3sPV+PAJU3IgLY/zNcFF4yhGtUKKxmIzRRf/Q5AN1zzz3YuHEj/P39cc8991z23E8//XTQhRHRyLP3fAMSw9n/M9wiA1UoZCM0Ub/0OQBpNBrHh5pGw9vbRNRdi7UTR8pN+Dn7f4ZdVyO0EYIgMHwS9VGfA9APh704BEZEP3ag5GL/D3sDh1tUsBrm1g4YLe3QcfkBoj4Z0DpAra2taGlpcXxdWlqK1atXY/v27U4rjIhGlj1n6xCo9saoAJXYpXic6O+2wSjkzvBEfTagADRnzhxs3rwZAGAymTB9+nS89tprmDNnDtauXevUAoloZPjmTB0mRGg4BCOCEF8F1HIZTlYxABH11YACUEFBAa677joAwMcffwydTofS0lJs3rwZr7/+ulMLJCLX19BsRWG1BeNHsT9QDBKJBNFBajZCE/XDgAJQS0sL/Pz8AADbt2/HPffcA6lUiquvvhqlpaVOLZCIXF/euXoIACZw/y/RRAX74ATvABH12YACUEJCArZt24by8nJ89dVXuOWWWwAANTU1XByRyAN9e64OowJUCPZViF2Kx4oOVqOkrhkt1k6xSyEaEQYUgJYuXYrFixcjJiYGqampmDFjBoCuu0GTJ092aoFE5Pp2n6nDeN79EVVMsA8EgDvDE/XRgALQvffei7KyMhw8eBDZ2dmO4zfddBP++Mc/Oq04InJ95Q0tKGtowYQI9v+IaVSACjKphI3QRH00oK0wAECn00Gn03U7Nn369EEXREQjy7dn6yCVAON4B0hUci8pRgWocJKN0ER9MqAA1NzcjOXLlyMnJwc1NTWw2+3dvn/+/HmnFEdEri+3qAYJYb7wVQz47ylykuggNU5UmcUug2hEGNAn1i9/+Uvs2rULP//5zxEeHs51P4g8VIfNjt1n63HrBN2VT6YhFx3sg4/zy2GzC5BJ+blMdDkDCkBffvkl/vWvf+Haa691ShFr1qzBq6++CoPBgOTkZLzxxhu9DqetX78emzdvxvHjxwEAU6ZMwf/93/91O18QBCxbtgzr16+HyWTCtddei7Vr12L06NFOqZeIuuSXXkBTeyeS9QFil0LomgnW1mlHcV0zEsJ8xS6HyKUNqAk6MDAQQUFBTilg69atyMzMxLJly1BQUIDk5GSkp6ejpqamx/Nzc3Nx3333YefOncjLy4Ner8ctt9yCyspKxzmvvPIKXn/9daxbtw779u2Dj48P0tPT0dbW5pSaiahLblEtNCpvxIb4iF0KoSsAAWAfEFEfSARBEPr7pPfeew+fffYZNm3aBLVaPagCUlNTMW3aNLz55psAALvdDr1ej8cffxzPPvvsFZ9vs9kQGBiIN998EwsWLIAgCIiIiMBTTz2FxYsXAwDMZjO0Wi02btyI+fPnX/GaFosFGo0GZrOZ6xoRXUb66q8R5qfAI7MSxC6FvvP4hwW4d4oez96aKHYpRMOuP7+/BzQE9tprr+HcuXPQarWIiYmBt7d3t+8XFBT06TpWqxX5+flYsmSJ45hUKkVaWhry8vL6dI2WlhZ0dHQ47kgVFxfDYDAgLS3NcY5Go0Fqairy8vJ6DEDt7e1ob293fG2x8K8noisxmNtQZGjELUlasUuhH4gO8mEjNFEfDCgAzZ071ykvXldXB5vNBq22+weoVqvFqVOn+nSN3/zmN4iIiHAEHoPB4LjGj6958Xs/lpWVhd///vf9LZ/Io+06XQOpBJg0KkDsUugHYkLU2HmqBoIgcIIK0WUMKAAtW7bM2XUMyPLly7Flyxbk5uZCqVQO+DpLlixBZmam42uLxQK9Xu+MEoncVk7hd9PflZz+7kpiQnzQ0NIBg6UN4RqV2OUQuawBNUEDgMlkwl/+8hcsWbIEDQ0NALqGvn7YjHwlISEhkMlkMBqN3Y4bjcZLFln8sZUrV2L58uXYvn07Jk2a5Dh+8Xn9uaZCoYC/v3+3BxH1rtVqw9enazEl2jmTIch5YoO7GtKPVXAYjOhyBhSAjh49ijFjxmDFihVYuXIlTCYTAODTTz/t1s9zJXK5HFOmTEFOTo7jmN1uR05OjmN/sZ688soreOmll5CdnY2pU6d2+15sbCx0Ol23a1osFuzbt++y1ySivtt1uhZtnXZMiwkUuxT6kSAfOTQqbxznlhhElzWgAJSZmYmFCxfizJkz3YaebrvtNnz99df9vtb69euxadMmFBYW4uGHH0ZzczMyMjIAAAsWLOgWqlasWIHnn38eGzZsQExMDAwGAwwGA5qamgAAEokETz75JF5++WX84x//wLFjx7BgwQJEREQ4rXeJyNN9dcKAqCA1h1hckEQiQUywGscreQeI6HIGNHh/4MABvP3225ccHzVqVK+Nxr2ZN28eamtrsXTpUhgMBqSkpCA7O9vRxFxWVgap9PuctnbtWlitVtx7773drrNs2TK88MILAIBnnnkGzc3NeOihh2AymTBz5kxkZ2cPqk+IiLp02Oz4d6ERN4/j7C9XFRvigz3n6sUug8ilDSgAKRSKHqeKnz59GqGhof2+3mOPPYbHHnusx+/l5uZ2+7qkpOSK15NIJHjxxRfx4osv9rsWIrq8vefr0djWiakx7P9xVTEhPth2uAo1jW0I8+MffkQ9GdAQ2F133YUXX3wRHR0dALoCR1lZGX7zm9/gpz/9qVMLJCLXkn3cgDA/BWKCB7cIKg2di43QJyrZB0TUmwEFoNdeew1NTU0IDQ1Fa2srbrjhBiQkJMDPzw9/+MMfnF0jEbmITpsdX50wYGp0INeYcWGhfgr4KGTsAyK6jAENgWk0GuzYsQPffvstjhw5gqamJlx11VXdVl8mIvfz7bl61DVZcU1CiNil0GVIJBLEBvvgOFeEJupVvwOQ3W7Hxo0b8emnn6KkpKTrP7Tvpp5z5VEi9/ZpQQUiAlSI4+anLi8mxAcFZRfELoPIZfVrCEwQBNx111345S9/icrKSkycOBHjx49HaWkpFi5ciLvvvnuo6iQikTW3d+KrEwbMTAjhHzojQGyID6pMbahvar/yyUQeqF93gDZu3Iivv/4aOTk5uPHGG7t97z//+Q/mzp2LzZs3Y8GCBU4tkojEl33cgLYOO2YmBItdCvVBXIgvAOBYpRmzxoaJXA2R6+nXHaAPP/wQv/3tby8JPwDwk5/8BM8++yzef/99pxVHRK7j00MVGBfuh1BOqx4RtP4K+Cq8cJRbYhD1qF8B6OjRo5g9e3av37/11ltx5MiRQRdFRK6l2tyKPWfrMTOh/+t8kTgkEgniQnxwpMIkdilELqlfAaihocGxQnNPtFotLlxg0x2Ru3l/bxmU3jJcHcfFD0eS2FAf3gEi6kW/ApDNZoOXV+9tQzKZDJ2dnYMuiohcR1uHDe/vK8X1Y0Khlg9o5QwSSXyIL2ob22Ewt4ldCpHL6denmSAIWLhwIRQKRY/fb2/nbAMid/P50WpcaOlAehL3/hpp4kK7lis4UmGCTqMTuRoi19KvAPTAAw9c8RzOACNyH4Ig4J1vi5GiD0B4AHd+H2mCfOQIVHvjaIUJ6eMZgIh+qF8B6J133hmqOojIBeWXXsCJKgt+M3us2KXQAEgkEsSGsA+IqCcD2guMiDzDn3LOIDJQhUmRAWKXQgMUF+qLoxVmCIIgdilELoUBiIh6tPd8Pb45U4d7p0RCypWfR6z4UB+YWztQ1tAidilELoUBiIguIQgCXv2qCHEhPpgew6nvI9nFFaEPl5vELYTIxTAAEdElck/XIr/0An42Vc99v0Y4f5U3tP5KBiCiH2EAIqJu2jpseOnzkxgX7ofkSI3Y5ZATJIT64FCZSewyiFwKAxARdfOnnDMob2hBxjWxvPvjJhLCfHGiyoz2TpvYpRC5DAYgInI4WmHCn3edx92TI6EPUotdDjlJQpgfOmwCCqsbxS6FyGUwABERAKCpvROZHx1BVLAadyaHi10OOVF0sBpeMgkOl3GvRqKLGICICB02Ox59vwBVplY8fEM8vKT8aHAn3jIpYoN9cIiN0EQO/JQj8nCCIOD5bcex+2wdnkwbw6EvNxUf5stGaKIfYAAi8mDtnTY8/fERbDlQjkXXxWLiKM76clcJob4oa2hBfRM3rSYCGICIPFZtYzv++y/78NnhKjwyKx43jAkTuyQaQqPDuhZEPFJhErcQIhfRr81QiWjk67TZ8f6+MqzcXgSpRILnbk/CGK2f2GXREAv1UyBA5Y1DZSb8JFErdjlEomMAIvIQrVYbPj1Ugb/uLkZxbTNuTAzDvGl6+Cu9xS6NhoFEIkFCmC8OlnAmGBHAAETk1gRBQH7pBXx6qBL/PFKF5vZOTI0Jwi+ujUV8qK/Y5dEwG631w7ZDlei02eElYwcEeTYGICI3VGNpw0cHy/G3gxUobWhBqK8cNyWG4caxYQjzV4pdHolkrNYPrR02nDI0YgIb3snDMQARuZHS+mas23UeH+eXQyqRIDU2CD+fEY1x4f6QclsLjxcb4gMvmQQHSxoYgMjjMQARuQFrpx1v7zqH1/9zBj5yL9x7VSTSkrRQy/mfOH1P7iVFXIgP8stMWHit2NUQiYufjkQjXFl9CxZtPogzNY24Y1IEfnpVJORe7O+gno3R+uFASYPYZRCJjp+SRCPYobILmPvWtzC3deAPd0/EfdOjGH7ossaE+cFgbkOVqVXsUohEJfon5Zo1axATEwOlUonU1FTs37+/13NPnDiBn/70p4iJiYFEIsHq1asvOeeFF16ARCLp9khMTBzCd0Akjj1n63Dfn/cixFeO3981HjHBPmKXRCPAaG3X7L/8Uk6HJ88magDaunUrMjMzsWzZMhQUFCA5ORnp6emoqanp8fyWlhbExcVh+fLl0Ol0vV53/PjxqK6udjx27949VG+BSBRFhkY89G4+xuj88LvbkriWD/VZgFqOcI2SAYg8nqgBaNWqVVi0aBEyMjKQlJSEdevWQa1WY8OGDT2eP23aNLz66quYP38+FApFr9f18vKCTqdzPEJCQobqLRANO6OlDQvf2Y8QXzmevGkMh7yo3xLCfNkHRB5PtE9Oq9WK/Px8pKWlfV+MVIq0tDTk5eUN6tpnzpxBREQE4uLicP/996OsrOyy57e3t8NisXR7ELkim13Ao+8XwNppx9PpiVDJZWKXRCNQos4fhdUWNLZ1iF0KkWhEC0B1dXWw2WzQarvvSaPVamEwGAZ83dTUVGzcuBHZ2dlYu3YtiouLcd1116GxsbHX52RlZUGj0Tgeer1+wK9PNJT+8s155JdewKM3JiDIRy52OTRCjdP5wS6wD4g8m9vdO7/11lvxs5/9DJMmTUJ6ejq++OILmEwmfPTRR70+Z8mSJTCbzY5HeXn5MFZM1DdFhkas3F6E2yeFY1y4v9jl0Aim0ygRoPbG/mIOg5HnEm0doJCQEMhkMhiNxm7HjUbjZRuc+ysgIABjxozB2bNnez1HoVBctqeISGw2u4CnPjoMnb8SP5vCO5Q0OBKJBIk6P+xjACIPJtodILlcjilTpiAnJ8dxzG63IycnBzNmzHDa6zQ1NeHcuXMIDw932jWJhtsnBRU4XmXBL6+LY9MzOUWizh9Hyk1o67CJXQqRKET9JM3MzMT69euxadMmFBYW4uGHH0ZzczMyMjIAAAsWLMCSJUsc51utVhw+fBiHDx+G1WpFZWUlDh8+3O3uzuLFi7Fr1y6UlJRgz549uPvuuyGTyXDfffcN+/sjcoYWayde/aoIM+KCMUbrJ3Y55CYSdX7otAs4VGYSuxQiUYi6Fca8efNQW1uLpUuXwmAwICUlBdnZ2Y7G6LKyMkil32e0qqoqTJ482fH1ypUrsXLlStxwww3Izc0FAFRUVOC+++5DfX09QkNDMXPmTOzduxehoaHD+t6InOXtXedharHivukc+iLn0Qeq4aOQYX9xA2bEB4tdDtGwkwiCIIhdhKuxWCzQaDQwm83w92ezKYmnprEN16/YiZuTtPh/qdFil0Nu5tWvTkHpLcMHi64WuxQip+jP7282ExC5sL/uLoZEIsGclFFil0JuKFHnj4LSC7B22sUuhWjYMQARuShzSwfezSvFzUla+ChEHa0mNzUu3B9tnXYcqzSJXQrRsGMAInJRm/JK0GkTcOsE5y0LQfRDsSE+UHnLsOdsvdilEA07BiAiF9Ri7cSGb4sxa2woAtRc8ZmGhkwqwbhwP+w5xwBEnocBiMgFfXSgHJbWDtwxietX0dBKCtegoOwC2ju5HhB5FgYgIhcjCAI27y3F9NgghPopxS6H3FxShD/aO+1cD4g8DgMQkYvJO1eP87XNuHmc9sonEw1SdLAavgovDoORx2EAInIxm/eWIjJQxQ1PaVhIJV19QHnn6sQuhWhYMQARuRCDuQ07Thhx8zgtJBKJ2OWQh0gK1+BQmQmtVvYBkedgACJyIR/uL4PcS4qZo0PELoU8yPgIf3TaBRws5e7w5DkYgIhchN0u4KOD5ZgRHwy1nAsf0vCJDFQhUO2N3Wc5DEaegwGIyEXkna9HtbkNN4zhxr00vCQSCcZHaPDNaQYg8hwMQEQu4pP8CoRrlBgd5it2KeSBJo7S4GS1BfVN7WKXQjQsGICIXEBTeye+OF6N60aHsvmZRDFhlAYAOB2ePAYDEJEL+OJYNdo77LiOzc8kkiAfOfSBKuw+w2Ew8gwMQEQu4JP8CkwYpUGIr0LsUsiDTRilwddnaiEIgtilEA05BiAikVWZWrGvuAEzE3j3h8Q1cZQG1eY2FNc1i10K0ZBjACIS2b+OVsNbJsHUmECxSyEPNy7cH15SCafDk0dgACIS2T+OVGGyPpBr/5DolN4yjNX5YdfpWrFLIRpyDEBEIiqpa8axSjNmxAeLXQoRAGDSKA32nK1Heye3xSD3xgBEJKLPj1ZB6S3F5KgAsUshAgAk6wPQ2mHDwZILYpdCNKQYgIhE9M8j1ZgSHQiFl0zsUogAAFFBagT5yJFbVCN2KURDigGISCRnjI0oMjbi6jgOf5HrkEgkmDRKg51F7AMi98YARCSSL48boPKWITkyQOxSiLpJ1gfgbE0TKk2tYpdCNGQYgIhE8sWxakyOCoC3jP8ZkmuZOEoDqQTYxbtA5Mb4yUskgtL6ZpwyNGJ6bJDYpRBdwkfhhTFaP/YBkVtjACISwZfHDVB4STn8RS4rOTIAu8/WcTo8uS0GICIRfHGsGpMiNVB6c/YXuabJUQFosdqw73yD2KUQDQkGIKJhVmVqxdEKM6bHcvYXua6oIDVCfRX4zykOg5F7YgAiGmZfnTDASyrBVVz8kFyYRCJBSlQAdpw0cnd4cksMQETD7KsTRkwY5c+9v8jlXRUViEpTK87UNIldCpHTMQARDSNTixUHihswJZqzv8j1JYX7Q+ktxb8LjWKXQuR0DEBEw+g/p2pgEwRMiQ4UuxSiK5J7STFxlAY5hewDIvcjegBas2YNYmJioFQqkZqaiv379/d67okTJ/DTn/4UMTExkEgkWL169aCvSTSctp8wYHSYLwLVcrFLIeqTyfpAHCq7gPqmdrFLIXIqUQPQ1q1bkZmZiWXLlqGgoADJyclIT09HTU3Pf220tLQgLi4Oy5cvh06nc8o1iYZLW4cNuadrefeHRpTJUQEQBCCHs8HIzYgagFatWoVFixYhIyMDSUlJWLduHdRqNTZs2NDj+dOmTcOrr76K+fPnQ6FQOOWaRMPl27N1aOuwY2oM+39o5AhQyzFG54evjhvELoXIqUQLQFarFfn5+UhLS/u+GKkUaWlpyMvLG9Zrtre3w2KxdHsQOdv2E0ZEBCgxKkAldilE/TI1OhDfnKlDc3un2KUQOY1oAaiurg42mw1arbbbca1WC4NhYH9pDPSaWVlZ0Gg0joderx/Q6xP1xmYX8O9CI66K4vAXjTzTYoJgtdnx9WlujkruQ/QmaFewZMkSmM1mx6O8vFzsksjNHC43ob7Ziqmc/k4jkNZfiaggNbaf5HR4ch+ircQWEhICmUwGo7H7f1BGo7HXBuehuqZCoei1p4jIGf5daIS/0gujw3zFLoVoQKbGBGLHSSM6bHZ4y/i3M418ov1bLJfLMWXKFOTk5DiO2e125OTkYMaMGS5zTSJn2HHSiMlRgZBKJWKXQjQg02KC0NjWib3n68UuhcgpRF2LPzMzEw888ACmTp2K6dOnY/Xq1WhubkZGRgYAYMGCBRg1ahSysrIAdDU5nzx50vH/KysrcfjwYfj6+iIhIaFP1yQabiV1zThb04S7JkWIXQrRgEUHqaH1V+DL4wZcNzpU7HKIBk3UADRv3jzU1tZi6dKlMBgMSElJQXZ2tqOJuaysDFLp9zepqqqqMHnyZMfXK1euxMqVK3HDDTcgNze3T9ckGm7/LjTCWybBxEiN2KUQDZhEIsG0mCB8eawaL941Hl4cBqMRTiJwm99LWCwWaDQamM1m+Pv7i10OjXDz3s5Dh82Op9MTxS6FaFDO1TbhuW3H8f4vU3FtQojY5RBdoj+/vxnhiYaQqcWKgyUXMJnT38kNxIX4IMxPgX8dqxa7FKJBYwAiGkK5RbWwCQLX/yG3IJFIMD02CNnHDei02cUuh2hQGICIhtD2kwYkhPogyIebn5J7SI0NRkOzFfuLG8QuhWhQGICIhkh7pw27imo5/EVuJT7UB6F+CnzOYTAa4RiAiIbIvvMNaLbauPs7uRWJRILU2CB8cawaHRwGoxGMAYhoiPy70IhQPwWigtRil0LkVNcmhMDU0oHdZ+vELoVowBiAiIaAIAjYcbJr81OJhKs/k3uJDlIjMlCFzw5Vil0K0YAxABENgZPVFlSb2zj8RW5JIpFgRlwwvjphRKvVJnY5RAPCAEQ0BHacNEItl2Gczk/sUoiGxDXxIWjtsCHnFHeIp5GJAYhoCHx1woAUfQC3CyC3pdMokRDmi88OV4ldCtGA8NOZyMkqLrSgsLoRUzn8RW7umvhg7DxVgwvNVrFLIeo3BiAiJ/v3SSO8pBIk6wPELoVoSF0THwK7IOCfR3kXiEYeBiAiJ9t+0ojxEf5Qy73ELoVoSGlU3kjRB+DjgxVil0LUbwxARE5kbunAvvMNnP1FHuP60aE4WmnG2ZpGsUsh6hcGICIn2llUA5sgYEp0kNilEA2Lq6ID4avwwsf5XBOIRhYGICInyj5hQDw3PyUP4i2T4uq4YHxaUAGbXRC7HKI+YwAicpJWqw25RTWYFsO7P+RZbhgTgprGdnxzplbsUoj6jAGIyEm+PlOLtg47pjMAkYeJD/VFVJAKWw6Ui10KUZ8xABE5yVfHDYgKUiE8QCV2KUTDSiKR4MaxYfj3SSNqG9vFLoeoTxiAiJzA2mnHjkIjpvLuD3momQmhkEiATwo4JZ5GBgYgIifIO1+PxrZODn+Rx/JVemF6bBA+3F8GQWAzNLk+BiAiJ8g+boDWX4GoILXYpRCJ5idjw1Ba34J9xQ1il0J0RQxARIPUabPjqxMGTIsJgkQiEbscItGMC/dHhEaJ9/aWil0K0RUxABEN0t7zDWhotuLquGCxSyESlUQiwc1JWnx53IAaS5vY5RBdFgMQ0SD961gVtP4KxIX4iF0KkeiuHxMKb5kEH+wvE7sUostiACIahA6bHV8eNyA1NpjDX0QA1HIvzEwIwXt7S2HttItdDlGvGICIBmHPuXqYWjo4/EX0A7ck6VDXZMVXJwxil0LUKwYgokH419EqhGuUiAnm7C+ii/RBaoyP8Mc73xaLXQpRrxiAiAbI2mlHNoe/iHo0e7wOBWUmFJRdELsUoh4xABEN0K7TtbC0dWJGPIe/iH7squhAhGuUWP/1ebFLIeoRAxDRAH1aUIGYYDUXPyTqgVQiwa0TwvHVCQNK65vFLofoEgxARANgbu1ATmENrk0IEbsUIpd1/ZgQ+Cq8sGE3e4HI9TAAEQ1A9vFqdNrtuCaeAYioNwovGW5O0mLrwXI0NFvFLoeoG5cIQGvWrEFMTAyUSiVSU1Oxf//+y57/t7/9DYmJiVAqlZg4cSK++OKLbt9fuHAhJBJJt8fs2bOH8i2Qh/m0oBITIjQI8pGLXQqRS7tlvA4A8Nfd7AUi1yJ6ANq6dSsyMzOxbNkyFBQUIDk5Genp6aipqenx/D179uC+++7Dgw8+iEOHDmHu3LmYO3cujh8/3u282bNno7q62vH48MMPh+PtkAeoNLViX3EDh7+I+sBf6Y20cVps/LYE5pYOscshchA9AK1atQqLFi1CRkYGkpKSsG7dOqjVamzYsKHH8//0pz9h9uzZePrppzFu3Di89NJLuOqqq/Dmm292O0+hUECn0zkegYGBw/F2yAN8kl8BpbcU02KCxC6FaES4fWI4OmwCNu4pEbsUIgdRA5DVakV+fj7S0tIcx6RSKdLS0pCXl9fjc/Ly8rqdDwDp6emXnJ+bm4uwsDCMHTsWDz/8MOrr63uto729HRaLpduDqCd2u4AtB8owIy4YKrlM7HKIRoQAtRw/SQzDX789j8Y23gUi1yBqAKqrq4PNZoNWq+12XKvVwmDoeQl1g8FwxfNnz56NzZs3IycnBytWrMCuXbtw6623wmaz9XjNrKwsaDQax0Ov1w/ynZG72n22DlWmNtw4NkzsUohGlDsmhaPVasOG3SVil0IEwAWGwIbC/Pnzcdddd2HixImYO3cuPv/8cxw4cAC5ubk9nr9kyRKYzWbHo7y8fHgLphFjy4EyRAWpkBDmK3YpRCNKsK8CNyfp8Oevz3FGGLkEUQNQSEgIZDIZjEZjt+NGoxE6na7H5+h0un6dDwBxcXEICQnB2bNne/y+QqGAv79/twfRj9U3tWP7CSNmjQ3j1hdEAzAnJQJ2AXhrZ8+fxUTDSdQAJJfLMWXKFOTk5DiO2e125OTkYMaMGT0+Z8aMGd3OB4AdO3b0ej4AVFRUoL6+HuHh4c4pnDzSJwUVkEiA6xJCxS6FaETyV3rj9knh2JxXiipTq9jlkIcTfQgsMzMT69evx6ZNm1BYWIiHH34Yzc3NyMjIAAAsWLAAS5YscZz/xBNPIDs7G6+99hpOnTqFF154AQcPHsRjjz0GAGhqasLTTz+NvXv3oqSkBDk5OZgzZw4SEhKQnp4uynukkc9mF7BpTylSY4Phq/QSuxyiEeu2CeFQyWVY+VWR2KWQhxM9AM2bNw8rV67E0qVLkZKSgsOHDyM7O9vR6FxWVobq6mrH+ddccw0++OAD/PnPf0ZycjI+/vhjbNu2DRMmTAAAyGQyHD16FHfddRfGjBmDBx98EFOmTME333wDhUIhynukke/fhUZUmloxe0LvQ61EdGUquQz3TonEp4cqcbjcJHY55MEkgiAIYhfhaiwWCzQaDcxmM/uBCAAw7+08mFo68MJd48UuhWjEs9sF/G7bMWjU3vj7w9dCKmVPHTlHf35/i34HiMjVnagyY19xA+/+EDmJVCrBghkxOFJuxrbDlWKXQx6KAYjoCjZ+W4IQXzlXfiZyonHh/rg6Lgh/+KKQW2SQKBiAiC6jytSKvx+qxC1JOsh4m57IqX5+dQxarTb84YuTYpdCHogBiOgy3t51DipvGdLGaa98MhH1S5CPHPdNj8JHByuQd6737YqIhgIDEFEvaixt+HB/OWZP0HHfL6Ih8pPEMIzT+eHZT4+i1drzdkVEQ4EBiKgXb399Ht4yCZufiYaQVCLBL6+Lg8HcxqEwGlYMQEQ9qGlsw/t7S5E+QQe1nAsfEg2liAAV/l9qFN7bW4b/nDJe+QlETsAARNSD1746DW8vKW6dwO1TiIbDzeO0mKwPwNN/O4raxnaxyyEPwABE9CMnqyz46GA5fnpVJHwVvPtDNBwkEgkeuj4ONruAxz4oQKfNLnZJ5OYYgIh+QBAEvPT5SUQEqHDTuDCxyyHyKAFqOR6/aTQOlDRgRfYpscshN8cARPQD208akXe+Hv8vNQpeUv7nQTTcksL9cX9qNNZ/U4x/HqkSuxxyY7y/T/Qdc0sHntt2HFdFBWCyPkDscog81q0TdCiua8ZTHx1BuEaJqVyFnYYA/8Ql+s6Ln59AS3snHpwZB4mEqz4TieViP1BCmC8e3HQQ52qbxC6J3BADEBGA/5wy4pOCStx/dTSCfORil0Pk8bxlUvzvzWPgp/TCz/+6DxUXWsQuidwMAxB5vIoLLXjqoyNI0Qdg1phQscshou/4Krzw7OxE2OwC5v95L6pMrWKXRG6EAYg8WluHDQ9tzoe3TIpHZsVz6IvIxQT7KvC725LQ3mnH/D/vRXkD7wSRczAAkcey2wX85pOjOFfb9N2tdm+xSyKiHoT6KfDcbeNg7bTj7rf24ESVWeySyA0wAJFHEgQBz392HP84XIVfXR+HmGAfsUsiossI81di2Z1J8Fd54b/W5WHX6VqxS6IRjgGIPI4gCHjhHyfw/r4yLLo+DjPiQ8QuiYj6IEAtx3O3JWGszg8LN+zHmp1nIQiC2GXRCMUARB6lqb0Tv3o3H5vySvHLmbG4cSxXeyYaSVRyGZ66ZSzunjwKr35VhF9uOoi6Ju4dRv3HAEQeo7Dagjlv7sY3Z+rw1C1jcNM4rdglEdEASCUS/GyqHotvGYuDpReQ/sevseMkd5Gn/mEAIrdnbu3AC/84gTte3w1rpx0vz52AqdFcWZZopJsSHYjl90xEdLAaizYfxK/ePcip8tRnEoEDqJewWCzQaDQwm83w9/cXuxwaoLM1TXhvbyn+drAcNkHAPZMjcesEHbxkzP1E7kQQBOw9X49395aitcOGh66Lw6Lr4ziz0wP15/c3A1APGIBGFrtdgLm1A5WmVpwyNOJYhQm5RbUobWiBv9ILP0nU4pbxWgSqucIzkTtrsXZi26FKfHXCCLVChkXXxeG/r46GRsUg5CkYgAaJAWhomVs7cLamESV1LSi/0AKjpQ21jVZcaLGiqa0TrR02dNjssNkFCAAkAKRSCWQSCaRSQIKuxQqtnXa0d9rQ2NaJTvv3/xpr/RWYEKFBij4AkyIDIPfiHR8iT9LQbMXfD1Vg1+laeMukmDdNj/tTo5AQ5id2aTTEGIAGiQHIeayddhytMGF/SQOOlJtwpMIMg7nN8f0AtTeCfeTwV3rDV+kFlbcMCi8pvGRSSCVAV/wRIAiAXRBwMecIggBvmRTeMil8FDL4q7quMypADZVcJsZbJSIXY2qxIvuEATtP1cDS1okp0YGYmxKB2yaGI9hXIXZ5NAQYgAaJAWjgBEHAudpm5BbV4OvTtdhf3IC2TjtU3jLEhfogLsQH0cE+iAxUQadRQuHFsEJEQ6vDZsfBkgbsOl2LY5Vdq0hPiQ7ETeO0uG50CMbp/CGVchscd8AANEgMQP1jtws4VH4B2ccN2H7SiNL6FnjLJEgK98eEURokhfsjJtiHHzBEJDpLWwcOlDTgUJkJxyvNaO+0w1/phWmxQZisD0CyPgBJ4f68QzRCMQANEgPQldntAg6WXsAXx6rxxbFq1DS2I0DtjauiAjElKhDjR/nz7g4RuTRrpx1na5twssqM08YmnKttQovVBgAI8ZVjdJgf4sN8EBPsg6ggNSID1YgIUEKj8ubGyS6qP7+/vYapJnIDgiCgoMyEz49W4V9Hu0JPsI8c02KDkBobhDFhfrzLQ0QjhtxLiqRwfySFd/2itAsCjOY2lDa0oKyhBVWmVnx9ug5/s1SgvdPueJ7SS4pQPwXC/JUI8ZUj2FeBILUcAWpvBKjl0Ki8oVF5w0/pBf/v/tdX7sXPRxfDAESXJQgCjlSY8cWxanx+pApV5jYEqr0xPTYYV8cFYYzWD1L+JUREbkAqkSA8QIXwABWujgt2HBeErqU2ahvbUd9sRUOzFaYWK0wtHTCY23Da2ISm9k40tnWgrcPe47UlAHwUXvBTdj0CVHJo1N4IUHkj0EeOQLUcwb5yhPjKEeqrRJi/AiG+CsgYmoYMAxBdwmYXkF/a1dOTfbwaVeY2aFTemBYThAdnxiKRDYNE5EEkEgkC1HIEqOUYfYVzO212NFttaG7vRIu1E83tNrRYbWjp6ESr1fbd151ottpQ29iOkrpmNLV3wtLWgeZ2W7drySQShPoroA9UQR+ohj5IjZgQNWKCfRAX6sv1jQaJAYgAdK3N8+3ZOuQU1mBnUQ0amq0IVHtjSnQQfjEzlrMkiIj6wEsmhUYlHVA46bTbYWnthKnFigstHWhotqKhuR21TVacqLLgP0U1MLV0OM4P9pVjTJgfxmh9MUbnh0SdH8bq/OGr4K/2vnCJf0pr1qzBq6++CoPBgOTkZLzxxhuYPn16r+f/7W9/w/PPP4+SkhKMHj0aK1aswG233eb4viAIWLZsGdavXw+TyYRrr70Wa9euxejRV8runqOtw4bD5SbsOVeP3WdqcaTcDJsgICpIjZkJIZgSHYiEMF8ObxERDRMvqRRBPnIE+fS+an2r1QaDpQ1VplZUmVtReaEVOadq8O7eUsc6aZGBKiSF+yMx3B9J4X4YF+4PfaCaf8T+iOgBaOvWrcjMzMS6deuQmpqK1atXIz09HUVFRQgLC7vk/D179uC+++5DVlYW7rjjDnzwwQeYO3cuCgoKMGHCBADAK6+8gtdffx2bNm1CbGwsnn/+eaSnp+PkyZNQKpXD/RZFZ7cLKGtowbFKM45WmHCw9AKOV5rRYRPgp/TCOJ0/Mq6NwaTIAIT6ceonEZGrUslliA3xQWyIT7fj1k47qsytKKvvauAua2jBvuIGmFu77hipvGUYrfXFOJ0/xuj8MDrMFwlhvgjXKD12Rpvo0+BTU1Mxbdo0vPnmmwAAu90OvV6Pxx9/HM8+++wl58+bNw/Nzc34/PPPHceuvvpqpKSkYN26dRAEAREREXjqqaewePFiAIDZbIZWq8XGjRsxf/78K9Y0UqfBN7Z1oOJCK8oaWlBc14zztU04bWxCkaERrR1dY8thfgrEh/pijNYPieF+iApS8y4PEZEbEgQBptYOlDd8H4oqLnTdNbLaupq1Vd7fB6ro4K4+o8hAFSICVND5K+EzwobTRsw0eKvVivz8fCxZssRxTCqVIi0tDXl5eT0+Jy8vD5mZmd2OpaenY9u2bQCA4uJiGAwGpKWlOb6v0WiQmpqKvLy8HgNQe3s72tvbHV+bzV0rhVoslgG/t74QBAE2u4DOiw+bHR2ddrR32mG129HeYUNbhx2tHV2Nc83tnWhq74C5pRPm1q4x4romK+qa2lDT2N6tgU7pLUW4Rgmdvwq3jdUgMlCF6GAf+HfbHdmO1qamIX2PREQkHjmAeI0M8Ro/ILZrLzS7IKCuqR2VplYYLW2oNrfhTGUN9hS2oaG5Az+8K+KrkCHYV4FQXwUCfb6b5q/0hr/KC35Kb/govOCjkEHl5QWlXAqltxQKLxnkMink321r5CWVwEsmgZe0a4ujobzjdPH3dl/u7YgagOrq6mCz2aDVarsd12q1OHXqVI/PMRgMPZ5vMBgc3794rLdzfiwrKwu///3vLzmu1+v79kZc1BmxCyAiIhJBY2MjNBrNZc8ZWfe2hsiSJUu63VWy2+1oaGhAcHCwx46NDobFYoFer0d5efmIGkJ0N/w5uAb+HFwDfw6uYah/DoIgoLGxEREREVc8V9QAFBISAplMBqPR2O240WiETqfr8Tk6ne6y51/8X6PRiPDw8G7npKSk9HhNhUIBhaJ7829AQEB/3gr1wN/fnx80LoA/B9fAn4Nr4M/BNQzlz+FKd34ukg7Jq/eRXC7HlClTkJOT4zhmt9uRk5ODGTNm9PicGTNmdDsfAHbs2OE4PzY2Fjqdrts5FosF+/bt6/WaRERE5FlEHwLLzMzEAw88gKlTp2L69OlYvXo1mpubkZGRAQBYsGABRo0ahaysLADAE088gRtuuAGvvfYabr/9dmzZsgUHDx7En//8ZwBdzVVPPvkkXn75ZYwePdoxDT4iIgJz584V620SERGRCxE9AM2bNw+1tbVYunQpDAYDUlJSkJ2d7WhiLisrg1T6/Y2qa665Bh988AGee+45/Pa3v8Xo0aOxbds2xxpAAPDMM8+gubkZDz30EEwmE2bOnIns7GyPXANIDAqFAsuWLbtkWJGGF38OroE/B9fAn4NrcKWfg+jrABERERENN1F7gIiIiIjEwABEREREHocBiIiIiDwOAxARERF5HAYgGhbt7e1ISUmBRCLB4cOHxS7Ho5SUlODBBx9EbGwsVCoV4uPjsWzZMlitVrFLc3tr1qxBTEwMlEolUlNTsX//frFL8jhZWVmYNm0a/Pz8EBYWhrlz56KoqEjssjza8uXLHUvWiIkBiIbFM88806elycn5Tp06BbvdjrfffhsnTpzAH//4R6xbtw6//e1vxS7NrW3duhWZmZlYtmwZCgoKkJycjPT0dNTU1IhdmkfZtWsXHn30Uezduxc7duxAR0cHbrnlFjQ3N4tdmkc6cOAA3n77bUyaNEnsUjgNnobel19+iczMTHzyyScYP348Dh061Ou2JDQ8Xn31Vaxduxbnz58XuxS3lZqaimnTpuHNN98E0LXKvV6vx+OPP45nn31W5Oo8V21tLcLCwrBr1y5cf/31YpfjUZqamnDVVVfhrbfewssvv4yUlBSsXr1atHp4B4iGlNFoxKJFi/Duu+9CrVaLXQ59x2w2IygoSOwy3JbVakV+fj7S0tIcx6RSKdLS0pCXlydiZWQ2mwGA//6L4NFHH8Xtt9/e7b8LMYm+EjS5L0EQsHDhQvz617/G1KlTUVJSInZJBODs2bN44403sHLlSrFLcVt1dXWw2WyOFe0v0mq1OHXqlEhVkd1ux5NPPolrr7222+4BNPS2bNmCgoICHDhwQOxSHHgHiPrt2WefhUQiuezj1KlTeOONN9DY2IglS5aIXbJb6uvP4YcqKysxe/Zs/OxnP8OiRYtEqpxIHI8++iiOHz+OLVu2iF2KRykvL8cTTzyB999/36W2pGIPEPVbbW0t6uvrL3tOXFwc/uu//gv//Oc/IZFIHMdtNhtkMhnuv/9+bNq0aahLdWt9/TnI5XIAQFVVFWbNmoWrr74aGzdu7LbHHjmX1WqFWq3Gxx9/3G0T5gceeAAmkwmfffaZeMV5qMceewyfffYZvv76a8TGxopdjkfZtm0b7r77bshkMscxm80GiUQCqVSK9vb2bt8bLgxANGTKyspgsVgcX1dVVSE9PR0ff/wxUlNTERkZKWJ1nqWyshI33ngjpkyZgvfee0+UDxtPk5qaiunTp+ONN94A0DX8EhUVhccee4xN0MNIEAQ8/vjj+Pvf/47c3FyMHj1a7JI8TmNjI0pLS7sdy8jIQGJiIn7zm9+INhzJHiAaMlFRUd2+9vX1BQDEx8cz/AyjyspKzJo1C9HR0Vi5ciVqa2sd39PpdCJW5t4yMzPxwAMPYOrUqZg+fTpWr16N5uZmZGRkiF2aR3n00UfxwQcf4LPPPoOfnx8MBgMAQKPRQKVSiVydZ/Dz87sk5Pj4+CA4OFjUXiwGICI3t2PHDpw9exZnz569JHjyBvDQmTdvHmpra7F06VIYDAakpKQgOzv7ksZoGlpr164FAMyaNavb8XfeeQcLFy4c/oLIZXAIjIiIiDwOuyCJiIjI4zAAERERkcdhACIiIiKPwwBEREREHocBiIiIiDwOAxARERF5HAYgIiIi8jgMQERERORxGICIiIjI4zAAERERkcdhACIiIiKPwwBEREREHuf/A3rbewsnbGP1AAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "\n",
    "df_pa = pd.DataFrame(ATTE_estimates, columns=['Estimate'])\n",
    "g = sns.kdeplot(df_pa, fill=True)\n",
    "plt.show()"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Repeated Cross-Sectional Data\n",
    "\n",
    "For repeated cross-sectional data, we assume that we observe i.i.d. data $W_i = (Y_{i}, D_i, X_i, T_i)$. \n",
    "\n",
    "Here $Y_i = T_i Y_{i1} + (1-T_i)Y_{i0}$ corresponds to the outcome of unit $i$ which is observed at time $T_i$.\n",
    "\n",
    "### Data\n",
    "\n",
    "As for panel data, we will use the implemented data generating process `make_did_SZ2020` to generate data according to the simulation in [Sant'Anna and Zhao (2020)](https://doi.org/10.1016/j.jeconom.2020.06.003) (Section 4.2). \n",
    "\n",
    "In this example, we will use `dgp_tpye=4`, which corresponds to the misspecified settings in [Sant'Anna and Zhao (2020)](https://doi.org/10.1016/j.jeconom.2020.06.003) (other data generating processes are also available via the `dgp_type` parameter). In all settings the true ATTE is zero."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In contrast to other `DoubleMLData` objects, we have to specify which column corresponds to our time variable $T$.\n",
    "\n",
    "The time variable can be simply set via the argument `t`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "================== DoubleMLData Object ==================\n",
      "\n",
      "------------------ Data summary      ------------------\n",
      "Outcome variable: y\n",
      "Treatment variable(s): ['d']\n",
      "Covariates: ['X1', 'X2', 'X3', 'X4']\n",
      "Instrument variable(s): None\n",
      "Time variable: t\n",
      "No. Observations: 1000\n",
      "\n",
      "------------------ DataFrame info    ------------------\n",
      "<class 'pandas.core.frame.DataFrame'>\n",
      "RangeIndex: 1000 entries, 0 to 999\n",
      "Columns: 7 entries, X1 to t\n",
      "dtypes: float64(7)\n",
      "memory usage: 54.8 KB\n",
      "\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "from doubleml.datasets import make_did_SZ2020\n",
    "from doubleml import DoubleMLData\n",
    "\n",
    "np.random.seed(42)\n",
    "n_obs = 1000\n",
    "x, y, d, t = make_did_SZ2020(n_obs=n_obs, dgp_type=4, cross_sectional_data=True, return_type='array')\n",
    "dml_data = DoubleMLData.from_arrays(x=x, y=y, d=d, t=t)\n",
    "print(dml_data)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### ATTE Estimation\n",
    "\n",
    "To estimate the ATTE with panel data, we will use the `DoubleMLDIDCS` class. \n",
    "\n",
    "As for all `DoubleML` classes, we have to specify learners, which have to be initialized first.\n",
    "Here, we will just rely on a tree based method. \n",
    "\n",
    "The learner `ml_g` is used to fit conditional expectations of the outcome $\\mathbb{E}[\\Delta Y_i| D_i=d, T_i =t, X_i]$ for all combinations of $d,t\\in\\{0,1\\}$, whereas the learner `ml_m` will be used to estimate the propensity score $P(D_i=1|X_i)$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "from lightgbm import LGBMClassifier, LGBMRegressor\n",
    "\n",
    "n_estimators = 30\n",
    "ml_g = LGBMRegressor(n_estimators=n_estimators)\n",
    "ml_m = LGBMClassifier(n_estimators=n_estimators)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `DoubleMLDIDCS` class can be used as any other `DoubleML` class. \n",
    "\n",
    "The score is set to `score='observational'`, since the we generated data where the treatment probability depends on the pretreatment covariates. Further, we will use `in_sample_normalization=True`, since normalization generally improved the results in our simulations (both `score='observational'` and `in_sample_normalization=True` are default values).\n",
    "\n",
    "After initialization, we have to call the `fit()` method to estimate the nuisance elements."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "================== DoubleMLDIDCS Object ==================\n",
      "\n",
      "------------------ Data summary      ------------------\n",
      "Outcome variable: y\n",
      "Treatment variable(s): ['d']\n",
      "Covariates: ['X1', 'X2', 'X3', 'X4']\n",
      "Instrument variable(s): None\n",
      "Time variable: t\n",
      "No. Observations: 1000\n",
      "\n",
      "------------------ Score & algorithm ------------------\n",
      "Score function: observational\n",
      "DML algorithm: dml2\n",
      "\n",
      "------------------ Machine learner   ------------------\n",
      "Learner ml_g: LGBMRegressor(n_estimators=30)\n",
      "Learner ml_m: LGBMClassifier(n_estimators=30)\n",
      "Out-of-sample Performance:\n",
      "Learner ml_g_d0_t0 RMSE: [[15.02897287]]\n",
      "Learner ml_g_d0_t1 RMSE: [[26.5602727]]\n",
      "Learner ml_g_d1_t0 RMSE: [[27.62403053]]\n",
      "Learner ml_g_d1_t1 RMSE: [[44.06834315]]\n",
      "Learner ml_m RMSE: [[0.47761563]]\n",
      "\n",
      "------------------ Resampling        ------------------\n",
      "No. folds: 5\n",
      "No. repeated sample splits: 1\n",
      "Apply cross-fitting: True\n",
      "\n",
      "------------------ Fit summary       ------------------\n",
      "       coef   std err         t     P>|t|     2.5 %     97.5 %\n",
      "d  5.096741  5.225034  0.975447  0.329339 -5.144137  15.337619\n"
     ]
    }
   ],
   "source": [
    "from doubleml import DoubleMLDIDCS\n",
    "dml_did = DoubleMLDIDCS(dml_data,\n",
    "                        ml_g=ml_g,\n",
    "                        ml_m=ml_m,\n",
    "                        score='observational',\n",
    "                        in_sample_normalization=True,\n",
    "                        n_folds=5)\n",
    "\n",
    "dml_did.fit()\n",
    "print(dml_did)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As usual, confidence intervals at different levels can be obtained via"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "      5.0 %     95.0 %\n",
      "d -3.497674  13.691157\n"
     ]
    }
   ],
   "source": [
    "print(dml_did.confint(level=0.90))"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Coverage Simulation\n",
    "\n",
    "Again, we add a small coverage simulation to highlight the difference to the linear implementation of [Sant'Anna and Zhao (2020)](https://doi.org/10.1016/j.jeconom.2020.06.003). We generate multiple datasets, estimate the ATTE and collect the results (this may take some time). "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Iteration: 0/200\n",
      "Iteration: 20/200\n",
      "Iteration: 40/200\n",
      "Iteration: 60/200\n",
      "Iteration: 80/200\n",
      "Iteration: 100/200\n",
      "Iteration: 120/200\n",
      "Iteration: 140/200\n",
      "Iteration: 160/200\n",
      "Iteration: 180/200\n"
     ]
    }
   ],
   "source": [
    "n_rep = 200\n",
    "ATTE = 0.0\n",
    "\n",
    "ATTE_estimates = np.full((n_rep), np.nan)\n",
    "coverage = np.full((n_rep), np.nan)\n",
    "ci_length = np.full((n_rep), np.nan)\n",
    "\n",
    "np.random.seed(42)\n",
    "for i_rep in range(n_rep):\n",
    "    if (i_rep % int(n_rep/10)) == 0:\n",
    "        print(f'Iteration: {i_rep}/{n_rep}')\n",
    "    dml_data = make_did_SZ2020(n_obs=n_obs, dgp_type=4, cross_sectional_data=True)\n",
    "\n",
    "    dml_did = DoubleMLDIDCS(dml_data, ml_g=ml_g, ml_m=ml_m, n_folds=5)\n",
    "    dml_did.fit()\n",
    "\n",
    "    ATTE_estimates[i_rep] = dml_did.coef.squeeze()\n",
    "    confint = dml_did.confint(level=0.95)\n",
    "    coverage [i_rep] = (confint['2.5 %'].iloc[0] <= ATTE) & (ATTE <= confint['97.5 %'].iloc[0])\n",
    "    ci_length[i_rep] = confint['97.5 %'].iloc[0] - confint['2.5 %'].iloc[0]"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us take a look at the corresponding coverage and the length of the confidence intervals."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Coverage: 0.94\n",
      "Average CI length: 23.34858240261807\n"
     ]
    }
   ],
   "source": [
    "print(f'Coverage: {coverage.mean()}')\n",
    "print(f'Average CI length: {ci_length.mean()}')"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As for panel data the coverage is still valid, since we did not rely on linear learners, so the setting is not misspecified in this example. \n",
    "\n",
    "If we know the conditional expectation is correctly specified (linear form), we can use this to obtain smaller confidence intervals but in many applications, we may want to safeguard against misspecification and use flexible models such as random forest or boosting.\n",
    "\n",
    "The distribution of the estimates takes the following form"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGdCAYAAAD60sxaAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAABZQUlEQVR4nO3dd3yTdeIH8E9Gk3Ql3UlbWrpbRmmVUVpwcFSK4MB1yHkyjvPcq+opnALq+at4h4cKJ+op6CnCcSoqYrVWGUoZHWwos7vpbtI9kuf3BxLtUaD7yfi8X6+8tE+/ST8JNPnwPN/n+0gEQRBARERE5ECkYgcgIiIiGmosQERERORwWICIiIjI4bAAERERkcNhASIiIiKHwwJEREREDocFiIiIiBwOCxARERE5HLnYAayR2WxGWVkZ3N3dIZFIxI5DREREPSAIAhoaGhAQEACp9NL7eFiAulFWVoagoCCxYxAREVEfFBcXY9iwYZccwwLUDXd3dwDnXkC1Wi1yGiIiIuoJo9GIoKAgy+f4pbAAdeP8YS+1Ws0CREREZGN6Mn2Fk6CJiIjI4bAAERERkcNhASIiIiKHwzlARETksEwmEzo6OsSOQT0kk8kgl8sHZIkaFiAiInJIjY2NKCkpgSAIYkehXnBxcYG/vz8UCkW/HocFiIiIHI7JZEJJSQlcXFzg6+vLRW9tgCAIaG9vR1VVFc6ePYvIyMjLLnZ4KSxARETkcDo6OiAIAnx9feHs7Cx2HOohZ2dnODk5obCwEO3t7VCpVH1+LE6CJiIih8U9P7anP3t9ujzOgDwKERERkQ3hITAiIqKflda3oK6pfch+nqerAoEePAQnBhYgIiIinCs/U1dsQ2uHech+pspJiswnrhW1BK1btw6PPfYY6uvrRcsgBhYgIiIiAHVN7WjtMOPBKRFDUkhK61uw+odTqGtq79XPmz9/Pt5///0LtqekpCA9Pf2S9w0JCcFjjz2Gxx57zLJt9uzZmDFjRo9/fl9ZW9FiASIiIvqVQA9nhPq4ih3jkqZPn461a9d22aZUKvv0WM7Ozg55JhwnQRMREdkYpVIJnU7X5ebp6QlBELBs2TIEBwdDqVQiICAAjzzyCADg2muvRWFhIR5//HFIJBLLGXDr1q2Dh4eH5bGXLVuG+Ph4vPfeewgODoabmxseeOABmEwmvPLKK9DpdPDz88NLL73UJdOrr76K2NhYuLq6IigoCA888AAaGxsBANu2bcOCBQtgMBgsP3vZsmUAgLa2Njz55JMIDAyEq6srEhISsG3btkF/DUUvQKtXr0ZISAhUKhUSEhKwd+/eS47ftGkTYmJioFKpEBsbi61bt14w5tixY7jpppug0Wjg6uqK8ePHo6ioaLCeAhH1kNks4IfjlVj2xRFM+8d2XPHCt7jihW8x5e/b8JfPDuGbI3p0moZu/gWRvfnkk0/wj3/8A2+99RZOnjyJzZs3IzY2FgDw6aefYtiwYXjhhRdQXl6O8vLyiz7O6dOn8fXXXyM9PR0ff/wx3n33XcycORMlJSXYvn07li9fjmeffRZ79uyx3EcqleL111/HkSNH8P777+P777/Hn//8ZwBAUlISVq5cCbVabfnZTz75JADgoYceQlZWFjZs2ICDBw/ijjvuwPTp03Hy5MlBfKVEPgS2ceNGpKamYs2aNUhISMDKlSuRkpKC/Px8+Pn5XTB+165dmDNnDtLS0nDDDTdg/fr1mDVrFnJzczF69GgA5/7QJk+ejIULF+L555+HWq3GkSNH+rVYEhH1jyAIyDhagVczTuC4vgFatRIj/TW4IsgTkAC1Te34Ib8SH+0pQpivK/6cEo2UUTqu0UJ0EVu2bIGbm1uXbYsXL4ZKpYJOp0NycjKcnJwQHByMCRMmAAC8vLwgk8ng7u4OnU53ycc3m81477334O7ujpEjR2LKlCnIz8/H1q1bIZVKER0djeXLl+OHH35AQkICAHSZVxQSEoK//vWvuO+++/DPf/4TCoUCGo0GEomky88uKirC2rVrUVRUhICAAADAk08+ifT0dKxduxb/93//NxAvV7dELUCvvvoq7rnnHixYsAAAsGbNGnz11Vd477338Mwzz1ww/rXXXsP06dPx1FNPAQBefPFFZGRkYNWqVVizZg0A4C9/+QtmzJiBV155xXK/8PDwIXg2RNQdY2sHnvnkILYe0mNUgBpLbxyJGJ2627Fnqhrxn+xi3PdhLq4bocWK2XFQq5yGODGR9ZsyZQrefPPNLtu8vLzQ1NSElStXIiwsDNOnT8eMGTNw4403Qi7v3cd9SEgI3N3dLV9rtVrIZLIuixBqtVpUVlZavv7uu++QlpaG48ePw2g0orOzE62trWhuboaLi0u3P+fQoUMwmUyIiorqsr2trQ3e3t69ytxboh0Ca29vR05ODpKTk38JI5UiOTkZWVlZ3d4nKyury3jg3Kz38+PNZjO++uorREVFISUlBX5+fkhISMDmzZsH7XkQ0cUd1xtx4xs/Ylt+FR6dGolnZ168/ABAmK8bnrl+BJ64Lgq7TlfjplU/4mRFwxAmJrINrq6uiIiI6HLz8vJCUFAQ8vPz8c9//hPOzs544IEHcPXVV/f6ivdOTl3/4SGRSLrdZjafO2RdUFCAG264AWPGjMEnn3yCnJwcrF69GsC5z/uLaWxshEwmQ05ODvbv32+5HTt2DK+99lqvMveWaAWouroaJpMJWq22y3atVgu9Xt/tffR6/SXHV1ZWorGxES+//DKmT5+Ob7/9FrfccgtuvfVWbN++/aJZ2traYDQau9yIqH8Olxow+63dkAD4v1tiMTGs5/+aGxfihRdnjYbJJOCOt7JwXM/fSaKecnZ2xo033ojXX38d27ZtQ1ZWFg4dOgQAUCgUMJlMA/4zc3JyYDabsWLFCkycOBFRUVEoKyvrMqa7n33FFVfAZDKhsrLygkJ3ucN0/WVXp8Gfb6I333wzHn/8cQBAfHw8du3ahTVr1uCaa67p9n5paWl4/vnnhywnkb07VGLAXe/uhq+bEouuHwFXZe/favw1zlhywyi8tPUo7npnDzbem4gIP7fL35Gon0rrW6z+57S1tV2ws0Aul2PLli0wmUxISEiAi4sLPvzwQzg7O2P48OEAzh3a2rFjB+68804olUr4+Pj06zmcFxERgY6ODrzxxhu48cYb8dNPP1mmppwXEhKCxsZGZGZmIi4uDi4uLoiKisJdd92FuXPnYsWKFbjiiitQVVWFzMxMjBkzBjNnzhyQfN0RrQD5+PhAJpOhoqKiy/aKioqLtj6dTnfJ8T4+PpDL5Rg5cmSXMSNGjMCPP/540SyLFi1Camqq5Wuj0YigoKBePR8iOqe4thlz39sDP3cVFl0fAxdF399m3FRyLJoxAi99dRS/e2c3vnx4MrRqntBAg8PTVQGVkxSrfzg1ZD9T5SSFp6ui1/dLT0+Hv79/l23R0dF4+eWX8fLLLyM1NRUmkwmxsbH48ssvLfNpXnjhBdx7770IDw9HW1sbBEEYkOcRFxeHV199FcuXL8eiRYtw9dVXIy0tDXPnzrWMSUpKwn333YfZs2ejpqYGS5cuxbJly7B27Vr89a9/xRNPPIHS0lL4+Phg4sSJuOGGGwYk28VIhIF69n2QkJCACRMm4I033gBwbg9OcHAwHnrooW4nQc+ePRvNzc348ssvLduSkpIwZswYS9NMSkpCeHg4/v3vf1vG3HLLLXB2dsb69et7lMtoNEKj0cBgMECtvvh8BSLqqqmtE7e+uQv1ze148ebRcB+gCcx1ze14bvNhBHm5YOO9E6GUywbkcclxtba24uzZswgNDe1yljCvBWb9LvZnB/Tu81vUQ2CpqamYN28exo0bhwkTJmDlypVoamqynBU2d+5cBAYGIi0tDQDw6KOP4pprrsGKFSswc+ZMbNiwAdnZ2Xj77bctj/nUU09h9uzZuPrqqzFlyhSkp6fjyy+/HJJFlYgcmSAIeOI/+1FU04znbxo1YOUHADxdFHgsOQovbDmC5zYfxvLbxvAUeRoUgR7OLCQOQtQCNHv2bFRVVWHJkiXQ6/WIj49Henq6ZaJzUVFRl1PukpKSsH79ejz77LNYvHgxIiMjsXnzZssaQMC5vT1r1qxBWloaHnnkEURHR+OTTz7B5MmTh/z5ETmSf+8uRPqRCjxxXRSCvLo/5bU/IvzcsHByKNZsP4MJod64feywAf8ZROQ4RD0EZq14CIyod05XNWLmaztxVZQv/jApdFB/1pvbTiGnsA7pj109KEWLHMOlDqOQdRuoQ2CiXwqDiGxbh8mMxzfsh5erAr+bEDzoP29eUghclHI8vnE/TGb++42I+oYFiIj65V87z+JwmQH3XxsBldPgT052UcjxwDXhyCmswzs7zwz6zyP7xoMgtmeg/sxYgIioz8rqW/Ba5glMH6Ub0jV6YvzVuD7WH//IOIHCmqYh+7lkP2Syc2X9UqsUk3Vqbm4GcOFq1b1lVwshEtHQenHLUTg7yXCbCBOS7xg7DPsKarHo00P46I8JPCuMekUul8PFxQVVVVVwcnLqcsINWSdBENDc3IzKykp4eHhYSmxfsQARUZ/sPFmFrw/r8eCUiH4tdthXKicZ/jApFMvTj+OT3FKeFUa9IpFI4O/vj7Nnz6KwsFDsONQLHh4eA3KZDBYgIuo1k1nAi1uOYoTOHZPCB/eKzZcSH+SBSRHeeOmro0ge4QcPl96vqEuOS6FQIDIykofBbIiTk1O/9/ycxwJERL32xYFSnKhoxAs3jRL90NNdCcPxxH8OYMW3J/DirNGXvwPRr0ilUp4G76B40JOIeqW904wV357A+BBPRGrdxY4DTxcFbh87DB/tKcThUoPYcYjIRrAAEVGvbNhXhNK6Ftwx1nouGDxtlBaBHs547vPDMHNtICLqARYgIuqx1g4TXs88hasifaxqFWa5VIr5SSHIK6rH5wdKxY5DRDaABYiIeuy/OSWobWrDrCsCxY5ygZEBGiSEeiFt63E0tXWKHYeIrBwLEBH1SKfJjDXbT2NCqBf8NdZ5tey7EoJR39yBf247JXYUIrJyLEBE1CNbDpajpK4FN8db396f83zdVbhhjD/e2XEWRTXNYschIivGAkREl2U2C/jntlOID/JAiLer2HEu6ca4ALir5Pi/rUfFjkJEVowFiIgua+epapyoaMRNcQFiR7kslZMMcyYEI/1IBXadrhY7DhFZKRYgIrqstT+dRaiPK2J04q/70xNJ4d6I0rrh+S+PotNkFjsOEVkhFiAiuqSz1U3Yll+FaSO1oq/63FMSiQRzE0OQr2/Ax/uKxY5DRFaIBYiILumDrAKoVXIkhfuIHaVXwn3dcG2UL/7+TT7qmnitJyLqigWIiC6qqa0Tm7JLMCXGDwq57b1dzB4fhA6TGX//Nl/sKERkZWzvHY2IhsxneaVobu/EdSO0YkfpE4+frxO2fk8RrxNGRF2wABHRRW3YW4Qrgz3h7aYUO0qfXTdSi2Geznh2M68TRkS/YAEiom4dKTPgcJkR10T7ih2lX+RSKRZMCsX+4nps4IRoIvoZCxARdes/+4rh6eKEK4I8xY7SbyP81bg22hdpXx9DVUOb2HGIyAqwABHRBVo7TPgsrxRXRfpCJrWNU98v53cTgiEB8NctXCGaiFiAiKgb3xzRw9jaiWujbPvw16+5q5xwV8JwfH6gDN8frxA7DhGJjAWIiC6wKbsEI/zd4e9hnVd976urIn0QH+SBpz85BENzh9hxiEhELEBE1EVlQyt2na7GpAjbWviwJyQSCe65KgzNbZ1Y9uURseMQkYhYgIioiy0HyiGVSJAQ6i12lEHh5arA3MQQfJZXivTD5WLHISKRsAARUReb80oRH+QBN6Vc7CiD5qpIH0wI9cKf/3sQxbXNYschIhGwABGRxdnqJhwsNdjl4a9fk0gk+NNVYVA5yfDQx7no4BXjiRwOCxARWWzOK4WzkwxXBtv+2j+X46qU4+HfROBwqREvf31c7DhENMRYgIgIACAIAjbvL8X4EE+bvPBpX0T4ueP3CcF498ez+A9XiSZyKI7xLkdEl3WsvAGFNc1IDLfvw1//K2WUDlNj/LD4s0PYc6ZG7DhENERYgIgIAJB+uByuShlGB6jFjjKkJBIJ5k8KQbTOHX/6dw7y9Q1iRyKiIcACREQAgK8OlePKYE/IZY73tiCXSvFYchQ8XZzwu3/txtnqJrEjEdEgc7x3OiK6wMmKBpyuarLbtX96wk0px6LrR8DZSYbfvbMbBSxBRHaNBYiI8PVhPZydZIgN1IgdRVRqZycsun4EJBIJbluzC0fKDGJHIqJBwgJERNh6qBxXBHs4zNlfl+LlqsDSG0bCw9kJs9/ajV2nqsWORESDgO92RA6usKYJx/UNmBDqJXYUq6F2dsJfZoxEmK8r7n53L/6dVSB2JCIaYCxARA4u42gFnGQSxA3zEDuKVXFWyPDnlBhcN0qL5z4/gkWfHkRrh0nsWEQ0QOz3Yj9E1CMZRyswOkADlZNM7ChWRyaVYF5iCIK9XLD2p7M4UGzAmt+PRbC3i9jRiKifuAeIyIEZmjuQXVCHK4fb/6Uv+mNKtB+ev2k0apvaMfP1nfj2iF7sSETUTyxARA5s24lKmATBIa791V+hPq546ZbRiPE/t2DiS18d5UVUiWwYCxCRA/vuaAXCfFzh5aoQO4pNcFHI8XhyFO6eOBzv/VSAue/uRV1Tu9ixiKgPWICIHFSHyYxt+VW4IthD7Cg2RSKRYEasPxZfH4PDZQbctOpHnKrk5TOIbA0LEJGD2ne2Fg1tnRg7nKe/98XIAA3+evNoAMDta7JwqISLJhLZEhYgIgf1/fFKeLkqEMIzmvrMT63CkhtGwcdNiTvfzsLes7ViRyKiHrKKArR69WqEhIRApVIhISEBe/fuveT4TZs2ISYmBiqVCrGxsdi6dWuX78+fPx8SiaTLbfr06YP5FIhszrYTVYgbpoFEIhE7ik1zU8mx+PoRCPFxxYK1e3n5DCIbIXoB2rhxI1JTU7F06VLk5uYiLi4OKSkpqKys7Hb8rl27MGfOHCxcuBB5eXmYNWsWZs2ahcOHD3cZN336dJSXl1tuH3/88VA8HSKbUFrfglOVjVz8cIA4K2R4clo0tGoV5q/dh5K6ZrEjEdFliF6AXn31Vdxzzz1YsGABRo4ciTVr1sDFxQXvvfdet+Nfe+01TJ8+HU899RRGjBiBF198EVdeeSVWrVrVZZxSqYROp7PcPD15mi/ReTtOVEEqAUY7+MVPB5LKSYanUqIhATB/7T40tXWKHYmILkHUAtTe3o6cnBwkJydbtkmlUiQnJyMrK6vb+2RlZXUZDwApKSkXjN+2bRv8/PwQHR2N+++/HzU1NRfN0dbWBqPR2OVGZM+25VciUusOVyUXgx9IHi4K/Hl6DErqmrHsiyNixyGiSxC1AFVXV8NkMkGr1XbZrtVqodd3v9KqXq+/7Pjp06fjgw8+QGZmJpYvX47t27fj+uuvh8nU/XV80tLSoNFoLLegoKB+PjMi69VhMuPHk9UYw70/gyLQwxnzk0KwKacEn+8vFTsOEV2EXf7z784777T8f2xsLMaMGYPw8HBs27YNU6dOvWD8okWLkJqaavnaaDSyBJHdyi2sQ1O7CXFBHmJHsVtXR/riYIkBiz89hLHDPTHMk2faEVkbUfcA+fj4QCaToaKiosv2iooK6HS6bu+j0+l6NR4AwsLC4OPjg1OnTnX7faVSCbVa3eVGZK+2n6iCxtkJoT6uYkexWxKJBAsnh0LpJMMLXx4VOw4RdUPUAqRQKDB27FhkZmZatpnNZmRmZiIxMbHb+yQmJnYZDwAZGRkXHQ8AJSUlqKmpgb+//8AEJ7Jh209UYXSgBlKe/j6oXBRy3JUQjG+PVmD7iSqx4xDR/xD9LLDU1FS88847eP/993Hs2DHcf//9aGpqwoIFCwAAc+fOxaJFiyzjH330UaSnp2PFihU4fvw4li1bhuzsbDz00EMAgMbGRjz11FPYvXs3CgoKkJmZiZtvvhkRERFISUkR5TkSWYvapnYcLTMilvN/hkRimDdG+qux9IvDaO/khVOJrInoBWj27Nn4+9//jiVLliA+Ph779+9Henq6ZaJzUVERysvLLeOTkpKwfv16vP3224iLi8N///tfbN68GaNHn1uSXiaT4eDBg7jpppsQFRWFhQsXYuzYsdi5cyeUSqUoz5HIWmSdroEAYHQAD/MOBYlEgnlJISiqacb7uwrEjkNEvyIRBEEQO4S1MRqN0Gg0MBgMnA9EdmXRp4ew40QV/n5HnNhRHMq/dp5BXnE9fnr6N3BWyMSOQ2S3evP5LfoeICIaOj+erMIo7v0ZcjfHB8DQ0oGP9hSKHYWIfsYCROQgimubUVzXwvk/IvB1V+GqCB+s2X4arR3dr0dGREOLBYjIQfx0qhpSCTCSe4BEMeuKQNQ2tWPD3iKxoxARWICIHMaPp6oR7ucGF4Vdrn9q9bRqFSZF+OCtHWdgMnPqJZHYWICIHIDZLOCnU9UYHcDDX2KaPkqHckMrvj9eKXYUIofHAkTkAPIrGlDX3IGR/jz8JaYwXzdE+Lnh31kFYkchcngsQEQOYPeZGshlEkRp3cWO4vCSR2ix42Q1CqqbxI5C5NBYgIgcwO4zNYjyc4dCzl95sSWGecNdKecp8UQi47shkZ0zmwXsPlOLEf7c+2MNFHIpron2xX+yS3hKPJGIWICI7NxxfQMMLZz/Y02mRPvB0NKBHzgZmkg0LEBEdm73mRo4ySSI8OMeIGsR4OGMMF9XbN5fKnYUIofFAkRk53afqUGUlvN/rM2kcB98f7wShuYOsaMQOSS+IxLZsXPzf2oQo+PhL2szMcwbnSYB6UfKxY5C5JBYgIjs2DG9EcbWTl7+wgp5uSowOlCDz/LKxI5C5JBYgIjs2N6ztZDLJIjwdRM7CnUjKdwbe87UQG9oFTsKkcNhASKyY3vO1iLS143zf6zUhFAvyGUSbD3Ew2BEQ43vikR2ShAE7D1bi2jO/7FaLgo5RgWo8e3RCrGjEDkcFiAiO3Wmugm1Te2I0fH0d2s2drgX9p2tRX1zu9hRiBwKCxCRndp7thZSCXj9Lyt3ZbAnTILAK8QTDTEWICI7te9sLUJ8XOGskIkdhS7By1WBCD83HgYjGmIsQER2avfZGsRw749NGDvcE9vzq3htMKIhxAJEZIdK61tQVt/KBRBtxLjhnmjpMGHX6WqxoxA5DBYgIju072wtACCaE6BtQqCHM/w1KmQc5TwgoqHCAkRkh/YW1GKYpzPUzk5iR6EekEgkiBvmgW35lRAEQew4RA6BBYjIDmUX1PLsLxszZpgG5YZWnKluEjsKkUNgASKyM4bmDpyoaGQBsjEj/NVwkkmw40SV2FGIHAILEJGdyS2qAwAugGhjVE4yRGvdsZ0FiGhIsAAR2Znswlp4uDjBz10pdhTqpdhhHth9pgZtnTwdnmiwsQAR2Zl9Z+sQ5ecOiUQidhTqpTHDNGjtMCOnoE7sKER2jwWIyI60d5pxoKSep7/bqGAvF3i4OGH7SR4GIxpsLEBEduRImQFtnWZOgLZRUokEsQEaToQmGgIsQER2JLugDkq5FCE+LmJHoT4aFajB8fIGXh2eaJCxABHZkezCWkT4uUEu5a+2rRrpr4YAYPeZWrGjENk1vksS2QlBEJBdWIdIPzexo1A/+LoroVUrsftMjdhRiOwaCxCRnSiubUFNYzvn/9iBGJ0aWSxARIOKBYjITuQUnTtkEunHAmTrRvqrka9vQF0T5wERDRYWICI7kVNYh2GeznBTycWOQv00MkANANhzlnuBiAYLCxCRncguqEOEL+f/2AMft/PzgDgRmmiwsAAR2YHGtk6cqGjg/B87MkKnxq7T1WLHILJbLEBEduBAcT3MAliA7MjIADVOVDSilvOAiAYFCxCRHcgprIObUg5/D5XYUWiAxOjOzQPKLuBhMKLBwAJEZAdyfl7/R8oLoNoNHzcFvN0UyC7khVGJBgMLEJGNM5sF5BbWIYILINoViUSCKD937OMeIKJBwQJEZONOVTWioa2T83/sULTOHYdKDGjtMIkdhcjusAAR2bi8ojpIJUA4T4G3O1Fad3SaBRwsMYgdhcjusAAR2bjcwnoEebnAWSETOwoNsGAvFzg7yXgYjGgQsAAR2bicQi6AaK9kUgki/dx4JhjRILCKArR69WqEhIRApVIhISEBe/fuveT4TZs2ISYmBiqVCrGxsdi6detFx953332QSCRYuXLlAKcmEp+hpQOnqhoRqWUBsldROndkF9bBbBbEjkJkV0QvQBs3bkRqaiqWLl2K3NxcxMXFISUlBZWVld2O37VrF+bMmYOFCxciLy8Ps2bNwqxZs3D48OELxn722WfYvXs3AgICBvtpEIniQHE9AF4A1Z5Fa93R0NqJU1WNYkchsiuiF6BXX30V99xzDxYsWICRI0dizZo1cHFxwXvvvdft+Ndeew3Tp0/HU089hREjRuDFF1/ElVdeiVWrVnUZV1paiocffhgfffQRnJychuKpEA253KKfF0DUcAFEexXh5wapBJwHRDTARC1A7e3tyMnJQXJysmWbVCpFcnIysrKyur1PVlZWl/EAkJKS0mW82WzG3XffjaeeegqjRo26bI62tjYYjcYuNyJbcH79HwkXQLRbKicZhnu7Yn9RvdhRiOyKqAWouroaJpMJWq22y3atVgu9Xt/tffR6/WXHL1++HHK5HI888kiPcqSlpUGj0VhuQUFBvXwmREPPbBaQV1SPSC6AaPfCfd2QW8QVoYkGkuiHwAZaTk4OXnvtNaxbt67H/ypetGgRDAaD5VZcXDzIKYn670z1uQUQI7kAot2L9HPD6aomGFo6xI5CZDdELUA+Pj6QyWSoqKjosr2iogI6na7b++h0ukuO37lzJyorKxEcHAy5XA65XI7CwkI88cQTCAkJ6fYxlUol1Gp1lxuRtcstrIcEQLivq9hRaJCd38t3ftI7EfWfqAVIoVBg7NixyMzMtGwzm83IzMxEYmJit/dJTEzsMh4AMjIyLOPvvvtuHDx4EPv377fcAgIC8NRTT+Gbb74ZvCdDNMTyiusQ5OUCF4Vc7Cg0yHQaFdxVcuRxHhDRgBH9nTM1NRXz5s3DuHHjMGHCBKxcuRJNTU1YsGABAGDu3LkIDAxEWloaAODRRx/FNddcgxUrVmDmzJnYsGEDsrOz8fbbbwMAvL294e3t3eVnODk5QafTITo6emifHNEgyims4+UvHIREIkEE5wERDSjRC9Ds2bNRVVWFJUuWQK/XIz4+Hunp6ZaJzkVFRZBKf9lRlZSUhPXr1+PZZ5/F4sWLERkZic2bN2P06NFiPQWiIdfY1omTFY24NtpP7Cg0RCL83PDNET0EQeBZf0QDQCIIApcX/R9GoxEajQYGg4Hzgcgq/XSqGnf9aw/+dvsYDPN0ETsODYGDJfVI+/o4vn/iGoRxzx9Rt3rz+W13Z4EROYK8ojq4KmQI8HAWOwoNkQg/N0gAzgMiGiAsQEQ2KLeoHuF+bpDyUIjDcFHIMczTmfOAiAYICxCRjREEAXlF51aAJscS5uuG/TwVnmhAsAAR2ZjCmmbUNXdwBWgHFO7rinx9A1o7TGJHIbJ5LEBENiav+NwhkAhfrgDtaMJ93dBpFnC0nNcrJOovFiAiG5NXVI8ADxXcVKKvYkFDLNjLBXKZhCtCEw0AFiAiG5NbWIcIngbtkOQyKUK9XXGwxCB2FCKbxwJEZENa2k04rm/gBVAdWKiPK/J4JhhRv7EAEdmQQ6UGdJoFToB2YBF+biioaYahmVeGJ+oPFiAiG5JXVAeVk5SrPzuw89d/O1haL24QIhvHAkRkQ/KKzl0AVSblAoiOSqdRwUUh4zwgon5iASKyEYIgIKeongsgOjipRIJwLohI1G8sQEQ2oszQiqqGNhYgQqiPKwsQUT+xABHZiPNn/kT68QwwRxfu64aqhjZUGFvFjkJks1iAiGxEXlE9tGoVNM5OYkchkYX5ugIA5wER9QMLEJGNyCmsQ8TPH3zk2LxdFfBwdsKhknqxoxDZLBYgIhvQ1mnCkTIDInj4iwBIJBKE+rjiAPcAEfVZnwrQmTNnBjoHEV3CkTIjOkwCIrWcAE3nhPq64mBJPQRBEDsKkU3qUwGKiIjAlClT8OGHH6K1lZPwiAZbbmEdFDIphntxAUQ6J8zHDXXNHSitbxE7CpFN6lMBys3NxZgxY5CamgqdTod7770Xe/fuHehsRPSz3KI6hPm6Qi7jUWs65/xE6EM8DEbUJ316N42Pj8drr72GsrIyvPfeeygvL8fkyZMxevRovPrqq6iqqhronEQOLbeQCyBSV54uCni7KjgPiKiP+vXPSblcjltvvRWbNm3C8uXLcerUKTz55JMICgrC3LlzUV5ePlA5iRxWuaEFemMrojgBmv5HqM+5eUBE1Hv9KkDZ2dl44IEH4O/vj1dffRVPPvkkTp8+jYyMDJSVleHmm28eqJxEDiuvqB4AEMEJ0PQ/Qn1ccajEwInQRH0g78udXn31Vaxduxb5+fmYMWMGPvjgA8yYMQNS6bk+FRoainXr1iEkJGQgsxI5pNzCOvi6K+HpohA7ClmZMF83NLR1oqCmGaE+XCOKqDf6VIDefPNN/OEPf8D8+fPh7+/f7Rg/Pz+8++67/QpHREBOUR3n/1C3LBOhSw0sQES91KdDYBkZGXj66acvKD+CIKCoqAgAoFAoMG/evP4nJHJgbZ0mHCk1IooFiLqhVjnB113JFaGJ+qBPBSg8PBzV1dUXbK+trUVoaGi/QxHROUfLjGg3mbkCNF3UuYnQPBOMqLf6VIAuNuGusbERKpWqX4GI6Bc5Py+AGOLNBRCpe6E+rjhcaoDZzInQRL3RqzlAqampAM5dh2bJkiVwcfnlTdlkMmHPnj2Ij48f0IBEjiynsA7hflwAkS4uzMcVTe0mFNQ0IcyXh0qJeqpXBSgvLw/AuT1Ahw4dgkLxy1kpCoUCcXFxePLJJwc2IZGDEgQBOYV1mBjmLXYUsmLnJz8fKjWwABH1Qq8K0A8//AAAWLBgAV577TWo1epBCUVEQJmhFZUNbYjkBGi6BHeVE7RqJQ6WGHBzfKDYcYhsRp9Og1+7du1A5yCi/5FTWAcAiNRyAjRdWoi3K68JRtRLPS5At956K9atWwe1Wo1bb731kmM//fTTfgcjcnS5hXXw16igcXYSOwpZuTAfV3x+oAxmswCpVCJ2HCKb0OMCpNFoIJFILP9PRIMrp7AOEZzTQT0Q5uuG5nYTzlQ3cdFMoh7qcQH69WEvHgIjGlzN7Z04WmbEvKQQsaOQDQj5eSL0wZJ6FiCiHurTubUtLS1obm62fF1YWIiVK1fi22+/HbBgRI7sYIkBJkFAFC+ASj3gppRDp1bhUCnnARH1VJ8K0M0334wPPvgAAFBfX48JEyZgxYoVuPnmm/Hmm28OaEAiR5RTWAdnJxmCPLkAIvXM+SvDE1HP9KkA5ebm4qqrrgIA/Pe//4VOp0NhYSE++OADvP766wMakMgR5RTUItLPjRNaqcfCfF1xuMwAE1eEJuqRPhWg5uZmuLufOzX322+/xa233gqpVIqJEyeisLBwQAMSORqzWUB2YR1Pf6deCfVxRWuHGaerGsWOQmQT+lSAIiIisHnzZhQXF+Obb77BtGnTAACVlZVcHJGon05VNcLY2okYHQsQ9VyoZSI0D4MR9USfCtCSJUvw5JNPIiQkBAkJCUhMTARwbm/QFVdcMaABiRxNdkEdpBLwbB7qFReFHAEeKhzmRGiiHunTStC33347Jk+ejPLycsTFxVm2T506FbfccsuAhSNyRNmFtQjxdoXKSSZ2FLIxod6uOFBSL3YMIpvQpwIEADqdDjqdrsu2CRMm9DsQkaPLLqjDCH8eSqbeC/Vxw6acYnSazJDL+rSDn8hh9KkANTU14eWXX0ZmZiYqKythNpu7fP/MmTMDEo7I0VQ1tKGothmz4gPEjkI2KMzXFW2dZpysbGSJJrqMPhWgP/7xj9i+fTvuvvtu+Pv7Wy6RQUT9k1NYCwCI4hlg1Ach3q6QADhUamABIrqMPhWgr7/+Gl999RUmTZo00HmIHFp2QR183ZTwdlOKHYVskLNChgAPZxwqMeC344LEjkNk1fp0kNjT0xNeXl4DnYXI4e0tqOXlL6hfQn1ccaC4XuwYRFavTwXoxRdfxJIlS7pcD6w/Vq9ejZCQEKhUKiQkJGDv3r2XHL9p0ybExMRApVIhNjYWW7du7fL9ZcuWISYmBq6urvD09ERycjL27NkzIFmJBsv5C6BG63jogvou3NcVx/RGtHeaLz+YyIH1qQCtWLEC33zzDbRaLWJjY3HllVd2ufXGxo0bkZqaiqVLlyI3NxdxcXFISUlBZWVlt+N37dqFOXPmYOHChcjLy8OsWbMwa9YsHD582DImKioKq1atwqFDh/Djjz8iJCQE06ZNQ1VVVV+eLtGQyCuqR6dZ4AKI1C9hvm7oMAnI1zeIHYXIqkkEQej1hWOef/75S35/6dKlPX6shIQEjB8/HqtWrQIAmM1mBAUF4eGHH8YzzzxzwfjZs2ejqakJW7ZssWybOHEi4uPjsWbNmm5/htFohEajwXfffYepU6deNtP58QaDgStb05D5R8YJvPfTWaz5/VhIeWIB9VFbpwl/WLcPL84ajbsShosdh2hI9ebzu0+ToHtTcC6lvb0dOTk5WLRokWWbVCpFcnIysrKyur1PVlYWUlNTu2xLSUnB5s2bL/oz3n77bWg0mi6LNv5aW1sb2traLF8bjcZePhOi/tt7thZRWneWH+oXpVyGYC8XHCw24K4EsdMQWa8+r5RVX1+Pf/3rX1i0aBFqa8+dupubm4vS0tIeP0Z1dTVMJhO0Wm2X7VqtFnq9vtv76PX6Ho3fsmUL3NzcoFKp8I9//AMZGRnw8fHp9jHT0tKg0Wgst6Agnj1BQ6vDZEZecR0Pf9GACPVxw36uCE10SX0qQAcPHkRUVBSWL1+Ov//976ivrwcAfPrpp1325ohpypQp2L9/P3bt2oXp06fjt7/97UXnFS1atAgGg8FyKy4uHuK05OgOlxrQ2mFmAaIBEebrilMVjWhpN4kdhchq9akApaamYv78+Th58iRUKpVl+4wZM7Bjx44eP46Pjw9kMhkqKiq6bK+oqLjgMhvn6XS6Ho13dXVFREQEJk6ciHfffRdyuRzvvvtut4+pVCqhVqu73IiG0r6CWijlUoT8fEVvov4I93WDSRBwtJwXRiW6mD4VoH379uHee++9YHtgYOBFD111R6FQYOzYscjMzLRsM5vNyMzMtFxh/n8lJiZ2GQ8AGRkZFx3/68f99TwfImuy92wtIv3cIJfy+k3Uf0FeznCSSXCgmAWI6GL6NAlaqVR2O1H4xIkT8PX17dVjpaamYt68eRg3bhwmTJiAlStXoqmpCQsWLAAAzJ07F4GBgUhLSwMAPProo7jmmmuwYsUKzJw5Exs2bEB2djbefvttAOeuU/bSSy/hpptugr+/P6qrq7F69WqUlpbijjvu6MvTJRpUZrOAfQV1SB7hJ3YUshNyqRQh3q44VMoCRHQxfSpAN910E1544QX85z//AQBIJBIUFRXh6aefxm233darx5o9ezaqqqqwZMkS6PV6xMfHIz093TLRuaioCNJf/as4KSkJ69evx7PPPovFixcjMjISmzdvxujRowEAMpkMx48fx/vvv4/q6mp4e3tj/Pjx2LlzJ0aNGtWXp0s0qE5UNsDQ0sFrN9GACvVxxX6uCE10UX1aB8hgMOD222/Hvn370NjYiICAAOj1eiQmJmLr1q1wdbXteQxcB4iG0vu7CvDXr47iX3PHQyHnITAaGNtPVGHN9tM4uGwa1ConseMQDYlBXwdIo9EgIyMDP/30Ew4cOIDGxkZceeWVSE5O7lNgIkeWdaYG4b5uLD80oCJ8z11T7lCJAZMiul8ChMiR9boAmc1mrFu3Dp9++ikKCgogkUgQGhoKnU4HQRAg4SJuRD0mCAL2nq3F1ZG9mztHdDn+Hiq4KGTYX1zPAkTUjV79k1MQBNx000344x//iNLSUsTGxmLUqFEoLCzE/PnzccsttwxWTiK7dKqyEbVN7Rjhz/V/aGBJJRKE+7pxHhDRRfRqD9C6deuwY8cOZGZmYsqUKV2+9/3332PWrFn44IMPMHfu3AENSWSvdp+pgUwqQZSWBYgGXpivK346Vc2980Td6NUeoI8//hiLFy++oPwAwG9+8xs888wz+OijjwYsHJG9232mFuG+rlA5ycSOQnYo3NcN1Y3t0BtbxY5CZHV6VYAOHjyI6dOnX/T7119/PQ4cONDvUESOQBAE7D5bgxgdzzSkwRHhd24i9P6ienGDEFmhXhWg2traCy5E+mtarRZ1dXX9DkXkCM5UN6GmsZ3r/9Cg8XRRwNtNwQujEnWjVwXIZDJBLr/4tCGZTIbOzs5+hyJyBLtOn5v/E835PzSIwn3cuAeIqBu9mgQtCALmz58PpVLZ7fd5rS2intt1qhoRvm5wVnD+Dw2ecF9XfH6gDCazAJmUE6GJzutVAZo3b95lx/AMMKLLM5sFZJ2pwZRoXv+LBleEnxua2004WdnA+WZEv9KrArR27drBykHkUI7rG1Df3IFRAfxAosEV6uMGqeTcRGgWIKJfcO19IhHsOl0NJ5kEkX6c/0ODy1khQ5CXC/I4D4ioCxYgIhHsOl2DaK07r/9FQyLC1w05hTxDl+jX+O5LNMQ6TWbsOVODkQEasaOQg4jUuuNUVSMMLR1iRyGyGixAREPsUKkBTe0mzv+hIRP584KIB3hdMCILFiCiIbbrdA2cnWQI83UVOwo5CH+NCu4qOXKLeBiM6DwWIKIhtvNkNUb4u0Mu5a8fDQ2JRIIIXzdOhCb6Fb4DEw2hlnYTcgprERvI+T80tCL83JBbVAezWRA7CpFVYAEiGkJ7ztagwyQgNtBD7CjkYCK17mho7cSZ6kaxoxBZBRYgoiH048lqeLsqEOChEjsKOZhwX1dIAOTyMBgRABYgoiG142QVRgdqIJHwmkw0tFwUcgR5uSCngBOhiQAWIKIhU9nQihMVjZz/Q6KJ0rpjX0Gt2DGIrAILENEQ+elUNQBw/R8STYzOHWeqm1DT2CZ2FCLRsQARDZGdJ6sx3NsFHi4KsaOQg4rWnbv2HC+LQcQCRDQkzGYB209U8fAXicrHTQkfNwWyWYCIWICIhsIxvRE1je2IG+YhdhRycFFad+w9y3lARCxARENgW34VVE5SxPx8CIJILNE6dxwuNaCl3SR2FCJRsQARDYHt+VUYFaCBXMZfORJXtNYdnWYBB0rqxY5CJCq+GxMNMmNrB3IK6xA3jPN/SHxBni5wVciQzdPhycGxABENsl2nqmESBM7/IasglUoQqXXjPCByeCxARINs+4kqBHo4w0/Ny1+QdYjWqZFdWIdOk1nsKESiYQEiGkSCIGBbfhViefiLrMhIfzWa2004XGYUOwqRaFiAiAZRfkUDyg2tiOfhL7IiYb6uUMql2H2mRuwoRKJhASIaRJnHKqFykmIkL39BVkQulSJa546s0yxA5LhYgIgGUeaxCsQGauDE09/Jyoz0V2NfQS3nAZHD4rsy0SCpbWpHXlE9rgj2FDsK0QU4D4gcHQsQ0SDZll8JAcAVQR5iRyG6QKivK1ROnAdEjosFiGiQZB6rRISvK6/+TlZJLpUiWst5QOS4WICIBkGHyYztJ6oQz8NfZMVG/DwPqIPzgMgBsQARDYJ9BbVobOvk4S+yaqMCNGhuN+EgrwtGDogFiGgQZBytgLerAqE+rmJHIbqoMB9XuCnl2HGiWuwoREOOBYhogAmCgG+O6HHlcE9IJBKx4xBdlFQqwcgANXaerBI7CtGQYwEiGmBHy40oq2/FuOGc/0PWb0ygBgeKDTC2dogdhWhIsQARDbCMoxVwVcgw0p+rP5P1iw3UwCQIPBuMHA4LENEA++aIHnFBHpBz9WeyAX5qFXRqFX48yXlA5Fj4Dk00gIprm3GsvAHjhnuJHYWox0YHarCD84DIwbAAEQ2g745VQC6VIC5II3YUoh4bE6hBYU0zimubxY5CNGSsogCtXr0aISEhUKlUSEhIwN69ey85ftOmTYiJiYFKpUJsbCy2bt1q+V5HRweefvppxMbGwtXVFQEBAZg7dy7KysoG+2kQ4etDeowO1MBFIRc7ClGPjQxQQyoBtp/gXiByHKIXoI0bNyI1NRVLly5Fbm4u4uLikJKSgsrKym7H79q1C3PmzMHChQuRl5eHWbNmYdasWTh8+DAAoLm5Gbm5uXjuueeQm5uLTz/9FPn5+bjpppuG8mmRA6psaMW+glpMCOXhL7Itrko5YnRqfH+8+/ddInskEQRBEDNAQkICxo8fj1WrVgEAzGYzgoKC8PDDD+OZZ565YPzs2bPR1NSELVu2WLZNnDgR8fHxWLNmTbc/Y9++fZgwYQIKCwsRHBx82UxGoxEajQYGgwFqNc/koZ75cHchlnx+GGt+PxbuKiex4xD1ypcHyvBpbgn2L50GlZNM7DhEfdKbz29R9wC1t7cjJycHycnJlm1SqRTJycnIysrq9j5ZWVldxgNASkrKRccDgMFggEQigYeHR7ffb2trg9Fo7HIj6q2th8oxKkDD8kM26cpgT7R2mnk6PDkMUQtQdXU1TCYTtFptl+1arRZ6vb7b++j1+l6Nb21txdNPP405c+ZctA2mpaVBo9FYbkFBQX14NuTIapvasedMLcaH8PAX2aYADxW0aiUyj1eIHYVoSIg+B2gwdXR04Le//S0EQcCbb7550XGLFi2CwWCw3IqLi4cwJdmDjKN6mAUB40O4+jPZJolEgiuCPJF5rBIiz4wgGhKiFiAfHx/IZDJUVHT9F0dFRQV0Ol2399HpdD0af778FBYWIiMj45LHApVKJdRqdZcbUW9sPaRHjL87PFwUYkch6rMrgj1QbmjFcX2D2FGIBp2oBUihUGDs2LHIzMy0bDObzcjMzERiYmK390lMTOwyHgAyMjK6jD9ffk6ePInvvvsO3t7eg/MEiADUNbXjp1PVSAjl3zOybSP81VA5SXk2GDkE0Q+Bpaam4p133sH777+PY8eO4f7770dTUxMWLFgAAJg7dy4WLVpkGf/oo48iPT0dK1aswPHjx7Fs2TJkZ2fjoYceAnCu/Nx+++3Izs7GRx99BJPJBL1eD71ej/b2dlGeI9m39CPnDn8l8PR3snFOMinGDPPAN4e7n1NJZE9EX61t9uzZqKqqwpIlS6DX6xEfH4/09HTLROeioiJIpb/0tKSkJKxfvx7PPvssFi9ejMjISGzevBmjR48GAJSWluKLL74AAMTHx3f5WT/88AOuvfbaIXle5Di+PFCGkf5qHv4iuzAhxAurfjiFkrpmDPN0ETsO0aARfR0ga8R1gKinqhrakPB/3+EPk0MxNUZ7+TsQWbnm9k7c92EOnp4egz9eFSZ2HKJesZl1gIhs3deHyyGRSDCBp7+TnXBRyDEm0ANbD5WLHYVoULEAEfXDlwfKEBvIxQ/JvowP9UJuUT0qjK1iRyEaNCxARH1UVt+CfQV1mBjGs7/Ivowd7gmZVIJvjnAyNNkvFiCiPvriQBkUMikXPyS746aUY3SAGl8d5GEwsl8sQER99FleKa4c7gEXhegnUxINuIlh3th7thblhhaxoxANChYgoj44rjciX9+ASRE+YkchGhQTQr3gJJPii/1lYkchGhQsQER9sDmvDO5KOeKHeYgdhWhQuCjkGDvcE5/mloodhWhQsAAR9ZLZLODz/aVICPOCXMZfIbJfkyN8kF/RgKNlRrGjEA04vnsT9dK+glqUG1p5+Ivs3pggDdQqOTbv514gsj8sQES99EluCbRqJaK07mJHIRpUcqkUE8O8sTmvFCYzLxpA9oUFiKgXWtpN+OpgOSZH+EAqkYgdh2jQXR3li8qGNuw4USV2FKIBxQJE1AvfHNGjqd2EqyJ9xY5CNCTCfFwR4u2Cj/YUih2FaECxABH1wn9zSjDC3x1atUrsKERDQiKRYOoILb4/Xsk1gciusAAR9VC5oQU/narGVRHc+0OOZVK4D5RyGTbsLRY7CtGAYQEi6qHP8kqhkEuREMYrv5NjcVbIkBTujQ37itBpMosdh2hAsAAR9YAgCNi4rxjjQ7x46QtySFNHaFFhbEPm8UqxoxANCBYgoh7Ye7YWhTXNmBLNw1/kmEJ9XBGtdce7O8+KHYVoQLAAEfXAxuxi6NQqjPBXix2FSDQzYv2xt6AWB4rrxY5C1G8sQESXYWztwNaD5bgmyhcSrv1DDmzccE/o1Cq8s/OM2FGI+o0FiOgyvjxQhnaTGVdH8fAXOTapVILrR+uw9VA5imubxY5D1C8sQESXsXFvMeKDPODlqhA7CpHoro7yhatCjnd/5Fwgsm0sQESXcLjUgIOlBkyJ9hM7CpFVUDnJMG2UDh/vLUKFsVXsOER9xgJEdAnr9xbB21WBK4I9xY5CZDVmxOogl0nw5rbTYkch6jMWIKKLaGzrxOa8UlwT7QuZlJOfic5zUcgxY7Q/1u8p4uUxyGaxABFdxBf7y9DaYcJvePiL6ALTR+uglEvxzx+4F4hsEwsQUTcEQcCHuwsRH+QBbzel2HGIrI6LQo6ZY/zx8d4inK1uEjsOUa+xABF140CJAUfLjZgaoxU7CpHVun60PzxcnJC29ZjYUYh6jQWIqBsf7CqAn7sS8UEeYkchsloKuRR3jg/Gt0crkHW6Ruw4RL3CAkT0P2oa27DlYDmmjtBCysnPRJeUFO6NSD83vLjlKExmQew4RD3GAkT0PzZmF0OAgGt54VOiy5JIJPj9xOE4Wm7E+r1FYsch6jEWIKJfMZnPTX5OCveBWuUkdhwimxCldceUaD8s//o4F0ckm8ECRPQrmccqUFbfiutGcvIzUW/8LiEYMqkEy744InYUoh5hASL6lfd+KkCk1g3hvm5iRyGyKW5KOe6eOBxfH9Yj42iF2HGILosFiOhnR8uM2H2mBteP0okdhcgmJYV748pgDzzz6UHUNLaJHYfokliAiH629qez8HZTYHyol9hRiGySRCLBPVeFob3TjMWfHYIg8Kwwsl4sQEQAqhvbsHl/KaaN0EIu5a8FUV95uCjwx8lh+OZIBf6bUyJ2HKKL4js9EYAPdxdCKpHgN1z5majfJoR64eooHyz5/AhOVTaKHYeoWyxA5PBaO0z4IKsQV0X6wk0lFzsOkV1YkBQKL1cFHvgoBy3tJrHjEF2ABYgc3n9zSlDf3I6Zsf5iRyGyGyonGR6dGonCmmYs/eKw2HGILsACRA7NZBbwzs4zGB/iBZ1GJXYcIrsS5OWCBZNC8J/sEny0p1DsOERdsACRQ8s4qkdhTTNujAsQOwqRXbomyg/XjdRi6edHsK+gVuw4RBYsQOSwBEHAm9tPY6S/mgsfEg2iuYnDEal1w30f5qCsvkXsOEQAWIDIgWWdrsGBYgNujOPcH6LBJJdK8ejUKEglEvxh3T40tXWKHYmIBYgc1+vfn0KYjyvihnmIHYXI7mmcnfDUtGgU1jTj4Y/zYDJzkUQSFwsQOaScwlrsPlODm+MDIZFIxI5D5BCCvFzwyNRIbMuvxItbjoodhxwcCxA5pDe+P4Vhns4YF+IpdhQihxIf5IH5SSFYt6sA/9p5Ruw45MC46hs5nEMlBmzLr8KDUyIg5d4foiF33Ugdqhvb8dJXx+CvccbMMZyHR0NP9D1Aq1evRkhICFQqFRISErB3795Ljt+0aRNiYmKgUqkQGxuLrVu3dvn+p59+imnTpsHb2xsSiQT79+8fxPRki1Zk5CPAwxlJYd5iRyFyWLPHByEpwhuPb9yPPWdqxI5DDkjUArRx40akpqZi6dKlyM3NRVxcHFJSUlBZWdnt+F27dmHOnDlYuHAh8vLyMGvWLMyaNQuHD/+yymhTUxMmT56M5cuXD9XTIBuSU1iHbflVuO3KQEil3PtDJBapRIJ7rw5HlNYNf3w/G/n6BrEjkYORCIIg2lT8hIQEjB8/HqtWrQIAmM1mBAUF4eGHH8YzzzxzwfjZs2ejqakJW7ZssWybOHEi4uPjsWbNmi5jCwoKEBoairy8PMTHx/cql9FohEajgcFggFqt7v0TI6t11792o6imGS/fNoaHv4isQHN7J17cchQtHSZ89sAkBHg4ix2JbFhvPr9F2wPU3t6OnJwcJCcn/xJGKkVycjKysrK6vU9WVlaX8QCQkpJy0fE91dbWBqPR2OVG9mfPmRr8dKoGt40dxvJDZCVcFHL8eXoMBAG4+909qG9uFzsSOQjRClB1dTVMJhO0Wm2X7VqtFnq9vtv76PX6Xo3vqbS0NGg0GsstKCioX49H1kcQBKR9fRxhPq4YH+Ildhwi+hVPFwWenh6DyoY2LHw/G60dvHo8DT7RJ0Fbg0WLFsFgMFhuxcXFYkeiAZZ+WI/9xfWYMyGYe3+IrFCAhzP+nBKNw6UGPMKFEmkIiFaAfHx8IJPJUFFR0WV7RUUFdDpdt/fR6XS9Gt9TSqUSarW6y43sR4fJjOXpxxEf5IHRgRqx4xDRRUT4ueORqZH47lgFln1xBCJOUSUHIFoBUigUGDt2LDIzMy3bzGYzMjMzkZiY2O19EhMTu4wHgIyMjIuOJwKADXuLUFjTjDkTgsWOQkSXcWWwJxZODsO/dxfin9tOix2H7JioCyGmpqZi3rx5GDduHCZMmICVK1eiqakJCxYsAADMnTsXgYGBSEtLAwA8+uijuOaaa7BixQrMnDkTGzZsQHZ2Nt5++23LY9bW1qKoqAhlZWUAgPz8fADn9h71d08R2Z765nasyDiBa6J8EezlInYcIuqB38T4obapDX/7Jh8BHirccsUwsSORHRK1AM2ePRtVVVVYsmQJ9Ho94uPjkZ6ebpnoXFRUBKn0l51USUlJWL9+PZ599lksXrwYkZGR2Lx5M0aPHm0Z88UXX1gKFADceeedAIClS5di2bJlQ/PEyGr8I+ME2jrMmD2eE9uJbMltVw5DTWM7ntp0EFp3FZIifMSORHZG1HWArBXXAbIPx8qNmPn6TsyZEIwbxgSIHYeIeqnTbMbfv8nH6aomfPpAEqK07mJHIitnE+sAEQ0mQRCw7Isj8Nc4Y/ooHvokskVyqRSPTI2Et5sC89/bi8qGVrEjkR1hASK79N+cEuw5W4u5icMhl/GvOZGtclHI8dS0aLR2mrFwXTaa2zvFjkR2gp8MZHeqGtrw4ldHMTnCB2OGeYgdh4j6ydtNiadSonGysgGpG/fDzDWCaACwAJHdeeHLI4AA3J04XOwoRDRAQrxd8dCUSHxzpAKvfJMvdhyyAyxAZFcyjlbgy4Pl+P3E4VCrnMSOQ0QDaOxwT9yVMBxrtp/Gf7K5Yj/1j6inwRMNpOrGNjz9yUGMDfbEZJ4yS2SXZsTqUG5oweJPDyHE2xUTQnltP+ob7gEiuyAIAp755CBMZgF/vCoUEl7vi8guSSQSzJ8UgmidO/7072wU1TSLHYlsFAsQ2YWP9xbju2OV+OPkUHi4KMSOQ0SDSC6V4tGpkVA5ybDw/X1oaO0QOxLZIBYgsnmHSw1Y9sURTI3xw7gQ7g4ncgTuKic8OS0apfUteGzDfl49nnqNBYhsmqGlA/d/mINhns6YmxgidhwiGkKBHs54+DcR+CG/En//lmeGUe+wAJHNMpkFpG7cj9rmdjw6NRIKOf86Ezma+CBPzJkQjDe3ncbn+0vFjkM2hJ8YZLOWpx/HD/mVePDaCPipVWLHISKRzIz1x9WRPvjzfw/iQHG92HHIRrAAkU36eG8R3t5xBr+fOBxXBHuKHYeIRCSRSLBwchiCvVxwzwfZqDTymmF0eSxAZHMyjlbg2c8OI3mElhc6JSIAgEIuxePXRcFkFnDPB9lo7TCJHYmsHAsQ2ZSfTlXjwY9yMXa4J+YnhXC9HyKy8HRR4PHronCsvAGLPj0EQeCZYXRxLEBkM/YV1OKP72djhL87HvpNBGRSlh8i6irc1w1/ujoMn+WV4q0dZ8SOQ1aMl8Igm/DjyWr88YN9CPd1w+PXRcFJxu5ORN2bFOGD4rpmLP/6OCJ83ZA8Uit2JLJC/BQhq/ftET3+sG4fRujU+HNKDJRymdiRiMjK/XZcEMaFeOKRDXk4rjeKHYesEAsQWbW1P53Fvf/OwRXBHki9Lopr/RBRj0glEjxwbQT83JX4w7p9qGpoEzsSWRl+mpBV6jCZsfTzw3j+y6OYOcYfj0yNhJyHvYioF1ROMjwxLRot7SaeGUYX4CcKWZ3Khlbc9c4efLinCH+YFIq7EoZDyrO9iKgPfNyUeGJaNI6VG/HkpgMw85ph9DMWILIqWadrcMPrP+JkZQOemzkS13HyIhH1U7ivG+6/NhxfHSzH33jNMPoZzwIjq9BpMuO1zJNY9f0pjAxQ48EpEfB0UYgdi4jsREKoN+5KaMeb204j0MMZv584XOxIJDIWIBLdyYoGpP7nAI6UGfDbcUG4KS4AUq7xQ0QDbEasDlWNbVjy+WFo1SruYXZwPARGojGZBby94zRmvv4japra8fxNozDrikCWHyIaFBKJBHMnDse4EC88tD4Xe87UiB2JRMQCRKLI1zfgltU/IW3rcSSP8EPaLbGI8HMXOxYR2TmpVIKHpkQg0s8NC9/PxpEyg9iRSCQsQDSkWjtMWPFtPma+vhO1ze1YdtMo3J0YwvV9iGjIOMnOXTjVT63E3e/uxcmKBrEjkQj4qUNDZtfpakx/bQfe3HYaN8cH4P9uiUWUlnt9iGjouSjkeGZ6DNyUcsx5Zw/OVjeJHYmGGAsQDbrapnY88Z/9+N07e6CSy/DyrWNw+9ggXs+LiETlrnLC4hkjoJRLMfutLJyqbBQ7Eg0hfgLRoBEEAf/NKcFvVmxD+hE97rkqDM/dMBKBns5iRyMiAgBonJ3w7MxfShCvG+Y4WIBoUJytbsLv3tmDJzcdwKgANf5+exx+E+PHFZ2JyOp4uCjw7A0j4a6SY/Zbu5FTWCd2JBoCLEA0oNo7zVj1/Umk/GMHzlQ34pnpMXhoSiQ8uKghEVkxtcoJf5k5EgEaFX73zm58e0QvdiQaZCxANGByi+pwwxs78WrGCaSM0mL5bWMQF+Qhdiwioh5xU8rxzPUjEB/kgfs+zMG/dp6BIPDaYfaKK0FTvzW1deJv3+Tj/V0FCPN1xUu3xCLE21XsWEREvaaQS/HIbyLx8b4i/PWrYzhaZsT/3RoLlZNM7Gg0wFiAqF9+yK/EXz49hJqmdtyVMBzTR+sg40rORGTDpFIJ7koYjuHernhnxxkc0xvxxpwrEeHnJnY0GkA8BEZ9UtPYhkc35GHB2n3wcVPildvGYOYYf5YfIrIbkyN88PzNo2Bo6cQNb+zEx3uLeEjMjkgE/mlewGg0QqPRwGAwQK1Wix3HqgiCgM/ySvHClqMwmQT8fuJwXBXpAwnP7iIiO9XaYcIHWQX4Ib8KkyK8kXbLGAR7u4gdi7rRm89vFqBusAB1r6imGX/ZfAg7T1YjKdwbd08czrO7iMhhHCiux7s/nUVDSwf+dE047rsmDC4KziSxJixA/cQC1FWHyYx/7TyL1747AXeVHPMnheLKYE+xYxERDbnWDhM+yyvF14fL4eGiwCO/icAd44I4SdpKsAD1EwvQL/acqcGzmw/jdFUjrh/tj9vHDuMvOhE5vEpjK/6TU4ys0zXwcVPij1eF4rfjgrhXXGQsQP3EAgRUGFuRtvUYNu8vQ5TWDQsmhfLUdiKi/1Fe34LNB0qRdboGMokEN8YF4NYrhyEh1AtSnhQy5FiA+smRC1Brhwnv/ngWq74/BYVcitnjgnBNtC8vYUFEdAmGlg58f7wS209UosLYBn+NCimjdLhupBbjQ7ygkPOk66HAAtRPjliATGYBm/NK8fdv81HZ0IbrRmpx+5XD4KrkBD8iop4SBAEnKxvx06lq5BTVoaaxHc5OMiSEeWFSuA/GhXhiVICGhWiQsAD1kyMVILNZQPoRPVZ8m4/TVU2YEOqFO8cHwV/DK7YTEfWHIAgoqGnGoZJ6HC4z4kRFA9o6zVDKpRgZoEZ8kAdGB2gwOlCDcF9XyGUsRf3FAtRPjlCAOkxmfHmgDP/cdhqnKhsRN0yD28cGcaVTIqJB0mk2o6C6GScqGnC6qhFnq5tQbmgFAChkUkRq3TDSX40YfzVG6NwR46+GlysnVfcGC1A/2XMBqm1qx8d7i/DvrELoja24MtgDN8UFIlrnLnY0IiKH09zeicKaZhTWNKGgphkldc0orm1Bu8kMAPB1V2KEzh0jAtQYoVNjhL+ae4sugQWon+ytAHWazNh5qhqbsouRcbQCEkgwKcIbKaN0GM4zu4iIrIrJLEBvaEVRbTOKaptQXNuCotpmVDW2Afhlb9GoADVG+p8rRSMC1FCrnEROLj6bK0CrV6/G3/72N+j1esTFxeGNN97AhAkTLjp+06ZNeO6551BQUIDIyEgsX74cM2bMsHxfEAQsXboU77zzDurr6zFp0iS8+eabiIyM7FEeeyhArR0mZJ2pwbdHKvD14XLUN3cgyNMZ10b7YXKkD39RiIhsTGNb57lSVNOEwprmc/9f14xO07mP8WGezr8UIn81YnTuCPZycajT8Xvz+S36KT4bN25Eamoq1qxZg4SEBKxcuRIpKSnIz8+Hn5/fBeN37dqFOXPmIC0tDTfccAPWr1+PWbNmITc3F6NHjwYAvPLKK3j99dfx/vvvIzQ0FM899xxSUlJw9OhRqFSqoX6KQ6K1w4QjZUbsPVuL3WdqsOdMDVo7zfBzV+KqCB9MDPNGqI8rr9lFRGSj3JRyjPQ/t9fnvE6zGWX1rZZDaMW1zdhzthaGlg4AgMpJiig/d0Tp3BHp54YIPzeE+bphmKcznBz8MJroe4ASEhIwfvx4rFq1CgBgNpsRFBSEhx9+GM8888wF42fPno2mpiZs2bLFsm3ixImIj4/HmjVrIAgCAgIC8MQTT+DJJ58EABgMBmi1Wqxbtw533nnnZTNZ8x6gTpMZpfUtOFvdhNNVTcjXG3G0zIjj+gZ0mgWonKSI1rpjVIAGVwR7INDDmaWHiMiBCIKA+pYOFNee20tUUteC0rpmlNa3oqXDBACQSyUI9HRGiLcrgr1cEOjpjEAPZwR4OEOnUcHPXWmTBclm9gC1t7cjJycHixYtsmyTSqVITk5GVlZWt/fJyspCampql20pKSnYvHkzAODs2bPQ6/VITk62fF+j0SAhIQFZWVndFqC2tja0tbVZvjYYDADOvZADTRAEGFs60W4yocMkoLXThPYOM1o6OtHSYUZTWycaWzvR0NYJQ3MH6prbUdPUhkpjOyqMrahubIP558rqJJMi0EOFAA9n3DHGG2G+rgjycoFccv4vrQnNjY0D/hyIiMi6KQCEa2QI17gDoedOchEEAXXN7dAbW6E3tqHS2IKqegNOllSiurENrR3mLo+hUcnh5aaAl4sSXm5O0Dg7wV3lBHeVHO4qOVyc5HBRyuCskEMll0LlJIVSLoNCLoVCJoXTzzeZTAK5VAKZVAKZRDKoh+TOf273ZN+OqAWouroaJpMJWq22y3atVovjx493ex+9Xt/teL1eb/n++W0XG/O/0tLS8Pzzz1+wPSgoqGdPRERnxA5ARER2qVjsAP3Q0NAAjUZzyTGizwGyBosWLeqyV8lsNqO2thbe3t52d/jIaDQiKCgIxcXFVnd4z97xtRcPX3tx8fUXj6O99oIgoKGhAQEBAZcdK2oB8vHxgUwmQ0VFRZftFRUV0Ol03d5Hp9Ndcvz5/1ZUVMDf37/LmPj4+G4fU6lUQqlUdtnm4eHRm6dic9RqtUP8Mlgjvvbi4WsvLr7+4nGk1/5ye37OE3WGk0KhwNixY5GZmWnZZjabkZmZicTExG7vk5iY2GU8AGRkZFjGh4aGQqfTdRljNBqxZ8+eiz4mERERORbRD4GlpqZi3rx5GDduHCZMmICVK1eiqakJCxYsAADMnTsXgYGBSEtLAwA8+uijuOaaa7BixQrMnDkTGzZsQHZ2Nt5++20AgEQiwWOPPYa//vWviIyMtJwGHxAQgFmzZon1NImIiMiKiF6AZs+ejaqqKixZsgR6vR7x8fFIT0+3TGIuKiqCVPrLjqqkpCSsX78ezz77LBYvXozIyEhs3rzZsgYQAPz5z39GU1MT/vSnP6G+vh6TJ09Genq63a4B1BtKpRJLly694JAfDT6+9uLhay8uvv7i4Wt/caKvA0REREQ01GxvlSMiIiKifmIBIiIiIofDAkREREQOhwWIiIiIHA4LkIMoKCjAwoULERoaCmdnZ4SHh2Pp0qVob2/vMu7gwYO46qqroFKpEBQUhFdeeUWkxPblpZdeQlJSElxcXC66yGZRURFmzpwJFxcX+Pn54amnnkJnZ+fQBrVTq1evRkhICFQqFRISErB3716xI9mdHTt24MYbb0RAQAAkEonl+oznCYKAJUuWwN/fH87OzkhOTsbJkyfFCWtn0tLSMH78eLi7u8PPzw+zZs1Cfn5+lzGtra148MEH4e3tDTc3N9x2220XLCrsaFiAHMTx48dhNpvx1ltv4ciRI/jHP/6BNWvWYPHixZYxRqMR06ZNw/Dhw5GTk4O//e1vWLZsmWWNJeq79vZ23HHHHbj//vu7/b7JZMLMmTPR3t6OXbt24f3338e6deuwZMmSIU5qfzZu3IjU1FQsXboUubm5iIuLQ0pKCiorK8WOZleampoQFxeH1atXd/v9V155Ba+//jrWrFmDPXv2wNXVFSkpKWhtbR3ipPZn+/btePDBB7F7925kZGSgo6MD06ZNQ1NTk2XM448/ji+//BKbNm3C9u3bUVZWhltvvVXE1FZAIIf1yiuvCKGhoZav//nPfwqenp5CW1ubZdvTTz8tREdHixHPLq1du1bQaDQXbN+6dasglUoFvV5v2fbmm28KarW6y58H9d6ECROEBx980PK1yWQSAgIChLS0NBFT2TcAwmeffWb52mw2CzqdTvjb3/5m2VZfXy8olUrh448/FiGhfausrBQACNu3bxcE4dxr7eTkJGzatMky5tixYwIAISsrS6yYouMeIAdmMBjg5eVl+TorKwtXX301FAqFZVtKSgry8/NRV1cnRkSHkZWVhdjYWMsCoMC5195oNOLIkSMiJrNt7e3tyMnJQXJysmWbVCpFcnIysrKyREzmWM6ePQu9Xt/lz0Gj0SAhIYF/DoPAYDAAgOX9PScnBx0dHV1e/5iYGAQHBzv0688C5KBOnTqFN954A/fee69lm16v7/IBDMDytV6vH9J8joav/eCorq6GyWTq9rXl6zp0zr/W/HMYfGazGY899hgmTZpkuUKCXq+HQqG4YP6ho7/+LEA27plnnoFEIrnk7fjx413uU1paiunTp+OOO+7APffcI1Jy29eX156IaDA9+OCDOHz4MDZs2CB2FKsn+rXAqH+eeOIJzJ8//5JjwsLCLP9fVlaGKVOmICkp6YLJzTqd7oKzAs5/rdPpBiawHenta38pOp3ugjOT+Nr3n4+PD2QyWbd/r/m6Dp3zr3VFRQX8/f0t2ysqKhAfHy9SKvvz0EMPYcuWLdixYweGDRtm2a7T6dDe3o76+voue4Ec/feABcjG+fr6wtfXt0djS0tLMWXKFIwdOxZr167tcpFZAEhMTMRf/vIXdHR0wMnJCQCQkZGB6OhoeHp6Dnh2W9eb1/5yEhMT8dJLL6GyshJ+fn4Azr32arUaI0eOHJCf4YgUCgXGjh2LzMxMzJo1C8C5QwSZmZl46KGHxA3nQEJDQ6HT6ZCZmWkpPEajEXv27LnomZHUc4Ig4OGHH8Znn32Gbdu2ITQ0tMv3x44dCycnJ2RmZuK2224DAOTn56OoqAiJiYliRLYOYs/CpqFRUlIiRERECFOnThVKSkqE8vJyy+28+vp6QavVCnfffbdw+PBhYcOGDYKLi4vw1ltviZjcPhQWFgp5eXnC888/L7i5uQl5eXlCXl6e0NDQIAiCIHR2dgqjR48Wpk2bJuzfv19IT08XfH19hUWLFomc3PZt2LBBUCqVwrp164SjR48Kf/rTnwQPD48uZ9xR/zU0NFj+XgMQXn31VSEvL08oLCwUBEEQXn75ZcHDw0P4/PPPhYMHDwo333yzEBoaKrS0tIic3Pbdf//9gkajEbZt29blvb25udky5r777hOCg4OF77//XsjOzhYSExOFxMREEVOLjwXIQaxdu1YA0O3t1w4cOCBMnjxZUCqVQmBgoPDyyy+LlNi+zJs3r9vX/ocffrCMKSgoEK6//nrB2dlZ8PHxEZ544gmho6NDvNB25I033hCCg4MFhUIhTJgwQdi9e7fYkezODz/80O3f8Xnz5gmCcO5U+Oeee07QarWCUqkUpk6dKuTn54sb2k5c7L197dq1ljEtLS3CAw88IHh6egouLi7CLbfc0uUfwI5IIgiCMIQ7nIiIiIhEx7PAiIiIyOGwABEREZHDYQEiIiIih8MCRERERA6HBYiIiIgcDgsQERERORwWICIiInI4LEBERETkcFiAiIiIyOGwABEREZHDYQEiIiIih8MCRERERA7n/wGjBe7p5BIvOwAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "\n",
    "df_pa = pd.DataFrame(ATTE_estimates, columns=['Estimate'])\n",
    "g = sns.kdeplot(df_pa, fill=True)\n",
    "plt.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "doubleml",
   "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.11.2"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "c17b3f7bd1a966b23df6bc421cb393eedb89185fca3be80230a2e6676a5e11a1"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}