Example.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "from sklearn.linear_model import LinearRegression, LogisticRegression\n",
    "\n",
    "\n",
    "def ope(a_num=500, x_dim=5, n_iter_big=15, n_iter=15, train_rounds=1000, test_rounds=100000, noise_std=0.1):\n",
    "    df_all = pd.DataFrame()\n",
    "    for _ in range(n_iter_big):\n",
    "        a = np.arange(a_num)\n",
    "        e_a = np.random.randn(a_num, x_dim + 1) # Initialize the model parameters\n",
    "\n",
    "        model_a = [LogisticRegression() for _ in a]\n",
    "        for i in a:\n",
    "            model_a[i].classes_ = np.array([0, 1])\n",
    "            model_a[i].intercept_ = e_a[i, -1:]\n",
    "            model_a[i].coef_ = e_a[i, :-1][:, np.newaxis].T\n",
    "\n",
    "        # On-policy evaluation\n",
    "        pi_e = np.ones(a_num) / a_num\n",
    "        x_test = np.random.randn(test_rounds, x_dim)\n",
    "        a_test = np.random.choice(a, test_rounds, p=pi_e)\n",
    "        r_test = np.array([model_a[a_test[i]].predict_proba([x_test[i]])[0] for i in range(test_rounds)])[:, 1]\n",
    "        value_test = r_test.mean() # ground truth policy value \n",
    "        # print(f\"Test value: {value_test}\")\n",
    "\n",
    "        pi_0 = (pi_0 := np.abs(np.random.exponential(size=a_num))) / pi_0.sum()\n",
    "        df = pd.DataFrame()\n",
    "\n",
    "        for _ in range(n_iter):\n",
    "            # Sample the logged dataset\n",
    "            x_train = np.random.randn(train_rounds, x_dim)\n",
    "            a_train = np.random.choice(a, train_rounds, p=pi_0)\n",
    "            r_train = np.array([model_a[a_train[i]].predict_proba([x_train[i]])[0] for i in range(train_rounds)])[:, 1] + (np.random.randn(train_rounds) * noise_std)\n",
    "\n",
    "            # IPS\n",
    "            acc_value = 0\n",
    "            for x_t, a_t, r_t in zip(x_train, a_train, r_train):\n",
    "                acc_value += pi_e[a_t] / pi_0[a_t] * r_t\n",
    "\n",
    "            result_ips = acc_value / len(x_train)\n",
    "            # print(f\"IPS: {result_ips}\")\n",
    "\n",
    "            # DM\n",
    "            a_u = np.unique(a_train)\n",
    "            dm_a = [LinearRegression() for _ in a]\n",
    "            for i in a_u:\n",
    "                mask = a_train == i\n",
    "                dm_a[i].fit(x_train[mask], r_train[mask])\n",
    "\n",
    "            acc_value = 0\n",
    "            for x_t, a_t, r_t in zip(x_train, a_train, r_train):\n",
    "                acc_value += sum([dm_a[i].predict([x_t]).sum() * pi_e[i] for i in a_u])\n",
    "\n",
    "            result_dm = acc_value / len(x_train)\n",
    "            # print(f\"DM: {result_dm}\")\n",
    "\n",
    "            # DR\n",
    "            acc_value = 0\n",
    "            for x_t, a_t, r_t in zip(x_train, a_train, r_train):\n",
    "                acc_value += sum([dm_a[i].predict([x_t]).sum() * pi_e[i] for i in a_u])\n",
    "                acc_value += pi_e[a_t] / pi_0[a_t] * (r_t - dm_a[a_t].predict([x_t]).sum())\n",
    "\n",
    "            result_dr = acc_value / len(x_train)\n",
    "            # print(f\"DR: {result_dr}\")\n",
    "\n",
    "            # Learned MIPS\n",
    "            w_x_a = pi_e / pi_0\n",
    "            pi_a_e = np.zeros((train_rounds, len(a)))\n",
    "            e_hat = np.zeros((a_num, x_dim + 1))\n",
    "            for i in a_u:\n",
    "                e_hat[i, :-1] = dm_a[i].coef_\n",
    "                e_hat[i, -1] = dm_a[i].intercept_\n",
    "            pi_a_e_estimator = LogisticRegression(max_iter=100000).fit(e_hat[a_train], a_train)\n",
    "            pi_a_e[:, np.unique(a_train)] = pi_a_e_estimator.predict_proba(e_hat[a_train])\n",
    "            w_x_e = (w_x_a * pi_a_e).sum(1)\n",
    "            result_learned_mips = w_x_e @ r_train / len(x_train)\n",
    "            # print(f\"Learned MIPS: {result_learned_mips}\")\n",
    "\n",
    "            # MIPS using true embeddings\n",
    "            w_x_a = pi_e / pi_0\n",
    "            pi_a_e = np.zeros((train_rounds, len(a)))\n",
    "            pi_a_e_estimator = LogisticRegression(max_iter=100000).fit(e_a[a_train], a_train)\n",
    "            pi_a_e[:, np.unique(a_train)] = pi_a_e_estimator.predict_proba(e_a[a_train])\n",
    "            w_x_e = (w_x_a * pi_a_e).sum(1)\n",
    "            result_mips = w_x_e @ r_train / len(x_train)\n",
    "            # print(f\"MIPS: {result_mips}\")\n",
    "\n",
    "            df = df.append({\n",
    "                'IPS': (result_ips - value_test) ** 2, 'DM': (result_dm - value_test) ** 2, 'DR': (result_dr - value_test) ** 2, 'Learned MIPS': (result_learned_mips - value_test) ** 2, 'MIPS': (result_mips - value_test) ** 2,\n",
    "                'n_actions': a_num, 'x_dim': x_dim, 'sample_size': train_rounds\n",
    "            }, ignore_index=True)\n",
    "\n",
    "        df_all = df_all.append(df.median(), ignore_index=True)\n",
    "    return df_all"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAVkAAADJCAYAAACJ40EsAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAABYnklEQVR4nO29eXwcR5n//66+pufSLVm2ZVu+z8S578NJIAfhSgIsxwJZEu5jWb5ssnyBxAb2F7MLYeHLtSxH2BD2IoGwywYSCCaEJJictknsOHF8y9Ytzd1X/f7omdHIlmxJlizJrvfrNZrpnu7q6lH3p5966qmnhJQShUKhUEwM2mRXQKFQKE5klMgqFArFBGJMdgVuvvlm2dLSMur9du7cSWtr6zEd+1jKGO2+o9l+pNuOx28w3ZnKv8HxrNtEHGsy77Gx7DeV7rF169Z9T0p5MwBSykl93X777XIsjHW/8SpjtPuOZvuRbjsev8F0Zyr/BsezbhNxrMm8x8ay31S6x4C1sqhx09ZdsGbNmkktY7T7jmb78Ti3k4Wp/Fsdz7pNxLEm8x4by35T9R4TcpKjC2688UbZ2trKmjVrpvQNMxXZsGGD+s0UiglkLPfYhg0buOyyy34opbwRpoDIrl27Vq5du3ZS66BQKBTjiRBinZRyLajoAoVCoZhQJl1kd+7cydq1a9mwYcNkV0WhUCiOmaKWtZaWJz2Eq7W1FeUuUCgUUwVf+jzzys94rutZrlt6E3VVraPav+jD3VlannSRVSgUiskmX+hnT+efaOzdwY/+/E9ssms4Q7N5JT6HuhXvO6ayJ11kS+4CFV2gUCiOJ+lCD/+z4z/Jv/gj/k14/J9MN91zr+XUWZdxXs1y6utX01izdNTlKnfBCNm4cSO33HILjuOwfv16brvtNhzH4corr6SjowNN0/jqV79KEAR84hOfIBqNEo1GeeaZZ3jb297GW9/61sk+BYVCcQi7+l5my3P/QKZ6Mdc+9UV+MG8NN2kmfz3zChqizdTULqO+ahHJWDOaNjZ5VO6CEXLOOeewZs0a0uk0l1xySflz6YFw6aWX8stf/pIgCNi7dy/33XcfAM8//zy//OUvJ7HmCoWiRBB4ZA5u5MH0Tpof/yTvaFrBN3p2k7YSvHD6J/h7q5ra6qXU1yymKjZrzMJ6JJTIjgHP8+js7KShoQEhBI8++ij33Xcfr3nNa1ixYgUrVqyY7CoqFCctuXwfj3f8kc5Nd/IdPeD73btItr6G3XOv4sc1y7Fa4iytXkx9zRKqYrPQdWtC6zMtRDb5rh+Pe5mpf337qPd57LHHWLt2LV1dXXz605/mnHPOAeAHP/gBX/va17jpppu4/PLLWb9+PYsXLx7vKisUiiEIAo/e7AH+/NJ/0hgU2PXiv/I/kRgL7SbeW30q+2deQ51m0Fq9iIaa5VTHJ15YK5l0kR1Jx9dYBHEiuOCCC4b0H1977bVce+21pFIp/v7v/563ve1tPPnkk8e/ggrFSUK+0E9/ehcP7vwp+9s20J87QM5M8v7G86ha8X7+UoKuaSSSrTTWraI6NhvDiByXuqmOr3HmrrvuQtd13vnOd5JMJnnzm9/M//7v/052tRSKE4og8OjPHWR39ws8tut+ZnRvYVVyAYnOJ3Fnr+Fy6xKSVi19mkZNbA6N9SupjrVgmtHjXlfV8TVCnnzySR555BEcx+HRRx8tf7733nu54YYbytstW7aML3zhC2zduhVd19m6dSvf+ta3JrHmCsWJQa7QRzp3gM6uzeR2/Td3uF18au/v+dX8q3lPrJnMrIuY3XAKs6UgnmihsW4VNfE5mGZssqs+CJUgRqFQTAn8wCWdPUhveg/b9zzA/xS6eMvL93Lf4hv4VPuf+dH8azhXi6MTIJBEYzNpqFlJbdU8LDMx2dUfRGWCGGXJKhSKSSNX6CWdO0hn9ws81bOJ8w/8kf+pWcB5L/6YwqLr6Fz8Vv6idhWvJJZyIRLbbqChbhW1yblErKrJrv6ImHSRVSO+FIqTB993SOfa6U3tprdvO3sOPMqv/BRX7H+UB2taWZJcxqtqT8G54Cvc5PQhZYBt2NQ3nk1dshU7Uj3Zp3BUVMeXQqE4ruQKvaRyB+jufZGe3m08m96D1v0c11QtYXfuIEsjSazTb+ETZjWemyLr5YjoUWY0X0hNspWYXTfZpzAqVMeXQqGYUHzfIZU7SF9qDz1928gXekhkD/BQ50ZedrpZhU5q1iWkq5YzL5KkOd8DgYuGpLnpXGqT84lHGyb7NMYNJbIKheKYyRV66U/vo7v/JbKZPXTkutjXs4nHghzvb38O+4xPsdrp5tzZl5FEUMj3kPHzmJ5Jc+PZ1Cbnk4g2ghCTfSrjjhJZhUIxakJr9UDZt+o6/XRk9uJHajn/pZ/w5qoGvtG3j9OWvgtj/l+Q0y2S9adDthPPjNHUcDq1VQtI2E0IbdLnDphQlMgOQ2UWriuvvJJcLgfApz71KTZt2sRtt93GK6+8wvbt27GsgSF6t956K3fffTef+9znuPnmmyer+grFuJPNd5PKtNHdt51Mdi+B79EmC4id9/PZ5Az+7ZWH+eXpf8OpC67j/tqlICWNhS78Qi/CiNBYdyp1NQtJ2s0nvLBWokR2GCqzcJU65h544AEuv/xyNm7cyJo1a8hms3z3u9/lQx/6EAAdHR1s3LiRWbNmKYFVTHtK1mpPajd9RWvVQ3Kw5wUe1ALW7nuMvvmvpc1KcnvjBcj5f8mrpUdXvgcvux9Ds6mvWUZd9VKqYs0ITZ/sU5oUJv1xMp3m+LrmmmuwLItf//rXANx2222sX7+eQqEAwNe//vWy4CoU05Eg8OhJ7eTF3Q/wzPP/zPYdP6Oj/Y9YQYGDex7g6zvu4ZWOjSwOPF469aO0NJ7DGUtvYp5RRTa9FyffRXX1YpbOv4HTlt/M/NmXU52YfVIJ7LQM4Vr70r+ybsePysu3L/hL1i56F+LBKwdtJ698cETbltaNhXnz5rF7924AVq1axfnnn893vvMd3vKWt6DrOo2NjWMqV6GYTNK5Drp6t9PVswnPK5DVBDt7n+eV1MuI9D7e03Au8Rnn88FILfWJFvzAx3F6SKf3oOkm1cnFNNQuDXOy6uZkn86kMi1DuNYueteQoiivfPCYth0Lu3btYu7cuezfvx+A22+/nauuuordu3dz6623smXLlnE5jkIx0RTcFD39uzjY9SxOvhMhdPKByxw3zfPbfsATTadzeWw2c5b8FRndpinwKTi9pFK70YRBdfVCGmqXURWbfVxTB043poXIThUeeugh8vk8r3rVq3jiiScAWLFiBZdccgmWZdHQcOLE9ilOTHzfoTezl46uzaTSuwCQus2uzif5rp/iH/Y+Rvair9J65m38daSaIAhwnF7S2Q6EECSTC2iYtZyaRAu6fnxSB053lMgOQ2UWri984Qtks1lc1+Xhhx/m2Wef5ZFHHiGdTnPHHXdwzz33AGHH1913301bWxv33HMP73jHOyb5LBQKkEFAKn+Qzt5t9PS8QBA4GEYMw25i9o57uVor8PedW3nzipswFrwTNI2C9HBTe0JhTcylftal1MTmYJj2ZJ/OtENl4VIoTlDyhX66+1+mvetZXKcPoZuYZi1/yu0hue1H/CHRzK3VqzhYu4xoNOxLcJwUhUIXdqSOmY3nUpOcN+VSB04HVBYuheIExfMK9KR30dG9iUxmLwIdy64ladVAoZsH/vz3FAwbZ86reU/9GfSbMaKE4uoUeohEalkw97XUVy08qSICJhIlsgrFNEcGAf3Z/XT0vEBv34sEgYdpJYjHWtjRv50d2+6n24jwIaOG5af8NXOtuvJggFBcu4lE6lkw73XUJVuVuI4zSmQVimlKNt9NV//LdHY9i+dm0HUbO9qIjiDt9NL++P/h000rWYfGmQvfTpcRZV5x35K4WnYDC1rfQF2i9aQahXU8USKrUEwjXDdHd/8O2rs3k88dBE3Htmux7Xq6+nfw3c5tfPb5H6Gdt57O1tdyz8xL0bWB27zg9OEWeonYjUpcjxOTLrIqabdCcWSCwKMvs5eOnhfo79uORGJZNSSScwBwOp5ifepFXrf/MRpb1rDl7NtpsWtpnX1FuYyC04tb6CNiN7Kw9TpqE3OVuE4Q03LEl0JxMpLOddDZu43uni34Xg7djBGNzUTTdPq8PE+9ci9NvsupdhNXGDXMv+BOVmqDR1sVnF6cfB92bAYLWy9X4nocmJYjviaDI2XhevHFFwd919PTw/79+7n77ruxbRVHqBg7BTdFd99O2rufw8l3omkWll1HNNoEQHd2P79s+y2ntD3BxpZLeFPVQjJNZ3P2IeXk8z24Tj/R+EwWL7icmrgS18lCiewwjCQLV+V3b3jDG/jpT3/K2972tsmrtGJaMtQoLCtSQyI5FwDHc3h+539D1Xyuev4ufrrsTcxtvpy/js88rKxCvhvX7ceOz6Z1zqtCcT0BE2FPJ5TIjoJrrrmGdevWlbNwVdLV1UVTU9Mk1EoxHakchdXd8zzSdzDMOPHYrLLF2eFl6Xn+O9xix/ha53a6609hz4X/wE3icIs0FNcUdnw28+ZcSU18jhLXKcL0ENmNa8e3vHPGXl5lFq6NGzfy+c9/ngceeID3vve9XHHFFUfZW3Gykyv00t3/Mh1dm3CcfnTdJBKpxyhmrpJBwB9TO7h4z8PoNUtoM2N8feaVtCx5H3OHKK+Q78J109jxFlrnXEV1vEWJ6xRjeojsMYjieFOZheucc87hs5/9LJdccgm33norN954o7rAFYcx7CisYnQAQDbfwyM9zxDs24ArDDbNey3La5Zz6axLDi9QSgqFbhwnRSw5j3lzr6E23nIcz0gxGqaHyE4RhsrCBXDppZdSXV3NT3/6U66//vpJrKFiqjBoFFbvVgIZlEdhVXZA7ep/mQW9OziIR2PnM4hlN3N6cZvDsopISb7QhetmiuL6Gmrjs4/reSlGjxLZYRhJFi7HcXjooYd49atfzbp163j3u99NoVBQnV8nMZlcJ519L9HdsxnPy6BrNnZsxqABAU7g8kjfVha9cBd/27CUH0qb+Uvfxfy51wxdqJTk8514XpZocj7z555DTXzWcTojxbGiRHYYzjrrLB5++OEhv6urqzvsu/POO49t27Ydj6opphiOm6Gnf+eQo7Aq6UnvJhFpZMZj/4cH5l7K+2Zfzo/nXImpDTOTgJQU8p24XpZY1QLmN55DzRARBYqpjRJZhWIMlEZhtfc8T6rvJSQBllVbHoVVQkrJKwf+wH35/Vy79xGSS99J5Lw7+HKkevjCi5ZrKK4LaW08W4nrNEaJrEIxClLZg3T2baenZwu+Vxg0CquSTODS2PU8j/Q8x+bcfi6sWcGc879MwogQDFe4lOTzHXhenmjVQhY0nUN1bMaEn5NiYlEiq1AchYKboqtvBx3dz+HkuytGYR0+r9XO9B6eaH+M5fv/QHrF+zmn4SzOrF1+xKgTGfgU8l14foFY9WIWNJ5FdUzFXJ8oTKjICiFeDywDTOBFKeV/TeTxFIrxojQKq71zE6n0ToTQBo3CqiQIPJLdf+bbbidvf+FH/HLVe5l72iepTcwlAIaTVxn4YYeWXyBevYSFTWdRFVXieqIxrMgKIZqBLwCrpZRnV6x/FXA90A5IKeW6I5T/lJTy50KIauB7gBJZxZRFBgH9uQN09mylp29reRRWIt4y5Lj/Di9D99Yf8OlIlA2pTl674mb8NZfwgeE6ssrHCcXV9wvEqpeyaMZZVNlqKvkTlSNZshcB9wOnlVYIIWLAt4GVUsqCEOJeIcQVwEFCQa7kfVLKfcXP1wFfGrdaKxTjSK7QS1ffdjq7tww5CmsQUiICjx1bvsq/CJ8rjSq+NOMyOhcvoPYoxxmwXB0S1UtoUeJ6UjCsyEopfyKEWHPI6vOBXVLKQnH5D8C1UspPAG8cqhwhxLXADmDfUN+X8skCUyanbGUGrjvvvJPzzjtvsqtU5tFHH+VjH/sYd95552G/1ec+9zluv/12vvWtb/GBD3ygvL6zs5PW1lZe//rX87nPfY6nn36aW265hQsuuID169fzxBNPlJcXLVrESy+9xBVXXMFNN91EW1sbf/u3f8uCBQvwfZ+nnnqKf/zHf+SUU045zmc+vpRGYbV3P0c2vT90BxwyCquSnJPi0a6n0doe4Y1Vy4m2XMVn4i3E7KNJK8jAo5DrwJUuyZrltDSeRfKQ8C7F9GfDhg2lXLJwDPlkm4BUxXJ/cd2QCCHeCNwKPAckgcPmyJ6K+WQrM3BNJYEFuOiiizj11FOH/O62227j5z//OXfccQd/9Vd/RSQSAeAb3/gG8+bN4+1vfzuLFi1i0aJFfPOb3+Ttb387c+fOZe7cueXl1772teTzeerq6njLW97Cl7/8Zc455xw+9rGPAXDXXXfR19d33M53PJFBQF92Lx09W+nr3TYwCis+e9g0gHqhh8d2/pwH8m28wXNxl7yD/ngrzSNIGxiKazuu9EnWrmBRwxlKXE9gKo3EdevW7SytH63IthOKZYmq4rohkVL+DPjZKI8xpdm3bx+f+tSnWLVqFS+99BLvf//7OfPMM/nZz37G/fffz9KlS9m8eTPf+ta3qKqq4i/+4i/YsWMHr371q/n973/P9ddfz69+9St0XefUU0/liSee4O1vfzvvfe97gVAoPc9D13WSySS33HILAB/72MdwXZcFCxawd+/eYeu3atUqstks3/3ud/nwhz9Me3s7QggaG0feLO3s7CQajRKJRGhubuZHP/oRZ599Nueddx433njjMf1+k8FIRmFV4gUef+jfzpJtdzNr9d+yKvCZs+LDzDKTQ25/KDLwwlCswCdZt4JF9WeQtOvG85QU04jRiuzjwDwhRKToMrgQ+OaxVGBE0898cy18q9i/9sHb4ew18J7L4Pu/hT9tGNt3Zw9zrKPwyU9+kte97nW8/e1vZ+fOnVx33XU888wz1NbW8k//9E9UV1dz5513cvfdd/PhD3+YL37xi1x44YWsW7eOfD5PW1sbp59+Op/61Ke444476Ojo4PLLL+e9730vv/rVr3jiiSd48MEHgfDJeOWVV7Jv3z62b9/OAw88AMD9999/xDrefvvtXHPNNdx000184xvf4KMf/eiwo9cq+fGPf8zGjRvp6urigQcewLIs/uZv/oZIJMJHPvIR9u7dy1vf+la+8IUvkEyOTHAmi4FRWM+Rz7WDZhC16w8bhVWJXujl++mXeMdz3+S+1iv4QNM5NFkJ4svfQ3wEx5R+UVylT7JuJXMaziQRqRm3c1JMD0Y8/YwQ4lLgncBMIcRngC9LKbNCiA8CXxNCdACbpJS/OZYKjchd8KG14auSzcX0GWevGdt3Y2TTpk00NTWxe/dupJQ0NTURBAGJRILPfe5zNDQ08PTTT7Ny5cryPosWLcI0TUzTJJlMsn//fpYsWQJAY2MjqVSqXHY2m2X9+vUAzJkzh46ODv785z+zePHicnkLFiw4Yh1XrlzJBRdcwLp164hGozQ0NIzo3Erugkp0XeejH/0oH/3oR9m7dy833XQT69at40tfmnr9mEHg0ZveS3vPFtL9OypGYQ2VJDBESsmu9j/ylJ/iYwefY2njKew689N8tXrhiI8bims7ngyoqjuVlobTlLiexIx4+hkp5e+A3w2x/iHgofGq0HSZSHHz5s20t7ezevVqrrjiCl7/+tcjpWT27NlomsbNN9/MV7/6VS655BK+853vsH///vK+QwWiD7Vu9erVPP744/zd3/0dAA8//DCLFi2iUCgMskR37Nhx1PrefvvtXHTRRbz44otjOd0yN954I5///OeZM2cOLS0tXHnllezateuYyhxvRjoKqxIpJU/u+QXP9W5FDwosrV7O3tWf4LxRTNESiutBPCmL4no6iSMNl1WcFKiJFEfAoRm4AF555RUuvvhivvSlL3HbbbexefNmDhw4wGWXXQbATTfdxOc//3kuu+wynnrqKXp6enjppZe466672LVrF9///vd5z3veQ6FQ4O6772bTpk08+eSTbNmyhb6+Pu69915uuOEGNm7cyKc+9SkMwyCfz7N+/XrmzJnDAw88wM0338ycOXOQUnL33Xdz5plnDmq233nnnWzatInf/OY3XHHFFXR1dQFw7733smvXLv7jP/6D1atXs3HjxmGXm5ubOeuss8plXnzxxXz84x9n1apVFAoFdu3axVe+8pXj+N8YGt936Op/hQMdf6KQ70DTIsOOwqpkb3Y/7b1bubZjCy/XLmBNw1ksmHkp2ijyAEvfLVquUF1/Ki31pxFX4qoocqglK6Q8tubzsbJ27Vo51URWMXXJF/rp6NtGR+dT+F4eK1KLZR3ZPxwEPtv3/5a7vC7uePkB/uWU9/OhxELcUQ5dDXyXQllcT6Ol/lQlroohEUKsk1KuhSlgyU4Xd4FicunPtnGwazO9vVsRCOxoA/pRhqB2eznsHT8lNu8aGnY/wKXL3k3h0jfwXsPGHcWxpe+EWbGkpKb+NFoaVhO3qo7thBQnLMpdoJg2+L5DT3o3be1Pks+1oetRYkfztQYBm/MH8F68m//S4HqzlvNlQP0FX+bqUR5f+g75XCeugJr605hTv5rYUaxmhWLEHV8KxWRRcFN09GylvfNpfC+HFaklkZx3lH2y9Gf2Mb9rM/t6NkHz+Xy6/gwS0Ua8UR4/FNcOXCGobTqTlrpTiJmJsZ+Q4qRGiaxiyjC0S+DIgygOZPayde9D/K7Qzrt9yJzyMa5Z8EYYYtrsoyF9h1yuA0/TqG06S4mrYlyYdJFVPtmTmyDw6E7tHJVLAAAvT3z7j/mI7nJHrp+/WvEB5lrVh08+OAJKlqunGdQ0nc2culOImiMZfqBQHI7yySqmBI6bob3neTo6n8HzspiR6qO6BKCYXnDz/+OXdcu4Pd7Cv828GNOMjqkO0iuQy7Xj6xa1TefSUreKqBkbU1kKRQnlkx0BjzzyCLfddhuvvPIK27dvx7IGYi9vvfVW7r77bj7+8Y/z/PPP8+CDD7J+/XquvvpqbrnlFh588EHe9773kc/n2b59O9/61rdoamri29/+Ns899xwzZsxgx44dzJo1qzyy62QilT3Iwe7N9Pa+AFISidZjR48+Is0JfJ7d9n1+nz/I0uQC/qrxAnpjTRw5c+vQhOLagaeb1M24gJa6lUpcFROGEtkhuOSSS1izZk050cqHPvQhADo6Oti4cSOzZs3illtuYcOGDWzZsoV3vetdQDg6asuWLWXL/Pbbb+dLX/oSn/nMZ/jsZz9bTtbieR4f+chHJuv0jjtB4IVRAh1PksvsR9MtotHmo7sEgIzv0Nn+R87p3sovaubynparqD3CMNkj1sPLh24B3aKu+XxaapW4Kiae0fcOjDMln2xFHsYpw2233cb69espFML0uV//+tfLgjsS2traaGhoIBKJIKXkzjvvpLOzE8Mw+Pa3vz1R1Z4yOG6G/R3P8Ny2H7Bj589x3RSJ5FxisaMLbMF3uK/9ccQjH+Tfgizdrddy3eyxCWzg5cmm95B2UlTPvIjVS/+KxU1nK4FVTAjT0if7yL1/4vc/fWrQuouvO3PQ8mi+v/i6M7nkhrM5GqtWreL888/nO9/5Dm95y1vQdf2oKQP379/P2rVryWaztLS08PGPfxzLsvjd737H+vXrWb9+PUuXLuUzn/kMV1892sjN6UE618HBrk309PwZCdjRI2e/GkTgI6WP9vuPsGXeqzjr9L/l72qWjWrwQLkoLx9mxdJt6psvoqV2JbZhj6EkhWLkTEuf7CU3nH1UUTzW74fj9ttv56qrrmL37t3ceuutbNmy5Yjbz5o1a8iHxsqVK7n77rvxfZ/77ruPG264gd27d1Nff2IkcS5lwDrQ+SSZ9F6EbmLHmofN2Xr4/j4b+v7MrO3/wZktV9N2zue5bYzTYQdujny+C8+I0Nh8KbNqlylxVUwa00JkJ5MVK1ZwySWXYFnWiFMGHsrOnTv5/Oc/z/e+9z10Xee6667j/e9//zjXdHJw3RwdfS9ysPNPeE4Kw6w6YmrBw5CSp7uf47m9D7Kq0Ev/8ptIx+dTNYpsWCVCce3EM2M0zl7D7JolRHQlrorJRYnsEJSycKXTae644w7uueceIOz4uvvuu2lra+OHP/whjzzyCG1tbdxzzz1ceeWV5e/uvPNOPvGJT5TLq66upquri7/5m7+hurqaV155hS9+8YvT2orN5Do52LOF7q5NSCQRux47Obrs/21dzzHXbuKcrf9Kz6kf4qz4gmGngTkSgZshn+/GM+M0zr6M2TVLieiRUZejUEwEKguXYsTIIKA3s5u2jqfIpPcgdAvbbhixS6DEvo4n+Up2F3+38yHsc79AIlIPo0g1WK6P75DLtoeWa9O5RctViati8lFZuBSjInQJbOdg55N4Th+GWUU83jJqq9NK7+W/ep4jse8R3jD7cqou+QaWOHoY12HIACffST5wqWk+n9b6U5W4KqYM0zK6QDE5DLgENodTudh1IxqVdSidqV08tv/XXBtfwLVYuOd8gag2lmEEELgpsvlu9KpFLGu+iOoRTMmtUBxPpmV0geL4MeASeIZMeheaZhGJNWGMQRS7Ct3s3vZD7rAsPu/65BZdSEK3xnTRha6Bg/hWFbPnvZ6ZVfPRxpAERqE43iiRVQDgeQU6+7ZzoONPuE4vhpkYk0sAoNfP07rjfh6qXcBM4J8XvJPasSa5lhI330kucKhpOofWxtNVxIBiWqFE9iQnm+/mYM/zdHU9h5Qull0/JpcAQCZw2bb9R9yfb+cD1cu5OrkE/9SxxSdDGDWQzXeiJ1tZOvMSakY6oEGhmEIokT0JkUFAX3YvBzqeIZXaidB0InYjhj42PynA7rZH2LTnl+yrX8WH572Rppql+GOtn+9RyB/E0W1mz72WWdWLlGtAMW1RInsS4XkFuvpe4kDnkziFHnQzTjw+e0wuAQAv8PlN11O86eWfU3fqRzjDquK19aeNvYJS4ha6yXtZkg1nsLzhLOwxpjFUKKYKky6yKoRr4skVemnvfp6O7meRvlOMEhhbJisIh8A+mH6Zi5/+Er+ZdxmrFt7A7FgziVjzmMv0vRy57EH0ZAuLW99A7VEmSVQopioqhOskQQYB/dn9tHU9S6r/JRAG9jG6BGQQUPCz1D7/XXbpglmnfpgv1pwyZks4LNOjkG3H0U1mzb2aWTVL0ccSO6tQTBFUCNcJjuc5dPe/TFvnn3Dy3ehmlHhsbFEClWzue5FtO+/nrMAnvvLDvM+IH3OZbr6bnJsmWb+aZU3nqtSDihMSJbInCLlCL+09W+noeprALxCJHJtLoFxueg9fOfg7Prj39+xe/dfMqVqMr5mMfhDsAL6XI5drR4vNZNHca6mPzzzmeioUUxUlstMYGQT05do40PUsqb7tCGEQsRswdOvoOx+FAz3Ps00I3v3n73HG8ncRb30L1x5ruTIgnz2Ai05zy6uZXbNs1HkPFIrphrrCpyGeV6C7fwdtnU9RyHdgmLFxcQkAdPXv4J6eZ7lw3yMw9xp2XvCPXD2G5C2H1bnQQ9ZJEa9dydIZ5xGzksdcpkIxHVAiO43I5rto73mBru7n8H2XSKSW5BgHDhyKyPdgpnfz33t+wdmxWSw594skxyHpSuA75LJtEG1g0YI3UZ9oGYfaKhTTByWyUxwZ+PRmdnOw42lS6T2gmdh2wzFFCVTSX+jlhR3/xQOGyf8jxjtOuxVzPHr3ZUAuexAPaJp5GS11qzB0dbkpTj7UVT9Fcd0MnT3baO96GsdNo5sJ4ok5iHFoupfY/udv8glD58vZTj6w6sN0RWrGNMX2obhOH/lCL9GaZSydcQHxSPU4lKpQTE8mXWTVYIQKpCSdPUB712a6+7YiAStaR8I+fMaBvBuQK3gITWAZGpauYehHF2An8Gnf+TMW165E6DH+edZlzBonl0PgO+SzB5GRaubPv56GxNxxfSgoFNMBNRhhChL4Dj39r3Cg809kc50IPYIdm4leMW22lFDwfNJ5l+60S8H1ESJcjwAkmLogFjGJRQwihoZlhuIrBHgy4Jk9v2J7x0bSyRbeb0RZtOzG8TkBKcnl2nGlR2PzhcytXz1u7gyFYrqhBiNMIQqFbjq7n6ejezNu4KBb1YNiW6WEnFsU1oyD6wUIwDZ0ktHDRcwPAnKOT3/OCfcn1N9aepF7vsW/N53KB2e8itZZFyC18bEwPTdFLtdFtGYxS5ovJhGpGZdyFYoTBSWyxxsZ0J/eQ3vH0/Skd4FmEonWk9DCGNSSsKayDj0ZB9eXCAG2qWPbR7YOdU1D1yBiasgg4ClvH6966Wvsbf0/VNW+nneby8lJjRf29WEaGnHLIGYZREwNy9AxdTHiqbYC3yWXOwhmktb519OYnKdcAwrFECiRPU74boau3u0c7HyKvJtGmDFiiTloQhBISBd80jmH7nQBX4IGREwd2xp97OsLbhtzdv4L2xPN1M5+MwvsJoJYM5Vpsz0/IJP36MsOWL0aEI3oxC2TaETH1EOXg14pnlJSyHfi+AUaGs9hTuMZmOMw+EGhOFFRIjuRFH2VnV2b6eh7AT+QmNF6EnYtvpRkCz59WZe+TCisuhBELB19jIMKdhXauKDvecz8Dh6f+TreEjtl2BFVhq5h6ACVfl+J60u6UgX8lCz7GyxDJ27p2HoBGfSQqJrP8gVrSKok2grFUVEiOwFIz6EvvZP29j/Rm29H6DZ2tBkTjVzBo6MvS2/WQUrQtWMTVoD2/F4eSf+Rdx78A39afAuzjUs5awzDVYUQWEYYrVCJ5zr09e6jQ0SIJi4gn51H30sO1fEuauMWCdskaulEj/E8FIoTESWy44hb6KWr53nauzZR8POISA3R6Fyyrsf+7gL9OQcJGEIQixhox+jDtAud/Dpo553bv87ji95H37JrmKXHkONzOgAETjeBnyNZfzq1VadhauEoMD+Q5Ao+vek0fiCRSDQhiNsGNbEINXGTqGUQi+hYhkpdqDh5USJ7rAQ+6cw+OjufoSu1k0DTEWYdvqilp88hle8rhldpxCPGMXcOBVLS1fcEv5Jd3N75Zy6ZeT1bVt3JG7TxnVww8HP4TifYs2lofA0xa7BrQNfCB0UsMnAJSSlxvID2vhz7ujPl8DLL0KiOWdTGTeJ2KL5RU0cbpwgHhWIqM21FtmNfD4HnY0YMLNtCNzR0Qy+/TzSBm6G37yUOdj5N2unDFVEC0Uh/xiOdzwNgGhrJiMGIu+yPhAzo6XqIR/PbqA4kZ8aX8eyCvyap28SPvfSKw/j4TjuBMEjWXU51bCGaNrLfUwhBxNSJmIO3L3WydaUKBDK0szUBcdukNm5RE7OwLZ2YZWAayt2gOLGYtiKb6c2G76k8MpBh7FNRy4QQmJaBZZuYEQMzYmKYelmENV0bm0UpJYV8O12dm+no20racynIKrKFWrIFD0Qey9BJ2uMjrH7g8HJ6MzuzW/hQNs3TTZdzYWQWtYlVExIu5Ts9BH4aM3EKdVVnYBnjM79W2MmmDWn1tvXk2N2RQSIRCCKmRk08Qm3cIhYxyr5eFR6mmK5MW5EFMCMGmn645SMDie8H5LMOmVSewA+AwbpnRUzMiFm0hE10Qy8K8RAi7Dv09++kq+sZ2lL7SXs6GSeO68cQgGlAMjp+YUy6n6Pv4E+50yiwtvsl9s97Jy80zKHKGE+bdYDAz+M7HRCZQV3DVcQjjRNynEqOZPX2ZR3a+/KADEPLhCAZNamNm1THLGzLIGbpGEP87xWKqcb0Fdmug6D5YJigG6DroYoKDaFrGEIDoYGugaGF7VOhgaYhEQRBQCHnkEuHIlwaHQXhB9MyMbQc6fxOOvq30eem6fcjSJFENwS2aWDb43eTu36WP3j7uHbHD+hd9Enm+D7vaXojQVMta8btKIMJXQMdBEInWbuG6vjiEbsGJgpD10joGokKF3MQhFbv3q4cOzsy5f+VbenUxixiloGug6Fp6LoW/quFQBMCTRNogsM/C4GmhZ+VlayYSKavyLbvB02G4kkQxnRW9quXB/ZTVM+BZQHomo6u66DpAyJdFOq828uB7l10pNvpdcDTEhiahWloCN0DAa7p4lkaWkQgTA3d0MDUELpA6Iz4xk1md/Pv/j7euvs/2DvjQl5ovobZeoKeOe9kIjOvem4vgZfCiq+gvvosTGPqzq+laQLb0rGtwQ8A1wvoTjt0BKGvt+TvBUCGz9TS/33YiAsZXiq6JoovDV0XGMXPhh6uNzRRFPGRibmAQUJeWq84+ZhQkRVCrAbOBqqAGinlbeN6gFgcjDGeggwgkOG7DAiyWXr799Pdu4PeXApPmpiaTaOmo2s94S6EAl1yASMFvhTIAMpZWihZz6AZoFkawhLopoY0NISpkzZ82pzHeUFL86FClItjjext+VtuNGtAaEgnQGoOCIEUAqlpRSt9HPy8voPvHkQzG6ifcTnxyPSdX8s0tHHpKJNSFv+nA++OFyBlgJRhRMeh24AcSM5TfJODmkOHHmR4Mdc0DUMrircmBoRdD5fLwi0GC/ZQVrkS86nHsAolhGgGvgCsllKeXbH+VcD1QDsgpZTrhitDSvmcECIFfBL46bjVejwQGoEWkM6k6OrZSVf/PhwZII0o0eQsosULNSi+hi2GQ+4rKcsdcX4AXlYSpCVC+kRkiscbnqRG38VWsYBz07XsKSxkpg6a1oOud4MQaCKgdJ9IAaJohgWaAF0n0DSkpiN1AykoukEEElFcFuWXFBpShPLvud1oQlBTfTpVkcWINJA5UBTvyn0qzqos7sX1pVaBGPheHrL/wLalX0YMqtNUQwiBXjqvo1GylqUc+MzAZxE+fYfeJgiQgSTwJdIL52gLJHhBQKFkiQdBWdRlEL5XHkOU61DxWcjwyR8Ujy0DdCEwhERHoolwNKGhhZ8NEY7z0zXQRSgCuhYOq9YIyxOAJmVZ2EXJOgcExfWAoELUB53vkD/ekX/P4b4bbr/hjiWHKY8jHGu476rrYeWZQx9nhBzJDLwIuB84rbRCCBEDvg2slFIWhBD3CiGuAA4SCnIl75NStkspdwghbgF+CDx0TLUdBwICMrkCPX0H6O19hYzXj9RNrEiS+HhM6icEopjXVQApkaIj8hRbLY1P9+zmWXEZc7JnsyioRUrI6qJkTIM3UIymS3QdNCPA0MNlTQQIEYQ3jpQI1wFkUYTlgMtEhjdC6boJvAxBkCEWmU1VbAlGj4Ce7UVBKNV0gEP906V1JWuMit3KQR2lL2RFGWLwjrL8uegjp+gnF1qYFUxoSKEVn1zF27nCipdFn3pZ2IU2WLSD8Icsi08QFIWpeEMWjy+CoFgxWW6WiKI4idJ2MiiL2cANXSF0lXkmKzRm0C9ZWigdWoiB/4sY+PHKP0v5B2bQNsMbyAOCMPCgKy4HEBR/81BzBIGU+AIcSha5KP9fgtL/KCwMURTO4SSz9K+VxWNqiKK3TRQ/CzShDbKw9ZL1rQ1Y4rrQ0PQB9wqaKAu3QBTdMaVLICy79Ln8ADjSg/uID/Wj7Od7kOk/wv4jY1hVkVL+RAix5pDV5wO7pJSF4vIfgGullJ8A3nh4PcVVUspfSSnTQoghZ84rJe0GJixxty8DMgWP/v4+uvv34eT34QkPzYxiR+uPeeTVYcgADw/N+lf+T+MCvt6xh3z+DTyjXcwp+eKxSq1cfejLWAahRvh5HeeQh3Z4k0scIYhYPrYBuiEROmhCohX7AAPpErg9YCWoqr6IqBVGDfjje7ajZ5BolVdWiB3hySMR/sCJa5WiN6gcDreSK6xoOcz68j6aIKi03Ct9AIdZ9NOLypbWRHZpykCWHltlt4onQSLD5xVF673Cdx4+APwBa72o3OXnTPHzYXdIsfzKh7yGQNPkgEALDR0QmgjviQrhDt0sINDCe6X8O4nwYSPD8gLPIyocRhprs2HDhlLCbjiGpN1NQKpiub+4bjgahRD/l7DFfddQG4wlabeUkl8/2kFVTYSmpjgNDVFqayKDfFFeEJB1PHrTefr7O/GdfbhBN4Guo9lxotr4J5VOa/t5IJHilv4HOCj/kr2s5Kv7T8fCZgWMqDVaQmhhv15/IGjLm7TlTNpyRvHd5EDeQACuDC0RS5NYmiSiSUwtIKL5WJpPzFpMzIoSMzVsoxfbEERNgW1q2LogYoCti+JngW2En21DENEHL1s64/NAOoJgjbKxqZgApJTF55cc1EAqv5U6Eg9rYcvDWv9lt0dxe1kU3lJftQxkudzScUsiOuAtKPnBGbSubE4z0FApNS5K5QYV24XP7mLLoFiPQBTdcWKgBVd6tkrPJ6b5XDbC363SSFy3bt3O0vrRimw7UGmRVhXXDYmU8kdHK3Cs08+cvrKaniwcPJhly5+7yOd9brhhIU8+3c6shQl6ujIYejfC349LDl+LYJi1mGL8YysjvMD3ag9wWf5FZgQrecx4OzO9BurchhHtn/EEB3Im+4viuT9nlsVUCMmsqEez7TIr6nJeQ5aZUZdm2yVmhBeSF4ATCJxAkHc9cm4OT6tHGC040qLgBRQ8SSGQFIKAQh76s9AlwZHgSihIcHzCbXxJwYe8L8l74XLekzg+WDoDglwhzBF9QIzDZcoCHRm0fkC4I2UB55DlsENIcewEQcm/W7QigwEB8/2AIJAEniQIJL4vi9tT4QqRA+ojxYAYVYpcpblZuX3RqhjUeCjXTAwsD9FwoNiZV/qmcmal4xVy5zkuTjo76v2OdfqZx4F5QohI0WVwIfDNUdeigrFYskII6mstmufE8IRGpuDR1Zfjhb29OH6O/o79PP+MRSAF0WQ9iWRAstojUusdvfARstfopiPyNIucLA2FxVyaqaXRvZnZcuhGWcEXHMgPWKJtOYP9eZMDOZO8L5gZdZkZ9ZgZdTmtJsfVM1PMirokzSN1u4WEYcAekaCfpB0h3rCciHn4vGCVlG+4IGy+yXKzfOBO0HWBYWrohsAwQv+aJwUO4ARQCCDvDRbivF/6THk540q68sHh2xXf8x7l9aV3YLAw64eId4UgRw1BzBTEDK34LireB9bFzXBbfRoKuCyKX1DuEBv4HwZ+gB+EAhn4xVcQfheacKFYSsSAOJZ9mxX+TU2gGypu+FgpGos7S8tHii64FHgnMFMI8Rngy1LKrBDig8DXhBAdwCYp5W8mtMbD0J936M345AKJJj2sIEU82Af1fbiawfJzo0jXJJMyyKYM+rpNNA0O7rVpWZClt8sknvSJJTxGGn+vS4eCvoOfVL/CmzLbeUa7jEV+PZ3aLOY44AawL1/RrM+b7M8ZHMiZ9HsaM2yPmXYopEuqClw6I81M26PW8o/J3ed7/UjpY8XmE7fnjGhAgRBhR8UgE6GCkrXjeQGuI0su0opoB0lMFyQNLRTimEDXdTRdDHRsHIOYecGAGA8Scj8U70LFd1lPknMlXXmfPSlJ1h1Yl/ECsq4k54Vin/MkpsaAIB9BlI+8rrjeEFgjmMDy0N81fLgxSDB9P0D6ofXplwQzCEIhDWQofmVH5YClWRJJUYz+ECIUS0McP6tPMTxH6vj6HfC7IdY/xDhGCYzVXdC3vw2SCeJWmkKhjZQsFF0CNVilCysisSIutQ1use4wd3EG3xMUcjrd7RFyWR075tPQXKCm3iEIBHZ0wHqUSF6yumjiQV7Uz+Wd6U083HcFm7pei5GL8V85k7aihdpdMKiPeMXmvMfcmMO59RlmRj3qI95wejZmgqBA4KUxrAbiscUY4zigoNK6Gf74ocWUz/vInCzGCw90Yggx2BIOreGiAOtHFmFDEyQsQWLczihEytAVknGDsvCWRDl7yLqUIzmYCch5bnndgFgH4bqi1R01iuJrCKK6IFp0gUSL7hVbo/yKGuH6qK6F22nF5aLo23pFXKyuoRtKLKcTJ8xstaf/fB3RVDdS03jwxk+w/OmNGG6BF9a8noUbH8Y3LJxYnLbFp2Cn+wl0g0KiCtOSmJZk7uLQ1xIEkE2Hll82ZZDP6ei6pOOgQVciz9m1D/D9ObVc29PMfS+v5it9l1FlBsyyXZqjLrOiHqtrc8yMujRGPMzjMJze9zPgF0CPkqhaTcQame93vNGK3bTD2c1hc1biFCSFnF/ujCj3ImthqI9hhCPmDFOrEGCKQjO+4iJE0Qdckamt7DKp8FeWfJOH+i1Lvsuy31II3CAg50Heh5wfkPMFeV+SC4p+7QByRSu805Pk8qFYZ31/wBL3IeeF690AbIOiWIeiHTPCesdK60pCXXKXlNcP7FcSfNtg8BRCJzBBsZPMlxTfZcXn0vowhj0A/CDcZvD34TrX9Zgxhu7WEbsLpjq/fu270ONRogg0K0LbstMwCmGKwUA3sDN9VHXsp23xqSx46hGsXJpnXvMOLv/u/4dn2TjROA+9/v00PbsRxzO5r+Vylr/yFH5tG1tmwBsaA3Z3Xchj6ddy3dY4NHq8v8qjYWYbiZiPrksM8/j1dwfSR/oZZOCim3VE48uwjLpyPONUJHRJFN0SQ1DyKbpuQMEpirCEI/mF9Up3hD5YhEsdPEFFk7zSbxkE4Yy+gV/czg/wy9Z32ASXoui3HKIpXno3rcHHjaCPq8Xty9BPHYpvKLzh+yHrPEl3QZLNhP7ubHk9g/bLe2BqA9a2rVMU7WHEuejrhqIgDRKpoQWpvC4IY25D8Rri+2HEr3ScIdeVBfHQfQ4vRxJ6wMIY3NJnUfE5fOAM/j5cp2sD+2lCoMmAa8chV5IYGFEyOdx4442ytbV11O6Cn63/JEYygWEeOftVqee+LW/SltXpS/k4qRx+OstTicW8LbORmUaW7Nm9XLPxt9SkPapzOo++5WNcfu+/4lkRnr7mHZz10++zN7EAPaKzcf7VdOy1iOgudp2gTnRT0+ghq+yxD/MdhkB6SC+NRGJGZhC1WzCNqqPveAIwyKoMDvELF2NydD3svSlZmUP6LZHhiCVBOCBgkGie+E3xoOgiqRTnbNG3XSnaWS+0vLNFfzccQZDEYEEa9vviSDJdE4eIn6jY//B1leWW12miWBbl96H2GS9K0QVXfOhto9pvw4YNXHbZZT+UUt4IU8CSHau7oJKheu5LsaWVPffNtsuiRp+Zc02aoiZXxR/iZXsn1/X28Ix5ATtefR213kIAEsDG629G8z3QBC9deCVWLoPu5Jk9P88FXQ/RJ6rZkTgNe/N+lv/p9zxS/xrsC5upe/wZZMwic+piTnvy57jROPlEFXtOPY/qA3twojFyVXWIwEfqQ/8LQn9rFiE0LHsedqQZY5zyu04XRuoXBlSv+BHQRNFSNQREJrs2Jz4nhLvgie0d/KBtDgf3xjlQsOh3NZoiXhhPWuq5b8owM+oe1nOvywJ91sP8Q53GF7pfZr/3erZq9VT5QwiYEARGOGihd9a8QV/tOutiAOoIoGUBz7OAOk+C5tOxeDl9qSj9W2O8YL2LuqCH8zOP099jsGLTRtxEkt2rz+PSu/4R37BI1zXxxze/n+Ubfk5/XR37Fixi3ovbkNVzkbXzSdfPx8j04QsDqY//IIrpjEqGMgyhvwTNc7FSXeiFLH4kipNsoG7bY6RnLUHzXOqffwTdyZGdMZ+eReew6Odfou3sNyACn9mP/SdetIqexefQuXINM57+BT1LzsOLxIn0HcSLVuMkagkiUzeD21RgWoqsqWu0RPKcW11gdhIajtJz75Lj97EeXle4j57gncxxJe/vvJDAv4o17vjVSyv6saoXmFTjAf14riCTirKn5mI6XrDZfMUNvPB0NdGdPrte9w9U2VmSRh++18/BltkE0VqSsWXU5Xdi9rxIsOcV0q2nseLHnybVspz959/A8n/7LG6shmxTK69c82Hm/PaHZJoX0t+6moYtv8W3ojiJOvoWnkm87SXceA1OohbdyeNbNiOOWVMcX6RE8xz0Qtgp60UTJPduJV8zAyEDal5+Gt3Jkq9ppm/+6bQ+9C90rroUKTTmbvhX9EKWvtZT2XfhWznza+9ix9Ufwo3XMO/h7+NbUXoWn0fb2a+jZsdTOFWNeHaCwIzgJmrJ1zTj2XH2n/cmso3zyuvNTB9uvBoReBi5NJqTJ95zgJZH/x0z10/34nM5eMZrWHn3Lbx87V+DEDQ99xBurJq++afT13oqNS8/Sf/cUwgMC90t4MaqkcbJYyyc0D7ZpNzPv1X18pbcA3w3cSYX5OqZk1+GLibvHywlOHmNTMogndLI9ms4BYNV5+Rp21XHslNmkO53SFRZQ1ppwveI9LRhZvsQgU9/62qanvkVhZomsg3zaPnDf6A5OXw7wc4r38fye/4v/a2r6V5yPqvu+gS6kyczcxGb3/NPrPrB39C99AK6VlxM64PfwY9EydW1sO/itzHzjz8lM2MBmeaFVO98Dt+K4saqyTYvxEx14UfiBNb4Tt44bQl8EIJI70H8SBykJLF/G7qTxYtVk2pZwezH/pPeBWcQGBazHv8JRiFLunkR+89/E6v/5cPsveiteNEki3/2j/iRGN1Lz2fXFTex7D/Wsv/c63CSDcx49pf4kRjpmYvpXnIeDVs2kJqzAs9OYPfsx7dieLFq3HhNmPxihA/TUsSE6/j0dhfIZz0K+TDDRX1jlEhUx3MDrIhOX08h7HzUBbquEY0bmLpEz/QjolFEOkVV+0tE8n04dbPIN8xm3m++T9s5byDavY+5v/0hRrafA+e8nvbVr2bhf3+Fna9+L7pToGbHU7ixKvrnriLX2Eq0Yye5xlYCwyomCjq+M2GMl0920kV27dq1ciw+2Z+t/yR6PIEZGSyyKc2l09zI1kiaD/ZGeM6qwwsWUh9MzNQtYyEIHAIvU/S3thCxZiKETU9HnmSNxXMbO8hlPRJJk6oai9oGm4YZoTvjmP2OMkDzHALTxu7aS2Da+FaUql2b0AtZAtOie/nFtDxyD5kZC8g2tTLvN99DL2Rxqhp4+XWf4LRvvY+u5RfRvfR8Vt59K74VIz17Gdve/BmW/fvtdC+7gL7W05j1+E/wI1HytbNpP/0qGjY/TLZxHvnamcQ6duFHYnh2EjdZFwrVcbSwhe+FCWkCH7unDd3JERgmmRkLaHru12SaF+BbUWY880v0QpZ83Szazn49K378GQ6edhVOVQNL7luP7mTpnX8G26+/ldX//CH2XvRWco1zafn9v+FHYqRmL+Xgmdcy+w//Se/CM3GSDST3Po9vxShUNZBrnEukrx03VkNgRiYkCY3vBXhegO9LutrzFHIe+ZxPPudRU28TSxhkUy61DTYH92exozoRO2zkxhMGQhN4boCmCdrbsuHwWz8cPNEwI4quazgFH8PU2LGtjyAIjyUDmLewCjum4xQC7KjO3p3pUKCLaZdbov1otXWI/l7qOl+iI2fhVTcg40lmbvk1zpkXUdv+Ik1//Dm6pbPrjOvIzl1Oy8b7aL/geiJOimT7DtxYFZnmxThV9Ri5NE6i9pivp7GKLIAQYp2Uci1MY5H95Vc/iiMaQIQiG9P7+UP1H2mUO9hsLuXcbCNN7nK00WRlmWB8Pwt+HvQIdnQ+ttWENkx6Rc8LSPU5pHodpIRkjcXeV1IsWFpNe1uWqhqLRJVFxJ7ESQalRHPz6IUcIvBwqpuo2vkcbqIWz05St/UP6E4Oz07QfsY1LPifr9LfuppMUysLf/G1sOlb18K2N3+Gc794HfvPu57uxeex+P5/wLeipGcv45WrP8TCn99Jz5LzSM9cROPmh/GtKIWaZnqWnEvNy0+Sr2kmMO3weIUsbqKW9tWvZvHPv0z3kvPJ18xgwf9+HcPJkpmxgO1vvIVzv3gduy97N6mWFSx44Ov4VpT+uSvZvebdLPzf/0fX8ovI1bdQt/UP+FaMfO1M+ltPJbnnefK1M/GiSfRCDj8SDf3kk/Q/kFLiexKn4COE4GBbhnzWLwtpVa1Fw4wovV0FmmbGOLAvgx01sKM6dtQgGjewIhPzcCt1SobRIWFd8zmvQqQl0Vh4/Qd+mLiluyM36Pu6RjvMOuh4GG6enXtd/ACk4+BhsLA2S1X3LnzXRzY1sy1VRzTdgVfdQMRJc9rBDViGoHPRuQT1M0i/tIv8rAUYvkPESVNTY+LVNuFFYhimRj4XFK10gQx89MA5eUV2wz9/gFQyQhtb2GSafKZzNz80r2CBE8UOatENWZFgYvKQUhL4aQgcNLOGaLQVy6hFjKFyvi8p5D0O7M3Q3+uQ7ndw3YCWeUnqm6JkMy41dRH6ex0ixZvIjuro02HCwaJViQyI9HWgF7JITSM7YwENW35Lrn4OXqyKpmd/Vbaq9593A8v+/XY6T7mc9KzFND/5P/hWlFx9C52nXE7Dlt+SmbEQp6oBu3s/fiSKZyfwYtUDYV5TnFI+glzWwzA19u9Ok6+wRBNJk3kLq2hvyzJrbqJCRAeE1DgeI2QmkcpQP9cJwgEkvkTmC8SDNHo2hWNXIZEUXtxBf9MijJ6DmJ37mJvbQc+Cs3DNKM1/vJ8/rHw35DLIQgFdg4UXLpj+IjsWn6zr5nj4R0u4uWkp3+7aw97Ya1nlVyEReK7AdQSFvEbghxeXrgflHKvHi0D6BF4KgcSwmrDtOVjm+Me3lp76YY4BnyCAtj3FGzHrU8h7tMxPUlUdwXUDEkmT/t4CkYob0bTGOEW64piQUuJ5YRLwVL9DNGaw++X+soDmcz7xpMnSVbXs25Vm7sIqDuzNYMd0bDv8/0WiuoqwGCeE7yKFhpXqwu5pQ+ayHJyx4uSMkzXNKH1yEd/tPY2EfQk1AeVBQqVhs7FEgO+B62o4BQ3XEeWEwBNp5Vb6WyP2HOzIzHHNKXAoerFpY1pauelVUzcQDDkwtDUcHuoWfHJZj56uQvlGbmlNlN0OEVunv9cZuJFjBhF7at3IQSDxitO+e0UjwRADEx9OlboGvsRxQl9lT2ee6toIL2/tHfQATFRZrDyjgfb9WeYvqcaOGdTU22UrtPQAXHpKmFVt/pLqST6rE5dSeKRT3YRT3YTnuDCGVIcnRJwsQJO7FC1y5Mhq3QDdCLCjxSGVnsBxBE5Ow3M1EBJdk+Ni5fp+DulnEbpNNLEE25oxrL/1eFIa2hqNFZ8qCZOa+sFRAVJKnLyPBPJZjyCQdHfkyxbV3PlJfF8SiRoIAel+Z1CTNBozwtl6x5GSkHp+MJBJv/hH1wS2qZOwrXIOgrznU/B8Co6PW87wDOWpUbTSpIXjI8JSSjw39IXaMZ32/VkammNs29RdfHh5uE5AosritHMb6enMU9tgh4nmZ8awY4NdOSURnbvg5BjNdzIx+SpwnNA00EpWbjwg8A+3cjWtOJfWCPWi0t8qjGriVauL+QSmlw9MCEEkGl4KdtQ4TIQBnIIPArLpUEAyKZeugznyeZ95C6vo73OorrVwnYBMyh3kE4wnzSH9gkEg8fxQSP1wtqhyOg5TE0QsnaRtYZs6hh7OTGvqAv0ooTwBQbFciRcEOJ5P3gnIuz55x8cv51QN0QlnhtW0sNdbFIfpOgWfQt4nWR12Ora0Jtn8VEf54QMQi5ucccEMUv0uM2YLGpqj5XOv7JRcsioU0dnzxjuvmGI8CVNOSnwpKRQ8jHFoFE1bkRUC8o6PkF44f48oZoQqjnU+ko8xdBmM3cod8LcGGNaMCfO3TiVKPdBWnT7IHVGiYUYUBGRSLlJCIefR31cgn/WZszDJwX1Z6mbYZNMu+ZyPZevEogaJpMmMBptYzMIyNEw9bPIfTUiPhIaGpYezOAyFF4QJafr6HVJph+o6kxf+3MOs1iTPPNlBIefjFnwMSyNi66w8qyFM54hk9rxkKKBRA7PiwbFkZS0AzbOnTqigYmgqhdTzZTivVxENgWVoRC2DWlvD9I5dIiddZMeaT7YhaePbNhgRXD/A9SW+74dTqHhDTBVYbDrqQiC00iRrDMxZX2Hl+n7YgebkNVx3wMoVogAyA2iYkRZi9qwJ9bdOB2QxbZxPGI+pRTVqbHvAR65rRC2deS1JIpZOJuWS6ndCazjj0b47y5z6OH/6wwFOWVXP7j1pMhmXRMIkHjdJxE1mzoyjaWAYI+ugkzLMcZvOuKTTLnNaEjz7XAerT23god/soaengOcF5fIvv6yFhmqbhTOSVJ8XhjSZtoEkwPMkec8nutzA8Xwi1QYSKPg+ed9HL87MGmZxElM6K9rJxEiENGYZRCMalqFjaAMPeFFs5gSuh5s9/tPPjDtjTRBj6RpEDAxzmDjTUj7QYq7QoPhje16AG4TNSdcLp+0YKpRW0wV2UmBLgefkcHMFPDeJEVmEHWnENK0wA9RJgCxm6vdLF27FBSsAy9CJWjq2ZWIbYdPe0ItN8EN+3LpEBGYefowrXz0XKUN3RV+fQzodCuSBg1mqayI88cc2zjpjBlv+3EWh4BNPhAIZjRrMb61i0+ZOVq6o53e/30t7ew7D0EgUt2meEaOqKoKUcN65zURtA/uQ+OJTVoU5eZtnHNkSDQhwvfD38Pwg9AW7AQXXJ+f45d+mVLJWdEXomig+qE+Oa+Z4UBJSrzSTxBBCGo8Y2NbwQjoRnDAdXwBaPo9AQ5qHD5MtTcRnHkUIg3IS32IiXxnG2IV+vB48r4BuVWMnV6OLBgI0PMcnm/PxXB8hBYYhsEwd3dDQBdPSmikJafmCrfBbaggihkYsYhCtuGANPfSRjscFW8q4VVdrU1d7uE/4mqtaATjj9EZSKZd0xiWTcenszNHSkqCpKYauCy6+cBaRyOCmPMCSxTUA1A3nQxghGhqRirvmUCeRFwR4QfgA9wKJ43rkiiLsuhUyIIszR2ihBayPU4fcicZh1+UQQpqwK4RU1zE1MeFCOhqmrcjm5rYQyUm0bB6tP4UoTVNsGEjLRBrGiEIGtGIOypIoe0GAG/RiGD5NNa3U1SwhaTciGehIKXeouGGTNJVy6Eu55LIufiARBsX8pWJgQs/ifO+6FgrxePVyj4bS3FGli1ZWXLA6AsvUSNomtqVjGmJKXrDV1RGqqw/3CSfi4YM2MkGjl0aKoWkYWuXMCwPDvoPyNSTxfB/Xk+Rdn4IbRkZ4wYAVLGUxL+tJ4Io4EYT0SExbkfWrEjgN1QS6Cb6PVnAQhTx6Oo2WSqP3p4Bib7WhI00LaQ4vvF5QwMn3IhBU1yyjvnohtjkw+7lAYOoCU9egwnCeUQxblFLiOAH5vEdfqkA25+F7EqmFsay+BMfzcf0AxwuKnwesxdKNBcVExEUR1sXompmVMaT+oaFPIuyxj9tG2Kw3BKauY+oCY5pFRExHNMJJFy0dMA9/GPgyqBDhgEIxKqLg+mQdf5D4CErGQeiKmOquq6MJacSc3kJ6JKatyA5C1wliUYhF8WvDXt4B4S0UhTc1pPA6OLheGt2waWo4k5rEXExj9NmlhBBEIjqRiE51dQTfD8gXfDIZj0zGCcdpRzQsU8coxpRKBizL8AIME2u4XoBTFGPXD/Dcwc13CD8KBIPt0dBNEjHCGNJoReiToSkhneroIpw0cajo79K14hZD3jy/GJLmhS2qXKUrgoHY4OPpijiZhfRInBgiOxSDhLcmXFcUXvI5vJ69iN6D2H6UpurlVMVmoOnRcDKhcRjXrusa8ZhGPGYiG2wcJyCX80inHbJZFwmhAJpaWXSPhC+D4gVcnNxPSlzfR9f0soV9rKFPiqmLQBQflAxpBZddEX5QjA0OKLgVscEM5YooWsGjaCkdKqS+lOVb5WQW0iMx6SI71hCuseALSUb2IE2f5OKzmNlwBtX2DMhnIZeFVA/0dkNvF8X5rcG0IGLDMaShq7Rya2oieF4Yp5lJh503QRCKunUEwdWFhl6K/RziJlOc3Ay4Ioa+fg5zRTh+OTIiWwgYbAeHw5SFEEN3giohPSInTAjXaPD8HPlcFwJBfcPpzKhdQcyuG9ggngxfDTPCZd8fWnhLF6Jugj124TWMUEzjMbPsy83lPFIph0w2nKqhZOWqHmfFeHCYK6JitqVKV4Tr+3iepOD5+IFUQjoGTqgQrqPhOCmcQg+GGWd288U01izDNEcwGaGuDy28hRxkM8MLb8QGa3TCO6SVW/BJZxzSGTd0XxzFylUojoVKV0RUtZLGnRNOZGUQkHe68N0cEbuB+fNeS12y9diTteg6xBLh61DhzWUh3Qs9XUXhLaYFG4Pwlq3cuEmTlBQKPrm8RyrlHubLVVauQjH1OWFE1g88CrlOgsAlWb2QmfWnUZ1omdiDVgpvfRPMA4JgwNWQ7g0t3r5ukMWhvro1YuEVQoTpBm2D2hobzwsjFtIph0zODadxEmLEnWcKheL4M+1FdsDfCnX1p9Jcd+pgf+vxRtNGIbwBYQCrCZHoUYXXMDQShkYiboZZopzQyk1XWLm6JrCsqZX/VaE4mZnWIpvL7MM0E8xuvpiGmiVY5hTNgDSs8OYglxlCeAnThFlFi3eI+FZNO9zKdRyfbNYjlQ7jcoUIOypMU818oKhABuF4cumH12EQhBNZBkFxRExxEIssjlU0TTAsMIzjPmPsicC0Fdna2pXEo43UVy2YEsmxR42mQSwevg4V3nwWUr2hf7e/Z0TCW/LlxmIm9fU2rhtQKISZrjLZMP2gQGBZyrVwwiCDsF9ASijOvlsWy0NCsgYl6xUaGGYompaN1A08TceXOr5mIDS9mHYOcJzQEMhnkb294TxsoligZoRlGGbF9Vh6mA83rVXF94duKgAZrpQizC+MPOTr4k6i9CCoKLGy/EE2hRj0Vvz+cKOjcpUQEHjBMOcwOqahOoW0zrxwsqsw/lQKb10jzFt8iPD2QV9JeP3i1SdK2caL7xpC07E0DSuqk4ybBMSUa2GqMkKrchAiTL+JboBphiJn20Vr0wzX6UZ4Teil68JAagI/CGNffdcvC6/QBJGYRTwWIRKzME0Dw9LRjYFIAyllWB/HQRZyUChAug+Z7odUP3hOuF1YQbAiyLL1K8pDxsOyGFwuw303UGJlYp3BnwdmwSiVVapquLq4LhgYGSkrtDOoWFnep3QMXWJZhyefGi3TVmRPGg4VXhYNFl7fA9cpvtzwYndd8FxwcuC5aL6PDdhAbaI4660bkM35pHp8Aj8cDGGYOqalD1gymlZMuqupZuLRGM6qlMFg5YAjWpVlkTTMcCCMphWF0qj4XLI0h384Bn6A5/r4no/0ZFlZLNskVmVjxyzMiIlp6ujm0aeVF0KEx7NthF0cdt7UPLCB54GTH3CBpfog0w+ZXhAl1RLhuZmRgXObyngeFEafT/ZQJl1kj+eIrxOGSuEdCVKGAhD44PvogU/M94kFPvWeh5t3KKSzZLozZPpzSM9FeD6W5mMIL7zYgopE6AOpxQYbWkIMFuXDhFowZPLeqcJYrEogDNcrNp1Na0RWJXpF6+MYCAKJ7/n4ro/vB+VfVzd17FgEO25h2SaGaWBaxsRl8jIMMIr9DjRWVhAK+QEBTvdDui/sh/Arrim9+NuZkbCsacxJOeLrpEeI4oV7+L9bECbjs4AkEAQBTt4lnymQ6s2SzTphy1QPc+ZqGuAHoRANeg9CMfbc0JLz3IHl0jrfG0llS1NWUJ77p/xeFOujCfVEW5WiKJgjtCrHAyklvhfguz6e5wMSIQVCLzb1q6LDNvUnFU2DaCx8VQMzZofrpQxbX4V8+MqkQss33QepwsD+QgxYvqY14b/zeHBSjfhSjB5N00ILKBahprEK3/Mp5Byy/XnSfVn84qy2pmViRsdgGflFsRtOqH2vQpy9AcH2PXAKEHiHC+UgRCh+JWE8jlbleOF7Ab4XJoWXstjVIwSWbRKvLorpKJr6UxIR+m2xIpCsHhjgAwOuh0J+oC8i0z94hCUU/6cRsKzw/zdFUSKrOCK6oRNLRoklo9TPqsEteBRyBTJ9OTL9+VAEhMCKGBjWCC4nXQf0Y7vyDm3KQ9HKNI6LVTleBMUOKM/zkf5Ab4xhhU39SNzCiliYlo5hTmBTf6oxyPXAwHRFQRA+aAv5cKRlpn/A/eBVtJJKD9gp4nqY/Boopg2iaE1ZtkmyNjHItZDuzZLtzxU7vTVMy0SfqFCxsttgely+lU193x9wWQhdYEcjJGqiRKIRDEvHtAy0YTJpnfRoGtjR8EUtMGvgu8NcD32hBTwFXA/T4ypVTEmO5looZEuuhQnudJlCDNvUj5okqqNE4hamZWJYBoZKxjJ+lIQzURXGnZfwim6mQm5414ModlxOkOtBiaxi3DjUteA5HvmsU3Qt5EbvWpjCBP6AmAZ+AMX53EpNfTseCf2mJ1tTf6phFKM+SpE4JdeDlAOWr5MvRjwU3Q+lDloZhMPdj7UKx1yCQjEEQohQZCImydo4MpAU8s7hrgVdw4xMoGvhGJHFECnP8wdGAAnQNJ1ILDw3y7ZUU3+6IUSF6wForJinvtL1cMRO1pGhRFZxXBCamPKuhZJl6nt+xXhNgR21SNba2LEBv+mUCZFSjD+VrodxQImsYlKodC00zK7FLbiDXAuledasiDHuYUqDmvpBcY4qCWbEJJYIxdS0jbLfdFqGSCmmDEpkFVOCI7kWcqn8mFwLR2rq23GLZF0CyzYxrTCAX5sicbKKEwslsoopx9CuBZdsKke693DXAmKgV993w6GakmLSc9XUV0wySmQVU57QtaATS9o0zBratWDaJrFksakfUU19xdThuIisEOK3wGellI8ej+MpTmwOdS1IpGrqK6YsEy6yQogrgcxEH0dxciI0NUW1YmozrMgKIZqBLwCrpZRnV6x/FXA90A5IKeW6I5QhgLOAJ8etxgqFQjGNOJIlexFwP3BaaYUQIgZ8G1gppSwIIe4VQlwBHCQU5EreB1wM/Ax4y3AHKeWTBVROWYVCMW3ZsGFDKZcsjCSfrJTyJ0KINYesPh/YJaUsZV34A3CtlPITwBsPLUMI0UqYwfcsIC6E2Cal7KjcRuWTHTsbNmxQDyWFYgIZzT1WaSSuW7duZ2n9aHsLmoBUxXJ/cd2QSCm/DPwSCAAf6Bvl8RRHoOKpqVAoJoDxuMdGK7LthAn0S1QV1w2LlHKnlPL1Usq/k1I6h35fcheM9mTG4+SPpYyJrK8Sz5EzlX+r41m3iTjWZN5jY9lvqtxjh04/Iw6dLbKSorvgS1LKs4rLMWATFT5Z4JtSyt+MtUJCiO8Ce8ewaysVUzyMkWMpY7T7jmb7kW472jqciLQydX+DVo5f3SbiWONR5ljLGMt+o9lnpNuOpR4ALVLKm+EIIiuEuBR4F3A18C3gy1LKnBDi1cCbgA7APVJ0gUKhUJzsHNGSVSgUCsWxoYbJKBQKxQSiRFahUCgmEJUg5gRCCPF6YBlgAi9KKf9rkqukUJxQCCFWA2cTRlbVSClvO+o+yic7tRnN8GYhxGwp5T4hRDXwPSnlmyan1grF9GG0KQSEEAuATwI/lVI+dLTylSU79Rnx8OaKULrrgC8d95oqFNOTUd1jUsodQohbgB8CRxVZ5ZOd4kgpf8LgUXYwzPBmACHEtcAOYN9xq6RCMY0ZzT0mhLiquE+awQOzhkVZstOTIYc3CyHeCNwKPEd4Abzj+FdNoTghGC6FQKMQ4v8Spgq4ayQFKZGdngw5vFlK+TPCrGcKheLYGO4e+9FoC1LugunJ48A8IUSkuHwh8ItJrI9CcaIxbveYii6Y4qjhzQrFxDLR95gSWYVCoZhAlLtAoVAoJhAlsgqFQjGBKJFVKBSKCUSJrEKhUEwgSmQVCoViAlEiq1AoFBOIElmFQqGYQJTIKk4KhBA3CiFqKpZ/JYQYdjp7hWK8UCKrOFm4EaipWL5aSnnE6ewVivFAJYhRTBmEEH8F3AF8BVgALAdeK6XsH2LbDwIrgYPAPOADUkpPCLEc+FvgBeAU4D8Bh3Bq548LIbYCzwNfE0J8XEq5QQhxIfBu4CXCmSU+A8SBHxBOB90GXAD8f1LKXxS3fw+wDTizeOye8f9FFCcEUkr1Uq8p8wI2AFcVP38DuGGY7V4HaMXPXwOuLX5+HDin+Lm5Yv0GoLVi/7uANYAA9gKNxfV/Afy4+PlG4N+Kn88G/rv4+Z+AzxK2BE8F4pP9u6nX1H0pd4FiKvJi8b2D4RMjZ4F/EEL8HbACaCyuP5XQIkVKeUBKebTMSQ1AlZSyo7j8ErD6KHX5e6Ae2Ai8lzC3qEIxJEpkFVORkWQt+gnwNSnlekKxK/EcsBBACDFbCPG64no/XCVOEULoFdt3An0VnWCLgWePUpfzpJQfJ7Rum4BrRlBfxUmK8skqpgzF1HLzgPcIIe4CLgFOEUL8osLSLPFt4BtCiEcJpwpZKoT4BXATcIsQYhvQQujjBfgl8HeADfwLocX7TuBPwFuBO4QQLwNLgU8KIWYQuiRqhRCLgL8kzC96BdAqhPgKYWLnHKErQqEYEpXqUKFQKCYQ5S5QKBSKCUS5CxRTGiHEaVRM1VzkWSnls8e9MgrFGFDuAoVCoZhAlLtAoVAoJpD/H0RE+5F/sJPYAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 360x216 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "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>IPS</th>\n",
       "      <th>DM</th>\n",
       "      <th>DR</th>\n",
       "      <th>Learned MIPS</th>\n",
       "      <th>MIPS</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>n_actions</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>50.0</th>\n",
       "      <td>0.62±0.12</td>\n",
       "      <td>0.48±0.15</td>\n",
       "      <td>0.48±0.15</td>\n",
       "      <td>0.76±0.24</td>\n",
       "      <td>0.92±0.24</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>100.0</th>\n",
       "      <td>0.93±0.53</td>\n",
       "      <td>1.41±0.48</td>\n",
       "      <td>1.41±0.48</td>\n",
       "      <td>0.60±0.25</td>\n",
       "      <td>0.52±0.22</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>200.0</th>\n",
       "      <td>0.84±0.21</td>\n",
       "      <td>5.93±1.09</td>\n",
       "      <td>5.93±1.09</td>\n",
       "      <td>0.52±0.25</td>\n",
       "      <td>0.37±0.17</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>500.0</th>\n",
       "      <td>0.70±0.21</td>\n",
       "      <td>27.01±1.41</td>\n",
       "      <td>27.01±1.41</td>\n",
       "      <td>0.63±0.25</td>\n",
       "      <td>0.57±0.20</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1000.0</th>\n",
       "      <td>0.73±0.15</td>\n",
       "      <td>64.23±0.65</td>\n",
       "      <td>64.23±0.65</td>\n",
       "      <td>0.51±0.28</td>\n",
       "      <td>0.52±0.22</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                 IPS          DM          DR Learned MIPS       MIPS\n",
       "n_actions                                                           \n",
       "50.0       0.62±0.12   0.48±0.15   0.48±0.15    0.76±0.24  0.92±0.24\n",
       "100.0      0.93±0.53   1.41±0.48   1.41±0.48    0.60±0.25  0.52±0.22\n",
       "200.0      0.84±0.21   5.93±1.09   5.93±1.09    0.52±0.25  0.37±0.17\n",
       "500.0      0.70±0.21  27.01±1.41  27.01±1.41    0.63±0.25  0.57±0.20\n",
       "1000.0     0.73±0.15  64.23±0.65  64.23±0.65    0.51±0.28  0.52±0.22"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import os\n",
    "import seaborn as sns\n",
    "import matplotlib.pyplot as plt\n",
    "from joblib import Parallel, delayed\n",
    "\n",
    "os.environ['PYTHONWARNINGS']='ignore::UserWarning'\n",
    "plt.style.use(['science', 'no-latex'])\n",
    "\n",
    "n_actions = [50, 100, 200, 500, 1000]\n",
    "train_rounds = 1000\n",
    "\n",
    "# change the arguments to n_iter=30, n_iter_big=50 to reproduce the results in the paper\n",
    "results = pd.concat(Parallel(n_jobs=-2)(delayed(ope)(a_num=a_num, train_rounds=train_rounds, n_iter=5, n_iter_big=5) for a_num in n_actions), ignore_index=True)\n",
    "plt.subplots(figsize=(5, 3), tight_layout=True)\n",
    "sns.lineplot(data=results.set_index('n_actions')[['IPS', 'DM', 'DR', 'Learned MIPS', 'MIPS']])\n",
    "plt.xscale('log')\n",
    "plt.yscale('log')\n",
    "plt.show()\n",
    "\n",
    "def mean_standard_error(serie):\n",
    "    return f\"{serie.mean():.2f}±{serie.std()/np.sqrt(len(serie)):.2f}\"\n",
    "\n",
    "results_a = results.copy()\n",
    "results_a[['IPS', 'DM', 'DR', 'Learned MIPS', 'MIPS']] = results_a[['IPS', 'DM', 'DR', 'Learned MIPS', 'MIPS']].multiply(results_a['sample_size'], axis=0)\n",
    "results_a.pivot_table(index='n_actions', aggfunc=mean_standard_error)[['IPS', 'DM', 'DR', 'Learned MIPS', 'MIPS']]"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3.9.12 ('obp')",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.8"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "9dd8cb6e9f7fdb75e7c84849fd3a485e3c0a281ba0b281a70d444f989a7c428c"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}