Chapter 8: Time Series Analysis and Forecasting

SWED

TextonSpedometerHPClogo.gif

Estimations of the nonnormalized spectral density of a stationary time series based on specified periodogram weights given the time series data.

Required Arguments

X — Vector of length NOBS containing the stationary time series.   (Input)

F — Vector of length NF containing the frequencies at which to evaluate the spectral density estimate.   (Input)
The units of F correspond to the scale specified by IFSCAL. The elements of F must be in the range (π/TINT, π/TINT), inclusive, for IFSCAL = 0 and (1/(2 * TINT),
1/(2 * TINT)), inclusive, for IFSCAL = 1.

WT — Vector of length NWT containing the weights used to smooth the periodogram.   (Input)
The actual weights are the values in WT normalized to sum to 1 with the current periodogram ordinate taking the middle weight for NWT odd or the weight to the right of the middle for NWT even.

PM — (N/2 + 1) by 5 matrix that contains a summarization of the periodogram analysis.   (Output)
For k = 0, 1, …, N/2, the (k + 1)-st element of the j-th column of PM is defined as

Col.

Description

1

Frequency, ωk where ωk = 2πk/N for IFSCAL = 0 or ωk = k/N for IFSCAL = 1.

2

Period, pk where pk = 2π/ωk for IFSCAL = 0 and pk = 1/ ωk for IFSCAL = 1. If
ωk = 0, pk is set to the missing value or NaN (not a number).

3

Periodogram ordinate, Ik).

4

Cosine transformation coefficient, Ak).

5

Sine transformation coefficient, Bk).

SM — NF by 3 matrix containing a summarization of the spectral analysis.   (Output)
The k-th element of the j-th column of SM is defined as

Col.

Description

1

Frequency, F(k).

2

Period, pk where pk = 2π/F(k) for IFSCAL = 0 and pk = 1/F(k) for IFSCAL = 1. If F(k) = 0, pk is set to missing.

3

Spectral density estimate at F(k) using the specified weights WT.

            where k = 1, …, NF.

Optional Arguments

NOBS — Number of observations in the stationary time series X.   (Input)
NOBS must be greater than or equal to two.
Default: NOBS = size (X,1).

IPRINT — Printing option.   (Input)
Default: IPRINT = 0.

IPRINT

Action

0

No printing is performed.

1

Prints the periodogram, cosine and sine transforms of the centered and padded time series, and the spectral density estimate based on a specified weight sequence.

XCNTR — Constant used to center the time series X.   (Input)
Default: XCNTR = the arithmetic mean.

NPAD — Number of zeroes used to pad the centered time series.   (Input)
NPAD must be greater than or equal to zero. The length of the centered and padded time series is N = NOBS + NPAD.
Default: NPAD = NOBS –1.

IFSCAL — Option for frequency scale.   (Input)
Default: IFSCAL = 0.

IFSCAL

Action

0

Frequency in radians per unit time.

1

Frequency in cycles per unit time.

NF — Number of frequencies at which to evaluate the spectral density estimate.   (Input)
NF must be greater than zero.
Default: NF = size (F,1).

TINT — Time interval at which the series is sampled.   (Input)
For a discrete parameter process, usually TINT = 1.0. For a continuous parameter process, TINT > 0.0. TINT is used to adjust the spectral density estimate.
Default: TINT = 1.0.

NWT — Number of weights.   (Input)
NWT must be greater than or equal to one.
Default: NWT = size (WT,1).

LDPM — Leading dimension of PM exactly as specified in the dimension statement in the calling program.   (Input)
LDPM must be greater than or equal to N/2 + 1.
Default: LDPM = size (PM,1).

LDSM — Leading dimension of SM exactly as specified in the dimension statement in the calling program.   (Input)
LDSM must be greater than or equal to NF.
Default: LDSM = size (SM,1).

FORTRAN 90 Interface

Generic:                              CALL SWED (X, F, WT, PM, SM [,…])

Specific:                             The specific interface names are S_SWED and D_SWED.

FORTRAN 77 Interface

Single:            CALL SWED (NOBS, X, IPRINT, XCNTR, NPAD, IFSCAL, NF, F, TINT, NWT, WT, PM, LDPM, SM, LDSM)

Double:                              The double precision name is DSWED.

Description

Routine SWED estimates the nonnormalized spectral density function of a stationary time series using a fixed sequence of weights, given a sample of n = NOBS observations {Xt}, for
 t = 1, 2, …, n.

Let

represent the centered and padded data where N = NOBS + NPAD,

and

is determined by

The modified periodogram of

is estimated by

where

and

represent the

cosine and sine transforms, respectively, and K is the scale factor equal to 1/(2πn). Since the periodogram is an even function of the frequency, it is sufficient to estimate the periodogram at the discrete set of nonnegative frequencies

(Here, a means the greatest integer less than or equal to a.) The routine PFFT is used to compute the modified periodogram of

Consider the sequence of m = NWT weights

{wj} for j = m/2, …, (m m/2 1)

where

Σjwj = 1

These weights are fixed in the sense that they do not depend on the frequency ω at which to estimate the nonnormalized spectral density hX(ω). The estimate of the nonnormalized spectral density is computed according to

where

and k(ω) is the integer such that ωk,0 is closest to ω. The weights specified by argument WT may be relative since they are normalized to sum to one in the actual computation of

Usually, m is odd with the weights symmetric about the middle weight w0. If m is even, the weight to the right of the middle is considered w0. Note that periodogram ordinate

is replaced by

and the sum reflects at each end. The nonnormalized spectral density is estimated over the set of frequencies

ω = fi,      i = 1, …, nf

where nf = NF. These frequencies are in the scale specified by the argument IFSCAL but are transformed to the scale of radians per unit time for computational purposes.

The above formula for

assumes the data {Xt} correspond to a realization of a discrete-parameter stationary process observed consecutively in time. In this case, the observations are equally spaced in time with interval Δt = TINT equivalent to one. However, if the data correspond to a realization of a continuous-parameter stationary process recorded at equal time intervals, then the estimate of the nonnormalized spectral density must be adjusted for the effect of aliasing. In general, the estimate of hX(ω) is given by

Note that the frequency ω of the desired spectral estimate is assumed to be input in a form already adjusted for the time interval Δt.

Approximate confidence intervals for h(ω) can be computed using formulas given in the introduction.

Comments

1.         Workspace may be explicitly provided, if desired, by use of S2ED/DS2ED. The reference is:

CALL S2ED (NOBS, X, IPRINT, XCNTR, NPAD, IFSCAL, NF, F, TINT,  NWT, WT, PM, LDPM, SM, LDSM,CX, COEF, WFFTC, CPY)

The additional arguments are as follows:

CX — Complex vector of length N containing the centered and padded time series X.   (Output)

COEF — Complex vector of length N containing the Fourier coefficients of the finite Fourier transform of CX.   (Output)
Note that COEF(k + 1) is the appropriately scaled Fourier coefficient at frequency ωk, k = 0, 1, …, N 1.

WFFTC — Work vector of length 4N + 15.

CPY — Work vector of length 2N.

2.         The centered and padded time series is defined by
CX(j) = X(j) XCNTR        for j = 1, …, NOBS
CX(j) = 0       for j = NOBS + 1, …, N
where N = NOBS + NPAD.

3.         The normalized spectral density estimate is obtained by dividing the nonnormalized spectral density estimate in matrix SM by an estimate of the variance of X.

Example

Consider the Wφlfer Sunspot Data (Anderson 1971, page 660) consisting of the number of sunspots observed each year from 1749 through 1924. The data set for this example consists of the number of sunspots observed from 1770 through 1869. Application of routine SWED to these data produces the following results:

 

      USE IMSL_LIBRARIES

 

      IMPLICIT   NONE

      INTEGER    LDPM, LDRDAT, LDSM, NDRDAT, NF, NOBS, NPAD, NWT

      PARAMETER  (LDRDAT=176, NDRDAT=2, NF=20, NOBS=100, NWT=7, &

                  LDSM=NF, NPAD=NOBS-1, LDPM=(NOBS+NPAD)/2+1)

!

      INTEGER    I, NROW, NVAR

      REAL       F(NF), PI, PM(LDPM,5), RDATA(LDRDAT,NDRDAT), &

                 REAL, SM(LDSM,3), WT(NWT), X(NOBS), IFSCAL, IPRINT, TINT

      CHARACTER  CLABEL(4)*20, FMT*20, RLABEL(1)*4, TITLE*28

      INTRINSIC  FLOAT

!

      EQUIVALENCE (X(1), RDATA(22,2))

!

      DATA WT/1.0, 2.0, 3.0, 4.0, 3.0, 2.0, 1.0/

      DATA IPRINT/0/, IFSCAL/0/, TINT/1.0/

      DATA FMT/'(F9.4, F6.2, F9.4)'/

      DATA RLABEL/'NONE'/

      DATA CLABEL/' ', '%/Frequency', '%/Period', 'Spectral%/Estimates' &

          /

      DATA TITLE/'Results of Spectral Analysis'/

!                                 Initializations

      PI = 2.0*ASIN(1.0)

      DO 10  I=1, NF

         F(I) = PI*FLOAT(I)/FLOAT(NF)

   10 CONTINUE

!                                 Wolfer Sunspot Data for years

!                                 1770 through 1869

      CALL GDATA (2, RDATA, NROW, NVAR)

!                                 Center on arithmetic mean

!                                 Spectral density

      CALL SWED (X, F, WT, PM, SM)

!                                 Print Results

      CALL WRRRL (TITLE, SM, RLABEL, CLABEL, FMT=FMT)

!

      END

Output

 

Results of Spectral Analysis
                    Spectral
Frequency  Period  Estimates
   0.1571   40.00   710.8386
   0.3142   20.00   116.3940
   0.4712   13.33   937.1508
   0.6283   10.00  1209.8268
   0.7854    8.00   538.9236
   0.9425    6.67    84.9561
   1.0996    5.71   128.0791
   1.2566    5.00    55.0304
   1.4137    4.44    40.2022
   1.5708    4.00    46.4240
   1.7279    3.64    21.0053
   1.8850    3.33    12.1449
   2.0420    3.08     8.8654
   2.1991    2.86     7.2589
   2.3562    2.67     6.8078
   2.5133    2.50     3.3873
   2.6704    2.35     3.9504
   2.8274    2.22     5.7418
   2.9845    2.11     4.4652
   3.1416    2.00     4.1216



http://www.vni.com/
PHONE: 713.784.3131
FAX:713.781.9260