ddribin/linear-fit.ipynb

## linear-fit.ipynb
{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Linear Fit using Python and NumPy\n",
    "\n",
    "This is a companion Jupyter notebook for my blog post, [Linear Fit using Python and NumPy](https://www.dribin.org/dave/blog/archives/2024/02/18/python-linear-fit/).\n",
    "\n",
    "## Create Some Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x=array([0, 1, 2, 3])\n",
      "y=array([-1. ,  0.2,  0.9,  2.1])\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "x = np.array([0, 1, 2, 3])\n",
    "y = np.array([-1, 0.2, 0.9, 2.1])\n",
    "print(f\"{x=}\")\n",
    "print(f\"{y=}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Using `Polynomial.fit`\n",
    "\n",
    "Here's how to use [`numpy.polynomial.Polynomial.fit()`][polynomial.fit] to get the coefficients of the linear equation. The third parameter to `fit()` is the degree and a `1` indicates it's a linear equation.\n",
    "\n",
    "[polynomial.fit]: https://numpy.org/doc/stable/reference/generated/numpy.polynomial.polynomial.Polynomial.fit.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$x \\mapsto \\text{0.55} + \\text{1.5}\\,\\left(\\text{-1.0} + \\text{0.66666667}x\\right)$"
      ],
      "text/plain": [
       "Polynomial([0.55, 1.5 ], domain=[0., 3.], window=[-1.,  1.], symbol='x')"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import numpy.polynomial.polynomial as poly\n",
    "\n",
    "linear = poly.Polynomial.fit(x, y, 1)\n",
    "linear"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$x \\mapsto \\text{-0.95} + \\text{1.0}\\,x$"
      ],
      "text/plain": [
       "Polynomial([-0.95,  1.  ], domain=[-1.,  1.], window=[-1.,  1.], symbol='x')"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "converted = linear.convert()\n",
    "converted"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[-0.95  1.  ]\n"
     ]
    }
   ],
   "source": [
    "coefs = linear.convert().coef\n",
    "print(coefs)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use destructuring to get the slope `m` and Y-intercept `b` from the coefficients:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "m=1.0 b=-0.95\n"
     ]
    }
   ],
   "source": [
    "b, m = coefs\n",
    "print(f\"{m=:.2} {b=:.2}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here's how to do this all in one line:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "m=1.0 b=-0.95\n"
     ]
    }
   ],
   "source": [
    "b, m = poly.Polynomial.fit(x, y, 1).convert().coef\n",
    "print(f\"{m=:.2} {b=:.2}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Using `polynomial.polyfit`\n",
    "\n",
    "Here's how to use [`numpy.polynomial.polynomial.polyfit()`][polynomial.polyfit] to get the coefficients of the linear equation. Again, the third parameter to `polyfit()` is the degree and a `1` indicates it's a linear equation.\n",
    "\n",
    "[polynomial.polyfit]: https://numpy.org/doc/stable/reference/generated/numpy.polynomial.polynomial.polyfit.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[-0.95  1.  ]\n"
     ]
    }
   ],
   "source": [
    "coefs = poly.polyfit(x, y, 1)\n",
    "print(coefs)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Again, use destructuring to get the slope `m` and Y-intercept `b` from the coefficients:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "m=1.0 b=-0.95\n"
     ]
    }
   ],
   "source": [
    "b, m = coefs\n",
    "print(f\"{m=:.2} {b=:.2}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And finally, here's how to do this all in one line:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "m=1.0 b=-0.95\n"
     ]
    }
   ],
   "source": [
    "b, m = poly.polyfit(x, y, 1)\n",
    "print(f\"{m=:.2} {b=:.2}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Plotting the Points and Line\n",
    "\n",
    "[Matplotlib](https://matplotlib.org) can plot the points and the line using only a few lines of code:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAi8AAAGdCAYAAADaPpOnAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAA9hAAAPYQGoP6dpAABXbUlEQVR4nO3deZyN5f/H8deZMYvBDGOZhcFIWZI9ogVlra9vql8LZcuSomwl2oQiRbTYSowsaWPShokQxm6KyBcRaYYszVhnjjn374+7OTVmcQ7nzJkz5/18POaR6z73fZ/PfBwzn67ruq/LYhiGgYiIiIiX8PN0ACIiIiLOUPEiIiIiXkXFi4iIiHgVFS8iIiLiVVS8iIiIiFdR8SIiIiJeRcWLiIiIeBUVLyIiIuJVink6AFez2Wz88ccflCpVCovF4ulwRERExAGGYXD69Gmio6Px88u/b6XIFS9//PEHMTExng5DRERErsDhw4epVKlSvucUueKlVKlSgPnNh4aGuvTeVquV5cuX07ZtWwICAlx676JGuXKccuU45cpxypVzlC/HuStXaWlpxMTE2H+P56fIFS9ZQ0WhoaFuKV5CQkIIDQ3Vh/sylCvHKVeOU64cp1w5R/lynLtz5ciUD03YFREREa+i4kVERES8iooXERER8SpFbs6LIwzD4OLFi2RmZjp1ndVqpVixYly4cMHpa32NcuU45cpxl8uVv78/xYoV0zIJIkWczxUvGRkZJCcnc+7cOaevNQyDyMhIDh8+rB+Ol6FcOU65cpwjuQoJCSEqKorAwMACjk5ECopPFS82m40DBw7g7+9PdHQ0gYGBTv2ysNlsnDlzhpIlS152AR1fp1w5TrlyXH65MgyDjIwM/vzzTw4cOMC1116rfIoUUT5VvGRkZGCz2YiJiSEkJMTp6202GxkZGQQHB+uH4mUoV45Trhx3uVwVL16cgIAAfvvtN/t5IlL0+ORPSv2CECm69O9bpOjTv3IRERFxSKbNYOOBk2w9bmHjgZNk2gyPxKHixUccPHgQi8VCUlKSw9fExcVRunRpj8cBULVqVSZPnuzSWERExHFLdyZzy/iVPDJrCx/u9eeRWVu4ZfxKlu5MLvBYVLxcoUybQeL+E3yRdITE/ScKpPo8fPgwjz76qH2ycZUqVRg4cCAnTpy47LUxMTEkJydTp04dh9/vwQcf5H//+9/VhOwx7ii8RER81dKdyTw+bxvJqReyHU9JvcDj87YVeAHjUxN2XWXpzhTGfL07219iVFgwIzvWpn2dKLe856+//kqzZs247rrr+Oijj4iNjeXnn3/mmWee4dtvv2XDhg2Eh4fnem1GRgaBgYFERkY69Z7FixenePHirghfRES8VKbNYNSXu8jtf9ENwAKM+nIXbWpH4u9XMMs9qOfFSSv2nKD/gu0FXn3279+fwMBAli9fTosWLahcuTIdOnTgu+++48iRIzz//PP2c6tWrcqYMWPo1q0boaGh9O3bN9fhmiVLlnDttdcSHBxMq1atmDNnDhaLhb/++gvI2Xvx8ssvU79+febOnUvVqlUJCwvjoYce4vTp0/Zzli5dyi233EJ4eDjVqlWjY8eO7N+/36nv9dixY3Ts2JHixYsTGxvL/Pnzc5zz5ptvcsMNN1CiRAliYmJ44oknOHPmDACrVq2iZ8+epKamYrFYsFgsvPzyywDMnTuXxo0bU6pUKSIjI+nSpQvHjh1zKj4REV+y6cDJHL/z/s0AklMvsOnAyQKLScWLEzJtBq9/92ue1SeY1aerh5BOnjzJsmXLeOKJJ3L0hERGRvLwww/z8ccfYxj/vO+ECROoV68e27dv58UXX8xxzwMHDvB///d/dOrUiR9//JHHHnssWwGUl/379xMfH89XX33FV199xerVq3nttdfsr589e5YhQ4awadMmvvjiC/z8/Ljnnnuw2WwOf789evTg8OHDfP/993z22WdMnTo1R4Hh5+fH22+/zc8//8ycOXNYuXIlw4YNA6B58+ZMnjyZ0NBQkpOTSU5O5umnnwbMFVrHjBnDjz/+SHx8PAcPHqRHjx4OxyYi4muOnc67cLmS81xBw0ZO2HzwJEdPZ+T5+r+rz2bXlHXZ++7duxfDMKhVq1aur9eqVYtTp07x559/UqFCBQBuv/12hg4daj/n4MGD2a6ZMWMGNWrU4I033gCgRo0a7Ny5k1dffTXfWGw2G3FxcZQqVQqArl27smLFCvt19913n/28ChUq8MEHHxAREcGuXbscmm/zv//9j2+//ZZNmzZx4403AvDBBx/k+N4HDRpk/3PVqlV55ZVX6NevH1OnTiUwMJCwsDAsFkuOobJHH33U/udq1arx9ttvc+ONN9oXPhMRkewqlHJsvSRHz3MF9bw44djpdAfPc0/1+e+elctp3Lhxvq/v2bPHXhxkadKkyWXvW7VqVXvhAhAVFZWtV2Tv3r107tyZ6tWrU7lyZapVqwbAoUOHHIp79+7dFCtWjEaNGtmP1axZM8fk2++++4477riDihUrUqpUKbp27cqJEycuu+3D1q1b6dixI5UrV6ZUqVK0aNHCqfhERHxNk9hwosKCyWs2iwVz3meT2NznXbqDihcnVCgV5OB5rq0+q1evjsViYffu3bm+vnv3bsqUKUP58uXtx0qUKOHSGLIEBARka1sslmxDQh07duTkyZPMmDGDhIQEEhMTAXPSsKscPHiQ//znP9StW5fPP/+crVu3MmXKlMu+z9mzZ2nXrh2hoaHMnz+fzZs3s3jxYpfHJyJSlPj7WRjZsTZAjgImqz2yY+0Cm6wLKl6ccmPVcCJKBRZ49Vm2bFnatGnD1KlTOX/+fLbXUlJSmD9/Pg8++KBT+zTVqFGDLVu2ZDu2efPmq4rzxIkT7NmzhxdeeIE77riDGjVqcOrUKafuUbNmTS5evMjWrVvtx/bs2WOfRAxm74nNZmPixIncdNNNXHfddfzxxx/Z7hMYGJhj1+FffvmFEydO8Nprr3HrrbdSs2ZNTdYVEXFA+zpRTHukIZFhwdRJ2UflU+bDKZFhwUx7pKHbnrTNi4oXJ/j7WRjW2hwGKejq89133yU9PZ127dqxZs0aDh8+zNKlS2nTpg0VK1a87FyVSz322GP88ssvPPvss/zvf//jk08+IS4uDuCKdzYuU6YMZcuW5b333mPfvn2sWbPGPlHWUTVq1KB9+/Y89thjbNy4ka1bt9K7d+9sE5WrV6+O1WrlnXfe4ddff2Xu3LlMnz49232qVq3KmTNnWLFiBcePH+fcuXNUrlyZwMBA+3VLlixhzJgxV/S9ioj4mva1I1hXbBtLFjzDgoTXmd+1Hmufvb3ACxdQ8eK0O2qUZUqXBkSGZR8acnf1ee2117JlyxaqVavGAw88wDXXXEPfvn1p1aoViYmJea7xkpfY2Fg+++wzFi1aRN26dZk2bZr9aaOgIMeGxy7l5+fHwoUL2bp1K3Xr1uW5555j/PjxTt9n9uzZREdH06JFC+6991769u1rn4gMUK9ePd58803Gjx9PnTp1mD9/PuPGjct2j+bNm9OvXz8efPBBypcvz+uvv0758uWJi4vj008/pXbt2rz22mtMmDDhir5XERGfkpwM7dvjN2wYflYrfhFlaRJZvECHirIx3Gjs2LFG48aNjZIlSxrly5c37r77buOXX3657HWffPKJUaNGDSMoKMioU6eO8fXXXzv8nqmpqQZgpKam5njt/Pnzxq5du4zz58879X1kyczMNE6dOmVkZmYaFzNtxvp9x4347b8b6/cdNy5m2q7onoXJK6+8YlSqVMkl9/p3riR/ypXjHMnV1f47LyoyMjKM+Ph4IyMjw9OheAXlKx9ffWUY5coZBhhG8eKGdepUI37xYpfnKr/f35dya8/L6tWr6d+/Pxs2bCAhIQGr1Urbtm05e/ZsntesX7+ezp0706tXL7Zv306nTp3o1KkTO3fudGeoTvP3s9DsmrLcXb8iza4p67nq8ypMnTqVzZs324de3njjDbp37+7psEREpDC4cAGeegr+8x84fhzq1oUtWzB694YrnF7gKm5d52Xp0qXZ2nFxcVSoUIGtW7dy22235XrNW2+9Rfv27XnmmWcAGDNmDAkJCbz77rs55jXI1dm7dy+vvPIKJ0+epHLlygwdOpQRI0Z4OiwREfG0XbvgoYdgxw6zPXAgvPYaBAeD1erZ2CjgRepSU1MB8p2fkZiYyJAhQ7Ida9euHfHx8bmen56eTnr6P+uvpKWlAeZKqtZLEmy1WjEMA5vN5tSKr1mMv9dZybqHt5s4cSITJ07McdwV31tRy5U7KVeOcyRXNpsNwzCwWq34+/sXZHiFStbPv0t/DkrulK+/GQZ+77+P39NPY7lwAaN8eTJnzsTo0MF8/V+/W12dK2fuV2DFi81mY9CgQdx88835rrSakpJCREREtmMRERGkpKTkev64ceMYNWpUjuPLly8nJCQk27FixYoRGRnJmTNnrmpdj3/v5SP5U64cp1w5Lr9cZWRkcP78edasWcPFixcLMKrCKSEhwdMheBVfzldAWhoNpkwhauNGAI7Vr8+2gQNJNwz45psc57s6V5dbZPTfCqx46d+/Pzt37mTt2rUuve+IESOy9dSkpaURExND27ZtCQ0NzXbuhQsXOHz4MCVLliQ42PmF5AzD4PTp05QqVeqKHyf2FcqV45QrxzmSqwsXLlC8eHFuu+22K/p3XlRYrVYSEhJo06ZNjsUlJSdfz5dl1Sr8+/fHcuQIRkAAtldfpcxTT3GHX86pse7KVdbIiSMKpHgZMGAAX331FWvWrKFSpUr5nhsZGcnRo0ezHTt69GiOPWqyBAUF5fpob0BAQI6kZmZmYrFY8PPzwy+Xv5DLyeqmzrqH5E25cpxy5ThHcuXn54fFYsn1Z4AvUh6c43P5slrh5Zdh3DgwDLjuOiwffYR/w4ZcbtDV1bly5l5u/UlpGAYDBgxg8eLFrFy5ktjY2Mte06xZM1asWJHtWEJCAs2aNXNXmCIiIr7n11/h1lth7FizcOnVC7ZuhYYNPR3ZZbm156V///4sWLCAL774glKlStnnrYSFhdlXTO3WrRsVK1a0LzI2cOBAWrRowcSJE7nrrrtYuHAhW7Zs4b333nNnqCIiIr5j3jx44gk4fRrCwuC99+CBBzwdlcPc2vMybdo0UlNTadmyJVFRUfavjz/+2H7OoUOHSE5OtrebN2/OggULeO+996hXrx6fffYZ8fHx+U7yFREREQekpUHXrubX6dNw883w449eVbhAAQwb5fbVo0cP+zmrVq2y76mT5f7772fPnj2kp6ezc+dO7rzzTneG6fVatmzJoEGDCuz94uLiKF26dJ6vHzx4EH9/f3b8vT7AqlWrsFgs2TZXFBGRArZxIzRoYPa6+PmZc11WrYIqVTwdmdM0O9BL9OjRA4vFkuNr3759LFq0KNsGg1WrVmXy5MnZrr9cweFOzZs3Jzk5mbCwMI+8f0FKTk6mS5cuXHfddfj5+TlcVB4+fJj//Oc/hISEUKFCBZ555plcH/Pt0aMHBw8edG3QLnLhwgX69+9P2bJlKVmyJPfdd1+OyfeXOnr0KD169CA6OpqQkBDat2/P3r17s53TsmVL++fd39+fMmXK8Pjjj7vzWxEpWjIzzQm5t9xiznOpXBnWrIGRI6FYgS735jIqXrxI+/btSU5OzvYVGxtLeHg4pUqV8nR4eQoMDCQyMtInHgNOT0+nfPnyvPDCC9SrV8+hazIzM3nwwQfJyMhg/fr1zJkzh7i4OF566SUATp48yZQpU+wLtAHs37+f+fPnuyzuqlWrsmrVqqu6x+DBg/nyyy/59NNPWb16NX/88Qf33ntvnucbhkGnTp349ddf+eKLL9i+fTtVqlShdevWObYQ6dOnD8nJyRw5coRffvnlijb8FPFJR45Amzbw3HNw8aI5PPTjj+ZwkRdT8eJFgoKCiIyMzPbl7++fbdioZcuW/PbbbwwePNj+f6urVq2iZ8+epKam2o+9/PLLgPnL9umnn6ZixYqUKFGCpk2b5vglFhcXR+XKlQkJCeGee+7hxIkTTsV96bBRVi/QsmXLqFWrFiVLlrQXZv82c+ZMatWqRXBwMDVr1mTq1KlOve/tt9/OgAEDsh37888/CQwMzPFEm6tUrVqVt956i27dujnc07R8+XL27NnD3LlzqV+/Ph06dGDMmDFMmTKFjIwMgoODOXLkCO3bt+f3339n+vTp9OjRg9jYWAzDoHXr1rRr185e3Jw8eZJKlSrZi5+CkJqaygcffMCbb77J7bffTqNGjZg9ezbr169nw4YNuV6zd+9eNmzYwLRp07jxxhupUaMG06ZN4/z583z00UfZzg0JCbF/5iMiInKs4SQiuYiPN/cj+v57KFECZs2ChQvBQ73wrqTixTDg7FnPfP3r/6RdZdGiRVSqVInRo0fbe2eaN2/O5MmTCQ0NtR97+umnAXMNnsTERBYuXMhPP/3E/fffn63rfuPGjfTq1YsBAwaQlJREq1ateOWVV646znPnzjFhwgTmzp3LmjVrOHTokD0mgPnz5/PSSy/x6quvsnv3bsaOHcuLL77InDlz7Oe0bNky2/ypS/Xu3ZsFCxZk2z5i3rx5VKxYkdtvvz3Xa3744QdKliyZ75crezwANmzYQO3atbOtLN2uXTvS0tL4+eefCQkJYezYsQwcOJBVq1axYcMGVq5cSfPmzbFYLMyZM4fNmzfz9ttvA9CvXz8qVqxYoMXL1q1bsVqttG7d2n6sZs2aVK5cmcTExFyvyfp7+fdCcn5+fgQFBeVYzHL+/PmUK1eOunXrMmrUKKdW4hTxOefOweOPwz33wMmT5qPP27ZBz54e31DRVbxzsMuVzp2DkiUdOtUPKO3K9z5zxqyGHfTVV19R8l+xdujQgU8//TTbOeHh4fj7+1OqVKlsC/uFhYVhsViyHTt06BCzZ8/m0KFDREdHA/D000+zdOlSZs+ezdixY+0bZQ4bNgyA6667jvXr1+fYdNNZVquV6dOnc8011wBmETV69Gj76yNHjmTixIn2YYfY2Fh27drFjBkz7DtfV65cmaioqDzf495772XAgAF88cUXPPD3TPq4uDj7/KHcNG7cmKSkpHxjv3T7iquVkpJChQoVcn2PlJQULly4wNixY9m4cSMtW7akcePGtG7dmjfeeIMmTZpQsWJFZsyYQbdu3UhJSeGbb75h+/btFCvAseyUlBQCAwNzzKvKb2uPrOJmxIgRzJgxgxIlSjBp0iR+//33bL1wXbp0oUqVKkRHR5OUlMTw4cM5ePAgixcvdue3JOKdfvoJOnc2N1YEePppePVVCAz0bFwupuLFi7Rq1Ypp06bZ2yWcKHxys2PHDjIzM7nuuuuyHU9PT6ds2bIA7N69m3vuuSfb682aNbvq4iUkJMReuABERUVx7NgxAM6ePcv+/fvp1asXffr0sZ9z8eLFbEMxH374Yb7vERwcTNeuXZk1axYPPPAA27ZtY+fOnSxZsiTPa4oXL0716tWv9Ntyi3PnzhEREcHSpUvp2bMn/fr1o0+fPiQmJtKkSRPAfEJv8eLFvPbaa0ybNo1rr70233v269ePefPmZXuPDh06ZNvI8MyZM7leO3bsWMaOHWtv78r6IemkgIAAFi1aRK9evexFd+vWrenQoUO2+T19+/a1//n6668nLCyMu+++m/3792f7DIn4NMOAd9+FZ56B9HSIiIAPP4S2bT0dmVuoeAkJMXtAHGCz2UhLSyM0NNQ1y7hfsnHk5ZQoUcKlv1jPnDmDv78/W7duzbH7bkkHe6Ou1KXLQFssFvsvrKxfmu+//z5NmzbNdp6zuwT37t2b+vXr8/vvvzN79mxuv/12quTzWOAPP/xAh6zdU/MwY8YMHn74YafiyE9kZGSOeSFZT+lERkYSHh5O//79s71+zTXXZPvFfe7cOfvf46VP6+Rm9OjR2YbpWrZsyfjx43PkOzf9+vWz92QBREdHExkZSUZGBn/99Ve23pf8tvYAaNSoEUlJSaSmppKRkUH58uVp2rQpjRs3zvcagH379ql4EQH4809zSOjrr832XXeZ81su6dEtSlS8WCyOD93YbOYjZyVKmM/IF1KBgYFkZmZe9liDBg3IzMzk2LFj3Hrrrbneq1atWmz8e4fRLHlNwHSViIgIoqOj+fXXX6+6SLjhhhto3Lgx77//PgsWLODdd9/N93xPDBvddNNNjB07lmPHjtl/0SckJBAaGkrt2rWznXvpmkhZhg4dip+fH99++y133nknd911V57zegAqVKiQbaiqWLFiVKxY0aHiODw8nPDw8GzHGjVqREBAACtWrOC+++4DYM+ePRw6dMihrT2yetT27t3Lli1bsj36f6ms9YPyGzIU8RkJCdCtG6SkQFAQvPEGDBhQZOa25EXFSxFUtWpV1qxZw0MPPURQUBDlypWjatWqnDlzhhUrVlCvXj1CQkK47rrrePjhh+nWrRsTJ06kQYMG/Pnnn6xYsYK6dety11138dRTT3HzzTczYcIE7r77bpYtW3bVQ0aOGDVqFE899RRhYWG0b9+e9PR0tmzZwqlTp+y7iF+6tUReevfuzYABAyhRokSOIbBLuWLYKKv4OXPmDH/++SdJSUkEBgbaC5HFixczYsQIfvnlFwDatm1LjRo16NatG2+88QYpKSm88MIL9O/fP9dNRy/19ddfM2vWLBITE2nYsCHPPPMM3bt356effqJMmTJX9b04KiwsjF69ejFkyBDCw8MJDQ3lySefpFmzZtx0003282rWrMm4cePsfw+ffvop5cuXp3LlyuzYsYOBAwfSqVMn2v7d1b1//34WLFjAnXfeSdmyZUlKSmLw4MHcdttt1K1bt0C+N5FCKSMDnn8eJkww27VqmU8S+cq/C6OISU1NNQAjNTU1x2vnz583du3aZZw/f/6K7p2ZmWmcOnXKyMzMvNownda9e3fj7rvvzvW1Fi1aGAMHDrS3ExMTjbp16xpBQUHGv/+K+/XrZ5QtW9YAjJEjRxqGYRgZGRnGSy+9ZFStWtUICAgwoqKijHvuucf46aef7Nd98MEHRqVKlYzixYsbHTt2NCZMmGCEhYXlGeuBAwcMwFizZo2RmZlpfP/99wZgnDp1yjAMw5g9e3aO6xcvXmxc+nGcP3++Ub9+fSMwMNAoU6aMcdtttxmLFi3K9n137949zziynD592ggJCTGeeOKJy57rCkCOrypVqthfnz17drbvNTMz0/jxxx+N9u3bG8WLFzfKlStnDB061LBarZd9r2PHjhkRERHG2LFj7ccyMjKMRo0aGQ888IDDMVepUsX4/vvvHT4/N+fPnzeeeOIJo0yZMkZISIhxzz33GMnJydnOAYzZs2fb22+99ZZRqVIlIyAgwKhcubLxwgsvGOnp6fbXDx06ZNx2221GeHi4ERQUZFSvXt148skn7Z+lvOK4mn/nRUVGRoYRHx9vZGRkeDoUr+BV+dqzxzAaNTIMc6aLYfTrZxhnzxbY27srV/n9/r6UxTDc8LyuB6WlpREWFkZqamqOtSAuXLjAgQMHiI2NzfZ4pqNcPuelCCtMuTp48CDXXHMNmzdvpmEh3C21MOWqsHMkV1f777yosFqtfPPNN9x555055phJTl6RL8OAuDh48klzuY3wcJg503wkugC5K1f5/f6+lIaNpMiyWq2cOHGCF154gZtuuqlQFi4iIg756y/o1w+yNjZu2RLmzoVKlTwZlcfof/OkyFq3bh1RUVFs3ryZ6dOnezocEZErs24d1K9vFi7+/ua6Ld9957OFC6jnRYqwli1bUsRGRUXEl1y8CGPHwqhR5tOusbHw0UfgwJIGRZ2KFxERkcLm0CF4+GHI2irj4Ydh6lTQvl6Aho1EREQKl88+g3r1zMKlZElzbsu8eSpc/sUne140lCBSdOnft3its2dh4ED44AOz3aQJLFgAWkk6B5/qecl6pEs70ooUXVn/vgvt464iudm2zdz9+YMPzNVxR4wwe15UuOTKp3pe/P39KV26tH0DwJCQkDx3F86NzWYjIyODCxcuaD2Oy1CuHKdcOS6/XBmGwblz5zh27BilS5d2eh8sEY+w2WDyZBg+HKxWiI42h4ny2d5DfKx4Aex7x2QVMM4wDIPz589TvHhxp4oeX6RcOU65cpwjuSpdunS+m0GKFBopKdCjByxbZrbvvtvseSlb1qNheQOfK14sFgtRUVFUqFABq9Xq1LVWq5U1a9Zw2223qUv6MpQrxylXjrtcrgICAtTjIt7h22/NwuXYMQgOhkmT4LHHivyGiq7ic8VLFn9/f6d/yPn7+3Px4kWCg4P1S+YylCvHKVeOU67E66Wnw7PPwltvme0bbjDXbrn+es/G5WV8tngREREpULt3Q+fO8OOPZvvJJ+H1182eF3GKZgeKiIi4k2HAe+9Bo0Zm4VKuHHz5Jbz9tgqXK6SeFxEREXc5eRL69IFFi8x269bw4YcQFeXZuLycel5ERETcYfVqc6XcRYugWDFziGjZMhUuLqCeFxEREVeyWmH0aHP3Z8OAa681V8pt3NjTkRUZKl5ERERc5cAB6NIFNmww2z17mnNbSpb0bFxFjIaNREREXGHBAqhf3yxcQkPNR6BnzVLh4gbqeREREbkap0/DgAHmRFyA5s1h/nyoWtWjYRVl6nkRERG5Ups3Q4MGZuHi5wcjR5oTdVW4uJV6XkRERJxls8Ebb8ALL8DFixATY/a23HqrpyPzCSpeREREnPHHH9C1K6xcabbvvx9mzIAyZTwblw/RsJGIiIijliyBunXNwiUkBGbOhI8/VuFSwNTzIiIicjnnz8PTT8PUqWa7QQPzaaIaNTwbl49Sz4uIiEh+du6EJk3+KVyGDoXERBUuHqSeFxERkdwYhlmwDB0K6ekQEQFz5kC7dp6OzOepeBEREbnU8ePw6KPm7s8AHTpAXBxUqODRsMSkYSMREZF/W7HCnJT75ZcQGAiTJ8PXX6twKUTU8yIiIgJYrFb8nnsOJk40h4xq1oSFC82doaVQUfEiIiKybx+3jhiB/759ZrtvX5g0yXwcWgodDRuJiIjvMgyYM4diTZpQZt8+jDJl4PPPzUXnVLgUWup5ERER35SaCo8/Dh99hAU4fv31hC1ZQkC1ap6OTC7DrT0va9asoWPHjkRHR2OxWIiPj8/3/FWrVmGxWHJ8paSkuDNMERHxNYmJUL++udCcvz+Zo0axbvRoc48iKfTcWrycPXuWevXqMWXKFKeu27NnD8nJyfavCprhLSIirpCZCa+8Ym6gePCgufvzDz9gGzEC/P09HZ04yK3DRh06dKBDhw5OX1ehQgVKly7t+oBERMR3HT4MjzwCa9aY7S5dzEXowsLAavVsbOKUQjlht379+kRFRdGmTRvWrVvn6XBERMTbLVpkPvK8Zg2ULGmulDtvnlm4iNcpVBN2o6KimD59Oo0bNyY9PZ2ZM2fSsmVLNm7cSMOGDXO9Jj09nfT0dHs7LS0NAKvVitXFlXTW/Vx936JIuXKccuU45cpxytXfzp7F/+mn8fvgAwBsjRuT+eGHUL06XLxoP035cpy7cuXM/SyGYRguffe83shiYfHixXTq1Mmp61q0aEHlypWZO3durq+//PLLjBo1KsfxBQsWEKLH3EREfFbor7/S+M03KfX77xgWC3vvvZdfHnoIIyDA06FJLs6dO0eXLl1ITU0lNDQ033MLVc9Lbpo0acLatWvzfH3EiBEMGTLE3k5LSyMmJoa2bdte9pt3ltVqJSEhgTZt2hCgD3++lCvHKVeOU64c59O5Mgz83nkHv+eew5KRgREVRebs2cTefjuxeVzi0/lykrtylTVy4ohCX7wkJSURFRWV5+tBQUEEBQXlOB4QEOC2D6A7713UKFeOU64cp1w5zudydfQo9OwJ335rtv/7XywffECxcuUcutzn8nUVXJ0rZ+7l1uLlzJkz7Mtaahk4cOAASUlJhIeHU7lyZUaMGMGRI0f48MMPAZg8eTKxsbFcf/31XLhwgZkzZ7Jy5UqWL1/uzjBFRKQoWLYMunc3C5jgYHOPoscfB4vF05GJi7m1eNmyZQutWrWyt7OGd7p3705cXBzJyckcOnTI/npGRgZDhw7lyJEjhISEULduXb777rts9xAREckmPR2eew7efNNs16ljLj5Xp45n4xK3cWvx0rJlS/KbDxwXF5etPWzYMIYNG+bOkEREpCjZswc6d4bt2812//7wxhtQvLhn4xK3KpTrvIiIiOTLMGDmTGjY0CxcypaFJUvg3XdVuPiAQj9hV0REJJtTp6BvX/jsM7N9xx3w4YcQHe3ZuKTAqOdFRES8xw8/mCvlfvYZFCsG48fD8uUqXHyMel5ERKTwu3gRxowxN1W02cwVchcsgBtv9HRk4gEqXkREpHA7eBAefhjWrzfb3bvDO+9AqVIeDUs8R8NGIiJSeH38MdSvbxYuoaFmb0tcnAoXH6eeFxERKXzOnIEnnzQLFYBmzWD+fIjNa4F/8SXqeRERkcJlyxbzEei4OPDzgxdfhDVrVLiInXpeRESkcLDZzCX9n38erFaoVMnsbbntNk9HJoWMihcREfG85GTo1g2++85s33svvP8+hId7Ni4plDRsJCIinvXVV1C3rlm4FC8O771nruOiwkXyoJ4XERHxjAsXYNgw87FnMJ8q+ugjqFnTo2FJ4aeeFxERKXg//wxNmvxTuAwaBBs2qHARh6jnRURECo5hwPTpMGSI2fNSoYL5VFGHDp6OTLyIihcRESkYJ05A794QH2+227WDOXMgIsKjYYn30bCRiIi43/ffm5Ny4+MhIADefBO++UaFi1wRFS8iIuI+Vis89xzccQf88QfUqAEbN8LgweYCdCJXQMNGIiLiHvv3Q5cusGmT2e7dGyZPhhIlPBqWeD+VvSIi4nrz5kGDBmbhUro0fPqpueicChdxAfW8iIiI66SlwRNPmMv6A9x6q1nIVK7s2bikSFHPi4iIuMaGDeZCc/Png78/jB5tTtRV4SIupp4XERG5OpmZMH48vPSS+ecqVWDBAmje3NORSRGl4kVERK7c779D166wapXZfughmDbNnOci4iYaNhIRkSsTHw/16pmFS4kSMHu22eOiwkXcTD0vIiLinHPnzOX9Z8ww240amRsqXnutZ+MSn6GeFxERcdxPP8GNN/5TuDzzDKxfr8JFCpR6XkRE5PIMw9wBetgwSE+HyEiYOxdat/Z0ZOKDVLyIiEj+/vwTevaEr7822//5D8yaBeXLezYu8VkaNhIRkbwtX25uqPj11xAUZPa+LFmiwkU8Sj0vIiKSU0YGPP88TJhgtmvXNifl1q3r2bhEUPEiIiKX+t//oHNn2LbNbD/+OEycCMWLezYukb9p2EhEREyGYa7V0rChWbiEh5truUydqsJFChX1vIiICPz1Fzz2GHzyidlu1cp8mqhiRY+GJZIb9byIiPi6devMlXI/+QSKFYNx4yAhQYWLFFrqeRER8VUXL8Krr5q7P9tsUK2aOSm3SRNPRyaSLxUvIiK+6Lff4JFHYO1as921K7z7LoSGejYuEQdo2EhExNd8+qk5TLR2LZQqBfPmwYcfqnARr6GeFxERX3H2LDz1lLk6LkDTpuYu0NWqeTYuESep50VExBds22Y+Aj1rFlgs5gJ0P/ygwkW8knpeRESKMpvNXGBuxAiwWs0niObNg5YtPR2ZyBVT8SIiUkQFnTqFf8eO5mPPAJ06wcyZULasR+MSuVoqXkREiiDLt9/SatAg/FJTzdVxJ02Cvn3NISMRL6fiRUSkKLlwAYYPp9hbb1EMMG64AcvChebGiiJFhIoXEZGiYvduc0PFH38EYP9//kPlBQsIKFXKw4GJuJZbnzZas2YNHTt2JDo6GovFQnx8/GWvWbVqFQ0bNiQoKIjq1asTFxfnzhBFRLyfYcB770GjRmbhUr48F+Pj2dm7NwQHezo6EZdza/Fy9uxZ6tWrx5QpUxw6/8CBA9x11120atWKpKQkBg0aRO/evVm2bJk7wxQR8V4nT8L//Z+5qeL589CmDfz4I8add3o6MhG3ceuwUYcOHejQoYPD50+fPp3Y2FgmTpwIQK1atVi7di2TJk2iXbt27gpTRMQ7rVplLvF/5AgEBJgbKg4eDH5+5mPRIkVUoZrzkpiYSOvWrbMda9euHYMGDcrzmvT0dNLT0+3ttLQ0AKxWK1YX/+PNup+r71sUKVeOU64cp1z9zWrFb8wY/MaPx2IYGNWrc3HePHMRusxMyMxUrpykfDnOXbly5n6FqnhJSUkhIiIi27GIiAjS0tI4f/48xYsXz3HNuHHjGDVqVI7jy5cvJyQkxC1xJmStmSCXpVw5TrlynC/nKiQlhUaTJhG+Zw8Av91xBzt69yYzJQW++SbH+b6cqyuhfDnO1bk6d+6cw+cWquLlSowYMYIhQ4bY22lpacTExNC2bVtCXbzJmNVqJSEhgTZt2hAQEODSexc1ypXjlCvH+XquLB99hP8zz2A5fRojLIzMqVOJvv9+onM519dz5Szly3HuylXWyIkjClXxEhkZydGjR7MdO3r0KKGhobn2ugAEBQURFBSU43hAQIDbPoDuvHdRo1w5TrlynM/l6vRp6N8f5s412zffjGX+fIpVqXLZS30uV1dJ+XKcq3PlzL0K1caMzZo1Y8WKFdmOJSQk0KxZMw9FJCLiYZs2QYMGZuHi5wcvv2xO1HWgcBEpqtxavJw5c4akpCSSkpIA81HopKQkDh06BJhDPt26dbOf369fP3799VeGDRvGL7/8wtSpU/nkk08YPHiwO8MUESl8bDZ47TW4+WbYvx8qV4bVq2HkSChWqDrNRQqcW/8FbNmyhVatWtnbWXNTunfvTlxcHMnJyfZCBiA2Npavv/6awYMH89Zbb1GpUiVmzpypx6RFxLccOQLdusHKlWb7/vthxgwoU8azcYkUEm4tXlq2bIlhGHm+ntvquS1btmT79u1ujEpEpBD74gvo1QtOnICQEHjnHejZ06kNFTNtBhsPnGTrcQtlD5ykWfUK+PtpQ0YpOtT3KCJSGJw/D08/DVOnmu2GDWHBAqhRw6nbLN2ZzKgvd5GcegHw58O9W4gKC2Zkx9q0rxPl+rhFPKBQTdgVEfFJO3bAjTf+U7gMHQrr119R4fL4vG1/Fy7/SEm9wOPztrF0Z7KrIhbxKBUvIiKeYhjw7rtm4fLzzxARAcuWwYQJkMsSEPnJtBmM+nIXuQ3UZx0b9eUuMm15D+WLeAsVLyIinnD8ONx9Nzz5JKSnw513wk8/Qdu2V3S7TQdO5uhx+TcDSE69wKYDJ68wYJHCQ8WLiEhB++47qFsXvvwSAgPhrbfgq6+gQoUrvuWx03kXLldynkhhpuJFRKSgZGTAs8+avSvJyVCrlrkI3VNPOfU0UW4qlAp26XkihZmKFxGRgrB3r7ng3Ouvm3NdHnsMtmyBevVccvsmseFEhQWTVwlkAaLCgmkSG+6S9xPxJBUvIiLuZBgwZ465xP+WLeZCc59/DtOnm+u4uIi/n4WRHWsD5ChgstojO9bWei9SJKh4ERFxl9RU6NIFevSAs2ehRQtzUu6997rl7drXiWLaIw2JDMs+NBQZFsy0RxpqnRcpMrRInYiIOyQmmoXLwYPg7w+jR5vzXfz93fq27etE0aZ2JIn7jrH8h420vbWpVtiVIkfFi4iIK2VmwtixMGqU+efYWHOl3JtuKrAQ/P0sNI0N58Rug6ax4SpcpMhR8SIi4iqHDsEjj8APP5jtLl3MVXPDwjwbl0gRozkvIiKu8Pnn5pNDP/wAJUvChx/C/PkqXETcQD0vIiJX4+xZGDwY3n/fbDdpYg4TXXONZ+MSKcLU8yIicqWSkqBxY7NwsVhg+HBYu1aFi4ibqedFRMRZNpu5pP/w4eaqudHRMHcu3H67pyMT8QkqXkREnHH0qLluy9KlZvu//4UPPoBy5Twalogv0bCRiIijli41N1RcuhSCg80nieLjVbiIFDD1vIiIXE56OowYAZMmme06deCjj8z/ikiBU/EiIpKfX36Bzp3NybkAAwaYmysWL+7RsER8mYoXEZHcGIY5l2XgQDh3DsqWhdmzoWNHT0cm4vNUvIiIXOrUKejbFz77zGy3bm3uDB0d7dm4RATQhF0Rkex++MFcKfezz6BYMXOIaNkyFS4ihYh6XkREAC5eNHd+fvVVcx2X6tXNSbmNG3s6MhG5hIoXEZGDB81NFBMTzXaPHvD221CqlCejEpE8aNhIRHzbwoXmMFFiIoSGmr0ts2ercBEpxNTzIiK+6fRpeOopiIsz282amRsqVq3qyahExAHqeRER37NlCzRsaBYufn7w4ouwZo0KFxEvoZ4XEfEdNhtMmADPP29O0I2JgXnz4LbbPB2ZiDhBxYuI+IbkZOjWDb77zmzfdx+8/z6UKePZuETEaRo2EpGi78svzQ0Vv/sOQkLMouXTT1W4iHgp9byISNF1/jwMGwbvvmu269c3nyaqWdOjYYnI1VHPi4gUTT//DE2a/FO4DB4MGzaocBEpAtTzIiJFi2HA9OkwZAhcuAAVKpj7ErVv7+nIRMRFVLyISNFx/Dj07g1ffGG227c3H4eOiPBoWCLiWho2EpGiYeVKc6XcL76AwECYNAm+/lqFi0gRpJ4XEfFuViu89BKMH28OGdWoYU7KbdDA05GJiJuoeBER77V/P3TuDJs3m+0+fcwelxIlPBuXiLiVho1ExDvNnWs++rx5M5Quba7b8t57KlxEfIB6XkTEqxQ7dw7/7t3NoSEwl/afN89c6l9EfIKKFxHxGpaNG2k5eDB+R4+Cvz+MHAnPPWf+WUR8hooXESn8MjPhtdfwHzmSEpmZGFWrYpk/H5o393RkIuIBKl5EpHD7/Xfo2hVWrcIC/H7rrUQsWkRAuXKejkxEPKRAJuxOmTKFqlWrEhwcTNOmTdm0aVOe58bFxWGxWLJ9BQcHF0SYIlLYLF5sbqi4ahWUKMHFmTPZOmQIhIV5OjIR8SC3Fy8ff/wxQ4YMYeTIkWzbto169erRrl07jh07luc1oaGhJCcn279+++03d4cpIoXJuXPQrx/cey+cOgWNG8P27RjduoHF4unoRMTD3F68vPnmm/Tp04eePXtSu3Ztpk+fTkhICLNmzcrzGovFQmRkpP0rQitkiviOH380i5UZM8z2sGGwbh1ce61n4xKRQsOtc14yMjLYunUrI0aMsB/z8/OjdevWJCYm5nndmTNnqFKlCjabjYYNGzJ27Fiuv/76XM9NT08nPT3d3k5LSwPAarVitVpd9J1gv+e//yt5U64cp1z9zTDwmzIFv+HDsWRkYERFkTlrFsYdd5iv/+vftM/nygHKlXOUL8e5K1fO3M9iGIbh0nf/lz/++IOKFSuyfv16mjVrZj8+bNgwVq9ezcaNG3Nck5iYyN69e6lbty6pqalMmDCBNWvW8PPPP1OpUqUc57/88suMGjUqx/EFCxYQEhLi2m9IRNwi8K+/aPDOO0Ru3QpA8o03kjRgABma2yLiM86dO0eXLl1ITU0lNDQ033MLXfFyKavVSq1atejcuTNjxozJ8XpuPS8xMTEcP378st+8s6xWKwkJCbRp04aAgACX3ruoUa4c5+u5siQk4P/oo1iOHsUICsL2+uvY+vXLdW6Lr+fKGcqVc5Qvx7krV2lpaZQrV86h4sWtw0blypXD39+fo0ePZjt+9OhRIiMjHbpHQEAADRo0YN++fbm+HhQURFBQUK7XuesD6M57FzXKleN8LlcZGeYCcxMnmu3rr8fy0Uf433ADl1tyzudydRWUK+coX45zda6cuZdbJ+wGBgbSqFEjVqxYYT9ms9lYsWJFtp6Y/GRmZrJjxw6ioqLcFaaIFLQ9e6BZs38KlyeeMPcouuEGz8YlIl7B7YvUDRkyhO7du9O4cWOaNGnC5MmTOXv2LD179gSgW7duVKxYkXHjxgEwevRobrrpJqpXr85ff/3FG2+8wW+//Ubv3r3dHaqIuJthwOzZ8OST5uPQ4eEwaxbcfbenIxMRL+L24uXBBx/kzz//5KWXXiIlJYX69euzdOlS++PPhw4dws/vnw6gU6dO0adPH1JSUihTpgyNGjVi/fr11K5d292hiog7/fUXPPYYfPKJ2b79dvjwQ6hY0aNhiYj3KZDtAQYMGMCAAQNyfW3VqlXZ2pMmTWLSpEkFEJWIFJi1a+Hhh+HQIShWDMaMgWee0YaKInJFtLeRiLjPxYvwyitmsWKzwTXXwIIF0KSJpyMTES+m4kVE3OO338zelnXrzHa3bvDuu1CqlGfjEhGvVyAbM4qIj/nkE6hXzyxcSpWC+fNhzhwVLiLiEup5ERHXOXMGBg40nyACaNrUHCaqVs2zcYlIkaKeFxFxjW3boFEjs3CxWOD55+GHH1S4iIjLqedFRK6OzQaTJsGIEWC1QqVKMG8etGjh6chEpIhS8SIiVy4lBbp3h+XLzfY998DMmebicyIibqJhIxG5Ml9/DXXrmoVL8eIwYwZ8/rkKFxFxO/W8iIhzLlyAZ5+Ft98223XrwkcfgVbBFpECop4XEXHcrl3mE0RZhcvAgbBxowoXESlQ6nkRkcszDHjvPRg8GM6fh/LlIS4O7rzT05GJiA9S8SIi+TtxAvr0gcWLzXbbtuaCc5GRno1LRHyWho1EJG+rVpkr5S5eDAEBMHEifPutChcR8SgVLyKSk9VqLjJ3++1w5Ahcdx1s2ABDhoCffmyIiGdp2EhEsvv1V+jSxZyIC9CrF0yeDCVLejQsEZEs+l8oEfnH/PlQv75ZuISFwccfm4vOqXARkUJEPS8iAmlpMGAAzJ1rtm++2SxkqlTxbFwiIrlQz4uIr9u0CRo0MAsXPz94+WVzoq4KFxEppNTzIuKrMjPh9dfhpZfg4kWoXBkWLDB7XURECjEVLyK+6MgR6NoVvv/ebD/wgLk3UenSHg1LRMQRGjYS8aBMm8HGAyfZetzCxgMnybQZ7n/TL74w9yP6/nsoUQJmzYKFC1W4iIjXUM+LiIcs3ZnMqC93kZx6AfDnw71biAoLZmTH2rSvE+X6Nzx/HoYOhWnTzHbDhuaGitdd5/r3EhFxI/W8iHjA0p3JPD5v29+Fyz9SUi/w+LxtLN2Z7No33LEDGjf+p3B5+mlITFThIiJeScWLSAHLtBmM+nIXuQ0QZR0b9eUu1wwhGQa8+y7ceKO5I3REBCxbBm+8AYGBV39/EREPUPEiUsA2HTiZo8fl3wwgOfUCmw6cvLo3+vNP+O9/4cknIT0d7roLfvrJ3FhRRMSLqXgRKWDHTudduFzJebn67jtzQ8WvvoKgIHj7bfjyS6hQ4crvKSJSSGjCrkgBq1Aq2KXnZZORAS+8YA4LAdSqZT5JVLeu8/cSESmk1PMiUsCaxIYTFRaMJY/XLUBUWDBNYsOdu/HevdC8+T+FS79+sGWLChcRKXJUvIgUMH8/CyM71gbIUcBktUd2rI2/X17lzSUMA+LizCX+t26F8HBYtMh8sigkxFVhi4gUGipeRDygfZ0opj3SkMiw7ENDkWHBTHukoePrvPz1F3TpAj17wtmz0LIl/Pgj3HOPy2MWESksNOdFxEPa14miTe1IEvcdY/kPG2l7a1OaVa/geI/L+vVm4fLbb+DvD6NHw7PPmn8WESnCVLyIeJC/n4WmseGc2G3QNDbcscIlMxNefdUsVjIzITbWXCm3aVP3BywiUgioeBHxJocOwSOPwA8/mO2HH4apUyE01LNxiYgUIM15EfEWn31mrt3yww9QsiTMnQvz5qlwERGfo54XkcLu7FkYNAhmzjTbTZrAggVwzTUeDUtExFPU8yJSmG3fDo0amYWLxQIjRsDatSpcRMSnqedFpDCy2eCtt2D4cHPV3Ohoc5jo9ts9HZmIiMepeBEpbI4ehR49YOlSs3333fDBB1C2rEfDEhEpLDRsJFKYLF1qLue/dCkEB5ur5C5erMJFRORf1PMiUgj4Wa34Pf20ufszwA03mGu3XH+9ZwMTESmEVLyIeNovv3DrsGH4Hzhgtp98El5/3ex5ERGRHDRsJOIphgHvv0+xpk0pfeAARrly8OWXZu+LChcRkTyp50XEE06ehL594fPPsQDH6tWjzJIlBFSu7OnIREQKvQLpeZkyZQpVq1YlODiYpk2bsmnTpnzP//TTT6lZsybBwcHccMMNfPPNNwURpkjBWLPGXCn388+hWDEyx40jceRIiHJwJ2kRER/n9uLl448/ZsiQIYwcOZJt27ZRr1492rVrx7Fjx3I9f/369XTu3JlevXqxfft2OnXqRKdOndi5c6e7QxVxr4sX4cUXoVUr+P13uPZaSEzENnQo+GkEV0TEUW7/ifnmm2/Sp08fevbsSe3atZk+fTohISHMmjUr1/Pfeust2rdvzzPPPEOtWrUYM2YMDRs25N1333V3qCLuc+AA3HYbvPKKuQBdz56wbRs0buzpyEREvI5b57xkZGSwdetWRowYYT/m5+dH69atSUxMzPWaxMREhgwZku1Yu3btiI+Pz/X89PR00tPT7e20tDQArFYrVqv1Kr+D7LLu5+r7FkXK1T8sCxfiP2AAlrQ0jNBQMqdMwXjwQfPFf31OlavLU64cp1w5R/lynLty5cz93Fq8HD9+nMzMTCIiIrIdj4iI4Jdffsn1mpSUlFzPT0lJyfX8cePGMWrUqBzHly9fTkhIyBVGnr+EhAS33Lco8uVcFTt/nhvee4/K338PwImaNdk6eDDnS5WCXOZx+XKunKVcOU65co7y5ThX5+rcuXMOn+v1TxuNGDEiW09NWloaMTExtG3bltDQUJe+l9VqJSEhgTZt2hAQEODSexc1vp4ry5Yt+HfrhmXfPgw/P2wjRhD6/PO0Kpbzn5yv58oZypXjlCvnKF+Oc1euskZOHOHW4qVcuXL4+/tz9OjRbMePHj1KZGRkrtdERkY6dX5QUBBBQUE5jgcEBLjtA+jOexc1Ppcrmw0mTIDnnzcn6MbEYJk/H/9bb8X/Mpf6XK6ugnLlOOXKOcqX41ydK2fu5dYJu4GBgTRq1IgVK1bYj9lsNlasWEGzZs1yvaZZs2bZzgezayqv80UKjT/+gLZt4dlnzcLl//4PfvwRbr3V05GJiBQpbh82GjJkCN27d6dx48Y0adKEyZMnc/bsWXr27AlAt27dqFixIuPGjQNg4MCBtGjRgokTJ3LXXXexcOFCtmzZwnvvvefuUEWu3Jdfmk8QnTgBISHmKrmPPgoWi6cjExEpctxevDz44IP8+eefvPTSS6SkpFC/fn2WLl1qn5R76NAh/P61xkXz5s1ZsGABL7zwAs899xzXXnst8fHx1KlTx92hijjv/Hl45hmYMsVs169vbqhYs6ZHwxIRKcoKZMLugAEDGDBgQK6vrVq1Ksex+++/n/vvv9/NUYlcpZ07oXNn878AQ4bA2LGQyxwsERFxHa9/2kikwBkGTJsGQ4fChQtQoQLMmQPt23s6MhERn6DiRcQZx49Dr16wZInZ7tABZs+GS9YmEhER99GGKiKOWrnS3FBxyRIIDITJk+Hrr1W4iIgUMPW8iFyO1WpuqPj66+aQUc2asHChWciIiEiBU/Eikp99+6BLF9i82Wz37QuTJpmPQ4uIiEdo2EgkN4YBH34IDRqYhUuZMvD55zBjhgoXEREPU8+LyKVSU+GJJ2DBArPdogXMnQsxMZ6NS0REAPW8iGS3YYPZ27JgAfj7wyuvwIoVKlxERAoR9byIAGRmwmuvwciR5p+rVjULGO2pJSJS6Kh4ETl8GLp2hdWrzXaXLjB1KoSFeTYuERHJlYaNxLctWmQ+8rx6NZQsaa6UO2+eChcRkUJMPS/im86dg8GDIWu38htvNIeJqlf3bFwiInJZ6nkR3/Pjj9CokVm4WCwwfDisXavCRUTES6jnRXyHYcDbb8OwYZCRAVFR5iPQd9zh6chERMQJKl7ENxw7Bj16wLffmu3//hc++ADKlfNoWCIi4jwNG0nRt2wZ1K1rFi7BwTBlCsTHq3AREfFS6nmRois9HZ57Dt5802zXqQMffWT+V0REvJaKFyma9uyBzp1h+3az3b8/vPEGFC/u2bhEROSqqXiRosUwYNYseOop83HosmVh9mzo2NHTkYmIiIuoeJGi49QpeOwx+PRTs33HHebO0NHRno1LRERcShN2pWhYuxbq1zcLl2LFYPx4WL5chYuISBGknhfxbhcvmjs/jxkDNpu50NyCBeaKuSIiUiSpeBHv9dtv8PDDsG6d2e7eHd55B0qV8mxcIiLiVho2Eu/08cfmhorr1kFoqNnbEhenwkVExAeo50W8y5kz5pNEs2eb7WbNYP58iI31bFwiIlJg1PMi3mPrVmjY0Cxc/PzgxRdhzRoVLiIiPkY9L1L42WzmKrnPPQdWK1SqZPa23HabpyMTEREPUPEihVtysjkRNyHBbN97L7z/PoSHezYuERHxGA0bSeH19dfmhooJCeay/u+9B599psJFRMTHqedFCp8LF2DYMPOxZzAXn/voI6hZ06NhiYhI4aCeFylcdu2CJk3+KVwGDYING1S4iIiInXpepHAwDJgxAwYPNnteKlQw123p0MHTkYmISCGj4kU878QJ6N0b4uPNdrt2MGcORER4NCwRESmcNGwknvX99+ZKufHxEBBgPhL9zTcqXEREJE8qXsQzrFZ4/nm44w44cgRq1ICNG81hIz99LEVEJG8aNpKC9+uv0KWLWayAOWQ0eTKUKOHRsERExDvof3GlYM2bZz76vHEjlC4Nn35qLjqnwkVERByknhcpGGlp0L+/WbwA3Hqr+efKlT0bl4iIeB31vIj7bdwIDRqYxYq/P4webU7UVeEiIiJXQD0v4j6ZmfiNHw+jRsHFi1ClCixYAM2bezoyERHxYipexD2OHKH5yy/jv2OH2X7oIZg2zZznIiIichU0bCSuFx9PsUaNKL9jB0aJEjB7ttnjosJFRERcQD0v4jrnzsHQoTB9Ohbgr2uuocSSJQTUru3pyEREpAhxa8/LyZMnefjhhwkNDaV06dL06tWLM2fO5HtNy5YtsVgs2b769evnzjDFFX76CW68EaZPByBzyBDWvPYaXHuthwMTEZGixq3Fy8MPP8zPP/9MQkICX331FWvWrKFv376Xva5Pnz4kJyfbv15//XV3hilXwzDMHaCbNDF3hI6MhIQEbK+9hhEQ4OnoRESkCHLbsNHu3btZunQpmzdvpnHjxgC888473HnnnUyYMIHo6Og8rw0JCSEyMtJdoYmr/Pkn9OwJX39ttv/zH5g1C8qXN5f/FxERcQO3FS+JiYmULl3aXrgAtG7dGj8/PzZu3Mg999yT57Xz589n3rx5REZG0rFjR1588UVCQkJyPTc9PZ309HR7Oy0tDQCr1YrVxb9As+7n6vt6I8t33+H/6KNYUlIwgoKwjR+P7fHHwWKBf+Veubo85cpxypXjlCvnKF+Oc1eunLmf24qXlJQUKlSokP3NihUjPDyclJSUPK/r0qULVapUITo6mp9++olnn32WPXv2sGjRolzPHzduHKNGjcpxfPny5XkWPFcrISHBLff1BharlVrz53NtfDwAaTExbB06lLSqVeHbb3Oc78u5cpZy5TjlynHKlXOUL8e5Olfnzp1z+Fyni5fhw4czfvz4fM/ZvXu3s7e1+/ecmBtuuIGoqCjuuOMO9u/fzzXXXJPj/BEjRjBkyBB7Oy0tjZiYGNq2bUtoaOgVx5Ebq9VKQkICbdq0IcAX53P873/4d+uG37ZtAGQ+9hjFX3+dW4oXz3Gqz+fKCcqV45QrxylXzlG+HOeuXGWNnDjC6eJl6NCh9OjRI99zqlWrRmRkJMeOHct2/OLFi5w8edKp+SxNmzYFYN++fbkWL0FBQQQFBeU4HhAQ4LYPoDvvXSgZBsTFwZNPwtmzEB4Os2bhf/fd+F/mUp/L1VVQrhynXDlOuXKO8uU4V+fKmXs5XbyUL1+e8uXLX/a8Zs2a8ddff7F161YaNWoEwMqVK7HZbPaCxBFJSUkAREVFORuquMJff0G/fvDxx2a7VSuYOxcqVvRoWCIi4rvc9qh0rVq1aN++PX369GHTpk2sW7eOAQMG8NBDD9mfNDpy5Ag1a9Zk06ZNAOzfv58xY8awdetWDh48yJIlS+jWrRu33XYbdevWdVeokpd166B+fbNwKVYMxo2DhAQVLiIi4lFuXWF3/vz5DBgwgDvuuAM/Pz/uu+8+3n77bfvrVquVPXv22CfpBAYG8t133zF58mTOnj1LTEwM9913Hy+88II7w5RLXbwIY8eaGyrabFCtGnz0kbmWi4iIiIe5tXgJDw9nwYIFeb5etWpVDMOwt2NiYli9erU7Q5LLOXQIHn4Y1q412127wrvvgosnP4uIiFwpbcwo//jsM6hXzyxcSpWCefPgww9VuIiISKGijRnFfIJo4ED44AOz3bSpuQt0tWqejUtERCQX6nnxddu2QcOGZuFiscDzz8MPP6hwERGRQks9L77KZoPJk2H4cHMfoooVzWGili09HZmIiEi+VLz4opQU6NEDli0z2506wcyZULasJ6MSERFxiIaNfM2335qTcpctg+LFYfp0WLRIhYuIiHgN9bz4ivR0ePZZeOsts123rrl2S+3ano1LRETESep58QW7d5tPEGUVLk89BRs3qnARERGvpJ6Xosww4P33YdAgOH8eypeH2bPhrrs8HZmIiMgVU/FSVJ08CX36mPNZANq0gTlzQBtcioiIl9OwUVG0erU5KXfRIggIgAkTYOlSFS4iIlIkqOelKLFaYfRoePVVc8jo2mvNSbmNGnk6MhEREZdR8VJUHDgAXbrAhg1m+9FHzQm6JUt6Ni4REREX07BRUbBgAdSvbxYuYWHw8cfmcv8qXEREpAhSz4s3O30aBgwwd34GuPlmmD8fqlTxbFwiIiJupJ4Xb7V5MzRoYBYufn7w8suwapUKFxERKfLU8+JtbDZ44w144QW4eBEqVzZ7W265xdORiYiIFAgVL97kjz+ga1dYudJs338/zJgBZcp4Ni4REZECpGEjb7Fkibkf0cqVEBJiTsj9+GMVLiIi4nPU81LYnT8PTz8NU6ea7YYNzaeLatTwbFwiIiIeop6XwmznTmjS5J/CZehQWL9ehYuIiPg09bwURoZhFixDh0J6OkREmE8VtW3r6chEREQ8TsVLYXP8uLk67pdfmu077zR3gq5QwbNxiYiIFBIaNipMVqwwJ+V++SUEBprL+3/1lQoXERGRf1HxUhhkZMDw4dCmDSQnQ61asGkTPPUUWCyejk5ERKRQ0bCRp+3bB507w5YtZvuxx+DNN83HoUVERCQH9bx4imHAnDnmEv9btpjrtXz+OUyfrsJFREQkH+p58YTUVHj8cfjoI7PdogXMmweVKnk2LhERES+gnpeClpgI9eubhYu/P7z6qjlRV4WLiIiIQ9TzUlAyM2HcOHP358xMiI01V8q96SZPRyYiIuJVVLwUhMOH4ZFHYM0as92li7kIXViYZ+MSERHxQho2crdFi6BePbNwKVnSXCl3/nwVLiIiIldIPS/ucvYsDB4M779vtps0MYeJrrnGs3GJiIh4OfW8uENSEjRubBYuFou5AN3atSpcREREXEA9L65kGOaS/s8+a66aGx0Nc+fC7bd7OjIREZEiQ8WLqxw9Cj17wrffmu3//hc++ADKlfNsXCIiIkWMho1cYdkyc1Lut99CcLD5JFF8vAoXERERN1DxcjXS02HoUGjf3ux5qVMHNm82V8/VhooiIiJuoWGjK7Vnj7mh4vbtZnvAAHj9dShe3LNxiYiIFHHqeXGWYWCZNQsaNjQLl7JlYckSeOcdFS4iIiIFQD0vzjh1isZvvEGx9evNduvW5s7Q0dGejUtERMSHqOfFUZs3U6xxYyquX49RrJg5RLRsmQoXERGRAqaeFwdllgnHOHmK1ApRHJ46mxvuaYu/nyblioiIFDS39by8+uqrNG/enJCQEEqXLu3QNYZh8NJLLxEVFUXx4sVp3bo1e/fudVeIDlu6M5lbPj3Ig51e4rYub9Fp80VuGb+SpTuTPR2aiIiIz3Fb8ZKRkcH999/P448/7vA1r7/+Om+//TbTp09n48aNlChRgnbt2nHhwgV3hXlZS3cm8/i8bSSnXmBrpdqcDQoBICX1Ao/P26YCRkREpIC5rXgZNWoUgwcP5oYbbnDofMMwmDx5Mi+88AJ33303devW5cMPP+SPP/4gPj7eXWHmK9NmMOrLXRi5vJZ1bNSXu8i05XaGiIiIuEOhmfNy4MABUlJSaN26tf1YWFgYTZs2JTExkYceeijX69LT00lPT7e309LSALBarVit1quKaeOBkySn5t3rYwDJqRdI3HeMprHhV/VeRU1W7q/278AXKFeOU64cp1w5R/lynLty5cz9Ck3xkpKSAkBERES24xEREfbXcjNu3DhGjRqV4/jy5csJCQm5qpi2HrcA/pc9b/kPGzmxW70vuUlISPB0CF5DuXKccuU45co5ypfjXJ2rc+fOOXyuU8XL8OHDGT9+fL7n7N69m5o1azpz26syYsQIhgwZYm+npaURExND27ZtCQ0Nvap7lz1wkg/3brnseW1vbaqel0tYrVYSEhJo06YNAQEBng6nUFOuHKdcOU65co7y5Th35Spr5MQRThUvQ4cOpUePHvmeU61aNWduaRcZGQnA0aNHiYqKsh8/evQo9evXz/O6oKAggoKCchwPCAi46qQ2q16BqLBgUlIv5DrvxQJEhgXTrHoFPTadB1f8PfgK5cpxypXjlCvnKF+Oc3WunLmXU8VL+fLlKV++vNMBOSI2NpbIyEhWrFhhL1bS0tLYuHGjU08suZK/n4WRHWvz+LxtWCBbAZNVqozsWFuFi4iISAFy29NGhw4dIikpiUOHDpGZmUlSUhJJSUmcOXPGfk7NmjVZvHgxABaLhUGDBvHKK6+wZMkSduzYQbdu3YiOjqZTp07uCvOy2teJYtojDYkMC852PDIsmGmPNKR9nag8rhQRERF3cNuE3Zdeeok5c+bY2w0aNADg+++/p2XLlgDs2bOH1NRU+znDhg3j7Nmz9O3bl7/++otbbrmFpUuXEhycvXAoaO3rRNGmdiSJ+46x/IeNtL21qYaKREREPMRtxUtcXBxxcXH5nmMY2WeSWCwWRo8ezejRo90V1hXz97PQNDacE7sNmsaGq3ARERHxEG3MKCIiIl5FxYuIiIh4FRUvIiIi4lVUvIiIiIhXUfEiIiIiXkXFi4iIiHgVFS8iIiLiVVS8iIiIiFdR8SIiIiJexW0r7HpK1qq9zmyt7Sir1cq5c+dIS0vTrqOXoVw5TrlynHLlOOXKOcqX49yVq6zf25euvp+bIle8nD59GoCYmBgPRyIiIiLOOn36NGFhYfmeYzEcKXG8iM1m448//qBUqVJYLK7dfygtLY2YmBgOHz5MaGioS+9d1ChXjlOuHKdcOU65co7y5Th35cowDE6fPk10dDR+fvnPailyPS9+fn5UqlTJre8RGhqqD7eDlCvHKVeOU64cp1w5R/lynDtydbkelyyasCsiIiJeRcWLiIiIeBUVL04ICgpi5MiRBAUFeTqUQk+5cpxy5TjlynHKlXOUL8cVhlwVuQm7IiIiUrSp50VERES8iooXERER8SoqXkRERMSrqHgRERERr6Li5RJTpkyhatWqBAcH07RpUzZt2pTv+Z9++ik1a9YkODiYG264gW+++aaAIvU8Z3IVFxeHxWLJ9hUcHFyA0XrOmjVr6NixI9HR0VgsFuLj4y97zapVq2jYsCFBQUFUr16duLg4t8dZGDibq1WrVuX4XFksFlJSUgomYA8aN24cN954I6VKlaJChQp06tSJPXv2XPY6X/yZdSW58tWfWdOmTaNu3br2BeiaNWvGt99+m+81nvhMqXj5l48//pghQ4YwcuRItm3bRr169WjXrh3Hjh3L9fz169fTuXNnevXqxfbt2+nUqROdOnVi586dBRx5wXM2V2CuxpicnGz/+u233wowYs85e/Ys9erVY8qUKQ6df+DAAe666y5atWpFUlISgwYNonfv3ixbtszNkXqes7nKsmfPnmyfrQoVKrgpwsJj9erV9O/fnw0bNpCQkIDVaqVt27acPXs2z2t89WfWleQKfPNnVqVKlXjttdfYunUrW7Zs4fbbb+fuu+/m559/zvV8j32mDLFr0qSJ0b9/f3s7MzPTiI6ONsaNG5fr+Q888IBx1113ZTvWtGlT47HHHnNrnIWBs7maPXu2ERYWVkDRFV6AsXjx4nzPGTZsmHH99ddnO/bggw8a7dq1c2NkhY8jufr+++8NwDh16lSBxFSYHTt2zACM1atX53mOL//M+jdHcqWfWf8oU6aMMXPmzFxf89RnSj0vf8vIyGDr1q20bt3afszPz4/WrVuTmJiY6zWJiYnZzgdo165dnucXFVeSK4AzZ85QpUoVYmJi8q3kfZ2vfq6uRv369YmKiqJNmzasW7fO0+F4RGpqKgDh4eF5nqPPlsmRXIF+ZmVmZrJw4ULOnj1Ls2bNcj3HU58pFS9/O378OJmZmURERGQ7HhERkef4eUpKilPnFxVXkqsaNWowa9YsvvjiC+bNm4fNZqN58+b8/vvvBRGyV8nrc5WWlsb58+c9FFXhFBUVxfTp0/n888/5/PPPiYmJoWXLlmzbts3ToRUom83GoEGDuPnmm6lTp06e5/nqz6x/czRXvvwza8eOHZQsWZKgoCD69evH4sWLqV27dq7neuozVeR2lZbCqVmzZtkq9+bNm1OrVi1mzJjBmDFjPBiZeLMaNWpQo0YNe7t58+bs37+fSZMmMXfuXA9GVrD69+/Pzp07Wbt2radDKfQczZUv/8yqUaMGSUlJpKam8tlnn9G9e3dWr16dZwHjCep5+Vu5cuXw9/fn6NGj2Y4fPXqUyMjIXK+JjIx06vyi4kpydamAgAAaNGjAvn373BGiV8vrcxUaGkrx4sU9FJX3aNKkiU99rgYMGMBXX33F999/T6VKlfI911d/ZmVxJleX8qWfWYGBgVSvXp1GjRoxbtw46tWrx1tvvZXruZ76TKl4+VtgYCCNGjVixYoV9mM2m40VK1bkOdbXrFmzbOcDJCQk5Hl+UXElubpUZmYmO3bsICoqyl1hei1f/Vy5SlJSkk98rgzDYMCAASxevJiVK1cSGxt72Wt89bN1Jbm6lC//zLLZbKSnp+f6msc+U26dDuxlFi5caAQFBRlxcXHGrl27jL59+xqlS5c2UlJSDMMwjK5duxrDhw+3n79u3TqjWLFixoQJE4zdu3cbI0eONAICAowdO3Z46lsoMM7matSoUcayZcuM/fv3G1u3bjUeeughIzg42Pj555899S0UmNOnTxvbt283tm/fbgDGm2++aWzfvt347bffDMMwjOHDhxtdu3a1n//rr78aISEhxjPPPGPs3r3bmDJliuHv728sXbrUU99CgXE2V5MmTTLi4+ONvXv3Gjt27DAGDhxo+Pn5Gd99952nvoUC8/jjjxthYWHGqlWrjOTkZPvXuXPn7OfoZ5bpSnLlqz+zhg8fbqxevdo4cOCA8dNPPxnDhw83LBaLsXz5csMwCs9nSsXLJd555x2jcuXKRmBgoNGkSRNjw4YN9tdatGhhdO/ePdv5n3zyiXHdddcZgYGBxvXXX298/fXXBRyx5ziTq0GDBtnPjYiIMO68805j27ZtHoi64GU9znvpV1Z+unfvbrRo0SLHNfXr1zcCAwONatWqGbNnzy7wuD3B2VyNHz/euOaaa4zg4GAjPDzcaNmypbFy5UrPBF/AcssTkO2zop9ZpivJla/+zHr00UeNKlWqGIGBgUb58uWNO+64w164GEbh+UxZDMMw3Nu3IyIiIuI6mvMiIiIiXkXFi4iIiHgVFS8iIiLiVVS8iIiIiFdR8SIiIiJeRcWLiIiIeBUVLyIiIuJVVLyIiIiIV1HxIiIiIl5FxYuIiIh4FRUvIiIi4lVUvIiIiIhX+X/QPz2R4UOKWwAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "\n",
    "plt.plot(x, y, \"o\", label=\"Original data\")\n",
    "plt.plot(x, m*x + b, \"red\",\n",
    "         label=f\"Fitted line: y = {m:.2}*x + {b:.2}\")\n",
    "plt.grid(True)\n",
    "plt.legend()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Using the \"Legacy\" `numpy.polyfit`\n",
    "\n",
    "For completeness, here's how to use the legacy [`numpy.polyfit()`][numpy.polyfit] API. Note that the coefficients are reversed from above, with slope being first and the Y-intercept being second:\n",
    "\n",
    "[numpy.polyfit]: https://numpy.org/doc/stable/reference/generated/numpy.polyfit.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 1.   -0.95]\n",
      "m=1.0 b=-0.95\n"
     ]
    }
   ],
   "source": [
    "coefs = np.polyfit(x, y, 1)\n",
    "print(coefs)\n",
    "m, b = coefs\n",
    "print(f\"{m=:.2} {b=:.2}\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "jupyter-notebook-3.12",
   "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.12.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}