Repository: vanvalenlab/deepcell-toolbox (language: python)

Issue:
Metrics package is slow
The metrics package is quite slow, and takes an order of magnitude longer to run than making the segmentation predictions themselves. This is likely because of scaling based on the pairwise interaction between the number of cells in each crop. Not pressing, but would be good to see if there's any optimization that can be done


Discussion / hints:
The metrics package also spams warnings during evaluation (using TensorFlow 2.x). The following message is spammed over and over:

> /usr/local/lib/python3.6/dist-packages/deepcell_toolbox/metrics.py:115: RuntimeWarning: invalid value encountered in long_scalars
>  precision = intersection.sum() / pred.sum()

Fix patch:
diff --git a/.coveragerc b/.coveragerc
index 47ec5a5..8bb9a66 100644
--- a/.coveragerc
+++ b/.coveragerc
@@ -14,3 +14,4 @@ show_missing = True
 
 omit =
     **/*_test.py
+    deepcell_toolbox/__version__.py
diff --git a/deepcell_toolbox/metrics.py b/deepcell_toolbox/metrics.py
index 99eefcc..860be97 100644
--- a/deepcell_toolbox/metrics.py
+++ b/deepcell_toolbox/metrics.py
@@ -42,7 +42,6 @@
 from __future__ import division
 
 import datetime
-import decimal
 import json
 import logging
 import operator
@@ -54,18 +53,144 @@
 import networkx as nx
 
 from scipy.optimize import linear_sum_assignment
+from scipy.stats import hmean
 from skimage.measure import regionprops
 from skimage.segmentation import relabel_sequential
 from sklearn.metrics import confusion_matrix
-from scipy.stats import hmean
+from tqdm import tqdm
 
 from deepcell_toolbox import erode_edges
 from deepcell_toolbox.compute_overlap import compute_overlap  # pylint: disable=E0401
 from deepcell_toolbox.compute_overlap_3D import compute_overlap_3D
 
 
-def stats_pixelbased(y_true, y_pred):
-    """Calculates pixel-based statistics
+def _cast_to_tuple(x):
+    try:
+        tup_x = tuple(x)
+    except TypeError:
+        tup_x = () if x is None else (x,)
+    return tup_x
+
+
+class Detection(object):  # pylint: disable=useless-object-inheritance
+    """Object to hold relevant information about a given detection."""
+
+    def __init__(self, true_index=None, pred_index=None):
+        # cast the indices as tuples if possible to make them immutable
+        try:
+            self.true_index = tuple(true_index)
+        except TypeError:
+            self.true_index = true_index
+        try:
+            self.pred_index = tuple(pred_index)
+        except TypeError:
+            self.pred_index = pred_index
+
+    def __eq__(self, other):
+        """Custom comparator. Detections with the same indices are the same."""
+        try:
+            is_true_same = self.true_index == other.true_index
+            is_pred_same = self.pred_index == other.pred_index
+            return is_true_same and is_pred_same
+        except AttributeError:
+            return False
+
+    def __hash__(self):
+        """Custom hasher, allow Detections to be hashable."""
+        return tuple((self.true_index, self.pred_index)).__hash__()
+
+    def __repr__(self):
+        return 'Detection({}, {})'.format(self.true_index, self.pred_index)
+
+    @property
+    def is_correct(self):
+        is_linked = self.true_index is not None and self.pred_index is not None
+        return is_linked and not self.is_split and not self.is_merge
+
+    @property
+    def is_gained(self):
+        return self.true_index is None and self.pred_index is not None
+
+    @property
+    def is_missed(self):
+        return self.true_index is not None and self.pred_index is None
+
+    @property
+    def is_split(self):
+        if self.is_gained or self.is_missed:
+            return False
+
+        try:
+            is_many_pred = len(self.pred_index) > 1
+        except TypeError:
+            is_many_pred = False
+
+        try:
+            is_single_true = len(tuple(self.true_index)) == 1
+        except TypeError:
+            is_single_true = isinstance(self.true_index, int)
+
+        return is_single_true and is_many_pred
+
+    @property
+    def is_merge(self):
+        if self.is_gained or self.is_missed:
+            return False
+
+        try:
+            is_many_true = len(self.true_index) > 1
+        except TypeError:
+            is_many_true = False
+
+        try:
+            is_single_pred = len(tuple(self.pred_index)) == 1
+        except TypeError:
+            is_single_pred = isinstance(self.pred_index, int)
+
+        return is_single_pred and is_many_true
+
+    @property
+    def is_catastrophe(self):
+        if self.is_gained or self.is_missed:
+            return False
+
+        try:
+            is_many_true = len(self.true_index) > 1
+        except TypeError:
+            is_many_true = False
+
+        try:
+            is_many_pred = len(self.pred_index) > 1
+        except TypeError:
+            is_many_pred = False
+
+        return is_many_true and is_many_pred
+
+
+class BaseMetrics(object):  # pylint: disable=useless-object-inheritance
+
+    """Base class for Metrics classes."""
+
+    def __init__(self, y_true, y_pred):
+        if y_pred.shape != y_true.shape:
+            raise ValueError('Input shapes must match. Shape of prediction '
+                             'is: {}.  Shape of y_true is: {}'.format(
+                                 y_pred.shape, y_true.shape))
+
+        if not np.issubdtype(y_true.dtype, np.integer):
+            warnings.warn('Casting y_true from {} to int'.format(y_true.dtype))
+            y_true = y_true.astype('int32')
+
+        if not np.issubdtype(y_pred.dtype, np.integer):
+            warnings.warn('Casting y_pred from {} to int'.format(y_pred.dtype))
+            y_pred = y_pred.astype('int32')
+
+        self.y_true = y_true
+        self.y_pred = y_pred
+
+
+class PixelMetrics(BaseMetrics):
+    """Calculates pixel-based statistics.
     (Dice, Jaccard, Precision, Recall, F-measure)
 
     Takes in raw prediction and truth data in order to calculate accuracy
@@ -80,9 +205,6 @@ def stats_pixelbased(y_true, y_pred):
         y_pred (numpy.array): Binary predictions for a single feature,
             (batch,x,y)
 
-    Returns:
-        dict: Containing a set of calculated statistics
-
     Raises:
         ValueError: Shapes of y_true and y_pred do not match.
 
@@ -91,39 +213,114 @@ def stats_pixelbased(y_true, y_pred):
         Make sure to input the same type of data for y_true and y_pred
     """
 
-    if y_pred.shape != y_true.shape:
-        raise ValueError('Shape of inputs need to match. Shape of prediction '
-                         'is: {}.  Shape of y_true is: {}'.format(
-                             y_pred.shape, y_true.shape))
+    def __init__(self, y_true, y_pred):
+        super(PixelMetrics, self).__init__(
+            y_true=(y_true != 0).astype('int'),
+            y_pred=(y_pred != 0).astype('int'))
+
+        self._y_true_sum = np.count_nonzero(self.y_true)
+        self._y_pred_sum = np.count_nonzero(self.y_pred)
+
+        # Calculations for IOU
+        self._intersection = np.count_nonzero(np.logical_and(self.y_true, self.y_pred))
+        self._union = np.count_nonzero(np.logical_or(self.y_true, self.y_pred))
+
+    @classmethod
+    def get_confusion_matrix(cls, y_true, y_pred, axis=-1):
+        """Calculate confusion matrix for pixel classification data.
+
+        Args:
+            y_true (numpy.array): Ground truth annotations after any
+                necessary transformations
+            y_pred (numpy.array): Prediction array
+            axis (int): The channel axis of the input arrays.
+
+        Returns:
+            numpy.array: nxn confusion matrix determined by number of features.
+        """
+        # Argmax collapses on feature dimension to assign class to each pixel
+        # Flatten is required for confusion matrix
+        y_true = y_true.argmax(axis=axis).flatten()
+        y_pred = y_pred.argmax(axis=axis).flatten()
+        return confusion_matrix(y_true, y_pred)
+
+    @property
+    def recall(self):
+        try:
+            _recall = self._intersection / self._y_true_sum
+        except ZeroDivisionError:
+            _recall = np.nan
+        return _recall
+
+    @property
+    def precision(self):
+        try:
+            _precision = self._intersection / self._y_pred_sum
+        except ZeroDivisionError:
+            _precision = 0
+        return _precision
+
+    @property
+    def f1(self):
+        _recall = self.recall
+        _precision = self.precision
+
+        # f1 is nan if recall is nan and no false negatives
+        if np.isnan(_recall) and _precision == 0:
+            return np.nan
+
+        f_measure = hmean([_recall, _precision])
+        # f_measure = (2 * _precision * _recall) / (_precision + _recall)
+        return f_measure
+
+    @property
+    def dice(self):
+        y_sum = self._y_true_sum + self._y_pred_sum
+        if y_sum == 0:
+            warnings.warn('DICE score is technically 1.0, '
+                          'but prediction and truth arrays are empty.')
+            return 1.0
+
+        return 2.0 * self._intersection / y_sum
+
+    @property
+    def jaccard(self):
+        try:
+            _jaccard = self._intersection / self._union
+        except ZeroDivisionError:
+            _jaccard = np.nan
+        return _jaccard
 
-    pred = y_pred
-    truth = y_true
+    def to_dict(self):
+        return {
+            'jaccard': self.jaccard,
+            'recall': self.recall,
+            'precision': self.precision,
+            'f1': self.f1,
+            'dice': self.dice,
+        }
 
-    if pred.sum() == 0 and truth.sum() == 0:
-        warnings.warn('DICE score is technically 1.0, '
-                      'but prediction and truth arrays are empty. ')
 
-    # Calculations for IOU
-    intersection = np.logical_and(pred, truth)
-    union = np.logical_or(pred, truth)
+def get_box_labels(arr):
+    """Get the bounding box and label for all objects in the image.
 
-    # Sum gets count of positive pixels
-    dice = (2 * intersection.sum() / (pred.sum() + truth.sum()))
-    jaccard = intersection.sum() / union.sum()
-    precision = intersection.sum() / pred.sum()
-    recall = intersection.sum() / truth.sum()
-    Fmeasure = (2 * precision * recall) / (precision + recall)
+    Args:
+        arr (np.array): integer label array of objects.
 
-    return {
-        'dice': dice,
-        'jaccard': jaccard,
-        'precision': precision,
-        'recall': recall,
-        'Fmeasure': Fmeasure
-    }
+    Returns:
+        tuple(list(np.array), list(int)): A tuple of bounding boxes and
+            the corresponding integer labels.
+    """
+    props = regionprops(np.squeeze(arr.astype('int')), cache=False)
+    boxes, labels = [], []
+    for prop in props:
+        boxes.append(np.array(prop.bbox))
+        labels.append(int(prop.label))
+    boxes = np.array(boxes).astype('double')
+    return boxes, labels
 
 
-class ObjectAccuracy(object):  # pylint: disable=useless-object-inheritance
+class ObjectMetrics(BaseMetrics):
     """Classifies object prediction errors as TP, FP, FN, merge or split
 
     The schema for this analysis was adopted from the description of
@@ -145,8 +342,6 @@ class ObjectAccuracy(object):  # pylint: disable=useless-object-inheritance
             smaller values are more conservative, default 0.4
         cutoff2 (:obj:`float`, optional): Threshold for overlap in unassigned
             cells, smaller values are better, default 0.1
-        test (:obj:`bool`, optional): Utility variable to control running
-            analysis during testing
         seg (:obj:`bool`, optional): Calculates SEG score for cell tracking
             competition
         force_event_links(:obj:'bool, optional): Flag that determines whether to modify IOU
@@ -165,216 +360,197 @@ def __init__(self,
                  y_pred,
                  cutoff1=0.4,
                  cutoff2=0.1,
-                 test=False,
-                 seg=False,
                  force_event_links=False,
                  is_3d=False):
-        self.cutoff1 = cutoff1
-        self.cutoff2 = cutoff2
-        self.seg = seg
-        self.is_3d = is_3d
 
-        if y_pred.shape != y_true.shape:
-            raise ValueError('Input shapes must match. Shape of prediction '
-                             'is: {}.  Shape of y_true is: {}'.format(
-                                 y_pred.shape, y_true.shape))
+        # If 2D, dimensions can be 3 or 4 (with or without channel dimension)
+        if not is_3d and y_true.ndim not in {2, 3}:
+            raise ValueError('Expected dimensions for y_true (2D data) are 2 '
+                             '(x, y) and 3 (x, y, chan). '
+                             'Got ndim: {}'.format(y_true.ndim))
 
-        if not np.issubdtype(y_true.dtype, np.integer):
-            warnings.warn('Casting y_true from {} to int'.format(y_true.dtype))
-            y_true = y_true.astype('int32')
+        elif is_3d and y_true.ndim != 3:
+            raise ValueError('Expected dimensions for y_true (3D data) is 3.'
+                             'Requires format is: (z, x, y)'
+                             'Got ndim: {}'.format(y_true.ndim))
 
-        if not np.issubdtype(y_pred.dtype, np.integer):
-            warnings.warn('Casting y_pred from {} to int'.format(y_pred.dtype))
-            y_pred = y_pred.astype('int32')
+        super(ObjectMetrics, self).__init__(y_true=y_true, y_pred=y_pred)
 
-        self.y_true = y_true
-        self.y_pred = y_pred
+        self.cutoff1 = cutoff1
+        self.cutoff2 = cutoff2
+        self.is_3d = is_3d
+
+        self.compute_overlap = compute_overlap_3D if is_3d else compute_overlap
 
         self.n_true = len(np.unique(self.y_true[np.nonzero(self.y_true)]))
         self.n_pred = len(np.unique(self.y_pred[np.nonzero(self.y_pred)]))
 
-        self.n_obj = self.n_true + self.n_pred
-
-        # Initialize error counters
-        self.correct_detections = 0
-        self.missed_detections = 0
-        self.gained_detections = 0
-
-        self.merge = 0
-        self.split = 0
-        self.catastrophe = 0
+        # keep track of every pair of objects through the detections dict
+        # using tuple(true_index, pred_index): Detection as a key/vaue pair
+        self._detections = set()
+
+        # store the keys of relevant Detections in a set for easy fetching
+        # types of detections
+        self._splits = set()
+        self._gained = set()
+        self._missed = set()
+        # types of errors
+        self._merges = set()
+        self._catastrophes = set()
+        self._correct = set()
+
+        # IoU: used to determine relative overlap of y_pred and y_true
+        self.iou = np.zeros((self.n_true, self.n_pred))
 
-        self.gained_det_from_split = 0
-        self.missed_det_from_merge = 0
-        self.true_det_in_catastrophe = 0
-        self.pred_det_in_catastrophe = 0
+        # used to determine seg score
+        self.seg_thresh = np.zeros((self.n_true, self.n_pred))
 
-        self.precision = 0
-        self.recall = 0
-        self.f1 = 0
+        # Check if either frame is empty before proceeding
+        if self.n_true == 0:
+            logging.info('Ground truth frame is empty')
 
-        # Initialize lists and dicts to store indices where errors occur
-        self.correct_indices = {}
-        self.correct_indices['y_true'] = []
-        self.correct_indices['y_pred'] = []
+        if self.n_pred == 0:
+            logging.info('Prediction frame is empty')
 
-        self.missed_indices = {}
-        self.missed_indices['y_true'] = []
+        self._calc_iou()  # set self.iou and update self.seg_thresh
 
-        self.gained_indices = {}
-        self.gained_indices['y_pred'] = []
+        self.iou_modified = self._get_modified_iou(force_event_links)
 
-        self.merge_indices = {
-            'y_true': [],
-            'y_pred': []
-        }
+        matrix = self._linear_assignment()
 
-        self.split_indices = {
-            'y_true': [],
-            'y_pred': []
-        }
-
-        self.catastrophe_indices = {
-            'y_true': []
-        }
-        self.catastrophe_indices['y_pred'] = []
-
-        # If 2D, dimensions can be 2 or 3 (with or without channel dimension)
-        if not self.is_3d:
-            if self.y_true.ndim not in {2, 3}:
-                raise ValueError('Expected dimensions for y_true (2D data) are 2 and 3'
-                                 'Accepts: (x, y), or (x, y, chan)'
-                                 'Got ndim: {}'.format(self.y_true.ndim))
+        # Identify direct matches as true positives
+        correct_index = np.nonzero(matrix[:self.n_true, :self.n_pred])
 
-        # If 3D, inputs must have 3 dimensions (batch, z, x, y) - cannot have channel dimension or
-        # _classify_graph breaks, as it expects input to be 2D or 3D
-        # TODO - add compatibility for multi-channel 3D-data
-        else:
-            if self.y_true.ndim != 3:
-                raise ValueError('Expected dimensions for y_true (3D data) is 3.'
-                                 'Requires format is: (z, x, y)'
-                                 'Got ndim: {}'.format(self.y_true.ndim))
+        for i, j in zip(correct_index[0], correct_index[1]):
+            self._add_detection(true_index=int(i), pred_index=int(j))
 
-        # Check if either frame is empty before proceeding
-        if self.n_true == 0:
-            logging.info('Ground truth frame is empty')
-            self.gained_detections += self.n_pred
-            self.empty_frame = 'n_true'
-        elif self.n_pred == 0:
-            logging.info('Prediction frame is empty')
-            self.missed_detections += self.n_true
-            self.empty_frame = 'n_pred'
-        elif test is False:
-            self.empty_frame = False
-            self._calc_iou()
-            self._modify_iou(force_event_links)
-            self._make_matrix()
-            self._linear_assignment()
-
-            # Check if there are loners before proceeding
-            if (self.loners_pred.shape[0] == 0) & (self.loners_true.shape[0] == 0):
-                pass
-            else:
-                self._assign_loners()
-                self._array_to_graph()
-                self._classify_graph()
-        else:
-            self.empty_frame = False
+        # Calc seg score for true positives if requested
+        iou_mask = np.where(self.seg_thresh == 0, self.iou, np.nan)
+
+        with warnings.catch_warnings():
+            warnings.simplefilter('ignore', category=RuntimeWarning)
+            # correct_index may be empty, suppress mean of empty slice warning
+            self.seg_score = np.nanmean(iou_mask[correct_index])
+
+        # Classify other errors using a graph
+        G = self._array_to_graph(matrix)
+        self._classify_graph(G)
+
+        # Calculate pixel-level stats
+        self.pixel_stats = PixelMetrics(y_true, y_pred)
+
+    def _add_detection(self, true_index=None, pred_index=None):
+        detection = Detection(true_index=true_index, pred_index=pred_index)
+
+        self._detections.add(detection)
+
+        # keep track of all error types
+        # TODO: better way to do this?
+        if detection.is_correct:
+            self._correct.add(detection)
+        if detection.is_gained:
+            self._gained.add(detection)
+        if detection.is_missed:
+            self._missed.add(detection)
+        if detection.is_split:
+            self._splits.add(detection)
+        if detection.is_merge:
+            self._merges.add(detection)
+        if detection.is_catastrophe:
+            self._catastrophes.add(detection)
 
     def _calc_iou(self):
         """Calculates IoU matrix for each pairwise comparison between true and
         predicted. Additionally, if seg is True, records a 1 for each pair of
         objects where $|Tbigcap P| > 0.5 * |T|$
         """
-
-        def get_box_labels(images):
-            props = regionprops(np.squeeze(images))
-            boxes, labels = [], []
-            for prop in props:
-                boxes.append(np.array(prop.bbox))
-                labels.append(int(prop.label))
-            boxes = np.array(boxes).astype('double')
-
-            return boxes, labels
-
-        self.iou = np.zeros((self.n_true, self.n_pred))
-
-        if self.seg:
-            self.seg_thresh = np.zeros((self.n_true, self.n_pred))
-
         # Use bounding boxes to find masks that are likely to overlap
-        y_true_boxes, y_true_labels = get_box_labels(self.y_true.astype('int'))
-        y_pred_boxes, y_pred_labels = get_box_labels(self.y_pred.astype('int'))
+        y_true_boxes, y_true_labels = get_box_labels(self.y_true)
+        y_pred_boxes, y_pred_labels = get_box_labels(self.y_pred)
+
+        if not y_true_boxes.shape[0] or not y_pred_boxes.shape[0]:
+            return  # cannot compute overlaps of nothing
 
         # has the form [gt_bbox, res_bbox]
-        if self.is_3d:
-            overlaps = compute_overlap_3D(y_true_boxes, y_pred_boxes)
-        else:
-            overlaps = compute_overlap(y_true_boxes, y_pred_boxes)
+        overlaps = self.compute_overlap(y_true_boxes, y_pred_boxes)
 
-        # Find the bboxes that have overlap at all
+        # Find the bboxes that have any overlap
         # (ind_ corresponds to box number - starting at 0)
         ind_true, ind_pred = np.nonzero(overlaps)
 
+        # TODO: this accounts for ~50+% of the time spent on calc_iou
         for index in range(ind_true.shape[0]):
-
             iou_y_true_idx = y_true_labels[ind_true[index]]
             iou_y_pred_idx = y_pred_labels[ind_pred[index]]
-            intersection = np.logical_and(self.y_true == iou_y_true_idx,
-                                          self.y_pred == iou_y_pred_idx)
-            union = np.logical_or(self.y_true == iou_y_true_idx,
-                                  self.y_pred == iou_y_pred_idx)
+
+            is_true = self.y_true == iou_y_true_idx
+            is_pred = self.y_pred == iou_y_pred_idx
+
+            intersection = np.count_nonzero(np.logical_and(is_true, is_pred))
+            union = np.count_nonzero(np.logical_or(is_true, is_pred))
+
+            iou = intersection / union
+
             # Subtract 1 from index to account for skipping 0
-            self.iou[iou_y_true_idx - 1, iou_y_pred_idx - 1] = intersection.sum() / union.sum()
+            self.iou[iou_y_true_idx - 1, iou_y_pred_idx - 1] = iou
 
-            if (self.seg) & \
-               (intersection.sum() > 0.5 * np.sum(self.y_true == index)):
+            if intersection > 0.5 * np.count_nonzero(self.y_true == index):
                 self.seg_thresh[iou_y_true_idx - 1, iou_y_pred_idx - 1] = 1
 
-    def _modify_iou(self, force_event_links):
-        """Modifies the IOU matrix to boost the value for small cells.
+    def _get_modified_iou(self, force_event_links):
+        """Modifies the IoU matrix to boost the value for small cells.
 
         Args:
-            force_event_links (:obj:`bool'): flag that determines whether to modify IOU values of
-             large cells if a small cell has been split or merged with them.
-        """
+            force_event_links (:obj:`bool'): Whether to modify IOU values of
+                large objects if they have been split or merged by
+                a small object.
 
+        Returns:
+            np.array: The modified IoU matrix.
+        """
         # identify cells that have matches in IOU but may be too small
-        true_labels, pred_labels = np.where(np.logical_and(self.iou > 0,
-                                                           self.iou < (1 - self.cutoff1)))
+        true_labels, pred_labels = np.nonzero(
+            np.logical_and(self.iou > 0, self.iou < 1 - self.cutoff1)
+        )
 
-        self.iou_modified = self.iou.copy()
+        iou_modified = self.iou.copy()
 
         for idx in range(len(true_labels)):
             # add 1 to get back to original label id
-            true_label, pred_label = true_labels[idx] + 1, pred_labels[idx] + 1
+            true_idx, pred_idx = true_labels[idx], pred_labels[idx]
+            true_label, pred_label = true_idx + 1, pred_idx + 1
             true_mask = self.y_true == true_label
             pred_mask = self.y_pred == pred_label
 
             # fraction of true cell that is contained within pred cell, vice versa
-            true_in_pred = np.sum(self.y_true[pred_mask] == true_label) / np.sum(true_mask)
-            pred_in_true = np.sum(self.y_pred[true_mask] == pred_label) / np.sum(pred_mask)
+            true_in_pred = np.count_nonzero(
+                self.y_true[pred_mask] == true_label) / np.sum(true_mask)
+            pred_in_true = np.count_nonzero(
+                self.y_pred[true_mask] == pred_label) / np.sum(pred_mask)
 
-            iou_val = self.iou[true_label - 1, pred_label - 1]
+            iou_val = self.iou[true_idx, pred_idx]
             max_val = np.max([true_in_pred, pred_in_true])
 
             # if this cell has a small IOU due to its small size,
             # but is at least half contained within the big cell,
             # we bump its IOU value up so it doesn't get dropped from the graph
             if iou_val <= self.cutoff1 and max_val > 0.5:
-                self.iou_modified[true_label - 1, pred_label - 1] = self.cutoff2
-
-                # optionally, we can also decrease the IOU value of the cell that
-                # swallowed up the small cell so that it doesn't directly match a different cell
-                if force_event_links:
-                    if true_in_pred > 0.5:
-                        fix_idx = np.where(self.iou[:, pred_label - 1] >= 1 - self.cutoff1)
-                        self.iou_modified[fix_idx, pred_label - 1] = 1 - self.cutoff1 - 0.01
-                    elif pred_in_true > 0.5:
-                        fix_idx = np.where(self.iou[true_label - 1, :] >= 1 - self.cutoff1)
-                        self.iou_modified[true_label - 1, fix_idx] = 1 - self.cutoff1 - 0.01
-
-    def _make_matrix(self):
+                iou_modified[true_idx, pred_idx] = self.cutoff2
+
+                # optionally, we can also decrease the IOU value of the cell
+                # that swallowed up the small cell so that it doesn't directly
+                # match a different cell
+                if force_event_links and true_in_pred > 0.5:
+                    fix_idx = np.nonzero(self.iou[:, pred_idx] >= 1 - self.cutoff1)
+                    iou_modified[fix_idx, pred_idx] = 1 - self.cutoff1 - 0.01
+
+                if force_event_links and pred_in_true > 0.5:
+                    fix_idx = np.nonzero(self.iou[true_idx, :] >= 1 - self.cutoff1)
+                    iou_modified[true_idx, fix_idx] = 1 - self.cutoff1 - 0.01
+
+        return iou_modified
+
+    def _get_cost_matrix(self):
         """Assembles cost matrix using the iou matrix and cutoff1
 
         The previously calculated iou matrix is cast into the top left and
@@ -383,109 +559,81 @@ def _make_matrix(self):
         value of cutoff1 the more likely it is for the linear sum assignment
         to pick unmatched assignments for objects.
         """
-
-        self.cm = np.ones((self.n_obj, self.n_obj))
+        n_obj = self.n_true + self.n_pred
+        matrix = np.ones((n_obj, n_obj))
 
         # Assign 1 - iou to top left and bottom right
-        self.cm[:self.n_true, :self.n_pred] = 1 - self.iou_modified
-        self.cm[-self.n_pred:, -self.n_true:] = 1 - self.iou_modified.T
+        cost = 1 - self.iou_modified
+        matrix[:self.n_true, :self.n_pred] = cost
+        matrix[n_obj - self.n_pred:, n_obj - self.n_true:] = cost.T
 
         # Calculate diagonal corners
-        bl = self.cutoff1 * \
-            np.eye(self.n_pred) + np.ones((self.n_pred, self.n_pred)) - \
-            np.eye(self.n_pred)
-        tr = self.cutoff1 * \
-            np.eye(self.n_true) + np.ones((self.n_true, self.n_true)) - \
-            np.eye(self.n_true)
+        bl = (self.cutoff1 * np.eye(self.n_pred)
+              + np.ones((self.n_pred, self.n_pred))
+              - np.eye(self.n_pred))
+        tr = (self.cutoff1 * np.eye(self.n_true)
+              + np.ones((self.n_true, self.n_true))
+              - np.eye(self.n_true))
 
         # Assign diagonals to cm
-        self.cm[-self.n_pred:, :self.n_pred] = bl
-        self.cm[:self.n_true, -self.n_true:] = tr
+        matrix[n_obj - self.n_pred:, :self.n_pred] = bl
+        matrix[:self.n_true, n_obj - self.n_true:] = tr
+        return matrix
 
     def _linear_assignment(self):
-        """Runs linear sun assignment on cost matrix, identifies true positives
-        and unassigned true and predicted cells.
+        """Runs linear sun assignment on cost matrix, identifies true
+        positives and unassigned true and predicted cells.
 
         True positives correspond to assignments in the top left or bottom
         right corner. There are two possible unassigned positions: true cell
         unassigned in bottom left or predicted cell unassigned in top right.
         """
+        cost_matrix = self._get_cost_matrix()
 
-        self.results = linear_sum_assignment(self.cm)
+        results = linear_sum_assignment(cost_matrix)
 
         # Map results onto cost matrix
-        self.cm_res = np.zeros(self.cm.shape)
-        self.cm_res[self.results[0], self.results[1]] = 1
+        assignment_matrix = np.zeros_like(cost_matrix)
+        assignment_matrix[results] = 1
+        return assignment_matrix
 
-        # Identify direct matches as true positives
-        correct_index = np.where(self.cm_res[:self.n_true, :self.n_pred] == 1)
-        self.correct_detections += len(correct_index[0])
-        self.correct_indices['y_true'].extend(list(correct_index[0] + 1))
-        self.correct_indices['y_pred'].extend(list(correct_index[1] + 1))
-
-        # Calc seg score for true positives if requested
-        if self.seg is True:
-            iou_mask = self.iou.copy()
-            iou_mask[self.seg_thresh == 0] = np.nan
-            self.seg_score = np.nanmean(iou_mask[correct_index[0], correct_index[1]])
-
-        # Collect unassigned cells
-        self.loners_pred, _ = np.where(
-            self.cm_res[-self.n_pred:, :self.n_pred] == 1)
-        self.loners_true, _ = np.where(
-            self.cm_res[:self.n_true, -self.n_true:] == 1)
-
-    def _assign_loners(self):
-        """Generate an iou matrix for the subset unassigned cells
-        """
-
-        self.n_pred2 = len(self.loners_pred)
-        self.n_true2 = len(self.loners_true)
-        self.n_obj2 = self.n_pred2 + self.n_true2
-
-        self.cost_l = np.zeros((self.n_true2, self.n_pred2))
-
-        for i, t in enumerate(self.loners_true):
-            for j, p in enumerate(self.loners_pred):
-                self.cost_l[i, j] = self.iou_modified[t, p]
-
-        self.cost_l_bin = self.cost_l >= self.cutoff2
-
-    def _array_to_graph(self):
+    def _array_to_graph(self, matrix):
         """Transform matrix for unassigned cells into a graph object
 
         In order to cast the iou matrix into a graph form, we treat each
         unassigned cell as a node. The iou values for each pair of cells is
         treated as an edge between nodes/cells. Any iou values equal to 0 are
         dropped because they indicate no overlap between cells.
-        """
 
-        # Use meshgrid to get true and predicted cell index for each val
-        tt, pp = np.meshgrid(self.loners_true, self.loners_pred, indexing='ij')
+        Args:
+            matrix (np.array): Assignment matrix.
+        """
+        # Collect unassigned objects
+        x, y = matrix.shape
+        gained, _ = np.nonzero(matrix[x - self.n_pred:, :self.n_pred])
+        missed, _ = np.nonzero(matrix[:self.n_true, y - self.n_true:])
 
-        df = pd.DataFrame({
-            'true': tt.flatten(),
-            'pred': pp.flatten(),
-            'weight': self.cost_l_bin.flatten()
-        })
+        # Use meshgrid to get true and predicted object index for each val
+        tt, pp = np.meshgrid(missed, gained, indexing='ij')
 
-        # Change cell index to str names
-        df['true'] = 'true_' + df['true'].astype('str')
-        df['pred'] = 'pred_' + df['pred'].astype('str')
+        true_nodes = tt.flatten()
+        pred_nodes = pp.flatten()
 
-        # Drop 0 weights to only retain overlapping cells
-        dfedge = df.drop(df[df['weight'] == 0].index)
+        # construct list of edges for networkx
+        G = nx.Graph()
 
-        # Create graph from edges
-        self.G = nx.from_pandas_edgelist(dfedge, source='true', target='pred')
+        for t, p in zip(true_nodes, pred_nodes):
+            # edges between overlapping objects only
+            if self.iou_modified[t, p] >= self.cutoff2:
+                G.add_edge('true_{}'.format(t), 'pred_{}'.format(p))
 
         # Add nodes to ensure all cells are included
-        nodes_true = ['true_' + str(node) for node in self.loners_true]
-        nodes_pred = ['pred_' + str(node) for node in self.loners_pred]
-        nodes = nodes_true + nodes_pred
-        self.G.add_nodes_from(nodes)
+        G.add_nodes_from(('true_{}'.format(n) for n in missed))
+        G.add_nodes_from(('pred_{}'.format(n) for n in gained))
+
+        return G
 
-    def _classify_graph(self):
+    def _classify_graph(self, G):
         """Assign each node in graph to an error type
 
         Nodes with a degree (connectivity) of 0 correspond to either false
@@ -497,121 +645,64 @@ def _classify_graph(self):
         error. If the top level node is a predicted cell, this indicates a merge
         event. If the top level node is a true cell, this indicates a split event.
         """
-
         # Find subgraphs, e.g. merge/split
-        for g in (self.G.subgraph(c) for c in nx.connected_components(self.G)):
+        for g in (G.subgraph(c) for c in nx.connected_components(G)):
             # Get the highest degree node
-            k = max(dict(g.degree).items(), key=operator.itemgetter(1))[0]
-
-            # Map index back to original cost matrix, adjust for 1-based indexing in labels
-            index = int(k.split('_')[-1]) + 1
-            # Process degree 0 nodes
-            if g.degree[k] == 0:
-                if 'pred' in k:
-                    self.gained_detections += 1
-                    self.gained_indices['y_pred'].append(index)
-                if 'true' in k:
-                    self.missed_detections += 1
-                    self.missed_indices['y_true'].append(index)
-
-            # Process degree 1 nodes
-            if g.degree[k] == 1:
-                for node in g.nodes:
-                    node_index = int(node.split('_')[-1]) + 1
-                    if 'pred' in node:
-                        self.gained_detections += 1
-                        self.gained_indices['y_pred'].append(node_index)
-                    if 'true' in node:
-                        self.missed_detections += 1
-                        self.missed_indices['y_true'].append(node_index)
-
-            # Process multi-degree nodes
-            elif g.degree[k] > 1:
-                node_type = k.split('_')[0]
-                nodes = g.nodes()
-                # Check whether the subgraph has multiple types of the
-                # highest degree node (true or pred)
-                n_node_type = np.sum([node_type in node for node in nodes])
-                # If there is only one of the high degree node type in the
-                # sub graph, then we have either a merge or a split
-                if n_node_type == 1:
-                    # Check for merges
-                    if 'pred' in node_type:
-                        self.merge += 1
-                        self.missed_det_from_merge += len(nodes) - 2
-                        true_merge_indices = [int(node.split('_')[-1]) + 1
-                                              for node in nodes if 'true' in node]
-                        self.merge_indices['y_true'] += true_merge_indices
-                        self.merge_indices['y_pred'].append(index)
-                    # Check for splits
-                    elif 'true' in node_type:
-                        self.split += 1
-                        self.gained_det_from_split += len(nodes) - 2
-                        self.split_indices['y_true'].append(index)
-                        pred_split_indices = [int(node.split('_')[-1]) + 1
-                                              for node in nodes if 'pred' in node]
-                        self.split_indices['y_pred'] += pred_split_indices
-
-                # If there are multiple types of the high degree node,
-                # then we have a catastrophe
-                else:
-                    self.catastrophe += 1
-                    true_indices = [int(node.split('_')[-1]) + 1
-                                    for node in nodes if 'true' in node]
-                    pred_indices = [int(node.split('_')[-1]) + 1
-                                    for node in nodes if 'pred' in node]
-
-                    self.true_det_in_catastrophe = len(true_indices)
-                    self.pred_det_in_catastrophe = len(pred_indices)
-
-                    self.catastrophe_indices['y_true'] += true_indices
-                    self.catastrophe_indices['y_pred'] += pred_indices
-
-            # Save information about the cells involved in the different error types
-            gained_label_image = np.zeros_like(self.y_pred)
-            for l in self.gained_indices['y_pred']:
-                gained_label_image[self.y_pred == l] = l
-            self.gained_props = regionprops(gained_label_image)
-
-            missed_label_image = np.zeros_like(self.y_true)
-            for l in self.missed_indices['y_true']:
-                missed_label_image[self.y_true == l] = l
-            self.missed_props = regionprops(missed_label_image)
-
-            merge_label_image = np.zeros_like(self.y_true)
-            for l in self.merge_indices['y_true']:
-                merge_label_image[self.y_true == l] = l
-            self.merge_props = regionprops(merge_label_image)
-
-            split_label_image = np.zeros_like(self.y_true)
-            for l in self.split_indices['y_true']:
-                split_label_image[self.y_true == l] = l
-            self.split_props = regionprops(split_label_image)
-
-    def print_report(self):
-        """Print report of error types and frequency
-        """
-        print(self.save_to_dataframe())
-
-    def save_to_dataframe(self):
-        """Save error results to a pandas dataframe
-
-        Returns:
-            pandas.DataFrame: Single row dataframe with error types as columns
-        """
-        try:
-            self.precision = self.correct_detections / self.n_pred
-        except ZeroDivisionError:
-            self.precision = 0
+            _, max_d = max(dict(g.degree).items(), key=operator.itemgetter(1))
+
+            true_indices, pred_indices = [], []
+
+            for node in g.nodes:
+                node_type, index = node.split('_')
+                index = int(index) + 1
+
+                if node_type == 'true':
+                    if max_d > 1:
+                        true_indices.append(index)
+                    else:
+                        self._add_detection(true_index=index)
+
+                if node_type == 'pred':
+                    if max_d > 1:
+                        pred_indices.append(index)
+                    else:
+                        self._add_detection(pred_index=index)
+
+            self._add_detection(
+                true_index=tuple(true_indices) if true_indices else None,
+                pred_index=tuple(pred_indices) if pred_indices else None,
+            )
+
+    def _get_props(self, detection_type):
+        prediction_types = {
+            'gained',
+        }
+        is_pred_type = detection_type in prediction_types
+        arr = self.y_pred if is_pred_type else self.y_true
+        label_image = np.zeros_like(arr)
+        attrname = '_{}'.format(detection_type)
 
         try:
-            self.recall = self.correct_detections / self.n_true
-        except ZeroDivisionError:
-            self.recall = 0
-
-        self.f1 = hmean([self.recall, self.precision])
-
-        D = {
+            detections = getattr(self, attrname)
+        except AttributeError:
+            raise ValueError('Invalid detection_type: {}'.format(
+                detection_type))
+
+        for det in detections:
+            idx = det.pred_index if is_pred_type else det.true_index
+            idx = idx if isinstance(idx, tuple) else (idx,)
+            for i in idx:
+                label_image[arr == i] = i
+
+        return regionprops(label_image)
+
+    def __repr__(self):
+        """Format the calculated statistics as a ``pd.DataFrame``."""
+        return json.dumps(self.to_dict())
+
+    def to_dict(self):
+        """Return a dictionary representation of the calclulated metrics."""
+        return {
             'n_pred': self.n_pred,
             'n_true': self.n_true,
             'correct_detections': self.correct_detections,
@@ -621,63 +712,162 @@ def save_to_dataframe(self):
             'gained_det_from_split': self.gained_det_from_split,
             'true_det_in_catastrophe': self.true_det_in_catastrophe,
             'pred_det_in_catastrophe': self.pred_det_in_catastrophe,
-            'merge': self.merge,
-            'split': self.split,
-            'catastrophe': self.catastrophe,
+            'merge': self.merges,
+            'split': self.splits,
+            'catastrophe': self.catastrophes,
             'precision': self.precision,
             'recall': self.recall,
-            'f1': self.f1
+            'f1': self.f1,
+            'seg': self.seg_score,
+            'jaccard': self.jaccard,
+            'dice': self.dice,
         }
 
-        if self.seg is True:
-            D['seg'] = self.seg_score
+    @property
+    def correct_detections(self):
+        return len(self._correct)
+
+    @property
+    def missed_detections(self):
+        return len(self._missed)
+
+    @property
+    def gained_detections(self):
+        return len(self._gained)
+
+    @property
+    def splits(self):
+        return len(self._splits)
+
+    @property
+    def merges(self):
+        return len(self._merges)
+
+    @property
+    def catastrophes(self):
+        return len(self._catastrophes)
+
+    @property
+    def gained_det_from_split(self):
+        gained_dets = 0
+        for det in self._splits:
+            true_idx = _cast_to_tuple(det.true_index)
+            pred_idx = _cast_to_tuple(det.pred_index)
+            gained_dets += len(true_idx) + len(pred_idx) - 2
+        return gained_dets
+
+    @property
+    def missed_det_from_merge(self):
+        missed_dets = 0
+        for det in self._merges:
+            true_idx = _cast_to_tuple(det.true_index)
+            pred_idx = _cast_to_tuple(det.pred_index)
+            missed_dets += len(true_idx) + len(pred_idx) - 2
+        return missed_dets
+
+    @property
+    def true_det_in_catastrophe(self):
+        return sum([len(d.true_index) for d in self._catastrophes])
+
+    @property
+    def pred_det_in_catastrophe(self):
+        return sum([len(d.pred_index) for d in self._catastrophes])
+
+    @property
+    def split_props(self):
+        return self._get_props('splits')
+
+    @property
+    def merge_props(self):
+        return self._get_props('merges')
+
+    @property
+    def missed_props(self):
+        return self._get_props('missed')
+
+    @property
+    def gained_props(self):
+        return self._get_props('gained')
+
+    @property
+    def recall(self):
+        try:
+            recall = self.correct_detections / self.n_true
+        except ZeroDivisionError:
+            recall = 0
+        return recall
 
-        # Calculate jaccard index for pixel classification
-        pixel_stats = stats_pixelbased(self.y_true != 0, self.y_pred != 0)
-        D['jaccard'] = pixel_stats['jaccard']
+    @property
+    def precision(self):
+        try:
+            precision = self.correct_detections / self.n_pred
+        except ZeroDivisionError:
+            precision = 0
+        return precision
 
-        df = pd.DataFrame(D, index=[0], dtype='float64')
+    @property
+    def f1(self):
+        return hmean([self.recall, self.precision])
 
-        # Change appropriate columns to int dtype
-        col = ['n_pred', 'n_true', 'correct_detections', 'missed_detections', 'gained_detections',
-               'missed_det_from_merge', 'gained_det_from_split', 'true_det_in_catastrophe',
-               'pred_det_in_catastrophe', 'merge', 'split', 'catastrophe']
-        df[col] = df[col].astype('int')
+    @property
+    def jaccard(self):
+        return self.pixel_stats.jaccard
 
-        return df
+    @property
+    def dice(self):
+        return self.pixel_stats.jaccard
 
-    def save_error_ids(self):
-        """Saves the ids of cells in each error category for subsequent visualization
+    def plot_errors(self):
+        """Plots the errors identified from linear assignment code.
 
-        Returns:
-            error_dict: dictionary containing {category_name: id list} pairs
+        This must be run with sequentially relabeled data.
+
+        TODO: this is not working!
         """
 
-        error_dict = {'splits': self.split_indices,
-                      'merges': self.merge_indices,
-                      'gains': self.gained_indices,
-                      'misses': self.missed_indices,
-                      'catastrophes': self.catastrophe_indices,
-                      'correct': self.correct_indices}
+        import matplotlib as mpl
+        import matplotlib.pyplot as plt
 
-        return error_dict, self.y_true, self.y_pred
+        # erode edges for easier visualization of adjacent cells
+        y_true = erode_edges(self.y_true.copy(), 1)
+        y_pred = erode_edges(self.y_pred.copy(), 1)
 
+        # semantic labels for each error
+        categories = ['Background', 'missed', 'splits', 'merges',
+                      'gained', 'catastrophes', 'correct']
 
-def to_precision(x, p):
-    """
-    returns a string representation of x formatted with a precision of p
+        # Background is set to zero
+        plotting_tif = np.zeros_like(y_true)
 
-    Based on the webkit javascript implementation taken from here:
-    https://code.google.com/p/webkit-mirror/source/browse/JavaScriptCore/kjs/number_object.cpp
-    """
-    decimal.getcontext().prec = p
-    dec = decimal.Decimal(x)
+        # missed detections are tracked with true labels
+        misses = [d.true_index for d in self._missed]
+        plotting_tif[np.isin(y_true, misses)] = 1
+
+        # skip background and misses, already done
+        for i, category in enumerate(categories[2:]):
+            # the rest are all on y_pred
+            labels = list(getattr(self, '_{}'.format(category)))
+            plotting_tif[np.isin(y_pred, labels)] = i + 2
 
-    return round(float(dec), p)
+        plotting_colors = ['Black', 'Pink', 'Blue', 'Green',
+                           'tan', 'Red', 'Grey']
+
+        cmap = mpl.colors.ListedColormap(plotting_colors)
+
+        fig, ax = plt.subplots(nrows=1, ncols=1)
+        mat = ax.imshow(plotting_tif, cmap=cmap,
+                        vmin=np.min(plotting_tif) - .5,
+                        vmax=np.max(plotting_tif) + .5)
+
+        # tell the colorbar to tick at integers
+        ticks = np.arange(len(categories))
+        cbar = fig.colorbar(mat, ticks=ticks)
+        cbar.ax.set_yticklabels(categories)
+        fig.tight_layout()
 
 
 class Metrics(object):
-    """Class to calculate and save various classification metrics
+    """Class to calculate and save various segmentation metrics.
 
     Args:
         model_name (str): Name of the model which determines output file names
@@ -693,8 +883,6 @@ class Metrics(object):
         feature_key (:obj:`list`, optional): List of strings, feature names
         json_notes (:obj:`str`, optional): Str providing any additional
             information about the model
-        seg (:obj:`bool`, optional): Calculates SEG score for
-            cell tracking competition
         force_event_links(:obj:`bool`, optional): Flag that determines whether to modify IOU
             calculation so that merge or split events with cells of very different sizes are
             never misclassified as misses/gains.
@@ -704,16 +892,9 @@ class Metrics(object):
     Examples:
         >>> from deepcell import metrics
         >>> m = metrics.Metrics('model_name')
-        >>> m.run_all(
-                y_true_lbl,
-                y_pred_lbl,
-                y_true_unlbl,
-                y_true_unlbl)
-        >>> m.all_pixel_stats(y_true_unlbl,y_pred_unlbl)
-        >>> m.calc_obj_stats(y_true_lbl,y_pred_lbl)
-        >>> m.save_to_json(m.output)
+        >>> all_metrics = m.run_all(y_true, y_pred)
+        >>> m.save_to_json(all_metrics)
     """
-
     def __init__(self, model_name,
                  outdir='',
                  cutoff1=0.4,
@@ -721,12 +902,11 @@ def __init__(self, model_name,
                  pixel_threshold=0.5,
                  ndigits=4,
                  crop_size=None,
-                 return_iou=False,
                  feature_key=[],
                  json_notes='',
-                 seg=False,
                  force_event_links=False,
-                 is_3d=False):
+                 is_3d=False,
+                 **kwargs):
         self.model_name = model_name
         self.outdir = outdir
         self.cutoff1 = cutoff1
@@ -734,71 +914,25 @@ def __init__(self, model_name,
         self.pixel_threshold = pixel_threshold
         self.ndigits = ndigits
         self.crop_size = crop_size
-        self.return_iou = return_iou
         self.feature_key = feature_key
         self.json_notes = json_notes
-        self.seg = seg
         self.force_event_links = force_event_links
         self.is_3d = is_3d
 
-        # Initialize output list to collect stats
-        self.output = []
-
-    def all_pixel_stats(self, y_true, y_pred):
-        """Collect pixel statistics for each feature.
-
-        y_true should have the appropriate transform applied to match y_pred.
-        Each channel is converted to binary using the threshold
-        'pixel_threshold' prior to calculation of accuracy metrics.
-
-        Args:
-            y_true (numpy.array): Ground truth annotations after transform
-            y_pred (numpy.array): Model predictions without labeling
-
-        Raises:
-            ValueError: If y_true and y_pred are not the same shape
-        """
-
-        if y_pred.shape != y_true.shape:
-            raise ValueError('Input shapes need to match. Shape of prediction '
-                             'is: {}.  Shape of y_true is: {}'.format(
-                                 y_pred.shape, y_true.shape))
-
-        n_features = y_pred.shape[-1]
-
-        # Intialize df to collect pixel stats
-        self.pixel_df = pd.DataFrame()
+        if 'seg' in kwargs:
+            warnings.warn('seg is deprecated and will be removed '
+                          'in a future release', DeprecationWarning)
 
-        # Set numeric feature key if existing key is not write length
-        if n_features != len(self.feature_key):
-            self.feature_key = range(n_features)
-
-        for i, k in enumerate(self.feature_key):
-            yt = y_true[:, :, :, i] > self.pixel_threshold
-            yp = y_pred[:, :, :, i] > self.pixel_threshold
-            stats = stats_pixelbased(yt, yp)
-            self.pixel_df = self.pixel_df.append(
-                pd.DataFrame(stats, index=[k]))
-
-        # Save stats to output dictionary
-        self.output = self.output + self.pixel_df_to_dict(self.pixel_df)
-
-        # Calculate confusion matrix
-        self.cm = self.calc_pixel_confusion_matrix(y_true, y_pred)
-        self.output.append(dict(
-            name='confusion_matrix',
-            value=self.cm.tolist(),
-            feature='all',
-            stat_type='pixel'
-        ))
-
-        self.print_pixel_report()
+        # Initialize output list to collect stats
+        self.object_metrics = []
+        self.pixel_metrics = []
 
-    def pixel_df_to_dict(self, df):
+    def df_to_dict(self, df, stat_type='pixel'):
         """Output pandas df as a list of dictionary objects
 
         Args:
             df (pandas.DataFrame): Dataframe of statistics for each channel
+            stat_type (str): Category of statistic.
 
         Returns:
             list: List of dictionaries
@@ -813,7 +947,7 @@ def pixel_df_to_dict(self, df):
                 name=k,
                 value=v,
                 feature='average',
-                stat_type='pixel'
+                stat_type=stat_type,
             ))
 
         # Save individual stats to list
@@ -823,40 +957,77 @@ def pixel_df_to_dict(self, df):
                     name=k,
                     value=v,
                     feature=i,
-                    stat_type='pixel'
+                    stat_type=stat_type,
                 ))
 
         return L
 
-    def calc_pixel_confusion_matrix(self, y_true, y_pred):
-        """Calculate confusion matrix for pixel classification data.
+    def calc_pixel_stats(self, y_true, y_pred, axis=-1):
+        """Calculate pixel statistics for each feature.
+
+        ``y_true`` should have the appropriate transform applied to match
+        ``y_pred``. Each channel is converted to binary using the threshold
+        ``pixel_threshold`` prior to calculation of accuracy metrics.
 
         Args:
-            y_true (numpy.array): Ground truth annotations after any
-                necessary transformations
-            y_pred (numpy.array): Prediction array
+            y_true (numpy.array): Ground truth annotations after transform
+            y_pred (numpy.array): Model predictions without labeling
 
         Returns:
-            numpy.array: nxn confusion matrix determined by number of features.
+            list: list of dictionaries with each stat being a key.
+
+        Raises:
+            ValueError: If y_true and y_pred are not the same shape
         """
+        n_features = y_pred.shape[axis]
 
-        # Argmax collapses on feature dimension to assign class to each pixel
-        # Flatten is requiremed for confusion matrix
-        y_true = y_true.argmax(axis=-1).flatten()
-        y_pred = y_pred.argmax(axis=-1).flatten()
+        pixel_metrics = []
 
-        return confusion_matrix(y_true, y_pred)
+        slc = [slice(None)] * y_pred.ndim
+        for i in range(n_features):
+            slc[axis] = slice(i, i + 1)
+            yt = y_true[slc] > self.pixel_threshold
+            yp = y_pred[slc] > self.pixel_threshold
+            pm = PixelMetrics(yt, yp)
+            pixel_metrics.append(pm.to_dict())
 
-    def print_pixel_report(self):
-        """Print report of pixel based statistics
-        """
+        pixel_df = pd.DataFrame.from_records(pixel_metrics)
+
+        # Calculate confusion matrix
+        cm = PixelMetrics.get_confusion_matrix(y_true, y_pred, axis=axis)
 
         print('\n____________Pixel-based statistics____________\n')
-        print(self.pixel_df)
+        print(pixel_df)
         print('\nConfusion Matrix')
-        print(self.cm)
+        print(cm)
+
+        output = self.df_to_dict(pixel_df)
+
+        output.append(dict(
+            name='confusion_matrix',
+            value=cm.tolist(),
+            feature='all',
+            stat_type='pixel'
+        ))
+        return output
+
+    def calc_pixel_confusion_matrix(self, y_true, y_pred, axis=-1):
+        """DEPRECATED: Use ``PixelMetrics.get_confusion_matrix``.
+
+        Calculate confusion matrix for pixel classification data.
+
+        Args:
+            y_true (numpy.array): Ground truth annotations after any
+                necessary transformations
+            y_pred (numpy.array): Prediction array
+            axis (int): The channel axis of the input arrays.
+
+        Returns:
+            numpy.array: nxn confusion matrix determined by number of features.
+        """
+        return PixelMetrics.get_confusion_matrix(y_true, y_pred, axis=axis)
 
-    def calc_object_stats(self, y_true, y_pred):
+    def calc_object_stats(self, y_true, y_pred, progbar=True):
         """Calculate object statistics and save to output
 
         Loops over each frame in the zeroth dimension, which should pass in
@@ -866,13 +1037,16 @@ def calc_object_stats(self, y_true, y_pred):
         Args:
             y_true (numpy.array): Labeled ground truth annotations
             y_pred (numpy.array): Labeled prediction mask
+            progbar (bool): Whether to show the progress tqdm progress bar
+
+        Returns:
+            list: list of dictionaries with each stat being a key.
 
         Raises:
             ValueError: If y_true and y_pred are not the same shape
             ValueError: If data_type is 2D, if input shape does not have ndim 3 or 4
             ValueError: If data_type is 3D, if input shape does not have ndim 4
         """
-
         if y_pred.shape != y_true.shape:
             raise ValueError('Input shapes need to match. Shape of prediction '
                              'is: {}.  Shape of y_true is: {}'.format(
@@ -890,143 +1064,169 @@ def calc_object_stats(self, y_true, y_pred):
         # TODO - add compatibility for multi-channel 3D-data
         else:
             if y_true.ndim != 4:
-                raise ValueError('Expected dimensions for y_true (3D data) is 4.'
-                                 'Required format is: (batch, z, x, y)'
+                raise ValueError('Expected dimensions for y_true (3D data) is 4. '
+                                 'Required format is: (batch, z, x, y) '
                                  'Got ndim: {}'.format(y_true.ndim))
 
-        self.stats = pd.DataFrame()
-        self.predictions = []
-
-        # boolean so that warning only gets displayed once
-        relabel_flag = False
-        for i in range(y_true.shape[0]):
+        all_object_metrics = []  # store all calculated metrics
+        is_batch_relabeled = False  # used to warn if batches were relabeled
 
+        for i in tqdm(range(y_true.shape[0]), disable=not progbar):
             # check if labels aren't sequential, raise warning on first occurence if so
             true_batch, pred_batch = y_true[i], y_pred[i]
             true_batch_relabel, _, _ = relabel_sequential(true_batch)
             pred_batch_relabel, _, _ = relabel_sequential(pred_batch)
 
-            if not (np.array_equal(true_batch, true_batch_relabel) and
-                    np.array_equal(pred_batch, pred_batch_relabel)):
-
-                    # segmentations were relabeled
-                    relabel_flag = True
-            o = ObjectAccuracy(true_batch_relabel,
-                               pred_batch_relabel,
-                               cutoff1=self.cutoff1,
-                               cutoff2=self.cutoff2,
-                               seg=self.seg,
-                               force_event_links=self.force_event_links,
-                               is_3d=self.is_3d)
-            self.stats = self.stats.append(o.save_to_dataframe())
-            predictions = o.save_error_ids()
-            self.predictions.append(predictions)
-            if i % 500 == 0:
-                logging.info('{} samples processed'.format(i))
-
-        if relabel_flag:
+            # check if segmentations were relabeled
+            if not is_batch_relabeled:  # only one True is required
+                is_batch_relabeled = not (
+                    np.array_equal(true_batch, true_batch_relabel)
+                    and np.array_equal(pred_batch, pred_batch_relabel)
+                )
+
+            o = ObjectMetrics(
+                true_batch_relabel,
+                pred_batch_relabel,
+                cutoff1=self.cutoff1,
+                cutoff2=self.cutoff2,
+                force_event_links=self.force_event_links,
+                is_3d=self.is_3d)
+
+            all_object_metrics.append(o)
+
+        if is_batch_relabeled:
             warnings.warn(
                 'Provided data is being relabeled. Cell ids from metrics will not match '
                 'cell ids in original data. Relabel your data prior to running the '
-                'metrics package if you wish to maintain cell ids')
+                'metrics package if you wish to maintain cell ids. ')
 
-        # Write out summed statistics
-        for k, v in self.stats.iteritems():
-            if k == 'seg':
-                self.output.append(dict(
-                    name=k,
-                    value=v.mean(),
-                    feature='mean',
-                    stat_type='object'
-                ))
-            else:
-                self.output.append(dict(
-                    name=k,
-                    value=v.sum().astype('float64'),
-                    feature='sum',
-                    stat_type='object'
-                ))
+        # print the object report
+        object_metrics = pd.DataFrame.from_records([
+            o.to_dict() for o in all_object_metrics
+        ])
+        self.print_object_report(object_metrics)
+        return object_metrics
+
+    def summarize_object_metrics_df(self, df):
+        correct_detections = int(df['correct_detections'].sum())
+        n_true = int(df['n_true'].sum())
+        n_pred = int(df['n_pred'].sum())
+
+        _round = lambda x: round(x, self.ndigits)
 
-        self.print_object_report()
+        seg = df['seg'].mean()
+        jaccard = df['jaccard'].mean()
 
-    def print_object_report(self):
+        try:
+            recall = correct_detections / n_true
+        except ZeroDivisionError:
+            recall = np.nan
+        try:
+            precision = correct_detections / n_pred
+        except ZeroDivisionError:
+            precision = 0
+
+        errors = [
+            'gained_detections',
+            'missed_detections',
+            'split',
+            'merge',
+            'catastrophe',
+        ]
+
+        bad_detections = [
+            'gained_det_from_split',
+            'missed_det_from_merge',
+            'true_det_in_catastrophe',
+            'pred_det_in_catastrophe',
+        ]
+
+        summary = {
+            'correct_detections': correct_detections,
+            'n_true': n_true,
+            'n_pred': n_pred,
+            'recall': _round(recall),
+            'precision': _round(precision * 100),
+            'seg': _round(seg * 100),
+            'jaccard': _round(jaccard),
+            'total_errors': 0,
+        }
+        # update bad detections
+        for k in bad_detections:
+            summary[k] = int(df[k].sum())
+        # update error counts
+        for k in errors:
+            count = int(df[k].sum())
+            summary[k] = count
+            summary['total_errors'] += count
+        return summary
+
+    def print_object_report(self, object_metrics):
         """Print neat report of object based statistics
+
+        Args:
+            object_metrics (pd.DataFrame): DataFrame of all calculated metrics
         """
+        summary = self.summarize_object_metrics_df(object_metrics)
+        errors = [
+            'gained_detections',
+            'missed_detections',
+            'split',
+            'merge',
+            'catastrophe'
+        ]
+
+        bad_detections = [
+            'gained_det_from_split',
+            'missed_det_from_merge',
+            'true_det_in_catastrophe',
+            'pred_det_in_catastrophe',
+        ]
 
         print('\n____________Object-based statistics____________\n')
-        print('Number of true cells:\t\t', self.stats['n_true'].sum())
-        print('Number of predicted cells:\t', self.stats['n_pred'].sum())
+        print('Number of true cells:\t\t', summary['n_true'])
+        print('Number of predicted cells:\t', summary['n_pred'])
 
         print('\nCorrect detections:  {}\tRecall: {}%'.format(
-            int(self.stats['correct_detections'].sum()),
-            to_precision(100 * self.stats['correct_detections'].sum() / self.stats['n_true'].sum(),
-                         self.ndigits)))
+            summary['correct_detections'], summary['recall']))
+
         print('Incorrect detections: {}\tPrecision: {}%'.format(
-            int(self.stats['n_pred'].sum() - self.stats['correct_detections'].sum()),
-            to_precision(100 * self.stats['correct_detections'].sum() / self.stats['n_pred'].sum(),
-                         self.ndigits)))
-
-        total_err = (self.stats['gained_detections'].sum()
-                     + self.stats['missed_detections'].sum()
-                     + self.stats['split'].sum()
-                     + self.stats['merge'].sum()
-                     + self.stats['catastrophe'].sum())
-
-        print('\nGained detections: {}\tPerc Error: {}%'.format(
-            int(self.stats['gained_detections'].sum()),
-            to_precision(100 * self.stats['gained_detections'].sum() / total_err, self.ndigits)))
-        print('Missed detections: {}\tPerc Error: {}%'.format(
-            int(self.stats['missed_detections'].sum()),
-            to_precision(100 * self.stats['missed_detections'].sum() / total_err, self.ndigits)))
-        print('Merges: {}\t\tPerc Error: {}%'.format(
-            int(self.stats['merge'].sum()),
-            to_precision(100 * self.stats['merge'].sum() / total_err, self.ndigits)))
-        print('Splits: {}\t\tPerc Error: {}%'.format(
-            int(self.stats['split'].sum()),
-            to_precision(100 * self.stats['split'].sum() / total_err, self.ndigits)))
-        print('Catastrophes: {}\t\tPerc Error: {}%\n'.format(
-            int(self.stats['catastrophe'].sum()),
-            to_precision(100 * self.stats['catastrophe'].sum() / total_err, self.ndigits)))
-
-        print('Gained detections from splits: {}'.format(
-            int(self.stats['gained_det_from_split'].sum())))
-        print('Missed detections from merges: {}'.format(
-            int(self.stats['missed_det_from_merge'].sum())))
-        print('True detections involved in catastrophes: {}'.format(
-            int(self.stats['true_det_in_catastrophe'].sum())))
-        print('Predicted detections involved in catastrophes: {}'.format(
-            int(self.stats['pred_det_in_catastrophe'].sum())), '\n')
-
-        if self.seg is True:
-            print('SEG:', to_precision(self.stats['seg'].mean(), self.ndigits), '\n')
+            summary['n_pred'] - summary['correct_detections'],
+            summary['precision']))
 
-        print('Average Pixel IOU (Jaccard Index):',
-              to_precision(self.stats['jaccard'].mean(), self.ndigits), '\n')
+        print('\n')
+        for k in errors:
+            v = summary[k]
+            name = k.replace('_', ' ').capitalize()
+            if not name.endswith('s'):
+                name += 's'
 
-    def run_all(self,
-                y_true_lbl,
-                y_pred_lbl,
-                y_true_unlbl,
-                y_pred_unlbl):
-        """Runs pixel and object base statistics and ouputs to file
+            try:
+                err_fraction = v / summary['total_errors']
+            except ZeroDivisionError:
+                err_fraction = 0
 
-        Args:
-            y_true_lbl (numpy.array): Labeled ground truth annotation,
-                (sample, x, y)
-            y_pred_lbl (numpy.array): Labeled prediction mask,
-                (sample, x, y)
-            y_true_unlbl (numpy.array): Ground truth annotation after necessary
-                transforms, (sample, x, y, feature)
-            y_pred_unlbl (numpy.array): Predictions, (sample, x, y, feature)
-        """
+            print('{name}: {val}{tab}Perc Error {percent}%'.format(
+                name=name, val=v,
+                percent=round(100 * err_fraction, self.ndigits),
+                tab='\t' * (1 if ' ' in name else 2)))
+
+        for k in bad_detections:
+            name = k.replace('_', ' ').capitalize().replace(' det ', ' detections')
+            print('{name}: {val}'.format(name=name, val=summary[k]))
 
-        logging.info('Starting pixel based statistics')
-        self.all_pixel_stats(y_true_unlbl, y_pred_unlbl)
+        print('SEG:', round(summary['seg'], self.ndigits), '\n')
 
-        logging.info('Starting object based statistics')
-        self.calc_object_stats(y_true_lbl, y_pred_lbl)
+        print('Average Pixel IOU (Jaccard Index):',
+              round(summary['jaccard'], self.ndigits), '\n')
 
-        self.save_to_json(self.output)
+    def run_all(self, y_true, y_pred, axis=-1):
+        object_metrics = self.calc_object_stats(y_true, y_pred)
+        pixel_metrics = self.calc_pixel_stats(y_true, y_pred, axis=axis)
+
+        object_list = self.df_to_dict(object_metrics, stat_type='object')
+        all_output = object_list + pixel_metrics
+        self.save_to_json(all_output)
 
     def save_to_json(self, L):
         """Save list of dictionaries to json file with file metadata
@@ -1036,7 +1236,7 @@ def save_to_json(self, L):
         """
         todays_date = datetime.datetime.now().strftime('%Y-%m-%d')
         outname = os.path.join(
-            self.outdir, self.model_name + '_' + todays_date + '.json')
+            self.outdir, '{}_{}.json'.format(self.model_name, todays_date))
 
         # Configure final output
         D = {}
@@ -1083,11 +1283,11 @@ def split_stack(arr, batch, n_split1, axis1, n_split2, axis2):
         >>> from deepcell import metrics
         >>> from numpy import np
         >>> arr = np.ones((10, 100, 100, 1))
-        >>> out = metrics.test_split_stack(arr, True, 10, 1, 10, 2)
+        >>> out = metrics.split_stack(arr, True, 10, 1, 10, 2)
         >>> out.shape
         (1000, 10, 10, 1)
         >>> arr = np.ones((100, 100, 1))
-        >>> out = metrics.test_split_stack(arr, False, 10, 1, 10, 2)
+        >>> out = metrics.split_stack(arr, False, 10, 1, 10, 2)
         >>> out.shape
         (100, 10, 10, 1)
     """
@@ -1112,18 +1312,19 @@ def split_stack(arr, batch, n_split1, axis1, n_split2, axis2):
     return split2con
 
 
-def match_nodes(gt, res):
+def match_nodes(y_true, y_pred):
     """Loads all data that matches each pattern and compares the graphs.
 
     Args:
-        gt (numpy.array): data array to match to unique.
-        res (numpy.array): ground truth array with all cells labeled uniquely.
+        y_true (numpy.array): ground truth array with all cells labeled uniquely.
+        y_pred (numpy.array): data array to match to unique.
 
     Returns:
         numpy.array: IoU of ground truth cells and predicted cells.
     """
-    num_frames = gt.shape[0]
-    iou = np.zeros((num_frames, np.max(gt) + 1, np.max(res) + 1))
+    num_frames = y_true.shape[0]
+    # TODO: does max make the shape bigger than necessary?
+    iou = np.zeros((num_frames, np.max(y_true) + 1, np.max(y_pred) + 1))
 
     # Compute IOUs only when neccesary
     # If bboxs for true and pred do not overlap with each other, the assignment
@@ -1131,8 +1332,8 @@ def match_nodes(gt, res):
 
     # Regionprops expects one frame at a time
     for frame in range(num_frames):
-        gt_frame = gt[frame]
-        res_frame = res[frame]
+        gt_frame = y_true[frame]
+        res_frame = y_pred[frame]
 
         gt_props = regionprops(np.squeeze(gt_frame.astype('int')))
         gt_boxes = [np.array(gt_prop.bbox) for gt_prop in gt_props]
@@ -1162,66 +1363,3 @@ def match_nodes(gt, res):
             iou[frame, iou_gt_idx, iou_res_idx] = intersection.sum() / union.sum()
 
     return iou
-
-
-def assign_plot_values(y_true, y_pred, error_dict):
-    """Generates a matrix with cells belong to error classes numbered for plotting
-
-    Args:
-        y_true: 2D matrix of true labels
-        y_pred 2D matrix of predicted labels
-        error_dict: dictionary produced by save_error_ids with IDs of all error cells
-
-    Returns:
-        plotting_tiff: 2D matrix with cells belonging to same error class having same value
-    """
-
-    plotting_tif = np.zeros_like(y_true)
-
-    # erode edges for easier visualization of adjacent cells
-    y_true = erode_edges(y_true, 1)
-    y_pred = erode_edges(y_pred, 1)
-
-    # missed detections are tracked with true labels
-    misses = error_dict.pop('misses')['y_true']
-    plotting_tif[np.isin(y_true, misses)] = 1
-
-    # all other events are tracked with predicted labels
-    category_id = 2
-    for key in error_dict.keys():
-        labels = error_dict[key]['y_pred']
-        plotting_tif[np.isin(y_pred, labels)] = category_id
-        category_id += 1
-
-    return plotting_tif
-
-
-def plot_errors(y_true, y_pred, error_dict):
-    """Plots the errors identified from linear assignment code
-
-    Due to sequential relabeling that occurs within the metrics code, only run
-    this plotting function on the outputs of save_error_ids so that values match up.
-
-    Args:
-        y_true: 2D matrix of true labels returned by save_error_ids
-        y_pred: 2D matrix of predicted labels returned by save_error_ids
-        error_dict: dictionary returned by save_error_ids with IDs of all error cells
-    """
-
-    import matplotlib as mpl
-    import matplotlib.pyplot as plt
-
-    plotting_tif = assign_plot_values(y_true, y_pred, error_dict)
-
-    plotting_colors = ['Black', 'Pink', 'Blue', 'Green', 'tan', 'Red', 'Grey']
-    cmap = mpl.colors.ListedColormap(plotting_colors)
-
-    fig, ax = plt.subplots(nrows=1, ncols=1)
-    mat = ax.imshow(plotting_tif, cmap=cmap, vmin=np.min(plotting_tif) - .5,
-                    vmax=np.max(plotting_tif) + .5)
-
-    # tell the colorbar to tick at integers
-    cbar = fig.colorbar(mat, ticks=np.arange(np.min(plotting_tif), np.max(plotting_tif) + 1))
-    cbar.ax.set_yticklabels(['Background', 'misses', 'splits', 'merges',
-                             'gains', 'catastrophes', 'correct'])
-    fig.tight_layout()
diff --git a/setup.py b/setup.py
index 43d848f..3e07eaf 100644
--- a/setup.py
+++ b/setup.py
@@ -105,7 +105,8 @@ def run(self, *args, **kwargs):
                         'numpy>=1.16.6,<1.20.0',
                         'scipy>=1.3.0',
                         'scikit-image',
-                        'scikit-learn'],
+                        'scikit-learn',
+                        'tqdm'],
       extras_require={
           'tests': ['pytest<6',
                     'pytest-pep8',


Test patch:
diff --git a/deepcell_toolbox/metrics_test.py b/deepcell_toolbox/metrics_test.py
index 273d680..473ed29 100644
--- a/deepcell_toolbox/metrics_test.py
+++ b/deepcell_toolbox/metrics_test.py
@@ -263,562 +263,471 @@ def _dense_sample():
     return true.astype('int'), pred.astype('int')
 
 
-class TestMetricFunctions():
+class TestDetection():
 
-    def test_pixelstats_output(self):
-        y_true = _get_image()
-        y_pred = _get_image()
+    def test_init(self):
+        # test correct detection
+        true_idx, pred_idx = 1, 1
+        detection = metrics.Detection(true_idx, pred_idx)
+        assert detection.is_correct
 
-        out2 = metrics.stats_pixelbased(y_true, y_pred)
-        assert isinstance(out2, dict)
+        # test missed detection
+        true_idx, pred_idx = 1, None
+        detection = metrics.Detection(true_idx, pred_idx)
+        assert detection.is_missed
 
-        # Test mistmatch size error
-        with pytest.raises(ValueError):
-            metrics.stats_pixelbased(np.ones((10, 10)), np.ones((20, 20)))
+        # test gained detection
+        true_idx, pred_idx = None, 1
+        detection = metrics.Detection(true_idx, pred_idx)
+        assert detection.is_gained
 
-        # Test empty input arguments
-        y_true = np.zeros_like(y_true)
-        y_pred = np.zeros_like(y_pred)
+        # test split detection
+        true_idx, pred_idx = 1, [1, 2]
+        detection = metrics.Detection(true_idx, pred_idx)
+        assert detection.is_split
 
-        with pytest.warns(UserWarning):
-            out3 = metrics.stats_pixelbased(y_true, y_pred)
-
-    def test_split_stack(self):
-        # Test batch True condition
-        arr = np.ones((10, 100, 100, 1))
-        out = metrics.split_stack(arr, True, 10, 1, 10, 2)
-        outshape = (10 * 10 * 10, 100 / 10, 100 / 10, 1)
-        testing.assert_equal(outshape, out.shape)
-
-        # Test batch False condition
-        arr = np.ones((100, 100, 1))
-        out = metrics.split_stack(arr, False, 10, 0, 10, 1)
-        outshape = (10 * 10, 100 / 10, 100 / 10, 1)
-        testing.assert_equal(outshape, out.shape)
-
-        # Test splitting in only one axis
-        out = metrics.split_stack(arr, False, 10, 0, 1, 1)
-        outshape = (10 * 1, 100 / 10, 100 / 1, 1)
-        testing.assert_equal(outshape, out.shape)
-
-        out = metrics.split_stack(arr, False, 1, 0, 10, 1)
-        outshape = (10 * 1, 100 / 1, 100 / 10, 1)
-        testing.assert_equal(outshape, out.shape)
-
-        # Raise errors for uneven division
-        with pytest.raises(ValueError):
-            metrics.split_stack(arr, False, 11, 0, 10, 1)
-        with pytest.raises(ValueError):
-            metrics.split_stack(arr, False, 10, 0, 11, 1)
+        # test merge detection
+        true_idx, pred_idx = [1, 2], 1
+        detection = metrics.Detection(true_idx, pred_idx)
+        assert detection.is_merge
 
+        # test catastrophe
+        true_idx, pred_idx = [1, 2], [2, 3]
+        detection = metrics.Detection(true_idx, pred_idx)
+        assert detection.is_catastrophe
 
-class TestMetricsObject():
+    def test_hash(self):
+        # test that Detections get hashed appropriately
+        detection_set = set()
+        d1 = metrics.Detection(1, 1)
+        detection_set.add(d1)
 
-    def test_Metrics_init(self):
-        m = metrics.Metrics('test')
+        d2 = metrics.Detection(1, 1)
+        assert d2 is not d1
+        # should be in the set since d2 == d1
+        assert d2 in detection_set
 
-        testing.assert_equal(hasattr(m, 'output'), True)
+    def test_eq(self):
+        # test Detection equality comparisons
+        detection_set = set()
+        d1 = metrics.Detection(1, 1)
+        detection_set.add(d1)
 
-    def test_all_pixel_stats(self):
-        m = metrics.Metrics('test')
+        d2 = metrics.Detection(1, 1)
+        assert d2 is not d1
+        assert d2 == d1
 
-        before = len(m.output)
+        assert d1 != 1
 
-        y_true = _generate_stack_4d()
-        y_pred = _generate_stack_4d()
+        print(d1)  # test that __repr__ is called
 
-        m.all_pixel_stats(y_true, y_pred)
 
-        # Check that items were added to output
-        assert before != len(m.output)
+class TestPixelMetrics():
 
-        # Check mismatch error
-        testing.assert_raises(ValueError, m.all_pixel_stats, np.ones(
-            (10, 10, 10, 1)), np.ones((5, 5, 5, 1)))
+    def test_init(self):
+        y_true, _ = _sample1(10, 10, 30, 30, True)
 
-    def test_df_to_dict(self):
-        m = metrics.Metrics('test')
-        df = _generate_df()
+        # Test basic initialization
+        o = metrics.PixelMetrics(y_true, y_true)
 
-        L = m.pixel_df_to_dict(df)
+        # Test mismatched input size
+        with pytest.raises(ValueError):
+            metrics.PixelMetrics(y_true, y_true[0])
 
-        # Check output types
-        assert len(L) != 0
-        assert isinstance(L, list)
-        assert isinstance(L[0], dict)
+        # using float dtype warns but still works
+        o = metrics.PixelMetrics(y_true.astype('float'), y_true.astype('float'))
 
-    def test_confusion_matrix(self):
-        y_true = _generate_stack_4d()
-        y_pred = _generate_stack_4d()
+    def test_y_true_equals_y_pred(self):
+        y_true, _ = _sample1(10, 10, 30, 30, True)
+        y_pred = y_true.copy()
 
-        m = metrics.Metrics('test')
+        # Test basic initialization
+        o = metrics.PixelMetrics(y_true, y_pred)
 
-        cm = m.calc_pixel_confusion_matrix(y_true, y_pred)
-        testing.assert_equal(cm.shape[0], y_true.shape[-1])
+        # metrics should be perfect since y_true == y_pred
+        assert o.recall == 1
+        assert o.precision == 1
+        assert o.f1 == 1
+        assert o.jaccard == 1
+        assert o.dice == 1
 
-    def test_metric_object_stats(self):
-        y_true = label(_generate_stack_3d())
-        y_pred = label(_generate_stack_3d())
+        # test both empty
+        o = metrics.PixelMetrics(np.zeros_like(y_true), np.zeros_like(y_pred))
+        assert o.dice == 1
 
-        m = metrics.Metrics('test')
-        before = len(m.output)
+    def test_y_pred_empty(self):
+        y_true, _ = _sample1(10, 10, 30, 30, True)
 
-        m.calc_object_stats(y_true, y_pred)
+        # Test basic initialization with empty array y_pred
+        o = metrics.PixelMetrics(y_true, np.zeros_like(y_true))
 
-        # Check data added to output
-        assert before != len(m.output)
+        assert o.recall == 0
+        assert o.precision == 0
+        assert o.f1 == 0
+        assert o.jaccard == 0
 
-        m.calc_object_stats(np.zeros_like(y_true), np.zeros_like(y_pred))
-        assert m.stats['precision'].sum() == 0
-        assert m.stats['recall'].sum() == 0
+    def test_y_true_empty(self):
+        y_pred, _ = _sample1(10, 10, 30, 30, True)
 
-        # Raise input size error
-        with testing.assert_raises(ValueError):
-            m.calc_object_stats(np.random.rand(10, 10), np.random.rand(10, 10))
+        # Test basic initialization with empty array y_pred
+        o = metrics.PixelMetrics(np.zeros_like(y_pred), y_pred)
 
-        # data that needs to be relabeled raises a warning
-        with pytest.warns(UserWarning):
-            y_pred[0, 0, 0] = 40
-            m.calc_object_stats(y_true, y_pred)
+        assert np.isnan(o.recall)
+        assert o.precision == 0
+        assert np.isnan(o.f1)
+        assert o.jaccard == 0
 
-    def test_save_to_json(self, tmpdir):
-        name = 'test'
-        tmpdir = str(tmpdir)
-        m = metrics.Metrics(name, outdir=tmpdir)
 
-        # Create test list to save
-        L = []
-        for i in range(10):
-            L.append(dict(
-                name=i,
-                value=i,
-                feature='test',
-                stat_type='output'
-            ))
+class TestObjectMetrics():
 
-        m.save_to_json(L)
-        todays_date = datetime.datetime.now().strftime('%Y-%m-%d')
-        outfilename = os.path.join(tmpdir, name + '_' + todays_date + '.json')
+    def test_init(self):
+        y_true, _ = _sample1(10, 10, 30, 30, True)
 
-        # Check that file exists
-        testing.assert_equal(os.path.isfile(outfilename), True)
+        # Test basic initialization
+        o = metrics.ObjectMetrics(y_true, y_true)
 
-        # Check that it can be opened
-        with open(outfilename) as json_file:
-            data = json.load(json_file)
+        # Test __repr__
+        print(o)
 
-        # Check data types from loaded data
-        assert isinstance(data, dict)
-        assert np.array_equal(sorted(list(data.keys())), ['metadata', 'metrics'])
-        assert isinstance(data['metrics'], list)
-        assert isinstance(data['metadata'], dict)
+        # Test using float dtype warns but still works
+        o = metrics.ObjectMetrics(
+            y_true.astype('float'),
+            y_true.astype('float'))
 
-    def test_run_all(self, tmpdir):
-        tmpdir = str(tmpdir)
-        y_true_lbl = label(_generate_stack_3d())
-        y_pred_lbl = label(_generate_stack_3d())
-        y_true_unlbl = _generate_stack_4d()
-        y_pred_unlbl = _generate_stack_4d()
+        # test errors thrown for improper ndim inputs
+        y_true = np.zeros(shape=(10))  # too few dimensions
+        with pytest.raises(ValueError):
+            metrics.ObjectMetrics(y_true, y_true)
 
-        name = 'test'
-        for seg in [True, False]:
-            m = metrics.Metrics(name, outdir=tmpdir, seg=seg)
+        y_true = np.zeros(shape=(10, 5, 5, 5))  # too many dimensions
+        with pytest.raises(ValueError):
+            metrics.ObjectMetrics(y_true, y_true)
 
-            before = len(m.output)
+        # test errors thrown for improper ndim inputs with 3d data
+        y_true = np.zeros(shape=(10, 15))  # too few dimensions
+        with pytest.raises(ValueError):
+            metrics.ObjectMetrics(y_true, y_true, is_3d=True)
 
-            m.run_all(y_true_lbl, y_pred_lbl, y_true_unlbl, y_pred_unlbl)
+        y_true = np.zeros(shape=(10, 15, 15, 10))  # too many dimensions
+        with pytest.raises(ValueError):
+            metrics.ObjectMetrics(y_true, y_true, is_3d=True)
 
-            # Assert that data was added to output
-            assert len(m.output) != before
+        # Test mismatched input size
+        with pytest.raises(ValueError):
+            metrics.ObjectMetrics(y_true, y_true[0])
 
-            # Check output file
-            todays_date = datetime.datetime.now().strftime('%Y-%m-%d')
-            outname = os.path.join(tmpdir, name + '_' + todays_date + '.json')
-            testing.assert_equal(os.path.isfile(outname), True)
+    def test__get_props(self):
+        y_true, y_pred = _sample1(10, 10, 30, 30, True)
+        o = metrics.ObjectMetrics(y_true, y_true)
 
+        props = o._get_props('correct')
+        # assert props != []
 
-class TestObjectAccuracy():
+        # Test _get_props with invalid detection type
+        with pytest.raises(ValueError):
+            o._get_props('invalid_type')
 
-    def test_init(self):
+    def test_y_true_equals_y_pred(self):
         y_true, _ = _sample1(10, 10, 30, 30, True)
+        y_pred = y_true.copy()
 
         # Test basic initialization
-        o = metrics.ObjectAccuracy(y_true, y_true, test=True)
+        o = metrics.ObjectMetrics(y_true, y_pred, force_event_links=True)
 
         # Check that object numbers are integers
         assert isinstance(o.n_true, int)
         assert isinstance(o.n_pred, int)
+        assert o.n_true == o.n_pred
+
+        # metrics should be perfect since y_true == y_pred
+        assert o.recall == 1
+        assert o.precision == 1
+        assert o.f1 == 1
+        assert o.jaccard == 1
+        assert o.dice == 1
+
+        # test other properties
+        assert o.correct_detections == o.n_pred
+        assert o.missed_detections == 0
+        assert o.gained_detections == 0
+        assert o.splits == 0
+        assert o.merges == 0
+        assert o.catastrophes == 0
+        assert o.gained_det_from_split == 0
+        assert o.missed_det_from_merge == 0
+        assert o.true_det_in_catastrophe == 0
+        assert o.pred_det_in_catastrophe == 0
+
+        assert o.missed_props == []
+        assert o.merge_props == []
+        assert o.split_props == []
+        assert o.gained_props == []
+
+    def test_y_pred_empty(self):
+        y_true, _ = _sample1(10, 10, 30, 30, True)
 
-        testing.assert_equal(o.empty_frame, False)
+        # Test basic initialization with empty array y_pred
+        o = metrics.ObjectMetrics(y_true, np.zeros_like(y_true),
+                                  force_event_links=True)
+
+        assert o.n_pred == 0
+        assert o.correct_detections == 0
+        assert o.missed_detections == o.n_true
+        assert o.recall == 0
+        assert o.precision == 0
+        assert o.f1 == 0
+        assert o.jaccard == 0
+        assert o.gained_det_from_split == 0
+        assert o.missed_det_from_merge == 0
+        assert o.true_det_in_catastrophe == 0
+        assert o.pred_det_in_catastrophe == 0
+
+    def test_y_true_empty(self):
+        y_pred, _ = _sample1(10, 10, 30, 30, True)
+
+        # Test basic initialization with empty array y_pred
+        o = metrics.ObjectMetrics(np.zeros_like(y_pred), y_pred,
+                                  force_event_links=True)
+
+        assert o.n_true == 0
+        assert o.correct_detections == 0
+        assert o.gained_detections == o.n_pred
+        assert o.recall == 0
+        assert o.precision == 0
+        assert o.f1 == 0
+        assert o.jaccard == 0
+        assert o.gained_det_from_split == 0
+        assert o.missed_det_from_merge == 0
+        assert o.true_det_in_catastrophe == 0
+        assert o.pred_det_in_catastrophe == 0
+
+    def test_merge_error(self):
+        y_true, y_pred = _sample1(10, 10, 30, 30, True)
+        o = metrics.ObjectMetrics(y_true, y_pred)
+        assert o.merges == 1
+        assert o.missed_det_from_merge == 1
+        assert o.merge_props != []
 
-        # test errors thrown for improper ndim inputs
-        y_true = np.zeros(shape=(10))
-        with pytest.raises(ValueError):
-            o = metrics.ObjectAccuracy(y_true, y_true, test=True)
+    def test_split_error(self):
+        y_true, y_pred = _sample1(10, 10, 30, 30, False)
+        o = metrics.ObjectMetrics(y_true, y_pred)
+        assert o.splits == 1
+        assert o.gained_det_from_split == 1
+        assert o.split_props != []
+
+    def test_multi_merge_error(self):
+        # 3 cells merged together
+        # forced event links to ensure accurate assignment
+        y_true, y_pred = _sample2_3(10, 10, 30, 30,
+                                    merge=True, similar_size=False)
+        o = metrics.ObjectMetrics(y_true, y_pred, force_event_links=True)
+        assert o.merges == 1
+        assert o.missed_det_from_merge == 2
+        assert o.merge_props != []
+
+    def test_multi_split_error(self):
+        # 1 cell split into 3
+        # forced event links to ensure accurate assignment
+        y_true, y_pred = _sample2_3(10, 10, 30, 30,
+                                    merge=False, similar_size=False)
+        o = metrics.ObjectMetrics(y_true, y_pred, force_event_links=True)
+        assert o.splits == 1
+        assert o.gained_det_from_split == 2
+        assert o.split_props != []
 
-        y_true = np.zeros(shape=(10, 5, 5, 5))
-        with pytest.raises(ValueError):
-            o = metrics.ObjectAccuracy(y_true, y_true, test=True)
+    def test_calc_iou(self):
+        # TODO: test correctness
+        y_true, y_pred = _sample1(10, 10, 30, 30, True)
+        o = metrics.ObjectMetrics(y_true, y_pred)
 
-        # test errors thrown for improper ndim inputs with 3d data
-        y_true = np.zeros(shape=(10, 15))
-        y_pred = y_true
-        with pytest.raises(ValueError):
-            o = metrics.ObjectAccuracy(y_true, y_pred, test=True, is_3d=True)
+        # Check that it is not equal to initial value
+        assert np.count_nonzero(o.iou) != 0
+        assert np.count_nonzero(o.seg_thresh) != 0
 
-        y_true = np.zeros(shape=(10, 15, 15, 10))
-        y_pred = y_true
-        with pytest.raises(ValueError):
-            o = metrics.ObjectAccuracy(y_true, y_pred, test=True, is_3d=True)
+        # check that image without any background passes
+        y_true, y_pred = _dense_sample()
+        o = metrics.ObjectMetrics(y_true=y_true, y_pred=y_pred)
 
-    def test_init_wrongsize(self):
-        # Test mismatched input size
-        y_true = label(_get_image())
-        y_wrong = label(_get_image(img_h=200, img_w=200))
-        testing.assert_raises(ValueError, metrics.ObjectAccuracy, y_true, y_wrong)
+        # Check that it is not equal to initial value
+        assert np.count_nonzero(o.iou) != 0
+        assert np.count_nonzero(o.seg_thresh) != 0
 
-    def test_init_emptyframe(self):
-        y_true, y_empty = _sample1(10, 10, 30, 30, True)
+    def test_calc_iou_3D(self):
+        # TODO: test correctness
+        y_true, y_pred = _sample1_3D(10, 10, 30, 30, True, 8)
+        o = metrics.ObjectMetrics(y_true, y_pred, is_3d=True)
 
-        # Check handling of empty frames
-        y_empty[:, :] = 0
-        y_empty = y_empty.astype('int')
+        # Check that it is not equal to initial value
+        assert np.count_nonzero(o.iou) != 0
+        assert np.count_nonzero(o.seg_thresh) != 0
 
-        oempty = metrics.ObjectAccuracy(y_true, y_empty)
-        testing.assert_equal(oempty.empty_frame, 'n_pred')
-        oempty = metrics.ObjectAccuracy(y_empty, y_true)
-        testing.assert_equal(oempty.empty_frame, 'n_true')
 
-    def test_init_noloners(self):
-        # Generate split error and use single cell in true for both pred/true
-        y_true, y_pred = _sample1(10, 10, 30, 30, False)
+class TestMetrics():
 
-        o = metrics.ObjectAccuracy(y_true, y_true)
+    def test_df_to_dict(self):
+        m = metrics.Metrics('test')
+        df = _generate_df()
 
-        # Look for creation of `self.n_pred2` which indicates assignloners ran inccorectly
-        testing.assert_equal(hasattr(o, 'n_pred2'), False)
+        L = m.df_to_dict(df)
 
-    def test_calc_iou(self):
-        y_true, y_pred = _sample1(10, 10, 30, 30, True)
-        o = metrics.ObjectAccuracy(y_true, y_pred, test=True)
+        # Check output types
+        assert len(L) != 0
+        assert isinstance(L, list)
+        assert isinstance(L[0], dict)
 
-        o._calc_iou()
+    def test_calc_pixel_stats(self):
+        m = metrics.Metrics('test')
 
-        # Check that iou was created
-        assert hasattr(o, 'iou')
+        y_true = _generate_stack_4d()
+        y_pred = _generate_stack_4d()
 
-        # Check that it is not equal to initial value
-        assert np.count_nonzero(o.iou) != 0
+        pixel_stats = m.calc_pixel_stats(y_true, y_pred)
 
-        # Test seg_thresh creation
-        o = metrics.ObjectAccuracy(y_true, y_pred, test=True, seg=True)
-        o._calc_iou()
+        for stat in pixel_stats:
+            assert 'name' in stat
+            assert 'stat_type' in stat
 
-        assert hasattr(o, 'seg_thresh')
+    def test_confusion_matrix(self):
+        y_true = _generate_stack_4d()
+        y_pred = _generate_stack_4d()
 
-        # check that image without any background passes
-        y_true, y_pred = _dense_sample()
-        o = metrics.ObjectAccuracy(y_true=y_true, y_pred=y_pred, test=False)
+        m = metrics.Metrics('test')
 
-        o._calc_iou()
+        cm = m.calc_pixel_confusion_matrix(y_true, y_pred)
+        testing.assert_equal(cm.shape[0], y_true.shape[-1])
 
-    def test_calc_iou_3D(self):
-        y_true, y_pred = _sample1_3D(10, 10, 30, 30, True, 8)
-        o = metrics.ObjectAccuracy(y_true, y_pred, test=True, is_3d=True)
+    def test_calc_object_stats(self):
+        y_true = label(_generate_stack_3d())
+        y_pred = label(_generate_stack_3d())
 
-        o._calc_iou()
+        m = metrics.Metrics('test')
 
-        # Check that iou was created
-        assert hasattr(o, 'iou')
+        # test that metrics are generated
+        object_metrics = m.calc_object_stats(y_true, y_pred)
+        # each row of metrics corresponds to a batch
+        assert len(object_metrics) == len(y_true)
 
-        # Check that it is not equal to initial value
-        assert np.count_nonzero(o.iou) != 0
+        object_metrics = m.calc_object_stats(
+            np.zeros_like(y_true), np.zeros_like(y_pred))
 
-        # Test seg_thresh creation
-        o = metrics.ObjectAccuracy(y_true, y_pred, test=True, seg=True)
-        o._calc_iou()
+        # test accuracy of metrics with blank predictions
+        assert object_metrics['precision'].sum() == 0
+        assert object_metrics['recall'].sum() == 0
 
-        assert hasattr(o, 'seg_thresh')
+        # Raise input size error
+        with testing.assert_raises(ValueError):
+            m.calc_object_stats(np.random.rand(10, 10), np.random.rand(10, 10))
 
-        m = metrics.Metrics('test', is_3d=True)
+        # Raise error if y_pred.shape != y_true.shape
+        with testing.assert_raises(ValueError):
+            m.calc_object_stats(np.random.rand(10, 10), np.random.rand(10,))
 
-        # test errors thrown for improper ndim inputs
-        y_true = np.zeros(shape=(10, 15, 11))
-        with pytest.raises(ValueError):
-            m.calc_object_stats(y_true, y_true)
+        # data that needs to be relabeled raises a warning
+        with pytest.warns(UserWarning):
+            y_pred[0, 0, 0] = 40
+            m.calc_object_stats(y_true, y_pred)
 
-        y_true = np.zeros(shape=(10, 15, 15, 10, 15))
-        with pytest.raises(ValueError):
-            m.calc_object_stats(y_true, y_true)
+        # seg is deprecated (TODO: this will be removed)
+        with pytest.warns(DeprecationWarning):
+            _ = metrics.Metrics('test', seg=True)
 
-        y_true = np.zeros(shape=(2, 3, 5, 2))
-        y_pred = np.zeros(shape=(1, 4, 11, 2))
-        with pytest.raises(ValueError):
-            m.calc_object_stats(y_true, y_pred)
+    def test_calc_object_stats_3d(self):
+        y_true = _generate_stack_4d()
+        y_pred = _generate_stack_4d()
 
-    def test_modify_iou(self):
-        y_true, y_pred = _sample1(10, 10, 30, 30, True)
-        o = metrics.ObjectAccuracy(y_true, y_pred, test=True)
+        m = metrics.Metrics('test', is_3d=True)
 
-        o._calc_iou()
-        o._modify_iou(force_event_links=False)
+        # test that metrics are generated
+        object_metrics = m.calc_object_stats(y_true, y_pred)
+        # each row of metrics corresponds to a batch
+        assert len(object_metrics) == len(y_true)
 
-        # Check that modified_iou was created
-        assert hasattr(o, 'iou_modified')
+        # test accuracy of metrics with blank predictions
+        object_metrics = m.calc_object_stats(
+            np.zeros_like(y_true), np.zeros_like(y_pred))
 
-    def test_make_matrix(self):
-        y_true, y_pred = _sample1(10, 10, 30, 30, True)
-        o = metrics.ObjectAccuracy(y_true, y_pred, test=True)
-        o._calc_iou()
-        o._modify_iou(force_event_links=False)
+        assert object_metrics['precision'].sum() == 0
+        assert object_metrics['recall'].sum() == 0
+
+        # Raise error if is_3d and ndim !=4
+        with testing.assert_raises(ValueError):
+            m3d = metrics.Metrics('test', is_3d=True)
+            m3d.calc_object_stats(np.random.random((32, 32, 1)),
+                                  np.random.random((32, 32, 1)))
+
+    def test_run_all(self, tmpdir):
+        tmpdir = str(tmpdir)
+        y_true = label(_generate_stack_3d())
+        y_pred = label(_generate_stack_3d())
 
-        o._make_matrix()
+        name = 'test'
+        m = metrics.Metrics(name, outdir=tmpdir)
 
-        assert hasattr(o, 'cm')
+        m.run_all(y_true, y_pred)
 
-        assert np.count_nonzero(o.cm) != 0
+    def test_save_to_json(self, tmpdir):
+        name = 'test'
+        tmpdir = str(tmpdir)
+        m = metrics.Metrics(name, outdir=tmpdir)
 
-    def test_linear_assignment(self):
-        y_true, y_pred = _sample1(10, 10, 30, 30, True)
-        o = metrics.ObjectAccuracy(y_true, y_pred, test=True)
-        o._calc_iou()
-        o._modify_iou(force_event_links=False)
-        o._make_matrix()
+        # Create test list to save
+        L = []
+        for i in range(10):
+            L.append(dict(
+                name=i,
+                value=i,
+                feature='test',
+                stat_type='output'
+            ))
 
-        o._linear_assignment()
+        m.save_to_json(L)
+        todays_date = datetime.datetime.now().strftime('%Y-%m-%d')
+        outfilename = os.path.join(tmpdir, name + '_' + todays_date + '.json')
 
-        cols = ['n_pred', 'n_true', 'correct_detections', 'missed_detections', 'gained_detections',
-                'missed_det_from_merge', 'gained_det_from_split', 'true_det_in_catastrophe',
-                'pred_det_in_catastrophe', 'merge', 'split', 'catastrophe']
+        # Check that file exists
+        testing.assert_equal(os.path.isfile(outfilename), True)
 
-        for obj in cols:
-            assert hasattr(o, obj)
+        # Check that it can be opened
+        with open(outfilename) as json_file:
+            data = json.load(json_file)
 
-        # Test condition where seg = True
-        o = metrics.ObjectAccuracy(y_true, y_pred, test=True, seg=True)
-        o._calc_iou()
-        o._modify_iou(force_event_links=False)
-        o._make_matrix()
-        o._linear_assignment()
+        # Check data types from loaded data
+        assert isinstance(data, dict)
+        assert np.array_equal(sorted(list(data.keys())), ['metadata', 'metrics'])
+        assert isinstance(data['metrics'], list)
+        assert isinstance(data['metadata'], dict)
 
-        for obj in ['results', 'cm_res', 'seg_score']:
-            assert hasattr(o, obj)
 
-    def test_assign_loners(self):
-        y_true, y_pred = _sample1(10, 10, 30, 30, True)
-        o = metrics.ObjectAccuracy(y_true, y_pred, test=True)
-        o._calc_iou()
-        o._modify_iou(force_event_links=False)
-        o._make_matrix()
-        o._linear_assignment()
+def test__cast_to_tuple():
+    assert metrics._cast_to_tuple(None) == ()
+    assert metrics._cast_to_tuple(1) == (1,)
+    assert metrics._cast_to_tuple((1,)) == (1,)
 
-        o._assign_loners()
-        assert hasattr(o, 'cost_l_bin')
 
-    def test_array_to_graph(self):
-        y_true, y_pred = _sample1(10, 10, 30, 30, True)
-        o = metrics.ObjectAccuracy(y_true, y_pred, test=True)
-        o._calc_iou()
-        o._modify_iou(force_event_links=False)
-        o._make_matrix()
-        o._linear_assignment()
-        o._assign_loners()
+def test_split_stack():
+    # Test batch True condition
+    arr = np.ones((10, 100, 100, 1))
+    out = metrics.split_stack(arr, True, 10, 1, 10, 2)
+    outshape = (10 * 10 * 10, 100 / 10, 100 / 10, 1)
+    testing.assert_equal(outshape, out.shape)
 
-        o._array_to_graph()
-        assert hasattr(o, 'G')
+    # Test batch False condition
+    arr = np.ones((100, 100, 1))
+    out = metrics.split_stack(arr, False, 10, 0, 10, 1)
+    outshape = (10 * 10, 100 / 10, 100 / 10, 1)
+    testing.assert_equal(outshape, out.shape)
 
-    def test_classify_graph(self):
-        y_true, y_pred = _sample1(10, 10, 30, 30, True)
-        # Test that complete run through is succesful
-        _ = metrics.ObjectAccuracy(y_true, y_pred)
+    # Test splitting in only one axis
+    out = metrics.split_stack(arr, False, 10, 0, 1, 1)
+    outshape = (10 * 1, 100 / 10, 100 / 1, 1)
+    testing.assert_equal(outshape, out.shape)
 
-        # Test for 0 degree graph
-        y_true, y_pred = _sample4_loner(10, 10, 30, 30, True)
-        _ = metrics.ObjectAccuracy(y_true, y_pred)
-        y_true, y_pred = _sample4_loner(10, 10, 30, 30, False)
-        _ = metrics.ObjectAccuracy(y_true, y_pred)
+    out = metrics.split_stack(arr, False, 1, 0, 10, 1)
+    outshape = (10 * 1, 100 / 1, 100 / 10, 1)
+    testing.assert_equal(outshape, out.shape)
 
-        # Test for splits in 1 degree graph
-        y_true, y_pred = _sample1(10, 10, 30, 30, False)
-        _ = metrics.ObjectAccuracy(y_true, y_pred)
-
-        # Test for catastrophic errors
-        y_true, y_pred = _sample3(10, 10, 30, 30)
-        _ = metrics.ObjectAccuracy(y_true, y_pred)
-
-    def test_save_error_ids(self):
-
-        # cell 1 in assigned correctly, cells 2 and 3 have been merged
-        y_true, y_pred = _sample1(10, 10, 30, 30, merge=True)
-        o = metrics.ObjectAccuracy(y_true, y_pred)
-        label_dict, _, _ = o.save_error_ids()
-        assert label_dict['correct']['y_true'] == [1]
-        assert label_dict['correct']['y_pred'] == [1]
-        assert set(label_dict['merges']['y_true']) == {2, 3}
-        assert label_dict['merges']['y_pred'] == [2]
-
-        # cell 1 in assigned correctly, cell 2 has been split
-        y_true, y_pred = _sample1(10, 10, 30, 30, merge=False)
-        o = metrics.ObjectAccuracy(y_true, y_pred)
-        label_dict, _, _ = o.save_error_ids()
-        assert label_dict['correct']['y_true'] == [1]
-        assert label_dict['correct']['y_pred'] == [1]
-        assert set(label_dict['splits']['y_pred']) == {2, 3}
-        assert label_dict['splits']['y_true'] == [2]
-
-        # gained cell in predictions
-        y_true, y_pred = _sample4_loner(10, 10, 30, 30, gain=True)
-        o = metrics.ObjectAccuracy(y_true, y_pred, cutoff1=0.2, cutoff2=0.1)
-        label_dict, _, _ = o.save_error_ids()
-        assert label_dict['correct']['y_true'] == [1]
-        assert label_dict['correct']['y_pred'] == [1]
-        assert label_dict['gains']['y_pred'] == [2]
-
-        # missed cell in true
-        y_true, y_pred = _sample4_loner(10, 10, 30, 30, gain=False)
-        o = metrics.ObjectAccuracy(y_true, y_pred, cutoff1=0.2, cutoff2=0.1)
-        label_dict, _, _ = o.save_error_ids()
-        assert label_dict['correct']['y_true'] == [1]
-        assert label_dict['correct']['y_pred'] == [1]
-        assert label_dict['misses']['y_true'] == [2]
-
-        # catastrophe between 3 cells
-        y_true, y_pred = _sample3(10, 10, 30, 30)
-        o = metrics.ObjectAccuracy(y_true, y_pred, cutoff1=0.2, cutoff2=0.1)
-        label_dict, _, _ = o.save_error_ids()
-        assert set(label_dict['catastrophes']['y_true']) == set(np.unique(y_true[y_true > 0]))
-        assert set(label_dict['catastrophes']['y_pred']) == set(np.unique(y_pred[y_pred > 0]))
-
-        # The tests below are more stochastic, and should be run multiple times
-        for _ in range(10):
-
-            # 3 cells merged together, with forced event links to ensure accurate assignment
-            y_true, y_pred = _sample2_3(10, 10, 30, 30, merge=True, similar_size=False)
-            o = metrics.ObjectAccuracy(y_true, y_pred, force_event_links=True,
-                                       cutoff1=0.2, cutoff2=0.1)
-            label_dict, _, _ = o.save_error_ids()
-            assert label_dict['correct']['y_true'] == [1]
-            assert label_dict['correct']['y_pred'] == [1]
-            assert set(label_dict['merges']['y_true']) == {2, 3, 4}
-            assert label_dict['merges']['y_pred'] == [2]
-
-            # 3 cells merged together, without forced event links. Cells must be similar size
-            y_true, y_pred = _sample2_3(10, 10, 30, 30, merge=True, similar_size=True)
-            o = metrics.ObjectAccuracy(y_true, y_pred, force_event_links=False,
-                                       cutoff1=0.2, cutoff2=0.1)
-            label_dict, _, _ = o.save_error_ids()
-            assert label_dict['correct']['y_true'] == [1]
-            assert label_dict['correct']['y_pred'] == [1]
-            assert set(label_dict['merges']['y_true']) == {2, 3, 4}
-            assert label_dict['merges']['y_pred'] == [2]
-
-            # 2 of 3 cells merged together, with forced event links to ensure accurate assignment
-            y_true, y_pred, y_true_merge, y_true_correct, y_pred_merge, y_pred_correct = \
-                _sample2_2(10, 10, 30, 30, similar_size=False)
-            o = metrics.ObjectAccuracy(y_true, y_pred, cutoff1=0.2, cutoff2=0.1,
-                                       force_event_links=True)
-            label_dict, _, _ = o.save_error_ids()
-            assert set(label_dict['correct']['y_true']) == y_true_correct
-            assert set(label_dict['correct']['y_pred']) == y_pred_correct
-            assert set(label_dict['merges']['y_true']) == y_true_merge
-            assert set(label_dict['merges']['y_pred']) == y_pred_merge
-
-            # 2 of 3 cells merged together, without forced event links. Cells must be similar size
-            y_true, y_pred, y_true_merge, y_true_correct, y_pred_merge, y_pred_correct = \
-                _sample2_2(10, 10, 30, 30, similar_size=True)
-            o = metrics.ObjectAccuracy(y_true, y_pred, cutoff1=0.2, cutoff2=0.1,
-                                       force_event_links=False)
-            label_dict, _, _ = o.save_error_ids()
-            assert set(label_dict['correct']['y_true']) == y_true_correct
-            assert set(label_dict['correct']['y_pred']) == y_pred_correct
-            assert set(label_dict['merges']['y_true']) == y_true_merge
-            assert set(label_dict['merges']['y_pred']) == y_pred_merge
-
-            # 1 cell split into three pieces, with forced event links to ensure accurate assignment
-            y_true, y_pred = _sample2_3(10, 10, 30, 30, merge=False, similar_size=False)
-            o = metrics.ObjectAccuracy(y_true, y_pred, cutoff1=0.2, cutoff2=0.1,
-                                       force_event_links=True)
-            label_dict, _, _ = o.save_error_ids()
-            assert label_dict['correct']['y_true'] == [1]
-            assert label_dict['correct']['y_pred'] == [1]
-            assert label_dict['splits']['y_true'] == [2]
-            assert set(label_dict['splits']['y_pred']) == {2, 3, 4}
-
-            # 1 cell split in three pieces, without forced event links. Cells must be similar size
-            y_true, y_pred = _sample2_3(10, 10, 30, 30, merge=False, similar_size=True)
-            o = metrics.ObjectAccuracy(y_true, y_pred, cutoff1=0.2, cutoff2=0.1,
-                                       force_event_links=False)
-            label_dict, _, _ = o.save_error_ids()
-            assert label_dict['correct']['y_true'] == [1]
-            assert label_dict['correct']['y_pred'] == [1]
-            assert label_dict['splits']['y_true'] == [2]
-            assert set(label_dict['splits']['y_pred']) == {2, 3, 4}
-
-            # 1 cell split into two pieces, one small accurate cell, with forced event links
-            y_true, y_pred, y_true_split, y_true_correct, y_pred_split, y_pred_correct = \
-                _sample2_2(10, 10, 30, 30, merge=False, similar_size=False)
-            o = metrics.ObjectAccuracy(y_true, y_pred, cutoff1=0.2, cutoff2=0.1,
-                                       force_event_links=True)
-            label_dict, _, _ = o.save_error_ids()
-            assert set(label_dict['correct']['y_true']) == y_true_correct
-            assert set(label_dict['correct']['y_pred']) == y_pred_correct
-            assert set(label_dict['splits']['y_true']) == y_true_split
-            assert set(label_dict['splits']['y_pred']) == y_pred_split
-
-            # 1 cell split into two pieces, one small accurate cell, without forced event links
-            y_true, y_pred, y_true_split, y_true_correct, y_pred_split, y_pred_correct = \
-                _sample2_2(10, 10, 30, 30, merge=False, similar_size=True)
-            o = metrics.ObjectAccuracy(y_true, y_pred, cutoff1=0.2, cutoff2=0.1,
-                                       force_event_links=False)
-            label_dict, _, _ = o.save_error_ids()
-            assert set(label_dict['correct']['y_true']) == y_true_correct
-            assert set(label_dict['correct']['y_pred']) == y_pred_correct
-            assert set(label_dict['splits']['y_true']) == y_true_split
-            assert set(label_dict['splits']['y_pred']) == y_pred_split
-
-    def test_optional_outputs(self):
-        y_true, y_pred = _sample1(10, 10, 30, 30, True)
-        o = metrics.ObjectAccuracy(y_true, y_pred)
-
-        o.print_report()
-
-        df = o.save_to_dataframe()
-        assert isinstance(df, pd.DataFrame)
-
-        columns = ['n_pred', 'n_true', 'correct_detections', 'missed_detections',
-                   'gained_detections', 'missed_det_from_merge', 'gained_det_from_split',
-                   'true_det_in_catastrophe', 'pred_det_in_catastrophe', 'merge', 'split',
-                   'catastrophe', 'jaccard', 'precision', 'recall', 'f1']
-        assert np.array_equal(sorted(columns), sorted(list(df.columns)))
-
-        # Check seg True case
-        o = metrics.ObjectAccuracy(y_true, y_pred, seg=True)
-        o.print_report()
-        df = o.save_to_dataframe()
-        columns = ['n_pred', 'n_true', 'correct_detections', 'missed_detections',
-                   'gained_detections', 'missed_det_from_merge', 'gained_det_from_split',
-                   'true_det_in_catastrophe', 'pred_det_in_catastrophe', 'merge', 'split',
-                   'catastrophe', 'seg', 'jaccard', 'precision', 'recall', 'f1']
-        assert np.array_equal(sorted(columns), sorted(list(df.columns)))
-
-    def test_assign_plot_values(self):
-        y_true, _ = random_shapes(image_shape=(200, 200), max_shapes=30, min_shapes=15,
-                                  min_size=10, multichannel=False)
-
-        # invert background
-        y_true[y_true == 255] = 0
-        y_true, _, _ = relabel_sequential(y_true)
-
-        error_dict = {'misses': {'y_true': [1, 2]}, 'splits': {'y_pred': [3, 4]},
-                      'merges': {'y_pred': [5, 6]}, 'gains': {'y_pred': [7, 8]},
-                      'catastrophes': {'y_pred': [9, 10]}, 'correct': {'y_pred': [11, 12]}}
-
-        plotting_tiff = metrics.assign_plot_values(y_true, y_true, error_dict)
-
-        # erode edges so that shape matches shape in plotting tiff
-        y_true = erode_edges(y_true, 1)
-
-        for error_type in error_dict.keys():
-            vals = list(error_dict[error_type].values())
-            mask = np.isin(y_true, vals)
-            assert len(np.unique(plotting_tiff[mask])) == 1
+    # Raise errors for uneven division
+    with pytest.raises(ValueError):
+        metrics.split_stack(arr, False, 11, 0, 10, 1)
+    with pytest.raises(ValueError):
+        metrics.split_stack(arr, False, 10, 0, 11, 1)