三重加权和答案

【问题标题】：Triple weighted sum三重加权和
【发布时间】：2016-12-30 18:14:09
【问题描述】：

我试图对某个加权和进行矢量化，但不知道该怎么做。我在下面创建了一个简单的最小工作示例。我想解决方案涉及 bsxfun 或 reshape 和 kronecker 产品，但我仍然没有设法让它工作。

rng(1);
N = 200;
T1 = 5;
T2 = 7;
T3 = 10;


A = rand(N,T1,T2,T3);
w1 = rand(T1,1);
w2 = rand(T2,1);
w3 = rand(T3,1);

B = zeros(N,1);

for i = 1:N
 for j1=1:T1
  for j2=1:T2
   for j3=1:T3
    B(i) = B(i) + w1(j1) * w2(j2) * w3(j3) * A(i,j1,j2,j3);
   end
  end
 end
end

A = B;

对于二维情况，有一个聪明的答案here。

【问题讨论】：

你需要概括吗？因为如果是这样，我会将您的 N、T1、T2、T3 改写为数组。
其实我只是想要三维外壳。但概括可能对其他人有用:)
下面的概括:)
认为您可能对基准测试感兴趣。对于在 Octave 中执行的 1000 次迭代，循环版本用了 20 分钟，我的版本用 repmat 用了 1.66 秒，没有用了 0.92 秒，而 Suever 用了 0.77 秒。 :)
@TasosPapastylianou 太棒了！很有用！

标签： matlab vectorization weighted-average

【解决方案1】：

如果我们还是要走具有功能的路线，并且偏爱性能而不是优雅/简洁，那么考虑一下：

function B = weightReduce(A, varargin)

    B = A;
    for i = length(varargin):-1:1
        N = length(varargin{i});
        B = reshape(B, [], N) * varargin{i};
    end

end

这是我看到的性能对比：

tic;
for i = 1:10000
    W = createWeights(w1,w2,w3);
    B = reshape(A, size(A,1), [])*W(:);
end
toc
Elapsed time is 0.920821 seconds.
tic; 
for i = 1:10000
    B2 = weightReduce(A, w1, w2, w3); 
end
toc
Elapsed time is 0.484470 seconds.

【讨论】：

【解决方案2】：

这就是它背后的逻辑：

ww1 = repmat (permute (w1, [4, 1, 2, 3]), [N, 1,  T2, T3]);
ww2 = repmat (permute (w2, [3, 4, 1, 2]), [N, T1, 1,  T3]);
ww3 = repmat (permute (w3, [2, 3, 4, 1]), [N, T1, T2, 1 ]);

B = ww1 .* ww2 .* ww3 .* A;
B = sum (B(:,:), 2)

您可以通过首先在适当的维度中创建w1、w2 和w3 来避免permute。你也可以使用bsxfun 而不是repmat 来获得额外的性能，我只是在这里展示逻辑，repmat 更容易理解。

编辑：任意输入维度的通用版本：

Dims   = {N, T1, T2, T3};  % add T4, T5, T6, etc as appropriate
Params = cell (1, length (Dims));

Params{1} = rand (Dims{:});
for n = 2 : length (Dims)
  DimSubscripts = ones (1, length (Dims));  DimSubscripts(n) = Dims{n};
  RepSubscripts = [Dims{:}];  RepSubscripts(n) = 1;
  Params{n} = repmat (rand (DimSubscripts), RepSubscripts);
end

B = times (Params{:});
B = sum (B(:,:), 2)

【讨论】：

【解决方案3】：

您可以使用额外的乘法来修改上一个答案中的w1 * w2' 网格，然后也乘以w3。然后，您可以再次使用矩阵乘法与A 的“扁平化”版本相乘。

W = reshape(w1 * w2.', [], 1) * w3.';
B = reshape(A, size(A, 1), []) * W(:);

您可以将权重的创建封装到它自己的函数中，并将其推广到N 权重。由于这使用递归，N 被限制为您当前的递归限制（默认为 500）。

function W = createWeights(W, varargin)
    if numel(varargin) > 0
        W = createWeights(W(:) * varargin{1}(:).', varargin{2:end});
    end
end

并将其用于：

W = createWeights(w1, w2, w3);
B = reshape(A, size(A, 1), []) * W(:);

更新

使用@CKT 非常好的建议使用kron 的一部分，我们可以稍微修改createWeights。

function W = createWeights(W, varargin)
    if numel(varargin) > 0
        W = createWeights(kron(varargin{1}, W), varargin{2:end});
    end
end

【讨论】：

【解决方案4】：

再一次，除非你创建一些函数来构造 Kronecker 乘积向量，否则你不能很好地对 N-D 进行概括，但是怎么样

  A = reshape(A, N, []) * kron(w3, kron(w2, w1));

【讨论】：