matlab - 基于长度向量的 MATLAB 重复数字

Question

是否有一种矢量化的方式来执行以下操作？（通过示例显示）：

input_lengths = [ 1 1 1 4       3     2   1 ]
result =        [ 1 2 3 4 4 4 4 5 5 5 6 6 7 ]

我已将 input_lengths 隔开，因此很容易理解如何获得结果

结果向量的长度为：sum(lengths)。我目前result使用以下循环进行计算：

result = ones(1, sum(input_lengths ));
counter = 1;
for i = 1:length(input_lengths)
    start_index = counter;
    end_index = counter + input_lengths (i) - 1;

    result(start_index:end_index) = i;
    counter = end_index + 1;
end

编辑：

我也可以使用 arrayfun 来做到这一点（虽然这不完全是一个矢量化函数）

cell_result = arrayfun(@(x) repmat(x, 1, input_lengths(x)), 1:length(input_lengths), 'UniformOutput', false);
cell_result : {[1], [2], [3], [4 4 4 4], [5 5 5], [6 6], [7]}

result = [cell_result{:}];
result : [ 1 2 3 4 4 4 4 5 5 5 6 6 7 ]

Answer 1

完全矢量化的版本：

selector=bsxfun(@le,[1:max(input_lengths)]',input_lengths);
V=repmat([1:size(selector,2)],size(selector,1),1);
result=V(selector);

缺点是，内存使用量为 O(numel(input_lengths)*max(input_lengths))

Answer 2

所有解决方案的基准

在之前的基准测试之后，我将这里给出的所有解决方案组合在一个脚本中，并运行几个小时进行基准测试。我这样做是因为我认为以输入长度为参数查看每个提议的解决方案的性能是很好的 - 我的目的不是要降低前一个解决方案的质量，它提供了有关效果的更多信息即时通讯。此外，每个参与者似乎都同意这一点，所有答案都做得很好，所以这篇很棒的帖子值得一个总结帖。

我不会在这里发布脚本的代码，这很长而且很无趣。基准测试的过程是针对一组不同长度的输入向量运行每个解决方案：10、20、50、100、200、500、1000、2000、5000、10000、20000、50000、100000、200000、500000、 1000000。对于每个输入长度，我根据泊松定律生成了一个随机输入向量，参数为 0.8（以避免大值）：

input_lengths = round(-log(1-rand(1,ILen(i)))/poisson_alpha)+1;

最后，我平均每个输入长度超过 100 次运行的计算时间。

我已经用 Matlab R2013b 在我的笔记本电脑（核心 I7）上运行了脚本；JIT 被激活。

这里是绘制的结果（对不起，彩色线），以对数刻度（x 轴：输入长度；y 轴：计算时间，以秒为单位）：

所以 Luis Mendo 是明显的赢家，恭喜！

对于任何想要数值结果和/或想要重新绘制它们的人，它们在这里（将表格分成 2 部分并近似为 3 位，以便更好地显示）：

N                   10          20          50          100         200         500         1e+03       2e+03
-------------------------------------------------------------------------------------------------------------
OP's for-loop       8.02e-05    0.000133    0.00029     0.00036     0.000581    0.00137     0.00248     0.00542 
OP's arrayfun       0.00072     0.00117     0.00255     0.00326     0.00514     0.0124      0.0222      0.047
Daniel              0.000132    0.000132    0.000148    0.000118    0.000126    0.000325    0.000397    0.000651
Divakar             0.00012     0.000114    0.000132    0.000106    0.000115    0.000292    0.000367    0.000641
David's for-loop    9.15e-05    0.000149    0.000322    0.00041     0.000654    0.00157     0.00275     0.00622
David's arrayfun    0.00052     0.000761    0.00152     0.00188     0.0029      0.00689     0.0122      0.0272
Luis Mendo          4.15e-05    4.37e-05    4.66e-05    3.49e-05    3.36e-05    4.37e-05    5.87e-05    0.000108
Bentoy13's cumsum   0.000104    0.000107    0.000111    7.9e-05     7.19e-05    8.69e-05    0.000102    0.000165
Bentoy13's sparse   8.9e-05     8.82e-05    9.23e-05    6.78e-05    6.44e-05    8.61e-05    0.000114    0.0002
Luis Mendo's optim. 3.99e-05    3.96e-05    4.08e-05    4.3e-05     4.61e-05    5.86e-05    7.66e-05    0.000111

N                   5e+03       1e+04       2e+04       5e+04       1e+05       2e+05       5e+05       1e+06
-------------------------------------------------------------------------------------------------------------
OP's for-loop       0.0138      0.0278      0.0588      0.16        0.264       0.525       1.35        2.73
OP's arrayfun       0.118       0.239       0.533       1.46        2.42        4.83        12.2        24.8
Daniel              0.00105     0.0021      0.00461     0.0138      0.0242      0.0504      0.126       0.264
Divakar             0.00127     0.00284     0.00655     0.0203      0.0335      0.0684      0.185       0.396
David's for-loop    0.015       0.0286      0.065       0.175       0.3         0.605       1.56        3.16
David's arrayfun    0.0668      0.129       0.299       0.803       1.33        2.64        6.76        13.6
Luis Mendo          0.000236    0.000446    0.000863    0.00221     0.0049      0.0118      0.0299      0.0637
Bentoy13's cumsum   0.000318    0.000638    0.00107     0.00261     0.00498     0.0114      0.0283      0.0526
Bentoy13's sparse   0.000414    0.000774    0.00148     0.00451     0.00814     0.0191      0.0441      0.0877
Luis Mendo's optim. 0.000224    0.000413    0.000754    0.00207     0.00353     0.00832     0.0216      0.0441

好的，我在列表中添加了另一个解决方案......我无法阻止自己优化 Luis Mendo 的最佳解决方案。不值得称赞，它只是 Luis Mendo 的变体，我稍后会解释。

显然，使用的解决方案arrayfun非常耗时。与其他解决方案相比，使用显式 for 循环的解决方案更快，但仍然很慢。所以是的，矢量化仍然是优化 Matlab 脚本的主要选项。

由于我看到最快解决方案的计算时间存在很大差异，尤其是输入长度在 100 到 10000 之间时，我决定更精确地进行基准测试。所以我把最慢的分开了（对不起），并在其他 6 个运行得更快的解决方案上重做基准测试。这个简化的解决方案列表的第二个基准测试是相同的，除了我的平均运行次数超过 1000 次。

（这里没有表，除非你真的想，否则和以前的数字完全一样）

正如所言，Daniel 的解决方案比 Divakar 的解决方案快一点，因为bsxfunwith @times 的使用似乎比 using 慢repmat。尽管如此，它们仍然比 for 循环解决方案快 10 倍：显然，在 Matlab 中进行矢量化是一件好事。

Bentoy13 和 Luis Mendo 的解决方案非常接近；第一个使用更多指令，但第二个在将 1 连接到cumsum(input_lengths(1:end-1)). 这就是为什么我们看到 Bentoy13 的解决方案在输入长度较大（大于 5.10^5）的情况下往往会更快一些，因为没有额外的分配。从这个考虑，我做了一个优化的解决方案，没有额外的分配；这是代码（如果他愿意，Luis Mendo 可以将这个代码放在他的答案中:)）：

result = zeros(1,sum(input_lengths));
result(1) = 1;
result(1+cumsum(input_lengths(1:end-1))) = 1;
result = cumsum(result);

欢迎任何改进意见。

Answer 3

更多的是评论，但我做了一些测试。我尝试了一个for循环和一个arrayfun，我测试了你的for循环和arrayfun版本。你的for循环是最快的。我认为这是因为它很简单，并且允许 JIT 编译进行最大程度的优化。我正在使用 Matlab，八度可能会有所不同。

和时间：

Solution:     With JIT   Without JIT  
Sam for       0.74       1.22    
Sam arrayfun  2.85       2.85    
My for        0.62       2.57    
My arrayfun   1.27       3.81    
Divakar       0.26       0.28    
Bentoy        0.07       0.06    
Daniel        0.15       0.16
Luis Mendo    0.07       0.06

所以 Bentoy 的代码真的很快，和 Luis Mendo 的速度几乎一模一样。而且我太依赖 JIT 了！

---

以及我尝试的代码

clc,clear
input_lengths = randi(20,[1 10000]);

% My for loop
tic()
C=cumsum(input_lengths);
D=diff(C);
results=zeros(1,C(end));
results(1,1:C(1))=1;
for i=2:length(input_lengths)
    results(1,C(i-1)+1:C(i))=i*ones(1,D(i-1));
end
toc()

tic()
A=arrayfun(@(i) i*ones(1,input_lengths(i)),1:length(input_lengths),'UniformOutput',false);
R=[A{:}];
toc()

Answer 4

result = zeros(1,sum(input_lengths));
result(cumsum([1 input_lengths(1:end-1)])) = 1;
result = cumsum(result);

这应该很快。并且内存使用量是最小的。

由于 Bentoy13，上述代码的优化版本（参见他非常详细的基准测试）：

result = zeros(1,sum(input_lengths));
result(1) = 1;
result(1+cumsum(input_lengths(1:end-1))) = 1;
result = cumsum(result);

Answer 5

这是@Daniel答案的一个轻微变体。该解决方案的关键是基于该解决方案。现在这个避免了repmat，所以这样它可能会更“矢量化”。这是代码 -

selector=bsxfun(@le,[1:max(input_lengths)]',input_lengths); %//'
V = bsxfun(@times,selector,1:numel(input_lengths)); 
result = V(V~=0)

对于所有绝望的单线搜索的人-

result = nonzeros(bsxfun(@times,bsxfun(@le,[1:max(input_lengths)]',input_lengths),1:numel(input_lengths)))

Answer 6

我搜索了一个优雅的解决方案，我认为大卫的解决方案是一个好的开始。我想到的是，可以生成从前一个元素添加一个的索引。

为此，如果我们计算cumsum输入向量的 ，我们得到：

cumsum(input_lengths)
ans = 1     2     3     7    10    12    13

这是相同数字序列末端的索引。这不是我们想要的，所以我们翻转向量两次以获得开始：

fliplr(sum(input_lengths)+1-cumsum(fliplr(input_lengths)))
ans = 1     2     3     4     8    11    13

这是诀窍。你翻转向量，cumsum它得到翻转向量的末端，然后翻转回来；但是您必须从输出向量的总长度中减去向量（+1，因为索引从 1 开始），因为 cumsum 适用于翻转的向量。

完成此操作后，它非常简单，您只需将 1 放在计算索引处，将 0 放在其他位置，然后对它求和：

idx_begs = fliplr(sum(input_lengths)+1-cumsum(fliplr(input_lengths)));
result = zeros(1,sum(input_lengths));
result(idx_begs) = 1;
result = cumsum(result);

编辑

首先，请看一下Luis Mendo 的解决方案，它与我的非常接近，但更简单、更快（即使非常接近，我也不会编辑我的）。我认为目前这是最快的解决方案。

其次，在查看其他解决方案时，我编造了另一个 one-liner，与我最初的解决方案和另一个 one-liner 略有不同。好的，这不是很可读，所以深呼吸：

result = cumsum( full(sparse(cumsum([1,input_lengths(1:end-1)]), ...
ones(1,length(input_lengths)), 1, sum(input_lengths),1)) );

我把它剪成两条线。好的，现在让我们解释一下。

类似的部分是构建索引数组，以增加当前元素的值。为此，我使用 Luis Mendo 的解决方案。为了在一行中构建解向量，我在这里使用了一个事实，即它实际上是二进制向量的稀疏表示，我们将在最后进行 cumsum。这个稀疏向量是使用我们计算的索引向量作为 x 位置、1 的向量作为 y 位置以及 1 作为放置在这些位置的值来构建的。第四个参数用于精确计算向量的总大小（如果 的最后一个元素input_lengths不是 1 则很重要）。然后我们得到这个稀疏向量的完整表示（否则结果是一个没有空元素的稀疏向量），我们可以求和。

除了为这个问题提供另一个解决方案之外，这个解决方案没有任何用处。由于内存负载较重，基准测试可以显示它比我的原始解决方案慢。