Physiological Features

This notebook demonstrates extraction of frequency-domain, nonlinear, and entropy-based physiological features from ECG and PPG signals using vitalDSP.

Setup

import numpy as np
import plotly.io as pio
pio.renderers.default = "sphinx_gallery"
from plotly import graph_objects as go
from vitalDSP.notebooks import load_sample_ecg_small, load_sample_ppg, plot_trace

fs = 128
signal_col, date_col = load_sample_ecg_small()
signal_col = np.array(signal_col)

ppg_col, ppg_date = load_sample_ppg()
ppg_col = np.array(ppg_col)

print(f"ECG signal length: {len(signal_col)} samples ({len(signal_col)/fs:.1f}s at {fs} Hz)")
print(f"PPG signal length: {len(ppg_col)} samples")

Nonlinear Features

from vitalDSP.physiological_features.nonlinear import NonlinearFeatures

nl = NonlinearFeatures(signal_col, fs=fs)

sample_entropy = nl.compute_sample_entropy(m=2, r=0.2)
app_entropy    = nl.compute_approximate_entropy(m=2, r=0.2)
dfa_alpha      = nl.compute_dfa()
lyapunov       = nl.compute_lyapunov_exponent()

print("Nonlinear Features:")
print(f"  Sample Entropy      : {sample_entropy:.4f}")
print(f"  Approximate Entropy : {app_entropy:.4f}")
print(f"  DFA Alpha           : {dfa_alpha:.4f}")
print(f"  Lyapunov Exponent   : {lyapunov:.4f}")

Poincaré Plot

poincare = nl.compute_poincare_features()
print("Poincaré Features:", poincare)

# Visualize Poincaré plot
fig = go.Figure()
fig.add_trace(go.Scatter(
    x=signal_col[:-1], y=signal_col[1:],
    mode='markers', marker=dict(size=3, opacity=0.5),
    name='Poincaré'
))
fig.update_layout(
    title='Poincaré Plot',
    xaxis_title='x(n)', yaxis_title='x(n+1)'
)
fig.show()

Symbolic Dynamics

from vitalDSP.physiological_features.symbolic_dynamics import SymbolicDynamics

sd = SymbolicDynamics(signal_col, n_symbols=4, word_length=3)

perm_entropy   = sd.compute_permutation_entropy()
shannon_entropy = sd.compute_shannon_entropy()
word_dist      = sd.compute_word_distribution()

print("Symbolic Dynamics:")
print(f"  Permutation Entropy : {perm_entropy:.4f}")
print(f"  Shannon Entropy     : {shannon_entropy:.4f}")

# Word distribution bar chart
words = list(word_dist.keys())
counts = list(word_dist.values())
fig = go.Figure(go.Bar(x=[str(w) for w in words], y=counts))
fig.update_layout(title='Symbolic Word Distribution', xaxis_title='Word', yaxis_title='Count')
fig.show()

Envelope Detection

from vitalDSP.physiological_features.envelope_detection import EnvelopeDetection

ed = EnvelopeDetection(signal_col)
hilbert_env = ed.hilbert_envelope()
rms_env     = ed.moving_average_envelope(window_size=32)
peak_env    = ed.peak_envelope()

fig = go.Figure()
fig.add_trace(go.Scatter(y=signal_col[:512], mode='lines', name='ECG Signal', opacity=0.4))
fig.add_trace(go.Scatter(y=hilbert_env[:512], mode='lines', name='Hilbert Envelope'))
fig.add_trace(go.Scatter(y=rms_env[:512], mode='lines', name='RMS Envelope'))
fig.add_trace(go.Scatter(y=peak_env[:512], mode='lines', name='Peak Envelope'))
fig.update_layout(title='Signal Envelopes', xaxis_title='Sample', yaxis_title='Amplitude')
fig.show()

Trend Analysis

from vitalDSP.physiological_features.trend_analysis import TrendAnalysis

ta = TrendAnalysis(signal_col)
ema_trend    = ta.compute_exponential_smoothing(alpha=0.1)
ma_trend     = ta.compute_moving_average(window_size=16)
linear_trend = ta.compute_linear_trend()

fig = go.Figure()
fig.add_trace(go.Scatter(y=signal_col[:512], mode='lines', name='ECG Signal', opacity=0.4))
fig.add_trace(go.Scatter(y=ema_trend[:512], mode='lines', name='EMA Trend'))
fig.add_trace(go.Scatter(y=ma_trend[:512], mode='lines', name='Moving Average'))
fig.add_trace(go.Scatter(y=linear_trend[:512], mode='lines', name='Linear Trend'))
fig.update_layout(title='Signal Trend Analysis', xaxis_title='Sample', yaxis_title='Amplitude')
fig.show()

Signal Change Detection

from vitalDSP.physiological_features.signal_change_detection import SignalChangeDetection

scd = SignalChangeDetection(signal_col)
zcr       = scd.zero_crossing_rate()
var_pts   = scd.variance_based_detection(window_size=64, threshold=2.0)

print(f"Zero Crossing Rate          : {zcr:.4f}")
print(f"Variance-based change points: {len(var_pts)}")

fig = go.Figure()
fig.add_trace(go.Scatter(y=signal_col[:512], mode='lines', name='ECG Signal'))
cp_in_range = [p for p in var_pts if p < 512]
if cp_in_range:
    fig.add_trace(go.Scatter(
        x=cp_in_range,
        y=[signal_col[p] for p in cp_in_range],
        mode='markers', marker=dict(color='red', size=8),
        name='Change Points'
    ))
fig.update_layout(title='Signal Change Detection', xaxis_title='Sample', yaxis_title='Amplitude')
fig.show()

Signal Segmentation

from vitalDSP.physiological_features.signal_segmentation import SignalSegmentation

seg = SignalSegmentation(signal_col)

# Fixed-size segments
fixed_segs = seg.fixed_size_segmentation(segment_size=128)
print(f"Fixed-size segments (128 samples): {len(fixed_segs)}")

# Variance-based adaptive segmentation
var_segs = seg.variance_based_segmentation(window_size=64, threshold_factor=1.5)
print(f"Variance-based segments: {len(var_segs)}")

Energy Analysis

from vitalDSP.physiological_features.energy_analysis import EnergyAnalysis

ea = EnergyAnalysis(signal_col, fs=fs)
total_energy  = ea.compute_total_energy()
band_energy   = ea.compute_band_energy(low_freq=0.5, high_freq=40.0)
spectral_energy = ea.compute_spectral_energy()

print("Energy Analysis:")
print(f"  Total Energy          : {total_energy:.4f}")
print(f"  Band Energy (0.5-40Hz): {band_energy:.4f}")
print(f"  Spectral Energy       : {spectral_energy:.4f}")

Signal Power Analysis

from vitalDSP.physiological_features.signal_power_analysis import SignalPowerAnalysis

spa = SignalPowerAnalysis(signal_col)
total_power = spa.compute_total_power()
rms_power   = spa.compute_mean_square()
snr         = spa.compute_snr()

print("Signal Power Analysis:")
print(f"  Total Power : {total_power:.4f}")
print(f"  RMS Power   : {rms_power:.4f}")
print(f"  SNR (dB)    : {snr:.4f}")

Nonlinear Features — Fractal Dimension and Recurrence

from vitalDSP.physiological_features.nonlinear import NonlinearFeatures

nl = NonlinearFeatures(signal_col, fs=fs)

fractal_dim = nl.compute_fractal_dimension()
recurrence  = nl.compute_recurrence_features(threshold=0.2, sample_size=200)

print(f'Fractal Dimension  : {fractal_dim:.4f}')
print('Recurrence Features:')
for k, v in recurrence.items():
    print(f'  {k}: {v:.4f}' if isinstance(v, float) else f'  {k}: {v}')

Symbolic Dynamics — Full Pipeline

from vitalDSP.physiological_features.symbolic_dynamics import SymbolicDynamics

sd = SymbolicDynamics(signal_col, n_symbols=4, word_length=3)

symbols      = sd.symbolize()
renyi        = sd.compute_renyi_entropy(alpha=2)
forbidden    = sd.detect_forbidden_words()
trans_matrix = sd.compute_transition_matrix()
all_features = sd.compute_symbolic_features()

print(f'Symbolized sequence (first 20): {symbols[:20]}')
print(f'Renyi Entropy (alpha=2)        : {renyi:.4f}')
print(f'Forbidden words count          : {len(forbidden)}')
print(f'Transition matrix shape        : {trans_matrix.shape}')
print('All symbolic features:', all_features)

fig = go.Figure(go.Heatmap(z=trans_matrix, colorscale='Blues', colorbar=dict(title='Prob')))
fig.update_layout(title='Symbol Transition Matrix',
                  xaxis_title='Next Symbol', yaxis_title='Current Symbol')
fig.show()

Envelope Detection — Extended Methods

from vitalDSP.physiological_features.envelope_detection import EnvelopeDetection

ed = EnvelopeDetection(signal_col)
abs_env     = ed.absolute_value_envelope(smoothing_factor=16)
wavelet_env = ed.wavelet_envelope(wavelet_name='db4', level=2)

fig = go.Figure()
fig.add_trace(go.Scatter(y=signal_col[:512], mode='lines', name='ECG Signal', opacity=0.4))
fig.add_trace(go.Scatter(y=abs_env[:512], mode='lines', name='Absolute Value Envelope'))
fig.add_trace(go.Scatter(y=wavelet_env[:512], mode='lines', name='Wavelet Envelope (db4)'))
fig.update_layout(title='Extended Envelope Methods',
                  xaxis_title='Sample', yaxis_title='Amplitude')
fig.show()

Trend Analysis — Extended Methods

from vitalDSP.physiological_features.trend_analysis import TrendAnalysis

ta = TrendAnalysis(signal_col)

poly_trend     = ta.compute_polynomial_trend(degree=3)
momentum       = ta.compute_momentum(window_size=16)
trend_strength = ta.compute_trend_strength()
reversals      = ta.detect_trend_reversal(window_size=32)

print(f'Trend Strength  : {trend_strength:.6f}')
print(f'Trend Reversals : {len(reversals)} detected')

fig = go.Figure()
fig.add_trace(go.Scatter(y=signal_col[:512], mode='lines', name='ECG Signal', opacity=0.4))
fig.add_trace(go.Scatter(y=poly_trend[:512], mode='lines', name='Polynomial Trend (deg=3)'))
fig.add_trace(go.Scatter(y=momentum[:512], mode='lines', name='Momentum'))
rev_in = [r for r in reversals if r < 512]
if rev_in:
    fig.add_trace(go.Scatter(x=rev_in, y=[signal_col[r] for r in rev_in],
                             mode='markers', marker=dict(color='red', size=8),
                             name='Trend Reversals'))
fig.update_layout(title='Trend Analysis Extended',
                  xaxis_title='Sample', yaxis_title='Amplitude')
fig.show()

Signal Change Detection — Extended Methods

from vitalDSP.physiological_features.signal_change_detection import SignalChangeDetection

scd = SignalChangeDetection(signal_col)

abs_diff   = scd.absolute_difference()
energy_pts = scd.energy_based_detection(window_size=64)
adapt_pts  = scd.adaptive_threshold_detection(threshold_factor=1.5, window_size=64)

print(f'Absolute difference mean  : {abs_diff.mean():.4f}')
print(f'Energy-based change points: {len(energy_pts)}')
print(f'Adaptive threshold points : {len(adapt_pts)}')

fig = go.Figure()
fig.add_trace(go.Scatter(y=abs_diff[:512], mode='lines', name='Absolute Difference'))
fig.update_layout(title='Absolute Difference Signal',
                  xaxis_title='Sample', yaxis_title='|delta|')
fig.show()

Signal Segmentation — Extended Methods

from vitalDSP.physiological_features.signal_segmentation import SignalSegmentation

seg = SignalSegmentation(signal_col)

thresh_segs = seg.threshold_based_segmentation(threshold=0.5)
peak_segs   = seg.peak_based_segmentation(min_distance=32, height=0.5)

print(f'Threshold-based segments : {len(thresh_segs)}')
print(f'Peak-based segments      : {len(peak_segs)}')

Energy Analysis — Extended Methods

from vitalDSP.physiological_features.energy_analysis import EnergyAnalysis

ea = EnergyAnalysis(signal_col, fs=fs)

seg_energy = ea.compute_segment_energy(start_idx=0, end_idx=256)
qrs_energy = ea.compute_qrs_energy(r_peaks=[100, 220, 340, 460])

print(f'Segment energy (first 256 samples) : {seg_energy:.4f}')
print(f'QRS energy (4 synthetic peaks)     : {qrs_energy:.4f}')

Signal Power Analysis — Extended Methods

from vitalDSP.physiological_features.signal_power_analysis import SignalPowerAnalysis

spa = SignalPowerAnalysis(signal_col)

rmse       = spa.compute_rmse()
peak_power = spa.compute_peak_power()
freqs, psd = spa.compute_psd(fs=fs, nperseg=256)
band_power = spa.compute_band_power(band=(0.5, 40.0), fs=fs, nperseg=256)

print(f'RMSE                  : {rmse:.4f}')
print(f'Peak Power            : {peak_power:.4f}')
print(f'Band Power (0.5-40Hz) : {band_power:.4f}')

fig = go.Figure()
fig.add_trace(go.Scatter(x=freqs[:60], y=psd[:60], mode='lines', name='PSD'))
fig.update_layout(title='Power Spectral Density',
                  xaxis_title='Frequency (Hz)', yaxis_title='Power')
fig.show()

Frequency Domain HRV Features

from vitalDSP.physiological_features.frequency_domain import FrequencyDomainFeatures
from vitalDSP.physiological_features.beat_to_beat import BeatToBeatAnalysis

bb = BeatToBeatAnalysis(signal_col, fs=fs, signal_type='ECG')
rr_intervals = bb.compute_rr_intervals()
print(f'RR intervals: {len(rr_intervals)} beats')

fdf = FrequencyDomainFeatures(rr_intervals)
psd_dict = fdf.compute_psd()

print('Frequency Domain HRV Features:')
print(f'  ULF power   : {fdf.compute_ulf():.4f} ms2')
print(f'  VLF power   : {fdf.compute_vlf():.4f} ms2')
print(f'  LF power    : {fdf.compute_lf():.4f} ms2')
print(f'  HF power    : {fdf.compute_hf():.4f} ms2')
print(f'  LF/HF ratio : {fdf.compute_lf_hf_ratio():.4f}')
print(f'  Total power : {fdf.compute_total_power():.4f} ms2')
print(f'  LFnu        : {fdf.compute_lfnu():.4f}')
print(f'  HFnu        : {fdf.compute_hfnu():.4f}')

fig = go.Figure()
fig.add_trace(go.Scatter(x=psd_dict['frequencies'], y=psd_dict['psd'],
                         mode='lines', name='HRV PSD'))
fig.update_layout(title='HRV Power Spectral Density',
                  xaxis_title='Frequency (Hz)', yaxis_title='Power (ms2)')
fig.show()