Merge pull request #44 from smoia/feat/hrv

Fix naming of metrics, clarify PULSE RATE vs HEART RATE

Author: smoia

phys2denoise/metrics/cardiac.py

@@ -15,6 +15,10 @@ def _cardiac_metrics(card, peaks, samplerate, metric, window=6, central_measure=
     Computes the average heart beats interval (HBI)
     or the average heart rate variability (HRV) in a sliding window.
 
+    We refer to HEART RATE (variability), however note that if using a PPG
+    recorded signal, it is more accurate to talk about PULSE RATE (variability).
+    See [3]_ for the differences and similarities between the two measures.
+
     Parameters
     ----------
     card : list or 1D numpy.ndarray
@@ -23,14 +27,15 @@ def _cardiac_metrics(card, peaks, samplerate, metric, window=6, central_measure=
         array of peak indexes for card.
     samplerate : float
         Sampling rate for card, in Hertz.
-    metrics : "hbi", "hrv", string
+    metrics : "hbi", "hr", "hrv", string
         Cardiac metric(s) to calculate.
     window : float, optional
         Size of the sliding window, in seconds.
         Default is 6.
-    central_measure : "mean","average", "avg", "median", "mdn",  string, optional
+    central_measure : "mean","average", "avg", "median", "mdn", "stdev", "std", string, optional
         Measure of the center used (mean or median).
         Default is "mean".
+
     Returns
     -------
     card_met : 2D numpy.ndarray
@@ -62,6 +67,9 @@ def _cardiac_metrics(card, peaks, samplerate, metric, window=6, central_measure=
         vol. 213, pp. 116707, 2020.
     .. [2] C. Chang, J. P. Cunningham, & G. H. Glover, "Influence of heart rate on the
         BOLD signal: The cardiac response function", NeuroImage, vol. 44, 2009
+    .. [3] N. Pinheiro et al., "Can PPG be used for HRV analysis?," 2016 38th
+    Annual International Conference of the IEEE Engineering in Medicine and
+    Biology Society (EMBC), doi: 10.1109/EMBC.2016.7591347.
     """
     # Convert window to samples, but halves it.
     halfwindow_samples = int(round(window * samplerate / 2))
@@ -70,6 +78,8 @@ def _cardiac_metrics(card, peaks, samplerate, metric, window=6, central_measure=
         central_measure_operator = np.mean
     elif central_measure in ["median", "mdn"]:
         central_measure_operator = np.median
+    elif central_measure in ["stdev", "std"]:
+        central_measure_operator = np.std
     else:
         raise NotImplementedError(
             f" {central_measure} is not a supported metric of centrality."
@@ -86,7 +96,10 @@ def _cardiac_metrics(card, peaks, samplerate, metric, window=6, central_measure=
         )
         if metric == "hbi":
             card_met[n] = central_measure_operator(diff) if diff.size > 0 else 0
+        elif metric == "hr":
+            card_met[n] = central_measure_operator(1 / diff) if diff.size > 0 else 0
         elif metric == "hrv":
+            central_measure_operator = np.std
             card_met[n] = central_measure_operator(1 / diff) if diff.size > 0 else 0
         else:
             raise NotImplementedError(
@@ -102,9 +115,13 @@ def _cardiac_metrics(card, peaks, samplerate, metric, window=6, central_measure=
 
 
 @due.dcite(references.CHANG_CUNNINGHAM_GLOVER_2009)
-def heart_rate_variability(card, peaks, samplerate, window=6, central_measure="mean"):
+def heart_rate(card, peaks, samplerate, window=6, central_measure="mean"):
     """
-    Compute average heart rate variability (HRV) in a sliding window.
+    Compute average heart rate (HR) in a sliding window.
+
+    We call this function HEART RATE, however note that if using a PPG
+    recorded signal, it is more accurate to talk about PULSE RATE.
+    See [2]_ for the differences and similarities between the two measures.
 
     Parameters
     ----------
@@ -130,7 +147,7 @@ def heart_rate_variability(card, peaks, samplerate, window=6, central_measure="m
 
     Notes
     -----
-    Heart rate variability (HRV) is taken from [1]_, and computed as the amounts of
+    Heart rate (HR) is taken from [1]_, and computed as the amounts of
     beats per minute.
     However, operationally, it is the average of the inverse of the time interval
     between two heart beats.
@@ -146,12 +163,71 @@ def heart_rate_variability(card, peaks, samplerate, window=6, central_measure="m
     ----------
     .. [1] C. Chang, J. P. Cunningham, & G. H. Glover, "Influence of heart rate on the
         BOLD signal: The cardiac response function", NeuroImage, vol. 44, 2009
+    .. [2] N. Pinheiro et al., "Can PPG be used for HRV analysis?," 2016 38th
+    Annual International Conference of the IEEE Engineering in Medicine and
+    Biology Society (EMBC), doi: 10.1109/EMBC.2016.7591347.
     """
     return _cardiac_metrics(
         card, peaks, samplerate, metric="hrv", window=6, central_measure="mean"
     )
 
 
+@due.dcite(references.PINHERO_ET_AL_2016)
+def heart_rate_variability(card, peaks, samplerate, window=6, central_measure="mean"):
+    """
+    Compute average heart rate variability (HRV) in a sliding window.
+
+    We call this function HEART RATE variability, however note that if using a PPG
+    recorded signal, it is more accurate to talk about PULSE RATE variability.
+    See [1]_ for the differences and similarities between the two measures.
+
+    Parameters
+    ----------
+    card : list or 1D numpy.ndarray
+        Timeseries of recorded cardiac signal
+    peaks : list or 1D numpy.ndarray
+        array of peak indexes for card.
+    samplerate : float
+        Sampling rate for card, in Hertz.
+    window : float, optional
+        Size of the sliding window, in seconds.
+        Default is 6.
+    central_measure : "mean","average", "avg", "median", "mdn",  string, optional
+        Measure of the center used (mean or median).
+        Default is "mean".
+    Returns
+    -------
+    card_met : 2D numpy.ndarray
+        Heart Beats Interval or Heart Rate Variability timeseries.
+        The first column is the raw metric, in Hertz.
+        The second column is the metric convolved with the CRF, cut to the length
+        of the raw metric.
+
+    Notes
+    -----
+    Heart rate variability (HRV) is taken from [1]_, and computed as the amounts of
+    beats per minute.
+    However, operationally, it is the average of the inverse of the time interval
+    between two heart beats.
+    This metric should be convolved with the cardiac response function
+    before being included in a GLM.
+
+    IMPORTANT : The unit of measure has a meaning, since they it's based on Hertz.
+    Hence, zscoring might remove important quantifiable information.
+
+    See `_cardiac_metrics` for full implementation.
+
+    References
+    ----------
+    .. [1] N. Pinheiro et al., "Can PPG be used for HRV analysis?," 2016 38th
+    Annual International Conference of the IEEE Engineering in Medicine and
+    Biology Society (EMBC), doi: 10.1109/EMBC.2016.7591347.
+    """
+    return _cardiac_metrics(
+        card, peaks, samplerate, metric="hrv", window=6, central_measure="std"
+    )
+
+
 @due.dcite(references.CHEN_2020)
 def heart_beat_interval(card, peaks, samplerate, window=6, central_measure="mean"):
     """

phys2denoise/references.py

@@ -11,6 +11,8 @@
 
 GLOVER_2000 = Doi("10.1002/1522-2594(200007)44:1<162::AID-MRM23>3.0.CO;2-E")
 
+PINHERO_ET_AL_2016 = Doi("10.1109/EMBC.2016.7591347")
+
 POWER_2018 = Doi("10.1073/pnas.1720985115")
 
 POWER_2020 = Doi("10.1016/j.neuroimage.2019.116234")