-2

假设我有这个工作示例：

lover_bound = 10;
upper_bound = 180;
steps = 10;
NumeroCestelli = 8; 

livello = [lover_bound:steps:upper_bound];
L = length(livello);
n_c = ceil((factorial(L+NumeroCestelli-1))/(factorial(NumeroCestelli)*factorial(L-1)));

randIdxs = randi([1,L],n_c,NumeroCestelli);
PianoSperimentale = single(livello(randIdxs));

我需要执行一个n_c x NumeroCestelli矩阵（称为PianoSperimentale），其中每一行都是唯一的。不允许任何形式的排列。使用 randi 我无法执行我的要求。

[10 20 30 40 50 60 70 80] is equal to [80 70 60 50 40 30 20 10]

PianoSperimentale应该是一个1081575x8矩阵。过去我使用的是Combinator ) 函数，但对于非常大的矩阵来说非常慢。

[PianoSperimentale] = combinator(L,NumeroCestelli,'c','r');

for i=1:L
    PianoSperimentale(PianoSperimentale==i)=livello(i);
end

那么，有一种方法可以快速执行相同的矩阵combinator吗randi？

编辑：我允许选择相同的数字两次（NumberOfCombinations = (NumeroCestelli+L-1)!/(NumeroCestelli!(L-1)!）

建议的编辑

我需要生成从 18 个元素的向量中选择任何 8 个数字时获得的完整组合（带有重复项）。这可以通过使用Combinator函数来完成，但对于非常大的矩阵来说非常慢。任何人都可以建议一种更快的方法来生成这个吗？

示例：使用“来自 4 的向量的样本 3”将产生以下结果：

1 1 1
1 1 2
1 1 3
1 1 4
1 2 2
1 2 3
1 2 4
1 3 3
1 3 4
1 4 4
2 2 2
2 2 3
2 2 4
2 3 3
2 3 4
2 4 4
3 3 3
3 3 4
3 4 4
4 4 4

我知道对于从 18 个元素中选择 8 个元素的向量，我将得到总共(18+8-1)!/8!*(18-1)!可能的组合，或 1081578 行 8 个值。谁能帮我找到一个快速的算法来做到这一点？

4

2 回答 2

2

出于“历史”的原因，我正在写一个新的答案，而不是删除我的旧答案（这个问题的很多来回都是理解这个问题的必要序言，所以这个答案甚至是有意义的）。TL;DR：最后的完整代码。

这是一个非常棘手的问题，但我想我已经明白了。关键的见解是，您对结果矩阵中正确数量的元素的表达式(L+H-1)!/(H!(L-1)!)强烈建议“从 L + H - 1 中选择 H”与您的问题的解决方案之间存在关系。Trick 正在寻找这种关系。我首先写出结果combnk(5, 3)（在这种大小下，您可以手动写出所有组合并寻找模式）：

1 2 3
1 2 4
1 2 5
1 3 4
1 3 5
1 4 5
2 3 4
2 3 5
2 4 5
3 4 5

我们如何将其转换为1 2 3（包括重复）的独特组合？我注意到有三组连续数字：

1 2 3
2 3 4
3 4 5

这给了我一个想法，即我需要对连续数字的差异做一些事情——不知何故，如果差异是 1，我需要重复这个数字。这种洞察力很快导致了以下代码：

L = 3; % pick three numbers
H = 3; % from three numbers: 1,2,3

a = combnk(1:L+H-1, L);  % generate all "combinations" of 1,2,3,4,5 without repeats

% the "magic" line: compress into "combinations with repeats"
b = cumsum([a(:,1)  diff(a,[],2) - 1],2);

对于上面的例子，这给出了

 1 1 1
 1 1 2
 1 1 3
 1 2 2
 1 2 3
 1 3 3
 2 2 2
 2 2 3
 2 3 3
 3 3 3

那是怎么发生的？嗯，diffa（沿着第二维）是

 1 1
 1 2
 1 3
 2 1
 2 2
 3 1
 1 1
 1 2
 2 1
 1 1

diff(a,[],2)-1)也是_

 0 0
 0 1
 0 2
 1 0
 1 1
 2 0
 0 0
 0 1
 1 0
 0 0

在这个表达式中，0表示“重复最后一个数字”，而1表示“加 1”，2表示“加 2”。我们可以通过使用cumsum(cumulative sum) 函数并从组合的第一个数字开始来完成所有这些相加。这导致表达式

b = cumsum([a(:, 1) diff(a, [], 2) - 1]);

作为最后一步，您必须将其转换回您正在使用的索引。您的完整代码将是

L = 8;
H = 18;
a = combnk(1:L+H-1, L);
b = cumsum([a(:,1)  diff(a,[],2) - 1],2);
livello = 10:10:180;
PianoSperimentale = livello(b);

这将创建一个大小为 的数组，其中b的值是livello。

我相信这对你有用（我无法测试这个，因为我的家用电脑上没有 Matlab），它会尽可能快地解决这个问题。

于 2013-11-08T23:03:23.487 回答
0

我正在尝试解析您的问题。从“无排列”和“10 20 30 40 50 60 70 80”等价于“80 70 60 50 40 30 20 10”，我推断您需要从 18 中抽取 8 个数字，而不是两次相同的样本。

这意味着您要生成所有可能的组合（来自 18 个可能值的 8 个样本），并从中进行选择。只有 43758 种方法可以做到这一点；之后，您将不得不包括排列。因此，上述问题（如果我理解正确的话）无法解决。

编辑现在问题已更新，我认为以下将是一个解决方案：

lover_bound = 10;
upper_bound = 180;
steps = 10;
NumeroCestelli = 8; 

livello = [lover_bound:steps:upper_bound];
livello = [livello livello];

PianoSperimentale = combnk(livello, 8);

因为每个数字出现两次，你可以让它重复。不幸的是，这将允许多个“双打”（例如 [10 10 20 20 30 30 40 40] 将被允许），这将比您计算的表达式大得多（即 36!/(28!8!) ~ 30M ）。一种可能的方法（最多允许一个双精度）是

livello = lover_bound:steps:upper_bound;

for ii = 1:numel(livello)
  PS(ii,:,:) = combnk([livello livello(ii)], 8);
end

PianoSperimentale = reshape(PS, [], 8);

这允许“每个循环一个重复的数字”，并且我相信这更接近您的想法，尽管组合的数量会18 * (19!/(19-8)!8!) = 1360476比您的表达式大一点。我现在无法测试这个，因为我这台电脑上没有 Matlab ......

于 2013-11-08T15:46:57.367 回答