Lieven Desmet
2007-Dec-12 15:04 UTC
[R] discrepancy between periodogram implementations ? per and spec.pgram
hello, I have been using the per function in package longmemo to obtain a simple raw periodogram. I am considering to switch to the function spec.pgram since I want to be able to do tapering. To compare both I used spec.pgram with the options as suggested in the documentation of per {longmemo} to make them correspond. Now I have found on a variety of examples that there is a shift between the log of the periodogram with per and that with spec.pgram. This vertical shift amounts to approx. 1.8 on the log scale (the graph of spec.pgram being above the one from per). Is there some explanation for this ? Is the one from spec.pgram the better one as suggested in the documentation of per {longmemo}? Finally how are these related to an estimate of the spectral density obtained from spec.arima ? Many thanks for help and clarification. Lieven Desmet Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm
Prof Brian Ripley
2007-Dec-12 15:44 UTC
[R] discrepancy between periodogram implementations ? per and spec.pgram
There are several definitions of a periodgram. Note that> log(2*pi)[1] 1.837877 See the comments in ?spectrum about scalings. I think the comments in ?per incorrectly ignore the scaling issues: per() does not take the base frequency into account and has an extra divisor of 2*pi. E.g.> x <- rnorm(64) > spec.pgram(x, taper=0, detrend=F)$spec/per(x)[-1][1] 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 [9] 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 [17] 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 [25] 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 On Wed, 12 Dec 2007, Lieven Desmet wrote:> hello, > > I have been using the per function in package longmemo to obtain a > simple raw periodogram. > I am considering to switch to the function spec.pgram since I want to be > able to do tapering. > To compare both I used spec.pgram with the options as suggested in the > documentation of per {longmemo} to make them correspond. > > Now I have found on a variety of examples that there is a shift between > the log of the periodogram with per and that with spec.pgram. This > vertical shift amounts to approx. 1.8 on the log scale (the graph of > spec.pgram being above the one from per). > > Is there some explanation for this ? Is the one from spec.pgram the > better one as suggested in the documentation of per {longmemo}? Finally > how are these related to an estimate of the spectral density obtained > from spec.arima ?What is spec.arima? If you meant spec.ar, that is on the same scale as spec.pgram for series with base frequency 1 (and for all series for R >= 2.7.0).> Many thanks for help and clarification. > > Lieven Desmet-- Brian D. Ripley, ripley at stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UK Fax: +44 1865 272595
Lieven Desmet
2007-Dec-14 14:31 UTC
[R] discrepancy between periodogram implementations ? per and spec.pgram
Prof Brian Ripley wrote:> There are several definitions of a periodgram. Note that > >> log(2*pi) > > [1] 1.837877 > > See the comments in ?spectrum about scalings. > > I think the comments in ?per incorrectly ignore the scaling issues: > per() does not take the base frequency into account and has an extra > divisor of 2*pi. E.g. > >> x <- rnorm(64) >> spec.pgram(x, taper=0, detrend=F)$spec/per(x)[-1] > > [1] 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 > 6.283185 > [9] 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 > 6.283185 > [17] 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 > 6.283185 > [25] 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 6.283185 > 6.283185 > > > On Wed, 12 Dec 2007, Lieven Desmet wrote: > >> hello, >> >> I have been using the per function in package longmemo to obtain a >> simple raw periodogram. >> I am considering to switch to the function spec.pgram since I want to be >> able to do tapering. >> To compare both I used spec.pgram with the options as suggested in the >> documentation of per {longmemo} to make them correspond. >> >> Now I have found on a variety of examples that there is a shift between >> the log of the periodogram with per and that with spec.pgram. This >> vertical shift amounts to approx. 1.8 on the log scale (the graph of >> spec.pgram being above the one from per). >> >> Is there some explanation for this ? Is the one from spec.pgram the >> better one as suggested in the documentation of per {longmemo}? Finally >> how are these related to an estimate of the spectral density obtained >> from spec.arima ? > > > What is spec.arima? If you meant spec.ar, that is on the same scale > as spec.pgram for series with base frequency 1 (and for all series for > R >= 2.7.0). > > >> Many thanks for help and clarification. >> >> Lieven Desmet > >Dear Prof. Ripley, thanks very much for a quick and helpful response. In the last question I wanted to hint at specARIMA which I am using to get the theoretical spectral density of an ARMA process. This works very well in general, however, in a simple example X_t=0.7*X_{t-1}+epsilon_t I obtain a value 1.768253 for funscaled[1] ( the first Fourier frequency 0.003141593) using str(f<-specARIMA(eta=c(H=0.5,phi=c(0.7),psi=c()),p=1,q=0,m=2000)) funscaled<-numeric(length(f$freq)) funscaled<-f$spec*f$theta1 where the theoretical value should be 0.901878 with b<-0.7 omega<-0.003141593 1/(2*pi)*(1-b^2)/(1+b^2-2*b*cos(omega)) [1] 0.9018088 using the formula (2.40) in Fan and Yao, Nonlinear Time Series ( Springer 2003 ), page 54-55 Is there also a simple explanation for this ? am I overlooking something ? Thanks and best regards, Lieven Desmet, maths dept - KULeuven - Belgium Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm