Standard deviation in 3D array from 1st to nth number array

Hello I am trying to calculate residual noise from auditory brainstem responses (ABR).

The 3D array that I'm working in is in the size [2 x 493 x 6000] double which is a reflection of [Nchannels x Ntime x Nepochs].

I have 2 recording channels and each channel consists of 493 time data points and there are 6000 trials (or epochs).

Will call our 3D array ABR.

To calculate residual noise I have to:

1. calculate the std of the data points (in rows) from the 1st data point of epoch no.1 (which is given by ABR(:,1,1)) to the last data point of epoch no. X (which is given by ABR(:, end, X)).
2. The issue that I'm facing is how to specify calculation of std from 1st epoch to nth epoch.
3. I've tried std(ABR,0,2,[ABR(:,1,1), ABR(:,end,1)]); but this is not the answer.
4. I need to calcuate 6000 of these in the following format...
5. std(ABR,0,2,[ABR(:,1,1), ABR(:,end,1)]);
6. std(ABR,0,2,[ABR(:,1,1), ABR(:,end,2)]);
7. std(ABR,0,2,[ABR(:,1,1), ABR(:,end,3)]);
8. std(ABR,0,2,[ABR(:,1,1), ABR(:,end,4)]);
9. std(ABR,0,2,[ABR(:,1,1), ABR(:,end,5)]); ... and so on until the last last number becomes 6000.
10. In this case ABR is the data set, 0 is the weighting, 2 means calculation in rows and [] specifies what arrays (or epochs) to calculate std from.

This is what I am tryin to do graphically - using excel.

As you can see each std data point is the std result of all the epochs upto and including the position of the std.

What is the best way to do this?

[Adam Danz]

It sounds like you want the std of the second dimension, right?

data = rand(2, 493, 6000); 
sd = squeeze(std(data, [], 2));
size(sd)

ans = 1×2

2 6000

[MinChul Park]

Thanks for that suggestion Adam!

I've attached a more clear diagram of what I need my code to do.

I hope this is enough information for you to understand what my issue is.

Based on the below diagram I think the code should look something like this.

I've already tried this code and it does not work - I write this to hopefully demonstrate to you what I have in mind.

% for data in 3D = 6 (No. of channels) X 439 (No. of time points) X 6000 (No. of epochs)

Nepochs = 6000

for i = 1 : Nepochs

data(i) = std(data(:, :, 1:i), 3); % dimension value 3 added for std calculation over pages

end

Hopefully, the eventual code gives me data in 2D = 6 (No. of channels) X 6000 (No. of epochs) with each element being the std output as shown by the diagram below.

Hope this makes sense. What are your thoughts?

[MinChul Park]

So further to the point above I realised that I needed a vector dimension to specify that the calculation must be row then pages. That is to say...

ESqrt = sqrt(linspace (1,(Nepochs),(Nepochs))); % Nepochs = 6000 epochs

RN = zeros(Nchannels, Nepochs); % Nchannels = 6 recording channels

for i = 1 : Nepochs

RN(i) = std(data(:, :, 1:i), 0, [2 3]); % data = 6 X 493 X 6000 3D matrix.

end

Here the std(A, 0, [2 3]) was used and [2 3] was specified as the calculation dimension is across the matrix first then over the pages.

However I keep getting the error "Unable to perform assignment because the left and right sides have a different number of elements."

What am I doing wrong?

I thought by making "RN = zeros(Nchannels, Nepochs);" which is a 6 X 6000 matrix, the "for" loop would then automatically insert the output data (6 X 1 matrix) 6000 times eventually completing the 6 X 6000 matrix??

[Adam Danz]

You're getting that error because you're not indexing R correctly. The output of your std line will be a column vector. Here's how to correctly index that:

data = rand(6, 493, 6000); 
RN = zeros(size(data,[1,3]));
for i = 1 : size(data,3)
    RN(:,i) = std(data(:, :, 1:i), 0, [2,3]);
end
size(RN)

BTW, great methods image! 🙂

[MinChul Park]

Thanks for that reply Adam! 😀

I've figured everything out and it works well.

[Adam Danz]

Glad it worked out! I'll move my comment to the answers section so your quesitons is moved from the unanswered questions queue.

Answer 1

[Adam Danz]

To cumulatively compute the std across dimension 3,

data = rand(6, 493, 6000); 
RN = zeros(size(data,[1,3]));
for i = 1 : size(data,3)
    RN(:,i) = std(data(:, :, 1:i), 0, [2,3]);
end
size(RN)