-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathfooof_for_eeg.py
165 lines (128 loc) · 5.68 KB
/
fooof_for_eeg.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
# -*- coding: utf-8 -*-
"""fooofas.ipynb
Automatically generated by Colaboratory.
Original file is located at
https://colab.research.google.com/drive/1umSeiaKx3GlOcE6MNccVSGuhuuCvybyH
"""
pip install fooof
from google.colab import drive
drive.mount('/content/drive')
import numpy as np
from scipy.io import loadmat, savemat
import glob
import pandas as pd
import os
import seaborn as sns
import matplotlib.pyplot as plt
from fooof import FOOOF
from fooof.bands import Bands
from fooof.sim.params import param_sampler
from fooof.sim.gen import gen_group_power_spectra
from fooof.sim.utils import set_random_seed
from fooof.analysis import get_band_peak_fm, get_band_peak_fg
from fooof import FOOOFGroup
from fooof.plts.annotate import plot_annotated_model
from fooof.plts.periodic import plot_peak_fits, plot_peak_params
from fooof.plts.aperiodic import plot_aperiodic_params, plot_aperiodic_fits
from fooof.utils.reports import methods_report_info
import numbers
# The first step is to load the power spectrum in two dimensions - one axis must contain the frequencies (for example from 0 to 30 every 0.5 Hz) and the other
# axis must contain the averaged powers of each electrode (obtained by applying the Fourier transform). The arrays must coincide. Lone PSD analysis.
data = loadmat('/content/drive/MyDrive/fooof/pirmas_powerintas.mat')
# Unpack data from dictionary, and squeeze numpy arrays
freqs = np.squeeze(data['freq'])
psd = np.squeeze(data['power'])
fm = FOOOF()
fm.report(freqs, psd, [0, 61])
# Extract FOOOF results from object
fooof_results = fm.get_results()
# Convert FOOOF results to a dictionary
# This is useful for saving out as a mat file
fooof_results_dict = fooof_results._asdict()
# Save FOOOF results out to a mat file
savemat('/content/drive/MyDrive/igb/fooof_results.mat', fooof_results_dict)
# Power spectral density parameters:
print('Frequency Range: \t', fm.freq_range)
print('Frequency Resolution: \t', fm.freq_res)
print('Frequency Values: \t', fm.freqs[0:5])
print('Power Values: \t\t', fm.power_spectrum[0:5])
# Print out model fit results
print('aperiodic params: \t', fm.aperiodic_params_)
#
# Aperiodic parameters: Offset (OFF), Knee (KNE) ir Exponent (EXP).
#
print('peak params: \t', fm.peak_params_)
print('r-squared: \t', fm.r_squared_)
print('fit error: \t', fm.error_)
print('fooofed spectrum: \t', fm.fooofed_spectrum_[0:5])
print('n peaks: \t', fm.n_peaks_)
# Power spectral density with aperiodic parameters plotted - OFF ir EXP.
model = plot_annotated_model(fm, plt_log=False)
model
# Representation of all the aperiodic components under study
aps1 = fg.get_params('aperiodic_params')
plot_aperiodic_params(aps1)
print(aps1)
aperiodic = pd.DataFrame(aps1)
print(aperiodic)
# Save aperiodic parameters of all patients
aperiodic.to_csv('/content/drive/MyDrive/fooof/off ir exp.csv')
# Save all alpha center frequencies of all patients
alphas = get_band_peak_fg(fg, bands.alpha)
alphasdf = pd.DataFrame(alphas)
alphasdf.to_csv('/content/drive/MyDrive/powerinta data/alpha center frequencies.csv')
# Representing aperiodic fits of all patients
freq_range = [1, 30]
plot_aperiodic_fits(aps1, freq_range, control_offset=True, log_freqs=False)
fm.save_report('FOOOF_report', '/content/drive/MyDrive/igb')
# Creating frequency bands
bands = Bands({'delta' : [0, 4],
'theta' : [4, 8],
'alpha' : [8, 15],
'beta' : [15, 30],
'gamma': [30, 100]})
# Center Frequency (CF)
# The peak frequency of identified peaks. The CF is a peak parameter, as part of the periodic component of the data.
#
# Power (PW)
# The power, over and above the aperiodic component, of identified peaks.
# The PW is a peak parameter, as part of the periodic component of the data.
#
# Bandwidth (BW)
# The bandwidth of identified peaks. The BW is a peak parameter, as part of the periodic component of the data.
#
alpha = get_band_peak_fm(fm, bands.alpha)
print(alpha)
# Multi power spectral density analysis - create an array with all patients average electrode powers and frequencies.
data = loadmat('/content/drive/MyDrive/fooof/visu poweriai powerinti.mat')
# Unpack data from dictionary, and squeeze numpy arrays
freqs = np.squeeze(data['freq']).astype('float')
psds = np.squeeze(data['power']).astype('float')
fg = FOOOFGroup()
fg.report(freqs, psds, [1, 100])
# It is possible to store a specific FOOOF measure that we are interested in, e.g. slope
exps = fg.get_params('aperiodic_params', 'exponent')
savemat('exps.mat', {'exps' : exps})
methods_report_info(fg)
r2s = pd.DataFrame(fg.get_params('r_squared'))
print(r2s)
r2s.to_csv('/content/drive/MyDrive/powerinta data/r_squared.csv')
expas = pd.DataFrame(fg.get_params('aperiodic_params', 'offset'))
print(expas)
expas.to_csv('/content/drive/MyDrive/powerinta data/offset.csv')
# Power points for individual subjects - csv files. Transpose, Lone PSD for analysis.
os.chdir('/content/drive/MyDrive/power exceliai')
extension = 'csv'
allfilenames = [i for i in glob.glob('*.{}'.format(extension))]
df = pd.concat([pd.read_csv(f) for f in allfilenames], axis = 1)
dff = df.loc[:, df.columns!='Hz\subj']
pirmas = pirmas[:-1]
antras = pirmas.groupby(by=pirmas.columns, axis=1).apply(lambda g: g.mean(axis=1) if isinstance(g.iloc[0,0], numbers.Number) else g.iloc[:,0])
trecias = np.square(antras)
ketvirtas = np.transpose(trecias)
ketvirtas.to_csv('/content/drive/MyDrive/fooof/visu_poweriai_powerinti.csv')
# Data frame for ARSQ questionnaire answers and periodic and aperiodic parameter values - a frame for a heatmap.
arsq = pd.read_excel('/content/drive/MyDrive/fooof/arsq ir parametrai be out powerinti.xlsx', sheet_name = "Sheet1")
df_arsq_subject = pd.DataFrame(arsq)
df_arsq = df_arsq_subject.drop(df_arsq_subject.iloc[:, 0:1],axis = 1)
df_arsq