5  RSP Algorithms Explained

A guide to the respiratory signal processing algorithms in Ripple, explaining each step in the analysis pipeline and how to use the Clojure implementation.

Ripple’s RSP module processes signals from respiration belts, thoracic impedance sensors, and similar devices. The pipeline extracts breathing rate, amplitude, phase, symmetry, and variability metrics.

References:

(ns ripple-book.rsp-algorithms-explained
  (:require
   ;; Ripple RSP public API:
   [scicloj.ripple.rsp :as rsp]
   ;; Table processing:
   [tablecloth.api :as tc]
   ;; High-performance array math:
   [tech.v3.datatype :as dtype]
   [tech.v3.datatype.functional :as dfn]
   ;; Interactive plotting via Plotly:
   [scicloj.tableplot.v1.plotly :as plotly]
   ;; Annotating kinds of visualizations:
   [scicloj.kindly.v4.kind :as kind]))

The Respiratory Signal Processing Pipeline

Respiration is a quasi-periodic signal. Inhalation produces a rising deflection (trough → peak) and exhalation a falling deflection (peak → trough). The processing pipeline extracts cycle-by-cycle features:

Step Function Purpose
1 Bandpass filtering Remove drift and high-frequency noise
2 Peak/trough detection Locate exhalation and inhalation onsets
3 Amplitude Measure depth of each breath
4 Phase Label inspiration vs expiration
5 Symmetry Characterize waveform shape
6 Breathing rate Instantaneous rate from inter-breath intervals
7 RRV Respiratory rate variability (time-domain)
8 RVT Respiratory volume per time

Synthetic RSP Signal

We start with a clean sinusoidal signal at 15 breaths per minute, sampled at 1000 Hz for 60 seconds. This gives us a known ground truth for validating each processing step.

(def sampling-rate 1000.0)
(def duration 60.0)
(def respiratory-rate 15.0)
(def rsp-clean-signal
  (rsp/synthetic-rsp {:sampling-rate sampling-rate
                      :duration duration
                      :respiratory-rate respiratory-rate
                      :seed 42}))
(def n-samples (count rsp-clean-signal))
n-samples
60000

Let’s also create a noisy version to demonstrate cleaning.

(def rsp-noisy
  (rsp/synthetic-rsp {:sampling-rate sampling-rate
                      :duration duration
                      :respiratory-rate respiratory-rate
                      :noise 0.3
                      :seed 42}))

Visualize both signals (first 10 seconds):

(let [n-show (int (* 10 sampling-rate))
      ds (tc/dataset {:t (dfn// (range n-show) sampling-rate)
                      :clean (dtype/sub-buffer rsp-clean-signal 0 n-show)
                      :noisy (dtype/sub-buffer rsp-noisy 0 n-show)})]
  (-> (tc/pivot->longer ds [:clean :noisy]
                        {:target-columns :series
                         :value-column-name :amplitude})
      (plotly/base {:=x :t
                    :=y :amplitude
                    :=color :series
                    :=title "Synthetic RSP — clean vs noisy"
                    :=x-title "Time (s)"
                    :=y-title "Amplitude"})
      (plotly/layer-line)
      plotly/plot))

Bandpass Filtering (Cleaning)

Raw respiratory signals contain baseline drift (< 0.05 Hz) and high-frequency noise. The Khodadad2018 method applies a 2nd-order Butterworth bandpass filter at 0.05–3 Hz, preserving the breathing band while removing artifacts.

(def rsp-cleaned
  (rsp/rsp-clean rsp-noisy sampling-rate {}))

Compare noisy vs cleaned (first 10 seconds):

(let [n-show (int (* 10 sampling-rate))
      ds (tc/dataset {:t (dfn// (range n-show) sampling-rate)
                      :noisy (dtype/sub-buffer rsp-noisy 0 n-show)
                      :cleaned (dtype/sub-buffer rsp-cleaned 0 n-show)})]
  (-> (tc/pivot->longer ds [:noisy :cleaned]
                        {:target-columns :series
                         :value-column-name :amplitude})
      (plotly/base {:=x :t
                    :=y :amplitude
                    :=color :series
                    :=title "RSP Cleaning — Khodadad2018 bandpass"
                    :=x-title "Time (s)"
                    :=y-title "Amplitude"})
      (plotly/layer-line)
      plotly/plot))

The cleaned signal closely tracks the original sinusoid. Let’s verify the correlation is high.

(let [n-check (int (* 50 sampling-rate))
      ;; Skip first 2s to avoid filter transient
      skip (int (* 2 sampling-rate))
      clean-seg (dtype/sub-buffer rsp-clean-signal skip n-check)
      cleaned-seg (dtype/sub-buffer rsp-cleaned skip n-check)
      ;; Normalize both to zero-mean, unit-variance for correlation
      norm (fn [s]
             (let [m (dfn/mean s)]
               (dfn/- s m)))
      a (norm clean-seg)
      b (norm cleaned-seg)
      corr (/ (double (dfn/sum (dfn/* a b)))
              (Math/sqrt (* (double (dfn/sum (dfn/* a a)))
                            (double (dfn/sum (dfn/* b b))))))]
  corr)
0.9844148948407483

Peak and Trough Detection

Peaks correspond to exhalation onsets (maximum lung volume) and troughs to inhalation onsets (minimum lung volume). The detector finds local extrema with a minimum distance constraint and removes amplitude outliers.

(def extrema
  (rsp/rsp-findpeaks rsp-cleaned sampling-rate {}))
(def peaks (:peaks extrema))
(def troughs (:troughs extrema))

We expect approximately respiratory-rate × duration / 60 cycles.

(count peaks)
18
(count troughs)
15

Visualize peaks and troughs on the signal (first 15 seconds):

(let [n-show (int (* 15 sampling-rate))
      ds-line (tc/dataset {:t (dfn// (range n-show) sampling-rate)
                           :amplitude (dtype/sub-buffer rsp-cleaned 0 n-show)})
      peak-in-range (filterv #(< % n-show) peaks)
      trough-in-range (filterv #(< % n-show) troughs)
      ds-peaks (tc/dataset {:t (dfn// peak-in-range sampling-rate)
                            :amplitude (mapv #(aget rsp-cleaned %) peak-in-range)
                            :type (repeat (count peak-in-range) "peak")})
      ds-troughs (tc/dataset {:t (dfn// trough-in-range sampling-rate)
                              :amplitude (mapv #(aget rsp-cleaned %) trough-in-range)
                              :type (repeat (count trough-in-range) "trough")})]
  (-> ds-line
      (plotly/base {:=x :t :=y :amplitude
                    :=title "Peak and trough detection"
                    :=x-title "Time (s)" :=y-title "Amplitude"})
      (plotly/layer-line {:=name "RSP"})
      (plotly/layer-point {:=dataset ds-peaks
                           :=x :t :=y :amplitude
                           :=mark-color "red" :=mark-size 8
                           :=name "peaks"})
      (plotly/layer-point {:=dataset ds-troughs
                           :=x :t :=y :amplitude
                           :=mark-color "blue" :=mark-size 8
                           :=name "troughs"})
      plotly/plot))

Respiratory Amplitude

Amplitude measures the depth of each breath — the difference between peak (maximum inflation) and the preceding trough (minimum). The per-cycle values are interpolated to a continuous signal using Akima monotone cubic spline.

Two methods:

  • Standard: peak minus preceding trough
  • Prepost: peak minus average of preceding and succeeding troughs
(def amplitude-standard
  (rsp/rsp-amplitude rsp-cleaned peaks troughs {}))
(def amplitude-prepost
  (rsp/rsp-amplitude rsp-cleaned peaks troughs {:method :prepost}))

For a pure sinusoid, amplitude should be close to 2.0 (peak-to-trough distance of a unit sine wave). The bandpass filter slightly attenuates, so we check it’s in [1.5, 2.5].

(let [mid (int (/ n-samples 2))]
  (aget amplitude-standard mid))
1.0766126954651902

Visualize (first 20 seconds):

(let [n-show (int (* 20 sampling-rate))
      ds (tc/dataset {:t (dfn// (range n-show) sampling-rate)
                      :standard (dtype/sub-buffer amplitude-standard 0 n-show)
                      :prepost (dtype/sub-buffer amplitude-prepost 0 n-show)})]
  (-> (tc/pivot->longer ds [:standard :prepost]
                        {:target-columns :method
                         :value-column-name :amplitude})
      (plotly/base {:=x :t :=y :amplitude :=color :method
                    :=title "Respiratory amplitude"
                    :=x-title "Time (s)" :=y-title "Amplitude"})
      (plotly/layer-line)
      plotly/plot))

Respiratory Phase

Phase labels each sample as either inspiration (1) or expiration (0). Phase completion runs from 0.0 to 1.0 within each phase segment.

(def phase-result
  (rsp/rsp-phase peaks troughs n-samples))
(def phase (:phase phase-result))
(def phase-completion (:phase-completion phase-result))

For a symmetric sinusoid, roughly 50% of samples should be in each phase.

(let [n-insp (long (reduce + (map #(aget phase %) (range n-samples))))
      pct (/ (double n-insp) n-samples 0.01)]
  pct)
41.44333333333333

Visualize phase and completion (first 10 seconds):

(let [n-show (int (* 10 sampling-rate))
      ds (tc/dataset {:t (dfn// (range n-show) sampling-rate)
                      :signal (dtype/sub-buffer rsp-cleaned 0 n-show)
                      :phase (dtype/sub-buffer phase 0 n-show)
                      :completion (dtype/sub-buffer phase-completion 0 n-show)})]
  (-> ds
      (plotly/base {:=x :t :=y :signal
                    :=title "RSP phase and completion"
                    :=x-title "Time (s)"})
      (plotly/layer-line {:=name "signal" :=y :signal})
      (plotly/layer-line {:=name "phase" :=y :phase
                          :=mark-color "red"})
      (plotly/layer-line {:=name "completion" :=y :completion
                          :=mark-color "green"})
      plotly/plot))

Respiratory Symmetry

Peak-trough symmetry characterizes the waveform shape of each breath cycle (Cole & Voytek, 2019): time_to_peak / cycle_duration. Values near 0.5 indicate symmetric breathing; lower values indicate faster inhalation than exhalation.

Cole & Voytek also define a rise-decay symmetry based on derivative zero-crossings, but for roughly sinusoidal respiratory signals the derivative zero-crossing coincides with the peak, making it redundant with peak-trough. We omit it here.

(def symmetry
  (rsp/rsp-symmetry rsp-cleaned peaks troughs))
(let [mid (int (/ n-samples 2))]
  (aget (:peak-trough symmetry) mid))
0.4127364746318107

Breathing Rate

Instantaneous breathing rate is computed from trough-to-trough intervals, converted to BPM, and interpolated to a continuous signal.

(def breathing-rate
  (rsp/rsp-rate troughs sampling-rate n-samples))

For our 15 BPM synthetic signal, the rate should be close to 15.

(let [mid (int (/ n-samples 2))]
  (aget breathing-rate mid))
15.138827628651578

Visualize (full signal):

(let [;; Downsample for plotting (every 100th sample)
      step 100
      indices (range 0 n-samples step)
      ds (tc/dataset {:t (mapv #(/ (double %) sampling-rate) indices)
                      :bpm (mapv #(aget breathing-rate %) indices)})]
  (-> ds
      (plotly/base {:=x :t :=y :bpm
                    :=title "Instantaneous breathing rate"
                    :=x-title "Time (s)" :=y-title "Breaths per minute"})
      (plotly/layer-line)
      plotly/plot))

Breath-to-Breath Intervals and RRV

Breath-to-breath (BB) intervals are the respiratory analog of RR intervals in ECG. Respiratory Rate Variability (RRV) quantifies the variation in breathing timing using the same time-domain metrics as HRV.

(def bb-intervals
  (rsp/peaks->bb-intervals troughs sampling-rate))

At 15 BPM, each cycle is 4000 ms.

(let [mid (int (/ (count bb-intervals) 2))]
  (aget bb-intervals mid))
4036.0

Compute RRV indices.

(def rrv (rsp/rrv-time bb-intervals))
rrv
{:mean-br 14.973948048292264,
 :cvsd 0.03225021044443647,
 :sd-br 0.28319942324707104,
 :mad-bb 68.9409,
 :sdsd 134.50641010366144,
 :median-bb 4002.0,
 :min-bb 3855.0,
 :cvbb 0.018848383692611756,
 :max-bb 4136.0,
 :sdbb 75.54970709247152,
 :max-br 15.56420233463035,
 :min-br 14.506769825918763,
 :rmssd 129.26805780714264,
 :mean-bb 4008.285714285714}

For a clean sinusoid with constant rate, variability should be very low.

(:sdbb rrv)
75.54970709247152
(:mean-br rrv)
14.973948048292264

Respiratory Volume per Time (RVT)

RVT is the ratio of respiratory amplitude to inter-peak interval. It captures the overall ventilatory drive and is used as an fMRI confound regressor (Power et al. 2020).

(def rvt
  (rsp/rsp-rvt rsp-cleaned peaks troughs sampling-rate))

RVT should be positive and relatively stable for a constant-rate signal.

(let [mid (int (/ n-samples 2))]
  (aget rvt mid))
0.8831688759174057

Visualize:

(let [step 100
      indices (range 0 n-samples step)
      ds (tc/dataset {:t (mapv #(/ (double %) sampling-rate) indices)
                      :rvt (mapv #(aget rvt %) indices)})]
  (-> ds
      (plotly/base {:=x :t :=y :rvt
                    :=title "Respiratory Volume per Time"
                    :=x-title "Time (s)" :=y-title "RVT"})
      (plotly/layer-line)
      plotly/plot))

Summary

Ripple’s RSP module provides a complete pipeline from raw respiratory signals to derived physiological metrics:

Function Output
rsp/synthetic-rsp Test signal generation
rsp/rsp-clean Bandpass-filtered signal
rsp/rsp-findpeaks Peak and trough indices
rsp/rsp-amplitude Continuous amplitude signal
rsp/rsp-phase Phase labels and completion
rsp/rsp-symmetry Cycle shape features
rsp/rsp-rate Instantaneous breathing rate
rsp/peaks->bb-intervals Inter-breath intervals
rsp/rrv-time Time-domain variability indices
rsp/rsp-rvt Respiratory volume per time

All functions work with double-array inputs and outputs, following the same conventions as the cardio module.

source: notebooks/ripple_book/rsp_algorithms_explained.clj