5.1-theano-and-neural-networks.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Necessary libraries:\n",
    "* `gcc` (if you're using windows, it's not that simple).\n",
    "* `theano`: `!pip install theano`\n",
    "* `keras`: `!pip install keras`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Low-level software to implement Neural networks:\n",
    "* Theano is a python library and an optimizing compiler, which allows to define, optimize and compute mathematical expressions effectively using multidimensional arrays.\n",
    "    * serves as a low-level backend for many libraries: keras, lasagne, blocks and many others\n",
    "    * very flexible, though slower then competitors\n",
    "    * awesome for research\n",
    "* Tensorflow\n",
    "    * also quite low-level, but has some of standard deep learning 'elements' built-in\n",
    "    * considered as a synonym for deep learning by people far from machine learning\n",
    "    * quite fast and nice if you want to use ML in production. However, mxnet / torch typically not worse\n",
    "    * better flexibility compared to mxnet / torch / caffe."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "from IPython.display import SVG\n",
    "import numpy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Import "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import theano\n",
    "# next line is a convention\n",
    "import theano.tensor as T"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Symbolic expessions for tensors\n",
    "Theano’s strength is in expressing symbolic calculations involving tensors. There are many types of symbolic expressions for tensors:\n",
    "\n",
    "* scalar\n",
    "* vector\n",
    "* matrix\n",
    "* tensor (3, 4 and more dimensions)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# create two vectors and one scalar. Those are abstract variables\n",
    "x_var = T.vector() \n",
    "y_var = T.vector()\n",
    "alpha_var = T.scalar()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### compute \n",
    "$z = x + \\alpha y + (\\sum x_i, ... ,\\sum x_i)^T$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# define mathematical expression (you can use any function, which you use for numpy arrays). \n",
    "z_var = x_var + alpha_var * y_var + T.sum(x_var)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compile defined expression\n",
    "So far we used \"symbolic\" variables and symbolic expressions defining the expression for computation, but not computing anything. To use the \"expression\", one should compile it.\n",
    "\n",
    "We started with defining mathematical function, process of compiling turns this into programmer's function (which is able to efficiently calculate outputs given the input values).\n",
    "\n",
    "`theano.function` returns a callable object (=function) that will calculate outputs from inputs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# input variables, output expessions\n",
    "inputs = [x_var, y_var, alpha_var]\n",
    "outputs = z_var\n",
    "# allow_input_downcast - automatic type casting for input parameters (e.g. float64 -> float32)\n",
    "compiled_expr = theano.function(inputs, outputs, allow_input_downcast=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compute compiled expression"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "alpha_val = 0.5\n",
    "\n",
    "# using function with lists:\n",
    "print \"using python lists:\"\n",
    "print compiled_expr([1, 2, 3], [4, 5, 6], alpha_val)\n",
    "print \"\\n\"\n",
    "\n",
    "# Or using numpy arrays (should be preferred):\n",
    "# BTW, that 'float' dtype is casted to second parameter dtype which is float32\n",
    "print \"using numpy arrays:\"\n",
    "print compiled_expr(numpy.arange(10),\n",
    "                    numpy.linspace(5, 6, 10, dtype='float'), alpha_val)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercises"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# create two vectors and two scalars\n",
    "x_var = T.vector()\n",
    "y_var = T.vector()\n",
    "alpha_var = T.scalar()\n",
    "beta_var = T.scalar()\n",
    "\n",
    "# define values for each variable \n",
    "x_val = numpy.arange(10)\n",
    "y_val = numpy.arange(10)\n",
    "alpha_val = 0.1\n",
    "beta_val = 0.3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* compute $z = (x_1 + y_1^2, x_2 + y_2^2, ...)^T$: define theano function and evaluate it"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* compute $||x||$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* compute $(<x, \\alpha y> + <\\beta x, y>)^2$, where $ <a , b >  $ is a scalar product"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Names for variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# define \n",
    "x_var = T.vector(name='x')\n",
    "y_var = T.vector(name='y')\n",
    "alpha_var = T.scalar(name='a')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "z_var = alpha_var * x_var * T.log(y_var) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### `theano.printing.pprint()`\n",
    "Theano provides the functions `theano.printing.pprint()` and `theano.printing.debugprint()` to print a graph to the terminal before or after compilation. `pprint()` is more compact and math-like, `debugprint()` is more verbose. Theano also provides `pydotprint()` that creates an image of the function. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "theano.pprint(x_var)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "theano.pprint(z_var)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "compiled_expr = theano.function(inputs=[x_var, y_var, alpha_var], outputs=[z_var], name='function')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "theano.printing.debugprint(compiled_expr)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "theano.printing.pydotprint(compiled_expr, outfile=\"graph.svg\", var_with_name_simple=True, format='svg')\n",
    "SVG('./graph.svg')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Gradient — why `theano` matters"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* `Theano` can compute derivatives and gradients automatically\n",
    "* Derivatives are computed symbolically, not numerically: `T.grad` computes symbolic gradients for one or more variables with respect to some variables.\n",
    "\n",
    "Limitations:\n",
    "* You can only compute a gradient of a scalar transformation over one or several scalar or vector (or tensor) transformations or inputs.\n",
    "* A transformation has to have `float32` or `float64` `dtype` throughout the whole computation graph, because derivative over an integer has no mathematical sense\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1D gradient (derivative)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "x_var = T.scalar(name='x')\n",
    "squaredx_var = x_var ** 2\n",
    "squaredx_derivative_var = T.grad(squaredx_var, x_var)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# result of analytical differentiation:\n",
    "theano.pprint(squaredx_derivative_var)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# let's compile it\n",
    "compiled_function = theano.function([x_var], squaredx_var)\n",
    "compiled_derivative = theano.function([x_var], squaredx_derivative_var)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# optimizations were done during compilation\n",
    "theano.printing.debugprint(compiled_derivative)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "x_val = numpy.linspace(-3, 3, 100)\n",
    "plt.plot(x_val, map(compiled_function, x_val), label='function')\n",
    "plt.plot(x_val, map(compiled_derivative, x_val), label='derivative')\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### NDimensional gradient"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "x_var = T.vector(name='x')\n",
    "function_var = T.sum(x_var * x_var)\n",
    "function_gradient_var = T.grad(function_var, x_var)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "compiled_gradient = theano.function([x_var], function_gradient_var)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "compiled_gradient([1, 2, 4, 0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "theano.printing.debugprint(compiled_gradient)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Why T.grad rocks"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def compute_volatility(my_scalar_var, my_vector_var):\n",
    "    wps = ((my_vector_var + my_scalar_var)**(1 + T.var(my_vector_var)) + \\\n",
    "                  1. / T.arcsinh(my_scalar_var)).mean() / (my_scalar_var**2 + 1) + \\\n",
    "                  0.01 * T.sin(2 * my_scalar_var**1.5) * (T.sum(my_vector_var) * my_scalar_var**2) \\\n",
    "                  * T.exp((my_scalar_var - 4)**2) / (1 + T.exp((my_scalar_var - 4)**2)) * \\\n",
    "                 (1 - (T.exp(-(my_scalar_var-4)**2)) / (1 + T.exp(-(my_scalar_var-4)**2)))**2\n",
    "    return wps.mean()\n",
    "\n",
    "# define varibales\n",
    "my_scalar_var = T.scalar(name='input', dtype='float64')\n",
    "my_vector_var = T.vector('float64')\n",
    "volatility_var = compute_volatility(my_scalar_var, my_vector_var)\n",
    "\n",
    "# define derivatives\n",
    "derivative_by_scalar_var = T.grad(volatility_var, my_scalar_var)\n",
    "derivative_by_vector_var = T.grad(volatility_var, my_vector_var)\n",
    "\n",
    "# compile the function and its derivatives\n",
    "compiled_fun_function = theano.function([my_scalar_var, my_vector_var], volatility_var)\n",
    "compiled_fun_derivative_scalar = theano.function([my_scalar_var, my_vector_var], derivative_by_scalar_var)\n",
    "\n",
    "# Optional exercise on calculus: compute derivatives"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# let's plot for different values of scalar\n",
    "my_scalar_vals = numpy.linspace(0, 7, 100)\n",
    "my_vector_val = [1, 2, 3]\n",
    "\n",
    "plt.plot(my_scalar_vals, [compiled_fun_function(val, my_vector_val) for val in my_scalar_vals], \n",
    "         label='function')\n",
    "plt.plot(my_scalar_vals, [compiled_fun_derivative_scalar(val, my_vector_val) for val in my_scalar_vals],\n",
    "         label='derivative')\n",
    "plt.legend(loc='best')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Shared variables\n",
    "* The inputs and transformations only exist when function is called\n",
    "* Shared variables always **stay in memory** like global variables\n",
    "    * Shared variables are shared between functions where they appear in.\n",
    "    * Shared variables can be included into a symbolic graph\n",
    "    * They can be set and evaluated using special methods\n",
    "        * but they can't change value arbitrarily during symbolic graph computation\n",
    "   \n",
    "* Hint: such variables are a perfect place to store network parameters\n",
    "    * e.g. weights or some metadata"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "w_shared = theano.shared(numpy.arange(10, dtype=float), name='weight')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "print w_shared, w_shared.get_value()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "x_var = T.vector('x')\n",
    "# dotproduct = theano.function([x_var], T.sum(x_var * w_shared))\n",
    "dotproduct_grad = theano.function([x_var], T.grad(T.sum(x_var * w_shared), x_var))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "dotproduct_grad(numpy.arange(10) * 0.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "dotproduct_grad(numpy.arange(10) * 10.)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Matrix"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Compute $||Ax||^2$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "x_var = T.vector('x')\n",
    "A_var = T.matrix('A')\n",
    "z_var = A_var.dot(x_var)\n",
    "normAx = theano.function([x_var, A_var], z_var.dot(z_var))\n",
    "normAx([0, 2], [[1, 1], [1, 1]])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercises"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* compile a function which takes an array $x$ with three elements and computes $x_0^3 + \\sin{x_1}*\\cos{x_2}$. Check it by evaluating at some points"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* compute gradient for previous function w.r.t. x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* compute gradient for $||Ax|| + \\alpha * ||x||$ w.r.t. x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Logistic regression "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# import function to create toy dataset for classification\n",
    "from sklearn.datasets import make_blobs, make_moons\n",
    "from sklearn.cross_validation import train_test_split\n",
    "from sklearn.metrics import roc_auc_score\n",
    "# logistic function (we introduce shortcut sigmoid)\n",
    "from scipy.special import expit as sigmoid"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### generate toy samples"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "n_features = 10 # number of features\n",
    "centers = 2 # number of classes\n",
    "X, y = make_blobs(n_samples=10000, centers=centers, n_features=n_features, random_state=42, cluster_std=10)\n",
    "trainX, testX, trainY, testY = train_test_split(X, y, train_size=0.5, random_state=42)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Logistic regression description:\n",
    "\n",
    "Implement the regular logistic regression training algorithm\n",
    "\n",
    "Tips:\n",
    "* Weights are represented as a `theano` vector\n",
    "* trainX and trainY are constants, stored in `numpy.arrays`\n",
    "\n",
    "Exercise:\n",
    "1. compute probabilities (**use** `T.nnet.sigmoid`):\n",
    "    $$p_i = \\sigma(\\sum_k X_{ik} w_k)$$\n",
    "2. compute logistic loss function\n",
    "    $$\\mathcal{L}=-\\sum_i y_i \\log{p_i} + (1-y_i)\\log{(1 - p_i)}\\qquad,$$ where $y \\in \\{0, 1\\}$\n",
    "\n",
    "3. compute loss function gradient using `theano` "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# write expression for probabilities\n",
    "# w_var = \n",
    "# p_var = ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# write expression for loss\n",
    "# loss_var = ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# compile loss expression, compile gradient expression for loss\n",
    "# loss_function = theano.function(...)\n",
    "# loss_grad_function = theano.function(...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# very simple check\n",
    "loss_function(numpy.random.random(n_features))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we have expessions for loss and its gradient and we need to use some optimization method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from scipy.optimize import minimize"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# minimize loss function using its gradient\n",
    "result = minimize(fun=loss_function, jac=loss_grad_function, x0=numpy.zeros(n_features))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "result"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "w_optimal = result['x']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Now predict output of logistic regression for the test sample and compute AUC"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# compute probabilities using numpy\n",
    "from scipy.special import expit\n",
    "proba = expit(testX.dot(w_optimal))\n",
    "print roc_auc_score(testY, proba)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Last element: updates\n",
    "* updates are a way of changing shared variables at after function call.\n",
    "* technically it's a dictionary {shared_variable : a recipe for new value} which is has to be provided when function is compiled\n",
    "\n",
    "That's how it works:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Multiply shared vector by a number and save the product back into shared vector\n",
    "\n",
    "ones_shared_vector = theano.shared(numpy.ones(10, dtype='float64'))\n",
    "input_scalar_var = T.scalar('coefficient', dtype='float32')\n",
    "expression_var = input_scalar_var * ones_shared_vector\n",
    "\n",
    "inputs = [input_scalar_var]\n",
    "outputs = [expression_var]  # return vector times scalar\n",
    "\n",
    "my_updates = {\n",
    "    ones_shared_vector: expression_var  # and write the result into ones_shared_vector\n",
    "}\n",
    "\n",
    "# updates appeared\n",
    "compute_and_save = theano.function(inputs, outputs, updates=my_updates)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# initial ones_shared_vector\n",
    "print \"initial shared value:\", ones_shared_vector.get_value()\n",
    "\n",
    "# evaluating the function (ones_shared_vector will be changed)\n",
    "print \"compute_and_save(2) returns\", compute_and_save(2)\n",
    "\n",
    "# evaluate new ones_shared_vector\n",
    "print \"new shared value:\", ones_shared_vector.get_value()\n",
    "\n",
    "# evaluating the function (ones_shared_vector will be changed)\n",
    "print \"compute_and_save(2) returns\", compute_and_save(2)\n",
    "\n",
    "# evaluate new ones_shared_vector\n",
    "print \"new shared value:\", ones_shared_vector.get_value()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Logistic regression with updates in `theano`:\n",
    "\n",
    "We did it the simplest way, now let's do it *the right way*!\n",
    "\n",
    "Tips:\n",
    "* Weights are represented as a shared variable\n",
    "* X and y are inputs of a function\n",
    "\n",
    "Compile 2 functions:\n",
    "* train_function(X, y) — returns an error and computes new values of weights (through updates)\n",
    "* predict_function(X) — just computes probabilities (\"y\") given data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# inputs and shareds\n",
    "# w_shared = ...\n",
    "# X_var = ...\n",
    "# y_var = ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# proba_var = ...\n",
    "# loss_var = ...\n",
    "# grad_var = ...\n",
    "\n",
    "# updates = ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# train_function = theano.function(...)\n",
    "# predict_function = theano.function(...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# training process\n",
    "for i in range(15):\n",
    "    print \"loss at iter %i:%.4f\" % (i, train_function(trainX, trainY)), \n",
    "    print \"train auc:\", roc_auc_score(trainY, predict_function(trainX)),\n",
    "    print \"test auc:\", roc_auc_score(testY, predict_function(testX))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Neural Network with theano\n",
    "\n",
    "## NN with one hidden layer\n",
    "This is a simple NN description with one hidden layer:\n",
    "\n",
    "Parameters: \n",
    "\n",
    "* $W$, $v$\n",
    "\n",
    "Calculations:\n",
    "\n",
    "* hidden activations: h = $\\sigma$(X.dot(W))\n",
    "* output = h.dot(v)\n",
    "* $p_{1}$ = $\\sigma$(output)\n",
    "* $p_{0} = 1 - p_{1}$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise: \n",
    "\n",
    "Write 1-hidden layer NN using theano"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# number of neurons in the hidden layer\n",
    "n_hidden = 10\n",
    "\n",
    "# define input variables\n",
    "# TODO\n",
    "\n",
    "# define NNs parameters\n",
    "# TODO"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# define train_function and predict_function\n",
    "# TODO"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# training process\n",
    "for i in range(15):\n",
    "    print \"loss at iter %i:%.4f\" % (i, train_function(trainX, trainY)), \n",
    "    print \"train auc:\", roc_auc_score(trainY, predict_function(trainX)),\n",
    "    print \"test auc:\", roc_auc_score(testY, predict_function(testX))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## NN with two hidden layers\n",
    "first let's take another dataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X, y = make_moons(n_samples=20000, noise=0.1)\n",
    "trainX, testX, trainY, testY = train_test_split(X, y, train_size=0.5, random_state=42)\n",
    "n_features = X.shape[1]\n",
    "plt.scatter(X[:, 0], X[:, 1], c=y, alpha=0.05, linewidths=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise\n",
    "\n",
    "Write 2 hidden layers NN (use code from 1-hidden layer NN)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# define train_function and predict_function\n",
    "# TODO"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# training process\n",
    "for i in range(100):\n",
    "    print \"loss at iter %i:%.4f\" % (i, train_function(trainX, trainY)), \n",
    "    print \"train auc:\", roc_auc_score(trainY, predict_function(trainX)),\n",
    "    print \"test auc:\", roc_auc_score(testY, predict_function(testX))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## So, what's good about theano?\n",
    "\n",
    "- symbolic computations with [rich set of operations](http://deeplearning.net/software/theano/library/tensor/index.html)\n",
    "- very nice correspondence with numpy\n",
    "- theano operates with symbolic graph to [optimize](http://deeplearning.net/software/theano/optimizations.html) for speed and stability\n",
    "- code we wrote can be evaluated on GPU \n",
    "- loops can be part of a function! Can you differentiate through loops? Theano [can](http://deeplearning.net/software/theano/tutorial/loop.html)\n",
    "- [optimization for gradients computation](http://deeplearning.net/software/theano/tutorial/gradients.html), like L-operators\n",
    "- and many other things"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Define special function for NNs\n",
    "\n",
    "This function simplifies experiments with activations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def define_NN(layers, activations, learning_rate=0.1):\n",
    "    X_var = T.matrix('float64')\n",
    "    y_var = T.vector('float64')\n",
    "    result = X_var\n",
    "    \n",
    "    shared_variables = []\n",
    "    for dim in range(1, len(layers)):\n",
    "        w_shared = theano.shared(numpy.random.random(size=(layers[dim - 1], layers[dim])) / 10.)\n",
    "        result = activations[dim - 1](result.dot(w_shared))\n",
    "        shared_variables.append(w_shared)\n",
    "        \n",
    "    v_shared = theano.shared(numpy.random.random(size=layers[-1]) / 10.)\n",
    "    shared_variables.append(v_shared)\n",
    "    \n",
    "    proba_var = T.nnet.sigmoid(result.dot(v_shared))\n",
    "    loss_var = -T.mean(y_var * T.log(proba_var) + (1 - y_var) * T.log(1 - proba_var))\n",
    "        \n",
    "    updates = {shared: shared - learning_rate * T.grad(loss_var, shared) for shared in shared_variables}\n",
    "\n",
    "    train_function = theano.function([X_var, y_var], loss_var, updates=updates)\n",
    "    predict_function = theano.function([X_var], proba_var)\n",
    "    return train_function, predict_function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "train_function, predict_function = define_NN((n_features, 20, 20), [T.nnet.sigmoid] * 2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# training process\n",
    "for i in range(100):\n",
    "    print \"loss at iter %i:%.4f\" % (i, train_function(trainX, trainY)), \n",
    "    print \"train auc:\", roc_auc_score(trainY, predict_function(trainX)),\n",
    "    print \"test auc:\", roc_auc_score(testY, predict_function(testX))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Popular activation functions\n",
    "\n",
    "* **Sigmoid**:\n",
    "\n",
    "    $f(x) = \\frac{1}  {1+e^{-x}}$\n",
    "\n",
    "\n",
    "* **ReLU - rectifier linear unit**\n",
    "\n",
    "    In the context of artificial neural networks, the rectifier is an activation function defined as\n",
    "\n",
    "    $f(x) = \\max(0, x)$\n",
    "\n",
    "    where x is the input to a neuron. This activation function has been argued to be more biologically plausible than the widely used logistic sigmoid (which is inspired by probability theory; see logistic regression) and its more practical counterpart, the hyperbolic tangent. The rectifier is the most popular activation function for deep neural networks.\n",
    "\n",
    "    A unit employing the rectifier is also called a rectified linear unit (ReLU).\n",
    "\n",
    "\n",
    "* **Softplus**\n",
    "    A smooth approximation to the rectifier is the analytic function\n",
    "\n",
    "    $f(x) = \\ln(1 + e^x)$\n",
    "\n",
    "    which is called softplus function."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Let's play with activations: \n",
    "\n",
    "Exercise: compare different intermediate activation functions:\n",
    "\n",
    "* sigmoid (which we used, `T.nnet.sigmoid`)\n",
    "* leaky ReLU (defined below)\n",
    "* softplus (`T.nnet.softplus`)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def LeakyReLU(x):\n",
    "    return T.switch(x > 0, x, 0.5 * x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Use above define_NN function to play with activations\n",
    "# TODO"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Feeling the power?\n",
    "\n",
    "Add multiclassification to the neural network!\n",
    "\n",
    "You'll need `softmax` activation function for the last layer to compute probabilities."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## NNs available in `hep_ml`\n",
    "library has some examples of neural networks"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X, y = make_moons(n_samples=2000, noise=0.25, random_state=42)\n",
    "trainX, testX, trainY, testY = train_test_split(X, y, train_size=0.5, random_state=42)\n",
    "plt.scatter(X[:, 0], X[:, 1], c=y, alpha=0.5, linewidths=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from hep_ml.nnet import PairwiseNeuralNetwork, RBFNeuralNetwork, MLPClassifier"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "for models in [PairwiseNeuralNetwork(random_state=42), \n",
    "               RBFNeuralNetwork(random_state=42), \n",
    "               MLPClassifier(random_state=42)]:\n",
    "    models.fit(trainX, trainY)\n",
    "    pred = models.predict_proba(testX)[:, 1]\n",
    "    print roc_auc_score(testY, pred)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# just for comparison ...\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "from sklearn.ensemble import RandomForestClassifier, GradientBoostingClassifier\n",
    "from sklearn.neighbors import KNeighborsClassifier\n",
    "\n",
    "knn_clf = KNeighborsClassifier(n_neighbors=15)\n",
    "knn_clf.fit(trainX, trainY)\n",
    "print 'KNN', roc_auc_score(testY, knn_clf.predict_proba(testX)[:, 1])\n",
    "\n",
    "log_clf = LogisticRegression(C=100)\n",
    "log_clf.fit(trainX, trainY)\n",
    "print 'logistic regression', roc_auc_score(testY, log_clf.predict_proba(testX)[:, 1])\n",
    "\n",
    "rf_clf = RandomForestClassifier(n_estimators=100, max_depth=20)\n",
    "rf_clf.fit(trainX, trainY)\n",
    "print 'random forest', roc_auc_score(testY, rf_clf.predict_proba(testX)[:, 1])\n",
    "\n",
    "gb_clf = GradientBoostingClassifier(n_estimators=200, max_depth=5, learning_rate=0.1)\n",
    "gb_clf.fit(trainX, trainY)\n",
    "print 'GB', roc_auc_score(testY, gb_clf.predict_proba(testX)[:, 1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Bagging over NN — using sklearn's meta algorithms"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from sklearn.ensemble import AdaBoostClassifier, BaggingClassifier"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "base = RBFNeuralNetwork(layers=(40,))\n",
    "meta_ada = AdaBoostClassifier(base_estimator=base, n_estimators=10, learning_rate=0.05, random_state=42)\n",
    "meta_ada.fit(trainX, trainY)\n",
    "\n",
    "print roc_auc_score(testY, meta_ada.predict_proba(testX)[:, 1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "base = RBFNeuralNetwork(layers=(40,))\n",
    "meta_bagging = BaggingClassifier(base_estimator=base, n_estimators=10, random_state=42)\n",
    "meta_bagging.fit(trainX, trainY)\n",
    "\n",
    "print roc_auc_score(testY, meta_bagging.predict_proba(testX)[:, 1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "** Nice, this works for simple datasets! Maybe try such tricks for some other data? **"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# References\n",
    "- [theano documentation](http://deeplearning.net/software/theano/)"
   ]
  }
 ],
 "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
}