Use the last RR interval in compute_hr

docs/changes.rst

@@ -29,6 +29,9 @@ Unreleased
 **Fix parameter handling in csv_to_wfdb**
   ``csv_to_wfdb`` now builds the same ``Record`` whether it writes files or returns the record (``record_only=True``), so all header parameters (such as ``dat_file_name``, ``counter_freq``, ``base_counter``, ``skew``, ``byte_offset``, ``adc_res`` and ``adc_zero``) are written instead of being ignored. Single values of any numeric type (including NumPy scalars and floats) are used for all channels, and an explicit ``0`` is no longer replaced by the default. The initial value and checksum are computed from the digital samples, and the default ADC resolution now follows ``fmt`` in both modes. The ``samps_per_frame``, ``init_value`` and ``checksum`` arguments are deprecated and ignored. The parameter descriptions in the docstring have been corrected.
 
+**Use the last RR interval in compute_hr**
+  ``compute_hr`` now uses the interval between the last two beats for the heart rate after the last beat, instead of repeating the previous rate. With exactly two beats it now returns their rate after the second beat instead of only NaN.
+
 **Use Ruff for formatting and linting**
   The ``dev`` extra now includes ``ruff`` instead of ``black`` and ``pylint``, and the code is formatted with a line length of 88 (the default).
 

src/wfdb/processing/hr.py

@@ -9,6 +9,10 @@
 def compute_hr(sig_len, qrs_inds, fs):
     """Compute instantaneous heart rate from peak indices.
 
+    The rate from each RR interval is assigned to the samples after its second
+    beat, up to and including the next beat. After the last beat, the rate
+    from the last RR interval is used up to the end of the signal.
+
     Parameters
     ----------
     sig_len : int
@@ -31,15 +35,15 @@ def compute_hr(sig_len, qrs_inds, fs):
     if len(qrs_inds) < 2:
         return heart_rate
 
-    for i in range(0, len(qrs_inds) - 2):
+    # The rate from each RR interval holds from just after its second beat up
+    # to the next beat, or to the end of the signal after the last one.
+    for i in range(0, len(qrs_inds) - 1):
         a = qrs_inds[i]
         b = qrs_inds[i + 1]
-        c = qrs_inds[i + 2]
+        end = qrs_inds[i + 2] + 1 if i + 2 < len(qrs_inds) else sig_len
         rr = (b - a) * (1.0 / fs) * 1000
         hr = 60000.0 / rr
-        heart_rate[b + 1 : c + 1] = hr
-
-    heart_rate[qrs_inds[-1] :] = heart_rate[qrs_inds[-1]]
+        heart_rate[b + 1 : end] = hr
 
     return heart_rate
 

tests/test_processing.py

@@ -174,6 +174,37 @@ def test_correct_peaks(self):
         assert np.array_equal(yz, expected_peaks)
 
 
+class TestComputeHr:
+    """Tests for compute_hr."""
+
+    def test_uses_last_rr_interval(self):
+        """The rate after the last beat comes from the last RR interval."""
+        hr = processing.compute_hr(sig_len=14, qrs_inds=[0, 2, 6, 12], fs=60)
+
+        expected = np.full(14, np.nan)
+        expected[3:7] = 60 * 60 / 2
+        expected[7:13] = 60 * 60 / 4
+        expected[13:] = 60 * 60 / 6
+        np.testing.assert_allclose(hr, expected)
+
+    def test_two_beats(self):
+        """Two beats give one RR interval, so the rate is known after them."""
+        hr = processing.compute_hr(sig_len=10, qrs_inds=[2, 6], fs=1)
+
+        expected = np.full(10, np.nan)
+        expected[7:] = 60 / 4
+        np.testing.assert_allclose(hr, expected)
+
+    def test_last_beat_past_end(self):
+        """A last beat beyond the signal doesn't raise an IndexError."""
+        hr = processing.compute_hr(sig_len=8, qrs_inds=[0, 2, 6, 12], fs=60)
+
+        expected = np.full(8, np.nan)
+        expected[3:7] = 60 * 60 / 2
+        expected[7:] = 60 * 60 / 4
+        np.testing.assert_allclose(hr, expected)
+
+
 class TestQrs:
     """Testing QRS detectors"""