Notebooks/Chap05/5_3_Multiclass_Cross_entropy_Loss.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/udlbook/udlbook/blob/main/Notebooks/Chap05/5_3_Multiclass_Cross_entropy_Loss.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "jSlFkICHwHQF"
      },
      "source": [
        "# **Notebook 5.3 Multiclass Cross-Entropy Loss**\n",
        "\n",
        "This notebook investigates the multi-class cross-entropy loss.  It follows from applying the formula in section 5.2 to a loss function based on the Categorical distribution.\n",
        "\n",
        "Work through the cells below, running each cell in turn. In various places you will see the words \"TO DO\". Follow the instructions at these places and make predictions about what is going to happen or write code to complete the functions.\n",
        "\n",
        "Contact me at udlbookmail@gmail.com if you find any mistakes or have any suggestions."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "PYMZ1x-Pv1ht"
      },
      "outputs": [],
      "source": [
        "# Imports math library\n",
        "import numpy as np\n",
        "# Used for repmat\n",
        "import numpy.matlib\n",
        "# Imports plotting library\n",
        "import matplotlib.pyplot as plt\n",
        "# Import math Library\n",
        "import math"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "Fv7SZR3tv7mV"
      },
      "outputs": [],
      "source": [
        "# Define the Rectified Linear Unit (ReLU) function\n",
        "def ReLU(preactivation):\n",
        "  activation = preactivation.clip(0.0)\n",
        "  return activation\n",
        "\n",
        "# Define a shallow neural network\n",
        "def shallow_nn(x, beta_0, omega_0, beta_1, omega_1):\n",
        "    # Make sure that input data is (1 x n_data) array\n",
        "    n_data = x.size\n",
        "    x = np.reshape(x,(1,n_data))\n",
        "\n",
        "    # This runs the network for ALL of the inputs, x at once so we can draw graph\n",
        "    h1 = ReLU(np.matmul(beta_0,np.ones((1,n_data))) + np.matmul(omega_0,x))\n",
        "    model_out = np.matmul(beta_1,np.ones((1,n_data))) + np.matmul(omega_1,h1)\n",
        "    return model_out"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "pUT9Ain_HRim"
      },
      "outputs": [],
      "source": [
        "# Get parameters for model -- we can call this function to easily reset them\n",
        "def get_parameters():\n",
        "  # And we'll create a network that approximately fits it\n",
        "  beta_0 = np.zeros((3,1));  # formerly theta_x0\n",
        "  omega_0 = np.zeros((3,1)); # formerly theta_x1\n",
        "  beta_1 = np.zeros((3,1));  # NOTE -- there are three outputs now (one for each class, so three output biases)\n",
        "  omega_1 = np.zeros((3,3)); # NOTE -- there are three outputs now (one for each class, so nine output weights, connecting 3 hidden units to 3 outputs)\n",
        "\n",
        "  beta_0[0,0] = 0.3; beta_0[1,0] = -1.0; beta_0[2,0] = -0.5\n",
        "  omega_0[0,0] = -1.0; omega_0[1,0] = 1.8; omega_0[2,0] = 0.65\n",
        "  beta_1[0,0] = 2.0; beta_1[1,0] = -2; beta_1[2,0] = 0.0\n",
        "  omega_1[0,0] = -24.0; omega_1[0,1] = -8.0; omega_1[0,2] = 50.0\n",
        "  omega_1[1,0] = -2.0; omega_1[1,1] = 8.0; omega_1[1,2] = -30.0\n",
        "  omega_1[2,0] = 16.0; omega_1[2,1] = -8.0; omega_1[2,2] =-8\n",
        "\n",
        "  return beta_0, omega_0, beta_1, omega_1"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "NRR67ri_1TzN"
      },
      "outputs": [],
      "source": [
        "# Utility function for plotting data\n",
        "def plot_multiclass_classification(x_model, out_model, lambda_model, x_data = None, y_data = None, title= None):\n",
        "  # Make sure model data are 1D arrays\n",
        "  n_data = len(x_model)\n",
        "  n_class = 3\n",
        "  x_model = np.squeeze(x_model)\n",
        "  out_model = np.reshape(out_model, (n_class,n_data))\n",
        "  lambda_model = np.reshape(lambda_model, (n_class,n_data))\n",
        "\n",
        "  fig, ax = plt.subplots(1,2)\n",
        "  fig.set_size_inches(7.0, 3.5)\n",
        "  fig.tight_layout(pad=3.0)\n",
        "  ax[0].plot(x_model,out_model[0,:],'r-')\n",
        "  ax[0].plot(x_model,out_model[1,:],'g-')\n",
        "  ax[0].plot(x_model,out_model[2,:],'b-')\n",
        "  ax[0].set_xlabel('Input, $x$'); ax[0].set_ylabel('Model outputs')\n",
        "  ax[0].set_xlim([0,1]);ax[0].set_ylim([-4,4])\n",
        "  if title is not None:\n",
        "    ax[0].set_title(title)\n",
        "  ax[1].plot(x_model,lambda_model[0,:],'r-')\n",
        "  ax[1].plot(x_model,lambda_model[1,:],'g-')\n",
        "  ax[1].plot(x_model,lambda_model[2,:],'b-')\n",
        "  ax[1].set_xlabel('Input, $x$'); ax[1].set_ylabel('$\\lambda$ or Pr(y=k|x)')\n",
        "  ax[1].set_xlim([0,1]);ax[1].set_ylim([-0.1,1.05])\n",
        "  if title is not None:\n",
        "    ax[1].set_title(title)\n",
        "  if x_data is not None:\n",
        "    for i in range(len(x_data)):\n",
        "      if y_data[i] ==0:\n",
        "        ax[1].plot(x_data[i],-0.05, 'r.')\n",
        "      if y_data[i] ==1:\n",
        "        ax[1].plot(x_data[i],-0.05, 'g.')\n",
        "      if y_data[i] ==2:\n",
        "        ax[1].plot(x_data[i],-0.05, 'b.')\n",
        "  plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "PsgLZwsPxauP"
      },
      "source": [
        "# Multiclass classification\n",
        "\n",
        "For multiclass classification, the network must predict the probability of $K$ classes, using $K$ outputs.  However, these probability must be non-negative and sum to one, and the network outputs can take arbitrary values.  Hence, we pass the outputs through a softmax function which maps $K$ arbitrary values to $K$ non-negative values that sum to one."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "uFb8h-9IXnIe"
      },
      "outputs": [],
      "source": [
        "# Softmax function that maps a vector of arbitrary values to a vector of values that are positive and sum to one.\n",
        "def softmax(model_out):\n",
        "  # This operation has to be done separately for every column of the input\n",
        "  # Compute exponentials of all the elements\n",
        "  # TODO: compute the softmax function (eq 5.22)\n",
        "  # Replace this skeleton code\n",
        "\n",
        "  # Compute the exponential of the model outputs\n",
        "  exp_model_out = np.zeros_like(model_out) ;\n",
        "  # Compute the sum of the exponentials (denominator of equation 5.22)\n",
        "  sum_exp_model_out = np.zeros_like(model_out) ;\n",
        "  # Normalize the exponentials (np.matlib.repmat might be useful here)\n",
        "  softmax_model_out = np.ones_like(model_out)/ exp_model_out.shape[0]\n",
        "\n",
        "  return softmax_model_out"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "VWzNOt1swFVd"
      },
      "outputs": [],
      "source": [
        "\n",
        "# Let's create some 1D training data\n",
        "x_train = np.array([0.09291784,0.46809093,0.93089486,0.67612654,0.73441752,0.86847339,\\\n",
        "                   0.49873225,0.51083168,0.18343972,0.99380898,0.27840809,0.38028817,\\\n",
        "                   0.12055708,0.56715537,0.92005746,0.77072270,0.85278176,0.05315950,\\\n",
        "                   0.87168699,0.58858043])\n",
        "y_train = np.array([2,0,1,2,1,0,\\\n",
        "                    0,2,2,0,2,0,\\\n",
        "                    2,0,1,2,1,2, \\\n",
        "                    1,0])\n",
        "\n",
        "# Get parameters for the model\n",
        "beta_0, omega_0, beta_1, omega_1 = get_parameters()\n",
        "\n",
        "# Define a range of input values\n",
        "x_model = np.arange(0,1,0.01)\n",
        "# Run the model to get values to plot and plot it.\n",
        "model_out= shallow_nn(x_model, beta_0, omega_0, beta_1, omega_1)\n",
        "lambda_model = softmax(model_out)\n",
        "plot_multiclass_classification(x_model, model_out, lambda_model, x_train, y_train)\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "MvVX6tl9AEXF"
      },
      "source": [
        "The left is model output and the right is the model output after the softmax has been applied, so it now lies in the range [0,1] and represents the probability, that y=0 (red), 1 (green) and 2 (blue).  The dots at the bottom show the training data with the same color scheme.  So we want the red curve to be high where there are red dots, the green curve to be high where there are green dots, and the blue curve to be high where there are blue dots  We'll compute the the likelihood and the negative log likelihood."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "YaLdRlEX0FkU"
      },
      "outputs": [],
      "source": [
        "# Return probability under categorical distribution for observed class y\n",
        "# Just take value from row k of lambda param where y =k,\n",
        "def categorical_distribution(y, lambda_param):\n",
        "    return np.array([lambda_param[row, i] for i, row in enumerate (y)])"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "4TSL14dqHHbV"
      },
      "outputs": [],
      "source": [
        "# Let's double check we get the right answer before proceeding\n",
        "print(\"Correct answer = %3.3f, Your answer = %3.3f\"%(0.2,categorical_distribution(np.array([[0]]),np.array([[0.2],[0.5],[0.3]]))))\n",
        "print(\"Correct answer = %3.3f, Your answer = %3.3f\"%(0.5,categorical_distribution(np.array([[1]]),np.array([[0.2],[0.5],[0.3]]))))\n",
        "print(\"Correct answer = %3.3f, Your answer = %3.3f\"%(0.3,categorical_distribution(np.array([[2]]),np.array([[0.2],[0.5],[0.3]]))))\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "R5z_0dzQMF35"
      },
      "source": [
        "Now let's compute the likelihood using this function"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "zpS7o6liCx7f"
      },
      "outputs": [],
      "source": [
        "# Return the likelihood of all of the data under the model\n",
        "def compute_likelihood(y_train, lambda_param):\n",
        "  # TODO -- compute the likelihood of the data -- the product of the categorical probabilities for each data point\n",
        "  # Top line of equation 5.3 in the notes\n",
        "  # You will need np.prod() and the categorical_distribution function you used above\n",
        "  # Replace the line below\n",
        "  likelihood = 0\n",
        "\n",
        "  return likelihood"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "1hQxBLoVNlr2"
      },
      "outputs": [],
      "source": [
        "# Let's test this\n",
        "beta_0, omega_0, beta_1, omega_1 = get_parameters()\n",
        "# Use our neural network to predict the parameters of the categorical distribution\n",
        "model_out = shallow_nn(x_train, beta_0, omega_0, beta_1, omega_1)\n",
        "lambda_train = softmax(model_out)\n",
        "# Compute the likelihood\n",
        "likelihood = compute_likelihood(y_train, lambda_train)\n",
        "# Let's double check we get the right answer before proceeding\n",
        "print(\"Correct answer = %9.9f, Your answer = %9.9f\"%(0.000000041,likelihood))"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "HzphKgPfOvlk"
      },
      "source": [
        "You can see that this gives a very small answer, even for this small 1D dataset, and with the model fitting quite well.  This is because it is the product of several probabilities, which are all quite small themselves.\n",
        "This will get out of hand pretty quickly with real datasets -- the likelihood will get so small that we can't represent it with normal finite-precision math\n",
        "\n",
        "This is why we use negative log likelihood"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "dsT0CWiKBmTV"
      },
      "outputs": [],
      "source": [
        "# Return the negative log likelihood of the data under the model\n",
        "def compute_negative_log_likelihood(y_train, lambda_param):\n",
        "  # TODO -- compute the negative log likelihood of the data -- don't use the likelihood function above -- compute the negative sum of the log probabilities\n",
        "  # You will need np.sum(), np.log()\n",
        "  # Replace the line below\n",
        "  nll = 0\n",
        "\n",
        "  return nll"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "nVxUXg9rQmwI"
      },
      "outputs": [],
      "source": [
        "# Let's test this\n",
        "beta_0, omega_0, beta_1, omega_1 = get_parameters()\n",
        "# Use our neural network to predict the parameters of the categorical distribution\n",
        "model_out = shallow_nn(x_train, beta_0, omega_0, beta_1, omega_1)\n",
        "# Pass the outputs through the softmax function\n",
        "lambda_train = softmax(model_out)\n",
        "# Compute the negative log likelihood\n",
        "nll = compute_negative_log_likelihood(y_train, lambda_train)\n",
        "# Let's double check we get the right answer before proceeding\n",
        "print(\"Correct answer = %9.9f, Your answer = %9.9f\"%(17.015457867,nll))"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "OgcRojvPWh4V"
      },
      "source": [
        "Now let's investigate finding the maximum likelihood / minimum negative log likelihood solution.  For simplicity, we'll assume that all the parameters are fixed except one and look at how the likelihood and negative log likelihood change as we manipulate the last parameter.  We'll start with overall y_offset, $\\beta_1$ (formerly $\\phi_0$)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "pFKtDaAeVU4U"
      },
      "outputs": [],
      "source": [
        "# Define a range of values for the parameter\n",
        "beta_1_vals = np.arange(-2,6.0,0.1)\n",
        "# Create some arrays to store the likelihoods, negative log likelihoods\n",
        "likelihoods = np.zeros_like(beta_1_vals)\n",
        "nlls = np.zeros_like(beta_1_vals)\n",
        "\n",
        "# Initialise the parameters\n",
        "beta_0, omega_0, beta_1, omega_1 = get_parameters()\n",
        "for count in range(len(beta_1_vals)):\n",
        "  # Set the value for the parameter\n",
        "  beta_1[0,0] = beta_1_vals[count]\n",
        "  # Run the network with new parameters\n",
        "  model_out = shallow_nn(x_train, beta_0, omega_0, beta_1, omega_1)\n",
        "  lambda_train = softmax(model_out)\n",
        "  # Compute and store the two values\n",
        "  likelihoods[count] = compute_likelihood(y_train,lambda_train)\n",
        "  nlls[count] = compute_negative_log_likelihood(y_train, lambda_train)\n",
        "  # Draw the model for every 20th parameter setting\n",
        "  if count % 20 == 0:\n",
        "    # Run the model to get values to plot and plot it.\n",
        "    model_out = shallow_nn(x_model, beta_0, omega_0, beta_1, omega_1)\n",
        "    lambda_model = softmax(model_out)\n",
        "    plot_multiclass_classification(x_model, model_out, lambda_model, x_train, y_train, title=\"beta1[0,0]=%3.3f\"%(beta_1[0,0]))\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "UHXeTa9MagO6"
      },
      "outputs": [],
      "source": [
        "# Now let's plot the likelihood and negative log likelihood as a function of the value of the offset beta1\n",
        "fig, ax = plt.subplots()\n",
        "fig.tight_layout(pad=5.0)\n",
        "likelihood_color = 'tab:red'\n",
        "nll_color = 'tab:blue'\n",
        "\n",
        "\n",
        "ax.set_xlabel('beta_1[0, 0]')\n",
        "ax.set_ylabel('likelihood', color = likelihood_color)\n",
        "ax.plot(beta_1_vals, likelihoods, color = likelihood_color)\n",
        "ax.tick_params(axis='y', labelcolor=likelihood_color)\n",
        "\n",
        "ax1 = ax.twinx()\n",
        "ax1.plot(beta_1_vals, nlls, color = nll_color)\n",
        "ax1.set_ylabel('negative log likelihood', color = nll_color)\n",
        "ax1.tick_params(axis='y', labelcolor = nll_color)\n",
        "\n",
        "plt.axvline(x = beta_1_vals[np.argmax(likelihoods)], linestyle='dotted')\n",
        "\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "aDEPhddNdN4u"
      },
      "outputs": [],
      "source": [
        "# Hopefully, you can see that the maximum of the likelihood fn is at the same position as the minimum negative log likelihood solution\n",
        "# Let's check that:\n",
        "print(\"Maximum likelihood = %f, at beta_1=%3.3f\"%( (likelihoods[np.argmax(likelihoods)],beta_1_vals[np.argmax(likelihoods)])))\n",
        "print(\"Minimum negative log likelihood = %f, at beta_1=%3.3f\"%( (nlls[np.argmin(nlls)],beta_1_vals[np.argmin(nlls)])))\n",
        "\n",
        "# Plot the best model\n",
        "beta_1[0,0] = beta_1_vals[np.argmin(nlls)]\n",
        "model_out = shallow_nn(x_model, beta_0, omega_0, beta_1, omega_1)\n",
        "lambda_model = softmax(model_out)\n",
        "plot_multiclass_classification(x_model, model_out, lambda_model, x_train, y_train, title=\"beta1[0,0]=%3.3f\"%(beta_1[0,0]))\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "771G8N1Vk5A2"
      },
      "source": [
        "They both give the same answer. But you can see from the likelihood above that the likelihood is very small unless the parameters are almost correct.  So in practice, we would work with the negative log likelihood.<br><br>\n",
        "\n",
        "Again, to fit the full neural model we would vary all of the 16 parameters of the network in the $\\boldsymbol\\beta_{0},\\boldsymbol\\Omega_{0},\\boldsymbol\\beta_{1},\\boldsymbol\\Omega_{1}$ until we find the combination that have the maximum likelihood / minimum negative log likelihood.<br><br>\n",
        "\n"
      ]
    }
  ],
  "metadata": {
    "colab": {
      "provenance": [],
      "include_colab_link": true
    },
    "kernelspec": {
      "display_name": "Python 3",
      "name": "python3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}