MATLAB 更高效地减去矩阵答案

【问题标题】：MATLAB Subtract Matrices More EfficientlyMATLAB 更高效地减去矩阵
【发布时间】：2013-06-25 14:25:37
【问题描述】：

我将两个矩阵相减。 dataClim 是 30 年期间每个月（12 个月）的平均数据。 dataAll 是超过 1257 天的每日数据。我需要从20100101到20130611的每个月的日数据中减去平均月数据（t = 1:31是一月，t = 32-57是二月，一直到十二月，然后363:393又是一月） .

此代码有效，但我想知道是否有任何方法可以使其更高效且不那么乏味。我不知道如何编写循环，因为月份的天数从 28 到 31 不等。

% Create new array in which data_Anom is the anomaly 
% dataAnom = dataAll - dataClim 

% January
dataAnom_1 = bsxfun(@minus, dataAll(:,:,[1:31, 363:393, 728:758, 1094:1124]), dataClim(:,:,1));

% February
dataAnom_2 = bsxfun(@minus, dataAll(:,:,[32:57, 394:421, 759:787, 1125:1152]), dataClim(:,:,2));

% March
dataAnom_3 = bsxfun(@minus, dataAll(:,:,[58:88, 422:452, 788:818, 1153:1183]), dataClim(:,:,3));

% April
dataAnom_4 = bsxfun(@minus, dataAll(:,:,[89:118, 453:482, 819:848, 1184:1213]), dataClim(:,:,4));

% May
dataAnom_5 = bsxfun(@minus, dataAll(:,:,[119:148, 483:513, 849:879, 1214:1244]), dataClim(:,:,5));

% June 
dataAnom_6 = bsxfun(@minus, dataAll(:,:,[149:178, 514:543, 880:909, 1245:1255]), dataClim(:,:,6));

% July
dataAnom_7 = bsxfun(@minus, dataAll(:,:,[179:209, 544:574, 910:940]), dataClim(:,:,7));

% August
dataAnom_8 = bsxfun(@minus, dataAll(:,:,[210:240, 575:605, 941:971]), dataClim(:,:,8));

% September
dataAnom_9 = bsxfun(@minus, dataAll(:,:,[241:270, 606:635, 972:1001]), dataClim(:,:,9));

% October
dataAnom_10 = bsxfun(@minus, dataAll(:,:,[271:301, 636:666, 1002:1032]), dataClim(:,:,10));

% November
dataAnom_11 = bsxfun(@minus, dataAll(:,:,[302:331, 667:696, 1033:1062]), dataClim(:,:,11));

% December
dataAnom_12 = bsxfun(@minus, dataAll(:,:,[332:362, 697:727, 1063:1093]), dataClim(:,:,12));

% Concatenate the seperate Anomalies
dataAnom = cat(3, dataAnom_1, dataAnom_2, dataAnom_3, dataAnom_4, dataAnom_5, dataAnom_6, dataAnom_7, dataAnom_8, dataAnom_9, dataAnom_10, dataAnom_11, dataAnom_12);

clear dataAnom_*

我的一个想法是先将每个月的天数连接在一起，然后为每个月创建 dataAnom。它可能更慢。

% Concatenation days below to each month in dataAll into dataMon so that each month is placed together. This
% makes it easier to do the anomaly subtraction later.

dataMon = cat(3, dataAll(:,:,1:31), dataAll(:,:,363:393), dataAll(:,:,728:758) , dataAll(:,:,1094:1124),... % January
    dataAll(:,:,32:57), dataAll(:,:,394:421), dataAll(:,:,759:787), dataAll(:,:,1125:1152),... % February
    dataAll(:,:,58:88), dataAll(:,:,422:452), dataAll(:,:,788:818), dataAll(:,:,1153:1183),... % March
    dataAll(:,:,89:118), dataAll(:,:,453:482), dataAll(:,:,819:848), dataAll(:,:,1184:1213),... % April
    dataAll(:,:,119:148), dataAll(:,:,483:513), dataAll(:,:,849:879), dataAll(:,:,1214:1244),... % May
    dataAll(:,:,149:178), dataAll(:,:,514:543), dataAll(:,:,880:909), dataAll(:,:,1245:1255),... % June. Last entry goes up to 20130611, not the full month
    dataAll(:,:,179:209), dataAll(:,:,544:574), dataAll(:,:,910:940),... % July
    dataAll(:,:,210:240), dataAll(:,:,575:605), dataAll(:,:,941:971),... % August
    dataAll(:,:,241:270), dataAll(:,:,606:635), dataAll(:,:,972:1001),... % Sept
    dataAll(:,:,271:301), dataAll(:,:,636:666), dataAll(:,:,1002:1032),... % Oct
    dataAll(:,:,302:331), dataAll(:,:,667:696), dataAll(:,:,1033:1062),... % Nov
    dataAll(:,:,332:362), dataAll(:,:,697:727), dataAll(:,:,1063:1093)); % Dec

% Create dataAnom
dataAnom1 = bsxfun(@minus, dataAll(:,:,1:124), dataClim(:,:,1);
dataAnom2 = bsxfun(@minus, dataAll(:,:,125:238:), dataClim(:,:,1);
.
.
.
dataAnom12 = ...

% Combine
dataAnom = cat(3, dataAnom1, dataAnom2, dataAnom3,....);

【问题讨论】：

顺便说一句，为什么不使用更好的矩阵结构呢？即12x31xNumOfYears?对于不是 31 天的月份，您可以使用 NaN 来表示缺失的日期。
为什么 2 月是 32-57？如果您的数据结构不正确，则必须对所有内容进行硬编码，这会很糟糕。
二月是 28 天或 29 天。数据只是简单地列出了一天又一天。矩阵的 x 和 y 是全球的 lat 和 lon，z 是从 20100101 到 20130611 的天数。我认为只添加 NaN 以使每个月有 31 天或其他东西不是一个好主意。
@user2258883 查看我的回答，以便更好地理解我的意思。
您的二月似乎既没有 28 天也没有 29 天，而是 26 天……另外请注意，第一年有 362 天。这怎么可能？不同年份的所有对应月份是否至少具有相同的天数？

标签： matlab matrix subtraction bsxfun

【解决方案1】：

我认为您需要重新考虑如何处理此处构建数据的方式。我不会创建大量的单行数组（dataAnom1 和dataAll），而是使用更好的结构化矩阵。

最简单的情况是，你可以使用这样的方案：31x12xNumYears，它会产生这样的结果：

Data = NaN(31x12xNumYears); % Blank Init
Data(1:31,1,1) = Rand(31,1); % January's populated values
Data(1:28,2,1) = Rand(28,1); % February's populated values
... % and so forth

这里的优点是矩阵运算更容易执行，并且您可以更好地了解数据实际代表什么。例如，假设 dataClim 是 30 年的月平均值（矩阵应为 12x30），dataAll 是每日读数（矩阵应为 31x12x30），您可以执行以下操作：

subValues = NaN(31,12,30);
for yr = 1:30
    for mn = 1:12
        subValues(1:31,mn,yr) = dataAll(1:31,mn,yr) - dataClim(mn,yr);
    end
end

根据您提供给我的附加信息，我认为这是您可能正在寻找的结构类型：days x months x years x 3 其中三个代表数据的纬度、经度和数据值。比如：

test = rand(31,12,30,3);
lat = test(1:end,1:end,1:end,1); 
long = test(1:end,1:end,1:end,2);
data = test(1:end,1:end,1:end,3);

【讨论】：

问题是我的 dataAll 包含全球每个纬度和经度点的数据，因此我无法更改它而不会丢失每个单元格对应的地球上的哪个点。也许我只是不明白你在说什么，但我不确定如何对其进行重组以使其更容易做到。
@user2258883 只是为了澄清一下，它不是每天的单个值，而是纬度和经度值？
纬度和经度值在不同的矩阵中。 X 是一个 1437x159 矩阵，1:360（长）沿第一列向下并在所有列中重复。 Y 的大小相同，-89.5 到 -50（南纬）穿过第一行并重复所有行。 dataAll (1437x159x1258) 是对应于 1258 天内这些点的数据（我上面的索引是错误的，因为您指出我在 2 月份错过了两天，而索引是基于没有这两天的）。所以 dataAll(1,1,1) 将是 1 (lon), -89.5 (lat), day 1 (20100101) 的海冰浓度。
@user2258883 好的，我对您的数据的外观有了更好的了解，但我仍然建议回到这些数据的来源，并对其进行重构以更好地满足您的需求。处理这个问题会很痛苦。
我认为顶部可能有错字。二月份的数据应该存储在Data(1:29,2,1)，对吧？

【解决方案2】：

虽然我仍然坚持您应该重构数据这一事实，但您也可以使用 datenum 的力量来帮助您。与尝试手动输入日期相比，以下内容将使您的生活更轻松：

clear all; clc;

% Init Data
dataAll = rand(1437,159,1258);

% Starting date of sampling.  Note that this assumes each day there was a sample, and only one sample
startDate = datenum('01-01-2010');
dateList = [0:1257] + startDate;
[yr, mn, ~, ~, ~, ~] = datevec(dateList);

% extract data depending on month that sample was taken
jan = dataAll(:,:,mn == 1);
feb = dataAll(:,:,mn == 2);
mar = dataAll(:,:,mn == 3);
... % and so forth

这应该会得到您想要的结果。从这里您可以进行结果计算：

dataAnom_1 = bsxfun(@minus, jan, dataClim(:,:,1));

使用来自 cmets 的更新信息，您可以执行以下操作以根据月份和年份分离数据：

jan2010 = dataAll(:,:,(mn == 1 & yr == 2010);
feb2010 = dataAll(:,:,(mn == 2 & yr == 2010);
... % and so forth

【讨论】：

谢谢。这行得通。但我发现我不能完全像这样这样做，因为这将所有一月的时间放在一起，然后是整个二月（我之前没有意识到这一点是我的错）。这破坏了日子的秩序。
@user2258883 “破坏了时代的秩序”是什么意思？
jan = dataAll(:,:,mn == 1); 将所有一月份的数据放在一起。所以当我稍后做dataAnom_1 = bsxfun(@minus, jan, dataClim(:,:,1)); 时，我最终会从dataClim 中减去所有一月份的数据。当我将 dataAnom_1、dataAnom_2 等组合回一个大矩阵时，它最终是错误的，因为我最终将 4 个不同的 1 月、2 月 4 日等的异常绘制在一起，而不是绘制 2010 年 1 月、2010 年 2 月、.. .. 2011 年 1 月，2011 年 2 月，...
@user2258883 我已经更新了答案，但这只是证实了我的原始答案，即您确实需要找出更好的结构来存储和操作数据。

【解决方案3】：

你不能这样做

% January
dataAnom1 = bsxfun(@minus, dataAll(:,:,[1:31 363:393 728:758 1094:1124]), dataClim(:,:,1));

？我认为这是一样的。

如果是这样，那么你可以这样做

dataAnom1=zeros(size(dataAll,1),size(dataAll,2), 128*12);
for v=1:12
    dataAnom1(:,:,1+((v-1)*128:v*128)) = bsxfun(@minus, dataAll(:,:,[1:31 363:393 728:758 1094:1124]+(v-1*32)), dataClim(:,:,v));
end

（索引可能有点偏离）

【讨论】：

哦，好像我之前尝试过这个时忘记了括号。谢谢。
问题是每个月有 31 或 30 天，而 2 月有 28 或 29 天。我不知道如何在代码中考虑到这一点。
啊，我明白了。是的，在这样的几个月里工作可能会很尴尬。您可以将每个月的信息放入一个 12 元素的单元格中，并对单元格中的数据使用向量表达式吗？例如。 MonthData = cell(12,1); MonthData{1}=dataAll(:,:,1:31);，然后使用 bsxfun(@minus,MonthData{1},dataClim(:,:,1)); 如果合适，您可以将它们包装在结构或单元格中多年：Years = struct(4); Years(1).year=2009; Years(1).Data = MonthData;
@user2258883 提醒一下，使用单元格和执行操作可能有点棘手。