使用Python-EEG工具库MNE-机器学习算法随机森林、SVM、逻辑回归、KNN、朴素贝叶斯、多层感知机、LSTM等方法判断睡眠类型
import numpy as np
import matplotlib.pyplot as plt
import mne
from sklearn.ensemble import RandomForestClassifier
from sklearn.metrics import accuracy_score
from sklearn.metrics import confusion_matrix
from sklearn.metrics import classification_report
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import FunctionTransformer
"""
可以通过
mne.datasets.sleep_physionet.age.fetch_data(subjects,recording,path)
来获取PhysioNet多导睡眠图数据集文件。
subjects:表示想要使用哪些受试者对象,可供选择的受试者对象范围为0-19。
recording:表示夜间记录的编号(索引),有效值为:[1]、[2]或[1、2]。
path:PhysioNet数据的存放地址,如果没有给定,则加载默认存放数据的地址;
如果默认存放数据集的地址不存在数据,则从网络中下载相关数据。
"""
ALICE, BOB = 0, 1
[alice_files, bob_files] = mne.datasets.sleep_physionet.age.fetch_data(subjects=[ALICE, BOB], recording=[1])
mapping = {'EOG horizontal': 'eog',
'Resp oro-nasal': 'misc',
'EMG submental': 'misc',
'Temp rectal': 'misc',
'Event marker': 'misc'}
raw_train = mne.io.read_raw_edf(alice_files[0])
annot_train = mne.read_annotations(alice_files[1])
raw_train.set_annotations(annot_train, emit_warning=False)
raw_train.set_channel_types(mapping)
raw_train.plot(duration=40, scalings='auto')
plt.show()
"""
睡眠表示与事件映射
"""
annotation_desc_2_event_id = {'Sleep stage W': 1,
'Sleep stage 1': 2,
'Sleep stage 2': 3,
'Sleep stage 3': 4,
'Sleep stage 4': 4,
'Sleep stage R': 5}
events_train, _ = mne.events_from_annotations(
raw_train, event_id=annotation_desc_2_event_id, chunk_duration=30.)
event_id = {'Sleep stage W': 1,
'Sleep stage 1': 2,
'Sleep stage 2': 3,
'Sleep stage 3/4': 4,
'Sleep stage R': 5}
mne.viz.plot_events(events_train, event_id=event_id,
sfreq=raw_train.info['sfreq'])
stage_colors = plt.rcParams['axes.prop_cycle'].by_key()['color']
Used Annotations descriptions: ['Sleep stage 1', 'Sleep stage 2', 'Sleep stage 3', 'Sleep stage 4', 'Sleep stage R', 'Sleep stage W']
tmax = 30. - 1. / raw_train.info['sfreq']
"""
所创建的是时间从tmin=0开始,到tmax为止的epochs
"""
epochs_train = mne.Epochs(raw=raw_train, events=events_train,
event_id=event_id, tmin=0., tmax=tmax, baseline=None)
print(epochs_train)
epochs_test = mne.Epochs(raw=raw_test, events=events_test, event_id=event_id,
tmin=0., tmax=tmax, baseline=None)
print(epochs_test)
fig, (ax1, ax2) = plt.subplots(ncols=2)
stages = sorted(event_id.keys())
for ax, title, epochs in zip([ax1, ax2],
['Alice', 'Bob'],
[epochs_train, epochs_test]):
for stage, color in zip(stages, stage_colors):
epochs[stage].plot_psd(area_mode=None, color=color, ax=ax,
fmin=0.1, fmax=20., show=False,
average=True, spatial_colors=False)
ax.set(title=title, xlabel='Frequency (Hz)')
ax2.set(ylabel='uV^2/hz (dB)')
ax2.legend(ax2.lines[2::3], stages)
plt.tight_layout()
plt.show()
from mne.time_frequency import psd_array_welch
def eeg_power_band(epochs):
"""脑电相对功率带特征提取
该函数接受一个 'mne.Epochs' 对象,
并基于与 scikit-learn 兼容的特定频带中的相对功率创建 EEG 特征。
Parameters
----------
epochs : mne.Epochs
The data.
Returns
-------
X : numpy array of shape [n_samples, 5]
Transformed data.
"""
FREQ_BANDS = {"delta": [0.5, 4.5],
"theta": [4.5, 8.5],
"alpha": [8.5, 11.5],
"sigma": [11.5, 15.5],
"beta": [15.5, 30]}
epochs_data = epochs.get_data()
psds, freqs = psd_array_welch(epochs_data, sfreq=epochs.info['sfreq'], fmin=0.5, fmax=30.)
psds /= np.sum(psds, axis=-1, keepdims=True)
X = []
for fmin, fmax in FREQ_BANDS.values():
psds_band = psds[:, :, (freqs >= fmin) & (freqs < fmax)].mean(axis=-1)
X.append(psds_band.reshape(len(psds), -1))
return np.concatenate(X, axis=1)
pipe = make_pipeline(FunctionTransformer(eeg_power_band, validate=False),
RandomForestClassifier(n_estimators=100, random_state=100))
epochs_train.load_data()
y_train = epochs_train.events[:, 2]
pipe.fit(epochs_train, y_train)
y_pred = pipe.predict(epochs_test)
y_test = epochs_test.events[:, 2]
acc = accuracy_score(y_test, y_pred)
print("Accuracy score: {}".format(acc))
print(classification_report(y_test, y_pred, target_names=event_id.keys()))
from sklearn.svm import SVC
pipe = make_pipeline(FunctionTransformer(eeg_power_band, validate=False),
SVC(kernel='linear'))
epochs_train.load_data()
y_train = epochs_train.events[:, 2]
pipe.fit(epochs_train, y_train)
y_pred = pipe.predict(epochs_test)
y_test = epochs_test.events[:, 2]
acc = accuracy_score(y_test, y_pred)
print("Accuracy score: {}".format(acc))
from sklearn.svm import SVC
pipe = make_pipeline(FunctionTransformer(eeg_power_band, validate=False),
SVC(kernel='poly'))
epochs_train.load_data()
y_train = epochs_train.events[:, 2]
pipe.fit(epochs_train, y_train)
y_pred = pipe.predict(epochs_test)
y_test = epochs_test.events[:, 2]
acc = accuracy_score(y_test, y_pred)
print("Accuracy score: {}".format(acc))
from sklearn.svm import SVC
pipe = make_pipeline(FunctionTransformer(eeg_power_band, validate=False),
SVC(kernel='rbf'))
epochs_train.load_data()
y_train = epochs_train.events[:, 2]
pipe.fit(epochs_train, y_train)
y_pred = pipe.predict(epochs_test)
y_test = epochs_test.events[:, 2]
acc = accuracy_score(y_test, y_pred)
print("Accuracy score: {}".format(acc))
from sklearn.svm import SVC
pipe = make_pipeline(FunctionTransformer(eeg_power_band, validate=False),
SVC(kernel='sigmoid'))
epochs_train.load_data()
y_train = epochs_train.events[:, 2]
pipe.fit(epochs_train, y_train)
y_pred = pipe.predict(epochs_test)
y_test = epochs_test.events[:, 2]
acc = accuracy_score(y_test, y_pred)
print("Accuracy score: {}".format(acc))
from sklearn.neighbors import KNeighborsClassifier
pipe = make_pipeline(FunctionTransformer(eeg_power_band, validate=False),
KNeighborsClassifier(n_neighbors=5))
epochs_train.load_data()
y_train = epochs_train.events[:, 2]
pipe.fit(epochs_train, y_train)
y_pred = pipe.predict(epochs_test)
y_test = epochs_test.events[:, 2]
acc = accuracy_score(y_test, y_pred)
print("Accuracy score: {}".format(acc))
from sklearn.linear_model import LogisticRegression
pipe = make_pipeline(FunctionTransformer(eeg_power_band, validate=False),
LogisticRegression(max_iter=1000))
epochs_train.load_data()
y_train = epochs_train.events[:, 2]
pipe.fit(epochs_train, y_train)
y_pred = pipe.predict(epochs_test)
y_test = epochs_test.events[:, 2]
acc = accuracy_score(y_test, y_pred)
print("Accuracy score: {}".format(acc))
from sklearn.neural_network import MLPClassifier
pipe = make_pipeline(FunctionTransformer(eeg_power_band, validate=False),
MLPClassifier(hidden_layer_sizes=(100, ), max_iter=1000))
epochs_train.load_data()
y_train = epochs_train.events[:, 2]
pipe.fit(epochs_train, y_train)
y_pred = pipe.predict(epochs_test)
y_test = epochs_test.events[:, 2]
acc = accuracy_score(y_test, y_pred)
print("Accuracy score: {}".format(acc))
from sklearn.naive_bayes import GaussianNB
pipe = make_pipeline(FunctionTransformer(eeg_power_band, validate=False),
GaussianNB())
epochs_train.load_data()
y_train = epochs_train.events[:, 2]
pipe.fit(epochs_train, y_train)
y_pred = pipe.predict(epochs_test)
y_test = epochs_test.events[:, 2]
acc = accuracy_score(y_test, y_pred)
print("Accuracy score: {}".format(acc))
import torch
import torch.nn as nn
import torch.optim as optim
from torch.utils.data import DataLoader, TensorDataset
from sklearn.metrics import accuracy_score
X_train_tensor = torch.Tensor(epochs_train.get_data())
y_train_tensor = torch.LongTensor(epochs_train.events[:, 2])
X_test_tensor = torch.Tensor(epochs_test.get_data())
y_test_tensor = torch.LongTensor(epochs_test.events[:, 2])
train_dataset = TensorDataset(X_train_tensor, y_train_tensor)
train_loader = DataLoader(train_dataset, batch_size=32, shuffle=True)
class EEGClassifier(nn.Module):
def __init__(self, input_size, num_classes):
super(EEGClassifier, self).__init__()
self.lstm = nn.LSTM(input_size, 256, batch_first=True)
self.dropout = nn.Dropout(0.5)
self.fc = nn.Linear(256, num_classes)
def forward(self, x):
out, _ = self.lstm(x)
out = self.dropout(out)
out = self.fc(out[:, -1, :])
return out
hidden_size1 = 200
hidden_size2 = 100
output_size = len(set(y_train_tensor.numpy()))
num_samples, num_channels, num_times = epochs_train.get_data().shape
num_features = num_channels * num_times
model = EEGClassifier(3000, output_size)
criterion = nn.CrossEntropyLoss()
optimizer = optim.Adam(model.parameters(), lr=0.001)
y_train_tensor -= 1
y_test_tensor -= 1
num_epochs = 10
for epoch in range(num_epochs):
for inputs, labels in train_loader:
optimizer.zero_grad()
outputs = model(inputs)
loss = criterion(outputs, labels)
loss.backward()
optimizer.step()
with torch.no_grad():
outputs = model(X_test_tensor)
_, predicted = torch.max(outputs, 1)
acc = accuracy_score(y_test_tensor.numpy(), predicted.numpy())
print("准确率: {}".format(acc))