pifeat

Zero-Dependency, Sub-Millisecond EEG & Biosignal Feature Extraction Engine in Rust + Python

pip install pifeat

Why pifeat? • Quickstart • 22 Canonical Features • Benchmarks • API Reference

Early development / experimental / no warranties. APIs, behavior, and packaging may change without notice. Research software provided as-is, without warranty of any kind. Not a medical device and not intended for clinical, diagnostic, or safety-critical use. See LICENSE for full terms.

⚡ Why pifeat?

Building real-time Brain-Computer Interfaces (BCI), neurofeedback loops, or edge ML pipelines on devices like a Raspberry Pi or embedded Linux usually forces a painful trade-off: pull in hundreds of megabytes of scipy / mne / numba / openblas dependencies—or write fragile DSP glue code by hand.

pifeat solves this in a single ~1.1 MB native wheel:

• Zero Runtime Dependencies: No numpy, scipy, mne, numba, or C/Fortran BLAS required at runtime. Works out of the box with Python's built-in array.array, memoryview, bytearray, or numpy.ndarray via the CPython Buffer Protocol.
• 100x–173x Faster than SciPy + NumPy: Extracts all 22 time, frequency (Welch PSD), and nonlinear features across 16 channels (500 samples @ 250 Hz) in ~146 µs with the GIL released.
• Zero Steady-State Allocations: Both the Rust core and extractor.extract_into(...) perform 0 heap allocations on the real-time hot path—eliminating GC pauses and latency spikes.
• Drop-In Scientific Equivalence: Matches scipy.signal.welch (window='hann', detrend='constant', scaling='density') to $< 10^{-12}$ relative error, with compensated two-pass variance to survive large uncalibrated ADC DC offsets ($10^5\text{--}10^6\ \mu\text{V}$).
• Hardware-Agnostic: Works with OpenBCI, BrainFlow, LSL (Lab Streaming Layer), MNE, PiEEG, or any raw ADC/SPI stream.

🚀 Developer Quickstart

30-Second Copy-Paste Example (Pure Python — Zero Dependencies)

import array, math
import pifeat as pf

# 1. Generate 2 seconds (500 samples @ 250 Hz) of synthetic 8-channel EEG (10 Hz Alpha + 20 Hz Beta)
n_ch, n_samples, fs = 8, 500, 250.0
raw = array.array("d", [
    15.0 * math.sin(2 * math.pi * 10.0 * (s / fs)) + 5.0 * math.cos(2 * math.pi * 20.0 * (s / fs))
    for _ in range(n_ch) for s in range(n_samples)
])

# 2. Extract all 22 features per channel in ~70 microseconds
extractor = pf.FeatureExtractor(num_channels=n_ch, window_samples=n_samples, fs=fs)
frame = extractor.extract(raw, layout="channels_first")

print("Features per channel:", len(pf.FEATURE_NAMES))        # 22
print("Ch0 Alpha Power (uV^2):", frame.to_dict()["alpha_power"][0])
print("Ch0 Alpha/Beta Ratio:  ", frame.to_dict()["alpha_beta_ratio"][0])

1. Batch Multichannel Extraction (numpy.ndarray, array.array, or memoryview)

Accepts contiguous float64 ('d') or float32 ('f') buffers in either channels_first ([channels, samples]) or samples_first ([samples, channels]) layout:

import array
import pifeat as pf

extractor = pf.FeatureExtractor(
    num_channels=16,
    window_samples=500,     # 2.0 s window at 250 Hz
    fs=250.0,
    nperseg=250,            # 1.0 s Welch segment (1.0 Hz bin spacing)
    noverlap=125,           # 50% overlap (3 Welch segments per window)
    pe_order=3,             # Permutation entropy embedding dimension
    pe_delay=1,             # Permutation entropy sample delay
    integration="trapezoid" # "trapezoid" (uV^2) or "bin_mean" (uV^2/Hz)
)

# Pass a 2D numpy.ndarray (shape [16, 500]), memoryview, or stdlib array.array('d'|'f')
raw_window = array.array("d", [0.0] * (16 * 500))

frame = extractor.extract(raw_window, layout="channels_first")
if frame.window_valid:
    matrix = frame.matrix()          # list[list[float]], shape [16][22]
    ch0 = frame.channel(0)           # list[float], length 22
    by_name = frame.to_dict()        # dict[str, list[float]]

# Zero-allocation hot path: writes directly into a preallocated writable f64 buffer (array or np.ndarray)
out_buf = array.array("d", [0.0] * (16 * pf.NUM_FEATURES_PER_CHANNEL))
valid = extractor.extract_into(raw_window, out_buf, layout="channels_first")

2. Continuous Real-Time Streaming (StreamingFeatureEngine)

StreamingFeatureEngine maintains per-channel circular buffers, tracks monotonic sequence numbers and timestamps, detects dropped upstream packets, and caches overlapping Welch segment periodograms so each hop only computes the newly completed segment FFT:

import pifeat as pf

stream = pf.StreamingFeatureEngine(
    num_channels=16,
    window_samples=500,     # 2.0 s sliding window
    hop_samples=125,        # Emit a FeatureFrame every 125 samples (500 ms at 250 Hz)
    fs=250.0,
    nperseg=250,
    noverlap=125,
    pe_order=3,
    pe_delay=1,
    integration="trapezoid",
)

# Push sample-by-sample from your hardware loop / LSL / BrainFlow / SPI reader
for seq, timestamp, channels in acquisition_stream:
    frame = stream.push_sample(seq, timestamp, channels)
    if frame is not None:
        # Inspect packet-loss & signal-quality flags before feeding your ML model
        if not frame.window_valid:
            print(f"Warning: {frame.missing_samples} dropped samples or invalid channel (flags={frame.channel_flags})")
        features = frame.matrix()

# Or push a contiguous chunk of samples at once ([samples, channels] or [channels, samples])
frames = stream.push_block(block_buf, start_seq=1, start_timestamp=0.0, layout="samples_first")

3. Standalone 1D DSP Functions

Need just one metric on a 1D signal (list[float], array.array, or 1D numpy.ndarray)?

import pifeat as pf

rms = pf.compute_rms(signal_1d)
var = pf.compute_variance(signal_1d)
ll  = pf.compute_line_length(signal_1d)
zcr = pf.compute_zero_crossing_rate(signal_1d)
activity, mobility, complexity = pf.compute_hjorth(signal_1d)
pe  = pf.compute_permutation_entropy(signal_1d, order=3, delay=1, normalize=True)

# Exact scipy.signal.welch equivalent (returns frequencies and one-sided PSD in uV^2/Hz)
freqs, psd = pf.compute_welch_psd(signal_1d, fs=250.0, nperseg=250, noverlap=125)

# Dict of absolute powers, relative powers, ratios, and spectral entropy
bands = pf.compute_band_powers(signal_1d, fs=250.0, nperseg=250, noverlap=125)

4. Using pifeat in Pure Rust (Cargo.toml)

pifeat can also be embedded directly in Rust firmware, Tauri apps, or native daemons without Python:

[dependencies]
pifeat = { version = "0.1", default-features = false }

🧠 The 22 Features Per Channel

Each channel produces NUM_FEATURES_PER_CHANNEL = 22 features in the canonical order given by pf.FEATURE_NAMES:

Numerical and Scientific Guarantees

• Welch PSD Equivalence (scipy.signal.welch): Uses a periodic Hann window $w[n] = 0.5 - 0.5\cos(2\pi n / L)$ (fftbins=True), per-segment constant mean detrending (detrend='constant'), and one-sided power spectral density scaling $U = f_s \sum w[n]^2$ (scaling='density'). Matches scipy.signal.welch to $< 10^{-12}$ relative error.
• Band Integration Modes:
  • "trapezoid" (default): Composite trapezoidal integration with linear interpolation at exact band boundaries ($[f_{\text{lo}}, f_{\text{hi}}]$), preserving Parseval's theorem in $\mu\text{V}^2$.
  • "bin_mean": Arithmetic mean of PSD bins in $[f_{\text{lo}}, f_{\text{hi}})$ in $\mu\text{V}^2/\text{Hz}$ (matches legacy pieeg_server.spectral).
• Compensated Variance and Single-Pass Hjorth: Two-pass compensated summation prevents catastrophic cancellation when large uncalibrated ADC DC offsets ($10^5\text{--}10^6\ \mu\text{V}$) sit atop microvolt EEG signals. Hjorth first- and second-difference means are evaluated in closed form from boundary samples, requiring zero intermediate allocations.
• Band Ratio Noise-Floor Guard: Denominator band power is guarded by max(1e-12, total_power * 1e-9) so band-limited signals do not produce arbitrary $10^{16}\text{--}10^{29}$ ratios from floating-point FFT floor differences between f32 and f64.
• Degenerate Signals and Validity Flags:
  • Flatlines ($\text{var}(x) \le 10^{-24}$) never produce NaN or Inf; Mobility, Complexity, Spectral Entropy, and Permutation Entropy return 0.0 and set FLAG_FLATLINE (1 << 1) in frame.channel_flags.
  • Non-finite inputs (NaN / Inf) set FLAG_NON_FINITE (1 << 0) on the affected channel and frame.window_valid = False.
  • Sequence number gaps in streaming mode set FLAG_SEQ_GAP (1 << 2), record the exact missing count in frame.missing_samples, and set frame.window_valid = False.

📖 API Reference

FeatureFrame

📊 Benchmarks (100x–170x Faster)

Measured at $f_s = 250\text{ Hz}$ extracting all 22 features per channel (time-domain, Welch PSD, band powers, ratios, spectral entropy, and permutation entropy). Steady-state heap allocations in Rust are verified via a custom #[global_allocator] counter in tests/bench_native.rs.

Native Rust Latency and Steady-State Allocations

Python API Overhead vs. Reference SciPy / NumPy (benches/bench_python_vs_native.py)

• State footprint: < 96 KB for a 16-channel, 500-sample StreamingFeatureEngine.
• Binary size: ~1.18 MB stripped native extension with zero external shared library dependencies.

Build and Test from Source

Requires Rust (rustup) and maturin:

pip install maturin
cd pifeat
maturin develop --release

Run the pure-Rust unit tests and zero-allocation benchmark suite:

cargo test --no-default-features --release
cargo test --no-default-features --release --test bench_native -- --nocapture

Run the Python scientific validation suite and benchmark (requires pytest, numpy, scipy for reference comparison only):

pytest tests/test_scientific_validation.py -v
python benches/bench_python_vs_native.py

Release and Wheel Distribution

.github/workflows/release.yml builds precompiled wheels via PyO3/maturin-action@v1 for:

• Linux: x86_64, aarch64 (Raspberry Pi 4/5 64-bit), armv7 (Raspberry Pi OS 32-bit)
• Windows: x64
• Python: CPython 3.10, 3.11, 3.12, 3.13 + source distribution (sdist)