crsspectrum(y, len [,reps]), y a REAL matrix with ncols(y) >= 2, len > 0 and reps > 0 integers. |

crsspectrum(y, len) computes smoothed periodogram and cross-periodogram estimates of the spectrum and cross-spectrum of the multivariate time series in the columns of REAL matrix y. Smoothing is accomplished by 4 successive convolutions of a "boxcar" smoother of length len. That is, the smoothing weights are rep(1/len,len). The total number of frequencies involved in each spectrum estimate is 4*len - 3. When len = 1 no smoothing is done. crsspectrum(y, len, reps) does the same, except that reps convolutions with a boxcar are used, where reps > 0 is an integer scalar. When reps is odd, then len must also be odd. See topic 'bandwidth' for information on the approximate bandwidth and equivalent degrees of freedom associated with these estimates. Suppose p = ncols(y) and N = nrows(y). Then, after Sy <- crsspectrum(y, len [,reps]), Sy is a N by Q matrix, where Q = p + p*(p-1)/2 = p*(p+1)/2. Sy[,j], j = 1,...,p are the estimated spectra of y[,j], in Real form. Sy[,p+1], Sy[,p+2],...,Sy[,Q] contain the estimated cross spectra of y[,i] and y[,j], i < j, in Hermitian form. The order is (i,j) = (1,2), (1,3), ..., (1,p), (2,3),..., (2,p), ... (p-1,p) Specifically, the estimated cross spectrum of y[,i] and y[,j] is in column i*(p - (i + 1)/2) + j. For example, when x and y are vectors, crsspectrum(hconcat(x,y), len) returns hconcat(Sxx, Syy, Sxy), where Sxx and Syy are estimated spectra and Sxy is the estimated cross spectrum. Sxx and Syy are in Real form and Sxy is in Hermitian form. The column means are subtracted before computing the periodograms and cross-periodograms, but no tapering or other detrending is done. The spectra and cross spectra are computed at the Nfreq Fourier frequencies 0, 1/Nfreq, 2/Nfreq, ..., (Nfreq-1)/Nfreq cycles per unit time, were Nfreq is determined as follows: When variable S is defined and is a positive integer, Nfreq = S. It is an error if S has a prime factor > 29. Otherwise, Nfreq = smallest integer >= 2*nrows(y) which has no prime factors > 29, that is Nfreq = goodfactors(2*nrows(y)). Most of the time, a better choice for estimating cross spectra is macro compfa() which includes tapering and detrending and for which the amount of smoothing is specified as a bandwidth or an equivalent degrees of freedom. See compfa() for details.

Gary Oehlert 2003-01-15