因子分析をExcelで理解する

主成分分析とは似て非なる手法として「因子分析」(Factor Analysis) があります。

主成分分析(PCA)では、説明変数に対して重み行列（固有ベクトル）a を線形結合した「主成分」 y<sub>PC1</sub>を合成しました。ここで、主成分は、説明変数と同じ数だけ定義します。

y<sub>PC1</sub> = a<sub>1,1</sub> x<sub>1</sub> + a<sub>1,2</sub>  x<sub>2</sub> +  a<sub>1,3</sub> x<sub>3</sub> + a<sub>1,4</sub> x<sub>4</sub> + a<sub>1,5</sub> + ...

因子分析では、説明変数（観測変数）x が「因子」(factor) という潜在変数から合成されるという考え方に基づき、その因子得点 f と重み行列（因子負荷） w 、そして独自因子 e を特定します（主成分分析に独自因子という考え方はありません）。

x<sub>1</sub> = w<sub>1,1</sub> f<sub>1</sub> + w<sub>1,2</sub> f<sub>2</sub> +  e<sub>1</sub>
x<sub>2</sub> = w<sub>2,1</sub> f<sub>1</sub> + w<sub>2,2</sub> f<sub>2</sub> +  e<sub>2</sub>
x<sub>3</sub> = w<sub>3,1</sub> f<sub>1</sub> + w<sub>3,2</sub> f<sub>2</sub> +  e<sub>3</sub>
x<sub>4</sub> = w<sub>4,1</sub> f<sub>1</sub> + w<sub>4,2</sub> f<sub>2</sub> +  e<sub>4</sub>
x<sub>5</sub> = w<sub>5,1</sub> f<sub>1</sub> + w<sub>5,2</sub> f<sub>2</sub> +  e<sub>5</sub>
x<sub>6</sub> = w<sub>6,1</sub> f<sub>1</sub> + w<sub>6,2</sub> f<sub>2</sub> +  e<sub>6</sub>

因子得点 f は、各個体（サンプル）が独自に持つ潜在変数です。因子得点と因子負荷の線形和（w<sub>1,1</sub> f<sub>1</sub> + w<sub>1,2</sub> f<sub>2</sub> など）を「共通因子」と呼び、観測変数独自の「独自因子」e と足し合わせることで「観測変数」として観測できるという考え方です。因子の個数は説明変数よりも小さな個数を使うのが普通で、あらかじめ決めておく必要があります。

（ただし、共通因子や因子などの用語は、調べてみた限り異なる人が異なる定義をしているように見えますので、大変紛らわしいです。私の説明のほうが間違っている可能性もありますのでその際はご容赦ください。）

因子分析を行う便利なツールとして、Pythonで利用可能なScikit-learnなどがありますが、ここでは、因子分析がどのような手法か概観するためにExcelで計算してみたいと思います。

# ワインのデータ

データは <a href="http://archive.ics.uci.edu/ml/datasets/Wine">UC Irvine Machine Learning Repository</a> から取得したものを少し改変しました。

https://github.com/maskot1977/ipython_notebook/blob/master/toydata/wine.xlsx
![01.png](https://qiita-image-store.s3.amazonaws.com/0/211162/86c6df6c-30a4-578e-9cf0-560842422873.png)

今回のデータは、20個のワインの銘柄の特徴が、６種類の「観測変数」(x<sub>1</sub>, x<sub>2</sub>, x<sub>3</sub>, x<sub>4</sub>, x<sub>5</sub>, x<sub>6</sub>) で表されています。この6変数を説明変数とします。

# 基準値

各列に対して、平均値を引いたものを標準偏差で割ります。

![02.png](https://qiita-image-store.s3.amazonaws.com/0/211162/d804b87e-4949-d4d3-940b-d4de9bc18286.png)

# 因子負荷量、因子得点、共通因子

因子負荷量 w は観測変数ごとに与えられる係数、因子得点は個体（サンプル）ごとに与えらえる係数で、それらの線形和が「共通因子」になります。因子分析の目的は、この因子負荷量と因子得点を求めることです。ここではまず、因子負荷量と因子得点を乱数で初期化しておき、その積を仮の共通因子とします。現時点ではランダムな数字で算出されたデタラメな数値が入ってます。

![30.png](https://qiita-image-store.s3.amazonaws.com/0/211162/0fe0e2d5-cbfc-3da5-f92c-b2166e6b3212.png)

# 独自因子

実際の観測値と共通因子との差を、共通因子では説明できない「独自因子」とします。現時点ではランダムな数字で算出されたデタラメな数値が入ってます。

![31.png](https://qiita-image-store.s3.amazonaws.com/0/211162/7c242f47-2a5b-3162-dff6-9facd6fbc1e5.png)

# 因子分析が推定する相関行列

まずは因子負荷量の自分自身に対する行列積を計算し、それを「もしも独自因子がまったくなく、共通因子だけで観測事象が説明できると仮定したときの相関行列」と考えます。現時点ではランダムな数字で算出されたデタラメな数値が入ってます。

![32.png](https://qiita-image-store.s3.amazonaws.com/0/211162/47cbecad-aded-5bea-5b2d-de3a31b85275.png)

次に、独自因子による影響として、対角成分に独自因子の分散を、それ以外の成分に０を入れた行列を作成します。現時点ではランダムな数字で算出されたデタラメな数値が入ってます。

![33.png](https://qiita-image-store.s3.amazonaws.com/0/211162/32b3cc3c-d28e-48a3-41c6-d25fbca5c0bc.png)

上記２つの行列の和を、「因子分析が推定する、共通因子と独自因子の和による相関行列」とします。現時点ではランダムな数字で算出されたデタラメな数値が入ってます。

![34.png](https://qiita-image-store.s3.amazonaws.com/0/211162/aa5fedcf-3de4-8171-05f7-614c69e6403a.png)

# 観測変数の相関行列

実際の観測変数の相関行列を作成します。

![35.png](https://qiita-image-store.s3.amazonaws.com/0/211162/45fd5945-53f7-ea15-9277-e745080de5a1.png)

# 2つの相関行列の差分を最小化する

「因子分析が推定する、共通因子と独自因子の和による相関行列」と、「実際の観測変数の相関行列」の差分の２乗和を計算します。これを最小化する因子負荷量と因子得点を計算するのが、因子分析の目的になります。

![36.png](https://qiita-image-store.s3.amazonaws.com/0/211162/74817f9e-b4f0-2402-e744-8f39bb4106bd.png)

目的変数と制約条件を整理します。目的変数は、２つの相関係数の差分の２乗和で、これの最小化が目的です。

![37.png](https://qiita-image-store.s3.amazonaws.com/0/211162/77383366-7278-b7f0-d7f6-a1b2d38e1a95.png)

制約条件その１としては、因子負荷量の積和をゼロにするというのがあります。ぴったりゼロにするのは難しいので、ここでは絶対値を 0.1 以下にすることにします。

![38.png](https://qiita-image-store.s3.amazonaws.com/0/211162/f4cfef98-322e-1cb3-43c7-8158b6627f1e.png)

制約条件その２として、各々の観測変数に対して、共通因子の分散と独自因子の分散の和を１にするというのがあります。ぴったり１にするのは難しいので、ここではその和から１を引いた数の絶対値を 0.1 以下にすることにします。

![39.png](https://qiita-image-store.s3.amazonaws.com/0/211162/867554bc-5bb0-cb8a-cb44-9039eccbe58a.png)

ソルバーを立ち上げ、以上の条件を間違えずに入力します。

![40.png](https://qiita-image-store.s3.amazonaws.com/0/211162/f2968378-42f9-ebc3-f5a8-f2816f9b87f5.png)

「解決」を押すと、けっこう時間がかかりますが、答えが返ってきます。収束しない場合は、ランダムな初期値を変更して計算し直してください。

![41.png](https://qiita-image-store.s3.amazonaws.com/0/211162/d61f6169-16ec-aaf6-63bf-92efc30a88d0.png)

収束後の目的変数と制約条件の数の例はこちらになります。

![42.png](https://qiita-image-store.s3.amazonaws.com/0/211162/8b0abdbd-70f3-5233-cb91-35ce3cdffb86.png)

# 最適化後の因子負荷量と因子得点

収束後（最適化後）の因子負荷量と因子得点の例はこちらになります。乱数による初期化を行なっているので、答えが必ずしも同じになるとは限りません。

![43.png](https://qiita-image-store.s3.amazonaws.com/0/211162/a1e4b99a-913a-81c2-e75d-6a43ea8f48a2.png)

# バリマックス回転

因子分析の肝の部分は以上で終了ですが、その後に直交回転を行うことがあります。以下に述べるバリマックス回転はそのひとつで、因子の解釈を容易にする目的で、全因子によって説明できる分散を最大にするように回転します。

これをソルバーで実行するには、以下のようにします。

* まず回転角を乱数で初期化し、回転行列を作成します。
* 回転行列との行列積により、回転前の因子負荷量から回転後の因子負荷量を算出します。
* 「共通性」を算出します。ここでいう共通性とは、因子負荷量の２乗和です。
* 因子負荷量を共通性で割って２乗したものを「寄与度」と言います。その分散を求めます。

※ 「観測値の２乗に対する、因子負荷量の２乗和の割合」のことを「共通性」と呼んでいる文献もあります。また、「共通因子」のことを「共通性」と呼んでいる文献もあります。非常に紛らわしい（単に私が誤解しているだけかもしれませんが）。

![44.png](https://qiita-image-store.s3.amazonaws.com/0/211162/80ac7faf-eb0a-cff4-c616-a366e25e7114.png)

* 寄与度の分散の和が最大になるような回転角をソルバーで求めます。

![45.png](https://qiita-image-store.s3.amazonaws.com/0/211162/f02a5812-84c7-52c2-a9ab-4abc5adc493e.png)

最適化後の値の例はこちらになります。

![46.png](https://qiita-image-store.s3.amazonaws.com/0/211162/ae9330aa-9eb4-a02c-4e99-1a9922119321.png)

# 回転後の因子得点と因子負荷量

因子得点を用いて、各個体（各サンプル、ここではワインの銘柄）をプロットできます。そのプロットの横軸・縦軸がいかほどの意味を持つかは、次の因子負荷量を用いた観測変数のプロットによります。

![48.png](https://qiita-image-store.s3.amazonaws.com/0/211162/2d50e06a-ad75-5460-3211-a7464169de7b.png)

因子負荷量を用いて、各観測変数をプロットできます。因子１は何を意味しそうか、因子２は何を意味しそうか、解釈してみましょう。

![47.png](https://qiita-image-store.s3.amazonaws.com/0/211162/93fe04a0-1952-abc4-9257-b864fd88940e.png)


# 課題

* ここでいう「共通性」は、回転前と回転後でどのように変わるか説明してください。

* 主成分分析による結果と比較してください。

* <a href="https://qiita.com/maskot1977/items/453eddf5aa154c04668c">実習用データ</a> から「都道府県別アルコール類の消費量」を取得し、同様に因子分析を行い、その結果について考察してください。また同様に、主成分分析による結果と比較してください。