{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Python: Difference-in-Differences"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Remark:**\n",
    "*This notebook is based on the deprecated version of the DiD implementation of DoubleML. Please check out the [DiD Section](https://docs.doubleml.org/dev/guide/models.html#difference-in-differences-models-did) to find out about the current implementation.*"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "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)."
   ]
  },
  {
   "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": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "================== DoubleMLDIDData Object ==================\n",
      "Time variable: None\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",
      "\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "from doubleml.did.datasets import make_did_SZ2020\n",
    "from doubleml import DoubleMLDIDData\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 = DoubleMLDIDData.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": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "from lightgbm import LGBMClassifier, LGBMRegressor\n",
    "\n",
    "n_estimators = 30\n",
    "ml_g = LGBMRegressor(n_estimators=n_estimators, verbose=-1)\n",
    "ml_m = LGBMClassifier(n_estimators=n_estimators, verbose=-1)"
   ]
  },
  {
   "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": 11,
   "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",
      "\n",
      "------------------ Score & Algorithm ------------------\n",
      "Score function: observational\n",
      "\n",
      "------------------ Machine Learner   ------------------\n",
      "Learner ml_g: LGBMRegressor(n_estimators=30, verbose=-1)\n",
      "Learner ml_m: LGBMClassifier(n_estimators=30, verbose=-1)\n",
      "Out-of-sample Performance:\n",
      "Regression:\n",
      "Learner ml_g0 RMSE: [[11.55855122]]\n",
      "Learner ml_g1 RMSE: [[11.38983785]]\n",
      "Classification:\n",
      "Learner ml_m Log Loss: [[0.67674909]]\n",
      "\n",
      "------------------ Resampling        ------------------\n",
      "No. folds: 5\n",
      "No. repeated sample splits: 1\n",
      "\n",
      "------------------ Fit Summary       ------------------\n",
      "       coef  std err         t     P>|t|     2.5 %    97.5 %\n",
      "d  0.305272  1.25855  0.242559  0.808347 -2.161441  2.771986\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": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "      5.0 %    95.0 %\n",
      "d -1.764859  2.375403\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": 13,
   "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": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Coverage: 0.965\n",
      "Average CI length: 5.517510812139451\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": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGdCAYAAAD60sxaAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjYsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvq6yFwwAAAAlwSFlzAAAPYQAAD2EBqD+naQAAWyNJREFUeJzt3XtclHWiP/DPzMDMMFxmuM4AgoDcvIKKomamLYXVVnbZY/32rMXp2G61nTqsp83dTdutDlbmupWbe9r11k27mG1WlJFYKgqCeBcFucPAcJuBQWZgZn5/mGwkKJeBZy6f9+s1r81nnnnmM7LCh+f7fb6PyGaz2UBERETkRsRCByAiIiIaayxARERE5HZYgIiIiMjtsAARERGR22EBIiIiIrfDAkRERERuhwWIiIiI3A4LEBEREbkdD6EDOCKr1Yq6ujr4+vpCJBIJHYeIiIgGwWazob29HWFhYRCLr36OhwWoH3V1dYiIiBA6BhEREQ1DdXU1xo0bd9V9WID64evrC+DSX6Cfn5/AaYiIiGgwDAYDIiIien+OXw0LUD8uD3v5+fmxABERETmZwUxf4SRoIiIicjssQEREROR2WICIiIjI7XAOEBERuS2LxYLu7m6hY9AgSSQSeHh42GWJGhYgIiJySx0dHaipqYHNZhM6Cg2BQqFAaGgopFLpiI7DAkRERG7HYrGgpqYGCoUCwcHBXPTWCdhsNpjNZuh0OpSXlyMuLu6aix1eDQsQERG5ne7ubthsNgQHB8PLy0voODRIXl5e8PT0RGVlJcxmM+Ry+bCPxUnQRETktnjmx/mM5KxPn+PY5ShEREREToRDYERERN+rbbuIVqN5zN7P31uKcBWH4ITAAkRERIRL5ecnr+Siq9s6Zu8p9xQj5zcLBS1BW7ZswZNPPom2tjbBMgiBBYiIiAhAq9GMrm4rHlsUOyaFpLbtIjbsLUWr0Tyk93vwwQexdevWK7anp6cjOzv7qq+NiorCk08+iSeffLJ329KlS3HrrbcO+v2Hy9GKFgsQERHRD4SrvBAd5C10jKtavHgxNm/e3GebTCYb1rG8vLzc8ko4ToImIiJyMjKZDBqNps/D398fNpsNzz77LCIjIyGTyRAWFob/+q//AgAsXLgQlZWV+O///m+IRKLeK+C2bNkClUrVe+xnn30WycnJ2LRpEyIjI+Hj44NHH30UFosFL730EjQaDUJCQvDCCy/0ybRu3TpMnToV3t7eiIiIwKOPPoqOjg4AQG5uLjIyMqDX63vf+9lnnwUAmEwmrFixAuHh4fD29kZqaipyc3NH/e+QZ4CISFA2mw0FFa34qLAG+0ubMM7fCxND/fD/UiMRr/YVOh6RU/noo4/w5z//Gdu3b8fkyZOh1Wpx7NgxAMDOnTuRlJSEhx9+GMuXL7/qccrKyvDFF18gOzsbZWVluPfee3HhwgXEx8dj3759OHjwIP7jP/4DaWlpSE1NBXDp8vRXX30V0dHRuHDhAh599FE89dRT+Otf/4p58+Zh/fr1WLVqFUpKSgAAPj4+AIBf//rXOH36NLZv346wsDB8/PHHWLx4MU6cOIG4uLhR+7tiASIiwXR1W/DkjmJkn9Qi2FeGGZH+aOs049PjdXgvvwq/v20ifjFnPNdqIfqR3bt39xaIy373u99BLpdDo9EgLS0Nnp6eiIyMxOzZswEAAQEBkEgk8PX1hUajuerxrVYrNm3aBF9fX0yaNAmLFi1CSUkJPv/8c4jFYiQkJODFF1/E3r17ewvQD+cVRUVF4fnnn8evfvUr/PWvf4VUKoVSqYRIJOrz3lVVVdi8eTOqqqoQFhYGAFixYgWys7OxefNm/O///q89/rr6xQJERILQX+zG8q1HUFzdhv+6MQ6pMQEQf190zD1WvHO4Eqs+OYXiqja88m9JLEFEP7Bo0SK88cYbfbYFBATAaDRi/fr1iImJweLFi3Hrrbfi9ttvh4fH0H7cR0VFwdf3X2dg1Wo1JBJJn0UI1Wo1Ghsbe//89ddfIysrC2fPnoXBYEBPTw+6urrQ2dkJhULR7/ucOHECFosF8fHxfbabTCYEBgYOKfNQsQAR0ZjrtljxwKZ8lDZ24Pe3TbxiqEvqIUbGddGIDfHBX3PLEKf2xSMLJwiUlsjxeHt7IzY29ortAQEBKCkpwddff409e/bg0Ucfxcsvv4x9+/bB09Nz0Mf/8b4ikajfbVbrpSUDKioq8NOf/hSPPPIIXnjhBQQEBGD//v146KGHYDabByxAHR0dkEgkKCwshEQi6fPcj89w2RsLEBGNuVdzzuN4TRuevX0y4q4yz+f6uGDUtXXhpeyzSAz1xaKEkDFMSeScvLy8cPvtt+P222/HY489hsTERJw4cQIzZsyAVCqFxWKx+3sWFhbCarXilVde6T1L9P777/fZp7/3nj59OiwWCxobG3H99dfbPdfVsAAR0ZgqqGjBhr2luGfGuKuWn8t+ljIOVS1GPPHeUexdsRCBPsO71JdosGrbLjr8+5hMJmi12j7bPDw8sHv3blgsFqSmpkKhUODtt9+Gl5cXxo8fD+DS0Na3336L++67DzKZDEFBQSP6DJfFxsaiu7sbr732Gm6//XYcOHAAGzdu7LNPVFQUOjo6kJOTg6SkJCgUCsTHx+PnP/85li1bhldeeQXTp0+HTqdDTk4Opk2bhttuu80u+frDAkREY6ar24IntxcjLsQXdyaHD+o1YpEIv7xhAjLfL8bar84h6+6po5yS3JW/txRyTzE27C0ds/eUe4rh7y0d8uuys7MRGhraZ1tCQgLWrFmDNWvWIDMzExaLBVOnTsWnn37aO5/mT3/6E375y19iwoQJMJlMsNlsdvkcSUlJWLduHV588UWsXLkSCxYsQFZWFpYtW9a7z7x58/CrX/0KS5cuRXNzM1avXo1nn30WmzdvxvPPP4/f/OY3qK2tRVBQEObMmYOf/vSndsk2EJHNXp/ehRgMBiiVSuj1evj5+Qkdh8hl/G1fGV76sgQv3zsNocqhLbyWfbIe2/Iq8enj8zElXDlKCclddHV1oby8HNHR0ZDL5b3beS8wxzfQ1w4Y2s9vngEiojGhv9iNDbmlWJQQMuTyAwBpk9TIOduIZ/95Ch/8ai6vCqNREa7yYiFxE1wJmojGxN/2lcHUbcXdMwY39PVjHmIxls2NwpHKVuw53WDndETkbliAiGjUNRq6sGl/OW6ZooG/YujzHS6bGq5EosYXf80ts9vcBSJyTyxARDTqtuZVQCQS4afTwkZ8rNuTwlBc3Yb88hY7JCMid8UCRESjqqvbgncOV2FhQjC8ZSOfdjg9QoXIAC/8NbfMDunI3fFMovOx19eMBYiIRtWuo7XQd3YjffLV7z00WJfPJO07p8OZeoNdjknu5/Kqw2bz2F3xRfbR2dkJ4MrVqoeKV4ER0aix2WzYdKAcM8f7Q+0nv/YLBmnuhEC8f6Qamw+U46V7k+x2XHIfHh4eUCgU0Ol08PT07HOPK3JMNpsNnZ2daGxshEqluuLWGUPFAkREo+ZgWTPONXTg97dOtOtxPcRiLEoIwT+P1eEPP50EP/nIfhMk9yMSiRAaGory8nJUVlYKHYeGQKVSXfNu9oPBAkREo+btQ5WI8PfC5DD7Lyi6MCEEHxXV4JOjtfjF3Ci7H59cn1QqRVxcHIfBnIinp+eIz/xcxgJERKOi1WjG12cacN+syFFZtDDAW4qZ4/3x9uEq/Puc8VwYkYZFLBZfsZowuQeHGPTcsGEDoqKiIJfLkZqaivz8/EG9bvv27RCJRFiyZEmf7TabDatWrUJoaCi8vLyQlpaG8+fPj0JyIhrIp8frYLUB18Xa52aL/bkxMQQl2nYUV7eN2nsQkWsSvADt2LEDmZmZWL16NYqKipCUlIT09HQ0NjZe9XUVFRVYsWIFrr/++iuee+mll/Dqq69i48aNOHz4MLy9vZGeno6urq7R+hhE9CMfHKnB9AgVlF6jNz9nWrgKwb4yvJdfNWrvQUSuSfACtG7dOixfvhwZGRmYNGkSNm7cCIVCgU2bNg34GovFgp///Of44x//iJiYmD7P2Ww2rF+/Hn/4wx9w5513Ytq0adi2bRvq6uqwa9euUf40RAQAJdp2nKjVY0Fc8Ki+j1gswoK4YHx2vB4XzZZRfS8ici2CFiCz2YzCwkKkpaX1bhOLxUhLS0NeXt6Ar/vTn/6EkJAQPPTQQ1c8V15eDq1W2+eYSqUSqampAx7TZDLBYDD0eRDR8H1UVANfuQemR6pG/b2ujwuC0WzB12d4fzAiGjxBC1BTUxMsFgvUanWf7Wq1Glqttt/X7N+/H//4xz/w5ptv9vv85dcN5ZhZWVlQKpW9j4iIiKF+FCL6ntVqw8dHazE3JhAektH/FqP2kyNe7YOPj9aO+nsRkesQfAhsKNrb2/GLX/wCb775JoKC7DexcuXKldDr9b2P6upqux2byN0cqWyFrt2EeRNGb/Lzj10XG4R953Ro7jCN2XsSkXMT9DL4oKAgSCQSNDT0PXXd0NDQ7yJHZWVlqKiowO233967zWq1Ari0qmdJSUnv6xoaGhAaGtrnmMnJyf3mkMlkkMlkI/04RATgs+N1CPSWIk7tM2bvOScmENvyKvHZiXos45pARDQIgp4BkkqlmDlzJnJycnq3Wa1W5OTkYO7cuVfsn5iYiBMnTqC4uLj3cccdd2DRokUoLi5GREQEoqOjodFo+hzTYDDg8OHD/R6TiOzHYrXh85NazI4OgHgM1+Xxk3siaZwSO4s4DEZEgyP4QoiZmZl44IEHkJKSgtmzZ2P9+vUwGo3IyMgAACxbtgzh4eHIysqCXC7HlClT+rxepVIBQJ/tTz75JJ5//nnExcUhOjoazzzzDMLCwq5YL4iI7OtIRQt07SbMiQkc8/e+LjYIr31TiqrmTkQGKsb8/YnIuQhegJYuXQqdTodVq1ZBq9UiOTkZ2dnZvZOYq6qqhnyTuqeeegpGoxEPP/ww2traMH/+fGRnZ3O1T6JR9tmJegT6SBEbMnbDX5fNiPSHzEOMz0/W41c3TBjz9yci5yKy2Ww2oUM4GoPBAKVSCb1eDz8/+9/DiMgVWaw2pP7v15gdFSDYvbnWf30ORlMPdv/XlQukEpHrG8rPb6e6CoyIHNfRqlY0dZiRKsDw12Wp0YE4WWdAVXOnYBmIyDmwABGRXew53QCll6cgw1+XTY9UQeYhxmcn6gXLQETOgQWIiOziy1NazIhUjenVXz8m95QgOUKF3cfrBMtARM6BBYiIRqy0sQMVzZ2YOT5A6CiYExOIU3UGVDYbhY5CRA6MBYiIRmzP6QbIPMSYGq4UOgqSIy4Ng2Wf7P/WN0REAAsQEdnBntNaTA1XQuoh/LcUuacE08Yp8QULEBFdhfDfrYjIqenaTTha1YaZ4/2FjtJrVlQAiqvboNV3CR2FiBwUCxARjcjekkYAlxYidBQzIv3hIRbhy1M8C0RE/WMBIqIR2Xu2EbFqH/h5eQodpZe3zANTwv04DEZEA2IBIqJh67ZY8d35JiSPUwkd5QopUQHIL29Gi9EsdBQickAsQEQ0bIWVregw9SApQiV0lCukjA+AzXZpgjYR0Y+xABHRsOWW6KD08kR0kLfQUa6g9PJEosYXX55qEDoKETkgFiAiGra9ZxuRNE4p6OrPV5MSFYD955vQYeoROgoRORgWICIalrq2iyhpaEeyAw5/XZYy3h9mixX7SnRCRyEiB8MCRETDsu+cDmIRMNUBJ0BfFuInR1SgAl/xcngi+hEWICIaltySRsSpfeEj8xA6ylXNHO+Pb842wtxjFToKETkQFiAiGrIeixUHSpsxzQHu/XUts6IC0G7qQd6FZqGjEJEDYQEioiE7VtOGDlMPpo1z/AIUGaBAiK+Mw2BE1AcLEBEN2bfnmuAj80BMkI/QUa5JJBIhZbw/vjylhdVqEzoOETkIFiAiGrJvz+swJdwPYrFjXv7+YylRAWjqMONodZvQUYjIQbAAEdGQ6Du7cay6DVPDVUJHGbQEtS+UXp74iqtCE9H3WICIaEgOljXBagOmOsEE6MvEYhFmRKqQfVILm43DYETEAkREQ/Tt+SaEq7wQ7CsTOsqQpEQFoLK5E6WNHUJHISIHwAJERINms9nw7TkdpjjR2Z/LpoQp4eUpwZe8GoyIwAJERENQ1dKJ2raLTjX8dZnUQ4ykCCWyT7IAERELEBENwYHSZohFwMRQX6GjDMusqACcrDOgtu2i0FGISGAsQEQ0aAdKmzAhxAcKqWPf/mIgyREqeEpEXBSRiFiAiGhwrFYbDpQ1YUqY8w1/XaaQemBquBJfnGABInJ3LEBENChntAa0dXZjSpif0FFGJCUqAAUVLdC1m4SOQkQCYgEiokE5UNoEmYcYcWrnnP9z2czx/hCJgK/PNAgdhYgE5BAFaMOGDYiKioJcLkdqairy8/MH3Hfnzp1ISUmBSqWCt7c3kpOT8dZbb/XZ58EHH4RIJOrzWLx48Wh/DCKXdqC0GQkaX3hKHOLbxrD5yT0xMdSPV4MRuTnBv5Pt2LEDmZmZWL16NYqKipCUlIT09HQ0Njb2u39AQAB+//vfIy8vD8ePH0dGRgYyMjLw5Zdf9tlv8eLFqK+v73289957Y/FxiFySuceKw+XNmOzE839+aFZUAA6UNkF/sVvoKEQkEMEL0Lp167B8+XJkZGRg0qRJ2LhxIxQKBTZt2tTv/gsXLsRdd92FiRMnYsKECXjiiScwbdo07N+/v89+MpkMGo2m9+Hv7z8WH4fIJRVXt6Gr2+r0838uSxnvjx6rDd+c5TAYkbsStACZzWYUFhYiLS2td5tYLEZaWhry8vKu+XqbzYacnByUlJRgwYIFfZ7Lzc1FSEgIEhIS8Mgjj6C5uXnA45hMJhgMhj4PIvqXvLJmeMskiAr0FjqKXQT6yBCn9sHnvBqMyG0JWoCamppgsVigVqv7bFer1dBqB/7GpNfr4ePjA6lUittuuw2vvfYabrrppt7nFy9ejG3btiEnJwcvvvgi9u3bh1tuuQUWi6Xf42VlZUGpVPY+IiIi7PMBiVzEwbImTAr1g1gsEjqK3aRGBWJfiQ7tXRwGI3JHgg+BDYevry+Ki4tRUFCAF154AZmZmcjNze19/r777sMdd9yBqVOnYsmSJdi9ezcKCgr67PNDK1euhF6v731UV1ePzQchcgJd3RYUVbViUqhrDH9dNjs6AGaLFd+c7X++IRG5NkGXcw0KCoJEIkFDQ99x+IaGBmg0mgFfJxaLERsbCwBITk7GmTNnkJWVhYULF/a7f0xMDIKCglBaWoqf/OQnVzwvk8kgkznXna2JxkphZSu6LTZMcpEJ0JcF+8oQG+KDz0/U487kcKHjENEYE/QMkFQqxcyZM5GTk9O7zWq1IicnB3Pnzh30caxWK0ymgRc1q6mpQXNzM0JDQ0eUl8gdHSxrgtLLE+P8vYSOYnezowKQW6KD0dQjdBQiGmOCD4FlZmbizTffxNatW3HmzBk88sgjMBqNyMjIAAAsW7YMK1eu7N0/KysLe/bswYULF3DmzBm88soreOutt/Dv//7vAICOjg78z//8Dw4dOoSKigrk5OTgzjvvRGxsLNLT0wX5jETO7GBZMyaG+kIscp35P5elRgfA1MNhMCJ3JPgdDZcuXQqdTodVq1ZBq9UiOTkZ2dnZvROjq6qqIBb/q6cZjUY8+uijqKmpgZeXFxITE/H2229j6dKlAACJRILjx49j69ataGtrQ1hYGG6++WY899xzHOYiGqIOUw+OV+vxwLwooaOMihA/OWKCvfHZiXrcnhQmdBwiGkMim81mEzqEozEYDFAqldDr9fDzc62Jn0RDsbekERmbC/DKz5IQpnK9ITAA+PRYHT4qqkHhMzfBRyb474RENAJD+fkt+BAYETmuQxeaEaDwRKhSLnSUUTMnJhCmHiu+Ps1FEYncCQsQEQ0or6wZE0P9IHLB+T+XBfvKEK/2wafH6oSOQkRjiAWIiPrV3tWNU7UGTHSR219czZyYQOw7p4O+k4siErkLFiAi6teRylZYbDaXWwCxP6nRgbBYbfjyNG+NQeQuWICIqF+X5/9o/Fx3/s9lAd5STAz1w24OgxG5DRYgIupXXlkzEl18/s8PzYkJxIHSZjR1DLyoKhG5DhYgIrrC5fk/7jD8dVlqTAAgAj47Xi90FCIaAyxARHQFd5r/c5mf3BNJ45T4+Git0FGIaAywABHRFXrn/7jw+j/9uS42CMXVbahsNgodhYhGGQsQEV3h0IVmJLjR/J/LZkT6Q+4pxifFnAxN5OpYgIioD6OpBydrDJiocZ/hr8vknhKkjA/Ax0drwbsEEbk2FiAi6qPQDef//NB1sUEobzLiZK1B6ChENIpYgIioj/zyFqi8PBGmcq/5P5dNDVdC5eWJnUdrhI5CRKOIBYiI+jh0oRkJGl+3m/9zmUQswrzYIHxSXIdui1XoOEQ0SliAiKhXV7cFx2rakOiG839+aEFcEFqMZnx7Tid0FCIaJSxARNTraFUbui02TAz1FTqKoMYHemN8oAI7i7gmEJGrYgEiol6Hy5vhK/dARIBC6CiCmx8bhD2nG6C/yDvEE7kiFiAi6nXoQgsS1L4Qu+n8nx+aNyEIPVYrPj/BW2MQuSIWICICAJh7rCiuanX7+T+XBXhLMW2cEh8cqRY6ChGNAhYgIgIAnKhtQ1ePFYluPv/nhxbEBaOoqg0XdB1CRyEiO2MBIiIAwOHyFnh5ShAV6C10FIcxc3wAvKUSfFTENYGIXA0LEBEBAA5faEG82gcSMef/XCb1EGPuhEB8VFgLi5W3xiByJSxARASL1YbCSs7/6c8N8cHQGrpwsKxJ6ChEZEcsQESEM/UGdJh6OP+nHxOCfRCu8uJkaCIXwwJERMgvb4GnRIQJwT5CR3E4IpEIC+KDkX2yAfpOrglE5CpYgIgI+eXNiAvxhaeE3xL6c33cpTWB/nm8TugoRGQn/G5H5OZsNhvyK1qRoOHw10D8FVIkR6jwQQGHwYhcBQsQkZsr0xnRYjQjkQXoqm6ID8HxWj1KtO1CRyEiO2ABInJz+eUtEIuAeDUL0NXMiFTBT+7BydBELoIFiMjNFVS0ICbIG3JPidBRHJqHRIzrYoOw82gtzD1WoeMQ0QixABG5uUMXmhHP9X8G5Yb4YLQYzdhb0ih0FCIaIYcoQBs2bEBUVBTkcjlSU1ORn58/4L47d+5ESkoKVCoVvL29kZycjLfeeqvPPjabDatWrUJoaCi8vLyQlpaG8+fPj/bHIHI6tW0XUa/v4vyfQRof6I2YYG8OgxG5AMEL0I4dO5CZmYnVq1ejqKgISUlJSE9PR2Nj/79hBQQE4Pe//z3y8vJw/PhxZGRkICMjA19++WXvPi+99BJeffVVbNy4EYcPH4a3tzfS09PR1dU1Vh+LyCkUlLcAAK8AG4Ib4oKx96wOje38fkLkzAQvQOvWrcPy5cuRkZGBSZMmYePGjVAoFNi0aVO/+y9cuBB33XUXJk6ciAkTJuCJJ57AtGnTsH//fgCXzv6sX78ef/jDH3DnnXdi2rRp2LZtG+rq6rBr164x/GREji+/ogXj/L3gJ/cUOorTmDchCGIxsOtordBRiGgEBC1AZrMZhYWFSEtL690mFouRlpaGvLy8a77eZrMhJycHJSUlWLBgAQCgvLwcWq22zzGVSiVSU1MHPKbJZILBYOjzIHIHhy80I4FXfw2Jj9wDKeMD8P6RGthsvEEqkbMStAA1NTXBYrFArVb32a5Wq6HVagd8nV6vh4+PD6RSKW677Ta89tpruOmmmwCg93VDOWZWVhaUSmXvIyIiYiQfi8gptBjNKNMZkRjKCdBDtTAhGKWNHSiubhM6ChENk+BDYMPh6+uL4uJiFBQU4IUXXkBmZiZyc3OHfbyVK1dCr9f3PqqrOcGRXF9BxaX5P5wAPXRTwpQI9JHio6IaoaMQ0TB5CPnmQUFBkEgkaGho6LO9oaEBGo1mwNeJxWLExsYCAJKTk3HmzBlkZWVh4cKFva9raGhAaGhon2MmJyf3ezyZTAaZTDbCT0PkXArKWxDsI0OQD/+/P1RisQjXxwbhk+I6/OG2SVxDicgJCXoGSCqVYubMmcjJyendZrVakZOTg7lz5w76OFarFSaTCQAQHR0NjUbT55gGgwGHDx8e0jGJXF1+eQviefZn2BbEBaO9qwd7Tjdce2cicjiCngECgMzMTDzwwANISUnB7NmzsX79ehiNRmRkZAAAli1bhvDwcGRlZQG4NF8nJSUFEyZMgMlkwueff4633noLb7zxBgBAJBLhySefxPPPP4+4uDhER0fjmWeeQVhYGJYsWSLUxyRyKEZTD07VGfDgdVFCR3FaoSovJKh98WFhDW5PChM6DhENkeAFaOnSpdDpdFi1ahW0Wi2Sk5ORnZ3dO4m5qqoKYvG/TlQZjUY8+uijqKmpgZeXFxITE/H2229j6dKlvfs89dRTMBqNePjhh9HW1ob58+cjOzsbcrl8zD8fkSM6WtUGi83GK8BG6Pr4IGzaXw6tvgsaJb+/EDkTkY3XcV7BYDBAqVRCr9fDz49XyJDr+fOec9h0oBwb/30mxCKR0HGcVqe5B4+8XYT/vikejyycIHQcIrc3lJ/fTnkVGBGNTH55C+LVviw/I6SQeiAlyh8fFXFNICJnwwJE5Ga6LVYcrW7l5e92cn1cEEobO3CqjguoEjkTFiAiN3OyVo+ubivn/9jJ1HAVVF6eXBOIyMmwABG5mSMVrZB5iBEd5C10FJcgEYswb0Ig/llch26LVeg4RDRILEBEbia/ogWxIT7wkPCfv73MjwtGs9GM787rhI5CRIPE74BEbsRqtaGgogUJnP9jV1GBCkQGeOGjQt4hnshZsAARuZELTR1o6+zm/B87E4lEuG5CEL4+04AOU4/QcYhoEFiAiNxIfnkrxCIgngXI7ubFBsHUY8WXJ7VCRyGiQWABInIjBRUtiA7y5s07R0GQjwwTNb74pJjDYETOgAWIyI1cXgCRRse82CAcKG2Grt0kdBQiugYWICI3Ua+/iNq2i0jU8PYuoyU1OgAQAZ8drxM6ChFdAwsQkZsoqGgFAMSrfQRO4rp85Z5IHqfCrmIWICJHxwJE5CYKylsQppJDpZAKHcWlzYsNRHF1G6pbOoWOQkRXwQJE5CYOlzfz8vcxMCPSHzIPMXYfrxc6ChFdBQsQkRvQX+zG+YYOLoA4BuSeEkyPVOHTYxwGI3JkLEBEbqCwsgU2gBOgx8ic6ECcrjegvMkodBQiGgALEJEbyC9vRYDCEyG+MqGjuIXkSBXknmLs5lkgIofFAkTkBgoqWhCv8YVIJBI6iluQeUgwI9Kf84CIHBgLEJGL6+q24HhNGxLUHP4aS3NjAlHS0I7SxnahoxBRP1iAiFzcseo2dFtsSAzlBOixNG2cCl6eEnx+gvcGI3JELEBELq6gogUKqQSR/gqho7gVqYcY0yNV+PwEh8GIHBELEJGLK6hoRbzaB2Ix5/+MtdnRATirbUcFrwYjcjgsQEQuzGK1obCylfN/BJI0TgWZhxhfnOQwGJGjYQEicmFn6g3oMPUgkQsgCkLuKUHSOBW+4DAYkcNhASJyYQUVLfCQiBATzBugCmV2dACO1+pR23ZR6ChE9AMsQEQuLL+8BbHBPpB68J+6UKZHquApESGbw2BEDoXfFYlclM1mQ35FC+//JTCF1ANTw5X4kgWIyKGwABG5qIrmTjR3mDn/xwGkjA/AkcoWNHeYhI5CRN9jASJyUfnlzRABiFezAAltxnh/2GxAztlGoaMQ0fdYgIhcVH55K6KCFFBIPYSO4vaUXp5I0Pjiq1McBiNyFA5RgDZs2ICoqCjI5XKkpqYiPz9/wH3ffPNNXH/99fD394e/vz/S0tKu2P/BBx+ESCTq81i8ePFofwwih5Jf3ox4rv/jMGaO98d355vQae4ROgoRwQEK0I4dO5CZmYnVq1ejqKgISUlJSE9PR2Nj/6eKc3Nzcf/992Pv3r3Iy8tDREQEbr75ZtTW1vbZb/Hixaivr+99vPfee2PxcYgcQoOhC9WtFzn/x4GkjA+AqceKb881CR2FiDDMAnThwgW7BVi3bh2WL1+OjIwMTJo0CRs3boRCocCmTZv63f+dd97Bo48+iuTkZCQmJuLvf/87rFYrcnJy+uwnk8mg0Wh6H/7+/nbLTOTo8stbAIAFyIFolHJEBig4DEbkIIZVgGJjY7Fo0SK8/fbb6OrqGvabm81mFBYWIi0t7V+BxGKkpaUhLy9vUMfo7OxEd3c3AgIC+mzPzc1FSEgIEhIS8Mgjj6C5uXnAY5hMJhgMhj4PImdWUNGCUKUcKoVU6Cj0AzMi/ZFzthE9FqvQUYjc3rAKUFFREaZNm4bMzExoNBr88pe/vOq8nYE0NTXBYrFArVb32a5Wq6HVDu63pN/+9rcICwvrU6IWL16Mbdu2IScnBy+++CL27duHW265BRaLpd9jZGVlQalU9j4iIiKG/FmIHMmhC81I4NVfDmfmeH/oL3ajsLJV6ChEbm9YBSg5ORl/+ctfUFdXh02bNqG+vh7z58/HlClTsG7dOuh0Onvn7NeaNWuwfft2fPzxx5DL5b3b77vvPtxxxx2YOnUqlixZgt27d6OgoAC5ubn9HmflypXQ6/W9j+rq6jHJTzQa2jrNONfQgcRQToB2NDHB3vBXePJyeCIHMKJJ0B4eHrj77rvxwQcf4MUXX0RpaSlWrFiBiIgILFu2DPX1V78BYFBQECQSCRoaGvpsb2hogEajuepr165dizVr1uCrr77CtGnTrrpvTEwMgoKCUFpa2u/zMpkMfn5+fR5Ezqqg4tLZhYmc/+NwxCIRkiP8sed0w7V3JqJRNaICdOTIETz66KMIDQ3FunXrsGLFCpSVlWHPnj2oq6vDnXfeedXXS6VSzJw5s88E5ssTmufOnTvg61566SU899xzyM7ORkpKyjVz1tTUoLm5GaGhoYP/cEROKr+8GYE+UgT7yoSOQv2YMV6F8iYjLug6hI5C5NaGtULaunXrsHnzZpSUlODWW2/Ftm3bcOutt0IsvtSnoqOjsWXLFkRFRV3zWJmZmXjggQeQkpKC2bNnY/369TAajcjIyAAALFu2DOHh4cjKygIAvPjii1i1ahXeffddREVF9c4V8vHxgY+PDzo6OvDHP/4R99xzDzQaDcrKyvDUU08hNjYW6enpw/m4RE7lcHkLEtW+EIlEQkehfkwNV0IqESPnTCNign2EjkPktoZVgN544w38x3/8Bx588MEBz6qEhITgH//4xzWPtXTpUuh0OqxatQparRbJycnIzs7unRhdVVXVW6wuv7fZbMa9997b5zirV6/Gs88+C4lEguPHj2Pr1q1oa2tDWFgYbr75Zjz33HOQyfgbMbk2o6kHp2oNePC6KKGj0ABkHhJMCffDntMNWL4gRug4RG5LZLPZbEN9UUVFBSIjI/sUE+DS3aerq6sRGRlpt4BCMBgMUCqV0Ov1nA9ETuW78zr84h/5WHtvEsL9vYSOQwPIOdOATQfKUfTMTVyqgMiOhvLze1hzgCZMmICmpitXM21paUF0dPRwDklEdpBf3gI/uQfCVPJr70yCmR7pD6sN2HdubK6YJaIrDasADXTSqKOjo8/l6EQ0tg5faEGChvN/HF2AtxTRQd74hpfDEwlmSHOAMjMzAQAikQirVq2CQqHofc5iseDw4cNITk62a0AiGpyubguKq9uwdBYX8nQGyREqfPP9qtAeEsFvy0jkdoZUgI4ePQrg0hmgEydOQCr919i1VCpFUlISVqxYYd+ERDQoxdVtMFusmMgFEJ3C9AgVPj5ai6PVbZgVFXDtFxCRXQ2pAO3duxcAkJGRgb/85S+cIEzkQA5faIG3VILxAYpr70yCmxDsA6WXJ74528gCRCSAYZ133bx5M8sPkYM5dKEZCRpfiMWc/+MMxGIRpo1TIucMV4UmEsKgzwDdfffd2LJlC/z8/HD33Xdfdd+dO3eOOBgRDZ65x4qjVa24Z+Y4oaPQEEyP8Mer35xHTWsnxvnzzB3RWBp0AVIqlb1XliiVylELRERDd7ymDV09nP/jbJIilJCIRdh7thG/mBsldBwitzLoArR58+Z+/5uIhHe4vAVenhJEBXoLHYWGQCH1QILaF3tLdCxARGNsWHOALl68iM7Ozt4/V1ZWYv369fjqq6/sFoyIBu/S/B8fSDj/x+kkRahwsLQJXd0WoaMQuZVhFaA777wT27ZtAwC0tbVh9uzZeOWVV3DnnXfijTfesGtAIrq6HosVRypakajh8Jczmh6hQlePFYfLW4SOQuRWhlWAioqKcP311wMAPvzwQ2g0GlRWVmLbtm149dVX7RqQiK7uRK0eF7stnP/jpMb5eyHIR4q9XBWaaEwNqwB1dnbC19cXAPDVV1/h7rvvhlgsxpw5c1BZWWnXgER0dYcutEDuKUZMMOf/OCORSITkCBULENEYG1YBio2Nxa5du1BdXY0vv/wSN998MwCgsbGR6wMRjbFDF5qRqPGFh5i3U3BWyRH+qGzpRHmTUegoRG5jWN8xV61ahRUrViAqKgqpqamYO3cugEtng6ZPn27XgEQ0sG6LFfnlLZjI+T9ObXKYHzwkIt4clWgMDasA3XvvvaiqqsKRI0eQnZ3du/0nP/kJ/vznP9stHBFd3eX5P5PCWICcmdxTgkmhfhwGIxpDQ7oX2A9pNBpoNJo+22bPnj3iQEQ0eIcuNMPLU4LoIB+ho9AIJY1TYXtBFTrNPVBIh/2tmYgGaVj/yoxGI9asWYOcnBw0NjbCarX2ef7ChQt2CUdEV3eo7NL9v7j+j/ObHqHCW4cqcbC0GWmT1ELHIXJ5wypA//mf/4l9+/bhF7/4BUJDQ3tvkUFEY6fbYkVBRSuWTA8XOgrZgUYph9pPhtxzjSxARGNgWAXoiy++wGeffYbrrrvO3nmIaJCO11ya/zOZ839cgkgkQtI4Fb452wibzcZfLIlG2bAmQfv7+yMgIMDeWYhoCA5daIZCyvt/uZLkCBXq2rpQpusQOgqRyxtWAXruueewatWqPvcDI6KxdaC0CRM1fpz/40ImhflBKhFj71md0FGIXN6whsBeeeUVlJWVQa1WIyoqCp6enn2eLyoqsks4IupfV7cFhZWt+LeUCKGjkB3JPCSYFOaLb842YvmCGKHjELm0YRWgJUuW2DkGEQ3F0ao2mHqsnP/jgpLG+eOdw5XoMPXAR8bL4YlGy7D+da1evdreOYhoCA6WNcFP7oGIAIXQUcjOkiNU2JpXgYOlTbh5subaLyCiYRn2zYPa2trw97//HStXrkRLSwuAS0NftbW1dgtHRP07UNqESWF+EPNKIZejUcoRppJjbwnnARGNpmGdATp+/DjS0tKgVCpRUVGB5cuXIyAgADt37kRVVRW2bdtm75xE9L0OUw+O1+ixbG6U0FFolEwbd+nu8Lwcnmj0DOsMUGZmJh588EGcP38ecrm8d/utt96Kb7/91m7hiOhKBRUt6LHaMIXzf1zW9AgVtIYunGvg5fBEo2VYBaigoAC//OUvr9geHh4OrVY74lBENLC8smYEekuhUcqvvTM5pUSNH2QeYuwt4c1RiUbLsAqQTCaDwWC4Yvu5c+cQHBw84lBENLDvzuswKdSPQyMuTOohxpQwJe8OTzSKhlWA7rjjDvzpT39Cd3c3gEtLuFdVVeG3v/0t7rnnniEfb8OGDYiKioJcLkdqairy8/MH3PfNN9/E9ddfD39/f/j7+yMtLe2K/W02G1atWoXQ0FB4eXkhLS0N58+fH3IuIkfT3GHCmfp2TB2nFDoKjbKkCCUKK1th6OoWOgqRSxpWAXrllVfQ0dGB4OBgXLx4ETfccANiY2Ph6+uLF154YUjH2rFjBzIzM7F69WoUFRUhKSkJ6enpaGzs/zef3Nxc3H///di7dy/y8vIQERGBm2++uc/VZy+99BJeffVVbNy4EYcPH4a3tzfS09PR1dU1nI9L5DAOljUDACaHsQC5uuQIf/RYbThwvknoKEQuSWSz2WzDffGBAwdw7NgxdHR0YMaMGUhLSxvyMVJTUzFr1iy8/vrrAACr1YqIiAg8/vjjePrpp6/5eovFAn9/f7z++utYtmwZbDYbwsLC8Jvf/AYrVqwAAOj1eqjVamzZsgX33XffNY9pMBigVCqh1+vh58eJpuQ4nv7oOPaXNuHle5OEjkJj4KkPj2FOTCBe/hm/3kSDMZSf30O+DN5qtWLLli3YuXMnKioqIBKJEB0dDY1GM+RLNs1mMwoLC7Fy5crebWKxGGlpacjLyxvUMTo7O9Hd3d17c9by8nJotdo+ZUypVCI1NRV5eXn9FiCTyQSTydT75/7mNxEJzWaz4dvzOkwLVwkdhcZIUoQKe0saYbXaIOY934jsakhDYDabDXfccQf+8z//E7W1tZg6dSomT56MyspKPPjgg7jrrruG9OZNTU2wWCxQq9V9tqvV6kFfTfbb3/4WYWFhvYXn8uuGcsysrCwolcreR0QE769EjqeyuRN1bV2YGs7hL3cxPUKFpg4zTtXxlzIiexvSGaAtW7bg22+/RU5ODhYtWtTnuW+++QZLlizBtm3bsGzZMruGHMiaNWuwfft25Obm9lmPaKhWrlyJzMzM3j8bDAaWIHI4+0ubIBYBE0M5LOsu4jW+UEgl2FvSyInvRHY2pDNA7733Hn73u99dUX4A4MYbb8TTTz+Nd955Z9DHCwoKgkQiQUNDQ5/tDQ0N0Giufg+ctWvXYs2aNfjqq68wbdq03u2XXzeUY8pkMvj5+fV5EDma/aVNiFP7wksqEToKjREPsRhTw5XIOdtw7Z2JaEiGVICOHz+OxYsXD/j8LbfcgmPHjg36eFKpFDNnzkROTk7vNqvVipycHMydO3fA17300kt47rnnkJ2djZSUlD7PXZ6P9MNjGgwGHD58+KrHJHJkFqsNB0ubuPqzG5oeqcLxaj2aO0zX3pmIBm1IBailpeWKuTU/pFar0draOqQAmZmZePPNN7F161acOXMGjzzyCIxGIzIyMgAAy5Yt6zNJ+sUXX8QzzzyDTZs2ISoqClqtFlqtFh0dl5aMF4lEePLJJ/H888/jn//8J06cOIFly5YhLCwMS5YsGVI2IkdxvKYNhq4eTBunEjoKjbGkcSrYAOw7x5ujEtnTkOYAWSwWeHgM/BKJRIKenp4hBVi6dCl0Oh1WrVoFrVaL5ORkZGdn9xatqqoqiMX/6mlvvPEGzGYz7r333j7HWb16NZ599lkAwFNPPQWj0YiHH34YbW1tmD9/PrKzs0c0T4hISN+ea4K3VIIJwT5CR6ExplJIMSHYG9+cbcTdM8YJHYfIZQxpHSCxWIxbbrkFMpms3+dNJhOys7NhsVjsFlAIXAeIHM09bxyERCTCf98UL3QUEsCHhTX46pQWRatugqdkWOvXErmFUVsH6IEHHrjmPmN1BRiRuzB0daO4qg0Z10UJHYUEMiNShY+KanCkohVzJwQKHYfIJQypAG3evHm0chDRAA6WNsNis3H9HzcWFeQNf4UnvjnbwAJEZCc8l0rk4L47r0OoUo4QP85hc1dikQjJEf74+gzvDk9kLyxARA7u23M6nv0hzIhUobzJiPImo9BRiFwCCxCRA6toMqK69SIvfydMCVfCQyLCN2d5FojIHliAiBxYbkkjPCUiTOYCiG5P7inB5FA/fH2aq0IT2QMLEJED21uiQ6LGD3JP3v6CgBnj/VFQ0QL9xW6hoxA5PRYgIgfV1W3BoQvNmMabYNL3Zkb6o8dq46rQRHbAAkTkoA5daIapx4rkCJXQUchBBPrIEB3kja9Pa4WOQuT0WICIHFRuiQ7BPjKEq7yEjkIOZEakCnvP6tBtsQodhcipsQAROajckkZMG6eESCQSOgo5kJnjA9Bu6kFBRYvQUYicGgsQkQOqau5ERXMnkjj8RT8SFahAoLcUX5/m5fBEI8ECROSA9pY0wkMswpQwToCmvkQiEaZH+uOr01oM4V7WRPQjLEBEDuibs42YGOoHLykvf6crpYz3R03rRZQ0tAsdhchpsQAROZhOcw/yypp59RcNaHKYHxRSCb46xUURiYaLBYjIwRwobYbZYsX0SJXQUchBeUjESI5QIfskL4cnGi4WICIH883ZRoQq5QhV8vJ3GljK+ACcrjegprVT6ChETokFiMiB2Gw2fHO2gcNfdE1JEUp4SkQcBiMaJhYgIgdyut6ABoMJ0yP9hY5CDk4h9cDkMCW+PMVhMKLhYAEiciB7zzbCy1OCiRpfoaOQE0j5/uaoLUaz0FGInA4LEJED+fpMI6aGK+Eh4T9NuraZ4/1hswFfn+YwGNFQ8bsskYNobO/Cseo2zBjP4S8aHJVCikSNLz4/WS90FCKnwwJE5CD2nm2ESARM5wRoGoJZ0QHYf74Jhq5uoaMQORUWICIH8fWZRsSrfeHn5Sl0FHIis6MC0GO14ZszvDcY0VCwABE5gK5uC747p8MMXv1FQxToI0NsiA++4DAY0ZCwABE5gINlTejqsXL+Dw3LrKgA7CvRodPcI3QUIqfBAkTkAPacvrT6c5hSLnQUckKp0QHo6rEit0QndBQip8ECRCQwq9WGr880YHqkP0QikdBxyAmp/eSIDvLG5yc4DEY0WCxARAI7VtMGXbsJszj8RSMwOzoAX59pwEWzRegoRE6BBYhIYF+dboCf3APxaq7+TMM3JzoQXd1W7C3h1WBEgyF4AdqwYQOioqIgl8uRmpqK/Pz8Afc9deoU7rnnHkRFRUEkEmH9+vVX7PPss89CJBL1eSQmJo7iJyAamS9PaTEj0h9iMYe/aPg0Sjligryx+3id0FGInIKgBWjHjh3IzMzE6tWrUVRUhKSkJKSnp6Oxsf/fYDo7OxETE4M1a9ZAo9EMeNzJkyejvr6+97F///7R+ghEI1La2IELOiNmRnH4i0YuNToA35xp5NVgRIMgaAFat24dli9fjoyMDEyaNAkbN26EQqHApk2b+t1/1qxZePnll3HfffdBJpMNeFwPDw9oNJreR1BQ0Gh9BKIR2XO6ATIPMaaFq4SOQi4gNSYQXT1WfHOWw2BE1yJYATKbzSgsLERaWtq/wojFSEtLQ15e3oiOff78eYSFhSEmJgY///nPUVVVNdK4RKPiy1NaTBunhNRD8NFocgFqPzligr2x+xiHwYiuRbDvuk1NTbBYLFCr1X22q9VqaLXaYR83NTUVW7ZsQXZ2Nt544w2Ul5fj+uuvR3t7+4CvMZlMMBgMfR5Eo63R0IXi6jakjA8QOgq5kDnRgfjmrA7tvDcY0VW53K+dt9xyC372s59h2rRpSE9Px+eff462tja8//77A74mKysLSqWy9xERETGGicldfXW6AWIRMD1SJXQUciFzJwTCbLHiq1MNQkchcmiCFaCgoCBIJBI0NPT9R9rQ0HDVCc5DpVKpEB8fj9LS0gH3WblyJfR6fe+jurrabu9PNJDsk1pMDlPCV86bn5L9BPnIkKjxxT85DEZ0VYIVIKlUipkzZyInJ6d3m9VqRU5ODubOnWu39+no6EBZWRlCQ0MH3Ecmk8HPz6/Pg2g0tXWacehCM1J49ReNgrkTArG/tAktRrPQUYgclqBDYJmZmXjzzTexdetWnDlzBo888giMRiMyMjIAAMuWLcPKlSt79zebzSguLkZxcTHMZjNqa2tRXFzc5+zOihUrsG/fPlRUVODgwYO46667IJFIcP/994/55yMaSM6ZRvRYbZz/Q6NiTnQgbDYbb41BdBUeQr750qVLodPpsGrVKmi1WiQnJyM7O7t3YnRVVRXE4n91tLq6OkyfPr33z2vXrsXatWtxww03IDc3FwBQU1OD+++/H83NzQgODsb8+fNx6NAhBAcHj+lnI7qaL09pEa/2QYC3VOgo5IL8vDwxNVyJT4/V4d/njBc6DpFDEtlsNpvQIRyNwWCAUqmEXq/ncBjZXae5B9P/tAf3zBiH25PChI5DLmrfOR3+tq8MB56+EWEqL6HjEI2Jofz8drmrwIgc3b4SHUw9VsyK4vAXjZ5ZUf7wlIg5GZpoACxARGPsi5P1GB+ogEYpFzoKuTCF1AMzxqvw8dFaoaMQOSQWIKIx1NVtwddnGjGbZ39oDMyPDUaJth1ntVzclejHWICIxtD+803oNFswO5oFiEZf0jglfGUe2HWUw2BEP8YCRDSGPj9Zj3H+XhjnrxA6CrkBD4kYqTGB+KS4FlYrr3ch+iEWIKIxYu6xYs/pBk5+pjE1PzYI9fouHC5vEToKkUNhASIaIwfLmtDe1YNUDn/RGIpX+0DtJ8MuToYm6oMFiGiMZJ/UIlQpR2QAh79o7IhEIlwXG4Tdx+vQ1W0ROg6Rw2ABIhoD3RYrsk9pMSsqACKRSOg45Gaujw2G0WzBV6d5h3iiy1iAiMbA4QstaOvs5vAXCUKjlCNe7YOdRTVCRyFyGCxARGPgsxP1UPvJEB3kLXQUclPzY4Px3bkmNLZ3CR2FyCGwABGNsh6LFdkn6zGbw18koLkxgRCLgX8Wc00gIoAFiGjU5Ze3oLWzG6kxgUJHITfmI/fA9Ah/fFjIYTAigAWIaNR9dqIeIb4yxHD4iwS2ID4YZ7XtOF3HW2MQsQARjSKL1Ybsk1rMjubwFwkvKUIJpZcnPuJkaCIWIKLRdLi8Gc1GM1KjOfxFwvMQizFvQiA+PlqLbotV6DhEgmIBIhpFnx2vR7CvDBOCOfxFjmFBfDBajGbsK9EJHYVIUCxARKOkx2LFFye1SOXwFzmQqEBvRAUq8CGHwcjNsQARjZLD5S1oMZoxh1d/kYO5Pi4YX59uQKvRLHQUIsGwABGNkt3HLy1+yKu/yNHMjw2CDcAnxbxBKrkvFiCiUXBp+KseqdGBHP4ih+Pn5YkZkSq8f4TDYOS+WICIRsHBsmbe+4sc2oL4YJyuN3BNIHJbLEBEo2D38Tpo/OS89xc5rOQIFZRenlwZmtwWCxCRnZl7Ll39NSeGw1/kuDzEYsyPDcLHR2tg7uGaQOR+WICI7Gx/qQ7tXT2YO4FXf5FjuyE+GK2d3fjmbKPQUYjGHAsQkZ3tPlaPcf5eiPD3EjoK0VVFBCgQG+KD949UCx2FaMyxABHZUVe3BV+e0vLqL3IaC+KCkVvSiAZDl9BRiMYUCxCRHeWW6GA0Wzj8RU5j3oRAeIjF2FnENYHIvbAAEdnR7uN1iApUIFzF4S9yDt4yD8yODsD7R6phs9mEjkM0ZliAiOzEaOrB16cbeOsLcjo3xAejvMmIwspWoaMQjRkWICI7+fpMA7p6rJjH4S9yMpPC/BDiK+NkaHIrghegDRs2ICoqCnK5HKmpqcjPzx9w31OnTuGee+5BVFQURCIR1q9fP+JjEtnLP4vrEKf2QbCvXOgoREMiFomwID4Ynx6rR4epR+g4RGNC0AK0Y8cOZGZmYvXq1SgqKkJSUhLS09PR2Nj/mhSdnZ2IiYnBmjVroNFo7HJMInvQd3Zj3zkd5nH4i5zUDfHB6Oq24PPj9UJHIRoTghagdevWYfny5cjIyMCkSZOwceNGKBQKbNq0qd/9Z82ahZdffhn33XcfZDKZXY5JZA/Zp+phtdmQygJETirIR4ap45TYXlAldBSiMSFYATKbzSgsLERaWtq/wojFSEtLQ15e3pge02QywWAw9HkQDcUnxXWYFOoHf4VU6ChEw7YwPhhFVW0obewQOgrRqBOsADU1NcFisUCtVvfZrlarodVqx/SYWVlZUCqVvY+IiIhhvT+5p0ZDFw5daMYcTn4mJ5cSFQBfmQcnQ5NbEHwStCNYuXIl9Hp976O6mv/4afA+PV4PiViE1GgWIHJunhIxrosNwkeFNei28Aap5NoEK0BBQUGQSCRoaGjos72hoWHACc6jdUyZTAY/P78+D6LB+uRoLZIjVPCReQgdhWjEFiYEo9loRs4ZXjhCrk2wAiSVSjFz5kzk5OT0brNarcjJycHcuXMd5phEV1PeZMTxWj3mTQgSOgqRXYwP9EZssDcnQ5PLE/RX1szMTDzwwANISUnB7NmzsX79ehiNRmRkZAAAli1bhvDwcGRlZQG4NMn59OnTvf9dW1uL4uJi+Pj4IDY2dlDHJLKnfxbXwctTghmR/kJHIbKbhQkh2HSgHPX6iwhV8rYu5JoELUBLly6FTqfDqlWroNVqkZycjOzs7N5JzFVVVRCL/3WSqq6uDtOnT+/989q1a7F27VrccMMNyM3NHdQxiezFZrNhV3EtZkX5Q+rB6XTkOuZOCMRbhyrx4ZEaPP6TOKHjEI0KkY13v7uCwWCAUqmEXq/nfCAa0PGaNtzx+gGsvCUR08aphI5DZFcb95WhtLED3z21CGKxSOg4RIMylJ/f/LWVaJh2FtXCX+GJKWFKoaMQ2d2NiSGobbuIg2XNQkchGhUsQETD0G2x4p/H6jBvQhB/OyaXFBfig3H+Xngvn5OhyTWxABENw3fndWgxmjE/jld/kWsSiURYlBCCL09p0dxhEjoOkd2xABENw86iWkQGeGF8gELoKESjZn5cEESiS/9/J3I1LEBEQ2To6sae0w24LjYYIhGHv8h1+ck9MSsqAO8VVIHXy5CrYQEiGqLsE1qYe6y4jvf+IjewKCEEF3RGFFS0Ch2FyK5YgIiG6P0j1ZgarkSgj0zoKESjblKYH0KVcrx7uFLoKER2xQJENAQVTUYcqWzFgvhgoaMQjQnx95OhPz+hRavRLHQcIrthASIago+KaqCQSjArKkDoKERjZkF8MKw2Gz4qqhE6CpHdsAARDZLFasOHhTWYGxPIW1+QW1F6XZoM/W4+J0OT6+B3caJByitrRr2+Czdw+Ivc0E8mXpoMnV/eInQUIrtgASIapPePVCFM5YXYEB+hoxCNuUmhfghTyvH2IU6GJtfAAkQ0CK1GM744qcWiBK79Q+5JJBLhJxPV+OKkFrp2rgxNzo8FiGgQdh6thc0GLIjj8Be5rwVxwRCLRHj/SLXQUYhGjAWI6BpsNhvey69CSpQ//Lw8hY5DJBgfuQfmTgjEO4crYbFyMjQ5NxYgomsoqmpFaWMHFiWECB2FSHA3TVKjrq0LuSWNQkchGhEWIKJreC+/Gmo/GaaEK4WOQiS4CcE+iAn2xrY8ToYm58YCRHQV+s5u7D5Wh4XxIRBz8jMRAOCmiWrsO6dDeZNR6ChEw8YCRHQVHxRWo8dqw8IETn4mumzehCD4yj2wLa9C6ChEw8YCRDQAq9WGbXmVmBMTCJVCKnQcIoch9RBjUUIIPjhSA6OpR+g4RMPCAkQ0gH3ndahq6cRNk9RCRyFyODdNUqPT3IOdR2uFjkI0LCxARAPYdrAC0UHeiOPKz0RXCPKRIWV8ALYcKOf9wcgpsQAR9aOy2YjcEh1umqjmys9EA1g8RYMynRH7zumEjkI0ZCxARP34x/5y+Hl54rrYIKGjEDmsRI0vYoK98ffvyoWOQjRkLEBEP9JqNOP9gmqkTVRD6sF/IkQDEYlEuHVKKPaXNuGs1iB0HKIh4Xd3oh9553AlrDbgZk5+Jrqm1JgABPpI8Q+eBSInwwJE9ANd3RZsPliBBfFBvO8X0SB4iMVIn6TBruJaNLZ3CR2HaNBYgIh+4OOjtWjpMOPWKaFCRyFyGjcmhsBDLMbmAxVCRyEaNBYgou91W6x4/ZtSpMYEIFTlJXQcIqfhLfNA2sQQvJVXCUNXt9BxiAaFBYjoex8frUVt20XcNX2c0FGInM4tU0Nh6rHgLd4klZyEQxSgDRs2ICoqCnK5HKmpqcjPz7/q/h988AESExMhl8sxdepUfP75532ef/DBByESifo8Fi9ePJofgZxcz/dnf2ZHBSAyQCF0HCKn46+QYkFcMDbtL0dXt0XoOETXJHgB2rFjBzIzM7F69WoUFRUhKSkJ6enpaGxs7Hf/gwcP4v7778dDDz2Eo0ePYsmSJViyZAlOnjzZZ7/Fixejvr6+9/Hee++NxcchJ/VJcR2qWjpx14xwoaMQOa3bk8LQ2mnGjoJqoaMQXZPIJvAa5qmpqZg1axZef/11AIDVakVERAQef/xxPP3001fsv3TpUhiNRuzevbt325w5c5CcnIyNGzcCuHQGqK2tDbt27RpWJoPBAKVSCb1eDz8/v2Edg5yHqceCG9fuQ5hKjsybEoSOQ+TUNuwtxbmGdnz71CLIPSVCxyE3M5Sf34KeATKbzSgsLERaWlrvNrFYjLS0NOTl5fX7mry8vD77A0B6evoV++fm5iIkJAQJCQl45JFH0NzcPGAOk8kEg8HQ50Hu4628StTrL2LprEihoxA5vbtnhKOpw4Tt+VVCRyG6KkELUFNTEywWC9TqvgvOqdVqaLXafl+j1Wqvuf/ixYuxbds25OTk4MUXX8S+fftwyy23wGLpf1w6KysLSqWy9xERETHCT0bOQt/ZjVe/OY9FCSEI55VfRCMWqvTC/LggvL63lHOByKEJPgdoNNx333244447MHXqVCxZsgS7d+9GQUEBcnNz+91/5cqV0Ov1vY/qao5fu4sNuaUw91hx70xe+UVkL3dPH4cWoxlvH+IVYeS4BC1AQUFBkEgkaGho6LO9oaEBGo2m39doNJoh7Q8AMTExCAoKQmlpab/Py2Qy+Pn59XmQ6yvTdWDzgXLcNjUMKoVU6DhELkPtJ8cN8cHYsLeU6wKRwxK0AEmlUsycORM5OTm926xWK3JycjB37tx+XzN37tw++wPAnj17BtwfAGpqatDc3IzQUK7uS5fYbDb8bucJBHrLcEdSmNBxiFzOvTMj0Gm2YGNumdBRiPol+BBYZmYm3nzzTWzduhVnzpzBI488AqPRiIyMDADAsmXLsHLlyt79n3jiCWRnZ+OVV17B2bNn8eyzz+LIkSP49a9/DQDo6OjA//zP/+DQoUOoqKhATk4O7rzzTsTGxiI9PV2Qz0iOZ2dRLQ6XtyDjuije8Z1oFAR4S3HLlFD8Y3856touCh2H6AqCf+dfunQp1q5di1WrViE5ORnFxcXIzs7unehcVVWF+vr63v3nzZuHd999F//3f/+HpKQkfPjhh9i1axemTJkCAJBIJDh+/DjuuOMOxMfH46GHHsLMmTPx3XffQSaTCfIZybG0GM14/rPTmDchENPGqYSOQ+Sybk8KhdxTgle+KhE6CtEVBF8HyBFxHSDXZbPZ8Mu3CpFX1owX750Gf879IRpVe05rsflABT5+7DokR6iEjkMuzmnWASIaazsKqvHV6QYsvz6G5YdoDNyYqEZUkDd+//EJWKz8fZscBwsQuY0Lug48++kp3JgYglnRAULHIXILErEIGfOicKrOgHcP87J4chwsQOQW2ru6sXzbEQQopPjFnPFCxyFyK3FqXyxKCMFL2SXQtZuEjkMEgAWI3IDFasN/vXcUWn0XfnNzAu9PRCSA+2dHQCwW4fcfnwCnnpIjYAEil2az2ZD1+RnsO6fD4zfGIYy3uyAShK/cEw9dF42vTjfg46O1QschYgEi1/aXnPP4+/5yLJsbhSRegUIkqFnRAZgfG4TVn5xCvZ5rA5GwWIDIZf01txTrvz6P+2ZFIH3ywLdKIaKx88C8S4uPPvFeMXosVqHjkBtjASKXY7Xa8MJnp/FSdgnumRGOO5PDhY5ERN/zkXng14ticaSyBS9/yQUSSTgsQORSLpot+PV7Rfj7d+V4YG4U7p0ZIXQkIvqRxFA//L/Z4/G3by/gy1NaoeOQm2IBIpdxvqEdd7y+HzlnGvHfN8Vj8RQOexE5qlunajA7OgBPbi/GyVq90HHIDbEAkdOzWG3YlleB21/bj4vdFjx35xTMiuJCh0SOTCQS4ZEbJiDc3wsPbM5HdUun0JHIzbAAkVM7XWfA3W8cwKpPTmF+XBCeXzIFEQEKoWMR0SDIPSVYcXMCPMViLNuUj8b2LqEjkRthASKnVNt2ESs+KMZtr36H5g4zVt8+CQ/Nj4HMg4scEjkTpZcnfrs4EfqL3fi3v+VBq2cJorHBu8H3g3eDd1zVLZ3YuK8M7x+phkLqgbunh+PGxBB4SNjliZxZg6ELL3x2GnJPCd75zzmIDOSZXBq6ofz8ZgHqBwuQY7HZbCiqasWm/RXIPqmFj9wDiydrkD5ZAy8pz/gQuQpduwn/+8UZdHVbsPHfZ2JOTKDQkcjJsACNEAuQYzCaevDPY3V4+1AlTtUZEKqUI32yBgsTgjnUReSiOrp68Jecczijbcfq2yfhF3PGQyQSCR2LnAQL0AixAAnHZrPhWI0eOwqq8ElxHbq6LUiOUOEnE9VIjlBBzG+ERC6vx2rF24eq8OUpLRYlBOOle5MQ7CsTOhY5ARagEWIBGnvNHSbsKq7DjoIqnGvoQJCPFDfEB2NhQgiCfPiNj8gdFVW24v++uwCRCHh6cSL+LeXSHeWJBsICNEIsQGPDYrXh2/M67CioxtenGwAAM8f7Y2FCCKaFK/mNjoigv9iNdw5X4rvzTUiOUOF3t07E7Giu80X9YwEaIRag0VXd0okPjlTj/SM10Bq6EBnghRviQzA/Lgh+ck+h4xGRAzpTb8BbhypR3mTEgrgg/NdP4pDCBU/pR1iARogFyP7MPVZ8dVqL9/KrcKC0GV6eEsybEIhFiSGICfLmJEciuiarzYaC8hZ8WFSDmtaLmB6pwvLrY3DTJDU8uRQGgQVoxFiA7OeCrgM7CqrxQWENWoxmJKh9sSgxGKnRgZB78kouIho6q82Go1Vt+OxEHc7UtyPYV4b7ZkXgZzMjuH6Qm2MBGiEWoJEx9Vjw5akGvHu4CocuNMNH5oH5cUG4MSGEt6kgIruqbDbi6zONOFDahIvdFqRGB+DuGeFYPCUUSi8OqbsbFqARYgEantLGDmzPr8KHRTVo6+zGxFBf3JioxuyoAEg9eHqaiEaPqceC/PIWfHe+CSdr9fCUiLEoMRh3JF1aLZ6LproHFqARYgEavItmCz4/UY/38qtwpLIVvnIPXB8XjBsTQxCu8hI6HhG5oRajGQfLmnDoQjPKdEZ4eUrwk4khuG1qKBYmsAy5MhagEWIBujqbzYbi6ja8f6QGnx6rQ4epB1PDlViYEIxZUQGcjEhEDkOr78Kh8mbklzejvKkTck8xFiaE4JYpGixMCOEwmYthARohFqD+1bVdxK7iWnxYWIMLOiMCfaS4IS4YC+KDofaTCx2PiOiq6vUXkV/egoKKFpTpjPAQi5AaE4C0iWr8JFHNCdQugAVohFiA/qWt04zsk1rsKq7D4QvN8JSIMSvaHwvigjEljIsVEpFzau4w4UhlK4oqW3G63oAeqw1RgQosTAjBnJhAzIkJgEohFTomDREL0Ai5ewFqNZqx53QDPjtRjwOlTbDabJgcpsR1sUGYFeUPhdRD6IhERHbTae7BqVoDjtW04WSdHg0GE0QAJoT4YGakP5IiVJgU5ocEte+YzB+yWm2wfP+jWSwSQcJfNAeNBWiE3K0A2Ww2lOk6kFuiw1enGnCksgU2G5AY6ovZUfxNiIjci67dhNP1BpxvaMf5xg7UtHbC+v1PSo2fHOMDFQhXeUGtlCNAIYXSyxMKmQQyDwk8JCLAdmmtInOPFaYeK4zmHhhNPejo6oGhqwftXT0wdHXDcLEb7V09aO/qRqfZgovdFph7rOix9v2xLBYBUg8xFFIP+Mk9oFJIEeQjRbCvDBo/L4Sq5IjwVyAyUIFQP7lbn5l3ugK0YcMGvPzyy9BqtUhKSsJrr72G2bNnD7j/Bx98gGeeeQYVFRWIi4vDiy++iFtvvbX3eZvNhtWrV+PNN99EW1sbrrvuOrzxxhuIi4sbVB53KECNhi7kXWjGwdJmfFeqQ11bFzwlIkwNV2LGeH/MjPRn6SEiwqWV7KtbO1HT2gmtvgsNBhNaOs1oNZrR3tWDi92Wq75eBMBLKoGXpwQKmQQKTw8opBJ4SSVQSC/9t8xDDLmnBJ4SETwkYoi/Xx3fZrOhx3qpTHV1XypJlwuUvrMbLUYz2i52976XVCLG+EAFYkN8EBfigwkhPohX+yI6yNstFp8dys9vwccyduzYgczMTGzcuBGpqalYv3490tPTUVJSgpCQkCv2P3jwIO6//35kZWXhpz/9Kd59910sWbIERUVFmDJlCgDgpZdewquvvoqtW7ciOjoazzzzDNLT03H69GnI5e43Wfei2YIzWgOOV7fheI0eBRUtqG69CACI8PfCtHAVfjFHiYmhfpB5uP4/ECKioZB6iDEh2AcTgn36fd5itcHUY0G3xQaL9fLQFeAhEUMqEcNTIhrV2/2Ye6xo6jChwdAFraEL9fouVLd04vCFZrR0XipHEpEIEQFeSND4Il7ti9iQS58nJtjbbac1CH4GKDU1FbNmzcLrr78OALBarYiIiMDjjz+Op59++or9ly5dCqPRiN27d/dumzNnDpKTk7Fx40bYbDaEhYXhN7/5DVasWAEA0Ov1UKvV2LJlC+67775rZnLGM0BWqw26DhOqWzpR2dyJC00dKGvswFltOyqbO2ED4CERISpAgVi1LxLVvkjQ+PIsDxGRC+sw9aCu7eL3Z7Auorb1IuraLqLZaO7dR+MnR1SQN6KDFIgIUCDCX4Fwfy+EKb0Q7CtzqjlITnMGyGw2o7CwECtXruzdJhaLkZaWhry8vH5fk5eXh8zMzD7b0tPTsWvXLgBAeXk5tFot0tLSep9XKpVITU1FXl5evwXIZDLBZDL1/lmv1wO49Bc5mi5PdOux2NBttV763x4rzBYrzBYLusxWXOy2oLPbgosmC9q7utFu6kZbZw/0F81oMZqhazdD19GFpnZzn3HjAIUnNEovRKvkWBAVgogAL4zzV8BT/IM1eiwmGNtN/SQjIiJXIAIQrgDCFYpL//G9zu4e1LddRL3BBK3+IhoNHdh/uhm69i50mq29+0lEIgR4eyLQV4ZAb+mlOU8KKVRenvCVS+Aj84S3/PshPU8xZB4SyDwlkEnEkHmI4ekhhkQsgqdEDA/xpQndo3k27PLP7cGc2xG0ADU1NcFisUCtVvfZrlarcfbs2X5fo9Vq+91fq9X2Pn9520D7/FhWVhb++Mc/XrE9IiJicB/EAVUDOCZ0CCIicnoVQgcYhvb2diiVyqvu454Dfz+ycuXKPmeVrFYrWlpaEBgYOKpNdbQYDAZERESgurraaYbw3B2/Zs6HXzPnwq+X8xnO18xms6G9vR1hYWHX3FfQAhQUFASJRIKGhoY+2xsaGqDRaPp9jUajuer+l/+3oaEBoaGhffZJTk7u95gymQwymazPNpVKNZSP4pD8/Pz4D93J8GvmfPg1cy78ejmfoX7NrnXm5zJBb9oklUoxc+ZM5OTk9G6zWq3IycnB3Llz+33N3Llz++wPAHv27OndPzo6GhqNps8+BoMBhw8fHvCYRERE5F4EHwLLzMzEAw88gJSUFMyePRvr16+H0WhERkYGAGDZsmUIDw9HVlYWAOCJJ57ADTfcgFdeeQW33XYbtm/fjiNHjuD//u//AAAikQhPPvkknn/+ecTFxfVeBh8WFoYlS5YI9TGJiIjIgQhegJYuXQqdTodVq1ZBq9UiOTkZ2dnZvZOYq6qqIP7BlUvz5s3Du+++iz/84Q/43e9+h7i4OOzatat3DSAAeOqpp2A0GvHwww+jra0N8+fPR3Z2ttusASSTybB69eorhvXIcfFr5nz4NXMu/Ho5n9H+mgm+DhARERHRWBN0DhARERGREFiAiIiIyO2wABEREZHbYQEiIiIit8MC5EZMJhOSk5MhEolQXFwsdBzqR0VFBR566CFER0fDy8sLEyZMwOrVq2E2m6/9YhozGzZsQFRUFORyOVJTU5Gfny90JBpAVlYWZs2aBV9fX4SEhGDJkiUoKSkROhYN0po1a3qXt7E3FiA38tRTTw1qeXASztmzZ2G1WvG3v/0Np06dwp///Gds3LgRv/vd74SORt/bsWMHMjMzsXr1ahQVFSEpKQnp6elobGwUOhr1Y9++fXjsscdw6NAh7NmzB93d3bj55pthNBqFjkbXUFBQgL/97W+YNm3aqByfl8G7iS+++AKZmZn46KOPMHnyZBw9enTAW4OQY3n55Zfxxhtv4MKFC0JHIQCpqamYNWsWXn/9dQCXVq+PiIjA448/jqefflrgdHQtOp0OISEh2LdvHxYsWCB0HBpAR0cHZsyYgb/+9a94/vnnkZycjPXr19v1PXgGyA00NDRg+fLleOutt6BQKISOQ0Ok1+sREBAgdAwCYDabUVhYiLS0tN5tYrEYaWlpyMvLEzAZDZZerwcA/ptycI899hhuu+22Pv/W7E3wlaBpdNlsNjz44IP41a9+hZSUFFRUVAgdiYagtLQUr732GtauXSt0FALQ1NQEi8XSu1L9ZWq1GmfPnhUoFQ2W1WrFk08+ieuuu67P3QPIsWzfvh1FRUUoKCgY1ffhGSAn9fTTT0MkEl31cfbsWbz22mtob2/HypUrhY7s1gb79fqh2tpaLF68GD/72c+wfPlygZITuY7HHnsMJ0+exPbt24WOQgOorq7GE088gXfeeWfUb1/FOUBOSqfTobm5+ar7xMTE4N/+7d/w6aefQiQS9W63WCyQSCT4+c9/jq1bt452VMLgv15SqRQAUFdXh4ULF2LOnDnYsmVLn/vhkXDMZjMUCgU+/PDDPjdXfuCBB9DW1oZPPvlEuHB0Vb/+9a/xySef4Ntvv0V0dLTQcWgAu3btwl133QWJRNK7zWKxQCQSQSwWw2Qy9XluJFiAXFxVVRUMBkPvn+vq6pCeno4PP/wQqampGDdunIDpqD+1tbVYtGgRZs6cibfffttu/9jJPlJTUzF79my89tprAC4Nq0RGRuLXv/41J0E7IJvNhscffxwff/wxcnNzERcXJ3Qkuor29nZUVlb22ZaRkYHExET89re/tevQJecAubjIyMg+f/bx8QEATJgwgeXHAdXW1mLhwoUYP3481q5dC51O1/ucRqMRMBldlpmZiQceeAApKSmYPXs21q9fD6PRiIyMDKGjUT8ee+wxvPvuu/jkk0/g6+sLrVYLAFAqlfDy8hI4Hf2Yr6/vFSXH29sbgYGBdp+3xQJE5ED27NmD0tJSlJaWXlFQebLWMSxduhQ6nQ6rVq2CVqtFcnIysrOzr5gYTY7hjTfeAAAsXLiwz/bNmzfjwQcfHPtA5DA4BEZERERuhzMriYiIyO2wABEREZHbYQEiIiIit8MCRERERG6HBYiIiIjcDgsQERERuR0WICIiInI7LEBERETkdliAiIiIyO2wABEREZHbYQEiIiIit8MCRERERG7n/wOTL3dENUMrOwAAAABJRU5ErkJggg==",
      "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": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "================== DoubleMLDIDData Object ==================\n",
      "Time variable: t\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",
      "\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "from doubleml.did.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 = DoubleMLDIDData.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": 17,
   "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": 18,
   "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",
      "No. Observations: 1000\n",
      "\n",
      "\n",
      "------------------ Score & Algorithm ------------------\n",
      "Score function: observational\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",
      "Regression:\n",
      "Learner ml_g_d0_t0 RMSE: [[14.54119805]]\n",
      "Learner ml_g_d0_t1 RMSE: [[26.40359107]]\n",
      "Learner ml_g_d1_t0 RMSE: [[27.98750578]]\n",
      "Learner ml_g_d1_t1 RMSE: [[45.69520523]]\n",
      "Classification:\n",
      "Learner ml_m Log Loss: [[0.67618978]]\n",
      "\n",
      "------------------ Resampling        ------------------\n",
      "No. folds: 5\n",
      "No. repeated sample splits: 1\n",
      "\n",
      "------------------ Fit Summary       ------------------\n",
      "       coef   std err         t     P>|t|     2.5 %     97.5 %\n",
      "d  7.667682  5.809278  1.319903  0.186868 -3.718294  19.053659\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": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "      5.0 %     95.0 %\n",
      "d -1.887731  17.223095\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": 20,
   "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": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Coverage: 0.96\n",
      "Average CI length: 23.357586986897548\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": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGhCAYAAACJaguMAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjYsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvq6yFwwAAAAlwSFlzAAAPYQAAD2EBqD+naQAAXCdJREFUeJzt3Xlc1HXiP/DXzMDMcA73DCDIDSoEHoiYZRaJZSVWm/nYTXPdtnsrqk39llptP7JW10rLakttyzQ7rEwpIrVMPDg8UMGLG4b7vgZmPr8/zNlYUQGBzxyv5+Mxj10/vAde0yi8+Hzen/dbIgiCACIiIiIrIhU7ABEREdFwYwEiIiIiq8MCRERERFaHBYiIiIisDgsQERERWR0WICIiIrI6LEBERERkdViAiIiIyOqwABEREZHVYQEiIiIiq2MSBWjt2rUICAiAUqlEXFwcDh48eNnxW7duRUREBJRKJaKiorBjx44eH5dIJL0+Xn/99aF8GURERGQmRC9AW7ZsQXJyMpYtW4bs7GxER0cjMTERVVVVvY7ft28f5s6di4ULFyInJwdJSUlISkpCbm6ucUxFRUWPx4cffgiJRIK77rpruF4WERERmTCJ2JuhxsXFITY2FmvWrAEAGAwG+Pn54fHHH8eiRYsuGj9nzhy0trZi+/btxmOTJk1CTEwM1q1b1+vXSEpKQnNzM9LT0/uUyWAwoLy8HE5OTpBIJAN4VURERDTcBEFAc3MzfHx8IJVe/hyPzTBl6pVOp0NWVhYWL15sPCaVSpGQkICMjIxen5ORkYHk5OQexxITE7Ft27Zex1dWVuK7777Dxo0bL5mjs7MTnZ2dxj+XlZVh9OjR/XglREREZCpKSkowYsSIy44RtQDV1NRAr9dDrVb3OK5Wq5GXl9frc7Raba/jtVptr+M3btwIJycn3HnnnZfMkZKSghdffPGi4yUlJXB2dr7SyyAiIiIT0NTUBD8/Pzg5OV1xrKgFaDh8+OGH+OMf/wilUnnJMYsXL+5xVunCf0BnZ2cWICIiIjPTl+krohYgDw8PyGQyVFZW9jheWVkJjUbT63M0Gk2fx//yyy/Iz8/Hli1bLptDoVBAoVD0Mz0RERGZK1HvApPL5Rg/fnyPyckGgwHp6emIj4/v9Tnx8fEXTWZOS0vrdfwHH3yA8ePHIzo6enCDExERkVkT/RJYcnIy5s+fjwkTJmDixIlYvXo1WltbsWDBAgDAvHnz4Ovri5SUFADAE088galTp2LlypWYOXMmNm/ejMzMTLz33ns9Pm9TUxO2bt2KlStXDvtrIiIiItMmegGaM2cOqqursXTpUmi1WsTExCA1NdU40bm4uLjHrWyTJ0/Gpk2b8Pzzz2PJkiUIDQ3Ftm3bEBkZ2ePzbt68GYIgYO7cucP6eoiIyHzo9Xp0dXWJHYP6SCaTwcbGZlCWqBF9HSBT1NTUBJVKhcbGRk6CJiKyUC0tLSgtLQV/DJoXe3t7eHt7Qy6XX/Sx/vz8Fv0MEBER0XDT6/UoLS2Fvb09PD09ueitGRAEATqdDtXV1SgoKEBoaOgVFzu8HBYgIiKyOl1dXRAEAZ6enrCzsxM7DvWRnZ0dbG1tUVRUBJ1Od9klbq5E9L3AiIiIxMIzP+bnas769Pg8g/JZiIiIiMwIL4ERERH9pqyhHfWtumH7eq4Ocvi68BKcGFiAiIiIcL783LRyNzq6DMP2NZW2UqQ/fYOoJWjDhg148skn0dDQIFoGMbAAERERAahv1aGjy4BHp4UMSyEpa2jH2l1nUN+q69fXu//++7Fx48aLjicmJiI1NfWyzw0ICMCTTz6JJ5980nhszpw5uPXWW/v89QfK1IoWCxAREdHv+LrYIdDDQewYlzVjxgysX7++x7GB7mlpZ2dnlXfCcRI0ERGRmVEoFNBoND0erq6uEAQBy5cvh7+/PxQKBXx8fPC3v/0NAHDDDTegqKgITz31FCQSifEOuA0bNsDFxcX4uZcvX46YmBh8+OGH8Pf3h6OjIx555BHo9Xq89tpr0Gg08PLywiuvvNIj06pVqxAVFQUHBwf4+fnhkUceQUtLCwBg9+7dWLBgARobG41fe/ny5QCAzs5OPPPMM/D19YWDgwPi4uKwe/fuIf9vyDNARGT12nV6fHesAntPV+NcTSsqGjvgpLCBu6Mc14xwwXWhHpgU5A6lrUzsqESX9cUXX+Bf//oXNm/ejDFjxkCr1eLIkSMAgC+//BLR0dH461//igceeOCyn+fs2bPYuXMnUlNTcfbsWdx99904d+4cwsLCsGfPHuzbtw9//vOfkZCQgLi4OADnb09/8803ERgYiHPnzuGRRx7B3//+d7z99tuYPHkyVq9ejaVLlyI/Px8A4OjoCAB47LHHcOLECWzevBk+Pj746quvMGPGDBw7dgyhoaFD9t+KBYiIrFa7To/V6afwyf5itHR2I8TLEb4udgjxckSHTo+G9i5sO1yGD/YWQGVniz/G+WP+5AConQe++BrRYNi+fbuxQFywZMkSKJVKaDQaJCQkwNbWFv7+/pg4cSIAwM3NDTKZDE5OTtBoNJf9/AaDAR9++CGcnJwwevRoTJs2Dfn5+dixYwekUinCw8OxYsUK7Nq1y1iAfj+vKCAgAP/4xz/w0EMP4e2334ZcLodKpYJEIunxtYuLi7F+/XoUFxfDx8cHAPDMM88gNTUV69evx//7f/9vMP5z9YoFiIis0oFztXj286PQNnZgRqQGN0Z49VpsBEFAaX07duVXYf2vhfj3LwV4cGoQHr4hGPZyfgslcUybNg3vvPNOj2Nubm5obW3F6tWrERQUhBkzZuDWW2/F7bffDhub/v1dDQgIgJOTk/HParUaMpmsxyKEarUaVVVVxj//+OOPSElJQV5eHpqamtDd3Y2Ojg60tbXB3t6+169z7Ngx6PV6hIWF9Tje2dkJd3f3fmXuL/7rJSKr81lmCRZ9cRRhaie8emcUvC9zB45EIoGfmz3mxQfg7vEj8O2RCqzbcxZbM0uRclcUpoV7DWNyovMcHBwQEhJy0XE3Nzfk5+fjxx9/RFpaGh555BG8/vrr2LNnD2xtbfv8+f93rEQi6fWYwXB+yYDCwkLcdtttePjhh/HKK6/Azc0Ne/fuxcKFC6HT6S5ZgFpaWiCTyZCVlQWZrOcl5v89wzXYWICIyKp8sLcAL28/gZsivPDnawMhlfZ9KwR7uQ3mxPrhhnBPrP+1AAvWH8J9k0Ziya2jYCfn/CAyDXZ2drj99ttx++2349FHH0VERASOHTuGcePGQS6XQ6/XD/rXzMrKgsFgwMqVK41niT777LMeY3r72mPHjoVer0dVVRWuu+66Qc91OSxARGQ1Pt5fhJe3n8Ad0T64N9ZvwPtAqZ2VeG5GBNJOVOKTA8U4VFiH9+dNgJ9b77/lknkpa2g3+a/T2dkJrVbb45iNjQ22b98OvV6PuLg42Nvb4+OPP4adnR1GjhwJ4PylrZ9//hn33nsvFAoFPDw8ruo1XBASEoKuri689dZbuP322/Hrr79i3bp1PcYEBASgpaUF6enpiI6Ohr29PcLCwvDHP/4R8+bNw8qVKzF27FhUV1cjPT0d11xzDWbOnDko+XrDAkREViHjbC2WfXMciWM0mDvR/6o/n0QiwfQxGoRrnPCvH0/htrf24u0/jsO1IYPzA4WGn6uDHEpbKdbuOjNsX1NpK4Wrg7zfz0tNTYW3t3ePY+Hh4Xj11Vfx6quvIjk5GXq9HlFRUfj222+N82leeuklPPjggwgODkZnZycEQRiU1xEdHY1Vq1ZhxYoVWLx4Ma6//nqkpKRg3rx5xjGTJ0/GQw89hDlz5qC2thbLli3D8uXLsX79evzjH//A008/jbKyMnh4eGDSpEm47bbbBiXbpUiEwXr1FqSpqQkqlQqNjY1wdnYWOw4RXaWSujbcvmYv/Fzt8dyMCMj6cdmrL1o6uvHWT6dxvKIJ//zDNZg9dsSgfn4afB0dHSgoKEBgYCCUyv9OfudeYKbvUu8d0L+f3zwDREQWTW8Q8Nin2VDYSPG3G0MHvfwAgKPSBs/OCMcHewvw1JYjqGzqxENTgwf969DQ83WxYyGxEixARGTR3v/lHI6WNGL5HWPgqBy6b3k2Uin+el0QXO3leHVnHtp0ejyVEDrgeUZENLRYgIjIYp2pasaqH05h5jXeCFM7XfkJV0kikeCeCX5Q2EjxZvppGAwCnp4exhJEZIJYgIjIIhkMAp79/Cg8nBT4w3i/Yf3as2J8IZNKsGbXGchtpPjbTUO3nD8RDQwLEBFZpK+PlCGnuAEv3DYacpvh3/f5tmt8oOs2YFXaKbg5yPGnSSOHPQNdGe8DMj+D9Z6xABGRxWnTdePVnXmIC3TDaG/x7uScPdYXzR3deGFbLtwd5LglyvvKT6JhcWHVYZ1OBzs7Tno2J21tbQAuXq26v1iAiMjivLvnHOpadVh8yyhRc0gkEtwXPxJNHV14asth+LnZI9JXJWomOs/Gxgb29vaorq6Gra1tjz2uyDQJgoC2tjZUVVXBxcXloq0z+ovrAPWC6wARma/Kpg5MfW0Xpg/SgoeDQddtwEvbj6O1U49vHrsWXtxN3iTodDoUFBQY97Mi8+Di4gKNRtPrzQVcB4iIrNY7u89CJpNgVoyP2FGM5DZSJN8cjhe+zsWDH2fhswfjYSvjGQexyeVyhIaGQqcbvoUP6erY2tpe9ZmfC1iAiMhiaBs7sOlAMWbF+MBeblrf3twc5HjyplC8uP0EVv5wCotuiRA7EgGQSqUXrSZM1oG/ghCRxXhn9/nbzmdEasSO0qtQtRPumeCHdXvOYs+parHjEFk1FiAisggVje3YdLAYt0Z5m9zZn9+77RpvRI9Q4akth1HT0il2HCKrxQJERBbh378UQG4jReIYtdhRLksqkeDhG0LQpTdg+TfHxY5DZLVYgIjI7DW2d+HTg8W4eZTapM/+XKCys8X8+ABsP1qB749rxY5DZJVYgIjI7G06UIwuvQHTx5jm3J/eTA52x3h/V/zfV7lobOsSOw6R1WEBIiKz1tmtx4e/FmBKiAdc7eVix+kziUSCP08JRGtnN1al5Ysdh8jqsAARkVn75nA5qps7MTPKdNb96Ss3BznuHOeLj/cX41Rls9hxiKwKCxARmS1BEPDB3gKM83eBr6t57uc0Y4wGXs4KvPTtCW7MSTSMWICIyGwdKqxHnrYZiWY09+d/2cik+GPcSOw9U4O0E5VixyGyGixARGS2Nu4rhI9KafYbjI7zd0GUrwqvfZ8PvYFngYiGAwsQEZklbWMHUo9rcfNoNaS9bIpoTiQSCe6ZMAJnqlqw/Wi52HGIrAILEBGZpU0HiiCXSXF9mKfYUQZFiJcTxvm7YPWPp3kWiGgYsAARkdnp0hvwycFiTAn1MIuFD/vqrnEjUFDTim+OlIkdhcjisQARkdlJO1GJ2hYdbh5l2tte9FeQpyMmjHTFGz+ehoFngYiGFAsQEZmdTw4UIVztBD83e7GjDLpZMT4orG3Djyd5RxjRUBK9AK1duxYBAQFQKpWIi4vDwYMHLzt+69atiIiIgFKpRFRUFHbs2HHRmJMnT+KOO+6ASqWCg4MDYmNjUVxcPFQvgYiGUVFtK349U4tpEV5iRxkSIV5OCFc74d+/FIgdhciiiVqAtmzZguTkZCxbtgzZ2dmIjo5GYmIiqqqqeh2/b98+zJ07FwsXLkROTg6SkpKQlJSE3Nxc45izZ89iypQpiIiIwO7du3H06FG88MILUCqVw/WyiGgIfXqwBA4KGeKD3MWOMmRujfLGwcI6HC1tEDsKkcWSCCIuPRoXF4fY2FisWbMGAGAwGODn54fHH38cixYtumj8nDlz0Nraiu3btxuPTZo0CTExMVi3bh0A4N5774WtrS3+85//DDhXU1MTVCoVGhsb4ezsPODPQ0SDS9dtQHxKOmID3DB/coDYcYaMwSAgeethxAW64c2548SOQ2Q2+vPzW7QzQDqdDllZWUhISPhvGKkUCQkJyMjI6PU5GRkZPcYDQGJionG8wWDAd999h7CwMCQmJsLLywtxcXHYtm3bZbN0dnaiqampx4OITM+PJytR26rDjRZ6+esCqVSCGWM0+O6oFhWN7WLHIbJIohWgmpoa6PV6qNU97+JQq9XQarW9Pker1V52fFVVFVpaWvDqq69ixowZ+OGHHzB79mzceeed2LNnzyWzpKSkQKVSGR9+fn5X+eqIaChsOVSCULWjRU5+/l/Xh3nCRibB1sxSsaMQWSTRJ0EPJoPBAACYNWsWnnrqKcTExGDRokW47bbbjJfIerN48WI0NjYaHyUlJcMVmYj6qLyhHT+fqsYNYZZ99ucCe7kN4oPcsflgMRdGJBoCohUgDw8PyGQyVFb2vNWzsrISGk3vGxtqNJrLjvfw8ICNjQ1Gjx7dY8yoUaMuexeYQqGAs7NzjwcRmZbPs0qhsJVa9OTn/3XTKC+UN3bgl9PVYkchsjiiFSC5XI7x48cjPT3deMxgMCA9PR3x8fG9Pic+Pr7HeABIS0szjpfL5YiNjUV+fn6PMadOncLIkSMH+RUQ0XAxGARsOVSCSYHusJPLxI4zbII9HeHvZo/NB7mMB9FgE3UN+eTkZMyfPx8TJkzAxIkTsXr1arS2tmLBggUAgHnz5sHX1xcpKSkAgCeeeAJTp07FypUrMXPmTGzevBmZmZl47733jJ/z2WefxZw5c3D99ddj2rRpSE1Nxbfffovdu3eL8RKJaBDsO1uLsoZ2/PX6ILGjDCuJRIIbI7zwn/1FqGrugJcTl/MgGiyizgGaM2cO/vnPf2Lp0qWIiYnB4cOHkZqaapzoXFxcjIqKCuP4yZMnY9OmTXjvvfcQHR2Nzz//HNu2bUNkZKRxzOzZs7Fu3Tq89tpriIqKwr///W988cUXmDJlyrC/PiIaHJ9lFmOEqx1CvRzFjjLspoR4QCaR4Kts7g9GNJhEXQfIVHEdICLT0djWhdhXfsTd40fg9mgfseOI4o30U2ho60Lqk9eLHYXIpJnFOkBERH3xzdFydBsMuC7UQ+woork22AN52macqmwWOwqRxWABIiKT9tmhEoz1c4WLvVzsKKKJ9nOBg0KGrw/zMhjRYGEBIiKTladtwrGyRkwN8xQ7iqhsZVJMDHDHtpxycNYC0eBgASIik7U1sxQqO1uM9XcRO4ropoS4o6yhHdnFDWJHIbIILEBEZJJ03QZ8mV2Ka4PdYSPjt6oIb2e4O8h5GYxokPC7ChGZpJ/yqlDf1oUbwq1j64srkUokmBjohp25Whi4NQbRVWMBIiKT9FlmCYI9Haxi49O+ig1wQ3VzJw6XNogdhcjssQARkcmpau7Anvxqq5/8/L/C1U5Q2dni+1yt2FGIzB4LEBGZnK+yyyCVAvHB1rv2T2+kUgnG+btiZ66Wd4MRXSUWICIyKYIg4LPMEsQGuMFRIep2hSZpYqAriuvakM9FEYmuCgsQEZmU7OIGnK1u5eWvSxjjo4KdrQypvAxGdFVYgIjIpGzNLIGnoxyRviqxo5gkW5kUY/1dWICIrhILEBGZjDZdN745Uo7rwzwhlUjEjmOyJox0RZ62GWUN7WJHITJbLEBEZDK+O1qBdp0eU8O49s/lRI1wgVQC7M6vEjsKkdliASIik/FZZgkifVXwdFKIHcWkOSpsEK5xwq48FiCigWIBIiKTcK66BYcK6zn5uY+iR7hg75kadHbrxY5CZJZYgIjIJGzJLIGjwgaxAW5iRzELMX4u6Ogy4MC5OrGjEJklFiAiEp2u24DPM0sxJdQDcht+W+oLfzd7uDvKsYvzgIgGhN9piEh06ScrUduqw43c+LTPJBIJYka44CfOAyIaEBYgIhLdpweLEap25Man/RTj54Ki2jYU1LSKHYXI7LAAEZGoSura8MvpGkzj2Z9+G+Ojgkwqwd7T1WJHITI7LEBEJKqtmSVQ2soQH+QudhSzYyeXIdTLEXvP1IgdhcjssAARkWi69AZ8eqgEU0I9oLSViR3HLI3xUWHfmVp06w1iRyEyKyxARCSaH09Uorq5EzdF8PLXQF0zQoXmzm4cK2sUOwqRWWEBIiLRfHygCOFqJ4x0dxA7itkK8nSAna0Me0/zMhhRf7AAEZEoCmpa8euZWtw0imd/roaNVIrRPs6cB0TUTyxARCSKTQeK4KSwQVwgJz9frUgfFbKK6tHa2S12FCKzwQJERMOuXafHlswSXB/myZWfB0HUCBW6DQIOFnJbDKK+4nceIhp2Xx8uQ3N7N24erRY7ikXwUSnh7iDHr5wHRNRnLEBENKwEQcD6fYUY5+8KtbNS7DgWQSKRYJS3M/afqxU7CpHZYAEiomF1sKAO+dpmTB/Dsz+DabS3M05UNKGxvUvsKERmgQWIiIbVhn2F8HGxQ6SvSuwoFmW0jzMMApDJeUBEfcICRETDpqyhHT8cr8T00WpIJRKx41gULycF3B3lOFDAAkTUFyxARDRsNu4rhFIuxdQwT7GjWByJRIJRGmdknOU8IKK+YAEiomHR0tmNTQeKcVOEmvt+DZHR3s44Xt6Ipg7OAyK6EhYgIhoWnx0qQXuXHtN56/uQ4Twgor5jASKiIdetN+CDvQWID3KHu6NC7DgW68I8oP3nWICIroQFiIiGXOpxLcoa2nFrlLfYUSyaRCJBBOcBEfUJCxARDSlBELB21xlE+aoQ6MFd34dahMYJJ8qbuC8Y0RWwABHRkNpzqhonK5pxR7SP2FGsQoTGCXpBwOGSBrGjEJk0FiAiGlJv7zqLEC9HjPFxFjuKVfBxsYOT0gYHuR4Q0WWZRAFau3YtAgICoFQqERcXh4MHD152/NatWxEREQGlUomoqCjs2LGjx8fvv/9+SCSSHo8ZM2YM5Usgol5kFdXhYGEd7oj2gYQLHw4LqUSCMLUTDvFOMKLLEr0AbdmyBcnJyVi2bBmys7MRHR2NxMREVFVV9Tp+3759mDt3LhYuXIicnBwkJSUhKSkJubm5PcbNmDEDFRUVxsenn346HC+HiH7nzfQzGOFqh/EjXcWOYlXC1U7ILq5Hl94gdhQikyV6AVq1ahUeeOABLFiwAKNHj8a6detgb2+PDz/8sNfxb7zxBmbMmIFnn30Wo0aNwssvv4xx48ZhzZo1PcYpFApoNBrjw9WV34CJhtPhkgbsOVWN2WN9ue3FMIvQOKGjy4Dj5U1iRyEyWaIWIJ1Oh6ysLCQkJBiPSaVSJCQkICMjo9fnZGRk9BgPAImJiReN3717N7y8vBAeHo6HH34YtbWXvi20s7MTTU1NPR5EdHVW/3gKI1ztMCnQXewoVifQwwEKGykXRCS6DFELUE1NDfR6PdTqnivDqtVqaLXaXp+j1WqvOH7GjBn46KOPkJ6ejhUrVmDPnj245ZZboNfre/2cKSkpUKlUxoefn99VvjIi63a4pAG78387+yPl2Z/hZiOTIsTLkROhiS7DRuwAQ+Hee+81/v+oqChcc801CA4Oxu7du3HTTTddNH7x4sVITk42/rmpqYkliOgq/CuNZ3/EFq52wq78KgiCwAnoRL0Q9QyQh4cHZDIZKisrexyvrKyERqPp9TkajaZf4wEgKCgIHh4eOHPmTK8fVygUcHZ27vEgooE5WFCHPaeqcefYETz7I6IwtRPq27pQUNMqdhQikyRqAZLL5Rg/fjzS09ONxwwGA9LT0xEfH9/rc+Lj43uMB4C0tLRLjgeA0tJS1NbWwtuby/ATDSVBELAiNQ+BHg6IC3ITO45VC/FyhARAdnGD2FGITJLod4ElJyfj/fffx8aNG3Hy5Ek8/PDDaG1txYIFCwAA8+bNw+LFi43jn3jiCaSmpmLlypXIy8vD8uXLkZmZicceewwA0NLSgmeffRb79+9HYWEh0tPTMWvWLISEhCAxMVGU10hkLXblVyGrqB5zJvjxzi+ROShs4Odmj6yierGjEJkk0ecAzZkzB9XV1Vi6dCm0Wi1iYmKQmppqnOhcXFwMqfS/PW3y5MnYtGkTnn/+eSxZsgShoaHYtm0bIiMjAQAymQxHjx7Fxo0b0dDQAB8fH0yfPh0vv/wyFAruQk00VAwGAStS8zHa2xnXjFCJHYcAhHo5IquIE6GJeiMRBEEQO4SpaWpqgkqlQmNjI+cDEfXR51mleGbrEbx4xxiEqZ3EjkM4vw/bu3vO4vCy6VDZ2Yodh2jI9efnt+iXwIjI/LXpuvFaah4mBbmx/JiQMC9HCAA3RiXqBQsQEV21938uQH2bDnNj/cWOQr+jUSnhpLThPCCiXrAAEdFVqWrqwLo9Z5E4RgMvZ6XYceh3JBIJwrycOA+IqBcsQER0VVak5sFGJkFSjK/YUagXoWpH5BQ3QG/gdE+i32MBIqIByy6uxxfZZZgzwQ8OCtFvKqVehKmd0KbTI1/bLHYUIpPCAkREA2IwCFj29XEEejhgWriX2HHoEoI8HSCVADklnAdE9HssQEQ0IFuzSnCsrBHz4wO45YUJU9jIEODugByuCE3UAwsQEfVbfasOKTvzMCXEA+Ea3vZu6oK9HJHNO8GIemABIqJ+e+37PHR1G/DHON72bg5CvRxxrqYVDW06saMQmQwWICLql+zienx6sAT3TPCDi71c7DjUByFejgC4ICLR77EAEVGfdesN+L+vjiHIwwEJo9Rix6E+0jifXxCR84CI/osFiIj6bMO+QuRVNOPPUwI58dmMSCQShHg5IqeY84CILmABIqI+Ka1vwz9/yMf0MRoEezqKHYf6KcTTETklDTBwQUQiACxARNQHgiBg6dfH4SC3wT0TRogdhwYgVO2E5o5unKtpFTsKkUlgASKiK9qZq8VPeVWYFx8AezlXfDZHwZ4OkOD8JHYiYgEioitobO/C0q+PY8JIV8QGuIodhwbIXm6DEa52nAhN9BsWICK6rBWpeWjt7Mb9kwMgkXDiszkL9uREaKILWICI6JIOFdZh04Fi3BvrB3dHhdhx6CqFeDniVGUz2nTdYkchEh0LEBH1qrNbj0VfHEWo2pFr/liIEC9HGATgWGmj2FGIRMcCRES9Wrf7HApr2/CXKUFc88dCjHC1h8JGyhWhicACRES9OFPVgjW7TuP2a7zh72YvdhwaJDKpBEGeDshhASJiASKingRBwP99dQzujgrMHss1fyxNsKcjDvNOMCIWICLq6fOsUhwoqMOfrw2E3IbfIixNiJcjtE0dqGzqEDsKkaj43Y2IjOpbdXhlx0lMCfFAlK9K7Dg0BEJ+28aE6wGRtWMBIiKjlJ0n0aU34I9x/mJHoSHi7qiAm4OcE6HJ6rEAEREAILOwDp9lluLeWH+42MvFjkNDKNjTgQsiktVjASIidOsNeH5bLoI9HXBjuJfYcWiIhXg64lhZI/TcGZ6sGAsQEeHj/UXI1zZjwbWBXPPHCgR7OaJNp8eZqhaxoxCJhgWIyMpVN3finz+cwo0RXgj+bYIsWbYgD0dIABzhPCCyYixARFbun9/nAwDmxPqJnISGi51chhGudjhc2iB2FCLRsAARWbHcskZ8llmCu8ePgJPSVuw4NIy4MzxZOxYgIislCAKWf3McI1ztuNmpFQr2csQpbQvadXqxoxCJggWIyEp9d6wCmUX1uC8+ADJOfLY6IV6O0AsCcsu5MzxZJxYgIivU2a3HqzvzMM7fhSs+Wym/33aG50RoslYsQERW6JP9xShvaMfciVzx2VrJpBIEenBneLJeLEBEVqapowtv/nQaU8O8MMLVXuw4JKIg7gxPVowFiMjKrNt9Fu06Pe4eP0LsKCSyEE9HlDW0o7q5U+woRMOOBYjIilQ1d+DDvQW4JVIDNwfu92XtQrwcAHBBRLJOLEBEVuSd3Wchk0ow8xofsaOQCfBwVMDFzhZHuCAiWSEWICIrUd7Qjo/3F+HWKG84KmzEjkMmQCKRIMjTATmcB0RWiAWIyEqs2XUGSlsZbon0FjsKmZBgT0ccKW2AgTvDk5VhASKyAiV1bdhyqAS3X+MDO7lM7DhkQkK8HNHc0Y2C2laxoxANK5MoQGvXrkVAQACUSiXi4uJw8ODBy47funUrIiIioFQqERUVhR07dlxy7EMPPQSJRILVq1cPcmoi8/Huz2fhIJfh5tHc8oJ6CvJ0BMCJ0GR9RC9AW7ZsQXJyMpYtW4bs7GxER0cjMTERVVVVvY7ft28f5s6di4ULFyInJwdJSUlISkpCbm7uRWO/+uor7N+/Hz4+nPBJ1quyqQNbDpXglkhvKG159od6clTYwMdFicMsQGRlRC9Aq1atwgMPPIAFCxZg9OjRWLduHezt7fHhhx/2Ov6NN97AjBkz8Oyzz2LUqFF4+eWXMW7cOKxZs6bHuLKyMjz++OP45JNPYGvLXa7Jer338znIZVJMH8OzP9S7YA8uiEjWR9QCpNPpkJWVhYSEBOMxqVSKhIQEZGRk9PqcjIyMHuMBIDExscd4g8GA++67D88++yzGjBlzxRydnZ1oamrq8SCyBLUtnfjkQBESx2hgL+edX9S7EC9HnKhoQkcXd4Yn6yFqAaqpqYFer4da3fM3U7VaDa1W2+tztFrtFcevWLECNjY2+Nvf/tanHCkpKVCpVMaHn59fP18JkWnamFEEQQBmRGrEjkImLNjLEd0GAScq+MsfWQ/RL4ENtqysLLzxxhvYsGEDJBJJn56zePFiNDY2Gh8lJSVDnJJo6LXr9PgooxA3hHvBScnLwHRpI93sYSuT8DIYWZUBFaBz584Nyhf38PCATCZDZWVlj+OVlZXQaHr/jVWj0Vx2/C+//IKqqir4+/vDxsYGNjY2KCoqwtNPP42AgIBeP6dCoYCzs3OPB5G5+zyrBE3tXbiVZ3/oCmxkUgR4OHAiNFmVARWgkJAQTJs2DR9//DE6OjoG/MXlcjnGjx+P9PR04zGDwYD09HTEx8f3+pz4+Pge4wEgLS3NOP6+++7D0aNHcfjwYePDx8cHzz77LL7//vsBZyUyJ3qDgPd/KcDEQDd4OSvFjkNmINjTETnF9WLHIBo2AypA2dnZuOaaa5CcnAyNRoMHH3zwimv3XEpycjLef/99bNy4ESdPnsTDDz+M1tZWLFiwAAAwb948LF682Dj+iSeeQGpqKlauXIm8vDwsX74cmZmZeOyxxwAA7u7uiIyM7PGwtbWFRqNBeHj4gDISmZsfjmtRXNeG27jnF/VRiKcjSurbUdvCneHJOgyoAMXExOCNN95AeXk5PvzwQ1RUVGDKlCmIjIzEqlWrUF1d3efPNWfOHPzzn//E0qVLERMTg8OHDyM1NdU40bm4uBgVFRXG8ZMnT8amTZvw3nvvITo6Gp9//jm2bduGyMjIgbwUIov0wd4CjNI4Ifi3Re6IriTE67cFEbkxKlkJiSAIV70BTGdnJ95++20sXrwYOp0Ocrkc99xzD1asWAFvb/Pbd6ipqQkqlQqNjY2cD0RmJ7esEbe9tRdPJoQiLtBd7DhkJgRBwEMfZ2H+5AA8PZ1ny8k89efn91XdBZaZmYlHHnkE3t7eWLVqFZ555hmcPXsWaWlpKC8vx6xZs67m0xPRAGzcVwgPRzkmjHQTOwqZEYlEghAvR+4MT1ZjQCujrVq1CuvXr0d+fj5uvfVWfPTRR7j11lshlZ7vU4GBgdiwYcMl77oioqFR29KJrw+X485xvpBJ+7YMBNEFwZ6OSM3VwmAQIOXfH7JwAypA77zzDv785z/j/vvvv+QlLi8vL3zwwQdXFY6I+mfzofNrWE2L8BI5CZmjEC9HNHd241xNq3FOEJGlGlABSktLg7+/v/GMzwWCIKCkpAT+/v6Qy+WYP3/+oIQkoivTGwR8vL8Ik4Pd4cyFD2kALkyaP1zSwAJEFm9Ac4CCg4NRU1Nz0fG6ujoEBgZedSgi6r/d+VWoaOzAzaO56SkNjIPCBr4udlwPiKzCgArQpW4ca2lpgVLJRdeIxPDx/iIEeTogiLe+01XgRGiyFv26BJacnAzg/N0CS5cuhb29vfFjer0eBw4cQExMzKAGJKIrK61vw+78avzluiCxo5CZC/VyxN7TNWjTdcNePqBZEkRmoV9/u3NycgCcPwN07NgxyOVy48fkcjmio6PxzDPPDG5CIrqiLYdKoLSVYXIw1/2hqxPi5Qi9IOBoaSMmBfHvE1mufhWgXbt2AQAWLFiAN954g4sEEpmAbr0Bmw+VYEqoB5S2MrHjkJnzc7WHna0MOcUNLEBk0QZ0fnP9+vWDnYOIBmhXfjWqmztxE299p0EglUoQ7OmAbE6EJgvX5wJ05513YsOGDXB2dsadd9552bFffvnlVQcjor757FAJgjwdMNLdQewoZCFCvBzxy5kaCIIAiYQLIpJl6nMBUqlUxn8IKpVqyAIRUd9VNXfgp7wqzJs8UuwoZEFCvZyw7XA5Suvb4edmf+UnEJmhPheg31/24iUwItPwVXYZpFJgcrCH2FHIglxYBDGnpIEFiCzWgNYBam9vR1tbm/HPRUVFWL16NX744YdBC0ZElycIArZkliA2wA2OCt6uTIPH2c4WGmclsos4D4gs14AK0KxZs/DRRx8BABoaGjBx4kSsXLkSs2bNwjvvvDOoAYmod9nFDThX3Ypp4Zz8TIMvxMuRBYgs2oAKUHZ2Nq677joAwOeffw6NRoOioiJ89NFHePPNNwc1IBH17ovsUng4yjHah8tR0OALVTviREUTOrr0YkchGhIDKkBtbW1wcnICAPzwww+48847IZVKMWnSJBQVFQ1qQCK6WEeXHt8eKceUEA9IeZcODYEwtRO6DQKOlTWKHYVoSAyoAIWEhGDbtm0oKSnB999/j+nTpwMAqqqquDgi0TD4Ka8KzR3dmBLqKXYUslB+rvZQ2kqRxctgZKEGVICWLl2KZ555BgEBAYiLi0N8fDyA82eDxo4dO6gBiehiX2aXIsTLEb4udmJHIQslk0oQ4unIAkQWa0C3jtx9992YMmUKKioqEB0dbTx+0003Yfbs2YMWjoguVtvSid351fjTJK79Q0MrVO2EPaequSAiWaQB3zur0Wig0Wh6HJs4ceJVByKiy9t+tAICgHhufEpDLEztiK9yylBU24YAD640TpZlQAWotbUVr776KtLT01FVVQWDwdDj4+fOnRuUcER0sS9zSjHWzwXOSluxo5CFC/E6f7NLVlE9CxBZnAEVoL/85S/Ys2cP7rvvPnh7e/PUKNEwKaptxZGSRvztxlCxo5AVcFTYwM/VDlnF9bhr/Aix4xANqgEVoJ07d+K7777DtddeO9h5iOgyvjlcDjtbGcaNdBE7ClmJEC9OhCbLNKC7wFxdXeHm5jbYWYjoMgRBwFeHyzAhwBUKG5nYcchKhKqdcErbjKaOLrGjEA2qARWgl19+GUuXLu2xHxgRDa0TFU04V93KjU9pWIWrnSAAyCluEDsK0aAa0CWwlStX4uzZs1Cr1QgICICtbc/JmNnZ2YMSjoj+6+vD5VDZ2SLSl4uN0vDxVinhrLRBVmEdpoZx4U2yHAMqQElJSYMcg4gux2AQ8M3hckwMdIONdEAnbokGRCKRIFzjhIOFnAdElmVABWjZsmWDnYOILiOruB7apg48GBwkdhSyQmFqJ3yRVYouvQG2MhZwsgwD/pvc0NCAf//731i8eDHq6uoAnL/0VVZWNmjhiOi87UfK4e4gR5jaSewoZIXC1U7o6DbgeHmT2FGIBs2AzgAdPXoUCQkJUKlUKCwsxAMPPAA3Nzd8+eWXKC4uxkcffTTYOYmslt4gYPuxCsQFunHndxJFoIcD5DIpMgvrEOPnInYcokExoDNAycnJuP/++3H69GkolUrj8VtvvRU///zzoIUjIuDAuVrUtui49QWJxkYmRYiXIw4V1okdhWjQDKgAHTp0CA8++OBFx319faHVaq86FBH917dHK+DlpECwp6PYUciKhamdcKiwHoIgiB2FaFAMqAApFAo0NV18LfjUqVPw9ORtkkSDpUtvwM7cCkwKcueWMySqcI0T6lp1KKzl+m9kGQZUgO644w689NJL6Oo6vzKoRCJBcXExnnvuOdx1112DGpDImmWcrUVDWxcmBfHyF4krTO0ICYBDBbwMRpZhQAVo5cqVaGlpgaenJ9rb2zF16lSEhITAyckJr7zyymBnJLJaO45VQOOsRIC7vdhRyMrZy20Q4GGPAyxAZCEGdBeYSqVCWloafv31Vxw5cgQtLS0YN24cEhISBjsfkdXq0huQelyL60M9efmLTEK4xhkHCmrFjkE0KPpdgAwGAzZs2IAvv/wShYWFkEgkCAwMhEajgSAI/EZNNEj2n+PlLzItozTOSM3VoryhHT4udmLHIboq/boEJggC7rjjDvzlL39BWVkZoqKiMGbMGBQVFeH+++/H7NmzhyonkdXZcawCamcFL3+RyYjQnF+I8yAvg5EF6NcZoA0bNuDnn39Geno6pk2b1uNjP/30E5KSkvDRRx9h3rx5gxqSyNp06w1IzdXiOl7+IhPibGeLEa52OFBQi6SxvmLHIboq/ToD9Omnn2LJkiUXlR8AuPHGG7Fo0SJ88skngxaOyFrtP1eH+rYuxAW6iR2FqIcIjTP2n+MZIDJ//SpAR48exYwZMy758VtuuQVHjhzpd4i1a9ciICAASqUScXFxOHjw4GXHb926FREREVAqlYiKisKOHTt6fHz58uWIiIiAg4MDXF1dkZCQgAMHDvQ7F5FYduSeX/ww0MNB7ChEPYzydkJBTSuqmjvEjkJ0VfpVgOrq6qBWqy/5cbVajfr6+n4F2LJlC5KTk7Fs2TJkZ2cjOjoaiYmJqKqq6nX8vn37MHfuXCxcuBA5OTlISkpCUlIScnNzjWPCwsKwZs0aHDt2DHv37kVAQACmT5+O6urqfmUjEoPeIOD7XC0mBrrx8heZnAiNMwDgUEH/vtcTmRqJ0I91zWUyGbRa7SVXe66srISPjw/0en2fA8TFxSE2NhZr1qwBcP4uMz8/Pzz++ONYtGjRRePnzJmD1tZWbN++3Xhs0qRJiImJwbp163r9Gk1NTVCpVPjxxx9x0003XTHThfGNjY1wdnbu82shGgz7z9Xi3vf246U7xiCUu7+TCUr+7DBuHq3GS7MixY5C1EN/fn73axK0IAi4//77oVAoev14Z2dnfz4ddDodsrKysHjxYuMxqVSKhIQEZGRk9PqcjIwMJCcn9ziWmJiIbdu2XfJrvPfee1CpVIiOjr5k7t9n722bD6LhsvNYBdwd5Qj24t5fZJoiNE7IOMv1gMi89asAzZ8//4pj+nMHWE1NDfR6/UWX1dRqNfLy8np9jlar7XX8/27Cun37dtx7771oa2uDt7c30tLS4OHh0evnTElJwYsvvtjn3ERDxWAQsCNXi9iRbpDy8heZqDE+KuzKr0Z1cyc8nXr/hZjI1PWrAK1fv36ocgy6adOm4fDhw6ipqcH777+Pe+65BwcOHICXl9dFYxcvXtzjrFJTUxP8/PyGMy4RACCnpB7VzZ28+4tM2mif85cW9p+rxe3RPiKnIRqYAe0FNlg8PDwgk8lQWVnZ43hlZSU0Gk2vz9FoNH0a7+DggJCQEEyaNAkffPABbGxs8MEHH/T6ORUKBZydnXs8iMSw45gWrva2COPcHzJhrvZy+LrYIeMcL4OR+RK1AMnlcowfPx7p6enGYwaDAenp6YiPj+/1OfHx8T3GA0BaWtolx//+8/Z3jhLRcBIEATtzKzB+pCukUl7+ItM2ytsZ+87UiB2DaMBELUAAkJycjPfffx8bN27EyZMn8fDDD6O1tRULFiwAcH5O0e8nST/xxBNITU3FypUrkZeXh+XLlyMzMxOPPfYYAKC1tRVLlizB/v37UVRUhKysLPz5z39GWVkZ/vCHP4jyGon6IresCeUNHZgYyL2/yPSN8XFGYW0btI1cD4jM04B2gx9Mc+bMQXV1NZYuXQqtVouYmBikpqYaJzoXFxdDKv1vT5s8eTI2bdqE559/HkuWLEFoaCi2bduGyMjzt2PKZDLk5eVh48aNqKmpgbu7O2JjY/HLL79gzJgxorxGor7YmVsBJ4UNRnnz8heZvtHe56cKZJyrweyxI0ROQ9R//VoHyFpwHSAaboIgYNrK3Rjp5oCHpgaLHYeoT5774igmBbnhtbt7X2KEaLj15+e36JfAiAg4XdWCwpo2TAzg3V9kPkZ7O+PXM5wITeaJBYjIBOw8poWdrQyRviqxoxD12RhfZ5Q1tKOotlXsKET9xgJEZAJ25lZgrL8L5Db8J0nmY7S3M6QS4JfTvBuMzA+/2xKJrKi2FXnaZsTy8heZGXu5DUK9nLCXt8OTGWIBIhJZaq4WcpkUMX4uYkch6rdI3/PrAekNvJ+GzAsLEJHIduZqEe2ngtJWJnYUon6L8nVBU0c3jpU1ih2FqF9YgIhEpG3swOGSBl7+IrMV7OUAO1sZ9p6uFjsKUb+wABGJ6PvjWsikEoz1dxU7CtGA2EilGO3tzInQZHZYgIhEtDNXi0gfZzgqRF+UnWjAIn1VyCqqR5uuW+woRH3GAkQkkrpWHQ4W1PLyF5m9a0ao0G0QsJ+7w5MZYQEiEknaCS0EARg/kpe/yLx5q5TwclJgTz7nAZH5YAEiEsnOXC0ivJ3gYi8XOwrRVZFIJLhmhAq7WYDIjLAAEYmgqaMLe0/X8PIXWYxoPxcU1bWhsIbbYpB5YAEiEsGuvCp0GwRufkoWY4y3CjZSCfac4lkgMg8sQEQi2HmsAiGeDnB3VIgdhWhQ2MlliNA4YXd+ldhRiPqEBYhomLXr9NidX40JPPtDFuaaES7IOFuLji692FGIrogFiGiY7TlVhY5uAy9/kcWJ8XNBR7cBhwrrxI5CdEUsQETDbOcxLfzd7ODtYid2FKJBNcLVDu6OcvyUx8tgZPpYgIiGUWe3Hul5VYgNcBc7CtGgk0gkGOvngh9PVEIQuDs8mTYWIKJh9OuZGrR0diMukJe/yDKN83dFSX07zla3iB2F6LJYgIiG0c5jWviolBjhystfZJnG+KigsJHix5O8DEamjQWIaJh06Q34/rgWEwPdIJFIxI5DNCTkNlJE+qrw44lKsaMQXRYLENEw2X+uFk0d3ZgYyPk/ZNnG+bsiu7ge9a06saMQXRILENEw2ZmrhZeTAgHu9mJHIRpSY/1dYBCA3ad4GYxMFwsQ0TDQGwR8n6tFbAAvf5Hlc7WXI9jTAWm8DEYmjAWIaBgcLKhDbasOk4J49xdZh3H+rtidX81VoclksQARDYOduRXwcJQj2NNR7ChEwyI2wA1tOj32na0ROwpRr1iAiIaYwSBgx7EKTOTlL7IiI1zt4K1S4vtcXgYj08QCRDTEsorrUdOiQ1wQ7/4i6yGRSDBhpCvSTlZCb+Cq0GR6WICIhth3Ryvg7iBHiBcvf5F1iQ1wQ12rDpncHJVMEAsQ0RAyGATszK1AbIAbpLz8RVYm2MsRrva2+P44L4OR6WEBIhpCOSX1qGzq5N5fZJWkEgnGj3RD6vEKbo5KJocFiGgIbT9aATcHOcI0TmJHIRJFXKAbyhs6cKysUewoRD2wABENEYNBwHdHz9/9xctfZK1GeTtDZWeL745WiB2FqAcWIKIhklVcj6rmTkzi3V9kxWRSCWIDXPHt0XJeBiOTwgJENEQu3P0VqubdX2TdJgW5o7yhA0dKeRmMTAcLENEQ0BsEfHesAhMDefmLaJTmwmWwcrGjEBmxABENgUOFdajm5S8iAIBUKsHEQDdsP8q7wch0sAARDYHtR8vh6ajg4odEv5kU5I6Kxg5kF9eLHYUIAAsQ0aDr1hvw3dEKxAXx8hfRBRFqJ7g5yPH1YV4GI9PAAkQ0yPadrUV9WxcmB3uIHYXIZEilEsQHuePbI+Xo0hvEjkPEAkQ02L49Ug5vlRIB7vZiRyEyKVNCPVDf1oW9p2vEjkJkGgVo7dq1CAgIgFKpRFxcHA4ePHjZ8Vu3bkVERASUSiWioqKwY8cO48e6urrw3HPPISoqCg4ODvDx8cG8efNQXs7TrjT0Orv1SM3VIj7IHRJe/iLqYaSbPUa42mFbTpnYUYjEL0BbtmxBcnIyli1bhuzsbERHRyMxMRFVVVW9jt+3bx/mzp2LhQsXIicnB0lJSUhKSkJubi4AoK2tDdnZ2XjhhReQnZ2NL7/8Evn5+bjjjjuG82WRlfr5VA2aO7t59xdRLyQSCa4N9sD3J7Ro7ewWOw5ZOYkg8j2JcXFxiI2NxZo1awAABoMBfn5+ePzxx7Fo0aKLxs+ZMwetra3Yvn278dikSZMQExODdevW9fo1Dh06hIkTJ6KoqAj+/v5XzNTU1ASVSoXGxkY4OzsP8JWRNXp8UzaOlDZgxV3RYkchMknVzR342+bD+NecaMweO0LsOGRh+vPzW9QzQDqdDllZWUhISDAek0qlSEhIQEZGRq/PycjI6DEeABITEy85HgAaGxshkUjg4uLS68c7OzvR1NTU40HUX62d3Ug7UcnJz0SX4emkxChvJ3yRxctgJC5RC1BNTQ30ej3UanWP42q1GlqtttfnaLXafo3v6OjAc889h7lz516yDaakpEClUhkffn5+A3g1ZO3STlSio9uAycG8/EV0OdeFeuLXMzUob2gXOwpZMdHnAA2lrq4u3HPPPRAEAe+8884lxy1evBiNjY3GR0lJyTCmJEux7XAZwjVO8HRSih2FyKRNCnSH3EaKL7NLxY5CVkzUAuTh4QGZTIbKysoexysrK6HRaHp9jkaj6dP4C+WnqKgIaWlpl70WqFAo4Ozs3ONB1B+1LZ345VQNz/4Q9YGdXIaJgW7YmlnKrTFINKIWILlcjvHjxyM9Pd14zGAwID09HfHx8b0+Jz4+vsd4AEhLS+sx/kL5OX36NH788Ue4u/OHEg2tHblaCBAwKZB/14j64oYwTxTVtSGziFtjkDhsxA6QnJyM+fPnY8KECZg4cSJWr16N1tZWLFiwAAAwb948+Pr6IiUlBQDwxBNPYOrUqVi5ciVmzpyJzZs3IzMzE++99x6A8+Xn7rvvRnZ2NrZv3w69Xm+cH+Tm5ga5XC7OCyWL9nVOGa4Z4QJnO1uxoxCZhQhvZ3g5KbA1swSxAW5ixyErJHoBmjNnDqqrq7F06VJotVrExMQgNTXVONG5uLgYUul/T1RNnjwZmzZtwvPPP48lS5YgNDQU27ZtQ2RkJACgrKwM33zzDQAgJiamx9fatWsXbrjhhmF5XWQ9Sn77LfaxaSFiRyEyG1KJBNeHeeLbIxV44bbRcFLylwcaXqKvA2SKuA4Q9ceb6afxzu6zeOdP46CwkYkdh8hs1LZ04m+bc/DSrEj8adJIseOQBTCbdYCIzJ0gCPgiuxQTAlxZfoj6yd1RgXH+rvjkQBEnQ9OwYwEiugqHSxpQVNuG60I9xY5CZJamRXjhZEUzjpU1ih2FrAwLENFV+CqnDG4Ocozx5qVSooGIGeECD0c5Nh0oFjsKWRkWIKIB0nUb8M2RckwOdodUyp3fiQZCKpXghnAvfH24HE0dXWLHISvCAkQ0QLvyq9DQ1sXLX0RXaVq4F3R6A77I4srQNHxYgIgG6IusUgR5OMDfzV7sKERmzc1BjokBbtiYUQiDgZOhaXiwABENQG1LJ37Kq+LZH6JBMn2MGoU1bfjlTI3YUchKsAARDcDXh8sBAJNDuPUF0WAIVzshwN0eG34tEDsKWQkWIKIB+DyrFOP8XeHM1WuJBoVEIsH00Rrszq9GUW2r2HHICrAAEfXTyYomnKhownVhHmJHIbIo14Z4wFFpg/W/FoodhawACxBRP205VAKVnS1i/FzEjkJkUeQ2Utw8So0th0rQ2MZb4mlosQAR9UNntx5f5ZThulAP2Ej5z4dosN08Wo1ugwGfHCwSOwpZOH4HJ+qHH45XorG9CzeEe4kdhcgiudjLMSXEAxt+LYSu2yB2HLJgLEBE/bDlUAnCNU7wdbETOwqRxbo1yhtVzZ349ki52FHIgrEAEfVRaX0bfj1TgxvCuPYP0VAa4WqPsX4uWLfnLBdGpCHDAkTUR1szS6G0lWFSENf+IRpqs2J8cbqqBT+erBQ7ClkoFiCiPujWG7DlUAnig92htJWJHYfI4oVrnDDK2wlrd52BIPAsEA0+FiCiPtidXw1tUwduiuDkZ6LhMivaF0dKG7HvbK3YUcgCsQAR9cEnB4oQ5OmAIE9HsaMQWY1rRqgQ5OmAN9NPix2FLBALENEVlDW0Y8+patzIW9+JhpVEIkFSjC8OFNRh/zmeBaLBxQJEdAVbDpVAYSPD5GBufUE03CaMdEWghwNW/XCKc4FoULEAEV3G+cnPxZgc7A47OSc/Ew03iUSCu8aNwMHCOmRwLhANIhYgostIO1GJyqZO3DxaLXYUIqs1zt8FwZ4OWJXGs0A0eFiAiC7jo4wihGucMNLdQewoRFbrwlmgzKJ67D5VLXYcshAsQESXcKaqGRnnanHzKJ79IRJbjJ8LRmmcsGJnHleHpkHBAkR0CR/vL4bKzhYTA93EjkJk9SQSCe6d6I88bTO+4R5hNAhYgIh60drZjc+zSjEt3BO2Mv4zITIFYWonxAa44vXv89HZrRc7Dpk5fmcn6sUX2aVo03XjJl7+IjIpcyb4o6KxHf/JKBI7Cpk5FiCi/2EwCPjw1wJMDHSDh6NC7DhE9Du+rna4McILb/x4GnWtOrHjkBljASL6H3tOV6Owpg0zxniLHYWIevGH8X7QCwL+lXZK7ChkxliAiP7Hhl8LEeTpgDA19/0iMkXOdraYPdYXnxwoQr62Wew4ZKZYgIh+52x1C/acqkbiaA0kEonYcYjoEmaM0UDtrMTyb45zcUQaEBYgot/5YG8BXOxsER/sLnYUIroMG5kU8+JHIuNcLb49WiF2HDJDLEBEv6lt6cQXWaWYPkbDW9+JzECMnysmBrrhpW9PoKmjS+w4ZGb4XZ7oN//Zf/622oRRXiInIaK+mjdpJFo6u7DqB06Ipv5hASIC0NGlx8Z9hZga5gknpa3YcYioj9wdFbhr3Ah8lFGIwyUNYschM8ICRATgy+wyNLR14dYo3vpOZG5uifRGgIcDnt16BLpug9hxyEywAJHV0xsEvPvzWcQGukHtrBQ7DhH1k0wqwV+vC8K5mla8vfuM2HHITLAAkdXbmVuBoto23BHtI3YUIhqgke4OuCPaB2t+OoOTFU1ixyEzwAJEVk0QBKzddQZRvioEe3LhQyJzNnusL3xc7PDk5sPcLJWuiAWIrNqeU9U4WdHMsz9EFsBWJsUjNwTjXE0L7wqjK2IBIqv29q6zCPFyxBgfZ7GjENEgGOnugLvH++G9n89h/7laseOQCRO9AK1duxYBAQFQKpWIi4vDwYMHLzt+69atiIiIgFKpRFRUFHbs2NHj419++SWmT58Od3d3SCQSHD58eAjTkznbf64WBwvrMCvah9teEFmQ26K8McrbGX/7NAe1LZ1ixyETJWoB2rJlC5KTk7Fs2TJkZ2cjOjoaiYmJqKqq6nX8vn37MHfuXCxcuBA5OTlISkpCUlIScnNzjWNaW1sxZcoUrFixYrheBpmp1T+eRoC7PcaPdBU7ChENIqlUgkenhaCjS4/kz47AYOBeYXQxiSDiLnJxcXGIjY3FmjVrAAAGgwF+fn54/PHHsWjRoovGz5kzB62trdi+fbvx2KRJkxATE4N169b1GFtYWIjAwEDk5OQgJiamX7mampqgUqnQ2NgIZ2deGrFEBwvqcM+7GUhOCENsoJvYcYhoCBwuqceK1HwsuiUCD00NFjsODYP+/PwW7QyQTqdDVlYWEhIS/htGKkVCQgIyMjJ6fU5GRkaP8QCQmJh4yfF91dnZiaamph4PsmyrfzyFke72GB/Asz9ElirGzxV3RPvgtdQ87DtTI3YcMjGiFaCamhro9Xqo1eoex9VqNbRaba/P0Wq1/RrfVykpKVCpVMaHn5/fVX0+Mm0HC+qw72wtZo/1hZRzf4gs2pwJfoj0UeGRTdkorW8TOw6ZENEnQZuCxYsXo7Gx0fgoKSkROxINEUEQ8Pr3eQhwt0dsAC99EVk6qVSCx24MgVwmxYP/yUKbrlvsSGQiRCtAHh4ekMlkqKys7HG8srISGo2m1+doNJp+je8rhUIBZ2fnHg+yTD+frsGhwnrcM8GPZ3+IrIST0hZP3RyGs9UteGrLYU6KJgAiFiC5XI7x48cjPT3deMxgMCA9PR3x8fG9Pic+Pr7HeABIS0u75Hii3xMEAa+l5iFc7YQYPxex4xDRMApwd8Cj00Lww/FKrPg+T+w4ZAJsxPziycnJmD9/PiZMmICJEydi9erVaG1txYIFCwAA8+bNg6+vL1JSUgAATzzxBKZOnYqVK1di5syZ2Lx5MzIzM/Hee+8ZP2ddXR2Ki4tRXl4OAMjPzwdw/uzR1Z4pIvOWmqvF8fImvHDbaK77Q2SFJox0w58mjcS7e85hhKs97ps0UuxIJCJRC9CcOXNQXV2NpUuXQqvVIiYmBqmpqcaJzsXFxZBK/3uSavLkydi0aROef/55LFmyBKGhodi2bRsiIyONY7755htjgQKAe++9FwCwbNkyLF++fHheGJkcXbcBr6bmIXqECqO9eYmTyFrdEqlBdUsnlm7LhbuDHLdGeYsdiUQi6jpAporrAFmeDb8W4KXtJ5By5zXwd7MXOw4RicggCHh79xkcLKjDhgUTcW2Ih9iRaJCYxTpARMOlqaMLq9NP4/pQT5YfIoJUIsFD1wdjjI8zFm48hIMFdWJHIhGwAJHFe3vXWbTr9PjDBK7vRETn2cikeCohHCGejrh//UFkFdWLHYmGGQsQWbTCmlZ8sPccZl7jDTcHudhxiMiEyG2keHp6OEa62WPeBwdwqJBngqwJCxBZtJe2n4CLnS3uiPYROwoRmSClrQx/nxGBAA8H3PfBAW6ZYUVYgMhi7cqrwk95Vfhj3EgobGRixyEiE6W0leHZxHCEqZ1w//pD+OH41W2vROaBBYgsUme3Hi9+exxjfJwxkbu9E9EVKGxkeGZ6OGL8XfDQx1n47BC3RLJ0LEBkkd7edRal9e2YHx/ARQ+JqE9sZVI8cWMopoV74e9fHMVb6afBlWIsl6gLIRINhTNVLXh79xncdo0P/HjbOxH1g1QqwcIpgXBzkGNl2imU1LfhldlRsJXxfIGlYQEiiyIIAv7vq2Nwd1Rg9lhfseMQkRmSSCS4c9wIeDop8N7P51BS1453/jQOLva8k9SSsNKSRdl8qAQHCurw52sDIbfhX28iGrjrQj2x+NZRyC1vxKy1v+JMVYvYkWgQ8ScEWYySuja8vP0EbozwQpSvSuw4RGQBRns74+VZkdAbBCSt/RXpJyvFjkSDhAWILILBIODvnx+FvVyGP8b5ix2HiCyI2lmJF+8YgwiNE/6yMRNvpZ+GwcDJ0eaOBYgswn/2FyHjXC3+en0w7OWc2kZEg8teboOnbg7DneNGYGXaKfz1P5lobO8SOxZdBRYgMnsnK5rwyncnMX20mpe+iGjISCUS3D1+BJ6dHo6Mc7W4/a29OFnRJHYsGiAWIDJr7To9HtuUDY1KiT/GjRQ7DhFZgXEjXfFKUhSkEiBp7a/YmslFE80RCxCZtRe/PY7S+nY8fmMI7/oiomFzfl5QJOKD3fHs50fx98+PoKNLL3Ys6gf+xCCzteVQMTYfKsH8+ACMcOWCh0Q0vOQ2Ujx4fTAevD4IXx8ux6w1v+JsNW+VNxcsQGSWjpQ04IVtx3FjhBemRXiJHYeIrNgN4V54aVYkmju7cftbe/H14TKxI1EfsACR2alq7sBDH2fB390e908OEDsOERH83ezxj1mRGOvvgic2H8biL4/xkpiJYwEis9Ku02Phhkx0dOnx5E2h3J+HiEyGnVyGR28IwQPXBeGLrFLMWsPVo00Zf3qQ2dAbBDyxOQenq5rxbGIE3B0VYkciIupBIpHgxggvvJwUiZbfLol9mV0qdizqBQsQmQVBELD061z8eLISj08LRaCHg9iRiIguyd/NHv9IikRsgCuSPzuCZ7YeRpuuW+xY9DssQGTyBEHAqzvz8MmBYvzluiCMG+kqdiQioitS2srw8A0heHhqMLYfrcDtb+1FnpYLJ5oKFiAyaYIgYPWPp/Huz+cwP34kpoXzji8iMi/Xh3niH0lR6DYImLXmV3y8vwiCwL3ExMYCRCZLEAS8mpqHN9JPY26sH2ZEeosdiYhoQHxd7PDSHZG4LtQTz2/LxcOfZKOxjXuJiYkFiEyS3iBg6dfH8e6ec5gXPxJ3xPiKHYmI6KrIbaRYOCUQTyaE4pfT1bjljZ+RWVgndiyrxQJEJqdN142H/pOFTw4U4S/XBeIWnvkhIgsSF+iOlNnXwMnOFve8m4E3fjyNbr1B7FhWhwWITIq2sQP3vrsfv5ypxtPTw3FThFrsSEREg87TSYEXZo7G7LG+eCP9FOa8tx8ldW1ix7IqLEBkMjLO1mLmm7+grKEdS28bg3H+vNuLiCyXTCrB3eP98MJto1Fc14YZq3/Gl9mlnCA9TFiASHTdegPW/HQaf/r3AXirlPh/s6O4zg8RWY0IjTNevTMKY/3Prxn08CfZqGvViR3L4kkEVs2LNDU1QaVSobGxEc7OzmLHsWjFtW146rPDyCmuxx3Rvrh7/AjIpBKxYxERiWL/uVp8uLcAclspUmZHYfoYjdiRzEp/fn7bDFMmoh669QZ8+GsBVv1wCs52tnjhttGI0LBsEpF1mxTkjnCNEz74pQB//U8WZsX4YPntY+DqIBc7msXhGaBe8AzQ0DpwrhbLvz2OfG0zEsdocM8EPyhtZWLHIiIyGYIgYO+ZGmzMKITCRoZlt4/GHdE+kEh4hvxy+vPzmwWoFyxAQ+NsdQte/z4fqblahHg54v7JAQj2dBQ7FhGRyWpo02FjRiH2n6vD9aEeeHFWJOdIXgYL0FViARpchTWtWLPrNL7MLoOrvRxzYv1wbYgHpPxNhoioT7KL6rExoxD1bTo8NDUYD98QDHs5Z7H8Lxagq8QCNDiOljbg3T3nsDO3As5KW8yK8cWNEV6Q2/DmQyKi/urs1mNbTjm+O1YOV3s5Ft0SgaQYX0h544gRC9BVYgEauI4uPVJztdiwrxCHSxqgdlZgZpQPpoZ5svgQEQ2CyqYOfHqwGAcK6hChccJzMyJwQ7gn5weBBeiqsQD1jyAIOFraiC+zS/FVThmaOroR5avC9NFqjPN35W8nRERD4FRlMzYfLMZJbTOi/VR44qZQTAv3suoixAJ0lViArkwQBJyqbMGOYxX4+kgZCmva4GJni+vDPHFDmCe8XezEjkhEZPEEQcCR0kZsyylDfmUzIjROeOC6INwe7WOVZ91ZgK4SC1DvdN0GZBXVI/1kJdJOVqKotg0OchnGjXTFlBAPjPFRcRFDIiIRCIKAExVN2H60AodLGuDheP6Gk3tj/eHnZi92vGHDAnSVWIDOu3CWZ9/ZGuw9U4OMs7Vo0+nham+Lsf6umDDSFVG+KtjIrO+3DCIiU1Va34a0E5XYe6YG7To9Jga64c5xvkgco4GLvWUvqGh2BWjt2rV4/fXXodVqER0djbfeegsTJ0685PitW7fihRdeQGFhIUJDQ7FixQrceuutxo8LgoBly5bh/fffR0NDA6699lq88847CA0N7VMeay1AbbpuHC1txOGSBmQW1iGzsB4N7V2wkUkQ5uWEqBEqRI9wwUh3e97CTkRk4jq69DhQUIdfz9Qgt6wRUokEcUFuuHm0GteHeSLIw8Hi5guZVQHasmUL5s2bh3Xr1iEuLg6rV6/G1q1bkZ+fDy8vr4vG79u3D9dffz1SUlJw2223YdOmTVixYgWys7MRGRkJAFixYgVSUlKwceNGBAYG4oUXXsCxY8dw4sQJKJXKK2ay9AIkCAIqGjtwqrIZ+dpmnKxoQm5ZE87VtMAgAEpbKYI9HRGucUKExhlhakcobLhSMxGRuapv053/xbaoHicrmtClF+CtUmJSkDtiA9wQ7adCmNoJtmZ+Rt+sClBcXBxiY2OxZs0aAIDBYICfnx8ef/xxLFq06KLxc+bMQWtrK7Zv3248NmnSJMTExGDdunUQBAE+Pj54+umn8cwzzwAAGhsboVarsWHDBtx7771XzGTuBUgQBDR1dKOyqQPlDe0oa2hHWX07imrbUFDTioKaVrR36QGcLzv+bvYY6e6AAHcHBHs6YISrPefyEBFZqI4uPU5UNOF4WSPyK5tRUNMKgwAobKQI9XJEuMYZIV6OCPQ4/7PB19UOzkpbsWP3idlshqrT6ZCVlYXFixcbj0mlUiQkJCAjI6PX52RkZCA5ObnHscTERGzbtg0AUFBQAK1Wi4SEBOPHVSoV4uLikJGR0WsB6uzsRGdnp/HPjY2NAM7/hxxKBoMAvSBAbxDQbRCg1wvoNhjQpTegSy9Ap9ejs8uAzm49OrsEtHfp0abrRptOj9bObrR06tHc0YXG9i40tOlQ39aF2tZO1LV2QddtMH4dqQRws5fD01kBTycFRkeooHFWwsfFDh6Oiv+5nGVAR2vLkL5uIiISV7irDcJd3YFId3R061FS14bC2laU1LUh63Qptmd1oF33358jjgoZPJ0U8HRUwt1RDld7W6js5HCys4Gj4vzDXmEDpY0UdrYyKG2lkNvIIJdJYSuTQHbhf6USyCTn/3coLr9d+Lndl3M7ohagmpoa6PV6qNXqHsfVajXy8vJ6fY5Wq+11vFarNX78wrFLjflfKSkpePHFFy867ufn17cXYgaKxA5ARERm7aTYAfqhubkZKpXqsmO4kQiAxYsX9zirZDAYUFdXB3d3d9EmiDU1NcHPzw8lJSVmeRnOWvF9Mz98z8wT3zfzMxzvmSAIaG5uho+PzxXHilqAPDw8IJPJUFlZ2eN4ZWUlNBpNr8/RaDSXHX/hfysrK+Ht7d1jTExMTK+fU6FQQKFQ9Djm4uLSn5cyZJydnfmP2wzxfTM/fM/ME9838zPU79mVzvxcIOp0b7lcjvHjxyM9Pd14zGAwID09HfHx8b0+Jz4+vsd4AEhLSzOODwwMhEaj6TGmqakJBw4cuOTnJCIiIusi+iWw5ORkzJ8/HxMmTMDEiROxevVqtLa2YsGCBQCAefPmwdfXFykpKQCAJ554AlOnTsXKlSsxc+ZMbN68GZmZmXjvvfcAABKJBE8++ST+8Y9/IDQ01HgbvI+PD5KSksR6mURERGRCRC9Ac+bMQXV1NZYuXQqtVouYmBikpqYaJzEXFxdDKv3viarJkydj06ZNeP7557FkyRKEhoZi27ZtxjWAAODvf/87Wltb8de//hUNDQ2YMmUKUlNT+7QGkKlQKBRYtmzZRZfmyLTxfTM/fM/ME98382Nq75no6wARERERDTfzXvKRiIiIaABYgIiIiMjqsAARERGR1WEBIiIiIqvDAkRERERWhwXIxBQWFmLhwoUIDAyEnZ0dgoODsWzZMuh0uh7jjh49iuuuuw5KpRJ+fn547bXXREpMF7zyyiuYPHky7O3tL7mSeHFxMWbOnAl7e3t4eXnh2WefRXd39/AGpR7Wrl2LgIAAKJVKxMXF4eDBg2JHot/8/PPPuP322+Hj4wOJRGLc9PoCQRCwdOlSeHt7w87ODgkJCTh9+rQ4YQnA+b01Y2Nj4eTkBC8vLyQlJSE/P7/HmI6ODjz66KNwd3eHo6Mj7rrrrot2eBgOLEAmJi8vDwaDAe+++y6OHz+Of/3rX1i3bh2WLFliHNPU1ITp06dj5MiRyMrKwuuvv47ly5cbF4Mkceh0OvzhD3/Aww8/3OvH9Xo9Zs6cCZ1Oh3379mHjxo3YsGEDli5dOsxJ6YItW7YgOTkZy5YtQ3Z2NqKjo5GYmIiqqiqxoxGA1tZWREdHY+3atb1+/LXXXsObb76JdevW4cCBA3BwcEBiYiI6OjqGOSldsGfPHjz66KPYv38/0tLS0NXVhenTp6O1tdU45qmnnsK3336LrVu3Ys+ePSgvL8edd945/GEFMnmvvfaaEBgYaPzz22+/Lbi6ugqdnZ3GY88995wQHh4uRjz6H+vXrxdUKtVFx3fs2CFIpVJBq9Uaj73zzjuCs7Nzj/eShs/EiROFRx991PhnvV4v+Pj4CCkpKSKmot4AEL766ivjnw0Gg6DRaITXX3/deKyhoUFQKBTCp59+KkJC6k1VVZUAQNizZ48gCOffI1tbW2Hr1q3GMSdPnhQACBkZGcOajWeAzEBjYyPc3NyMf87IyMD1118PuVxuPJaYmIj8/HzU19eLEZH6ICMjA1FRUcZVzoHz71tTUxOOHz8uYjLrpNPpkJWVhYSEBOMxqVSKhIQEZGRkiJiM+qKgoABarbbH+6dSqRAXF8f3z4Q0NjYCgPFnWFZWFrq6unq8bxEREfD39x/2940FyMSdOXMGb731Fh588EHjMa1W2+OHKADjn7Va7bDmo77j+2ZaampqoNfre31P+H6YvgvvEd8/02UwGPDkk0/i2muvNW5XpdVqIZfLL5onKcb7xgI0TBYtWgSJRHLZR15eXo/nlJWVYcaMGfjDH/6ABx54QKTk1m0g7xsREQGPPvoocnNzsXnzZrGj9Er0zVCtxdNPP43777//smOCgoKM/7+8vBzTpk3D5MmTL5rcrNFoLpoxf+HPGo1mcAITgP6/b5ej0WguusOI75t4PDw8IJPJev23xPfD9F14jyorK+Ht7W08XllZiZiYGJFS0QWPPfYYtm/fjp9//hkjRowwHtdoNNDpdGhoaOhxFkiMf3csQMPE09MTnp6efRpbVlaGadOmYfz48Vi/fj2k0p4n6uLj4/F///d/6Orqgq2tLQAgLS0N4eHhcHV1HfTs1qw/79uVxMfH45VXXkFVVRW8vLwAnH/fnJ2dMXr06EH5GtR3crkc48ePR3p6OpKSkgCcP2Wfnp6Oxx57TNxwdEWBgYHQaDRIT083Fp6mpiYcOHDgkndi0tATBAGPP/44vvrqK+zevRuBgYE9Pj5+/HjY2toiPT0dd911FwAgPz8fxcXFiI+PH/awZEJKS0uFkJAQ4aabbhJKS0uFiooK4+OChoYGQa1WC/fdd5+Qm5srbN68WbC3txfeffddEZNTUVGRkJOTI7z44ouCo6OjkJOTI+Tk5AjNzc2CIAhCd3e3EBkZKUyfPl04fPiwkJqaKnh6egqLFy8WObn12rx5s6BQKIQNGzYIJ06cEP76178KLi4uPe7UI/E0Nzcb/x0BEFatWiXk5OQIRUVFgiAIwquvviq4uLgIX3/9tXD06FFh1qxZQmBgoNDe3i5ycuv18MMPCyqVSti9e3ePn19tbW3GMQ899JDg7+8v/PTTT0JmZqYQHx8vxMfHD3tWFiATs379egFAr4/fO3LkiDBlyhRBoVAIvr6+wquvvipSYrpg/vz5vb5vu3btMo4pLCwUbrnlFsHOzk7w8PAQnn76aaGrq0u80CS89dZbgr+/vyCXy4WJEycK+/fvFzsS/WbXrl29/puaP3++IAjnb4V/4YUXBLVaLSgUCuGmm24S8vPzxQ1t5S7182v9+vXGMe3t7cIjjzwiuLq6Cvb29sLs2bN7/JI/XCS/BSYiIiKyGrwLjIiIiKwOCxARERFZHRYgIiIisjosQERERGR1WICIiIjI6rAAERERkdVhASIiIiKrwwJEREREVocFiIiIiKwOCxARERFZHRYgIiIisjr/H1AxGawBQ8YCAAAAAElFTkSuQmCC",
      "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": ".venv",
   "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.12.3"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}