重复数组元素的副本:在MATLAB中运行长度解码

我试图插入多个值到一个数组使用'值'数组和'计数器'数组。 例如,如果:

a=[1,3,2,5] b=[2,2,1,3] 

我想要一些函数的输出

 c=somefunction(a,b) 

成为

 c=[1,1,3,3,2,5,5,5] 

(1)重复b(1)次,a(2)重复b(2)次,等等

MATLAB中有内置函数吗? 如果可能,我想避免使用for循环。 我已经试过了'repmat()'和'kron()'的变种无济于事。

这基本上是Run-length encoding

    问题陈述

    我们有一系列值, valsrunlensrunlens

     vals = [1,3,2,5] runlens = [2,2,1,3] 

    我们需要在runlens每个对应元素的vals次重复每个元素。 因此,最终的输出是:

     output = [1,1,3,3,2,5,5,5] 

    前瞻性的方法

    MATLAB中速度最快的工具之一是cumsum ,在处理不规则模式的矢量化问题时非常有用。 在陈述的问题中, runlens的不同元素会带来不规则性。

    现在,为了利用cumsum ,我们需要在这里做两件事:初始化一个zeros数组,并将“合适的”值放置在零数组上的“关键”位置,使得在应用“ cumsum ”之后, runlens时间的重复vals的最后阵列。

    步骤:让我们对上述步骤进行编号,使预期的方法更容易:

    1)初始化零数组:长度是多少? 由于我们重复runlens时间,零数组的长度必须是所有runlens的总和。

    2)查找关键位置/索引:现在这些关键位置是沿着零点阵列的位置,其中来自vals每个元素开始重复。 因此,对于runlens = [2,2,1,3] ,映射到零数组上的关键位置将是:

     [X 0 X 0 XX 0 0], where X's are those key positions. 

    3)找到合适的价值:在使用cumsum之前最后打的钉子是将“合适”的数值放到这些关键位置。 现在,既然我们会在不久之后做cumsum ,如果你仔细考虑,你需要differentiateddiff values版本,这样cumsum就会带回我们的values 。 由于这些差异值将被放置在由runlens距离分隔的位置处的零点阵列上,因此在使用cumsum我们将每个vals元素重复runlens次数作为最终输出。

    解决方案代码

    这是实施上述所有步骤的实现 –

     %// Calculate cumsumed values of runLengths. %// We would need this to initialize zeros array and find key positions later on. clens = cumsum(runlens) %// Initalize zeros array array = zeros(1,(clens(end))) %// Find key positions/indices key_pos = [1 clens(1:end-1)+1] %// Find appropriate values app_vals = diff([0 vals]) %// Map app_values at key_pos on array array(pos) = app_vals %// cumsum array for final output output = cumsum(array) 

    预先分配Hack

    可以看出,上面列出的代码使用了预分配零。 现在,根据这个UNDOCUMENTED MATLAB blog on faster pre-allocation ,可以实现更快的预分配,

     `array(clens(end)) = 0` instead of `array = zeros(1,(clens(end)))` 

    包装:功能代码

    总结一切,我们将有一个紧凑的功能代码来实现这种运行长度解码,

     function out = rle_cumsum_diff(vals,runlens) clens = cumsum(runlens); idx(clens(end))=0; idx([1 clens(1:end-1)+1]) = diff([0 vals]); out = cumsum(idx); return; 

    标杆

    基准代码

    下面列出的是基准代码,用于比较本文中所述的cumsum+diff方法的运行时间和加速比在MATLAB 2014B上的其他cumsum-only方法 –

     datasizes = [reshape(linspace(10,70,4).'*10.^(0:4),1,[]) 10^6 2*10^6]; %//' fcns = {'rld_cumsum','rld_cumsum_diff'}; %// approaches to be benchmarked for k1 = 1:numel(datasizes) n = datasizes(k1); %// Create random inputs vals = randi(200,1,n); runs = [5000 randi(200,1,n-1)]; %// 5000 acts as an aberration for k2 = 1:numel(fcns) %// Time approaches tsec(k2,k1) = timeit(@() feval(fcns{k2}, vals,runs), 1); end end figure, %// Plot runtimes loglog(datasizes,tsec(1,:),'-bo'), hold on loglog(datasizes,tsec(2,:),'-k+') set(gca,'xgrid','on'),set(gca,'ygrid','on'), xlabel('Datasize ->'), ylabel('Runtimes (s)') legend(upper(strrep(fcns,'_',' '))),title('Runtime Plot') figure, %// Plot speedups semilogx(datasizes,tsec(1,:)./tsec(2,:),'-rx') set(gca,'ygrid','on'), xlabel('Datasize ->') legend('Speedup(x) with cumsum+diff over cumsum-only'),title('Speedup Plot') 

    rld_cumsum.m相关功能代码:

     function out = rld_cumsum(vals,runlens) index = zeros(1,sum(runlens)); index([1 cumsum(runlens(1:end-1))+1]) = 1; out = vals(cumsum(index)); return; 

    运行时间和加速计划

    在这里输入图像描述

    在这里输入图像描述

    结论

    提出的方法似乎给了我们一个明显的加速超过cumsum-only方法,这是约3x

    为什么这个新的基于cumsum+diff的方法比以前的cumsum-only方法更好?

    那么,理由的精髓在于需要将“累积”值映射到valscumsum-only方法的最后一步。 在新的基于cumsum+diff的方法中,我们正在执行diff(vals)而对于其中的MATLAB只处理n元素(其中n是runLengths的数量)与sum(runLengths)元素的映射相比, cumsum-only这个方法,这个数字必须比n多很多倍,因此这个新方法的显着加速!

    基准

    针对R2015b进行了更新 :对于所有数据大小而言, repelem速度最快。


    测试功能:

    1. repelem中增加了MATLAB的内置repelem功能
    2. gnovice的cumsum解决方案( rld_cumsum
    3. Divakar的cumsum + diff解决方案( rld_cumsum_diff
    4. knedlsepp的accumarray解决方案( knedlsepp5cumsumaccumarray )从这个职位
    5. 天真的基于循环的实现( naive_jit_test.m )来测试即时编译器

    test_rld.mtest_rld.m结果b

    重复时间

    旧的时间情节使用R2015 这里 。

    调查结果

    • repelem总是最快的大约是2倍。
    • rld_cumsum_diff始终比rld_cumsum快。
    • 对于小数据量(小于大约300-500个元素), repelem是最快的,
    • rld_cumsum_diff变得比repelem 5000个元素快得多
    • 在三万到三十万个元素之间, repelem变得比rld_cumsum
    • rld_cumsumrld_cumsum具有大致相同的性能
    • naive_jit_test.m具有几乎恒定的速度,与rld_cumsumknedlsepp5cumsumaccumarray对于较小的尺寸,对于较大的尺寸稍快

    在这里输入图像描述

    老利率情节使用R2015 这里 。

    结论

    在下面大约5000个元素和cumsum + diff解上面使用repelem

    没有我知道的内置函数,但是这里有一个解决方案:

     index = zeros(1,sum(b)); index([1 cumsum(b(1:end-1))+1]) = 1; c = a(cumsum(index)); 

    说明:

    首先创建与输出数组长度相同的矢量(即b中所有复制的总和)。 然后将一个元素放置在第一个元素中,每个后续元素代表输出中新值序列的开始位置。 然后可以使用矢量index的累加和来索引a ,将每个值复制为期望的次数。

    为了清楚起见,这是问题中给出的ab的值的各种向量的样子:

      index = [1 0 1 0 1 1 0 0] cumsum(index) = [1 1 2 2 3 4 4 4] c = [1 1 3 3 2 5 5 5] 

    编辑:为了完整起见,有另一种使用ARRAYFUN的替代方案,但这似乎需要20-100倍以上的时间运行比上述解决方案向量多达10,000个元素长:

     c = arrayfun(@(x,y) x.*ones(1,y),a,b,'UniformOutput',false); c = [c{:}]; 

    最后有(作为R2015a )内置和记录的功能来做到这一点, repelem 。 下面的语法,其中第二个参数是一个向量,在这里是相关的:

    W = repelem(V,N)与向量V和向量N创建一个向量W ,其中元素V(i)被重复N(i)次。

    换句话说,“ N每个元素指定重复V的对应元素的次数”。

    例:

     >> a=[1,3,2,5] a = 1 3 2 5 >> b=[2,2,1,3] b = 2 2 1 3 >> repelem(a,b) ans = 1 1 3 3 2 5 5 5 

    MATLAB的内置repelem的性能问题在repelem已经被固定。 我已经从chappjc的R2015b文章中运行了test_rld.m程序,现在repelem比其他算法快了大约2倍:

    比较不同方法的重复执行时间的更新情节在cumsum + diff上重复加速(约2倍)