Acquisition quality control

Before cleaning a recording, check how well it was made: a flat channel, a loose EOG electrode, a head that moved, or a sidecar that disagrees with the data are better found the day the data were recorded than months later.

meeg_utils.qc answers these questions without changing the data. This tutorial runs it on MNE’s sample recording (Neuromag MEG and EEG).

import mne

import meeg_utils as meu
from meeg_utils import qc

mne.set_log_level("ERROR")

data_path = mne.datasets.sample.data_path()
meg_dir = data_path / "MEG" / "sample"
raw = mne.io.read_raw_fif(meg_dir / "sample_audvis_raw.fif").load_data()
raw.info["line_freq"] = 60.0

Run every check

inspect() picks the checks that apply to the data (MEG, EEG, the acquisition system, BIDS) and lists the others with the reason they were skipped:

report = qc.inspect(raw)
print(report)

<QCReport /home/runner/mne_data/MNE-sample-data/MEG/sample/sample_audvis_raw.fif level='warn'>
  ok   amplitude
  ok   bridging
  warn outlier_channels
         - outlying EEG channels (LOF > 1.5): 1 (> 0); EEG 007
  warn narrowband_noise
         - GRAD channels with unusually strong line noise: 2 (> 0); check their contact: MEG 1033, MEG 1042
         - EEG channels with unusually strong line noise: 1 (> 0); check their contact: EEG 025
  info muscle
  info blinks
  info heart_rate
  warn squid_jumps
         - SQUID jumps: 8 (> 0); in ['MEG 2313']; mark the channels bad or annotate the jumps
  ok   digitization
  ok   events
  skip impedance: found no impedances in the recording (pass impedances=...)
  skip head_movement: found no continuous head localization (cHPI) in the recording
  skip chpi_snr: found no continuous head localization (cHPI) in the recording
  skip empty_room: needs an empty-room recording (the recording is not in a BIDS dataset)
  skip bids_metadata: the recording is not in a BIDS dataset

Each finding has a level:

fail / warn

A problem by definition (a flat channel, a SQUID jump, a bridged electrode, a sidecar that disagrees with the data) or a channel that stands out from the others in this recording.

ok

The criterion was checked and passed.

info

A measurement without a verdict: head movement, blink rate, heart rate, muscle activity, … Their normal range depends on the participants and the setup (children move more than adults, passive electrodes have higher impedances than active ones), so meeg-utils does not pretend to know it.

report.flags lists the problems, worst first:

for finding in report.flags:
    print(f"{finding.level:4}  {finding.check}: {finding.message}")

warn  outlier_channels: outlying EEG channels (LOF > 1.5): 1 (> 0); EEG 007
warn  narrowband_noise: GRAD channels with unusually strong line noise: 2 (> 0); check their contact: MEG 1033, MEG 1042
warn  narrowband_noise: EEG channels with unusually strong line noise: 1 (> 0); check their contact: EEG 025
warn  squid_jumps: SQUID jumps: 8 (> 0); in ['MEG 2313']; mark the channels bad or annotate the jumps

Look at what was measured

Every check draws what it measured. Channels that stand out from the others (Local Outlier Factor, per channel type):

figs = report.checks["outlier_channels"].plot()

Blinks show whether the EOG electrodes work; the average blink should look like one:

figs = report.checks["blinks"].plot()

Without an ECG channel, MNE builds a synthetic ECG from the magnetometers. meeg-utils judges the heartbeat detection by how regular the intervals between beats are, and does not report a rate it cannot trust:

print(report.checks["heart_rate"].metrics_)
figs = report.checks["heart_rate"].plot()

{'heart_rate_bpm': 62.0, 'reliable': True, 'n_beats': 286, 'rr_sd_ms': 104.6, 'irregular_fraction': 0.007, 'channel': 'synthetic'}

The digitized head shape, the sphere fitted to it, and the MEG sensors:

figs = report.checks["digitization"].plot()

The events found in the recording:

figs = report.checks["events"].plot()

The empty room

The empty-room recording of the session shows the noise of the sensors and the room on the day. In a BIDS dataset it is found automatically; here we pass the file:

er = qc.EmptyRoom(empty_room=meg_dir / "ernoise_raw.fif").compute(raw)
print(er.metrics_["mag"], er.metrics_["grad"], sep="\n")
figs = er.plot()

{'noise_floor': 3.989, 'unit': 'fT/√Hz', 'recording_above_db': 8.23, 'outlier_channels': ['MEG 0321', 'MEG 0411', 'MEG 0541', 'MEG 0631', 'MEG 2211', 'MEG 2441']}
{'noise_floor': 2.718, 'unit': 'fT/cm/√Hz', 'recording_above_db': 1.5, 'outlier_channels': ['MEG 0613', 'MEG 1332', 'MEG 2422']}

Your own thresholds

The measurements reported as info become verdicts when you give them a threshold that fits your participants and equipment. Every threshold is a parameter of its check:

report = qc.inspect(
    raw,
    checks=[
        qc.Blinks(min_rate=5),  # blinks per minute
        qc.Muscle(warn_fraction=0.05),  # share of the recording
        qc.Events(expected={"1": 72, "2": 73, "3": 73, "4": 71}),
    ],
)
print(report)

<QCReport /home/runner/mne_data/MNE-sample-data/MEG/sample/sample_audvis_raw.fif level='ok'>
  ok   blinks
  ok   muscle
  ok   events

A whole dataset

Across many recordings the info measurements get a reference: the rest of the dataset. inspect_dataset() runs the checks on every recording of a BIDS dataset and flags the recordings that stand out on any metric (modified z-score above 3.5):

result = meu.qc.inspect_dataset("bids/", n_jobs=8)
result.outliers          # recording, metric, value, z
result.to_csv("qc/summary.csv")
result.plot()["levels"]  # recordings x checks

or from the command line: meu qc bids/ --out qc/, which also writes an HTML report per recording with meeg_utils.report.build():

html = meu.report.build(qc=qc.inspect(raw), title="QC of sample_audvis_raw")
html

<Report | 1 title | 13 items | QC of sample_audvis_raw  | 0.4 MB
 Quality checks                                         | 0.4 MB
>