非官方中文翻译
本页译自 MNE-Python 官方文档：
查看英文原文 ↗。
文档版权归 MNE Developers 所有（BSD-3-Clause）。
说明
转到末尾下载完整示例代码。
标注肌肉伪迹#
肌肉收缩产生的高频活动会掩盖感兴趣的脑信号。咬紧牙关、吞咽或抽动颅部肌肉都可能产生肌肉伪迹，这类伪迹在 110–140 Hz 范围最明显。
本示例使用 annotate_muscle_zscore() 标注可能存在肌肉活动的时段：先在 110–140 Hz 范围对数据带通滤波，再用 Hilbert 解析信号取包络，只考虑高频信号的绝对幅度而非相位。包络做 z 分数化后跨通道求和，并除以通道数的平方根。由于肌肉伪迹持续数百毫秒，再对平均 z 分数施加 4 Hz 低通滤波以去除瞬态峰。超过设定阈值的时段标注为 BAD_muscle。此外，min_length_good 参数决定肌肉伪迹之间短暂的“良好数据”是否并入周围的“BAD”注释。
# Authors: Adonay Nunes <adonay.s.nunes@gmail.com>
# Luke Bloy <luke.bloy@gmail.com>
# License: BSD-3-Clause
# Copyright the MNE-Python contributors.
import matplotlib.pyplot as plt
import numpy as np
from mne.datasets.brainstorm import bst_auditory
from mne.io import read_raw_ctf
from mne.preprocessing import annotate_muscle_zscore
# Load data
data_path = bst_auditory.data_path()
raw_fname = data_path / "MEG" / "bst_auditory" / "S01_AEF_20131218_01.ds"
raw = read_raw_ctf(raw_fname, preload=False)
raw.crop(130, 160).load_data() # just use a fraction of data for speed here
raw.resample(300, npad="auto")
ds directory : /home/circleci/mne_data/MNE-brainstorm-data/bst_auditory/MEG/bst_auditory/S01_AEF_20131218_01.ds
res4 data read.
hc data read.
Separate EEG position data file read.
Quaternion matching (desired vs. transformed):
2.51 74.26 0.00 mm <-> 2.51 74.26 -0.00 mm (orig : -56.69 50.20 -264.38 mm) diff = 0.000 mm
-2.51 -74.26 0.00 mm <-> -2.51 -74.26 -0.00 mm (orig : 50.89 -52.31 -265.88 mm) diff = 0.000 mm
108.63 0.00 0.00 mm <-> 108.63 0.00 0.00 mm (orig : 67.41 77.68 -239.53 mm) diff = 0.000 mm
Coordinate transformations established.
Reading digitizer points from ['/home/circleci/mne_data/MNE-brainstorm-data/bst_auditory/MEG/bst_auditory/S01_AEF_20131218_01.ds/S01_20131218_01.pos']...
Polhemus data for 3 HPI coils added
Device coordinate locations for 3 HPI coils added
5 extra points added to Polhemus data.
Measurement info composed.
Finding samples for /home/circleci/mne_data/MNE-brainstorm-data/bst_auditory/MEG/bst_auditory/S01_AEF_20131218_01.ds/S01_AEF_20131218_01.meg4:
System clock channel is available, checking which samples are valid.
360 x 2400 = 864000 samples from 340 chs
Current compensation grade : 3
Reading 0 ... 72000 = 0.000 ... 30.000 secs...
Finding events on: UDIO001, UPPT001, UTRG001
4 events found on stim channel UDIO001
Event IDs: [64]
23 events found on stim channel UPPT001
Event IDs: [1 2]
Finding events on: UDIO001, UPPT001, UTRG001
4 events found on stim channel UDIO001
Event IDs: [64]
23 events found on stim channel UPPT001
Event IDs: [1 2]
对数据做陷波滤波：
说明
若存在工频噪声，应在检测肌肉伪迹前先做陷波滤波。示例参见 Power line noise。
raw.notch_filter([60, 120])
Filtering raw data in 1 contiguous segment
Setting up band-stop filter
FIR filter parameters
---------------------
Designing a one-pass, zero-phase, non-causal bandstop filter:
- Windowed time-domain design (firwin) method
- Hamming window with 0.0194 passband ripple and 53 dB stopband attenuation
- Lower transition bandwidth: 0.50 Hz
- Upper transition bandwidth: 0.50 Hz
- Filter length: 1981 samples (6.603 s)
# The threshold is data dependent, check the optimal threshold by plotting
# ``scores_muscle``.
threshold_muscle = 5 # z-score
# Choose one channel type, if there are axial gradiometers and magnetometers,
# select magnetometers as they are more sensitive to muscle activity.
annot_muscle, scores_muscle = annotate_muscle_zscore(
raw,
ch_type="mag",
threshold=threshold_muscle,
min_length_good=0.2,
filter_freq=[110, 140],
)
Removing 5 compensators from info because not all compensation channels were picked.
Filtering raw data in 1 contiguous segment
Setting up band-pass filter from 1.1e+02 - 1.4e+02 Hz
FIR filter parameters
---------------------
Designing a one-pass, zero-phase, non-causal bandpass filter:
- Windowed time-domain design (firwin) method
- Hamming window with 0.0194 passband ripple and 53 dB stopband attenuation
- Lower passband edge: 110.00
- Lower transition bandwidth: 27.50 Hz (-6 dB cutoff frequency: 96.25 Hz)
- Upper passband edge: 140.00 Hz
- Upper transition bandwidth: 10.00 Hz (-6 dB cutoff frequency: 145.00 Hz)
- Filter length: 99 samples (0.330 s)
Setting up low-pass filter at 4 Hz
FIR filter parameters
---------------------
Designing a one-pass, zero-phase, non-causal lowpass filter:
- Windowed time-domain design (firwin) method
- Hamming window with 0.0194 passband ripple and 53 dB stopband attenuation
- Upper passband edge: 4.00 Hz
- Upper transition bandwidth: 2.00 Hz (-6 dB cutoff frequency: 5.00 Hz)
- Filter length: 495 samples (1.650 s)
绘制整个记录期间的肌肉 z 分数#
fig, ax = plt.subplots()
ax.plot(raw.times, scores_muscle)
ax.axhline(y=threshold_muscle, color="r")
ax.set(xlabel="time, (s)", ylabel="zscore", title="Muscle activity")
查看注释#
order = np.arange(144, 164)
raw.set_annotations(annot_muscle)
raw.plot(start=5, duration=20, order=order)
脚本总运行时间：（0 分 1.290 秒）