5.2 主成分分析
5.3.1 因子分析 factanal()
共分散行列もしくはデータ行列に対し，最尤法による因子分析
(factor analysis)
を実行 する．書式：
factanal(x, factors, data = NULL, covmat = NULL, n.obs = NA, subset, na.action, start = NULL,
scores = c("none", "regression", "Bartlett"), rotation = "varimax", control = NULL, ...)
引数：
x
モデル式，数値行列，または数値行列に変換できるオブジェクトfactors
当てはめられる因子の数data
オプションのデータフレーム(
もしくは類似物，model.frame()
参照)
でx
が 公式のときだけ使われる．既定では変数はenvironment(formula)
から取ら れるcovmat
共分散行列か，cov.wt()
が返すような共分散のリスト．相関行列でも良いn.obs
観測値数で，covmat
が共分散行列のとき使われるsubset
もしx
が行列か公式の時，使われる部分データセットを指定するna.action
もしx
が公式の時，欠損値の処理法を与えるstart NULL
か，初期値の行列で，各列が独自性(uniqueness)
の初期値を与えるscores
もしあれば，出力スコアのタイプ．既定では無し．"regression"
はThomp-son
スコアを与える．"Bartlett"
はBartlett
の重みつき最小2
乗スコアを与 える．文字列の一部を与えるだけで良いrotation
文字列．"none"
か，因子を回転するために使われる関数の名前．この関 数は最初の引き数が負荷行列として呼び出され，回転された負荷を与える成分loadings
を持つ引き数か，もしくは単に回転された負荷そのものを返さなければならない
control
次の制御変数のリスト：nstart start=NULL
の時，試みられる初期値の数．既定では1
trace
論理値．実行過程の情報を出力するか？既定ではFALSE
lower
最適化実行中の独自性に対する下限で正値(
既定値0.005)
opt
関数optim()
に引き渡される制御変数のリストrotate
回転関数のための追加引き数リスト... control
の成分をfactanal()
への名前付き変数として引き渡すことができる 返り値： クラス"factanal"
のリストで，以下の成分を持つ：loadings
因子負荷量の行列で，一つの列が各因子に対応する．因子は負荷量の2
乗和に関して降順に並べられ，負荷量の和が正になるように符号を付けられる
uniquenesses
計算された独自性correlation
使われた相関行列criteria
最適化の結果．負対数尤度の値と繰り返しに関する情報factors
引き数factors
dof
因子分析モデルに対する自由度の数method
メソッド，常に"mle"
scores
もし必要とされたなら，得点の行列n.obs
もし得られるなら観測値の数，さもなければNA
call
マッチした呼び出し式na.action
もし関係するならSTATISTIC,PVAL
もし計算可能なら，有意性検定統計量とp
値因子分析モデルは，
n
個体に関するp
個の観測値ベクトルを表すp × n
観測値行列x
，p × k
因子負荷(loading)
行列Λ
，p × k
個の因子得点(score)
行列f
，そしてp × n
個の 誤差行列(
独自因子) ε
からなる次のようなモデルである：x = Λf + ε
x
以外全て観測されないが，主な制約として，因子得点は無相関で単位分散，p
個の誤差 ベクトルは互いに独立で分散ベクトル(
独自性，uniquness) Φ
を持つ，と仮定する．この ように，因子分析は本質的にx
の共分散行列Σ = Λ
tΛ + Ψ
に対するモデルである
(Ψ
はΦ
を対角成分に持つ対角行列)
．もしΛ
をp × p
直交行列G
を用いてGΛ
と変更してもモデルは変わらないから，依然として未確定さが残る．こ うした行列G
は回転(rotation)
と呼ばれる(
この言葉は必ずしも直交でない正則行列に 対しても使われることがある)
．もし
covmat
が与えられれば，それが使われる．さもなければ，もしそれが行列であれば，
x
から共分散行列が作られる．もしくは，モデル式x
が「データ」に適応され，それ から共分散行列が作られる．モデル式における目的変数は意味が無い，全ての変数は数値 である必要がある．共分散が与えられるか，x
から計算されると，それは解析のために相 関行列に変換される．この相関行列は結果のcorrelation
成分として返される．当てはめは，独自因子が多変量正規分布に従うとして，対数尤度を最適化することによ り行われる．与えられた独自因子に対する因子負荷の最大化は解析的に求めることがで きる
(Lawley & Maxwell)
．start
で与えられた全ての初期値が順に試めされ，得られた 最良のものが使われる．もしstart = NULL
ならば，最初の当てはめはJl‘eskog
が提案し，
Lawley & Maxwell
で与えられている初期値を用い，それから，独自因子が全て等しいとしてランダムに選んだ
control$nstart-1
個の他の値が試される．独自性は技術的に区間
[0, 1]
にあるとされるが，0
に近い値は問題を引き起こしやすい ので，最適化は既定値が0.005
である 下限control$lower
を用いて行われる(Lawley
& Maxwell)
．因子得点はデータ行列が与えられ，そして使われたときだけ計算できる．最初の方法は
Thomson
の回帰法であり，二つ目はBartlett
の重みつき最小自乗法である．両者は未知の得点
f
の推定値である．Thomson
の方法は未知のf
を(
母集団で) s
に回帰しf b = Λ
tΣ
−1x
とし，それから右辺の量の標本推定値を代入する．
Bartlett
の方法は標準誤差の2
乗和 を，(
当てはめられた) Λ
を与えた上で，f
に付いて最小化する．もし
x
がモデル式なら，標準的な欠損値処理が得点(
必要とされるなら)
に対して適用 される，napredict()
を参照せよ．注意：因子分析には無数の変種があり，他のプログラムの出力と比較するのは困難であ る．更に，因子分析における最大尤度の計算は困難であり，我々が比較した多くの他の例 は，この関数による当てはめより劣る結果を示した．特に，
Heywood
のケース(
一つも しくは複数の独自性が本質的に0)
である解は，多くのテキストや他のプログラムがそう 思わせるよりも，はるかに普通に起こる．関 連：
print.loadings(), varimax(),princomp(),ability.cov, Harman23.cor, Harman74.cor.’
# テストデータベクトル
> v1 <- c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 3, 3, 3, 3, 3, 4, 5, 6)
> v2 <- c(1, 2, 1, 1, 1, 1, 2, 1, 2, 1, 3, 4, 3, 3, 3, 4, 6, 5)
> v3 <- c(3, 3, 3, 3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 5, 4, 6)
> v4 <- c(3, 3, 4, 3, 3, 1, 1, 2, 1, 1, 1, 1, 2, 1, 1, 5, 6, 4)
> v5 <- c(1, 1, 1, 1, 1, 3, 3, 3, 3, 3, 1, 1, 1, 1, 1, 6, 4, 5)
> v6 <- c(1, 1, 1, 2, 1, 3, 3, 3, 4, 3, 1, 1, 1, 2, 1, 6, 5, 4)
> m1 <- cbind(v1, v2, v3, v4, v5, v6) # 列ベクトルとして行列に結合
> cor(m1) # 相関行列
v1 v2 v3 v4 v5 v6
v1 1.0000000 0.9393083 0.5128866 0.4320310 0.4664948 0.4086076 v2 0.9393083 1.0000000 0.4124441 0.4084281 0.4363925 0.4326113 v3 0.5128866 0.4124441 1.0000000 0.8770750 0.5128866 0.4320310 v4 0.4320310 0.4084281 0.8770750 1.0000000 0.4320310 0.4323259 v5 0.4664948 0.4363925 0.5128866 0.4320310 1.0000000 0.9473451 v6 0.4086076 0.4326113 0.4320310 0.4323259 0.9473451 1.0000000
> factanal(m1, factors = 3) # 因子数3で因子分析実行
Call:
factanal(x = m1, factors = 3)
Uniquenesses: # 独自性
v1 v2 v3 v4 v5 v6
0.005 0.101 0.005 0.224 0.084 0.005
Loadings: # 因子負荷量
Factor1 Factor2 Factor3
v1 0.944 0.182 0.267 v2 0.905 0.235 0.159 v3 0.236 0.210 0.946 v4 0.180 0.242 0.828 v5 0.242 0.881 0.286 v6 0.193 0.959 0.196
Factor1 Factor2 Factor3 SS loadings 1.893 1.886 1.797 Proportion Var 0.316 0.314 0.300 Cumulative Var 0.316 0.630 0.929
# 因子モデルの自由度と，負の最大尤度値
The degrees of freedom for the model is 0 and the fit was 0.4755
> factanal(m1, factors = 3, rotation = "promax") # promax回転を指定して再解析 Call:
factanal(x = m1, factors = 3, rotation = "promax")
Uniquenesses: # 独自性(同じ値)
v1 v2 v3 v4 v5 v6
0.005 0.101 0.005 0.224 0.084 0.005
Loadings: # 負荷量(大きく異なる)
Factor1 Factor2 Factor3
v1 0.985
v2 0.951
v3 1.003
v4 0.867
v5 0.910 v6 1.033
Factor1 Factor2 Factor3 # ほぼ同じ値
SS loadings 1.903 1.876 1.772 Proportion Var 0.317 0.313 0.295 Cumulative Var 0.317 0.630 0.925
The degrees of freedom for the model is 0 and the fit was 0.4755
# モデル式を用いた因子分析(目的変数は指定する必要は無い)．得点だけ出力
> factanal(~v1+v2+v3+v4+v5+v6,factors=3,scores="Bartlett")$scores Factor1 Factor2 Factor3
1 -0.9039949 -0.9308984 0.9475392 2 -0.8685952 -0.9328721 0.9352330 (途中省略)
18 1.8822320 0.3086244 1.9547752
# 知的能力検査結果データability.covを使用
> ability.FA <- factanal(factors = 1, covmat = ability.cov)
# 検定結果(1因子帰無仮説は棄却される)
Test of the hypothesis that 1 factor is sufficient.
The chi square statistic is 75.18 on 9 degrees of freedom.
The p-value is 1.46e-12
> update(ability.FA, factors = 2) (途中省略)
# 検定結果(2因子帰無仮説は棄却されない)
Test of the hypothesis that 2 factors are sufficient.
The chi square statistic is 6.11 on 4 degrees of freedom.
The p-value is 0.191
# promax回転を指定して再実行(update関数の使用に注意)
> update(ability.FA, factors = 2, rotation = "promax") (途中省略)
# 検定結果(やはり2因子帰無仮説は棄却されない)
Test of the hypothesis that 2 factors are sufficient.
The chi square statistic is 6.11 on 4 degrees of freedom.
The p-value is 0.191