非官方中文翻译 本页译自 MNE-Python 官方文档: 查看英文原文 ↗。 文档版权归 MNE Developers 所有(BSD-3-Clause)。

从 EEG 数据使用 XDAWN 解码#

使用 Xdawn 进行 ERP 解码 [1] [2]。对每种事件类型训练一组空间 Xdawn 滤波器并应用到信号上;随后拼接通道并重新缩放,构造特征向量输入逻辑回归。

# Authors: Alexandre Barachant <alexandre.barachant@gmail.com>
#
# License: BSD-3-Clause
# Copyright the MNE-Python contributors.
import matplotlib.pyplot as plt
import numpy as np
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import classification_report, confusion_matrix
from sklearn.model_selection import StratifiedKFold
from sklearn.multiclass import OneVsRestClassifier
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import MinMaxScaler

from mne import Epochs, io, pick_types, read_events
from mne.datasets import sample
from mne.decoding import Vectorizer, XdawnTransformer, get_spatial_filter_from_estimator
from mne.utils import check_version

print(__doc__)

data_path = sample.data_path()

设置参数并读取数据

meg_path = data_path / "MEG" / "sample"
raw_fname = meg_path / "sample_audvis_filt-0-40_raw.fif"
event_fname = meg_path / "sample_audvis_filt-0-40_raw-eve.fif"
tmin, tmax = -0.1, 0.3
event_id = {
    "Auditory/Left": 1,
    "Auditory/Right": 2,
    "Visual/Left": 3,
    "Visual/Right": 4,
}
n_filter = 3

# Setup for reading the raw data
raw = io.read_raw_fif(raw_fname, preload=True)
raw.filter(1, 20, fir_design="firwin")
events = read_events(event_fname)

picks = pick_types(raw.info, meg=False, eeg=True, stim=False, eog=False, exclude="bads")

epochs = Epochs(
    raw,
    events,
    event_id,
    tmin,
    tmax,
    proj=False,
    picks=picks,
    baseline=None,
    preload=True,
    verbose=False,
)

# Create classification pipeline
kwargs = dict()
if check_version("sklearn", "1.8"):  # TODO VERSION remove on sklearn 1.8+
    kwargs["l1_ratio"] = 1
else:
    kwargs["penalty"] = "l1"
clf = make_pipeline(
    XdawnTransformer(n_components=n_filter),
    Vectorizer(),
    MinMaxScaler(),
    OneVsRestClassifier(
        LogisticRegression(solver="liblinear", random_state=157, **kwargs)
    ),
)

# Get the data and labels
# X is of shape (n_epochs, n_channels, n_times)
X = epochs.get_data(copy=False)
y = epochs.events[:, -1]

# Cross validator
cv = StratifiedKFold(n_splits=10, shuffle=True, random_state=163)

# Do cross-validation
preds = np.empty(len(y))
for train, test in cv.split(epochs, y):
    clf.fit(X[train], y[train])
    preds[test] = clf.predict(X[test])

# Classification report
target_names = ["aud_l", "aud_r", "vis_l", "vis_r"]
report = classification_report(y, preds, target_names=target_names)
print(report)

# Normalized confusion matrix
cm = confusion_matrix(y, preds)
cm_normalized = cm.astype(float) / cm.sum(axis=1)[:, np.newaxis]

# Plot confusion matrix
fig, ax = plt.subplots(1, layout="constrained")
im = ax.imshow(cm_normalized, interpolation="nearest", cmap=plt.cm.Blues)
ax.set(title="Normalized Confusion matrix")
fig.colorbar(im)
tick_marks = np.arange(len(target_names))
plt.xticks(tick_marks, target_names, rotation=45)
plt.yticks(tick_marks, target_names)
ax.set(ylabel="True label", xlabel="Predicted label")
Normalized Confusion matrix
Opening raw data file /home/circleci/mne_data/MNE-sample-data/MEG/sample/sample_audvis_filt-0-40_raw.fif...
    Read a total of 4 projection items:
        PCA-v1 (1 x 102)  idle
        PCA-v2 (1 x 102)  idle
        PCA-v3 (1 x 102)  idle
        Average EEG reference (1 x 60)  idle
    Range : 6450 ... 48149 =     42.956 ...   320.665 secs
Ready.
Reading 0 ... 41699  =      0.000 ...   277.709 secs...
Filtering raw data in 1 contiguous segment
Setting up band-pass filter from 1 - 20 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: 1.00
- Lower transition bandwidth: 1.00 Hz (-6 dB cutoff frequency: 0.50 Hz)
- Upper passband edge: 20.00 Hz
- Upper transition bandwidth: 5.00 Hz (-6 dB cutoff frequency: 22.50 Hz)
- Filter length: 497 samples (3.310 s)

Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Computing rank from data with rank='full'
    data: rank 59 from info
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Computing rank from data with rank='full'
    data: rank 59 from info
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Computing rank from data with rank='full'
    data: rank 59 from info
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Computing rank from data with rank='full'
    data: rank 59 from info
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Computing rank from data with rank='full'
    data: rank 59 from info
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Computing rank from data with rank='full'
    data: rank 59 from info
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Computing rank from data with rank='full'
    data: rank 59 from info
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Computing rank from data with rank='full'
    data: rank 59 from info
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Computing rank from data with rank='full'
    data: rank 59 from info
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Estimating covariance using EMPIRICAL
Done.
Computing rank from data with rank='full'
    data: rank 59 from info
              precision    recall  f1-score   support

       aud_l       0.78      0.68      0.73        72
       aud_r       0.70      0.78      0.74        73
       vis_l       0.99      0.97      0.98        73
       vis_r       0.94      0.97      0.96        70

    accuracy                           0.85       288
   macro avg       0.85      0.85      0.85       288
weighted avg       0.85      0.85      0.85       288

已拟合的 XdawnTransformer 实例(这里取最后一个交叉验证折)的模式可通过 SpatialFilter 容器可视化。

# Instantiate SpatialFilter
spf = get_spatial_filter_from_estimator(
    clf, info=epochs.info, step_name="xdawntransformer"
)

# Let's first examine the scree plot of generalized eigenvalues
# for each class.
spf.plot_scree(title="")

# We can see that for all four classes ~five largest components
# capture most of the variance, let's plot their patterns.
# Each class will now return its own figure
components_to_plot = np.arange(5)
figs = spf.plot_patterns(
    # Indices of patterns to plot,
    # we will plot the first three for each class
    components=components_to_plot,
    show=False,  # to set the titles below
)

# Set the class titles
event_id_reversed = {v: k for k, v in event_id.items()}
for fig, class_idx in zip(figs, clf[0].classes_):
    class_name = event_id_reversed[class_idx]
    fig.suptitle(class_name, fontsize=16)
  • decoding xdawn eeg
  • decoding xdawn eeg
  • decoding xdawn eeg
  • decoding xdawn eeg
  • Auditory/Left, Pattern0, Pattern1, Pattern2, Pattern3, Pattern4, AU
  • Auditory/Right, Pattern0, Pattern1, Pattern2, Pattern3, Pattern4, AU
  • Visual/Left, Pattern0, Pattern1, Pattern2, Pattern3, Pattern4, AU
  • Visual/Right, Pattern0, Pattern1, Pattern2, Pattern3, Pattern4, AU

参考文献#

脚本总运行时间:(0 分 3.969 秒)

图库由 Sphinx-Gallery 生成