pipeline_sandbox.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Testing framework \"pipeline\"\n",
    "\n",
    "_Alex Malz (NYU)_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import string\n",
    "import itertools\n",
    "import random\n",
    "import numpy as np\n",
    "import scipy.stats as sct\n",
    "import scipy.integrate as spi\n",
    "import sklearn as skl\n",
    "from sklearn import metrics\n",
    "from pycm import ConfusionMatrix\n",
    "import pandas as pd\n",
    "import os\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import proclam\n",
    "# from proclam import *"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. True classes"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Import the `proclam.Simulator` superclass and the simulator subclass you want to test.  In this notebook, I'm going to use an unbalanced distribution of true classes such that the probability of an object being in class $m$ (with $0 \\leq m \\leq M$) is proportional to $10^{y}$, where $y$ is a draw from a uniform distribution $U(0,M)$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "from proclam.simulators import simulator\n",
    "from proclam.simulators import logunbalanced"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create an instantiation of the simulator for the true dataset class distribution scheme.  If you use the base superclass instead of a subclass, the default scheme will be to assign all objects the same class."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = proclam.simulators.logunbalanced.LogUnbalanced()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then, simulate a truth catalog.  In this case, there are 3 true classes and 100 objects in the catalog.  The output will be a `numpy.ndarray` with 100 entries, each of which is the index of the class for that catalog member."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = proclam.simulators.logunbalanced.LogUnbalanced(seed=None)\n",
    "M_classes = 5\n",
    "N_objects = 1000\n",
    "names = [''.join(random.sample(string.ascii_lowercase, 2)) for i in range(M_classes)]\n",
    "truth = A.simulate(M_classes, N_objects, base=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can check that the class distribution is as expected with a histogram of the true classes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0.5, 0, 'class')"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZkAAAEKCAYAAADAVygjAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAElpJREFUeJzt3XuwXWV9xvHvQwARgaiAbQXGoDBYpEoEqdWhTvEWpFzGAnJRqjDQzAjqTO2YVrx10FZr6whSaChIVe4qMZRYYJRyK8VcRCRE2og4hEG5FQQ0XH/9Y++M2+M5YR/Oec8+Z5/vZ2ZP9nrXe9b6rdlJnvOud+21UlVIktTCJoMuQJI0vAwZSVIzhowkqRlDRpLUjCEjSWrGkJEkNWPISJKaMWQkSc0YMpKkZjYddAGDtt1229W8efMGXYYkzSgrV668v6q2f7Z+sz5k5s2bx4oVKwZdhiTNKEl+2k8/T5dJkpoxZCRJzRgykqRmDBlJUjOGjCSpGUNGktSMISNJasaQkSQ1M+u/jClJAPMWXT5q+51bHNV2v+vPb7r9sdz59wdMyX4cyUiSmjFkJEnNGDKSpGack5H0W8aanxivYZ3PUP8cyUiSmjFkJEnNGDKSpGYMGUlSM4aMJKkZQ0aS1IwhI0lqxpCRJDUzlF/GTHIIcACwDXB2VV054JIkaVZqNpJJskWS7yX5QZLVST41gW2dk+TeJLeOsm5BktuTrE2yCKCqllTV8cBC4F3P/SgkSRPR8nTZ48B+VfUaYE9gQZLX93ZI8pIkW49o22WUbZ0LLBjZmGQOcDqwP7A7cGSS3Xu6nNxdL0kagGYhUx2Pdhc3675qRLc3AUuSPA8gyfHAaaNs61rgwVF2sw+wtqruqKongAuBg9PxWeDbVbVqco5IkjReTSf+k8xJcjNwL3BVVd3Uu76qLgGuAC5KcjRwLHDYOHaxA3BXz/K6bttJwFuAQ5MsHKO2A5Msfvjhh8exO0nSeDQNmap6uqr2BHYE9kmyxyh9PgesB84ADuoZ/Uxkv6dW1V5VtbCqzhyjz2VVdcLcuXMnujtJ0him5BLmqnoIuJrR51X2BfYALgU+Mc5N3w3s1LO8Y7dNkjQNtLy6bPskL+y+fz7wVuBHI/rMBxYDBwPvA7ZNcso4drMc2DXJzkk2B44Alk5G/ZKkiWs5kvk94Ookt9AJg6uq6t9H9NkSOLyqflxVzwDHAD8duaEkFwA3ArslWZfkOICqego4kc68zhrg4qpa3eyIJEnj0uzLmFV1CzD/WfrcMGL5SeCsUfoduZFtLAOWPccyJUkNeVsZSVIzhowkqRlDRpLUjCEjSWrGkJEkNWPISJKaMWQkSc0YMpKkZgwZSVIzhowkqRlDRpLUjCEjSWrGkJEkNWPISJKaMWQkSc0YMpKkZgwZSVIzhowkqRlDRpLUjCEjSWrGkJEkNWPISJKaMWQkSc0YMpKkZgwZSVIzhowkqRlDRpLUjCEjSWrGkJEkNWPISJKaMWQkSc0YMpKkZjYddAEtJDkEOADYBji7qq4ccEmSNCs1G8kk2SnJ1UluS7I6yQcnsK1zktyb5NZR1i1IcnuStUkWAVTVkqo6HlgIvOu5H4UkaSJajmSeAv6yqlYl2RpYmeSqqrptQ4ckLwF+VVWP9LTtUlVrR2zrXOBLwFd6G5PMAU4H3gqsA5YnWdqzj5O766XnbN6iywey3zu3OKr5PuatP7/5PjS7NRvJVNU9VbWq+/4RYA2ww4hubwKWJHkeQJLjgdNG2da1wIOj7GYfYG1V3VFVTwAXAgen47PAtzfUIEmaelMyJ5NkHjAfuKm3vaouSbIzcFGSS4Bj6YxK+rUDcFfP8jrgD4GTgLcAc7sjozNHqelA4MBddtllHLuTJI1H86vLkmwFfAP4UFX9YuT6qvocsB44Azioqh6d6D6r6tSq2quqFo4WMN0+l1XVCXPnzp3o7iRJY2gaMkk2oxMw51XVN8fosy+wB3Ap8Ilx7uJuYKee5R27bZKkaaDl1WUBzgbWVNU/jdFnPrAYOBh4H7BtklPGsZvlwK5Jdk6yOXAEsHRilUuSJkvLkcwbgfcA+yW5uft6x4g+WwKHV9WPq+oZ4BjgpyM3lOQC4EZgtyTrkhwHUFVPAScCV9C5sODiqlrd7pAkSePRbOK/qq4H8ix9bhix/CRw1ij9jtzINpYBy55jmZKkhrytjCSpGUNGktSMISNJasaQkSQ1Y8hIkpoxZCRJzRgykqRmDBlJUjOGjCSpGUNGktSMISNJasaQkSQ1Y8hIkpoxZCRJzRgykqRmDBlJUjOGjCSpGUNGktSMISNJasaQkSQ101fIJPlgkm3ScXaSVUne1ro4SdLM1u9I5tiq+gXwNuBFwHuAv29WlSRpKPQbMun++Q7gq1W1uqdNkqRR9RsyK5NcSSdkrkiyNfBMu7IkScNg0z77HQfsCdxRVb9Msi3wvnZlSZKGQb8jmauqalVVPQRQVQ8AX2hXliRpGGx0JJNkC2BLYLskL+LX8zDbADs0rk2SNMM92+myvwA+BLwUWMmvQ+YXwJca1iVJGgIbDZmq+iLwxSQnVdVpU1STJGlI9DXxX1WnJXkDMK/3Z6rqK43qkiQNgb5CJslXgVcANwNPd5sLMGQkSWPq9xLmvYHdq6paFiNJGi79XsJ8K/C7LQuRJA2ffkcy2wG3Jfke8PiGxqo6qElVkqSh0G/IfLJlEZKk4dTv1WXXtC5EkjR8+r267BE6V5MBbA5sBjxWVdu0KkySNPP1O5LZesP7JAEOBl7fqihJ0nAY9+OXq2MJ8PYG9UiShki/p8ve2bO4CZ3vzaxvUpEkaWj0e3XZgT3vnwLupHPKTJKkMfU7J+MDyiRJ49bXnEySHZNcmuTe7usbSXZsXZwkaWbrd+L/y8BSOs+VeSlwWbdNkqQx9Rsy21fVl6vqqe7rXGD7hnVJkoZAvyHzQJJ3J5nTfb0beKBlYZKkma/fkDkWOBz4GXAPcCjw3kY1SZKGRL+XMP8t8OdV9X8ASV4MfJ5O+EiSNKp+RzKv3hAwAFX1IDC/TUmSpGHRb8hskuRFGxa6I5l+R0GSpFmq36D4R+DGJJd0lw8DPt2mJEnSsOj3G/9fSbIC2K/b9M6quq1dWZKkYdD3Ka9uqBgskqS+jftW/5Ik9cuQkSQ1Y8hIkpoxZCRJzRgykqRmDBlJUjOGjCSpGUNGktSMISNJasaQkSQ1M1R3Uk5yCHAAsA1wdlVdOeCSJGlWm/YjmSTnJLk3ya0j2hckuT3J2iSLAKpqSVUdDywE3jWIeiVJvzbtQwY4F1jQ25BkDnA6sD+wO3Bkkt17upzcXS9JGqBpHzJVdS3w4IjmfYC1VXVHVT0BXAgcnI7PAt+uqlVTXask6TdN+5AZww7AXT3L67ptJwFvAQ5NsnCsH05yQpIVSVbcd999bSuVpFlsqCb+q+pU4NQ++i0GFgPsvffe1bouSZqtZupI5m5gp57lHbttkqRpZKaGzHJg1yQ7J9kcOAJYOuCaJEkjTPuQSXIBcCOwW5J1SY6rqqeAE4ErgDXAxVW1epB1SpJ+27Sfk6mqI8doXwYsm+JyJEnjMO1HMpKkmcuQkSQ1Y8hIkpoxZCRJzRgykqRmZm3IJDkwyeKHH3540KVI0tCatSFTVZdV1Qlz584ddCmSNLRmbchIktozZCRJzRgykqRmDBlJUjOGjCSpGUNGktSMISNJasaQkSQ1Y8hIkpoxZCRJzczakPHeZZLU3qwNGe9dJkntzdqQkSS1Z8hIkpoxZCRJzRgykqRmDBlJUjOGjCSpGUNGktSMISNJasaQkSQ1Y8hIkpoxZCRJzRgykqRmZm3IeBdmSWpv1oaMd2GWpPZmbchIktozZCRJzRgykqRmDBlJUjOGjCSpGUNGktSMISNJasaQkSQ1Y8hIkpoxZCRJzQxlyCR5eZKzk3x90LVI0mzWNGSSvDDJ15P8KMmaJH/0HLdzTpJ7k9w6yroFSW5PsjbJIoCquqOqjpto/ZKkiWk9kvki8B9V9UrgNcCa3pVJXpJk6xFtu4yynXOBBSMbk8wBTgf2B3YHjkyy++SULkmaqGYhk2Qu8MfA2QBV9URVPTSi25uAJUme1/2Z44HTRm6rqq4FHhxlN/sAa7sjlyeAC4GDJ+8oJEkT0XIkszNwH/DlJN9P8q9JXtDboaouAa4ALkpyNHAscNg49rEDcFfP8jpghyTbJjkTmJ/kr0f7QZ8nI0nttQyZTYHXAmdU1XzgMWDRyE5V9TlgPXAGcFBVPTrRHVfVA1W1sKpeUVV/N0YfnycjSY21DJl1wLqquqm7/HU6ofMbkuwL7AFcCnxinPu4G9ipZ3nHbpskaRpoFjJV9TPgriS7dZveDNzW2yfJfGAxnXmU9wHbJjllHLtZDuyaZOckmwNHAEsnXLwkaVK0vrrsJOC8JLcAewKfGbF+S+DwqvpxVT0DHAP8dORGklwA3AjslmRdkuMAquop4EQ68zprgIuranWzo5EkjcumLTdeVTcDe29k/Q0jlp8Ezhql35Eb2cYyYNkEypQkNTKU3/iXJE0PhowkqRlDRpLUTNM5GQ2feYsun5Tt3LnFUZOynY2Zt/785vuQtHGOZCRJzRgykqRmDBlJUjPOyUzAZM1PPBet5zScz5A0GRzJSJKaMWQkSc0YMpKkZgwZSVIzhowkqZlZGzI+flmS2pu1IePjlyWpvVkbMpKk9gwZSVIzqapB1zBQSe5jlEc+zwDbAfcPuogpNNuOFzzm2WKmHvPLqmr7Z+s060NmpkqyoqrGfLT1sJltxwse82wx7Mfs6TJJUjOGjCSpGUNm5lo86AKm2Gw7XvCYZ4uhPmbnZCRJzTiSkSQ1Y8jMUEkeHXQNmhrD/lkn+WSSDw+6jqmQZF6SWwddx1QyZCRJzRgyM0CSdyf5XpKbk/xLkjnd9i8kWZ3kO0me9UtRM0GSjyW5Pcn1SS5I8uEk/9k91hVJ1iR5XZJvJvnfJKcMuuaJSvJXST7Qff+FJN/tvt8vyXk9/bZLcmOSAwZV62RJ8tEk/5PkemC3bttQf8495iQ5q/tv98okz0/ygSS3JbklyYWDLnAyGTLTXJLfB94FvLGq9gSeBo4GXgCsqKpXAdcAnxhclZMjyeuAPwNeA+wP9H5B7YnuF9bOBL4FvB/YA3hvkm2nutZJdh2wb/f93sBWSTbrtl0LkOR3gMuBj1fV5QOpcpIk2Qs4AtgTeAfwup7Vw/w5b7ArcHr33+5DdP7OLwLmV9WrgYWDLG6ybTroAvSs3gzsBSxPAvB84F7gGeCibp+vAd8cSHWT643At6pqPbA+yWU965Z2//whsLqq7gFIcgewE/DAlFY6uVYCeyXZBngcWEUnbPYFPgBsBnwHeH9VXTOwKifPvsClVfVLgCRLe9YN8+e8wU+q6ubu+5XAPOAW4LwkS4AlgyqsBUcy01+Af6uqPbuv3arqk6P0G/Zr0R/v/vlMz/sNyzP6l6WqehL4CfBe4L/ojGz+BNgFWAM8Rec/o7cPqMSpNLSfc4/e43qaznEdAJwOvJbOL5TDcqyGzAzwHeDQJC8BSPLiJC+j89kd2u1zFHD9gOqbTDcABybZIslWwJ8OuqApdB3wYTqnx66jc8rk+9X5IlsBxwKvTPKRwZU4aa4FDunORWwNHDjoggZsE2Cnqroa+AgwF9hqsCVNnqFJy2FVVbclORm4MskmwJN0zlM/BuzTXXcvnXmbGa2qlndPndwC/JzOKZPZ8ujS64CPAjdW1WNJ1nfbAKiqp5McCSxN8khV/fOgCp2oqlqV5CLgB3T+7i4fcEmDNgf4WpK5dM5cnFpVDw24pknjN/41rSTZqqoeTbIlnd94T6iqVYOuS9Jz40hG083iJLsDW9CZizJgpBnMkYwkqRkn/iVJzRgykqRmDBlJUjOGjDRAs+kOxJqdDBlJUjOGjDSFkhzTvdPuD5J8dcS645Ms7677Rve7QiQ5LMmt3fYNN8x8Vc+duW9Jsusgjkd6Nl7CLE2RJK8CLgXeUFX3J3kxnRtgPlpVn0+ybVU90O17CvDzqjotyQ+BBVV1d5IXVtVDSU4D/ruqzkuyOTCnqn41qGOTxuJIRpo6+wGXVNX9AFX14Ij1eyS5rhsqRwOv6rbfAJyb5Hg6tyABuBH4m+69zF5mwGi6MmSk6eNc4MSq+gPgU3TuekBVLQROpnOr+5XdEc/5wEHAr4BlSfYbTMnSxhky0tT5LnDYhodvdU+X9doauKf7wLKjNzQmeUVV3VRVHwfuA3ZK8nLgjqo6lc7DvV49JUcgjZP3LpOmSFWtTvJp4JokTwPfB+7s6fIx4CY6QXITndAB+IfuxH7oPPrhB3RuCf+eJE8CPwM+MyUHIY2TE/+SpGY8XSZJasaQkSQ1Y8hIkpoxZCRJzRgykqRmDBlJUjOGjCSpGUNGktTM/wPlywg9cu2opQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "d = np.diff(np.unique(truth)).min()\n",
    "left_of_first_bin = truth.min() - float(d)/2\n",
    "right_of_last_bin = truth.max() + float(d)/2\n",
    "plt.hist(truth, np.arange(left_of_first_bin, right_of_last_bin + d, d), log=True)\n",
    "plt.xticks(range(max(truth)+1), names)\n",
    "\n",
    "plt.hist(truth, log=True)\n",
    "plt.ylabel('counts')\n",
    "plt.xlabel('class')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "def save_info(truth, df, name, seed=0):\n",
    "    pathname = \"examples/ProClaM/\"+name+\"/\"\n",
    "    if os.path.isdir(pathname) == False:\n",
    "        os.system('mkdir '+ pathname)\n",
    "    df.to_csv(pathname + \"predicted_prob_\"+name+\".csv\", sep=' ', index=False, header=False)\n",
    "    truth_mask = proclam.metrics.util.det_to_prob(truth)\n",
    "    fileloc = pathname + 'truth_table_'+name+'.csv'\n",
    "    np.savetxt(fileloc, truth_mask, delimiter=' ')\n",
    "    return"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Mock classifications"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Import the `proclam.Classifier` superclass and the classifier subclass you want to test.  In this notebook, I'm going to use a very stupid classifier that takes a random guess."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "from proclam.classifiers import classifier\n",
    "from proclam.classifiers import guess\n",
    "from proclam.classifiers import from_cm"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create an instantiation of the classifier of a particular scheme.  If you use the base superclass instead of a subclass, the default classification scheme will return the true classes.\n",
    "\n",
    "Then, \"classify\" the \"data.\"  By default, classification results will also include an extra column for \"other\" classes beyond the number in the training set, but in this example let's assume it knows of the M classes in the training set without leaving room for additional classes.  The output will be a `numpy.ndarray` with N rows and column entries representing each catalog member's posterior probability for being of each class."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "B = proclam.classifiers.guess.Guess(seed=None)\n",
    "predictionB = B.classify(M_classes, truth, other=False)\n",
    "save_info(truth, pd.DataFrame(predictionB), 'Guess', seed=0)\n",
    "# print(prediction)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's compare this to a smarter classifier based on a confusion matrix."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[0.77459445 0.11209536 0.06052605 0.00362192 0.07104932]\n",
      " [0.11003839 0.80009249 0.01965666 0.02884268 0.01506941]\n",
      " [0.09383029 0.07687496 0.65189991 0.07650355 0.12909478]\n",
      " [0.12409745 0.0805176  0.00174395 0.75464134 0.06014014]\n",
      " [0.06169615 0.08333079 0.04914953 0.02018145 0.77150318]]\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'true class')"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARYAAAEMCAYAAAABAJmyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAEr5JREFUeJzt3XuQXnV9x/H3ZzeJuRJu0REbCXIVKESzUdo0HS5W5KYooXLTUtRUTQ1tJ4qO3NrRqY7UTJlBS6gdqdyiFCGC3CYQciGWbALkwkUqgdqWEkzkkoQku8m3f5yT+mTZhfMkv3POPruf18zOc57znOf8PtlsPvmd8zz7HEUEZmYptdUdwMwGHheLmSXnYjGz5FwsZpaci8XMknOxmFlyLpYGkjbWnWEw6O/fZ0lXSppVd44iJE2QtLruHD25WMwsuUFbLJIukPSIpMckXSupPV8/W9IaSfMljas7J4CkyyQ9LWmxpJslzZK0IM/aKelJSZMl3SbpGUnfqDnvlyXNzJdnS3ogXz5R0o0N2+0vaamk0+rK2pDl65J+KWkxcHi+rt9+j3tol3Rd/nN7n6QRkmZKekLSSkm3VB1oUBaLpPcCnwSmRMREYDtwPjAK6IyIo4CHgCvqS5mRNBk4CzgWOAXoaHh4W0R0AP8E3AHMAI4GLpS0X9VZGywCpubLHcBoSUPzdQsBJL0DuAu4PCLuqiVlTtIk4BxgInAqMLnh4f76PW50KHBN/nP7MtnPy1eB90XEMcDnqw40pOoB+4mTgEnAMkkAI4B1wA5gbr7NDcBttaTb1RTgjojYAmyR9LOGx+blt6uANRHxAoCkZ4HxwPpKk/7OcmCSpL2ArcAKsoKZCswEhgLzgRkR8VBNGRtNBX4aEZsBJM1reKy/fo8brY2Ix/Ll5cAEYCVwo6TbgdurDjQoZyyAgOsjYmL+dXhEXNnLdv39F6m25rc7GpZ33q/tP42I6ALWAhcCD5PNYE4ADgGeBLrJ/gGcXFPEZvTL73EPjbm2k+U6DbgGeD/Zf6CVZh2sxTIfmCbp7QCS9pV0INn3Y1q+zXnA4pryNVoCnCFpuKTRwOl1BypoETCL7NBnEdl0/NHIfus1gIuAIyRdUl/E/7cQODM/NzEGOKPuQHuoDRgfEQ8ClwBjgdFVBugvjVupiHhC0qXAfZLagC6yY+dNwAfyx9aRnYepVUQsy6fmK4EXyabkr9SbqpBFwNeBpRGxSdKWfB0AEbFd0rnAPEmvRcT36goaESskzQUeJ/t7X1ZXlkTagRskjSWbnV8dES9XGUD+2IT+T9LoiNgoaSTZ/67TI2JF3bnM+jIoZywtaI6kI4HhZOeGXCrWr3nGYmbJDdaTt2ZWIheLmSXnYjGz5FwsDSRNrztDM1otL7Re5lbLC/0js4tlV7X/hTSp1fJC62VutbzQDzK7WMwsuZZ7uXm/fdti/Phy3n6zfv0O9tsvfdeufWb/5PsE2Na9mWFDRqbfcXd3+n3mtu14nWFtI5LvN7rKydzFVobytlL2rWFDS9nvtu2vM6w9/ff49e5X2bb9dRXZtuXeIDd+/BDu+3k5/1DL8qnTP1t3hKa0rftt3RGa1v3C/9YdoWlDDhhfd4SmPPw/N771RjkfCplZci4WM0vOxWJmyblYzCw5F4uZJediMbPkXCxmlpyLxcySc7GYWXIuFjNLzsViZsm5WMwsOReLmSXnYjGz5FwsZpaci8XMknOxmFlytRaLpI11jm9m5fCMxcySq6xYJF0g6RFJj0m6VlJ7vn62pDWS5ksaV1UeMytPJcUi6b3AJ4EpETER2A6cD4wCOiPiKOAh4Io+nj9dUqekzvXrd1QR2cz2QFWf0n8SMAlYJglgBLAO2AHMzbe5AbittydHxBxgDsDEY4e11vVKzAahqopFwPUR8bVdVkqX9djOpWE2AFR1jmU+ME3S2wEk7SvpwHz8afk25wGLK8pjZiWqZMYSEU9IuhS4T1Ib0AXMADYBH8gfW0d2HsbMWlxlV0KMiLn87nzKTqOrGt/MquP3sZhZci4WM0vOxWJmyblYzCw5F4uZJediMbPkXCxmlpyLxcySc7GYWXIuFjNLzsViZsm5WMwsOReLmSXnYjGz5FwsZpaci8XMklNEa33M7NiRB8Rxh32m7hhNufueW+qO0JSPHPTBuiM0TUMq+8yyZGJbV90RmvKLrnt4dcd6FdnWMxYzS87FYmbJuVjMLDkXi5kl52Ixs+RcLGaWnIvFzJJzsZhZci4WM0vOxWJmyblYzCw5F4uZJediMbPkXCxmlpyLxcySc7GYWXIuFjNLzsViZsm5WMwsOReLmSVX6icQS7oMuAB4Cfg1sBw4HXgUmAqMAj4NfA34fWBuRFxaZiYzK19pxSJpMnAWcCwwFFhBViwA2yKiQ9LFwB3AJGAD8CtJsyNifVm5zKx8ZR4KTQHuiIgtEfEa8LOGx+blt6uANRHxQkRsBZ4FxvfckaTpkjoldW7r3lRiZDNLoa5zLFvz2x0Nyzvvv2EWFRFzIqIjIjqGDRlVRT4z2wNlFssS4AxJwyWNJju3YmaDQGnnWCJimaR5wErgRbLDnlfKGs/M+o+yr0t5VURcKWkksBBYHhHX7XwwIhYACxruH19yHjOrQNnFMkfSkcBw4PqIWFHyeGbWD5RaLBFxXpn7N7P+ye+8NbPkXCxmlpyLxcySc7GYWXIuFjNLzsViZsm5WMwsOReLmSXnYjGz5FwsZpaci8XMknOxmFlyLhYzS87FYmbJNVUsktok7VVWGDMbGN7y81gk3QR8HtgOLAP2kvSPEfGdssP1KgJt7a5l6N112pSP1R2hKX+87Km6IzRtyUcPrztC07aPG1t3hOasXlB40yIzliMj4lXgTOBu4CDgU7sVzMwGhSLFMlTSULJimRcRXUCUG8vMWlmRYrkWeI7scqgLJR0IvFpmKDNrbW95jiUirgaublj1vKQTyotkZq3uLWcski6WtJcyP5C0Ajixgmxm1qKKHApdlJ+8/TCwD9mJ22+VmsrMWlqRYlF+eyrwo4hY07DOzOwNihTLckn3kRXLvZLGkF283cysV0UuWPYZYCLwbERslrQf8OflxjKzVlbkVaEdktYCh0kaXkEmM2txRd7S/1ngYuD3gMeA44Cl+JUhM+tDkXMsFwOTgecj4gTgfcDLpaYys5ZWpFi2RMQWAElvi4ingNb7jS8zq0yRk7f/JWlv4Hbgfkm/BZ4vN5aZtbIiJ28/ni9eKelBYCxwT6mpzKyl9VkskvbtZfWq/HY0sKGURGbW8t5sxrKc7OMRGt9lu/N+AO8pMZeZtbA+iyUiDqoyiJkNHEV+u/njksY23N9b0pmpAkjamGpfZtY/FHm5+YqIeGXnnYh4GbiivEhm1uqKFEtv2xR5mRoASV+WNDNfni3pgXz5REk3Nmy3v6Slkk4rum8z65+KFEunpO9KOjj/+i7Zid2iFgFT8+UOYHT+GbpTgYUAkt4B3AVcHhF39dyBpOmSOiV1buve3MTQZlaHIsXyJWAbMBe4BdgCzGhijOXApPx6RFvJfs+og6xYFgFDgfnAVyLi/t52EBFzIqIjIjqGDRnZxNBmVocib5DbBHx1dweIiK78t6MvBB4GVgInAIcATwLdZOVzMvDQ7o5jZv1HVZdYXQTMIjv0WUR2AbRHIyLI3hNzEXCEpEsqymNmJaqyWN4JLI2IF8kOpxbtfDAitgPnAidK+mJFmcysJIVf3dkTETGf7FzKzvuHNSyPzm+3kh0OmVmLK/IGucMkzZe0Or9/jKRLy49mZq2qyKHQdcDXgC6AiFgJnFNmKDNrbUWKZWREPNJjXXcZYcxsYChSLL+RdDD5heAlTQNeKDWVmbW0IidvZwBzyF4O/m9gLXBBqanMrKUVeYPcs8CHJI0C2iLitfJjmVkrK3L5j8t73AcgIv6upExm1uKKHAptalgeDpxO9lZ8M7NeFTkU+ofG+5KuAu4tLZGZtbzdeUv/SLKrIpqZ9arIOZZV5C81A+3AOMDnV8ysT0XOsZzesNwNvBgRfoOcmfXpTYtFUjtwb0QcUVEeMxsA3vQcS/5xBk9LendFecxsAChyKLQPsEbSIzS89BwRHy0tlZm1tCLFclnpKcxsQClSLKdGxC4fGSnp29T0+bTRJnaMGV7H0LutfUNrnet+6JgRdUdo2veev6nuCE2bedy0uiM0RV1dhbct8j6WP+ll3SmFRzCzQafPGYukLwBfBN4jaWXDQ2OAJWUHM7PW9WaHQjcBdwN/z66X/3gtIjaUmsrMWlqfxZJfr/kVsk/PNzMrrKrLf5jZIOJiMbPkXCxmlpyLxcySc7GYWXIuFjNLzsViZsm5WMwsOReLmSXnYjGz5FwsZpaci8XMknOxmFlyLhYzS66yYpF0paRZVY1nZvXxjMXMkiu1WCR9XdIvJS0GDs/XLZA0W1KnpCclTZZ0m6RnJH2jzDxmVo0in9K/WyRNAs4BJubjrACW5w9vi4gOSRcDdwCTgA3AryTNjoj1ZeUys/KVOWOZCvw0IjZHxKvAvIbHdi6vAtZExAsRsRV4Fhjfc0eSpucznM6u7s0lRjazFOo6x7I1v93RsLzz/htmURExJyI6IqJj6JCRVeQzsz1QZrEsBM6UNELSGOCMEscys36ktHMsEbFC0lzgcWAdsKysscysfymtWAAi4pvAN3usvqrh8QXAgob7x5eZx8yq4fexmFlyLhYzS87FYmbJuVjMLDkXi5kl52Ixs+RcLGaWnIvFzJJzsZhZci4WM0vOxWJmyblYzCw5F4uZJediMbPkXCxmlpyLxcySc7GYWXKlfoJcGdS1nfYXNtQdoykx4m11R2hK+7hxdUdo2ozDTqo7QtP+9T/+re4ITTn51JcLb+sZi5kl52Ixs+RcLGaWnIvFzJJzsZhZci4WM0vOxWJmyblYzCw5F4uZJediMbPkXCxmlpyLxcySc7GYWXIuFjNLzsViZsm5WMwsOReLmSXnYjGz5CorFkkTJK2uajwzq49nLGaWXNXF0i7pOklrJN0naYSkmZKekLRS0i0V5zGzElT9Kf2HAudGxOck/Rg4C/gqcFBEbJW0d29PkjQdmA4wvH1MZWHNbPdUPWNZGxGP5cvLgQnASuBGSRcA3b09KSLmRERHRHQMaxtRTVIz221VF8vWhuXtZDOm04BrgPcDyyS13LWOzGxXdZ+8bQPGR8SDwCXAWGB0vZHMbE/VPTtoB26QNBYQcHVEFL/cmpn1S5UVS0Q8BxzdcP+qqsY2s2rVfShkZgOQi8XMknOxmFlyLhYzS87FYmbJuVjMLDkXi5kl52Ixs+RcLGaWnIvFzJJzsZhZci4WM0vOxWJmyblYzCw5F4uZJediMbPkXCxmlpwiou4MTZH0EvB8SbvfH/hNSfsuQ6vlhdbL3Gp5obzMB0bEuCIbtlyxlElSZ0R01J2jqFbLC62XudXyQv/I7EMhM0vOxWJmyblYdjWn7gBNqjSvpI357QGSbn2Lbf9K0sheHuozs6TjJd3ZRJ4Fksqe8rfazwT0g8w+xzLISWqPiO0Ft90YEYUuKCfpOaAjIgqfRJR0PDArIk4vuP2CfPvOomNYNTxjGaAkTZD0lKQbJT0p6dadMwhJz0n6tqQVwNmSDpZ0j6TlkhZJOiLf7iBJSyWtkvSNHvtenS+3S7pK0mpJKyV9SdJM4ADgQUkP5tt9ON/XCkk/kTQ6X/+RPOcK4BN9/FneMEYv23xfUqekNZL+tmH9tyQ9kT/vqnzd2fm+Hpe0MM133HYREf4agF/ABCCAKfn9fyH73x3gOeArDdvOBw7Nlz8IPJAvzwM+nS/PADY27Ht1vvwF4FZgSH5/34Yx9s+X9wcWAqPy+5cAlwPDgV8Dh5JdCfPHwJ29/Fn6GmMB2ayocV17vv4YYD/gaX43M987v10FvKtxnb/SfnnGMrD9OiKW5Ms3AH/U8NhcgHzm8IfATyQ9BlwLvDPfZgpwc778oz7G+BBwbUR0A0TEhl62OQ44EliSj/FnwIHAEcDaiHgmsn/lN+zBGH+az3oeBY7Kx3sF2AL8QNIngM35tkuAH0r6HFkRWWJ1X7vZytXzBFrj/U35bRvwckRMLLiP3SHg/og4d5eVUl9jNrdz6SBgFjA5In4r6YfA8IjolvQB4CRgGvCXwIkR8XlJHwROA5ZLmhQR61NksYxnLAPbuyX9Qb58HrC45wYR8SqwVtLZAMocmz+8BDgnXz6/jzHuB/5C0pD8+fvm618DxuTLvwCmSDok32aUpMOAp4AJkg7Ot9uleAqMsdNeZEX5iqR3AKfk240GxkbEz4G/Bo7N1x8cEf8eEZcDLwHj+xjXdpOLZWB7Gpgh6UlgH+D7fWx3PvAZSY8Da4CP5esvzp+/CnhXH8/9Z+A/gZX588/L188B7pH0YES8BFwI3CxpJbAUOCIitgDTgbvyw5h1TY4BQEQ8TnYI9BRwE1khQlZsd+ZjLgb+Jl//nfyE9GrgYeDxPsa13eSXmwcoSRPIToQeXXMUG4Q8YzGz5DxjMbPkPGMxs+RcLGaWnIvFzJJzsZhZci4WM0vOxWJmyf0fFeFaBj1BNuAAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 288x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "cm = np.eye(M_classes) + 0.2 * np.random.uniform(size=(M_classes, M_classes))\n",
    "cm /= np.sum(cm, axis=1)\n",
    "print(cm)\n",
    "np.flip(cm,0)\n",
    "plt.matshow(cm)#, vmin=0., vmax=1.)\n",
    "plt.xticks(range(max(truth)+1), names)\n",
    "plt.yticks(range(max(truth)+1), names)\n",
    "plt.xlabel('predicted class')\n",
    "plt.ylabel('true class')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We want to support classifiers with an extra class for classes not represented in the training set, but the infrastructure isn't there yet."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[0.11284647 0.80101817 0.03185446 0.03936846 0.01491243]\n",
      " [0.09973749 0.07725771 0.60755048 0.08775271 0.12770162]\n",
      " [0.08582221 0.08195349 0.05213019 0.02797892 0.75211519]\n",
      " ...\n",
      " [0.07873401 0.08532457 0.05993838 0.03413578 0.74186726]\n",
      " [0.05587737 0.46613713 0.36092296 0.04309275 0.0739698 ]\n",
      " [0.09950859 0.07750803 0.61382232 0.07556721 0.13359384]]\n"
     ]
    }
   ],
   "source": [
    "C = proclam.classifiers.from_cm.FromCM(seed=None)\n",
    "predictionC = C.classify(cm, truth, other=False)\n",
    "print(predictionC)\n",
    "save_info(truth, pd.DataFrame(predictionC), 'Idealized', seed=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also try this with real data from `SNPhotCC`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "# old_snphotcc_names = []\n",
    "# for prefix in ['templates_', 'wavelets_']:\n",
    "#     for suffix in ['boost_forest', 'knn', 'nb', 'neural_network', 'svm']:\n",
    "#         old_snphotcc_names.append(prefix+suffix+'.dat')\n",
    "\n",
    "# snphotcc = {}\n",
    "# snphotcc['label'] = 'SNPhotCC'\n",
    "# prefixes = ['Templates', 'Wavelets']\n",
    "# suffixes = ['BoostForest', 'KNN', 'NB', 'NeuralNetwork', 'SVM']\n",
    "# snphotcc['names'] = []\n",
    "# for prefix in prefixes:\n",
    "#     for suffix in suffixes:\n",
    "#         snphotcc['names'].append(prefix+suffix)\n",
    "# snphotcc['dirname'] = 'examples/' + snphotcc['label'] + '/'\n",
    "\n",
    "# truthloc = 'examples/SNPhotCC/key.txt'\n",
    "# truth_snphotcc = pd.read_csv(truthloc, sep=' ')\n",
    "# # truth_snphotcc -= 1.\n",
    "# # truth_snphotcc = truth_snphotcc.astype(int)[:N_objects]\n",
    "# # truth_snphotcc = pd.DataFrame(truth_snphotcc)\n",
    "\n",
    "# for i in range(len(snphotcc['names'])):\n",
    "#     name = old_snphotcc_names[i]\n",
    "#     fileloc = 'examples/SNPhotCC/classifications/'+name\n",
    "#     snphotcc_info = pd.read_csv(fileloc, sep=' ')\n",
    "#     full = snphotcc_info.set_index('Object').join(truth_snphotcc.set_index('Object'))\n",
    "#     name = snphotcc['names'][i]\n",
    "    \n",
    "#     snphotcc_truth = full['Type'] - 1\n",
    "#     snphotcc_truth_table = proclam.metrics.util.det_to_prob(snphotcc_truth)\n",
    "#     fileloc = snphotcc['dirname']+name+'/truth_table_'+name+'.csv'\n",
    "#     np.savetxt(fileloc, snphotcc_truth_table, delimiter=' ')\n",
    "    \n",
    "#     probs = full[['1', '2', '3']]\n",
    "#     fileloc = snphotcc['dirname']+name+'/predicted_prob_'+name+'.csv'\n",
    "#     probs.to_csv(fileloc, sep=' ', index=False, header=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "# print(snphotcc_info)\n",
    "# print(truth_snphotcc)\n",
    "# print(full)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "# probs = np.loadtxt('examples/classifications/templates_knn.dat', usecols=[1, 2, 3], skiprows=1)\n",
    "# # print(np.shape(probs))\n",
    "# probs = probs[:N_objects]\n",
    "# # print(np.shape(probs))\n",
    "\n",
    "# # truth_snphotcc = np.loadtxt('examples/key.txt', usecols=[1], skiprows=1)\n",
    "# # truth_snphotcc -= 1.\n",
    "# # truth_snphotcc = truth_snphotcc.astype(int)[:N_objects]\n",
    "# # truth_snphotcc = pd.DataFrame(truth_snphotcc)\n",
    "\n",
    "# CM = proclam.metrics.util.prob_to_cm(probs, truth_snphotcc, vb=True)\n",
    "# # print(CM)\n",
    "# plt.matshow(CM)#, vmin=0., vmax=1.)\n",
    "# plt.xticks(range(max(truth)+1)[:3], names[:3])\n",
    "# plt.yticks(range(max(truth)+1)[:3], names[:3])\n",
    "# plt.xlabel('predicted class')\n",
    "# plt.ylabel('true class')\n",
    "\n",
    "# F = proclam.classifiers.from_cm.FromCM(seed=None)\n",
    "# predictionF = F.classify(CM, truth_snphotcc, other=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Below is tunnel vision classifier which means that it performs well on one class only and classifies other classes agnostically so it might classify other classes correctly by chance. The class that is classified correctly is chosen in the code randomly each time."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There are different scenarios for the confusion matrix which are discussed and shown below:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**1) Tunnel classifier CM:** \n",
    "\n",
    "This is where the confusion matrix has one class that has a high value on the diagnoal and small value on the off-diagonal elements for that specific class and has small values for the diagonal component of other classes or 1/M value for all entries (classifies y chance)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[0.15842418 0.18857113 0.16013153 0.17989965 0.01515237]\n",
      " [0.17176994 0.13661988 0.13941695 0.09105078 0.10789923]\n",
      " [0.17019976 0.06886277 0.77397623 0.01511753 0.06713813]\n",
      " [0.04351717 0.12580389 0.14106496 0.06285216 0.11175521]\n",
      " [0.06891721 0.0174588  0.12293988 0.07999228 0.12880625]]\n"
     ]
    }
   ],
   "source": [
    "# N = 3 #len(truth)\n",
    "# M = len(cm)\n",
    "#if other: M += 1\n",
    "CM = np.zeros((M_classes, M_classes))\n",
    "class_corr = np.array([2]) #np.random.randint(0, M_classes, size=1) # randomly choose which class to work well on\n",
    "\n",
    "for i in range(M_classes):\n",
    "    for j in range(M_classes):\n",
    "        CM[i][j] = np.random.uniform(0., 1./M_classes, 1)\n",
    "        CM[class_corr[0]][class_corr[0]] = np.random.uniform(1./M_classes, 1., 1)\n",
    "        \n",
    "print(CM)\n",
    "CM = CM/np.sum(CM, axis=1)[:, np.newaxis] #normalizing to make sure each row sums to 1\n",
    "CMtunnel = CM\n",
    "\n",
    "# saving probabilities\n",
    "G           = proclam.classifiers.from_cm.FromCM(seed=None)\n",
    "predictionG = G.classify(CMtunnel, truth, other=False)\n",
    "df          = pd.DataFrame(predictionG)\n",
    "save_info(truth, df, 'Tunnel', seed=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x7f4c5f729630>"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARMAAAD0CAYAAAC4n8I2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAGuBJREFUeJzt3X2UXVWd5vHvkxAISQDBqK2AktYgIo0oAe1hmAZ8i6LgC3YDOtOOthlHELpZOOISo8PYq1ubkWnXopXYutr2DdS2NS3R6CDIi4h5EQIJojGoJM0IAUFezFvVM3+cU/FWUVX3VO65de+pej5rnVXnnLvr7B83yY999tlnb9kmIqJTM3odQERMDUkmEVGLJJOIqEWSSUTUIskkImqRZBIRtUgyiYhaJJlERC2STCKmIUmLJd0laaOki0b5/DJJt5bbTyU91PaaGQEb0f9eefJcP/DgQKWya9ZtX2l78VifS5oJ/BR4ObAZWAWcZXvDGOXfDbzQ9tvGq3evStFFRE9tfXCAW1YeUqnsrKf/fH6bIscDG21vApB0JXA6MGoyAc4CPtiu3iSTiEYwAx6sWni+pNUtx8tsL2s5Phi4p+V4M/Di0S4k6VnAAuB77SpNMoloAAODVO6S2Gp7UU1Vnwl81Xbbe6wkk4iGGKRyy6SdLcChLceHlOdGcyZwTpWLJplENIAxA/U9LFkFLJS0gCKJnAmcPbKQpCOAA4Gbq1x02j4alvRor2Pod738jiR9SNKFvaq/JY7DJN3R6ziguM2psrVjexdwLrASuBP4su31ki6RdFpL0TOBK13xkW9aJhENYGCgep9J++vZK4AVI84tHXH8oYlcc1q0TCS9RdKPygE4V5TP2YcG5qyXdI2kp0xCHB8oBwrdKOlLki6UdF0Zx2pJd0o6TtLXJP1M0oe7HM97JJ1X7l8m6Xvl/imSvtBSbr6kmyWd2uV43l8OkLoReG55rmffT4uZkj5V/l35jqR9JZ0naYOkdeWj1a4ysNODlbZemfLJRNLzgD8DTrB9DDAAvBmYC6y2/Xzg+1R4jt5hHMcBbwReALwKaO1t31H2vn8S+AZFh9dRwFslPbmLYd0AnFjuLwLmSZpVnru+jPtpwNXAUttXdysQScdSNKuPAV4NHNfyca++nyELgcvLvysPUfw5XkQxkOto4J2TEAODFbdemQ63OS8FjgVWSQLYF7iP4nu/qizzeeBrXY7jBOAbtrcB2yT9W8tny8uftwPrbd8LIGkTRa/7A12KaQ1wrKT9ge3AWoqkciJwHjALuAY4x/b3uxTDkBOBf7X9OICk5S2f9er7GXK37VvL/TXAYcA64AuSvg58vcv1Fx2wNd7mdMOUb5kAAj5r+5hye+4Y94K9/JPaXv4cbNkfOu5awre9E7gbeCvwA4qWysnAcyg65nZR/ON5ZbdiqKgn388o9UPRst0LOBW4HHgRxf+ouhuHYaDi1ivTIZlcA5wh6akAkg4qR/XNAM4oy5wN3NjlOG4CXitptqR5wGu6XF9VNwAXUtzW3EDRZP9x2YNv4G3AEZLe2+U4rgdeV/ZH7Ae8tsv1dWIGcKjta4H3AgcA87pZYTFoLbc5PWV7g6SLge9ImgHspLjnfgw4vvzsPop+lW7Gsapsuq8Dfk3RZH+4m3VWdAPwfuBm249J2laeA8D2gKSzgOWSHrH9D90IwvZaSVcBt1H8eazqRj01mQl8XtIBFC3fj9tu+1ZtZ8QA6m4VHcpbw5NI0jzbj0qaQ/F/4iW21/Y6ruh/Rx29t//l6nbv7xWOeOa9a2ocTl/ZlG+Z9Jllko4EZlP04ySRRCUGdvR5r0SSySSy/YQhyxFVDbq/b3OSTCIaoBgBm2QSER0yYqDPb3P6O7oukrSk1zG06qd4+ikWSDxDBq1KW69M22QC9NVfUPornn6KBRLP7tucKluv5DYnohHEgPv7//2NSyZ7zZ7rfeYd1PF19p57IHPnH9rxIJvBmR2HAsCs/Q5k36d1Hs+MahOYj6uu76auv/uz5h3InKd2Hs9e9z9WRzjMZg7766CO43mE32y1XeltdQM7qekvW5c0LpnsM+8gnvfav+p1GLttf1J/9bDv83D/DELcOae/vpunfvKWXocwzP8duOqXVcvaaZlERE0G82g4IjpVdMCmZRIRHcttTkTUoJiCIMkkIjpkxA7naU5E1GAwtzkR0akmdMD2d3QRAZQv+rnaVoWkxeWyKxslXTRGmT8tl/NYL+mL7a6ZlklEQ9TVAVuuG3U58HJgM8WE2Mttb2gpsxB4H8USMb8ZmkN5PEkmEQ1gU+ej4eOBjbY3AZSLiJ0ObGgp8w6KtYJ+U9Tv+9pdNLc5EY0gBituwPxyBcShbeRbzgcD97Qcby7PtTocOFzSTZJ+KGlxuwjTMoloAAM7XPmf69YaJpTei2Ilw5OAQ4DrJf3ReLPwJ5lENICpdeKjLRQrIQ45pDzXajNwy9BCbZJ+SpFcxlyCpKe3OZIe7WX9EU0ywIxKWwWrgIWSFkjam2KN5+UjynydolWCpPkUtz2bxrtoWiYRDWDqG7Rme5ekc4GVFAuKfcb2ekmXAKttLy8/e4WkDRRLor7H9rhrOk9aMpH0ForFsPcGbgHeVZ6/DHgF8P+AM23fP1kxRTRHvVMy2l4BrBhxbmnLvoELyq2SSbnNkfQ8iuU3T7B9DEWmezMwlyITPh/4PvDBMX5/yVDP9K5t9cyWFdEkQy2TKluvTFbL5KXAsRSDYwD2pVhPdhC4qizzeeBro/2y7WXAMqCW6QQjmijr5hREsRzm+4adlD4wolwSRcQobLFzsL+7OCerTXQNcMbQkFxJB0l6Vln/GWWZs4EbJymeiEYp5jOpPGitJyYl1dneIOli4DuSZgA7gXOAx4Djy8/uo+hXiYgnyExru9m+it/3jwyZN1n1RzRZ0QGbPpOIqEG/z2eSZBLRADUPp++KJJOIhsiE0hHRMRt2DiaZRESHitucJJOIqEFGwEZEx/JoOCJqktuciKhJL4fKV5FkEtEAxez0SSYR0SEjdg1mreGIqEFuc2rmmbDjgP75UufcN9jrEIbppz66+et+1+sQhtnxihf1OoThvjXyvdex5WlORNQmT3MionPOi34RUYOhmdb6WZJJREP0e8ukv2/CIgIoWia7BmdU2qqQtFjSXZI2SrpolM/fKul+SbeW21+0u2ZaJhENUOfkSJJmApcDL6dYU3iVpOW2N4woepXtc6teNy2TiIaocXb644GNtjfZ3gFcCZzeaXxJJhFN4KLPpMpWwcHAPS3Hm8tzI71R0jpJX5V0aLuLJplENMDQoLWKyWT+0HK65bZkD6r8N+Aw20cD3wU+2+4X0mcS0RAT6DPZanvROJ9vAVpbGoeU53az/UDL4T8CH21XaZJJRAMYMVDfHLCrgIWSFlAkkTMpVtTcTdLTbd9bHp4G3NnuokkmEQ1R16A127sknQusBGYCn7G9XtIlwGrby4HzJJ0G7AIeBN7a7rpJJhENYNc7aM32CmDFiHNLW/bfB7xvItdMMoloCPf5CNgkk4hGyIt+EVGTad0ykfQB4C3A/RSDZNYArwF+DJwIzAX+C8W92R9RDN+9uJsxRTTRtJ4cSdJxwBuBFwCzgLUUyQRgh+1Fks4HvgEcS9Fj/HNJl414xk056GYJwKz9DuxWyBH9qwETSndzBOwJwDdsb7P9CMWIuiHLy5+3A+tt32t7O7CJ4YNpALC9zPYi24v22nduF0OO6E+muM2psvVKr/pMtpc/B1v2h47TjxPxBP3fAdvNlslNwGslzZY0j6KvJCL2kF1t65WutQJsr5K0HFgH/JrilubhbtUXMdVN66c5wKW2PyRpDnA9sMb2p4Y+tH0dcF3L8UldjieikYpWx/ROJsskHQnMBj5re22X64uYsvq9z6SrycT22e1LRUQVg4PTOJlERD1Mbx/7VpFkEtEQPXxQU0mSSUQTpAM2ImrT502TJJOIhphSLRNJM4B5tn/bpXgiYgy9HN1aRdvh9JK+KGl/SXOBO4ANkt7T/dAiYogNHpxRaeuVKjUfWbZEXgd8C1gA/OeuRhURT9Dv7+ZUSSazJM2iSCbLbe+k77uCIqYgV9x6pEoyuQL4BcWsaNdLehaQPpOISVVtLpO+ns/E9seBj7ec+qWkk7sX0vhm/s48ef329gUnyewNW9oXmkRXr/l2r0PY7ZXPOKbXIQwz5w+e1usQOtPn9wNVOmDPLztgJenTktYCp0xCbBExxPXOtCZpsaS7JG2UdNE45d4oyZLGW24UqHab87ayA/YVwIEUna9/WyniiKhPTX0mkmYClwOvAo4Ezirf7h9Zbj/gfOCWKuFVSSZDqe7VwOdsr285FxGTxaq2tXc8sNH2Jts7gCuB00cp97+AjwDbqly0SjJZI+k7FMlkZZmtBqtcPCJqVN/TnIMplp4Zsrk8t5ukFwGH2r66anhVRsC+HTgG2GT7cUlPBv5r1QoiogamaqsDYL6k1S3Hy2wvq/rL5Uj3j1FhsfJWVZ7mDEq6Gzhc0uyJXDwi6jOBAWlbbY/XYbqF4UvKHFKeG7IfcBRwnSSAPwCWSzrNdmuSGqZtMpH0FxSdMIcAtwIvAW4mT3QiJld9j4ZXAQslLaBIImcCu2dFtP0wMH/oWNJ1wIXjJRKo1mdyPnAc8EvbJwMvBB6aaPQR0aGaOmBt7wLOBVYCdwJftr1e0iWSTtvT8Kr0mWyzvU0Skvax/RNJz93TCiNiDxhU42MP2yuAFSPOLR2j7ElVrlklmWyW9CTg68B3Jf0G+GWVi0dEXSo/9u2ZKh2wry93PyTpWuAAoH/GbEdMF30+nH7MZCLpoFFO317+nAc82JWIImJ0TU0mwBqK8FvbVkPHBv6wi3FFxEhNTSa2F0xGAJIetT1vMuqKaKyJDVrriSpvDb9e0gEtx0+S9LruhhURI8nVtl6pMs7kg+UgFgBsPwR8sGoFkt4j6bxy/zJJ3yv3T5H0hZZy8yXdLOnU6uFHTCNTYKa10cpMZFb7G4ATy/1FwLxyGsgTgesBJD0NuBpYOpEXiyKmk6nQMlkt6WOSnl1uH6PonK1qDXCspP2B7RRD8RdRJJMbgFnANcD/sP3d0S4gaYmk1ZJW79z52ASqjphC6puCoCuqJJN3AzuAqyjmPdgGnFO1gnIC6rsp3kD8AUUCORl4DsVQ3l0UCeeV41xjme1FthfNmjW3atURU0fVW5wetkyqDFp7DBhzWreKbgAuBN5GMVblY8Aa25bk8vxXJL3X9kc6rCtiaurzR8OTtWLPDcDTgZtt/5qidXPD0Ie2B4CzgFMkvWuSYopolH7vM5mUtYZtX0PRNzJ0fHjL/rzy53bGudWJmPb6vGWShcsjGkA1vzXcDVUGrR0u6RpJd5THR0u6uPuhRcQwU+BpzqeA9wE7AWyvo5iZKSImU9Of5gBzbP+onAtyyK4uxRMRY+hl52oVVZLJVknPpsx5ks4A7u1qVBHxRFMgmZwDLAOOkLSFYgDaW7oaVUQM1+PHvlVUGbS2CXiZpLnADNuPdD+siHiCpicTSUtHHANg+5IuxRQRo+j3R8NVbnNa36ybDbyG4p2aiIjdqtzm/O/WY0mXUqy3ERGTqem3OaOYQ7G6X0RMlgZ0wFYZAXu7pHXlth64C/g/3Q8tIoapcdCapMWS7pK0UdITZgWQ9M7y3/6tkm6UdGS7a1ZpmbymZX8X8OtyecGemLFjgH3u6Z/VSf/9Df01Sf9Rf98/L10f/JJHex3CMPc9r8/mwvnMBMvX1DKRNBO4HHg5sBlYJWm57Q0txb5o+5Nl+dMopg1ZPN51x00mZaUrbR/RSfAR0RlR623O8cDGctgHkq4ETgd2JxPbv20pP5cKqWzcZGJ7oGwKPdP2r/Yo7Ijo3MTeGp4vaXXL8TLby1qODwbuaTneDLx45EUknQNcAOwNnNKu0iq3OQcC6yX9iJbHxLb3eLX0iNgD1VsmW20v6rg6+3LgcklnAxcDfz5e+SrJ5AOdBhURNajvNmcLcGjL8SHlubFcCXyi3UWrTEHwatvfb92AV1f4vYioUY3TNq4CFkpaIGlviilFlg+rS1rYcngq8LN2F62STF4+yrlXVfi9iKhTTY+Gy6ex51IMPr0T+LLt9ZIuKZ/cAJwrab2kWyn6Tca9xYFxbnMk/XfgXcAfSlrX8tF+wE3tQ46I2tQ88ZHtFcCKEeeWtuyfP9Frjtdn8kXgW8DfMHypi0dsPzjRiiKiM4190a9cX/hhiiUoIqLH+n04fWanj2iKJJOI6FiPJ4uuIskkogFUbv0sySSiKdIyiYg6pAO2JOlDwKO2L52sOiOmlKY+Go6IPjIVZlrrhKT3S/qppBuB55bnrpN0maTVku6UdJykr0n6maQPdzOeiEabAsuD7hFJx1K8QHRMWc9aYE358Q7biySdD3wDOBZ4EPi5pMtsPzDiWkuAJQCz99q/WyFH9LXp3DI5EfhX24+Xsza1vpU4tH87sN72vba3A5sY/mo0ALaX2V5ke9HeM+d0MeSIPjZdWyZtbC9/DrbsDx2nHydiFNO5ZXI98DpJ+0raD3htF+uKmNqqtkqmYsvE9lpJVwG3AfdRTMgSEXtANPit4TrY/mvgr0ecvrTl8+uA61qOT+pmPBGN1ue3OemfiGgIub+zSZJJRBPkreGIqEu/P81JMoloiiSTiKhDWiYR0bmJLQ/aE0kmEU3R5y2Trr41HBH1ELWu6IekxZLukrRR0kWjfH6BpA2S1km6RtKz2l0zySSiKexqWxuSZgKXU6zMeSRwlqQjRxT7MbDI9tHAV4GPtrtukklEQ9TYMjke2Gh7k+0dFAuTn95awPa1th8vD39Isbj5uJJMIppgYi/6zS8nHxvaloy42sHAPS3Hm8tzY3k7xeqe40oHbERDTOBpzlbbi2qpU3oLsAj4k3ZlG5dMBveZye8WHNjrMHZ7+ufu6HUIw/zqXUf1OoTdHn/Gvr0OYZjHDu73lWfGV+Oj4S0Mn4TskPLc8PqklwHvB/6knLxsXLnNiWgCU1sHLMV0IAslLZC0N8X0qq0zISLphcAVwGm276ty0ca1TCKmq7pGwNreJelcYCUwE/iM7fWSLgFW214O/B0wD/iKJIBf2T5tvOsmmUQ0RY2D1myvAFaMOLe0Zf9lE71mkklEAwwNWutnSSYRTVC9P6RnkkwiGiIv+kVELXKbExGdMzDY39kkySSiKfo7lySZRDRFbnMioh55mhMRdUjLJCI6JoPSARsRtejzcSaT9tawpMMk9df7+hENIrvS1itpmUQ0QQOWB53s+UxmSvqUpPWSviNpX0nntcyCfeUkxxPREBXnMplGLZOFwFm23yHpy8AbgYuABba3S3rSaL9UzmG5BGCf2aMWiZjy+v1pzmS3TO62fWu5vwY4DFgHfKGca3LXaL9ke5ntRbYXzdp77uREGtFv+rxlMtnJpHUeyQGKltGpFGt4vAhYJSn9OBEjGTTgSluv9HoO2BnAobavBd4LHEAxVVxEjFR9qYue6HUrYCbweUkHUEwm9XHbD/U4poi+1MvHvlVMWjKx/QvgqJbjSyer7ogpIckkIjpm+n4EbJJJRAOI3o5urSLJJKIp+jyZ9PppTkRUYWDA1bYKJC2WdJekjZIuGuXz/yRpraRdks6ocs0kk4iGqOtFP0kzKcZ2vQo4EjhL0pEjiv0KeCvwxarx5TYnoinqu805HthoexNA+U7c6cCG31flX5SfVe72TcskohFqfdHvYOCeluPN5bmOpGUS0QRmIi2T+ZJWtxwvs72s/qCGSzKJaIrq40y22l40zudbgENbjg8pz3UktzkRDVHjTGurgIWSFkjaGzgTWN5pfEkmEU1gYGCw2tbuUvYu4FxgJXAn8GXb6yVdIuk0AEnHSdoMvAm4QtL6dtfNbU5EI9Q7V4ntFcCKEeeWtuyvorj9qaxxyeTR327Zev23L/plDZeaD2yt4Tp1qSeev+k8EKbqd/MvnQdSquv7edaESvf5CNjGJRPbT6njOpJWt+mkmlT9FE8/xQKJZ7ckk4jomIEswhURnTO4v+cgmM7JpOuDeCaon+Lpp1gg8fz+aU4fm7aPhidjROBE1BmPpEfLn8+Q9NU2Zf9S0pyJxCLpJEnfnEA810na4z6GqfxnNcGKMzt9dK5803NCbP+77Xavj/8lMKdNmegHSSYxnnIN5p9I+oKkOyV9dailIOkXkj4iaS3wJknPlvRtSWsk3SDpiLLcAkk3S7pd0odHXPuOcn+mpEsl3VGunvhuSecBzwCulXRtWe4V5bXWSvqKpHnl+cVlnGuBN4zx3/KEOkYp8wlJq8tVHf9ny/m/bVnZ8dLy3JvKa90m6fp6vvGmyop+Uc1zgbfbvknSZ4B3AUMTbj9g+0UAkq4B3mn7Z5JeDPwDcArw98AnbP+zpHPGqGMJxaJnx9jeJekg2w9KugA42fZWSfOBi4GX2X5M0nuBCyR9FPhUWddG4KqqdYxS5v1lvTOBayQdTfFeyOuBI2y7ZWXHpcArbW8Za7XHacPAYPpMor17bN9U7n8e+I8tn10FULYQ/gPwFUm3AlcATy/LnAB8qdz/3Bh1vAy4ohxKje0HRynzEorJcm4q6/hzioFVR1Csxvgz2y5j3NM6/rRs3fwYeH5Z38PANuDTkt4APF6WvQn4J0nvoFgWZXpLyyQqGPk3oPX4sfLnDOAh28dUvMaeEPBd22cNOymNVefELi4tAC4EjrP9G0n/BMwuWzHHAy8FzqB4b+QU2+8sW2CnAmskHWv7gTpiaaQ+H7SWlkl/eKakPy73zwZuHFnA9m+BuyW9CUCFF5Qf30Tx5ifAm8eo47vAfxtafrXlFuQRYL9y/4fACZKeU5aZK+lw4CfAYZKeXZYblmwq1DFkf4rk+LCkp1FMGzjU6jqgfF/kr4AXlOefbfuW8p2R+xn+2vz0YuOBgUpbrySZ9Ie7gHMk3QkcCHxijHJvBt4u6TZgPcVUewDnl79/O2PPmPWPFPN6rit//+zy/DLg25KutX0/xbyfX5K0DriZoh9jG0V/yNXlLcp9E6wDANu3Udze/IRibtGhW7v9gG+Wdd4IXFCe/7uyU/kO4AfAbWPUOz0MutrWI3KfN52mOkmHAd+0fVSbojGNHbDXU/zH+53eviCw8qFPr+nFu0PpM4loArvvn+YkmfTYyDWYI8bU53cRSSYRDeG0TCKic70dQ1JFkklEExjo4WPfKpJMIhrAgDM5UkR0zJkcKSJq0u8tkwxai2gASd+mmBW/iq22F3czntEkmURELfJuTkTUIskkImqRZBIRtUgyiYhaJJlERC2STCKiFkkmEVGLJJOIqEWSSUTU4v8DXq5Ucvg2FXoAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 288x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.matshow(CMtunnel)\n",
    "plt.xticks(range(max(truth)+1), names)\n",
    "plt.yticks(range(max(truth)+1), names)\n",
    "plt.xlabel('predicted class')\n",
    "plt.ylabel('true class')\n",
    "plt.colorbar()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**2) Broadbrush classifier CM:** \n",
    "\n",
    "This is where the confusion matrix takes thre most populous class and classifies that one correctly. This is very similar to tunnel classifier except that we kow which class it will classify correctly and it is not randomly chosen."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "CM = np.zeros((M_classes, M_classes))\n",
    "class_corr = sct.mode(truth)[0] #takes the most common class\n",
    "\n",
    "for i in range(M_classes):\n",
    "    for j in range(M_classes):\n",
    "        CM[i][j] = np.random.uniform(0., 1./M_classes, 1)\n",
    "        CM[class_corr[0]][class_corr[0]] = np.random.uniform(1./M_classes, 1., 1)\n",
    "        \n",
    "CM = CM/np.sum(CM, axis=1)[:, np.newaxis] #normalizing to make sure each row sums to 1\n",
    "CMbroadbrush = CM\n",
    "\n",
    "# saving probabilities\n",
    "H           = proclam.classifiers.from_cm.FromCM(seed=None)\n",
    "predictionH = H.classify(CMbroadbrush, truth, other=False)\n",
    "df          = pd.DataFrame(predictionH)\n",
    "\n",
    "save_info(truth, df, 'Broadbrush', seed=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x7f4c5f679cf8>"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARMAAAD0CAYAAAC4n8I2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAGgtJREFUeJzt3X+YXmV95/H3h5AYSCBCoxQBJSI/iggoA+pStoKoKApaoQvIblmpKSsUul644KpoWXu1WlZ2vS5qiT9at4qg1koUSmApSEDU/BACCSAxSAlrxYAgYPNr5rN/nDPxmcnMPGfynOfXzOd1Xeeac85zP+d8M5P5zn3f537uW7aJiGjVTt0OICKmhiSTiKhFkklE1CLJJCJqkWQSEbVIMomIWiSZREQtkkwiohZJJhFRi527HUBENPeW4+f4yacGK5VdsWrTEtsntTmk7SSZRPSBDU8N8oMl+1YqO3Pvn8xvczhjSjKJ6Atm0EPdDmJCSSYRfcDAEL39odwkk4g+MURqJhHRImMGe3y6kGn7aFjSc92Oodd183sk6eOSLu7W/Rvi2F/S/d2OA4pmTpWtCkknSXpI0lpJl45T5g8krZG0WtI1za6ZmklEHzAwWFOfiaQZwFXAm4D1wDJJi22vaShzIPAh4Fjbv5T04mbXnRY1E0lnS/qhpHskXV1+M5F0ZZl1b5X0og7E8dHyr8Gdkr4q6WJJt5dxLJf0gKSjJX1T0sOSPtHmeD4o6cJy/0pJ/1zunyDpKw3l5ku6W9LJbY7nw5J+LOlO4ODyXNe+Pw1mSPpc+X/lZkm7SLqw/Ku9StK17Q7AwBYPVdoqOAZYa3ud7c3AtcCpo8q8D7jK9i8BbD/R7KJTPplI+h3gP1Bk2COBQeA9wBxgue1XAt8FPtbmOI4G3g0cAbwVGGh4ebPtAeBvgOuB84HDgHMk/VYbw1oKHFfuDwBzJc0sz91Rxr0XcANwme0b2hWIpKOAM4AjgbcBRze83K3vz7ADKX6xXgk8TfFzvBR4te3DgfM6EANDFTdgfpl8h7eFoy61D/BYw/H68lyjg4CDJN0l6fuSmg6Cmw7NnDcCR1FU5QB2AZ6g+L5fV5b5MvDNNsdxLHC97Y3ARknfbnhtcfn1PmC17Z8BSFoH7Ac82aaYVgBHSdod2ASspEgqxwEXAjOBW4HzbX+3TTEMOw74R9u/BpC0uOG1bn1/hj1i+55yfwWwP7AK+IqkbwHfavP9iw7Y6s2cDWXybcXOFEn0DcC+wB2SXmX76fHeMOVrJoCAL9k+stwOtv3xMcp1s6t8U/l1qGF/+LhtCd/2FuAR4BzgexQ1leOBVwAPAFspfnne0q4YKurK92eM+0NRs90ZOJmi3+E1FH+o2huHYbDiVsHjFEl42L7luUbrgcW2t9h+BPgxRXIZ13RIJrcCpw13IEnaU9LLKP7tp5VlzgLubHMcdwHvkDRb0lzg7W2+X1VLgYspmjVLKarsP3KxbIGB9wKHSLqkzXHcAbyz7I/YDXhHm+/Xip2A/WzfBlwCzAPmtvOGxaC1ys2cZpYBB0paIGkWRfNy8agy36KolSBpPkWzZ91EF53yzRzbayR9BLhZ0k7AFoo29/PAMeVrT1D0q7QzjmVl1X0V8HOKKvsz7bxnRUuBDwN3235e0sbyHAC2ByWdCSyW9Kztv25HELZXSroOuJfi57GsHfepyQzgy5LmUdR8PzNR9b8eYhDVciXbWyVdACyh+Ld80fZqSZdT9CMuLl97s6Q1FLWxD9qesDmprJvTOZLm2n5O0q4Uf4kX2l7Z7bii9x12+Cz/ww3VPr93yEt/tqKGPpNJm/I1kx6zSNKhwGyKfpwkkqjEwOYe75VIMukg22d1O4boX0Oup5nTLkkmEX2gGAGbZBIRLTJisMebOb0dXRuNMSqwq3opnl6KBRLPsCGr0tYt0zaZAD31H5TeiqeXYoHEs62ZU2XrljRzIvqCGHRv/+3vu2QySy/wbOa0fJ3Z7Mru2rP1QTa7zm75EgCzZ81j9zkvaTmerXNmtBzLzLl7sOuL92s5lrr+78/cbQ92+e3W45k5b3Md4TB7r93Y/eC9Wo7n2R8/scF2pU+rG9hC6z/bduq7ZDKbObxWb+x2GNvo0Fd2O4QRfv66ed0OYZvBF3Q7gpH2Pvlfuh3CCP/3hP/1aNWydmomEVGToTwajohWFR2wqZlERMvSzImIGhRTECSZRESLjNjsPM2JiBoMpZkTEa1KB2xE1MKIwUxBEBF1SAdsRLTMJo+GI6IOygjYiGidgc3u7V/X3o4uIoCiAzZzwE5A0nO227p4UcRUkUfDEdEy0/uD1joWnaSzJf1Q0j2SrpY0ozx/paTVkm6VVGmimIjpp9qUjd2ctrEjyUTS71Asv3ms7SMplht8DzCHYjnCVwLfBT42zvsXSlouafmWEWtIR0wPwzWTKlu3dOrObwSOolgt/p7y+OUU6yxfV5b5MvC7Y73Z9iLbA7YHZtJj03dFdEidNRNJJ0l6SNJaSZeO8fo5kn5RtiTukfRHza7ZqT4TUSyH+aERJ6WPjiqXhY8jxmCLLUP1/LqWXQxXAW8C1lP8kV9se82ootfZvqDqdTtVM7kVOE3SiwEk7SnpZeX9TyvLnAXc2aF4IvpKMZ+JKm0VHAOstb3O9mbgWuDUVmPsSDIpM95HgJslrQJuAfYGngeOkXQ/cAJweSfiieg/xUxrVbYK9gEeazheX54b7d2SVkn6hqT9ml20Y4+GbV/Hb/pHhmWMSUQFRQds5Sc18yUtbzheZHvRJG/5beCrtjdJ+mPgSxR/8MeVcSYRfWISg9Y22B6Y4PXHgcaaxr7luW1sP9lw+HngU81u2tujYCIC+M1w+prWGl4GHChpgaRZwBnA4sYCkvZuODwFeKDZRVMziegTdc1nYnurpAuAJcAM4Iu2V0u6nGLc12LgQkmnAFuBp4Bzml03ySSiD9iwZai+hoTtG4EbR527rGH/Q8CHRr9vIkkmEX2gaOb0dq9EkklEn+jm526qSDKJ6AOTfDTcFUkmEX0hzZyIqEnmgI2IlhWz0yeZRESLjNg6lLWGI6IGaebUTDvvzIz5L+52GNs8cmlvdYrtu+djzQt1yKxze+s/v2+Y3e0Qdlie5kREbfI0JyJaV/1DfF2TZBLRB4ZnWutlSSYRfSI1k4homYGtNX5quB2STCL6QNYajojapM8kIlrn9JlERA0yaC0iapNkEhEtM2IwT3Miog7pgI2IljkdsBFRFyeZRETrMmgtImoyrWsmkj4KnA38AngMWAG8HfgRcBwwB/hPFCuHvQq4zvZH2hlTRD+a1uNMJB0NvBs4ApgJrKRIJgCbbQ9Iugi4HjiKYj3Tn0i6ctQK7EhaCCwEmL3T3HaFHNG7+mBC6XY+uD4WuN72RtvPAt9ueG14xfX7gNW2f2Z7E7AO2G/0hWwvsj1ge2DWTru0MeSI3mSKZk6VrQpJJ0l6SNJaSZdOUO7dkixpoNk1u9Vnsqn8OtSwP3ycfpyI7dTXAStpBnAV8CZgPbBM0mLba0aV2w24CPhBleu2s2ZyF/AOSbMlzaXoK4mIHWRX2yo4Blhre53tzcC1wKljlPsfwCeBjVUu2rZkYnsZRXNmFfBPFE2aZ9p1v4iprsZmzj4UD0SGrS/PbSPpNcB+tm+oGl+7mxRX2P64pF2BO4AVtj83/KLt24HbG47f0OZ4IvpSUeuo3MyZL2l5w/Ei24uqvlnSTsCngXOqR9j+ZLJI0qHAbOBLtle2+X4RU9Yk+kw22J6ow/RxRj7o2Lc8N2w34DDgdkkAvw0slnSK7cYkNUJbk4nts9p5/YjpZGiotkfDy4ADJS2gSCJnANt+V20/A8wfPpZ0O3DxRIkE2tsBGxE1MdX6S6o0hWxvBS4AlgAPAF+zvVrS5ZJO2dEY8xg2ok9Ue1BT8Vr2jcCNo85dNk7ZN1S5ZpJJRD+YXAdsVySZRPSLOqsmbZBkEtEnplTNpHz+PNf2r9oUT0SMo+Lo1q5p+jRH0jWSdpc0B7gfWCPpg+0PLSKG2eChnSpt3VLlzoeWNZF3UgyLXwD8x7ZGFRHbqfGzOW1RJZnMlDSTIpkstr2Fnu8KipiCXHHrkirJ5GrgpxSzot0h6WVA+kwiOqq+QWvt0rQD1vZngM80nHpU0vHtC2liW/Z4Af/6rgO6dfvt/NkR13Q7hBH+9qjDux3CNoObN3c7hBE2nP2abocw0oOTLN/j7YEqHbAXlR2wkvQFSSuBEzoQW0QMc70zrbVDlWbOe8sO2DcDe1B0vv5lW6OKiO31eJ9JlXEmw6nubcDflx8I6u3RMxFTUY8PWqtSM1kh6WaKZLKknBdyqL1hRcR2pkDN5FzgSGCd7V9L+i3gP7c3rIgYwfR8zaTK05whSY8AB0ma3YGYImIMvT6cvmkykfRHFNPd7wvcA7wOuJs80YnorB5PJlX6TC4CjgYetX088Grg6bZGFRHbs6ptXVKlz2Sj7Y2SkPQC2w9KOrjtkUXEbxjU4489qiST9ZJeCHwLuEXSL4FH2xtWRIzU3VpHFVU6YN9V7n5c0m3APOCmtkYVEdvr8T6TcZOJpD3HOH1f+XUu8FRbIoqIsfVrMgFWUITfWLcaPjbw8jbGFRGj9Wsysb2gEwFIes723E7cK6Jv9cGgtSqfGn6XpHkNxy+U9M72hhURo8nVtm6pMs7kY+VygQDYfhr4WNUbSPqgpAvL/Ssl/XO5f4KkrzSUmy/pbkknVw8/Yhrp8c/mVEkmY5WZzKz2S4Hjyv0BYG45DeRxwB0AkvYCbgAus33DJK4dMW1MhZrJckmflnRAuX2aonO2qhXAUZJ2BzZRDMUfoEgmS4GZwK3Af7N9y1gXkLRQ0nJJy7f+2/OTuHXEFNLjI2CrJJM/ATYD1wHXAhuB86veoJyA+hHgHOB7FAnkeOAVFIsmb6VIOG+Z4BqLbA/YHth5lzlVbx0xdVRt4lSsmUg6SdJDktZKunSM18+TdJ+keyTdKenQZtdsmkxsP2/70vKX+Wjb/932ZKsHS4GLKZo1S4HzgB/ZHv7nvxc4RNIlk7xuxPRRUzKRNAO4CngrcChw5hjJ4hrbr7J9JPAp4NPNrtupFXuWAnsDd9v+OUXtZunwi7YHgTOBEyS9v0MxRfSVGvtMjgHW2l5nezNFi+PUxgKjVu2cQ4U01ZG1hm3fStE3Mnx8UMP+3PLrJiZo6kRMe/V1ru4DPNZwvB547ehCks4HPgDMosKUI91bSzAiKlP5qeEqGzB/+IFFuS3ckXvavsr2AcAlwEeala8yOdJBwGeBvWwfJulw4BTbn9iRACNiB1V/UrPB9sAErz8O7NdwvG95bjzXUuSACVWpmXwO+BCwBcD2KuCMCu+LiDrV9zRnGXCgpAWSZlH8Pi9uLCDpwIbDk4GHm120Sp/JrrZ/OGp1i60V3hcRNaprQJrtrZIuAJYAM4AvlkvYXA4st70YuEDSiRSViF8Cf9jsulWSyQZJB1DmPEmnAT/bwX9HROyoGke32r4RuHHUucsa9i+a7DWrJJPzgUUU40AepxiAdvZkbxQRLejyUPkqqsy0tg44UdIcYCfbz7Y/rIjYTr8nE0mXjToGwPblbYopIsYwFSaUbhw6Pxt4O8VnaiIitqnSzPmfjceSrqDoBY6ITur3Zs4YdqUY5BIRnTIVOmAl3cdvcuIM4EVA+ksiOq3fkwlFH8mwrcDPbXd10Jp36p2Jdf/mgtO7HcIIz542s3mhDtn8wt75OQHMfK7Hfxub6fHwJ0wm5bwHS2wf0qF4ImIMovebORN+NqecZ+QhSS/tUDwRMZbJfWq4K6o0c/YAVkv6IQ2PiW2f0raoImJ7PV4zqZJMPtr2KCKiuSmQTN5me8TcrJI+CXy3PSFFxFj6us+k9KYxzr217kAiookeX4Rr3JqJpP8CvB94uaRVDS/tBtzV7sAiokGXE0UVEzVzrgH+CfgLoHFdjWdtP9XWqCJiO337Qb9yfeFnKJagiIgu6/U+k44sdRERNUgyiYiW9XmfSUT0CJVbL0syiegXqZlERB3SAVuS9HHgOdtXdOqeEVNKvz4ajoge0gczrbV14XJJH5b0Y0l3AgeX526XdGW5oPIDko6W9E1JD0vK+sUR4+nX4fStknQUxRqmR5b3WQmsKF/ebHtA0kXA9cBRwFPATyRdafvJUddaCCwEmDl3j3aFHNHTpnPN5DjgH23/2vavGLkw8vD+fcBq2z+zvQlYx8jV2QGwvcj2gO2BnXeZ08aQI3pYj9dM2trMmcCm8utQw/7wcfpxIsYgV9sqXUs6SdJDktZKunSM1z8gaY2kVZJulfSyZtdsZzK5A3inpF0k7Qa8o433ipjaqtZKKiSTcm7nqyimEjkUOFPSoaOK/QgYsH048A3gU82u27ZkYnslcB1wL8Wnj5e1614RU52odQ7YY4C1ttfZ3gxcC5zaWMD2bbZ/XR5+nwprZbW1SWH7z4E/H3X6iobXbwdubzh+Qzvjiehr1ftD5kta3nC8yPaihuN9gMcajtcDr53geudSVAgmlP6JiD4hV84mG2wP1HJP6WxgAPi9ZmWTTCL6Qb1Pah5n5FPTfctzI0g6Efgw8Hvl09YJdetpTkRMUo1Pc5YBB0paIGkWxXiwxqEbSHo1cDVwiu0nqlw0ySSiX9T0NKdc3vcCYAnwAPA126slXS5peD2svwLmAl+XdI+kxeNcbps0cyL6RJ0jYG3fCNw46txlDfsnTvaaSSYR/cB9PKF0RPSYHv9sTpJJRB8Qvf9BvySTiH5RfZxJVySZRPSJ1EwionVZ6iIi6pKnOTWb+exW9r7tF90OY5uN+83rdggjvOgfVnc7hG3+7fUHdTuEEW772893O4QRZixqXqZRkklEtM6kAzYi6pEO2IioR5JJRLQqg9Yioh52+kwioh55mhMRtUgzJyJaZ2Cot7NJkklEv+jtXJJkEtEv0syJiHrkaU5E1CE1k4homQxKB2xE1KLHx5l0bN0cSftLur9T94uYamRX2rolNZOIftAHM611ekW/GZI+J2m1pJsl7SLpQklrJK2SdG2H44noE/7N53OabV3S6ZrJgcCZtt8n6WvAu4FLgQW2N0l64VhvkrQQWAgwe+buHQs2opf0+tOcTtdMHrF9T7m/AtgfWAV8RdLZwNax3mR7ke0B2wOzZuzamUgjek2P10w6nUw2NewPUtSMTgauAl4DLJOUfpyI0QwadKWtCkknSXpI0lpJl47x+r+XtFLSVkmnVblmp5PJWPffz/ZtwCXAPIqV1yNiNFfcmpA0g+IP+FuBQ4EzJR06qti/AOcA11QNr9u1gBnAlyXNo5hM6jO2n+5yTBE9qcbHvscAa22vAygffJwKrBkuYPun5WuVR7d0LJmUwR3WcHxFp+4dMSVUTybzJS1vOF5ku3FhjX2AxxqO1wOvbTG6rtdMIqIKM5kRsBtsD7QvmLElmUT0AVHr6NbHgf0ajvctz7Wk2x2wEVFVfY+GlwEHSlogaRZwBrC41fCSTCL6gYFBV9uaXcreClwALAEeAL5me7WkyyWdAiDpaEnrgdOBqyU1XXc2zZyIPlHnh/hs3wjcOOrcZQ37yyiaP5UlmUT0i8y0FhGtyyJcEVEHk2QSETXp8ZnWkkwi+kQ3Z1GrIskkoh8YGOztqkmSSURfSAds7X618V83LFnzF4/WcKn5wIaWr7KmeZGK6omnHvXEclPrgZRqiWfG3jVEUqjrZ/WySZVOMqmX7RfVcR1Jy7vxYajx9FI8vRQLJJ5tkkwiomUGsghXRLTO4HTA9qpFzYt0VC/F00uxQOLpi6c50/ZTw6Nmnuq6OuOR9Fz59SWSvtGk7J9KGjHlf7NYJL1B0ncmEc/tkna4j2Eq/6wmeePMTh+tKycBnhTb/892s5nF/xTI+iH9IMkkJlKuwfygpK9IekDSN4ZrCpJ+KumTklYCp0s6QNJNklZIWirpkLLcAkl3S7pP0idGXfv+cn+GpCsk3V+unvgnki4EXgLcJum2styby2utlPR1SXPL8yeVca4Efn+cf8t29xijzGclLS9XdfyzhvN/2bCy4xXludPLa90r6Y56vuP9Kiv6RTUHA+favkvSF4H3A8MTbj9p+zUAkm4FzrP9sKTXAn8NnAD8b+Cztv+PpPPHucdCikXPjrS9VdKetp+S9AHgeNsbJM0HPgKcaPt5SZcAH5D0KeBz5b3WAtdVvccYZT5c3ncGcKukwymmDHwXcIhtN6zseBnwFtuPj7fa47RhYCh9JtHcY7bvKve/DPxuw2vXAZQ1hH8HfF3SPcDVwPAwrGOBr5b7fz/OPU4Eri5n2cL2U2OUeR3FOip3lff4Q4qBVYdQrMb4sG2XMe7oPf6grN38CHhleb9ngI3AFyT9PvDrsuxdwN9Jeh/FsijTW2omUcHo/wGNx8+XX3cCnrZ9ZMVr7AgBt9g+c8RJabx7Tu7i0gLgYuBo27+U9HfA7LIWcwzwRuA0iikFT7B9XlkDOxlYIeko20/WEUtf6vFBa6mZ9IaXSnp9uX8WcOfoArZ/BTwi6XQAFY4oX76LYlJggPeMc49bgD8eXn61oQnyLLBbuf994FhJryjLzJF0EPAgsL+kA8pyI5JNhXsM250iOT4jaS+KFeWGa13zyqkE/ytwRHn+ANs/KKcT/AUjZ1SfXmw8OFhp65Ykk97wEHC+pAeAPYDPjlPuPcC5ku4FVlOswgZwUfn++ygWWBrL5ymWfFxVvv+s8vwi4CZJt9n+BcWSkF+VtAq4m6IfYyNFf8gNZRPliUneAwDb91I0bx6kWHZyuGm3G/Cd8p53Ah8oz/9V2al8P/A94N5x7js9DLna1iVyj1edpjpJ+wPfsX1Yk6Ixjc3b+UV+/W6nNi8ILHn6CyuyCFdEjM3u+ac5SSZdNnoN5ohx9XgrIskkok84NZOIaF1mWouIOhjo4mPfKpJMIvqAAWdypIhomTM5UkTUpNdrJhm0FtEHJN1EMSt+FRtsn9TOeMaSZBIRtchncyKiFkkmEVGLJJOIqEWSSUTUIskkImqRZBIRtUgyiYhaJJlERC2STCKiFv8fbW6dzibKLuUAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 288x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.matshow(CMbroadbrush)\n",
    "plt.xticks(range(max(truth)+1), names)\n",
    "plt.yticks(range(max(truth)+1), names)\n",
    "plt.xlabel('predicted class')\n",
    "plt.ylabel('true class')\n",
    "plt.colorbar()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**3) Cruise control classifier CM:** \n",
    "\n",
    "This is where the confusion matrix has high values on the column of one specific class which means that the classifier constantly classifies all entries as one specific class"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "CM = np.zeros((M_classes, M_classes))\n",
    "class_corr = 2 #np.random.randint(0, M_classes, size=1) # randomly choose which class to classify all class on\n",
    "\n",
    "for i in range(M_classes):\n",
    "    for j in range(M_classes):\n",
    "        CM[i][j] = np.random.uniform(0., 1./M_classes, 1)\n",
    "        CM[i][class_corr] = np.random.uniform(1./M_classes, 1., 1)\n",
    "        \n",
    "CM = CM/np.sum(CM, axis=1)[:, np.newaxis] #normalizing to make sure each row sums to 1\n",
    "CMcruise = CM\n",
    "\n",
    "# saving probabilities\n",
    "I = proclam.classifiers.from_cm.FromCM(seed=None)\n",
    "predictionI = I.classify(CMcruise, truth, other=False)\n",
    "df          = pd.DataFrame(predictionI)\n",
    "\n",
    "save_info(truth, df, 'Cruise', seed=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x7f4c5f70ca58>"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARMAAAD0CAYAAAC4n8I2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAGjxJREFUeJzt3X2UHXWd5/H3JyEhkISghnEUUDIIYnQQJaAuww4gKoqCjjgr6O6wPmRdYcDl6AJHRYZ1dtVhZcZzGIf4cMYZH0AdRzISRU4G5EHUkAiBBJEYVMJ6hPAkD4Yk3Z/9o+rC7aa7b3Vu3b630p/XOXW6qm51/b7pdL751a9+9S3ZJiKiWzP6HUBE7BySTCKiFkkmEVGLJJOIqEWSSUTUIskkImqRZBIRtUgyiYhaJJlERC126XcAEdHZ646e6/sfGKp07Oq1T1xp+7geh/Q0SSYRDbD5gSF+fOU+lY6d9ZxfLOxxOGNKMoloBDPk4X4HMaEkk4gGMDDMYD+Um2QS0RDDpGcSEV0yZmjAy4VM21vDkh7tdwyDrp8/I0nnS/pgv9pvi2M/Sbf1Ow4oLnOqLP2SnklEAxgYGvAxk2nRM5H0Tkk/kXSzpEskzSz3XyRpnaSVkvaagjg+KukOSddL+pqkD0q6pozjJkm3SzpM0rck3Snp4z2O50OSzijXL5L07+X6MZK+0nbcQkk3Sjq+x/F8WNLPJV0PvLDc17efT5uZkj5X/q58X9Juks6QtF7SWkmX9joAA9s8XGnpl50+mUh6EfCfgCNsHwIMAe8A5gI32X4x8APgYz2O4zDgrcBLgdcDS9o+3mp7CfAPwOXAacBLgFMlPauHYV0HHFmuLwHmSZpV7ru2jPvZwBXAebav6FUgkg4F3g4cArwBOKzt4379fFoOAC4uf1ceovh7PAd4me2DgfdNQQwMV1z6ZTpc5rwaOBRYJQlgN+Beip/7ZeUxXwa+1eM4jgAut70F2CLp39o+W15+vRVYZ/s3AJI2AvsC9/coptXAoZL2AJ4A1lAklSOBM4BZwErgNNs/6FEMLUcC/2r7cQBJy9s+69fPp+Uu2zeX66uB/YC1wFckfRv4do/bLwZgc5nTdwK+ZPuQcnmh7fPHOK6ff1NPlF+H29Zb2z1L+La3AXcBpwI/pOipHA28ALgd2E7xj+d1vYqhor78fMZoH4qe7S7A8cDFwMsp/qPqbRyGoYpLv0yHZLISOEnSHwBIeqak51P82U8qjzkFuL7HcdwAvEnSHEnzgDf2uL2qrgM+SHFZcx1Fl/2nLl5bYOBdwEGSzu5xHNcCby7HI+YDb+pxe92YAexr+2rgbGABMK+XDRaT1nKZ01e210v6CPB9STOAbRTX3I8Bh5ef3UsxrtLLOFaVXfe1wG8puuwP97LNiq4DPgzcaPsxSVvKfQDYHpJ0MrBc0iO2/74XQdheI+ky4BaKv49VvWinJjOBL0taQNHz/Yzth3rbpBhCvW2iS8p7c6aOpHm2H5W0O8X/xEttr+l3XDH4XnLwbP/LFdWe3zvoeb9ZXQ5YT6mdvmcyYJZJWgzMoRjHSSKJSgxsHfBRiSSTKWT7lH7HEM017MG+zEkyiWiAYgZskklEdMmIoQG/zBns6HpI0tJ+x9BukOIZpFgg8bQMW5WWfpm2yQQYqF9QBiueQYoFEs+TlzlVln7JZU5EI4ghD/b//Y1LJrN32d1zdt2z6/PMmb2APeY+t+tJNsO7zuw6FoDZu+/JvGfu23U8u//h413HMu8Pd+cPFj+r61i23D2n61gA5uy6gD3m7911PNqytY5wmDNjHgtm7dV1PL/bvnmz7UpPqxvYRj2/a73SuGQyZ9c9eeXiwen1Prqop7OoJ+2Qs2/ufNAU+fkHXtTvEEbY5Y67+x3CCFfed8mvqh5r19szkXQc8HcUs3k/b/sTYxzz58D5FLnslk5TGxqXTCKmq+GaxkPKej4XA68BNlE8qLjc9vq2Yw4AzqUo3fFg69m2iSSZRDRAMQBbW8/kcGCD7Y0AZXGnE4H1bce8l6KGy4MAtu/tdNLBHtGJiFJxmVNlARaWlelay+hxgb2B9mu+TeW+dgcCB0q6QdKPysuiCaVnEtEARQmCyv/3b67hQb9dKCrMHQXsA1wr6Y8nejo6ySSiAYzY6tru5txDUaGuZZ9yX7tNwI9bBbQk/ZwiuYxbGiKXORENMewZlZYKVgEHSFokaTZF7d3lo475NkWvBEkLKS57Nk500vRMIhqgzgFY29slnQ5cSXFr+Iu210m6gKLI+vLys9dKWk9RqvJDtiestZtkEtEARgzV+NyN7RXAilH7zmtbN3BWuVSSZBLREJMYgO2LJJOIBrDJszkRUQfVNgO2V5JMIhrAwFYP9j/XwY4uIoBiADY1YCcg6VHbg/XYbcSAGvSyjemZRDSAoeqEtL6ZsugkvVPSTyTdLOmS8jFoJF0kaZ2klZIqFYqJmH6qlWzsZ9nGKUkmkl5E8frNI2wfQjGj7h3AXIoZdy8GfgB8bJzvX9p6AnLb9u4riUU0TatnUtN0+p6YqsucVwOHUhRhAdiN4n2yw8Bl5TFfBr411jfbXgYsA2optRjRRHlvTkEUr8M8d8RO6aOjjkuiiBiDLbYND/YQ51T1iVYCJ7VKv0l6pqTnl+2fVB5zCnD9FMUT0ShFPRNVWvplSlKd7fWSPgJ8X9IMYBtwGvAYcHj52b0U4yoR8TR51cWTbF/GU+MjLZljElFBMQCbMZOIqEEmrUVE1zKdPiJqk3omEdE1G7YNJ5lERJeKy5wkk4ioQWbARkTXcms4ImqSy5yIqElqwEZE14rq9IOdTAa73xQRQHE3Z/vwzEpLFZKOk3SHpA2Szhnj81Ml3VcWM7tZ0ns6nTM9k4iGqOsyp6xyeDHwGooXlK+StNz2+lGHXmb79KrnbV4ymSGGdpvV7yieNH/DI/0OYYRX7bGh3yE86c4bt/c7hBF+/9qX9zuEkb5b/dCa7+YcDmywvRFA0qXAicDoZDIpucyJaIgayzbuDdzdtr2p3DfaWyWtlfRNSft2OmmSSUQTuHjQr8oCLGzVTC6XpTvQ4r8B+9k+GLgK+FKnb2jeZU7ENNSqtFbRZttLJvj8HqC9p7FPue+p9uz72zY/D3yqU6PpmUQ0xCR6Jp2sAg6QtEjSbODtwPL2AyQ9p23zBOD2TidNzySiAQxsr+mpYdvbJZ0OXAnMBL5oe52kCyhePbMcOEPSCcB24AHg1E7nTTKJaIC6iyPZXgGsGLXvvLb1c4FzR3/fRJJMIhoi0+kjonvOU8MRUYOUIIiI2iSZRETXjBhKDdiIqEMGYCOia84AbETUxUkmEdG9vNEvImoyrXsmkj4KvBO4j6J+wmrgjcBPgSOBucB/oZi2+8cUlZ0+0suYIppoWs8zkXQY8FbgpcAsYA1FMgHYanuJpDOBy4FDKR4m+oWki0Y9/kxZj2EpwK67LuhVyBGDa5oXlD4CuNz2FtuPUBRbaWk97nwrsM72b2w/AWxkZJ0FAGwvs73E9pLZs+b2MOSIwWSKy5wqS7/0a8zkifLrcNt6azvjOBFPM/gDsL3smdwAvEnSHEnzKMZKImIH2dWWfulZL8D2KknLgbXAbykuaR7uVXsRO7tpfTcHuND2+ZJ2B64FVtv+XOtD29cA17RtH9XjeCIaqeh1TO9kskzSYmAO8CXba3rcXsROa9DHTHqaTGyf0svzR0wnw8PTOJlERD1Mf2/7VpFkEtEQfbxRU0mSSUQTZAA2Imoz4F2TJJOIhtipeiaSZgDzbP+uR/FExDj6Obu1io7T6SV9VdIekuYCtwHrJX2o96FFRIsNHp5RaalC0nGS7pC0QdI5Exz3VkmWNNGL0IFqz+YsLnsibwa+CywC/nOliCOiNnU9myNpJnAx8HpgMXByObl09HHzgTOBH1eJr0oymSVpFkUyWW57GwM/FBSxE3LFpbPDgQ22N9reClwKnDjGcf8L+CSwpcpJqySTS4BfUlRFu1bS84GMmURMqWq1TMpB2oWSbmpblo462d4UlQ9bNpX7nmpNejmwr+0rqkbYcQDW9meAz7Tt+pWko6s20BMzB2dUe+uzdut3CCNcft8h/Q7hSTP3H6yqeLv8fnu/Q+hO9euBzbY7jnGMp7zR8mng1Ml8X5UB2DPLAVhJ+oKkNcAxOxZmROwQ11pp7R5GVjTcp9zXMh94CXCNpF8CrwSWdxqErXKZ865yAPa1wDMoBl8/USXiiKhRfWMmq4ADJC2SNBt4O0+VUsX2w7YX2t7P9n7Aj4ATbN800UmrJJNWqnsD8M+217Xti4ipYlVbOp3G3g6cDlwJ3A583fY6SRdIOmFHw6syaW21pO9T3BI+t7xdNLyjDUbEDqrxHqrtFcCKUfvOG+fYo6qcs0oyeTdwCLDR9uOSngX81yonj4iamEq9jn6qcjdnWNJdwIGS5kxBTBExhkGfTt8xmUh6D8UsuH2AmylGdm8kd3QiptaAJ5MqA7BnAocBv7J9NPAy4KGeRhURT1fTAGyvVBkz2WJ7iyQk7Wr7Z5Je2PPIIuIpBg34bY8qyWSTpD2BbwNXSXoQ+FVvw4qIkfrb66iiygDsW8rV8yVdDSwAvtfTqCLi6QZ8zGTcZCLpmWPsvrX8Og94oCcRRcTYmppMgNUU4bf3rVrbBv6oh3FFxGhNTSa2F01FAJIetT1vKtqKaKwGTFqr8tTwWyQtaNveU9KbextWRIwmV1v6pco8k4/Zfri1Yfsh4GNVG5D0IUlnlOsXSfr3cv0YSV9pO26hpBslHV89/IhppL6nhnuiSjIZ65jJVLW/DjiyXF8CzCvLQB4JXAsg6dnAFcB5k6nsFDGd7Aw9k5skfVrS/uXyaYrB2apWA4dK2gN4gmIq/hKKZHIdMAtYCfxP21eNdQJJS1sl6LZue2wSTUfsRAZ8BmyVZPKXwFbgMorCs1uA06o2UBagvouiBNwPKRLI0cALKGopbKdIOK+b4BzLbC+xvWT2rLlVm47YeVS9xOljz6TKpLXHgHHfq1HRdcAHgXdRzFX5NLDatiW53P8NSWfb/mSXbUXsnAb81nC1N/Z07zrgOcCNtn9L0bu5rvWh7SHgZOAYSe+fopgiGmXQx0ym5F3DtldSjI20tg9sW59Xfn2CCS51Iqa9Ae+Z5MXlEQ2gBjw1XGXS2oGSVkq6rdw+WNJHeh9aRIywE9zN+RxwLrANwPZaitL4ETGVmn43B9jd9k+kERmv4a9Gi2iefg6uVlElmWyWtD9lzpN0EvCbnkYVEU+3EyST04BlwEGS7qGYgPbOnkYVESP1+bZvFR3HTGxvtH0ssBdwkO0/sf3LnkcWESPVOGYi6ThJd0jaIOlpk1IlvU/SrZJulnS9pMWdzlnlVRfnjdoGwPYF1cKOiDrUdWtY0kzgYuA1wCZglaTltte3HfZV2/9QHn8Cxaz14yY6b5W7OY+1LUPA64H9JvsHiIiBcTiwobzq2ErxzN2J7QfY/l3b5lwq9HmqPJvzf9u3JV1I8cLjiJhK9Y2Z7A3c3ba9CXjF6IMknQacBcymwkv3duTZnN0p3u4XEVOl4nM55SDtwlbJjnJZukNN2hfb3h84G+g4UbXKmMmtPJUTZ1IMxGa8JGKqVe+ZbLa9ZILP7wH2bdvep9w3nkuBz3ZqtMqt4Te2rW8Hfmu7v5PWhgbnHtkuP5hMnajee88l6/odwpP+ljf1O4QRts6f1fmgQVbfr/0q4ABJiyiSyNuBU9oPkHSA7TvLzeOBO+lgwmRSjvpeafugHQo5Imoh6ptnYnu7pNMpxj5nAl+0vU7SBcBNtpcDp0s6luIxmgeBv+h03gmTie2h8l7082z/uvs/RkTskJqfGra9Algxat95betnTvacVS5zngGsk/QTitvDrcZOmGxjEdGFwbm6H1OVZPLRnkcREZ3tBMnkDbbPbt8h6ZPAD3oTUkSMpfHP5lBMuR3t9XUHEhEdNLWeiaT/Drwf+CNJa9s+mg/c0OvAIqJNnxNFFRNd5nwV+C7wfxj5qotHbD/Q06gi4mkGvQbsuMmkfL/wwxSvoIiIPhv0MZNUp49oiiSTiOhaw8dMImJAqFwGWZJJRFOkZxIRdcgAbEnS+cCjti+cqjYjdipNvTUcEQNkZ3jVRTckfVjSzyVdD7yw3HeNpIvKcnK3SzpM0rck3Snp472MJ6LRmjqdvluSDqWo4HRI2c4aoFWWbKvtJZLOBC4HDgUeAH4h6SLb948611JgKcCuuy7oVcgRA20690yOBP7V9uNl2fzlbZ+11m8F1tn+je0ngI2MrE0JgO1ltpfYXjJ71twehhwxwKZrz6SDJ8qvw23rre2M40SMYTr3TK4F3ixpN0nzYcCqC0c0SdVeyc7YM7G9RtJlwC3AvRQVsSNiB4gGPzVcB9t/Dfz1qN0Xtn1+DXBN2/ZRvYwnotEG/DIn4xMRDSEPdjZJMologjw1HBF1GfS7OUkmEU0x4Mmkp9PpI6I+crWl0rmk48q3dW6QdM4Yn58lab2ktZJWSnp+p3MmmUQ0Qfl60CpLJ+U7xC+meGXNYuBkSYtHHfZTYIntg4FvAp/qdN4kk4imqG/S2uHABtsbbW8FLgVOHNGUfbXtx8vNHwH7dDppkklEA4hJXeYsLJ/Kby1LR51ub+Dutu1N5b7xvJvitTcTygBsRFNUn2ey2faSOpqU9E5gCfCnnY5NMoloiBpvDd/DyKfz9yn3jWxPOhb4MPCn5VP9E8plTkQT1Pug3yrgAEmLJM2mqDvUXiIESS8DLgFOsH1vlZOmZxLREHU96Gd7u6TTgSuBmcAXba+TdAFwk+3lwN8A84BvSAL4te0TJjpv85LJY79n5o9u63cUT5qx1179DmGE1+6+rd8hPOlvZw3Wr9eM7QM+66uDOp8atr0CWDFq33lt68dO9pyD9bcdEWMzkxmA7Yskk4iGyLM5EVGPJJOI6FZr0togSzKJaAI7YyYRUY9pXQM2IuqTy5yI6J6B4cHOJkkmEU0x2LkkySSiKXKZExH1yN2ciKhDeiYR0TUZlAHYiKjFgM8zmbLiSJL2kzQ4tQMiGkZ2paVf0jOJaIIGvB50qss2zpT0OUnrJH1f0m6Szmh72c+lUxxPREP4qedzOi19MtU9kwOAk22/V9LXgbcC5wCLbD8hac+xvqks1b8UYA67T1mwEYNk0O/mTHXP5C7bN5frq4H9gLXAV8qS+tvH+ibby2wvsb1klnadmkgjBs2A90ymOpm0l8sfougZHU/xqsKXA6skZRwnYjSDhlxp6Zd+v+piBrCv7auBs4EFFBWxI2K0+l510RP97gXMBL4saQFFManP2H6ozzFFDKR+3vatYsqSie1fAi9p275wqtqO2CkkmURE10xmwEZE90S12a9VL4UkHSfpDkkbJJ0zxuf/UdIaSdslnVTlnEkmEU1R061hSTMp7qC+HlgMnCxp8ajDfg2cCny1ani5zIloAgP13fY9HNhgeyNAOfP8RGD9k80VY5xI1ctYJ5lENMQk7uYslHRT2/Yy28vatvcG7m7b3gS8osvwkkwiGqN6Mtlse0kvQxlLkklEI9Q6Vf4eYN+27X3KfV3JAGxEE5g6n81ZBRwgaZGk2cDbgeXdhphkEtEUwxWXDmxvB04HrgRuB75ue52kCySdACDpMEmbgLcBl0ha1+m8ucyJaIg6p9PbXgGsGLXvvLb1VRSXP5UlmUQ0gYGhwZ4Cm2QS0Qj9rVVSReOSySN+cPNV2y79VQ2nWghs7vos93YfSKmWeGY+p4ZI6vrZ8L+7P0Whnng6XvVXVtPPh+dP6ugkk3rZ3quO80i6qR/34sczSPEMUiyQeJ6UZBIRXTOQl3BFRPcMzgDsoFrW+ZApNUjxDFIskHgacTdn2k5aG/XgU9/VGY+kR8uvz5X0zQ7HfkDSiPeHdIpF0lGSvjOJeK6RtMNjDDvz39UkG051+uheWYNiUmz/P9udCtt8APIyokZIMomJlO9g/pmkr0i6XdI3Wz0FSb+U9ElJa4C3Sdpf0vckrZZ0naSDyuMWSbpR0q2SPj7q3LeV6zMlXSjptvLtiX8p6QzgucDVkq4uj3ttea41kr4haV65/7gyzjXAn43zZ3laG2Mc81lJN5Vvdfyrtv2faHuz44XlvreV57pF0rX1/MSbKm/0i2peCLzb9g2Svgi8H2gV3L7f9ssBJK0E3mf7TkmvAP4eOAb4O+Cztv9J0mnjtLGU4qVnh9jeLumZth+QdBZwtO3NkhYCHwGOtf2YpLOBsyR9Cvhc2dYG4LKqbYxxzIfLdmcCKyUdTPHE6luAg2y77c2O5wGvs33PeG97nDYMDGfMJDq72/YN5fqXgT9p++wygLKH8B+Ab0i6GbgEaE1ROwL4Wrn+z+O0cSxwSfmQF7YfGOOYV1KU8buhbOMvKCZWHUTxNsY7bbuMcUfb+POyd/NT4MVlew8DW4AvSPoz4PHy2BuAf5T0XorXokxv6ZlEBaN/A9q3Hyu/zgAesn1IxXPsCAFX2T55xE5pvDYnd3JpEfBB4DDbD0r6R2BO2Ys5HHg1cBLFE63H2H5f2QM7Hlgt6VDb99cRSyMN+KS19EwGw/MkvapcPwW4fvQBtn8H3CXpbQAqvLT8+AaKmhQA7xinjauA/9Z6/WrbJcgjwPxy/UfAEZJeUB4zV9KBwM+A/STtXx43ItlUaKNlD4rk+LCkZ1MUNG71uhaUT7L+D+Cl5f79bf+4fJr1PkYW9JlebDw0VGnplySTwXAHcJqk24FnAJ8d57h3AO+WdAvFkyYnlvvPLL//Vor6nmP5PEXF8bXl959S7l8GfE/S1bbvo6hI/jVJa4EbKcYxtlCMh1xRXqKM90TSeG0AYPsWisubn1FUPW9d2s0HvlO2eT1wVrn/b8pB5duAHwK3jNPu9DDsakufyAPeddrZSdoP+I7tl3Q4NKaxBbvs5VfNP7HzgcCVD31hdWrARsTY7IG/m5Nk0mej38EcMa4Bv4pIMoloCKdnEhHdS6W1iKiDgT7e9q0iySSiAQw4xZEiomtOcaSIqMmg90wyaS2iASR9j6IqfhWbbR/Xy3jGkmQSEbXIszkRUYskk4ioRZJJRNQiySQiapFkEhG1SDKJiFokmURELZJMIqIWSSYRUYv/D76H5qQMUnT9AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 288x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.matshow(CMcruise)\n",
    "plt.xticks(range(max(truth)+1), names)\n",
    "plt.yticks(range(max(truth)+1), names)\n",
    "plt.xlabel('predicted class')\n",
    "plt.ylabel('true class')\n",
    "plt.colorbar()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**4) Subsumed classifier CM:** \n",
    "\n",
    "This is where the classifier consistently classifies class M and class M'. So, the confusion matrix will have a higher value on the element on row M and column M' compared to other values in that row so that it is misclassified most of the times. Other classes may be classified correctly randomly as I assign them random uniform values between 0 and 1. There could be more restriction on the other classes if we want."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "CM = np.zeros((M_classes, M_classes))\n",
    "class_M      = 1 #np.random.randint(0., M_classes, size=1)[0]\n",
    "class_Mprime = 2 #np.random.randint(0., M_classes, size=1)[0]\n",
    "\n",
    "#This is to make sure we do not misclassify class M with the same class (which makes no sense)\n",
    "n = 0\n",
    "while class_M == class_Mprime:\n",
    "    n+=1\n",
    "    print((\"going in to the loop\", n, \"times!\"))\n",
    "    class_Mprime = np.random.randint(0., M_classes, size=1)\n",
    "\n",
    "CM = cm.copy()#np.random.uniform(0., 1, 1)\n",
    "CM[class_M] = CM[class_Mprime]\n",
    "# for i in range(M_classes):\n",
    "#     for j in range(M_classes):\n",
    "#         CM[i][j] = np.random.uniform(0, 1., 1)\n",
    "#         CM[class_M][class_Mprime]  = np.random.uniform(1./M_classes, 1., 1)\n",
    "\n",
    "CM = CM/np.sum(CM, axis=1)[:, np.newaxis] #normalizing to make sure each row sums to 1\n",
    "CMsubsumedto2 = CM\n",
    "\n",
    "# saving probabilities\n",
    "J = proclam.classifiers.from_cm.FromCM(seed=None)\n",
    "predictionJ = J.classify(CMsubsumedto2, truth, other=False)\n",
    "df          = pd.DataFrame(predictionJ)\n",
    "\n",
    "save_info(truth, df, 'SubsumedTo', seed=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2 3\n"
     ]
    }
   ],
   "source": [
    "CM = np.zeros((M_classes, M_classes))\n",
    "class_M      = 2 #np.random.randint(0., M_classes, size=1)[0]\n",
    "class_Mprime = 3 #np.random.randint(0., M_classes, size=1)[0]\n",
    "\n",
    "#This is to make sure we do not misclassify class M with the same class (which makes no sense)\n",
    "n = 0\n",
    "while class_M == class_Mprime:\n",
    "    n+=1\n",
    "    print((\"going in to the loop\", n, \"times!\"))\n",
    "    class_Mprime = np.random.randint(0, M_classes, size=1)[0]\n",
    "print(class_M, class_Mprime)\n",
    "\n",
    "# for i in range(M_classes):\n",
    "#     for j in range(M_classes):\n",
    "CM = cm.copy()#np.random.uniform(0., 1, 1)\n",
    "CM[class_M] = CM[class_Mprime]#np.random.uniform(1./M_classes, 1., 1)\n",
    "\n",
    "# print(CM,np.sum(CM, axis=1)[:, np.newaxis])\n",
    "CM = CM/np.sum(CM, axis=1)[:, np.newaxis] #normalizing to make sure each row sums to 1\n",
    "CMsubsumedfrom2 = CM\n",
    "\n",
    "\n",
    "# saving probabilities\n",
    "J1           = proclam.classifiers.from_cm.FromCM(seed=None)\n",
    "predictionJ1 = J1.classify(CMsubsumedfrom2, truth, other=False)\n",
    "df           = pd.DataFrame(predictionJ1)\n",
    "\n",
    "save_info(truth, df, 'SubsumedFrom', seed=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x7f4c5f56da58>"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARMAAAD0CAYAAAC4n8I2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAGVZJREFUeJzt3XuUXlV5x/HvL5PEEBKiGHShRkKRSxEhmuCllBYQBUUpFmwBUamX1IpCS7HoApG2drVaKqushZZYXbUFJEqtRFEuK4IkGJUkQiBBRQm0WGoM98vKZWae/nHOG96ZzMy7J+857zln5vdZ66y8lzNnPwzhYe999tmPIgIzs25NqToAM5sYnEzMrBBOJmZWCCcTMyuEk4mZFcLJxMwK4WRiZoVwMjGzQjiZmFkhplYdgJl1dtzRu8cjjw4knbtm3dYbI+L4kkPaiZOJWQNsfnSAH934sqRzp+39y7klhzMiJxOzRggGYrDqIMbkZGLWAAEMUu+Hcp1MzBpiEPdMzKxLQTBQ8+1CJu2tYUlPVx1D3VX5O5J0saTzqmq/LY75ku6pOg7IhjkpR1XcMzFrgAAGaj5nMil6JpLOkPRjSXdKukJSX/75pZLWS1ouaa8exPFJST+TtFLSVyWdJ+nWPI7Vku6VdLikb0i6T9KnS47nY5LOzl9fKul7+etjJF3Vdt5cSasknVByPBdI+rmklcCB+WeV/X7a9En6Yv535SZJu0k6W9IGSeskXVN2AAFsj8GkoyoTPplI+m3gj4EjImIBMAC8C9gdWB0RrwS+D3yq5DgOB04GDgPeAixq+3pbRCwC/gW4DjgLOAQ4U9ILSwxrBXBk/noRMEvStPyz2/K4XwxcD1wUEdeXFYikhcCpwALgrcDhbV9X9ftp2R+4PP+78jjZv8ePA6+OiEOBD/UgBgYTj6pMhmHOG4GFwB2SAHYDNpH93pfm51wJfKPkOI4ArouILcAWSd9q+25Z/ufdwPqIeBhA0v3APOCRkmJaAyyUtAewFVhLllSOBM4GpgHLgbMi4vslxdByJPBfEfEsgKRlbd9V9ftp2RgRd+av1wDzgXXAVZK+CXyz5PazCVgPcyon4CsRsSA/DoyIi0c4r8p/U1vzPwfbXrfel5bwI2I7sBE4E/gBWU/laOAVwL1AP9l/PMeVFUOiSn4/I7QPWc92KnACcDnwGrL/UZUbR8BA4lGVyZBMlgOnSHoRgKQ9Je1D9s9+Sn7O6cDKkuO4HXi7pBmSZgFvK7m9VCuA88iGNSvIuuw/iaxsQQDvAw6SdH7JcdwGnJTPR8wG3l5ye92YAsyLiFuA84E5wKwyG8wWrXmYU6mI2CDpQuAmSVOA7WRj7meA1+bfbSKbVykzjjvyrvs64NdkXfYnymwz0QrgAmBVRDwjaUv+GQARMSDpNGCZpKci4vNlBBERayUtBe4i+/dxRxntFKQPuFLSHLKe72UR8Xi5TYoBVG4TXZLr5vSOpFkR8bSkmWT/J14cEWurjsvq75BDp8d/Xp/2/N5BL394TT5h3VMTvmdSM0skHQzMIJvHcSKxJAFsq/mshJNJD0XE6VXHYM01GPUe5jiZmDVAtgLWycTMuhSIgZoPc+odXYkkLa46hnZ1iqdOsYDjaRkMJR1VmbTJBKjVX1DqFU+dYgHHs2OYk3KkkHR8/pzYLyR9fITvL82fZbszf16q461vD3PMGkEMRDH/788fdL0ceBPwENkK3mURsaF1TkT8Rdv5HwVe3em6jUsmL9xzSsyb133YL3tpHwsOm971IpuN9xWzd++MaXOYM/Ml3S/66e/vPpa+WcyZ/qKuY4nt3ccCMIOZ7KE9u45H06cVEQ4z+mYz53kv7jqeJ7dt2hwRSU+rB7Cdvm6bbHkt8IuIuB8gf+r5D4ANo5x/GgkPwjYumcybN5WbvlPJ5tsjevfbPlB1CENM2fRY1SHs0P/w/1UdwhBTXzKv6hCGuOGBSx9MPTeiuJ4J8FLgf9rePwS8bqQT80dP9gW+1+mijUsmZpPVYPqt4bmSVre9XxIRS3ax2VOBayOiY9EeJxOzBsgmYJN7Jps7LKf/FdnWDS0vyz8byalkz7J15GRi1giFDnPuAPaXtC9ZEjmV7Mn5oS1KBwEvAFalXNTJxKwBsi0IikkmEdEv6SPAjWRPQH85ItZL+huy3Qdbm1GdClwTiU8DO5mYNUAgtkVhd3OIiO8A3xn22UXD3l88nms6mZg1xGBxw5xSOJmYNcA4J2Ar4WRi1gCBGPAWBGZWhKImYMviZGLWABEUeWu4FE4mZo2g8ayArYSTiVkDBLAt6v2fa72jMzMgm4D1HrBjkPR0RJRavMhsovCtYTPrWlD/RWs9i07SGZJ+nG8Dd0W+21Nre7j1kpZLStooxmzySduyscod7HuSTCT9Nln5zSMiYgFZ8ed3AbuTPVj0SuD7jLKbk6TFklZLWv3II1VWUzWrRqtnknJUpVfDnDcCC8n2mgTYjaye7CCwND/nSuAbI/1wvrHLEqCQrRbNmsh1czIiK4f5iSEfSp8cdp4ThdkIIsT2wXpPcfaqT7QcOEXSiwAk7ZnvLTkFOCU/53RgZY/iMWuUbD8TJR1V6Umqi4gNki4EbpI0BdhOthXcM8Br8+82kc2rmNlOCt1prRQ96zdFxFKemx9p8RoTswTZBKznTMysAF60ZmZd83J6MytM3fczqXd0ZgZk+5lsH5ySdKToVLg8P+ePJG3IV6hf3ema7pmYNUA2zOld4XJJ+wOfIFu1/lhrWcdY3DMxa4gCn83ZUbg8IrYBrcLl7T4IXB4RjwFExKZOF3UyMWuA1q3hlIO81nDbsXjY5UYqXP7SYeccABwg6XZJP5R0fKcYPcwxa4RxDXM61RpOMRXYHziKrBbxbZJeFRGPj/YD7pmYNUSBy+lTCpc/BCyLiO0RsRH4OVlyGZWTiVkDZLvTK+lIsKNwuaTpZDWFlw0755tkvRIkzSUb9tw/1kU9zDFrgED0DxZTazixcPmNwJslbSDbf+hjEfHIWNd1MjFriCKfCO5UuDwiAjg3P5I0Lpls/PkLec+x7606jB2mbH2m6hCGeMOND1Qdwg6rTjyg6hCGGNhrTtUhDPVA+ql+0M/MClP3DaWdTMyaIPygn5kVoLXTWp05mZg1hHsmZta1APoTnwiuipOJWQN4cyQzK4znTMyse+E5EzMrgBetmVlhnEzMrGuBGPDdHDMrgidgzaxr4QlYMytKOJmYWfe8aM3MCjKpeyaSPgmcAfyGbGv9NcDbgJ8ARwK7A+8hK/bzKmBpRFxYZkxmTTSp15lIOhw4GTgMmAasJUsmANsiYpGkc4DrgIXAo8AvJV06fK/JvO7HYoAZU/coK2Sz+so3lK6zMm9cHwFcFxFbIuIp4Ftt37V2wr4bWB8RD0fEVrLdr+cNuw4RsSQiFkXEoulTZ5YYslk9BdkwJ+VI0anWsKQzJf1G0p358YFO16xqzmRr/udg2+vWe8/jmO2kuAnYlFrDuaUR8ZHU65bZM7kdeLukGZJmkc2VmNkuikg7EqTUGh630pJJRNxBNpxZB3yXbEjzRFntmU10BQ5zUmoNA5wsaZ2kayXtNP0wXNmL/S+JiAOA44B9gDURcVRErAaIiFsjYkePpf07M3tO1utITiadCpen+BYwPyIOBW4GvtLpB8qen1gi6WBgBvCViFhbcntmE9Y45kw6FS7vWGt42B3VfwU+26nRUpNJRJxe5vXNJpPBwcJuDe+oNUyWRE4Fhvy3KmnviHg4f3sicG+ni/rOiVkDBOm3fTteK63W8NmSTgT6ydaAndnpuk4mZg2RdqMm8Vqdaw1/gmxlejInE7MmiEn+bI6ZFajIrkkJnEzMGmJC9UwkTQFmRcSTJcVjZqNIXN1amY6L1iRdLWkPSbsD9wAbJH2s/NDMrCUCYnBK0lGVlJYPznsiJ5Eti98XeHepUZnZTgp8NqcUKclkmqRpZMlkWURsp/ZTQWYTUCQeFUlJJlcAD5DtinabpH0Az5mY9VTaczlVTtJ2nICNiMuAy9o+elDS0eWFZGYjqvl4IGUC9px8AlaSviRpLXBMD2Izs5bxPTVciZRhzvvyCdg3Ay8gm3z9h1KjMrOd1XzOJGWdSSvVvRX4j/yBoHqvnjGbiGq+aC2lZ7JG0k1kyeRGSbPJ9mo1s16aAD2T9wMLgPsj4llJLwT+pNywzGyIoPY9k5S7OYOSNgIHSJrRg5jMbAR1X07fMZnk9TLOIdva7U7g9cAqfEfHrLdqnkxS5kzOAQ4HHoyIo4FXA4+XGpWZ7SyUdlQkZc5kS0RskYSk50XETyUdWHpkZvacANX8tkdKMnlI0vOBbwI3S3oMeLDcsMxsqGp7HSlSJmDfkb+8WNItwBzghlKjMrOdNXXORNKeww+yqnwrgVk9i9DMMgWuM+lUuLztvJMlhaSx6vAAY/dM1uShtfetWu8D+K20sM2sEAX1TFILl+cLVM8BfpRy3VGTSUTsu+vhppP0dES4p2M2lmIXre0oXA4gqVW4fMOw8/4W+AyQtLNiylPD75A0p+398yWdlBq1mRVDkXbQudZwx8Llkl4DzIuI61PjS1ln8qmIeKL1JiIeBz6V2oCkj0k6O399qaTv5a+PkXRV23lzJa2SdELqtc0mlfQ5k80RsajtWDKeZvKN4z8H/OV4fi4lmYx0znh2tV8BHJm/XgTMyreBPBK4DUDSi4HrgYvGkwnNJpNx9Ew66VS4fDZwCHCrpAfIVr0v6zQJm5JMVkv6nKT98uNzZJOzqdYACyXtAWwlW4q/iCyZrACmAcuBv4qIm0e6gKTFrS7btv5nx9G02QRS3ArYHYXLJU0nK1y+bEczEU9ExNyImB8R84EfAidGxOqxLpqSTD4KbAOWAtcAW4CzUiLOA9sObCQrfPwDsgRyNPAKssrq/WQJ57gxrrGk1WWbPnVmatNmE0fqECehZxIR/UCrcPm9wNdahcvzYuW7JGXR2jPAqPehE60AzgPeR7ZW5XPAmogISZF//nVJ50fEZ7psy2xiKnDRWqfC5cM+Pyrlmr2q2LMC2BtYFRG/JuvdrGh9GREDwGnAMZI+3KOYzBqlwDmTUvSk1nBELCebG2m9P6Dt9az8z62MMdQxm/RqvpzehcvNGkANeGo4ZdHaAZKWS7onf3+opAvLD83Mhqj5fiYpcyZfBD4BbAeIiHVkt5LMrJcmwIbSMyPix8OqW/SXFI+ZjaLKydUUKclks6T9yHOepFOAh0uNysx2NgGSyVnAEuAgSb8iW4B2RqlRmdlQFd/2TZGyaO1+4FhJuwNTIuKp8sMys500PZlIumjYewAi4m9KisnMRlD3W8Mpw5xn2l7PAN5Gtp7fzGyHlGHOP7W/l3QJ2QNCZtZLTR/mjGAm2f4HZtYrE2ECVtLdPJcT+4C9AM+XmPVa05MJ2RxJSz/w63w/hErEFDE4uz710/u212v93opD6/O7+fyDV1cdwhBnv/6UqkPoTpOTSb4l/o0RcVCP4jGzEYj6D3PGfDYn32fkZ5Je3qN4zGwk+VPDKUdVUoY5LwDWS/oxbbeJI2KXt3czs11Q855JSjL5ZOlRmFlnEyCZvDUizm//QNJngO+XE5KZjaTRcya5N43w2VuKDsTMOuhh4XJJH5J0t6Q7Ja2UdHCna46aTCT9Wb7G5EBJ69qOjcC6tJDNrBAFlrpoK1z+FuBg4LQRksXVEfGqiFgAfJasosSYxhrmXA18F/h7hpa6eCoiHu0cspkVqcA7NR0Ll0fEk23n705Cmho1meT1hZ8gK0FhZhUbx5zJXEnt1feWDKs3PFLh8tft1J50FnAuMB04plOj3p3erCnSk8nmiBizLnBScxGXA5dLOh24EHjvWOf3qgiXmXWjwDkTOhcuH+4a4KROF3UyMWsAjeNIMGbhcgBJ+7e9PQG4r9NFPcwxa4qC1plERL+kVuHyPuDLrcLlwOqIWAZ8RNKxZCVuHqPDEAecTMwao8hFa50Kl0fEOeO9Zs+SiaSLgacj4pJetWk2oUyAPWDNrGoN2Gmt1AlYSRdI+rmklcCB+We3SrpU0mpJ90o6XNI3JN0n6dNlxmPWaBOgPOgukbSQbJZ4Qd7OWmBN/vW2iFgk6RzgOmAh8CjwS0mXRsQjw661GFgMMGP6nLJCNqu1ydwzORL4r4h4Nl+a237rqfX6bmB9RDwcEVuB+xl6/xuAiFgSEYsiYtG0qTNLDNmsxiZrz6SDrfmfg22vW+89j2M2gsncM7kNOEnSbpJmA28vsS2zia3YFbClKK0XEBFrJS0F7gI2ka26M7NdICZGedBdFhF/B/zdsI8vafv+VuDWtvdHlRmPWaPVfJjj+QmzhlDUO5s4mZg1QcXzISmcTMwaou53c5xMzJrCycTMiuCeiZl1Lyb5rWEzK5B7JmbWLeFhjpkVxetMzKwIde+ZeHd6syYo+EG/hFrD50rakJcEXi5pn07XdDIxawgNph0dr5NWa/gnwKKIOBS4lqze8JgaN8zR9gH6Hq5PqeOYOaPqEIbo22uvqkPY4SMHHlt1CEN89/4bqg5hiL69x3d+j2sN39J2/g+BMzpd1D0TsyYIsgnYlKOzkWoNv3SM898PfLfTRRvXMzGbrAosXJ7epnQGsAj4/U7nOpmYNUVxhcuTag3nFf0uAH4/36N5TB7mmDVAa9FaypEgpdbwq4ErgBMjYlPKRd0zMWuC9PmQhEsl1Rr+R2AW8HVJAP8dESeOdV0nE7OGKPJBv4Raw+O+FedkYtYQdV8B62Ri1gQBDNY7mziZmDVFvXOJk4lZU3iYY2bF8BYEZlYE90zMrGsKkCdgzawQNd9QumfL6SXNl3RPr9ozm2gUkXRUxT0TsyZoQHnQXj/o1yfpi5LWS7pJ0m6Szm7bHu6aHsdj1hCJe5lMop7J/sBpEfFBSV8DTgY+DuwbEVslPX+kH5K0GFgMMKNvds+CNauTut/N6XXPZGNE3Jm/XgPMB9YBV+WbsPSP9EMRsSQiFkXEoulTdutNpGZ1U/OeSa+TSfsGKwNkPaMTyDa3fQ1whyTP45gNF6CBSDqqUvXmSFOAefnmtecDc8j2UDCz4QosdVGGqnsBfcCVkuaQbSZ1WUQ8XnFMZrVU5W3fFD1LJhHxAHBI2/tLetW22YTgZGJmXQtqvwLWycSsAUS1q1tTOJmYNUXNk0nVd3PMLEUAA5F2JEgoXP57ktZK6pd0Sso1nUzMGqKoB/0SC5f/N3AmcHVqfB7mmDVFccOclMLlD+TfJU/7umdi1gjjetBvrqTVbcfiYRcbb+HyJO6ZmDVBMJ6eSadaw6VwMjFriuLWmSQVLh8vD3PMGqLAndY6Fi7fFU4mZk0QwMBg2tHpUhH9QKtw+b3A11qFyyWdCCDpcEkPAe8ErpC0vtN1Pcwxa4Ri9ypJKFx+B9nwJ1njksmT2zdtvuGhyx4s4FJzgc0FXKcodYqnTrFAQfH07V1AJJmifj/7jOvsmq+AbVwyiYi9iriOpNVVzHiPpk7x1CkWcDw7OJmYWdcCcBEuM+teQNR7D4LJnEyWVB3AMHWKp06xgON57m5OjU3aW8MRUau/oEXGI+np/M+XSLq2w7l/LmnmeGKRdJSkb48jnlsl7fIcw0T+dzXOhr07vXUvf9JzXCLifyOi0+Pjfw7M7HCO1YGTiY0lr8H8U0lXSbpX0rWtnoKkByR9RtJa4J2S9pN0g6Q1klZIOig/b19JqyTdLenTw659T/66T9Ilku7Jqyd+VNLZwEuAWyTdkp/35vxaayV9XdKs/PPj8zjXAn84yj/LTm2McM4X8ofP1kv667bP/6GtsuMl+WfvzK91l6TbivmNN5Ur+lmaA4H3R8Ttkr4MfBhobbj9SES8BkDScuBDEXGfpNcBnweOAf4Z+EJE/Luks0ZpYzFZ0bMFEdEvac+IeFTSucDREbFZ0lzgQuDYiHhG0vnAuZI+C3wxb+sXwNLUNkY454K83T5guaRDyZ4LeQdwUEREW2XHi4DjIuJXo1V7nDQCGPSciXX2PxFxe/76SuB3275bCpD3EH4H+LqkO4ErgNYyrCOAr+av/2OUNo4FrsiXUhMRj45wzuvJNsu5PW/jvWQLqw4iq8Z4X0REHuOutvFHee/mJ8Ar8/aeALYAX5L0h8Cz+bm3A/8m6YNkZVEmN/dMLMHwvwHt75/J/5wCPB4RCxKvsSsE3BwRpw35UBqtzfFdXNoXOA84PCIek/RvwIy8F/Na4I3AKWTPjRwTER/Ke2AnAGskLYyIR4qIpZFqvmjNPZN6eLmkN+SvTwdWDj8hIp4ENkp6J4Ayh+Vf30725CfAu0Zp42bgT1vlV9uGIE8BrWrwPwSOkPSK/JzdJR0A/BSYL2m//LwhySahjZY9yJLjE5JeTLZtYKvXNSd/XuQvgMPyz/eLiB/lz4z8hqGPzU8uEcTAQNJRFSeTevgZcJake4EXAF8Y5bx3Ae+XdBewnmyrPYBz8p+/m9F3zPpXsn091+U/f3r++RLgBkm3RMRvyPb9/KqkdcAqsnmMLWTzIdfnQ5RN42wDgIi4i2x481OyvUVbQ7vZwLfzNlcC5+af/2M+qXwP8APgrlHanRwGI+2oiKLmXaeJTtJ84NsRcUiHU20SmzN1r3jD7D/ofCJw4+NfWuOd1sxsZBG1v5vjZFKx4TWYzUZV81GEk4lZQ4R7JmbWvWrXkKRwMjFrggAqvO2bwsnErAECCG+OZGZdC2+OZGYFqXvPxIvWzBpA0g1ku+Kn2BwRx5cZz0icTMysEH42x8wK4WRiZoVwMjGzQjiZmFkhnEzMrBBOJmZWCCcTMyuEk4mZFcLJxMwK8f81O2FPmxFvcwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 288x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.matshow(CMsubsumedto2)\n",
    "plt.xticks(range(max(truth)+1), names)\n",
    "plt.yticks(range(max(truth)+1), names)\n",
    "plt.xlabel('predicted class')\n",
    "plt.ylabel('true class')\n",
    "plt.colorbar()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x7f4c5f7aca90>"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARMAAAD0CAYAAAC4n8I2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAGiZJREFUeJzt3X2UXXV97/H3J0NiyANRGhRRBIogRUQwidZyaQFRUQSxYC8gt1KpuVxRqBQucHmQ22vX9YHKravUEqpLKyBBaiWK8rAikRBR80AIJKggiIVSAgSQh5Wnme/9Y+8znJnMzPmdnL3PPnvm82LtNWef2fPbXyYz3/nt3/7t31cRgZlZpyZVHYCZjQ9OJmZWCCcTMyuEk4mZFcLJxMwK4WRiZoVwMjGzQjiZmFkhnEzMrBA7VB2AmbX23sOnx9Mb+pOOXblm0y0RcVTJIW3DycSsBp7a0M/Pbnl90rGTX/vr2SWHMyInE7NaCPpjoOogxuQxE7MaCGCASNpSSDpK0i8lPSjp/BE+/wZJt0u6W9IaSe9v1aaTiVlNDCT+14qkPuAK4H3A/sBJkvYfdthFwPURcTBwIvCPrdr1ZY5ZDQRBf3HLhbwdeDAiHgKQdB3wQWDdkFPCTvnrWcB/tGp0wvZMJL1QdQy9rsrvkaRLJZ1T1fmb4thT0n1VxwFtXebMlrSiaZs/rKnXAf/etP9o/l6zS4FTJD0K/AD4VKv43DMxq4EA+hPHQ4CnImJuh6c8Cfh6RPydpHcC35R0QMToo8ATomci6RRJP5e0WtKV+TUjki6XtFbSYkm7dCGOi/NBrzslfUvSOZKW5HGskHS/pHmSviPpAUmfLTmecyWdmb++XNKP8tdHSLqm6bjZku6SdHTJ8Vwo6VeS7gTelL9X2fenSZ+kq/KflVsl7SjpTEnr8sHJ68oOIIAtMZC0JXgM2L1p//X5e81OA64HiIi7gKnAmLecx30ykfQHwH8FDomIg4B+4CPAdGBFRLwZ+DHwmZLjmAccD7yVbOCr+S/H5vwvyT8BNwJnAAcAp0r6vRLDWgocmr+eC8yQNDl/74487tcANwGXRMRNZQUiaQ7ZQN9BwPuBeU2frur707APcEX+s/Is2b/j+cDBEXEgcHoXYmAgcUuwHNhH0l6SppB93xcNO+a3wLtg8HdoKvDkWI1OhMucdwFzgOWSAHYE1pN93xfmx1wNfKfkOA4BboyIjcBGSd9r+lzjH/JeYG1EPA4g6SGyvyBPlxTTSmCOpJ2ATcAqsqRyKHAmMBlYDJwRET8uKYaGQ4F/i4iXACQ1/3BX9f1peDgiVuevVwJ7AmuAayR9F/huyefPBmDTL3PGbitiq6RPArcAfcDXImKtpL8h+wO7CPhr4CpJnybrGJ0aLRaMngjJRMA3IuKCIW9KFw87rsqVtTflHweaXjf2S/s3iogtkh4GTgV+QvYLcjjwRuB+YCvZL897yXpvVank+zPC+SHr2e4IHA38MXAMcKGkt0TE1tIiCOgv8Cc0In5ANrDa/N4lTa/Xkf0BTDbuL3PI/rKeIOnVAJJ2lrQH2f/7CfkxJwN3lhzHMuAYSVMlzQA+UPL5Ui0FziG7rFlK1mW/O/8rFMDHgP0knVdyHHcAx+XjETPJfkl71SRg94i4HTiP7NbpjDJPmE1aK+wypxTjvmcSEeskXQTcKmkSsIXsmvtF4O3559aTjauUGcfyvOu+BniCrMv+XJnnTLQUuBC4KyJelLQxfw+AiOiXdBKwSNLzEdFy8tL2iIhVkhYC95D9eywv4zwF6QOuljSLrOf75Yh4ttxTin5U7ik6JNfN6R5JMyLiBUnTyP4Sz4+IVVXHZb3vgAOnxL/elPb83n5veHxlAbeG2zbueyY9ZkE+bXkq2TiOE4klCWBzj49KOJl0UUScXHUMVl8D0duXOU4mZjWQzYB1MjGzDgWiv8cvc3o7uhKN8PBTpXopnl6KBRxPw0AoaavKhE0mQE/9gNJb8fRSLOB4Bi9zUraq+DLHrBZEf/T23/7aJZNX7twXu72+r+N2dn1dH/sfOKXjSTaPPVjM2r1TJ+/ErGm7dT7pZ2vnM7qn9s1g1pRXdxxLbClmdvlUprGTdu44Hk2ZXEQ4TO2byaxXvKbjeH63ef1TEZH0tHoAW+j8575MtUsmu72+j6u/t2vVYQy64JiPVh3CEJOeeqbqEAZt/c8nqg5hiB12TVvdvVtu/u3/eyT12Aj3TMysIAO+NWxmncoGYN0zMbOO+TLHzAqQLUHgZGJmHQrE5vDdHDMrwECPX+b0dnRmBrw8AJuypUgoD3p5Xs1hdV4xoOXiT+6ZmNVAIPoLeu6mqTzou8kKcC2XtChf9zU7X8Snm47/FHBwq3bdMzGriQEmJW0JBsuDRsRmoFEedDQnAd9q1ah7JmY1EEE7t4ZnS1rRtL8gIhY07Y9UHvQdIzWUL76+F/CjVid1MjGrBbUzA7aI8qANJwI3RER/qwOdTMxqIIDNUdiva0p50IYTyao5tORkYlYDQaELHw2WByVLIieS1Y4aQtJ+wKuAu1IarTSZSHohIkotXmQ2XhT1bE5ieVDIksx1rcqCNrhnYlYDQbGT1lqVB833L22nza7dGpZ0iqSf55NgrszvdTcmx6yVtFhS0kIxZhNP2pKNVS7b2JVkIukPyMpvHhIRB5EVf/4IMJ2sW/VmssLYnxnl6+dLWiFpxTMbqqymalaNRs8kZatKty5z3gXMIZtpB1kV+fVkdZYX5sdcDXxnpC/O75EvAApZatGsjlw3JyOycpgXDHlTunjYcU4UZiOIEFsGenuIs1t9osXACZJeDSBp53xm3STghPyYk4E7uxSPWa1k65koaatKV1JdRKyTdBFwq6RJwBayiTAvAm/PP7eebFzFzLbhldYGRcRCXh4fafAcE7ME2QCsx0zMrABeUNrMOlbwdPpSOJmY1YQXlDazjkXAlgEnEzPrUHaZ42RiZgXwDFgz65hvDZtZQXyZY2YFqXKqfAonE7MayFandzIxsw4FYutAb9ca7u2LMDMbVORTw63Kg+bH/JmkdflKiNe2arN2PZPHHpzNBR/8aNVhDPrhrddVHcIQR+01Yi2lSkyaPr3qEIbo/8/1VYew3Yq8m5NSHlTSPsAFZKsjPtNYPmQstUsmZhNVgXdzBsuDAkhqlAdd13TMx4ErIuIZgIhomYl9mWNWB5E96JeykZcHbdrmD2ttpPKgrxt2zL7AvpKWSfqppKNaheieiVkNNFZaS1REedAdgH2Aw8gq/t0h6S0R8exYX2BmNVDgDNiU8qCPAj+LiC3Aw5J+RZZclo/WqC9zzGoggK0Dk5K2BIPlQSVNIavct2jYMd8l65UgaTbZZc9DYzXqnolZDRS5OFJiedBbgPdIWkdW5+rciHh6rHadTMxqosjp9K3Kg+b1hc/OtyROJmZ1EH5q2MwK4CUIzKwwTiZm1rFA9HsNWDMrgtczMbOOhQdgzawo4WRiZp1zRT8zK8iE7plIuhg4BXiS7JHnlcAHgLuBQ4HpwJ+TLcLyFmBhRFxUZkxmdTSh55lImgccD7wVmAysIksmAJsjYq6ks4AbgTnABuDXki4f/gxAvh7DfICpk3cqK2Sz3lWDBaXLvHF9CHBjRGyMiOeB7zV9rvGE4r3A2oh4PCI2kT2VuPuwdoiIBRExNyLmTtmht5YCNOuGILvMSdmqUtWYyab840DT68a+x3HMttH7A7Bl9kyWAcdImippBtlYiZltp4i0rSql9QIiYrmkRcAa4AmyS5rnyjqf2Xg3oe/mAJdFxKWSpgF3ACsj4qrGJyNiCbCkaf+wkuMxq6Ws1zGxk8kCSfsDU4FvRMSqks9nNm71+phJqckkIk4us32ziWRgoLeTSW8/02xmQLYEQZG3hluVB5V0qqQnJa3Ot79s1aZvw5rVRFE3alLKg+YWRsQnU9t1z8SsDqLQSWuD5UEjYjPQKA/aEScTs7qIxK21lPKgAMdLWiPpBknbzEwfzsnErCba6Jm0qjWc4nvAnhFxIHAb8I1WX9DWmImkScCMiPjddgRnZh1oY3Zrq1rDLcuDDnvY9p+BL7Q6acueiaRrJe0kaTpwH7BO0rmtvs7MihMBMTApaUvQsjyopNc27R4L3N+q0ZQz75/3RI4DfgjsBfy3lIjNrDhFPZsTEVuBRnnQ+4HrG+VBJR2bH3ampLWS7gHOBE5t1W7KZc5kSZPJksk/RMQWSRU+TmQ2QRX4W5dQHvQCskXLkqX0TK4EfkO2KtodkvYAPGZi1lXFTlorQ8ueSUR8Gfhy01uPSDq8vJBaxCOIvr6qTr+Now/p+PZ8oWLTI1WHMOj8dcurDmGILx7y7qpDGOrxNo/v8euBlAHYs/IBWEn6qqRVwBFdiM3MGoqdtFaKlMucj+UDsO8BXkU2+Pq5UqMys20VN2mtFCkDsI1U937gm/mob28/vmg2HvX4EgQpPZOVkm4lSya3SJpJtlarmXXTOOiZnAYcBDwUES9J+j3gL8oNy8yGCHq+Z5JyN2dA0sPAvpKmdiEmMxtBlYtFp2iZTPJFUc4im7+/GvhD4C58R8esu3o8maSMmZwFzAMeiYjDgYOBZ0uNysy2FUrbKpIyZrIxIjZKQtIrIuIXkt5UemRm9rIA9fhtj5Rk8qikVwLfBW6T9AzQO9MszSaEansdKVIGYD+Uv7xU0u3ALODmUqMys231+JjJqMlE0s4jvH1v/nEGsKGUiMxsZHVNJsBKsvCb+1aN/QB+v8S4zGy4uiaTiNirGwFIeiEiZnTjXGa1VYNJaylPDX9I0qym/VdKOq7csMxsOEXaVpWUeSafiYjnGjsR8SzwmdQTSDpX0pn568sl/Sh/fYSka5qOmy3pLklHp4dvNoH0+LM5KclkpGPaWdV+KXBo/nouMCNfBvJQ4A4ASa8BbgIuiYib2mjbbMIosmfSqjxo03HHSwpJY612D6QlkxWSviRp73z7EtngbKqVwBxJOwGbyKbizyVLJkuBycBi4H9GxG0jNSBpfqMGyJatL7VxarNxpKAZsE3lQd8H7A+cJGn/EY6bSTYD/mcp4aUkk08Bm4GFZGUENwJnpDQOEBFbgIfJVrf+CVkCORx4I9nK2FvJEs57x2hjQUTMjYi5k3eYlnpqs/Ej9RInrWeSWh70/wCfJ/udb6llMomIFyPi/PyXeV5E/K+IeDEp5JctBc4hu6xZCpwO3B0Rjf/9jwH7STqvzXbNJo4ulgeV9DZg93aGHbpVHnQp8Frgroh4gizTLW18MiL6gZOAIyR9oksxmdVKG2MmHZUHzSt3fgn463a+rq3yoNsrIhaTjY009vdtej0j/7iJMS51zCa87pUHnQkcACzJV2jdFVgk6diIWDFao11JJmbWGRX71PBgeVCyJHIicHLjk/lUkNmD55aWAOeMlUggbdLavpIWS7ov3z9Q0kXb9b9gZtuvoLs5ieVB25bSM7kKOJessh8RsUbStcBnt/ekZrYdulgedNj7h6W0mZJMpkXEz4dVt9ia0riZFafXK3ynJJOnJO1NnhclnUD7hQ3NrFPjIJmcASwgmwfyGNkEtFNKjcrMhqr4Ib4UKSutPQQcKWk6MCkini8/LDPbRt2TiaRLhu0DEBF/U1JMZjaC8bCgdPPU+anAB8huJ5mZDUq5zPm75n1Jl5Hdnzazbqr7Zc4IppFNvzWzbhkPA7CS7uXlnNgH7AJ4vMSs2+qeTMjGSBq2Ak/k03HNrJvqnEzyFZluiYj9uhSPmY1A9P5lzpgP+uXrjPxS0hu6FI+ZjSR/ajhlq0rKZc6rgLWSfk7TbeKI2O6nC81sO/R4zyQlmVxcehRm1to4SCbvj4gha7NK+jzw43JCMrOR1HrMJPfuEd57X9GBmFkLPV6Ea9SeiaT/AXwC+H1Ja5o+NRNYVnZgZtak4kSRYqzLnGuBHwL/F2iu+PV8RGwoNSoz20ZtH/TLF5V9jqwEhZlVbDyMmZhZLyhwzKRVrWFJp0u6V9JqSXeOVD50OCcTszoosDxoYq3hayPiLRFxEPAFsqJcY3IyMasBtbElaFlrOCJ+17Q7nYQ05SJcZnWRPmYyW1JzwawFEbGgaX+kWsPvGN6IpDOAs4EpwBGtTupkYlYTbQzAtioPmiQirgCukHQycBHw0bGO79pljqRLJZ3TrfOZjTsDiVtrrWoND3cdcFyrRj1mYlYH+UprKVuCwVrDkqaQ1Rpe1HyApH2ado8GHmjVaKnJRNKFkn4l6U7gTfl7SyRdLmmFpPslzZP0HUkPSHLJUbPRFHQ3J7HW8CclrZW0mmzcZMxLHChxzETSHLKMd1B+nlXAyvzTmyNirqSzgBuBOcAG4NeSLo+Ip4e1NR+YDzB1yqyyQjbraUVOWmtVazgizmq3zTJ7JocC/xYRL+W3mZq7UY3X9wJrI+LxiNgEPMTQazkAImJBRMyNiLmTd5hWYshmPayuD/qVbFP+caDpdWPfd5jMRjCRp9PfARwnaUdJM4FjSjyX2fhW4AzYspTWC4iIVZIWAvcA68lGkM1sO4gaPzVchIj4W+Bvh719WdPnlwBLmvYPKzMes1rr8cscj0+Y1YSit7OJk4lZHdR8pTUz6yG9fjfHycSsLpxMzKwI7pmYWedigt8aNrMCuWdiZp0Svswxs6J4nomZFcE9EzPrnCetmVlRfDenYNrST9/6Z6oO42WvmFJ1BEP07bJL1SEM+tz+86oOYYgFD/xr1SEMsdc2y4CNrdeTiReUNquDIBuATdkSJJQHPVvSOklrJC2WtEerNp1MzGqiqNXpE8uD3g3MjYgDgRvISoSOycnErC6KW2ktpTzo7RHxUr77U7LaOmNyMjGrgcaktcSeyey8lExjmz+suZHKg75ujNOfBvywVYy1G4A1m5DaGA+hoPKgAJJOAeYCf9LqWCcTs5oo8G5OUnlQSUcCFwJ/kpeiGZMvc8xqosvlQQ8GrgSOjYj1KY26Z2JWBwEMFDMFNiK2SmqUB+0DvtYoDwqsiIhFwBeBGcC3JQH8NiKOHbVRnEzM6qO75UGPbLdNJxOzmvCDfmZWDC9BYGZFcM/EzDqmABU0AFsWJxOzuvBTwxlJe0q6r1vnMxtvFJG0VcU9E7M6qMFKa92eAdsn6SpJayXdKmlHSWc2rZtwXZfjMauJxLVMJlDPZB/gpIj4uKTrgeOB84G9ImKTpFeO9EX5U4/zAab2zexasGa9pNfv5nS7Z/JwRKzOX68E9gTWANfkTyduHemLImJBRMyNiLlTJu3YnUjNek2P90y6nUyanzzsJ+sZHU226tPbgOWSPI5jNlyA+iNpq0rVTw1PAnaPiNuB84BZZA8Xmdlwxa20VoqqewF9wNWSZpEtJvXliHi24pjMelKVt31TdC2ZRMRvgAOa9i/r1rnNxgUnEzPrWNDzM2CdTMxqQFQ7uzWFk4lZXTiZmFnHAqjwtm+Kqm8Nm1miIh/0SygP+seSVknaKumElDadTMzqoqAZsInlQX8LnApcmxqeL3PMaqHQqfKD5UEB8gdsPwisGzxbNpUDKb1aj3smZnUQFPlsTrvlQZO4Z2JWF+nzTGZLWtG0vyAiFhQf0FBOJmY10cY8k1a1hpPKg7bLycSsDgLoL2wK7GB5ULIkciJwcqeNeszErBaKW2ktIrYCjfKg9wPXN8qDSjoWQNI8SY8CHwaulLS2VbuKHp9VN5ykJ4FHCmhqNvBUAe0UpZfi6aVYYPzGs0dE7JJy4Kypu8Yf7f7nSY3e/OAXV7a4zClF7S5zUr/5rUhaUcU3fDS9FE8vxQKOZ1CP/+GvXTIxm5ACcBEuM+tcQPT2GgQTOZmUft+9Tb0UTy/FAo6n6Ls5pZiwd3O6MYmnHUXGI+mF/ONukm5ocexfSZrWTiySDpP0/TbiWSJpu8cYxvO/VZsn9ur01rn84ay2RMR/RESrJz7/CpjW4hjrBU4mNpa8BvMvJF0j6X5JNzR6CpJ+I+nzklYBH5a0t6SbJa2UtFTSfvlxe0m6S9K9kj47rO378td9ki6TdF9ePfFTks4EdgNul3R7ftx78rZWSfq2pBn5+0flca4C/nSU/5dtzjHCMV+RtCKv6vi/m97/XFNlx8vy9z6ct3WPpDuK+Y7XlSv6WZo3AadFxDJJXwM+ATQW3H46It4GIGkxcHpEPCDpHcA/AkcAfw98JSL+RdIZo5xjPlnRs4MiYquknSNig6SzgcMj4ilJs4GLgCMj4kVJ5wFnS/oCcFV+rgeBhannGOGYC/Pz9gGLJR1INgvzQ8B+ERFNlR0vAd4bEY+NVu1xwghgwGMm1tq/R8Sy/PXVwH9p+txCgLyH8EfAtyWtBq4EXpsfcwjwrfz1N0c5x5HAlfnsRyJiwwjH/CHZ+hbL8nN8FNgD2I+sGuMDkc1yvLqDc/xZ3ru5G3hzfr7ngI3AVyX9KfBSfuwy4OuSPk5WFmVic8/EEgz/CWjefzH/OAl4NiIOSmxjewi4LSJOGvKmNNo522s8exbkHGBeRDwj6evA1LwX83bgXcAJZFO9j4iI0/Me2NHASklzIuLpImKppR6ftOaeSW94g6R35q9PBu4cfkBE/A54WNKHAZR5a/7pZWQPawF8ZJRz3Ab890b51aZLkOeBRjX4nwKHSHpjfsx0SfsCvwD2lLR3ftyQZJNwjoadyJLjc5JeQ7bSV6PXNSsifgB8Gnhr/v7eEfGziLgEeJKhT7pOLBFEf3/SVhUnk97wS+AMSfcDrwK+MspxHwFOk3QPsJZsdSyAs/Kvv5fRF7n5Z7Kl+NbkX994SnQBcLOk2yPiSbKl+r4laQ1wF9k4xkay8ZCb8kuU9W2eA4CIuIfs8uYXZMsBNi7tZgLfz895J3B2/v4X80Hl+4CfAPeMct6JYSDStorU7kG/8UbSnsD3I+KAFofaBDZrh13inTM/2PpA4JZnv+oH/cxsFBE9fzfHyaRiw2swm42qx68inEzMaiLcMzGzzlU7hySFk4lZHQRQ4W3fFE4mZjUQQHhxJDPrWHhxJDMrSK/3TDxpzawGJN1Mtip+iqci4qgy4xmJk4mZFcLP5phZIZxMzKwQTiZmVggnEzMrhJOJmRXCycTMCuFkYmaFcDIxs0I4mZhZIf4/FVuyrg+xoKQAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 288x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.matshow(CMsubsumedfrom2)\n",
    "plt.xticks(range(max(truth)+1), names)\n",
    "plt.yticks(range(max(truth)+1), names)\n",
    "plt.xlabel('predicted class')\n",
    "plt.ylabel('true class')\n",
    "plt.colorbar()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Metric values"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Import the `proclam.Metric` superclass and the metric subclass you want to test."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "from proclam.metrics import metric\n",
    "from proclam.metrics import logloss\n",
    "from proclam.metrics import brier"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Probabilistic metrics"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's try out the [logloss metric](https://en.wikipedia.org/wiki/Loss_functions_for_classification#Cross_entropy_loss).  The logloss metric is a sum over $LL_{ij} = -y_{ij}\\ln[p_{ij}]$ for predicted probabilities $p$ and true class indicators $y$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [],
   "source": [
    "# weights = np.array([0.3,0.15,0.1,0.15,0.3])\n",
    "weights = np.array([0.0,0.0,1.0,0.0,0.0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "proclam implementation of log-loss: 1.9467692362663755\n",
      "proclam implementation of log-loss: 0.5873053661481069\n",
      "proclam implementation of log-loss: 0.4985535342418646\n",
      "proclam implementation of log-loss: 3.2716260697307087\n",
      "proclam implementation of log-loss: 0.8391912085925096\n",
      "proclam implementation of log-loss: 0.5795189847083857\n",
      "proclam implementation of log-loss: 4.44547996411847\n"
     ]
    }
   ],
   "source": [
    "for candidate in [predictionB, predictionC, predictionG, predictionH, predictionI, predictionJ, predictionJ1]:\n",
    "    D = proclam.metrics.LogLoss()\n",
    "    performance = D.evaluate(candidate, truth, averaging=weights)\n",
    "    print('proclam implementation of log-loss: '+str(performance))\n",
    "    #alternative = skl.metrics.log_loss(truth, candidate, normalize=True)\n",
    "    #print('scikit-learn implementation of log-loss: '+str(alternative))\n",
    "\n",
    "#D = proclam.metrics.LogLoss()\n",
    "#performance = D.evaluate(predictionF, truth_snphotcc, averaging=weights)\n",
    "#print('proclam implementation of log-loss: '+str(performance))\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "tunnel logloss = 1.3831099142162573\n",
      "broadbrush logloss = 2.147514136377367\n",
      "cruise logloss = 1.8925640255539373\n",
      "subsumedto2 logloss = 0.8697035708801626\n",
      "subsumedfrom logloss = 1.2313487087147357\n"
     ]
    }
   ],
   "source": [
    "test_cases = {'tunnel': CMtunnel, 'broadbrush': CMbroadbrush, 'cruise': CMcruise, 'subsumedto2': CMsubsumedto2, 'subsumedfrom': CMsubsumedfrom2}\n",
    "LL_metric = proclam.metrics.LogLoss()\n",
    "for candidate_cm in test_cases.keys():\n",
    "    competitor = proclam.classifiers.FromCM(seed=None, scheme=candidate_cm)\n",
    "    pred_competitor = competitor.classify(test_cases[candidate_cm], truth, other=False)\n",
    "    metval = LL_metric.evaluate(pred_competitor, truth, averaging='per_class')\n",
    "    print(candidate_cm+' logloss = '+str(metval))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "tunnel logloss = 0.45996840473970635\n",
      "broadbrush logloss = 3.3173356388856186\n",
      "cruise logloss = 0.7884312373010852\n",
      "subsumedto logloss = 0.5772789861655064\n",
      "subsumedfrom logloss = 4.401321122622816\n"
     ]
    }
   ],
   "source": [
    "test_cases = {'tunnel': CMtunnel, 'broadbrush': CMbroadbrush, 'cruise': CMcruise, 'subsumedto': CMsubsumedto2,'subsumedfrom': CMsubsumedfrom2}\n",
    "LL_metric = proclam.metrics.LogLoss()\n",
    "for candidate_cm in test_cases.keys():\n",
    "    competitor = proclam.classifiers.FromCM(seed=None, scheme=candidate_cm)\n",
    "    pred_competitor = competitor.classify(test_cases[candidate_cm], truth, other=False)\n",
    "    metval = LL_metric.evaluate(pred_competitor, truth, averaging=weights)\n",
    "    print(candidate_cm+' logloss = '+str(metval))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's try the [Brier score](https://en.wikipedia.org/wiki/Brier_score#Original_definition_by_Brier) this time.  The multi-class Brier score is $BS = \\frac{1}{N}\\sum\\limits _{t=1}^{N}\\sum\\limits _{i=1}^{R}(f_{ti}-o_{ti})^2$ for $N$ objects, $R$ classes, predicted probabilities $f$, and $o_{i}$ of 1 for true class $i$ and 0 for other true classes.\n",
    "\n",
    "First we create an instantiation of the metric.  Then, we calculate the metric value.  For binary classes, we can compare to the implementation in `scikit-learn`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "proclam implementation of Brier score: 0.17481523141779132\n",
      "proclam implementation of Brier score: 0.03349662221027182\n",
      "proclam implementation of Brier score: 0.13994060524323776\n",
      "proclam implementation of Brier score: 0.16442263732364013\n",
      "proclam implementation of Brier score: 0.19760600980918508\n",
      "proclam implementation of Brier score: 0.08131350128672857\n",
      "proclam implementation of Brier score: 0.08047203028278058\n"
     ]
    }
   ],
   "source": [
    "for candidate in [predictionB, predictionC, predictionG, predictionH, predictionI, predictionJ, predictionJ1]:\n",
    "    E = proclam.metrics.Brier()\n",
    "    performance = E.evaluate(candidate, truth)\n",
    "    print('proclam implementation of Brier score: '+str(performance))\n",
    "    \n",
    "    if M_classes == 2:\n",
    "        skl_truth = proclam.metrics.util.truth_reformatter(truth).T[0]\n",
    "        alternative = skl.metrics.brier_score_loss(skl_truth, prediction.T[0])\n",
    "        print('scikit-learn implementation: '+str(alternative))\n",
    "        \n",
    "# E = proclam.metrics.brier.Brier()\n",
    "# performance = E.evaluate(predictionF, truth_snphotcc)\n",
    "# print('proclam implementation of Brier score: '+str(performance))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "tunnel Brier = 0.13835879922919994\n",
      "broadbrush Brier = 0.16647095421593264\n",
      "cruise Brier = 0.19722890971661755\n",
      "subsumedto Brier = 0.07719957610858942\n",
      "subsumedfrom Brier = 0.0838596231599401\n"
     ]
    }
   ],
   "source": [
    "test_cases = {'tunnel': CMtunnel, 'broadbrush': CMbroadbrush, 'cruise': CMcruise, 'subsumedto': CMsubsumedto2, 'subsumedfrom': CMsubsumedfrom2}\n",
    "B_metric = proclam.metrics.Brier()\n",
    "for candidate_cm in test_cases.keys():\n",
    "    competitor = proclam.classifiers.FromCM(seed=None, scheme=candidate_cm)\n",
    "    pred_competitor = competitor.classify(test_cases[candidate_cm], truth, other=False)\n",
    "    metval = B_metric.evaluate(pred_competitor, truth, averaging='per_class')\n",
    "    print(candidate_cm+' Brier = '+str(metval))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "tunnel Brier = 0.03829109430194667\n",
      "broadbrush Brier = 0.2519580848801754\n",
      "cruise Brier = 0.07399778110647401\n",
      "subsumedto Brier = 0.0487604318197624\n",
      "subsumedfrom Brier = 0.28825253455457156\n"
     ]
    }
   ],
   "source": [
    "test_cases = {'tunnel': CMtunnel, 'broadbrush': CMbroadbrush, 'cruise': CMcruise, 'subsumedto': CMsubsumedto2, 'subsumedfrom': CMsubsumedfrom2}\n",
    "B_metric = proclam.metrics.Brier()\n",
    "for candidate_cm in test_cases.keys():\n",
    "    competitor = proclam.classifiers.FromCM(seed=None, scheme=candidate_cm)\n",
    "    pred_competitor = competitor.classify(test_cases[candidate_cm], truth, other=False)\n",
    "    metval = B_metric.evaluate(pred_competitor, truth, averaging=weights)\n",
    "    print(candidate_cm+' Brier = '+str(metval))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Deterministic metrics\n",
    "\n",
    "Let's check that reducing the probabilities to class point estimates actually does what we want; the one based on a good confusion matrix should do better than the random guesser.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### True/False positive/negative rates\n",
    "\n",
    "Let's compare `proclam`'s calculation of the standard rates to that of `pycm`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "rates(TPR=array([0.98113208, 1.        , 0.94117647, 0.95548961, 0.96296296]), FPR=array([0.00527983, 0.00740741, 0.01256281, 0.01357466, 0.01540832]), FNR=array([0.01886792, 0.        , 0.05882353, 0.04451039, 0.03703704]), TNR=array([0.99472017, 0.99259259, 0.98743719, 0.98642534, 0.98459168]), TP=array([ 52.,  55., 192., 322., 338.]), FP=array([ 5.,  7., 10.,  9., 10.]), FN=array([ 1.,  0., 12., 15., 13.]), TN=array([942., 938., 786., 654., 639.]))\n"
     ]
    }
   ],
   "source": [
    "detC = proclam.metrics.util.prob_to_det(predictionC)\n",
    "cmC = proclam.metrics.util.det_to_cm(detC, truth)\n",
    "print(proclam.metrics.util.cm_to_rate(cmC, vb=True))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "rates(TPR={0: 0.9811320754716981, 1: 1.0, 2: 0.9411764705882353, 3: 0.9554896142433235, 4: 0.9629629629629629}, FPR={0: 0.0052798310454065245, 1: 0.007407407407407418, 2: 0.012562814070351758, 3: 0.013574660633484115, 4: 0.015408320493066285}, FNR={0: 0.018867924528301883, 1: 0.0, 2: 0.05882352941176472, 3: 0.04451038575667654, 4: 0.03703703703703709}, TNR={0: 0.9947201689545935, 1: 0.9925925925925926, 2: 0.9874371859296482, 3: 0.9864253393665159, 4: 0.9845916795069337}, TP={0: 52, 1: 55, 2: 192, 3: 322, 4: 338}, FP={0: 5, 1: 7, 2: 10, 3: 9, 4: 10}, FN={0: 1, 1: 0, 2: 12, 3: 15, 4: 13}, TN={0: 942, 1: 938, 2: 786, 3: 654, 4: 639})\n"
     ]
    }
   ],
   "source": [
    "compare = ConfusionMatrix(truth, detC)\n",
    "print(proclam.metrics.util.RateMatrix(TPR=compare.TPR, FPR=compare.FPR, FNR=compare.FNR, TNR=compare.TNR, TP=compare.TP, FP=compare.FP, FN=compare.FN, TN=compare.TN))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### ROC AUC\n",
    "\n",
    "The ROC is the true positive rate as a function of the false positive rate, where the rates are calculated from a confusion matrix derived by deterministically assigning classes based on a series of threshold values in probability."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ROC AUC for prediction B = 0.493\n",
      "ROC AUC for prediction C = 0.995\n",
      "ROC AUC for prediction B with weird weights = 0.498\n"
     ]
    }
   ],
   "source": [
    "ROC_metric = proclam.metrics.ROC()\n",
    "rocB = ROC_metric.evaluate(predictionB,truth, 0.1, averaging='per_class', vb=False)\n",
    "rocC = ROC_metric.evaluate(predictionC,truth, 0.1)\n",
    "print('ROC AUC for prediction B = %.3f'%rocB)\n",
    "print('ROC AUC for prediction C = %.3f'%rocC)\n",
    "rocB = ROC_metric.evaluate(predictionB,truth,0.1, averaging=[0.1,0.3,0.2,0.2,0.2])\n",
    "print('ROC AUC for prediction B with weird weights = %.3f'%rocB)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's compare the ROC AUC calculated by `proclam` to that of `scikit-learn`.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [],
   "source": [
    "# this runs if you change line 78 in roc.py to return auc_class instead of auc_allclass\n",
    "# truth_score = proclam.metrics.util.det_to_prob(truth)\n",
    "# from_proclam, proclam_auc, from_skl = [],[],[]\n",
    "# for m in range(M_classes):\n",
    "#     fpr, tpr, thresholds = skl.metrics.roc_curve(truth_score.T[m].T.astype(int), predictionB.T[m].T)\n",
    "#     proclam_says = ROC_metric.evaluate(predictionB, truth, thresholds, vb=True)[m]\n",
    "#     proclam_auc.append(ROC_metric.evaluate(predictionB, truth, thresholds, vb=False)[m])\n",
    "#     i = np.argsort(proclam_says[0])\n",
    "#     new_auc = proclam.metrics.util.auc(proclam_says[0][i], proclam_says[1][i])\n",
    "#     print('proclam: '+str(new_auc))\n",
    "#     from_proclam.append(proclam.metrics.util.auc(proclam_says[0], proclam_says[1]))\n",
    "#     skl_says = proclam.metrics.util.auc(fpr, tpr)\n",
    "#     print('skl: '+str(skl_says))\n",
    "#     from_skl.append(skl_says)\n",
    "# print('proclam says %.3f'%np.mean(proclam_auc)+' scikit-learn says %.3f'%np.mean(from_skl))\n",
    "# print('proclam says '+str(proclam_auc)+' scikit-learn says '+str(from_skl))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### PRC AUC\n",
    "\n",
    "The precision is the number of correctly classified positives divided by the number all items classified as positive, whereas the recall is the number of correctly classified positives divided by the number of items whose true class was positive.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/aimalz/Code/proclam/proclam/metrics/util.py:455: RuntimeWarning: invalid value encountered in double_scalars\n",
      "  p = np.asarray(TP / (TP + FP))\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Precision/Recall AUC for prediction B = 0.203\n",
      "Precision/Recall AUC for prediction C = 0.937\n"
     ]
    }
   ],
   "source": [
    "metric = proclam.metrics.PRC()\n",
    "prB = metric.evaluate(predictionB,truth, 0.01)\n",
    "prC = metric.evaluate(predictionC,truth, 0.01)\n",
    "\n",
    "print('Precision/Recall AUC for prediction B = %.3f'%prB)\n",
    "print('Precision/Recall AUC for prediction C = %.3f'%prC)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [],
   "source": [
    "# compare = ConfusionMatrix(truth, proclam.metrics.util.prob_to_det(predictionC))\n",
    "# print((compare.TPR, compare.PPV))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### F-score"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [],
   "source": [
    "prC = metric.evaluate(predictionC,truth, 0.01, vb=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "f1 score for prediction B = 0.184\n",
      "f1 score for prediction C = 0.953\n"
     ]
    }
   ],
   "source": [
    "metric = proclam.metrics.F1()\n",
    "fB = metric.evaluate(predictionB,truth)\n",
    "fC = metric.evaluate(predictionC,truth)\n",
    "\n",
    "print('f1 score for prediction B = %.3f'%fB)\n",
    "print('f1 score for prediction C = %.3f'%fC)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [],
   "source": [
    "# compare = ConfusionMatrix(truth, proclam.metrics.util.prob_to_det(predictionB))\n",
    "# print((compare.F1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Matthews Correlation Coefficient"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "MCC for prediction B = 0.009\n",
      "MCC for prediction C = 0.942\n"
     ]
    }
   ],
   "source": [
    "metric = proclam.metrics.MCC()\n",
    "mccB = metric.evaluate(predictionB, truth)\n",
    "mccC = metric.evaluate(predictionC,truth)\n",
    "\n",
    "print('MCC for prediction B = %.3f'%mccB)\n",
    "print('MCC for prediction C = %.3f'%mccC)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [],
   "source": [
    "# compare = ConfusionMatrix(truth, proclam.metrics.util.prob_to_det(predictionC))\n",
    "# print((compare.MCC))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Accuracy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Accuracy for prediction B = 0.207\n",
      "Accuracy for prediction C = 0.968\n"
     ]
    }
   ],
   "source": [
    "metric = proclam.metrics.Accuracy()\n",
    "accB = metric.evaluate(predictionB, truth)\n",
    "accC = metric.evaluate(predictionC,truth)\n",
    "\n",
    "print('Accuracy for prediction B = %.3f'%accB)\n",
    "print('Accuracy for prediction C = %.3f'%accC)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{0: 0.785, 1: 0.769, 2: 0.675, 3: 0.601, 4: 0.596}\n"
     ]
    }
   ],
   "source": [
    "# compare = ConfusionMatrix(truth, proclam.metrics.util.prob_to_det(predictionB))\n",
    "# print((compare.ACC))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Next steps"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Once there are more simulators, classifiers, and metrics, we'll loop over tests and plot comparisons.  Stay tuned for more!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Z = proclam.classifiers.guess.Guess(seed=None)\n",
    "# predictionZ = Z.classify(M_classes, truth, other=False)\n",
    "# # print(prediction)\n",
    "# predictionZ\n",
    "\n",
    "# for candidate in [predictionZ]:\n",
    "#     D = proclam.metrics.LogLoss()\n",
    "#     performance = D.evaluate(candidate, truth, averaging='per_item')\n",
    "#     print('proclam implementation of log-loss: '+str(performance))\n",
    "# #     alternative = skl.metrics.log_loss(truth, candidate, normalize=True)\n",
    "# #     print('scikit-learn implementation of log-loss: '+str(alternative))\n",
    "\n",
    "# D = proclam.metrics.LogLoss()\n",
    "# performance = D.evaluate(predictionF, truth_snphotcc, averaging='per_item')\n",
    "# print('proclam implementation of log-loss: '+str(performance))\n",
    "\n",
    "# print(type(predictionZ))\n",
    "# print(predictionZ.shape)\n",
    "\n",
    "# rows = [0, 1, 30, 50, 62, 60]\n",
    "# cols = [0, 2, 4, 1, 0, 0]\n",
    "\n",
    "# ## fill with fake zeros\n",
    "# predictionZ_0 = predictionZ\n",
    "# predictionZ_0[rows,cols] = 0\n",
    "# predictionZ_0\n",
    "\n",
    "# for candidate in [predictionZ_0]:\n",
    "#     D = proclam.metrics.LogLoss()\n",
    "#     performance = D.evaluate(candidate, truth, averaging='per_item')\n",
    "#     print('proclam implementation of log-loss: '+str(performance))\n",
    "# #     alternative = skl.metrics.log_loss(truth, candidate, normalize=True)\n",
    "# #     print('scikit-learn implementation of log-loss: '+str(alternative))\n",
    "\n",
    "# D = proclam.metrics.LogLoss()\n",
    "# performance = D.evaluate(predictionF, truth_snphotcc, averaging='per_item')\n",
    "# print('proclam implementation of log-loss: '+str(performance))\n",
    "\n",
    "# print(np.where(predictionZ_0 < 0.000003))\n",
    "\n",
    "# predictionZ_1 = predictionZ\n",
    "# predictionZ_1[rows,cols] = 1\n",
    "# predictionZ_1\n",
    "\n",
    "# for candidate in [predictionZ_1]:\n",
    "#     D = proclam.metrics.LogLoss()\n",
    "#     performance = D.evaluate(candidate, truth, averaging='per_item')\n",
    "#     print('proclam implementation of log-loss: '+str(performance))\n",
    "# #     alternative = skl.metrics.log_loss(truth, candidate, normalize=True)\n",
    "# #     print('scikit-learn implementation of log-loss: '+str(alternative))\n",
    "\n",
    "# D = proclam.metrics.LogLoss()\n",
    "# performance = D.evaluate(predictionF, truth_snphotcc, averaging='per_item')\n",
    "# print('proclam implementation of log-loss: '+str(performance))\n",
    "\n",
    "# np.where(predictionZ_1 > 0.9995)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# for i in range(1, 100000):\n",
    "#     print(i)\n",
    "#     Z = proclam.classifiers.guess.Guess(seed=i)\n",
    "#     predictionZ = Z.classify(M_classes, truth, other=False)\n",
    "#     # print(prediction)\n",
    "#     predictionZ\n",
    "\n",
    "#     for candidate in [predictionZ]:\n",
    "#         D = proclam.metrics.LogLoss()\n",
    "#         performance = D.evaluate(candidate, truth, averaging='per_item')\n",
    "#         # print('proclam implementation of log-loss: '+str(performance))\n",
    "#     #     alternative = skl.metrics.log_loss(truth, candidate, normalize=True)\n",
    "#     #     print('scikit-learn implementation of log-loss: '+str(alternative))\n",
    "\n",
    "#     D = proclam.metrics.LogLoss()\n",
    "#     performance = D.evaluate(predictionF, truth_snphotcc, averaging='per_item')\n",
    "#     # print('proclam implementation of log-loss: '+str(performance))\n",
    "\n",
    "#     # print(type(predictionZ))\n",
    "#     # print(predictionZ.shape)\n",
    "\n",
    "#     rows = [0, 1, 30, 50, 62, 60]\n",
    "#     cols = [0, 2, 4, 1, 0, 0]\n",
    "\n",
    "#     ## fill with fake zeros\n",
    "#     predictionZ_0 = predictionZ\n",
    "#     predictionZ_0[rows,cols] = 0\n",
    "#     predictionZ_0\n",
    "\n",
    "#     for candidate in [predictionZ_0]:\n",
    "#         D = proclam.metrics.LogLoss()\n",
    "#         performance = D.evaluate(candidate, truth, averaging='per_item')\n",
    "#         # print('proclam implementation of log-loss: '+str(performance))\n",
    "#     #     alternative = skl.metrics.log_loss(truth, candidate, normalize=True)\n",
    "#     #     print('scikit-learn implementation of log-loss: '+str(alternative))\n",
    "\n",
    "#     D = proclam.metrics.LogLoss()\n",
    "#     performance = D.evaluate(predictionF, truth_snphotcc, averaging='per_item')\n",
    "#     # print('proclam implementation of log-loss: '+str(performance))\n",
    "\n",
    "#     np.where(predictionZ_0 < 0.0005)\n",
    "\n",
    "#     predictionZ_1 = predictionZ\n",
    "#     predictionZ_1[rows,cols] = 1\n",
    "#     predictionZ_1\n",
    "\n",
    "#     for candidate in [predictionZ_1]:\n",
    "#         D = proclam.metrics.LogLoss()\n",
    "#         performance = D.evaluate(candidate, truth, averaging='per_item')\n",
    "#         # print('proclam implementation of log-loss: '+str(performance))\n",
    "#     #     alternative = skl.metrics.log_loss(truth, candidate, normalize=True)\n",
    "#     #     print('scikit-learn implementation of log-loss: '+str(alternative))\n",
    "\n",
    "#     D = proclam.metrics.LogLoss()\n",
    "#     performance = D.evaluate(predictionF, truth_snphotcc, averaging='per_item')\n",
    "#     # print('proclam implementation of log-loss: '+str(performance))\n",
    "\n",
    "#     np.where(predictionZ_1 == 1)\n",
    "#     print(\"ITER :\\n\",i)\n",
    "#     i += 1\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "proclam (Python 3)",
   "language": "python",
   "name": "proclam_3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}