Match each test annotation at most once in compare_annotations (#606)

`compare_annotations` can match one test annotation to two reference
annotations:

```python
c = processing.compare_annotations(np.array([0, 2, 4, 5]), np.array([1, 5]), window_width=6)
c.matching_sample_nums  # [0, -1, 0, 1]: test annotation 0 is used twice
c.fp, c.positive_predictivity  # -1, 1.5
```

When a reference annotation loses a contested test annotation to the
next reference, `compare` falls back to the previous test annotation. It
only checks that the *previous reference* didn't take it, so if that
reference was unmatched, a test annotation matched earlier on can be
matched again.

**Fix:** keep track of the last matched test annotation. Matches are
made in increasing order, so the fallback is used only if it comes after
that one. This replaces the old check (which it covers).

**Checking that nothing else changes:**
- On 200,000 random small cases, the new code gives exactly the same
matches as `main` whenever `main` doesn't match a test annotation twice
(about 6% of the cases did). With the fix, none do.
- On record 100 (gqrs and xqrs against the reference beats, two window
sizes), the matches are identical.

**Tests:** the example above, plus a randomized check that matches are
unique, in order and within the window, and that `fp` is never negative.
Both fail on `main`. `pytest -m "not network"` and the docs build pass.

This is separate from #566, which handles empty inputs in `_calc_stats`.

Author: tompollard

docs/changes.rst

@@ -23,6 +23,9 @@ Unreleased
 **Speed up XQRS backsearch**
   ``XQRS`` no longer rechecks peaks that an earlier backsearch already rejected. Long stretches without QRS complexes (for example a disconnected lead) made detection take quadratic time, up to hours on long records. The detected QRS locations are unchanged.
 
+**Fix test annotations matched twice in compare_annotations**
+  ``compare_annotations`` (``Comparitor.compare``) could match one test annotation to two reference annotations, which gave a negative number of false positives and a positive predictivity above 1. Each test annotation is now matched at most once. Results are unchanged when no annotation was matched twice.
+
 **Fix record name and add write directory in csv2ann**
   ``csv2ann`` now takes the record name from the CSV file name without its directory and extension, as ``csv_to_wfdb`` does. Paths with a directory (such as ``data/x.csv``) used to fail when writing, and ``./x.csv`` gave an empty name. The annotation file is written to the new ``write_dir`` argument, which defaults to the current directory. ``csv2ann`` is now included in the API documentation.
 

src/wfdb/processing/evaluate.py

@@ -155,6 +155,9 @@ def compare(self):
         """
         test_samp_num = 0
         ref_samp_num = 0
+        # The last test annotation that was matched. Test annotations are
+        # matched in increasing order, so none up to this one can be reused.
+        last_matched_test_num = -1
 
         # Iterate through the reference sample numbers
         while ref_samp_num < self.n_ref and test_samp_num < self.n_test:
@@ -186,11 +189,7 @@ def compare(self):
                 # if not already assigned to a previous sample.
                 # It will be the previous testing sample number in any
                 # possible case (scenario D below), or nothing.
-                if closest_samp_num and (
-                    not ref_samp_num
-                    or closest_samp_num - 1
-                    != self.matching_sample_nums[ref_samp_num - 1]
-                ):
+                if closest_samp_num - 1 > last_matched_test_num:
                     # The previous test annotation is inspected
                     closest_samp_num = closest_samp_num - 1
                     smallest_samp_diff = abs(
@@ -200,6 +199,7 @@ def compare(self):
                     # Assign the reference-test pair if close enough
                     if smallest_samp_diff < self.window_width:
                         self.matching_sample_nums[ref_samp_num] = closest_samp_num
+                        last_matched_test_num = closest_samp_num
                     # Set the starting test sample number to inspect
                     # for the next reference sample.
                     test_samp_num = closest_samp_num + 1
@@ -213,6 +213,7 @@ def compare(self):
                 # Assign the reference-test pair if close enough
                 if smallest_samp_diff < self.window_width:
                     self.matching_sample_nums[ref_samp_num] = closest_samp_num
+                    last_matched_test_num = closest_samp_num
                 # Increment the starting test sample number to inspect
                 # for the next reference sample.
                 test_samp_num = closest_samp_num + 1

tests/test_processing.py

@@ -244,6 +244,35 @@ def test_rejects_unordered_unsigned(self):
             processing.compute_hr(sig_len=20, qrs_inds=qrs_inds, fs=1)
 
 
+class TestCompareAnnotations:
+    """Tests for compare_annotations."""
+
+    def test_test_annotation_matched_once(self):
+        """A test annotation isn't matched to two reference annotations."""
+        comparitor = processing.compare_annotations(
+            np.array([0, 2, 4, 5]), np.array([1, 5]), window_width=6
+        )
+
+        assert comparitor.matching_sample_nums.tolist() == [0, -1, -1, 1]
+        assert (comparitor.tp, comparitor.fp, comparitor.fn) == (2, 0, 2)
+        assert comparitor.positive_predictivity == 1.0
+
+    def test_random_annotations(self):
+        """Matches are unique, in order and within the window."""
+        rng = np.random.default_rng(0)
+        for _ in range(2000):
+            ref = np.sort(rng.choice(40, rng.integers(1, 9), replace=False))
+            test = np.sort(rng.choice(40, rng.integers(1, 9), replace=False))
+            window_width = int(rng.integers(1, 10))
+            comparitor = processing.compare_annotations(ref, test, window_width)
+
+            matches = comparitor.matching_sample_nums
+            matched = matches[matches != -1]
+            assert np.all(np.diff(matched) > 0)
+            assert np.all(np.abs(ref[matches != -1] - test[matched]) < window_width)
+            assert 0 <= comparitor.fp <= len(test)
+
+
 class TestQrs:
     """Testing QRS detectors"""