代码之家  ›  专栏  ›  技术社区  ›  Tommaso Belluzzo

迭代模拟的严重性能问题

  •  5
  • Tommaso Belluzzo  · 技术社区  · 8 年前

    我最近在实现一个模拟算法时偶然发现了一个性能问题。我找到了瓶颈函数(很明显,它是对 arrayfun 这会减慢速度):

    function sim = simulate_frequency(the_f,k,n)
    
        r = rand(1,n); % 
        x = arrayfun(@(x) find(x <= the_f,1,'first'),r);
        sim = (histcounts(x,[1:k Inf]) ./ n).';
    
    end
    

    它被用于代码的其他部分,如下所示:

    h0 = zeros(1,sims);
    
    for i = 1:sims
        p = simulate_frequency(the_f,k,n);
        h0(i) = max(abs(p - the_p));
    end
    

    以下是一些可能的值:

    % Test Case 1
    sims = 10000;
    the_f = [0.3010; 0.4771; 0.6021; 0.6990; 0.7782; 0.8451; 0.9031; 0.9542; 1.0000];
    k = 9;
    n = 95;
    
    % Test Case 2
    sims = 10000;
    the_f = [0.0413; 0.0791; 0.1139; 0.1461; 0.1760; 0.2041; 0.2304; 0.2552; 0.2787; 0.3010; 0.3222; 0.3424; 0.3617; 0.3802; 0.3979; 0.4149; 0.4313; 0.4471; 0.4623; 0.4771; 0.4913; 0.5051; 0.5185; 0.5314; 0.5440; 0.5563; 0.5682; 0.5797; 0.5910; 0.6020; 0.6127; 0.6232; 0.6334; 0.6434; 0.6532; 0.6627; 0.6720; 0.6812; 0.6901; 0.6989; 0.7075; 0.7160; 0.7242; 0.7323; 0.7403; 0.7481; 0.7558; 0.7634; 0.7708; 0.7781; 0.7853; 0.7923; 0.7993; 0.8061; 0.8129; 0.8195; 0.8260; 0.8325; 0.8388; 0.8450; 0.8512; 0.8573; 0.8633; 0.8692; 0.8750; 0.8808; 0.8864; 0.8920; 0.8976; 0.9030; 0.9084; 0.9138; 0.9190; 0.9242; 0.9294; 0.9344; 0.9395; 0.9444; 0.9493; 0.9542; 0.9590; 0.9637; 0.9684; 0.9731; 0.9777; 0.9822; 0.9867; 0.9912; 0.9956; 1.000];
    k = 90;
    n = 95;
    

    标量 sims 必须在范围内 1000 1000000 .累积频率矢量 the_f 从不包含超过 100 元素。标量 k 表示中的元素数 阿芙 .最后,标量 n 表示经验样本向量中的元素数,甚至可以非常大(最多 10000 元素,据我所知)。

    关于如何提高这个过程的计算时间有什么线索吗?

    3 回复  |  直到 8 年前
        1
  •  5
  •   Anthony    8 年前

    如果我想要表演,我会避免 arrayfun 是的。甚至这个 for 循环速度更快:

    function sim = simulate_frequency(the_f,k,n)
    
        r = rand(1,n); % 
        for i = 1:numel(r)
            x(i) = find(r(i)<the_f,1,'first');
        end
        sim = (histcounts(x,[1:k Inf]) ./ n).';
    
    end
    

    使用第一组样本数据运行10000个sims会给出以下计时。

    你的 阵列风 功能:

    >Elapsed time is 2.848206 seconds.
    

    这个 对于 循环函数:

    >Elapsed time is 0.938479 seconds.
    

    受Cris Luengo回答的启发,我建议如下:

    function sim = simulate_frequency(the_f,k,n)
    
        r = rand(1,n); % 
        x = sum(r > the_f)+1;
        sim = (histcounts(x,[1:k Inf]) ./ n)';
    
    end
    

    时间:

    >Elapsed time is 0.264146 seconds.
    
        2
  •  6
  •   Cris Luengo    8 年前

    在第二个测试用例中,而不是第一个测试用例中,这对我来说似乎要快一点。时间差可能会更长 the_f 以及更大的 n 是的。

    function sim = simulate_frequency(the_f,k,n)
        r = rand(1,n); % 
        [row,col] = find(r <= the_f); % Implicit singleton expansion going on here!
        [~,ind] = unique(col,'first');
        x = row(ind);
        sim = (histcounts(x,[1:k Inf]) ./ n).';
    end
    

    我使用隐式单重展开 r <= the_f ,使用 bsxfun 如果你有一个旧版本的Matlab(但是你知道这个练习)。

    find然后将行和列返回到 r 大于 阿芙 是的。 unique 在每个列的第一个元素的结果中查找索引。

    学分: Andrei Bobrov over on MATLAB Answers


    另一个选项(源自 this other answer )是有点短,但也有点模糊我:

    mask = r <= the_f;
    [x,~] = find(mask & (cumsum(mask,1)==1));
    
        3
  •  2
  •   rahnema1    8 年前

    你可以用 histcounts 具有 r 作为输入:

    r = rand(1,n);
    sim = (histcounts(r,[-inf ;the_f]) ./ n).';
    

    如果 histc 是用来代替 历史计数 整个模拟可以矢量化:

    r = rand(n,sims);
    p = histc(r, [-inf; the_f],1);
    p = [p(1:end-2,:) ;sum(p(end-1:end,:))]./n;
    h0 = max(abs(p-the_p(:)));    %h0 = max(abs(bsxfun(@minus,p,the_p(:))));