notebooks/tensorflow_linear_regression.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "tensorflow_linear_regression.ipynb",
      "provenance": [],
      "collapsed_sections": [],
      "authorship_tag": "ABX9TyOn9SRNlj/aTpmBusIyZmOQ",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/letianzj/QuantResearch/blob/master/notebooks/tensorflow_linear_regression.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "AH_tdJZQQqRC",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "import numpy as np \n",
        "import pandas as pd\n",
        "from matplotlib import pyplot as plt \n",
        "import tensorflow as tf"
      ],
      "execution_count": 1,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "BMkHwIeqThuR",
        "colab_type": "text"
      },
      "source": [
        "### 0. Data Preparation"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "E77lkB1URQ_0",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 102
        },
        "outputId": "761372b2-5c8b-43cb-bacc-e088029e20e4"
      },
      "source": [
        "tf.random.set_seed(0)\n",
        "np.random.seed(0)\n",
        "\n",
        "\n",
        "sample_size = 500\n",
        "batch_size = 10\n",
        "sigma_e = 3.0             # true value of parameter error sigma\n",
        "xs = 100 * np.random.rand(sample_size)\n",
        "es = np.random.normal(0, sigma_e, sample_size)\n",
        "ys = 1.0 + 2.0 * xs + es\n",
        "\n",
        "ds = tf.data.Dataset.from_tensor_slices((xs, ys)).shuffle(buffer_size=sample_size).batch(batch_size)\n",
        "for x, y in ds.take(1):         # take first batch\n",
        "    print(x, y)"
      ],
      "execution_count": 119,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "tf.Tensor(\n",
            "[94.53015335 42.40889884 42.40322519 84.90383084 43.44166256  2.01075462\n",
            " 55.20782767 40.71832972 49.04588086 10.20448107], shape=(10,), dtype=float64) tf.Tensor(\n",
            "[189.93610426  85.85229439  84.8538211  170.60625922  86.47141214\n",
            "   3.74855637 115.39504027  81.28272318  97.11208653  21.64374268], shape=(10,), dtype=float64)\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "_am7GSBUTG-X",
        "colab_type": "text"
      },
      "source": [
        "### 1. Low-Level Implementation\n",
        "\n",
        "The convergence depends on initial states of w and b, batch_size, learning rate, and loss function (mse or rmse), etc."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "8WiG4duHSjmp",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 102
        },
        "outputId": "742717a9-5d06-4868-9de6-7f55a30126d8"
      },
      "source": [
        "w = tf.Variable(tf.random.normal(shape=[1], dtype=tf.float64))      # scaler, shape=[] or shape=[1,]\n",
        "b = tf.Variable(0,  dtype = tf.float64)          # scaler\n",
        "\n",
        "epochs = 100\n",
        "learning_rate = 0.001\n",
        "ws = []\n",
        "bs = []\n",
        "for epoch in tf.range(1,epochs+1):      # for each epoch\n",
        "  for x, y in ds:           # for each batch\n",
        "    with tf.GradientTape() as tape:\n",
        "      y_bar = x*w + b                    # broadcasting\n",
        "      loss = tf.sqrt(tf.reduce_mean( (y-y_bar)**2 ))\n",
        "      # loss = tf.reduce_mean( (y-y_bar)**2 )\n",
        "    # Back propagation to calculate gradient\n",
        "    dloss_dw, dloss_db = tape.gradient(loss, [w, b])\n",
        "    # apply_gradients\n",
        "    w.assign(w - learning_rate * dloss_dw)         # assign, keeps tf.Variable; operator ==> tf.Tensor\n",
        "    b.assign(b - learning_rate * dloss_db)\n",
        "    ws.append(w.numpy()[0])\n",
        "    bs.append(b.numpy())\n",
        "\n",
        "  if epoch % 20 == 0:\n",
        "    tf.print('w=', w[0], ' ,b=', b, ' ,mse=', loss)"
      ],
      "execution_count": 120,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "w= 2.0092433160037424  ,b= 0.083804579058864792  ,mse= 2.5152901410985193\n",
            "w= 2.0133854164832647  ,b= 0.12503606117564978  ,mse= 3.4116469321528675\n",
            "w= 1.9745402878296359  ,b= 0.16570829832847597  ,mse= 2.181400506321292\n",
            "w= 2.0332908978757693  ,b= 0.20492633276977928  ,mse= 1.9839655160821144\n",
            "w= 1.9998285047327671  ,b= 0.23774513013090734  ,mse= 4.39657509939333\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Z7cukP17ZVRu",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 382
        },
        "outputId": "c6c8b4e7-c780-4573-ceac-cfaef0c3c146"
      },
      "source": [
        "fig, axes = plt.subplots(nrows=1, ncols=2, figsize=(12, 5))\n",
        "fig.tight_layout()\n",
        "axes[0].plot(ws)\n",
        "axes[1].plot(bs)"
      ],
      "execution_count": 121,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "[<matplotlib.lines.Line2D at 0x7ff96d9d2a58>]"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 121
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1kAAAFcCAYAAAA6U3SpAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nOzdd5xU1f3/8ffZTu9N2tIRUBGWpogg0o0YS77WiImS2JNYAjYUC8TEJPrTaIixxd5FQRFFFBGRrhRRytKRvpTtc8/vj50ddnZnYded3Tt77+v5ePBg7rl3Zj57KbPvPed+rrHWCgAAAAAQHXFuFwAAAAAAXkLIAgAAAIAoImQBAAAAQBQRsgAAAAAgighZAAAAABBFCW4XUJrGjRvb1NRUt8sAALhgyZIle6y1Tdyu43j4rAIA/zrWZ1XMhqzU1FQtXrzY7TIAAC4wxmxyu4ay4LMKAPzrWJ9VLBcEAAAAgCgiZAEAAABAFBGyAAAAACCKCFkAAAAAEEWELAAAAACIIkIWAAAAAEQRIQsAAAAAooiQBQAAAABRRMgCAAAAgCgiZAEAAABAFBGyAAAAACCKCFkAAAAAEEWELJf98NMhnXj3R9p2IMvtUgAAAADPOpyTrz+8ukypE2bos+93Vep7EbLKYO3OQ8rMza+U13554WZl5QU0a+XOcj0vOy+gf81dp/yAUyl1VVc7M7J1OCf6f1YHMnOV5/Nzvfdwjp6et0HWWrdLqVK5+Y6+WrfH7TIAAMDPlJGVp7QHPlGPSbP07vLtkqRGtZMq9T19E7Icx2pHRpZWbc/Qi19vUkZmnnZmZCvgHPsbxj2HczTin1/oupeWKuBYfb/zYKnHHs7J1z8/+aFcwScjK0+SlJPvKCs3IElauS1DOfmB4HhA63cfDh3/5Nz1Sp0wQ/+Y/YMe/mitXliwSU/P26D7P1itI8FwYa0NfSOcF3C0blfB87PzAnpryVblBRxl5uZryN/m6puN+8LqsdZq9uqf5ATPy8zvdigjM++YX0NewDnm1/ze8m26+dVlxz0XUz5co9mrfwpt7zmcE/rzWbZ5f6imY+k/5VOd89i8qAWBXYey9e6ybeo5ebb++NryqLxmoR0ZWcrOCxzzmAOZuWX6ussjN7/kn9XB7DzN+3G3np63QX/56Hut3XlIUvif7Z9eX6EHZqzRG4u36rO1R3/6cyQnX+t2HSr1/fYdydVH5fwhQmm27s8sMev73vJtWpRe8Pd4wfq9en3xljK/Xm6+o637M495zF8++l6XPr1QK7YciLg/fc8RHcw+9r+RsrjupSV6b/m2UvcX/XcdDYdz8rW9DDPos1f/pMnvr47a+wIAUJWy8wI65b6PtedwTtj4ya3qV+r7mlj9qXRaWppdvHhx1F7vic/W6a+z1kbc9+Jv+2lgp8ah7S37MnXGw5/p2av66KEZa/TjrsNKSYzT7wZ10KOf/qgPbhyo5IQ4PfHZOr27fLsGtG+kO8ecqPeWb9N/5m2UJF03uINuG9FFs1b9pHW7DunawR0VH2ckSWt2HNTctbs1flB7dbhjZlgtX08cqv5TPtWv0lrp4QtP0S2vr9BbS7fqkYtO0fm9WqrdxILjx/Y8Qe8Fk3ihGonxevHqvrrgyQVq26imPr9tiCa9t1LPL9ikyWO76573VoWOfevaAbrgyQWSpPSpY5SZm6+UhHi9vWybbn1jhSSpVlK8jgSD34aHRmv34RzlBRz9+NNhXfXcIt015kSd3rGxRj06Tz1a1lXbRrW0ZsdBvXpNfzWunaxDOfk65b6PQ+/5zZ1D1ffBT3X1wHa665xukqSFG/bq1DYNlBBn1D54LtKnjtHOjGz1n/KpkhPiNLBjY336/S71b99Qp7Sur4mjTpRU8E3nTwdz1LxeSug9UifMCD1edOfZWrJpn47kBJSW2kDN6qYoJTE+7M951faD6teuoU69f3bo78IpreupTkqiHpyxOvTnWSh96hjl5Ad0JCeghrWO/gRkcfo+XfjUAl11eqom/aK7JGn97sNq37iWRj06Tye3qqeHLzwldPw9763UCws2qV3jWnrlmv5qXi9Fy7cc0HlPzNeoHs3VplFN/SqttYY+8rluPKujujSvo0/X7NIjF52iuDijgGNDf5+K2ns4R5+s+UmPfbpOBzJztWrySEkFQeRgdr7O6NhYg/82V1cPbKc7x5yoJz9fr/N6ttR5T8zXrkPh//m0bVRTm/ZmqkOTWnr72tP162e/CQsa6VPH6NM1P+mmV5bpSG5AG6eM1sOz1qpmYrwGd2mq+Dijwzn5unf6Kq3ecTD0nEj+PvsH7T6Uo98OTNWugzm67c1v9dglp+qCJ7/Sa+P7q1/7RmF/vmsmj9TnP+yWY62ue2mpJGncaal67qv0gnP/0GjtyMjSmh2HNKxbM416dJ6GdGmiX57aUi0b1FDNpAQtTt+ny55eqJx8RyvvG6HXF23Rr/q0Vu3khNCf30crd2r+uj36av1e3TKss64b0lFrdhxUreQEHcnJV8emtdX17o/UpmFNXdS7la4Z1D7s79juQzl6b/k29W/fSD1a1ov4tS/ZtF+vLdqs1xdvLfUcZWTmafDfPtORnIB+eHCUMnPz9dqiLbpyQKri4kwoiAes1axVOzW4S1Nl5Qa0PzNX327N0K1vrNCyu4fpp0PZ6ty0juLijM7+++dat+tw2PvtOZwjI6lR7WRJ0gsL0kP/bxQeZ63VviO5oWMqizFmibU2rVLfJAqi/VkFAIiuot8bFvritiFq06hmhV/7WJ9VvglZv/zXfC3bHPkn0ZKU2qimkhLi9MNPR2eN6qQk6FB22ZaeNa2TrDEnt9Cz89NLPeb2kV1UOzkhLOwcS7vGtbQzI1tZEWY7zj+1pd5eVvpPvSVpzi1n6toXl2rtT6XPMkjSb05vp2fmbzzmMcXr2rjnyDGPaV43RTsPZpe6/9J+bfTyws0R97VtVFMjujfXtC82RNw//YbTdXKr+qF/NI9cdIpqJSdoUOfG6nbPrGPW9c0dQ9X3oU+PeYwkvXv96TrvifklxpffM0w9JxcEsoV3DNXM73bovmI/5X/ysl5avvWA/v35Bg3r1iw0O9evXUMtLDZzWKhn6/rq2bp+KCSUpn7NRF3Wr42e+Gy9Rp/UXA+cd5IOZ+frpW82affBHK3fcyQsCJ3cqp5evLqfTr63IOz+X1prvVaOmZ7y6Ni0dmjWtDRvXTtAHZrUVkJ8nA5n56v/lNL/LM7o1FjzfixYpvfq+P6qmRSvcx8v+DM51rksrkfLulq5LXwG+tcD2uqFBZtC29cP6aAnPlsvSWpZv4auG9JBd76zssRr1U5OOOZy1MK6Zt50hto3qaWud38Utv/y/m2051CuNu45orU/HSpRhyT9/VenyBjpj6+t0N9/dYpOblVPZ//9i9D+MSe30Ixvd0iSOjWtrbyAo/S94bNxI7o306xVPymSE+qlaPafzlT3SQX/Vto2qqlGtZK0tMj/jxunjNY1LyzRJ2siv0ahm4Z2UkZmru49t7u+3ZqhfUdydWbnJjJGMqbkDwHKg5AFAKiIXz/zjb74YXdo+9lxfTSka9OovgchS9L5/5of9k0Eqre0tg20eNN+t8sAIKlx7STtOZwb2r5+SAfdNqJrhV6TkAUAKK9tB7J03UtLVSc5QV8Wu566tNU0FXGsz6qEKLx4a0kvSGomyUqaZq19tNgxRtKjkkZLypQ0zlq7tKLvXR7Hm3lB9ULAAmJH0YAlSc9/tanCIQsAgLL68sc92p+ZqxtfKdkD4PPbBqtto1pVXlOFQ5akfEm3WGuXGmPqSFpijJltrS26hmqUpE7BX/0kPRn8vcrsP07zBgBAdFRGh08AACLZvDdTl/93YcR9lTF7VVYV7i5ord1ROCtlrT0kaY2klsUOGyvpBVvga0n1jTEtKvre5dG0TuVepA0AAACgamTnBXTjK8s06K+fhY2P7XmCJGnVfSPcKCskGjNZIcaYVEmnSioeJ1tKKnq1/dbg2I5izx8vabwktWnTJpqlKSFCJzYAAAAA1ceKLQf0zPyNJbps/25Qe90yvIuSEuL06MWnulTdUVG7T5YxpraktyT9wVpb+s2kjsFaO81am2atTWvSpEm0SpMkndax8fEPqmS3jehS5mNn3DTwuMc0rp2k0Sc1172/6BYae/KyXj+rtmhbevcwt0uIeQ1qJkYcP//U4hPBKKuW9Wv87Od2aVYnipWEq0hdZTVxVOxcA7UmeOsAAACiadoX6zX2ifklAta0K3pr4ugTlZQQO7cAjkolxphEFQSsl6y1b0c4ZJuk1kW2WwXHqkyN4L1rLujVSpI0qHMTrX1gpN69/nTdfU43vX3daTqxRV1J0imt6+v8U1uqUfA+SLP+MCj0OulTxyh96hiNOy017PW/vXe4nr2qj84LTlFGkpwQp/SpY/ToxT119onNVDMpXs+O66NzTm6hZnULljM+dXkv3fuLbup+Qj1teGi0JOmaM9pp1X0jdP95PdSjZUGN405L1Td3nK1/XdZb405vp9fG99dzV/XRqJNa6L3rT4/4/k9d3kvpU8fokz8Nirg/kppJBectJTFOLYL3o5o/4Sytum+Enrg0PNB1bX70m9Si95Aq6sazOoYe331ON00qEhAl6YMbj4bL2X8cpLvGnKj0qWP05u8HhMY7N6td4j3bNa6lSb/ops9uHSxJqpOcoG/vHa7//bZviTolaVSP5rphSEclxBk9/5u+6tq8jn47sJ2uH9JBj17cUwODobxuSsnJ3q8mnKUFE88KG3v/hoFKjA+fLU2fOkb/l1bw135ghJBftK71wT9rSfr7//XU6snln+Ju36SWfnxwlNKnjtGz4/qU2H9R71YlxvqmNgzbLrxHVHHJCXH6z6/TtPTuYXp2XB+9fHU/rXtwVGhKvlZSfMTnFbp9ZMEPGM7o1Fh3jC5bGLipyN8VSVow8SzdNLRTaLtxsfs0/e+3fTV/wlnq0KTg4tY2DWvq31f0liT94exOob9bU84/SfNuH6IbhoS//sybz9DLV/dTckKcLuvXRuNOS1X9molKii/4b3LDQ6PVvnH4hbPpU8eodcOSAeqxS47+BO2B83poxk0DVafIuf323uER/6+ItKx5xaThuumsjqF9NYud64a1knT/2O763ZkdQv8/FGrdsIauOaNd6DyUVaNaSXrq8t56/4aB+nriUM27fYhWFlt2UTs5ocS/g0I1jvP3AQCA8li9/aBSJ8zQQzO/Dxt/bXx/rbxvhIZ3b+5SZaWrcAv3YOfA5yXts9b+oZRjxki6QQXdBftJesxa2/dYrxvttrgT3/5On6z5SYvuPLvUY15euFl3vPOdVtwzXPWKzTKM+McXWvvTodAFdI5jlRtwtHbnITWqnaRWDQpuaPbfLzfq/g9W661rT1OzusmqWyNRk99frTeXbNWkX3TTVae3K1fdjmPD7jkT6Qa8kTzz5Ub1bFNfJ7Wsp053figp/OK/KTPXqG+7hsrOc3T9y0v16MU9NbZny9C9p/51WS8N6txESfFxWrvzkLqfUFcfr96pu95dqfkTzlJyQnyonodnrdWTc9froV+epAt6t1R+wKpW8BtKa62e/ypdw7o3V3JCnBrXTg69R/rUMQo4Vi9+vUmTph+94WlewFFmTqDEn8GWfZn6f3N+1APnnaTLn16ocaen6szOTbTvSK5aNzx6Q7ntB7LUvG6K4oosEb3plWXq1LS2Hpn9g64b3EG3jzz+N/rzftytvu0a6u8f/6Bxp6fq5leW67L+bTS2Z8FM01fr9+ilrzfrvFNbali3ZpKkHpNmhS76L36x5bhnv9Hctbt1/9juuvu9VVp4x1AdzslXUnycWjesGXZeJOmGl5eqUa0k/XFYZxkZ1auZqPW7D2voI5/r5Wv66eNVP4XdV+ui3q3014sKbnh8OCdfPSYdvWdY4Wuu3JahBev3akjXpurYtCCsXv70Qq3cnqG/XniKerWpr7eWbtXFfdvonndXatIvuisn31FKYpzq14wcnAs9N3+j7n1/tV74TV/1bttAR3LyNfqxeXr5mv7q3KyONu45orYNayouzmjhhr36av1e/apPa50+dU7Y64zo3kzXnNFep7ZpEPGGy0UtWL9Xy7cckGOtri8WmsriSE6+5q/bo417juh3Z3Y47vF5AUeLNu7TpU8vVPvGtTTn1sHadShbfR/8VL87s71SEuLVskENXdCrle585zsNPbFZ6O/Gk3PX6y8ffa/bRnQJ1bpyW4baN6ml1xZt0YW9W6l2coK+2bhPDWolqXOxmbVD2Xl6aeFmXT2wnTLzAkpJiFd8nAk7Rzn5AXW5q+DeXB/cODB0A+TsvIBueX2FbhvRRamNa+l/X2/Syws3a03wJtF1UxJ0MHhPwPUPjY543h3Hhm4YLkn3nNNNvxnYLvT39penttQ7y7ZpwcSz1KJexWfuaOEOAAg4Vtl5AQ2Y8mnoc0qSnr2qj87s1CTsez03VOp9sowxAyXNk/SdJCc4fIekNpJkrX0qGMQelzRSBS3cr7LWHvNTKdofXLe/uULzftyjBROH/qznZ+UGlJMfOO43mtZardt1WJ2KfIO0MyNbf3p9uZ68rHeJ4FAVin/zfiyHsvN0MDu/3Mub1uw4qK7N65TpBqRb92eqfs2ksFmTaV+sV9fmdTWoc3SXiRb3/c6D6hi8GW5lyM4L6PE563R5/7YlgnBWbkDbDmSpQ5NaynesEovVkDphhhrWSirzUsvsvIDeXLJV2XkBHc7J1+/P7KCUxKMzCD/8dEhb9mWqc7M6YSG0slhrtWVfVrnvoL4zI1sLNuzRmh2HtHJbhl6+pn8lVRgdGVl5OuW+j3X/eT10Rf+2ZX5efsDRJ2t+0ojuzSt8o95j2X8kVzn5znF/EGOtVU6+o8JSNu45oq/X79W4Y/wg6PqXlur8Xi019MRmobFbXl+ht5ZujXoHJ0IWAPjbviO56nX/7LCxuikJmnPr4BKrWdzCzYgl/en15fpm4z59+efIy1u8bGdGtg5m55X4yThiy7dbD6h5vRQ1rXPsb44BPyBkAYD/FK76KI2bLdkjqdSbEVcXAcf6tsNg83opx/2pNtx3cqv6bpcAAABQ5f63IF13v7cq4r6W9Wto75EcLb6rejVV803Iynfsca/vAAAAAFB1Pvh2e6kBa83kkdW2mZJvQlYgYJUQFzttHQEAAAA/27jniG54eVnY2NoHRoYarFVnvglZzGQBAAAA7jqUnadrXlisrzfsCxt/5Zr+GtChkUtVRZ+PQpajhHhCFgAAAOCGJZv26YInF5QY/+RPg9SxqbcatPkmZAWYyQIAAACqnONYTZu3QVM/LNk58OELT/ZcwJJ8FLLyA/7tLggAAABUtVvfWKE3l2wtMd6lWR39d1yaWjWo/Pt4usU3IaughTuNLwAAAIDK9PS8DXpgxpqI+1ZMGq56NRKruKKq55uQle84Sk70zZcLAAAAVLqAY9XhjpmSpGsHd9CQLk0jBqxL+rbRlPNPquryXOOb1ME1WQAAAEB0FQYsSXpy7no9OXd9aPv8Xi1101mdlNq4lhuluco3ISvf4ZosAAAAIFpmfLuj1H3L7h6mBrWSqrCa2OKbi5SYyQIAAACiw3Gsrn95qSRpwqiu2jhltK4d3EGS9MGNA30dsCTfzWT5JlMCAAAAUXMwO09PfLZOt4/oGrZEUJJ+f2ZBuPrzyK7688iubpQXc3wTspjJAgAAAMrv9cVbdPub30qS/v35hrB9vx7Q1o2SYp5vpnbyHYeQBQAAAJTD5z/sDgWs4i7u01qTx/ao4oqqB9/MZDmOFGcIWQAAAMCxWGvVbuLMEuNTzj9JE9/+TlPOP0lndW2qZnVTXKiuevBNyCpYLuh2FQAAAEBsG/PYlyXGvr9/pFIS43VJ3zYuVFT9+CdkWa7JAgAAACLJzXf0xQ+7dfULi0NjZ3ZuolXbM/Ti1f2UkhjvYnXVj29CluNYlgsCAAAAxVhr1fmuD8PGWtavoed/09eliqo//4QsZrIAAACAkLyAo7U7D+mJz9aV2PfOdae5UJF3+CZkBZjJAgAAACRJE976Vq8u2lJi/P0bBqpTs9osD6wg34Qsx9JdEAAAAHhv+bYSASs5IU5rHxjlUkXe45uQRXdBAAAA+N3XG/bq5leXh7b7pDbQGZ2a6KahnVysynv8E7KsVRzXZAEAAMBHjuTk679fblTdlARt2Z+l/365MbRv/UOj6VlQSXwTshzHKp7lggAAAPC4fUdyNef7XVq+Zb9e/HpzxGPm3HImAasS+SZkcZ8sAAAAeN38dXt02dMLj3nMugdHKYHraCqVL0KWtVaWxhcAAADwoG+3HtC5j89Xz9b1tXzLgbB9J7eqp8v7tdW32w5o8rk9ZIxk+J640vkiZAUcK0nMZAEAAMBzzn18viSFBaxnx/VR1xZ11KJeDUnSr/q0dqU2v/JHyLKELAAAAHiLtVbnP/lVifFrzminIV2bulARCvkiZAUzlpgZBQDvM8aMlPSopHhJT1trpxbb/ydJV0vKl7Rb0m+stZuC+66UdFfw0Aestc9XWeEAUA55AUed7vwwtN2+SS09f1VfHcrOV7cT6rpYGSTJF1e8hZYLkrIAwNOMMfGSnpA0SlI3SZcYY7oVO2yZpDRr7cmS3pT0cPC5DSVNktRPUl9Jk4wxDaqqdgA4nl2HsnXnO98pMzdfve+fHbZvzi2D1bphTQJWjPDFTBbLBQHAN/pKWmet3SBJxphXJY2VtLrwAGvtZ0WO/1rS5cHHIyTNttbuCz53tqSRkl6pgroB4JjeW74tdBPhlxaGt2X/7t7hbpSEY/DFTJYTnMmiuyAAeF5LSVuKbG8NjpXmt5IK19uU6bnGmPHGmMXGmMW7d++uYLkAcHzW2lDAKm7DQ6NVJyWxiivC8fgiZNFdEABQnDHmcklpkv5anudZa6dZa9OstWlNmjSpnOIAIGjP4Ry1mzgztH3POQUroMedlqpldw9THN/fxiRfLRfkLyEAeN42SUX7FLcKjoUxxpwt6U5JZ1prc4o8d3Cx586tlCoBoIzSHvgk9PgPZ3fSbwa2028GtnOxIpSFL2ayHKfgdxpfAIDnLZLUyRjTzhiTJOliSdOLHmCMOVXSvyWda63dVWTXLEnDjTENgg0vhgfHAKBK5QccfbRyp1InzAgb/8PZnV2qCOXlq5mseF9ESgDwL2ttvjHmBhWEo3hJz1hrVxljJktabK2droLlgbUlvWEKfvi22Vp7rrV2nzHmfhUENUmaXNgEAwCqyrpdh3X23z8vMb7hodEuVIOfyxchi8YXAOAf1tqZkmYWG7unyOOzj/HcZyQ9U3nVAUBkby/dqi9/3KO3l4WvcL7xrI66ZXgXl6rCz+WPkGUJWQAAAIhNv3lukeZ8vytsrFWDGpp+w0A1rJXkUlWoCF+ELLoLAgAAINYMmPKpdmRklxi/on9b3X9eDxcqQrT4ImQ5dBcEAABADHh98Rad1LKe3l+xvUTAeuP3A7Ru12Fd0reNS9UhWnwRsgJ0FwQAAIDLincLLGrN5JGqkRSvPqkNq7AiVBafhCy6CwIAAMA9eYU/9S8mfeqYKq4EVcEXsYPGFwAAAHCDtVZfrdujTnd+WGLf4rtKbXaKas5nM1mELAAAAFSNP7y6TO8u3x42NueWM9WucS1JkmECwLOiMpNljHnGGLPLGLOylP2DjTEZxpjlwV/3RDqusgRofAEAAIAq9PeP15YIWOf1PEHtm9SWMYaA5XHRmsl6TtLjkl44xjHzrLXnROn9yqXwZsQ0vgAAAEBlK97g4oR6KXrrutPUol4NlypCVYtKyLLWfmGMSY3Ga1WGYMZiuSAAAAAqTVZuQL/81/ywsY1TRjNr5UNV2fhigDFmhTHmQ2NM90gHGGPGG2MWG2MW7969O2pvXHhNFn+/AQAAUFmue2mJvt95KLS94SECll9VVchaKqmttfYUSf9P0ruRDrLWTrPWpllr05o0aRK1Ny/sLshyQQAAAETLkZx87TmcI0l65ZvN+mzt0UmCf13Wi34APlYl3QWttQeLPJ5pjPmXMaaxtXZPVbw/3QUBAAAQTQez83TyvR9H3PfdvcNVJyWxiitCLKmSmSxjTHMTnCs1xvQNvu/eqnhvie6CAAAAiJ6AY0sNWPNuH0LAQnRmsowxr0gaLKmxMWarpEmSEiXJWvuUpAslXWuMyZeUJelia4PJpwrQXRAAAAAVlR9w9NX6vfrja8sj7v/xwVFKjK/KlgeIVdHqLnjJcfY/roIW765guSAAAAAqYvBfP1P63sywsbm3Dlb9mokyMqqVHK8EAhaCquSaLLcVNr6IYyYLAAAA5fT5D7tLBCxJSm1cy4VqUB34ImQFnILfmckCAABAWa3deUgj/vlFifFh3ZrpiUt7uVARqgt/hKzCFu7M4AIAAKAM9hzOKRGwuLEwysoXIauwxwb/KAAAAFAax7Ea9o/PtX73kRL7RvVozveSKDNfhKwA3QUBAABwDEdy8tV90qwS4yO6N9O/r0hzoSJUZ75YQEd3QQAAAJQmKzcQMWBJImDhZ/HFTJbDzYgBAABQive/3R56/J9fp2lYt2YuVgMv8EXICnUXZLkgAAAAgk6d/LH2Z+aFtts3qUXAQlT4Y7lgaCbL5UIAAAAQE75atycsYEnSnFsGu1MMPMcXM1kOjS8AAAAgadaqnVq6eb/+/fmGsPGR3Zu7VBG8yBchi8YXAAAAeHreBj0wY03YWGK8Ub0aiXrqit4uVQUv8kXIovEFAACA/ziOVZ7j6L3l21U3JbFEwJpx00B1P6GeS9XBy/wVslguCAAA4Bvt75gZcfyL24aodcMa3FwYlcYXIYvuggAAAP5y4t0fRRz/z6/T1KZRzSquBn7ji5Dl0F0QAADA89L3HNGKrQd086vLw8ZTG9XUqW0aaGzPEzS4S1OXqoOf+CJkBeguCAAA4GnWWg3+29wS4+lTx1R9MfA9X8zt0F0QAADA29pNDL/+qnHtJH1//0iXqoHf+WImy7FWxoiLGwEAADxkz+EcpT3wSdjYiS3q6sObz3CpIqCAL0JWwLEsFQQAAPAIa22JmStJ6tq8DgELMSWvgDgAACAASURBVMEfywWt5R5ZAAAAHrF08/6I4wQsxApfzGQ5zGQBAAB4wvItB3TBkwtC2/ed213GSFf0b8ulIYgZ/ghZlqYXAAAA1dGq7RlqUjtZF0/7Whv2HAnbR+dAxCpfhKyAU9D4AgAAANWDtVaTP1itZ+enR9w/oH2jqi0IKAdfhCzHWmayAAAAqokvftitXz/zTan7F94xVM3qplRhRUD5+CJk0V0QAAAg9m3YfVhnPfJ5xH2X9G2jK09rq67N61ZxVUD5+SJkOXQXBAAAiHmRAtb3949USmK8C9UAP58vQhYzWQAAALGtwx3h971aMWm46tVIdKkaoGJ8ErLoLggAABBLDmTmqk5KopZvOaDvth5QwLGhfXQNRHXni5BVsFzQ7SoAAAAgSVv3Z2rgXz4rMd6oVpKW3D3MhYqA6PJF9GC5IAAAQGzYczgnYsCSpA9vPqOKqwEqhz9CFo0vAAAAXDd9xXalPfBJxH292tRXU9qywyN8sVzQWqs4ZrIAAABcEXCsdh3K1k2vLAuNTTn/JP3ilBNUO9kX347CZ3zxt5rlggAAAO4p3jnwyct6adRJLVyqBqh8PglZYrkgAACAC4Y+Mjdse83kkaqRxH2v4G2+CFmOtYr3xdVnAAAAseFAZq56Tp4dNjbnljMJWPAFX4QslgsCAABUreIBa+EdQ9WMxhbwCV/M7zh0FwQAXzDGjDTGrDXGrDPGTIiwf5AxZqkxJt8Yc2GxfQFjzPLgr+lVVzXgLdZapU6YEdpuWb+GNjw0moAFX2EmCwDgCcaYeElPSBomaaukRcaY6dba1UUO2yxpnKRbI7xElrW2Z6UXCnhU+4kz5NjwsdtGdNH1Qzq6UxDgIt+ELGayAMDz+kpaZ63dIEnGmFcljZUUClnW2vTgPseNAgEvchyrW99cUSJgNamTrOsGd3CnKMBlvghZjrVKiPPFykgA8LOWkrYU2d4qqV85np9ijFksKV/SVGvtu5EOMsaMlzRektq0afMzSwWqv4Bj9egnP+ixOetK7EuMN1p059kuVAXEBp+ELImMBQA4jrbW2m3GmPaS5hhjvrPWri9+kLV2mqRpkpSWlmaL7wf8wHFsiXtfSdKq+0aoFjcXBvzR+CLgWMVxTRYAeN02Sa2LbLcKjpWJtXZb8PcNkuZKOjWaxQFekZGZpyuf/abE+IKJZxGwgCBf/EsouE8WIQsAPG6RpE7GmHYqCFcXS7q0LE80xjSQlGmtzTHGNJZ0uqSHK61SoBpyHKuAtTpl8sdh4yO6N9Pjl/ZSIjclBUJ8EbLoLggA3metzTfG3CBplqR4Sc9Ya1cZYyZLWmytnW6M6SPpHUkNJP3CGHOftba7pBMl/TvYECNOBddkrS7lrQBfah9heeDGKaNl+B4LKME3IYvuggDgfdbamZJmFhu7p8jjRSpYRlj8eV9JOqnSCwSqmT2Hc5T2wCclxv920Sk6sUUdAhZQiqjM6xpjnjHG7DLGrCxlvzHGPBa8OeS3xphe0XjfsnIsM1kAAADlsf9IbsSA9ey4Prqwdyt1P6GeC1UB1UO0Fs8+J2nkMfaPktQp+Gu8pCej9L5lEnC4JgsAAKA8Tr1/domx7+4driFdm7pQDVC9RGW5oLX2C2NM6jEOGSvpBWutlfS1Maa+MaaFtXZHNN7/eApauBOyAAAASvP20q3aeTBbJ9SroT+8tjw0PvX8k3RxX+4JB5RHVV2TFekGkS0lVUnIKmh8URXvBAAAUP08PudH/e3jH0qMP/jLHgQs4GeIqV6bxpjxxpjFxpjFu3fvjtrrOpbGFwAAAKWJFLAGd2miy/q1daEaoPqrqpBVphtEWmunWWvTrLVpTZo0idqbO9yMGAAAIKKXFm4KPe7Woq4kacZNA/XcVX3dKgmo9qpqueB0STcYY16V1E9SRlVdjyVJAboLAgAAhMkLOHpufroenLlGknT3Od3024HtXK4K8IaohCxjzCuSBktqbIzZKmmSpERJstY+pYJ7loyWtE5SpqSrovG+ZRVwaHwBAABQaPPeTA3662dhYwQsIHqi1V3wkuPst5Kuj8Z7/RyOtYqPqavPAAAA3JGdFygRsD68+QyXqgG8qaqWC7qqoLsgM1kAAMDfsvMC6nr3R6Ht53/TV2d2jt518AAK+GJ+x3HoLggAAPzNWhsWsCQRsIBK4o+ZLBpfAAAAn+rz4CfKDzjan5kXNj7jpoEuVQR4nz9ClmMVz0wWAADwmT2Hc7T7UE7YWMv6NTR/wlkuVQT4gz+WC3IzYgAA4DMZWXlKe+CTsLH6NRP15Z+HuFQR4B++mMlyrETGAgAAfpCTH9C901fplW+2hMam33C6Tm5V38WqAH/xxUwW3QUBAIBf/PWjtWEB69mr+hCwgCrm+ZDlOFYSNyMGAADeZ63V019uDG0/dXlvDenS1MWKAH/y/HLBgC0IWcxkAQAAr7p3+io991V62Fj61DHuFAPAByGLmSwAAOBRjmPV/o6ZJcafvaqPC9UAKOT95YKFM1mELAAA4CHWRg5Yd405kSWCgMt8M5PFckEAAOAVby3ZqlveWBE2xvJAIHZ4PmQ5TsHvLBcEAADVnbVW7SaGz1798ezOuvnsTi5VBCASzy8XPNr4wuVCAAAAKujXz3wTtn3OyS0IWEAM8v5MlqXxBQAAqP5SJ8wIPU5OiNPMm89Qhya1XawIQGm8H7IKuwtyTRYAAKiGcvID6nLXR2Fjax8Y5VI1AMrC8yErQHdBAABQDW3ae0R7Dufqgie/ChtfMWm4SxUBKCvvhyy6CwIAgGpky75MnfHwZxH3bZwyWobvaYCY5/mQRXdBAABQXWTm5pcasGjRDlQfng9ZR5cLulwIAADAMfz5zW/12uItoe0B7RupZlK8Fm7cp7m3DXavMADl5v2QReMLAAAQw3YfylGfBz8JGzv7xGZ6+so0lyoCUFGeD1kOjS8AAEAMG/HPL8K2X7mmvwZ0aORSNQCiwfMhi8YXAAAgVt3z3krtO5Ib2l5+zzDVr5nkYkUAosHzIYubEQMAgFizdPN+nf+v8NbsNLYAvMPz7SBC3QWZyQIAADHgfwvSSwSs2X8c5E4xACqF52ey6C4IAABiwfYDWTpt6pywsb7tGur5q/qqRlK8S1UBqAyejx50FwQAAG5ZvuWAhj4yVyu3ZZQIWNOu6K3XfzeAgAV4kOdnsuguCAAA3BBwrM57Yr4k6Zz/92XYvov7tNbw7s3dKAtAFfB8yKK7IAAAqErWWt317kodzM4vse/GszrqluFdXKgKQFXyfMhyHLoLAgCAqpGTH1CXuz4KGzu9YyMlJ8Trnxf3VN2URJcqA1CVPB+yAiwXBAAAVaR4wJKkF3/bT4YVNYCveD9k0fgCAABUImut7njnO73yzZaw8Td/P0BpqQ1dqgqAmzwfsmh8AQAAKtN9768OC1hdm9fRR3/gvleAn3m+hfvRmxG7WwcAAPCm575KD9smYAHw/ExW4TVZLBcEAADRlB9w9OaSraHt9KljXKwGQCzxfMgq7C7IckEAABBNHe/8MPS4f3uuvQJwlOdDFt0FAQBANOTkB/TBih265Y0VJfY9O66vCxUBiFXeD1l0FwQAAFEQqT378G7N9PilvZSU4PnL3AGUg+dDFt0FAQBARY1+dF7E8acu7604vscAUIznQ1Yg2F0wnpksAABQTlM//F5Pfb6+xHj61DGy1nKTYQAReT5kFTa+iGMWHwAAlMOTc9eXCFhf/nmImtRJliQCFoBSeT5k0fgCAACUx0sLN+nOd1aWGKdFO4Cy8vz8jsN9sgDAV4wxI40xa40x64wxEyLsH2SMWWqMyTfGXFhs35XGmB+Dv66suqoRK15euLlEwLrmjHYELADl4vmZLIfuggDgG8aYeElPSBomaaukRcaY6dba1UUO2yxpnKRbiz23oaRJktIkWUlLgs/dXxW1w10Bx+rLdXt0xzvfhY3PvOkMndiijktVAaiuPB+yAtyMGAD8pK+kddbaDZJkjHlV0lhJoZBlrU0P7nOKPXeEpNnW2n3B/bMljZT0SuWXDTelTphRYuz+83roiv5tXagGgBd4P2QVZCy6CwKAP7SUtKXI9lZJ/Srw3JZRqgsxKDsvoK53l7z31Q8PjOK+VwAqJCohyxgzUtKjkuIlPW2tnVps/zhJf5W0LTj0uLX26Wi89/HQXRAAEE3GmPGSxktSmzZtXK4G5fX+iu268ZVlEfc1rZOsV8f3J2ABqLAKh6wyrn+XpNestTdU9P3Ki+6CAOAr2yS1LrLdSkd/wFeW5w4u9ty5xQ+y1k6TNE2S0tLS7M8pEu7YfSgnYsBKbVRTc28b4kJFALwqGj+qCa1/t9bmSipc/x4TAjS+AAA/WSSpkzGmnTEmSdLFkqaX8bmzJA03xjQwxjSQNDw4Bg8IOFZ9HvykxPgV/dsSsABEXTSWC5Z1/fsFxphBkn6Q9Edr7ZbiB1TGEgyHxhcA4BvW2nxjzA0qCEfxkp6x1q4yxkyWtNhaO90Y00fSO5IaSPqFMeY+a213a+0+Y8z9KghqkjS5sAkGqjdrrf75yQ+h7aV3D1PDWkkuVgTA66qq8cX7kl6x1uYYY34n6XlJZxU/qDKWYISWCzKTBQC+YK2dKWlmsbF7ijxepIKlgJGe+4ykZyq1QFSpp+dt0AMz1oS2xw9qT8ACUOmisVzwuOvfrbV7rbU5wc2nJfWOwvuWSXAiS3HMZAEA4Bv7juTqw+92hAUsSfrzyK4uVQTAT6IxkxVa/66CcHWxpEuLHmCMaWGt3RHcPFdS+P94lchxrMhXAAD4x6HsPPW6f3aJ8Zb1a3D5AIAqUeGQVZb175JuMsacKylf0j5J4yr6vmUVsJb/UAEA8LiMrDzVSIxXdn5AJ9/7cdi+t687TbWSEtSleR2XqgPgN1G5JqsM698nSpoYjfcqr4KZLEIWAABeZa3VKfd9XGL8V2mtdPvIrmpcO9mFqgD4WVU1vnBNwGEmCwAAL+ty90clxj67dbDaNa7lQjUAEJ3GFzEtYC2dBQEA8KhdB7OVm++EjZ1/aksCFgBXeX4my3EsnQUBAPCgnPyA+j70aWg7feoYF6sBgKP8MZNFyAIAwFNWbstQl7uOLhP87NbB7hUDAMV4fiYr4IjGFwAAVHPf7zyoO99ZqSWb9pfY99a1p7E8EEBM8XzIchyreM/P1wEA4G0j/zkv4njvtg3Uu22DKq4GAI7N+yHL0sIdAIDq6o3FW3Tbm99G3Pfu9aerW4u6VVwRAByf50NWgJAFAEC1dNMryzR9xfawse/vH6mUxHiXKgKAsvF8yHK4TxYAANXKrFU7teNAFgELQLXl+ZAVsCJkAQBQTVwy7Wst2LA3bOyN3w/Qya3qKTmBgAWgevB8yHIcKzIWAACxzVqrdhNnlhife+tgpdI5EEA14/m+ewGWCwIAEPNOmzonbLtx7WR9/MdBBCwA1ZLnZ7JofAEAQGw7nJOvHRnZoW2uvQJQ3Xk+ZNH4AgCA2HT/B6v13y83ho2lTx3jUjUAED3eXy5oCVkAAMSanw5mlwhYvx7Q1qVqACC6PB+yHCsZlgsCABAzFm7Yq34PfVpi/L5zu7tQDQBEnz+WC5KxAABwVXZeQL9/cYnmrt0dNr568ght3pepLs3q8ENRAJ7h+ZBFd0EAANx1xX8Xat6Pe0qMP3FpL9VMSlDX5nVdqAoAKo/3QxbdBQEAcEVuvqPcgBMxYL3+uwHq266hC1UBQOXzfMhyHKukBM9fegYAQMw59/Ev9f3OQ6HtubcOVuuGNVlhAsDzPB+y6C4IAEDVeurz9Zr64fdhYzNuGsiNhQH4hudDluOwXBAAgKqwYP1eXfKfryPu635CvSquBgDc4/mQxUwWAACVb/ehnBIBa8xJLTSoc2ONPqmFS1UBgDu8H7IcMZMFAEAl+mjlTv3+xSVhYwvvGKpmdVNcqggA3OX5kGWtVTx9LwAAqDRFA9YNQzrq1hFdXKwGANzn+fgR4JosAAAqzfrdh8O2CVgA4IOZrIC1iuOaLAAAKsXQRz6XJP3j/07RL09t5XI1ABAbPB+yHMcqnpksAACi5rn5G7VsywFd0b9taIyABQBHeT5k0V0QAIDo+MfsH/TjrkOa+d1OSdJ7y7dLki7r18bNsgAg5ng+ZDl0FwQAoELeXLJVt76xotT9d5/TrQqrAYDY5/mQFXDoLggAwM+1/UBWxIBVOzlBV57WVrcO7yLDDzMBIIz3QxbLBQEA+FnS9xzR4L/NLTE+99bBSm1cq+oLAoBqwvMhy6GFOwAA5XbRU19pUfr+0HackdY9OJqOvQBQBp4PWcxkAQBQdv9bkK6731sVNtamYU19cfsQdwoCgGrI81crMZMFAEDZfLV+T4mANbbnCZpzy5kuVQQA1ZPnZ7IcS3dBAACOZ8/hHF36n4VhY1/+eYhaNajpUkUAUH15PmTRXRAAgGMrvkQwfeoY94oBAA/wfPwIWMtFugAAFLN6+0F9vWGvrLVhAWvJXWe7WBUAeIPnZ7Icxyqe5YIAAIRk5uZr9GPzSoz/Kq2VGtVOdqEiAPAWz4csugsCAHDU9gNZ+tW/F5QYf+LSXhpzcgsXKgIA7/F0yLLWytL4AgCAkNOmzok4TsACgOjxdMgKOFaSmMkCAEDS64u3hG2nTx2j1dsPqkPTWi5VBADe5O2QZQlZAAAczM7Tyfd+HDa2ccpoSVK3E+q6URIAeJqnQ1YwY7FcEADgW2c9Mlcbdh8JG6NFOwBULk+HrMLlgkxkAQD85kBmrnpOnl1i/Lt7h7tQDQD4S1Tuk2WMGWmMWWuMWWeMmRBhf7Ix5rXg/oXGmNRovO/xsFwQAOBXxQPWa+P7K33qGNVJSXSpIgDwjwqHLGNMvKQnJI2S1E3SJcaYbsUO+62k/dbajpL+IekvFX3fsnBCM1mELACAf7z6zebQ4ysHtFX61DHq176RixUBgL9EY7lgX0nrrLUbJMkY86qksZJWFzlmrKR7g4/flPS4McZYW3jVVOWguyAAwA+stTqSG1Dt5ASlTpgRtu++sT1cqgoA/CsaywVbSiraE3ZrcCziMdbafEkZkkr8SM0YM94Ys9gYs3j37t0VLqxwuWAcIQsAfOPnLmE3xqQaY7KMMcuDv56q6tp/jtx8R099vkE9Js0qEbAKOwgCAKpWTDW+sNZOkzRNktLS0io8y+U4Bb/Hs1wQAHyhyBL2YSr4od8iY8x0a23R1RWhJezGmItVsIT9/4L71ltre1Zp0RWw7UCWTi/l5sIL7xgqw+cfALgiGjNZ2yS1LrLdKjgW8RhjTIKkepL2RuG9j+lo44vKficAQIwILWG31uZKKlzCXtRYSc8HH78paaippmkkUsD6v7TWSp86Rs3qprhQEQBAis5M1iJJnYwx7VQQpi6WdGmxY6ZLulLSAkkXSppT2ddjSTS+AAAfirSEvV9px1hr840xRZewtzPGLJN0UNJd1tp5xd/AGDNe0nhJatOmTXSrL4fnv0oP276odytdM6i9Ojer405BAICQCoes4AfUDZJmSYqX9Iy1dpUxZrKkxdba6ZL+K+l/xph1kvapIIhVOocW7gCAstshqY21dq8xprekd40x3a21B4seFO2l7T/HPe+t1AsLNoW2ubkwAMSWqFyTZa2dKWlmsbF7ijzOlnRRNN6rPALMZAGA35RnCfvWokvYgyssciTJWrvEGLNeUmdJiyu96nIo3tziyz8PcakSAEBpPH21kkN3QQDwm9ASdmNMkgpWTkwvdkzhEnapyBJ2Y0yTYOMMGWPaS+okaUMV1V0m2w9khW0vu3uYWjWo6VI1AIDSxFR3wWgL0F0QAHylgkvYB0mabIzJk+RI+r21dl/VfxWlO61Io4uNU0bTPRAAYpTHQxbdBQHAb37uEnZr7VuS3qr0An+mossEv79/JAELAGKYp+NHaLkgH0QAgGrqcE5+WMB65KJTlJIY72JFAIDj8XTIOjqTRcgCAFQ/2w5kqcekWaHtKwe01QW9W7lYEQCgLLwdsmh8AQCopjJz88NuNhxnpPvG9nCxIgBAWXk6ZBXejJjGFwCA6mTLvkx1u2dW2NiGKdwLCwCqC580viBkAQCqh/Q9RzT4b3ND2zcN7aQ/DevsXkEAgHLzdMgKZiwxkQUAqA4uePIrLdm0P2yMgAUA1Y+3lwtalgsCAKqHgGPDAtb1QzoofSpLBAGgOvL0TBbLBQEA1cWh7LzQ45k3naFuJ9R1sRoAQEV4eiaL7oIAgOoiI6sgZD1y0SkELACo5jwdsuguCACoLgpDVr0aiS5XAgCoKE+HLJYLAgCqiwOZwZBVk5AFANWdp0NWYeOLOGayAAAxjpksAPAOT4esgFPwOzNZAIBYR8gCAO/wdsgqbOHu6a8SAOAFhCwA8A5Pxw/LckEAQDVxMCtPyQlxSkmMd7sUAEAFeTpkFTa+IGQBAGJdRlYes1gA4BG+CFlckwUAiHWELADwDk+HLIebEQMAqokDmYQsAPAKT4esUHdBlgsCAGIcM1kA4B3eDlmhmSyXCwEA4DgIWQDgHZ6OH07hNVnMZAEAYtzBrDzVq0nIAgAv8HTIovEFAKA6CDhWh3LymckCAI/wdMii8QUAoDo4yI2IAcBTfBGyWC4IAIhlBwhZAOApng5Zhd0FuRkxACCWZRCyAMBTPB2yHLoLAgCqAUIWAHiLp+NHgO6CAIBqgJAFAN7ij5BF4wsAQAwLhSxauAOAJ3g6ZDnWyhjJMJMFAIhhdBcEAG/xdMgKOJalggCAmJeRlaeUxDglJ8S7XQoAIAq8HbKs5R5ZAICYl5GZxywWAHiIp0OWw0wWAKAaOJCVS8gCAA/xdMgKODS9AADEvowsZrIAwEs8HbIKG18AABDLMrLyCVkA4CGeD1nMZAEAYt3BrDzVJWQBgGd4OmTRXRAAUB2wXBAAvMXTIcuhuyAAIMblBxwdzslX/RpJbpcCAIgST4csZrIAALHuYHa+JKlejQSXKwEARIvHQxbdBQEAsS0jK0+SVK8mywUBwCs8HbIKlgu6XQUAAKULhSyuyQIAz/B0BGG5IAAg1h3IzJVEyAIAL/F2yKLxBQAgxjGTBQDeU6GQZYxpaIyZbYz5Mfh7g1KOCxhjlgd/Ta/Ie5aHtcxkAQBi28FgyOI+WQDgHRWdyZog6VNrbSdJnwa3I8my1vYM/jq3gu9ZZgHHKo6QBQCIYcxkAYD3VDRkjZX0fPDx85LOq+DrRVXAEcsFAQAxLSMrTymJcUpOiHe7FABAlFQ0ZDWz1u4IPt4pqVkpx6UYYxYbY742xlRZEHOsVbynrzoDAFR3GVl53IgYADzmuHc+NMZ8Iql5hF13Ft2w1lpjjC3lZdpaa7cZY9pLmmOM+c5auz7Ce42XNF6S2rRpc9zij4fuggCAWJeRlcdSQQDwmOOGLGvt2aXtM8b8ZIxpYa3dYYxpIWlXKa+xLfj7BmPMXEmnSioRsqy10yRNk6S0tLTSAluZOXQXBADEOEIWAHhPRRfTTZd0ZfDxlZLeK36AMaaBMSY5+LixpNMlra7g+5YJM1kAgFiXkZVPZ0EA8JiKhqypkoYZY36UdHZwW8aYNGPM08FjTpS02BizQtJnkqZaa6ssZDGTBQCIZRmZucxkAYDHHHe54LFYa/dKGhphfLGkq4OPv5J0UkXe5+dyrFVCHJ0vAACxi+WCAOA9nk4gjpXimckCAMSovICjI7kBQhYAeIynQ1bAseKSLADwF2PMSGPMWmPMOmPMhAj7k40xrwX3LzTGpBbZNzE4vtYYM6Kyaz0YuhFxhRaWAABijKdDVsF9skhZAOAXxph4SU9IGiWpm6RLjDHdih32W0n7rbUdJf1D0l+Cz+0m6WJJ3SWNlPSv4OtVmozCkFWTmSwA8BJPhyy6CwKA7/SVtM5au8FamyvpVUljix0zVtLzwcdvShpqjDHB8VettTnW2o2S1gVfr9IUhixuRgwA3uLp9QnPjusjQ8gCAD9pKWlLke2tkvqVdoy1Nt8YkyGpUXD862LPbVn8DYwx4yWNl6Q2bdpUqNgeLevpyz8PUcNahCwA8BJPz2Q1rZuiJnWS3S4DAOAh1tpp1to0a21akyZNKvRaifFxatWgpmomefpnngDgO54OWQAA39kmqXWR7VbBsYjHGGMSJNWTtLeMzwUA4LgIWQAAL1kkqZMxpp0xJkkFjSymFztmuqQrg48vlDTHWmuD4xcHuw+2k9RJ0jdVVDcAwENYnwAA8IzgNVY3SJolKV7SM9baVcaYyZIWW2unS/qvpP8ZY9ZJ2qeCIKbgca9LWi0pX9L11tqAK18IAKBaI2QBADzFWjtT0sxiY/cUeZwt6aJSnvugpAcrtUAAgOexXBAAAAAAooiQBQAAAABRRMgCAAAAgCgiZAEAAABAFBGyAAAAACCKCFkAAAAAEEWELAAAAACIIkIWAAAAAESRsda6XUNExpjdkjZF4aUaS9oThdfxGs5LZJyXyDgvkXFeIovGeWlrrW0SjWIqE59VlY7zEhnnJTLOS2Scl8gq9bMqZkNWtBhjFltr09yuI9ZwXiLjvETGeYmM8xIZ56X8OGeRcV4i47xExnmJjPMSWWWfF5YLAgAAAEAUEbIAAAAAIIr8ELKmuV1AjOK8RMZ5iYzzEhnnJTLOS/lxziLjvETGeYmM8xIZ5yWySj0vnr8mCwAAAACqkh9msgAAAACgyhCyAAAAACCKPBuyjDEjjTFrjTHrjDET3K6nKhhjnjHG7DLGrCwy1tAYM9sY82Pw9wbBcWOMeSx4fr41xvQq8pwrg8f/aIy50o2vJVqMMa2NMZ8ZY1YbY1YZY24Ojvv9vKQYY74xxqwInpf7guPtjDELg1//a8aYpOB4cnB7XXB/apHXmhgcX2uMGeHOVxRdxph4Y8wyY8wHwW3fnxdjTLox5jtjzHJjzOLgmK//HUWD3z6r+JyKjM+qyPisOjY+hKXEHwAABIZJREFUq0qKqc8qa63nfkmKl7ReUntJSZJWSOrmdl1V8HUPktRL0soiYw9LmhB8PEHSX4KPR0v6UJKR1F/SwuB4Q0kbgr83CD5u4PbXVoFz0kJSr+DjOpJ+kNSN8yIjqXbwcaKkhcGv93VJFwfHn5J0bfDxdZKeCj6+WNJrwcfdgv++kiW1C/67i3f764vC+fmTpJclfRDc9v15kZQuqXGxMV//O4rCOfXdZxWfU6WeFz6rIp8XPquOfX74rCp5TmLms8qrM1l9Ja2z1m6w1uZKelXSWJdrqnTW2i8k7Ss2PFbS88HHz0s6r8j4C7bA15LqG2NaSBohaba1dp+1dr+k2ZJGVn71lcNau8NauzT4+JCkNZJaivNirbWHg5uJwV9W0lmS3gyOFz8vhefrTUlDjTEmOP6qtTbHWrtR0joV/PurtowxrSSNkfT0/2/n7l3lqMI4jn8f8BUVNUFTGEEDAStRCBIxRRAMKJIqhSAoKlhbCRLwTxAtLC2DhRgxnW9JrwSjRuLLFQS9RC8Iia0vj8V5NpmYXUO4Y2az5/uB5e6c2bvseZiZ387M2VPLgXVZpOv9aATdZZU5NZ9ZNZ9ZtZhZdVkm2Y9W9STrLuCnwfLP1dajbZl5up7/Amyr54tqtLK1q9vjD9KuhHVflxpmcALYoB1AfgDOZOaf9ZJhH8/1v9afBbaygnUBXgdeBv6u5a1YF2hfbD6MiOMR8WK1db8fbZL1aNyOBsyqC5lVC5lV8y1NVl1zuf+gq1dmZkR0OWd/RNwMvAu8lJm/tws4Ta91ycy/gAci4jbgPeC+iT/S5CLiSWAjM49HxN6pP8+S2ZOZ6xFxJ/BRRHwzXNnrfqRx9b4dmVUXM6suZlb9p6XJqlW9k7UO3D1Y3l5tPfq1bn1SfzeqfVGNVq52EXEtLbQOZebhau6+LjOZeQY4BjxMu1U+u/gy7OO5/tf6W4HfWL26PALsj4gfaUO3HgXewLqQmev1d4P2Rech3I82y3o0bkeYVZdiVl3ArFpgmbJqVU+yPgN21iwr19F+5Hdk4s80lSPAbFaUZ4H3B+3P1Mwqu4GzdSv1A2BfRNxes6/sq7arUo05fgs4lZmvDVb1Xpc76qogEXEj8BjtNwDHgAP1sn/XZVavA8DRzMxqf6pmLroX2Al8emV6Mb7MfCUzt2fmPbTjxtHMfJrO6xIRN0XELbPntO3/JJ3vRyMwq5rutyOzaj6zaj6zar6ly6pcgplA/o8HbcaQ72hjdw9O/XmuUJ/fBk4Df9DGj75AG3P7CfA98DGwpV4bwJtVn6+AXYP3eZ7248c14Lmp+7XJmuyhjc/9EjhRjyesC/cDn1ddTgKvVvsO2gF2DXgHuL7ab6jltVq/Y/BeB6te3wKPT923EWu0l/MzNnVdl+r/F/X4enZM7X0/Gqm2XWWVObWwLmbV/LqYVZeukVl1vi9LlVVRbyRJkiRJGsGqDheUJEmSpEl4kiVJkiRJI/IkS5IkSZJG5EmWJEmSJI3IkyxJkiRJGpEnWZIkSZI0Ik+yJEmSJGlE/wBKz8NrlY5OhwAAAABJRU5ErkJggg==\n",
            "text/plain": [
              "<Figure size 864x360 with 2 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "bcSvUnxhfqKm",
        "colab_type": "text"
      },
      "source": [
        "### 2. High-Level Implementation"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "v62AgNvnfwLV",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 187
        },
        "outputId": "c437311c-0549-47a9-e1d4-6f7b1190817e"
      },
      "source": [
        "tf.keras.backend.clear_session()\n",
        "\n",
        "model = tf.keras.models.Sequential()\n",
        "model.add(tf.keras.layers.Dense(1,input_shape =(1,)))        # total param # = 2\n",
        "model.summary()"
      ],
      "execution_count": 123,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Model: \"sequential\"\n",
            "_________________________________________________________________\n",
            "Layer (type)                 Output Shape              Param #   \n",
            "=================================================================\n",
            "dense (Dense)                (None, 1)                 2         \n",
            "=================================================================\n",
            "Total params: 2\n",
            "Trainable params: 2\n",
            "Non-trainable params: 0\n",
            "_________________________________________________________________\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "iMMxiROvgYKW",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "loss_fn = tf.keras.losses.mean_squared_error\n",
        "optimizer = tf.keras.optimizers.Adam(learning_rate=0.001)\n",
        "model.compile(optimizer=optimizer, loss=loss_fn, metrics=[tf.keras.metrics.RootMeanSquaredError()])"
      ],
      "execution_count": 125,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "nvUCh47ti2bb",
        "colab_type": "text"
      },
      "source": [
        "We can either mix it with low-level api"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "ufaKd-4pire9",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 187
        },
        "outputId": "78450502-9945-4894-800c-2ce7a518a509"
      },
      "source": [
        "epochs = 200\n",
        "\n",
        "for epoch in tf.range(1,epochs+1):      # for each epoch\n",
        "  loss = tf.constant(0.0)\n",
        "  for x, y in ds:           # for each batch\n",
        "    with tf.GradientTape() as tape:\n",
        "      y_bar = model(y)                # __call__\n",
        "      loss = loss_fn(tf.reshape(y, [-1]), tf.reshape(y_bar, [-1]))\n",
        "    # Back propagation to calculate gradient\n",
        "    grads = tape.gradient(loss, model.variables)\n",
        "    # apply_gradients\n",
        "    model.optimizer.apply_gradients(zip(grads, model.variables))\n",
        "  \n",
        "  if epoch%20==0:\n",
        "    tf.print('w=', w, ' ,b=', b, ' ,mse=', loss)"
      ],
      "execution_count": 126,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "w= [1.9998285047327671]  ,b= 0.23774513013090734  ,mse= 32856.3203\n",
            "w= [1.9998285047327671]  ,b= 0.23774513013090734  ,mse= 6936.47\n",
            "w= [1.9998285047327671]  ,b= 0.23774513013090734  ,mse= 790.176208\n",
            "w= [1.9998285047327671]  ,b= 0.23774513013090734  ,mse= 7.72515583\n",
            "w= [1.9998285047327671]  ,b= 0.23774513013090734  ,mse= 0.698684156\n",
            "w= [1.9998285047327671]  ,b= 0.23774513013090734  ,mse= 1.83272243\n",
            "w= [1.9998285047327671]  ,b= 0.23774513013090734  ,mse= 1.08908737\n",
            "w= [1.9998285047327671]  ,b= 0.23774513013090734  ,mse= 0.85788995\n",
            "w= [1.9998285047327671]  ,b= 0.23774513013090734  ,mse= 0.855531335\n",
            "w= [1.9998285047327671]  ,b= 0.23774513013090734  ,mse= 0.84222728\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "pIo0f3WMjA9t",
        "colab_type": "text"
      },
      "source": [
        "Or just call model.fit"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "iWVw5VxcjDJi",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        },
        "outputId": "fed7b0e2-42ca-4862-bf61-122d98875081"
      },
      "source": [
        "tf.keras.backend.clear_session()\n",
        "model.reset_states()\n",
        "\n",
        "model.fit(tf.convert_to_tensor(xs), tf.convert_to_tensor(ys), batch_size = 10, epochs = 200, verbose=0) \n",
        "tf.print('w=', model.layers[0].kernel, ' b=', model.layers[0].bias)"
      ],
      "execution_count": 128,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "w= [[1.99535918]]  b= [0.967212737]\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "APG0pQJOl5cS",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        ""
      ],
      "execution_count": null,
      "outputs": []
    }
  ]
}