{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Communities and Crime Cross-Validation\n",
    "Data is taken from [UCI Machine Learning Repository](https://archive.ics.uci.edu/ml/datasets/Communities+and+Crime).\n",
    "> U. S. Department of Commerce, Bureau of the Census, Census Of Population And Housing 1990 United States: Summary Tape File 1a & 3a (Computer Files),\n",
    ">\n",
    "> U.S. Department Of Commerce, Bureau Of The Census Producer, Washington, DC and Inter-university Consortium for Political and Social Research Ann Arbor, Michigan. (1992)\n",
    ">\n",
    "> U.S. Department of Justice, Bureau of Justice Statistics, Law Enforcement Management And Administrative Statistics (Computer File) U.S. Department Of Commerce, Bureau Of The Census Producer, Washington, DC and Inter-university Consortium for Political and Social Research Ann Arbor, Michigan. (1992)\n",
    ">\n",
    "> U.S. Department of Justice, Federal Bureau of Investigation, Crime in the United States (Computer File) (1995)\n",
    ">\n",
    "> Redmond, M. A. and A. Baveja: A Data-Driven Software Tool for Enabling Cooperative Information Sharing Among Police Departments. European Journal of Operational Research 141 (2002) 660-678.\n",
    "\n",
    "*This notebook compares the performance of warped linear regression to ordinary least squares on a cross-validation of the data.*"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Import Dependencies"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import math\n",
    "from scipy.stats import norm\n",
    "from sklearn.linear_model import LinearRegression\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "from sklearn.model_selection import KFold\n",
    "np.random.seed(1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def drop_missing(df):\n",
    "    drop_columns = []\n",
    "    for column in df.columns:\n",
    "        if '?' in list(df[column].values):\n",
    "            drop_columns.append(column)\n",
    "    return df.drop(columns=drop_columns)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def load_dataset(filename):\n",
    "    df = pd.read_csv(filename, header=None)\n",
    "    df = df.drop(range(5), axis=1)\n",
    "    df = drop_missing(df)\n",
    "    X = np.array(df.iloc[:,:-1].values, dtype=float)\n",
    "    y = np.array(df.iloc[:,-1].values, dtype=float)\n",
    "    return X, y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "X, y = load_dataset(\"communities.data\")\n",
    "X = StandardScaler().fit_transform(X)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Set up Cross-Validation Splits"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "cv_splits = list(KFold(10, shuffle=True).split(X))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Set Warping Parameters\n",
    "The warping parameters were fit to maximize the log-likelihood of the training data for each cross-validation using a second order optimizer."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "warping_parameters = np.array([\n",
    "    [2.45539788e+02, 8.21200550e+00, 1.29687713e-01], \n",
    "    [2.20807880e+02, 8.18124658e+00, 1.29775771e-01], \n",
    "    [2.47853763e+02, 8.33731224e+00, 1.26915050e-01], \n",
    "    [2.37965788e+02, 8.22328782e+00, 1.28942838e-01], \n",
    "    [2.26799485e+02, 8.24759809e+00, 1.27894342e-01], \n",
    "    [2.50334254e+02, 8.11527398e+00, 1.32680489e-01], \n",
    "    [2.30428371e+02, 8.34087828e+00, 1.25444850e-01], \n",
    "    [-2.48123866e+02,  8.17254228e+00,  1.31001290e-01], \n",
    "    [2.29218765e+02, 8.16332304e+00, 1.30279786e-01], \n",
    "    [2.21048428e+02, 8.34043339e+00, 1.25391233e-01]])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Decorator for Vector Functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "def vectorizable(f):\n",
    "    def f_(x):\n",
    "        if hasattr(x, '__iter__'):\n",
    "            return np.array([f(xi) for xi in x])\n",
    "        return f(x)\n",
    "    return f_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Build Warping Function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "def make_warper(parameters, y):\n",
    "    y_mean = np.mean(y)\n",
    "    y_std = np.std(y)\n",
    "    def f(t):\n",
    "        t = (t - y_mean) / (y_std * np.sqrt(len(y)))\n",
    "        result = t\n",
    "        for i in range(0, len(parameters), 3):\n",
    "            a, b, c = parameters[i:(i+3)]\n",
    "            result += a**2*math.tanh(b**2*(t + c))\n",
    "        return result\n",
    "    mean = np.mean([f(yi) for yi in y])\n",
    "    @vectorizable\n",
    "    def f_(t):\n",
    "        return f(t) - mean\n",
    "    @vectorizable\n",
    "    def fp(t):\n",
    "        t = (t - y_mean) / (y_std * np.sqrt(len(y)))\n",
    "        result = 1.0\n",
    "        for i in range(0, len(parameters), 3):\n",
    "            a, b, c = parameters[i:(i+3)]\n",
    "            u = math.tanh(b**2*(t + c))\n",
    "            result += a**2*b**2*(1 - u**2)\n",
    "        result /= y_std * np.sqrt(len(y))\n",
    "        return result\n",
    "    return f_, fp"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compute Log-Likelihood Proxy of Training Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "def compute_log_likelihood_proxy(X, y, phi):\n",
    "    f, fp = make_warper(phi, y)\n",
    "    z = f(y)\n",
    "    model = LinearRegression(fit_intercept=False)\n",
    "    model.fit(X, z)\n",
    "    z_pred = model.predict(X)\n",
    "    rss = sum((z-z_pred)**2)\n",
    "    return -len(y)/2*np.log(rss) + sum(np.log(fp(y)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Verify Maximums\n",
    "Tweak the warping parameters to verify they're a local maximum of the log-likelihood."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "def verify_log_likelihood_opt(X, y, phi):\n",
    "    delta_x = 1.0e-3\n",
    "    for i, phi_i in enumerate(phi):\n",
    "        def f(x):\n",
    "            phi_copy = np.array(phi)\n",
    "            phi_copy[i] = x\n",
    "            return compute_log_likelihood_proxy(X, y, phi_copy)\n",
    "        f0 = f(phi_i)\n",
    "        for x in [phi_i - delta_x, phi_i + delta_x]:\n",
    "            delta_f = f(x) - f0\n",
    "            relative_delta_f = delta_f / delta_x\n",
    "            if relative_delta_f > 0 and np.abs(relative_delta_f) > 1.0e-3:\n",
    "                print(i, x, \"\\t\", delta_f, relative_delta_f)\n",
    "                assert False, \"Can't verify optimum\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "for index, (train_indexes, _) in enumerate(cv_splits):\n",
    "    X_train = X[train_indexes, :]\n",
    "    y_train = y[train_indexes]\n",
    "    verify_log_likelihood_opt(X_train, y_train, warping_parameters[index])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Make Predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "def compute_prediction_dist(model, xi):\n",
    "    noise_stddev, beta, A_inv = model\n",
    "    xi = np.array(list(xi) + [1])\n",
    "    std_err = noise_stddev * np.sqrt(1 + np.dot(xi, np.dot(A_inv, xi)))\n",
    "    pred = np.dot(xi, beta)\n",
    "    return pred, std_err"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "def compute_prediction_pdf(model, xi):\n",
    "    pred, std_err = compute_prediction_dist(model, xi)\n",
    "    @vectorizable\n",
    "    def pdf(y):\n",
    "        return norm.pdf(y, loc=pred, scale=std_err)\n",
    "    return pdf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "def compute_warped_prediction_pdf(model, f, fp, xi):\n",
    "    pred, std_err = compute_prediction_dist(model, xi)\n",
    "    @vectorizable\n",
    "    def pdf(y):\n",
    "        return norm.pdf(f(y), loc=pred, scale=std_err)*fp(y)\n",
    "    return pdf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "def fit_linear_regression_model(X, y):\n",
    "    X = np.hstack((X, np.ones((X.shape[0], 1))))\n",
    "    A = np.dot(X.T, X)\n",
    "    A_inv = np.linalg.inv(A)\n",
    "    beta = np.dot(A_inv, np.dot(X.T, y))\n",
    "    y_pred = np.dot(X, beta)\n",
    "    rss = sum((y-y_pred)**2)\n",
    "    noise_stddev = np.sqrt(rss / (X.shape[0] - X.shape[1]))\n",
    "    return noise_stddev, beta, A_inv"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "def compute_linear_model_log_likelihood(X_train, y_train, X_test, y_test):\n",
    "    model = fit_linear_regression_model(X_train, y_train)\n",
    "    result = 0.0\n",
    "    for xi, yi in zip(X_test, y_test):\n",
    "        pred, std_err = compute_prediction_dist(model, xi)\n",
    "        result += norm.logpdf(yi, loc=pred, scale=std_err)\n",
    "    return result"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "def compute_warped_linear_model_log_likelihood(phi, X_train, y_train,\n",
    "                                               X_test, y_test):\n",
    "    f, fp = make_warper(phi, y_train)\n",
    "    z_train = f(y_train)\n",
    "    z_test = f(y_test)\n",
    "    result = compute_linear_model_log_likelihood(X_train, z_train,\n",
    "                                                 X_test, z_test)\n",
    "    result += np.sum(np.log(fp(y_test)))\n",
    "    return result"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compute Cross-Validation Scores"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "linear_cv_log_likelihoods = []\n",
    "warped_cv_log_likelihoods = []\n",
    "for index, (train_indexes, test_indexes) in enumerate(cv_splits):\n",
    "    X_train = X[train_indexes, :]\n",
    "    y_train = y[train_indexes]\n",
    "    X_test = X[test_indexes, :]\n",
    "    y_test = y[test_indexes]\n",
    "    log_likelihood1 = \\\n",
    "        compute_linear_model_log_likelihood(X_train, y_train,\n",
    "                                            X_test, y_test)\n",
    "    log_likelihood2 = \\\n",
    "        compute_warped_linear_model_log_likelihood(\n",
    "            warping_parameters[index], \n",
    "            X_train, y_train,\n",
    "            X_test, y_test)\n",
    "    linear_cv_log_likelihoods.append(log_likelihood1)\n",
    "    warped_cv_log_likelihoods.append(log_likelihood2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.608949\t0.959150\n",
      "0.450565\t1.100966\n",
      "0.624678\t1.028367\n",
      "0.443151\t0.817843\n",
      "0.706899\t1.126734\n",
      "0.589558\t1.017814\n",
      "0.734414\t1.043960\n",
      "0.602487\t0.992717\n",
      "0.467184\t0.992180\n",
      "0.464309\t0.965338\n"
     ]
    }
   ],
   "source": [
    "for index, (log_likelihood1, log_likelihood2) in \\\n",
    "    enumerate(zip(linear_cv_log_likelihoods, warped_cv_log_likelihoods)):\n",
    "        n = len(cv_splits[index][1])\n",
    "        print(\"%f\\t%f\" % (log_likelihood1 / n, log_likelihood2 / n))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.569145\t1.004451\n"
     ]
    }
   ],
   "source": [
    "linear_cv_log_likelihood_avg = np.sum(linear_cv_log_likelihoods) / len(y)\n",
    "warped_cv_log_likelihood_avg = np.sum(warped_cv_log_likelihoods) / len(y)\n",
    "print(\"%f\\t%f\" % (linear_cv_log_likelihood_avg, warped_cv_log_likelihood_avg))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Plot Error Distribution"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "fold = 0\n",
    "y_range = np.arange(np.min(y), np.max(y), 0.0001)\n",
    "train_indexes, test_indexes = cv_splits[fold]\n",
    "X_train = X[train_indexes, :]\n",
    "y_train = y[train_indexes]\n",
    "phi = warping_parameters[fold]\n",
    "f, fp = make_warper(phi, y_train)\n",
    "z_train = f(y_train)\n",
    "linear_model = fit_linear_regression_model(X_train, y_train)\n",
    "warped_model = fit_linear_regression_model(X_train, z_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_error_distributions(index, plot):\n",
    "    xi = X_test[index,:]\n",
    "    yi = y_test[index]\n",
    "    p1 = compute_prediction_pdf(linear_model, xi)\n",
    "    p2 = compute_warped_prediction_pdf(warped_model, f, fp, xi)\n",
    "    plot.plot(y_range, p1(y_range), label = 'OLS')\n",
    "    plot.plot(y_range, p2(y_range), label = 'WLR')\n",
    "    plot.axvline(yi, c='tab:green', label='y_test')    \n",
    "    plot.legend()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaQAAAEYCAYAAAATRII7AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOzdd3xUVfr48c+ZyaQ30iCQhAQIHRJ6EwRRsaAi2FARFER3XcWyqLu4P/W7u66LveyiiIqsHayAIooiaOgQSoIEhAAJAZIQSEjPzPn9cSchgfRpdybn/XqNk7lz586TMI/PnHPPPUdIKVEURVEUVzO4OgBFURRFAVWQFEVRFJ1QBUlRFEXRBVWQFEVRFF1QBUlRFEXRBS9nvllERISMj49v8PnMwkwA4oMb3kdRnG3btm15UspIV8cBTedQNZVLip40N4ecWpDi4+PZunVrg8/fuepOAN694l1nhaQoTRJCHHZ1DNWayqFqKpcUPWluDqkuO0VRFEUXVEFSFEVRdEG/BakwB07+5uooFMX9VJVDeZGro1CUFnPqOaRmM1fC25dB0XH402YI6+LqiDxeZWUlWVlZlJWVuToUl/H19SUmJgaTyeTqUGzz3V8hZydE9XR1JG2KyiHbc0ifBenwr3DmqPbzns9hzJ9dG08bkJWVRVBQEPHx8QghXB2O00kpyc/PJysri4SEBFeH03pSQvrXEAQU57s6mjZF5ZDtOaTPLrtjO7T7kDjI/MW1sbQRZWVlhIeHt8lEAhBCEB4e7v7fbotzofik9rPqtnMqlUO255A+C1LOTgjtDF3GwPFd2rc+xeHaaiJV84jf/9RB7d4nCKrKoLLUtfG0MR7xGbKBrb+/PgtS/u8Q2QM6JEFJPhTluDoiRXEPpw5p9/5h2n3+AdfFoigtpM+CVJgNITEQkag9zv/dtfEoTpOVlcV1111HYmIiXbt2Zc6cOVRUVLB27VomTpx4wf4rVqxgwIABJCUl0bt3b958800XRN04IcQ7QoiTQog9tbY9JYTIFkKkWm9X2eXNCjIBAb6htR4rbYk755D+CpK0aK2i4E4QZj0xVnDItTEpTiGlZPLkyUyaNIn9+/eTkZHB2bNnmTdvXr37V1ZWMnv2bJYvX87OnTvZsWMHY8eOdW7QzbMYuKKe7S9JKZOtt2/s8k7FJ7XWkZev9vhMtl0Oq7gHd88h/Y2yqyrX7kNiIDgGDKZz/eKKR/vxxx/x9fXlzju1aW+MRiMvvfQSCQkJjBs37oL9i4qKqKqqIjw8HAAfHx969Ojh1JibQ0q5TggR75Q3K84F/wgwmkAYoDDLKW+r6IO755C+C5LRC0LjzvWLK07x9PI00o8V2vWYvTsG8+Q1fRrdJy0tjUGDBtXZFhwcTFxcHAcOXHguJCwsjGuvvZbOnTszfvx4Jk6cyNSpUzEY9Nfwb8CfhBB3AFuBR6SUBfXtJISYDcwGiIuLa/yIxfkQYJ3D0ssbCo/ZMVyluVQOtY7+Mtdcod0HRWv3YQmqH7yNkFLWO0qnoe0AixYtYs2aNQwdOpTnn3+eu+66y9Fh2ssCoCuQDOQALzS0o5RyoZRysJRycGRkExMml+RBgPZtF6OP6rJrY9w9h/TXQrJUavf+1qQK7gg5u1wXTxvU1LcwR+nTpw+fffZZnW2FhYUcPXqUrl27Nvi6fv360a9fP6ZNm0ZCQgKLFy92cKS2k1KeqP5ZCPEWsMIuBy7OhYAxwCnw8oEC1UJyBZVDraPDFlIlCCP4hmiPgztpJ2qrKlwbl+Jw48ePp6SkhCVLlgBgNpt55JFHmDFjBv7+/hfsf/bsWdauXVvzODU1lc6dOzsrXJsIIaJrPbwe2NPQvs1mroLSAu0cEmjnX4tz1XV8bYi755D+CpKlSmsdVTcvgztq92ePuy4mxSmEEHzxxRcsXbqUxMREunfvjq+vL8888wwAa9asISYmpua2Y8cO5s+fT48ePUhOTubJJ5/UZetICPERsAHoIYTIEkLMBOYLIXYLIXYB44CHbH6jEutUQQHWgmQ0QVUpVJy1+dCKe3D3HNJfl5258lx3HZwrSIXHtAEOikeLjY1l+fLlF2wfO3YspaUXzjowevRoZ4RlEynl1Ho2v233Nyq1jonwD4NTgNFbe1ycq83coLQJ7pxDOmwhnVeQgqoLkjo5qyiNKreO6vKxdncbrTMun811TTyK0kI6LEhV56Y9gbotJEVRGlZmLUi+wdp9dUGqnmxVUXROfwXp/C473xAwBWgL9imK0rCaFpK1e652l52iuAH9FaTqQQ3VhIDgaHXFuaI0paYgWVtIBtVlp7gXfRUkS5U2RLV2QQKt265IjbJTlEad32UnrJOsqhaS4ib0V5AA/NrV3R4UrZagUJSmlBcBQuvirhYQqc4hKW5DnwWp+htetaAOWgtJXeDnsR566CFefvnlmscTJkxg1qxZNY8feeQRXnzxRfr27XvBa2fMmEFCQgLJyckkJSWxZs0ap8SsO+WFWndd7XnIAqNUl10b4Qk5pLOCZNbuz79mIiham+OutN65JxUPMHLkSFJSUgCwWCzk5eWRlpZW83xKSgqjRo1q8PXPPfccqampvPzyy9x7770Oj1eXygov/DIXEKHNb6d4PE/IITcpSB20ezX022ONGjWqJpnS0tLo27cvQUFBFBQUUF5ezt69e2nXrl0TR4ERI0aQnd1Gr1krL7wwd/wjoFgVpLbAE3JIXzM1yOqCdH6XnfVapKLj0OHC5qZiZ98+Dsd32/eYHfrBlc82+HTHjh3x8vLiyJEjpKSk1CTFhg0bCAkJoX///nh7ezf5NqtWrWLSpEn2jNx9VHfZ1RYQAaWntHnujPpKd4+mcqhVmvyECiHeASYCJ6WUfa3bwoBPgHggE7ipobVcWqT6HFJDLSQ1sMGjVX/DS0lJ4eGHHyY7O5uUlBRCQkIYOXJko6+dO3cujz76KCdPnmTjxo1Oilhnygq1c0a1Va+NVHrqwucUj+PuOdScr0yLgdeBJbW2PQ6skVI+K4R43Pr4MZujaarLTg39do5GvoU5UnUf+O7du+nbty+xsbG88MILBAcHN7lGy3PPPcfkyZN59dVXmT59Otu2bXNS1DpSXgjh5y0xUH0JRXGuKkjOpHKoVZo8hySlXIc2VWNt1wHvWX9+D7BP+85SpV074eVbd7uXD/iFqRaShxs1ahQrVqwgLCwMo9FIWFgYp0+fZsOGDYwYMaLJ1xsMBubMmYPFYuG7775zQsQ6U15UT5edtYWkziO1Ce6eQ60d1NBeSpkDYL1v8KuXEGK2EGKrEGJrbm4Tw08tZjAYzy09UVtQtGohebh+/fqRl5fH8OHD62wLCQkhIkJbUmHfvn11ps9funRpnWMIIXjiiSeYP3++U2PXhYZG2YEaaddGuHsOOfwsp5RyIbAQYPDgwY1fSCTNIBoIKagDFKlRdp7MaDRSWFhYZ1vttVni4+OprKy84HU33nhjncdTpkxhypQpDolRt6rKwVyuWkhtnLvnUGtbSCeqV7y03tvnUvDqFlJ9glULSVEaJmDSAug+oe5mv3bac6ogKW6gtQXpa2C69efpwFd2icZS1XBBCoqGsyfODXxQFOUcL29IvlUbGlybwagNbFDz2SluoMmC1MDyy88Clwkh9gOXWR/brqkuO2lRiaUoLaVma1DcRJPnkBpYfhlgvJ1j0VpIXn71PxcUrd0X5ZwbBq4oStPUbA2Km9Df1EENdtmpa5EUpVUCVEFS3IMbFSRrC0nNZ6coLaO67BQ3oZ+CVFWhnSMyNNCLGBAFwqBaSG1cZmYmH374Yatf/8wzz9gxGjcREKnNlG++cLiv0vboOYf0U5Aqzmr3ooEWktFLK0pqtoY2Tc/J1BAhxDtCiJNCiD21toUJIb4XQuy33jc9DXNrVU8fVHL+hCtKW6TnHNJPQSq3XszVUAsJzi3Up3icv/3tb7zyyis1j+fNm8err756wX6PP/4469evJzk5mZdeegmz2czcuXMZMmQI/fv358033wQgJyeHMWPGkJycTN++fVm/fj2PP/44paWlJCcnc9tttzntd0ObD/KK838VtPkgE4E11seOUXNxrBqh6sk8IYf0Mx99eZF239A5JNDOI53Jck48bdi/N/+b3079Ztdj9gzryWNDG55/d+bMmUyePLlmHq2PP/6YzZs3X7Dfs88+y/PPP8+KFSsAWLhwISEhIWzZsoXy8nJGjRrF5Zdfzueff86ECROYN28eZrOZkpISRo8ezeuvv05qaqpdf7emSCnXCSHiz9t8HTDW+vN7wFrsMUFxfdT0QU6ncqh13KwgdYCsLc6JR3Gq+Ph4wsPD2bFjBydOnGDAgAGEh4c3+brVq1eza9culi1bBsCZM2fYv38/Q4YM4a677qKyspJJkyaRnJzs6F+hperMBymEaHQ+SGA2QFxcXMvfSU0f1CZ4Qg7ppyCVWbvsGjqHBFoLqSRPGwDh1fRCU0rrNPYtzJFmzZrF4sWLOX78eJNT5VeTUvLaa68xYcKEC55bt24dK1euZNq0acydO5c77rjD3iE7RYvmg6yPv7WFpAqS06gcah0dnUOqbiE1UiODrUO/z6rzSJ7o+uuvZ9WqVWzZsqXe5AAICgqiqKio5vGECRNYsGBBzYSRGRkZFBcXc/jwYaKiorj77ruZOXMm27dvB8BkMtU7uaQLOGY+yPr4tdNGqKouO4/n7jmknxZSzaCGJlpIoA1sCG1F14Wia97e3owbN47Q0FCMxvo/B/3798fLy4ukpCRmzJjBnDlzyMzMZODAgUgpiYyM5Msvv2Tt2rU899xzmEwmAgMDWbJEW19y9uzZ9O/fn4EDB/LBBx8489c7X/V8kM9iz/kg62MwqPns2gi3zyEppdNugwYNkg1a/5KcsaCbnPHNHQ3vk7NLyieDpUz7suF9lFZJT093dQjSbDbLpKQkmZGR4bIY6vs7AFulDZ974CMgB6gEsoCZQDja6Lr91vuw5hyr0RyqZca3M+SMb2ec2/D6MCk/urWlfw6lBVQOaWzJIR21kIq0hflEI72ItVtIikdJT09n4sSJXH/99SQmJro6HLuSzpwPsiFq+iCP5wk5pLOC1Eh3HWjLmBtM6uJYD9S7d28OHjxY83j37t1Mmzatzj4+Pj5s2rTJ2aF5hoAIOL7b1VEoDuQJOaSvgtTYgAbQ+sLVxbFtQr9+/Zx+vZBHC4hU55DaGHfMIR2NsitsfEBDtaBoNcGqg2hdvW2XR//+/hFQdkbNZ+dgHv0ZagZbf383LEiqheQIvr6+5Ofnt9mEklKSn5+Pr6+vq0NxjJrZGvJdG4cHUzlkew7pq8vOp5ktpIM/Oz6eNiYmJoasrCxyc9tut46vry8xMTGuDsMxqgtSca5a4NJBVA7ZnkP6Kkh+Pk3vF9QBys9ARTF4Bzg+rjbCZDKRkJDg6jAUR1GzNTicyiHb6ajLrqj555BAddspSkuo+ewUN+CGBUktZa4oLaZm/FbcgD4KkrkKKktANKMHMbijdq+uRVKU5vMN1a7zU0O/FR3TR0EyGOHBPecmT21MTQtJFSRFabaa+exUC0nRL30UJCEgNLbpC2MBfILB5K+67BSlpQIi1bBvRdf0UZBaQgjrtUiqhaQoLRKgZvxW9M39ChJoI+1UC0lRWiYgUnXZKbrmpgVJtZAUpcX81Yzfir65aUGKhsIcaKNTdChKqwREaBeVV1W4OhJFqZf7FqSqUm2ySEVRmkddi6TonE1TBwkhMoEiwAxUSSkH2yOoJlUPDy88Bn6hTnlLRXF7AVHa/dkT567nUxQdsUcLaZyUMtlpxQigXbx2X3DIaW+pKG6vugidyXZtHIrSAPfssgvrot2fOtj4foqinBMSq90XqoKk6JOtBUkCq4UQ24QQs+vbQQgxWwixVQix1W7Tsvu105YzVwVJUZovIAKMPnAmy9WRKEq9bC1Io6SUA4ErgfuEEGPO30FKuVBKOVhKOTgyMtLGt6slrIsqSIrSEkJo3XaqhaTolE0FSUp5zHp/EvgCGGqPoJpFFSRFabmQGHUOSdGtVhckIUSAECKo+mfgcmCPvQJrUlgXreuhqtxpb6ko9iaEyBRC7BZCpAohtjr8DYM7qRaSolu2DPtuD3whhKg+zodSylV2iao5wrqAtEDBYYjs7rS3VRQHGCeldM7FQSGdtMslLObmrT+mKE7U6oIkpTwIJNkxlpapGWn3uypIitJcwZ1AmrW5IEM6uToaRanDPYd9A0R00+7zMlwbh6LYxrkjVdXQb0XH3Lcg+bWDoI5wIt3VkSiKLZw7UjU0TrsvOGzbcRTFAdy3IAFE9YKTaa6OQlFazekjVdvFA0Lr6lYUnXHvgtS+N+RmgLnK1ZEoSou5ZKSqyVcb+p2vCpKiPzZNrupyUX3AXK5dj6QGNijuxzUjVdU1fIpOuXlB6qXdn0xXBUlxOy4bqRrWBdK/dPrbKkpT3LvLLrIHCCOccN71uIri9sK7QmkBlJxydSSKUod7FySTH0T1hizHX+CuKB4jrKt2f0ot36Loi3sXJICYwZC9HSwWV0eiKO4hIlG7z9vn2jgU5Ty6KUir9uRQUmFGyha+MGYIlJ+B/P0OiUtRPE5YF/Dyg+Oqq1vRF10MaiitMPOHD7bjG3sagxBc+/ovjE6MYExiJIPjwzAaRMMvjrEuVJu1VTunpChK4wxGbUDQid2ujkRR6tBFC8nHy8D3D40hMSqIDiG++HgZeOPng9y8cCMj/rWGf65M57fjhfW/ODwRfELg6CbnBq0o7qx9H62F1OIuCUVxHF20kAwGQbeoIMIDvQnHm3dvHUlhWSXrM/L4MjWbxSmZvLX+ECO7hjPzogTG9YjCUN1qMhgg/iI4+JOWXKKR1pSieLAXV++jc3gAPToEYZESQ2O50KEf7PgfFOVoi/Ypig7ooiDVJ9jXxNX9o7m6fzSniitYtu0o7/6aycz3ttIlMoA54xO5pn9HrTB1HQf7VmoX+4V3dXXoiuJ0Z0oqeWPdQSqqtME9AZ0LCPL1YsHa37m0VxSJ7YPqvqBDP+0+Z5cqSIpu6KLLrilhAd7MHtOVdY+O45VbkjEZDMz5OJUrX1nPqj3HkV3GaTv+/qNrA1UUFwnxN5H+9AR+eHgMr00dQPtgXyrNkn+v+o3LXlrHFS+v479rD3CisEx7QXQyGLzg6EbXBq4otbhFQapmMhq4LrkT384ZzatTB1BpsXDv+9uY9FEO5YExcOAHV4eoKC7jZTTQLSqIa5I60jncn/4xIWz+63ieuqY3ft5G5q/ax6hnf+RPH25ny7EyZHQyHFEFSdEP3XbZNcZgEFyb1JGr+nbgix3ZvLA6gyUl/bnz7GqOHTtGXEfVBaEoAFHBvswYlcCMUQlk5hXz/sbDfLr1KCt25fBiaByTKpZDRSkGbz9Xh6oo7tVCOp+X0cCNg2P56c9jCR58E15U8Z//vsI/V6ZzprTS1eEpiq7ERwTwxMTebPzreP55fV+2yR4YLJXMfeUdvtyRTZVZXVyuuJZbF6Rqft5Gbr5uElUhccwK2cKiXw4x9rmfWPzroZqTvIqiaPy9vbhtWGeennMvFuHF8MotPPhJKuNf/JlPthxROaO4jEcUJACEwGvQHSQWb2P17e3pFR3MU8vTufyln1m5KweprrdQlDq8/EMxdBnDDUG7WThtECF+Jh77bDfjnl/L+xsPU15ldnWIShvjOQUJYNCdYPQh8dAHfDBrGO/OGIK3l4H7PtzO9f9NYfMhNbuxotTR4ypE/gEujyzgq/tG8e6dQ4gM8uGJL/dw8fy1vJeSSVmlKkyKc3hWQQqIgKRbYMcHiDNHGdczim/njGH+Df3JOVPKTW9uYNZ7W9l/osjVkSqKPvS+Thv+veN9hBCM6xHFF38cyf9mDiU2zI8nv05j9PyfWLT+IKUVqjApjuVZBQng4se02RrW/B8ARoPgpsGxrP3zOOZO6MHGg/lc/vI6/vThdvbmNDAdkaK0FYFR0PNqSP0QqsoBEEIwOjGST+8ZwUd3D6dbZCD/WLmX0fN/5I2ff6e4vMrFQSueyvMKUkgnGPEn2L0UMr6r2eznbeS+cd34ee5Y/nBxV9buy+XKV9Zz95KtpB497cKAFcXFhsyC0lOw7b06m4UQjOgazkezh7P03hH0ig7m2W9/Y9S/f+S1Nfs5VVzhooAVT+V5BQng4kehfV/46j44faTOU+GBPjx6RU9+fewSHrw0kc2HTjHpP78y6T+/8vn2LNVfrrQ98aOh80Ww7jkor787e0h8GP+bOYwv/jiSgXHteOH7DIb/aw0PfZLKtsMFatCQYheeWZC8fGDK21BVAf+bDEUnLtglxN/Eg5d255fHxvHkNb0pLKvk4U93MvLZH/nnynT2ZJ9RSaa0DULAZU9DSR6s+kujuw6Ia8c7M4aw+qEx3DIklu/TTzBlQQpXvfoLb607SM6ZUicFrXgizyxIAFE94daPoTAb3rpEW1W2HkG+Ju4clcCahy/m/ZnDGBLfjsUpmUx87RfGv/gzL/+QoYqT4vliBsOoB7UZwLcsanL37u2D+L/r+rLJepGtySj45zd7Gfnsj9z05gaWbMjk6KkSx8eteBS3nDqo2TqPhDu/hY9vhUWXwrB7tKQLan/BrkIILkqM4KLECE6XVPDtnuN8lZrNK2v28/IP+4kM8uHi7pGMToxgYFw7Ytr5IdRSF4onGTcPTu6FlX+G0gK46BFteZdGBPhoF9neNqwzh/KKWb7zGF/vPMb/+yoNSKNLZABju0cxoms4A+NCCQ/0cc7vorglmwqSEOIK4BXACCySUj5rl6jsqWMy/OFX+P5J2LgAtrwNva7RhrsmjAG/0AteEurvzdShcUwdGsfJojLWZeSxdt9Jvk8/wbJtWQBEBvkwMC6Uvh1DSGwfSGL7IDqH+eNl9NxGp2J/usohoxfc+C58/QD8+A/4bSWM/Qt0u1RbZbYJCREBPDA+kfsv6cahvGLW7stlbUYu7286zDu/HgKgc7g/A+Pa0Ss6iMT2QSRGBdIpVH25UzStLkhCCCPwH+AyIAvYIoT4WkqZbq/g7MavHVz7KoyaAxv+A2lfwJ5l2nORPbW1YcK6aLegDuAfXnOLCvLlhkEx3DAohiqzhd+OF7HjSAHbj5xm+5ECvks7d37K22ggpp0f0aG+RIf40THUj+gQX9r5mwj19ybU30Son3bv42VQSdjG6TKHTH4weaFWhH54Cj68CQI7QJeLIXaotkJzeFcIiNTO1dZDCEGXyEC6RAZy10UJlFWa2Z19hu2HC9h+pIBfD+TxxY7smv0DvI3EhvkTHeJLdKgfnUL9aB9cnTcmQqw5E+JnwqS+8Hk0W1pIQ4EDUsqDAEKIj4HrAP0VpGrhXWHii3DlfDiyQZt6P2uztvz5ns9A1jOHl8GkJanJDy+TH31N/vT18mWawQvCjZjDjZRUSUoqJUWVUFIJJSclxdlQUiWxYKAcOG69SQQSgUCbtdxoEBiFAaNBYDAYtJsAhEAIgUC7R1DncU0xqylqtYpbA3Wu/s3ntrbV+ihNAQy/97+ueGt95pAQkHQz9LkeMr7VvsAdWAO7Pqm7n3cQ+LcDU4BWnLx8z90bTdpxhAFfYWCI9Ya/gF4Gys2SwnILZ8qqOFNaRXGFhZLjVZQcMVNWaaYCOGG91WYwCLwMAqMQGI1a3hitj4UQNbljsP4aQoBBSx4Mwvppb+CDLhp6VO/uKm9qs/gEM2L2azYfx5aC1Ak4WutxFjDs/J2EELOB2QBxcXE2vJ0dGb0gYbR2q1ZVrg0RP3sCSvLP3SqKobIMKkugynpfWQaWKpBmjBYzQUYzQaKK9iYzWKw3aUaaqzBbLJgtFixSYrFILNafpZRICRLtZ6QEs8RiBiElEhBIQIIE63/QSpp1tXbrNkHLBlyo4Rl1nTUENb2TY+g7h7y8ta7t3tdpH7jCbMjbDwWHrPlxSruvLNHyp6pMu5WdBnMV2mfXUs9N4iMlkdJCZPU2AJN2k4DZmiMWSc29xZozAFJatJ/NIKu0PKpWO13qfNblBVuaRdQ5evNf05acMYTZ5Ti2FKT6vhdc8K8gpVwILAQYPHiwfv+VvHwgIlG72YlA+wN79sgR9xfuurd2nxwSAkJitBvjHPtWqJxxN/YpR7YN+84CYms9jgGO2RaOorQpKocUpRZbCtIWIFEIkSCE8AZuAb62T1iK0iaoHFKUWlrdMpZSVgkh/gR8hzZk9R0pZZrdIlMUD6dySFHqEs6cgUAIkQscbmSXCCDPSeG0lJ5jAxWfLZqKrbOUMtJZwTRG5ZBDqfhazy455NSC1BQhxFYp5WBXx1EfPccGKj5b6Dm2ltLz76Ln2EDFZwt7xaauMlMURVF0QRUkRVEURRf0VpAWujqARug5NlDx2ULPsbWUnn8XPccGKj5b2CU2XZ1DUhRFUdouvbWQFEVRlDZKFSRFURRFF5xSkIQQVwgh9gkhDgghHq/neR8hxCfW5zcJIeJrPfcX6/Z9QogJLorvYSFEuhBilxBijRCic63nzEKIVOvNIVfZNyO+GUKI3FpxzKr13HQhxH7rbboLYnupVlwZQojTtZ5z6N9OCPGOEOKkEGJPA88LIcSr1th3CSEG1nrOoX+3llI55PD4VA7VH5tzc0jWzDztmBvaFei/A10Ab2An0Pu8ff4IvGH9+RbgE+vPva37+wAJ1uMYXRDfOMDf+vMfquOzPj6rg7/fDOD1el4bBhy03rez/tzOmbGdt//9aLMROOtvNwYYCOxp4PmrgG/R5vMcDmxyxt/NQZ8BlUO2xadyqP73c2oOOaOFVLPmi5SyAqhe86W264D3rD8vA8YLIYR1+8dSynIp5SHggPV4To1PSvmTlLLE+nAj2iSYztKcv19DJgDfSylPSSkLgO+BK1wY21TgIzu+f6OklOuAU43sch2wRGo2AqFCiGgc/3drKZVDDo6vESqHnJhDzihI9a350qmhfaSUVcAZtFUBmvNaZ8RX20y0bwTVfIUQW4UQG4UQk+wcW0vim2JtMi8TQlTPIO3ov1+zj2/tokkAfqy12dF/u6Y0FL8zPnctoXLIOfGpHGo5u+aQM5Ydac6aLw3t06z1YmzU7B5/2Y8AACAASURBVPcQQtwODAYurrU5Tkp5TAjRBfhRCLFbSvm7k+NbDnwkpSwXQtyL9k35kma+1tGxVbsFWCalNNfa5ui/XVNc+blrCZVDjo9P5VDr2PVz54wWUnPWfKnZRwjhBYSgNROdsV5Ms95DCHEpMA+4VkpZXr1dSnnMen8QWAsMcHZ8Usr8WjG9BQxq7msdHVstt3BeV4MT/nZNaSh+va1TpHLIwfGpHGo1++aQI0+IWU9ueaGd0Erg3Em7Puftcx91T8h+av25D3VPyB7E/idkmxPfALQTj4nnbW8H+Fh/jgD208gJSQfGF13r5+uBjfLcicVD1jjbWX8Oc2Zs1v16AJlYL8R21t/Oeux4Gj4hezV1T8hudsbfzUGfAZVDtsWncqjhGJ2WQ85KqKuADOsHcp512/+hfVMC8AWWop1w3Qx0qfXaedbX7QOudFF8PwAngFTr7Wvr9pHAbuuHaDcw00Xx/QtIs8bxE9Cz1mvvsv5dDwB3Ojs26+OngGfPe53D/3Zo3yZzgEq0b2wzgXuBe63PC+A/1th3A4Od9XdzwGdA5ZBt8akcqj82p+aQmjpIURRF0QU1U4OiKIqiC6ogKYqiKLqgCpKiKIqiC6ogKYqiKLqgCpKiKIqiC6ogKYqiKLqgCpKiKIqiC6ogKYqiKLqgCpKiKIqiCzYVJCHEQ0KINCHEHiHER0IIX3sFpiiKorQtrS5IQohOwANocxf1RVv58BZ7BaYoiqK0Lbauh+QF+AkhKgF/mphePCIiQsbHx9v4lnVlFmYCEB9s3+MqSrVt27blSSkjXR0HOCaHXEXlbtvR3BxqdUGSUmYLIZ4HjgClwGop5erz9xNCzAZmA8TFxbF169bWvmW97lx1JwDvXvGuXY+rKNWEEIddHUO1+Ph4u+eQq6jcbTuam0O2dNm1Q1tPPQHoCARYV4OsQ0q5UEo5WEo5ODJSF18yFUVRFB2yZVDDpcAhKWWulLIS+BxtfQ5FURRFaTFbziEdAYYLIfzRuuzGA87pSzj5G6x/HnpPcsrbKYqib2WVZr5LO07KgXzyiyuIDPJhTGIEl/Vuj5dRXd3iLmw5h7RJCLEM2A5UATuAhfYKrEEWMyydDrm/QdqXMOBSMPk5/G09XWVlJVlZWZSVlbk6FJfx9fUlJiYGk8nk6lCUFliddpwnv04j50wZ7fxNtA/2ZfOhfD7afIT4cH+endKf4V3CHR6HyiHbc8imUXZSyieBJ205RotlrteK0RXPwg9PQWE2hHdzagieKCsri6CgIOLj4xFCuDocp5NSkp+fT1ZWFgkJCa4OR2kGKSUv/bCfV9fsp1d0MPNv6M+orhEYDAKzRbJm7wme+WYvU9/ayJMTezNjlGP/XVUO2Z5D7teW3bcKjD4w8A6ty644D9Qy7DYrKysjPDy8TSYSgBCC8PDwNv3t1t28sDqDV9fs58ZBMXx13yhGJ0ZiMGifX6NBcHmfDqx8YDSX9mrPU8vTWbT+oEPjUTlkew65X0HK/AU6jwTvAOh1DViqoOyMq6PyCG01kaq19d/fnfxv42Fe/+kAU4fGMv+G/nh71f+/sgAfLxbcNpCr+0Xzz2/2smpPjkPjauufIVt/f/cqSJVlkLsXOg3UHne9BISAstOujUtRFKfZk32Gvy9PZ1yPSP4xqV+T/xP0Mhp44aYkkmNDefjTnWTmFTspUqWl3KsgnUzTWkTRSdpjb3/wCVItJA+SlZXFddddR2JiIl27dmXOnDlUVFSwdu1aJk6ceMH+K1asYMCAASQlJdG7d2/efPNNF0StOEtphZk/fbid8EBvXrwpGaOhed/IfU1G/nPrQLwMgoc+TaXKbHFwpK7jzjnkXgUpZ6d236H/uW0+wVBxFirUtx53J6Vk8uTJTJo0if3795ORkcHZs2eZN29evftXVlYye/Zsli9fzs6dO9mxYwdjx451btCKU72yZj+Z+SW8cFMS7QK8W/TajqF+/PP6fuw4cpq3fznkoAhdy91zyL0KUt4B8PKDdvHntvmGaIMasre5LCzFPn788Ud8fX25805tShmj0chLL73EO++8Q0lJyQX7FxUVUVVVRXi4NqTXx8eHHj16ODVmxXn2HS9i0fqD3DgohpFdI1p1jGuSOnJpr/a8smY/x8943gAWd88hWydXda6CQ1oxqt1n7B2o3R/bAQljXBKWp3l6eRrpxwrteszeHYN58po+je6TlpbGoEGD6mwLDg4mLi6OAwcOXLB/WFgY1157LZ07d2b8+PFMnDiRqVOnYjC41/cspWlSSp78eg9Bvl785apeNh3ryWt6c+mLP/PPb/by2tQBdoqwLpVDreNemXvqEISdN77daAIvHziW6pqYFLuRUtZ7grqh7QCLFi1izZo1DB06lOeff5677rrL0WEqLvBzRi4bD57iocu6E9bCrrrzxYb5c8/FXVm+8xg7j3rWgCh3zyH3aSFJCQWZ2si683kHQo4qSPbS1LcwR+nTpw+fffZZnW2FhYUcPXqUrl27Nvi6fv360a9fP6ZNm0ZCQgKLFy92cKSKM1kskvmr9hEb5sctQ+LscszZY7rw/sbDPL96H/+bOcwux6xN5VDruE8Lqeg4VJVe2EIC8AmEUwfVaDs3N378eEpKSliyZAkAZrOZRx55hBkzZuDv73/B/mfPnmXt2rU1j1NTU+ncubOzwlWcZMXuHNJzCnnksh4NXm/UUoE+Xvzh4q6s35/HpoP5djmmHrh7DrlPQSqwjoppV09Bqj6PVD0KT3FLQgi++OILli5dSmJiIt27d8fX15dnnnkGgDVr1hATE1Nz27FjB/Pnz6dHjx4kJyfz5JNPqtaRh7FYJK//uJ/u7QO5NqmjXY89bURn2gf78OL3GXY9riu5ew65T5ddoXUx2pCYC5+rLkjHd6uBDW4uNjaW5cuXX7B97NixlJaWXrB99OjRzghLcZGf9p0k48RZXro5qWZaIHvxNRmZPaYrf1+Rzo4jBQyIa2fX47uKO+eQ+7SQqgtScPSFzxlNEBAJJ9KdG5OiKA61YO3vdAr1Y2J/+7aOqt0yJJYQPxNv/uzYee6U5nGfglR0HEz+2oWw9Ynqrc3koCiKR9iaeYqthwuYNToBk4PWNArw8WLa8M58l36cg7lnHfIeSvO5UUHKgaAOda9Bqq19H23hPovZuXEpiuIQb647SDt/EzcPiXXo+0wfGY/JaOAtB88GrjTNzQpSI832qF7aKLyCTKeFpCiKYxw9VcIPe09w27DO+Hs79lR3ZJAPUwZ24vPt2ZwpqXToeymNc7OC1KHh56Os4/5PqvNIiuLuPth0BAHcOsw+1x01ZdrweMqrLCzddtQp76fUz6aCJIQIFUIsE0L8JoTYK4QYYa/A6pBSO4fUaEHqCQg1sEFR3FxZpZlPthzhst7t6Rjq55T37N0xmMGd2/G/jYexWNSCn65iawvpFWCVlLInkATstT2kepSdhqoyCG6ky847QJvnTg1scEsPPfQQL7/8cs3jCRMmMGvWrJrHjzzyCC+++CJ9+/a94LUzZswgISGB5ORkkpKSWLNmjVNiVhzjm905FJRUMm14vFPfd9qIzhzOL2H9gTynvq+9eEIOtbogCSGCgTHA2wBSygoppWMmhiq0rvLYWAsJtIENqoXklkaOHElKSgoAFouFvLw80tLOfblISUlh1KhRDb7+ueeeIzU1lZdffpl7773X4fEqjrNkw2G6RAYwqlu4U9/3yr7RRAT68L8NmU59X3vxhByypYXUBcgF3hVC7BBCLBJCBJy/kxBithBiqxBia25ubuve6exx7T6wiYIU1RtO/Q6VF178pejbqFGjapIpLS2Nvn37EhQUREFBAeXl5ezdu5d27Zq+cHHEiBFkZ2c7OlzFQXZnnSH16GmmDe/s9OXAvb0MTB0ay5rfTpJVcOFSDXrnCTlky/AVL2AgcL+UcpMQ4hXgceBvtXeSUi4EFgIMHjy4dZ2zxdYmdEBk4/u17w3SArn7oGNyq95KAb59XJv1wp469IMrn23w6Y4dO+Ll5cWRI0dISUmpSYoNGzYQEhJC//798fZuepbnVatWMWnSJHtGrjjRx1uO4GsyMHlgPTOyOMEtQ+N4/acDLN2axUOXdW/9gVQOtYotBSkLyJJSbrI+XoZWkOyvpiA1sShXzUi7vaoguaHqb3gpKSk8/PDDZGdnk5KSQkhICCNHjmz0tXPnzuXRRx/l5MmTbNy40UkRK/ZUVmlm+c5jTOjTgRA/k0ti6BTqx0XdIli2LYs54xPtPl2Ro7l7DrW6IEkpjwshjgohekgp9wHjAcecwCnJA2EE39DG9wvrAkYfOLHHIWG0GY18C3Ok6j7w3bt307dvX2JjY3nhhRcIDg5uco2W5557jsmTJ/Pqq68yffp0tm1TKwi7mx/2nqCwrIobBrmmdVTtpsGx3P/RDn79PY/RiU30yjRE5VCr2DrK7n7gAyHELiAZeMb2kOpRkg/+YdDUKoZGL4jsoa5FclOjRo1ixYoVhIWFYTQaCQsL4/Tp02zYsIERI5q+osBgMDBnzhwsFgvfffedEyK2jRAiVgjxk/WSiTQhxBxXx+RKn23LIjrEt9XLk9vL5X3aE+pv4tOtWS6NozXcPYdsKkhSylQp5WApZX8p5SQpZYG9AqujOA/8m/khVSPt3Fa/fv3Iy8tj+PDhdbaFhIQQEaH9++/bt6/O9PlLly6tcwwhBE888QTz5893auytVAU8IqXsBQwH7hNC9HZxTC5xsrCMnzNyuX5AJ4wu7ibz8TIyKbkT36Ud53RJhUtjaSl3zyH3WH6iJB/8mzkEtH0f2PkRlJzSWlWK2zAajRQWFtbZVnttlvj4eCorL5za5cYbb6zzeMqUKUyZMsUhMdqTlDIHyLH+XCSE2At0wlFd3zr2ZWo2FglTXNxdV+3GwTEsTsnkq9RjTB8Z7+pwms3dc8g9pg4qzoOAZhakKOsXzBPqAlnFfQgh4oEBwKbG9/Q8UkqWbctiQFwoXSMDXR0OAH06htCnYzCfblVTCTmTexSkkvyWddmBOo+kuA0hRCDwGfCglLKwnudtv5ZPx/ZkF5Jx4qzLBzOc7+YhsaQdK2RP9hlXh9Jm6L8gWcxQWtD0kO9qge3BL0y1kBS3IIQwoRWjD6SUn9e3j5RyofVc7eDIyFaO+tKxZduO4u1lcNgifK11bVJHvI0Glm1zv8EN7kr/BankFCCbfw5JCOvABlWQFH0T2lQEbwN7pZQvujoeV7BIyVc7j3F57/Yuu/aoIaH+3lzWpz1fpmZTXqXWWXMGNyhI1otim1uQwLp67F6wWBwTk6LYxyhgGnCJECLVervK1UE50+mSSk6XVOpmMMP5bhgUw+mSSn7ce9LVobQJ+h9lV5Kv3Te3yw60KYQqi+H0YQhLcExcimIjKeUvgHtNBWBnuWfLiQryYXQ311571JAxiZG0D/Zh2bYsruwX7epwPJ7+W0jV0wY1d1ADQHvr9OpqYIPHyczM5MMPP2z16595xjHXbistV2m2cLqkkusHdMLLqM//FRkNgskDY1ibkcvJojJXh2MXes4hfX4Kaitp5jx2tUX21O7VBbIeR8/JpLRM3tkKpJSO664rOwNZ22DvctjxPmxaCNsWw65P4eDPcPqINmiqCTcMisFskXy5wzNmkddzDum/y67Y2mXn14KLXH0C1WJ9buZvf/sbERERzJmjzZ4zb9482rdvzwMPPFBnv8cff5y9e/eSnJzM9OnTeeCBB3j88cdZu3Yt5eXl3Hfffdxzzz3k5ORw8803U1hYSFVVFQsWLGDlypWUlpaSnJxMnz59+OCDD1zxqypWuWfLCfTxonv7IPsc8PQR2P89HP4VjmyEwmYUEO9AiBkMnUdBj6u0AVHnLXvRNTKQgXGhLN2axd2juzh9WYzm8oQc0n9BKskDnxDwanra9Dqi1BRCrfXvzf/mt1O/2fWYPcN68tjQxxp8fubMmUyePLlmHq2PP/6YzZs3X7Dfs88+y/PPP8+KFSsAWLhwISEhIWzZsoXy8nJGjRrF5Zdfzueff86ECROYN28eZrOZkpISRo8ezeuvv05qaqpdfzel5dKOnaGkvIqEiAuWUGuZ4nzY9Qns+Qyyt2rbgqKh80htuYaI7hASC77BWvExV0BFiVasCg7B8T1a8frpGfjpnxDeDZJvg0Ez6sz0csOgWP76xW52Zp0hObaJSZ5ROdRablCQ8ps/S0Nt7XtDxiqoLAOTr/3jUuwqPj6e8PBwduzYwYkTJxgwYADh4U3/u69evZpdu3axbNkyAM6cOcP+/fsZMmQId911F5WVlUyaNInkZLUciZ58ti0bgxCEB/i07gDHUmHzQti9DMzl0KE/jH8Sel0L4V0vaOVcIKIbcPG5x2dPal17ez6HNU/Dz/+GpFtgzFwIiWFiUjRPL09j2bajzSpIruAJOaT/glSc17Ih39WieoM0Q94+iE6yf1werLFvYY40a9YsFi9ezPHjx5ucKr+alJLXXnuNCRMmXPDcunXrWLlyJdOmTWPu3Lnccccd9g5ZaYVKs4WvUrMJjTPhZWxh91fWNlj7DBz4AUwBMOB2GHo3RPWyLajAKBgyU7udSIfNb0Lqh5D6EQyZRfCYP3NF3w58nXqMJ67uja/J2OjhVA61jv4HNZSeal1Bqp5CSHXbuY3rr7+eVatWsWXLlnqTAyAoKIiioqKaxxMmTGDBggU1E0ZmZGRQXFzM4cOHiYqK4u6772bmzJls374dAJPJVO/kkorzrN2XS35xBZFBLWgdnUiHD2+GRZdA9na49Cl4OB0mvmh7MTpf+95wzStw/zbofyNsWgCvD+GPYdsoLKvk+/QT9n0/O3L3HNJ/C6nklNYcb6mwruDlZ11GeKrdw1Lsz9vbm3HjxhEaGorRWP830P79++Pl5UVSUhIzZsxgzpw5ZGZmMnDgQKSUREZG8uWXX7J27Vqee+45TCYTgYGBLFmyBIDZs2fTv39/Bg4cqAY1uMhn27KICPQm1K8Z54VLC+Cnf8GWRdpgpUuegGH3go+dBkI0JjQOrvsPDPsDLJ9Dj5RH+NSvPx9ueoxrkvQ1zVE1t88hKaXTboMGDZIt9vf2Uq76a4NPz/h2hpzx7Yz6n3zrUinfubLl79kGpaenuzoEaTabZVJSkszIyHBZDPX9HYCt0ol50titVTmkI6fOlstuf10p/748rfHcNZul3PKOlM/GS/lUqJQrHpayON+5wZ4fz+a3ZPnTHWTB/4uWpzZ/esEuKoc0tuSQvrvsKkuhqrT16xp1TIacnc261kBxrfT0dLp168b48eNJTEx0dTiKg3y98xiV5iauPco7AIuvghUPat1x96yDq19w7fpmBgMMmcXJW78nU7an3cpZ8PUDUFXuupjO4wk5ZHOXnRDCCGwFsqWUE20PqZaSU9p9S65Bqq3jAG0kTv4BbWlzRbd69+7NwYMHax7v3r2badOm1dnHx8eHTZva3HJBHmXZtiz6dAymV3Qw7DzvSXMVbPyPNgTbywcmLYCkqU2PmHOimG59mRv9Ktecepdbt7+nTeJ88/sQ7PpphTwhh+xxDmkOsBcItsOx6iq1FqTWfjOKtg5TPJaqCpKb6devn7peyMP8dryQ3dlneOqaelZpz82AL+6BY9uh50StRRTUwflBNsP1g+N59LMpDLvyUrr+8mdYeLFWlHDCea0WcMccsqnLTggRA1wNLLJPOOcpacUsDbVFdAeTPxzbYb+YPJjW1dt2tfXf39GWbc3CZBRcm9yp7hPbFsObY6AgE254V/ufu06LEcBV/aPxMxl5K68vzPoBTH6w+GqoLGnznyFbf39bzyG9DDwKNLjOg02rXVZ32bVm2DeA0UsboZfjXt8SXMHX15f8/Pw2m1BSSvLz8/H1VRdRO0Kl2cKXqdmM79mesADr6DpLlbZMzPI5EDcM/pACfSfrqouuPoE+XlzVL5oVu3IobdcD7v4JOg7E99hm8rMPqhyyIYda3WUnhJgInJRSbhNCjG1oPynlQmAhwODBg1v2L2Vrlx1oAxu2L9EGNhgav5itLYuJiSErKwtPXCK7uXx9fYmJ0ee6PO5u7b5c8s5WnFum/PAGrefCXAGX/R+MuF8bOOAmbhgUw2fbs1iVlsP1A2Lgji+J+fJ+svZIck8mg2+o7gurI9iaQ7acQxoFXGtdUMwXCBZCvC+lvN2GY9ZVUqDdt7bLDrSBDZvegLwM+19A50FMJhMJCWrtKMUxlm07SkSgDxd3j4BNb8J3f4VOHaFDEoya4+rwWmxYQhixYX4s25alFSSTH6Ypb5LwzVz4+i+QfDtc+6r6EtxCrf5KIqX8i5QyRkoZD9wC/GjXYgRaC8k7sOUTq9bWcaB2n7XVPjEpitIi+WfLWbP3JDf2D8P09R/g20ch8XJt0JFPoKvDaxWDQTBlYAwpv+eTVVBi3WjUBmNc/Bikvq8N0jBXuTZQN6PvNnLJKdtaR6DN3uvXDo5utE9MiqK0yFepx4i05DHn8J+0tYjGPQE3fwAG/U8U05gpA2OQEj7fXmuZCyFg3F+1iV53L4Vld0JVheuCdDN2KUhSyrV2vwYJtFF2tl4MZzBA7HA4ot+x94riybZvWscKv6fwPXsUbv0ULp7rVueLGhIb5s+ILuEs25aFxXLe6fHRD8OEf8Her+HTO7RVB5Qm6ftTUXrKPldnxw2D/P3nlkNXFMUpDm/6imcLH8XHxwR3rYLul7s6JLu6cXAMR06VsCXz1IVPjvgjXP0iZHwLn07T1awOeqXvgmSPLjuAuBHa/RHVbacoTrN7GTHf3skROmC+8/tzM/B7kCv6diDQx4ul27Lq32HITLjmVdi/GpbeCWY103xj9F2Q7NVCik4Go7c6j6QozpL6IfLzu9kue7CkxwJC2se5OiKH8Pf24up+0XyzO4fi8gYGMAyaDlc9D/tWwmcz1UCHRui3IJmroOyMfVpIJl9ttJ1qISmK423/H3z5R06ED2Va+VymjOzp6ogc6obBMZRUmPlmd07DOw29GyY8A+lfwZf3qgmfG6DfglRqvQbJXjP8xg3X5rQrP2uf4ymKcqH0r+Dr+6HrJfxJPkZsVDiDO7dzdVQONbhzO+LD/Rvutqs24r5zo+++vh8sDU5w02bpuCDZOG3Q+RLGgKUSDqfY53iKotR1aD18NgtihpA25r9szS7l1mFxCA+fsUAIwQ2DYth86BSH84sb33n0wzD2L5D6Aax8CNroNEMN0W9Bqll6wk7frjqPBKMPHPzJPsdTFOWc3Az4+FZolwC3fsKH23Px8TIweUDbmIpp8sAYhNBWw23SxY/BRQ9rk8p++5gqSrXo98o0G+exO1texdbMUxzMLeZkUTlCwG3ByYSmf0/5RU8RHuhjx2AVpQ0rK9SKkdEbbv+MYmMwX6Vu4er+0YT4m1wdnVN0DPXjom4RfLo1i/vHJ2IyNvJdXwgY//+0efw2vK7NRHPZ39vk3Hfn029BasXifFVmC9+nn+DDzUfY8Hs+VdaL1UxG7R/aQhf+YtrE+H8upWNsArcMieXapE74eav5phTXEEK8A1RPVNzX1fG0mMUCX9wLpw7C9K8hNJbPNx7mbHkVtw3zzJF1DbljRDx3L9nK6rQTXN2/iQX7hIDL/wFVZZDyGnj5wSXznBOojum3ILWwhfTt7hye+24fB/OK6Rjiy91junBRtwh6RQfTzt+EEIKizPaw+CP+3i+XF07E8dhnu3nuu308dFl3bh4ci1dj32oUxTEWA68DS1wcR+tsWqANZ77iWYi/CItFsvjXQ/SPCWFgnGcPZjjfJT2jiA3z491fDzVdkEArSlc+p10wu26+tkrumD87PlAd029BKskHg0mbXLURFVUWDuYV84e120mMCuSN2wdyWe8OGA0XNn+D4gZAQCSXe+/msgcfZNOhU7y4OoN5X+zh/Y1HePGmJG1pZUVxEinlOiFEvKvjaJUT6fDD09Djahh2LwDr9ufye24xL9+c7PGDGc5nNAimj4jnHyv3sjvrDP1iQpp+kcEA17yidd/9+Hfw8oWRf3J8sDql3yZBySlthF0jH+qVu3LYlXWaorJKnr62D6seHMMVfaPrLUaA9o/f/QrIWI0wVzC8Szif3DOc/942kNyiMq59/RcWrW+7C2wp+mTTIpeOUlUOn88G32Dtf6jWPH3310yigny4ql8zWgge6KYhsQR4G3k35VDzX2QwwnX/hd6TYPU82PyW4wLUOf0WpNKCBrvrzBbJ/FW/cd+H2/H1NtKvUwjTR8Y3XIhq6zkRKoq0IapoQzav6hfN6ocu5pKeUfxj5V4e/CSV0gp14ZqiD1LKhVLKwVLKwZGRka4OR7P+RTixG659HQK1mA6cLOLnjFymDe+Mt5d+/9fiSMG+Jm4YFMOKnTnkFrVg7jqjF0xZpLU2v/mztqhoG6TfT01xLgRcmHylFWZmL9nKf9f+ztShcfSJDsHX1IJBCV3Gat2Avy2vszkswJs3bh/E3Ak9+HrnMW58M6VlHyhFaSvyf4dfXoK+N0CPK2o2v/trJt5eBm5tY4MZzjd9ZDwVZgsfbDrcshcaTXDju9DtUvj6Adj5iWMC1DG3KkhnSiuZ9vYmftx3kr9f14d/Te7X8pGSJl/tH/y3by6YvkMIwX3jurHojsH8frKYG99I4eipEht/EUXxIFJq3+C9fGDCP2s2550t11ZPTe7U5i+p6BIZyKW9ongvJbPh+e0a4uUDN78PCaO1KYbSvnBMkDql44KUV6cg5RaVc8vCjezMOs3rUwcybUR864/d6xooPtngrA3je7Xn/VnDKCip5IY3Uth/oqj176UojRBCfARsAHoIIbKEEDNdHVOj0r+C33+ES56AoA41m9/+5RCVZgv3XNzFhcHpxx/HdaOgpJKPNh9p+YtNfjD1Y4gdps188dtK+weoU/osSJVlUF4IAREAnCquYOpbG8nMq0lpzgAAFzFJREFUK2bR9CHNG1LZmB5Xat12Oz9ucJdBndvx6T0jkBKmvrWJAyfVHHiK/Ukpp0opo6WUJilljJTybVfH1KCqCvjhKYjqDYPP1c0zpZX8b8NhruoXTZdI91yS3N4GxrVjRJdw3lp/kPKqVpyP9g7QFjOMToJPp8P+7+0fpA61uiAJIWKFED8JIfYKIdKEEHPsFlWJdSG9gEgKyyq5451NHD1Vwrt3DuHi7nY4qesdoI1oSf8SKhqee6pHhyA+mj0cgFvf2sihvCbmqVIUT7b9PSg4BJc+rZ2Et1qSksnZ8ir+OLabC4PTn/vGdeNEYXndJc5bwjcYbv8MonrBJ7fDwbV2jU+PbGkhVQGPSCl7AcOB+4QQve0SVbE2tLXcJ5yZi7fwW04Rb9w+iOFd7DTRKkDyVKg422RzuGtkIB/dPQyzRXLrWxs5kq/OKSltUHkR/Pxv6HwRJF5Ws7mkoop3fj3E+J5R9O6oruGrbVS3cJJiQliw9neqzK2c2duvHUz7EsK6woc3e3xLqdUFSUqZI6Xcbv25CNgLdLJLVNalxv+9Po+thwt46eZkxvWMssuha8SNhNC4Zg2vTGwfxPuzhlFaaWbqWxvJPl1q31gURe82/Ef7onjZ03WuDVyckklBSSV/HKdaR+erHiR15FRJ61tJAAHhMH05RHSHj6ZC+tf2C1Jn7HIOyXql+QBgkz2OZy46CcD3Ryz86/p+XJPU0R6HrctggEF3QuZ6OJHW5O69ooN5f+YwCssquX3RJk4Wldk/JkXRo7JC2Phf7Rq+mME1m8+UVPLG2t+5tFcUgzx8zaPWuqx3e5JiQ3nphwzKKm24trG6KHUcAEtneOyQcJsLkhAiEPgMeFBKWVjP8y26ytxikXy7aTcAMy8fyi1DHXhNw6AZ2qSGm95o1u59O4Ww+M4hnCgs4/ZFmygornBcbIqiF9sWa6s3X/Rwnc0Lfv6dovIq/jyhh2vicgNCCB67ogc5Z8p4f2MLr0s6n18oTPtCW0rni3tg67v2CVJHbCpIQggTWjH6QEr5eX37tOQqcykl/1i5l+zsI1QafJgxto8t4TXNPwz63wS7PoXi/Ga9ZFDnMBbdMZjM/BLueGczhWWVjo1RUVypqlzrrksYAzGDajafKCxjccohrkvqSM8O6txRY0Z2jWB0YgT/+emA7f+/8AmE25Zq5/FWPKjNmOFBU53ZMspOAG8De6WUL9ojmJd/2M87vx5iWJQFr6Ao56wPMvyP1qR7vdkvGdktgjdvH8Rvxwu5690tlFS08OI3RXEXuz6Bs8dh1IN1Nr+weh9VZsnDl6nWUXM8OqEnBSWVLFj7u+0HM/nBzR9A3ymw5mlY+cgFF/m7K1taSKOAacAlQohU6+2q1h5s0fqDvLJmPzcOiiEptBxRz7RBDhHVE/pcD5verBlM0Rzjekbxyi0D2H6kgLuXbLWtf1hR9Mhihl9fgQ79oeslNZt3HCng061ZzLwogbhwfxcG6D76xYQweUAnFq0/yMFcO1zT6OUNkxfBqDmw9W1tWHiF+48AtmWU3S9SSiGl7C+lTLbevmnNsU4UlvHC6gyu6teBZ6f0RxTlQLADBjI0ZOxfoKoUfn25RS+7ql80z9+YxK8H8rnvg+1UVLVyaKei6NFvKyD/AFz0UE1vhcUiefLrNKKCfLh/fKKLA3Qvj1/VEx8vI08tT7fPigIGA1z2f3DV85CxCt67Bs7qZDb4VtLFTA3tg31Z9ocRvHzzAG3G7sJjzi1Ikd2h/y2waaE2cWQLTB4Ywz8m9WXNbyd56JNUzBbP6c9V2jAp4ZeXoV0C9L6uZvMnW4+yK+sM867uRaCPfpdT06OoIF8evDSRdRm5rE4/Yb8DD71bm//uRBq8NQ5ydtrv2E6mi4IE/7+9Mw+Pokj/+KfmyoSQQA4IGI4EDDdyyC2wXC5RdPF6lBUVjMohiuIq6ur+5Ke/XZFHRVHQxV1RUUEEXMEFDw4FdLmRUwghXCHIGQIkIZNk6vdHdWDIJjCTuTpSn+fpp3u6q7u/qfQ771RX1ftC66tqqZD1rnwoyoPoEOdTGfACWB2w+GmfOwnv6daY5we15N9bDzN+7hbc2ilpqjt7V0DORrhurMrXA/yad46/LfqFrilx/CEYUzGuAIb1SKZ5YjQvfLmdvMIADohqMQjSF6vvrn8OhC2fB+7aIcQ0Duk8pw+rdUxg5th6TXQ96PtnyPxOBZD0kQd7NWHcgGbM25jNCwu26yR/murNqskQVRfa3Q2oEbDPzN9CcambV26/5orLBhso7FYLk+64hmNni3hx4Y7AXvyqDjDie0jqCPMfhG+fh9LqNeDKhA7JmNEcE4aMk11GQP328NW4C47RB8b2v5qRv2vCzNX7mbh4p3ZKmupJzs+QtRy6P6zStQCfb8jm+13HeDqtBckJUWEWWL1p17A2D/dpyryN2XwXyFd3oJIl3vcldH4IfnpL9SudOhjYewQR8zmkM4YjiA7DK4GyrI3FhSoXidu3QQpCCJ5Ja8F93Rvz9xVZTFmaGSShGk0Q+fENiIiBTukAZB49y4QF2+mSEscwf9K+aM7zaL9UWtaP4dn5WzhyOsBRX6x2GPQq3Dodft0C7/asNuGGzOeQTueodThaSAAJqXDDRBVZd+n/+ny6EIIJN7fmjmsbMHlJBu+tyAq8Ro0mWJzYo15Zd0oHZy0KXaWM+WQjTruVN4e0x2LRr+oCgcNmYcqQ9hS4Snnk040UVzX46qVodxeMXAFxKTDnXlj4uAqSa2LM55BOHYDIOJUiIlx0HKYM8sc3YONMn0+3WASv3H4Ng66pz18X/cLby3br13ea6sGPb4LFDt1GI6XkuS+2knH0DG/c1Z76tSLDre43RWpiNC/f1pZ1+3KZuHhncG4S3xTSv4Uej6oQUNO6Q+aS4NwrAJjPIeXuVR49nAgBN0yCJn1h4VjYNs/nS1gtgjfuas+tHZJ49dsM/vzFtqqHoNdoQsHpHNg8CzrcA9H1mLI0k/mbDjFuQDN6ByIPmea/GNw+iWHdG/PPVXv5dE0Vsst6g80Bv/8/SP9aRXn4+Hb4YjQUnAzO/fzAfA7p5F419yHcWO1qbH/DbiqN8JY5Pl/CbrXw+p3teLhPU2atPcCImRvIL6peo140VxD/maqiM1w3lnkbspm8JIPbOibxaD+dWiKYPH9TK/o2r8Pz/9rKkkAPcvCkUTcYuRJ6/UmFhHqrI6x9z1Qj8czlkEpckHcw/C2kMsoCGTa+DuY/BMtf9nmOkhCC8WkteOmWNny/6yiDp/6o06FrzEfBSRU9uu0dLDjg4Km5m+nRNJ6Jt+kh3sHGbrUwdWhH2ibV4pFZG1m12/sQZr7fzAn9/wdGrYTENrDoSTXoYc+y4N3TB8zlkPIOgnRDXJNwK7lARE2VRrj9UPhhIsy+26eYd2Xc260xMx/oSm6+i8Fvr2LxVt+HlWs0QeM/U6E4n2UJd/P47E10So7jvfs6qcnqmqBTw2Hj/eGdSY6PIv3DdSzbGcSWEkBia5Vf6a6PobgAZt4KMwbBvlXBve9lMNfTdjxDreOahldHeWwRMHgqpE1UHYLTusEvC31uLV13dQJfje1JamI0oz/ZyPi5m3X6Ck34OXsUuXoamXWuJ31RPp2T45gxvDNROjRQSImvGcHsEd1onhjNyJkbmLMuyPOHhICWN8OYtZD2CpzYDR8Mgg9uUqOMwzAQy1wO6VeVmI/EVuHVURFCQLfRaiZ0zUQVXXfmLXDEt9nW9WtFMmdkd8b0bcrcDdmkTV7BiozqHRBRU70pXj4Rd3ERD2ancXO7q/gwvYt2RmGidg0HnzzUlW5N4hk/bwsTFmwPzpBwT+xO6DYKHtsMA1+GY7vgo8HwTg/Y8KGalxkizOeQ4ppARHS4lVROYmvllG6YpGa0v9MDPrsXcjZ5fQmHzcJTA1swb3QPnA4r972/lhEfrWf/ifygydZoKmL3tvWIDR8wq7QPg/v1ZsqQ9jjt1nDLuqKJcdqZMbwzD/RM4YOf9nHbtJ/IPBqC+UP2SBWd4/Gt6o2QsKpRxq+3gkXj1XdckFtN5nJIhzerjjazY7VD15EwdpMasZL1A0zvo4IabvoYirwbtNChUSyLxvZifFpzVmUe5/rXVzBhwXYO54XuF4nmyuT0uWJeWriN3DljyCeS5nf9jXHXN9MDGEyCzWrhLze14p2hHcnOLWDQlFVMWbqbQlcI8q7ZnWro/6iVMOwrlS14wwz1HTetG6x8TbWiguCczOOQTh2AU/tVvvjqQo046P8XGLdV5SUpOA5fjoHXmsOcYWqoeGHuJS/htFt5uM/VLH+yD7d0uIqPV++n96TlPDt/CztyTofoD9FcKRS4SvjHyiz6vfoDxavfo4tlJ7aBL9G5TYtwS9NUwA1t6/PNuN70b1mX17/LoN9r3/PZugMUlYTAMQkBKb3gzg/hyQy4abIKKbX0RZjaBd66VgVw3bdKjZAOAOZ5Ubx3hVon9wqvjqrgrKUyN/YYCwfXqMmFOxfBjn+BxQZXdVSONrknNOyiypcjMcbJpDva8Wi/VP6+Yg9z1mcza+1B2jWoxZ2dG/L7VvWoEx0Rhj9O81tg3/F8Pt9wkFlrD3Iy38XdDY4zwf0pNBlAVNfh4ZanuQR1o51MG3ota7JO8NdFv/D0vK289m0Gw3okc1vHpNBE0IiMVdFrOqVD3iHYtUgtq99VQVyj6sKfdp5PVVJVRChD2nTq1EmuX7++4oMzb4XjmapjzeJ9w+3+r+8HYEbajEBIDBxuNxzaoP5p+1ap3DJuYwJabArUa6uWxDYQfzXENlaj+Qxy8118sekQs9cdIOPIWYSAaxvF0q9lXbqmxNE2qbYekhsihBAbpJSdwq0DLmNDHkgpyThyluW7jrJkxxHW78/FIqBv87o80UHS+ps/qj6DkStUSz8MmNZ2TYyUklWZx5m+IouVu4+rsVYp8aS1qUev1ARSEqJC+9r13GnYt1JF+ejyUKXFvLUhv1pIQog04E3ACvxDSjmxShfKP6H6YXo94ZMzMjUWCzTsrBZQ+e4PrYcDa+DIVvh1G/ziGYFXQK0GEJsMscnERtcnPTqR+9MS2VcUw/JsCwv2FDLp610ARNgstGtQm5b1o2lWL5rmidEkJ0QRH+XQ/QDViEDZkKvEzZq9J9h88BSbs/PYfPAUR88UAdCiXjRPDWzO7R2SqHd4iXqtbI2Ae78ImzPSVA0hBL1S69ArtQ77jufz5c85fPnzIV5YsB2ApNqRtG9YmzZJtWiTFEOTOjWpF+NUmbiDgTNGJQcMEFV2SEIIKzAVuB7IBtYJIRZIKX3POhUVrwYI2JxVlWN+HDVU52BK7wv7is7A0Z1wMkstuXvVOuMbyD8GSASQYizpgIyKoMgewxlqcPJ4JEcOO8l11yBDOtmCA5clAoczCocziojIKGwRNbAb25FOJxGOCGwOOzabHbsjArvdjsPuwO5wYLHasVhtCKsDi9WGxSKwWCwIiwWrsCIsAqvFgsUiEAj1jrlsLSznt4UodwzjeNkxDRBYG3KVurnv/bUg3TRLcNI/JYpOjerSMyWGRPdRyFkGcz+H7LWqZX7nTPNERNFUieSEKB4bkMpjA1I5cKKAlZnH+CnzBFsOneLfHhPv7VZBg9gaJNWOJC7KQVyUg/goB7FRDqKdNpx2K5F2K5EOtXbarditAqvFWITHtsdiEQKBcpIWQUBs258WUhcgU0qZhRI1GxgMVC0NYmxjP6RUUyKiL25FeVJaopzS2V/hzBG1zj+OOJeH81weznOnqHMuj2bn8ijNz6G06CyiuBBb6TksrlJwAb/xMRFu6Z0BXO6l9CkRQ/yEIAW2vDQBs6GawkWW816ELIWzQIaxeJLQDG58VUWztzn81a4xEY3iazA0vjFDu6rv0VMFLnbknGbfiQIOnCzgYG4Bh3ILOZhbwMmzLs4EOKZmYkwEa/48wO/r+OOQkgDPqcTZQNfyhYQQI4ARAI0aNfLjdlcYVpvKCXWZvFAC9U+86B9ZWqwms5Wcg+IC3K5CCs6e4UxhIUVFLopLXBS7iikpcVFcXExpsYvSkmJwF6vgmu4ShLsEt9ttpM1wIyVIdylIaeyTIEGgtgUSYZRVlH2+cMy73spKSpXbLTx2XPq6XtzVHkn3y5cKBoGzIVsEoufjYHWoaQkWu7Ftg5gGUKe5eh2sW6hXBLVrOOhxdQI9KomL6ypxk1vgIr+ohMLiUgpdpRetS0olpVJS6r54cUtJSdm2W9m0lBAVEZi5a/44pIqe7P+yfinldGA6qA5ZP+6n8RarXS3EAGpsf81EqBlWUZoKCJwNWawqaKZG4wUOm4XEGPN1kfgzgiAbaOjxuQGQ458cjeaKQtuQRuOBPw5pHZAqhEgRQjiAIUD1SNyu0ZgDbUMajQdVfmUnpSwRQjwCfIMasvq+lHJ7wJRpNL9xtA1pNBcT0omxQohjwP5LFEkAgpidyi/MrA20Pn+4nLbGUkpT5PDWNhRUtL6qExAbCqlDuhxCiPVmmRFfHjNrA63PH8yszVfM/LeYWRtoff4QKG2/kbAIGo1Go6nuaIek0Wg0GlNgNoc0PdwCLoGZtYHW5w9m1uYrZv5bzKwNtD5/CIg2U/UhaTQajebKxWwtJI1Go9FcoWiHpNFoNBpTEBKHJIRIE0LsEkJkCiGeqeB4hBDiM+P4GiFEssexZ439u4QQA8Ok7wkhxA4hxBYhxFIhRGOPY6VCiJ+NJSiz7L3QN1wIccxDx4Mex4YJIXYby7AwaJvsoStDCHHK41hQ604I8b4Q4qgQYlslx4UQYoqhfYsQoqPHsaDWm69oGwq6Pm1DFWsLrQ1JI3pzsBbUDPQ9QBPAAWwGWpUr8zDwrrE9BPjM2G5llI9ApQTaA1jDoK8vUMPYHl2mz/h81gT1Nxx4u4Jz44AsYx1rbMeGUlu58o+iohGEqu56Ax2BbZUcvxFYjApy2g1YE4p6C9IzoG3IP33ahiq+X0htKBQtpPM5X6SULqAs54sng4EPje25QH8hhDD2z5ZSFkkp9wKZxvVCqk9KuVxKWWB8XI0KghkqvKm/yhgIfCelPCmlzAW+A9LCqO2PwKwA3v+SSClXACcvUWQw8JFUrAZqCyHqE/x68xVtQ0HWdwm0DYXQhkLhkCrK+ZJUWRkpZQmQB8R7eW4o9HnyAOoXQRlOIcR6IcRqIcQtAdbmi77bjSbzXCFEWQTpYNef19c3XtGkAMs8dge77i5HZfpD8dz5grah0OjTNuQ7AbUhf/IheYs3OV8qK+NVvhg/8foeQoh7gE7A7zx2N5JS5gghmgDLhBBbpZR7QqxvITBLSlkkhBiF+qXcz8tzg62tjCHAXCllqce+YNfd5Qjnc+cL2oaCr0/bUNUI6HMXihaSNzlfzpcRQtiAWqhmYijyxXh1DyHEAOA54A9SyqKy/VLKHGOdBXwPdAi1PinlCQ9N7wHXentusLV5MIRyrxpCUHeXozL9ZstTpG0oyPq0DVWZwNpQMDvEjM4tG6pDK4ULnXaty5UZw8UdsnOM7dZc3CGbReA7ZL3R1wHV8Zhabn8sEGFsJwC7uUSHZBD11ffYvhVYLS90LO41dMYa23Gh1GaUaw7sw5iIHaq6M66dTOUdsoO4uEN2bSjqLUjPgLYh//RpG6pcY8hsKFQGdSOQYTyQzxn7XkT9UgJwAp+jOlzXAk08zn3OOG8XcEOY9C0BjgA/G8sCY38PYKvxEG0FHgiTvpeB7YaO5UALj3PTjXrNBO4PtTbj8wRgYrnzgl53qF+Th4Fi1C+2B4BRwCjjuACmGtq3Ap1CVW9BeAa0DfmnT9tQxdpCakM6dJBGo9FoTIGO1KDRaDQaU6Adkkaj0WhMgXZIGo1GozEF2iFpNBqNxhRoh6TRaDQaU6Adkkaj0WhMgXZIGo1GozEF/w/6yu3BaoH/UgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 4 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "nrows = 2\n",
    "ncols = 2\n",
    "indexes = np.random.choice(range(len(cv_splits[0][1])), nrows*ncols)\n",
    "fig, axes = plt.subplots(nrows=nrows, ncols=ncols)\n",
    "index = 0\n",
    "for i in range(nrows):\n",
    "    for j in range(ncols):\n",
    "        plot_error_distributions(indexes[index], axes[i, j])\n",
    "        index += 1\n",
    "fig.tight_layout()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "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.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}