An efficient way to calculate the "product" of discrete convolution

Question

I am looking for an elegant way to calculate the "product" of a discrete convolution, not a sum.

Here is the discrete convolution formula :

In this case, we can use: conv(x,y)

Now I would like to implement these operations

Of course, I can use a loop, but I'm looking for a trick to linearize this operation.

Example:

f = [2 4 3 9 7 1]
g = [3 2 1]

dist = length(g)-1;

for ii = 1:length(f)-dist
    x(ii) = prod(f(ii:ii+dist).*g)
end

x =

144    648   1134    378

Answer 1

Another solution, partially inspired by Dev-iL's answer to a relatively same question

exp(sum(log(g))+conv(log(f),[1 1 1],'valid'))

or

exp(sum(log(g))+movsum(log(f),numel(g),'Endpoints', 'discard'))

as exp(sum(log(x))) = prod(x)

f g.

:

:

f= rand(1,1000000);
g= rand(1,100);

disp('----------EXP(LOG)------')
tic
    exp(sum(log(g))+conv(log(f),ones(1,numel(g))));
toc

disp('----------BSXFUN------')
tic
    ind = bsxfun(@plus, 0:numel(f)-numel(g), (1:numel(g)).');
    x = prod(bsxfun(@times, f(ind), g(:)),1);
toc
disp('----------HANKEL------')
tic
    prod(g)*prod(hankel(f(1:numel(g)), f(numel(g):end)));
toc

disp('----------CUMPROD-----')
tic
    pf = cumprod(f);
    x = prod(g)*pf(numel(g):end)./[1 pf(1:(end-numel(g)))];
toc

:

----------EXP(LOG)------%rahnema1
Elapsed time is 0.211445 seconds.
----------BSXFUN--------%Luis Mendo
Elapsed time is 1.94182 seconds.
----------HANKEL--------%gnovice
Elapsed time is 1.46593 seconds.
----------CUMPROD-------%gnovice
Elapsed time is 0.00748992 seconds.

Answer 2

cumprod : ( )

>> pf = cumprod(f);
>> x = prod(g).*pf(numel(g):end)./[1 pf(1:(end-numel(g)))]

x =

         144         648        1134         378

f, cumprod. 3 , numel(g) f. g.

.. f ( ), / . f , :

c = ...set a scaling factor...
pf = cumprod(f./c);
x = prod(c.*g).*pf(numel(g):end)./[1 pf(1:(end-numel(g)))];

c - mean(abs(f)) max(abs(f)), f , (.. 1). .

---

hankel : ( , )

>> x = prod(g).*prod(hankel(f(1:numel(g)), f(numel(g):end)))

x =

         144         648        1134         378

hankel , numel(g) . , g . f / g .

---

:

6 ( , 2 rahnema (conv(log) movsum(log)), bsxfun Luis cumprod hankel) f = rand(1,1000000); g = rand(1,100); 40 . ( Windows 7 x64, 16 , MATLAB R2016b):

solution    | avg. time (s)
------------+---------------
loop        |       1.10671
conv(log)   |       0.04984
movsum(log) |       0.03736
bsxfun      |       1.20472
cumprod     |       0.01469
hankel      |       1.17704

Answer 3

:

ind = bsxfun(@plus, 0:numel(f)-numel(g), (1:numel(g)).');
x = prod(bsxfun(@times, f(ind), g(:)), 1);

This works as follows. indrepresents a sliding window of indices, each of which corresponds to an offset. For example, if it ghas a size 3, then the matrix indwill be

1 2 3 4 ...
2 3 4 5 ...
3 4 5 6 ...

Used for indexing in f. The result is multiplied (with translation) by gas a column. Finally, the product of the elements of each column is calculated.