NAG Toolbox: nag_tsa_uni_spectrum_daniell (g13cb)

Purpose

Syntax

Description

References

Parameters

Compulsory Input Parameters

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

$n$, the length of the time series.

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

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

Whether the data are to be initially mean or trend corrected.

$\mathbf{mtx}=0$

For no correction.

$\mathbf{mtx}=1$

For mean correction.

$\mathbf{mtx}=2$

For trend correction.

Constraint: $0\le \mathbf{mtx}\le 2$.

3:     $\mathrm{px}$ – double scalar

The proportion of the data (totalled over both ends) to be initially tapered by the split cosine bell taper. (A value of $0.0$ implies no tapering.)

Constraint: $0.0\le \mathbf{px}\le 1.0$.

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

The value of $M$ which determines the frequency width of the smoothing window as $2\pi /M$. A value of $n$ implies no smoothing is to be carried out.

Constraint: $1\le \mathbf{mw}\le \mathbf{nx}$.

5:     $\mathrm{pw}$ – double scalar

$p$, the shape parameter of the trapezium frequency window.

A value of $0.0$ gives a triangular window, and a value of $1.0$ a rectangular window.

If $\mathbf{mw}=\mathbf{nx}$ (i.e., no smoothing is carried out), pw is not used.

Constraint: $0.0\le \mathbf{pw}\le 1.0$.

6:     $\mathrm{l}$ – int64int32nag_int scalar

$L$, the frequency division of smoothed spectral estimates as $2\pi /L$.

Constraints:

• $\mathbf{l}\ge 1$;
• l must be a factor of kc.

7:     $\lg$ – int64int32nag_int scalar

Indicates whether unlogged or logged spectral estimates and confidence limits are required.

$\mathbf{lg}=0$

For unlogged.

$\mathbf{lg}\ne 0$

For logged.

8:     $\mathrm{xg}\left(\mathbf{kc}\right)$ – double array

The $n$ data points.

Optional Input Parameters

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

Default: the dimension of the array xg.

$K$, the order of the fast Fourier transform (FFT) used to calculate the spectral estimates. kc should be a multiple of small primes such as $2^m$ where $m$ is the smallest integer such that $2^m\ge 2n$, provided $m\le 20$.

Constraints:

• $\mathbf{kc}\ge 2\times \mathbf{nx}$;
• kc must be a multiple of l. The largest prime factor of kc must not exceed $19$, and the total number of prime factors of kc, counting repetitions, must not exceed $20$. These two restrictions are imposed by the internal FFT algorithm used.

Output Parameters

1:     $\mathrm{xg}\left(\mathbf{kc}\right)$ – double array

Contains the ng spectral estimates $\hat{f}\left({\omega }_{\mathit{i}}\right)$, for $\mathit{i}=0,1,\dots ,\left[L/2\right]$, in $\mathbf{xg}\left(1\right)$ to $\mathbf{xg}\left(\mathbf{ng}\right)$ (logged if $\mathbf{lg}\ne 0$). The elements $\mathbf{xg}\left(\mathit{i}\right)$, for $\mathit{i}=\mathbf{ng}+1,\dots ,\mathbf{kc}$, contain $0.0$.

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

The number of spectral estimates, $\left[L/2\right]+1$, in xg.

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

Four associated statistics. These are the degrees of freedom in $\mathbf{stats}\left(1\right)$, the lower and upper $95\%$ confidence limit factors in $\mathbf{stats}\left(2\right)$ and $\mathbf{stats}\left(3\right)$ respectively (logged if $\mathbf{lg}\ne 0$), and the bandwidth in $\mathbf{stats}\left(4\right)$.

4:     $\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}=1$

On entry,	$\mathbf{nx}<1$,
or	$\mathbf{mtx}<0$,
or	$\mathbf{mtx}>2$,
or	$\mathbf{px}<0.0$,
or	$\mathbf{px}>1.0$,
or	$\mathbf{mw}<1$,
or	$\mathbf{mw}>\mathbf{nx}$,
or	$\mathbf{pw}<0.0$ and $\mathbf{mw}\ne \mathbf{nx}$,
or	$\mathbf{pw}>1.0$ and $\mathbf{mw}\ne \mathbf{nx}$,
or	$\mathbf{l}<1$.

   $\mathbf{ifail}=2$

On entry,	$\mathbf{kc}<2\times \mathbf{nx}$,
or	kc is not a multiple of l,
or	kc has a prime factor exceeding $19$,
or	kc has more than $20$ prime factors, counting repetitions.

   $\mathbf{ifail}=3$

This indicates that a serious error has occurred. Check all array subscripts and function argument lists in calls to nag_tsa_uni_spectrum_daniell (g13cb). Seek expert help.

W  $\mathbf{ifail}=4$

One or more spectral estimates are negative. Unlogged spectral estimates are returned in xg, and the degrees of freedom, unlogged confidence limit factors and bandwidth in stats.

W  $\mathbf{ifail}=5$

The calculation of confidence limit factors has failed. This error will not normally occur. Spectral estimates (logged if requested) are returned in xg, and degrees of freedom and bandwidth in stats.

   $\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: g13cb_example

function g13cb_example


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

% Smoothing parameters
nx  = int64(131);
mtx = int64(1);
px  = 0.2;
mw  = int64(30);
pw = 0.5;
l  = int64(100);
lg = int64(1);

% Series
xg = zeros(400, 1);
xg(1:131) = ...
  [11.500;  9.890;  8.728;  8.400;  8.230;  8.365;  8.383;  8.243;  8.080;
    8.244;  8.490;  8.867;  9.469;  9.786; 10.100; 10.714; 11.320; 11.900;
   12.390; 12.095; 11.800; 12.400; 11.833; 12.200; 12.242; 11.687; 10.883;
   10.138;  8.952;  8.443;  8.231;  8.067;  7.871;  7.962;  8.217;  8.689;
    8.989;  9.450;  9.883; 10.150; 10.787; 11.000; 11.133; 11.100; 11.800;
   12.250; 11.350; 11.575; 11.800; 11.100; 10.300;  9.725;  9.025;  8.048;
    7.294;  7.070;  6.933;  7.208;  7.617;  7.867;  8.309;  8.640;  9.179;
    9.570; 10.063; 10.803; 11.547; 11.550; 11.800; 12.200; 12.400; 12.367;
   12.350; 12.400; 12.270; 12.300; 11.800; 10.794;  9.675;  8.900;  8.208;
    8.087;  7.763;  7.917;  8.030;  8.212;  8.669;  9.175;  9.683; 10.290;
   10.400; 10.850; 11.700; 11.900; 12.500; 12.500; 12.800; 12.950; 13.050;
   12.800; 12.800; 12.800; 12.600; 11.917; 10.805;  9.240;  8.777;  8.683;
    8.649;  8.547;  8.625;  8.750;  9.110;  9.392;  9.787; 10.340; 10.500;
   11.233; 12.033; 12.200; 12.300; 12.600; 12.800; 12.650; 12.733; 12.700;
   12.259; 11.817; 10.767;  9.825;  9.150];

% Calculate smooth spectrum
[xg, ng, stats, ifail] = ...
  g13cb( ...
	 nx, mtx, px, mw, pw, l, lg, xg);

% Display results
fprintf('Frequency width of smoothing window = 1/%d\n',mw);
fprintf('Degrees of freedom = %4.1f      Bandwidth = %7.4f\n\n', ...
	stats(1), stats(4));
fprintf('95 percent confidence limits -  Lower = %7.4f  Upper = %7.4f\n', ...
	stats(2:3));
fprintf('\n      Spectrum       Spectrum       Spectrum       Spectrum\n');
fprintf('      estimate       estimate       estimate       estimate\n');
result = [double([1:ng]); xg(1:ng)'];
for j = 1:4:ng
  fprintf('%4d%10.4f', result(:,j:min(j+3,ng)));
  fprintf('\n');
end

g13cb example results

Frequency width of smoothing window = 1/30
Degrees of freedom =  7.0      Bandwidth =  0.1767

95 percent confidence limits -  Lower = -0.8275  Upper =  1.4213

      Spectrum       Spectrum       Spectrum       Spectrum
      estimate       estimate       estimate       estimate
   1   -0.1776   2   -0.4561   3   -0.1784   4    1.9042
   5    2.1094   6    1.7061   7   -0.7659   8   -1.4734
   9   -1.5939  10   -2.1157  11   -2.9151  12   -2.7055
  13   -2.8200  14   -3.4077  15   -3.8813  16   -3.6607
  17   -4.0601  18   -4.4756  19   -4.2700  20   -4.3092
  21   -4.5711  22   -4.8111  23   -4.5658  24   -4.7285
  25   -5.4386  26   -5.5081  27   -5.2325  28   -5.0262
  29   -4.4539  30   -4.4764  31   -4.9152  32   -5.8492
  33   -5.5872  34   -4.9804  35   -4.8904  36   -5.2666
  37   -5.7643  38   -5.8620  39   -5.5011  40   -5.7129
  41   -6.3894  42   -6.4027  43   -6.1352  44   -6.5766
  45   -7.3676  46   -7.1405  47   -6.1674  48   -5.8600
  49   -6.1036  50   -6.2673  51   -6.4321