Documentation for Box.py

class ArrayGrid(Grid):
    """
    A dense, lossy implementation of the grid.
    """

    @staticmethod
    def init_boxes(rtboxes, mzboxes):
        return np.array([[False for mz in mzboxes] for rt in rtboxes])

    def approx_non_overlap(self, box):
        rt_box_range, mz_box_range, total_boxes = self.get_box_ranges(box)
        boxes = self.boxes[rt_box_range[0]: rt_box_range[1], mz_box_range[0]: mz_box_range[1]]
        return (total_boxes - np.sum(boxes)) / total_boxes

    def register_box(self, box):
        rt_box_range, mz_box_range, _ = self.get_box_ranges(box)
        self.boxes[rt_box_range[0]: rt_box_range[1], mz_box_range[0]: mz_box_range[1]] = True

class BoxExact(BoxGeometry):
    @abstractmethod
    def _get_overlapping_boxes(self, box):
        pass

    @staticmethod
    def non_overlap_boxes(box, other_boxes):
        new_boxes = [box]
        for b in other_boxes:
            if box.overlaps_with_box(b):  # quickly exits any box not overlapping new box
                updated_boxes = []
                for b2 in new_boxes:
                    # if your box is contained within a
                    # previous box area is 0 and box is not carried over
                    if not b.contains_box(b2):
                        split_boxes = b2.non_overlap_split(b)
                        if split_boxes is not None:
                            updated_boxes.extend(split_boxes)
                        else:
                            updated_boxes.append(b2)
                if not updated_boxes:
                    return []
                new_boxes = updated_boxes
        return new_boxes

    @staticmethod
    def splitting_non_overlap(box, other_boxes):
        boxes = BoxExact.non_overlap_boxes(box, other_boxes)
        if boxes == []:
            return 0.0
        else:
            return sum(b.area() for b in boxes) / box.area()

    def non_overlap(self, box):
        return float(self.splitting_non_overlap(box, self._get_overlapping_boxes(box)))

    @staticmethod
    def splitting_intensity_non_overlap(box, other_boxes, current_intensity, scoring_params):
        """Will give nonsense results if other_boxes overlap each other."""
        areas = [box.overlap_raw(b) for b in other_boxes]
        non_overlapped = max(0.0, 1.0 - sum(areas) / box.area())
        non_overlap = current_intensity**non_overlapped
        refragment = sum(
            max(0.0, current_intensity - b.intensity) ** (area / box.area())
            for b, area in zip(other_boxes, areas)
        )
        return non_overlap + scoring_params["theta1"] * refragment

    def intensity_non_overlap(self, box, current_intensity, scoring_params):
        return BoxExact.splitting_intensity_non_overlap(
            box=box,
            other_boxes=self._get_overlapping_boxes(box),
            current_intensity=current_intensity,
            scoring_params=scoring_params,
        )

    def flexible_non_overlap(self, box, current_intensity, scoring_params):
        other_boxes = self._get_overlapping_boxes(box)
        areas = [box.overlap_raw(b) for b in other_boxes]
        non_overlapped = max(0.0, 1.0 - sum(areas) / box.area())
        norm_areas = [ar / box.area() for ar in areas]
        # NB: this seems wrong, but i tried to preserve it as it was - vinny?
        intensity_diff = [
            np.log(current_intensity) - np.log(max(1.0, b.intensity)) * na
            for b, na in zip(other_boxes, norm_areas)
        ]

        non_overlap = np.log(current_intensity) * non_overlapped
        refragment = sum(max(0.0, int_diff) for int_diff in intensity_diff)
        refragment2 = sum(intensity_diff)

        new_peak_score = sum(
            np.log(current_intensity) * na
            for b, na in zip(other_boxes, norm_areas)
            if math.isclose(0.0, b.intensity)
        )

        return (
            non_overlap
            + scoring_params["theta1"] * refragment
            + scoring_params["theta2"] * refragment2
            + scoring_params["theta3"] * new_peak_score
        )

    def case_control_non_overlap(self, box, current_intensity, scoring_params):
        other_boxes = self._get_overlapping_boxes(box)
        areas = [box.overlap_raw(b) for b in other_boxes]
        non_overlapped = max(0.0, 1.0 - sum(areas) / box.area())
        norm_areas = [ar / box.area() for ar in areas]
        # NB: this seems wrong, but i tried to preserve it as it was - vinny?
        intensity_diff = [
            np.log(current_intensity) - np.log(max(1.0, b.intensity)) * na
            for b, na in zip(other_boxes, norm_areas)
        ]

        non_overlap = np.log(current_intensity) * non_overlapped
        refragment = sum(max(0.0, int_diff) for int_diff in intensity_diff)
        refragment2 = sum(intensity_diff)

        new_peak_score = sum(
            np.log(current_intensity) * na
            for b, na in zip(other_boxes, norm_areas)
            if math.isclose(0.0, b.intensity)
        )

        if box.pvalue is None:
            model_score = 0.0
        else:
            model_score = sum(
                max(0.0, np.log(current_intensity) - max(0.0, np.log(b.intensity)) * na)
                for b, na in zip(other_boxes, norm_areas)
            ) * (1 - box.pvalue)

        return (
            non_overlap
            + scoring_params["theta1"] * refragment
            + scoring_params["theta2"] * refragment2
            + scoring_params["theta3"] * new_peak_score
            + scoring_params["theta4"] * model_score
        )

class BoxGeometry(metaclass=ABCMeta):
    """
    Describes the interface for an abstract class which can do geometric
    operations on points, intervals and rectangles. Different subclasses use
    different data structures, and hence the choice of data structure matters
    for performance.
    """

    def set_active_boxes(self, *args):
        pass

    @abstractmethod
    def get_all_boxes(self):
        pass

    @abstractmethod
    def point_in_box(self, pt):
        pass

    @abstractmethod
    def point_in_which_boxes(self, pt):
        pass

    @abstractmethod
    def interval_overlaps_which_boxes(self, inv):
        pass

    @abstractmethod
    def interval_covers_which_boxes(self, inv):
        pass

    @abstractmethod
    def non_overlap(self, box):
        pass

    @abstractmethod
    def intensity_non_overlap(self, box):
        pass

    @abstractmethod
    def register_boxes(self, boxes):
        pass

    @abstractmethod
    def clear(self):
        pass

class BoxLineSweeper(BoxExact):
    def __init__(self):
        self.lswp = LineSweeper()
        self.running_scores = dict()

    def get_all_boxes(self):
        return self.lswp.get_all_boxes()

    def set_active_boxes(self, current_loc):
        self.lswp.set_active_boxes(current_loc)

    def point_in_box(self, pt):
        return self.lswp.point_in_box(pt)

    def point_in_which_boxes(self, pt):
        return self.lswp.point_in_which_boxes(pt)

    def interval_overlaps_which_boxes(self, inv):
        return self.lswp.interval_overlaps_which_boxes(inv)

    def interval_covers_which_boxes(self, inv):
        return self.lswp.interval_covers_which_boxes(inv)

    def _get_overlapping_boxes(self, box):
        raise NotImplementedError("Not used in this class: this line shouldn't be reached.")

    def non_overlap(self, box):
        """NB: This won't work if the boxes are capable of moving between time updates."""
        sliced = box.copy()
        sliced.pt1.x = self.lswp.previous_loc

        other_boxes = (
            self.lswp.interval_overlaps_which_boxes(
                Interval(
                    self.lswp.current_loc,
                    self.lswp.current_loc,
                    sliced.pt1.y,
                    sliced.pt2.y,
                )
            )
            + self.lswp.was_active
        )  # TODO: If the manual intersection check is removed from this
        # non-overlap, then filter to intersecting boxes
        sliced_uncovered = sum(b.area() for b in BoxExact.non_overlap_boxes(sliced, other_boxes))

        running_uncovered, running_total = self.running_scores.get(sliced.id, (0, 0))
        running_uncovered += sliced_uncovered
        running_total += sliced.area()
        self.running_scores[sliced.id] = (running_uncovered, running_total)

        return running_uncovered / running_total

    def intensity_non_overlap(self, box):
        raise NotImplementedError("BoxLineSweeper doesn't implement this yet!")

    def register_boxes(self, boxes):
        self.lswp.register_boxes(boxes)

    def clear(self):
        self.__init__()

class DictGrid(Grid):
    """
    A sparse, lossless implementation of the grid.
    """

    @staticmethod
    def init_boxes(rtboxes, mzboxes):
        return defaultdict(list)

    def approx_non_overlap(self, box):
        rt_box_range, mz_box_range, total_boxes = self.get_box_ranges(box)
        active = sum(
            float(not self.boxes[(rt, mz)])
            for rt in range(*rt_box_range)
            for mz in range(*mz_box_range)
        )
        return active / total_boxes

    def register_box(self, box):
        rt_box_range, mz_box_range, _ = self.get_box_ranges(box)
        for rt in range(*rt_box_range):
            for mz in range(*mz_box_range):
                self.boxes[(rt, mz)].append(box)

class GenericBox(Box):
    """Makes no particular assumptions about bounding boxes."""

    def __repr__(self):
        return "Generic{}".format(super().__repr__())

    def contains_point(self, pt):
        return (
            self.pt1.x <= pt.x and self.pt1.y <= pt.y and self.pt2.x >= pt.x and self.pt2.y >= pt.y
        )

    def interval_contains(self, inv):
        if inv.is_vertical():
            return (
                self.pt1.x <= inv.pt1.x
                and self.pt2.x >= inv.pt1.x
                and self.pt1.y >= inv.pt1.y
                and self.pt2.y <= inv.pt2.y
            )
        else:
            return (
                self.pt1.x >= inv.pt1.x
                and self.pt2.x <= inv.pt2.x
                and self.pt1.y <= inv.pt1.y
                and self.pt2.y >= inv.pt1.y
            )

    def overlaps_with_box(self, other_box):
        return (
            self.pt1.x < other_box.pt2.x
            and self.pt2.x > other_box.pt1.x
            and self.pt1.y < other_box.pt2.y
            and self.pt2.y > other_box.pt1.y
        )

    def contains_box(self, other_box):
        return (
            self.pt1.x <= other_box.pt1.x
            and self.pt1.y <= other_box.pt1.y
            and self.pt2.x >= other_box.pt2.x
            and self.pt2.y >= other_box.pt2.y
        )

    def combine_max(self, other_box):
        return GenericBox(
            min(self.pt1.x, other_box.pt1.x),
            max(self.pt2.x, other_box.pt2.x),
            min(self.pt1.y, other_box.pt1.y),
            max(self.pt2.y, other_box.pt2.y),
            parents=(self.parents + other_box.parents),
            intensity=max(self.intensity, other_box.intensity),
        )

    def apply_min_box_ppm(self, xwidth=None, ywidth=None, round_digits=8):
        x1, y1 = self.pt1
        x2, y2 = self.pt2

        if xwidth is not None:
            mid = (x1 + x2) / 2
            dist = (xwidth / 1e6) * mid
            if dist > x2 - x1:
                x1 = mid - dist / 2
                x2 = mid + dist / 2

        if ywidth is not None:
            mid = (y1 + y2) / 2
            dist = (ywidth / 1e6) * mid
            if dist > y2 - y1:
                y1 = mid - dist / 2
                y2 = mid + dist / 2

        new_box = self.copy()
        new_box.pt1.x, new_box.pt1.y = x1, y1
        new_box.pt2.x, new_box.pt2.y = x2, y2
        new_box.round(ndigits=round_digits)
        return new_box

    def overlap_raw(self, other_box):
        if not self.overlaps_with_box(other_box):
            return 0.0
        return (min(self.pt2.x, other_box.pt2.x) - max(self.pt1.x, other_box.pt1.x)) * (
            min(self.pt2.y, other_box.pt2.y) - max(self.pt1.y, other_box.pt1.y)
        )

    def overlap_2(self, other_box):
        if not self.overlaps_with_box(other_box):
            return 0.0
        b = type(self)(
            max(self.pt1.x, other_box.pt1.x),
            min(self.pt2.x, other_box.pt2.x),
            max(self.pt1.y, other_box.pt1.y),
            min(self.pt2.y, other_box.pt2.y),
        )
        return b.area() / (self.area() + other_box.area() - b.area())

    def overlap_3(self, other_box):
        if not self.overlaps_with_box(other_box):
            return 0.0
        b = type(self)(
            max(self.pt1.x, other_box.pt1.x),
            min(self.pt2.x, other_box.pt2.x),
            max(self.pt1.y, other_box.pt1.y),
            min(self.pt2.y, other_box.pt2.y),
        )
        return b.area() / self.area()

    def non_overlap_split(self, other_box):
        """Finds 1 to 4 boxes describing the polygon of area of this box
        not overlapped by other_box. If one box is found, crops this box to
        dimensions of that box, and returns None. Otherwise, returns list of
        2 to 4 boxes. Number of boxes found is equal to number of edges
        overlapping area does NOT share with this box."""
        if not self.overlaps_with_box(other_box):
            return None
        x1, x2, y1, y2 = self.pt1.x, self.pt2.x, self.pt1.y, self.pt2.y
        split_boxes = []
        if other_box.pt1.x > self.pt1.x:
            x1 = other_box.pt1.x
            split_boxes.append(
                type(self)(self.pt1.x, x1, y1, y2, parents=self.parents, intensity=self.intensity)
            )
        if other_box.pt2.x < self.pt2.x:
            x2 = other_box.pt2.x
            split_boxes.append(
                type(self)(x2, self.pt2.x, y1, y2, parents=self.parents, intensity=self.intensity)
            )
        if other_box.pt1.y > self.pt1.y:
            y1 = other_box.pt1.y
            split_boxes.append(
                type(self)(x1, x2, self.pt1.y, y1, parents=self.parents, intensity=self.intensity)
            )
        if other_box.pt2.y < self.pt2.y:
            y2 = other_box.pt2.y
            split_boxes.append(
                type(self)(x1, x2, y2, self.pt2.y, parents=self.parents, intensity=self.intensity)
            )
        return split_boxes

    def split_all(self, other_box):
        if not self.overlaps_with_box(other_box):
            return None, None, None
        both_parents = self.parents + other_box.parents
        both_box = type(self)(
            max(self.pt1.x, other_box.pt1.x),
            min(self.pt2.x, other_box.pt2.x),
            max(self.pt1.y, other_box.pt1.y),
            min(self.pt2.y, other_box.pt2.y),
            parents=both_parents,
            intensity=max(self.intensity, other_box.intensity),
        )
        b1_boxes = self.non_overlap_split(other_box)
        b2_boxes = other_box.non_overlap_split(self)
        return b1_boxes, b2_boxes, both_box

class Grid(metaclass=ABCMeta):
    """
    Partitions a 2D space into a number of rectangles of fixed size for faster
    lookup. If a query object and a saved object touch the same rectangle, then
    the saved object should be factored into the query.
    """

    @staticmethod
    @abstractmethod
    def init_boxes():
        pass

    def __init__(self, min_rt, max_rt, rt_box_size, min_mz, max_mz, mz_box_size):
        self.min_rt, self.max_rt = min_rt, max_rt
        self.min_mz, self.max_mz = min_mz, max_mz
        self.rt_box_size, self.mz_box_size = rt_box_size, mz_box_size
        self.box_area = float(Decimal(rt_box_size) * Decimal(mz_box_size))

        self.rtboxes = range(0, int((self.max_rt - self.min_rt) / rt_box_size) + 2)
        self.mzboxes = range(0, int((self.max_mz - self.min_mz) / mz_box_size) + 2)
        self.boxes = self.init_boxes(self.rtboxes, self.mzboxes)

    def get_box_ranges(self, box):
        rt_box_range = (
            int((box.pt1.x - self.min_rt) / self.rt_box_size),
            int((box.pt2.x - self.min_rt) / self.rt_box_size) + 1,
        )
        mz_box_range = (
            int((box.pt1.y - self.min_mz) / self.mz_box_size),
            int((box.pt2.y - self.min_mz) / self.mz_box_size) + 1,
        )
        total_boxes = (rt_box_range[1] - rt_box_range[0]) * (mz_box_range[1] - mz_box_range[0])
        return rt_box_range, mz_box_range, total_boxes

    @abstractmethod
    def register_box(self, box):
        pass

    def clear(self):
        self.__init__(
            self.min_rt, self.max_rt, self.rt_box_size, self.min_mz, self.max_mz, self.mz_box_size
        )

class LocatorGrid(Grid):
    """
    A dense, lossless implementation of the grid.
    """

    @staticmethod
    def init_boxes(rtboxes, mzboxes):
        arr = np.empty((max(rtboxes), max(mzboxes)), dtype=object)
        for i, row in enumerate(arr):
            for j, _ in enumerate(row):
                arr[i, j] = set()
        return arr

    def get_boxes(self, box):
        rt_box_range, mz_box_range, _ = self.get_box_ranges(box)
        boxes = set()
        for row in self.boxes[
            rt_box_range[0]: rt_box_range[1], mz_box_range[0]: mz_box_range[1]
        ]:
            for s in row:
                boxes |= s
        return boxes

    def get_all_boxes(self):
        return reduce(or_, (s for row in self.boxes for s in row), set())

    def register_box(self, box):
        rt_box_range, mz_box_range, _ = self.get_box_ranges(box)
        for row in self.boxes[
            rt_box_range[0]: rt_box_range[1], mz_box_range[0]: mz_box_range[1]
        ]:
            for s in row:
                s.add(box)

    def clear(self):
        for i, row in enumerate(self.boxes):
            for j, _ in enumerate(row):
                self.boxes[i, j] = set()