ブログ

「統計検定1級テキスト」 練習問題解答例 Ch.3 「統計的推定」

当記事は「統計検定$1$級対応 統計学(東京図書)」の読解サポートにあたって第$3$章の「統計的推定」に関して演習問題を中心に解説を行います。十分統計量、分解定理、尤度、不偏性・一致性・有効性、クラメル・ラオの不等式、AICなど重要トピックが多いので、抑えておくと良いと思います。

本章のまとめ

練習問題解説

問$3$.$1$

$[1]$
$y = \log{(1+x)}$は下記のように$x$に関して解くことができる。
$$
\large
\begin{align}
y &= \log{(1+x)} \\
e^y &= 1 + x \\
x &= e^y – 1
\end{align}
$$

また、$x$を$y$で微分を行うと下記が得られる。
$$
\large
\begin{align}
\frac{dx}{dy} &= \frac{d}{dy} (e^y – 1) \\
&= e^y
\end{align}
$$

上記を元に確率変数$Y$に関する確率密度関数$g(y)$は下記のように得られる。
$$
\large
\begin{align}
g(y) &= f(e^{y}-1) \frac{dx}{dy} \\
&= \theta (1+e^{y}-1)^{-(1+\theta)} \times e^y \\
&= \theta e^{-\theta y}
\end{align}
$$

上記の期待値を$E[Y]$、分散を$V[Y]$とおき、以下それぞれの計算を行う。
・$E[Y]$
$$
\large
\begin{align}
E[Y] &= \int_{0}^{\infty} y g(y) dy \\
&= \int_{0}^{\infty} y \theta e^{-\theta y} dy \\
&= -\frac{1}{\theta} \left[ y \theta e^{-\theta y} \right]_{0}^{\infty} + \int_{0}^{\infty} e^{-\theta y} dy \\
&= 0 – \frac{1}{\theta} \left[ e^{-\theta y} \right]_{0}^{\infty} \\
&= \frac{1}{\theta}
\end{align}
$$

・$V[Y]$
$V[Y]=E[Y^2]-E[Y]^2$を用いるにあたって$E[Y^2]$を先に計算する。
$$
\large
\begin{align}
E[Y^2] &= \int_{0}^{\infty} y^2 g(y) dy \\
&= \int_{0}^{\infty} y^2 \theta e^{-\theta y} dy \\
&= -\frac{1}{\theta} \left[ y^2 \theta e^{-\theta y} \right]_{0}^{\infty} + \frac{2}{\theta} \int_{0}^{\infty} y \theta e^{-\theta y} dy \\
&= 0 + \frac{2}{\theta}E[Y] \\
&= \frac{2}{\theta^2}
\end{align}
$$

よって$V[Y]$は下記のように計算できる。
$$
\large
\begin{align}
V[Y] &= E[Y^2] – E[Y]^2 \\
&= \frac{2}{\theta^2} – \left( \frac{1}{\theta} \right)^{2} \\
&= \frac{2}{\theta^2} – \frac{1}{\theta^2} = \frac{1}{\theta^2}
\end{align}
$$

$[2]$

$[3]$
下記のように$T$の期待値$E[T]$を計算することができる。
$$
\large
\begin{align}
E[T] &= E \left[ \frac{1}{n} \sum_{i=1}^{n} \log{(X_i+1)} \right] \\
&= E \left[ \frac{1}{n} \sum_{i=1}^{n} Y_i \right] \\
&= \frac{1}{n} E \left[ \sum_{i=1}^{n} Y_i \right] \\
&= \frac{1}{n} \times \frac{n}{\theta} \\
&= \frac{1}{\theta}
\end{align}
$$

上記より$T$は$\displaystyle \frac{1}{\theta}$の不偏推定量であることがわかる。

また、$T$の分散$V[T]$は下記のように得られる。
$$
\large
\begin{align}
V[T] &= V \left[ \frac{1}{n} \sum_{i=1}^{n} Y_i \right] \\
&= \frac{1}{n^2} V \left[ \sum_{i=1}^{n} Y_i \right] \\
&= \frac{1}{n^2} \times n V[Y_i] \\
&= \frac{1}{n} \times \frac{1}{\theta^2} \\
&= \frac{1}{n \theta^2}
\end{align}
$$

ここで$T$が$\displaystyle \frac{1}{\theta} = \alpha$の推定量であると考えるとき、$y_1,…,y_n$に基づくパラメータ$\alpha$に関する尤度関数$L(\alpha)$は下記のように得られる。
$$
\large
\begin{align}
L(\alpha) &= g(y_1,y_2,…,y_n) \\
&= \prod_{i=1}^{n} g(y_i) \\
&= \prod_{i=1}^{n} \frac{1}{\alpha} e^{-(y/\alpha)}
\end{align}
$$

上記より対数尤度$\log{L(\alpha)}$は下記のように考えることができる。
$$
\large
\begin{align}
\log{L(\alpha)} &= \log \left[ \prod_{i=1}^{n} \frac{1}{\alpha} e^{-(y_i/\alpha)} \right] \\
&= – n \log{\alpha} – \frac{\displaystyle \sum_{i=1}^{n} y_i}{\alpha} \\
&= – n \log{\alpha} – \frac{n t}{\alpha}
\end{align}
$$

上記に対し、対数尤度$\log{L(\alpha)}$の二階微分を考える。
$$
\large
\begin{align}
\frac{\partial \log{L(\alpha)}}{\partial \alpha} &= -\frac{n}{\alpha} + \frac{n t}{\alpha^2} \\
\frac{\partial^2 \log{L(\alpha)}}{\partial \alpha^2} &= \frac{n}{\alpha^2} – \frac{2nt}{\alpha^3}
\end{align}
$$

よってパラメータ$\alpha$に関するフィッシャー情報量$\mathit{I}_{n}(\alpha)$は下記のように考えることができる。
$$
\large
\begin{align}
\mathit{I}_{n}(\alpha) &= – E \left[ \frac{\partial^2 \log{L(\alpha)}}{\partial \alpha^2} \right] \\
&= – E \left[ \frac{n}{\alpha^2} – \frac{2nT}{\alpha^3} \right] \\
&= – \frac{n}{\alpha^2} + \frac{2n}{\alpha^3} E[T] \\
&= – \frac{n}{\alpha^2} + \frac{2n}{\alpha^3} \times \alpha \\
&= \frac{n}{\alpha^2}
\end{align}
$$

ここで$\displaystyle V[T]=\frac{1}{n \theta^2}=\frac{\alpha^2}{n}=\mathit{I}_{n}(\alpha)^{-1}$より$V[T]$は$\displaystyle \alpha = \frac{1}{\theta}$に関するクラメル・ラオの不等式の下限と一致することが確認できる。また、$T$は有効推定量であると考えることができる。

問$3$.$2$

$[1]$
標本平均の期待値$E[\bar{X}]$に関して下記が成立する。
$$
\large
\begin{align}
E[\bar{X}] &= E \left[ \frac{T}{n} \right] \\
&= E \left[ \frac{X_1+X_2+…+X_n}{n} \right] \\
&= E[X] \\
&= \int_{0}^{\infty} x \times \lambda e^{-\lambda x} dx \\
&= \left[ x \frac{\lambda}{-\lambda} e^{-\lambda x} \right]_{0}^{\infty} + \int_{0}^{\infty} e^{-\lambda x} dx \\
&= 0 – \frac{1}{\lambda} \left[ e^{-\lambda x} \right]_{0}^{\infty} \\
&= \frac{1}{\lambda}
\end{align}
$$

また、標本平均の分散$V[\bar{X}]$に関して下記が成立する。
$$
\large
\begin{align}
V[\bar{X}] &= V \left[ \frac{X_1+X_2+…+X_n}{n} \right] \\
&= \frac{1}{n^2} \times n V[X_1] = \frac{1}{n} V[X] \\
&= \frac{1}{n} (E[X^2] – E[X]^2) \quad (1)
\end{align}
$$

$E[X]$はすでに計算を行ったので、$E[X^2]$に関して計算することを考える。
$$
\large
\begin{align}
E[X^2] &= \int_{0}^{\infty} x^2 \times \lambda e^{-\lambda x} dx \\
&= \frac{1}{\lambda} \left[ x^2 \lambda e^{-\lambda x} \right]_{0}^{\infty} + \frac{2}{\lambda} \int_{0}^{\infty} x \times \lambda e^{-\lambda x} dx \\
&= 0 + \frac{2}{\lambda} E[X] = \frac{2}{\lambda^2} \quad (2)
\end{align}
$$

$(2)$式を$(1)$式を代入することで$V[\bar{X}]$は下記のように得られる。
$$
\large
\begin{align}
V[\bar{X}] &= \frac{1}{n} (E[X^2] – E[X]^2) \\
&= \frac{1}{n} \left( \frac{2}{\lambda^2} – \left( \frac{1}{\lambda} \right)^2 \right) \\
&= \frac{1}{n \lambda^2}
\end{align}
$$

$[2]$
観測値$x_1,…,x_n$に基づいてパラメータ$\lambda$の尤度関数$L(\lambda)$は下記のように得られる。
$$
\large
\begin{align}
L(\lambda) &= P(x_1,x_2,…,x_n|\lambda) \\
&= \prod_{i=1}^{n} P(x_i|\lambda) = \prod_{i=1}^{n} \lambda e^{-\lambda x_i} \\
&= \prod_{i=1}^{n} \exp(\log{\lambda}) \exp(-\lambda x_i) \\
&= \prod_{i=1}^{n} \exp(\log{\lambda} – \lambda x_i) \\
&= \exp \left[ \sum_{i=1}^{n} (\log{\lambda} – \lambda x_i) \right] \\
&= \exp \left( n \log{\lambda} – \lambda \sum_{i=1}^{n} x_i \right) \\
&= \exp \left( n \log{\lambda} – \lambda t \right)
\end{align}
$$

上記より対数尤度$\log{L(\lambda)}$は下記のように得られる。
$$
\large
\begin{align}
\log{L(\lambda)} &= \log \left[ \exp \left( n \log{\lambda} – \lambda t \right) \right] \\
&= n \log{\lambda} – \lambda t
\end{align}
$$

ここで$\log{L(\lambda)}$を$\lambda$で偏微分すると下記が得られる。
$$
\large
\begin{align}
\frac{\partial \log{L(\lambda)}}{\partial \lambda} &= \frac{\partial}{\partial \lambda} \left( n \log{\lambda} – \lambda t \right) \\
&= \frac{n}{\lambda} – t
\end{align}
$$

上記は$\lambda$に関して単調減少であるので、$\lambda$の最尤推定量$\hat{\lambda}$に関して$\displaystyle \frac{n}{\hat{\lambda}} – t = 0$が成立する。
$$
\large
\begin{align}
\frac{n}{\hat{\lambda}} – t &= 0 \\
\frac{n}{\hat{\lambda}} &= t \\
\hat{\lambda} &= \frac{n}{t}
\end{align}
$$

$[3]$
フィッシャー情報量$\mathit{I}_{n}(\lambda)$は下記のように得られる。
$$
\large
\begin{align}
\mathit{I}_{n}(\lambda) &= – E \left[ \frac{\partial^2 \log{L(\lambda)}}{\partial \lambda^2} \right] \\
&= – E \left[ \frac{\partial}{\partial \lambda} \left( \frac{n}{\lambda} – t \right) \right] \\
&= E \left[ \frac{n}{\lambda^2} \right] = \frac{n}{\lambda^2}
\end{align}
$$

$[4]$

参考

・統計検定1級 統計数理 関連まとめ
https://www.hello-statisticians.com/stat_certifi_1_math

ディリクレ分布(Dirichlet distribution)と多項分布の事後分布

ディリクレ分布(Dirichlet distribution)は多項分布のパラメータの事後確率を考えるときに用いる共役事前分布です。当記事では多項分布の式を確認し、その共役事前分布のディリクレ分布の導入やディリクレ分布を用いた事後分布の計算に関して取り扱います。
「パターン認識と機械学習」の$2$.$2$節の「多項分布」や$2$.$2$.$1$節の「ディリクレ分布」を参考に作成を行いました。

多項分布

多項分布の定義と確率関数

多項分布(multinomial distribution)は二項分布の拡張であると考えることができる。二項分布ではコインの「表/裏」のような$2$値を取り扱ったが、多項分布ではサイコロの出目のような$K$値を取り扱う。

二項分布ではコインの表を$p$、裏を$1-p$のようにパラメータ$p$で考えたが、多項分布でも同様に下記のように$K$次元パラメータベクトル$\mu$を定義する。
$$
\large
\begin{align}
\mu &= \left(\begin{array}{c} \mu_{1} \\ \vdots \\ \mu_{K} \end{array} \right) \quad (1) \\
\mu_{k} & \geq 0 \quad (2) \\
\sum_{k=1}^{K} \mu_{k} &= 1 \quad (3)
\end{align}
$$

$(1)$式のように$K$次元ベクトルの形式でパラメータベクトルの$\mu$を定義することができるが、確率であることより$(2)$式のように$\mu$の各要素は$0$以上であるかつ$(3)$式のように和が$1$にならなければならない。

二項分布(Binomial distribution)と対応させて考えると、$(1)$式は$\displaystyle \left(\begin{array}{c} p \\ 1-p \end{array} \right)$のベクトルに対応し、$1-p$を用いることで$(3)$式も同時に内包することも確認できる。また、$(2)$式に関しては$0 \leq p \leq 1$を設定することでこれも成立する。このように多項分布は二項分布の拡張であると考えることができる。

さて、上記のように定義した$\mu$に基づいて多項分布$\mathrm{Multi}(m_1,…,m_K|\mu,N)$の確率関数を考える。確率変数のベクトルを$M$とおき、$N$回の試行に対して$M$に対応する観測回数のベクトルを$m$とおく。ここでベクトル$m$は下記を表すと考える。
$$
\large
\begin{align}
m &= \left(\begin{array}{c} m_{1} \\ \vdots \\ m_{K} \end{array} \right)
\end{align}
$$

上記は式だけを考えるとわかりにくいが、「$K$面体のサイコロを$N$回投げた時のそれぞれの出目が観測される回数」と対応づけて理解すればそれほど難しくないことがわかる。この多項分布に対し。確率関数$p(m_1,…,m_K|\mu,N)$は下記のように表すことができる。
$$
\large
\begin{align}
p(m_1,…,m_K|\mu,N) &= \left(\begin{array}{c} N \\ m_{1} m_{2} … m_{K} \end{array} \right) \prod_{k=1}^{K} \mu_{k}^{m_{k}} \\
\left(\begin{array}{c} N \\ m_{1} m_{2} … m_{K} \end{array} \right) &= \frac{N!}{m_{1}! m_{2}! … m_{K}!} \\
\sum_{k=1}^{K} m_{k} &= N
\end{align}
$$

多項分布と尤度関数

観測値$\mathcal{D} = (x_1,…,x_N) \sim \mathrm{Multi}(m_1,…,m_K|\mu,N) \quad \mathrm{i.i.d}$に対して尤度関数$L(m_1,…,m_K)$は下記のように考えることができる。
$$
\large
\begin{align}
L(m_1,…,m_K) &= p(x_1,…,x_N|\mu,N) \\
&= p(m_1,…,m_K|\mu,N) \\
&= \left(\begin{array}{c} N \\ m_{1} m_{2} … m_{K} \end{array} \right) \prod_{k=1}^{K} \mu_{k}^{m_{k}} \\
m_{k} &= \sum_{n=1}^{N} {x_{nk}} \\
x_{i} &= \left(\begin{array}{c} x_{i1} \\ \vdots \\ x_{iK} \end{array} \right)
\end{align}
$$

ここで$x_{i}$は$1$hotベクトルである。

ディリクレ分布

多項分布の尤度関数の関数形への着目

共役分布の見つけ方」でも取り扱ったが、「尤度関数をパラメータの関数と見たときに同じ関数形の確率分布を考える」ことで共役事前分布を導出することができる。

多項分布の尤度関数をパラメータベクトル$\mu$に関して確認した際に、$\mu$の事前分布に下記のような関数形を考えることができる。
$$
\large
\begin{align}
p(\mu|\alpha) & \propto \prod_{k=1}^{K} \mu_{k}^{\alpha_{k}-1} \\
\alpha &= \left(\begin{array}{c} \alpha_{1} \\ \vdots \\ \alpha_{K} \end{array} \right)
\end{align}
$$

$\alpha$は事前分布のパラメータを表すベクトルであると考えれば良い。上記を確率分布に対応させるにあたって正規化を行うと下記のようなディリクレ分布$\mathrm{Dir}(\mu|\alpha)$の確率密度関数$p(\mu|\alpha)$が得られる。
$$
\large
\begin{align}
p(\mu|\alpha) = \frac{\Gamma(\alpha_{1}+…+\alpha_{K})}{\Gamma(\alpha_{1}) … \Gamma(\alpha_K)} \prod_{k=1}^{K} \mu_{k}^{\alpha_{k}-1}
\end{align}
$$

ディリクレ分布の正規化定数の導出

多項分布の事後分布

多項分布のパラメータベクトル$\mu$に関する事後分布を$p(\mu|\mathcal{D},\alpha)$のようにおくとき、ベイズの定理に基づいて事後分布は下記のように考えることができる。
$$
\large
\begin{align}
p(\mu|\mathcal{D},\alpha) & \propto p(\mathcal{D}|\mu) p(\mu|\alpha) \\
& \propto \prod_{k=1}^{K} \mu_{k}^{\alpha_{k}+m_{k}-1}
\end{align}
$$

参考

・共役事前分布
・統計学実践ワークブック ベイズ法 演習
・(YouTube)【ゲーム×統計@桃鉄】# 04 勝率を推論しよう(ベイズ推論を使う)

「統計検定1級テキスト」 練習問題解答例 Ch.2 「種々の確率分布」

当記事は「統計検定$1$級対応 統計学(東京図書)」の読解サポートにあたって第$2$章の「種々の確率分布」に関して演習問題を中心に解説を行います。正規分布や指数分布、$\chi^2$分布のような基本的な分布からガンマ分布やベータ分布など様々な分布について取り扱われているので、抑えておくと良いと思います。

本章のまとめ

練習問題解説

問$2$.$1$

問$2$.$2$

ポアソン分布$\mathrm{Poisson}(\lambda)$を$\lambda$に関してガンマ分布$\mathrm{Ga}(\alpha,\beta)$で混合した分布を$f(x)$とおくと、$f(x)$は下記のように考えることができる。
$$
\large
\begin{align}
f(x) &= \int_{0}^{\infty} \exp(-\lambda) \frac{\lambda^{x}}{x!} \times \frac{\beta^{\alpha}}{\Gamma(\alpha)} \lambda^{\alpha-1} \exp(-\beta \lambda) d \lambda \\
&= \frac{\beta^{\alpha}}{x! \Gamma(\alpha)} \int_{0}^{\infty} \lambda^{x+\alpha-1} \exp(-(\beta+1)x) d \lambda \quad (1)
\end{align}
$$

$(1)$式の積分はガンマ分布に関する式を元に考えると良い。
$$
\large
\begin{align}
\int_{0}^{\infty} x^{\alpha} \exp(-\beta x) dx = \frac{\Gamma(\alpha)}{\beta^{\alpha}} \quad (2)
\end{align}
$$

$(1)$式に対して$(2)$式を適用すると、下記のように変形を行うことができる。
$$
\large
\begin{align}
f(x) &= \frac{\beta^{\alpha}}{x! \Gamma(\alpha)} \int_{0}^{\infty} \lambda^{x+\alpha-1} \exp(-(\beta+1)x) d \lambda \\
&= \frac{\beta^{\alpha}}{x! \Gamma(\alpha)} \times \frac{\Gamma(x + \alpha)}{(\beta+1)^{x+\alpha}} \\
&= \frac{\Gamma(x + \alpha)}{x! \Gamma(\alpha)}\frac{\beta^{\alpha}}{(\beta+1)^{x+\alpha}} \quad (3)
\end{align}
$$

ここで$(3)$式に対し、$\displaystyle \alpha = r, \beta = \frac{p}{1-p}$の代入を行う。
$$
\large
\begin{align}
f(x) &= \frac{\Gamma(x + \alpha)}{x! \Gamma(\alpha)}\frac{\beta^{\alpha}}{(\beta+1)^{x+\alpha}} \\
&= \frac{\Gamma(x + r)}{x! \Gamma(r)} \times \left( \frac{p}{1-p} \right)^{r} \times (1-p)^{x+r} \\
&= \frac{(x+r-1)!}{x! (r-1)!} p^{r}(1-p)^{x}
\end{align}
$$

上記が負の二項分布の確率関数に一致する。なお、途中計算では$\displaystyle \frac{1}{\beta+1}$に関して下記が成立することを用いた。
$$
\large
\begin{align}
\frac{1}{\beta+1} &= \frac{1}{\displaystyle \frac{p}{1-p}+1} \\
&= \frac{1}{\displaystyle \frac{p+(1-p)}{1-p}} \\
&= \frac{1}{\displaystyle \frac{1}{1-p}} = 1-p
\end{align}
$$

問$2$.$3$

$Y \sim U(0,1)$より、$Y$の累積分布関数に関して$P(Y \leq y) = y$が成立する。また、$X = F^{-1}(Y)$かつ$F$が狭義単調増加であることより下記が成立する。
$$
\large
\begin{align}
X &= F^{-1}(Y) \\
Y &= F(X)
\end{align}
$$

このとき$X$に関する累積分布関数$P(X \leq x)$を考えると、$P(X \leq x)$は下記のように変形できる。
$$
\large
\begin{align}
P(X \leq x) &= P(F^{-1}(Y) \leq x) \\
&= P(Y \leq F(x)) = F(x)
\end{align}
$$

上記より、$F$が$X$の累積分布関数を表すことが確認できる。

・考察
「$X = F^{-1}(Y)$」かつ$「F$が狭義単調増加」は、「$X$を$F$によって$0$〜$1$の値に対応させる」かつ「$X$が大きくなるにつれて$F(X)$も大きくなる」に対応すると考えることができるので、これは累積分布関数に意味すると直感的に解釈することができる。

問$2$.$4$

$$
\large
\begin{align}
X_1 &= \frac{Y_1}{Y_1 + Y_2 + Y_3} \\
X_2 &= \frac{Y_2}{Y_1 + Y_2 + Y_3} \\
X_3 &= Y_1 + Y_2 + Y_3
\end{align}
$$

$Y_1, Y_2, Y_3$に対応する同時確率密度関数$g(y_1,y_2,y_3)$に対し、上記に基づいて変数変換を行うことを考える。同時確率密度関数$g(y_1,y_2,y_3)$は下記のように表せる。
$$
\large
\begin{align}
g(y_1,y_2,y_3) &= g(y_1)g(y_2)g(y_3) \\
&= \prod_{i=1}^{3} \frac{\beta^{\alpha_i}}{\Gamma(\alpha_i)} y_i^{\alpha_i-1} \exp(-\beta y_i) \quad (1)
\end{align}
$$

このとき$x_3=y_1+y_2+y_3$に着目することで、$y_1, y_2, y_3$は$x_1, x_2, x_3$を用いて下記のように表すことができる。
$$
\large
\begin{align}
y_1 &= x_1 x_3 \quad (2) \\
y_2 &= x_2 x_3 \quad (3) \\
y_3 &= x_3(1-x_1-x_2) \quad (4)
\end{align}
$$

ここで下記のようにヤコビ行列$\displaystyle J = \left( \frac{\partial y_i}{\partial x_j} \right)$を考える。
$$
\large
\begin{align}
J &= \left( \frac{\partial y_i}{\partial x_j} \right) \\
&= \left(\begin{array}{ccc} \frac{\partial y_1}{\partial x_1} & \frac{\partial y_1}{\partial x_2} & \frac{\partial y_1}{\partial x_3} \\ \frac{\partial y_2}{\partial x_1} & \frac{\partial y_2}{\partial x_2} & \frac{\partial y_2}{\partial x_3} \\ \frac{\partial y_3}{\partial x_1} & \frac{\partial y_3}{\partial x_2} & \frac{\partial y_3}{\partial x_3} \end{array} \right) \\
&= \left(\begin{array}{ccc} x_3 & 0 & x_1 \\ 0 & x_3 & x_2 \\ -x_3 & -x_3 & 1-x_1-x_2 \end{array} \right)
\end{align}
$$

上記よりヤコビアン$|\det J| = x_3^2(1-x_1-x_2) – (x_3^2x_1 + x_3^2x_2) = x_3^2$が計算できる。

$(1)$式に対し$(2)$〜$(4)$式に基づいて変数変換を行うことを考える。$X_1,X_2,X_3$に関する確率密度関数を$f(x_1,x_2,x_3)$とおくと、$f(x_1,x_2,x_3)$は下記のように考えることができる。
$$
\large
\begin{align}
f(x_1,x_2,x_3) = & g(y_1,y_2,y_3) |\det J| \\
= & \frac{\beta^{\alpha_1}}{\Gamma(\alpha_1)} (x_1 x_3)^{\alpha_1-1} \exp(-\beta x_1 x_3) \times \frac{\beta^{\alpha_2}}{\Gamma(\alpha_2)} (x_2 x_3)^{\alpha_2-1} \exp(-\beta x_2 x_3) \\
& \times \frac{\beta^{\alpha_3}}{\Gamma(\alpha_3)} (x_3(1-x_1-x_2))^{\alpha_3-1} \exp(-\beta x_3(1-x_1-x_2)) \times x_3^2 \\
= & \frac{\beta^{\alpha_1+\alpha_2+\alpha_3}}{\Gamma(\alpha_1)\Gamma(\alpha_2)\Gamma(\alpha_3)} x_1^{\alpha_1-1} x_2^{\alpha_2-1} (1-x_1-x_2)^{\alpha_3-1} x_3^{\alpha_1+\alpha_2+\alpha_3-1} \exp(-\beta x_3)
\end{align}
$$

ここで$(5)$式を$x_3$に関して$0$から$\infty$まで積分することを考えると、$f(x_1,x_2)$に対応する。
$$
\large
\begin{align}
f(x_1,x_2) &= \int_{0}^{\infty} f(x_1,x_2,x_3) d x_3 \\
&= \int_{0}^{\infty} \frac{\beta^{\alpha_1+\alpha_2+\alpha_3}}{\Gamma(\alpha_1)\Gamma(\alpha_2)\Gamma(\alpha_3)} x_1^{\alpha_1-1} x_2^{\alpha_2-1} (1-x_1-x_2)^{\alpha_3-1} x_3^{\alpha_1+\alpha_2+\alpha_3-1} \exp(-\beta x_3) d x_3 \\
&= \frac{\beta^{\alpha_1+\alpha_2+\alpha_3}}{\Gamma(\alpha_1)\Gamma(\alpha_2)\Gamma(\alpha_3)} x_1^{\alpha_1-1} x_2^{\alpha_2-1} (1-x_1-x_2)^{\alpha_3-1} \int_{0}^{\infty} x_3^{\alpha_1+\alpha_2+\alpha_3-1} \exp(-\beta x_3) d x_3
\end{align}
$$

上記に対して$\displaystyle \int_{0}^{\infty} x^{\alpha} \exp(-\beta x) dx = \frac{\Gamma(\alpha)}{\beta^{\alpha}}$を適用することを考える。
$$
\large
\begin{align}
f(x_1,x_2) &= \frac{\beta^{\alpha_1+\alpha_2+\alpha_3}}{\Gamma(\alpha_1)\Gamma(\alpha_2)\Gamma(\alpha_3)} x_1^{\alpha_1-1} x_2^{\alpha_2-1} (1-x_1-x_2)^{\alpha_3-1} \int_{0}^{\infty} x_3^{\alpha_1+\alpha_2+\alpha_3-1} \exp(-\beta x_3) d x_3 \\
&= \frac{\beta^{\alpha_1+\alpha_2+\alpha_3}}{\Gamma(\alpha_1)\Gamma(\alpha_2)\Gamma(\alpha_3)} x_1^{\alpha_1-1} x_2^{\alpha_2-1} (1-x_1-x_2)^{\alpha_3-1} \times \frac{\Gamma(\alpha_1+\alpha_2+\alpha_3)}{\beta^{\alpha_1+\alpha_2+\alpha_3}} \\
&= \frac{\Gamma(\alpha_1+\alpha_2+\alpha_3)}{\Gamma(\alpha_1)\Gamma(\alpha_2)\Gamma(\alpha_3)} x_1^{\alpha_1-1} x_2^{\alpha_2-1} (1-x_1-x_2)^{\alpha_3-1}
\end{align}
$$

上記は、$\displaystyle f(x_1,x_2) = \frac{\Gamma(\alpha_1+\alpha_2+\alpha_3)}{\Gamma(\alpha_1)\Gamma(\alpha_2)\Gamma(\alpha_3)} x_1^{\alpha_1-1} x_2^{\alpha_2-1} (1-x_1-x_2)^{\alpha_3-1}$が成立することを示す。

・参考
$n \times n$正方行列の行列式
https://www.hello-statisticians.com/explain-books-cat/matrix_determinants1.html

問$2$.$5$

問$2$.$6$

確率変数$X$が正規分布$\displaystyle \mathcal{N} \left( 0,\frac{1}{\theta} \right)$に従うとき、この確率密度関数は下記のように表すことができる。
$$
\large
\begin{align}
\sqrt{\frac{\theta}{2 \pi}} \exp \left( – \frac{x^2 \theta}{2} \right)
\end{align}
$$

上記の分布をパラメータの$\theta$に関してガンマ分布$\mathrm{Ga}(\alpha,\alpha)$を元に混合した際の確率密度関数を$f(x)$とおくと、$f(x)$は下記のように考えることができる。
$$
\large
\begin{align}
f(x) &= \int_{0}^{\theta} \sqrt{\frac{\theta}{2 \pi}} \exp \left( – \frac{x^2 \theta}{2} \right) \times \frac{\alpha^{\alpha}}{\Gamma(\alpha)} \theta^{\alpha-1} \exp{(-\alpha \theta)} d \theta \\
&= \frac{\alpha^{\alpha}}{\sqrt{2 \pi} \Gamma(\alpha)} \int_{0}^{\infty} \theta^{\alpha + \frac{1}{2} – 1} \exp \left[ – \left( \alpha + \frac{x^2}{2} \right) \theta \right] d \theta \quad (1)
\end{align}
$$

$(1)$式に対して下記の式を適用することを考える。
$$
\large
\begin{align}
\int_{0}^{\infty} x^{\alpha} \exp(-\beta x) dx &= \frac{\Gamma(\alpha)}{\beta^{\alpha}} \quad (2) \\
B(a,b) &= \frac{\Gamma(a)\Gamma(b)}{\Gamma(a+b)} \quad (3) \\
\Gamma \left( \frac{1}{2} \right) &= \sqrt{\pi} \quad (4)
\end{align}
$$

$(1)$式に対して$(2), (3), (4)$式を適用すると、下記のように変形を行うことができる。
$$
\large
\begin{align}
f(x) &= \frac{\alpha^{\alpha}}{\sqrt{2 \pi} \Gamma(\alpha)} \int_{0}^{\infty} \theta^{\alpha + \frac{1}{2} – 1} \exp \left[ – \left( \alpha + \frac{x^2}{2} \right) \theta \right] d \theta \\
&= \frac{\alpha^{\alpha}}{\sqrt{2 \pi} \Gamma(\alpha)} \times \frac{\displaystyle \Gamma \left( \alpha + \frac{1}{2} \right)}{\displaystyle \left( \alpha + \frac{x^2}{2} \right)^{\alpha + \frac{1}{2}}} \\
&= \frac{\displaystyle \Gamma \left( \alpha + \frac{1}{2} \right)}{\sqrt{2 \pi} \Gamma(\alpha) \alpha^{\frac{1}{2}}} \times \frac{\alpha^{\alpha + \frac{1}{2}}}{\displaystyle \left( \alpha + \frac{x^2}{2} \right)^{\alpha + \frac{1}{2}}} \\
&= \frac{\displaystyle \Gamma \left( \frac{1}{2} \right)}{\sqrt{2 \pi} B \left( \alpha,\frac{1}{2} \right) \alpha^{\frac{1}{2}}} \times \frac{1}{\displaystyle \left( 1 + \frac{x^2}{2 \alpha} \right)^{\alpha + \frac{1}{2}}} \\
&= \frac{\sqrt{\pi}}{\sqrt{2 \pi} B \left( \alpha,\frac{1}{2} \right) \alpha^{\frac{1}{2}}} \times \frac{1}{\displaystyle \left( 1 + \frac{x^2}{2 \alpha} \right)^{\alpha + \frac{1}{2}}} \\
&= \frac{1}{\sqrt{2} B \left( \alpha,\frac{1}{2} \right) \alpha^{\frac{1}{2}}} \times \left( 1 + \frac{x^2}{2 \alpha} \right)^{-\left( \alpha + \frac{1}{2} \right)}
\end{align}
$$

ここで$(5)$式に$\displaystyle \alpha = \frac{p}{2}$を代入すると下記のように変形できる。
$$
\large
\begin{align}
f(x) &= \frac{1}{\sqrt{2} B \left( \alpha,\frac{1}{2} \right) \alpha^{\frac{1}{2}}} \times \left( 1 + \frac{x^2}{2 \alpha} \right)^{-\left( \alpha + \frac{1}{2} \right)} \\
&= \frac{1}{\sqrt{2} B \left( \frac{p}{2},\frac{1}{2} \right) \left( \frac{p}{2} \right)^{\frac{1}{2}}} \times \left( 1 + \frac{x^2}{p} \right)^{-\left( \frac{p}{2} + \frac{1}{2} \right)} \\
&= \frac{1}{\sqrt{p} B \left( \frac{p}{2},\frac{1}{2} \right)} \left( 1 + \frac{x^2}{p} \right)^{-\frac{p+1}{2}}
\end{align}
$$

上記は自由度$p$の$t$分布の確率密度関数の式に一致する。

参考

・統計検定1級 統計数理 関連まとめ
https://www.hello-statisticians.com/stat_certifi_1_math

オッズ(odds)の概要と一般的な使用例・統計学における使用例

ロジスティック回帰などに関連してオッズ(odds)という考え方が用いられますが、教科書などの解説ではオッズそのものの取り扱いがあまり多くありません。そこで当記事ではオッズという考え方を中心に取り扱い、一般的な使用例も確認しながら具体的に理解できるように取りまとめを行いました。

作成にあたってはWikipediaのオッズのまとめや「統計学実践ワークブック」の$18$章「質的回帰」、$28$章「分割表」などを参考にしました。
・オッズ Wikipedia

オッズの概要

オッズの定義

ある事象が起こる確率を$p$、起こらない確率を$1-p$とするとき、下記のようにオッズ$\mathrm{Odds}$を定義する。
$$
\large
\begin{align}
\mathrm{Odds} = \frac{p}{1-p}
\end{align}
$$

定義域の$0 \leq p < 1$において上記は$p$に関する単調増加関数であり、式定義より$0 \leq \mathrm{Odds} < \infty$が成立する。式だけだと少々わかりにくいので、下記のPythonを実行することでオッズのグラフを描くことができる。

import numpy as np
import matplotlib.pyplot as plt

p = np.arange(0.,1.,0.05)
odds = p/(1-p)

plt.plot(p,odds)
plt.show()

・実行結果

オッズの一般的な使用例

オッズは競馬などの賭けを行う際の儲けの計算に伝統的に用いられてきた。たとえば$5$回に$1$回に$A$、$5$回に$4$回に$B$が起こる場合を考える。このとき$A$が起こるオッズの$\mathrm{Odds}_{A}$と$B$が起こるオッズ$\mathrm{Odds}_{B}$はそれぞれ下記のように計算できる。
$$
\large
\begin{align}
\mathrm{Odds}_{A} &= \frac{0.2}{1-0.2} = 0.25 \\
\mathrm{Odds}_{B} &= \frac{0.8}{1-0.8} = 4
\end{align}
$$

上記のようにオッズを計算する場合、オッズの逆数$+1$が取り分の倍率となる。

一方で競馬などでのオッズは「払戻金」を「的中券$100$円分に対する金額」で表現される。これをパリミュチュエル方式という。

$\displaystyle \mathrm{Odds} = \frac{p}{1-p}$で定義されるオッズは的中時の倍率を重視する一般的な用法と少々異なるので注意が必要であると考えられる。

・参考
払戻金 Wikipedia

統計学におけるオッズの取り扱い

ロジスティック回帰とオッズ

「ロジスティック回帰」は「事象が生じる確率を予測する回帰に対応する」と解釈することができる。事象を表す確率変数を$Y$と定義し、$y_1,…,y_n$が観測されると考える。
$$
\large
\begin{align}
\pi_{i} &= \hat{y}_i = \frac{1}{1+\exp{(f(x_i))}} \quad (1) \\
f(x_i) &= \beta_0 + \beta_1 x_{i1} + … + \beta_{k} x_{ik}
\end{align}
$$

ここで上記の$(1)$式を$\exp{(f(x_i))}$に関して解くことを考える。
$$
\large
\begin{align}
\pi_{i} &= \frac{1}{1+\exp{(-f(x_i))}} \\
\pi_{i} &= \frac{\exp{(f(x_i))}}{\exp{(f(x_i))}+1} \\
\pi_{i}(1+\exp{(f(x_i))}) &= \exp{(f(x_i))} \\
\pi_{i} + \pi_{i} \exp{(f(x_i))} &= \exp{(f(x_i))} \\
\pi_{i} &= (1-\pi_{i}) \exp{(f(x_i))} \\
\frac{\pi_{i}}{1-\pi_{i}} &= \exp{(f(x_i))} \quad (2)
\end{align}
$$

$(2)$式より$\exp{(f(x_i))} = \exp{(\beta_0 + \beta_1 x_{i1} + … + \beta_{k} x_{ik})} = \exp{\beta_0} \exp{\beta_1 x_{i1}}…\exp{\beta_{k} x_{ik})}$は確率のオッズに対応することが読み取れる。

また、$\exp{(f(x_i))}$は下記のように表すことができる。
$$
\large
\begin{align}
\exp{(f(x_i))} &= \exp{(\beta_0 + \beta_1 x_{i1} + … + \beta_{k} x_{ik})} \\
&= e^{\beta_0} e^{\beta_1 x_{i1}}…e^{\beta_{k} x_{ik}} \quad (3)
\end{align}
$$

$(2), (3)$式より$\exp{(f(x_i))} = e^{\beta_0} e^{\beta_1 x_{i1}}…e^{\beta_{k} x_{ik}}$は「確率変数$Y$に関する事象が観測される確率」のオッズに対応することが読み取れる。

分割表とオッズ比

分割表における$2$群の一様性などを考えるにあたって、オッズの比率のオッズ比を考える。

参考

・質的回帰 演習 ロジスティック回帰とオッズ
・オッズ Wikipedia

数列の和の公式とその導出 〜等差数列の和、等比数列の和、Σの公式 etc〜

数列の和の計算は基本的な数学のトピックである一方で、単に公式を抑えるだけではわからなくなりがちです。そこで当記事では数列の和の公式とその導出に関して取り扱いました。
等差数列・等比数列の和の公式の導出は考え方の理解が中心である一方で、$\displaystyle \sum$に関連する和の導出は二項定理を元に考えるなど少々トリッキーですが、簡単な流れだけでも抑えておくと良いと思います。

・数学まとめ
https://www.hello-statisticians.com/math_basic

等差数列・等比数列の和

等差数列の和

$1,3,5,…$のように、値が同じ数だけ変化する数列を等差数列という。

ここで$n$番目の項を$a_{n}$のように考えたとき、等差数列では$a_1$を初項、$a_{n+1}-a_{n}=d$を公差とそれぞれ考える。$a_1$ではなく$a_0$を初項と考える場合もあるが、これは定義次第なので状況に合わせて考えれば良い。

数列の表記にあたっては「要素を書き下す定義」、「一般項を用いた定義」、「漸化式を用いた定義」の主に$3$通りがあるが、「要素を書き下す定義」が$1,3,5,…$であるとき、「一般項を用いた定義」は$a_{n} = 1 + 2(n-1)$、「漸化式を用いた定義」は$a_{n+1}=a_{n}+2, a_{1}=1$のようにそれぞれ表すことができる。特に「漸化式」が難しいように見えがちだが、「一般項」と同様に単に数列の表記を取り扱うと考えればそれほど難しくない。

ここまでで等差数列の概要に関して確認ができたので、$a_{n} = a + d(n-1)$で表される数列の和$a_{1}+a_{2}+…+a_{n}$を計算することを考える。これを考えるにあたっては$\displaystyle \frac{a_{1}+a_{n}}{2}$が$a_{1},a_{2},…,a_{n}$の平均であることに基づいて考えることでシンプルに考えることができる。$\displaystyle \frac{a_{1}+a_{n}}{2}$が$a_{1},a_{2},…,a_{n}$の平均であるので、$\displaystyle \frac{a_{1}+a_{n}}{2}$に個数をかけることで下記のように等差数列の和を計算することができる。
$$
\large
\begin{align}
a_{1} + a_{2} + … + a_{n} &= \frac{a_{1}+a_{n}}{2} \times n \\
&= \frac{a + a + d(n-1)}{2} \times n \\
&= \frac{n(2a + d(n-1))}{2}
\end{align}
$$

上記が初項$a$、公差$d$の等差数列の和の公式である。

等比数列の和

$1,2,4,8,16,…$のように、次の数が$1$つ前の数の定数倍になる数列を等比数列という。初項$a$、定数倍に対応する公比が$r$の等比数列は$a, ar, ar^2,…$のように表される。

また、「一般項を用いた等比数列の定義」は$a_{n} = ar^{n-1}$、「漸化式を用いた等比数列の定義」は$a_{n+1}=ra_{n}, a_{1}=a$のように表される。以下、この等比数列の和の導出を考える。

等比数列の和の導出にあたっては$S = a_{1}+…+a_{n}$とおき、$S-rS = (1-r)S$を考えればよい。詳細の計算を下記に示す。
$$
\large
\begin{align}
S-rS &= (a_{1}+…+a_{n}) – r(a_{1}+…+a_{n}) \\
&= (a+ar+ar^2+…+ar^{n-1}) – r(a+ar+ar^2+…+ar^{n-1}) \\
&= (a+ar+ar^2+…+ar^{n-1}) – (ar+ar^2+…+ar^{n}) \\
&= a + (ar-ar) + (ar^2-ar^2) + … + (ar^{n-1}-ar^{n-1}) – ar^n \\
&= a(1-r^n)
\end{align}
$$

上記より、下記のように等比数列の公式が導出できる。
$$
\large
\begin{align}
S-rS &= (1-r)S = a(1-r^n) \\
S &= \frac{a(1-r^n)}{1-r}
\end{align}
$$

$\displaystyle \sum$に関する公式とその導出

$\displaystyle \sum_{k=1}^{n} k = \frac{1}{2}n(n+1)$

$\displaystyle \sum$は数列の和を表すので、$\displaystyle \sum_{k=1}^{n} k = 1+2+…+n$のように考えればよい。また、要素の和の形式で表すことで初項$1$、公差$1$の等差数列であることも確認できる。これより、$\displaystyle \sum_{k=1}^{n} k$は下記のように計算できる。
$$
\large
\begin{align}
\sum_{k=1}^{n} k &= 1 + 2 + … + n \\
&= \frac{1+n}{2} \times n \\
&= \frac{1}{2}n(n+1)
\end{align}
$$

上記は$\displaystyle \sum_{k=1}^{n} k = \frac{1}{2}n(n+1)$が成立することを示す。

$\displaystyle \sum_{k=1}^{n} k^2 = \frac{1}{6}n(n+1)(2n+1)$

二項定理より、$(k+1)^3=k^3+3k^2+3k+1$が成立することを元に、$\displaystyle \sum_{k=1}^{n} k^2 = \frac{1}{6}n(n+1)(2n+1)$の導出を行う。

まず$(k+1)^3=k^3+3k^2+3k+1$は下記のように変形することができる。
$$
\large
\begin{align}
(k+1)^3 &= k^3 + 3k^2 + 3k + 1 \\
(k+1)^3 – k^3 &= 3k^2 + 3k + 1
\end{align}
$$

上記は任意の$k$に関して成立するので、$k=1$から$k=n$の和を計算しても同様に成立する。
$$
\large
\begin{align}
\sum_{i=1}^{n} ((k+1)^3 – k^3) &= \sum_{i=1}^{n} (3k^2 + 3k + 1) \\
(n+1)^3 – 1 &= \sum_{i=1}^{n} (3k^2 + 3k + 1) \\
3 \sum_{i=1}^{n} k^2 &= (n^3+3n^2+3n+1-1) – \sum_{i=1}^{n} (3k + 1) \\
&= n(n^2+3n+3) – \frac{3n(n+1)}{2} – n \\
&= \frac{n}{2}(2(n^2+3n+3) – (3n+3) – 2) \\
&= \frac{n}{2}(2n^2+3n+1) \\
&= \frac{1}{2}n(n+1)(2n+1)
\end{align}
$$

上記より下記の公式が成立することが示せる。
$$
\large
\begin{align}
3 \sum_{i=1}^{n} k^2 &= \frac{1}{2}n(n+1)(2n+1) \\
\sum_{i=1}^{n} k^2 &= \frac{1}{6}n(n+1)(2n+1)
\end{align}
$$

$\displaystyle \sum_{k=1}^{n} k^3 = \left[ \frac{1}{2}n(n+1) \right]^3$

二項定理より、$(k+1)^4=k^4+4k^3+6k^2+4k+1$が成立することを元に、$\displaystyle \sum_{k=1}^{n} k^3 = \left[ \frac{1}{2}n(n+1) \right]^3$の導出を行う。

まず$(k+1)^4=k^4+4k^3+6k^2+4k+1$は下記のように変形することができる。
$$
\large
\begin{align}
(k+1)^4 &= k^4 + 4k^3 + 6k^2 + 4k + 1 \\
(k+1)^4 – k^4 &= 4k^3 + 6k^2 + 4k + 1
\end{align}
$$

上記は任意の$k$に関して成立するので、$k=1$から$k=n$の和を計算しても同様に成立する。
$$
\large
\begin{align}
\sum_{i=1}^{n} ((k+1)^4 – k^4) &= \sum_{i=1}^{n} (4k^3 + 6k^2 + 4k + 1) \\
(n+1)^4 – 1 &= \sum_{i=1}^{n} (4k^3 + 6k^2 + 4k + 1) \\
4 \sum_{i=1}^{n} k^3 &= (n^4+4n^3+6n^2+4n+1-1) – \sum_{i=1}^{n} (6k^2 + 4k + 1) \\
&= n(n^3+4n^2+6n+4) – n(n+1)(2n+1) – 2n(n+1) – n \\
&= n((n^3+4n^2+6n+4) – (n+1)(2n+1) – 2(n+1) – 1) \\
&= n(n^3 + (4-2)n^2 + (6-3-2)n 4-1-2-1) \\
&= n(n^3 + 2n^2 + n) \\
&= n^2(n+1)^2
\end{align}
$$

上記より下記の公式が成立することが示せる。
$$
\large
\begin{align}
4 \sum_{i=1}^{n} k^3 &= n^2(n+1)^2 \\
\sum_{i=1}^{n} k^3 &= \left[ \frac{1}{2}n(n+1) \right]^3
\end{align}
$$

参考

・スピアマンの順位相関係数の導出
https://www.hello-statisticians.com/explain-terms-cat/spearman_coef1.html

スピアマンの順位相関係数(Spearman correlation coefficient)の導出

統計検定の準1級ワークブックなどに出てくるスピアマンの順位相関係数(Spearman correlation coefficient)は通常の相関係数の式から導出できるとされる一方で詳しい導出がないので、当記事ではスピアマンの順位相関係数の式の導出に関して取り扱いました。

スピアマンの順位相関係数の概要

ピアソンの積率相関係数

ピアソンの積率相関係数は一般的に用いられる相関係数であり、標本$(x_1,y_1), …, (x_n,y_n)$に対して相関係数$r$は下記のように定義される。
$$
\large
\begin{align}
r = \frac{\displaystyle \sum_{i=1}^{n}(x_i-\bar{x})(y_i-\bar{y})}{\displaystyle \sqrt{\sum_{i=1}^{n}(x_i-\bar{x})^2\sum_{i=1}^{n}(x_i-\bar{x})^2}}
\end{align}
$$

また、標本$(x_i,y_i)$に対応する確率変数を$X, Y$のようにおくとき、$r$は下記のように定義することもできる。
$$
\large
\begin{align}
r = \frac{\mathrm{Cov}[X,Y]}{\sqrt{V[X]V[Y]}}
\end{align}
$$

ここで$V[X], V[Y], \mathrm{Cov}[X,Y]$はそれぞれ$X$の分散、$Y$の分散、$X,Y$の共分散を表す。$V[X], V[Y], \mathrm{Cov}[X,Y]$に関連する公式は下記で取りまとめた。

・期待値、分散の公式
https://www.hello-statisticians.com/explain-terms-cat/expectation-variance-covariance.html

スピアマンの順位相関係数の式

前項のように$(x_i,y_i)$を考えるとき、スピアマンの順位相関係数(Spearman correlation coefficient)を$r_s$のようにおくと、$r_s$は下記のような式で表される。
$$
\large
\begin{align}
r = 1 – \frac{\displaystyle 6 \sum_{i=1}^{n} (x_i-y_i)^2}{n(n^2-1)}
\end{align}
$$

スピアマンの順位相関係数の導出

$V[X], \mathrm{Cov}[X,Y]$に関する公式

$V[X], \mathrm{Cov}[X,Y]$に関して下記の式を公式のように用いることができる。
$$
\large
\begin{align}
V[X] &= E[X^2] – E[X]^2 \\
\mathrm{Cov}[X,Y] &= E[XY] – E[X]E[Y]
\end{align}
$$

詳しい導出は下記で取り扱った。
https://www.hello-statisticians.com/explain-terms-cat/expectation-variance-covariance.html

数列の和の公式

$\displaystyle \sum_{k=1}^{n} k, \sum_{k=1}^{n} k^2$に関して下記のような式が成立する。
$$
\large
\begin{align}
\sum_{k=1}^{n} k &= \frac{n(n+1)}{2} \\
\sum_{k=1}^{n} k^2 &= \frac{1}{6}n(n+1)(2n+1)
\end{align}
$$

詳しい導出に関しては下記で取り扱った。
・数列の和の公式とその導出
https://www.hello-statisticians.com/explain-terms-cat/sum_formula1.html

スピアマンの順位相関係数の導出

以下、ピアソンの相関係数の式に対し順位を取り扱うことによる制約を考慮することでスピアマンの順位相関係数の式を導出する。
$$
\large
\begin{align}
r &= \frac{\mathrm{Cov}[X,Y]}{\sqrt{V[X]V[Y]}} \\
&= \frac{E[XY]-E[X]E[Y]}{\sqrt{V[X]V[Y]}} \quad (1)
\end{align}
$$

上記がピアソンの積率相関係数の式だが、$x_i,y_i$がどちらも順序を表すことから期待値に関して$E[X]=E[Y]$や分散に関して$V[X]=V[Y]$が成立する。これより$(1)$式は下記のように表すことができる。
$$
\large
\begin{align}
r &= \frac{E[XY]-E[X]E[Y]}{\sqrt{V[X]V[Y]}} \quad (1) \\
&= \frac{E[XY]-E[X]^2}{\sqrt{V[X]^2}} \\
&= \frac{E[XY]-E[X]^2}{V[X]} \quad (2)
\end{align}
$$

ここで$E[X], V[X]$は$x_1,x_2,…,x_n$では$1$から$n$が$1$回ずつ出現することから下記のように計算することができる。
$$
\large
\begin{align}
E[X] &= \frac{1}{n} \sum_{i=1}^{n} i \\
&= \frac{n(n+1)}{2n} = \frac{n+1}{2} \quad (3) \\
V[X] &= E[X^2] – E[X]^2 \\
&= \frac{1}{n} \sum_{i=1}^{n} i^2 – \left( \frac{n+1}{2} \right)^2 \\
&= \frac{n(n+1)(2n+1)}{6n} – \frac{(n+1)^2}{4} \\
&= \frac{n+1}{12}(4n+2) – \frac{n+1}{12}(3n+3) \\
&= \frac{(n+1)(n-1)}{12} = \frac{n^2-1}{12} \quad (4)
\end{align}
$$

また、$E[XY]$に関して下記のような変形を行うことができる。
$$
\large
\begin{align}
E[XY] &= E[XY] – E[X^2] + E[X^2] \\
&= \frac{1}{2}(2E[XY] – 2E[X^2]) + E[X^2] \\
&= E[X^2] – \frac{1}{2}(E[X^2] – 2E[XY] + E[Y^2]) \\
&= E[X^2] – \frac{1}{2}E[(X-Y)^2] \\
&= E[X^2] – \frac{1}{2n} \sum_{i=1}^{n}(x_i-y_i)^2 \quad (5)
\end{align}
$$

ここで$(3)$〜$(5)$式を$(2)$式を代入することで下記のように変形を行うことができる。
$$
\large
\begin{align}
r &= \frac{E[XY]-E[X]^2}{V[X]} \\
&= \frac{\displaystyle E[X^2] – \frac{1}{2n} \sum_{i=1}^{n}(x_i-y_i)^2 – E[X]^2}{V[X]} \\
&= \frac{E[X^2]-E[X]^2}{V[X]} – \frac{\displaystyle \frac{1}{2n} \sum_{i=1}^{n}(x_i-y_i)^2}{V[X]} \\
&= \frac{V[X]}{V[X]} – \frac{1}{2n} \sum_{i=1}^{n}(x_i-y_i)^2 \times \frac{12}{n^2-1} \\
&= 1 – \frac{\displaystyle 6 \sum_{i=1}^{n}(x_i-y_i)^2}{n(n^2-1)}
\end{align}
$$

参考

・統計学実践ワークブック 例$13$.$5$

正方行列(Square matrix)のトレース(trace)の定義と定義から成立する公式

行列を考える際に「トレース(trace)」は正方行列(Square matrix)の対角成分の和を表す概念です。機械学習や統計学を学ぶにあたってトレースの記号が出てくるときがあるので当記事ではトレース(trace)の定義と定義から成立する公式に関して取りまとめました。
「統計のための行列代数(Matrix Algebra From a Statistician’s Perspective)」のCh.$5$や「パターン認識と機械学習」の「Appendix C」を参考に作成を行いました。

基本事項のまとめ

トレース(trace)の定義

トレース(trace)は正方行列(Square matrix)の対角成分の和で定義され、$n \times n$正方行列$A=(a_{ij})$のトレース$\mathrm{tr}(A)$を下記のように表記する。
$$
\large
\begin{align}
\mathrm{tr}(A) &= a_{11} + a_{22} + \cdots + a_{nn} \\
&= \sum_{i=1}^{n} a_{ii}
\end{align}
$$

また、定数$k$を元に$[k]$を$1 \times 1$正方行列と考えるとき、$\mathrm{tr}([k])=k$のように考える。

単位行列(Identity matrices)$I_{n}$のトレース

$$
\large
\begin{align}
I_{n} = \left(\begin{array}{cccc} 1 & 0 & \cdots & 0 \\ 0 & 1 & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & 1 \end{array} \right)
\end{align}
$$

上記のように$n \times n$単位行列$I_{n}$を考えるとき、$\mathrm{tr}(I_{n})$は下記のように計算できる。
$$
\large
\begin{align}
\mathrm{tr}(I_{n}) &= \sum_{i=1}^{n} 1 \\
&= n
\end{align}
$$

Matrices of ones$J_{n}$のトレース

$$
\large
\begin{align}
J_{n} = \left(\begin{array}{cccc} 1 & 1 & \cdots & 1 \\ 1 & 1 & \cdots & 1 \\ \vdots & \vdots & \ddots & \vdots \\ 1 & 1 & \cdots & 1 \end{array} \right)
\end{align}
$$

上記のように$n \times n$行列$J_{n}$を考えるとき、$\mathrm{tr}(J_{n})$は下記のように計算できる。
$$
\large
\begin{align}
\mathrm{tr}(I_{n}) &= \sum_{i=1}^{n} 1 \\
&= n
\end{align}
$$

トレース(trace)の公式

基本的な公式

$\mathrm{tr}(kA) = k \mathrm{tr}(A)$

定数$k$と$n \times n$正方行列$A=(a_{ij})$のトレースに関して$\mathrm{tr}(kA) = k \mathrm{tr}(A)$が成立することを以下に示す。
$$
\large
\begin{align}
\mathrm{tr}(kA) &= \sum_{i=1}^{n} k a_{ii} \\
&= k \sum_{i=1}^{n} a_{ii} \\
&= k \mathrm{tr}(A)
\end{align}
$$

$\mathrm{tr}(A+B) = \mathrm{tr}(A)+\mathrm{tr}(B)$

$n \times n$正方行列$A=(a_{ij}), B=(b_{ij})$のトレースに関して$\mathrm{tr}(A+B) = \mathrm{tr}(A)+\mathrm{tr}(B)$が成立することを以下に示す。
$$
\large
\begin{align}
\mathrm{tr}(A+B) &= \sum_{i=1}^{n} (a_{ii} + b_{ii}) \\
&= \sum_{i=1}^{n} a_{ii} + \sum_{i=1}^{n} b_{ii} \\
&= \mathrm{tr}(A) + \mathrm{tr}(B)
\end{align}
$$

$\mathrm{tr}(A^{\mathrm{T}}) = \mathrm{tr}(A)$

$n \times n$正方行列$A=(a_{ij})$のトレースに関して$\mathrm{tr}(A^{\mathrm{T}}) = \mathrm{tr}(A)$が成立することを以下に示す。
$$
\large
\begin{align}
\mathrm{tr}(A^{\mathrm{T}}) &= \sum_{i=1}^{n} a_{ii} \\
&= \mathrm{tr}(A)
\end{align}
$$

行列の積のトレースの公式

$\mathrm{tr}(AB) = \mathrm{tr}(B^{\mathrm{T}}A^{\mathrm{T}})$

$\mathrm{tr}(B^{\mathrm{T}}A^{\mathrm{T}}) = \mathrm{tr}((AB)^{\mathrm{T}})$なので、正方行列$X$に対して$\mathrm{tr}(X^{\mathrm{T}}) = \mathrm{tr}(X)$が成立することより$\mathrm{tr}(AB) = \mathrm{tr}(B^{\mathrm{T}}A^{\mathrm{T}})$が成立することがわかる。ここで$m \times n$行列$A$と$n \times m$行列$B$の積である$AB$は$m \times m$正方行列であることは注意しておくとよい。

$\mathrm{tr}(AB) = \mathrm{tr}(BA)$

$m \times n$行列$A=(a_{ij})$と$n \times m$行列$B=(b_{ij})$に関して$\mathrm{tr}(AB) = \mathrm{tr}(B^{\mathrm{T}}A^{\mathrm{T}})$が成立することを以下に示す。

$AB$の$i,j$成分を$(AB)_{ij}$、$BA$の$i,j$成分を$(BA)_{ij}$のように表す。このとき、$(AB)_{ij}, (BA)_{ij}$は下記のように表すことができる。
$$
\large
\begin{align}
(AB)_{ij} &= a_{i1}b_{1j} + a_{i2}b_{2j} + … + a_{in}b_{nj} \\
&= \sum_{k=1}^{n} a_{ik}b_{kj} \\
(BA)_{ij} &= a_{i1}b_{1j} + a_{i2}b_{2j} + … + a_{im}b_{mj} \\
&= \sum_{l=1}^{m} a_{il}b_{kl}
\end{align}
$$

上記より、$(AB)_{ii}, (BA)_{ii}$は下記のように表すことができる。
$$
\large
\begin{align}
(AB)_{ii} &= \sum_{k=1}^{n} a_{ik}b_{ki} \\
(BA)_{ii} &= \sum_{l=1}^{m} a_{il}b_{li}
\end{align}
$$

よって、$\mathrm{tr}(AB)$は下記のように表せる。
$$
\large
\begin{align}
\mathrm{tr}(AB) &= \sum_{l=1}^{m} \sum_{k=1}^{n} a_{lk}b_{kl} \\
&= \sum_{k=1}^{n} \sum_{l=1}^{m} a_{lk}b_{kl} = \mathrm{tr}(BA)
\end{align}
$$

上記より、$\mathrm{tr}(AB) = \mathrm{tr}(B^{\mathrm{T}}A^{\mathrm{T}})$が成立する。

$\mathrm{tr}(ABC) = \mathrm{tr}(CAB) = \mathrm{tr}(BCA)$

$ABC$を$A$と$BC$の積や$AB$と$C$の積のように見ることで、前項の$\mathrm{tr}(AB) = \mathrm{tr}(BA)$を用いて$\mathrm{tr}(ABC) = \mathrm{tr}(CAB) = \mathrm{tr}(BCA)$が成立することを示すことができる。

固有値・固有ベクトルに関する公式

$\displaystyle \mathrm{tr}(A) = \sum_{i=1}^{M} \lambda_{i}$

$A$を固有ベクトルに基づいて列を構成する直交行列$U$と、固有値を対角に並べた対角行列$\Lambda$を用いて下記のように固有値分解を行うことを考える。
$$
\large
\begin{align}
A = U \Lambda U^{\mathrm{T}}
\end{align}
$$

このとき$\mathrm{tr}(A)$は下記のように変形できる。
$$
\large
\begin{align}
\mathrm{tr}(A) &= \mathrm{tr}(U \Lambda U^{\mathrm{T}}) \\
&= \mathrm{tr}(U^{\mathrm{T}} U \Lambda) \\
&= \mathrm{tr}(\Lambda) = \sum_{i=1}^{M} \lambda_{i}
\end{align}
$$

参考

・行列の定義まとめ
https://www.hello-statisticians.com/explain-terms-cat/matrix_def1.html
・統計のための行列代数 Ch.$1$ Matrices

・パターン認識と機械学習 Appendix C

転置行列(transposed matrix)の定義とよく出てくる公式

行列$A$の転置行列$A^{\mathrm{T}}$は様々な場面でよく出てきますが、$(AB)^{\mathrm{T}}=B^{\mathrm{T}}A^{\mathrm{T}}$などは詳しく確認していないとわからなくなりがちです。そこで当記事では転置行列の定義とよく出てくる公式に関して取りまとめを行いました。
統計のための行列代数(Matrix Algebra From a Statistician’s Perspective)」のCh.$1$を参考に作成を行いました。

基本的な考え方

転置行列の定義

「統計のための行列代数(Matrix Algebra From a Statistician’s Perspective)」のSection$1$.$2$.$d$「Transposition」の例を元に$2 \times 3$行列の$A$を下記のように考える。
$$
\large
\begin{align}
A = \left(\begin{array}{ccc} 2 & -3 & 0 \\ 1 & 4 & 2 \end{array} \right)
\end{align}
$$

このとき、$A$の転置行列(transposed matrix)$A^{\mathrm{T}}$は下記のように表される。
$$
\large
\begin{align}
A^{\mathrm{T}} = \left(\begin{array}{cc} 2 & 1 \\ -3 & 4 \\ 0 & 2 \end{array} \right)
\end{align}
$$

このように$A$の転置行列$A^{\mathrm{T}}$は、「$A$の$(i,j)$成分を$(j,i)$成分に持つ行列である」と定義される。上記では$2 \times 3$行列で考えたが、$m \times n$行列でも同様に考えることができる。

行列の成分表示

$m \times n$の行列$A$の$i$行$j$列の要素(element)を$a_{ij}$のように定義し、行列$A$を下記のように表すことを考える。
$$
\large
\begin{align}
A = ( a_{ij} )
\end{align}
$$

上記のような行列の表記を考えることで、行列の式を示す際に「要素全てに対してではなく、$1$つの要素に関してのみ示す」ことで等号の成立を示すことができる。このことに基づくと行列に関する式の導出などがシンプルになり行いやすくなる。

ここで、行列$A$の$(i,j)$成分を考えるにあたって$A = ( a_{ij} )$のようにおき、$a_{ij}$のように取り扱う場合、$A-B+C$の成分を考える場合などに表記が複雑になると思われる。よって、行列$A$の$(i,j)$成分を$(A)_{ij}$のように表す場合がある。

以下、表記の簡略化にあたって、基本的には行列$A$の$(i,j)$成分は$(A)_{ij}$のように表す。また、「統計の森」では表記の簡略化を行うにあたって$(A)_{ij}$の表記をメインで用いる。

抑えておきたい公式とその導出

$(AB)^{\mathrm{T}}=B^{\mathrm{T}}A^{\mathrm{T}}$の導出

行列に関する導出を行うにあたって、$m \times k$行列の$A$と$k \times n$行列の$B$に対して$(AB)^{\mathrm{T}}=B^{\mathrm{T}}A^{\mathrm{T}}$は用いることがかなり多い。よって、単に式表記を抑えるだけでなく、式の意図も合わせて掴んでおく必要がある。

前項の「行列の成分表示」の内容を元に、$(AB)_{ij} = (B^{\mathrm{T}}A^{\mathrm{T}})_{ji}$であることを示せばこの式が成立することが確認できる。
$$
\large
\begin{align}
(AB)_{ij} &= \left(\begin{array}{cccc} a_{i1} & a_{i2} & \cdots & a_{ik} \end{array} \right) \left(\begin{array}{c} b_{1j} \\ b_{2j} \\ \vdots \\ b_{kn} \end{array} \right) \\
&= \sum_{k} a_{ik} b_{kj} \\
(B^{\mathrm{T}}A^{\mathrm{T}})_{ji} &= \left(\begin{array}{cccc} b_{1j} & b_{2j} & \cdots & b_{kj} \end{array} \right) \left(\begin{array}{c} a_{i1} \\ a_{i2} \\ \vdots \\ a_{ik} \end{array} \right) \\
&= \sum_{k} b_{kj} a_{ik} \\
&= \sum_{k} a_{ik} b_{kj}
\end{align}
$$

上記は任意の$(i,j)$成分に関して成立するので、$(AB)^{\mathrm{T}}=B^{\mathrm{T}}A^{\mathrm{T}}$が成立することを示すことができる。

参考

・行列の定義まとめ
https://www.hello-statisticians.com/explain-terms-cat/matrix_def1.html

・統計のための行列代数 Ch.$1$ Matrices

多次元尺度構成法(MDS; Multi-Dimensional Scaling)の導出まとめ

多次元尺度構成法(MDS; Multi-Dimensional Scaling)は個体間の類似度が与えられているときに、類似度に基づいてそれぞれの個体の位置を表現する手法の一つです。当記事では計量MDSの基本的な導出に関して取りまとめを行いました。
「統計学実践ワークブック」の$26$章の「その他の多変量解析手法」の内容を参考に作成を行いました。

 

理論的な導出

問題設定

多次元尺度構成法(MDS; Multi-Dimensional Scaling)は$n$個の個体間の類似度が与えられているときに、類似度に基づいてそれぞれの個体の位置を表現する手法の一つである。

下記の行列$D$を用いて$n$個の個体の類似度が得られている場合を考える。
$$
\large
\begin{align}
D = \left(\begin{array}{ccccc} 0 & d_{12}^2 & d_{13}^2 & \cdots & d_{1n}^2 \\ d_{21}^2 & 0 & d_{23}^2 & \cdots & d_{2n}^2 \\ d_{31}^2 & d_{32}^2 & 0 & \cdots & d_{3n}^2 \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ d_{n1}^2 & d_{n2}^2 & d_{n3}^2 & \cdots & 0 \end{array} \right)
\end{align}
$$

上記のように類似度の行列が与えられたとき、行列を元に個体$i$の位置を$x_{i} = (x_{i1},x_{i2},…,x_{ik})^{\mathrm{T}}$のように表すことを考える。

このとき、下記が成立すると考えても一般性を失わないので、下記が成立すると考える。
$$
\large
\begin{align}
\sum_{i=1}^{n} x_{i} &= \sum_{i=1}^{n} \left(\begin{array}{c} x_{i1} \\ \vdots \\ x_{ik} \end{array} \right) \\
&= \left(\begin{array}{c} 0 \\ \vdots \\ 0 \end{array} \right)
\end{align}
$$

また、表記の簡略化にあたって、$\displaystyle a = \sum_{i}^{n} ||x_i||^{2}$とおく。ここで$||x||^{2}$はフロベニウスノルムであると考える。

・参考
特異値分解 フロベニウスノルム(Frobenius norm)

 

二重中心化(double centering)

・二重中心化の概要
前項の問題設定に基づいて、$x_{i}$を並べた行列$X$に関して対称行列$X^{\mathrm{T}}X$は下記のように表すことができる。
$$
\large
\begin{align}
X X^{\mathrm{T}} &= – \frac{1}{2} \left( I_{n} – \frac{1}{n} J_{n} \right) D \left( I_{n} – \frac{1}{n} J_{n} \right) \quad (1) \\
I_{n} &= \left(\begin{array}{c} 1 & 0 & 0 & \cdots & 0 \\ 0 & 1 & 0 & \cdots & 0 \\ 0 & 0 & 1 & \cdots & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & 0 & \cdots & 1 \end{array} \right) \\
J_{n} &= \left(\begin{array}{c} 1 & 1 & 1 & \cdots & 1 \\ 1 & 1 & 1 & \cdots & 1 \\ 1 & 1 & 1 & \cdots & 1 \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 1 & 1 & 1 & \cdots & 1 \end{array} \right)
\end{align}
$$

・二重中心化の導出
以下、$(1)$式の導出を行う。導出にあたって、$\displaystyle d_{ij}^{2}, \sum_{i=1}^{n} d_{ij}^{2}, \sum_{j=1}^{n} d_{ij}^{2}$を$x_i, x_j, a$を用いて表すことを考える。
$$
\large
\begin{align}
d_{ij}^{2} &= ||x_i-x_j||^2 = ||x_i||^2 + ||x_j||^2 – 2 x_i^{\mathrm{T}} x_j \\
\sum_{i=1}^{n} d_{ij}^{2} &= \sum_{i=1}^{n} (||x_i||^2 + ||x_j||^2 + – x_j^{\mathrm{T}} x_i) \\
&= \sum_{i=1}^{n} ||x_i||^2 + n||x_j||^2 + – x_j^{\mathrm{T}} \sum_{i=1}^{n} x_i = a + n||x_j||^2 \\
\sum_{j=1}^{n} d_{ij}^{2} &= \sum_{j=1}^{n} (||x_i||^2 + ||x_j||^2 + – x_i^{\mathrm{T}} x_j) \\
&= n||x_i||^2 + \sum_{j=1}^{n} ||x_j||^2 – x_i^{\mathrm{T}} \sum_{j=1}^{n} x_j = n||x_i||^2 + a
\end{align}
$$

上記の計算にあたっては$x_i^{\mathrm{T}} x_j = x_j^{\mathrm{T}} x_i$が成立することを用いた。また、$\displaystyle \sum_{i=1}^{n} \sum_{j=1}^{n} d_{ij}^{2}$は下記のように表すことができる。
$$
\large
\begin{align}
\sum_{i=1}^{n} \sum_{j=1}^{n} d_{ij}^{2} &= \sum_{i=1}^{n} \sum_{j=1}^{n} (||x_i||^2 + ||x_j||^2 + – x_i^{\mathrm{T}} x_j) \\
&= 2na
\end{align}
$$

ここまでの内容を元に、$x_i^{\mathrm{T}} x_j$は下記のように表すことができる。
$$
\large
\begin{align}
d_{ij}^{2} &= ||x_i||^2 + ||x_j||^2 – 2 x_i^{\mathrm{T}} x_j \\
x_i^{\mathrm{T}} x_j &= – \frac{1}{2} (d_{ij}^2 – (||x_i||^2 + ||x_j||^2)) \\
&= – \frac{1}{2} \left( d_{ij}^2 – \frac{1}{n} \left( \sum_{i=1}^{n} d_{ij}^{2} – a \right) – \frac{1}{n} \left( \sum_{j=1}^{n} d_{ij}^{2} – a \right) \right) \\
&= – \frac{1}{2} \left( d_{ij}^2 – \frac{1}{n} \sum_{i=1}^{n} d_{ij}^{2} – \frac{1}{n} \sum_{j=1}^{n} d_{ij}^{2} – \frac{2a}{n} \right) \\
&= – \frac{1}{2} \left( d_{ij}^2 – \frac{1}{n} \sum_{i=1}^{n} d_{ij}^{2} – \frac{1}{n} \sum_{j=1}^{n} d_{ij}^{2} – \frac{2na}{n^2} \right) \\
&= – \frac{1}{2} \left( d_{ij}^2 – \frac{1}{n} \sum_{i=1}^{n} d_{ij}^{2} – \frac{1}{n} \sum_{j=1}^{n} d_{ij}^{2} – \frac{1}{n^2} \sum_{i=1}^{n} \sum_{j=1}^{n} d_{ij}^{2} \right)
\end{align}
$$

ここで「行列の成分表示」の表記に基づいて行列$A$の$(i,j)$成分を$(A)_{ij}$のように表すと考えるとき、$(1)$式に関して下記のような式がそれぞれ成立する。
$$
\large
\begin{align}
(D)_{ij} &= d_{ij}^2 \\
\left( \frac{1}{n} J_{n} D \right)_{ij} &= \frac{1}{n} \left(\begin{array}{cccc} 1 & 1 & \cdots & 1 \end{array} \right) \left(\begin{array}{c} d_{1j}^2 \\ d_{2j}^2 \\ \vdots \\ d_{nj}^2 \end{array} \right) \\
&= \frac{1}{n} \sum_{i=1}^{n} d_{ij}^2 \\
\left( \frac{1}{n} D J_{n} \right)_{ij} &= \frac{1}{n} \left(\begin{array}{cccc} d_{i1}^2 & d_{i2}^2 & \cdots & d_{in}^2 \end{array} \right) \left(\begin{array}{c} 1 \\ 1 \\ \vdots \\ 1 \end{array} \right) \\
&= \frac{1}{n} \sum_{j=1}^{n} d_{ij}^2 \\
\left( \frac{1}{n^2} J_{n} D J_{n} \right)_{ij} &= \frac{1}{n^2} \left(\begin{array}{cccc} 1 & 1 & \cdots & 1 \end{array} \right) \left(\begin{array}{cccc} d_{11}^2 & d_{12}^2 & \cdots & d_{1n}^2 \\ d_{21}^2 & d_{22}^2 & \cdots & d_{2n}^2 \\ \vdots & \vdots & \ddots & \vdots \\ d_{n1}^2 & d_{n2}^2 & \cdots & d_{nn}^2 \end{array} \right) \left(\begin{array}{c} 1 \\ 1 \\ \vdots \\ 1 \end{array} \right) \\
&= \frac{1}{n^2} \left(\begin{array}{cccc} 1 & 1 & \cdots & 1 \end{array} \right) \left(\begin{array}{c} \sum_{j=1}^{n} d_{1j}^2 \\ \sum_{j=1}^{n} d_{2j}^2 \\ \vdots \\ \sum_{j=1}^{n} d_{nj}^2 \end{array} \right) \\
&= \frac{1}{n^2} \sum_{i=1}^{n} \sum_{j=1}^{n} d_{ij}^{2}
\end{align}
$$

上記には$d_{ii}^2$が出てくるが$d_{ii}^2=0$と考えれば、$D$の行列の定義に反しない。上記を元に下記のように考えることができる。
$$
\large
\begin{align}
(X X^{\mathrm{T}})_{ij} &= x_i^{\mathrm{T}} x_j \\
&= – \frac{1}{2} \left( d_{ij}^2 – \frac{1}{n} \sum_{i=1}^{n} d_{ij}^{2} – \frac{1}{n} \sum_{j=1}^{n} d_{ij}^{2} – \frac{1}{n^2} \sum_{i=1}^{n} \sum_{j=1}^{n} d_{ij}^{2} \right) \\
&= – \frac{1}{2} \left( D_{ij} – \left( \frac{1}{n} J_{n} D \right)_{ij} – \left( \frac{1}{n} D J_{n} \right)_{ij} – \left( \frac{1}{n^2} J_{n} D J_{n} \right)_{ij} \right) \\
&= – \frac{1}{2} \left( I_{n} – \frac{1}{n} J_{n} \right) D \left( I_{n} – \frac{1}{n} J_{n} \right)_{ij} \\
&= \left( – \frac{1}{2} \left( I_{n} – \frac{1}{n} J_{n} \right) D \left( I_{n} – \frac{1}{n} J_{n} \right) \right)_{ij}
\end{align}
$$

上記より$(1)$式が成立することがわかる。

・参考
$\displaystyle \left( I_{n}-\frac{1}{n}J_{n} \right)A\left( I_{n}-\frac{1}{n}J_{n} \right)$の計算による二重中心化の計算の詳細の確認
https://www.hello-statisticians.com/explain-terms-cat/double_centering1.html

 

エッカート・ヤング分解と寄与度

$n \times n$の対称行列$B=X X^{\mathrm{T}}$の固有値$\lambda_1,…,\lambda_n$を元に下記のように対角行列$\Lambda$を考える。
$$
\large
\begin{align}
\Lambda = \left(\begin{array}{cccc} \lambda_1 & 0 & \cdots & 0 \\ 0 & \lambda_2 & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & \lambda_n \end{array} \right)
\end{align}
$$

このとき$\Lambda^{1/2}$は下記のように表すことができる。
$$
\large
\begin{align}
\Lambda^{\frac{1}{2}} = \left(\begin{array}{cccc} \sqrt{\lambda_1} & 0 & \cdots & 0 \\ 0 & \sqrt{\lambda_2} & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & \sqrt{\lambda_n} \end{array} \right)
\end{align}
$$

ここで$B=X X^{\mathrm{T}}$の固有ベクトルを元に直交行列$U$を考えると$B=U \Lambda U^{\mathrm{T}}$のように表すことができる。このとき$X=U \Lambda^{1/2}$のように考えると$X X^{\mathrm{T}} = (U \Lambda^{1/2})(U \Lambda^{1/2})^{\mathrm{T}} = U \Lambda^{1/2} (\Lambda^{1/2})^{\mathrm{T}} U^{\mathrm{T}} = U \Lambda U^{\mathrm{T}} = B$が成立する。これをエッカート・ヤング分解(Eckart-Young decomposition)という。

このように考えることで、$B$がわかっているときに$X$を計算することができる。ここで計算した$X$は$B$が非負定値行列であれば$D$と対応する位置を表す。

 

具体例

行列分解まとめ② 〜Factorization Machines・DeepFM〜

行列分解まとめ①では特異値分解(SVD; singular value decomposition)などの基本的な行列分解の手法に関してまとめました。②では近年用いられているFactorization MachinesやDeepFMに関してご紹介します。

Factorization Machines

Factorization Machinesの問題設定

特異値分解(SVD; singular value decomposition)などの行列分解の手法では行列$X$に対して$X = UV^{\mathrm{T}}$のような分解を行うことのみを考えたが、Factorization Machinesでは回帰のように$(X,y)$に対して主に取り扱いが行われる。Factorization Machinesの論文のFigure.$1$では下記のように$(X,y)$の具体例が紹介される。

Factorization Machines論文Fig.$1$

ここで上記の例では「青枠のUserがオレンジの枠のMovieを$y$のように評価した」ということに基づいて$X$と$y$の具体例が作成されている。このとき、特定のUserに対する特定のMovieの評価がないときに、評価を予測するというのがここでの問題設定である。

上記の$X$は行列の要素の多くが$0$を持つスパース(sparse)な行列になるが、このようなスパースな行列の取り扱いにあたってFactorization Machinesが考案された。

・参考
Factorization Machines論文 Section.$2$ Prediction Under Sparsity

Factorization Machinesの予測式

Factorization Machinesの予測式は下記のように定義される。
$$
\large
\begin{align}
\hat{y}(x) &= w_{0} + \sum_{i=1}^{n} w_i x_i + \sum_{i=1}^{n} \sum_{j=i+1}^{n} \langle v_i, v_j \rangle x_i x_j \quad (1) \\
w_{0} &\in \mathbb{R}, w \in \mathbb{R}^{n}, V \in \mathbb{R}^{n \times k} \quad (2) \\
\langle v_i, v_j \rangle &= \sum_{f=1}^{k} v_{i,f} \cdot v_{j,f} \quad (3)
\end{align}
$$

上記はFactorization Machines論文の$(1)$〜$(3)$式に対応する。ここで$n$は$X$の行数や$y$の個数に対応し、$k$は$V$で表される潜在空間の次元に対応すると考えることができる。式だけだと少々分かりにくいので概略を下記のように図式化を行った。

Factorization Machinesの概略 左$)$ $\displaystyle \sum_{i=1}^{n} w_i x_i$、右$)$ $\displaystyle \sum_{i=1}^{n} \sum_{j=i+1}^{n} \langle v_i, v_j \rangle x_i x_j$

・参考
Factorization Machines論文 Section.$3$-A Factorization Machine Model

Factorization Machinesのパラメータ推定

Factorization Machinesのパラメータ推定では前項の$(2)$式で表した、$w_{0} \in \mathbb{R}, w \in \mathbb{R}^{n}, V \in \mathbb{R}^{n \times k} $の推定を行う。

ここでパラメータの推定にあたっては確率的勾配降下法(SGD)を用いる。$(1)$式のそれぞれのパラメータでの偏微分はそれぞれ下記のように表される。
$$
\large
\begin{align}
\frac{\partial}{\partial w_0} \hat{y}(x) &= 1 \\
\frac{\partial}{\partial w_i} \hat{y}(x) &= x_i \\
\frac{\partial}{\partial v_{i,f}} \hat{y}(x) &= x_i \sum_{j=1}^{n} v_{j,f} x_{j} – v_{i,f} x_{i}^2
\end{align}
$$

ここで上記の$\hat{y}(x)$はloss functionそのものではないが、論文の数式をそのまま用いた。一方で、loss functionが予測と実測の差を数値化したものである以上、合成関数の微分的に上記の式が出てくることは理解しておくと良い。たとえば二乗和誤差の場合は$((y-\hat{y}(x))^2)’ = 2(y-\hat{y}(x))\hat{y}(x)’$のような計算になるので、$\hat{y}(x)$の偏微分がわかっていればスムーズに各勾配の要素の計算を行うことができる。

また、$\displaystyle \frac{\partial}{\partial v_{i,f}} \hat{y}(x)$の式がやや複雑であるが、$\displaystyle \hat{y}(x) = w_{0} + \sum_{i=1}^{n} w_i x_i + \sum_{i=1}^{n} \sum_{j=i+1}^{n} \langle v_i, v_j \rangle x_i x_j$の式に出てくる$\displaystyle \sum_{i=1}^{n} \sum_{j=i+1}^{n}$より、$i < j$であり、これは上三角行列のように対角成分と左下の要素が$0$である行列を考え、これに基づいて式を解釈するとわかりやすいと思われる。

・参考
Factorization Machines論文 Section.$3$-C Learning Factorization Machines

DeepFM

DeepFMの概要

DeepFMはFactorization-Machineに基づくニューラルネットワークであり、FM ComponentとDeep Componentの二つの要素から構成される。このことはDeepFM論文の$(1)$式で表されるので下記に表した。
$$
\large
\begin{align}
\hat{y} &= sigmoid(y_{FM} + y_{DNN}) \\
sigmoid(x) &= \frac{1}{1+\exp(-x)}, \quad \hat{y} \in (0,1)
\end{align}
$$

また、パラメータ学習にあたってはAdamなどが用いられる。以下、FM ComponentとDeep Componentの詳細に関してそれぞれ確認を行う。

FM Component

FM Componentは前節で取り扱ったFactorization Machinesの考え方に基づく。”Besides a linear (order-1) interactions among features, FM models pairwise (order-2) feature interactions as inner product of respective feature latent vectors”と表現されるように、特徴量の$1$次式に加えて、$2$つの特徴量の相互作用を潜在ベクトルの内積に基づいて計算すると考えると良い。Factorization Machinesで定義される$V$は潜在空間(latent space)を取り扱っており、これは通常の行列分解を潜在空間と考えることと同様に理解することができる。

DeepFM論文の$(2)$式ではFM Componentが下記のように表記される。
$$
\large
\begin{align}
y_{FM} &= \langle w, x \rangle + \sum_{i=1}^{d} \sum_{j=i+1}^{d} \langle V_i, V_j \rangle x_i \cdot x_j
\end{align}
$$

上記はFactorization Machinesの式を行列表記を用いて整理を行った式であると考えることができる。これを図式化を行ったのが論文のFigure.$2$である。

DeepFM論文 Figure 2: The architecture of FM

論文のFigure.$2$の解釈にあたっては、Dense Embeddingsが$V_i, V_j$に対応すると考えればよい。また、FM layerではAdditionとInner Productの二つの演算が表されていることにも着目すると良い。

Deep Component

Deep Componentの概要は論文のFigure.$3$で図示される。

DeepFM論文 Figure 3: The architecture of DNN

上記のHidden Layerより上に関しては基本的には通常のMulti Layer Perceptronと同様であるが、ここで着目する必要があるのがDense Embeddingを作成する処理である。Dense Embeddingの作成は論文のFigure.$4$で表される。

DeepFM論文 Figure 4: The structure of the embedding layer

上記で表されるembedding layerの構造の理解にあたっては、「①$0$以外の値を持つ特徴量(Field)を固定長$k$のembeddingに変換する」、「②Input LayerからEmbeding Layerの計算にあたっては潜在パラメータの$V$を用いる」をそれぞれ抑えておくとよい。

・参考
DeepFM論文 Section.$2$.$1$ DeepFM

参考

・Factorization Machines論文
・DeepFM論文