根据有放回的组合计算列的乘积
Calculate products of columns according to combinations with replacement
问题
解释起来有点困难，但我会尽力而为。我知道找到替换组合数的方程式。假设我有 6 个向量：A、B、C、D、E、F。如果我想找到这 6 个变量的所有可能的三次乘积，它将是 (6+3-1)!/3!(6- 1）！ = 56 种组合（见结尾）。同样，如果我想要每个二次积，它是 21。对于线性，当然是 6（只是每个变量本身）。我想计算所有 6+21+56 = 83 种组合。我正在考虑 3 个循环，每个内部循环都从其外部循环开始迭代，如
for i1=1:6
X(:,?) = X.*X(:,i1)
for i2=i1:6
X(:,?) = X.*X(:,i2)
for i3=i2:6
X(:,?) = X.*X(:,i3)
但是左侧存储所有数据的83列矩阵的索引让我感到困惑。如您所见，它们标有问号。
PS：可能也需要对第 5 阶执行此操作，因此它会添加另外 126 和 252 列，总共 461 列。因此，不对三阶进行硬编码的更通用的代码更好。但如果硬编码到第 5 位也没关系，因为我绝对不会超过那个。
MATLAB 或 Python 都可以，因为我可以在两者之间轻松切换。
计算出的二次组合实例
这是我期望的 6 个变量（A 到 F）的二次组合的 21 列示例。在 Excel 中完成。我为每个向量取了 3 个样本。
立方组合列表
这里是我需要计算的56个组合：
一个，一个，一个
A,A,B
A,A,C
A,A,D
A,A,E
A,A,F
A,B,B
A、B、C
A,B,D
A,B,E
A,B,F
A,C,C
A、C、D
A、C、E
A,C,F
A,D,D
A,D,E
A,D,F
A,E,E
A,E,F
A,F,F
B,B,B
B,B,C
B,B,D
B,B,E
B,B,F
B,C,C
B,C,D
B,C,E
B,C,F
B,D,D
B,D,E
B,D,F
B,E,E
B,E,F
B,F,F
C,C,C
C,C,D
C,C,E
C,C,F
C,D,D
C,D,E
C,D,F
C,E,E
C,E,F
C,F,F
D,D,D
D,D,E
D,D,F
D,E,E
D,E,F
D,F,F
E,E,E
E,E,F
E,F,F
F,F,F
您可以通过使用计数器来避免混淆索引：
clear all; close all
% Original matrix
M = [
2 2 3 2 8 8;
5 1 7 9 4 4;
4 1 2 7 2 9
];
% Number of combinations
order = 3;
sizeX = nchoosek(size(M,2)+order-1,order);
% Combinations
imat = ones(sizeX,order);
for c=2:sizeX
imat(c,:) = imat(c-1,:);
for o=order:-1:1
if (imat(c-1,o)<size(M,2))
imat(c,o:end) = imat(c-1,o)+1;
break
end
end
end
% Transpose & display combinations
imat = transpose(imat)
% Computations of products
X = ones(size(M,1),sizeX);
for o=1:order
X = X.*M(:,imat(o,:));
end
% Display result
X
当您执行脚本时，您会得到：
>> test_script
imat =
Columns 1 through 16
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 2 2 2 2 2 3 3 3 3 4
1 2 3 4 5 6 2 3 4 5 6 3 4 5 6 4
Columns 17 through 32
1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2
4 4 5 5 6 2 2 2 2 2 3 3 3 3 4 4
5 6 5 6 6 2 3 4 5 6 3 4 5 6 4 5
Columns 33 through 48
2 2 2 2 3 3 3 3 3 3 3 3 3 3 4 4
4 5 5 6 3 3 3 3 4 4 4 5 5 6 4 4
6 5 6 6 3 4 5 6 4 5 6 5 6 6 4 5
Columns 49 through 56
4 4 4 4 5 5 5 6
4 5 5 6 5 5 6 6
6 5 6 6 5 6 6 6
X =
Columns 1 through 16
8 8 12 8 32 32 8 12 8 32 32 18 12 48 48 8
125 25 175 225 100 100 5 35 45 20 20 245 315 140 140 405
64 16 32 112 32 144 4 8 28 8 36 16 56 16 72 196
Columns 17 through 32
32 32 128 128 128 8 12 8 32 32 18 12 48 48 8 32
180 180 80 80 80 1 7 9 4 4 49 63 28 28 81 36
56 252 16 72 324 1 2 7 2 9 4 14 4 18 49 14
Columns 33 through 48
32 128 128 128 27 18 72 72 12 48 48 192 192 192 8 32
36 16 16 16 343 441 196 196 567 252 252 112 112 112 729 324
63 4 18 81 8 28 8 36 98 28 126 8 36 162 343 98
Columns 49 through 56
32 128 128 128 512 512 512 512
324 144 144 144 64 64 64 64
441 28 126 567 8 36 162 729
我针对 order=4 测试了它，它应该可以工作。
这是 Matlab 中的矢量化方法。它应该很快，但内存效率不高，因为它会生成所有列索引的笛卡尔元组，然后只保留那些非递减的。
x = [2 2 3 2 8 8; 5 1 7 9 4 4; 4 1 2 7 2 9]; % data
P = 2; % product order
ind = cell(1,P);
[ind{end:-1:1}] = ndgrid(1:size(x,2)); % Cartesian power of column indices with order P
ind = reshape(cat(P+1, ind{:}), [], P); % 2D array where each Cartesian tuple is a row
ind = ind(all(diff(ind, [], 2)>=0, 2), :); % keep only non-decreasing rows
result = prod(reshape(x(:,ind.'), size(x,1), P, []), 2); % apply index into data. This
% creates an intermediate 3D array. Compute products
result = permute(result, [1 3 2]); % convert to 2D array
问题
解释起来有点困难，但我会尽力而为。我知道找到替换组合数的方程式。假设我有 6 个向量：A、B、C、D、E、F。如果我想找到这 6 个变量的所有可能的三次乘积，它将是 (6+3-1)!/3!(6- 1）！ = 56 种组合（见结尾）。同样，如果我想要每个二次积，它是 21。对于线性，当然是 6（只是每个变量本身）。我想计算所有 6+21+56 = 83 种组合。我正在考虑 3 个循环，每个内部循环都从其外部循环开始迭代，如
for i1=1:6
X(:,?) = X.*X(:,i1)
for i2=i1:6
X(:,?) = X.*X(:,i2)
for i3=i2:6
X(:,?) = X.*X(:,i3)
但是左侧存储所有数据的83列矩阵的索引让我感到困惑。如您所见，它们标有问号。
PS：可能也需要对第 5 阶执行此操作，因此它会添加另外 126 和 252 列，总共 461 列。因此，不对三阶进行硬编码的更通用的代码更好。但如果硬编码到第 5 位也没关系，因为我绝对不会超过那个。
MATLAB 或 Python 都可以，因为我可以在两者之间轻松切换。
计算出的二次组合实例
这是我期望的 6 个变量（A 到 F）的二次组合的 21 列示例。在 Excel 中完成。我为每个向量取了 3 个样本。
立方组合列表
这里是我需要计算的56个组合：
一个，一个，一个
A,A,B
A,A,C
A,A,D
A,A,E
A,A,F
A,B,B
A、B、C
A,B,D
A,B,E
A,B,F
A,C,C
A、C、D
A、C、E
A,C,F
A,D,D
A,D,E
A,D,F
A,E,E
A,E,F
A,F,F
B,B,B
B,B,C
B,B,D
B,B,E
B,B,F
B,C,C
B,C,D
B,C,E
B,C,F
B,D,D
B,D,E
B,D,F
B,E,E
B,E,F
B,F,F
C,C,C
C,C,D
C,C,E
C,C,F
C,D,D
C,D,E
C,D,F
C,E,E
C,E,F
C,F,F
D,D,D
D,D,E
D,D,F
D,E,E
D,E,F
D,F,F
E,E,E
E,E,F
E,F,F
F,F,F
您可以通过使用计数器来避免混淆索引：
clear all; close all
% Original matrix
M = [
2 2 3 2 8 8;
5 1 7 9 4 4;
4 1 2 7 2 9
];
% Number of combinations
order = 3;
sizeX = nchoosek(size(M,2)+order-1,order);
% Combinations
imat = ones(sizeX,order);
for c=2:sizeX
imat(c,:) = imat(c-1,:);
for o=order:-1:1
if (imat(c-1,o)<size(M,2))
imat(c,o:end) = imat(c-1,o)+1;
break
end
end
end
% Transpose & display combinations
imat = transpose(imat)
% Computations of products
X = ones(size(M,1),sizeX);
for o=1:order
X = X.*M(:,imat(o,:));
end
% Display result
X
当您执行脚本时，您会得到：
>> test_script
imat =
Columns 1 through 16
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 2 2 2 2 2 3 3 3 3 4
1 2 3 4 5 6 2 3 4 5 6 3 4 5 6 4
Columns 17 through 32
1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2
4 4 5 5 6 2 2 2 2 2 3 3 3 3 4 4
5 6 5 6 6 2 3 4 5 6 3 4 5 6 4 5
Columns 33 through 48
2 2 2 2 3 3 3 3 3 3 3 3 3 3 4 4
4 5 5 6 3 3 3 3 4 4 4 5 5 6 4 4
6 5 6 6 3 4 5 6 4 5 6 5 6 6 4 5
Columns 49 through 56
4 4 4 4 5 5 5 6
4 5 5 6 5 5 6 6
6 5 6 6 5 6 6 6
X =
Columns 1 through 16
8 8 12 8 32 32 8 12 8 32 32 18 12 48 48 8
125 25 175 225 100 100 5 35 45 20 20 245 315 140 140 405
64 16 32 112 32 144 4 8 28 8 36 16 56 16 72 196
Columns 17 through 32
32 32 128 128 128 8 12 8 32 32 18 12 48 48 8 32
180 180 80 80 80 1 7 9 4 4 49 63 28 28 81 36
56 252 16 72 324 1 2 7 2 9 4 14 4 18 49 14
Columns 33 through 48
32 128 128 128 27 18 72 72 12 48 48 192 192 192 8 32
36 16 16 16 343 441 196 196 567 252 252 112 112 112 729 324
63 4 18 81 8 28 8 36 98 28 126 8 36 162 343 98
Columns 49 through 56
32 128 128 128 512 512 512 512
324 144 144 144 64 64 64 64
441 28 126 567 8 36 162 729
我针对 order=4 测试了它，它应该可以工作。
这是 Matlab 中的矢量化方法。它应该很快，但内存效率不高，因为它会生成所有列索引的笛卡尔元组，然后只保留那些非递减的。
x = [2 2 3 2 8 8; 5 1 7 9 4 4; 4 1 2 7 2 9]; % data
P = 2; % product order
ind = cell(1,P);
[ind{end:-1:1}] = ndgrid(1:size(x,2)); % Cartesian power of column indices with order P
ind = reshape(cat(P+1, ind{:}), [], P); % 2D array where each Cartesian tuple is a row
ind = ind(all(diff(ind, [], 2)>=0, 2), :); % keep only non-decreasing rows
result = prod(reshape(x(:,ind.'), size(x,1), P, []), 2); % apply index into data. This
% creates an intermediate 3D array. Compute products
result = permute(result, [1 3 2]); % convert to 2D array