Устранить моргание от сигнала ЭЭГ с помощью ICA

Question

Я новичок в scikit-learn, но я пытаюсь убрать моргание глаз (шумовые пики) внутри одного канала ЭЭГ. Я искал по Интернету, но вижу только более сложные чтения с MNE, PyEEG или другими модулями Python. Я просто хочу что-то простое и зависимое только от склеарн. Вот мой код:

#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()

То, что я ожидал, было разделением на два компонента, более чистый сигнал ЭЭГ и моргание глаз. Я не могу понять, в чем проблема. Может кто-нибудь протянуть руку помощи?

Answer 1

В случае, если вы работаете с одноканальной ЭЭГ, простой способ реализации может быть следующим:

1) Обнаружение миганий в сигнале x с использованием простого определения пиков на основе порогов (вам может потребоваться выяснить это, посмотрев на сигналы в течение нескольких миганий). Устройства Neurosky, Muse и т. Д. Поставляются с API для обнаружения морганий. Вы можете использовать это при необходимости.

2) Возьмите кадр, соответствующий миганию ( xb ). Установите на нем плавную линию ( xbs ).

3) Вычтите плавную линию ( xbs ) из мигания ( xb ) и добавьте к этому среднее значение сигнала. Давайте назовем это как xbc .

4) Заменить xbc вместо xb в исходном сигнале x .

Answer 2

Вы заметили, что ваши "компоненты" в точности соответствуют исходному сигналу и перевернуты? Это потому, что вы не можете получить больше компонентов, чем сигналов.

Вам необходимо выполнить следующие шаги:

1. подача всех каналов ЭЭГ в 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.

---

Если у вас есть записи EOG, есть надежда! Существует Полностью автоматизированный метод коррекции артефактов EOG в записях ЭЭГ . Этот подход основан на канонической корреляции или регрессии (я не помню деталей), но вам необходимо записать сигналы EOG вместе с EEG.

---

Я создал рабочий пример с имитацией данных "ЭЭГ". Данные состоят из трех каналов: фронтального, центрального и теменного. Альфа-активность 10 Гц является наиболее сильной в задней части, а несколько похожих на мигание шипов являются наиболее сильными в передней части.

Надеюсь, этот пример лучше иллюстрирует, как удалять компоненты из многоканальных данных.

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')