コレポンは何をやっているのか②検討の流れ - 統計学といくつかのよしなしごと

前回

　前回の記事でマーケティングリサーチにおいて多用されるコレスポンデンス分析について、その可視化までの数理を概観した。今回は出力されたマップやその他の統計量を用いた検討のポイントを見ていく。前回同様主な参考文献は以下。

作者: Sten-Erik Clausen,藤本一男
出版社/メーカー: オーム社
発売日: 2015/12/01
メディア: 単行本（ソフトカバー）
この商品を含むブログ (1件) を見る

「対応分析法・数量化Ⅲ類の考え方」
http://www.wordminer.org/wp-content/uploads/2013/04/63_18.pdf

f:id:yanbow221:20180717203905j:plain — 前回出力した同時布置図

固有値と寄与率、累積寄与率

　コレスポンデンス分析によってサンプル、銘柄ともに5次元の成分スコアが得られ、上記の同時布置図ではそのうち第1成分、第2成分を軸に取った。このように取っ掛かりとしては出力された成分の順に沿って2軸を取ることが多いが、これは対応する固有値の大きさに由来している。5成分の固有値を確認してみると、

> result_Xi$values[1:5] # 行方向に主成分分析を行って得た固有値
[1] 0.62603096 0.18766678 0.13451054 0.04519841 0.02603776
> result_Xj$values[1:5] # 列行方向に主成分分析を行って得た固有値
[1] 0.62603096 0.18766678 0.13451054 0.04519841 0.02603776

行方向、列方向それぞれで得た固有値の値が一致している(6成分以降の固有値はほとんど0になっている)。解いた固有値問題が元データに対していくつか操作を加えて得たデータの分散共分散行列に対するものだったことを思い出すと、この固有値の総和はデータが持っていた分散の総和に等しいことが分かる。実際それぞれ計算してみると、

> sum(diag(Vxi)) # 元のデータの分散の総和（行方向）
[1] 1.019444
> sum(diag(Vxj)) # 元のデータの分散の総和（列方向）
[1] 1.019444
> sum(result_Xi$values) # 固有値の総和（行方向）
[1] 1.019444
> sum(result_Xj$values) # 固有値の総和（列方向）
[1] 1.019444

となっていてこれが確認できる。なお、分散の総和と固有値の総和が等しくなるメカニズムについては、例によって以下で丁寧に解説してある。

多変量解析法入門 (ライブラリ新数学大系)

作者: 永田靖,棟近雅彦
出版社/メーカー: サイエンス社
発売日: 2001/04/01
メディア: 単行本
購入: 2人クリック: 9回
この商品を含むブログ (1件) を見る

　このように、元のデータが持つ分散をできるだけ表現できる軸を順に求めなおしてやって作った空間上に再配置するのが主成分分析である。というわけでこの固有値を使って、元のデータをどのくらい表現できているかを評価するための統計量が定義される。

${\displaystyle v_k = \frac {λ_k}{\sum^K_{k=1} λ_k} × 100(\%) \ \ \ k = 1,2,...K \ \ \ (第k成分の寄与率) \\ \\ v_{k'} = \frac {\sum^{k'}_{k=1} λ_k}{\sum^K_{k=1} λ_k} × 100(\%) \ \ \ 1 \leq k' \leq K \ \ \ (第k'成分までの累積寄与率) }$

　今回は第1，2成分までで分散の約80％を説明できていることになるので、これらで軸をとってやってデータの意味付けを検討していくのが基本になる。このように元のデータの持つデータ項目が多く考えにくい時に、データの持つ情報をなるべく担保しつつ考えやすいより少ない次元に落とし込んでやろう*1というのが主成分分析とかコレスポンデンス分析の意義である。軸の選択については「累積寄与率がX%以上(領域によって慣例的に囁かれている)」とか「スクリープロット*2」の利用等いろいろと成分の採否に関する考え方があるが、コレスポンデンス分析の場合はマッピングに基づく解釈がより重要になるので、寄与率の大きな軸だけで見ずにいろいろと見てみるのが良いと思う*3。

同時布置図の検討

　そんな形で作成した同時布置図を見て、あれこれ「軸」と「点」について解釈していくのが次の段階になる。検討に際し有用な統計量は大きく2つで、「ある軸の形成に対して各々の行項目*4や列項目*5がどの程度寄与しているか」を示す絶対寄与度と、「各々の行項目/列項目が各軸(成分)によってどの程度表現されているか」を示す平方相関(相対寄与度)がある。

(1)絶対寄与度
　各行項目への絶対寄与度 $C_k(i)$ 、各列項目への絶対寄与度 $C_k(j)$ はそれぞれ以下で計算できる。
${\displaystyle C_k(i) = \frac {p_{i+}(z_{ik})^2}{λ_k}, \ \ \ \sum^{m}_{i=1} C_k(i) = 1 }$

${\displaystyle C_k(j) = \frac {p_{+j}(z_{jk}^*)^2}{λ_k}, \ \ \ \sum^{n}_{j=1} C_k(j) = 1 }$

${\displaystyle \left( i \in I,\ \ j \in J,\ \ k = 1,2,...,K,\ \ K = min\{m, n\} - 1 \right) }$

　それぞれ求めてみるとこんな感じ。

> round((Pi %*% (rscore2 ^ 2) %*% diag(1 / result_Xj$values[1:5])) * 100, 2)
       [,1]  [,2]  [,3]  [,4]  [,5]
 [1,]  0.47  4.85 34.84  3.34 10.87
 [2,]  9.26 12.15 30.48 18.59  7.57
 [3,] 14.67  5.67  8.17  1.56  0.85
 [4,]  3.16  3.87  0.11  8.16 23.86
 [5,]  8.58  1.36  1.06 17.63 33.29
 [6,]  1.30 13.29 15.36  9.15  0.05
 [7,] 18.25  4.65  0.21 14.78 11.34
 [8,] 14.67  5.67  8.17  1.56  0.85
 [9,]  7.22 17.15  0.99 24.64  4.36
[10,] 22.43 31.35  0.60  0.61  6.95
> round((Pj %*% (cscore2 ^ 2) %*% diag(1 / result_Xj$values[1:5])) * 100, 2)
      [,1]  [,2]  [,3]  [,4]  [,5]
[1,] 19.53  0.01  4.77 46.71 15.20
[2,]  0.44 32.65  4.12 22.15 19.95
[3,] 17.61  0.01 61.52  1.47  2.14
[4,] 18.94 49.50  0.18 23.97  0.51
[5,]  8.68  4.14  4.82  0.66 61.02
[6,] 34.80 13.69 24.60  5.04  1.18

　例えば第1成分の軸については、縦方向に見れば順にサンプル10, 7, 3等が、横方向に見れば順に銘柄F, D, C等が軸の形成に寄与していることが分かる。同様にして第2成分の軸については、縦方向に見れば順にサンプル10, 9, 6, 2等が、横方向に見れば順に銘柄D, B, F等が寄与している。一方でサンプル1や2は第1成分、2成分には寄与しておらず、むしろ第3成分の形成に寄与していることが分かる。銘柄Cも第3成分の形成において大きく寄与している。こういったことを確認しながら、必要に応じて軸の解釈も行いつつデータ項目の関係性を確認していく。

(2)平方相関(相対寄与度)
　各行項目に対する平方相関 $C^*_k(i)$ 、各列項目に対する平方相関 $C^*_k(j)$ はそれぞれ以下で計算できる。
${\displaystyle C^*_k(i) = \frac {z_{ik}^2}{ \sum^{n}_{j=1} p_{+j} \left( \frac {p_{ij} - p_{i+} p_{+j}}{p_{i+}p_{+j}} \right) } }$

${\displaystyle C^*_k(j) = \frac {(z_{ik}^*)^2}{ \sum^{m}_{i=1} p_{i+} \left( \frac {p_{ij} - p_{i+} p_{+j}}{p_{i+}p_{+j}} \right) } }$

${\displaystyle \left( i \in I,\ \ j \in J,\ \ k = 1,2,...,K,\ \ K = min\{m, n\} - 1 \right) }$

　こんな感じで求まる。

> round(diag(1 / apply((rscore2 ^ 2), 1, sum)) %*% (rscore2 ^ 2) * 100, 2)
       成分スコア
        第1成分 第2成分 第3成分 第4成分 第5成分
   [1,]    4.61   14.40   74.13    2.39    4.48
   [2,]   43.88   17.26   31.02    6.36    1.49
   [3,]   80.28    9.30    9.61    0.61    0.19
   [4,]   53.31   19.59    0.39    9.95   16.76
   [5,]   72.27    3.43    1.92   10.72   11.66
   [6,]   14.08   43.08   35.68    7.14    0.02
   [7,]   85.97    6.57    0.21    5.02    2.22
   [8,]   80.28    9.30    9.61    0.61    0.19
   [9,]   49.68   35.36    1.47   12.24    1.25
  [10,]   69.46   29.11    0.40    0.14    0.89
> round(diag(1 / apply((cscore2 ^ 2), 1, sum)) %*% (cscore2 ^ 2) * 100, 2)
      成分スコア
       第1成分 第2成分 第3成分 第4成分 第5成分
  [1,]   79.51    0.02    4.17   13.73    2.57
  [2,]    3.25   72.28    6.54   11.81    6.13
  [3,]   56.76    0.01   42.60    0.34    0.29
  [4,]   53.25   41.72    0.11    4.87    0.06
  [5,]   64.11    9.16    7.64    0.35   18.74
  [6,]   78.02    9.20   11.85    0.82    0.11

　サンプル3, 7, 8、銘柄1や6等は第1成分によってよく表現されていることがわかる。ここで注意しないといけないのは、例えばサンプル1を同時布置上で見てみると真ん中あたりにプロットされているが、これを第1成分、第2成分で定義された空間上での「平均」と解釈するのは誤りである。サンプル1の第1成分、第2成分への絶対寄与度は0.47%, 4.85%となっておりそもそも軸の形成に寄与していない。また平方相関を見ても第1成分、第2成分それぞれ4.61%, 14.40%でありこれらの軸はサンプル1にそれ程強く寄与していないことが分かる。サンプル1においてはむしろ第3成分との関連が強くあり、これはつまり第1成分、第2成分で描かれた空間ではサンプル1は「表現されていない」と解釈した方が良さそうということになる。
　このように、同時布置図でプロットされた全ての点の位置関係をもって素朴に解釈すると思わぬ誤りを犯してしまうことがあるので、求めた成分軸それぞれで同時布置図と統計量を確認することが重要となる。以下は「対応分析法・数量化Ⅲ類の考え方」から引用。

ここで留意すべきことは、散布図の情報は"そこでいま眺めている成分軸の組み合わせの中での射影図"であり多次元データとしての全情報ではない、ということである。よって、絶対寄与度、相対寄与度も図と併せて観察し、どの成分でどう寄与するかを多次元的に知ることが必要である。

その他いくつかの補足

(1)サプリメンタリーポイント
　対応分析は行方向の分析と列方向の分析を遷移方程式を使って関係づけるが、この式によって追加処理を行いデータの解釈を助けることができる。例えば先の10サンプルについて実は「性別」という属性があったとする。つまりデータセットがこんな感じであったとする。

> sex_sample <- matrix(c(1, 0, 1, 1, 1, 0, 0, 1, 0, 0,
+                        0, 1, 0, 0, 0, 1, 1, 0, 1, 1), nrow = 10, ncol = 2, byrow = F)
> colnames(sex_sample) <- c("男性", "女性")
> (F_col_added <- cbind(F, sex_sample))
           銘柄A 銘柄B 銘柄C 銘柄D 銘柄E 銘柄F 男性 女性
サンプル1      0     1     0     0     1     1    1    0
サンプル2      0     0     0     0     0     1    0    1
サンプル3      0     0     1     0     0     1    1    0
サンプル4      0     1     1     0     1     1    1    0
サンプル5      0     1     1     0     0     1    1    0
サンプル6      1     1     1     0     1     0    0    1
サンプル7      1     1     0     1     1     0    0    1
サンプル8      0     0     1     0     0     1    1    0
サンプル9      1     1     0     0     1     0    0    1
サンプル10     1     0     0     1     1     0    0    1

　遷移方程式を使うことで、銘柄項目に対して算出した成分スコアを担保した状態で新たに性別項目の成分スコアを求めることができる。図と併せてこんな感じ。

> # 追加処理（サプリメンタリーポイント）
> P_col_added <- F_col_added / sum(F_col_added)
> pj_col_added <- matrix(apply(P_col_added, 2, sum))
> C_col_added <- diag(1 / pj_col_added[, 1]) %*% t(P_col_added)
> Dm1 <- diag(1 / sqrt(result_Xi$values[1:5]))
> cscore_col_added <- C_col_added %*% rscore2 %*% Dm1 # 遷移方程式
> 
> # 銘柄の成分スコアが担保されていることの確認
> cscore2 # 元々算出した成分スコア
       成分スコア
銘柄       第1成分     第2成分     第3成分     第4成分     第5成分
  銘柄A -0.9413812 -0.01306953 -0.21557326  0.39122449  0.16936819
  銘柄B -0.1153204  0.54418081  0.16369085 -0.21996465  0.15846620
  銘柄C  0.7996653  0.01204477 -0.69280014 -0.06201891 -0.05690933
  銘柄D -1.3113409 -1.16064830 -0.05848639 -0.39637626  0.04373046
  銘柄E -0.5125064  0.19370276  0.17693374  0.03799750 -0.27710881
  銘柄F  1.0261401 -0.35232510  0.39991984  0.10495867  0.03857811
> cscore_col_added # 性別を付与したデータセットについて算出した成分スコア
           [,1]        [,2]        [,3]        [,4]        [,5]
[1,] -0.9413812 -0.01306953 -0.21557326  0.39122449  0.16936819
[2,] -0.1153204  0.54418081  0.16369085 -0.21996465  0.15846620
[3,]  0.7996653  0.01204477 -0.69280014 -0.06201891 -0.05690933
[4,] -1.3113409 -1.16064830 -0.05848639 -0.39637626  0.04373046
[5,] -0.5125064  0.19370276  0.17693374  0.03799750 -0.27710881
[6,]  1.0261401 -0.35232510  0.39991984  0.10495867  0.03857811
[7,]  0.9035441 -0.04731066 -0.11472595 -0.33848483 -0.25003063
[8,] -0.4252808 -0.38593508  0.42217115  0.77741483  0.43181891
> 
> # 追加処理後の同時布置図
> library(gridExtra)
> 
> dimnames(cscore_col_added) <- list(銘柄 = c("銘柄A", "銘柄B", "銘柄C", "銘柄D", "銘柄E", "銘柄F", "男性", "女性"),
+                             成分スコア = c("第1成分", "第2成分", "第3成分", "第4成分", "第5成分"))
> dat_map_col_added <- rbind(rscore2, cscore_col_added)
> dat_map_col_added <- transform(dat_map_col_added, label = c(rep("サンプル", 10), rep("銘柄", 6), rep("性別", 2)))
> 
> map2 <- ggplot(dat_map_col_added, aes(x = 第1成分, y = 第2成分, color = label)) + 
+   xlim(-1.5, 1.5) + ylim(-1.5, 1.5) +
+   geom_hline(yintercept = 0, color = "black") +
+   geom_vline(xintercept = 0, color = "black") +
+   geom_point(size = 5) +
+   scale_color_manual(
+     values = c("サンプル" = "tomato", "銘柄" = "skyblue", "性別" = "mediumseagreen")
+   ) +
+   geom_text(aes(label = rownames(dat_map_col_added)), vjust = -1, size = 5) +
+   theme_bw()
> 
> grid.arrange(map1, map2, ncol = 2) # 左図:元の同時布置図/右図:性別を追加処理した同時布置図