Lecture_Huber.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Тема лекции: \n",
    "## Подход Хьюбера к определению оптимальной функции потерь при $ \\alpha $ -загрязненных помехах"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Выбросы и их причины\n",
    "\n",
    "* Ошибки измерений и регистрации данных: ошибки точности измерений, некорректно проведенные наблюдения, некорректная регистрация данных или ошибки при вводе данных в БД(опечатки, не там стоит запятая).\n",
    "\n",
    "<img src=\"./img/i_dpz_crop.jpg\" />"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Ошибочное предположение о распределении данных. Например, истинное распределение может быть сильно ассиметричным или симметричным, но иметь тяжелые хвосты.\n",
    "\n",
    "\n",
    "* Редкие наблюдения. В выборках, отобранных(как предполагается) из заданных распределений маловероятные наблюдения могут появиться в очень редких случаях. Экстремальные наблюдения в этом случае обычно принимаются за выбросы, но они ими не являются. \n",
    "\n",
    "\n",
    "* Загрязнение данных: данные принадлежат двум и более распределениям, т.е. имеется одно основное распределение и одно или несколько дополнительных распределений (примесей), загрязняющий данные. Экстремальные наблюдения, принадлежащие распределениям-примесям, могут появиться как выбросы основного распределения.\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. $ \\alpha$-загрязненная помеха"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from scipy.stats import rv_continuous\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "\n",
    "class gaussian_gen(rv_continuous):\n",
    "    \n",
    "    def _pdf(self, x, alpha, sigma_1, sigma_2):\n",
    "        \"Gaussian distribution\"\n",
    "        def gaussian(self, x, sigma):\n",
    "            return np.exp(-x**2 / 2. / sigma**2) / (np.sqrt(2.0 * np.pi)*sigma)\n",
    "        \n",
    "        return (1-alpha)*gaussian(self, x, sigma_1)+alpha*gaussian(self, x, sigma_2)\n",
    "\n",
    "gaussian = gaussian_gen(name='dirty_distrib')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Приближенно нормальное распределение"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAD8CAYAAABpcuN4AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJztnWd4G9eVsN87ANhFqpLq1ZJlVcuSa1xoO64pjp1mp9nrbBzvZrMt2Y0d50uyyWbTnLKbsomz9sYbO3Zc4u64SrTcJFm9d6pQpEQVSuwog/v9mBkApEAQAAEMSJz3efgQGMxgDmbu3HNPuecqrTWCIAiC0BeG2wIIgiAI+Y0oCkEQBCEhoigEQRCEhIiiEARBEBIiikIQBEFIiCgKQRAEISGiKARBEISEiKIQBEEQEiKKQhAEQUiI120BMsHo0aP11KlT0z6+o6OD8vLyzAmUIUSu1BC5UidfZRO5UiNdudasWXNMaz2m3x211oP+b/HixXogLFu2bEDHZwuRKzVErtTJV9lErtRIVy5gtU6ijxXXkyAIgpAQURSCIAhCQkRRCIIgCAlxTVEopSYppZYppbYppbYopf7B3j5SKfWqUmqX/X+EWzIKgiAI7loUIeArWuuzgAuALyml5gB3Aa9rrWcCr9vvBUEQBJdwTVForZu01mvt123ANmACcAPwoL3bg8BH3JFQEARBAFA6D1a4U0pNBZYD84ADWuvhMZ+1aK1Pcz8ppe4A7gCoqalZ/Oijj6Z9/vb2dioqKtI+PluIXKkhcqVOvsomcqVGunJdfvnla7TWS/rdMZkc2mz+ARXAGuAm+/3JXp+39PcdMo8it4hcWi/ddkQfaulMat98vV5a569sIldqDOl5FEopH/Ak8LDW+s/25iNKqXH25+OAZrfkE4S++NuH1/KHFfvdFkMQcoKbWU8KuB/YprX+acxHzwK32q9vBZ7JtWyC0B/+kEkgFHZbDEHICW7Wenof8Flgk1Jqvb3t68APgMeUUp8HDgAfd0k+QYhLOKwJazDD7sf3BCEXuKYotNZvAaqPj6/MpSyCkAqmnQASCotFIRQGMjNbEFLEsSRCplgUQmEgikIQUiTkKApxPQkFgigKQUgR07YkJEYhFAqiKAQhRZzYhFgUQqEgikIQUsSxJEwJZgsFgigKQUiRkASzhQJDFIUgpIgpwWyhwBBFIQgpIllPQqEhikIQUsSJTUiMQigURFEIQopIjEIoNERRCEKKhGQehVBgiKIQhBRxFERQFIVQIIiiEIQUCck8CqHAEEUhCCkiRQGFQkMUhSCkSCiS9SSKQigM3F4K9QGlVLNSanPMtm8rpQ4ppdbbf9e7KaMg9CZawkMUhVAYuG1R/B64Ns72n2mtz7b/XsyxTIKQkFAkmC0xCqEwcFVRaK2XAyfclEEQUiVSZlxiFEKB4LZF0Rd/p5TaaLumRrgtjCDEIiU8hEJDae1uY1dKTQWe11rPs9/XAMcADXwXGKe1vj3OcXcAdwDU1NQsfvTRR9OWob29nYqKirSPzxYiV2rkSq5Vh0P8er2fYUXwiyvK80audMhX2USu1EhXrssvv3yN1npJvztqrV39A6YCm1P9LPZv8eLFeiAsW7ZsQMdnC5ErNXIl19PrGvSUrz2vF3z75aT2z9frpXX6sv3781v0N57alFlhYsjXazbU5AJW6yT6aW/KKijLKKXGaa2b7Lc3ApsT7S8IuSY6j6Jwg9kbG07RHSrc319ouKoolFKPALXAaKVUA/AtoFYpdTaW62kf8EXXBBSEOEiMwvrthawoCw1XFYXW+pY4m+/PuSCCkAIyj8KypoKiKAqGfM16EoS8Jdai0C4ng7hF0NRSwqSAEEUhCClixoykC9WoCIXDMuGwgBBFIQgpEhubKFT3S8jUBEMFqiULEFEUgpAisbGJQo1TBMPhSHHEVNh86BQnOwNZkEjIJqIoBCFFYi2KQs18CpmaQBrpsbf8bgUPvL0v8wIJWUUUhSCkiFgUdjA7xd+utabdH6LTH8qSVEK2EEUhCCkSO38gHffLUCAUDqec9WSGNVoXrhU2mBFFkSQr9x7n7j9vclsMIQ/o4Xoq0BTRkKkJmOGU0oODpsw/GayIokiS5buO8siqAzIbVRDXE9Fsr1R+f8A+RiyKwYcoiiRxRkPBAh1B5iNLtx/hG0/n3sqTYHZ6ZUwc5SKDrcGHKIokcTI80sn0ELLDGzuO8ue1h3J+3p4WReG1B6115BoEUuj0Q+J6GrSIokgS54FI5cEQsos/FHal04kNYBeiRRFrVacSowmK62nQIooiSSIWhSiKvCEQChN2odaSmWYwOxAK8/0Xt9HaHcyGWDkjVlGmMjM9kEZcQ8gPRFEkifNABMX1lDf4zbAro9NY5ZDK+bc2tfLb5XtZsed4NsTKGbEWRSqKImpRyDM02BBFkSRiUeQf/mAYrcl5Bdd0YxROGxrsrpce80hScT2FJEYxWBFFkSQSzM4/3HJlpDuPImKVDvLBRrpFESU9dvDiqqJQSj2glGpWSm2O2TZSKfWqUmqX/X+EmzI6SDA7/wiETCD3HU+68yictuMf5IONWOWQSrp4SGIUgxa3LYrfA9f22nYX8LrWeibwuv3edRxLQmIU+YPT4eY6oJ1u1pPTdga9RdEjRpNKjMJZa1wUxWDDVUWhtV4OnOi1+QbgQfv1g8BHcipUHwTFosg7HOWd6xFqj6ynVGIUQyQhIt2sp3Rmcwv5gatrZvdBjda6CUBr3aSUqo63k1LqDuAOgJqaGurq6tI+YXt7e7/HHz/ZBcCadRswD+XmsiUjlxvki1wtrZ0AvLH8LSqKVM7kaj7ajUeBqWH9hk0Yh7cl3N+Ra2OjVTV1287d1AX3Z13OZEjnmh1siyqH99aso63ek9Rxa49Yv//EyZP9njNf2lhvClWufFQUSaG1vg+4D2DJkiW6trY27e+qq6ujv+NL1i2H1jZmz5lL7bxxaZ8r03K5Qb7I5V21DOjkwosuYlRFcc7kun/PSkpaW+gImJw1Zy618xO3B0eu5vcOwsaNTJ46jdraM7IuZzKkc802NZyCt98CYO78BVwyc0xSx7VvbIR16yivGEZt7cUZlysXFKpcbsco4nFEKTUOwP7f7LI8QGwwW8zmfCHiesp1jMLUFPusUXQqMQr/EMl6il0rO72sL3mGBhv5qCieBW61X98KPOOiLBEkPTb/cJR3rudvmWFNsdeIvE6WoRjMTi1GIfMoBitup8c+ArwLnKmUalBKfR74AXCVUmoXcJX93nUCQ2Q0OJTwB5302Nzek1A4HFEUvdvDK1sOs/7gybjHDZURdSjN9FiZmT14cTVGobW+pY+PrsypIEkgFkX+4a5F4Ym8juV7L25j3oQqfvWpc047LpI5N8jbUDDNrK+gS1lqwsDJR9dTXjJUZtUOFcJhHXVl5HwehabEZ0Rex9IVMCOWTm8CQ8b1lK5FkfoaFkJ+IIoiSZyHfLDPqh0qxM5nyfWaEIksiu6g2WcbCUQWvxrcbajnzGypHlsIiKJIgnBYR0ZBg/0hHyrEdsa5viWhsKa4D4uiOxTGH4wv0JBxPfVYjyKd6rGiKAYboiiSIHb0Otgf8qFCoIeiyP3MbCeYHdtRhsOaQCiMPxTf9TRkgtkxFlwq6eKywt3gRRRFEqRragvZI7YzdqPWk+N6ih0d+/txTw6VUvUDtigG+e8vRERRJEHs6FUsivwg9j7kvHqsGXU9xY6Ou+0gdt8xiqESzE5v4SaJUQxeRFEkQQ/X0yB/yIcKPYPZuc96KvKcHqPotq2cvrKegkMkmN3D9ZTCwEliFIMXURRJ4KzMBRAISSPPB2IDxrl2PZlhjdej8BqqR8ZVd7A/15OlQIKDvA0F0y0zLivcDVpEUSRBwDRjXg/u0eBQIfY+5Hp9g1BY4zUMPIbqce7+XE9OBzvY21DaS6HGWBS5Xr5WGBiiKJIg1ooY7GsJDBViXR5uWBQew7IoQnFjFP1lPQ3uNuT85iKvkZLSi53RLVbF4EIURRJIjCL/iO2Mcx+jCOM1FF6P0SuYHU1/jSfTUJmZ7chfVuRJzaJwMQFBGBiiKJJA0mPzD7fnUUQtipgYRYzyihfkDQyVeRSmRiko8hhprXAHYlEMNkRRJIHz0BtKSnjkC34XFYUVo1B4DNXj3LHZTvHcT0NmZnY4jM8w8HmMlJRej7iSKIpBhSiKJHAe7PJir1gUeUIPRZHDGEU4rNEaPIaB11A9OsrumEyseAMKJ+snGfflsh3N/PWD7+Vl0DdkWllfPo9KLetJLIpBiyiKJHAe7Ipi76AfDQ4VegSzc9jpOCNhr0fh8fS0KLpiLIruOHMpUglmr9h7nNe2NeelBRsyozGaVAZOoTTTagX3yVtFoZTap5TapJRar5Ra7aYsYlHkH27NzHYUgxWjMOJmPUF8i8LZlkzmXKff+q4Of2hA8maDYFjj86Tueso3i+LvH1nH/3t6c0rHaK25+8+bWHugJUtS5SeuLlyUBJdrrY+5LYTTwMuLvZzqDLgsjQA9O+Jcpsc6I2Gv0feEOyBuBdlUigI6CqIzYDJqQBJnnpAZjrqeUiozHlsjyn1FsedoO+VFqXWBXSF4ZNUBRlcUcc7kEVmSLP/IW4sin3BGrxXFnkGfsTJUcCvrKdai6GvCHfQTzDbD/cYe2m1F0RHIP4siZFoTDge7RWGGNW0pWmyddpyprTv/7ks2yWdFoYFXlFJrlFJ3uClIxKIo8rrmMw6ZYbYfbnXl3PlI7Gz5XLqeIjEKQ+H1qLi1nqCPYHYKxfQ6A47rKf7kPTexXE/KDuanFsz2GArIj6ynoBmm3R9M6ZjOoCV3a3dqxw128tn19D6tdaNSqhp4VSm1XWu93PnQVh53ANTU1FBXV5f2idrb2xMev3Wf1SjaTx6nyx8a0LnSlWtFU4jfbvDzs8tLGV7srn7v73rlgt31/sjrbdu2U9e+JydytXRbHeOe3bvobA/R7O+InHNPjEzvrV1PsMF6vBy5ugMhDAVhDUuXvUGxV/V5nsajXQC8+95a2uo9Wfo16d3LpsPdBLrDtNON3yTp49s6uvApjQmsWLmKhmF9t+Nc3Mu2jk46gjql8xxv6wIU9Q2HXX8GYsn29cpbRaG1brT/NyulngLOA5bHfH4fcB/AkiVLdG1tbdrnqqurI9HxW+t2w/YdzJg8gXVHGxLum0li5dq5fA96w3ZmL1zC7LGVOTl/MnK5xdJTm2HffgDOmDWL2vOn5ESuhpZOqFvGnNmz2dR+kBKfQW3tBQC8dHwjHDgIwJlnzaV2/jjAul6XXXYZoZdeZFiJl7buEBdcdDFVZb4+z/O9tW8A7Zwxey6188Zm7fekc80e2r+adtVJdVUJxzsC1NZenNRxnrdfo8Kr6W4PsGjxYuaOr8qoXKniW7EUf3c3l112GUr1rbRjWfun1wA/xRVV1NZemFX5UiHb1ysvXU9KqXKl1DDnNXA1kFp6QgZx8t/dzHpq7bJ91i5lwdQf62DJv7/GwROdrpy/N4GQlaIJ7sUovAljFD3bieNqqSi2xmb9zaWIBrPzzxceCofxeQw7PTa1GEWJz17wKQ9ifaFwmFBYp+ROlhhFflEDvKWU2gCsAl7QWr/kljAB08RjKEp8VjpkLvP2HRyfqFsNdHtTK8fa/exqbnPl/L3xh8KUFlmdTi4VRew8Cm+veRTdwXBEEfQOZjvB9zJb5n4VhROjCORfjCJ2wl1KMYpQmFLf6SsDuoWjrFKJN3Tau7bnYdpyNslL15PWei+w0G05HIKmtVCNz16sJmCGKTGy5zeOx6kuq4W6Fdw8aZ+/pSOYF+maAbvTaesOuZj1ZBAKx0yyC5lUlfpo94dOG6XGplhD4rkUWuuoRZGHHVLQjJbwSG0pVB1RlPmQ9eQoq/buENXDkjtGLAqhTwKhMD6PotgbVRS5ptXuqFPN0sgUJ+2hVEuezCPxh8KRTien8yjMmKyn3kUBg5aigNPnUQRiMucg8exsfygc6cTy0qKILNyUvOtJa20NsCIWhfsTVx0ll4p14GQ9tXUH87K8SrYQRZEEATNMkdcTsSjcWJOi1R7BuDWSOdllKQjHsnEbf8ik1O503ZmZHW/honBUUfRyPTkdanlx/66nzhjlkEmLoitgcumPlvHGzqMD+h5rwp1BkTd515NzjxxFkQ8WRTDGokiWTnvXoJlabGOwI4oiCQKhMMVegyIXLQq3XU+n8syisFxP1v3Iba2n6Mxs32kxCpOKEm/cKsPRGIVjUfQtc2zCQiYtikMnOzlwopPNh04N6HuCpsYXp4RJIhyFmk8xCufetaaiKIJRuQtpLoUoiiQImpbrKWpRuBDMzhvXU348HAEzHOl0c6m3e8coepQZD1mulWKvJ0GMwtPjfTxiXSGZzHI72mYp+RMdA1P2oXA4EsxP1rp2BleRBASXs560ji4ulZLrKebZL6Q4hSiKJAiEwhT1sChyP6p3Ri/trgWzbddTnigKfzBMiW1RmDn0d/eYmR1nKdQSr0Gxz+ixNgXEFJYs6j+YHZsSm8n02GPt1oTAloEqCtMqCljkMQgmee0dxViSJxZFrEXXnmLWk5OWLYpC6IEVzDYo8ij7fW4buT9kRgrOuZWWl2/BbCtuZMUJcrkexem1nnoGs4t9BsVeo0eBQIhZPjSJeRTOYMBjqIy6Go+2WYri+AAVRTAyj0IlPR/C+f2leRKjiD1/qhbFuOElgBXQLhREUSSB0ym5FaNwJtuBexPuHEVxMk8sCitu5MGjVE5dTz3mUZxmUYQp8Tqup74sCsf11H+MYnRFUXYsigEqe6soYDRGkUz2j+OudaxAt7OeYi2hVCyDzqBm4vCylI8b7IiiSALHoojMo8hxtkNs0CyVDI1M4rieTuaJReEPmRR5LIsil+mxjpvLyXpyRqZaa7pDph2jMOLEKKKz+633fbchR1FUDyvJSDB755E2wmEdsSgGGqMImtrOejIi7/sjkGcWRawllGwF2XBY0xWCCSNKrePEohBiCZp21pMTzM65RWE1yCKPkXJZ5EzQHbRcX2VFHjoCpuv+ZbAtCt/pKarZJnYehc8TzfqxSodbI+ZiXzxFEV0l0ZG/L6KKonjA6bEHT3Ryzc+X88KmpszFKMLhSPVY530yx0A0mO12G4qVOdnBV3sghAYmDHcUhVgUQgwBM2zNzPa6Y1E4qbFjq0pccT05558yqhyA9kB+KIoij2FXY3UvRuG8d2IS0aynnpaAP4USHo4VMWZY8YAtir3HOtAatja1cqw9EPn+eEu1JkvsehSQXBZg1PWUH7WeYs+fbIzCGbCNH16CUqml1Q52RFEkQTAUzfIAF2IUdoMcP7zElWC2E5eYOsryzXYE3Xcd+O1MNK/HcKfWk9FzPQYny6nYcT31EcxOxvXU7g/hNRTDywYeozjUYpUrrz/awdE2P7YRMKBYUzRd3PqyZDKfnGcmWsLD3RhFD9dTki4kJ1ZYVVpERZFXXE99YVdyLTicYHaxSxaFM5KZMLzMJUVhjUSnjrZu/+5TJnO++RLrD57MuSxg+YpDYU2x14OhlEszs/uwKOx20pfryekoE6bH+kOUF3sjKyoOpL0dOmlV+91ztJ3jHX6m2lbhQOIUkRIeKbhie2c9Ofes/lgHX318Q9wVAbNJKI1gthMrrCz1RsrFFwpJKQql1EVKqa3ANvv9QqXUr7MqWR7RO5idSoxi+c6jLNvePKDzn4oxeQOhcM4fKqcg4DS7k3nvsIk/FGZV/fEBf3cyGTO/Wrabbzy9KfLeGZ1a6bG5npntWBRGJOvJCWRD366n3jGKRAHgdr9JRbE3MqFwIFZFg21R7D7aTtDUzKypANLPfHImqsW6npJxI/WeR+Eo2Bc3NfHEmga2NuZ29UbnPhZ7jZRdT5UlPoaV+MSiiMPPgGuA4wBa6w3ApdkSKt84LT3WHuEt297MvmMdCY/99nNb+M7zWwd0/tbuIEVeg1HlRUDuy3g4k+wci2L7Cev8Ow63p/Q97+45zi33rejhH//U71byo5e293lMIBTmd2/u5ck1hyKdi+PWKfIaVnqsC1lPsSPqsI6uRVHi81ASJ5gdKeGRxDyKDn+IsiJPZBb3QOIUjuvJuUSzaqwyqelaFI6C6+F6SsGi6D3hbsdhq2z9riOptaWB4sgzoqwoeUXR7biefGJR9IXW+mCvTflX1jJLWIFT1cOiWLH3OLc/+B4/e21nn8c1t3Wz92gH9cc6BlQXprUrRFWpL+LfznVAuyXierJiFE4fuPNIamtT1O1o5t29x3l16xHAGlWu3n+CfcdPV7ZPrzvEI6sO8PbuY5zsDNIVNCP7vbL1MADjqkrw9Kq3lG1iYxSemKyfaDDbsCyK06rHWseV+fov4dERsFxPEYtiAPf70MkuptixJYgqinQtikitK4+B10g+PdbZp/caIk4bytQ6J0EzTHNrd7/7OecfXuajrTuUlGUba1FUDEBRaK15aXOTa4ugpUOyiuKgUuoiQCulipRSX8V2QxUCwV4WxYmOIF99fANak7DA2qr6E5HXqRZiazzZFelcWruCVJZYflFI7FPtDprc/1Y9j6/urdfT52RXEJ9HMaaiOBLQL/Ia7DzSllIn7XT0T607BFij3Xg++A5/iG88vZm7/7yJH7+8IxKA3dLYypHWbr77/FbOmzqSa+eOtSfcpa8olm4/kpJrsPcKd2C5XmItCis9Nr7rqdhn9AiCx6PDH6Ki2HuaRaF1cpPbYs95pLWbS2aOjmw7o9pyPSVrUSzb3sy6Ay0x3xmbHpzYogiHdeQ6nBajMDVBM8yeo5YlsTNDFsWvl+2h9t66flOAnd8xvMyHGdanzaQ/fKqbG3/9NpsaTuEPmXz+9+/x8Epr6d2KEm9Srqd1B1q48id1p13r1ftbuPOhtTy+uiHh8fXHOuLGp9r9IQ6e6MxpJedkFcWdwJeACUADcLb9Pmsopa5VSu1QSu1WSt2VzXO9vC/Iv/dyD21tbI007mgJD+ty/c+be2k82UXtmWPYe6yjTx/yqvoTkQD4pobkFcWpriBX/fQNnthpNbDW7iCVsRZFH+er29HMlT95g+8+v5VvP7slYyOWk51BqkqLUEpF1nm+bt5Y/KEwB1JYGnX/cWvfN3Ye5Xi7n3pbcQR6jUifWd9Iuz/EiDIfW5taueHsCfg8iq2NrfzklR34Q2F++LEFGIbCGEAJj62Nrdz5h7X823NbADje7rfWxE5AdB5FNLmh3R+KKgpvHxPu7PfOAliJZ2ablBV5elgUb+8+xoXfX8qvlu1O+vcdPtVNWMP8CVUMt+/b2MoSqkp9nOgI0NDS2W8iwHef38pn71/FLnvk75QsiY3ZxX6HGY4q/p+/tpMr7n0DM6yjitJroJTlwqs/1kHQ1JT6POxuTqwozLBm37EOjttzQfriL5ub6AyYPLP+UML9nN8x0nbntvUqtvndF7ay7sBJnt/YyNr9J3l9ezMHW7oYU2pZksNKvP26rH65dDd7jnawv5fF7Awa/7K5qc9jj7b5ueZny/nO81t6bDfDmmt/vpxLfrSMc//9NQ6d7EooQ6ZISlForY9prT+tta7RWldrrT+jtR54JLMPlFIe4FfAdcAc4Bal1JxsnEtrzV/qg/xhxf5IA995pI0P/OJNfrF0dyTDJtaiaPOH+NyFU/n0+VPQGrY1xQ/Erdx7gvOmjWTiiFI2xrEotNb8ceUB6nY09+jUn93QSEfA5O3GEP6QSWtXkKpSXyQQ2t4d6jGaOdUV5F8e38Bt//sepUUe7rxsBh0Bk7X7W047p0MgFGbXkbakAuOnugKRjmZEmQ8FfHzxJCDqY+4PrTX7jndw8RmjMcOa5zY0Um+PJmMzgLTWPLRiP7PHDuN3n1tCRbGXT58/mZnVw1i7v4UXNx3mwwvHM82Ol3gNdVowe1tTK9//yzbMsGbnkTa+/Mi604oZdgVM/u6RtQTMMPuOd3KyM8DXntzEh3/5Nic6Aphh3WMAoLXmsdUHI4rR41HMn1gFwJr9LXSHermeQmH+4dF1fP0pKwgfMMMoZVkiPo+KtLWugMmRXq6SjoBtUdiK4rmNjXzm/pUcbu1m2Y7oWhLNbd384C/b6QqYtHUH+c0be3q4OJ1A9oThZUwfXY7XUFSV+hhZXsT2pjau+Mkb/GlH4pF3d9Ck3R/i8w+u5lRXsEcJk3jJHd95bgvX/nw5/pDJn1Yf5NDJLtYfbInMo4gtfeK0nSvOqubQya4+O969R9s573uvUXtvHRf/cBkPr9wf17JqaOlku/2dj/UxWt9xuI2mU10xridbUcRY6W/uOsoLG5vwGIoV9SdYWX8cpWDF3VfyH5dYk+2GlXhp7Qr1mUix92g7r9uWalevGJMTuH9nz/E+LZ8XNjYSMMM8uupgD0Xz3r4TNLR08anzJxMww7y4sW9lk0mSzXr6kVKqUinlU0q9rpQ6ppT6TBblOg/YrbXeq7UOAI8CN2TjRFsaWznptxYh2dxodeZ/XHkAreHBd/axxja7K4q9eAyFoaCmspivXD2LeRMqI9/RHTQjjSYQCrP2QAs7jrRxwfRRLJhYFbEomlu7efCdfZzqCrKy/gRff2oTt/3ve1zzs+WRzuzx1QcpL/LQEYSl25o51RW0/KK2onh+YxML/+0VXt16hC2Np7jmZ8t5cm0Df1s7g+e/fDF/e/kMPIZi+a5op9LuD/H9F7fxlcc2oLXml0t3cdXPljPvWy9z3X++yTef2RyZueuw+dApPvGbd1lV38Jwe0GeUeXFjKtQnDNlOEr1VBThsOaRVQc4fCra8a0/eJIH3qqnuc1PdzDMNXNrOKO6gqU7jrLPtjCcwG530OTeV3awtamVz1wwhSVTR7Lp21ezZOpI5o6vZNW+E7T7Q3xk0YTI9/dOj+0KmHzp4bX89o29rKw/zgNv1fPchkb+9ckNPTqX5zY0svdoB1+8bDoAaw+08O6eY5zoCPC1Jzfy4V++xUf/+93I/st2NPOvT2zk9+/sAywFtWDicMqKPKzYezySQlxaZFkUZljzwsYmHnvvIKf81upuPo+BUooiryfym7/8yDquuLeuR2fQYafHOq6nJ9Y0ML6qlE8umcQqzCgBAAAgAElEQVTWxtZIJ/cfL2zjN2/s4U/vHeD+t+r5wV+28/0Xo4kBjnU0cUQp8yZUMWVUGYahGFHmY9W+EwRCYZYesNwYfREwNfMnVNHQ0smPXtoeUQo+wyoKCFFFcaoryJ9WH2TvsQ6+/exWjrRa7Wnp9uZopponWvpk55E2PIbi2rljAXpYFauaQtT+eBlrD7Twy6W76QyY/OCm+SyZOoJ7ntrM/W/Vnybr69usjvm2i6aytamVd/ccpzMQYv/xjogb7qZfv82PX9oRWbTIadfO7Gx/yOSbz2xh6qgyPn/xNDYfOsXS7c3MGVfJyPIifLa7cd74KgJmmB+8tJ3vPr+Vy++to7kt2u4feDsqX2dvRdHUytjKEsyw5n/ftu5b7DMD8PT6RqaOKsPrUfzopR2Re/78xkZKfR6+8YGzmDehkhc25UZRJLtm9tVa639VSt2I5Xr6OLAMeChLck0AYp3sDcD5sTsope4A7gCoqamhrq4urRM9szuAQqNR/On11Ryb7OWxVZ1MGmZwsC3Ip+57l6piRXXXfurqDvCB6T7mjlKsWfE2WmuG+eDFVdv59atbUQqunuLj+b1BjnfbwbvW/ZT7TQ6cCPLZX7zMu40hQhqWr9+B34QSD3x2ThEPbO7g9t++zjVTfWxs6OaW2UW8uDfEN55cx4luzdzKIJvWWUrr+Q0NhDV8+eHVeA0o8ii+cUEJ00sOs+JtK9A7o0rx4pp6zi0+THNnmB+s6uaELdOCkmP8ab2fKZUG80Z5ONDWwR9XtPLU6v18aVEJs0d68Juab73TxSm/Rmuo0iZ1dXVcO9akfViYVe+8xZhSxVub97DQa5n5bx0K8j+bAkyvMvj6+SUcbAvzw1XddJvQuN9ymZxs2M2UEpO39rRz9JhllJ5oOUVdXR2/Xt/NqsMmF473UN2xl7q66MNW1GEp0eHFCv/BTdQ1WA9sV2cXzWYnDz77Oi/u7uKnq19l7zETrwG/fGENG46GqCxSvLzlCF/7/atcP80aQf5hdTdjShULvdb1+vkL6+gImEyrNCLB9mFFUFdXR1hrvvl2FyUe6Laf+bfeXI7XUEyvhNc2HWCFT1FTpti1fiUNB6xOx1FgS+s76NIHMHTY+r5QgIMNjTz83FFe22aN+v/qvuV8/fwSvIairSvI8SONbFhjdXxBU7NgRIhK/xG6giaPvriMgKl5en03HgW/em0b3abGZ8Ajqw4wXTUzc4SHt3dbymvXxlW8rxyWzLN/T5fVKU2rMmhoM7nroTf5woLiuM9Hl99Ptcdk/BQvD688QGmndb127dzOqTJrnLli9Xr8Bzy8tj9EdzBMZZElR7EHxpcbPLe6nksmWl3NqhXvQjjMvv0HONqlqS6FjgZLuT33xnucnOhj2YEg/7fVj0bxhQfepcWvef9kL2M79/JX0zVtJz386C/bqGzbR3VZdKz7+HvdjC1TLClp5mEDbvndishnk4cZjCpVdARM6g8dZp06BsCxxgMAvLNqDS17PDy7J0D9sSBfWVyM6mzEDGs2Npzimile6urqaG9vp66ujgqtuXySl/uW7wVAAV998A1un1fM3lMmf1zZzZkjDHa0hFm9YRPe5m2RNrGjqZOrpvpYHVT811Lrudi+Zz+3zSvmeFeY1oBm/cFuPjHLh3+4h2c2NbGu/jAfn1XEs5v9zB/lYdU7bzG7LMATu4I8+ZelFJudafeByZCsovDZ/68HHtFan1BKZUkkwLruvelh42mt7wPuA1iyZImura1N60Q/3fwWM4abBIwSTnqH0TZ8LJ2hDdx/83n8/LWdrKw/wX99+lwumzUGgN6nOXvvSt7cZTW60RXFPLjVzxnVFfy/G85g4aThTBtdzjt7jvH4zpWsPBLmk+dNpqUzwGtbm/EYihsXT+KemxYwZvke/uPF7aw5YlJe5OFfP34ZrQ/V8UJ9kBvOHs9/3DjfKlVR9wp+ExZNHs6B4514PYo/3XFhJHXVYaO5i5+9tpOyKQv47yc3Eja8PHDbQu74vzUsP1HJkc5mvvuRuXz2gikA7DrSxm3/+x6vN5dw500X8a1nNnO4Yz9//ML5XDh9FLH3u66ujtraWs49tJblO48ydvY5TB5ZxtfuraN6WDF7T/m5b1cJWw6doqK0mO52P6tPlgEBPnzFhcxsauX1h9ayo8X2W5eVU1t7KT9Yv5wrZ5dy/23nnnafSvce5+HtK/jYuVO54vKoF3L45rcYUV7EsZIqVjbvBkxuf980TnUFeXKt5X64/9bFPLGmgcc2H2bcxKncdtFUtr3yGn99yXSuv2o2MzbWsfGoNaL/45dqeXJNA+sPnmTF3hPU1tby1LoGGto38KtPncNDK/bz3r4TXFFbi2EotrGHH760HdB87drZXF47gwPv7IMdW/B5FHPGVfJucytXzR9P2bHD1NbWMuy9ZYwaM5x13Qalvia+8cGzuOepzXxzZZgZ1RWEdAdzZ03n/RdOgbpXAPibD5xHsdfD7zYtp2zCLF5Y3cCo8jB3XTebf3liIwAP3LaEbzy1me+v6mbehCqOtplUDzO46orLe1zLF49tYP3RBu7+8CKefmsjL9YHWXTmVL5y9Sx6P9d66UtMmzyJf7pqFlf99A2e2289hvPnzWVWTQWsfJP/Wudn0shSgiGDhZOGc8u5k7jrz5u4bv54zhxbyQ9f2k7JqAlAPbWXXULJiqWMmzCBHTuaOXtaJR+7dhHffPdlThVXc7JqNA9uXc/CMV6++uHFfO6BVXgNxXc+dRljq6zy3med08VVP13Od1cFMcOaX35qERfNGM3OV1/hcxdO5YNXz2Hq3FNsbDjFya4ARR6Dn7yykwNtlpavGjGSs+ZOhrVrmHPmTNi5lQULF3Lm2GF88bWlXD9/LF/++GI6AyF+vvYVQmHNxy47m9o5NZG2D3DJpWF+/MoOFkwYzroDLdz/dj3XnnsWD63eS01lCT//3Llc959vMu2MWdSeOxmA7YdbCb3yJteeP5frlWLN/haOtvl5bdsRbq6dzb/+cW3EevjHmy5hbGUJV21u4mev7uQ/11pt9Pb3L6R23jimzOvgiXvraKmYyhnmAdLtA5MhWUXxnFJqO9AF/K1SagzQfw5a+jQAk2LeTwQaM32S5tZuNjac4qMzfeiKEdTtPMrmQ63MrK7ggukj+fnNZ7P9cFtEScRj3oQq3tx1jA8vHM93PzKP5TuPcvXcGoq9nsg+F04fxa8/fQ5LpoygurKEgyc6eWXLEQLBMB+zff1/fbHlAin1ebhsVjUjyou44Qwff3XNEhZPGYFSPX3xHzl7AtfNG4vPYzDCDsjFUnvmGH766k4+8dt3KfYa/PEL57N4ykgumTk64ju9ek5NZP+ZNcO45bxJ3PvKTtbsb+HhlQf4zAWTuWjG6NO+2+Gu62azev8JPvM/q6ipLOZIq58n/+ZCHl11kGfWN/LBBeP456tn8YnfvMv6gyfxGooJw0upLPGhlDX/AKKup0AoHEmf7M2iySP4wiXTuP3iaT22G7YbI2CG8Rmw6d+updhr8NbuYzy5toERZT4unTWGy2aN4e4/b+I/X9/FK1uPEAprPjB/HAALJw1nz9EOZo8dxriqUv7uipn88KXtkbWl1x84ybBiL9fPH8v500eypbEVw3ZBXDhjFGDFHj662HKJOUHuRZNG8KnzJ/OPf1rPu3uOR7KEfB6Dw63drNnfwqfOm8ynz59C9bAS/rBiP4daOvnCJdP45LmTIsHsySPLmD+hCjOsKfEZPPZeA+/uPc7Xr5/NTedM5FfLdlNV6uPyM6t55I4LeGJNA6v3teCtUlw2q/q0a7l4ygj2HO3gitnVGEd8VI0Zyy+X7aa6spjPXTiV5zY0MmNMBXPGV0YW7Sov9vKRRRP4dd0e6zcYijNrhvHAbUvY3dzOa9uaWVV/gruvn801c8eyqv4Et188DZ/H4IcvbedpO9PNKSYYNMMcbw9QU1mC12Nw2awxPLGmgSfWNHDetJH89RndXDJzDD/86ALMsI4oCYBxVaXc+/GFPL+xkVe3HuHNXceoHlZCwAxz9uThgPVczptQFTnmvGkjeWZ9Iyv2HicQCkeSEopj5nU0tHThD4W5adFEwFq2dv7EKtYfPMl5U0eedh29HoO7rzsLgItnjub5jU18/alNeAzFH24/j7GVlsyx856c+MTc8ZWcUT2MDy0cz47DbbywqYm/eXgNU0eVc/O5k6gq9THeLj74wQXjuXrOWP67bg9rDrRQe6Z1T6eNLuescZWsPdDCGRPIKkkpCq31XUqpHwKtWmtTKdVBlmIGNu8BM5VS04BDwM3ApzJ9kqPtfs6eNJyzq/14qkfy53WH8BiKJ+68EKUU46pKGVdVmvA7rpxdzfKdR7nnA2dRVerjQwvHn7aPUorr7U4JYNLIMm6/eBpr9rdwjt2wDUNxx6UzehxX5FEsiWmghqGs2EXA5OKZo6muLKEvFkwczv23LqEzYDJnfCUzxlhpkR9aOJ5lO45yzuTh1PQ6/oMLxnPvKzv524fXENaaL/aSpzeTRpbxf7efz1cf30Cx1+CbH5zD4ikjWTRpBN+5YV6k07945mgeW93AxBGleG3FNntsJduaWvF5Yuol2Z1SPIq8Bvd84PR8Bic9NhAK4zWiE7oumjGaKaPKeP9ZNZGg648+toBRFcX85o09TB5ZFokxLZw4nD+vPRTp9MHKXgqaOpI6WVbsQSnF6IriHgOHeeMrqSr1cf60kVQPs65nsb3mwkVnjOL9c2rwKqs430S7PLXPY7BmfwtBU3ODHW+5ak4NV8Uo7ug1LuWT505CKWvp0bPGVfLu3uOU+jx8cslkPIbisS9eiMdQKKWYMqqcr1x9ZsL79slzJ/NJe4TrNRT/ceN8/rL5cGTS2zee3sz188fxvY/MIxTWket3xezqiKLw2vGWK2bXcMXsGu64dAYtHYHIoOWnnzwbsJIAvnDJNB5eeYDKEi8+I1rx1wnaA/zmM4tZur2ZFXuP8/dXzGTdqrcB+MSS2PFilGvnjeXaeWP5yK/eZsfhtkisbPbYYXH3XzBxOAsmDueW+1ZYiiIczcICKwvLmUzpi2mDn7twCvMnVEWy/fqiqtTHq/98KYdOdjGirIiaypJIJlxX0OTQyS5ufWBV5JxOKRWAM8cO45KZo1lZf4JffmoRc8dXnfb9RV6Df3j/zNO2P/T58xhZXsQbb7yRUL6BklBRKKWu0FovVUrdFLMtdpc/Z0MorXVIKfV3wMuAB3hAa72ln8NSZu74Kp7+0vuoq6tj8vSRKAVfqp3Boskjkv6OJVNH8sLfX5Lyub9+/VkpHwNWUbmqUh/TR/dfduvKs07veK6aU8Oo8iI+unjiaZ9NHV3OgolVbGw4xQfmj2PSyLLT9unNmWOH8dyXL+6xzTBUD8vg4pljeGx1A5NjHo6LZoxiW1Mr00dXRNa6CNjl3FMhYlHYisLBYyhe/sdLI50cWG33rutmM3d8JaMqiiJt+Xz73l8xOzr6djp7p2SKo4B64/UYPH7nhYyuiPr4q+wA6SUzx1BR7GX2SA+bj5vROSgeRdC0rIN5cTqFWJZ+pRZPzDM3f0IV6w6c5KZzJkQ6r0QDhmRQSlFZ4qO1O0g4rGnrDuIPmZFif47yXjR5BMPLfJzsDEYC2bHEs2yVUtzzgTn801WzaPeHMOwFj9r9IcI6WiTRY6g+lWUiZo8dxitbj7DjSBtFHiNS4bgvirwGHYFQxKKIlhTpWULe4cZFE7lx0enPSjyGlfiYPTaqUIrtFRg7AyG2N7VGgvXnTR0ZmdXv8MtbzuFoezdnVMdXdH0xqiJ+bCnT9GdRXAYsBT4U5zNNlhQFgNb6ReDFbH1/b6aPqaDuq7VMTqJzdJOFk4Yzq6biNF9ysgwr8bHqnvdj9HH4hxeOZ2PDKT5/ybT4O6TB++yR+tSYGcIfWzyRfcc6GFFexFLbFeaUDk8Fr2GlmgbNcCQjxaGvzr231Td7bCUr776yR4frKCxnLY5ECsyZ7exw6cwxPH7nhSyeYg04zq62FYX9HY7yWjRpRJ8WlIOv1/U4b9pIHll1gNsumprwuFRxSlJ0BKwO3Lqm0XIdYHXml80awzPrG/EZqd2nsqLoTHOPoSKTxSqK49+jZDlz7DAefe8gb+46xozqitOuV2+KvEYCiyI6mTITKKUo83noDJiR1N9HvnABZ8axeqrKfP1aLW6SUFForb9l//+r3IjjLv2NRvKB331uyYC/I9GDcOtFUzlv2kgWTBw+4PM4jKoo5he3LGLBxOjo+axxldx/27l857mtkTkFgQSup75w1swOhMLEGeQmTe9RuRNj8vdjUcTD6zE4N8ZluHCMh4e2RTt95/+5U5O3XB0+MH8cF0wf1cOCyQSWoghG5hMEQuHIfYlV3lfMruaZ9Y2UFqW/QoE3RlE4FkW6zB5ruQ+3NbVy46L+HfWOogj2sihC9nwpR75MUVrkoStGUcwYUx6Z5DeYSOouKaWqgG8RLQT4BvAdrXVqdSmEvMfnMTKqJBzixW4AfF4VDWabqSsKww7yB+300EzhjDT9of4tiv4YU2Ywe+ywiDvO8YEviRMg7Q8nTpJpKkt8HDjRGZmwFzDD0TkTMb/9gwvG4/MYLJqUupJz8GRUUURH5/FG6r0p9hgEzKj14KzhbYZ1ZIZ/piwKsMrKdwbMyDyNipKB/V63SFbqB4DNwCfs958F/he4qc8jBCEJij0GQfvBNcOaIk9qrgivbVH4Q+GMPuBOjMIfCtMdMiNB13T5r1sWRVbiK/JYEzcXTc68Qk4Xq3ZRKK5FEevO8Rg9EzPSoafraWDXdUR5EdXDimlu8yelKKIWheN6itaeMs3MK4rSIm/E9WSoaK2rwUayd2mG1vqjMe//TSm1PhsCCYWFz2OgdXTNhZQtCjuDxkmPzRQR11MwjD8YZlT5wB7w2DjGjOoKgqZmWEn++KStkhTBSIXU2M401bhRf3g9KlKCpqyPdOhUmD2ukua2o5xZ07+i8NkWRex6FGBZFKEMxyjA+n1dQUsBlxd7044tuk2yiqJLKXWx1votAKXU+7DmVAjCgHAUg+PDTTlGoRRhrQn2ynoaKCW+GNdTyIxYGJng7uvOSqkKbC6oLPHSHghFRvoBM9xjgahM4jGMyByagVoUABdMH0n9sXbGVfWf/eVYFFHXUzRGYUZiFJn7vWVFHtr9Idr9IYZl4Le6RbKS/w3woB2rAGgBbs2OSEIh4bg1HB9uqrEAZz2KgJlZRdEjmB0MU+LNrMsg30aWlaU+tIYmu+ZQIBSOFPLrL5MoVWKDxQONUQDceekM/vri6Uld09NdT9GsJycTKqOuJ5+Ho21+2rtDgzY+Ackrim3Aj4AZwHDgFPARYGOW5BIKBCdQ2jYAi8IpY53JbJXYYLY/wxZFPuKsdeKUrQ6EohaFbyDpZHHItKIwDEVRkve+yGMQstuLoaJK0OxhUWQ+mB07uXAwkqzkzwAngbVYM6UFISMU2w+qs2pfyhZFTHpsWQb7M0cxdAet1esybVHkG068xFk6NTbrKRsxCofyDMQoUsEZiHQGTLweA4/HWaUwOzEKJ5jd1h2KKOPBSLKST9RaX5tVSYSCxOe1HkrH9ZRqp2Slx1oj4MoMPodR11MBWxSh7MUowIoD9Z6hnG2cgUhXwMQbs0phD4sigxZUWZGHroAVoxg/fGAz6N0k2bv0jlJqflYlEQoSJx023WC2tQhO2I5RZN711BkwrXIbQ9yiqLQtisaTp1sU2YpRuOGKcX5Lh60oDJVdi6KsyENn0FpYqhBcTxcDtyml6gE/VhlwrbVekDXJhILA8X93pKkorFpPoHWmg9nWl7V2WXKVFIhF4SyyExvwzbSicDriTMQnUqUoYlGE8HmMnhaF/XszmfVUWuRBa2uN8ori/EmHTpVk79R1WZVCKFh8vdNjU+yUPAaEteU2yGRf7qRNOumiA5mZPRjoPafDCmbbEwS92QlmO7WfcklRJCZm4jFURGmFwhpnGfOMWhR2OwqaeuhnPWmt92dbEKEwcYLZ7XbN/tRdT9ayoyEzjDeDKacRi8KeGJZKrafBSO9AayisI2WyU50t3x+eiOsp99c0EswOmpGlaa2lWaNlxjOa9RRjNbnxezPF0B4mCXlP1KKwOuR0aj2Zdq2nTD7gXntt54hFMcRdTyU+z2nX3nEH+rJkUbjtevLGVMU1w2QtRuEwmF1PQ7v1C3lPrCsA0kmPJSaYnVnZir1GpKTFUA9mgzU7G6IdeURRZDxGYX2fm4rCcT2BXS8sHM5KraceimIQu57yTlEopb6tlDqklFpv/13vtkxC9ojMzI7EKFLrkA1D4Q85LoPMyharKIa6RQHROIVTBjtdd2B/RLKeXIxRdAZCkTU1PIbqmfWUQRdmqS/6GwuhhEeu+ZnW+l63hRCyT1HveRRppMc6ZZMy6XoCay7FqQKyKJw4xegKqxprR5oJBv3hTHJz06KwJtzFup6shAhDEVkPPROIRSEIGcCxIDoCac7Mjhn9Zdyi8BkxMYqhryicuRSjKhyLIjuuJ5+bwWxPtHy8N8b15FgUmUyNhZ6KotwFCypT5Kui+Dul1Eal1ANKqfRXSBHyntNmZqcxjyLyXRluzSVeD61pFiscjDgWxRh7YaQOew2FTPrsIT9iFEBkVrjHUJimxgxndk0ToMfa8YVQwiOjKKVeA8bG+ege4L+B72Ktyf1d4CfA7XG+4w7gDoCamhrq6urSlqe9vX1Ax2eLQpCrNWD5jU60dQLw7ttvpuRCOrA/EHkdCvgzer0CXV2Rsg4b163h2K70lEW+3kfoKVt7ix+AzpYjADQcOYpXkXHZGw9Z9+zQ/r3U1R3sV65M0tQejp6j9RR1dXWEAgEaGpso8wHaTHjeVOVqC0TLyW9YvZI9RdmpGpztNuaKotBavz+Z/ZRSvwOe7+M77gPuA1iyZImura1NW566ujoGcny2KAS5WruDsPQVAtpAKZMrL69NqQT3VnbDrh0AVJSWZPR6/Wr7O9S3tgBwyUUXMGlkWVrfk6/3EXrK9lb7Vt48VM+is87gpX3b8ZVVUtzRnnHZV3Zvh317WDT/LGoXTexXrkxy8EQnvLUMgNGjRlBbewHlq5ZRXTOCyhIvxUcbE543Vbm6gyYsfQmAq6+4LOOJAenKlSp5Z08rpWLXWbwRawlWYYgSzUIxKbInQKVCVmMUMQHsoT7hDqJZT6Ns11N7dzDjgWyImUfhRtZTrOvJdoH1jFFkOiHCQCnrf7aURC7IR6fZj5RSZ2O5nvYBX3RXHCGbxAZK03mQYn3K2UiPjbwugPTYaWPKqSzxRlaK6/CbGQ9kQ+zMbPfSYyFaZyw6M1tnPEahlKLM5xn0A428UxRa68+6LYOQO5x6O2ZYpxUw7qkoMjwajFEOhZAe+6EF47jqrBp2N7cDVjB7hD2nIpPkw8xsiLYdj73uejaynsBak6J8EJfvgDx0PQmFhzPKS8fNkV2Lwnq4rZXQ8mvp0myglKK0KFrKoz0QysrvzsusJ3seRaYtCrBSZAdziXEQRSHkAU5nlI7rychqjMKw/3vybo3rbOLcB60zP4cC3F2PwmsonFsZfx6FKIp4DG7phSGB0zGloyhiH+wM166L+JWH+loUvYm1IrIRgJ1ZU8GMMeWMKM99kTylFD6PQSAUjriZohZFOKOzsh2uPKv6tDLugw1RFILrRFxP6VgUPSbcZT5jxfo/uP3LqRJ7H7KR9VR7ZjW1Z1Zn/HuTpdhWFI5CjJaqz45F8S/XzM74d+aawhoqCXmJU2o8rRhFDlxPhWZRFMcUZsyG68ltHEUYG8zOZoxiKDD0WoEw6HAURDojd68ni4rCdj0VskXhG8S5/33h/D5HCXo91rrr2YpRDAWGXisQBh2+gbieVBbTYwvUoujpehp6Hafz+5z24ix+FdZiUfRFYT0BQl4ykGB2LibcFULl2Fhi15IezLOJ+8KxYD2eXllPZnbmUQwF5KoIrjOQYHYu5lEUQuXY3jj3ZCjGKJzf5Dst60ksir4Yeq1AGHQ4pcaLBxjMzrSHyJmZPdjLL6RDbz/+UCLienIsCo8zjyLcI+YlRBl6rUAYdGTCojBUz3hFJihoi6IQFEUk68kgLBZFQoZeKxAGHQMKZmfRl17QFkUkE23odRHFEYsi+9VjhwpDrxUIg46iAcyjcB7sbIx8oxPuCu8xKY5YFEOv4yzy9LQonKwnsSj6pvCeACHvGNDMbNvdlI3O3HE9FaRFUUCuJ8uiCGeteuxQQK6K4DoDcT1F0jizaFEU2jwKGNqKIjrRLpomKxZFYoZeKxAGHZmYR5GNGcSOgii0mdkwMOWd7xT1cqtFYxRhURR94EorUEp9XCm1RSkVVkot6fXZ3Uqp3UqpHUqpa9yQT8gtvgGU8MiuRVGY1WNhYGuE5DvRWk8x8yhMTTiMKIo+cKsVbAZuApbHblRKzQFuBuYC1wK/VkoV3nCuwBiQRaGyF8x2FtapKB7cJaLTofeoeygRnUx4ukUhWU/xcaXMuNZ6GxBvMZgbgEe11n6gXim1GzgPeDe3Egq5xKknlM6EOyf2aHVsZgalgpHlRfz+r87l3KkjM/q9g4GIohiCrqdIeqwRtSxMLTGKROTbehQTgBUx7xvsbaehlLoDuAOgpqaGurq6tE/a3t4+oOOzRaHI1XAwAMDuXTuo69iT2rFtYQA621tpbzezcr3eaxrY8fl6H6Fv2U6d6AZg766d1HXV51iq7F6zxgarvW3ftoXyEztoOBggZIbp8gc43NRIXd1xV+QaCNmWK2uKQin1GjA2zkf3aK2f6euwONt0vB211vcB9wEsWbJE19bWpiMmAHV1dQzk+GxRKHLtUHtg13bOnj+X2vnjUjp2d3MbvL2c6lEjqajoKojrlUn6ku3ZI+tZdfgQC+bNofbsuGM1V+TKBBvNXbB3J2cvmCqq7JgAAAs9SURBVE/tWTWsDe4kvHcXyvAyZdIkamvnuCLXQMi2XFlTFFrr96dxWAMwKeb9RKAxMxIJ+YpvAIFTZx7FUMzOcZOhnB5bFGdmNmAtjzoEYzKZIN9awbPAzUqpYqXUNGAmsMplmYQsM7A1s4du0NVNhrSi6KUgnLhEwJT02L5wKz32RqVUA3Ah8IJS6mUArfUW4DFgK/AS8CWtdWYjlELeMbA1s+3vKMC5DtlkIPck34k3M9tBsp7i41bW01PAU3189j3ge7mVSHATp8z4gCbciUWRUYZ0emwv11OsFSEWRXyG3nBBGHQMpJy382AXYuG+bDKQQo35Tm/Xk1gU/TP0WoEw6Lhk5mi+ctUsZo+tTPlYZ8LdUOzQ3GQoxyhKi3oWe+xpUQy935sJ8m0ehVCADCvx8eUrZ6Z1rCeLZcYLmaEco6g9cwz3fnwhs2oqgJ7KQZpRfERRCIOabC5cVMgMZYui2OvhY4snRt57xaLoF7kqwqAmmwsXFTIlBbQMrEdiFP0iFoUwqPF5DMqLPIyuKAK/29IMHd4/p4a7O2czcUSp26JkndhJdpL1FJ+hP1wQhjQ+j8HL/3Qpnzh3Uv87C0kzsryIL142I17hziGHWBT9IxaFMOiZOKLMbRGEQYxX5lH0i1gUgiAUNEaM1SS1nuIjikIQhIKmZ4xCusR4yFURBKGgiVUOEqOIjygKQRAKGolR9I8oCkEQChrJeuofURSCIBQ0YlH0jygKQRAKGikz3j9uLVz0caXUFqVUWCm1JGb7VKVUl1Jqvf33GzfkEwShcBBF0T9uTbjbDNwE/DbOZ3u01mfnWB5BEAqUnjEKcbLEw60V7rYBBVEeQBCE/Mbbo8y49EnxyEf1OU0ptU4p9YZS6hK3hREEYWgjWU/9o7TW2flipV4Dxsb56B6t9TP2PnXAV7XWq+33xUCF1vq4Umox8DQwV2vdGuf77wDuAKipqVn86KOPpi1re3s7FRUVaR+fLUSu1BC5UidfZculXEc6wnztzS4A/u2iEqZUevJCrlRIV67LL798jdZ6Sb87aq1d+wPqgCXpfu78LV68WA+EZcuWDej4bCFypYbIlTr5Klsu5TpwvENP+drzesrXntfbmk4l3HeoXS9gtU6ir84r15NSaoxSymO/ng7MBPa6K5UgCEMZcT31j1vpsTcqpRqAC4EXlFIv2x9dCmxUSm0AngDu1FqfcENGQRAKA1kKtX/cynp6CngqzvYngSdzL5EgCIWKWBT9I+pTEISCRtJj+0cUhSAIBY3HIxZFf4iiEAShoIlVDoYoiriIohAEoaCRGEX/iKIQBKGg8SgpCtgfoigEQShoDEPh6AopChgfuSqCIBQ8jstJLIr4iKIQBKHgcRSExCjiI4pCEISCx2sYKCVZT30hikIQhILHYyixJhIgikIQhILHayiJTyRAFIUgCAWPYSjJeEqAXBlBEAoesSgSI4pCEISCR2IUiRFFIQhCweM1lGQ8JUAUhSAIBY9YFIlxa4W7HyultiulNiqlnlJKDY/57G6l1G6l1A6l1DVuyCcIQmHhNQyJUSTALYviVWCe1noBsBO4G0ApNQe4GZgLXAv82llDWxAEIVuIRZEYVxSF1voVrXXIfrsCmGi/vgF4VGvt11rXA7uB89yQURCEwsEjWU8JyYcYxe3AX+zXE4CDMZ812NsEQRCyhkfmUSREaa2z88VKvQaMjfPRPVrrZ+x97gGWADdprbVS6lfAu1rrh+zP7wde1Fo/Gef77wDuAKipqVn86KOPpi1re3s7FRUVaR+fLUSu1BC5UidfZcu1XN9b0YXfhO+8rzThfkPtel1++eVrtNZL+t1Ra+3KH3Ar8C5QFrPtbuDumPcvAxf2912LFy/WA2HZsmUDOj5biFypIXKlTr7Klmu5PvGbd/SHfvFmv/sNtesFrNZJ9NduZT1dC3wN+LDWujPmo2eBm5VSxUqpacBMYJUbMgqCUDh4PQpDSYyiL7wunfeXQDHwqrJuzgqt9Z1a6y1KqceArUAI+JLW2nRJRkEQCoTPXzyNQCjsthh5iyuKQmt9RoLPvgd8L4fiCIJQ4Fwxu8ZtEfIaCfMLgiAICRFFIQiCICREFIUgCIKQEFEUgiAIQkJEUQiCIAgJEUUhCIIgJEQUhSAIgpAQURSCIAhCQrJWFDCXKKWOAvsH8BWjgWMZEieTiFypIXKlTr7KJnKlRrpyTdFaj+lvpyGhKAaKUmq1TqaCYo4RuVJD5EqdfJVN5EqNbMslridBEAQhIaIoBEEQhISIorC4z20B+kDkSg2RK3XyVTaRKzWyKpfEKARBEISEiEUhCIIgJKSgFYVS6lql1A6l1G6l1F0uyjFJKbVMKbVNKbVFKfUP9vZvK6UOKaXW23/XuyTfPqXUJluG1fa2kUqpV5VSu+z/I3Is05kx12W9UqpVKfWPblwzpdQDSqlmpdTmmG1xr4+y+C+7zW1USp2TY7l+rJTabp/7KaXUcHv7VKVUV8x1+0225EogW5/3Til1t33NdiilrsmxXH+KkWmfUmq9vT1n1yxBH5GbdpbMeqlD8Q/wAHuA6UARsAGY45Is44Bz7NfDgJ3AHODbwFfz4FrtA0b32vYj4C779V3AD12+l4eBKW5cM+BS4Bxgc3/XB7ge+AuggAuAlTmW62rAa7/+YYxcU2P3c+maxb139rOwAWtVzGn2c+vJlVy9Pv8J8M1cX7MEfURO2lkhWxTnAbu11nu11gHgUeAGNwTRWjdprdfar9uAbcAEN2RJgRuAB+3XDwIfcVGWK4E9WuuBTLpMG631cuBEr819XZ8bgP/TFiuA4UqpcbmSS2v9itY6ZL9dAUzMxrn7o49r1hc3AI9qrf1a63pgN9bzm1O5lLVu8yeAR7Jx7kQk6CNy0s4KWVFMAA7GvG8gDzpnpdRUYBGw0t70d7bp+ECu3TsxaOAVpdQapdQd9rYarXUTWI0YqHZJNoCb6fnw5sM16+v65FO7ux1r1OkwTSm1Tin1hlLqEpdkinfv8uWaXQIc0VrvitmW82vWq4/ISTsrZEWh4mxzNQVMKVUBPAn8o9a6FfhvYAZwNtCEZfa6wfu01ucA1wFfUkpd6pIcp6GUKgI+DDxub8qXa9YXedHulFL3ACHgYXtTEzBZa70I+Gfgj0qpyhyL1de9y4trBtxCzwFJzq9ZnD6iz13jbEv7mhWyomgAJsW8nwg0uiQLSikfVgN4WGv9ZwCt9RGttam1DgO/I0vmdn9orRvt/83AU7YcRxxT1v7f7IZsWMprrdb6iC1jXlwz+r4+rrc7pdStwAeBT2vboW27dY7br9dgxQFm5VKuBPcuH66ZF7gJ+JOzLdfXLF4fQY7aWSEriveAmUqpafao9GbgWTcEsX2f9wPbtNY/jdke61O8Edjc+9gcyFaulBrmvMYKhm7Gula32rvdCjyTa9lseozy8uGa2fR1fZ4FPmdnpVwAnHJcB7lAKXUt8DXgw1rrzpjtY5RSHvv1dGAmsDdXctnn7evePQvcrJQqVkpNs2VblUvZgPcD27XWDc6GXF6zvvoIctXOchGxz9c/rMyAnVgjgXtclONiLLNwI7De/rse+AOwyd7+LDDOBdmmY2WcbAC2ONcJGAW8Duyy/490QbYy4DhQFbMt59cMS1E1AUGskdzn+7o+WC6BX9ltbhOwJMdy7cbyXTvt7Df2vh+17+8GYC3wIReuWZ/3DrjHvmY7gOtyKZe9/ffAnb32zdk1S9BH5KSdycxsQRAEISGF7HoSBEEQkkAUhSAIgpAQURSCIAhCQkRRCIIgCAkRRSEIgiAkRBSFIAiCkBBRFIIgCEJCRFEIgiAICfn/yTk+P3Li6W4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "alpha = 0.1 # вероятность появления наблюдения из примеси \n",
    "sigma_1 = 0.5 # разброс основного распределения\n",
    "sigma_2 = 10 # разброс распределения-примеси\n",
    "N=200 # размер выборки\n",
    "noise = gaussian.rvs(alpha, sigma_1, sigma_2, size=N)\n",
    "\n",
    "plt.plot(noise)\n",
    "plt.ylabel('noise')\n",
    "plt.grid()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Подход Хьюбера"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Рассмотрим задачу идентификации коэффициентов линейного регрессионного объекта: <br> \n",
    "    <br> $$ \\large y(i) = b_0 + b_1 u_1(i)+...+b_n u_n(i) + \\eta(i) $$ <br> в предположении, что\n",
    "распределение $ \\large f (\\eta)$ помех принадлежит классу $ \\large \\alpha$-загрязненных\n",
    "распределений:<br><br>\n",
    "$$ \\large Ф = \\{ f (\\eta) : f (\\eta) = (1 - \\alpha)h(\\eta) + \\alpha g(\\eta)\\} $$<br>\n",
    "где $ \\large h(\\eta) $ — некоторая известная плотность распределения,<br>\n",
    "$ \\large g(\\eta) $ — произвольная неизвестная плотность распределения и <br> $\\large \\alpha $ —\n",
    "вероятность появления «выброса» с распределением $ \\large g(\\eta) $, $\\large \\alpha $ удовлетворяет условию\n",
    "$\\large 0  \\leq  \\alpha  \\leq 1 $"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### __Теорема Хьюбера (Хубера)__ .\n",
    "Пусть $\\large h(\\eta)$ — дважды непрерывно дифференцируемая плотность распределения, такая, что $ \\large (- ln (h(\\eta)) $ — выпуклая вниз функция. <br>\n",
    "Тогда АМКО линейных регрессионных объектов — $ \\large V(F, f ) $ имеет седловую точку, т.е. существует плотность распределения $ \\large f^* (\\eta) = (1 - \\alpha)h(\\eta) + \\alpha g^* (\\eta) $ и функция $ \\large F^* (\\eta) = - ln f^* (\\eta) $, такие, что:<br><br>\n",
    "$$\\large V(F^*, f ) \\leq V(F^*, f^*) \\leq V(F, f^*) $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Далее, пусть  $\\large \\eta_0$ и $\\large \\eta_1$ ( $ \\large \\eta_0 < \\eta_1 $) — концы интервала (один или\n",
    "оба конца могут быть бесконечными), где $$ \\large |\\frac{h'(\\eta)}{h(\\eta)}| \\leq k, \\eta [\\eta_0, \\eta_1], $$\n",
    "\n",
    "и $ \\large k, \\alpha, \\eta_0,  \\eta_1 $ связаны соотношением:\n",
    "$$ \\large (1-\\alpha)^{-1} = \\int_{\\eta_0}^{\\eta_1}h(\\eta)d\\eta + \\frac{h(\\eta_0)+h(\\eta_1)}{k} $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Тогда плотность $ \\large f^*(\\eta) $ имеет вид:\n",
    "$$ \\large f^*(\\eta) = \\left\\{\n",
    "    \\begin{array}\\\\\n",
    "        (1-\\alpha)h(\\eta_0)exp(k(\\eta - \\eta_0)), & \\eta \\leq \\eta_0 \\\\\n",
    "        (1-\\alpha)h(\\eta), & \\eta_0 \\leq \\eta \\leq \\eta_1\\\\\n",
    "        (1-\\alpha)h(\\eta_1)exp(-k(\\eta - \\eta_1)), &  \\eta \\geq \\eta_1\n",
    "    \\end{array}  \n",
    "\\right. $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Соответствующая оптимальная функция потерь:\n",
    "    $$ \\large F^*(\\varepsilon) = \\left\\{\n",
    "    \\begin{array}\\\\\n",
    "        -ln[(1-\\alpha)h(\\eta_0)] - k(\\varepsilon - \\eta_0), & \\varepsilon \\leq \\eta_0 \\\\\n",
    "        -ln[(1-\\alpha)h(\\varepsilon)], & \\eta_0 \\leq \\varepsilon \\leq \\eta_1\\\\\n",
    "        -ln[(1-\\alpha)h(\\eta_1)]+k(\\varepsilon - \\eta_1), &  \\varepsilon \\geq \\eta_1\n",
    "    \\end{array}  \n",
    "\\right. $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Случай, когда $\\large h(\\eta) $ - нормальное распределение"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$\\large h(\\eta) = f_N(\\eta) = \\frac{1}{\\sigma_\\eta\\sqrt{2\\pi}}exp(-\\frac{\\eta^2}{2\\sigma_\\eta})$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$ \\large F^*(\\varepsilon) = \\left\\{\n",
    "    \\begin{array}\\\\\n",
    "        -\\frac{\\Delta}{\\sigma^2_\\eta}\\varepsilon - \\frac{1}{2}\\frac{\\Delta^2}{\\sigma^2_\\eta}, & \\varepsilon \\leq -\\Delta \\\\\n",
    "        \\frac{1}{2\\sigma^2_\\eta}\\varepsilon^2, & | \\varepsilon | < \\Delta \\\\\n",
    "        \\frac{\\Delta}{\\sigma^2_\\eta}\\varepsilon - \\frac{1}{2}\\frac{\\Delta^2}{\\sigma^2_\\eta}, &  \\varepsilon \\geq \\Delta\n",
    "    \\end{array}  \n",
    "\\right. $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "где $ \\large \\Delta = \\xi \\sigma_\\eta $"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Связь вероятности выброса $\\large \\alpha$ и множителя $\\large \\xi$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\\large \\frac{1}{1-\\alpha} = \\int_{-\\xi}^{\\xi}\\frac{1}{\\sqrt{2\\pi}}exp(-\\frac{x^2}{2})dx + \\frac{2}{\\xi\\sqrt{2\\pi}}exp(-\\frac{\\xi^2}{2})$$\n",
    "\n",
    "<img src=\"./img/alpha_xi.png\" />"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Чувствительность МНК к выбросам. Пример МНК и HuberRegressor"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$ MSE : \\sum (y_i - \\hat{y}_i)^2 - квадратичная\\_\\_функция\\_\\_потерь$\n",
    "<br>$ MAE : \\sum |y_i - \\hat{y}_i| - модульная\\_\\_функция\\_\\_потерь$ \n",
    "<br> $ Huber: \\begin{equation} L_\\delta(y_i, \\hat{y}_i))= \\left\\{\n",
    "    \\begin{array}\\\\\n",
    "        \\frac{1}{2}(y_i - \\hat{y}_i)^2, & |y_i - \\hat{y}_i| < \\delta \\\\\n",
    "        \\delta|y_i - \\hat{y}_i| - \\frac{1}{2}\\delta^2, & иначе\\\\\n",
    "    \\end{array}  \n",
    "\\right.\\end{equation}- функция\\_\\_потерь\\_\\_Хьюбера$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"./img/losses.png\" />"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.linear_model import LinearRegression, HuberRegressor\n",
    "import random"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 50\n",
    "x = np.linspace(0, 10, n)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "y = 4 + 0.5*x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "noise = np.random.normal(0, 0.5, n)\n",
    "y_new = y + noise"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x7f2ea71aa828>"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAW4AAAD8CAYAAABXe05zAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xt81OWZ9/HPDUgOchQVUKrBylMFlEPAqghIoR4K4oJuRdTWVsuiVtGu62PVigmLrfvQLtjCuD4qWlsLraArtg9SW7I0fYoYxEU5tTVQTQ2IgJAAwUCu/eOXYAiTOWV+M/Ob+b5fr7wkk9/MXHfAK3fu33VftzMzREQkONqlOwAREYmPEreISMAocYuIBIwSt4hIwChxi4gEjBK3iEjAKHGLiASMEreISMAocYuIBEwHP1705JNPtqKiooSeu3//fk488cTkBpThNObsl2vjBY05XmvXrv3YzE6J5VpfEndRUREVFRUJPbesrIxLL700uQFlOI05++XaeEFjjpdz7m+xXqulEhGRgFHiFhEJGCVuEZGA8WWNO5z6+nqqqqqoq6uLeF3Xrl3ZtGlTiqLKDLGOOT8/nz59+nDCCSekICoRyVQpS9xVVVV07tyZoqIinHOtXldTU0Pnzp1TFVZGiGXMZsauXbuoqqqib9++KYpMRDJRypZK6urq6NGjR8SkLa1zztGjR4+ov7GIBE11TTWjnx3N9trt6Q4lMFK6xq2k3Tb6/kk2mrVqFuXvl1P6X6XpDiUwUrZUIiLSXMHsAuoOf/YbZKgiRKgiRH6HfA4+eDCNkWU+VZUk6Ctf+QqffPJJxGsefvhhXn/99YRev6ysjAkTJiT0XJEgqLyrkqkDp1LYoRCAwg6F3HDeDWydsTXNkWU+zbjjZGaYGb/5zW+iXltaql/9RFrTu3NvuuR1oe5IHfkd8qk7UkeXvC706tQr3aFlvJyZcX/ve99j3rx5Rz9/8MEHefzxx4+77kc/+hEDBw5k4MCBzJ07F4Bt27Zx7rnncvvttzN06FA++OADioqK+PjjjwGYNWsW55xzDl/+8pe5/vrrmTNnDgA333wzL774IuC1AZg5cyZDhw7lvPPOY/PmzQCsWbOGcePGMWTIEC6++GK2bNni6/dBJJPs2L+D6cXTWX3LaqYXT9cNyhilZ8Z9N/B2+C8VHCmA9gm85mBgbutfvuWWW5g8eTIzZsygoaGBRYsWsWbNmmOuWbt2LQsXLuSNN97AzPjiF7/I6NGj6d69O1u2bGHhwoUsWLDgmOdUVFSwZMkS1q1bx+HDhxk6dCjFxcVhYzj55JN56623WLBgAXPmzOGpp57inHPOYfny5XTv3p3XX3+dBx54gCVLliTwDRAJnqXXLT365/nj56cxkmDJmaWSoqIievTowbp169ixYwdDhgyhR48ex1xTXl7OpEmTjnb3mjx5Mn/4wx+YOHEiZ555JhdeeOFxr1teXs7VV19NQUEBAFdddVWrMUyePBmA4uJili71/sHu3buX22+/na1bt+Kco76+PinjFZHslZ7EHWFmfLDmoG8bcG699VaeffZZtm/fzje/+c3jvm5mrT63tVaNkZ7TUl5eHgDt27fn8OHDgLeEM3LkSJYtW8a2bdtyrpuaiMQvZ9a4ASZNmsTy5ct58803ufzyy4/7+qhRo3j55Zc5cOAA+/fv56WXXmLkyJERX/OSSy5h2bJl1NXVUVtby69//eu4Ytq7dy+nnXYaAM8++2xczxWR3JQzSyUAHTt2ZMyYMXTr1o327Y9fSB86dCg333wzF1xwAeDN0IcMGcK2bdtafc3hw4czceJEBg0axJlnnsmwYcPo2rVrzDHdd9993HTTTYRCIb70pS/FPSYRyUFN5W3J/CguLraWNm7ceNxj4ezbty+m6xJx5MgRGzRokP35z39O6uvW1NSYmdn+/futuLjY1q5dG9fz4xlzrN/HTLdy5cp0h5BSuTZeM405XkCFxZhjc2apZOPGjZx99tmMHTuWfv36JfW1p02bxuDBgxk6dCjXXHMNQ4cOTerri4g0lzNLJf3796eystKX137hhRd8eV0RkXByZsYtIpItlLhFRAJGiVtEJGCUuEVEAkaJu43a0rq1ueZNq0REIsnoxB2EI41KS0sZN25cusMQkRyS0Yk7mUcaxdLWtal967e+9S0GDBjAZZddxsGD3kkcb7/9NhdeeCHnn38+kyZNYs+ePcCxrVvvv/9++vfvz/nnn8+9994LwM6dO7nmmmsYPnw4w4cP549//CMAu3bt4rLLLmPIkCHMmDEjrp4nIpLbYkrczrkZzrl3nXMbnHN3+x1UwewCXIkjVBGiwRoIVYRwJY6C2QUJv+Ytt9zCc889B3C0resNN9xw3HV/+ctfuOOOO9iwYQPdunU72mL1a1/7Go899hjr16/nvPPOo6Sk5Jjn7d69m5deeokNGzawfv16HnroIQBmzJjBPffcw5tvvsmSJUu49dZbASgpKeGSSy5h3bp1XHnllbz//vsJj01Eki+Tf+OPmridcwOBbwEXAIOACc655G49bMGPI42at3VdsWJF2LauAH379mXw4MGA135127Zt7N27l08++YTRo0cD8PWvf51Vq1Yd87wuXbqQn5/PrbfeytKlSyks9GJ//fXX+fa3v83gwYOZOHEi+/bto6amhlWrVnHjjTcCcMUVV9C9e/eExyYiyZfJhxjHsnPyXGC1mR0AcM79FzAJ+De/gvLrSKNobV3hs9ar4LVfbVoqiaZDhw6sWbOG3/3udyxatIif/OQn/P73v6ehoYE//elPR/t1N6dT20UyT8KHGK+B0146DS71P8ZYEve7wGznXA/gIPAVoKLlRc65acA0gJ49e1JWVnbM17t27UpNTU3UNzty5Ag1NTX8fe/f+eb53+Qb53+DhesXUvVJVUzPj2TcuHE89NBDHD58mP/4j/847vVqa2tpaGg4+vihQ4c4dOgQ7dq1o2vXrrz22mtcfPHFPPXUU1x00UXU1NRQX1/PwYMHqa6u5uDBg4wcOZIBAwYwePBgampqGDNmDD/84Q+ZMWMGAOvXr+f888/nwgsv5JlnnuG+++5j+fLl7Nmzh9ra2mN+cIRTV1d33Pc2iGpra7NiHLHKtfFCcMf8s2E/I1QZovzjcg41HCKvXR4jTx7JbZ+/Lex4Ct8vpO/TfTll1SmccdIZrLpyFQ35Db7GGDVxm9km59xjwG+BWuC/gcNhrnsSeBJg2LBh1vJAgE2bNsV0QEJNTQ2dO3fmlRteOfrYxWddHPV5sRo7dizdunWjW7dux32tU6dOtGvX7miceXl51NfX07lzZ55//nmmT5/OgQMHOOuss1i4cCGdO3fmhBNOODqbnjJlCnV1dZgZc+fOpXPnzoRCIe644w5GjBjB4cOHGTVqFE888QSzZ8/m+uuvZ/To0Vx00UWcccYZdOrUKer3KD8/nyFDhiTt+5EuZWVlOXVoRK6NF4I95tdffZ2VO1eS3yGfT498Sr8z+jH5ssnHXlQFlADPAIXAI/Dm8DcZdcUo3+OLqcmUmT0NPA3gnHsUL+TAaWhoYPXq1fzqV78K+/WioiLefffdo583VYYADB48mNWrVx/3nOaHH7Q8wxK8cyYXL1583OM9evRgxYoVgPfDav58nbcnkimaDjGeVjyNJ9c+SXVt9Wdf3A38APgxcAS4E3gQOAWOlB1JSXwxJW7n3Klm9pFz7gxgMnCRv2El38aNG5kwYQKTJk1KeltXEckuYQ8x3g88DtU/qWbKZVNYfP1iej3cC4pSH1+sbV2XNK5x1wN3mNkeH2PyhZ9tXUUki9XjrTeUANth1u2zKD+1nNJhpSwoWpCWkGJdKol88GKMzEyVFG2gTToiKdQA/BJ4CHgPCr5XQF37BKpNfJCynZP5+fns2rVLySdBZsauXbvIz89PdyiShTJ5s0nKGfAaMAy4HigAlkHlvcnfX5KolJ2A06dPH6qqqti5c2fE6+rq6nIuOcU65vz8fPr06ZOCiCTXNN9ssmB8en79zwhrgPuBlXhr18/jJe/20Bt/9pckImWJ+4QTTqBv375RrysrK8uKcrd45OKYJTMkvNkk22zGWxJZApwCPI63K6XFtoqI1SYplDNnTorI8SrvquTeFffy8uaXOXD4AIUdCpl07iTmXDYn3aGlRsta7BLgHqCV7RRhq03SQIlbJIf51V4i4zWvxW4A7gIewJttB0BGt3UVEf81/fq/+pbVTC+enlE3KJN+03Q/8H3gLGAOcB1Ur6tm9KDRbC/InHFHoxm3SI7LlF//w0naTdMWtdhMBGYDA2HWr4N3Y1aJW0QyTtJumjYAv8K78fhX4BLgRWBE43ssCeaNWS2ViEjGaXNPfgNWAMOBKXi12K8Cq4ARSXqPNFLiFpGM06abpmuAccDleDchnwfWAeOBZhu3g3xjVolbRDJS3DdNNwPXAl8E3gHmNT52I9A+Se+RIbTGLSIZKeabpmH6YvMdWq3FTug9Moxm3CISTLuB+4B+wE+BbwPvATOJKWknW3VNNTPenpGSWbsSt4gES8ta7K8CW4B5UF2QvmZZs1bN4p2976TkcGElbhFJSMo7CtYDTwBn4+1yHA2sB57j6GEG6TiZvWB2Aa7EEaoIYRihihCuxFEw+/gDwpNFiVtEEhJvkkw40TcAi4H+wG14ibsc+E9goHdJ8+TZYA0pSZ5N0lFWqMQtInFJNEnGPRuOoRa7SbTk6edvB83LCju265iSskIlbhGJS7wzzIQS/RpgLJ/VYv+UsLXYTaLVZPu9hNJUVjh/yPyUlBUqcYtIXKIlyabZ7e5PdwNxJvrmtdjv8lkt9k20WovdJFxNdqqWUJZet5T54+dzdqezmT9+/jFlhn5QHbeIxC3SgQJNs9tun3ZjMpNj26HYhlrsJuFqsqtrqrOy37gSt4jELVySbNkY6pXqV3AljvwO+Vx59pXhE33Lvth3Ag+StL7YQd7WHokSt4gkRcvTdPLa5XHtgGuZc9mcYxLl/PHzP6vFfgzYh7cUUsLRsr5kypTjxpJJiVtEkqLl7PbQ4UPHz25b9sW+CniUo2V9fgjqtvZIYro56Zy7xzm3wTn3rnPuF8653DqGXURi0vwG4cTTJn5WXdFaLfYr+Jq0s1XUGbdz7nS8E9n6m9lB59wv8aoqn/U5NhEJmOaz27v73c2loy+F14Dv4pXzDQSW0WpZn8Qm1nLADkCBc64D3j3fD/0LSUSyQeeNnb1a7CuAPXi12G8DE1DSbqOoidvM/o7XyuV9oBrYa2Yr/A5MRAJqM3ANFN9RHHcttsTGmVnkC5zrDiwBrgM+wTvB7UUz+1mL66YB0wB69uxZvGjRooQCqq2tpVOnTgk9N6g05uwX5PHuOrSL0k2lzOw/k5M6ntTqdXk78yh6tohey3txJO8I7016j49u+IgjhUdSGG16teXvecyYMWvNbFhMF5tZxA/gH4Gnm33+NWBBpOcUFxdbolauXJnwc4NKY85+QR7vba/eZu1K2tltr94W/oJdZnavmeWZWUczm2FmHwV7zIlqy5iBCouSj5s+YikHfB+40DlXCBzEW7WqSOAHiogESNST1vfjLYP8G77XYsuxYlnjfgPvQPu38E5yawc86XNcIuKDeLrktdpj5PatEMIr6XsQGMVxfbHFXzFVlZjZTDM7x8wGmtlNZnbI78BEJLp425XG0yUv7Hbx97vQa3gvuB34PPAHVIudBuoOKBJgsSbiRLvkHd1Qc+5qpr83ne1/2u71xV6Gl7QvSdpQJA7a8i4SQFHXn1to2Uck1i55S/suhfuBlTC/aD6UAi+gsr4004xbJIDiPcwg7i55jbXYifTFFv8pcYsEQMu17ETalYY7aOA4VcCtwAC8Y8MeAd7Da3qRl+xRSaK0VCISAM3XsheMXwDE3640Ype83XhtVn+Md9ZjkvtiS3IpcYtksFjXshNuV6pa7EDSUolIBot3LTtm9Rxbiz0a1WIHiGbcIhks6UdvNQC/BB7CW7u+BG973YhkRSypoBm3SIaL6aZiNIZ3s3E4cD1ec+ZXgVUoaQeQZtwiPqmuqWbKkiksvnZxm16nzUdvreFoLTZFeH2xp6KyvgDTjFvEJ/FsL/fFZuBaPqvFfhzVYmcJzbhFkixsJQgh8v8Yfldj0lXhVYYsxNueXgLcA3T2/60lNTTjFkmycJUg404d1/ZKkGh2A/cB/aB6STWj7x/N9ne3w8MoaWcZJW6RJAtXCVLYvjDxSpBo9uNtnjkL75DB62DWglmUdyyn9N00LdP4LN6uiNlGiVvEBy0rQfbU70n+m9QDT+DVYj8AjIaCkgJcX0doS3xdAIMm7fcP0kxr3CIxal4lEm323LISpKysLHmBNOCd/PoQ8FeOqcWurEmsC2BQxNsVMVtpxi0So7TP8prXYk/Bu/HYohY76Rt2MoxvO0kDRjNukSgyYpbXshb7ebyNNGHK+uJtPhUk2f6DKVZK3CJRJHoIQVJsxlsSWYLXqW8e8E9EbLHa5g07GS6bfzDFSolbJIq0zPKaarGfwduePhP4Z1TWR/b/YIqFErdIDFI2y9sN/ACvL/YRvL7YDwCn+vN2EkxK3CIx8H2Wtx+q51UzZesUFi9aTK+v9lJfbGmVqkpE0qlZLfas/z+L8tPLKZ1fqr7YEpESt0g6NACLgf5QUFWAm+4IDQ/R4BoIbc3OTTOSPFETt3PuC865t5t97HPO3Z2K4ESyTsta7HyoLFZtssQnauI2sy1mNtjMBgPFwAHgJd8jE8k2a4CxwOXALry+2G9D70mqTZb4xLtUMhZ4z8z+5kcwIllpMwyYOeCzvtjzgC0c0xc7KafcSM5wZhb7xc49A7xlZj8J87VpwDSAnj17Fi9atCihgGpra+nUqVNCzw0qjTk75e3M48znzqT3/+vNkbwjfHDdB1T9YxVHCo+kO7SUyIW/45baMuYxY8asNbNhMV1sZjF9AB2Bj4Ge0a4tLi62RK1cuTLh5waVxpxldpnZv5hZvpmdYGYzzMqXlqc3pjTI6r/jVrRlzECFxZiP41kquRJvtr0j3p8kItnsaG/oj7Yf2xf7q8CfgblQ370+rTFKdokncV8P/MKvQESCalbZLMr/Vk7pHaXeLsdRwHpUiy2+iWnnpHOuEPgyXnsbESFM18CBIUIDG7sGDsyd3tCSejHNuM3sgJn1MLO9fgckkgptOvqqsRa78teVTF0/lcLDqr+W1NLOSclJCR+K0KwWu3d1b7qM6ELdCaq/ltRSkynJKQkfitBKX+wdL+9geqfc7g0tqafELTkl7kMRWvbFfgT4Dkf7Yqs3tKSDlkokK8S6Zh3zoQi7gfuAfnhb0+8EKvEONNBhBpJmStySFVpbsw6X0CNuL98PPMqxtdhbgLl4SyQiGUBLJRJo0dasmyf0BeMXAK0sb9QDTwGlwHbgKrwEPjBFAxGJg2bcEmiVd4VviWpmuBJHqCJEgzUQqmilx3UDsAg4F7gdOBsoB15BSVsylhK3BFpra9ZbZ2yN3OPagNeAYXh7gguAZcAqYEQ6RiISOyVuCbxwa9YRb0I21WJfAezhaF9sJgAujQMRiZHWuCXwWivJO+5k9upquAZYyjG12OSlOGCRNlLilqx1NKFXwfyX5sNCwtZiiwSNlkokbdrULyQWu4F/wbvh+DyqxZasocQtaZNwv5Bomtdi/xC4DtViS1bRUomkXML9QqJRLbbkCM24JeUq76pk7KljWy/Vi5dqsSXHKHFLyvXu3JsT258YvV9INI19sVWLLblGiVvSYk/9ntb7hbTimJuZzfpiqxZbco3WuCUtSgeUcumllwKxt0OdtarxbMcHSlnw4wWqxZacpcQtGe+4m5k9QoQeCZHfPp+Dd+lsR8k9WiqRzLYbKvdVMnXDVArrW9zMvFtnO0puUuKWzHQA+D5wFvT+P73p8rnUnO3o+6YgkSRQ4pbMUg88gVfS9wAwGlgPOy7YwfRh8d3MTIRvm4JEkkhr3JI01TXVTFkyhcXXLo5/NtwA/ArvQN6/Apc0ft5Y1rd0oL9nO/q2KSiMNn2fRIhxxu2c6+ace9E5t9k5t8k5d5HfgUnwJDRbbarFHg5MwavFfpWU12K3diBDwpuCItCsXtoq1hn3PGC5mV3rnOuI12NNclC42WLCs9U3gfuB3wNFeI2grgfa+xZ+q2I+RLgNUjmrl+wWdcbtnOsCjAKeBjCzT83sE78Dk8wUbrYY92x1MwyYOQAuAN7BmxZsBm4kLUm7ScRDhJMglbN6yW6xzLjPAnYCC51zg4C1wAwz2+9rZJJRos0WY5qtVgElwDPQPa97xvXFbu1AhmRJxaxecoMzs8gXODcMWA2MMLM3nHPzgH1m9r0W100DpgH07NmzeNGiRQkFVFtbS6dOnRJ6blAFYcy7Du0iVBmi/ONyDjUcIq9dHiNPHsltn7+NkzqexMMbHuakjicxofcEXq1+ld2f7qZ0gDcr77CvA2e8cAanv3Q6rsHx4dUfsvEfNtKxT8c0jyp1mv6OI32fsk0Q/l0nW1vGPGbMmLVmNiymi80s4gfQC9jW7PORwK8jPae4uNgStXLlyoSfG1RBGfP0ZdOtXUk7y//XfGtX0s5ue/W2yE+oNbNHzayrmTkz+5qZbfW+FJQxJ0uujddMY44XUGFR8nHTR9Q1bjPbDnzgnPtC40NjgY0J/ECRFPBzA0nMa8Ata7FHAf8NPId3E1JE2iTWqpI7gZ83VpRUAt/wLyRpi+Y3DxeMX5DU1466BtyyFntE4+eXJDUMkZwXU+I2s7fxuh5LhkprqZkBvwW+C7yFd3jBMmA8arEq4gNtec8SaSs1a94Xexfqiy2SAkrcWSLlpWabgWuBLwLv4tVibwFuIq212CK5QIk7i/i9gQTwarG/BQwAXsOrxX4PuAsdZiCSImoylUV83UCyG/gB8GPgCN7t6geAU5P7NiISnRK3RLYfeBx4DNiHtxRSgsr6RNJIiVvCq8frTlMCbAeuAmYD56UzKBEBJW5pSbXYIhlPNyfF07Ivdj7wCvAHlLRFMowSt7Rei30VSanF1jmOIsmlxJ3LUlSLrRNfRJJLiTsXtazFnokvtdgFswtwJY5QRYgGayBUEcKVOApmF0R8nmboIpEpcadQOhNSdU01o//vaLb/7+3QD2855E68lmGP4MthBoluw9cMXSQyJe4USltC2g+z5s6ivKqc0g2l8FW8JZG5wCn+vW282/ATnaGL5Bol7hRIW0Kqh4KSAtwcR6ghREO7BkLDQ7izHAU/T00yjGcbvs5kFImNEncKpDwhNQCLgf5QOaeSqdVTKWyXnmS49LqlzB8/n0G9BjF//PxjtuW3pDMZRWKjxJ0C0RLSrkO7krP2HaYWu/cvetNlQhfqLPx7Z9qNwJQ0yhIJOCXuFImUkH76t5+2fe27eS32bo7pix3pvTPtRmA8M3SRXKUt7wmqrqlmypIpLL52cUy/yofr3JeUU2s2Aw8CS/FuNM4D/oljyvp8e28RSQvNuBOUjJlq09p3Xjsvy8a1/lwF3IpXi72CuPti60agSHBpxh2nZM5Um9a+P234NPabcbuB7+P1xTa8WuwHibusTzcCRYJLM+44JXumumP/DiaeNjH6zbj9wKPAWcAPgetocy22bgSKBJNm3HFK9kx16XVLKSsrO3oz7jj1wFNAKZ/1xX4U7yT1NvL1xBwR8Y1m3AlIyUy1AVgEnAvcDpwNlOO1Wk1C0haR4Ippxu2c2wbU4J02eNjMhvkZVKbzdaZqwG+B+4F1eEl6GTCepLRYFZHgi2epZIyZfexbJOLVYt8PrMQ70/GnwFSS2mJVRIJPa9wZoOD9Aq8v9hJarcUWEWkS6xq3ASucc2udc9P8DCinNPbFvuAbF3h9sR8hplrsTNumLiKp5cws+kXOnWZmHzrnTsVbgb3TzFa1uGYaMA2gZ8+exYsWLUoooNraWjp16pTQc4Oiw74OnPHCGZz+0uk4c2y7chvV36imvlt9TM//9z//O8uql3FV76u453/d43O0/siFv+fmcm28oDHHa8yYMWtjvn9oZnF94M0L7410TXFxsSVq5cqVxz324b4PbdTCUVZdU53w62aE/Wb2qJl1NTNnZl8zs63hx2x2/Ljz/zXfeITjPvL/NT8l4SdTa2POVrk2XjONOV5AhcWYh6MulTjnTnTOdW76M3AZ3gmFKZNpjZDiVg88gVfS9wAwClgPPId3E7IVLcetbeoiArHdnOwJvOSca7r+BTNb7mtUjYLYCOmY5lOFveBXwEPAX4FL8D4fEfk1Io375kE3a5u6SI6LOuM2s0ozG9T4McDMZqciMAjmDPPoLPlnpcf0xWYZsIqoSRsij1vb1EUko8sBg9QI6bhZ8o4QoYkh8q9u/O0gjlrsSOPWNnURyfgt70GZYVaOr2TqnqkUfto4S6aQG/rfwNbvbE1oA01Qxi0iqZfRM24IQCOkKqAEei/sTZerulA3uI789vnUNdTRpbBtzaeaZOS4RSRtMj5xZ6zdwA/w+mI3AN+GHUN3ML3HdKYVT+PJtU9SXVud3hhFJCspccfrAN6W9MeAfcBNQAlQBEvJnllyvEeziUjqZPwad8YIV4v930StxQ6qwNfOi2QxzbijaSChWuygCmLtvEiu0Yy7NYZ3CG+CtdhBFcTaeZFco8QdzhpgLHA53k3InwJvAxPI+sMMglQ7L5KrlLib24zXF/uLeN1Y5jU+dhM5dZiBashFMpvWuOFoLTbPAIV4/Q+/A3ROY0xppBpykcyW24m7ZS32ncCDeKfQiIhkqNxM3PuBxwlbiy0ikulyK3HXA0/jJentwFXAo3gnqYuIBERuJO5wtdgvktVlfSKSvbK7qiSOWmwdwCsiQZG9iTvOWmxt8RaRoMi+pZLNeEsiS/CqQ+YB/wTkhb9cW7xFJGiyZ8ZdBXwLGADV5dWMnjma7eu3w120mrRBW7xFJHiCn7h3A/cB/fCWQ+6EWT+ZRXm7ckrXRl/20BZvEQma4C6VhKnFLuhXQN2ROtjgXRLrskfTFm8dgCAiQZBRibu6ppoZb8/gtWGvtT7jbVmLPQGvFvs8qKyp5N4V9/Ly5pc5cPgAhR0KmXTuJOaT4XCdAAAF4klEQVRcNifi+2qLt4gESUYtlcxaNYt39r4TvrKjAVgM9AduwzvQoByvvO887xIte4hILsiIGXfEyo4HDsJvge8Cb+HtclwGjCdsWZ+WPUQk28WcuJ1z7YEK4O9mNiGZQVTe1coSxylzvFrslcCZeDcfpxKxxaqWPUQk28WzVDID2ORHEM2XODq260jd4Tq6lHWh16hen/XF3kLO9cUWEQknpsTtnOuDtzjxlF+B7Ni/g+lfmM6r619l+pvT2b5ru9cX+z2i1mKLiOSSWJdK5uJVS/t2tMDSy5fC56ChroEv3/Fl7yT1U/16NxGR4HJmFvkC5yYAXzGz251zlwL3hlvjds5NA6YB9OzZs3jRokVxB9NzeU8+7Pch7T+fW+shtbW1dOrUKd1hpFSujTnXxgsac7zGjBmz1syGxXSxmUX8AL6Pt6F8G17l9AHgZ5GeU1xcbIlauXJlws8NKo05++XaeM005ngBFRYlHzd9RF3jNrPvmlkfMyvCa476ezO7MaEfKSIi0mYZtQFHRESii2sDjpmVAWW+RCIiIjHRjFtEJGCUuEVEAkaJW0QkYHIucetQYBEJupxL3PEeCqxELyKZJmcSd8HsAlyJI1QRosEaCFWEcCWOgtkFEZ+n099FJNMEOnHHMxuO91DgRBO9iIjfAp2445kNx3s6jk5/F5FMlREn4MQr4ok5SToUWMegiUimCmTibvXEnCQfCqxj0EQkEwUycadqNqxj0EQkEwV2jbtpNrz6ltVML56ucj0RyRmBnHGDZsMikrsCO+MWEclVStwiIgGjxC0iEjBK3CIiAaPELSISMErcIiIB47xT4ZP8os7tBP6W4NNPBj5OYjhBoDFnv1wbL2jM8TrTzE6J5UJfEndbOOcqzGxYuuNIJY05++XaeEFj9pOWSkREAkaJW0QkYDIxcT+Z7gDSQGPOfrk2XtCYfZNxa9wiIhJZJs64RUQkgoxJ3M65K5xzW5xzf3XO3Z/uePzmnPucc26lc26Tc26Dc25GumNKFedce+fcOufcq+mOJRWcc92ccy865zY3/n1flO6Y/Oacu6fx3/W7zrlfOOfy0x1TsjnnnnHOfeSce7fZYyc5537rnPtL43+7+/HeGZG4nXPtgfnAlUB/4HrnXP/0RuW7w8A/m9m5wIXAHTkw5iYzgE3pDiKF5gHLzewcYBBZPnbn3OnAXcAwMxsItAempDcqXzwLXNHisfuB35lZP+B3jZ8nXUYkbuAC4K9mVmlmnwKLgKvTHJOvzKzazN5q/HMN3v/Mp6c3Kv855/oA44Gn0h1LKjjnugCjgKcBzOxTM/skvVGlRAegwDnXASgEPkxzPElnZquA3S0evhp4rvHPzwH/4Md7Z0riPh34oNnnVeRAEmvinCsChgBvpDeSlJgL3Ac0pDuQFDkL2AksbFweeso5d2K6g/KTmf0dmAO8D1QDe81sRXqjSpmeZlYN3uQMONWPN8mUxO3CPJYT5S7OuU7AEuBuM9uX7nj85JybAHxkZmvTHUsKdQCGAiEzGwLsx6dfnzNF47ru1UBf4DTgROfcjemNKrtkSuKuAj7X7PM+ZOGvVi05507AS9o/N7Ol0a7PAiOAic65bXjLYV9yzv0svSH5rgqoMrOm36ZexEvk2WwcsNXMdppZPbAUuDjNMaXKDudcb4DG/37kx5tkSuJ+E+jnnOvrnOuIdyPjlTTH5CvnnMNb99xkZj9KdzypYGbfNbM+ZlaE93f8ezPL6pmYmW0HPnDOfaHxobHAxjSGlArvAxc65wob/52PJctvyDbzCvD1xj9/HfhPP94kIw4LNrPDzrlvA6/h3YF+xsw2pDksv40AbgLecc693fjYA2b2mzTGJP64E/h546SkEvhGmuPxlZm94Zx7EXgLr3pqHVm4i9I59wvgUuBk51wVMBP4AfBL59wteD/A/tGX99bOSRGRYMmUpRIREYmREreISMAocYuIBIwSt4hIwChxi4gEjBK3iEjAKHGLiASMEreISMD8DyGWKQK6acySAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(x, y, 'magenta', label='y original')\n",
    "plt.plot(x, y_new, 'g*', label='y noised')\n",
    "plt.grid()\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Добавим выбросов в данные"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "outliers = [25.6, 34.2, 44, 29, 50, 100, 90]\n",
    "outliers_ind = random.sample(range(len(x)), len(outliers))\n",
    "\n",
    "for i in range(len(outliers)): \n",
    "    y_new[outliers_ind[i]] = outliers[i]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x7f2ea7124fd0>"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD8CAYAAAB5Pm/hAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xt81NWd//HXh4QwM4YRxS6J0hVcWZHF5ZJ42bqtpLRuFS8N2B+YdH/g0h9EtoWtjb25XZTUbWvzcK3bJA+zWqUtLVTEaLVeWkqq6aPYJoDFghd2goqGgKhIIMMlOb8/MsSQZHKZSyb55v18PHhk5jvn+/2eM9F3zpw53/M15xwiIuJdI1JdARERSS4FvYiIxynoRUQ8TkEvIuJxCnoREY9T0IuIeJyCXkTE4xT0IiIep6AXEfG49FRXAOCss85yEyZMiGnfw4cPc9pppyW2QoOc2jw8qM3DQzxtrqure8c595Heyg2KoJ8wYQK1tbUx7VtdXc2sWbMSW6FBTm0eHtTm4SGeNpvZ630pp6EbERGPU9CLiHicgl5ExOMGxRi9iHjT8ePH2bNnD+FwuE/lTz/9dHbu3JnkWg0ufWmzz+dj/PjxjBw5MqZzKOhFJGn27NnD6NGjmTBhAmbWa/lDhw4xevToAajZ4NFbm51zHDhwgD179jBx4sSYztHr0I2Z/cjM9pnZSx22nWlmvzaz1yI/z4hsNzO718x2mdmfzWxmTLUS6aeGQw1c8dAV7G3am+qqSAfhcJixY8f2KeSle2bG2LFj+/ypqDt9GaN/CPhMp21fBzY65yYBGyPPAa4CJkX+LQEqYq6ZSD+UPFdCzRs1rPrdqlRXRTpRyMcv3vew16Eb59xzZjah0+brgVmRx6uBauBrke0/dm33J9xsZmPMLNs51xBXLUWi8N/pJ3ziw55ORW0FFbUV+NJ9NN/WnMKaiQwe1pd7xkaC/gnn3NTI8/edc2M6vP6ec+4MM3sC+K5zriayfSPwNedcl6uhzGwJbb1+xo0bl7N27dqYGtDU1ERmZmZM+w5VavOHDhw9QEWogpp3ajjaepRRI0bx8bM+zs1/czNnZpyZgpomjhd+z6effjrnn39+n8u3tLSQlpaWxBrFZ968eTzwwAOMGTMmaplvf/vbXH755eTl5fXpmB3b/Pzzz3Pvvffy8MMPdym3a9cuDh48eMq2vLy8Oudcbm/nSPSXsd19vuj2L4lzrhKoBMjNzXWxXhmmK+mGh57a/JsnfsOm/Zvwpfs41nKMSX89iblXzh3YCiaBF37PO3fu7NeXq4P1y1jnHM45nn322V7Lfu973+vXsTu2ORAIkJ6e3u174PP5mDFjRr+OfVKs8+gbzSwbIPJzX2T7HuCjHcqNB96O8RwifdJ4uJGinCI2L95MUU6RvpCVdt/61rf4wQ9+0P78tttu49577+1S7u6772bq1KlMnTqVe+65B4Ddu3dz4YUXsmzZMmbOnMmbb77JhAkTeOeddwAoKSlh8uTJfPrTn+bGG2+ktLQUgEWLFrF+/XqgbXmXlStXMnPmTC666CJefvllAP74xz/ysY99jBkzZvCpT32KV155JanvQ6w9+seBhcB3Iz8f67D9i2a2FrgUOKjxeUm2DfM3tD8um1OWwppIj/4N2NZzEX+LH/ozcjMduCf6y4sXL2bu3LmsWLGC1tZW1q5dyx//+MdTytTV1fHggw/ywgsv4Jzj0ksv5YorruCMM87glVde4cEHH6S8vPyUfWpra3nkkUfYunUrJ06cYObMmeTk5HRbh7POOostW7ZQXl5OaWkp999/P5MnT+a5554jPT2dxx9/nG9+85s88sgj/Wh4//Qa9Gb2c9q+eD3LzPYAK2kL+F+Y2WLgDeBzkeK/Aq4GdgFHgJuSUGcRkT6ZMGECY8eOZevWrTQ2NjJjxgzGjh17Spmamhry8/PbV5CcO3cuzz//PNdddx3nnnsul112WZfj1tTUcP311+P3+wG49tpro9Zh7ty2YcScnBw2bGjrlBw8eJCFCxfy2muv4ZyjpaUlIe2Npi+zbm6M8tLsbso64F/jrZSIeFAPPe+Tmg81J3yM/gtf+AIPPfQQe/fu5V/+5V+6vN7ThJRoywf3ZRLLSaNGjQIgLS2NEydOAG1DSnl5eTz66KO89NJLXHPNNX0+Xiy01o2IeFp+fj5PP/00f/rTn/inf/qnLq9/4hOfoKqqiiNHjnD48GEeffRRPv7xj/d4zH/8x3/kl7/8JeFwmKamJp588sl+1engwYOcc845AKxZs6Zf+8ZCSyCIiKdlZGSQl5fHmDFjup26OXPmTBYtWsQll1wCtH0CmDFjBrt37456zIsvvpjrrruOadOmce6555Kbm8vpp5/e5zp99atfZeHChdx9991cfvnl/W5Tv52cNpTKfzk5OS5WmzZtinnfoUptHh680OYdO3b0q/wHH3yQ8Dq0tLS4adOmuVdffTWhxz106JBzzrnDhw+7nJwcV1dXF9Nx+trm7t5LoNb1IWM1dCMinrVjxw7OP/98Zs+ezaRJkxJ67CVLljB9+nRmzpzJvHnzmDlz8C7tpaEbEfGsKVOmEAqFknLsn/3sZ0k5bjKoRy8i4nEKehERj1PQi4h4nIJeRMTjFPQiIv300EMP8fbbH67XOGvWLGpr21Zjv/rqq3n//fdTVbVuKehFZFAZCreF7Bz0Hf3qV7/qcb36VFDQi8igksjbQsa7TPHUqVPby5SWlnL77bezfv16amtrKSwsZPr06TQ3n3ons45LGf/0pz/lkksuYfr06SxdupSWlhZaWlpYtGgRU6dO5aKLLuKHP/xh3O3sjYJ+AA2FnopIqvjv9BO8O0hFbQWtrpWK2grsDsN/pz/mYy5evJjVq1cDtC9TXFhYeEqZjssUb968mf/5n/9h69atUY95ww03kJuby5o1a9i2bVv7Cpad7dy5k3Xr1vH73/+ebdu2kZaW1r7PW2+9xUsvvcT27dv5/Oc/H3P7+kpBP4B0A2uR6ELLQ3xu8ucIpAcACKQHKLyokPoV9TEfs+Myxc8++2yvyxRnZma2L1Mcr40bN1JXV8fFF1/M9OnT2bhxI6FQiPPOO49QKMSXvvQlnn76aYLBYNzn6o2ujB0AuoG1SO+yR2czOmM04ZYwvnQf4ZYwwVFBsjKz4jpurMsUp6en09ra2v48HA53Wy4a5xwLFy7kO9/5TpfXXnzxRZ555hnKyspYs2YNP/nJT/p17P5Sj34AhJaHKJhakNCeiogX7T+yP+G3hYx1meJx48axb98+Dhw4wNGjR3niiSfa9xk9ejSHDh3q8byzZ89m/fr17NvXdqfVd999l9dff5133nmH1tZW5s2bR0lJCS+++GLcbeyNevQDIHt0NsFRwYT3VES8Zs11a9pvPJKo20LGukwxwH/8x39w6aWXMnHiRCZPnty+z6JFiygqKsLv9/OHP/yh2/NOmTKFb3/721x55ZW0trYycuRIysrK8Pv93HTTTe2fFlauXJmQdvZEQT9ATt7AeknOEirrKmlo0q10RQZCa2srmzdv5uGHH45a5pZbbuGWW27psn358uUsX768y/Z58+Yxb9689ufV1dXtjzuuYz9//nzmz5/fZf8tW7a0P+7tk0EiKOgHiG5gLTLwduzYwTXXXEN+fn7ClykeShT0IuJZyVymeCjRl7EiklTRZrVI38X7HiroRSRpfD4fBw4cUNjHwTnHgQMH8Pl8MR9DQzcikjTjx49nz5497N+/v0/lw+FwXIE2FPWlzT6fj/Hjx8d8DgW9iCTNyJEjmThxYp/LV1dXt09tHC4Gos0auhER8TgFvYiIxynoRUQ8TkEvIuJxCnoREY9T0IuIeJyCXkTE4xT0IiIeF1fQm9mXzewvZvaSmf3czHxmNtHMXjCz18xsnZllJKqyIiLSfzEHvZmdAywHcp1zU4E0YAHwPeC/nHOTgPeAxYmoqIiIxCbeoZt0wG9m6UAAaAA+CayPvL4a+Gyc5xARkThYPKvKmdkK4E6gGXgWWAFsds6dH3n9o8BTkR5/532XAEsAxo0bl7N27dqY6tDU1ERmZmZsDRii1ObhQW0eHuJpc15eXp1zLre3cjEvamZmZwDXAxOB94GHgau6KdrtXxLnXCVQCZCbm+tmzZoVUz2qq6uJdd+hSm0eHtTm4WEg2hzP0M2ngHrn3H7n3HFgA/AxYExkKAdgPPB2nHUUEZE4xBP0bwCXmVnAzAyYDewANgE3RMosBB6Lr4oiIhKPmIPeOfcCbV+6bgG2R45VCXwNuMXMdgFjgQcSUE8REYlRXDcecc6tBFZ22hwCLonnuCIikji6MlZExOMU9CIiHqegFxHxOAW9iIjHKehFRDxOQS8i4nEKehERj1PQi4h4nIJeRMTjFPQiIh6noBcR8TgFvYiIxynoRUQ8TkEvIuJxCnoREY9T0IuIeJyCXkTE4xT0IiIep6AXEfE4Bb2IiMcp6EVEPE5BLyLicQp6ERGPU9CLiHicgl5ExOMU9CIiHqegFxHxOAW9iIjHKehFRDxOQS8i4nEKehERj1PQi4h4XFxBb2ZjzGy9mb1sZjvN7B/M7Ewz+7WZvRb5eUaiKisiIv0Xb4/+B8DTzrnJwDRgJ/B1YKNzbhKwMfJcRERSJOagN7Mg8AngAQDn3DHn3PvA9cDqSLHVwGfjraSIiMTOnHOx7Wg2HagEdtDWm68DVgBvOefGdCj3nnOuy/CNmS0BlgCMGzcuZ+3atTHVo6mpiczMzJj2HarU5uFBbR4e4mlzXl5enXMut9eCzrmY/gG5wAng0sjzHwAlwPudyr3X27FycnJcrDZt2hTzvkOV2jw8qM3DQzxtBmpdH/I6njH6PcAe59wLkefrgZlAo5llA0R+7ovjHCIiEqeYg945txd408wuiGyaTdswzuPAwsi2hcBjcdVQRETikh7n/l8C1phZBhACbqLtj8cvzGwx8AbwuTjPISIicYgr6J1z22gbq+9sdjzHFRGRxNGVsSIiHqegFxHxOAW9iIjHKehFRDxOQS8i4nEKehERj1PQi4h4nIJeRMTjFPQiIh6noBcR8TgFvYiIxynoRUQ8TkEf0XCogSseuoK9TXtTXRURkYRS0EeUPFdCzRs1rPrdqlRXRUQkoeJdj37I89/pJ3wi3P68oraCitoKfOk+mm9rTmHNREQSY9j36EPLQxRMLSCQHgAgkB6g8KJC6lfUp7hmIiKJMeyDPnt0NsFRQcItYXzpPsItYYKjgmRlZqW6aiIiCTHsgx6g8XAjRTlFbF68maKcIn0hKyKeMuzH6AE2zN/Q/rhsTlkKayIiknjq0YuIeJyCXkTE4xT0IiIep6AXEfE4Bb2IiMcp6EU8SGs3SUcKehEP0tpN0pHm0Yt4iNZuku6oRy/iIVq7SbqjoBfxEK3dJN1R0It4jNZuks40Ri/iMVq7STpTj15ExOPiDnozSzOzrWb2ROT5RDN7wcxeM7N1ZpYRfzVFRCRWiejRrwB2dnj+PeC/nHOTgPeAxQk4h4iIxCiuoDez8cAc4P7IcwM+CayPFFkNfDaec4iISHzMORf7zmbrge8Ao4FiYBGw2Tl3fuT1jwJPOeemdrPvEmAJwLhx43LWrl0bUx2amprIzMyMad+hSm0eHtTm4SGeNufl5dU553J7KxfzrBszuwbY55yrM7NZJzd3U7TbvyTOuUqgEiA3N9fNmjWru2K9qq6uJtZ9hyq1eXhQm4eHgWhzPNMrLweuM7OrAR8QBO4BxphZunPuBDAeeDv+aoqISKxiHqN3zn3DOTfeOTcBWAD81jlXCGwCbogUWwg8FnctRUQkZsmYR/814BYz2wWMBR5IwjlERKSPEnJlrHOuGqiOPA4BlyTiuCIiEj9dGSsi4nEKehERj1PQi4h4nIJekk73LxVJLQW9JJ3uXyqSWlqPXpJG9y8VGRzUo5ek0f1LRQYHBb0kje5fKjI4KOglqXT/UpHU0xi9JJXuXyqSeurRi4h4nIJeRMTjFPQiIh6noBcR8TgFvYgklJa8GHwU9CKSUFryYvDR9EoRSQgteTF4qUcvIgmhJS8GLwW9iCSElrwYvBT0IpIwWvJicNIYvYgkjJa8GJzUoxcR6YaXpokq6EVEuuGlaaIauhER6cCL00TVoxcR6cCL00QV9CIiHXhxmqiCXkSkE69NE9UYvYhIJ16bJqoevYiIxynoRUQ8TkEvIuJxCnoREY+LOejN7KNmtsnMdprZX8xsRWT7mWb2azN7LfLzjMRVV0RE+iueHv0J4CvOuQuBy4B/NbMpwNeBjc65ScDGyHMREUmRmIPeOdfgnNsSeXwI2AmcA1wPrI4UWw18Nt5KiohI7Mw5F/9BzCYAzwFTgTecc2M6vPaec67L8I2ZLQGWAIwbNy5n7dq1MZ27qamJzMzMmPYdqtTm4UFtHh7iaXNeXl6dcy6314LOubj+AZlAHTA38vz9Tq+/19sxcnJyXKw2bdoU875Dldo8PKjNw0M8bQZqXR9yOq5ZN2Y2EngEWOOcO3kpWaOZZUdezwb2xXMOERGJTzyzbgx4ANjpnLu7w0uPAwsjjxcCj8VePRERiVc8a91cDvwzsN3MtkW2fRP4LvALM1sMvAF8Lr4qiohIPGIOeudcDWBRXp4d63FFRCSxdGWsiIjHKehFRDxOQS8i4nEKehERj1PQi4h4nIJeRMTjFPQiIh6noBcR8TgFvYiIxynoRUQ8TkEvIuJxCnoREY9T0IvIkNNwqIErHrqCvU17U12VIUFBLyJDTslzJdS8UcOq361KdVWGBAW9iMQkFb1q/51+7A6joraCVtdKRW0Fdofhv9M/YHWIJtr7MRg+fSjoRSQm/e1V9xR4J19799i7Pe4TWh6iYGoBgfQAAIH0AIUXFVK/or7Xc8RSr/6Uj/Z+9PQ+NRxqYMW2FUn/I6CgFxlGEhGEsfaqewq8k6+t3r26x32yR2cTHBUk3BLGl+4j3BImOCpIVmZWr+fotV7Vq+AIcADYA+wCtgN/Ap4DngUeg5IHSqh5vYZVZavgbvDf0f37Ybdb9+/Tt/zwV8BoKCkoYfv725M+BGVtNxJPrdzcXFdbWxvTvtXV1cyaNSuxFRrk1OahoeFQAwseWcC6G9a1B1F/yiejzcueXMZ9dfexNGcp5XPK+1TXzvs0HGqg+Nliql6u4siJIwTSA+RfmE/plaXdttN/p5/wiXCX7b50H0C3r0XjG+HjqjFXkT0imyWjl1D5TiUNxxp46uhThOnmHK0+mv/STENLAwsmLGDdtnVkHcyCZvDf4Cec3s0+x30039nctR23+QmP7Fp+1PFRzNs5j6oLqzgy8giB4wHy38zn1l23ctff3UVVVhVH0o4QaAmQ35RP6QelTBw/kbB1/54039b13NGYWZ1zLre3curRi/RRf3vD0XqX/f3on4g69dYL7+7c0fY5757zCL4fJHwijM98hE+ECb4eJGttFvwQ+D5QQtsdpL8MobdDFDQVEGiNDLe0BChsKKT+qXpCj4YoCBUQON72mv+4n8I/F7KtYhsFfy4gcCyyz7EAhX8upP6uejYs30DZF8uYtnAaZV8pY8M3NhAqDZ1a/niAwlcLqX+4HnZCyekl1Jxew6q/XQXnABdBaE+IgiMFBFxkHxegcEQh9WfVQyXwY+Bh4JfAbyD0yRAF5xQQSIuUTwtQeEEhu7+8m+BNQcIZkU8ZGWGC/yfItN9PI3hdkHB6ZHt6mODsIFn3ZxH6Ss9DUIkWz83BRQZcvL3kWLfDqWHYsTfcWecebEVtBRW1Fe29tc7HiVY+Y0QGR2cdbdvYAoSB5lN/lvyphJo3a1j1k1WUZ5ef+nqHx6GjIYrTi6nKqOLIiCMEWgPkH8hn/dj12B3Wta4nfIT+O0Tx7GKqJldxJOMIgWMB8l/Op/TZUpZdvYyipiKW1C2hMqeShswG+EWnN2Ik4INsfzbBvCDhyWF8LT7CI8IEW4NkBbNgHATPjIShi4Th1CDTcqcRdEHCx8P4iIRnXpCsoizwAf7Iv8jjbF82weeDhF8M40vzEbYwwYIgE6dMPPW9PbuCirM//F0EnwgS3hLZpyVMcGaQrDnd/3eVTTbBI0HCb3cYNhodJGtsFo2HGynKKWJJzhIq6yppaGoAiLq94xBUxoiMLkNQiaagH+b6G5ypPndfw7a38v3ZHjW403w0FzTT8G4DC/6wgHVT1pHVmkXowhDFrxdTdaiKI+4IAQLkWz7rj3cfqqNaR1Gwv4CqsyIf8U8EyH89n//87X9CKdAMDaMaWHDDAtatX0dWU1aXYYSKfRVU7KvoftjBImF7VZDw1DC+Ez7CaWGCh4PUv1JP8QXFVJ0ZOXdrgPyj+ZSeKCXri1kEM4KEXYew/XiQrM9nscG/oT1ky/xlbY9L+TCAR3FKujSua6Qos1Pgze/62spfrmRvcC/Mj7LPp6P/rhvDXUM1tDwUdZgJogdx1HNEKb9h/ob2MmVzytofR9ve8Vgz3Ay22tZezx0PjdEPYtGC8JFnHuHehnsTEs79HbNN5B+G/px7VMkojrUe63IMX5qP5v/X3KUn6/+Nn7Dr+9gvDrCum30nfIR+HqL44mKq/iYyBtuhZ5vVlMWyOcu4L+c+ltYtpfzJtnbcPOdmKnMqyWjJ4FjaMZbWLeVbv/sWxVcVU/W3keOcCJD/Vj6lr5Ryx5Q7qDynkgyXwTE7xtIPlvLvu/6ds887G/ywzLeM+7iPpaOWUj62nIaRDRS/X0zVB5E/JiMC5GfnU5pTStYZWaf0dhnZ1ra56+aSnZl9SkhtmL+Bm5+4mcotlWSkZXCs5dgpv49o+yRLMv5/7ql9g0E8bR72Y/SDYe5qX/V3zPbHr/+4x+lafZnLG8uYbdTtLUATsB94E3gVeBF4AagGngIeBX4G/Cj6LAX/Sj/cCCXfKKFmdw2rvrwKLgH+Hnb8eAcFuz4cy20fs/1ePWQBE4ELgZnA5RD6fqjbMd72sd/jHcZy/7eQbU9vo+DNAgItHcaR3yuk/s16svOyCU7oMLyQESaYG2Ri8UTsdqPi4gpaR7RScXEFdrvhX+Wn8fpGiiYXsXn+ZoqmFbG3cC/ZH2QT/L8dxnJHhgnOC5L1+ywaP9VI0cVFbL55M0UXF7H30r28esur+D/ix3xGBRW00krF0QrsbeO8t84jeH6QMJFjuTDB7CBZM7JgAjAOOB3IoP0P2Ib5GyibU8a0rGmUzSlrD+yTPcvNizdTlFN0yn870fYZSnpq33Dh2R59tN5ioiSzZxttlkI07d/UO1j2+DLu23YfSy9cSnlOeXsPd9m2Zdz39n0sHbuU8o+U09DUQPH+YqqOVnGEIwRcgPwT+axPX89RO9r1JNF6vB2GChoyTx1eiKYhs4HiKzuM/R4PkB/KZ/3frOdoetdz+1p97K7Zze3Tb6dydCUZZHCMYywdsZTyzPJux2zxwc3/ezOVb1aSMSKDY63HWDp1KeWfKefm395M5dauPbz+9mzLri7r1+yTaMeJFp7V1dVckHNB1HMse3LZgPa2B4KXP6FHMxA9+iE9Rn/yYoNncp9p/x+rty/CEqXHsWIHHKPLF2ddhhe2nDq80D5m60ZRcKyAqpEdvjh7L59bd9/KXWffxaMfeZTm9Oa2j/6hfEqrS/GHT50qVrGzgoqdFV0CuuJABRUH2sZyF21bRDgn8gVZWpjg9iD1W+opvqKYqolVHEmPDC/sz+fWvbdy17l3UTUmUicXIP94PqWuFO4A/FDSUkLN0RpWrVpF+V+Xdw3fyONsfzbBmiDhlzp8cVYYpP4T9VFDbWfeThobo4/zdqdxbyNFuZ3KnwaNR/r3xRlEH2vtaU53d3oas+1OT/PG+3ssGb6GdNCXPFfC9oPbTwnbbr98OT+f0stK4Q1oeK+BBc8vYN1F68hyWV1mJ3T+2RBuYIFvAesaI1+CXeAnPKKbPyQnfDSXt40VN6Q3sGBe7z3bUGbo1J7tybHfTaXcMesOwn/f4Yuzg0GmHZhGMDvI0bSj+Foj288OkjUvi9CxEMVpxVRZpIdOgPzT8rl13K3c9d5dVB2s4khrZCx3fD6ll5aybNIyioKRYNteScOUBrIfyz51JoKFCV4bZNqcaV1nKHysbYaC/04/4SMd3pMPKqh4qec/ro3Hu4ZqT6H2Mi/3O9j6+wVZLMHZ3y/zYjEQ5xBvG5JBH7XX3uKj+aFmgpcECf9dpKfaGib4kyBZCyJXzs0poSanhlXfX9X+xVlU6ZHy02tY1bSK8q3lhOpDFE8r/vAiiNYPL4LgRtp6toESakbUsOobqyg/o/zDXm3Hnu3JKWHbg4Rfi4RnRpjgTUGyHsnqOuPgwob2mQjXfnAtd1x7xym92myyuwbxhZ0C+mR4/lWQrClZbJjSIdj++sNg6+9Usd5mNnQnWqgOtVAbiF61eu4SryE5Rt9+Zd5fIjMOWgLkH8yn9J1SsjKymJs1l2zLZsmIJVSOrKQhrYGnWqJcOTfCR/N1zV2GF/z3Rr+ab9G0Rd2O5fZ0BWC0nm0ssxqijelFO9ZAzJxI9swGjd0OD2pz/3h6jL79Iz6Riw0sTPDKDy902ECHHhBtPaCeLtsms+s5euqlLntyWdJ7trFI5JBEfw21nrjIcDIkgx76f7FBb4sh9ad8tODs7zm8RMMLIoPXkA36k8FSXV3NF2Z9oU/7JOoquETvIyKSTEM26GORqFkbid5HRCSZknJlrJl9xsxeMbNdZvb1ZJxDRET6JuFBb2ZpQBlwFTAFuNHMpiT6PCIi0jfJ6NFfAuxyzoWcc8eAtcD1STiPiIj0QcLn0ZvZDcBnnHNfiDz/Z+BS59wXO5VbAiwBGDduXM7atWtjOl9TUxOZmd3Mj/QwtXl4UJuHh3janJeXl7J59N0sfUWXvybOuUra7uNCbm6ui/WCAV1gMTyozcOD2pwcyRi62QN8tMPz8cDbSTiPiIj0QTKGbtJpW5F8NvAWbfdQL3DO/aWHffYDr8d4yrOAd2Lcd6hSm4cHtXl4iKfN5zrnPtJboYQP3TjnTpjZF4FngDTgRz2FfGSfXisajZnV9mWMykvU5uFBbR4eBqLNSblgyjn3K+BXyTi2iIj0j2dvJSgiIm28EPQ+vKwRAAAC50lEQVSVqa5ACqjNw4PaPDwkvc2DYj16ERFJHi/06EVEpAdDOuiH2+JpZvZRM9tkZjvN7C9mtiLVdRoIZpZmZlvN7IlU12UgmNkYM1tvZi9Hftf/kOo6JZuZfTny3/RLZvZzM/Oluk6JZmY/MrN9ZvZSh21nmtmvzey1yM8zknHuIRv0w3TxtBPAV5xzFwKXAf86DNoMsALYmepKDKAfAE875yYD0/B4283sHGA5kOucm0rbtOwFqa1VUjwEfKbTtq8DG51zk4CNkecJN2SDnmG4eJpzrsE5tyXy+BBtAXBOamuVXGY2HpgD3J/qugwEMwsCnwAeAHDOHXPOvZ/aWg2IdMAfueAygAevpnfOPQe822nz9cDqyOPVwGeTce6hHPTnAG92eL4Hj4deR2Y2AZgBvJDamiTdPcBXgdZUV2SAnAfsBx6MDFfdb2anpbpSyeScewsoBd4AGoCDzrlnU1urATPOOdcAbR054K+ScZKhHPR9WjzNi8wsE3gE+Dfn3Aeprk+ymNk1wD7nXF2q6zKA0oGZQIVzbgZwmCR9nB8sIuPS1wMTgbOB08zs86mtlbcM5aAflounmdlI2kJ+jXNuQ2/lh7jLgevMbDdtQ3OfNLOfprZKSbcH2OOcO/lJbT1twe9lnwLqnXP7nXPHgQ3Ax1Jcp4HSaGbZAJGf+5JxkqEc9H8CJpnZRDPLoO3Lm8dTXKekMjOjbex2p3Pu7lTXJ9mcc99wzo13zk2g7ff7W+ecp3t6zrm9wJtmdkFk02xgRwqrNBDeAC4zs0Dkv/HZePwL6A4eBxZGHi8EHkvGSYbszcFjWTzNAy4H/hnYbmbbItu+GVlbSLzjS8CaSAcmBNyU4voklXPuBTNbD2yhbWbZVjx4hayZ/RyYBZxlZnuAlcB3gV+Y2WLa/uB9Linn1pWxIiLeNpSHbkREpA8U9CIiHqegFxHxOAW9iIjHKehFRDxOQS8i4nEKehERj1PQi4h43P8Hkg9im6Jhy58AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(x, y, 'magenta', label='y original')\n",
    "plt.plot(x, y_new, 'g*', label='y outlies')\n",
    "plt.grid()\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Определим параметры $ y=kx+b $ с помощью МНК"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "k= -0.478 b= 15.464\n"
     ]
    }
   ],
   "source": [
    "reg_mse = LinearRegression()\n",
    "reg_mse.fit(x[:, np.newaxis], y_new)\n",
    "print('k=', round(reg_mse.coef_[0], 3), 'b=', round(reg_mse.intercept_, 3))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x7f2ea709ae48>"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD8CAYAAAB5Pm/hAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xt4VNW9//H3SkIyE8OAYkui9AhWjogglwS1pVVSWltFsQH7A5P2Bx78QeS00FKsWk+PSmovNo9HbZM85niBVixUxGjVKpUSMT5FmgAWCiiYoKAhIHJJJENIsn5/ZIi5zCSZWy47n9fzzDOz96y991oT+K611157bWOtRUREnCumpzMgIiLRpUAvIuJwCvQiIg6nQC8i4nAK9CIiDqdALyLicAr0IiIOp0AvIuJwCvQiIg4X19MZADj33HPt8OHDQ9r2008/5ayzzopshno5lbl/UJn7h3DKXFZW9rG19nOdpesVgX748OGUlpaGtG1xcTFTpkyJbIZ6OZW5f1CZ+4dwymyMeb8r6dR1IyLicAr0IiIOp0AvIuJwvaKPXkSc6fTp0xw4cACv19ul9IMGDWLXrl1RzlXv0pUyu1wuhg0bxoABA0I6hgK9iETNgQMHGDhwIMOHD8cY02n66upqBg4c2A056z06K7O1liNHjnDgwAFGjBgR0jE67boxxjxhjDlkjNnRYt05xpi/GmP2+N7P9q03xphHjDF7jTH/NMZMDClXIkGqrK7k6uVXc7DmYE9nRVrwer0MGTKkS0Fe/DPGMGTIkC6fFfnTlT765cC32qy7E1hvrR0JrPctA1wLjPS95gMFIedMJAg5G3Mo+aCEZa8v6+msSBsK8uEL9zfstOvGWrvRGDO8zeobgSm+zyuAYuAO3/rf26bnE24yxgw2xqRYayvDyqVIAO773XjrP2vpFJQWUFBagCvORe3dtT2YM5Hew3TlmbG+QP+itXaMb/mYtXZwi++PWmvPNsa8CPzKWlviW78euMNa2+5uKGPMfJpa/QwdOjR11apVIRWgpqaGpKSkkLbtq1Tmzxw5dYSC8gJKPi7hVOMpEmIS+Oq5X+W2L97GOfHn9EBOI8cJf+dBgwZx0UUXdTl9Q0MDsbGxUcxReGbOnMnjjz/O4MGDA6b5+c9/zuTJk0lPT+/SPluW+Y033uCRRx7hmWeeaZdu7969HD9+vNW69PT0MmttWmfHiPTFWH/nF35rEmttIVAIkJaWZkO9M0x30vUPHZX5tRdfY8PhDbjiXNQ11DHy30Yy45oZ3ZvBKHDC33nXrl1BXVztrRdjrbVYa1m3bl2naX/9618Hte+WZU5MTCQuLs7vb+ByuZgwYUJQ+z4j1HH0VcaYFADf+yHf+gPAF1qkGwZ8FOIxRLqk6tMqslOz2TRvE9mp2bogK81+9rOf8fDDDzcv33333TzyyCPt0j344IOMGTOGMWPG8NBDDwGwb98+LrnkEhYuXMjEiRPZv38/w4cP5+OPPwYgJyeHUaNG8Y1vfIObb76Z3NxcAObOncuaNWuApuld7rnnHiZOnMjYsWPZvXs3AJs3b+bLX/4yEyZM4Otf/zrvvPNOVH+HUFv0LwBzgF/53p9vsf77xphVwBXAcfXPS7StnbW2+XPetLwezIl06IfAto6TuBvcEEzPzXjgocBfz5s3jxkzZrB48WIaGxtZtWoVmzdvbpWmrKyMJ598krfeegtrLVdccQVXX301Z599Nu+88w5PPvkk+fn5rbYpLS3l2WefZevWrdTX1zNx4kRSU1P95uHcc89ly5Yt5Ofnk5uby2OPPcaoUaPYuHEjcXFxvPDCC/z0pz/l2WefDaLgwek00Btj/kjThddzjTEHgHtoCvB/MsbMAz4AvuNL/jJwHbAXOAncEoU8i4h0yfDhwxkyZAhbt26lqqqKCRMmMGTIkFZpSkpKyMjIaJ5BcsaMGbzxxhtMnz6dCy64gCuvvLLdfktKSrjxxhtxu90A3HDDDQHzMGNGUzdiamoqa9c2NUqOHz/OnDlz2LNnD9ZaGhoaIlLeQLoy6ubmAF9N9ZPWAv8ZbqZExIE6aHmfUVtdG/E++ltvvZXly5dz8OBB/uM//qPd9x0NSAk0fXBXBrGckZCQAEBsbCz19fVAU5dSeno6zz33HDt27OD666/v8v5CobluRMTRMjIyeOWVV/jHP/7BN7/5zXbfX3XVVRQVFXHy5Ek+/fRTnnvuOb761a92uM+vfOUr/PnPf8br9VJTU8NLL70UVJ6OHz/O+eefD8DKlSuD2jYUmgJBRBwtPj6e9PR0Bg8e7Hfo5sSJE5k7dy6XX3450HQGMGHCBPbt2xdwn5MmTWL69OmMGzeOCy64gLS0NAYNGtTlPP3kJz9hzpw5PPjgg0yePDnoMgXtzLChnnylpqbaUG3YsCHkbfsqlbl/cEKZd+7cGVT6EydORDwPDQ0Ndty4cfbdd9+N6H6rq6uttdZ++umnNjU11ZaVlYW0n66W2d9vCZTaLsRYdd2IiGPt3LmTiy66iKlTpzJy5MiI7nv+/PmMHz+eiRMnMnPmTCZO7L1Te6nrRkQca/To0ZSXl0dl308//XRU9hsNatGLiDicAr2IiMMp0IuIOJwCvYiIwynQi4gEafny5Xz00WfzNU6ZMoXS0qbZ2K+77jqOHTvWU1nzS4FeRHqVvvBYyLaBvqWXX365w/nqe4ICvYj0KpF8LGS40xSPGTOmOU1ubi733nsva9asobS0lKysLMaPH09tbesnmbWcyvipp57i8ssvZ/z48SxYsICGhgYaGhqYO3cuY8aMYezYsfzud78Lu5ydUaDvRn2hpSLSU9z3u/E86KGgtIBG20hBaQHmPoP7fnfI+5w3bx4rVqwAaJ6mOCsrq1WaltMUb9q0if/93/9l69atAfd50003kZaWxsqVK9m2bVvzDJZt7dq1i9WrV/Pmm2+ybds2YmNjm7f58MMP2bFjB9u3b+e73/1uyOXrKgX6bqQHWIsEVr6onO+M+g6JcYkAJMYlkjU2i4rFFSHvs+U0xevWret0muKkpKTmaYrDtX79esrKypg0aRLjx49n/fr1lJeXc+GFF1JeXs4PfvADXnnlFTweT9jH6ozujO0GeoC1SOdSBqYwMH4g3gYvrjgX3gYvngQPyUnJYe031GmK4+LiaGxsbF72er1+0wVirWXOnDn88pe/bPfd22+/zauvvkpeXh4rV67kD3/4Q1D7DpZa9N2gfFE5mWMyI9pSEXGiwycPR/yxkKFOUzx06FAOHTrEkSNHOHXqFC+++GLzNgMHDqS6urrD406dOpU1a9Zw6FDTk1Y/+eQT3n//fT7++GMaGxuZOXMmOTk5vP3222GXsTNq0XeDlIEpeBI8EW+piDjNyukrmx88EqnHQoY6TTHAf//3f3PFFVcwYsQIRo0a1bzN3Llzyc7Oxu128/e//93vcUePHs3Pf/5zrrnmGhobGxkwYAB5eXm43W5uueWW5rOFe+65JyLl7IgCfTc58wDr+anzKSwrpLJGj9IV6Q6NjY1s2rSJZ555JmCaJUuWsGTJknbrFy1axKJFi9qtnzlzJjNnzmxeLi4ubv7cch77WbNmMWvWrHbbb9mypflzZ2cGkaBA3030AGuR7rdz506uv/56MjIyIj5NcV+iQC8ijhXNaYr7El2MFRFxOAV6ERGHU6AXEXE4BXoREYdToBcRcTgFehERh1OgFxHH6so0xfv27WPUqFHceuutjBkzhqysLF577TUmT57MyJEj2bx5MwCvv/4648ePZ/z48UyYMKH5Rqff/OY3TJo0icsuu6xb7nINhcbRi0i3+OEPYdu2jtM0NLjxM0tBQOPHg2/6eL/mzZvHjBkzWLx4cfM0xWcCd0t79+7lmWeeobCwkEmTJvH0009TUlLCCy+8wC9+8QuKiorIzc0lLy+PyZMnU1NTg8vlYt26dezZs4fNmzdjrWX69Ols3LiRq666quuF6AZq0YuIY3VlmmKAESNGMHbsWGJiYrj00kuZOnUqxhjGjh3bPKXB5MmTWbJkCY888gjHjh0jLi6OdevWNe934sSJ7N69mz179nRzKTunFr2IdIuOWt5nVFfXNk9qFimdTVMMkJCQ0Pw5JiameTkmJob6+noA7rzzTqZNm8bLL7/MlVdeyWuvvYa1lrvuuosFCxZENM+Rpha9iDhaZ9MUd9V7773H2LFjueOOO0hLS2P37t1885vf5IknnqCmpgaADz/8sHla4t5ELXoRcbTOpinuqoceeogNGzYQGxvL6NGjufbaa0lISGDXrl186UtfAiApKYmnnnqKz3/+85HKfkSEFeiNMT8CbgUssB24BUgBVgHnAFuA71lr68LMp4hISDqbpnj48OHs2LGjeXn58uV+v/vtb3/rd/vFixezePHiyGU4CkLuujHGnA8sAtKstWOAWGA28Gvgf6y1I4GjwLxIZFREJFg7d+7koosuYurUqZqmOMzt3caY00AiUAl8Dcj0fb8CuBcoCPM4IiJB0zTFTUIO9NbaD40xucAHQC2wDigDjllr633JDgDn+9veGDMfmA8wdOjQVk9oCUZNTU3I2/ZVKnP/4IQyDxo0iBMnTmCM6VL6hoaGbnniUm/SlTJba/F6vSH/ewg50BtjzgZuBEYAx4BngGv95dHf9tbaQqAQIC0tzU6ZMiWkfBQXFxPqtn2Vytw/OKHMFRUV1NXVMWTIkC4F++rq6ogPr+ztOiuztZYjR44wePDg5mfZBiucrpuvAxXW2sMAxpi1wJeBwcaYOF+rfhjwURjHEJE+bNiwYRw4cIDDhw93Kb3X68XlckU5V71LV8rscrkYNmxYyMcIJ9B/AFxpjEmkqetmKlAKbABuomnkzRzg+TCOISJ92IABAxgxYkSX0xcXF4fcau2ruqPMIY+6sda+BayhaQjldt++CoE7gCXGmL3AEODxCORTRERCFNaoG2vtPUDb6drKgcvD2a+IiESOpkAQEXE4BXoREYdToBcRcTgFehERh1OgFxFxOAV6ERGHU6AXEXE4BXoREYdToBcRcTgFehERh1OgFxFxOAV6ERGHU6AXEXE4BXoREYdToBcRcTgFehERh1OgFxFxOAV6ERGHU6AXEXE4BXoREYdToBcRcTgFehERh1OgFxFxOAV6ERGHU6AXEXE4BXoREYdToBcRcTgFehERh1OgFxFxOAV6ERGHU6AXEXE4BXoREYcLK9AbYwYbY9YYY3YbY3YZY75kjDnHGPNXY8we3/vZkcqsiIgEL9wW/cPAK9baUcA4YBdwJ7DeWjsSWO9bFhGRHhJyoDfGeICrgMcBrLV11tpjwI3ACl+yFcC3w82kiIiEzlhrQ9vQmPFAIbCTptZ8GbAY+NBaO7hFuqPW2nbdN8aY+cB8gKFDh6auWrUqpHzU1NSQlJQU0rZ9lcrcP6jM/UM4ZU5PTy+z1qZ1mtBaG9ILSAPqgSt8yw8DOcCxNumOdrav1NRUG6oNGzaEvG1fpTL3Dypz/xBOmYFS24V4HU4f/QHggLX2Ld/yGmAiUGWMSQHwvR8K4xgiIhKmkAO9tfYgsN8Yc7Fv1VSaunFeAOb41s0Bng8rhyIiEpa4MLf/AbDSGBMPlAO30FR5/MkYMw/4APhOmMcQEZEwhBXorbXbaOqrb2tqOPsVEZHI0Z2xIiIOp0AvIuJwCvQiIg6nQC8i4nAK9CIiDqdALyLicAr0IiIOp0AvIuJwCvQiIg6nQC8i4nAK9CIiDqdALyLicAr0PpXVlVy9/GoO1hzs6ayIiESUAr1PzsYcSj4oYdnry3o6KyIiERXufPR9nvt+N956b/NyQWkBBaUFuOJc1N5d24M5ExGJjH7foi9fVE7mmEwS4xIBSIxLJGtsFhWLK3o4ZyIikdHvA33KwBQ8CR68DV5ccS68DV48CR6Sk5J7OmsiIhHR7wM9QNWnVWSnZrNp3iayU7N1QVZEHKXf99EDrJ21tvlz3rS8HsyJiEjkqUUvIuJwCvQiIg6nQC8i4nAK9CIiDqdALyLicAr0Ig6kuZukJQV6EQfS3E3SksbRiziI5m4Sf9SiF3EQzd0k/ijQiziI5m4SfxToRRxGczdJW+qjF3EYzd0kbalFLyLicGEHemNMrDFmqzHmRd/yCGPMW8aYPcaY1caY+PCzKSIioYpEi34xsKvF8q+B/7HWjgSOAvMicAwREQlRWIHeGDMMmAY85ls2wNeANb4kK4Bvh3MMEREJj7HWhr6xMWuAXwIDgaXAXGCTtfYi3/dfAP5irR3jZ9v5wHyAoUOHpq5atSqkPNTU1JCUlBTStn2Vytw/qMz9QzhlTk9PL7PWpnWWLuRRN8aY64FD1toyY8yUM6v9JPVbk1hrC4FCgLS0NDtlyhR/yTpVXFxMqNv2VSpz/6Ay9w/dUeZwhldOBqYbY64DXIAHeAgYbIyJs9bWA8OAj8LPpoiIhCrkPnpr7V3W2mHW2uHAbOBv1tosYANwky/ZHOD5sHMpIiIhi8Y4+juAJcaYvcAQ4PEoHENERLooInfGWmuLgWLf53Lg8kjsV0REwqc7Y0VEHE6BXkTE4RToRUQcToFeok7PLxXpWQr0EnV6fqlIz9J89BI1en6pSO+gFr1EjZ5fKtI7KNBL1Oj5pSK9gwK9RJWeXyrS89RHL1Gl55eK9Dy16EVEHE6BXkTE4RToRUQcToFeRMThFOhFJKI05UXvo0AvIhGlKS96Hw2vFJGI0JQXvZda9CISEZryovdSoBeRiNCUF72XAr2IRIymvOid1EcvIhGjKS96J7XoRUT8cNIwUQV6ERE/nDRMVF03IiItOHGYqFr0IiItOHGYqAK9iEgLThwmqkAvItKG04aJqo9eRKQNpw0TVYteRMThFOhFRBxOgV5ExOEU6EVEHC7kQG+M+YIxZoMxZpcx5l/GmMW+9ecYY/5qjNnjez87ctkVEZFghdOirwd+bK29BLgS+E9jzGjgTmC9tXYksN63LCIiPSTkQG+trbTWbvF9rgZ2AecDNwIrfMlWAN8ON5MiIhI6Y60NfyfGDAc2AmOAD6y1g1t8d9Ra2677xhgzH5gPMHTo0NRVq1aFdOyamhqSkpJC2ravUpn7B5W5fwinzOnp6WXW2rROE1prw3oBSUAZMMO3fKzN90c720dqaqoN1YYNG0Letq9SmfsHlbl/CKfMQKntQpwOa9SNMWYA8Cyw0lp75layKmNMiu/7FOBQOMcQEZHwhDPqxgCPA7ustQ+2+OoFYI7v8xzg+dCzJyIi4QpnrpvJwPeA7caYbb51PwV+BfzJGDMP+AD4TnhZFBGRcIQc6K21JYAJ8PXUUPcrIiKRpTtjRUQcToFeRMTh+vR89Pv3w/btHgYOBJcL3O6m95avGFVlItLP9elA/8c/wh13TOwwTXx8+0rAX4XQdl1X0gTar9sNCQmqZESkd+jTgX7WLDDmbS6+eBxeL82v2tqm16lTrdedeT916rPlI0daf9/yVV8fXv78VTLhVCRnlt9552zi4jquwFTJiMgZfTrQX3ABTJp0lClTorP/+vrWlUOgSqNtRdG20miZvuX6Tz4JXCl1XMmM6zTvZyqZcM9aAlUkgdIlJEBsbMT+BCISAX060EdbXBwkJTW9ult9feszj5avv/99K6NGTfBbkXS1Ujp6FCor/VdAp0+Hl/cBA0I7a+moItm7dwj19Z3vS5WMSHsK9L1UXFzT66yz2n9XU3M8amcxAA0NnZ/JdOWMJdAZzrFjgbvSAlcyY7uU9zOVTGdnHuGe3fhbp0pGeisFemknNrapgvFXyURbQ0P7M5naWnjzzTIuvTTVb2XS9kym7bYt93fsWOAKKtwzGX/XTTq6YN9ZxfLee5+jurprFZAqGemIAr30KrGxkJjY9Grp44+r+cpXontsf2cyHVUagc5Y/K2vrYXjx+HgwWDOZC7tct47uzgfzPWXYLdRJdP7KdCL+PS2M5mNGzdz2WWXd9pd5q9Sart8/DhUVflPV1cXXt4DncmE0v1VUTGUQ4e6vn2cIliX6GcS6QX8ncns33+S1NToH7tlJRPqmUxHo87OVDL+9t++krkkqLyfqWTCHUkWynDnvlTJ9KGsikg0BOou6w6Nja0rkeLiTYwff2WnZyxVx07why3PcOMXZxHXmBSwoqmuhsOH/XevnToVXt5jY8PvJnO7YeDAs6I6uAIU6EWkB8XENAU7t7tp+bzzvIwe3fl2C1+6k/2Nj3I6tYyHp+WHdOy2lUxH3WL+1retXNpWSmcqGX9nOS0rmSVLPCHlPxgK9CISksrqSmY/O5vVN60mOSm5W47pvt+Nt97bvFxQWkBBaQGuOBe1d9cGta+2lUy4Av0e/tafqWS8Xti8uQq4ODKZCEA3yotISHI25lDyQQnLXl/WpfSV1ZVcvfxqDtYcDPjdJ3WfdLhN+aJyMsdkkhjX1M+UGJdI1tgsKhZXdHqMUPIVTPpAv4e/9WcqGW9cJT/Z9YOg8hsKBXqRfiQSgdB9vxtzn6GgtIBG20hBaQHmPoP7/o6bxh1VDGe+W7FvRYfbpAxMwZPgwdvgxRXnwtvgxZPgaW4pB1v5tNqmeBmcBI4AB4C9wHbgH8BGYB3wPOQ8nkPJ+yUsy1sGD4L7Pv+/h7nX+P+dfuaGzwMDISczh+3HtgeV31CYpgeJ96y0tDRbWloa0rbFxcVMifaVjF5GZe4bgu3aaJs+GmVe+NJCHi17lAWpC8hv0bfdUV7bblNZXcnSdUsp2l3EyfqTJMYlknFJBrnX5PotZ9vuljNccS4Av98F4opxce3ga0mJSWH+wPkUflxIZV0lfzn1F7z4OUaji9p/1VLZUMns4bNZvW01yceToRbcN7nxxvnZ5rSL2vvbdwO573bjHdA+fcLpBGbumknRJUWcHHCSxNOJZOzP4Pa9t/PApQ9QlFzEydiTJDYkklGTQe6JXEYMG4HX+P9NgumCMsaUWWvTOkunFr1IFwXbGg7Uugz21D8SeeqsFe7v2IG2ufChC/Ec8+Ct9+IyLrz1Xjzve0helQy/A34D5ND0BOkfQflH5WTWZJLY6OtuaUgkqzKLir9UUP5cOZnlmSSebvrOfdpN1j+z2Fawjcx/ZpJY59umLpGsf2ZR8UAFaxetJe/7eYybM468H+ex9q61lOeWt05/OpGsd7OoeKYCdkHOoBxKBpWw7N+XwfnAWCg/UE7myUwSrW8bm0hWTBYV51ZAIfB74Bngz8BrUP61cjLPzyQx1pc+NpGsi7PY96N9eG7x4I33nWXEe/H8Hw/j3hyHZ7oHb5xvfZwXz1QPyY8lU/7jjrugIk0XY6VPCbeVHOp6aB0M8zsY6dHZBcO2+wmUPj4mnlNTfMMzGgAvUNv6PecfOZTsL2HZH5aRn5Lf+vsWn8tPlbM0bilF8UWcjDlJYmMiGUcyWDNkDeY+0z6v9S7Kf1vO0qlLKRpVxMn4kyTWJZKxO4PcdbksvG4h2TXZzC+bT2FqIZVJlfCnNj/EAMAFKe4UPOkevKO8uBpceGO8eBo9JHuSYSh4zvEFQ+sLhmM8jEsbh8d68J724sIXPNM9JGcngwtw+16+zymuFDxvePC+7cUV68JrvHgyPYwYPaL1b3teAQXnffa38LzowbvFt02DF89ED8nT/P+7SiEFz0kP3o9adBsN9JA8JJmqT6vITs1mfup8CssKqaypBAi4vmUXVHxMfLsuqEhToO/nemLkRDjH7mqw7Sx9MOsDBu5YF7WZtVR+Usnsv89m9ejVJDcmU35JOUvfX0pRdREn7UkSSSTDZLDmtP+gmtCYQObhTIrO9Z3i1yeS8X4Gv/jbLyAXqIXKhEpm3zSb1WtWk1yT3K4boeBQAQWHCvx3OxhfsL3Wg3eMF1e9C2+sF8+nHireqWDpxUspOsd37MZEMk5lkFufS/L3k/HEe/DaFsH2qx6Sv5vMWvfa5iCb585r+pzLZwE4gVbRpWp1FdlJbQLerPbf3fPnezjoOQizAmzzjcB/6ypv+6Bavqg8YDcTBA7EAY8RIP3aWWub0+RNy2v+HGh9y31NsBPYarZ2euxwqI++FwsUCJ999VkeqXwkIsE52D7bSFYMwRw7ISeBusb29+q7Yl3U/r/adi1Z92tuvLbrfb9YwLRf7ap3Uf7HcpZOWkrRF319sC1atsk1ySyctpBHUx9lQdkC8l9qKsdt026jMLWQ+IZ46mLrWFC2gJ+9/jOWXruUon/37ac+kYwPM8h9J5f7Rt9H4fmFxNt46kwdC04s4L/2/hfnXXgeuGGhayGP8igLEhaQPySfygGVLD22lKITvsokJpGMlAxyU3NJPju5VWuXAU1lm7F6BilJKa2C1NpZa7ntxdso3FJIfGw8dQ11rf4egbaJlmj8f+6ofL1BOGXu9330oYwu6CnB9tn+/v3fB+zLDbSvYEdOBDNUjAagBjgM7AfeBd4G3gKKgb8AzwFPA08EHqXgvscNN0POXTmU7Cth2Y+WweXAZbDz9zvJ3PtZX25zn+2vKyAZGEHT3fMTgclQ/ptyv328zX2/p1v05b6XxbZXtpG5P5PEhhb9yEezqNhfQUp6Cp7hLboX4r140jyMWDoCc6+hYFIBjTGNFEwqwNxrcC9zU3VjFdmjstk0axPZ47I5mHWQlBMpeP5vi77cAV48Mz0kv5lM1deryJ6UzabbNpE9KZuDVxzk3SXv4v6cG+MyFFBAI40UnCrAfGS48MML8VzkwYtvX9aLJ8VD8oRkGA4MBQYB8TRXYGtnrSVvWh7jkseRNy2vOWCfaVlumreJ7NTsVv92Am3Tl3RUvv7CsS36QK3FSIlmyzbQKIVAmq/UW1j4wkIe3fYoCy5ZQH5qfnMLd+G2hTz60aMsGLKA/M/lU1lTydLDSyk6VcRJTpJoE8moz2BN3BpOGT/3hgdq8bboKqhMat29EEhlUiVLr2nR93s6kYzyDNZ8cQ2n4tof29XoYl/JPu4dfy+FAwuJJ5466lgQs4D8pHy/fba44Lb3bqNwfyHxMfHUNdaxYMwC8r+Vz21/u43Cre1beMG2bPOuywtq9Emg/QQKnsXFxVycenHAYyx8aWG3tra7g5PP0APpjhZ9n+6jr6yuZPG2xbya9mrzf6xI3jnXkQ4jA9hfAAAHeElEQVT7ii1QR7sLZ+26F7a07l5o7rO1CWTWZVI0oMWFs6MZ3L7vdh447wGe+9xz1MbVNp36l2eQW5yL29t6qFjBrgIKdhW0C9AFRwooONLUlzt321y8qb4LZLFePNs9VGypYOnVSykaUcTJOF/3wuEMbj94Ow9c8ABFg315solknM4g1+bCfYAbchpyKDlVwrJly8j/t/z2wdf3OcWdgqfEg3dHiwtnWR4qrqoIGNR2pe+iqipwP68/VQeryE5rk/4sqDoZ3IUzCNzX2tGYbn866rP1p6Nx48HuS/qvPh3oczbmsP349lbB1u/Fl4syyL0yFz6AyqOVzH5jNqvHribZJrcbndD2vdJbyWzXbFZX+S6CXezGG+OnIql3UZvf1FdcGVfJ7Jmdt2zLk8pbt2zP9P1uyOW+KffhvazFhbPjHsYdGYcnxcOp2FO4Gn3rz/OQPDOZ8rpylsYupcj4WugkknFWBrcPvZ0Hjj5A0fEiTjb6+nKHZZB7RS4LRy4k2+MLbNsLqRxdScrzKa1HIhgvnhs8jJs2rv0IhS83jVBw3+/Ge7LFb3KigIIdHVeuVafbB9WOgtpudgcd2IK9QBZK4Az2Yl4ouuMY4mx9MtAHbLU3uKhdXovncg/eS30t1UYvnj94SJ7tu3NuWg4lqSUs+82y5gtnAcX50o8vYVnNMvK35lNeUc7ScUs/uwmi8bObILiZppZtYg4lMSUsu2sZ+Wfnf9aqbdmyPTMkbLsH7x5f8Iz34rnFQ/Kzye1HHFxS2TwS4YYTN3DfDfe1atWmkNI+EF/SJkCfCZ6f95A8Opm1o1sEtn/7LLAFO1Sss5EN/gQKqn0tqHVHq1otdwlXn+yjb74z71++EQcNiWQczyD341yS45OZkTyDFJPC/Jj5FA4opDK2kr80BLhzLsZF7fTadt0L7kcC3803d9xcv325Hd0BGKhlG8qohkB9eoH21R0jJ6I9skF9t/2DyhwcR/fRN5/i47vZwHjxXPPZjQ5radECoqkF1NFt2yS1P0ZHrdSFLy2Mess2FJHskghWX2uJi/QnfTLQQ/A3G3Q2GVIw6QMFzmCP4STqXhDpvfpsoD8TWIqLi7l1yq1d2iZSd8FFehsRkWjqs4E+FJEatRHpbUREoikqd8YaY75ljHnHGLPXGHNnNI4hIiJdE/FAb4yJBfKAa4HRwM3GmC48BVJERKIhGi36y4G91tpya20dsAq4MQrHERGRLoj4OHpjzE3At6y1t/qWvwdcYa39fpt084H5AEOHDk1dtWpVSMerqakhKcnP+EgHU5n7B5W5fwinzOnp6T02jt7P1Fe0q02stYU0PceFtLQ0G+oNA7rBon9QmfsHlTk6otF1cwD4QovlYcBHUTiOiIh0QTS6buJompF8KvAhTc9Qz7TW/quDbQ4D74d4yHOBj0Pctq9SmfsHlbl/CKfMF1hrP9dZooh33Vhr640x3wdeBWKBJzoK8r5tOs1oIMaY0q70UTmJytw/qMz9Q3eUOSo3TFlrXwZejsa+RUQkOI59lKCIiDRxQqAv7OkM9ACVuX9QmfuHqJe5V8xHLyIi0eOEFr2IiHSgTwf6/jZ5mjHmC8aYDcaYXcaYfxljFvd0nrqDMSbWGLPVGPNiT+elOxhjBhtj1hhjdvv+1l/q6TxFmzHmR75/0zuMMX80xrh6Ok+RZox5whhzyBizo8W6c4wxfzXG7PG9nx2NY/fZQN9PJ0+rB35srb0EuBL4z35QZoDFwK6ezkQ3ehh4xVo7ChiHw8tujDkfWASkWWvH0DQse3bP5ioqlgPfarPuTmC9tXYksN63HHF9NtDTDydPs9ZWWmu3+D5X0xQAzu/ZXEWXMWYYMA14rKfz0h2MMR7gKuBxAGttnbX2WM/mqlvEAW7fDZeJOPBuemvtRuCTNqtvBFb4Pq8Avh2NY/flQH8+sL/F8gEcHvRaMsYMByYAb/VsTqLuIeAnQGNPZ6SbXAgcBp70dVc9Zow5q6czFU3W2g+BXOADoBI4bq1d17O56jZDrbWV0NSQAz4fjYP05UDfpcnTnMgYkwQ8C/zQWnuip/MTLcaY64FD1tqyns5LN4oDJgIF1toJwKdE6XS+t/D1S98IjADOA84yxny3Z3PlLH050PfLydOMMQNoCvIrrbVrO0vfx00Gphtj9tHUNfc1Y8xTPZulqDsAHLDWnjlTW0NT4HeyrwMV1trD1trTwFrgyz2cp+5SZYxJAfC9H4rGQfpyoP8HMNIYM8IYE0/TxZsXejhPUWWMMTT13e6y1j7Y0/mJNmvtXdbaYdba4TT9ff9mrXV0S89aexDYb4y52LdqKrCzB7PUHT4ArjTGJPr+jU/F4RegW3gBmOP7PAd4PhoH6bMPBw9l8jQHmAx8D9hujNnmW/dT39xC4hw/AFb6GjDlwC09nJ+osta+ZYxZA2yhaWTZVhx4h6wx5o/AFOBcY8wB4B7gV8CfjDHzaKrwvhOVY+vOWBERZ+vLXTciItIFCvQiIg6nQC8i4nAK9CIiDqdALyLicAr0IiIOp0AvIuJwCvQiIg73/wF/LNBeZX8cjQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(x, y, 'magenta', label='y original')\n",
    "plt.plot(x, y_new, 'g*', label='y outlies')\n",
    "plt.plot(x, reg_mse.coef_[0]*x + reg_mse.intercept_, 'b', label='y mse')\n",
    "plt.grid()\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Попробуем определить параметры $ y = kx+b $ с помощью функции потерь Хьюбера"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Huber:  k= 0.515 b= 4.068\n",
      "MSE:    k= -0.478 b= 15.464\n",
      "Original: k= 0.5 b= 4\n"
     ]
    }
   ],
   "source": [
    "reg_huber = HuberRegressor()\n",
    "reg_huber.fit(x[:, np.newaxis], y_new)\n",
    "print('Huber: ', 'k=', round(reg_huber.coef_[0], 3), 'b=', round(reg_huber.intercept_, 3))\n",
    "print('MSE:   ', 'k=', round(reg_mse.coef_[0], 3), 'b=', round(reg_mse.intercept_, 3))\n",
    "print('Original:', 'k=', round(0.5, 3), 'b=', round(4, 3))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x7f2ea6f93dd8>"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD8CAYAAAB5Pm/hAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xt4VNWh9/HvykySmRDCtQKCh2DliBg0kKBWrIIoXrBaQKUVLXjwAHIsWOtbbX1bL2hbTzlWaYHTeMOjKFagiujxRomIr0gTwIqAgnhDAwgIJGQmyUzW+8dMYi6T21xy2fl9nmeezOzZl7US+K211+xZ21hrERER50pq6wKIiEhiKehFRBxOQS8i4nAKehERh1PQi4g4nIJeRMThFPQiIg6noBcRcTgFvYiIw7nbugAAvXv3tpmZmVFte+zYMbp06RLfArVzqnPnoDp3DrHUubCw8IC19jtNrdcugj4zM5OCgoKots3Pz2f06NHxLVA7pzp3Dqpz5xBLnY0xnzVnPQ3diIg4nIJeRMThFPQiIg7XLsboRcSZKioq2LNnD36/v1nrd+vWje3btye4VO1Lc+rs8XgYMGAAycnJUR1DQS8iCbNnzx66du1KZmYmxpgm1y8uLqZr166tULL2o6k6W2s5ePAge/bsYdCgQVEdo8mhG2PMY8aY/caYrTWW9TTGvG6M2Rn+2SO83BhjFhhjdhlj/mmMGRFVqURaqKi4iPOWnMfekr1tXRSpwe/306tXr2aFvERmjKFXr17NPiuKpDlj9EuAi+ssux1YY60dDKwJvwa4BBgcfswAFkddMpEWmLduHus/X889b97T1kWROhTysYv1d9jk0I21dp0xJrPO4iuA0eHnTwD5wG3h5f9jQ/cn3GCM6W6M6WetLYqplCIN8N7nxR/4tqezuGAxiwsW43F78N3ha8OSibQfpjn3jA0H/WprbVb49WFrbfca739jre1hjFkN/N5auz68fA1wm7W23rehjDEzCPX66dOnT86yZcuiqkBJSQnp6elRbdtRqc7fOlh2kMW7F7P+wHrKKstITUrl+72/z43fvZGeKT3boKTx44S/c7du3TjppJOavX4wGMTlciWwRLGZNGkSjz76KN27d29wnXvvvZdRo0YxZsyYZu2zZp3feustFixYwHPPPVdvvV27dnHkyJFay8aMGVNorc1t6hjx/jA20vlFxJbEWpsH5AHk5ubaaL8Zpm/SdQ6N1fmN1W+w9uu1eNweyoPlDP6XwUwcN7F1C5gATvg7b9++vUUfrrbXD2OttVhree2115pc9/7772/RvmvWOS0tDbfbHfF34PF4GD58eIv2XSXa6+j3GWP6AYR/7g8v3wOcUGO9AcBXUR5DpFn2HdvHrJxZbJi+gVk5s/SBrFT79a9/zUMPPVT9+o477mDBggX11nvggQfIysoiKyuLBx98EIBPP/2UU045hdmzZzNixAi++OILMjMzOXDgAADz5s1jyJAhXHjhhfz4xz9m/vz5AEybNo3ly5cDoeld7rzzTkaMGMGwYcPYsWMHABs3buTss89m+PDhXHDBBXz44YcJ/T1E26NfBUwFfh/++UKN5TcZY5YBZwJHND4vibZy8srq5wvHL2zDkkijbga2NL6KN+iFlozcZAMPNvz29OnTmThxInPnzqWyspJly5axcePGWusUFhby+OOP8+6772Kt5cwzz+S8886jR48efPjhhzz++OMsWrSo1jYFBQWsWLGCzZs3EwgEGDFiBDk5ORHL0Lt3bzZt2sSiRYuYP38+jzzyCEOGDGHdunW43W5WrVrFr371K1asWNGCirdMk0FvjHmG0AevvY0xe4A7CQX8X40x04HPgavCq78MXArsAkqB6xNQZhGRZsnMzKRXr15s3ryZffv2MXz4cHr16lVrnfXr1zNhwoTqGSQnTpzIW2+9xeWXX87AgQM566yz6u13/fr1XHHFFXi9XgB+8IMfNFiGiRNDw4g5OTmsXBnqlBw5coSpU6eyc+dOrLUEg8G41Lchzbnq5scNvDU2wroW+I9YCyUiDtRIz7uKr9gX9zH6G264gSVLlrB3717+7d/+rd77jV2Q0tD0wc25iKVKamoqAC6Xi0AgAISGlMaMGcPf/vY3tm7dymWXXdbs/UVDc92IiKNNmDCBV155hX/84x9cdNFF9d4/99xzef755yktLeXYsWP87W9/4/vf/36j+zznnHN48cUX8fv9lJSU8NJLL7WoTEeOHKF///4ALF26tEXbRkNTIIiIo6WkpDBmzBi6d+8e8dLNESNGMG3aNM444wwgdAYwfPhwPv300wb3OXLkSC6//HJOP/10Bg4cSG5uLt26dWt2mX7xi18wdepUHnjgAUaNGtXiOrVY1WVDbfnIycmx0Vq7dm3U23ZUqnPn4IQ6b9u2rUXrHz16NO5lCAaD9vTTT7cfffRRXPdbXFxsrbX22LFjNicnxxYWFka1n+bWOdLvEiiwzchYDd2IiGNt27aNk046ibFjxzJ48OC47nvGjBlkZ2czYsQIJk2axIgR7XdqLw3diIhjDR06lN27dydk308//XRC9psI6tGLiDicgl5ExOEU9CIiDqegFxFxOAW9iEgLLVmyhK+++na+xtGjR1NQEJqN/dJLL+Xw4cNtVbSIFPQi0q50hNtC1g36ml5++eVG56tvCwp6EWlX4nlbyFinKc7KyqpeZ/78+dx1110sX76cgoICpkyZQnZ2Nj5f7TuZ1ZzK+KmnnuKMM84gOzubmTNnEgwGCQaDTJs2jaysLIYNG8af//znmOvZFAV9K+oIPRWRtuK9z0vGAxksLlhMpa1kccFizN0G733eqPc5ffp0nnjiCYDqaYqnTJlSa52a0xRv2LCBhx9+mM2bNze4zyuvvJLc3FyWLl3Kli1bqmewrGv79u08++yzvP3222zZsgWXy1W9zZdffsnWrVt5//33ufbaa6OuX3Mp6FuRbmAt0rDdc3Zz1ZCrSHOnAZDmTmPKsCl8MveTqPdZc5ri1157rclpitPT06unKY7VmjVrKCwsZOTIkWRnZ7NmzRp2797NiSeeyO7du/npT3/KK6+8QkZGRszHaoq+GdsKdANrkab169qPrild8Qf9eNwe/EE/GakZ9E3vG9N+o52m2O12U1lZWf3a7/dHXK8h1lqmTp3K7373u3rvvffee7z66qssXLiQpUuX8uSTT7Zo3y2lHn0r2D1nN9dkXRPXnoqIE31d+nXcbwsZ7TTFffr0Yf/+/Rw8eJCysjJWr15dvU3Xrl0pLi5u9Lhjx45l+fLl7N8futPqoUOH+Oyzzzhw4ACVlZVMmjSJefPm8d5778Vcx6aoR98K+nXtR0ZqRtx7KiJOs/TypdU3HonXbSGjnaYY4De/+Q1nnnkmgwYNYsiQIdXbTJs2jVmzZuH1ennnnXciHnfo0KHce++9jBs3jsrKSpKTk1m4cCFer5frr7+++mzhzjvvjEs9G6OgbyVVN7CekTODvMI8ikp0K12R1lBZWcmGDRt47rnnGlznlltu4ZZbbqm3fM6cOcyZM6fe8kmTJjFp0qTq1/n5+dXPa85jP3nyZCZPnlxv+02bNlU/b+rMIB4U9K1EN7AWaX3btm3jsssuY8KECXGfprgjUdCLiGMlcprijkQfxoqIOJyCXkTE4RT0IiIOp6AXEXE4Bb2IiMMp6EVEHE5BLyKO1Zxpij/99FOGDBnCDTfcQFZWFlOmTOGNN95g1KhRDB48mI0bNwLw5ptvkp2dTXZ2NsOHD6/+otMf/vAHRo4cyWmnndYq33KNhq6jF5FWcfPNsGVL4+sEg14izFLQoOxsCE8fH9H06dOZOHEic+fOrZ6muCq4a9q1axfPPfcceXl5jBw5kqeffpr169ezatUqfvvb3/L8888zf/58Fi5cyKhRoygpKcHj8fDaa6+xc+dONm7ciLWWyy+/nHXr1nHuuec2vxKtQD16EXGs5kxTDDBo0CCGDRtGUlISp556KmPHjsUYw7Bhw6qnNBg1ahS33HILCxYs4PDhw7jdbl577bXq/Y4YMYIdO3awc+fOVq5l09SjF5FW0VjPu0pxsa96UrN4aWqaYoDU1NTq50lJSdWvk5KSCAQCANx+++2MHz+el19+mbPOOos33ngDay2//OUvmTlzZlzLHG/q0YuIozU1TXFzffzxxwwbNozbbruN3NxcduzYwUUXXcRjjz1GSUkJAF9++WX1tMTtiXr0IuJoTU1T3FwPPvgga9euxeVyMXToUC655BJSU1PZvn073/ve9wBIT0/nqaee4rjjjotX8eMipqA3xvwMuAGwwPvA9UA/YBnQE9gEXGetLY+xnCIiUWlqmuLMzEy2bt1a/XrJkiUR3/vTn/4Ucfu5c+cyd+7c+BU4AaIeujHG9AfmALnW2izABfwIuB/4o7V2MPANMD0eBRURaalt27Zx0kknMXbsWE1THOP2XmNMBZAGFAHnA9eE338CuAtYHONxRERaTNMUh0Qd9NbaL40x84HPAR/wGlAIHLbWBsKr7QH6R9reGDMDmAHQp0+fWndoaYmSkpKot+2oVOfOwQl17tatW4vuoBQMBlvljkvtSXPr7Pf7o/73EHXQG2N6AFcAg4DDwHPAJRFWjXiLdWttHpAHkJuba0ePHh1VOfLz84l2245Kde4cnFDn7du3t+hyyeLi4rhfXtneNbfOHo+n+l62LRXL5ZUXAJ9Ya7+21lYAK4Gzge7GmKoGZADwVQzHEBGRGMUS9J8DZxlj0owxBhgLbAPWAleG15kKvBBbEUVEJBZRB7219l1gOaFLKN8P7ysPuA24xRizC+gFPBqHcoqIJMSnn35KVlZWi7ZZsmQJN910U4JKFH8xXXVjrb0TqDtd227gjFj2KyLiZMFgMKYvb7WUpkAQEcdqzjTFEAref//3f+fUU09l3Lhx+Hw+AEaPHk1BQQEABw4cIDMzs3qbL774gosvvpiTTz6Zu+++u3r5U089xRlnnEF2djYzZ84kGAwCoW/N/uY3v+HMM8/knXfeSUR1G6QpEESkdTRjnmJvMEg85ylu7jTFO3fu5JlnnuHhhx/m6quvZsWKFVx77bWNHnrjxo1s3bqVtLQ0Ro4cyfjx4+nSpQvPPvssb7/9NsnJycyePZulS5fyk5/8hGPHjpGVlcU999zT/PrFiYJeRByr5jTF+/bta3Sa4uzsbABycnKqpyZuzIUXXli9r4kTJ7J+/XrcbjeFhYWMHDkSAJ/PVz3vjcvlYtKkSXGqWcso6EWkdTRjnmJfAq6jb+k0xS6Xq3roxu12U1lZCYS+sFRT6GLD2q+ttUydOpXf/e539Y7h8XhadVy+Jo3Ri4ijxTJNcWZmJoWFhQAsX7681nuvv/46hw4dwufz8fzzzzNq1CjGjh3L8uXLq6cqPnToEJ999ll8KhID9ehFxNFimab41ltv5eqrr+bJJ5/k/PPPr/XeOeecw3XXXceuXbu45ppryM3NBeDee+9l3LhxVFZWkpyczMKFCxk4cGDc6hMNBb2IOFpLpym+9dZbq58PGTKEf/7zn9Wv7733XgCmTZvGtGnTIu5v8uTJTJ48ud7yqpuTtAUN3YiIY2ma4hD16EXEsTRNcYh69CKSUNZGnMBWWiDW36GCXkQSxuPxcPDgQYV9DKy1HDx4EI/HE/U+NHQjIgkzYMAA9uzZw9dff92s9f1+f0yB1hE1p84ej4cBAwZEfQwFvYgkTHJyMoMGDWr2+vn5+VHfXKOjao06a+hGRMThFPQiIg6noBcRcTgFvYiIwynoRUQcTkEvIuJwCnoREYdT0IuIOJyCXkTE4RT0IiIOp6AXEXE4Bb2IiMMp6EVEHE5BLyLicAp6ERGHU9CLiDicgl5ExOEU9CIiDqegFxFxuJiC3hjT3Riz3Bizwxiz3RjzPWNMT2PM68aYneGfPeJVWBERablYe/QPAa9Ya4cApwPbgduBNdbawcCa8GsREWkjUQe9MSYDOBd4FMBaW26tPQxcATwRXu0J4IexFlJERKJnrLXRbWhMNpAHbCPUmy8E5gJfWmu711jvG2ttveEbY8wMYAZAnz59cpYtWxZVOUpKSkhPT49q245Kde4cVOfOIZY6jxkzptBam9vkitbaqB5ALhAAzgy/fgiYBxyus943Te0rJyfHRmvt2rVRb9tRqc6dg+rcOcRSZ6DANiOvYxmj3wPssda+G369HBgB7DPG9AMI/9wfwzFERCRGUQe9tXYv8IUx5uTworGEhnFWAVPDy6YCL8RUQhERiYk7xu1/Ciw1xqQAu4HrCTUefzXGTAc+B66K8RgiIhKDmILeWruF0Fh9XWNj2a+IiMSPvhkrIuJwCnoREYdT0IuIOJyCXkTE4RT0IiIOp6AXEXE4Bb2IiMMp6EVEHE5BLyLicAp6ERGHU9CLiDicgl5ExOEU9GFFxUWct+Q89pbsbeuiiIjElYI+bN66eaz/fD33vHlPWxdFRCSuYp2PvsPz3ufFH/BXv15csJjFBYvxuD347vC1YclEROKj0/fod8/ZzTVZ15DmTgMgzZ3GlGFT+GTuJ21cMhGR+Oj0Qd+vaz8yUjPwB/143B78QT8ZqRn0Te/b1kUTEYmLTh/0APuO7WNWziw2TN/ArJxZ+kBWRByl04/RA6ycvLL6+cLxC9uwJCIi8acevYiIwynoRUQcTkEvIuJwCnoREYdT0IuIOJyCXsSBNHeT1KSgF3Egzd0kNek6ehEH0dxNEol69CIOormbJBIFvYiDaO4miURBL+IwmrtJ6tIYvYjDaO4mqUs9ehERh4s56I0xLmPMZmPM6vDrQcaYd40xO40xzxpjUmIvpoiIRCsePfq5wPYar+8H/mitHQx8A0yPwzFERCRKMQW9MWYAMB54JPzaAOcDy8OrPAH8MJZjiIhIbIy1NvqNjVkO/A7oCtwKTAM2WGtPCr9/AvC/1tqsCNvOAGYA9OnTJ2fZsmVRlaGkpIT09PSotu2oVOfOQXXuHGKp85gxYwqttblNrRf1VTfGmMuA/dbaQmPM6KrFEVaN2JJYa/OAPIDc3Fw7evToSKs1KT8/n2i37ahU585Bde4cWqPOsVxeOQq43BhzKeABMoAHge7GGLe1NgAMAL6KvZgiIhKtqMforbW/tNYOsNZmAj8C/m6tnQKsBa4MrzYVeCHmUoqISNQScR39bcAtxphdQC/g0QQcQ0REmiku34y11uYD+eHnu4Ez4rFfERGJnb4ZKyLicAp6ERGHU9CLiDicgl4STvcvFWlbCnpJON2/VKRtaT56SRjdv1SkfVCPXhJG9y8VaR8U9JIwun+pSPugoJeE0v1LRdqexugloXT/UpG2px69iIjDKehFRBxOQS8i4nAKehERh1PQi0hcacqL9kdBLyJxpSkv2h9dXikicaEpL9ov9ehFJC405UX7paAXkbjQlBftl4JeROJGU160TxqjF5G40ZQX7ZN69CIiETjpMlEFvYhIBE66TFRDNyIiNTjxMlH16EVEanDiZaIKehGRGpx4maiCXkSkDqddJqoxehGROpx2mah69CIiDqegFxFxOAW9iIjDKehFRBwu6qA3xpxgjFlrjNlujPnAGDM3vLynMeZ1Y8zO8M8e8SuuiIi0VCw9+gDwc2vtKcBZwH8YY4YCtwNrrLWDgTXh1yIi0kaiDnprbZG1dlP4eTGwHegPXAE8EV7tCeCHsRZSRESiZ6y1se/EmExgHZAFfG6t7V7jvW+stfWGb4wxM4AZAH369MlZtmxZVMcuKSkhPT09qm07KtW5c1CdO4dY6jxmzJhCa21ukytaa2N6AOlAITAx/Ppwnfe/aWofOTk5Nlpr166NetuOSnXuHFTnziGWOgMFthk5HdNVN8aYZGAFsNRaW/VVsn3GmH7h9/sB+2M5hoiIxCaWq24M8Ciw3Vr7QI23VgFTw8+nAi9EXzwREYlVLHPdjAKuA943xmwJL/sV8Hvgr8aY6cDnwFWxFVFERGIRddBba9cDpoG3x0a7XxERiS99M1ZExOEU9CIiDteh56P/4gt4//0MunYFjwe83tDPmo8kNWUi0sl16KB/5hm47bYRja6TklK/EYjUINRd1px1Gtqv1wupqWpkRKR96NBBP3kyGPMeJ598On4/1Q+fL/QoK6u9rOpnWdm3rw8erP1+zUcgEFv5IjUysTQkVa8//LAHbnfjDZgaGRGp0qGDfuBAGDnyG0aPTsz+A4HajUNDjUbdhqJuo1Fz/ZrLDx1quFFqvJE5vcmyVzUysZ61NNSQNLReaiq4XHH7E4hIHHTooE80txvS00OP1hYI1D7zqPl4553NDBkyPGJD0txG6ZtvoKgocgNUURFb2ZOToztraawh2bWrF4FA0/tSIyNSn4K+nXK7Q48uXeq/V1JyJGFnMQDBYNNnMs05Y2noDOfw4YaH0hpuZIY1q+xVjUxTZx6xnt1EWqZGRtorBb3U43KFGphIjUyiBYP1z2R8Pnj77UJOPTUnYmNS90ym7rY193f4cMMNVKxnMpE+N2nsA/umGpaPP/4OxcXNa4DUyEhjFPTSrrhckJYWetR04EAx55yT2GNHOpNprNFo6Iwl0nKfD44cgb17W3Imc2qzy97Uh/Mt+fylpduokWn/FPQiYe3tTGbduo2cdtoZTQ6XRWqU6r4+cgT27Yu8Xnl5bGVv6EwmmuGvTz7pw/79zd/erQRrFv2aRNqBSGcyX3xRSk5O4o9ds5GJ9kymsavOqhqZSPuv38ic0qKyVzUysV5JFs3lzh2pkelARRWRRGhouKw1VFbWbkTy8zeQnX1Wk2cs+w4f5clNz3HFdyfjrkxvsKEpLoavv448vFZWFlvZXa7Yh8m8XujatUtCL64ABb2ItKGkpFDYeb2h18cf72fo0Ka3m/3S7XxR+Rcqcgp5aPyiqI5dt5FpbFgs0vK6jUvdRqmqkYl0llOzkbnlloyoyt8SCnoRiUpRcRE/WvEjnr3yWfqm922VY3rv8+IP+KtfLy5YzOKCxXjcHnx3+Fq0r7qNTKwa+n1EWl7VyPj9sHHjPuDk+BSiAfqivIhEZd66eaz/fD33vHlPs9YvKi7ivCXnsbdkb4PvHSo/1Og2u+fs5pqsa0hzh8aZ0txpTBk2hU/mftLkMaIpV0vWn7duHhs+eYv7//f/wv798NlnsGMHjz8+h+D6t3j2oRnw0kuwfDlJS5/E+2QeSXnzeHjNDS0qbzTUoxfpRKLphdfdJtpedc2GYVGd4Zaq97qXd2ciExvcpl/XfmSkZuAP+vG4PfiDfjJSM6rrUmv9i/9U/2vjdR9+P6/8v0WcvPMt1r1zFVef8AMo9kGJD4754Zgv9Cit2qYU39FtPFi2D35xCgQy2F/2Od4ALKiARRZCd1h9tLoOvwo/4MXw41vdgD6XEvF3Ek8mdCPxtpWbm2sLCgqi2jY/P5/Rif4ko51RnTuGloZq3fUTUefZL83mL4V/YWbOzFrB0lhZ625TVFzEra/dyvM7nqc0UEqaO40Jp0xg/rj5Eevpvc9LWbkfTwC8AfBWgCcA3W0K3gAk+cvxVoTe84Tf9zbwMz2YxCD6kRFwk1nZl/2ln2HKfZSVH8UTsLXWTamM5TeVAnhrPQr7bseXbPElg89N9c9yVxInHB3M1uM+5mhqgGBSMv9afDpjD17GC//yGut7F3A4pRzrSuXMitHM8d/GeQMu5nByOT43FKdCIPxdhJYOQRljCq21uU2tpx69SDO1NLgb6sE2tJ/GeryxlqlWL9zCY+8s5pm3FtPDprJ7xgc88uod+N57i6X7ZvHz4bPB5+MnT1+Fq6wCbwXcHABv/mLuu38xXQMuZqaczLjDpaQFk0itKGWwfZO+t10FZT4oDz8q/FDhozRgMcFIpWreBfwVSaFAtcZDenlvXDadUPi66MpQwIvPZSjov5WPeu6hODVAwOVmUPEpjNt3BV2TjuOpQSt4te86vl8yjhlH5oLXw4Eux/iv4xbxcvpaDiX7MW4vF6Zezn09/ou+XfuDJ3yY8M/jTRG37r6V5/c+T2mwlDRXGhNOmsD8C+dz9zt3k7dpJykuD+XBcmbmjGTa+Dv5YPVentu0oXr5d3NOpPf4MbxR/Gl1Yxmo01gmgoJeOpRYe8nRLofmB3FTQxu19nPpQrrd48X4yvAGYGAFvPnyYnJfWEzXSjejK14ODRkc88HR8KPEHxpeOObjg09fZPI377E7bwx9Xdnhyzp8oUeN0C2p6E5p4BCuQDmeQM0P58rgP0/i18CvAXgh/ID/iVC3oDEY66E4dRdZJp20il4cSCvhWMoxOJgMZFCdjMYLbi8m1cvqQWt4u+8mKlwuSlKCnFUyhmmHbwSvlweP/wvPdVtNIDmZI+5yJrmv5b4uD3Cz+5csso/hcqVSTjkze13PopMXfRvANULY64Gn37qRvPfyaoTtOVy7ZX6tv8VTvMpMXq3+Wxxe/RJbN5WHtykjNacnfcf3j/h37Uc/Mkoz8H9VY9ioawZ9e/Vl37F9zMqZxYycGeQV5lFUUgTQ4PKaQ1ApSSn1hqDiTUM3HVA869wWV05Ec+yqOjc09NCQhtaPuNxablo1i6fefZiZQ3/C/efcBT4f3/tzNkn+8nrDChmVbv582h84emg/K3cv5cqMy0gvS6K05CCbv/l/HPDtIaUiSFogif6BHpRWHMJTYesNR0R7RUQl3w4f+MNDCX63Ydj+kdQddsDlZW3muxT02Uq5y0VJSiUjis9m3NHLebr/33i7+6bq4YUz7Hn8h7mV3l0GcIf3Xh5zP0PQncKR5HKm95zFooGLaoVsvedeIJVa3ciJz06kX3q/WoG3cvLKeu/d+eKdJGUksXLyyka3iSTS+gsvXdjoMFM8jtHY+o2p2tdwO5zNZnNU+2ru0I2Cvh1rKAhXvLqCBUUL4hLOLR2zjWfD0OSxf7iUvq5u4Pdz0u/74yqvqDdum1GZzDPnPBbu6X77mP/RvSRXBBsd74300xXDf4fyJBcplRmAl/1pJexLP0qZy3As2dLbdwLfPXQqm47fzq4eeyhODRJwuRnoG8LYby7jlX5vsrb3hlCPNznIuRUXcfX+qfQaMADSvcxL+wNPuP/KlV2u4/d9HqIo7TADHlDaAAAKnklEQVS3lvyK54tfoNSWkpaUxoR+E5ifM5++PfrWDt9kwDQcUjeuvpG8TXmkuFLCPeFv/x7xDLbmSMT/58bq1x7EUudOH/Rt2VNtqYbK2lAQXpF3Bav3ro74D7a5AV13eKFK1SltvWNXVoLPx/954SaeLXyC60/+EXefdXtoWKHEB0fCj+IajxJfnasWQlc/rDyylJSKyloB6wlAWgUM9B+Hr/Ib3MEKkmP4MM2SSmlyJUdTA5QmW8pcBk+gB8cdG8gnPYr4rNvXlKQEqXC56Oc/kdOPnsnbfTbzXo+POJpSQdCVzNCKHK4un0r3lL4s6PkIK9NeJuBO5mhyBT9M+zFP2RUcdpeFe9KACf8Okzxc0veSUEBmzSBvZx5FZUWs/NFKbnwpcuhECtQ5feZwyduXNPh3mnb6tLgEWGuHeWMSEfTtqX6RKOiboaFfUktP8VsqkT3bhkK4iqmk1hUM3Ujhg+sLoNTHH9/4La9vX8UVfccxs/9PwoHrZ+UnS9lycANnu7O52H0ex0oPUnh0PfvLPyelIkh6IIn+5b04FjyAp8LW2r83AKkRP0xrrmRqntdXmGT2dDtAUfoRSpMrKXcl0atsADu776EkpbI6OKuubAi43Nz16X+yYtBLvNT171Qkuyl2B7g4ZQK3dbsLunghI/zo6oEMD6QlcePHN5L3RR4pSSmUV5YzM2smiy5exI1/v5G8zfUDsqU926aGBSJpSejk5+dzcs7JDR5j9kuz23WARcPJZ+gNaY2g79AfxhYVFzF3y1xezX21+j9WPL8515jqD9Ty72bRBX+sfW1uqQ+K/XC4xgdotXq5fijx8fvP7yS5IkhWBTwSAO9fF/P8zxezKmA4ofw7HAseJLWikrSAoUe5h67lbiorj5ESrNvVLYc/ngbAz8IPeDX8CJkYflQkbabYvRmfG0616RyXFsTvNpQmV2L9KQwtOZ9/fucj3uv2FcUpoeGFEyr+lZG+c1nbewMFXT/gSHIFle4UTjNn8hP3THqm9YcuHu5zPcCTlcuZ1ONa7hv0AHTzhkK3i6vWWG6yF/5z/Y3kba35wdl4fn3urxsMtc0FO1i97y36pd9YHWzvlhTB5IZvSLJv7z5m5db5IKwL7Ctt2QdnQK0AXTh+YfXzxq7pjqSh/TSksevGW7ov6bw6dNDPWzeP94+8/+1VENaye9YOfv3Srfx920tYv4+e1sMlfUfzs8Gz4Ok3OHTgSxZsuZ+b+8yie1lKKHhL/aFwLq0R1uGrFvzlR9heuYlTfSeSUh7g46SteAKWeytgfgA8gf8G/juq8t9GEmUuN8UpAfzhD9JSgt3oV3Iin/TYy/7ulZS5kziWXMlAXybnlFzAq33eYkP3LVS43RQnBzkjeC5TmM43aX4eSVvK2uQNfJNcBm4P52RcwPX9ZvOn0kdY6ftfDib5SE1OY8KACcw/cz6z/jGbfhmhHuHS97/tET5e3bOtCuHz+MH4RaxbfSN5m/5ZY3kWPcdPqde4/pYn+e3nTzbauO6rqB+qjYXaDna0ONgaWr+lyxvTWOMQL61xDHG2Dhn0VcEyowA+WwfeisUcCyzGE4B+Fh6ptbYfeCX8gJ7AXQDMjbTneo+i7vs40u0IO10HODWYy/Huk9hw/D/5sOtnFKeErtc90Z7KJYGJZKT2Aa+HR7os5cWUNzjfcxlze94G6V7o6oV0T2h4oZsXunsx6cn8bOts8nbWHC74MYsuW8Qv65zirygp4pzJC/jLsxOpPDqQu39wd/XyKZOvowewe/UmXt30VnUQD885gaHjL4HVqziwqYwUVzg8j8ug79C+rBxaI9j+5dtga+mlYrvn7G6wJ96QhkK1o4Vaa/Sq1XOXWHXIoK8Klv27VvDGiWVUuNwMqDyJc8svIj3lOzzeYxne5J6cnTyW1z3r2Z9ylLeTNnEkuaLWN9p8yUBKKkU/PAwZqZBmqocXvAvqjpN/AXxR50Owqp7t2Uwe/5t6PdtVvMjNvBjq2f68gZ7t9gjDCzTe68zPz+f0vqfX+0/f0oBuSEt7vE19Lb0lFGoi8dchg74qWJYNruDlf00hYAPMzBnDpeMfBOB67qhed3r4Z2Nf2ybdU+8YjfVSZ780O+E922jEc0iipTpaT1ykM+mQQQ/fBkvNLxs0pqW9zmg+BItnz7ajUU9cpP3qsEFfFSz5+fncMPqGZm3T0l5nNL1U9WxFpL3psEEfjXhdtRHvbUREEikhNx4xxlxsjPnQGLPLGHN7Io4hIiLNE/egN8a4gIXAJcBQ4MfGmGbcBVJERBIhET36M4Bd1trd1tpyYBlwRQKOIyIizRD3uW6MMVcCF1trbwi/vg4401p7U531ZgAzAPr06ZOzbNmyqI5XUlJCenp6bIXuYFTnzkF17hxiqfOYMWPabK4bE2FZvdbEWpsH5EFoUrNoJ/XRJEidg+rcOajOiZGIoZs9wAk1Xg8AvkrAcUREpBkSMXTjBj4CxgJfAv8ArrHWftDINl8Dn0V5yN7AgSi37ahU585Bde4cYqnzQGvtd5paKe5DN9bagDHmJkJz5LqAxxoL+fA2TRa0IcaYguaMUTmJ6tw5qM6dQ2vUOSFfmLLWvgy8nIh9i4hIyyTkC1MiItJ+OCHo89q6AG1Ade4cVOfOIeF1bhf3jBURkcRxQo9eREQa0aGDvrNNnmaMOcEYs9YYs90Y84ExJtL9EB3HGOMyxmw2xqxu67K0BmNMd2PMcmPMjvDf+nttXaZEM8b8LPxveqsx5hljTP27AXVwxpjHjDH7jTFbayzraYx53RizM/yzRyKO3WGDvpNOnhYAfm6tPQU4C/iPTlBnCN3gd3tbF6IVPQS8Yq0dApyOw+tujOkPzAFyrbVZhC7L/lHbliohlgAX11l2O7DGWjsYWBN+HXcdNujphJOnWWuLrLWbws+LCQVA/7YtVWIZYwYA46l7z3eHMsZkAOcCjwJYa8uttYfbtlStwg14w1+4TMOB36a31q4DDtVZfAXwRPj5E8APE3Hsjhz0/QndsbvKHhweejUZYzKB4cC7bVuShHsQ+AVQ2dYFaSUnAl8Dj4eHqx4xxnRp60IlkrX2S2A+8DlQBByx1r7WtqVqNX2stUUQ6sgBxyXiIB056Js1eZoTGWPSgRXAzdbao21dnkQxxlwG7LfWFrZ1WVqRGxgBLLbWDgeOkaDT+fYiPC59BTAIOB7oYoy5tm1L5SwdOeg75eRpxphkQiG/1Fq7sqn1O7hRwOXGmE8JDc2db4x5qm2LlHB7gD3W2qozteWEgt/JLgA+sdZ+ba2tAFYCZ7dxmVrLPmNMP4Dwz/2JOEhHDvp/AIONMYOMMSmEPrxZ1cZlSihjjCE0drvdWvtAW5cn0ay1v7TWDrDWZhL6+/7dWuvonp61di/whTHm5PCiscC2NixSa/gcOMsYkxb+Nz4Wh38AXcMqYGr4+VTghUQcpMPeHDyaydMcYBRwHfC+MWZLeNmvwnMLiXP8FFga7sDsBq5v4/IklLX2XWPMcmAToSvLNuPAb8gaY54BRgO9jTF7gDuB3wN/NcZMJ9TgXZWQY+ubsSIiztaRh25ERKQZFPQiIg6noBcRcTgFvYiIwynoRUQcTkEvIuJwCnoREYdT0IuIONz/B5fr+PFCVKruAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(x, y, 'magenta', label='y original')\n",
    "plt.plot(x, y_new, 'g*', label='y outlies')\n",
    "plt.plot(x, reg_mse.coef_[0]*x + reg_mse.intercept_, 'b', label='y mse')\n",
    "plt.plot(x, reg_huber.coef_[0]*x + reg_huber.intercept_, 'r', label='y huber')\n",
    "plt.grid()\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Реализация HuberRegressor в scikit-learn"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "https://github.com/scikit-learn/scikit-learn/blob/95d4f0841d57e8b5f6b2a570312e9d832e69debc/sklearn/linear_model/_huber.py#L17"
   ]
  },
  {
   "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.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}