1

I have a 3 dimensional (or higher) array that I want to aggregate by another vector. The specific application is to take daily observations of spatial data and average them to get monthly values. So, I have an array with dimensions <Lat, Lon, Day> and I want to create an array with dimensions <Lat, Lon, Month>.

Here is a mock example of what I want. Currently, I can get the correct output using a loop, but in practice, my data is very large, so I was hoping for a more efficient solution than the second loop:

% Make the mock data
A = [1 2 3; 4 5 6];
X = zeros(2, 3, 9);
for j = 1:9
  X(:, :, j) = A;
  A = A + 1;
end

% Aggregate the X values in groups of 3 -- This is the part I would like help on
T = [1 1 1 2 2 2 3 3 3];
X_agg = zeros(2, 3, 3);
for i = 1:3
  X_agg(:,:,i) = mean(X(:,:,T==i),3);
end

In 2 dimensions, I would use accumarray, but that does not accept higher dimension inputs.

4

1 回答 1

0

在得到你的答案之前，让我们首先以更一般的方式重写你的代码：

ag = 3; % or agg_size

X_agg = zeros(size(X)./[1 1 ag]);

for i = 1:ag
  X_agg(:,:,i) = mean(X(:,:,(i-1)*ag+1:i*ag), 3);
end

为了避免使用for循环，一个想法是将reshape您的 X 矩阵转换为可以mean直接使用该函数的东西。

splited_X = reshape(X(:), [size(X_agg), ag]);

所以现在splited_X(:,:,:,i)是包含所有应该聚合的矩阵的第 i 个部分X(:,:,(i-1)*ag+1:i*ag))（如上）

现在您只需要在 的第三维中找到平均值splited_X：

temp = mean(splited_X, 3);

然而，这会产生一个 4D 矩阵（其第 3 维大小为 1）。您可以使用以下函数再次将其转换为 3D 矩阵reshape：

X_agg = reshape(temp, size(X_agg))

我没有尝试过它的效率有多高，但它应该对大型矩阵做得更好，因为它不使用for循环。

于 2013-05-22T04:03:18.643 回答