はじめに
『トピックモデルによる統計的潜在意味解析』の学習時のメモです。基本的な内容は、数式の行間を読んで埋めたものになります。本と併せて読んでください。
この記事では、3.6.4項のLDAのハイパーパラメータをニュートン・ラフソン法を用いて推定する方法について書いています。
数学よく解らない自分が理解できるレベルまで落として数式を書き下していますので、分かる人にはかなりくどいです。
【前節の内容】
www.anarchive-beta.com
【他の節一覧】
www.anarchive-beta.com
【この節の内容】
3.6.4 ニュートン・ラフソン法
ニュートン・ラフソン法を用いてLDAのハイパーパラメータを推定する。
・変分下限の導出
3.3.6項(LDAの変分ベイズ法(2))の変分下限(3.102)を用いる。
$$
\begin{align}
F[q(\boldsymbol{z}, \boldsymbol{\theta}, \boldsymbol{\phi} | \boldsymbol{\xi}^{\theta}, \boldsymbol{\xi}^{\phi})]
&= \sum_{k=1}^K \left[
\log \frac{
\Gamma(\sum_{v=1}^V \beta_v)
}{
\prod_{v=1}^V \Gamma(\beta_v)
}
- \log \frac{
\Gamma(\sum_{v=1}^V \xi_{k,v}^{\phi})
}{
\prod_{v=1}^V \Gamma(\xi_{k,v}^{\phi})
}
\right] \\
&\qquad
+ \sum_{k=1}^K \sum_{v=1}^V (
\mathbb{E}_{q(\boldsymbol{z})}[
n_{k,v}
]
+ \beta_v - \xi_{k,v}^{\phi}
)
\mathbb{E}_{q(\boldsymbol{\phi}_k | \boldsymbol{\xi}_k^{\phi})}[
\log \phi_{k,v}
] \\
&\qquad
+ \sum_{d=1}^M \left[
\log \frac{
\Gamma(\sum_{k=1}^K \alpha_k)
}{
\prod_{k=1}^K \Gamma(\alpha_k)
}
- \log \frac{
\Gamma(\sum_{k=1}^K \xi_{d,k}^{\theta})
}{
\prod_{k=1}^K \Gamma(\xi_{d,k}^{\theta})
}
\right] \\
&\qquad
+ \sum_{d=1}^M \sum_{k=1}^K (
\mathbb{E}_{q(\boldsymbol{z}_d)}[
n_{d,k}
]
+ \alpha_k - \xi_{d,k}^{\theta}
)
\mathbb{E}_{q(\boldsymbol{\theta}_d | \boldsymbol{\xi}_d^{\theta})}[
\log \theta_{d,k}
] \\
&\qquad
- \sum_{d=1}^M \sum_{i=1}^{n_d} \sum_{k=1}^K
q(z_{d,i} = k)
\log q(z_{d,i} = k)
\tag{3.102}
\end{align}
$$
(次の微分で他の項は消えるので)ここから、$\boldsymbol{\alpha}$に関係のある項を取り出して
$$
F[\boldsymbol{\alpha}]
= \sum_{d=1}^M \left[
\log \Gamma \left(
\sum_{k=1}^K \alpha_k
\right)
- \sum_{k=1}^K
\log \Gamma(\alpha_k)
+ \sum_{k=1}^K
\alpha_k
\mathbb{E}_{q(\boldsymbol{\theta}_d | \boldsymbol{\xi}_d^{\theta})} [
\log \theta_{d,k}
]
\right]
\tag{3.195}
$$
とおく。
よって、変分下限(3.102)を$\boldsymbol{\alpha}$で微分すると、プサイ関数を用いたDirichlet分布の期待値計算(3.74)を用いて
$$
\begin{align}
\frac{\partial F[\boldsymbol{\alpha}]}{\partial \alpha_k}
&= \sum_{d=1}^M \left[
\frac{
\partial \log \Gamma \left(
\sum_{k=1}^K \alpha_k
\right)
}{
\partial \alpha_k
}
- \frac{
\partial \log \Gamma(\alpha_k)
}{
\partial \alpha_k
}
+ \mathbb{E}_{q(\boldsymbol{\theta}_d | \boldsymbol{\xi}_d^{\theta})} [
\log \theta_{d,k}
]
\right]
\\
&= \sum_{d=1}^M \left[
\Psi \left(
\sum_{k=1}^K \alpha_k
\right)
- \Psi(\alpha_k)
+ \Psi(\xi_{d,k}^{\theta})
- \Psi \left(
\sum_{k=1}^K \xi_{d,k}^{\theta}
\right)
\right]
\\
&= M \left[
\Psi \left(
\sum_{k=1}^K \alpha_{k}
\right)
- \Psi(\alpha_k)
\right]
+ \sum_{d=1}^M \left[
\Psi(\xi_{d,k}^{\theta})
- \Psi \left(
\sum_{k'=1}^K \xi_{d,k}^{\theta}
\right)
\right]
\tag{3.196}
\end{align}
$$
となる。
最適化問題を解くために、$\frac{\partial F[\boldsymbol{\alpha}]}{\partial \alpha_k}$が0となる$\boldsymbol{\alpha}$を求めたいが、この式を解析的に求めるのは困難である。そこで、ニュートン・ラフソン法と呼ばれる最適化手法を用いる。
・テイラー展開による近似
通常の勾配法では、1次のテイラー展開により勾配の導出を行ったが、ニュートン・ラフソン法では、2次のテイラー展開を行う。そのため、目的関数が持つ2次情報を勾配に利用することができる。
・ヘッセ行列
関数$f(\boldsymbol{x}),\ \boldsymbol{x} = (x_1, x_2, \cdots, x_N)$を$\boldsymbol{x} = \boldsymbol{a}$で2次までテイラー展開すると
$$
f(\boldsymbol{x})
= f(\boldsymbol{a})
+ (\boldsymbol{x} - \boldsymbol{a})^{\top}
\frac{\partial}{\partial \boldsymbol{x}} f(\boldsymbol{x})
+ \frac{1}{2}
(\boldsymbol{x} - \boldsymbol{a})^{\top}
\frac{\partial^2}{\partial \boldsymbol{x}^2} f(\boldsymbol{a})
(\boldsymbol{x} - \boldsymbol{a})
$$
となる(たぶん)。ここで、2次偏導関数を並べた行列のことをヘッセ行列$H(\boldsymbol{x})$呼ぶ。
$$
H(\boldsymbol{x})
= \frac{\partial^2}{\partial \boldsymbol{x}^2} f(\boldsymbol{a})
= \begin{bmatrix}
\frac{\partial^2 f(\boldsymbol{x})}{\partial x_1^2}
& \cdots
& \frac{\partial^2 f(\boldsymbol{x})}{\partial x_1 \partial x_N}
\\
\vdots
& \ddots
& \vdots
\\
\frac{\partial^2 f(\boldsymbol{x})}{\partial x_n \partial x_1}
& \cdots
& \frac{\partial^2 f(\boldsymbol{x})}{\partial x_N^2}
\end{bmatrix}
$$
ヘッセ行列の$k'k$成分は
$$
H(\boldsymbol{x})_{k'k}
= \frac{\partial^2}{\partial x_{k'} \partial x_k}
f(\boldsymbol{x})
$$
である。
$F[\boldsymbol{\alpha}]$を$\boldsymbol{\alpha}$の近傍点$\hat{\boldsymbol{\alpha}}$を中心として2次のテイラー展開したものを$\tilde{F}[\boldsymbol{\alpha}]$とおく。
$$
\begin{align}
\tilde{F}[\boldsymbol{\alpha}]
&= F[\hat{\boldsymbol{\alpha}}]
+ (\boldsymbol{\alpha} - \hat{\boldsymbol{\alpha}})^{\top}
\left.
\frac{\partial}{\partial \boldsymbol{\alpha}}
F[\boldsymbol{\alpha}]
\right|_{\boldsymbol{\alpha}=\hat{\boldsymbol{\alpha}}}
+ \frac{1}{2}
(\boldsymbol{\alpha} - \hat{\boldsymbol{\alpha}})^{\top}
\left.
\frac{\partial^2}{\partial \boldsymbol{\alpha}^2}
F[\boldsymbol{\alpha}]
\right|_{\boldsymbol{\alpha}=\hat{\boldsymbol{\alpha}}}
(\boldsymbol{\alpha} - \hat{\boldsymbol{\alpha}})
\\
&= F[\hat{\boldsymbol{\alpha}}]
+ (\boldsymbol{\alpha} - \hat{\boldsymbol{\alpha}})^{\top}
\left.
\frac{\partial}{\partial \boldsymbol{\alpha}}
F[\boldsymbol{\alpha}]
\right|_{\boldsymbol{\alpha}=\hat{\boldsymbol{\alpha}}}
+ \frac{1}{2}
(\boldsymbol{\alpha} - \hat{\boldsymbol{\alpha}})^{\top}
H(\hat{\boldsymbol{\alpha}})
(\boldsymbol{\alpha} - \hat{\boldsymbol{\alpha}})
\tag{3.197}
\end{align}
$$
($|_{\boldsymbol{\alpha}=\hat{\boldsymbol{\alpha}}}$は$\hat{\boldsymbol{\alpha}}$の周りでテイラー展開したって示してるだけで深い意味はないですよね?検索の仕方が分からないのでよく分かってないまま)
ヘッセ行列の$k'k$成分は
$$
H(\hat{\boldsymbol{\alpha}})_{k'k}
= \left.
\frac{\partial^2 F[\boldsymbol{\alpha}]}{\partial \alpha_{k'} \partial \alpha_k}
\right|_{\boldsymbol{\alpha}=\hat{\boldsymbol{\alpha}}}
\tag{3.98}
$$
である。
$F[\boldsymbol{\alpha}]$の近似式$\tilde{F}[\boldsymbol{\alpha}]$も$\boldsymbol{\alpha}$で偏微分して$\boldsymbol{0}$となる$\boldsymbol{\alpha}$を求める。
$$
\begin{align}
\frac{\partial \tilde{F}[\boldsymbol{\alpha}]}{\partial \boldsymbol{\alpha}}
= \cdots
&= 0
\\
\boldsymbol{\alpha}
&= \hat{\boldsymbol{\alpha}}
- H^{-1} (\hat{\boldsymbol{\alpha}})
g(\hat{\boldsymbol{\alpha}})
\tag{3.199}
\end{align}
$$
(…線形代数と微積勉強しますちゃんとします)
ここで
$$
g(\hat{\boldsymbol{\alpha}})
= \left.
\frac{\partial}{\partial \boldsymbol{\alpha}}
F[\boldsymbol{\alpha}]
\right|_{\boldsymbol{\alpha}=\hat{\boldsymbol{\alpha}}}
\tag{3.200}
$$
とおく。
式(3.199)がハイパーパラメータの更新式になる。この計算を行うために、式に含まれるヘッセ行列の逆行列$H^{-1} (\hat{\boldsymbol{\alpha}})$を求めていく。
・ヘッセ行列の逆行列の導出
トリガンマ関数
$$
\Psi^{(1)}(x)
= \frac{d}{dx} \Psi(x)
= \frac{d^2}{dx^2} \log \Gamma(x)
$$
を用いて、ヘッセ行列の$k'k$成分は、式(3.196)より
$$
\begin{align}
H(\boldsymbol{\alpha})_{k'k}
= \frac{\partial^2 F[\boldsymbol{\alpha}]}{\partial \alpha_{k'} \partial \alpha_k}
&= M \left[
\frac{\partial}{\partial \alpha_{k'}}
\Psi \left(\sum_{k=1}^K \alpha_k \right)
- \frac{\partial}{\partial \alpha_{k'}}
\Psi(\alpha_k)
\right]
\\
&= M \left[
\Psi^{(1)} \left(\sum_{k=1}^K \alpha_k \right)
- \delta(k' = k)
\Psi^{(1)} (\alpha_k)
\right]
\tag{3.201}
\end{align}
$$
となる。後の因子は
$$
\Psi^{(1)} (\alpha_k)
= \frac{\partial}{\partial \alpha_{k}}
\Psi(\alpha_k)
= \frac{\partial}{\partial \alpha_{k}^2}
\Gamma(\alpha_k)
= \frac{\partial}{\partial \alpha_{k}^2}
\alpha_k \Gamma(\alpha_k - 1)
$$
である。つまり、$k' \neq k$のとき($\alpha_k$以外の$\boldsymbol{\alpha}$で偏微分したとき)0になる。そこで、$k' = k$のとき1となり$k' \neq k$のとき0となるデルタ関数$\delta(k' = k)$を用いている。
式(3.201)の前の因子を
$$
y = M
\Psi^{(1)} \left(\sum_{k=1}^K \hat{\alpha}_k \right)
\tag{3.203}
$$
とおく。また、後ろの因子を
$$
h_k = - M
\Psi^{(1)} (\hat{\alpha}_k)
\tag{3.202}
$$
とおき、更に、${\rm diag}(\boldsymbol{h})$を$\boldsymbol{h} = (h_1, h_2, \cdots, h_K)$を対角要素とする対角行列とする。
これらと、K個の1を要素とする行ベクトル
$$
\boldsymbol{1}^{\top}
= (1, 1, \cdots, 1)
\tag{3.204}
$$
を用いて、式(3.201)を参考に(k成分以外を含めて)、ヘッセ行列$H(\hat{\boldsymbol{\alpha}})$を
$$
H(\hat{\boldsymbol{\alpha}})
= {\rm diag}(\boldsymbol{h})
+ y
\boldsymbol{1}
\boldsymbol{1}^{\top}
\tag{3.205}
$$
とする。
次に
$$
Y = {\rm diag}(\boldsymbol{h})^{-1}
- \frac{
{\rm diag}(\boldsymbol{h})^{-1}
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
}{
\frac{1}{y}
+ \boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
\boldsymbol{1}
}
\tag{A.4}
$$
とおく。(何故こういう形なのかは深く考えない…と言うかこれが逆行列であることを求めていく。)
これを式(3.205)と掛けると、単位行列を$E$として
$$
\begin{aligned}
H(\hat{\boldsymbol{\alpha}})
Y
&= \Bigl(
{\rm diag}(\boldsymbol{h})
+ y
\boldsymbol{1}
\boldsymbol{1}^{\top}
\Bigr)
\left(
{\rm diag}(\boldsymbol{h})^{-1}
- \frac{
{\rm diag}(\boldsymbol{h})^{-1}
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
}{
\frac{1}{y}
+ \boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
\boldsymbol{1}
}
\right)
\\
&= E
- \frac{
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
}{
\frac{1}{y}
+ \boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
\boldsymbol{1}
}
+ y
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
- \frac{
y
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
}{
\frac{1}{y}
+ \boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
\boldsymbol{1}
}
\\
&= E
+ y
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
- \frac{
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
+ y
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
}{
\frac{1}{y}
+ \boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
\boldsymbol{1}
}
\end{aligned}
$$
となる。
ここで
$$
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
\boldsymbol{1}
= \sum_{k=1}^K \frac{1}{h_k}
$$
より、置き換えると
$$
\begin{align}
H(\hat{\boldsymbol{\alpha}})
Y
&= E
+ y
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
- \frac{
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
+ y
\boldsymbol{1}
\left(
\sum_{k=1}^K \frac{1}{h_k}
\right)
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
}{
\frac{1}{y}
+ \sum_{k=1}^K \frac{1}{h_k}
}
\\
&= E
+ y
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
- \frac{
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
\left\{
1
+ y \left(
\sum_{k=1}^K \frac{1}{h_k}
\right)
\right\}
}{
\frac{1}{y}
+ \sum_{k=1}^K \frac{1}{h_k}
}
\\
&= E
+ y
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
- \frac{
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
y
\left\{
\frac{1}{y}
+ \sum_{k=1}^K \frac{1}{h_k}
\right\}
}{
\frac{1}{y}
+ \sum_{k=1}^K \frac{1}{h_k}
}
\\
&= E
+ y
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
- \boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
y
\\
&= E
\tag{A.5}
\end{align}
$$
となる。
右逆行列の存在定理より
$$
\begin{align}
Y
H(\hat{\boldsymbol{\alpha}})
&= E
\\
Y &= H^{-1} (\hat{\boldsymbol{\alpha}})
\end{align}
$$
となる。
従って、式(A.4)がヘッセ行列$H(\hat{\boldsymbol{\alpha}})$の逆行列である。(逆行列の導出は難しい(あるいは紙面をとる)からこういう形で一応証明したということ?)
・ハイパーパラメータの更新式の導出
ヘッセ行列の逆行列
$$
H^{-1} (\hat{\boldsymbol{\alpha}})
= {\rm diag}(\boldsymbol{h})^{-1}
- \frac{
{\rm diag}(\boldsymbol{h})^{-1}
\boldsymbol{1}
\boldsymbol{1}^{\top}
{\rm diag}(\boldsymbol{h})^{-1}
}{
y^{-1}
+ \sum_{k=1}^K h_k^{-1}
}
\tag{3.206}
$$
の$kk$成分は
$$
H^{-1} (\hat{\boldsymbol{\alpha}})_{kk}
= h_k^{-1}
- \frac{
h_k^{-1}
\sum_{k=1}^K h_k^{-1}
}{
y^{-1}
+ \sum_{k=1}^K h_k^{-1}
}
$$
である(?)。
また、$g(\hat{\boldsymbol{\alpha}})$のk番目の要素$g(\hat{\boldsymbol{\alpha}})_k$は$\frac{\partial F[\boldsymbol{\alpha}]}{\partial \alpha_k}$なので、式(3.196)である。
これらを用いて、$H^{-1} (\hat{\boldsymbol{\alpha}}) g(\hat{\boldsymbol{\alpha}})$のk番目の要素$\Bigl( H^{-1} (\hat{\boldsymbol{\alpha}}) g(\hat{\boldsymbol{\alpha}}) \Bigr)_k$は
$$
\begin{align}
\Bigl(
H^{-1} (\hat{\boldsymbol{\alpha}})
g(\hat{\boldsymbol{\alpha}})
\Bigr)_k
&= \left(
h_k^{-1}
- \frac{
h_k^{-1}
\sum_{k=1}^K h_k^{-1}
}{
y^{-1}
+ \sum_{k=1}^K h_k^{-1}
}
\right)
g(\hat{\boldsymbol{\alpha}})_k
\\
&= \left(
\frac{1}{h_k}
- \frac{1}{h_k}
\frac{
1
}{
y^{-1}
+ \sum_{k=1}^K h_k^{-1}
}
\sum_{k=1}^K \frac{1}{h_k}
\right)
g(\hat{\boldsymbol{\alpha}})_k
\\
&= \frac{
g(\hat{\boldsymbol{\alpha}})_k
}{
h_k
}
- \frac{1}{h_k}
\frac{
1
}{
y^{-1}
+ \sum_{k=1}^K h_k^{-1}
}
\sum_{k=1}^K
\frac{
g(\hat{\boldsymbol{\alpha}})_k
}{
h_k
}
\tag{3.207}
\end{align}
$$
となる(???式が入れ子過ぎて整理できない)。(あと$g(\hat{\boldsymbol{\alpha}})_k$が$\sum_k$の影響下になるのは、$\boldsymbol{1} \boldsymbol{1}^{\top} {\rm diag}(\boldsymbol{h})^{-1} * g(\hat{\boldsymbol{\alpha}})$の計算を理解しないといけないんでしょうね)
従って、式(3.199)より$\hat{\boldsymbol{\alpha}}$を1ステップ前の$\boldsymbol{\alpha}$とすると、$\boldsymbol{\alpha}$の各要素の更新式は
$$
\alpha_k^{(s+1)}
= \alpha_k^{(s)}
- \Bigl(
H^{-1} (\boldsymbol{\alpha}^{(s)})
g(\boldsymbol{\alpha}^{(s)})
\Bigr)_k
\tag{3.209}
$$
が得られる。
単語分布のパラメータ$\boldsymbol{\beta}$も同様に求めることができる。
参考文献
- 佐藤一誠『トピックモデルによる統計的潜在意味解析』(自然言語処理シリーズ 8)奥村学監修,コロナ社,2015年.
おわりに
3章が終わったらしっかりと時間をとって線形代数と微分積分をやります。
【次節の内容】
www.anarchive-beta.com
