如何在较大矩阵的子集上矢量化运行函数的代码？

Question

我们假设我有以下9 x 5矩阵:

myArray = [
   54.7    8.1   81.7   55.0   22.5
   29.6   92.9   79.4   62.2   17.0
   74.4   77.5   64.4   58.7   22.7
   18.8   48.6   37.8   20.7   43.5
   68.6   43.5   81.1   30.1   31.1
   18.3   44.6   53.2   47.0   92.3
   36.8   30.6   35.0   23.0   43.0
   62.5   50.8   93.9   84.4   18.4
   78.0   51.0   87.5   19.4   90.4
];

我有11个这个矩阵的"子集",我需要max在每个子集上运行一个函数(比方说).可以使用以下逻辑逻辑来识别子集(按列标识,而不是按行标识):

myLogicals = logical([
    0 1 0 1 1
    1 1 0 1 1
    1 1 0 0 0
    0 1 0 1 1
    1 0 1 1 1
    1 1 1 1 0
    0 1 1 0 1
    1 1 0 0 1
    1 1 0 0 1
]);

或通过线性索引:

starts = [2 5 8 10 15 23 28 31 37 40 43]; #%index start of each subset
ends =   [3 6 9 13 18 25 29 33 38 41 45]; #%index end of each subset

这样第一个子集是2:3,第二个子集是5:6,依此类推.

我可以找到max每个子集的内容并将其存储在向量中,如下所示:

finalAnswers = NaN(11,1); 
for n=1:length(starts) #%i.e. 1 through the number of subsets
    finalAnswers(n) = max(myArray(starts(n):ends(n)));
end

循环运行后,finalAnswers包含每个数据子集的最大值:

74.4  68.6  78.0  92.9  51.0  81.1  62.2  47.0  22.5  43.5  90.4

是否可以在不使用for循环的情况下获得相同的结果？换句话说,这个代码可以被矢量化吗？这种方法会比现在更有效吗？

---

编辑:我做了一些测试建议的解决方案.我使用的数据是1,510 x 2,185矩阵,其中10,103个子集的长度从2到916不等,子集长度的标准偏差为101.92.

我将每个解决方案tic;for k=1:1000 [code here] end; toc;都包含在内,结果如下:

• for 循环方式--- Elapsed time is 16.237400 seconds.
• Shai的方法--- Elapsed time is 153.707076 seconds.
• 丹的方法--- Elapsed time is 44.774121 seconds.
• 迪瓦卡的方法#2 --- Elapsed time is 127.621515 seconds.

笔记:

• 我还尝试通过将k=1:1000 for循环包裹在线周围来对Dan的方法进行基准测试accumarray(因为其余的理论上只运行一次).在这种情况下,时间是28.2​​9秒.
• 对Shai的方法进行基准测试,同时将lb = ...线路排除在外k,时间为113.48秒.
• 当我运行Divakar的代码时,我遇到了Non-singleton dimensions of the two input arrays must match each other.错误bsxfun.我通过使用共轭转置(撇号操作者"固定的"本'上)trade_starts(1:starts_extent)和 intv(1:starts_extent)在代码调用的线bsxfun.我不确定为什么会出现这个错误......

我不确定我的基准设置是否正确,但for在这种情况下,循环实际上运行速度最快.

Answer 1

bwlabeln与垂直连接一起使用：

lb = bwlabeln( myLogicals, [0 1 0; 0 1 0; 0 1 0] );

现在每个区域都有一个标签 1..11。

要获得最大值，您可以使用regionprops

props = regionprops( lb, myArray, 'MaxIntensity' );
finalAnswers = [props.MaxIntensity];

您可以使用它regionprops来获取每个子集的一些其他属性，但它不是太通用。
如果您希望对每个区域应用更通用的函数，例如 ，median您可以使用accumarray：

finalAnswer = accumarray( lb( myLogicals ), myArray( myLogicals ), [], @median );