lab2/lab2_CNN.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import matplotlib\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import sklearn.datasets\n",
    "import theano\n",
    "import theano.tensor as T\n",
    "import lasagne\n",
    " "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Convolutional Neural networks 101\n",
    "\n",
    "Convolution neural networks are one of the most succesfull types of neural networks for image recognition and an integral part of reigniting the interest in neural networks. \n",
    "\n",
    "In this lab we'll experiment with inserting 2D-convolution layers in the fully connected neural networks introduced in LAB1. We'll furhter experiment with stacking of convolution layers, max pooling and strided convolutions which are all important techniques in current convolution neural network architectures. Lastly we'll try to visualize the learned convolution filters and try to understand what kind of features they learn to recognize.\n",
    "\n",
    "\n",
    "If you are unfamilar with the the convolution operation  https://github.com/vdumoulin/conv_arithmetic have a nice visualization of different convolution variants. For a more indept tutorial please see http://cs231n.github.io/convolutional-networks/ or http://neuralnetworksanddeeplearning.com/chap6.html. Lastly if you are ambitious and want implement a convolution neural network from scratch please see an exercise for our Deep Learning summer school last year https://github.com/DTU-deeplearning/day2-Conv"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#LOAD the mnist data. To speed up training we'll only work on a subset of the data.\n",
    "#Note that we reshape the data from (nsamples, num_features)= (nsamples, nchannels*rows*cols)  -> (nsamples, nchannels, rows, cols)\n",
    "# in order to retain the spatial arrangements of the pixels\n",
    "data = np.load('mnist.npz')\n",
    "num_classes = 10\n",
    "nchannels,rows,cols = 1,28,28\n",
    "x_train = data['X_train'][:10000].astype('float32')\n",
    "x_train = x_train.reshape((-1,nchannels,rows,cols))\n",
    "targets_train = data['y_train'][:10000].astype('int32')\n",
    "\n",
    "x_valid = data['X_valid'][:500].astype('float32')\n",
    "x_valid = x_valid.reshape((-1,nchannels,rows,cols))\n",
    "targets_valid = data['y_valid'][:500].astype('int32')\n",
    "\n",
    "x_test = data['X_test'][:500].astype('float32')\n",
    "x_test = x_test.reshape((-1,nchannels,rows,cols))\n",
    "targets_test = data['y_test'][:500].astype('int32')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#plot a few MNIST examples\n",
    "idx = 0\n",
    "canvas = np.zeros((28*10, 10*28))\n",
    "for i in range(10):\n",
    "    for j in range(10):\n",
    "        canvas[i*28:(i+1)*28, j*28:(j+1)*28] = x_train[idx].reshape((28, 28))\n",
    "        idx += 1\n",
    "plt.figure(figsize=(7, 7))\n",
    "plt.imshow(canvas, cmap='gray')\n",
    "plt.title('MNIST handwritten digits')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#Define a simple feed forward neural network\n",
    "\n",
    "from lasagne.nonlinearities import leaky_rectify, softmax, tanh, elu\n",
    "from lasagne.layers import InputLayer, DenseLayer, Conv2DLayer, batch_norm, DropoutLayer, MaxPool2DLayer\n",
    "\n",
    "#defined the model\n",
    "num_class = 10\n",
    "num_features = x_train.shape[1]\n",
    "\n",
    "l_in = InputLayer(shape=(None,nchannels,rows,cols)) #note that we use a 4D input since we need to retain the spatial arrangement of the pixels when working with convolutions.\n",
    "#l_conv = Conv2DLayer(l_in,num_filters=16,filter_size=5)\n",
    "l_hid = DenseLayer(l_in, num_units=100, nonlinearity=elu) #remember to connect the new conv-layer here\n",
    "l_out = DenseLayer(l_hid, num_units=num_class, nonlinearity=softmax)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#Setting up the graph in theano\n",
    "sym_x = T.tensor4('sym_x') # a symbolic variable, this is now a 4-D tensor.\n",
    "sym_t = T.ivector('sym_t') # a symbolic variable taking on the value of the target batch.\n",
    "\n",
    "# Get network output\n",
    "train_out = lasagne.layers.get_output(l_out, sym_x, deterministic=False)\n",
    "eval_out = lasagne.layers.get_output(l_out, sym_x, deterministic=True)\n",
    "\n",
    "\n",
    "# Get list of all trainable parameters in the network.\n",
    "all_params = lasagne.layers.get_all_params(l_out, trainable=True)\n",
    "\n",
    "cost = T.nnet.categorical_crossentropy(train_out+1e-8, sym_t).mean()\n",
    "# Let Theano do its magic and get all the gradients we need for training\n",
    "all_grads = T.grad(cost, all_params)\n",
    "\n",
    "\n",
    "# Set the update function for parameters \n",
    "# you might wan't to experiment with more advanded update schemes like rmsprob, adadelta etc.\n",
    "updates = lasagne.updates.adam(all_grads, all_params, learning_rate=0.001)\n",
    "\n",
    "\n",
    "f_eval = theano.function([sym_x],\n",
    "                     eval_out, on_unused_input='warn')\n",
    "\n",
    "f_train = theano.function([sym_x, sym_t],\n",
    "                          [cost],\n",
    "                          updates=updates, on_unused_input='warn')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#Test the forward pass\n",
    "x = np.random.normal(0,1, (45, 1,28,28)).astype('float32') #dummy data\n",
    "\n",
    "model = lasagne.layers.get_output(l_out, sym_x)\n",
    "out = model.eval({sym_x:x}) #this could also include mask etc if used\n",
    "print(\"l_out\", out.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#Training Loop\n",
    "from confusionmatrix import ConfusionMatrix\n",
    "batch_size = 100\n",
    "num_epochs = 10\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",
    "train_acc, train_loss = [], []\n",
    "valid_acc, valid_loss = [], []\n",
    "test_acc, test_loss = [], []\n",
    "cur_loss = 0\n",
    "loss = []\n",
    "try:\n",
    "    for epoch in range(num_epochs):\n",
    "        #Forward->Backprob->Update params\n",
    "        cur_loss = 0\n",
    "        for i in range(num_batches_train):\n",
    "            idx = range(i*batch_size, (i+1)*batch_size)\n",
    "            x_batch = x_train[idx]\n",
    "            target_batch = targets_train[idx]    \n",
    "            batch_loss = f_train(x_batch,target_batch) #this will do the complete backprob pass\n",
    "            cur_loss += batch_loss[0]\n",
    "        loss += [cur_loss/batch_size]\n",
    "\n",
    "        confusion_valid = ConfusionMatrix(num_classes)\n",
    "        confusion_train = ConfusionMatrix(num_classes)\n",
    "\n",
    "        for i in range(num_batches_train):\n",
    "            idx = range(i*batch_size, (i+1)*batch_size)\n",
    "            x_batch = x_train[idx]\n",
    "            targets_batch = targets_train[idx]\n",
    "            net_out = f_eval(x_batch)   \n",
    "            preds = np.argmax(net_out, axis=-1) \n",
    "            confusion_train.batch_add(targets_batch, preds)\n",
    "\n",
    "        confusion_valid = ConfusionMatrix(num_classes)\n",
    "        for i in range(num_batches_valid):\n",
    "            idx = range(i*batch_size, (i+1)*batch_size)\n",
    "            x_batch = x_valid[idx]\n",
    "            targets_batch = targets_valid[idx]\n",
    "            net_out = f_eval(x_batch)   \n",
    "            preds = np.argmax(net_out, axis=-1) \n",
    "\n",
    "            confusion_valid.batch_add(targets_batch, preds)\n",
    "\n",
    "        train_acc_cur = confusion_train.accuracy()\n",
    "        valid_acc_cur = confusion_valid.accuracy()\n",
    "\n",
    "        train_acc += [train_acc_cur]\n",
    "        valid_acc += [valid_acc_cur]\n",
    "        print \"Epoch %i : Train Loss %e , Train acc %f,  Valid acc %f \" \\\n",
    "        % (epoch+1, loss[-1], train_acc_cur, valid_acc_cur)\n",
    "except KeyboardInterrupt:\n",
    "    pass\n",
    "    \n",
    "\n",
    "#get test set score\n",
    "confusion_test = ConfusionMatrix(num_classes)\n",
    "net_out = f_eval(x_test)    \n",
    "preds = np.argmax(net_out, axis=-1) \n",
    "confusion_test.batch_add(targets_test, preds)\n",
    "print \"\\nTest set Acc:  %f\" %(confusion_test.accuracy())\n",
    "\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 Acc','Val Acc'])\n",
    "plt.xlabel('Epochs'), plt.ylabel('Acc'), plt.ylim([0.75,1.03])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Assignments 1\n",
    "\n",
    " 1) Note the performance of the standard feedforward neural network. Add a 2D convolution layer before the dense hidden layer and confirm that it increases the generalization performance of the network (try num_filters=16 and filter_size=5 as a starting point). \n",
    " \n",
    " 2) Can the performance be increases even further by stacking more convolution layers ?\n",
    " \n",
    " 3) Maxpooling is a technique for decreasing the spatial resolution of an image while retaining the important features. Effectively this gives a local translational invariance and reduces the computation by a factor of four. In the classification algorithm which is usually desirable. Try to either: \n",
    " \n",
    "     a) add a maxpool layer(add arguement pool_size=2)  after the convolution layer or\n",
    "     b) set add stride=2 to the arguments of the convolution layer. \n",
    "  Verify that this decreases spatial dimension of the image. (print l_conv.output_shape or print   l_maxpool.output_shape). Does this increase the performance of the network (you may need to stack multiple layers or increase the number of filters to increase performance) ?\n",
    "  \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Visualization of filters\n",
    "Convolution filters can be interpreted as spatial feature detectors picking up different image features such as edges, corners etc. Below we provide code for visualization of the filters. The best results are obtained with fairly large filters of size 9 and either 16 or 36 filters. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "### If you get an error with l_conv not being defined you need define l_conv when the network is defined!\n",
    "np_W = l_conv.W.get_value() #get the filter values from the conv layer\n",
    "print np_W.shape, \"i.e. the shape is num_filters, num_channels, filter_size, filter_size\"\n",
    "num_filters,num_channels,filter_size,_= np_W.shape\n",
    "n = int(num_filters**0.5)\n",
    "\n",
    "np_W_res = np_W.reshape(n,n,num_channels,filter_size,filter_size)\n",
    "fig, ax = plt.subplots(n,n)\n",
    "print \"learned filter values\"\n",
    "for i in range(n):\n",
    "    for j in range(n):\n",
    "        ax[i,j].imshow(np_W_res[i,j,0], cmap='gray',interpolation='none')\n",
    "        ax[i,j].xaxis.set_major_formatter(plt.NullFormatter())\n",
    "        ax[i,j].yaxis.set_major_formatter(plt.NullFormatter())\n",
    "\n",
    "\n",
    "idx = 1\n",
    "plt.figure()\n",
    "plt.imshow(x_train[idx,0],cmap='gray',interpolation='none')\n",
    "plt.title('Inut Image')\n",
    "plt.show()\n",
    "\n",
    "#visalize the filters convolved with an input image\n",
    "from scipy.signal import convolve2d\n",
    "np_W_res = np_W.reshape(n,n,num_channels,filter_size,filter_size)\n",
    "fig, ax = plt.subplots(n,n,figsize=(9,9))\n",
    "print \"Response from input image convolved with the filters\"\n",
    "for i in range(n):\n",
    "    for j in range(n):\n",
    "        ax[i,j].imshow(convolve2d(x_train[1,0],np_W_res[i,j,0],mode='same'), cmap='gray',interpolation='none')\n",
    "        ax[i,j].xaxis.set_major_formatter(plt.NullFormatter())\n",
    "        ax[i,j].yaxis.set_major_formatter(plt.NullFormatter())\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Assignment 2\n",
    "\n",
    "The visualized filters will likely look most like noise due to the small amount of training data.\n",
    "\n",
    " 1) Try to use 10000 traning examples instead and visualise the filters again\n",
    " \n",
    " 2) Dropout is a very usefull technique for preventing overfitting. Try to add a DropoutLayer after the convolution layer and hidden layer. This should increase both performance and the \"visual appeal\" of the filters\n",
    " \n",
    " 3) Batch normalization is a recent innovation for improving generalization performance. Try to insert batch normalization layers into the network to improve performance. \n",
    " \n",
    " \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# More Fun with convolutional networks\n",
    "### Get the data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "!wget -N https://s3.amazonaws.com/lasagne/recipes/datasets/mnist_cluttered_60x60_6distortions.npz"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the data the each mnist digit (20x20 pixels) has been placed randomly in a 60x60 canvas. To make the task harder each canvas has then been cluttered with small pieces of digits. In this task it is helpfull for a network if it can focus only on the digit and ignore the rest.\n",
    "\n",
    "The ``TransformerLayer`` lets us do this. The transformer layer learns an affine transformation which lets the network zoom, rotate and skew. If you are interested you should read the paper, but the main idea is that you can let a small convolutional network determine the the parameters of the affine transformation. You then apply the affine transformation to the input data. Usually this also involves downsampling which forces the model to zoom in on the relevant parts of the data. After the affine transformation we can use a larger conv net to do the classification. \n",
    "This is possible because you can backprop through a an affine transformation if you use bilinear interpolation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import os\n",
    "import matplotlib\n",
    "import numpy as np\n",
    "np.random.seed(123)\n",
    "import matplotlib.pyplot as plt\n",
    "import lasagne\n",
    "import theano\n",
    "import theano.tensor as T\n",
    "conv = lasagne.layers.Conv2DLayer\n",
    "pool = lasagne.layers.MaxPool2DLayer\n",
    "NUM_EPOCHS = 500\n",
    "BATCH_SIZE = 256\n",
    "LEARNING_RATE = 0.001\n",
    "DIM = 60\n",
    "NUM_CLASSES = 10\n",
    "mnist_cluttered = \"mnist_cluttered_60x60_6distortions.npz\"\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def load_data():\n",
    "    data = np.load(mnist_cluttered)\n",
    "    X_train, y_train = data['x_train'], np.argmax(data['y_train'], axis=-1)\n",
    "    X_valid, y_valid = data['x_valid'], np.argmax(data['y_valid'], axis=-1)\n",
    "    X_test, y_test = data['x_test'], np.argmax(data['y_test'], axis=-1)\n",
    "\n",
    "    # reshape for convolutions\n",
    "    X_train = X_train.reshape((X_train.shape[0], 1, DIM, DIM))\n",
    "    X_valid = X_valid.reshape((X_valid.shape[0], 1, DIM, DIM))\n",
    "    X_test = X_test.reshape((X_test.shape[0], 1, DIM, DIM))\n",
    "    \n",
    "    print \"Train samples:\", X_train.shape\n",
    "    print \"Validation samples:\", X_valid.shape\n",
    "    print \"Test samples:\", X_test.shape\n",
    "\n",
    "    return dict(\n",
    "        X_train=lasagne.utils.floatX(X_train),\n",
    "        y_train=y_train.astype('int32'),\n",
    "        X_valid=lasagne.utils.floatX(X_valid),\n",
    "        y_valid=y_valid.astype('int32'),\n",
    "        X_test=lasagne.utils.floatX(X_test),\n",
    "        y_test=y_test.astype('int32'),\n",
    "        num_examples_train=X_train.shape[0],\n",
    "        num_examples_valid=X_valid.shape[0],\n",
    "        num_examples_test=X_test.shape[0],\n",
    "        input_height=X_train.shape[2],\n",
    "        input_width=X_train.shape[3],\n",
    "        output_dim=10,)\n",
    "data = load_data()\n",
    "\n",
    "idx = 0\n",
    "canvas = np.zeros((DIM*10, 10*DIM))\n",
    "for i in range(10):\n",
    "    for j in range(10):\n",
    "        canvas[i*DIM:(i+1)*DIM, j*DIM:(j+1)*DIM] = data['X_train'][idx].reshape((DIM, DIM))\n",
    "        idx += 1\n",
    "plt.figure(figsize=(10, 10))\n",
    "plt.imshow(canvas, cmap='gray')\n",
    "plt.title('Cluttered handwritten digits')\n",
    "plt.axis('off')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Building the model\n",
    "\n",
    "We use a model where the localization network is a two layer convolution network which operates directly on the image input. The output from the localization network is a 6 dimensional vector specifying the parameters in the affine transformation.\n",
    "\n",
    "We set up the transformer layer to initially do the identity transform, similarly to [1]. If the output from the localization networks is [t1, t2, t3, t4, t5, t6] then t1 and t5 determines zoom, t2 and t4 determines skewness, and t3 and t6 move the center position. By setting the initial values of the bias vector to \n",
    "\n",
    "```\n",
    "|1, 0, 0|\n",
    "|0, 1, 0|\n",
    "```\n",
    "and the final W of the localization network to all zeros we ensure that in the beginning of training the network works as a pooling layer. \n",
    "\n",
    "The output of the localization layer feeds into the transformer layer which applies the transformation to the image input. In our setup the transformer layer downsamples the input by a factor 3.\n",
    "\n",
    "Finally a 2 layer convolution layer and 2 fully connected layers calculates the output probabilities.\n",
    "\n",
    "\n",
    "### The model\n",
    "```\n",
    "Input -> localization_network -> TransformerLayer -> output_network -> predictions\n",
    "   |                                |\n",
    "   >--------------------------------^\n",
    "```\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def build_model(input_width, input_height, output_dim,\n",
    "                batch_size=BATCH_SIZE):\n",
    "    ini = lasagne.init.HeUniform()\n",
    "    l_in = lasagne.layers.InputLayer(shape=(None, 1, input_width, input_height),)\n",
    "\n",
    "    # Localization network\n",
    "    b = np.zeros((2, 3), dtype=theano.config.floatX)\n",
    "    b[0, 0] = 1\n",
    "    b[1, 1] = 1\n",
    "    b = b.flatten()\n",
    "    loc_l1 = pool(l_in, pool_size=(2, 2))\n",
    "    loc_l2 = conv(\n",
    "        loc_l1, num_filters=8, filter_size=(5, 5), W=ini)\n",
    "    loc_l3 = pool(loc_l2, pool_size=(2, 2))\n",
    "    loc_l4 = conv(loc_l3, num_filters=8, filter_size=(5, 5), W=ini)\n",
    "    loc_l5 = lasagne.layers.DenseLayer(\n",
    "        loc_l4, num_units=50, W=lasagne.init.HeUniform('relu'))\n",
    "    loc_out = lasagne.layers.DenseLayer(\n",
    "        loc_l5, num_units=6, b=b, W=lasagne.init.Constant(0.0), \n",
    "        nonlinearity=lasagne.nonlinearities.identity)\n",
    "    \n",
    "    # Transformer network\n",
    "    l_trans1 = lasagne.layers.TransformerLayer(l_in, loc_out, downsample_factor=3.0)\n",
    "    print \"Transformer network output shape: \", l_trans1.output_shape\n",
    "    \n",
    "    # Classification network\n",
    "    class_l1 = conv(\n",
    "        l_trans1,\n",
    "        num_filters=16,\n",
    "        filter_size=(3, 3),\n",
    "        nonlinearity=lasagne.nonlinearities.rectify,\n",
    "        W=ini,\n",
    "    )\n",
    "    class_l2 = pool(class_l1, pool_size=(2, 2))\n",
    "    class_l3 = conv(\n",
    "        class_l2,\n",
    "        num_filters=16,\n",
    "        filter_size=(3, 3),\n",
    "        nonlinearity=lasagne.nonlinearities.rectify,\n",
    "        W=ini,\n",
    "    )\n",
    "    class_l4 = pool(class_l3, pool_size=(2, 2))\n",
    "    class_l5 = lasagne.layers.DenseLayer(\n",
    "        class_l4,\n",
    "        num_units=256,\n",
    "        nonlinearity=lasagne.nonlinearities.rectify,\n",
    "        W=ini,\n",
    "    )\n",
    "\n",
    "    l_out = lasagne.layers.DenseLayer(\n",
    "        class_l5,\n",
    "        num_units=output_dim,\n",
    "        nonlinearity=lasagne.nonlinearities.softmax,\n",
    "        W=ini,\n",
    "    )\n",
    "\n",
    "    return l_out, l_trans1\n",
    "\n",
    "model, l_transform = build_model(DIM, DIM, NUM_CLASSES)\n",
    "model_params = lasagne.layers.get_all_params(model, trainable=True)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X = T.tensor4()\n",
    "y = T.ivector()\n",
    "\n",
    "# training output\n",
    "output_train = lasagne.layers.get_output(model, X, deterministic=False)\n",
    "\n",
    "# evaluation output. Also includes output of transform for plotting\n",
    "output_eval, transform_eval = lasagne.layers.get_output([model, l_transform], X, deterministic=True)\n",
    "\n",
    "sh_lr = theano.shared(lasagne.utils.floatX(LEARNING_RATE))\n",
    "cost = T.mean(T.nnet.categorical_crossentropy(output_train, y))\n",
    "updates = lasagne.updates.adam(cost, model_params, learning_rate=sh_lr)\n",
    "\n",
    "train = theano.function([X, y], [cost, output_train], updates=updates)\n",
    "eval = theano.function([X], [output_eval, transform_eval])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Training the model\n",
    "Unfortunately NVIDIA has yet to squeeze a TitanX into a labtop and training convnets on CPU is painfully slow. After 10 epochs you should see that model starts to zoom in on the digits. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def train_epoch(X, y):\n",
    "    num_samples = X.shape[0]\n",
    "    num_batches = int(np.ceil(num_samples / float(BATCH_SIZE)))\n",
    "    costs = []\n",
    "    correct = 0\n",
    "    for i in range(num_batches):\n",
    "        if i % 10 == 0:\n",
    "            print i,\n",
    "        idx = range(i*BATCH_SIZE, np.minimum((i+1)*BATCH_SIZE, num_samples))\n",
    "        X_batch = X[idx]\n",
    "        y_batch = y[idx]\n",
    "        cost_batch, output_train = train(X_batch, y_batch)\n",
    "        costs += [cost_batch]\n",
    "        preds = np.argmax(output_train, axis=-1)\n",
    "        correct += np.sum(y_batch == preds)\n",
    "    print \"\"\n",
    "    return np.mean(costs), correct / float(num_samples)\n",
    "\n",
    "\n",
    "def eval_epoch(X, y):\n",
    "    output_eval, transform_eval = eval(X)\n",
    "    preds = np.argmax(output_eval, axis=-1)\n",
    "    acc = np.mean(preds == y)\n",
    "    return acc, transform_eval"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "valid_accs, train_accs, test_accs = [], [], []\n",
    "try:\n",
    "    for n in range(NUM_EPOCHS):\n",
    "        train_cost, train_acc = train_epoch(data['X_train'], data['y_train'])\n",
    "        valid_acc, valid_trainsform = eval_epoch(data['X_valid'], data['y_valid'])\n",
    "        test_acc, test_transform = eval_epoch(data['X_test'], data['y_test'])\n",
    "        valid_accs += [valid_acc]\n",
    "        test_accs += [test_acc]\n",
    "        train_accs += [train_acc]\n",
    "\n",
    "        if (n+1) % 20 == 0:\n",
    "            new_lr = sh_lr.get_value() * 0.7\n",
    "            print \"New LR:\", new_lr\n",
    "            sh_lr.set_value(lasagne.utils.floatX(new_lr))\n",
    "\n",
    "        print \"Epoch {0}: Train cost {1}, Train acc {2}, val acc {3}, test acc {4}\".format(\n",
    "                n, train_cost, train_acc, valid_acc, test_acc)\n",
    "except KeyboardInterrupt:\n",
    "    pass"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Plot errors and zoom"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "plt.figure(figsize=(9,9))\n",
    "plt.plot(1-np.array(train_accs), label='Training Error')\n",
    "plt.plot(1-np.array(valid_accs), label='Validation Error')\n",
    "plt.legend(fontsize=20)\n",
    "plt.xlabel('Epoch', fontsize=20)\n",
    "plt.ylabel('Error', fontsize=20)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "plt.figure(figsize=(7,14))\n",
    "for i in range(3):\n",
    "    plt.subplot(321+i*2)\n",
    "    plt.imshow(data['X_test'][i].reshape(DIM, DIM), cmap='gray', interpolation='none')\n",
    "    if i == 0:\n",
    "        plt.title('Original 60x60', fontsize=20)\n",
    "    plt.axis('off')\n",
    "    plt.subplot(322+i*2)\n",
    "    plt.imshow(test_transform[i].reshape(DIM//3, DIM//3), cmap='gray', interpolation='none')\n",
    "    if i == 0:\n",
    "        plt.title('Transformed 20x20', fontsize=20)\n",
    "    plt.axis('off')\n",
    "    \n",
    "    \n",
    "plt.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "# A few pointers for image classification\n",
    "If you want do image classification using a pretrained model is often a good choice, especially if you have limited amounts of labeled data. \n",
    "\n",
    "An often used pretrained network is the VGG16 and VGG19. Lasagne has pretrained models in the [modelzoo](https://github.com/Lasagne/Recipes/tree/master/modelzoo). Torch7 and Tensorflow have similar pretrained models that you can find with google. \n",
    "\n",
    "Currently the best performing image networks is [ResNet](https://arxiv.org/pdf/1512.03385v1.pdf). Torch7 has an interesting blog post about Residual nets. http://torch.ch/blog/2016/02/04/resnets.html\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}