复制向量将它们向右移动

Question

在 Matlab 中，我在 2x249 矩阵中有两个单行 (1x249) 向量，我必须通过多次复制它们来创建矩阵 A，每次将 2 个位置的向量向右移动。我想用零填充左侧的条目。有没有聪明的方法来做到这一点？目前，我正在使用 for 循环和 circshift，并且在每次迭代中添加新行到 A，但这可能效率非常低。

代码（myMat是我要移位的矩阵）：

A = [];
myMat = [1 0 -1 zeros(1,246); 0 2 0 -2 zeros(1,245)];
N = 20;
for i=1:N-1
    aux = circshift(myMat,[0,2*(i-1)]);
    aux(:,1:2*(i-1)) = 0;
    A =[A; aux];
end

Answer 1

您可能知道，Matlab 中的循环效率不高。 我知道 Mathworks 一直说 JIT 不再如此 编译，但我还没有体验过快速循环。

我把你构建矩阵 A 的方法放在一个函数中：

function A = replvector1(myMat,shift_right,width,N)

    pre_alloc = true; % make implementation faster using pre-allocation yes/no

    % Pad myMat with zeros to make it wide enough
    myMat(1,width)=0;

    % initialize A
    if pre_alloc
       A = zeros(size(myMat,1)*(N-1),width);
    else
       A = [];
    end

    % Fill A
    for i=1:N-1
        aux = circshift(myMat,[0,shift_right*(i-1)]);
        aux(:,1:min(width,shift_right*(i-1))) = 0;

        A(size(myMat,1)*(i-1)+1:size(myMat,1)*i,:) =aux;
    end

您的矩阵运算看起来很像克罗内克积，但是 块矩阵具有重叠的列范围，因此是直接克罗内克积 不管用。相反，我构造了以下函数：

function A = replvector2(myMat,shift_right,width,N)
    [i,j,a] = find(myMat);
    i = kron(ones(N-1,1),i) + kron([0:N-2]',ones(size(i))) * size(myMat,1);
    j = kron(ones(N-1,1),j) + kron([0:N-2]',ones(size(j))) * shift_right;
    a = kron(ones(N-1,1),a);
    ok = j<=width;
    A = full(sparse(i(ok),j(ok),a(ok),(N-1)*size(myMat,1),width));

您可以通过删除分号并查看中间值来遵循该算法 结果。

以下主程序运行您的示例，并且可以轻松修改为 运行类似的例子：

% inputs (you may vary them to see that it always works)
shift_right = 2;
width       = 249;
myMat1      =  [ 1 0 -1  0 ;
                 0 2  0 -2 ];
N           = 20;

% Run your implementation
tic;
A = replvector1(myMat,shift_right,width,N);
disp(sprintf('\n   original implementation took %e sec',toc))

% Run the new implementation
tic;
B = replvector2(myMat,shift_right,width,N);
disp(sprintf('   new implementation took %e sec',toc))

disp(sprintf('\n   norm(B-A)=%e\n',norm(B-A)))

Answer 2

我采用了 Nathan 的代码 (see his answer to this question)，并添加了另一个可能的实现 (replvector3)。

我的想法源于您并不真正需要循环移位。您需要右移并在左侧添加零。如果您从一个预先分配的数组开始（这确实是您在时间上取得重大胜利的地方，剩下的就是花生），那么您已经有了零。现在您只需将myMat 复制到正确的位置即可。

这些是我看到的时间（MATLAB R2017a）：

OP's, with pre-allocation: 1.1730e-04
Nathan's:                   5.1992e-05
Mine:                       3.5426e-05
                            ^ shift by one on purpose, to make comparison of times easier

这是完整的副本，复制粘贴到 M 文件中并运行：

function so

shift_right = 2;
width       = 249;
myMat       =  [ 1 0 -1  0 ;
                 0 2  0 -2 ];
N           = 20;

A = replvector1(myMat,shift_right,width,N);
B = replvector2(myMat,shift_right,width,N);
norm(B(:)-A(:))
C = replvector3(myMat,shift_right,width,N);
norm(C(:)-A(:))

timeit(@()replvector1(myMat,shift_right,width,N))
timeit(@()replvector2(myMat,shift_right,width,N))
timeit(@()replvector3(myMat,shift_right,width,N))

% Original version, modified to pre-allocate
function A = replvector1(myMat,shift_right,width,N)
    % Assuming width > shift_right * (N-1) + size(myMat,2)
    myMat(1,width) = 0;
    M = size(myMat,1);
    A = zeros(M*(N-1),width);
    for i = 1:N-1
        aux = circshift(myMat,[0,shift_right*(i-1)]);
        aux(:,1:shift_right*(i-1)) = 0;
        A(M*(i-1)+(1:M),:) = aux;
    end

% Nathan's version
function A = replvector2(myMat,shift_right,width,N)
    [i,j,a] = find(myMat);
    i = kron(ones(N-1,1),i) + kron((0:N-2)',ones(size(i))) * size(myMat,1);
    j = kron(ones(N-1,1),j) + kron((0:N-2)',ones(size(j))) * shift_right;
    a = kron(ones(N-1,1),a);
    ok = j<=width;
    A = full(sparse(i(ok),j(ok),a(ok),(N-1)*size(myMat,1),width));

% My trivial version with loops
function A = replvector3(myMat,shift_right,width,N)
    % Assuming width > shift_right * (N-1) + size(myMat,2)
    [M,K] = size(myMat);
    A = zeros(M*(N-1),width);
    for i = 1:N-1
        A(M*(i-1)+(1:M),shift_right*(i-1)+(1:K)) = myMat;
    end