Lesson5.Factor Analysis
Kosugitti
2015年11月20日
毎回やる儀式
- インターネットに接続しておいてください（このサイトが見えていればOK）。
- Rstudioを起動し，この授業のプロジェクトを開いてください。(File>Open Project)
- File>New File>R scriptで新しいコードファイルを作ってください。
この度ファイルをちょっと更新しました。2015がついていないbaseball,batter,pitcherデータを用意してください。
ピッチャーのデータ，pitcher.csvをb.dataというオブジェクトに入れているとします。 できてないひとは，次の通りやってみよう。
p.data <- read.csv("pitcher.csv")今日使うパッケージ
- psych
- lavaan
目に見えないものを見えるように
目に見えないものを見えるようにというとなんだか魔法のような，怪しげな術があるのかと思うかもしれませんが，今回のテーマはまさにそういうことなのです。
目に見えないもの＝心とまでは言いませんが，ここで言わんとすることを専門的に言うと「潜在変数」とか「合成変数」ということになります。
例えば「身長」と「体重」と「年齢」を組み合わせた変数でプロ野球選手を表現してみたいと思います。プロ野球選手の特徴を際立たせるのが目的ですから（みんな同じになっても仕方がないので），この三つの変数に重みをつけて，なにかその特徴を表現する変数を作ることをイメージをしてみてください。
数学的にはこうなります。 \[ P = w_1 * height + w_2 * weight + w_3 * age\]
ここで，作られる変数\(P\)の分散を最大にするように，三つの変数\(w_1,w_2,w_3\)の重みを決めます。ここで\(P\)を主成分と呼び，この方法が主成分分析(PCA,Principle Component Analysis)と言われるものです。数学的に解を得るために\(\sum w_i^2 = 1\)という上限の制約をかけます。
こうして作られた合成変数\(P\)はプロ野球選手の特徴を表現したスコアです。これを「選手力」と呼べばほら，目に見えなかったものが見えるようになったでしょ？
今日のお話はそういう感じのものです。
因子分析と主成分分析
今日のテーマは因子分析で，主成分分析ではありません。 この二つのどこが違うかというと，
- 相関係数から計算を始めるのが因子分析，分散共分散行列から計算を始めるのが主成分分析
- 項目(変数)に誤差を想定するのが因子分析，想定しないのが主成分分析
の2点が挙げられます。
主成分分析は標準化してないスコアでも計算できます。因子分析は普通，変数間関係を標準化した相関行列から分析を始めます。
主成分分析は変数のスコアをそのまま使いますが，因子分析は知能テストや性格検査を背景に生まれてきたこともあって，得られた変数の数字には誤差が含まれているのではないかという過程をおきます。この辺の詳しいところは，心理教育測定法で説明しますね。
とりあえずやってみましょう
とりあえず因子分析をやってみましょう。psychパッケージを使いますので，まだの方はインストールしてください。
使う変数を一部に限定します。量的変数だけにし，内容的に重複するものは避けました。
# ここで使う変数は次の通りです。
fa.dat <- subset(p.data, select=c("win","lose","save","hold","hits","HR","BB",
"DB","K","W.pitch","R"))
# 勝利 win # 敗北 lose # セーブ save
# 被安打 hits # 被本塁打 HR # 与四球 BB # 与死球 DB
# 奪三振 K # 暴投 W.pitch # 失点 R因子数を決める手法
因子分析を始める場合，まず説明する潜在変数＝因子の数をいくつにするか，決めなければなりません。因子数の決め方は幾つかあります。
- スクリープロットを見て考える
- 固有値が1.0を超えるところまで採用する
- 累積寄与率が5割をこえるところまで採用する
- 適合度指標をみて最適値を得る
- 乱数に基づくデータと比べる（平行分析）
- 因子で相関データを綺麗に説明できるように (最小平平均偏相関法，MAP)
1-3は古典な方法，4-6はより客観的で現代的な方法です。おすすめは5と6です。
実際にやってみましょう。
因子数を決める
# psychパッケージを使います
# install.packages("psych")
library(psych)## Warning: package 'psych' was built under R version 3.2.2fa.parallel(fa.dat,SMC=TRUE)
## Parallel analysis suggests that the number of factors = 5 and the number of components = 1VSS(fa.dat)
##
## Very Simple Structure
## Call: vss(x = x, n = n, rotate = rotate, diagonal = diagonal, fm = fm,
## n.obs = n.obs, plot = plot, title = title, use = use, cor = cor)
## VSS complexity 1 achieves a maximimum of 0.91 with 1 factors
## VSS complexity 2 achieves a maximimum of 0.93 with 2 factors
##
## The Velicer MAP achieves a minimum of 0.06 with 1 factors
## BIC achieves a minimum of NA with 2 factors
## Sample Size adjusted BIC achieves a minimum of NA with 6 factors
##
## Statistics by number of factors
## vss1 vss2 map dof chisq prob sqresid fit RMSEA BIC SABIC complex
## 1 0.91 0.00 0.056 44 136.192 2.3e-11 3.73 0.91 0.191 -48 90.37 1.0
## 2 0.50 0.93 0.069 34 73.655 9.5e-05 3.11 0.93 0.146 -69 38.25 1.7
## 3 0.44 0.80 0.088 25 49.591 2.4e-03 2.42 0.94 0.136 -55 23.55 2.0
## 4 0.40 0.68 0.115 17 28.847 3.6e-02 1.88 0.96 0.118 -42 11.14 2.3
## 5 0.30 0.63 0.134 10 11.922 2.9e-01 1.38 0.97 0.074 -30 1.51 2.6
## 6 0.39 0.65 0.186 4 4.897 3.0e-01 0.97 0.98 0.079 -12 0.73 2.5
## 7 0.32 0.58 0.266 -1 1.059 NA 0.94 0.98 NA NA NA 2.7
## 8 0.33 0.59 0.411 -5 0.046 NA 0.81 0.98 NA NA NA 2.6
## eChisq SRMR eCRMS eBIC
## 1 46.451 0.0800 0.089 -138
## 2 33.595 0.0680 0.087 -109
## 3 21.565 0.0545 0.081 -83
## 4 12.420 0.0414 0.074 -59
## 5 3.676 0.0225 0.053 -38
## 6 0.800 0.0105 0.039 -16
## 7 0.325 0.0067 NA NA
## 8 0.022 0.0017 NA NA結果として，平行分析では2因子，MAP法では1因子，適合度基準(BIC)では2因子が良いということでした。ここでは仮に2因子だとしましょう。
因子分析を行う
fa.result <- fa(fa.dat,nfactors =2,fm="ML")## Loading required namespace: GPArotationprint(fa.result,sort=T)## Factor Analysis using method = ml
## Call: fa(r = fa.dat, nfactors = 2, fm = "ML")
## Standardized loadings (pattern matrix) based upon correlation matrix
## item ML1 ML2 h2 u2 com
## R 11 0.97 0.04 1.00 0.0049 1.0
## HR 6 0.93 -0.14 0.70 0.3023 1.0
## lose 2 0.74 0.13 0.69 0.3074 1.1
## hits 5 0.57 0.49 0.97 0.0285 2.0
## save 3 -0.53 0.13 0.20 0.8000 1.1
## BB 7 0.52 0.38 0.70 0.3042 1.8
## hold 4 -0.43 -0.19 0.34 0.6604 1.4
## DB 8 0.36 0.27 0.34 0.6600 1.8
## W.pitch 10 0.25 0.22 0.19 0.8058 2.0
## K 9 -0.02 0.86 0.71 0.2932 1.0
## win 1 0.06 0.83 0.76 0.2389 1.0
##
## ML1 ML2
## SS loadings 4.10 2.49
## Proportion Var 0.37 0.23
## Cumulative Var 0.37 0.60
## Proportion Explained 0.62 0.38
## Cumulative Proportion 0.62 1.00
##
## With factor correlations of
## ML1 ML2
## ML1 1.00 0.72
## ML2 0.72 1.00
##
## Mean item complexity = 1.4
## Test of the hypothesis that 2 factors are sufficient.
##
## The degrees of freedom for the null model are 55 and the objective function was 10.25 with Chi Square of 620.17
## The degrees of freedom for the model are 34 and the objective function was 1.24
##
## The root mean square of the residuals (RMSR) is 0.07
## The df corrected root mean square of the residuals is 0.09
##
## The harmonic number of observations is 66 with the empirical chi square 32.65 with prob < 0.53
## The total number of observations was 66 with MLE Chi Square = 73.39 with prob < 1e-04
##
## Tucker Lewis Index of factoring reliability = 0.884
## RMSEA index = 0.145 and the 90 % confidence intervals are 0.091 0.174
## BIC = -69.06
## Fit based upon off diagonal values = 0.98
## Measures of factor score adequacy
## ML1 ML2
## Correlation of scores with factors 1.00 0.97
## Multiple R square of scores with factors 0.99 0.94
## Minimum correlation of possible factor scores 0.99 0.87因子分析の結果の見方は次の通りです。
- 因子負荷量を見て，まとまりの意味を考えます。
- 累積寄与率もチェックしておきましょう
- 因子間相関も見ておきましょう
- その下には適合度指標が出ていますので，参考にしましょう
ちょっとわかりにくいところもあるので，因子負荷量が小さなところを表示しないようにしてみます。
print(fa.result,sort=T,cut=0.3)## Factor Analysis using method = ml
## Call: fa(r = fa.dat, nfactors = 2, fm = "ML")
## Standardized loadings (pattern matrix) based upon correlation matrix
## item ML1 ML2 h2 u2 com
## R 11 0.97 1.00 0.0049 1.0
## HR 6 0.93 0.70 0.3023 1.0
## lose 2 0.74 0.69 0.3074 1.1
## hits 5 0.57 0.49 0.97 0.0285 2.0
## save 3 -0.53 0.20 0.8000 1.1
## BB 7 0.52 0.38 0.70 0.3042 1.8
## hold 4 -0.43 0.34 0.6604 1.4
## DB 8 0.36 0.34 0.6600 1.8
## W.pitch 10 0.19 0.8058 2.0
## K 9 0.86 0.71 0.2932 1.0
## win 1 0.83 0.76 0.2389 1.0
##
## ML1 ML2
## SS loadings 4.10 2.49
## Proportion Var 0.37 0.23
## Cumulative Var 0.37 0.60
## Proportion Explained 0.62 0.38
## Cumulative Proportion 0.62 1.00
##
## With factor correlations of
## ML1 ML2
## ML1 1.00 0.72
## ML2 0.72 1.00
##
## Mean item complexity = 1.4
## Test of the hypothesis that 2 factors are sufficient.
##
## The degrees of freedom for the null model are 55 and the objective function was 10.25 with Chi Square of 620.17
## The degrees of freedom for the model are 34 and the objective function was 1.24
##
## The root mean square of the residuals (RMSR) is 0.07
## The df corrected root mean square of the residuals is 0.09
##
## The harmonic number of observations is 66 with the empirical chi square 32.65 with prob < 0.53
## The total number of observations was 66 with MLE Chi Square = 73.39 with prob < 1e-04
##
## Tucker Lewis Index of factoring reliability = 0.884
## RMSEA index = 0.145 and the 90 % confidence intervals are 0.091 0.174
## BIC = -69.06
## Fit based upon off diagonal values = 0.98
## Measures of factor score adequacy
## ML1 ML2
## Correlation of scores with factors 1.00 0.97
## Multiple R square of scores with factors 0.99 0.94
## Minimum correlation of possible factor scores 0.99 0.87これを見て因子の意味を考え，因子に命名することが一般的です。
落穂拾い
因子分析の留意点を幾つか挙げておきます。
- 推定方法，今回はMLを使いましたが他にもあります。デフォルトではminres法です。
- 回転方法，今回はデフォルトのoblimin回転を使っていますが，回転方には他にもあります。
- 直交回転はバリマックスが有名です。
- 斜交回転はプロマックス法，独立クラスター法，ジオミン法などが有名です。
- GPArotationパッケージを使うと他にもいろいろな回転ができます。
本日の課題
- バッターのデータをつかって因子分析を行ってください。使う変数は，次の通りです。
- 身長
- 体重
- 打数(AB)
- 安打(Hit)
- 得点(R)
- 二塁打(Double)
- 三塁打(Three)
- 本塁打(HR)
- 盗塁(steal)
- 三振(K)
- 併殺打(DP)
- 四球(BB)
- 因子数を決めた基準と回転法を明記してください。因子分析の結果の解釈，命名も行ってください。