「入門はじめての多変量解析」をPythonで写経 Vol.19 ～ ４章「はじめての因子分析」④因子の回転

４章「はじめての因子分析」

書籍の著者　石村貞夫 先生、石村光資郎 先生

書籍「入門はじめての多変量解析」４章「はじめての因子分析」の Python写経活動記録 です。 

多変量解析の入門を Python と一緒に学ぶ写経シリーズです。
因子分析のうち「探索的因子分析」を取り扱います。

この記事は、因子分析における 因子の回転の一種「直交回転」 を実践します。
軸を回転することで因子の解釈をしやすくなるようです。

ChatGPT 活用型学習で進めてまいります！
では書籍を開いて多変量解析の旅に出かけましょう🚀

はじめに

このブログシリーズは、書籍「入門はじめての多変量解析」（東京図書、「テキスト」と呼びます）の Python 写経を通じて得た「多変量解析の楽しさ」をご紹介します。

書籍の紹介と引用表記はリンク先の記事に掲載しています。

４章 はじめての因子分析

この記事は４章の以下のSectionを取り扱います。

記事に用いるデータは、テキストに掲載されたデータそのものを引用しています。
データ件数の少ないものはコード上でデータを登録し、データ件数の多いものはCSVファイル化してデータを読み込みしています。

この記事で用いるライブラリをインポートします。

### インポート

# 数値計算
import numpy as np
import pandas as pd
import sympy

# 因子分析
import statsmodels.api as sm

# LaTeX表示
from IPython.display import Math

# 描画
import matplotlib.pyplot as plt
import matplotlib.patches as patches
plt.rcParams['font.family'] = 'Meiryo'  # または import japanize_matplotlib

因子分析と回転

突然、「回転」と言われても… という感じですよね。
ざっくり回転は「共通因子の解釈をしやすくする目的で、変数と共通因子の関係性にメリハリをつける手法の一つ」だと思っています。
ChatGPTにやさしくざっくり解説してもらいましょう。

🔄 回転ってなに？
因子分析で見つけた「因子」は、まるでグラフ上に立った棒みたいなもの。
そこから

🔄 なんで回転するの？

1. ぴったりハマるパズルにするため

   • 回転前：あちこちのピース（変数）が、いくつもの穴（因子）に少しずつ入っちゃってる

   • 回転後：「このピースはここ！」と１つの穴にだけスッと入るようになる
     → どの変数がどの因子にグッと効いているかが、一目でわかりやすくなります。

2. 説明力はそのまま、見た目だけスッキリ

   • 回転前後で「どれだけデータを説明できるか」は変わりません。

   • ちょうど、写真の向きを変えても景色そのものは同じ。

   • 見せたい部分をくっきり見える角度に変えるだけです。

3. ２つの回し方

   • 直交回転（例: Varimax）：因子同士はずっと直角同士。まずはコレ！

   • 斜交回転（例: Promax）：因子同士がちょっとだけ仲良く（相関）なるように回す。
     「因子同士が完全に独立じゃなくてもいいかな？」というときに使います。

🔄 まとめ

• 回転 は「因子の向きをちょっとだけ変えて、解釈しやすくする技」

• データの本質（説明力）はそのままに、

• 変数と因子の“くっつき具合”をパッと見でわかりやすく整えます😊

最初は Varimax 回転で遊んでみて、
「どの質問がどの因子にパチッとはまるかな？」を楽しんでみてくださいね！

この記事は「直交回転」を扱います。
直角に交差する軸を維持しています。

回転を体感

テキスト p.154 表 4.6.1～表 4.6.3 の因子負荷行列と因子変換行列を用いて、因子の回転をやってみます。

### 因子の回転 p.154

# 回転前の因子負荷行列Λf
Λf = np.array([[-0.354, 0.724],
               [0.616, 0.587],
               [0.910, -0.116]])
Λf_df = pd.DataFrame(Λf, columns=['因子1', '因子2'], index=['x1', 'x2', 'x3'])
print('【回転前の因子負荷行列】')
display(Λf_df)

# 因子変換行列 T の設定
T = np.array([[0.854, -0.520],
              [0.520, 0.854]])

# 回転後の因子負荷行列 Λg = Λf @ T
Λg = Λf @ T
Λg_df = pd.DataFrame(Λg, columns=['因子1', '因子2'], index=['x1', 'x2', 'x3'])

# ハイライト付きで回転後の因子負荷行列を表示
print('【回転後の因子負荷行列】')
color = 'background-color: lightpink'
display(Λg_df
 .style.map(lambda x: color if abs(x) > 0.7  else '')
 .format(precision=3))

【実行結果】
回転前・後の因子負荷行列を並べています。
回転後の因子負荷行列には目印のハイライトを入れています。

【回転前後の変化点と効果】
変数 $${x_1, x_2}$$ に対応する因子負荷量は、１つの共通因子の値が大きくなり、他の共通因子の値が小さくなっています。
⇒【効果】変数が影響を受けている共通因子を特定しやすくなった！
⇒ 第１因子は変数 $${x_2, x_3}$$ を用いて解釈し、第２因子は変数 $${x_1}$$ を用いて解釈すればよさそう！

続いて「これがなぜ回転なの？」に応えましょう。
回転前・後の因子負荷行列を可視化します。

### 座標を固定したときの回転前・後の因子負荷の可視化

# 設定と準備
loagings = [Λf, Λg]
titles = ['回転前', '回転後']
var_names = ['$x_1$', '$x_2$', '$x_3$']

# 描画領域の設定
fig, axes = plt.subplots(1, 2, figsize=(7, 4), tight_layout=True)
# 回転前・回転後の散布図描画を繰り返し処理
for loading, title, ax in zip(loagings, titles, axes.flat):
    # 因子負荷の散布図の描画
    ax.plot(loading[:, 0], loading[:, 1], 'o')
    # 3つの変数の変数名描画を繰り返し処理
    for x, y, s in zip(loading[:, 0], loading[:, 1], var_names):
        # 変数名の描画
        ax.text(x=x, y=y-0.15, s=s, ha='center')
    # x=0,y=0の直線の描画
    ax.axhline(0, color='black', ls='--', lw=0.5)
    ax.axvline(0, color='black', ls='--', lw=0.5)
    # 修飾
    ax.set(xlim=(-1.1, 1.1), ylim=(-1.1, 1.1), yticks=[-1.0, -0.5, 0, 0.5, 1.0],
           xlabel='第１因子 $f_1$', ylabel='第２因子 $f_2$', title=title,
           aspect='equal')

【実行結果】

【補足説明】
夜空の星座が時間とともに移動するように、因子負荷量の各点は原点を中心にして時計回りに移動＝回転しています！
因子負荷量が時計回りに回転したのは、グラフの軸が反時計回りに回転したからです。
「因子の回転は軸の回転」なのです。

📕 ここまでのまとめ
・回転は「共通因子の解釈をしやすくする」手段です。
・実際に軸が回転しています。

回転の計算

因子負荷行列 $${\Lambda_f}$$ に「因子変換行列」と呼ばれる回転行列 $${T}$$ を掛けて（行列積）、回転後の因子負荷行列 $${\Lambda_g}$$ を得ます。

冒頭の体感で用いた以下の因子変換行列 $${\Lambda_f}$$ と因子変換行列 $${T}$$ で因子負荷行列を回転させます。
Python で計算しましょう。

# 因子負荷行列と因子変換行列の行列積
Λf @ T

【実行結果】
冒頭で用いた回転後の因子負荷行列が求まりました。

因子変換行列 $${T}$$ による回転は、$${T}$$ が三角関数を用いた「回転行列」になっていることで実現できています。

$${\cos \theta = 0.854}$$ を用いて、角度 $${\theta}$$ を求めてみましょう。
テキスト p.157 上部の計算に相当します。

### cosθ=0.854 を解く p.157
cosθ = 0.854
print(f'θ = {np.degrees(np.arccos(cosθ)):.2f}°')

【実行結果】
$${\theta = 31.15^{\circ}}$$ です。
軸が反時計回りに $${31.15^{\circ}}$$ 回転していることを示しています。

回転をもう一度図示します。
「軸の回転」を強調する可視化です。
テキスト p.157 図 4.6.2 に相当します。

### 因子の回転の可視化 p.157 図4.6.2

## 描画領域の設定
fig, ax = plt.subplots(figsize=(7, 7))

## y=0の水平線の描画
# y=0の水平線（黒）の描画
ax.plot([-1.1, 1.1], [0, 0], color='black', lw=1)
# y=0の水平線上の軸目盛り点の描画
ax.plot([-1.0, -0.5, 0.5, 1.0], [0, 0, 0, 0], 'o', color='black', ms=5)
# y=0の水平線上の軸目盛りラベルの描画
for x in [-1.0, -0.5, 0.5, 1.0]:
    ax.text(x=x, y=-0.08, s=x, ha='center')

## x=0の垂直線の描画
# x=0の垂直線（黒）の描画
ax.plot([0, 0], [-1.1, 1.1], color='black', lw=1)
# x=0の垂直線上の軸目盛り点の描画
ax.plot([0, 0, 0, 0], [-1.0, -0.5, 0.5, 1.0], 'o', color='black', ms=5)
# x=0の垂直線上の軸目盛りラベルの描画
for y in [-1.0, -0.5, 0.5, 1.0]:
    ax.text(x=0.03, y=y, s=y, ha='left', va='center')

## 回転後のy軸の描画
# y=0の水平線がθ回転した斜線（赤）の描画
x_h, y_h = np.array([1.3, 0]) @ T
ax.plot([-x_h, x_h], [y_h, -y_h], color='tab:red')
# y=0の水平線がθ回転した斜線上の軸目盛り点(-0.5, 0.5)の描画
x1, y1 = np.array([0.5, 0]) @ T * np.array([1, -1])
ax.plot(x1, y1, 'o', color='tab:red', ms=5)
ax.plot(-x1, -y1, 'o', color='tab:red', ms=5)
# y=0の水平線がθ回転した斜線上の軸目盛りラベル(-0.5, 0.5)の描画
ax.text(x=x1-0.07, y=y1+0.03, s=0.5, color='tab:red')
ax.text(x=-x1-0.1, y=-y1+0.03, s=-0.5, color='tab:red')
# y=0の水平線がθ回転した斜線上の軸目盛り点(-1.0, 1.0)の描画
x2, y2 = np.array([1, 0]) @ T * np.array([1, -1])
ax.plot(x2, y2, 'o', color='tab:red', ms=5)
ax.plot(-x2, -y2, 'o', color='tab:red', ms=5)
# y=0の水平線がθ回転した斜線上の軸目盛りラベル(-1.0, 1.0)の描画
ax.text(x=x2-0.07, y=y2+0.03, s=1.0, color='tab:red')
ax.text(x=-x2-0.1, y=-y2+0.03, s=-1.0, color='tab:red')

## 回転後のx軸の描画
# x=0の垂直線がθ回転した斜線（赤）の描画
x_v, y_v = np.array([0, 1.3]) @ T
ax.plot([-x_v, x_v], [y_v, -y_v], color='tab:red')
# x=0の垂直線がθ回転した斜線上の軸目盛り点(-0.5, 0.5)の描画
x3, y3 = np.array([0, 0.5]) @ T * np.array([1, -1])
ax.plot(x3, y3, 'o', color='tab:red', ms=5)
ax.plot(-x3, -y3, 'o', color='tab:red', ms=5)
# x=0の垂直線がθ回転した斜線上の軸目盛りラベル(-0.5, 0.5)の描画
ax.text(x=x3-0.13, y=y3-0.05, s=-0.5, color='tab:red')
ax.text(x=-x3-0.1, y=-y3-0.05, s=0.5, color='tab:red')
# x=0の垂直線がθ回転した斜線上の軸目盛り点(-1.0, 1.0)の描画
x4, y4 = np.array([0, 1]) @ T * np.array([1, -1])
ax.plot(x4, y4, 'o', color='tab:red', ms=5)
ax.plot(-x4, -y4, 'o', color='tab:red', ms=5)
# x=0の垂直線がθ回転した斜線上の軸目盛りラベル(-1.0, 1.0)の描画
ax.text(x=x4-0.13, y=y4-0.05, s=-1.0, color='tab:red')
ax.text(x=-x4-0.1, y=-y4-0.05, s=1.0, color='tab:red')

## 因子負荷の散布図（青）の描画
ax.scatter(Λf[:, 0], Λf[:, 1], s=80)

## 因子負荷の座標(a_i1, a_i2)の描画
coords = ['($a_{11}$, $a_{12}$)', '($a_{21}$, $a_{22}$)', '($a_{31}$, $a_{32}$)']
adjs = [0.00, 0.16, -0.05]
for i in range(3):
    # 因子負荷の座標の取得
    x, y = Λf[i]
    # 座標(a1*, a2*)（青い字）の表示
    ax.text(x=x+0.03, y=y+adjs[i], s=coords[i], color='tab:blue')
    # 回転前の座標（黒い字）の表示
    ax.text(x=x+0.03, y=y+adjs[i]-0.07, s=f'({x:.3f}, {y:.3f})', color='black')
    # 回転後の座標（赤い字）の表示
    xr, yr = Λg[i]
    ax.text(x=x+0.03, y=y+adjs[i]-0.14, s=f'({xr:.3f}, {yr:.3f})',
            color='tab:red')

## 角度 θ の描画
# 0度からθの弧（緑）の描画
arc1 = patches.Arc(
    xy=(0, 0), width=0.4, height=0.4, 
    theta1=0, theta2=31.35, edgecolor="green", linewidth=1)
ax.add_patch(arc1)
# 0度からθの弧のテキストの描画
ax.text(x=0.2, y=0.045, s='$\\theta=31.35$°', color='tab:green')
# 90度からθの弧（緑）の描画
arc2 = patches.Arc(
    xy=(0, 0), width=0.4, height=0.4, 
    theta1=90, theta2=90+31.35, edgecolor="green", linewidth=1)
ax.add_patch(arc2)
# 90度からθの弧のテキストの描画
ax.text(x=-0.08, y=0.2, s='$\\theta$', color='tab:green')
# 0のテキストの描画
ax.text(x=-0.06, y=-0.1, s='$0$', color='black', fontsize=14)

## 修飾
ax.set(xlim=(-1.1, 1.1), ylim=(-1.1, 1.1), xticks=[], yticks=[])
ax.set_title('因子の回転の可視化')
ax.set_xlabel('第１因子 $f_1$', fontsize=12)
ax.set_ylabel('第２因子 $f_2$', fontsize=12);

【実行結果】
回転前の黒い軸、回転後の赤い軸に注目しましょう。
反時計回りに $${31.15^{\circ}}$$ 回転しています！
各因子負荷量のデータ点には「回転前（黒い字）・回転後（赤い字）の座標＝因子負荷量」を付記しています。

■ ここまでの分析を statsmodels でやってみる

テキスト Section 4.6 にはデータそのものが掲載されていません。
見えないデータに対する SPSS の因子分析結果のみ、具体的には回転前・後の「因子負荷行列」と「因子変換行列」が開示されています。

ところで statsmodels の因子分析はデータを与えるほかに「データの相関行列」を与えることも可能です。
相関行列さえ推定できれば、statsmodels で因子分析できるのです！

データの分散共分散行列（相関行列）$${\Sigma}$$ と因子負荷行列 $${\Lambda_f}$$ の関係を思い出しましょう。

💡 因子負荷行列からデータの相関行列を推定できるかも！？
回転前の因子負荷行列の行列積を算出して、対角成分に１を埋めれば、データの相関行列に近づけそうです！

ということで、statsmodels で相関行列をもとに因子分析を行ってみます。

### statsmodelsによる因子分析の実行 p.154 ※SPSSと比べて第２因子が正負逆転している

# データの相関行列の推定 ※因子負荷行列の行列積→対角成分に1.0を設定
corr_matrix = Λf @ Λf.T             # 因子負荷行列の行列積
np.fill_diagonal(corr_matrix, 1.0)  # 対角成分に1.0を設定

# 因子分析の実行 ※相関行列を元にして因子分析を行う
result_sm = sm.multivariate.Factor(
    corr=corr_matrix,
    n_factor=3, 
    method='pa',
    smc=True,
    endog_names=['x1', 'x2', 'x3']
).fit()

# 因子の回転：バリマックス
result_sm.rotate(method='varimax')

# 結果表示
result_sm.summary()

【実行結果】
回転前の因子負荷行列が Pre-rotated loadings、回転後の因子負荷行列が varimax rotated loadings です。
テキストと比べて第２因子（factor 1）が正負逆転していますが、気にしない気にしない。

回転後の因子負荷行列を見てみましょう。
因子分析の結果 result_sm に対して get_loadings_frame() を実行すると…

# 因子負荷行列（ハイライト付き）
result_sm.get_loadings_frame(threshold=0)

【実行結果】
共通因子と変数の関係が自動でハイライトされます！見やすい！
第１因子は変数 $${x_2, x_3}$$、第２因子は変数 $${x_1}$$ から解釈します。

直交回転に用いた因子変換行列を表示しましょう。
因子分析の結果 result_sm の rotation_matrix 属性で取り出します。

# 因子変換行列の表示
result_sm.rotation_matrix

【実行結果】
テキストの SPSS の因子変換行列と近い値です（ちょっとずれてる）。

回転前・後の因子負荷行列を可視化しましょう。

### 座標を固定したときの回転前・後の因子負荷の可視化

# 設定と準備
loagings = [result_sm.loadings_no_rot, result_sm.loadings]
titles = ['回転前', '回転後']
var_names = result_sm.endog_names

# 描画領域の設定
fig, axes = plt.subplots(1, 2, figsize=(7, 4), tight_layout=True)
# 回転前・回転後の散布図描画を繰り返し処理
for loading, title, ax in zip(loagings, titles, axes.flat):
    # 因子負荷の散布図の描画
    ax.plot(loading[:, 0], loading[:, 1], 'o')
    # 3つの変数の変数名描画を繰り返し処理
    for x, y, s in zip(loading[:, 0], loading[:, 1], var_names):
        # 変数名の描画
        ax.text(x=x, y=y-0.15, s=s, ha='center')
    # x=0,y=0の直線の描画
    ax.axhline(0, color='black', ls='--', lw=0.5)
    ax.axvline(0, color='black', ls='--', lw=0.5)
    # 修飾
    ax.set(xlim=(-1.1, 1.1), ylim=(-1.1, 1.1), yticks=[-1.0, -0.5, 0, 0.5, 1.0],
           xlabel='第１因子 $f_1$', ylabel='第２因子 $f_2$', title=title,
           aspect='equal')

【実行結果】
先程のプロットと比べると第２因子の軸で正負反転していますが、しっかりと反時計回りに回転している様子がつかめると思います！

データの相関行列を推定すること、相関行列から因子分析ができることを学びました。

因子の回転に関する不定性

テキスト p.158 の疑問１に、因子負荷行列を回転したら別の因子に（別の意味合いに）変化してしまうのでは？があります。

答えは「因子の回転に関する不定性」の性質により、因子負荷行列を回転させても、共通性・独自性は変わらないから、だそうです。

■ 回転で共通性が変わらないことの確認（具体的な数値版）
共通性は因子負荷行列 $${\Lambda_f}$$ の行列積 $${\Lambda_f \bm \cdot \Lambda_f^{\top}}$$ であり、独自性 $${D}$$＝１－共通性、ですので、共通性が回転前・後で変わらないことを Python で確認してみましょう。
テキスト p.160 の計算に相当します。

最初に回転前の因子負荷行列の行列積を計算します。

### 回転前の因子負荷の行列積 p.160
(Λf @ Λf.T).round(4)

【実行結果】

続いて回転後の因子負荷行列の行列積を計算します。

### 回転後の因子負荷の行列積 p.160 ※回転前とほぼ同じ
(Λg @ Λg.T).round(4)

【実行結果】

共通性は回転前・後でほぼ同じ値になることが分かりました。

■ 回転で共通性が変わらないことの確認（文字式版）
テキストは「このことを、文字式で確認してみましょう」と続きます。
sympy ライブラリを使って文字式で確認してみましょう！
テキスト p.161 の計算に相当します。

以下の式で回転前の因子負荷行列 $${\Lambda_{f2}}$$ と 因子変換行列 $${T_2}$$ を設定して、回転後の因子負荷行列 $${\Lambda_{g2}}$$ を計算します。

コードの出力は回転後の因子負荷行列の行列積 $${\Lambda_{g2} \bm \cdot \Lambda_{g2}^{\top}}$$ です。

### 因子の回転に関する不定性を文字式で確認 p.161

# sympyの変数の定義
a11, a12, a21, a22, a31, a32, θ = sympy.symbols('a11 a12 a21 a22 a31 a32 θ')

# 因子負荷行列の定義
Λf2 = sympy.Matrix([[a11, a12],
                    [a21, a22],
                    [a31, a32]])

# 因子変換行列の定義
T2 = sympy.Matrix([[sympy.cos(θ), -sympy.sin(θ)],
                   [sympy.sin(θ), sympy.cos(θ)]])

# 因子の回転
Λg2 = Λf2 @ T2

# 回転後の因子負荷行列の行列積の算出
Λg2_Λg2_T2_expand = (Λg2 * Λg2.T).expand().applyfunc(sympy.simplify)

# 結果の表示
display(Math(f'{sympy.latex(Λg2)} \\bm \cdot {sympy.latex(Λg2.T)}'))
display(Math(f'= {sympy.latex(Λg2_Λg2_T2_expand)}'))

【実行結果】

行列積の結果だけにフォーカスすると：

回転前の因子負荷行列の行列積  $${\Lambda_{f2} \bm \cdot \Lambda_{f2}^{\top}}$$ を計算しましょう。

# 回転前の因子負荷行列の行列積の算出
display(Math(f'{sympy.latex(Λf2)} \\bm \cdot {sympy.latex(Λf2.T)}'))
display(Math(f'= {sympy.latex(Λf2 @ Λf2.T)}'))

【実行結果】

回転前の因子負荷行列の行列積  $${\Lambda_{f2} \bm \cdot \Lambda_{f2}^{\top}}$$ と回転後の因子負荷行列の行列積 $${\Lambda_{g2} \bm \cdot \Lambda_{g2}^{\top}}$$ は一致しています。

回転を行っても共通性は変化しないのですね！

解釈しやすい因子負荷

回転すると共通因子の解釈がしやすくなることは、ここまでで眺めてきたとことです。
テキストは p.158 疑問２で「なぜ、因子を回転するのですか？」とあり、回答は「共通因子の解釈をしやすくするため」となります。

テキスト p.162 表 4.6.4 の因子負荷行列等をお借りして、変数の多いケースを見ていきましょう。
こちらは回転前の因子負荷行列です。

### 因子負荷の解釈 p.162 表4.6.4
Λf3 = pd.DataFrame(
    {'因子1': [-0.719, 0.473, 0.369, 0.917, 0.754, 0.674, 0.697],
     '因子2': [0.293, -0.638, -0.464, 0.192, 0.430, 0.044, 0.219]},
     index=['ストレス', '運動量', '健康', '仕事', '地域活動', '趣味', '家庭生活'])
Λf3

【実行結果】

テキストの SPSS の因子変換行列を用いて、回転を実行します。

# 因子変換行列の設定
T3 = np.array([[0.871, 0.492],
               [0.492, -0.871]])

# 回転後の因子負荷行列の算出
Λg3 = pd.DataFrame(
    Λf3.values @ T3, index=Λf3.index, columns=['因子1', '因子2']).round(3)

# 回転後の因子負荷行列の表示（ハイライト付き）
color = 'background-color: lightpink'
(Λg3
 .style.map(lambda x: color if (abs(x)>0.5)  else '')
 .format(precision=3))

【実行結果】
回転前よりも「変数と共通因子の関係が際立つ」ように変化しています。
関連の強い共通因子の因子負荷が大きくなり、関連の薄い共通因子の因子負荷が小さくなっている様子＝「解釈しやすくなったこと」を確認できます。

【共通因子の解釈】
テキストは回転後の因子負荷行列を参照して「第１因子は外面的充実度」と解釈しています。
そうすると第２因子は「内面的な心身状態の良さ」でしょうか？

回転前・後の因子負荷行列を可視化しましょう。

### 図 4.6.3, 4.6.4 回転前・後の因子プロット ※裏返しになっていて、回転ではない

# 設定と準備
loadings = [Λf3, Λg3]
titles = ['回転前', '回転後']

# 描画領域の設定
fig, axes = plt.subplots(1, 2, figsize=(8, 5), tight_layout=True)
# 回転前・回転後の散布図描画を繰り返し処理
for loading, title, ax in zip(loadings, titles, axes.flat):
    # 因子負荷行列の散布図の描画
    ax.plot(loading.iloc[:, 0], loading.iloc[:, 1], 'o')
    # 3つの変数の変数名描画を繰り返し処理
    for x, y, s in zip(loading.iloc[:, 0], loading.iloc[:, 1], loading.index):
        # 変数名の描画
        ax.text(x=x, y=y-0.1, s=s, ha='right')
    # x=0,y=0の直線の描画
    ax.axhline(0, color='black', ls='--', lw=0.5)
    ax.axvline(0, color='black', ls='--', lw=0.5)
    # 修飾
    ax.set(xlim=(-1.1, 1.1), ylim=(-1.1, 1.1), yticks=[-1.0, -0.5, 0, 0.5, 1.0],
           xlabel='第１因子 $f_1$', ylabel='第２因子 $f_2$', title=title,
           aspect='equal')

【実行結果】
第１因子の大きい変数グループ、第２因子の大きい因子グループは、関連の薄い軸に近づいていることが分かります。
「関連の薄い軸に寄っていく」のが回転の特徴です。

ちなみに、回転前から「回転」しても回転後にはなりません。
「反転」されているからです。
テキストは p161で「因子の回転よりも変換といった方がより正確な表現です。実は直交変換です」と説明しています（うさぎさん談）。

■ ここまでの分析を statsmodels でやってみる
因子負荷行列からデータの相関行列を推定して、statsmodels で因子分析を実行しましょう。

データの相関行列の推定 ⇒ 因子分析の実行 ⇒ 因子の回転 ⇒ 結果の表示の順で処理しています。
statsmodels は因子分析の実行と因子の回転が別々の処理になる感じです。

### statsmodelsによる因子分析の実行 p.162 ※SPSSと比べて第２因子が正負逆転している

# データの相関行列の推定 ※因子負荷行列の行列積→対角成分に1.0を設定
corr_matrix3 = Λf3.values @ Λf3.values.T
np.fill_diagonal(corr_matrix3, 1.0)

# 因子分析の実行 ※相関行列を元にして因子分析を行う
result_sm3 = sm.multivariate.Factor(
    corr=corr_matrix3,
    n_factor=2, 
    method='pa',
    smc=True,
    endog_names=Λf3.index
).fit()

# 因子の回転：バリマックス
result_sm3.rotate(method='varimax')

# 結果表示
result_sm3.summary()

【実行結果】

回転後の因子負荷行列をハイライト付きで表示しましょう。

# 因子負荷行列（ハイライト付き）
result_sm3.get_loadings_frame(threshold=0)

【実行結果】
テキストのグループと同じになっています。

因子変換行列を表示します。

# 因子変換行列の表示
result_sm3.rotation_matrix

【実行結果】
テキストと符号が異なっているようです。

回転前・後の因子負荷行列を可視化しましょう。

### 座標を固定したときの回転前・後の因子負荷の可視化 ※こちらは回転になっている

# 設定と準備
loagings = [result_sm3.loadings_no_rot, result_sm3.loadings]
titles = ['回転前', '回転後']
var_names = result_sm3.endog_names

# 描画領域の設定
fig, axes = plt.subplots(1, 2, figsize=(8, 5), tight_layout=True)
# 回転前・回転後の散布図描画を繰り返し処理
for loading, title, ax in zip(loagings, titles, axes.flat):
    # 因子負荷の散布図の描画
    ax.plot(loading[:, 0], loading[:, 1], 'o')
    # 3つの変数の変数名描画を繰り返し処理
    for x, y, s in zip(loading[:, 0], loading[:, 1], var_names):
        # 変数名の描画
        ax.text(x=x, y=y-0.15, s=s, ha='center')
    # x=0,y=0の直線の描画
    ax.axhline(0, color='black', ls='--', lw=0.5)
    ax.axvline(0, color='black', ls='--', lw=0.5)
    # 修飾
    ax.set(xlim=(-1.1, 1.1), ylim=(-1.1, 1.1), yticks=[-1.0, -0.5, 0, 0.5, 1.0],
           xlabel='第１因子 $f_1$', ylabel='第２因子 $f_2$', title=title,
           aspect='equal')

【実行結果】
じつはこちらは「回転」が成立しています（反転をしていないです）。
テキストと比べて、回転前の第２因子が正負反転しているからでしょう。

今回もテキストのデータを活用して Python による因子分析を実践できました。
楽しかったですね！

あれ？
回転行列（因子変換行列）ってどうやって算出するのだろう…

記事の最後はChatGPTが締めくくります。
今回は清く澄みわたる詩的な雰囲気で。

📘 ChatGPTのひとこと：

今回の写経は以上です。

シリーズの記事

次の記事

前の記事

目次

ブログの紹介

note で７つのシリーズ記事を書いています。
ぜひ覗いていってくださいね！

１．のんびり統計

統計検定２級の問題集を手がかりにして、確率・統計をざっくり掘り下げるブログです。
雑談感覚で大丈夫です。ぜひ覗いていってくださいね。
統計検定２級公式問題集CBT対応版に対応しています。
Python、EXCELのサンプルコードの配布もあります。

２．実験！たのしいベイズモデリング１＆２をPyMC Ver.5で

書籍「たのしいベイズモデリング」・「たのしいベイズモデリング２」の心理学研究に用いられたベイズモデルを PyMC Ver.5で描いて分析します。
この書籍をはじめ、多くのベイズモデルは Ｒ言語＋Stanで書かれています。
PyMCの可能性を探り出し、手軽にベイズモデリングを実践できるように努めます。
身近なテーマ、イメージしやすいテーマですので、ぜひぜひPyMCで動かして、一緒に楽しみましょう！

３．実験！岩波データサイエンス1のベイズモデリングをPyMC Ver.5で

書籍「実験！岩波データサイエンスvol.1」の４人のベイジアンによるベイズモデルを PyMC Ver.5で描いて分析します。
この書籍はベイズプログラミングのイロハをざっくりと学ぶことができる良書です。
楽しくPyMCモデルを動かして、ベイズと仲良しになれた気がします。
みなさんもぜひぜひPyMCで動かして、一緒に遊んで学びましょう！

４．楽しい写経 ベイズ・Python等

ベイズ、Python、その他の「書籍の写経活動」の成果をブログにします。
主にPythonへの翻訳に取り組んでいます。
写経に取り組むお仲間さんのサンプルコードになれば幸いです🍀

５．RとStanではじめる心理学のための時系列分析入門 を PythonとPyMC Ver.5 で

書籍「RとStanではじめる心理学のための時系列分析入門」の時系列分析をPythonとPyMC Ver.5 で実践します。
この書籍には時系列分析のテーマが盛りだくさん！
時系列分析の懐の深さを実感いたしました。
大好きなPythonで楽しく時系列分析を学びます。

６．データサイエンスっぽいことを綴る

統計、データ分析、AI、機械学習、Pythonのコラムを不定期に綴っています。
統計・データサイエンス書籍にまつわる記事が多いです。
「統計」「Python」「数学とPython」「R」のシリーズが生まれています。

７．Python機械学習プログラミング実践記

書籍「Python機械学習プログラミング PyTorch & scikit-learn編」を学んだときのさまざまな思いを記事にしました。
この書籍は、scikit-learnとPyTorchの教科書です。
よかったらぜひ、お試しくださいませ。

最後までお読みいただきまして、ありがとうございました。