algorithm - 向量化帕累托前沿算法

Question

首先，这是我的设置:

• x是一个 n x 1包含第一个成本函数值的向量。
• y是另一个n x 1包含第二个成本函数值的向量。
• a是一个 m x 1包含 x 索引的向量和 y要检查，这用于有选择地从算法中排除值。除非需要，这可以替换为 1:n .
• 假设所有组合都是安全的 (x,y)是独一无二的。

任务是找到帕累托最优值对集 (x,y) ，即所有未被支配的对。如果存在另一对 (u,v)，则称一对被支配这样 u <= x && v <= y其中一个比较是严格的:u < x || v < y .换句话说，如果另一对在一个值上更好并且在另一个值上不差，则一对是主导的。

到目前为止，我的研究已经产生了三种可行的算法，不幸的是，它们都依赖于循环。这是它们的工作方式以及我用 x 运行它们的时间, y和 a长度1e8 :

1. Sort x in ascending order. Add the first pair to the pareto set.
2. Loop through x. Add each pair to the pareto set where y is lower than the previous pareto pair's y.

Elapsed time is 80.204052 seconds.

1. Find min(x). Add that pair to the pareto set.
2. Select all pairs where y is lower than the previously added pair's y.
3. Go to step 1 unless step 2 resulted in the empty set.

Elapsed time is 2.993350 seconds.

1. Loop through all pairs (x,y).
2. Remove all pairs (u,v) with x >= u && y >= v.

Elapsed time is 105.924814 seconds.

现在我要做的是创建一个矢量化算法。它不必基于上述之一，但我无法找到任何其他工作算法。我能做的最好的是:

ap = a(y < min(y(x == min(x))) | x < min(x(y == min(y))));

通常会找到所有帕累托最优对，但包括所有不被 min(x) 处的对支配的对或 min(y) ，即使一个人支配另一个人。我说通常是因为如果只有一个全局最优对支配所有其他对，它就会完全失败。替换 <与 <=解决了第二个问题，但找到了更多的支配对(那些只有一个更差值的对)。我还通过与上面相同的计时器运行了它:

Elapsed time is 0.800385 seconds.

---

这是我用来检查算法效果的测试脚本，请随意使用它

for i=1:25
    x = randi(8,10,1);
    y = randi(8,10,1);
    a = 1:10;
    ap = a(y < min(y(x == min(x))) | x < min(x(y == min(y)))); %// algorithm here
    figure(1);
    subplot(5,5,i);
    plot(a,x,'b',a,y,'r',ap,x(ap),'b.',ap,y(ap),'r.','MarkerSize',20);
    axis([0,11,0,9]);
    set(gca,'XGrid','on','YGrid','on','XTick',1:10,'YTick',0:8);
    figure(2);
    subplot(5,5,i);
    plot(x,y,'b.',x(ap),y(ap),'ro','MarkerSize',10);
    axis([0,9,0,9]);
end

Answer 1

因此，如果速度是主要特征(在正确性之后)，那么我发现更快循环版本的递归版本要快 30% 以上:

>> testPareto(1e8);
Recursive:
Elapsed time is 4.507267 seconds.
Loop:
Elapsed time is 6.136856 seconds.
Vector:
Elapsed time is 7.246806 seconds.

同样，时间取决于机器，甚至可能取决于 matlab 的版本。这是代码:

function testPareto(dim)

x = rand(dim, 1);
y = rand(dim, 1);

tic;
rR = paretoRecursive(x, y);
disp('Recursive:');
toc;

tic;
rL = paretoLoop(x, y);
disp('Loop:');
toc;

tic;
rV = paretoVector(x, y);
disp('Vector:');
toc;

end

function result = paretoLoop(x, y)
    result = zeros(numel(x), 2);
    off = 1;
    loop = true;
    while loop
        xmin = min(x);
        ymin = min(y(x == xmin));
        yfilter = y < ymin;
        result(off, :) = [xmin ymin];
        off = off + 1;
        if any(yfilter)
            x = x(yfilter);
            y = y(yfilter);
        else
            loop = false;
            result(off:end, :) = [];
        end
    end
end

function result = paretoRecursive(x, y)
    xmin = min(x);
    ymin = min(y(x == xmin));
    yfilter = y < ymin;
    if any(yfilter)
        result = [xmin ymin; paretoRecursive(x(yfilter), y(yfilter))];
    else
        result = [xmin ymin];
    end
end

function result = paretoVector(x, y)
    xmin = min(x);
    xfilter = x == xmin;
    ymin = min(y(xfilter));
    yfilter = y < ymin;
    if any(yfilter)
        [x, ind] = sort(x(yfilter));
        y = y(yfilter);
        y = y(ind);
        yfilter = [true; y(2:end) < cummin(y(1:end-1))];
        result = [xmin x(yfilter)'; ymin y(yfilter)']';
    else
        result = [xmin ymin];
    end
end