{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# LASSO when conditioning on less: ROSI\n", "\n", "Instead of conditioning on the active set and signs, \n", "one can work in the full model and for each feature $j$ selected\n", "construct p-values and confidence intervals\n", "conditional only on the event $j$ was selected.\n", "This is the approach of [Liu et al.](https://arxiv.org/abs/1801.09037), which\n", "can be extended as ROSI (Relevant One-step Selective Inference)\n", "beyond squared-error loss (described in forthcoming work, though\n", "code is already available).\n" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import numpy as np, pandas as pd\n", "import matplotlib.pyplot as plt\n", "import statsmodels.api as sm\n", "%matplotlib inline\n", "\n", "from selectinf.tests.instance import gaussian_instance # to generate the data\n", "from selectinf.algorithms.api import ROSI\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We will know generate some data from an OLS regression model and fit the LASSO\n", "with a fixed value of $\\lambda$. In the simulation world, we know the\n", "true parameters, hence we can then return\n", "pivots for each variable selected by the LASSO. These pivots should look\n", "(marginally) like a draw from `np.random.sample`. This is the plot below." ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": true }, "outputs": [], "source": [ "np.random.seed(0) # for replicability\n", "\n", "def simulate(n=500, \n", " p=100, \n", " s=5, \n", " signal=(5, 10), \n", " sigma=1): \n", "\n", " # description of statistical problem\n", "\n", " X, y, truth = gaussian_instance(n=n,\n", " p=p, \n", " s=s,\n", " equicorrelated=False,\n", " rho=0., \n", " sigma=sigma,\n", " signal=signal,\n", " random_signs=True,\n", " scale=False)[:3]\n", "\n", " sigma_hat = np.linalg.norm(y - X.dot(np.linalg.pinv(X).dot(y))) / np.sqrt(n - p)\n", " L = ROSI.gaussian(X, y, 2 * np.sqrt(n), sigma=sigma_hat)\n", " soln = L.fit()\n", " active_vars = soln != 0\n", " \n", " if active_vars.sum() > 0:\n", " projected_truth = np.linalg.pinv(X[:, active_vars]).dot(X.dot(truth))\n", " S = L.summary(truth=projected_truth)\n", " S0 = L.summary()\n", "\n", " pivot = S['pval'] # these should be pivotal\n", " pvalue = S0['pval']\n", " return pd.DataFrame({'pivot':pivot,\n", " 'pvalue':pvalue})" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's take a look at what we get as a return value:" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "Index(['pivot', 'pvalue'], dtype='object')" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "while True:\n", " df = simulate()\n", " if df is not None:\n", " break\n", "df.columns" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "dfs = []\n", "for i in range(200):\n", " df = simulate()\n", " if df is not None:\n", " dfs.append(df)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeMAAAHSCAYAAADfUaMwAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0\ndHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOzde3yO9R/H8de12cz5fByico7CnI/Z\nHEKI2CgSRUgnogMRUTlVkgoVyWEIySlnOklUv4SSYznkfDY7Xr8/vubebM7bfd27934+HnvY9b2v\nts8ye+/7vb4Hy7ZtRERExDk+ThcgIiKS3imMRUREHKYwFhERcZjCWERExGEKYxEREYcpjEVERByW\nwalPnDdvXrt48eJOfXoRERG32rx58zHbtvMl95pjYVy8eHE2bdrk1KcXERFxK8uy9l3tNQ1Ti4iI\nOExhLCIi4jCFsYiIiMMUxiIiIg5TGIuIiDjMsdnUN+LMmTMcOXKE6Ohop0uRa/Dz8yN//vxkz57d\n6VJERNIkjw3jM2fOcPjwYQIDA8mUKROWZTldkiTDtm0iIiI4cOAAgAJZROQWeOww9ZEjRwgMDCRz\n5swKYg9mWRaZM2cmMDCQI0eOOF2OiEia5LFhHB0dTaZMmZwuQ25QpkyZ9DhBROQWeWwYA+oRpyH6\nuxIRuXUeHcYiIiLpgcJYRETEYQrjVDRkyBAsy7r8VrhwYdq2bcuuXbsA6NKlC0FBQSn+eXfs2MGQ\nIUM4depUin9sERFJeR67tMlb5MiRg2XLlgGwe/duBg0aRHBwMFu3bmXQoEFERESk+OfcsWMHr7/+\nOl26dCFnzpwp/vFFRCRlXTeMLcv6FGgBHLFt+55kXreA94BmwAWgi23bv6R0oWlVhgwZqFGjBgA1\natSgWLFi1K1blyVLltCuXTuHqxMREU9wI8PUU4Cm13j9AaDkpbfuwIe3X5b3qlKlCgB79+5NNEy9\nZ88eLMti8eLFie6PjY2lYMGCDBw48HLb6tWrqV69OgEBARQoUIBevXpx7tw5ANauXcuDDz4IQIkS\nJbAsi+LFi7vhKxMR8R4XL8KuRdvd9vmuG8a2ba8HTlzjllbA57axAchpWVahlCrQ2+zduxeAggUL\nJmovUaIE1apVY/bs2Yna161bx+HDhwkLCwNg69atNG3alLx58/Lll1/y+uuvM2PGDB5++GEAKleu\nzOjRowGYN28eP/74I/Pnz0/lr0pExHtE7D7EhhJhlHiwPH9N3eCWz5kSz4wDgX8TXO+/1HYoBT52\nIp6wlNW2b/6/iYmJAcwz4169epEtWzZCQkJYtWpVovvCwsJ4/fXXiYyMJGPGjACEh4dTvnx57rnH\nPCEYNmwYd9xxBwsXLsTX1xeA3LlzExoayo8//kjNmjUpXbo0AJUqVVKvWETkRsXFEfX+x1x8YQBz\n4s5SHIju1oP/Gm6iYFG/VP3Ubp1NbVlWd8uyNlmWteno0aPu/NSOOX78OH5+fvj5+VG6dGl2795N\neHg4hQolHTxo3749Z86cuTzhKyYmhnnz5hEaGnr5no0bN/LQQw9dDmKAtm3bkiFDBr777rvU/4JE\nRLzR//5HbI1a+D/Xi4FxZ5kArAFiy99LwRwpP9H2SinRMz4AFE1wXeRSWxK2bU8EJgIEBQXdQh8z\n7cmRIwcrV67EsiwKFixI4cKFr7pbVWBgIHXq1CE8PJxWrVqxatUqjh07dnmIGuDQoUMUKFAg0X/n\n6+tLnjx5OHHiWk8TREQkiV27YMgQ7OnT8b009NkfuIMC+D0+g3s/beiWMlKiZ7wQ6GwZNYDTtm2n\n+BA1mCFip99uVoYMGQgKCqJKlSoEBgZed9vI0NBQvv76ayIiIggPD6dSpUqULFny8uuFChVKciBD\nbGwsx48fJ3fu3DdfoIhIerR/P/ToAWXKwBdfcMG2GQ1E4MdnDMYaupdH3RTEcANhbFnWTOBHoLRl\nWfsty+pmWdZTlmU9demWJcBuYCcwCeiVatWmA+3atSMiIoL58+czf/78RL1igOrVqzN//nxiY2Mv\nt82bN4+YmBjq1KkDgL+/PwAXL150X+EiImnBjh3w5JNw110wcSLExHAWsyxoAFCaz8kycggvDgpw\na1nXHaa2bbvDdV63gd4pVlE6lz9/fho0aEC/fv04deoU7du3T/T6wIEDqVSpEq1bt6Znz57s37+f\nAQMG0KRJE2rWrAlweQLXxx9/TFhYGJkzZ6ZChQpu/1pERDyCbcOGDTBmDMybl2iY8xRQm2xs4wIw\nnefGhPLCC+4vUdtheqCwsDAOHTpEjRo1ksyGLl++PEuXLuXIkSO0adOGgQMH0qFDB+bOnXv5njvu\nuIPRo0czb948ateufXndsYhIuhEdDStXQp8+cMcdUKsWfPlloiD+o0gQd1CKbVwE5jBggDNBDGDZ\nt/IgNAUEBQXZmzZtuurr27dvp2zZsm6sSG6X/s5ExFG2DT/+CNOmQXg4nDyZ/H3NmrGq6kuEvG4B\nDwLTCAtrwfTp4JOKXVTLsjbbtp3sgQTam1pERNK2Cxdg7FiYMsXMjk5OrlzQqhU8/zxLD5SmdeuM\nl17YQ716OZkyJXWD+HoUxiIiknbt3w8tW8KvvyZ9rWhReOghaN0a6tQBPz8WLDhImzaVsO0Xgccp\nUyYn8+dDxoxJ/3N3UhiLiEja9NNPJmj/+8/VliMHtGsHnTtD7dqJurtff/0vbdo0xLb/A+4mMBCW\nLgVPWBWqMBYRkbRnxgzo2hUiI811hgwwerRZOxyQdFnS4sV7aN26IbZ9AlhOgQI1Wb0aPGXHYIWx\niIikDUePwqxZZoLWzz+72nPnhrlz4f77k/3PVq48TcuW9YmLOwesIk+eIFauhFKl3FP2jVAYi4iI\n5zp2DBYtMsuSli2DSwfvXFa2LCxcCHffneQ/tW2YMAGefz4HcXEvAA3ImfM+VqyAS2fveAyFsYiI\neJbdu+Grr2DBAvjuO4iLS3qPnx906ADjxpnnxFc4fx7CwrayaNEFoCrwHLlymWfElSql+ldw0xTG\nIiLiLNuG336D+fNNAG/ZcvV7a9WCTp2gffurzrzauxcaNfofO3eGAPmB36lc2Ze5c6FEidT4Am6f\nwlhERJxz8CB07w6LFyf/umVBzZpm1nSbNmZP6WvYvBkaN97EiRONgSzAArp182X8+GTndXkMbYeZ\nioYMGYJlWZffChcuTNu2bdl1tUXpKcyyLMaPH++WzyUiclNsG6ZPNw9vrwzijBmheXOYNAkOHYLv\nv4cXX7xuEC9dCnXqbODEiWAgB35+65k0qSSTJ3t2EIN6xqkuR44cLFu2DIDdu3czaNAggoOD2bp1\nK1myZHG4OhERBxw+DD17mmHphEJDzRrhJk0ga9ab+pCffmo62LGxHwD5yJ59NV9/XYx69VKu7NSk\nME5lGTJkoEaNGgDUqFGDYsWKUbduXZYsWUK7du0crk5ExM3mzDFBfPy4q614cfjsM2jQ4JY+5Oef\nQ7duNmABkwkMPMXy5QUoVy4F6nUTDVO7WZUqVQDYu3dvsq+XKFGCF198MUl7u3btLp9XfP78eZ5+\n+mlKly5N5syZKVGiBL179+bMmTPX/NzFixenX79+idqmTJmCZVmcO3fuctuJEyfo3r07BQoUICAg\ngFq1avHTTz/dzJcpIpLY8eMQFmYmXiUM4h494PffbzmIDx6Enj2XAzWAY9x3X0Y2bkxbQQwKY7eL\nD+GCBQsm+3r79u2ZM2dOorZz586xePFiwsLCALhw4QKxsbEMHz6cpUuXMmzYMFavXp0iPe3IyEhC\nQkJYuXIlo0aNYsGCBeTLl4+QkBD+S7jlnIjIjTh2DD74AMqXNycpxStSBL75Bj76CLJlu+UP//DD\ni7lw4UEgkhIlbNatg8KFb79sd9MwtRvEXFqkvnv3bnr16kW2bNkICQlJ9t6wsDBGjhzJhg0bLg9v\nf/3110RFRV0O23z58vHhhx8m+vglSpSgTp06/PPPPxQrVuyWa/3iiy/4448/2Lp1KyVLlgQgJCSE\n0qVLM2bMGEaNGnXLH1tE0onISPj6a7NT1pIlSTfqePxxeOedZNcH34wBAxbw44/tgYrAcj75JDfZ\ns9/Wh3RM2gpjy3K6gkQHU9+I48eP4+fnd/m6WLFihIeHU7BgwcshDWbms6+vL5UqVaJUqVKEh4df\nDuPw8HDq169PgQIFLt8/bdo0xo4dy99//8358+cvt+/YseO2wnjlypVUqVKFEiVKJKqvfv36XOv8\naRER9u0zPd1PPjFbV16pUCGYOBFatLjtTzVr1iJGjmwHBAFLeeKJnFfbDTNN0DB1KsuRIwc///wz\nmzZtYv/+/ezdu5cHHniAdevW4efnd/ktODj48n8TGhrKnDlzsG2bM2fOsGzZsstD1ADz58+nc+fO\n1KxZkzlz5rBhwwbmX5qVePHixduq99ixY2zYsCFRbX5+fnz22Wf8+++/t/WxRcQLnTsH8+aZYwzv\nvBPeeitpENeoYYaqt29PkSAG+PrrykBH4BsKFsxJWh+0S1s94zQoQ4YMBAUFJWmvUqUKPyfY6Dxb\ngmcmoaGhDBs2jO+++449e/YQFxdHmzZtLr8+Z84cqlevzoQJEy63rVu37rq1BAQEEBUVlajt5MmT\nia5z585NUFBQomHweBmdPvBTRDzDiRMmgL/6ClascJ2clFDRomY4+tFH4dIjr5SwZs0azp2ry4wZ\nhYGpgMn5nDlT7FM4Im2F8U0OEXuybNmyJRvSAOXLl+eee+4hPDycPXv2EBISQp48eS6/HhERkSQY\np0+fft3PWaRIEbZv356obfny5Ymug4ODWb58OcWKFSN//vw3+uWISHrw778wdqwZar5wIfl7GjeG\nXr3Mph0ZUjZiJk6cSI8ePfDxGQWYlSEPPWQ25krr0lYYpyOhoaG89957nD59mkmTJiV6rVGjRvTu\n3Zvhw4dTvXp1lixZwqpVq677MR966CH69OnDiBEjqFq1Kl9++SVbt25NdE/nzp356KOPaNCgAf36\n9ePOO+/k+PHjbNy4kYIFC/L888+n6NcpIh7Ots0ek+PHmx2zrpyMBVChArRqZfaMTqVzCd9/fzzP\nPNMHaEZc3NMABAaaXrE3UBh7qLCwMAYNGkTGjBlp3bp1otd69OjB7t27ee+997h48SKNGjVixowZ\nlyd8XU337t3ZtWsX48aNIzIyks6dOzNw4EB69Ohx+Z6AgADWrFnDa6+9xuDBgzl8+DD58+enWrVq\ntGzZMlW+VhHxMNHRsG6d6+Sk/fuT3lOhAnTpYkL4OttU3q5Ro8bSv39foBUQDmSkfHmz/WWhQqn6\nqd3Gsh0a+g0KCrKvNTt3+/btlC1b1o0Vye3S35lIGmbbsHGjWY40a1bijTkSqlcPXnoJmjZ1ywqX\n2bMP0KFDaeLimgHTAT8aNDA7aaa158SWZW22bTvZ55PqGYuIpGfnzsGECWY50o4dyd+TK5eZBd2z\npzlByQ127jRnQyxYEAj8CJQFMhAWBlOmmLMkvInCWEQkPbp4ET7+GEaMgCNHkr4eGAht25ph6Lp1\nIcF+CakpOhqGD7cZNmwgcXGFgKeBCmTODIMGQf/+4OOFi3IVxiIi6YVtw59/wrJlZlb0lc+Cs2aF\nhx82E7EaNHB76m3ZAp072/z2Wz9gLPAUYNO5s8WIEeb3A2+lMBYR8WbR0TB3rjkzePVqcz7wlYoU\nMd3ORx+FzJndXmJMDIwaBa+9ZhMT8wwwHuhD9erv8f77FlWrur0kt1MYi4h4o/PnYfJkGDPGrA9O\nToEC8Mor5iDggAD31neJbcNjj8GMGTbQE/gYX9++vP32KJ57zsLX15Gy3M6jw9i2bSxP2I9arsup\nWfkicoXTp+Hdd2HcOLNT1pVy5oT77zebc3TqBFmyuL/GBKZOhRkzwJxFXJpChV5m1arhlC2bvn72\ne2wY+/n5ERERQWYHhkzk5kVERCQ6EENE3Oz8eRPAo0bBFdvcki8f9O5tZkTfdx+e0t3cuRN6944B\ndgDl6NLleSZNSvGNu9IEj/2S8+fPz4EDBwgMDCRTpkzqIXso27aJiIjgwIEDiU6VEhE3OX4cPv/c\nHNBw5azoEiXM+qAuXSBTJkfKu5roaAgLi+bChUeAZdx111988EGhdBnE4MFhnP3SoZQHDx4kOjra\n4WrkWvz8/ChQoMDlvzMRSWWRkWZC1rRp5s8rf0bedRcMGQJhYR7bzRw4MJLNm0OBr/DxGcOcOYWc\nmDvmMTzzb+mS7Nmz6we8iKRvMTGwZAn8+qtZlvTnn/DXXxARkfTeokXhtdfMjCgPfmy0YsVFRo5s\nCywB3uftt5+mUiWnq3KWR4exiEi6FRcH4eEweDD8/fe1761WzRxX2KWLY7Oib9ShQ9C27ThgKfAx\nwcHdeeEFp6tynsJYRMST2DYsXGjW/W7ZcvX7SpSARx4xa4NLl3ZffbchKsrsKXL27PNAFXLnDmbq\nVO/cUetmKYxFRDyBbcOKFTBwIPz8c+LXcuY0Q88VKkCZMiZ88+Z1ps5bdObMGapWfZYdO94CCuDj\nE0x4uHfvqnUzFMYiIk777jt49VVYvz5xe5Ys8Nxz0K9f2juiKIFTp05RuXJT9uzZDDwMNGfkSAgJ\ncboyz6EwFhFxwsWLMGeOOTFpw4bEr2XMCL16maMK8+d3pr4UcuLECWrVasyePb8Dc4HmdOyInhNf\nQWEsIuIutg2//w4zZ5ojC48dS/x6hgzwxBNmqNoLxm+PHj1K9eqN2LPnT2AB0Iz77oNJk9xyFHKa\nojAWEUlNcXFm+HnBAvO2b1/Se/z9oWNHM2nrzjvdX2MqsG0YP95mz54MwNdAI/LkgfnzHTmLwuMp\njEVEUkNkJEyfDiNHmnXBySlWDJ56Crp1S/PD0Qn9889hBg7MzbRp+YGNgA933GEmiRcv7nBxHkph\nLCKSks6fh48+MucFHzyY9PUcOaBZM7M7VvPmHrNPdEr5449/qFq1IRcv1gM+BXyoUwe+/NKrft9I\ncQpjEZGUsnat2Xxj797E7dmzmzXBDz0E9eubYWkvtGnTHmrXbkhU1EmgB2D+d3z4oZmTJlenMBYR\nuV0XLpiZz++/n7i9YEF4/nno0cP0iL3Y+vV/ExzckJiYC8AqoApvv23OqdBkretTGIuI3IqoKLND\n1s8/myHphFtW5soFw4ebbqGHb0+ZErZvjyU4uBUxMReB1fj43MukSdC1q9OVpR0KYxGRGxUdDR9/\nDF98Ab/9ZiZpXal5c5g4EQoXdn99Dli/Htq29SUm5lMgG35+5Zkxw2x7KTdOYSwicj22DUuXmp0q\nrjYzOnt2ePddc1hDOhmXHTz4N9544zvi4p4GapA5s1m61Lix05WlPQpjEZGrOXkSNm6Ed96Bb75J\n+nqJElC1qjk1qUOHdNMbjomBzp03MXNmYyAr0JkCBbLz1VdQvbrT1aVNCmMRETBLkn791TwDjn/b\nuTPpfdmzmx2yHn88zR3WkBIiIiAk5Ed++KEpkBtYw333mSAuVszp6tIuhbGIpF8XLsCoUWYR7Nat\nZresq/HxgSefhKFD0+2C2QsXoF699Wze3BwoCKzmoYeKMm2aOdNCbp3CWETSp6VLoXdv2LPn6vf4\n+UHFimbstUcP8346df48tGwJmzf/BRQBVvHSS4UZPlznEacEhbGIpC8HD5q1v7NnJ263LHNWcNWq\nrufAFSumi6VJ13PuHDRteobvv88OPAl0YtiwAAYOdLoy76EwFhHvFxkJixbBtGmwZIlZohQvd254\n6y2zPWW2bM7V6KFOn4ZatRaxbVtnYAlQgzffDOCll5yuzLsojEXEO128CKtXm7U2c+fCqVNJ7+nc\nGUaPhnz53F9fGnD0KFSvPp89e0KBe4FSjBoF/fo5XZn3URiLiPewbRO+M2eaZ8Lnzyd/X40aZoes\nhg3dW18a8u+/UL36bA4d6ghUBZbxzjs5eO45pyvzTgpjEfEOf/4JPXuawxqSU7w4PPoodOoEpUq5\ns7I0Z8cOqFfvew4f7gDUxrIWM3lyNm1vmYoUxiKStl28CG++aZ77RkUlfq1kSWjd2rzVqKFpvzfg\nhx/M/66jR2sAI8iQ4Wlmzsyi7S1TmcJYRNKubdvMsYQ7drjafH3hmWfMmuAyZdLN1pQpYepU6NZt\nGrGx9wNFyJx5APPmQZMmTlfm/fRrooikTevXQ+3aiYO4enXYvNmcolS2rIL4BsXGQv/+0KXLOGJj\nOwNvkzcvrFihIHYX9YxFJO2ZPds8+40fls6SxcyK7t5dQ9E3KSICQkPh669HAf2BhyhXbgyLF5vH\n7OIeCmMRSTtiY82hDS++6GorUMCsHa5c2bm60qjISGjbFpYuHQ4MBEJp3nwaM2f6acm1mymMRcTz\nHTsGn34KH32UePvK0qVh2TJ14W5BdLQ5aGrp0ghgNtCJfv0+5a23MuDr63R16Y/CWEQ804UL5qHl\n3LkwZ47pxiVUuzYsXGh20JKbEhsLnTvbzJ8fA2QC1vHKK9kYPlwp7BSFsYh4jqgoCA83pygtX24e\naF4pVy5zaMPgwdo3+hZERUH37jazZvUFdgJf0rdvTt54w+nK0jeFsYg4z7Zh8WLo2zfx7OiEqlaF\nXr3MbKNMmdxbn5f4/Xfo3DmO//3vGeAD4Bl69szAqFGaeO40hbGIOGvrVnOK0ooVSV8rW9bsQNG2\nLVSp4v7avERMDLz9NgwZEkdMTA9gMtCPLl1GMn68pSD2AApjEXHGuXPw2mswbpx5iBkvRw7TQw4N\n1baVKeDYMWjeHDZuBHgemIyv76uMHDmM556ztBLMQyiMRcT9vvoKnn4a9u93tfn4mHXCQ4fqFKUU\nEhlpNigzQQzQmaJFA1mxoj+lSztZmVxJvxOJiPv8+adrr+iEQXz//fDbb/DhhwriFGLb8NRT8N13\n0cBcLAtGjKjC7t0KYk90Q2FsWVZTy7L+sixrp2VZSY6UtiyrmGVZayzL+tWyrN8ty2qW8qWKSJq1\ncSO0aQPlyplecbx8+eCLL2DVKqhQwbn6vNCoUTBlSiTwMNCOp5/exMsvQwaNh3qk6/61WJbli5l2\n1wjYD/xsWdZC27a3JbhtIDDbtu0PLcsqBywBiqdCvSLiyWJi4KefYN8+OHAADh40e0V/+23Se598\n0py0pHXCKe6rr2DAgAigLbCUGjU+4L33gpwuS67hRn5HqgbstG17N4BlWbOAVkDCMLaB7JfezwEc\nTMkiRcTDXbgAn31m9ofeu/fa97ZoAa++ao40lBT3yy/QseMFzI/pVZQsOYl1657QjGkPdyNhHAj8\nm+B6P1D9inuGAMsty+oDZAFCUqQ6EfFcZ87AX3/BN9+YGdFHj179Xl9f6NjRHA10zz3uqzGd+ftv\naNoULlxYC6wlb97P+OGHx/D3d7oyuZ6UenrQAZhi2/YYy7JqAtMsy7rHtu24hDdZltUd6A5QrFix\nFPrUIuI2CxbABx+Yc4QPXmUALE8eMyErMNC8FS4M9epB0aLurTWdOXgQGjWyOXrUApqRLdtfrFt3\nJ3nzOl2Z3IgbCeMDQMJ/RUUutSXUDWgKYNv2j5ZlBQB5gSMJb7JteyIwESAoKMi+xZpFxN2iosxJ\nSePGXf2eYsWgXz/o2tUcaShuc/IkhIScZN++1sBAMmVqxLJld1KunNOVyY26kTD+GShpWVYJTAiH\nAR2vuOcfIBiYYllWWSAAuMaYlYikGQcPQrt28MMPidv9/aFkSShTBlq1grAw8PNzpsZ07OxZaNLk\nONu3NwL+wMcngi+/hFq1nK5MbsZ1w9i27RjLsp4GvgF8gU9t295qWdZQYJNt2wuBvsAky7Kex0zm\n6mLbtnq+ImndunVmJ6zDh11tDz1kZkHfeafWyTjszz+hZcsj/P13CLAD+IrPP3+ABx5wujK5WTf0\nL8m27SWY5UoJ215L8P42oHbKliYijrFtGDsWBgxwbVXp42NCuF8/nSrgAebOhS5dTnL+fANgL7CI\nd98N4ZFHnK1Lbo1+rRWRxM6eNc995851teXLB7NmQcOGztUlAERHw8svw5gxYFaShuDv35bJk+vT\nqZPDxcktUxiLiMu2beaEpD//dLXVqAFz5kCRIs7VJQBs2QKPPQa//roP80SwOHfdNY4vv4R773W6\nOrkd2ptaRCAuDt55xxxTmDCIe/c2z40VxI6KiYE33zR/Pb/+uhuoB7ShRQubTZsUxN5APWOR9G7X\nLnj88cRbVmbKBJMmoQeQztu710xU/+knMJO0GgIRPPfcfMaM0RGI3kJhLJJe2TZ89JGZkHXhgqv9\n3nth2jQd3OABtm2DRo3i91fZBgSTIUMsc+asoXXrig5XJylJYSySHh0/biZpLVzoavP1hVdegYED\n0f6Jztu0yWxtefy4ubas/mTJYvP992upWFG7eXgbhbFIerN2LTz6qDlVKV758jB1qnkoKY5buxZa\ntjQT28FsaDZ9+jTKlTtGyZIlHa1NUofCWCQ9OHXKTMRavBgmTzZD1PGeeQbefhsCApyrTy5bssQc\n/RwZCbARf/8xLF06lbp1cwG5HK5OUovCWMQbXbgA338Pq1fDqlXmTOG4uMT35MkDU6aYIw3FI3z9\ntVlZFh0N8AOW1ZR8+fJy553HMQfoibdSGIt4i3/+gZkzYelS+PFHc7jD1TRsaCZpFS7svvrkmubP\nNzuPmiBej2U14447CvPtt6sJDFQQezuFsUhaduaM2Slr2jQzDH21LeEtyzwPDg4203Pvvx+tifEc\nc+dChw5mPTGsxrJacNddxVm/fhWFChVyujxxA4WxSFoTEwPLl5sAXrAALl5M/r5y5UwPODgY6teH\nXHre6ImmTze7asVvAV6sWF4CA6uzYEE4+fPnd7Y4cRuFsUhasW+fOU/4iy/gyJGkr/v4mOB95BFo\n0gQKFnR/jXLDbBtGjTJncRjbKVWqDGvWVKRQodVYOowjXVEYi3i6P/6AkSPN82AzjplYxYrQqRN0\n7KhnwGlEbCw8+yx88EF8y5dAGE8+OYHChZ8EFMTpjcJYxBNFRsKiRfDpp2aty5UKFzbh26mTCWNJ\nMy5cMH91X30V3zILeJRq1arz5JPtHaxMnKQwFvEUtm1mQU+dCrNnm7XBV2rQAF580QxD+/q6vUS5\nPXv2mKVLv/4a3zINy+pC7SxTf7IAACAASURBVNp1WLJkEdmyZXOyPHGQwljEaXFxppv09tvxpwEk\nZlnQurV5uFi9uvvrkxSxZInZ+OzkyfiWffj4dKVBgwYsXLiQLFmyOFmeOExhLOKUqCgzGWvUqMTH\nFsYrUcL89O7cGe6+2/31SYqIjYWhQ2HYMNfKMz8/GDfuDu66awl16tQhU6ZMzhYpjlMYi7jb2bMw\ncaI5Pzjh/tBgDmjo1Am6dIHatU2vWNKsixfN+uEFC1xtOXO+z8CBJXjqqRZAI8dqE8+iMBZxlyNH\n4P33Yfz4pM+Ds2eHnj3NFFtt8uAVzp41TxdWr3a13X33KHbu7M8vv3QEtA2puCiMRVLbnj0wZgx8\n8knSDToKFIDnnjNBnCOHM/VJijt+HJo1g40bXW21ar3BDz8MIiwsjKlTpzpXnHgkhbFIatm+Hd54\nA8LDXdsrxbvrLujf3zwP1mlJXuXgQWjcGLZujW+xadhwMKtXD6NTp0589tln+GomvFxBYSyS0nbv\nhtdfN5OzrjwpqXJlMyu6bVstTfJC//5rtv3etctcW5Z5KrFt2wlKlOjGxx9/rCCWZCmMRVLKrl1m\nZvQnnyTdKSs4GF56yfypSVlead8+E8R79phrX1+b8eOP8NRTBYiLGweAjw7nkKvQd4bI7YiNNYfQ\nPvCAWX708ceJg7hpU/PgcOVKCAlREHupPXvMWRzxQZwhQxyNG/fmjTeqcPToUXx8fBTEck367hC5\nFbGxMGmSefbbsiUsW5b49Xr1YP16c7Zw1arO1ChusWuX2Rht3z5z7ecXS8OG3Vm69EMeeeQR8ubN\n62h9kjYojEVu1po15tlv9+6un8Bger3Nm5vjDdeuhbp1HStRUl9EhJmfd++98M8/ps3fP4Z69R5n\n+fJPGDRoEG+99ZZOX5IbomfGIjdqyxZ47bXEOzgA5MkD3brBU0+ZXbPEq9m22Tq8f39XCIOZFB8a\n+hZTp05j2LBhDBw40LkiJc1RGItcS1QUzJsHEybAt98mfi1zZjMpq29f8754vehoaNcu4YlLxj33\nmE3Vypd/huDgO+jUqZMzBUqapTAWSU5kJIwbZzbrOHw46eudOsGbb0JgoPtrE8c891ziIM6bFwYP\njuS//4ZTseIAsmTJriCWW6JnxiIJ2bbpCZcrZ8YhEwZxhgzQvr2ZHf355wridGb8eDNAEq9HD9iy\nJYJFi1oxfPgwVq5c6VxxkuapZywS7+efTQCvXZu4PTDQ/OR94gntG51OLV9utg2PFxYGo0efp1Wr\nlqxZs4bJkyfTqlUr5wqUNE9hLOmbbZs1wG+9lXhHf4BcucxOWk89Zc68k3Rp+3bznDh+M7Vq1eC9\n987SrFlzvv/+e6ZOnaqhabltCmNJf86fh19+MT3h6dPN+wn5+kLv3jB4MOTO7UyN4riICDMs/eab\ncOaMaStSxEymP3PmMLt372bGjBmEhoY6W6h4BYWxpA//+x9Mm2bGG7duTbpnNJgQDg2FgQOhbFn3\n1ygeISrK7Gj6xhvm0Id4mTNDePh5ChbMjGXdzY4dO8isWfSSQhTG4r0OHTKHNUybZtYIX01AgHke\n3LcvFC/utvLE8+zdCy1aJDxxybjjDhg79hi9ezeiWbNmDB8+XEEsKUphLN7Fts0ErAkTYP78pEcX\ngtkpq1w5s01ltWrmBKX8+d1eqniWX381ZxD/95+rrVAhM1Dy4IOHadYshJ07d/Lmm286V6R4LYWx\npG3nz8OOHfDnn2amzdy55s8rZcoErVvDo4+afaOzZnV/reKxvvkGHn4Yzp0z1/7+MGwYPP00nDp1\nkODgYPbt28eiRYsIDg52tljxSgpjSXtiYkzojhkDmzZd+966daFrV9P7zZbNPfVJmjJ1qnlKEX/Y\nVs6cZpJW/foQHR1NSEgI+/fvZ9myZdSrV8/ZYsVrKYwl7YiIgClTzJnB8WfVJSdrVujcGXr2NPsU\nilxFeDh06eK6LlrUHLRVvry59vPzY+jQoRQuXJhatWo5UqOkDwpjSRtWrYLHHoMDBxK3+/qac4TL\nlIHSpaFCBXOkYfbsztQpaca335rf2eLdey8sWQKFC8OuXbvYunUrLVu25OGHH3auSEk3FMbi2aKj\nzUlJb79tJmfFy50b+vQxD/V0XqzcpO3boVUrs4wJzO9yq1ebb6u//vqLhg0bEhcXR3BwMFmyZHG2\nWEkXFMbiuXbvhg4dzF7Q8fLlg1deMQ/5NAlLbsF//5lZ0ydPmusCBczQdO7csG3bNho2bIht26xc\nuVJBLG6jMBbPNHOm2YYyfusjgMaNzWybggWdq0vStB07zL4ue/ea68yZYfFis7z8999/JyQkBF9f\nX9asWUNZbfwibqRTm8SznDsHjz8OHTu6gjhDBjNpa+lSBbHckshIs814hQrw22+mzccHZs+GKlXM\n9bx58/D392fdunUKYnE7y074HM6NgoKC7E3XW5Yi6csvv5jjcP7+29V2112ml1y1qnN1SZq2dq05\ndGvHDlebjw98/HH8kqYYMmTIgG3bHDt2jHz58jlWq3g3y7I227YdlNxr6hmL87ZuNTOlq1dPHMSd\nOpltkRTEcotmz4bg4MRBXK0abN5sgvj777+nXLly/PXXX1iWpSAWxyiMxRkXL8KaNWYZ0j33wOef\nu3ZdyJrVXH/+uTbqkFu2YoXZcC3+TJDs2WH8ePjhB7jvPli7di1NmjTBsiyyajKgOEwTuMQ9Ll6E\nOXPgu+/M0YVbtrjCN6GGDc344d13u79G8RobN8JDD5mVcWCWLq1cCYGB5nrlypW0bNmSEiVKsHLl\nSgoVKuRcsSIojCW1RUebXbOGDYN//03+Hssy+0YPGGCGqkVuw/btZunS+fPmumhRc3JmfBD/8MMP\ntGjRglKlSrFy5Ury65AQ8QAappbUERsL06eb05G6d08+iEuXNsuXtm2DefMUxHLb9u0zK+COHzfX\nefKYIC5a1HXPvffeS9euXVmzZo2CWDyGesaSsmzb7LI/aFDSQ2Hz54fevaF2bbOeJGdOZ2oUr3Tw\noJmstX+/uc6SxWxvWaaMuV6xYgU1atQgW7ZsTJgwwblCRZKhnrGkDNs259BVqwZt2iQO4ly54M03\nzY5ar71mfmIqiCUFHTlivq127TLX/v5msKVaNXM9c+ZMHnjgAQYNGuRckSLXoJ6x3L69e80e0YsX\nJ27PmhWeew769lX4Sqo5cQIaNTJHWoPZI2bOHDNcDTB16lS6du1K3bp1eeONN5wrVOQaFMZy66Kj\n4d13YcgQuHDB1R4QYIajBwwwe0mLpJKzZ6FJE/j9d3Pt42OmKrRsaa4nT55M9+7dCQ4O5quvviJz\n5szOFStyDQpjuTWrV5te75YtrjbLgiefhMGDzTl0IqnsxRchfiM/y4LPPoP27c31uXPneP3112na\ntCnz5s0jICDAuUJFrkNhLDfnxx/h1VfNhh0JVawIEydqRrS4zfr1Zkl6vA8+cJ1PbNs2WbNm5dtv\nv6VQoUJkzJjRmSJFbpAmcMmN2bIFmjeHWrUSB3HmzOYQh02bFMTiNhcvmkGYeC1bmlVyAG+99RbP\nPvsstm1TvHhxBbGkCQpjubYjR8xPufvuM+tE4vn6mp+Gf/4J/fqBn59zNUq6M2yYa7/pbNlMrxhs\nhg4dyssvv8yxY8eIjY11skSRm6JhaknemTMwaRIMHZr4TGHLMscbDhmiLSvFEf/7H4wc6boeORIC\nA20GDhzIiBEjeOyxx/jkk0/w9fV1rkiRm6QwFiM62jyEW7XKTM7atMnsopVQ48YwZow52EHEATEx\n8ccemuu6dc0Gb6+++ipvvvkmTz75JB999BE+Phr0k7RFYZzeHThgesATJ8KhQ8nfU6oUjB1rNvy1\nLPfWJ3JJXBw884xr9rS/v/nW9fGBmjVr8uyzzzJ27FgFsaRJCuP06pdfYMQIs3Vlcs/WLMs8J+7S\nBXr21DNhcVRcnPk2nDjR1TZwYBynTv0MVOfBBx/kwQcfdKw+kdulME5voqLM7Jc330wawoUKmXPn\ngoOhQQPInduREkUSiouDHj1g8mRXW2hoLLt2PUmtWlP59ddfqVixonMFiqQAhXF68ttv8Nhjru2K\n4tWvb3bMat1aPWDxKDEx5pnwZ5+52jp2jAG6MHXqdAYPHkyFChUcq08kpSiMvd1//5l1wStWwLRp\nrpkvAPXqwfjxoB9m4oE2bTJB/OuvrrZOnaKJiHiEuXPnMHz4cF555RXnChRJQQpjbxQRYfaM/uIL\nc1bwlTJlMsPUffqY2S8iHuTMGRg40KwdjotztXfrBvffP5dHH53D6NGj6du3r3NFiqQwhbE3sW2Y\nPRv694d//kn+npo1YcoUM0NaxMP88AO0a2fOJo4XEGC2O+/fHywrjOLFi1G7dm3nihRJBeoWeYvN\nm6FOHQgLSxzE/v7mmfDQofD99/Dddwpi8UjLlkFISOIgbtIENm2K4PffO7Jt2x9YlqUgFq+knnFa\nFxsLb78Nr72WeHZ03rxm1nTnzmb/aBEPNmcOPPKI2XsGzLfv++9DixbnadnyQdauXUvz5s25RxvO\niJe6oZ6xZVlNLcv6y7KsnZZlvXSVe9pblrXNsqytlmXNSNkyJVkHDphT1V991RXEfn7Qty/s3Gn2\nlFYQi4f75BMzoBMfxMWKmeHq5s3P0qzZA6xbt47PP/+cRx55xNlCRVLRdXvGlmX5Ah8AjYD9wM+W\nZS20bXtbgntKAi8DtW3bPmlZVv7UKlgwz4bnzTOLL48fd7XXrAlTp0LJks7VJnIT3nvPHIsdr0wZ\nM/E/W7bTNGnyABs3bmTGjBmEhoY6V6SIG9xIz7gasNO27d22bUcBs4BWV9zzJPCBbdsnAWzbPpKy\nZQpgliVNnw733gsPP+wKYsuCQYPM3tIKYkkjxoxJHMSVK5tv4SJFwN/fn1y5cjFnzhwFsaQLN/LM\nOBD4N8H1fuDKg2tLAViW9T3gCwyxbXtZilQoZn3H5MlmOdLevYlfCww0S5gaNHCiMpFb8vbb8FKC\nB161a8PixRAdfYyTJ33JlSsXixYtwtJe6JJOpNQErgxASaABUARYb1lWBdu2TyW8ybKs7kB3gGLF\niqXQp/Zyu3bB44/Dt98mbs+SxZwnPHAg5MnjTG0it2D4cPNtG69ePRPE588fJjg4mLx587JmzRoF\nsaQrNzJMfQAomuC6yKW2hPYDC23bjrZtew+wAxPOidi2PdG27SDbtoPy5ct3qzWnD3FxMGGCGZJO\nGMR58phlSv/8A++8oyCWNGXo0MRBfP/9sGQJnDlzkAYNGrBnzx5ee+01BbGkOzcSxj8DJS3LKmFZ\nlj8QBiy84p4FmF4xlmXlxQxb707BOtOXXbvM2cG9e8P586bN19fMmt63zzwf1iEOkobYtll9N3iw\nqy0kBBYtghMn/qV+/frs37+fZcuW0bBhQ+cKFXHIdYepbduOsSzraeAbzPPgT23b3mpZ1lBgk23b\nCy+91tiyrG1ALPCibdvHr/5RJVlRUTBqFLzxBly86GovW9bMkq5a1bnaRG6RbZve8IgRrrbGjc3p\nnZkyQYsWj3HkyBGWL19OzZo1nStUxEGWbduOfOKgoCB7U/wp4eldXJw5zKFPH9i+3dXu42PWDA8d\navYEFEljbBteftlM2Ir3wANmZV78t/Tu3bs5ceIEQUFBzhQp4iaWZW22bTvZb3Rth+mUqCj45htz\nYnrRombMLmEQV64MP/0EI0cqiCXNsW2zTOmBBxIHcfPmMH8+7N37Jy+99BJxcXHceeedCmJJ97Qd\nprv9/Td89JE5oPXkyaSvZ81qhql794YM+uuRtCUuDr7+Gt56CzZsSPxaq1YQHg5///0HISEh2LZN\n7969KVq0aPIfTCQd0U97d7Bts3bj/fdh+fLk78mTB9q0MbNcihRxb30iKeDcOejQwUzKSsiyoGtX\nszhg+/b/ERISgp+fH6tXr1YQi1yiME5tf/8NvXrBypVJXytWDNq2hdatoVYt9YQlzTp8GFq0gITT\nQPz94bHHoF8/c1DY5s2badSoEVmyZGH16tWU1G5xIpfpp39qiYw0z3uHDzfvx7Ms81OrVy8zpdRH\nj+0lbduxA5o2hT17XG29e5uVeIUKudqOHz9OgQIFWLJkCSVKlHB/oSIeTGGcki5cMGcGr14Nc+ea\nk5Pi+fiYAO7bF4oXd6xEkZS0ahWEhrq2SffxgQ8/hO7dXfccO3aMvHnz0rhxY7Zs2UIGjQCJJKFu\nWUrYvNn0cnPlMn++9VbiIA4Kgp9/Ns+MFcTiBXbtMk9YQkJcQZw5M3z1VeIgXrt2LXfeeSdffvkl\ngIJY5CoUxrcjLs4cPVOzpjn3LSoq8evZssG4cWZaaeXKztQokoLOnoUBA6BcObNWOF6+fLB2rXkC\nE2/FihU0a9aMokWLUrt2bbfXKpKW6NfUW3X4MHTpAsuuOJyqXDkIDoaGDc1b9uyOlCeS0iIizOFg\nv/ySuP3RR81gUGCgq23JkiW0adOGMmXKsGLFCrQXvci1KYxv1qlTZmvKN980gRwvKMgcZVi6tHO1\niaSi559PHMQ1asC770L1Kw5U/fvvv2ndujUVK1Zk+fLl5NY+6iLXpTC+Ub/+ahZKzphhJmol1L8/\nDBtm1nKIeKHZs+Hjj13XI0eaJUvJHa5UsmRJPvzwQ9q2bUvOnDndV6RIGqa9qa/l4kUzK3rCBPjx\nx6SvFyxoesmNG7u/NhE32b0bKlWCM2fMdfv2MGtW0iCeNWsWpUqVorLmR4gk61p7U6tnnJy//4ZP\nP4XJk+HYsaSvV6hgFlI+8ojZvlLES0VFmaVL8UF8550wcWLSIJ4yZQpdu3alTZs2zJ071/2FiqRx\nCuN4x46ZjXOnTTMHNFzJzw8eftisFa5dO/nxOREv07+/a1ctPz/TI86RI/E9EydOpEePHjRq1IjP\nP//c/UWKeIH0HcbxR8t88IE5SiYmJuk9RYvCU09Bt25QoID7axRxyLvvwnvvua5Hjkx6pPb48ePp\n06cPzZs3Z+7cuQTohDGRW5I+w/jCBXNq0oQJsG1b0tf9/KBZM3j8cXPmmzYqkHRm2jQzezpe69bw\n7LOJ74mLi+Obb76hdevWhIeH468JjCK3LP2lzNat5nSkHTuSvlajBnTqZB6S5cnj/tpEPMDixeb3\n0Hi1a8P06YmfzERERJApUybmzJmDr68vfn5+7i9UxIukrx24Zs6EatUSB3HWrNCzJ/z+u5kx3auX\ngljSre+/h3btIDbWXFeoYM4nzpzZXNu2zZAhQ6hVqxanT58mICBAQSySAtJHGEdFmTG2jh1da4Qz\nZzYPxQ4cMMPVFSo4W6OIw/75B1q2NDttgdlGfdkys+U6mCB+9dVXef3117nvvvvIqpUEIinG+4ep\nL16EBx9MfJ5wqVJmY93y5Z2rS8SDREdDhw5w4oS5zp8fli+HwoXNtW3b9OvXj7Fjx9KjRw8mTJiA\nj47/FEkx3v2vKSbG9IYTBvFDD5kTlBTEIpcNHgw//GDe9/WFL7+EkiVdr48YMYKxY8fSp08fPvzw\nQwWxSArz3p6xbZslSfPnu9peew2GDNEaYZEEli83W63HGzYM6tRJfE/nzp3x9fVlwIABWPr3I5Li\nvPfX21degU8+cV2/8IKCWOQKhw6ZBQTxGjUyRyQCxMbGMnnyZGJjYylatCgvvfSSglgklXhnGI8b\nZ850i9e5M4wapSAWSSA21hx/eOSIuS5Y0Kwv9vGBmJgYOnfuzJNPPsmiRYucLVQkHfC+YerISBg0\nyHXdooXZY1rPuEQSGTECVq8271uWWUtcoABER0fzyCOPMGfOHEaMGEGrVq2cLVQkHfC+MD550rWr\nfY4c5uw3rYMUSWT9evPUJt7AgdCwIURGRhIWFsaCBQsYM2YML7zwgmM1iqQn3hfG0dGu97NmhUyZ\nnKtFxAMdPWqWMcXFmet69czcRoBt27axfPly3n//fZ5++mnnihRJZ7w7jNUjFkkkLg66dIGDB811\nnjwwYwZYVizgS6VKldi5cyeFChVyskyRdMf7HqQqjEWu6p13YMkS1/Xnn0OOHOdo1KgRH3/8MYCC\nWMQBCmORdCA62qz2e/FFV1u/flCnzhmaNm3KunXryJYtm3MFiqRzGqYW8XK7d5uN6H76ydVWvTq8\n+OIpGjduyubNm5k1axbt2rVzrkiRdE49YxEvFh4OlSolDuJGjWDu3EiaNQvhl19+Ye7cuQpiEYcp\njEW81McfQ1iYa6VfhgwwcqQ5ialIkYx06NCBBQsWaB2xiAfQMLWIF5o61WzNHu+uu8xx3kWL/sfm\nzf9StWpV+vbt61yBIpKIwljEy4SHQ9euruuqVWHFCjh37gD16zfk/Pnz7Ny5k4CAAOeKFJFEFMYi\nXuSrr+CRR1wbetx7rxmWPn36Hxo2bMiRI0dYunSpgljEw+iZsYiXmDUL2rUzB0AAlC1resSnT++h\nXr16HDt2jBUrVlC7dm1nCxWRJBTGImmcbZtDyTp0cH373303rFoF+fLB6NGjOXv2LKtWraJ69erO\nFisiyVIYi6RhsbHwzDPQv7+rrWxZcxpT/EZa77zzDhs2bKBKlSrOFCki16UwFkmjTp2Chx+G8eNd\nbfXqwfffw+nTfxAcHMzRo0fx9/enZMmSzhUqItelMBZJg5Ytg3vugQULXG3t28M338C+fb/RoEED\n/vzzT06ePOlckSJywxTGImnI2bPQvTs88AAcOOBqf+EFs474jz820bBhQzJnzsy6desoVaqUc8WK\nyA3T0iaRNGL7dmjWDPbudbXlywcffght28KmTZsIDg4md+7crFmzhuLFiztVqojcJPWMRdKArVuh\nQYPEQdy2Lfzxh/kToEiRItStW5f169criEXSGIWxiIfbssUE8ZEj5jprVpgxA+bMgfz54bfffiMm\nJoaCBQuyaNEiihYt6mi9InLzFMYiHux//4P774djx8x1tmxmklaHDmBZ8M0331CzZk0GDx7sbKEi\nclsUxiIe6vffoWFDOH7cXGfPDsuXQ61a5nrRokW0bNmSMmXK8PzzzztXqIjcNoWxiAfavRuaNIET\nJ8x1jhxma8saNcz1/PnzadOmDRUrVmTVqlXkzZvXuWJF5LZpNrWIh/nvP2jc2PwJpke8ciUEBZnr\n06dP061bN6pUqcKyZcvIkSOHc8WKSIpQGIt4kNOnoWlT2LXLXAcEwNdfu4IYIEeOHCxfvpzSpUuT\nLVs2ZwoVkRSlYWoRDxERAS1bmklbAL6+5mzievXM9aeffsp7770HQFBQkIJYxIsojEU8wPHj0KgR\nrF/vaps82YQzwEcffUS3bt1YunQpsfFnJIqI11AYizhszx6oXdsc8BBv1Cjo0sW8P27cOHr27Enz\n5s1ZsGABvr6+jtQpIqlHYSzioM2boWZN+Osvc21ZMHYs9OtnrkePHs2zzz7LQw89xLx58wgICHCu\nWBFJNZrAJeKQ1avNMPT58+Y6Y0aYNg3atXPdkzFjRkJDQ5k2bRp++n4W8VrqGYs44Ntv4cEHXUGc\nK5dZR9yuHdi2zd5Lm1D36dOHmTNnKohFvJzCWMTNNm6E5s3hwgVzHRhonhfXrWuC+OWXX+aee+7h\nr0tj15ZlOVitiLiDwljEjX77zeysdfasuS5Y0AxXly1rgrhv3768/fbbdOrUiZIlSzpbrIi4jcJY\nxE22bzfLl06dMtd585qdtUqVgri4OPr06cM777zDM888w4QJE/Dx8b5/niKSPO/7164wFg909Cg0\na+Y6fSlnTnPoQ/ny5nrKlCl88MEH9OvXj3fffVdD0yLpjGZTi6SyyEho0wYuzckiSxZYtgwqVXLd\n07lzZzJnzkxoaKiCWCQdUs9YJBXZNjz1FHz3nbm2LJg5E6pXh5iYGPr378+hQ4fIkCEDYWFhCmKR\ndEphLJKKRo+GKVNc1yNHmiVNUVFRhIWFMWrUKBYvXuxYfSLiGTRMLZJKFiyAAQNc1127Qt++EBkZ\nSfv27Vm4cCFjx47liSeecK5IEfEICmORVLB+PYSFmWFqMCcvffghXLwYQdu2bVm6dCkffPABvXr1\ncrZQEfEIGqYWSWG//WaGoiMjzfVdd8GXX4K/P1y4cIH9+/czadIkBbGIXKaesUgK2rULmjaFM2fM\ndcGCZglTQMA5oqL8yZMnD5s2bcLf39/ZQkXEo6hnLJJCDh2Cxo3h8GFznSMHfPMN5MlzmiZNmvDo\no49i27aCWESSUBiLpIB//oEGDWD3bnMdEACLFkHRoidp3LgxGzdupH379lq6JCLJ0jC1yG366y+z\nzeW//5prX1+YOxfKlj1OcHAj/vjjD+bOnUurVq2cLVREPNYN9Ywty2pqWdZflmXttCzrpWvc19ay\nLNuyrKCUK/EmKYzFjX75xZy2FB/E/v4wezY0a2bTpk0btm3bxldffaUgFpFrum7P2LIsX+ADoBGw\nH/jZsqyFtm1vu+K+bMCzwE+pUegNUxiLm3z/vdlvOn6yVubMZm1xo0YAFiNHjuTs2bOEhIQ4WaaI\npAE30jOuBuy0bXu3bdtRwCwguV/zhwFvAxdTsL6bExvrWthpWWa8UCQV7Nhhli/FB3HOnOYEpnLl\nDjBx4kQAqlevriAWkRtyI2EcCPyb4Hr/pbbLLMuqDBS1bdvZff2iolzvq1csqeTkSRPEJ0+a6/z5\nYd06CAz8h/r169OvXz8OHTrkbJEikqbc9gQuy7J8gLFAlxu4tzvQHaBYsWK3+6mT0hC1pLLoaGjX\nzvSMwcyaXrwYsmbdTb16DTl16hQrVqygUKFCzhYqImnKjfSMDwBFE1wXudQWLxtwD7DWsqy9QA1g\nYXKTuGzbnmjbdpBt20H58uW79aqvRmEsqci24ZlnYNUqV9vnn0OOHH9Tv359zp49y+rVq6levbpz\nRYpImnQjYfwzUNKyrBKWZfkDYcDC+Bdt2z5t23Ze27aL27ZdHNgAtLRte1OqVHwtCmNJRWPHwkcf\nua6HDjW95A0bNhAViMGxQgAAG29JREFUFcWaNWuoXLmycwWKSJp13TC2bTsGeBr4BtgOzLZte6tl\nWUMty2qZ2gXeFIWxpIK4OHj5ZejXz9XWoQO8+KLZfLpTp07s2LGDihUrOlShiKR1N7TO2LbtJbZt\nl7Jt+y7btodfanvNtu2FydzbwJFeMSiMJcVFRkKnTvDWW662WrWgT59fKVWqJOvWrQMgR44cDlUo\nIt7Au7bDVBhLCjp1Ch54AGbMcLW1aAEjRmykWbOGWJZFkSJFnCtQRLyGd22HqTCWFHLhgtm8Y1OC\nMZ4ePaBjxx9o0aIpefPmZc2aNdxxxx3OFSkiXkNhLHKFuDh47LHEQTxiBLRqtZ1q1RpTuHBhVq9e\nrV6xiKQYDVOLXGHwYHPQQ7z33zcTuMqUKc3zzz/PunXrFMQikqK8N4x1Zqzcghkz4I03XNd9+kCZ\nMiv5559/8PHxYdiwYdrQQ0RSnPeGsXrGcpM2bICuXV3XTZpAw4Zf07x5c1544QXnChMRr6cwFgG+\n/dacwBRplg5Ttix07Pgl7dq14d5772XSpEnOFigiXk1hLOnenDlm5nT8wQ958kD37jPp2jWUatWq\nsWLFCnLlyuVskSLi1RTGkq698w6Ehrp6xPnzw5IlscyY8Q61a9dm2bJl2tBDRFKdljZJuhQdDX37\nmpnS8UqXhsWL47jrLl+WLl1KQEAAWbJkca5IEUk31DOWdOfQIQgOThzEtWtDt24TePbZlkRGRpIn\nTx4FsYi4jcJY0pVvv4XKlc2f8dq2hVat3qV//974+HjXPwkRSRs0TC3pgm3DuHFmaDo21rT5+Jg1\nxT4+I+nffwBt27ZlxowZ+GuNuoi4mXd1AxTGkozISHjiCXjuOVcQ580Ly5eDv/8YXnppAGFhYcya\nNUtBLCKOUM9YvNrRo9CmDXz3nautWjWz3WXRopAr1/307t2b9957D19fX+cKFZF0TT1j8VpbtkDV\nqomDuHNnWLfO5s8/VwBQuXJlxo8fryAWEUcpjMUrzZ4NNWvCvn3m+v/t3X+cjXXex/HXxzCspfzM\nbSkkkojcQm4bZvzKz70li7RItyJZbKLftkj2fpRtd0PK3EJkExqWGb8GtcWqtG201pAkNon1m/n1\nvf+4ZnaGjDnMmXOdc+b9fDzm8TjXda7HOW9fY96+3+vMdZnBb34D//d/jgkTxtCpUyeSkpL8DSki\nkk3L1BJV0tPh0Ufht7/N3VeunHcDiG7dshg5ciQzZsxg9OjRdO7c2b+gIiJ5qIwlahw8CH37nr8s\nfcMNsHQpNGyYxbBhw5g9ezaPPvooL7zwAmbmX1gRkTy0TC0R74svYNQoaNDg/CLu1Qs++ggaNYIP\nPviAhIQEnnrqKRWxiIQdzYwlYiUleeeBU1LO31+iBEye7C1X51zDo02bNnzyySc0bdo09EFFRAqg\nmbFEpNdegzvv/GER33ij9/vDEyZARkYaAwYM+PcHtVTEIhKuVMYScZYsgQcfzN2OifF+l3jtWm/J\nOj4ezp07R58+fVi4cCG7du3yL6yISAC0TC0RJSUF+veHrCxvu1kzWLbMu4BHjjNnztC7d2+SkpKY\nPn06w4cP9yesiEiAVMYSMT75xPtQVlqat12/Pqxa5d2DOMfZs2fp0aMH69ev5/XXX2fo0KH+hBUR\nuQwqY4kI+/ZBly5w4oS3/ZOfeOeG8xYxQGxsLPXr12fQoEHce++9oQ8qInIFVMYS9pzzbvTw3Xfe\ndsWKXhHXqpV7zLFjxzh69Ci1a9dm+vTp/gQVEblCKmMJe7NnwxrvUtKYeeeIb7459/mjR4/SuXNn\njh49yvbt23XnJRGJOCpjCWv79sHYsbnbY8fCHXfkbh8+fJhOnTqxfft2Fi9erCIWkYikMpaw5RwM\nG5Z7nrh+fXjuudznDx06RIcOHdi1axfvvvsuXbp08SeoiEghqYwlbCUkQHKy99jM2/7Rj3KfHzdu\nHKmpqaxYsYL4+Hh/QoqIBIEu+iFhaf/+85enR4+G//qv8495+eWXSUlJURGLSMRTGUtYGj8ejh/3\nHterB5MmeY/37t3Lfffdx5kzZ6hQoQItW7b0L6SISJCojCXsbNni3X84x6xZULYs7N69m7Zt27J0\n6VL27NnjX0ARkSBTGUtYcQ7GjMnd7t0b2rWDnTt3cscdd3Dq1CnWr1/PzXl/t0lEJMLpA1wSVhYt\ngg8/9B7Hxnq3SNyxYwdxcXE450hJSaFx48b+hhQRCTLNjCVsnDnjnSvOMWoU1K0LWVlZVKtWjQ0b\nNqiIRSQqaWYsYWPaNO8iHwBVqsCQIftxrgaNGjVi27ZtlCgRXf93FBHJEV0/3VTGEevAAZgyJXd7\nyJAttG7diJdeeglARSwiUS26fsKpjCPSiRPQsyecPOlt16nzZ2bO7EjlypXp06ePv+FEREJAZSy+\nOnsWfvYz+Phjb9tsAwcPdqZ69eps2rSJWnlvzSQiEqVUxuKbzEwYOBDWr8/Z8z2lS/fk+utrsXHj\nRmrUqOFnPBGRkNEHuMQXzsGIEfDOO7n7Jk+uTOPGb9KqVSuqVq3qXzgRkRBTGYsvnnrKu7KW5116\n9IDHHuuFWQ8/Y4mI+CJ6lqkzM73pFni3+ImJ8TeP5Ovll2Hy5JytxZj14V//ehFwPqYSEfFP9JSx\nZsURYf587w5MnoWY9eP221uyYsUKzMzPaCIivlEZS8isXAlDhuRsvQEMpE2bNiQnJ3HVVVf5mExE\nxF/Rc85YZRzWUlKgTx/IyPC2K1f+lEaN4li58l3Kli3rbzgREZ+pjKXILVgAgwfn/BUdo3btq3n/\n/ZeoUiWN0qVL+5xORMR/WqaWIuOcd9ele+7J+euZRokSDZkzZx81apiKWEQkm8pYikRmJjz8cN67\nME0FxtKlS2tat67uYzIRkfCjMpagcw4eeABeeSVnz3PABO66qz/vvruQUvr7ERE5j8pYgm7SJJg9\nO2drNvA0AwcOYtGieZQsGT0fUxARCZbo+cmoMg4Lb7wBTz+du92/f1+aNj3CI4/8SrdBFBHJh8pY\ngmbtWrj/fvCupPUH2rUbwpw55YmNHedzMhGR8BY9UxWVsa/+9je46y7IyMgCRgCj6NlzPrGxficT\nEQl/mhlLoZ05A337wvHjmcAwIIGHHprA6NEP+B1NRCQiaGYshfb44/D3v2cAg4EEHnzwGX7/++d1\nrWkRkQBpZiyFkpICv/0twEFgLb16TWbGjMd9TiUiElmip4zT0nIfq4xD4tgxGDQoHe/b6Fo6dNjO\n0qWV/I4lIhJxtEwtV2zkyLN8/XVvYByVKsHcuZXQyrSIyOVTGcsVWbToDPPn9wJWAPWYOROq6yqX\nIiJXJHqWqVXGIbNz5ykGDuwBbAASGDBgCHff7XMoEZEIpjKWy5Ke7mjZshcZGRuBudSsOZA//MHv\nVCIikU3L1HJZJk40jh0bASwkJmYgCxdCxYp+pxIRiWyaGUtAjhw5wsyZf2HKlC5AbwCefRbatPE3\nl4hINFAZS4EOHz5M+/Yd2b59F859CVSlY0eYMMHvZCIi0UFlLJf07bff0r59PDt37sa5d4GqVKsG\n8+aBbsIkIhIcKmPJ14EDB4iLi2fXrn1kZf0JiMMM3nwTqlXzO52ISPRQGUu+5s1bSGrqfrKykoCf\nAt6lL+Pj/c0lIhJtAlpoNLMuZrbTzFLN7AdnCs1srJntMLPPzGydmdUKftQCqIyDxjlHejp88MFY\nMjM/I6eIX3gBRo3yN5uISDQqsIzNLAZ4BbgTaAj0N7OGFxy2DWjunLsFWAz8JthBC6QyDorU1FSa\nNGnBT3/6BYmJBtQB4JlnYPx4f7OJiESrQJapWwCpzrk9AGb2FtAL2JFzgHMuJc/xm4GBwQwZEJVx\noX322d9p0yaOEyfSgXP/3j9unFfGIiJSNAJZpq4BfJ1ne3/2vvwMBVYVJtQVyVvGsbEhf/tIt2DB\n5zRr1o4TJzKBFKApJUrAk0/C1KnoBhAiIkUoqB/gMrOBQHOgbT7PDwOGAVx33XXBfGvNjAth5sy/\nM3x4e6AUsB5owH/+J8yaBc2a+RxORKQYCGRm/A1wbZ7tmtn7zmNmHYAngJ7OuXMXPg/gnJvlnGvu\nnGtetWrVK8mbP5XxFVm4EEaOvA7oDGykfPkG/O53sGWLilhEJFQCmRlvBeqZWR28Eu4HDMh7gJnd\nCrwKdHHOHQp6ykCojC/bE09s4/nnrweuBuZTpw6sWQN16/qdTESkeClwZuycywBGAsnAF8AfnXPb\nzexZM+uZfdj/AuWAt83sUzNLLLLE+VEZX5YRI97n+efvAB4CoGFDeP99FbGIiB8COmfsnFsJrLxg\n39N5HncIcq7LpzIO2NixG5gxoxve2Yep3HYbrFoFlSv7nUxEpHiKnqsLq4wD8tRTa5g2rStQG9hA\n27Y1WLdORSwi4ieVcTGSlJTOpEkjgHrABm677T9YvhzKl/c7mYhI8aZrUxcTW7dCnz6lgCSgAjfe\nWJmVK1XEIiLhQDPjKJeVBUOHvk3r1qM4dcoBdalZszKrV0OVKn6nExER0Mw4qu3dC926vcmOHb8A\nWgNnqFSpLKtXQ7CvuSIiIldOM+MoNXcuNGgwhx077sW7INoqbr65LCkpcNNNfqcTEZG8VMZRaMYM\nGDTodc6dGwJ0wGwFEyaU4+OP4ZZb/E4nIiIX0jJ1lJkzB0aMAKgO/Df16i1g7twytGrlby4REcmf\nyjiKvPUW3HffF8BNQDdatOjGmjVw1VV+JxMRkUvRMnWUWLYMBgyYgnONgE00bQpJSSpiEZFIoDKO\nAps3O/r0+TXOPQ7056abWrN6NVSs6HcyEREJhJapI9zBg46OHZ8gM3MKMJi6dV9n3boYgn2HShER\nKTqaGUewtDTo0GENJ09OAR6gQoXZJCfHUL2638lERORyaGYcwX75S9ixoyOwHLNuLFpkugWiiEgE\n0sw4AmVlZdG58zhmzvwrYEB3pk41OnXyO5mIiFwJzYwjTGZmJl273s/q1XOAikAT+vWDRx7xOZiI\niFyx6CjjzExwzntsBjEx/uYpIhkZGXTvPojVqxcAE4HHaNIEZs/2/tgiIhKZoqOMi8GsOD09nZ49\n7yE5+W3geeAxataExEQoW9bvdCIiUhgq4wjx1VdZbNhwAngRGEulSujuSyIiUUJlHObOnj3L11+f\noVevipw9uwKI4cc/hpUrdfclEZFooTIOY6dPn6Zbt5/xl78c4/TpPwMlKVUKliyBli39TiciIsGi\nMg5TJ0+epGvXHrz33kYgASiJGcybh36FSUQkyqiMw9Dx48fp0qUrH374ITAfGABAQgL8/Oe+RhMR\nkSKgMg5D99//AJs3bwHeAu4GYPp0GDzYz1QiIlJUouMKXFFUxllZcObMCzi3jJwifuklGD7c31wi\nIlJ0NDMOE9999x0zZszg5MknWbGiFlALgEmTYMwYf7OJiEjRUhmHgX/+85/Ex8eza9ce0tN7A40A\neOghePxxf7OJiEjRUxn77JtvviEuLo6vvtpPRsZKcoq4Z094+WVd5lJEpDhQGfto3759xMXFceDA\nISAZ59oA0KIFLFwYtZfYFhGRC6iMfbR79x4OHTrF2bOrca4VAHXqwPLlut60iEhxojL2wcmTJ4mN\nLcf8+e04cWI34DXvtddCUhJcc42/+UREJLT0q00h9sUXX1C/fgMaN55HQgLkFPHtt8PWrVC/vq/x\nRETEByrjEPr8889p27Ydhw9n8I9/NPv3/l/8AlJSoFo1H8OJiIhvVMYh8umnn9KuXTuOHy9JevpG\n4GYApk6FOXOgdGlf44mIiI90zjgEvv32W+Li4jh3rhznzq0HbgBg2jQYPdrfbCIi4j/NjEPgmmuq\n0aTJs5w+vYmcIn78cRWxiIh4NDMuQu+99x5mZXjjjdvYsGHkv/cPG+Zd5lJERARUxkVm3bp19OjR\nkxIlmnLq1PuAdymt3r29OzDpyloiIpJDy9RFIDk5ma5du5OWdj2nTi0hp4gHDIA339SVtURE5Hwq\n4yBbtmwF3br1JC2tAZmZKUA1SpaE3/8e5s+HMmV8jSciImFIy9RBtHcvDB36BpmZtwDJQCV+8hNY\nvNi7qIeIiMjFaGYcJAsWZNCkCRw5Mh9YC1QiPh4++URFLCIil6YyDoLevedzzz0tOH78e6A0MTFX\nM2UKrF6tq2qJiEjBtExdSH36JLB06f1AO6AMderAggXQqlVIY4iISATTzLgQ+vZ9lXfeGQp0BFbQ\nqdOP2bZNRSwiIpcnOmbGaWm5j0NUxkOGzOXttx8EugGLad++DMuWwY9+FJK3FxGRKKKZ8RV47TWY\nMyceGAUsoXXrMiQmqohFROTKqIwvg3Pw858vY9iwTKAG8DLNm8eyciWUK1dkbysiIlFOZRygtDTH\nrbdO5I9//G8gAYBbb4XkZLj66iJ5SxERKSZUxgE4ftzRoMHj/PWvvwaGAPfRuTNs3AiVKgX97URE\npJhRGRdgxw5HnTq/4ssvXwAeBF5nyJAYli+H8uWD+lYiIlJMqYwvYfFiuO22VI4ceRXvw1rTmTix\nBLNn+34JbBERiSLR8atNQS7j9HQYP94xbZoB9YBPKVPmBl57zRg4sNAvLyIich6V8QXOnoXevTNZ\ntep+4DZgBHXr1mPJErjllkK9tIiIyEVpmTqPM2egR48MVq26F5gDfE+PHvDRRypiEREpOirjbKdP\nQ/fuaaxd2w9YCLzAY489xbJlUKFCUFKKiIhclJapgVOnoHv3LDZsuBtIBF5i4sQxPPNM0BKKiIjk\nq1jPjJ2DFSugWTPYsKEE0Bb4A889pyIWEZHQKbYz4+3bYcwYWLPmNPAPoCkwlilTYMKEoggpIiJy\nccVyZjx5sveBrDVrTgJdgfaUL3+UV19VEYuISOgVu5nx9Onw5JMAx/GKeDPx8fNYuLAiVasWYUYR\nEZF8FKuZ8fLl8PDDAEeBjpht4cUXF7F2bX8VsYiI+Cb6ZsaxsRc9ZOtW6NcPsrIAXsRsG2+99Q59\n+/YMSUQREZH8FIuZ8ZdfQvfu3u8SA9Su/Qx/+tP7KmIREQkLUV/G3sU84NChg0BvKlT4lqSkUtx5\nZ4vQZhQREclH9C1TX1DGjzwCO3bsB+KAAzz//G5uvLFaSOOJiIhcSlTPjBMTYcaMr/Au5vFPxo9P\nZvjw1iGPJyIicimRPzPOzPQupQVgBjExABw8CIMHfwm0A47Tvv1apkzR0rSIiISfyJ8ZX2RWnJUF\ngwfD0aNlgWupWnUdixe3wMyXhCIiIpcU+TPji5Txk09+yerVNYFqwHssWmRUquRLOhERkQIFNDM2\nsy5mttPMUs3sBxeMNLPSZrYo+/ktZlY72EHzlaeMXalSDB36N6ZMaQmMBmD8eKN9+5ClERERuWwF\nlrGZxQCvAHcCDYH+ZtbwgsOGAkedczcA04CpwQ6arzxlvOkkJCS0B2KBX9K8OTz7bMiSiIiIXJFA\nZsYtgFTn3B7nXBrwFtDrgmN6AW9kP14MxJuF6AxtdhlvBXplHAN+DGyka9f6JCfne0EuERGRsBFI\nGdcAvs6zvT9730WPcc5lAMeAysEIWJAli9I5B/QGriIG2MSvf12X5cvReWIREYkIIf0Al5kNA4YB\nXHfddUF5zXq10ykNvA1kWE1OrqxFly5BeWkREZGQCKSMvwGuzbNdM3vfxY7Zb2YlgauB7y98Iefc\nLGAWQPPmzd2VBL5Q4zsqsq3rE+z8OJ0egytRSUUsIiIRJpAy3grUM7M6eKXbDxhwwTGJwCDgQ6AP\nsN45F5SyLVDVqjRdMYlGGQXeylhERCQsFVjGzrkMMxsJJAMxQIJzbruZPQt85JxLBGYD88wsFTiC\nV9ghY6YiFhGRyBXQOWPn3Epg5QX7ns7z+Cxwd3CjiYiIFA+RfzlMERGRCKcyFhER8ZnKWERExGcq\nYxEREZ+pjEVERHymMhYREfGZylhERMRnKmMRERGfqYxFRER8pjIWERHxmcpYRETEZypjERERn6mM\nRUREfKYyFhER8ZnKWERExGfmnPPnjc2+A74K4ktWAQ4H8fWKK41j4WkMC09jWHgaw8IL9hjWcs5V\nvdgTvpVxsJnZR8655n7niHQax8LTGBaexrDwNIaFF8ox1DK1iIiIz1TGIiIiPoumMp7ld4AooXEs\nPI1h4WkMC09jWHghG8OoOWcsIiISqaJpZiwiIhKRIq6MzayLme00s1Qzm3CR50ub2aLs57eYWe3Q\npwxvAYzhWDPbYWafmdk6M6vlR85wVtAY5jnuLjNzZqZPtV5EIONoZn2zvx+3m9mCUGcMdwH8e77O\nzFLMbFv2v+mufuQMV2aWYGaHzOzzfJ43M/td9vh+ZmbNiiSIcy5ivoAYYDdwPRAL/BVoeMExI4CZ\n2Y/7AYv8zh1OXwGOYXugbPbj4RrDyx/D7OPKA5uAzUBzv3OH21eA34v1gG1Axezta/zOHU5fAY7h\nLGB49uOGwF6/c4fTF3AH0Az4PJ/nuwKrAANaAVuKIkekzYxbAKnOuT3OuTTgLaDXBcf0At7IfrwY\niDczC2HGcFfgGDrnUpxzp7M3NwM1Q5wx3AXyfQjwHDAVOBvKcBEkkHH8H+AV59xRAOfcoRBnDHeB\njKEDrsp+fDVwIIT5wp5zbhNw5BKH9ALmOs9moIKZVQ92jkgr4xrA13m292fvu+gxzrkM4BhQOSTp\nIkMgY5jXULz/FUquAscweynrWufcn0IZLMIE8r1YH6hvZn82s81m1iVk6SJDIGM4ERhoZvuBlcDD\noYkWNS73Z+YVKRnsF5ToYWYDgeZAW7+zRBIzKwG8BAz2OUo0KIm3VN0Ob4Vmk5k1ds79y9dUkaU/\nMMc596KZ3Q7MM7NGzrksv4NJrkibGX8DXJtnu2b2voseY2Yl8ZZlvg9JusgQyBhiZh2AJ4Cezrlz\nIcoWKQoaw/JAI2CDme3FO8+UqA9x/UAg34v7gUTnXLpz7kvgH3jlLJ5AxnAo8EcA59yHQBm8ay5L\nYAL6mVlYkVbGW4F6ZlbHzGLxPqCVeMExicCg7Md9gPUu+yy8AAGMoZndCryKV8Q6R/dDlxxD59wx\n51wV51xt51xtvPPuPZ1zH/kTN2wF8u95Gd6sGDOrgrdsvSeUIcNcIGO4D4gHMLOb8Mr4u5CmjGyJ\nwC+yP1XdCjjmnDsY7DeJqGVq51yGmY0EkvE+RZjgnNtuZs8CHznnEoHZeMswqXgn5fv5lzj8BDiG\n/wuUA97O/uzbPudcT99Ch5kAx1AKEOA4JgOdzGwHkAmMc85ppStbgGP4K+A1MxuD92GuwZqg5DKz\nhXj/4auSfV79GaAUgHNuJt559q5AKnAaGFIkOfR3IiIi4q9IW6YWERGJOipjERERn6mMRUREfKYy\nFhER8ZnKWERExGcqYxEREZ+pjEVERHymMhYREfHZ/wPNSlaT2CEdcgAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "results = pd.concat(dfs)\n", "import statsmodels.api as sm\n", "grid = np.linspace(0, 1, 101)\n", "fig = plt.figure(figsize=(8, 8))\n", "plt.plot(grid, sm.distributions.ECDF(results['pivot'])(grid), 'b-', linewidth=3, label='Pivot')\n", "plt.plot(grid, sm.distributions.ECDF(results['pvalue'])(grid), 'r-', linewidth=3, label='P-value')\n", "plt.plot([0, 1], [0, 1], 'k--')\n", "plt.legend(fontsize=15);" ] } ], "metadata": { "jupytext": { "cell_metadata_filter": "all,-slideshow", "formats": "ipynb,Rmd" }, "kernelspec": { "display_name": "Python 3", "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.6.2" } }, "nbformat": 4, "nbformat_minor": 2 }