2_intermediate/1.2-FFN.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Credits\n",
    "\n",
    "See main readme for credits."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import torch\n",
    "from torch.autograd import Variable\n",
    "from torch.nn.parameter import Parameter\n",
    "import torch.nn as nn\n",
    "import torch.nn.functional as F\n",
    "import torch.optim as optim\n",
    "import torch.nn.init as init\n",
    "\n",
    "%matplotlib inline\n",
    "import matplotlib\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# MNIST dataset\n",
    "MNIST is a dataset that is often used for benchmarking. The MNIST dataset consists of 70,000 images of handwritten digits from 0-9. The dataset is split into a 50,000 images training set, 10,000 images validation set and 10,000 images test set. The images are 28x28 pixels, where each pixel represents a normalised value between 0-255 (0=black and 255=white).\n",
    "\n",
    "## Primer\n",
    "We use a feedforward neural network to classify the 28x28 mnist images. `num_features` is therefore $28 * 28=784$, i.e. we represent each image as a vector. The ordering of the pixels in the vector does not matter, so we could permutate all images using the same permutation and still get the same performance. (You are of course encouraged to try this using ``numpy.random.permutation`` to get a random permutation. This task is therefore called the _permutation invariant_ MNIST. Obviously this throws away a lot of structure in the data. In the next module we'll fix this with the convolutional neural network wich encodes prior knowledgde about data that has either spatial or temporal structure.  \n",
    "\n",
    "## Ballpark estimates of hyperparameters\n",
    "__Optimizers:__\n",
    "    1. SGD + Momentum: learning rate 1.0 - 0.1 \n",
    "    2. ADAM: learning rate 3*1e-4 - 1e-5\n",
    "    3. RMSPROP: somewhere between SGD and ADAM\n",
    "\n",
    "__Regularization:__\n",
    "1. [Dropout](http://pytorch.org/docs/master/nn.html?highlight=dropout#torch.nn.Dropout). Dropout rate 0.1-0.5\n",
    "  - Remember to pick the correct version according to the input dimensionality\n",
    "  - **NOTE** call `net.train()` before training to activate random dropout, and call `net.eval()` to deactivate dropout while validating or running inference with model.\n",
    "2. L2 (weight decay of optimization functions, e.g. [Adam](http://pytorch.org/docs/master/optim.html#torch.optim.Adam)) and [L1 regularization](http://pytorch.org/docs/master/nn.html#torch.nn.L1Loss).\n",
    "  - Not used that often in deep learning, but 1e-4  -  1e-8\n",
    "3. [Batchnorm](http://pytorch.org/docs/master/nn.html#torch.nn.BatchNorm1d): Batchnorm also acts as a regularizer - Often very useful (faster and better convergence)\n",
    "  - Remember to pick the correct version according to the input dimensionality\n",
    "  - **NOTE** call `net.train()` before training to activate, and call `net.eval()` to have a non-stochastic variant while validating or running inference with model.\n",
    "    \n",
    "    \n",
    "__Parameter initialization:__\n",
    "Parameter initialization is extremely important. PyTorch has a lot of different initializers, check the [PyTorch API](http://pytorch.org/docs/master/nn.html#torch-nn-init). Often used initializer are\n",
    "1. Kaming He\n",
    "2. Xavier Glorot\n",
    "3. Uniform or Normal with small scale (0.1 - 0.01)\n",
    "4. Orthogonal (this usually works very well for RNNs)\n",
    "\n",
    "Bias is nearly always initialized to zero using the [torch.nn.init.constant(tensor, val)](http://pytorch.org/docs/master/nn.html#torch.nn.init.constant)\n",
    "\n",
    "__Number of hidden units and network structure:__\n",
    "Probably as big network as possible and then apply regularization. You'll have to experiment. One rarely goes below 512 units for feedforward networks unless your are training on CPU...\n",
    "There's some research into stochastic depth networks: https://arxiv.org/pdf/1603.09382v2.pdf, but in general this is trial and error.\n",
    "\n",
    "__Nonlinearity:__ [The most commonly used nonliearities are](http://pytorch.org/docs/master/nn.html#non-linear-activations)\n",
    "    \n",
    "1. ReLU\n",
    "2. Leaky ReLU\n",
    "3. Elu\n",
    "3. Sigmoids squash the output [-1, 1], and are used if your output is binary (not used in the hidden layers)\n",
    "4. Softmax normalizes the the output to 1, and is used as output if you have a classification problem\n",
    "\n",
    "See the plot below.\n",
    "\n",
    "__Mini-batch size:__\n",
    "Usually people use 16-256. Bigger is not allways better. With smaller mini-batch size you get more updates and your model might converge faster. Also small batch sizes use less memory, which means you can train a model with more parameters.\n",
    "\n",
    "Hyperparameters can be found by experience (guessing) or some search procedure. Random search is easy to implement and performs decent: http://www.jmlr.org/papers/volume13/bergstra12a/bergstra12a.pdf . \n",
    "More advanced search procedures include [Spearmint](https://github.com/JasperSnoek/spearmint) and many others."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Illustrate different output units\n",
    "x = np.linspace(-6, 6, 100)\n",
    "units = {\n",
    "    \"ReLU\": lambda x: np.maximum(0, x),\n",
    "    \"Leaky ReLU\": lambda x: np.maximum(0, x) + 0.1 * np.minimum(0, x),\n",
    "    \"Elu\": lambda x: (x > 0) * x + (1 - (x > 0)) * (np.exp(x) - 1),\n",
    "    \"Sigmoid\": lambda x: (1 + np.exp(-x))**(-1),\n",
    "    \"tanh\": lambda x: (np.exp(x) - np.exp(-x)) / (np.exp(x) + np.exp(-x))\n",
    "}\n",
    "\n",
    "plt.figure(figsize=(5, 5))\n",
    "[plt.plot(x, unit(x), label=unit_name, lw=2) for unit_name, unit in units.items()]\n",
    "plt.legend(loc=2, fontsize=16)\n",
    "plt.title('Non-linearities', fontsize=20)\n",
    "plt.ylim([-2, 5])\n",
    "plt.xlim([-6, 6])\n",
    "\n",
    "# assert that all class probablities sum to one\n",
    "softmax = lambda x: np.exp(x) / np.sum(np.exp(x))\n",
    "print(\"softmax should sum to one (approxiamtely):\", np.sum(softmax(x)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "## MNIST\n",
    "First let's load the MNIST dataset and plot a few examples:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!if [ ! -f mnist.npz ]; then wget -N https://www.dropbox.com/s/qxywaq7nx19z72p/mnist.npz; else echo \"mnist.npz already downloaded\"; fi"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#To speed up training we'll only work on a subset of the data\n",
    "data = np.load('mnist.npz')\n",
    "num_classes = 10\n",
    "x_train = data['X_train'][:1000].astype('float32')\n",
    "targets_train = data['y_train'][:1000].astype('int32')\n",
    "\n",
    "x_valid = data['X_valid'][:500].astype('float32')\n",
    "targets_valid = data['y_valid'][:500].astype('int32')\n",
    "\n",
    "x_test = data['X_test'][:500].astype('float32')\n",
    "targets_test = data['y_test'][:500].astype('int32')\n",
    "\n",
    "print(\"Information on dataset\")\n",
    "print(\"x_train\", x_train.shape)\n",
    "print(\"targets_train\", targets_train.shape)\n",
    "print(\"x_valid\", x_valid.shape)\n",
    "print(\"targets_valid\", targets_valid.shape)\n",
    "print(\"x_test\", x_test.shape)\n",
    "print(\"targets_test\", targets_test.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "#plot a few MNIST examples\n",
    "idx, dim, classes = 0, 28, 10\n",
    "# create empty canvas\n",
    "canvas = np.zeros((dim*classes, classes*dim))\n",
    "\n",
    "# fill with tensors\n",
    "for i in range(classes):\n",
    "    for j in range(classes):\n",
    "        canvas[i*dim:(i+1)*dim, j*dim:(j+1)*dim] = x_train[idx].reshape((dim, dim))\n",
    "        idx += 1\n",
    "\n",
    "# visualize matrix of tensors as gray scale image\n",
    "plt.figure(figsize=(6, 6))\n",
    "plt.axis('off')\n",
    "plt.imshow(canvas, cmap='gray')\n",
    "plt.title('MNIST handwritten digits')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Hyperparameters\n",
    "num_classes = 10\n",
    "num_l1 = 512\n",
    "num_features = x_train.shape[1]\n",
    "\n",
    "# define network\n",
    "class Net(nn.Module):\n",
    "\n",
    "    def __init__(self, num_features, num_hidden, num_output):\n",
    "        super(Net, self).__init__()  \n",
    "        # input layer\n",
    "        self.W_1 = Parameter(init.xavier_normal(torch.Tensor(num_hidden, num_features)))\n",
    "        self.b_1 = Parameter(init.constant(torch.Tensor(num_hidden), 0))\n",
    "        # hidden layer\n",
    "        self.W_2 = Parameter(init.xavier_normal(torch.Tensor(num_output, num_hidden)))\n",
    "        self.b_2 = Parameter(init.constant(torch.Tensor(num_output), 0))\n",
    "        # define activation function in constructor\n",
    "        self.activation = torch.nn.ELU()\n",
    "\n",
    "    def forward(self, x):\n",
    "        x = F.linear(x, self.W_1, self.b_1)\n",
    "        x = self.activation(x)\n",
    "        x = F.linear(x, self.W_2, self.b_2)\n",
    "        return F.softmax(x, dim=1)\n",
    "\n",
    "\n",
    "net = Net(num_features, num_l1, num_classes)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "optimizer = optim.SGD(net.parameters(), lr=0.1)\n",
    "criterion = nn.CrossEntropyLoss()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Test the forward pass with dummy data\n",
    "x = np.random.normal(0, 1, (45, dim*dim)).astype('float32')\n",
    "\n",
    "print(net(Variable(torch.from_numpy(x))).size())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Build the training loop\n",
    "\n",
    "We train the network by calculating the gradient w.r.t the cost function and update the parameters in direction of the negative gradient. \n",
    "\n",
    "\n",
    "When training neural network you always use mini batches. Instead of calculating the average gradient using the entire dataset you approximate the gradient using a mini-batch of typically 16 to 256 samples. The paramters are updated after each mini batch. Networks converge much faster using mini batches because the parameters are updated more often.\n",
    "\n",
    "We build a loop that iterates over the training data. Remember that the parameters are updated each time ``optimizer.step()`` is called."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# we could have done this ourselves,\n",
    "# but we should be aware of sklearn and it's tools\n",
    "from sklearn.metrics import accuracy_score\n",
    "\n",
    "# setting hyperparameters and gettings epoch sizes\n",
    "batch_size = 100\n",
    "num_epochs = 100\n",
    "num_samples_train = x_train.shape[0]\n",
    "num_batches_train = num_samples_train // batch_size\n",
    "num_samples_valid = x_valid.shape[0]\n",
    "num_batches_valid = num_samples_valid // batch_size\n",
    "\n",
    "# setting up lists for handling loss/accuracy\n",
    "train_acc, train_loss = [], []\n",
    "valid_acc, valid_loss = [], []\n",
    "test_acc, test_loss = [], []\n",
    "cur_loss = 0\n",
    "losses = []\n",
    "\n",
    "get_slice = lambda i, size: range(i * size, (i + 1) * size)\n",
    "\n",
    "for epoch in range(num_epochs):\n",
    "    # Forward -> Backprob -> Update params\n",
    "    ## Train\n",
    "    cur_loss = 0\n",
    "    net.train()\n",
    "    for i in range(num_batches_train):\n",
    "        slce = get_slice(i, batch_size)\n",
    "        x_batch = Variable(torch.from_numpy(x_train[slce]))\n",
    "        output = net(x_batch)\n",
    "        \n",
    "        # compute gradients given loss\n",
    "        target_batch = Variable(torch.from_numpy(targets_train[slce]).long())\n",
    "        batch_loss = criterion(output, target_batch)\n",
    "        optimizer.zero_grad()\n",
    "        batch_loss.backward()\n",
    "        optimizer.step()\n",
    "        \n",
    "        cur_loss += batch_loss   \n",
    "    losses.append(cur_loss / batch_size)\n",
    "\n",
    "    net.eval()\n",
    "    ### Evaluate training\n",
    "    train_preds, train_targs = [], []\n",
    "    for i in range(num_batches_train):\n",
    "        slce = get_slice(i, batch_size)\n",
    "        x_batch = Variable(torch.from_numpy(x_train[slce]))\n",
    "        \n",
    "        output = net(x_batch)\n",
    "        preds = torch.max(output, 1)[1]\n",
    "        \n",
    "        train_targs += list(targets_train[slce])\n",
    "        train_preds += list(preds.data.numpy())\n",
    "    \n",
    "    ### Evaluate validation\n",
    "    val_preds, val_targs = [], []\n",
    "    for i in range(num_batches_valid):\n",
    "        slce = get_slice(i, batch_size)\n",
    "        x_batch = Variable(torch.from_numpy(x_valid[slce]))\n",
    "        \n",
    "        output = net(x_batch)\n",
    "        preds = torch.max(output, 1)[1]\n",
    "        val_preds += list(preds.data.numpy())\n",
    "        val_targs += list(targets_valid[slce])\n",
    "\n",
    "    train_acc_cur = accuracy_score(train_targs, train_preds)\n",
    "    valid_acc_cur = accuracy_score(val_targs, val_preds)\n",
    "    \n",
    "    train_acc.append(train_acc_cur)\n",
    "    valid_acc.append(valid_acc_cur)\n",
    "    \n",
    "    if epoch % 10 == 0:\n",
    "        print(\"Epoch %2i : Train Loss %f , Train acc %f, Valid acc %f\" % (\n",
    "                epoch+1, losses[-1], train_acc_cur, valid_acc_cur))\n",
    "\n",
    "epoch = np.arange(len(train_acc))\n",
    "plt.figure()\n",
    "plt.plot(epoch, train_acc, 'r', epoch, valid_acc, 'b')\n",
    "plt.legend(['Train Accucary','Validation Accuracy'])\n",
    "plt.xlabel('Updates'), plt.ylabel('Acc')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "# Questions\n",
    "\n",
    "Try and add these modifications:\n",
    "- Kaiming He initialization instead of Xavier Glorot\n",
    "- add an extra layer\n",
    "- use the relu activation function\n",
    "- add dropout to the network (**note** the `net.train()` and `net.eval()` already in the code)\n",
    "- add momentum to the optimizer\n",
    "- use the ADAM optimizer instead of stochastic gradient descent\n",
    "- add L2/weight decay to the optimizer\n",
    "- add L1 regularization"
   ]
  }
 ],
 "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.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}