首页 文章

使用向量索引没有线性索引的矩阵

提问于
浏览 次
5

G 'day, I' m试图找到一种方法来使用[x，y]点的矢量从MATLAB中的大矩阵中索引 . 通常，我会将下标点转换为矩阵的线性索引 . （例如.Use a vector as an index to a matrix）但是，矩阵是4维的，我想要采用具有相同的第3和第4维的所有元素第一和第二维度 . 让我希望通过一个例子来证明：

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）等进行索引，并且不仅返回 Matrix(:,:,1,1) 中存在的999和888，而且返回 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

我还想补充一点，由于我正在处理的矩阵的大小，速度是一个问题 - 因此我想使用下标索引来索引数据 .

我还应该提一下（不像这个问题：Accessing values using subscripts without using sub2ind），因为我想要存储在i和jth索引的额外维度3和4中的所有信息，我不认为 sub2ind 稍微快一点的版本仍然不会剪了它..

matlab vector matrix indexing subscript

2 回答

  • 8

    我可以想到三种方法

    简单循环

    只需遍历您拥有的所有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阵列） . 为4D制作专用版本无疑会更快 .

    然而，由于JIT，使用简单循环的索引最简单，最简单，也很快 .

    回复于
  • 1

    所以，这是一个可能的答案......但它很混乱 . 我怀疑它的计算成本会比更直接的方法更昂贵......而这绝对不是我的首选答案 . 如果我们能够在没有任何for循环的情况下得到答案那就太棒了！

    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 不是正方形 . 但是，正如我上面所说，它并不理想 .

    回复于

相关问题