diff --git a/metrics_evaluation.ipynb b/metrics_evaluation.ipynb index 61d8c61..1e5857c 100644 --- a/metrics_evaluation.ipynb +++ b/metrics_evaluation.ipynb @@ -16,8 +16,10 @@ }, { "cell_type": "code", - "execution_count": null, - "metadata": {}, + "execution_count": 1, + "metadata": { + "scrolled": true + }, "outputs": [], "source": [ "# import string\n", @@ -36,8 +38,10 @@ }, { "cell_type": "code", - "execution_count": null, - "metadata": {}, + "execution_count": 2, + "metadata": { + "scrolled": true + }, "outputs": [], "source": [ "import matplotlib as mpl\n", @@ -64,8 +68,10 @@ }, { "cell_type": "code", - "execution_count": null, - "metadata": {}, + "execution_count": 3, + "metadata": { + "scrolled": true + }, "outputs": [], "source": [ "metricslist = ['Brier', 'LogLoss']\n", @@ -90,8 +96,10 @@ }, { "cell_type": "code", - "execution_count": null, - "metadata": {}, + "execution_count": 4, + "metadata": { + "scrolled": true + }, "outputs": [], "source": [ "mystery = {}\n", @@ -108,8 +116,10 @@ }, { "cell_type": "code", - "execution_count": null, - "metadata": {}, + "execution_count": 5, + "metadata": { + "scrolled": true + }, "outputs": [], "source": [ "snphotcc = {}\n", @@ -131,8 +141,10 @@ }, { "cell_type": "code", - "execution_count": null, - "metadata": {}, + "execution_count": 6, + "metadata": { + "scrolled": true + }, "outputs": [], "source": [ "plasticc = {}\n", @@ -141,38 +153,16 @@ ] }, { - "cell_type": "code", - "execution_count": null, + "cell_type": "markdown", "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", - "# for i in range(len(snphotcc_names)):\n", - "# name = old_snphotcc_names[i]\n", - "# fileloc = dirname+'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", - "# truth = full['Type'] - 1\n", - "# snphotcc_truth_table = proclam.metrics.util.det_to_prob(truth)\n", - "# fileloc = 'examples/'+name+'/truth_table_'+name+'.csv'\n", - "# with open(fileloc, 'wb') as truth_place:\n", - "# np.savetxt(fileloc, snphotcc_truth_table, delimiter=' ')\n", - " \n", - "# probs = full[['1', '2', '3']]\n", - "# fileloc = 'examples/'+name+'/predicted_prob_'+name+'.csv'\n", - "# probs.to_csv(fileloc, sep=' ', index=False, header=True)" - ] + "source": [] }, { "cell_type": "code", - "execution_count": null, - "metadata": {}, + "execution_count": 7, + "metadata": { + "scrolled": true + }, "outputs": [], "source": [ "# more_names = snphotcc_names\n", @@ -183,8 +173,10 @@ }, { "cell_type": "code", - "execution_count": null, - "metadata": {}, + "execution_count": 8, + "metadata": { + "scrolled": true + }, "outputs": [], "source": [ "def make_class_pairs(data_info_dict):\n", @@ -201,11 +193,21 @@ }, { "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], + "execution_count": 9, + "metadata": { + "scrolled": true + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "{'label': 'ProClaM', 'names': ['Idealized', 'Guess', 'Tunnel', 'Broadbrush', 'Cruise', 'SubsumedTo', 'SubsumedFrom'], 'dirname': 'examples/ProClaM/', 'classifications': ['Idealized/predicted_prob_Idealized.csv', 'Guess/predicted_prob_Guess.csv', 'Tunnel/predicted_prob_Tunnel.csv', 'Broadbrush/predicted_prob_Broadbrush.csv', 'Cruise/predicted_prob_Cruise.csv', 'SubsumedTo/predicted_prob_SubsumedTo.csv', 'SubsumedFrom/predicted_prob_SubsumedFrom.csv'], 'truth_tables': ['Idealized/truth_table_Idealized.csv', 'Guess/truth_table_Guess.csv', 'Tunnel/truth_table_Tunnel.csv', 'Broadbrush/truth_table_Broadbrush.csv', 'Cruise/truth_table_Cruise.csv', 'SubsumedTo/truth_table_SubsumedTo.csv', 'SubsumedFrom/truth_table_SubsumedFrom.csv']}\n" + ] + } + ], "source": [ - "for dataset in [mystery, snphotcc, plasticc]:\n", + "for dataset in [ plasticc]: #mystery, snphotcc,\n", " dataset = make_file_locs(dataset)\n", " dataset['class_pairs'] = make_class_pairs(dataset)" ] @@ -221,8 +223,10 @@ }, { "cell_type": "code", - "execution_count": null, - "metadata": {}, + "execution_count": 10, + "metadata": { + "scrolled": true + }, "outputs": [], "source": [ "def plot_cm(probs, truth, name, loc=''):\n", @@ -234,14 +238,17 @@ " plt.ylabel('true class')\n", " plt.colorbar()\n", " plt.title(name)\n", - " plt.savefig(loc+name+'_cm.png')\n", - " plt.close()" + " #plt.savefig(loc+name+'_cm.png')\n", + " plt.show()\n", + " #plt.close()" ] }, { "cell_type": "code", - "execution_count": null, - "metadata": {}, + "execution_count": 11, + "metadata": { + "scrolled": true + }, "outputs": [], "source": [] }, @@ -263,6 +270,7 @@ " nobj_truth = np.shape(truth_values)[0]\n", " nclass_truth = np.shape(truth_values)[1]\n", " tvec = np.where(truth_values==1)[1]\n", + " print(tvec)\n", "# if nclass_truth!= nclass:\n", "# print('Truth table of size %i x %i and prob matrix of size %i x %i do not match up in size'%(nobj,nclass,nobj_truth,nclass_truth))\n", "# else:\n", @@ -274,8 +282,10 @@ }, { "cell_type": "code", - "execution_count": null, - "metadata": {}, + "execution_count": 12, + "metadata": { + "scrolled": true + }, "outputs": [], "source": [ "def make_patch_spines_invisible(ax):\n", @@ -324,6 +334,7 @@ " plt.legend(handles, metric_names)\n", " plt.suptitle(title)\n", " plt.savefig(fileloc)\n", + " plt.show()\n", " return" ] }, @@ -336,13 +347,28 @@ }, { "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], + "execution_count": 13, + "metadata": { + "scrolled": false + }, + "outputs": [ + { + "ename": "KeyError", + "evalue": "'class_pairs'", + "output_type": "error", + "traceback": [ + "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", + "\u001b[0;31mKeyError\u001b[0m Traceback (most recent call last)", + "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0mdataset\u001b[0m \u001b[0;32min\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0mmystery\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0msnphotcc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mplasticc\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2\u001b[0m \u001b[0mdata\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mempty\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mlen\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmetricslist\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mlen\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdataset\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'names'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 3\u001b[0;31m \u001b[0;32mfor\u001b[0m \u001b[0mcc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpair\u001b[0m \u001b[0;32min\u001b[0m \u001b[0menumerate\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdataset\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'class_pairs'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 4\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mpair\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 5\u001b[0m \u001b[0mprobm\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtruthv\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mread_class_pairs\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mpair\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mdataset\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mcc\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;31m#loc=dataset['dirname'], title=dataset['label']+' '+dataset['names'][cc])\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;31mKeyError\u001b[0m: 'class_pairs'" + ] + } + ], "source": [ "for dataset in [mystery, snphotcc, plasticc]:\n", " data = np.empty((len(metricslist), len(dataset['names'])))\n", " for cc, pair in enumerate(dataset['class_pairs']):\n", + " print(pair)\n", " probm, truthv = read_class_pairs(pair, dataset, cc)#loc=dataset['dirname'], title=dataset['label']+' '+dataset['names'][cc])\n", "# plot_cm(probm, truthv, str(cc), loc='./sandbox/')\n", " det = proclam.metrics.util.prob_to_det(probm)\n", @@ -361,45 +387,6 @@ "# metric_plot(dataset, metricslist, markerlist, colors)" ] }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# more_data = np.empty((len(metricslist), len(more_names)))\n", - "# for cc, pair in enumerate(more_class_pairs):\n", - "# probm, truthv = read_class_pairs(pair, dirname)\n", - "# for count, metric in enumerate(metricslist):\n", - "# D = getattr(proclam.metrics, metric)()\n", - "# hm = D.evaluate(probm, truthv)\n", - "# more_data[count][cc] = hm" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# metric_plot(more_names, metricslist, more_data, markerlist, colors, title='SNPhotCC', fileloc=dirname+'snphotccdata.png')" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# data = np.empty((len(metricslist), len(names)))\n", - "# for cc, pair in enumerate(class_pairs):\n", - "# probm, truthv = read_class_pairs(pair, dirname)\n", - "# for count, metric in enumerate(metricslist):\n", - "# D = getattr(proclam.metrics, metric)()\n", - "# hm = D.evaluate(probm, truthv)\n", - "# data[count][cc] = hm" - ] - }, { "cell_type": "code", "execution_count": null, @@ -407,9 +394,7 @@ "scrolled": true }, "outputs": [], - "source": [ - "# metric_plot(names, metricslist, data, markerlist, colors, title='Mystery Dataset', fileloc=dirname+'mysterydata.png')" - ] + "source": [] }, { "cell_type": "code", @@ -443,7 +428,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.6.5" + "version": "3.6.8" } }, "nbformat": 4, diff --git a/metrics_plots.py b/metrics_plots.py new file mode 100644 index 0000000..be1d7be --- /dev/null +++ b/metrics_plots.py @@ -0,0 +1,207 @@ +# import string +# import itertools +# import random +# import os +# import csv + +import numpy as np +import pandas as pd + +import proclam +from proclam import * +import matplotlib as mpl +import pylab as pl +mpl.use('Agg') +mpl.rcParams['text.usetex'] = False +mpl.rcParams['mathtext.rm'] = 'serif' +mpl.rcParams['font.family'] = 'serif' +mpl.rcParams['font.serif'] = 'Times New Roman' +mpl.rcParams['axes.titlesize'] = 16 +mpl.rcParams['axes.labelsize'] = 14 +mpl.rcParams['savefig.dpi'] = 250 +mpl.rcParams['savefig.format'] = 'pdf' +mpl.rcParams['savefig.bbox'] = 'tight' +import matplotlib.pyplot as plt +metricslist = ['Brier', 'LogLoss'] +colors = ['teal', 'magenta'] +dirname = 'examples/' +markerlist = ['d', 'o', 's', '*'] +plasticc = {} +plasticc['label'] = 'plasticc' +#plasticc['names'] = ['Submission_alpha_0.5_190516_1756', 'submission_40_avocado', 'submission_probe99_40_avocado'] +plasticc['names'] = ['3_MajorTom'] #'2_MikeSilogram' ] #'1_Kyle'] + +#, '2_MikeSilogram', '3_MajorTom'] + + + +list = [6,15,16,42,52,53,62,64,65,67,88,90,92,95,99] +itemlist=['uLens-Point', 'TDE', 'EBE', 'SNCC-II', 'MIRA', 'SNCC-Ibc', 'KN', 'Mdwarf', 'SNIa-91bg', 'AGN', 'SNIa-normal', 'RRlyrae', 'SLSN-I', 'Other'] + +choices=['All'] + +#, 'uLens-Point', 'TDE', 'EBE', 'SNCC-II', 'MIRA', 'SNCC-Ibc', 'KN', 'Mdwarf', 'SNIa-91bg', 'AGN', 'SNIa-normal', 'RRlyrae', 'SLSN-I', 'Other'] + +# 90 SNIa-normal +#62 SNCC-Ibc +#42 SNCC-II +#67 SNIa-91bg +#52 SNIa-x +#64 KN +#95 SLSN-I +#99 Other +#15 TDE +#88 AGN +#92 RRlyrae +#65 Mdwarf +#16 EBE +#53 MIRA +#6 uLens-Point + + +#Idealized', 'Guess', 'Tunnel', 'Broadbrush', 'Cruise', 'SubsumedTo', 'SubsumedFrom'] +def make_class_pairs(data_info_dict): + return zip(data_info_dict['classifications'], data_info_dict['truth_tables']) + +def make_file_locs(data_info_dict): + names = data_info_dict['names'] + data_info_dict['dirname'] = dirname + data_info_dict['label'] + '/' +# data_info_dict['classifications'] = ['%s/predicted_prob_%s.csv'%(name, name) for name in names] +# data_info_dict['truth_tables'] = ['%s/truth_table_%s.csv'%(name, name) for name in names] + data_info_dict['classifications'] = ['%s/%s.csv'%(name, name) for name in names] + data_info_dict['truth_tables'] = ['%s/%s_truth.csv'%(name, name) for name in names] + print(data_info_dict) + return data_info_dict + +def plot_cm(probs, truth, name, loc=''): + print(np.shape(probs), np.shape(truth), 'checking sizes of probs and truth') + cm = proclam.metrics.util.prob_to_cm(probs, truth) + pl.clf() + plt.matshow(cm.T, vmin=0., vmax=1.) +# plt.xticks(range(max(truth)+1), names) +# plt.yticks(range(max(truth)+1), names) + plt.xlabel('predicted class') + plt.ylabel('true class') + plt.colorbar() + plt.title(name) + plt.savefig(loc+name+'_cm.png') + #plt.show() + #plt.close() + + + +def read_class_pairs(pair, dataset, cc):#loc='', title=''): + loc=dataset['dirname'] + title=dataset['label']+' '+ dataset['names'][cc] + clfile = pair[0] + truthfile = pair[1] + print(clfile, truthfile) + prob_mat = pd.read_csv(loc+clfile)#, delim_whitespace=True) + nobj = np.shape(prob_mat)[0] + nclass = np.shape(prob_mat)[1]-1 #since they have object ID as an element + cols=prob_mat.columns.tolist() + + objid = prob_mat[cols[0]] + pmat=np.array(prob_mat[cols[1:]]) + + truth_values = pd.read_csv(loc+truthfile) #, delim_whitespace=True) + nobj_truth = np.shape(truth_values)[0] + nclass_truth = np.shape(truth_values)[1]-1 + truvals = np.array(truth_values[cols[1:]]) + tvec = np.where(truvals==1)[1] + + +# pmat = prob_mat[:,1:] + plot_cm(pmat, tvec, title, loc=loc+dataset['names'][cc]+'/') + return pmat, tvec + + +def make_patch_spines_invisible(ax): + ax.set_frame_on(True) + ax.patch.set_visible(False) + for sp in ax.spines.values(): + sp.set_visible(False) + +def per_metric_helper(ax, n, data, metric_names, codes, shapes, colors): + plot_n = n+1 + in_x = np.arange(len(codes)) + ax_n = ax + n_factor = 0.1 * (plot_n - 2) + if plot_n>1: + ax_n = ax.twinx() + rot_ang = 270 + label_space = 15. + else: + rot_ang = 90 + label_space = 0. + if plot_n>2: + ax_n.spines["right"].set_position(("axes", 1. + 0.1 * (plot_n-1))) + make_patch_spines_invisible(ax_n) + ax_n.spines["right"].set_visible(True) + handle = ax_n.scatter(in_x+n_factor*np.ones_like(data[n]), data[n], marker=shapes[n], s=10, color=colors[n], label=metric_names[n]) + ax_n.set_ylabel(metric_names[n], rotation=rot_ang, fontsize=14, labelpad=label_space) +# ax_n.set_ylim(0.9 * min(data[n]), 1.1 * max(data[n])) + return(ax, ax_n, handle) + +def metric_plot(dataset, metric_names, shapes, colors, choice): + codes = dataset['names'] + data = dataset['results'] + title = dataset['label']+' results focusing on: '+str(choice) + fileloc = dataset['dirname']+dataset['label']+'_'+str(choice)+'_results.png' + xs = np.arange(len(codes)) + pl.clf() + fig, ax = plt.subplots() + fig.subplots_adjust(right=1.) + handles = [] + for n in range(len(metric_names)): + (ax, ax_n, handle) = per_metric_helper(ax, n, data, metric_names, codes, shapes, colors) + handles.append(handle) + plt.xticks(xs, codes) + for tick in ax.get_xticklabels(): + tick.set_rotation(90) + plt.xlabel('Classifiers', fontsize=14) + leg=plt.legend(handles, metric_names, numpoints=1, loc='lower right') + #leg.draw_frame(False) + plt.suptitle(title) + plt.savefig(fileloc) + #plt.show() + return + + +for dataset in [ plasticc]: #mystery, snphotcc, + dataset = make_file_locs(dataset) + dataset['class_pairs'] = make_class_pairs(dataset) + + +for choice in choices: + pl.clf() + + if choice=='All': + print('ignoring weights for %s'%choice) + #weights= None + weights=np.ones(len(list)) + print(len(list), 'length of list') + else: + weights= np.zeros(len(list)) #1e-5*np.ones(len(list)) + ind = itemlist.index(choice) + print(itemlist) + print(ind, 'check ind', choice) + print(itemlist[ind], choice, 'checking choice') + weights[ind]=1.0 + for dataset in [plasticc]: + data = np.empty((len(metricslist), len(dataset['names']))) + + for cc, pair in enumerate(dataset['class_pairs']): + probm, truthv = read_class_pairs(pair, dataset, cc) + + for count, metric in enumerate(metricslist): + print(weights, 'checking huh') + print(len(weights), 'how many weights?') + D = getattr(proclam.metrics, metric, weights)() + hm = D.evaluate(probm, truthv, weights) + data[count][cc] = hm + dataset['results'] = data + + metric_plot(dataset, metricslist, markerlist, colors, choice) + +#----------------------------------------------------- diff --git a/organize_csv.py b/organize_csv.py new file mode 100644 index 0000000..a98a4f8 --- /dev/null +++ b/organize_csv.py @@ -0,0 +1,48 @@ +import pylab as pl +import pandas as pd +import numpy as np +dir = 'examples/plasticc/' +list = ['2_MikeSilogram', '3_MajorTom'] + +kyle = pd.read_csv('examples/plasticc/1_Kyle/1_Kyle.csv') +kylecols = kyle.columns.tolist() + + +matstruc = kyle.copy() +matstrucdat = np.array(matstruc) +indexlist = np.zeros(len(kylecols)) +print(kylecols) +print(matstrucdat[0,:]) + + +truth = pd.read_csv('examples/plasticc/1_Kyle/1_Kyle_truth.csv') + +for file in list: + name = dir+file+'/'+file+'.csv' + print(name) + mat = pd.read_csv(name) + matdat = np.array(mat) + + cols = mat.columns.tolist() + index = mat.index + #print(np.shape(mat)) + print(cols) + print(matdat[0,:], 'test before') + + + for i in range(len(kylecols)): +# df[col] = df[col].replace(findL, replaceL) + for j in range(len(cols)): + if cols[j]==kylecols[i]: + indexlist[i]=j + matstrucdat[:,i] = matdat[:,j] + print(i, j, cols[j], kylecols[i]) + + matstruc.rename(columns={"A": "a", "B": "c"}) + print(matstrucdat[0,:], 'test after') + + newname = dir+file+'/'+file+'_reordered.csv' + newnametruth = dir+file+'/'+file+'_reordered_truth.csv' + + matstruc.to_csv(newname, index=False) + truth.to_csv(newnametruth) diff --git a/paper/kaggle-run.ipynb b/paper/kaggle-run.ipynb deleted file mode 100644 index 26948f9..0000000 --- a/paper/kaggle-run.ipynb +++ /dev/null @@ -1,695 +0,0 @@ -{ - "cells": [ - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "# The Photometric LSST Astronomical Time-series Classification Challenge (PLAsTiCC): code that runs the performance metric\n", - "\n", - "*Alex Malz (NYU)*, *Renee Hlozek (U. Toronto)*, *Tarek Alam (UCL)*, *Anita Bahmanyar (U. Toronto)*, *Rahul Biswas (U. Stockholm)*, *Rafael Martinez-Galarza (Harvard)*, *Gautham Narayan (STScI)*" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "import numpy as np\n", - "import pandas as pd\n", - "\n", - "# This is the code available on GitHub for calculating metrics,\n", - "# as well as performing other diagnostics on probability tables." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "\n", - "def make_class_pairs(data_info_dict):\n", - " \"\"\"\n", - " Paris the paths to classifier output and truth tables for each classifier.\n", - "\n", - " Parameters\n", - " ----------\n", - " data_info_dict: dictionary\n", - " \n", - " Returns\n", - " -------\n", - " data_info_dict: dictionary\n", - " updated keywords: class_pairs, dict - classifier: [path_to_class_output, path_to_truth_tables] \n", - " \"\"\"\n", - " \n", - " for name in data_info_dict['names']:\n", - " data_info_dict['class_pairs'][name] = [data_info_dict['classifications'][name], data_info_dict['truth_tables'][name]]\n", - " \n", - " return(data_info_dict['class_pairs'])\n", - " \n", - "def make_file_locs(data_info_dict):\n", - " \"\"\"\n", - " Set paths to data directory, classifier output and truth tables.\n", - "\n", - " Parameters\n", - " ----------\n", - " data_info_dict: dictionary \n", - " \n", - " Returns\n", - " -------\n", - " data_info_dict: dictionary\n", - " updated keywords: dirname - data directory, str\n", - " classifications, dict - classifier: path to classifier output, str\n", - " truth_tables, dict - classifier: path to truth tables - str\n", - " \"\"\"\n", - " \n", - " # get the names of classifiers to be considered\n", - " names = data_info_dict['names']\n", - " \n", - " # set data directory\n", - " data_info_dict['dirname'] = topdir + data_info_dict['label'] + '/'\n", - "\n", - " for name in names:\n", - " # get the path to classifier output\n", - " data_info_dict['classifications'][name] = '%s/predicted_prob_%s.csv'%(name, name)\n", - " \n", - " # get the path to truth table\n", - " data_info_dict['truth_tables'][name] = '%s/truth_table_%s.csv'%(name, name)\n", - " \n", - " return data_info_dict\n", - "\n", - "def process_strings(dataset, cc):\n", - " \"\"\"\n", - " Get info on directory name and classifier.\n", - "\n", - " Parameters\n", - " ----------\n", - " dataset: dictionary \n", - " cc: classifier name, str\n", - " \n", - " Returns\n", - " -------\n", - " loc: data directory, str\n", - " text: version label, str\n", - " \"\"\"\n", - " \n", - " loc = dataset['dirname']\n", - " text = dataset['label'] + ' ' + dataset['names'][cc]\n", - " \n", - " return loc, text\n", - "\n", - "def just_read_class_pairs(pair, dataset, cc):\n", - " \"\"\"\n", - " Reads predicted probabilities and truth table.\n", - "\n", - " Parameters\n", - " ----------\n", - " pair: list of str - [path_to_classifier_output, path_to_truth_table]\n", - " dataset: dictionary\n", - " cc: classifier name, str\n", - " \n", - " Returns\n", - " -------\n", - " prob_mat: probability matrix (output from classifier)\n", - " tvec: truth vector\n", - " \"\"\"\n", - " \n", - " loc, text = process_strings(dataset, cc)\n", - " clfile = pair[0]\n", - " truthfile = pair[1]\n", - " \n", - " # read classifier output\n", - " prob_mat = pd.read_csv(loc + clfile, delim_whitespace=True).values\n", - " \n", - " # read truth table\n", - " truth_values = pd.read_csv(loc + truthfile, delim_whitespace=True).values\n", - " tvec = np.where(truth_values==1)[1]\n", - " \n", - " return prob_mat, tvec" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# Build dictionary to store classification results\n", - "mystery = {}\n", - "mystery['label'] = 'Unknown'\n", - "mystery['names'] = ['RandomForest', 'KNeighbors', 'MLPNeuralNet']\n", - "mystery['classifications'] = {}\n", - "mystery['truth_tables'] = {}\n", - "mystery['class_pairs'] = {}\n", - "mystery['probs'] = {}\n", - "mystery['truth'] = {}" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# Read classifier output and truth tables\n", - "topdir = '../examples/'\n", - "mystery = make_file_locs(mystery)\n", - "mystery['class_pairs'] = make_class_pairs(mystery)\n", - "for nm, name in enumerate(mystery['names']):\n", - " probm, truthv = just_read_class_pairs(mystery['class_pairs'][name], mystery, nm)\n", - " mystery['probs'][name] = probm\n", - " mystery['truth'][name] = truthv\n", - "M_classes = np.shape(probm)[-1]\n", - "\n", - "# we need the class labels in the dataset in a consistently sorted order \n", - "# and will assume the weights of the weightvec correspond to this order\n", - "class_labels = sorted(np.unique(mystery['truth']['RandomForest']))" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "topdir = '../examples/'\n", - "mystery = make_file_locs(mystery)\n", - "mystery['class_pairs'] = make_class_pairs(mystery)\n", - "mystery['class_pairs']" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Method (Metric)\n", - "======\n", - "\n", - "The log-loss is defined as\n", - "\\begin{eqnarray*}\n", - "L &=& -\\sum_{m=1}^{M}\\frac{w_{m}}{N_{m}}\\sum_{n=1}^{N_{m}}\\ln[p_{n}(m | m)]\n", - "\\end{eqnarray*}\n", - "\n", - "We calculate the metric within each class $m$ by taking an average of its value $-\\ln[p_{n}(m | m)]$ for each true member $n$ of the class. Then we weight the metrics for each class by an arbitrary weight $w_{m}$ and take a weighted average of the per-class metrics to produce a global scalar metric $L$." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "\"\"\"\n", - "This is all from proclam but copied here so no one has to install it.\n", - "\"\"\"\n", - "\n", - "import numpy as np\n", - "import sys\n", - "import collections\n", - "\n", - "\"\"\"\n", - "Utility functions for PLAsTiCC metrics\n", - "\"\"\"\n", - "\n", - "# from __future__ import absolute_import, division\n", - "# __all__ = ['sanitize_predictions',\n", - "# 'weight_sum', 'averager', 'check_weights',\n", - "# 'det_to_prob',\n", - "# 'prob_to_det',\n", - "# 'det_to_cm', 'prob_to_cm',\n", - "# 'cm_to_rate', 'det_to_rate', 'prob_to_rate']\n", - "\n", - "RateMatrix = collections.namedtuple('rates', 'TPR FPR FNR TNR')\n", - "\n", - "def sanitize_predictions(predictions, epsilon=sys.float_info.epsilon):\n", - " \"\"\"\n", - " Replaces 0 and 1 with 0+epsilon, 1-epsilon\n", - "\n", - " Parameters\n", - " ----------\n", - " predictions: numpy.ndarray, float\n", - " N*M matrix of probabilities per object, may have 0 or 1 values\n", - " epsilon: float\n", - " small placeholder number, defaults to floating point precision\n", - "\n", - " Returns\n", - " -------\n", - " predictions: numpy.ndarray, float\n", - " N*M matrix of probabilities per object, no 0 or 1 values\n", - " \"\"\"\n", - " assert epsilon > 0. and epsilon < 0.0005\n", - " mask1 = (predictions < epsilon)\n", - " mask2 = (predictions > 1.0 - epsilon)\n", - "\n", - " predictions[mask1] = epsilon\n", - " predictions[mask2] = 1.0 - epsilon\n", - " predictions = predictions / np.sum(predictions, axis=1)[:, np.newaxis]\n", - " return predictions\n", - "\n", - "def det_to_prob(dets, prediction=None):\n", - " \"\"\"\n", - " Reformats vector of class assignments into matrix with 1 at true/assigned class and zero elsewhere\n", - "\n", - " Parameters\n", - " ----------\n", - " dets: numpy.ndarray, int\n", - " vector of classes\n", - " prediction: numpy.ndarray, float, optional\n", - " predicted class probabilities\n", - "\n", - " Returns\n", - " -------\n", - " probs: numpy.ndarray, float\n", - " matrix with 1 at input classes and 0 elsewhere\n", - "\n", - " Notes\n", - " -----\n", - " det_to_prob formerly truth_reformatter\n", - " Does not yet handle number of classes in truth not matching number of classes in prediction, i.e. for having \"other\" class or secret classes not in training set. The prediction keyword is a kludge to enable this but should be replaced.\n", - " \"\"\"\n", - " N = len(dets)\n", - " indices = range(N)\n", - "\n", - " if prediction is None:\n", - " prediction_shape = (N, int(np.max(dets) + 1))\n", - " else:\n", - " prediction, dets = np.asarray(prediction), np.asarray(dets)\n", - " prediction_shape = np.shape(prediction)\n", - "\n", - " probs = np.zeros(prediction_shape)\n", - " probs[indices, dets] = 1.\n", - "\n", - " return probs\n", - "\n", - "def prob_to_det(probs):\n", - " \"\"\"\n", - " Converts probabilistic classifications to deterministic classifications by assigning the class with highest probability\n", - "\n", - " Parameters\n", - " ----------\n", - " probs: numpy.ndarray, float\n", - " N * M matrix of class probabilities\n", - "\n", - " Returns\n", - " -------\n", - " dets: numpy.ndarray, int\n", - " maximum probability classes\n", - " \"\"\"\n", - " dets = np.argmax(probs, axis=1)\n", - "\n", - " return dets\n", - "\n", - "def det_to_cm(dets, truth, per_class_norm=True, vb=False):\n", - " \"\"\"\n", - " Converts deterministic classifications and truth into confusion matrix\n", - "\n", - " Parameters\n", - " ----------\n", - " dets: numpy.ndarray, int\n", - " assigned classes\n", - " truth: numpy.ndarray, int\n", - " true classes\n", - " per_class_norm: boolean, optional\n", - " equal weight per class if True, equal weight per object if False\n", - " vb: boolean, optional\n", - " if True, print cm\n", - "\n", - " Returns\n", - " -------\n", - " cm: numpy.ndarray, int\n", - " confusion matrix\n", - "\n", - " Notes\n", - " -----\n", - " I need to fix the norm keyword all around to enable more options, like normed output vs. not.\n", - " \"\"\"\n", - " pred_classes, pred_counts = np.unique(dets, return_counts=True)\n", - " true_classes, true_counts = np.unique(truth, return_counts=True)\n", - " if vb: print((pred_classes, pred_counts), (true_classes, true_counts))\n", - "\n", - " M = max(max(pred_classes), max(true_classes)) + 1\n", - "\n", - " cm = np.zeros((M, M), dtype=float)\n", - " # print((np.shape(dets), np.shape(truth)))\n", - " coords = np.array(list(zip(dets, truth)))\n", - " indices, index_counts = np.unique(coords, axis=0, return_counts=True)\n", - " # if vb: print(indices, index_counts)\n", - " indices = indices.T\n", - " # if vb: print(np.shape(indices))\n", - " cm[indices[0], indices[1]] = index_counts\n", - " if vb: print(cm)\n", - "\n", - " if per_class_norm:\n", - " # print(type(cm))\n", - " # print(type(true_counts))\n", - " # cm = cm / true_counts\n", - " # cm /= true_counts[:, np.newaxis] #\n", - " cm = cm / true_counts[np.newaxis, :]\n", - "\n", - " if vb: print(cm)\n", - "\n", - " return cm\n", - "\n", - "def prob_to_cm(probs, truth, per_class_norm=True, vb=False):\n", - " \"\"\"\n", - " Turns probabilistic classifications into confusion matrix by taking maximum probability as deterministic class\n", - "\n", - " Parameters\n", - " ----------\n", - " probs: numpy.ndarray, float\n", - " N * M matrix of class probabilities\n", - " truth: numpy.ndarray, int\n", - " N-dimensional vector of true classes\n", - " per_class_norm: boolean, optional\n", - " equal weight per class if True, equal weight per object if False\n", - " vb: boolean, optional\n", - " if True, print cm\n", - "\n", - " Returns\n", - " -------\n", - " cm: numpy.ndarray, int\n", - " confusion matrix\n", - " \"\"\"\n", - " dets = prob_to_det(probs)\n", - "\n", - " cm = det_to_cm(dets, truth, per_class_norm=per_class_norm, vb=vb)\n", - "\n", - " return cm\n", - "\n", - "def cm_to_rate(cm, vb=False):\n", - " \"\"\"\n", - " Turns a confusion matrix into true/false positive/negative rates\n", - "\n", - " Parameters\n", - " ----------\n", - " cm: numpy.ndarray, int or float\n", - " confusion matrix, first axis is predictions, second axis is truth\n", - " vb: boolean, optional\n", - " print progress to stdout?\n", - "\n", - " Returns\n", - " -------\n", - " rates: named tuple, float\n", - " RateMatrix named tuple\n", - "\n", - " Notes\n", - " -----\n", - " BROKEN!\n", - " This can be done with a mask to weight the classes differently here.\n", - " \"\"\"\n", - " if vb: print(cm)\n", - " diag = np.diag(cm)\n", - " if vb: print(diag)\n", - "\n", - " TP = np.sum(diag)\n", - " FN = np.sum(np.sum(cm, axis=0) - diag)\n", - " FP = np.sum(np.sum(cm, axis=1) - diag)\n", - " TN = np.sum(cm) - TP\n", - " if vb: print((TP, FN, FP, TN))\n", - "\n", - " T = TP + TN\n", - " F = FP + FN\n", - " P = TP + FP\n", - " N = TN + FN\n", - " if vb: print((T, F, P, N))\n", - "\n", - " TPR = TP / P\n", - " FPR = FP / N\n", - " FNR = FN / P\n", - " TNR = TN / N\n", - "\n", - " rates = RateMatrix(TPR=TPR, FPR=FPR, FNR=FNR, TNR=TNR)\n", - " if vb: print(rates)\n", - "\n", - " return rates\n", - "\n", - "def det_to_rate(dets, truth, per_class_norm=True, vb=False):\n", - " cm = det_to_cm(dets, truth, per_class_norm=per_class_norm, vb=vb)\n", - " rates = cm_to_rate(cm, vb=vb)\n", - " return rates\n", - "\n", - "def prob_to_rate(probs, truth, per_class_norm=True, vb=False):\n", - " cm = prob_to_cm(probs, truth, per_class_norm=per_class_norm, vb=vb)\n", - " rates = cm_to_rate(cm, vb=vb)\n", - " return rates\n", - "\n", - "def weight_sum(per_class_metrics, weight_vector, norm=True):\n", - " \"\"\"\n", - " Calculates the weighted metric\n", - "\n", - " Parameters\n", - " ----------\n", - " per_class_metrics: numpy.float\n", - " the scores separated by class (a list of arrays)\n", - " weight_vector: numpy.ndarray floar\n", - " The array of weights per class\n", - " norm: boolean, optional\n", - "\n", - " Returns\n", - " -------\n", - " weight_sum: np.float\n", - " The weighted metric\n", - " \"\"\"\n", - " weight_sum = np.dot(weight_vector, per_class_metrics)\n", - "\n", - " return weight_sum\n", - "\n", - "def check_weights(avg_info, M, truth=None):\n", - " \"\"\"\n", - " Converts standard weighting schemes to weight vectors for weight_sum\n", - "\n", - " Parameters\n", - " ----------\n", - " avg_info: str or numpy.ndarray, float\n", - " keyword about how to calculate weighted average metric\n", - " M: int\n", - " number of classes\n", - " truth: numpy.ndarray, int, optional\n", - " true class assignments\n", - "\n", - " Returns\n", - " -------\n", - " weights: numpy.ndarray, float\n", - " relative weights per class\n", - " \"\"\"\n", - " if type(avg_info) != str:\n", - " avg_info = np.asarray(avg_info)\n", - " weights = avg_info / np.sum(avg_info)\n", - " assert(np.isclose(sum(weights), 1.))\n", - " elif avg_info == 'per_class':\n", - " weights = np.ones(M) / float(M)\n", - " elif avg_info == 'per_item':\n", - " classes, counts = np.unique(truth, return_counts=True)\n", - " weights = np.zeros(M)\n", - " weights[classes] = counts / float(len(truth))\n", - " assert len(weights) == M\n", - " return weights\n", - "\n", - "def averager(per_object_metrics, truth, M):\n", - " \"\"\"\n", - " Creates a list with the metrics per object, separated by class\n", - " \"\"\"\n", - " group_metric = per_object_metrics\n", - " class_metric = np.empty(M)\n", - " for m in range(M):\n", - " true_indices = np.where(truth == m)[0]\n", - " how_many_in_class = len(true_indices)\n", - " try:\n", - " assert(how_many_in_class > 0)\n", - " per_class_metric = group_metric[true_indices]\n", - " # assert(~np.all(np.isnan(per_class_metric)))\n", - " class_metric[m] = np.average(per_class_metric)\n", - " except AssertionError:\n", - " class_metric[m] = 0.\n", - " # print((m, how_many_in_class, class_metric[m]))\n", - " return class_metric\n", - "\n", - "\"\"\"\n", - "A superclass for metrics\n", - "\"\"\"\n", - "class Metric(object):\n", - "\n", - " def __init__(self, scheme=None, **kwargs):\n", - " \"\"\"\n", - " An object that evaluates a function of the true classes and class probabilities\n", - "\n", - " Parameters\n", - " ----------\n", - " scheme: string\n", - " the name of the metric\n", - " \"\"\"\n", - " self.scheme = scheme\n", - "\n", - " def evaluate(self, prediction, truth, weights=None, **kwds):\n", - " \"\"\"\n", - " Evaluates a function of the truth and prediction\n", - "\n", - " Parameters\n", - " ----------\n", - " prediction: numpy.ndarray, float\n", - " predicted class probabilities\n", - " truth: numpy.ndarray, int\n", - " true classes\n", - " weights: numpy.ndarray, float\n", - " per-class weights\n", - "\n", - " Returns\n", - " -------\n", - " metric: float\n", - " value of the metric\n", - " \"\"\"\n", - " print('No metric specified: returning true positive rate based on maximum value')\n", - "\n", - " return # metric\n", - "\n", - "\"\"\"\n", - "A metric subclass for the log-loss\n", - "\"\"\"\n", - "class LogLoss(Metric):\n", - " def __init__(self, scheme=None):\n", - " \"\"\"\n", - " An object that evaluates the log-loss metric\n", - "\n", - " Parameters\n", - " ----------\n", - " scheme: string\n", - " the name of the metric\n", - " \"\"\"\n", - " super(LogLoss, self).__init__(scheme)\n", - " self.scheme = scheme\n", - "\n", - " def evaluate(self, prediction, truth, averaging='per_class'):\n", - " \"\"\"\n", - " Evaluates the log-loss\n", - "\n", - " Parameters\n", - " ----------\n", - " prediction: numpy.ndarray, float\n", - " predicted class probabilities\n", - " truth: numpy.ndarray, int\n", - " true classes\n", - " averaging: string or numpy.ndarray, float\n", - " 'per_class' weights classes equally, other keywords possible\n", - " vector assumed to be class weights\n", - "\n", - " Returns\n", - " -------\n", - " logloss: float\n", - " value of the metric\n", - "\n", - " Notes\n", - " -----\n", - " This uses the natural log.\n", - " \"\"\"\n", - " prediction, truth = np.asarray(prediction), np.asarray(truth)\n", - " prediction_shape = np.shape(prediction)\n", - " (N, M) = prediction_shape\n", - "\n", - " weights = check_weights(averaging, M, truth=truth)\n", - " truth_mask = det_to_prob(truth, prediction)\n", - "\n", - " prediction = sanitize_predictions(prediction)\n", - "\n", - " log_prob = np.log(prediction)\n", - " logloss_each = -1. * np.sum(truth_mask * log_prob, axis=1)[:, np.newaxis]\n", - "\n", - " # use a better structure for checking keyword support\n", - " class_logloss = averager(logloss_each, truth, M)\n", - "\n", - " logloss = weight_sum(class_logloss, weight_vector=weights)\n", - "\n", - " assert(~np.isnan(logloss))\n", - "\n", - " return logloss\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# This is how you run the metric with a random weight vector.\n", - "\n", - "metric = 'LogLoss'\n", - "weightvec = np.ones(M_classes)\n", - "\n", - "# dummy example for SNPhotCC demo data\n", - "special_classes = (1, 10, 11, 12)\n", - "\n", - "# we should be using this for the PLAsTiCC data\n", - "# special_clases = (51, 99)\n", - "\n", - "mask = np.array([True if classname in special_classes else False for classname in class_labels])\n", - "weightvec[mask] = 2\n", - "weightvec = weightvec / sum(weightvec)\n", - "name = np.random.choice(mystery['names'])\n", - "probm = mystery['probs'][name]\n", - "truthv = mystery['truth'][name]\n", - "LL = LogLoss()\n", - "val = LL.evaluate(prediction=probm, truth=truthv, averaging=weightvec)\n", - "print(name+' with weights '+str(weightvec)+' has '+metric+' = '+str(val))" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Acknowledgments\n", - "===============\n", - "\n", - "The DESC acknowledges ongoing support from the Institut National de Physique Nucleaire et de Physique des Particules in France; the Science & Technology Facilities Council in the United Kingdom; and the Department of Energy, the National Science Foundation, and the LSST Corporation in the United States.\n", - "\n", - "DESC uses resources of the IN2P3 Computing Center (CC-IN2P3--Lyon/Villeurbanne - France) funded by the Centre National de la Recherche Scientifique; the National Energy Research Scientific Computing Center, a DOE Office of Science User Facility supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC02-05CH11231; STFC DiRAC HPC Facilities, funded by UK BIS National E-infrastructure capital grants; and the UK particle physics grid, supported by the GridPP Collaboration.\n", - "\n", - "This work was performed in part under DOE Contract DE-AC02-76SF00515." - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Contributions\n", - "=======\n", - "\n", - "Alex Malz: conceptualization, data curation, formal analysis, investigation, methodology, project administration, software, supervision, validation, visualization, writing - original draft\n", - "\n", - "Renee Hlozek: data curation, formal analysis, funding acquisition, investigation, project administration, software, supervision, validation, visualization, writing - original draft\n", - "\n", - "Tarek Alam: investigation, software, validation\n", - "\n", - "Anita Bahmanyar: formal analysis, investigation, methodology, software, writing - original draft\n", - "\n", - "Rahul Biswas: conceptualization, methodology, software\n", - "\n", - "Rafael Martinez-Galarza: data curation, software, visualization\n", - "\n", - "Gautham Narayan: data curation, formal analysis" - ] - } - ], - "metadata": { - "anaconda-cloud": {}, - "celltoolbar": "Tags", - "kernelspec": { - "display_name": "Python 3", - "language": "python3", - "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.8" - } - }, - "nbformat": 4, - "nbformat_minor": 1 -} diff --git a/paper/main_paperfigs.ipynb b/paper/main_paperfigs.ipynb index 476ae20..75f2287 100644 --- a/paper/main_paperfigs.ipynb +++ b/paper/main_paperfigs.ipynb @@ -2,7 +2,7 @@ "cells": [ { "cell_type": "code", - "execution_count": null, + "execution_count": 1, "metadata": { "tags": [ "hideme" @@ -18,7 +18,7 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 2, "metadata": { "tags": [ "hideme" @@ -51,7 +51,7 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 3, "metadata": {}, "outputs": [], "source": [ @@ -121,14 +121,33 @@ { "cell_type": "code", "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": 5, "metadata": { "tags": [ "hideme" ] }, - "outputs": [], - "source": [ - "import proclam\n", + "outputs": [ + { + "ename": "ModuleNotFoundError", + "evalue": "No module named 'proclam'", + "output_type": "error", + "traceback": [ + "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", + "\u001b[0;31mModuleNotFoundError\u001b[0m Traceback (most recent call last)", + "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0;31m#import proclam\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0;32mfrom\u001b[0m \u001b[0mproclam\u001b[0m \u001b[0;32mimport\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m", + "\u001b[0;31mModuleNotFoundError\u001b[0m: No module named 'proclam'" + ] + } + ], + "source": [ + "#import proclam\n", "from proclam import *" ] }, @@ -1913,9 +1932,9 @@ "anaconda-cloud": {}, "celltoolbar": "Tags", "kernelspec": { - "display_name": "ProClaM (Python 3)", + "display_name": "Python 3", "language": "python", - "name": "proclam_3" + "name": "python3" }, "language_info": { "codemirror_mode": { @@ -1927,7 +1946,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.6.5" + "version": "3.6.8" } }, "nbformat": 4, diff --git a/proclam/metrics/brier.py b/proclam/metrics/brier.py index b06f923..f7b8f6b 100644 --- a/proclam/metrics/brier.py +++ b/proclam/metrics/brier.py @@ -27,7 +27,7 @@ def __init__(self, scheme=None): super(Brier, self).__init__(scheme) - def evaluate(self, prediction, truth, averaging='per_class'): + def evaluate(self, prediction, truth, weightvector, averaging='per_class'): """ Evaluates the Brier score @@ -37,6 +37,8 @@ def evaluate(self, prediction, truth, averaging='per_class'): predicted class probabilities truth: numpy.ndarray, int true classes + wieghtvec: numpy.ndarray, int + relative weights averaging: string, optional 'per_class' weights classes equally, other keywords possible vector assumed to be class weights @@ -60,7 +62,9 @@ def evaluate(self, prediction, truth, averaging='per_class'): q_each = (prediction - truth_mask) ** 2 class_brier = averager(q_each, truth, M) - metric = weight_sum(class_brier, weight_vector=weights) + weight_vector=weights*weightvector + + metric = weight_sum(class_brier, weight_vector) assert(~np.isnan(metric)) diff --git a/proclam/metrics/logloss.py b/proclam/metrics/logloss.py index a3fdc29..bed540b 100644 --- a/proclam/metrics/logloss.py +++ b/proclam/metrics/logloss.py @@ -30,7 +30,7 @@ def __init__(self, scheme=None): super(LogLoss, self).__init__(scheme) self.scheme = scheme - def evaluate(self, prediction, truth, averaging='per_class'): + def evaluate(self, prediction, truth, weightvector, averaging='per_class'): """ Evaluates the log-loss @@ -40,6 +40,8 @@ def evaluate(self, prediction, truth, averaging='per_class'): predicted class probabilities truth: numpy.ndarray, int true classes + weightvector: numpy.ndarray, float + per class weights averaging: string or numpy.ndarray, float 'per_class' weights classes equally, other keywords possible vector assumed to be class weights @@ -53,11 +55,14 @@ def evaluate(self, prediction, truth, averaging='per_class'): ----- This uses the natural log. """ + print(weightvector, 'checking') + prediction, truth = np.asarray(prediction), np.asarray(truth) prediction_shape = np.shape(prediction) (N, M) = prediction_shape weights = check_weights(averaging, M, truth=truth) + print('average weights', weights) truth_mask = truth_reformatter(truth, prediction) prediction = sanitize_predictions(prediction) @@ -67,8 +72,9 @@ def evaluate(self, prediction, truth, averaging='per_class'): # use a better structure for checking keyword support class_logloss = averager(logloss_each, truth, M) - - logloss = weight_sum(class_logloss, weight_vector=weights) + weight_vector = weights*weightvector + print('ok ready to go', weight_vector) + logloss = weight_sum(class_logloss, weight_vector=weight_vector) #=weights) assert(~np.isnan(logloss)) diff --git a/proclam/metrics/util.py b/proclam/metrics/util.py index 98c4671..9570978 100644 --- a/proclam/metrics/util.py +++ b/proclam/metrics/util.py @@ -349,19 +349,33 @@ def det_to_cm(dets, truth, per_class_norm=False, vb=False): ----- I need to fix the norm keyword all around to enable more options, like normed output vs. not. """ + pred_classes, pred_counts = np.unique(dets, return_counts=True) true_classes, true_counts = np.unique(truth, return_counts=True) - # if vb: print('by request '+str(((pred_classes, pred_counts), (true_classes, true_counts)))) + if vb: print('by request '+str((pred_classes, pred_counts), (true_classes, true_counts))) + +# print(pred_classes, true_classes, 'huh') M = np.int(max(max(pred_classes), max(true_classes)) + 1) # if vb: print('by request '+str((np.shape(dets), np.shape(truth)), M)) cm = np.zeros((M, M), dtype=int) coords = np.array(list(zip(dets, truth))) +# print(coords, 'huzzah') + #print(np.shape(coords), 'shape coords') indices, index_counts = np.unique(coords, axis=0, return_counts=True) if vb: print(indices.T, index_counts) index_counts = index_counts.astype(int) + + if per_class_norm: + # print(type(cm)) + # print(type(true_counts)) + print(np.shape(cm), np.shape(true_counts), 'shapes') + # cm = cm / true_counts + # cm /= true_counts[:, np.newaxis] # + cm = cm / true_counts[np.newaxis, :] + indices = indices.T.astype(int) cm[indices[0], indices[1]] = index_counts @@ -372,31 +386,238 @@ def det_to_cm(dets, truth, per_class_norm=False, vb=False): return cm -# def prob_to_cm(probs, truth, per_class_norm=True, vb=False): -# """ -# Turns probabilistic classifications into confusion matrix by taking maximum probability as deterministic class -# -# Parameters -# ---------- -# probs: numpy.ndarray, float -# N * M matrix of class probabilities -# truth: numpy.ndarray, int -# N-dimensional vector of true classes -# per_class_norm: boolean, optional -# equal weight per class if True, equal weight per object if False -# vb: boolean, optional -# if True, print cm -# -# Returns -# ------- -# cm: numpy.ndarray, int -# confusion matrix -# """ -# dets = prob_to_det(probs) -# -# cm = det_to_cm(dets, truth, per_class_norm=per_class_norm, vb=vb) -# -# return cm +def prob_to_cm(probs, truth, per_class_norm=True, vb=False): + """ + Turns probabilistic classifications into confusion matrix by taking maximum probability as deterministic class + + Parameters + ---------- + probs: numpy.ndarray, float + N * M matrix of class probabilities + truth: numpy.ndarray, int + N-dimensional vector of true classes + per_class_norm: boolean, optional + equal weight per class if True, equal weight per object if False + vb: boolean, optional + if True, print cm + + Returns + ------- + cm: numpy.ndarray, int + confusion matrix + """ + dets = prob_to_det(probs) + +# print(truth, 'huh 1') + cm = det_to_cm(dets, truth, per_class_norm=per_class_norm, vb=vb) + + return cm + +def cm_to_rate(cm, vb=False): + """ + Turns a confusion matrix into true/false positive/negative rates + + Parameters + ---------- + cm: numpy.ndarray, int or float + confusion matrix, first axis is predictions, second axis is truth + vb: boolean, optional + print progress to stdout? + + Returns + ------- + rates: named tuple, float + RateMatrix named tuple + + Notes + ----- + BROKEN! + This can be done with a mask to weight the classes differently here. + """ + if vb: print('by request '+str(cm)) + diag = np.diag(cm) + if vb: print('by request '+str(diag)) + + TP = np.sum(diag) + FN = np.sum(np.sum(cm, axis=0) - diag) + FP = np.sum(np.sum(cm, axis=1) - diag) + TN = np.sum(cm) - TP + if vb: print('by request '+str((TP, FN, FP, TN))) + + T = TP + TN + F = FP + FN + P = TP + FP + N = TN + FN + if vb: print('by request '+str((T, F, P, N))) + + TPR = TP / P + FPR = FP / N + FNR = FN / P + TNR = TN / N + + rates = RateMatrix(TPR=TPR, FPR=FPR, FNR=FNR, TNR=TNR) + if vb: print('by request '+str(rates)) + + return rates + +def det_to_rate(dets, truth, per_class_norm=True, vb=False): + cm = det_to_cm(dets, truth, per_class_norm=per_class_norm, vb=vb) + rates = cm_to_rate(cm, vb=vb) + return rates + +def prob_to_rate(probs, truth, per_class_norm=True, vb=False): + cm = prob_to_cm(probs, truth, per_class_norm=per_class_norm, vb=vb) + rates = cm_to_rate(cm, vb=vb) + return rates + +def weight_sum(per_class_metrics, weight_vector, norm=True): + """ + Calculates the weighted metric + + Parameters + ---------- + per_class_metrics: numpy.float + the scores separated by class (a list of arrays) + weight_vector: numpy.ndarray float + The array of weights per class + norm: boolean, optional + + Returns + ------- + weight_sum: np.float + The weighted metric + """ + weight_sum = np.dot(weight_vector, per_class_metrics) + + #print(weight_sum, 'weight_sum') + return weight_sum + +def check_weights(avg_info, M, chosen=None, truth=None): + """ + Converts standard weighting schemes to weight vectors for weight_sum + + Parameters + ---------- + avg_info: str or numpy.ndarray, float + keyword about how to calculate weighted average metric + M: int + number of classes + chosen: int, optional + which class is to be singled out for down/up-weighting + truth: numpy.ndarray, int, optional + true class assignments + + Returns + ------- + weights: numpy.ndarray, float + relative weights per class + + Notes + ----- + Assumes a random class + """ + if type(avg_info) != str: + avg_info = np.asarray(avg_info) + weights = avg_info / np.sum(avg_info) + assert(np.isclose(sum(weights), 1.)) + elif avg_info == 'per_class': + weights = np.ones(M) / float(M) + elif avg_info == 'per_item': + classes, counts = np.unique(truth, return_counts=True) + weights = np.zeros(M) + weights[classes] = counts / float(len(truth)) + assert len(weights) == M + elif avg_info == 'flat': + weights = np.ones(M) + elif avg_info == 'up' or avg_info == 'down': + if chosen is None: + chosen = np.random.randint(M) + if avg_info == 'up': + weights = np.ones(M) / np.float(M) + weights[chosen] = 1. + elif avg_info == 'down': + weights = np.ones(M) + weights[chosen] = 1./np.float(M) + else: + print('something has gone wrong with avg_info '+str(avg_info)) + return weights + +def averager(per_object_metrics, truth, M, vb=False): + """ + Creates a list with the metrics per object, separated by class + + Notes + ----- + There is currently a kludge for when there are no true class members, causing an improvement when that class is upweighted due to increasing the weight of 0. + """ + group_metric = per_object_metrics + class_metric = np.empty(M) + for m in range(M): + true_indices = np.where(truth == m)[0] + how_many_in_class = len(true_indices) + try: + assert(how_many_in_class > 0) + per_class_metric = group_metric[true_indices] + # assert(~np.all(np.isnan(per_class_metric))) + class_metric[m] = np.average(per_class_metric) + except AssertionError: + class_metric[m] = 0. + if vb: print('by request '+str((m, how_many_in_class, class_metric[m]))) + return class_metric + +def binary_rates(dets, truth, m): + + tp = np.sum(dets[truth == m]) + fp = np.sum(dets[truth != m]) + tpr = tp/len(dets[truth == m]) + fpr = fp/len(dets[truth != m]) + + return tpr,fpr + +def precision(classifications,truth,class_idx): + + tp = np.sum(classifications[truth == class_idx]) + fp = np.sum(classifications[truth != class_idx]) + + precision = tp/(tp+fp) + if precision != precision: + import pdb; pdb.set_trace() + + return tp/(tp+fp) + +def recall(classifications,truth,class_idx): + + tp = np.sum(classifications[truth == class_idx]) + fp = np.sum(classifications[truth != class_idx]) + fn = len(np.where((classifications == 0) & (truth == class_idx))[0]) + #import pdb; pdb.set_trace() + #print(fn) + + return tp/(tp+fn) + +def auc(x, y): + """ + Computes the area under curve (just a wrapper for trapezoid rule) + + Parameters + ---------- + x: numpy.ndarray, int or float + + y: numpy.ndarray, int or float + + Returns + ------- + rates: named tuple, float + RateMatrix named tuple + """ + + x = np.concatenate(([0.], x, [1.]),) + y = np.concatenate(([0.], y, [1.]),) + + i = np.argsort(x) + auc = trapz(y[i], x[i]) + + return auc # def det_to_rate(dets, truth, per_class_norm=True, vb=False): # cm = det_to_cm(dets, truth, per_class_norm=per_class_norm, vb=vb) @@ -466,3 +687,4 @@ def precision(TP, FP): # #print(fn) # # return tp/(tp+fn) + diff --git a/proclam/metrics/util.py~ b/proclam/metrics/util.py~ new file mode 100644 index 0000000..d595388 --- /dev/null +++ b/proclam/metrics/util.py~ @@ -0,0 +1,781 @@ +""" +Utility functions for PLAsTiCC metrics +""" + +from __future__ import absolute_import, division +__all__ = ['sanitize_predictions', + 'weight_sum', 'check_weights', 'averager', + 'cm_to_rate', 'precision', + 'auc', 'check_auc_grid', 'prep_curve', + 'det_to_prob', 'prob_to_det', + 'det_to_cm'] + +import collections +import numpy as np +# import pycm +import sys +from scipy.integrate import trapz + +RateMatrix = collections.namedtuple('rates', 'TPR FPR FNR TNR TP FP FN TN') + +def sanitize_predictions(predictions, epsilon=1.e-8): + """ + Replaces 0 and 1 with 0+epsilon, 1-epsilon + + Parameters + ---------- + predictions: numpy.ndarray, float + N*M matrix of probabilities per object, may have 0 or 1 values + epsilon: float + small placeholder number, defaults to floating point precision + + Returns + ------- + predictions: numpy.ndarray, float + N*M matrix of probabilities per object, no 0 or 1 values + """ + assert epsilon > 0. and epsilon < 0.0005 + mask1 = (predictions < epsilon) + mask2 = (predictions > 1.0 - epsilon) + + predictions[mask1] = epsilon + predictions[mask2] = 1.0 - epsilon + predictions = predictions / np.sum(predictions, axis=1)[:, np.newaxis] + return predictions + +def weight_sum(per_class_metrics, weight_vector): + """ + Calculates the weighted metric + + Parameters + ---------- + per_class_metrics: numpy.ndarray, float + vector of per-class scores + weight_vector: numpy.ndarray, float + vector of per-class weights + + Returns + ------- + weight_sum: np.float + The weighted metric + """ + weight_sum = np.dot(weight_vector, per_class_metrics) + return weight_sum + +def check_weights(avg_info, M, chosen=None, truth=None): + """ + Converts standard weighting schemes to weight vectors for weight_sum + + Parameters + ---------- + avg_info: str or numpy.ndarray, float + keyword about how to calculate weighted average metric + M: int + number of classes + chosen: int, optional + which class is to be singled out for down/up-weighting + truth: numpy.ndarray, int, optional + true class assignments + + Returns + ------- + weights: numpy.ndarray, float + relative weights per class + + Notes + ----- + Assumes a random class + """ + if type(avg_info) != str: + avg_info = np.asarray(avg_info) + weights = avg_info / np.sum(avg_info) + assert(np.isclose(sum(weights), 1.)) + elif avg_info == 'per_class': + weights = np.ones(M) / float(M) + elif avg_info == 'per_item': + classes, counts = np.unique(truth, return_counts=True) + weights = np.zeros(M) + weights[classes] = counts / float(len(truth)) + assert len(weights) == M + elif avg_info == 'flat': + weights = np.ones(M) + elif avg_info == 'up' or avg_info == 'down': + if chosen is None: + chosen = np.random.randint(M) + if avg_info == 'up': + weights = np.ones(M) / np.float(M) + weights[chosen] = 1. + elif avg_info == 'down': + weights = np.ones(M) + weights[chosen] = 1./np.float(M) + else: + print('something has gone wrong with avg_info '+str(avg_info)) + weights = None + return weights + +def averager(per_object_metrics, truth, M, vb=False): + """ + Creates a list with the metrics per object, separated by class + + Notes + ----- + There is currently a kludge for when there are no true class members, causing an improvement when that class is upweighted due to increasing the weight of 0. + """ + group_metric = per_object_metrics + class_metric = np.empty(M) + for m in range(M): + true_indices = np.where(truth == m)[0] + how_many_in_class = len(true_indices) + try: + assert(how_many_in_class > 0) + per_class_metric = group_metric[true_indices] + # assert(~np.all(np.isnan(per_class_metric))) + class_metric[m] = np.average(per_class_metric) + except AssertionError: + class_metric[m] = 0. + if vb: print('by request '+str((m, how_many_in_class, class_metric[m]))) + return class_metric + +def cm_to_rate(cm, vb=False): + """ + Turns a confusion matrix into true/false positive/negative rates + + Parameters + ---------- + cm: numpy.ndarray, int or float + confusion matrix, first axis is predictions, second axis is truth + vb: boolean, optional + print progress to stdout? + + Returns + ------- + rates: named tuple, float + RateMatrix named tuple + + Notes + ----- + This can be done with a mask to weight the classes differently here. + """ + cm = cm.astype(float) + # if vb: print('by request cm '+str(cm)) + tot = np.sum(cm) + # mask = range(len(cm)) + # if vb: print('by request sum '+str(tot)) + + T = np.sum(cm, axis=1) + F = tot[np.newaxis] - T + P = np.sum(cm, axis=0) + N = tot[np.newaxis] - P + # if vb: print('by request T, F, P, N'+str((T, F, P, N))) + + TP = np.diag(cm) + FN = P - TP + TN = F - FN + FP = T - TP + # if vb: print('by request TP, FP, FN, TN'+str((TP, FP, FN, TN))) + + TPR = TP / P + FPR = FP / N + FNR = FN / P + TNR = TN / N + # if vb: print('by request TPR, FPR, FNR, TNR'+str((TPR, FPR, FNR, TNR))) + + rates = RateMatrix(TPR=TPR, FPR=FPR, FNR=FNR, TNR=TNR, TP=TP, FN=FN, TN=TN, FP=FP) + # if vb: print('by request TPR, FPR, FNR, TNR '+str(rates)) + + return rates + +def prep_curve(x, y): + """ + Makes a curve for AUC + + Parameters + ---------- + x: numpy.ndarray, float + x-axis + y: numpy.ndarray, float + y-axis + + Returns + ------- + x: numpy.ndarray, float + x-axis + y: numpy.ndarray, float + y-axis + """ + x = np.concatenate(([0.], x, [1.]),) + y = np.concatenate(([0.], y, [1.]),) + return (x, y) + +def auc(x, y): + """ + Computes the area under curve (just a wrapper for trapezoid rule) + + Parameters + ---------- + x: numpy.ndarray, float + x-axis + y: numpy.ndarray, float + y-axis + + Returns + ------- + auc: float + the area under the curve + """ + i = np.argsort(x) + auc = trapz(y[i], x[i]) + return auc + +def check_auc_grid(grid): + """ + Checks if a grid for an AUC metric is valid + + Parameters + ---------- + grid: numpy.ndarray, float or float or int + array of values between 0 and 1 at which to evaluate AUC or grid spacing or number of grid points + + Returns + ------- + thresholds_grid: numpy.ndarray, float + grid of thresholds + """ + if type(grid) == list or type(grid) == np.ndarray: + thresholds_grid = np.concatenate((np.zeros(1), np.array(grid), np.ones(1))) + elif type(grid) == float: + if grid > 0. and grid < 1.: + thresholds_grid = np.arange(0., 1., grid) + else: + thresholds_grid = None + elif type(grid) == int: + if grid > 0: + thresholds_grid = np.linspace(0., 1., grid) + else: + thresholds_grid = None + try: + assert thresholds_grid is not None + return np.sort(thresholds_grid) + except AssertionError: + print('Please specify a grid, spacing, or density for this AUC metric.') + return + +def det_to_prob(dets, prediction=None): + """ + Reformats vector of class assignments into matrix with 1 at true/assigned class and zero elsewhere + + Parameters + ---------- + dets: numpy.ndarray, int + vector of classes + prediction: numpy.ndarray, float, optional + predicted class probabilities + + Returns + ------- + probs: numpy.ndarray, float + matrix with 1 at input classes and 0 elsewhere + + Notes + ----- + formerly truth_reformatter + Does not yet handle number of classes in truth not matching number of classes in prediction, i.e. for having "other" class or secret classes not in training set. The prediction keyword is a kludge to enable this but should be replaced. + """ + N = len(dets) + indices = range(N) + + if prediction is None: + prediction_shape = (N, int(np.max(dets) + 1)) + else: + prediction, dets = np.asarray(prediction), np.asarray(dets) + prediction_shape = np.shape(prediction) + + probs = np.zeros(prediction_shape) + probs[indices, dets] = 1. + + return probs + +def prob_to_det(probs, m=None, threshold=None): + """ + Converts probabilistic classifications to deterministic classifications by assigning the class with highest probability + + Parameters + ---------- + probs: numpy.ndarray, float + N * M matrix of class probabilities + m: int + class relative to binary decision + threshold: float, optional + value between 0 and 1 at which binary decision is made + + Returns + ------- + dets: numpy.ndarray, int + maximum probability classes + """ + if m == None and threshold == None: + dets = np.argmax(probs, axis=1) + else: + try: + assert(type(m) == int and type(threshold) == np.float64) + except AssertionError: + print(str(m)+' is '+str(type(m))+' and must be int; '+str(threshold)+' is '+str(type(threshold))+' and must be float') + dets = np.zeros(np.shape(probs)[0]).astype(int) + dets[probs[:, m] >= threshold] = 1 + + return dets + +def det_to_cm(dets, truth, per_class_norm=False, vb=False): + """ + Converts deterministic classifications and truth into confusion matrix + + Parameters + ---------- + dets: numpy.ndarray, int + assigned classes + truth: numpy.ndarray, int + true classes + per_class_norm: boolean, optional + equal weight per class if True, equal weight per object if False + vb: boolean, optional + if True, print cm + + Returns + ------- + cm: numpy.ndarray, int + confusion matrix + + Notes + ----- + I need to fix the norm keyword all around to enable more options, like normed output vs. not. + """ + + pred_classes, pred_counts = np.unique(dets, return_counts=True) + true_classes, true_counts = np.unique(truth, return_counts=True) +<<<<<<< HEAD + + if vb: print('by request '+str((pred_classes, pred_counts), (true_classes, true_counts))) +======= + # if vb: print('by request '+str(((pred_classes, pred_counts), (true_classes, true_counts)))) +>>>>>>> 86c514163e7143486f8927aaf93f16f4346579b3 + +# print(pred_classes, true_classes, 'huh') + M = np.int(max(max(pred_classes), max(true_classes)) + 1) + + # if vb: print('by request '+str((np.shape(dets), np.shape(truth)), M)) + cm = np.zeros((M, M), dtype=int) + + coords = np.array(list(zip(dets, truth))) +# print(coords, 'huzzah') + #print(np.shape(coords), 'shape coords') + indices, index_counts = np.unique(coords, axis=0, return_counts=True) + if vb: print(indices.T, index_counts) + index_counts = index_counts.astype(int) +<<<<<<< HEAD + #print(index_counts, 'index_counts') + #print(indices, 'indices') + if vb: print('by request '+str(index_counts)) + # if vb: print(indices, index_counts) + indices = indices.T + # if vb: print(np.shape(indices)) + #print(cm, 'yo') + cm[indices[0], indices[1]] = index_counts + #if vb: print(cm) + #print(cm, 'hi') + + if per_class_norm: + # print(type(cm)) + # print(type(true_counts)) + print(np.shape(cm), np.shape(true_counts), 'shapes') + # cm = cm / true_counts + # cm /= true_counts[:, np.newaxis] # + cm = cm / true_counts[np.newaxis, :] +======= + indices = indices.T.astype(int) + cm[indices[0], indices[1]] = index_counts + + if per_class_norm: + cm = cm.astype(float) / true_counts[np.newaxis, :].astype(float) +>>>>>>> 86c514163e7143486f8927aaf93f16f4346579b3 + + if vb: print('by request '+str(cm)) + + return cm + +<<<<<<< HEAD +def prob_to_cm(probs, truth, per_class_norm=True, vb=False): + """ + Turns probabilistic classifications into confusion matrix by taking maximum probability as deterministic class + + Parameters + ---------- + probs: numpy.ndarray, float + N * M matrix of class probabilities + truth: numpy.ndarray, int + N-dimensional vector of true classes + per_class_norm: boolean, optional + equal weight per class if True, equal weight per object if False + vb: boolean, optional + if True, print cm + + Returns + ------- + cm: numpy.ndarray, int + confusion matrix + """ + dets = prob_to_det(probs) + +# print(truth, 'huh 1') + cm = det_to_cm(dets, truth, per_class_norm=per_class_norm, vb=vb) + + return cm + +def cm_to_rate(cm, vb=False): + """ + Turns a confusion matrix into true/false positive/negative rates + + Parameters + ---------- + cm: numpy.ndarray, int or float + confusion matrix, first axis is predictions, second axis is truth + vb: boolean, optional + print progress to stdout? + + Returns + ------- + rates: named tuple, float + RateMatrix named tuple + + Notes + ----- + BROKEN! + This can be done with a mask to weight the classes differently here. + """ + if vb: print('by request '+str(cm)) + diag = np.diag(cm) + if vb: print('by request '+str(diag)) + + TP = np.sum(diag) + FN = np.sum(np.sum(cm, axis=0) - diag) + FP = np.sum(np.sum(cm, axis=1) - diag) + TN = np.sum(cm) - TP + if vb: print('by request '+str((TP, FN, FP, TN))) + + T = TP + TN + F = FP + FN + P = TP + FP + N = TN + FN + if vb: print('by request '+str((T, F, P, N))) + + TPR = TP / P + FPR = FP / N + FNR = FN / P + TNR = TN / N + + rates = RateMatrix(TPR=TPR, FPR=FPR, FNR=FNR, TNR=TNR) + if vb: print('by request '+str(rates)) + + return rates +======= +# def prob_to_cm(probs, truth, per_class_norm=True, vb=False): +# """ +# Turns probabilistic classifications into confusion matrix by taking maximum probability as deterministic class +# +# Parameters +# ---------- +# probs: numpy.ndarray, float +# N * M matrix of class probabilities +# truth: numpy.ndarray, int +# N-dimensional vector of true classes +# per_class_norm: boolean, optional +# equal weight per class if True, equal weight per object if False +# vb: boolean, optional +# if True, print cm +# +# Returns +# ------- +# cm: numpy.ndarray, int +# confusion matrix +# """ +# dets = prob_to_det(probs) +# +# cm = det_to_cm(dets, truth, per_class_norm=per_class_norm, vb=vb) +# +# return cm + +#def cm_to_rate(cm, vb=False): +# """ +# Turns a confusion matrix into true/false positive/negative rates + +# Parameters +# ---------- +# cm: numpy.ndarray, int or float +# confusion matrix, first axis is predictions, second axis is truth +# vb: boolean, optional +# print progress to stdout? + +# Returns +# ------- +# rates: named tuple, float +# RateMatrix named tuple + +# Notes +# ----- + # BROKEN! + # This can be done with a mask to weight the classes differently here. + # """ + # if vb: print('by request '+str(cm)) + # diag = np.diag(cm) + # if vb: print('by request '+str(diag))# +# +# TP = np.sum(diag) +# FN = np.sum(np.sum(cm, axis=0) - diag) +# FP = np.sum(np.sum(cm, axis=1) - diag) +# TN = np.sum(cm) - TP +# if vb: print('by request '+str((TP, FN, FP, TN))) + +# T = TP + TN +# F = FP + FN +# P = TP + FP +# N = TN + FN +# if vb: print('by request '+str((T, F, P, N))) + +# TPR = TP / P +# FPR = FP / N +# FNR = FN / P +# TNR = TN / N + +# rates = RateMatrix(TPR=TPR, FPR=FPR, FNR=FNR, TNR=TNR) +# if vb: print('by request '+str(rates)) + +# return rates +>>>>>>> 86c514163e7143486f8927aaf93f16f4346579b3 + +def det_to_rate(dets, truth, per_class_norm=True, vb=False): + cm = det_to_cm(dets, truth, per_class_norm=per_class_norm, vb=vb) + rates = cm_to_rate(cm, vb=vb) + return rates + +def prob_to_rate(probs, truth, per_class_norm=True, vb=False): + cm = prob_to_cm(probs, truth, per_class_norm=per_class_norm, vb=vb) + rates = cm_to_rate(cm, vb=vb) + return rates + +def weight_sum(per_class_metrics, weight_vector, norm=True): + """ + Calculates the weighted metric + + Parameters + ---------- + per_class_metrics: numpy.float + the scores separated by class (a list of arrays) + weight_vector: numpy.ndarray float + The array of weights per class + norm: boolean, optional + + Returns + ------- + weight_sum: np.float + The weighted metric + """ + weight_sum = np.dot(weight_vector, per_class_metrics) + + #print(weight_sum, 'weight_sum') + return weight_sum + +def check_weights(avg_info, M, chosen=None, truth=None): + """ + Converts standard weighting schemes to weight vectors for weight_sum + + Parameters + ---------- + avg_info: str or numpy.ndarray, float + keyword about how to calculate weighted average metric + M: int + number of classes + chosen: int, optional + which class is to be singled out for down/up-weighting + truth: numpy.ndarray, int, optional + true class assignments + + Returns + ------- + weights: numpy.ndarray, float + relative weights per class + + Notes + ----- + Assumes a random class + """ + if type(avg_info) != str: + avg_info = np.asarray(avg_info) + weights = avg_info / np.sum(avg_info) + assert(np.isclose(sum(weights), 1.)) + elif avg_info == 'per_class': + weights = np.ones(M) / float(M) + elif avg_info == 'per_item': + classes, counts = np.unique(truth, return_counts=True) + weights = np.zeros(M) + weights[classes] = counts / float(len(truth)) + assert len(weights) == M + elif avg_info == 'flat': + weights = np.ones(M) + elif avg_info == 'up' or avg_info == 'down': + if chosen is None: + chosen = np.random.randint(M) + if avg_info == 'up': + weights = np.ones(M) / np.float(M) + weights[chosen] = 1. + elif avg_info == 'down': + weights = np.ones(M) + weights[chosen] = 1./np.float(M) + else: + print('something has gone wrong with avg_info '+str(avg_info)) + return weights + +def averager(per_object_metrics, truth, M, vb=False): + """ + Creates a list with the metrics per object, separated by class + + Notes + ----- + There is currently a kludge for when there are no true class members, causing an improvement when that class is upweighted due to increasing the weight of 0. + """ + group_metric = per_object_metrics + class_metric = np.empty(M) + for m in range(M): + true_indices = np.where(truth == m)[0] + how_many_in_class = len(true_indices) + try: + assert(how_many_in_class > 0) + per_class_metric = group_metric[true_indices] + # assert(~np.all(np.isnan(per_class_metric))) + class_metric[m] = np.average(per_class_metric) + except AssertionError: + class_metric[m] = 0. + if vb: print('by request '+str((m, how_many_in_class, class_metric[m]))) + return class_metric + +def binary_rates(dets, truth, m): + + tp = np.sum(dets[truth == m]) + fp = np.sum(dets[truth != m]) + tpr = tp/len(dets[truth == m]) + fpr = fp/len(dets[truth != m]) + + return tpr,fpr + +def precision(classifications,truth,class_idx): + + tp = np.sum(classifications[truth == class_idx]) + fp = np.sum(classifications[truth != class_idx]) + + precision = tp/(tp+fp) + if precision != precision: + import pdb; pdb.set_trace() + + return tp/(tp+fp) + +def recall(classifications,truth,class_idx): + + tp = np.sum(classifications[truth == class_idx]) + fp = np.sum(classifications[truth != class_idx]) + fn = len(np.where((classifications == 0) & (truth == class_idx))[0]) + #import pdb; pdb.set_trace() + #print(fn) + + return tp/(tp+fn) + +def auc(x, y): + """ + Computes the area under curve (just a wrapper for trapezoid rule) + + Parameters + ---------- + x: numpy.ndarray, int or float + + y: numpy.ndarray, int or float + + Returns + ------- + rates: named tuple, float + RateMatrix named tuple + """ + + x = np.concatenate(([0.], x, [1.]),) + y = np.concatenate(([0.], y, [1.]),) + + i = np.argsort(x) + auc = trapz(y[i], x[i]) + + return auc + +# def det_to_rate(dets, truth, per_class_norm=True, vb=False): +# cm = det_to_cm(dets, truth, per_class_norm=per_class_norm, vb=vb) +# rates = cm_to_rate(cm, vb=vb) +# return rates +# +# def prob_to_rate(probs, truth, per_class_norm=True, vb=False): +# cm = prob_to_cm(probs, truth, per_class_norm=per_class_norm, vb=vb) +# rates = cm_to_rate(cm, vb=vb) +# return rates + +# def binary_rates(dets, truth, m): +# +# tp = np.sum(dets[truth == m]) +# fp = np.sum(dets[truth != m]) +# tpr = tp/len(dets[truth == m]) +# fpr = fp/len(dets[truth != m]) +# +# return tpr,fpr + +# def recall(rates): +# return 1. - rates.FNR + +# def recall(rates): +# """ +# Calculates recall from rates +# +# Parameters +# ---------- +# rates: namedtuple +# named tuple of 'TPR FPR FNR TNR' +# +# Returns +# ------- +# recall: float +# recall +# """ +# return 1. - rates.FNR + +def precision(TP, FP): + """ + Calculates precision from rates + + Parameters + ---------- + TP: float + number of true positives + FP: float + number of false positives + + Returns + ------- + p: float + precision + """ + p = np.asarray(TP / (TP + FP)) + if np.any(np.isnan(p)): + p[np.isnan(p)] = 0. + return p +# +# def recall(classifications,truth,class_idx): +# +# tp = np.sum(classifications[truth == class_idx]) +# fp = np.sum(classifications[truth != class_idx]) +# fn = len(np.where((classifications == 0) & (truth == class_idx))[0]) +# #import pdb; pdb.set_trace() +# #print(fn) +# +# return tp/(tp+fn) +