如何在不创建子矩阵的情况下对 Matlab 中的部分矩阵求和？

Question

所以在 Matlab 中，假设我有一个大小为 N 乘以 N 的矩阵 X，并且 i 是一个大小为 1 乘以 N 的逻辑索引向量。那么我可以做

sum(X(i,i))

问题是相当于先为

分配内存

Y=X(i,i),

然后计算 Y 的和，然后删除 Y。我说得对吗？（Hoki的回答证明是对的。）

有没有更快的方法来计算总和而不（隐含地）创建 Y？如果 Y 很大，内存操作会消耗大量时间。换句话说，是否可以执行以下操作：

S=zeros(1,nnz(i));

for k=find(i)
    for j=find(i)
        S(k)=S(k)+X(j,k);
    end
end

这样，除了 X 之外，我们需要的所有内存都是向量 S - 我们不需要为大 Y 分配内存。当然，循环可能很慢，但你明白我的意思。

Answer 1

有两个答案，如果您一直在寻找完整的专栏，答案很简单

t=sum(X);

是一行所有列的总和

然后

ans=sum(t(i))

就是你想要的。

如果您正在寻找奇怪的形状，线性索引可能就是您要找的。

见 sub2ind

首先在矩阵中创建一个线性索引（一维索引）然后直接使用该索引

第 i 列中六个项目（5 到 10）的使用总和

ind = sub2ind(size(X) , ones(6,1)* i , (5:10)'*ones(1,N))  ; 
sum(X(ind))

Answer 2

您对内存管理的运作方式假设过多。

时间：

我运行 timeit 的基准测试。从 N=10 到 N=20000，两种形式的执行时间绝对没有明显差异。

此外，如果我关闭 JIT 加速，我仍然会得到完全相同的结果...所以优化可能只是 Matlab lazy-copy 行为的结果。

内存使用：

在内存方面，似乎有所不同。间接方法（使用临时变量）似乎为此临时变量分配内存（分配的大小与临时变量的大小完全对应）。另一方面，直接方法不需要为 return 结果分配任何额外的内存。

这已经达到了我对这些东西的把握的极限。我不够专业，无法假装解释为什么这种内存使用差异不会导致时序差异。我知道内存很快，但对于 N 的高阶，我认为它会有所作为。显然不是...

更多信息：

关于Matlab内存管理的更多细节，我邀请你阅读这篇来自Loren at Matlab的文章：
Memory Management for Functions and Variables

或者如果您想阅读更深入的机制测试：
Internal Matlab memory optimizations

时间基准：

基准测试结果：

基准代码：

function ExecTimes = benchmark_sumcol

%// prepare logarithmic progression (up to what my 16GB RAM can take)
nOrder = (1:9).' * 10.^(1:3) ; nOrder = [nOrder(:) ; 10000 ; 20000] ; %'
npt = numel(nOrder) ;

ExecTimes = zeros( npt , 2 ) ;

for k = 1:npt
    %// Sample data
    N = nOrder(k) ;
    X = rand(N) ;
    ci = logical(randi([0 1],1,N)) ;

    %// Benchmark
    f1 = @() direct_sum(X,ci) ;
    f2 = @() indirect_sum(X,ci) ;
    ExecTimes(k,1) = timeit( f1 ) ;
    ExecTimes(k,2) = timeit( f2 ) ;

    clear X ci
    disp(N)
end

function R = direct_sum(X,ci)
R = sum(X(:,ci)) ;

function R = indirect_sum(X,ci)
Y = X(:,ci) ;
R = sum(Y) ;

内存基准：

两个函数的摘要

间接求和的详细信息，带有临时变量。我突出显示了内存分配：

直接求和的详细信息：

内存基准代码

%% // set profiler options
clear all
profile('-memory','on');
setpref('profiler','showJitLines',1);
profile on

%% // sample data
N = 1000 ;
X = rand(N) ;
ci = logical(randi([0 1],1,N)) ;

%% // Benchmark
R2 = bench_indirect_sum(X,ci) ;
R1 = bench_direct_sum(X,ci) ;

%% // result
profile viewer
p = profile('info');
profsave(p,'profile_results')

上次编辑：

我将你的 loop 版本添加到测试中，尽管我不得不对其进行一些修改以使其实际工作（并给出与其他版本相同的结果）：

function R = bench_loop_sum(X,ci)
    R = zeros(1,nnz(ci));
    idxRes=1 ;
    for k=find(ci)
         for j=1:size(X,1)
             R(idxRes)=R(idxRes)+X(j,k);
         end
        idxRes = idxRes+1 ;
    end

结果在内存方面还可以（即没有为临时数组分配额外的内存），但在速度方面是灾难性的：

正如我们对循环的预期，关闭 JIT 会更糟：

现在一个简单的改变来抑制内部循环使事情变得更好，但仍然有点落后于直接方式（注意这个版本不为临时分配内存做总和的列）：

function R = bench_loop_sum(X,ci)
    R = zeros(1,nnz(ci));
    idxRes=1 ;
    for k=find(ci)
        R(idxRes) = sum(X(:,k));
        idxRes = idxRes+1 ;
    end

启用 JIT。

如何在不创建子矩阵的情况下对 Matlab 中的部分矩阵求和？

How to sum part of a matrix in Matlab without creating a submatrix?

indexing

matlab

sum

matrix

时间基准：

内存基准：

上次编辑：