1
2
3
4
5
6
7
8
9 """Utility class to compute the transfer error of classifiers."""
10
11 __docformat__ = 'restructuredtext'
12
13 import mvpa.misc.copy as copy
14
15 import numpy as N
16
17 from sets import Set
18 from StringIO import StringIO
19 from math import log10, ceil
20
21 from mvpa.base import externals
22
23 from mvpa.misc.errorfx import meanPowerFx, rootMeanPowerFx, RMSErrorFx, \
24 CorrErrorFx, CorrErrorPFx, RelativeRMSErrorFx, MeanMismatchErrorFx, \
25 AUCErrorFx
26 from mvpa.base import warning
27 from mvpa.misc.state import StateVariable, Stateful
28 from mvpa.base.dochelpers import enhancedDocString, table2string
29 from mvpa.clfs.stats import autoNullDist
30
31 if __debug__:
32 from mvpa.base import debug
33
34 if externals.exists('scipy'):
35 from scipy.stats.stats import nanmean
36 else:
37 from mvpa.clfs.stats import nanmean
38
39 -def _p2(x, prec=2):
40 """Helper to print depending on the type nicely. For some
41 reason %.2g for 100 prints exponential form which is ugly
42 """
43 if isinstance(x, int):
44 return "%d" % x
45 elif isinstance(x, float):
46 s = ("%%.%df" % prec % x).rstrip('0').rstrip('.').lstrip()
47 if s == '':
48 s = '0'
49 return s
50 else:
51 return "%s" % x
52
56 """Basic class to collect targets/predictions and report summary statistics
57
58 It takes care about collecting the sets, which are just tuples
59 (targets, predictions, values). While 'computing' the matrix, all
60 sets are considered together. Children of the class are
61 responsible for computation and display.
62 """
63
64 _STATS_DESCRIPTION = (
65 ('# of sets',
66 'number of target/prediction sets which were provided',
67 None), )
68
69
70 - def __init__(self, targets=None, predictions=None, values=None, sets=None):
71 """Initialize SummaryStatistics
72
73 targets or predictions cannot be provided alone (ie targets
74 without predictions)
75
76 :Parameters:
77 targets
78 Optional set of targets
79 predictions
80 Optional set of predictions
81 values
82 Optional set of values (which served for prediction)
83 sets
84 Optional list of sets
85 """
86 self._computed = False
87 """Flag either it was computed for a given set of data"""
88
89 self.__sets = (sets, [])[int(sets is None)]
90 """Datasets (target, prediction) to compute confusion matrix on"""
91
92 if not targets is None or not predictions is None:
93 if not targets is None and not predictions is None:
94 self.add(targets=targets, predictions=predictions,
95 values=values)
96 else:
97 raise ValueError, \
98 "Please provide none or both targets and predictions"
99
100
101 - def add(self, targets, predictions, values=None):
102 """Add new results to the set of known results"""
103 if len(targets) != len(predictions):
104 raise ValueError, \
105 "Targets[%d] and predictions[%d]" % (len(targets),
106 len(predictions)) + \
107 " have different number of samples"
108
109 if values is not None and len(targets) != len(values):
110 raise ValueError, \
111 "Targets[%d] and values[%d]" % (len(targets),
112 len(values)) + \
113 " have different number of samples"
114
115
116
117
118 nonetype = type(None)
119 for i in xrange(len(targets)):
120 t1, t2 = type(targets[i]), type(predictions[i])
121
122
123 if t1 != t2 and t2 != nonetype:
124
125
126 if isinstance(predictions, tuple):
127 predictions = list(predictions)
128 predictions[i] = t1(predictions[i])
129
130 self.__sets.append( (targets, predictions, values) )
131 self._computed = False
132
133
134 - def asstring(self, short=False, header=True, summary=True,
135 description=False):
136 """'Pretty print' the matrix
137
138 :Parameters:
139 short : bool
140 if True, ignores the rest of the parameters and provides consise
141 1 line summary
142 header : bool
143 print header of the table
144 summary : bool
145 print summary (accuracy)
146 description : bool
147 print verbose description of presented statistics
148 """
149 raise NotImplementedError
150
151
153 """String summary over the `SummaryStatistics`
154
155 It would print description of the summary statistics if 'CM'
156 debug target is active
157 """
158 if __debug__:
159 description = ('CM' in debug.active)
160 else:
161 description = False
162 return self.asstring(short=False, header=True, summary=True,
163 description=description)
164
165
167 """Add the sets from `other` s `SummaryStatistics` to current one
168 """
169
170
171
172 othersets = copy.copy(other.__sets)
173 for set in othersets:
174 self.add(*set)
175 return self
176
177
179 """Add two `SummaryStatistics`s
180 """
181 result = copy.copy(self)
182 result += other
183 return result
184
185
187 """Actually compute the confusion matrix based on all the sets"""
188 if self._computed:
189 return
190
191 self._compute()
192 self._computed = True
193
194
196 self._stats = {'# of sets' : len(self.sets)}
197
198
199 @property
201 """Return a list of separate summaries per each stored set"""
202 return [ self.__class__(sets=[x]) for x in self.sets ]
203
204
205 @property
207 raise NotImplementedError
208
209
210 @property
212 self.compute()
213 return self._stats
214
215
217 """Cleans summary -- all data/sets are wiped out
218 """
219 self.__sets = []
220 self._computed = False
221
222
223 sets = property(lambda self:self.__sets)
224
227 """Generic class for ROC curve computation and plotting
228 """
229
231 """
232 :Parameters:
233 labels : list
234 labels which were used (in order of values if multiclass,
235 or 1 per class for binary problems (e.g. in SMLR))
236 sets : list of tuples
237 list of sets for the analysis
238 """
239 self._labels = labels
240 self._sets = sets
241 self.__computed = False
242
243
245 """Lazy computation if needed
246 """
247 if self.__computed:
248 return
249
250 labels = self._labels
251 Nlabels = len(labels)
252 sets = self._sets
253
254
255 def _checkValues(set_):
256 """Check if values are 'acceptable'"""
257 x = set_[2]
258
259 if (x is None) or len(x) == 0: return False
260 for v in x:
261 try:
262 if Nlabels <= 2 and N.isscalar(v):
263 continue
264 if (isinstance(v, dict) or
265 ((Nlabels>=2) and len(v)!=Nlabels)
266 ): return False
267 except Exception, e:
268
269
270
271 if __debug__:
272 debug('ROC', "Exception %s while checking "
273 "either %s are valid labels" % (str(e), x))
274 return False
275 return True
276
277 sets_wv = filter(_checkValues, sets)
278
279 Nsets_wv = len(sets_wv)
280 if Nsets_wv > 0 and len(sets) != Nsets_wv:
281 warning("Only %d sets have values assigned from %d sets" %
282 (Nsets_wv, len(sets)))
283
284
285
286
287 for iset,s in enumerate(sets_wv):
288
289 values = s[2]
290
291 if isinstance(values, N.ndarray) and len(values.shape)==1:
292 values = list(values)
293 rangev = None
294 for i in xrange(len(values)):
295 v = values[i]
296 if N.isscalar(v):
297 if Nlabels == 2:
298 def last_el(x):
299 """Helper function. Returns x if x is scalar, and
300 last element if x is not (ie list/tuple)"""
301 if N.isscalar(x): return x
302 else: return x[-1]
303 if rangev is None:
304
305
306 values_ = [last_el(x) for x in values]
307 rangev = N.min(values_) + N.max(values_)
308 values[i] = [rangev - v, v]
309 else:
310 raise ValueError, \
311 "Cannot have a single 'value' for multiclass" \
312 " classification. Got %s" % (v)
313 elif len(v) != Nlabels:
314 raise ValueError, \
315 "Got %d values whenever there is %d labels" % \
316 (len(v), Nlabels)
317
318 sets_wv[iset] = (s[0], s[1], N.asarray(values))
319
320
321
322
323
324 ROCs, aucs = [], []
325 for i,label in enumerate(labels):
326 aucs_pl = []
327 ROCs_pl = []
328 for s in sets_wv:
329 targets_pl = (s[0] == label).astype(int)
330
331 ROC = AUCErrorFx()
332 aucs_pl += [ROC([x[i] for x in s[2]], targets_pl)]
333 ROCs_pl.append(ROC)
334 if len(aucs_pl)>0:
335 ROCs += [ROCs_pl]
336 aucs += [nanmean(aucs_pl)]
337
338
339
340 self._ROCs = ROCs
341 self._aucs = aucs
342 self.__computed = True
343
344
345 @property
347 """Compute and return set of AUC values 1 per label
348 """
349 self._compute()
350 return self._aucs
351
352
353 @property
355 self._compute()
356 return self._ROCs
357
358
359 - def plot(self, label_index=0):
360 """
361
362 TODO: make it friendly to labels given by values?
363 should we also treat labels_map?
364 """
365 externals.exists("pylab", raiseException=True)
366 import pylab as P
367
368 self._compute()
369
370 labels = self._labels
371
372 ROCs = self.ROCs[label_index]
373
374 fig = P.gcf()
375 ax = P.gca()
376
377 P.plot([0, 1], [0, 1], 'k:')
378
379 for ROC in ROCs:
380 P.plot(ROC.fp, ROC.tp, linewidth=1)
381
382 P.axis((0.0, 1.0, 0.0, 1.0))
383 P.axis('scaled')
384 P.title('Label %s. Mean AUC=%.2f' % (label_index, self.aucs[label_index]))
385
386 P.xlabel('False positive rate')
387 P.ylabel('True positive rate')
388
391 """Class to contain information and display confusion matrix.
392
393 Implementation of the `SummaryStatistics` in the case of
394 classification problem. Actual computation of confusion matrix is
395 delayed until all data is acquired (to figure out complete set of
396 labels). If testing data doesn't have a complete set of labels,
397 but you like to include all labels, provide them as a parameter to
398 the constructor.
399
400 Confusion matrix provides a set of performance statistics (use
401 asstring(description=True) for the description of abbreviations),
402 as well ROC curve (http://en.wikipedia.org/wiki/ROC_curve)
403 plotting and analysis (AUC) in the limited set of problems:
404 binary, multiclass 1-vs-all.
405 """
406
407 _STATS_DESCRIPTION = (
408 ('TP', 'true positive (AKA hit)', None),
409 ('TN', 'true negative (AKA correct rejection)', None),
410 ('FP', 'false positive (AKA false alarm, Type I error)', None),
411 ('FN', 'false negative (AKA miss, Type II error)', None),
412 ('TPR', 'true positive rate (AKA hit rate, recall, sensitivity)',
413 'TPR = TP / P = TP / (TP + FN)'),
414 ('FPR', 'false positive rate (AKA false alarm rate, fall-out)',
415 'FPR = FP / N = FP / (FP + TN)'),
416 ('ACC', 'accuracy', 'ACC = (TP + TN) / (P + N)'),
417 ('SPC', 'specificity', 'SPC = TN / (FP + TN) = 1 - FPR'),
418 ('PPV', 'positive predictive value (AKA precision)',
419 'PPV = TP / (TP + FP)'),
420 ('NPV', 'negative predictive value', 'NPV = TN / (TN + FN)'),
421 ('FDR', 'false discovery rate', 'FDR = FP / (FP + TP)'),
422 ('MCC', "Matthews Correlation Coefficient",
423 "MCC = (TP*TN - FP*FN)/sqrt(P N P' N')"),
424 ('AUC', "Area under (AUC) curve", None),
425 ) + SummaryStatistics._STATS_DESCRIPTION
426
427
428 - def __init__(self, labels=None, labels_map=None, **kwargs):
429 """Initialize ConfusionMatrix with optional list of `labels`
430
431 :Parameters:
432 labels : list
433 Optional set of labels to include in the matrix
434 labels_map : None or dict
435 Dictionary from original dataset to show mapping into
436 numerical labels
437 targets
438 Optional set of targets
439 predictions
440 Optional set of predictions
441 """
442
443 SummaryStatistics.__init__(self, **kwargs)
444
445 if labels == None:
446 labels = []
447 self.__labels = labels
448 """List of known labels"""
449 self.__labels_map = labels_map
450 """Mapping from original into given labels"""
451 self.__matrix = None
452 """Resultant confusion matrix"""
453
454
455
456
457 @property
463
464
466 """Actually compute the confusion matrix based on all the sets"""
467
468 super(ConfusionMatrix, self)._compute()
469
470 if __debug__:
471 if not self.__matrix is None:
472 debug("LAZY",
473 "Have to recompute %s#%s" \
474 % (self.__class__.__name__, id(self)))
475
476
477
478
479 try:
480
481 labels = \
482 list(reduce(lambda x, y: x.union(Set(y[0]).union(Set(y[1]))),
483 self.sets,
484 Set(self.__labels)))
485 except:
486 labels = self.__labels
487
488
489 labels_map = self.__labels_map
490 if labels_map is not None:
491 labels_set = Set(labels)
492 map_labels_set = Set(labels_map.values())
493
494 if not map_labels_set.issuperset(labels_set):
495 warning("Provided labels_map %s is not coherent with labels "
496 "provided to ConfusionMatrix. No reverse mapping "
497 "will be provided" % labels_map)
498 labels_map = None
499
500
501 labels_map_rev = None
502 if labels_map is not None:
503 labels_map_rev = {}
504 for k,v in labels_map.iteritems():
505 v_mapping = labels_map_rev.get(v, [])
506 v_mapping.append(k)
507 labels_map_rev[v] = v_mapping
508 self.__labels_map_rev = labels_map_rev
509
510 labels.sort()
511 self.__labels = labels
512
513 Nlabels, Nsets = len(labels), len(self.sets)
514
515 if __debug__:
516 debug("CM", "Got labels %s" % labels)
517
518
519 mat_all = N.zeros( (Nsets, Nlabels, Nlabels), dtype=int )
520
521
522
523
524 counts_all = N.zeros( (Nsets, Nlabels) )
525
526
527 rev_map = dict([ (x[1], x[0]) for x in enumerate(labels)])
528 for iset, (targets, predictions, values) in enumerate(self.sets):
529 for t,p in zip(targets, predictions):
530 mat_all[iset, rev_map[p], rev_map[t]] += 1
531
532
533
534
535
536 self.__matrix = N.sum(mat_all, axis=0)
537 self.__Nsamples = N.sum(self.__matrix, axis=0)
538 self.__Ncorrect = sum(N.diag(self.__matrix))
539
540 TP = N.diag(self.__matrix)
541 offdiag = self.__matrix - N.diag(TP)
542 stats = {
543 '# of labels' : Nlabels,
544 'TP' : TP,
545 'FP' : N.sum(offdiag, axis=1),
546 'FN' : N.sum(offdiag, axis=0)}
547
548 stats['CORR'] = N.sum(TP)
549 stats['TN'] = stats['CORR'] - stats['TP']
550 stats['P'] = stats['TP'] + stats['FN']
551 stats['N'] = N.sum(stats['P']) - stats['P']
552 stats["P'"] = stats['TP'] + stats['FP']
553 stats["N'"] = stats['TN'] + stats['FN']
554 stats['TPR'] = stats['TP'] / (1.0*stats['P'])
555
556
557 stats['TPR'][stats['P'] == 0] = 0
558 stats['PPV'] = stats['TP'] / (1.0*stats["P'"])
559 stats['NPV'] = stats['TN'] / (1.0*stats["N'"])
560 stats['FDR'] = stats['FP'] / (1.0*stats["P'"])
561 stats['SPC'] = (stats['TN']) / (1.0*stats['FP'] + stats['TN'])
562
563 MCC_denom = N.sqrt(1.0*stats['P']*stats['N']*stats["P'"]*stats["N'"])
564 nz = MCC_denom!=0.0
565 stats['MCC'] = N.zeros(stats['TP'].shape)
566 stats['MCC'][nz] = \
567 (stats['TP'] * stats['TN'] - stats['FP'] * stats['FN'])[nz] \
568 / MCC_denom[nz]
569
570 stats['ACC'] = N.sum(TP)/(1.0*N.sum(stats['P']))
571 stats['ACC%'] = stats['ACC'] * 100.0
572
573
574
575 ROC = ROCCurve(labels=labels, sets=self.sets)
576 aucs = ROC.aucs
577 if len(aucs)>0:
578 stats['AUC'] = aucs
579 if len(aucs) != Nlabels:
580 raise RuntimeError, \
581 "We must got a AUC per label. Got %d instead of %d" % \
582 (len(aucs), Nlabels)
583 self.ROC = ROC
584 else:
585
586 stats['AUC'] = [N.nan] * Nlabels
587 self.ROC = None
588
589
590
591 for k,v in stats.items():
592 stats['mean(%s)' % k] = N.mean(v)
593
594 self._stats.update(stats)
595
596
597 - def asstring(self, short=False, header=True, summary=True,
598 description=False):
599 """'Pretty print' the matrix
600
601 :Parameters:
602 short : bool
603 if True, ignores the rest of the parameters and provides consise
604 1 line summary
605 header : bool
606 print header of the table
607 summary : bool
608 print summary (accuracy)
609 description : bool
610 print verbose description of presented statistics
611 """
612 if len(self.sets) == 0:
613 return "Empty"
614
615 self.compute()
616
617
618 labels = self.__labels
619 labels_map_rev = self.__labels_map_rev
620 matrix = self.__matrix
621
622 labels_rev = []
623 if labels_map_rev is not None:
624 labels_rev = [','.join([str(x) for x in labels_map_rev[l]])
625 for l in labels]
626
627 out = StringIO()
628
629 Nlabels = len(labels)
630 Nsamples = self.__Nsamples.astype(int)
631
632 stats = self._stats
633 if short:
634 return "%(# of sets)d sets %(# of labels)d labels " \
635 " ACC:%(ACC).2f" \
636 % stats
637
638 Ndigitsmax = int(ceil(log10(max(Nsamples))))
639 Nlabelsmax = max( [len(str(x)) for x in labels] )
640
641
642 L = max(Ndigitsmax+2, Nlabelsmax)
643 res = ""
644
645 stats_perpredict = ["P'", "N'", 'FP', 'FN', 'PPV', 'NPV', 'TPR',
646 'SPC', 'FDR', 'MCC']
647
648 if self.ROC is not None: stats_perpredict += [ 'AUC' ]
649 stats_pertarget = ['P', 'N', 'TP', 'TN']
650 stats_summary = ['ACC', 'ACC%', '# of sets']
651
652
653
654 prefixlen = Nlabelsmax + 1
655 pref = ' '*(prefixlen)
656
657 if matrix.shape != (Nlabels, Nlabels):
658 raise ValueError, \
659 "Number of labels %d doesn't correspond the size" + \
660 " of a confusion matrix %s" % (Nlabels, matrix.shape)
661
662
663 printed = []
664 underscores = [" %s" % ("-" * L)] * Nlabels
665 if header:
666
667 printed.append(['@l----------. '] + labels_rev)
668 printed.append(['@lpredictions\\targets'] + labels)
669
670 printed.append(['@l `------'] \
671 + underscores + stats_perpredict)
672
673
674 for i, line in enumerate(matrix):
675 l = labels[i]
676 if labels_rev != []:
677 l = '@r%10s / %s' % (labels_rev[i], l)
678 printed.append(
679 [l] +
680 [ str(x) for x in line ] +
681 [ _p2(stats[x][i]) for x in stats_perpredict])
682
683 if summary:
684
685
686
687
688
689
690 printed.append(['@lPer target:'] + underscores)
691 for stat in stats_pertarget:
692 printed.append([stat] + [
693 _p2(stats[stat][i]) for i in xrange(Nlabels)])
694
695
696
697
698 mean_stats = N.mean(N.array([stats[k] for k in stats_perpredict]),
699 axis=1)
700 printed.append(['@lSummary \ Means:'] + underscores
701 + [_p2(stats['mean(%s)' % x])
702 for x in stats_perpredict])
703
704 for stat in stats_summary:
705 printed.append([stat] + [_p2(stats[stat])])
706
707 table2string(printed, out)
708
709 if description:
710 out.write("\nStatistics computed in 1-vs-rest fashion per each " \
711 "target.\n")
712 out.write("Abbreviations (for details see " \
713 "http://en.wikipedia.org/wiki/ROC_curve):\n")
714 for d, val, eq in self._STATS_DESCRIPTION:
715 out.write(" %-3s: %s\n" % (d, val))
716 if eq is not None:
717 out.write(" " + eq + "\n")
718
719
720 result = out.getvalue()
721 out.close()
722 return result
723
724
725 - def plot(self, labels=None, numbers=False, origin='upper',
726 numbers_alpha=None, xlabels_vertical=True, numbers_kwargs={},
727 **kwargs):
728 """Provide presentation of confusion matrix in image
729
730 :Parameters:
731 labels : list of int or basestring
732 Optionally provided labels guarantee the order of
733 presentation. Also value of None places empty column/row,
734 thus provides visual groupping of labels (Thanks Ingo)
735 numbers : bool
736 Place values inside of confusion matrix elements
737 numbers_alpha : None or float
738 Controls textual output of numbers. If None -- all numbers
739 are plotted in the same intensity. If some float -- it controls
740 alpha level -- higher value would give higher contrast. (good
741 value is 2)
742 origin : basestring
743 Which left corner diagonal should start
744 xlabels_vertical : bool
745 Either to plot xlabels vertical (benefitial if number of labels
746 is large)
747 numbers_kwargs : dict
748 Additional keyword parameters to be added to numbers (if numbers
749 is True)
750 **kwargs
751 Additional arguments given to imshow (\eg me cmap)
752
753 :Returns:
754 (fig, im, cb) -- figure, imshow, colorbar
755 """
756
757 externals.exists("pylab", raiseException=True)
758 import pylab as P
759
760 self.compute()
761 labels_order = labels
762
763
764 labels = self.__labels
765 labels_map = self.__labels_map
766 labels_map_rev = self.__labels_map_rev
767 matrix = self.__matrix
768
769
770 labels_indexes = dict([(x,i) for i,x in enumerate(labels)])
771
772 labels_rev = []
773 if labels_map_rev is not None:
774 labels_rev = [','.join([str(x) for x in labels_map_rev[l]])
775 for l in labels]
776
777 if labels_order is not None:
778 labels_order_filtered = filter(lambda x:x is not None, labels_order)
779 labels_order_filtered_set = Set(labels_order_filtered)
780
781 if Set(labels) == labels_order_filtered_set:
782
783 labels_plot = labels_order
784 elif len(labels_rev) \
785 and Set(labels_map.keys()) == labels_order_filtered_set:
786
787
788 labels_plot = []
789 for l in labels_order:
790 v = None
791 if l is not None: v = labels_map[l]
792 labels_plot += [v]
793 else:
794 raise ValueError, \
795 "Provided labels %s do not match set of known " \
796 "original labels (%s) or mapped labels (%s)" % \
797 (labels_order, labels, labels_rev)
798 else:
799 labels_plot = labels
800
801
802 isempty = N.array([l is None for l in labels_plot])
803 non_empty = N.where(N.logical_not(isempty))[0]
804
805 NlabelsNN = len(non_empty)
806 Nlabels = len(labels_plot)
807
808 if matrix.shape != (NlabelsNN, NlabelsNN):
809 raise ValueError, \
810 "Number of labels %d doesn't correspond the size" + \
811 " of a confusion matrix %s" % (NlabelsNN, matrix.shape)
812
813 confusionmatrix = N.zeros((Nlabels, Nlabels))
814 mask = confusionmatrix.copy()
815 ticks = []
816 tick_labels = []
817
818
819 for i, l in enumerate(labels_plot):
820 if l is not None:
821 j = labels_indexes[l]
822 confusionmatrix[i, non_empty] = matrix[j, :]
823 confusionmatrix[non_empty, i] = matrix[:, j]
824 ticks += [i + 0.5]
825 if labels_map_rev is not None:
826 tick_labels += ['/'.join(labels_map_rev[l])]
827 else:
828 tick_labels += [str(l)]
829 else:
830 mask[i, :] = mask[:, i] = 1
831
832 confusionmatrix = N.ma.MaskedArray(confusionmatrix, mask=mask)
833
834
835 if P.matplotlib.get_backend() == 'TkAgg':
836 P.ioff()
837
838 fig = P.gcf()
839 ax = P.gca()
840 ax.axis('off')
841
842
843 xticks_position, yticks, ybottom = {
844 'upper': ('top', ticks[::-1], 0.1),
845 'lower': ('bottom', ticks, 0.2)
846 }[origin]
847
848
849
850 axi = fig.add_axes([0.15, ybottom, 0.7, 0.7])
851 im = axi.imshow(confusionmatrix, interpolation="nearest", origin=origin,
852 aspect='equal', **kwargs)
853
854
855 if numbers:
856 numbers_kwargs_ = {'fontsize': 10,
857 'horizontalalignment': 'center',
858 'verticalalignment': 'center'}
859 maxv = float(N.max(confusionmatrix))
860 colors = [im.to_rgba(0), im.to_rgba(maxv)]
861 for i,j in zip(*N.logical_not(mask).nonzero()):
862 v = confusionmatrix[j, i]
863
864 if numbers_alpha is None:
865 alpha = 1.0
866 else:
867
868 alpha = 1 - N.array(1 - v / maxv) ** numbers_alpha
869 y = {'lower':j, 'upper':Nlabels-j-1}[origin]
870 numbers_kwargs_['color'] = colors[int(v<maxv/2)]
871 numbers_kwargs_.update(numbers_kwargs)
872 P.text(i+0.5, y+0.5, '%d' % v, alpha=alpha, **numbers_kwargs_)
873
874 maxv = N.max(confusionmatrix)
875 boundaries = N.linspace(0, maxv, N.min(maxv, 10), True)
876
877
878 P.xlabel("targets")
879 P.ylabel("predictions")
880
881 P.setp(axi, xticks=ticks, yticks=yticks,
882 xticklabels=tick_labels, yticklabels=tick_labels)
883
884 axi.xaxis.set_ticks_position(xticks_position)
885 axi.xaxis.set_label_position(xticks_position)
886
887 if xlabels_vertical:
888 P.setp(P.getp(axi, 'xticklabels'), rotation='vertical')
889
890 axcb = fig.add_axes([0.8, ybottom, 0.02, 0.7])
891 cb = P.colorbar(im, cax=axcb, format='%d', ticks = boundaries)
892
893 if P.matplotlib.get_backend() == 'TkAgg':
894 P.ion()
895 P.draw()
896 return fig, im, cb
897
898
899 @property
901 self.compute()
902 return 1.0-self.__Ncorrect*1.0/sum(self.__Nsamples)
903
904
905 @property
907 self.compute()
908 return self.__labels
909
910
912 return self.__labels_map
913
914
916 if val is None or isinstance(val, dict):
917 self.__labels_map = val
918 else:
919 raise ValueError, "Cannot set labels_map to %s" % val
920
921 self.__labels_map_rev = None
922 self._computed = False
923
924
925 @property
927 self.compute()
928 return self.__matrix
929
930
931 @property
933 self.compute()
934 return 100.0*self.__Ncorrect/sum(self.__Nsamples)
935
936 labels_map = property(fget=getLabels_map, fset=setLabels_map)
937
940 """Class to contain information and display on regression results.
941
942 """
943
944 _STATS_DESCRIPTION = (
945 ('CCe', 'Error based on correlation coefficient',
946 '1 - corr_coef'),
947 ('CCp', 'Correlation coefficient (p-value)', None),
948 ('RMSE', 'Root mean squared error', None),
949 ('STD', 'Standard deviation', None),
950 ('RMP', 'Root mean power (compare to RMSE of results)',
951 'sqrt(mean( data**2 ))'),
952 ) + SummaryStatistics._STATS_DESCRIPTION
953
954
956 """Initialize RegressionStatistics
957
958 :Parameters:
959 targets
960 Optional set of targets
961 predictions
962 Optional set of predictions
963 """
964
965 SummaryStatistics.__init__(self, **kwargs)
966
967
969 """Actually compute the confusion matrix based on all the sets"""
970
971 super(RegressionStatistics, self)._compute()
972 sets = self.sets
973 Nsets = len(sets)
974
975 stats = {}
976
977 funcs = {
978 'RMP_t': lambda p,t:rootMeanPowerFx(t),
979 'STD_t': lambda p,t:N.std(t),
980 'RMP_p': lambda p,t:rootMeanPowerFx(p),
981 'STD_p': lambda p,t:N.std(p),
982 'CCe': CorrErrorFx(),
983 'CCp': CorrErrorPFx(),
984 'RMSE': RMSErrorFx(),
985 'RMSE/RMP_t': RelativeRMSErrorFx()
986 }
987
988 for funcname, func in funcs.iteritems():
989 funcname_all = funcname + '_all'
990 stats[funcname_all] = []
991 for i, (targets, predictions, values) in enumerate(sets):
992 stats[funcname_all] += [func(predictions, targets)]
993 stats[funcname_all] = N.array(stats[funcname_all])
994 stats[funcname] = N.mean(stats[funcname_all])
995 stats[funcname+'_std'] = N.std(stats[funcname_all])
996 stats[funcname+'_max'] = N.max(stats[funcname_all])
997 stats[funcname+'_min'] = N.min(stats[funcname_all])
998
999 self._stats.update(stats)
1000
1001
1002 - def asstring(self, short=False, header=True, summary=True,
1003 description=False):
1004 """'Pretty print' the statistics"""
1005
1006 if len(self.sets) == 0:
1007 return "Empty"
1008
1009 self.compute()
1010
1011 stats = self.stats
1012
1013 if short:
1014 if short == 'very':
1015 return "%(# of sets)d sets CCe:%(CCe).2f CCp:%(CCp).2g" \
1016 " RMSE:%(RMSE).2f" \
1017 " RMSE/RMP_t:%(RMSE/RMP_t).2f" \
1018 % stats
1019 else:
1020 return "%(# of sets)d sets CCe:%(CCe).2f+-%(CCe_std).3f" \
1021 " RMSE:%(RMSE).2f+-%(RMSE_std).3f" \
1022 " RMSE/RMP_t:%(RMSE/RMP_t).2f+-%(RMSE/RMP_t_std).3f" \
1023 % stats
1024
1025 stats_data = ['RMP_t', 'STD_t', 'RMP_p', 'STD_p']
1026
1027 stats_ = ['CCe', 'RMSE', 'RMSE/RMP_t']
1028 stats_summary = ['# of sets']
1029
1030 out = StringIO()
1031
1032 printed = []
1033 if header:
1034
1035 printed.append(['Statistics', 'Mean', 'Std', 'Min', 'Max'])
1036
1037 printed.append(['----------', '-----', '-----', '-----', '-----'])
1038
1039 def print_stats(printed, stats_):
1040
1041 for stat in stats_:
1042 s = [stat]
1043 for suffix in ['', '_std', '_min', '_max']:
1044 s += [ _p2(stats[stat+suffix], 3) ]
1045 printed.append(s)
1046
1047 printed.append(["Data: "])
1048 print_stats(printed, stats_data)
1049 printed.append(["Results: "])
1050 print_stats(printed, stats_)
1051
1052 if summary:
1053 for stat in stats_summary:
1054 printed.append([stat] + [_p2(stats[stat])])
1055
1056 table2string(printed, out)
1057
1058 if description:
1059 out.write("\nDescription of printed statistics.\n"
1060 " Suffixes: _t - targets, _p - predictions\n")
1061
1062 for d, val, eq in self._STATS_DESCRIPTION:
1063 out.write(" %-3s: %s\n" % (d, val))
1064 if eq is not None:
1065 out.write(" " + eq + "\n")
1066
1067 result = out.getvalue()
1068 out.close()
1069 return result
1070
1071
1072 @property
1076
1080 """Compute (or return) some error of a (trained) classifier on a dataset.
1081 """
1082
1083 confusion = StateVariable(enabled=False)
1084 """TODO Think that labels might be also symbolic thus can't directly
1085 be indicies of the array
1086 """
1087
1088 training_confusion = StateVariable(enabled=False,
1089 doc="Proxy training_confusion from underlying classifier.")
1090
1091
1092 - def __init__(self, clf, labels=None, train=True, **kwargs):
1093 """Initialization.
1094
1095 :Parameters:
1096 clf : Classifier
1097 Either trained or untrained classifier
1098 labels : list
1099 if provided, should be a set of labels to add on top of the
1100 ones present in testdata
1101 train : bool
1102 unless train=False, classifier gets trained if
1103 trainingdata provided to __call__
1104 """
1105 Stateful.__init__(self, **kwargs)
1106 self.__clf = clf
1107
1108 self._labels = labels
1109 """Labels to add on top to existing in testing data"""
1110
1111 self.__train = train
1112 """Either to train classifier if trainingdata is provided"""
1113
1114
1115 __doc__ = enhancedDocString('ClassifierError', locals(), Stateful)
1116
1117
1123
1124
1125 - def _precall(self, testdataset, trainingdataset=None):
1126 """Generic part which trains the classifier if necessary
1127 """
1128 if not trainingdataset is None:
1129 if self.__train:
1130
1131
1132
1133
1134
1135
1136
1137
1138
1139
1140
1141 if self.states.isEnabled('training_confusion'):
1142 self.__clf.states._changeTemporarily(
1143 enable_states=['training_confusion'])
1144 self.__clf.train(trainingdataset)
1145 if self.states.isEnabled('training_confusion'):
1146 self.training_confusion = self.__clf.training_confusion
1147 self.__clf.states._resetEnabledTemporarily()
1148
1149 if self.__clf.states.isEnabled('trained_labels') and \
1150 not testdataset is None:
1151 newlabels = Set(testdataset.uniquelabels) \
1152 - Set(self.__clf.trained_labels)
1153 if len(newlabels)>0:
1154 warning("Classifier %s wasn't trained to classify labels %s" %
1155 (`self.__clf`, `newlabels`) +
1156 " present in testing dataset. Make sure that you have" +
1157 " not mixed order/names of the arguments anywhere")
1158
1159
1160
1161
1162
1163
1164
1165
1166
1167
1168
1169
1170 - def _call(self, testdataset, trainingdataset=None):
1171 raise NotImplementedError
1172
1173
1174 - def _postcall(self, testdataset, trainingdataset=None, error=None):
1176
1177
1178 - def __call__(self, testdataset, trainingdataset=None):
1179 """Compute the transfer error for a certain test dataset.
1180
1181 If `trainingdataset` is not `None` the classifier is trained using the
1182 provided dataset before computing the transfer error. Otherwise the
1183 classifier is used in it's current state to make the predictions on
1184 the test dataset.
1185
1186 Returns a scalar value of the transfer error.
1187 """
1188 self._precall(testdataset, trainingdataset)
1189 error = self._call(testdataset, trainingdataset)
1190 self._postcall(testdataset, trainingdataset, error)
1191 return error
1192
1193
1194 @property
1197
1198
1199 @property
1202
1206 """Compute the transfer error of a (trained) classifier on a dataset.
1207
1208 The actual error value is computed using a customizable error function.
1209 Optionally the classifier can be trained by passing an additional
1210 training dataset to the __call__() method.
1211 """
1212
1213 null_prob = StateVariable(enabled=True,
1214 doc="Stores the probability of an error result under "
1215 "the NULL hypothesis")
1216 samples_error = StateVariable(enabled=False,
1217 doc="Per sample errors computed by invoking the "
1218 "error function for each sample individually. "
1219 "Errors are available in a dictionary with each "
1220 "samples origid as key.")
1221
1224 """Initialization.
1225
1226 :Parameters:
1227 clf : Classifier
1228 Either trained or untrained classifier
1229 errorfx
1230 Functor that computes a scalar error value from the vectors of
1231 desired and predicted values (e.g. subclass of `ErrorFunction`)
1232 labels : list
1233 if provided, should be a set of labels to add on top of the
1234 ones present in testdata
1235 null_dist : instance of distribution estimator
1236 """
1237 ClassifierError.__init__(self, clf, labels, **kwargs)
1238 self.__errorfx = errorfx
1239 self.__null_dist = autoNullDist(null_dist)
1240
1241
1242 __doc__ = enhancedDocString('TransferError', locals(), ClassifierError)
1243
1244
1253
1254
1255 - def _call(self, testdataset, trainingdataset=None):
1256 """Compute the transfer error for a certain test dataset.
1257
1258 If `trainingdataset` is not `None` the classifier is trained using the
1259 provided dataset before computing the transfer error. Otherwise the
1260 classifier is used in it's current state to make the predictions on
1261 the test dataset.
1262
1263 Returns a scalar value of the transfer error.
1264 """
1265
1266 clf = self.clf
1267 if testdataset is None:
1268
1269
1270 import traceback as tb
1271 filenames = [x[0] for x in tb.extract_stack(limit=100)]
1272 rfe_matches = [f for f in filenames if f.endswith('/rfe.py')]
1273 cv_matches = [f for f in filenames if
1274 f.endswith('cvtranserror.py')]
1275 msg = ""
1276 if len(rfe_matches) > 0 and len(cv_matches):
1277 msg = " It is possible that you used RFE with stopping " \
1278 "criterion based on the TransferError and directly" \
1279 " from CrossValidatedTransferError, such approach" \
1280 " would require exposing testing dataset " \
1281 " to the classifier which might heavily bias " \
1282 " generalization performance estimate. If you are " \
1283 " sure to use it that way, create CVTE with " \
1284 " parameter expose_testdataset=True"
1285 raise ValueError, "Transfer error call obtained None " \
1286 "as a dataset for testing.%s" % msg
1287 predictions = clf.predict(testdataset.samples)
1288
1289
1290
1291
1292
1293
1294
1295 states = self.states
1296 if states.isEnabled('confusion'):
1297 confusion = clf._summaryClass(
1298
1299 targets=testdataset.labels,
1300 predictions=predictions,
1301 values=clf.states.get('values', None))
1302 try:
1303 confusion.labels_map = testdataset.labels_map
1304 except:
1305 pass
1306 states.confusion = confusion
1307
1308 if states.isEnabled('samples_error'):
1309 samples_error = []
1310 for i, p in enumerate(predictions):
1311 samples_error.append(self.__errorfx(p, testdataset.labels[i]))
1312
1313 states.samples_error = dict(zip(testdataset.origids, samples_error))
1314
1315
1316 error = self.__errorfx(predictions, testdataset.labels)
1317
1318 return error
1319
1320
1321 - def _postcall(self, vdata, wdata=None, error=None):
1322 """
1323 """
1324
1325
1326 if not self.__null_dist is None and not wdata is None:
1327
1328
1329
1330 null_terr = copy.copy(self)
1331 null_terr.__null_dist = None
1332 self.__null_dist.fit(null_terr, wdata, vdata)
1333
1334
1335
1336 if not error is None and not self.__null_dist is None:
1337 self.null_prob = self.__null_dist.p(error)
1338
1339
1340 @property
1341 - def errorfx(self): return self.__errorfx
1342
1343 @property
1344 - def null_dist(self): return self.__null_dist
1345
1349 """For a given classifier report an error based on internally
1350 computed error measure (given by some `ConfusionMatrix` stored in
1351 some state variable of `Classifier`).
1352
1353 This way we can perform feature selection taking as the error
1354 criterion either learning error, or transfer to splits error in
1355 the case of SplitClassifier
1356 """
1357
1358 - def __init__(self, clf, labels=None, confusion_state="training_confusion",
1359 **kwargs):
1360 """Initialization.
1361
1362 :Parameters:
1363 clf : Classifier
1364 Either trained or untrained classifier
1365 confusion_state
1366 Id of the state variable which stores `ConfusionMatrix`
1367 labels : list
1368 if provided, should be a set of labels to add on top of the
1369 ones present in testdata
1370 """
1371 ClassifierError.__init__(self, clf, labels, **kwargs)
1372
1373 self.__confusion_state = confusion_state
1374 """What state to extract from"""
1375
1376 if not clf.states.isKnown(confusion_state):
1377 raise ValueError, \
1378 "State variable %s is not defined for classifier %s" % \
1379 (confusion_state, `clf`)
1380 if not clf.states.isEnabled(confusion_state):
1381 if __debug__:
1382 debug('CERR', "Forcing state %s to be enabled for %s" %
1383 (confusion_state, `clf`))
1384 clf.states.enable(confusion_state)
1385
1386
1387 __doc__ = enhancedDocString('ConfusionBasedError', locals(),
1388 ClassifierError)
1389
1390
1391 - def _call(self, testdata, trainingdata=None):
1397