Documentation of rsp_clean() seems inconsistent with implementation

[danibene]

Hello! I am a bit confused about the default method for cleaning respiration signals.

In the docstring, it is written that there is linear detrending followed by a 5th order lowpass Butterworth filter with a cutoff of 2Hz:

NeuroKit/neurokit2/rsp/rsp_clean.py

Lines 18 to 19 in 3d004b4

	* **khodadad2018**: Linear detrending followed by a fifth order 2Hz low-pass IIR Butterworth
	filter)

But the implementation seems to just be a 2nd order bandpass Butterworth filter with cutoffs of 0.5 Hz and 3 Hz:

NeuroKit/neurokit2/rsp/rsp_clean.py

Lines 121 to 143 in 3d004b4

	def _rsp_clean_khodadad2018(rsp_signal, sampling_rate=1000):
	"""The algorithm is based on (but not an exact implementation of) the "Zero-crossing algorithm with amplitude
	threshold" by `Khodadad et al. (2018)
	<https://iopscience.iop.org/article/10.1088/1361-6579/aad7e6/meta>`_.
	"""
	# Slow baseline drifts / fluctuations must be removed from the raw
	# breathing signal (i.e., the signal must be centered around zero) in order
	# to be able to reliable detect zero-crossings.
	# Remove baseline by applying a lowcut at .05Hz (preserves breathing rates
	# higher than 3 breath per minute) and high frequency noise by applying a
	# highcut at 3 Hz (preserves breathing rates slower than 180 breath per
	# minute).
	clean = signal_filter(
	rsp_signal,
	sampling_rate=sampling_rate,
	lowcut=0.05,
	highcut=3,
	order=2,
	method="butterworth",
	)

Also, in the cited paper, Khodadad et al. (2018) wrote "However, further improvement is also obtained by pre-processing the data using a digital high-pass filter to remove the dominating low frequency contents. Here, a second order high-pass Butterworth filter has been used with a cut-off frequency of 15 breaths/min." Wouldn't that correspond to 0.25 Hz?

Should we update the documentation to reflect that it is a 2nd order bandpass Butterworth filter with cutoffs of 0.5 Hz and 3 Hz?

[DominiqueMakowski]

I think we should adjust the method to match the paper if there's a discrepancy no?

[danibene]

Seems like the discrepancy between the docstring and the implementation started here:

d947009

Any idea why that was done? Maybe another option could be to have a default "neurokit" method (in case this implementation does work better) as well as the original implementation based on the one from the paper

[DominiqueMakowski]

Any idea why that was done?

indeed that's strange and unfortunately the committer disappeared from GH it seems.
I think the idea was to replace the detrending by the highpass, but...

in case this implementation does work better

The problem is how to benchmark that, I ran a quick search for annotated RSP data just to have some empirical evidence and found this dataset:

https://www.physionet.org/content/bidmc/1.0.0/

That means we could maybe try to compare the time of peaks&troughts after various cleaning procedures vs. their annotations?

And then, if the current version is better, we move it to a new neurokit default method. If not, we drop it. What do you think?

---

Also they describe another method here based on Lu2006 that I don't think we implement:

https://physiodatatoolbox.leidenuniv.nl/docs/user-guide/physioanalyzer-modules/resp-module.html

[danibene]

Making the decision based on its performance on that dataset sounds good! I just don't know when I would be able to do that (though you or any kind stranger lurking here is welcome to), so my proposal would be:

1. PR changing the docstring to reflect the current implementation, mentioning this issue & that it currently doesn't reflect the implementation in the paper (so that the default functionality isn't changed before we check its performance)

2. PR downloading and processing the open-access dataset

3. PR adding study comparing methods

4. PR changing or adding implementation depending on results

What do you think?

[DominiqueMakowski]

sounds good! 👌