とあるデータのスペクトル解析

結論から言うと、今回は「やりたいことが実現できなーい。どうすればいいのか分からん!」という話。
なので、読んでもいいことないです。いつも以上に単なる備忘録です。

さてさて、今回はこんなデータ。

CSVを読み込んでtsに変換。

jniv<-read.csv("データ.csv")
jniv.close<-jniv[2]
jniv.close.ts<-ts(jniv.close)
jniv.close.ts<-ts(as.numeric(jniv.close.ts))

スペクトル解析

spec.pgram(jniv.close.ts,log="no", xlim=c(0,0.002))

グラフ化せずに生データを見たいときは、「plot=FALSE」とする。

> spec.pgram(jniv10Year.close.ts,log="no",xlim=c(0,0.002), plot=FALSE)
$freq
   [1] 0.0002222222 0.0004444444 0.0006666667 0.0008888889 0.0011111111
   [6] 0.0013333333 0.0015555556 0.0017777778 0.0020000000 0.0022222222
  [11] 0.0024444444 0.0026666667 0.0028888889 0.0031111111 0.0033333333
  [16] 0.0035555556 0.0037777778 0.0040000000 0.0042222222 0.0044444444
  
・・・・・・・・・・・・・・・・・
[2236] 0.4968888889 0.4971111111 0.4973333333 0.4975555556 0.4977777778
[2241] 0.4980000000 0.4982222222 0.4984444444 0.4986666667 0.4988888889
[2246] 0.4991111111 0.4993333333 0.4995555556 0.4997777778 0.5000000000

$spec
   [1] 3.063381e+03 2.232736e+04 3.555620e+04 4.828382e+03 1.111915e+04
   [6] 3.021952e+03 9.457354e+03 1.109082e+04 8.995889e+02 2.844646e+02
  [11] 9.599256e+02 5.158237e+03 4.985057e+03 2.780245e+03 5.226908e+03
  [16] 8.107157e+02 1.879867e+03 1.869383e+03 2.184182e+02 1.788623e+03
  [21] 5.462029e+03 1.010075e+03 7.125709e+02 3.086387e+03 2.357322e+03
  [26] 8.773852e+02 2.078958e+02 8.100285e+02 9.715714e+01 4.294267e+02
  [31] 2.259698e+03 1.571764e+03 1.130904e+03 2.759914e+02 8.386062e+01
 ・・・・・・・・・・・・・・・・・・
[2236] 2.092352e+00 4.042008e-01 1.339844e+00 1.866943e+00 7.919000e-01
[2241] 2.731446e+00 9.306263e-01 1.085729e+00 3.989569e-01 1.309799e-01
[2246] 8.265553e-01 2.611334e-01 9.262818e-02 2.484483e+00 3.953657e+00

$coh
NULL

$phase
NULL

$kernel
NULL

$df
[1] 1.75974

$bandwidth
[1] 6.415003e-05

$n.used
[1] 4500

$orig.n
[1] 4420

$series
[1] "jniv10Year.close.ts"

$snames
NULL

$method
[1] "Raw Periodogram"

$taper
[1] 0.1

$pad
[1] 0

$detrend
[1] TRUE

$demean
[1] FALSE

attr(,"class")
[1] "spec"

◎スペクトルの最大値を求める

max(spec.pgram(jniv10Year.close.ts,log="no",xlim=c(0,0.002), plot=FALSE)$spec)

ここで漸く本題。実は、この最大値に対応した周波数を知りたい。だけど求め方が分からないって話。
周波数は、$freqがnumeric型なので、which関数を使えば次のように(何番目の値かが)求められるかと思っていたが駄目だった。

> vec<-c(1,3,6,3,5)
> class(vec)
[1] "numeric"
> which(vec==1)
[1] 1
> which(vec==3)
[1] 2 4
> which(vec==6)
[1] 3

ところが、which(spec.pgram(jniv10Year.close.ts,log="no",xlim=c(0,0.002), plot=FALSE)$freq==スペクトルの最大値)としても戻り値はinteger(0)。
これは対応する添字が存在しない場合に返却される。
ここでお手上げ。

恐らく、数値計算の誤差の問題だろうから、自分で数値計算して求めないといけないのかなと思っているが、関数無いのかな????