5.2K Views
July 25, 23
スライド概要
神戸大学経営学研究科で2022年度より担当している「統計的方法論特殊研究(多変量解析)」の講義資料「06_因子分析」です。
神戸大学経営学研究科准教授 分寺杏介(ぶんじ・きょうすけ)です。 主に心理学的な測定・教育測定に関する研究を行っています。 講義資料や学会発表のスライドを公開していきます。 ※スライドに誤りを見つけた方は,炎上させずにこっそりお伝えいただけると幸いです。
Chapter 6 因子分析 Chapter 6 因子分析 本資料は,神戸大学経営学研究科で 2022 年度より担当している「統計的方法論特殊研究(多 変量解析)」の講義資料です。CC BY-NC 4.0 ライセンスの下に提供されています。 作成者連絡先 神戸大学大学院経営学研究科 分寺杏介(ぶんじ・きょうすけ) mail: [email protected] website: https://www2.kobe-u.ac.jp/~bunji/ おまたせしました。今回からいよいよ本編です。まずは因子分析 (factor analysis) です。は じめに因子分析の基本的&数理的な考え方を紹介した後で,R でやってみます。その後,因 子分析を実施する上でのプラクティカルなポイントを抑えていく予定です。 Contents 6.1 6.2 6.3 6.4 6.5 6.6 6.7 6.8 6.9 6.10 6.11 因子分析とは . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 因子分析モデルの一般化 . . . . . . . . . . . . . . . . . . . . . . . . . . パラメータの推定方法 . . . . . . . . . . . . . . . . . . . . . . . . . . . とりあえずやってみよう . . . . . . . . . . . . . . . . . . . . . . . . . . カテゴリカル因子分析 因子分析をする前に 回転の方法 . . . . . 因子の推定法 . . . . 因子数の決め方 . . . 𝜔 係数 . . . . . . . 主成分分析との違い . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 読み込み 1 1 dat <- readRDS("chapter04.rds") . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 6 7 15 28 30 34 37 40 47 55
6.1 因子分析とは 6.1 因子分析とは 今回使用している dat の Q1_11 から Q1_15 はいずれも外向性 (extraversion) を測定してい る項目です。実際に項目を見てみると*1 11. Don’t talk a lot. (【逆転】あまり喋らない。) 12. Find it difficult to approach others. (【逆転】人に近づきづらいと思う。) 13. Know how to captivate people. (人を魅了できる。) 14. Make friends easily. (友達を作りやすい。) 15. Take charge. (リーダータイプである。) となっており,たしかに他人との関わりに関する内容を尋ねているようです。改めて項目間 の相関を確認してみると,逆転項目も処理済みなので全ての項目間に正の相関があります。 この相関構造を一つの指標にまとめたものの一つが 𝛼 係数でした。 項目間相関 1 # round()関数を使って表示桁数を減らしています 2 round(cor(dat[, paste0("Q1_", 11:15)]), 2) 1 Q1_11 Q1_12 Q1_13 Q1_14 Q1_15 2 Q1_11 1.00 0.47 0.32 0.42 0.31 3 Q1_12 0.47 1.00 0.39 0.53 0.39 4 Q1_13 0.32 0.39 1.00 0.42 0.39 5 Q1_14 0.42 0.53 0.42 1.00 0.32 6 Q1_15 0.31 0.39 0.39 0.32 1.00 因子分析では,こうした一連の項目に関して共通する潜在的な特性値によって回答が決まっ ていると考えます。これをグラフィカルモデルで表すと以下のようになります。 図6.1 因子分析モデル 1 グラフィカルモデルでは,潜在変数は楕円で表すのが一般的です。図6.1に関して,項目ごと に見るとこれは回帰分析モデルのように表すことができそうです。例えば項目 Q1_11 は,説 明変数を「外向性」𝑓𝑝 として,対応する回帰係数を 𝑏11 (添字は”11” の一つだけ)とすると 𝑦𝑝,11 = 𝑓𝑝 𝑏11 + 𝜀𝑝,11 *1 2 日本語訳はこのサイトから拝借しました。 (6.1)
6.1 因子分析とは と書くことが出来ます。この表現で言えば,回帰分析の説明変数に当たるものが 𝑥𝑝 から 𝑓𝑝 に変わっているだけです。同様にして,他の項目も 𝑦𝑝,12 = 𝑓𝑝 𝑏12 + 𝜀𝑝,12 𝑦𝑝,13 = 𝑓𝑝 𝑏13 + 𝜀𝑝,13 (6.2) 𝑦𝑝,14 = 𝑓𝑝 𝑏14 + 𝜀𝑝,14 𝑦𝑝,15 = 𝑓𝑝 𝑏15 + 𝜀𝑝,15 と表すことが出来ます。つまり,因子分析は項目レベルで見ると説明変数が潜在変数になっ た回帰分析なのです。なお,因子分析では本来回帰分析と同じように結果変数 𝑦 は連続量 (間隔尺度か比率尺度)を考えています。実際に心理尺度で得られる回答(項目レベル)はカ テゴリカル変数であり間隔尺度である保証はないのですが,とりあえず今はポリシリアル相 関のときと同じように 𝑦 は回答の背後の潜在変数だと考えるか,間隔尺度だと仮定されたも のと考えておいてください。 ここまでは実際のデータに合わせて説明をしてきましたが,以後の説明ではもう少し一般化 して説明していきます。ということで,ここで今後登場する添字を説明しておきます。なお, 以下の記号について小文字は一つの値を,大文字は総数を表しているものとします。 𝑝 回答者の番号 (person) を表します。 (𝑝 = 1, 2, ⋯ , 𝑃 − 1, 𝑃 ) 𝑖 項目の番号 (item) を表します。 (𝑖 = 1, 2, ⋯ , 𝐼 − 1, 𝐼 ) (6.1-6.2) 式では,項目 𝑖 は 11, 12, 13, 14, 15 と対応しています。この表記を用いると,回 答者 𝑝 の項目 𝑖 に対する回答は以下のように表されます。 (6.3) 𝑦𝑝𝑖 = 𝑓𝑝 𝑏𝑖 + 𝜀𝑝𝑖 ただ,いちいち 𝑝 を書いているとちょっと面倒なので,Chapter5 のときと同じように,回 答者の要素についてはベクトル表記にしておきます。𝐲𝑖 = [𝑦1𝑖 𝐟 = [𝑓1 ⊤ 𝑓2 ⋯ 𝑓𝑃 −1 𝑓𝑃 ] ,𝛆𝑖 = [𝜀1𝑖 ⊤ 𝑦2𝑖 ⋯ 𝑦𝑃 −1,𝑖 𝑦𝑃 𝑖 ] , ⊤ 𝜀2𝑖 ⋯ 𝜀𝑃 −1,𝑖 𝑦1𝑖 𝑓1 𝜀1𝑖 ⎡𝜀 ⎤ ⎡𝑦 ⎤ ⎡𝑓 ⎤ ⎢ 2𝑖 ⎥ = ⎢ 2 ⎥ 𝑏𝑖 + ⎢ 2𝑖 ⎥ ⎢ ⋮ ⎥ ⎢ ⋮ ⎥ ⎢ ⋮ ⎥ 𝑦 𝑓 ⎣𝜀𝑃 𝑖 ⎦ ⎣ 𝑃 𝑖⎦ ⎣ 𝑃 ⎦ ⟶ 𝐲 𝑖 = 𝐟 𝑏 𝑖 + 𝛆𝑖 𝜀𝑃 𝑖 ] として, (6.4) としておきます。 因子分析では,回帰分析と違って説明変数(𝐟 )が直接観測できない潜在変数 (latent variable) になっています。説明変数が潜在変数だと,ちょっと困ったことが起こります。それは変数の スケールの問題です。いま,因子 𝐟 には身長や体重のときのような「単位」というものがあり ません。例えば,100 点満点のテストの点数を 500 点満点に換算して扱うように,𝐟 の値を 5 倍したものを 𝐟 ′ とします。さらにこの得点にゲタを履かせる要領で,𝐟 ″ = 𝐟 ′ + 10 = 5𝐟 + 10 と置いてみます。これを使うと,先程の回帰式は 𝐲𝑖 = 𝐟 𝑏𝑖 + 𝛆𝑖 = (5𝐟 + 10)(0.2𝑏𝑖 ) + (𝛆𝑖 − 2𝑏𝑖 ) = 3 𝐟 ″ 𝑏𝑖′ + 𝛆′𝑖 (6.5)
6.1 因子分析とは と変形が可能です。つまり 𝐟 のスケールが変わることで 𝑏𝑖′ = 0.2𝑏𝑖 および 𝛆′𝑖 = 𝛆𝑖 − 2𝑏𝑖 とい う新しい解が得られてしまいます。このままでは回帰係数を一意に求めることが出来ないの で,通常因子分析モデルでは,因子の分散 𝜎𝐟2 = 1 かつ平均 0 という制約を置いています*2 。 これに加えて,通常因子分析で考える潜在変数は正規分布に従う連続量であると仮定されま す。これらをまとめると,因子分析では 𝐟 ∼ 𝑁 (0, 1) という仮定が置かれていることになり ます。 そしてこの問題は観測変数 𝑦 の側についても,心理尺度のように値のスケールに単位が存在 しない限り同様のことがいえます。これも 100 点満点のテストの点数を 500 点満点に換算し て扱うように,5 件法のリッカート尺度の得点を 5, 10, 15, 20, 25 と 5 刻みにしても問題は 無いはずです。ということで 𝐲′𝑖 = 5𝐲𝑖 として回帰式を変形させてみると, 𝐲′𝑖 = 5𝐲𝑖 = 5𝐟 𝑏𝑖 + 5𝛆𝑖 (6.6) = 𝐟 × 5𝑏𝑖 + 5𝛆𝑖 = 𝐟 𝑏𝑖″ + 𝛆′𝑖 ということで,𝑏𝑖″ = 5𝑏𝑖 および 𝛆𝑖 = 5𝛆𝑖 という新しい解が得られてしまいます。このよう に,因子分析の解は回答側のスケールの影響も受けてしまうため,通常は各項目の回答は全 て標準化された値を使用します。説明変数(𝐟 )も結果変数(𝐲𝑖 )も標準化されているという 意味では,因子負荷は標準偏回帰係数のように解釈することができる,というわけです。 因子分析モデルでは,各個人の特性値の強さ (𝑓𝑝 ) を因子得点 (factor score),回帰係数 (𝑏𝑖 ) を因子負荷 (factor loading) と呼びます。因子負荷が高い項目ほど,その人の回答が因子得 点の高低と強く関係しています。これは,合計点を X 軸に,グループごとの平均値を Y 軸 に置いたトレースラインの右上がり度の強さと近いものです。図6.2でいえば,赤い項目の方 が強い右上がりになっているため,青い項目よりも因子負荷が高い項目だと考えることが出 4 3 1 2 項目平均値 5 6 来ます。 10 15 20 25 30 グループ平均値 図6.2 トレースライン *2 4 因子分析モデルの制約は,正確に言うと「パラメータが一つ固定されていれば良い」ので,別の制約として第 1 項目の因子負荷を 1 に固定するという方法を取ることもあります。
6.1 因子分析とは なおグラフィカルモデルでは図6.3のように,因子負荷の強さを矢印の横に書くのが一般的 です。 図6.3 因子分析モデル 2(因子負荷を追加) また,因子分析では誤差 𝜀 に関しても「因子」として取り扱います。複数の項目に共通する 成分(𝐟 )を共通因子 (common factor),各項目に特有の成分(あるいは共通因子で説明さ れなかった残り)を独自因子 (unique factor) と呼び,各項目への回答はこの 2 種類の因子 の和のみによって表されていると考えるわけです。したがって,先程示した因子分析モデル の基本の式(6.1)式における 𝛆𝑖 は,実際には「誤差」ではなく「もう一つの説明変数」であ ると考えるわけです(とはいえ実用上は誤差と同じように考えても差し支えないのですが)。 ということで,(6.1)式は厳密には以下のような形をしています。 error 𝐲𝑖 = 𝐟 𝑏𝑖 + 𝛆𝑖 𝑢𝑖 ⏞ +0 (6.7) ここで 𝑢𝑖 は,独自因子に対する因子負荷を表しています。このように,因子分析モデルでは 回帰分析で言うところの「誤差」は 0 と考えることが多いようです。 ということでグラフィカルモデルにもこの独自因子を書き込むと,図6.4のようになるわけで す。ただし独自因子の因子負荷 𝑢𝑖 を自由に推定しようとすると,独自因子得点 𝛆𝑖 との間に またスケールの問題が生じてしまい解が一つに決まらないため,通常は 𝑢𝑖 = 1 に固定され ます。 図6.4 因子分析モデル 3(誤差を追加) 実際に論文などに載せる場合,独自因子の因子負荷はどうせ固定されるのだから,矢印の横 の「1」は省略されることが多いです。更にいうと,因子分析ではそもそも共通因子にのみ関 心があることが圧倒的に多いので,独自因子自体を省略する形(つまり図6.3)になることも 多いです。 5
6.2 因子分析モデルの一般化 6.2 因子分析モデルの一般化 回帰分析的に因子分析を見ると,独自因子は誤差のようなものなのでなるべくその影響は少 ないほうが良い,といえます。極端な例を考えてみましょう。因子負荷 𝑏𝑖 = 0 の項目がある とすると,その項目の因子分析モデル式は 𝐲𝑖 = 𝐟 × 0 + 𝛆 𝑖 = 𝛆𝑖 (6.8) となります。共通因子の影響を全く受けること無く独自因子のみで回答が決定している状態 です。因子得点 𝑓𝑝 を推定する上では何の情報も持ち合わせていないこの項目は,あってもな くても同じです。この項目に高い点数をつけていても低い点数をつけていても,𝑓𝑝 には何も 関与していないのです。ということで,基本的には因子分析では共通因子の影響が多い方が 望ましいわけです。 回帰分析では,𝑦 に占める誤差分散の割合を減らし予測の精度を上げるために,複数の説明 変数を用いることができました。因子分析でも同様に,複数の説明変数=共通因子がある状 態を考えることができます。ここでは,5 項目が 2 つの共通因子の影響を受けている,と考 えてみましょう。グラフィカルモデルでは図6.5のような状態です。なお,因子番号を表す添 字は trait の頭文字 𝑡 (𝑡 = 1, 2) としておきます。 図6.5 2 因子モデル 当然ながら,因子 1 の項目 1 に対する因子負荷 𝑏11(添字は 2 つ:1 と 1)と因子 2 の因子負 荷 𝑏21 は異なります。したがって,𝐲𝑖 に関する因子分析モデルは 𝐲𝑖 = 𝐟1 𝑏1𝑖 + 𝐟2 𝑏2𝑖 + 𝛆𝑖 (6.9) となります。 もちろん因子の数はいくつでも良いので,更に一般化して 𝑇 個の因子がある場合を考えると, 𝐲 に関する因子分析モデルは 𝐲𝑖 = 𝐟1 𝑏1𝑖 + 𝐟2 𝑏2𝑖 + ⋯ + 𝐟𝑇 𝑏𝑇 𝑖 + 𝛆𝑖 = [𝐟1 𝐟2 = 𝐅𝐛𝑖 + 𝛆𝑖 と書くことが出来ます。 6 𝑏1𝑖 ⎡𝑏 ⎤ ⋯ 𝐟𝑇 ] ⎢ 2𝑖 ⎥ + 𝛆𝑖 ⎢ ⋮ ⎥ ⎣𝑏𝑇 𝑖 ⎦ (6.10)
6.3 パラメータの推定方法 6.3 パラメータの推定方法 実際に因子分析のパラメータを推定する際には,観測変数の相関行列 𝚺 をもとにして計算し ます。まずは因子数が 1 の場合でその仕組みを確認してみましょう。 6.3.1 1 因子の場合 相関行列 𝚺 と各項目の因子負荷 (𝑏1 , 𝑏2 , ⋯, 𝑏𝐼 ) の関係性を理解するために,同じ因子 の影響を受ける 2 つの変数の関係を見てみます。ここまでに見てきた因子分析モデルでは, 「全ての項目が共通因子の影響を受けている」という一因子構造を仮定しています。理想的な 一因子構造では,共通成分は共通因子にすべて吸収された結果,残った独自因子同士の相関 はゼロであることが期待されています。各項目の回答は標準化されているため,変数間の共 分散と相関係数は一致します。いま,ある項目への回答 𝐲𝑖 と別の項目への回答 𝐲𝑗 (𝑖 ≠ 𝑗) の 共分散 𝜎𝐲𝑖 ,𝐲𝑗 および相関係数 𝑟𝐲𝑖 ,𝐲𝑗 は, 𝜎𝐲𝑖 ,𝐲𝑗 = 𝑟𝐲𝑖 ,𝐲𝑗 = 𝜎(𝑏𝑖 𝐟 +𝛆𝑖 ),(𝑏𝑗 𝐟 +𝛆𝑗 ) = 𝜎(𝑏𝑖 𝐟 ),(𝑏𝑗 𝐟 ) + 𝜎(𝑏𝑖 𝐟 ),(𝛆𝑗 ) + 𝜎(𝑏𝑗 𝐟 ),(𝛆𝑖 ) + 𝜎𝛆𝑖 ,𝛆𝑗 = 𝑏𝑖 𝑏𝑗 𝜎𝐟2 (6.11) + 𝜎(𝑏𝑖 𝐟 ),(𝛆𝑗 ) + 𝜎(𝑏𝑗 𝐟 ),(𝛆𝑖 ) + 𝜎𝛆𝑖 ,𝛆𝑗 となります。ここで,ある項目の共通因子成分と別の項目の独自因子成分は,通常の回帰分 析と同じように無関係(共分散ゼロ)と見なすことができるため, 𝜎(𝑏𝑖 𝐟 ),(𝛆𝑗 ) = 𝜎(𝑏𝑗 𝐟 ),(𝛆𝑖 ) = 0 (6.12) が成り立ちます。ということで,2 項目の相関係数は 𝑟𝐲𝑖 ,𝐲𝑗 = 𝑏𝑖 𝑏𝑗 𝜎𝐟2 + 𝜎𝛆𝑖 ,𝛆𝑗 となりました。 先程説明したように,因子の分散 𝜎𝐟2 は 1 に固定するという制約が置かれています。したが って,2 項目の相関係数は最終的に 𝑟𝐲𝑖 ,𝐲𝑗 = 𝑏𝑖 𝑏𝑗 + 𝜎𝛆𝑖 ,𝛆𝑗 (6.13) つまり「因子負荷の積」と「独自因子同士の共分散」の和によって表すことが出来ます。こ こで,独自因子同士の共分散 𝜎𝛆𝑖 ,𝛆𝑗 はゼロとは限りません。図6.4のモデルで言えば,共通因 子は「5 項目全てに共通する成分」ですが,もしかしたら図6.6のように「5 項目全てでは無 いけれど 2 項目だけに共通する成分」などがまだ残っているかもしれません。 図6.6 誤差共分散がゼロにならないケース 7
6.3 パラメータの推定方法 この状態では,𝜎𝐲𝑖 ,𝐲𝑗 のときと同じ要領で 𝜎𝛆𝑖 ,𝛆𝑗 も 2 つの独自因子に対する因子負荷の積 (𝑏𝜀1 𝑏𝜀2 ≠ 0)になってしまいます。2 項目の誤差の相関がゼロではない限り,少なからずこ のような「別の因子」が存在しているわけですが,この「別の因子」があるということは一因 子構造からは逸脱していることになります。因子分析モデルにおいて重要なポイントは,こ のような「別の因子」は存在していないという仮定を置いている,ということです。その仮 定のもとでは,2 項目間の相関は因子負荷の積になっているため,因子負荷が強い項目間の 相関は必然的に高くなります。 同様にして,項目 𝑖 の分散(=項目 𝑖 と項目 𝑖 の共分散)は 𝜎𝐲2 𝑖 = 𝜎𝐲𝑖 ,𝐲𝑖 = 𝑏𝑖2 + 𝜎𝛆2 𝑖 (6.14) となります。つまり因子負荷の 2 乗と誤差分散の和になっています。また 𝜎𝐲𝑖 ,𝐲𝑖 = 𝑟𝐲𝑖 ,𝐲𝑖 = 1 でもあります。回帰分析では,結果変数 𝑦 の分散のうち説明変数で説明可能な分散の割合を 決定係数 𝑅2 と呼びました。因子分析では,この決定係数にあたるもの,すなわち観測変数 𝐲𝑖 の分散 𝜎𝐲2 𝑖 = 1 のうち共通因子で説明可能な分散の割合を共通性 (communality) と呼びま す。また,誤差の部分も(独自)因子として考えるため,観測変数 𝐲𝑖 の分散 𝜎𝐲2 𝑖 = 1 のうち独 自因子で説明可能な分散(=共通因子では説明できない分散)の割合を独自性 (uniqueness) と呼びます。上の式では独自性は 𝜎𝛆2 𝑖 の部分になるため,共通性は𝑏𝑖2 ,つまり因子負荷の二 乗に合致します。そして共通性と独自性の和は必ず 1 になります。 以上をまとめると,観測変数の相関行列 𝚺 は以下のように因子分析のパラメータの関数に変 換することが出来ます。 𝜎𝑦2 𝜎𝑦1 ,𝑦2 𝜎𝑦1 ,𝑦3 𝜎𝑦1 ,𝑦4 ⎡𝜎 1 𝜎𝑦22 𝜎𝑦2 ,𝑦3 𝜎𝑦2 ,𝑦4 𝑦2 ,𝑦1 ⎢ 𝜎𝑦3 ,𝑦4 𝜎𝑦23 𝚺 = ⎢𝜎𝑦3 ,𝑦1 𝜎𝑦3 ,𝑦2 ⎢𝜎𝑦4 ,𝑦1 𝜎𝑦4 ,𝑦2 𝜎𝑦4 ,𝑦3 𝜎𝑦24 ⎣𝜎𝑦5 ,𝑦1 𝜎𝑦5 ,𝑦2 𝜎𝑦5 ,𝑦3 𝜎𝑦5 ,𝑦4 𝑏12 + 𝜎𝛆2 𝑏1 𝑏2 𝑏1 𝑏3 ⎡ 𝑏 𝑏 1 𝑏2 + 𝜎 2 𝑏2 𝑏3 𝛆2 2 ⎢ 2 1 𝑏3 𝑏2 𝑏32 + 𝜎𝛆2 3 = ⎢ 𝑏3 𝑏1 ⎢ 𝑏4 𝑏1 𝑏4 𝑏2 𝑏4 𝑏3 𝑏5 𝑏2 𝑏5 𝑏3 ⎣ 𝑏5 𝑏1 𝑏12 ⎡𝑏 𝑏 ⎢ 2 1 = ⎢𝑏3 𝑏1 ⎢𝑏4 𝑏1 ⎣𝑏5 𝑏1 𝑏1 𝑏2 𝑏22 𝑏3 𝑏2 𝑏4 𝑏2 𝑏5 𝑏2 𝑏1 ⎡𝑏 ⎤ ⎢ 2⎥ = ⎢𝑏3 ⎥ [𝑏1 𝑏2 ⎢𝑏4 ⎥ ⎣𝑏5 ⎦ = 𝐛⊤ 𝐛 + 𝚿 ここで登場した 𝐛 = [𝑏1 𝑏1 𝑏3 𝑏2 𝑏3 𝑏32 𝑏4 𝑏3 𝑏5 𝑏3 𝑏3 𝑏2 𝑏1 𝑏4 𝑏2 𝑏4 𝑏3 𝑏4 𝑏42 𝑏5 𝑏4 𝑏4 ⋯ 𝜎𝑦1 ,𝑦5 𝜎𝑦2 ,𝑦5 ⎤ ⎥ 𝜎𝑦3 ,𝑦5 ⎥ 𝜎𝑦4 ,𝑦5 ⎥ 𝜎𝑦25 ⎦ 𝑏1 𝑏4 𝑏2 𝑏4 𝑏3 𝑏4 𝑏42 + 𝜎𝛆2 4 𝑏5 𝑏4 2 𝜎𝛆1 𝑏1 𝑏5 ⎡ ⎤ 0 𝑏2 𝑏5 ⎥ ⎢ 𝑏3 𝑏5 ⎥ + ⎢ 0 𝑏4 𝑏5 ⎥ ⎢ 0 𝑏52 ⎦ ⎣ 0 𝑏1 𝑏5 𝑏2 𝑏5 ⎤ ⎥ 𝑏3 𝑏5 ⎥ 𝑏4 𝑏5 ⎥ 𝑏52 + 𝜎𝛆2 5 ⎦ 0 𝜎𝛆2 2 0 0 0 0 0 𝜎𝛆2 3 0 0 0 0 0 𝜎𝛆2 4 0 0 0 ⎤ ⎥ 0 ⎥ 0 ⎥ 𝜎𝛆2 5 ⎦ (6.15) 𝑏5 ] + 𝚿 𝑏𝐼 ] は(6.10)式の 𝐛 とは異なり,すべての項目に対する因 子負荷を横に並べたベクトルです。また 𝚿 は観測変数の誤差分散を対角成分にもつ対角行列 8
6.3 パラメータの推定方法 です。実際には誤差共分散行列なわけですが,一因子構造の仮定では,誤差共分散=非対角 成分はゼロになっているわけですね。 実際のところ,計算のイメージとしては次のような感じです。 1. 観測変数の相関行列 𝚺 の非対角成分について,𝐛⊤ 𝐛 がなるべく近い値になるような 因子負荷ベクトル 𝐛 を求める 2. 対角成分の値が 1 になるように 𝚿 で調整する こうして得られた因子負荷および誤差分散を固定した上で,𝐲 とのズレが最も小さくなるよ うに因子得点 𝐟 の値が計算されます。 6.3.2 多因子の場合 因子の数が増えても,基本的には同じような形で観測変数の相関行列 𝚺 を分解していくこと になります。ちょっと長くなりますが,2 因子モデルでの 2 項目の共分散=相関係数を展開 してみましょう。異なる項目の共通因子成分 𝑏𝑡𝑖 𝐟𝑡 と誤差成分 𝜀𝑗 が無相関,かつ誤差項目間 の相関も共通因子間の相関もゼロ,すなわち 𝜎(𝑏𝑡𝑖 𝐟𝑡 ),(𝛆𝑗 ) = 𝜎(𝑏𝑡𝑗 𝐟𝑡 ),(𝛆𝑖 ) = 0 𝜎(𝛆𝑖 ),(𝛆𝑗 ) = 0 (6.16) 𝜎(𝑏1𝑖 𝐟1 ),(𝑏2𝑗 𝐟2 ) = 𝜎(𝑏2𝑖 𝐟2 ),(𝑏1𝑗 𝐟1 ) = 0 だと仮定すると,項目 𝐲𝑖 , 𝐲𝑗 の共分散=相関行列 𝚺 の (𝑖, 𝑗) 成分は 𝜎𝐲𝑖 ,𝐲𝑗 =𝑟𝐲𝑖 ,𝐲𝑗 =𝜎(𝑏1𝑖 𝐟1 +𝑏2𝑖 𝐟2 +𝛆𝑖 ),(𝑏1𝑗 𝐟1 +𝑏2𝑗 𝐟2 +𝛆𝑗 ) =𝜎(𝑏1𝑖 𝐟1 ),(𝑏1𝑗 𝐟1 ) + 𝜎(𝑏1𝑖 𝐟1 ),(𝑏2𝑗 𝐟2 ) + 𝜎(𝑏1𝑖 𝐟1 ),(𝛆𝑗 ) + 𝜎(𝑏2𝑖 𝐟2 ),(𝑏1𝑗 𝐟1 ) + 𝜎(𝑏2𝑖 𝐟2 ),(𝑏2𝑗 𝐟2 ) + 𝜎(𝑏2𝑖 𝐟2 ),(𝛆𝑗 ) + 𝜎(𝛆𝑖 ),(𝑏1𝑗 𝐟1 ) + 𝜎(𝛆𝑖 ),(𝑏2𝑗 𝐟2 ) + 𝜎(𝛆𝑖 ),(𝛆𝑗 ) =𝑏1𝑖 𝑏1𝑗 𝜎𝐟21 + 𝑏2𝑖 𝑏2𝑗 𝜎𝐟22 (6.17) =𝑏1𝑖 𝑏1𝑗 + 𝑏2𝑖 𝑏2𝑗 = [𝑏1𝑖 𝑏2𝑖 ] [ 𝑏1𝑗 ] 𝑏2𝑗 =𝐛⊤ 𝑖 𝐛𝑗 となります。同様にして,項目 𝐲𝑖 の分散= 𝚺 の (𝑖, 𝑖) 成分は 2 2 𝜎𝐲𝑖 ,𝐲𝑖 = 𝑏1𝑖 + 𝑏2𝑖 + 𝜎𝛆2 𝑖 = [𝑏1𝑖 𝑏2𝑖 ] [ 𝑏1𝑖 ] + 𝜎𝛆2 𝑖 𝑏2𝑖 (6.18) 2 = 𝐛⊤ 𝑖 𝐛𝑖 + 𝜎𝛆𝑖 となります。ここから,多因子の場合でも共通性は因子負荷の二乗に合致するということが わかります。 あとはこれを複数の項目・因子に拡張するだけです。各項目の因子負荷ベクトルを横に並べ 9
6.3 パラメータの推定方法 ていけば,全項目での全因子の負荷行列は 𝐛2 𝐁 = [𝐛1 𝑏11 ⎡ 𝑏 ⎢ 21 =⎢ ⋮ ⎢𝑏𝑇 −1,1 ⎣ 𝑏𝑇 ,1 ⋯ 𝐛𝐼 ] 𝑏12 𝑏22 ⋮ ⋯ 𝑏1𝐼 ⋯ 𝑏2𝐼 ⎤ ⎥ ⋱ ⋮ ⎥ ⋯ 𝑏𝑇 −1,𝐼 ⎥ ⋯ 𝑏𝑇 ,𝐼 ⎦ 𝑏𝑇 −1,2 𝑏𝑇 ,2 (6.19) となります。すると観測変数の相関行列 𝚺 は結局 𝚺 = 𝐁⊤ 𝐁 + 𝚿 𝐛⊤ 1 ⎡𝐛⊤ ⎤ = ⎢ 2 ⎥ [𝐛1 𝐛2 ⋯ ⎢ ⋮ ⎥ ⎣𝐛⊤ 𝐼⎦ 𝐛⊤ 𝐛⊤ ⋯ 1 𝐛1 1 𝐛2 ⎡𝐛⊤ 𝐛 𝐛⊤ 𝐛 ⋯ 1 2 2 2 =⎢ ⎢ ⋮ ⋮ ⋱ 𝐛⊤ ⋯ ⎣𝐛⊤ 𝐼 𝐛1 𝐼 𝐛2 𝐛𝐼 ] + 𝚿 (6.20) 𝐛⊤ 1 𝐛𝐼 ⎤ 𝐛⊤ 2 𝐛𝐼 ⎥ ⋮ ⎥ 𝐛⊤ 𝐼 𝐛𝐼 ⎦ +𝚿 と分解でき,1 因子のときと同じようにして推定が可能となります。 6.3.3 因子の数 因子分析も回帰分析モデルと同じように,共通因子の成分 𝐅𝐛𝑖 によって回答 𝐲 を予測(ある いは近似)することが目的です。回帰分析では決定係数が高いほど予測の精度が良かったよ うに,因子分析でも共通性の割合が高いほうが基本的には嬉しいものです。当然ながら,共 通性の割合は因子の数が増えるほど高くなるわけですが,因子分析では因子の数は自由に設 定できます。ここで,因子の数に関する数理的な性質に少し触れておきます。実際に因子数 を決めるための方法に関しては6.9節にて紹介します。 相関行列をベクトルに分解するためには,固有値分解を利用します。細かいことは省略しま すが,固有値分解は対称行列(ここでは 𝐼 × 𝐼 相関行列 𝚺)を 𝚺 = 𝐗𝚲𝐗⊤ (6.21) という要素に分解する計算です。ここで,𝐗 は 𝐼 本の線形独立な長さ 𝐼 の縦ベクトル(固有 ベクトル)が横に並んだ形 𝐗 = [𝐱1 𝐱2 𝑥11 ⎡𝑥 ⋯ 𝐱𝐼 ] = ⎢ 21 ⎢ ⋮ ⎣𝑥𝐼1 をしています。また,𝚲 は固有値(𝜆1 , 𝜆2 , ⋯ , 𝜆1 ⎡0 𝚲=⎢ ⎢ ⋮ ⎣0 0 𝜆2 ⋮ 0 𝑥12 𝑥22 ⋮ 𝑥𝐼2 ⋯ 𝑥1𝐼 ⋯ 𝑥2𝐼 ⎤ ⎥ ⋱ ⋮ ⎥ ⋯ 𝑥𝐼𝐼 ⎦ (6.22) 𝜆𝐼 ) を対角成分にもつ対角行列 ⋯ 0 ⋯ 0⎤ ⎥ ⋱ ⋮ ⎥ ⋯ 𝜆𝐼 ⎦ (6.23) になっています。R では eigen() という関数で固有値分解を行うことができるので一応見 てみましょう。 10
6.3 パラメータの推定方法 相関行列を固有値分解 1 eig <- eigen(cor(dat[, 2:6])) 2 eig$values # 固有値 1 [1] 2.4143113 0.8819440 0.6987614 0.5517027 0.4532806 1 eig$vectors # 固有ベクトル 1 [,1] 2 [1,] -0.3333328 0.8503498 -0.132707026 -0.3843622 -0.02118902 3 [2,] -0.4936057 0.1622930 0.002823959 0.7651945 4 [3,] -0.5164665 -0.1419179 0.239729691 0.0954805 -0.80407389 5 [4,] -0.4058417 -0.3777438 -0.802813913 -0.2100389 0.06305300 6 [5,] -0.4623825 -0.2962012 0.45227951 [,2] [,3] [,4] [,5] 0.38011021 0.529528437 -0.4620716 ここで,固有ベクトルに関してはそれぞれ要素の二乗和が 1,つまり 𝐼 ∑ 𝑥2𝑖𝑗 = 𝐱𝑗⊤ 𝐱𝑗 = 1 ∀𝑗 (6.24) 𝑖=1 となるように標準化された値が算出されています。 固有値分解から相関行列を復元 1 # 復元した値 2 eig$vectors %*% diag(eig$values) %*% t(eig$vectors) 1 [,2] [,3] [,4] [,5] [1,] 1.0000000 0.3527771 0.2744489 0.1616950 0.1945084 3 [2,] 0.3527771 1.0000000 0.4974111 0.3501907 0.3925362 4 [3,] 0.2744489 0.4974111 1.0000000 0.3848005 0.5131434 5 [4,] 0.1616950 0.3501907 0.3848005 1.0000000 0.3211529 6 [5,] 0.1945084 0.3925362 0.5131434 0.3211529 1.0000000 1 # もとの相関係数 2 cor(dat[, 2:6]) 1 11 [,1] 2 Q1_1 Q1_2 Q1_3 Q1_4 Q1_5 2 Q1_1 1.0000000 0.3527771 0.2744489 0.1616950 0.1945084 3 Q1_2 0.3527771 1.0000000 0.4974111 0.3501907 0.3925362 4 Q1_3 0.2744489 0.4974111 1.0000000 0.3848005 0.5131434 5 Q1_4 0.1616950 0.3501907 0.3848005 1.0000000 0.3211529 6 Q1_5 0.1945084 0.3925362 0.5131434 0.3211529 1.0000000
6.3 パラメータの推定方法 確かに完全に一致しています。ここでは固有値分解の意味までは踏み込みません。とりあえ ず固有値分解を使えば行列をベクトル(とスカラー)に分解できるということを理解してお いてください。 そして固有値分解の式を,𝚲 の対角成分の要素(固有値)ごとに分けて見てみると, 𝚺 = 𝜆1 𝐱1 𝐱1⊤ + 𝜆2 𝐱2 𝐱2⊤ + ⋯ + 𝜆𝐼 𝐱𝐼 𝐱𝐼⊤ (6.25) という形になっています。つまり各固有ベクトルに,対応する固有値をかけたものの和にな っているということです。これの何が嬉しいかというと,固有値の大きいところだけを使っ た近似ができるということです。例えば上の式に関して,固有値の一番大きい 𝜆1 に関する ところだけを使った 𝚺 ≅ 𝜆1 𝐱1 𝐱1⊤ −0.333 ⎡−0.494⎤ ⎢ ⎥ = 2.414 ⎢−0.516⎥ [−0.333 −0.494 −0.516 ⎢−0.406⎥ ⎣−0.462⎦ 0.268 0.397 0.416 0.327 0.372 ⎡0.397 0.588 0.615 0.484 0.551⎤ ⎢ ⎥ = ⎢0.416 0.615 0.644 0.506 0.577⎥ ⎢0.327 0.484 0.506 0.398 0.453⎥ ⎣0.372 0.551 0.577 0.453 0.516⎦ −0.406 −0.462] は,単一のベクトルと定数による近似の中で最も良い近似になっている*3 のです。そしてこ の形は,先程因子分析の解を出した時に出てきた 𝐛⊤ 𝐛 に対応しているので, 𝐛⊤ 𝐛 = 𝜆1 𝐱1 𝐱1⊤ (6.26) という等式を考えれば,因子負荷ベクトルを 𝐛 = √𝜆1 𝐱1⊤ としたときに,観測変数の相関行 列を最も良く近似できるということがわかります。 同様に, 𝚺 ≅ 𝜆1 𝐱1 𝐱1⊤ + 𝜆2 𝐱2 𝐱2⊤ (6.27) は 2 つのベクトルを使った近似の中では最も良い近似になっており,因子 1,2 に関する因子 負荷の横ベクトルをそれぞれ 𝐛1 , 𝐛2 とおけば, 𝐁⊤ 𝐁 = [𝐛⊤ 1 = [𝐱1 𝐛⊤ 2][ 𝐱2 ] [ = [√𝜆1 𝐱1 𝐛1 ] 𝐛2 𝜆1 0 0 𝐱1⊤ ][ ] 𝜆2 𝐱2⊤ √𝜆2 𝐱2 ] [ (6.28) √𝜆1 𝐱1⊤ ] √𝜆2 𝐱2⊤ という形で,因子の数がいくつであっても一因子のときと同様に,相関行列を最も良く近似 できる因子負荷行列 𝐁 を求めることができるわけです。 *3 厳密には,ここで近似しているのは 𝚺 − 𝚿 = 𝐛⊤ 𝐛 です。しかし実際には 𝚿 はわからないので,計算上は 「対角成分を 𝚿 に置き換えた 𝚺」について固有値分解を繰り返すことで,「非対角成分だけ近似された状態」 を求めます。なので上の近似も対角成分(全て 1 になるはず)はあまりうまく近似されていません。 12
6.3 パラメータの推定方法 ここまでの議論からは • 因子の数は最大でも 𝐼 個である(固有値分解するため) • 因子の数が増えるほどデータの相関行列をより良く近似できる • 固有ベクトルはすべて標準化(𝐱𝑖⊤ 𝐱𝑖 = 1)されており固有値は大きい順に並んでいる ため,1 個目の因子が一番強い といったことがわかります。 ではここで,固有ベクトルによる相関行列の近似のイメージを掴むための例を見てみます*4 。 (仮想的な)48 個の因子の相関行列をヒートマップに落とし込んだものです。 以下の図6.7は, 無相関のところを白色として,黄色が濃いほど正の相関が強く,黒が濃いほど負の相関が強 い場所を表すように色をつけたところ,偶然ニッコリ笑顔にみえるような相関行列が得られ たとします。 0 10 cor 20 1.0 y 0.5 0.0 -0.5 30 40 50 0 10 20 x 30 40 50 図6.7 相関行列 では,この相関行列を固有値分解して,最初の固有ベクトルのみを使った近似 𝜆1 𝐱1 𝐱1⊤ を 見てみましょう。うまく近似できたら図6.7のようなニッコリ笑顔になっているはずですが, 図6.8にはその面影は全くありません…。 0 10 cor 20 y 0.50 0.25 0.00 30 40 50 0 10 20 30 40 50 x 図6.8 *4 13 第一固有ベクトルのみでの近似は失敗 ここで紹介しているのは画像処理における次元圧縮のマネです。
6.3 パラメータの推定方法 ということで,もう少し多くの固有ベクトルを使って近似してみましょう。図6.9ではそれぞ れ c(2:5, 10, 48) 個の固有ベクトルを使った近似に基づくヒートマップが描かれていま す。当然 48 個全てを使えば完璧な復元ができる((6.25)式)わけですが,この図を見ると大 体 5 個くらいでもニッコリ笑顔の感じは割と復元されているように見えます。固有値分解を 用いることで,本来の変数の数(この例では 48)よりもかなり少ない数(この例では 5-10 く らい)で,もとの相関行列の大体の感じを説明可能となります。これこそが因子分析が目指 している「なるべく少ない共通因子で説明する」の本質なのです。観測された 𝐼 × 𝐼 相関行列 𝚺 を固有値分解し,そのうち 𝑇 (𝑇 < 𝐼) 個の固有値と固有ベクトルを使って,相関行列を 𝚺 ≅ 𝜆1 𝐱1 𝐱1⊤ + 𝜆2 𝐱2 𝐱2⊤ + ⋯ + 𝜆𝑇 𝐱𝑇 𝐱𝑇⊤ 𝑇 (6.29) = ∑ 𝜆𝑡 𝐱𝑡 𝐱𝑡⊤ 𝑡=1 という形で近似しているのです。 2個で近似 3個で近似 cor 20 10 1.0 y 0.5 30 0 10 20 x 30 40 40 50 50 5個で近似 0 10 20 x 30 40 40 50 50 40 0 10 20 x 30 40 50 x 30 40 50 10 cor 1.5 1.0 0.5 0.0 -0.5 20 30 40 50 20 0 10 20 x 30 40 cor 20 50 0.0 30 -0.5 40 50 1.0 0.5 y y 30 10 48個全部 10 1.5 1.0 0.5 0.0 -0.5 0 0 y cor 20 1.0 0.5 0.0 -0.5 30 0 10 cor 20 10個で近似 0 50 1.0 0.5 0.0 -0.5 30 -0.5 10 cor 20 0.0 40 0 y 10 50 4個で近似 0 y 0 0 10 20 x 30 40 50 図6.9 もっといっぱい使って近似 今回の例では 1,2 個の固有ベクトルによる近似はほぼ失敗していましたが,これはこの48 個の変数に共通する要因を1,2個だけの共通因子で説明するのは無理があることを意味して います。例えば最初の固有ベクトル 1 個だけを使った近似は,因子分析で言えばすべての 項目が概ね単一の共通因子のみの影響を受けている状態(図6.3)を考えているわけですが, 図6.7の相関行列では負の相関を持つ(黒く塗りつぶされた)ペアもそれなりに存在している ため,一因子ではさすがに無理がある,ということを表しています。 ちなみに,すべて項目の因子負荷の二乗和を固有値分解的に表すと, 𝑇 𝐼 tr(𝐁⊤ 𝐁) = ∑ ∑ 𝜆𝑡 𝑥2𝑖𝑡 𝑡=1 𝑖=1 𝑇 𝐼 = ∑ 𝜆𝑡 ∑ 𝑥2𝑖𝑡 𝑡=1 14 𝑖=1 (6.30)
6.4 とりあえずやってみよう という形になります。固有値分解で得られる固有ベクトルは標準化されているため, 𝐼 ∑ 𝑥2𝑖𝑡 = 1 ∀𝑡 (6.31) 𝑖=1 となっています。したがって, 𝑇 tr(𝐁⊤ 𝐁) = ∑ 𝜆𝑡 (6.32) 𝑡=1 となります。これは,固有値こそが因子負荷の合計を決めており,ひいては分散説明率(共 通性)の合計を規定していることを表しています。ここからも,因子の数を決める際に何を 意識したら良いかが見えてきます。 • 因子の数が増えるほど分散説明率が増える=共通性の値が大きくなる • 後ろの固有値になるほど影響が小さくなるので,無くても良くなってくる 実際にどういった指標を用いて因子数を決めたら良いかは6.9節までお待ち下さい。 6.4 とりあえずやってみよう このあたりでひとまず R で因子分析をやってみて,結果を見ながら更に解説していきたいと 思います。R にはデフォルトで factanal() という関数が用意されているのですが,ちょっ とモダンな方法に関して不十分なので,もっと色々できる関数を使用します。psych パッケ ージにある fa() という関数を使用しましょう。 とりあえず因子分析 1 library(psych) 2 # まずは2因子でやってみます 3 res_fa <- fa(dat[, paste0("Q1_", 1:10)], nfactors = 2) 4 res_fa 1 Factor Analysis using method = 2 Call: fa(r = dat[, paste0("Q1_", 1:10)], nfactors = 2) 3 Standardized loadings (pattern matrix) based upon correlation matrix 4 MR1 MR2 h2 u2 com 5 Q1_1 -0.08 0.42 0.16 0.84 1.1 6 Q1_2 0.00 0.68 0.47 0.53 1.0 7 Q1_3 -0.03 0.77 0.58 0.42 1.0 8 Q1_4 0.14 0.46 0.27 0.73 1.2 0.61 0.38 0.62 1.0 9 Q1_5 0.03 10 Q1_6 0.57 -0.05 0.31 0.69 1.0 11 Q1_7 0.64 -0.02 0.40 0.60 1.0 12 Q1_8 0.55 0.02 0.31 0.69 1.0 13 Q1_9 0.68 0.00 0.46 0.54 1.0 14 Q1_10 0.57 0.06 0.35 0.65 1.0 15 16 15 MR1 MR2 minres
6.4 とりあえずやってみよう 17 SS loadings 1.87 1.83 18 Proportion Var 0.19 0.18 19 Cumulative Var 0.19 0.37 20 Proportion Explained 0.51 0.49 21 Cumulative Proportion 0.51 1.00 22 23 With factor correlations of 24 MR1 25 MR1 1.00 0.32 26 MR2 0.32 1.00 MR2 27 28 Mean item complexity = 29 Test of the hypothesis that 2 factors are sufficient. 1 30 31 df null model = 32 df of Square = 45 with the objective function = 2.14 with Chi 5190.53 the model are 26 and the objective function was 0.16 33 34 The root mean square of the residuals (RMSR) is 35 The df corrected root mean square of the residuals is 0.04 0.05 36 37 The harmonic n.obs is 38 The total n.obs was prob < prob < 2432 with the empirical chi square 311.44 with 1.1e-50 2432 with Likelihood Chi Square = 379.19 with 2.2e-64 39 40 Tucker Lewis Index of factoring reliability = 41 RMSEA index = 42 BIC = 43 Fit based upon off diagonal values = 0.98 44 Measures of factor score adequacy 0.075 0.881 and the 90 % confidence intervals are 0.068 0.082 176.48 45 MR1 MR2 46 Correlation of (regression) scores with factors 0.87 0.88 47 Multiple R square of scores with factors 0.76 0.78 48 Minimum correlation of possible factor scores 0.51 0.56 fa() 関数では,因子数を nfactors という引数で与えます。今回は 2 因子にしてみました。 早速結果を見ていきましょう。 6.4.1 因子負荷量 3 Standardized loadings (pattern matrix) based upon correlation matrix 4 16 MR1 MR2 h2 u2 com 5 Q1_1 -0.08 0.42 0.16 0.84 1.1 6 Q1_2 0.00 0.68 0.47 0.53 1.0
6.4 とりあえずやってみよう 7 Q1_3 -0.03 0.77 0.58 0.42 1.0 8 Q1_4 0.14 0.46 0.27 0.73 1.2 9 Q1_5 0.03 0.61 0.38 0.62 1.0 10 Q1_6 0.57 -0.05 0.31 0.69 1.0 11 Q1_7 0.64 -0.02 0.40 0.60 1.0 12 Q1_8 0.55 0.02 0.31 0.69 1.0 13 Q1_9 0.68 0.00 0.46 0.54 1.0 14 Q1_10 0.57 0.06 0.35 0.65 1.0 3-14 行 目 の 部 分 で す。Standardized loadings (pattern matrix) based upon correlation matrix とあるように,データの相関行列から因子負荷行列を出しています。 左の MR1 および MR2 という部分が因子負荷(𝐁⊤ )を表しています*5 。論文では,良く「値 が一定以上のところを太字にする」や「一定以下のところを空白にする」という形で記載さ れています。試しに値が 0.4 以上のところを太字にすると,表6.1のようになります。 表6.1 因子負荷の表 MR1 MR2 Q1_1 -0.081 0.417 Q1_2 0.003 0.684 Q1_3 -0.027 0.770 Q1_4 0.140 0.463 Q1_5 0.028 0.605 Q1_6 0.574 -0.047 Q1_7 0.638 -0.017 Q1_8 0.553 0.021 Q1_9 0.680 0.001 Q1_10 0.575 0.057 また,print() 関数に引数 sort=TRUE を与えると,因子負荷行列がいい感じになるように 項目を並び替えてくれます。 1 print(res_fa, sort = TRUE) 1 Factor Analysis using method = 2 Call: fa(r = dat[, paste0("Q1_", 1:10)], nfactors = 2) 3 Standardized loadings (pattern matrix) based upon correlation matrix 4 5 Q1_9 *5 17 item MR1 9 0.68 MR2 h2 minres u2 com 0.00 0.46 0.54 1.0 MR というのは,因子負荷行列を最小残差法 (Minimum Residual) で計算したことを表しています。詳細 は6.8節にて説明します。
6.4 とりあえずやってみよう 6 Q1_7 7 7 Q1_10 10 8 Q1_6 6 0.57 -0.05 0.31 0.69 1.0 0.64 -0.02 0.40 0.60 1.0 0.57 0.06 0.35 0.65 1.0 9 Q1_8 8 0.55 0.02 0.31 0.69 1.0 10 Q1_3 3 -0.03 0.77 0.58 0.42 1.0 11 Q1_2 2 0.00 0.68 0.47 0.53 1.0 12 Q1_5 5 0.03 0.61 0.38 0.62 1.0 13 Q1_4 4 0.14 0.46 0.27 0.73 1.2 14 Q1_1 1 -0.08 0.42 0.16 0.84 1.1 15 16 (以下略) 出力を見ると,確かに因子負荷行列の行の並び順が変わっています。一番上には「因子1へ の負荷が最も高い項目のうち,最も MR1 が高い項目」がおかれ,後もそれに続いています。 つまり全ての項目を「どれか一つの因子のグループに属するとしたら」という感じで並び替 えてくれるわけです。論文に載せる際なども,これで表示される順にすると見やすいのでは ないかと思います。 表6.2 因子負荷の表(ソート後) MR1 MR2 Q1_9 0.680 0.001 Q1_7 0.638 -0.017 Q1_6 0.574 -0.047 Q1_10 0.575 0.057 Q1_8 0.553 0.021 Q1_3 -0.027 0.770 Q1_2 0.003 0.684 Q1_5 0.028 0.605 Q1_4 0.140 0.463 Q1_1 -0.081 0.417 また,psych パッケージの fa.diagram() 関数を使うと,図6.10のようにグラフィカルモデ ル的なものを描いてくれます。もちろんそのまま論文に載せたりできるレベルではないです が,項目と因子の関係をざっと眺めたいときには役に立ちそうです。 1 fa.diagram(res_fa, cut = 0.3) デフォルトでは各項目には「最も高い因子負荷を持つ因子」のみが矢印を引くようですが, 引数 simple=FALSE とすることで複数因子からの矢印を引くことも可能です。ただ数字が重 なったりして見にくいかもしれません。また,引数 cut は「その値以下の因子負荷や因子間 相関の矢印は引かない」というものです。なので先程のコードでは,因子負荷や因子間相関 18
6.4 とりあえずやってみよう Factor Analysis Q1_9 0.7 Q1_7 0.6 Q1_10 0.6 MR1 0.6 Q1_6 0.6 Q1_8 0.3 Q1_3 0.8 Q1_2 0.7 MR2 0.6 Q1_5 0.5 Q1_4 0.4 Q1_1 図6.10 グラフィカルモデル が 0.3 以下のものは省略されている…と言っても「最も高い因子負荷」に 0.3 以下のものが無 いので違いが分かりませんが…。 因子分析は相関行列を対象にした分析法ですが,その目的の一つは項目をグルーピングする とともに,その特性を表す軸を決定することです。図6.11は X 軸に MR1 を,Y 軸に MR2 の 値をとった各項目の因子負荷量のプロットです。psych パッケージにある fa.plot() 関数 で簡単に描くことが出来ます*6 。 1 fa.plot(res_fa) 0.8 Factor Analysis 0.4 2 5 4 1 0.2 MR2 0.6 3 0.0 10 8 6 0.0 0.2 0.4 7 9 0.6 MR1 図6.11 因子負荷のプロット たしかに Q1_1 から Q1_5 のグループと,Q1_6 から Q1_10 のグループに分かれているように 見えますね。そしてこのプロットにおいて,それぞれの軸は因子を表しています。Q1_1 から Q1_5 のグループでは MR1 の値がゼロに近く,反対に Q1_6 から Q1_10 のグループでは MR2 の値がゼロに近くなっています。この状態ならば,Q1_1 から Q1_5 のグループの共通項(因 子)を解釈する上では MR1 は無視して良く,Q1_1 から Q1_5 のグループの共通項には MR2 *6 19 ちなみに 3 因子以上の場合,fa.plot() は各因子ペアごとの散布図を並べてくれます。
6.4 とりあえずやってみよう は無関係とみなして良さそうなので解釈が楽そうです。もしも「複数の因子に高い負荷量を もつ項目」があると,共通項の解釈はかなり大変そうです。 ということで,因子分析では一つの項目はなるべく一つの因子だけに高い負荷を持つという 状態が望ましいとされています。これを単純構造 (simple structure) と呼びます。因子負荷 量の出力の一番右の列にあった com はこの単純構造の程度を表した指標です。これが 1 に近 いほど,その項目は 1 つの因子にのみ高い負荷を持っているという判断ができます。今回の 結果では,Q1_4 が 1.2 と最も高い値になっています。たしかに,他の項目と比べると MR1 と MR2 の差が小さいようです。が,実際のところ 1.2 は相当マシな値です。どれくらい大きか ったら気をつけたら良いか,という基準は特に無いですが… ちなみに,固有値分解によって因子負荷を推定すると解は一つ(√𝜆𝑡 𝐱𝑡 )に決まっているの で,普通に考えると単純構造になるかどうかは運任せ(たまたま単純構造に近い解が固有値 分解で得られるか)になってしまいそうです。ただ実際にはそんなことはありません。これ には因子分析の解の不定性が大きく関係しています。 6.4.2 因子の回転 まずは数理的に説明しましょう。(6.10) 式では,単一の項目 𝑖 に関する因子分析モデルを,回 帰分析と同様に 𝐲𝑖 = 𝐟1 𝑏1𝑖 + 𝐟2 𝑏2𝑖 + ⋯ + 𝐟𝑇 𝑏𝑇 𝑖 + 𝛆𝑖 = [𝐟1 𝑏1𝑖 ⎡𝑏 ⎤ ⋯ 𝐟𝑇 ] ⎢ 2𝑖 ⎥ + 𝛆𝑖 ⎢ ⋮ ⎥ ⎣𝑏𝑇 𝑖 ⎦ 𝐟2 (6.33) = 𝐅𝐛𝑖 + 𝛆𝑖 と説明しました。これを全項目についてまとめて線形代数で表すと,以下のようになります。 𝐘 = [𝐲1 𝐲2 ⋯ 𝑓11 𝑓12 ⎡𝑓 𝑓22 = ⎢ 21 ⎢ ⋮ ⋮ ⎣𝑓𝑃 1 𝑓𝑃 2 = 𝐅𝐁 + 𝐄 𝐲𝐼 ] ⋯ ⋯ ⋱ ⋯ 𝑓1𝑇 𝑏11 𝑓2𝑇 ⎤ ⎡ 𝑏21 ⎥⎢ ⋮ ⎥⎢ ⋮ 𝑓𝑃 𝑇 ⎦ ⎣𝑏𝑇 1 𝑏12 𝑏22 ⋮ 𝑏𝑇 2 ⋯ 𝑏1𝐼 𝜀11 ⋯ 𝑏2𝐼 ⎤ ⎡ 𝜀21 ⎥+⎢ ⋱ ⋮ ⎥ ⎢ ⋮ ⋯ 𝑏𝑇 𝐼 ⎦ ⎣𝜀𝑃 1 𝜀12 𝜀22 ⋮ 𝜀𝑃 2 ⋯ 𝜀1𝐼 ⋯ 𝜀2𝐼 ⎤ ⎥ ⋱ ⋮ ⎥ ⋯ 𝜀𝑃 𝐼 ⎦ (6.34) 表記がこんがらがってきたのでおさらいしておきますが, • 因子得点 𝑓𝑝𝑡 にある 2 つの添字は回答者(𝑃 )と因子(𝑇 ) • 因子負荷 𝑏𝑡𝑖 にある 2 つの添字は因子(𝑇 )と項目(𝐼 ) • 残差 𝜀𝑝𝑖 にある 2 つの添字は回答者(𝑃 )と項目(𝐼 ) です。ここで,共通因子に関する部分 𝐅𝐁 について考えてみます。 • 因子負荷行列 𝐁 に対して,正則な(=逆行列が存在する)対称行列 𝐀 を左からかけ たものを 𝐁̃ = 𝐀𝐁 とする • 因子得点行列 𝐅 に対して,その逆行列 𝐀−1 を右にかけたものを 𝐅̃ = 𝐅𝐀−1 とする 20
6.4 とりあえずやってみよう 表6.3 因子負荷行列 (𝐁⊤ ) F1 F2 0.702 0.238 0.683 0.500 0.787 0.532 0.811 -0.492 0.694 -0.499 0.665 -0.273 すると,𝐅̃ 𝐁̃ = 𝐅𝐀−1 𝐀𝐁 = 𝐅𝐁 となります。そして正則な対称行列 𝐀 というのは無限に 存在しています。したがって,𝐅𝐁 は割と好きなように変換することができるわけです。 𝐅𝐁 を自由に変換できると何が嬉しいのでしょうか。もう少しイメージを交えて説明しまし ょう。ある 6 項目の相関行列に対して固有値分解を行い,得られた 2 つの因子に対する因子 負荷量(表6.3)を二次元にプロットしたものが図6.12です。どうやら 3 項目ずつのグループ に分かれそうですね。ですがこのままでは全ての項目が因子 1 に対して高い負荷を持ってい るため,因子 1 は「総合的な何か」としか解釈できません。 1 0.75 0.5 0.25 -1 -0.75 -0.5 -0.25 0.25 0.5 0.75 1 -0.25 -0.5 -0.75 -1 図6.12 仮想的な 2 因子 図6.12のようなプロットに対して,正則な対称行列 𝐀 による変換は,プロットの軸を回転さ せることに相当します。例えば 2 次元の回転ならば, 𝐀=[ cos 𝜃 sin 𝜃 − sin 𝜃 ] cos 𝜃 (6.35) という行列を左からかけることは,原点を中心に時計回りに 𝜃 度回転することに対応してい ます。また,この行列には逆行列 𝐀−1 = [ 21 cos 𝜃 − sin 𝜃 sin 𝜃 ] cos 𝜃 (6.36)
6.4 とりあえずやってみよう が存在しているため,𝐅̃ = 𝐅𝐀−1 の計算も可能です。試しに図6.13のように軸を時計回りに 45 度回転させると,図6.12よりは単純構造に近づけることができそうです。 25 1. 25 . -1 1 -1 75 0. 5 .7 -0 0. 5 .5 -0 25 0. 25 . -0 25 . -0 25 0. .5 -0 5 0. 5 .7 -0 75 0. -1 1 5 .2 -1 25 1. 図6.13 45 度回転させると 実際に因子負荷行列を変換させてみると, 𝐀𝐁 = [ cos 45∘ sin 45∘ − sin 45∘ 0.702 ][ 0.238 cos 45∘ −0.707 0.702 ][ 0.707 0.238 0.683 0.787 0.500 0.532 =[ 0.707 0.707 =[ 0.328 0.129 0.180 0.921 0.665 0.836 0.933 0.226 0.683 0.787 0.500 0.532 0.811 −0.492 0.811 −0.492 0.694 −0.499 0.694 −0.499 0.665 ] −0.273 0.665 ] −0.273 0.844 0.663 ] 0.138 0.277 となり,表6.3よりはいい感じに 3 項目ずつに分かれました。 ここまでは,2 つの軸が直交のままに回転をさせてきました。このような回転のことを直交 回転 (orthogonal rotation) と呼びます。固有値分解で得られる固有ベクトルは,互いに直 交するようになっています。直交するということは,言い換えると相関がないということで す。つまり,直交回転は因子間相関がゼロという制約のもとで単純構造を目指す回転なので す。ですが,ほとんどの場合ではわざわざ因子間相関をゼロに制約する必要はないでしょう。 むしろ,因子の間には少なからず相関関係が見られるはずです。 2 つの因子の相関係数は,2 つの軸が作る角度 𝜃 に対して cos 𝜃 になります。直交回転では必 ず相関が cos 90∘ = 0 になっているわけですが,もし 2 つの因子の軸が作る角度が 90∘ 以外 になっても良いならば,図6.13よりも更に良い単純構造が目指せそうです。実際に,例えば 回転行列 𝐀 を 𝐀=[ 0.613 −0.902 ] 0.582 0.922 と設定すると,変換後の因子負荷行列 𝐁̃ は 𝐀𝐁 = [ =[ 22 0.613 −0.902 0.702 ][ 0.582 0.922 0.238 0.216 −0.032 0.628 0.859 0.683 0.787 0.500 0.532 0.811 −0.492 0.002 0.941 0.876 0.949 0.018 −0.057 0.694 −0.499 0.654 ] 0.135 0.665 ] −0.273
6.4 とりあえずやってみよう となり,直交回転のときよりもさらに単純構造になりました。 図6.14 斜交回転 図6.14のように,因子間の相関を認めながら回転を行うことを斜交回転 (oblique rotation) と呼びます。なお,回転行列 𝐀 によって初期解を回転させた時,因子得点行列 𝐅 は 𝐅̃ = 𝐅𝐀−1 と変換されているため,因子間相関行列 𝚽 は (𝐀𝐀⊤ )−1 によって求めることができます。例 えば先程の数値例の場合は, 𝚽 = (𝐀𝐀⊤ )−1 = ([ ≃[ 0.613 −0.902 0.613 0.582 ][ ]) 0.582 0.922 −0.902 0.922 −1 1 0.400 ] 0.400 1 となり,回転後の因子間相関は 0.4 と計算できます。cos 𝜃 = 0.4 ということは,図6.14の 2 軸が作る角度がおよそ 𝜃 = 66.4 度になっている,ということです。 回転後の因子間相関 あるデータ行列 𝐗 は,各行が個人を,各列が変数を表している 𝑃 × 𝐼 の行列としま す(dat などと同じ形式)。これに対する分散共分散行列 𝚺 の (𝑠, 𝑡) 成分は変数 𝑠 と 変数 𝑡 の共分散なので, 𝜎𝑠𝑡 = 1 𝑃 ∑(𝑥 − 𝑥𝑠̄ )(𝑥𝑝𝑡 − 𝑥𝑡̄ ) 𝑃 𝑝=1 𝑝𝑠 = 𝐸 [(𝑥𝑝𝑠 − 𝑥𝑠̄ )(𝑥𝑝𝑡 − 𝑥𝑡̄ )] と単純に因子得点の積の期待値になります。ただし 𝑥𝑠̄ は変数 𝑠 の平均値です。 これを因子得点行列 𝐅 に置き換えて考えてみると,因子分析では基本的な設定の一 つとして,各行(各因子の因子得点)の平均値(𝑥𝑠̄ や 𝑥𝑡̄ に相当するもの)はそれぞ 23
6.4 とりあえずやってみよう れ 0 という制約を置いているため,因子 𝑠 と因子 𝑡 の共分散は 𝜎𝑠𝑡 = 1 𝑃 ∑𝑓 𝑓 𝑃 𝑝=1 𝑝𝑠 𝑝𝑡 = 𝐸 [𝑓𝑝𝑠 𝑓𝑝𝑡 ] となります。同様に,各列(因子得点)の標準偏差は 1 に制約されているため,分散 共分散行列 𝚺 はそのまま相関行列として見ることができます。 ということで因子得点の相関行列は,上の式を行列全体に拡張することで 𝐸[𝐅⊤ 𝐅] と表すことができます。回転前の初期解(固有ベクトル同士)は必ず無相関なので, ⊤ 𝐸[𝐅⊤ 𝐅] = 𝐈 ということです。𝐅̃ ⊤ = (𝐅𝐀−1 )⊤ = 𝐀−1 𝐅⊤ であることを利用すると, 回転後の因子間相関は ⊤ 𝐸[𝐅̃ ⊤ 𝐅]̃ = 𝐸[𝐀−1 𝐅⊤ 𝐅𝐀−1 ] ⊤ = 𝐀−1 𝐸[𝐅⊤ 𝐅]𝐀−1 ⊤ = 𝐀−1 𝐈𝐀−1 ⊤ = 𝐀−1 𝐀−1 = (𝐀𝐀⊤ )−1 となるわけです。 直交回転にせよ斜交回転にせよ,図6.12-6.14における各項目のプロットの位置関係は変わり ません。また,共通性と独自性の値も変わりません。回転の不定性を利用することによって, 「相関行列から項目のグルーピングを行う」という基本姿勢は変わること無く,より解釈しや すいように因子負荷量をいい感じにしてあげるわけです。 3,4 列目の h2 と u2 はそれぞれ共通性と独自性を表しています。Q1_1 の場合,共通性が 0.16, 独自性が 0.84 ということで,観測変数の分散(標準化されているので 1)のうち 0.16 だけ がこの 2 つの因子によって説明されていることになります。なお h2 の値は単純な因子負荷 ̃ の値になっています*7 。 の二乗和ではなく,回転後の因子間相関を考慮した diag(𝐁̃ ⊤ 𝚽𝐁) 6.4.3 因子の解釈と項目の選定 因子分析では,得られた因子負荷行列をもとに,因子の解釈を行います。表6.1や(2 因子な ら)図6.11をもとに,それぞれの因子に負荷が高い項目のグループを見て,共通の要因を見 つけてあげます。ここは領域の知識などが要求されるところです。なので同じデータであっ ても人によっては異なる名前をつける可能性があります。そして解釈した因子の名前が真っ 当かどうかは,査読者ないし読者が判断することになります。 また,当然ながら因子分析の結果は手元のデータから得られるものです。そのため,先行研 究で作成された尺度をそのまま用いた場合でも,因子構造が同じものが得られるとは限りま *7 24 回転前後で共通性の値は変わらないため,h2 の値は回転前の初期解の因子負荷の二乗和には一致します。回 転して因子間相関を持ったことで,各因子の影響のうち僅かな部分が別の因子を経由したものに置き換わっ たというイメージです。
6.4 とりあえずやってみよう せん。特に先行研究と属性の異なる集団でデータを集めた場合や,古い尺度を使用する場合 には改めて因子分析を行い,同じ因子構造が得られることを確認してから分析を進めるのが 良いでしょう。 今回のデータのように各項目がほぼ一つの因子のみに高い負荷を持っている場合には解釈は 容易なのですが,多くの実データではきっと複数の因子に高い負荷を持っている項目や,ど の因子にも高い負荷を持っていない項目が出現するでしょう。どの因子にも高い負荷を持っ ていない項目は,言い換えると独自性が高い項目であり,因子得点を算出する上ではあって もなくてもあまり変わらないと予想されます。したがって,項目数を減らしたい場合などに は優先的に除外対象としたら良いでしょう。 複数の因子に高い負荷を持っている項目の扱いは結構悩ましいところです。尺度作成論文で は, 「複数の因子に 0.3 以上の負荷を持つ項目は削除した」的な手続きが書かれていることが よくあります。この手続きが正しいかどうかは,時と場合によると思われます。 一般的に心理尺度を作成する場合には,なるべく「みんなが使いやすいモノサシ」を開発し ようとしています。そして,共通のモノサシを用いることで異なる研究間での結果を比較可 能にしようと考えたりしています。そのため,尺度得点としては因子得点ではなく和得点を 用いることが多いです。因子得点を用いた場合,データ毎に因子構造が変わってしまうので 異なるデータ間での因子得点の比較などは難しくなってしまいます。さて,下位概念ごとの 和得点を計算する場合,直感的には「その因子に高い負荷を持っている項目の和」を使えば 良さそうです。この時に,複数の因子に高い負荷を持つ項目があると,その項目の回答は複 数の下位概念でカウントされてしまうため,実質的に他の項目の倍の価値を持ってしまうこ とになります。こうした問題点を避けるという目的であれば,複数の因子に高い負荷を持つ 項目は使わずに,もっと使いやすい(単純構造な)項目を使って尺度を構成するのは悪くな い方法だと思います。 一方で,尺度開発が目的ではない(例えば人に使ってもらうというよりは,単純に項目をグ ルーピングして共通項を取り出したい)場合には,項目数を減らす必要も無ければ和得点を 使う必要もありません。因子得点を計算する上では,複数の因子に高い負荷を持っている項 目も重要(むしろきっと共通性が高いのでかなり役に立つ)です。このような場合に, 「他の 研究でもやっているから」という理由で安直に「複数の因子に高い負荷を持っている項目は 削除」としてしまうのは悪手と言えるでしょう。 また,構造的な問題として「1 項目だけが高い負荷を持っている因子」は通常使用しません。 因子分析の目的は「複数の項目に共通して影響する潜在変数を抽出する」ことなので,他の どの項目とも独立している一匹狼に対して「因子」という見方はあまりしないのです。統計 的にも,1 項目だけの因子の場合,その項目への回答と因子得点は完全に一致するため,標準 化してしまえば回答値をそのまま使えば良いということがいえます。加えて,実は SEM 的 には 1 因子には 3 項目以上あることが望ましかったりもします。 以上の話は,全て統計的な視点からのことです。実際に項目を削除するかは,構成概念と項 目の内容との関係を改めて精査したり,因子の解釈においてノイズになっていないかを検討 したりといった,妥当性の内容的な側面からのチェックも忘れずに行ってください。 25
6.4 とりあえずやってみよう 6.4.4 因子寄与率・分散説明率 16 MR1 MR2 17 SS loadings 1.87 1.83 18 Proportion Var 0.19 0.18 19 Cumulative Var 0.19 0.37 20 Proportion Explained 0.51 0.49 21 Cumulative Proportion 0.51 1.00 一番上の SS loadings は因子負荷の二乗和 (sum of squared loadings) を表していると思 われます。実際には単純な二乗和 diag(𝐁⊤ 𝐁) ではなく,これに回転後の因子間相関を考慮 ̃ が出ています*8 。例えば因子 1 を例にとって見ると,因子 1 の得点 𝐟1 した値 diag(𝐁̃ ⊤ 𝚽𝐁) は,それぞれの項目に対する因子負荷量 𝑏1𝑖 を通じて各項目の分散をある程度 (𝑏1𝑖 ) 説明しま す。ただ,項目 1 の影響はこれにとどまりません。因子 2 と相関があるということは,𝐟1 の 値とある程度連動して 𝐟2 の値も変化します。仮にこの相関関係を回帰分析的に言えば,𝐟1 は 𝐟2 に対して僅かな分散説明率を持っているわけです。ということで,𝐟1 が 𝐟2 を介して各項 目に対してもつ僅かな分散説明率 𝑏2𝑖 も合算されて SS loadings が計算されています。 2 つ目の Proportion Var は,観測変数全体に対する各因子の分散説明率を表しています。 今回は 10 項目のデータで因子分析を行っていますが,各観測変数は標準化されているため, 分散の合計は 10 です*9 。この 10 に対して,因子 1 の寄与の合計 (SS loadings) は 1.87 な ので,Proportion Var は 0.187 となります。つづく Cumulative Var は,その分散説明率 の累積(累積寄与率)です。今回の場合,因子 1 も因子 2 もそれぞれ 0.187 および 0.183 とい う分散説明率(因子寄与率)なので,第二因子の Cumulative Var には,その合計で 0.370 という値が入るわけです。言い換えると,今回の分析では 2 因子でこの 10 項目の分散のう ちおよそ 37% が説明される,ということになります。この Cumulative Var の値があまり にも低い場合は,まだ各項目の独自性が高すぎるということです。せっかく因子得点を算出 しても,結局各項目への回答はその因子以外の要因で決まっている割合が大きいということ になり,因子得点の意味が問われてしまいます。もうお分かりかもしれませんが,因子分析 では因子の数を増やすほど累積寄与率は上がるので,累積寄与率が一定の値になるまで因子 数を増やすというアプローチも無いことは無いでしょう*10 。 一番下の 2 つ Proportion Explained と Cumulative Proportion は,共通因子で説明可 能な部分のうち,各因子が占める割合(および累積割合)を表しています。ここまで解釈す ることはあまり無い気がしますが,例えば Proportion Explained が低すぎる因子がある 場合には,相対的に「なくても良い」という判断に使えるかもしれませんね。 ̃ 𝐁̃ ⊤ ) でした。それに対して 項目の共通性は,因子負荷行列 𝐁 を列(項目)ごとに二乗和したもの diag(𝐁𝚽 ̃ という順番になっている 因子寄与率は因子負荷行列を行(因子)ごとに二乗和しているため,diag(𝐁̃ ⊤ 𝚽𝐁) のです。 *9 「合計点の分散」ではなく「分散の合計」です。なので観測変数間の相関などは気にする必要はありません。 *10 一説には,50% くらいの累積寄与率はほしいという声もあったりするそうです。 *8 26
6.4 とりあえずやってみよう 6.4.5 因子間相関 23 With factor correlations of 24 MR1 25 MR1 1.00 0.32 26 MR2 0.32 1.00 MR2 見ての通り,回転後の因子間の相関行列 𝚽 = (𝐀𝐀⊤ )−1 です。相関があまりに高い場合には, もしかしたら因子を統合してしまっても良いのかもしれません。もちろん内容次第ですが。 6.4.6 因子得点 因子得点 𝐅 は,fa の出力の中に $scores という名前で入っています。 1 head(res_fa$scores) 1 MR1 2 1 -1.38258689 -0.8922837 3 2 -0.34481092 -0.2707195 4 3 -0.11683020 -0.4538879 5 4 -1.08523045 6 5 -0.06184927 -0.8948835 7 6 1.35556644 MR2 0.1082199 0.4285821 観測された因子得点の相関を見てみると,モデル上の回転後の因子間相関 𝚽 とは少し異なり ます。報告する際には気をつけてください。 因子得点の相関と因子間相関は別もの 1 cor(res_fa$scores) 1 MR1 2 MR1 1.000000 0.392582 3 MR2 0.392582 1.000000 MR2 図6.15は因子得点と合計点の散布図です。確かに合計点が高い人ほど因子得点も高くなって いるようです。 1 # MR2はQ1_1からQ1_5に対応 == Q1_A 2 plot(dat$Q1_A, res_fa$scores[, 2]) 尺度得点として,合計点の代わりに因子得点を用いることのメリットは大きく分けて 2 つあ ると考えます。1 つ目は,得点のばらつきです。図6.15を見るとわかるように,合計点が同じ 人であっても因子得点はかなりばらつきがあります。この理由は,項目ごとに因子負荷が異 なるためです。因子負荷が高いほど,因子得点が高い人と低い人を明確に識別できる項目で 27
0 -1 -2 -3 res_fa$scores[, 2] 1 6.5 カテゴリカル因子分析 5 10 15 20 25 30 dat$Q1_A 図6.15 因子得点と合計点の散布図 す。したがって,因子負荷が高い項目に「とても当てはまる」と回答した人はかなりの確度 で因子得点が高いといえます。一方で,因子負荷がほぼゼロの項目では,因子得点が高いか ら「当てはまる」と回答しているとはいえません。したがって,この項目に「とても当ては まる」と回答した人の因子得点が高いという保証はありません。尺度得点として合計点を使 うということは,このように各項目が因子得点の推定に与える情報の強さを無視して,一律 同じパワーを持っていると仮定していることと同じです。そのため,特に因子負荷がばらつ くような状況下では,単純な合計点がその人の潜在的な特性の強さを正確に表しているとは あまり言えないかもしれないのです*11 。 因子得点を用いるもう一つのメリットは,モデル上は誤差が分離されるという点です。図6.4に あるように,因子分析モデルでは各項目への回答は「共通因子」と「独自因子」の影響を受 けていると考えます。このとき,誤差(特に偶然誤差)は,他の項目とは無関係に生じるも のなので独自因子の中に吸収されているという見方ができます。したがって,共通因子の部 分には誤差の影響は含まれていないと考え,因子得点はある意味では純粋な値として扱うこ とができる,とされているのです。 (ここからは因子分析を行うにあたってのプラクティカルな話題をいくつか紹介します。) 6.5 カテゴリカル因子分析 因子分析では,相関行列をもとに因子負荷行列を求めています。その相関行列に関して, Chapter 4 で紹介した I-T 相関では,カテゴリカル変数の場合の特別な相関(ポリシリアル 相関)の話をしました。因子分析の時にも同じことが言えます。因子分析での項目間の相関 関係は,本当は背後にある潜在的な値に対して求められるべきです。そこで,因子分析に用 いる相関をピアソンの積率相関係数ではなくポリコリック相関に変更することを考えます。 このように,各観測変数を連続変数(間隔尺度)ではなく順序尺度とみなして因子分析を行 うことをカテゴリカル因子分析と呼んだりします。 *11 28 尺度としての使いやすさや計算の簡便性を考えると合計点を使うのはありなのですが,個人的には例えば尺 度得点として重み付け和を使う,みたいな方法はもっと広まっても良いのでは,と思っています。
6.5 カテゴリカル因子分析 psych::fa() では,使用する相関係数を引数 cor によって指定することが出来ます。 カテゴリカル因子分析 1 # 毎回列名を指定するのが面倒なので先に作っておく 2 cols <- paste0("Q1_", 1:10) 3 res_fa_cat <- fa(dat[, cols], nfactors = 2, cor = "poly") 4 res_fa_cat 1 Factor Analysis using method = 2 Call: fa(r = dat[, cols], nfactors = 2, cor = "poly") 3 Standardized loadings (pattern matrix) based upon correlation matrix 4 MR2 5 Q1_1 6 7 MR1 h2 minres u2 com -0.09 0.48 0.21 0.79 1.1 Q1_2 0.00 0.74 0.55 0.45 1.0 Q1_3 -0.03 0.81 0.65 0.35 1.0 8 Q1_4 0.16 0.49 0.32 0.68 1.2 0.64 0.43 0.57 1.0 9 Q1_5 0.04 10 Q1_6 0.62 -0.04 0.36 0.64 1.0 11 Q1_7 0.68 -0.03 0.45 0.55 1.0 12 Q1_8 0.58 0.02 0.35 0.65 1.0 13 Q1_9 0.72 0.01 0.53 0.47 1.0 14 Q1_10 0.60 0.06 0.39 0.61 1.0 15 16 MR2 MR1 17 SS loadings 2.11 2.11 18 Proportion Var 0.21 0.21 19 Cumulative Var 0.21 0.42 20 Proportion Explained 0.50 0.50 21 Cumulative Proportion 0.50 1.00 22 23 With factor correlations of 24 MR2 25 MR2 1.00 0.34 26 MR1 0.34 1.00 MR1 27 28 Mean item complexity = 29 Test of the hypothesis that 2 factors are sufficient. 1 30 31 df null model = Square = 32 df of 45 with the objective function = 2.78 with Chi 6749.62 the model are 26 and the objective function was 0.22 33 34 The root mean square of the residuals (RMSR) is 35 The df corrected root mean square of the residuals is 0.04 0.05 36 37 29 The harmonic n.obs is 2432 with the empirical chi square 367.22 with
6.6 因子分析をする前に prob < 38 5.9e-62 The total n.obs was prob < 2432 with Likelihood Chi Square = 540.84 with 1.2e-97 39 40 Tucker Lewis Index of factoring reliability = 41 RMSEA index = 42 BIC = 43 Fit based upon off diagonal values = 0.98 44 Measures of factor score adequacy 0.09 0.867 and the 90 % confidence intervals are 0.084 0.097 338.13 45 MR2 MR1 46 Correlation of (regression) scores with factors 47 Multiple R square of scores with factors 0.79 0.82 48 Minimum correlation of possible factor scores 0.58 0.64 0.89 0.91 第 9 回では,ピアソンの積率相関係数を使って因子分析を行いましたが,それと比べると因 子負荷の値が基本的には大きくなっています。もともと相関係数が過小推定されていたため に,因子負荷も過小推定されていたと考えるのが妥当でしょう。今回カテゴリカル因子分析 を行ったことによって,より正しい因子負荷量を求めることができたといえます。 ただ,萩生田・繁桝 (1996) では,シミュレーション研究をもとに,多くの場合はピアソン の積率相関係数で十分(ポリコリック相関と大差ない)という結論を出していたりもします。 数理的に厳密なのはポリコリック相関ですが,実際には不適解を出してしまったり,回答者 のいないカテゴリがあると面倒だったりと,実用上は問題点もあったりします。また,ポリ コリック相関係数を安定的に推定するためにはデータが数百はあることが望ましいです。そ のため,まずはカテゴリカル因子分析をやってみて,解が得られない or 不安定な場合には 普通にピアソンの積率相関を用いた方法にしても良いのではないかと思います。 6.6 因子分析をする前に 因子分析は,複数の観測変数の背後に共通の要因が存在し,また観測変数間の相関関係がそ の共有の要因のみによって説明される,というモデルを仮定した分析です。したがって,因 子分析を行う前には • そもそも観測変数の間には相関関係が存在しており, • その相関関係が, (特定の2変数間の関係ではなく)概ね全ての観測変数群の関係性に よって説明可能である ことを確認すると良いかもしれません。2点目がちょっと分かりづらいかもしれませんが (言葉でどう説明したら良いのか迷いました) ,これは図6.6に示されるように,全ての観測変 数に共通する相関の要因(共通因子)とは異なる「共通因子では説明できない要因」によっ て2変数の相関が規定されてほしくない,ということを意味しています。そこで,この節で は上記の2つを検証するための方法*12 を簡単に紹介します。 *12 どうやら SPSS ではこの後紹介する分析の結果がデフォルトで出てくるようです。一方でここで紹介してい る分析を(SPSS 以外で分析をする人が)論文にも書かなかったとしても咎められる可能性は高くないので 「出来そうならやってみてもいいんじゃないかなぁ」くらいの位置づけとして紹介します。 30
6.6 因子分析をする前に 6.6.1 Bartlett の球面性検定 Bartlett の球面性検定 (Bartlett, 1950, 1951) では,共通因子の影響を除いたあとの「独自 因子得点の相関行列」𝚿 が単位行列である,という帰無仮説を立てて検定を行っています。 具体的には,以下の検定統計量 (Browne, 1968) 𝜒2𝑚 = − (𝑃 − が自由度 1 2 2𝐼 + 11 2𝑇 − ) ln|𝚿| 6 3 [(𝐼 − 𝑇 )2 − 𝐼 − 𝑇 ] の 𝜒2 分布に漸近的に従うことを利用して統計的仮説検定を 行います。念のため再度示しておきますが, 𝑃 : 回答者 (Person) の数 𝐼 : 項目 (Item) の数 𝑇 : 共通因子 (Trait) の数 を表しています。また,|𝚿| は 𝚿 の行列式を表しています。細かい式の意味はともかく,重 要なのは 𝜒2𝑚 の中の ln|𝚿| の部分です。𝚿 に限らず,相関行列(全ての対角成分が 1)の行 列式の値は必ず 0 から 1 の間の値になること,そして行列式の値が 1 になるのは単位行列 (全ての変数間が無相関)の場合であることが知られています*13 。したがって,その対数を取 った ln|𝚿| は,変数間の相関が高いほど大きなマイナスの値を取るようになっています。と いうことで,これの符号をひっくり返すことで「全体としての相関の強さ」を表しているわ けです。 そんな Bartlett 検定ですが,検定の目的からすると,帰無仮説「𝚿 が単位行列である」が棄 却されなければ嬉しいような気がします。しかし実際には,独自因子間の相関を検証する目 的,つまり因子分析の結果得られた𝚿に対してBartlettの球面性検定を使用することは現状 ほぼ無いようです。 実際の Bartlett の球面性検定(SPSS で出力される結果)では,共通因子の数が0(𝑇 = 0) であると設定して検定を行います。このとき「独自因子得点の相関行列」𝚿 イコール「観測 変数の相関行列」𝚺 となるため,球面性検定の帰無仮説は「観測変数の相関行列が単位行列 である」 ,すなわち全ての観測変数は無相関であると見ることができます。この場合,帰無仮 説が棄却されてくれれば「観測変数間に相関関係がある」ことが言えるため,共通因子を置 いた因子分析を行う意味があると言えるわけです。 R では,psych パッケージに cortest.bartlett() という関数が用意されています。 Bartlett の球面性検定 31 1 # psych::polychoric()関数も使える 2 cor_poly <- polychoric(dat[, cols]) 3 # polychoric()は相関行列以外にもlistで色々返してくる *13 これは,アダマールの不等式というものによって証明される性質です。
6.6 因子分析をする前に 4 # listの中のrhoがポリコリック相関行列の本体 5 cortest.bartlett(cor_poly$rho, n = nrow(dat)) 1 $chisq 2 [1] 6749.621 3 4 $p.value 5 [1] 0 6 7 $df 8 [1] 45 出力された p.value が 0.05 より小さければ帰無仮説が棄却されるのですが,上の 𝜒2𝑚 の式 からもわかるように,この方法も一般的な統計的仮説検定の例に漏れず,サンプルサイズが ある程度大きくなったらほぼ確実に有意になってしまいます。その意味では,回帰分析にお ける F 検定のように,Bartlett の球面性検定はしなかったところで大した問題ではないと言 えるでしょう*14 。 6.6.2 KMO (Kaiser-Meyer-Olkin) 続いては,独自因子に起因する相関があるかを確認する方法を紹介します。ここでは KMO (Kaiser-Meyer-Olkin) または MSA (measure of sampling adequacy) と呼ばれる指 標 (Kaiser, 1970) を紹介します。 KMO はイメージ分析 (Guttman, 1953) と呼ばれる方法に基づく考え方です。これは,複数 の観測変数間の相関関係に関する多変量解析法として,主成分分析・因子分析に続く第3の 方法として提案されたものらしく,因子分析のような不定性が無いことや, (多変量正規)分 布の仮定が必要ないという点で因子分析よりも general な方法とされています(が,現状あ まり使われていないと思います…)。 イメージ分析では,各観測変数 𝑥𝑖 について「それ以外の全ての観測変数」を用いた重回帰分 析を行います。 𝑥𝑖 = 𝑥𝑖̂ + 𝑒𝑖 𝑥𝑖̂ = 𝑤1 𝑥1 + 𝑤2 𝑥2 + ⋯ + 𝑤𝑖−1 𝑥𝑖−1 + 𝑤𝑖+1 𝑥𝑖+1 + ⋯ + 𝑤𝐼 𝑥𝐼 (6.37) すると,観測変数 𝑥𝑖 は「いずれかの観測変数と関係がある部分」𝑥𝑖̂ と「いずれの観測変数 とも関係がない部分」𝑒𝑖 に分けることができます。そしてイメージ分析では,𝑥𝑖̂ のことを イメージスコア,𝑒𝑖 のことを反イメージスコアと呼びます。式からもなんとなく分かるよう に,イメージスコア 𝑥𝑖̂ は,他の観測変数と(相関)関係がある部分を表しているため,これ は因子分析の文脈では共通因子得点に相当しています。同様に,反イメージスコアは独自因 子得点に相当するわけです。したがって,反イメージスコアの相関行列を見ることによって, *14 32 実際に私は SPSS を使わないので,ここに来て修士論文の審査をするまでこの検定の存在を知りませんでし た。
6.6 因子分析をする前に 図6.6に示されるような共通因子では取り切れない「独自因子間の相関」を見ることができる, という考えに至ります。 具体的なやり方としては,まず反イメージスコアの相関行列を以下の式 (Kaiser, 1970) によ って求めます。 𝐐 = 𝑆𝚺−1 𝑆 (6.38) 𝑆 = diag(𝚺−1 )−1/2 なんだかややこしい式ですが,最終的に得られる行列 𝐐 の各要素は,対応する2変数の反イ メージスコアの相関係数=「それ以外の全ての観測変数」の影響を除いた偏相関係数に等し くなっています。つまり,𝐐 の非対角成分が0に近いほど「観測変数間の相関の大部分が他 の観測変数との共通成分で説明できる」ことを意味するため,因子分析を行うことが妥当で あると示されるわけです。 ここまでの背景を説明して,ようやく KMO の登場です。KMO は上で説明したものをまと めて一つの指標としたもので, KMO = 2 ∑ ∑𝑖≠𝑗 𝑟𝑖𝑗 (6.39) 2 + ∑∑ 𝑞2 ∑ ∑𝑖≠𝑗 𝑟𝑖𝑗 𝑖≠𝑗 𝑖𝑗 という式で表されます (Kaiser and Rice, 1974)。ここで 𝑟𝑖𝑗 , 𝑞𝑖𝑗 はそれぞれデータから計算 された相関行列 𝚺 と反イメージ相関行列 𝐐 における変数 𝑖, 𝑗 の相関係数を表しています。 つまり,KMO は「そもそも観測変数間には十分な相関 𝑟𝑖𝑗 がある」かつ「他の変数と無関係 な独自成分同士の相関 𝑞𝑖𝑗 は小さい」という程度を表した指標と言えるでしょう。 なお,KMO は上述のようにデータ全体に対してだけでなく,各変数ごとに計算することも できます。といっても考え方は簡単で,上の式を特定の変数の列のみで計算したら良いだけ です。 KMO𝑖 = 2 ∑𝑗=1,⋯𝑖−1,𝑖+1,⋯,𝐼 𝑟𝑖𝑗 2 +∑ 𝑞2 ∑𝑗=1,⋯𝑖−1,𝑖+1,⋯,𝐼 𝑟𝑖𝑗 𝑗=1,⋯𝑖−1,𝑖+1,⋯,𝐼 𝑖𝑗 (6.40) KMO𝑖 が低い場合,その変数は他の変数との相関がそもそも強くない可能性が高いため,因 子分析から除外することを検討すると良いかもしれません。 R では,psych パッケージに KMO() という関数が用意されています。 KMO(またの名を MSA) 33 1 # cor_polyは先程計算済みのポリコリック相関 2 KMO(cor_poly$rho) 1 Kaiser-Meyer-Olkin factor adequacy 2 Call: KMO(r = cor_poly$rho) 3 Overall MSA = 4 MSA for each item = 0.8 5 Q1_1 Q1_2 Q1_3 Q1_4 Q1_5 Q1_6 Q1_7 Q1_8 Q1_9 Q1_10 6 0.76 0.80 0.78 0.84 0.82 0.80 0.79 0.85 0.79 0.79
6.7 回転の方法 出力のうち,Overall MSA がデータ全体に計算した KMO の値です。また,MSA for each item = 以下に示されているのが,変数ごとに計算した KMO𝑖 の値です。具体的にどの程度 あれば良いか,というのはあまり明確ではないような気もしますが,Kaiser (1974) では,以 下のような基準を示しています。 値 評価 .9 以上 marvelous .8 以上 meritorious .7 以上 middling .6 以上 mediocre .5 以上 miserable .5 未満 unacceptable ということで,Overall MSA の値が 0.5 以上であればとりあえず “acceptable” と言えるよ うにも見えますが,式の定義からすると,KMO が 0.5 になるのは「データにおける相関係 数」と「反イメージ相関係数」が同じ大きさの場合なので,実質的に共通因子が存在しない ようなケース(各変数を独立に乱数生成した場合とか;Shirkey and Dziuban, 1976)に相 当すると考えられます。そういうわけで,実際にデータ分析をする際には,なんとなくです が 0.8 くらいはあったほうが良いのではないか,という気がします*15 。 6.7 回転の方法 6.4.2項では,回転のメカニズムおよびその目的を説明しました。回転するためには,回転後 の因子負荷行列がいい感じに単純構造に近づくように回転行列を指定する必要があります。 そのためには,「いい感じ」の基準を決めなければいけません。 そもそも「単純構造」とは何なのでしょうか。Thurstone (1947) による具体的 (にして理想 的) な単純構造の条件は以下のとおりです。 1. 各項目は少なくとも 1 つの因子に対して因子負荷が 0 であること 2. 各因子は少なくとも因子数 𝑇 個の項目に対して因子負荷が 0 であること 3. 2 つの因子に着目した場合に,一方の因子の要素は 0 で他方の因子の要素は 0 でない 項目があること(=異なる因子を識別可能な項目が存在すること) 4. 因子数が 4 以上の時に,そのうち 2 つの因子に着目した場合に,どちらの要素も 0 で ある項目があること 5. 2 つの因子に着目した場合に,どちらの要素も 0 である項目が少ないこと *15 34 こんなことになっている原因として考えられるのは,KMO には複数の定義があり,取りうる値の範囲が異 なること (Kaiser, 1981),表に示された基準が厳密には KMO に関するものではない(単純構造の程度に関 するものである)こと,などがありそうですがよく分かりません。
6.7 回転の方法 この「単純構造」の条件になるべく近づけるための回転の方法について,数多くの基準が提 案されてきました。fa() 関数で利用可能なもののうち,有名と思われるものについていくつ か紹介します。 【直交回転系】 直交回転で最も有名なのはバリマックス回転でしょう。この方法では,因子 負荷行列に対して 2 𝐼 ⎧ ⎫ { 𝐼 } 1 2 2 2 ∑ ⎨∑ (𝑏𝑡𝑖 ) − (∑ 𝑏𝑡𝑖 ) ⎬ 𝐼 𝑖=1 } 𝑡=1 { ⎩ 𝑖=1 ⎭ 𝑇 (6.41) という計算を行います。最初の ∑ 記号によって,この計算が「列ごとに何かを計算したも のの総和を求めている」ということがわかります。列ごとに求めているものは,因子負荷量 の二乗和の分散のようなものです。ある変数 𝑥 の分散を求める方法として, 𝜎𝑥2 = 2 1 1 ∑ 𝑥2 − ( ∑ 𝑥) 𝑛 𝑛 (6.42) というものがありました。すこし調整の項があったりしますが,基本的にはこれと同じよう なことをしている(𝑥 のところに 𝑏𝑡𝑖 が入っている)わけです。variance を max にするので vari-max という名前なわけですね。直感的には,因子負荷の分散の二乗和が最大になるよう にしてやると,それぞれの因子について「負荷が高い観測変数」と「負荷が低い観測変数」が あることになるため,単純構造の条件 2 に近づきそうです。 【斜交回転系】 斜交回転で最も有名なのはプロマックス回転だと思います。プロマックス回 転ではプロクラステス回転をベースに回転を行うため,まずはプロクラステス回転の説明を しておきます。プロクラステス回転では,回転のゴールをただの単純構造ではなく,仮説や 理論から導き出されるターゲット行列に設定し,ターゲット行列とのズレが最小になるよう に回転を行います。ターゲット行列には,例えば先行研究で観測された因子負荷行列などを 設定することができるのですが,そのような行列を探してきて設定するのは結構大変なこと です。そこでプロマックス法では,まずバリマックス法などの直交回転による回転を行い, 得られた因子負荷行列を調整することでターゲット行列をデータから作成します。具体的に は,バリマックス回転で得られた因子負荷行列について,値が大きい要素はより大きく,値 が小さい要素はより小さくすることで,バリマックス解よりもさらに「理想に近い」行列を 作成しています。プロマックス回転は,それ自体が「単純構造」に向かって何らかの基準を 最適化する方法ではありません。バリマックス回転をベースに解を求めているということは, 例えばそもそもバリマックス回転でうまく単純構造が出なかった場合には失敗する可能性が あります。計算量が少ないことなどから昔はよく使われていた方法のようです。 fa() 関数のデフォルトはオブリミン回転と呼ばれる方法です。実際には,コーティミン回転 と呼ばれる方法を行っているような気がします*16 。コーティミン回転では,以下の基準を最 小化します。 𝑇 𝑇 𝐼 2 2 ∑ ∑ {∑ 𝑏𝑡𝑖 𝑏𝑠𝑖 } 𝑡=1 𝑠≠𝑡 *16 35 (6.43) 𝑖=1 オブリミン回転は,コーティミン回転やコバリミン回転などの回転法を一般化したものなので,コーティミ ン回転はオブリミン回転の特殊なケースなのです。
6.7 回転の方法 最初の 2 つの ∑ 記号によって,因子のペアごとに何かしらの計算を行い総和を求めている ことがわかります。その中身は同じ項目 𝑖 に対する因子負荷の二乗の積です。この式は,同 じ項目に因子 𝑡 と因子 𝑠 が両方とも高い因子負荷を持っていた場合にのみ大きな値になりま す。したがって,コーティミン基準を最小化するということは特に単純構造の条件 3 を目指 している,と言えるでしょう。 ここで紹介した以外にも,因子の回転法はかなりたくさんあります。その中でどれを使えば よいかについては,特に決まっていない(個人の好みな)気がします。強いて言うなら,大 抵の場合は斜交回転をしておけば良いでしょう。方法によって,因子間相関が高く出やすい ものがあれば低く出やすいものもあったりするので,いろいろな方法を試してみていい感じ の因子負荷行列になる回転法を選べば良いのではないかと思います。 fa() 関数では,引数 rotate によって回転法を指定することができます。これまでは指定し ていなかったので rotate="oblimin" になっていたわけですが,試しに直交回転であるバ リマックス回転で因子負荷を求めてみましょう。 直交回転 1 res_fa_varimax <- fa(dat[, cols], nfactors = 2, cor = "poly", rotate = "varimax") 2 res_fa_varimax 1 Factor Analysis using method = 2 Call: fa(r = dat[, cols], nfactors = 2, rotate = "varimax", cor = minres "poly") 3 Standardized loadings (pattern matrix) based upon correlation matrix 4 MR1 MR2 h2 u2 com 5 Q1_1 -0.01 0.45 0.21 0.79 1.0 6 Q1_2 0.13 0.73 0.55 0.45 1.1 7 Q1_3 0.12 0.80 0.65 0.35 1.0 8 Q1_4 0.24 0.51 0.32 0.68 1.4 9 Q1_5 0.15 0.64 0.43 0.57 1.1 10 Q1_6 0.60 0.06 0.36 0.64 1.0 11 Q1_7 0.66 0.09 0.45 0.55 1.0 12 Q1_8 0.58 0.12 0.35 0.65 1.1 13 Q1_9 0.71 0.14 0.53 0.47 1.1 14 Q1_10 0.61 0.16 0.39 0.61 1.1 15 16 MR1 17 SS loadings 2.12 2.11 18 Proportion Var 0.21 0.21 19 Cumulative Var 0.21 0.42 20 Proportion Explained 0.50 0.50 21 Cumulative Proportion 0.50 1.00 22 23 36 MR2 Mean item complexity = 1.1
6.8 因子の推定法 24 Test of the hypothesis that 2 factors are sufficient. 25 26 df null model = 27 df of Square = 45 with the objective function = 2.78 with Chi 6749.62 the model are 26 and the objective function was 0.22 28 29 The root mean square of the residuals (RMSR) is 30 The df corrected root mean square of the residuals is 0.04 0.05 31 32 The harmonic n.obs is prob < 33 The total n.obs was prob < 2432 with the empirical chi square 367.22 with 5.9e-62 2432 with Likelihood Chi Square = 540.84 with 1.2e-97 34 35 Tucker Lewis Index of factoring reliability = 36 RMSEA index = 37 BIC = 38 Fit based upon off diagonal values = 0.98 39 Measures of factor score adequacy 0.09 0.867 and the 90 % confidence intervals are 0.084 0.097 338.13 40 MR1 MR2 41 Correlation of (regression) scores with factors 0.88 0.9 42 Multiple R square of scores with factors 0.77 0.8 43 Minimum correlation of possible factor scores 0.54 0.6 出力を見ると,直交回転なので因子間相関が出力されていません。また res_fa_cat(オブ リミン解)と比べると,各項目について負荷量の低い側の因子の負荷量がやや大きめになっ ていることがわかります。このように,基本的には斜交回転のほうが単純構造に近い結果を 得ることができるため,よっぽどの理由が無ければ斜交回転にしておくのがおすすめです。 6.8 因子の推定法 因子分析では,固有値分解によって因子負荷行列(の初期解)を求める,と言うのはすでに 説明したとおりです。この際分解されるものは,観測変数の相関行列 𝚺 そのものではなく, 𝚺 = 𝐁⊤ 𝐁 + 𝚿 (6.20) と分解した上での 𝐁⊤ 𝐁 の部分でした。ただし,誤差共分散行列 𝚿 の中身はもちろんわかり ません。当然 𝚿 の値次第で 𝐁⊤ 𝐁 の値も変わるため,得られる因子負荷も変わってきます。 このままでは解が一つに決まりません。では実際には,どのようにしての推定を行っている のでしょうか。 𝐁 の推定に際しては,まず大前提として 𝚿 の非対角成分はゼロ,すなわち共通因子の影響 を取り除けば誤差共分散は全て消えるという仮定をおきます。もし本当に誤差共分散が全て ゼロなのであれば,𝚺 = 𝐁⊤ 𝐁 + 𝚿 という分解は完全に一致します。ですが,実際には少な からず誤差共分散が残っています。その分だけ 𝚺 と 𝐁⊤ 𝐁 + 𝚿 の間にはズレが生じること になります。そこで,何らかの基準のもとで最適な分解のかたちを推定する必要があります。 37
6.8 因子の推定法 ここでもいくつかの選択肢があり,場合によって使い分けが必要になります。 6.8.1 主因子法 一昔前のスタンダードだったのが主因子法です。主因子法では,以下の手順によって各変数 の共通性の値を求めます。 1. 共通性の初期値を適当に決める 2. データの相関行列の対角成分に共通性の推定値を入れた行列を固有値分解する 3. 2. の結果をもとに共通性を推定する 4. 共通性の値が変わらなくなる(計算が収束する)まで 2-3. を繰り返す 実は,数学的には主因子法は後述する最小二乗法と同じらしいです。 6.8.2 最小二乗法 fa() 関数のデフォルトになっている方法です。各変数の独自性(𝚺 の対角成分)は「共通性 の残りカス」で調整してしまえば良いので無視できます。それ以外の非対角成分の部分につ いて,𝐁⊤ 𝐁 による近似と 𝚺 の差の二乗和が最小になるように「共通性の推定」と「固有値 分解」を繰り返すのが最小二乗法です。 6.8.3 最尤法 共通因子と独自因子が多変量正規分布 𝑁 (𝛍, 𝐒) = ([ 𝟎𝐶 𝚽 ] , [ 𝐶𝐶 𝟎𝑈 𝟎𝑈𝐶 𝟎𝐶𝑈 ]) 𝐈𝑈𝑈 (6.44) に従うと仮定します。なお,添字の 𝐶 は共通因子に関する部分,𝑈 は独自因子に関する部分 を表しています。したがってここでは • 共通因子,独自因子ともに平均はゼロ • 共通因子同士の相関は認める • 独自因子同士の相関は認めないので,対角成分が 1 の単位行列 • 共通因子と独自因子の間は無相関 ということが仮定されています。 この状態で,標準化された観測変数 𝐲 も多変量正規分布に従うと仮定します。 𝐲 ∼ 𝑁 (𝟎, 𝐁⊤ 𝐁 + 𝚿) (6.45) あとは普通の最尤法に従って,データと最も整合的な 𝐁 および 𝚿 を求めたら良いわけです。 最尤法は,最小二乗法などと比べると,データの多変量正規性を仮定しているなどの理由か ら,統計的に高度な方法となっています。そのため,計算がうまくいく場合には最小二乗法 などよりも最尤法のほうが望ましいとされています。一方で,多変量正規性の仮定などが推 38
6.8 因子の推定法 定に悪影響を及ぼす場合があります。例えば「データが少ない」「因子数が必要以上に多い」 といったときに共通性の推定値が 1 を超える (Heywood case と呼ばれる) といった不適解 が出てしまうこともしばしばあります。そのような場合には,最小二乗法などの別の方法に よって推定を行うのが良いかもしれません。 fa() 関数では,fm という引数によって因子の推定法を指定することが出来ます。試しに最 尤法でやってみましょう。 因子の推定法を指定 1 res_fa_ml <- fa(dat[, cols], nfactors = 2, cor = "poly", fm = "ml") 2 res_fa_ml 1 Factor Analysis using method = 2 Call: fa(r = dat[, cols], nfactors = 2, fm = "ml", cor = "poly") 3 Standardized loadings (pattern matrix) based upon correlation matrix 4 ML2 ML1 h2 ml u2 com 5 Q1_1 -0.08 0.47 0.20 0.80 1.1 6 Q1_2 0.01 0.71 0.51 0.49 1.0 7 Q1_3 -0.03 0.82 0.66 0.34 1.0 8 Q1_4 0.16 0.49 0.32 0.68 1.2 0.67 0.46 0.54 1.0 9 Q1_5 0.03 10 Q1_6 0.61 -0.03 0.36 0.64 1.0 11 Q1_7 0.66 -0.02 0.43 0.57 1.0 12 Q1_8 0.58 0.02 0.34 0.66 1.0 13 Q1_9 0.73 0.00 0.54 0.46 1.0 14 Q1_10 0.62 0.05 0.41 0.59 1.0 15 16 ML2 ML1 17 SS loadings 2.12 2.11 18 Proportion Var 0.21 0.21 19 Cumulative Var 0.21 0.42 20 Proportion Explained 0.50 0.50 21 Cumulative Proportion 0.50 1.00 22 23 With factor correlations of 24 ML2 25 ML2 1.00 0.34 26 ML1 0.34 1.00 ML1 27 28 Mean item complexity = 29 Test of the hypothesis that 2 factors are sufficient. 1 30 31 df null model = Square = 32 33 39 df of 45 with the objective function = 2.78 with Chi 6749.62 the model are 26 and the objective function was 0.22
6.9 因子数の決め方 34 The root mean square of the residuals (RMSR) is 35 The df corrected root mean square of the residuals is 0.04 0.05 36 37 The harmonic n.obs is prob < 38 The total n.obs was prob < 2432 with the empirical chi square 377.19 with 5.6e-64 2432 with Likelihood Chi Square = 532.87 with 5.4e-96 39 40 Tucker Lewis Index of factoring reliability = 41 RMSEA index = 42 BIC = 43 Fit based upon off diagonal values = 0.98 44 Measures of factor score adequacy 0.09 0.869 and the 90 % confidence intervals are 0.083 0.096 330.16 45 ML2 ML1 46 Correlation of (regression) scores with factors 47 Multiple R square of scores with factors 0.79 0.82 48 Minimum correlation of possible factor scores 0.59 0.64 0.89 0.91 今回のようにデータがきれいな場合には,最尤法でも問題なく解を得ることが出来ます。こ のような場合には,最尤法での結果を報告するほうが良いかもしれません。 6.9 因子数の決め方 因子分析では因子数を分析者が自由に選ぶことが出来ます。これまでは,dat の最初の 10 項目を使って,因子数を 2 に決め打ちして因子分析を行っていました。あらかじめ因子数が 分かっていたり,理論的に目処が立っている場合にはその因子数で分析を行えば良いのです が,尺度をボトムアップに作成するような場合や,自分が集めたデータでも先行研究と同様 の因子構造が得られるかを検証する場合には,因子数はわからないことが多いです。ここで は,統計的に因子数を決定するためのいくつかの基準を紹介します。 6.9.1 カイザー(・ガットマン)基準 因子分析は,データの相関行列(正確には独自性を除いた 𝐁⊤ 𝐁)を固有値分解することで 因子負荷行列を求めていました。このとき,各因子に対する因子負荷は固有ベクトルを固有 値倍したもの(√𝜆𝑡 𝐱𝑡⊤ ) として算出されていました。したがって,固有値こそが因子負荷の合 計を決めており,ひいては分散説明率(共通性)の合計を規定しているということが言えま す。実は,𝐼 × 𝐼 相関行列を固有値分解して得られる 𝐼 個の固有値は,合計が必ず𝐼 になると いう性質があります。つまり,固有値は 𝐼 個の項目に対する分散説明率を再配分した値であ る,とも言えそうです。イメージ的には,例えば第一因子に対する固有値(固有値の中の最 大値)が 3 だとすると,その因子は一つで合計 3 項目分の分散を説明しているということで す。1 個で 3 項目分の説明力を持っているということは,結果的に次元数を減らすことに役 立っていると言えそうです。 このように考えると,固有値が 1 より小さい因子は 1 項目分の情報も持っていないことにな るため,因子としてはザコと言えます。また,独自因子は 1 項目のみに影響を与えている因 40
6.9 因子数の決め方 子なので,固有値を計算するとしたら 1 になります。この意味でも,固有値が 1 以下の因子 はやはり要らなさそうです。この「固有値が 1 以上の因子のみを使う」という基準をカイザ ー(・ガットマン)基準と呼びます。カイザー基準は(未だに!)よく用いられている方法で すが,その問題点も長らく指摘されています。言ってしまえば「カイザー基準は大雑把であ り,それほどあてにならない」(堀, 2005) のです。 6.9.2 スクリーテスト 固有値を順に並べたものをスクリープロット (scree plot) と呼びます。psych パッケージで は scree() 関数によって図6.16のように簡単にプロット可能です。 スクリープロット scree(dat[, cols]) PC FA 0.5 1.0 1.5 2.0 2.5 3.0 Scree plot 0.0 Eigen values of factors and components 1 2 4 6 8 10 factor or component number 図6.16 スクリープロット なお,結果は FA の方を見てください。PC という方は主成分分析での結果です。ほとんどの 場合主成分分析のプロットの方が上に来ます。なぜそうなるかは6.11節を見てください。ま た,デフォルトでは Y 軸が 1 のところに実線が引かれています。これは先程のカイザー基準 を表しています。つまり今回の場合カイザー基準的には 2 因子が妥当ということになります。 scree() 関数にデータフレームを与えた場合には,そのデータフレームに対してピアソンの 積率相関を求めて固有値を計算します。先程紹介したポリコリック相関に基づく固有値を見 たい場合には,事前にポリコリック相関を計算しておき,それを scree() 関数に与える必要 があります。 ポリコリック相関でスクリープロット 41 1 # psych::polychoric()関数も使える 2 cor_poly <- polychoric(dat[, cols]) 3 # polychoric()は相関行列以外にもlistで色々返してくる 4 # listの中のrhoがポリコリック相関行列の本体 5 scree(cor_poly$rho)
6.9 因子数の決め方 0.5 1.0 1.5 2.0 2.5 3.0 PC FA 0.0 Eigen values of factors and components 3.5 Scree plot 2 4 6 8 10 factor or component number 図6.17 ポリコリック相関でスクリープロット あまり違いはないように見えますが,最初の因子を比べると,わずかに図6.17のほうが 図6.16よりも大きな固有値になっているはずです。 このスクリープロットを用いて行うスクリーテスト (scree test) と呼ばれる因子数選択では, プロット的に傾きがなだらかになる手前の因子数を選択する,という方法です。図6.16や 図6.17では,いずれも最初の 2 因子は傾きが強く,第 3 因子以降はなだらかになっています。 したがってスクリーテスト的にも因子数は 2 が妥当,という判断ができます。 イメージ的には固有値分解では,まず第一固有値が「全てに共通する成分」をごっそり持っ ていきます。第 2 固有値以降は,その残りの中から共通成分を取り出していきます。したが って,あるところからはもう絞りカスしか残っていない状態になり,固有値的には下の方で 大差ない状態になってしまいます。スクリーテストでは,この「下の方で大差ない状態」と くらべて如実に大きな固有値にはなにか意味があると考えて,その因子数を採用しているわ けです。 ただ,今回のデータのように因子数がハッキリとしたきれいなデータではない場合には,こ れほどハッキリと変化点を見出すことが出来ないケースが結構多いです。そのような場合に 客観的に因子数を決めづらかったりもするので,やはりまだアバウトな方法と言えそうです。 6.9.3 平行分析 カイザー基準では一律で 1 以上の固有値のある因子には意味がある,という判断をしていま した。平行分析 (parallel analysis) では,代わりに「乱数での固有値」との比較を行います。 完全に乱数のデータでは本来相関は無いはずですが,標本誤差の関係などで相関が完全にゼ ロになることはありません。相関が完全にゼロではない限りは,固有値分解を行うと第一固 有値の値は 1 よりも少しだけ大きくなります。試しにやってみましょう。 乱数で固有値を出してみる 42 1 # 1000行10列の一様乱数行列 2 dat_rnd <- matrix(runif(10 * 1000), ncol = 10) 3 eigen(cor(dat_rnd))
6.9 因子数の決め方 1 eigen() decomposition 2 $values 3 [1] 1.1531292 1.0889700 1.0765136 1.0248099 1.0075078 0.9926984 4 [7] 0.9533219 0.9367722 0.8834080 0.8828690 5 6 $vectors 7 [,1] [,2] 8 [1,] 0.112436984 -0.37121765 9 [2,] 0.484608727 [,3] [,4] [,5] 0.47341966 -0.20738992 0.103139623 0.05575061 -0.01770245 0.32996851 0.067515555 0.28649476 0.500656646 10 [3,] -0.228692064 -0.48815828 0.20308039 11 [4,] -0.232490192 0.40499060 0.29262518 -0.48025212 -0.123849547 12 [5,] -0.031949135 0.16118059 -0.61178243 0.08781265 13 [6,] -0.557884289 0.24802394 0.25838357 -0.001868273 14 [7,] 15 [8,] -0.240054398 -0.05764926 -0.32196831 -0.47149177 16 17 0.15782882 0.161289155 0.007055743 -0.48064664 -0.17069328 -0.39982071 -0.395600998 0.598441643 [9,] -0.015405949 -0.27818657 -0.32376559 -0.02988766 -0.222218221 [10,] 0.525349112 18 0.23988743 0.10898482 -0.27903649 0.358854999 [,6] [,7] [,8] [,9] [,10] 0.07380739 0.68971166 0.07746406 0.28930807 0.03051726 0.09815149 -0.72312285 0.13720316 0.15336519 19 [1,] 20 [2,] -0.27094157 21 [3,] -0.12276860 -0.09821372 -0.09707071 -0.54777825 -0.03495240 22 [4,] -0.32906603 0.16989052 -0.25314902 -0.44589225 23 [5,] 0.62789616 -0.01514297 -0.32785652 -0.02723450 24 [6,] -0.01486085 25 [7,] 26 [8,] -0.09714458 -0.14690604 -0.19076435 0.40467304 [9,] -0.79986315 0.04312236 -0.23139465 27 28 0.24956474 0.14087695 -0.29134682 0.22420762 0.27164195 -0.59735607 0.29098840 -0.10511832 -0.46338991 -0.15375669 -0.29889905 0.15216131 [10,] -0.03729309 -0.07919072 0.22122246 0.15420830 0.11141185 -0.16749641 -0.63231066 確かに第一固有値の値は 1.153 となり,1 よりも大きな値になりました。というか,もしこ のデータにカイザー基準を当てはめると,5 因子もあることになってしまいます。 実際には乱数には因子など無いはずなので,平行分析では乱数で得られた固有値よりも小さ い値は因子として無意味と判断します。R では fa.parallel() という関数があるのでやっ てみましょう。 平行分析 1 # この関数には引数corが用意されている 2 fa.parallel(dat[, cols], cor = "poly") 1 Parallel analysis suggests that the number of factors = number of components = 2 and the 2 青い線が実際のデータでの固有値を,赤い点線が乱数での固有値を表しています。 (乱数を使 うのでやるたびに微妙に結果が変わるかもしれませんが)Parallel analysis suggests 43
3.5 Parallel Analysis Scree Plots 0.5 1.0 1.5 2.0 2.5 3.0 PC Actual Data PC Simulated Data PC Resampled Data FA Actual Data FA Simulated Data FA Resampled Data 0.0 eigenvalues of principal components and factor analysis 6.9 因子数の決め方 2 4 6 8 10 Factor/Component Number 図6.18 平行分析 that the number of factors = 2 ということで,因子数 2 が提案されました*17 。 6.9.4 累積寄与率 分散説明率の観点から言えば,共通因子によってある程度の割合は説明されていないと, 「結 局独自因子(よくわからないもの)で回答が決まってるんでしょ」となってしまいます。因 子の数を増やせば累積寄与率は必ず上昇するので,因子数の最低条件として,累積寄与率が 一定の値以上になっていない場合には因子の数を増やすというアプローチを取ることもある ようです。必要な累積寄与率の基準は一律で決まっているわけではないですが,一説による と最低 50% はほしいと言われています (小杉, 2018)。 6.9.5 VSS 因子分析では単純構造が望ましい状態であり,なるべくこれを目指して回転が行われます。 VSS (Very Simple Structure) 基準 (Revelle and Rocklin, 1979) では,得られた結果を強 制的に単純構造にしたときのデータとの整合性を見ます。 完全な単純構造では,各項目は 1 つの因子にのみ負荷を持っているはずです,そこで VSS で は最も高い負荷を持っている因子以外の負荷をゼロにした因子負荷行列を考えます。先程出 した res_fa_cat で言えば,表6.5のような状態にする,ということです。 ̃ そして,この因子負荷行列(𝐁̃ VSS とします) に対して,𝐁̃ ⊤ VSS 𝚽𝐁VSS によって相関行列を復 元します。これをデータの相関行列 𝚺 と比較して,差が小さければ「単純構造とみなした時に ̃ データと整合性がある」ということになります。そこで,まず要素ごとの差 𝐁̃ ⊤ VSS 𝚽𝐁VSS − 𝚺 をとり,この行列の要素の二乗和を 𝑀 𝑆VSS とします。同様に,データの相関行列 𝚺 の要素 の二乗和を 𝑀 𝑆 とすると,𝑉 𝑆𝑆 は 𝑉 𝑆𝑆 = 1 − *17 44 𝑀 𝑆VSS 𝑀𝑆 ここでも number of components は主成分分析の話なので無視して OK です。 (6.46)
6.9 因子数の決め方 表6.5 very simple structure MR2 MR1 Q1_1 0.000 0.478 Q1_2 0.000 0.743 Q1_3 0.000 0.813 Q1_4 0.000 0.486 Q1_5 0.000 0.642 Q1_6 0.617 0.000 Q1_7 0.680 0.000 Q1_8 0.580 0.000 Q1_9 0.720 0.000 Q1_10 0.603 0.000 となります。この計算を,様々な因子数で行い,最も 𝑉 𝑆𝑆 の値が小さくなったところが「単 純構造として解釈した時に最もデータと整合的な因子数」ということになるわけです。 psych パッケージには,VSS() という関数が用意されています。 VSS 1 # 引数corが用意されています 2 VSS(dat[, cols], cor = "poly") 0.8 2 1 4 3 2 1 4 3 2 4 3 2 2 1 1 1 6 7 8 4 3 2 4 3 1 0.2 0.4 0.6 1 1 3 2 0.0 Very Simple Structure Fit 1.0 Very Simple Structure 1 2 3 4 5 Number of Factors 図6.19 VSS のプロット 1 2 Very Simple Structure 3 Call: vss(x = x, n = n, rotate = rotate, diagonal = diagonal, fm = fm, 4 VSS complexity 1 achieves a maximimum of 0.72 with 2 factors 6 VSS complexity 2 achieves a maximimum of 0.81 with 4 factors 7 45 n.obs = n.obs, plot = plot, title = title, use = use, cor = cor) 5
6.9 因子数の決め方 8 The Velicer MAP achieves a minimum of 0.04 9 BIC achieves a minimum of 10 10.65 with 5 with 2 factors factors Sample Size adjusted BIC achieves a minimum of 26.53 with 5 factors 11 12 Statistics by number of factors 13 vss1 vss2 map dof chisq prob sqresid fit RMSEA BIC SABIC 14 1 0.61 0.00 0.059 35 2779.5 0.0e+00 7.1 0.61 0.180 2507 15 2 0.72 0.80 0.035 26 540.8 1.2e-97 3.6 0.80 0.090 338 421 16 3 0.66 0.80 0.060 18 308.0 1.1e-54 3.1 0.83 0.081 168 225 17 4 0.63 0.81 0.084 11 130.7 1.3e-22 2.5 0.86 0.067 45 80 18 5 0.60 0.79 0.127 5 49.6 1.7e-09 2.3 0.88 0.061 11 27 19 6 0.47 0.73 0.186 1.9 0.89 NA NA 20 complex eChisq 0 1.1 SRMR eCRMS NA eBIC 21 1 1.0 4210.1 0.1387 0.157 3937.2 22 2 1.1 367.2 0.0410 0.054 164.5 23 3 1.2 205.1 0.0306 0.048 64.8 24 4 1.4 78.2 0.0189 0.038 -7.6 25 5 1.7 28.3 0.0114 0.034 -10.7 26 6 27 1.8 0.4 0.0014 NA 2618 NA NA [ reached 'max' / getOption("max.print") -- omitted 2 rows ] complexity 1 というのは,𝐁̃ VSS として各項目が 1 つの因子のみに負荷している状態,つ まり完全な単純構造を仮定した状態での結果です。同様に,各項目が上位 2 つの因子にの み負荷している状態を 𝐁̃ VSS に仮定した場合の結果が complexity 2 となっています*18 。 図6.19は各 complexity における VSS の値です。1 の線を見ると,Number of Factors が 2 のところで最も高い値になっています。これが出力されているわけです。ということで VSS(complexity 1) 的には 2 因子が妥当だという判断がされました。 VSS() では他にも,次に示す MAP を始めとした様々な指標が出力されます。論文などで報 告されることがあるものをざっと確認しましょう。 6.9.6 MAP MAP (Minimum average partial; 最小平均偏相関) では,データの相関行列のうち共通因 子の影響を取り除いた偏相関の部分に着目します。いわばこれは誤差どうしの相関であり, もしこれが全てゼロであれば,独自因子の間に別の共通因子が存在していないことになる (図6.6を参照)ため,因子数として妥当だと判断できます。計算自体はさほど難しくなく,ま ずはふつうに 𝑘 因子で因子分析の初期解を求めます。つまり,データの相関行列 𝚺 を固有値 分解して,得られた固有ベクトル [𝐱1 𝑘 個のみを用いた 𝐗𝑘 = [𝐱1 𝐱2 𝐱2 ⋯ 𝐱𝑘−1 𝐱𝑘 𝐱𝑘+1 ⋯ 𝐱𝐼 ] のうち最初の ⋯ 𝐱𝑘 ] と,対応する 𝑘 個の固有値を対角成分においた 行列 𝚲𝑘 を用いて 𝚺 を近似します。 *18 46 一応内部的には complexity 8 まで計算されていますが,それはもう単純構造では無いのでは?という気が するので大抵の場合は complexity 2 くらいまで見ておけば十分でしょう。
6.10 𝜔 係数 𝚺 ≅ 𝐗 𝑘 𝚲𝑘 𝐗 ⊤ 𝑘 (6.47) 𝐏 = 𝚺 − 𝐗 𝑘 𝚲𝑘 𝐗 ⊤ 𝑘 (6.48) そして,𝚺 との差を取ります。 1 1 このとき,diag(𝐏)− 2 𝐏diag(𝐏)− 2 の非対角成分が,𝑘 個の共通因子の影響を取り除いた後 の各変数間の偏相関係数を表します。したがって,MAP ではこの非対角成分の要素の二乗 の平均値を計算します。これを 𝑘 の数を変えて行い,最も小さい MAP を出した 𝑘 が妥当な 因子数 𝑇 だという判断ができるわけです。 6.9.7 その他 他にも,VSS() 関数では次のような基準による評価を行っています。 • 情報量規準 BIC が最も小さくなる因子数 • 共通因子の影響を取り除いた残差相関行列に対する 𝜒2 統計量 • モデルの適合度指標 (RMSEA, SRMR など:SEM のときに説明します) 6.9.8 最後は解釈可能性 ここまで様々な基準によって最適な因子数を推定しました。ですがその結果は一貫していま せん。今回の場合,多くの指標は(データの設計通り)因子数 2 を提示しましたが,一方で 3 や 5 といった数を出したものもあります。それぞれの指標・基準は異なる視点からの評価 なので,このように結果が変わるのは当然のことです。その上で,実際に因子数を決める場 合にはどうしたら良いのでしょうか。 因子分析で重要なのは,解釈可能性だと考えています。因子の名前を決める際には,分析者 が自身の知識と経験をフル活用して項目の共通点をラベリングしなければいけません。そし て,統計的に最適と思われる因子数が提案されたとしても,その結果は標本誤差の影響を受 けている可能性がゼロとは言えません。したがって,探索的に因子数を決定する際には, 1. まず VSS() や fa.parallel() によって因子数にアタリをつける 2. 提案された因子数(多くの場合は VSS() で出力される MAP と BIC の間)に関して実 際に因子分析を行う 3. 得られた結果をもとに因子の解釈を行う 4. 最もしっくり来る解釈が得られた因子数および因子構造を採用する という手順で行われることが多いでしょう。もしも全ての統計的な基準が指している因子数 において解釈が難しいという事態に直面したら,そのときは泣いてください。 6.10 𝜔 係数 信頼性係数のところでは,𝛼 係数はあまり良くない指標であることをお話しました。ここで は,因子分析の考え方に基づく,より良い(とされている)信頼性係数の 𝜔 係数を紹介し ます。 47
6.10 𝜔 係数 古典的テスト理論 (𝑥 = 𝑡 + 𝑒) における信頼性係数は 𝜌= 𝜎𝑡2 𝜎2 = 1 − 𝑒2 2 𝜎𝑥 𝜎𝑥 という式で定義されていました。ここで古典的テスト理論のモデルを因子分析 (𝐲𝑖 = 𝐟 𝑏𝑖 +𝛆𝑖 ) に当てはめて考えてみます。真値 𝑡 にあたる部分が 𝐟 𝑏𝑖 ,誤差 𝑒 に当たる部分が 𝛆𝑖 だと考え ると,信頼性係数 𝜌 は,観測変数 𝑦 の分散のうち 𝐟 𝑏𝑖 の分散説明率になります。また,共通 の因子から影響を受けている 2 変数の相関係数は,各変数への因子負荷の積 (𝑏𝑖 𝑏𝑗 ) で表され ました。さらに,観測変数 𝑦 の共通性(=分散説明率)は,その変数に対する因子負荷の二 乗和で表されていました。これらを総合すると,各項目の因子負荷が大きいほど,内的一貫 性が高いということが言えそうです。 𝜔 係数はこの考え方に基づき,一因子モデルであれば以下の式によって係数を求めます(𝑢 は unidimensional の頭文字)。 𝐼 𝜔𝑢 = 2 (∑𝑖=1 𝑏𝑖 ) 𝜎𝑦2𝑡𝑜𝑡𝑎𝑙 𝐼 2 (∑𝑖=1 𝑏𝑖 ) = 2 𝐼 𝐼 (6.49) (∑𝑖=1 𝑏𝑖 ) + ∑𝑖=1 𝜎𝜀2𝑖 なお,𝜎𝑦2𝑡𝑜𝑡𝑎𝑙 は合計点の分散を表しています。 𝜔 係数の基本的な考え方はこれだけなのですが,実際にはいろいろな 𝜔 を考えることができ ます。 Composite Reliability (CR) と Average Variance Extracted (AVE) (主にマーケティングの分野では?)妥当性の収束的証拠の一つとして合成信頼性 (Composite Reliability; CR) や Average Variance Extraced (AVE) と呼 ばれる指標 (Fornell and Larcker, 1981; Bagozzi and Yi, 1988) が用いられているよ うです。これらの指標の考え方は,結局のところ 𝜔 係数と同じで,各観測変数の分散 に対して,それらの変数から抽出された共通因子による説明力の割合が大きければ, それらの変数は確かに同じ構成概念を測定しているだろう,と考えているわけです。 実際に,CR の式は上で紹介した 𝜔𝑢 と全く同じです (Padilla and Divers, 2016)。 これに対して AVE は,もう少し直接的に「観測変数の分散の合計のうちどれだけの 割合が共通因子で説明できるか」を求めているようです。したがって,AVE の計算 式は CR(𝜔𝑢 )と似て非なるものになっています。 𝐼 AVE = ∑𝑖=1 𝑏𝑖2 𝐼 𝐼 ∑𝑖=1 𝑏𝑖2 + ∑𝑖=1 𝜎𝜀2𝑖 考え方の違いは名前にも表れているのですが, CR(𝜔𝑢 ) 複数の観測変数の集合体 (composite) としての和得点における真値の分 散説明率を考えている。複数の観測変数が集まった場合,誤差が無相関であれ ばそれぞれ打ち消し合うことが期待されるため,変数の数が増えるほどその値 は大きくなりやすい(これが因子負荷を先に全て足してから 2 乗するあたりに 現れている)。 AVE 観測変数一つ一つにおける平均的な真値の分散説明率を考えている。そのため 48
6.10 𝜔 係数 観測変数が増えるほどその値は大きくなりやすいとは限らず,基本的に CR よ りは低い値が出る。しかし観測変数そのものにおける信頼性を評価していると いう点で,CR とは異なる価値を持っている。 という感じです。たぶん。 6.10.1 Bifactor な 𝜔 例えば学力を考えてみると,ある科目で良い点が取れる生徒は他の科目でも良い点を取るこ とが多いでしょう。因子分析的に言えば,これは全ての科目が一つの因子の影響を受けてい る,という見方が出来ます。ですが,果たして国語と数学と英語と…様々な科目の学力をた ったひとつの因子で説明できるのでしょうか。図6.20のように,特に知能などを扱う領域で は「総合的なもの」である G 因子と,特定の項目にのみ影響を与える因子がそれぞれ存在す るという双因子 (Bifactor) モデルというものを考えることがあります。 図6.20 双因子モデル 双因子モデルに基づけば,観測変数は G 因子による部分が追加されて 𝐲𝑖 = 𝐠𝑏𝑖𝑔 + 𝐟 𝑏𝑖 + 𝛆𝑖 (6.50) となります。こうなると,「信頼性」に複数の解釈が生まれます。 まず,図6.20の 5 科目のテストに共通している部分のみ(全体としての内的一貫性)を考え ると,これは G 因子による分散説明率のみを取り上げる必要があります。したがって, 𝐼 𝜔ℎ = (∑𝑖=1 𝑏𝑖𝑔 ) 2 (6.51) 𝜎𝑦2𝑡𝑜𝑡𝑎𝑙 という 𝜔 係数が考えられます (ℎ は hierarchical の頭文字)。 一方で,特定の項目にのみ影響を与えるものを含めた「共通因子によって説明可能な部分」 (どちらかというと元の信頼性の定義に近い)を考えると, 𝐼 𝜔𝑡 = 49 2 𝐼 (∑𝑖=1 𝑏𝑖𝑔 ) + (∑𝑖=1 𝑏𝑖 ) 𝜎𝑦2𝑡𝑜𝑡𝑎𝑙 2 (6.52)
6.10 𝜔 係数 という 𝜔 係数が考えられます(𝑡 は total の頭文字)。 どの 𝜔 係数も間違いではありません。ただ視点が異なるだけなのです。その上でどのように 使い分けるかですが,これは仮説での因子構造によります (Flora, 2020)。心理尺度を作成す る場合には,尺度全体として一つの構成概念を反映していると考えるでしょう。もし一因子 性が満たされている場合には,最初に紹介した 𝜔𝑢 を用いるのが良いでしょう。 もし • 因子構造として多因子を仮定しているが,尺度得点としては全項目の合計点を利用し ようと考えている場合 • もともと一因子で尺度を作ろうとしていたけれども結果的に一因子性が満たされなく なってしまった場合 には,適切な因子構造のもとで 𝜔 を計算する必要があります。つまり特定の項目にのみ影響 を与える因子に加えて G 因子に相当するものを用意し,この G 因子の分散説明率をもって 尺度全体の内的一貫性を求めるべき,ということで 𝜔ℎ を使用すべきかもしれません。 なお,尺度の内的一貫性という観点で言えば,𝜔𝑡 を使うことは基本的には無いと思っていて よい気がします。 psych::omega() という関数は,𝜔 係数を求めることができる関数です。 omega_u を求める 1 # 引数poly=TRUEにするとポリコリック相関でやってくれる 2 omega(dat[, cols], nfactors = 1, poly = TRUE) 1 Warning in schmid(m, nfactors, fm, digits, rotate = rotate, n.obs = 2 n.obs, : Omega_h and Omega_asymptotic are not meaningful with one factor 1 Omega 2 Call: omegah(m = m, nfactors = nfactors, fm = fm, key = key, flip = flip, 3 digits = digits, title = title, sl = sl, labels = labels, 4 plot = plot, n.obs = n.obs, rotate = rotate, Phi = Phi, option = 5 covar = covar) option, 6 Alpha: 0.78 7 G.6: 0.81 8 Omega Hierarchical: 0.79 9 Omega H asymptotic: 1 Omega Total 0.78 10 11 12 Schmid Leiman Factor loadings greater than 13 14 50 g Q1_1 0.31 F1* h2 u2 p2 0.09 0.91 1 0.2
6.10 𝜔 係数 15 Q1_2 0.57 0.33 0.67 1 16 Q1_3 0.59 0.35 0.65 1 17 Q1_4 0.54 0.29 0.71 1 18 Q1_5 0.54 0.29 0.71 1 19 Q1_6 0.46 0.21 0.79 1 20 Q1_7 0.52 0.27 0.73 1 21 Q1_8 0.49 0.24 0.76 1 22 Q1_9 0.58 0.34 0.66 1 23 Q1_10 0.54 0.30 0.70 1 24 25 With Sums of squares 26 g F1* 27 2.7 0.0 of: 28 29 general/max 30 mean percent general = 31 Explained Common Variance of the general factor = 3.262578e+16 max/min = 1 1 with sd = 0 and cv of 0 1 32 33 The degrees of freedom are 35 34 The number of observations was 35 The root mean square of the residuals is 36 The df corrected root mean square of the residuals is 37 RMSEA index = 38 BIC = prob < and the fit is 2432 1.15 with Chi Square = 2779.53 with 0 0.18 0.14 0.16 and the 10 % confidence intervals are 0.174 0.185 2506.65 39 40 Compare this with the adequacy of just a general factor and no group factors 41 The degrees of freedom for just the general factor are 35 is 42 and the fit 1.15 The number of observations was prob < 2432 with Chi Square = 2779.53 with 0 43 The root mean square of the residuals is 44 The df corrected root mean square of the residuals is 0.14 0.16 45 46 RMSEA index = 47 BIC = 0.18 and the 10 % confidence intervals are 2506.65 48 49 Measures of factor score adequacy 50 g F1* 51 Correlation of scores with factors 0.89 0 52 Multiple R square of scores with factors 0.79 0 53 Minimum correlation of factor score estimates 0.59 -1 54 55 56 51 Total, General and Subset omega for each subset g F1* 0.174 0.185
6.10 𝜔 係数 57 Omega total for total scores and subscales 0.78 0.79 58 Omega general for total scores and subscales 0.79 0.79 59 Omega group for total scores and subscales 0.00 0.00 引数 nfactors は,G因子以外の共通因子の数を指定します。nfactors=1 にした場合,G 因子との識別性を失うため,実質的に 1 因子(G 因子のみ)での因子分析を行っているのと 同じ状態になります。デフォルトでは nfactors=3 となっていますが,この設定で計算した 場合,実質 4 因子(G 因子に加えて 3 つの共通因子)で因子分析を行うことになります。そ の結果,G 因子のみを考慮する 𝜔ℎ では過小推定してしまい,一方 𝜔𝑡 では因子数が多すぎて 過大推定になってしまいがちなようです (岡田, 2011)。 出力を見てみましょう。 6 Alpha: 0.79 7 G.6: 0.82 8 Omega Hierarchical: 0.79 9 Omega H asymptotic: 1 Omega Total 0.79 10 一番上には alpha() 関数と同じように 𝛼 係数および G6 を出してくれています。その下に 続いているのが 𝜔ℎ (omega hierarchical) です。報告する際にはこの値を示せば OK です。そ の下にある Omega H asymptotic は,「項目数が無限にあったら」,つまり「独自因子の影 響を無視できたら」という状況での 𝜔ℎ 係数なので報告の必要がありません。最後には 𝜔𝑡 係 数が表示されています。nfactors=1 のときには 𝜔ℎ と同じ値になります*19 。 12 Schmid Leiman Factor loadings greater than 13 g F1* h2 0.2 u2 p2 14 Q1_1 0.32 0.10 0.90 1 15 Q1_2 0.59 0.34 0.66 1 16 Q1_3 0.60 0.36 0.64 1 17 18 (以下略) 因子負荷行列です。nfactors=1 の場合,G 因子のみでの因子分析の結果になります。 fa(dat[,cols],nfactors=1) と同じ値が出ているはずです。 今回は,最初に作成した cols に含まれる 10 項目で omega() を実行しましたが,この 10 項 目は 2 因子ぶんの項目になっています。続いて 2 因子にして結果を見比べてみましょう。 2 因子だとどうなるの 52 1 omega(dat[, cols], nfactors = 2, poly = TRUE) *19 nfactors=1 のときには G 因子と共通因子を識別できないので,𝜔𝑡 のほうも G 因子のみで推定を行ってい ます。そのために 2 つは同じ値になるのです。
6.10 𝜔 係数 Omega Q1_9 0.4 Q1_7 0.4 Q1_6 0.3 Q1_10 0.4 g 0.4 Q1_8 0.5 Q1_3 0.4 Q1_2 0.4 0.4 Q1_5 0.2 Q1_4 Q1_1 図6.21 0.6 0.6 0.5 0.5 0.5 F1* 0.7 0.6 F2* 0.5 0.4 0.4 omega() のプロット 1 Omega 2 Call: omegah(m = m, nfactors = nfactors, fm = fm, key = key, flip = flip, 3 digits = digits, title = title, sl = sl, labels = labels, 4 plot = plot, n.obs = n.obs, rotate = rotate, Phi = Phi, option = option, 5 covar = covar) 6 Alpha: 0.78 7 G.6: 0.81 8 Omega Hierarchical: 0.42 9 Omega H asymptotic: 0.51 Omega Total 0.83 10 11 12 Schmid Leiman Factor loadings greater than 13 g F1* F2* h2 u2 p2 14 Q1_1 0.22 0.39 0.21 0.79 0.24 15 Q1_2 0.43 0.60 0.55 0.45 0.34 16 Q1_3 0.46 0.66 0.65 0.35 0.33 17 Q1_4 0.38 0.39 0.32 0.68 0.45 18 Q1_5 0.40 19 Q1_6 0.33 0.50 0.36 0.64 0.31 20 Q1_7 0.38 0.55 0.45 0.55 0.32 21 Q1_8 0.35 0.47 0.35 0.65 0.36 22 Q1_9 0.43 0.58 0.53 0.47 0.35 23 Q1_10 0.39 0.49 0.39 0.61 0.38 0.52 0.43 0.57 0.37 24 25 With Sums of squares 26 g F1* F2* 27 1.5 1.4 1.4 of: 28 29 53 general/max 1.06 max/min = 1 0.2
6.10 𝜔 係数 30 mean percent general = 31 Explained Common Variance of the general factor = 0.34 with sd = 0.05 and cv of 0.16 0.35 32 33 The degrees of freedom are 26 34 The number of observations was 35 The root mean square of the residuals is 36 The df corrected root mean square of the residuals is 37 RMSEA index = 38 BIC = prob < and the fit is 2432 0.22 with Chi Square = 540.84 with 1.2e-97 0.09 0.04 0.05 and the 10 % confidence intervals are 0.084 0.097 338.13 39 40 Compare this with the adequacy of just a general factor and no group factors 41 The degrees of freedom for just the general factor are 35 is 42 and the fit 1.35 The number of observations was prob < 2432 with Chi Square = 3282.92 with 0 43 The root mean square of the residuals is 44 The df corrected root mean square of the residuals is 0.19 0.21 45 46 RMSEA index = 47 BIC = 0.195 and the 10 % confidence intervals are 0.19 0.201 3010.05 48 49 Measures of factor score adequacy 50 g F1* F2* 51 Correlation of scores with factors 0.66 0.75 0.76 52 Multiple R square of scores with factors 0.43 0.56 0.58 53 Minimum correlation of factor score estimates -0.14 0.12 0.16 54 55 Total, General and Subset omega for each subset 56 g F1* F2* 57 Omega total for total scores and subscales 0.83 0.78 0.78 58 Omega general for total scores and subscales 0.42 0.27 0.28 59 Omega group for total scores and subscales 0.39 0.51 0.51 先ほどとは打って変わって,Omega Hierarchical の値がかなり低下しました。データがか なりきれいに 2 因子に分かれているため,「すべてに共通の 1 因子 (G 因子)」の影響がかな り弱いことがわかります。このような状態では,10 項目の合計点の内的一貫性は低い,とい う判断になるわけです。このように,𝜔ℎ では内的一貫性だけでなく「一因子性」も評価して あげることが出来ている感じがします。 omega() 関数では,デフォルトで図6.21のようなプロットを出してくれます。F1 および F2 に関しては,因子分析を行った結果をもとに,勝手に最適なグルーピングをしてくれている わけです。ボトムアップに尺度を作成するような場面では,fa() の結果に基づくこのやり 方で 𝜔ℎ を求めるのが良いでしょう。一方で,もともと因子構造の仮定があるような場合は, omega() 関数は勝手にグループを決定してしまうので使えません。そのような場合には,次 54
6.11 主成分分析との違い 回とりあげる SEM の枠組み(lavaan パッケージ)で 𝜔ℎ を計算しましょう。 6.11 主成分分析との違い 因子分析は,よく主成分分析と対比される手法です。ここで少しだけその違いを説明してお きましょう。 主成分分析の主目的は,次元圧縮です。これは,複数の変数を合体させて最も高い説明力を 持つ「最強の変数」を作るようなイメージです。グラフィカルモデル的には図6.22のように なります。因子分析モデル(図6.4)と比べると,矢印の向きが逆になっているわけです。 図6.22 主成分分析モデル そして,主成分分析ではとにかく「なるべく分散説明率の高い主成分を作る」ことで「なる べく少ない数の変数にまとめる」ことが重要となります。そのため,全ての項目間に正の相 関がある場合には,基本的に第一主成分に対しては全ての項目が高い負荷量を示すと考えら れます*20 。これと比べると,因子分析では項目をいくつかのグループに分ける意識が強いと 言えます。そのため,因子数は最小限というよりも「メリハリのある分かれ方(=単純構造) をする数」になるように設定されることが多いです。 また,数理的には主成分分析では誤差を想定しないという違いがあります。なので図6.22に おいても,これに誤差項(独自因子)を追加することはありません。もう少し具体的に説明 しましょう。因子分析では,観測変数の相関行列 𝚺 を 𝐁⊤ 𝐁 + 𝚿 と分解した上で,𝐁⊤ 𝐁 の 部分に対して固有値分解を行うことで因子負荷量を計算していました。これに対して主成分 分析では 𝚺 そのものを固有値分解したものが主成分負荷量になります。主成分分析のほうが 𝚿 のぶんだけ大きな値の行列を固有値分解するため,得られる固有値も大きな値になります。 R の psych パッケージでは,pca() という関数で行うことが出来ます。固有値分解の結果と 見比べてみましょう。 固有値分解 55 1 eig <- eigen(cor(dat[, cols])) 2 # 主成分負荷量は固有ベクトルに固有値の平方根をかけた値になる 3 eig$vectors[, 1] * sqrt(eig$values[1]) *20 ただ,主成分の解釈をしようと考えるケースもあるらしく,その場合には主成分軸の回転をすることもあり ます (林他, 2008)。ただしこの場合得られる得点は「なるべく分散説明率の高い」を追求した得点ではない ため,厳密には主成分分析ではなくなる,という見方もあるようです。
6.11 主成分分析との違い 1 [1] -0.3240019 -0.5916129 -0.6105479 -0.5661146 -0.5645371 -0.5047118 2 [7] -0.5737352 -0.5470221 -0.6143070 -0.5920611 主成分分析 1 # デフォルトではバリマックス回転してしまうので,引数rotateを指定する 2 pca(dat[, cols], nfactors = 2, rotate = "none") 1 Principal Components Analysis 2 Call: principal(r = r, nfactors = nfactors, residuals = residuals, 3 rotate = rotate, n.obs = n.obs, covar = covar, scores = scores, 4 missing = missing, impute = impute, oblique.scores = oblique.scores, 5 6 method = method, use = use, cor = cor, correct = 0.5, weight = NULL) Standardized loadings (pattern matrix) based upon correlation matrix 7 PC1 PC2 h2 u2 com 8 Q1_1 0.32 0.46 0.32 0.68 1.8 9 Q1_2 0.59 0.49 0.59 0.41 1.9 10 Q1_3 0.61 0.52 0.64 0.36 1.9 11 Q1_4 0.57 0.28 0.40 0.60 1.5 12 Q1_5 0.56 0.44 0.51 0.49 1.9 13 Q1_6 0.50 -0.45 0.46 0.54 2.0 14 Q1_7 0.57 -0.44 0.53 0.47 1.9 15 Q1_8 0.55 -0.38 0.45 0.55 1.8 16 Q1_9 0.61 -0.42 0.56 0.44 1.8 17 Q1_10 0.59 -0.35 0.47 0.53 1.6 18 19 PC1 PC2 20 SS loadings 3.08 1.83 21 Proportion Var 0.31 0.18 22 Cumulative Var 0.31 0.49 23 Proportion Explained 0.63 0.37 24 Cumulative Proportion 0.63 1.00 25 26 Mean item complexity = 27 Test of the hypothesis that 2 components are sufficient. 1.8 28 29 30 The root mean square of the residuals (RMSR) is with the empirical chi square 1871.22 0.09 with prob < 0 31 32 Fit based upon off diagonal values = 0.88 PC1 列は符号こそ異なりますが,データの相関行列をそのまま固有値分解した結果と一致す る値が得られています。 56
参考文献 参考文献 Bagozzi, R. P. and Yi, Y. (1988). On the evaluation of structural equation models. Journal of the Academy of Marketing Science, 16(1), 74–94. https://doi.org/10.100 7/BF02723327 Bartlett, M. S. (1950). Tests of significance in factor analysis. British Journal of Statistical Psychology, 3(2), 77–85. https://doi.org/10.1111/j.2044-8317.1950.tb00285.x Bartlett, M. S. (1951). The effect of standardization on a 𝜒2 approximation in factor analysis. Biometrika, 38(3/4), 337–344. https://doi.org/10.2307/2332580 Browne, M. W. (1968). A comparison of factor analytic techniques. Psychometrika, 33(3), 267–334. https://doi.org/10.1007/BF02289327 Flora, D. B. (2020). Your coefficient alpha is probably wrong, but which coefficient omega is right? A tutorial on using R to obtain better reliability estimates. Advances in Methods and Practices in Psychological Science, 3(4), 484–501. https://doi.org/ 10.1177/2515245920951747 Fornell, C. and Larcker, D. F. (1981). Evaluating structural equation models with unobservable variables and measurement error. Journal of Marketing Research, 18(1), 39–50. https://doi.org/10.2307/3151312 Guttman, L. (1953). Image theory for the structure of quantitative variates. Psychometrika, 18(4), 277–296. https://doi.org/10.1007/BF02289264 萩生田伸子・繁桝算男 (1996). 順序付きカテゴリカルデータへの因子分析の適用に関するい くつかの注意点 心理学研究,67 (1), 1–8. https://doi.org/10.4992/jjpsy.67.1 林邦好・冨田誠・田中豊 (2008). 主成分分析における軸の回転について 計算機統計学,19(2), 89–101. https://doi.org/10.20551/jscswabun.19.2_89 堀啓造 (2005). 因子分析における因子数決定法――並行分析を中心にして―― 香川大学経済 論叢,77 (4), 35–70. URL:https://www.ec.kagawa-u.ac.jp/~hori/yomimono/factor2 005.pdf Kaiser, H. F. (1970). A second generation little jiffy. Psychometrika, 35(4), 401–415. https://doi.org/10.1007/BF02291817 57
参考文献 Kaiser, H. F. (1974). An index of factorial simplicity. Psychometrika, 39(1), 31–36. https://doi.org/10.1007/BF02291575 Kaiser, H. F. (1981). A revised measure of sampling adequacy for factor-analytic data matrices. Educational and Psychological Measurement, 41(2), 379–381. https://doi. org/10.1177/001316448104100216 Kaiser, H. F. and Rice, J. (1974). Little jiffy, mark IV. Educational and Psychological Measurement, 34(1), 111–117. https://doi.org/10.1177/001316447403400115 小杉考司 (2018). 言葉と数式で理解する多変量解析入門 北大路書房 岡田謙介 (2011). クロンバックの 𝛼 に代わる信頼性の推定法について――構造方程式 モ デ リ ン グ に よ る方法・McDonald の 𝜔 の 比較―― 日本 テ ス ト 学 会 誌,7, 37–50. https://doi.org/10.24690/jart.7.1_37 Padilla, M. A. and Divers, J. (2016). A comparison of composite reliability estimators: Coefficient omega confidence intervals in the current literature. Educational and Psychological Measurement, 76(3), 436–453. https://doi.org/10.1177/0013164415593776 Revelle, W. and Rocklin, T. (1979). Very simple structure: An alternative procedure for estimating the optimal number of interpretable factors. Multivariate Behavioral Research, 14(4), 403–414. https://doi.org/10.1207/s15327906mbr1404_2 Shirkey, E. C. and Dziuban, C. D. (1976). A note on some sampling characteristics of the measure of sampling adequacy (MSA). Multivariate Behavioral Research, 11(1), 125–128. https://doi.org/10.1207/s15327906mbr1101_9 Thurstone, L. L. (1947). Multiple factor analysis: A development and expansion of vectors of the mind, Chicago. University of Chicago Press. 58