performance - 通过数组操作高效搜索包含子排列的排列？

Question

我有一组整数，比如说 S = {1,...,10}，还有两个矩阵 N 和 M，它们的行是 S 阶元素的一些(但不一定是所有可能的)排列，比如说, 3 和 5 分别，例如N = [1 2 3; 2 5 3;...], M = [1 2 3 4 5; 2 4 7 8 1;...].

排列 P 的子排列 Q 只是 P 的索引子集，使得 Q 元素的索引顺序与它们在 P 中的索引顺序相同。示例:[2,4, 7]是[2,3,4,6,7,1]的子排列，但[1,2,3]不是后者的子排列。

我需要一种有效的方法(例如尽可能矢量化和尽可能小的 for 循环)来查找

(1) M 的所有排列具有 N 的子排列

和

(2) N 的每个子排列要在 M 中找到多少次。

到目前为止，我拥有的是一个矢量化代码，用于检查给定的单个子排列是否包含在 M 中(以及多少次)，但随后我必须使用 N 的 parfor 循环，这会变慢非常大的 N-s。请注意，如果 N 不是太大，也可以通过简单地从给定的 3 元组构造可接受的 5 元组并将结果与​​ M 进行比较来解决问题，但是如果 N足够大。

查看问题的另一种方法如下:检查其行的 N 模排列是否是一般意义上的 M 的子矩阵，即是否有可能通过删除 N 的行来获得排列来自 M 的元素。

抱歉，如果我的问题太初级，我的背景是算术代数几何和表示论，而且我对 MATLAB 还很陌生。

编辑: 这是我检查 M 中是否存在单个 k 元组的代码: [代码]

function [A,f] = my_function(x,M)
%// returns all rows in M that contain x and the absolute frequency of x in M
%// suboptimal for checking combinations rather than permutations byy at least ~ 50%
k = size(x,2);
m = size(M,1);
R = zeros(m,k);
I = R;
Z = I;
    for j = 1:k   
        [R(:,j),I(:,j)] = max((M == x(j)),[],2); 
        Z(:,j) = R(:,j).*I(:,j);
    end
z = zeros(m,k-1);
    for j = 1:(k-1)
        z(:,j) = (Z(:,j) > 0 & Z(:,j) < Z(:,j+1)); 
    end
[v,~] = find(sum(z,2) == k-1);    
A = M(v,:);
f = length(v);
end

使用此函数，检查 N 只是一个简单的 (par)for 循环问题，我希望避免这种情况以支持更快的矢量化解决方案。

Answer 1

方法 #1

[val,ind] = max(bsxfun(@eq,permute(M,[4 2 1 3]),permute(N,[2 3 4 1])),[],2)
matches = squeeze(all(diff(ind,1)>0,1).*all(val,1))
out1 = any(matches,2) %// Solution - 1
out2 = sum(matches,1) %// Solution - 2

---

方法 #2

另一种避免排列 N 的方法可能更适合较长的 N -

[val,ind] = max(bsxfun(@eq,N,permute(M,[3 4 1 2])),[],4)
matches = squeeze(all(diff(ind,[],2)>0,2).*all(val,2))
out1 = any(matches,1) %// Solution - 1
out2 = sum(matches,2) %// Solution - 2

---

方法 #3

大数据量的内存 scroogey 方法 -

out1 = false(size(M,1),1);  %// Storage for Solution - 1
out2 = zeros(size(N,1),1);  %// Storage for Solution - 2
for k=1:size(N,1)
    [val3,ind3] = max(bsxfun(@eq,N(k,:),permute(M,[1 3 2])),[],3);
    matches = all(diff(ind3,[],2)>0,2).*all(val3,2);
    out1 = or(out1,matches);
    out2(k) = sum(matches);
end

---

方法 #4

GPU 的内存节约方法 -

gM = gpuArray(M);
gN = gpuArray(N);

gout1 = false(size(gM,1),1,'gpuArray');  %// GPU Storage for Solution - 1
gout2 = zeros(size(gN,1),1,'gpuArray');  %// GPU Storage for Solution - 2
for k=1:size(gN,1)
    [val3,ind3] = max(bsxfun(@eq,gN(k,:),permute(gM,[1 3 2])),[],3);
    matches = all(diff(ind3,[],2)>0,2).*all(val3,2);
    gout1 = or(gout1,matches);
    gout2(k) = sum(matches);
end
out1 = gather(gout1);  %// Solution - 1
out2 = gather(gout2);  %// Solution - 2

现在，这种 GPU 方法已经超越了所有其他方法。它是用 M : 320000X5 和 N : 2100X3 (与您的输入大小相同)填充随机整数运行的。使用 GTX 750 Ti，只用了 13.867873 秒!!因此，如果您的 GPU 具有足够的内存，这也可能是您的成功方法。

---

方法 #5

GPU 的极端内存消耗方法 -

gM = gpuArray(M);
gN = gpuArray(N);

gout1 = false(size(gM,1),1,'gpuArray');  %// GPU Storage for Solution - 1
gout2 = zeros(size(gN,1),1,'gpuArray');  %// GPU Storage for Solution - 2
for k=1:size(gN,1)
    [val2,ind2] = max(bsxfun(@eq,gM,permute(gN(k,:),[1 3 2])),[],2);
    matches = all(diff(ind2,[],3)>0,3).*all(val2,3);
    gout1 = or(gout1,matches);
    gout2(k) = sum(matches);
end
out1 = gather(gout1);  %// Solution - 1
out2 = gather(gout2);  %// Solution - 2