linear/LinearRegressionUnivariate-TensorFlow.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Linear regression using Tensor Flow 2.0\n",
    "\n",
    "Reuse the same univariate polynomial model as in the [univariate linear regression workbook](LinearRegressionUnivariate.html)  ([Jupyter Notebook](LinearRegressionUnivariate.ipynb)) but with a TensorFlow implementation\n",
    "\n",
    "From TensorFlow 2.0, recommended API is clearly the one of Keras. Building and fitting a model is performed in most cases with few lines of code. This is hidding a lot of the details on how the fit is done.\n",
    "\n",
    "This tutorial is a mid step between the home made gradient descent of above mentionned tutorial and the fully wrapped Keras model as in [Bivariate linear regression with Keras](LinearRegressionBivariate-Keras.html) ([Notebook](LinearRegressionBivariate-Keras.ipynb))\n",
    "\n",
    "The model is mostly hand designed, even if doing the same as the Keras Dense layer. But the optimizers are the powerful ones of Keras and TensorFlow.\n",
    "\n",
    "Learning goals:\n",
    "- Design a model in TensorFlow 2.0\n",
    "- Use TensorFlow to perform gradient descents. \n",
    "- Compare several optimizers\n",
    "\n",
    "References:\n",
    "- https://databricks.com/tensorflow/training-and-convergence"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import tensorflow as tf # TF 2.0 required\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import pandas\n",
    "from sklearn import metrics "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Data model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Univariate\n",
    "numFeatures = 1\n",
    "\n",
    "def generateBatch(N, stochastic = False):\n",
    "    #\n",
    "    xMin = 0\n",
    "    xMax = 0.5\n",
    "    #\n",
    "    b = 0.35\n",
    "    std = 0.01\n",
    "    #\n",
    "    if stochastic:\n",
    "        x = np.random.uniform(xMin, xMax, N)\n",
    "    else:\n",
    "        x = np.linspace(xMin, xMax, N)\n",
    "    yClean = x**4 + (x-0.3)**3 + b\n",
    "    y =  yClean + np.random.normal(0, std, N) \n",
    "    return (x, y, yClean)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Reference values for linear regression\n",
    "\n",
    "From [LinearRegressionUnivariate.ipynb](LinearRegressionUnivariate.ipynb)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "wRef, bRef = 0.145, 0.323"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Test data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Ntest = 100000\n",
    "(xTest, yTest1, yTestClean1) = generateBatch(Ntest)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Shared parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# (Mini) Batch size\n",
    "nBatch = 128\n",
    "# Number of batches per Epoch\n",
    "nBatchPerEpoch = 10\n",
    "# Stop threshold on the Epoch MSE\n",
    "threshold = 1e-4\n",
    "# Safe guard to stop on number of epochs\n",
    "nEpochMax = 2000"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Helpers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot helper to show target reference value\n",
    "def plotToRef(label, subNum, subRow, subCol, ref, values):\n",
    "    nIter = len(values)\n",
    "    r = range(nIter)\n",
    "    plt.subplot(subNum, subRow, subCol, alpha=1)\n",
    "    plt.title(label)\n",
    "    plt.plot(r, values, r, np.ones((nIter))*ref, alpha=0.5)\n",
    "    plt.grid();"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Training model\n",
    "\n",
    "From TensorFlow 2.0, the symbolic variables are no longer explicit, they look alike any Python variables.\n",
    "\n",
    "The model is generated at \"compile\" time, rewritting and checking the code to act as a symbolic model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Simple custom layer exposing the linear regression model\n",
    "class MyDenseLayer(tf.keras.layers.Layer):\n",
    "    def __init__(self, *args, **kwargs):\n",
    "        super(MyDenseLayer, self).__init__(*args, **kwargs)\n",
    "    \n",
    "    def build(self, input_shape):\n",
    "        self.w = self.add_weight(\n",
    "            shape=1.0,\n",
    "            dtype=self.dtype,\n",
    "            initializer=tf.keras.initializers.ones(),\n",
    "            #regularizer=tf.keras.regularizers.l2(0.02),\n",
    "            trainable=True)\n",
    "        self.b = self.add_weight(\n",
    "            shape=1.0,\n",
    "            dtype=self.dtype,\n",
    "            initializer=tf.keras.initializers.ones(),\n",
    "            #regularizer=tf.keras.regularizers.l2(0.02),\n",
    "            trainable=True)\n",
    "\n",
    "    @tf.function\n",
    "    def call(self, x, training=None):\n",
    "        return x * self.w + self.b"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Gradient descent optimizer"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "# Model 1, instantiate the custom layer\n",
    "model1 = tf.keras.Sequential([MyDenseLayer(input_shape=[numFeatures], dtype=\"float64\")])\n",
    "\n",
    "# Stochastic Gradient Descent Optimizer\n",
    "optim1 = tf.keras.optimizers.SGD(0.01)\n",
    "\n",
    "# Perform a train step on a mini-batch\n",
    "#  This function's code is rewritten by TensorFlow 2.0 and shall be compiled at every execution of the optimizer\n",
    "@tf.function\n",
    "def trainStep1(x, y):\n",
    "    with tf.GradientTape() as tape:\n",
    "        predictions = model1(x, training=True)\n",
    "        loss = tf.keras.losses.mean_squared_error(y, predictions)\n",
    "        \n",
    "        gradients = tape.gradient(loss, model1.trainable_variables)\n",
    "        optim1.apply_gradients(zip(gradients, model1.trainable_variables))\n",
    "        return loss\n",
    "    \n",
    "# Initialize values and loop on epochs and mini batch\n",
    "epoch = 0\n",
    "mse_epoch = 1\n",
    "hist = []\n",
    "while mse_epoch > threshold and epoch < nEpochMax:\n",
    "    mse_cumul = 0\n",
    "    for b in range(0, nBatchPerEpoch):  \n",
    "        (xTrain, yTrain, yTrainClean) = generateBatch(nBatch, True)\n",
    "        mse_cumul += trainStep1(xTrain, yTrain)\n",
    "        \n",
    "    W = model1.get_weights()\n",
    "    mse_epoch = mse_cumul / nBatchPerEpoch\n",
    "    hist.append((W[1][0], W[0][0], mse_epoch))\n",
    "    epoch += 1\n",
    "        \n",
    "print(\"Predicted model: {a:.3f} x + {b:.3f}, num epochs={c}\".format(a=w, b=b, c=len(wLearn)))\n",
    "df1 = pandas.DataFrame(hist, columns = ('b', 'w', 'MSE'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(16,10))\n",
    "plotToRef('b', 2, 2, 1, bRef, df1['b'])\n",
    "plotToRef('w', 2, 2, 2, wRef, df1['w'])\n",
    "plt.subplot(2,2,3)\n",
    "plt.semilogy(df1['MSE'])\n",
    "plt.grid()\n",
    "plt.title(('Loss (MSE)'));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Test model "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "yEst1 = w * xTest + b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.plot(xTest, yTestClean1, xTest, yEst1);\n",
    "plt.legend(('Test data (clean)', 'SGD'), loc='upper left')\n",
    "mse1 = metrics.mean_squared_error(yTest1, yEst1)\n",
    "print('Gradient Optimizer MSE = {:.3e}'.format(mse1));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Moment optimizer"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Model, instantiate the custom layer\n",
    "model2 = tf.keras.Sequential([MyDenseLayer(input_shape=[numFeatures], dtype=\"float64\")])\n",
    "\n",
    "# Gradient Descent Optimizer\n",
    "optim2 = tf.keras.optimizers.SGD(0.01, momentum=0.0001) # <---\n",
    "\n",
    "# Perform a train step on a mini-batch\n",
    "#  This function's code is rewritten by TensorFlow 2.0 and shall be compiled at every execution of the optimizer\n",
    "@tf.function\n",
    "def trainStep2(x, y):\n",
    "    with tf.GradientTape() as tape:\n",
    "        predictions = model2(x, training=True)\n",
    "        loss = tf.keras.losses.mean_squared_error(y, predictions)\n",
    "        \n",
    "        gradients = tape.gradient(loss, model2.trainable_variables)\n",
    "        optim2.apply_gradients(zip(gradients, model2.trainable_variables))\n",
    "        return loss\n",
    "    \n",
    "# Initialize values and loop on epochs and mini batch\n",
    "epoch = 0\n",
    "mse_epoch = 1\n",
    "wLearn = []\n",
    "while mse_epoch > threshold and epoch < nEpochMax:\n",
    "    mse_cumul = 0\n",
    "    for b in range(0, nBatchPerEpoch):  \n",
    "        (xTrain, yTrain, yTrainClean) = generateBatch(nBatch, True)\n",
    "        mse_cumul += trainStep2(xTrain, yTrain) # <---\n",
    "    \n",
    "    W = model1.get_weights()\n",
    "    mse_epoch = mse_cumul / nBatchPerEpoch\n",
    "    hist.append((W[1][0], W[0][0], mse_epoch))\n",
    "    epoch += 1\n",
    "        \n",
    "print(\"Predicted model: {a:.3f} x + {b:.3f}, num iterations={c}\".format(a=w, b=b, c=len(wLearn)))\n",
    "df2 = pandas.DataFrame(wLearn, columns = ('b', 'w', 'MSE'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(16,10))\n",
    "plotToRef('b', 2, 2, 1, bRef, df2['b'])\n",
    "plotToRef('w', 2, 2, 2, wRef, df2['w'])\n",
    "plt.subplot(2,2,3)\n",
    "plt.semilogy(df2['MSE'])\n",
    "plt.grid()\n",
    "plt.title(('Loss (MSE)'));"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "yEst2 = w * xTest + b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "plt.plot(xTest, yTestClean1, xTest, yEst1, xTest, yEst2);\n",
    "plt.legend(('Test data (clean)', 'SGD', 'Momentum'), loc='upper left')\n",
    "mse2 = metrics.mean_squared_error(yTest1, yEst2)\n",
    "print('Moment Optimizer MSE = {:.3e}'.format(mse2));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Adam optimizer\n",
    "\n",
    "https://www.tensorflow.org/api_docs/python/tf/train/AdamOptimizer\n",
    "\n",
    "http://arxiv.org/pdf/1412.6980.pdf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Model, instantiate the custom layer\n",
    "model3 = tf.keras.Sequential([MyDenseLayer(input_shape=[numFeatures], dtype=\"float64\")])\n",
    "\n",
    "# Gradient Descent Optimizer\n",
    "optim3 = tf.keras.optimizers.Adam(0.01) # <---\n",
    "\n",
    "# Perform a train step on a mini-batch\n",
    "#  This function's code is rewritten by TensorFlow 2.0 and shall be compiled at every execution of the optimizer\n",
    "@tf.function\n",
    "def trainStep3(x, y):\n",
    "    with tf.GradientTape() as tape:\n",
    "        predictions = model3(x, training=True)\n",
    "        loss = tf.keras.losses.mean_squared_error(y, predictions)\n",
    "        \n",
    "        gradients = tape.gradient(loss, model3.trainable_variables)\n",
    "        optim3.apply_gradients(zip(gradients, model3.trainable_variables))\n",
    "        return loss\n",
    "    \n",
    "# Initialize values and loop on epochs and mini batch\n",
    "epoch = 0\n",
    "mse_epoch = 1\n",
    "wLearn = []\n",
    "while mse_epoch > threshold and epoch < nEpochMax:\n",
    "    mse_cumul = 0\n",
    "    for b in range(0, nBatchPerEpoch):  \n",
    "        (xTrain, yTrain, yTrainClean) = generateBatch(nBatch, True)\n",
    "        mse_cumul += trainStep3(xTrain, yTrain)\n",
    "        \n",
    "    W = model1.get_weights()\n",
    "    mse_epoch = mse_cumul / nBatchPerEpoch\n",
    "    hist.append((W[1][0], W[0][0], mse_epoch))\n",
    "    epoch += 1\n",
    "        \n",
    "print(\"Predicted model: {a:.3f} x + {b:.3f}, num iterations={c}\".format(a=w, b=b, c=len(wLearn)))\n",
    "df3 = pandas.DataFrame(wLearn, columns = ('b', 'w', 'MSE'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(16,10))\n",
    "plotToRef('$b$', 2, 2, 1, bRef, df3['b'])\n",
    "plotToRef('$w$', 2, 2, 2, wRef, df3['w'])\n",
    "plt.subplot(2,2,3)\n",
    "plt.semilogy(df3['MSE'])\n",
    "plt.grid()\n",
    "plt.title(('Loss (MSE)'));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There is a clear gain compared to the standard and momentum gradient descent : \n",
    "- Less iterations\n",
    "- Less undershoot on $b$\n",
    "- Clear convergence of the MSE to the noise floor\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "yEst3 = w * xTest + b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "plt.plot(xTest, yTestClean1, xTest, yEst1, xTest, yEst3);\n",
    "plt.legend(('Test data (clean)', 'SGD', 'Adam'), loc='upper left')\n",
    "mse3 = metrics.mean_squared_error(yTest1, yEst3)\n",
    "print('Adam Optimizer MSE = {:.3e}'.format(mse3));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Where to go from here ?\n",
    "\n",
    "__Other single feature linear implementation__ [using closed form, Scipy, SKLearn or \"raw\" Python](LinearRegressionUnivariate.html) ([Notebook](LinearRegressionUnivariate.ipynb))\n",
    "\n",
    "__More complex bivariation models__ [using \"raw\" Python](LinearRegressionBivariate.html) ([Notebook](LinearRegressionBivariate.ipynb)) up to the gradient descent with regularizer, or [using Keras](LinearRegressionBivariate-Keras.html) ([Notebook](LinearRegressionBivariate-Keras.ipynb))\n",
    "\n",
    "__Compare with the single feature binary classification using logistic regression__ [using \"raw\" Python or libraries](../classification/ClassificationContinuousSingleFeature.html) ([Notebook](../classification/ClassificationContinuousSingleFeature.ipynb]))"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}