使用载体来索引的矩阵而不线性指数(Use a vector to index a matrix wi

2019-07-04 03:21发布

天儿真好,我试图找到一种方法,从MATLAB中的大型矩阵使用[X,Y]点的矢量指数。 通常,我将下标点转换为基体的线性索引。(对于例如。 使用一个向量作为索引的矩阵 )。然而,该矩阵是4维,和我想利用所有的元素的第三,并且具有相同的第1和第2维第四维。 让我希望用一个例子证明:

Matrix = nan(4,4,2,2); % where the dimensions are (x,y,depth,time)
Matrix(1,2,:,:) = 999; % note that this value could change in depth (3rd dim) and time (4th time) 
Matrix(3,4,:,:) = 888; % note that this value could change in depth (3rd dim) and time (4th time) 
Matrix(4,4,:,:) = 124;

现在,我希望能够与标的指数(1,2)和(3,4)等,并返回不仅是999和888存在于Matrix(:,:,1,1)但其内容存在于Matrix(:,:,1,2) Matrix(:,:,2,1)Matrix(:,:,2,2)等等(IRL,的尺寸Matrix可能更喜欢size(Matrix) = (300 250 30 200)

我不希望使用线性指标,因为我想的结果是类似的矢量时尚。 例如,我想结果是这样的:

ans(time=1)
999 888 124
999 888 124
ans(time=2)
etc etc etc
etc etc etc

我还想补充一点,由于我负责的,速度是这里的一个问题矩阵的大小 - 因此,为什么我想使用下标索引来索引到的数据。

我还要提到的(不像这个问题: 访问使用标值,而无需使用sub2ind ),因为我想存储在额外维度,3和4,i和第j个指标的全部信息,我不认为一个稍快版本sub2ind仍然不会削减它..

Answer 1:

我能想到的三种方式去了解这个

简单的循环

就在所有的你有,并用冒号访问其余尺寸的2D指数环:

for jj = 1:size(twoDinds,1)
    M(twoDinds(jj,1),twoDinds(jj,2),:,:) = rand;
end

线性指数的矢量计算

跳过sub2ind和矢量化线性指数的计算:

% generalized for arbitrary dimensions of M

sz = size(M);
nd = ndims(M);

arg = arrayfun(@(x)1:x, sz(3:nd), 'UniformOutput', false);

[argout{1:nd-2}] = ndgrid(arg{:});

argout = cellfun(...
    @(x) repmat(x(:), size(twoDinds,1),1), ...
    argout, 'Uniformoutput', false);

twoDinds = kron(twoDinds, ones(prod(sz(3:nd)),1));

% the linear indices
inds = twoDinds(:,1) + ([twoDinds(:,2) [argout{:}]]-1) * cumprod(sz(1:3)).';

Sub2ind

只需使用现成的工具,船舶用Matlab:

inds = sub2ind(size(M), twoDinds(:,1), twoDinds(:,2), argout{:});

速度

那么,哪一个是最快的? 让我们来看看:

clc

M = nan(4,4,2,2);

sz = size(M);
nd = ndims(M);

twoDinds = [...
    1 2
    4 3
    3 4
    4 4
    2 1];

tic
for ii = 1:1e3
    for jj = 1:size(twoDinds,1)
        M(twoDinds(jj,1),twoDinds(jj,2),:,:) = rand;
    end
end
toc


tic
twoDinds_prev = twoDinds;
for ii = 1:1e3

    twoDinds = twoDinds_prev;

    arg = arrayfun(@(x)1:x, sz(3:nd), 'UniformOutput', false);

    [argout{1:nd-2}] = ndgrid(arg{:});

    argout = cellfun(...
        @(x) repmat(x(:), size(twoDinds,1),1), ...
        argout, 'Uniformoutput', false);

    twoDinds = kron(twoDinds, ones(prod(sz(3:nd)),1));
    inds = twoDinds(:,1) + ([twoDinds(:,2) [argout{:}]]-1) * cumprod(sz(1:3)).';

    M(inds) = rand;

end
toc


tic
for ii = 1:1e3

  twoDinds = twoDinds_prev;

    arg = arrayfun(@(x)1:x, sz(3:nd), 'UniformOutput', false);

    [argout{1:nd-2}] = ndgrid(arg{:});

    argout = cellfun(...
        @(x) repmat(x(:), size(twoDinds,1),1), ...
        argout, 'Uniformoutput', false);

    twoDinds = kron(twoDinds, ones(prod(sz(3:nd)),1));

    inds = sub2ind(size(M), twoDinds(:,1), twoDinds(:,2), argout{:});

    M(inds) = rand;
end
toc

结果:

Elapsed time is 0.004778 seconds.  % loop
Elapsed time is 0.807236 seconds.  % vectorized linear inds
Elapsed time is 0.839970 seconds.  % linear inds with sub2ind

结论:使用循环。

诚然,上述测试是由JIT未能编译最后两个循环在很大程度上影响和非特异性4D阵列(最后两个方法也ND阵列工作)。 制作四维专用版本无疑会快很多。

然而,与简单的循环索引,没错,是最简单的事,最容易对眼睛非常快过,由于JIT。



Answer 2:

所以,这里是一个可能的答案...但它是凌乱。 我怀疑它会更昂贵的计算则更为直接的方法......这绝对不会是我的首选答案。 这将是巨大的,如果我们能找到答案,而无需任何的循环!

Matrix = rand(100,200,30,400);
grabthese_x = (1 30 50 90);
grabthese_y = (61 9 180 189);
result=nan(size(length(grabthese_x),size(Matrix,3),size(Matrix,4));
for tt = 1:size(Matrix,4)
subset = squeeze(Matrix(grabthese_x,grabthese_y,:,tt));
 for NN=1:size(Matrix,3)
 result(:,NN,tt) = diag(subset(:,:,NN));
 end
end

将得到的矩阵, result应具有尺寸size(result) = (4 N tt)

我认为这应该工作,即使Matrix不是方形。 然而,这不是理想的,正如我上面所说。



文章来源: Use a vector to index a matrix without linear index