Skip to content

INTERICTAL_PIPELINE_STEPS.md

Latest commit

 

History

History
282 lines (212 loc) · 12.9 KB

INTERICTAL_PIPELINE_STEPS.md

File metadata and controls

282 lines (212 loc) · 12.9 KB
Raw

Interictal (IED) pipeline – user guide

This document describes the interictal epileptiform discharge (IED) detection pipeline: how to run it, what it does, and what it outputs. The method is based on Janca et al. (2015).

Reference: Janca, R., et al. (2015). Detection of interictal epileptiform discharges using signal envelope distribution modelling. Brain Topography, 28(6), 770–778. https://link.springer.com/article/10.1007/s10548-014-0379-1

Quick start

Requirements: a bidsified session — an _ieeg.edf with a _channels.tsv next to it, under <bids_data_dir>/sub-<sub_id>/ses-<ses_id>/. The _channels.tsv decides which channels are sEEG; channels typed anything else (EKG, TRIG, DC…) are read but not analysed, and are not carried into the output. Run the bidsify pipeline first if a session does not have one.

Option 1: Config file

  1. Copy the example config and edit the session entities:

    cp configs/example_configs/example_ied_config.json my_config.json
# Edit my_config.json: set bids_data_dir, sub_id, ses_id

  2. Run the pipeline:

    python -m ieeg_prep.ied.analysis --config my_config.json

Option 2: Command-line only

python -m ieeg_prep.ied.analysis \
  --bids-data-dir /path/to/bids \
  --sub-id EMOP0649 --ses-id ieeg20260310 \
  --dataset-description configs/ied_dataset_description.json

Option 3: HPC (SLURM)

sbatch batch_script/run_ied_analysis.batch

Config file format

config.json can contain:

Key Required Description
bids_data_dir Yes Root of the BIDS dataset. Output goes to its derivatives/ directory
sub_id Yes Subject label, with or without the sub- prefix
ses_id Yes Session label, with or without the ses- prefix
task_name Only if the session holds more than one task task- label
run_id Only if the session holds more than one run run- index, zero-padded or not
dataset_description Only when the derivatives directory does not exist yet Path to the dataset_description.json to copy in
derivatives_name No (default: ied-detection) Directory under derivatives/ to write to
channels_above_rate_thresh No (default: 6.5) Spike rate threshold (spikes/min) for the channels_above_rate boolean
interictal_config No Object with pipeline parameters (see below)
metadata No Extra keys merged into _metadata.json, for anything the pipeline cannot derive

notch_freq defaults to the recording's own PowerLineFrequency from its _ieeg.json, not to 60 Hz; setting it in the config or on the command line still wins.

Any option can be overridden from the command line. Example config:

{
  "bids_data_dir": "/path/to/bids",
  "sub_id": "EMOP0649",
  "ses_id": "ieeg20260310",
  "task_name": "pilotbattery",
  "run_id": "01",
  "derivatives_name": "ied-detection",
  "dataset_description": "configs/ied_dataset_description.json",
  "channels_above_rate_thresh": 6.5,
  "interictal_config": {
    "sfreq_target": 200,
    "bandpass_low": 10,
    "bandpass_high": 60,
    "notch_method": "custom_biquad",
    "notch_freq": 60,
    "win_coeff": 5,
    "overlap_coeff": 4,
    "threshold_k": 3.65,
    "polyspike_union_sec": 0.12,
    "trim_edges_sec": 1
  }
}

Pipeline parameters (interictal_config)

Parameter Default Description
sfreq_target 200 Resampling target (Hz)
bandpass_low 10 Bandpass low cutoff (Hz)
bandpass_high 60 Bandpass high cutoff (Hz)
bandpass_order 8 IIR filter order
bandpass_ftype "cheby2" IIR type ("cheby2", etc.)
bandpass_rs 60 Cheby2 stopband ripple (dB)
notch_freq 60 Line-noise notch frequency (Hz)
notch_radius 0.985 Pole radius for custom biquad notch
notch_method "custom_biquad" "custom_biquad" or "mne"
win_coeff 5 Envelope segment window (seconds)
overlap_coeff 4 Envelope segment overlap (seconds)
smooth_window_len 5 Moving-average length for μ/σ smoothing
threshold_k 3.65 Threshold multiplier: k × (mode + median)
polyspike_union_sec 0.12 Merge spikes within this window (seconds)
trim_edges_sec 1 Drop first/last N seconds

Outputs

Results are written as a BIDS derivatives dataset:

<bids_data_dir>/derivatives/ied-detection/
    dataset_description.json
    sub-EMOP0004/ses-litieeg20260213a/ieeg/
        sub-EMOP0004_ses-litieeg20260213a_task-pilotbattery_run-01_desc-ied_events.tsv
                                                                 ..._desc-ied_events.json
                                                                 ..._desc-ied_channels.tsv
                                                                 ..._desc-ied_ieeg.fif
                                                                 ..._desc-ied_ieeg.json
                                                                 ..._desc-ied_metadata.json
                                                                 ..._desc-ied_results.npz

Every file carries the source recording's entities plus a desc-ied label, so a derivative can always be traced back to what produced it. Re-running overwrites in place.

The derivatives/ied-detection/ directory is created on first use, but only if a dataset_description.json is supplied to copy into it — a derivatives directory without one is not a BIDS dataset, and its contents (dataset name, generating software) cannot be worked out from the data. The run stops with an error if one is neither already present nor provided. See configs/ied_dataset_description.json.

_desc-ied_events.tsv

One row per detection, ordered by onset and then by channel:

Column Description
onset Seconds from the start of the source recording
duration Always 0 — detections are point events
trial_type Always ied
channel sEEG channel the detection is on

_desc-ied_events.json is the accompanying column dictionary.

_desc-ied_channels.tsv

The source recording's _channels.tsv for the sEEG channels only, in pipeline order, with two changes:

  • sampling_frequency, low_cutoff, high_cutoff and notch are rewritten to describe the .fif next to it, not the source EDF.
  • Three columns are appended: ied_count, ied_rate (spikes/min), and above_ied_rate_threshold (true/false against channels_above_rate_thresh).

_desc-ied_ieeg.fif

The resampled, bandpassed and notched signal the detector actually ran on, sEEG channels only. In pipeline terms:

Raw EDF → resample / bandpass / notch → THIS FILE → envelope → lognormal μ/σ
        → threshold → candidates → polyspike merge → detections (_events.tsv)

So it is the detector's input, not its output — the last point at which the data is still a continuous voltage trace. Load it with mne.io.read_raw_fif to eyeball a detection against the signal without recomputing the filtering.

Two caveats. It is bandpassed to 10–60 Hz, so it is unusable for anything outside IED work. And onsets in _events.tsv sit exactly on this file's 200 Hz sample grid, so raw.time_as_index(onset) is exact here, but only ±2.5 ms against the original 2048 Hz EDF.

_desc-ied_ieeg.json

The BIDS sidecar for that .fif — required, since BIDS expects an _ieeg.json beside every iEEG data file. Descriptive fields (TaskName, PowerLineFrequency, iEEGReference, Manufacturer, …) are carried over from the source recording's own sidecar; everything the filtering changed is recomputed:

  • SamplingFrequency and RecordingDuration describe the .fif, not the source.
  • SEEGChannelCount is the channels actually kept; the other count fields are explicitly 0, since dropped is different from unknown.
  • SoftwareFilters records the bandpass and notch as a structured object. The source usually says "n/a"; this is the only machine-readable record of how the .fif was produced.

_desc-ied_results.npz

NumPy archive (load with np.load(...)):

Key Shape / type Description
spike_markers (n_channels, n_samples) bool 1 where a spike was detected
channel_spike_counts (n_channels,) int Total spikes per channel
channel_spike_rates (n_channels,) float Spikes per minute per channel
channels_above_rate (n_channels,) bool Channels with rate > channels_above_rate_thresh
ch_names 1D array of str sEEG channel names (same order as above)
sfreq scalar Sampling frequency (Hz) after resampling

_desc-ied_metadata.json

What the detection run found, and nothing else:

Key Description
Description Fixed text naming what this dataset is
Sources bids:: URIs of the EDF and sidecars that were read
RecordingDurationMinutes Duration the rates are computed over
SEEGChannelCount Channels analysed
IEDCount Total detections across all channels
IEDRateThreshold Rate in spikes/min a channel is flagged above
ChannelsAboveIEDRateThreshold How many channels exceeded it
Parameters Full interictal_config used for the run

Deliberately small. Anything another file in the run already holds is not repeated: per-channel counts and rates are columns of _channels.tsv, channel names are its name column, the arrays are in the .npz, and the recording's own properties are in _ieeg.json. The two JSON files divide by subject — _ieeg.json describes the recording, _metadata.json describes the detection — and state nothing in common.

Anything the pipeline cannot derive goes in the config's metadata object and is merged in last, so it can also override any key above.

Pipeline steps (what the code does)

  1. Locate the session from its BIDS entities, load the EDF and set channel types from its _channels.tsv (bidsify.query.apply_bids_channel_types).
  2. Resample to sfreq_target (default 200 Hz).
  3. Bandpass 10–60 Hz (Cheby2 IIR, order 8) on sEEG only.
  4. Notch 60 Hz line noise (custom biquad or MNE IIR).
  5. Hilbert envelope of the filtered sEEG.
  6. Lognormal fit per sliding segment (window win_coeff s, overlap overlap_coeff s) → μ, σ per segment.
  7. Smooth μ and σ along segments (moving average).
  8. Interpolate μ and σ to every sample (cubic spline).
  9. Threshold = k × (lognormal mode + lognormal median).
  10. Detection: envelope > threshold → crossings → one peak per crossing → merge nearby spikes (polyspike union) → one local max per merged event → trim first/last N seconds.
  11. Rates: spike count and rate per channel; boolean channels_above_rate using the chosen rate threshold.

Using results in Python

import json
from pathlib import Path

import mne
import numpy as np
import pandas as pd

bids = Path("/path/to/bids")
sub, ses = "sub-EMOP0004", "ses-litieeg20260213a"
out = bids / "derivatives/ied-detection" / sub / ses / "ieeg"
stem = f"{sub}_{ses}_task-pilotbattery_run-01_desc-ied"

# Per-channel rates, and which channels are above threshold. pandas may parse
# the true/false column as bool or leave it as text, so normalise before testing.
channels = pd.read_csv(out / f"{stem}_channels.tsv", sep="\t")
above = channels["above_ied_rate_threshold"].astype(str).str.lower() == "true"
high_rate = channels.loc[above, "name"]

# Individual detections
events = pd.read_csv(out / f"{stem}_events.tsv", sep="\t")
print(events.groupby("channel").size())

# The filtered signal the detector ran on
raw = mne.io.read_raw_fif(out / f"{stem}_ieeg.fif")

# Or the arrays, in the layout earlier versions wrote
data = np.load(out / f"{stem}_results.npz", allow_pickle=True)
spike_markers = data["spike_markers"]           # (n_channels, n_samples) bool

with open(out / f"{stem}_metadata.json") as f:
    meta = json.load(f)
print(meta["IEDCount"], meta["RecordingDurationMinutes"], meta["IEDRateThreshold"])

Command-line reference

python -m ieeg_prep.ied.analysis [OPTIONS]

Required (if not in config): --bids-data-dir, --sub-id, --ses-id.

Common options:

  • --config PATH – load options from JSON
  • --bids-data-dir PATH – root of the BIDS dataset
  • --sub-id LABEL, --ses-id LABEL – which session to run on
  • --task-name LABEL, --run-id INDEX – only needed when the session holds more than one recording
  • --dataset-description PATH – dataset_description.json to install if the derivatives directory is new
  • --derivatives-name NAME – directory under derivatives/ (default: ied-detection)
  • --channels-above-rate-thresh FLOAT – rate threshold for boolean mask (default: 6.5)
  • --notch-method {custom_biquad,mne} – notch filter method
  • --sfreq-target, --bandpass-low, --bandpass-high, etc. – override config

Run with --help for the full list.