Preprocessing_and_ICA_all_subject¶
Preprocessing and ICA
import numpy as np
import mne
import os
from mne.channels import make_standard_montage
from mne.io import read_raw_edf
from mne.viz import plot_events
from mne import find_events
from mne.preprocessing import (ICA, create_eog_epochs, create_ecg_epochs,
corrmap)
# define filename
x = "/mnt/c/Users/anna-/Desktop/HC/s01.edf"
os.path.basename(x).split(".")[0]
# loading files
files = ["/mnt/c/Users/anna-/Desktop/HC/s01.edf",
"/mnt/c/Users/anna-/Desktop/HC/s02.edf",
"/mnt/c/Users/anna-/Desktop/HC/s03.edf",
"/mnt/c/Users/anna-/Desktop/HC/s04.edf",
"/mnt/c/Users/anna-/Desktop/HC/s05.edf",
"/mnt/c/Users/anna-/Desktop/HC/s06.edf",
"/mnt/c/Users/anna-/Desktop/HC/s07.edf",
"/mnt/c/Users/anna-/Desktop/HC/s08.edf",
"/mnt/c/Users/anna-/Desktop/HC/s09.edf",
"/mnt/c/Users/anna-/Desktop/HC/s10.edf",
"/mnt/c/Users/anna-/Desktop/HC/s11.edf",
"/mnt/c/Users/anna-/Desktop/HC/s12.edf",
"/mnt/c/Users/anna-/Desktop/HC/s13.edf",
"/mnt/c/Users/anna-/Desktop/HC/s14.edf",
"/mnt/c/Users/anna-/Desktop/HC/h01.edf",
"/mnt/c/Users/anna-/Desktop/HC/h02.edf",
"/mnt/c/Users/anna-/Desktop/HC/h03.edf",
"/mnt/c/Users/anna-/Desktop/HC/h04.edf",
"/mnt/c/Users/anna-/Desktop/HC/h05.edf",
"/mnt/c/Users/anna-/Desktop/HC/h06.edf",
"/mnt/c/Users/anna-/Desktop/HC/h07.edf",
"/mnt/c/Users/anna-/Desktop/HC/h08.edf",
"/mnt/c/Users/anna-/Desktop/HC/h09.edf",
"/mnt/c/Users/anna-/Desktop/HC/h10.edf",
"/mnt/c/Users/anna-/Desktop/HC/h11.edf",
"/mnt/c/Users/anna-/Desktop/HC/h12.edf",
"/mnt/c/Users/anna-/Desktop/HC/h13.edf",
"/mnt/c/Users/anna-/Desktop/HC/h14.edf"]
# creating an own template as reference for ICA
template = np.array([ 0.94491604, 0.2894043 , -0.01358259, -0.10380207, -0.12087584,
1.06496116, 0.4053908 , -0.02073163, -0.07473682, -0.10292878,
0.2996572 , 0.04038614, -0.04146407, 0.31626096, 0.09097505,
-0.02723849, 0.3196638 , 0.04634076, -0.00174722])
# first ICA with "Fp1" and "Fp2" as reference (not so nice output -> using a template)
import numpy as np
import mne
import os
from mne.channels import make_standard_montage
from mne.io import read_raw_edf
from mne.viz import plot_events
from mne import find_events
from mne.preprocessing import (ICA, create_eog_epochs, create_ecg_epochs,
corrmap)
mne.sys_info()
for file in files:
data = read_raw_edf(file)
data.info
montage = make_standard_montage(kind="standard_1020")
data.set_montage(montage)
#data.plot_sensors();
#data.plot_psd(average=True);
#data.plot_psd(average=True,fmax=50);
#data.plot_psd(average=False, fmax=50);
#%matplotlib notebook
data.plot(n_channels=19, scalings=dict(eeg=1e-4));
raw = data
raw.load_data()
raw = raw.filter(l_freq=0.5,
h_freq=45.0,
filter_length='auto',
l_trans_bandwidth='auto',
h_trans_bandwidth='auto',
method='fir',
phase='zero',
fir_window='hamming',
fir_design='firwin')
raw = raw.set_eeg_reference(ref_channels="average")
eog_evoked = create_eog_epochs(raw, ch_name = ["Fp1","Fp2"]).average()
eog_evoked.apply_baseline(baseline=(None, -0.2))
eog_evoked.plot_joint();
filt_raw = raw.copy()
filt_raw = filt_raw.filter(l_freq=1.0, h_freq=None)
ica = ICA(n_components=9, max_iter='auto', method = "fastica", random_state=97)
ica.fit(filt_raw)
ica.plot_components();
#ica.plot_sources(raw);
ica.exclude = []
# find which ICs match the EOG pattern
eog_indices, eog_scores = ica.find_bads_eog(raw, ch_name="Fp1", measure="correlation", threshold=0.75)
ica.exclude = eog_indices
ica.plot_scores(eog_scores);
ica.plot_sources(eog_evoked);
# ica.apply() changes the Raw object in-place, so let's make a copy first:
reconst_raw = raw.copy()
ica.apply(reconst_raw)
raw.plot(n_channels=len(raw),
show_scrollbars=False);
reconst_raw.plot(n_channels=len(raw),
show_scrollbars=False);
reconst_raw.save("/mnt/c/Users/anna-/Desktop/reconst_raw/"+os.path.basename(file).split(".")[0]+"-raw.fif")
del reconst_raw;
# second ICA with template as reference
# position of the eeg sensors on the head
montage = make_standard_montage(kind="standard_1020")
#
for file in files:
subject = os.path.basename(file)
raw = read_raw_edf(file, preload=True)
raw.set_montage(montage)
# interpolate bad channels
if subject == "s12.edf":
raw.info["bads"] = ["T5"]
elif subject == "h10.edf":
raw.info["bads"] = ["Fp1", "Fz"]
if raw.info["bads"]:
raw.interpolate_bads(mode="accurate")
raw = raw.filter(l_freq=0.5,
h_freq=45.0,
filter_length='auto',
l_trans_bandwidth='auto',
h_trans_bandwidth='auto',
method='fir',
phase='zero',
fir_window='hamming',
fir_design='firwin')
raw.plot(n_channels=19, scalings=dict(eeg=1e-4), start=150.0);
filt_raw = raw.copy()
filt_raw = filt_raw.filter(l_freq=1.0, h_freq=None)
ica = ICA(n_components=7, max_iter='auto', method = "fastica", random_state=97)
ica.fit(filt_raw, reject=dict(eeg=250e-6))
ica.plot_components();
corrmap(icas=[ica], template=template, threshold=0.85, label="blink")
ica.exclude = ica.labels_["blink"]
reconst_raw = raw.copy()
ica.apply(reconst_raw)
reconst_raw.plot(n_channels=19, scalings=dict(eeg=1e-4), start=150.0);
reconst_raw.save("/mnt/c/Users/anna-/Desktop/reconst_raw/"+os.path.basename(file).split(".")[0]+"-raw.fif", overwrite=True)
del reconst_raw;
#check epochs of one example subject
file = "/mnt/c/Users/anna-/Desktop/reconst_raw/s06-raw.fif"
from mne.io import read_raw_fif
from mne import Epochs
from mne import make_fixed_length_epochs
data = read_raw_fif(file, preload=True)
epochs = make_fixed_length_epochs(data, duration=2.0, preload=True)
epochs.plot_psd_topomap(ch_type="eeg", normalize=True, dB=True, vlim=(0.1,0.8),);