使用ICA去除32通道腦電數(shù)據(jù)中的眼電

% 去除聽音樂1的腦電數(shù)據(jù)眼電

% Method: ICA

% 1.導入數(shù)據(jù) 32通道 1000hz采樣率

% 2. 數(shù)據(jù)預(yù)處理(50hz陷波+0.3hz-250hz帶通，可以自己改參數(shù)調(diào)整)

% 3. ICA(下載的FastICA_25)

% 4. 計算ICA的各ICs分量與原水平、垂直眼電信號相關(guān)性

% 5. 將眼電噪聲的ICs置零后重建成沒有眼電的腦電信號

% 注意代碼中設(shè)置斷點

% copyright：seu-ww

% 2019-05-07

clc;clear all;close all;

MUSIC1PATH=‘D:\Share\music_breath_2018\1-EEG\MUSIC2’;

SAVEPATH1 = ‘D:\Share\music_breath_2018\1-EEG\CLEAN_DATA\MUSIC2’;

data_list = dir(MUSIC1PATH);

eeg_number = length(data_list)-2;

fs = 1000;

for i =3:48 %45名被試的數(shù)據(jù)

tic;

filename = data_list(i).name;

datapath = [MUSIC1PATH,filename];

num_med = regexp(data_list(i).name,‘D’,‘split’);

subjectID = str2num(cell2mat(num_med(1)));

load(datapath);

data_org = data_music2([1:18,20:23,25:36]?;% M1、M2

%% Preprocessing(Notch + bandpass)

Wo = 50/(fs/2); BW = Wo/35;

[b,a] = iirnotch(Wo,BW);

data_notch= filtfilt(b,a,(double(data_org(:?))’)’;

[b,a] = butter(3,[0.3/(fs/2) 250/(fs/2)],‘bandpass’);

data_bp= filtfilt(b,a,data_notch’)’;%使用零相移數(shù)字濾波器(IIR)

L=length(data_bp);t=(0:(L-1))/fs;

%% ICA

level = 32; %ICA分層數(shù)，最高設(shè)為腦電的通道數(shù)

[Zica,A,W] = fastica(data_bp,‘numOfIC’, level);

plot_ICAs(Zica,fs,10);

%% Correlation with VEOG and HEOG

for k=1:(length(Zica(:,1)))

s=corrcoef(data_bp(33,:)’,Zica(k,:)’);

ss1(k) = s(2,1);

s=corrcoef(data_bp(34,:)’,Zica(k,:)’);

ss2(k) = s(2,1);

end

figure;plot(abs(ss1));title(‘各ICs與原水平眼電信號相關(guān)性’);

figure;plot(abs(ss2));title(‘各ICs與原垂直眼電信號相關(guān)性’);

%% Reconstruction

rec_no = [3，7] %這里輸入的是與原垂直眼電信號相關(guān)性最高的IC分量

Zica_del = Zica; Zica_eye = zeros(size(Zica,1),size(Zica,2));

Zica_del(rec_no,:) = 0; Zica_eye(rec_no,:) = Zica(rec_no,:);

data_rec = AZica;

data_clean = AZica_del;

data_eog = A*Zica_eye;

%去眼電前后對比(以FP1通道的數(shù)據(jù)為例)

figure;plot(t,data_bp(1,:),‘k’);title('32通道ICA去眼電前EEG信號 ');

figure;plot(t,data_clean(1,:),‘k’);title('32通道ICA去眼電后EEG信號 ');

data_music2=data_clean;

save([SAVEPATH1,filename],‘data_music2’); %保存數(shù)據(jù)

end

%畫ICA分解后ICs分量的函數(shù)plot_ICAs

% 畫出ICA 獨立分量

% Zica : ICs

% fs : fs

% d : offset

function plot_ICAs(Zica,fs,d)

level = length(Zica(:,1));

L = length(Zica(1,:));

t=(0:(L-1))/fs;

figure;

offs= 1:1:level;

d= 10;

for ics = 1:level

plot(t,Zica(ics,:) + offs(ics)d,‘k’);hold on;

end

ylim([0 d(level+1)]); title(‘Independent Components’);

IC = regexp(sprintf('IC-%02d ‘,[1:level]),’ ',‘split’);

set(gca,‘YTick’,[d:d:d*(level)],‘YTickLabel’,IC(1:end));

end

ICA分解后的獨立分量

3明顯為垂直眼電，7明顯為水平眼電。

去眼電前后的腦電信號對比：

總結(jié)

以上是生活随笔為你收集整理的matlab去eog,使用ICA去除32通道脑电数据中的眼电的全部內(nèi)容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯，歡迎將生活随笔推薦給好友。