NAG Toolbox: nag_stat_moving_average (g01wa)

Purpose

Syntax

Description

References

Parameters

Compulsory Input Parameters

1:     $\mathrm{m}$ – int64int32nag_int scalar

$m$, the length of the rolling window.

If $\mathbf{pn}\ne 0$, m must be unchanged since the last call to nag_stat_moving_average (g01wa).

Constraint: $\mathbf{m}\ge 1$.

2:     $\mathrm{x}\left(\mathbf{nb}\right)$ – double array

The current block of observations, corresponding to $x_{\mathit{i}}$, for $\mathit{i}=k+1,\dots ,k+b$, where $k$ is the number of observations processed so far and $b$ is the size of the current block of data.

Optional Input Parameters

1:     $\mathrm{nb}$ – int64int32nag_int scalar

Default: the dimension of the array x.

$b$, the number of observations in the current block of data. The size of the block of data supplied in x (and when $\mathbf{iwt}=1$, wt) can vary; therefore nb can change between calls to nag_stat_moving_average (g01wa).

Constraints:

• $\mathbf{nb}\ge 0$;
• if rcomm is not provided, $\mathbf{nb}\ge \mathbf{m}$.

2:     $\mathrm{iwt}$ – int64int32nag_int scalar

Default: $0$

The type of weighting to use.

$\mathbf{iwt}=0$

No weights are used.

$\mathbf{iwt}=1$

Each observation has its own weight.

$\mathbf{iwt}=2$

Each position in the window has its own weight.

$\mathbf{iwt}=3$

Each position in the window has a weight equal to its position number.

If $\mathbf{pn}\ne 0$, iwt must be unchanged since the last call to nag_stat_moving_average (g01wa).

Constraint: $\mathbf{iwt}=0$, $1$, $2$ or $3$.

3:     $\mathrm{wt}\left(:\right)$ – double array

The dimension of the array wt must be at least $\mathbf{nb}$ if $\mathbf{iwt}=1$ and at least $\mathbf{m}$ if $\mathbf{iwt}=2$

The user-supplied weights.

If $\mathbf{iwt}=1$, $\mathbf{wt}\left(\mathit{i}\right)={\nu }_{\mathit{i}+k}$, for $\mathit{i}=1,2,\dots ,b$.

If $\mathbf{iwt}=2$, $\mathbf{wt}\left(\mathit{j}\right)=w_{\mathit{j}}$, for $\mathit{j}=1,2,\dots ,m$.

Otherwise, wt is not referenced.

Constraints:

• if $\mathbf{iwt}=1$, $\mathbf{wt}\left(\mathit{i}\right)\ge 0$, for $\mathit{i}=1,2,\dots ,\mathbf{nb}$;
• if $\mathbf{iwt}=2$, $\mathbf{wt}\left(1\right)\ne 0$ and ${\sum }_{\mathit{j}=1}^m\mathbf{wt}\left(\mathit{j}\right)>0$;
• if $\mathbf{iwt}=2$ and $\mathbf{lrsd}\ne 0$, $\mathbf{wt}\left(\mathit{j}\right)\ge 0$, for $\mathit{j}=1,2,\dots ,\mathbf{m}$.

4:     $\mathrm{pn}$ – int64int32nag_int scalar

Default: $0$

$k$, the number of observations processed so far. On the first call to nag_stat_moving_average (g01wa), or when starting to summarise a new dataset, pn must be set to $0$.

If $\mathbf{pn}\ne 0$, it must be the same value as returned by the last call to nag_stat_moving_average (g01wa).

Constraint: $\mathbf{pn}\ge 0$.

5:     $\mathrm{wantsd}$ – logical scalar

Default: $\mathit{false}$

If the standard deviations are required then wantsd should be set to true.

6:     $\mathrm{rcomm}\left(2\mathbf{m}+20\right)$ – double array

Communication array, used to store information between calls to nag_stat_moving_average (g01wa). If $\mathbf{lrcomm}=0$, rcomm is not referenced and all the data must be supplied in one go.

Output Parameters

1:     $\mathrm{pn}$ – int64int32nag_int scalar

Default: $0$

$k+b$, the updated number of observations processed so far.

2:     $\mathrm{rmean}\left(\max \left(0,\mathbf{nb}+\min \left(0,\mathbf{pn}-\mathbf{m}+1\right)\right)\right)$ – double array

${\mu }_{\mathit{l}}$, the (weighted) moving averages, for $\mathit{l}=1,2,\dots ,b+\min \left(0,k-m+1\right)$. Therefore, ${\mu }_l$ is the mean of the data in the window that ends on $\mathbf{x}\left(l+m-\min \left(k,m-1\right)-1\right)$.

If, on entry, $\mathbf{pn}\ge \mathbf{m}-1$, i.e., at least one windows worth of data has been previously processed, then $\mathbf{rmean}\left(l\right)$ is the summary corresponding to the window that ends on $\mathbf{x}\left(l\right)$. On the other hand, if, on entry, $\mathbf{pn}=0$, i.e., no data has been previously processed, then $\mathbf{rmean}\left(l\right)$ is the summary corresponding to the window that ends on $\mathbf{x}\left(\mathbf{m}+l-1\right)$ (or, equivalently, starts on $\mathbf{x}\left(l\right)$).

3:     $\mathrm{rsd}\left(\mathbf{lrsd}\right)$ – double array

If sdp was set to true, then ${\sigma }_l$, the (weighted) standard deviation. The ordering of rsd is the same as the ordering of rmean.

4:     $\mathrm{rcomm}\left(2\mathbf{m}+20\right)$ – double array

Communication array, used to store information between calls to nag_stat_moving_average (g01wa).

5:     $\mathrm{ifail}$ – int64int32nag_int scalar

$\mathbf{ifail}=\mathbf{0}$ unless the function detects an error (see Error Indicators and Warnings).

Error Indicators and Warnings

Cases prefixed with W are classified as warnings and do not generate an error of type NAG:error_n. See nag_issue_warnings.

   $\mathbf{ifail}=11$

Constraint: $\mathbf{m}\ge 1$.

   $\mathbf{ifail}=12$

Constraint: if $\mathbf{pn}>0$, m must be unchanged since previous call.

   $\mathbf{ifail}=21$

Constraint: $\mathbf{nb}\ge 0$.

   $\mathbf{ifail}=22$

Constraint: if $\mathbf{lrcomm}=0$, $\mathbf{nb}\ge \mathbf{m}$.

   $\mathbf{ifail}=41$

Constraint: $\mathbf{iwt}=0$, $1$, $2$ or $3$.

   $\mathbf{ifail}=42$

Constraint: if $\mathbf{pn}>0$, iwt must be unchanged since previous call.

   $\mathbf{ifail}=51$

Constraint: $\mathbf{wt}\left(i\right)\ge 0$.

   $\mathbf{ifail}=52$

Constraint: if $\mathbf{iwt}=2$, $\mathbf{wt}\left(1\right)>0$.

W  $\mathbf{ifail}=53$

On entry, at least one window had all zero weights.

W  $\mathbf{ifail}=54$

On entry, unable to calculate at least one standard deviation due to the weights supplied.

   $\mathbf{ifail}=55$

Constraint: if $\mathbf{iwt}=2$, the sum of the weights $>0$.

   $\mathbf{ifail}=61$

Constraint: $\mathbf{pn}\ge 0$.

   $\mathbf{ifail}=62$

Constraint: if $\mathbf{pn}>0$, pn must be unchanged since previous call.

   $\mathbf{ifail}=91$

Constraint: $\mathbf{lrsd}=0$ or .

   $\mathbf{ifail}=101$

rcomm has been corrupted between calls.

   $\mathbf{ifail}=111$

lrcomm is too small.

   $\mathbf{ifail}=-99$

An unexpected error has been triggered by this routine. Please contact NAG.

   $\mathbf{ifail}=-399$

Your licence key may have expired or may not have been installed correctly.

   $\mathbf{ifail}=-999$

Dynamic memory allocation failed.

Accuracy

Further Comments

Example

Open in the MATLAB editor: g01wa_example

function g01wa_example


fprintf('g01wa example results\n\n');

m   = int64([5,10,15]);
iwt = int64([0, 0, 2]);
wt{1}  = [1:double(m(1))];
wt{2}  = [1:double(m(2))];
wt{3}  = [ -3.0; -6.0; -5.0;  3.0; 21.0; 46.0; 67.0; 74.0;
           67.0; 46.0; 21.0;  3.0; -5.0; -6.0; -3.0];
data{1} = [-2170; -1770; -1660; -1360; -1100];
data{2} = [ -950;  -640;  -370;  -140;  -250; 
	    -510;  -620;  -730;  -880; -1130];
data{3} = [-1200;  -830;  -330;  -190;   210;
             170;   440;   440;   780;   880;
	    1220;  1260;  1140;   850;   640];
n = 30;

fprintf('\n Interval          Mean\n');
fprintf('------------------------\n');

% Loop over window lengths
for im = 1:numel(m);
  pn = int64(0);
  x{im} = 1820 + m(im)/2 + [1:n-m(im)+1];
  for b=1:3
    % Calculate the number of summaries we can produce
    nb = numel(data{b});
    nsummaries = max(0, nb + min(0, pn-m(im)+1));
    
    % Calculate summary statistics for this block of data
    if b==1
      [pn, rmean, rsd, rcomm, ifail] = ...
      g01wa( ...
	     m(im), data{b}, 'iwt', iwt(im), 'wt', wt{im}, 'wantsd', im<3);
    else
      [pn, rmean, rsd, rcomm, ifail] =  ...
      g01wa( ...
	     m(im), data{b}, 'iwt', iwt(im), 'wt', wt{im}, 'pn', pn, ...
             'rcomm', rcomm, 'wantsd', im<3);
    end

    % Number of results printed so far
    offset = max(0, pn-nb-m(im)+1);

    % Display the results for this block of data
    for i=1:nsummaries
      if im==3
	fprintf(' [%2d,%3d]%14.1f\n', i+offset, i+m(im)-1+offset, rmean(i));
      else
	y{im}(i+offset) = rmean(i);
	z{im}(i+offset) = rsd(i);
      end
    end
  end
end

fprintf('\nTotal number of observations : %d\n', pn);
fprintf('Length of window             : %d\n', m(3));

g01wa_plot(x,y,z,data);



function g01wa_plot(x, y, z, data)

  d = [data{1}; data{2}; data{3}];
  fig1 = figure;
  plot(x{1},y{1},x{2},y{2},[1821:1850],d,'*');
  title({'Raw data and mean from a rolling window of', ...
          'changes in rate of Earth''s rotation (μs)'});
  xlabel('Year');
  ylabel('Change in day length');
  legend('m = 5', 'm = 10', 'Location', 'SouthEast');
  set(gca,'xtick',[1820:5:1850])
  set(gca,'ytick',[-2500:500:1500])
  set(gca, 'PlotBoxAspectRatio', [1,1/2,1])
  fig2 = figure;
  plot(x{1},z{1},x{2},z{2});
  title({'Standard deviation from a rolling window of', ...
	 'changes in rate of Earth''s rotation (μs)'});
  ylabel('Standard deviation');
  xlabel('Year');
  set(gca,'ytick',[100:100:700])
  set(gca, 'PlotBoxAspectRatio', [1,1/3,1])

g01wa example results


 Interval          Mean
------------------------
 [ 1, 15]        -427.6
 [ 2, 16]        -332.5
 [ 3, 17]        -337.1
 [ 4, 18]        -438.2
 [ 5, 19]        -604.4
 [ 6, 20]        -789.4
 [ 7, 21]        -935.4
 [ 8, 22]        -990.6
 [ 9, 23]        -927.1
 [10, 24]        -752.1
 [11, 25]        -501.3
 [12, 26]        -227.2
 [13, 27]          23.2
 [14, 28]         236.2
 [15, 29]         422.4
 [16, 30]         604.2

Total number of observations : 30
Length of window             : 15