「同時度数分布と相関係数」の編集履歴(バックアップ)一覧に戻る

同時度数分布と相関係数 - (2006/06/09 (金) 21:04:46) のソース

とりあえず東証平均の収益率と同時期の新日鉄の収益率を入力する。
(ここではエクセルに入力してあったデータをコピペした)

 tosho<-scan()
 1: 1.9 0.0 -2.1 2.9 0.0 1.1 -0.3 1.6 2.7 1.5 -0.6 0.3 
 13: 3.3 -0.4 5.3 6.0 -0.1 4.7 1.0 -0.4 -6.1 0.4 0.4 2.9 
 25: 2.2 -4.9 -3.0 2.3 0.1 -1.3 -1.5 -0.1 -0.5 3.6 6.6 2.7 
 37: -0.9 0.7 5.0 2.4 0.9 3.1 1.5 2.0 2.4 -0.8 0.6 6.5 
 49: 6.1 -0.2 12.7 -1.0 -9.9 2.6 -4.1 7.0 1.4 4.6 2.0 4.3 
 61: 
 Read 60 items

 shinnittetu<-scan()
 1: 9.2 2.3 -6.5 9.0 5.3 -4.3 -3.7 7.0 7.6 1.4 -3.4 0.7 
 13: 2.8 -1.4 17.6 17.8 5.5 -1.9 1.9 9.0 -10.3 -10.3 -7.7 6.5 
 25: -0.6 -11.8 3.5 1.9 -5.5 -9.1 -5.7 2.3 -4.9 -0.8  8.0 6.7 
 37: -2.8 9.3 11.4 3.0 -7.5 2.5 -0.6 1.8 5.1 -2.3 -6.0 10.6 
 49: 0.0 -5.7 10.6 -0.6 -11.2 -3.8 -5.2 6.2 -4.2 2.1 0.6 4.7 
 61: 
 Read 60 items

コピペして変数に格納したデータをts()で、時系列データにしておく。

 tosho<-ts(tosho,frequency=12,start=1980)
 tosho
       Jan  Feb  Mar  Apr  May  Jun  Jul  Aug  Sep  Oct  Nov  Dec
 1980  1.9  0.0 -2.1  2.9  0.0  1.1 -0.3  1.6  2.7  1.5 -0.6  0.3
 1981  3.3 -0.4  5.3  6.0 -0.1  4.7  1.0 -0.4 -6.1  0.4  0.4  2.9
 1982  2.2 -4.9 -3.0  2.3  0.1 -1.3 -1.5 -0.1 -0.5  3.6  6.6  2.7
 1983 -0.9  0.7  5.0  2.4  0.9  3.1  1.5  2.0  2.4 -0.8  0.6  6.5
 1984  6.1 -0.2 12.7 -1.0 -9.9  2.6 -4.1  7.0  1.4  4.6  2.0  4.3

 shinnittetu<-ts(shinnittetu,frequency=12,start=1980)
 shinnittetu
        Jan   Feb   Mar   Apr   May   Jun   Jul   Aug   Sep   Oct   Nov   Dec
 1980   9.2   2.3  -6.5   9.0   5.3  -4.3  -3.7   7.0   7.6   1.4  -3.4   0.7
 1981   2.8  -1.4  17.6  17.8   5.5  -1.9   1.9   9.0 -10.3 -10.3  -7.7   6.5
 1982  -0.6 -11.8   3.5   1.9  -5.5  -9.1  -5.7   2.3  -4.9  -0.8   8.0   6.7
 1983  -2.8   9.3  11.4   3.0  -7.5   2.5  -0.6   1.8   5.1  -2.3  -6.0  10.6
 1984   0.0  -5.7  10.6  -0.6 -11.2  -3.8  -5.2   6.2  -4.2   2.1   0.6   4.7

東証平均と新日鉄のデータについて、データの範囲と大体の分布をみるために、5 数要約 fivenum() : 最小値・下側ヒンジ・中央値・上側ヒンジ・最大値を見る。

 fivenum(tosho)
 [1] -9.90 -0.35  1.25  2.90 12.70
 fivenum(shinnittetu)
 [1] -11.80  -4.25   1.05   5.85  17.80

さて、この2つのデータから同時度数分布を作ろう。
まず連続変数をカテゴリー変数に変換するために使用するのは cut 関数である。
区切り(-20~20、2刻み)を変数breaksにいれて


 breaks<-seq(-10,14,by=2)
 breaks2<-seq(-12,18,by=2)

cut()を使う。「○○以上□□未満」の場合は right = FALSE にする必要がある。

 c_tosho<-cut(tosho,breaks,right=F)
 c_shinnittetu<-cut(sinnittetu,breaks2,right=F)

カテゴリー化できたら、それらをつかってクロス表をつくる。
 table(c_tosho,c_shinnittetu)

これすなわち同時度数分布となっている。
          c_shinnittetu
 c_tosho    [-12,-10) [-10,-8) [-8,-6) [-6,-4) [-4,-2) [-2,0) [0,2) [2,4)
  [-10,-8)         1        0       0       0       0      0     0     0
  [-8,-6)          1        0       0       0       0      0     0     0
  [-6,-4)          1        0       0       1       0      0     0     0
  [-4,-2)          0        0       1       0       0      0     0     1
  [-2,0)           0        1       0       3       4      2     0     1
  [0,2)            1        0       2       4       0      1     3     1
  [2,4)            0        0       0       0       1      2     3     3
  [4,6)            0        0       0       0       0      1     0     1
  [6,8)            0        0       0       0       0      0     1     0
  [8,10)           0        0       0       0       0      0     0     0
  [10,12)          0        0       0       0       0      0     0     0
  [12,14)          0        0       0       0       0      0     0     0
          c_shinnittetu
 c_tosho    [4,6) [6,8) [8,10) [10,12) [12,14) [14,16) [16,18)
  [-10,-8)     0     0      0       0       0       0       0
  [-8,-6)      0     0      0       0       0       0       0
  [-6,-4)      0     0      0       0       0       0       0
  [-4,-2)      0     0      0       0       0       0       0
  [-2,0)       1     0      1       0       0       0       0
  [0,2)        1     1      2       0       0       0       0
  [2,4)        1     3      1       0       0       0       0
  [4,6)        1     0      0       1       0       0       1
  [6,8)        0     1      1       1       0       0       1
  [8,10)       0     0      0       0       0       0       0
  [10,12)      0     0      0       0       0       0       0
  [12,14)      0     0      0       1       0       0       0


これを3次元のグラフにしよう。
同時度数分布なので、同時度数グラフとか3次元ヒストグラムとか呼ばれる。

まず3次元データとして使いやすいようにデータを加工する。変数matにクロス集計を格納して、
 mat<-table(c_tosho,c_shinnittetu)

これを次のようなデータ・フレーム形式にする。
 s3d.dat <- data.frame(cols=as.vector(col(mat)),rows=as.vector(row(mat)),value=as.vector(mat))

何をやっているかは、結果を見た方が理解が早い。
 s3d.dat
     cols rows value
 1      1    1     1
 2      1    2     1
 3      1    3     1
 4      1    4     0
 5      1    5     0
 6      1    6     1
 7      1    7     0
 ……(中略)……
 176   15    8     1
 177   15    9     1
 178   15   10     0
 179   15   11     0
 180   15   12     0
つまりm行目のn列目の数値はいくつであるか、を並べているのである。

これをscatterplot3d()に渡してやると、まずは三次元にぶつぶつが浮かぶ。
 scatterplot3d(s3d.dat)
#ref(scatterplot1.gif)
これを棒グラフっぽくするには、type="h"とし、太さを lwd=5で少し太らせ、 pch=" "で、棒の先端についたキャラクタを外す。
 scatterplot3d(s3d.dat,type="h", lwd=5, pch=" ")
#ref(scatterplot2.gif)

**散布図
二つの量の関係を確かめるには、散布図が一番だ。
散布図は普通、
 plot(x,y)
の形で描けるが、いま扱っているデータは、toshoもshinnitetuも時系列データで、そのまま
 plot(tosho,shinnitetu)
では、時系列の順に結ばれる折れ線グラフとなる。
普通の散布図が欲しい場合には、as.vector()関数で、データをベクトル・データに変換してやればいい。

ここでは、少し工夫して、それぞれのヒストグラムをつけた散布図を書いてみる。
 xhist <- hist(tosho,  plot=FALSE) # 棒グラフ情報を記録
 yhist <- hist(shinnittetu, plot=FALSE) # 棒グラフ情報を記録
 top <- max(c(xhist$counts, yhist$counts))
 xrange <- c(-15,15)
 yrange <- c(-15,15)
 # layout 関数で画面を横に 3:1、縦に1:3 の4画面に分割
 # 副画面の番号を [1,1] -> 2, [1,2] -> 0, [2,1] -> 1, [2,2] -> 3 とする
 # 番号 0 の副画面には何も描かない(?) 
 nf <- layout(matrix(c(2,0,1,3),2,2,byrow=TRUE), c(3,1), c(1,3), TRUE)
 par(mar=c(3,3,1,1)) #番号1の副画面の余白指定
 plot(as.vector(tosho), as.vector(shinnittetu), xlim=xrange, ylim=yrange, xlab="", ylab="") # 番号1の副画面にプロット
 par(mar=c(0,3,1,1)) #番号2の副画面の余白指定
 barplot(xhist$counts, axes=FALSE, ylim=c(0, top), space=0) # 番号2の副画面にプロット
 par(mar=c(3,0,1,1)) #番号3の副画面の余白指定
 barplot(yhist$counts, axes=FALSE, xlim=c(0, top), space=0, horiz=TRUE) # 番号3の副画面にプロット

#ref(plot.gif)
目安箱バナー