spec.pgram関数の使い方とtsオブジェクトのfrequencyの設定について

なかなか私を混乱させる仕様になっていたのでメモ。
ここでは例題としてcos関数のパワースペクトルを計算(表示)することにする。
以下のように条件設定。

#cosを何回繰り返すのか
REP <- 10
#cosの周期(この逆数が周波数f = 1/T = 1/5 = 0.2)
T <- 5
#サンプリングする個数
N <- 100

これから計算したいcos関数を
y = \cos(\frac{2\pi x}{T})
みたいな形だとしてまずは時系列データの生成。

#データ生成
x <- (REP * T) * (1:N) / N
y <- ts(cos(2 * pi * x / T),start=1,frequency=(N/(REP*T)))

ここで頭を悩ませたのがfrequency引数の設定でこれは結局単位時間あたりのデータ個数を入力するものらしい。今の場合は周期のオーダーで物をはかりたいのでN/(REP*T)としている。つまり周期が1のとき、その中に何個のデータ点があるのか?という値を指定している。

プロットすると

となってて、山から山(1周期)の感覚がだいたい5になってるのがわかる。

次にこれのパワースペクトルを計算してやる。
この関数は周波数fが0.2の所にしか成分を持っていないのでそういう結果が出ることを期待して実行すると・・・

spec.pgram(y, log = "no", )


となって、確かに周波数fが0.2(= 1/T = 0.2)となっている、正しい位置にしか成分を持たないことがわかる。



結局、大事なのはtsオブジェクトを作るときのfrequencyの設定を正しく単位時間あたりのデータ個数にするってこと。