離散カントロビッチ問題のエントロピー正則化とその双対問題

離散カントロビッチ問題とは？

次のような問題を考えます

商品を $n$ 個の工場から, $m$ 個のお店に輸送することを考えます.
商品の輸送に掛かる労力が「重量」と「始点と終点の距離」の積で決まっているとき,
トータルの労力が最小になるような運び方はどんなものでしょうか？

数学的に定式化すると

それぞれの工場での生産量を $a \in \mathbb{R}_+^n$ として表し, それぞれのお店での需要量を $b \in \mathbb{R}_+^m$ で表します.
また, 行列 $C = (c_{ij})_{ij}$ により, 単位体積を工場 $i$ からお店 $j$ へ輸送したときの距離（労力）を表します.

さらに商品全体の総重量を $1$ で正規化し,

\begin{align*}
\sum_{i}a_i = 1,\quad \sum_{j}b_j = 1.
\end{align*}

と仮定して一般性を失いません*1.

これに対して, 工場 $i$ からお店 $j$ への $nm$ 通りの輸送量を決めていきます.
工場 $i$ からお店 $j$ への輸送量を未知変数を $p_{ij}$ とすると, $i,j$ 間の輸送コストは

\begin{equation}
c_{ij}p_{ij}
\end{equation}

となります.

さらに, $P = (p_{ij})_{ij}$ と確率行列を用いて表すことにすると, 全工場から全お店への総輸送コストは

\begin{equation}
C(P) = \langle C,P \rangle = \sum_{ij}c_{ij}p_{ij}
\end{equation}

と表せることがわかります.

こうして, 離散カントロビッチ問題は, 以下の最小化問題により定式化されます*2:

\begin{equation}\tag{P}\label{P}
\mathop{\rm minimize\,}_{P\in H} \langle C,P \rangle
\end{equation}

したがって, この最小化問題は線形計画問題となっています.
ここで, 実行可能領域 $H$ は

\begin{align*}
H &= H_1 \cap H_2 \cap H_3,\\
H_1 &= \left\{p \in \mathbb{R}^{n\times m} \, |\, p_{ij} \ge 0\right\},\\
H_2 &= \left\{p \in \mathbb{R}^{n\times m} \, |\, P 1_m = a \right\},\\
H_3 &= \left\{p \in \mathbb{R}^{n\times m} \, |\, P^T1_n = b\right\}.
\end{align*}

により記述されます.

1つ目の条件は, 輸送量はマイナス値を取らないこと
2つ目の条件は, それぞれの輸送元 $i$ について, 送り出した商品の総量が $a_i$ であること
3つ目の条件は, それぞれの輸送先 $j$ について, 受けとった商品の総量が $b_j$ であること

をそれぞれ表します.

双対問題を導出しよう

ここで, $H_2, H_3$ を等式制約と見なすと,
対応する Lagrange 関数 $L: H_1 \times \mathbb{R}^n\times \mathbb{R}^m\rightarrow \mathbb{R}$ を,

\begin{equation}
L(P,f,g)
= \langle{C},{P}\rangle - \langle{f},{P1_m -a}\rangle - \langle{g},{P^T1_n- b}\rangle
\end{equation}

と定めることにより, 問題 (\ref{P}) を

\begin{align}\tag{P}
\mathop{\rm minimize\,}_{P\in H_1} \max_{(f,g)\in \mathbb{R}^n\times \mathbb{R}^m}L(P,f,g)
\end{align}

と書き直すことができます*3.
以後, 最小化問題と呼んだときには, こちらの表現を念頭に置いて考察しましょう.

したがってその双対問題は, min と max をひっくり返して計算することにより,

\begin{align}\tag{Q}\label{Q}
&\mathop{\rm maximize\,}_{(f,g)\in \mathbb{R}^n\times \mathbb{R}^m}
\min_{P\in H_1}L(P,f,g)\\
&\quad = \mathop{\rm maximize\,}_{(f,g)\in \mathbb{R}^n\times \mathbb{R}^m} \left\{
\langle{f},{a}\rangle +\langle{g},{b}\rangle -\min_{P\in H_1} \langle{C-f1_m^T-1_ng^T},{P}\rangle\right\}\\
&\quad = \mathop{\rm maximize\,}_{(f,g)\in R(C)} \left\{
\langle{f},{a}\rangle +\langle{g},{b}\rangle\right\}, \quad R(C) = \left\{(f,g)\in \mathbb{R}^n\times \mathbb{R}^m| f_i + g_j \le c_{ij} \right\}
\end{align}

となります.
すなわち, 制約条件付き最大化問題として表すことができました *4.

また, 主問題が線形計画問題なので双対問題との双対ギャップは発生しません.

エントロピー 正則化

エントロピー関数を

\begin{equation}
H(P) = -\sum_{i,j} p_{ij}\left(\log{(p_{ij})}-1\right)
\end{equation}

により定義します*5.

さらにパラメータ $\varepsilon >0$ を固定します.

最小化問題 \ref{P} におけるコスト関数にエントロピー関数 $H(P)$ を加えて,

\begin{equation}\tag{P$_\varepsilon$}\label{P_eps}
\mathop{\rm minimize\,}_{P\in H'} \left\{\langle{C},{P}\rangle -\varepsilon H(P)\right\}
\end{equation}

により最小化問題を緩和することをエントロピー正則化と呼びます*6.

このとき, $H' =H_2 \cap H_3$ で, 不等号制約 $P \in H_1$ は自動的に外せます.
なぜならば $P\not\in H_1$ のとき, 不等号 $p_{ij}<0$ が成り立つことにより $-\varepsilon H(P) = +\infty$ であり, 最小化問題の解の候補になりえないからです.

さらに, 目的関数が線形関数と強凸関数の和なので強凸計画問題に変形*7されたことがわかります.

対応する Lagrange 関数 $\mathcal{E}: \mathbb{R}^{n\times m} \times \mathbb{R}^n\times \mathbb{R}^m\rightarrow \mathbb{R}$ は,

\begin{equation}
\mathcal{E}(P,f,g)
= \langle{C},{P}\rangle - \varepsilon H(P)- \langle{f},{P1_m -a}\rangle - \langle{g},{P^T1_n- b}\rangle
\end{equation}

と定めることができます.

このとき, 双対問題は

\begin{align}\tag{Q$_\varepsilon$}\label{Q_eps}
&\mathop{\rm maximize\,}_{(f,g)\in \mathbb{R}^n\times \mathbb{R}^m}
\min_{P\in \mathbb{R}^{n\times m}}\mathcal{E}(P,f,g)\\
&\quad = \mathop{\rm maximize\,}_{(f,g)\in \mathbb{R}^n\times\mathbb{R}^m} \left\{
\langle{f},{a}\rangle +\langle{g},{b}\rangle - \varepsilon\langle{\exp{(f/\varepsilon)}},{K\exp{(g/\varepsilon)}}\rangle\right\}
\end{align}

であり, こちらは制約条件のない最大化問題となっています.

ただしここで, $K_{ij} = \exp{(-c_{ij}/\varepsilon)}$ であり, 式変形には停留条件:

\begin{align*}
&\dfrac{\partial \mathcal{E}(P,f,g)}{\partial p_{ij}} = c_{ij} +\varepsilon \log(p_{ij}) -f_i -g_j = 0\\
&\Leftrightarrow p_{ij} = \exp{\left(\frac{f_i+g_j - c_{ij}}{\varepsilon}\right)}
\end{align*}

を用いました.
エントロピー関数の形状が凸であり, また境界近傍の傾きが $\log$ により発散していることから, 微分により停留条件を求めることができています*8 *9.