bugGNU Octave - Bugs: bug #65245, movfun/movmean: scope for speed up...

 
 

bug #65245: movfun/movmean: scope for speed up with conv

Submitter:  A.R. Burgers <arb>
Submitted:  Sat 03 Feb 2024 12:18:16 PM UTC
   
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Performance
Status:  Confirmed Assigned to:  None
Originator Name:  Open/Closed:  * Open
Release:  * dev Operating System:  * Any
Fixed Release:  None Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup
   

Fri 09 Feb 2024 02:31:07 PM UTC, comment #2: 

No, not so easy to integrate this in the current framework, just wanted to have a record of the idea for the time being.

movstd and movvar also yield to this approach


np = 100;
nw = 13;
window = ones(1, nw);
a1 = rand(1, np);
a2 = a1 .^2;
a0 = ones(1, np);
count  = conv(a0, window, 'same');
count1 = count - 1;

fc = @(x) conv(x, window, 'same');

conv_sum  = fc(a1);
conv_mean = fc(a1) ./ count;
conv_var0 = (fc(a2) - conv_mean .* (2 * fc(a1)  - nw * conv_mean)) ./ count1;
conv_var1 = (fc(a2) - conv_mean .* (2 * fc(a1)  - nw * conv_mean)) ./ count;
conv_std0 = sqrt(conv_var0);
conv_std1 = sqrt(conv_var1);

mov_sum  = movsum(a1, nw);
mov_mean = movmean(a1, nw);
mov_var0 = movvar(a1, nw);
mov_var1 = movvar(a1, nw, 1);
mov_std0 = movstd(a1, nw);
mov_std1 = movstd(a1, nw, 1);

assert(all(abs(mov_sum - conv_sum)) < 1.e-10);
assert(all(abs(mov_mean - conv_mean)) < 1.e-10);
assert(all(abs(mov_var0 - conv_var0)) < 1.e-10);
assert(all(abs(mov_var1 - conv_var1)) < 1.e-10);
assert(all(abs(mov_std0 - conv_std0)) < 1.e-10);
assert(all(abs(mov_std1 - conv_std1)) < 1.e-10);


A.R. Burgers <arb>
Sat 03 Feb 2024 06:31:33 PM UTC, comment #1: 

Yes, this is true.  It might be a good addition, but there is also a bigger strategic issue that is worth thinking about.  The current architecture is that the movXXX functions are simple in that they do just a little bit of input validation and then pass on the function to movfun().  This structure makes maintenance easy and any improvements to movfun() are automatically reflected in the derivative movXXX functions.  movmean() is a special case in that there is a convenient kernel that can be used with conv() (which is written in C++) to accelerate that function.  But other functions, such as movmad(), do not seem to have any easy mapping on to conv().  Perhaps it would be better to re-write movfun() in C++ where loops would not be a problem.

Rik <rik5>
Group administrator
Sat 03 Feb 2024 12:18:16 PM UTC, original submission:  

see also bug #55389

with a a vector, window size wlen odd, conv can give much faster results because of the linearity of the operator in the movsum/movmean cases.

The following two snippets show how conv can calculate the same result for movsum and movmean


result = movsum(a, wlen);
result = conv(a, ones(1, wlen), 'same');



result = movmean(a, wlen)
count = conv(ones(size(a)), ones(1, wlen), 'same');
result = conv(a, ones(1, wlen), 'same') ./ count;


Speed-ups can be a factor ~30


movsum time 0.321551
movsum with conv time 0.0108559
movmean time 0.520523
movmean with conv time 0.0155551


Timings obtained with this script


np = 10000; ns = 13;
nrep = 100;
a = rand(1, np);
a1 = ones(1, np);

tic;
for i = 1 : nrep
  as = movsum(a, ns);
end
fprintf('movsum time %g\n', toc);

tic;
for i = 1 : nrep
  as_c = conv(a, ones(1, ns), 'same');
end
fprintf('movsum with conv time %g\n', toc);

assert(all(abs(as - as_c)) < 1.e-10);

tic;
for i = 1 : nrep
  as = movmean(a, ns);
end
fprintf('movmean time %g\n', toc);

tic;
for i = 1 : nrep
  as_c = conv(a, ones(1, ns), 'same') ./ conv(a1, ones(1, ns), 'same');
end
fprintf('movmean with conv time %g\n', toc);


A.R. Burgers <arb>

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:
   
   
Comment:
   

No files currently attached

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by rik5 (Posted a comment)
  • -email is unavailable- added by arb (Submitted the item)
  •  

    There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

    Only group members can vote.

     

    Follows 1 latest change.

    Date Changed by Updated Field Previous Value => Replaced by
    2024-02-03 rik5 StatusNone Confirmed

    Back to the top

    Powered by Savane 3.12.
    Corresponding source code