How can I condense this block of code?

2 views (last 30 days)
CONTEXT: I have EEG data in a matrix. I have parsed through the data and separated it into several cell arrays of smaller matrices based on behavioral time stamps. For example, EEG data is 500,000 points long and 4 channels. One cell array might have 30 cells; in one of these cells, I might have a matrix of 4 channels of data, each only 1000 points long. I am then looping over these arrays separately and performing the same operations on them. These seems bulky and inefficient. Can someone let me know if this code can be condensed and if it can, advice? For clarity, one of the cell arrays is called restLFPs, each loop is exactly the same.
CODE:
% Rest
thetaPowerRest = NaN(length(restLFPs),1);
betaPowerRest = NaN(length(restLFPs),1);
for i = 1:length(restLFPs)
% Obtain Fourier coefficients and normalize them for the signal length
% and consequently, number of FFT points
ft = fft(restLFPs);
nft = 2.*(ft(1:floor(length(ft)/2),:).*conj(ft(1:floor(length(ft)/2),:)))./length(ft)^2;
% Create frequency vector using sampling frequency and signal
% length/NFFT
frex = [0:length(ft)-1]*(downsamprate/length(ft));
frex = frex(1:floor(length(ft)/2));
% Theta
Hz3 = dsearchn(frex',3);
Hz8 = dsearchn(frex',8);
thetaPower = mean(mean(nft(Hz3:Hz8,:),2));
thetaPowerRest(i,1) = thetaPower;
% Beta
Hz10 = dsearchn(frex',10);
Hz30 = dsearchn(frex',30);
betaPower = mean(mean(nft(Hz10:Hz30,:),2));
betaPowerRest(i,1) = betaPower;
end
clear thetaPower betaPower ft nft frex Hz3 Hz8 Hz10 Hz30 i
% Control(All Congruent)
thetaPowerControl = NaN(length(controlLFPs),1);
betaPowerControl = NaN(length(controlLFPs),1);
for i = 1:length(controlLFPs)
% Obtain Fourier coefficients and normalize them for the signal length
% and consequently, number of FFT points
ft = fft(controlLFPs{i});
nft = 2.*(ft(1:floor(length(ft)/2),:).*conj(ft(1:floor(length(ft)/2),:)))./length(ft)^2;
% Create frequency vector using sampling frequency and signal
% length/NFFT
frex = [0:length(ft)-1]*(downsamprate/length(ft));
frex = frex(1:floor(length(ft)/2));
% Theta
Hz3 = dsearchn(frex',3);
Hz8 = dsearchn(frex',8);
thetaPower = mean(mean(nft(Hz3:Hz8,:),2));
thetaPowerControl(i,1) = thetaPower;
% Beta
Hz10 = dsearchn(frex',10);
Hz30 = dsearchn(frex',30);
betaPower = mean(mean(nft(Hz10:Hz30,:),2));
betaPowerControl(i,1) = betaPower;
end
clear thetaPower betaPower ft nft frex Hz3 Hz8 Hz10 Hz30 i
% Congruent
thetaPowerCongruent = NaN(length(congruentLFPs),1);
betaPowerCongruent = NaN(length(congruentLFPs),1);
for i = 1:length(congruentLFPs)
% Obtain Fourier coefficients and normalize them for the signal length
% and consequently, number of FFT points
ft = fft(congruentLFPs{i});
nft = 2.*(ft(1:floor(length(ft)/2),:).*conj(ft(1:floor(length(ft)/2),:)))./length(ft)^2;
% Create frequency vector using sampling frequency and signal
% length/NFFT
frex = [0:length(ft)-1]*(downsamprate/length(ft));
frex = frex(1:floor(length(ft)/2));
% Theta
Hz3 = dsearchn(frex',3);
Hz8 = dsearchn(frex',8);
thetaPower = mean(mean(nft(Hz3:Hz8,:),2));
thetaPowerCongruent(i,1) = thetaPower;
% Beta
Hz10 = dsearchn(frex',10);
Hz30 = dsearchn(frex',30);
betaPower = mean(mean(nft(Hz10:Hz30,:),2));
betaPowerCongruent(i,1) = betaPower;
end
clear thetaPower betaPower ft nft frex Hz3 Hz8 Hz10 Hz30 i
% Incongruent
thetaPowerIncongruent = NaN(length(incongruentLFPs),1);
betaPowerIncongruent = NaN(length(incongruentLFPs),1);
for i = 1:length(incongruentLFPs)
% Obtain Fourier coefficients and normalize them for the signal length
% and consequently, number of FFT points
ft = fft(incongruentLFPs{i});
nft = 2.*(ft(1:floor(length(ft)/2),:).*conj(ft(1:floor(length(ft)/2),:)))./length(ft)^2;
% Create frequency vector using sampling frequency and signal
% length/NFFT
frex = [0:length(ft)-1]*(downsamprate/length(ft));
frex = frex(1:floor(length(ft)/2));
% Theta
Hz3 = dsearchn(frex',3);
Hz8 = dsearchn(frex',8);
thetaPower = mean(mean(nft(Hz3:Hz8,:),2));
thetaPowerIncongruent(i,1) = thetaPower;
% Beta
Hz10 = dsearchn(frex',10);
Hz30 = dsearchn(frex',30);
betaPower = mean(mean(nft(Hz10:Hz30,:),2));
betaPowerIncongruent(i,1) = betaPower;
end
clear thetaPower betaPower ft nft frex Hz3 Hz8 Hz10 Hz30 i

Accepted Answer

Chunru
Chunru on 18 Aug 2021
Consider this:
% Rest/Control/Congruent/Incongruent
thetaPower = NaN(length(restLFPs), 4);
betaPower = NaN(length(restLFPs), 4);
for j = 1:4
for i = 1:length(restLFPs)
% Obtain Fourier coefficients and normalize them for the signal length
% and consequently, number of FFT points
% The switch structure can be simplified if input is a 2D cell
% array
switch j
case 1
ft = fft(restLFPs);
case 2
ft = fft(controlLFPs{i});
case 3
ft = fft(congruentLFPs{i});
case 4
ft = fft(incongruentLFPs{i});
end
nft = 2.*(ft(1:floor(length(ft)/2),:).*conj(ft(1:floor(length(ft)/2),:)))./length(ft)^2;
% Create frequency vector using sampling frequency and signal
% length/NFFT
frex = [0:length(ft)-1]*(downsamprate/length(ft));
frex = frex(1:floor(length(ft)/2));
% Theta
Hz3 = dsearchn(frex',3);
Hz8 = dsearchn(frex',8);
thetaPower = mean(mean(nft(Hz3:Hz8,:),2));
thetaPower(i,j) = thetaPower;
% Beta
Hz10 = dsearchn(frex',10);
Hz30 = dsearchn(frex',30);
betaPower = mean(mean(nft(Hz10:Hz30,:),2));
betaPower(i,j) = betaPower;
end
end
  2 Comments
Anas Khan
Anas Khan on 18 Aug 2021
Dude.. thank you so much! My first instinct was that switch-case command statements would be very useful here but I couldn’t wrap my head around how to implement it. This is great!
Anas Khan
Anas Khan on 18 Aug 2021
Edited: Anas Khan on 18 Aug 2021
This won’t work perfectly as is because the length of the “-LFPs” cell arrays are not all equal(in fact restLFPs is not a cell array at all-just a matrix), but I think with some adjustments, this is what I needed.

Sign in to comment.

More Answers (0)

Categories

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

Products


Release

R2021a

Community Treasure Hunt

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

Start Hunting!