spectral features extraction from eeg

5 views (last 30 days)
Tarek Lajnef
Tarek Lajnef on 11 Jun 2012
Answered: BOMMALA SILPA on 15 Nov 2021
hi evrey body, i am novice in signal proccessing and matlab scripts, i try to write a function to extract relative power of the eeg frequencies bands from the PSD, the code i used is the next but i am not sure if it is rigth, can someone help me to correct it, thanks
function [p_tot Prdelta Prtheta Pralpha Prsigma Prbeta]=freq_feat(data,debut,fin,Fs,win)
epoch=data(debut:fin);
%%%%%calcul de la densité spéctrale avec la methode de welch
[Pxx,f]= pwelch(epoch,win,10,1024,250);
figure;
Hy=spectrum.welch;
psd(Hy,epoch,'Fs',Fs,'NFFT',1024);
%%%%%%%%%extraction des caractéristiques
%%%%puissance total
p_tot=sum(Pxx);
%%%%%puissnance relative bande delta
%Filtrage bande delta
n=10;
Wn1=[0.5 3]/50;
[b1,a1] = butter(n,Wn1,'bandpass');
delta = filter(b1,a1,epoch);
%%%%DSP delta
[Pdelta,fd] =pwelch(delta);
%%%puissance relative
Ptot_delta =sum(Pdelta);
Prdelta= Ptot_delta./p_tot;
%%%%puissnance relative bande theta
%Filtrage bande theta
n=10;
Wn2=[4 7]/50;
[b2,a2] = butter(n,Wn2,'bandpass');
theta = filter(b2,a2,epoch);
%%%%DSP delta
[Ptheta,ft] = pwelch(theta);
%%%puissance relative
Ptot_theta =sum(Ptheta);
Prtheta= Ptot_theta./p_tot;
%%%%%puissnance relative bande alpha
%Filtrage bande alpha
n=10;
Wn3=[8 12]/50;
[b3,a3] = butter(n,Wn3,'bandpass');
alpha = filter(b3,a3,epoch);
%%%%DSP alpha
[Palpha,fa] = pwelch(alpha);
%%%puissance relative
Ptot_alpha =sum(Palpha);
Pralpha= Ptot_alpha./p_tot;
%%%%%puissnance relative bande sigma
%Filtrage bande sigma
n=10;
Wn4=[12 14]/50;
[b4,a4] = butter(n,Wn4,'bandpass');
sigma = filter(b4,a4,epoch);
%%%%DSP sigma
[Psigma,fsig] = pwelch(sigma);
%%%puissance relative sigma
Ptot_sigma =sum(Psigma);
Prsigma= Ptot_sigma./p_tot;
%%%%%puissnance relative bande beta
%Filtrage bande beta
n=10;
Wn5=[13 22]/50;
[b5,a5] = butter(n,Wn5,'bandpass');
beta = filter(b5,a5,epoch);
%%%%DSP beta
[Pbeta,fb] = pwelch(beta);
%%%puissance relative
Ptot_beta =sum(Pbeta);
Prbeta= Ptot_beta./p_tot;

Answers (2)

Kosai
Kosai on 28 Jun 2012
you can do the same JOB easily with wavelet transformation look at wavelet toolbox and Wavelet - Decomposition Level ...

BOMMALA SILPA
BOMMALA SILPA on 15 Nov 2021
Hello iam silpa. Iam also searching for code of PSD of EEG frequency rythms can you help me

Categories

Find more on EEG/MEG/ECoG in Help Center and File Exchange

Tags

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!