Reputation: 139
I'm trying to create a matrix that contains the averages of the kxk submatrices of a larger nxn matrix, where n is divisible by k. I can accomplish this fairly efficiently with something like this
mat = mat2cell(mat, k*ones(1,n/k), k*ones(1,n/k))
mat = cellfun(@mean,mat,'UniformOutput',false);
mat = cellfun(@mean,mat,'UniformOutput',false); %repeated to collapse cells to 1x1
mat = cell2mat(mat)
However, since I have a very large amount of data all in very large matrices, repeating this process can still take a long time even on a cluster, and combining the matrices is not an option due to memory limitations. I'm wondering if it's possible to rewrite this code using arrayfun
instead so I can utilize its GPU capability (as GPU cannot process cell arrays), but I'm running into problems since the helper function
function avg = blockavg(mat,i,j,k)
i1 = (i-1)*k+1;
i2 = i*k;
j1 = (j-1)*k+1;
j2 = j*k
avg = mean(mean(mat(i1:i2,j1:j2)));
end
calls for two broadcast variables mat
, an nxn matrix, and k
, a scalar, and these are not array inputs when plugged into arrayfun
. When running
ii = 1:(n/k);
jj = 1:(n/k);
matavg = arrayfun(@blockavg,mat,ii,jj,k)
An error message returns stating that input arguments must be arrays of the same size and shape, as only ii
and jj
are the array inputs, while mat
and k
do not vary across elements. I'm not quite sure how to work around this issue, or even if this type of operation is possible at all for arrayfun
. Any suggestion is appreciated.
Thanks!
Upvotes: 2
Views: 1173
Reputation: 36710
Another nice possibility without using the image processing toolbox:
squeeze(mean(mean(reshape(X,k,n/k,k,n/k),1),3))
This code reshapes the matrix in a way, that each block can be accessed by Y(a,:,b,:), which means they have identical first and third dimension index. Then mean is applied to first and third dimension, finally singleton dimensions are removed to get a 2d-result.
Upvotes: 2
Reputation: 36710
What you are trying to do is exactly what blockproc
from the image processing toolbox provides:
blockproc(X,[k k],@(x)(mean(mean(x.data))))
To speed up, you may try setting the UseParallel
-parameter.
Upvotes: 2