{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 一.简介\n",
    "\n",
    "Dropout技术在深度学习中对防止过拟合起到了很好的作用，google甚至为其申请了专利，论文《Dart:Dropouts Meet Multiple Additive Regression Trees》将dropout应用到了gbdt中，这种技术称作DART。简单来说就是在训练过程中暂时丢弃部分已生成的树，使得模型中树的贡献更加均衡（一般最先生成的树的贡献更大），防止过拟合。\n",
    "\n",
    "### 二.流程\n",
    "分两步：  \n",
    "\n",
    "（1）在进行每一轮训练时，对当前已经生成好的$n$颗树随机丢弃掉$k$颗，对对剩下的$n-k$颗树计算其负梯度，并训练一颗新的回归树去拟合该负梯度；  \n",
    "\n",
    "（2）执行标准化操作，由于丢掉了部分的树，所以新训练的树的预测结果其实是超出了拟合目标的，需要对其做标准化操作，对丢弃的树乘以$\\frac{k}{k_+1}$的权重，对新训练的树乘以$\\frac{1}{k+1}$的权重\n",
    "\n",
    "### 三.代码实现\n",
    "代码实现很简单，就直接在GBDTRegressor和GBDTClassifier上面微调即可"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "os.chdir('../')\n",
    "from ml_models.tree import CARTRegressor\n",
    "import copy\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "\n",
    "\n",
    "\"\"\"\n",
    "DART回归模型，封装到ml_models.ensemble\n",
    "\"\"\"\n",
    "\n",
    "\n",
    "class DARTRegressor(object):\n",
    "    def __init__(self, base_estimator=None, n_estimators=10, loss='ls', huber_threshold=1e-1,\n",
    "                 quantile_threshold=0.5, dropout=0.5):\n",
    "        \"\"\"\n",
    "        :param base_estimator: 基学习器，允许异质；异质的情况下使用列表传入比如[estimator1,estimator2,...,estimator10],这时n_estimators会失效；\n",
    "                                同质的情况，单个estimator会被copy成n_estimators份\n",
    "        :param n_estimators: 基学习器迭代数量\n",
    "        :param loss:表示损失函数ls表示平方误差,lae表示绝对误差,huber表示huber损失,quantile表示分位数损失\n",
    "        :param huber_threshold:huber损失阈值，只有在loss=huber时生效\n",
    "        :param quantile_threshold损失阈值，只有在loss=quantile时生效\n",
    "        :param dropout:每个模型被dropout的概率\n",
    "        \"\"\"\n",
    "        self.base_estimator = base_estimator\n",
    "        self.n_estimators = n_estimators\n",
    "        if self.base_estimator is None:\n",
    "            # 默认使用决策树桩\n",
    "            self.base_estimator = CARTRegressor(max_depth=2)\n",
    "        # 同质分类器\n",
    "        if type(base_estimator) != list:\n",
    "            estimator = self.base_estimator\n",
    "            self.base_estimator = [copy.deepcopy(estimator) for _ in range(0, self.n_estimators)]\n",
    "        # 异质分类器\n",
    "        else:\n",
    "            self.n_estimators = len(self.base_estimator)\n",
    "        self.loss = loss\n",
    "        self.huber_threshold = huber_threshold\n",
    "        self.quantile_threshold = quantile_threshold\n",
    "        self.dropout = dropout\n",
    "        # 记录模型权重\n",
    "        self.weights = []\n",
    "\n",
    "    def _get_gradient(self, y, y_pred):\n",
    "        if self.loss == 'ls':\n",
    "            return y - y_pred\n",
    "        elif self.loss == 'lae':\n",
    "            return (y - y_pred > 0).astype(int) * 2 - 1\n",
    "        elif self.loss == 'huber':\n",
    "            return np.where(np.abs(y - y_pred) > self.huber_threshold,\n",
    "                            self.huber_threshold * ((y - y_pred > 0).astype(int) * 2 - 1), y - y_pred)\n",
    "        elif self.loss == \"quantile\":\n",
    "            return np.where(y - y_pred > 0, self.quantile_threshold, self.quantile_threshold - 1)\n",
    "\n",
    "    def _dropout(self, y_pred):\n",
    "        # 选择需要被dropout掉的indices\n",
    "        dropout_indices = []\n",
    "        no_dropout_indices = []\n",
    "        for index in range(0, len(y_pred)):\n",
    "            if np.random.random() <= self.dropout:\n",
    "                dropout_indices.append(index)\n",
    "            else:\n",
    "                no_dropout_indices.append(index)\n",
    "        if len(dropout_indices) == 0:\n",
    "            np.random.shuffle(no_dropout_indices)\n",
    "            dropout_indices.append(no_dropout_indices.pop())\n",
    "        k = len(dropout_indices)\n",
    "        # 调整对应的weights\n",
    "        for index in dropout_indices:\n",
    "            self.weights[index] *= (1.0 * k / (k + 1))\n",
    "        # 返回新的pred结果以及dropout掉的数量\n",
    "        y_pred_result = np.zeros_like(y_pred[0])\n",
    "        for no_dropout_index in no_dropout_indices:\n",
    "            y_pred_result += y_pred[no_dropout_index] * self.weights[no_dropout_index]\n",
    "        return y_pred_result, k\n",
    "\n",
    "    def fit(self, x, y):\n",
    "        # 拟合第一个模型\n",
    "        self.base_estimator[0].fit(x, y)\n",
    "        self.weights.append(1.0)\n",
    "        y_pred = [self.base_estimator[0].predict(x)]\n",
    "        new_y_pred, k = self._dropout(y_pred)\n",
    "        new_y = self._get_gradient(y, new_y_pred)\n",
    "        for index in range(1, self.n_estimators):\n",
    "            self.base_estimator[index].fit(x, new_y)\n",
    "            self.weights.append(1.0 * (1 / (k + 1)))\n",
    "            y_pred.append(self.base_estimator[index].predict(x))\n",
    "            new_y_pred, k = self._dropout(y_pred)\n",
    "            new_y = self._get_gradient(y, new_y_pred)\n",
    "\n",
    "    def predict(self, x):\n",
    "        return np.sum(\n",
    "            [self.base_estimator[0].predict(x) * self.weights[0]] +\n",
    "            [self.base_estimator[i].predict(x) * self.weights[i] for i in\n",
    "             range(1, self.n_estimators - 1)] +\n",
    "            [self.base_estimator[self.n_estimators - 1].predict(x) * self.weights[-1]]\n",
    "            , axis=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = np.linspace(1, 10, num=100)\n",
    "target = np.sin(data) + np.random.random(size=100)  # 添加噪声\n",
    "data = data.reshape((-1, 1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX8AAAD8CAYAAACfF6SlAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzt3X2UXHWd5/H3tzudpOkm3YSEPDTBJMiDQJSWHhFzDqsRFh+YEBgdxZlZ9LCLjOPooMNO2Nl1PB7PkBk8urq6uBFZ0XFABjFGieJDYJgFB+nQIITnJ0l3EpIQOmmSTtLp/PaPW9Wprr63q6rr3rpPn9c5fbqr6tL32lZ987vf3/f3/ZlzDhERyZemuC9AREQaT8FfRCSHFPxFRHJIwV9EJIcU/EVEckjBX0QkhxT8RURySMFfRCSHFPxFRHJoWtwXEGTOnDlu8eLFcV+GiEiqbNq0aZdzbm6l4xIb/BcvXkxvb2/clyEikipm9vtqjlPaR0QkhxT8RURySMFfRCSHFPxFRHJIwV9EJIcU/EVEciixpZ5Jtq5vgBvufpqtg8Ms7Gzl2otOY1V3V9yXJSJSNQX/Gq3rG+C6Ox9jeGQUgIHBYa678zEA/QMgIqmhtE+Nbrj76bHAXzQ8MsoNdz8d0xWJiNROI/8abR0crul5ERmfKu1obcEMBvePKG0ao9wH/1rflAs7WxnwCfQLO1sbdckiqVKeKh0cHhl7TWnT+OQ67VN8Uw4MDuPw3pSv7R/BcfRNua5vYNx/c+1Fp9Ha0jzuudaWZq696LTGXbhIivilSkspbRqPXI/8q31Tlo5Iij+r2kekOtWkRGtJm6raLhy5Dv5TfVOu6u7Sm02kSkGp0vJjqqFqu/DkOu1TzRtOuXyR+vilSkvVkjZVtV14ch38w3xTioi/Vd1dXH/ZMro6WzGgs7WF445pwYCuzlauv2xZ1aN2VduFJ9dpn/L8vUrQRKIRVqpU1XbhyXXwB+XvRdLk2otOG5fzB92hT1Xug7+IpIeq7cKj4F8nlZ2JNJbu1sOR6eAfdWBW2ZnIeBoMpUdmq33KV+8Grdith8rORI5qxGdOwpPZ4N+IwKyyM5GjNBhKl8ymfRoRmMMsO9PtsqSdBkPpktmRf1AADrMeOKwmb7pdlixoxGdOwpPZ4N+I7pvlKxdrXa1YpNtlyQJ1vE2XzKZ9GlUPHEbZmW6XJQtUg58umQ3+0Ph64Knm7bVkXbJCNfjpkdm0T6PVk7fX7bKINFoowd/MbjazHWb2eMDr7zSzPWb2SOHrc2GcN0nqyduHNXcgIlKtsNI+3wG+Dnx3kmP+zTl3cUjnS5xq8/ZBqSHdLotII4Uy8nfO3QfsDuN3pVU1ZW4q6RSRpGhkzv88M3vUzH5mZmc28LwTrOsbYPmajSxZfRfL12wMJfhWk7dXSaeIJEWjqn0eBt7gnHvdzN4HrANOKT/IzK4CrgI46aSTIrmQqJqxVVPmppJOSTOtQs+WhgR/59zekp83mNn/NrM5zrldZcetBdYC9PT0uCiuZbLRd71v5Ep5e5V0Slqpg232NCT4m9l84BXnnDOzt+Glm15txLnLTWX0HdaIR7sQSRr4vd+jHDRJPEIJ/mZ2K/BOYI6Z9QN/B7QAOOe+CXwA+HMzOwwMAx92zkUysg9SfEMHnTRo9B3miEcrICXpgt7v5YG/SCnL9Aol+DvnLq/w+tfxSkFjUf6GLjfZ6DvsEY9KOiXJgt7vzWaM+ozXlLJMr0y3dyjye0MXdVUYfUc5SasJNEmaoPf1qHO0tjSnMmWpz5m/XLR3CHpDG3D/6hUVJ2lreb5aqvmXJAp6XxdXnadtFbo+Z8GyOfK/5RbYunXs4d888jx7h0cmHDartQWu/92kv2rt1r3cvXk7h0ePjD03rbmJi86cD9c/6D0xYwZcdRW0t1d9iZpAkySarCghjSlLfc6CZTP433gjPPjg2MOrJzv27sl/1ZmFrwnuKXvc0QFXXlnV5YFq/iU69aQ5slaUoM9ZsGwG//vug7LJqfWPDPDlXz7DtsEDLOicyWcuPJWVZ4fwhnYOTjgBHn64puCvmn+JQhjVaWkc4QfR5yxYNoP/9OkTnlp57lJWnrs0mvOdfTb09dX0n6jmX6KQxTRHPXcy+pwFy2bwb7TubrjpJhgdhebmyseTvdtrSYaspTnqvZPR5yyYgn8Yurth/3549lk4/fSq/7Nqbq9Vpia1yFqaI4w7mSylscKk4B+G7m7ve19fTcG/EvVTkVpNluYoHUh0tLZgBoP7RxI9qMjanUyS5KLOP3JnnOHNM9SY969ELaClVkG7wgHj6t0Hh0d4bf9I4mvfo1pnIxr5h6OlBc46K7TgXxyh+d2+g0Y94qllV7jlazYGrnKH5E4Ka8I2Ogr+YenuhnXrvNJPsyn/mkp9iMAb9WguIN9qTQlWM2BI4qBCE7bRUfAPS3c3fPvb0N8PixZN+ddM1ocIvFHPu06fq7mAnKt1IjRoIrj8mIoOH4aDB2u61rpMn64J24go5x+W0knfOkw2+irmb+95aqfmAnKu1olQv21GS1WVSjl8GJYu9dqYNOpr7lx47bWq/y5SPY38w/LmN3vpnr4+WLlyymmZoBFaV2cr969eAcA1P3jE979N4m27RKPWks7y9MmUqn2efBK2bIErroAz69+G+/GBPfzs8e2MlPTNamlu4r1nzeesrg6vdPpb34KnnoLzzqv7fDKegn9Y2tvh1FOhr6+uEs1qJriyVssttat2IjTUuaHeXu/7ddfBaUfPE3SOSuf++JqNDMyc+D5e19nK/deugM2bveD/4osK/hFQ8A9Tdzc88EBdC1OqmeBSBYRU8z4JfZ3Ipk1w7LFwyikVz9H7+938cNPApOeumLpavNj7/tJLtV+rVKTgH6bubrjtNvZvewVaZ014udq0TKUJLlVACFR+n4Te56e3F976Vmg6OlUYdI5bH9wyYeev8nNXvINta/OaJr74Yu3XKhUp+AeY0u1yYdL3ut5/4emZx094uaO1Bb7yTCjXtwpYNQN4x8mwckUov1OyJdTVsSMj8Oij8IlPVPW7/LZ8LD++mtXI32juZOSehxnoG9DgJmQK/j6mfLvc0wMdHfzxAz8KPuauMK8UbxQ2POzbyVTyLdS5oSeegAMHvPc4RwdH/iGeqvb8DbqDBcY+f1s65rFs+3P8J5Uyh07B38eUb5ePOw5eeQUOHOCnj27laxufZfueA8zvmMmnVpzCxW9ZGO6FfvObsHo1vP46zJ4d7u+W1At1bmjTJu/7OedUXIjY2tLMH53TNS7nH3TuSquRt3TO46JnfsPBg4cSuQI5zRT8fdR1uzxjBsyYwcXnd3Dx+W8K+crKzJkDwN0PPs8XHn1E+X8ZJ9S5od5emDUL3vhGbvjHewMDf1fJOXreMHtK5y79nPV3zGP6kcPMe303W5uqa5cu1VHw95GaUsq2NgC+tr6PgQ7vQ6XVvlIqtNWxmzaNTfYGDYIMxtai1HPu0s/flo55ACza8wpNJ51U+3VLIK3w9eG3GjKRpZSF4N+0f/+4p7XaV0J16JA32VvI90fVaXNd3wDL12xkYHCYYnes/kLwXzq0I3mfv5TTyN9Hakop29sBOGbkwISXtNpXajFpddvmzV4/n3POAaJZZ1I+j+Dw7iQGZp3AETOuWAhvStrnL+UU/AOkoplUYeTvF/wTl6KSxJqsug3gqS/+E6uBD/Ue4vLTBiIZHPkVWThg7pxZNC1cyJsOvDrl3y3+Qgn+ZnYzcDGwwzl3ls/rBnwVeB+wH/ioc+7hMM6da4WRf+eRQ+OeTmSKShIrqLrt8+s3c/DwEf77C0+wd0Ybv20+jt+VzCeFOTiatMhi8WIt9IpAWDn/7wDvmeT19wKnFL6uAm4M6bwTFPOGS1bfxfI1GxO5O1FoCiP/Pzvr+Ak7NyX+rkUSIyjwDg6PMDwyylnbn+fxeSfjrCmy+aRJ5xGWLFGLhwiEMvJ3zt1nZosnOeQS4LvOOQf8u5l1mtkC59y2MM5flLs9bwvB/61zpnP/X2qVb27cfLO3xiMkd23by6GSzprlznzleW76g1Vjj6OYT5p0HmH/Evjnf/ZWGbe0hH7uvGpUzr8L2FLyuL/wXKjBP/ReJklXSPuwb1+81yGNdfvt8Mwz8I53hPLrZk9v55ntQ+NW5Dab0dRkjIwe4Z6Te/jRWUcHF1HMJ006j9C3GI4cgZdfhpNPDv3cedWo4O+3r+GEtd9mdhVeWoiTplDTG2ovkzSYMcNr7/D663FfiTTS0JBXc79hQyi/bj7w7z7VPkBDu8cGziMsWeJ9f+klBf8QNSr49wOlexueCGwtP8g5txZYC9DT0xPUNiRQUhdnRbbfrpk3+tfIP1+GhsZWd4dlsgnc2Euei8Ffk76halTwXw980sxuA84F9oSd74dk9rmPfB6irU0j/7wZGvJaLTRAIkqeTzwRmpsV/EMWVqnnrcA7gTlm1g/8HdAC4Jz7JrABr8zzObxSz4+Fcd5ySVycFfk8RFubRv55MzTkbaqSF9OmwaJFqvgJWVjVPpdXeN0BfxHGuSpJxEilROTzEEr75E/egj94qR+N/EOl3j4Ri6oPyhilffLl0CHvK2/BXwu9QqfgH7HIm8Rp5J8vQ0Pe97wF/yVLYPt2b+OiKuRqsecUqbdPxCKfh2hrg/7+cH6XJF9eg//Spd73Y47xqtwKXOGr2Amu+MpKBytL//vrwJl/zXmgP/kT+N73xh5GVrUXEwX/Boh0HkITvvmS1+C/ciX8/d+PG/k/tX0vG5/cwciRo1XhLU3GtGZjeGTiiuVZM6fxseVLqjvfnXfCQw+NPcxi9wAF/7Rrb1fOP0/yGvyPPRauu27cU1eu2cjA7OoLJwz42BfeX93B27fD+vVjD7PYPUA5/7TTyD9f8hr8fdRaMVdTkUVnJwwOVjxXmrsHKPinXVubdys86r+nqmTM3r3e9wYt8kqyWoJ5zUUWnZ3eBjYHDkx6rri7B9RDwT/tis3dyrZylIzSyH+MXyWdnym1Oe/s9L7v2RN4rri7B9RLOf+0K7R1Zt8+BYQ8UPAfU1pJ59fTCyZuKl9q0uqdYvAfHIR58xLZPaBeCv5pVxz5a9I39aoqJVTwH6dYSVfc+L1cUFqmYvVOafAvO1dWKO2TdqUjf0mtYjAaGBzGcTQYTVicNDQE06d7XzKm1rTMZNU76/oGuHr9cwBc83/uzewCMY38U8R3ZKjgnwlVlxLmsa9PFWpNywRV6RT/0T1xxAuNI6/uTn09fxAF/5QIuk2dvfQA54PSPilXdSmhgn+gWtIyQXt/NJsxPDLK3pneoGrWwX2pr+cPorRPSgSNDG/u2+E90Mg/Nfz6zlRdSqjgH4qgNFFxK8u9M7y5tFkHvM9Vmuv5gyj4p0TQm++lg4VuJRr5p0JQbv9dp8+tLmet4B+KVd1dXH/ZMro6WzGOloN2Ff6xHW6ZwUhTM8ce9IJ/muv5gyjtkxJBt6mzji9UJWjknwpBd3D3PLWT6y9bVjFn/dorr/L0oelcvvquTJQbxikoTVRMr+6d0casg/tSX88fRME/JYK2qLzqPy6D61HwT4nJcvuVctbr+gZYtn03u+a8YdxdQ1GWatDjUjpxvHdmG/OOHKh9gVhKKPinRFA1w8Vvnu8doLRPKgTdwVWTVrjh7qf54cF97Jt+9NjhkVE+v34zBw8fyVTHyTiN/SN8x0KWnDADMvo3VPBPkcCR4cyZGvmnRNAdXDVpha2Dw7QdGh4X/AEGh0cmHJvVCpWo+JZRlzV3yxoF/yxQW+fUqKdNwMKOmbQdOsDr06ubfMxihUq9/II84FtG/QfMoGtwa5yXGykF/yxQW+dUmWqbgNXnL6IJx74ZR4N/a0szM1uaeG3/xNF/FitU6hG0VmZmS5PvJPzDexxdhcZuWaTgnwUK/rnwh0u9Es9pHR0YBI5cIf0dJ6MQVGlV/lzRdpuptI8knNI++VBo6nbtB3q49k8n7kilap/J1ZoGc52dXqv0Q4cy2UtJwT8LNPLPh0k6emat42QUgiqtOltbxlVLgXfndG73UvglXk//uXMbeKWNoRW+WaCRfz6onXNdglo6fH7lmb6rfd+ybLF3UEZTPxr5Z4FG/vmg4F+XSpVWE+6c+ju87wr+wczsPcBXgWbgJufcmrLXPwrcABQbY3/dOXdTGOcWvOCvkX/2FYO/9u+dsprSYz4bukShqk18IlB38DezZuAbwIVAP/CQma13zj1RdugPnHOfrPd84qO9XSP/PNDIv7HK9vGNQsUdxSIURs7/bcBzzrkXnHOHgNuAS0L4vVKtYtqn0I5WMkrBv7EKwX/Nrb8Z1347TJNt4hO1MNI+XcCWksf9wLk+x/2RmZ0PPANc45zb4nOMTEV7O4yOwsGDXqsHyaa9e73vxd3bJFI//f1+LgZGdw9G1kiv6k18IhDGyN98nisfgv4EWOycezPwK+AW319kdpWZ9ZpZ786dO0O4tJzQVo75MDTk/UPfpCK9Rrj+3/oZtSZmHTz6uSo20qtqv+UqVL2JTwTCeBf1A4tKHp8IjGuI4Zx71Tl3sPDwW8A5fr/IObfWOdfjnOuZm8G62sgUg78mfVPHb1evQNrIpaG27jlQ6Ok//nM1ODwSWqqm1o3nwxRG2uch4BQzW4JXzfNh4COlB5jZAufctsLDlcCTIZxXitq9Lec08k+Xmif7FPwbamFnK0MzjhnbyrGSqaRq6mn0V6+6g79z7rCZfRK4G6/U82bn3GYz+wLQ65xbD3zKzFYCh4HdwEfrPa+UUNonlSpN9k0ICAr+DXXtRafx+lfbx7ZyhGga6cW1OjuUOn/n3AZgQ9lznyv5+TrgujDOJT6KI3+lfVIlaKRYvAMovSO45gePsPCx3zNtegsv9w2olUMDrOruYueJ8xjZtTeTjfS0wjcLAkb+cS0ekeoE9ZppNptwR+CA9kPDDMw8Vrt0NdDcRfOYe2AvL67JXiM9Bf8s8JnwjXPxiFQnaFevoBbDbYeGeX16q3bpaqSA3byy0EhPNWNZ4DPhG+fiEanOqu4u34ZiXQG549ItHLVLV4NkeCtHjfyzwCftE+fiEale0Aiy/I4AvLRPcQtH7dLVIJ2d3h314cMwLVvhMlv/a/LKJ+0TlE9W0Gicqc65lJb/DQwOY0Dz6GFmHj7EvumtqZxcTK3S/j7HHx/vtYRMaZ8smDEDmpvHjfzjXDwiR+dcproKdFV3F/evXsFLa97PVz50Nm9s9RbNT+vo4PrLlqU+35waHYW2zhncy1fBPwvMJrR1DsonK2g0RphzLqu6u/j5ld0AfPrSt+r/w0ZqUFvnOCjtkxU+bZ2zUJGQVqHPuaijZzwU/CXxCm2dq80zaw1AtEKfc1Hwj4eCvyReWxvbB3ZNWttfDPjFScRi61WtAQhfUA3/lOdcFPzjoeAvidfezsDArknzzKXBqLznthYOhSv0hl3awjEeCv6SeG1tTBt+xfelrYPDvhOQfsfVQ6mk8UKdc9HIPx7HHusVVGSw2kfBPyva25k1esj3pYWdrVUF9nrWAKidRMQU/OPR1OTdbWVw5K9Sz6xoa2Ne0+HA2v5Kgb3eNQBqJxExBf/4ZLTFg4J/VrS1cczIgcDafr9FX8X9N8NYA6B2EhEbGoKWFm9BnzRWRoO/0j5ZUajzD8ozR71jkNpJREwbucRHwV8Sra0NDhyA0VGv1YOPKBd9hV7aKOMp+MensxN6e+HLX/Z9+bGBPdzz1A72Do8wq7WFd51+Asu6OqZ8HADLl8O554b5v2ICBf+sKLZ1vvZamD69Mee85BI47zwg3r1Ic2HvXgX/uJx5Jvz4x/DZz/q+vKzwNWaD72FVHwfAF78YefA358orvpOhp6fH9fb2xn0Z6fGrX8Gll8LIxL1FI3HoEFxwAfziF405X9I8+SRcfTUcPNiY8z31FJxxBjzwQGPOJ0c5d3TCvcyFX/5Xtu05MOH5BR0z+eVn/kPNx42ZMWPK8ztmtsk511PpOI38s+KCCwLfoJG49FJ4/vnGnS9p7r0X7rsPVqzwJmKj9va3w4c+FP15ZCKzwMV1zx1ows04xud5WPf80NidsKMJAo6La+Gegr9UpXwB1/daZrF027a4Lys+u3Z533/+88YEf0mkoEKHjtYW3w15/P77uKjUUyry601/1ytHvADYqDRT0uza5fV6V+DPtaB9M8yoGPjjLohQ8JeK/BZwbW0t9DzZsSOGK0qAXbtgzpy4r0JiFrRvxuD+4EFRUvbXUNpHKvJbqLWz7Tjvh+3boSuHFT0K/lLgV0Jd7J5brquzlftXr2jUpU1KI3+pyC8vubOtMPLPa95fwV8mkYZtVBX8pSK/N/LezrneD9u3R37+dX0DLF+zkSWr72L5mo1V74MbKQV/mUQatlENJe1jZu8Bvgo0Azc559aUvT4D+C5wDvAq8CHn3EthnFui57eA65pLl8P/IvLgn9huoQr+UkHSt1GtO/ibWTPwDeBCoB94yMzWO+eeKDnsSuA159wbzezDwD8AKlpOoKCe/L5v5NmzIw/+k3ULje2DtX+/96XgLykWRtrnbcBzzrkXnHOHgNuAS8qOuQS4pfDzHcC7zcyQRPEr6bzuzseC0yzz50ee809kt9Bijb+Cv6RYGMG/C9hS8ri/8JzvMc65w8Ae4PgQzi0hqrkn/4IFkY/8gxbBxNotVMFfMiCM4O83gi9vGFTNMZjZVWbWa2a9O3fuDOHSpBY1j7Lnzw8M/mFN0iayakLBXzIgjAnffmBRyeMTga0Bx/Sb2TSgA9hd/oucc2uBteA1dgvh2qQGNffkLwZ/57z+JwVhTtKG2S00tD2GFfwlQJr2sQ4j+D8EnGJmS4AB4MPAR8qOWQ9cAfwG+ACw0SW1nWgOFd+wA4PDGONvySYdZc+f7018Dg2Na04V9iRtGFUToVYNKfiLj8RWpgWoO+1TyOF/ErgbeBK43Tm32cy+YGYrC4d9GzjezJ4DPgOsrve8Eo7SSV7wAn/V2zsuWOB9L0v9JHGSNtQ9hnft8u50jjsupKuTLEjbPtah1Pk75zZQtjWBc+5zJT8fAD4YxrkkXH5vWEeVy9Dnz/e+b98Op5469nQSt3QM9R+kXbu8Mtfm5lTd5ku0kjjomYxW+OZcXW/Y0uBfImiz+IHB4dhW6IZaNVRY4FVzaaxkWiIr0yah4J9zdb1hi8G/rNa/dGk7MG4eIa4AGWrV0K5dMHdu6m7zJVqJrEybhIJ/ztX1hp092+tn71Puuaq7i/tXr6Crs3VCTW8cATLUXiuFkX/abvMlWmno51NKLZ1zrq5SSrNJa/0hWXnQ0Hqt7NoF556byLkNiVfS+/mUUvCX+t6wFYJ/mgOk72Tu2QvHRv7XXnTahK36knybL1JKwV/qM38+vPxy4MtpC5BBax6KcxXTXh/i4pERmDMn1AVoIo2m4C/1WbAAfvvbwJfTFCDLF+n4zVXc8pNeLoaxBV5pus0XKaXgLzUrTYf8j+eG+djOndjoKDQ3+x6flgDpV71TbmR7Yc9ire6VOiRhfYiqfaQm5bXtL0w7FjtyhJ9t/F3cl1a3aiah39h00PtBwV+mKCnrQzTyl0B+o5Py0XFxL98f/OS3vPfC7rguNRRBk9NFrS3NXH7yMd6DGoN/EkZ6kgxJ2aBII3/xFTQ6KQ+OO9pmez9si34v36gFrUyGozXb57QVPrQ1BP+kjPQkGZJS/qyRv/gKGp00mzFa0pB1Z7vX3OwUt68h1xXlCLqqyel/2QXTpo3rYlpJUkZ6kgxJKX9W8BdfQaOQUedobWkeC2Y7j/HSPqsW+E/2hqnalrn1/ANRcXK6uHF7DbuQJmWkJ8mQlPJnpX3EV9AopJj+KC5hnzO3k5H2WZxp0Y/8q+mlE3mKpRj8a5C2hl8SraS0gdDIX3xNNjqZMDq+ZSE8/7z3FaFpL77AScDh5ma2Hjt3bPRdOoKOKsVSvJv4ykPPMG36NF7uG6j69yVlpCfJkYTyZwV/8VXT4qxFi2DDBu8rQv9a8vNfv++vuGPZBcD4EXQUKZbSdNNxw3t5pu2kmnZoStNCN8kPBX8JVPXo5MYb4YEHIr+eTS/t5vbefv7L/7uNKzb9hDuWXTBhBB3FZFrp3cRxw3t5rXVWzXcTSRjpiZRS8Jf6nXyy9xWxc4AtfQOs+/wh/nr913j36y/zh1euHBdUo0ixFO8azB3huOEhdrfOGve8SBop+EuqrOrugu9+ARbexLcPPgzdfz7xdcJNsRTvJmYd2EezO8Jrx8wae14krRT8pS6xrFzt6IDLL4dbb4Uvfcl7XCLsFEvxbmL27r0A7G6dpQlbST2VesqU1VtWua5vgOVrNrJk9V217+179dWwbx98//tTu/gaFEvzTm32+vo0zZ2b6B2aRKphzpU3rk2Gnp4e19vbG/dlyCSWr9noO7na1dnK/atXTPrfli/YAi83X1NQPeccOHwYHnmkpkVXU7Z+PVxyCfT2eucWSSAz2+Sc66l0nNI+MmX1lFWGUo//8Y97X5/4RE3tFqZs82bvuzp6SgYo+MuU1VNWGUo9/kc+AjfcAN/5TvX/Tb1OPZX120f5h1s3qmZfUk3BX6asnrLKUOrx29vh2WerPz5ALZPW1fYXEkk6TfjKlNXTo8SvfXIcFTS1TlpX019IJA3qGvmb2WzgB8Bi4CXgj51zr/kcNwo8Vnj4snNuZT3nleSYalllUloe1Dr3oA6dkhX1pn1WA792zq0xs9WFx3/jc9ywc+7sOs8lGZOElge1BvOk9GIXqVe9wf8S4J2Fn28B7sU/+IskSjHPH1ToHBTMK81zaLtGqUcj3z/1Bv95zrltAM65bWZ2QsBxM82sFzgMrHHOravzvCJT5rfGoNRkcw+Tpas0GSz1aPT7p2LwN7NfAfN9XvrbGs5zknNuq5ktBTaa2WPOuQnN383sKuAqgJNOOqmGXy9SPb88f1FXFaOtoHSVtmuUejT6/VMx+DvnLgh6zcxeMbMFhVH/AmBHwO/YWvj+gpndC3QDE4K/c24tsBa8Fb5V/S8QqVFQPt+g4srkqfxeTQZLNRr9/qm31HM9cEXh5yuAH5cfYGbHmdmMws9zgOXAE3WBLfKoAAAF6ElEQVSeV2TKotpWUds1Sj0a/f6pN/ivAS40s2eBCwuPMbMeM7upcMybgF4zexS4By/nr+Av49TV5K1GUa0xSMraBUmnRr9/6prwdc69Crzb5/le4D8Xfn4AWFbPeSTbGj3RFdUag6SsXZB0avT7R109JXb1dAcVkfGq7eqp9g4SO02UijSegr/EThOlIo2n4C+xC2Oiq5ETxiJZoJbOErt6J7q0slakdgr+kgj1NHnTylqR2intI6mnCWOR2in4S+ppwlikdgr+knpaWStSO+X8JfW0slakdgr+kglJ2BVMJE2U9hERySGN/KWhtM2hSDIo+EvDaDGWSHIo7SMNM9liLBFpLI38pWHiWIylNJOIP438pWEavRirmGYaGBzGcTTNpKZvIgr+0kCNXoylNJNIMKV9pGEavRhLPX9Egin4S0M1cjHWws5W3+0h1fNHRGkfyTD1/BEJppG/ZJZ6/ogEU/CXxAmzPFM9f0T8KfhLomgVsEhjKPhLotSyJaMWcIlMnYK/JEq15Zm6QxCpj6p9JFGqXQWsBVwi9akr+JvZB81ss5kdMbOeSY57j5k9bWbPmdnqes4p2VZteaYWcInUp96R/+PAZcB9QQeYWTPwDeC9wBnA5WZ2Rp3nlYxa1d3F9Zcto6uzFQO6Olu5/rJlE1I52rRdpD515fydc08CmNlkh70NeM4590Lh2NuAS4An6jm3ZFc15ZnXXnTauJw/aAGXSC0akfPvAraUPO4vPDeBmV1lZr1m1rtz584GXJqkVbV3CCLir+LI38x+Bcz3eelvnXM/ruIcfrcFzu9A59xaYC1AT0+P7zEiRVrAJTJ1FYO/c+6COs/RDywqeXwisLXO3ykiInVoRNrnIeAUM1tiZtOBDwPrG3BeEREJUG+p56Vm1g+cB9xlZncXnl9oZhsAnHOHgU8CdwNPArc75zbXd9kiIlKPeqt9fgT8yOf5rcD7Sh5vADbUcy4REQmPVviKiOSQgr+ISA6Zc8msqDSzncDv476OEMwBdsV9EQmhv8V4+nscpb/FePX8Pd7gnJtb6aDEBv+sMLNe51xg36M80d9iPP09jtLfYrxG/D2U9hERySEFfxGRHFLwj97auC8gQfS3GE9/j6P0txgv8r+Hcv4iIjmkkb+ISA4p+EfAzBaZ2T1m9mRhp7NPx31NSWBmzWbWZ2Y/jfta4mRmnWZ2h5k9VXiPnBf3NcXJzK4pfE4eN7NbzWxm3NfUSGZ2s5ntMLPHS56bbWa/NLNnC9+PC/u8Cv7ROAx81jn3JuDtwF9o9zIAPo3X3ynvvgr83Dl3OvAWcvw3MbMu4FNAj3PuLKAZr/ljnnwHeE/Zc6uBXzvnTgF+XXgcKgX/CDjntjnnHi78PIT34c5143kzOxF4P3BT3NcSJzObBZwPfBvAOXfIOTcY71XFbhrQambTgGPIWct359x9wO6ypy8Bbin8fAuwKuzzKvhHzMwWA93Ag/FeSez+J/BfgSNxX0jMlgI7gf9bSIHdZGZtcV9UXJxzA8CXgJeBbcAe59wv4r2qRJjnnNsG3mASOCHsEyj4R8jM2oEfAn/lnNsb9/XExcwuBnY45zbFfS0JMA14K3Cjc64b2EcEt/RpUchlXwIsARYCbWb2p/FeVT4o+EfEzFrwAv/3nXN3xn09MVsOrDSzl4DbgBVm9k/xXlJs+oF+51zxTvAOvH8M8uoC4EXn3E7n3AhwJ/COmK8pCV4xswUAhe87wj6Bgn8EzMzwcrpPOue+HPf1xM05d51z7kTn3GK8ybyNzrlcju6cc9uBLWZ2WuGpdwNPxHhJcXsZeLuZHVP43LybHE+Al1gPXFH4+Qqgmv3Sa1LXZi4SaDnwZ8BjZvZI4bn/VtjURuQvge8XtjV9AfhYzNcTG+fcg2Z2B/AwXpVcHzlb7WtmtwLvBOYUdkb8O2ANcLuZXYn3D+QHQz+vVviKiOSP0j4iIjmk4C8ikkMK/iIiOaTgLyKSQwr+IiI5pOAvIpJDCv4iIjmk4C8ikkP/H1GXj68TCqa+AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x2098aa527b8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "model = DARTRegressor(base_estimator=CARTRegressor())\n",
    "model.fit(data, target)\n",
    "plt.scatter(data, target)\n",
    "plt.plot(data, model.predict(data), color='r')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ml_models import utils\n",
    "\"\"\"\n",
    "DART分类模型，封装到ml_models.ensemble\n",
    "\"\"\"\n",
    "\n",
    "class DARTClassifier(object):\n",
    "    def __init__(self, base_estimator=None, n_estimators=10, dropout=0.5):\n",
    "        \"\"\"\n",
    "        :param base_estimator: 基学习器，允许异质；异质的情况下使用列表传入比如[estimator1,estimator2,...,estimator10],这时n_estimators会失效；\n",
    "                                同质的情况，单个estimator会被copy成n_estimators份\n",
    "        :param n_estimators: 基学习器迭代数量\n",
    "        :param dropout: dropout概率\n",
    "        \"\"\"\n",
    "        self.base_estimator = base_estimator\n",
    "        self.n_estimators = n_estimators\n",
    "        self.dropout = dropout\n",
    "        if self.base_estimator is None:\n",
    "            # 默认使用决策树桩\n",
    "            self.base_estimator = CARTRegressor(max_depth=2)\n",
    "        # 同质分类器\n",
    "        if type(base_estimator) != list:\n",
    "            estimator = self.base_estimator\n",
    "            self.base_estimator = [copy.deepcopy(estimator) for _ in range(0, self.n_estimators)]\n",
    "        # 异质分类器\n",
    "        else:\n",
    "            self.n_estimators = len(self.base_estimator)\n",
    "\n",
    "        # 扩展class_num组分类器\n",
    "        self.expand_base_estimators = []\n",
    "\n",
    "        # 记录权重\n",
    "        self.weights = None\n",
    "\n",
    "    def _dropout(self, y_pred_score_):\n",
    "        y_pred_score_results = []\n",
    "        ks = []\n",
    "        for class_index in range(0, self.class_num):\n",
    "            dropout_indices = []\n",
    "            no_dropout_indices = []\n",
    "            for index in range(0, len(y_pred_score_[class_index])):\n",
    "                if np.random.random() <= self.dropout:\n",
    "                    dropout_indices.append(index)\n",
    "                else:\n",
    "                    no_dropout_indices.append(index)\n",
    "            if len(dropout_indices) == 0:\n",
    "                np.random.shuffle(no_dropout_indices)\n",
    "                dropout_indices.append(no_dropout_indices.pop())\n",
    "            k = len(dropout_indices)\n",
    "            # 调整对应的weights\n",
    "            for index in dropout_indices:\n",
    "                self.weights[class_index][index] *= (1.0 * k / (k + 1))\n",
    "            # 返回新的pred结果以及dropout掉的数量\n",
    "            y_pred_result = np.zeros_like(y_pred_score_[class_index][0])\n",
    "            for no_dropout_index in no_dropout_indices:\n",
    "                y_pred_result += y_pred_score_[class_index][no_dropout_index] * self.weights[class_index][\n",
    "                    no_dropout_index]\n",
    "            y_pred_score_results.append(y_pred_result)\n",
    "            ks.append(k)\n",
    "        return y_pred_score_results, ks\n",
    "\n",
    "    def fit(self, x, y):\n",
    "        # 将y转one-hot编码\n",
    "        class_num = np.amax(y) + 1\n",
    "        self.class_num = class_num\n",
    "        y_cate = np.zeros(shape=(len(y), class_num))\n",
    "        y_cate[np.arange(len(y)), y] = 1\n",
    "\n",
    "        self.weights = [[] for _ in range(0, class_num)]\n",
    "\n",
    "        # 扩展分类器\n",
    "        self.expand_base_estimators = [copy.deepcopy(self.base_estimator) for _ in range(class_num)]\n",
    "\n",
    "        # 拟合第一个模型\n",
    "        y_pred_score_ = [[] for _ in range(0, self.class_num)]\n",
    "        # TODO:并行优化\n",
    "        for class_index in range(0, class_num):\n",
    "            self.expand_base_estimators[class_index][0].fit(x, y_cate[:, class_index])\n",
    "            y_pred_score_[class_index].append(self.expand_base_estimators[class_index][0].predict(x))\n",
    "            self.weights[class_index].append(1.0)\n",
    "        y_pred_result, ks = self._dropout(y_pred_score_)\n",
    "        y_pred_result = np.c_[y_pred_result].T\n",
    "        # 计算负梯度\n",
    "        new_y = y_cate - utils.softmax(y_pred_result)\n",
    "        # 训练后续模型\n",
    "        for index in range(1, self.n_estimators):\n",
    "            for class_index in range(0, class_num):\n",
    "                self.expand_base_estimators[class_index][index].fit(x, new_y[:, class_index])\n",
    "                y_pred_score_[class_index].append(self.expand_base_estimators[class_index][index].predict(x))\n",
    "                self.weights[class_index].append(1.0 / (ks[class_index] + 1))\n",
    "            y_pred_result, ks = self._dropout(y_pred_score_)\n",
    "            y_pred_result = np.c_[y_pred_result].T\n",
    "            new_y = y_cate - utils.softmax(y_pred_result)\n",
    "\n",
    "    def predict_proba(self, x):\n",
    "        # TODO:并行优化\n",
    "        y_pred_score = []\n",
    "        for class_index in range(0, len(self.expand_base_estimators)):\n",
    "            estimator_of_index = self.expand_base_estimators[class_index]\n",
    "            y_pred_score.append(\n",
    "                np.sum(\n",
    "                    [estimator_of_index[0].predict(x)* self.weights[class_index][0]] +\n",
    "                    [self.weights[class_index][i] * estimator_of_index[i].predict(x) for i in\n",
    "                     range(1, self.n_estimators - 1)] +\n",
    "                    [estimator_of_index[self.n_estimators - 1].predict(x) * self.weights[class_index][-1]]\n",
    "                    , axis=0)\n",
    "            )\n",
    "        return utils.softmax(np.c_[y_pred_score].T)\n",
    "\n",
    "    def predict(self, x):\n",
    "        return np.argmax(self.predict_proba(x), axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "#造伪数据\n",
    "from sklearn.datasets import make_classification\n",
    "data, target = make_classification(n_samples=100, n_features=2, n_classes=2, n_informative=1, n_redundant=0,\n",
    "                                   n_repeated=0, n_clusters_per_class=1, class_sep=.5,random_state=21)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzs3Xd0VNXawOHfmT6T3hNS6TUh9N6bSlGqigUromLlE722a7v23kEFARFFQKRKld5b6ISWkAQS0uv0c74/ck2cG5AAkwzE/azlWubkzD7vDJN39rxnF0lRFARBEIS6Q+XpAARBEAT3EoldEAShjhGJXRAEoY4RiV0QBKGOEYldEAShjhGJXRAEoY4RiV0QBKGOEYldEAShjhGJXRAEoY7ReOKiQcE+SnRMiCcuLQiCcN1K2ns6R1GUSyZPjyT26JgQ1mx43ROXFgRBuG4F+9yZWp3zRClGEAShjhGJXRAEoY4RiV0QBKGOEYldEAShjhGJXRAEoY4RiV0QBKGOEYldEAShjhGJXRAEoY4RiV0QBKGOEYldEAShjhGJXRAEoY4RiV0QBKGOEYldEAShjhGJXRAEoY4RiV0QBKGOEYldEAShjhGJXRAEoY4RiV0QBKGOEYldEAShjhGJXRAEoY4RiV0QBKGOEYldEAShjhGJXRAEoY4RiV0QBKGOEYldEAShjhGJXRAEoY4RiV0QBKGOEYldEAShjhGJXRAEoY4RiV0QBKGOEYldEAShjhGJXRAEoY4RiV0QBKGOEYldEAShjhGJXRAEoY7RXG0DkiRFAzOBcEAGpiqK8snVtit4TlFhGSuW7yEjPYv6DSLpPygRLy+Dp8MSBKGa3NFjdwCTFEVpDnQGHpUkqYUb2hU8IDXlPA+Oe4e0o6tpHnmUgzuW8dA975OdXejp0ARBqKarTuyKopxTFGXPf/+/GDgCRF5tu4JnTPl8IfeM1vHvyfUYdXMIb78cxY29ZaZPXe7p0ARBqCa31tglSYoD2gDb3dmuUDvsdgd7diVz801BLsdHDAlk25YDHoqqZsmyTMrp82Rm5ns6FEFwm6uusf9JkiRvYD7wpKIoRRf4/XhgPEBUdND//lq4BqjVKtRqFRarjF5f+ZlfZpbR67QejKxm7NiWzKcf/IxGMlNmlomKjWbyi2MJDw/wdGiCcFXc0mOXJElLeVKfrSjKggudoyjKVEVR2iuK0j4o2NcdlxXcTKVS0atvW6bOOI+iKADIssI3s7LpM6DjVbVtNtvYtuUYO7YlY7M53BHuVclIz+XdN6bz8lMm5n8fx+LZ9enRNp+Xnv0WWZY9HZ4gXJWrTuySJEnAd8ARRVE+vPqQBE96+LFhHE4JZOxDqbz23llG3nOaImsMd9/X/4rb3LDuIGNHvsrCObOZO2MWd4x6jT27Trox6su3fMkOhg4w0i7RBwCNRuKuW0NQywUc2J/q0dgE4Wq5oxTTDbgLOCBJ0r7/HnteUZRlbmhbqGU+viY++fIxDuxPJSM9l5tuDadpsyu/F56VVcDH787mi7fDaNrIBMCepGKe/fc0Zvz0Ij4+RneFflnycgvo0Mz17S9JEtGRWvJySzwSkyC4iztGxWxSFEVSFCVBUZTE//4nkvp1TJIkElrHcePgdleV1AH+WLWfAb30FUkdoG1rH9rFa9i04fDVhnrFWrRqyB+bzRUlJ4CSUid79ptp3jLKY3EJgjuImadCjSors+DvK1U5HuAvUVZm9UBE5foPSiQr359X3s1gT1Ix6zYVMPHZNPre0E3cPBWueyKxCzWqfacmrFpvwWKpvCFZVOxg3WYzHTo29lhcBoOODz57lNC4Xnz2vYqfl3kz7NYxPPLYUI/FJAju4rbhjoJwIfEJsTRPSOSBJ5MYMdgLhwN+WVxC/xt7EhMb4tHYvL0N3H1vP+6+t59H4xAEdxOJXahRkiQx6bnRbN2SwOb1SajVKh6d1Ja27Rt6OjRBqLNEYhdqnEqlolv35nTr3tzToQjCP4KosQuCINQxIrELgiDUMSKxC4Ig1DGixi7UKYcPpfHD9OUkHztDWFggI8b0o9/A1p4OSxBqleixC3XG0SPpvPzslwzsmsvsLyN49G6ZH6fPYeH8LZ4OTRBqleixC3XGT7NWMeFuH4beUL4sdHCQlrdf0vLov35nyM2d0GjUFefKsszunSc5cjidkFBfevVphcmk91ToguBWoscu1BknjqfRvo23y7EGcUYkxUZ+fuXCXlarneeensr0L6ejLlvPrvW/cd/Ytzh1MrO2QxaEGiF67EKdUS8yhKPHC4mOrNx4O+u8DZtDhZ+fV8WxeT9vwkebwRdfxKJSla9js2h5Lh+98xOfTX2y1uMWBHcTPXahzhh1Wz8+/aaAvfuLURSFtAwLr7x7jqHDe6HTVfZhNq3fzR2jAyqSOsCQQYGcyzgnNu0W6gTRYxfqjI6dm/Dgo2N545Ml5OZkotPrGTaiN3fd47oWjAT8ZbVeQahzRGIXPOr0qSy2bj6KVqehV59WhIb6XVV7ffon0LtfPKWlVoxGHWp11S+lPfp0YNYvq4lv4YVaXVmKiYiqR0jI1V1fEK4FIrELHjNt6gpWLv2DAb0M5FkUJsxYzCNP3Eb/QYlX1a4kSXh7Gy76+5FjuvHynmPc+XAK3TroOZ3m5NhJNW99ePtVXVcQrhUisQsecfhQGmt//4Mfp8Tg51v+Nrz1FjPjJ/1Mx85N8PUzXaKFK6fTaXjrgwfZu/sURw6n07WpL8++3hKjUVdj1xSE2iQSu3BF8vNKWLl8D2fPnqdho2j6D0q8rHHgG9ftZ9ggU0VSh/Khie1b69i29RgDb2hTE2FXkCSJtu0biuWDhTpJjIoRLtuJ4+cYf8+7ZJ1eS6uYY+zftpSH7/uA3Jziy2rnYjcwJanqVnqCIFSfSOzCZfvyk/k8fLeRF56ux4ihIbzz7yh6d7Izc/rKarfRs08Ci1aUkV9grzh2MsXMriQbnbs0qYmwBeEfQ5RihMtSVmYl+Wgqg99s5HJ8+JAAHp58ABhZrXaat4hmwJC+jH1oLf17GjBbFNZvtfH4pNvx8a25+rog/BOIxC5clvLhgxIWq4z3X9ZeKS2T0Ru0l9XWPfcPpE+/RLZtOUawTsPUCS0JDvF1c8SC8M8jEnstk2WZndtPkLT3JP7+3vQbmEhQsI+nw6o2vV5L524JfDvrJE88FI4kSTgcCt/OyqHPgN6X3V5sXCixcaHuD1QQ/sFEjb0W2e0OXpj8HTO+/p4g3U7Op6xh/Li32LP7pKdDuywTnxzO3qN+jH0olVfeOcuIcadRdI24/c7eng5NEAREj71WLV+yG8mayvefx1bMeOzfs4jX35zNrF9eRKW6Pj5n/QO8+GzKE+xPSuHc2XxuvjOCRo0jPB2WIAj/JRJ7LdqycS9jhvlWJHWADm198TLmczz5HE2bRXowussjSRKtE+vTOrG+p0OpFlmWMZttmEx6MZxSqPNEYq9Fao0au8N18LaiKNjt8gXXNBGunizLzPlhHb/+8gdWi4WAQD/G3T9UbJcn1Gkisdei3v06MHvuXLp19MVoLB9RsmpdPpLGj4aNwj0cXd00e+Za9m5ZyzcfhBMdqefA4VJeemsORpOert2beTo8QagRIrHXon4DEti/7zhjHthN904GzmXJHD+t4vV3x4vyQDUoisLK3/ex+vetWMwW2ndOYOSY7hdd8MvhcLJw3jqmfRxBVL3y5Q4SWnrz1AQHs39aXSuJXVEUnE7ZZVs+QahpIrHXIpVKxaRnR3PyRA/2J6XQvIMXL3Zvhl5/eeO//6m++HQRx5K2c+/t/vj6qPlt+QYmPZbEJ189hlqt4uCBMyiKQnxCLFqthpJiC4rTVpHU/9SsiYlzGdk1GuufJaCF89ZRVFhKXINw7ht/M53ErFqhFojE7gENG4WL0stlOncunz9WbGHBjDi8vcp7v4nx3jz9UjpTv1rBpnW7iAiRkVA4m61m8gt30bZ9A/QGE8dOlNG0UeVs1p17iqnfKKpG4/3+25Uc2r2Br94JIzY6km27inj9zWm8+MYEElrH1ei1BcEtd+wkSZomSdJ5SZIOuqM9QfhfRw+n0661sSKpQ/nInI5tdMyb8ztvPOvLtE+j+e7TGN78ly9vvjqdokIzd957Ey/8J5OtOwvJzbOzbFUuX0wr4o5xg2osVovFxpKFG3jj+XrUjzWiUkl07ejHI/f68cucNTV2XUH4k7t67N8DnwMz3dSeILgIDPLmTIYdRVFc7kf8sbGIXl1MtEusnL3bJsGHbu2LWPfHAW4Z0RmTl5GvZq0iK+scjRpH8/KbY2nZKqbGYs3PK8HLqBAa4rq+e8tmJmbMy6yx6wrCn9yS2BVF2SBJUpw72hKEC4lPiMUp+TPzp2zuGB2CWg2795WwbY+Fe24NqHJ+cKCKkmILAH36xdOnX3ytxRoU7EupRUXGOSuREZX1/T1JJcQ1qNkSkCCAWFLgH0W5jndwVqlUvPHug2zaG8Dg208yYtxpXv+4jKcm3876rTbMZmfFuRaLzJqNZtp3bPQ3LdYcnU7D6Nv789xrZ9l/qISSUifLV+fx7ewSbruzv0diEv5Zau3mqSRJ44HxAFHRQbV12X+88tEZ61m8YD25ecW0aBnLfQ8NvW5mjP6V3ebAy9tAQaGD4hIVA25syqAbEkk9ncH4p/cyepg3kiTxy6JiEtq19ehM3tvu6I1Wq+PJF5eSm1NAeL0QHn/mjutqdnFdc+xoBmtW7cFus9OlWys6dGpcZ4cZS+7qxf23FLNEUZRWlzo3sW0DZc2G191yXeHvTfliCScPbeWZiaHERBn4Y1M+739ZyFsfPEbjpvU8HV61FReVMX7ce4wdrmb44CDMFpmpM85zPD2EDz57mM2bjrJh7R4AevZtS9t2DVi9IomDScfwC/TnxsEda3U9m4L8Up6e+Bn1I8vo1tHAidN2Vq638erbD9VofV+4sLk/bmDhL8sYfqMJo0Fi0coyGjVvzTPP33pdJfdgnzt3K4rS/lLnieGOdVhJiYXfl2xi7rfRBAaUj5Uf0DuQ3FwH8+eu47mXxno4wupb+fs+2sbLjB1VnpyNRjXPPlGP28encvhQOj17taRnr5ZA+fOe9NgX1AsqoH8vL86eO8PzT2/j0afG0qtv7dTaf/xhLR0TLDzzWGVNPbFVPp9+MJevp026rpLJ9S47u5CfZi1j9tdRFTe0hw+Rueex/ezd3aFO7nvrruGOc4CtQFNJktIlSbrfHe0KV+d8VgEhQeqKpP6nhJZenEk566Gorkx6WiYJzV1HmUiSRKtmetJSXScb/bZgKw2jCnnv1Shu6BfIfXeG8eFroXzxyS/Y7Y5aiXfn1gPcfJO/y7E+PfzJzjxPXm5JrcQglNu14wRdOhhcRikZDCpu6mdk29bDHoys5rhrVMzt7mhHcK+wcH+yc2Xy8u0uyX3fwVJi4i7cS1EUhT27TrJu7V4UWaZbz9Z07trU4z3MuPr12LMriVE3Vx47frKMZavyuDOmBJvNgU5X/nbes+Mg94/xdYm5eVMvAnxyOH0yiyZ/qXPbbA5+X7qb7VuTMBgM9B/U0S3PV2/QUlrqdDlmsyk4nKDTieUFapPBoKOkVK5yvKRUwaDXX+AR1z8xKqYO8/IycNOwHjz3WgbJJ8uwWGSWr85jxtxSRt3W+4KP+earZXzxwXc0j0wmPu4E3381kw/fnefxETX9ByVy6LiWb2ZmkpNr47lXTnDvo4fp0k7NoV1ruPf2t0g5fR4Ak7eJ/ELXnrnTqVBY7ESr17Bp42EWLdzB0SPp/GvSN+xYv5jh/fLp3jqdbz+fwfRvq78p98Xj7cq3s/KwWssTiqIoTP/xPK3bNhN7utayTl2acOCozN79xRXH0s9aWbyyjL4DEj0YWc1x283TyyFuntYeWZb5Zc5GFv26ntycYlrGx3Hv+CG0io+tcm7K6fNMfvwDfv42Fl+f8t6v2ezkjglnmPzvhz1+0y8rq4Dvvl7KsiU7iQ638tX7cURG+qNSSfy2LJe5y3R89d3TbNp4hJlTZvHlu1EEBmhRFIVZP2ezfL2BsjILESFWYqPUrFlfiMViY93ieDSa8j5OQaGD0fel8vWM5wkN9bviWJ1Omfff+pnd2/fRrrWRkyk2NPoQ3nj3AQKDrp+tEOuKPbtO8p9XptOkgQqTUWJ3kpUHHxnB4GEdPR3aZanuzVOR2K8DxcVmVv2+l9SUs0THRDDoxjY10uubN3cL2adXMvkx19EyX08/i03fnXsfuLIx2CUlFnKyiwgL98dorKxzFuSXsmf3SbQ6DR06NsJg0P1NK5UmP/kVtw0uple3yhq2LCvcfNdp3vnk/4iKDmLm9NUs/GUNrZrpOZdlR60PwulQGDPYxshhwQCcOZPLi/9JY8TQSG4ZElLR1r9ez6BT3xH0G3D1a7afSc3m2NEMwsL9iU+I9XhJ65/MYrGxc8cJ7DYH7Ts0wtfv+vvmJEbF1BHnzuUz6dHPSGzpJLGljv2HDzD+p9W8/9lEIqPcOx/AZNRRWFT1g76gCIJjLr8W6XA4+fqLJaxevpXAADUFhQojbu3HHXf3ZfHCHbz/+hJ8CcCpOLAbynj387to2+7SIxScDgc6nWuClCTQ6VTY7Q4kSWLcfQMYdksXjhxOJyDQCz9/L56a8A7DhzQA4Hy2jaSDFnp2MbLij1yXxH4+x4mvmz44Y2JDiIkNufSJQo0zGHT06NnC02HUCpHYr3HffLmYkTdJ3DO2/IbfyGHww9zzTP1iEa++da9br9WjVwu++WoB+w6UkBjvDcCxE2Ws2WhlyoyEy25v5rTVnDu5k/nTY/D305B53sbkV9ZgtTqZ8eUm2hg746XzAiC7NJvJE2exeO3zLr36C+ncvQ3zFq+gUztfVKryBL9tVxFOxURc/dCK8wICvSvWXE87k4NWo0JRFN78IJUlK7Jp3EDLqVQbRcUKhUUOfH3ULFmRR06BoU4OgRP+OURiv0ZkpOeSk1NEgwZhLmWW7VsO8dJE13r4iCFBfPX94SoLYl0tH18TL7xyH8++PoP60floNBJHTzh5avIdhIRcXr1ZlmUWL9zAzM/D8fcrf5uFh+p4+uFgHvvXaoKc9SuSOkCIVwgZRd7s2JZMrz5/P8dt2PBObNt8gAefPEOf7gbSzzn5Y7OdF19/4KIbgkdFB6Ez+vHca6c4dqyA2V+GEhykobDIwUdfF9D35gM0qB+ASuvHG+/eL7YqFK5rIrF7WHGxmbdem82Jo8eJjtRzKtXGLaP7cfe9/ZEkCZ1eQ2mZs+JmJkBpmRO9XlMj9dr2HRsxe97L7NtzGqdT5uV2DS7Zg74Qq9WB1WIlPNT1sTFRekpKzJiUqm89NWrMZtsl29brtbz78UNs2XyUA0mnCI7xY8qDiQQFX/ympCRJPPXs7Yy79VW+eieYuBgdsqxgMuqY9GgwK9af47FnHiCxbX1RBxeueyKxe9jH780jKiidD3+oj1arIjfPzuP/WktkVCj9B7am36COTPl+Ny8/E4lKJSHLClNmnKfvwA41FpNer73qnX4MBi1RMWFs311Elw6Vvf0NWwpp0aoBybszccoNUKvKx3Sb7WYKlFw6dGxcrfbVahU9era4YM1UlmW2b01m965kfHxMDBjUhnqRQcQnxKJWq4iJ0qMoKtRqCZVaRb0wGUlyEhkdKJK6UCeI75seVFxUxq7tB3jsgTC02vJ/iqBALRPuCWD5oo0A3PfgDWQV1mPM/Sm8+u45bn0glTPnw3lgwk2eDP2SJEni/odu5rX3c/l1SQ5Hksv4Ye55vppRwpPPjKbXkCbsLNzCibwTHM07xp7SbTw2+Ya/7XVXh9Mp8+/nZzBr6kyi/ZJwFmzmsfHvsWFd+R4wOr2RpatL0WrVqNUqJGDj9jIcDjXBwb5ueOaC4Hmix+5BJaVWvEwqTCbXmYihIToKC8snU5hMet756CGOHkkn5fR5Bo0KoUXL6OuiZ9mpSxP+/dYjzPvpD+Ytz6J+w4a892lf6jcI4+U3xrBz6HHWrzmM3qBh0OAhbln5cO3q/ZTmn2D6Z7FoNOWv0aB+ZTz2/Bw6dXmF4aN68cP8lRQXK3TrZOREio1vZhXSZ1CXi9bnBeF6IxK7B4WF+aHSeLmMQgFYubaQhLaVY6glSaJ5i2iat4j2RJhXpWWrGFq+Ma7KcUmS6Ni5CR07u3dz560b9zFisE9FUgdo1thEXJSKQwfOMPGpYWRl5bJ49UE27igmv0ChVWI8/37tDrfGIQieJBK7B6lUKiY8NpJ/vTGTO0eVUT9Wz8ZtpWzeqeLjr/t4OrzrklanxWJxVjlutshotGoMBh3vffwQJ46fIzXlPDExIdfV8sWCUB0isXtY954tCA19jMULN7M1KZcmLVrz+TddCQj0vvSDhSr6DujA1E8P0q9XQMVIonWbCigqNdAqvnJJhEaNI2p1fXZBqE0isV8DmjSLZNJzYzwdRp3QsXNjDu7vyej71tG1g5GcPJkTqRKvvzNe1NCFfwyR2IU6pXw0zo3cNLQTe3adpL2vkU5dmqLXay/9YDdyOmWOHklHURSaNY9CoxFL9Qq1RyR2D9qflMKMb5dy7MgZgkN8uWVUX24e0fm6GPFyrYuoF8jgYYEeuXbSvtO889pM/HxsSJJEfqGWZ168s1rr4FyJdWsPlO9pm1tIy/hG3H5XP6Kig2vkWsL1QXw39ZBjRzN4/YUpjBxUzPKf4njt/0ysWbqY2TPXejo04SoUF5Xx2ovf8fwTJn74KpZZX8bw70levPHyNAoLSt1+vXk/b2LWNz9y5y1lvPeSNw1Ck3nq0U84dzbP7dcSrh8isXvI3B/X8MAd3gzqG4jJpKZlcy/eeimC+T+vwWq1ezq8a5qiKBTkl16Tr9OGdYfo0FrtMtu2Q1tfurXX8MfaA269ltVq58eZy/nw9Xr07OpPgzgj94wNY/gNWubOWefWawnXF1GK8ZDU0xmMv9V15EtEmB4fE+RkF13xkrxnUrNZsmgbOVm5NG3egJuGtr+qtdtlWWbunM3Mmb6ZvNwSWreL5bH/u8ktk4muxLYtx/jmy1/Jy81DllX0GdCBCROHVnst95pWXGwmJKhqKS00WKKo0OzWa2Wk5xIcAFH1XJdU7t7Jh3e+Ou3WawnXF9Fj95Co6HAOHnb9an4+20ZJGVe8w86uHSeY9OhH+Kv3MrBrJmnJa3j0wY/Jz7vyzZO//nQF3769kXpFTens1Yu8nTom3DW1Yhu62nTsaAYfvvU9k8ZrWDWvIfOmRWMtSOLDd36p9Vgupk27hqzfYsViqdxj02aTWbvJStv2Ddx6raAgH7JzHVX2Vj1x2kxomGvH4Ny5fJYt2c2G9YeuyW86gnuJxO4ho27vy5SZRWzcWoAsK5xKMfPim2cZckvPK1pNUZZlPv9oLq9ODuShe8IZ0DuQV5+NpHt7Gz/OWktxURmyXHVD379TXGzm55lbSPBtS4AxAK1aS5x/LEGWSGZ/v+Gy2lIUhdOnsjh+7Oxlx/GnhfM2MG6MF53al29UHRig5YWnIti1bT852UVVzs/KKmDGtNV89N48Vizfi83muECr7tW0WSTxbROZ8H9nWL46j+Wr85gw6QyNW7Ry+9aCfv5edOnehjc/PkthUflzO3i4lG9nF3HzqF5A+es+beoKJt7/Nkd2Lub3+T9z15g3OHok3a2xCNcWUYrxkFbxsTzz4v18881innk1mYAAL24ZNYDb7uh1Re2dP1+IpayIju3quxzv3knL/Y8vZdXyTZhMXowddwNDbu5UrTYzz+Wjw4BBY3A5HmQI4tjBjGrHdupkJm+9OhNrWT56vYoyi55J/7rjgptZnD6VxYK560k/c47YBtGMHNOT6JjyER7nM7Np3Nfocr7RqCayno7s7EKCQyoX8dqz+yRv/vs7BvXW0yJKw/rlSSxasJ53P56Al5fr83G3Sc+NZv0fLVizdjeKojD01hvp2z++RkY7PT5pBF9+qmH4uJ0YDRJqjYkHHrmdxDbl74NdO06wae06fpkWi59v+Z/7hi0FvP7idGbOfUGsO19HicTuQR06NaZDp6eRZfmqJ8+YjHrMFhmrVcFgKE8g+XklnDtbRGIrAz9MbcTR42W89NavGAx6+g+69O7s4eH+WBULNqcNnbryW0S+JZ+WTcOqFZfN5uD5SV/TKdFOxhkZh0OmTTuJ11/6lqkzn3PZwOPA/lRee34KY0cYubm3iT1Jh3j6kd385/1HaNIskoZN49i2ay/tEitLVTm5dtLPOlyG98myzCfv/cQrzwRU3MQcOUzhpbcyWPDLZu66p98l43Y6ZRbO38Yfq7bhcDjo1C2R0bf1xNvbUHGNvNwSfP1M6HSuf0YqlYo+/eLp0y++Wq/R1TAYdDw9eRQTJg4tr++H+Lq8l9au2sVtt3hXJHWAnl39+Xb2GQ4eSKV1Yv0LNStc58TH9TXAHTMiff1MJLZtztQZWchy+b6l6WkFLFhawuhbyreLa9bYxOSJwcz7ebXLY202B6Wllipt+viauHlMe5IK91JiK0FWZDKKMjivO8PYe3pUK65tW4+Rn5PPsnkOsnfHULQ/jkU/KpQVlrBy+V6Xc7/9ciGTHvHl7tvCaJPgw/13hTFhnBfff7sMgOGjerB4tYMZc7JIP2tl554iJr2cwbCRffDxqezJp6flojhK6Ny+sgcvSRIjh/ixY0tSteJ+67Uf2bVhKRPHyfzrUS15aZt49smvsdkcLFu8m8G932RE//cZ2PU1vvhoGU7nlZWX3MVk0hMW5l/lvWSz2TEYqr6/jAYJu63qmjpC3SB67HXIU5PH8OqL3zPynhRio3WsXJPB+LsjGHpDZW+2aWMTmWfTgPIx11988hubN+xDUWQaNorikSdHuYx4eXLyUPwDvJgzYxPFuWZatIri4+fupXGT6i2cdexoBplnZG6MaoPmv5tqhCtBrD27naR9p7jj7t5AeQ/4yJEz9H3PdbXH/r0D+GjKCQAiIgL46IsnmDV9JfMmJ+Mf4MOQEcO5cYjrpu06vQazRUaWQf2XCZ8lpU50+kuPEDpx/ByH9x9k/vQ4dLrypNiyuYmJz6Yx5cvfmTdtFy2NifguegjSAAAgAElEQVT7+2O2m/lt2n4AHn3q2lsjv3O3BH6df4xBfQMrVrxMPlnGqTMyrRLcW/MXrh0eSew5ZQV8s/83T1y6zms+PpKcND+Kc814Hy2idQcNedbCit+v31KMFCIxNWkhy97fSdfGNqZ/E4jJpGLdhmwenvgOQ17shk/gX+rQXWFMl0QUGVRqiZ3sZ+f+/dWKZ0v6cXRmf2RFwaH8t4eogK40hEPnUyveB4qiYNPYOZiSTWS9yrLPidMWbHqby/sl+BYv+tzSBoAMzvLtgUVVrmvzU/PNT6mMGuEPgNks8+mMTPw6Naloq1FYEH3Culd57NHD6XRuZ6xI6lDe4+/eUc+HX2ynvqY5/obydo1aIy19Evhl9lYefHRglbKMp/Xtn8DGdXu57/HjDOpjJDdPZulqM48/PfaaGSIquJ9H3oVKUQDyils9cel/hKD//meP3csrr8zk6QleNG3oxb5DxXz6jYMO7R8ha7YRzu1g4r+boJIlKIG+bYM40vksR6aE0b7tDRds+3ILDsHZSThUP5J2RiYwQI1KgvxCJxbZQTNlEPKKvhXnNgr14dMP1/HSpGC8TBqKih189nERDQNHu5xXHd0bDmDmjM9ZsayY2CgtO/aaCQ3pQ5PS0cgrVJy0luC8cRHH0hdUeWyYjzenz1QdEng6zYm5zIGvt+tOS0atEbkUigrLXG7gXgvUahWv/GdcxVaB3r4mPpnS5ornSQjXB48kdpNBS2LjcE9cutaUmctYunIB23dtRlGgQ7vODLthFCajV63FkNj4Rg7GRTDjp/nk5J0lPDSOu0Y+QeP6zdh7cCetmnnhbXLttbVuYSItvdht/z4t4wLZvmURJQU2yPdBksCsciAFKYwadCN+vgEV58Y3GMevyxXufHgdkRE6Ms7ZSWw5iGEDb7uC+xDhdE/8nORTRygqLmDiPQ0JDa58TokAaU9XedTbhdsZ0HwTh86a+WLWaYYP80eths1bSli12UyTduFk7Eonzj8Wvbq8rl9gKcDLT0dgkDfZ2YUcO5pBSIgfTZrWuybW/VGpVHTp1owu3Zp5OhShlkiKotT6RRs3aKZ88tp3tX7d2iLLMu988gqlaWXEGBojIXHGcgJdhJrnn34dtdrzX9fP52TyzY//x4JpDV1KDq9/mIZKGkXfbhfusV+JjHNnmPL9pxTlFqGSVGhMGu6/+xGaNqy6ETVAcUkRufk5hASF4mWq/XXp7Q47ew/uZO3m+Vgs59BqVZyzWyjMU+EssWEvttLIqwnhgaHYFQupzuOMeK4FuemFHPzjBJ0S/DmTbsfkG8qrb95/1fu4CsKfgn3u3K0oSvtLnef5DFMHHT1xkNyMHBK9ulX02Jp6tSYpaysHjybRumU7D0cIocHh1I/uwnNv7GT83SH4+mhYtiqXHXv0PPlg9Ua8VFdkRAyvPvceZzPTcMpOIiNiUKsuvoytj7cvPt6eKWkcPLqP72Z9hdauRVYUFL0vnTp2I2vdejrrEvH3CyRTm0ZS4XZS8zPxiYiiXqcxrN1kISB7ETOnhOPrW/7cps/K5OHn3mPwpI7EBPpfsJ4vCDVBJPYakJGZjpfs7/I1XJIkvJ1+nM1Kc2tit1gtHDt5GFl20rRhS0zG6q8Lc+uwh1i7OYbn/7MSq81C4/rteXjc6BopF0mSRGTEtT0KI78wj6nTPqOpKhF/Q3kNOttyjt8W/0Jb/2746wIBhXBTNN5aX46r9/Pus++hUqn4aubLTBgXTpg+gPyC8qGjI2/w4cefTrP7h+6YR6zlRNZvBPgYGVV/oAefpfBPIBJ7DQgJCsWsrro+i1ldQnBg9Sb2VMfh5P1Mnf4ZBocXEhKlqiLuHvsg7Vt3rtbj1WoNA3oOZUDPoW6L6Xq2a99W/J3BFUkdIMQQgcHuTX5JLvY8GVmW0el0+PsFUGItwe6wo9fpsdrK8P/vJKAAY/mIIsWgEOpnYED9cELTnmb+wRM0HrbUI89N+GcRE5RqQKumiej9NZwsPYJDtuOQHZwuPYrkC4lu6q2XmUuZ8t2nNFYSiDd2pJWxAy1U7Znxw1TyCnLcco0LKS4twmJx7yqF14oycykauepOSzpJT07peYyY8FH7obKrSck+iZfJC522/OZzw9j2LF2d7/K4vQeKkGVfggPLJ4g11It9bIXa4ZYeuyRJNwCfAGrgW0VR3nZHu9crjUbDpIkv8dOCGWw7uAaAVi1a8/CIx9Fq3TN2OOnQbryd/vibKnuXPlo/Amwh7N6/nQE9B7vlOn9KSTvFopVTySs4g1OWaBTXjhE3PuCxWnhNaNa4FetWr8WpOFFL5XVyq8NKnnwerUpHrnIePwIpoZAzSjJNAptWlNv6dB3MlzO289LbqfTsYiI1zcaCZRZG3vSM2GtVqHVXndglSVIDXwADgHRgpyRJixRFOXy1bV8rnLKTg0f2cfzUUfz9AuiQ2MVlmN6F+PsFMOHeJ3E4HCgoaDXu3XPTareiUqregFTJaqxWq1uvVVCYz8xfXueZiT707toEi1Vmxs9Hmf7z2zx233/cPqTPYjGTnXeeAL9AvL1qb0RJkwbNaZ7QgqSkrYRJUciKkwxnCiaTF218u3Ki8DBptuOYNN408mqBbK9cLdLby4fH73uTbXs3snDZYXy9Q3nozr4uQywFoba4o8feETihKMopAEmSfgJuBupEYrfbbXwy9R0yT2fi5wzCprKyZNmvPPrQJBrXv/S4YI2mZm5jtGgcz3xpDjanFZ26fKMFh2wnX51Nq2at3Xqt7XvXM7CPhr7dy78dmIxqJoyrx6btpzl95gQNYhu75TqKorB01QJWrF6KHgMWxUznjt25bfi4q3odz53PoKi4kMjw6L/9oJAkiQfunMjexJ3s2rMNtVrNgMRBfD/rK7w0PnQK7V1xbkrpMSJjYl0ebzAY6d1lIFD15qjdYeds5mm0p4tQ4pVrYny7UHe5I+tEAml/+TkdqLIurCRJ44HxACFB7ruBWNM2bltLzqlcEk1dK/4Ysy3nmP7D17zxwoce+5odGhzOwAGDWb3qd0LkekhIZEtn6dS1K7FR7t3QoaAok/ZtXEtIkiTRsL6OvMJcGuCexL5pxx+s/X01icZuGNRG7LKNA1v3YjTMZeTQsZfdXnFpEVOmf8KZlFSMKhNlSgmDBgxh8IDhF02sKpWKdgmdaJdQ+RbO6HuGdSvX0FDXAi+ND1nmDM5rMriv7/hqxXHgyF5+Xf4FgUE2zKcLuH9uBi+8eg8NG4nevFAz3JHYL/QXUmXWk6IoU4GpUD5ByQ3XrRW79mwnXB3jkgiC9eGkFB4jK/scEWGe2SIOYOigkbRoGs/OvVuRZScjW4+iacOWbu8N1gtrxLZdOxkyoPKYzSaz94CZ+2513xDG1WuXU1/bDMN/Z3RqVTqaGONZv3kNt9w05rInds34cSrFp8ro6NUHSZKwOi2sXrGCehFRtI3vWO12hgwcga+vH6vXLqewKJ9GTZrx9ODnq/Vvn5N3ngXLP+LD18KIi9GDrxe7tii8NHkqM35+Hq1WDEwT3M8d76p0IPovP0cBZ93Q7jVBrdYgK1WXN5WRUasvPsmmtjSMa0LDuCaXPvESLFYLh5P3Y7fbaNa4FX4+/hW/65DYnU++Xcyn36Qx7IYgSkqdfDs7m7ioLoSHVm+Vx+ooKi4kRuNa3tKrjNjNdux2+2Ul9sLiAo4eO0RHr74VH3R6tYFoVUPWbVh1WYldpVLRu+sAencdcOmT/8fOfZsYOlBPiybelFnsSJLEjf0D+O33NHZuP0HX7lc3zb+4qAzgkvva2u0Odmw/Tn5uCS3jY6jf4Pr51ixcPnck9p1AY0mS6gMZwG3A5X9vvkZ169qLX07/SLAchlpV/nJlmFMICQ+9rkpKf+fYycPMWfgeCS3UmIzw8TdWenW5g95dypcVMOgNPDLudVZt+JUnXtiOXmegdYvb6NXVvRNtmjRqTtahDGK9GlUcy7FmEh4WgV5/ebsemc1laCVtxeiWPxnURnJL8twSb7XisBQRFlK1AxAWoqawsPQCj6iejPRcPnl/LkcPnwYUmjSrzxP/N6Zit6m/OpOazfP/N4WIECtREWp+mGahfee2PP3sKDFip4666sSuKIpDkqSJwArKhztOUxTl0FVHdo3omNiV5OOH2blzPf4EYZUsqHwUnrj7uTpxA8xqszJn4fu8+3IgrVuWD13MzrVx/5M/0CCmGTGRcUD5NP8RN40DxtVYLMNuGsV7x1/DWWonQBtCkb2ATFUqE4Y/edmvdUhwGFqTlnxrDgG6ymSXaUunbas27g79ohrGxrNy3TpGDK6sPhYVO9i2y8wdD13Z7kVWq53JT37J2JtVfPpq+f2UhUtzefbJL/lu9nMue+YqisI7r//APaNVDB9S/sXaYpF59Nl9rFjekBsHe355C8H93FLgUxRlGbDMHW1da1QqFXffOp7+vW/iVOpxfLz9aNkkocZGu9SUw8kHWLxsPhmZaYQFRzDkxuG0btmOI8cP0LKpVJHUAUKCdIwYbGLvwS0Vib02RIZH8/yk11m1bhkpqacID6vHnX3uIiby8hOgWqVm7Jh7mfb9V4SWRGFSe5MrZ0Ggg/69am9DjFbNE9m2tzFPv3ycmwZ4YVYVsXhRPgNu6kFEvcAranPj+sM0iLZz28ioimOjbwlh2550Nqw7xKAbKz+4MtJzycvJ4uabKl9Dg0HFuFv9mbNku0jsddT1lZ08qF5YFPXCoi594jXocPIBvp76MbFSU9rou1N4Po/vpn3F3Xc/ACgYjVV7w0aDhNNpc8v1FUXh+OmjJB3cg0ajoUObLkRdZN2Y0OBw7hh1n1uum9iyPc889TIbNq8hNy+HXk17071Tn1pdOlmtUnPfbZPZsXczU2atQR1bysQJY+nWvfkVt5mVWUCjuKrlncb1VWSec539arc70etU/O8XHoNehd1Wdc15oW4Qib0GlJQWc+T4QQBaNk2o1URyIb8t/YU4qRlhxvJRHCGGCNRWDQsXz+XZJ19h8Wor57IsRISV17EtVieLV5rp2618ddAycxk79mwmPSOVehFRdGrXvdrL6SqKwpwF09mxdRtBcjgKMn+sXcktN4+mbw/3LQ18MdH1YrljtHs+KK6UVqOlW4feePk3QzXoZ7onXHi54upq0rQe07608YisoFKVZ2xZVti2284dD7iO1ImNC0GWvNi2q6hiY29ZVvhlUSGdu13+zWDh+iASu5vt3LeZJau/pkMbPYoMS1bbuOWGibRpVf1RGO6WcTaNDvreLscCdMEczNmB0WBiYM97ePDp6dxyowkvk8SSVWaCA7rSrFFLcvNzePeTV9CU6PFW/DkqHWX5qkU88/jL1ZpVeTI1mR1bt9HG0A2Nqnz2baQzjgW/zaVtQif8/f5+Bm9d47Q7+Wb/bxfdlq862nVoyByfSF56K507RwchSTB7Xh4qfTgdO7vOKVCpVDz93Fheeelb+nQrJSZSzdpNVtBGcPPI6i0WJ1x/RGJ3o7yCHJat/ZpvPogkNrp8LPbJlFIeefZzGsZ+hq+Pn0fiCgkOpSAvj2B95SieInsB/n4BqNVqunXoR4OYZuzev5nikgJUFHPkyCk+zHgTi60Mr6IAGnpXlg5SSpOZu/AHxt02nt+W/8KefTtQq9V06dSDm/oPx/CXESz7D+0hUA4rT+qKQqm5FLPFjNqpZcP21QwbOLpWXwun04HD6USv01frfFmWycw+i4REeOjV7YiU2DicfWvHsl6dg7PXyirb8g1s2pwG3pcu0ahUKv7z3gP8NHsdL7+3C0WBHr27MfFfvS84yiWxTX2+nv4sq1bsJTW3kFtub0C3Hs3RaDw/XFeoGWIHJTdas2k5RuN8Jj0c7XL8jY/S0Glup0enfh6Ja1fSNmbN/I4m2gT8tIEUOwo5Zk1ixOgx9OxSGVNRcSH/+eBFDEXehOkiMTvL2Ju9lQa+TWnqn1BxnkN2sKV0BaEh4ahytEQZGiIrTlKtxwlqFMjTDz9fkQAXr5zH7t/30MirJedzMrFb7WjQkiwn4fCxct/dj9CtQ68afw0sVguLV84i6cgGZNlOvbCGDOl/H3HRDS/6mJS0k3w74wtKCktQFAW/QF8eGDfxim7mXsrbhdsZ1GrzBX/XNCrU5edYr6BqfQAIdY/YQckDHA4HxgsMtzYZJcxmR9Vf1JL2rTvjGOvgt6W/cKAgDx9vX4YPHUWPzq4bRK/bvAp9sYkm3vEA+BFIE1Vrjhfvp4Fvc7T/LaU4FQd2pw17noN478pRFS00bdl9agOnUo9XTJpqn9iFlSuXEVASit1qx6TyplguwKoy09GrFz/Pm0m7+I4YDMYafQ1mL/iEBrHJzP8uFh9vDX9szuX9L19n4r3vExRQdex3mbmUT79+lyhbI5oZyidhZean88lX7/Cflz52+VbiDs/5dYK0Kitx8HbhdlYcdD02qNVm1NpkAnyqvmZiEw8BRGJ3q1ZNWzNrwc/cOdqBr3f5S5tfYGfNRjPj73DvwlyXq3O77nRq2w273YZWq7tgSSH5+BGC1K6TrgJ9glEXaii2FxCoD0FRFE6Zj1AvPApDtuuSvZIk4UMAZ7PSKhJ7RGgkY0bdydRpn+Et+6NSVJRKRSQGd8ZPF4jR4s2pMydo0SS+xp575vmzZOcd4Ov3GqLRlJcq+vcM5vgpK1t2rmLowNsByM3PYe2G3zmVchKH04q21ECYb+XNyAhjNHllWew7tIvObWtnm7vn/Kome9I68dmpg1UOd++1khmWJVWOi235/nlEYnejyIgY4psO4b4nlnLzDUacTvjtdzMdWo+o8eVbi0uKOJmajEFvoHH9Zhecfi9JErq/1JZPpiSzfNVvnD2XQWRkNGqtilJnMcFUxurj44uzzM5R214C5RBKlCKiGkaT2Lo9axaucmlfURRKKaqyS1SPzn05evwQx7YfI9wQRYghAo1Ki6Io2BUrxhrurefknadhnKEiqf+peRMD+w9nAJCVfY53Pn4FP3MwAdoQThQfgjIJi96MQV8Zn9app7iksEbjrY7HGrSqcmzf2mAOWV137jprstD9Kur5wvVJJHY3G9z/Nlo0ac/+I9sBiVuHdXb7aov/a/WG5SxcPBdfArBjQ+0t8dj4Z/52j9Ejxw/y5ZSPiFIaEKtrRv7hHFLlZJyyAz9bEP66QByygxPmQ3Ts2JWbbxrFuawMwkIiiI1qgMVqZsWqxaSUJBNtKq+xp5iT8Q/3o2nDqsP5+vQYyIGkJAL0wRWjY86aUzEFmGr89YkIi2LBcjNmixOjofKG4c69ZsKCy5cvWLR8HoGWcOp7N634/d6yreTkZVeMuXcqTorUuTSMvfq1eWpCYuNwEi/0izTXD4G3C7cDm4EjVU4d2LRqshcfANcfcfP0OncyJZlPP3+XBH1nDOryhaDOlZ0h1/cc/3nxo4uuBfLG+y/glRVAqLFyEa9Mcxo5PmdxOBzYzTYcioOE+DbcOeaBC26SfT4nkznzZ3D02EEkSUWb1u25bcS4i+6qtGbj7yxcNBdv/LBjRe+nY+L4Z4gIrfkVMn9eNAWkLUwYF0qgv5YVf+Qy51d48oH38PH25dHJ92As9UECQowRhBui2Ju7lYKyPFqFtEOSJM46U2mc0Jjx4x6vE8tJ/K/vDIurHGsbchy1Vk2jsKAqvxPlndonbp7+Q2zdsYFQJaoiqQNEmGLILD7DqTPHaRTXtMpjFEUhLSOF7v/TEwvRR3Ci8CD/fvZt1mxagdVsISG+bcW+nv8rNDicJx56FrvDjkqSLrn6Yr8eN9CpbTeOpxxDccq0bJpw2Yt7XalRQx7gj831+NcbK7FY82gU15aH7x6Nj7cv6zavojC3AJMzAINk4FRZMun6FFoFtGerZhVyrBWVpGJ4x1F0bt+jTiZ1gPstF9jUPO3PHr6rQa02cyZvCV4G1/dGkMkkEv41QCT261yZ2YxGqpp4NWixWi0XfIwkSfj5+lNiK8JXW7k8b7GjCL3OwJvvv0ygMwydYuDI3h9YF7eKJx569qL7tV7Otn8Hj+5j/m9zsJRZkNQSPbv1ZfjgWy97rfXLcS4rg0XL55F8/Ai+vn707zOSrh16IUkSpWUlzF/4Ix18emEpsmKUTIRQj8PWXSQVbWVQvyHcNqLmFj67HlzoBu6+tbGsUlfdNL17r5WcyPrN5Vi/Bk1EOaeWicR+nUts3Y5f9v9IpBJb0ZMsdRRTpir+2y3rBvYbwrKFi2gutcWk8aLMUcIJ20FK5RISDV3xN5R/9Y5RGnLg1A4271x/ReuR/9XBY0nM+XEmTbWt8TMFYnGa2bGuvDc4atgdV9X2xWTnZvHuJ68SbI2kuaE95rwS5v/0E/kFuQwZOJKTqcl4S/4E+4ZRQD7FRYVIqPBW/JEDLIwcervbY5JlmaMnDpF88jDeXr60b935upuBW17PrzogYN/aYOS//Px7aCpO+4Xr+RPajai5AP/hRGK/zrWL78jWJhvYl7yFYMKxKTZyVGe59dZxGA0X33yhb/dBWK0WVq5ZCmUgaSU6dOvMnq278NdV1lMlSSJcHc3efTuvOrGvWLWYGFUj/HTlqxoa1EaaGluzYfMaht0wymXEjrusXr+cAEsocd7lH3JGtQmTxocVq5fSr+eNGPRG7Er5Ymf+vgH4evvhcNixWEpo0rbNRb+lXCmn08GU7z/hxJHj+DtDsKtsLF46j4cffIpmjaqOdLneJDZ2TfaJhF9wfP53hsV8vXtBleNqrZp+DarenBY9/ssjEvt1Tq3WMPGB/yPp8B4OHNyLyeRF5w4PEF0v9m8fp1KpGDxgOAP7DKG4pBAfbz/O52Syc9t2FMV1s2W7YsdbXzXplplL+WPTSvYm7cJoNNGre1/aJXS+aA06Jzeb+hrXETMGtRGsKkrKSgisgcR++vRJArWuMzcNaiN6u4Hs3PM0imuK1kfD2eIz1DPGoFKpcKjs5GrO0bWj+xcP27F3C6cPnaaNV3dUUvmN7TxrNt/N/Iq3//1JjZakriX3W4a67JSsKAonUo4xO+83fg9OIbBxJNJ/FzhrG3Ico+kkMYH+VdoR9fwL+2e8i+o4tVpD2/iOVbZ7M1vKOHH6GFqtlkZxzdBoNNhsVsrMZfj6+KFSqdBqtAT6l8+8rBcWRUBIIBnnU4gylU+bt8s2zikp3NPZdeNmi9XCe5++hi3LSbg2GptsZfapGaT2SWHkkAuXLxrUb8T5pEy8tD4Vx4rsBeiMOvx8q/7R/ikl7RT7Du5CrVbRNr4jkRExKIrCoeT9bNuxEZvNTvu2HWmX0KlKYgwPjyDzbDYB+srZpQ7ZjkW2EOAXiEqlYuL4/+Pzqe+TVXwGjaQjz36eyIhofpo/i2ZNmtO7+0CXrQKvxs7d2whXR1ckdYBAfQgp5mOcyUihfkyjv3l03WSxWpj20zsoyiluaK3j0BY7RStCeWDsi/h4+7JvbSar1DmUBrnuOPXnDdz/Na75kNoK/ZolEnsdJMsyv/0+l8W/L8BH7Y9WrQWjk8aNmnHw0H6QJYxeBkbdcgcd23SteJwkSTx0z+N8OuVdcorOocNAoZJL376DSGjR1uUaO/ZsxnLeRkuv9hU99CA5lD/WraRfjxsuWDO+aeAtvHv4VSiFEH04xfZCzsjJjBlxF2rVhRekWrjsZ9auXUWQMwwkWLVyOUOGjKCktJiNa/8gjGjUkppfDs9hV4ttTLj3KZchnv1738j7+97AZPEmWB+OVbZw3HyAjp06VwzLjAyP5j8vfsTJ1GT2HdzFhnVr0Wf5YNR4sfvMbrZs38jzT7+Gn+/V18HVajVWxXWde0VRkHFeE3voesLKdQto0jCVF56sj0oloSgKU2ae5bcV07lz5BMXrefPX151CObx6FzMZQtQa11fywcTbq6x+K9FIrHXMcWlRbz36WucPnaKQCmUEgrx0vtAicKutJ10rtcXg8ZIoTWPH36Yho+3L80bV9Z2w0Pr8foLH3DsxGFKy0poENv4gmupJB8/QoAU4lJ20ap0+KoCSE0/hb9f1Z15IsOjmfzEyyxduZBTpw8RGB7E/QMeIaF52yrnAqSdTWXt2lW0NXRHqyqvdUc5G7Dwt5+RUejs1a/ieLgSzd4jmzl8/ACtmlYu3xATWZ8JDz7JT/NnkJy9H7VWTY9efRg++DaXa6lUKupHN+Lr7z6mha59xWihIH0YyUX7WbV+GaOGXv0N3i6devDDkemEyvUqJmplWTLw8vMiKuLvy2d11YGj6/jq3dCKteUlSeKOkWEMXbgdp9Nx0fLUyFYX+nbTiH1rY11u4O6O34nTfuF6fl1N+CKx1zFzf52F+YyVeKkzRrUJRVE4ZkkiT8kkXtUZxaGABvx0gUQ5GrJyzRKXxA7lu/5cau2WgIBAMuRMl2OKomCWS/92eeLIiBjGj3u8Ws/l4JF9BDpDKpI3lNfHvex+lFHiclwlqQhwhnDk2EGXxA7Qokk8rz73HmaLGZ1Wd9FtDbPzspCtCr5G17JLiK4eh48chAsM875cbVp14GjXg2zdsp4AQrCprNgNVp6499ka2Vg6IzON/Yd3o1arSWzZocaXtrgSTtmJRuN6X0arlQCFK5k/WeUG7v/U8/+0N/rDC97ANZp0Ver519uKmiKx1yFOp4Pd+3bQWJ2AAycK5TdBQ1WR5Nqz0Kh1yHJlX8ZH68fZ3JNXdK3unfuwfuNaAq0hBOpDkBWZlLJkgiKC/nYp3MuhUquRpap/2YoKZMVZ5bhDsmMyXXgkkCRJF5w9+1deRm/sig2n7ECtqvzTMDvK8PN1z1r6kiQxduR99Ok+iBMpRzEZvYlvllgjI4J+XzefvQd+5ab+BhwOmDp7Lr06302PTtfWzknNG3dh/pJtTLy/crnr337PplFc6xrdW7hN2tNVju07nsnC0NQqxwe12kxwYNVPh2t1NU2R2OsQm91OcXEhpdZSkMHmtI+5pRAAABVISURBVKHXGDBiwoGNUqWIcH1Exfm51kwaxF98rPvfCQuJYPz9jzFrzrecLD2EU3FQv0FD7rvzKbfNzGzTqj2LlyygzFGCSVO+FV+RvQCrvhSjzkCmOZ1wY/k+tIX/396dx0dZ3Xsc//xmMjPZExKSELIQ9lXWsIkCrlC1qEXrhlIXqFVbb7W112trW6335a222tbWrXqLC3qLQkVcWMqqAoISIBBiwhIIJAQSSEKWySzn/pGIxASSkEkmM/zer1deL57J8zzznZD85sx5znNOXRllIYcZ345ZF6OjYhg2dAR52TvoHz4Mq1ipdp+g0OxmztR7ffKavpaclEJyUsdNpXCwaD9bti/itecyiI2p7/K57rtOfvCTeQwbNIZuMWe3kHZHmD71+7zw+k72FBQwbrSN7TlutmbbmDvr9k7Pcqb+/E/CG9/wFz5kL0fLFhIW3nhIbFe4eKtzxQSRBYvfYMniRSR7Mkgyabi9buqo4bAc4HBIIdGOWAZHjCLSFs0RZxFH7Af5xU9/0665WrxeL4ePFuGwhxIX2/RiVnt9+vlq3l4wj2gTD2KotBzj9lk/Ij6uO8+/8ix1lXVYJQSXrY4f3DyXEUOb9u23RXVNFf/75vPk7NpBqCUclzi5+qrrOmV9Vl/6aOW/6Jm8hLtnN1705ck/H8Ai/lv05XRcrjq27NhM0eECusf17JQ5+n3h3ez8RtvHM3MYnZDXZD9f9efrXDHnGJerjnWfriIzfjJbjq6n2ltJuCWKMm8J5SGlPPzAYwjC0hVL2HPsAP0GD+COaXPaPQGXxWIhPDSc0mNHsdvsREZEtXxQG0waN5XzBo9ix1fbsFosDB044uRC2k/88hkKCvfgcrvond6vTVMbnE54WAT33vUzyo4fpaKygh6JPX2+qEZnsFgseDxNPzm53BDm6Hqjb2w2O+NGng+c3+K+XUmTC7i1/drcn3/ht7oufdGXr4U9SNTU1mA8hrjwBCYnT+dgdQHVrkoSJAmbw0Kv1D50i4lrd4v2VC63izcXvMLmLzYSbo2kwnWcAQMGcenk7zBkwHCf9Y9GR8UwccyFTR63WCwdNu47Lrb7yfH9gWjk0LG8PP+fXPddJ0kJ9f33BQdq+GSjkwfmjvJzunPP6frzPz1vE58WNF5dbXTCwibLIULbbsbSwh4kIiOiiIiM5HhtKbH2eFLDM8gq3cCRmmKiHNH86vEHmTj+Qm6cObvZMeNut5udeds4dryMXql96JXau8W+8sUfLWDXphwyw6aw90QuxRWHyN2Qy57sPYTHh3H/3Q+dcU541XGSEpKZPP5WZv/4dS6aFIrLDes2OLnykrv9tqi6amxk/x71I3a+5d2P8lmaVtrosYSkEnITmrb4T0cLe5CwWCxcd+3NvPbay6S6+1FcfYC66jqGWjNJjk9FLJC1fhNJST24dPIVjY49WnaEZ/72BK5yD6HeCCqkjAFDBjH3tp+cttXt9XpZ8+lKhodN4JirlMKKAkZaJmHBQp3TiVQZnn/1WR57+OkOGcanWjZ5wjSGDRpD9q6tOGwWfjp3lM/uoFUdp75759tdPPUFH1pX3PUvLohkjpjAPT96AHtfOOQsoHfEIJJ7pGILsRFisZFhH8TqtSuaHDdv/otEHItlRNhEBkYMZ0zYZPZl72P1Z8tO+1wej5u6OiehljAOnthHMunYxI4FC16vl+SwNKqPV7P/4N6OfMkBp6BwD0uWL2TpqvcpOVrc8gHtFBfbnckTLuGCcRdpUQ9wzd+Q1TxtsQeZQf2GkpHWl115O0mKSG40J4nDGkp1deP5NipPVLBnbz4Twr8ZJWERC2m2vny2YW2T1v3XbDY7aSm9OFx8EI/xYJX67h2XcREaGoqIYCUEl8vVAa8y8BhjeOf9N1m3djXRrm5UeU7w7r/eYvasOUwae5G/46kgoy32IBTqCCU9NYPDNYWNHi+qKWDYsMarYno8HgRBvvWrYBELHk/Tm4BOdcPMWymQXMRqKPLup9ZTjUvqiI2No7zuGG5bHRnpvrlZqaMdKy9j4Qdv8fs//ZZ5b73I/oP7fHr+PQV5rFuzmpi67hyo2Iur2oWnwvCXF56moHCPT59LKS3sQeqm62ZzICSf/KodFNcUklu1lYrIUmZMn9lov9iYbiT37Mmhmv0nHzPGcLBuH2MzJ57xOfplDOThBx9j5IWjMbFu8h3bcUXVsM+ZS47nS267eY5PhiB2tKNlJTzx9CNsXbkNx6FoCjcX8fSzj5Odu/W0x1TXVFFccgiXu3WfSLK2b8bhDKWocj/DLRMZaBnJCMv5pHj68NzLT+OP+0lU8NKumCDVO70fv/z5E6xdv5LiooMM6j2BC8df3OyIiNtumsszf/1vKqpLcXjCqbCWEZcWd9pumFMlJ6Vw2w1zuOX6O8jOyWJn7nYiIyMZP/qCLjkvSXOWLF1EdFV3+jaMH+7uSCLSGcPbC+bx+CN/aDw3vdvFPxe9zvrP12HDjgnxMuOKmS3ewGSxWilzltKT3tilfvihiJAoKewpy+bwkSJ6JPY84zmUai0t7EEsIT6p0dzobrebTVnr2blrO5GRUZw/djLJSSmk9ezF7375RzZlrafsWCkZ6X0YPnhUmxZ9sFqsjBg6xqfj5DtLTm42fUMbT4QWb08k79g2Kk9UNHozXPj+W2xfv43M8CnYLHaq3JUs/tdCYmPjmsyHf6oxI8bzzsL5WE/5k/MYD8ZiiHBEUnua9WmVOhvtKuwicj3wG2AwMM4Ys9kXoZTvudwu/vTikxzec5h4kqgz+1mzZgW33TKHzJETCA+LYMrES/0d0y+iI6OpKakmIuSbu2Zdpg4s4DjlrtO6OiefrF/NqLBvphGOCImil3sAy1Z80KiwG2PYkbuVdZ+twumsZfTIcYweM5ZdG3cR7okAEbzixh5tB4chVcf7Kx9qb4s9G/ge8KIPsqgOtPHLTyjZfYQRERNPdi0kuHryxv+9yoiho32+tqe/GWP4ak8OO3dtIzQ0jMyRE0iIT2p234unTmPB/PlEeWJwWEPxGA/5NdmMHz8JxymzLtY4a8ADjtDGUwxEhERRXN54RsDFS99h9fIV9CAdm8XOB/mL6ZYeS89ByRws3E289MDYPByy7uYHN/6wQ2cxVOeedv02GWNyAJ/N5hfM3G43n25axc68NRhjGDJgCheMvbjT/qC3ZG0m0ZrS6P8q2haLrcbOvsI99O89qFNydAav18s/3nqB7Vu2EudJxG1x8+HH7zF71lwyR0xosv+EMRdytPQIy1YsIVQiqPFUMXz4KL5/za2N9ouKiCYiOpJjNUfpZv9muoESZxH9Bw88uX2svIzlKz5kTOhk7Nb6N4ZE05Ot+9fzvZtuwOv1kJObTUx0LJMm/IieSakd9JNQ56pOayaIyFxgLnDallOwMsbw2jvPEBm5g/vurF9e7e1F85m3YAt33PhQp7wxhoWFUeU9dnLbYzwcd5ZS5aoMiJErbZGdm0X2l9sYFX7ByfH1la50Xp//d4YNGtlkUi8R4bvTZnLJ5OkcPlJEbExcs9PaWiwWrr/mFua9/jIprt5E2WIpdR2mzF7MndPmnNxvT0EeMZa4k0X96+eII4m83TncdsNcJo2b2jEvXilaMdxRRFaISHYzX22ag9IY85IxJtMYk3mmhYuDUf6+XKprtvPUo70YNyqWcaNi+f2jvah17iRv765OyXD++CkUsx+np5ai6gOsOvg+W0s2cqLyBK+8/jcOHS5s+SQBYkvWJhJIOVnUoX5RkTATSd6enNMeFx4WQe/0fmecq3zMiPHcd8+DRAx2UBJbQMb4dB5+8LFGs2RGhEfiNE0vhtZRS5TO06I6QYstdmPMuXlFzYf27s9n8kQHISHfvI+GhFiYPNHG3v35DOjT8UtuDRlwHpdPv5LFH75D+bFyBshwImzRJHRP4uixIv78wv/wxC+fadNImK4qJMTW7ApLXuPxyesb0GfwGf/P+vcZhD3aRmH5XlLCMhARyuvKKLUUc/7Y+9r9/Eq1RG9Q6gTRUbEUFHqbPF5QaIiO7LwW3BWXXcPFUy+nd9QA0hP7kJKcht1mp2dYLzwnDLm7d3Zalo40bsz5lMhB6jzOk48ddR7G43DRv0/HX0uwWqz85IcPUZtQweaaNWTVfkq+dTu3z76bpITklk+gVDu1d7jjtcBfgATgAxHJMsZM80myIDJyaCa//9s8lq85yiUXxiMC/15XypdbhYfuGdvu8xtjyNqxmVVrlnHiRCXDho3gsilXEhUZ3WRfV52bSFs0oY7Gq9PYCaXqW/PIdFUutwuLxdLs9MNQ32K+bNp0Pl66hG6SgFtc1NqquPeuBzvtekKPxJ48+tCTHCo+gLPOSVpKRtBdy1BdV3tHxSwCFvkoS9AKdYRyx42/4sXX/sRfX61fPDokJJHbb/yZT5b/+njlYpZ9+CGplr7EW1PI+ncWm7/cyH898HiTFY2GDDqPLRu+IN30O3nRts7rpNyU0jdjQLuzdKSDxQd4+9155O3OxWq1MD5zEtfNmNXsItVXXT6TCZmTyd29A4cjlGEDm1407WgiovPRK78I/A7VAJGanM4Dc5/mSOlhoH5kkC9Gw1TXVPHh0vcYGTqJUGv9m0SsPZ6c41tYt3EV37l4RqP9hw8Zxep+y9mat4Ee1lTcXjdFFHDZZVd0yJqlvlJRWc4f/vI7Ep1pTIq4HLfXTd7GHJ4/+kceuOeRZn+W3eMS6B43tfPDKuVn2sfeiUSExO49SOzew2dDHAuL9hNO5Mmi/rV4SxK5uTua7G+1hnDfnJ9zzY0zsfcTug2P4q459zBj2nU+ydNRPtu0hojaWNLC+2ARK3arg4ERIziwdz8HDu3zdzyluhRtsQe4mKhYarxVeI230dzr1d4TpMU1v1C1LcTGpHFTA2osdXFxEZE0vmYgIkRaYigpPUx6Sm8/JVOq69HCHuCSEpLp3acv+bt30Dd8CFaxcryulGI5wK0XzPZ3vDYzxrAp6zP+vXop5eXHGTRwCFdefi3p6RnkfZEHfFPAvcZDhfcYKUlp/gusVBekXTFBYM7sH5M4OJ7Pa1ayuWY1++y7uGP23QHZiv145WLeev01wopi6OMcyoFNB3nymV8zsO9gPNF17D6xk1pPNZWucrKrNjNk2DCSk5r/ZKLUuUpb7EEgMiKKe+/6GRWV5VTXVJHQPem0QwG7stramvoLwWHnE2qtH+nSxzaYvBNePtu0lofuf5T3PnyHbdmfY3c4uOCiqUz/1sVhpZQW9qASHRXT7EIagaKktBgHoSeL+tfibYnk53/F9TNmcfvNd/spnVKBQ7tiVJcRE92NWm8Nbq+70eOVnnISE8+tieOUag8t7KrLiImKZdTITHKrs6jzOjHGUOYsoVgKuGTKmZeeU0p9Q7tiVJcy6/t3scDxBhs+X4t4hajYaO6ceS8ZaX39HU2pgKGFXXUpDruDWdffyfVXz6K2tpqoyBgsFv1gqVRbaGFXXZLD7mi0LJ1SqvW0KaSUUkFGC7tSSgUZLexKKRVktLArpVSQ0cKulFJBRgu7UkoFGS3sSikVZLSwK6VUkNHCrpRSQUYLu1JKBRkt7EopFWS0sCulVJDRwq6UUkFGC7tSSgUZLexKKRVktLArpVSQ0cKulFJBRgu7UkoFGS3sSikVZLSwK6VUkGlXYReRp0Rkl4hsE5FFIhLrq2BKKaXOTntb7MuBYcaY4cBXwMPtj6SUUqo92lXYjTHLjDHuhs0NQGr7IymllGoPX/ax3wF85MPzKaWUOgshLe0gIiuAHs186xFjzHsN+zwCuIE3z3CeucBcgIT4pLMKq5RSqmUtFnZjzKVn+r6IzAauAi4xxpgznOcl4CWA/n0GnXY/pZRS7dNiYT8TEZkO/AKYYoyp9k0kpZRS7dHePvbngChguYhkicgLPsiklFKqHdrVYjfG9PNVEKWUUr6hd54qpVSQ0cKulFJBRgu7UkoFGS3sSikVZLSwK6VUkNHCrpRSQUYLu1JKBRk5wywAHfekIkeAgobN7sDRTg/hG5q98wVqbtDs/hKo2ZvL3csYk9DSgX4p7I0CiGw2xmT6NcRZ0uydL1Bzg2b3l0DN3p7c2hWjlFJBRgu7UkoFma5Q2F/yd4B20OydL1Bzg2b3l0DNfta5/d7HrpRSyre6QotdKaWUD3WJwi4ij4vItoY53ZeJSE9/Z2otEXlKRHY15F8kIrH+ztQaInK9iOwQEa+IBMSIARGZLiK5IpIvIv/p7zytJSKvikiJiGT7O0tbiEiaiKwSkZyG35X7/Z2ptUQkVEQ+F5GtDdl/6+9MbSUiVhHZIiJL2npslyjswFPGmOHGmJHAEuBRfwdqg+XAMGPMcOAr4GE/52mtbOB7wFp/B2kNEbECfwW+AwwBbhKRIf5N1Wr/AKb7O8RZcAMPGmMGAxOAewPoZ+4ELjbGjABGAtNFZIKfM7XV/UDO2RzYJQq7MabilM0IIGA6/o0xy4wx7obNDUCqP/O0ljEmxxiT6+8cbTAOyDfG7DHG1AFvA1f7OVOrGGPWAmX+ztFWxpgiY8yXDf+upL7IpPg3VeuYeicaNm0NXwFTV0QkFbgS+PvZHN8lCjuAiDwhIgeAWwisFvup7gA+8neIIJUCHDhlu5AAKTLBQEQygFHARv8mab2GrowsoARYbowJmOzAs8BDgPdsDu60wi4iK0Qku5mvqwGMMY8YY9KAN4H7OitXa7SUvWGfR6j/6Pqm/5I21prcAUSaeSxgWmCBTEQigXeB//jWp+suzRjjaejeTQXGicgwf2dqDRG5Cigxxnxxtudo15qnbWGMubSVu84HPgB+3YFx2qSl7CIyG7gKuMR0ofGjbfiZB4JCIO2U7VTgkJ+ynDNExEZ9UX/TGLPQ33nOhjHmuIispv46RyBcwJ4EzBCRK4BQIFpE3jDGzGrtCbpEV4yI9D9lcwawy19Z2kpEpgO/AGYYY6r9nSeIbQL6i0hvEbEDNwKL/ZwpqImIAK8AOcaYP/o7T1uISMLXI9REJAy4lACpK8aYh40xqcaYDOp/z1e2pahDFynswJMNXQTbgMupvxocKJ4DooDlDcM1X/B3oNYQkWtFpBCYCHwgIkv9nelMGi5Q3wcspf4i3j+NMTv8m6p1ROQtYD0wUEQKReROf2dqpUnArcDFDb/bWQ2tyECQDKxqqCmbqO9jb/OwwUCld54qpVSQ6SotdqWUUj6ihV0ppYKMFnallAoyWtiVUirIaGFXSqkgo4VdKaWCjBZ2pZQKMlrYlVIqyPw/PEyHyNShgzIAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x2099de17ac8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from ml_models.linear_model import LinearRegression\n",
    "classifier = DARTClassifier(base_estimator=[LinearRegression(),LinearRegression(),LinearRegression(),CARTRegressor(max_depth=2)])\n",
    "classifier.fit(data, target)\n",
    "utils.plot_decision_function(data, target, classifier)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 四.讨论\n",
    "\n",
    "（1）DART其实可以看做介于随机森林和GBDT之间的一种树，当`dropout=0`时，等价于GBDT，当`dropout=1`时，等价于randomforest；   \n",
    "\n",
    "（2）另外需要注意一下的是，当xgboost使用dart时，由于进入了随机性，会使得early stopping操作变得不稳定"
   ]
  },
  {
   "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.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}