PhysioCrossModal provides a comprehensive toolkit for analyzing coupling and synchrony between physiological signals from different modalities (e.g., EEG, EMG, ECG, EDA, MoCap). It extends the PhysioCore ecosystem with:
MultiPhysioExperiment container for multi-modal
recordings at different sampling ratesThis vignette walks through a complete cross-modal coupling analysis workflow.
PhysioCrossModal includes three data generators for quick prototyping and testing. We start by creating a simulated EEG-EMG dataset with corticomuscular coherence (CMC) at 20 Hz.
When you have signals from multiple modalities recorded
simultaneously, the MultiPhysioExperiment class provides a
unified container. Each modality is stored as a named
PhysioExperiment with its own sampling rate.
# Wrap modalities into a MultiPhysioExperiment
mpe <- MultiPhysioExperiment(EEG = pe_eeg, EMG = pe_emg)
mpe
# Inspect
modalities(mpe) # c("EEG", "EMG")
samplingRates(mpe) # EEG=500, EMG=1000
nModalities(mpe) # 2
# Extract a single modality
eeg <- mpe[["EEG"]]When coupling functions receive signals at different sampling rates, they automatically resample to a common rate. You can also align signals explicitly:
Coherence measures the linear frequency-domain relationship between two signals. Values range from 0 (no linear relationship) to 1 (perfect linear relationship) at each frequency.
# Coherence between first EEG and first EMG channels
coh <- coherence(pe_eeg, pe_emg,
channels_x = 1, channels_y = 1,
nperseg = 256L)
# Peak coherence and its frequency
peak_idx <- which.max(coh$coherence)
cat("Peak coherence:", round(coh$coherence[peak_idx], 3),
"at", round(coh$frequencies[peak_idx], 1), "Hz\n")
cat("95% confidence limit:", round(coh$confidence_limit, 3), "\n")For short signals or when frequency resolution is critical, multitaper estimation provides lower variance than Welch’s method:
The cross-spectrum captures both magnitude and phase relationships:
Phase-based measures quantify the consistency of the phase relationship between signals, independent of amplitude.
PLV measures phase consistency over time. It requires specifying a frequency band for bandpass filtering:
PLI is insensitive to zero-lag (volume conduction) effects:
wPLI reduces noise sensitivity by weighting phase differences by their cross-spectral magnitude:
Granger causality tests whether one signal’s past helps predict another signal beyond its own history, indicating directional information flow.
# Create directed signals: x drives y with a 10-sample lag
dir_data <- make_directed_signals(n = 5000, sr = 500,
lag_samples = 10, coupling = 0.7)
gc <- grangerCausality(dir_data$x, dir_data$y, sr = 500, order = 15)
cat("GC x->y:", round(gc$gc_xy, 4), "\n")
cat("GC y->x:", round(gc$gc_yx, 4), "\n")
cat("Net GC:", round(gc$net_gc, 4), "(positive = x drives y)\n")Cross-correlation identifies time delays between signals:
Wavelet methods reveal how coupling varies jointly in time and frequency.
signals <- make_coupled_signals(sr1 = 200, sr2 = 200,
coupling_freq = 10, duration = 5)
wcoh <- waveletCoherence(signals$x, signals$y, sr = 200,
frequencies = seq(2, 30),
n_cycles = 7)
# Visualize time-frequency coherence
if (requireNamespace("ggplot2", quietly = TRUE)) {
plotWaveletCoherence(wcoh)
}The couplingAnalysis() wrapper dispatches to any
coupling method with a uniform API:
# Same call structure for any method
res_coh <- couplingAnalysis(pe_eeg, pe_emg,
method = "coherence",
channels_x = 1, channels_y = 1,
nperseg = 256L)
res_plv <- couplingAnalysis(pe_eeg, pe_emg,
method = "plv",
channels_x = 1, channels_y = 1,
freq_band = c(18, 22))
res_gc <- couplingAnalysis(dir_data$x, dir_data$y,
method = "granger",
sr = 500, order = 15)
# Also works directly with MultiPhysioExperiment
res_mpe <- couplingAnalysis(mpe,
modality_x = "EEG", modality_y = "EMG",
method = "coherence",
channels_x = 1, channels_y = 1)Compute coupling between all pairs of channels across two modalities:
# Coherence matrix
coh_mat <- coherenceMatrix(pe_eeg, pe_emg, nperseg = 128L)
coh_mat$matrix
# Generic coupling matrix (any method)
pli_mat <- couplingMatrix(pe_eeg, pe_emg,
method = "coherence",
nperseg = 128L)
# Visualize as heatmap
if (requireNamespace("ggplot2", quietly = TRUE)) {
plotCouplingMatrix(coh_mat$matrix, title = "EEG-EMG Coherence")
}Test whether observed coupling exceeds what would be expected by chance, using phase-randomized or time-shifted surrogates:
sig_test <- surrogateTest(pe_eeg, pe_emg,
method = "coherence",
channels_x = 1, channels_y = 1,
n_surrogates = 199,
surrogate_type = "phase",
nperseg = 128L)
cat("Observed statistic:", round(sig_test$statistic, 3), "\n")
cat("p-value:", round(sig_test$p_value, 4), "\n")
cat("95th percentile threshold:", round(sig_test$threshold_95, 3), "\n")Estimate the precision of the coupling statistic using block bootstrap:
Test significance for every channel pair simultaneously, with FDR or Bonferroni correction:
Putting it all together – a minimal script for analyzing EEG-EMG corticomuscular coherence:
library(PhysioCrossModal)
# 1. Create or load data
data <- make_eeg_emg(n_sec = 10, sr_eeg = 500, sr_emg = 1000)
# 2. Wrap into MultiPhysioExperiment
mpe <- MultiPhysioExperiment(EEG = data$eeg, EMG = data$emg)
# 3. Compute coherence matrix
coh_mat <- coherenceMatrix(mpe,
modality_x = "EEG", modality_y = "EMG",
nperseg = 256L,
freq_range = c(15, 25))
# 4. Test significance with FDR correction
sig <- surrogateMatrixTest(mpe,
modality_x = "EEG", modality_y = "EMG",
method = "coherence",
n_surrogates = 199,
correction = "fdr",
nperseg = 256L,
freq_range = c(15, 25))
# 5. Visualize
if (requireNamespace("ggplot2", quietly = TRUE)) {
plotCouplingMatrix(sig$matrix, title = "CMC (15-25 Hz)")
}
# 6. Detailed spectrum for the most coupled pair
best <- which(sig$matrix == max(sig$matrix, na.rm = TRUE), arr.ind = TRUE)
coh_detail <- coherence(mpe,
modality_x = "EEG", modality_y = "EMG",
channels_x = best[1, 1],
channels_y = best[1, 2])
if (requireNamespace("ggplot2", quietly = TRUE)) {
plotCoherenceSpectrum(coh_detail)
}