examples_from_paper.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Look Elsewhere Effect in 2-d\n",
    "\n",
    "Kyle Cranmer, Nov 19, 2015\n",
    "\n",
    "Based on\n",
    "*Estimating the significance of a signal in a multi-dimensional search* by  Ofer Vitells and Eilam Gross http://arxiv.org/pdf/1105.4355v1.pdf\n",
    "\n",
    "This is for the special case of a likelihood function of the form \n",
    "$L(\\mu, \\nu_1, \\nu_2)$ where $\\mu$ is a single parameter of interest and\n",
    "$\\nu_1,\\nu_2$ are two nuisance parameters that are not identified under the null.\n",
    "For example, $\\mu$ is the signal strength of a new particle and $\\nu_1,\\nu_2$ are the\n",
    "unknown mass and width of the new particle. Under the null hypothesis, those parameters \n",
    "don't mean anything... aka they \"are not identified under the null\" in the statistics jargon.\n",
    "This introduces a 2-d look elsewhere effect.\n",
    "\n",
    "The LEE correction in this case is based on \n",
    "\n",
    "\\begin{equation}\n",
    "E[ \\phi(A_u) ] = P(\\chi^2_1 > u) + e^{-u/2} (N_1 + \\sqrt{u} N_2) \\,\n",
    "\\end{equation}\n",
    "where \n",
    "   * $A_u$ is the 'excursion set above level $u$ (eg. the set of parameter points in $(\\nu_1,\\nu_2)$ that have a -2 log-likelihood ratio greater than $u$ )\n",
    "   * $\\phi(A_u)$ is the Euler characteristic of the excursion set\n",
    "   * $E[ \\phi(A_u) ]$ is the expectation of the Euler characteristic of those excursion sets under the null\n",
    "   * $P(\\chi^2_1 > u)$ is the standard chi-square probability \n",
    "   * and $N_1$ and $N_2$ are two coefficients that characterize the chi-square random field.\n",
    "   \n",
    "### structure of the notebook\n",
    "\n",
    "The notebook is broken into two parts.\n",
    "  * calculation of $N_1$ and $N_2$ based on $E[ \\phi(A_u) ]$ at two different levels $u_1$ and $u_2$\n",
    "  * calculation of LEE-corrected 'global p-value' given $N_1,N_2$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Populating the interactive namespace from numpy and matplotlib\n"
     ]
    }
   ],
   "source": [
    "%pylab inline --no-import-all\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from lee2d import *"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Test numerical solution to $N_1, N_2$  from the example in the paper"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Usage: calculate n1,n2 based on expected value of Euler characteristic (calculated from toy Monte Carlo) at two different levels u1, u2. For example: \n",
    "   * $u_1=0$ with $E[ \\phi(A_{u=u_1})]=33.5 $ \n",
    "   * $u_2=1$ with $E[ \\phi(A_{u=u_2})]=94.6 $\n",
    "\n",
    "would lead to a call like this"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "32.5 122.945875625\n"
     ]
    }
   ],
   "source": [
    "# An example from the paper\n",
    "n1, n2 = get_coefficients(u1=0., u2=1., exp_phi_1=33.5, exp_phi_2=94.6)\n",
    "print n1, n2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[]"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": [
       "iVBORw0KGgoAAAANSUhEUgAAAYsAAAESCAYAAAAMifkAAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\n",
       "AAALEgAACxIB0t1+/AAAHQ5JREFUeJzt3XmYXFWZx/FvpbMQshCDAhICFRYxKg4BWSJEStEhUQQZ\n",
       "eAQVZIIyMsPIKA6GjM7Y47gCDmrUKMQAigYEZXGQVYnGhbCGPZAEIlmgg+JCYCYk9Dt/nBu6qqnq\n",
       "vlVdt8659/w+z1NP171dXfWeS3hPnXs2EBERERERERERERERERERERERERERERERERERkVS6fAcw\n",
       "gKOBTwAfAP4CPOY3HBERCdkEYIHvIEREpHMWAj3A/f3OzwSWAyuAOf1+dx6wb/ahiYhIKGYA06it\n",
       "LLqAlUAZGAEsA6YCJeDLwOGdDVFERPob3uHPW4KrFKodiKssVifHl+H6K96OqyjGA3sC3+lIhCIi\n",
       "8jKdrizqmQSsqTpeCxwEfBSY5yUiERGpEUJlYUP425XAHu0KREQkEqtwd2yCVqa2z+Jg4Iaq47m8\n",
       "vJO7kaFUNJ3U7TuAlLp9B5BSt+8AUur2HUBK3b4DSKnbdwApdPsOIKWmc+ewLKJo0p3AXrhKZCRw\n",
       "PHCtz4BERKRWpyuLRcBvgdfg+ilmA1uAfwZuBB4CLgcebuI9u4FKO4MUESmoCvlp/bRVXm5DVXwH\n",
       "kFLFdwApVXwHkFLFdwApVXwHkFLFdwApVHwHkFJecmfbRFdgEZE2yGWfhYiIBK4IlUU3+Wn6iYj4\n",
       "VEF9FiIikpJuQ4mISPupshARkUGpshARkUEVobLoRh3cIiJpVFAHt4iIpKQObhERaT9VFiIiMihV\n",
       "FiIiMihVFh1hJd8RiIgMhSqLzNlrgdvA3uA7EhGRWOVgNJSVwD4E9jTYJ8G6fEckItHLQe5srxwV\n",
       "2Mpgi8F+Dba772hEJGo5yp3tkbMC2zCwM5NWxqnqyxART3KWO4fOyOUMbns92D1g14Dt4DsaEYlG\n",
       "BZczo6wscspGgn0RbD3Ykb6jEZGo5Dh3tqYABbYZYI+DzQcb4zsaEYlCAXJncwpSYNsO7Htgj4C9\n",
       "yXc0IlJ4BcmdzhRgAXDFAK8pVIHBjgfbADZXQ2xFJEMFy51ORJUFgO2aDLFdDDbZdzQiUkjBrzq7\n",
       "EOgB7u93fiawHFgBzOlwTIEpPQEcDtwE3AV2nOeAREQ6bgYwjdrKogtYCZSBEcAyYGrV7yNrWVSz\n",
       "A8FWgC1Q57eItFHwLYslwJ/6nTsQV1msBjYDlwFHAxOBbwP7Em1ro3Q7sB+uQr0bbH/PAYlIpIb7\n",
       "DgCYBKypOl4LHAQ8A5yW4u+7q54vTh4FUnoWmA12AnA92JeB86HU6zkwEcmPCrmbvOxuN1XfhjoW\n",
       "uLDq+ERgXsr3KvhtqP5sCthvwW4E28l3NCKSW8HfhqpnHVA96mcyrnUhL1N6HHgLcDvuttQRngMS\n",
       "EclMmdqWxXBgVXJ+JC/v4B5IZC2LavZWsLVg57qlQ0REUgs+dy4C1gObcP0Us5Pzs4BHcB3dc5t4\n",
       "v5wuJNgu9kqwn4LdDraH72hEJHgVtJBgrKwEdkay7Pn7fEcjIrkQXe6MrsCN2TSwR8G+qzkZIjKI\n",
       "6HJndAUemI1LFiR8CGwf39GISLCiy52R91k0Yh9Mbkt9RLvxiUiVCuqzkFq2N9gysB+5JdBFRF4S\n",
       "Xe6MrsDNsW3Avgm2CuwA39GISDCiy53RFbg1dlyyT8bHdFtKRIgwd0ZX4NbZ7sl8jGvAJvqORkS8\n",
       "ii53qoO7KTYS7L/BVoNN9x2NiHRcBXVwS3p2FFgP2FlgIawPJiKdFV3ujK7A7WO7gf0O7Dq3bIiI\n",
       "RCS63BldgdvLRiQLET4BdojvaESkY6LLndEVOBt2ZHJb6pO6LSUShehyZ3QFzo7tmtyW+h+w7X1H\n",
       "IyKZii53RlfgbNkIsPPAfq/RUiKFFl3ujK7AnWFHJZP4/lWT+EQKKbrcGV2BO8fKYEuTSXyv8B2N\n",
       "iLRVdLkzugJ3lo0E+yrY41pbSqRQosudmsHdEfZ3yW2pj+q2lEiuVdAMbsmW7QF2F9gVWvJcJPei\n",
       "y53RFdgv2wbsW2ArwPb1HY2ItCy63BldgcNg70924jtVt6VEcqlQuXMMcAlwAfD+Bq8pVIHzxV4L\n",
       "9kCy5/cY39GISFMKlTtPAt6VPL+swWsKVeD8sTFglySVxmt9RyMiqTWdOzu9DtBCoAe4v9/5mcBy\n",
       "YAUwJzk3CViTPH+xI9FJk0rPAX8PfA1YAvY+v/GISFHMAKZRW1l0ASuBMjACWAZMBU6kr2WxqMH7\n",
       "qWURDJsGtjLZ83uU72hEZEC5yJ1laiuL6cANVcdnJ49tcS2RbwGNvrHmosDxsO3AfgJ2h5sBLiKB\n",
       "ajp3Ds8iiiZV324CWAscBDwPnJLi77urni9OHuJF6S9gxwIfB5aCnQKl63xHJSJUyOHk5TK1LYtj\n",
       "gQurjk8E5qV8L7UsgmWHgK0B+wJYCF9KRKRP8B3c9awDJlcdT8a1LiTXSr8B9gcOAG4G28lzQCKS\n",
       "M2VqWxbDgVXJ+ZH0dXCnobWhgmddYP8JthbsLb6jEYlchZysDbUIWA9swvVTzE7OzwIewY2KmtvE\n",
       "+wVfYNnKZoI9BXaWZn2LeBdd7oyuwPlmuyZ7ZFwFNsF3NCIRiy53Rlfg/LNRYN9I5mT8je9oRCIV\n",
       "Xe5Un0VuvbQY4ezBXysibVIhJ30W7RZdgYvFXgf2MNgCsNG+oxGJSHS5M7oCF4+NA7sM7G6w3X1H\n",
       "IxKJ6HJndAUuJislW7b2gL3bdzQiEYgud0ZX4GKz6WBPaNa3SOaiy53RFbj47FVgt4D9HGwH39GI\n",
       "FFR0uVOjoQrJusA+l6wt9Wbf0YgUSAWNhpLisSPBNoCdkcz6PgQmXAkTfoy+IIgMRXS5M7oCx8d2\n",
       "dyOlnrgFdnwOvmEw32Dc88DbfUcnklPR5c7oChwnGw0/ewKeMbDkcbHBxFt8RyaSU21fonw/4Fxg\n",
       "KW7v7KeS5+fitkcV6YDS/8KJD8JDVee2AUoaMSUSgJ8BPwCOB6bg/u8cDewOnJD8zvcuaGpZxGMW\n",
       "vOI5uMlgo8F3NsN2x/kOSiSn2po7d0zxGt9DG1VZxOUomPgb2PsO+P1dYEvAXu07KJEcii53Rldg\n",
       "2cqGgf0H2DptqiTStExy50bg2eSxCegF/prFB7VAlUX07IhkmZBPaFMlkdQyz53DgPcAX8r6g1JS\n",
       "ZSGA7QZ2B9gVbmFCERlEx3Lnsk590CA0g1sStg3YBWDL3dLnIlJHhRZncKdpth9b9XwYsD9wGDC9\n",
       "2Q/LgJGuDBINmw2cA5wOpR/5jkYkUE3nzjQvvpi+WmgLsBq4ENjQzAdlRJWF1GH7AVcCVwNzoLTZ\n",
       "c0AioYkud6rPQhqwiWDXgf1Kw2tFXqZQuXMKsAC4YoDXFKrA0m41w2sP9R2NSEAKmTtVWcgQ2cxk\n",
       "eO2/aHitCBBo7lyIW1fq/n7nZwLLgRXAnAH+XpWFtIFNSfb5vgxsrO9oRDzLNHceA7Qyhn0GbtHB\n",
       "6sqiC1gJlIERuKG4U4GTgPOBnateq8pC2sRGgy0EewDsNb6jEfEos9y5B7AZ+McW/75MbWUxHbih\n",
       "6vjs5FFtIvBtBm55qLKQJlkJ7B+STZXe4zsaEU+azp1pl3g+BTd2fTYwv9kPqWMSsKbqeC1wUL/X\n",
       "PAOcluK9uqueL04eIg2UDLgAbBlwBdhBwKeh9KLnwESyVKEDk5eHA4/ibkHdDPxNC+9RprZlcSxu\n",
       "rsZWJwLzWnhftSxkCOxVYLeA3eyei0Sj7ZsfAcwCfodbSPAi4EPNfkgd64DJVceTca0LkQ4qPQ0c\n",
       "AdwB3Al2oOeARHLtGtzyHuA2P3oMGNnke5SpbVkMB1Yl50fS18HdLK0NJW1ix4A9nfRnaHitFFWF\n",
       "FteGGswreHkfwDm41kZai4D1uOXN1+D6PUje4xHcqKi5Lcan21DSRrY32INg33ULE4oUVnS5M7oC\n",
       "S9ZsLNiPwO4EK/uORiQj0eXO6AosnWAlsDPBngL7W9/RiGQgutypPgvJkB0Gth7sU26dKZHcq5BR\n",
       "n0XooiuwdJpNAvsd2NVg2/mORqRNMsmdEwd5+KTKQjrARoF9C+xRsNf7jkakDTLJnauBXuCPyaM3\n",
       "Ofc4bhitT6ospIPs5GR47Xt9RyIyRJnkzguBd1YdzwIuyOKDWqDKQjrMpoE9DnYeWNrlckRCk0nu\n",
       "fCDlOR9UWYgHtj3YDWC3gu3gOxqRFmSy3Md64NO42dZTgE/hlusIRTcaDSUdVfoj8C7g18BdYAd7\n",
       "DkgkrQq1i6+21fbA14F7ksfX8N+xvZVaFuKZHZUsd/4RLRMiORJd7oyuwBIie02yoZKWCZG8iC53\n",
       "RldgCZWNBbsc7A5YvDduh8j9SXerV6TTosud0RVYQmYlWPl56HkRTt0IO2+E8T/HbR0sEpKO5c5m\n",
       "lyjPiioLCcyEW+HSF8EMthhUnoPSGb6jEuknk9z5S9woqK0OBO7L4oNaoMpCAjN+HTxgrrIwg0cN\n",
       "dr3Id1Qi/WSSO48AlgOnA1/AjYjaL4sPaoEqCwnMhOvgzBeg12CjwVWb4Q/r3V4ZIsHILHe+FdgC\n",
       "PAnslNWHtECVhYRmRxi3HF61Ecb8H4z7AWz+h2R47Xt8ByeSyCR3/jtuxvZ04CO43e2OzOKDWqDK\n",
       "QkLUBewFTOo7ZQeCPQH2ObAuX4GJJDLJnV/F7b291W7AzVl8UAu0n4XkiO0A9otkqZBQJrZKXCpo\n",
       "PwuRPLDhYOeCPQa2r+9oJFqZ5M4dgPOA64Fbk8cvsvigFqiykJyy45Plzk/yHYlEKZPceTPwYdyI\n",
       "qMOAi4BzsvigFqiykByzN4CtAPs6mCbuSSdlkjvvTn5Wz624M4sPquNo3N4ZlwHvqPN7VRaSczYB\n",
       "7H/AloCFNNJQii2T3Hlb8vMm3Cio/YBVWXzQACYAC+qcV2UhBWDDwLrB1oBN9x2NRCGT3PluXLLe\n",
       "B1iMa2kc1eR7LAR6gPv7nZ+Ju721ApgzwN+fB9TrDFRlIQViRybzMU7TcueSsWBz5wxgGrWVRRew\n",
       "Erep0ghgGTAVOAk4H9gZKAFfBg5v8L7BFlikNbaXljuXDsgkd+6OS95XAT9NHte28D5laiuL6cAN\n",
       "VcdnJ49qZ+D6R+bjJgT2p8pCCqh6uXPb1Xc0UkhN5840G85fjesv+CnQ2+oH1TEJWFN1vBY4qN9r\n",
       "vp48BtJd9Xxx8hDJsdJGsBOATwBLwd4PpVt9RyW5VqEDk5dvb9P7lKltWRwLXFh1fCIwr8n3VMtC\n",
       "Cs4OB3sK7Ez1Y0gbZdKymIf79n4jsKnq/N11X53eOmBy1fFkXOtCRF5S+jnYQcBPgAPAPgyl53xH\n",
       "JVLPl3CJ/Zf0zeBupUlcprZlMRw3BLeM20xpawd3M7Q2lETCRoNdDHYf2J6+o5HcqpDh2lCrGPrO\n",
       "eIuA9biWyRpgdnJ+Fm4V25XA3BbeV7ehJCJWAjsdrAdslu9oJNcyyZ1XAztm8cZtoMpCImSHgq0D\n",
       "+3c3oU+kaU3nzjQdZr8E3gjcQV+fhdH8xLwsGOnKIFIwtjNwJbAB+CCU/uo5IMmXpnNnmg7uzzT4\n",
       "oFB0oyGzEp3SerAK8DXgdrBjoPSw56AkfBUy6ONNU+v4/lYfUqUl4omdkix3fozvSCQ32po7fwmc\n",
       "Bbymzu/2xq3l9Kt2fmALVFmIAGAHJNu2fl7btkoKbc2do4BTcPtZPAk8ilvw78nk3N8z9FFSQ6XK\n",
       "QuQltgPYYrDrtW2rDCKz3NmFGxG1Y/I8FKosRGrYCLDzwVaBvdF3NBKstubO0cDHgW/iFvFL0xne\n",
       "aaosROqyDyT9GCf4jkSC1Nbc+SPgUuA04BrcqIvQaAa3SEO2L9hjYOeBhfhlTzqvQgYzuPsvzXFP\n",
       "O9+8TdSyEBmQTQS7EeznYK/0HY0Eo+ncOdDszy0NnotIbpSeAd6JWz36TrD9PAckBfQi8GzVY0vV\n",
       "81Bmi6plIZKaHZf0Y3wwObEnbmfKmQz8xVGKJ7rcGV2BRYbGXg+2Au69FsY/B0c/C3s9C+OvRxVG\n",
       "TKLLndEVWGTobALc/AL82cAMXjB4/bOAZoDHI5PNj0SkUEp/ha6uvnVBRwAHd8GDk3xGJWFTs1Mk\n",
       "Pr0w5l747IvuC+YWYEwXcJvnuEQyo9tQIq3ZBcY9AMO3wD6b4Jn1YPPBfC/hI50RXe7UpDyRoRkD\n",
       "DAMbD3YV2G+SvTKkmCpkuK1qyKIrsEh2bBjYp5Nd+N7sOxrJVHS5M7oCi2TP3gW2Aew0t++3FFB0\n",
       "uTO6Aot0hu0F9iDYArBtfEcjbVeo3PlaYD5uQcMPNXhNoQosEhYbB3Yl2FKwXXxHI21VyNw5DFdh\n",
       "1FPIAouEw0pgZ4OtB5vhOxppmyBz50Kgh9pVbMGtR7Mct/venAZ/+27geuDvGvw+yAKLFI8dAdYD\n",
       "drr6MQohyNw5A5hGbWXRBawEyrjpo8uAqbhFzc4H+g/du6bBewdZYJFisj3A7gO7SP0YuRds7ixT\n",
       "W1lMB26oOj47eVQ7DLfh0neAjzV432ALLFJMNgbscrA7wCb7jkZaFmzuLFNbWRwHXFh1fCIwr4X3\n",
       "DbbAIsVlJbCzwJ4EO8x3NNKS3Cwk2M4k3131fHHyEJHMlAw4F2wZcDnYF4B5yXkJU4WcrHRRprZl\n",
       "cTC1t6Hm0riTeyD6xynilU1xlYZ9D2y072gktWBzZ5mX7+m9Kjk/kr4O7mYFW2CReNi2YD8Euwts\n",
       "N9/RSCpB5s5FwHrc4vlrgNnJ+VnAI7hRUXNbfG8tJCgSBCuBnQn2FNjbfEcjDVXQQoIi4p8dnlQY\n",
       "H9d8jKBFlzujK7BI+KwMdg/Ype4WlQSo6dypnfJEpM1Kq4FDkoNfqx9DQqA+C5FgWSm5HaV+jHBU\n",
       "UJ+FiITJ3pZM4DtT/RjBiC53RldgkXyy3cDuBvuB+jGCoD4LEQlR6fe4foxe4DeuE1ykc9SyEMkV\n",
       "K4F9LOnHONx3NBGLLndGV2CRYrC3qh/Dq+hyp0ZDieSW+jE8qKDRUCKSPzYa7PvJJL6y72giEl3u\n",
       "jK7AIsWjfgwPosud0RVYpLjUj9FB0eXO6AosUmzqx+gQzbMQkTzTfAzJhloWIoVU3Y9x8T8BHwbe\n",
       "Buj2VHtElzujK7BIXD7xbejphUs3wS4bYfx83xEVRHS5M7oCi0RkAozaBOsMzOAFg0nPAfv4DqwA\n",
       "ouyz6EaT8kSKaHsYvxl2Tg5HALeOgve+0WdQOVfB5czoqGUhUlwjYNseWNgLvQY3Gpz9f7Blg/bH\n",
       "GLLocmd0BRaJzBtg7GoY1gvb/gGoJPtjPJV0gKvDuzXR5c7oCiwSqZG1hy/t8/19t2SINKlwuXMM\n",
       "cAfwrga/L1yBRSQt2xbsh2B3ge3qO5qcKVwH9yeBy30HISIhKj0PfAD4IbAUrOI3HhmqhUAPcH+/\n",
       "8zOB5cAKYE6dv3sHcDxwMmpZiMiA7O1JP8YZ6sdIJcjcOQOYRm1l0QWsBMq48XDLgKnAScD5uLFy\n",
       "n0ue3whcTf2Zm0EWWER8sClg94Jdon6MQQWbO8vUVhbTgRuqjs9OHvWcDLyzwe+CLbCI+GBjwC4D\n",
       "uxNssu9oApabPotJwJqq47XJuXouAX6WeUQiUgCl54D34fo6l4K9xXNAhTHc0+e2s0XQXfV8cfIQ\n",
       "kWiVDDjX3ZLiCrD/Ar6ZnI9VhZysdFGm9jbUwdTehppL/U7uwcT8H19EBmV7gN0PthBsG9/RBCTY\n",
       "3FmmtrIYDqxKzo+kr4O7WcEWWERCYWPBrgC7HWwX39EEIsjcuQhYD2zC9VPMTs7PAh7BjYqa2+J7\n",
       "G1pIUEQGZSWwuWDrwA71HY1HFVzODLKyyFJ0BRaRobBZYD1gp0U+HyO63BldgUVkqGxPsAfALgAb\n",
       "5TsaT6LLndEVWETawcaB/Rjst2A7D/76wokud6rPQkRaZMPAPg22Fmy672g6pIL6LEREWmFHgm0A\n",
       "O9V3JB0UXe6MrsAikgXbG+xhsG+BjRz89bkXXe6MrsAikhUbD3YN2BKwnXxHk7Hocmd0BRaRLNkw\n",
       "sM+APQF2gO9oMhRd7oyuwCLSCXZ00o9xsu9IMhJd7tRoKBHJiL0O7FGwr4ON8B1Nm1TQaCgRkXaz\n",
       "CWA/A1sM9irf0bRRdLkzugKLSKdZF9jnwFaD7ec7mjaJLndGV2AR8cWOA3sa7AO+I2mD6HJndAUW\n",
       "EZ9sH7BVYF8B87V5XDtElzujK7CI+GYTwW4Cuxlse9/RtCi63BldgUUkBNYFdg7YY2BvxG3oNhtK\n",
       "nwXeA4S+/Hl0uTO6AotISOwE6H0aPr4M3rQR/sOgvBHGnOM7skFElzujK7CIhKb7ZFjTC1sMzOBp\n",
       "g5EvABN8RzaA6HKnJuWJiG8z4Yi/uIrCDHoNtnsemOw7sDoqaFKeiIgX28O2f4KLe+EPBp/aDGMf\n",
       "BYb5DmwA0eXO6AosIkHaF7Z7AEb/Fbb7FTDJd0CDiC53RldgEZE2aDp3htxMqgBLgPnAYX5DERGJ\n",
       "W8iVRS/wLDAKWOs5FhERydhCoAe4v9/5mcByYAUwp87fbZ3UsgNwaYP3zsttqIrvAFKq+A4gpYrv\n",
       "AFKq+A4gpYrvAFKq+A4ghYrvAFIK8jbURbiKoVoX8I3k/OuA9wFTgZOA84Gd6SvMn3Gtizyr+A4g\n",
       "pYrvAFKq+A4gpYrvAFKq+A4gpYrvAFKo+A4gK51YCGsJUO537kBgJbA6Ob4MOBr4EvD95NwxwBG4\n",
       "iS3zsg5SREQa87Vq4iRgTdXxWuCgfq+5KnmIiIhnnVrsqgz8FNgnOT4Wdwvq1OT4RFxl8dEm33cl\n",
       "sEcb4hMRickqYM9m/sBXy2IdtVPhJ9PaiKemCisiImErUzsaajiuZisDI4FluA5uERGJ1CJgPbAJ\n",
       "108xOzk/C3gEdytprp/QRESk6AabpxGK1cB9wD3A7X5DqVFv/stE4GbgUeAmwlhiuV6c3bjblvck\n",
       "j/5DszttMnAr8CDwAHBGcj6069kozm7Cup7bAEtxdxweAr6YnA/tejaKs5uwrie46Qr34PqOIbxr\n",
       "mZkuXIukDIwg7NtYj+P+w4RmBjCN2iR8DvDJ5Pkc3FBm3+rF+RngTD/h1LUTsG/yfCyuxTyV8K5n\n",
       "ozhDu54A2yY/hwO3AYcS3vWE+nGGeD3PBH4AXJscN30tQ17uYyDV8zQ20zdPI1QhbrG4BPhTv3NH\n",
       "AZckzy/BbQ/pW704Iaxr+hTuCwvARuBh3PDw0K5nozghrOsJ8HzycyTuy+GfCO96Qv04IazruQvw\n",
       "TmABfXE1fS3zWlnUm6cR6pLABtwC3EnfUOFQ7Yi75UPyc0ePsQzmo8C9wHcJqwldxrWElhL29Szj\n",
       "4rwtOQ7teg7DVWw99N06C/F61osTwrqe5wNn4dbb26rpa5nXyiIva0IBHIL7n3IWcDrutkoeGOFe\n",
       "5/nAFNwtlSeBr/gN5yVjgR8D/4JbBLNaSNdzLHAlLs6NhHk9e3Hx7AK8BXhrv9+Hcj37x1khrOt5\n",
       "JLAB11/RqLWT6lrmtbJo1zyNTngy+fk0bkb6gR5jGUwP7r42wKtx/8hCtIG+f+ALCOOajsBVFN8H\n",
       "rk7OhXg9t8Z5KX1xhng9t/oLcB2wP2Fez622xvkmwrqeb8bdcnocNzL1bbh/o01fy7xWFncCe9E3\n",
       "T+N4+jpuQrItMC55Pgb4W16++m5IrgVOTp6fTF8yCc2rq54fg/9rWsLdbngI+GrV+dCuZ6M4Q7ue\n",
       "r6Tv1s1o4B24b8ahXc9Gce5U9Rrf1/PfcF+mpwAnAL/ALdga2rXMVB7maUzB3c9chhuqGFKcW+e/\n",
       "vEDf/JeJuP6VkIbT9Y/zFOB7uOHI9+L+kfu+d30o7nbEMmqHS4Z2PevFOYvwruc+wN24OO/D3W+H\n",
       "8K5nozhDu55bHUbfl+rQrqWIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiISGNlapd4+FfcHgci\n",
       "uZXXtaFE8iSE1VFFhkSVhYiIDEqVhUj7baH2/63RvgIRaRdVFiLt1wPsgFvZcxRuAxqRXBvuOwCR\n",
       "AtoMfBa4HbdR10Oo30JEREREREREREREREREREREREREREREREREREQkDv8PVyLAEq1UyTIAAAAA\n",
       "SUVORK5CYII=\n"
      ],
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10cd61050>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# reproduce Fig 5 from paper (the markers are read by eye)\n",
    "u = np.linspace(5,35,100)\n",
    "global_p = global_pvalue(u,n1,n2)\n",
    "plt.plot(u, global_p)\n",
    "plt.scatter(35,2.E-5) #from Fig5\n",
    "plt.scatter(30,2.E-4) #from Fig5\n",
    "plt.scatter(25,2.5E-3) #from Fig5\n",
    "plt.scatter(20,2.5E-2) #from Fig5\n",
    "plt.scatter(15,.3) #from Fig5\n",
    "plt.xlabel('u')\n",
    "plt.ylabel('P(max q > u)')\n",
    "plt.semilogy()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Check Euler characteristic from Fig 3 example in the paper"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x10ca5cb90>"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": [
       "iVBORw0KGgoAAAANSUhEUgAAAPYAAAD7CAYAAABZjGkWAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\n",
       "AAALEgAACxIB0t1+/AAACZBJREFUeJzt3U+MXVUdwPHvm5k2QEskBqNAJ7YbDCbEQrAxtNARlBSC\n",
       "sHEhidFgwkql0aRB3dQdcVVNjBuhBMSURRFig0ZBO6VoLCIdLLQllLRJC+1AUlKDjabTeS7OHZlO\n",
       "Onl35v457/3m+0le5m/fOdPMd+679937DkiSJEmSJEmSJEkDq1PDfRwDPl3D/UhauD3A2NxPDtVw\n",
       "x0atbDqdTqVb1fvodruVblu3bq3y73cBGy/2/1JH2JL6jGFLARm2lNHY2Fgj91vHwbNuDfchLcrM\n",
       "fnIu09PTOYff1el0vsJFOnaLLQVk2FJAhi0FZNhSQIYtBVQm7E3AYeAt4KFmpyOpDr3CHgZ+Tor7\n",
       "s8B9wHVNT0pSNb3CXgccIV3ocQ54Cri34TlJqqhX2NcAx2d9fKL4nKQ+1itszyqTBtBIj6+/A4zO\n",
       "+niUtNWWlMH4+Djj4+MATE5Ofma+7+t1ou0I8CZwO/Au8DLpANqhWd/jVl3ZeK74xc8V77XFngK+\n",
       "A/yBdIT8US6MWlIf6hU2wO+Lm6QB4ZlnUkCGLQVk2FJAhi0FZNhSQIYtBWTYUkCGLQVk2FJAhi0F\n",
       "ZNhSQGXOFe9rua/uOX/+fNbxc+t2817cNzIy8L/CjXCLLQVk2FJAhi0FZNhSQIYtBWTYUkCGLQVk\n",
       "2FJAhi0FZNhSQIYtBVQm7O3AJHCg4blIqkmZsB8jrY8taUCUCXsv8EHTE5FUH/expYC8mFUaIHUt\n",
       "oztjNbALuP4iX8t6pb0vtJDXUn+hhX5dRteH4lJAZcLeAfwVuBY4Dtzf6IwkVVbmccx9jc9CUq18\n",
       "KC4FZNhSQIYtBWTYUkCGLQVk2FJAhi0FZNhSQIYtBWTYUkCGLQVUyzVvOS+dzH3ZZO7LRpe63JeN\n",
       "Dg3157axP2clqRLDlgIybCkgw5YCMmwpIMOWAjJsKSDDlgIybCkgw5YCMmwpoDJhjwK7gTeA14EH\n",
       "G52RpMrKXARyDvgeMAGsBP4BPA8canBekioos8U+RYoa4ENS0Fc3NiNJlS10H3s1cAOwr/6pSKrL\n",
       "Qq7HXgnsBDaTttz/N/eaWK9Rlpoxp7V518cuu8VeBjwNPAk8O/eLnU7ngpukZsxp7c35vq9M2B3g\n",
       "UeAg8NO6JiipOWXCXg98HfgisL+4bWpyUpKqKbOP/RKeyCINFIOVAjJsKSDDlgIybCkgw5YCMmwp\n",
       "IMOWAjJsKSDDlgIybCkgw5YCqmV97Jy8TDSv3P//uddHHx4ezjr+fNxiSwEZthSQYUsBGbYUkGFL\n",
       "ARm2FJBhSwEZthSQYUsBGbYUkGFLAZUJ+xLSInwTpNVAHm50RpIqK3MRyH9Iq4CcLb7/JWBD8VZS\n",
       "Hyr7UPxs8XY5MAycbmY6kupQNuwh0kPxSWA36SG5pD5VNuxpYC2wCrgVGJv9xW63e8FNUisqr489\n",
       "4wzwHHDT7E+6PraURaX1sa8ErijevxT4MmkpXUl9qsxR8auAx0l/BIaAXwF/anJSkqopE/YB4Mam\n",
       "JyKpPp55JgVk2FJAhi0FZNhSQIYtBWTYUkCGLQVk2FJAhi0FZNhSQIYtBVTL+tg5r8HOvT5x7vWZ\n",
       "c8t9/f3IyMAv8d4It9hSQIYtBWTYUkCGLQVk2FJAhi0FZNhSQIYtBWTYUkCGLQVk2FJAZcMeJq3+\n",
       "savBuUiqSdmwN5NW2HTFPWkAlAl7FXAX8AjginvSACgT9jZgC2kpXUkDoFfYdwPvkfav3VpL/WXR\n",
       "62PfDNwDHAV2ALcBT9Q3L0kVLHp97B8Bo8Aa4GvAn4Fv1DcvSU1Y6PPYHhWXBsBCXjBqT3GT1Oc8\n",
       "80wKyLClgAxbCsiwpYAMWwrIsKWADFsKyLClgAxbCsiwpYAMWwpo4BcXnp7O+/oPudfnnpqayjp+\n",
       "7p8/9/rc/cotthSQYUsBGbYUkGFLARm2FJBhSwEZthSQYUsBGbYUkGFLARm2FFDZc8WPAf8CzgPn\n",
       "gHVNTUhSdWXD7gJjwOnmpiKpLgt5KO5qm9KAKBt2F3gBeAV4oLnpSKpD2Yfi64GTwCeA54HDwN6m\n",
       "JiWplEWvjz3jZPH2feAZPHgm9YNFr48NcBlwefH+CuAO4EANk5LUkDIPxT9J2krPfP+vgT82NiNJ\n",
       "lZUJ+yiwtumJSKqPZ55JARm2FJBhSwEZthSQYUsBGbYUkGFLARm2FJBhSwEZthSQYUsBDfz62Lkt\n",
       "9fW5XZ+6P7nFlgIybCkgw5YCMmwpIMOWAjJsKSDDlgIybCkgw5YCMmwpoDJhXwHsBA4BB4EvNDoj\n",
       "SZWVOVf8Z8DvgK8W37+i0RlJqqxX2B8DbgG+WXw8BZxpdEaSKuv1UHwNaSG+x4BXgV+S1vKS1Md6\n",
       "hT0C3Aj8onj7b+AHTU9KUjW9wj5R3P5efLyTFLik/Ba9PvYp4DhwbfHxl4A3apqUpGrmXR+7zFHx\n",
       "75KWzl0OvA3cX9OkJDWkTNivAZ9veiKS6uOZZ1JAhi0FZNhSQIYtBWTYUkCGLQVk2FJAhi0FZNhL\n",
       "nGtvxWTYUkCGLQVU1zK6r1b4t1cBJ2uah+M7/lIa/+35vtBZ5B3ONg5srOF+JC3cHmAs9yQkSZIk\n",
       "9YVNwGHgLeChlsfeDkwCB1oed8YosJv0UlOvAw+2PP4lwD5ggrQQxMMtjw8wDOwHdmUYG+AY8M9i\n",
       "Di+3PHbYhTiGgSPAamAZ6RfsuhbHvwW4gXxhfwpYW7y/kvT6VW3+/PDRS0mPAH8DNrQ8/vdJL7v1\n",
       "25bHnXEU+HimsR8HvlW8P0J6Df/a5Hweex0p7GPAOeAp4N4Wx98LfNDieHOdIv0xA/iQ9Jf76pbn\n",
       "cLZ4u5z0h/Z0i2OvAu4CHqGeZ2cWK8fYMwtxbC8+rn0hjpxhX0N6BdQZJ4rPLUWrSY8e9rU87hDp\n",
       "j8skabfgYItjbwO2ANMtjjlXF3gBeAV4oMVxG1+II2fYnqScrCTta20mbbnbNE3aHVgF3Ep7z4fe\n",
       "DbxH2rfNubVeT/qDeifwbdJWtA2NL8SRM+x3SAeQZoySttpLyTLgaeBJ4NmM8zgDPAfc1NJ4NwP3\n",
       "kPZxdwC3AU+0NPZsM2d8vQ88Q9o9bEPohThGSKfErSbt47V98Ixi7FwHzzqkX+Ztmca/knRkFuBS\n",
       "4EXg9gzz2Eieo+KXAZcX768A/gLc0eL4L/LRQhw/Bn7S4tiNu5N0NPgI8MOWx94BvAv8l7Sv3/ZC\n",
       "CBtID4UnSA9J95Oe/mvL9aT9uwnSUz5bWhx7to3kOSq+hvSzT5Cebmz79+9zpC32a8BvqPmouCRJ\n",
       "kiRJkiRJkiRJktSo/wGyjceY1YEC0AAAAABJRU5ErkJggg==\n"
      ],
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10ca5a610>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#create Fig 3 of http://arxiv.org/pdf/1105.4355v1.pdf\n",
    "a = np.zeros((7,7))\n",
    "a[1,2]=a[1,3]=a[2,1]=a[2,2]=a[2,3]=a[2,4]=1\n",
    "a[3,1]=a[3,2]=a[3,3]=a[3,4]=a[3,5]=1\n",
    "a[4,1]=a[4,2]=a[4,3]=a[4,4]=1\n",
    "a[5,3]=1\n",
    "a[6,0]=a[6,1]=1\n",
    "a=a.T\n",
    "plt.imshow(a,cmap='gray',interpolation='none')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#should be 2\n",
    "calculate_euler_characteristic(a) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Try a big matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#Fully filled, should be 1\n",
    "randMatrix = np.zeros((100,100))+1\n",
    "calculate_euler_characteristic(randMatrix)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": [
       "iVBORw0KGgoAAAANSUhEUgAAAP0AAAD/CAYAAAA6/dD3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\n",
       "AAALEgAACxIB0t1+/AAACktJREFUeJzt3EGMXHUdwPHv2K0IVLv00lYoKRIaOChCjCAh2QcpCRKs\n",
       "nAATSIOpNwPECG09DTcgIcULF0rIhoiKxTTUiG5FXjmYIEirQFsLalMK6YKQKnKCMB7+f3Y3053d\n",
       "t519nTf9fT/JZOe9fbP70+53/m/ekAFJkiRJkiRJkiRJQ+x64CDwBrB5wLNIqtkS4E1gLbAU2Adc\n",
       "MsiBJM3vc3089puk6A8DHwO/AL67CDNJqlE/0Z8LvDVj+2jeJ6nB+om+M98BY2NjnXycN2/eTuEt\n",
       "tzerkV7fqOBtYM2M7TWk1X7Knj17GBsboygKAIqimLrfZO12m3a7PegxKhu2eWH4Zm76vGVZUpbl\n",
       "1PZ9993X89h+on8ZuIh0Ie8d4Bbge90HFUXR6P+zpNNB94JaV/SfAD8Efk+6kv8YcKCPnyfpFOgn\n",
       "eoBn862nYTid7zZsMw/bvDB8Mw/bvHNp1fzzO51Oz+sJkmrSarWgR9/9XL2XNISMXgrG6KVgjF4K\n",
       "xuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG\n",
       "6KVgjF4KxuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG6KVgjF4KxuilYKpEvwZ4\n",
       "HngdeA24M+9fAewGDgETwGgdA0paXK0Kx6zKt33AMuAvwE3AHcC/gQeBzcA5wJaux3Y6nc6iDSup\n",
       "mlarBT36rrLSHyMFD/A/4ABwLrABGM/7x0lPBJIabqGv6dcClwEvAiuBybx/Mm9LariRBRy7DHga\n",
       "uAv4sOt7nXw7QbvdnrpfFAVFUSxoQEnzK8uSsiwrHVvlNT3AUuA3wLPAw3nfQaAgnf6vJl3su7jr\n",
       "cb6mlwag39f0LeAxYD/TwQM8A2zM9zcCO09+REmnSpWV/mrgBeBvTJ/CbwX+DDwFnA8cBm4Gjnc9\n",
       "1pVeGoC5Vvqqp/cny+ilAej39F7SacTopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG6KVgjF4K\n",
       "xuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG\n",
       "6KVgjF4KxuilYIxeCsbopWCMXgrG6KVgjF4Kpmr0S4C9wK68vQLYDRwCJoDRxR9NUh2qRn8XsB/o\n",
       "5O0tpOjXAc/lbUlDoEr05wE3ANuBVt63ARjP98eBmxZ/NEl1qBL9NuAe4NMZ+1YCk/n+ZN6WNATm\n",
       "i/5G4F3S6/lWj2M6TJ/2S2q4kXm+fxXpVP4G4AvAl4AnSKv7KuAYsJr0xDCrdrs9db8oCoqi6Gde\n",
       "SbMoy5KyLCsd22v1ns0Y8GPgO8CDwPvAA6SLeKPMfjGv0+l4EiCdaq1WC3r0vdD36T8r+H7gOtJb\n",
       "dtfmbUlDYCEr/clwpZcGYDFXeklDzuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG\n",
       "6KVgjF4KxuilYIxeCsbopWDm+4y8vk1MTNT9KyQtQO3Rb9u2re5fIWkBao/+yJEjdf8KSQtQ+2fk\n",
       "bd++veZfIanbpk2boEfffjCmdBrygzElTTF6KRijl4IxeikYo5eCMXopGKOXgjF6KRijl4IxeikY\n",
       "o5eCMXopGKOXgjF6KRijl4KpGv0osAM4AOwHrgBWALuBQ8BEPkZSw1WN/qfAb4FLgK8BB4EtpOjX\n",
       "Ac/lbUkNV+WTc5YDe4GvdO0/CIwBk8AqoAQu7jrGT86RBqDfT865AHgPeBx4BXgUOBtYSQqe/HVl\n",
       "v4NKql+V6EeAy4FH8tePOPFUvpNvkhquykdgH823l/L2DmArcIx0Wn8MWA28O9uD2+321P2iKCiK\n",
       "4qSHlTS7siwpy7LSsVU/DfcFYBPpSn0bOCvvfx94gLTyjzLLGYCv6aVTb67X9FWjvxTYDnwe+Adw\n",
       "B7AEeAo4HzgM3Awc73qc0UsDsBjRnyyjlwbAz72XNMXopWCMXgrG6KVgjF4KxuilYIxeCsbopWCM\n",
       "XgrG6KVgjF4KxuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxe\n",
       "CsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG6KVgjF4KxuilYKpEvxV4HXgVeBI4A1gB7AYO\n",
       "ARPAaF0DSlpc80W/FvgBcDnwVWAJcCuwhRT9OuC5vC1pCMwX/X+Bj4GzgJH89R1gAzCejxkHbqpr\n",
       "QEmLa77oPwAeAo6QYj9OWuFXApP5mMm8LWkIjMzz/QuBu0mn+f8BfgXc1nVMJ99m1W63p+4XRUFR\n",
       "FAufUtKcyrKkLMtKx7bm+f4twHXAprx9O3AlcC1wDXAMWA08D1w8y+M7nU7P5wNJNWm1WtCj7/lO\n",
       "7w+SIj8z/4D1wH5gF7AxH7MR2LkYg0qq33wrPcC9pLA/BV4hrfpfBJ4CzgcOAzeTXu93c6WXBmCu\n",
       "lb5K9P0wemkA+jm9l3SaMXopGKOXgjF6KRijl4IxeikYo5eCMXopGKOXgjF6KRijl4IxeikYo5eC\n",
       "MXopGKOXgjF6KRijl4IxeikYo5eCMXopGKOXgjF6KRijl4IxeikYo5eCMXopGKOXgjF6KRijl4Ix\n",
       "eikYo5eCMXopGKOXgjF6KRijl4KpPfqyLOv+FYtu2GYetnlh+GYetnnnYvSzGLaZh21eGL6Zh23e\n",
       "uXh6LwVj9FIwrZp/fgmM1fw7JJ1oD1AMeghJkiRJQ+t64CDwBrB5wLP0sgZ4HngdeA24M+9fAewG\n",
       "DgETwOhAputtCbAX2JW3mz7vKLADOADsB66g+TNvJf1dvAo8CZxB82ceqCXAm8BaYCmwD7hkkAP1\n",
       "sAr4er6/DPg7ac4HgXvz/s3A/ad+tDn9CPgZ8Ezebvq848D38/0RYDnNnnkt8E9S6AC/BDbS7JkH\n",
       "7lvA72Zsb8m3ptsJrCedoazM+1bl7aY4D/gDcA3TK32T511OCqhbk2deQVoAziE9Se0CrqPZM1dW\n",
       "1/v05wJvzdg+mvc12VrgMuBF0j/sZN4/yfQ/dBNsA+4BPp2xr8nzXgC8BzwOvAI8CpxNs2f+AHgI\n",
       "OAK8AxwnndY3eebK6oq+U9PPrcsy4GngLuDDru91aM7/nhuBd0mv53v9NxZNmhfSSnk58Ej++hEn\n",
       "nvU1beYLgbtJC8GXSX8ft3Ud07SZK6sr+rdJF8k+s4a02jfRUlLwT5BO7yE9i6/K91eTQmuCq4AN\n",
       "wL+AnwPXkuZu6ryQ/t2PAi/l7R2k+I/R3Jm/AfwJeB/4BPg16SVrk2eurK7oXwYuIj1Tfh64hemL\n",
       "Tk3SAh4jXVF+eMb+Z0gXbshfd9IMPyE9gV4A3Ar8Ebid5s4LKZS3gHV5ez3pqvgumjvzQeBK4EzS\n",
       "38h60t9Ik2duhG+TLoa8SXr7o4muJr023kc6Zd5LeqtxBeliWZPfmhlj+om06fNeSlrp/0paNZfT\n",
       "/JnvZfotu3HSGWHTZ5YkSZIkSZIkSZIkSZI0SP8HzxrHmksDOHQAAAAASUVORK5CYII=\n"
      ],
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10ca18f10>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# split in half vertically, should be 2\n",
    "randMatrix[50,:]=0\n",
    "plt.imshow(randMatrix,cmap='gray')\n",
    "calculate_euler_characteristic(randMatrix)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "6"
      ]
     },
     "execution_count": 38,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": [
       "iVBORw0KGgoAAAANSUhEUgAAAP0AAAD/CAYAAAA6/dD3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\n",
       "AAALEgAACxIB0t1+/AAACvhJREFUeJzt3G+MFPUdgPHnCmittlDeABWSM0YSUvtHY+sfaG5qTkKN\n",
       "obxSm2hIwb5qq7apCr653XdqYrRvfFGh5mq0rUVDJGnhkN5e0pgYDWD1gKJtCaDhUAhWjCYSri/m\n",
       "J1yOW3aP22Fn+T6fZLI7c7O3X9h7dmZ3LweSJEmSJEmSJEmSOtgyYA/wDvBQm2eRVLBpwLtANzAD\n",
       "2AksaudAkhr70hRu+33y6PcBnwN/An7cgpkkFWgq0V8OHBizfjBtk1RiU4l+tNEOPT09o2k/FxeX\n",
       "87ik9iY0lejfAxaMWV9AfrQ/ZWhoiJ6eHvr6+ujr62NwcJDR0dHSL319fQ332bx5M0uXLuXqq69m\n",
       "/fr1pZ+3bEu7Z57s49fueRstg4ODpzrr6+tjaGiobrjTpxD9G8BV5G/kvQ/cAfxk/E5ZllGpVKZw\n",
       "N5IaybKMLMtOrVer1br7TiX6E8AvgC3k7+SvB3ZP4ftJOg+mEj3A39JS19hnn07RaTN32rzQeTN3\n",
       "2rxnM5XX9E3pxP+sTpu50+aFzpu50+Y9m8Kjl1QuRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RS\n",
       "MEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIw\n",
       "Ri8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjDNRL8AGASGgbeBe9P22cBWYC8wAMwq\n",
       "YkBJrdVM9J8DvwK+CdwA/BxYBKwhj34hsC2tSyq5ZqI/BOxM148Du4HLgeVAf9reD6xo+XSSWm6y\n",
       "r+m7gWuA14A5wEjaPpLWJZXc9EnsexnwInAf8PG4r42m5QyVSuXU9SzLyLJsUgNKaqxWq1Gr1Zra\n",
       "t9noZ5AH/yywMW0bAeaSn/7PAw5PdMOx0UsqxvgDarVarbtvM6f3XcB6YBfw5JjtLwMr0/WVnH4y\n",
       "kFRizRzpFwN3Af8EdqRta4FHgBeA1cA+4PYC5pPUYs1E/w/qnxH0tnAWSeeBv5EnBWP0UjBGLwVj\n",
       "9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0\n",
       "UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8F02z0\n",
       "04AdwKa0PhvYCuwFBoBZrR9NUhGajf4+YBcwmtbXkEe/ENiW1iV1gGainw/cCqwDutK25UB/ut4P\n",
       "rGj9aJKK0Ez0TwAPACfHbJsDjKTrI2ldUgdoFP1twGHy1/NddfYZ5fRpv6SSm97g6zeRn8rfCnwZ\n",
       "+BrwLPnRfS5wCJhH/sQwoUqlcup6lmVkWTaVeSVNoFarUavVmtq3UfQPpwWgB/gNcDfwGLASeDRd\n",
       "bqz3DcZGL6kY4w+o1Wq17r6T/Zz+i9P4R4BbyD+yuzmtS+oAjY70Yw2lBeAo0Nv6cSQVzd/Ik4Ix\n",
       "eikYo5eCMXopGKOXgjF6KRijl4IxeikYo5eCMXopGKOXgjF6KRijl4IxeikYo5eCMXopmMn8EY1z\n",
       "MjAwUPRdtMX27ds5cuQIx48fZ3h4+IL9d16oIj9+9f7CbauMLlu2rOC7aI8PP/yQffv28dlnn9Hd\n",
       "3c38+fPbPZIm4UJ//DZv3gx1+i78SL9///6i76Itjh8/zqeffsqJEyc4evQoJ0+ebHwjlUbkx6/w\n",
       "I/26desKvov2+OKU8NixYyxdupTFixe3eyRNwoX++N1zzz3QriP96tWri76LttiyZQvDw8N0dXWx\n",
       "ZMkSVq1a1e6RNAkX+uOXop+Q795LwRi9FIzRS8EYvRSM0UvBGL0UjNFLwRi9FIzRS8EYvRSM0UvB\n",
       "GL0UjNFLwRi9FIzRS8E0G/0sYAOwG9gFXA/MBrYCe4GBtI+kkms2+t8CfwUWAd8G9gBryKNfCGxL\n",
       "65JKrpnoZwI/AH6f1k8AHwHLgf60rR9Y0fLpJLVcM9FfAXwAPANsB54GLgXmACNpn5G0Lqnkmol+\n",
       "OnAt8FS6/IQzT+VH0yKp5Jr5w5gH0/J6Wt8ArAUOAXPT5Tzg8EQ3rlQqp65nWUaWZec8rKSJ1Wo1\n",
       "arVaU/s2E/0h4AD5G3Z7gV5gOC0rgUfT5caJbjw2eknFGH9ArVardfdt9k9g/xJ4DrgI+DfwU2Aa\n",
       "8AKwGtgH3H4uw0o6v5qN/k3gexNs723hLJLOA38jTwrG6KVgjF4KxuilYIxeCsbopWCMXgrG6KVg\n",
       "jF4KxuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCM\n",
       "XgrG6KVgjF4KxuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCaiX4tMAy8BTwPXAzMBrYC\n",
       "e4EBYFZRA0pqrUbRdwM/A64FvgVMA+4E1pBHvxDYltYldYBG0f8P+Bz4CjA9Xb4PLAf60z79wIqi\n",
       "BpTUWo2iPwo8Duwnj/0Y+RF+DjCS9hlJ65I6wPQGX78SuJ/8NP8j4C/AXeP2GU3LhCqVyqnrWZaR\n",
       "Zdnkp5R0VrVajVqt1tS+jaK/DngVOJLWXwJuBA4Bc9PlPOBwvW8wNnpJxRh/QK1Wq3X3bXR6vwe4\n",
       "AbgE6AJ6gV3AJmBl2mclsPGcp5V0XjU60r8J/AF4AzgJbAd+B3wVeAFYDewDbi9uREmt1Ch6gMfS\n",
       "MtZR8qO+pA7jb+RJwRi9FIzRS8EYvRSM0UvBGL0UjNFLwRi9FIzRS8EYvRSM0UvBGL0UjNFLwRi9\n",
       "FIzRS8EYvRSM0UvBGL0UjNFLwRi9FIzRS8EYvRSM0UvBGL0UjNFLwRi9FIzRS8EYvRSM0UvBGL0U\n",
       "jNFLwRi9FIzRS8EYvRSM0UvBFB59rVYr+i5artNm7rR5ofNm7rR5z8boJ9BpM3favNB5M3favGfj\n",
       "6b0UjNFLwXQV/P1rQE/B9yHpTENA1u4hJEmSJHWsZcAe4B3goTbPUs8CYBAYBt4G7k3bZwNbgb3A\n",
       "ADCrLdPVNw3YAWxK62WfdxawAdgN7AKup/wzryX/uXgLeB64mPLP3FbTgHeBbmAGsBNY1M6B6pgL\n",
       "fDddvwz4F/mcjwEPpu0PAY+c/9HO6tfAc8DLab3s8/YDq9L16cBMyj1zN/Af8tAB/gyspNwzt92N\n",
       "wOYx62vSUnYbgV7yM5Q5advctF4W84FXgB9y+khf5nlnkgc0Xplnnk1+APg6+ZPUJuAWyj1z04r6\n",
       "nP5y4MCY9YNpW5l1A9cAr5E/sCNp+winH+gyeAJ4ADg5ZluZ570C+AB4BtgOPA1cSrlnPgo8DuwH\n",
       "3geOkZ/Wl3nmphUV/WhB37colwEvAvcBH4/72ijl+ffcBhwmfz1f73csyjQv5EfKa4Gn0uUnnHnW\n",
       "V7aZrwTuJz8QfIP85+OucfuUbeamFRX9e+Rvkn1hAfnRvoxmkAf/LPnpPeTP4nPT9XnkoZXBTcBy\n",
       "4L/AH4Gbyecu67yQP+4HgdfT+gby+A9R3pmvA14FjgAngJfIX7KWeeamFRX9G8BV5M+UFwF3cPpN\n",
       "pzLpAtaTv6P85JjtL5O/cUO63Eg5PEz+BHoFcCfwd+Buyjsv5KEcABam9V7yd8U3Ud6Z9wA3AJeQ\n",
       "/4z0kv+MlHnmUvgR+Zsh75J//FFGS8hfG+8kP2XeQf5R42zyN8vK/NFMD6efSMs+73fIj/Rvkh81\n",
       "Z1L+mR/k9Ed2/eRnhGWfWZIkSZIkSZIkSZIkSVI7/R+XkH7e+jgBRwAAAABJRU5ErkJggg==\n"
      ],
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10de96410>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#split in half horizontally twice, should be 6\n",
    "randMatrix[:,25]=0\n",
    "randMatrix[:,75]=0\n",
    "plt.imshow(randMatrix,cmap='gray')\n",
    "calculate_euler_characteristic(randMatrix)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "5"
      ]
     },
     "execution_count": 39,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": [
       "iVBORw0KGgoAAAANSUhEUgAAAP0AAAD/CAYAAAA6/dD3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\n",
       "AAALEgAACxIB0t1+/AAAC4tJREFUeJzt3X9onPUdwPF31mpn1dn2n/4wxeiwIHPrFDerdeSZxFJF\n",
       "OlFQB0qZdX9tUzemVv9J8p8KovvHP6adZKLbnJZiYavpXC8wBLHEOk3tqttCq9LUptRZUbA0++P5\n",
       "2oYk17uYe3rP9fN+wcPd89yT3Kde3vc8dxciSJIkSZIkSZIkSWphq4FdwLvA/U2eRVLBZgHvAR3A\n",
       "acAO4KJmDiSptq/N4Gu/Tx79MPAF8EfgRw2YSVKBZhL9ucDecevvp22SSmwm0Y/V2qGzs3Ms7efi\n",
       "4nISl9TelGYS/QfA0nHrS8mP9scMDAzQ2dlJd3c33d3dbNu2jbGxsdIv3d3dNffZsmULq1at4uKL\n",
       "L2bDhg2ln7dsS7Nnnu7j1+x5ay3btm071ll3dzcDAwNVw509g+i3AxeSv5H3IXAL8OOJO2VZRk9P\n",
       "zwzuRlItWZaRZdmx9d7e3qr7ziT6I8DPgZfJ38nfALwzg+8n6SSYSfQAf01LVeOffVpFq83cavNC\n",
       "683cavOeyExe09elFf9jtdrMrTYvtN7MrTbviRQevaRyMXopGKOXgjF6KRijl4IxeikYo5eCMXop\n",
       "GKOXgjF6KRijl4IxeikYo5eCMXopGKOXgjF6KRijl4IxeikYo5eCMXopGKOXgjF6KRijl4IxeikY\n",
       "o5eCMXopGKOXgjF6KRijl4IxeikYo5eCMXopGKOXgjF6KRijl4KpJ/qlwDZgCHgbuCttXwBsBXYD\n",
       "/cC8IgaU1Fj1RP8F8EvgW8AK4GfARcB68uiXAa+kdUklV0/0+4Ad6fph4B3gXGAN0Je29wE3NHw6\n",
       "SQ033df0HcAlwGvAQmAkbR9J65JKbvY09j0LeBG4G/hkwm1jaZmkp6fn2PUsy8iybFoDSqqtUqlQ\n",
       "qVTq2rfe6E8jD/4ZYFPaNgIsIj/9Xwzsn+oLx0cvqRgTD6i9vb1V963n9L4N2ADsBB4ft/0lYG26\n",
       "vpbjTwaSSqyeI/1K4Dbgn8AbadsDwEPA88A6YBi4uYD5JDVYPdH/g+pnBF0NnEXSSeBv5EnBGL0U\n",
       "jNFLwRi9FMx0fjlHgYyOjjI8PMzo6Oik29rb2+no6GDu3LlNmEwzZfSa0vDwMBs3bmRwcHDSbdde\n",
       "ey033XST0bcoo9eUDhw4wPbt2+nv759023nnncfq1aubMJUawdf0UjBGLwVj9FIwRi8FY/RSMEYv\n",
       "BWP0UjBGLwXjL+doSmeffTYXXHABy5cvn3Rbe3s7c+bMacJUagSj15Q6Ojq48cYbWbFixaTbli1b\n",
       "xvz585swlRrB6DWlJUuWsGTJkmaPoQL4ml4KxuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbo\n",
       "pWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgqm3uhnAW8Am9P6AmArsBvoB+Y1fjRJRag3+ruBncBY\n",
       "Wl9PHv0y4JW0LqkF1BN9O3Ad8BTQlratAfrS9T7ghsaPJqkI9UT/GHAvcHTctoXASLo+ktYltYBa\n",
       "0V8P7Cd/Pd9WZZ8xjp/2Syq5Wn8C+0ryU/nrgK8D3wCeIT+6LwL2AYvJnxim1NPTc+x6lmVkWTaT\n",
       "eSVNoVKpUKlU6tq3VvQPpgWgE/g1cDvwCLAWeDhdbqr2DcZHL6kYEw+ovb29Vfed7uf0X57GPwRc\n",
       "Q/6R3dVpXVILmM7/4WYgLQAHga7GjyOpaP5GnhSM0UvBGL0UjNFLwRi9FIzRS8EYvRSM0UvBGL0U\n",
       "jNFLwRi9FIzRS8EYvRSM0UvBGL0UjNFLwUznj2h8Jf39/UXfRVMMDg4yOjrK4cOHGRoaOmX/naeq\n",
       "yI9ftb9w2yhjq1evLvgumuPAgQMMDw/z+eef09HRQXt7e7NH0jSc6o/fli1boErfhR/p9+zZU/Rd\n",
       "NMXhw4f57LPPOHLkCAcPHuTo0aO1v0ilEfnxK/xI/9RTTxV8F83x5SnhoUOHWLVqFStXrmz2SJqG\n",
       "U/3xu/POO6FZR/p169YVfRdN8fLLLzM0NERbWxtXXXUVd9xxR7NH0jSc6o9fin5KvnsvBWP0UjBG\n",
       "LwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVTb/TzgBeAd4CdwOXA\n",
       "AmArsBvoT/tIKrl6o/8N8BfgIuA7wC5gPXn0y4BX0rqkkqsn+nOAHwC/S+tHgI+BNUBf2tYH3NDw\n",
       "6SQ1XD3Rnw98BDwNDAJPAmcCC4GRtM9IWpdUcvVEPxu4FHgiXX7K5FP5sbRIKrl6/jDm+2l5Pa2/\n",
       "ADwA7AMWpcvFwP6pvrinp+fY9SzLyLLsKw8raWqVSoVKpVLXvvVEvw/YS/6G3W6gCxhKy1rg4XS5\n",
       "aaovHh+9pGJMPKD29vZW3bfeP4H9C+BZ4HTg38BPgFnA88A6YBi4+asMK+nkqjf6N4HvTbG9q4Gz\n",
       "SDoJ/I08KRijl4IxeikYo5eCMXopGKOXgjF6KRijl4IxeikYo5eCMXopGKOXgjF6KRijl4IxeikY\n",
       "o5eCMXopGKOXgjF6KRijl4IxeikYo5eCMXopGKOXgjF6KRijl4IxeikYo5eCMXopGKOXgjF6KRij\n",
       "l4IxeikYo5eCMXopGKOXgqkn+geAIeAt4DlgDrAA2ArsBvqBeUUNKKmxakXfAfwUuBT4NjALuBVY\n",
       "Tx79MuCVtC6pBdSK/n/AF8BcYHa6/BBYA/SlffqAG4oaUFJj1Yr+IPAosIc89kPkR/iFwEjaZySt\n",
       "S2oBs2vc/k3gHvLT/I+BPwO3TdhnLC1T6unpOXY9yzKyLJv+lJJOqFKpUKlU6tq3VvSXAa8Co2l9\n",
       "I3AFsA9YlC4XA/urfYPx0UsqxsQDam9vb9V9a53e7wJWAGcAbUAXsBPYDKxN+6wFNn3laSWdVLWO\n",
       "9G8Cvwe2A0eBQeC3wNnA88A6YBi4ubgRJTVSregBHknLeAfJj/qSWoy/kScFY/RSMEYvBWP0UjBG\n",
       "LwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYv\n",
       "BWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVTePSVSqXo\n",
       "u2i4Vpu51eaF1pu51eY9EaOfQqvN3GrzQuvN3Grznoin91IwRi8F01bw968AnQXfh6TJBoCs2UNI\n",
       "kiRJalmrgV3Au8D9TZ6lmqXANmAIeBu4K21fAGwFdgP9wLymTFfdLOANYHNaL/u884AXgHeAncDl\n",
       "lH/mB8h/Lt4CngPmUP6Zm2oW8B7QAZwG7AAuauZAVSwCvpuunwX8i3zOR4D70vb7gYdO/mgn9Cvg\n",
       "WeCltF72efuAO9L12cA5lHvmDuA/5KED/AlYS7lnbrorgC3j1tenpew2AV3kZygL07ZFab0s2oG/\n",
       "AT/k+JG+zPOeQx7QRGWeeQH5AWA++ZPUZuAayj1z3Yr6nP5cYO+49ffTtjLrAC4BXiN/YEfS9hGO\n",
       "P9Bl8BhwL3B03LYyz3s+8BHwNDAIPAmcSblnPgg8CuwBPgQOkZ/Wl3nmuhUV/VhB37coZwEvAncD\n",
       "n0y4bYzy/HuuB/aTv56v9jsWZZoX8iPlpcAT6fJTJp/1lW3mbwL3kB8IlpD/fNw2YZ+yzVy3oqL/\n",
       "gPxNsi8tJT/al9Fp5ME/Q356D/mz+KJ0fTF5aGVwJbAG+C/wB+Bq8rnLOi/kj/v7wOtp/QXy+PdR\n",
       "3pkvA14FRoEjwEbyl6xlnrluRUW/HbiQ/JnydOAWjr/pVCZtwAbyd5QfH7f9JfI3bkiXmyiHB8mf\n",
       "QM8HbgX+DtxOeeeFPJS9wLK03kX+rvhmyjvzLmAFcAb5z0gX+c9ImWcuhWvJ3wx5j/zjjzK6ivy1\n",
       "8Q7yU+Y3yD9qXED+ZlmZP5rp5PgTadnnXU5+pH+T/Kh5DuWf+T6Of2TXR35GWPaZJUmSJEmSJEmS\n",
       "JEmSJDXT/wEdb549J7fEZwAAAABJRU5ErkJggg==\n"
      ],
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10daeb3d0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#remove a hole from middle of one, should be 5\n",
    "randMatrix[25:30,50:53]=0\n",
    "plt.imshow(randMatrix,cmap='gray')\n",
    "calculate_euler_characteristic(randMatrix)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x10e3065d0>"
      ]
     },
     "execution_count": 46,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": [
       "iVBORw0KGgoAAAANSUhEUgAAAP0AAAD/CAYAAAA6/dD3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\n",
       "AAALEgAACxIB0t1+/AAAC+NJREFUeJzt3VuMVIUdgPFvi0qLWpEXLruE0aYkprRWYytemj01lKAx\n",
       "1EBibaIhFfvUVm1TFX3Z3Tc1MdoXH6rUbI29WCVEkhaXWoakMTGSVWsXKbXtBC9hV4bQuipGwvbh\n",
       "HGGzO8PMunOYM/y/X3Kyc86c3fnj8M05c3EBSZIkSZIkSZIkSR1sDbAX+CdwT5tnkZSzOcCbQAk4\n",
       "E3gVuKidA0lq7HOz+N5vkkZfAT4Bfgd8twUzScrRbKLvBt6atP52tk1Sgc0m+olGO/T29k5k+7m4\n",
       "uJzCJWuvptlE/w6wdNL6UtKj/XG7du2it7eXvr4++vr62LlzJxMTE4Vf+vr6Gu6zfft2Vq9ezYoV\n",
       "K9i8eXPh5y3a0u6ZZ3r/tXveRsvOnTuPd9bX18euXbvqhnvGLKLfDXyZ9IW8d4HvAd+fulOSJPT3\n",
       "98/iZiQ1kiQJSZIcXx8YGKi772yiPwr8GHie9JX8zcAbs/h5kk6B2UQP8KdsqWvyo0+n6LSZO21e\n",
       "6LyZO23ek5nNc/qmdOJ/rE6budPmhc6budPmPZnco5dULEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYv\n",
       "BWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8F\n",
       "Y/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwzUS/FNgJjAB/B27Pti8AdgD7\n",
       "gCFgfh4DSmqtZqL/BPgp8BVgJfAj4CJgE2n0y4EXsnVJBddM9AeAV7PL48AbQDewFhjMtg8CN7R8\n",
       "OkktN9Pn9CXgEuAlYCEwmm0fzdYlFdwZM9j3HOBZ4A7g/SnXTWTLNP39/ccvJ0lCkiQzGlBSY+Vy\n",
       "mXK53NS+zUZ/JmnwTwJbs22jwCLS0//FwFitb5wcvaR8TD2gDgwM1N23mdP7LmAzsAd4ZNL254AN\n",
       "2eUNnHgwkFRgzRzprwJuBv4GvJJtuxe4H3ga2AhUgBtzmE9SizUT/V+pf0awqoWzSDoF/ESeFIzR\n",
       "S8EYvRSM0UvBzOTDOQqkWq1SqVSoVqvTruvp6aFUKjFv3rw2TKbZMnrVVKlU2LJlC8PDw9Ouu/ba\n",
       "a1m/fr3RdyijV00HDx5k9+7dDA0NTbtu2bJlrFmzpg1TqRV8Ti8FY/RSMEYvBWP0UjBGLwVj9FIw\n",
       "Ri8FY/RSMH44RzWde+65XHjhhVx88cXTruvp6WHu3LltmEqtYPSqqVQqsW7dOlauXDntuuXLl3P+\n",
       "+ee3YSq1gtGrpiVLlrBkyZJ2j6Ec+JxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG6KVg\n",
       "jF4KxuilYIxeCsbopWCMXgrG6KVgjF4Kptno5wCvANuy9QXADmAfMATMb/1okvLQbPR3AHuAiWx9\n",
       "E2n0y4EXsnVJHaCZ6HuA64DHga5s21pgMLs8CNzQ+tEk5aGZ6B8G7gKOTdq2EBjNLo9m65I6QKPo\n",
       "rwfGSJ/Pd9XZZ4ITp/2SCq7Rr8C+kvRU/jrg88AXgSdJj+6LgAPAYtIHhpr6+/uPX06ShCRJZjOv\n",
       "pBrK5TLlcrmpfRtFf1+2APQCPwduAR4ENgAPZF+31vsBk6OXlI+pB9SBgYG6+870ffpPT+PvB75D\n",
       "+pbdNdm6pA4wk3/hZle2ABwCVrV+HEl58xN5UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIw\n",
       "Ri8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FM5NfovGZDA0N5X0TbTE8PEy1WmV8fJyRkZHT9s95uop8\n",
       "/9X7DbetMrFmzZqcb6I9Dh48SKVS4ciRI5RKJXp6eto9kmbgdL//tm/fDnX6zv1Iv3///rxvoi3G\n",
       "x8f56KOPOHr0KIcOHeLYsWONv0mFEfn+y/1I//jjj+d8E+3x6Snh4cOHWb16NVdddVW7R9IMnO73\n",
       "32233QbtOtJv3Lgx75toi+eff56RkRG6urq4+uqrufXWW9s9kmbgdL//suhr8tV7KRijl4IxeikY\n",
       "o5eCMXopGKOXgjF6KRijl4IxeikYo5eCMXopGKOXgjF6KRijl4IxeimYZqOfDzwDvAHsAS4HFgA7\n",
       "gH3AULaPpIJrNvpfAH8ELgK+BuwFNpFGvxx4IVuXVHDNRH8e8C3gV9n6UeC/wFpgMNs2CNzQ8ukk\n",
       "tVwz0V8AvAc8AQwDjwFnAwuB0Wyf0WxdUsE1E/0ZwKXAo9nXD5h+Kj+RLZIKrplfjPl2trycrT8D\n",
       "3AscABZlXxcDY7W+ub+///jlJElIkuQzDyuptnK5TLlcbmrfZqI/ALxF+oLdPmAVMJItG4AHsq9b\n",
       "a33z5Ogl5WPqAXVgYKDuvs3+CuyfAE8BZwH/An4AzAGeBjYCFeDGzzKspFOr2ehfA75RY/uqFs4i\n",
       "6RTwE3lSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjBG\n",
       "LwVj9FIwRi8FY/RSMEYvBWP0UjBGLwVj9FIwRi8FY/RSMEYvBWP0UjDN/lt2CqZarVKpVPj4449Z\n",
       "tmwZ3d3d7R5JLWL0qqlSqbBlyxaq1Srr1q0z+tOI0aumI0eOMDY2xtjYGOPj4+0eRy1k9KqpVCqx\n",
       "fv16PvzwQ1asWNHucdRCRq+auru7PaU/TfnqvRSM0UvBGL0UjNFLwTQT/b3ACPA68BtgLrAA2AHs\n",
       "A4aA+XkNKKm1GkVfAn4IXAp8FZgD3ARsIo1+OfBCti6pAzSK/n/AJ8A80rf35gHvAmuBwWyfQeCG\n",
       "vAaU1FqNoj8EPATsJ439MOkRfiEwmu0zmq1L6gCNPpzzJeBO0tP8/wJ/AG6ess9EttTU399//HKS\n",
       "JCRJMvMpJZ1UuVymXC43tW+j6C8DXgSq2foW4ArgALAo+7oYGKv3AyZHLykfUw+oAwMDdfdtdHq/\n",
       "F1gJfAHoAlYBe4BtwIZsnw3A1s88raRTqtGR/jXg18Bu4BgwDPwSOBd4GtgIVIAb8xtRUis18z/c\n",
       "PJgtkx0iPepL6jB+Ik8KxuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG6KVgjF4K\n",
       "xuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG6KVgjF4KxuilYIxeCsbopWCMXgrG\n",
       "6KVgjF4KxuilYIxeCsbopWCMXgom9+jL5XLeN9FynTZzp80LnTdzp817MkZfQ6fN3GnzQufN3Gnz\n",
       "noyn91IwRi8F05Xzzy8DvTnfhqTpdgFJu4eQJEmS1LHWAHuBfwL3tHmWepYCO4ER4O/A7dn2BcAO\n",
       "YB8wBMxvy3T1zQFeAbZl60Wfdz7wDPAGsAe4nOLPfC/p34vXgd8Acyn+zG01B3gTKAFnAq8CF7Vz\n",
       "oDoWAV/PLp8D/IN0zgeBu7Pt9wD3n/rRTupnwFPAc9l60ecdBG7NLp8BnEexZy4B/yYNHeD3wAaK\n",
       "PXPbXQFsn7S+KVuKbiuwivQMZWG2bVG2XhQ9wJ+Bb3PiSF/kec8jDWiqIs+8gPQAcD7pg9Q24DsU\n",
       "e+am5fU+fTfw1qT1t7NtRVYCLgFeIr1jR7Pto5y4o4vgYeAu4NikbUWe9wLgPeAJYBh4DDibYs98\n",
       "CHgI2A+8CxwmPa0v8sxNyyv6iZx+bl7OAZ4F7gDen3LdBMX581wPjJE+n6/3GYsizQvpkfJS4NHs\n",
       "6wdMP+sr2sxfAu4kPRAsIf37cfOUfYo2c9Pyiv4d0hfJPrWU9GhfRGeSBv8k6ek9pI/ii7LLi0lD\n",
       "K4IrgbXAf4DfAteQzl3UeSG9398GXs7WnyGN/wDFnfky4EWgChwFtpA+ZS3yzE3LK/rdwJdJHynP\n",
       "Ar7HiRediqQL2Ez6ivIjk7Y/R/rCDdnXrRTDfaQPoBcANwF/AW6huPNCGspbwPJsfRXpq+LbKO7M\n",
       "e4GVwBdI/46sIv07UuSZC+Fa0hdD3iR9+6OIriZ9bvwq6SnzK6RvNS4gfbGsyG/N9HLigbTo815M\n",
       "eqR/jfSoeR7Fn/luTrxlN0h6Rlj0mSVJkiRJkiRJkiRJkiS10/8BwgbAoCs6o+8AAAAASUVORK5C\n",
       "YII=\n"
      ],
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10dad7e50>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#remove a single pixel hole\n",
    "randMatrix[75,50]=0\n",
    "plt.imshow(randMatrix,cmap='gray')"
   ]
  }
 ],
 "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.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}