Created
March 19, 2020 07:20
-
-
Save kiwamizamurai/7b59f2b051b8740431c1b95c3fe46f8f to your computer and use it in GitHub Desktop.
DoublyRobust-AIPTW.ipynb
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| { | |
| "nbformat": 4, | |
| "nbformat_minor": 0, | |
| "metadata": { | |
| "colab": { | |
| "name": "DoublyRobust-AIPTW.ipynb", | |
| "provenance": [], | |
| "authorship_tag": "ABX9TyPPJq6/3zkjoSSADOMfJWCT", | |
| "include_colab_link": true | |
| }, | |
| "kernelspec": { | |
| "name": "python3", | |
| "display_name": "Python 3" | |
| } | |
| }, | |
| "cells": [ | |
| { | |
| "cell_type": "markdown", | |
| "metadata": { | |
| "id": "view-in-github", | |
| "colab_type": "text" | |
| }, | |
| "source": [ | |
| "<a href=\"https://colab.research.google.com/gist/kiwamizamurai/7b59f2b051b8740431c1b95c3fe46f8f/doublyrobust-aiptw.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>" | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "metadata": { | |
| "id": "Ktx7_leBce2J", | |
| "colab_type": "text" | |
| }, | |
| "source": [ | |
| "## データ\n", | |
| "\n", | |
| "- https://qiita.com/deaikei/items/df3626486986566cb65c\n", | |
| "\n", | |
| "前回と同じです\n", | |
| "\n", | |
| "- https://gist.github.com/kiwamizamurai/f3aba3e0deddc294776883da70bdd143" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "ydqx-ZMXbuJG", | |
| "colab_type": "code", | |
| "colab": {} | |
| }, | |
| "source": [ | |
| "import pandas as pd\n", | |
| "import numpy as np\n", | |
| "from sklearn.linear_model import LogisticRegression\n", | |
| "from sklearn.linear_model import LinearRegression\n", | |
| "\n", | |
| "import matplotlib.pyplot as plt\n", | |
| "plt.rcParams[\"figure.figsize\"] = (18, 9)\n", | |
| "plt.rcParams[\"font.size\"] = 25\n", | |
| "import seaborn as sns\n", | |
| "sns.set()" | |
| ], | |
| "execution_count": 0, | |
| "outputs": [] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "Y5jAf06ob43I", | |
| "colab_type": "code", | |
| "colab": {} | |
| }, | |
| "source": [ | |
| "age = [\"10y\", \"20y\", \"30y\", \"40y\", \"50y\", \"60y\", \"70y\"]\n", | |
| "gender = [\"female\", \"male\"]" | |
| ], | |
| "execution_count": 0, | |
| "outputs": [] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "3LxWLsHKb451", | |
| "colab_type": "code", | |
| "colab": {} | |
| }, | |
| "source": [ | |
| "def generate_users(n_users, **features):\n", | |
| " \"\"\"make user list where columns are 'user_id' and feature columns ('age', 'gender',... etc)\n", | |
| " values of the feature columns are randomly chosen from the given list\n", | |
| " \"\"\"\n", | |
| " user_id = [i + 1 for i in range(n_users)]\n", | |
| " df = pd.DataFrame({\"user_id\":user_id})\n", | |
| "\n", | |
| " np.random.seed(100)\n", | |
| " for col in features:\n", | |
| " feature = features[col]\n", | |
| " values = np.random.choice(feature, n_users)\n", | |
| " df[col] = values\n", | |
| "\n", | |
| " return df" | |
| ], | |
| "execution_count": 0, | |
| "outputs": [] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "6pig1dy4b49J", | |
| "colab_type": "code", | |
| "colab": {} | |
| }, | |
| "source": [ | |
| "n_users = 100000\n", | |
| "df = generate_users(n_users, age=age, gender=gender)" | |
| ], | |
| "execution_count": 0, | |
| "outputs": [] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "8xDbvs9Nb4_g", | |
| "colab_type": "code", | |
| "colab": {} | |
| }, | |
| "source": [ | |
| "def add_flag(df, **targets):\n", | |
| " \"\"\"add flag column named as 'x_flag'\n", | |
| " bias is added to the users specified in the variables '**targets'\n", | |
| " \"\"\"\n", | |
| " n_users = df.shape[0]\n", | |
| "\n", | |
| " # specify users to whom bias is added\n", | |
| " is_target = [True for _ in range(n_users)]\n", | |
| " for col in targets:\n", | |
| " target_values = targets[col]\n", | |
| " is_value = [False for _ in range(n_users)] \n", | |
| " for value in target_values:\n", | |
| " tmp_bool = (df[col] == value)\n", | |
| " is_value = np.logical_or(is_value, tmp_bool)\n", | |
| " is_target = np.logical_and(is_target, is_value)\n", | |
| "\n", | |
| " #assign flag\n", | |
| " flag = [0, 1]\n", | |
| " #for ordinary users, 10% of them are assigned to 'flag=1'\n", | |
| " weight = [0.9, 0.1]\n", | |
| " np.random.seed(200)\n", | |
| " df[\"treatment\"] = np.random.choice(flag, n_users, p=weight)\n", | |
| " #for users specified above, 90% of them are assigned to 'flag=1' \n", | |
| " weight = [0.1, 0.9]\n", | |
| " np.random.seed(300)\n", | |
| " df[\"treatment\"][is_target] = np.random.choice(flag, sum(is_target), p=weight)\n", | |
| " return df" | |
| ], | |
| "execution_count": 0, | |
| "outputs": [] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "-S7gbJ2cb5Fw", | |
| "colab_type": "code", | |
| "colab": {} | |
| }, | |
| "source": [ | |
| "def add_rate(df, **targets):\n", | |
| " \"\"\"add real valued column named as 'x_rate', where the values are normalized between 0 and 1\n", | |
| " for the users with 'flag=1', values in 'x_rate' are increased\n", | |
| " for the users specified by '**targets', values in 'x_rate' are much higher, regardless of 'flag=1 or 0'\n", | |
| " \"\"\"\n", | |
| " n_users = df.shape[0]\n", | |
| "\n", | |
| " #specify users to whom high 'x_rate' values are set\n", | |
| " is_target = [True for _ in range(n_users)]\n", | |
| " for col in targets:\n", | |
| " target_values = targets[col]\n", | |
| " is_value = [False for _ in range(n_users)] \n", | |
| " for value in target_values:\n", | |
| " tmp_bool = (df[col] == value)\n", | |
| " is_value = np.logical_or(is_value, tmp_bool)\n", | |
| " is_target = np.logical_and(is_target, is_value)\n", | |
| "\n", | |
| " #for users with 'flag=0', 'x_rate' values are sampled from Normal(0.2, 0.1)\n", | |
| " np.random.seed(200)\n", | |
| " df[\"outcome\"] = np.random.normal(0.2, 0.1, n_users)\n", | |
| "\n", | |
| " #for users with 'flag=1', 'x_rate' values are increased by 0.2 point\n", | |
| " is_flag = (df[\"treatment\"] == 1)\n", | |
| " df[\"outcome\"][is_flag] += 0.2\n", | |
| "\n", | |
| " #for users specified by '**targets', 'x_rate' values are sampled from Normal(0.8, 0.2)\n", | |
| " np.random.seed(300)\n", | |
| " df[\"outcome\"][is_target] = np.random.normal(0.8, 0.2, sum(is_target))\n", | |
| "\n", | |
| " #for convenience normalize to [0 ,1]\n", | |
| " df[\"outcome\"][df[\"outcome\"] > 1] = 1 \n", | |
| " df[\"outcome\"][df[\"outcome\"] < 0] = 0\n", | |
| "\n", | |
| " return df" | |
| ], | |
| "execution_count": 0, | |
| "outputs": [] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "3ubMHBphb5KK", | |
| "colab_type": "code", | |
| "colab": { | |
| "base_uri": "https://localhost:8080/", | |
| "height": 364 | |
| }, | |
| "outputId": "c463cb2e-e4a5-45f1-caf4-4abefcbd425f" | |
| }, | |
| "source": [ | |
| "df = add_flag(df, age=[\"20y\", \"30y\"], gender=[\"male\"])\n", | |
| "df = add_rate(df, age=[\"20y\", \"30y\"], gender=[\"male\"])" | |
| ], | |
| "execution_count": 7, | |
| "outputs": [ | |
| { | |
| "output_type": "stream", | |
| "text": [ | |
| "/usr/local/lib/python3.6/dist-packages/ipykernel_launcher.py:26: SettingWithCopyWarning: \n", | |
| "A value is trying to be set on a copy of a slice from a DataFrame\n", | |
| "\n", | |
| "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n", | |
| "/usr/local/lib/python3.6/dist-packages/ipykernel_launcher.py:24: SettingWithCopyWarning: \n", | |
| "A value is trying to be set on a copy of a slice from a DataFrame\n", | |
| "\n", | |
| "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n", | |
| "/usr/local/lib/python3.6/dist-packages/ipykernel_launcher.py:28: SettingWithCopyWarning: \n", | |
| "A value is trying to be set on a copy of a slice from a DataFrame\n", | |
| "\n", | |
| "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n", | |
| "/usr/local/lib/python3.6/dist-packages/ipykernel_launcher.py:31: SettingWithCopyWarning: \n", | |
| "A value is trying to be set on a copy of a slice from a DataFrame\n", | |
| "\n", | |
| "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n", | |
| "/usr/local/lib/python3.6/dist-packages/ipykernel_launcher.py:32: SettingWithCopyWarning: \n", | |
| "A value is trying to be set on a copy of a slice from a DataFrame\n", | |
| "\n", | |
| "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n" | |
| ], | |
| "name": "stderr" | |
| } | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "metadata": { | |
| "id": "rbNvEzsRpubB", | |
| "colab_type": "text" | |
| }, | |
| "source": [ | |
| "0.2ポイントが広告効果\n", | |
| "\n", | |
| "ATE = 0.20" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "MeJtlZXfb5MU", | |
| "colab_type": "code", | |
| "colab": { | |
| "base_uri": "https://localhost:8080/", | |
| "height": 362 | |
| }, | |
| "outputId": "42c784ba-11d2-4c07-bfde-6279d8725729" | |
| }, | |
| "source": [ | |
| "df.head(10)" | |
| ], | |
| "execution_count": 8, | |
| "outputs": [ | |
| { | |
| "output_type": "execute_result", | |
| "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>user_id</th>\n", | |
| " <th>age</th>\n", | |
| " <th>gender</th>\n", | |
| " <th>treatment</th>\n", | |
| " <th>outcome</th>\n", | |
| " </tr>\n", | |
| " </thead>\n", | |
| " <tbody>\n", | |
| " <tr>\n", | |
| " <th>0</th>\n", | |
| " <td>1</td>\n", | |
| " <td>10y</td>\n", | |
| " <td>male</td>\n", | |
| " <td>1</td>\n", | |
| " <td>0.254905</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>1</th>\n", | |
| " <td>2</td>\n", | |
| " <td>10y</td>\n", | |
| " <td>female</td>\n", | |
| " <td>0</td>\n", | |
| " <td>0.391095</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>2</th>\n", | |
| " <td>3</td>\n", | |
| " <td>40y</td>\n", | |
| " <td>female</td>\n", | |
| " <td>0</td>\n", | |
| " <td>0.271188</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>3</th>\n", | |
| " <td>4</td>\n", | |
| " <td>10y</td>\n", | |
| " <td>male</td>\n", | |
| " <td>0</td>\n", | |
| " <td>0.175226</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>4</th>\n", | |
| " <td>5</td>\n", | |
| " <td>30y</td>\n", | |
| " <td>male</td>\n", | |
| " <td>1</td>\n", | |
| " <td>0.502966</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>5</th>\n", | |
| " <td>6</td>\n", | |
| " <td>70y</td>\n", | |
| " <td>male</td>\n", | |
| " <td>0</td>\n", | |
| " <td>0.196705</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>6</th>\n", | |
| " <td>7</td>\n", | |
| " <td>50y</td>\n", | |
| " <td>male</td>\n", | |
| " <td>0</td>\n", | |
| " <td>0.177865</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>7</th>\n", | |
| " <td>8</td>\n", | |
| " <td>30y</td>\n", | |
| " <td>female</td>\n", | |
| " <td>1</td>\n", | |
| " <td>0.447726</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>8</th>\n", | |
| " <td>9</td>\n", | |
| " <td>60y</td>\n", | |
| " <td>female</td>\n", | |
| " <td>0</td>\n", | |
| " <td>0.130806</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>9</th>\n", | |
| " <td>10</td>\n", | |
| " <td>30y</td>\n", | |
| " <td>female</td>\n", | |
| " <td>1</td>\n", | |
| " <td>0.479201</td>\n", | |
| " </tr>\n", | |
| " </tbody>\n", | |
| "</table>\n", | |
| "</div>" | |
| ], | |
| "text/plain": [ | |
| " user_id age gender treatment outcome\n", | |
| "0 1 10y male 1 0.254905\n", | |
| "1 2 10y female 0 0.391095\n", | |
| "2 3 40y female 0 0.271188\n", | |
| "3 4 10y male 0 0.175226\n", | |
| "4 5 30y male 1 0.502966\n", | |
| "5 6 70y male 0 0.196705\n", | |
| "6 7 50y male 0 0.177865\n", | |
| "7 8 30y female 1 0.447726\n", | |
| "8 9 60y female 0 0.130806\n", | |
| "9 10 30y female 1 0.479201" | |
| ] | |
| }, | |
| "metadata": { | |
| "tags": [] | |
| }, | |
| "execution_count": 8 | |
| } | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "Xnecbg9hdDBe", | |
| "colab_type": "code", | |
| "colab": { | |
| "base_uri": "https://localhost:8080/", | |
| "height": 34 | |
| }, | |
| "outputId": "aae1bb2d-9dde-4208-b401-32e4654b69a3" | |
| }, | |
| "source": [ | |
| "X = df[['age', 'gender']]\n", | |
| "X = pd.get_dummies(X).values\n", | |
| "\n", | |
| "z = df['treatment']\n", | |
| "z = z.values\n", | |
| "\n", | |
| "y = df['outcome']\n", | |
| "y = y.values\n", | |
| "\n", | |
| "X.shape" | |
| ], | |
| "execution_count": 14, | |
| "outputs": [ | |
| { | |
| "output_type": "execute_result", | |
| "data": { | |
| "text/plain": [ | |
| "(100000, 9)" | |
| ] | |
| }, | |
| "metadata": { | |
| "tags": [] | |
| }, | |
| "execution_count": 14 | |
| } | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "metadata": { | |
| "id": "NwjuP7SRc3Ns", | |
| "colab_type": "text" | |
| }, | |
| "source": [ | |
| "## Linear model for outcome" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "KlqesOyqb5Or", | |
| "colab_type": "code", | |
| "colab": { | |
| "base_uri": "https://localhost:8080/", | |
| "height": 34 | |
| }, | |
| "outputId": "1d7f86e0-124e-4f6e-a7d4-71776dbbdd9c" | |
| }, | |
| "source": [ | |
| "model_linear = LinearRegression().fit(X, y)\n", | |
| "\n", | |
| "print('score: ', model_linear.score(X, y))" | |
| ], | |
| "execution_count": 28, | |
| "outputs": [ | |
| { | |
| "output_type": "stream", | |
| "text": [ | |
| "score: 0.4145992303231701\n" | |
| ], | |
| "name": "stdout" | |
| } | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "metadata": { | |
| "id": "yLwJFE8bc6ND", | |
| "colab_type": "text" | |
| }, | |
| "source": [ | |
| "## IPTW model for z" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "KQkKxKvsb5Sg", | |
| "colab_type": "code", | |
| "colab": { | |
| "base_uri": "https://localhost:8080/", | |
| "height": 34 | |
| }, | |
| "outputId": "92019333-392f-47f5-d289-ed9be3cbab9e" | |
| }, | |
| "source": [ | |
| "model_iptw = LogisticRegression().fit(X, z)\n", | |
| "\n", | |
| "print('accuracy: ', sum(model_iptw.predict(X)==z)/len(z))" | |
| ], | |
| "execution_count": 30, | |
| "outputs": [ | |
| { | |
| "output_type": "stream", | |
| "text": [ | |
| "accuracy: 0.90047\n" | |
| ], | |
| "name": "stdout" | |
| } | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "d2shzhtcd2Qs", | |
| "colab_type": "code", | |
| "colab": { | |
| "base_uri": "https://localhost:8080/", | |
| "height": 362 | |
| }, | |
| "outputId": "ab11b484-eb55-429f-df10-bf4897066bb1" | |
| }, | |
| "source": [ | |
| "df_ps_m = df\n", | |
| "df_ps_m['ps'] = model_iptw.predict_proba(X)[:, np.where(model_iptw.classes_ == 1)].flatten()\n", | |
| "df_ps_m['m'] = model_linear.predict(X).flatten()\n", | |
| "df_ps_m.head(10)" | |
| ], | |
| "execution_count": 39, | |
| "outputs": [ | |
| { | |
| "output_type": "execute_result", | |
| "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>user_id</th>\n", | |
| " <th>age</th>\n", | |
| " <th>gender</th>\n", | |
| " <th>treatment</th>\n", | |
| " <th>outcome</th>\n", | |
| " <th>ps</th>\n", | |
| " <th>m</th>\n", | |
| " </tr>\n", | |
| " </thead>\n", | |
| " <tbody>\n", | |
| " <tr>\n", | |
| " <th>0</th>\n", | |
| " <td>1</td>\n", | |
| " <td>10y</td>\n", | |
| " <td>male</td>\n", | |
| " <td>1</td>\n", | |
| " <td>0.254905</td>\n", | |
| " <td>0.171632</td>\n", | |
| " <td>0.297852</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>1</th>\n", | |
| " <td>2</td>\n", | |
| " <td>10y</td>\n", | |
| " <td>female</td>\n", | |
| " <td>0</td>\n", | |
| " <td>0.391095</td>\n", | |
| " <td>0.029579</td>\n", | |
| " <td>0.146484</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>2</th>\n", | |
| " <td>3</td>\n", | |
| " <td>40y</td>\n", | |
| " <td>female</td>\n", | |
| " <td>0</td>\n", | |
| " <td>0.271188</td>\n", | |
| " <td>0.029094</td>\n", | |
| " <td>0.145020</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>3</th>\n", | |
| " <td>4</td>\n", | |
| " <td>10y</td>\n", | |
| " <td>male</td>\n", | |
| " <td>0</td>\n", | |
| " <td>0.175226</td>\n", | |
| " <td>0.171632</td>\n", | |
| " <td>0.297852</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>4</th>\n", | |
| " <td>5</td>\n", | |
| " <td>30y</td>\n", | |
| " <td>male</td>\n", | |
| " <td>1</td>\n", | |
| " <td>0.502966</td>\n", | |
| " <td>0.726206</td>\n", | |
| " <td>0.576660</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>5</th>\n", | |
| " <td>6</td>\n", | |
| " <td>70y</td>\n", | |
| " <td>male</td>\n", | |
| " <td>0</td>\n", | |
| " <td>0.196705</td>\n", | |
| " <td>0.167537</td>\n", | |
| " <td>0.295410</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>6</th>\n", | |
| " <td>7</td>\n", | |
| " <td>50y</td>\n", | |
| " <td>male</td>\n", | |
| " <td>0</td>\n", | |
| " <td>0.177865</td>\n", | |
| " <td>0.167257</td>\n", | |
| " <td>0.297363</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>7</th>\n", | |
| " <td>8</td>\n", | |
| " <td>30y</td>\n", | |
| " <td>female</td>\n", | |
| " <td>1</td>\n", | |
| " <td>0.447726</td>\n", | |
| " <td>0.280675</td>\n", | |
| " <td>0.425293</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>8</th>\n", | |
| " <td>9</td>\n", | |
| " <td>60y</td>\n", | |
| " <td>female</td>\n", | |
| " <td>0</td>\n", | |
| " <td>0.130806</td>\n", | |
| " <td>0.030219</td>\n", | |
| " <td>0.145508</td>\n", | |
| " </tr>\n", | |
| " <tr>\n", | |
| " <th>9</th>\n", | |
| " <td>10</td>\n", | |
| " <td>30y</td>\n", | |
| " <td>female</td>\n", | |
| " <td>1</td>\n", | |
| " <td>0.479201</td>\n", | |
| " <td>0.280675</td>\n", | |
| " <td>0.425293</td>\n", | |
| " </tr>\n", | |
| " </tbody>\n", | |
| "</table>\n", | |
| "</div>" | |
| ], | |
| "text/plain": [ | |
| " user_id age gender treatment outcome ps m\n", | |
| "0 1 10y male 1 0.254905 0.171632 0.297852\n", | |
| "1 2 10y female 0 0.391095 0.029579 0.146484\n", | |
| "2 3 40y female 0 0.271188 0.029094 0.145020\n", | |
| "3 4 10y male 0 0.175226 0.171632 0.297852\n", | |
| "4 5 30y male 1 0.502966 0.726206 0.576660\n", | |
| "5 6 70y male 0 0.196705 0.167537 0.295410\n", | |
| "6 7 50y male 0 0.177865 0.167257 0.297363\n", | |
| "7 8 30y female 1 0.447726 0.280675 0.425293\n", | |
| "8 9 60y female 0 0.130806 0.030219 0.145508\n", | |
| "9 10 30y female 1 0.479201 0.280675 0.425293" | |
| ] | |
| }, | |
| "metadata": { | |
| "tags": [] | |
| }, | |
| "execution_count": 39 | |
| } | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "rjD0qFsxfF64", | |
| "colab_type": "code", | |
| "colab": { | |
| "base_uri": "https://localhost:8080/", | |
| "height": 574 | |
| }, | |
| "outputId": "7a4bbc2d-e041-403b-a714-e4967bf54ecd" | |
| }, | |
| "source": [ | |
| "t_is_1 = df_ps_m[df_ps_m['treatment'] == 1].sort_values(by=['ps'])[['user_id', 'ps', 'outcome']].round(4)\n", | |
| "t_is_0 = df_ps_m[df_ps_m['treatment'] == 0].sort_values(by=['ps'])[['user_id', 'ps', 'outcome']].round(4)\n", | |
| "\n", | |
| "bins = 25\n", | |
| "\n", | |
| "sns.distplot(t_is_1.ps, hist=True, kde=True, bins=bins, kde_kws={'linewidth': 4}, label='1')\n", | |
| "sns.distplot(t_is_0.ps, hist=True, kde=True, bins=bins, kde_kws={'linewidth': 4}, label='0')\n", | |
| "plt.legend(frameon=False, loc=0, ncol=1, prop={'size': 20})" | |
| ], | |
| "execution_count": 40, | |
| "outputs": [ | |
| { | |
| "output_type": "execute_result", | |
| "data": { | |
| "text/plain": [ | |
| "<matplotlib.legend.Legend at 0x7fe0095d2898>" | |
| ] | |
| }, | |
| "metadata": { | |
| "tags": [] | |
| }, | |
| "execution_count": 40 | |
| }, | |
| { | |
| "output_type": "display_data", | |
| "data": { | |
| "image/png": "iVBORw0KGgoAAAANSUhEUgAABBIAAAIcCAYAAABLgo2MAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjAsIGh0\ndHA6Ly9tYXRwbG90bGliLm9yZy8GearUAAAgAElEQVR4nOzdeZRkd3kf/OdWVa8zo9k0Go00WkAL\nCIQAIYyNLUwkjPBrMDZewMJ2gokdk9h+Ewef88YvB5wYcoJNfHLIIcE5JnkdgwGT5IWXxYjFMlZY\nZLEjJCQNQttoJI1mn+mluqru+8dopvve2z1T3V3rrc/nnJz43q6q+5O4M+r7ref3PEmapmkAAAAA\ntKHS7wUAAAAAw0OQAAAAALRNkAAAAAC0TZAAAAAAtE2QAAAAALRNkAAAAAC0TZAAAAAAtK3W7wUc\nOnQiWq2038tgSGzfvjEOHDje72VAz7n3GVXufUaZ+59R5d7vv0olia1bN6z4874HCa1WKkhgVdwv\njCr3PqPKvc8oc/8zqtz7g83WBgAAAKBtggQAAACgbYIEAAAAoG2CBAAAAKBtggQAAACgbYIEAAAA\noG2CBAAAAKBtggQAAACgbYIEAAAAoG2CBAAAAKBtggQAAACgbYIEAAAAoG2CBAAAAKBtggQAAACg\nbYIEAAAAoG2CBAAAAKBttX4vAAAAANai0YqYX2j0exkrmhirRa0DX9/feuvn4pvf/Hrcd9+9sWfP\nfTEzcyJe/vKfjLe+9Q/X/+FrIEgAAABgKM0vNOKOux/v9zJW9MKrdkZtYv2P3X/+5/819uy5N6am\npuO8886LBx880YHVrZ0gAQAAAAbY7/zO78aOHefF7t0XxTe+8bX4nd/5zb6uR5AAAAAAA+zaa6/r\n9xIyBAkAAABP6dWe+07tnYd+ECQAAAA8pVd77ju1dx76QQYGAAAAtE2QAAAAALRNkAAAAAC0TZAA\nAAAAtE2QAAAAALRNkAAAAAC0zbwRAAAAGGB/93d/G7fd9rcREXHgwIGIiLjzzm/HO97xBxERsXnz\nlvit3/rnPVuPIAEAAAAG2H333RN//defyJx79NG98eijeyMi4vzzdwkSAAAA4Gwmxmrxwqt29nsZ\nK5oY68wj9xvf+E/ijW/8Jx35rE4QJAAAADCUapWI2oTH2l7TbBEAAABomyABAAAAaJsgAQAAAGib\nIAEAAABomyABAAAAaJsgAQAAAGibIAEAAABomyABAAAAaJsgAQAAAGibIAEAAABomyABAAAAaJsg\nAQAAAGibIAEAAABoW63fCwAAAIC1GI/5iIW5fi9jZWOTUY+Jjn3cE088Hn/2Z++N22//chw9eiS2\nbz83rr/+pfGGN/x6nHPOOR27ztkIEgAAABhOC3NxfM83+r2KFW28/PkRY50JEvbufSR+8zd/LQ4d\nOhjXX//jcfHFl8bdd383PvKRD8btt38p/vN/fl9s3rylI9c6G0ECAAAADLh//+//XRw6dDD++T9/\nc/z8z7/u9Pn/+B//JD784b+M//Jf/lP83u/9fk/WctYeCYcOHYpf//Vfj5tuuile9apXxW/91m/F\nwYMHIyLim9/8Zvz0T/903HTTTfFrv/ZrceDAga4vGAAAAEbJ3r2PxN///Vdi164L4jWv+cXMz974\nxn8SU1NTccstn4rZ2dmerOesQUKSJPGP//E/jltuuSU+/vGPx0UXXRTvete7otVqxe/93u/FW9/6\n1rjlllviuuuui3e96129WDMAAACMjK9//asREfHCF74oKpXsY/z09IZ4znOeG3Nzc/Hd736nJ+s5\na5CwZcuWeNGLXnT6+HnPe148+uijceedd8bExERcd911ERHxute9Lj796U93b6UAAAAwgh566MGI\niLjookuW/fnu3RdFRMTDDz/Uk/Wsavxjq9WKD37wg3HDDTfEvn374oILLjj9s23btkWr1YrDhw93\nfJEAAAAwqo4fPx4RERs3blz25xs2bHzqdcd6sp5VNVv8wz/8w5ieno5f/uVfjs9+9rMdWcD27cv/\ni4CV7Nixqd9LgL5w7zOq3PuMMvd/76UHZ2LTxsmuX2d6eiJ2bJvu+nWGVbv3/sLh2Ug2dm68YqdN\nTY/H5i3r/3M8NTUWERGbNk0u++9mw4aJ0/9/L/7eaDtIeOc73xkPPvhgvPe9741KpRK7du2KRx99\n9PTPDx48GJVKJbZsWd24iQMHjkerla7qPYyuHTs2xf79vUnZYJC49xlV7n1Gmfu/P2bmG3Hs+Fz3\nrzMzH/ubza5fZxit5t4fX6jH8ePzXV7R2qUz9agvrP/PcbV6MijYt+/JZf/dPPHEycEHSTLekb83\nKpXkjF/6t7W14U/+5E/izjvvjPe85z0xPj4eERFXX311zM3NxVe/erLpw4c+9KF4xStese4FAwAA\nAIsuvvhkb4SHH35w2Z8/8sjDERFx0UUX92Q9Z61IuO++++JP//RP49JLL43Xve7krMrdu3fHe97z\nnvijP/qjeNvb3hbz8/Nx4YUXxh//8R93fcEAAAAwSq699uSQgzvuuD1arVZmcsPMzIn4zne+FZOT\nk/HsZz+nJ+s5a5BwxRVXxD333LPsz6699tr4+Mc/3vFFAQAAACddeOHu+KEf+uH4+7//Svyv//VX\n8fM//7rTP3vf+/40Zmdn49Wvfk1MTU31ZD2rarYIAAAA9N6//Jf/V/zmb/5a/If/8K742tfuiEsu\neVrcdded8fWvfzUuuuji+I3f+Kc9W4sgAQAAAAbchRfujj/7s/8e73vfn8btt38pvvzlL8b27efG\nL/zCL8Ub3vDrcc455/RsLYIEAAAAhtPYZGy8/Pn9XsXKxjo7SnTnzvPj93//bR39zLUQJAAAADCU\n6jERMTbR72WMnLbGPwIAAABECBIAAACAVRAkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAA\nAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAA\nAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAA\nbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABt\nEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0T\nJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMk\nAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQA\nAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAA\nAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAA\nAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAA\nbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABt\nEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG0TJAAAAABtEyQAAAAAbRMkAAAAAG2r\ntfOid77znXHLLbfE3r174+Mf/3hceeWVERFxww03xPj4eExMTERExJvf/Oa4/vrru7daAAAAoK/a\nChJuvPHG+NVf/dV4/etfX/jZu9/97tPBAgAAAFBubQUJ1113XbfXAQAAAAyBtoKEM3nzm98caZrG\nC17wgvjd3/3dOOecczqxLgAAAGAAJWmapu2++IYbboj3vve9p7cy7Nu3L3bt2hX1ej3e8Y53xIkT\nJ+Jd73pX1xYLAADQTU8cnImv3/NE169z7TPOi/O2TXf9OtAN66pI2LVrV0REjI+Px8033xxvetOb\nVv0ZBw4cj1ar7SyDEbdjx6bYv/9Yv5cBPefeZ1S59xll7v/+mJlvxLHjc92/zsx87G82u36dYeTe\n779KJYnt2zeu/PO1fvDMzEwcO3byf9w0TeNTn/pUXHXVVWv9OAAAAGAItFWR8Pa3vz0+85nPxJNP\nPhlveMMbYsuWLfHe9743fvu3fzuazWa0Wq247LLL4m1ve1u31wsAAAD00ap6JHSDrQ2sRi/LnMZj\nPmKh+2VtERExNhn1mOjNtRhKSvwYVe59Rpn7vz9OzDfijrsf7/p1XnjVztgwse7e96Xk3u+/s21t\ncOfCShbm4vieb/TkUhsvf37EmCABAAAYfGvukQAAAACMHkECAAAA0DZBAgAAANA2QQIAAADQNkEC\nAAAA0DZBAgAAANA2QQIAAADQNkECAAAA0DZBAgAAANA2QQIAAADQNkECAAAA0DZBAgAAANA2QQIA\nAADQNkECAAAA0DZBAgAAANA2QQIAAADQNkECAAAA0LZavxcAg6h54OE4/u1PRGNuLqoXXR1Jdazf\nSwIAABgIggTISRvzMfvJP4p07tjJE8161C57UX8XBQAAMCBsbYCc5hP3L4YIEdE8+HAfVwMAADBY\nBAmQ0zq8L3ti7nikzYX+LAYAAGDACBIgpxAkREQ6e7QPKwEAABg8ggTIaR15rHAunTnch5UAAAAM\nHkEC5CxbkSBIAAAAiAhBAmSkjXqkxw4UzwsSAAAAIkKQABmto49HRFo4n84c6f1iAAAABpAgAZZY\nbltDREQ6eyzSVqPHqwEAABg8ggRYonW42GjxpNTkBgAAgBAkQMZKFQkREekJfRIAAAAECbDEcqMf\nT9EnAQAAQJAAp6VpeuaKBJMbAAAABAlwSjp7JGJhbuWfCxIAAAAECXBKvhqhsunczHE6dzTSVrOX\nSwIAABg4tX4vAAZFPkio7bgk6vMzEfWZkyfSNNLZY5Fs2NKH1dFvrZnDUf/mJyPSNMaveUUhaAIA\ngFGhIgGekh/9WN18XiTT2dAgnbW9YVTN3/bnsXDnZ2Phu5+L2c+8O9K01e8lAQBAXwgS4CmtI9mK\nhOo550VlenPmnMkNoyltNaLx0DdPH7cOPBTNfff0cUUAANA/ggR4SqEiYcvOYkWChosjKT1xKCJN\nM+ca9325T6sBAID+EiRARKSNeqTHnsycq56zQ5BARES0jh0onFu4/45IG/U+rAYAAPpLkAAR0Tr6\neEQsfuOcbNweSW08kmW2NtgbP3rS48UgIRZmo/HQt3q/GAAA6DNBAkRxW0Nly66IiEjGJiPGJhd/\nkLYinTvey6UxAFq5apVTGntsbwAAYPQIEiCKox8rm88//X/b3kB6fIUg4aFvCZYAABg5ggSIZYKE\nLYIEFrWW29oQEdFqxsL9d/R2MQAA0GeCBIiI1pH81oYLFv/vQp8EQcKoWa7Z4im2NwAAMGoECYy8\nNE2LPRLOuLXhSE/WxWBI09byzRaf0nzs3mgd29/DFQEAQH8JEhh56eyRiIXZxRO1iUg2bD19uNzW\nhjRNg9GQzh6NaDUWT4xPR2Xn5ZnXLOz5So9XBQAA/SNIYOQt1x8hSZLFE2OTEbWJJW9oRsxrsDcq\n0tzEhsqmc2Ps8h/JnGvc9yXhEgAAI0OQwMgrbmvYlTlOkiSSXJ+Elj4JIyM/+rGycXvULvuhiKS6\n+JrD+6J14MFeLw0AAPpCkMDIO9PEhlP0SRhd+YkNyaZzozK5KaoXPSdzfuE+TRcBABgNggRGXnFi\nw67Ca4yAHF35RouVjdsjImLsihdnzjf2fCXSVrNn6wIAgH4RJDDyChUJm4sVCRVBwsjKb21IngoS\napc872T/jKeks0ei+ejdPV0bAAD0gyCBkZY26sVmessECfkeCSY3jI70eLHZYkREUhuP2tNemPmZ\n7Q0AAIwCQQIjrXX0iYhYDASSjdsjGZsovnB8OqI6tnjcbETUZ7q/QPoqTdNoHSv2SDhl7Irs9Ibm\nE9/vyboAAKCfBAmMtHa2NUScmtxge8PImT8R0ZhfPK6NRzKx8fRhZdvuzMvT2aO9WhkAAPSNIIGR\n1s7EhlOMgBw9rfy2ho3nRpIkp49PhgqLx1GfibTV6NHqAACgPwQJjLR2JjacUqhI8O1z6RW3NWzP\nHlcqkUxuzJxL5453fV0AANBPggRGWutwLkjYfIYgYXw6e2JhfvkXUhqFRpwbtxdek0xtyr5HwAQA\nQMkJEhhprSP5rQ1nCBJq45njtCFIKLv81oaljRZPn5s8J3Oczh7r6poAAKDfBAmMrLRRj6jPLp5I\nKpFs2LryG2q5aQ6NencWxsBIj2e3NlQ2LhMk5CsS5gQJAACUmyCBkZXmxjcmExsyjfQKcmMhVSSU\nX75HwrJbGyZtbQAAYLQIEhhZ6dyJ7ImJDWd8fX5rg4qE8mtra8NUbmuDigQAAEpOkMDISuvZICE5\nS5AQywQJadrq8KoYFGl9NmJ+yT1SqRZGgEaoSAAAYPQIEhhd86sLEpKkElEdy55sLHR6VQyIVq4/\nQrJh28l7IEePBAAARo0ggZGVrjJIiIhlGi7qk1BWaW5bQ2WZbQ0RpjYAADB6BAmMrGKQMH3W9xRH\nQOqTUFb5RovJMhMbIooVCa05WxsAACg3QQIja00VCWMqEkZFYfTjpuLEhohlmi2qSAAAoOQECYys\ntQQJhYqEBUFCWbWO5bY2LDP6MeKp+2bp2ND6TKTNRjeXBgAAfSVIYGR1pkeCrQ1l1c7ox4iTTTiT\niY2ZcxouAgBQZoIERlY6P5M5XlNFgq0NpZXmeiSsVJEQscz2BkECAAAlJkhgZOUrEkJFAk9JG/VI\nZ48snkiSSDZsW/H1yWRuBOSshosAAJSXIIGRtbYeCdkgQUVCOaUnDmaOk+mtkVRrK74+P7lBRQIA\nAGUmSGB0ralHQnZrg4qEcsqPfjzTtoaIiGTS5AYAAEaHIIGRlKZpsSJhfPqs70ty4x9NbSinYqPF\nswQJ+YoEWxsAACgxQQKjaWEuIm0tHlfHC40Ul6UiYSSkhdGPy09sOKXQI8HWBgAASkyQwEgq9kc4\nezVChB4Jo6J1PLu1ITnb1gZTGwAAGCGCBEZSMUjY2N4bVSSMhDQXJFQ2ra4ioWVrAwAAJSZIYCSl\n9ZnMcbsVCVEdi4hk8bjViLTV7NzCGAitY6vtkaAiAQCA0SFIYCSlc8czx21NbIiIJEmWqUqwvaFM\n0lYz0hOHMufOOrVBs0UAAEaIIIGRlK9IiDaDhIjlJjfY3lAm6YlDmUacyeSmQm+MvGRiQ0SypFKl\nPhtps9GtJQIAQF8JEhhJ6Vy+R0L7QYKKhHIrNFo8S3+EiIgkqZjcAADAyBAkMJrqaw8SipMbVCSU\nSXH045m3NZxSCBJsbwAAoKQECYykYkVCm80WIyLyZe4qEkplLRUJEcsECSoSAAAoKUECIyktVCS0\nOf4xorC1QUVCuay5IiE/uUFFAgAAJSVIYCSl82sc/xjFrQ0qEsolX5FQ2agiAQAAlhIkMJLS+fz4\nx1VUJIzlmy2qSCiTwmjQ6c1tva9YkSBIAACgnAQJjKROViSkCyoSyqQQMk22FzIlU/mKBFsbAAAo\nJ0ECIyn/sBjrGv+oIqFM0vm1TfQoTm1QkQAAQDkJEhg5aasVUZ/NnEvG11GRoEdCaaTNRsTC3OKJ\nJIkYn2rrvfmtDS09EgAAKClBAqOnnt3WEGNTkVSq7b+/UJEgSCiLQjXC+IZIkvb+mlSRAADAqBAk\nMHIKD4uTq9jWEBHJWL4iwdaGssjfG7GKe0OPBAAARoUggZGz3LfOq7JMRUKaputcFYNgrf0RTr82\nSRZP1GcjbS50amkAADAwBAmMnHVXJFRqEUu3QqRpRLPRiaXRb+sYC5okleL2htwoSQAAKANBAiOn\nWJHQfqPF03INF/VJKIf1VCRELNcnwfYGAADKR5DAyFnvw2JERJLb3mByQzmkc+sMEnKTG1KTGwAA\nKCFBAiOnE0FCsSJBw8UySAtbG1QkAABAniCBkVPozK8igacU+2e03yMhYrnJDSoSAAAoH0ECIyed\nn8kcr6kiYUxFQhmtv0dCbmvDrCABAIDyESQweurd6JEgSCiDYpCw3ooEWxsAACgfQQIjpys9EhZs\nbSiD/LjGVY8GLfRIUJEAAED5nDVIeOc73xk33HBDPOMZz4h777339Pkf/OAH8drXvjZuuummeO1r\nXxsPPPBAN9cJHWNqAytZ99aG3NSGlh4JAACU0FmDhBtvvDE+8IEPxIUXXpg5/7a3vS1uvvnmuOWW\nW+Lmm2+Ot771rV1bJHRS8WFxevUfYmpDKRWnNqxya4OKBAAARsBZg4Trrrsudu3alTl34MCBuOuu\nu+KVr3xlRES88pWvjLvuuisOHjzYnVVCB613H3yEioQySlvNiPps9uT46kKmylS+2aIeCQAAlM+a\neiTs27cvdu7cGdVqNSIiqtVqnHfeebFv376OLg46LW0uZKsHkkrE2OTqP8jUhtJJ69lpHjE+HUll\nlX9FTkyfvKdOWZg9ec8BAECJ1Pq9gO3bV/9tMKNtx45NZ3/RChrHD8XS4vXK1MY477xzln3twuHZ\nSDZOLPuzZrIpDi05Tpr12LTCa9sxNT0em7es/Z+L9asfOBZLa1Vq05vWdK/NTp8TzROHTx9vm06j\ndk5n/rddz70Pw8y9zyhz//deenAmNm1cwxdNqzQ9PRE7tq1hi+2IcO8PtjUFCbt27YrHH388ms1m\nVKvVaDab8cQTTxS2QLTjwIHj0Wqla1kGI2jHjk2xf//a9503Dz2ROU7Hplb8vPGFehw/vvyWhTT3\nJXOrPhfHVnhtO9KZetQX7Kfvp+bjj2eOW2PTa7rX0vGNEUuChCf37ovq/PgZ3tGe9d77MKzc+4wy\n939/zMw34tjxue5fZ2Y+9jebXb/OMHLv91+lkpzxS/81bW3Yvn17XHXVVfGJT3wiIiI+8YlPxFVX\nXRXbtm1b2yqhRzoy+jEiItcjIZoLkaatNa6KQdCpeyOZyjVcNLkBAICSOWuQ8Pa3vz1e8pKXxGOP\nPRZveMMb4qd+6qciIuIP/uAP4v3vf3/cdNNN8f73vz/+9b/+111fLKxbpx4Wk0pEdSx7Up+Eodax\nIKEwuUHDRQAAyuWsWxve8pa3xFve8pbC+csuuyw+8pGPdGVR0C2Fh8XxNVYkRJwcAbm0kV6jvrbG\njQyETkzziFimIsEISAAASmZNWxtgWHVsa0NEJLnJDUZADrd07njmOJlca0VCbgSkrQ0AAJSMIIGR\nUgwS1tEpN98nwdaGoZbO54KETvVIsLUBAICSESQwUjpVvh4RkeSChHRBRcIw69jWhnyPBBUJAACU\njCCBkdLZioTs1gYVCcOtc1MbslsbWioSAAAoGUECI6WrFQl6JAy1dM74RwAAaIcggZGS1meyJ1Qk\n8JR8yBRrbLZYyTdbNLUBAICSESQwWvKd+dcztaFmakOZFJstrrFaZWI6IlnyV+vCbKRLx4QCAMCQ\nEyQwUvIVCesJEmIsP7VBkDCs0rQVMZ+/N9ZWrZIklUgmsyGEqgQAAMpEkMDISNO0Y/vgI5apSFiw\ntWFo1WcjIl08HpuMpFJb88cl+e0N+iQAAFAiggRGR2M+Im0uHlfHCg0TVyX/XhUJQ6vQhHNy7U04\nI5ZpuGhyAwAAJSJIYGR0arzf6fcXeiSoSBhWaQd7Z0REJJMmNwAAUF6CBEZGp4OE4tQGFQnDqmON\nFk+9X0UCAAAlJkhgZHQ8SKjWIiJZPG41I2011veZ9EXHq1Wm9EgAAKC8BAmMjPzDYoyvrSv/KUmS\nRIzlqxJsbxhGnWzCGbHM1gYVCQAAlIgggZFRbKi3zoqEiEKzRpMbhlPHKxLyQUI+xAIAgCEmSGB0\n5B8Wx9cfJJjcUA6FHgnrDJnyQYQgAQCAMhEkMDLS+ZnMcWcqEkxuKINiRcI6my0KEgAAKDFBAiOj\n8LDYkYoEkxvKoPCgv+6tDdkgQpAAAECZCBIYGd3okZDf2pAKEoZSx3sk5CsS5gQJAACUhyCBkdGN\nioT81gZTG4bUXL5Hwvq2NkRtIiKpLh4367a9AABQGoIERkbxW+f1jX+MiIix/NQGFQnDqOMVCUlS\nqHixvQEAgLIQJDAyOt1QL0JFQhmkadrxIGG5zxAkAABQFrV+LwB6pdhQrwMVCZotDr+F2Yi0tXhc\nm4ikOrbiy8djPmJh7qwfOzs+GUs+NWonHo+xTZvaX9fYZNRj4uyvAwCAHhMkMBLSViuiPps515Fv\nnQvNFlUkDJt8I8Sz3hcLc3F8zzfO+rmtxkLmeOaBO6N6/Fjb69p4+fMjxgQJAAAMHlsbGA31mYhI\nF4/HJiOpVFd8edtUJAy9rkzziFjm3hAyAQBQDoIERkJan8kcd6IaISIiyTdb9LA4dNL53MSGDvTO\niCj2zzAaFACAshAkMBK60UwvIiJyWxtiYT7SNF3+tQykrt0b+W0JQiYAAEpCkMBI6NbDYlKpRWS2\nSKQRzUZHPpve6Nq9ke+fYTQoAAAlIUhgJHTtW+cIfRKGXDqX29ow2ZmtDe4LAADKSpDASCgGCR0Y\n/XjqswqTGzwwDpOeVSS4LwAAKAlBAiOh+LDYoW+dI3TnH3L5eyP0SAAAgDMSJDASig+LHaxIKExu\n8M3zUOlSyFSY2qBHAgAAJSFIYDT0siJhwTfPw6TQI6FjEz30SAAAoJwECYyEdH4mc9zZHgm5b549\nMA6VwraXyU4FCWPZ4+ZCpGmrM58NAAB9JEhgJKTzuW+dxzs5tSG7tcFe+OFSuDc6tbUhqbg3AAAo\nJUECI6FQkdCpb51DRcIwS9O0p6NB9UkAAKAMBAmMhMLDoooEIk72LWg1F4+rY4WxjeuRD5n0SQAA\noAwECYyEru2Dj4hkTEXCsOpqNUJERGGih5AJAIDhJ0ig9NJmI/tNcJJEjE127gIqEoZWYWLDZAen\neYSKBAAAykmQQOktt60hSTp36xd6JNgHPzS6XpHg3gAAoIQECZReWs8+LEbHHxbzFQkeFodFMUjo\ndEWCewMAgPIRJFB+c/mHxZTqUFAAACAASURBVOnOfn7+YbG5EGna6uw16Iru90jIb22w7QUAgOEn\nSKD08hUJnX5YTJJKRHUse9ID41BI57M9EjpdrWLbCwAAZSRIoPTSQkVCh791jijuhRckDIVis8Xu\n9kiwtQEAgDIQJFB6aX0mc9yNICE/AtID45DocY8EARMAAGUgSKD0Ct86d6UiIdcnQQn7UOh9jwT3\nBQAAw0+QQOn1pCLB1oah1O0gQY8EAADKSJBA6XX9W+cIIyCHVNf7ZxTui3qkadrZawAAQI8JEii9\nfJAQnR7/GCoShlV+akMy2eEeCZVqRKW29IoRzYWOXgMAAHpNkEDpFSsSOvuwGBHLfvPM4OvJvaFP\nAgAAJSNIoPwKD4tdqEgYy1ckeFgcdGmjnq0OqNSKgVAHFCY36JMAAMCQEyRQev2pSPCwOOiW652R\nJEnnL1TLVySoVgEAYLgJEii1NE2XeWDsQY+EBQ+Lg67YH6ELTThjuf4ZQiYAAIabIIFya9QjWs3F\n42ototr58nUVCcOnOLGhC5UqEcUeCbY2AAAw5AQJlFqhGmG8O+XrpjYMn0JFQjfGgsYyPRLcGwAA\nDDlBAqW23D74rijsg/et86ArjgV1bwAAQDsECZRaz4KEai0illQ6tJqRthrduRYdUdza0KWKBBM9\nAAAoGUECpVb81rnzjRYj4uR2ifxeeCXsg61HWxv0zwAAoGwECZRbL0Y/nvrs/F54TfUGWjqXn9rQ\nnXvDRA8AAMpGkECp9WL042mFvfAeGAdZcfzjpu5cSI8EAABKRpBAqRWDhB5WJHhgHGg9q0jQIwEA\ngJIRJFBqaX0mc6wigVPSuWOZ4+5VJOR6JCzMR5qm3bkWAAD0gCCBUit869ythnoRhQdG3zwPtl5V\nJESlFpEs+as2bUW0mt25FgAA9IAggVIrViR0L0jIl7CHpnoDK01bxR4JXdr2kiSJPgkAAJSKIIFS\nK/ZIUJFARMzPRCzdXjA2FUm11rXLFfskCJkAABheggRKrZdBQn7Mnx4Jg6s4saF7TTgjYtk+CQAA\nMKwECZRaPkiIrlYkKF8fFj3rj3Dq82smNwAAUB6CBEorTVsnS9iX6ObUBuMfh0fPJjackq9IcG8A\nADDEBAmUV302Ipbug5+MpNK9ffDGPw6Pnlck6JEAAECJCBIorUJ/hPHuVSNERCRjuYoE++AHVnEs\naLd7JOQnerg3AAAYXoIESivNb2uY7GJ/hIhlytfrkS6dDMDAKG5t0CMBAADaJUigtAqd+ce7GyQk\nlVpEpbp0BRHNha5ek7Upbm3QIwEAANolSKC0ChUJ3ZzYcIo+CUOh1+MfCz0SFtwXAAAML0ECpVV4\nWOxBkGByw3DodbNFo0EBACgTQQKlla9ICBUJPKXX4x/1SAAAoEwECZRWYWpDLyoSTG4YCj2vSBgr\nNuIEAIBhJUigvApBQnfHP0aEioQhkKatZba9dDlIqOaChOZCpK1Wd68JAABdIkigtIoVCV1+WAwl\n7ENhfiZi6VjOsalIqrWuXjJJEn0SAAAoDUECpVUMEnpRkWDM36Dr+baGU9fRiBMAgJIQJFBaaX0Q\nKhJsbRg0vR79eNqYbS8AAJSDIIHSSudUJFDU64kNp+VDJo04AQAYUoIESiutZ8c/9mZqg4qEQde/\nrQ16JAAAUA7d7TAGfZK2GhELc0vOJBHjU6v6jFZE1BvN1V0390cqXZhv6zP07++d/lUk6JEAAEA5\nCBIopXQ+W40QE9ORJKsrwGk2W/HAvqOrek9toR4XLjluzM/F3jY+4zlXtCLGVnUp1qhQkdCDSpWI\nYrVKLKhWAQBgONnaQDkVJjb05mGxVcl+61xpelgcNMWtDX3qkaAiAQCAIaUigVIqjn7sVZCQLSuo\npI2ItBWxymoIVmc85nNbWVY2N3soczw2VonxhSNtX6saq9vuckp+/KOpDQAADCtBAqXUryAhkiRa\nSe1kgPCUSmshWtWJM7yJdVuYi+N7vtHWSxuHH88czz/5aCw00rYvtfnpz1rV0k5TkQAAQEn4mpRS\nKgQJ4z0Y/fiUZjW3vaG10LNr04b82MWxyZ5cttAjQZAAAMCQEiRQSoUgoUcj/iKW2d7QUsI+SNLc\nFojCA3635Lc25AMNAAAYEoIESik/taGXFQnFhosqEgZFmqbF3gS13gQJSWFrg4AJAIDhpEcCpZTO\n92fEX8RyFQmChIHRqEfEkn4I1bFIKtXeXDsfWDTmI03TSJKkN9cHAPqq0WzFXL0Z8/VmzC8046vf\neyI2TtbiGRdvjY1T5oAzXAQJlFKhIqGnQUK+R4JvngdFfltDr/ojREQklUpEdSxiaYVKox7Rq60V\nAEBP7d1/Iu564GAcOVGP+Xozmq3lmzuP1yrxomftjBuu3R2XnN+jsdSwToIESilfkRC9DBKq2US5\nqiJhcOQaHPasP8IptfFckDAvSACAkjl8fD6+9r39sffJE2d/cUTUG6247dv74rZv74vLd2+OG669\nMK57xnlRq9qFzuASJFBOfa1IsLVhUKW5Bof5vgXdltQmMo1A00Y9bGwAgHKYqzfjW3uejHsfPhxp\n+5OlM/Y8ciT2PHIkPrRhT/z4cy+IG1+wO87ZMH72N0KPCRIopcLUhn5ubWja2jAw+ri14eT1jIAE\ngLJpttK456FD8e09B6LeaK34ukqSxMR4NSbGKjE5Xovzt03HQ08ci4NHi78PHD1Rj49/6YG49Rt7\n47d/7jlxxe4t3fxHgFUTJFBK/QwSmioSBlba560NSW18aavHQoUE65M25iOdPRrJxu2RJMpBAei+\nxw/OxJfufCyOzSz/+955W6fi2it3xNZNE1GrJpkmyy+8amdMjlXiW3sOxOe/9kjc/eChwvuPzy7E\nH3/wm/Ebr3pWXPfM87r2zwGrJUigdNI07XNFgiBhYOUrEnq8tWG5yQ10RvPg3pj9xL+LdO5YVC96\nTkz9xG9HUlMKCkD37N1/Im79xt5oLdNEcePUWLzgGTvi4p0bzzihqVqpxLVX7ohrr9wRjz55Im79\n+t744p37Yq7ePP2aRrMV//mjd8Zrb7wiXv7Ci7ryzwKr5SsbyqdZj2g1Fo8r1ZNN7nrE1IbBVeiR\n0OOtDfnrpfW5FV7JatW/+j8jnTsWERHNh78T81/+YJ9XBECZPfrk8iHCWLUS1z5jR7z6+kvjkvM3\nrWrM8wXnbojXv/zK+Pf/7Efjp37kkszP0oj40Ofviw9+7r5orbUBA3SQigRKZ7nRj6v5S3y98lMb\nVCQMkEa+R0KPtzaMT2WO04XZnl6/rNLGfDQevjNzbuHuW6N6wTNj7LIX9WlVAJTVqcqBfIhwxe7N\n8bwrzo2pifU9Yk1N1OLnfvyyOG/rVPz3T9+TGRv52a8+HIeOzcWvv+pZMVarrus6sB4qEiid/OjH\nXm5riFhua4OKhEFRnNrQ62aL2SAh6oKETmg88t2TlUg5c3/336J15PE+rAiAstp34GSI0MyFCD/0\nrPPiR64+f90hwlLXX3NB/J+/cE1MjGcDg6/esz/++EPfjOOzvqyif9YdJNxwww3xile8Il796lfH\nq1/96rjttts6sS5Ys3xFQvQ4SEiTWqRLhvpV0lYkreYZ3kHP5JsbqkgohcYDX1/+BwtzMfu5/xRp\nQ5gHwPrtO3Ai/uZry4QIV50Xz7x4a1euefXTtse/ev21sXljduvsnkeOxL/9i6/FkeP6LdEfHalI\nePe73x0f+9jH4mMf+1hcf/31nfhIWLN+VyREkqhKGFBprtliz6c25IMEFQnrlraa0Xzwmyv+vHXg\nwZj/yod7uCIAyuixAzPLhggvfOZ58cxLuhMinHLxzk3xll+5Li44N/s77WMHZ+I9H70zGs2VR05C\nt9jaQPnkeySMT/d8CcWGi0rP+i1NWxH5b6YHYGtDqmHSujQf35MND8eno3bpCzKvWbjr87Fw/x09\nXhkAZfHYgZn4/NceKYQI1z1zR1x1aXdDhFO2b56Mf/XL18YzLtqSOb/nkSPx4c/v6ckaYKmOBAlv\nfvOb41WvelX8wR/8QRw9erQTHwlrVqhImNzY8zVouDiAGvU42fP4KdWxSCo9zlKrtYjKkr2TaWvZ\nvf20L7+toXbxc2PypW+MZNOOzPm5L/zXaB19opdLA6AEjhyfj7/5+vIhwrMu3dbTtWyYHIvffe3z\n4tm58OLzX38kvvidfT1dC6y7G8gHPvCB2LVrV9Tr9XjHO94R/+bf/Jt417ve1fb7t2/v/UMew23H\njk1n/PnBaiOW7hbbsHVrbDvLe5Zz+LFjMTkxdvYXLqc2EUsXMVltRXKGz6qN1WLLGtbISQuHZyPZ\neOZtCs0Ts7H0kb0yPhWbzvKe5YyNVdf0vlMOTkxHa3YxcJ2uNaO2zOdNTY/H5i3Ze+Js9/6oSdM0\nHn44u61h2zU/Ghsv3BnzP//m2Pvn//fiKNiF2Wh+8b/Fzl99ex9Wynq59xll7v/eSw/OxKaNk9Fo\ntuKTX34wGs1siPDia3bF8688b93XmZ6eiB3bVl85+/u/9sPxL/7DF+KJg4tVuH9xyz1x9ZXnxeW7\nt5zhncPFvT/Y1h0k7Nq1KyIixsfH4+abb443velNq3r/gQPHC6NTYCU7dmyK/fuPnfE1c4cOZY5n\nm2Nnfc9ykoVGzM2vrZJgIWqx9NGwMTcbc2Mrf1ZjobGmNXLS+EI9jp+l2VDraPbfb1obj2NraFBU\nWWiu6X2L182GBicOH41KFH+JSGfqUV9YXHM79/6oaR58OBqHl0xlqNZi5pzLYnb/sYixnTHxw6+N\n+S994PSP5x6+Ox7//vejcs76f/mjd9z7jDL3f3/MzDfi2PG5uP2ux+PAkWx/pWuvPDcuv+CcOHZ8\nboV3r+I6M/Oxv7m2htxv+ulnx799/9dioXGyP0K90Yq3v+/2eOs/ui42TY+f5d2Dz73ff5VKcsYv\n/ddV1zszMxPHjp38HzhN0/jUpz4VV1111Xo+EtYtnT+ROe5PjwRbGwZNodFir/sjnLquyQ0dk9/W\nUL3w2Zl/v2PPfllUznt65jXN/Q/0YmkADLmHHj8W9zx0OHPu0l2b4tlP6+12hpVccv6m+EeveGbm\n3IGjc/Gn/993fUlLT6wrSDhw4ED8yq/8SrzqVa+KV77ylfGDH/wg3va2t3VqbbAmhSBhssdTG2K5\nZov2wfddn0c/nmZyQ8c0HvhG5rh2yfMzx0mSRPX8KzPnmvt/0PV1ATDcDh6diy/d+Vjm3Mapsfjh\nZ++MJElWeFfv/cjV58eNL9idOXfXA4fif/7d9/u0IkbJurY2XHTRRfHRj360U2uBjihWJPQhSNBs\nceCkjWyQ0OvRj4vXzU1uUJGwJq3jB6L15ANLziRRu+R5hddVz700lv7py74HALKarVb8+V9/L+oL\niyMVkyTiJc/bFeO1ah9XtrzX3nB5PPT4sbjvkSOnz/31Vx6Kp51/Tlz3TFv56B7jHymdNDf+MSb6\nUZGQDRKqKhL6L7e1oeejH59S2NqgImFNGg9mqxEqOy+LynSxwVR1x6WZ4+aTDxq5CcCKPva/H4j7\nH81Oobv2yh1x7uapFd7RX7VqJf7pz1wdmzdmq2Hf96m7Y++TJ1Z4F6yfIIHSKYx/7EuQkN/aoCKh\n39KFwahIyG9tCEHCmhS3NVy77OuSc86LWFoFUp+J9Nj+bi4NgCF19wMH45NfeiBz7sJzN8SzcuMW\nB83mjRPxz372OVGtLG67mK83488+flc0W60zvBPWTpBAqaRpK6KerUjoT5CQ29rQFCT0XSNXkdC3\nrQ3ZSoh8E0jOLp0/Ec1Hv5c5N3bp85d9bZJUilUJGi4CkHP0RD3+y8fviqU1a1MT1Xjxc84fqL4I\nK7n8ws1x88uuyJx78PFj8Zk7Hu7Tiig7QQLlsjAXsbRsuTYRSXXdU05XrTi1wdaGfitUJNjaMLQa\nD387Il0cl1XZsisqW3at+PrKuZdkjvVJAGCpVprG+z55dxw5kf197ceu2RVTE73/PXKtXvr8Cwt9\nET562w/i8YMzK7wD1k6QQKmkc7lGi32oRoiIaNraMHjy3/z3a2tDvtlifdae/VXKj32sXbr8toZT\nqudemjluChIAWOIL33w0vnP/gcy55zx9W+za3p/fI9cqSZJ4/U9cGRsmF8OPhUYr/p+//l60/K5B\nhwkSKJW0PhhBwrJTG/wF3leD0iMhqdYiMvdHGpGbKMHK0uZCNB7+TuZcfuxj3nJbG4Q3AEREHD4+\nH//jb7PjEp9+wTnx3MvP7dOK1mfzhvF43Y3ZLQ73PHw4/u5bj/ZpRZSVIIFSKVYkTPdnIUk1Wsni\niKAk0kjSRn/WwsneGY3c9pI+bW2IsL1hPZp7785UlyRTm6Ny3tPP+B4NFwFYyYc+f1/Mzi/+jjYx\nXo1/+JPPjEpl8PsirOTFV58fVz9tW+bcR27dEweP6stE5wgSKJVBqUiI0HBxoDTqEUvbJ1XHIqn0\n8a+/ZbY30J7mY/dmjmuXPD+S5Mz/WyZJJaq5PgkaLgJw5/0H4u/vfiJz7mevf3psO6d/XzZ0QpIk\n8auveEZMjC1+qTU734z3f+ZeFXl0jCCBUhmUHgkRxRGQ1Zby9X7Jb2uIsf7+glCoSFgQJLSrefCR\nzHH1/CtWeGVWJbe9QcNFgNFWX2jGX3zmnsy5i3dujBtfcGGfVtRZ526eip/78WzF3jf3PBl3fO+J\nFd4BqyNIoFTyFQnRxyChWc3uwa82BQl90xiM/ginr29rw5q1Du3NHFe2tfcLn4aLACz18S89EPsP\nL9kqFxH/8BXPjGo/KxY77IZrd8flF27OnPvAZ++NYzOmibF+5fmTAhER89nxNv2sSGjm9uBXGx4W\n+yXNTWzo1+jH02xtWJO0PpvtbZAkUdlyQVvvLTRcfPJB5Z0AI2rv/uPx6dsfypy74QW742m7zunT\nirqjUkniH/3kM6NWXez3cGxmIT70+fv6uCrKQpBAqbRmj2SOk4mNfVpJRLOafVisNjW46ZvC1oZ+\nVyRkg4x80MHyCtUI5+yMpDa+wquzCg0X509ouAgwglppGv/9lnui2VoMk7dsHI/XvOTMjXuH1QXn\nbohX/ejTMue+/N3H49vfP7DCO6A9ggRKJT1xOHNc2bC1TyuJaFSzD4s1QULfpAO2tSFsbViTZmFb\nw+6236vhIgAREf/72/vivkeyXzzd/LIrY2qi1qcVdd9Pvuji2L0j++XaX9xyT9QXmn1aEWUgSKBU\n0hOHMsdJH4OEZi5IUJHQR/lv/Pu8tSGxtWFNWrlGi6sJEiI0XAQYdUdP1OMjt+7JnLvmsu3xgmfs\n6NOKeqNWrcQb/o9nRrJkouWBo3OF7R2wGoIESqWVDxKmt/RpJcv1SBAk9Et+akO/KxJMbVibQpCw\ndXWdtTVcBBhtH/6bPXFirnH6eLxWiV/+iSsjWfqEXVJP23VOvOwFF2XOffIrD8aTR/wOwtoIEiiN\ndGEuYukDWaUaydSmvq1HRcIAyYc4/d7akB8/uTAXadrqz1qGSD5IqG67aIVXLk/DRYDRdfeDh+LL\n330sc+7VP/a0OHfL1ArvKJ9X/9ilsWl67PTxQqMVf3Xr9/u4IoaZIIHSKGxrmN4SSdK/W7xZnYil\njyjVVj0itRetHwoVCf3e2lCpRtRyYYaGi2fUmjkS6dyxxRPVsZMNFFdBw0WA0dRsteIvP3tv5tzu\nHRviJ164ukB62E1PjsXP/fhlmXNf/d4TcfeDh1Z4B6xMkEBpFLY19LE/wskFVKJZzT4sVnNN/+iR\n/EN6vysSYpntDfoknFFhYsPWCyJZ5azvZRsu2t4AUHp/+41HY++TJzLnfvWmZ0atOnqPQj92za64\n9Pxsxe4HP3dvNFsqI1md0fvTQ2mlM7mJDX3sj3BKfntDrelhsR8GrUdCRGS/GY/QcPEsWgcfzhyv\nttHi6fflGy6a3ABQasdnF+Kjt92fOfejV58fl+/e3KcV9VclSeLmn7gyc+6R/Sfib7/xaJ9WxLAS\nJFAaA1eREBHNavZhUZ+E3kvTVkSjnj3Z560NERourlbrYLYiobp1bUGChosAo+Vj//sHmQaLE2PV\neE2uvH/UXH7h5njx1ednzn30tvvj2Ex9hXdAkSCB0sj3SKgMQpBgckP/NeoRS7tVVMdWXRLfDcl4\n9t5I6+6NM2muc/TjKcUgQcNFgLLa++SJuPXr2SD6lS++JLZuGoDKxD77+ZdeFhPj1dPHJ+Ya8f/e\n9oM+rohh0//fpqFDCs0WByFIMLmh7/LbGgoTE/rF1oa2pWmr2CNhjUFCsvm87D2g4SJAKaVpGh/+\n/H3RWhIWn7t5Ml4+Yg0WV7Jl40T89IsvzZz7wjf2xoOPHVv+DZAjSKA0WjP5qQ39DxIagoT+awxg\nf4SwtWE10mNPZv93nNgQyRp7oJxsuHhp5pztDQDl8+3vH4g7f3Awc+4X/8HlMVarrvCO0fOy6y6K\nnVsXfx9JI+IvP3evSj3aIkigNNITuWaLA1iRUBMk9Fyam9jQ79GPp5ja0L78tobqtt2RJMmaP0/D\nRYByazRb8aG/2ZM594yLtsQLnrGjTysaTGO1SvzSy67InLvvkSNx+92P92lFDBNBAqWQtlqFqQ0D\nsbWhlmu22PCw2HPzM9nj8cEIEmxtaF8r3x9h64Xr+rzl+iQAUB5/87VH4vGDi//9TyLil152xbpC\n6LK65rJz45rLtmfOfeTW78d8vdmnFTEsBAmUQjp3NCJdMv92YkMktfH+Legpy/ZIUC7WU+n88cxx\nMrGxTyvJsrWhfYUgYY39EU5ZbmuDMk6Acjg6U4+PffGBzLnrn3tBXLxzU38WNAR+6cYrolpZDFkO\nHZuPv75dyM6ZCRIohcLEhgHojxAR0aqMRbrkj1klbUaSNs7wDjotncsFCZPrCxLSiKg3muv/f0kt\nTn5H8pRGPebr9dM/b620gBHUOpQPEtbXKCvZfF7E0qBx/kSks0fW9ZkADIaP/t39MTu/+LvW1EQ1\nXvOSp/dxRYNv57bpQhPKT9/+UBw8aksuK6v1ewHQCa3CxIa1NWLruCSJRm0yxhqL5XXVxlw0xsf6\nuKjRUqxI2LC+z0vTeGDf0XV9xim7K+NRbS02EXxk75Ont8M854pWhNsk0uZCtA5n92pWt16wrs9M\nkkpUtl4Yrf2LY65aBx+JyhobOAIwGB5+4nh84VuPZs696sVPi3M29L9KddC98sWXxhe/sy+OzixE\nRES90Yr/8YXvx2+86tl9XhmDSkUCpVCoSBiA/ginaLjYX52uSOikZjU7QcJUj6LW4cci0sV9msmG\nbesOgyIiKluz2yNaB/eu8EoAhkGapvHBz92b2UG6c+tUvOy69W2HGxVTE7X42Vzlxle++3h8f6+K\nPZYnSKAU8kHCIDRaPGXZPgn0RNpsRGSmNiQR4+t/CO2UYpAwv8IrR1fr4MOZ4/X2Rzilui3bsDG/\nfQKA4fL1e5+M7z2Ubbz92huuiFrV4067rr/mgti9I/uFywc/f58+QizLnyxKoTWTCxIGpEdCREQz\nN27Q5IbeSedPZE9MTEdSGZy/9oohkyAhr3UoWylQ7VCQkA8k8iMmARgeC41W/NWt92XOPfvSrfHc\ny7ev8A6WU6kkhXGQ9z96NG6/yzhIigbnN2pYh/RENoEerK0NuRGQKhJ6Z0AnNpxia8PZ5R/w1zv6\n8fTn5IKE1qG9kaZaXAIMo89+9eHYf3jxv6GVJInX3Wjc41pcdcnWeP4V52bOfeRvvx/zC8ZBkiVI\noBTSfEXCoDRbjIiGrQ19k85lKxKSycHZ1hBRDBIqKhIKOj368ZRkanPE0l4LjXqkx57syGcD0DuH\nj8/Hx7/0QObcS59/QVy4Y7C+PBgmv3jD5YVxkLfc/lAfV8QgEiRQCoWpDYO0tUGzxb4pTmwYrF8q\nbG04s7Q+G+nxA4snkkpUtuzqyGcnSVLYJmF7A8Dw+V9fuD/m64vflm+YrMXPXG/c43rs3DodP3Fd\ndhzkp25/0DhIMgQJDL10YT6ivqTvQKUaydSm/i0op9gjwV/CvTLIExsiIlq2NpxRoRph885Iap0b\n4VWc3CBIABgmDzx2NL74nX2Zc6/+safFxinzk9frlS++NDZNL/57rC+04n9+4f4+rohBI0hg6BUm\nNkxviSQZnFt72akNut/2xOBXJJjacCbNXKPFTm1rWPy83OQGQQLA0EjTNP7yc/fF0t+odm2fjpc+\nvzO9dEbd9GQtfjZX2fHl7z4W9z96tE8rYtAMztMWrFFhYsMANVqMiEgrtWgltdPHSaRRadX7uKLR\nka9IyOyJHwC2NpxZYfTj1k4HCdmyzfyECAAG1x3feyL2PHIkc+6XbjTusZOuf+6u2L0j+7vTBz9/\nr3GQRIQggRLIVyRUpgen0eIpDdsbei5ttbJbXmLwmi22KmORxmIzo0raiKTV6OOKBkvrYHcrEqpb\nL8he7/BjkTb9+wcYdPMLzfjIrXsy5665bHtc/XTjHjupWqnE627MjoP8/t6j8RXjIAlBAiVQaLQ4\nYBUJEcs1XJxd4ZV0TP1ExNKCx7GpSCq1FV/eF0lie8MK0jQtbDXIN0dcr2RiQyQbti25aDNaR/at\n/AYABsIttz8UB44u/veyWknitTdc3scVldezLt0Wz7s8Ow7yr27dE7PzgvdRJ0hg6BUqEgYySJjK\nHGuq133FRouDVY1wihGQy0tnj2R7XFTHI9m0o+PXyVc55KsgABgsB4/OxadufzBz7sYX7I5d2wfz\nv/Nl8NobL49adbGC8sjxenwiN3KT0TNgX8/B6hWaLQ5ikJDf2iBI6LpBb7R4yslqlcU9nmW/N8Zj\nPmLh7P+M9SfuyRxXt+6MieaxiOYKb1hGdawWzYUzf2OysPncaC5pxZA8+f0Yv+Sq9i8SETE2GfWY\nOPvrAFi3//GF70d9oXX6eOPUWPz0j17avwWNgJ1bp+OmH7o4PvnlxQDnM3c8HD92zS4BzggTJDD0\nCs0WpwcwSMg31dMjTA3JdwAAIABJREFUoevSuROZ40Eb/XhKcQRkySsSFubi+J5vnPVljYe/kzlO\nq5NtvW+pzU9/Vhy//64zvqZZzzY+nX/ke9HavLqO3xsvf37EmCABoNvuffhwfOW72f35r3nJ02N6\n0rjHbvupH7kkvnTnY3Ho2MnfU5qtND74+fviX/zCcyNJkrO8mzKytYGhl544nDkexK0NjeVGQNJV\nw1ORMGJBQpvSEwczx5WN21Z45frkK5jywSQAg6HZasX7P5OtVtu9Y2O85LkXrPAOOmlyvBa/8A8u\ny5y78/6D8c09T/ZpRfSbIIGhlrZakc5kg4SB3NpQaLYoSOi2wujHge2RYATkctLj2SAh0xSxg5Kp\nzRFLJmfE3PFImwtduRYAa/f5r+2NR/Znqw1/+eVXRqXi2/BeedFVO+PKi7LT0T70+ftiobGKfYeU\nhiCBoZbOHY1IF/fJxcSGSGrj/VvQCpq1XLNFWxu6Lp3PbW0YloqElnsjbS5EOpudDd6tgDCp1iKZ\n2pS9fq7KCYD+Onx8Pj562/2Zcz/y7PMLD7V0V5IkcfPLroilOxn2H56LT9/+UP8WRd8IEhhqhYkN\nA9gfIWK5h8X5bABCR6VpKyK/tWFAeyTY2lBUeJCf3NTVgDDfVyW1vQFgoPzVrXtirr74rffURDV+\nMVdmT29cvHNT/IPnZ3sJffLLD8aBI74IGTWCBIZa/oEj2TCgyXRSiWYl/8DoL9yuqc9GpOnicW0i\nkupgNmKytaGoV/0RTsn/vaEiAWBw3PPQoUKDxZ+5/umxeaMmt/3yM9c/PTZOLf5eVW+04q9u3dPH\nFdEPggSGWr4x2iA2WjylkR8BaXtD1+T7IwxqNUJEsSKh0pzPhiAjqHX8QOa4W/0RTn/+tIaLAIOo\n0WzF+z9zb+bcRedtjBuuXd10HTpr49RYvObHn545d8f3noi7H/Tfz1EiSGCo5bc2DGKjxVM0XOyd\n4sSGwWy0GBGRJrVoJYt/FVfSZiRpo48r6r98RULS9YqE3NaGE34RAhgEn/vqI7H3yWKDxWrFI0y/\nveSaC+KSndkeQ3/5uXuj2bJ1d1T4U8hQa+WDhAHtkRCxXAm7IKFb0rnhaLQYERFJEi3bG05L01ax\n98mG7V29ZjK1KWJJmBMLc5HW/fkE6KdDx+bjY1/8Qebcjz7n/Lhi94BuYx0xlUoSr3/5lZlze/ef\niM/e8UifVkSvCRIYasUHjsH9j0szv7VBkNA1+YqEGOCtDRERjdz2hlpjpk8r6b909mhEa8kYqbHJ\niPGpld/QAUlSiWQ61yfB9gaAvvrw39wX85kGi7X4hZde3scVkXf5hZvjxVefnzn30dvu///bu/P4\nOM4q4fe/qupV3erWvlmWbcn7GsdLnMQBYickBJtAIAmTQJgJH8Jk7tyZYTYyHy5kGAbum5mX+w7D\nzgwTdggQIMRZwUkwTuw4jh07trzKWqx9a229d1XdP1pWq7VYbVtSd8vn+0mnu56u6j6ySt1Vp57n\nPHT2BdMUkZhNkkgQWW3swX42DW3QYvIhO10MIBLTR25GcDDped3qTHr+Sm4zUb0gZk3uGmiNDMzA\nu2QHc2jMsAZXAYoy83OEjx/eIAUXhRAiXU409HLgRGdS213vqMbjyrwpvq92d7+rhhy7ZWQ5EjP4\nwfMnMa/yek9XA8vUqwiRubJpaENMS76qKj0Spo+uGzS0JU6+K/yDjJ6joXUAouHpOTnPr5n+L8aI\n1ZO0bIsOTrLm3DfbMzZcMLZHghHwoc3KOwshhBgtphv86HfJBRarSt3jphwUmcHrtnPPtsV877mT\nI221DT5eO9bOjWvK0xiZmGnSI0FkLTMajk/zd4Gqxcc6ZygptjhLTBNNTx4aoFty0hRMaqI26ZFw\ngTFBj4TZIAUXhRAiMzz9agNtPcnf4x959zJUdeZ7p4nLc9PacpZXJSfkf7b7DP3+SJoiErNBEgki\na42bsSEnD0XJ3F16XI0Emf5xRqhGGNVMVAw2FAuGar3IFukXHTu0ITqIaV59VY9N05z1GRsuGDt1\nrBnok26ZQggxyxrbB3lmX2NS29a15Sye501TRCIViqLwsfcsx2pJHIf7QzF++vvTF9lKZLvMPesS\nYgpj53rP5PoIAIZqwySRTVfNGIpxdU/zNxMsY2pPxCxOmIUx9ldC1xzoo5IdqqljDl2FV8QjQYiO\nSrCpGorTM/n608mWA9qohJMehbB/8vWFEEJMq5hu8N1nTmCMSuJ6XTbuuVkKLGaD0vwc3r91UVLb\ngROdvHW2O00RiZkmiQSRtcbN2JCTuTM2AKAoMgXkLBg740Esw4c1AKAoRMfUSdD72tIUTPoYY3sj\nuPJnrZeRoijjkpFjk5VCCCFmzq7XGmjuSp516YHbl+F2ZnavQpHw7s3zqSpNninrhy+cIhiWC2dz\nkSQSRNYaV2gxw3skAMRkeMOMy8pEAuPrJBh97WmKJH0mmrFhNo2bAlJmbhBCiFkx0ZCGLatKWb+k\nOE0RicuhqSp/9p4VqKN6gvoGwzz5h7o0RiVmiiQSRNYa1yMhCxIJ0iNh5o2dVlO3OCdZM7NErJJI\nMP09ScuzNWPDyPtJwUUhhJh1Md3gf549gW4khjR4XDbuu2VpGqMSl2tBWS63bZ6f1PbyoRbONEty\nfq6RRILIWuOKLWZFIiH5pNaiBydZU1yu7O2RIEMbxvdIKJzV9x/bA8IY6JxkTSGEENPlmX2NnO9M\nHtLwsdtkSEM2u3PrIkryEse8JvC9504SjV19haTnMkkkiKxlBJIzm0pOFiQSxg5tkB4J0y5rEwlj\neyT0d2AaepqimX1mLIoZGhzVMr5mwUxT3IUwuiZDeAhTCi4KIcSMaeoYZNdrDUltW1aWsn6pDGnI\nZjarxsfeszypra0nwFN769MUkZgJkkgQWWtODG2QGgnTyzQnnrUhCxianZhmH9WgYwx0pC+gWTZu\n2kenB0WzzGoMimZByS1KajP6r57fgRBCzKYLszSMG9JwqwxpmAtWLMjnprXlSW3P7W+ktqF3ki1E\ntpFEgshKpmlgBvqT2rJjaIP0SJhJqhFFNROVgQ1FxVDtF9kis4yducHobUlTJLNvXCJhlusjXKB6\nSpKWZXiDEELMjGcnGNLwgAxpmFPu2baYPLdtZNkE/uvpWgb8kfQFJaaNJBJEVjKDA2CO6vZtd6FY\nbJNvkCFikkiYUdqYYQ26lgOjKgdnunHDG3xXTyLBSPOMDSPv6ylNWjavol4hQggxWxraB3h6zJCG\n61aWcq0MaZhTXA4rn9i5itFHYv3+CP/9TC2GaU66ncgOkkgQWWnstGxqFtRHgPE1EiyxIJhSeGa6\nZOuwhgvGFlw0epvTFMnsG9sjYbZnbBh53zE9Eky/DzMWTkssQggxF/lDUb7x62PJQxpyrNwvQxrm\npBUL8tlxw8KktmPnennxwPn0BCSmjSQSRFYaP2ND3iRrZhZTtSb1SlAwsUX6L7KFuBTZWmjxgnFT\nQF4lPRJMw5jgbzpNPRKsdpSc5M8TY6ArLbEIIcRcY5gm//10Ld39yT0yP3rbchnSMIe9b+tCllR6\nk9qe/EMd9W0DaYpITAdJJIispPcmZzFV9+xOE3clwvbkEyR7WIrOTBeLnt2JhKhtgpkbYnN/HKEZ\n7E/umWPLQbGlrzeJDG8QQoiZ8dz+Ro7U9SS1bb+2kg3LZEjDXKapKg/tXIXLkSiirBsm33rqGMFw\n7CJbikw2uyWxhZgmesfZpGW1pPqKX9M0TfyhGL7BMAOBCGp0iNYhDU010RTQFLAMP3ZYTNTLHHof\ntufjCrSOLNvDvousLS7F2KEN+hUMbQjGFIaiCqGYQkgfvo8pnH61g8ZGJ7oJVtXEpplYVYbvTWzD\nj/PsBl67cUn7SbzHihOLPvxzmAZGfztaYdVl/xzZwBxTH0FNU2+Ekff3lmK0nxpZNvql4KIQQlyp\nE40+frXnXFJbdYWHe7cvTlNEYjYVeh382R0r+Nqv3h5p6+oL8YMXTvHQzpUoWVTTSsRJIkFkHdM0\nxiUStNIlKW5r0tkX5Mz5flq6h/ANhukbDOMbCtM3FCEaG1uvwD3h66hK/EQx325Q4NApcBjkOwwK\n7DqOKf6qIvbkeg42SSRMm7HFFlPpkWCY4AupdAU1uoIX7jX80Uk6bLX2AKkV9lQVc3gfie8n+XaD\nQqdOocPAMsnLR2weLMFEQsTwtcz9RII/+epUumZsuGBcnYTBLkxDR1G1NEUkhBDZzTcY5ttPHWN0\nfT2308rDd67GokkH6avFtUuL2XbtPF46lBi6+XptBysX5nPT2oo0RiYuhyQSRNYxfG0QGXXCaHeh\n5pVNuG5MN2jqGOJMcx9nm/s509zHQCB65TGYCr0hjd6QRl1/8pg+t9Wg0h2jyhOjyh3DY0+uShux\n5WGioBBvt8YCGMFBFEd2FIzMZKnUSNANaPFrNA5YaBq00B3U0M2ZyYIbpkJPSKMnpAGJ/URTTMpd\nOpXuGJW5McpdOtbh46ioNReCia70V8MUkJkyY8PI+zvcYHdB2B9vMA3MwW4Ub+nFNxRCCDFOTDf4\n5lPHko6/FOChnSsp9Dom31DMSfduW8zp8/00dyWm/vzx705TU+GlosiVxsjEpZJEgsg6eseZpGWt\npAZFSWSzh4JRDp3u4o2TnZw530dkXC+DmTUUVTnps3HSF79q7bXrVLl1qnJjzM+NkWPViNi82COJ\nmSf07kYs+XP7qvNMU4wYmpE4SDFR0DUHpgl9YZWGAQuNAxbOD1mIGuntPqebCs1DFpqHLNAeTyyU\n5ehU5sbYYMtjdDkifY7P3GCaZsbM2JAUg6cUoyvRBdcY6ECVRIIQQlyyX75Sx9nm5MLSO29cyOrq\n7KlvJaaP1aLx53eu4l++/waRaPwYPRI1+Oqv3uYzH90gRTeziCQSRNYZP6xhMYFQjMNn4smD4/W9\nSVMKXQq7TSPfbSfPbQMjRr+vD91Q0E2ImaAbCjEDQnrq3fD6wxpvhzXe7oknFua5YnzQVUQNoxIJ\nXY1Yltx0WTGLOFs4eUrQsJbDH1qcnO2z0B+59C7pmmLitRk4LSYOi4lDi98vqqnE19yEpppEdYWo\noRDRid8bClEdAjGV3pBKIJbafqKbCi1+Cy1+C81aKf84KpMw5xMJfh+MLiipWcCRO/kGs0T1lCQl\nEswBqZMghBCX6uDJTl58I7lA9qpFBbzvxkVpikhkgooiF/ffspTHnzs50tbRG+CrTx7l7z98DVaL\nDCXMBpJIEFlnbCLht6dVfrd7LzE9tZ4HdqtGdYWHmnleSvOd5OXayXfbyc+147Qn/iSUkI+397wy\n4WuEdegNafhC8RPG3nD8cV9YnbKbfIvfwqvRMmrciZ8j3FGPPaXoxWScoeQp+o4EynjTn9q/qsti\nUJSjU+wwKM7RKXbq5DsMtAl+les3FXM4fHb8ExMIxeL7SW9YHR4Ko9IZ0BiarP4C0Kl7MUwFVRlO\nhg11859PHOCGa6upKcslP3du7SnG2BlYvOUZUXBp7DAGY6AT0zQzIjYhhMgGbT1+/ufZE0ltBR47\nD+1ciXq5FavFnLF1bTknm/rYd7x9pO1Mcz/ffeYED71vVRojE6mSRILIKnpgALM/8YGjmwq7G6zE\nmDyJ4HXbWFKZx5J5XpbM9zK/xI2mXllhH7sG5S6dcpeeHJ8BbQGN84Px8fdtfg1jgsRCQyx5mqNI\nZyP/64dvsHlVOZtXlEq3rkswFIzyxtEuKnw9eEf9Wk9HyyfdxmU1WJgbY4EnPtzEZb28HixTcVig\nwq1T4daB+LAL04SBSHxow/nB+PCGgUgi8CgWuoxcSrXE3Mq9TQ18oz4+lnDxPC+blpewYVkxBZ7s\nH1tq9IxJJBTOT1MkyZScPLDYEr0lYhHMQB+KS2qZCCHEVHoHQnz5ibcIRRLHSZqq8PD7V5Obk1rB\nYjG3KYrCn75nOT0DIU6fT/QqPXCik0Kvg7+4e30aoxOpkESCyAqmaXL6fB91B/cyegBAq55PhPEn\n3eWFOWxaXsKmFaVUFObM2lVETYVKt06lW+f68jBRPd4DoWnQQtOAhc5gvKtWj+Fm0HCQq4YAsCsx\nhtqb+FHLID/9/RnW1hRyw+py1tYUYp2svP9VLBzVOXK2m33H2jlW34vNDPOlvORx9qejiQKcmmIy\nzx1joSeePChyGKTrwrKigNdu4rVHWVUYTy4MRBSaBy00Dlqo77fQruclJRLKtT6a9CIAzrb0c7al\nn5/uPpP1SQUzEsAc6k5qUwsq0xRNMkVR4sMbRg0tMfs7QBIJQghxUQOBCF9+4i16B8JJ7R/evoSa\nCu8kW4mrkdWi8pd3reFLP3yT9t5Ewezn9jexqDKfjYuljkYmk0SCyGiRqM6+4+3sfrOF5q4hdjhP\ngTPxfP2oK/vFeQ42ryhl84pSKotdGdEF2arBQk/8BJZ50B9WON1n5bTPSkOsiDW2xEnKAksXrXo+\numFy+Ew3h89043JY2LyilOtXl1FT4cmInylddMPgRKOP/cc7ePN0F+FRVzkWW9sTQwGAtlgeQ6aT\nak+U5QVRarxRrBk83M5jM1lZGGVlYRTDBK0zBxIzQFKu9U243dikwsblJWzMoqSCMab+g+IuQrFN\nPWXnbFE8pTAqRmOgA61ieRojEkKIzBYMx/g/Pz9CW0/yLErvWFfOtmvnpSkqkcncTiufumcdX/zB\nwaSZPb715BH+6kNrWVtTlMboxMVIIkFkpGA4xsuHW3jxQFPSh8oiS/I4+BbKuGVDJdevLmNhWW7G\nn2h77SabSiNsKo1g68mFwcRziyxd7AsvTVrfH4r/O7x8uIWSfCc3rCpjy+oySvKcXA1M06ShfZD9\nxzs4cKKDfn9kwvWWWtuSlgOOYj65cBCnZWaGLMwkVQGn252USFhXHOa4mceZ8xMnFCCRVPhZFiUV\nZnpYgwlEYvqU6026vTv54MXo75jw9WZ3XhghhMhM0ZjOV588SmP7YFL7hmXFPHDb8ow/RhPpU5zn\n5K/vXsdjPz40MtuaYcI3f3OcT9+/noVlnjRHKCYiiQSRUYaCUX5/8Dy/P9hMIBxLek7FoMqS3A36\nox+5A2dBdk7JpubkJyUSVnv6KbXn0NEbmHD9Tl+Q3+yt5zd761lS6eWG1WVsWl5CjmPu1VPo9AXY\nf7yD/bUdSV3dJrPS3pG0nFtQQDALkwgXRK3JsxYUGj38f596JyfOdHLwVHx2kvq2gUm2njipcO2S\nIooyIAFlED+5N/UY9LUmPad7KzCu4MR/3HuZJg0X+XeakmmlChXlQqogEuD8+XZ0S3KviTVLDCYY\nYSWEEFcN3TD41lPHOdmUnPBeuTCfh3aukuKKYkqLyj188s5VfO1Xb2MOH8KFozpf+cVRPvPABoq8\n6T+GEckkkSAyQt9QmBcPnOflwy2EoxOfSMy39GJTEs8prnwc+SWzFeK0C9vzMIELX62ucDdffGA1\n9b06+46183ptB/5QbMJtzzT3c6a5nx//7gzXLCni+pWlrFpUgC2T++9PYSAQ4Y0Tnew/3k5d69Qn\nf163jetWlHJDtZ383yUOXEwUQo7sHlMXtbpB1cCI7+9moA89OERRnpPbr6vi9uuq6O4PcvDkpSUV\n5hW7uGZxEetqiqiu8KTlwE7XDRraBnAG2ikxEn/PMc1JS78GA1dw4j9Gfs0VJpMUjbA9H0e4Z6TJ\nEerB786c4RdCCJFuhmnyvWdPcvhM8sWe6goPf3nXGqn1JFK2fkkx992ylB//7vRIW78/wv/+2Vv8\n/b3XZMQFEZEgiQSRVt39QZ57vYk/HmmbdPpGu03j5vXzuNUVhLcS7Vrp4qzuJmeqVqJWD7Zo4sTJ\n6KqnZv4aaiq8fHj7Et6u6+G1Y+0cqesmpo8/KYrpBgdPdnLwZCc2q8qaRYWsX1rEusVFuLKgp0JH\nb4DDZ7p562w3Z5r7RjLQk3HYNDYsK2bLqjJWVOWjqgrR03sJjVonbM/HVDP/Z78oRY0X+utLDNmI\ndp8HR6IQYZF3fFLh4KlOzl0kCdPS5aely88z+xrJzbGytrqQdYuLWLmwgBzH7H4dOAPJvUiCOaWk\nrQLmRYQdBUmJBHuoF787M2aWEEKIdDNNk5+/dJZXj7UntVcUufibu9fhsMmphrg02zdU0tUX5MU3\nEsMfO31B/t8fH+Jv772GeUWuNEYnRpO/bpEWbT1+nt3fyP7jHejGxGePLoeFWzfOZ/vGSlwOK8Hf\nv8Do6/Na6eLZCXYGhe35SYkEvbMOy/w1AFg0lfVLi1m/tDg+xeHJTvYda+dsS/+ErxWJGrx5uos3\nT3ehqQrLqvK4dmkx1ywuyphx8oZhcq51gMNnu3jrTPe4YkwT0VSFtTWFbFlVxrqawnG9LmLNx5OW\nQ47kqTWzlZpXlpRIiHQ2QdXEMxpcTlJhMBDl1WPtvHqsHUWBqtJcllflsawqn6WVeTObWDBNnMHk\ng86As2ySldMrZC/Ey5mRZfuopIIQQlzNDNPkly/XJZ3wARR5HfzdvdfIVNbist2zbTH+sM6rRxND\nIH2DYf7Xj97kb+5ZJ7N/ZAhJJIhZ1dQxyK59jbx5spPJLj57XDZu31zFO6+pwGmP76KmaaK3n0la\nTytdMsPRzrywPZ/cocaRZb2zbsL13E4rN6+fx83r59HhC7DvWDuvHWunuz804fq6YVLb4KO2wceP\nXjxNkdfBkkovSyrzWFLppbzIhToLV3/DEZ2G9gHqWgeoa4kPxxgKRqfeEFha6WXLqjI2Li+Z9GDE\nNE30ltqktpBzbiQStLxyYhweWY50NUHVDVNuNzqp0NMf4uCpTg6f6eZscz/GJF0+TBMa2wdpbB/k\nhQPnURRYUJrL8qp8auZ5WVDmptDjmLYeQLZIPxY9se8aikbIkZlVmcOOgqQhSLboIKoewdBkHnQh\nxNUrEtX5r121vHkquQi2x2Xj7z58Dfm59jRFJuYCVVH4u/s3EI3GOHCic6TdH4rxv3/6Fn/5wTWs\nWliQxggFSCJBzJKzLf3seq2Bo3WTX80r9Ni5/boF3LS2fNxVZ9PfixkYVcBHs6IWVs1UuLMmbE/+\nENQ7z2GaBooy+XjC0vwc3n9TNXduXcTZln4O1HZy6EwXvsHwpNt094fo7g+x73i8O7nLYaFmnpea\nCg/F+U6KvU6K8px4cqyXdbKoGwZ9gxF6B0N0+oKcaxvgXMsA5zuHJj15nUhlsZvrVpZw3crSlIrq\nGL5WzGCih4YxPKZ9LlDzkq/QR7rOX3I9v0Kvg9s2V3Hb5iqGglGOnevhrbPdvH2ul2B44vobEE8s\nNLQP0jCq8rbLYWFBWS5VpbksKM1lQVkuJfnOy0pIOQPJvRFCzpJ4TYgMZKpWojYvtkhiP7OHewjm\nlKcxKiGESJ9+f4SvPnl0XK83p93C396zjtJ8qSMjrpzVovLQzlXkOKy8crhlpD1egPEID+1cxcbl\n2VsrbS6QRIKYMaZpcqLRx67XGsZV8R2ttCCHO7ZUcf2qMizaxCfQ43ojFC9C0bJ/941Z3RiKBdUc\nPqkL+zH7O1Dypj5JURRluIdBHvfduoSG9kEOne7i0OmuKYcM+EMxjtb1jEvs2CwqhV4HRV4n+bl2\nNE1BRUFR4u+nKPEsccww6BsM4xsM0zsYpm8oPGV9g4loqsLS+XmsX1LENYsvfVYBvWXssIYiuEgS\nJpto+cn7QKSrEYuho1zmCbfbaWXLqjK2rCojphucbe7nrbPdHG/opaXLP+X2/lBspJfLBVaLSnGe\nk5I8JyX5TkrznZTk51Cc76TQY0dTJ/5djB/WkNkzr4TsBUmJBNdQsyQShBBXpZZuP1/5xZFxPSI9\nLht//aG1VJXmTrKlEJdOVRU++u6luJ1Wdr3WMNIe002++dQxHggt453XzEtfgFe57D8TExnHNE2O\nnO1h176Gi47Rrix2s+OGBWxcVjJl9Xi942zS8lyojwCAohC25+MMJboG6p11qCkkEpJfRmFRuYdF\n5R4++M4a2nr8I0UM61sHJq1DMVYkZtDWE0ipdsHlcto11lQXsn5JMWuqC65o+srYuETC3BjWAKC4\nC8DqgGj8YM0IDhFregvrwg1X/NoWTWX5gnyWL4j33hgIRDjd1MfJJh+nmvpo6Z46sQAQjRm0dvtp\nnWB9RYHcHBt5Lhtetx2vy4bXbaNY7WP9qJNyExiwl6GaGVlrEYBgTjmewfqR5ZxAK5aon5hVCj4J\nIa4etQ29fP3Xx8b1aJtX5OKv714r0/OJGaEoCne9oxq308rPdicuLJomfP/5Uwz4I7z3hoWzMmRX\nJJNEgpg2hmFy8FQnu15rpLlraNL1aio8vPeGhayrKUy5G/34REL210e4IGwvSE4kdJzFunTrFb1m\neaGL8kIXd2xZQCSq09A+yJnmPs4093O2uZ/ARbq1T7fSfCfVFV4Wz/NQXeGlssQ16ZXqS2EaMfS2\nU0ltc6U+AoCiqFirNxM9tWekLXp897QkEsby5NjYuLxkpIvggD/C6fN9nD7fR2PHIE2dQ4QjE0/L\nOhnTjL/OgD8CnYnPgxvsp1k/6vy7IVrMV46VoCkmDs3EYZngPqnNwGExcQ4vW9WZT0CEHEVERs2w\nogC5A3X4CtfO7BsLIUSG+OORVn7wwqlxFyZWLSrg4TtXz/rMP+Lq8+5N83E5LDz+7MmkYbO//mM9\np8/38fEdK8lzS22O2SR/9eKKhSM6rx5r43dvnKfDF5x0vRUL8tlx/QKWL8i/pHH4ZjSM0dOU1KaW\n1lx2vJlm7Jj+yQouXi6bVWPp/DyWzs8D4lWWW7v9nGnup70nQHd/kK6+EN39QUKXeLI4mifHSr7H\nQUGunXnFLmoqvFRXeMjNmZmidEZn/cjVegBdtRO1zq0uldZV25MSCXpLLbqvBS1/ZrvxeVzJiQXD\nNOn0BWlsH6SpYzCeXOgYSrlw5mirrcnVvY9F4zNR6KaCP6bgv8Qcl8rkSYfGA134umw4NJMcq0Gu\n1STXZnDJU5rkswbwAAAgAElEQVQrCgPeGoq6E8Uv3UNN9Octl6KLQog5LRzR+eUrdew+1DzuuXde\nU8H9ty6ddFiqENPtxjXl5DgsfPM3x5OmjT/e4ONz3z3Ag3es4JolmVm8eS6SRIK4bL0DIXYfambP\nW634Q5Mf/a+rKWTHDQupmXd5U7XoXefATHxYKN5SVKfnsl4rE0XGJBKM3mbMaAjFOjNTNqqKQmWx\nm8pid1K7aZr4QzG6+4N094UYDEQwzHi7OXxvmGAOz7fhddkoyHVQ4LGTn2vHapndYnmxiWZrmGPd\n2rSiBWhlS9HbT4+0RY/vRtv6wKzGoSoKZQU5lBXkcN3KRD2DoWCUrr4gHb4Anb4gXb4gHX1BOnsD\nDATGJxlsRFlqbUtqOx6deErLVBkoBGIKgRgwpt7ooc4uYHxXW5fFINdmkGsz8dgMPDaDPLtBvkPH\nYzOZaKSV31VJnu/EyGwTqqnjHmxgIG/pFcUvhBCZ6kRDL48/d3LCGaLuuXkxt22eP22z+QiRqvVL\nivm7e9fxtV+9nXT+MRSM8p9PHuXma+dx782LxxVuF9NPEgniktW3DfDiG+c5eLJz0rH3CrBxeQnv\nvX7BFRfeGT/t4xypjzDM0GyonmKMgeHhDaaJ3lWPpWLFrMahKApupxW308rCssxP1IwttBjM0OkD\nr5R11fbkRMLpV7Fv/hCKLf1VsS/sL4vKx+8v0ZjBYCBC31CEfn+Y/qEI1vajWM8nkoI+w8Wg6kEz\nTHRz9g5G/TEVf0ylfYJSIJpikme/kFgwKLDrFDoNCh0wmLuI/L4TI+vmDp5jwDt3ekcJIQRAMBzj\nFy+f5ZW3Wsc9Z7WofGLHSqmWf4VUPYQ97MMcVDFtZZKQuUTLqvL5/IOb+e9dteMKur98qIVTTX18\n8n2rmF/inuQVxHSQRIJISTiq8+apTl55q5Wzzf2TrqepCltWlnLH9QsoL5yeQmRzuT7CBVrxwkQi\ngfjwhtlOJGQTMxpC70geAjKX6iOMZlm0ASUnLzH9aSxM9PSr2Fbfmt7ApmC1qBR4HBR4Ej1rgq88\nx+i+S5q3lIeq4/UTogaEYgohXYnfj36sK4RiKkF9fHvMmN6DL91U6Alp9IQ0GPNRV2pbxd+7zmBT\n4j+FRQ+TM9QywasIIUR2Onauh+89f5LegfFTShfnOfjk+1ZTXZH5FxuuVI4WZVnRZUxHNZZpovXW\no/qaUP09qIEeFH8PaiRepNg8BoH8CiyLr8e6eAtq7tw8lpkJBR4Hf//h9Tz3eiO/+WN90sXN1m4/\nX/j+Qe56RzW3bKyU4TczRBIJ4qIa2wfZc7SV/cc7LjrvvMth4V3r57Ht2kryc6ev0IlpGuNqBszJ\nRELRAqJ1b4wsxxoOYVt3B8ocmcpwuultp8BM1HNQc4vRLem/Qj8TFNWCdcXNRN789Uhb5PhurKu2\nZ9X+oXc1EDvzWlJbMCcxTMKqgtVmksulHbhFDQiPSS4Eh++9ZZU0NrYSiqkMRRUGI/H7eJ+pS9cR\nyWGfZTHvdJwcaYt11vPFn51iwbwCFpblsqAsl4pC15Qz0QghRCbxh6I8sfsse99uG/ecAmzfWMkH\n31GD3XZ1dBdXYyF6at+YesWLsIV6yfcdxxHuveh6hq+VyBtPEnnjSbTSJViWXI+1ejOKQ66mT0VV\nFd57/UJWLCjgO789TmdfolZbTDf4+ctneelQMx94RzXXrSyVmR2mmSQSxDiBUJT9tR3sOdJKU8fk\nsy8AlBXkcOum+dywugz7DIxFip54BcKJqeVUew5q/tybv10rrU5aNjrPET25B9uKd6UnoAw3tj6C\nVr4UImkKZhZYV7yTyOGnwYgn88z+dvSWWiyVq9McWWpMPUboD/+dVOskpjkITcNwlAsJCPcECYj1\nN5Zy2DyX1KabMBSJJxUGhm/9ERVfSMUXVgnGLp6ceSW0gpvsp1CV+PuVa31Yes6wuz1RANNmVakq\niScVFg7fyiW5IITIQEPBKC++0cTuN5sJhscXXC4tyOHBO5azpDIvDdFlJ0t0iDxfLa7A+KTMVPSO\nM+gdZwjv+xm2De/DtvZ2FFVO16ZSXeHh0T/bxE9+f5pX325Peq67P8R/PV3L86838cF31rCmukCG\nkkwT2TMFEB8Pd6SumzdPdnH0XA/RmHHR9VcuzOfdm+azurpwxrJ7xmAX4dd/ntTmXHxtVl2FTZWW\nV45WtQ696chIW/j1n2NZcA1qjnx5j2aG/cTOHUhqs5QvhcaLZ/yzmZqTh3vF9Qwd/+NIW+TY77Im\nkRA5/DRGb3LF797CdZCGv2VNAa/dxGvXgfEHzaEY+MLaSGKhN6TRHVTpC6uYKPQauRyJVLHe3jiy\nzc2OWk5GE4mESNTgbEs/Z1sSYyMkuSCEyCT9/ggvHGji5UMthKPjPwsVBW7bVMX7b1okRetSpOph\nvH2nyB1sQJmid52JQsySg1UPJiXZR+gRIgd+Sezsfhw3/emcqw82E5x2Cx9/70rWVBfy/edPjetJ\nfb5ziP/4xRGWzc/jQ++quewi8CJBEglXsUAoxpGz3Rw81cnb53qTplGZSG6OlRtXl3PTuvJpq38w\nGdM0Ce35XtL0flgdFN78EXxz9Mqz48aP4m89AbHhHzASILzvpzi3P5zewDKIaRoEX/4vTL8v0ahq\nWMpq5nQiAcCz8T1JiQS96SjGQCeqJ7MLXundjUQO70pqG3JVEswpS1NEF+ewQLlFp9yVfGAdNaA3\npNIV1GgNLGG9nkgkLLe2UaH10qoXTPq6EyYXLCrzS90sLPVQVeZmYZmH8sIcGcsphJgxvsEwz+1v\n5A9HWie9aFRR5OLP7lhOTYWcaKXENHEPNpDvq0U1Jx4GHHCWEnKWELW6iVlcxCxOUFRWbd5ItOU0\n0bP7MMbUBIP4TF6Bp76IdeW24ULL42ciEsk2ryhlWVU+T79azx/eah1XGP7U+T6++MM3WVdTyLvW\nz2N1dQGaKt+7l0MSCVcRc3gu+NqGXo7W9XC8oZeYfvGMqQKsri7kHevKWbe4aNYOcKMn/zCuKr/9\nunuweIuha3BWYphtam4R9o13Ed7/s5G2WN3rxJbeiGX+2jRGljkiR55Hb3orqc266hYU+8wmtjKB\nfd5S1KKFGN0Nwy0mkdqXcGz5cDrDuijTuDCkIXFSrjhz8RWsSWNUl8eqQmmOQWmOAYUuQm0FSeNe\nPzq/iWfUxTS2DzI4wdSXE4nEDOpaBqhrGRhps2gq80vcIz0XFpTmMq/YJckFIcRli+kGtQ0+Xq/t\n4I2THZMe+7mdVm7bPJ93b6rCapHPnFRYIwMU9ByZtA5C2JaHL38l4UkKQqsON7ZV27Gt2o4x0En0\n7D4ib7+YNKwXTKK1u4k1HsJ+w/1YFm6QrvlT8LpsfOTdy7h103x+veccB050jlvnSF0PR+p6yHPb\nuHFNOVvXlFNaMDfrbc0USSTMcQP+CLWNvdQ2+DjR0EvPBFV4J1LkdbB1TTlb15YnVV2fDcZQT9LJ\nNIBWsQLrVVAvwLr61nhWujtxtTO09we4PvRFFOv0FbHMRrG2U0Te+GVSm1pcjX3z3RCdm8ml0RRF\nwbZqO6E/fHekLXrqj9g3fgDFkpn7RuTwMxg955PaHNfdjdHQk6aIps+AZzGOrsQQm4rB4/zfN1+P\npWYrvQNhGtoHaewYiN9fQnIhphvUtw1Q35ZILmiqQmVxPLmwoCyX+cVuKopc5DjkK1wIMTHDNDlz\nvo/XT3Ry8GQnQ8HJP4M8Lhu3b67iXesrcNjkcyUlho63/wze/tMTDmOIWnLoy1tBwDUvPk4kBaqn\nBPu1d2JdcTPh/T8bV6DY9PsI/e5raPNWYb/+PrSCeZO8krigND+HP79zNe+5bpAn/1DHsfrxCZ++\noQjP7GvkmX2NLK30ctO6CjYuK7lqCoteCfm0mEMMw6Stx0992yD17QOcbe7nfOfFiyWOVuR1sGl5\nCRuXl7CwLDct2c74kIbHk4c0WOw43vHgnKyNMJaiajhu+jMCv/k8mPEvJnOwm/Cbv8Gx5d40R5c+\nRqCf0O5vJo8jtLtw3vIXKJoFUjtHy3qWmutQ9j+BGR7+uw77iRx9Adv6nRl3dULvOU/k0G+T2iyL\nt2CtWgMNr6QnqGkUzCkjanFhjQ1fNTJ1Qi99G9tAFwXrd1LodbBhWfwKlGmaycmFtkEaO1JPLuiG\nSWNHfBsSZVQo8NipLHYzr8jFvGIX84rclBXmzEjhWyFE5ovG4onIw2e6OHCiE9/gxS8e5efaec91\nVbxjXYXUQUiVaeIIdVHQ8zbW2PhjbEPR6MtbzqBnESiX92+qOj04b36I2JIbCe39PuZA8tV0veU4\ngSc/Gx/usOH9MrtDChaU5fK3917DiUYfv3ylLilZP9rp5n5ON/fzwxdOsaTSy8qFBaxcWMD8UrfM\n+DABSSRkqZhu0OkLcr5ziPq2ARraBmjsGJqwYM7FlOQ748mDZSVUlbrTfjISO70XvflYUpv9urtR\nPVfPvLpa8UKsq24leuzFkbbo2y9gXbwFrWhBGiNLD9MwCL38bcxAX1K78+ZPoOZeedX/bKJYbPEZ\nHN56ZqQtcvBXGH1tOLY+kDFjJ81YeIIhDR4cN3yEOZP1URR8+Ssp6UqeHixy8FeYg13Yb/rYSKVt\nRVEo9DrGJRd8g2Ea2+MJgobh+/6h1IvA9A6E6R0Ic7QuuYdHgcdOaX4OZQXxW2lBDmUFTgq9DhkH\nKsQcEgzHqGvp53RzH6fP93OudWDKelcQv3B0x5YF3LimXIYwpMo0yQm04uk/gz3SP+EqAWcpvYVr\np206akvlKlwf+lcih35L5MhzSd+pmAbR478ndnY/tk13YV3+ThRVkkFTWbEgn//ngQ3UtQyw52gr\nb5zonPDcKRIzON7g43iDD6jD5bCwYkE+KxYWsLwqj9L8HCmWjCQSMlpMNxjwR+j0BWnvDSRuPQG6\n+oMXLlhfEoumsnR+PMO2elEB80vSnzy4wPD7CO37SVKbVr4M68ptaYoofeyb7iJWfxDTP9wFyzQI\n/fF7OG/7G9Scq6v4UeTQU+hjpnu0XbMDS9U1aYoovawrt8XHT+qJE/LY2X34O87i3P4wWkn1Rbae\nWaYeJXriFSKHd2EGkw+07FsfiF81Cfkm2Tr7BF0VdJvXUtR7BIzEgUj01B8xhnpx3vp/odgmPqBU\nFIUCj4MCj4P1SxOJ0r6h4eTCqATDVFcVx7qQYDjRmPxvrSoK+bl2Cr0OioZvhZ74fYHHgddtk27N\nQmQg0zTpG4rQ1uOnrSdAa4+fc60DNHUMpnws6HJY2LCsmOuGC9HJSVCKDB23/zye/rOJHmhj6Kqd\n3sI1BHIqUh7GkCrFYsO++UNYFm8h/NqP0VtPJD1vhocI7/0BkcO7sC7egmXJDWgFldMaw1yjKAqL\nK70srvRy3y1LeONkJ3uPtnGmeeIEEYA/FOPgqS4OnuoC4sWSy4tcVA73BqwsdjGv2E2e25Yx51Wz\n4YqPGOrr63nkkUfo6+sjLy+Pxx57jIULF05DaHOPYZgEIzH8oRjBUAx/KEogFCMQjtE/FMY3FKFv\nMIxvMIxvKMygPzLF5DFTU4h354l3zcln8TxvxnZfCx/4JUSCiQbNdtUMaRhLsTpwbP0owRe+MtJm\ndNXj/8nfYqnehG3ldtTSxXP2w8qMBIjVv0n07P5xRTe18mXYNn4gTZGln+ouxHHzQ4Re/i/QE1ev\nzcEuAk99EdumD2Jbd/us/t2YRozoqb1EDv02kfwaxVK9GeuijbMWz2zyu+dTde31BP7wvaTiWHrL\ncQJPfQn75g+hlS+dNKEwVp7bTt5iO+sWJ3rbDPgjNA0nFc53DtHS7ae9J4BxidlkwzTpGQjRMxDi\n9PmJ13HYNLxuO3kuG3m5drwuG163DbfTOu7mcljlZESIaWCaJv5QjL7BMH1D8WPAvsEwHb4gXf0h\nzncMEopcWo9TALtVY/2SIjavKGV1dYEUbU2BaZqYg93o5w9R2H0YZ6ADzZg8mTvoXkBf/koMzTaj\ncWkFlTjf+4/EGg4R3v8zzMGu5Lj9vUSOPEvkyLOohfOxLrkBS80WVFf+jMaV7Rw2CzetreCmtRW0\n9fjZ+3Ybr9d20DtFTblIzBhJ+o/mtFso8NjJc9vJd9vJy7WT745/nxZ6HBl1AXc6XHEi4dFHH+W+\n++7jzjvv5KmnnuJzn/scP/jBD6YjtoyiGwa/2nOO2nofhmmiKgqKAqo6fK8oKEDMMInFDKK6QUw3\niOkm0Vj88eV8CVwqt9PKonIPC8tyWVTuYXGlF7fTOuPvOx309lNJy/bNH0L1lqYpmvSzLFiPZdFG\nYvUHE42GTuzsfmJn96MWVmFdtR2taCGK04PiyI3XC8gypmlAJIgZ9qN3N8Z/vvNHQB8/hZLi9ODY\n/vBV333PWr0JNb+C0O5vYvQ2J54wdSIHfk7szKuoxYvQCipRCypR8+eh5ORNy5eXacQwh3oxBjox\nBrowB7uInntj3EHNBYq3FPuNH7ni981klrLF5Nz5GYLP/Z+kfwfD10zwhf8ARUEtWoRl3gq0ihWo\nhVUo9pyRoQ9T8bhsrK4uZHV14UhbNGbQ3hugpSueWGjp8tPSPUR3f+iyeqtdEIrohHoDdPQGplxX\nAXIcFpx2Cw6bBaddw2m3jNwcNg2bRcVm1bBqKlaris2iYrVoWC0XHqvYhpcTbRoWTUFVFUlUiIxm\nGObIsV4kphOJGUSiOpHo8H0s/jgQjg1fOIqOXEwKhOMXlPqHIvQNRVIakpCKAo+dpZV5XLOkiHU1\nRVIwbhJmNIwZ7McI9GEG+jH9PvSuevS2UyMJ8YtVHgg4y+jPW0rEPnsn6oqiYF20Acv8NUTefiE+\ntXJs/Amv0XOecM8ThPc/geL0oOaVo3rL4re8MhRHLqahx3vS6TFMIwamiVZSfVUnHsoLXdz9rsV8\n6J01tPcGONHoo7bBx8lGH4HwxNN6jhUMx2jpitHSNXHvlZoKD/9437VzZkjRFZ119PT0UFtby+OP\nPw7Ajh07+MIXvkBvby8FBZPPp52Ndr/ZwnP7m9IdRpL8XDvlhTksKM1lYbmHRWW5FHodWZvpsiy4\ndqQugDZ/DdbVt6Q5ovSz3/gRDF8rRl/ruOeMnibCex5PbrTloDhz41c/FRUUJX5lWlEAJX6f9FhN\ndMNLOvswk5dTes4cv+6FgpGMaTNNiAQww37MsD+5iOJkFAXH9odRc/KmXvcqoOXPI+f9nyO8/wmi\ntbuTnjN8LRi+FpK+9mxOFKsTLFYUzQra8L2qTdwV09DjBxp6NHGgEYvGa1Wk8vuy2LCtugXbNe+9\nKqbn1PIqyHn/Zwm+8B8YneeSnzRNjK5zRLrOwaj6Fmg2FJszXtvClgOqinLhbxNG/V7GtsUfFyoK\nhcBaRQGXAq74X2E4ahCJ6ISjw7fhx9N1ojKODgSHbymurgOhqVa8QIELh1zxjzAl8RiF4f+G/4kU\nEv8f/SDp4fCCwkTfliMV2JUxrzVBXNn0bTtRZfnRrrQHZPJ7zcKbTPGi5qinkr7BRn1HJa8zqn3s\nV96FNczhK9ZjXmcqOcO3i64wRqr7ls2qkWO3kOOwkGO3JE5Q6kCvg6lTgrPsSjKdl/5mmMPfYeix\n+GMjhhkJJhf1TvnVFPzu+fR7FhOz5c5AvKlRLDbs63diXbqVyMFfE63bD7GJ6+uYwQH04AB626kJ\nnx/zyjjv+HsslaumN+AsoygK5YUuygtdbLu2EmO44HFtQy8nGn00dQxddBaUi6lrjRdD3bxiblwo\nvaJEQltbG6WlpWhaPNupaRolJSW0tbWlnEjIlqsNumFSkj/7hczcTitet53iPCfFeQ6Kvc7hx05s\n1rmRzbrAeeN92BauAyOGVrkWZZKiYLO1zyiahsM1OydAqqZhTvBzqe583B98lNi5N4ieeRVzsDuF\nVzMgkvpsHbNh0t+Y3RG/TUHNr8C2fida6ZKJX39Wf1eWWXmvifaJcfu+zU7OOx4gVn0t4dd/DpFU\nDhkN0MPx2xTG/d4sGngKJ1o1QbNgqbke68ptqM7xB1pz+Xelury47/wnwq/9FL357Ut4BR0i0zOF\nqQI4h2+ogH34JoSYo3Rg+CQyPHwTI5K+x1RAtYAlF0g9EWAoGn5XJUPuKgyLAwvTV2RusuO/lLbN\nLcBy88cxb/oIevMxovWHMDpOX1Gyxmh+C7Vq9WVvP90y4TxRVRVq5nmpmedl542LME2ToWB0uHZd\nkLYePx29QTp8ASIpFL0vK8yeQo1TxamYl5JSHePYsWN8+tOf5plnEldY7rjjDv793/+dVauu7myW\nEEIIIYQQQggxF13RJe3y8nI6OjrQ9Xj2Rdd1Ojs7KS8vn5bghBBCCCGEEEIIkVmuKJFQWFjIihUr\n2LVrFwC7du1ixYoVc64+ghBCCCGEEEIIIeKuaGgDQF1dHY888ggDAwN4PB4ee+wxqqvTN4+5EEII\nIYQQQgghZs4VJxKEEEIIIYQQQghx9ZhbZf+FEEIIIYQQQggxoySRIIQQQgghhBBCiJRJIkEIIYQQ\nQgghhBApk0SCEEIIIYQQQgghUiaJBCGEEEIIIYQQQqRMEgki49TX13Pvvfdy2223ce+999LQ0DBu\nHV3X+fznP88tt9zCrbfeyi9+8YvZD1SIaZbKvv/1r3+d9773vezcuZO77rqLP/7xj7MfqBAzIJX9\n/4Jz586xbt06HnvssdkLUIgZkuq+/+yzz7Jz50527NjBzp076e7unt1AhZhmqez7PT09PPTQQ+zc\nuZP3vOc9/PM//zOxWGz2gxXjSCJBZJxHH32U++67jxdeeIH77ruPz33uc+PWefrpp2lqauLFF1/k\niSee4Ktf/SrNzc1piFaI6ZPKvr927Vp++ctf8vTTT/OlL32JT33qU4RCoTREK8T0SmX/h3gi+dFH\nH+WWW26Z5QiFmBmp7Ptvv/02X/va1/if//kfdu3axU9+8hNyc3PTEK0Q0yeVff9b3/oWNTU1PP30\n0/z2t7/l+PHjvPjii2mIVowliQSRUXp6eqitrWXHjh0A7Nixg9raWnp7e5PWe/bZZ7n77rtRVZWC\nggJuueUWnn/++XSELMS0SHXfv+mmm3A6nQAsW7YM0zTp6+ub9XiFmE6p7v8A3/nOd3jXu97FwoUL\nZzlKIaZfqvv+9773PR588EGKi4sByM3NxW63z3q8QkyXVPd9RVHw+/0YhkEkEiEajVJaWpqOkMUY\nkkgQGaWtrY3S0lI0TQNA0zRKSkpoa2sbt15FRcXIcnl5Oe3t7bMaqxDTKdV9f7Tf/OY3VFVVUVZW\nNlthCjEjUt3/T548yd69e/nTP/3TNEQpxPRLdd+vq6vj/Pnz3H///XzgAx/gG9/4BqZppiNkIaZF\nqvv+X/zFX1BfX8/WrVtHbhs2bEhHyGIMSSQIIUQWOnDgAF/5ylf48pe/nO5QhJgV0WiUz372s3z+\n858fOfAU4mqh6zqnTp3i8ccf54c//CF79uzhqaeeSndYQsy4559/nmXLlrF371727NnDwYMHpRdy\nhpBEgsgo5eXldHR0oOs6EP/i7OzspLy8fNx6ra2tI8ttbW1yVVZktVT3fYDDhw/zD//wD3z961+n\nurp6tkMVYtqlsv93dXXR1NTEQw89xLZt2/j+97/Pz3/+cz772c+mK2whrliqn/0VFRXcfvvt2Gw2\n3G4327dv5+jRo+kIWYhpkeq+/6Mf/Yj3ve99qKpKbm4u27Zt4/XXX09HyGIMSSSIjFJYWMiKFSvY\ntWsXALt27WLFihUUFBQkrXf77bfzi1/8AsMw6O3t5fe//z233XZbOkIWYlqkuu8fPXqUT33qU/zn\nf/4nq1atSkeoQky7VPb/iooKXn/9dV566SVeeuklPvaxj3HPPffwhS98IV1hC3HFUv3s37FjB3v3\n7sU0TaLRKPv372f58uXpCFmIaZHqvl9ZWcmePXsAiEQi7Nu3jyVLlsx6vGI8xZQBViLD1NXV8cgj\njzAwMIDH4+Gxxx6jurqaT3ziE/zVX/0Va9asQdd1/uVf/oVXX30VgE984hPce++9aY5ciCuTyr7/\nwQ9+kJaWlqRCQ//2b//GsmXL0hi5EFculf1/tK9+9asEAgE+/elPpyliIaZHKvu+YRg89thj7Nmz\nB1VV2bp1K5/+9KdRVbkmKLJXKvt+U1MTjz76KN3d3ei6znXXXcdnPvMZLBZLusO/6kkiQQghhBBC\nCCGEECmTNKYQQgghhBBCCCFSJokEIYQQQgghhBBCpEwSCUIIIYQQQgghhEiZJBKEEEIIIYQQQgiR\nMkkkCCGEEEIIIYQQImWSSBBCCCGEEEIIIUTKJJEghBBCCCGEEEKIlEkiQQghhBBCCCGEECmTRIIQ\nQgghLsm2bdv49re/zR133MGmTZv4p3/6J8LhML29vXzyk59k48aNbN68mfvuuw/DMNIdrhBCCCGm\nmSXdAQghhBAi+zz99NN897vfxel08ud//ud84xvfwDAMSktL2bdvHwBHjhxBUZQ0RyqEEEKI6SY9\nEoQQQghxye6//37Ky8vJy8vj4Ycf5plnnsFisdDV1UVraytWq5WNGzdKIkEIIYSYgySRIIQQQohL\nVl5ePvK4oqKCzs5OPv7xj7NgwQIefPBBtm/fzne+8500RiiEEEKImSKJBCGEEEJcsra2tpHHra2t\nlJSU4Ha7eeSRR9i9ezff/OY3efzxx0eGOQghhBBi7pBEghBCCCEu2U9+8hPa29vp6+vjW9/6Fnfc\ncQcvv/wyjY2NmKZJbm4umqbJ0AYhhBBiDpJii0IIIYS4ZDt27ODBBx+ks7OT7du38/DDD/PEE0/w\nhS98gd7eXjweD3/yJ3/Cli1b0h2qEEIIIaaZYpqmme4ghBBCCJE9tm3bxr/+679yww03pDsUIYQQ\nQqSBDG0QQgghhBBCCCFEyiSRIIQQQgghhBBCiJTJ0AYhhBBCCCGEEEKkTHokCCGEEEIIIYQQImWS\nSBBCCDc2krUAAABDSURBVCGEEEIIIUTKJJEghBBCCCGEEEKIlEkiQQghhBBCCCGEECmTRIIQQggh\nhBBCCCFSJokEIYQQQgghhBBCpOz/B3aIfDdB/o16AAAAAElFTkSuQmCC\n", | |
| "text/plain": [ | |
| "<Figure size 1296x648 with 1 Axes>" | |
| ] | |
| }, | |
| "metadata": { | |
| "tags": [] | |
| } | |
| } | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "metadata": { | |
| "id": "s1HUcRGWf2JX", | |
| "colab_type": "text" | |
| }, | |
| "source": [ | |
| "## IPTW\n", | |
| "\n", | |
| "" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "SV3GdAFcftbl", | |
| "colab_type": "code", | |
| "colab": { | |
| "base_uri": "https://localhost:8080/", | |
| "height": 34 | |
| }, | |
| "outputId": "6c331618-0eb0-4fb6-d7f2-e61969676f74" | |
| }, | |
| "source": [ | |
| "w1 = df_ps_m[\"treatment\"] / df_ps_m[\"ps\"]\n", | |
| "w0 = (1 - df_ps_m[\"treatment\"]) / (1 - df_ps_m[\"ps\"]) \n", | |
| "IPW_E1 = (np.sum(df_ps_m[\"outcome\"] * w1)) / np.sum(w1) # E(y1)\n", | |
| "IPW_E0 = (np.sum(df_ps_m[\"outcome\"] * w0)) / np.sum(w0) # E(y0)\n", | |
| "\n", | |
| "IPTW = IPW_E1-IPW_E0\n", | |
| "print('caucal effect ', IPTW)" | |
| ], | |
| "execution_count": 75, | |
| "outputs": [ | |
| { | |
| "output_type": "stream", | |
| "text": [ | |
| "caucal effect 0.20887929345544462\n" | |
| ], | |
| "name": "stdout" | |
| } | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "hEPYtbo7p7ue", | |
| "colab_type": "code", | |
| "colab": { | |
| "base_uri": "https://localhost:8080/", | |
| "height": 34 | |
| }, | |
| "outputId": "3ba0399d-3119-43e5-a261-b5622913bd97" | |
| }, | |
| "source": [ | |
| "IPW_E1 = (np.sum(df_ps_m[\"outcome\"] * w1)) / len(df_ps_m) # E(y1)\n", | |
| "IPW_E0 = (np.sum(df_ps_m[\"outcome\"] * w0)) / len(df_ps_m) # E(y0)\n", | |
| "\n", | |
| "IPTW = IPW_E1-IPW_E0\n", | |
| "print('caucal effect ', IPTW)" | |
| ], | |
| "execution_count": 74, | |
| "outputs": [ | |
| { | |
| "output_type": "stream", | |
| "text": [ | |
| "caucal effect 0.5088342732661757\n" | |
| ], | |
| "name": "stdout" | |
| } | |
| ] | |
| }, | |
| { | |
| "cell_type": "markdown", | |
| "metadata": { | |
| "id": "JpGly8kmf8yG", | |
| "colab_type": "text" | |
| }, | |
| "source": [ | |
| "## AIPTW\n", | |
| "\n", | |
| "" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "Kid-sQthf7E7", | |
| "colab_type": "code", | |
| "colab": { | |
| "base_uri": "https://localhost:8080/", | |
| "height": 34 | |
| }, | |
| "outputId": "a865b5a4-e6f9-49e4-b42e-c1734701def7" | |
| }, | |
| "source": [ | |
| "delta_dr_first_iptw = np.sum(df_ps_m[\"outcome\"] * w1)\n", | |
| "delta_dr_first_augmentation = -np.sum(w1*df_ps_m[\"m\"]) +np.sum(df_ps_m[\"m\"])\n", | |
| "\n", | |
| "delta_dr_first = (delta_dr_first_iptw - delta_dr_first_augmentation) / np.sum(w1)\n", | |
| "\n", | |
| "\n", | |
| "delta_dr_second_iptw = np.sum(df_ps_m[\"outcome\"] * w0)\n", | |
| "delta_dr_second_augmentation = -np.sum(w0*df_ps_m[\"m\"]) +np.sum(df_ps_m[\"m\"])\n", | |
| "\n", | |
| "delta_dr_second = (delta_dr_second_iptw + delta_dr_second_augmentation) / np.sum(w0)\n", | |
| "\n", | |
| "delta_dr = delta_dr_first - delta_dr_second\n", | |
| "\n", | |
| "print('caucal effect ', delta_dr)" | |
| ], | |
| "execution_count": 76, | |
| "outputs": [ | |
| { | |
| "output_type": "stream", | |
| "text": [ | |
| "caucal effect 0.21218653258073766\n" | |
| ], | |
| "name": "stdout" | |
| } | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "AepbyW3hmBJr", | |
| "colab_type": "code", | |
| "colab": { | |
| "base_uri": "https://localhost:8080/", | |
| "height": 34 | |
| }, | |
| "outputId": "f8bdb3e7-d0ca-4991-8d1d-5985f389c69f" | |
| }, | |
| "source": [ | |
| "size = len(df_ps_m)\n", | |
| "size" | |
| ], | |
| "execution_count": 71, | |
| "outputs": [ | |
| { | |
| "output_type": "execute_result", | |
| "data": { | |
| "text/plain": [ | |
| "100000" | |
| ] | |
| }, | |
| "metadata": { | |
| "tags": [] | |
| }, | |
| "execution_count": 71 | |
| } | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "I6QMivyDn5LK", | |
| "colab_type": "code", | |
| "colab": { | |
| "base_uri": "https://localhost:8080/", | |
| "height": 34 | |
| }, | |
| "outputId": "bd6cde46-7093-4768-d257-c926b67a9b11" | |
| }, | |
| "source": [ | |
| "delta_dr_first_iptw = np.sum(df_ps_m[\"outcome\"] * w1)\n", | |
| "delta_dr_first_augmentation = -np.sum(w1*df_ps_m[\"m\"]) +np.sum(df_ps_m[\"m\"])\n", | |
| "\n", | |
| "delta_dr_first = (delta_dr_first_iptw - delta_dr_first_augmentation) / size\n", | |
| "\n", | |
| "\n", | |
| "delta_dr_second_iptw = np.sum(df_ps_m[\"outcome\"] * w0)\n", | |
| "delta_dr_second_augmentation = -np.sum(w0*df_ps_m[\"m\"]) +np.sum(df_ps_m[\"m\"])\n", | |
| "\n", | |
| "delta_dr_second = (delta_dr_second_iptw + delta_dr_second_augmentation) / size\n", | |
| "\n", | |
| "delta_dr = delta_dr_first - delta_dr_second\n", | |
| "\n", | |
| "print('caucal effect ', delta_dr)" | |
| ], | |
| "execution_count": 72, | |
| "outputs": [ | |
| { | |
| "output_type": "stream", | |
| "text": [ | |
| "caucal effect 0.5385921600164205\n" | |
| ], | |
| "name": "stdout" | |
| } | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "1K5fc8LaoLsE", | |
| "colab_type": "code", | |
| "colab": { | |
| "base_uri": "https://localhost:8080/", | |
| "height": 69 | |
| }, | |
| "outputId": "e9bec871-2ee1-4e59-e71c-8db811732751" | |
| }, | |
| "source": [ | |
| "print('Actual Caucal Effect ', 0.20)\n", | |
| "print('IPTW Caucal Effect ', IPTW)\n", | |
| "print('AIPTW Caucal Fffect ', delta_dr)" | |
| ], | |
| "execution_count": 80, | |
| "outputs": [ | |
| { | |
| "output_type": "stream", | |
| "text": [ | |
| "Actual Caucal Effect 0.2\n", | |
| "IPTW Caucal Effect 0.20887929345544462\n", | |
| "AIPTW Caucal Fffect 0.21218653258073766\n" | |
| ], | |
| "name": "stdout" | |
| } | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "metadata": { | |
| "id": "o2qKFWTsqO9r", | |
| "colab_type": "code", | |
| "colab": {} | |
| }, | |
| "source": [ | |
| "" | |
| ], | |
| "execution_count": 0, | |
| "outputs": [] | |
| } | |
| ] | |
| } |
丁寧な紹介記事ありがとうございます！
model_linear にはXだけではなく、Zも説明変数に入れる必要があるのではないでしょうか？
X_and_z = np.hstack([X, z[:, None]])
model_linear = LinearRegression().fit(X_and_z, y)その上で、m_1とm_0は別個に
df_ps_m['m_1'] = model_linear.predict(np.hstack([X, np.repeat([1], X.shape[0])[:, None]])).flatten()
df_ps_m['m_0'] = model_linear.predict(np.hstack([X, np.repeat([0], X.shape[0])[:, None]])).flatten()と計算し、
delta_dr_first_iptw = np.sum(df_ps_m["outcome"] * w1)
delta_dr_first_augmentation = np.sum(w1*df_ps_m["m_1"]) - np.sum(df_ps_m["m_1"])
delta_dr_first = (delta_dr_first_iptw - delta_dr_first_augmentation) / size
delta_dr_second_iptw = np.sum(df_ps_m["outcome"] * w0)
delta_dr_second_augmentation = np.sum(w0*df_ps_m["m_0"]) - np.sum(df_ps_m["m_0"])
delta_dr_second = (delta_dr_second_iptw - delta_dr_second_augmentation) / size
delta_dr = delta_dr_first - delta_dr_second
print('caucal effect ', delta_dr)とすると、 delta_dr ~ 0.203となりました。
PythonでAIPWやっている例はあまりない（Rが多い）ので助かりました。
この辺、z を説明変数に入れたり、z毎に異なる線形モデルを立てたりいろいろ流儀があるみたいでややこしいですよね。。
@tohtsky ですね！ここらへんは少しややこしいですね。僕がやってたときもrがほとんどでしたね、そもそも因果推論も全然少なかった気がします。いまはやっとCFMLが少しずつ人気が出ていますね。
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
基本的にAIPTWはIPTWよりも精度は良いはずですがなぜでしょうか、、、