トップページ -> データの扱い -> R システムを用いた統計,主成分分析,データマイニング -> R を用いた主成分分析
[サイトマップへ], [サイト内検索へ]

R を用いた主成分分析

サイト構成 データベース関連技術 データの扱い インストール,設定,利用 プログラミング 情報工学の講義実習資料 サポートページ 連絡先,業績など

この Web ページでは,R を用いた主成分分析のための操作手順の一例を図解で説明する ここでは,R の princomp() 関数を使う.

【補足】princomp() 関数では,特に指定しない限り, 不偏分散共分散行列 あるいは 不偏相関係数行列 が使われます (このことは R の princomp() 関数のマニュアルに書いてあります).


あらかじめ決めておく事項

主成分分析を行う CSV ファイル名と属性名を決めておいてください.

この Web ページでは,次のように書きます.


前準備

R システムのインストール


princomp() 関数を使って主成分分析を行う手順

princomp() 関数では,

を選ぶことができる.どちらが「良い」とは簡単に言えません.この Web ページでは,どちらか片方を選ぶことはせず, 両方の手順を併記します.


CSVファイルを読み込み,データテーブルに格納

  1. (前準備) 使用する CSV ファイルの作成

    ※ ここでは Book1.csv をダウンロードし,分かりやすいディレクトリに置く

    (参考: 「外国為替データ(時系列データ)の情報源の紹介」の Web ページ

    以下の説明では、

    として説明を続ける.

    自前の CSV ファイルを使うときの注意: read.table() 関数を使うので, 属性名は英語になっていること.属性名は,CSV ファイルの第一行目に書いていること.

  2. 使用する CSV ファイルの確認

    属性名が CSV ファイルの1行目に書かれていることを確認する.

    [image]
  3. R の起動

    ◆ Windows での動作画面例

    [image]

    ◆ Ubuntu での動作画面例

    [image]
  4. read.table」を用いて,ファイルを変数に読み込み

    端末で,次のコマンドを実行.

    ◆ Windows での動作手順例

    X <- read.table("C:/R/Book1.csv", header=TRUE, sep=",", na.strings="NA", dec=".", strip.white=TRUE);
    

    [image]

    ◆ Linux での動作手順例

    X <- read.table("/tmp/Book1.csv", header=TRUE, sep=",", na.strings="NA", dec=".", strip.white=TRUE);
    

    [image]

    オプション

  5. 変数 X の中身の確認

    端末で,次のコマンドを実行.。

    edit(X);
    

    [image]

    端末で,次のコマンドを実行.。

    str(X)
    

    [image]

主成分分析不偏相関係数行列から主成分を求める場合)

  1. 主成分分析不偏相関係数行列から主成分を求める場合)

    PC <- princomp(~USD+EUR+AUD, cor=TRUE, data=X)
    PC
    

    [image]
  2. (オプション)確認のため PC の構造の表示

    [image]
  3. 固有ベクトルの表示

    ※ "loadings" という名前が付いているので,「主成分負荷量」のことでは無いだろうかと勘違いしやすいので注意してください.princomp() 関数の場合には,固有ベクトルが求まります.

    unclass(loadings(PC))
    

    [image]

    結果が3行3列の行列になっているのは間違いではない. princomp() 関数が扱うデータフレームは, 列数を変数の数,行数を標本数とみなし,i 行 j 列目の要素は,第 j 変数の第 i 番目の標本という意味になる.その主成分分析の結果である主成分は,変数の数を行数及び列数とする行列をなす.

    ※ 第1列(一番左の列)が第1主成分,第2列が第2主成分, 第3列が第3主成分.主成分の数は元の変数の数と一致する.

  4. 固有値の表示

    PC$sd^2
    

    [image]

    固有値の総和は3になる.

  5. 主成分負荷量

    主成分負荷量は,構想係数とも呼ばれます. 相関係数行列から主成分を求めた場合,「主成分負荷量=(固有値)の平方根×固有ベクトル」が成り立ちます.次の式で求めることができる.

    t( PC$sd * t( PC$loadings ) )[, drop = FALSE]
    

    [image]
  6. 主成分分析不偏分散共分散行列から主成分を求める場合)

    princomp()関数では,「cor=FALSE」と指定する.これは,不偏分散共分散行列を使うという意味.

    PC <- princomp(~USD+EUR+AUD, cor=FALSE, data=X)
    PC
    

    [image]


主成分分析の結果のグラフ表示

  1. 相関係数行列の固有ベクトルがなす行列から2列を選ぶ

    2次元の散布図を書きたいので、そのために、主成分分析で得られた主成分の第1列と第2列を選ぶ方法を練習しておきます.

    ※ [,c(1, 2)] と書いて,第1列と第2列を選ぶ.

    PC <- princomp(~USD+EUR+AUD, cor=TRUE, data=X)
    unclass(loadings(PC))[,c(1, 2)]
    

    [image]
  2. 主成分スコアを求める

    「unclass(loadings(PC))[,c(1, 2)]」で2列を選んでいるので,スコアは2個になる.

    ※ 「data.matrix( norm[,3:5] ) %*% eigen( cor( norm[,3:5] ) )$vectors[,c(1,2)]」と書いても同じ結果になる.

    D <- data.matrix( X[,c(3,5,11)] ) %*% unclass(loadings(PC))[,c(1, 2)]
    edit(D); 
    

    [image]

    [image]
  3. 第一主成分と第二主成分から2次元の散布図を描く.

    plot( data.matrix( X[,c(3,5,11)] ) %*% unclass(loadings(PC))[,c(1, 2)] )
    

    [image]
  4. biplot() 関数を用いた主成分分析結果のバイプロット

    R の biplot() 関数を使うと, 第一主成分と第二主成分から2次元の散布図を描くことが簡単にできる.

    ※ 上の plot() 関数での結果と同様のグラフになっているので,確かに,biplot() 関数は, 第一主成分と第二主成分から2次元の散布図を描いているということが確認できる.

    biplot( PC )
    

    [image]

    [image]
  5. screeplot() 関数

    screeplot( PC )
    

    [image]

    [image]

(オプション)より基本的な関数を用いた主成分分析

主成分分析を行う手法としては,大きく分けて,相関係数行列による方法と,分散共分散行列による方法が知られています.

princomp() 関数を使う方法だと、中身がブラックボックスです. より理解を深めるために,固有値分解を行うeigen() 関数と, 不偏相関係数行列を求める cor() 関数と, 不偏分散共分散行列を求める cor() 関数を組み合わせて,主成分を求めてみます. つまり,以下の通りです.

R による分散共分散行列 の Web ページで,不偏相関係数行列と 不偏分散共分散行列を説明している.

  1. データの読み込み

    下記のデータを変数 X に読み込む

    X <-read.table("C:/R/Book1.csv", header=TRUE, sep=",", na.strings="NA", strip.white=TRUE )
    X
    

    [image]
  2. 不偏相関係数行列の固有値と固有ベクトルを求める

    固有値の総和は3 $vectorsの各列が主成分.つまり,1番左の列が第一主成分

    eigen( cor( X[,c(3,5,11)] ) )
    

    [image]

    princomp() 関数と比べてみましょう.
    同じ結果が得られていることが分かります.
    「unclass(loadings(PC))」のところは固有ベクトル, 「PC$sd^2」のところは固有値, 「t( sqrt ( PC$sd ...」のところは主成分負荷量です.

    PC <- princomp(~USD+EUR+AUD, cor=TRUE, data=X)
    unclass(loadings(PC))
    PC$sd^2
    t( PC$sd * t( PC$loadings ) )[, drop = FALSE]
    

    書きかけ

  3. 不偏分散共分散行列の固有値と固有ベクトルを求める

    $vectorsの各列が主成分.つまり,1番左の列が第一主成分

    eigen( cov( X[,c(3,5,11)] ) )
    

    [image]

    princomp() 関数と比べてみましょう.
    固有ベクトルは同じ値になっていますが,固有値は値が少し違うということが見て取れます.
    「unclass(loadings(PC))」のところは固有ベクトル, 「PC$sd^2」のところは固有値, 「t( sqrt ( PC$sd ...」のところは主成分負荷量です.

    PC <- princomp(~USD+EUR+AUD, cor=FALSE, data=X)
    unclass(loadings(PC))
    PC$sd^2
    t( PC$sd / sqrt( c( sd( X[,3] ), sd( X[,4] ), sd( X[,5] ) ) ) * t( PC$loadings ) )[, drop = FALSE]
    

    書きかけ


princomp2() 関数を用いた主成分分析

http://aoki2.si.gunma-u.ac.jp/R/princomp2.html で公開されている princomp2() 関数が便利です.実行例を以下に示しておきます.


本サイトは金子邦彦研究室のWebページです.サイトマップは,サイトマップのページをご覧下さい. 本サイト内の検索は,サイト内検索のページをご利用下さい.

問い合わせ先: 金子邦彦(かねこ くにひこ) [image]