# Spectral density estimation for limited data

We need an idea of the frequencies in the time series shown.

What can we expect from nonparametric spectrum estimation? Observing that much power is in frequency 4 years, of which we have only one period.

May we improve the result by proper window design, and what recommendations then?

We have modeled the logseries with an AR1 process. Far from perfect fit, but it serves some practical risk management purposes. How would you decide between nonparametric spectrum estimation or using the model's spectrum?

We are doing our research on spectral methods, but hope for some advice from practitioners. Given the limited data, how best go about obtaining important frequencies? When attempting to determine what frequencies contribute the most to the variance of a process my research group relies heavily on multitaper spectral estimation techniques. This involves estimating multiple spectra with different tapers / windows and using the weighted average as the overall estimate.

Specifically, we use Discrete Prolate Spheriodal Sequences (DPSS's or Slepian sequences (named after their "discoverer")) as our windows because they are doubly orthogonal and provide an optimal trade off between bias and variance (Spectrum Estimation and Harmonic Analysis - DJ Thomson - 1982).

One method to determine which frequencies are important is the harmonic F-test which determines "how sinusoidal" a particular frequency is in the data. There are other methods that can be used, the harmonic F-test is one of the more simple.

R has a software package "multitaper" and Matlab has the function "pmtm()". This method of estimation is parametric in the sense that in order to compute the windows, one parameters must be assigned: NW (the time bandwidth parameter) which is normally 3 < NW < 20 (for most purposes I have seen it is defaulted to 4). This parameter sets the width of the window (W) based on the amount of data that is available (N). Lastly, you must decide how many windows you would like to use (K) where K <= floor(2*NW-1). The tapers after 2*NW actually concentrate power outside of -W to W.

If you use the multitaper package in R, the defaults are all set and you can get a general idea of the spectrum with very little tinkering.

• Thanks! You're last phrase is great news for us :) Also, you nowhere say that it's foolish to expect results from such small dataset, am I interpreting you correctly? – user3817704 Feb 27 '15 at 21:04
• My understanding is that you require at least 2 periods worth of data in order to obtain a reliable result at a frequency of f=1/P (if 'f' is the frequency of interest with period, P). In other words - the lowest frequency that you would expect to get reliable results would be at a frequency of 1/<2 years> = 0.5 cycles/year. How many data points do you have out of curiosity? – driegert Feb 27 '15 at 21:16
• Yes sorry I haven't specified much. Time series is 250 weekly observations. Ok now I understand we shouldn't expect spectral info for periods above 2 years. The time series shown is the observed excitation of our system. We want to estimate its spectral content so as to predict spectrum of our system's output H(w)^2 S(w). – user3817704 Feb 27 '15 at 21:26