いきなり因子分析（その5）：高次因子分析モデルと双因子モデル）

因子分析についてのシリーズものの記事。

その1（とりあえず試した編）
その2（因子数の決定編）
その3（様々な推定法編）
その4（因子軸の回転編）

今回はもう少し複雑なモデルということで高次因子分析モデルと階層因子分析モデルのうち双因子モデルを試してみる。

高次因子分析って？

高次因子分析というのは、複数の因子をさらに少数の因子で説明するモデルである。たとえば、因子分析の結果因子間に強い相関が仮定される場合に、それらに影響を与える上位の因子を想定するモデルである。

こうした場合に、観測を変数を説明する因子（通常の因子分析で想定する因子）を1次因子と言い、それらの因子を説明する上位の因子を2次因子というそうだ。

パス図で表すとこんな感じ。

f:id:nekomosyakushimo:20180714232805p:plain

双因子モデルって？

双因子モデルというのは、階層因子分析の一部とされる。階層因子分析というのは、一般因子とグループ因子で階層を分けて観測変数を説明するモデルである。一般因子は、全ての観測変数に影響を与えている因子で、グループ因子は通常の因子分析と同じように、一部の観測辺境に影響を与えるモデルらしい。

パス図で表すとこんな感じ。

f:id:nekomosyakushimo:20180714232832p:plain

Rで高次因子分析を実践

psychを使って、高次因子分析モデルや双因子モデルで因子分析する方法は次のpdfを参考にした。 http://personality-project.org/r/psych/HowTo/factor.pdf

データは前の記事と変わらずHolzingerSwineford1939を使うことにする。

どちらのモデルもpsychに含まれるomega()という関数を使ってできるらしい。

ちなみに日本語でこの関数についてググると大体がω係数の算出に関連したものがヒットする。ω係数とは、信頼性の指標であるαと同じように、内的一貫性を表す指標だとされている。因子負荷と誤差分散を用いて算出される値らしく、この関数を使うとそれが算出されるのでそうした記事がヒットするようだ。

それはさておき、さっそく使ってみる。

library(psych)
model_higherorder <- omega(m=dat[7:15],
                           nfactors = 3,
                           sl=FALSE
                           )

引数のm=でデータの指定。ローデータでも相関行列でもどちらでも良いっぽい。

nfactorsには因子数を指定する。省略すると、何だかよく分からないメカニズムで因子数を推定して結果を返すけれどもきっと指定した方が良いでしょう。（因子数の決め方についてはその2 を参照）

slというのは、Schmid Leiman変換というもので、これをTRUEにすると双因子モデルになり、FALSEだと高次因子分析モデルになる。

他にもfmで、前記事で扱った推定法を選べたりもする。（デフォルトは"minres"なようだ）

さて、結果を見てみる。summary(model_higherorder)を打つと結果が表示される。

> summary(model_higherorder)
Omega
Alpha:                 0.76
G.6:                   0.81
Omega Hierarchical:    0.45
Omega H asymptotic:    0.53
Omega Total            0.85

まず最初のブロックはω係数とか信頼係数の値が表示されている。Aplhaはクロンバックのα係数、G.6はガットマンの提案した信頼性係数の6つ目の下界である（詳しく知りたい人は岡田先生のこの論文を読むと良い）。

そこから3つがω係数について。調べた感じだとOmega Hierarchicalは一般因子の因子負荷係数を基に算出したもので、Omega Totalの方は一般因子だけでなく、グループ因子も足して計算したものらしい。omega H asymptoticというのは、計算式をみる限りだと全ての因子負荷に占める一般因子部分の割合みたいに読めるのだけれどちょっとどう解釈して良いかは分からない。信頼性係数として報告する際に使うのはOmega Totalの値である。

With eigenvalues of:
   g  F1*  F2*  F3*
1.44 1.62 0.77 1.03

次には一般因子とグループ因子の固有値が出てくる。まぁこれはそのまま固有値ですね。

The degrees of freedom for the model is 12  and the fit was  0.08
The number of observations was  301  with Chi Square =  22.56  with prob <  0.03

The root mean square of the residuals is  0.02
The df corrected root mean square of the residuals is  0.05

RMSEA and the  0.9 confidence intervals are  0.055 0.016 0.088
BIC =  -45.93Explained Common Variance of the general factor =  0.3

そこから先は、モデルの評価に関するあれやこれやの指標。普通に因子分析をしたときと同様のものが出てくる。一つ違うのは、一般因子のECVというものが出力されていることでしょうか。これはomega()のhelpに次のように書いてあった。

Following suggestions by Steve Reise, the Explained Common Variance (ECV) is also reported. This is the ratio of the general factor eigen value to the sum of all of the eigen values. As such, it is a better indicator of unidimensionality than of the amount of test variance accounted for by a general factor.

Total, General and Subset omega for each subset
                                                g  F1*  F2*  F3*
Omega total for total scores and subscales    0.85 0.89 0.65 0.72
Omega general for total scores and subscales  0.45 0.24 0.27 0.19
Omega group for total scores and subscales    0.35 0.65 0.37 0.53

最後には、それぞれの因子が合計得点および下位尺度のばらつきをどの程度説明されるかが表示される。１行目はOmega totalとは、一般因子Omega generalと（2行目）とグループ因子Omega group（３行目）を足したものである。一番左の列は、結果の最初に表示されたω係数と同じである。右の3列のサブセットというのがそれぞれの一次因子が負荷していると想定される観測変数のω係数ということだろうか（よく分からん）。

summaryだと、これしか出ないのでが、高次の因子から下位の因子への負荷ぐらい知っておきたいだろう。そういう場合は、次のようにするとその値を取り出せる。

print(model_higherorder$schmid$gloading, digits=2)

Loadings:
   MR1  
F1 0.514
F2 0.628
F3 0.415

                 MR1
SS loadings    0.831
Proportion Var 0.277

双因子モデルを実践

今度は、双因子モデルを試してみる。といっても、引数sl=TRUEと指定するのみだけれど。

model_higherorder <- omega(m=dat[7:15],
                           nfactors = 3,
                           sl=TRUE
)

summary(model_bifactor)

結果は高次因子分析モデルと同様に出てくる。ただ、一般因子から観測因子への負荷が知りたいこともあるでしょう。その場合には、print(model_bifactor$schmid$sl, digit = 2)とすると値が取り出せる。

> print(model_bifactor$schmid$sl, digit = 2)
      g    F1*    F2*     F3*   h2   u2   p2
x1 0.49  0.169  0.460  0.0282 0.48 0.52 0.49
x2 0.29  0.037  0.396 -0.1113 0.26 0.74 0.33
x3 0.41 -0.053  0.534  0.0173 0.45 0.55 0.37
x4 0.45  0.726  0.013  0.0070 0.73 0.27 0.28
x5 0.42  0.760 -0.050  0.0061 0.75 0.25 0.23
x6 0.46  0.691  0.062 -0.0118 0.69 0.31 0.30
x7 0.23  0.038 -0.119  0.6701 0.52 0.48 0.10
x8 0.35 -0.029  0.097  0.6246 0.52 0.48 0.23
x9 0.45  0.027  0.297  0.4152 0.46 0.54 0.43