hide long namesshow long names
hide short namesshow short names
Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

NAG Toolbox: nag_tsa_uni_spectrum_lag (g13ca)

Purpose

nag_tsa_uni_spectrum_lag (g13ca) calculates the smoothed sample spectrum of a univariate time series using one of four lag windows – rectangular, Bartlett, Tukey or Parzen window.

Syntax

[c, xg, ng, stats, ifail] = g13ca(nx, mtx, px, iw, mw, ic, c, kc, l, lg, xg, 'nc', nc, 'nxg', nxg)
[c, xg, ng, stats, ifail] = nag_tsa_uni_spectrum_lag(nx, mtx, px, iw, mw, ic, c, kc, l, lg, xg, 'nc', nc, 'nxg', nxg)

Description

The smoothed sample spectrum is defined as
(ω) = 1/(2π)
( M1 )
C0 + 2wkCkcos(ωk)
k = 1
,
f^(ω)=12π (C0+2k=1 M-1wkCkcos(ωk)) ,
where MM is the window width, and is calculated for frequency values
ωi = (2πi)/L,  i = 0,1,,[L / 2],
ωi=2πiL,  i=0,1,,[L/2],
where [][] denotes the integer part.
The autocovariances CkCk may be supplied by you, or constructed from a time series x1,x2,,xnx1,x2,,xn, as
nk
Ck = 1/nxtxt + k,
t = 1
Ck=1nt=1 n-kxtxt+k,
the fast Fourier transform (FFT) being used to carry out the convolution in this formula.
The time series may be mean or trend corrected (by classical least squares), and tapered before calculation of the covariances, the tapering factors being those of the split cosine bell:
(1/2)(1cos(π(t(1/2)) / T)), 1tT
(1/2)(1cos(π(nt + (1/2)) / T)), n + 1Ttn
1, otherwise,
12(1-cos(π (t-12)/T)), 1tT 12(1-cos(π (n-t+12)/T)), n+ 1-Ttn 1, otherwise,
where T = [(np)/2] T=[ np2]  and pp is the tapering proportion.
The smoothing window is defined by
wk = W (k/M) ,  kM1,
wk=W (kM) ,  kM-1,
which for the various windows is defined over 0α < 10α<1 by
rectangular:
W(α) = 1
W(α)=1
Bartlett:
W(α) = 1α
W(α )= 1-α
Tukey:
W(α) = (1/2)(1 + cos(πα))
W(α)=12(1+cos(πα))
Parzen:
W(α) = 1 6α2 + 6α3, 0α(1/2)
W(α) = 2(1α)3, (1/2) < α < 1.
W(α)= 1- 6α2+ 6α3, 0α12 W(α)= 2 (1-α) 3, 12<α< 1.
The sampling distribution of (ω)f^(ω) is approximately that of a scaled χd2χd2 variate, whose degrees of freedom dd is provided by the function, together with multiplying limits mumu, mlml from which approximate 95%95% confidence intervals for the true spectrum f(ω)f(ω) may be constructed as [ ml × (ω) , mu × (ω) ] [ ml × f ^ ( ω ) , mu × f ^ ( ω ) ] . Alternatively, log (ω)f^(ω) may be returned, with additive limits.
The bandwidth bb of the corresponding smoothing window in the frequency domain is also provided. Spectrum estimates separated by (angular) frequencies much greater than bb may be assumed to be independent.

References

Bloomfield P (1976) Fourier Analysis of Time Series: An Introduction Wiley
Jenkins G M and Watts D G (1968) Spectral Analysis and its Applications Holden–Day

Parameters

Compulsory Input Parameters

1:     nx – int64int32nag_int scalar
nn, the length of the time series.
Constraint: nx1nx1.
2:     mtx – int64int32nag_int scalar
If covariances are to be calculated by the function (ic = 0ic=0), mtx must specify whether the data are to be initially mean or trend corrected.
mtx = 0mtx=0
For no correction.
mtx = 1mtx=1
For mean correction.
mtx = 2mtx=2
For trend correction.
Constraint: if ic = 0ic=0, 0mtx20mtx2
If covariances are supplied (ic0ic0), mtx is not used.
3:     px – double scalar
If covariances are to be calculated by the function (ic = 0ic=0), px must specify the proportion of the data (totalled over both ends) to be initially tapered by the split cosine bell taper.
If covariances are supplied (ic0)(ic0), px must specify the proportion of data tapered before the supplied covariances were calculated and after any mean or trend correction. px is required for the calculation of output statistics. A value of 0.00.0 implies no tapering.
Constraint: 0.0px1.00.0px1.0.
4:     iw – int64int32nag_int scalar
The choice of lag window.
iw = 1iw=1
Rectangular.
iw = 2iw=2
Bartlett.
iw = 3iw=3
Tukey.
iw = 4iw=4
Parzen.
Constraint: 1iw41iw4.
5:     mw – int64int32nag_int scalar
MM, the ‘cut-off’ point of the lag window. Windowed covariances at lag MM or greater are zero.
Constraint: 1mwnx1mwnx.
6:     ic – int64int32nag_int scalar
Indicates whether covariances are to be calculated in the function or supplied in the call to the function.
ic = 0ic=0
Covariances are to be calculated.
ic0ic0
Covariances are to be supplied.
7:     c(nc) – double array
nc, the dimension of the array, must satisfy the constraint mwncnxmwncnx.
If ic0ic0, c must contain the nc covariances for lags from 00 to (nc1)(nc-1), otherwise c need not be set.
8:     kc – int64int32nag_int scalar
If ic = 0ic=0, kc must specify the order of the fast Fourier transform (FFT) used to calculate the covariances. kc should be a product of small primes such as 2m2m where mm is the smallest integer such that 2mnx + nc2mnx+nc, provided m20m20.
If ic0ic0, that is covariances are supplied, kc is not used.
Constraint: kcnx + nckcnx+nc. The largest prime factor of kc must not exceed 1919, and the total number of prime factors of kc, counting repetitions, must not exceed 2020. These two restrictions are imposed by the internal FFT algorithm used.
9:     l – int64int32nag_int scalar
LL, the frequency division of the spectral estimates as (2π)/L 2πL . Therefore it is also the order of the FFT used to construct the sample spectrum from the covariances. l should be a product of small primes such as 2m2m where mm is the smallest integer such that 2m2M12m2M-1, provided m20m20.
Constraint: l2 × mw1l2×mw-1. The largest prime factor of l must not exceed 1919, and the total number of prime factors of l, counting repetitions, must not exceed 2020. These two restrictions are imposed by the internal FFT algorithm used.
10:   lg – int64int32nag_int scalar
Indicates whether unlogged or logged spectral estimates and confidence limits are required.
lg = 0lg=0
Unlogged.
lg0lg0
Logged.
11:   xg(nxg) – double array
nxg, the dimension of the array, must satisfy the constraint
  • if ic = 0ic=0, nxgmax (kc,l)nxgmax(kc,l);
  • if ic0ic0, nxglnxgl.
If the covariances are to be calculated, then xg must contain the nx data points. If covariances are supplied, xg may contain any values.

Optional Input Parameters

1:     nc – int64int32nag_int scalar
Default: The dimension of the array c.
The number of covariances to be calculated in the function or supplied in the call to the function.
Constraint: mwncnxmwncnx.
2:     nxg – int64int32nag_int scalar
Default: The dimension of the array xg.
The dimension of the array xg as declared in the (sub)program from which nag_tsa_uni_spectrum_lag (g13ca) is called.
Constraints:
  • if ic = 0ic=0, nxgmax (kc,l)nxgmax(kc,l);
  • if ic0ic0, nxglnxgl.

Input Parameters Omitted from the MATLAB Interface

None.

Output Parameters

1:     c(nc) – double array
If ic = 0ic=0, c will contain the nc calculated covariances.
If ic0ic0, the contents of c will be unchanged.
2:     xg(nxg) – double array
Contains the ng spectral estimates, (ωi)f^(ωi), for i = 0,1,,[L / 2]i=0,1,,[L/2] in xg(1)xg1 to xg(ng)xgng respectively (logged if lg = 1lg=1). The elements xg(i)xgi, for i = ng + 1,,nxgi=ng+1,,nxg contain 0.00.0.
3:     ng – int64int32nag_int scalar
The number of spectral estimates, [L / 2] + 1[L/2]+1, in xg.
4:     stats(44) – double array
Four associated statistics. These are the degrees of freedom in stats(1)stats1, the lower and upper 95%95% confidence limit factors in stats(2)stats2 and stats(3)stats3 respectively (logged if lg = 1lg=1), and the bandwidth in stats(4)stats4.
5:     ifail – int64int32nag_int scalar
ifail = 0ifail=0 unless the function detects an error (see [Error Indicators and Warnings]).

Error Indicators and Warnings

Errors or warnings detected by the function:

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

  ifail = 1ifail=1
On entry,nx < 1nx<1,
ormtx < 0mtx<0 and ic = 0ic=0,
ormtx > 2mtx>2 and ic = 0ic=0,
orpx < 0.0px<0.0,
orpx > 1.0px>1.0,
oriw < 1iw<1,
oriw > 4iw>4,
ormw < 1mw<1,
ormw > nxmw>nx,
ornc < mwnc<mw,
ornc > nxnc>nx,
ornxg < max (kc,l)nxg<max(kc,l) and ic = 0ic=0,
ornxg < lnxg<l and ic0ic0.
  ifail = 2ifail=2
On entry,kc < nx + nckc<nx+nc,
orkc has a prime factor exceeding 1919,
orkc has more than 2020 prime factors, counting repetitions.
This error only occurs when ic = 0ic=0.
  ifail = 3ifail=3
On entry,l < 2 × mw1l<2×mw-1,
orl has a prime factor exceeding 1919,
orl has more than 2020 prime factors, counting repetitions.
W ifail = 4ifail=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 ifail = 5ifail=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.

Accuracy

The FFT is a numerically stable process, and any errors introduced during the computation will normally be insignificant compared with uncertainty in the data.

Further Comments

nag_tsa_uni_spectrum_lag (g13ca) carries out two FFTs of length kc to calculate the covariances and one FFT of length l to calculate the sample spectrum. The time taken by the function for an FFT of length nn is approximately proportional to nlog(n)nlog(n) (but see Section [Further Comments] in (c06pa) for further details).

Example

function nag_tsa_uni_spectrum_lag_example
nx = int64(256);
mtx = int64(1);
px = 0.1;
iw = int64(4);
mw = int64(100);
ic = int64(0);
c = zeros(100, 1);
kc = int64(360);
l = int64(200);
lg = int64(0);
xg = zeros(360,1);
xg(1:256) = [5;
     11;
     16;
     23;
     36;
     58;
     29;
     20;
     10;
     8;
     3;
     0;
     0;
     2;
     11;
     27;
     47;
     63;
     60;
     39;
     28;
     26;
     22;
     11;
     21;
     40;
     78;
     122;
     103;
     73;
     47;
     35;
     11;
     5;
     16;
     34;
     70;
     81;
     111;
     101;
     73;
     40;
     20;
     16;
     5;
     11;
     22;
     40;
     60;
     80.9;
     83.4;
     47.7;
     47.8;
     30.7;
     12.2;
     9.6;
     10.2;
     32.4;
     47.6;
     54;
     62.9;
     85.9;
     61.2;
     45.1;
     36.4;
     20.9;
     11.4;
     37.8;
     69.8;
     106.1;
     100.8;
     81.6;
     66.5;
     34.8;
     30.6;
     7;
     19.8;
     92.5;
     154.4;
     125.9;
     84.8;
     68.1;
     38.5;
     22.8;
     10.2;
     24.1;
     82.9;
     132;
     130.9;
     118.1;
     89.9;
     66.6;
     60;
     46.9;
     41;
     21.3;
     16;
     6.4;
     4.1;
     6.8;
     14.5;
     34;
     45;
     43.1;
     47.5;
     42.2;
     28.1;
     10.1;
     8.1;
     2.5;
     0;
     1.4;
     5;
     12.2;
     13.9;
     35.4;
     45.8;
     41.1;
     30.1;
     23.9;
     15.6;
     6.6;
     4;
     1.8;
     8.5;
     16.6;
     36.3;
     49.6;
     64.2;
     67;
     70.9;
     47.8;
     27.5;
     8.5;
     13.2;
     56.9;
     121.5;
     138.3;
     103.2;
     85.7;
     64.6;
     36.7;
     24.2;
     10.7;
     15;
     40.1;
     61.5;
     98.5;
     124.7;
     96.3;
     66.6;
     64.5;
     54.1;
     39;
     20.6;
     6.7;
     4.3;
     22.7;
     54.8;
     93.8;
     95.8;
     77.2;
     59.1;
     44;
     47;
     30.5;
     16.3;
     7.3;
     37.6;
     74;
     139;
     111.2;
     101.6;
     66.2;
     44.7;
     17;
     11.3;
     12.4;
     3.4;
     6;
     32.3;
     54.3;
     59.7;
     63.7;
     63.5;
     52.2;
     25.4;
     13.1;
     6.8;
     6.3;
     7.1;
     35.6;
     73;
     85.1;
     78;
     64;
     41.8;
     26.2;
     26.7;
     12.1;
     9.5;
     2.7;
     5;
     24.4;
     42;
     63.5;
     53.8;
     62;
     48.5;
     43.9;
     18.6;
     5.7;
     3.6;
     1.4;
     9.6;
     47.4;
     57.1;
     103.9;
     80.6;
     63.6;
     37.6;
     26.1;
     14.2;
     5.8;
     16.7;
     44.3;
     63.9;
     69;
     77.8;
     64.9;
     35.7;
     21.2;
     11.1;
     5.7;
     8.7;
     36.1;
     79.7;
     114.4;
     109.6;
     88.8;
     67.8;
     47.5;
     30.6;
     16.3;
     9.6;
     33.2;
     92.6;
     151.6;
     136.3;
     134.7;
     83.9;
     69.4;
     31.5;
     13.9;
     4.4;
     38];
[cOut, xgOut, ng, stats, ifail] = ...
    nag_tsa_uni_spectrum_lag(nx, mtx, px, iw, mw, ic, c, kc, l, lg, xg);
 ng, stats, ifail
 

ng =

                  101


stats =

    9.0000
    0.4731
    3.3329
    0.1165


ifail =

                    0


function g13ca_example
nx = int64(256);
mtx = int64(1);
px = 0.1;
iw = int64(4);
mw = int64(100);
ic = int64(0);
c = zeros(100, 1);
kc = int64(360);
l = int64(200);
lg = int64(0);
xg = zeros(360,1);
xg(1:256) = [5;
     11;
     16;
     23;
     36;
     58;
     29;
     20;
     10;
     8;
     3;
     0;
     0;
     2;
     11;
     27;
     47;
     63;
     60;
     39;
     28;
     26;
     22;
     11;
     21;
     40;
     78;
     122;
     103;
     73;
     47;
     35;
     11;
     5;
     16;
     34;
     70;
     81;
     111;
     101;
     73;
     40;
     20;
     16;
     5;
     11;
     22;
     40;
     60;
     80.9;
     83.4;
     47.7;
     47.8;
     30.7;
     12.2;
     9.6;
     10.2;
     32.4;
     47.6;
     54;
     62.9;
     85.9;
     61.2;
     45.1;
     36.4;
     20.9;
     11.4;
     37.8;
     69.8;
     106.1;
     100.8;
     81.6;
     66.5;
     34.8;
     30.6;
     7;
     19.8;
     92.5;
     154.4;
     125.9;
     84.8;
     68.1;
     38.5;
     22.8;
     10.2;
     24.1;
     82.9;
     132;
     130.9;
     118.1;
     89.9;
     66.6;
     60;
     46.9;
     41;
     21.3;
     16;
     6.4;
     4.1;
     6.8;
     14.5;
     34;
     45;
     43.1;
     47.5;
     42.2;
     28.1;
     10.1;
     8.1;
     2.5;
     0;
     1.4;
     5;
     12.2;
     13.9;
     35.4;
     45.8;
     41.1;
     30.1;
     23.9;
     15.6;
     6.6;
     4;
     1.8;
     8.5;
     16.6;
     36.3;
     49.6;
     64.2;
     67;
     70.9;
     47.8;
     27.5;
     8.5;
     13.2;
     56.9;
     121.5;
     138.3;
     103.2;
     85.7;
     64.6;
     36.7;
     24.2;
     10.7;
     15;
     40.1;
     61.5;
     98.5;
     124.7;
     96.3;
     66.6;
     64.5;
     54.1;
     39;
     20.6;
     6.7;
     4.3;
     22.7;
     54.8;
     93.8;
     95.8;
     77.2;
     59.1;
     44;
     47;
     30.5;
     16.3;
     7.3;
     37.6;
     74;
     139;
     111.2;
     101.6;
     66.2;
     44.7;
     17;
     11.3;
     12.4;
     3.4;
     6;
     32.3;
     54.3;
     59.7;
     63.7;
     63.5;
     52.2;
     25.4;
     13.1;
     6.8;
     6.3;
     7.1;
     35.6;
     73;
     85.1;
     78;
     64;
     41.8;
     26.2;
     26.7;
     12.1;
     9.5;
     2.7;
     5;
     24.4;
     42;
     63.5;
     53.8;
     62;
     48.5;
     43.9;
     18.6;
     5.7;
     3.6;
     1.4;
     9.6;
     47.4;
     57.1;
     103.9;
     80.6;
     63.6;
     37.6;
     26.1;
     14.2;
     5.8;
     16.7;
     44.3;
     63.9;
     69;
     77.8;
     64.9;
     35.7;
     21.2;
     11.1;
     5.7;
     8.7;
     36.1;
     79.7;
     114.4;
     109.6;
     88.8;
     67.8;
     47.5;
     30.6;
     16.3;
     9.6;
     33.2;
     92.6;
     151.6;
     136.3;
     134.7;
     83.9;
     69.4;
     31.5;
     13.9;
     4.4;
     38];
[cOut, xgOut, ng, stats, ifail] = g13ca(nx, mtx, px, iw, mw, ic, c, kc, l, lg, xg);
 ng, stats, ifail
 

ng =

                  101


stats =

    9.0000
    0.4731
    3.3329
    0.1165


ifail =

                    0



PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

© The Numerical Algorithms Group Ltd, Oxford, UK. 2009–2013