Parameter exploration - psyrun.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Parameter exploration\n",
    "\n",
    "Purpose: Run the simulation with varying parameters and characterize the effects of those parameters\n",
    "    \n",
    "Parameter exploration can either be done as grid search where a multidimensional \"regular\" grid of parameter values is explored or as random search where random parameters values are picked.\n",
    "\n",
    "# psyrun tutorial"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from __future__ import print_function\n",
    "from pprint import pprint"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This tutorial will introduce you to the [\"psyrun\"](https://github.com/jgosmann/psyrun) tool for parameter space exploration and serial farming step by step. It also integrates well with [\"ctn_benchmarks\"](https://github.com/ctn-waterloo/ctn_benchmarks). We will use the `CommunicationChannel` as an example."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from ctn_benchmark.nengo.communication import CommunicationChannel"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Parameter space exploration\n",
    "\n",
    "Let us try to find out how the RMSE of the communication channel varies with the number of dimensions and number of neurons. Usually you would write a bunch of nested `for` loops like this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "running CommunicationChannel#20160328-152630-11e3af7f\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 2\n",
      "_L = 2\n",
      "_N = 10\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.16572583026875315\n",
      "sim_speed = 7.043639184918233\n",
      "running CommunicationChannel#20160328-152630-4d65aacb\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 2\n",
      "_L = 2\n",
      "_N = 50\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.070963177680149134\n",
      "sim_speed = 8.502611008400637\n",
      "running CommunicationChannel#20160328-152630-4d65aacb\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 2\n",
      "_L = 2\n",
      "_N = 100\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.028531004425303052\n",
      "sim_speed = 8.638787279617729\n",
      "running CommunicationChannel#20160328-152630-4d65aacb\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 3\n",
      "_L = 2\n",
      "_N = 10\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.012725007735396695\n",
      "sim_speed = 9.304754087450362\n",
      "running CommunicationChannel#20160328-152630-09efe4b2\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 3\n",
      "_L = 2\n",
      "_N = 50\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.017750297895821786\n",
      "sim_speed = 8.73140287694902\n",
      "running CommunicationChannel#20160328-152630-09efe4b2\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 3\n",
      "_L = 2\n",
      "_N = 100\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.027605168028911279\n",
      "sim_speed = 8.48161741659033\n",
      "running CommunicationChannel#20160328-152631-09efe4b2\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 4\n",
      "_L = 2\n",
      "_N = 10\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.16662160881813956\n",
      "sim_speed = 10.729041004783465\n",
      "running CommunicationChannel#20160328-152631-09efe4b2\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 4\n",
      "_L = 2\n",
      "_N = 50\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.087142299358396944\n",
      "sim_speed = 8.984726643331777\n",
      "running CommunicationChannel#20160328-152631-09efe4b2\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 4\n",
      "_L = 2\n",
      "_N = 100\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.023944215909647099\n",
      "sim_speed = 8.744710084042374\n"
     ]
    }
   ],
   "source": [
    "rmses = []\n",
    "for D in np.arange(2, 5):\n",
    "    for N in [10, 50, 100]:\n",
    "        rmses.append(CommunicationChannel().run(D=D, N=N)['rmse'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[0.16572583026875315,\n",
      " 0.070963177680149134,\n",
      " 0.028531004425303052,\n",
      " 0.012725007735396695,\n",
      " 0.017750297895821786,\n",
      " 0.027605168028911279,\n",
      " 0.16662160881813956,\n",
      " 0.087142299358396944,\n",
      " 0.023944215909647099]\n"
     ]
    }
   ],
   "source": [
    "pprint(rmses)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "However there are a few annoyances with this approach: We have to add another `for` loop for each dimensions in the parameter space; and we have to handle storing the results ourselves. Here, especially, the assignment to input parameters in only implicit.\n",
    "\n",
    "With psyrun, we can define a parameter space in a more natural way:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Param(D=[2, 2, 2, 3, 3, 3, 4, 4, 4], N=[10, 50, 100, 10, 50, 100, 10, 50, 100])\n"
     ]
    }
   ],
   "source": [
    "from psyrun import Param\n",
    "\n",
    "pspace = Param(D=np.arange(2, 5)) * Param(N=[10, 50, 100])\n",
    "print(pspace)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We use the `Param` class to define sets of parameters with keyword arguments to the constructor. As long as all sequences have the same length it is possible to give multiple arguments to the constructor:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Param(a=[1, 2], b=[3, 4])\n"
     ]
    }
   ],
   "source": [
    "print(Param(a=[1, 2], b=[3, 4]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If at least one argument is a sequence, non-sequence arguments will generate lists of the required length:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Param(a=[1, 2], b=[3, 3])\n"
     ]
    }
   ],
   "source": [
    "print(Param(a=[1, 2], b=3))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Such parameters definitions support basic operations like a cartesian product (as shown above), concatenation with `+`, and two subtraction like operations with `-` and the `psyrun.pspace.missing` function.\n",
    "\n",
    "Once we defined a parameter space we want to explore, it is easy to apply a function (i.e. run our model) to each set of parameters."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "running CommunicationChannel#20160328-152631-09efe4b2\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 2\n",
      "_L = 2\n",
      "_N = 10\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.16572583026875315\n",
      "sim_speed = 9.76084634413682\n",
      "running CommunicationChannel#20160328-152631-4d65aacb\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 2\n",
      "_L = 2\n",
      "_N = 50\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.070963177680149134\n",
      "sim_speed = 9.590115169072405\n",
      "running CommunicationChannel#20160328-152631-4d65aacb\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 2\n",
      "_L = 2\n",
      "_N = 100\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.028531004425303052\n",
      "sim_speed = 9.921217898528957\n",
      "running CommunicationChannel#20160328-152632-4d65aacb\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 3\n",
      "_L = 2\n",
      "_N = 10\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.012725007735396695\n",
      "sim_speed = 10.659266205662675\n",
      "running CommunicationChannel#20160328-152632-09efe4b2\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 3\n",
      "_L = 2\n",
      "_N = 50\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.017750297895821786\n",
      "sim_speed = 9.910996533529302\n",
      "running CommunicationChannel#20160328-152632-09efe4b2\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 3\n",
      "_L = 2\n",
      "_N = 100\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.027605168028911279\n",
      "sim_speed = 8.251111478325189\n",
      "running CommunicationChannel#20160328-152632-09efe4b2\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 4\n",
      "_L = 2\n",
      "_N = 10\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.16662160881813956\n",
      "sim_speed = 10.275069695885861\n",
      "running CommunicationChannel#20160328-152632-09efe4b2\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 4\n",
      "_L = 2\n",
      "_N = 50\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.087142299358396944\n",
      "sim_speed = 8.868312774603872\n",
      "running CommunicationChannel#20160328-152632-09efe4b2\n",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 4\n",
      "_L = 2\n",
      "_N = 100\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.023944215909647099\n",
      "sim_speed = 7.7395821223813455\n"
     ]
    }
   ],
   "source": [
    "import psyrun\n",
    "\n",
    "result = psyrun.map_pspace(CommunicationChannel().run, pspace)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'D': [2, 2, 2, 3, 3, 3, 4, 4, 4],\n",
      " 'N': [10, 50, 100, 10, 50, 100, 10, 50, 100],\n",
      " 'rmse': [0.16572583026875315,\n",
      "          0.070963177680149134,\n",
      "          0.028531004425303052,\n",
      "          0.012725007735396695,\n",
      "          0.017750297895821786,\n",
      "          0.027605168028911279,\n",
      "          0.16662160881813956,\n",
      "          0.087142299358396944,\n",
      "          0.023944215909647099],\n",
      " 'sim_speed': [9.76084634413682,\n",
      "               9.590115169072405,\n",
      "               9.921217898528957,\n",
      "               10.659266205662675,\n",
      "               9.910996533529302,\n",
      "               8.251111478325189,\n",
      "               10.275069695885861,\n",
      "               8.868312774603872,\n",
      "               7.7395821223813455]}\n"
     ]
    }
   ],
   "source": [
    "pprint(result)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `map_pspace` function conveniently returns a dictionary with all the input parameters and output values. This makes it easy to convert it to other data structurs like a Pandas data frame for further analysis:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>D</th>\n",
       "      <th>N</th>\n",
       "      <th>rmse</th>\n",
       "      <th>sim_speed</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2</td>\n",
       "      <td>10</td>\n",
       "      <td>0.165726</td>\n",
       "      <td>9.760846</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2</td>\n",
       "      <td>50</td>\n",
       "      <td>0.070963</td>\n",
       "      <td>9.590115</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2</td>\n",
       "      <td>100</td>\n",
       "      <td>0.028531</td>\n",
       "      <td>9.921218</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>3</td>\n",
       "      <td>10</td>\n",
       "      <td>0.012725</td>\n",
       "      <td>10.659266</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>3</td>\n",
       "      <td>50</td>\n",
       "      <td>0.017750</td>\n",
       "      <td>9.910997</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>3</td>\n",
       "      <td>100</td>\n",
       "      <td>0.027605</td>\n",
       "      <td>8.251111</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>4</td>\n",
       "      <td>10</td>\n",
       "      <td>0.166622</td>\n",
       "      <td>10.275070</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>4</td>\n",
       "      <td>50</td>\n",
       "      <td>0.087142</td>\n",
       "      <td>8.868313</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>4</td>\n",
       "      <td>100</td>\n",
       "      <td>0.023944</td>\n",
       "      <td>7.739582</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   D    N      rmse  sim_speed\n",
       "0  2   10  0.165726   9.760846\n",
       "1  2   50  0.070963   9.590115\n",
       "2  2  100  0.028531   9.921218\n",
       "3  3   10  0.012725  10.659266\n",
       "4  3   50  0.017750   9.910997\n",
       "5  3  100  0.027605   8.251111\n",
       "6  4   10  0.166622  10.275070\n",
       "7  4   50  0.087142   8.868313\n",
       "8  4  100  0.023944   7.739582"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import pandas as pd\n",
    "pd.DataFrame(result)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Parallelizaton\n",
    "\n",
    "It would be nice to parallelize all these simulations to make the best use of multicore processors (though, your BLAS library might be parallelized already). This is easy with psyrun as we just have to use a different map function. The only problem is that the function we want to apply needs to be pickleable which requires us to put it into an importable Python module. (It is also possible to switch the parallelization backend from `multiprocessing` to `threading`, but that can introduce other issues.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%%sh\n",
    "cat << EOF > model.py\n",
    "from ctn_benchmark.nengo.communication import CommunicationChannel\n",
    "def evaluate(**kwargs):\n",
    "    return CommunicationChannel().run(**kwargs)\n",
    "EOF"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "running CommunicationChannel#20160328-152633-09efe4b2\n",
      "running CommunicationChannel#20160328-152633-09efe4b2\n",
      "running CommunicationChannel#20160328-152633-09efe4b2\n",
      "running CommunicationChannel#20160328-152633-09efe4b2\n",
      "running CommunicationChannel#20160328-152633-09efe4b2\n",
      "running CommunicationChannel#20160328-152633-09efe4b2\n",
      "running CommunicationChannel#20160328-152633-09efe4b2\n",
      "running CommunicationChannel#20160328-152633-09efe4b2\n",
      "running CommunicationChannel#20160328-152633-09efe4b2\n",
      "\r",
      "Simulation finished in 0:00:01.                                                 \n",
      "\r",
      "Simulation finished in 0:00:01.                                                 \n",
      "\r",
      "Simulation finished in 0:00:01.                                                 \n",
      "\r",
      "Simulation finished in 0:00:01.                                                 \n",
      "\r",
      "Simulation finished in 0:00:01.                                                 \n",
      "\r",
      "Simulation finished in 0:00:01.                                                 \n",
      "\r",
      "Simulation finished in 0:00:01.                                                 \n",
      "\r",
      "Simulation finished in 0:00:01.                                                 \n",
      "\r",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 2\n",
      "_L = 2\n",
      "_N = 50\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.070963177680149134\n",
      "sim_speed = 4.957144022134159_D = 2\n",
      "_L = 2\n",
      "_N = 10\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.16572583026875315\n",
      "sim_speed = 4.187779439633052_D = 3\n",
      "_L = 2\n",
      "_N = 50\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.017750297895821786\n",
      "sim_speed = 6.460906186554772_D = 2\n",
      "_L = 2\n",
      "_N = 100\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.028531004425303052\n",
      "sim_speed = 4.271733440339595_D = 3\n",
      "_L = 2\n",
      "_N = 10\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.012725007735396695\n",
      "sim_speed = 4.16351067504735_D = 4\n",
      "_L = 2\n",
      "_N = 50\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.087142299358396944\n",
      "sim_speed = 5.09186174545693_D = 4\n",
      "_L = 2\n",
      "_N = 100\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.023944215909647099\n",
      "sim_speed = 4.6014242098114915_D = 3\n",
      "_L = 2\n",
      "_N = 100\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.027605168028911279\n",
      "sim_speed = 4.232498062516146_D = 4\n",
      "_L = 2\n",
      "_N = 10\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.16662160881813956\n",
      "sim_speed = 4.321293862744472\n",
      "\n",
      "\n",
      "\n",
      "\n",
      "\n",
      "\n",
      "\n",
      "\n"
     ]
    }
   ],
   "source": [
    "from model import evaluate\n",
    "result2 = psyrun.map_pspace_parallel(evaluate, pspace)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'D': [2, 2, 2, 3, 3, 3, 4, 4, 4],\n",
      " 'N': [10, 50, 100, 10, 50, 100, 10, 50, 100],\n",
      " 'rmse': [0.16572583026875315,\n",
      "          0.070963177680149134,\n",
      "          0.028531004425303052,\n",
      "          0.012725007735396695,\n",
      "          0.017750297895821786,\n",
      "          0.027605168028911279,\n",
      "          0.16662160881813956,\n",
      "          0.087142299358396944,\n",
      "          0.023944215909647099],\n",
      " 'sim_speed': [4.187779439633052,\n",
      "               4.957144022134159,\n",
      "               4.271733440339595,\n",
      "               4.16351067504735,\n",
      "               6.460906186554772,\n",
      "               4.232498062516146,\n",
      "               4.321293862744472,\n",
      "               5.09186174545693,\n",
      "               4.6014242098114915]}\n"
     ]
    }
   ],
   "source": [
    "pprint(result2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Serial farming"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "After doing some initial experiments, you might want to run a larger number of simulations. Because all the simulations are independent from one another, they could run one after another on your computer, but that would take a very long time. It would be faster to distribute the individual simulations to a high performance cluster (like Sharcnet) with many CPU cores that can run the simulations in parallel. This is called serial farming.\n",
    "\n",
    "To do serial farming, first make sure that it is a Python module that can be installed, is installed, and can be imported. (It is possible to change `sys.path` if you don't want to install it as a Python module, but this approach is more complicated and error prone.) Then you need to define a \"task\" for what you intend to run. To define a task you have to create a file `task_<name>.py` in a directory called `psy-tasks`. Here is an example task for our communication channel:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# task_cc1.py\n",
    "\n",
    "from ctn_benchmark.nengo.communication import CommunicationChannel\n",
    "import numpy as np\n",
    "import psyrun\n",
    "from psyrun import Param\n",
    "\n",
    "pspace = Param(D=np.arange(2, 5)) * Param(N=[10, 50, 100])\n",
    "min_items = 1\n",
    "max_jobs = 6\n",
    "\n",
    "def execute(**kwargs):\n",
    "    return CommunicationChannel().run(**kwargs)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There are a number of variables with a special meaning in these kind of task files. `pspace` defines the parameter space to explore which should be familar by now. `execute` is the function that gets called for each assignment of parameters and returns a dictionary with the results. `max_jobs` says that a maximum of 6 processing jobs should be created (a few other jobs not doing the core processing might be created in addition).\n",
    "\n",
    "When running this task (we will get to that in a moment), the parameter space will be split into parts of the same size and assigned to the processing jobs to be processed. At the end all the results from the processing jobs will be merged into a single file. `min_items` defines how many parameter assignments should at least be processd by each job.\n",
    "\n",
    "To run tasks defined in this way, we use the `psy-doit` command. As the name suggests, this tool is based on the excellent [doit](http://pydoit.org/) automation tool. Be sure that your working directory contains the `psy-tasks` directory when invoking this command.\n",
    "\n",
    "First we can list the available tasks:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "cc1   \n",
      "cc2   \n"
     ]
    }
   ],
   "source": [
    "%%sh\n",
    "psy-doit list"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us test our task first:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "cc1\n",
      "running CommunicationChannel#20160328-152634-1b7add91\n",
      "\r",
      "Simulation finished in 0:00:01.                                                 \n",
      "_D = 2\n",
      "_L = 2\n",
      "_N = 10\n",
      "_pstc = 0.01\n",
      "_T = 1.0\n",
      "_backend = 'nengo'\n",
      "_dt = 0.001\n",
      "_seed = 1\n",
      "_hide_overlay = False\n",
      "_gui = False\n",
      "rmse = 0.16572583026875315\n",
      "sim_speed = 10.075970288371915\n"
     ]
    }
   ],
   "source": [
    "%%sh\n",
    "psy-doit test cc1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This immediatly runs the first (and only the first) parameter assignment as a way to check that nothing crashes. Now let us run the full task."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      ".  cc1:split\n",
      ".  cc1:process:0\n",
      ".  cc1:process:1\n",
      ".  cc1:process:2\n",
      ".  cc1:process:3\n",
      ".  cc1:process:4\n",
      ".  cc1:process:5\n",
      ".  cc1:merge\n"
     ]
    }
   ],
   "source": [
    "%%sh\n",
    "psy-doit cc1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This produces standard doit output with dots (`.`) indicating that a task was executed. As you can see the first task is splitting the parameter space and then several processing jobs are started and in the end the results are merged.\n",
    "\n",
    "All intermediary files for this will be written to the `psy-work/<task name>` directory by default. This is also where we find the result file. psyrun provides some helper functions to load this file:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'D': array([2, 2, 2, 3, 4, 3, 3, 4, 4]),\n",
      " 'N': array([ 10, 100,  50, 100,  10,  10,  50,  50, 100]),\n",
      " 'rmse': array([ 0.16572583,  0.028531  ,  0.07096318,  0.02760517,  0.16662161,\n",
      "        0.01272501,  0.0177503 ,  0.0871423 ,  0.02394422]),\n",
      " 'sim_speed': array([  9.86115448,   9.05615532,   9.2535156 ,   8.22138503,\n",
      "        10.00549618,  10.20564068,   9.70149143,   9.1584289 ,   9.24154737])}\n"
     ]
    }
   ],
   "source": [
    "data = psyrun.NpzStore().load('psy-work/cc1/result.npz')\n",
    "pprint(data)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Per default NumPy `.npz` files are used to store this data. However, this format has a few disadvantages. For example, to append to such a file (which happens in the merge step), the whole file has to be loaded into memory. Thus, it is possible to switch to HDF5 which is better in this regard. Just add the line `store = psyrun.H5Store()` to your task file.\n",
    "\n",
    "Note that this will require [pytables](http://www.pytables.org/) to be installed and that the Sharcnet pytables seems to be broken at the moment.\n",
    "\n",
    "Speaking of Sharcnet: Currently psy-doit runs the jobs sequantially and immediatly. To actually get serial farming we have to tell psyrun to invoke the right job scheduler. For Sharcnet this would be `sqsub`. Here is the updated task file:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# task_cc2.py\n",
    "\n",
    "import platform\n",
    "\n",
    "from ctn_benchmark.nengo.communication import CommunicationChannel\n",
    "import numpy as np\n",
    "from psyrun import Param, Sqsub\n",
    "\n",
    "pspace = Param(D=np.arange(2, 5)) * Param(N=[10, 50, 100])\n",
    "min_items = 1\n",
    "max_jobs = 6\n",
    "\n",
    "sharcnet_nodes = ['narwhal', 'bul', 'kraken', 'saw']\n",
    "if any(platform.node().startswith(x) for x in sharcnet_nodes):\n",
    "    workdir = '/work/jgosmann/tcm'\n",
    "    scheduler = Sqsub(workdir)\n",
    "    scheduler_args = {\n",
    "        'timelimit': '10m',\n",
    "        'memory': '512M'\n",
    "    }\n",
    "\n",
    "def execute(**kwargs):\n",
    "    return CommunicationChannel().run(**kwargs)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First of all, I check whether I am actually on a Sharcnet node before I do the Sharcnet specific settings. This is not required, but allows me to run the same file on my local computer. (Note that the list of Sharcnet nodes contains only the ones I typically use. You might want to use a different list).\n",
    "\n",
    "`workdir` is a special variable to set the location for intermediary and result files. I'm changing the location because on Sharcnet `/work` offers more storage space and a better performance. The essential line to change to the sqsub scheduler is `scheduler = Sqsub(workdir)` and then we finally add some scheduler specific arguments. For sqsub you have to specify a runtime and memory usage limit.\n",
    "\n",
    "Now you can run `psy-doit cc2` on Sharcnet and it will submit a number of jobs and return. Use `sqjobs` to check the status of your jobs and see when they are finished.\n",
    "\n",
    "Sometimes jobs crash. If that the case run `psy-doit cc2` again and it should figure out by itself what jobs have to be resubmitted."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Load balancing"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sometimes the parameters can significantly influence the runtime of your models. That can be problematic on Sharcnet because if you estimate the time limit too low your jobs will get killed before completion, but a high time limit will push your job back in the queue. In those cases it might be good to do some simple load balancing.\n",
    "\n",
    "If you insert `backend = psyrun.LoadBalancingBackend`, psyrun will submit `max_jobs` jobs and each job will process as many parameters as it can, always grabbing the next in queue, until it gets killed. Unfortunately, the parameters assignment of a killed job will not be reentered into the queue. Thus, after all jobs have completed, you have to check the results file for completeness. If results are missing reexecute `psy-doit <task name>`. This will only process parameters assignments for missing results.\n",
    "\n",
    "It is best to use the load balancing backend with the `H5Store`, because it does a lot more appending to files which is inefficient with the `NpzStore`. Also note, that the load balancing backend is less tested as it is newer."
   ]
  },
  {
   "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
}