This repository provides the implementation, data processing scripts, and analysis tools for:
The codebase enables full reproduction of all main figures and tables in the paper, including CRC phase synchronization, BPRSA coupling analysis, and rPPG-based non-contact measurement.
- Required Data: The raw data is from SUMS. Batch processing is performed to extract relevant features and generate the data used for reproducing results. Place the processed input files in the
data/directory (see Repository Structure for details).
- Reproduce Results: Run
scripts/draw_update.ipynbto generate all main figures and metrics for the paper.
├── data/ # Data for reproducing results
├── scripts/ # Analysis and plotting scripts
├── figures/ # Output figures and images
├── synchrogram.py # Synchronization Degree Calculation
├── README.md # Project documentation
└── BPRSA.py # Pulse Wave–Respiration Coupling
Corresponds to Figure 2, Figure 3, Table 1 in the paper.
- Find the phase values of the slower wave at the time points where the faster wave peaks occur. When two signals are n:m phase synchronized, you will see n nearly horizontal lines in the m-normalized relative phase plot.
- Example of a synchrogram:
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import hilbert, find_peaks
# Parameter settings
omega_x = 3
omega_y = 5
t = np.arange(0, 50, 0.01) # Time vector
x = np.sin(omega_x * np.pi * t) # Slower wave (x)
y = np.sin(omega_y * np.pi * t) # Faster wave (y)
# Calculate the phase of the slower wave (x) using Hilbert transform
phi_x = np.angle(hilbert(x))
new_phi_x = adjust_phase_to_0_mpi(phi_x) # Adjust phase to [0, 2π] range
phi_x = np.mod(new_phi_x, omega_x * 2 * np.pi) # Adjust phase to [0, 2π] range
# Find the peak positions of the faster wave (y)
peaks, _ = find_peaks(y)
# Get the phase values of the slower wave at the peaks of the faster wave
sync = phi_x[peaks]
sync_1 = new_phi_x[peaks] # Adjusted phase values
# Plotting
fig, axes = plt.subplots(2, 1, figsize=(10, 6))
# Plot the signals
axes[0].plot(t, x, 'b', label='x (slow)')
axes[0].plot(t, y, 'r', label='y (fast)')
axes[0].legend()
# Plot the phase synchrogram
axes[1].plot(peaks * 0.01, sync / (2 * np.pi), 'o') # peaks are indices, multiply by 0.01 to convert to time
axes[1].set_title(r'$\mathrm{\psi}_2(t_k)$', fontsize=13)
plt.tight_layout()
plt.show()- Refer to Cardiorespiratory Phase Synchronization in Elderly Patients with Periodic and non-Periodic Breathing Patterns
- To calculate the synchronization degree of a segment of cardiac and respiratory signals, first obtain the phase of the respiratory signal and the timing of heartbeats.
- Normalize the relative phase: calculate the respiratory phase at heartbeat times, constrain it to the 0~2πm interval, and then normalize by multiplying by 1/(2π).
- Merge the n horizontal lines into a single line, use a sliding window, and calculate the synchronization degree within the window using a method similar to circular variance.
-
Use the Hilbert transform to obtain the respiratory phase (range: -π to π), perform phase unwrapping, and then use the mod operation to constrain the phase within [0, 2πm] to observe heartbeats within m respiratory cycles.
-
Signal processing: apply bandpass filtering to the signals.
-
Use a 20-second window to calculate and record the initial guess of n:1 for each heartbeat.
-
For each point (heartbeat time), use a window of size N=50 to calculate the synchronization degree for all possible n:m ratios, select the maximum synchronization degree as the value for that point, and record the corresponding synchronization ratio.
-
Determine the synchronization threshold and calculate the relevant evaluation metrics:
-
%Sync: the percentage of synchronization defined as the sum of all the synchronization periods over the recording.
-
AvgDurSync: the average duration of the synchronization as the arithmetic mean of all periods of synchronization over the duration of the segment.
-
NumSync: the number of synchronized epochs defined as the sum of synchronization periods over the duration of the segment.
-
FreqRat: the frequency ratio between the respiratory and cardiac oscillators, as an approximation to the number of heartbeats per each respiration.
- Refer to the method from the paper Mechanical ventilatory modes and cardioventilatory phase synchronization in acute respiratory failure patients: First, determine an initial n1:1 ratio, which can be obtained by the ratio of heart rate to respiratory rate. For example, if the ratio is 3.27:1, the first guess for n:m is 3:1.
- When m = 1, the process selects the two ratios closest to the initial guess, i.e., (n1−1):1 and (n1+1):1. These are called the leftmost and rightmost ratios. For example, when m = 1, the leftmost and rightmost ratios are 2:1 and 4:1. Then, 2:1 and 4:1 propagate to m = 2 as 4:2 and 8:2, and to m = 3 as 6:3 and 12:3. For a given m, all propagated ratios (including the extremes) are accepted as candidate values.
Corresponds to Figure 4, Table 2 in the paper.
- Refer to the BPRSA method used in the paper Cardiorespiratory coupling in preterm infants:
- PRSA (Phase-Rectified Signal Averaging) is a signal analysis technique that can detect and quantify quasi-periodic oscillations masked by non-stationary characteristics, artifacts, and noise in composite signals. PRSA quantifies the coherent time of each quasi-periodicity when the signal is at a specific phase. Phase identification is based on detecting anchor points (AP), and averaging the signal around these APs yields a compressed, phase-aligned version of the original signal. This PRSA signal contains information about the physiological mechanisms (possibly quasi-periodic) that generate the corresponding APs. The following AP criteria are used to identify the rising edge of the pulse wave; replacing '>' with '<' identifies the falling edge.
- BPRSA (Bivariate PRSA) is an extension of the PRSA technique, where APs are defined in one (trigger) signal and averaging is performed in another (target) signal, enabling the identification and study of coupling from the trigger to the target signal.
- In this context, the pulse wave is used as the trigger signal and the respiratory wave as the target signal to examine the coupling from pulse wave to respiration.
-
MRA (Maximum Respiratory Amplitude):
The peak-to-peak amplitude difference of the respiratory BPRSA curve within the window, reflecting the maximum tidal volume change. -
SAP (Slope at Anchor Point):
The instantaneous slope of the respiratory curve at the anchor point; a positive value corresponds to inspiration, and a negative value corresponds to expiration.
Title: Contact Sensors to Remote Cameras: Quantifying Cardiorespiratory Coupling in High-Altitude Exercise Recovery
Jiankai Tang, Meng Kang, Yiru Zhang, Kegang Wang, Daniel Mcduff, Xin Liu, Yuanchun Shi, Yuntao Wang, "Contact Sensors to Remote Cameras: Quantifying Cardiorespiratory Coupling in High-Altitude Exercise Recovery", UbiComp Companion, 2025
@article{tang2025contact,
title={Contact Sensors to Remote Cameras: Quantifying Cardiorespiratory Coupling in High-Altitude Exercise Recovery},
author={Jiankai Tang and Meng Kang and Yiru Zhang and Kegang Wang and Daniel Mcduff and Xin Liu and Yuanchun Shi and Yuntao Wang},
year={2025},
eprint={2508.00773},
archivePrefix={arXiv},
primaryClass={cs.CE},
url={https://arxiv.org/abs/2508.00773},
}