使用 ICA 消除 EEG 信号中的眨眼现象

我是 scikit-learn 的新手，但我正在尝试消除单个脑电图通道内的眨眼（噪声峰值）。 我在互联网上搜索过，但只看到了 MNE、PyEEG 或其他 Python 模块的更复杂的读物。我只想要一些简单的东西并且只依赖于sklearn。这是我的代码：

#The channel containing some eye-blinks
X = f1ep1_data[:,[4]]

#run ICA on signal
ica = FastICA(n_components=2)
ica.fit(X)

#reconstruct signal with independent components
components = ica.fit_transform(X)
X_restored = ica.inverse_transform(components)

fig1 = plt.figure()
plt.subplot(3,1,1)
plt.title("Original signal")
plt.plot(f1ep1_timescale, X)

plt.subplot(3,1,2)
plt.title("Components")
plt.plot(f1ep1_timescale, components)

plt.subplot(3,1,3)
plt.title("Signal Reconstructed")
plt.plot(f1ep1_timescale, X_restored)
plt.draw()

我期待的是分成两个部分，更清晰的脑电图信号和眨眼。我不知道出了什么问题。有人可以伸出援助之手吗？

您是否注意到您的“组件”正是按比例缩放并颠倒的原始信号？这是因为你无法获得比信号更多的组件。

您需要执行以下步骤：

1. feed allEEG 通道进入 ICA
2. 手动删除包含眨眼或其他伪影的组件
3. 使用逆变换重建

我们来详细看看第2步：为什么要手动删除组件？ICA 对眨眼一无所知。它根据统计测量将信号分离成多个分量。如果幸运的话，其中一些组件看起来就像眨眼一样。

到目前为止还好，但真正的问题是组件的顺序没有定义。运行 ICA，您可能会发现组件 1 包含眨眼。再次运行，它们位于组件 3 中。再次运行，它们位于组件 2 和 5 中...

无法提前知道要删除哪些组件以及要删除多少组件。这就是为什么您需要在每次运行时手动将其告知算法。

代码看起来像这样：

# Use all channels - they will contain eye blinks to varying degrees
X = f1ep1_data[:, :]

# run ICA on signal
ica = FastICA(n_components=x.shape[1])  # we want *all* the components
ica.fit(X)

# decompose signal into components
components = ica.fit_transform(X)

# plot components and ask user which components to remove
# ...
remove_indices = [0, 1, 3]  # pretend the user selected components 0, 1, and 3

# "remove" unwanted components by setting them to 0 - simplistic but gets the job done
components[:, remove_indices] = 0

#reconstruct signal
X_restored = ica.inverse_transform(components)

您很可能对手动步骤不满意。运气不好:) 2013 年，还没有强大的自动算法来标记眨眼组件。我不认为这已经改变，但如果有什么东西你会发现它是领域特定的软件包之一，如 MNE 或 PyEEG。

不过，如果您碰巧有眼电图记录，那就还有希望！那里存在脑电图记录中脑电图伪影的全自动校正方法 https://www.ncbi.nlm.nih.gov/pubmed/17088100。该方法基于典型相关或回归（我不记得细节），但您需要将 EOG 信号与 EEG 一起记录。

我使用模拟“EEG”数据创建了一个工作示例。 数据由三个通道组成：额叶通道、中央通道和顶叶通道。 10 Hz 阿尔法活动在后面最强，一些类似眨眼的尖峰在前面最强。

希望这个示例能够更好地说明如何从多通道数据中删除组件。

import numpy as np
import scipy.signal as sps
from sklearn.decomposition import FastICA
import matplotlib.pyplot as plt

np.random.seed(42)

n = 1000
fs = 100
noise = 3

# simulate EEG data with eye blinks
t = np.arange(n)
alpha = np.abs(np.sin(10 * t / fs)) - 0.5
alpha[n//2:] = 0
blink = np.zeros(n)
blink[n//2::200] += -1
blink = sps.lfilter(*sps.butter(2, [1*2/fs, 10*2/fs], 'bandpass'), blink)

frontal = blink * 200 + alpha * 10 + np.random.randn(n) * noise
central = blink * 100 + alpha * 15 + np.random.randn(n) * noise
parietal = blink * 10 + alpha * 25 + np.random.randn(n) * noise

eeg = np.stack([frontal, central, parietal]).T  # shape = (100, 3)

# plot original data
plt.subplot(3, 1, 1)
plt.plot(frontal + 50)
plt.plot(central + 100)
plt.plot(parietal + 150)
plt.yticks([50, 100, 150], ['Fz', 'Cz', 'Pz'])
plt.ylabel('original data')

# decompose EEG and plot components
ica = FastICA(n_components=3)
ica.fit(eeg)
components = ica.transform(eeg)

plt.subplot(3, 1, 2)
plt.plot([[np.nan, np.nan, np.nan]])  # advance the color cycler to give the components a different color :)
plt.plot(components + [0.5, 1.0, 1.5])
plt.yticks([0.5, 1.0, 1.5], ['0', '1', '2'])
plt.ylabel('components')

# looks like component #2 (brown) contains the eye blinks
# let's remove them (hard coded)!
components[:, 2] = 0

# reconstruct EEG without blinks
restored = ica.inverse_transform(components)

plt.subplot(3, 1, 3)
plt.plot(restored + [50, 100, 150])
plt.yticks([50, 100, 150], ['Fz', 'Cz', 'Pz'])
plt.ylabel('cleaned data')