Chapter 8: Time Series Analysis and Forecasting

CSWEP

Estimates the nonnormalized cross-spectral density of two stationary time series using a weighted cross periodogram given the spectral densities and cross periodogram.

Required Arguments

N — Number of observations in each of the appropriately centered and padded
time series X and Y.   (Input)
N must be greater than or equal to two.

SX — Vector of length NF containing the estimate of the spectral density of the first time series X.   (Input)

SY — Vector of length NF containing the estimate of the spectral density of the second time series Y.   (Input)

CPREAL — Vector of length N/2 +1 containing the real part of the cross periodogram between X and Y.   (Input)
The real part of the cross periodogram evaluated at (angular) frequency
ωk = 2πk/N is given by CPREAL(k + 1), k = 0, 1, …, N/2.

CPIMAG — Vector of length N/2. + 1 containing the imaginary part of the cross periodogram between X and Y.   (Input)
The imaginary part of the cross periodogram evaluated at (angular) frequency
ωk = 2πk/N is given by CPIMAG(k + 1), k = 0, 1, …, N/2.

F — Vector of length NF containing the (angular) frequencies at which the spectral density is estimated.   (Input)

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.

COSPEC — Vector of length NF containing the estimate of the cospectrum.   (Output)

QUADRA — Vector of length NF containing the estimate of the quadrature spectrum.   (Output)

CRAMPL — Vector of length NF containing the estimate of the cross-amplitude spectrum.   (Output)

PHASE — Vector of length NF containing the estimate of the phase spectrum.   (Output)

COHERE — Vector of length NF containing the estimate of the coherence.   (Output)

Optional Arguments

NF — Number of (angular) frequencies.   (Input)
F must be greater than or equal to one.
Default: NF = size (F,1).

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

FORTRAN 90 Interface

Generic:          CALL CSWEP (N, SX, SY, CPREAL, CPIMAG, F, WT, COSPEC, QUADRA, CRAMPL, PHASE, COHERE [,…])

Specific:                             The specific interface names are S_CSWEP and D_CSWEP.

FORTRAN 77 Interface

Single:            CALL CSWEP (N, SX, SY, CPREAL, CPIMAG, NF, F, NWT, WT, COSPEC, QUADRA, CRAMPL, PHASE, COHERE)

Double:                              The double precision name is DCSWEP.

Description

Routine CSWEP estimates the nonnormalized cross-spectral density function of two jointly stationary time series using a fixed sequence of weights given the modified cross-periodogram and spectral densities of the appropriately centered and padded data

and

for t = 1, …, N. The routine CPFFT may be used to compute the modified periodograms

and cross-periodogram

over the discrete set of nonnegative frequencies

(Here, a means the greatest integer less than or equal to a.) Either routine SSWP or routine SWEP may be applied to the periodograms to obtain nonnormalized spectral density estimates

over the set of frequencies

ω = ƒi,     i = 1, …, nƒ

where nf = NF. These frequencies are in the scale of radians per unit time. The time sampling interval Δt is assumed to be equal to one. Note that the spectral window or weight sequence used to compute

may differ from that used to compute

The cross-spectral density function is complex-valued in general and may be written as

The cospectrum is estimated by

and the quadrature spectrum is estimated by

where

and k(ω) is the integer such that ωk,0 is closest to ω. The sequence of m = NWT weights {wj} for
j = m/2, …, (m m/2 1) satisfies Σjwj = 1. These weights are fixed in the sense that they do not depend on the frequency ω at which to estimate hXY(ω). 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. The argument WT may contain relative weights since they are normalized to sum to one in the actual computations. The nonnormalized cross-spectral density estimate is computed over the same set of frequencies as the nonnormalized spectral density estimates. However, the particular weight sequence used to compute

need not correspond to either the weight sequence or spectral window used to compute either

An equivalent representation of hXY(ω) is the polar form defined by

The cross-amplitude spectrum is estimated by

and the phase spectrum is estimated by

Finally, the coherency spectrum is estimated by

The coherence or squared coherency is output.

Comments

1.         The periodograms of X and Y and cross periodogram between X and Y may be computed via the routine CPFFT. The spectral densities of X and Y may then be estimated using any of the routines SSWD, SWED, SSWP, or SWEP. Thus, different window types and/or weight sequences may be used to estimate the spectral and cross-spectral densities given either the series or their periodograms. Note that use of the modified periodograms and modified cross periodogram ensures that the scales of the spectral and cross-spectral densities and their estimates are equivalent.

2          The time sampling interval, TINT, is assumed to be equal to one. This assumption is appropriate for discrete parameter processes. The adjustment for continuous parameter processes (TINT > 0.0) involves multipication of the frequency vector F by 1/TINT and multiplication of the spectral and cross-spectral density estimates by TINT.

3.         To convert the frequency scale from radians per unit time to cycles per unit time, multiply F by 1/(2π).

Example

Consider the Robinson Multichannel Time Series Data (Robinson 1967, page 204) where X is the Wφlfer sunspot number and Y is the northern light activity for the years 1770 through 1869. Application of routine CSWEP to these data produces the following results.

 

      USE IMSL_LIBRARIES

 

      IMPLICIT   NONE

      INTEGER    LDCPM, LDCSM, LDRDAT, N, NDRDAT, NF, NOBS, &

                NPAD, NWT

      PARAMETER  (LDRDAT=100, NDRDAT=4, NF=10, NOBS=100, &

                NWT=7, LDCSM=NF, NPAD=NOBS-1, N=NOBS+NPAD, &

                LDCPM=N/2+1)

!

      INTEGER    I, IPVER, NROW, NVAR

      REAL       COHERE(NF), COSPEC(NF), CPIMAG(LDCPM), &

                 CPM(LDCPM,10), CPREAL(LDCPM), CRAMPL(NF), &

                 CSM(LDCSM,9), F(NF), FLOAT, PHASE(NF), PI, PX(LDCPM), &

                 PY(LDCPM), QUADRA(NF), RDATA(LDRDAT,NDRDAT), &

                 SX(NF), SY(NF), WT(NWT), X(NOBS), Y(NOBS)

      CHARACTER  CLABEL1(5)*24, CLABEL2(6)*16, FMT*8, RLABEL(1)*6, &

                 TITLE1*32, TITLE2*40

      INTRINSIC  FLOAT

 

!

      EQUIVALENCE (X(1), RDATA(1,2)), (Y(1), RDATA(1,3))

      EQUIVALENCE (PX(1), CPM(1,3)), (PY(1), CPM(1,6))

      EQUIVALENCE (CPREAL(1), CPM(1,9)), (CPIMAG(1), CPM(1,10))

!

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

      DATA FMT/'(F12.4)'/

      DATA CLABEL1/'%/%/ k', '%/%/Frequency', '%/%/Period', &

          'Spectral%/Estimate%/of X', 'Spectral%/Estimate%/of Y'/

      DATA CLABEL2/'%/ k', '%/Cospectrum', '%/Quadrature', &

          'Cross%/Amplitude', '%/Phase', '%/Coherence'/

      DATA RLABEL/'NUMBER'/

      DATA TITLE1/'Results of the Spectral Analyses'/

      DATA TITLE2/'%/Results of the Cross-Spectral Analysis'/

!                                 Initialization

      PI = 2.0*ASIN(1.0)

      DO 10  I=1, NF

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

         CALL SCOPY (NF, F, 1, CSM(1:,1), 1)

         CSM(I,2) = 2.0*FLOAT(NF)/FLOAT(I)

   10 CONTINUE

!                                 Robinson data

      CALL GDATA (8, RDATA, NROW, NVAR)

!                                 Center on arithmetic means

!                                 Frequency in radians per unit time

!                                 Modified periodogram version

      IPVER = 1

!                                 Compute the cross periodogram

      CALL CPFFT (X, Y, CPM, IPVER=IPVER)

!                                 Estimate the spectral densities

      CALL SWEP (N, PX, F, WT, SX)

      CALL SWEP (N, PY, F, WT, SY)

!                                 Estimate the cross-spectral density

      CALL CSWEP (N, SX, SY, CPREAL, CPIMAG, F, WT, COSPEC, &

                 QUADRA, CRAMPL, PHASE, COHERE)

!                                 Print results

!

!                                 Copy results to output matrices

      CALL SCOPY (NF, SX, 1, CSM(1:,3), 1)

      CALL SCOPY (NF, SY, 1, CSM(1:,4), 1)

      CALL SCOPY (NF, COSPEC, 1, CSM(1:,5), 1)

      CALL SCOPY (NF, QUADRA, 1, CSM(1:,6), 1)

      CALL SCOPY (NF, CRAMPL, 1, CSM(1:,7), 1)

      CALL SCOPY (NF, PHASE, 1, CSM(1:,8), 1)

      CALL SCOPY (NF, COHERE, 1, CSM(1:,9), 1)

!                                 Call printing routines

      CALL WRRRL (TITLE1, CSM, RLABEL, CLABEL1, NF, 4, FMT=FMT)

      CALL WRRRL (TITLE2, CSM(1:,5), RLABEL, CLABEL2, NF, 5, FMT=FMT)

!

      END

Output

 

             Results of the Spectral Analyses
                                    Spectral      Spectral
                                    Estimate      Estimate
 k     Frequency        Period          of X          of Y
 1        0.3142       20.0000      116.9550     1315.8370
 2        0.6283       10.0000     1206.6086     1005.1219
 3        0.9425        6.6667       84.8369      317.2589
 4        1.2566        5.0000       55.2120      270.2111
 5        1.5708        4.0000       46.5748      115.6768
 6        1.8850        3.3333       12.4050      250.0125
 7        2.1991        2.8571        7.0934       82.6773
 8        2.5133        2.5000        3.4091       62.3267
 9        2.8274        2.2222        5.6828       12.8970
10        3.1416        2.0000        4.0346       17.6441

                Results of the Cross-Spectral Analysis
                                      Cross
 k    Cospectrum    Quadrature     Amplitude         Phase     Coherence
 1       94.0531     -254.0125      270.8659        1.2162        0.4767
 2      702.5118       21.9823      702.8557       -0.0313        0.4073
 3       70.2379      -31.4431       76.9547        0.4209        0.2200
 4       -1.8715      -36.1639       36.2123        1.6225        0.0879
 5       36.6366      -18.5925       41.0843        0.4696        0.3133
 6       32.7071       -6.6569       33.3776        0.2008        0.3592
 7        9.4887       -9.1692       13.1950        0.7683        0.2969
 8        9.2534       -0.3000        9.2583        0.0324        0.4034
 9       -0.5568        2.9455        2.9977       -1.7576        0.1226
10        1.7640       -1.8321        2.5433        0.8043        0.0909

 



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