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

所有解决方案的基准

按照之前的基准测试,我将脚本中给出的所有解决方案分组并运行几个小时以获得基准测试.我之所以这样做是因为我觉得很高兴看到每个提出的解决方案的性能与输入长度作为参数 - 我的意图不是要降低前一个解决方案的质量,这提供了有关其效果的其他信息. JIT.而且,每个参与者似乎都同意这一点,在所有答案中都做了相当不错的工作,所以这个伟大的帖子值得结论.

我不会在这里发布脚本代码,这很长很无趣.基准测试的过程是针对一组不同长度的输入向量运行每个解决方案: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的最佳解决方案.没有信用,它只是路易斯·门多的一个变种,我稍后会解释.

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

由于我看到最快解决方案的计算时间有很大差异,特别是输入长度在100到10000之间,我决定更准确地进行基准测试.所以我把最慢的分开(对不起),并重新运行其他6个运行速度更快的解决方案.这个减少的解决方案列表的第二个基准是相同的,除了我平均超过1000次运行.

(这里没有表格,除非你真的想要,它和以前的数字完全相同)

正如有人所说,Daniel的解决方案比Divakar的解决方案快一点,因为似乎使用bsxfun@times比使用慢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

这是诀窍.你翻转矢量,收集它以得到翻转矢量的末端,然后翻转; 但是你必须从输出向量的总长度中减去向量(+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的解决方案,它非常接近我的但更简单,更快一点(即使它非常接近,我也不会编辑我的).我想在这个日期,这是所有人中最快的解决方案.

其次,虽然看着别人的解决方案,我做了一个又一个内胆,一点点不同,我最初的解决方案,并从其他的单行.好吧,这不会是非常易读的,所以请一试:

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

我把它切成两条线.好的,现在让我们解释一下.

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

除了给出这个问题的另一种解决方案之外,没有使用这种解决方案.基准测试可以显示它比我原来的解决方案慢,因为内存负载较重.