Calculating a mean matrix from a large number of matrices

I have been trying to calculate an average matrix from 50 matrices. Using squeeze and cat, i have calculated the mean to two matrices. My question was; how is it possible to get the mean of 50 matrices. I have a code to open the files from a folder and reformat and bin them.

myFolder = uigetdir('C:\Users\c13459232\Documents\MATLAB'); 
if ~isdir(myFolder) 
errorMessage = sprintf('Error: the following folder does not exist: \n%s', myFolder); 
uiwait(warndlg(errorMessage)); 
return;
end
outFolder = fullfile(myFolder, 'output'); 
mkdir(outFolder); 
filePattern = fullfile(myFolder, '*.asc'); 
Files = dir(filePattern);  
for k = 1 : length(Files) 
baseFileName = Files(k).name;
FileName = fullfile(myFolder, baseFileName); 
fprintf(1, 'Now reading %s\n', FileName); 
fid = fopen(FileName); 
Cell = textscan( fid, '%d', 'delimiter', ';')
fclose(fid);
Data = cell2mat(Cell); 
N = 1024; 
skip = 2; 
Finish0 = reshape(Data, N, [])'; 
Finish1 = Finish0(1:skip:end, 1:skip:end);
newFileName = fullfile(outFolder, ['finish_' sprintf('%03d', k) '.txt']); 
dlmwrite(newFileName, Finish0, ';');
eval(['finish_' sprintf('%03d', k) ' = Finish1;']); 
disp([':: Now writing ' newFileName]); 
end

Where could i insert squeeze(cat()) to average all of these matrices. The example I found to get the mean of two matrices is:

c2 = squeeze(mean(cat(3,finish_002,finish_004), 3));

Would you suggest a variation on this code or is there a better way to get an average matrix?

[Stephen23]

Why not read and learn from the answers to your earlier questions, which told you that using awful eval will make your code slow, buggy, difficult to understand, and difficult to debug:

https://www.mathworks.com/matlabcentral/answers/323811-creating-multiple-equally-sized-matrices-from-a-single-numerical-cell

A much simpler, faster, and easier alternative is to use indexing. If you cannot load all of the data at once, then here are some alternatives:

• load each file in turn, use indexing to extract a subset of the loaded matrix, store this in an ND array, once all have been extracted take the mean. Then repeat for the next subset.
• use tall arrays to access data directly from the harddrive, use indexing to extract a subset of the loaded matrix, store this in an ND array, once all have been extracted take the mean. Then repeat for the next subset.

Whatever you do avoid eval: it will not make your life better or easier trying to solve this problem.

Rather than using cat, you can use indexing into a preallocated array. Then calculate the mean. Rinse and repeat.

[Aaron Smith]

Thanks. I only have the 2012 version of matlab, it's the last time my uni got matlab. Aren't Tall arrays only on versions since 2016?

[Stephen23]

@Aaron Smith: you did not much information in your question, but I assumed that:

1. you cannot load all matrices simultaneously.
2. you wish to calculate the mean along the third dimension (i.e. as if the matrices were concatenated along the third dimension) giving you one single matrix result with the same size as all of your input matrices.
3. all input matrices have the same size.

Are these correct?

Tall arrays were added in 2016b. If you are using an earlier version then you can use indexing, as my earlier comment suggests. Of course this will be slow.

[Aaron Smith]

Using eval I have all the matrices in the workspace and using code to create a cell array they are all saved in a cell array in the workspace. I wish to calculate the mean across all 50 matrices which is essentially calculating the mean across a certain stretch to time so yes, it would be along a third dimension. All the input matrices do have the same size.

[Stephen23]

If you can load all of the matrices into memory at once then you can load them into one ND array or cell array (using indexing). Then your task would be quite simple to solve (perhaps just one command).

Why make your life more difficult by using eval when simple indexing would make this trivially easy?

You picked a slow, buggy, obfuscated way to write code (and obviously did not read any of the links I have given you advising to avoid this awful way to write code), and you are now finding out how difficult it is to do anything with fifty separate variables in your workspace.

Two options:

1. load your data into one array using indexing, then your question can be resolved with perhaps just one command.
2. Continue writing pointlessly complicated code using eval, which makes accessing data difficult, and your code slow and hard to debug.

Which would you prefer?

[Aaron Smith]

You see, when each full matrix is saved as a single cell in a cell array, I'm not sure if this will give me the freedom i need to plot each one individually and to bin each one individually. Basically will i be able to get at the files as easily as i need. I tried extracting them from the cell array and looking at them individually and my difficulty with this is why i decided to go with the eval function. All the files being open is a little slow and it does clutter the workspace but it made the files easier to access and work with. That's why I chose to use it. I guess what i'm looking for is a way to write a line that'll basically get

mean( finish_001... finish_050)

If there isn't a way to write something like this I will have to change the way i've been working with the data.

[Aaron Smith]

I guess, When plotting I could create a loop so that it plots the first

contourf(finishCell{1,1})

And then for each successive navigation it will plot the next.

[Stephen23]

@Aaron Smith: using indexing is going to be much easier to work with.

If you want to know why, read that link I gave you.

[Guillaume]

"You see, when each full matrix is saved as a single cell in a cell array, I'm not sure if this will give me the freedom i need to plot each one individually and to bin each one individually"

It's the exact opposite. Having all the matrices as cells of a cell array gives you the freedom to treat each them individually any way you want. Use a fixed index in your cell array, and the freedom to operate on all them at once. Having individually named matrices you sacrifice the latter and don't gain anything in return.

As a rule, if you're numbering variables you're doing something wrong. These numbered variables are obviously related so they should be stored together in a container of some sort (matrix, cell array, containers.map, whatever).

This would make your question dead easy to solve:

%finish: cell array of matrices
meanofallmatrices = mean(cat(3, finish{:}), 3);  %no need for squeeze

And since all the matrices are obviously the same size, you could store them directly as a 3D array instead of cell array (with again no loss of freedom)

Stephen is right, eval is completely the wrong approach. Here be dragons!

Answer 1

[Thorsten]

In the for-loop, you can just use

if k == 1, 
  M = 1/length(Files)*Data;
else
  M = M + 1/length(Files)*Data;
end

Then after the loop, M is the mean of the matrices.

[Aaron Smith]

Where exactly in the for loop should I insert this if loop? I've tried placing it right before the end but the error that arises is where i'm a little hazy. You haven't used the i in your if loop and I also haven't used it in my original code so what is the i in your for statement for? I expect the answer will be a lot more simple than I am making it. Thanks

[Thorsten]

The i is the running index of the for loop, so here if should be k. I corrected it above. Place the code in the loop after you have determined "Data".

[Stephen23]

@Aaron Smith: this is a simple method that avoids needing to load all of the matrices at once (or use eval). You should use this method.

[Aaron Smith]

myFolder = uigetdir('C:\Users\c13459232\Documents\MATLAB'); % Generate command window to choose a folder
if ~isdir(myFolder) % if the directory is not a valid path
errorMessage = sprintf('Error: the following folder does not exist: \n%s', myFolder); % print this error message
uiwait(warndlg(errorMessage)); % block the execution of program and wait to resume
return;
end
outFolder = fullfile(myFolder, 'output'); % build full file name from parts in folder 'Output'
mkdir(outFolder); % create folder
filePattern = fullfile(myFolder, '*.asc'); % Call all files with '.asc' from the chosen folder
Files = dir(filePattern);  % list folder contents
finishCell = cell(length(Files), k);
for k = 1 : length(Files) % for all files files in the folder
baseFileName = Files(k).name;
FileName = fullfile(myFolder, baseFileName); %
fprintf(1, 'Now reading %s\n', FileName); % generates statements to indicate files are being read
fid = fopen(FileName); % open the file from chosen folder
Cell = textscan( fid, '%d', 'delimiter', ';'); % scanning data from files
fclose(fid); % close file from chosen folder
Data = cell2mat(Cell); % convert the cell data to matrix
N = 1024; % Number of numbers per row
skip = 2; % skip every second value to bin the matrix in half
Finish0 = reshape(Data, N, [])'; % reshape the data into the correct format
Finish1 = Finish0(1:skip:end, 1:skip:end);
finishCell{k} = Finish0;
newFileName = fullfile(outFolder, ['finish_' sprintf('%03d', k) '.txt']); % finish_001.asc, finish_002.asc, ...
dlmwrite(newFileName, Finish0, ';'); % create new files
eval(['finish_' sprintf('%03d', k) ' = Finish1;']); % create matrices in the workspace
disp([':: Now writing ' newFileName]); % generates statement to say the data is being written to file
end

This code is the original I had written using eval. I added the cell array lines to simply test how it works. This works completely fine, creates the 50 separate files (eval) and creates the cell array. Now, below is a new code I've written with eval removed and the if loop above added.

myFolder = uigetdir('C:\Users\c13459232\Documents\MATLAB'); % Generate command window to choose a folder
if ~isdir(myFolder) % if the directory is not a valid path
errorMessage = sprintf('Error: the following folder does not exist: \n%s', myFolder); % print this error message
uiwait(warndlg(errorMessage)); % block the execution of program and wait to resume
return;
end
 filePattern = fullfile(myFolder, '*.asc'); % Call all files with '.asc' from the chosen folder
Files = dir(filePattern);  % list folder contents
finishCell = cell(length(Files), k);
for k = 1 : length(Files) % for all files files in the folder
baseFileName = Files(k).name;
FileName = fullfile(myFolder, baseFileName); %
fid = fopen(FileName); % open the file from chosen folder
Cell = textscan( fid, '%d', 'delimiter', ';'); % scanning data from files
fclose(fid); % close file from chosen folder
Data = cell2mat(Cell); % convert the cell data to matrix
N = 1024; % Number of numbers per row
Finish0 = reshape(Data, N, [])'; % reshape the data into the correct format
finishCell{k} = Finish0;
if k == 1,
M = 1/length(Files)*Finish0;
else
M = M + 1/length(Files)*Finish0;
end
end

I'm not sure if I have written the new code correctly but I can fix that later. The problem is an error occurred which said "undefined function or variable k". After trying this new code, this even began happening with the original code which worked fine previously. Where might this problem with k have cone from?

[Guillaume]

finishCell = cell(length(Files), k);
for k = ...

You use the variable k before you create it. If your code is a script, it would have worked before if you had a leftover k from previous runs. If the code is a function, it would never have worked.

In any case, I don't see the purpose of the k in the offending line. You only want one column anyway, so:

finishCell = cell(length(Files), 1);  %And I would use numel everywhere instead of length

%or possibly better:
finishCell = cell(size(Files));

And since you're now storing the matrices in a cell array, I wouldn't bother calculating the mean iteratively in the loop. I would calculate it after the loop. This will be more accurate, you'll get accumulation errors with your loop approach. As per my comment to the question, getting the mean after the loop is a one-liner:

M = mean(cat(3, finishCell{:}), 3);

[Aaron Smith]

The problem must be the k before the loop. I must have never tried the code to create a cell array with an empty workspace and k still had a value assigned from the previous code. That's about the greenest of rookie mistakes.

myFolder = uigetdir('C:\Users\c13459232\Documents\MATLAB'); % Generate command window to choose a folder
if ~isdir(myFolder) % if the directory is not a valid path
errorMessage = sprintf('Error: the following folder does not exist: \n%s', myFolder); % print this error message
uiwait(warndlg(errorMessage)); % block the execution of program and wait to resume
return;
end
 filePattern = fullfile(myFolder, '*.asc'); % Call all files with '.asc' from the chosen folder
Files = dir(filePattern);  % list folder contents
finishCell = cell(length(Files));
for k = 1 : length(Files) % for all files files in the folder
baseFileName = Files(k).name;
FileName = fullfile(myFolder, baseFileName); %
fid = fopen(FileName); % open the file from chosen folder
Cell = textscan( fid, '%d', 'delimiter', ';'); % scanning data from files
fclose(fid); % close file from chosen folder
Data = cell2mat(Cell); % convert the cell data to matrix
N = 1024; % Number of numbers per row
Finish0 = reshape(Data, N, [])'; % reshape the data into the correct format
finishCell{k} = Finish0;
if k == 1,
M = 1/length(Files)*Finish0;
else
M = M + 1/length(Files)*Finish0;
end
end

This code doesn't present any errors but I'm fairly sure the way i'm using finish0 inside my if loop is wrong though.

[Guillaume]

if you have N files, cell(length(Files)) will create a NxN cell array, of which you will only fill the first row. It's either:

finishCell = cell(length(Files), 1);
%or
finishCell = cell(numel(Files), 1); 
%or
finishCell = cell(size(Files));

to create a Nx1 cell array.

I can't see anything wrong with the if block (which is not a loop) other than it's inefficient and introduces accumulation error as previously stated. It will produce more or less the correct result. If you really insist on calculating the mean iteratively, I would at least perform the division only once, at the end.

To see what I mean about accumulation error:

>>v = ones(1, 30); %the mean of this is obviously 1
>>meanyourway = sum(v/numel(v));  %each element is divided by the count and then summed
>>meantheproperway = sum(v)/numel(v);  %sum each element then divide at the end, or simply use the mean function
>>1 - meantheproperway
ans = 
     0
>>1 - meanyourway
ans = 
     1.1102e-16

[Aaron Smith]

Thanks Guillaume. I had realized the problem with my k before the loop and added my second comment before I saw your comment was posted. The following error occurs when i use the one line mean calculation:

Mean = mean(cat(3, finishCell{:}, 3));
Error using cat
CAT arguments dimensions are not consistent.

From what I can tell from the mathworks page on cat, the first 3 determines the dimension along which the matrices are concatenated. So what is the second 3 for? Is that a specification of the number of matrices? I have tried replacing it with k which generates the same error and I've tried removing it which gives a result, but only the mean values of the columns. I need a 1024 x 1024 mean of all 50 matrices. Essentially (cell1 + cell2 +... cell50)/50

Thanks again

[Guillaume]

You've got a closing bracket in the wrong location, the second 3 is an argument to mean, not cat.

mean(cat(3, finishCell{:}), 3);
%bracket goes here ------^

That second 3 tells mean to average across pages instead of the default across rows as you saw.

[Aaron Smith]

Ah, I see.

Thanks a million Guillaume