3 つの変数を n 回観測して得られる n 個の標本から構成される長さ n のベクトル に対する不偏分散共分散行列と不偏相関係数行列の定義は次の通り。
R システムの CRAN の URL: https://cran.r-project.org/
Book1.csv をダウンロード (参考: 「外国為替データ(時系列データ)の情報源の紹介」の Web ページ)
以下の説明では、
として説明を続ける.
※ 自前の CSV ファイルを使うときの注意: read.table() 関数を使うので, 属性名は英語になっていること.属性名は,CSV ファイルの第一行目に書いていること.
属性名が CSV ファイルの1行目に書かれていることを確認する.
次のコマンドを実行.
◆ Windows での動作手順例
X <- read.table("C:/R/Book1.csv", header=TRUE, sep=",", na.strings="NA", dec=".", strip.white=TRUE);
◆ Linux での動作手順例
X <- read.table("/tmp/Book1.csv", header=TRUE, sep=",", na.strings="NA", dec=".", strip.white=TRUE);
R の read.table のオプション
次のコマンドを実行.
edit(X);
次のコマンドを実行.
str(X)
cov() 関数を使って不偏分散共分散行列を求める
次のコマンドを実行. データはUSDとEURとAUDの3列(それぞれ3列目と4列目と5列目)なので、 「3:5」と書く.これは「c(3,4,5)」と同じ意味.
不偏分散共分散行列が結果として得られる. 結果が3行3列の配列になっている.
cov( X[,3:5] )
R の cov() 関数を使う方法だと、中身がブラックボックスです. 実は,以下の 2 つの R の式が,同じ結果になります. ここでは、確かに同じ結果になることを確認する手順を,図解などで説明する.
cov( X[,3:5] )
crossprod( D ) / ( nrow(D) - 1 ) ※ この式は「t( D ) %*% t( D ) / ( nrow(D) - 1 )」と同じ意味
但し、D は次の手順で作る. D は,行列 X の各列の平均値を X の各要素から引き,転置した行列になっている. ここで「転置」になっているのは R のリサイクル規則を使って簡単に作れるから。
D <- t( t( X[,3:5] ) - colMeans( X[,3:5] ) )
この Web ページの先頭で説明したように、不偏分散共分散行列は,「不偏分散を拡張した形の分散共分散行列」になっていることが分かります.
上では「crossprod( D ) / ( nrow(D) - 1 )」, 「D <- t( t( X[,3:5] ) - colMeans( X[,3:5] ) )」のように書きました。 式が複雑なのでこれを一歩一歩確かめて、R の理解を深めます.
次のコマンドを実行.
◆ Windows での動作手順例
X <- read.table("/tmp/Book1.csv", header=TRUE, sep=",", na.strings="NA", dec=".", strip.white=TRUE);
◆ Linux での動作手順例
X <- read.table("/tmp/Book1.csv", header=TRUE, sep=",", na.strings="NA", dec=".", strip.white=TRUE);
read.table() のオプション
次のコマンドを実行.
edit(X);
データはUSDとEURとAUDの3列(それぞれ3列目と4列目と5列目)あることに注意して下さい.この 3 列について, 各列ごとに平均を求めます. 「3:5」は「c(3,4,5)」と同じ意味.平均を求めた結果はベクトルになる.
colMeans( X[,3:5] )
つまり、次のことを行いたいのです.
※ ここで行っていることは, 変数 X 内の3列目と4列目と5列目がなす「3次元のベクトル集合」の平均が原点になるように,平行移動させる操作とみなすことができる.
行列からベクトルの引き算を行うとき、リサイクル規則が使われます. ここで、リサイクル規則を、例を使って、簡単に説明する. 次のような,要素 1, 2, 3, 4, 5, 6 が入った3行2列の行列があるとします.
1 4 2 5 3 6
この行列に対して、ベクトル c(2, 5) の引き算を行うと,リサイクル規則から、 ベクトル c(2, 5) が、次のような行列に引き伸ばされる.
2 5 5 2 2 5
matrix( c( 1, 2, 3, 4, 5, 6 ), nrow = 3, ncol = 2 ) matrix( c( 1, 2, 3, 4, 5, 6 ), nrow = 3, ncol = 2 ) - c( 2, 5 )
実行結果の例は次の通りです.
今度は,要素 1, 4, 2, 5, 3, 6 が入った2行3列の行列があるとします.
1 2 3 4 5 6
この行列に対して、ベクトル c(2, 5) の引き算を行うと,リサイクル規則から、 ベクトル c(2, 5) が、次のような行列に引き伸ばされる.
2 2 2 5 5 5
matrix( c( 1, 4, 2, 5, 3, 6 ), nrow = 2, ncol = 3 ) matrix( c( 1, 4, 2, 5, 3, 6 ), nrow = 2, ncol = 3 ) - c( 2, 5 )
実行結果の例は次の通りです.
行列の行数と,ベクトルの長さが等しいときは, 行列の各要素から,ベクトルの値を引く(行列の行番号と、ベクトルの要素番号が等しい)ことが見て取れます.
「行列の各要素から,列の平均値を引く」ことをしたいので, 行列 X[,3:5] を転置させて,3行610列の行列を作り,列の平均値のベクトル(ベクトルの長さは3)で引き算します. リサイクル規則が使えるので、次のような簡単な式になる.(結果は3行610列の行列です).
※ データを転置したいので関数 t() を使う.
D <- t( t( X[,3:5] ) - colMeans( X[,3:5] ) ) edit(D);
上のように書く代わりに、次のように書いても同じ意味ですが、リサイクル規則を使わないと、書くのが面倒です ※ 「V[,3] - colMeans(V[,3])」は,「3列目」の全ての要素から「3列目の平均値」を引くという意味.
hoge <- matrix( c( ( X[,3] - colMeans(X[,3]) ), ( X[,4] - colMeans(X[,4]) ), ( X[,5] - colMeans(X[,5]) ) ), nrow = 2834, ncol = 3 );
crossprod() 関数はクロス積(この場合 t(D) %*% D)を求める関数です.D は、610行3列の行列になっていることに注意.
crossprod( D ) / ( nrow(D) - 1 )