sparsity_causal_simulations.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "03e44f3e",
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1f0d15cd",
   "metadata": {},
   "source": [
    "In this post, I'll introduce the econml python package and use it to compare double machine learning and doubly robust learning. I'll look at whether bootstrapping standard errors improves the coverage of confidence intervals, and I'll also look at whether sample size influences estimation accuracy. \n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7f488f07",
   "metadata": {},
   "source": [
    "I. Introduction to econml\n",
    "\n",
    "The econml package was developed by the [ALICE team](https://www.microsoft.com/en-us/research/group/alice/) at microsoft research to facilitate estimation of causal effects. This package builds upon [scikit-learn](https://scikit-learn.org/stable/) and presents several different algorithms for causal inference in one unified framework. \n",
    "\n",
    "One of the main reasons to use econml is its integration with scikit-learn. Because econml uses scikit-learn style estimators under the hood, it is easy to evaluate how a model would do if a different estimation style were used by simply swapping out one line in the code.\n",
    "\n",
    "Another big reason to use econml is that estimators defined in econml all inherit from the same base class. In practice, this means that the methods used to fit models and present treatment effects are the same, regardless of which estimator you are using. So, in this post, switching between double machine learning and doubly robust machine learning is very simple--no additional documentation reading is needed. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a0675303",
   "metadata": {},
   "source": [
    "II. Brief introduction to double machine and doubly robust learning. \n",
    "\n",
    "Doubly robust machine learning involves fitting two models: one model to predict the treatment X from possible confounders, and another model to predict the outcome Y from X as well as possible confounders. One of the parameters estimated in the second model is used to estimate the average treatment effect of X on Y. This method was developed by [Robins et al.](https://www.tandfonline.com/doi/pdf/10.1080/01621459.1994.10476818?needAccess=true&role=button) in 1994, and has been frequently been used under the name inverse probability weighting. For more on doubly robust machine learning, see [this](https://towardsdatascience.com/doubly-robust-estimators-for-causal-inference-in-statistical-estimation-3c00847e9db) blog post.\n",
    "\n",
    "Double machine learning involves fitting two models: one model to predict X from confounders and another model to predict Y from the confounders. The residuals of the two models are then used to generate estimates of the causal effect of X on Y. I wrote [this post](https://medium.com/@clementzach_38631/double-machine-learning-for-causal-inference-from-a-partially-linear-model-ada4c39914e3) on double machine learning for more of an overview. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "74adbbdd",
   "metadata": {},
   "source": [
    "III. Data generation. \n",
    "\n",
    "For this simulation, I generated data in a similar method as my previous [blog post](https://medium.com/@clementzach_38631/causal-inference-with-high-dimensional-data-an-evaluation-of-the-double-machine-learning-algorithm-a882d2135bd5). However, I converted x to a binary variable because doubly robust learning only works with categorical variables. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "5128afda",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "def generate_data(n = 400, n_confounders = 10, \n",
    "                n_possible_confounders = 300, \n",
    "                causal_effect_size = 2, \n",
    "                beta_c_x = None, beta_c_y = 1, \n",
    "                sd_x_error = None, \n",
    "                sd_confounders = 1):\n",
    "    \n",
    "    if beta_c_x is None:\n",
    "        beta_c_x = np.sqrt(1/(n_confounders + 1)) #this lets our x have a variance of 1\n",
    "    \n",
    "    if sd_x_error is None:\n",
    "        sd_x_error = np.sqrt(1/(n_confounders + 1)) #same as above\n",
    "    \n",
    "    zero_sequence = np.zeros(n) # a sequence of zeroes the length of our model\n",
    "    data_dict = dict()\n",
    "    x = m1 = m2 = zero_sequence\n",
    "    \n",
    "    for i in range(n_confounders): # We will have some variables that are true confounders\n",
    "        c = np.random.normal(size = n, scale = sd_confounders) # Create some variance in the confounder\n",
    "        m1 = m1 + c*beta_c_x\n",
    "        m2 = m2 + c*beta_c_y\n",
    "\n",
    "        data_dict[\"c\" + str(i)] = c\n",
    "        \n",
    "    for i in range(n_confounders,n_possible_confounders): #because we are testing sparsity, the rest of the variables will have zero effect. \n",
    "        c = np.random.normal(size = n, scale = sd_confounders)\n",
    "\n",
    "        data_dict[\"fc\" + str(i)] = c #label them as fake\n",
    "    x = ((np.random.normal(size = n, scale = sd_x_error) + m1 ) > 0) * 1 #convert to binary to make more estimation methods\n",
    "\n",
    "\n",
    "    y = x * causal_effect_size  + m2 + np.random.normal(size = n, scale = np.std(m2) /10)\n",
    "\n",
    "    data_dict[\"x\"] = x #replaces the column of zeroes we made before\n",
    "    data_dict[\"y\"] = y\n",
    "    \n",
    "    output_df = pd.DataFrame(data_dict)\n",
    "    \n",
    "    return output_df\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4145b693",
   "metadata": {},
   "source": [
    "Next, I defined a function to generate an estimate of the average treatment effect from a simulated dataset. As you can see, I can pass in either a LinearDRLearner or a LinearDML object to the function and use the same methods on the object to get simulated effects. \n",
    "\n",
    "Both methods generate confidence intervals using traditional statistical techniques (\"statsmodels\") or with bootstrapping. Here, I allow my function to conduct estimation under both techniques. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "17504b61",
   "metadata": {},
   "outputs": [],
   "source": [
    "from econml.dr import LinearDRLearner\n",
    "from econml.dml import LinearDML\n",
    "\n",
    "\n",
    "\n",
    "estimators_to_test = {\n",
    "    \"dml_linear_bootstrap\": [LinearDML(), \"bootstrap\"],\n",
    "    \"dml_linear_statsmodels\": [LinearDML(), \"statsmodels\"],\n",
    "    \"dr_linear_bootstrap\": [LinearDRLearner(), \"bootstrap\"],\n",
    "    \"dr_linear_statsmodels\": [LinearDRLearner(), \"statsmodels\"]\n",
    "}\n",
    "\n",
    "def get_econml_estimate(sim_df, est, inference_method = \"bootstrap\", sample_size = 0, estimator_name = \"\", alpha = 0.05):\n",
    "    est.fit(sim_df[\"y\"], T = sim_df[\"x\"], W=sim_df.drop([\"x\",\"y\"], axis = 1), inference = inference_method)\n",
    "    lb, ub = est.ate_interval(alpha=alpha) # OLS confidence intervals\n",
    "    point_estimate = est.ate()\n",
    "    return {\"lower\": lb, \"estimate\": point_estimate, \"upper\": ub, \"n\": sample_size, \"estimator_name\": estimator_name}\n",
    "    \n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fd7464f0",
   "metadata": {},
   "source": [
    "I'll also define a function to get effects using a standard linear model. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "9d5b79a1",
   "metadata": {},
   "outputs": [],
   "source": [
    "import statsmodels.api as sm\n",
    "def get_lm_interval(sim_df, sample_size = 0, alpha = 0.05):\n",
    "    y = sim_df[\"y\"]\n",
    "    #sim_df[\"interaction\"] = sim_df[\"moderator\"] * sim_df[\"x\"]\n",
    "    x = sim_df.drop(\"y\", axis = 1)\n",
    "    ols_model = sm.OLS(y, x)\n",
    "    ols_results = ols_model.fit()\n",
    "    results_df = pd.DataFrame(ols_results.summary(alpha = alpha).tables[1].data)\n",
    "    results_df.index = results_df[0]\n",
    "    \n",
    "    return {\"lower\": float(results_df.loc[\"x\",5]), \n",
    "            \"upper\": float(results_df.loc[\"x\",6]),\n",
    "            \"estimate\": float(results_df.loc[\"x\",1]),\n",
    "            \"n\": sample_size, \n",
    "            \"estimator_name\": \"standard_linear\"\n",
    "           }\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "fae9268e",
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_lm_interval_small(sim_df, sample_size = 0, alpha = 0.05):\n",
    "    y = sim_df[\"y\"]\n",
    "    x = sim_df.drop([\"y\"] + [col for col in sim_df.columns if col.startswith(\"f\")], axis = 1)\n",
    "    ols_model = sm.OLS(y, x)\n",
    "    ols_results = ols_model.fit()\n",
    "    results_df = pd.DataFrame(ols_results.summary(alpha = alpha).tables[1].data)\n",
    "    results_df.index = results_df[0]\n",
    "    \n",
    "    return {\"lower\": float(results_df.loc[\"x\",5]), \n",
    "            \"upper\": float(results_df.loc[\"x\",6]),\n",
    "            \"estimate\": float(results_df.loc[\"x\",1]),\n",
    "            \"n\": sample_size, \n",
    "            \"estimator_name\": \"standard_linear_small\"\n",
    "           }"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "04bc7d6a",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "400\n",
      "1600\n",
      "6400\n"
     ]
    }
   ],
   "source": [
    "import warnings\n",
    "warnings.filterwarnings('ignore')\n",
    "import pickle\n",
    "\n",
    "all_estimates = []\n",
    "num_reps = 50\n",
    "for sample_size in [400, 1600, 6400]:\n",
    "    print(sample_size)\n",
    "    for i in range(num_reps):\n",
    "        sim_df = generate_data(sample_size)\n",
    "        for estimator_name in estimators_to_test.keys():\n",
    "            estimate_dict = get_econml_estimate(sim_df, \n",
    "                                                estimators_to_test[estimator_name][0], \n",
    "                                                estimators_to_test[estimator_name][1], \n",
    "                                               sample_size = sample_size,\n",
    "                                               estimator_name = estimator_name)\n",
    "            all_estimates.append(estimate_dict)\n",
    "        all_estimates.append(get_lm_interval_small(sim_df, sample_size = sample_size))\n",
    "        all_estimates.append(get_lm_interval(sim_df, sample_size = sample_size))\n",
    "    with open('simulations_4_feb.pickle', 'wb') as handle:\n",
    "        pickle.dump(all_estimates, handle)  \n",
    "        \n",
    "        \n",
    "            \n",
    "            \n",
    "            \n",
    "        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "9bf3d2c9",
   "metadata": {},
   "outputs": [],
   "source": [
    "import pickle\n",
    "\n",
    "with open('all_estimates_4_feb.pickle', 'wb') as handle:\n",
    "    pickle.dump(all_estimates, handle, protocol=pickle.HIGHEST_PROTOCOL)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "81ddcbee",
   "metadata": {},
   "outputs": [],
   "source": [
    "all_estimates_df = pd.DataFrame(all_estimates)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "a637ba84",
   "metadata": {},
   "outputs": [],
   "source": [
    "all_estimates_df[\"width\"] = all_estimates_df[\"upper\"] - all_estimates_df[\"lower\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "8384a4f1",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>lower</th>\n",
       "      <th>estimate</th>\n",
       "      <th>upper</th>\n",
       "      <th>n</th>\n",
       "      <th>estimator_name</th>\n",
       "      <th>width</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1.773386</td>\n",
       "      <td>1.936886</td>\n",
       "      <td>2.100386</td>\n",
       "      <td>400</td>\n",
       "      <td>dml_linear_bootstrap</td>\n",
       "      <td>0.326999</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>1.813699</td>\n",
       "      <td>1.900463</td>\n",
       "      <td>1.987227</td>\n",
       "      <td>400</td>\n",
       "      <td>dml_linear_statsmodels</td>\n",
       "      <td>0.173528</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>-699.627111</td>\n",
       "      <td>843.452201</td>\n",
       "      <td>2386.531513</td>\n",
       "      <td>400</td>\n",
       "      <td>dr_linear_bootstrap</td>\n",
       "      <td>3086.158625</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>-4117.111658</td>\n",
       "      <td>745.121911</td>\n",
       "      <td>5607.355481</td>\n",
       "      <td>400</td>\n",
       "      <td>dr_linear_statsmodels</td>\n",
       "      <td>9724.467139</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>1.902000</td>\n",
       "      <td>1.951100</td>\n",
       "      <td>2.000000</td>\n",
       "      <td>400</td>\n",
       "      <td>standard_linear_small</td>\n",
       "      <td>0.098000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>1.841000</td>\n",
       "      <td>1.940600</td>\n",
       "      <td>2.041000</td>\n",
       "      <td>400</td>\n",
       "      <td>standard_linear</td>\n",
       "      <td>0.200000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>1.827129</td>\n",
       "      <td>2.009646</td>\n",
       "      <td>2.192164</td>\n",
       "      <td>400</td>\n",
       "      <td>dml_linear_bootstrap</td>\n",
       "      <td>0.365035</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>1.964629</td>\n",
       "      <td>2.053412</td>\n",
       "      <td>2.142196</td>\n",
       "      <td>400</td>\n",
       "      <td>dml_linear_statsmodels</td>\n",
       "      <td>0.177567</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>-2318.921907</td>\n",
       "      <td>2.247979</td>\n",
       "      <td>2323.417865</td>\n",
       "      <td>400</td>\n",
       "      <td>dr_linear_bootstrap</td>\n",
       "      <td>4642.339772</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>2.198664</td>\n",
       "      <td>2.270068</td>\n",
       "      <td>2.341472</td>\n",
       "      <td>400</td>\n",
       "      <td>dr_linear_statsmodels</td>\n",
       "      <td>0.142808</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "         lower    estimate        upper    n          estimator_name  \\\n",
       "0     1.773386    1.936886     2.100386  400    dml_linear_bootstrap   \n",
       "1     1.813699    1.900463     1.987227  400  dml_linear_statsmodels   \n",
       "2  -699.627111  843.452201  2386.531513  400     dr_linear_bootstrap   \n",
       "3 -4117.111658  745.121911  5607.355481  400   dr_linear_statsmodels   \n",
       "4     1.902000    1.951100     2.000000  400   standard_linear_small   \n",
       "5     1.841000    1.940600     2.041000  400         standard_linear   \n",
       "6     1.827129    2.009646     2.192164  400    dml_linear_bootstrap   \n",
       "7     1.964629    2.053412     2.142196  400  dml_linear_statsmodels   \n",
       "8 -2318.921907    2.247979  2323.417865  400     dr_linear_bootstrap   \n",
       "9     2.198664    2.270068     2.341472  400   dr_linear_statsmodels   \n",
       "\n",
       "         width  \n",
       "0     0.326999  \n",
       "1     0.173528  \n",
       "2  3086.158625  \n",
       "3  9724.467139  \n",
       "4     0.098000  \n",
       "5     0.200000  \n",
       "6     0.365035  \n",
       "7     0.177567  \n",
       "8  4642.339772  \n",
       "9     0.142808  "
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "all_estimates_df.head(10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "id": "78578a1a",
   "metadata": {},
   "outputs": [],
   "source": [
    "all_estimates_df.to_pickle(\"all_estimates_5_feb.pkl\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "43188df6",
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "all_estimates_df = pd.read_pickle(\"all_estimates_5_feb.pkl\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "43024732",
   "metadata": {},
   "source": [
    "First, let's look at whether different estimators are biased. \n",
    "\n",
    "The bias of an estimator d is defined as the expected value of the estimator minus the true parameter value: \n",
    "\n",
    "$$E(\\hat{d} - d)$$\n",
    "\n",
    "For an estimator to be useful, the estimator must have zero bias. If an estimator systematically delivers high or low estimates, that estimator is useless. In the presence of confounding (such as in this situation), it is especially important for the estimator to eliminate bias which originates from confounders.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 52,
   "id": "c3ce6677",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "## Bias\n",
    "\n",
    "def bias(x):\n",
    "    return np.mean(x) - 2\n",
    "def mse(x):\n",
    "    return np.mean((np.array(x) - 2) * (np.array(x) - 2))\n",
    "\n",
    "def sd(x):\n",
    "    return np.std(x)\n",
    "    \n",
    "bias_df = all_estimates_df.groupby([\"estimator_name\", \"n\"])[\"estimate\"].agg(\n",
    "    [bias, sd]\n",
    ")\n",
    "\n",
    "bias_df[\"bias_sd\"] = np.round(bias_df[\"bias\"]).astype(str) + np.round(bias_df[\"sd\"]).astype(str)\n",
    "\n",
    "#"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "id": "76c9b471",
   "metadata": {},
   "outputs": [],
   "source": [
    "bias_df[\"bias_sd\"] = np.round(bias_df[\"bias\"], 3).astype(str) + \" (\" + np.round(bias_df[\"sd\"]/ np.sqrt(50), 3).astype(str) + \")\""
   ]
  },
  {
   "cell_type": "markdown",
   "id": "aad6f756",
   "metadata": {},
   "source": [
    "|    |   A |   B |\n",
    "|:---|----:|----:|\n",
    "| a  |   1 |   1 |\n",
    "| a  |   2 |   2 |\n",
    "| b  |   3 |   3 |"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "id": "ba45c0f4",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th>n</th>\n",
       "      <th>400</th>\n",
       "      <th>1600</th>\n",
       "      <th>6400</th>\n",
       "      <th>Name</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0.017 (0.009)</td>\n",
       "      <td>0.002 (0.003)</td>\n",
       "      <td>0.003 (0.002)</td>\n",
       "      <td>Double Machine Learning</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>161.905 (108.094)</td>\n",
       "      <td>-188.097 (119.696)</td>\n",
       "      <td>0.048 (0.027)</td>\n",
       "      <td>Doubly Robust Estimation</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>-0.004 (0.008)</td>\n",
       "      <td>-0.001 (0.002)</td>\n",
       "      <td>0.001 (0.001)</td>\n",
       "      <td>Standard Linear Model</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "n                400                1600           6400  \\\n",
       "0      0.017 (0.009)       0.002 (0.003)  0.003 (0.002)   \n",
       "2  161.905 (108.094)  -188.097 (119.696)  0.048 (0.027)   \n",
       "4     -0.004 (0.008)      -0.001 (0.002)  0.001 (0.001)   \n",
       "\n",
       "n                      Name  \n",
       "0   Double Machine Learning  \n",
       "2  Doubly Robust Estimation  \n",
       "4     Standard Linear Model  "
      ]
     },
     "execution_count": 54,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "bias_df = bias_df.reset_index().pivot(index = \"estimator_name\", columns = \"n\", values = \"bias_sd\").reset_index().loc[[0, 2, 4],:]\n",
    "\n",
    "\n",
    "bias_df[\"Name\"] = [\"Double Machine Learning\", \"Doubly Robust Estimation\", \"Standard Linear Model\"]\n",
    "\n",
    "bias_df.drop(\"estimator_name\", axis = 1, inplace = True)\n",
    "bias_df"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "79978063",
   "metadata": {},
   "source": [
    "Next, let's look at mean squared error, or how close estimates tend to be to the actual parameter value. An estimator with a lower mean squared error is more useful than an estimator with a higher mean squared error because it will estimate the parameter of interest with more precision. \n",
    "\n",
    "Mean squared error is defined as:\n",
    "\n",
    "$$E[(\\hat{d} - d)^2]$$\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 168,
   "id": "e4510bb6",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "## mean squared error\n",
    "mse_df = all_estimates_df.groupby([\"estimator_name\", \"n\"])[\"estimate\"].agg([mse]).reset_index()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 169,
   "id": "b6952924",
   "metadata": {},
   "outputs": [],
   "source": [
    "translation_dict = {\"dml_linear_bootstrap\": \"Double Machine Learning (Bootstrap)\",\n",
    "                    \"dml_linear_statsmodels\": \"Double Machine Learning (Bootstrap)\",\n",
    "                    \"dr_linear_bootstrap\": \"Doubly Robust Estimation (Bootstrap)\",\n",
    "                    \"dr_linear_statsmodels\": \"Doubly Robust Estimation (Bootstrap)\",\n",
    "                    \"standard_linear\": \"Standard Linear Model\"\n",
    "                   }\n",
    "\n",
    "translation_dict = {\"dml_linear_bootstrap\": \"Double Machine Learning\",\n",
    "                    \"dml_linear_statsmodels\": \"Double Machine Learning\",\n",
    "                    \"dr_linear_bootstrap\": \"Doubly Robust Estimation\",\n",
    "                    \"dr_linear_statsmodels\": \"Doubly Robust Estimation\",\n",
    "                    \"standard_linear\": \"Standard Linear Model\"\n",
    "                   }"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 170,
   "id": "94a5b294",
   "metadata": {},
   "outputs": [],
   "source": [
    "mse_df = mse_df.loc[mse_df[\"estimator_name\"].isin(\n",
    "    [\"dml_linear_statsmodels\", \"dr_linear_statsmodels\", \"standard_linear\"]),:]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 171,
   "id": "f71b5692",
   "metadata": {},
   "outputs": [],
   "source": [
    "mse_df[\"Estimator\"] = mse_df[\"estimator_name\"].apply(lambda x: translation_dict[x])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 172,
   "id": "f97b8a07",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAs4AAAHUCAYAAAAqSa5MAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAABlWklEQVR4nO3deVhV5f7//9eWGUEUkMkQ0dQkMxOKoJwyFU3D4SSdDMcsslKk0pzSrKMH+6Yec6zjnCl1zKEyc0ZNMucsyaxwyCCHFBxBYP3+8Of+tGVwYehWez6ua1+17vVe932vxWbx9t73vpfFMAxDAAAAAEpVwd4dAAAAAG4FJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgA0e/ZsWSwWWSwWrV+/vsh+wzB05513ymKxqFmzZje8f2Vx8eJFTZ8+Xffff7+8vb3l7u6ukJAQxcbGavHixfbu3nW3fv36En+Of/bnn3lxr6sdfytISUnR3XffLTc3N1ksFu3atavYuMvXrKTX7NmzTbd57tw5jRw5stjrd/maHzhw4JrO56/47bffNHLkyBKvAQBzHO3dAQA3D09PT82YMaNIcpyamqqff/5Znp6e9ulYGcTHx+uTTz5RYmKi3njjDbm4uOiXX37RihUr9OWXX6pjx4727uJNZdasWbrrrruKlIeFhdmhN+Xn2LFjio+PV0xMjKZMmSIXFxfVqVOn1GNGjx6t5s2bFymvVauW6XbPnTunN954Q5KK/B499thjSktLU2BgoOn6ystvv/2mN954QzVq1FDDhg1vePvA7YLEGYBVXFyc5s+fr8mTJ6tSpUrW8hkzZigqKko5OTl27N3VZWRkKCUlRa+//ro1eZGkFi1aqE+fPiosLLRj70pnGIYuXLggNze3G9pu/fr1FRERUaZjSuvr+fPn5erqKovFcs19OnfunNzd3a/5eEn68ccfdfHiRT399NNq2rSpqWNq166tBx988C+1W5qqVauqatWq161+eyiPnxVwK2GqBgCrf/7zn5KkBQsWWMuys7O1aNEi9erVq9hj8vLy9NZbb+muu+6Si4uLqlatqp49e+rYsWM2cSkpKWrVqpUCAwPl5uamevXq6bXXXtPZs2dt4nr06CEPDw/99NNPatu2rTw8PBQcHKyXX35Zubm5pfb/xIkTklTiiF6FCra3vB9++EExMTFyd3eXr6+vEhIS9OmnnxaZqlCjRg316NGjSH3NmjWzGVW8cOGCXn75ZTVs2FBeXl7y9vZWVFSUli5dWuRYi8WiF198UdOmTVO9evXk4uKiOXPmSJL279+vp556Sn5+fnJxcVG9evU0efLkInUU1//Tp0+Xeo2uRUl9vTz1YOXKlerVq5eqVq0qd3d35ebmqrCwUGPHjrW+L/z8/NStWzf9+uuvNnU3a9ZM9evX14YNGxQdHS13d/cS32uXLVu2TFFRUXJ3d5enp6datmyptLQ06/4ePXro4YcflnTpH4PlOcVo7dq1atasmXx8fOTm5qbq1aurc+fOOnfunA4cOGBNjN944w3rVI/L753ipmpcPv+0tDRFR0fLzc1NNWrU0KxZsyRJn3/+uRo1aiR3d3fdc889WrFihU1/fvrpJ/Xs2VO1a9eWu7u7qlWrpvbt22vPnj3WmPXr1+v++++XJPXs2dPar5EjR5q+ppI0cuRIWSwW7dixQ//4xz9UpUoV62j8L7/8oieffFJBQUFycXGRv7+/WrRowdQQ3HYYcQZgValSJf3jH//QzJkz9dxzz0m6lERXqFBBcXFxmjBhgk18YWGhYmNjtXHjRg0cOFDR0dE6ePCgRowYoWbNmmnbtm3WUcn9+/erbdu2SkxMVMWKFfXDDz8oOTlZ33zzjdauXWtT78WLF/X444+rd+/eevnll7Vhwwa9+eab8vLy0uuvv15i/+vVq6fKlSvrjTfeUIUKFdSqVSvVqFGj2Njff/9dTZs2lZOTk6ZMmSJ/f3/Nnz9fL7744jVfv9zcXP3xxx965ZVXVK1aNeXl5Wn16tXq1KmTZs2apW7dutnEL1myRBs3btTrr7+ugIAA+fn5ae/evYqOjlb16tX1zjvvKCAgQF9++aX69eun48ePa8SIEeXa/4KCAuXn59uUWSwWOTg4XLWvW7dulST16tVLjz32mObNm6ezZ8/KyclJzz//vN577z29+OKLateunQ4cOKDhw4dr/fr12rFjh3x9fa11Z2Zm6umnn9bAgQM1evToIv/A+bMPP/xQXbt2VatWrbRgwQLl5uZq7NixatasmdasWaOHH35Yw4cP1wMPPKAXXnjBOv3iz5+glKSwsLDItZAkR8dLfyoPHDigxx57TI0bN9bMmTNVuXJlHTlyRCtWrFBeXp4CAwO1YsUKxcTEqHfv3nrmmWck6aqjzFlZWerZs6cGDhyoO+64Q++++6569eqlw4cP63//+5+GDBkiLy8vjRo1Sh06dNAvv/yioKAgSZemYPj4+Ojf//63qlatqj/++ENz5sxRZGSkdu7cqbp166pRo0aaNWuWevbsqWHDhumxxx6TJN1xxx2mr+mfderUSU8++aQSEhKs//Bt27atCgoKNHbsWFWvXl3Hjx/X5s2bderUqated+CWYgD425s1a5Yhydi6dauxbt06Q5Lx3XffGYZhGPfff7/Ro0cPwzAM4+677zaaNm1qPW7BggWGJGPRokU29W3dutWQZEyZMqXY9goLC42LFy8aqamphiRj9+7d1n3du3c3JBkfffSRzTFt27Y16tate9Vz+fzzzw1fX19DkiHJ8PHxMZ544glj2bJlNnGDBg0yLBaLsWvXLpvyli1bGpKMdevWWctCQkKM7t27F2mradOmNtfjSvn5+cbFixeN3r17G/fdd5/NPkmGl5eX8ccff9iUt27d2rjjjjuM7Oxsm/IXX3zRcHV1tcaXpf/FufwzL+7l4OBgqq+X6+jWrZtNeXp6uiHJ6Nu3r035li1bDEnGkCFDrGVNmzY1JBlr1qwptb+GYRgFBQVGUFCQcc899xgFBQXW8tOnTxt+fn5GdHS0tezy+/jjjz++ar2XY0t6HT582DAMw/jf//5nSCpyzf/s2LFjhiRjxIgRRfZdvl4ZGRlFzn/btm3WshMnThgODg6Gm5ubceTIEWv5rl27DEnGxIkTS2w/Pz/fyMvLM2rXrm0MGDDAWn75d3LWrFk28WW5piNGjDAkGa+//rpNHcePHzckGRMmTCixX8DtgqkaAGw0bdpUtWrV0syZM7Vnzx5t3bq1xI/OP/vsM1WuXFnt27dXfn6+9dWwYUMFBATYTHf45Zdf9NRTTykgIEAODg5ycnKyzj1NT0+3qddisah9+/Y2ZQ0aNNDBgwev2v+2bdvq0KFDWrx4sV555RXdfffdWrJkiR5//HGb0dh169bp7rvv1r333mtz/FNPPXXVNkrz8ccf66GHHpKHh4ccHR3l5OSkGTNmFDlHSXrkkUdUpUoV6/aFCxe0Zs0adezYUe7u7jbXtG3btrpw4YK+/vrrcu3/3LlztXXrVpvXli1brtrXP+vcubPN9rp16ySpyPSWBx54QPXq1dOaNWtsyqtUqaJHHnnkqn3dt2+ffvvtN8XHx9uMSnt4eKhz5876+uuvde7cuavWU5Lk5OQi12Lr1q3y9/eXJDVs2FDOzs569tlnNWfOHP3yyy/X3NafBQYGKjw83Lrt7e0tPz8/NWzY0DqyLF36REWSze9Bfn6+Ro8erbCwMDk7O8vR0VHOzs7av39/se+5K13LNb3y5+3t7a1atWrp7bff1rhx47Rz586b+vsEwF/BVA0ANiwWi3r27KmJEyfqwoULqlOnjho3blxs7O+//65Tp07J2dm52P3Hjx+XJJ05c0aNGzeWq6ur3nrrLdWpU0fu7u46fPiwOnXqpPPnz9sc5+7uLldXV5syFxcXXbhwwdQ5uLm5qUOHDurQoYMk6dChQ2rTpo0mT56s559/XnfffbdOnDih0NDQIscGBASYaqM4n3zyibp06aInnnhCr776qgICAuTo6KipU6dq5syZReKvnIt94sQJ5efn691339W7775bbBuXr2l59b9evXqmvhxY2koQxZ1HSccEBQUV+QeQ2VUmrlZvYWGhTp48ec1fVqtZs2ap16JWrVpavXq1xo4dqxdeeEFnz55VzZo11a9fP/Xv3/+a2pQuJZ5XcnZ2LlJ++ffsz78HSUlJmjx5sgYNGqSmTZuqSpUqqlChgp555pkiv1fFuZZremWsxWLRmjVrNGrUKI0dO1Yvv/yyvL291bVrV/3rX/+6JVbjAcwicQZQRI8ePfT6669r2rRp+te//lVinK+vr3x8fIp8Yemyy38w165dq99++03r16+3WeHgRs1/rF69up599lklJibq+++/19133y0fHx9lZWUViS2uzNXVtdgvJh4/ftxmru4HH3yg0NBQpaSk2KwqUdKXGq9ceaJKlSpycHBQfHy8XnjhhWKPuZwsl6X/5aG0VTKu3Ofj4yPp0tzly/NoL/vtt99srtnV6i6p3iv99ttvqlChQomj4uWlcePGaty4sQoKCrRt2za9++67SkxMlL+/v5588snr2nZxPvjgA3Xr1k2jR4+2KT9+/LgqV6581eOv5ZoW9/MKCQnRjBkzJF1a0eSjjz7SyJEjlZeXp2nTppk9HeCmx1QNAEVUq1ZNr776qtq3b6/u3buXGNeuXTudOHFCBQUFioiIKPKqW7eupP/7Q+vi4mJz/PTp08u136dPn9aZM2eK3Xf5Y+vLH303b95c33//vXbv3m0T9+GHHxY5tkaNGvr2229tyn788Uft27fPpsxiscjZ2dkmscjKyip2VY3iuLu7q3nz5tq5c6caNGhQ7DW9nOiUpf832uVpFx988IFN+datW5Wenq4WLVpcU71169ZVtWrV9OGHH8owDGv52bNntWjRIuuqEDeCg4ODIiMjraud7NixQ9L/vcfNjPaWB4vFUuT36vPPP9eRI0dsykrq1/W4pnXq1NGwYcN0zz33WK8LcLtgxBlAsf79739fNebJJ5/U/Pnz1bZtW/Xv318PPPCAnJyc9Ouvv2rdunWKjY1Vx44dFR0drSpVqighIUEjRoyQk5OT5s+fXyTp+6v27dun1q1b68knn1TTpk0VGBiokydP6vPPP9d7772nZs2aKTo6WpKUmJiomTNn6rHHHtNbb71lXZXihx9+KFJvfHy8nn76afXt21edO3fWwYMHNXbs2CKrJbRr106ffPKJ+vbtq3/84x86fPiw3nzzTQUGBmr//v2mzuE///mPHn74YTVu3FjPP/+8atSoodOnT+unn37Sp59+al2BpCz9L813331X7EoStWrVuuY1h+vWratnn31W7777ripUqKA2bdpYV9UIDg7WgAEDrqneChUqaOzYseratavatWun5557Trm5uXr77bd16tQpU+/Z0uzfv986h/zP7rjjDt1xxx2aNm2a1q5dq8cee0zVq1fXhQsXrFNwHn30UUmXPmUJCQnR0qVL1aJFC3l7e8vX17fE1V3+qnbt2mn27Nm666671KBBA23fvl1vv/12kZH+WrVqyc3NTfPnz1e9evXk4eGhoKAgBQUF/eVr+u233+rFF1/UE088odq1a8vZ2Vlr167Vt99+q9dee+26nDdgN/b+diIA+/vzqhqluXJVDcMwjIsXLxr/7//9P+Pee+81XF1dDQ8PD+Ouu+4ynnvuOWP//v3WuM2bNxtRUVGGu7u7UbVqVeOZZ54xduzYUeSb/t27dzcqVqxYpO3L3+gvzcmTJ4233nrLeOSRR4xq1aoZzs7ORsWKFY2GDRsab731lnHu3Dmb+L179xotW7Y0XF1dDW9vb6N3797G0qVLi6xKUVhYaIwdO9aoWbOm4erqakRERBhr164tdlWNf//730aNGjUMFxcXo169esb7779fbN8lGS+88EKx55GRkWH06tXLqFatmuHk5GRUrVrViI6ONt56661r6n9xSltVQ5Lx/vvvX7Wvpb1vCgoKjOTkZKNOnTqGk5OT4evrazz99NPWFSoua9q0qXH33XeX2tcrLVmyxIiMjDRcXV2NihUrGi1atDC++uorm5jyXFVj6NChhmEYRlpamtGxY0cjJCTEcHFxMXx8fIymTZsWWbFl9erVxn333We4uLgYkqwrspS0qkZx5x8SEmI89thjRcqv/FmcPHnS6N27t+Hn52e4u7sbDz/8sLFx48Zi35sLFiww7rrrLsPJyanIyh9mrunl9/GxY8dsyn///XejR48exl133WVUrFjR8PDwMBo0aGCMHz/eyM/PL/G6A7cii2H86bMZAPibW79+vZo3b65169aV20MzAAC3B+Y4AwAAACaQOAMAAAAmMFUDAAAAMIERZwAAAMAEEmcAAADABBJnAAAAwAQegHIdFRYW6rfffpOnp6fpR8oCAADgxjEMQ6dPn1ZQUJAqVCh9TJnE+Tr67bffFBwcbO9uAAAA4CoOHz5c5KmbVyJxvo48PT0lXfpBVKpUyc69AQAAwJVycnIUHBxszdtKQ+J8HV2enlGpUiUSZwAAgJuYmWm1fDkQAAAAMIHEGQAAADCBxBkAAAAwgTnOAAAAN4hhGMrPz1dBQYG9u/K34eDgIEdHx3JZGpjEGQAA4AbIy8tTZmamzp07Z++u/O24u7srMDBQzs7Of6keEmcAAIDrrLCwUBkZGXJwcFBQUJCcnZ15ONoNYBiG8vLydOzYMWVkZKh27dpXfchJaUicAQAArrO8vDwVFhYqODhY7u7u9u7O34qbm5ucnJx08OBB5eXlydXV9Zrr4suBAAAAN8hfGe3EtSuv685PDwAAADCBxBkAAAAwgcQZAAAApbJYLFqyZMl1bWPkyJFq2LDhdW3jryJxBgAAsLOjR4/queeeU/Xq1eXi4qKAgAC1bt1aaWlp9u5auVm0aJEiIyPl5eUlT09P3X333Xr55Zet+1955RWtWbPGjj28OlbVAAAAsLPOnTvr4sWLmjNnjmrWrKnff/9da9as0R9//GHvrpWL1atX68knn9To0aP1+OOPy2KxaO/evTaJsoeHhzw8POzYy6tjxBkAAMCOTp06pU2bNik5OVnNmzdXSEiIHnjgAQ0ePFiPPfaYNW7cuHG65557VLFiRQUHB6tv3746c+aMdf/s2bNVuXJlffbZZ6pbt67c3d31j3/8Q2fPntWcOXNUo0YNValSRS+99JLNkwtr1KihN998U0899ZQ8PDwUFBSkd999t9Q+HzlyRHFxcapSpYp8fHwUGxurAwcOlBj/2Wef6eGHH9arr76qunXrqk6dOurQoYNNO1dO1bBYLEVeNWrUsO7fu3ev2rZtKw8PD/n7+ys+Pl7Hjx83ccWvHYkzAACAHV0eaV2yZIlyc3NLjKtQoYImTpyo7777TnPmzNHatWs1cOBAm5hz585p4sSJWrhwoVasWKH169erU6dOWr58uZYvX6558+bpvffe0//+9z+b495++201aNBAO3bs0ODBgzVgwACtWrWq2H6cO3dOzZs3l4eHhzZs2KBNmzbJw8NDMTExysvLK/aYgIAAff/99/ruu+9MX5fMzEzr66efftKdd96pJk2aWPc1bdpUDRs21LZt27RixQr9/vvv6tKli+n6rwVTNQAAt4yQLuPs3YXbxsGPkuzdBfz/HB0dNXv2bPXp00fTpk1To0aN1LRpUz355JNq0KCBNS4xMdH6/6GhoXrzzTf1/PPPa8qUKdbyixcvaurUqapVq5Yk6R//+IfmzZun33//XR4eHgoLC1Pz5s21bt06xcXFWY976KGH9Nprr0mS6tSpo6+++krjx49Xy5Yti/R34cKFqlChgv773/9an344a9YsVa5cWevXr1erVq2KHPPSSy9p48aNuueeexQSEqIHH3xQrVq1UteuXeXi4lLsdQkICJB06el/nTt3lpeXl6ZPny5Jmjp1qho1aqTRo0db42fOnKng4GD9+OOPqlOnTukX/Rox4gwAAGBnnTt31m+//aZly5apdevWWr9+vRo1aqTZs2dbY9atW6eWLVuqWrVq8vT0VLdu3XTixAmdPXvWGuPu7m5NmiXJ399fNWrUsJk77O/vr6NHj9q0HxUVVWQ7PT292L5u375dP/30kzw9Pa2j5d7e3rpw4YJ+/vnnYo+pWLGiPv/8c/30008aNmyYPDw89PLLL+uBBx7QuXPnSr02Q4YMUVpampYsWSI3NzdrH9atW2dt38PDQ3fddZckldiH8sCIMwAAwE3A1dVVLVu2VMuWLfX666/rmWee0YgRI9SjRw8dPHhQbdu2VUJCgt588015e3tr06ZN6t27ty5evGitw8nJyaZOi8VSbFlhYeFV+3N5NPlKhYWFCg8P1/z584vsq1q1aql11qpVS7Vq1dIzzzyjoUOHqk6dOkpJSVHPnj2Ljf/ggw80fvx4rV+/XnfccYdNH9q3b6/k5OQixwQGBpbah7+CxBkAAOAmFBYWZl07edu2bcrPz9c777xjfXz0Rx99VG5tff3110W2L4/gXqlRo0ZKSUmRn5+fKlWqdM1t1qhRQ+7u7jYj5n+WlpamZ555RtOnT9eDDz5YpA+LFi1SjRo15Oh449JZpmoAAADY0YkTJ/TII4/ogw8+0LfffquMjAx9/PHHGjt2rGJjYyVdGqnNz8/Xu+++q19++UXz5s3TtGnTyq0PX331lcaOHasff/xRkydP1scff6z+/fsXG9u1a1f5+voqNjZWGzduVEZGhlJTU9W/f3/9+uuvxR4zcuRIDRw4UOvXr1dGRoZ27typXr166eLFi8XOo87KylLHjh315JNPqnXr1srKylJWVpaOHTsmSXrhhRf0xx9/6J///Ke++eYb/fLLL1q5cqV69epls2JIeSNxBgAAsCMPDw9FRkZq/PjxatKkierXr6/hw4erT58+mjRpkiSpYcOGGjdunJKTk1W/fn3Nnz9fY8aMKbc+vPzyy9q+fbvuu+8+vfnmm3rnnXfUunXrYmPd3d21YcMGVa9eXZ06dVK9evXUq1cvnT9/vsQR6KZNm+qXX35Rt27ddNddd6lNmzbKysrSypUrVbdu3SLxP/zwg37//XfNmTNHgYGB1tf9998vSQoKCtJXX32lgoICtW7dWvXr11f//v3l5eVlHZG/HiyGYRjXrfa/uZycHHl5eSk7O/svfZQBALiEVTXKD6tq3FgXLlxQRkaGQkND5erqau/u2KhRo4YSExNtVu243ZR2/cuSrzHiDAAAAJhA4gwAAACYwKoaAAAAf2OlPSobtuw+4jxlyhTrfJPw8HBt3Lix1PjU1FSFh4fL1dVVNWvWLPYbpYsWLVJYWJhcXFwUFhamxYsX2+yfOnWqGjRooEqVKqlSpUqKiorSF198YRNjGIZGjhypoKAgubm5qVmzZvr+++//+gkDAADglmTXxDklJUWJiYkaOnSodu7cqcaNG6tNmzY6dOhQsfEZGRlq27atGjdurJ07d2rIkCHq16+fFi1aZI1JS0tTXFyc4uPjtXv3bsXHx6tLly7asmWLNeaOO+7Qv//9b23btk3btm3TI488otjYWJvEeOzYsRo3bpwmTZqkrVu3KiAgQC1bttTp06ev3wUBAADATcuuq2pERkaqUaNGmjp1qrWsXr166tChQ7FLrAwaNEjLli2zeQRkQkKCdu/erbS0NElSXFyccnJybEaQY2JiVKVKFS1YsKDEvnh7e+vtt99W7969ZRiGgoKClJiYqEGDBkmScnNz5e/vr+TkZD333HPF1pGbm6vc3Fzrdk5OjoKDg1lVAwDKCatqlB9W1bixbuZVNf4ObvlVNfLy8rR9+3a1atXKprxVq1bavHlzscekpaUViW/durW2bdtmfdxkSTEl1VlQUKCFCxfq7Nmz1ue0Z2RkKCsry6YeFxcXNW3atMR6JGnMmDHy8vKyvoKDg0uMBQAAwK3Fbonz8ePHVVBQIH9/f5tyf39/ZWVlFXtMVlZWsfH5+fk6fvx4qTFX1rlnzx55eHjIxcVFCQkJWrx4scLCwqx1XD7ObN8kafDgwcrOzra+Dh8+XGIsAAAAbi12X1XDYrHYbBuGUaTsavFXlpups27dutq1a5dOnTqlRYsWqXv37kpNTbUmz9fSNxcXF7m4uJS4HwAAALcuuyXOvr6+cnBwKDKCe/To0SIjvZcFBAQUG+/o6CgfH59SY66s09nZWXfeeackKSIiQlu3btV//vMfTZ8+XQEBAZIujTwHBgaa6hsAAEB5udHz+Znzbo7dpmo4OzsrPDxcq1atsilftWqVoqOjiz0mKiqqSPzKlSsVEREhJyenUmNKqvMywzCsX+wLDQ1VQECATT15eXlKTU29aj0AAAB/J2PGjJHFYrF5ZLeZZX1zc3P10ksvydfXVxUrVtTjjz+uX3/99Qb3vmzsuhxdUlKS/vvf/2rmzJlKT0/XgAEDdOjQISUkJEi6NGe4W7du1viEhAQdPHhQSUlJSk9P18yZMzVjxgy98sor1pj+/ftr5cqVSk5O1g8//KDk5GStXr3a5oc5ZMgQbdy4UQcOHNCePXs0dOhQrV+/Xl27dpUk6w9/9OjRWrx4sb777jv16NFD7u7ueuqpp27MxQEAALjJbd26Ve+9954aNGhgU25mWd/ExEQtXrxYCxcu1KZNm3TmzBm1a9dOBQUFN/o0TLPrHOe4uDidOHFCo0aNUmZmpurXr6/ly5crJCREkpSZmWmzpnNoaKiWL1+uAQMGaPLkyQoKCtLEiRPVuXNna0x0dLQWLlyoYcOGafjw4apVq5ZSUlIUGRlpjfn9998VHx+vzMxMeXl5qUGDBlqxYoVatmxpjRk4cKDOnz+vvn376uTJk4qMjNTKlSvl6el5A64MAADAze3MmTPq2rWr3n//fb311lvWcsMwNGHCBA0dOlSdOnWSJM2ZM0f+/v768MMP9dxzzyk7O1szZszQvHnz9Oijj0qSPvjgAwUHB2v16tVq3bq1Xc7pauy6jvPtrizrAgIAro51nMsPc1pvrLKu43wrzHHu3r27vL29NX78eDVr1kwNGzbUhAkT9Msvv6hWrVrasWOH7rvvPmt8bGysKleurDlz5mjt2rVq0aKF/vjjD1WpUsUac++996pDhw564403yuW8LiuvdZztvqoGAAAAbi0LFy7Ujh07tHXr1iL7SlvW9+DBg9YYZ2dnm6T5ckxpS//aG4kzAAAATDt8+LD1O2WljZ6XdVlfszH2ZNcvBwIAAODWsn37dh09elTh4eFydHSUo6OjUlNTNXHiRDk6OlpHmktbHjggIEB5eXk6efJkiTE3IxJnAAAAmNaiRQvt2bNHu3btsr4iIiLUtWtX7dq1SzVr1rzqsr7h4eFycnKyicnMzNR33313Uy/9y1QNAAAAmObp6an69evblFWsWFE+Pj7W8svL+tauXVu1a9fW6NGjbZb19fLyUu/evfXyyy/Lx8dH3t7eeuWVV3TPPfdYV9m4GZE4AwAA3GRu9VVPzCzrO378eDk6OqpLly46f/68WrRoodmzZ8vBwcGOPS8dy9FdRyxHBwDli+Xoys+tnpjdasq6HB3KV3ktR8ccZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEHoByG2O90/LDeqcAAIARZwAAAMAERpwBAABuMo0Hp9zQ9jaOibuh7d2qGHEGAABAmWzYsEHt27dXUFCQLBaLlixZUiQmPT1djz/+uLy8vOTp6akHH3xQhw4dsu7Pzc3VSy+9JF9fX1WsWFGPP/64fv31V5s6Tp48qfj4eHl5ecnLy0vx8fE6derUdT67kpE4AwAAoEzOnj2re++9V5MmTSp2/88//6yHH35Yd911l9avX6/du3dr+PDhcnV1tcYkJiZq8eLFWrhwoTZt2qQzZ86oXbt2KigosMY89dRT2rVrl1asWKEVK1Zo165dio+Pv+7nVxKmagAAAKBM2rRpozZt2pS4f+jQoWrbtq3Gjh1rLatZs6b1/7OzszVjxgzNmzdPjz76qCTpgw8+UHBwsFavXq3WrVsrPT1dK1as0Ndff63IyEhJ0vvvv6+oqCjt27dPdevWvU5nVzJGnAEAAFBuCgsL9fnnn6tOnTpq3bq1/Pz8FBkZaTOdY/v27bp48aJatWplLQsKClL9+vW1efNmSVJaWpq8vLysSbMkPfjgg/Ly8rLG3GgkzgAAACg3R48e1ZkzZ/Tvf/9bMTExWrlypTp27KhOnTopNTVVkpSVlSVnZ2dVqVLF5lh/f39lZWVZY/z8/IrU7+fnZ4250ZiqAQAAgHJTWFgoSYqNjdWAAQMkSQ0bNtTmzZs1bdo0NW3atMRjDcOQxWKxbv/5/0uKuZEYcQYAAEC58fX1laOjo8LCwmzK69WrZ11VIyAgQHl5eTp58qRNzNGjR+Xv72+N+f3334vUf+zYMWvMjUbiDAAAgHLj7Oys+++/X/v27bMp//HHHxUSEiJJCg8Pl5OTk1atWmXdn5mZqe+++07R0dGSpKioKGVnZ+ubb76xxmzZskXZ2dnWmBuNqRoAAAAokzNnzuinn36ybmdkZGjXrl3y9vZW9erV9eqrryouLk5NmjRR8+bNtWLFCn366adav369JMnLy0u9e/fWyy+/LB8fH3l7e+uVV17RPffcY11lo169eoqJiVGfPn00ffp0SdKzzz6rdu3a2WVFDYnEGQAA4KZzsz/Jb9u2bWrevLl1OykpSZLUvXt3zZ49Wx07dtS0adM0ZswY9evXT3Xr1tWiRYv08MMPW48ZP368HB0d1aVLF50/f14tWrTQ7Nmz5eDgYI2ZP3+++vXrZ1194/HHHy9x7egbwWIYhmG31m9zOTk58vLyUnZ2tipVqnTD2w/pMu6Gt3m7OvhRkr27AEDc18oT97Ub68KFC8rIyFBoaKjNQ0BwY5R2/cuSrzHHGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAwgcQZAAAAMIFHbgMAANxker+/7oa2N6NP86sHgRFnAAAAlN2RI0f09NNPy8fHR+7u7mrYsKG2b99ebOxzzz0ni8WiCRMm2JTn5ubqpZdekq+vrypWrKjHH39cv/76q03MyZMnFR8fLy8vL3l5eSk+Pl6nTp26TmdVOhJnAAAAlMnJkyf10EMPycnJSV988YX27t2rd955R5UrVy4Su2TJEm3ZskVBQUFF9iUmJmrx4sVauHChNm3apDNnzqhdu3YqKCiwxjz11FPatWuXVqxYoRUrVmjXrl2Kj4+/nqdXIqZqAAAAoEySk5MVHBysWbNmWctq1KhRJO7IkSN68cUX9eWXX+qxxx6z2Zedna0ZM2Zo3rx5evTRRyVJH3zwgYKDg7V69Wq1bt1a6enpWrFihb7++mtFRkZKkt5//31FRUVp3759qlu37vU7yWIw4gwAAIAyWbZsmSIiIvTEE0/Iz89P9913n95//32bmMLCQsXHx+vVV1/V3XffXaSO7du36+LFi2rVqpW1LCgoSPXr19fmzZslSWlpafLy8rImzZL04IMPysvLyxpzI5E4AwAAoEx++eUXTZ06VbVr19aXX36phIQE9evXT3PnzrXGJCcny9HRUf369Su2jqysLDk7O6tKlSo25f7+/srKyrLG+Pn5FTnWz8/PGnMjMVUDAAAAZVJYWKiIiAiNHj1aknTffffp+++/19SpU9WtWzdt375d//nPf7Rjxw5ZLJYy1W0Yhs0xxR1/ZcyNwogzAAAAyiQwMFBhYWE2ZfXq1dOhQ4ckSRs3btTRo0dVvXp1OTo6ytHRUQcPHtTLL79snQsdEBCgvLw8nTx50qaeo0ePyt/f3xrz+++/F2n/2LFj1pgbicQZAAAAZfLQQw9p3759NmU//vijQkJCJEnx8fH69ttvtWvXLusrKChIr776qr788ktJUnh4uJycnLRq1SprHZmZmfruu+8UHR0tSYqKilJ2dra++eYba8yWLVuUnZ1tjbmRmKoBAACAMhkwYICio6M1evRodenSRd98843ee+89vffee5IkHx8f+fj42Bzj5OSkgIAA60oYXl5e6t27t15++WX5+PjI29tbr7zyiu655x7rKhv16tVTTEyM+vTpo+nTp0uSnn32WbVr1+6Gr6gh3QSJ85QpU/T2228rMzNTd999tyZMmKDGjRuXGJ+amqqkpCR9//33CgoK0sCBA5WQkGATs2jRIg0fPlw///yzatWqpX/961/q2LGjdf+YMWP0ySef6IcffpCbm5uio6OVnJxs8wPo0aOH5syZY1NvZGSkvv7663I6cwAAgOLd7E/yu//++7V48WINHjxYo0aNUmhoqCZMmKCuXbuWqZ7x48fL0dFRXbp00fnz59WiRQvNnj1bDg4O1pj58+erX79+1tU3Hn/8cU2aNKlcz8csuybOKSkpSkxM1JQpU/TQQw9p+vTpatOmjfbu3avq1asXic/IyFDbtm3Vp08fffDBB/rqq6/Ut29fVa1aVZ07d5Z0admSuLg4vfnmm+rYsaMWL16sLl26aNOmTdalTFJTU/XCCy/o/vvvV35+voYOHapWrVpp7969qlixorW9mJgYm/UJnZ2dr/MVAQAAuDW0a9dO7dq1Mx1/4MCBImWurq5699139e6775Z4nLe3tz744INr6WK5sxiGYdir8cjISDVq1EhTp061ltWrV08dOnTQmDFjisQPGjRIy5YtU3p6urUsISFBu3fvVlpamiQpLi5OOTk5+uKLL6wxMTExqlKlihYsWFBsP44dOyY/Pz+lpqaqSZMmki6NOJ86dUpLliy55vPLycmRl5eXsrOzValSpWuu51qFdBl3w9u8XR38KMneXQAg7mvlifvajXXhwgVlZGQoNDRUrq6u9u7O305p178s+ZrdvhyYl5en7du32yx6LUmtWrUqcUHrtLS0IvGtW7fWtm3bdPHixVJjSlskOzs7W9Klf9H82fr16+Xn56c6deqoT58+Onr0aKnnlJubq5ycHJsXAAAAbg92S5yPHz+ugoKCIkuJ/HnR6ytlZWUVG5+fn6/jx4+XGlNSnYZhKCkpSQ8//LDq169vLW/Tpo3mz5+vtWvX6p133tHWrVv1yCOPKDc3t8RzGjNmjLy8vKyv4ODgki8AAAAAbil2/3LglYtXX21B6+LirywvS50vvviivv32W23atMmmPC4uzvr/9evXV0REhEJCQvT555+rU6dOxdY1ePBgJSX930dfOTk5JM8AAAC3Cbslzr6+vnJwcCgyEvznRa+vFBAQUGy8o6OjdcmTkmKKq/Oll17SsmXLtGHDBt1xxx2l9jcwMFAhISHav39/iTEuLi5ycXEptR4AAPD3Zcevlv2tldd1t9tUDWdnZ4WHh9ssei1Jq1atKnFB66ioqCLxK1euVEREhJycnEqN+XOdhmHoxRdf1CeffKK1a9cqNDT0qv09ceKEDh8+rMDAQFPnBwAAcNnlPOXcuXN27snf0+XrfvnncK3sOlUjKSlJ8fHxioiIUFRUlN577z0dOnTIui7z4MGDdeTIEc2dO1fSpRU0Jk2apKSkJPXp00dpaWmaMWOGzWoZ/fv3V5MmTZScnKzY2FgtXbpUq1evtpmK8cILL+jDDz/U0qVL5enpaR2h9vLykpubm86cOaORI0eqc+fOCgwM1IEDBzRkyBD5+vrarAcNAABghoODgypXrmxdaMDd3b3UqakoH4Zh6Ny5czp69KgqV65ssz70tbBr4hwXF6cTJ05o1KhRyszMVP369bV8+XLr4xozMzOtzzyXpNDQUC1fvlwDBgzQ5MmTFRQUpIkTJ1rXcJak6OhoLVy4UMOGDdPw4cNVq1YtpaSkWNdwlmRd/q5Zs2Y2/Zk1a5Z69OghBwcH7dmzR3PnztWpU6cUGBio5s2bKyUlRZ6entfxigAAgNtVQECAJF11lS6Uv8qVK1uv/19h13Wcb3es43z7YL1T4ObAfa38cF+zn4KCAusyurj+nJycSh1pLku+ZvdVNQAAAP5OHBwc/vKUAdiH3b4cCAAAANxKSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE+yeOE+ZMkWhoaFydXVVeHi4Nm7cWGp8amqqwsPD5erqqpo1a2ratGlFYhYtWqSwsDC5uLgoLCxMixcvttk/ZswY3X///fL09JSfn586dOigffv22cQYhqGRI0cqKChIbm5uatasmb7//vu/fsIAAAC4Jdk1cU5JSVFiYqKGDh2qnTt3qnHjxmrTpo0OHTpUbHxGRobatm2rxo0ba+fOnRoyZIj69eunRYsWWWPS0tIUFxen+Ph47d69W/Hx8erSpYu2bNlijUlNTdULL7ygr7/+WqtWrVJ+fr5atWqls2fPWmPGjh2rcePGadKkSdq6dasCAgLUsmVLnT59+vpdEAAAANy0LIZhGPZqPDIyUo0aNdLUqVOtZfXq1VOHDh00ZsyYIvGDBg3SsmXLlJ6ebi1LSEjQ7t27lZaWJkmKi4tTTk6OvvjiC2tMTEyMqlSpogULFhTbj2PHjsnPz0+pqalq0qSJDMNQUFCQEhMTNWjQIElSbm6u/P39lZycrOeee87U+eXk5MjLy0vZ2dmqVKmSqWPKU0iXcTe8zdvVwY+S7N0FAOK+Vp64rwGXlCVfs9uIc15enrZv365WrVrZlLdq1UqbN28u9pi0tLQi8a1bt9a2bdt08eLFUmNKqlOSsrOzJUne3t6SLo1sZ2Vl2dTj4uKipk2bllpPbm6ucnJybF4AAAC4PdgtcT5+/LgKCgrk7+9vU+7v76+srKxij8nKyio2Pj8/X8ePHy81pqQ6DcNQUlKSHn74YdWvX99ax+XjzNYjXZo77eXlZX0FBweXGAsAAIBbi92/HGixWGy2DcMoUna1+CvLy1Lniy++qG+//bbYaRxl7dvgwYOVnZ1tfR0+fLjEWAAAANxaHO3VsK+vrxwcHIqM4B49erTISO9lAQEBxcY7OjrKx8en1Jji6nzppZe0bNkybdiwQXfccYdNO9KlkefAwEBTfZMuTedwcXEpcT8AAABuXXYbcXZ2dlZ4eLhWrVplU75q1SpFR0cXe0xUVFSR+JUrVyoiIkJOTk6lxvy5TsMw9OKLL+qTTz7R2rVrFRoaahMfGhqqgIAAm3ry8vKUmppaYt8AAABwe7PbiLMkJSUlKT4+XhEREYqKitJ7772nQ4cOKSEhQdKlqQ9HjhzR3LlzJV1aQWPSpElKSkpSnz59lJaWphkzZthMs+jfv7+aNGmi5ORkxcbGaunSpVq9erU2bdpkjXnhhRf04YcfaunSpfL09LSOUHt5ecnNzU0Wi0WJiYkaPXq0ateurdq1a2v06NFyd3fXU089dQOvEAAAAG4Wdk2c4+LidOLECY0aNUqZmZmqX7++li9frpCQEElSZmamzZrOoaGhWr58uQYMGKDJkycrKChIEydOVOfOna0x0dHRWrhwoYYNG6bhw4erVq1aSklJUWRkpDXm8vJ3zZo1s+nPrFmz1KNHD0nSwIEDdf78efXt21cnT55UZGSkVq5cKU9Pz+t0NQAAAHAzs+s6zrc71nG+fbDeKXBz4L5WfrivAZfcEus4AwAAALcSEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADAhDIlzmPHjtX58+et2xs2bFBubq51+/Tp0+rbt2/59Q4AAAC4SZQpcR48eLBOnz5t3W7Xrp2OHDli3T537pymT59efr0DAAAAbhJlSpwNwyh1GwAAALhdMccZAAAAMIHEGQAAADDBsawH/Pe//5WHh4ckKT8/X7Nnz5avr68k2cx/BgAAAG4nZUqcq1evrvfff9+6HRAQoHnz5hWJAQAAAG43ZUqcDxw4cJ26AQAAANzcmOMMAAAAmFCmxHnLli364osvbMrmzp2r0NBQ+fn56dlnn7V5IAoAAABwuyhT4jxy5Eh9++231u09e/aod+/eevTRR/Xaa6/p008/1ZgxY8q9kwAAAIC9lSlx3rVrl1q0aGHdXrhwoSIjI/X+++8rKSlJEydO1EcffVTunQQAAADsrUyJ88mTJ+Xv72/dTk1NVUxMjHX7/vvv1+HDh8uvdwAAAMBNokyJs7+/vzIyMiRJeXl52rFjh6Kioqz7T58+LScnp/LtIQAAAHATKFPiHBMTo9dee00bN27U4MGD5e7ursaNG1v3f/vtt6pVq1a5dxIAAACwtzKt4/zWW2+pU6dOatq0qTw8PDR79mw5Oztb98+cOVOtWrUq904CAAAA9lamxLlq1arauHGjsrOz5eHhIQcHB5v9H3/8sTw9Pcu1gwAAAMDNoEyJc69evUzFzZw585o6AwAAANysypQ4z549WyEhIbrvvvtkGMb16hMAAABw0ylT4pyQkKCFCxfql19+Ua9evfT000/L29v7evUNAAAAuGmUaVWNKVOmKDMzU4MGDdKnn36q4OBgdenSRV9++SUj0AAAALitlSlxliQXFxf985//1KpVq7R3717dfffd6tu3r0JCQnTmzJnr0UcAAADA7sqcOP+ZxWKRxWKRYRgqLCy8pjqmTJmi0NBQubq6Kjw8XBs3biw1PjU1VeHh4XJ1dVXNmjU1bdq0IjGLFi1SWFiYXFxcFBYWpsWLF9vs37Bhg9q3b6+goCBZLBYtWbKkSB09evSwnt/l14MPPnhN5wgAAIBbX5kT59zcXC1YsEAtW7ZU3bp1tWfPHk2aNEmHDh2Sh4dHmepKSUlRYmKihg4dqp07d6px48Zq06aNDh06VGx8RkaG2rZtq8aNG2vnzp0aMmSI+vXrp0WLFllj0tLSFBcXp/j4eO3evVvx8fHq0qWLtmzZYo05e/as7r33Xk2aNKnU/sXExCgzM9P6Wr58eZnODwAAALcPi1GGycl9+/bVwoULVb16dfXs2VNPP/20fHx8rrnxyMhINWrUSFOnTrWW1atXTx06dNCYMWOKxA8aNEjLli1Tenq6tSwhIUG7d+9WWlqaJCkuLk45OTn64osvrDExMTGqUqWKFixYUKROi8WixYsXq0OHDjblPXr00KlTp4odjTYrJydHXl5eys7OVqVKla65nmsV0mXcDW/zdnXwoyR7dwGAuK+VJ+5rwCVlydfKtKrGtGnTVL16dYWGhio1NVWpqanFxn3yySdXrSsvL0/bt2/Xa6+9ZlPeqlUrbd68udhj0tLSijyZsHXr1poxY4YuXrwoJycnpaWlacCAAUViJkyYcNU+XWn9+vXy8/NT5cqV1bRpU/3rX/+Sn59fifG5ubnKzc21bufk5JS5TQAAANycypQ4d+vWTRaLpVwaPn78uAoKCuTv729T7u/vr6ysrGKPycrKKjY+Pz9fx48fV2BgYIkxJdVZkjZt2uiJJ55QSEiIMjIyNHz4cD3yyCPavn27XFxcij1mzJgxeuONN8rUDgAAAG4NZX4ASnm7MhE3DKPU5Ly4+CvLy1pnceLi4qz/X79+fUVERCgkJESff/65OnXqVOwxgwcPVlLS/330lZOTo+Dg4DK1CwAAgJtTmRLn8uTr6ysHB4ciI8FHjx4tMmJ8WUBAQLHxjo6O1rnWJcWUVKdZgYGBCgkJ0f79+0uMcXFxKXE0GgAAALe2v7Qc3V/h7Oys8PBwrVq1yqZ81apVio6OLvaYqKioIvErV65URESEnJycSo0pqU6zTpw4ocOHDyswMPAv1QMAAIBbk91GnCUpKSlJ8fHxioiIUFRUlN577z0dOnRICQkJki5NfThy5Ijmzp0r6dIKGpMmTVJSUpL69OmjtLQ0zZgxw2a1jP79+6tJkyZKTk5WbGysli5dqtWrV2vTpk3WmDNnzuinn36ybmdkZGjXrl3y9vZW9erVdebMGY0cOVKdO3dWYGCgDhw4oCFDhsjX11cdO3a8QVcHAAAANxO7Js5xcXE6ceKERo0apczMTNWvX1/Lly9XSEiIJCkzM9NmTefQ0FAtX75cAwYM0OTJkxUUFKSJEyeqc+fO1pjo6GgtXLhQw4YN0/Dhw1WrVi2lpKQoMjLSGrNt2zY1b97cun15XnL37t01e/ZsOTg4aM+ePZo7d65OnTqlwMBANW/eXCkpKfL09LzelwUAAAA3oTKt44yyYR3n2wfrnQI3B+5r5Yf7GnBJWfI1u81xBgAAAG4lJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAl2T5ynTJmi0NBQubq6Kjw8XBs3biw1PjU1VeHh4XJ1dVXNmjU1bdq0IjGLFi1SWFiYXFxcFBYWpsWLF9vs37Bhg9q3b6+goCBZLBYtWbKkSB2GYWjkyJEKCgqSm5ubmjVrpu+///4vnSsAAABuXXZNnFNSUpSYmKihQ4dq586daty4sdq0aaNDhw4VG5+RkaG2bduqcePG2rlzp4YMGaJ+/fpp0aJF1pi0tDTFxcUpPj5eu3fvVnx8vLp06aItW7ZYY86ePat7771XkyZNKrFvY8eO1bhx4zRp0iRt3bpVAQEBatmypU6fPl1+FwAAAAC3DIthGIa9Go+MjFSjRo00depUa1m9evXUoUMHjRkzpkj8oEGDtGzZMqWnp1vLEhIStHv3bqWlpUmS4uLilJOToy+++MIaExMToypVqmjBggVF6rRYLFq8eLE6dOhgLTMMQ0FBQUpMTNSgQYMkSbm5ufL391dycrKee+45U+eXk5MjLy8vZWdnq1KlSqaOKU8hXcbd8DZvVwc/SrJ3FwCI+1p54r4GXFKWfM1uI855eXnavn27WrVqZVPeqlUrbd68udhj0tLSisS3bt1a27Zt08WLF0uNKanO4mRkZCgrK8umHhcXFzVt2rTUenJzc5WTk2PzAgAAwO3Bbonz8ePHVVBQIH9/f5tyf39/ZWVlFXtMVlZWsfH5+fk6fvx4qTEl1VlSO5ePK0s9Y8aMkZeXl/UVHBxsuk0AAADc3Oz+5UCLxWKzbRhGkbKrxV9ZXtY6y6tvgwcPVnZ2tvV1+PDhMrcJAACAm5OjvRr29fWVg4NDkRHco0ePFhnpvSwgIKDYeEdHR/n4+JQaU1KdJbUjXRp5DgwMNF2Pi4uLXFxcTLcDAACAW4fdRpydnZ0VHh6uVatW2ZSvWrVK0dHRxR4TFRVVJH7lypWKiIiQk5NTqTEl1Vmc0NBQBQQE2NSTl5en1NTUMtUDAACA24fdRpwlKSkpSfHx8YqIiFBUVJTee+89HTp0SAkJCZIuTX04cuSI5s6dK+nSChqTJk1SUlKS+vTpo7S0NM2YMcNmtYz+/furSZMmSk5OVmxsrJYuXarVq1dr06ZN1pgzZ87op59+sm5nZGRo165d8vb2VvXq1WWxWJSYmKjRo0erdu3aql27tkaPHi13d3c99dRTN+jqAAAA4GZi18Q5Li5OJ06c0KhRo5SZman69etr+fLlCgkJkSRlZmbarOkcGhqq5cuXa8CAAZo8ebKCgoI0ceJEde7c2RoTHR2thQsXatiwYRo+fLhq1aqllJQURUZGWmO2bdum5s2bW7eTki4tydO9e3fNnj1bkjRw4ECdP39effv21cmTJxUZGamVK1fK09Pzel4SAAAA3KTsuo7z7Y51nG8frHcK3By4r5Uf7mvAJbfEOs4AAADArYTEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADDB0d4dAG4FjQen2LsLt4WNY+Ls3QUAAK4ZI84AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJJM4AAACACSTOAAAAgAkkzgAAAIAJdk+cp0yZotDQULm6uio8PFwbN24sNT41NVXh4eFydXVVzZo1NW3atCIxixYtUlhYmFxcXBQWFqbFixeXud0ePXrIYrHYvB588MG/drIAAAC4Zdk1cU5JSVFiYqKGDh2qnTt3qnHjxmrTpo0OHTpUbHxGRobatm2rxo0ba+fOnRoyZIj69eunRYsWWWPS0tIUFxen+Ph47d69W/Hx8erSpYu2bNlS5nZjYmKUmZlpfS1fvvz6XAgAAADc9CyGYRj2ajwyMlKNGjXS1KlTrWX16tVThw4dNGbMmCLxgwYN0rJly5Senm4tS0hI0O7du5WWliZJiouLU05Ojr744gtrTExMjKpUqaIFCxaYbrdHjx46deqUlixZcs3nl5OTIy8vL2VnZ6tSpUrXXM+1Cuky7oa3ebuqXquavbtwW9g4Js7eXcAtjvta+Tn4UZK9uwDcFMqSr9ltxDkvL0/bt29Xq1atbMpbtWqlzZs3F3tMWlpakfjWrVtr27ZtunjxYqkxl+ssS7vr16+Xn5+f6tSpoz59+ujo0aOlnlNubq5ycnJsXgAAALg92C1xPn78uAoKCuTv729T7u/vr6ysrGKPycrKKjY+Pz9fx48fLzXmcp1m223Tpo3mz5+vtWvX6p133tHWrVv1yCOPKDc3t8RzGjNmjLy8vKyv4ODgq1wFAAAA3Coc7d0Bi8Vis20YRpGyq8VfWW6mzqvFxMX930fK9evXV0REhEJCQvT555+rU6dOxfZt8ODBSkr6v4++cnJySJ4BAABuE3ZLnH19feXg4FBkdPno0aNFRoMvCwgIKDbe0dFRPj4+pcZcrvNa2pWkwMBAhYSEaP/+/SXGuLi4yMXFpcT9AAAAuHXZbaqGs7OzwsPDtWrVKpvyVatWKTo6uthjoqKiisSvXLlSERERcnJyKjXmcp3X0q4knThxQocPH1ZgYKC5EwQAAMBtxa5TNZKSkhQfH6+IiAhFRUXpvffe06FDh5SQkCDp0tSHI0eOaO7cuZIuraAxadIkJSUlqU+fPkpLS9OMGTOsq2VIUv/+/dWkSRMlJycrNjZWS5cu1erVq7Vp0ybT7Z45c0YjR45U586dFRgYqAMHDmjIkCHy9fVVx44db+AVAgAAwM3CrolzXFycTpw4oVGjRikzM1P169fX8uXLFRISIknKzMy0WVs5NDRUy5cv14ABAzR58mQFBQVp4sSJ6ty5szUmOjpaCxcu1LBhwzR8+HDVqlVLKSkpioyMNN2ug4OD9uzZo7lz5+rUqVMKDAxU8+bNlZKSIk9Pzxt0dQAAAHAzses6zrc71nG+fbCOc/lgHWf8VdzXyg/rOAOX3BLrOAMAAAC3EhJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAEEmcAAADABBJnAAAAwAQSZwAAAMAER3t3AAAA3HiNB6fYuwu3hY1j4uzdBdxAjDgDAAAAJpA4AwAAACaQOAMAAAAmkDgDAAAAJpA4AwAAACaQOAMAAAAmkDgDAAAAJpA4AwAAACaQOAMAAAAmkDgDAAAAJpA4AwAAACY42rsDAP4+er+/zt5duG3M6NPc3l0AgL8dRpwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE0icAQAAABNInAEAAAATSJwBAAAAE3jkNgAAwDXq/f46e3fhtjGjT3N7d+Gq7D7iPGXKFIWGhsrV1VXh4eHauHFjqfGpqakKDw+Xq6uratasqWnTphWJWbRokcLCwuTi4qKwsDAtXry4zO0ahqGRI0cqKChIbm5uatasmb7//vu/drIAAAC4Zdk1cU5JSVFiYqKGDh2qnTt3qnHjxmrTpo0OHTpUbHxGRobatm2rxo0ba+fOnRoyZIj69eunRYsWWWPS0tIUFxen+Ph47d69W/Hx8erSpYu2bNlSpnbHjh2rcePGadKkSdq6dasCAgLUsmVLnT59+vpdEAAAANy0LIZhGPZqPDIyUo0aNdLUqVOtZfXq1VOHDh00ZsyYIvGDBg3SsmXLlJ6ebi1LSEjQ7t27lZaWJkmKi4tTTk6OvvjiC2tMTEyMqlSpogULFphq1zAMBQUFKTExUYMGDZIk5ebmyt/fX8nJyXruuedMnV9OTo68vLyUnZ2tSpUqleHKlI+QLuNueJu3q+q1qtm7C7eFOjX97N2F28at8JHm9cB9rfxwXysf3NfKj73ua2XJ1+w2xzkvL0/bt2/Xa6+9ZlPeqlUrbd68udhj0tLS1KpVK5uy1q1ba8aMGbp48aKcnJyUlpamAQMGFImZMGGC6XYzMjKUlZVl05aLi4uaNm2qzZs3l5g45+bmKjc317qdnZ0t6dIPxB4KL16wS7u3o/zcc/buwm0h7/xZe3fhtmGv+4q9cV8rP9zXygf3tfJjr/va5XbNjCXbLXE+fvy4CgoK5O/vb1Pu7++vrKysYo/JysoqNj4/P1/Hjx9XYGBgiTGX6zTT7uX/Fhdz8ODBEs9pzJgxeuONN4qUBwcHl3gMbg2/2rsDt4mv7d2B28gH/e3dA9zquK+VD+5r5cfe97XTp0/Ly8ur1Bi7r6phsVhstg3DKFJ2tfgry83UWV4xfzZ48GAlJSVZtwsLC/XHH3/Ix8en1OOAvyonJ0fBwcE6fPiwXaYFAUB5476GG8UwDJ0+fVpBQUFXjbVb4uzr6ysHB4cio8tHjx4tMtJ7WUBAQLHxjo6O8vHxKTXmcp1m2g0ICJB0aeQ5MDDQVN+kS9M5XFxcbMoqV65cYjxQ3ipVqsQfGAC3Fe5ruBGuNtJ8md1W1XB2dlZ4eLhWrVplU75q1SpFR0cXe0xUVFSR+JUrVyoiIkJOTk6lxlyu00y7oaGhCggIsInJy8tTampqiX0DAADA7c2uUzWSkpIUHx+viIgIRUVF6b333tOhQ4eUkJAg6dLUhyNHjmju3LmSLq2gMWnSJCUlJalPnz5KS0vTjBkzrKtlSFL//v3VpEkTJScnKzY2VkuXLtXq1au1adMm0+1aLBYlJiZq9OjRql27tmrXrq3Ro0fL3d1dTz311A28QgAAALhZ2DVxjouL04kTJzRq1ChlZmaqfv36Wr58uUJCQiRJmZmZNmsrh4aGavny5RowYIAmT56soKAgTZw4UZ07d7bGREdHa+HChRo2bJiGDx+uWrVqKSUlRZGRkabblaSBAwfq/Pnz6tu3r06ePKnIyEitXLlSnp6eN+DKAGXj4uKiESNGFJkqBAC3Ku5ruBnZdR1nAAAA4FZh90duAwAAALcCEmcAAADABBJnAAAAwAQSZ9zWmjVrpsTExFJjatSoYX0k+63ian0+cOCALBaLdu3adcP6dLMx87MH/q6u5z3ir9bN/euSst7DZs+ezbMjbgASZ9xQPXr0kMVikcVikZOTk/z9/dWyZUvNnDlThYWF9u7eNbt8Tl9/bfvw1dzcXOuTI9evX3/D+hMcHGxdMeZ6GjlypBo2bHhd27hWn3zyid588017dwOQdOkBWs8995yqV68uFxcXBQQEqHXr1kpLS7PGWCwWLVmyxH6dvIFKSwpv1P3rWl1O7B0dHXXkyBGbfZmZmXJ0dJTFYtGBAwfs00FcVyTOuOFiYmKUmZmpAwcO6IsvvlDz5s3Vv39/tWvXTvn5+fbu3jULDg7WrFmzbMoWL14sDw+PG94XBwcHBQQEyNHRritOXhcXL140Feft7c3ykbhpdO7cWbt379acOXP0448/atmyZWrWrJn++OMPe3ftmuXl5V2Xem+W+1dBQUGpAzpBQUHW50xcNmfOHFWrVu16dw12ROKMG+7yaEu1atXUqFEjDRkyREuXLtUXX3yh2bNnW+MOHTqk2NhYeXh4qFKlSurSpYt+//136/4ePXqoQ4cONnUnJiaqWbNmNmX5+fl68cUXVblyZfn4+GjYsGEqbRXG7OxsPfvss/Lz81OlSpX0yCOPaPfu3Vc9r+7du2vhwoU6f/68tWzmzJnq3r17kdhBgwapTp06cnd3V82aNTV8+PAiCeGyZcsUEREhV1dX+fr6qlOnTjb7z507p169esnT01PVq1fXe++9Z9135Ued69evl8Vi0Zo1axQRESF3d3dFR0dr3759NnV++umnCg8Pl6urq2rWrKk33njjL/1j5siRI4qLi1OVKlXk4+Oj2NhYm1GYrVu3qmXLlvL19ZWXl5eaNm2qHTt22NRhsVg0bdo0xcbGqmLFinrrrbesI93z5s1TjRo15OXlpSeffFKnT5+2HnfliFaNGjU0evToEq+ZJG3evFkNGzaUq6urIiIitGTJEj4yxl926tQpbdq0ScnJyWrevLlCQkL0wAMPaPDgwXrsscckXXp/SlLHjh1lsVis2z///LNiY2Pl7+8vDw8P3X///Vq9erVN/Wbe2998843uu+8+63t7586dNvsLCgrUu3dvhYaGys3NTXXr1tV//vMfm5jL99wxY8YoKChIderUMVV3WV2v+9e4ceN0zz33qGLFigoODlbfvn115swZ6/7LUx0+++wzhYWFycXFRQcPHiyxn927dy8yWDJ79uxi7/mpqal64IEH5OLiosDAQL322ms2fTt79qy6desmDw8PBQYG6p133ilSR15engYOHKhq1aqpYsWKioyMvKGfZOISEmfcFB555BHde++9+uSTTyRJhmGoQ4cO+uOPP5SamqpVq1bp559/VlxcXJnrnjNnjhwdHbVlyxZNnDhR48eP13//+99iYw3D0GOPPaasrCwtX75c27dvV6NGjdSiRYurjgyFh4crNDRUixYtkiQdPnxYGzZsUHx8fJFYT09PzZ49W3v37tV//vMfvf/++xo/frx1/+eff65OnTrpscce086dO61/MP7snXfesf6R6tu3r55//nn98MMPpfZx6NCheuedd7Rt2zY5OjqqV69e1n1ffvmlnn76afXr10979+7V9OnTNXv2bP3rX/8qtc6SnDt3Ts2bN5eHh4c2bNigTZs2ycPDQzExMdaRqtOnT6t79+7auHGjvv76a9WuXVtt27a1SYAlacSIEYqNjdWePXusff7555+1ZMkSffbZZ/rss8+Umpqqf//736X2qbRrdvr0abVv31733HOPduzYoTfffFODBg26pnMH/szDw0MeHh5asmSJcnNzi43ZunWrJGnWrFnKzMy0bp85c0Zt27bV6tWrtXPnTrVu3Vrt27e3eTiYVPp7++zZs2rXrp3q1q2r7du3a+TIkXrllVdsji8sLNQdd9yhjz76SHv37tXrr7+uIUOG6KOPPrKJW7NmjdLT07Vq1Sp99tlnpuouL3/1/lWhQgVNnDhR3333nebMmaO1a9dq4MCBNm2cO3dOY8aM0X//+199//338vPzK7E/jz/+uE6ePGl9MvGmTZv0xx9/qH379jZxR44cUdu2bXX//fdr9+7dmjp1qmbMmKG33nrLGvPqq69q3bp1Wrx4sVauXKn169dr+/btNvX07NlTX331lRYuXKhvv/1WTzzxhGJiYrR///6yX0xcOwO4gbp3727ExsYWuy8uLs6oV6+eYRiGsXLlSsPBwcE4dOiQdf/3339vSDK++eabEuvq37+/0bRpU+t206ZNjXr16hmFhYXWskGDBlnbMQzDCAkJMcaPH28YhmGsWbPGqFSpknHhwgWbemvVqmVMnz69xPOSZCxevNiYMGGC0bx5c8MwDOONN94wOnbsaJw8edKQZKxbt67E48eOHWuEh4dbt6OiooyuXbuWGB8SEmI8/fTT1u3CwkLDz8/PmDp1qmEYhpGRkWFIMnbu3GkYhmGsW7fOkGSsXr3aesznn39uSDLOnz9vGIZhNG7c2Bg9erRNO/PmzTMCAwNL7MeIESOMe++9t9h9M2bMMOrWrWtz7XNzcw03Nzfjyy+/LPaY/Px8w9PT0/j000+tZZKMxMTEIu26u7sbOTk51rJXX33ViIyMtG43bdrU6N+/v3X7atds6tSpho+Pj/V6GIZhvP/++zbXEbhW//vf/4wqVaoYrq6uRnR0tDF48GBj9+7dNjGX7yNXExYWZrz77rvW7au9t6dPn254e3sbZ8+etcZMnTr1qu/tvn37Gp07d7Zud+/e3fD39zdyc3OtZdda95W/n392o+5fH330keHj42PdnjVrliHJ2LVrV4nHXNm/xMREo2fPnoZhGEbPnj2NAQMGGDt37jQkGRkZGYZhGMaQIUOK3AsnT55seHh4GAUFBcbp06cNZ2dnY+HChdb9J06cMNzc3KzX6KeffjIsFotx5MgRm760aNHCGDx4sLX/Xl5epfYdfx0jzrhpGIYhi8UiSUpPT1dwcLCCg4Ot+8PCwlS5cmWlp6eXqd4HH3zQWq8kRUVFaf/+/SooKCgSu337dp05c0Y+Pj7WUSIPDw9lZGTo559/vmpbTz/9tNLS0vTLL79o9uzZNiMif/a///1PDz/8sAICAuTh4aHhw4fbjCDt2rVLLVq0KLWtBg0aWP/fYrEoICBAR48eNX1MYGCgJFmP2b59u0aNGmVz3n369FFmZqbOnTtX+okXY/v27frpp5/k6elprc/b21sXLlywXsujR48qISFBderUkZeXl7y8vHTmzJkio2lXjrZLlz6e/vMc5sDAwDKd/5XXbN++fWrQoIFcXV2tMQ888ECZzxsoTufOnfXbb79p2bJlat26tdavX69GjRrZTE8rztmzZzVw4EDr/c/Dw0M//PBDkd+R0t7b6enpuvfee+Xu7m6NiYqKKtLWtGnTFBERoapVq8rDw0Pvv/9+kXbuueceOTs7W7fN1l0e/ur9a926dWrZsqWqVasmT09PdevWTSdOnNDZs2et9To7O9u0czW9e/fWxx9/rKysLH388cfF3vPT09MVFRVl83fooYce0pkzZ/Trr7/q559/Vl5ens118/b2Vt26da3bO3bskGEYqlOnjs05pqammvrbhPJz+31zCLes9PR0hYaGSrJNov/sz+UVKlQoMlfZ7BfHSlJYWKjAwMBi542ZWebHx8dH7dq1U+/evXXhwgW1adOmyLSDr7/+Wk8++aTeeOMNtW7dWl5eXlq4cKHNnDY3N7ertuXk5GSzbbFYrroyyZ+PuXwdLx9TWFioN954o8hcakk2yaRZhYWFCg8P1/z584vsq1q1qqRLcyaPHTumCRMmKCQkRC4uLoqKiirypaOKFSuWei6Xz6cs53/lMcW95658fwF/haurq1q2bKmWLVvq9ddf1zPPPKMRI0aoR48eJR7z6quv6ssvv9T/+3//T3feeafc3Nz0j3/8o8jvyNXe21fz0UcfacCAAXrnnXcUFRUlT09Pvf3229qyZYtN3JW/izfyd+Sv3L8OHjyotm3bKiEhQW+++aa8vb21adMm9e7d2+bvhpubW7F/e0pSv3593XXXXfrnP/+pevXqqX79+kW+E1HavcVisZi6hoWFhXJwcND27dvl4OBgs88eX0D/OyNxxk1h7dq12rNnjwYMGCDp0ujyoUOHdPjwYeuo8969e5Wdna169epJupR8fffddzb17Nq1q8gfkCuXiLs8l/bKm48kNWrUSFlZWXJ0dLR+OaesevXqpbZt22rQoEHFtvHVV18pJCREQ4cOtZZd+QWUBg0aaM2aNerZs+c19eFaNGrUSPv27dOdd95ZbvWlpKRYv2RZnI0bN2rKlClq27atpEvzwo8fP14u7ZfVXXfdpfnz5ys3N1cuLi6SpG3bttmlL/h7CAsLs1l+zsnJqcgnYRs3blSPHj3UsWNHSZfmPJd1mbOwsDDNmzdP58+ft/6j/Mr74saNGxUdHa2+fftay8yMZJqp+0a42v1r27Ztys/P1zvvvKMKFS592H7l/O1r1atXL/Xt21dTp04tdn9YWJgWLVpkk0Bv3rxZnp6eqlatmqpUqSInJyd9/fXXql69uiTp5MmT+vHHH9W0aVNJ0n333aeCggIdPXpUjRs3Lpd+49owVQM3XG5urrKysnTkyBHt2LFDo0ePVmxsrNq1a6du3bpJkh599FE1aNBAXbt21Y4dO/TNN9+oW7duatq0qfVj+0ceeUTbtm3T3LlztX//fo0YMaJIIi1dSsaSkpK0b98+LViwQO+++6769+9fbN8effRRRUVFqUOHDvryyy914MABbd68WcOGDTOdRMXExOjYsWMaNWpUsfvvvPNOHTp0SAsXLtTPP/+siRMnavHixTYxI0aM0IIFCzRixAilp6drz549Gjt2rKn2r9Xrr7+uuXPnauTIkfr++++Vnp6ulJQUDRs2rNTjzp8/r127dtm8fvrpJ3Xt2lW+vr6KjY3Vxo0blZGRodTUVPXv31+//vqrpEvXYt68eUpPT9eWLVvUtWtXU6Pt18NTTz2lwsJCPfvss0pPT7eO8kkq0wgUcKUTJ07okUce0QcffKBvv/1WGRkZ+vjjjzV27FjFxsZa42rUqKE1a9YoKytLJ0+elHTpd+STTz7Rrl27tHv3buv7tCyeeuopVahQQb1799bevXu1fPly63v7sjvvvFPbtm3Tl19+qR9//FHDhw+3fkHxr9ZdkmPHjhW5d2RlZZXp3C672v2rVq1ays/P17vvvqtffvlF8+bN07Rp066prSv16dNHx44d0zPPPFPs/r59++rw4cN66aWX9MMPP2jp0qUaMWKEkpKSVKFCBXl4eKh379569dVXtWbNGn333Xfq0aOHNcGXpDp16qhr167q1q2bPvnkE2VkZGjr1q1KTk7W8uXLy+U8YA6JM264FStWKDAwUDVq1FBMTIzWrVuniRMnaunSpdYR2ssPAqhSpYqaNGmiRx99VDVr1lRKSoq1ntatW2v48OEaOHCg7r//fp0+fdqaeP9Zt27ddP78eT3wwAN64YUX9NJLL+nZZ58ttm8Wi0XLly9XkyZN1KtXL9WpU0dPPvmkDhw4IH9/f1PnZ7FY5OvrazMP8M9iY2M1YMAAvfjii2rYsKE2b96s4cOH28Q0a9ZMH3/8sZYtW6aGDRvqkUceKfKRaXlr3bq1PvvsM61atUr333+/HnzwQY0bN04hISGlHvfjjz/qvvvus3k988wzcnd314YNG1S9enV16tRJ9erVU69evXT+/HnrCPTMmTN18uRJ3XfffYqPj1e/fv1K/Rb79VSpUiV9+umn2rVrlxo2bKihQ4fq9ddfl3RtU1WAyzw8PBQZGanx48erSZMmql+/voYPH64+ffpo0qRJ1rh33nlHq1atUnBwsO677z5J0vjx41WlShVFR0erffv2at26tRo1alTm9j/99FPt3btX9913n4YOHark5GSbmISEBHXq1ElxcXGKjIzUiRMnbEaf/0rdJfnwww+L3DuuNZm92v2rYcOGGjdunJKTk1W/fn3Nnz9fY8aMuaa2ruTo6ChfX98S152uVq2ali9frm+++Ub33nuvEhIS1Lt3b5tBibfffltNmjTR448/rkcffVQPP/ywwsPDbeqZNWuWunXrppdffll169bV448/ri1btth8FwjXn8VgEh8AFGv+/Pnq2bOnsrOz7TYSDgC4eTDHGQD+f3PnzlXNmjVVrVo17d69W4MGDVKXLl1ImgEAkkicAcAqKytLr7/+urKyshQYGKgnnnjimh8AAwC4/TBVAwAAADCBLwcCAAAAJpA4AwAAACaQOAMAAAAmkDgDAAAAJpA4AwAAACaQOAPA39zs2bNVuXJle3cDAG56JM4AcIvp0aOHLBZLkVdMTMxVj61Ro4YmTJhgUxYXF6cff/zxOvX2/5CgA7jV8QAUALgFxcTEaNasWTZlLi4u11SXm5vbLfV0xIKCAlksFlWowNgPgBuLuw4A3IJcXFwUEBBg86pSpYokaeTIkapevbpcXFwUFBSkfv36SZKaNWumgwcPasCAAdZRaqnoSPDIkSPVsGFDzZw5U9WrV5eHh4eef/55FRQUaOzYsQoICJCfn1+RpyqOGzdO99xzjypWrKjg4GD17dtXZ86ckSStX79ePXv2VHZ2trXtkSNHSpJOnjypbt26qUqVKnJ3d1ebNm20f/9+a72X+/fZZ58pLCxMLi4uOnjw4PW6tABQIhJnALiN/O9//9P48eM1ffp07d+/X0uWLNE999wjSfrkk090xx13aNSoUcrMzFRmZmaJ9fz888/64osvtGLFCi1YsEAzZ87UY489pl9//VWpqalKTk7WsGHD9PXXX1uPqVChgiZOnKjvvvtOc+bM0dq1azVw4EBJUnR0tCZMmKBKlSpZ237llVckXZp6sm3bNi1btkxpaWkyDENt27bVxYsXrXWfO3dOY8aM0X//+199//338vPzux6XDwBKxVQNALgFffbZZ/Lw8LApGzRokCpWrKiAgAA9+uijcnJyUvXq1fXAAw9Ikry9veXg4CBPT08FBASUWn9hYaFmzpwpT09PhYWFqXnz5tq3b5+WL1+uChUqqG7dukpOTtb69ev14IMPSpISExOtx4eGhurNN9/U888/rylTpsjZ2VleXl6yWCw2be/fv1/Lli3TV199pejoaEnS/PnzFRwcrCVLluiJJ56QJF28eFFTpkzRvffe+5evHQBcKxJnALgFNW/eXFOnTrUp8/b21tmzZzVhwgTVrFlTMTExatu2rdq3by9Hx7Ld7mvUqCFPT0/rtr+/vxwcHGzmFfv7++vo0aPW7XXr1mn06NHau3evcnJylJ+frwsXLujs2bOqWLFise2kp6fL0dFRkZGR1jIfHx/VrVtX6enp1jJnZ2c1aNCgTOcAAOWNqRoAcAuqWLGi7rzzTpuXt7e3goODtW/fPk2ePFlubm7q27evmjRpYjPtwQwnJyebbYvFUmxZYWGhJOngwYNq27at6tevr0WLFmn79u2aPHmyJJXatmEYJZZfnoMtXfoC45+3AcAeSJwB4Dbj5uamxx9/XBMnTtT69euVlpamPXv2SLo0cltQUFDubW7btk35+fl655139OCDD6pOnTr67bffbGKKazssLEz5+fnasmWLtezEiRP68ccfVa9evXLvJwD8FUzVAIBbUG5urrKysmzKHB0d9dlnn6mgoECRkZFyd3fXvHnz5ObmppCQEEmXpmBs2LBBTz75pFxcXOTr61su/alVq5by8/P17rvvqn379vrqq680bdo0m5gaNWrozJkzWrNmje699165u7urdu3aio2NVZ8+fTR9+nR5enrqtddeU7Vq1RQbG1sufQOA8sKIMwDcglasWKHAwECb18MPP6zKlSvr/fff10MPPaQGDRpozZo1+vTTT+Xj4yNJGjVqlA4cOKBatWqpatWq5dafhg0baty4cUpOTlb9+vU1f/58jRkzxiYmOjpaCQkJiouLU9WqVTV27FhJ0qxZsxQeHq527dopKipKhmFo+fLlRaaGAIC9WYySJpgBAAAAsGLEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAwgcQZAAAAMIHEGQAAADCBxBkAAAAw4f8Da/9w3fg0xPgAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 800x500 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import seaborn as sns\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "mse_df[\"MSE\"] = mse_df[\"mse\"]\n",
    "\n",
    "mse_df[\"Sample Size\"] = mse_df[\"n\"]\n",
    "sns.barplot(data=mse_df.loc[ (mse_df[\"estimator_name\"] != \"dr_linear_statsmodels\"), :], \n",
    "            x=\"Estimator\", y=\"MSE\", hue=\"Sample Size\",\n",
    "           palette = sns.color_palette(\"Blues_r\"))\n",
    "plt.title(\"Mean Squared Error of Estimators\")\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 124,
   "id": "f874e15b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th>n</th>\n",
       "      <th>Name</th>\n",
       "      <th>400</th>\n",
       "      <th>1600</th>\n",
       "      <th>6400</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>Double Machine Learning</td>\n",
       "      <td>0.002816</td>\n",
       "      <td>0.000671</td>\n",
       "      <td>0.000133</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>Doubly Robust Estimation</td>\n",
       "      <td>1778149.937155</td>\n",
       "      <td>376383.085348</td>\n",
       "      <td>0.016019</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>Standard Linear Model</td>\n",
       "      <td>0.003114</td>\n",
       "      <td>0.000235</td>\n",
       "      <td>4.2e-05</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "n                      Name             400           1600      6400\n",
       "0   Double Machine Learning        0.002816       0.000671  0.000133\n",
       "1  Doubly Robust Estimation  1778149.937155  376383.085348  0.016019\n",
       "2     Standard Linear Model        0.003114       0.000235   4.2e-05"
      ]
     },
     "execution_count": 124,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mse_df[\"mse\"]  = mse_df[\"mse\"].astype(float).round(6).astype(str)\n",
    "mse_df.pivot(index = \"Name\", columns = \"n\", values = \"mse\").reset_index()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9b06000f",
   "metadata": {},
   "source": [
    "Finally, let's look at coverage probabilities. We calculated 95% confidence intervals for our estimates using both bootstrap estimation and estimation based on a known distribution. In bootstrap variance estimation, the dataset is duplicated many times by sampling with replacement and the quantiles from the distribution of these estimations are used to determine confidence intervals. In estimation based on a known distribution, statistical theory is used to determine what the variance of the estimate should be, and confidence intervals are generated from the normal or t distribution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 143,
   "id": "8c41c8e5",
   "metadata": {},
   "outputs": [],
   "source": [
    "translation_dict = {\"dml_linear_bootstrap\": \"Double ML\\n (Bootstrap)\",\n",
    "                    \"dml_linear_statsmodels\": \"Double ML\\n (Non-Bootstrap)\",\n",
    "                    \"dr_linear_bootstrap\": \"Doubly Robust\\n(Bootstrap)\",\n",
    "                    \"dr_linear_statsmodels\": \"Doubly Robust\\n(Non-Bootstrap)\",\n",
    "                    \"standard_linear\": \"Standard Linear\\nModel\"\n",
    "                   }"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 144,
   "id": "267b5177",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "## coverage\n",
    "\n",
    "all_estimates_df[\"covered\"] = (all_estimates_df[\"upper\"] > 2) & (all_estimates_df[\"lower\"] < 2)\n",
    "coverage_df = all_estimates_df.groupby([\"estimator_name\", \"n\"])[\"covered\"].agg(\"mean\").reset_index()\n",
    "\n",
    "coverage_df = coverage_df.loc[coverage_df[\"estimator_name\"] != \"standard_linear_small\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 145,
   "id": "656b9b78",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>estimator_name</th>\n",
       "      <th>n</th>\n",
       "      <th>covered</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>dml_linear_bootstrap</td>\n",
       "      <td>400</td>\n",
       "      <td>0.98</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>dml_linear_bootstrap</td>\n",
       "      <td>1600</td>\n",
       "      <td>0.96</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>dml_linear_bootstrap</td>\n",
       "      <td>6400</td>\n",
       "      <td>1.00</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>dml_linear_statsmodels</td>\n",
       "      <td>400</td>\n",
       "      <td>0.90</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>dml_linear_statsmodels</td>\n",
       "      <td>1600</td>\n",
       "      <td>0.88</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>dml_linear_statsmodels</td>\n",
       "      <td>6400</td>\n",
       "      <td>0.96</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>dr_linear_bootstrap</td>\n",
       "      <td>400</td>\n",
       "      <td>0.98</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>dr_linear_bootstrap</td>\n",
       "      <td>1600</td>\n",
       "      <td>0.96</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>dr_linear_bootstrap</td>\n",
       "      <td>6400</td>\n",
       "      <td>1.00</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>dr_linear_statsmodels</td>\n",
       "      <td>400</td>\n",
       "      <td>0.64</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>dr_linear_statsmodels</td>\n",
       "      <td>1600</td>\n",
       "      <td>0.66</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>dr_linear_statsmodels</td>\n",
       "      <td>6400</td>\n",
       "      <td>0.62</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>standard_linear</td>\n",
       "      <td>400</td>\n",
       "      <td>0.92</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13</th>\n",
       "      <td>standard_linear</td>\n",
       "      <td>1600</td>\n",
       "      <td>0.96</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14</th>\n",
       "      <td>standard_linear</td>\n",
       "      <td>6400</td>\n",
       "      <td>0.96</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "            estimator_name     n  covered\n",
       "0     dml_linear_bootstrap   400     0.98\n",
       "1     dml_linear_bootstrap  1600     0.96\n",
       "2     dml_linear_bootstrap  6400     1.00\n",
       "3   dml_linear_statsmodels   400     0.90\n",
       "4   dml_linear_statsmodels  1600     0.88\n",
       "5   dml_linear_statsmodels  6400     0.96\n",
       "6      dr_linear_bootstrap   400     0.98\n",
       "7      dr_linear_bootstrap  1600     0.96\n",
       "8      dr_linear_bootstrap  6400     1.00\n",
       "9    dr_linear_statsmodels   400     0.64\n",
       "10   dr_linear_statsmodels  1600     0.66\n",
       "11   dr_linear_statsmodels  6400     0.62\n",
       "12         standard_linear   400     0.92\n",
       "13         standard_linear  1600     0.96\n",
       "14         standard_linear  6400     0.96"
      ]
     },
     "execution_count": 145,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "coverage_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 161,
   "id": "99ca067f",
   "metadata": {},
   "outputs": [],
   "source": [
    "coverage_df[\"Estimator\"] = coverage_df[\"estimator_name\"].apply(lambda x: translation_dict[x])\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4c44fc62",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 147,
   "id": "879e0046",
   "metadata": {},
   "outputs": [],
   "source": [
    "coverage_df[\"Sample Size\"] = coverage_df[\"n\"]\n",
    "\n",
    "coverage_df[\"Coverage Probability\"] = coverage_df[\"covered\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 154,
   "id": "b04ca322",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['ArtistList',\n",
       " '_AxesBase__clear',\n",
       " '_PROPERTIES_EXCLUDED_FROM_SET',\n",
       " '__class__',\n",
       " '__delattr__',\n",
       " '__dict__',\n",
       " '__dir__',\n",
       " '__doc__',\n",
       " '__eq__',\n",
       " '__format__',\n",
       " '__ge__',\n",
       " '__getattribute__',\n",
       " '__getstate__',\n",
       " '__gt__',\n",
       " '__hash__',\n",
       " '__init__',\n",
       " '__init_subclass__',\n",
       " '__le__',\n",
       " '__lt__',\n",
       " '__module__',\n",
       " '__ne__',\n",
       " '__new__',\n",
       " '__reduce__',\n",
       " '__reduce_ex__',\n",
       " '__repr__',\n",
       " '__setattr__',\n",
       " '__setstate__',\n",
       " '__sizeof__',\n",
       " '__str__',\n",
       " '__subclasshook__',\n",
       " '__weakref__',\n",
       " '_add_text',\n",
       " '_adjustable',\n",
       " '_agg_filter',\n",
       " '_alias_map',\n",
       " '_alpha',\n",
       " '_anchor',\n",
       " '_animated',\n",
       " '_aspect',\n",
       " '_autotitlepos',\n",
       " '_axes',\n",
       " '_axes_class',\n",
       " '_axes_locator',\n",
       " '_axis_map',\n",
       " '_axis_names',\n",
       " '_axisbelow',\n",
       " '_box_aspect',\n",
       " '_callbacks',\n",
       " '_check_no_units',\n",
       " '_children',\n",
       " '_clipon',\n",
       " '_clippath',\n",
       " '_cm_set',\n",
       " '_colorbars',\n",
       " '_convert_dx',\n",
       " '_current_image',\n",
       " '_default_contains',\n",
       " '_deprecate_noninstance',\n",
       " '_errorevery_to_mask',\n",
       " '_facecolor',\n",
       " '_fill_between_x_or_y',\n",
       " '_frameon',\n",
       " '_fully_clipped_to_axes',\n",
       " '_gci',\n",
       " '_gen_axes_patch',\n",
       " '_gen_axes_spines',\n",
       " '_get_aspect_ratio',\n",
       " '_get_lines',\n",
       " '_get_pan_points',\n",
       " '_get_patches_for_fill',\n",
       " '_get_view',\n",
       " '_gid',\n",
       " '_gridOn',\n",
       " '_in_layout',\n",
       " '_init_axis',\n",
       " '_internal_update',\n",
       " '_label',\n",
       " '_label_outer_xaxis',\n",
       " '_label_outer_yaxis',\n",
       " '_left_title',\n",
       " '_make_twin_axes',\n",
       " '_mouseover',\n",
       " '_mouseover_set',\n",
       " '_navigate',\n",
       " '_navigate_mode',\n",
       " '_originalPosition',\n",
       " '_parse_scatter_color_args',\n",
       " '_path_effects',\n",
       " '_pcolor_grid_deprecation_helper',\n",
       " '_pcolorargs',\n",
       " '_picker',\n",
       " '_position',\n",
       " '_prepare_view_from_bbox',\n",
       " '_process_unit_info',\n",
       " '_projection_init',\n",
       " '_quiver_units',\n",
       " '_rasterization_zorder',\n",
       " '_rasterized',\n",
       " '_remove_legend',\n",
       " '_remove_method',\n",
       " '_request_autoscale_view',\n",
       " '_right_title',\n",
       " '_sci',\n",
       " '_set_alpha_for_array',\n",
       " '_set_artist_props',\n",
       " '_set_gc_clip',\n",
       " '_set_lim_and_transforms',\n",
       " '_set_position',\n",
       " '_set_title_offset_trans',\n",
       " '_set_view',\n",
       " '_set_view_from_bbox',\n",
       " '_shared_axes',\n",
       " '_sharex',\n",
       " '_sharey',\n",
       " '_sketch',\n",
       " '_snap',\n",
       " '_stale',\n",
       " '_stale_viewlims',\n",
       " '_sticky_edges',\n",
       " '_subclass_uses_cla',\n",
       " '_subplotspec',\n",
       " '_tight',\n",
       " '_transform',\n",
       " '_transformSet',\n",
       " '_twinned_axes',\n",
       " '_unit_change_handler',\n",
       " '_unstale_viewLim',\n",
       " '_update_image_limits',\n",
       " '_update_line_limits',\n",
       " '_update_patch_limits',\n",
       " '_update_props',\n",
       " '_update_set_signature_and_docstring',\n",
       " '_update_title_position',\n",
       " '_update_transScale',\n",
       " '_url',\n",
       " '_use_sticky_edges',\n",
       " '_validate_converted_limits',\n",
       " '_viewLim',\n",
       " '_visible',\n",
       " '_xaxis_transform',\n",
       " '_xmargin',\n",
       " '_yaxis_transform',\n",
       " '_ymargin',\n",
       " 'acorr',\n",
       " 'add_artist',\n",
       " 'add_callback',\n",
       " 'add_child_axes',\n",
       " 'add_collection',\n",
       " 'add_container',\n",
       " 'add_image',\n",
       " 'add_line',\n",
       " 'add_patch',\n",
       " 'add_table',\n",
       " 'angle_spectrum',\n",
       " 'annotate',\n",
       " 'apply_aspect',\n",
       " 'arrow',\n",
       " 'artists',\n",
       " 'autoscale',\n",
       " 'autoscale_view',\n",
       " 'axes',\n",
       " 'axhline',\n",
       " 'axhspan',\n",
       " 'axis',\n",
       " 'axison',\n",
       " 'axline',\n",
       " 'axvline',\n",
       " 'axvspan',\n",
       " 'bar',\n",
       " 'bar_label',\n",
       " 'barbs',\n",
       " 'barh',\n",
       " 'bbox',\n",
       " 'boxplot',\n",
       " 'broken_barh',\n",
       " 'bxp',\n",
       " 'callbacks',\n",
       " 'can_pan',\n",
       " 'can_zoom',\n",
       " 'child_axes',\n",
       " 'cla',\n",
       " 'clabel',\n",
       " 'clear',\n",
       " 'clipbox',\n",
       " 'cohere',\n",
       " 'collections',\n",
       " 'containers',\n",
       " 'contains',\n",
       " 'contains_point',\n",
       " 'contour',\n",
       " 'contourf',\n",
       " 'convert_xunits',\n",
       " 'convert_yunits',\n",
       " 'csd',\n",
       " 'dataLim',\n",
       " 'drag_pan',\n",
       " 'draw',\n",
       " 'draw_artist',\n",
       " 'end_pan',\n",
       " 'errorbar',\n",
       " 'eventplot',\n",
       " 'figure',\n",
       " 'fill',\n",
       " 'fill_between',\n",
       " 'fill_betweenx',\n",
       " 'findobj',\n",
       " 'fmt_xdata',\n",
       " 'fmt_ydata',\n",
       " 'format_coord',\n",
       " 'format_cursor_data',\n",
       " 'format_xdata',\n",
       " 'format_ydata',\n",
       " 'get_adjustable',\n",
       " 'get_agg_filter',\n",
       " 'get_alpha',\n",
       " 'get_anchor',\n",
       " 'get_animated',\n",
       " 'get_aspect',\n",
       " 'get_autoscale_on',\n",
       " 'get_autoscalex_on',\n",
       " 'get_autoscaley_on',\n",
       " 'get_axes_locator',\n",
       " 'get_axisbelow',\n",
       " 'get_box_aspect',\n",
       " 'get_children',\n",
       " 'get_clip_box',\n",
       " 'get_clip_on',\n",
       " 'get_clip_path',\n",
       " 'get_cursor_data',\n",
       " 'get_data_ratio',\n",
       " 'get_default_bbox_extra_artists',\n",
       " 'get_facecolor',\n",
       " 'get_fc',\n",
       " 'get_figure',\n",
       " 'get_frame_on',\n",
       " 'get_gid',\n",
       " 'get_gridspec',\n",
       " 'get_images',\n",
       " 'get_in_layout',\n",
       " 'get_label',\n",
       " 'get_legend',\n",
       " 'get_legend_handles_labels',\n",
       " 'get_lines',\n",
       " 'get_mouseover',\n",
       " 'get_navigate',\n",
       " 'get_navigate_mode',\n",
       " 'get_path_effects',\n",
       " 'get_picker',\n",
       " 'get_position',\n",
       " 'get_rasterization_zorder',\n",
       " 'get_rasterized',\n",
       " 'get_renderer_cache',\n",
       " 'get_shared_x_axes',\n",
       " 'get_shared_y_axes',\n",
       " 'get_sketch_params',\n",
       " 'get_snap',\n",
       " 'get_subplotspec',\n",
       " 'get_tightbbox',\n",
       " 'get_title',\n",
       " 'get_transform',\n",
       " 'get_transformed_clip_path_and_affine',\n",
       " 'get_url',\n",
       " 'get_visible',\n",
       " 'get_window_extent',\n",
       " 'get_xaxis',\n",
       " 'get_xaxis_text1_transform',\n",
       " 'get_xaxis_text2_transform',\n",
       " 'get_xaxis_transform',\n",
       " 'get_xbound',\n",
       " 'get_xgridlines',\n",
       " 'get_xlabel',\n",
       " 'get_xlim',\n",
       " 'get_xmajorticklabels',\n",
       " 'get_xminorticklabels',\n",
       " 'get_xscale',\n",
       " 'get_xticklabels',\n",
       " 'get_xticklines',\n",
       " 'get_xticks',\n",
       " 'get_yaxis',\n",
       " 'get_yaxis_text1_transform',\n",
       " 'get_yaxis_text2_transform',\n",
       " 'get_yaxis_transform',\n",
       " 'get_ybound',\n",
       " 'get_ygridlines',\n",
       " 'get_ylabel',\n",
       " 'get_ylim',\n",
       " 'get_ymajorticklabels',\n",
       " 'get_yminorticklabels',\n",
       " 'get_yscale',\n",
       " 'get_yticklabels',\n",
       " 'get_yticklines',\n",
       " 'get_yticks',\n",
       " 'get_zorder',\n",
       " 'grid',\n",
       " 'has_data',\n",
       " 'have_units',\n",
       " 'hexbin',\n",
       " 'hist',\n",
       " 'hist2d',\n",
       " 'hlines',\n",
       " 'ignore_existing_data_limits',\n",
       " 'images',\n",
       " 'imshow',\n",
       " 'in_axes',\n",
       " 'indicate_inset',\n",
       " 'indicate_inset_zoom',\n",
       " 'inset_axes',\n",
       " 'invert_xaxis',\n",
       " 'invert_yaxis',\n",
       " 'is_transform_set',\n",
       " 'label_outer',\n",
       " 'legend',\n",
       " 'legend_',\n",
       " 'lines',\n",
       " 'locator_params',\n",
       " 'loglog',\n",
       " 'magnitude_spectrum',\n",
       " 'margins',\n",
       " 'matshow',\n",
       " 'minorticks_off',\n",
       " 'minorticks_on',\n",
       " 'mouseover',\n",
       " 'name',\n",
       " 'patch',\n",
       " 'patches',\n",
       " 'pchanged',\n",
       " 'pcolor',\n",
       " 'pcolorfast',\n",
       " 'pcolormesh',\n",
       " 'phase_spectrum',\n",
       " 'pick',\n",
       " 'pickable',\n",
       " 'pie',\n",
       " 'plot',\n",
       " 'plot_date',\n",
       " 'properties',\n",
       " 'psd',\n",
       " 'quiver',\n",
       " 'quiverkey',\n",
       " 'redraw_in_frame',\n",
       " 'relim',\n",
       " 'remove',\n",
       " 'remove_callback',\n",
       " 'reset_position',\n",
       " 'scatter',\n",
       " 'secondary_xaxis',\n",
       " 'secondary_yaxis',\n",
       " 'semilogx',\n",
       " 'semilogy',\n",
       " 'set',\n",
       " 'set_adjustable',\n",
       " 'set_agg_filter',\n",
       " 'set_alpha',\n",
       " 'set_anchor',\n",
       " 'set_animated',\n",
       " 'set_aspect',\n",
       " 'set_autoscale_on',\n",
       " 'set_autoscalex_on',\n",
       " 'set_autoscaley_on',\n",
       " 'set_axes_locator',\n",
       " 'set_axis_off',\n",
       " 'set_axis_on',\n",
       " 'set_axisbelow',\n",
       " 'set_box_aspect',\n",
       " 'set_clip_box',\n",
       " 'set_clip_on',\n",
       " 'set_clip_path',\n",
       " 'set_facecolor',\n",
       " 'set_fc',\n",
       " 'set_figure',\n",
       " 'set_frame_on',\n",
       " 'set_gid',\n",
       " 'set_in_layout',\n",
       " 'set_label',\n",
       " 'set_mouseover',\n",
       " 'set_navigate',\n",
       " 'set_navigate_mode',\n",
       " 'set_path_effects',\n",
       " 'set_picker',\n",
       " 'set_position',\n",
       " 'set_prop_cycle',\n",
       " 'set_rasterization_zorder',\n",
       " 'set_rasterized',\n",
       " 'set_sketch_params',\n",
       " 'set_snap',\n",
       " 'set_subplotspec',\n",
       " 'set_title',\n",
       " 'set_transform',\n",
       " 'set_url',\n",
       " 'set_visible',\n",
       " 'set_xbound',\n",
       " 'set_xlabel',\n",
       " 'set_xlim',\n",
       " 'set_xmargin',\n",
       " 'set_xscale',\n",
       " 'set_xticklabels',\n",
       " 'set_xticks',\n",
       " 'set_ybound',\n",
       " 'set_ylabel',\n",
       " 'set_ylim',\n",
       " 'set_ymargin',\n",
       " 'set_yscale',\n",
       " 'set_yticklabels',\n",
       " 'set_yticks',\n",
       " 'set_zorder',\n",
       " 'sharex',\n",
       " 'sharey',\n",
       " 'specgram',\n",
       " 'spines',\n",
       " 'spy',\n",
       " 'stackplot',\n",
       " 'stairs',\n",
       " 'stale',\n",
       " 'stale_callback',\n",
       " 'start_pan',\n",
       " 'stem',\n",
       " 'step',\n",
       " 'sticky_edges',\n",
       " 'streamplot',\n",
       " 'table',\n",
       " 'tables',\n",
       " 'text',\n",
       " 'texts',\n",
       " 'tick_params',\n",
       " 'ticklabel_format',\n",
       " 'title',\n",
       " 'titleOffsetTrans',\n",
       " 'transAxes',\n",
       " 'transData',\n",
       " 'transLimits',\n",
       " 'transScale',\n",
       " 'tricontour',\n",
       " 'tricontourf',\n",
       " 'tripcolor',\n",
       " 'triplot',\n",
       " 'twinx',\n",
       " 'twiny',\n",
       " 'update',\n",
       " 'update_datalim',\n",
       " 'update_from',\n",
       " 'use_sticky_edges',\n",
       " 'viewLim',\n",
       " 'violin',\n",
       " 'violinplot',\n",
       " 'vlines',\n",
       " 'xaxis',\n",
       " 'xaxis_date',\n",
       " 'xaxis_inverted',\n",
       " 'xcorr',\n",
       " 'yaxis',\n",
       " 'yaxis_date',\n",
       " 'yaxis_inverted',\n",
       " 'zorder']"
      ]
     },
     "execution_count": 154,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "dir(ax)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 156,
   "id": "bf022ebb",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[Text(0, 0, 'Double ML\\n (Bootstrap)'),\n",
       " Text(1, 0, 'Double ML\\n (Non-Bootstrap)'),\n",
       " Text(2, 0, 'Doubly Robust\\n(Bootstrap)'),\n",
       " Text(3, 0, 'Doubly Robust\\n(Non-Bootstrap)'),\n",
       " Text(4, 0, 'Standard Linear\\nModel')]"
      ]
     },
     "execution_count": 156,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ax.get_xticklabels()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 164,
   "id": "94ecee10",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAArMAAAHkCAYAAAAthKfJAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAABz8klEQVR4nO3dd1QU19sH8O/Sll4EaYqAIAixa1Q0RhEVazRNk9iwJYol0RQ1mlhiiRoRTSyJjZiowcQau7HFWBJF7FhiwwKxRREL9Xn/8GV+jiy4a8Bl9fs5Z89x79yZee7cu+uzw50ZjYgIiIiIiIhMkJmxAyAiIiIielJMZomIiIjIZDGZJSIiIiKTxWSWiIiIiEwWk1kiIiIiMllMZomIiIjIZDGZJSIiIiKTxWSWiIiIiEwWk1kiIiIiMllMZum5c+jQIXTr1g3+/v6wtraGvb09atSogYkTJ+LGjRvGDs+kxMXFQaPRKC8LCwuULVsW3bp1w6VLl4p0XxqNBv369Suy7Z07dw4ajQZfffXVY+vmtfPcuXNKWVRUFPz8/FT1/Pz8EBUVpby/fPkyRo4ciQMHDuTb5siRI6HRaJ4w+uIXHx+PF154ATY2NtBoNDrbMHDgQGg0Ghw/frzA7QwbNgwajQb79+8vkrg0Gg1GjhxZJNsqalFRUarPw6MvQ+zatQsjR47EzZs38y1r1KgRGjVqVDRBG2jt2rUl9vjT88vC2AEQPU2zZ89GdHQ0goOD8fHHHyM0NBRZWVnYt28fZs2ahd27d2P58uXGDtPkzJ8/HxUrVsS9e/fw+++/Y/z48di+fTsOHz4MOzs7Y4f3n7Vq1Qq7d++Gl5dXofWWL18OR0dH5f3ly5cxatQo+Pn5oVq1aqq6PXv2RPPmzYsj3P/s6tWr6Ny5M5o3b44ZM2ZAq9UiKCgoX70ePXogNjYW8+bNw8SJE/Mtz83NxYIFC1CtWjXUqFGjSGLbvXs3ypYtWyTbKg42NjbYsmXLf97Orl27MGrUKERFRcHZ2Vm1bMaMGf95+09q7dq1mD59OhNaKlGYzNJzY/fu3ejTpw+aNm2KFStWQKvVKsuaNm2KDz/8EOvXrzdihAW7e/cubG1tjR1GgSpVqoRatWoBAMLDw5GTk4MvvvgCK1asQMeOHXWuU9Lb9LDSpUujdOnSj61XvXp1vbdZtmzZEpuUnTx5EllZWejUqRMaNmxYYL1KlSqhdu3a+OGHHzBu3DhYWKj/S9m4cSMuXryIwYMH/6d4RAT379+HjY0N6tat+5+2VdzMzMyKPcbQ0NBi3b4x3Lt3DzY2NsYOg0wUpxnQc2PcuHHQaDT47rvvVIlsHisrK7zyyivK+9zcXEycOBEVK1aEVquFu7s7unTpgosXLyp1PvjgA9jZ2SEtLS3f9jp06AAPDw9kZWUpZfHx8QgLC4OdnR3s7e0RGRmJxMRE1XpRUVGwt7fH4cOH0axZMzg4OCAiIgIAsGnTJrRt2xZly5aFtbU1AgMD8d577+HatWv59r9y5UpUqVIFWq0W5cuXx9SpU3X+aVtEMGPGDFSrVg02NjZwcXHBG2+8gTNnzuh5ZPPL+8/8/Pnzj23TjRs3EB0djTJlysDKygrly5fHsGHDkJGRoXPb3377LYKCgqDVahEaGoqffvpJtfzq1auIjo5GaGgo7O3t4e7ujsaNG2PHjh06t5ebm4uxY8eiXLlysLa2Rq1atbB582ZVHV3TDHR5eJrBtm3b8OKLLwIAunXrpvypOe+MVkHTDPQZI2fOnMFbb70Fb29vaLVaeHh4ICIiQudUgEetWrUKYWFhsLW1hYODA5o2bYrdu3cry6OiovDSSy8BeDCGNRpNoX/S7tGjB1JTU7Fu3bp8y+bPnw+tVouOHTvi/v37+PDDD1GtWjU4OTmhVKlSCAsLw8qVK/OtlzelZNasWQgJCYFWq8X333+vLHv0rOCRI0fQtm1buLi4wNraGtWqVVPq5ymoD7dt2waNRoNt27YpZYmJiWjdujXc3d2h1Wrh7e2NVq1aqT77/0Vubi7GjBmD4OBg2NjYwNnZGVWqVMHUqVMBPBgbH3/8MQDA399fGTt5MT46zSBvysykSZMwYcIE+Pn5wcbGBo0aNVJ+mAwZMgTe3t5wcnLCq6++iitXrqhiio+PR7NmzeDl5QUbGxuEhIRgyJAhuHPnjlInKioK06dPBwDV9Im8Y3r//n0MHToU/v7+sLKyQpkyZdC3b998UyX8/PzQunVrLFu2DNWrV4e1tTVGjRoFAPj5559Rp04dODk5wdbWFuXLl0f37t2L5LjTM0yIngPZ2dlia2srderU0Xudd999VwBIv379ZP369TJr1iwpXbq0+Pj4yNWrV0VE5ODBgwJAZs+erVr333//Fa1WK4MGDVLKxo4dKxqNRrp37y6rV6+WZcuWSVhYmNjZ2cnRo0eVel27dhVLS0vx8/OT8ePHy+bNm2XDhg0iIjJz5kwZP368rFq1SrZv3y7ff/+9VK1aVYKDgyUzM1PZxrp168TMzEwaNWoky5cvl59//lnq1Kkjfn5+8ujHvlevXmJpaSkffvihrF+/XhYtWiQVK1YUDw8PSU1NLfQYzZ8/XwDI3r17VeVTp04VAPLdd98V2qZ79+5JlSpVxM7OTr766ivZuHGjfPbZZ2JhYSEtW7ZUbROA+Pj4SGhoqCxevFhWrVolzZs3FwDy888/K/WOHz8uffr0kZ9++km2bdsmq1evlh49eoiZmZls3bpVqXf27Fllmy+99JIsXbpUfv75Z3nxxRfF0tJSdu3ala+dZ8+eVfWTr6+vKkZfX1/p2rWriIjcunVLWW/48OGye/du2b17t1y4cEFEREaMGJGvL/QdI8HBwRIYGCg//PCDbN++XZYuXSoffvihqn26LFy4UABIs2bNZMWKFRIfHy81a9YUKysr2bFjh4iI/P333zJ9+nQBIOPGjZPdu3er9v2otLQ0sbW1lXbt2qnKb9y4IVqtVt566y0REbl586ZERUXJDz/8IFu2bJH169fLRx99JGZmZvL999+r1gUgZcqUkSpVqsiiRYtky5YtcuTIEWXZiBEjlLrHjx8XBwcHCQgIkAULFsiaNWvk7bffFgAyYcIEpZ6uPhQR2bp1qwBQjl16erq4urpKrVq1ZMmSJbJ9+3aJj4+X3r17y7Fjxwo9vl27dhU7OzvJysrK98rJyVHqjR8/XszNzWXEiBGyefNmWb9+vcTGxsrIkSNFROTChQvSv39/ASDLli1Txs6tW7dERKRhw4bSsGFDZXt5Y9nX11fatGkjq1evlh9//FE8PDwkKChIOnfuLN27d5d169bJrFmzxN7eXtq0aaOK/YsvvpApU6bImjVrZNu2bTJr1izx9/eX8PBwpc7ff/8tb7zxhgBQYtq9e7fcv39fcnNzJTIyUiwsLOSzzz6TjRs3yldffSV2dnZSvXp1uX//vrIdX19f8fLykvLly8u8efNk69at8tdff8muXbtEo9HIW2+9JWvXrpUtW7bI/PnzpXPnzoUedyIms/RcSE1NFQDKf6yPk5SUJAAkOjpaVf7nn38KAPn000+Vsho1aki9evVU9WbMmCEA5PDhwyIikpycLBYWFtK/f39Vvdu3b4unp6e0b99eKevatasAkHnz5hUaY25urmRlZcn58+cFgKxcuVJZ9uKLL4qPj49kZGSo9uXq6qpKoHbv3i0AZPLkyaptX7hwQWxsbOSTTz4pNIa8BGHPnj2SlZUlt2/fltWrV0vp0qXFwcFBSYYLatOsWbMEgCxZskRVPmHCBAEgGzduVMoAiI2NjSrBzs7OlooVK0pgYGCBMWZnZ0tWVpZERETIq6++qpTnJQDe3t5y7949pTwtLU1KlSolTZo0yddOQ5JZEZG9e/cKAJk/f36+uB5NZvUdI9euXRMAEhsbW2CbdcnJyRFvb2+pXLmyKrG6ffu2uLu7q8ZwXoL38I+EwuT9WPnnn3+Usq+//loAyKZNm3Suk9cvPXr0kOrVq6uWARAnJye5ceNGvvUeTWbfeust0Wq1kpycrKrXokULsbW1lZs3b4qI/snsvn37BICsWLFCr7Y/LG+c63pFREQo9Vq3bi3VqlUrdFuTJk3SGa9Iwcls1apVVX0bGxsrAOSVV15Rrf/BBx8IACU5flTed8v27dsFgBw8eFBZ1rdv33w/wkRE1q9fLwBk4sSJqvL4+HjVD1uRB58Tc3NzOXHihKruV199JQCUPiPSF6cZEOmwdetWAFBdmQ4AtWvXRkhIiOrP0N26dcOuXbtw4sQJpWz+/Pl48cUXUalSJQDAhg0bkJ2djS5duiA7O1t5WVtbo2HDhqo/ceZ5/fXX85VduXIFvXv3ho+PDywsLGBpaQlfX18AQFJSEgDgzp072LdvH9q1awcrKytlXXt7e7Rp00a1vdWrV0Oj0aBTp06quDw9PVG1alWdcelSt25dWFpawsHBAa1bt4anpyfWrVsHDw+PQtu0ZcsW2NnZ4Y033lCV5x33R//cHxERodqmubk5OnTogL///lv1J+BZs2ahRo0asLa2Vo7T5s2blWP0sNdeew3W1tbKewcHB7Rp0wa///47cnJy9Gp/UdB3jJQqVQoBAQGYNGkSYmJikJiYiNzc3Mdu/8SJE7h8+TI6d+4MM7P/ffXb29vj9ddfx549e3D37t0nir1Hjx7IysrCDz/8oJTNnz8fvr6+ynQS4MGfkOvXrw97e3ulX+bOnauzXxo3bgwXF5fH7nvLli2IiIiAj4+PqjwqKgp3795VTaHQR2BgIFxcXDB48GDMmjULx44dM2h9Gxsb7N27N9/r4Yu2ateujYMHDyI6OhobNmzQOU3pSbRs2VLVtyEhIQAeXMD4sLzy5ORkpezMmTN455134OnpCXNzc1haWirzpXX1z6PyLnp79DvzzTffhJ2dXb7PcpUqVfJdVJg3Lad9+/ZYsmRJkd8RhZ5dTGbpueDm5gZbW1ucPXtWr/rXr18HAJ1Xr3t7eyvLAaBjx47QarWIi4sDABw7dgx79+5Ft27dlDr//PMPgAdf1paWlqpXfHx8vjmvtra2qqvigQfz7Jo1a4Zly5bhk08+webNm/HXX39hz549AB5cQAEA//77L0QkXyIJIF/ZP//8o9R9NK49e/bonIury4IFC7B3714kJibi8uXLOHToEOrXr//YNl2/fh2enp755o66u7vDwsJCdZwBwNPTM9++88ry6sbExKBPnz6oU6cOli5dij179mDv3r1o3ry5coz02WZmZibS09P1aH3R0HeMaDQabN68GZGRkZg4cSJq1KiB0qVLY8CAAbh9+3aB23/cmM7NzcW///77RLE3aNAAQUFBmD9/PoAHt7/bv3+/MlcYAJYtW4b27dujTJky+PHHH7F7927s3bsX3bt3x/379/Nt83F3jni4XQW1KW+5IZycnLB9+3ZUq1YNn376KV544QV4e3tjxIgRqvnvBTEzM0OtWrXyvR5O3IYOHYqvvvoKe/bsQYsWLeDq6oqIiAjs27fPoFgfVapUKdX7vB+zBZXnHff09HQ0aNAAf/75J8aMGYNt27Zh7969WLZsGQDo/Nw86vr167CwsMh3oaRGo4Gnp2e+ftDVZy+//DJWrFih/KgrW7YsKlWqhMWLFz92//R8490M6Llgbm6OiIgIrFu3DhcvXnzsVeSurq4AgJSUlHx1L1++DDc3N+W9i4sL2rZtiwULFmDMmDGYP38+rK2t8fbbbyt18ur/8ssvypnUwui6MOjIkSM4ePAg4uLi0LVrV6X877//VtVzcXGBRqNRkqOHpaamqt67ublBo9Fgx44dOi+K01WmS0hIiHI3g4LoapOrqyv+/PNPiIhq+ZUrV5Cdna06zrrif7gsr89+/PFHNGrUCDNnzlTVKyjRK2ibVlZWsLe3L7RNRcmQMeLr64u5c+cCeHDngSVLlmDkyJHIzMzErFmzdK7z8Jh+1OXLl2FmZqbXmdCCdO/eHUOGDMFff/2FRYsWwczMTHWW7scff4S/vz/i4+NVfV3QhX763pfV1dW1wDYB/zuueWffH92frh9slStXxk8//QQRwaFDhxAXF4fRo0fDxsYGQ4YM0SuuwlhYWGDQoEEYNGgQbt68id9++w2ffvopIiMjceHChad+l48tW7bg8uXL2LZtm+ruFbrucVsQV1dXZGdn4+rVq6qEVkSQmpqqnHXNU1D/tm3bFm3btkVGRgb27NmD8ePH45133oGfnx/CwsIMaxg9N3hmlp4bQ4cOhYigV69eyMzMzLc8KysLv/76K4AHf+IEHvwH/LC9e/ciKSlJ9adT4MFUg8uXL2Pt2rX48ccf8eqrr6ruDRkZGQkLCwucPn1a51mbxyWCwP++/B9NML/99lvVezs7O9SqVQsrVqxQtTM9PR2rV69W1W3dujVEBJcuXdIZU+XKlR8b138RERGB9PR0rFixQlW+YMECZfnDNm/erErSc3JyEB8fj4CAAOVHh0ajyXeMDh06VOCfm5ctW6Y6M3j79m38+uuvaNCgAczNzZ+4bcD/+kqfM1tPOkaCgoIwfPhwVK5cudAHEwQHB6NMmTJYtGgRREQpv3PnDpYuXarc4eBJde3aFRYWFvj222+xcOFCREREqJJyjUYDKysrVRKTmpqq824GhoiIiFCSsYctWLAAtra2yp018h5wcejQIVW9VatWFbhtjUaDqlWrYsqUKXB2di6yBz88zNnZGW+88Qb69u2LGzduKHcGMGTs/Ff6frcUFlfeZ/XR78ylS5fizp07+T7Lj6PVatGwYUNMmDABAPLd0YPoYTwzS8+NsLAwzJw5E9HR0ahZsyb69OmDF154AVlZWUhMTMR3332HSpUqoU2bNggODsa7776Lr7/+GmZmZmjRogXOnTuHzz77DD4+Phg4cKBq282aNUPZsmURHR2N1NRU1RQD4MF/pKNHj8awYcNw5swZNG/eHC4uLvjnn3/w119/wc7OTrk1TUEqVqyIgIAADBkyBCKCUqVK4ddff8WmTZvy1R09ejRatWqFyMhIvP/++8jJycGkSZNgb2+vespZ/fr18e6776Jbt27Yt28fXn75ZdjZ2SElJQV//PEHKleujD59+vyHo164Ll26YPr06ejatSvOnTuHypUr448//sC4cePQsmVLNGnSRFXfzc0NjRs3xmeffQY7OzvMmDEDx48fV92eq3Xr1vjiiy8wYsQINGzYECdOnMDo0aPh7++P7OzsfDGYm5ujadOmGDRoEHJzczFhwgSkpaU9tj/0ERAQABsbGyxcuBAhISGwt7eHt7e38ifwh+k7Rg4dOoR+/frhzTffRIUKFWBlZYUtW7bg0KFDhZ41NDMzw8SJE9GxY0e0bt0a7733HjIyMjBp0iTcvHkTX3755X9qq6enJ1q2bIn58+dDRNCjRw/V8rxbMUVHR+ONN97AhQsX8MUXX8DLywunTp164v2OGDECq1evRnh4OD7//HOUKlUKCxcuxJo1azBx4kQ4OTkBeDB9Izg4GB999BGys7Ph4uKC5cuX448//lBtb/Xq1ZgxYwbatWuH8uXLQ0SwbNky3Lx5E02bNn1sPLm5ucrUn0dVr14dWq0Wbdq0Ue7NXLp0aZw/fx6xsbHw9fVFhQoVAED5ITl16lR07doVlpaWCA4OhoODwxMfq4LUq1cPLi4u6N27N0aMGAFLS0ssXLgQBw8ezFc3L64JEyagRYsWMDc3R5UqVdC0aVNERkZi8ODBSEtLQ/369XHo0CGMGDEC1atXR+fOnR8bx+eff46LFy8iIiICZcuWxc2bNzF16lTV/F0inYx04RmR0Rw4cEC6du0q5cqVEysrK+XWMZ9//rlcuXJFqZeTkyMTJkyQoKAgsbS0FDc3N+nUqZNya6VHffrpp8qtnh6+ovhhK1askPDwcHF0dBStViu+vr7yxhtvyG+//abUybu9jy7Hjh2Tpk2bioODg7i4uMibb74pycnJ+a7wFhFZvny5VK5cWaysrKRcuXLy5ZdfyoABA8TFxSXfdufNmyd16tQROzs7sbGxkYCAAOnSpYvs27ev0GNZ0K25HlVYm65fvy69e/cWLy8vsbCwEF9fXxk6dKjqVj4iD65i79u3r8yYMUMCAgLE0tJSKlasKAsXLlTVy8jIkI8++kjKlCkj1tbWUqNGDVmxYkW+uw/kXQE+YcIEGTVqlJQtW1asrKykevXqyq3QHm2noXczEBFZvHixVKxYUSwtLVX9pOvWXCKPHyP//POPREVFScWKFcXOzk7s7e2lSpUqMmXKFMnOztZ5jB/dfp06dcTa2lrs7OwkIiJCdu7cqapj6N0M8qxcuVIASKlSpfL1n4jIl19+KX5+fqLVaiUkJERmz56t8zjk9bUuusb64cOHpU2bNuLk5CRWVlZStWpVnXeQOHnypDRr1kwcHR2ldOnS0r9/f1mzZo3qbgbHjx+Xt99+WwICAsTGxkacnJykdu3aEhcX99j2F3Y3AwBy6tQpERGZPHmy1KtXT9zc3JTPZ48ePeTcuXOq7Q0dOlS8vb3FzMxMFWNBdzOYNGmSav2C+lHX53bXrl0SFhYmtra2Urp0aenZs6fs378/3904MjIypGfPnlK6dGnRaDSqz8W9e/dk8ODB4uvrK5aWluLl5SV9+vSRf//9V7V/X19fadWqVb7jt3r1amnRooWUKVNGrKysxN3dXVq2bKncNo6oIBqRh/7eRETPrKysLFSrVg1lypTBxo0bjR0OERFRkeA0A6JnVI8ePdC0aVN4eXkhNTUVs2bNQlJSkvKUISIiomcBk1miZ9Tt27fx0Ucf4erVq7C0tESNGjWwdu3afPNQiYiITBmnGRARERGRyeKtuYiIiIjIZDGZJSIiIiKTxWSWiIiIiEzWc3cBWG5uLi5fvgwHBwe9H5dIRERERE+PiOD27dvw9vaGmVnh516fu2T28uXL8PHxMXYYRERERPQYFy5cUB5XXpDnLpnNexTghQsX4OjoaORoiIiIiOhRaWlp8PHx0esRzs9dMps3tcDR0ZHJLBEREVEJps+UUF4ARkREREQmi8ksEREREZksJrNEREREZLKYzBIRERGRyWIyS0REREQmi8ksEREREZksJrNEREREZLKYzBIRERGRyWIyS0REREQmi8ksEREREZksJrNEREREZLKMmsz+/vvvaNOmDby9vaHRaLBixYrHrrN9+3bUrFkT1tbWKF++PGbNmlX8gRIRERFRiWTUZPbOnTuoWrUqvvnmG73qnz17Fi1btkSDBg2QmJiITz/9FAMGDMDSpUuLOVIiIiIiKoksjLnzFi1aoEWLFnrXnzVrFsqVK4fY2FgAQEhICPbt24evvvoKr7/+ejFFSUREREQllUnNmd29ezeaNWumKouMjMS+ffuQlZWlc52MjAykpaWpXkRERET0bDDqmVlDpaamwsPDQ1Xm4eGB7OxsXLt2DV5eXvnWGT9+PEaNGvW0Qnzm9Zi91dghYG6vcGOHQPRM4uebiEyRSZ2ZBQCNRqN6LyI6y/MMHToUt27dUl4XLlwo9hiJiIiI6OkwqTOznp6eSE1NVZVduXIFFhYWcHV11bmOVquFVqt9GuERERER0VNmUmdmw8LCsGnTJlXZxo0bUatWLVhaWhopKiIiIiIyFqMms+np6Thw4AAOHDgA4MGttw4cOIDk5GQAD6YIdOnSRanfu3dvnD9/HoMGDUJSUhLmzZuHuXPn4qOPPjJG+ERERERkZEadZrBv3z6Eh/9vsv+gQYMAAF27dkVcXBxSUlKUxBYA/P39sXbtWgwcOBDTp0+Ht7c3pk2bxttyERERET2njJrMNmrUSLmAS5e4uLh8ZQ0bNsT+/fuLMSoiIiIiMhUmdQGYMfm2jzF2CDi/ZJCxQyB6JvHzTURkukzqAjAiIiIioocxmSUiIiIik8VkloiIiIhMFpNZIiIiIjJZTGaJiIiIyGQxmSUiIiIik8VkloiIiIhMFu8za0IaDI03dggIKu9u7BDoKeoxe6uxQ8DcXuGPr/QM4Oeb6NnFz/cDxfV9zjOzRERERGSymMwSERERkcliMktEREREJovJLBERERGZLCazRERERGSymMwSERERkcliMktEREREJovJLBERERGZLCazRERERGSymMwSERERkcliMktEREREJovJLBERERGZLCazRERERGSyLIwdABEREVFx8W0fY+wQUC6gjLFDeKbxzCwRERERmSwms0RERERkspjMEhEREZHJYjJLRERERCaLySwRERERmSwms0RERERkspjMEhEREZHJ4n1miXQoCfclPL9kkLFDICIiKvF4ZpaIiIiITBaTWSIiIiIyWUxmiYiIiMhkMZklIiIiIpPFZJaIiIiITBaTWSIiIiIyWbw1F1EJ1WBovLFDQFB5d2OHQEREVCiemSUiIiIik8VkloiIiIhMFpNZIiIiIjJZTGaJiIiIyGQxmSUiIiIik8VkloiIiIhMFpNZIiIiIjJZTGaJiIiIyGQxmSUiIiIik8VkloiIiIhMFpNZIiIiIjJZTGaJiIiIyGQxmSUiIiIik8VkloiIiIhMFpNZIiIiIjJZTGaJiIiIyGQxmSUiIiIik8VkloiIiIhMFpNZIiIiIjJZTGaJiIiIyGQxmSUiIiIik8VkloiIiIhMFpNZIiIiIjJZTGaJiIiIyGQxmSUiIiIik2X0ZHbGjBnw9/eHtbU1atasiR07dhRaf+HChahatSpsbW3h5eWFbt264fr1608pWiIiIiIqSYyazMbHx+ODDz7AsGHDkJiYiAYNGqBFixZITk7WWf+PP/5Aly5d0KNHDxw9ehQ///wz9u7di549ez7lyImIiIioJDBqMhsTE4MePXqgZ8+eCAkJQWxsLHx8fDBz5kyd9ffs2QM/Pz8MGDAA/v7+eOmll/Dee+9h3759Be4jIyMDaWlpqhcRERERPRuMlsxmZmYiISEBzZo1U5U3a9YMu3bt0rlOvXr1cPHiRaxduxYign/++Qe//PILWrVqVeB+xo8fDycnJ+Xl4+NTpO0gIiIiIuMxWjJ77do15OTkwMPDQ1Xu4eGB1NRUnevUq1cPCxcuRIcOHWBlZQVPT084Ozvj66+/LnA/Q4cOxa1bt5TXhQsXirQdRERERGQ8Rr8ATKPRqN6LSL6yPMeOHcOAAQPw+eefIyEhAevXr8fZs2fRu3fvArev1Wrh6OioehERERHRs8HCWDt2c3ODubl5vrOwV65cyXe2Ns/48eNRv359fPzxxwCAKlWqwM7ODg0aNMCYMWPg5eVV7HETERERUclhtDOzVlZWqFmzJjZt2qQq37RpE+rVq6dznbt378LMTB2yubk5gAdndImIiIjo+WLUaQaDBg3CnDlzMG/ePCQlJWHgwIFITk5Wpg0MHToUXbp0Ueq3adMGy5Ytw8yZM3HmzBns3LkTAwYMQO3ateHt7W2sZhARERGRkRhtmgEAdOjQAdevX8fo0aORkpKCSpUqYe3atfD19QUApKSkqO45GxUVhdu3b+Obb77Bhx9+CGdnZzRu3BgTJkwwVhOIiIiIyIiMmswCQHR0NKKjo3Uui4uLy1fWv39/9O/fv5ijIiIiIiJTYPS7GRARERERPSkms0RERERkspjMEhEREZHJYjJLRERERCaLySwRERERmSwms0RERERkspjMEhEREZHJYjJLRERERCaLySwRERERmSwms0RERERkspjMEhEREZHJYjJLRERERCaLySwRERERmSwms0RERERkspjMEhEREZHJYjJLRERERCaLySwRERERmSwms0RERERkspjMEhEREZHJYjJLRERERCaLySwRERERmSwms0RERERksgxOZkeOHInz588XRyxERERERAYxOJn99ddfERAQgIiICCxatAj3798vjriIiIiIiB7L4GQ2ISEB+/fvR5UqVTBw4EB4eXmhT58+2Lt3b3HER0RERERUoCeaM1ulShVMmTIFly5dwrx583Dp0iXUr18flStXxtSpU3Hr1q2ijpOIiIiIKJ//dAFYbm4uMjMzkZGRARFBqVKlMHPmTPj4+CA+Pr6oYiQiIiIi0umJktmEhAT069cPXl5eGDhwIKpXr46kpCRs374dx48fx4gRIzBgwICijpWIiIiISMXgZLZKlSqoW7cuzp49i7lz5+LChQv48ssvERgYqNTp0qULrl69WqSBEhERERE9ysLQFd588010794dZcqUKbBO6dKlkZub+58CIyIiIiJ6HIPPzIoIXFxc8pXfu3cPo0ePLpKgiIiIiIj0YXAyO2rUKKSnp+crv3v3LkaNGlUkQRERERER6eOJzsxqNJp85QcPHkSpUqWKJCgiIiIiIn3oPWfWxcUFGo0GGo0GQUFBqoQ2JycH6enp6N27d7EESURERESki97JbGxsLEQE3bt3x6hRo+Dk5KQss7Kygp+fH8LCwoolSCIiIiIiXfROZrt27QoA8Pf3R7169WBpaVlsQRERERER6UOvZDYtLQ2Ojo4AgOrVq+PevXu4d++ezrp59YiIiIiIipteyayLiwtSUlLg7u4OZ2dnnReA5V0YlpOTU+RBEhERERHpolcyu2XLFuVOBVu3bi3WgIiIiIiI9KVXMtuwYUOd/yYiIiIiMia9ktlDhw7pvcEqVao8cTBERERERIbQK5mtVq0aNBoNRKTQepwzS0RERERPk17J7NmzZ4s7DiIiIiIig+mVzPr6+hZ3HEREREREBtMrmV21ahVatGgBS0tLrFq1qtC6r7zySpEERkRE9KxqMDTe2CFgx/gOxg6BqEjolcy2a9cOqampcHd3R7t27QqsxzmzRERERPQ06ZXM5ubm6vw3EREREZExmRk7ACIiIiKiJ/VEyezmzZvRunVrBAQEIDAwEK1bt8Zvv/1W1LERERERERXK4GT2m2++QfPmzeHg4ID3338fAwYMgKOjI1q2bIlvvvmmOGIkIiIiItJJrzmzDxs/fjymTJmCfv36KWUDBgxA/fr1MXbsWFU5EREREVFxMjiZTUtLQ/PmzfOVN2vWDIMHDy6SoIiIiIqLb/sYY4eAcgFljB0C0TPD4GkGr7zyCpYvX56vfOXKlWjTpk2RBEVEREREpA+9zsxOmzZN+XdISAjGjh2Lbdu2ISwsDACwZ88e7Ny5Ex9++GHxRElEREREpINeyeyUKVNU711cXHDs2DEcO3ZMKXN2dsa8efMwfPjwoo2QiIiIiKgAeiWzZ8+eLe44iIiIiIgMxocmEBEREZHJMvhuBgBw8eJFrFq1CsnJycjMzFQti4kx/lWiRERERPR8MDiZ3bx5M1555RX4+/vjxIkTqFSpEs6dOwcRQY0aNYojRiIiIipiPWZvNXYImNsr3Ngh0DPA4GkGQ4cOxYcffogjR47A2toaS5cuxYULF9CwYUO8+eabxREjEREREZFOBiezSUlJ6Nq1KwDAwsIC9+7dg729PUaPHo0JEyYUeYBERERERAUxOJm1s7NDRkYGAMDb2xunT59Wll27dq3oIiMiIiIiegyDk9m6deti586dAIBWrVrhww8/xNixY9G9e3fUrVvX4ABmzJgBf39/WFtbo2bNmtixY0eh9TMyMjBs2DD4+vpCq9UiICAA8+bNM3i/RERERGT6DL4ALCYmBunp6QCAkSNHIj09HfHx8QgMDMz3cIXHiY+PxwcffIAZM2agfv36+Pbbb9GiRQscO3YM5cqV07lO+/bt8c8//2Du3LkIDAzElStXkJ2dbWgziIiIiOgZYHAyW758eeXftra2mDFjxhPvPCYmBj169EDPnj0BALGxsdiwYQNmzpyJ8ePH56u/fv16bN++HWfOnEGpUqUAAH5+fk+8fyIiIiIybU/80IR9+/bhhx9+wI8//oiEhASD18/MzERCQgKaNWumKm/WrBl27dqlc51Vq1ahVq1amDhxIsqUKYOgoCB89NFHuHfvXoH7ycjIQFpamupFRERERM8Gg8/MXrx4EW+//TZ27twJZ2dnAMDNmzdRr149LF68GD4+Pnpt59q1a8jJyYGHh4eq3MPDA6mpqTrXOXPmDP744w9YW1tj+fLluHbtGqKjo3Hjxo0C582OHz8eo0aN0r+BRERERGQyDD4z2717d2RlZSEpKQk3btzAjRs3kJSUBBFBjx49DA5Ao9Go3otIvrI8ubm50Gg0WLhwIWrXro2WLVsiJiYGcXFxBZ6dHTp0KG7duqW8Lly4YHCMRERERFQyGXxmdseOHdi1axeCg4OVsuDgYHz99deoX7++3ttxc3ODubl5vrOwV65cyXe2No+XlxfKlCkDJycnpSwkJAQigosXL6JChQr51tFqtdBqtXrHRURERESmw+Azs+XKlUNWVla+8uzsbJQpU0bv7VhZWaFmzZrYtGmTqnzTpk2oV6+eznXq16+Py5cvK3dTAICTJ0/CzMwMZcuW1XvfRERERPRsMDiZnThxIvr37499+/ZBRAA8uBjs/fffx1dffWXQtgYNGoQ5c+Zg3rx5SEpKwsCBA5GcnIzevXsDeDBFoEuXLkr9d955B66urujWrRuOHTuG33//HR9//DG6d+8OGxsbQ5tCRERERCZOr2kGLi4uqnmsd+7cQZ06dWBh8WD17OxsWFhYoHv37mjXrp3eO+/QoQOuX7+O0aNHIyUlBZUqVcLatWvh6+sLAEhJSUFycrJS397eHps2bUL//v1Rq1YtuLq6on379hgzZoze+yQiIiKiZ4deyWxsbGyxBRAdHY3o6Gidy+Li4vKVVaxYMd/UBCIiIiJ6PumVzHbt2rW44yAiIiIiMpjBdzMAgJycHKxYsQJJSUnQaDQIDQ3FK6+8AnNz86KOj4iIiIioQAYns3///TdatmyJS5cuITg4GCKCkydPwsfHB2vWrEFAQEBxxElERERElI/BdzMYMGAAAgICcOHCBezfvx+JiYlITk6Gv78/BgwYUBwxEhERERHpZPCZ2e3bt2PPnj0oVaqUUubq6oovv/zSoIcmEBERERH9VwafmdVqtbh9+3a+8vT0dFhZWRVJUERERERE+jA4mW3dujXeffdd/PnnnxARiAj27NmD3r1745VXXimOGImIiIiIdDI4mZ02bRoCAgIQFhYGa2trWFtbo379+ggMDMTUqVOLI0YiIiIiIp0MmjMrIrh16xYWL16My5cvIykpCSKC0NBQBAYGFleMREREREQ6GZzMVqhQAUePHkWFChWYwBIRERGRURk0zcDMzAwVKlTA9evXiyseIiIiIiK9GTxnduLEifj4449x5MiR4oiHiIiIiEhvBt9ntlOnTrh79y6qVq0KKysr2NjYqJbfuHGjyIIjIiIiIiqMwcnslClToNFoiiMWIiIiIiKDGJzMvv3228jOzoadnV1xxENEREREpDe958xeu3YNrVq1gr29PRwdHVGvXj2cOXOmOGMjIiIiIiqU3sns0KFDkZCQgFGjRmHSpEm4du0a3nvvveKMjYiIiIioUHpPM9iwYQPmzZuHli1bAgBatmyJSpUqISsrC5aWlsUWIBERERFRQfQ+M3v58mVUr15deV+xYkVYWVnh8uXLxRIYEREREdHj6J3MiggsLNQnci0sLJCbm1vkQRERERER6UPvaQYigoiICFVCe/fuXbRp0wZWVlZK2f79+4s2QiIiIiKiAuidzI4YMSJfWdu2bYs0GCIiIiIiQ/ynZJaIiIiIyJj0njNLRERERFTSMJklIiIiIpPFZJaIiIiITBaTWSIiIiIyWf8pmb1//35RxUFEREREZDCDk9nc3Fx88cUXKFOmDOzt7XHmzBkAwGeffYa5c+cWeYBERERERAUxOJkdM2YM4uLiMHHiRNXDEipXrow5c+YUaXBERERERIUxOJldsGABvvvuO3Ts2BHm5uZKeZUqVXD8+PEiDY6IiIiIqDAGJ7OXLl1CYGBgvvLc3FxkZWUVSVBERERERPowOJl94YUXsGPHjnzlP//8M6pXr14kQRERERER6UPvx9nmGTFiBDp37oxLly4hNzcXy5Ytw4kTJ7BgwQKsXr26OGIkIiIiItLJ4DOzbdq0QXx8PNauXQuNRoPPP/8cSUlJ+PXXX9G0adPiiJGIiIiISCeDz8wCQGRkJCIjI4s6FiIiIiIig/AJYERERERksgw+M+vi4gKNRpOvXKPRwNraGoGBgYiKikK3bt2KJEAiIiIiooIYnMx+/vnnGDt2LFq0aIHatWtDRLB3716sX78effv2xdmzZ9GnTx9kZ2ejV69exREzERERERGAJ0hm//jjD4wZMwa9e/dWlX/77bfYuHEjli5diipVqmDatGlMZomIiIioWBk8Z3bDhg1o0qRJvvKIiAhs2LABANCyZUucOXPmv0dHRERERFQIg5PZUqVK4ddff81X/uuvv6JUqVIAgDt37sDBweG/R0dEREREVAiDpxl89tln6NOnD7Zu3YratWtDo9Hgr7/+wtq1azFr1iwAwKZNm9CwYcMiD5aIiIiI6GEGJ7O9evVCaGgovvnmGyxbtgwigooVK2L79u2oV68eAODDDz8s8kCJiIiIiB71RA9NqF+/PurXr1/UsRARERERGeSJktk89+7dQ1ZWlqrM0dHxPwVERERERKQvgy8Au3v3Lvr16wd3d3fY29vDxcVF9SIiIiIieloMTmY//vhjbNmyBTNmzIBWq8WcOXMwatQoeHt7Y8GCBcURIxERERGRTgZPM/j111+xYMECNGrUCN27d0eDBg0QGBgIX19fLFy4EB07diyOOImIiIiI8jH4zOyNGzfg7+8P4MH82Bs3bgAAXnrpJfz+++9FGx0RERERUSEMTmbLly+Pc+fOAQBCQ0OxZMkSAA/O2Do7OxdlbEREREREhTI4me3WrRsOHjwIABg6dKgyd3bgwIH4+OOPizxAIiIiIqKCGDxnduDAgcq/w8PDcfz4cezbtw8BAQGoWrVqkQZHRERERFQYg87MZmVlITw8HCdPnlTKypUrh9dee42JLBERERE9dQYls5aWljhy5Ag0Gk1xxUNEREREpDeD58x26dIFc+fOLY5YiIiIiIgMYvCc2czMTMyZMwebNm1CrVq1YGdnp1oeExNTZMERERERERXG4GT2yJEjqFGjBgCo5s4C4PQDIiIiInqqDE5mt27dWhxxEBEREREZzOA5s3n+/vtvbNiwAffu3QMAiEiRBUVEREREpA+Dk9nr168jIiICQUFBaNmyJVJSUgAAPXv2xIcffljkARIRERERFcTgZHbgwIGwtLREcnIybG1tlfIOHTpg/fr1RRocEREREVFhDE5mN27ciAkTJqBs2bKq8goVKuD8+fMGBzBjxgz4+/vD2toaNWvWxI4dO/Rab+fOnbCwsEC1atUM3icRERERPRsMTmbv3LmjOiOb59q1a9BqtQZtKz4+Hh988AGGDRuGxMRENGjQAC1atEBycnKh6926dQtdunRBRESEQfsjIiIiomeLwcnsyy+/jAULFijvNRoNcnNzMWnSJISHhxu0rZiYGPTo0QM9e/ZESEgIYmNj4ePjg5kzZxa63nvvvYd33nkHYWFhhoZPRERERM8Qg2/NNWnSJDRq1Aj79u1DZmYmPvnkExw9ehQ3btzAzp079d5OZmYmEhISMGTIEFV5s2bNsGvXrgLXmz9/Pk6fPo0ff/wRY8aMeex+MjIykJGRobxPS0vTO0YiIiIiKtkMPjMbGhqKQ4cOoXbt2mjatCnu3LmD1157DYmJiQgICNB7O9euXUNOTg48PDxU5R4eHkhNTdW5zqlTpzBkyBAsXLgQFhb65eHjx4+Hk5OT8vLx8dE7RiIiIiIq2Qw+MwsAnp6eGDVqVJEE8OhTw0RE55PEcnJy8M4772DUqFEICgrSe/tDhw7FoEGDlPdpaWlMaImIiIieEQYns/7+/ujUqRM6deqE4ODgJ96xm5sbzM3N852FvXLlSr6ztQBw+/Zt7Nu3D4mJiejXrx8AIDc3FyICCwsLbNy4EY0bN863nlarNfjCNCIiIiIyDQZPM+jfvz/Wr1+PkJAQ1KxZE7GxscqDEwxhZWWFmjVrYtOmTaryTZs2oV69evnqOzo64vDhwzhw4IDy6t27N4KDg3HgwAHUqVPH4BiIiIiIyLQZnMwOGjQIe/fuxfHjx9G6dWvMnDkT5cqVQ7NmzVR3OdB3W3PmzMG8efOQlJSEgQMHIjk5Gb179wbwYIpAly5dHgRqZoZKlSqpXu7u7rC2tkalSpVgZ2dnaFOIiIiIyMQZnMzmCQoKwqhRo3DixAns2LEDV69eRbdu3QzaRocOHRAbG4vRo0ejWrVq+P3337F27Vr4+voCAFJSUh57z1kiIiIien490QVgef766y8sWrQI8fHxuHXrFt544w2DtxEdHY3o6Gidy+Li4gpdd+TIkRg5cqTB+yQiIiKiZ4PByezJkyexcOFCLFq0COfOnUN4eDi+/PJLvPbaa3BwcCiOGImIiIiIdDI4ma1YsSJq1aqFvn374q233oKnp2dxxEVERERE9FgGJ7PHjx836D6vRERERETFxeBkNi+RTUhIQFJSEjQaDUJCQlCjRo0iD46IiIiIqDAGJ7NXrlzBW2+9hW3btsHZ2Rkiglu3biE8PBw//fQTSpcuXRxxEhERERHl80QPTUhLS8PRo0dx48YN/Pvvvzhy5AjS0tIwYMCA4oiRiIiIiEgng8/Mrl+/Hr/99htCQkKUstDQUEyfPh3NmjUr0uCIiIiIiApj8JnZ3NxcWFpa5iu3tLREbm5ukQRFRERERKQPg5PZxo0b4/3338fly5eVskuXLmHgwIGIiIgo0uCIiIiIiApjcDL7zTff4Pbt2/Dz80NAQAACAwPh7++P27dv4+uvvy6OGImIiIiIdDJ4zqyPjw/279+PTZs24fjx4xARhIaGokmTJsURHxERERFRgQxOZvM0bdoUTZs2LcpYiIiIiIgMovc0gy1btiA0NBRpaWn5lt26dQsvvPACduzYUaTBEREREREVRu9kNjY2Fr169YKjo2O+ZU5OTnjvvfcQExNTpMERERERERVG72T24MGDaN68eYHLmzVrhoSEhCIJioiIiIhIH3ons//884/O+8vmsbCwwNWrV4skKCIiIiIifeidzJYpUwaHDx8ucPmhQ4fg5eVVJEEREREREelD72S2ZcuW+Pzzz3H//v18y+7du4cRI0agdevWRRocEREREVFh9L411/Dhw7Fs2TIEBQWhX79+CA4OhkajQVJSEqZPn46cnBwMGzasOGMlIiIiIlLRO5n18PDArl270KdPHwwdOhQiAgDQaDSIjIzEjBkz4OHhUWyBEhERERE9yqCHJvj6+mLt2rX4999/8ffff0NEUKFCBbi4uBRXfEREREREBXqiJ4C5uLjgxRdfLOpYiIiIiIgMovcFYEREREREJQ2TWSIiIiIyWUxmiYiIiMhkMZklIiIiIpPFZJaIiIiITBaTWSIiIiIyWUxmiYiIiMhkMZklIiIiIpPFZJaIiIiITBaTWSIiIiIyWUxmiYiIiMhkMZklIiIiIpPFZJaIiIiITBaTWSIiIiIyWUxmiYiIiMhkMZklIiIiIpPFZJaIiIiITBaTWSIiIiIyWUxmiYiIiMhkMZklIiIiIpPFZJaIiIiITBaTWSIiIiIyWUxmiYiIiMhkMZklIiIiIpPFZJaIiIiITBaTWSIiIiIyWUxmiYiIiMhkMZklIiIiIpPFZJaIiIiITBaTWSIiIiIyWUxmiYiIiMhkMZklIiIiIpPFZJaIiIiITBaTWSIiIiIyWUxmiYiIiMhkMZklIiIiIpPFZJaIiIiITBaTWSIiIiIyWUxmiYiIiMhkMZklIiIiIpPFZJaIiIiITJbRk9kZM2bA398f1tbWqFmzJnbs2FFg3WXLlqFp06YoXbo0HB0dERYWhg0bNjzFaImIiIioJDFqMhsfH48PPvgAw4YNQ2JiIho0aIAWLVogOTlZZ/3ff/8dTZs2xdq1a5GQkIDw8HC0adMGiYmJTzlyIiIiIioJLIy585iYGPTo0QM9e/YEAMTGxmLDhg2YOXMmxo8fn69+bGys6v24ceOwcuVK/Prrr6hevbrOfWRkZCAjI0N5n5aW9ti4RATZ2dnIyclRysqUstGnScXK3cHS2CHAWav/7x8BkJkjuJctxRcQERERPdeMlsxmZmYiISEBQ4YMUZU3a9YMu3bt0msbubm5uH37NkqVKlVgnfHjx2PUqFEGxZWSkoK7d++qyke9U0PvbRQXCwtzY4cASwNiEAC5IrhwKxN/Xb6PO1lMaomIiKhoGS2ZvXbtGnJycuDh4aEq9/DwQGpqql7bmDx5Mu7cuYP27dsXWGfo0KEYNGiQ8j4tLQ0+Pj466+bm5uLs2bMwNzeHt7c3rKysoNFoAACZ5lf1iqk4WVoZ/8ys1kr/ISMiyM3Jhr3tv3CzNcfS4+nIZT5LRERERcio0wwAKMliHhHJV6bL4sWLMXLkSKxcuRLu7u4F1tNqtdBqtXrFkpmZidzcXPj4+MDW1lYdp7nxE0kzCytjhwBzS0OHjBZO5ha4ey8D9lYapGUwmyUiIqKiY7Rk1s3NDebm5vnOwl65ciXf2dpHxcfHo0ePHvj555/RpEmTIo/NzMzoN3l4pmg0GmgAmEGDB5MPiIiIiIqG0bI2Kysr1KxZE5s2bVKVb9q0CfXq1StwvcWLFyMqKgqLFi1Cq1atijtMIiIiIirBjDrNYNCgQejcuTNq1aqFsLAwfPfdd0hOTkbv3r0BPJjveunSJSxYsADAg0S2S5cumDp1KurWrauc1bWxsYGTk5PR2kFERERExmHUZLZDhw64fv06Ro8ejZSUFFSqVAlr166Fr68vACAlJUV1z9lvv/0W2dnZ6Nu3L/r27auUd+3aFXFxccUaa5uhi4p1+4/6dfw7T3V/RERERKbI6JNDo6Ojce7cOWRkZCAhIQEvv/yysiwuLg7btm1T3m/btg0iku9V3ImsqfnumykI8XHFuJGfKmUigm9iJuDlmqGoFlgGXd58BadOHFetl5mRgTGfDUZYlQqoEeSD6G4dkZpy6WmHT0RERKQ3oyezVLQOH9iPJYsWIDjkBVX5nJnTEDd7BoaPmYAlq3+DW2l39HjnNdxJv63UGTfyU/y2fg0mT5+NH5etwd276egT9Y7q4RFEREREJQmT2WfInTvp+HhAb4yeMAWOTs5KuYhgwdxv8V7/QWjWog2CKobgyynTcf/+PaxesRQAcDstDcviF+KTz75AvQaNEFqpCiZOnYWTx49h947tRmoRERERUeGYzD5Dvhj+CRo2bop6DRqpyi8mn8e1K/+g/svhSpmVVosX69RDYsJfAICjhw8gKytLVcfd0wsVgkOUOkREREQljdEfmkBFY83KZTh2+BB+Xv1bvmXXrl4BALi5lVaVu5Z2x+WLFx7UuXIFllZWcHJ2VtdxK62sT0RERFTSMJl9BqRcvoTxIz/FnIW/QGttXXDFJ3jamkC/J7IRERERGQOnGTwDjh46gOvXruKNlo1Ryc8dlfzcsXfPTvw47ztU8nOH6/+fkX30DOuNa1fhWvrBMjd3d2RlZuLWzZuP1LmmrE9ERERU0jCZfQaEvfQyVm76A8vWb1delapUQ+tX38Cy9dvh4+sHN3cP7NqxTVknMzMTe//cheo1awMAXqhcDZaWlqo6V/5JxakTSUodIiIiopKG0wyeAXb2DgiqGKIqs7G1g7NLKaW8S4/38N03U+DrVx6+/gH47pspsLa2Qet2rwMAHBwd8VqHjpj4xWdwdnGBk7MLJo35HEEVQxHWoOFTbxMRERGRPpjM6snUn8jVs88AZNy/j9HDP0HarZuoUq0m5ixcCjt7B6XO0BFjYWFhgYF9eiDj/n3UfellzIiZDnNzcyNGTkRERFQwJrPPqAU/r1K912g06DdoMPoNGlzgOlprawz/YgKGfzGhuMMjIiIiKhKcM0tEREREJovJLBERERGZLCazRERERGSymMwSERERkcliMktEREREJovJLBERERGZLCazRERERGSymMwSERERkcliMktEREREJotPANNT3znbn+r+pvds+FT3R0RERGSKeGb2GbF3zy706fYOXq4ZihAfV/y2fk2+OqdPnUB0t454MdQPNSuWQ4dXmuHypYvK8syMDIz5bDDCqlRAjSAfRHfriNSUS6pt3Lr5LwZGv4sqAT6oEuCDgdHvIu3WzeJuHhEREZFOTGafEffu3UVwyAsYPmaCzuXJ586i42ut4B9YAd8vWYUVG35Hn/c/glarVeqMG/kpflu/BpOnz8aPy9bg7t109Il6Bzk5OUqd93v3RNLRw4j7aSniflqKpKOHMbDve8XePiIiIiJdOM3gGfFyeBO8HN6kwOWxE8fi5cZN8PGwkUqZj6+f8u/baWlYFr8QX8bORL0GjQAAE6fOQnidKti9YzteatQYp0+dwPYtv2HZus2oXrMWAGD85Gl4rWUTnP77FAICKxRH04iIiIgKxDOzz4Hc3Fxs37IRfv6B6NnxDdSvFowObZqqpiIcPXwAWVlZqP9yuFLm7umFCsEhSEz4CwBwIGEfHBydlEQWAKrXehEOjk7Yv/fPp9cgIiIiov/HZPY5cP3aVdy9cwdzZkzFS40iMGfhL2jSvBUGvNsVf+3eCQC4duUKLK2s4OTsrFrX1a00rl298qDO1X/g6uaWb/uubm64euWfYm8HERER0aM4zeA5ILm5AIDGzVogqlcfAEDIC5WRuG8v4n+MQ+2w+gWvC4FGo1HeP/zv/1US3eVERERExYxnZp8DzqVcYWFhgYAKQary8hUqIOXyg7sZuLm7IyszE7du3lTVuXHtGlzdSj+oU9oD165ezbf969evw620e/EET0RERFQIJrPPASsrK1SqWh1nz/ytKj935jS8y/gAAF6oXA2WlpbYtWObsvzKP6k4dSIJ1WvWBgBUq1kLt9Nu4cD+BKVOYsI+3E67hRov1in2dhARERE9itMMnhF37qQj+dxZ5f3FC8lIOnoYTs4u8C5TFt3f64cP+/ZErTr1UCfsJfyxfTO2/bYB3y9ZBQBwcHTEax06YuIXn8HZxQVOzi6YNOZzBFUMRViDBw9wCKgQjIaNm2DooAEY91UsAGDoh++jcbPmvJMBERERGQWTWT2V9CdyHT10AF3bt1XeTxg9HADQ7o23MH7KdDRt0Rojxk3Gd9NjMe7zofAPCMTUb+NQs3ZdZZ2hI8bCwsICA/v0QMb9+6j70suYETMd5ubmSp3YmbMx8tPB6NL+VQBARGQLjP5y0lNqJREREZEak9lnRO2wl5B04XqhdV5/qyNef6tjgcu11tYY/sUEDP9C94MXAMDZpRRiZ85+4jiJiIiIihLnzBIRERGRyWIyS0REREQmi8ksEREREZksJrNEREREZLKYzBIRERGRyWIyS0REREQmi8ksEREREZksJrNEREREZLKYzBIRERGRyeITwPQ0dUvSU93f+41Dnur+iIiIiEwRz8w+Q/5JuYxPBryHupUDUb1CWbwa2RBHDx3QWXfEkEEI8XHF93NmqcozMzIw5rPBCKtSATWCfBDdrSNSUy6p6ty6+S8GRr+LKgE+qBLgg4HR7yLt1s1iahURERFRwZjMPiNu3byJd15rCQsLS3y3IB6rt+zCJ5+NhoOjU766v61fg0OJCXD38My3bNzIT/Hb+jWYPH02fly2BnfvpqNP1DvIyclR6rzfuyeSjh5G3E9LEffTUiQdPYyBfd8r1vYRERER6cJpBs+IOTOnwsurDMbFfKOUlfEpl6/ePymXMeazwZj94y/oHfWWatnttDQsi1+IL2Nnol6DRgCAiVNnIbxOFezesR0vNWqMv0+ewPYtv2HZus2oXrMWAGD85Gl4rWUTnP77FAICKxRfI4mIiIgewTOzz4itm9bjhSrV8EHvbqhfLRivNW+EJYsWqOrk5uZi8Ad90L13f1QIrphvG0cPH0BWVhbqvxyulLl7eqFCcAgSE/4CAOzf9xccHJ2URBYAqtd6EQ6OTti/989iah0RERGRbkxmnxEXks/jpx/nw9evPGb/+DM6dIrCuM+HYsUvPyl15syYCnNzC3Tu/q7ObVy7cgWWVlZwcnZWlbu6lca1q1cAAFev/ANXN7d867q6ueHqlX+KrkFEREREeuA0g2eE5ObihSrVMHDIZwCA0EpV8PfJ4/jph/lo98ZbOHroAH6Y9x2Wrt0CjUZj2LYhqnV0ri9i8HaJiIiI/iuemX1GuLl7IKBCsKqsfGAQUi5dBADs+2sPrl+7isZ1q6KSnzsq+bnj8sULmPjFZ4gIq/b/23BHVmYmbt28qdrOjWvX4OpWGgBQ2t0D165ezbf/69evw620e9E3jIiIiKgQTGafETVq1cG503+rys6dOQ3vsj4AgFdeb48VG3dg2frtysvdwxPde/fDnB9/BgC8ULkaLC0tsWvHNmUbV/5JxakTSahes/b/76c2bqfdwoH9CUqdxIR9uJ12CzVerFO8jSQiIiJ6BKcZPCO69uyNd15tgW+/jkHz1u1w+MB+/LxoAUZNiAEAuLiUgotLKdU6FpaWcCvtAf+AB3cgcHB0xGsdOmLiF5/B2cUFTs4umDTmcwRVDEVYg4YAgMCgYDRs3ARDBw3AuK9iAQBDP3wfjZs1550MiIiI6KljMqunkv5ErsrVamDa7AWY8uUXmDH1K5T1KYchI8eizatvGrSdoSPGwsLCAgP79EDG/fuo+9LLmBEzHebm5kqd2JmzMfLTwejS/lUAQERkC4z+clKRtoeIiIhIH0xmnyHhTSIR3iRS7/qbdx/IV6a1tsbwLyZg+BcTClzP2aUUYmfOfpIQiYiIiIoU58wSERERkcliMktEREREJovJLBERERGZLCazRERERGSymMzqICLGDuGZIiIQADyqREREVNSYzD7E0tISAHD37l0jR/JsycnMQE5uLu5m5Ro7FCIiInrG8NZcDzE3N4ezszOuXLkCALC1tYVGowEASE6WMUMDAORmGzsCIMdM/4RURJCTmYEb16/h2NX7YC5LRERERY3J7CM8PT0BQElo81y5mmaMcFQsLMwfX6mYWRoQgwDIyc3Fsav3kZCaWXxBERER0XOLyewjNBoNvLy84O7ujqys/52N7TF1vhGjesCrnIexQ4BfWVe96wqAu1m5PCNLRERExYbJbAHMzc1Vj3C9dOOeEaN5wNzF+FMdnDOYmRIREVHJYfQLwGbMmAF/f39YW1ujZs2a2LFjR6H1t2/fjpo1a8La2hrly5fHrFmznlKkRERERFTSGDWZjY+PxwcffIBhw4YhMTERDRo0QIsWLZCcnKyz/tmzZ9GyZUs0aNAAiYmJ+PTTTzFgwAAsXbr0KUdORERERCWBUZPZmJgY9OjRAz179kRISAhiY2Ph4+ODmTNn6qw/a9YslCtXDrGxsQgJCUHPnj3RvXt3fPXVV085ciIiIiIqCYw2ZzYzMxMJCQkYMmSIqrxZs2bYtWuXznV2796NZs2aqcoiIyMxd+5cZGVlKfeJfVhGRgYyMjKU97du3QIApKUZdneC3Kz7BtUvDtkZxr//bea9O8YOweC+exLs7wfY308P+/sB9vfTw/5+etjfDxjS33l19XqQlRjJpUuXBIDs3LlTVT527FgJCgrSuU6FChVk7NixqrKdO3cKALl8+bLOdUaMGJH38Cm++OKLL7744osvvkzodeHChcfmlEa/m0HeQwnyiEi+ssfV11WeZ+jQoRg0aJDyPjc3Fzdu3ICrq2uh+3nWpKWlwcfHBxcuXICjo6Oxw6Fixv5+vrC/ny/s7+fL89rfIoLbt2/D29v7sXWNlsy6ubnB3NwcqampqvIrV67Aw0P3/VQ9PT111rewsICrq+77n2q1Wmi1WlWZs7Pzkwdu4hwdHZ+rD8Pzjv39fGF/P1/Y38+X57G/nZyc9KpntAvArKysULNmTWzatElVvmnTJtSrV0/nOmFhYfnqb9y4EbVq1dI5X5aIiIiInm1GvZvBoEGDMGfOHMybNw9JSUkYOHAgkpOT0bt3bwAPpgh06dJFqd+7d2+cP38egwYNQlJSEubNm4e5c+fio48+MlYTiIiIiMiIjDpntkOHDrh+/TpGjx6NlJQUVKpUCWvXroWvry8AICUlRXXPWX9/f6xduxYDBw7E9OnT4e3tjWnTpuH11183VhNMhlarxYgRI/JNuaBnE/v7+cL+fr6wv58v7O/H04joc88DIiIiIqKSx+iPsyUiIiIielJMZomIiIjIZDGZJSIiIiKTxWT2GdOoUSN88MEHhdbx8/NDbGzsU4mHihf72zSZUr/pEysVjv39/Dp37hw0Gg0OHDhQ4rZdnLE9bUxmi1lUVBQ0Gg00Gg0sLS3h4eGBpk2bYt68ecjNzTV2eE8sr0179uxRlWdkZChPV9u2bZuq/ooVK55ukEbA/v5ffVPq72ex37Zt26a0SaPRwNXVFY0bN8bOnTuNHVo+T3u8sL+Nqzj7+8qVK3jvvfdQrlw5aLVaeHp6IjIyErt3734q+y9pCvtx4uPjo9xJytQxmX0KmjdvjpSUFJw7dw7r1q1DeHg43n//fbRu3RrZ2dnGDu+J+fj4YP78+aqy5cuXw97e3kgRlQzsb9P0rPbbiRMnkJKSgm3btqF06dJo1aoVrly5YuywjI79/Wx6/fXXcfDgQXz//fc4efIkVq1ahUaNGuHGjRvGDu2JZWZmFst2zc3N4enpCQsLo96lFTk5Of/5RyST2acg79dhmTJlUKNGDXz66adYuXIl1q1bh7i4OKVecnIy2rZtC3t7ezg6OqJ9+/b4559/lOVRUVFo166datsffPABGjVqpCrLzs5Gv3794OzsDFdXVwwfPhyF3YHt1q1bePfdd+Hu7g5HR0c0btwYBw8efGy7unbtip9++gn37t1TyubNm4euXbs+dt1nGfvbNJXkfuvevTtat26db31PT0/Mmzev0Ha5u7vD09MTlStXxvDhw3Hr1i38+eefyvLt27ejdu3a0Gq18PLywpAhQ/Ilc4+LVdeZLmdnZ+W4ZWZmol+/fvDy8oK1tTX8/Pwwfvx4AA/+vA4Ar776KjQajfK+uLG/n73+vnnzJv744w9MmDAB4eHh8PX1Re3atTF06FC0atWq0P2fPn0abdu2hYeHB+zt7fHiiy/it99+U23fz88P48aNQ/fu3eHg4IBy5crhu+++U9X566+/UL16dVhbW6NWrVpITExULc/JyUGPHj3g7+8PGxsbBAcHY+rUqao6eWNq/Pjx8Pb2RlBQkF7bNtSj0wzyzu5v3rwZtWrVgq2tLerVq4cTJ06o1vv1119Rs2ZNWFtbo3z58hg1apRqDMXExKBy5cqws7ODj48PoqOjkZ6eriyPi4uDs7MzVq9ejdDQUGi1Wpw/f/4/tYXJrJE0btwYVatWxbJlywAAIoJ27drhxo0b2L59OzZt2oTTp0+jQ4cOBm/7+++/h4WFBf78809MmzYNU6ZMwZw5c3TWFRG0atUKqampWLt2LRISElCjRg1EREQ89pdszZo14e/vj6VLlwIALly4gN9//x2dO3c2OOZnHfvbNJWUfuvZsyfWr1+PlJQUpWzt2rVIT09H+/bt9drf3bt3lTPreY//vnTpElq2bIkXX3wRBw8exMyZMzF37lyMGTPmiWPVZdq0aVi1ahWWLFmCEydO4Mcff1SSiL179wIA5s+fj5SUFOW9MbC/DY9VF2P1t729Pezt7bFixQpkZGTorFPQ/tPT09GyZUv89ttvSExMRGRkJNq0aaN6cBMATJ48WUkko6Oj0adPHxw/fhwAcOfOHbRu3RrBwcFISEjAyJEj8z2hNDc3F2XLlsWSJUtw7NgxfP755/j000+xZMkSVb3NmzcjKSkJmzZtwurVq/XadlEZNmwYJk+ejH379sHCwgLdu3dXlm3YsAGdOnXCgAEDcOzYMXz77beIi4vD2LFjlTpmZmaYNm0ajhw5gu+//x5btmzBJ598otrH3bt3MX78eMyZMwdHjx6Fu7v7fwtaqFh17dpV2rZtq3NZhw4dJCQkRERENm7cKObm5pKcnKwsP3r0qACQv/76q8Btvf/++9KwYUPlfcOGDSUkJERyc3OVssGDByv7ERHx9fWVKVOmiIjI5s2bxdHRUe7fv6/abkBAgHz77bcFtguALF++XGJjYyU8PFxEREaNGiWvvvqq/PvvvwJAtm7dmq/+s479ra5vKkp6v4mIhIaGyoQJE5T37dq1k6ioqALbtHXrVgEgdnZ2YmdnJxqNRgBIzZo1JTMzU0REPv30UwkODlbFMX36dLG3t5ecnBy9Y9XV305OTjJ//nwREenfv780btxYtY2HPe3xwv5+dvv7l19+ERcXF7G2tpZ69erJ0KFD5eDBg0+0/9DQUPn666+V976+vtKpUyflfW5urri7u8vMmTNFROTbb7+VUqVKyZ07d5Q6M2fOFACSmJhY4H6io6Pl9ddfV9537dpVPDw8JCMjQyl70m03bNhQ3n//fZ3Lzp49q1o/bwz99ttvSp01a9YIALl3756IiDRo0EDGjRun2s4PP/wgXl5eBcawZMkScXV1Vd7Pnz9fAMiBAwcKXMdQPDNrRCICjUYDAEhKSoKPjw98fHyU5aGhoXB2dkZSUpJB261bt66yXQAICwvDqVOnkJOTk69uQkIC0tPT4erqqvyqtbe3x9mzZ3H69OnH7qtTp07YvXs3zpw5g7i4ONUvOFJjf5umktBvwIOzdXln2q5cuYI1a9bodfx37NiB/fv3Y/HixfD19UVcXJxypi4pKQlhYWGqOOrXr4/09HRcvHjxiWN9VFRUFA4cOIDg4GAMGDAAGzdu1Gs9Y2B/m3Z/v/7667h8+TJWrVqFyMhIbNu2DTVq1FBNHdHlzp07+OSTT5T+tbe3x/Hjx/Odma1SpYryb41GA09PT2VOclJSEqpWrQpbW1ulTlhYWL59zZo1C7Vq1ULp0qVhb2+P2bNn59tP5cqVYWVlpbzXd9tF4eE2enl5AYDSxoSEBIwePVr1/0evXr2QkpKCu3fvAgC2bt2Kpk2bokyZMnBwcECXLl1w/fp13LlzR9mulZWVaj//lXFn/T7nkpKS4O/vD0D9Bfqwh8vNzMzyzbHKysr6TzHk5ubCy8tLdSV6Hmdn58eu7+rqitatW6NHjx64f/8+WrRogdu3b/+nmJ5V7G/TVBL6DQC6dOmCIUOGYPfu3di9ezf8/PzQoEGDx67n7+8PZ2dnBAUF4f79+3j11Vdx5MgRaLVane3Ji11XOwui0WgKbXONGjVw9uxZrFu3Dr/99hvat2+PJk2a4JdfftF7H08L+/vxSnp/W1tbo2nTpmjatCk+//xz9OzZEyNGjEBUVFSB63z88cfYsGEDvvrqKwQGBsLGxgZvvPFGvouv8n4Y5NFoNMrFS48eE12WLFmCgQMHYvLkyQgLC4ODgwMmTZqkmtcMAHZ2dqr3+my7qDzcxrxxkdfG3NxcjBo1Cq+99lq+9aytrXH+/Hm0bNkSvXv3xhdffIFSpUrhjz/+QI8ePVRjxMbGxqAx9zg8M2skW7ZsweHDh/H6668DePBrPzk5GRcuXFDqHDt2DLdu3UJISAgAoHTp0qo5VAB03h/u0dsn7dmzBxUqVIC5uXm+ujVq1EBqaiosLCwQGBioerm5uenVlu7du2Pbtm3o0qWLzn0Q+9tUlZR+Ax78kGjXrh3mz5+P+fPno1u3bga3p3PnzsjNzcWMGTOU9uzatUv1H+WuXbvg4OCAMmXK6B3ro20+deqUcpYmj6OjIzp06IDZs2cjPj4eS5cuVeZpW1pa6n3Wrzixv/WL1dT6OzQ0VHVWUNf+d+zYgaioKLz66quoXLkyPD09ce7cOYP3c/DgQdVFso8eyx07dqBevXqIjo5G9erVERgYqNdfxfTZ9tNQo0YNnDhxIt//H4GBgTAzM8O+ffuQnZ2NyZMno27duggKCsLly5eLP7Aim7BAOnXt2lWaN28uKSkpcvHiRUlISJCxY8eKvb29tG7dWrKzs0Xkwdyb6tWrS4MGDSQhIUH+/PNPqVmzpmre1fr160Wj0cj3338vJ0+elM8//1wcHR3zzc2yt7eXgQMHyvHjx2XRokViZ2cns2bNUuo8PDcrNzdXXnrpJalataqsX79ezp49Kzt37pRhw4bJ3r17C2wXHppzlJubK1evXlXm9xQ0hzImJkYSExNVr9u3b/+3A1zCsL//V9+U+ruk91uejRs3ipWVlZibm8ulS5cKbVPe/Ld///1XVT5t2jRxd3eXO3fuyMWLF8XW1lb69u0rSUlJsmLFCnFzc5MRI0YYFOtbb70lISEhkpCQIHv37pXGjRuLpaWlMocyJiZGFi9eLElJSXLixAnp0aOHeHp6KvM0K1SoIH369JGUlBS5cePGY3rrv2N/P5v9fe3aNQkPD5cffvhBDh48KGfOnJElS5aIh4eHdO/eXamna//t2rWTatWqSWJiohw4cEDatGkjDg4OqvmmuvqoatWqyvG7ffu2uLm5ydtvvy1Hjx6VNWvWSGBgoGpeamxsrDg6Osr69evlxIkTMnz4cHF0dJSqVasq29Q1D1ufbevSsGFDeeedd/J9F6ekpBQ4Z/bhMZSYmCgA5OzZsyLyYLxbWFjIiBEj5MiRI3Ls2DH56aefZNiwYar6sbGxcvr0aVmwYIGUKVNGtd358+eLk5NTgTE/CSazxaxr164CQACIhYWFlC5dWpo0aSLz5s1TPth5zp8/L6+88orY2dmJg4ODvPnmm5Kamqqq8/nnn4uHh4c4OTnJwIEDpV+/fvm+NKOjo6V3797i6OgoLi4uMmTIENVE/Ec/kGlpadK/f3/x9vYWS0tL8fHxkY4dO6ouenjUw8nNowpKbnS9Hq7zLGB//6++KfW3KfSbyIPkytfXV1q2bPnYNhWU3KSnp4uLi4tycdG2bdvkxRdfFCsrK/H09JTBgwdLVlaWQbFeunRJmjVrJnZ2dlKhQgVZu3at6oKg7777TqpVqyZ2dnbi6OgoERERsn//fmX9VatWSWBgoFhYWIivr+9j2/Zfsb+fzf6+f/++DBkyRGrUqCFOTk5ia2srwcHBMnz4cLl7926h+z979qyEh4eLjY2N+Pj4yDfffJPv4qnHJbMiIrt375aqVauKlZWVVKtWTZYuXapKGO/fvy9RUVHi5OQkzs7O0qdPHxkyZMhjk1l9tq1Lw4YNdX4Xjxgx4omSWZEHCW29evXExsZGHB0dpXbt2vLdd98py2NiYsTLy0tsbGwkMjJSFixYUOzJrEbkKU7EICKiJ3b37l14e3tj3rx5Oues0bOF/U2kH14ARkRUwuXm5iI1NRWTJ0+Gk5MTXnnlFWOHRMWI/U1kGCazREQlXHJyMvz9/VG2bFnExcUZ/fGTVLzY30SG4TQDIiIiIjJZvDUXEREREZksJrNPQefOnTFu3Dhjh1Fivfjii8qz0EsC9lfJVhLGy/Xr1+Hu7m7wfSifB4cPH0bZsmVV9/U0NvZXyVUSx4up2LZtGzQaDW7evKn3On5+foiNjS22mIyFyWwxO3ToENasWYP+/fsrZY0aNYJGo4FGo4GZmRk8PDzw5ptv4vz580W67ycZ6FFRUWjXrl2RxvE4n332GYYMGaI8YcSYCuuvn376SVU3NjYWfn5+xR4Tx4taSRgv48ePR5s2beDn54dz584p/aPRaGBlZYXAwECMGTOmyJ/a8yTHW6PRYMWKFUUaR2EqV66M2rVrY8qUKU9tn4+jq7/c3d3zPb2uWrVqGDlyZLHGwvGiVhLHS1GJioqCRqNB79698y2Ljo6GRqMp9KlkpD8ms8Xsm2++wZtvvgkHBwdVed6zjC9duoSVK1fiwoUL6NSpk5GiNFxRPK4xT6tWrXDr1i1s2LChyLb5pArqL2trawwfPrxI220Ijpf/MfZ4uXfvHubOnYuePXuqyn/77TekpKTg1KlTGDVqFMaOHYt58+YZJUZDFfW47tatG2bOnFkinupVUH/dvn0bX331lZGi4nh5WEkaL0XNx8cHP/30k+rJXffv38fixYtRrlw5I0b2bGEyW4xyc3Px888/67ytiq2tLTw9PeHl5YW6deuib9++2L9/v6rO9u3bUbt2bWi1Wnh5eWHIkCHIzs5WlmdkZGDAgAFwd3eHtbU1XnrpJezduxfAg1//4eHhAAAXFxfVL8BffvkFlStXho2NDVxdXdGkSRPcuXMHI0eOxPfff4+VK1cqZw22bdumnElYsmQJGjVqBGtra/z444+4fv063n77bZQtWxa2traoXLkyFi9erGpDo0aN0K9fP/Tr1w/Ozs5wdXXF8OHDVWcgzM3N0bJly3zrPm2F9dfbb7+NW7duYfbs2YVuY+bMmQgICICVlRWCg4Pxww8/qJZrNBrMmTMHr776KmxtbVGhQgWsWrXqsbFxvJSc8bJu3TpYWFggLCxMVe7q6gpPT0/4+vqiY8eOqFevnqqPcnNzMXr0aJQtWxZarRbVqlXD+vXrVds4fPgwGjdurBzrd999F+np6QBQ4PHOzMxEv3794OXlBWtra/j5+WH8+PEAoPzl4NVXX4VGo1Hejxw5EtWqVcO8efNQvnx5aLVaiAjWr1+Pl156STn2rVu3Vj1qM69vf/rpJ9SrVw/W1tZ44YUXsG3bNlU7IiMjcf36dWzfvr0oDvl/UlB/9e/fHzExMbhy5UqB6/7777/o0qULXFxcYGtrixYtWuDUqVPK8ri4ODg7O2PDhg0ICQmBvb09mjdvnu/xtrpwvPxPSRovRa1GjRooV66camrUsmXL4OPjg+rVqytlhX0/51m7di2CgoJgY2OD8PBwndNmdu3ahZdffhk2Njbw8fHBgAEDno8pHEX6CAZSyXtyxqNPi3n0qSLXr1+XNm3aSHh4uFKW99jB6OhoSUpKkuXLl+d77OCAAQPE29tb1q5dK0ePHpWuXbuKi4uLXL9+XbKzs5Wng5w4cUJSUlLk5s2bcvnyZbGwsJCYmBg5e/asHDp0SKZPny63b9+W27dvS/v27ZXHPKakpEhGRobylBA/Pz9ZunSpnDlzRi5duiQXL16USZMmSWJiopw+fVqmTZsm5ubmsmfPHlVb7e3t5f3335fjx4/Ljz/+KLa2tqqnhYiIzJgxQ/z8/Iq2Awz0uP6KiYkRDw8PSU9PFxGRKVOmqJ5es2zZMrG0tJTp06fLiRMnZPLkyWJubi5btmxR6gCQsmXLyqJFi+TUqVMyYMAAsbe3l+vXrxcYF8dLyRov77//vjRv3lx5/+hTdERE9u7dK87OzvL9998rZTExMeLo6CiLFy+W48ePyyeffCKWlpZy8uRJERG5c+eOeHt7y2uvvSaHDx+WzZs3i7+/v3Tt2lVEpMDjPWnSJPHx8ZHff/9dzp07Jzt27JBFixaJiMiVK1cEgMyfP19SUlLkypUrIiIyYsQIsbOzk8jISNm/f78cPHhQcnNz5ZdffpGlS5fKyZMnJTExUdq0aSOVK1dWnoqV19ayZcvKL7/8IseOHZOePXuKg4ODXLt2TXWcateuLSNHjizy42+ogvpr//79Uq1aNenbt6+y7NGnOb3yyisSEhIiv//+uxw4cEAiIyMlMDBQMjMzReTBk4wsLS2lSZMmsnfvXklISJCQkBB55513CoyH46Vkj5eilPckr5iYGImIiFDKIyIiZMqUKdK2bVulvwr7fhYRSU5OFq1Wq/pu9PDwUD1Z69ChQ2Jvby9TpkyRkydPys6dO6V69eoSFRWl7FvXU8yeBUxmi9Hy5cvF3Nxc9RhAkQf/YVtaWoqdnZ3Y2toKAAkKClI9Lu7TTz+V4OBg1brTp08Xe3t7ycnJkfT0dLG0tJSFCxcqyzMzM8Xb21smTpwoIrofTZeQkCAA5Ny5czpj1vUYvbwvpNjY2Me2uWXLlvLhhx+q2hoSEqJqx+DBgyUkJES13sqVK8XMzCzfoySfpsL66/3335f79++Lr6+vjB49WkTyJ7P16tWTXr16qdZ98803VY+iBCDDhw9X3qenp4tGo5F169YVGBfHS8kaL23btlU95z2vvTY2NmJnZyeWlpYCQN59913Vet7e3jJ27FhV2YsvvijR0dEi8uARoC4uLsqPJRGRNWvWiJmZmfIDS9fx7t+/vzRu3DjfuM0DHY8iHjFihFhaWirJSkHykpvDhw+r2vrll18qdbKysqRs2bLK41LzvPrqq6r/RI2loP5KTEyU9evXi6Wlpfz9998iok5mT548KQBk586dyrrXrl0TGxsbWbJkiYg8SGYBKOuLPPjceXh4FBgPx0vJHi9FKe/4X716VbRarZw9e1bOnTsn1tbWcvXqVSWZ1ef7eejQoTq/Gx/+zu7cuXO+cbRjxw4xMzOTe/fuicizm8xymkExunfvHrRaLTQaTb5lHTt2xIEDB3Dw4EH88ccfCAwMRLNmzZQLEpKSkhAWFqZat379+khPT8fFixdx+vRpZGVloX79+spyS0tL1K5dG0lJSQXGVLVqVURERKBy5cp48803MXv2bPz77796tadWrVqq9zk5ORg7diyqVKkCV1dX2NvbY+PGjUhOTlbVq1u3rqodYWFhOHXqlGp+lI2NDXJzc5GRkaFXLMWhsP4CAK1Wi9GjR2PSpEm4du1avuVJSUmq/gAe9Nmj/VGlShXl33Z2dnBwcFD+1PnCCy/A3t4e9vb2aNGihVKP46XkjJd79+7B2to6X3l8fLzSR/Hx8Vi5ciWGDBkCAEhLS8Ply5cLHR9JSUmoWrUq7OzsVMtzc3Nx4sSJAuOJiorCgQMHEBwcjAEDBmDjxo16tcPX1xelS5dWlZ0+fRrvvPMOypcvD0dHR/j7+wNAvj56+E/2FhYWqFWrVr5xZGNjg7t37+oVS3EqqL+AB3/efumll/DZZ5/lW5aUlAQLCwvUqVNHKXN1dUVwcLCqrba2tggICFDee3l5KZ/nHTt2KJ9ne3t7LFy4UKnH8VIyx0txcHNzQ6tWrfD9999j/vz5aNWqFdzc3JTl+nw/JyUl6fxufFhCQgLi4uJUYy4yMhK5ubk4e/ZsMbfSuJjMFiM3NzfcvXsXmZmZ+ZY5OTkhMDAQgYGBqF+/PubOnYtTp04hPj4eACAi+ZIq+f95gxqNRvXvR+sUlIwBD+Ybbtq0CevWrUNoaCi+/vprBAcH6zXQH/7SBIDJkydjypQp+OSTT7BlyxYcOHAAkZGROtv7ODdu3ICtrS1sbGwMXreoFNZfeTp16gQ/Pz+MGTNG53J9+sPS0jLfOnlX5q9duxYHDhzAgQMHMGfOHKUOx4uaMceLm5ubzoTex8cHgYGBCAkJQfv27fHBBx9g8uTJuH//vlKnsONfWF8U1kc1atTA2bNn8cUXX+DevXto37493njjjce249H+AYA2bdrg+vXrmD17Nv7880/8+eefAKBXHz0a440bN/IlP8ZQUH/l+fLLLxEfH4/ExERVed5n5lGP9pOuz3PeurVq1VI+zwcOHFDNx+d4KZnjpbh0794dcXFx+P7779G9e3fVMn2+nwsajw/Lzc3Fe++9pxpzBw8exKlTp1Q/uJ5FTGaLUbVq1QAAx44de2xdc3NzAFCueAwNDcWuXbtUA3jXrl1wcHBAmTJlEBgYCCsrK/zxxx/K8qysLOzbtw8hISEAACsrKwDId4WoRqNB/fr1MWrUKCQmJsLKygrLly9X1tH3itIdO3agbdu26NSpE6pWrYry5curLo7Is2fPnnzvK1SooLQZAI4cOYIaNWrotd/iok9/mZmZYfz48Zg5c2a+yfchISGq/gAe9Flef+jD19dXSVrLlClTYD2OF+ONl+rVq+v9mc7OzkZmZiYcHR3h7e1d6PgIDQ3FgQMHVBdr7Ny5E2ZmZggKCgJQ8PF2dHREhw4dMHv2bMTHx2Pp0qW4ceMGgAfJlj59dP36dSQlJWH48OGIiIhASEhIgUngw32UnZ2NhIQEVKxYUVXnyJEjqgtcjOVx/VW7dm289tprylnRPKGhocjOzlYSNODBMTp58qTen2kbGxvl8xwYGJjvLikP43gpGeOluDRv3hyZmZnIzMxEZGSkapk+38+hoaE6vxsfVqNGDRw9elQ15vJeed/vzyoms8WodOnSqFGjRr4vJAC4e/cuUlNTkZqaioMHDyI6OhrW1tZo1qwZgAf3oLtw4QL69++P48ePY+XKlRgxYgQGDRoEMzMz2NnZoU+fPvj444+xfv16HDt2DL169cLdu3fRo0cPAA8SI41Gg9WrV+Pq1atIT0/Hn3/+iXHjxmHfvn1ITk7GsmXLcPXqVeUD4+fnh0OHDuHEiRO4du1aobdgCQwMxKZNm7Br1y4kJSXhvffeQ2pqar56Fy5cwKBBg3DixAksXrwYX3/9Nd5//31VnR07dihtN5bC+uthrVq1Qp06dfDtt9+qyj/++GPExcVh1qxZOHXqFGJiYrBs2TJ89NFH/zk2jpeSM14iIyNx9OjRfP9xX79+Hampqbh48SLWrVuHqVOnIjw8HI6OjgAejI8JEyYgPj4eJ06cwJAhQ3DgwAGlbR07doS1tTW6du2KI0eOYOvWrejfvz86d+4MDw8PALqP95QpU/DTTz/h+PHjOHnyJH7++Wd4enrC2dlZWWfz5s1ITU0t9Ayli4sLXF1d8d133+Hvv//Gli1bMGjQIJ11p0+fjuXLl+P48ePo27cv/v33X9XZpnPnzuHSpUto0qTJEx/nolJQfz1s7Nix2LJli+rP8xUqVEDbtm3Rq1cv/PHHHzh48CA6deqEMmXKoG3btv85Lo6Xkjleiou5uTmSkpKQlJSk+mEOQK/v5969e+P06dPKd+OiRYsQFxen2s7gwYOxe/du9O3bFwcOHMCpU6ewatUq1X3Tn1lPdYbuc2jWrFlSt25dVVnDhg0FgPJycXGRhg0bqq56FxHZtm2bvPjii2JlZSWenp4yePBgycrKUpbfu3dP+vfvL25ubqLVaqV+/fry119/qbYxevRo8fT0FI1GI127dpVjx45JZGSklC5dWrRarQQFBcnXX3+t1L9y5Yo0bdpU7O3tBYBs3bpV59W3Ig+uqm/btq3Y29uLu7u7DB8+XLp06aK64KBhw4YSHR0tvXv3FkdHR3FxcZEhQ4aoJrFfvHhRLC0t5cKFC096mItMQf318N0ERER27dolAFQXgIk8uMq+fPnyYmlpKUFBQbJgwQLVcui4uMLJyUnmz59fYEwcLyVvvNStW1dmzZolIv+7yCXvZW5uLmXLlpVevXqpLpjJycmRUaNGSZkyZcTS0lKqVq2a78K/Q4cOSXh4uFhbW0upUqWkV69ecvv2bWW5ruP93XffSbVq1cTOzk4cHR0lIiJC9u/fr6yzatUqCQwMFAsLC2W8jhgxQqpWrZqvXZs2bZKQkBDRarVSpUoV2bZtm2rM5rV10aJFUqdOHbGyspKQkBDZvHmzajvjxo2TyMjI/3KIi5Su/np0fL777rsCQHU3gxs3bkjnzp3FyclJbGxsJDIyUrmbgMiDC8CcnJxU21m+fLkU9l8rx0vJHy9FRdcFeA97+G4G+nw///rrrxIYGCharVYaNGgg8+bNy3fR7l9//aX0uZ2dnVSpUkV1IeGzegEYk9lidu/ePSlXrpzs2rXL2KEYha5E8FEfffRRvrsAGMvz3l/GZirjZc2aNRISEmLUu28YQ0GJ4MPu378vPj4+8scffzy9wB7jee0vYzPV8UKmx+JpnP19nllbW2PBggU6r36nB9zd3YvkT/FFgf1V8pWE8dKyZUucOnUKly5dgo+Pj1FjKWnOnz+PYcOG5bsS35jYXyVXSRwvZHqYzD4FDRs2NHYIJdrHH39s7BBU2F8lW0kZL4/O46UHgoKClAuQShL2V8lUUscLmRaNiB73eyAiIiIiKoF4NwMiIiIiMllMZomIiIjIZDGZJSIiIiKTxWSWiIiIiEwWk1kiIiIiMllMZomISqC4uDjlEaNERFQwJrNEREUgKioKGo0m36t58+aPXdfPzw+xsbGqsg4dOuDkyZPFFO3/MGkmIlPHhyYQERWR5s2bY/78+aoyrVb7RNuysbGBjY1NUYT1VOTk5ECj0cDMjOdIiOjp4rcOEVER0Wq18PT0VL1cXFwAACNHjkS5cuWg1Wrh7e2NAQMGAAAaNWqE8+fPY+DAgcrZXCD/GdORI0eiWrVqmDdvHsqVKwd7e3v06dMHOTk5mDhxIjw9PeHu7o6xY8eqYoqJiUHlypVhZ2cHHx8fREdHIz09HQCwbds2dOvWDbdu3VL2PXLkSADAv//+iy5dusDFxQW2trZo0aIFTp06pWw3L77Vq1cjNDQUWq0W58+fL65DS0RUICazRETF7JdffsGUKVPw7bff4tSpU1ixYgUqV64MAFi2bBnKli2L0aNHIyUlBSkpKQVu5/Tp01i3bh3Wr1+PxYsXY968eWjVqhUuXryI7du3Y8KECRg+fDj27NmjrGNmZoZp06bhyJEj+P7777FlyxZ88sknAIB69eohNjYWjo6Oyr4/+ugjAA+mTezbtw+rVq3C7t27ISJo2bIlsrKylG3fvXsX48ePx5w5c3D06FG4u7sXx+EjIioUpxkQERWR1atXw97eXlU2ePBg2NnZwdPTE02aNIGlpSXKlSuH2rVrAwBKlSoFc3NzODg4wNPTs9Dt5+bmYt68eXBwcEBoaCjCw8Nx4sQJrF27FmZmZggODsaECROwbds21K1bFwDwwQcfKOv7+/vjiy++QJ8+fTBjxgxYWVnByckJGo1Gte9Tp05h1apV2LlzJ+rVqwcAWLhwIXx8fLBixQq8+eabAICsrCzMmDEDVatW/c/HjojoSTGZJSIqIuHh4Zg5c6aqrFSpUrhz5w5iY2NRvnx5NG/eHC1btkSbNm1gYWHYV7Cfnx8cHByU9x4eHjA3N1fNU/Xw8MCVK1eU91u3bsW4ceNw7NgxpKWlITs7G/fv38edO3dgZ2encz9JSUmwsLBAnTp1lDJXV1cEBwcjKSlJKbOyskKVKlUMagMRUVHjNAMioiJiZ2eHwMBA1atUqVLw8fHBiRMnMH36dNjY2CA6Ohovv/yy6k/2+rC0tFS912g0Ostyc3MBAOfPn0fLli1RqVIlLF26FAkJCZg+fToAFLpvESmwPG9OL/DgIrWH3xMRGQOTWSKip8DGxgavvPIKpk2bhm3btmH37t04fPgwgAdnOHNycop8n/v27UN2djYmT56MunXrIigoCJcvX1bV0bXv0NBQZGdn488//1TKrl+/jpMnTyIkJKTI4yQi+i84zYCIqIhkZGQgNTVVVWZhYYHVq1cjJycHderUga2tLX744QfY2NjA19cXwIPpA7///jveeustaLVauLm5FUk8AQEByM7Oxtdff402bdpg586dmDVrlqqOn58f0tPTsXnzZlStWhW2traoUKEC2rZti169euHbb7+Fg4MDhgwZgjJlyqBt27ZFEhsRUVHhmVkioiKyfv16eHl5qV4vvfQSnJ2dMXv2bNSvXx9VqlTB5s2b8euvv8LV1RUAMHr0aJw7dw4BAQEoXbp0kcVTrVo1xMTEYMKECahUqRIWLlyI8ePHq+rUq1cPvXv3RocOHVC6dGlMnDgRADB//nzUrFkTrVu3RlhYGEQEa9euzTetgYjI2DRS0OQoIiIiIqISjmdmiYiIiMhkMZklIiIiIpPFZJaIiIiITBaTWSIiIiIyWUxmiYiIiMhkMZklIiIiIpPFZJaIiIiITBaTWSIiIiIyWUxmiYiIiMhkMZklIiIiIpPFZJaIiIiITNb/Aa8mEWtnG8fiAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 800x500 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import seaborn as sns\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib\n",
    "matplotlib.rcParams['figure.figsize'] = [8, 5]\n",
    "\n",
    "ax = sns.barplot(data=coverage_df, \n",
    "            x=\"Estimator\", y=\"Coverage Probability\", hue=\"Sample Size\",\n",
    "           palette = sns.color_palette(\"Blues_r\"))\n",
    "\n",
    "plt.legend(loc='lower left')\n",
    "plt.title(\"Coverage Probabilities of Various Estimators\")\n",
    "#ax.set_xticklabels(ax rotation=0, fontsize=8)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "193c17db",
   "metadata": {},
   "source": [
    "To summarize, here are my main conclusions from this simulation: \n",
    "\n",
    "1. econml is a great package for causal inference  \n",
    "2. Bootstrapping can be used to generate accurate coverage probabilities even when traditional statistical estimation fails\n",
    "3. Double Machine Learning out-performs doubly robust estimation in conditions of sparsity and low sample size"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "id": "27fd682a",
   "metadata": {},
   "outputs": [],
   "source": [
    "from scipy.stats import binom"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "id": "8b5b7a25",
   "metadata": {},
   "outputs": [],
   "source": [
    "observed_list = []\n",
    "prob_list = []\n",
    "for num_observed in range(30, 51):\n",
    "    observed_list.append(num_observed)\n",
    "    if num_observed < 0.95 * 50:\n",
    "    \n",
    "        prob_list.append(binom.cdf(num_observed, 50, .95))\n",
    "        \n",
    "    else:\n",
    "        prob_list.append(1 - binom.cdf(num_observed -1, 50, .95))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "id": "988809a5",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>probs</th>\n",
       "      <th>observed</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1.668198e-12</td>\n",
       "      <td>30</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2.058903e-11</td>\n",
       "      <td>31</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2.340397e-10</td>\n",
       "      <td>32</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>2.446165e-09</td>\n",
       "      <td>33</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>2.346135e-08</td>\n",
       "      <td>34</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>2.059933e-07</td>\n",
       "      <td>35</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>1.651038e-06</td>\n",
       "      <td>36</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>1.203974e-05</td>\n",
       "      <td>37</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>7.956627e-05</td>\n",
       "      <td>38</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>4.743368e-04</td>\n",
       "      <td>39</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>2.537013e-03</td>\n",
       "      <td>40</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>1.209576e-02</td>\n",
       "      <td>41</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>5.101349e-02</td>\n",
       "      <td>42</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13</th>\n",
       "      <td>1.885832e-01</td>\n",
       "      <td>43</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14</th>\n",
       "      <td>6.044188e-01</td>\n",
       "      <td>44</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>15</th>\n",
       "      <td>1.657869e+00</td>\n",
       "      <td>45</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>16</th>\n",
       "      <td>3.833473e+00</td>\n",
       "      <td>46</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17</th>\n",
       "      <td>7.351470e+00</td>\n",
       "      <td>47</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>18</th>\n",
       "      <td>8.648530e+00</td>\n",
       "      <td>48</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>19</th>\n",
       "      <td>4.470908e+00</td>\n",
       "      <td>49</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>20</th>\n",
       "      <td>1.231120e+00</td>\n",
       "      <td>50</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "           probs  observed\n",
       "0   1.668198e-12        30\n",
       "1   2.058903e-11        31\n",
       "2   2.340397e-10        32\n",
       "3   2.446165e-09        33\n",
       "4   2.346135e-08        34\n",
       "5   2.059933e-07        35\n",
       "6   1.651038e-06        36\n",
       "7   1.203974e-05        37\n",
       "8   7.956627e-05        38\n",
       "9   4.743368e-04        39\n",
       "10  2.537013e-03        40\n",
       "11  1.209576e-02        41\n",
       "12  5.101349e-02        42\n",
       "13  1.885832e-01        43\n",
       "14  6.044188e-01        44\n",
       "15  1.657869e+00        45\n",
       "16  3.833473e+00        46\n",
       "17  7.351470e+00        47\n",
       "18  8.648530e+00        48\n",
       "19  4.470908e+00        49\n",
       "20  1.231120e+00        50"
      ]
     },
     "execution_count": 46,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pd.DataFrame({\"probs\": np.array(prob_list) * 16, \"observed\": observed_list})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3373fa3a",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "dml_env",
   "language": "python",
   "name": "dml_env"
  },
  "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.9.16"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}