以下では、$\Theta$を統計モデル統計モデル$\{p(\cdot|\theta)\}_{\theta\in \Theta}$のパラメータ集合、$\Theta_*$をKL情報量を最小にするという意味で最適なパラメータ集合とします。$\Theta_*$に含まれる任意の$\theta,\theta’$が
$$p(x|\theta)=p(x|\theta’),\hspace{5mm} x\in {\cal X}$$を満足するとき、その統計モデル$\{p(\cdot|\theta)\}_{\theta\in \Theta}$は同質であるといいます。また、
$${\mathbb E}_X\left[\left\{\log \frac{p(X|\theta_*)}{p(X|\theta)}\right\}^2\right]\leq
c{\mathbb E}_X\left[\log \frac{p(X|\theta_*)}{p(X|\theta)}\right]$$なる定数$c>0$が存在するとき、「相対的に有限な分散をもつ」といいます。命題2の1「相対的に有限な分散をもてば、同質」の証明(第3章の付録)ので、
\begin{eqnarray*}
0&=&D(q\| p(\cdot|\theta_2))-D(q\| p(\cdot|\theta_1))=\int_{\cal X}q(x)f(x,\theta_1,\theta_2)dx\\
&\geq &\gamma \int_{\cal X} q(x)f(x,\theta_1,\theta_2)^2dx\geq 0
\end{eqnarray*}なる定数$\gamma>0$が存在するので、$f(\cdot,\theta_1,\theta_2)=\log \frac{p(x|\theta_+)}{p(x|\theta)}$は関数として0となり、$\theta_1,\theta_2$が同じ分布になるとしています。

質問をいただきましたが、「関数として0」は正しくありません。正確には、「確率$q(\cdot)$に関していたるところ0」であるが正しいです(命題1を用いています)。そして、同質の定義も、確率$q(\cdot)$に関していたるところ$p(x|\theta)=p(x|\theta’)$というようにすべきです。

同様に、実現可能である$p(x|\theta)=q(x)$も「確率$q(\cdot)$に関していたるところ0」で等号が成立するというようにすべきです。

次回の改定では、修正します。

確率空間を$(E,{\cal F},\mu)$とします。機械学習のためのカーネル第5章で、$k: E\times E\rightarrow {\mathbb R}$が可測であることを仮定しています。この意味について確認したいと思います。

任意の$\mathbb R$のBorel集合$B$について、$$\{\omega\in E\mid f(\omega)\in B\}\in {\cal F}$$であるとき、$f: E\rightarrow {\mathbb R}$は可測であるといいます。

$k: E\times E\rightarrow {\mathbb R}$が可測であるというのは、$\{(x,y)\in E\times E\mid k(x,y)\in B\}\in {\cal F}\times {\cal F}$を意味しています。ここで、${\cal F}\times {\cal F}$は1変数の可測集合の直積集合$C:=\{A\times B|A\in {\cal F},B\in {\cal F}\}$をすべて含む最小の$\sigma$集合体$\sigma(C)$になります。この仮定のもとで、以下のふたつが成立します。

• 各$x\in E$について、$E\ni y\mapsto k(x,y)\in {\mathbb R}$は可測
• $E\ni x\mapsto k(x,x)\in {\mathbb R}$は可測

いずれの証明でも、$k$の可測性から、$B$を任意の${\mathbb R}$のBorel集合としたときに、$A=k^{-1}(B) \in {\cal F}\times {\cal F}$であることを用います。

1つめは、$k_x(\cdot)=k(x,\cdot)$, $f_x: E\ni y\mapsto (x,y)\in E\times E$とおくと、$k_x=k\circ f_x$とできます。
このとき、$$k_x^{-1}(B)=f_x^{-1}\circ k^{-1}(B)=f_x^{-1}(A)$$となり、$A$の第2成分で$\cal F$の要素となります。したがって、$k_x$は可測になります。

２つめは、$g: E\ni x\rightarrow (x,x)$, $h: E\ni x\mapsto k(x,x)$として、$h=k\circ g$とできます。このとき、$$h^{-1}(B)=g^{-1}\circ k^{-1}(B)=g^{-1}(A)$$が得られます。このとき、$A$は同じ事象の対からなる事象となり、$h$は可測になります。

いろいろな教科書にあると思いますが、講義でよく用いるので、私の好きな書き方で証明してみました。「$H_0$を内積空間とする。Hilbert空間$H$があって、$H_0$は$H$の稠密な部分空間と等長的同型になる。このような$H$($H_0$の完備化)は同型を除いて一意に決まる。」。下記の証明では、実数全体$H={\mathbb R}$は有理数からなるCauchy列の集合の同値類(同じ値に収束するCauchy列どうしを同値とする)と読み替えることができます。

$H$とその内積を定義

体$K$上の内積空間$H_0$のCauchy列全体を$X$、$\{x_n\}_{n=1}^\infty,\{y_n\}_{n=1}^\infty\in X$の同値関係$\{x_n\}_{n=1}^\infty\sim\{y_n\}_{n=1}^\infty$を$\|x_n-y_n\|_{H_0}\rightarrow 0$ ($n\rightarrow \infty$)とし、その商を$H:=X/\sim$と書くものとします。そして、$\{x_n\}_{n=1}^\infty\in X$のクラスを$[\{x_n\}_{n=1}^\infty]\in H$と書くと、$\{x_n\}_{n=1}^\infty,\{y_n\}_{n=1}^\infty\in X$, $\alpha,\beta\in {K}$について、$[\{x_n\}_{n=1}^\infty],[\{y_n\}_{n=1}^\infty]\in H$, $\alpha x_n+\beta y_n\in H_0$, $$\alpha[\{x_n\}_{n=1}^\infty]+\beta[\{y_n\}_{n=1}^\infty]=[\{\alpha x_n+\beta y_n\}_{n=1}^\infty]\in H$$が成立し、$H$は線形空間となります(クラスの代表元によらないことはお確かめください)。また、$H$における内積を$$\langle [\{x_n\}_{n=1}^\infty], [\{y_n\}_{n=1}^\infty]\rangle_H:=\lim_{n\rightarrow\infty}\langle x_n, y_n\rangle_{H_0}$$で定義します。右辺の極限が存在し、その代表元によらないことは、以下のように示せます。$\{x_n\}_{n=1}^\infty\in X$に対して、$|\|x_n\|_{H_0}-\|x_m\|_{H_0}|\leq \|x_n-x_m\|_{H_0}$より$\{\|x_n\|_{H_0}\}$はCauchyで、したがって有界です。そして、Cauchy列$\{x_n\}_{n=1}^\infty,\{y_n\}_{n=1}^\infty \in X$に対して、\begin{eqnarray*}|\langle x_n,y_n\rangle_{H_0} – \langle x_m,y_m\rangle_{H_0}|&\leq& |\langle x_n-x_m,y_n\rangle_{H_0}| + |\langle x_m,y_n-y_m\rangle_{H_0}|\\&\leq& \|y_n\|\cdot \|x_n-x_m\|_{H_0} + \|x_m\|\cdot \|y_n-y_m\|_{H_0}\rightarrow 0\end{eqnarray*}が成立し、実数列$\langle x_n,y_n\rangle_{H_0}$は収束します。また、$\{x_n\}_{n=1}^\infty\sim\{x_n’\}_{n=1}^\infty$, $\{y_n\}_{n=1}^\infty\sim\{y_n’\}_{n=1}^\infty$のとき、
$$|\langle x_n,y_n\rangle_{H_0}-\langle x_n’,y_n’\rangle_{H_0}|\leq \|y_n’\|\cdot \|x_n’-x_n\|_{H_0} + \|x_n\|\cdot \|y_n-y_n’\|_{H_0}\rightarrow 0$$より、代表元のとりかたによらないことがわかります。

$H_0$は$H$で稠密

次に、$T: H_0\rightarrow H$を$x\mapsto [\{x,x,\ldots\}]$で定義すると、$\langle Tx,Tx\rangle_{H}=\langle x,x\rangle_{H_0}$を満たすので等長写像です。$T$の像は、$H$で稠密でもあります。実際、$\{x_n\}_{n=1}^\infty\in X$について、$Tx_l=[\{x_l,x_l,\ldots\}]$は、\begin{eqnarray*}\|Tx_l-[\{x_n\}_{n=1}^\infty]\|_H^2&=&
\langle Tx_l,Tx_l\rangle_{H}-2\langle Tx_l,[\{x_n\}_{n=1}^\infty]\rangle_{H}+\langle [\{x_n\}_{n=1}^\infty],[\{x_n\}_{n=1}^\infty]\rangle_{H}\\
&=&\langle x_l,x_l\rangle_{H_0}-2\lim_{n\rightarrow \infty}\langle x_l,x_n\rangle_{H_0}+\lim_{n\rightarrow \infty}\langle x_n,x_n\rangle_{H_0}=\lim_{n\rightarrow \infty}\|x_l-x_n\|^2_{H_0}
\end{eqnarray*}さらに、$\{x_n\}_{n=1}^\infty$はCauchyであるので、$l\rightarrow \infty$で右辺は0となります。すなわち、$x_1\leftrightarrow [\{x_1,\ldots\}]$, $x_2\leftrightarrow [\{x_2,\ldots\}]$, $\ldots$, が$[\{x_1,x_2,\ldots\}]$に収束します。$[\{x_n\}_{n=1}^\infty]$は$H$の任意の要素ですから、$H_0$は$H$の稠密な部分空間であることが示せました。

$H$は完備

最後に$H$の完備性を示します。$\{x_n^{(1)}\}_{n=1}^\infty, \{x_n^{(2)}\}_{n=1}^\infty,\ldots \in X$について、$[\{x_n^{(1)}\}_{n=1}^\infty], [\{x_n^{(2)}\}_{n=1}^\infty],\ldots \in H$がCauchyであると仮定します。内積の定義から、$$\|[\{x_n^{(l)}\}_{n=1}^\infty]-[\{x_n^{(m)}\}_{n=1}^\infty]\|_H=\lim_{n\rightarrow\infty}\|x_n^{(l)}-x_n^{(m)}\|_{H_0}$$であり、$[\{x_n^{(1)}\}_{n=1}^\infty], [\{x_n^{(2)}\}_{n=1}^\infty],\ldots$がCauchyであるので、任意の$\epsilon>0$について、ある正整数$N$があって、$$l,m,n\geq N \Longrightarrow \|x_n^{(l)}-x_n^{(m)}\|_{H_0}<\epsilon$$とできます。また、$$m,n\geq N \Longrightarrow \|x_n^{(n)}-x_m^{(m)}\|_{H_0}<\epsilon$$とできるので、$\{x_n^{(n)}\}_{n=1}^\infty$は$H_0$のCauchy列となります。特に、$\{x_n^{(n)}\}_{n=1}^\infty\in X$, $[\{x_n^{(n)}\}_{n=1}^\infty]\in H$が成立します。したがって、$$\|[\{x_n^{(n)}\}_{n=1}^\infty]-[\{x_n^{(l)}\}_{n=1}^\infty]\|_H^2=\lim_{n\rightarrow \infty}\|x_n^{(n)}-x_n^{(l)}\|^2_{H_0}$$において$l$を十分大きくとると、右辺を十分に小さくできます。したがって、$l\rightarrow \infty$で、$$H\ni [\{x_n^{(l)}\}_{n=1}^\infty]\rightarrow [\{x_n^{(n)}\}_{n=1}^\infty]\in H$$が示されました。

平滑化スプラインの問題は、
$$\sum_{i=1}^N\{y_i-f(x_i)\}^2+\lambda \int_0^1 \{f^{(q)}(x)\}^2dx \rightarrow {\rm min}$$を最小にする$f \in H:=W_q[0,1]$を求める問題として定式化できます。$q=2$の場合は、3次の自然なスプラインの範囲で最適解が存在することが知られています (「統計的機械学習の数理100問」の第6章)。また、Sobolev空間は$H=H_0+H_1$という直交するHilbert空間の直和で書けて、$P$を$H$から$H_1$への写像とすると、上記の問題は$$\sum_{i=1}^N\{y_i-f(x_i)\}^2+\lambda \|Pf\|^2dx \rightarrow {\rm min}$$とかけます。さらに、$H_0,H_1$もRKHSであって、その再生核を$k_0,k_1$と書くと、$H$の再生核は$k=k_0+k_1$となります。

まず、$M={\rm span}(k(x_1,\cdot),\ldots,k(x_N,\cdot))$、$M^\perp$をその直交補空間とし、$f=f_*+f_{\perp}$, $f_*\in M$, $f_{\perp}\in M^\perp$とかくと、$$\sum_{i=1}^N\{y_i-f_*(x_i)\}^2+\lambda \|Pf_*\|^2dx \rightarrow {\rm min}$$の最小化に帰着できます。というのも、$f(x_i)=\langle f,k(x_i,\cdot)\rangle_H=f_*(x_i)$とでき、
$$\|Pf\|^2= \|Pf_*+Pf_\perp\|^2=\|Pf_*|^2+\|Pf_\perp\|^2\geq \|Pf_*|^2$$ここで、$\langle Pf_*,Pf_\perp\rangle_H=0$を用いました。したがって、一般性を失うことなく、$f=\sum_{i=1}^N\beta_ik(x_i,\cdot)$とおいて、その最小化をはかればよいことになります。その場合、$$Pf=\sum_{i=1}^N\beta_iPk(x_i,\cdot)=\sum_{i=1}^N\beta_ik_1(x_i,\cdot)$$となるので、
$$\sum_{i=1}^N\{y_i-\sum_{j=1}^N\beta_j k(x_i,x_j)\}^2+\lambda \sum_{i=1}^N\sum_{j=1}^N\beta_i\beta_jk_1(x_i,x_j)$$
$K=(k(x_i,x_j))\in {\mathbb R}^{N\times N}$, $K_1=(k_1(x_i,x_j))\in {\mathbb R}^{N\times N}$, $y=[y_1,\ldots,y_N]^{\top}$とおけば、$\beta=[\beta_1,\ldots,\beta_N]^{\top}$は
$$\beta=(K^2+\lambda K_1)^{-1}Ky$$となります。

正定値対称行列$A\in {\mathbb R}^{N\times N}$を下三角行列$R\in {\mathbb R}^{N\times N}$を用いて$RR^\top$と書くことを、$A$をCholesky分解するといいます。以下の手順をふみます。
$$ A_{j,i}=\sum_{h=1}^NR_{j,h}R_{i,h}$$という関係式があるので、最初$R=O$(零行列)とし、第$i-1$列まで$R$の値が確定したとして、
1. $R_{i,i}=\sqrt{A_{i,i}-\sum_{h=1}^{i-1}R_{i,h}^2}$
2. $j=i+1,\ldots,N$に対して、$R_{j,i}=(A_{j,i}-\sum_{h=1}^{i-1}R_{j,h}R_{i,h})/R_{i,i}$
によって第$i$列の$i+1$行目以降を設定します。この操作を$i=1,\ldots,N$に対して行います。

そして、$R$の第$r+1$列以降を0とおいた台形の形をした$\tilde{R}$を用いて$A\approx \tilde{R}\tilde{R}^\top$というように(階数$r$に)近似することを$A$を不完全Cholesky分解するといいます。ただ、上記の手順を途中の$r$列目までで止めるだけでは芸がありません。近似の精度が良くなるように工夫します。

まず、最初に$B:=A$とし、上記の手順を以下のように修正します。$$B_{j,i}=\sum_{h=1}^NR_{j,h}R_{i,h}$$という関係式があるので、最初$R=O$とし、第$i-1$列まで$R$の値が確定したとして、$R_{k,k}=\sqrt{B_{k,k}-\sum_{h=1}^{i-1}R_{k,h}^2}$ ($i\leq k\leq N$)を最大にする$k$について、$B$の$i,k$行と$i,k$列、$R$の$i,k$行をいれかえます。その上で、$j=i+1,\ldots,N$に対して、$R_{j,i}=(B_{j,i}-\sum_{h=1}^{i-1}R_{j,h}R_{i,h})/R_{i,i}$によって第$i$列の$i+1$行目以降を設定します。また、大きさ$N$の単位行列で、$Q_i(i,k)=Q_{i}(k,i)=1$, $Q_i(i,i)=Q_{i}(k,k)=0$とした行列$Q_i$を設定します。そのようにして、$r$列目まですすめ、$P:=Q_{1}\cdots Q_r$を用いて以下のようにできます。$$A= PBP^\top\approx P\tilde{R}\tilde{R}^\top P^\top$$

このようにすると、誤差$B-\tilde{R} \tilde{R}^\top$をそのトレースで評価するとして、実行終了時の$\sum_{i=r+1}^NR_{i,i}^2$に等しくなります。実際、$r+1\leq i\leq N$について、$R^\top R$の$(i,i)$成分が$\sum_{j=1}^r R_{i,j}^2$となり、これと$B_{i,i}$との差が$R_{i,i}^2$となります。つまり、$R_{i,i}^2$の値の大きな$i$から順に選択していきます。事前に用意した$r$を用いて$r$列で止めてもよいし、$R_{i,i}^2$が事前に用意したある値にしてもよいと思います。また、上記の手順でいくと、$R^2_{i,i}$が単調減少になることが保証されます。実際、$R_{k,k}^2=B_{k,k}-\sum_{h=1}^{i-1}R_{k,h}^2$が第$i$列として選択されない($R^2_{i,i}より小さい)$と、$i+1$列目では、$R_{k,k}^2=B_{k,k}-\sum_{h=1}^{i}R_{k,h}^2$となって、さらに$R_{k,i}^2$だけ小さくなります。

「渡辺澄夫ベイズ理論 with R/Stan」に関しては、2023年9月の出版以来多くの方にお読みいただき、大変感謝しています。

2024年7月初旬に「渡辺澄夫ベイズ理論 with Python/Stan」を出版するはこびとなりました。理論面も若干加筆していますが、Pythonでの実装の手順を詳しく乗せています。WAIC/WBICは、マルコフ連鎖モンテカルロ(MCMC)の手続きが必要です。私の研究室では、MCMCの実現手段であるStanをGoogle Colab(の有料版の高い方)で走らせる人が多く、Rが得意な人でもPythonで実装する人の方が多いようです。

本書では、PythonでStanを動かす手順を詳しく乗せています。皆様のお役にたてればと考えています。疑問点などございましたら連絡してください。このサイトの記事として取り上げたいと思っています。

命題17に関して若干の修正があります(最終的な結論には影響しません)。
まず、$h(\theta)=\log p(x|\theta)-\log p(x|\theta_*)$, $$\displaystyle d\mu(\theta)=\frac{\displaystyle p(x|\theta)^\alpha p(\theta|x_1,\ldots,x_n)d\theta}{
\displaystyle \int_\Theta p(x|\theta)^\alpha p(\theta|x_1,\ldots,x_n)d\theta}
$$ $$s_k^*(x,\alpha):=\int_\Theta |h(\theta)|^kd\mu(\theta)$$とおきます。また、ヘルダーの不等式(Cauchy-Schwarzの不等式)より、$1\leq k\leq l$について、$$s_l^*(x,\alpha)^{1/l} \leq s_k^*(x,\alpha)^{1/k}$$が成立します。したがって、
$$|s^{(k)}(x,\alpha)|\leq C_k |s_k^*(x,\alpha)|$$が成立します。テキストでは、$s_l(x,\alpha)^{1/l} \leq s_k(x,\alpha)^{1/k}$(偽)を用いて、$s_k(x,\alpha):=\int_\Theta h(\theta)^kd\mu(\theta)$に対して、$|s^{(k)}(x,\alpha)|\leq C_k |s_k(x,\alpha)|$を主張していますが、$|s^{(k)}(x,\alpha)|\leq C_k |s_k^*(x,\alpha)|$が正しいです。

ただ、そのように修正しても7章の結論には影響しません。

行列$A\in {\mathbb R}^{m\times n}$を特異値分解して$A=U\Sigma V^\top$と書くものとします。ただし、 $U\in {\mathbb R}^{m\times n}$, $\Sigma \in {\mathbb R}^{n\times n}$, $V\in {\mathbb R}^{n\times n}$とし、$U^\top U=V^\top V=I$が成立し、$\Sigma$は対角行列でその各成分(特異値)は非負の値をとります。このとき、$$A^\top A= V\Sigma^2V^\top=(V \Sigma V^\top)^2$$とできます。また、正方行列$|A|:=V \Sigma V^\top$のトレースは、$\Sigma$の特異値$\lambda_1,\ldots,n$の和で常に非負になります。$|A|$の$(i,j)$成分は、$|A|_{i,j}=\sum_{k=1}^nv_{i,k}\lambda_{k}v_{j,k}$とかけるので、$|A|$のトレース(対角成分の和)は$$\sum_{i=1}^n |A|_{i,i}=\sum_{i=1}^n\sum_{k=1}^nv_{i,k}\lambda_{k}v_{j,i}=\sum_{k=1}^n \lambda_k \sum_{i=1}^n v_{i,k}^2=\sum_{i=1}^n \lambda_k$$となります。トレースというと$m=n$のときの$\sum_{i=1}^nA_{i,i}$をさすことが多いですが、関数解析では、$\sum_{i=1}^n |A|_{i,i}$のことをトレースノルムといいます。両者は非負定値対称のときに一致し、特異値$\lambda_1,\ldots,\lambda_n$は固有値になります。

同様の概念を一般の線型作用素にも適用できます。線型作用素の$T$について、$|T|:=\{T^*T\}^{1/2}$と書いて、$\sum_{i=1}^\infty\langle |T|e_i,e_i\rangle$と書きます。$\{\cdot\}^{1/2}$は作用素の平方根でと呼ばれるものです。トレースノルムがが有限のときその作用素はトレースクラスであるといいます。そして、自己共役のときに$\sum_{i=1}^\infty\langle Te_i,e_i\rangle$をトレースといい、さらに非負定値のときにトレース、トレースクラスは一致します。「機械学習のためのカーネル」の第2章の最後のトレース作用素の箇所では、自己共役かつ非負定値の線型作用素のみを扱っていて、両者は一致します。

「統計的機械学習の数理」の8章(サポートベクトルマシン)3節と同じ議論は、“The Elements of Statistical Learning”やBishop本(“Pattern Recognition and Machine Learning”)でもされていますが、$y_i(\beta_0+x_i\beta)=1$なる$i$が少なくとも1個存在する理由について書かれていません。拙書で記載したのですが、よく検討したら誤っていました。直感的には、「サポートベクトルが存在しないとマージンをもっと大きくできる」が正解ですが、それを数学的にしかも短く書けないかと思っていました。

まず、$y_1,\ldots,y_N$がすべて等しいということはないと仮定します。この仮定は一般的です。すべて等しい場合は、マージンが存在しません。

その上で、$A=\{i|y_i(\beta_0+x_i\beta)>1\}$, $B=\{i|y_i(\beta_0+x_i\beta)<1\}$の和集合が$\{1,\ldots,N\}$であるとして、矛盾を導きます。ただし、$B$が空集合で$A$ばかりの場合は$\beta=0$となり、$y_i\beta_0>1$がすべての$i\in A$で成立することを意味し、仮定と矛盾します。そこで、以下では、$B$が空集合ではないと仮定します。

まず、(8.18), (8.19)がそれぞれ$$\beta=C\sum_{i\in B}y_ix_i^\top$$ $$\sum_{i\in B}y_i=0$$ (8.16)で$i=1,\ldots, N$ に関して和を取ると、$i\in A$について$\alpha_i=0$であるので、$$C\sum_{i\in B}\{ y_i(\beta_0+x_i\beta)-(1-\epsilon_i)\}=0$$であり、上の2式を代入すると、$$\beta^\top \beta=C\sum_{i\in B}(1-\epsilon_i)$$が成立します。それらを(8.3)に代入すると、$$L_P=C\{-\frac{1}{2}\sum_{i\in B}(1-\epsilon_i)+\sum_{i\in A}\epsilon_i+N\}$$となります。ここで、$i\in A$では(8.17)より$\epsilon_i=0$、$i\in B$では(8.14)より$\epsilon_i>0$となります。このことは、集合$A,B$の要素を変えない範囲で、$i\in B$の$\epsilon_i$のどれかを小さくすれば、$L_P$の値をさらに小さくでき、最適性と矛盾します。したがって、$A\cup B\not=\{1,\ldots,N\}$が成立します。

いずれの場合でも、$y_i(\beta_0+x_i\beta)=1$なる$i$が1個は存在します。