全ては収束する〜ボルツマン方程式とResNetにみるtrend to equilibrium~

1081

自己紹介と概要

はじめまして。ShinKunと申します。今回は東京大学航空宇宙工学科／専攻 Advent Calendar 2020一日目、ならびに物工/計数 Advent Calendar 2020二日目の記事として、流体最強の方程式として知られるBoltzmann方程式と、深層学習において非常によく知られるResNetの、ある意味での共通点、equillibrium trendについてお話していきたいと思います。全く異なるように思えるこれら２つの対象は、実は平均場解析と呼ばれる物理的視点により繋がっています。つまり、ある法則に従う多変数の系を直接相手にするのではなく、それらの変数が従う確率分布の発展を考えることで系全体の漸近的な振る舞いを解析できるという発想でこれら２つを理解できます。

本記事ではれより少し進んだ知見を得たいと考えています。結論を先にいうと、適切に学習された（推定対象から適切にパラメータが定められた）ResNetの層方向への情報伝達は、とある運動則に基づくBoltzmann方程式のtrend to equilibriumだと理解できるということになります。

この結論に向かっていきたいのですが、確率分布ないし確率測度の発展方程式を扱い、さらにその漸近挙動を調べることはそれほど簡単なことではありません。その最も単純な理由の一つは確率測度の空間が線形空間でないためであるといえるでしょう。

そこで役に立つのがWasserstein距離$ W_2$です。この距離を用いることによって確率測度を解析的に扱うことができ、特にBoltzmann方程式の漸近挙動、特に平衡解(equilibrium)への収束(trend to equilibrium)を(a priori評価なしに)解析することができます。

以上の流れを踏まえ、今回の記事では

Wasserstein距離の定義・性質を速習（証明はしない）
Boltzmann方程式の漸近挙動を調べる
ResNetはBoltzmann方程式(の"線形化"であるFokker-Planck方程式)とみなして理解してみる

という構成で以下お話していきたいと思います。よろしくお願いします！

Wasserstein距離とその性質速習

この節では一般的な状況として$(X,d)$を可分完備距離空間とします。Boltzmann方程式を述べるだけであれば$X=\mathbb{R}^3,d(x,y)=|x-y|(\text{Euclidian norm})(x,y\in X)$と思っていれば問題ありません。また$P(X)$を$X$上の確率測度全体の空間とし、$P_2(X)$を$P(X)$の元で二次モーメントが有界なもの全体の集合とします。つまり、$\mu\in P_2(X)$について、ある$x_0\in X$が存在して $\int_X d(x_0,x)\dd\mu(x)<\infty$が成り立ちます。

二つの確率測度$\mu,\nu\in P_2(X)$の間のWasserstein距離$W_2$とは、直感的には"砂の山$\mu$を穴$\nu$に移すときの最小コスト"を表します。しかし次の定義からではそれが理解しにくいかもしれません。

$W_2$

確率測度$\mu,\nu\in P_2(X)$に対し
\begin{equation} W_2(\mu,\nu)=\sqrt{\inf_{\pi\in\Pi(\mu,\nu)}\int_{X^2}d(x,y)^2\dd\pi(x,y)}\label{aa} \end{equation}

とする。ここで$\Pi(\mu,\nu)=\left\{\pi\in P(X^2) \mid\forall A\subset X:可測,\pi[A\times X]=\mu[A],\pi[X\times A]=\nu[A]\right\}$

定義に$\inf$が入っており複雑ですが、この$\inf$はきちんと下に有界であり、距離の定義を満たすこと、特に三角不等式を満たすことが知られています。

$W_2$は距離である

上の定義の下限$\inf$を達成する$\pi$は存在し、かつ$W_2$は$P_2(X)$上の距離となる。

証明は Villaniの本の定理7.3を参照してください。

$W_2$は定義は複雑ですが、通常目にする距離が満たしてほしい性質は結構満たしてくれています。例えば次の凸性は基本的ですが有用です。

$W_2$の凸性

任意の$\mu_1,\mu_2,\nu_1,\nu_2\in P_2(X),\alpha\in[0,1]$について
$$ W_2(\alpha\mu_1+(1-\alpha)\mu_2,\alpha\nu_1+(1-\alpha)\nu_2)^2\leq\alpha W_2(\mu_1,\nu_1)^2+(1-\alpha)W_2(\mu_2,\nu_2)^2 $$
が成り立つ。

しかし最も重要な性質と言えうるのは次の位相的な性質です。

$W_2$によって確率測度の弱*収束は距離化できる

$P_2(X)$内の点列$(\mu_k)_{k\in\mathbb{N}}$と$\mu\in P(X)$に対し、以下二つは同値。

$W_2(\mu_k,\mu)\xrightarrow[k\to\infty]{}0$。
$\mu_k\xrightarrow[k\to\infty]{\text{w}\ast}\mu$であり、かつ、ある$x_0\in X$に対し$$\lim_{R\to\infty}\limsup_{k\to\infty}\int_{d(x_0,x)\geq R}d(x_0,x)^2\dd\mu_k(x)=0$$が成り立つ。

ここでいう弱*収束とは、任意の$\varphi\in C_0^\infty(X)$に対して
$$ \lim_{k\to\infty}\int_X\varphi(x)\dd\mu_k(x)=\int_X\varphi(x)\dd\mu(x) $$
となることを指します。したがって通常確率測度の収束を述べる際にはすべての関数$\varphi$について上記が成り立つかチェックしなければいけないので面倒なのですが、Wasserstein距離により上の定理が成り立つのでその必要がなくなります。やったね。

Boltzmann方程式とそのtrend to equilibrium〜田中の定理〜

前節の性質を使ってBoltzmann方程式のtrend to equilibriumを示すことができます。その前に、航空Advent calendar2018内の kanepleaseさんのブログ記事を参考にしながらBoltzmann方程式について復習しましょう。端的に言うとボルツマン方程式とは時刻$t$、位置$x$、速度$v$の粒子の存在確率$f(t,x,v)$が従う方程式のことで、次のような形をしているのでした。
$$ \frac{\partial f}{\partial t}+\boldsymbol{v} \cdot \frac{\partial f}{\partial \boldsymbol{r}}+\frac{\boldsymbol{F}}{m} \cdot \frac{\partial f}{\partial \boldsymbol{v}}=\int_{\boldsymbol{v}_{1}} \int_{\Omega^{\prime}}\left(f f_{1}^{\prime}-f f_{1}\right) \boldsymbol{V} \sigma d \Omega^{\prime} d^{3} \boldsymbol{v}_{1} $$

今回は簡単のため、位置$x$に$f$が依存しない、つまり空間斉次性の仮定を置き、さらに次のような簡略化されたBoltzmann方程式を考えます。
$$ \frac{\partial f}{\partial t}=Q(f,f)=\int_{\mathbb{R}^3}\dd v_\ast\int_{S^2}\dd\sigma\ b(\cos\theta)(f'f'_\ast-ff_\ast)\tag{1} $$
上式の正式な意味は、任意の$\varphi\in C_0^\infty(\mathbb{R}^3)$に対して
$$ \frac{\mathrm{d}}{\mathrm{d}t}\int_{\mathbb{R}^3}\varphi \dd f(t)=\int_{\mathbb{R}^3}\varphi\dd Q(f,f) $$
が成り立つ、ということです。ここで$Q$はBoltzmann collison operaterと呼ばれる$P_2(X)$上の作用素であり、$f,f',f'_\ast,f_\ast$はそれぞれ$f(t,v),f(t,v'),f(t,v'_\ast),f(t,v_\ast)$を略記したものである。ここでさらに$v_\ast$は速度$v$の粒子と衝突した粒子の速度、$v_\ast,v'_\ast$はそれぞれ衝突前の二つの粒子の速度を意味します。何故こういった文字が登場するかはブログ記事を参照することにして、ここでは端的に$v',v'_\ast$は次式で与えられるという結果のみ述べます。
$$ \left\{\begin{array}{l} v^{\prime}=\dfrac{v+v_{*}}{2}+\dfrac{\left|v-v_{*}\right|}{2} \sigma \\ v_{*}^{\prime}=\dfrac{v+v_{*}}{2}-\dfrac{\left|v-v_{*}\right|}{2} \sigma \end{array}\right. \tag{2} $$
端的には運動量保存則とエネルギー保存則(と散乱方向を表すパラメータ$\sigma\in S^2$)を用いることによって導かれるといえます。
また、$b$は(Maxwellian) Boltzmann collison kernelと呼ばれ、散乱方向$\sigma\in S^2$を制御するための非負値関数であり、
$$ \cos\theta=\frac{v-v_\ast}{|v-v_{\ast}|}\cdot\sigma $$
です。

事実として、任意の$k\in S^2$に対し$\int_{S^2} b(k\cdot\sigma)(1-k\cdot\sigma)\dd\sigma<\infty$ならBoltzmann方程式の絵画存在し、一意であること、さらに初期値$f_0(v)$に解$f(t,v)$が(弱*収束の意味で)安定に依存することがわかっています[ DiPerna+,89 ]。以下では常にこの仮定が満たされているものとします。しかしこの仮定のみで解の大域的性質を解析することには技術的困難を伴います。そこで63年に編み出されたのが次のGrad's angular cut-off
$$ \int_{S^2}b(k\cdot\sigma)\dd\sigma=1<+\infty $$
を課した解を考え、その解で一般の解を近似するテクニックです。この仮定を課すと、Boltzmann方程式は次式のように簡潔になることがわかります。
$$ \frac{\partial f}{\partial t}=Q^+(f,f)-f,Q^+(f,f)=\int_{\mathbb{R}^3}\dd v_\ast\int_{S^2}\dd\sigma\ b(\cos\theta)f'f'_\ast\label{simple}\tag{3} $$

上ではBoltzmann方程式の形式的な説明をしてきました。一方、このBoltzmann方程式は運動量保存$\int_{\mathbb{R}^3} f(t,v)v\dd v=\text{const.}$とエネルギー保存$\int_{\mathbb{R}^3}f(t,v)v^2=\text{const.}$を満たすことが変数変換、つまり式(2)に注目して$(v,v_\ast)$と$(v',v'_\ast)$を入れ替えることによりわかります（演習問題）。

これらの保存則の他にも、$H$定理と呼ばれるものが成立するのでした。$H$定理については航空Advent calendar2019内の kanepleaseさんのブログ記事が非常に参考になります。主張を簡単にまとめると"Boltzmann方程式の解は負エントロピーを単調減少させる。さらにエントロピーが一定となる場合は$f\equiv M$という平衡状態しかない"ということです。ここで$M$はMaxwell分布、特にBoltzmann方程式の初期値$f_0$として
$$ \int_{\mathbb{R}^3}f_0(v)v\dd v=0,\int_{\mathbb{R}^3}f_0(v)\frac{v^2}{2}\dd v=\frac{3}{2}\tag{4} $$
というものを選んだ場合は
$$ M(v)=\frac{1}{(2\pi)^{3/2}}\exp\left(-\frac{|v|^2}{2}\right) $$
となります。この$M$はBoltzmann方程式の平衡解、つまり$Q(M,M)=0$となっています。

$H$定理とtrend to equilibrium

実はこの$H$定理を用いてtrend to equilibrium、つまり$f\to M,t\to\infty$を証明することもできるのですが、時刻$t=0$で負エントロピーが上に有界であるという仮定を一度置く必要があります。以下で述べる田中の定理は、その仮定を課さずにtrend to equilibriumを証明することができるという意味で優れているといえます。

以上でtrend to equilibriumを述べるのに必要な準備が終わりました。事実の羅列が多くて申し訳ありません。しかしこれで記号等にも慣れたと思うので、この記事の主定理を述べたいと思います。この定理は主張自体は[ Tanaka,78 ]で述べられたものですが、後にVillaniが整理し直しました[ Villani,03 ]。以下の証明は基本的にVillaniの証明を倣っていますが、部分的に"形式的には成り立つが証明はそれほど自明ではない"ところがあるので、そのような箇所が登場したら逐次remarkします。

Tanaka's theorem

$f_0,g_0\in P_2(\mathbb{R}^3)$を式(4)を満たすような確率測度とする。また、$f_0,g_0$を初期値とするようなBoltzmann方程式(1)の解をそれぞれ$f(t,\bullet),g(t,\bullet)$とする。このとき、以下が成り立つ。

任意の$t\geq0$について
$$ W_2(f(t,\bullet),g(t,\bullet))\leq W_2(f_0,g_0) $$
が成り立つ。

2.$$W_2(f(t,\bullet),M)\to0,\text{as}\ t\to\infty$$特に三角不等式から$W_2(f(t,\bullet),g(t,\bullet))\xrightarrow[t\to\infty]{}0$が成り立つ。

Step 1 上で述べた注意から、Grad's cut-offの仮定を課した解を考えれば十分である（[ DiPerna+,89 ]等の手法によりcut-offされない場合でも成立することがわかる）。そこで以下では簡潔になったBoltzmann方程式(3)を考えればよい。

Step 2 まず主張1を示す。命題2により、$\int_{\mathbb{R}^3}f(v)v\dd v=\int_{\mathbb{R}^3}g(v)v\dd v$を満たす確率測度$f,g\in P_2()$について次が成り立つことを示せば十分である。
\begin{equation} W_2(Q^+(f,f),Q^+(g,g))\leq W_2(f,g)\tag{5} \end{equation}

形式的な$\because$)
$f(t+\Delta t)=f(t)+\Delta t(Q^+(f(t),f(t))-f(t))+\varepsilon$程度なので、$f(t),g(t)$を$f,g$と略記すれば
\begin{align} W_2(f(t+\Delta t),g(t+\Delta t))^2&= W_2(f+\Delta t(Q^+(f,f)-f),g+\Delta t(Q^+(g,g)-g))^2+\varepsilon\\ &=W_2((1-\Delta t)f+\Delta tQ^+(f,f),(1-\Delta t)g+\Delta tQ^+(g,g))^2+\varepsilon\\ &\leq (1-\Delta t)W_2(f,g)^2++\Delta tW_2(Q^+(f,f),Q^+(g,g))^2+\varepsilon\\ &\leq (1-\Delta t)W_2(f,g)^2++\Delta tW_2(f,g)^2+\varepsilon\\ &=W_2(f,g)^2+\varepsilon \end{align}
となる。この議論を区間$[0,t]$を$\Delta t$ごとに分割して行えば大筋は示せる。

Step 3 命題2によりJensenの不等式が従うが、実は式(5)をある$\theta_0\in [0,\pi]$が存在して$b(\cos\theta)\sin\theta=\delta_{\theta_0}$となる場合について証明すれば十分である。

形式的な$\because$)
形式的なDirac関数を用いて$$\int_{S^2}\dd\sigma b(\cos\theta)=\frac{1}{2\pi}\int_0^\pi b(\cos\theta)\sin\theta\dd\theta=\int_0^\pi\lambda_\theta\delta\theta$$と表わせる。これとJensenの不等式から
\begin{align} W_2(Q^+(f,f),Q^+(g,g))^2&=W_2\left(\int\dd v_\ast\int\lambda\delta_{\theta}f'f'_\ast,\int\dd v_\ast\int\lambda\delta_{\theta}g'g'_\ast\right)^2\\ &\leq\int\lambda_\theta W_2\left(\int\dd v_{\ast}\delta_{\theta}f'f'_\ast,\int\dd v_{\ast}\delta_{\theta}g'g'_\ast \right)^2\\ &\leq\int\lambda W_2(f,g)^2\\ &=W_2(f,g)^2 \end{align}
が従う。

Step 4 ここで、任意の$\varphi\in C_0^\infty(\mathbb{R}^3)$について、変数変換により$Q^+(f,f)$は
\begin{align} \int_{\mathbb{R}^3}\varphi Q^+(f,f)&=\int_{\mathbb{R}^3}\dd v\varphi(v)\int_{\mathbb{R}^3}\dd v_{\ast}\int_0^{2\pi}\dd\phi \left. f'f'_{\ast}\right|_{\theta=\theta_0}\\ &= \int_{\mathbb{R}^3}\dd v'\varphi(v)\int_{\mathbb{R}^3}\dd v'_{\ast}\int_0^{2\pi}\dd\phi \left. f'f'_{\ast}\right|_{\theta=\theta_0}\\ &= \int_{\mathbb{R}^3}\dd v\varphi(v')\int_{\mathbb{R}^3}\dd v_{\ast}\int_0^{2\pi}\dd\phi \left. f'f'_{\ast}\right|_{\theta=\theta_0}\\ \end{align}
を満たすので、$Q^+(f,f)$は次のような確率測度とみなせる(これをTanaka's representation of the gainと呼ぶ)。
$$ Q^+(f,f)=\int_{\mathbb{R}^6}\dd v\dd v_{\ast} ff_\ast\Pi_{v,v,\theta_0} = \mathrm{E}\left[\Pi_{V,V_\ast,\theta_0}\right] $$
ここで、$V,V_\ast$は互いに独立な、$f$に従う確率変数であり、$\mathrm{E}$はそれらの確率変数に関する期待値である。また、$\Pi_{v,v_\ast,\theta_0}$は円$C_{c,r,k}$上の一様分布であり、制限によって$\mathbb{R}^3$上の確率分布とみなす。ここで円$C_{c,r,k}$は、中心$c$、半径$r$、そして円の法線ベクトル$k$が
$$ c=\frac{v+v_\ast+(v-v_\ast)\cos\theta_0}{2},r=\frac{|v-v_\ast|}{2}\sin\theta_0,k=\frac{v-v_\ast}{|v-v_\ast|} $$
となるようなものである。

Step 5 ここで、次の補題を示す; $U_{c,r,k}$で$C_{c,r,k}$上を表すこととする。このとき、以下の評価が従う。
$$ W_2(U_{c_1,r_1,k_1},U_{c_2,r_2,k_2})^2\leq |c_1-c_2|^2+r_1^2+r_2^2-r_1r_2(1+|k_1\cdot k_2|) $$

$\because$)
適当なUnitary変換と平行移動により、$c_1=0,k_1=(0,0,1),k_2=(0,-\sin\gamma,\cos\gamma)(\gamma\in[0,\pi/2])$の場合に次が成り立つことを示せば十分である。
$$ W_2(U_{c_1,r_1,k_1},U_{c_2,r_2,k_2})^2\leq |c_2|^2+r_1^2+r_2^2-r_1r_2(1+\cos\gamma) $$
そこで$\Omega= S^1\subset\mathbb{R}^2$を一様分布によって確率空間と見做す。また$\omega\in\Omega$に対し、
\begin{align} V_1(\omega)&=r_1(\cos\omega,\sin\omega,0),\\ V_2(\omega)&=c_2+r_2(\cos\omega,\sin\omega\cos\gamma,\sin\omega\sin\gamma) \end{align}
と置くと、これらはそれぞれ$U_{c_1,r_1,k_1}, U_{c_2,r_2,k_2}$に従う確率変数である。
$ U_{c_1,r_1,k_1}\otimes U_{c_2,r_2,k_2} \in\Pi( U_{c_1,r_1,k_1},U_{c_2,r_2,k_2} )$あるから、$W_2$の定義より
\begin{align} W_2(U_{c_1,r_1,k_1},U_{c_2,r_2,k_2})^2&\leq\mathrm{E}_{V_1\sim U_{c_1,r_1,k_1},V_2\sim U_{c_2,r_2,k_2}}\left[|V_1-V_2|^2\right]\\ &=\frac{1}{2\pi}\int_0^{2\pi}\left(r_1\cos\omega-{c_2}_x-r_2\cos\omega\right)^2+\left(r_1\sin\omega-{c_2}_y-r_2\sin\omega\cos\gamma\right)^2+\left({c_2}_z+r_2\sin\omega\sin\gamma\right)^2\dd\omega\\ &=|c_2|^2+r_1^2+r_2^2-r_1r_2(1+\cos\gamma) \end{align}
であり、示せた。

Step 6 Step5で示した補題を用いてStep3を式(5)を示す。$(V,V_\ast)$を周辺分布がそれぞれ$f$に従う確率変数（独立とは限らない）、$(W,W_\ast)$を周辺分布が$g$に従う確率変数で、$(V,W)$と$(V_\ast,W_\ast)$は独立であるようなものとする。Step4とJensenの不等式から
\begin{align} W_2(Q^+(f,f),Q^+(g,g))^2&=W_2(\mathrm{E}\Pi_{V,V_\ast,\theta_0},\mathrm{E}\Pi_{W,W_\ast,\theta_0})^2\\ &\leq\mathrm{E}W_2(\Pi_{V,V_\ast,\theta_0},\Pi_{W,W_\ast,\theta_0})^2\\ &\leq\mathrm{E}\left|\frac{V+V_\ast+(V-V_\ast)\cos\theta_0}{2}-\frac{W+W_\ast+(W-W_\ast)\cos\theta_0}{2}\right|^2\\ &\quad\quad+\frac{\sin^2\theta_0}{4}\left(|V-V_\ast|^2+|W-W_\ast|^2-|V-V_\ast||W-W_\ast|-|(V-V_\ast)\cdot(W-W_\ast)|\right)\\ &\leq\mathrm{E}\left|\frac{1+\cos\theta_0}{2}(V-W)+\frac{1-\cos\theta_0}{2}(V_\ast-W_\ast)\right|^2+\frac{\sin^2\theta_0}{4}\left|(V-W)-(V_\ast-W_\ast)\right|^2\\ &=\frac{1+\cos\theta_0}{2}\mathrm{E}|V-W|^2+\frac{1-\cos\theta_0}{2}\mathrm{E}|V_\ast-W_\ast|^2\\ &=\mathrm{E}|V-W|^2 \end{align}
が成り立つ。最左辺と最右辺を比較し、最右辺の$(V,W)$は$\Pi(f,g)$の中で任意に動かせるから、$\inf$を取れば主張1を得る。

Step 7 ここで主張1の等号成立が$f=g$のときのみに限ることを示す。等号が成立すると仮定すればStep6での等号成立が必要であり、$(V-V_\ast)\cdot(W-W_\ast)=|V-V_\ast||W-W_\ast|\quad\text{a.s.}$つまり
$$ \frac{V-V_\ast}{|V-V_\ast|}=\frac{W-W_\ast}{|W-W_\ast|}\quad\text{a.s.} $$
が従う。ここで、Brenierの定理[ Villani,03 ]により、輸送写像と呼ばれる可測写像$\psi:\mathbb{R}^3\to\mathbb{R}^3$が一意に存在して$f=\psi\#g$が成り立つ。よって$V=\psi(W),V_\ast=\psi(W_\ast)$が成り立ち、ほとんど全ての$w,w_\ast\in\mathbb{R}^3$に対し
$$ \frac{\psi(w)-\psi(w_\ast)}{|\psi(w)-\psi(w_\ast)|}=\frac{w-w_\ast}{|w-w_\ast|} $$
が従う。これから$\psi$が$\psi(y)=ay+b$のような形式になることがわかるが、式(4)、つまり平均ゼロ、二次モーメントが単位的であるという仮定から$a=1,b=0$となることがわかる。よって$\psi=\mathrm{Id}$であり、これは$f=g$を表す。

他の等号成立も$f=g$で従うので、示したいことは言えた。

Step 8 次に主張2を示す。（中心極限定理の証明と似ているが）$\int f_0(v)|v|^4\dd v<\infty$のときに示せば十分である（$\because$そうでない$f_0$に関しては$1_{|f|v|^4|\leq n}$のようなcut-off functionを掛け、単調収束定理を用いて$n\to\infty$として今から考える場合に帰着させればよい。）。このとき少々の議論によって$\sup_{t\geq0}f(t,v)|v|^4<\infty$であることがいえる。これとCauchy-Schwartzの不等式から
$$ \lim_{K\to\infty}\limsup_{t\geq0}\int_{|v|\geq K}f(t,v)|v|^2\dd v=0 $$
を得る。また主張1から分かる有界性から(Banach-Alaogruの定理の要領で)、部分列に降りて、ある$\mu_0 \in P(X)$が存在し$f(t_k)\xrightarrow[k\to\infty]{\text{w}\ast}\mu_0$である。故に定理3から$W_2(f(t_k),M)\xrightarrow[k\to\infty]{}W_2(\mu_0,M)$が従う。一方、$\mu(t)$を初期値を$\mu_0$とするようなBoltzmann方程式の解とすると、以前remarkした解の初期値に対する安定性[ DiPerna+,89 ]から$f(t_k+1)\xrightarrow[k\to\infty]{\text{w}\ast}\mu(1)$であり、上と同様の議論から$W_2(f(t_k+1),M)\xrightarrow[k\to\infty]{}W_2(\mu(1),M)$が従う。完備性から$W_2(f(t),M)\xrightarrow[t\to\infty]{}W_2(\mu_0,M)$と収束先は一意に決まるから、結局$W_2(\mu_0,M)=W_2(\mu(1),M)$である。今までの議論で$1$とした部分は任意の$t\in[0,1]$で成り立つから、$W_2(\mu(t),M)$は$t$に依らず一定である。これとStep 7から$\mu=M$が得られ、証明が終わる。

お疲れ様でした。初等的な主張の積み重ねと解析学の基礎的な議論が両方登場するので勉強になりますね。

Boltzmann方程式を"線形化"する

前節でBoltzmann方程式がMaxwell分布に収束していくことがわかりました。しかしその証明はcollison operaterが（特異）積分作用素だったのでやや面倒でしたね。そこで苦労した意味があまりなくなってしまう気もしますが、以下では次のように"線形化"した$Q$を考えましょう、すなわちBoltzmann方程式
$$ \frac{\mathrm{d}}{\mathrm{d}t}\int\varphi \dd f(t)=\int\varphi\dd Q(f,f) $$
の右辺を取り替えて
$$ \frac{\mathrm{d}}{\mathrm{d}t}\int\varphi \dd f(t)=\int\frac{\varphi(v')-\varphi(v)}{\tau}\dd f\tag{6} $$
のようになっているものを考えます。ここで$\tau$は適当な係数です。$v'$はBoltzmann方程式のときは衝突後の粒子の速度を表しました。しかし以下では衝突に限らない適当な粒子の運動則に従って$v$から$v'$が定まるとして、表現の幅を広げておきましょう。このようなある種の"線形化"は例えば[ Mathiaud+,16 ]等によってもう少し複雑なモデルが提案されていますが、今回参照している[ Herty+,20 ]では上記のように非常にシンプルなモデルを考えています。

このような線形化を施しても田中の定理同様にtrend to equilibriumの性質があることがわかっています[ Villani,06 ]。

ResNetの情報伝達をtrend to equilibriumとして理解する

さて、前節の内容をまとめて言い換えると、粒子の運動則が定まるとそれら粒子の分布の発展がBoltzmann方程式(6)が記述でき、そのBoltzmann方程式はtrend to equilibriumという漸近的な性質を満たすのでした。

ここで、粒子というのを素子ないし層と読み替えると、Neural Networkと関連しそうな気がしてきます。Neural Network、特にResidual Neural Networkの概要については Alicia Solidさんの動画がわかりやすいですが、一番簡単で基本的な構成は前の層からの隠れ変数$x$を入力として、出力$x'$を重み$W$、バイアス$b$、そして活性化関数$\sigma$を用いて
$$ x'=x+\epsilon\sigma(Wx+b) $$
で表されます。ここで$x$は位置を表すものではなく、$\epsilon>0$は十分小さいパラメータとして便宜的に導入します。今回重みやバイアスはすべての層で共通としますが、実際は重みやバイアスは確率的勾配法によって学習されます。そこでさらに、上式を確率微分方程式とした
$$ x'=x+\epsilon\sigma(Wx+b)+\sqrt{\epsilon}K(x)\eta\tag{7} $$
をResNetのモデル化として採用することにします。ここで$\eta$は標準正規分布に従う確率変数で、$K(x)$を拡散効果を表す適当な関数とします。
　
そこで式(7)をResNetの粒子が従う運動則だとして、それら粒子が従う分布$f$が従う発展方程式をBoltzmann方程式(5)を用いておおらかに（細かいことを気にせず）求めてみます。すなわち"伊藤の公式"ないしEuler-丸山スキームを用いて$\varphi(x')= \varphi(x)+\epsilon\sigma(Wx+b)\partial_x\varphi(x)+\sqrt{\epsilon}K(x)\eta\partial_x\varphi(x)+\dfrac{\epsilon}{2}K^2(x)\partial_{xx}\varphi(x)$とし、式(6)に$\tau=\epsilon$として代入して$\eta$に関しては期待値を取ることにすると、次のBoltzmann方程式が導出できます。
$$ \frac{\partial f}{\partial t}=\frac{\partial}{\partial x}\left(\sigma(Wx+b)f+\frac{\partial}{\partial x}\left(\frac{1}{2}K^2(x)f\right)\right)\tag{8} $$

このような方程式はFokker-Planck方程式と呼ばれ、やはりtrend to equilibriumを有します。特に今回のResNet(8)の場合は、平衡解$f_\infty$が式(8)から計算でき

$$ f_\infty(x)=\frac{C}{K^2(x)}\exp\left(-\int\frac{2\sigma(Wx+b)}{K^2(x)}\dd x\right) $$

となります。$C$は正規化定数です。ここでさらに簡単のため、隠れ層の次元を$1$、活性化関数をReLU(もどき)$\sigma_I(x)=x$、$K(x)=1$とすると、$f_\infty$はMaxwell分布
$$ f_\infty(x)=C\exp(-Wx^2-2bx) $$
となります。

この結果は何を示唆しているのでしょうか。

equilibrium、つまり最終的なResNetの出力というのは、確率分布を学習したいようなResNetにとっては推定したい分布であって欲しいはずです。その推定したい分布に応じてパラメータ$W,b$を学習していくわけですが、例えば標準正規分布を学習したい場合なら上の計算結果よりパラメータは$W=1/2,b=0$と学習するのがベストです。そこで、ResNetのパラメータを$W=1/2,b=0$とし、適当な入力（例えば一様分布）を入れて各層の様子を観察すると、次図のように推定したい分布に入力が近づいていくことがわかります。