pipeline_sandbox.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Testing framework \"pipeline\"\n",
    "\n",
    "_Alex Malz (NYU)_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import string\n",
    "import itertools\n",
    "import random\n",
    "import numpy as np\n",
    "import scipy.stats as sct\n",
    "import sklearn as skl\n",
    "from sklearn import metrics\n",
    "import pandas as pd\n",
    "import os\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "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": null,
   "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": null,
   "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": null,
   "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)"
   ]
  },
  {
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": null,
   "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": null,
   "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": null,
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": null,
   "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": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# print(snphotcc_info)\n",
    "# print(truth_snphotcc)\n",
    "# print(full)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": null,
   "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": null,
   "metadata": {},
   "outputs": [],
   "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 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": null,
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": null,
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": null,
   "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": null,
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": null,
   "metadata": {},
   "outputs": [],
   "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": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# from proclam.metrics import util as pmu\n",
    "\n",
    "# for candidate in [predictionB, predictionC]:\n",
    "#     intermediate = pmu.prob_to_det(candidate)\n",
    "#     rates = pmu.det_to_rate(intermediate, truth, vb=False)\n",
    "\n",
    "#     print('true positives '+str(rates.TPR))\n",
    "#     print('false negatives'+str(rates.FNR))\n",
    "    \n",
    "# intermediate = pmu.prob_to_det(predictionF)\n",
    "# rates = pmu.det_to_rate(intermediate, truth_snphotcc, vb=False)\n",
    "# print('true positives '+str(rates.TPR))\n",
    "# print('false negatives'+str(rates.FNR))"
   ]
  },
  {
   "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": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "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
}