{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "2f8e19e4",
   "metadata": {},
   "source": [
    "# Lineare Regression mit mehreren Features ($d>1$)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "643861b2",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "# plotting settings\n",
    "pd.plotting.register_matplotlib_converters()\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "import seaborn as sns\n",
    "from tqdm.notebook import tqdm"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "315bd31f",
   "metadata": {},
   "source": [
    "Wir verwenden hier beispielhaft den Datensatz [Melbourne Housing Snapshot](https://www.kaggle.com/datasets/dansbecker/melbourne-housing-snapshot). Diesen finden Sie auch im Moodle unter `data/kaggle/melb_data.csv`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "e3381ac0",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Index(['Suburb', 'Address', 'Rooms', 'Type', 'Price', 'Method', 'SellerG',\n",
       "       'Date', 'Distance', 'Postcode', 'Bedroom2', 'Bathroom', 'Car',\n",
       "       'Landsize', 'BuildingArea', 'YearBuilt', 'CouncilArea', 'Lattitude',\n",
       "       'Longtitude', 'Regionname', 'Propertycount'],\n",
       "      dtype='object')"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "melbourne_file_path = 'data/melb_data.csv'\n",
    "melbourne_data = pd.read_csv(melbourne_file_path)\n",
    "melbourne_data = melbourne_data.dropna(axis=0)  # entfernen von Daten mit fehlenden Werten\n",
    "melbourne_data.columns  # Spaltennamen der Tabelle (potentielle Features)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "0f80237c",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Suburb</th>\n",
       "      <th>Address</th>\n",
       "      <th>Rooms</th>\n",
       "      <th>Type</th>\n",
       "      <th>Price</th>\n",
       "      <th>Method</th>\n",
       "      <th>SellerG</th>\n",
       "      <th>Date</th>\n",
       "      <th>Distance</th>\n",
       "      <th>Postcode</th>\n",
       "      <th>...</th>\n",
       "      <th>Bathroom</th>\n",
       "      <th>Car</th>\n",
       "      <th>Landsize</th>\n",
       "      <th>BuildingArea</th>\n",
       "      <th>YearBuilt</th>\n",
       "      <th>CouncilArea</th>\n",
       "      <th>Lattitude</th>\n",
       "      <th>Longtitude</th>\n",
       "      <th>Regionname</th>\n",
       "      <th>Propertycount</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>Abbotsford</td>\n",
       "      <td>25 Bloomburg St</td>\n",
       "      <td>2</td>\n",
       "      <td>h</td>\n",
       "      <td>1035000.0</td>\n",
       "      <td>S</td>\n",
       "      <td>Biggin</td>\n",
       "      <td>4/02/2016</td>\n",
       "      <td>2.5</td>\n",
       "      <td>3067.0</td>\n",
       "      <td>...</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>156.0</td>\n",
       "      <td>79.0</td>\n",
       "      <td>1900.0</td>\n",
       "      <td>Yarra</td>\n",
       "      <td>-37.8079</td>\n",
       "      <td>144.9934</td>\n",
       "      <td>Northern Metropolitan</td>\n",
       "      <td>4019.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>Abbotsford</td>\n",
       "      <td>5 Charles St</td>\n",
       "      <td>3</td>\n",
       "      <td>h</td>\n",
       "      <td>1465000.0</td>\n",
       "      <td>SP</td>\n",
       "      <td>Biggin</td>\n",
       "      <td>4/03/2017</td>\n",
       "      <td>2.5</td>\n",
       "      <td>3067.0</td>\n",
       "      <td>...</td>\n",
       "      <td>2.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>134.0</td>\n",
       "      <td>150.0</td>\n",
       "      <td>1900.0</td>\n",
       "      <td>Yarra</td>\n",
       "      <td>-37.8093</td>\n",
       "      <td>144.9944</td>\n",
       "      <td>Northern Metropolitan</td>\n",
       "      <td>4019.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>Abbotsford</td>\n",
       "      <td>55a Park St</td>\n",
       "      <td>4</td>\n",
       "      <td>h</td>\n",
       "      <td>1600000.0</td>\n",
       "      <td>VB</td>\n",
       "      <td>Nelson</td>\n",
       "      <td>4/06/2016</td>\n",
       "      <td>2.5</td>\n",
       "      <td>3067.0</td>\n",
       "      <td>...</td>\n",
       "      <td>1.0</td>\n",
       "      <td>2.0</td>\n",
       "      <td>120.0</td>\n",
       "      <td>142.0</td>\n",
       "      <td>2014.0</td>\n",
       "      <td>Yarra</td>\n",
       "      <td>-37.8072</td>\n",
       "      <td>144.9941</td>\n",
       "      <td>Northern Metropolitan</td>\n",
       "      <td>4019.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>Abbotsford</td>\n",
       "      <td>124 Yarra St</td>\n",
       "      <td>3</td>\n",
       "      <td>h</td>\n",
       "      <td>1876000.0</td>\n",
       "      <td>S</td>\n",
       "      <td>Nelson</td>\n",
       "      <td>7/05/2016</td>\n",
       "      <td>2.5</td>\n",
       "      <td>3067.0</td>\n",
       "      <td>...</td>\n",
       "      <td>2.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>245.0</td>\n",
       "      <td>210.0</td>\n",
       "      <td>1910.0</td>\n",
       "      <td>Yarra</td>\n",
       "      <td>-37.8024</td>\n",
       "      <td>144.9993</td>\n",
       "      <td>Northern Metropolitan</td>\n",
       "      <td>4019.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>Abbotsford</td>\n",
       "      <td>98 Charles St</td>\n",
       "      <td>2</td>\n",
       "      <td>h</td>\n",
       "      <td>1636000.0</td>\n",
       "      <td>S</td>\n",
       "      <td>Nelson</td>\n",
       "      <td>8/10/2016</td>\n",
       "      <td>2.5</td>\n",
       "      <td>3067.0</td>\n",
       "      <td>...</td>\n",
       "      <td>1.0</td>\n",
       "      <td>2.0</td>\n",
       "      <td>256.0</td>\n",
       "      <td>107.0</td>\n",
       "      <td>1890.0</td>\n",
       "      <td>Yarra</td>\n",
       "      <td>-37.8060</td>\n",
       "      <td>144.9954</td>\n",
       "      <td>Northern Metropolitan</td>\n",
       "      <td>4019.0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>5 rows × 21 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "       Suburb          Address  Rooms Type      Price Method SellerG  \\\n",
       "1  Abbotsford  25 Bloomburg St      2    h  1035000.0      S  Biggin   \n",
       "2  Abbotsford     5 Charles St      3    h  1465000.0     SP  Biggin   \n",
       "4  Abbotsford      55a Park St      4    h  1600000.0     VB  Nelson   \n",
       "6  Abbotsford     124 Yarra St      3    h  1876000.0      S  Nelson   \n",
       "7  Abbotsford    98 Charles St      2    h  1636000.0      S  Nelson   \n",
       "\n",
       "        Date  Distance  Postcode  ...  Bathroom  Car  Landsize  BuildingArea  \\\n",
       "1  4/02/2016       2.5    3067.0  ...       1.0  0.0     156.0          79.0   \n",
       "2  4/03/2017       2.5    3067.0  ...       2.0  0.0     134.0         150.0   \n",
       "4  4/06/2016       2.5    3067.0  ...       1.0  2.0     120.0         142.0   \n",
       "6  7/05/2016       2.5    3067.0  ...       2.0  0.0     245.0         210.0   \n",
       "7  8/10/2016       2.5    3067.0  ...       1.0  2.0     256.0         107.0   \n",
       "\n",
       "   YearBuilt  CouncilArea Lattitude  Longtitude             Regionname  \\\n",
       "1     1900.0        Yarra  -37.8079    144.9934  Northern Metropolitan   \n",
       "2     1900.0        Yarra  -37.8093    144.9944  Northern Metropolitan   \n",
       "4     2014.0        Yarra  -37.8072    144.9941  Northern Metropolitan   \n",
       "6     1910.0        Yarra  -37.8024    144.9993  Northern Metropolitan   \n",
       "7     1890.0        Yarra  -37.8060    144.9954  Northern Metropolitan   \n",
       "\n",
       "  Propertycount  \n",
       "1        4019.0  \n",
       "2        4019.0  \n",
       "4        4019.0  \n",
       "6        4019.0  \n",
       "7        4019.0  \n",
       "\n",
       "[5 rows x 21 columns]"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "melbourne_data.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "b4939e52",
   "metadata": {},
   "outputs": [],
   "source": [
    "#features = ['BuildingArea', Rooms', 'Bathroom', 'Landsize', 'Lattitude', 'Longtitude', 'YearBuilt', 'Distance']\n",
    "features = ['Rooms', 'BuildingArea']\n",
    "data = melbourne_data[features + ['Price']]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "47f35849",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Rooms</th>\n",
       "      <th>BuildingArea</th>\n",
       "      <th>Price</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>count</th>\n",
       "      <td>6196.000000</td>\n",
       "      <td>6196.000000</td>\n",
       "      <td>6.196000e+03</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>mean</th>\n",
       "      <td>2.931407</td>\n",
       "      <td>141.568645</td>\n",
       "      <td>1.068828e+06</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>std</th>\n",
       "      <td>0.971079</td>\n",
       "      <td>90.834824</td>\n",
       "      <td>6.751564e+05</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>min</th>\n",
       "      <td>1.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>1.310000e+05</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25%</th>\n",
       "      <td>2.000000</td>\n",
       "      <td>91.000000</td>\n",
       "      <td>6.200000e+05</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50%</th>\n",
       "      <td>3.000000</td>\n",
       "      <td>124.000000</td>\n",
       "      <td>8.800000e+05</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>75%</th>\n",
       "      <td>4.000000</td>\n",
       "      <td>170.000000</td>\n",
       "      <td>1.325000e+06</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>max</th>\n",
       "      <td>8.000000</td>\n",
       "      <td>3112.000000</td>\n",
       "      <td>9.000000e+06</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "             Rooms  BuildingArea         Price\n",
       "count  6196.000000   6196.000000  6.196000e+03\n",
       "mean      2.931407    141.568645  1.068828e+06\n",
       "std       0.971079     90.834824  6.751564e+05\n",
       "min       1.000000      0.000000  1.310000e+05\n",
       "25%       2.000000     91.000000  6.200000e+05\n",
       "50%       3.000000    124.000000  8.800000e+05\n",
       "75%       4.000000    170.000000  1.325000e+06\n",
       "max       8.000000   3112.000000  9.000000e+06"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data.describe()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "ed0fdea0",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Rooms</th>\n",
       "      <th>BuildingArea</th>\n",
       "      <th>Price</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2</td>\n",
       "      <td>79.0</td>\n",
       "      <td>1035000.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3</td>\n",
       "      <td>150.0</td>\n",
       "      <td>1465000.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>4</td>\n",
       "      <td>142.0</td>\n",
       "      <td>1600000.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>3</td>\n",
       "      <td>210.0</td>\n",
       "      <td>1876000.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>2</td>\n",
       "      <td>107.0</td>\n",
       "      <td>1636000.0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   Rooms  BuildingArea      Price\n",
       "1      2          79.0  1035000.0\n",
       "2      3         150.0  1465000.0\n",
       "4      4         142.0  1600000.0\n",
       "6      3         210.0  1876000.0\n",
       "7      2         107.0  1636000.0"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b5126919",
   "metadata": {},
   "source": [
    "## Definition der Funktionen für die Lineare Regression"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "dacabf66",
   "metadata": {},
   "source": [
    "Aus der Vorlesung:\n",
    "\n",
    "$$ h(x, w) = w^T x . $$\n",
    "\n",
    "In der Vorlesung haben wir $\\theta$ statt $w$ verwendet.\n",
    "\n",
    "**Wichtig:** Diese Definition von $h$ nimmt an, dass die erste Komponente von $x$, also in Python code `x[0]`, immer 1 ist.\n",
    "\n",
    "Wir können auch eine vektorisierte Form von $h(x, w)$ definieren, bei der der Input $X$ mehrere (oder alle) Trainingsbeispiele umfasst und der Output ein Vektor aus den zugehörigen Werten von h zu jedem der Trainingsbeispiele ist. In Matrixschreibweise:\n",
    "\n",
    "$$ h(X, w) = X w , $$\n",
    "\n",
    "wobei die Zeilen von $X$ aus je einem Trainingsbeispiel (inkl. der \"1\" in der ersten Komponente) bestehen.\n",
    "\n",
    "Aufgrund der Art wie `numpy` den Spezialfall der Multiplikation zweier Vektoren handhabt können wir den Code für beide oben erwähnten Varianten von $h$ vereinheitlichen und eine Funktion $h(x, w)$ definieren, die sowohl mit einer Inputzeile als auch mit mehreren Inputzeilen umgehen kann.\n",
    "\n",
    "Bei der Multiplikation zweier numpy arrays (also zweier Vektoren) mittels `@`-Operator bildet numpy stets das Skalarprodukt der Vektoren, ohne dass man einen der Vektoren transponieren müsste. D.h., wenn wir zwei Spaltenvektoren $w, x$ haben, lautet die korrekte Schreibweise eigentlich:\n",
    "$$w^T x$$\n",
    "numpy erlaubt es uns aber einfach `w @ x` oder auch `x @ w` zu schreiben anstelle (des ebenfalls möglichen) `w.T @ x`.\n",
    "\n",
    "Dies ermöglicht es uns eine vektorisierte Form von $h(x, w)$ leicht aufzuschreiben, die sowohl mit einem Parameter `x` bestehend aus einer Zeile an Inputdaten (also z.B. einem einzelnen Trainingsbeispiel) funktioniert als auch mit der gesamten Feature-Matrix `X`, bestehend aus allen (oder mehreren) Trainingsdaten auf einmal."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "14116a52",
   "metadata": {},
   "outputs": [],
   "source": [
    "def h(x, w):\n",
    "    \"\"\"x und w sind numpy arrays; x kann auch die komplette Feature-Matrix sein\"\"\"\n",
    "    # Diese Form erlaubt es für x eine ganze (Feature-)Matrix zu übergeben. Die Matrix enthält\n",
    "    # zeilenweise je einen Datenpunkt, für den h berechnet werden soll.\n",
    "    # w @ x.T ist dann ein Vektor mit je einem Ergebnis in den Komponenten des Vektors pro Zeile\n",
    "    # der übergebenen (Feature-)Matrix.\n",
    "    return x @ w"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "82129a25",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Definition der Kostenfunktion\n",
    "def J(w, X, y):\n",
    "    \"\"\"\n",
    "    w, X, y müssen numpy arrays sein\n",
    "    X: Feature-Matrix aller Trainingsdaten inkl. Spalte mit 1; Dimension: n x (d+1)\n",
    "    y: Vektor aller Targets zu X\n",
    "    \"\"\"\n",
    "    errors = y - h(x=X, w=w)\n",
    "    mse = 1.0/(2.0*len(y)) * ( errors @ errors )\n",
    "    return mse"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "4209dc9c",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(6196, 3)"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "8b34a5c7",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1.],\n",
       "       [1.],\n",
       "       [1.],\n",
       "       ...,\n",
       "       [1.],\n",
       "       [1.],\n",
       "       [1.]])"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.ones((len(data),1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "6632cabe",
   "metadata": {},
   "outputs": [],
   "source": [
    "def feature_matrix_from_data(data):\n",
    "    # hier erzeugen wir die Matrix mit unseren Input-Daten (Features) inklusive der Spalte mit \"1\"\n",
    "    return np.hstack((np.ones((len(data),1)), data.to_numpy(copy=True)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "74556a69",
   "metadata": {},
   "outputs": [],
   "source": [
    "# hier erzeugen wir die Matrix mit unseren Input-Daten (Features) inklusive der Spalte mit \"1\"\n",
    "#X = np.hstack((np.ones((len(data),1)), data[features].to_numpy(copy=True)))\n",
    "X = feature_matrix_from_data(data[features])\n",
    "# und ausserdem den Vektor der Targets\n",
    "y = data.Price.to_numpy(copy=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "79f5e3e0",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(6196, 3)"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "8f9724c3",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[  1.  ,   2.  ,  79.  ],\n",
       "       [  1.  ,   3.  , 150.  ],\n",
       "       [  1.  ,   4.  , 142.  ],\n",
       "       ...,\n",
       "       [  1.  ,   1.  ,  35.64],\n",
       "       [  1.  ,   2.  ,  61.6 ],\n",
       "       [  1.  ,   6.  , 388.5 ]])"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "85733d6b",
   "metadata": {},
   "source": [
    "**Hinweis:** Die Matrix $X$ hat zwar die gleiche Dimension wie `data`, allerdings enthält data eine Spalte `Price`, die in $X$ nicht enthalten ist. Dafür hat $X$ als erste Spalte die \"1er\"."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1d8e64e6",
   "metadata": {},
   "source": [
    "## Analytische Lösung der linearen Regression\n",
    "\n",
    "Die analytische Lösung verläuft identisch zum Fall mit nur einem Feature.\n",
    "\n",
    "`np.linalg.solve(A, b)` berechnet $w$ im linearen Gleichungssystem\n",
    "\n",
    "$ A w = b $\n",
    "\n",
    "$A$ - Matrix,\n",
    "$w$ - Vektor (unsere unbekannten),\n",
    "$b$ - Vektor.\n",
    "\n",
    "Wir suchen die Lösung $w$ im folgenden Gleichungssystem:\n",
    "\n",
    "$$ X^T X w = X^T Y $$\n",
    "\n",
    "Mit $A = X^TX$ und $b = X^T Y$ berechnet `np.linalg.solve(A, b)` unsere gesuchten Paramter für die lineare Regression."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "fc1d2c0a",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Die 3 Parameter der linearen Regression:\n",
      "[ 42769.88494072 232612.86504788   2431.15453776]\n",
      "Kostenfunktion J(w_ana): 147658829426.14856\n",
      "CPU times: user 11.3 ms, sys: 2.13 ms, total: 13.4 ms\n",
      "Wall time: 1.7 ms\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "w_ana = np.linalg.solve(X.T @ X, X.T @ y)\n",
    "print('Die {} Parameter der linearen Regression:\\n{}'.format(len(w_ana), w_ana))\n",
    "J_ana = J(w=w_ana, X=X, y=y)\n",
    "print('Kostenfunktion J(w_ana): {}'.format(J_ana))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "daab4572",
   "metadata": {},
   "source": [
    "## Numerische Lösung mit Gradient Descent"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "b314f36a",
   "metadata": {},
   "outputs": [],
   "source": [
    "## Numerische Lösung mit Gradient Descent\n",
    "def grad_desc_upd(w, alpha, x, y):\n",
    "    \"\"\"y, x sind Vektoren (numpy-arrays)\"\"\"\n",
    "    errors = y - h(x=x, w=w)\n",
    "    w = w + alpha / len(y) * (x.T @ errors)\n",
    "    return w"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "3dc2775c",
   "metadata": {},
   "outputs": [],
   "source": [
    "def grad_desc(w, alpha, x, y, n_iterations):\n",
    "    J_all = [[0], [J(w=w, X=x, y=y)]]\n",
    "    for it in tqdm(range(n_iterations)):\n",
    "        w = grad_desc_upd(w=w, alpha=alpha, x=x, y=y)\n",
    "        if it % 100 == 0:\n",
    "            J_all[1].append(J(w=w, X=x, y=y))\n",
    "            J_all[0].append(it)\n",
    "    return w, J_all"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "a801cac3",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([  2.43686014,   5.00088371, 206.19316114])"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "grad_desc_upd(w=np.ones(X.shape[1]), alpha=1e-6, x=X[:7], y=y[:7])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "dc6f778a",
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "e9c9403f9b08472294edb52bb2c10c1d",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "  0%|          | 0/10000 [00:00<?, ?it/s]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 1.94 s, sys: 1.87 s, total: 3.81 s\n",
      "Wall time: 417 ms\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "w_init = np.ones(X.shape[1])\n",
    "alpha = 3.1e-10  # verschiedene alpha ausprobieren\n",
    "n_iterations = 10000\n",
    "_, J_tmp = grad_desc(w=w_init, alpha=alpha, x=X, y=y, n_iterations=n_iterations)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "c04ebb9f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "57bac370a44f48e9951dc5dba56b29ef",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "  0%|          | 0/100000 [00:00<?, ?it/s]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Die 3 Parameter der linearen Regression:\n",
      "[  10.49970832   31.89447067 1601.95876825]\n",
      "Kostenfunktion J: 540959857400.77966\n",
      "J relativ zu Startkosten: 0.6771395257663181\n",
      "Vergleich Kostenfunktion zu analytischer Lösung: 3.66*J_ana\n",
      "Relative Abweichung der Parameter zu analytischer Lösung: [2.45493022e-04 1.37113958e-04 6.58929222e-01]*w_ana\n",
      "CPU times: user 21.5 s, sys: 11 s, total: 32.6 s\n",
      "Wall time: 3.53 s\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "w_init = np.ones(X.shape[1])\n",
    "alpha = 1e-10  # verschiedene alpha ausprobieren\n",
    "n_iterations = 100000\n",
    "w_gd, J_all = grad_desc(w=w_init, alpha=alpha, x=X, y=y, n_iterations=n_iterations)\n",
    "print('Die {} Parameter der linearen Regression:\\n{}'.format(len(w_gd), w_gd))\n",
    "print('Kostenfunktion J: {}'.format(J_all[1][-1]))\n",
    "print('J relativ zu Startkosten: {}'.format(J_all[1][-1]/J_all[1][0]))\n",
    "print('Vergleich Kostenfunktion zu analytischer Lösung: {:.2f}*J_ana'.format(J_all[1][-1]/J_ana))\n",
    "print('Relative Abweichung der Parameter zu analytischer Lösung: {}*w_ana'.format((w_gd)/w_ana))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "8b4db3ee",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<Axes: >"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAicAAAGsCAYAAAAGzwdbAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAA9hAAAPYQGoP6dpAABHWUlEQVR4nO3deVhU9eIG8PfMDAwgzIDIqoC4i6KiuAAuebXILc3UNBT3LU2xsvLX9bZY2aZl5r6XW+7mHmqugIqKiiiKKKCyuMGA7Mz5/WHNjevGIHBmeT/Pc57ncuZ7Zt45VLx35pzvVxBFUQQRERGRgZBJHYCIiIjon1hOiIiIyKCwnBAREZFBYTkhIiIig8JyQkRERAaF5YSIiIgMCssJERERGRSWEyIiIjIoLCdERERkUFhOiIiIyKAYVTk5cuQIevXqBXd3dwiCgG3btul1fH5+PoYNGwZfX18oFAr06dPnsTGpqal466230KBBA8hkMoSFhVVIdiIiIioboyonDx8+RPPmzTFv3rxyHV9SUgJra2tMmjQJXbt2feKYgoICODk54d///jeaN2/+InGJiIioHBRSB9BHt27d0K1bt6c+XlBQgI8//hjr1q1DZmYmmjZtim+++QYvvfQSAKBatWpYsGABAOD48ePIzMx87Dlq166NOXPmAACWL19e4e+BiIiIns2oPjl5nokTJyIyMhLr16/H+fPn0b9/f7z66qu4evWq1NGIiIiojEymnCQnJ2PFihXYuHEjOnTogLp16+L9999H+/btsWLFCqnjERERURkZ1dc6z3LhwgWUlJSgQYMGpfYXFBTA0dFRolRERESkL5MpJzk5OZDL5Th9+jTkcnmpx2xtbSVKRURERPoymXLi5+eHkpISZGRkoEOHDlLHISIionIyqnKSk5ODhIQE3c/Xr19HTEwMqlevjgYNGiAkJAShoaGYNWsW/Pz8cOfOHRw4cADNmjVDjx49AABxcXEoLCzE/fv3kZ2djZiYGABAixYtdM/7976cnBzcuXMHMTExsLS0hI+PT1W9VSIiIrMliKIoSh2irA4dOoTOnTs/tn/o0KFYuXIlioqK8MUXX+CXX37BrVu3UKNGDbRr1w6fffYZfH19ATy6VTgpKemx5/jnaRAE4bHHvby8cOPGjYp7M0RERPRERlVOiIiIyPSZzK3EREREZBpYToiIiMigGMUFsVqtFrdv34adnd0TrwchIiIiwyOKIrKzs+Hu7g6ZrOyfhxhFObl9+zY8PDykjkFERETlkJKSglq1apV5vFGUEzs7OwCP3pxKpZI4DREREZWFRqOBh4eH7u94WRlFOfn7qxyVSsVyQkREZGT0vSSDF8QSERGRQWE5ISIiIoPCckJEREQGheWEiIiIDArLCRERERkUlhMiIiIyKCwnREREZFBYToiIiMig6FVOSkpKMH36dHh7e8Pa2hp169bFjBkzIIriM487dOgQWrZsCaVSiXr16mHlypUvkpmIiIhMmF4zxH7zzTdYsGABVq1ahSZNmiA6OhrDhw+HWq3GpEmTnnjM9evX0aNHD4wbNw5r1qzBgQMHMGrUKLi5uSE4OLhC3gQRERGZDkF83sce/9CzZ0+4uLhg2bJlun1vvPEGrK2tsXr16ice8+GHH2LXrl2IjY3V7Rs4cCAyMzOxd+/eMr2uRqOBWq1GVlYWp68nIiIyEuX9+63X1zqBgYE4cOAArly5AgA4d+4cjh07hm7duj31mMjISHTt2rXUvuDgYERGRj71mIKCAmg0mlIbERERmQe9yslHH32EgQMHolGjRrCwsICfnx/CwsIQEhLy1GPS0tLg4uJSap+Liws0Gg3y8vKeeMzMmTOhVqt1m4eHhz4xy+zIlTsY/Us08otKKuX5iYiISH96lZMNGzZgzZo1WLt2Lc6cOYNVq1bh+++/x6pVqyo01LRp05CVlaXbUlJSKvT5AeBhQTHCfotBeFw6hq84hdzC4gp/DSIiItKfXuVk6tSpuk9PfH19MWTIEEyZMgUzZ8586jGurq5IT08vtS89PR0qlQrW1tZPPEapVEKlUpXaKlo1pQKLhrSCrVKByMR7CF12Etn5RRX+OkRERKQfvcpJbm4uZLLSh8jlcmi12qceExAQgAMHDpTaFx4ejoCAAH1eulK0rl0dv45sA5WVAtFJDzB42Ulk5bKgEBERSUmvctKrVy98+eWX2LVrF27cuIGtW7di9uzZeP3113Vjpk2bhtDQUN3P48aNQ2JiIj744ANcvnwZ8+fPx4YNGzBlypSKexcvwM/TAWtHt4ODjQXOpWRi0JIo3H9YKHUsIiIis6VXOZk7dy769euHt99+G40bN8b777+PsWPHYsaMGboxqampSE5O1v3s7e2NXbt2ITw8HM2bN8esWbOwdOlSg5rjpGlNNdaPCUANW0vEpWowcHEk7mQXSB2LiIjILOk1z4lUqmqek4SMHIQsjUK6pgB1nKph7ah2cFVbVdrrERERmbIqmefE1NVztsWGsQGoaW+NxDsPMWBRJG4+yJU6FhERkVlhOfkfXo7V8NvYdvCsboPk+7l4c1EUbtx9KHUsIiIis8Fy8gS1HGywYWwA6jhVw63MPLy5OBIJGTlSxyIiIjILLCdP4aq2wm9jAtDQxQ7pmgIMXByJy2mcRp+IiKiysZw8g5OdEuvGtEMTdxXu5hRi4OIoxN7KkjoWERGRSWM5eY7q1SyxdlQ7NPewR2ZuEQYticLZ5AdSxyIiIjJZLCdloLaxwOqRbdC6tgOy84sxeOkJRCXekzoWERGRSWI5KSM7KwusGtEGgXUd8bCwBEOXn8Sh+AypYxEREZkclhM92FgqsHxYa3Rp5IyCYi1G/xKNPRdSpY5FRERkUlhO9GRlIcfCIa3Qo5kbikpETFh7BlvO3JQ6FhERkclgOSkHC7kMPw30wwD/WtCKwLsbzuHXqCSpYxEREZkElpNykssEfN23GYYF1gYATN8Wi0WHr0kbioiIyASwnLwAmUzAJ718MKFzXQDAzD2XMTv8CoxgLUUiIiKDxXLyggRBwNTgRpga3BAA8NOBq/hi1yUWFCIionJiOakgEzrXw2evNQEALDt2Hf+39QJKtCwoRERE+mI5qUBDA2vj237NIBOAdSdT8O6GGBSVaKWORUREZFRYTirYAH8P/DTIDwqZgO0xt/H2mjMoKC6ROhYREZHRYDmpBD2buWNxaCtYKmQIj0vHqFXRyC0sljoWERGRUWA5qST/auSClcNaw8ZSjqNX72Lo8pPQ5BdJHYuIiMjgsZxUosB6NfDryLaws1Lg1I0HCFlyAg8eFkodi4iIyKCxnFSyVl4OWDe6HapXs8SFW1l4c3EkMjT5UsciIiIyWCwnVaBpTTU2jG0HF5USV9JzMGBRJFLu50odi4iIyCCxnFSRes522Dg2EB7VrXHjXi76L4zE1fRsqWMREREZHJaTKuTpaIONYwNR39kWaZp8DFgUiXMpmVLHIiIiMigsJ1XMVW2FDWMD0LyWGg9yi/DWkihEXrsndSwiIiKDwXIiAYdqllgzuh0C6zriYWEJhq44ifC4dKljERERGQSWE4nYKhVYPqw1XvFxQWGxFuNWn8aWMzeljkVERCQ5lhMJWVnIMT+kJfq2rIkSrYh3N5zDyuPXpY5FREQkKZYTiSnkMnzfrzmGBdYGAHy6Iw5z9l+FKHJFYyIiMk8sJwZAJhPwSS8fTOnaAADww/4r+HxnHLRaFhQiIjI/LCcGQhAETO5aH5/08gEArDh+A1M3nUdxiVbiZERERFWL5cTADA/yxuwBzSGXCdh85ibeXnMG+UUlUsciIiKqMiwnBqhvy1pYOLgVLBUy/BGXjhErTyGnoFjqWERERFWC5cRAvezjgpXDW6OapRwR1+4hZEkUVzQmIiKzwHJiwALr1sC6Me3gYGOBczezMGBRJNKyuKIxERGZNpYTA9eslj02jguAq8oKVzNy0G9hBG7cfSh1LCIiokrDcmIE6jnbYdP4ANR2tMHNB3notzASl1I1UsciIiKqFCwnRqKWgw02jgtEYzcV7uYUYMCiSJxI5IKBRERkelhOjIiTnRLrx7RDm9rVkZ1fjCHLT+KPi2lSxyIiIqpQLCdGRm1tgV9GtkHXxv9dMHDDqRSpYxEREVUYlhMjZGUhx8LBLdG/VS1oReCDzeex4NA1rsdDREQmgeXESCnkMnzbrxnGdaoLAPhm72V8sesS1+MhIiKjx3JixARBwEfdGuHfPRoDAJYdu473Np5DEdfjISIiI8ZyYgJGdaiD2QOaQyETsPXsLYz+JRq5hZzunoiIjBPLiYno27IWloT6w8pChkPxdxCy9ASnuyciIqPEcmJCOjdyxppR7aC2tsDZ5Ez0XxSJ25l5UsciIiLSC8uJiWnl5aCb7j4hIwf9FkQgISNb6lhERERlxnJighq42GHz24Go41QNt7Py0W9hJM4mP5A6FhERUZmwnJiomvbW2DQuEM1rqZGZW4S3lpzA4St3pI5FRET0XCwnJqx6NUusHd0OHerXQF5RCUauPIXtMbekjkVERPRMLCcmrppSgWVDW6NXc3cUa0VMXh+DFcevSx2LiIjoqVhOzIClQoY5b7bAsMDaAIDPdsTh272XOd09EREZJJYTMyGTCfiklw/ef6UBAGD+oWt4f+N5ziZLREQGh+XEjAiCgIn/qo+v+/pCLhOw+cxNjFoVjYcFnE2WiIgMh17lpHbt2hAE4bFtwoQJTxy/cuXKx8ZaWVlVSHAqv4FtPLF4SCtYWchw+ModDFoShbs5BVLHIiIiAqBnOTl16hRSU1N1W3h4OACgf//+Tz1GpVKVOiYpKenFElOF6NLYBetGt4ODjQXO38zCGwsikHTvodSxiIiI9CsnTk5OcHV11W07d+5E3bp10alTp6ceIwhCqWNcXFxeODRVDD9PB2weH4haDtZIupeLvvMjcP5mptSxiIjIzJX7mpPCwkKsXr0aI0aMgCAITx2Xk5MDLy8veHh4oHfv3rh48eJzn7ugoAAajabURpWjjpMttrwdiCbuKtx7WIiBi6NwKD5D6lhERGTGyl1Otm3bhszMTAwbNuypYxo2bIjly5dj+/btWL16NbRaLQIDA3Hz5s1nPvfMmTOhVqt1m4eHR3ljUhk421lh/Zh2aF+vBnILSzBqVTQ2nX7274iIiKiyCGI5J7sIDg6GpaUlduzYUeZjioqK0LhxYwwaNAgzZsx46riCggIUFPz3Ak2NRgMPDw9kZWVBpVKVJy6VQWGxFh9sOodtMbcBAFODG+Ltl+o+85MxIiKip9FoNFCr1Xr//VaU58WSkpKwf/9+bNmyRa/jLCws4Ofnh4SEhGeOUyqVUCqV5YlGL8BSIcPsAS3gorLCoiOJ+G5fPNKy8vHpa00gl7GgEBFR1SjX1zorVqyAs7MzevTooddxJSUluHDhAtzc3MrzslQFZDIB07o3xn96+kAQgF+jkjBhzRnkF5VIHY2IiMyE3uVEq9VixYoVGDp0KBSK0h+8hIaGYtq0abqfP//8c/zxxx9ITEzEmTNnMHjwYCQlJWHUqFEvnpwq1Yj23pg7yA+Wchn2XkxD6LKTyMotkjoWERGZAb3Lyf79+5GcnIwRI0Y89lhycjJSU1N1Pz948ACjR49G48aN0b17d2g0GkRERMDHx+fFUlOV6NnMHatGtIGdUoGTN+6j38II3M7MkzoWERGZuHJfEFuVyntBDVWMy2kaDF1+EumaAriqrLBqRBs0dLWTOhYRERm48v795to69FyNXFXY8nYQ6jnbIk2Tj34LIxCVeE/qWEREZKJYTqhMatpbY9O4APh7OSA7vxihy07i93O3pY5FREQmiOWEyszexhKrR7XFq01cUViixaR1Z7Hw8DUYwTeDRERkRFhOSC9WFnLMC2mJ4UG1AQBf77mM6dtjUVyilTYYERGZDJYT0ptcJuCTXk0w/a+5UFZHJWPc6tPILSyWOhoREZkAlhMqt5HtvTH/rZZQKmTYfykDgxZH4U52wfMPJCIiegaWE3oh3XzdsHZ0WzjYWODczSz0XXAc1+7kSB2LiIiMGMsJvbBWXtWx5e0geDnaIOV+Ht5YEIFTN+5LHYuIiIwUywlVCO8a1bBlfCBaeNgjM7cIIUtPYOd53mpMRET6YzmhCuNoq8S60e3wio8LCou1mLj2LBYf4a3GRESkH5YTqlDWlnIsGNwKwwJrAwC+2n0Zn/5+ESVaFhQiIioblhOqcI9uNfbBv3s0BgCsikzCuNWnkVdYInEyIiIyBiwnVCkEQcCoDnUw762WsFTIEB6XjoFLonA3h7caExHRs7GcUKXq0cwNa0e1hb2NBc6lZKLv/Agk8lZjIiJ6BpYTqnT+tatjy/hAeFa3QfL9XPRdEIFo3mpMRERPwXJCVaKOky22vB2I5rXUyMwtwltLT2AHVzUmIqInYDmhKlPDVol1Y9qha+NHtxq/s+4s5v2ZwFuNiYioFJYTqlI2lgosGtIKI4K8AQDf7YvHh5vPo4irGhMR0V9YTqjKyWUC/tPLB5/3bgKZAGyIvomhy08iK69I6mhERGQAWE5IMqEBtbFsaGtUs5Qj4to9vLEgAin3c6WORUREEmM5IUl1buSMDeMC4KqyQkJGDvrMO44zyQ+kjkVERBJiOSHJNXFXY9uEIDRxV+Hew0IMWhyF3RdSpY5FREQSYTkhg+CqtsKGsQHo0sgZBcVavL3mDBYc4qKBRETmiOWEDEY1pQKLQ/11iwZ+s/cypm25wDt5iIjMDMsJGRS5TMCnrzXBJ718IBOA9adSMHzFKWjyeScPEZG5YDkhgzQ8yBtLQv1hYynHsYS7eGM+7+QhIjIXLCdksLo0dsGGsQFwUSlxNSMHr88/jpiUTKljERFRJWM5IYPWtOajO3kau6lwN6cQby6KxB7eyUNEZNJYTsjguamtsXFcADo3dEJBsRbj15zBosO8k4eIyFSxnJBRsFUqsCTUH0MDvAAAM/dcxv9t5Z08RESmiOWEjIZCLsNnvZvik14+EARg3ckUDF1+Epm5hVJHIyKiCsRyQkZneJA3lob669bkeX1+BBLv5Egdi4iIKgjLCRmlLo1dsGl8IGraW+P63Yd4fX4EIhLuSh2LiIgqAMsJGa3GbipsmxAEP097ZOUVIXT5Saw9kSx1LCIiekEsJ2TUnOyUWDe6HV5r7o5irYj/23oBn++IQ4mWd/IQERkrlhMyelYWcswZ2ALvvtwAALD8+HWM/iUa2ZzynojIKLGckEkQBAGTutTHz2/5QamQ4eDlDPRbEMkp74mIjBDLCZmUns3csWFsAJzslIhPz0afecdxOum+1LGIiEgPLCdkcpp72OP3iUHwcVPh3sNCDFp8AtvO3pI6FhERlRHLCZkkN7U1No0PwCs+Ligs0SLstxh8vy8eWl4oS0Rk8FhOyGTZWCqwcHArjH+pLgDg5z8TMHHdGeQVlkicjIiInoXlhEyaTCbgw1cb4bt+zWAhF7D7QhreXByJdE2+1NGIiOgpWE7ILPT398CaUe3gYGOB8zez0Pvn44i9lSV1LCIiegKWEzIbbbyrY/uE9qjnbIs0TT76L4zEngupUsciIqL/wXJCZsXT0QZb3g5ExwZOyCsqwfg1Z/Dj/iu8UJaIyICwnJDZUVlZYPlQf4wI8gYA/Lj/KiauO4PcwmKJkxEREcByQmZKIZfhP7188O0b/71Qtt+CSNzKzJM6GhGR2WM5IbM2oLUH1o5uB8dqlohL1aD3z8cQfYMzyhIRSYnlhMxe69rV8fs77dHYTYW7OYUYtCQKG6JTpI5FRGS2WE6IANS0t8bm8QHo1tQVRSUiPth0HjN2xqG4RCt1NCIis8NyQvQXG0sF5r3VEmFd6wMAlh27jhGropGVVyRxMiIi88JyQvQPMpmAsK4NMD+kJawt5Dhy5Q5en3cc1+7kSB2NiMhssJwQPUF3XzdsGh8Ad7UVEu8+RJ95x3H4yh2pYxERmQWWE6KnaOKuxvaJ7dHKywHZ+cUYvuIklh5NhChywjYiosqkVzmpXbs2BEF4bJswYcJTj9m4cSMaNWoEKysr+Pr6Yvfu3S8cmqiqONkpsXZ0WwzwrwWtCHyx6xKmbjqPgmKubExEVFn0KienTp1CamqqbgsPDwcA9O/f/4njIyIiMGjQIIwcORJnz55Fnz590KdPH8TGxr54cqIqolTI8c0bzTC9pw9kArDp9E28teQE7mQXSB2NiMgkCeILfEYdFhaGnTt34urVqxAE4bHH33zzTTx8+BA7d+7U7WvXrh1atGiBhQsXlvl1NBoN1Go1srKyoFKpyhuX6IUdvnIHE9eeQXZ+MdzUVlgS6o+mNdVSxyIiMkjl/ftd7mtOCgsLsXr1aowYMeKJxQQAIiMj0bVr11L7goODERkZ+cznLigogEajKbURGYJODZywfUIQ6tSohtSsfPRbGIHfz92WOhYRkUkpdznZtm0bMjMzMWzYsKeOSUtLg4uLS6l9Li4uSEtLe+Zzz5w5E2q1Wrd5eHiUNyZRhavjZIutE4LQsYET8ou0mLTuLGbuuYQSrmxMRFQhyl1Oli1bhm7dusHd3b0i8wAApk2bhqysLN2WksKpxMmwqK0tsGJYa4ztVAcAsOhwIoavPIWsXE7YRkT0ospVTpKSkrB//36MGjXqmeNcXV2Rnp5eal96ejpcXV2feZxSqYRKpSq1ERkauUzAtG6N8dMgP1hZyHDkyh30nncMV9KzpY5GRGTUylVOVqxYAWdnZ/To0eOZ4wICAnDgwIFS+8LDwxEQEFCelyUySK81d8emcYGoaW+NG/dy8fq849h38dlfXRIR0dPpXU60Wi1WrFiBoUOHQqFQlHosNDQU06ZN0/08efJk7N27F7NmzcLly5fx6aefIjo6GhMnTnzx5EQGpGlNNX6fGIR2darjYWEJxv56GrPDr0DL61CIiPSmdznZv38/kpOTMWLEiMceS05ORmpqqu7nwMBArF27FosXL0bz5s2xadMmbNu2DU2bNn2x1EQGyNFWiV9HtsXwoNoAgJ8OXMWYX08jO5/XoRAR6eOF5jmpKpznhIzNxugUfLwtFoXFWtRztsXiIa1Qx8lW6lhERFWqyuc5IaKn6+/vgQ1jA+CqskJCRg56zzuOPy9nSB2LiMgosJwQVZIWHvb4/Z0g3cKBI1adwrw/E7hwIBHRc7CcEFUiZzsrrBvdDoPaeEIUge/2xWPi2rPILSyWOhoRkcFiOSGqZJYKGWb29cWXrzeFhVzArgup6Ds/Ain3c6WORkRkkFhOiKpISFsvrB3dDjVslbiclo1ePx/D8YS7UsciIjI4LCdEVah17erY8U4QmtdSIzO3CEOWncDSo4m8DoWI6B9YToiqmJvaGr+NDcAbLWtBKwJf7LqEyetjeB0KEdFfWE6IJGBlIcf3/Zvh014+UMgE/H7uNvrOj8CNuw+ljkZEJDmWEyKJCIKAYUHej12HcvBy+vMPJiIyYSwnRBJr410duya1R0tPe2TnF2PkqmjM2X+V6/IQkdliOSEyAC4qK6wfE4Ah7bwgisAP+69g9C/RyMrjujxEZH5YTogMhKVChhl9muK7fs1gqZDhwOUM9P75GOLTsqWORkRUpVhOiAxMf38PbB4XiJr21rhxLxd95h3HjnO3pY5FRFRlWE6IDJBvLTV2vNMe7evVQF5RCd5ZdxZf7opDcYlW6mhERJWO5YTIQFWvZolVI9pgXKe6AIAlR69j8LITuJtTIHEyIqLKxXJCZMDkMgEfdWuEBSEtUc1SjqjE++g19xhiUjKljkZEVGlYToiMQDdfN2yfGIQ6TtWQmpWPAQsj8dupZKljERFVCpYTIiNRz9kO2ycE4WUfFxSWaPHh5guYtuUCCopLpI5GRFShWE6IjIidlQUWDW6FqcENIQjAupPJGLAoCrcz86SORkRUYVhOiIyMTCZgQud6WDm8DdTWFjiXkomec4/h2NW7UkcjIqoQLCdERqpTAyfsfKc9mrircP9hIYYsP4GfD3LaeyIyfiwnREbMo7oNNo8PxJv+HhBF4Ps/rmDUL9HIyuW090RkvFhOiIyclYUc3/Rrhm/faAalQoaDlzPQ8+ejiL2VJXU0IqJyYTkhMhEDWntg8/hAeFa3Qcr9PPRdEMHbjYnIKLGcEJmQpjXV2DGxPbo2dkZh8aPbjaduPIf8It5uTETGg+WEyMSobSyweIg/pgY3hEwANp6+ib7zI5B076HU0YiIyoTlhMgE/X278eqRbeFYzRJxqRr0nHsM4XHpUkcjInoulhMiExZYrwZ2TeqAlp72yM4vxuhfovHt3stc3ZiIDBrLCZGJc1VbYf2YAAwPqg0AmH/oGoYsO4k72VzdmIgME8sJkRmwVMjwSa8mmDvIDzaWckQm3kPPuUcRfeO+1NGIiB7DckJkRno1d8fvE4NQz9kW6ZoCDFwcheXHrkMUOassERkOlhMiM/P36sY9m7mhWCvi851xmLjuLHIKiqWORkQEgOWEyCxVUyowd5AfPu3lA4VMwK7zqXjt52O4nKaROhoREcsJkbkSBAHDgrzx29h2cFVZIfHOQ/SZdxwbolOkjkZEZo7lhMjMtfKqjl2T2qNjAyfkF2nxwabzeH/jOeQVclZZIpIGywkRwdFWiZXDWuP9VxpAJgCbTt9E73nHkJCRLXU0IjJDLCdEBODRrLIT/1Ufa0a1g5OdElfSc/Daz8ex7ewtqaMRkZlhOSGiUgLqOmLXpPYIqOOI3MIShP0Wg2lbLnDxQCKqMiwnRPQYZzsrrB7VFpO61IcgAOtOJqPv/AjcuMvFA4mo8rGcENETyWUC3n25AVYNb1Nq8cBd51OljkZEJo7lhIieqWMDJ+ya1AGtazsgp6AYE9aewae/X0RBMb/mIaLKwXJCRM/lqrbCutHtMK5TXQDAyogbGLAwEin3cyVORkSmiOWEiMpEIZfho26NsHyYP+xtLHDuZhZ6/HQU4XHpUkcjIhPDckJEevlXIxfsmtQBfp720OQXY/Qv0fhyVxyKSrRSRyMiE8FyQkR6q2lvjd/GBGBke28AwJKj1/HmokjczsyTOBkRmQKWEyIqF0uFDNN7+mDh4Faws1LgTHImuv90FAcu8WseInoxLCdE9EJebeqKXe90gG9NNTJzizByVTS+2BmHwmJ+zUNE5cNyQkQvzNPRBpvGB2BE0KOveZYeu47+CyOQfI938xCR/lhOiKhCKBVy/KeXD5aE+kNt/d+7eXaevy11NCIyMiwnRFShXvZxwe7JHeDv5YDsgmJMXHsW/7eVa/MQUdmxnBBRhatpb431Y9phQue6EARg7Ylk9Jl3HAkZOVJHIyIjwHJCRJVCIZdhanAj/DKiDWrYWuJyWjZ6zT2GTadvSh2NiAwcywkRVaoO9Z2we3IHBNVzRF5RCd7feA7v/haDhwXFUkcjIgPFckJElc7Zzgq/jGiL919pAJkAbDl7C71+Poa42xqpoxGRAdK7nNy6dQuDBw+Go6MjrK2t4evri+jo6KeOP3ToEARBeGxLS0t7oeBEZFzkMgET/1Uf68cEwFVlhcQ7D9Fn/nH8GpUEURSljkdEBkSvcvLgwQMEBQXBwsICe/bsQVxcHGbNmgUHB4fnHhsfH4/U1FTd5uzsXO7QRGS82nhXx+7JHdClkTMKi7WYvi0Wb685g6y8IqmjEZGBUOgz+JtvvoGHhwdWrFih2+ft7V2mY52dnWFvb69XOCIyTdWrWWLpUH8sO3Yd3+y9jD2xabhwKws/v9USLTzspY5HRBLT65OT33//Hf7+/ujfvz+cnZ3h5+eHJUuWlOnYFi1awM3NDS+//DKOHz/+zLEFBQXQaDSlNiIyLYIgYFSHOtg0LhAe1a1x80Ee+i2IwOIj16DV8mseInOmVzlJTEzEggULUL9+fezbtw/jx4/HpEmTsGrVqqce4+bmhoULF2Lz5s3YvHkzPDw88NJLL+HMmTNPPWbmzJlQq9W6zcPDQ5+YRGREmnvYY9ekDujRzA3FWhFf7b6MYStP4U52gdTRiEgigqjHlWiWlpbw9/dHRESEbt+kSZNw6tQpREZGlvlFO3XqBE9PT/z6669PfLygoAAFBf/9D5NGo4GHhweysrKgUqnK/DpEZDxEUcS6kyn4fOdF5BdpUcPWEt/3b46XGvL6NCJjpdFooFar9f77rdcnJ25ubvDx8Sm1r3HjxkhOTtbnadCmTRskJCQ89XGlUgmVSlVqIyLTJggC3mrriR0T26ORqx3u5hRi2IpT+GJnHAqKOfU9kTnRq5wEBQUhPj6+1L4rV67Ay8tLrxeNiYmBm5ubXscQkXmo72KHbROCMDTg0X9Xlh67jr7zI5B4h1PfE5kLvcrJlClTEBUVha+++goJCQlYu3YtFi9ejAkTJujGTJs2DaGhobqff/zxR2zfvh0JCQmIjY1FWFgYDh48WOoYIqJ/srKQ47PeTbE01B8ONha4eFuDnnOPYUN0CudEITIDepWT1q1bY+vWrVi3bh2aNm2KGTNm4Mcff0RISIhuTGpqaqmveQoLC/Hee+/B19cXnTp1wrlz57B//3506dKl4t4FEZmkrj4u2DO5IwLrOiK3sAQfbDqPSetjoMnnnChEpkyvC2KlUt4LaojINJRoRSw6cg2z/riCEq2IWg7WmDPQD628nj8BJBFJp0ouiCUikoJcJuDtl+ph07gA3ZwoAxZF4ueDV1HCOVGITA7LCREZDT9PB+ye1AG9W7ijRCvi+z+uIGRpFFKz8qSORkQViOWEiIyKnZUFfnyzBWb1b45qlnJEJd5HtzlHse8iFxMlMhUsJ0RkdARBwButamHnpA7wralGZm4Rxv56GtO3xSK/iHOiEBk7lhMiMlreNaph8/hAjO1YBwDwa1QSev98HPFp2RInI6IXwXJCREbNUiHDtO6N8cuINqhhq0R8ejZe+/kYfo28wTlRiIwUywkRmYSODZywN6wDOjd0QkGxFtO3X8TIVdFcQJDICLGcEJHJqGGrxPJhrfFJLx9YKmQ4eDkD3eYcwcHL6VJHIyI9sJwQkUkRBAHDg7xLLSA4YmU0pm+LRV4hL5YlMgYsJ0Rkkhq6PlpAcGR7bwCPLpbt9fMxxN7KkjgZET0PywkRmSwrCzmm9/TBryPbwNlOiYSMHLw+/zgWHb4GLWeWJTJYLCdEZPI61HfC3rCOCG7igqISETP3XMbgZSc4syyRgWI5ISKzUL2aJRYOboWv+/rC2kKOiGv38OqPR7HrfKrU0Yjof7CcEJHZEAQBA9t4YvfkDmheS42svCJMWHsG7204h5yCYqnjEdFfWE6IyOx416iGTeMDMbFzPcgEYPOZm+g+5yhOJz2QOhoRgeWEiMyUhVyG94MbYv2YANS0t0by/VwMWBSJH8KvoLhEK3U8IrPGckJEZq2Nd3XsCeuAPi3cUaIVMefAVfRfFImkew+ljkZktlhOiMjsqaws8ONAP8wZ2AJ2SgXOJmei+5yj2BidwvV5iCTAckJE9JfeLWpiT1gHtKldHQ8LSzB103mMX30G9x8WSh2NyKywnBAR/UMtBxusG9MOU4MbQiETsPdiGl75gevzEFUllhMiov8hlwmY0Lketk0IQn1nW9zNKcCIldGYtuUCHvKWY6JKx3JCRPQUTWuqseOd9rr1edadTEb3n3jLMVFlYzkhInqGv9fnWTuqLdzVVki6l4v+CyPw3b7LKCzmLcdElYHlhIioDALr1cCesI7o61cTWhGY9+c1vD7/OK6mZ0sdjcjksJwQEZWR2toCs99sgfkhLWFvY4GLtzXoMfcYlh27zlWOiSoQywkRkZ66+7rhj7COeKmhEwqLtZixMw6Dl53ArUyuckxUEVhOiIjKwVllhRXDWuOLPk3/scrxEWw9e5MTtxG9IJYTIqJyEgQBg9t5YffkDvDztEd2fjGm/HYOE9aewQNO3EZUbiwnREQvyLtGNWwcG4D3Xm4AhUzA7gtpCP7xCP6Mz5A6GpFRYjkhIqoACrkM73Spj61vB6Gesy0ysgswfMUpfLz1AnILOXEbkT5YToiIKpBvLTV2vtMew4NqAwDWnEhG9zlHcTrpvrTBiIwIywkRUQWzspDjk15NsGZUW7iprXDjXi76LYzEzN2XkF9UInU8IoPHckJEVEmC6tXA3rCOeKNlLYgisOhIInrNPYbzNzOljkZk0FhOiIgqkdraArMGNMeSUH/UsFXiakYOXp8fgdnhVzj9PdFTsJwQEVWBl31cED6lI3o2c0OJVsRPB66iz7zjuJymkToakcFhOSEiqiIO1Szx81st8fNbfnCwsUBcqga95h7DvD8TUFzCT1GI/sZyQkRUxXo2c8e+KR3RtbELikpEfLcvHm8sjERCRo7U0YgMAssJEZEEnO2ssCS0FWb1bw47KwXOpWSix09HsfRoIhcRJLPHckJEJBFBEPBGq1r4Y0pHdKhfAwXFWnyx6xIGLolC8r1cqeMRSYblhIhIYm5qa/wyog2+fL0pbCzlOHn9Pl6dcwSro5K4iCCZJZYTIiIDIAgCQtp6Ye/kjmjrXR25hSX497ZYhC4/iduZeVLHI6pSLCdERAbE09EG60a3w396+kCpkOHo1bsI/uEINkan8FMUMhssJ0REBkYmEzCivTd2T+4AP097ZBcUY+qm8xi5KhppWflSxyOqdCwnREQGqq6TLTaODcAHrzaEpVyGg5cz8PIPh/kpCpk8lhMiIgOmkMvw9kv1sHNSezSvpUZ2/qNPUYavPIXULF6LQqaJ5YSIyAg0cLHD5vGB+KhbI1gqZDgUfwevzD6C304l81MUMjksJ0RERkIhl2Fcp7rYPam97lqUDzdfQOjyk7jFO3rIhLCcEBEZmXrOdtg0LhAfd29c6o6etSf4KQqZBpYTIiIjJJcJGN2xDnZP7oBWXg7IKSjG/229gCHLTiLlPmeXJePGckJEZMTqOtliw9gATO/pAysLGY4l3MWrPx7Br1FJXKOHjBbLCRGRkZPLBIxs7409kzuiTe3qeFhYgunbYhGy9AQ/RSGjxHJCRGQivGtUw/ox7fBpLx9YW8gRmXgPwT8ewaqIG/wUhYwKywkRkQmRyQQMC/LG3rAOujV6Pvn9IgYuiULSvYdSxyMqE5YTIiIT5OVYDetGt8OM3k10Kx0H/3gEy49d56coZPBYToiITJRMJmBIQG3sC+uIgDqOyC/S4vOdcRiwKBIJGTlSxyN6Kr3Lya1btzB48GA4OjrC2toavr6+iI6OfuYxhw4dQsuWLaFUKlGvXj2sXLmyvHmJiEhPHtVtsGZUW3zRpymqWcoRnfQA3eccxbw/E1BUopU6HtFj9ConDx48QFBQECwsLLBnzx7ExcVh1qxZcHBweOox169fR48ePdC5c2fExMQgLCwMo0aNwr59+144PBERlY1MJmBwOy/88W4nvNTQCYUlWny3Lx6v/XwcF25mSR2PqBRB1GM6wY8++gjHjx/H0aNHy/wCH374IXbt2oXY2FjdvoEDByIzMxN79+4t03NoNBqo1WpkZWVBpVKV+bWJiOhxoihie8xtfLbjIh7kFkEmAKM71sGUrg1gZSGXOh6ZkPL+/dbrk5Pff/8d/v7+6N+/P5ydneHn54clS5Y885jIyEh07dq11L7g4GBERkY+9ZiCggJoNJpSGxERVQxBENDHrybC3+2EXs3doRWBRYcT8eqPRxCVeE/qeET6lZPExEQsWLAA9evXx759+zB+/HhMmjQJq1ateuoxaWlpcHFxKbXPxcUFGo0GeXlPXqhq5syZUKvVus3Dw0OfmEREVAY1bJWYO8gPS0L94aJS4sa9XAxcHIX/23oBmvwiqeORGdOrnGi1WrRs2RJfffUV/Pz8MGbMGIwePRoLFy6s0FDTpk1DVlaWbktJSanQ5yciov962ccF4e92wlttPQEAa08k45XZR3DgUrrEychc6VVO3Nzc4OPjU2pf48aNkZyc/NRjXF1dkZ5e+h/w9PR0qFQqWFtbP/EYpVIJlUpVaiMiosqjsrLAV6/7Yt3odqjtaIM0TT5GrorGpHVncS+nQOp4ZGb0KidBQUGIj48vte/KlSvw8vJ66jEBAQE4cOBAqX3h4eEICAjQ56WJiKgKBNR1xJ7JHTG2Yx3IBOD3c7fRdfZhbDt7C3rcP0H0QvQqJ1OmTEFUVBS++uorJCQkYO3atVi8eDEmTJigGzNt2jSEhobqfh43bhwSExPxwQcf4PLly5g/fz42bNiAKVOmVNy7ICKiCmNtKce07o2xbUIQGrna4UFuEcJ+i8GIladwO/PJ1woSVSS9yknr1q2xdetWrFu3Dk2bNsWMGTPw448/IiQkRDcmNTW11Nc83t7e2LVrF8LDw9G8eXPMmjULS5cuRXBwcMW9CyIiqnDNatljxzvt8f4rDWApl+HP+Dt45Ycj+DUqiVPgU6XSa54TqXCeEyIiaSVkZOODTedxJjkTANCmdnXMfMMXdZ1spQ1GBq1K5jkhIiLzVM/ZDhvHBeLTXj6PFhK8cR/d5hzF3ANXUVjMKfCpYrGcEBFRmchlAoYFeWNfWEd0bOCEwmItZoVfQY+fjuJ00n2p45EJYTkhIiK9eFS3warhrTFnYAs4VrPE1Ywc9FsYiX9v4+RtVDFYToiISG+CIKB3i5rY/24n9G9VC6IIrI5KxsuzD2NvbJrU8cjIsZwQEVG5OVSzxHf9m2Pt6Lao7WiDdE0Bxq0+jTG/RCMtK1/qeGSkWE6IiOiFBdatgb1hHTGxcz0oZAL+iEtH19mH8UvkDZTwtmPSE8sJERFVCCsLOd4PbohdkzrAz9MeOQXF+M/2i+i3MAKX07i6PJUdywkREVWohq522DQuEJ/3bgJbpQJnkzPR86dj+H5fPPKLSqSOR0aA5YSIiCqcXCYgNKA2wt/tiJd9XFCsFfHznwnoNucoIq7dlToeGTiWEyIiqjRuamssCfXHwsGt4GynxPW7D/HWkhOYuvEcHjwslDoeGSiWEyIiqnSvNnXF/vc6YUg7LwgCsPH0TXSdfRjbY7jaMT2O5YSIiKqEysoCM/o0xaZxAWjgYot7DwsxeX0Mhq44haR7D6WORwaE5YSIiKpUK6/q2PlOB7z3cgNYKmQ4cuXRasdzD1xFQTEvmCWWEyIikoClQoZ3utTHvrCOCKrniIK/1unpPucoohLvSR2PJMZyQkREkvGuUQ2rR7bFnIEtUMPWEtfuPMTAxVF4b8M53MspkDoeSYTlhIiIJPX3Oj0H3n0JIW09IQjA5jM30WX2Yfx2KhlazjBrdgTRCC6T1mg0UKvVyMrKgkqlkjoOERFVojPJD/Dx1lhcSn00q2zr2g748nVfNHCxkzgZ6au8f7/5yQkRERmUlp4O2DExCB93bwwbSzlO3XiA7nOO4pu9l5FXyAtmzQHLCRERGRyFXIbRHesg/N1OeOWvGWYXHLqGl384jIOX06WOR5WM5YSIiAxWTXtrLA71x5JQf7irrXDzQR5GrIzG+NWnkZqVJ3U8qiQsJ0REZPBe9nFB+LudMKZjHchlAvbEpqHrrMNYfuw6iku0UsejCsZyQkRERqGaUoH/694YOya2h5+nPR4WluDznXHoM/84zqVkSh2PKhDLCRERGRUfdxU2jwvEl683hcpKgdhbGvSZfxz/2R6LrLwiqeNRBWA5ISIioyOTCQhp64UD772EPi3cIYrAL5FJ6DLrMLacucnFBI0cywkRERktJzslfhzohzWj2qKOUzXczSnAuxvO4c1FUYhPy5Y6HpUTywkRERm9oHo1sHdyR3zwakNYW8hx8sZ9dP/pKL7YGYecgmKp45GeWE6IiMgkWCpkePuletj/XicEN3FBiVbE0mPX0WXWIew4d5tf9RgRlhMiIjIpNe2tsWiIP1YMbw0vRxukawrwzrqzGLLsJBIycqSOR2XAckJERCapc0Nn7AvriLCu9WGpkOFYwl10m3ME3+y9jNxCftVjyFhOiIjIZFlZyBHWtQH2T+mEfzVyRlHJX9Pgzz6CvbFp/KrHQLGcEBGRyfN0tMGyof5YPKQVatpb41ZmHsatPo3hK08h6d5DqePR/2A5ISIisyAIAl5p4or973bChM51YSEXcCj+Dl7+4Qh+CL+C/CKueGwoWE6IiMisWFvKMTW4EfaGdUT7ejVQWKzFnANX8coPR7jisYFgOSEiIrNU18kWv45sg3lvtYSrygrJ93MxYmU0Rv8SjZT7uVLHM2ssJ0REZLYEQUCPZm7Y/96jFY8VMgHhcenoOvsw5uy/yq96JCKIRnCpskajgVqtRlZWFlQqldRxiIjIRF1Jz8b0bbE4cf0+AKCWgzWm9/TBKz4uEARB4nTGp7x/v1lOiIiI/kEURew8n4qvdl9CalY+AKBD/Rr49LUmqOtkK3E648JyQkREVIFyC4sx788ELDlyHYUlWljIBYwI8sY7XerDVqmQOp5RYDkhIiKqBDfuPsTnO+Nw8HIGAMDZTon/694YvVu486ue52A5ISIiqkQHLqXj851xSLr36E6e1rUd8OlrTdDEXS1xMsPFckJERFTJ8otKsOzYdfx8MAF5RSWQCUBIWy+890oD2NtYSh3P4LCcEBERVZHbmXn4avcl7DyfCgBwsLHA+8ENMbC1J+QyftXzN5YTIiKiKhZx7S4++z0O8enZAICmNVX47LWmaOXlIHEyw8ByQkREJIHiEi1+jUrC7PAryM4vBgC80bIWPuzWEM52VhKnkxbLCRERkYTu5hTgu73x+C06BQBgq1QgrGt9hAbUhqXCPCdkZzkhIiIyADEpmfhkeyzO3cwCANRxqobpPX3QuaGzxMmqHssJERGRgdBqRWw6fRPf7ruMuzmFAIDODZ3w754+ZjXLLMsJERGRgdHkF+HngwlYcfw6ikpEKGQChgbWxqQu9aG2tpA6XqVjOSEiIjJQ1+8+xJe74rD/0qNZZqtXs8T7rzTEm609TPrWY5YTIiIiA3f4yh3M2BmHhIwcAEBjNxU+6eWDdnUcJU5WOVhOiIiIjEBRiRZr/rr1WPPXrcfdfV0xrVtjeFS3kThdxWI5ISIiMiL3Hxbih/ArWHMiCVoRsFTIMLZjHYx/qS5sLE1j1WOWEyIiIiN0OU2Dz3fEIeLaPQCAq8oKH3ZriN7Na0Jm5NejsJwQEREZKVEU8UdcOr7YFYeU+3kAAD9Pe3zSqwlaeNhLG+4FlPfvt15T1n366acQBKHU1qhRo6eOX7ly5WPjrazMeypfIiKi/yUIAoKbuCJ8Sid88GpD2FjKcTY5E33mHcd7G84hXZMvdcQqpfeXWk2aNMH+/fv/+wSKZz+FSqVCfHy87mdBMO6PqIiIiCqLlYUcb79UD2+0rIVv98Zj85mb2HzmJvbEpmJC53oY2d4bVhZyqWNWOr3LiUKhgKura5nHC4Kg13giIiJz56KywqwBzTEkwAuf7biIs8mZ+G5fPNaeSMaH3RqhVzM3k/4/+3qvRHT16lW4u7ujTp06CAkJQXJy8jPH5+TkwMvLCx4eHujduzcuXrz43NcoKCiARqMptREREZmbFh722DI+ED++2QJuaivcyszDpHVn0XdBBM4kP5A6XqXR64LYPXv2ICcnBw0bNkRqaio+++wz3Lp1C7GxsbCzs3tsfGRkJK5evYpmzZohKysL33//PY4cOYKLFy+iVq1aT32dTz/9FJ999tlj+3lBLBERmau8whIsPZqIBYevIbewBADwWnN3fPBqQ9RyMMz5USS5WyczMxNeXl6YPXs2Ro4c+dzxRUVFaNy4MQYNGoQZM2Y8dVxBQQEKCgp0P2s0Gnh4eLCcEBGR2UvX5GPWH/HYePomxL/mRxnV3hvjX6oLOyvDWq+nSu7W+V/29vZo0KABEhISyjTewsICfn5+zx2vVCqhUqlKbURERPToepRv+zXHznfaI6COIwqLtZh/6Bo6f38Ia08ko0Rr8DOEPNcLlZOcnBxcu3YNbm5uZRpfUlKCCxculHk8ERERPVkTdzXWjm6LJaH+8K5RDXdzCvF/Wy+gx09HcfTqHanjvRC9ysn777+Pw4cP48aNG4iIiMDrr78OuVyOQYMGAQBCQ0Mxbdo03fjPP/8cf/zxBxITE3HmzBkMHjwYSUlJGDVqVMW+CyIiIjMkCAJe9nHBvrCO+E9PH6itLXA5LRtDlp3E8BUnkZCRLXXEctHrVuKbN29i0KBBuHfvHpycnNC+fXtERUXByckJAJCcnAyZ7L9958GDBxg9ejTS0tLg4OCAVq1aISIiAj4+PhX7LoiIiMyYpUKGEe290bdlTfx0IAG/RN7An/F3cOTqXYS09URY1waoXs1S6phlxunriYiITEzinRzM3HMZ4XHpAAA7KwUm/as+QgO9oFRU3SRuXFuHiIiISom4dhdf7LyEuNRH84V5VrfBtG6N8GpT1yqZxI3lhIiIiB5TohWx+cxNfL8vHhnZj6bpaFO7Oj7u0RjNK3lRQZYTIiIieqqHBcVYdCQRi49cQ36RFsCjSdymBjeER/XKmcSN5YSIiIieKzUrD9/vu4ItZ/87idvwoNp4+6V6UFtX7CRukkzCRkRERMbFTW2NWQOaY8fE9gis+2gSt0WHE3EoPkPqaDp6r0pMRERExq9pTTXWjGqLP+MzsD3mNno1c5c6kg7LCRERkZkSBAH/auSCfzVykTpKKfxah4iIiAwKywkREREZFJYTIiIiMigsJ0RERGRQWE6IiIjIoLCcEBERkUFhOSEiIiKDwnJCREREBoXlhIiIiAwKywkREREZFJYTIiIiMigsJ0RERGRQWE6IiIjIoBjFqsSiKAIANBqNxEmIiIiorP7+u/333/GyMopykp2dDQDw8PCQOAkRERHpKzs7G2q1uszjBVHfOiMBrVaL27dvw87ODoIgVNjzajQaeHh4ICUlBSqVqsKelx7Hc101eJ6rBs9z1eB5rjqVda5FUUR2djbc3d0hk5X9ShKj+OREJpOhVq1alfb8KpWK/+BXEZ7rqsHzXDV4nqsGz3PVqYxzrc8nJn/jBbFERERkUFhOiIiIyKCYdTlRKpX45JNPoFQqpY5i8niuqwbPc9Xgea4aPM9Vx9DOtVFcEEtERETmw6w/OSEiIiLDw3JCREREBoXlhIiIiAwKywkREREZFLMuJ/PmzUPt2rVhZWWFtm3b4uTJk1JHMhgzZ85E69atYWdnB2dnZ/Tp0wfx8fGlxuTn52PChAlwdHSEra0t3njjDaSnp5cak5ycjB49esDGxgbOzs6YOnUqiouLS405dOgQWrZsCaVSiXr16mHlypWP5TGX39XXX38NQRAQFham28fzXDFu3bqFwYMHw9HREdbW1vD19UV0dLTucVEU8Z///Adubm6wtrZG165dcfXq1VLPcf/+fYSEhEClUsHe3h4jR45ETk5OqTHnz59Hhw4dYGVlBQ8PD3z77bePZdm4cSMaNWoEKysr+Pr6Yvfu3ZXzpqtYSUkJpk+fDm9vb1hbW6Nu3bqYMWNGqXVVeJ7L58iRI+jVqxfc3d0hCAK2bdtW6nFDOq9lyfJcoplav369aGlpKS5fvly8ePGiOHr0aNHe3l5MT0+XOppBCA4OFlesWCHGxsaKMTExYvfu3UVPT08xJydHN2bcuHGih4eHeODAATE6Olps166dGBgYqHu8uLhYbNq0qdi1a1fx7Nmz4u7du8UaNWqI06ZN041JTEwUbWxsxHfffVeMi4sT586dK8rlcnHv3r26Mebyuzp58qRYu3ZtsVmzZuLkyZN1+3meX9z9+/dFLy8vcdiwYeKJEyfExMREcd++fWJCQoJuzNdffy2q1Wpx27Zt4rlz58TXXntN9Pb2FvPy8nRjXn31VbF58+ZiVFSUePToUbFevXrioEGDdI9nZWWJLi4uYkhIiBgbGyuuW7dOtLa2FhctWqQbc/z4cVEul4vffvutGBcXJ/773/8WLSwsxAsXLlTNyahEX375pejo6Cju3LlTvH79urhx40bR1tZWnDNnjm4Mz3P57N69W/z444/FLVu2iADErVu3lnrckM5rWbI8j9mWkzZt2ogTJkzQ/VxSUiK6u7uLM2fOlDCV4crIyBABiIcPHxZFURQzMzNFCwsLcePGjboxly5dEgGIkZGRoig++pdJJpOJaWlpujELFiwQVSqVWFBQIIqiKH7wwQdikyZNSr3Wm2++KQYHB+t+NoffVXZ2tli/fn0xPDxc7NSpk66c8DxXjA8//FBs3779Ux/XarWiq6ur+N133+n2ZWZmikqlUly3bp0oiqIYFxcnAhBPnTqlG7Nnzx5REATx1q1boiiK4vz580UHBwfdef/7tRs2bKj7ecCAAWKPHj1KvX7btm3FsWPHvtibNAA9evQQR4wYUWpf3759xZCQEFEUeZ4ryv+WE0M6r2XJUhZm+bVOYWEhTp8+ja5du+r2yWQydO3aFZGRkRImM1xZWVkAgOrVqwMATp8+jaKiolLnsFGjRvD09NSdw8jISPj6+sLFxUU3Jjg4GBqNBhcvXtSN+edz/D3m7+cwl9/VhAkT0KNHj8fOBc9zxfj999/h7++P/v37w9nZGX5+fliyZInu8evXryMtLa3U+1er1Wjbtm2p82xvbw9/f3/dmK5du0Imk+HEiRO6MR07doSlpaVuTHBwMOLj4/HgwQPdmGf9LoxZYGAgDhw4gCtXrgAAzp07h2PHjqFbt24AeJ4riyGd17JkKQuzLCd3795FSUlJqf+YA4CLiwvS0tIkSmW4tFotwsLCEBQUhKZNmwIA0tLSYGlpCXt7+1Jj/3kO09LSnniO/37sWWM0Gg3y8vLM4ne1fv16nDlzBjNnznzsMZ7nipGYmIgFCxagfv362LdvH8aPH49JkyZh1apVAP57np71/tPS0uDs7FzqcYVCgerVq1fI78IUzvNHH32EgQMHolGjRrCwsICfnx/CwsIQEhICgOe5shjSeS1LlrIwilWJSVoTJkxAbGwsjh07JnUUk5OSkoLJkycjPDwcVlZWUscxWVqtFv7+/vjqq68AAH5+foiNjcXChQsxdOhQidOZjg0bNmDNmjVYu3YtmjRpgpiYGISFhcHd3Z3nmfRilp+c1KhRA3K5/LE7HtLT0+Hq6ipRKsM0ceJE7Ny5E3/++Sdq1aql2+/q6orCwkJkZmaWGv/Pc+jq6vrEc/z3Y88ao1KpYG1tbfK/q9OnTyMjIwMtW7aEQqGAQqHA4cOH8dNPP0GhUMDFxYXnuQK4ubnBx8en1L7GjRsjOTkZwH/P07Pev6urKzIyMko9XlxcjPv371fI78IUzvPUqVN1n574+vpiyJAhmDJliu5TQZ7nymFI57UsWcrCLMuJpaUlWrVqhQMHDuj2abVaHDhwAAEBARImMxyiKGLixInYunUrDh48CG9v71KPt2rVChYWFqXOYXx8PJKTk3XnMCAgABcuXCj1L0R4eDhUKpXuD0VAQECp5/h7zN/PYeq/qy5duuDChQuIiYnRbf7+/ggJCdH9b57nFxcUFPTYrfBXrlyBl5cXAMDb2xuurq6l3r9Go8GJEydKnefMzEycPn1aN+bgwYPQarVo27atbsyRI0dQVFSkGxMeHo6GDRvCwcFBN+ZZvwtjlpubC5ms9J8VuVwOrVYLgOe5shjSeS1LljIp86WzJmb9+vWiUqkUV65cKcbFxYljxowR7e3tS93xYM7Gjx8vqtVq8dChQ2Jqaqpuy83N1Y0ZN26c6OnpKR48eFCMjo4WAwICxICAAN3jf9/i+sorr4gxMTHi3r17RScnpyfe4jp16lTx0qVL4rx58554i6s5/a7+ebeOKPI8V4STJ0+KCoVC/PLLL8WrV6+Ka9asEW1sbMTVq1frxnz99deivb29uH37dvH8+fNi7969n3grpp+fn3jixAnx2LFjYv369UvdipmZmSm6uLiIQ4YMEWNjY8X169eLNjY2j92KqVAoxO+//168dOmS+Mknnxj1La7/NHToULFmzZq6W4m3bNki1qhRQ/zggw90Y3ieyyc7O1s8e/asePbsWRGAOHv2bPHs2bNiUlKSKIqGdV7LkuV5zLaciKIozp07V/T09BQtLS3FNm3aiFFRUVJHMhgAnritWLFCNyYvL098++23RQcHB9HGxkZ8/fXXxdTU1FLPc+PGDbFbt26itbW1WKNGDfG9994Ti4qKSo35888/xRYtWoiWlpZinTp1Sr3G38zpd/W/5YTnuWLs2LFDbNq0qahUKsVGjRqJixcvLvW4VqsVp0+fLrq4uIhKpVLs0qWLGB8fX2rMvXv3xEGDBom2traiSqUShw8fLmZnZ5cac+7cObF9+/aiUqkUa9asKX799dePZdmwYYPYoEED0dLSUmzSpIm4a9euin/DEtBoNOLkyZNFT09P0crKSqxTp4748ccfl7o1lee5fP78888n/jd56NChoiga1nktS5bnEUTxH1P3EREREUnMLK85ISIiIsPFckJEREQGheWEiIiIDArLCRERERkUlhMiIiIyKCwnREREZFBYToiIiMigsJwQERGRQWE5ISIiIoPCckJEREQGheWEiIiIDArLCRERERmU/wf2G5JrgT8YPAAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sns.lineplot(x=J_all[0], y=J_all[1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "id": "31574b9a",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "16a65185cb4541ffbdb5a02b33027697",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "  0%|          | 0/1000000 [00:00<?, ?it/s]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Die 3 Parameter der linearen Regression:\n",
      "[  81.08915614  248.45370379 6493.32860783]\n",
      "Kostenfunktion J: 201611248738.63248\n",
      "J relativ zu Startkosten: 0.37282268754553793\n",
      "Vergleich Kostenfunktion zu analytischer Lösung: 1.36539*J_ana\n",
      "Relative Abweichung der Parameter zu analytischer Lösung: [1.89594048e-03 1.06809958e-03 2.67088270e+00]*w_ana\n",
      "CPU times: user 3min 26s, sys: 1min 26s, total: 4min 53s\n",
      "Wall time: 31.6 s\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "alpha = 3.1e-10  # verschiedene alpha ausprobieren\n",
    "n_iterations = 1000000\n",
    "w_gd2, J_all2 = grad_desc(w=w_gd, alpha=alpha, x=X, y=y, n_iterations=n_iterations)\n",
    "print('Die {} Parameter der linearen Regression:\\n{}'.format(len(w_gd2), w_gd2))\n",
    "print('Kostenfunktion J: {}'.format(J_all2[1][-1]))\n",
    "print('J relativ zu Startkosten: {}'.format(J_all2[1][-1]/J_all2[1][0]))\n",
    "print('Vergleich Kostenfunktion zu analytischer Lösung: {:.5f}*J_ana'.format(J_all2[1][-1]/J_ana))\n",
    "print('Relative Abweichung der Parameter zu analytischer Lösung: {}*w_ana'.format((w_gd2)/w_ana))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "4434e050",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<Axes: >"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiMAAAG+CAYAAABBOgSxAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA4xElEQVR4nO3deXxU5f3//ffMJDPZBwLZgAiyr2FzISh1o0WlFtp+rUUU2y/a2kKL7U+/LdXeLmiDrXjXb7UobtRbMRZvwRYRtShQDVhWjYDIngBZEMhMEmCSzJzfH0kGAknIJDNzMsnr+XicR5yT68z5zIGH8+a6rnMdi2EYhgAAAExiNbsAAADQuRFGAACAqQgjAADAVIQRAABgKsIIAAAwFWEEAACYijACAABMRRgBAACmIowAAABTEUYAAICpIiqMrFu3TjfddJN69Oghi8Wi5cuXB3T86dOn9aMf/UgjRoxQVFSUpk6del6boqIi3XrrrRo4cKCsVqvuueeeoNQOAAAaF1FhpLKyUiNHjtQzzzzTquO9Xq9iY2P1y1/+UhMnTmy0jcfjUUpKih544AGNHDmyLeUCAIAWiDK7gEDccMMNuuGGG5r8vcfj0f3336/XX39dZWVlGj58uB5//HFdffXVkqT4+HgtXLhQkvTJJ5+orKzsvPfo06ePnnrqKUnSSy+9FPTPAAAAGoqonpELmT17ttavX6/c3Fx9/vnnuvnmm3X99ddr9+7dZpcGAACa0GHCSEFBgV5++WUtXbpUEyZMUL9+/XTvvffqyiuv1Msvv2x2eQAAoAkRNUzTnPz8fHm9Xg0cOLDBfo/Ho27duplUFQAAuJAOE0YqKipks9m0efNm2Wy2Br9LSEgwqSoAAHAhHSaMjB49Wl6vV6WlpZowYYLZ5QAAgBaKqDBSUVGhPXv2+F/v379f27ZtU3JysgYOHKjp06drxowZWrBggUaPHq2jR49q9erVysrK0uTJkyVJO3bsUFVVlY4fP67y8nJt27ZNkjRq1Cj/+9bvq6io0NGjR7Vt2zbZ7XYNHTo0XB8VAIBOw2IYhmF2ES21Zs0aXXPNNeftv+OOO7R48WJVV1fr0Ucf1SuvvKLDhw+re/fuGjdunB5++GGNGDFCUu2tuwcPHjzvPc6+DBaL5bzf9+7dWwcOHAjehwEAAJIiLIwAAICOp8Pc2gsAACITYQQAAJgqIiaw+nw+HTlyRImJiY3O5wAAAO2PYRgqLy9Xjx49ZLU23f8REWHkyJEjyszMNLsMAADQCoWFherVq1eTv4+IMJKYmCip9sMkJSWZXA0AAGgJt9utzMxM//d4UyIijNQPzSQlJRFGAACIMBeaYhHQBNaHHnpIFoulwTZ48OAm2y9evPi89jExMYGcEgAAdHAB94wMGzZM//rXv868QVTzb5GUlKRdu3b5XzMBFQAAnC3gMBIVFaX09PQWt7dYLAG1BwAAnUvA64zs3r1bPXr0UN++fTV9+nQVFBQ0276iokK9e/dWZmampkyZou3bt1/wHB6PR263u8EGAAA6poDCyOWXX67Fixdr1apVWrhwofbv368JEyaovLy80faDBg3SSy+9pLfffluvvvqqfD6fxo8fr0OHDjV7npycHDmdTv/Gbb0AAHRcbXo2TVlZmXr37q0nn3xSM2fOvGD76upqDRkyRNOmTdO8efOabOfxeOTxePyv628Ncrlc3E0DAECEcLvdcjqdF/z+btOtvV26dNHAgQO1Z8+eFrWPjo7W6NGjL9je4XDI4XC0pTQAABAh2vRsmoqKCu3du1cZGRktau/1epWfn9/i9gAAoOMLKIzce++9Wrt2rQ4cOKC8vDx997vflc1m07Rp0yRJM2bM0Ny5c/3tH3nkEb3//vvat2+ftmzZottuu00HDx7UnXfeGdxPAQAAIlZAwzSHDh3StGnTdOzYMaWkpOjKK6/Uhg0blJKSIkkqKCho8CCcEydO6K677lJxcbG6du2qsWPHKi8vT0OHDg3upwAAABGrTRNYw6WlE2AAAED70dLv7zbNGQEAAGirTh1GfvH6Vt38bJ6KXafNLgUAgE4rIp7aGyob9x9Xsfu0SstPK93JA/wAADBDp+4Z6RpvlyQdr6wyuRIAADqvTh1GutWFkRMnCSMAAJilU4eR+p6RYxWEEQAAzNKpwwg9IwAAmK9Th5GuccwZAQDAbJ06jCQnEEYAADBb5w4j9IwAAGC6zh1GuLUXAADTEUZEGAEAwEyEEUllp6rl9bX75wUCANAhdeow0iUuWpJkGFIZt/cCAGCKTh1Gom1WOWNrAwlrjQAAYI5OHUakM0M1rMIKAIA5CCOswgoAgKk6fRipX4X1GHfUAABgik4fRvzPpyGMAABgik4fRvxP7iWMAABgik4fRugZAQDAXJ0+jNAzAgCAuTp9GOnG3TQAAJiq04eRrv5hmmqTKwEAoHPq9GGkm3+YxmNyJQAAdE6dPozU94ycrvbpVJXX5GoAAOh8On0YibfbZI+qvQz0jgAAEH6dPoxYLBYlxzFvBAAAs3T6MCKd9bA8ekYAAAg7woh4WB4AAGYijOisnpEKwggAAOFGGNGZMHKcVVgBAAg7woik7gm1YeTrCuaMAAAQboQRSd0THJIYpgEAwAyEEUnd6sIIPSMAAIQfYURnD9PQMwIAQLgRRnTWME2lR4ZhmFwNAACdC2FEUreEM8+nqeT5NAAAhBVhRFKcPUpxdpsk6RjzRgAACKuAwshDDz0ki8XSYBs8eHCzxyxdulSDBw9WTEyMRowYoZUrV7ap4FDpxu29AACYIuCekWHDhqmoqMi/ffzxx022zcvL07Rp0zRz5kxt3bpVU6dO1dSpU/XFF1+0qehQ6O6/o4ZJrAAAhFPAYSQqKkrp6en+rXv37k22feqpp3T99dfrvvvu05AhQzRv3jyNGTNGTz/9dJuKDoVu8dzeCwCAGQIOI7t371aPHj3Ut29fTZ8+XQUFBU22Xb9+vSZOnNhg36RJk7R+/fpmz+HxeOR2uxtsoZaSyPNpAAAwQ0Bh5PLLL9fixYu1atUqLVy4UPv379eECRNUXl7eaPvi4mKlpaU12JeWlqbi4uJmz5OTkyOn0+nfMjMzAymzVegZAQDAHAGFkRtuuEE333yzsrKyNGnSJK1cuVJlZWX6+9//HtSi5s6dK5fL5d8KCwuD+v6NqZ/ASs8IAADhFdWWg7t06aKBAwdqz549jf4+PT1dJSUlDfaVlJQoPT292fd1OBxyOBxtKS1g3VkSHgAAU7RpnZGKigrt3btXGRkZjf4+Oztbq1evbrDvgw8+UHZ2dltOGxLc2gsAgDkCCiP33nuv1q5dqwMHDigvL0/f/e53ZbPZNG3aNEnSjBkzNHfuXH/7OXPmaNWqVVqwYIG+/PJLPfTQQ9q0aZNmz54d3E8RBCn+JeEZpgEAIJwCGqY5dOiQpk2bpmPHjiklJUVXXnmlNmzYoJSUFElSQUGBrNYz+Wb8+PFasmSJHnjgAf3ud7/TgAEDtHz5cg0fPjy4nyII6p/cW3ayWtVen6JtLE4LAEA4WIwIeDKc2+2W0+mUy+VSUlJSSM7h8xka8MC78voMffq765SWFBOS8wAA0Fm09Pubf/7XsVotSo6vnTdytJx5IwAAhAth5Czd6sII80YAAAgfwshZUhLrbu+lZwQAgLAhjJzlTM8IYQQAgHAhjJylfuEzVmEFACB8CCNnqb+99ygLnwEAEDaEkbN05/k0AACEHWHkLDyfBgCA8COMnIXn0wAAEH6EkbPU39p7rKJKPl+7X5gWAIAOgTBylvphmhqfoRMnmTcCAEA4EEbOEm2z+peEL2XhMwAAwoIwco7UuqEank8DAEB4EEbOUT9vhJ4RAADCgzByjtTEGElSaflpkysBAKBzIIycw98z4qZnBACAcCCMnMM/Z4S1RgAACAvCyDlSk+rCCD0jAACEBWHkHCkJ9RNYmTMCAEA4EEbOkZpUO4GVW3sBAAgPwsg56ueMVFZ5VempMbkaAAA6PsLIOeIdUYqz2ySx1ggAAOFAGGkEq7ACABA+hJFGsPAZAADhQxhpREoSC58BABAuhJFG1N/ey8JnAACEHmGkEan0jAAAEDaEkUYwZwQAgPAhjDQihbtpAAAIG8JII7i1FwCA8CGMNKI+jByrrFK112dyNQAAdGyEkUZ0jbMrymqRJB2rqDK5GgAAOjbCSCOsVou68/ReAADCgjDSBG7vBQAgPAgjTaifN1JCzwgAACFFGGlCWlLtWiMlLsIIAAChRBhpQnpdGCl2E0YAAAglwkgT0pz1YYQ5IwAAhFKbwsj8+fNlsVh0zz33NNlm8eLFslgsDbaYmJi2nDYsMpwM0wAAEA5RrT1w48aNeu6555SVlXXBtklJSdq1a5f/tcViae1pw6Z+mKbIdcrkSgAA6Nha1TNSUVGh6dOn6/nnn1fXrl0v2N5isSg9Pd2/paWltea0YVU/TOM+XaNTVV6TqwEAoONqVRiZNWuWJk+erIkTJ7aofUVFhXr37q3MzExNmTJF27dvb7a9x+OR2+1usIVboiNKcXabJCaxAgAQSgGHkdzcXG3ZskU5OTktaj9o0CC99NJLevvtt/Xqq6/K5/Np/PjxOnToUJPH5OTkyOl0+rfMzMxAy2wzi8Vy5o4a5o0AABAyAYWRwsJCzZkzR6+99lqLJ6FmZ2drxowZGjVqlK666iq99dZbSklJ0XPPPdfkMXPnzpXL5fJvhYWFgZQZNP61RugZAQAgZAKawLp582aVlpZqzJgx/n1er1fr1q3T008/LY/HI5vN1ux7REdHa/To0dqzZ0+TbRwOhxwORyClhUS6k7VGAAAItYDCyHXXXaf8/PwG+3784x9r8ODB+s1vfnPBICLVhpf8/HzdeOONgVVqgjSGaQAACLmAwkhiYqKGDx/eYF98fLy6devm3z9jxgz17NnTP6fkkUce0bhx49S/f3+VlZXpT3/6kw4ePKg777wzSB8hdNLrHpZHGAEAIHRavc5IUwoKCmS1npmKcuLECd11110qLi5W165dNXbsWOXl5Wno0KHBPnXQMUwDAEDoWQzDMMwu4kLcbrecTqdcLpeSkpLCdt5thWWa+swnynDGaP3c68J2XgAAOoKWfn/zbJpm1N/aW1rukdfX7jMbAAARiTDSjO4Jdlktktdn6FgFD8wDACAUCCPNiLJZlZJYO4m1iEmsAACEBGHkAvyrsDKJFQCAkCCMXACrsAIAEFqEkQvw397LMA0AACFBGLmANIZpAAAIKcLIBfDkXgAAQoswcgEZXQgjAACEEmHkAnp2iZUkHS47pQhYrBYAgIhDGLmA+gmsnhqfTpysNrkaAAA6HsLIBTiibOqeULvw2ZGyUyZXAwBAx0MYaYGedfNGDhNGAAAIOsJIC/SomzdCzwgAAMFHGGkBwggAAKFDGGmBM2GE23sBAAg2wkgL9Ki7o+aIi54RAACCjTDSAgzTAAAQOoSRFqgPI6XlHlXV+EyuBgCAjoUw0gLd4u2yR1llGFIJD8wDACCoCCMtYLVa/PNGWGsEAIDgIoy0EPNGAAAIDcJIC9WHkSKe3gsAQFARRlqIYRoAAEKDMNJCDNMAABAahJEWIowAABAahJEWqg8jh0+ckmEYJlcDAEDHQRhpoR5daueMVFZ55T5dY3I1AAB0HISRFoqzR6lrXLQkhmoAAAgmwkgAzh6qAQAAwUEYCUBm1zhJUuGJkyZXAgBAx0EYCUBmcm3PSOFxekYAAAgWwkgAMpPpGQEAINgIIwHwD9McJ4wAABAshJEA9OrKWiMAAAQbYSQAvep6Rso9NXKdqja5GgAAOgbCSABi7TZ1T3BIYhIrAADBQhgJkP+OGiaxAgAQFISRADGJFQCA4GpTGJk/f74sFovuueeeZtstXbpUgwcPVkxMjEaMGKGVK1e25bSmqp/EeohVWAEACIpWh5GNGzfqueeeU1ZWVrPt8vLyNG3aNM2cOVNbt27V1KlTNXXqVH3xxRetPbWpWGsEAIDgalUYqaio0PTp0/X888+ra9euzbZ96qmndP311+u+++7TkCFDNG/ePI0ZM0ZPP/10qwo2G8M0AAAEV6vCyKxZszR58mRNnDjxgm3Xr19/XrtJkyZp/fr1TR7j8XjkdrsbbO1F/QTWQ6w1AgBAUEQFekBubq62bNmijRs3tqh9cXGx0tLSGuxLS0tTcXFxk8fk5OTo4YcfDrS0sMhwxspqkTw1Ph0t9yg1KcbskgAAiGgB9YwUFhZqzpw5eu211xQTE7ov4blz58rlcvm3wsLCkJ0rUPYoq9LrAkghk1gBAGizgHpGNm/erNLSUo0ZM8a/z+v1at26dXr66afl8Xhks9kaHJOenq6SkpIG+0pKSpSent7keRwOhxwORyClhVWv5DgdcZ3WoRMnNbZ383NmAABA8wLqGbnuuuuUn5+vbdu2+bdLLrlE06dP17Zt284LIpKUnZ2t1atXN9j3wQcfKDs7u22Vm4hJrAAABE9APSOJiYkaPnx4g33x8fHq1q2bf/+MGTPUs2dP5eTkSJLmzJmjq666SgsWLNDkyZOVm5urTZs2adGiRUH6COHnX4WVJeEBAGizoK/AWlBQoKKiIv/r8ePHa8mSJVq0aJFGjhypN998U8uXLz8v1ESSi+rWGimgZwQAgDazGBFwf6rb7ZbT6ZTL5VJSUpLZ5WjzwRP6/sI89XDGKG/udWaXAwBAu9TS72+eTdMKfbrV9owccZ3W6WqvydUAABDZCCOtkBxvV6KjdroNk1gBAGgbwkgrWCwW9ekeL0k6cIwwAgBAWxBGWql33VDNwWOVJlcCAEBkI4y0Up9u9T0jhBEAANqCMNJKZ3pGGKYBAKAtCCOtdGbOCD0jAAC0BWGklep7Rg6fOKWqGp/J1QAAELkII62UkuBQnN0mnyEdOsFQDQAArUUYaSWLxaLedZNYmTcCAEDrEUbaoH4lVuaNAADQeoSRNqBnBACAtiOMtAE9IwAAtB1hpA3oGQEAoO0II23Qp3ttz0jh8ZOq8XJ7LwAArUEYaYO0xBjFRFtV4zN06MQps8sBACAiEUbawGq16OLuCZKkfV9XmFwNAACRiTDSRv1SaueN7C1lEisAAK1BGGmjvim1PSN7j9IzAgBAaxBG2qi+Z2TfUXpGAABoDcJIG/WjZwQAgDYhjLRR37qekWOVVSo7WWVyNQAARB7CSBvF2aPUwxkjSdrLUA0AAAEjjARBv1SGagAAaC3CSBD07V53ey9hBACAgBFGgqC+Z4Q7agAACBxhJAj6dmeYBgCA1iKMBEG/1NphmoJjJ1XNA/MAAAgIYSQI0pNiFGe3qcZnqOD4SbPLAQAgohBGgsBisfjXG9lbylANAACBIIwEyZmVWJnECgBAIAgjQdK/LozsLi03uRIAACILYSRIBqQlSpJ2lzBMAwBAIAgjQTIovS6MlJbL6zNMrgYAgMhBGAmSi5Lj5Iiy6nS1T4XcUQMAQIsRRoLEZrWof91KrF+VMG8EAICWIowE0aC6eSOEEQAAWo4wEkQD/GGESawAALRUQGFk4cKFysrKUlJSkpKSkpSdna133323yfaLFy+WxWJpsMXExLS56PZqUDrDNAAABCoqkMa9evXS/PnzNWDAABmGob/97W+aMmWKtm7dqmHDhjV6TFJSknbt2uV/bbFY2lZxOzYgtbZnZO/RClV7fYq20fEEAMCFBBRGbrrppgavH3vsMS1cuFAbNmxoMoxYLBalp6e3vsII0rNLrOLtNlVWeXXwWKX614UTAADQtFb/093r9So3N1eVlZXKzs5usl1FRYV69+6tzMxMTZkyRdu3b7/ge3s8Hrnd7gZbJLBaLf55I7uKmTcCAEBLBBxG8vPzlZCQIIfDobvvvlvLli3T0KFDG207aNAgvfTSS3r77bf16quvyufzafz48Tp06FCz58jJyZHT6fRvmZmZgZZpmoFpzBsBACAQFsMwAloutKqqSgUFBXK5XHrzzTf1wgsvaO3atU0GkrNVV1dryJAhmjZtmubNm9dkO4/HI4/H43/tdruVmZkpl8ulpKSkQMoNuxf+vU+PvrNTNwxP18LbxppdDgAApnG73XI6nRf8/g5ozogk2e129e/fX5I0duxYbdy4UU899ZSee+65Cx4bHR2t0aNHa8+ePc22czgccjgcgZbWLgysH6ahZwQAgBZp8+0ePp+vQS9Gc7xer/Lz85WRkdHW07Zb9c+oOfB1pU5Xe02uBgCA9i+gnpG5c+fqhhtu0EUXXaTy8nItWbJEa9as0XvvvSdJmjFjhnr27KmcnBxJ0iOPPKJx48apf//+Kisr05/+9CcdPHhQd955Z/A/STuRmuhQt3i7jlVWaVdxuUZmdjG7JAAA2rWAwkhpaalmzJihoqIiOZ1OZWVl6b333tM3v/lNSVJBQYGs1jOdLSdOnNBdd92l4uJide3aVWPHjlVeXl6L5pdEKovFoqE9kvTv3V9rR5GbMAIAwAUEPIHVDC2dANNe5KzcqefW7dPt43pr3tThZpcDAIApWvr9zRKhITC0R+0F31kUGeujAABgJsJICAzNOBNGfL523/EEAICpCCMhcHH3eNmjrKqs8qrg+EmzywEAoF0jjIRAlM2qwXW3+O5gqAYAgGYRRkKkfqhmxxHCCAAAzSGMhEj9JFZ6RgAAaB5hJEToGQEAoGUIIyEyuC6MFLtP61hFy5bLBwCgMyKMhEiCI0p9usVJknYW8dA8AACaQhgJofp5I18ccZlcCQAA7RdhJIRG9OwiSfr8UJmpdQAA0J4RRkJoZC+nJOnzQ/SMAADQFMJICA3rWRtGDp04xSRWAACaQBgJIWdstPp2j5ckfX6Y3hEAABpDGAmxrLqhmnyGagAAaBRhJMSyenWRxCRWAACaQhgJsSwmsQIA0CzCSIgN6+GUzWpRablHxa7TZpcDAEC7QxgJsVi7TQNSEyRJnzFUAwDAeQgjYTCybt4Ik1gBADgfYSQMRtTNG6FnBACA8xFGwmBUZhdJ0meFZfL5DHOLAQCgnSGMhMHg9ETFRtvkPl2jvUcrzC4HAIB2hTASBlE2q793ZPPBE+YWAwBAO0MYCZMxvbtIIowAAHAuwkiYjO3dVZK0uYAwAgDA2QgjYTI6szaM7DtaqeOVVSZXAwBA+0EYCZOu8Xb1S6l9gu9WekcAAPAjjISRf6iGeSMAAPgRRsKoPoxsoWcEAAA/wkgYjbmoNox8VuhStddncjUAALQPhJEw6peSoKSYKJ2q9mpnkdvscgAAaBcII2FktVp0SZ9kSdJ/9h83uRoAANoHwkiYjetbG0Y27DtmciUAALQPhJEwG9e3m6TanhEvD80DAIAwEm5DM5KU4IiS+3QN80YAABBhJOyibFZd2qf2rhqGagAAIIyYon6oZsM+JrECAEAYMcGZeSPHmDcCAOj0AgojCxcuVFZWlpKSkpSUlKTs7Gy9++67zR6zdOlSDR48WDExMRoxYoRWrlzZpoI7gmE9zswb+bKYeSMAgM4toDDSq1cvzZ8/X5s3b9amTZt07bXXasqUKdq+fXuj7fPy8jRt2jTNnDlTW7du1dSpUzV16lR98cUXQSk+UkXZrLrEP2+EoRoAQOdmMQyjTeMEycnJ+tOf/qSZM2ee97tbbrlFlZWVWrFihX/fuHHjNGrUKD377LMtPofb7ZbT6ZTL5VJSUlJbym03nl27V/Pf/VITh6TphTsuMbscAACCrqXf362eM+L1epWbm6vKykplZ2c32mb9+vWaOHFig32TJk3S+vXrm31vj8cjt9vdYOtoxvernTfy6b5jPKcGANCpBRxG8vPzlZCQIIfDobvvvlvLli3T0KFDG21bXFystLS0BvvS0tJUXFzc7DlycnLkdDr9W2ZmZqBltnvDejjVNS5a5Z4afVZYZnY5AACYJuAwMmjQIG3btk2ffvqpfvazn+mOO+7Qjh07glrU3Llz5XK5/FthYWFQ3789sFktuqJ/d0nSuq+OmlwNAADmCTiM2O129e/fX2PHjlVOTo5Gjhypp556qtG26enpKikpabCvpKRE6enpzZ7D4XD479ip3zqibwxIkSSt2/21yZUAAGCeNq8z4vP55PF4Gv1ddna2Vq9e3WDfBx980OQck85mwsDanpHPD5Wp7GSVydUAAGCOgMLI3LlztW7dOh04cED5+fmaO3eu1qxZo+nTp0uSZsyYoblz5/rbz5kzR6tWrdKCBQv05Zdf6qGHHtKmTZs0e/bs4H6KCJXhjNWA1AT5DClvL0vDAwA6p4DCSGlpqWbMmKFBgwbpuuuu08aNG/Xee+/pm9/8piSpoKBARUVF/vbjx4/XkiVLtGjRIo0cOVJvvvmmli9fruHDhwf3U0SwCXVDNf/ezbwRAEDn1OZ1RsKhI64zUu+jXaX68csb1bNLrD7+zTWyWCxmlwQAQFCEfJ0RBMflFyfLbrPqcNkp7fu60uxyAAAIO8KIyeLsUbrs4mRJ0oc7S02uBgCA8COMtAPXDUmVJP1rZ8kFWgIA0PEQRtqBiUNqV6nddPAEt/gCADodwkg7kJkcp4FpCfL6DK3ZxV01AIDOhTDSTtT3jjBUAwDobAgj7cR1dWFk7VdHVVXDU3wBAJ0HYaSdGJXZRd3i7So/XaONB46bXQ4AAGFDGGknbFaLrh3MXTUAgM6HMNKOTBxaO1Tz/vYSRcDCuAAABAVhpB25amCK4u02HS47pW2FZWaXAwBAWBBG2pGYaJuurZvIujK/6AKtAQDoGAgj7czkERmSpJX5xQzVAAA6BcJIO3P1IIZqAACdC2GknWGoBgDQ2RBG2qHJI9IlMVQDAOgcCCPt0NWDUhVXN1SzlaEaAEAHRxhph2KibZo0rLZ3ZNmWwyZXAwBAaBFG2qnvjekpSfrn50d4Vg0AoEMjjLRT4/t1V1qSQ2Unq/XRrlKzywEAIGQII+2UzWrR1FG1vSNvbTlkcjUAAIQOYaQd+96YXpKkD78s1YnKKpOrAQAgNAgj7dig9EQNzUhStdfQis+PmF0OAAAhQRhp5+onsr6xqdDkSgAACA3CSDv3vTG9ZLdZ9cVhtz4/VGZ2OQAABB1hpJ1LjrfrhroVWZd8WmByNQAABB9hJALcetlFkqR/fHZE5aerTa4GAIDgIoxEgMsuTlb/1ASdrPJq+TYmsgIAOhbCSASwWCz+3pHXNhzk4XkAgA6FMBIhvj+mlxxRVn1ZXK7NB0+YXQ4AAEFDGIkQzrhofXd07W2+L3683+RqAAAIHsJIBPnvKy+WJL23vVgFx06aXA0AAMFBGIkgA9MSddXAFPkM6aVP6B0BAHQMhJEIc+eE2t6Rv28qlOsUt/kCACIfYSTCXNm/uwalJepklVev/4dF0AAAkY8wEmEsFotm1vWOvPjxfp2u9ppcEQAAbUMYiUBTR/VUzy6xOlru0RsbeYAeACCyEUYikD3Kqruv7idJWrhmrzw19I4AACIXYSRC/eCSXkpPilGx+7SWbjpkdjkAALRaQGEkJydHl156qRITE5WamqqpU6dq165dzR6zePFiWSyWBltMTEybiobkiLLp7qv6SqrtHamq8ZlcEQAArRNQGFm7dq1mzZqlDRs26IMPPlB1dbW+9a1vqbKystnjkpKSVFRU5N8OHjzYpqJR64eXXaSURIcOl53SG5uYOwIAiExRgTRetWpVg9eLFy9WamqqNm/erG984xtNHmexWJSent7i83g8Hnk8Hv9rt9sdSJmdRky0TbOv6a8H/7FdT/1rt743uqfiHQH9kQIAYLo2zRlxuVySpOTk5GbbVVRUqHfv3srMzNSUKVO0ffv2Ztvn5OTI6XT6t8zMzLaU2aFNu+wi9ekWp68rPHr+3/vMLgcAgIBZjFY+j97n8+k73/mOysrK9PHHHzfZbv369dq9e7eysrLkcrn0xBNPaN26ddq+fbt69erV6DGN9YxkZmbK5XIpKSmpNeV2aO98XqRZS7Yozm7T2vuuUUqiw+ySAACQ2+2W0+m84Pd3q8PIz372M7377rv6+OOPmwwVjamurtaQIUM0bdo0zZs3r0XHtPTDdFaGYWjqX/P0WWGZbht3kR6dOsLskgAAaPH3d6uGaWbPnq0VK1boo48+CiiISFJ0dLRGjx6tPXv2tObUaITFYtHcGwZLkpZ8WqAdR5hjAwCIHAGFEcMwNHv2bC1btkwffvihLr744oBP6PV6lZ+fr4yMjICPRdPG9e2myVkZ8hnS//P2F/L5WtXhBQBA2AUURmbNmqVXX31VS5YsUWJiooqLi1VcXKxTp07528yYMUNz5871v37kkUf0/vvva9++fdqyZYtuu+02HTx4UHfeeWfwPgUkSQ9MHqI4u02bDp7QW1sPm10OAAAtElAYWbhwoVwul66++mplZGT4tzfeeMPfpqCgQEVFRf7XJ06c0F133aUhQ4boxhtvlNvtVl5enoYOHRq8TwFJUoYzVnOuGyBJylm5U65T1SZXBADAhbV6Ams4MYG15apqfLrxf/+tPaUV+uGlmZr//SyzSwIAdFIhncCK9sseZdVjU4dLknI3FmrdV0dNrggAgOYRRjqgy/t204/G95Ek/fb//1zlpxmuAQC0X4SRDup/rh+ki5LjdMR1Wn9Y+aXZ5QAA0CTCSAcVZ4/SH/+rdr7I6/8p0Ac7SkyuCACAxhFGOrBxfbvpzitr14K5783PdKTs1AWOAAAg/AgjHdz/XD9YWb2cKjtZrV++vlU1Xp/ZJQEA0ABhpIOzR1n1l2mjleiI0qaDJ/TkB1+ZXRIAAA0QRjqB3t3ilfP92ofn/XXNXr3zedEFjgAAIHwII53Et7N6+OeP/J+l2/TFYZfJFQEAUIsw0onMvXGIrhqYotPVPt31yiaVlp82uyQAAAgjnYnNatH/ThutvinxKnKd1szFm1ThqTG7LABAJ0cY6WScsdF68Y5LlRxvV/5hl37yyiZ5arxmlwUA6MQII53Qxd3jtfjHlyreblPe3mP61Rvb5PW1++clAgA6KMJIJ5XVq4sWzbhEdptVK/OL9X/+vo01SAAApiCMdGJX9O+u/502SlFWi5ZvO6J73timagIJACDMCCOd3PXDM/TM9DGKtlm04vMi/WLJVuaQAADCijACTRqWruduHyu7zapV24s148X/yHWy2uyyAACdBGEEkqRrB6fppR9dqgRHlD7df1zffzZPhcdPml0WAKATIIzA78oB3bX07mylJ8VoT2mFvvvXPG08cNzssgAAHRxhBA0MyUjS8llXaEhGkr6u8OiHizbohX/vk2Fw6y8AIDQIIzhPujNGb96dre+M7CGvz9Cj7+zU7CVb5T7NPBIAQPARRtCoeEeUnvrhKD1001BFWS16J79I1/+/6/TJnq/NLg0A0MEQRtAki8WiH11xsd74abZ6d4vTEddpTX/hUz38z+06WcUzbQAAwUEYwQWN7d1VK385QdMvv0iS9PInB/TNJ9fpve3FzCUBALQZYQQtEu+I0mPfHaHFP75UPbvE6nDZKf30/9usHy/eqP1fV5pdHgAgglmMCPinrdvtltPplMvlUlJSktnldHqnqrx6+qPdWrRun6q9hqKsFt1yaabmXDdAqUkxZpcHAGgnWvr9TRhBq+09WqF5K3Zoza6jkqSYaKt+fMXFmnnlxeqe4DC5OgCA2QgjCJsN+47pj6u+1JaCMkmSI8qqmy/ppbsm9FXvbvHmFgcAMA1hBGFlGIb+tbNUT3+0R58VlkmSrBbpm0PTdOvlvTWhf3dZrRZziwQAhBVhBKYwDEOf7j+u59bu1Ud1wzeS1KtrrH54aaa+P7aXMpyxJlYIAAgXwghM91VJuZZ8WqC3thyS+/SZdUku65Osb4/M0A3DM5SSyNwSAOioCCNoN05VebUyv0i5Gwu08cAJ/36rRbq0T7KuHZyqawanakBqgiwWhnIAoKMgjKBdOlJ2Sivzi/TPz47os0OuBr/r2SVWVw1K0bi+3XT5xclK4zZhAIhohBG0e4XHT+rDL0v10a5Srd97TJ4aX4Pf9+4Wp0v7JOvSPl01omcXDUhLULSNdfoAIFIQRhBRTlV5tWHfMa396qg2HjiunUVu+c75m2mPsmpIeqKG93RqeE+nBqYlqn9qgpyx0eYUDQBoFmEEEc19ulqbD57Qxv3HtaXghLYfdqvc0/jD+VISHeqXEq/+qQnql5Kg3t3i1LNLnHp1jVW8IyrMlQMA6hFG0KH4fIYKjp/UF0dcyj/s0vbDbu0prVCx+3Szx3WNi1avrnHq2SVWPbvGKjXRoZSztwSHusbZWQMFAEKAMIJOocJTo72lFdp7tEJ76n4eOnFKh06ckutUdYvew2a1qFu8Xcnxdjljo9UlLlpdYu1yxkXLGRvt3+eMjVZiTLTi7TbF2m2Kt0cpzmGT3WblLiAAaERLv78D6sPOycnRW2+9pS+//FKxsbEaP368Hn/8cQ0aNKjZ45YuXarf//73OnDggAYMGKDHH39cN954YyCnBhqV4IjSyMwuGpnZ5bzflZ+u1uGyUzp0/JQOnTipI67T+rrco6MVHh0tr92OVVbJ6zNUWu5RabmnVTVEWS2Ks9sU74jy/4yNtikm2iZ7lFWOKGvdT5scda/P3nd2m2ibVVFWi6LqftqsljM/bRbZrPW/r9/fsF2UzSqb1SKrRbJaLLJaLLLU/bfFv0+EJwDtSkBhZO3atZo1a5YuvfRS1dTU6He/+52+9a1vaceOHYqPb/wZJHl5eZo2bZpycnL07W9/W0uWLNHUqVO1ZcsWDR8+PCgfAmhMYky0BqdHa3B602m82uvT8coqlbo9OnGySmWnquU6VS3XySqVnaz97zP7qlXhqdHJqhpVVnlVVXf3T43PkPt0TYOF3SJBfSg5+6fVYpFFZ4UXq8UfYKSzQ07tMRaLVJ9rao88+3Xdz7OCj+Wc/zi3zZnXjb+nmmjf4JgL1KMmzxUcwQh6waslSO8TjIqCVksQ3qMdXZeg1RKE93nm1jHqEmdv+xu1QpuGaY4eParU1FStXbtW3/jGNxptc8stt6iyslIrVqzw7xs3bpxGjRqlZ599tkXnYZgG7VG116eTVV6drKqp/enxqrKqLqx4vPLU+FRV45Onxlv3s6nXtfs8NT5Ve33y+gzV+Izan966nz6fahq8NuSt21f/usbrO+8OJABoqf/cf51SE4O7vlNIhmnO5XLVLlqVnJzcZJv169fr17/+dYN9kyZN0vLly5s8xuPxyOM502XudrvbUiYQEtE2q5yx1nZ1a7GvLpj46v6N4TMM+Yzan4ZPMnTmtc8wZNT/7pyfPqP2OUNn/2y8Xe156jPQmX/aGA1en/t7o4njDDU8oKnfn/++Dd9PTbW/0HFBEtyZeMF7s2DPEAzm2wWzNiPIf6LBrS2I7xXkP9CkGPP+X9bqMOLz+XTPPffoiiuuaHa4pbi4WGlpaQ32paWlqbi4uMljcnJy9PDDD7e2NKDTslotsnNnEIAI0+rlLGfNmqUvvvhCubm5waxHkjR37ly5XC7/VlhYGPRzAACA9qFVPSOzZ8/WihUrtG7dOvXq1avZtunp6SopKWmwr6SkROnp6U0e43A45HDwNFcAADqDgHpGDMPQ7NmztWzZMn344Ye6+OKLL3hMdna2Vq9e3WDfBx98oOzs7MAqBQAAHVJAPSOzZs3SkiVL9PbbbysxMdE/78PpdCo2NlaSNGPGDPXs2VM5OTmSpDlz5uiqq67SggULNHnyZOXm5mrTpk1atGhRkD8KAACIRAH1jCxcuFAul0tXX321MjIy/Nsbb7zhb1NQUKCioiL/6/Hjx2vJkiVatGiRRo4cqTfffFPLly9njREAACCJ5eABAECItPT7u9V30wAAAAQDYQQAAJiKMAIAAExFGAEAAKYijAAAAFMRRgAAgKkIIwAAwFStfmpvONUvheJ2u02uBAAAtFT99/aFljSLiDBSXl4uScrMzDS5EgAAEKjy8nI5nc4mfx8RK7D6fD4dOXJEiYmJslgsQXtft9utzMxMFRYWsrJrCHGdw4drHR5c5/DgOodHKK+zYRgqLy9Xjx49ZLU2PTMkInpGrFarevXqFbL3T0pK4i96GHCdw4drHR5c5/DgOodHqK5zcz0i9ZjACgAATEUYAQAApurUYcThcOjBBx+Uw+Ewu5QOjescPlzr8OA6hwfXOTzaw3WOiAmsAACg4+rUPSMAAMB8hBEAAGAqwggAADAVYQQAAJiqw4eRZ555Rn369FFMTIwuv/xy/ec//2m2/dKlSzV48GDFxMRoxIgRWrlyZZgqjWyBXOfnn39eEyZMUNeuXdW1a1dNnDjxgn8uqBXo3+d6ubm5slgsmjp1amgL7EACvdZlZWWaNWuWMjIy5HA4NHDgQP7/0QKBXuc///nPGjRokGJjY5WZmalf/epXOn36dJiqjUzr1q3TTTfdpB49eshisWj58uUXPGbNmjUaM2aMHA6H+vfvr8WLF4e2SKMDy83NNex2u/HSSy8Z27dvN+666y6jS5cuRklJSaPtP/nkE8Nmsxl//OMfjR07dhgPPPCAER0dbeTn54e58sgS6HW+9dZbjWeeecbYunWrsXPnTuNHP/qR4XQ6jUOHDoW58sgS6HWut3//fqNnz57GhAkTjClTpoSn2AgX6LX2eDzGJZdcYtx4443Gxx9/bOzfv99Ys2aNsW3btjBXHlkCvc6vvfaa4XA4jNdee83Yv3+/8d577xkZGRnGr371qzBXHllWrlxp3H///cZbb71lSDKWLVvWbPt9+/YZcXFxxq9//Wtjx44dxl/+8hfDZrMZq1atClmNHTqMXHbZZcasWbP8r71er9GjRw8jJyen0fY/+MEPjMmTJzfYd/nllxs//elPQ1pnpAv0Op+rpqbGSExMNP72t7+FqsQOoTXXuaamxhg/frzxwgsvGHfccQdhpIUCvdYLFy40+vbta1RVVYWrxA4h0Os8a9Ys49prr22w79e//rVxxRVXhLTOjqQlYeR//ud/jGHDhjXYd8sttxiTJk0KWV0ddpimqqpKmzdv1sSJE/37rFarJk6cqPXr1zd6zPr16xu0l6RJkyY12R6tu87nOnnypKqrq5WcnByqMiNea6/zI488otTUVM2cOTMcZXYIrbnW//jHP5Sdna1Zs2YpLS1Nw4cP1x/+8Ad5vd5wlR1xWnOdx48fr82bN/uHcvbt26eVK1fqxhtvDEvNnYUZ34UR8aC81vj666/l9XqVlpbWYH9aWpq+/PLLRo8pLi5utH1xcXHI6ox0rbnO5/rNb36jHj16nPeXH2e05jp//PHHevHFF7Vt27YwVNhxtOZa79u3Tx9++KGmT5+ulStXas+ePfr5z3+u6upqPfjgg+EoO+K05jrfeuut+vrrr3XllVfKMAzV1NTo7rvv1u9+97twlNxpNPVd6Ha7derUKcXGxgb9nB22ZwSRYf78+crNzdWyZcsUExNjdjkdRnl5uW6//XY9//zz6t69u9nldHg+n0+pqalatGiRxo4dq1tuuUX333+/nn32WbNL61DWrFmjP/zhD/rrX/+qLVu26K233tI777yjefPmmV0a2qjD9ox0795dNptNJSUlDfaXlJQoPT290WPS09MDao/WXed6TzzxhObPn69//etfysrKCmWZES/Q67x3714dOHBAN910k3+fz+eTJEVFRWnXrl3q169faIuOUK35O52RkaHo6GjZbDb/viFDhqi4uFhVVVWy2+0hrTkSteY6//73v9ftt9+uO++8U5I0YsQIVVZW6ic/+Ynuv/9+Wa38+zoYmvouTEpKCkmviNSBe0bsdrvGjh2r1atX+/f5fD6tXr1a2dnZjR6TnZ3doL0kffDBB022R+uusyT98Y9/1Lx587Rq1Spdcskl4Sg1ogV6nQcPHqz8/Hxt27bNv33nO9/RNddco23btikzMzOc5UeU1vydvuKKK7Rnzx5/4JOkr776ShkZGQSRJrTmOp88efK8wFEfAA0esxY0pnwXhmxqbDuQm5trOBwOY/HixcaOHTuMn/zkJ0aXLl2M4uJiwzAM4/bbbzd++9vf+tt/8sknRlRUlPHEE08YO3fuNB588EFu7W2BQK/z/PnzDbvdbrz55ptGUVGRfysvLzfrI0SEQK/zubibpuUCvdYFBQVGYmKiMXv2bGPXrl3GihUrjNTUVOPRRx816yNEhECv84MPPmgkJiYar7/+urFv3z7j/fffN/r162f84Ac/MOsjRITy8nJj69atxtatWw1JxpNPPmls3brVOHjwoGEYhvHb3/7WuP322/3t62/tve+++4ydO3cazzzzDLf2ttVf/vIX46KLLjLsdrtx2WWXGRs2bPD/7qqrrjLuuOOOBu3//ve/GwMHDjTsdrsxbNgw45133glzxZEpkOvcu3dvQ9J524MPPhj+wiNMoH+fz0YYCUyg1zovL8+4/PLLDYfDYfTt29d47LHHjJqamjBXHXkCuc7V1dXGQw89ZPTr18+IiYkxMjMzjZ///OfGiRMnwl94BPnoo48a/X9u/bW94447jKuuuuq8Y0aNGmXY7Xajb9++xssvvxzSGi2GQd8WAAAwT4edMwIAACIDYQQAAJiKMAIAAExFGAEAAKYijAAAAFMRRgAAgKkIIwAAwFSEEQAAYCrCCAAAndS6det00003qUePHrJYLFq+fHnA72EYhp544gkNHDhQDodDPXv21GOPPRbQe3TYp/YCAIDmVVZWauTIkfrv//5vfe9732vVe8yZM0fvv/++nnjiCY0YMULHjx/X8ePHA3oPloMHAACyWCxatmyZpk6d6t/n8Xh0//336/XXX1dZWZmGDx+uxx9/XFdffbUkaefOncrKytIXX3yhQYMGtfrcDNMAAIBGzZ49W+vXr1dubq4+//xz3Xzzzbr++uu1e/duSdI///lP9e3bVytWrNDFF1+sPn366M477wy4Z4QwAgAAzlNQUKCXX35ZS5cu1YQJE9SvXz/de++9uvLKK/Xyyy9Lkvbt26eDBw9q6dKleuWVV7R48WJt3rxZ//Vf/xXQuZgzAgAAzpOfny+v16uBAwc22O/xeNStWzdJks/nk8fj0SuvvOJv9+KLL2rs2LHatWtXi4duCCMAAOA8FRUVstls2rx5s2w2W4PfJSQkSJIyMjIUFRXVILAMGTJEUm3PCmEEAAC02ujRo+X1elVaWqoJEyY02uaKK65QTU2N9u7dq379+kmSvvrqK0lS7969W3wu7qYBAKCTqqio0J49eyTVho8nn3xS11xzjZKTk3XRRRfptttu0yeffKIFCxZo9OjROnr0qFavXq2srCxNnjxZPp9Pl156qRISEvTnP/9ZPp9Ps2bNUlJSkt5///0W10EYAQCgk1qzZo2uueaa8/bfcccdWrx4saqrq/Xoo4/qlVde0eHDh9W9e3eNGzdODz/8sEaMGCFJOnLkiH7xi1/o/fffV3x8vG644QYtWLBAycnJLa6DMAIAAEzFrb0AAMBUhBEAAGAqwggAADAVYQQAAJiKMAIAAExFGAEAAKYijAAAAFMRRgAAgKkIIwAAwFSEEQAAYCrCCAAAMNX/BS7/tMjM8pZgAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sns.lineplot(x=J_all2[0], y=J_all2[1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "id": "4d0fbfee",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "d680d4cc18984adab5920af97da76e2e",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "  0%|          | 0/10000000 [00:00<?, ?it/s]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Die 3 Parameter der linearen Regression:\n",
      "[ 540.5020477  1598.80804052 6469.41806027]\n",
      "Kostenfunktion J: 200954758401.09796\n",
      "J relativ zu Startkosten: 0.9967438136028319\n",
      "Vergleich Kostenfunktion zu analytischer Lösung: 1.36*J_ana\n",
      "Relative Abweichung der Parameter zu analytischer Lösung: [0.01263744 0.00687326 2.66104765]*w_ana\n",
      "CPU times: user 37min 33s, sys: 9min 27s, total: 47min 1s\n",
      "Wall time: 5min 1s\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "alpha = 3.1e-10  # verschiedene alpha ausprobieren\n",
    "n_iterations = 10000000\n",
    "w_gd3, J_all3 = grad_desc(w=w_gd2, alpha=alpha, x=X, y=y, n_iterations=n_iterations)\n",
    "\n",
    "print('Die {} Parameter der linearen Regression:\\n{}'.format(len(w_gd3), w_gd3))\n",
    "print('Kostenfunktion J: {}'.format(J_all3[1][-1]))\n",
    "print('J relativ zu Startkosten: {}'.format(J_all3[1][-1]/J_all3[1][0]))\n",
    "print('Vergleich Kostenfunktion zu analytischer Lösung: {:.2f}*J_ana'.format(J_all3[1][-1]/J_ana))\n",
    "print('Relative Abweichung der Parameter zu analytischer Lösung: {}*w_ana'.format((w_gd3)/w_ana))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "id": "252656f1",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<Axes: >"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjUAAAG+CAYAAABrivUeAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAA9hAAAPYQGoP6dpAABPtElEQVR4nO3deVxU9f4/8NcZdmQYBWVRUHFBXFDcWdXKiwoutKipuWUuCYL1bbO6V28b2c1ubrmVaCmRmKgh2SUXAgFNBZUUXBEUQRQdcJB1zu8Pb9wfCcKwHWbm9Xw8zuNxPfM557w4t5xX854ZBFEURRARERFpOZnUAYiIiIiaAksNERER6QSWGiIiItIJLDVERESkE1hqiIiISCew1BAREZFOYKkhIiIincBSQ0RERDqBpYaIiIh0AksNERER6QS9LDW//fYbJkyYgI4dO0IQBOzdu1ej40tKSjBnzhy4urrC0NAQAQEBj625desWpk+fDmdnZ8hkMixdurRJshMREVHN9LLUqFQqDBgwAOvXr2/Q8ZWVlTAzM0NwcDBGjx5d45rS0lJ06NAB77//PgYMGNCYuERERFQPhlIHkMK4ceMwbty4Wh8vLS3Fe++9h++//x73799Hv379sHLlSowaNQoA0KZNG2zYsAEAcOzYMdy/f/+xc3Tt2hWrV68GAGzdurXJfwYiIiKqTi9fqalLUFAQkpKSEBERgbNnz2Ly5MkYO3YsLl26JHU0IiIiqgVLzV9kZWUhLCwMkZGR8PHxQffu3fHGG2/A29sbYWFhUscjIiKiWujl+OlJzp07h8rKSjg7O1fbX1paCmtra4lSERERUV1Yav7iwYMHMDAwwKlTp2BgYFDtMQsLC4lSERERUV1Yav5i4MCBqKysxO3bt+Hj4yN1HCIiIqonvSw1Dx48wOXLl6v+fO3aNaSmpsLKygrOzs6YMWMGZs2ahVWrVmHgwIHIz8/HoUOH0L9/f/j7+wMAzp8/j7KyMhQUFKCoqAipqakAADc3t6rz/rnvwYMHyM/PR2pqKoyNjdGnT5+W+lGJiIj0hiCKoih1iJZ29OhRPPXUU4/tnz17NrZt24by8nJ89NFH+Pbbb3Hz5k20b98e7u7u+Oc//wlXV1cAjz6yff369cfO8f/fTkEQHnu8S5cuyMzMbLofhoiIiADoaakhIiIi3cOPdBMREZFOYKkhIiIinaA3bxRWq9XIycmBXC6v8b0uRERE1PqIooiioiJ07NgRMtmTX4vRm1KTk5MDR0dHqWMQERFRA2RnZ8PBweGJa/Sm1MjlcgCPboqlpaXEaYiIiKg+CgsL4ejoWPU8/iR6U2r+HDlZWlqy1BAREWmZ+rx1hG8UJiIiIp3AUkNEREQ6gaWGiIiIdAJLDREREekElhoiIiLSCSw1REREpBNYaoiIiEgnsNQQERGRTmCpISIiIp3AUkNEREQ6gaWGiIiIdAJLTRMoq1BLHYGIiEjvsdQ0UlmFGlM2JSE05gLKK1luiIiIpKI3v6W7uRxOv43U7PtIzb6P3zMLsG76IHRsayZ1LCIiIr3DV2oaaWw/O2x8aRDkpoY4nXUffmvicehCntSxiIiI9I5GpSY0NBRDhw6FXC6HjY0NAgICkJGRUedxkZGRcHFxgampKVxdXRETE1Pt8T179sDX1xfW1tYQBAGpqak1nicpKQlPP/002rRpA0tLS4wYMQIPHz7U5EdoFmP72ePAEh/0d1DgfnE55m0/yXEUERFRC9Oo1MTFxSEwMBDJycmIjY1FeXk5fH19oVKpaj0mMTER06ZNw7x585CSkoKAgAAEBAQgLS2tao1KpYK3tzdWrlxZ63mSkpIwduxY+Pr64sSJE/j9998RFBQEmax1vNjU2dockYs8MMezKwBg029XMXVTEm7el750ERER6QNBFEWxoQfn5+fDxsYGcXFxGDFiRI1rpk6dCpVKhejo6Kp97u7ucHNzw8aNG6utzczMhJOTE1JSUuDm5lbtMXd3d/ztb3/Dhx9+2KCshYWFUCgUUCqVsLS0bNA56utg2i28ufssikoq0NbcCKsmD8AzvW2b9ZpERES6SJPn70a9zKFUKgEAVlZWta5JSkrC6NGjq+0bM2YMkpKS6n2d27dv4/jx47CxsYGnpydsbW0xcuRIJCQk1HpMaWkpCgsLq20theMoIiKiltfgUqNWq7F06VJ4eXmhX79+ta7Lzc2FrW31VylsbW2Rm5tb72tdvXoVALBixQrMnz8fBw8exKBBg/DMM8/g0qVLNR4TGhoKhUJRtTk6Otb7ek2B4ygiIqKW1eBSExgYiLS0NERERDRlnhqp1Y9e4Vi4cCHmzp2LgQMH4t///jd69eqFrVu31njMsmXLoFQqq7bs7Oxmz/lXJoYGWDGxb7VPR/nz01FERETNokGlJigoCNHR0Thy5AgcHByeuNbOzg55edWfxPPy8mBnZ1fv69nb2wMA+vTpU21/7969kZWVVeMxJiYmsLS0rLZJheMoIiKi5qdRqRFFEUFBQYiKisLhw4fh5ORU5zEeHh44dOhQtX2xsbHw8PCo93W7du2Kjh07Pvbx8YsXL6JLly71Po+U/hxHzfXqCoDjKCIioqam0TcKBwYGIjw8HPv27YNcLq96X4xCoYCZ2aNv0Z01axY6deqE0NBQAEBISAhGjhyJVatWwd/fHxERETh58iQ2b95cdd6CggJkZWUhJycHAKrKi52dHezs7CAIAt58800sX74cAwYMgJubG7Zv34709HTs3r278XehhZgYGmD5hL4Y7mSFN3efrRpH8dNRRERETUDUAIAat7CwsKo1I0eOFGfPnl3tuF27donOzs6isbGx2LdvX/HAgQPVHg8LC6vxvMuXL6+2LjQ0VHRwcBDNzc1FDw8PMT4+vt7ZlUqlCEBUKpWa/MjNJuuuSpywNl7s8na02OXtaPHjA+fFsopKqWMRERG1Kpo8fzfqe2q0SUt+T019lVZU4tOf0xF2LBMAMKhzW6ydPgid+LujiIiIALTg99RQ4/w5jtr40uD//e6o1fx0FBERUUOw1LQCY/vZISbYBwMcFFA+fPTpqE/46SgiIiKNsNS0Eo5W5ohc5Fn16ajNv13FFH46ioiIqN5YaloRY0NZtXFUyn/HUb+e5ziKiIioLiw1rdBfx1GvfMtxFBERUV1YalqpP8dRL3s9+oLDP8dRN+4VS5yMiIiodWKpacWMDWX4x4Q+2DTzf+Mo/zUJHEcRERHVgKVGC4zp+99xlGPbqnHUR9HnUVbBcRQREdGfWGq0hKOVOSIXemCe96Nx1NcJ1zBlUxKyCziOIiIiAlhqtIqxoQx/H98Hm2cOhqWpIVKzH/3uqP/8kSt1NCIiIsmx1Ggh3752OBDsAzfHtigsqcCC707hg584jiIiIv3GUqOlHK3MsWuhB+b7PBpHbT12DZM3JnIcRUREeoulRosZG8rwnn8ffD1rCBRmRjhzQwm/NfE4mMZxFBER6R+WGh0wuo8tYkJ8MLBzWxSVVGDRjlNYsf8PlFZUSh2NiIioxbDU6IhObc2wa6EHFo7oBgDYlpiJFzYkIesux1FERKQfWGp0iJGBDMv8emPrnCFoa26EczeV8F8Tj5/P3ZI6GhERUbNjqdFBT7vYIibYB4O7tENRaQVe3Xkay/elcRxFREQ6jaVGR3Vsa4aIBe5YNLI7AGB70nU8vyER1++qJE5GRETUPFhqdJiRgQzvjHNB2JyhaGduhLSbhRi/JgEHznIcRUREuoelRg885WKDmBAfDPnvOCow/DT+vjcNJeUcRxERke5gqdET9opH46jFox6No75LfjSOunaH4ygiItINLDV6xNBAhrfGumDb3KGwamOMP3IKMWFtAn46kyN1NCIiokZjqdFDo3rZICbYB8O6WuFBaQWWfJ+C96LOcRxFRERajaVGT9kpTBE+fziCnuoBQQB2Hs/Cs18l4mr+A6mjERERNQhLjR4zNJDhjTG9sH3uMFi3McaFW4/GUftSb0odjYiISGMsNYQRzh0QE+KD4U5WUJVVIiQiFcv2cBxFRETahaWGAAC2lqbY+cpwBD/9aBz1/YksBKw/hiscRxERkZZgqaEqhgYyvO7bC9+9PBztLYyRnluECWsTsDeF4ygiImr9WGroMd492yMm2Afu3axQXFaJpT+k4p0fz3IcRURErRpLDdXIxtIUO19xR8gzPSEIQMTv2Zi07hgu3+Y4ioiIWieWGqqVgUzAa39zxo55w9HewgQZeUWYuC4Be07fkDoaERHRY1hqqE5ePdojJsQbnt2tUVxWidd3ncFbu8/gYRnHUURE1Hqw1FC92MhN8d284XhttDMEAdh18gYmrU/ApbwiqaMREREBYKkhDRjIBISM7omdrwxHB7kJLuY9wMR1x7D7FMdRREQkPZYa0phn90efjvLu0R4PyyvxRuQZvBF5BsVlFVJHIyIiPcZSQw3SQW6C7S8Pw//9zRkyAdh96gYmrTuGixxHERGRRFhqqMEMZAKWPNMT4fPdYSM3waXbDzBxXQIiT2ZLHY2IiPQQSw01mns3a8SE+MCnZ3uUlKvx5u6zeH1XKlSlHEcREVHLYamhJtHewgTb5w7Dm2N6QSYAe07fxMR1CcjI5TiKiIhaBksNNRmZTEDgUz3w/Xx32Fqa4Eq+CpPWJ+CH37MgiqLU8YiISMdpVGpCQ0MxdOhQyOVy2NjYICAgABkZGXUeFxkZCRcXF5iamsLV1RUxMTHVHt+zZw98fX1hbW0NQRCQmpr62DlGjRoFQRCqbYsWLdIkPrWQ4d2sERPsgxHOHVBSrsbbP57D67vOcBxFRETNSqNSExcXh8DAQCQnJyM2Nhbl5eXw9fWFSqWq9ZjExERMmzYN8+bNQ0pKCgICAhAQEIC0tLSqNSqVCt7e3li5cuUTrz9//nzcunWravvss880iU8tyNrCBNvmDMVbY3vBQCYgKuUmJqxLwIVbhVJHIyIiHSWIjZgL5Ofnw8bGBnFxcRgxYkSNa6ZOnQqVSoXo6Oiqfe7u7nBzc8PGjRurrc3MzISTkxNSUlLg5uZW7bFRo0bBzc0NX375ZYOyFhYWQqFQQKlUwtLSskHnoIb5PbMAS8JTkFtYAhNDGVZM7IsXhzpCEASpoxERUSunyfN3o95To1QqAQBWVla1rklKSsLo0aOr7RszZgySkpI0vt7OnTvRvn179OvXD8uWLUNxcXGta0tLS1FYWFhtI2kM7WqFmBAfjOrVAaUVaizbcw5Lf0jFA46jiIioCTW41KjVaixduhReXl7o169fretyc3Nha2tbbZ+trS1yc3M1ut706dOxY8cOHDlyBMuWLcN3332Hl156qdb1oaGhUCgUVZujo6NG16OmZdXGGFtnD8U741xgIBOwLzUHE9cm4HwOyyYRETUNw4YeGBgYiLS0NCQkJDRlnlotWLCg6n+7urrC3t4ezzzzDK5cuYLu3bs/tn7ZsmV4/fXXq/5cWFjIYiMxmUzAopHdMaRLOyz5PgVX76gQ8NUxLJ/QB9OHdeY4ioiIGqVBr9QEBQUhOjoaR44cgYODwxPX2tnZIS8vr9q+vLw82NnZNeTSVYYPHw4AuHz5co2Pm5iYwNLSstpGrcOQrlaICfbB0y42KKtQ472oNCz5PgVFJeVSRyMiIi2mUakRRRFBQUGIiorC4cOH4eTkVOcxHh4eOHToULV9sbGx8PDw0CzpX/z5sW97e/tGnYek0a6NMb6eNQTv+rnAUCYg+uwtTFibgLSbSqmjERGRltJo/BQYGIjw8HDs27cPcrm86n0xCoUCZmZmAIBZs2ahU6dOCA0NBQCEhIRg5MiRWLVqFfz9/REREYGTJ09i8+bNVectKChAVlYWcnJyAKDqu2/s7OxgZ2eHK1euIDw8HH5+frC2tsbZs2fx2muvYcSIEejfv3/j7wJJQiYTsGBEdwzuYoUl4aeRebcYz21IxN/H98FLwzmOIiIizWj0ke7anmTCwsIwZ84cAI8+et21a1ds27at6vHIyEi8//77yMzMRM+ePfHZZ5/Bz8+v6vFt27Zh7ty5j513+fLlWLFiBbKzs/HSSy8hLS0NKpUKjo6OePbZZ/H+++/Xe6zEj3S3bveLy/BG5Bn8euE2AMC/vz1Cn3OFpamRxMmIiEhKmjx/N+p7arQJS03rJ4oivkm4hk9/TkeFWkQXa3Osnz4I/ToppI5GREQSabHvqSFqSoIg4BWfbti1yAOd2prh+t1iPPdVIr5NyuTvjiIiojqx1FCrM6hzOxwI9sbo3rYoq1TjH/v+QGD4aRTy01FERPQELDXUKrU1N8aWWYPx9/F9YGQgIOZcLsavScC5G/x0FBER1YylhlotQRAwz9sJkYs84dDODFkFxXh+QyK2HbvGcRQRET2GpYZaPTfHtjiwxAe+fR6No1b8dB6v7jgN5UOOo4iI6H9YakgrKMyNsGnmYCyf8GgcdfCPXIxfG48z2feljkZERK0ESw1pDUEQMNfLCbsXecLRygzZBQ/xwsZEbE3gOIqIiFhqSAsNcGyL6CU+GNvXDuWVIj6IPo+F352CspjjKCIifcZSQ1pJYWaEDS8Nwj8n9oWxgQz/OZ8HvzXxSOU4iohIb7HUkNYSBAGzPbvix1c90dnKHDfvP8QLGxLxdfxVjqOIiPQQSw1pPVcHBaKDveHnaocKtYiPDlzA/G9P4X5xmdTRiIioBbHUkE6wNDXC+umD8OGkR+OoXy/kwX9NAk5n3ZM6GhERtRCWGtIZgiBgpkdX7FnsiS7Wj8ZRUzYmYctvHEcREekDlhrSOf06KRC9xBv+/e1RoRbxccwFvLL9JO6pOI4iItJlLDWkk+SmRlg3bSA+CugHY0MZDqXfhv+aeJy6znEUEZGuYqkhnSUIAl5y74KoxZ5wat8GOcoSTN2UhE1xV6BWcxxFRKRrWGpI5/XtqMD+IC9MGNARFWoRoT+nY97231HAcRQRkU5hqSG9IDc1wpoX3fDJs64wNpThSEY+/NfE42RmgdTRiIioibDUkN4QBAHTh3fG3sVe6Na+DW4pSzB1czI2HOU4iohIF7DUkN7p09ES+5d4Y5JbR1SqRaw8mI6Xt/+Ouw9KpY5GRESNwFJDesnCxBBfTnXDp8+5wsRQhqMZ+fBfk4AT1ziOIiLSViw1pLcEQcCLwzpjb6AXunVog9zCEkzbkoz1Ry5zHEVEpIVYakjv9ba3xE9B3nh2YCdUqkX865cMzNnGcRQRkbZhqSEC0MbEEF9MGYDPnu8PUyMZfruYD7818Th+9a7U0YiIqJ5Yaoj+SxAETBnqiH2B3ujeoQ3yCksxbUsy1h2+xHEUEZEWYKkh+otednLsD/LGc4M6QS0Cn//nImaHncAdjqOIiFo1lhqiGjwaR7nhXy88GkfFX7oDv9XxSLrCcRQRUWvFUkP0BJOHOGJ/kDd62ljgdlEpZnydjNW/XkIlx1FERK0OSw1RHZxt5dgX5IXJgx2gFoF//3oRs7YeR34Rx1FERK0JSw1RPZgbG+Jfkwfg88kDYGZkgGOX78JvTTwSL9+ROhoREf0XSw2RBl4Y7ID9QV5wtrVAflEpZnxzHF/+epHjKCKiVoClhkhDPW3l2BfojSlDHCCKwJe/XsLMb47jdlGJ1NGIiPQaSw1RA5gZG+CzFwbgiykDYG5sgMQrd+G3OgHHOI4iIpIMSw1RIzw3yAH7g7zRy1aOOw9K8dI3x/FFLMdRRERSYKkhaqQeNhbYG+iFF4c6QhSBNYcuYcbXybhdyHEUEVFLYqkhagJmxgb49Pn++HKqG8yNDZB8tQB+a+IRfylf6mhERHqDpYaoCQUM7ISflnjDxU6OOw/KMGvrCXz+SwYqKtVSRyMi0nksNURNrHuHR+Oo6cM7QxSBdUcuY/rXx5HHcRQRUbNiqSFqBqZGBvjkWVesftENbYwNcOJaAfxWxyPuIsdRRETNhaWGqBlNcns0juptb4m7qjLM3noC//olneMoIqJmwFJD1My6dbBA1GJPzBjeGQCw/sgVTN9yHLeUDyVORkSkWzQqNaGhoRg6dCjkcjlsbGwQEBCAjIyMOo+LjIyEi4sLTE1N4erqipiYmGqP79mzB76+vrC2toYgCEhNTa31XKIoYty4cRAEAXv37tUkPpFkTI0M8PGzrlg7bSAsTAxxIvPROOpIxm2poxER6QyNSk1cXBwCAwORnJyM2NhYlJeXw9fXFyqVqtZjEhMTMW3aNMybNw8pKSkICAhAQEAA0tLSqtaoVCp4e3tj5cqVdWb48ssvIQiCJrGJWo0JAzoieok3+na0xL3icswN+x2f/pyOco6jiIgaTRBFscFffZqfnw8bGxvExcVhxIgRNa6ZOnUqVCoVoqOjq/a5u7vDzc0NGzdurLY2MzMTTk5OSElJgZub22PnSk1Nxfjx43Hy5EnY29sjKioKAQEBNV63tLQUpaWlVX8uLCyEo6MjlEolLC0tNf9hiZpQSXklPj5wAd8lXwcADO7SDmunDUTHtmYSJyMial0KCwuhUCjq9fzdqPfUKJVKAICVlVWta5KSkjB69Ohq+8aMGYOkpCSNrlVcXIzp06dj/fr1sLOzq3N9aGgoFApF1ebo6KjR9Yiak6mRAT4M6IevZgyC3MQQp67fg9+aePx6Pk/qaEREWqvBpUatVmPp0qXw8vJCv379al2Xm5sLW1vbavtsbW2Rm5ur0fVee+01eHp6YtKkSfVav2zZMiiVyqotOztbo+sRtQQ/V3tEB3vDtZMC94vL8cq3J/FR9HmUVXAcRUSkKcOGHhgYGIi0tDQkJCQ0ZZ4a7d+/H4cPH0ZKSkq9jzExMYGJiUkzpiJqGl2s22D3qx749Od0hB3LxNcJ13Dy+j2snTYQjlbmUscjItIaDXqlJigoCNHR0Thy5AgcHByeuNbOzg55edVfUs/Ly6vXCOlPhw8fxpUrV9C2bVsYGhrC0PBRF3v++ecxatQojfMTtTYmhgZYPqEvNs0cDEtTQ6Rm34f/mnj88odmr2gSEekzjUqNKIoICgpCVFQUDh8+DCcnpzqP8fDwwKFDh6rti42NhYeHR72v+8477+Ds2bNITU2t2gDg3//+N8LCwjT5EYhatTF97XAg2Adujm1RWFKBhd+dwor9f6C0olLqaERErZ5G46fAwECEh4dj3759kMvlVe+LUSgUMDN79KmNWbNmoVOnTggNDQUAhISEYOTIkVi1ahX8/f0RERGBkydPYvPmzVXnLSgoQFZWFnJycgCg6rtv7Ozsqm1/1blz53oVKyJt4mhljl0LPfCvX9KxJf4atiVm4tT1e1g/fRA6W3McRURUG41eqdmwYQOUSiVGjRoFe3v7qu2HH36oWpOVlYVbt25V/dnT0xPh4eHYvHkzBgwYgN27d2Pv3r3V3ly8f/9+DBw4EP7+/gCAF198EQMHDnzsI99E+sLYUIb3/Pvgm9lD0NbcCOduKuG/Jh4x527VfTARkZ5q1PfUaBNNPudO1Jrk3H+IJd+n4NT1ewCAme5d8J5/b5gaGUicjIio+bXY99QQUfPr2NYMEQvc8eqo7gCA75Kv47mvEnHtTu3f5E1EpI9Yaoi0gJGBDG+PdcG2uUNh1cYY528VYvyaeOxLvSl1NCKiVoOlhkiLjOplg5hgHwxzsoKqrBIhEalYtucsSsr56SgiIpYaIi1jpzBF+CvDseTpHhAE4PsT2QhYfwyXbz+QOhoRkaRYaoi0kKGBDP/n2wvfvTwc7S2MkZ5bhInrErDn9A2poxERSYalhkiLefdsj5hgH3h0s0ZxWSVe33UGb0aeQXFZhdTRiIhaHEsNkZazsTTFjleG47XRzpAJQOSpG5i07hgu5hVJHY2IqEWx1BDpAAOZgJDRPbHzFXd0kJvg0u0HmLguAbtOZkNPvoqKiIilhkiXeHS3xs8hPvDp2R4l5Wq8tfssXt91BqpSjqOISPex1BDpmPYWJtg+dxjeHNMLMgGISrmJCesScOFWodTRiIiaFUsNkQ6SyQQEPtUDEQs8YGdpiqv5KgSsP4bw41kcRxGRzmKpIdJhw5ysEBPig1G9OqC0Qo13o84hOCIVRSXlUkcjImpyLDVEOs6qjTG2zh6KZeNcYCAT8NOZHExYm4C0m0qpoxERNSmWGiI9IJMJWDiyO3Yt9EBHhSky7xbjua8S8V1SJsdRRKQzWGqI9MjgLu0QE+KD0b1tUVapxt/3/YHA8NMo5DiKiHQASw2Rnmlrbowtswbjff/eMDIQEHMuF/5r4nH2xn2poxERNQpLDZEeEgQBr/h0Q+QiTzi0M0N2wUM8vyERWxOucRxFRFqLpYZIj7k5tsWBYB+M7WuH8koRH0Sfx8LvTkFZzHEUEWkflhoiPacwM8KGlwbhnxP7wthAhv+cz4PfmnikZN2TOhoRkUZYaogIgiBgtmdX/PiqJ7pYm+Pm/YeYvDEJW367ynEUEWkNlhoiquLqoMBPS7zh398eFWoRH8dcwCvbT+KeqkzqaEREdWKpIaJqLE2NsG7aQHwU0A/GhjIcSr8NvzXxOJlZIHU0IqInYqkhoscIgoCX3LsgarEnnNq3wS1lCaZuTsZXRy9DreY4iohaJ5YaIqpV346PxlGT3DqiUi3is4MZmLvtd9x9UCp1NCKix7DUENETWZgY4supblj5vCtMDGWIu5gPvzXxOH71rtTRiIiqYakhojoJgoCpQztjf5A3undog7zCUkzbkoy1hy6hkuMoImolWGqIqN562cnx0xJvPD/IAWoRWBV7EbO3nkB+EcdRRCQ9lhoi0oi5sSFWTRmAzycPgJmRARIu38G41fFIvHxH6mhEpOdYaoioQV4Y7ID9QV5wtrXAnQelmPHNcXwRe5HjKCKSDEsNETVYT1s59gV648WhjhBFYM2hS5jxdTLyCkukjkZEeoilhogaxczYAJ8+3x+rX3RDG2MDJF8tgN/qePx2MV/qaESkZ1hqiKhJTHLrhJ+WeKO3vSXuqsowO+wE/vVLOioq1VJHIyI9wVJDRE2mWwcLRC32xIzhnSGKwPojVzB9y3HcUj6UOhoR6QGWGiJqUqZGBvj4WVesmz4QFiaGOJH5aBx1JP221NGISMex1BBRsxjfvyOil3ijXydL3Csux9xtvyM05gLKOY4iombCUkNEzaZr+zb48VVPzPHsCgDY9NtVTN2UhJv3OY4ioqbHUkNEzcrE0AArJvbFxpcGQW5qiNNZ9+G3Oh6x5/OkjkZEOoalhohaxNh+9ogJ9sEABwWUD8sx/9uT+DD6PMoqOI4ioqbBUkNELcbRyhyRizwxz9sJAPBNwjVM3pSE7IJiiZMRkS5gqSGiFmVsKMPfx/fBlllDoDAzwpns+/BbE4+DabekjkZEWk6jUhMaGoqhQ4dCLpfDxsYGAQEByMjIqPO4yMhIuLi4wNTUFK6uroiJian2+J49e+Dr6wtra2sIgoDU1NTHzrFw4UJ0794dZmZm6NChAyZNmoT09HRN4hNRK/K3PrY4EOyNgZ3boqikAot2nMbyfWkoraiUOhoRaSmNSk1cXBwCAwORnJyM2NhYlJeXw9fXFyqVqtZjEhMTMW3aNMybNw8pKSkICAhAQEAA0tLSqtaoVCp4e3tj5cqVtZ5n8ODBCAsLw4ULF/DLL79AFEX4+vqispJ/ARJpK4d25ti10AMLR3YDAGxPuo4XNiTh+t3a/04hIqqNIIpig3+lbn5+PmxsbBAXF4cRI0bUuGbq1KlQqVSIjo6u2ufu7g43Nzds3Lix2trMzEw4OTkhJSUFbm5uT7z22bNnMWDAAFy+fBndu3evM2thYSEUCgWUSiUsLS3r/uGIqEUdSb+N13el4l5xOSxMDPHp864Y37+j1LGISGKaPH836j01SqUSAGBlZVXrmqSkJIwePbravjFjxiApKanB11WpVAgLC4OTkxMcHR1rXFNaWorCwsJqGxG1Xk+52CAmxAdDu7bDg9IKBIWn4L2ocygp56uxRFQ/DS41arUaS5cuhZeXF/r161frutzcXNja2lbbZ2tri9zcXI2v+dVXX8HCwgIWFhb4+eefERsbC2Nj4xrXhoaGQqFQVG21lR8iaj3sFWb4fr47Ap/qDkEAdh7PwrNfJeJq/gOpoxGRFmhwqQkMDERaWhoiIiKaMs8TzZgxAykpKYiLi4OzszOmTJmCkpKSGtcuW7YMSqWyasvOzm6xnETUcIYGMrw5xgXb5w6DdRtjXLhViAlrE7Av9abU0YiolWtQqQkKCkJ0dDSOHDkCBweHJ661s7NDXl71bw7Ny8uDnZ2dxtdVKBTo2bMnRowYgd27dyM9PR1RUVE1rjUxMYGlpWW1jYi0xwjnDogJ8YF7NyuoyioREpGKd348i4dlHEcRUc00KjWiKCIoKAhRUVE4fPgwnJyc6jzGw8MDhw4dqrYvNjYWHh4emiWtIYsoiigtLW3UeYio9bK1NMXOV9wR8kxPCAIQ8Xs2AtYfw+XbRVJHI6JWSKNSExgYiB07diA8PBxyuRy5ubnIzc3Fw4f/++V0s2bNwrJly6r+HBISgoMHD2LVqlVIT0/HihUrcPLkSQQFBVWtKSgoQGpqKs6fPw8AyMjIQGpqatX7bq5evYrQ0FCcOnUKWVlZSExMxOTJk2FmZgY/P79G3QAiat0MZAJe+5szds4bjvYWJsjIK8KEtcew+9QNqaMRUSujUanZsGEDlEolRo0aBXt7+6rthx9+qFqTlZWFW7f+982gnp6eCA8Px+bNmzFgwADs3r0be/furfbm4v3792PgwIHw9/cHALz44osYOHBg1Ue+TU1NER8fDz8/P/To0QNTp06FXC5HYmIibGxsGnUDiEg7ePZoj59DfODdoz0ellfijcgz+L9dZ1BcViF1NCJqJRr1PTXahN9TQ6QbKtUivjpyGf/+9SLUItDDxgLrpw9CLzu51NGIqBm02PfUEBG1NAOZgCXP9ET4fHfYWprg8u0HmLQ+AT/8ngU9+W80IqoFSw0RaSX3btaICfbBSOcOKClX4+0fz+G1H1LxoJTjKCJ9xVJDRFrL2sIEYXOG4u2xLjCQCdibmoOJaxNwPoffIE6kj1hqiEiryWQCXh3VHT8scIe9whRX76gQ8NUx7Dx+neMoIj3DUkNEOmFIVyvEBPvgGRcblFWo8V5UGoK+T0FRSbnU0YiohbDUEJHOaNfGGF/PHoL3/HrDUCbgwNlbGL82AWk3lVJHI6IWwFJDRDpFEATMH9ENuxZ5oFNbM1y/W4znvkrE9sRMjqOIdBxLDRHppEGd2yEm2Ae+fWxRVqnG8v1/YPHO01A+5DiKSFex1BCRzlKYG2HTzMFYPqEPjAwE/JyWi/Fr43Em+77U0YioGbDUEJFOEwQBc72csHuRJxytzJBd8BAvbEzENwnXOI4i0jEsNUSkFwY4tsWBYB/4udqhvFLEh9HnMf/bU7hfXCZ1NCJqIiw1RKQ3LE2NsH76IHw4qS+MDWT49UIe/Nck4NT1e1JHI6ImwFJDRHpFEATM9OiKPYs90dXaHDfvP8TUTUnYFHcFajXHUUTajKWGiPRSv04K/LTEGxMGdESFWkToz+l45duTKFBxHEWkrVhqiEhvyU2NsOZFN4Q+5woTQxkOp9+G3+p4/J5ZIHU0ImoAlhoi0muCIGDasM7YG+iFbh3aILewBC9uTsb6I5c5jiLSMiw1REQAettb4qcgbzw3sBMq1SL+9UsGZoedwJ0HpVJHI6J6YqkhIvqvNiaGWDVlAD57oT9MjWSIv3QHfqvjkXTlrtTRiKgeWGqIiP4/giBgyhBH7A/yRk8bC9wuKsWMr5Ox+tdLqOQ4iqhVY6khIqqBs60c+4K8MHmwA9Qi8O9fL2LW1uO4XVQidTQiqgVLDRFRLcyNDfGvyQPwxZQBMDc2wLHLd+G3OgHHLt+ROhoR1YClhoioDs8NcsD+IG+42Mlx50EpXvrmOL74TwbHUUStDEsNEVE99LCxwN5AL0wb1hmiCKw5fBnTtyQjr5DjKKLWgqWGiKieTI0MEPqcK9ZMG4g2xgY4fq0A41bHI+5ivtTRiAgsNUREGps4oCOig33Qx94SBaoyzN56AisPpqOiUi11NCK9xlJDRNQATu3bYM9iT8x07wIA2HD0Cl7cnIyc+w8lTkakv1hqiIgayNTIAB8G9MNXMwZBbmKIk9fvwW9NPA6n50kdjUgvsdQQETWSn6s9DgT7oL+DAveLy/HytpP4+MB5lFVwHEXUklhqiIiaQGdrc0Qu8sDLXk4AgC3x1zBlUxKyC4olTkakP1hqiIiaiImhAf4xoQ82zxwMS1NDpGbfh9+aeBxMuyV1NCK9wFJDRNTEfPvaISbEB4M6t0VRSQUW7TiN5fvSUFpRKXU0Ip3GUkNE1Awc2pnjh4UeWDiyGwBge9J1PL8hEZl3VBInI9JdLDVERM3EyECGZeN6I2zuUFi1MUbazUKMX5uA/WdypI5GpJNYaoiImtlTvWwQE+yDYU5WeFBageDvU7BszzmUlHMcRdSUWGqIiFqAncIU4a8MR/DTPSAIwPcnsjBp3TFcvl0kdTQincFSQ0TUQgwNZHjdtxd2zBuO9hYmyMgrwoS1x7D71A2poxHpBJYaIqIW5tWjPX4O8YF3j/Z4WF6JNyLP4PVdqVCVVkgdjUirsdQQEUmgg9wE218ehjd8nSETgD2nb2LiugSk5xZKHY1Ia7HUEBFJxEAmIOjpnohY4AE7S1NcyVdh0rpj+P5EFkRRlDoekdZhqSEiktgwJyvEhPhgVK8OKK1QY9mecwiOSEVRSbnU0Yi0ikalJjQ0FEOHDoVcLoeNjQ0CAgKQkZFR53GRkZFwcXGBqakpXF1dERMTU+3xPXv2wNfXF9bW1hAEAampqdUeLygowJIlS9CrVy+YmZmhc+fOCA4OhlKp1CQ+EVGrZdXGGFtnD8W7fi4wlAn46UwOxq9NQNpN/j1HVF8alZq4uDgEBgYiOTkZsbGxKC8vh6+vL1Sq2r8hMzExEdOmTcO8efOQkpKCgIAABAQEIC0trWqNSqWCt7c3Vq5cWeM5cnJykJOTg88//xxpaWnYtm0bDh48iHnz5mkSn4ioVZPJBCwY0R27FnmgU1szXL9bjOe+SsS2Y9c4jiKqB0FsxL8p+fn5sLGxQVxcHEaMGFHjmqlTp0KlUiE6Orpqn7u7O9zc3LBx48ZqazMzM+Hk5ISUlBS4ubk98dqRkZF46aWXoFKpYGhoWGfWwsJCKBQKKJVKWFpa1v3DERFJSFlcjjd3n8F/zucBAMb0tcVnzw+AwtxI4mRELUuT5+9Gvafmz/GPlZVVrWuSkpIwevToavvGjBmDpKSkxly66oerrdCUlpaisLCw2kZEpC0U5kbYNHMwVkzoA2MDGX75Iw9+a+KRknVP6mhErVaDS41arcbSpUvh5eWFfv361bouNzcXtra21fbZ2toiNze3oZfGnTt38OGHH2LBggW1rgkNDYVCoajaHB0dG3w9IiIpCIKAOV5O+PFVT3S2MsfN+w8xeWMStvx2FWo1x1FEf9XgUhMYGIi0tDREREQ0ZZ46FRYWwt/fH3369MGKFStqXbds2TIolcqqLTs7u+VCEhE1IVcHBaKDveHf3x4VahEfx1zAK9+eRIGqTOpoRK1Kg0pNUFAQoqOjceTIETg4ODxxrZ2dHfLy8qrty8vLg52dncbXLSoqwtixYyGXyxEVFQUjo9pnyyYmJrC0tKy2ERFpK0tTI6ybNhCfPOsKY0MZDqffht/qeJy4ViB1NKJWQ6NSI4oigoKCEBUVhcOHD8PJyanOYzw8PHDo0KFq+2JjY+Hh4aFR0MLCQvj6+sLY2Bj79++HqampRscTEWk7QRAwfXhn7Av0QrcObZBbWIJpW5Kx/shljqOIoGGpCQwMxI4dOxAeHg65XI7c3Fzk5ubi4cOHVWtmzZqFZcuWVf05JCQEBw8exKpVq5Ceno4VK1bg5MmTCAoKqlpTUFCA1NRUnD9/HgCQkZGB1NTUqvfd/FloVCoVvvnmGxQWFlZdu7KyslE3gIhI2/S2t8RPQd54bmAnVKpF/OuXDMwOO4H8olKpoxFJSqOPdAuCUOP+sLAwzJkzBwAwatQodO3aFdu2bat6PDIyEu+//z4yMzPRs2dPfPbZZ/Dz86t6fNu2bZg7d+5j512+fDlWrFiBo0eP4qmnnqrx2teuXUPXrl3rzM6PdBORLoo8mY1/7PsDD8sr0UFugi+nusGrR3upYxE1GU2evxv1PTXahKWGiHTVpbwiBIWnICOvCIIALHm6J0Ke6QkDWc3/IUqkTVrse2qIiEh6PW3l2BvohReHOkIUgTWHLmH6lmTkFZZIHY2oRbHUEBHpADNjA3z6fH+sftENbYwNcPxaAcatjsfRjNtSRyNqMSw1REQ6ZJJbJ0QH+6CPvSUKVGWYE/Y7Pv05HeWVaqmjETU7lhoiIh3j1L4N9iz2xCyPLgCAjXFXMHVTEm7ef1jHkUTajaWGiEgHmRoZ4INJ/bBhxiDITQ1xOus+/FbHI/Z8Xt0HE2kplhoiIh02ztUeMcE+GOCggPJhOeZ/exIf/HQeZRUcR5HuYakhItJxjlbmiFzkiVe8H30L/NZj1/DCxkRk3S2WOBlR02KpISLSA8aGMrw/vg++njUEbc2NcPaGEv5r4hFz7pbU0YiaDEsNEZEeGd3HFjHBPhjSpR2KSiuweOdpvL/3HErK+StnSPux1BAR6ZmObc3w/QJ3LB7VHQCwIzkLz36ViKv5DyRORtQ4LDVERHrIyECGt8a6YPvLw2DdxhgXbhVi/NoE7E25KXU0ogZjqSEi0mMjnTsgJsQHHt2sUVxWiaU/pOKt3WfwsIzjKNI+LDVERHrO1tIUO14ZjqWje0IQgF0nb2DiugRczCuSOhqRRlhqiIgIBjIBS0c7Y+crw9FBboJLtx9g4roE7Po9G6IoSh2PqF5YaoiIqIpn9/b4OcQHPj3bo6Rcjbd+PIvXfkjFg9IKqaMR1YmlhoiIqmlvYYLtc4fhrbG9YCATsDc1BxPXJuCPHKXU0YieiKWGiIgeI5MJWDyqB35Y4A57hSmu3lHh2a8S8V3ydY6jqNViqSEioloN6WqFmGAfPONig7IKNf6+Nw1B4SkoLCmXOhrRY1hqiIjoidq1McbXs4fgff/eMJQJOHDuFvzXxONM9n2poxFVw1JDRER1EgQBr/h0w+5XPeHQzgzZBQ/xwsZEfJNwjeMoajVYaoiIqN7cHNviQLAPxva1Q3mliA+jz2P+t6dwv7hM6mhELDVERKQZhZkRNrw0CB9O6gtjAxl+vZAHv9XxOHW9QOpopOdYaoiISGOCIGCmR1fsWeyJrtbmyFGWYMqmZGw4egVqNcdRJA2WGiIiarB+nRSIDvbBxAEdUakWsfJgOuZu+x13H5RKHY30EEsNERE1ioWJIVa/6IZPn3OFiaEMcRfz4bcmHslX70odjfQMSw0RETWaIAh4cVhn7A/yRg8bC+QVlmL6lmSs/vUSKjmOohbCUkNERE2ml50c+4O8MHmwA9Qi8O9fL2LmN8dxu7BE6mikB1hqiIioSZkbG+JfkwfgiykDYG5sgMQrd+G3Jh7xl/KljkY6jqWGiIiaxXODHLA/yBsudnLceVCGWVtP4PNfMlBRqZY6GukolhoiImo2PWwssDfQC9OHd4YoAuuOXMa0Lcm4pXwodTTSQSw1RETUrEyNDPDJs65YO20gLEwM8XvmPfitjsfh9Dypo5GOYakhIqIWMWFAR0Qv8YZrJwXuFZfj5W0n8UnMBZRVcBxFTYOlhoiIWkzX9m2w+1UPzPHsCgDY/NtVTNmUhOyCYmmDkU5gqSEiohZlYmiAFRP7YtPMwbA0NURq9n34r4nHwbRbUkcjLcdSQ0REkhjT1w4Hgn3g5tgWhSUVWLTjNJbvS0NpRaXU0UhLsdQQEZFkHK3MEbnIAwtHdAMAbE+6juc3JCLzjkriZKSNWGqIiEhSRgYyLPPrjbA5Q9HO3AhpNwsxfm0C9p/JkToaaRmWGiIiahWecrFBTIgPhnW1woPSCgR/n4Jle86hpJzjKKoflhoiImo17BVmCJ8/HEue7gFBAL4/kYWA9cdw+XaR1NFIC7DUEBFRq2JoIMP/+fbCty8PQ3sLY6TnFmHC2mPYfeqG1NGoldOo1ISGhmLo0KGQy+WwsbFBQEAAMjIy6jwuMjISLi4uMDU1haurK2JiYqo9vmfPHvj6+sLa2hqCICA1NfWxc2zevBmjRo2CpaUlBEHA/fv3NYlORERaxqdnB8SE+MCzuzUellfijcgz+L9dZ6AqrZA6GrVSGpWauLg4BAYGIjk5GbGxsSgvL4evry9UqtrfpZ6YmIhp06Zh3rx5SElJQUBAAAICApCWlla1RqVSwdvbGytXrqz1PMXFxRg7dizeffddTSITEZEWs5Gb4rt5w/F/f3OGTAB+PH0DE9clID23UOpo1AoJoiiKDT04Pz8fNjY2iIuLw4gRI2pcM3XqVKhUKkRHR1ftc3d3h5ubGzZu3FhtbWZmJpycnJCSkgI3N7caz3f06FE89dRTuHfvHtq2bVvvrIWFhVAoFFAqlbC0tKz3cURE1DokX72LkIgU5BWWwsRQhhUT++LFoY4QBEHqaNSMNHn+btR7apRKJQDAysqq1jVJSUkYPXp0tX1jxoxBUlJSYy5dp9LSUhQWFlbbiIhIe7l3s0ZMsA9GOndAaYUay/acQ3BEKopKyqWORq1Eg0uNWq3G0qVL4eXlhX79+tW6Ljc3F7a2ttX22draIjc3t6GXrpfQ0FAoFIqqzdHRsVmvR0REzc/awgRhc4binXEuMJAJ+OlMDiasTUDaTaXU0agVaHCpCQwMRFpaGiIiIpoyT5NZtmwZlEpl1ZadnS11JCIiagIymYBFI7tj10IPdGprhsy7xXjuq0RsT8xEI95RQTqgQaUmKCgI0dHROHLkCBwcHJ641s7ODnl5edX25eXlwc7OriGXrjcTExNYWlpW24iISHcM7tIOB4K98bc+tiirVGP5/j+waMcpKIs5jtJXGpUaURQRFBSEqKgoHD58GE5OTnUe4+HhgUOHDlXbFxsbCw8PD82SEhER/UVbc2NsnjkY/xjfB0YGAn75Iw/+a+ORknVP6mgkAUNNFgcGBiI8PBz79u2DXC6vel+MQqGAmZkZAGDWrFno1KkTQkNDAQAhISEYOXIkVq1aBX9/f0RERODkyZPYvHlz1XkLCgqQlZWFnJxHv+fjz+++sbOzq3pFJzc3F7m5ubh8+TIA4Ny5c5DL5ejcufMT36hMRES6TRAEvOzthCFd2yEoPAVZBcWYvDEJb491wTxvJ8hk/HSU3hA1AKDGLSwsrGrNyJEjxdmzZ1c7bteuXaKzs7NobGws9u3bVzxw4EC1x8PCwmo87/Lly6vWLF++vM5rP4lSqRQBiEqlUpMfmYiItIjyYZm4eMcpscvb0WKXt6PFuWEnxIIHpVLHokbQ5Pm7Ud9To034PTVERPpBFEXsPJ6FD6LPo6xCDXuFKdZMG4ihXfmqvjZqse+pISIiam0EQcBL7l2wd7EXurVvg1vKEry4ORnrj1yGWq0X/x2vt1hqiIhIJ/XpaImflnjj2YGdUKkW8a9fMjA77ATyi0qljkbNhKWGiIh0VhsTQ3wxZQA+e6E/TI1kiL90B35r4pF4+Y7U0agZsNQQEZFOEwQBU4Y4Yn+QN5xtLZBfVIoZ3xzHF7EXUclxlE5hqSEiIr3gbCvHvkBvTB3iCFEE1hy6hOlbkpFXWCJ1NGoiLDVERKQ3zIwNsPKF/vhyqhvMjQ1w/FoBxq2Ox9GM21JHoybAUkNERHonYGAnRC/xRm97SxSoyjAn7Hd8+nM6yivVUkejRmCpISIivdStgwWiFntipnsXAMDGuCt4cXMybt5/KHEyaiiWGiIi0lumRgb4MKAfvpoxCHITQ5y6fg9+q+MRez6v7oOp1WGpISIivefnao8DwT7o76CA8mE55n97Eh/89OgbiUl7sNQQEREB6Gxtjt2LPPGylxMAYOuxa3hhYyKy7hZLnIzqi6WGiIjov4wNZfjHhD7YMmsIFGZGOHtDCf818Thw9pbU0ageWGqIiIj+4m99bBET4oPBXdqhqLQCgeGn8f7ecygpr5Q6Gj0BSw0REVENOrU1Q8QCd7w6qjsAYEdyFgLWH8OV/AcSJ6PasNQQERHVwshAhrfHumD7y8Ng3cYY6blFmLA2AXtO35A6GtWApYaIiKgOI507ICbEBx7drFFcVonXd53BG5FnUFxWIXU0+v+w1BAREdWDraUpdrwyHK+NdoZMAHafuoEJaxOQnlsodTT6L5YaIiKiejKQCQgZ3RM7X3GHjdwEV/JVmLTuGMKPZ0EU+Ru/pcZSQ0REpCGP7tb4OcQHI507oLRCjXejzmHJ9ykoKimXOppeY6khIiJqAGsLE4TNGYpl41xgKBMQffYWxq9NwLkbSqmj6S2WGiIiogaSyQQsHNkdPyz0QKe2Zrh+txjPbTiGrQnXOI6SAEsNERFRIw3u0g4xwT7w7WOL8koRH0Sfx/xvT+F+cZnU0fQKSw0REVETUJgbYdPMwVgxoQ+MDWT49UIe/FbH49T1Aqmj6Q2WGiIioiYiCALmeDlhz2JPdLU2R46yBFM2JeOro5ehVnMc1dxYaoiIiJpYv04K/LTEGxMGdESlWsRnBzMwO+wE7jwolTqaTmOpISIiagZyUyOsedENnz7nClMjGeIv3cG41fFIvHJH6mg6i6WGiIiomQiCgBeHdca+QG/0tLFAflEpZnx9HF/EXkQlx1FNjqWGiIiomfWyk2NfkBemDHGAKAJrDl3C9C3JyFWWSB1Np7DUEBERtQBzY0N89sIAfDnVDW2MDXD8WgH81sTjSMZtqaPpDJYaIiKiFhQwsBN+WuKNPvaWKFCVYW7Y7wj9+QLKK9VSR9N6LDVEREQtrFsHC+xZ7IlZHl0AAJvirmLKpiRkFxRLnEy7sdQQERFJwNTIAB9M6ocNMwZBbmqIlKz78F8Tj4NpuVJH01osNURERBIa52qPmGAfDHBsi8KSCizacQor9v+B0opKqaNpHZYaIiIiiTlamSNyoQcWjOgGANiWmInnNyTi2h2VxMm0C0sNERFRK2BsKMO7fr2xdc4QtDM3QtrNQoxfE499qTeljqY1WGqIiIhakaddbBET4oNhXa2gKqtESEQq3vnxLB6WcRxVF5YaIiKiVsZeYYbw+cMR/HQPCAIQ8Xs2Jq1PwKW8IqmjtWosNURERK2QoYEMr/v2wo55w9HewgQX8x5gwroE7Po9G6LIX7FQE5YaIiKiVsyrR3v8HOIDn57tUVKuxls/nsVrP6TiQWmF1NFaHZYaIiKiVq6D3ATb5w7Dm2N6wUAmYG9qDiasTcAfOUqpo7UqGpWa0NBQDB06FHK5HDY2NggICEBGRkadx0VGRsLFxQWmpqZwdXVFTExMtcf37NkDX19fWFtbQxAEpKamPnaOkpISBAYGwtraGhYWFnj++eeRl5enSXwiIiKtJZMJCHyqByIWuMNeYYprd1R4dn0ivk3K5DjqvzQqNXFxcQgMDERycjJiY2NRXl4OX19fqFS1f44+MTER06ZNw7x585CSkoKAgAAEBAQgLS2tao1KpYK3tzdWrlxZ63lee+01/PTTT4iMjERcXBxycnLw3HPPaRKfiIhI6w3taoWYYB+M7m2Dsko1/rHvD7y64zSUD8uljiY5QWxEvcvPz4eNjQ3i4uIwYsSIGtdMnToVKpUK0dHRVfvc3d3h5uaGjRs3VlubmZkJJycnpKSkwM3NrWq/UqlEhw4dEB4ejhdeeAEAkJ6ejt69eyMpKQnu7u51Zi0sLIRCoYBSqYSlpWUDfloiIqLWQxRFbD2WiU9/voDyShEO7cywdtpADOzcTupoTUqT5+9GvadGqXw0y7Oysqp1TVJSEkaPHl1t35gxY5CUlFTv65w6dQrl5eXVzuPi4oLOnTvXep7S0lIUFhZW24iIiHSFIAiY5+2E3Ys84Whlhhv3HmLyxiRs/u0K1Gr9HEc1uNSo1WosXboUXl5e6NevX63rcnNzYWtrW22fra0tcnPr/wu7cnNzYWxsjLZt29b7PKGhoVAoFFWbo6Njva9HRESkLQY4tsWBYB/4u9qjQi3ik5h0zNv+OwpUZVJHa3ENLjWBgYFIS0tDREREU+ZpMsuWLYNSqazasrOzpY5ERETULCxNjbBu+kB8/Gw/GBvKcCQjH36r43H86l2po7WoBpWaoKAgREdH48iRI3BwcHjiWjs7u8c+pZSXlwc7O7t6X8/Ozg5lZWW4f/9+vc9jYmICS0vLahsREZGuEgQBM4Z3wd7FXujWoQ1yC0swbUsy1h66hEo9GUdpVGpEUURQUBCioqJw+PBhODk51XmMh4cHDh06VG1fbGwsPDw86n3dwYMHw8jIqNp5MjIykJWVpdF5iIiIdF2fjpb4Kcgbzw3qBLUIrIq9iFlbj+N2UYnU0ZqdoSaLAwMDER4ejn379kEul1e9n0WhUMDMzAwAMGvWLHTq1AmhoaEAgJCQEIwcORKrVq2Cv78/IiIicPLkSWzevLnqvAUFBcjKykJOTg4AVH33jZ2dHezs7KBQKDBv3jy8/vrrsLKygqWlJZYsWQIPD496ffKJiIhIn7QxMcQXU9zg2b09/r43Dccu34Xf6nj8e6obfHp2kDpe8xE1AKDGLSwsrGrNyJEjxdmzZ1c7bteuXaKzs7NobGws9u3bVzxw4EC1x8PCwmo87/Lly6vWPHz4UFy8eLHYrl070dzcXHz22WfFW7du1Tu7UqkUAYhKpVKTH5mIiEirXcorFMf8O07s8na02PWdaPGzgxfE8opKqWPVmybP3436nhptwu+pISIifVVSXokPos8j/HgWAGBIl3ZYM20gOrY1kzhZ3Vrse2qIiIio9TM1MsAnz7pi7bSBsDAxxMnr9+C3Jh6/ntetXzfEUkNERKQnJgzoiAPB3nDtpMD94nK88u1JfBh9HmUVaqmjNQmWGiIiIj3SxboNdr/qgZe9Hn2C+ZuEa3hhYyKy7hZLnKzxWGqIiIj0jImhAf4xoQ+2zBoChZkRzt5Qwn9NPA6cvSV1tEZhqSEiItJTf+tji5gQHwzu0g5FpRUIDD+N96LOoaS8UupoDcJSQ0REpMc6tTVDxAJ3LB7VHQCw83gWAtYfw+XbDyROpjmWGiIiIj1nZCDDW2NdsP3lYbBuY4z03CJMXJeAPadvSB1NIyw1REREBAAY6dwBP4f4wLO7NYrLKvH6rjP4v11noCqtkDpavbDUEBERURUbS1N8N284Xv+bM2QC8OPpG5i4LgEXbhVKHa1OLDVERERUjYFMQPAzPRE+3x22lia4kq9CwPpjCD+ehdb8iwhYaoiIiKhG7t2sERPsg1G9OqC0Qo13o84h6PsUFJaUSx2tRiw1REREVCtrCxNsnT0U7/q5wFAm4MDZWxi/JgFnb9yXOtpjWGqIiIjoiWQyAQtGdMeuRR7o1NYMWQXFeH5DIrYmXGtV4yiWGiIiIqqXQZ3bISbYB2P62qK8UsQH0ecx/9tTuF9cJnU0ACw1REREpAGFuRE2vjQYH0zqC2MDGX69kAe/1fE4db1A6mgsNURERKQZQRAwy6Mr9iz2RFdrc+QoSzBlUzK+OnoZarV04yiWGiIiImqQfp0UiA72wSS3jqhUi/gu6TqKSqT7oj5Dya5MREREWs/CxBBfTnWDV/f2cOrQBgpzI8mysNQQERFRowiCgClDHaWOwfETERER6QaWGiIiItIJLDVERESkE1hqiIiISCew1BAREZFOYKkhIiIincBSQ0RERDqBpYaIiIh0AksNERER6QSWGiIiItIJLDVERESkE1hqiIiISCew1BAREZFO0Jvf0i2KIgCgsLBQ4iRERERUX38+b//5PP4kelNqioqKAACOjtL/anQiIiLSTFFRERQKxRPXCGJ9qo8OUKvVyMnJgVwuhyAITXruwsJCODo6Ijs7G5aWlk16bvof3ueWwfvcMnifWw7vdctorvssiiKKiorQsWNHyGRPfteM3rxSI5PJ4ODg0KzXsLS05L8wLYD3uWXwPrcM3ueWw3vdMprjPtf1Cs2f+EZhIiIi0gksNURERKQTWGqagImJCZYvXw4TExOpo+g03ueWwfvcMnifWw7vdctoDfdZb94oTERERLqNr9QQERGRTmCpISIiIp3AUkNEREQ6gaWGiIiIdAJLTT2tX78eXbt2hampKYYPH44TJ048cX1kZCRcXFxgamoKV1dXxMTEtFBS7abJfd6yZQt8fHzQrl07tGvXDqNHj67z/xd6RNN/nv8UEREBQRAQEBDQvAF1hKb3+f79+wgMDIS9vT1MTEzg7OzMvzvqQdP7/OWXX6JXr14wMzODo6MjXnvtNZSUlLRQWu3022+/YcKECejYsSMEQcDevXvrPObo0aMYNGgQTExM0KNHD2zbtq3Zc0KkOkVERIjGxsbi1q1bxT/++EOcP3++2LZtWzEvL6/G9ceOHRMNDAzEzz77TDx//rz4/vvvi0ZGRuK5c+daOLl20fQ+T58+XVy/fr2YkpIiXrhwQZwzZ46oUCjEGzdutHBy7aLpff7TtWvXxE6dOok+Pj7ipEmTWiasFtP0PpeWlopDhgwR/fz8xISEBPHatWvi0aNHxdTU1BZOrl00vc87d+4UTUxMxJ07d4rXrl0Tf/nlF9He3l587bXXWji5domJiRHfe+89cc+ePSIAMSoq6onrr169Kpqbm4uvv/66eP78eXHt2rWigYGBePDgwWbNyVJTD8OGDRMDAwOr/lxZWSl27NhRDA0NrXH9lClTRH9//2r7hg8fLi5cuLBZc2o7Te/zX1VUVIhyuVzcvn17c0XUCQ25zxUVFaKnp6f49ddfi7Nnz2apqQdN7/OGDRvEbt26iWVlZS0VUSdoep8DAwPFp59+utq+119/XfTy8mrWnLqkPqXmrbfeEvv27Vtt39SpU8UxY8Y0YzJR5PipDmVlZTh16hRGjx5dtU8mk2H06NFISkqq8ZikpKRq6wFgzJgxta6nht3nvyouLkZ5eTmsrKyaK6bWa+h9/uCDD2BjY4N58+a1REyt15D7vH//fnh4eCAwMBC2trbo168fPvnkE1RWVrZUbK3TkPvs6emJU6dOVY2orl69ipiYGPj5+bVIZn0h1fOg3vxCy4a6c+cOKisrYWtrW22/ra0t0tPTazwmNze3xvW5ubnNllPbNeQ+/9Xbb7+Njh07PvYvEv1PQ+5zQkICvvnmG6SmprZAQt3QkPt89epVHD58GDNmzEBMTAwuX76MxYsXo7y8HMuXL2+J2FqnIfd5+vTpuHPnDry9vSGKIioqKrBo0SK8++67LRFZb9T2PFhYWIiHDx/CzMysWa7LV2pIJ3z66aeIiIhAVFQUTE1NpY6jM4qKijBz5kxs2bIF7du3lzqOTlOr1bCxscHmzZsxePBgTJ06Fe+99x42btwodTSdcvToUXzyySf46quvcPr0aezZswcHDhzAhx9+KHU0agJ8paYO7du3h4GBAfLy8qrtz8vLg52dXY3H2NnZabSeGnaf//T555/j008/xa+//or+/fs3Z0ytp+l9vnLlCjIzMzFhwoSqfWq1GgBgaGiIjIwMdO/evXlDa6GG/PNsb28PIyMjGBgYVO3r3bs3cnNzUVZWBmNj42bNrI0acp///ve/Y+bMmXjllVcAAK6urlCpVFiwYAHee+89yGT8b/2mUNvzoKWlZbO9SgPwlZo6GRsbY/DgwTh06FDVPrVajUOHDsHDw6PGYzw8PKqtB4DY2Nha11PD7jMAfPbZZ/jwww9x8OBBDBkypCWiajVN77OLiwvOnTuH1NTUqm3ixIl46qmnkJqaCkdHx5aMrzUa8s+zl5cXLl++XFUaAeDixYuwt7dnoalFQ+5zcXHxY8XlzyIp8lchNhnJngeb9W3IOiIiIkI0MTERt23bJp4/f15csGCB2LZtWzE3N1cURVGcOXOm+M4771StP3bsmGhoaCh+/vnn4oULF8Tly5fzI931oOl9/vTTT0VjY2Nx9+7d4q1bt6q2oqIiqX4EraDpff4rfvqpfjS9z1lZWaJcLheDgoLEjIwMMTo6WrSxsRE/+ugjqX4EraDpfV6+fLkol8vF77//Xrx69ar4n//8R+zevbs4ZcoUqX4ErVBUVCSmpKSIKSkpIgDxiy++EFNSUsTr16+LoiiK77zzjjhz5syq9X9+pPvNN98UL1y4IK5fv54f6W5N1q5dK3bu3Fk0NjYWhw0bJiYnJ1c9NnLkSHH27NnV1u/atUt0dnYWjY2Nxb59+4oHDhxo4cTaSZP73KVLFxHAY9vy5ctbPriW0fSf5/8fS039aXqfExMTxeHDh4smJiZit27dxI8//lisqKho4dTaR5P7XF5eLq5YsULs3r27aGpqKjo6OoqLFy8W79271/LBtciRI0dq/Pv2z3s7e/ZsceTIkY8d4+bmJhobG4vdunUTw8LCmj2nIIp8vY2IiIi0H99TQ0RERDqBpYaIiIh0AksNERER6QSWGiIiItIJLDVERESkE1hqiIiISCew1BAREZFOYKkhIiIincBSQ0RERI3y22+/YcKECejYsSMEQcDevXs1On7FihUQBOGxrU2bNhqdh6WGiIiIGkWlUmHAgAFYv359g45/4403cOvWrWpbnz59MHnyZI3Ow1JDREREjTJu3Dh89NFHePbZZ2t8vLS0FG+88QY6deqENm3aYPjw4Th69GjV4xYWFrCzs6va8vLycP78ecybN0+jHCw1RERE1KyCgoKQlJSEiIgInD17FpMnT8bYsWNx6dKlGtd//fXXcHZ2ho+Pj0bXYakhIiKiZpOVlYWwsDBERkbCx8cH3bt3xxtvvAFvb2+EhYU9tr6kpAQ7d+7U+FUaADBsisBERERENTl37hwqKyvh7OxcbX9paSmsra0fWx8VFYWioiLMnj1b42ux1BAREVGzefDgAQwMDHDq1CkYGBhUe8zCwuKx9V9//TXGjx8PW1tbja/FUkNERETNZuDAgaisrMTt27frfI/MtWvXcOTIEezfv79B12KpISIiokZ58OABLl++XPXna9euITU1FVZWVnB2dsaMGTMwa9YsrFq1CgMHDkR+fj4OHTqE/v37w9/fv+q4rVu3wt7eHuPGjWtQDkEURbHRPw0RERHpraNHj+Kpp556bP/s2bOxbds2lJeX46OPPsK3336Lmzdvon379nB3d8c///lPuLq6AgDUajW6dOmCWbNm4eOPP25QDpYaIiIi0gn8SDcRERHpBJYaIiIi0gksNURERKQTWGqIiIhIJ7DUEBERkU5gqSEiIiKdwFJDREREOoGlhoiIiHQCSw0RERHpBJYaIiIi0gksNURERKQT/h/o03QTlOeeQgAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sns.lineplot(x=J_all3[0], y=J_all3[1])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e8b1f648",
   "metadata": {},
   "source": [
    "## $R^2$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "id": "50022cc2",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "erklärte Varianz (R^2): 0.3520362618371272\n"
     ]
    }
   ],
   "source": [
    "X = feature_matrix_from_data(data[features])\n",
    "y = data.Price.to_numpy(copy=True)\n",
    "J_ana = J(w=w_ana, X=X, y=y)\n",
    "MSE = 2*J_ana\n",
    "mu_y = sum(y)/len(y)\n",
    "sigma_y_quadrat = ( (y - mu_y) @ (y - mu_y) ) / len(y)\n",
    "R2 = 1 - MSE/sigma_y_quadrat\n",
    "print('erklärte Varianz (R^2): {}'.format(R2))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "104ad3d4",
   "metadata": {},
   "source": [
    "$R^2$ ist größer als beim Modell mit nur 1 Feature (BuildingArea)."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.11.4"
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}