5

こんにちは、[x,y] ポイントのベクトルを使用して、MATLAB の大きな行列からインデックスを作成する方法を見つけようとしています。通常、添字の点を行列の線形インデックスに変換します(たとえば、ベクトルを行列のインデックスとして使用します)ただし、行列は4次元であり、のすべての要素を取得したい同じ 1 次元と 2 次元を持つ 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) などでインデックスを作成し、 に存在する 999 および 888 だけでなく、 、および などMatrix(:,:,1,1)に存在するコンテンツを返すことができるようにしたいと考えています(IRL 、の次元はもっと似ているかもしれませんMatrix(:,:,1,2)Matrix(:,:,2,1)Matrix(:,:,2,2)Matrixsize(Matrix) = (300 250 30 200)

結果を同様のベクトル形式にしたいので、線形インデックスを使用したくありません。たとえば、次のような結果が必要です。

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

また、私が扱っている行列のサイズのために、ここでは速度が問題になることを追加したいと思います.

また、(この質問とは異なり: sub2ind を使用せずに添え字を使用して値にアクセスする) i および j 番目のインデックスの余分な次元 3 および 4 にすべての情報を格納する必要があるため、それが少し速いとは思わないことにも言及する必要があります。のバージョンはsub2indまだそれをカットしません..

4

2 に答える 2

8

これについて3つの方法を考えることができます

単純なループ

持っているすべての 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 が最後の 2 つのループをコンパイルできなかったことと、4D 配列に対する非特異性 (最後の 2 つの方法は ND 配列でも機能します) の影響を大きく受けています。4D に特化したバージョンを作成すると、間違いなくはるかに高速になります。

それにもかかわらず、単純なループによるインデックス作成は、JIT のおかげで、最も簡単に行うことができ、目に優しく、非常に高速でもあります。

于 2012-11-26T06:07:44.230 に答える
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

結果の行列は、resultsize を持つ必要がありますsize(result) = (4 N tt)

Matrix正方形でなくても、これはうまくいくはずです。ただし、上で述べたように、それは理想的ではありません。

于 2012-11-26T05:41:24.090 に答える