色々と考えてみる: 文系のための「多次元データの要約」（１）

2012/09/01

文系のための「多次元データの要約」（１）

今回は、分散と共分散を楽して「行列」を使って計算するという話。

まずは、復習から。分散の計算は、以下の式で表すことができた。
本来は、二乗の形で表すが、共分散の式との比較のため、このように表す。

$\sigma^2=\frac{1}{n-1}\sum_{i=1}^{n}(\bar%20x%20-%20x_i)^2=\frac{1}{n-1}\sum_{i=1}^{n}(\bar%20x%20-%20x_i)(\bar%20x%20-%20x_i)$

一方、共分散の計算は、以下のように表した。
共分散の場合は、二つの変数間のバラツキを表すので、
$x_{1}$ と $x_{2}$ の二つの変数が使われている。

$\textup{cov}_{(x_1,x_2)}=\frac{1}{n-1}\sum_{i=1}^{n}(\bar x_1 - x_{1i})(\bar x_2 - x_{2i})$

さて、ここからが今回の話の重要なポイント。
今回は、ちょっとした、行列の掛け算の魔法を垣間見ることになる。

まずは、「 $x_{1}$ 」の分散の計算から考えてみる。
以下は、「 $x_{1}$ 」の偏差であり、「 $\tilde x_1$ 」という記号で表した。
「〜」の記号は、チルダと呼ぶのだった。

$\boldsymbol{\tilde{x_1}}= \begin{vmatrix} \bar x_1 - x_{11} \\ \bar x_1 - x_{12} \\ \vdots \\ \bar x_1 - x_{1n} \end{vmatrix}$

ここで、 $\boldsymbol{\tilde{x_1}}^t\boldsymbol{\tilde{x_1}}$ というベクトルの計算を考えてみる。これを展開するとどのようになるか？
tの記号は、転置（transpose）を意味したのだから、
$\boldsymbol{\tilde{x_1}}^t\boldsymbol{\tilde{x_1}}= \begin{vmatrix} \bar x_1 - x_{11} && \bar x_1 - x_{12} && \hdots && \bar x_1 - x_{1n} \end{vmatrix} \begin{vmatrix} \bar x_1 - x_{11} \\ \bar x_1 - x_{12} \\ \vdots \\ \bar x_1 - x_{1n} \end{vmatrix}$

このようになる。これは、どのような計算だったか？
忘れた人は、文系のための内積（１）を参照。

確か、「横✕縦＋横✕縦＋横✕縦・・・となるから・・・。」

おぉ、なるほど。二乗したのを足していく訳か。
あるベクトルの転置と元のベクトルの掛け算は「二乗和」になるのだった。

したがって、以下のようになる。

$\boldsymbol{\tilde x}^t \boldsymbol{\tilde x}= (\bar x_1 - x_{11})^2+ (\bar x_1 - x_{12})^2+ \hdots + (\bar x_1 - x_{1n})^2$

うん？これは、確か分散の式のΣの内側と同じ計算となっている。
何が言いたいかと言うと、分散の計算をベクトルで表すと、以下のようになる。

$\sigma^2_{(x_1)}=\frac{1}{n-1}\boldsymbol{\tilde x_1}^t \boldsymbol{\tilde x_1}$

$x_{1}$ の分散は、見事にベクトルの式で表すことができる。
確か、共分散のΣの中身は、添字が異なるだけで、分散と同じだった。
つまり、掛ける側のベクトルを $x_{2}$ とすると、

$\boldsymbol{\tilde x_1}^t \boldsymbol{\tilde x_2}= \begin{vmatrix} \bar x_1 - x_{11} && \bar x_1 - x_{12} && \hdots && \bar x_1 - x_{1n} \end{vmatrix} \begin{vmatrix} \bar x_2 - x_{21} \\ \bar x_2 - x_{22} \\ \vdots \\ \bar x_2 - x_{2n} \end{vmatrix}$

となり、したがって、

$\boldsymbol{\tilde x_1}^t \boldsymbol{\tilde x_2}= (\bar x_1 - x_{11})(\bar x_2 - x_{21})+(\bar x_1 - x_{12})(\bar x_2 - x_{22})+\hdots+(\bar x_1 - x_{1n})(\bar x_2 - x_{2n})$

という計算になる。素晴らしい。分散と同様に見事に計算できる。
分散と同様に、ベクトルの式で表すと、次のようになる。

$\textup{cov}_{(x_1,x_2)}=\frac{1}{n-1}\boldsymbol{\tilde x_1}^t \boldsymbol{\tilde x_2}$

ここまでが理解できれば大丈夫。ここまでは、ベクトルの掛け算であったが、
次のような行列を考えてみるとどうなるか？

$\frac{1}{n-1}\boldsymbol{\tilde X}^t\boldsymbol{\tilde X} & = \frac{1}{n-1}\begin{vmatrix} \bar x_1 - x_{11} & \bar x_1 - x_{12} & \hdots & \bar x_1 - x_{1n}\\ \bar x_2 - x_{21} & \bar x_2 - x_{22} & \hdots & \bar x_2 - x_{2n}\\ \bar x_3 - x_{31} & \bar x_3 - x_{32} & \hdots & \bar x_3 - x_{3n}\\ \end{vmatrix} \begin{vmatrix} \bar x_1 - x_{11} & \bar x_2 - x_{21} & \bar x_3 - x_{31} \\ \bar x_1 - x_{12} & \bar x_2 - x_{22} & \bar x_3 - x_{31} \\ \vdots & \vdots & \vdots \\ \bar x_1 - x_{1n} & \bar x_1 - x_{2n} & \bar x_3 - x_{31}\\ \end{vmatrix} \\$

ちょっと、難しいかもしれないが、落ち着いて考えれば大丈夫。
頑張れる人は、手計算で。自信が無い人は、行列の計算の話を参考に。
頑張れない人は、あとでRを使って練習。

計算の順番と解となる行列での位置は次のようになる。

１行目の横ベクトル×１列目の縦ベクトル（ $\boldsymbol{\tilde{x_1}}^t\boldsymbol{\tilde{x_1}}$ ） ⇒ Answer[1, 1]
１行目の横ベクトル×２列目の縦ベクトル（ $\boldsymbol{\tilde{x_1}}^t\boldsymbol{\tilde{x_2}}$ ） ⇒ Answer[1, 2]
１行目の横ベクトル×３列目の縦ベクトル（ $\boldsymbol{\tilde{x_1}}^t\boldsymbol{\tilde{x_3}}$ ） ⇒ Answer[1, 3]
２行目の横ベクトル×１列目の縦ベクトル（ $\boldsymbol{\tilde{x_2}}^t\boldsymbol{\tilde{x_1}}$ ） ⇒ Answer[2, 1]
２行目の横ベクトル×２列目の縦ベクトル（ $\boldsymbol{\tilde{x_2}}^t\boldsymbol{\tilde{x_2}}$ ） ⇒ Answer[2, 2]
２行目の横ベクトル×３列目の縦ベクトル（ $\boldsymbol{\tilde{x_2}}^t\boldsymbol{\tilde{x_3}}$ ） ⇒ Answer[2, 3]
３行目の横ベクトル×１列目の縦ベクトル（ $\boldsymbol{\tilde{x_3}}^t\boldsymbol{\tilde{x_1}}$ ） ⇒ Answer[3, 1]
３行目の横ベクトル×２列目の縦ベクトル（ $\boldsymbol{\tilde{x_3}}^t\boldsymbol{\tilde{x_2}}$ ） ⇒ Answer[3, 2]
３行目の横ベクトル×３列目の縦ベクトル（ $\boldsymbol{\tilde{x_3}}^t\boldsymbol{\tilde{x_3}}$ ） ⇒ Answer[3, 3]

つまり、以下のように書き換えることが可能。

$\frac{1}{n-1}\boldsymbol{\tilde X}^t\boldsymbol{\tilde X} = \frac{1}{n-1}\begin{vmatrix} \boldsymbol{\tilde x_1}^t\boldsymbol{\tilde x_1} & \boldsymbol{\tilde x_1}^t\boldsymbol{\tilde x_2} & \boldsymbol{\tilde x_1}^t\boldsymbol{\tilde x_3} \\ \boldsymbol{\tilde x_2}^t\boldsymbol{\tilde x_1} & \boldsymbol{\tilde x_2}^t\boldsymbol{\tilde x_2} & \boldsymbol{\tilde x_2}^t\boldsymbol{\tilde x_3} \\ \boldsymbol{\tilde x_3}^t\boldsymbol{\tilde x_1} & \boldsymbol{\tilde x_3}^t\boldsymbol{\tilde x_2} & \boldsymbol{\tilde x_3}^t\boldsymbol{\tilde x_3} \\ \end{vmatrix} = \begin{vmatrix} \textup{cov}_{(x_1,x_1)} & \textup{cov}_{(x_1,x_2)} & \textup{cov}_{(x_1,x_3)} \\ \textup{cov}_{(x_2,x_1)} & \textup{cov}_{(x_2,x_2)} & \textup{cov}_{(x_2,x_3)} \\ \textup{cov}_{(x_3,x_1)} & \textup{cov}_{(x_3,x_2)} & \textup{cov}_{(x_3,x_3)} & \end{vmatrix}$

ここで、注意して見てもらいたいのが、対角成分に何が入っているか？ということ。
自分同士の共分散ということは、何を意味するのだったか？答えは、分散。

もうひとつのポイントは、二つの変数の共分散の計算は、
順番を変えても答えは同じ（スカラーでは、A×B＝B×A）なので、
対角成分を挟んで上側と下側が同じになる、ということも重要。

したがって、この行列の計算を通して計算された結果は、
対角成分に分散が並び、それ以外に共分散が入った対象行列である。

これを、「分散共分散行列」と呼ぶ。

さらに、元の行列が標準化されていた場合、すなわち、
平均が「0」で分散が「1」となるように基準化されていた場合には、
対角成分が全て「1」で、それ以外の成分が相関係数となる。

例えば、先ほどの行列 X を標準化した行列を Z としたとき、

$\frac{1}{n-1}\boldsymbol{Z}^t\boldsymbol{Z}=\begin{vmatrix} 1 & \textup{cor}_{(x_1,x_2)} & \textup{cor}_{(x_1,x_3)}\\ \textup{cor}_{(x_2,x_1)} & 1 & \textup{cor}_{(x_2,x_3)}\\ \textup{cor}_{(x_3,x_1)} & \textup{cor}_{(x_3,x_2)} & 1 \end{vmatrix}$

という関係が成り立つ。これを、「相関係数行列」と呼ぶ。

分散共分散行列と相関係数行列は、データを要約した指標として重要であると同時に、
この行列を用いることで、主成分分析の計算が容易になる。この話は、いずれ。

では、今度は実際のデータを用いて計算の練習をしてみる。
サンプルデータは、分散と標準偏差の説明で用いたものを再び使う。

# Windows と Linux の人は次のコマンド
X <- read.table("clipboard", sep=",", header=TRUE)

# Mac の人は次のコマンド
X <- read.table(pipe("pbpaste") , sep=",", header=TRUE)

以下が、そのデータ。前回のものを再掲しただけ。
データの説明に関しては、平均の話を参照。

Oroshi_Kg, Taka, Naka, Yasu
Tai, 1336, 3150, 1217, 53
Chinu, 226, 630, 513, 105
Sawara, 212, 1575, 1259, 840
Yazu, 4834, 840, 315, 179
Suzuki, 618, 1575, 1146, 525
Akou, 21, 4725, 3129, 1575
Kochi, 28, 2100, 874, 210
Okoze, 28, 5250, 2635, 210
Ainame, 29, 3150, 621, 315
Hage, 1205, 1890, 383, 210
Konoshiro, 21, 2100, 1100, 105
Sayori, 12, 5250, 3916, 1260
Anago, 1637, 2625, 1721, 630
Mebaru, 330, 4200, 1680, 315
Tachiuo, 1211, 2310, 930, 105
Hamo, 1622, 3938, 592, 53
Karei, 41, 6300, 3178, 525
Hirame, 540, 2100, 1496, 210
Aji, 5149, 3150, 789, 210
Saba, 5413, 3675, 625, 53
Ika, 2414, 2888, 1083, 105
Tako, 1844, 1890, 1180, 525
Ebi, 560, 7350, 1420, 525
KurumaEbi, 222, 8925, 6508, 2625
Koiwashi, 1316, 1155, 617, 328
ObaIwashi, 2716, 499, 267, 175
Mentai, 678, 1155, 859, 394
Hamachi, 4772, 998, 632, 105
Isaki, 268, 2625, 1465, 840

まずは、行列 X と同じサイズで、平均が並んでいる行列が必要だから、

# ステップ1：データの総数（n）を求める。
n <- nrow(X)
# ステップ2：卸売数量の平均を求める。転置が必要。
m <- t(colMeans(X))
# ステップ3：行列の引き算ができるように、平均値を総数個複製する。
M <- matrix(rep(m, each=n), nrow=n, ncol=4)
# ステップ4：元の行列から、平均値が複製された行列Mを引き算する。
tilde_X <- as.matrix(X-M)

最終的に、「tilde_X」という変数に偏差行列が入った。
ここで、分散共分散の式をそのまま当てはめると、

(t(tilde_X) %*% tilde_X)/(n-1)

となり、実行結果は以下のようになる。

> (t(tilde_X) %*% tilde_X)/(n-1)
Oroshi_Kg Taka Naka Yasu
Oroshi_Kg 2814792.6 -1026368.1 -920404.3 -335784.2
Taka -1026368.1 4210888.7 2159783.0 685439.7
Naka -920404.3 2159783.0 1752709.8 637497.8
Yasu -335784.2 685439.7 637497.8 306021.7

疑い深い人は、分散と共分散の話の結果と照合してみよう。

間違いなく、対角成分には分散が並び、それ以外には共分散が入っている。

さて、今度は、元のデータを標準化して計算をしてみる。

# 標準化した行列を作る

Z <- scale(X)

# 共分散の計算。平均は0なので偏差行列は不要

(t(Z) %*% Z)/(n-1)

実行結果は以下の通り。

> cov(Z)

Oroshi_Kg Taka Naka Yasu

Oroshi_Kg 1.0000000 -0.2981213 -0.4143816 -0.3617937

Taka -0.2981213 1.0000000 0.7950020 0.6038183

Naka -0.4143816 0.7950020 1.0000000 0.8704575

Yasu -0.3617937 0.6038183 0.8704575 1.0000000

対角成分には「1」が並んでいて、それ以外の所は相関係数が入っている。

このようにして、相関係数行列を計算することができる。

ちなみに、前回、登場したcov()関数とcor()関数は、

両方とも行列に対応しているので、共分散行列と相関係数行列を出せる。

cov(X) # 分散共分散行列を直接計算する

cor(X) # 相関係数行列を直接計算する

以下は、実行結果。

> cov(X) # 分散共分散行列を直接計算する

Oroshi_Kg Taka Naka Yasu

Oroshi_Kg 2814792.6 -1026368.1 -920404.3 -335784.2

Taka -1026368.1 4210888.7 2159783.0 685439.7

Naka -920404.3 2159783.0 1752709.8 637497.8

Yasu -335784.2 685439.7 637497.8 306021.7

> cor(X) # 相関係数行列を直接計算する

Oroshi_Kg Taka Naka Yasu

Oroshi_Kg 1.0000000 -0.2981213 -0.4143816 -0.3617937

Taka -0.2981213 1.0000000 0.7950020 0.6038183

Naka -0.4143816 0.7950020 1.0000000 0.8704575

Yasu -0.3617937 0.6038183 0.8704575 1.0000000

手計算の場合には、結局、展開して計算しないといけないので、

手間は同じであるが、Rのように行列の計算ができるソフトウェアを使えば、

分散、共分散、相関係数を一度に計算することができ、便利なのである。

10 件のコメント:

文系A2013/03/26 16:28:00
「まずは、復習から。分散の計算は、以下の式で表すことができた。
本来は、二乗の形で表すが、共分散の式との比較のため、このように表す。」

の下の式は、シグマの下が"i=i"になってますが？

\sigma^2=\frac{1}{n-1}\sum_{i=i}^{n}(\bar x - x_i)(\bar x - x_i)
返信削除
返信
文系A2013/03/26 16:45:00
さらに同じ式で、"(\bar x - x_i"が２回あるのは正しいのでしょうか？また、この部分は"x_i - \bar x"が正しいのではないでしょうか？
返信削除
返信
Yu Fujimoto2013/03/26 20:08:00
文系Aさんから上記のような質問があったので、
誤解が少しでも軽減するように、
最初の分散の式を少し書き換えました。
返信削除
返信
文系A2013/03/27 11:22:00
そうでした。二乗でしたね。書き方が少し変わっただけで混乱してしまうとは、ちゃんと理解できていない証拠ですね。お恥ずかしい。ご説明ありがとうございます。

ということは、"i=i"も正しいのですね。
返信削除
返信
文系A2013/03/27 15:27:00
以下のAnswer[2,1]は正しくはAnswer[2,2]でしょうか？

２行目の横ベクトル×２列目の縦ベクトル（） ⇒ Answer[2, 1]
返信削除
返信
文系A2013/03/27 15:47:00
以下の部分は、実行結果のコマンドが見本と異なりますが、同じ事ができるという

cov(Z)
は
Z <- scale(X)
(t(Z) %*% Z)/(n-1)
と同じ結果を返す、ということでしょうか？

qt.
# 共分散の計算。平均は0なので偏差行列は不要
(t(Z) %*% Z)/(n-1)

実行結果は以下の通り。

> cov(Z)
Oroshi_Kg Taka Naka Yasu
Oroshi_Kg 1.0000000 -0.2981213 -0.4143816 -0.3617937
Taka -0.2981213 1.0000000 0.7950020 0.6038183
Naka -0.4143816 0.7950020 1.0000000 0.8704575
Yasu -0.3617937 0.6038183 0.8704575 1.0000000
unqt.
返信削除
返信

コメントを追加