physiopy / prep4phys

A toolbox for physiological peak detection analyses. Formerly peakdet.
Apache License 2.0
1 stars 1 forks source link

Interpolation of physio objects after peak detection does not interpolate peaks/troughs #3

Open smoia opened 6 months ago

smoia commented 6 months ago

Quick check inspired by PR#5 in physioqc. When resampling a physio object post peak detection, the indexes of peaks and troughs are not interpolated.

Expected Behavior

Say p is a physio object with fs=4, on which peak detection was performed. g is the downsampled version of p to fs=2, h is the upsampled version of p to fs=8.

Either:

In [25]: len(p)
Out[25]: 100

In [26]: p.peaks
Out[26]: array([14, 17, 29, 31, 38, 40, 53, 56, 58, 65, 81, 87, 91, 96])

In [27]: len(g)
Out[27]: 50

In [28]: g.peaks
Out[28]: array([7, 9, 15, 16, 19, 20, 27, 28, 29, 33, 41, 44, 46, 48])
# Eventually, cases like that 15-16, 19-20, 27-28-29 should be solved by running a quick local maxima detection - in fact the maxima detection should be run on all neighbours of the index, especially when the modulo of the partition is not 0.

In [29]: len(h)
Out[29]: 200

In [30]: h.peaks
Out[30]: array([28,  34,  58,  62,  76,  80, 106, 112, 116, 130, 162, 174, 182, 192])

Or raise a warning on resampling (specifically downsampling) that peak data will be removed.

Actual Behavior

In [25]: len(p)
Out[25]: 100

In [26]: p.peaks
Out[26]: array([14, 17, 29, 31, 38, 40, 53, 56, 58, 65, 81, 87, 91, 96])

In [27]: len(g)
Out[27]: 50

In [28]: g.peaks
Out[28]: array([14, 17, 29, 31, 38, 40, 53, 56, 58, 65, 81, 87, 91, 96])

In [29]: len(h)
Out[29]: 200

In [30]: h.peaks
Out[30]: array([14, 17, 29, 31, 38, 40, 53, 56, 58, 65, 81, 87, 91, 96])

(same for troughs btw).

Steps to Reproduce the Problem

import numpy as np
import peakdet as pk
import matplotlib.pyplot as plt
a = np.random.randn(100)
p = pk.Physio(a, fs=4)
pk.plot_physio(p)
plt.show()
p = pk.peakfind_physio(p)
pk.plot_physio(p)
plt.show()
g = pk.interpolate_physio(p, 2)
h = pk.interpolate_physio(p, 8)
plt.show(g)
pk.plot_physio(g)
len(p)
p.peaks
len(g)
g.peaks
len(h)
h.peaks

Specifications

- Python version: 3.8.10
- peakdet version: 0.3.0
- Platform: Ubuntu

Possible solution

Best:

  1. compute the new indexes of peaks and troughs based on the resampling ratio (target fs / current fs)
  2. understand if rounding the computed indexes (up/down/both) constantly hits the new true peaks after resampling
  3. on downsampling, run local maxima detection (two neighbouring indexes, for a total of 2-3) especially when peak indexes are clustered together (see expected behaviour output)
  4. on upsampling, if (2) does not hit the right index all the time, implement a maxima detection on neighbouring indexes (5 adjacent or based on resampling ratio)

Less good but fair: remove peak and trough metadata on down-/re-sampling, raise a warning, invite to run peak detection again (and remind that the history can help repeat steps so far).