How to correctly compute the EEG Frequency Bands with Python?

How to correctly compute the EEG Frequency Bands with Python?

so I am trying to compute the EEG25 channels, 512 sampling rate, 248832channel bands alpha, beta, gamma, etc. with Python. I managed to do so by: firstly filtering the signal with a butterworth filter that looks like this: def butter_bandpass_filterdata, lowcut, highcut, fs, order2: nyq 0.5 fs low lowcut nyq high highcutnyq b, a butterorder, low, high, btypeband printb,a y lfilterb, a, data return y I used the filter with the low set on 0.1 and high on 80. Than I compute the fft of the signal and store it in fft1, on which I use again the butterwort filter to extract the frequencies of each band and it looks like this: for i in np.arangen: alpha1 butter_bandpass_filterfft1i, :, 8.1, 12.0, 256 beta1 butter_bandpass_filterfft1i, :, 16.0, 36.0, 256 gamma1 butter_bandpass_filterfft1i, :, 36.1, 80, 256 delta1 butter_bandpass_filterfft1i, :, 0.0, 4.0, 256 sigma1 butter_bandpass_filterfft1i, :, 12.1, 16.0, 256 theta1 butter_bandpass_filterfft1i, :, 4.1, 8.0, 256 sumalpha1 sumabsalpha1 sumbeta1 sumabsbeta1 sumgamma1 sumabsgamma1 sumdelta1 sumabsdelta1 sumsigma1 sumabssigma1 sumtheta1 sumabstheta1 objects sumalpha1, sumbeta1, sumgamma1, sumdelta1, sumsigma1, sumtheta1 N lenobjects ra rangeN plt.titlesignal_labelsi plt.autoscale somestuffneeded np.arange6 ticks alpha,beta,gamma,delta,sigma,theta plt.xtickssomestuffneeded, ticks plt.barra, objects plt.show The problem is that I get a high value for gamma which represents a high cognitive activity, which the person clearly didnt have. The result for the first channel is: Can anyone point out a better solution has any ideas about what is wrong can tell me what is wrong Am I using the filter incorect Thanks

Ok, so for those interested, Ive computed the frequency bands of an eeg by using the butterworth filter described in the problem description. During the eeg analysis class I came to the conclusion that the frequency bands were computed from the fft of the eeg which was not enough because the fft should have been multiplied with its conjugate so here is the code in python which computes the total power, the relative and the absolute frequency bands. alpha1 np.zerosn, g.getNSamples0 beta1 np.zerosn, g.getNSamples0 gamma1 np.zerosn, g.getNSamples0 delta1 np.zerosn, g.getNSamples0 sigma1 np.zerosn, g.getNSamples0 theta1 np.zerosn, g.getNSamples0 I only used the channels with the following indexes for i in 1,5,10,15,18: note that fft1 is precomputed in another part of the code i used ccFFT1 fft1i,:np.conjugatefft1i,: the butterworth filter is used to compute each frequency band according to its featuring freq. alpha1 butter_bandpass_filterccFFT1:, 8.1, 12.0, 512 beta1 butter_bandpass_filterccFFT1:, 16.0, 36.0, 512 gamma1 butter_bandpass_filterccFFT1:, 36.1, 80, 512 delta1 butter_bandpass_filterccFFT1:, 0.0, 4.0, 512 sigma1 butter_bandpass_filterccFFT1:, 12.1, 16.0, 512 theta1 butter_bandpass_filterccFFT1:, 4.1, 8.0, 512 this sums frequencies for each band sumalpha1 sumabsalpha1 sumbeta1 sumabsbeta1 sumgamma1 sumabsgamma1 sumdelta1 sumabsdelta1 sumsigma1 sumabssigma1 sumtheta1 sumabstheta1 compute the total power totalPower sumalpha1sumbeta1sumgamma1sumdelta1sumsigma1sumtheta1 storing them in objects to plot this computes the relativ power objects sumalpha1totalPower, sumbeta1totalPower, sumgamma1totalPower, sumdelta1totalPower, sumsigma1totalPower, sumtheta1totalPower to plot the absolute power, dont divide each sum through the total power N lenobjects ra rangeN plt.titlestrsignal_labelsi relativ plt.autoscale somestuffneeded np.arange6 ticks alpha,beta,gamma,delta,sigma,theta plt.xtickssomestuffneeded, ticks plt.barra, objects plt.show Hope this helps I am glad to answer questions if there are any regarding this task of computing the frequency bands with python

I believe there is a much simpler way to do this with numpy.fft.rfft and numpy.fft.rfftfreq. In the below example, I have two seconds of random data between 0.0 and 100.0 sampled at 512 Hz. Im plotting the band amplitude below, but if you wanted power, it would just be the square of the values. I think the reason you saw overly large values for Gamma, was that your gamma range is much larger than the others, and youre taking the sum of all the values in that range. This is fine if you are going to compare one EEG sets gamma value to another relative amppower, but not if youre just looking at one set and want to compare gamma to say alpha. If you just want to know what EEG state the person is in, you could use np.max instead of np.mean. Then whichever band has the max will tell you about the persons state. import numpy as np fs 512 Sampling rate 512 Hz data np.random.uniform0, 100, 1024 2 sec of data bw 0.0-100.0 Get real amplitudes of FFT only in postive frequencies fft_vals np.absolutenp.fft.rfftdata Get frequencies for amplitudes in Hz fft_freq np.fft.rfftfreqlendata, 1.0fs Define EEG bands eeg_bands Delta: 0, 4, Theta: 4, 8, Alpha: 8, 12, Beta: 12, 30, Gamma: 30, 45 Take the mean of the fft amplitude for each EEG band eeg_band_fft dict for band in eeg_bands: freq_ix np.wherefft_freq eeg_bandsband0 fft_freq eeg_bandsband10 eeg_band_fftband np.meanfft_valsfreq_ix Plot the data using pandas here cause its easy import pandas as pd df pd.DataFramecolumnsband, val dfband eeg_bands.keys dfval eeg_band_fftband for band in eeg_bands ax df.plot.barxband, yval, legendFalse ax.set_xlabelEEG band ax.set_ylabelMean band Amplitude

I would use scipy.signal.welch for that: def calc_bands_powerx, dt, bands: from scipy.signal import welch f, psd welchx, fs1. dt power band: np.meanpsdnp.wheref lf f hf for band, lf, hf in bands.items return power Also, mne-python is a great package for EEGMEG analysis, its worth taking a look

Комментарии

Популярные сообщения из этого блога

Skipping acquire of configured file 'contrib/binary-i386/Packages' as repository … doesn't support architecture 'i386'

Connection string for MariaDB using ODBC

Celery like system based on django channels