Featured image of post Quadratic Pseudo-Boolean Optimization

Quadratic Pseudo-Boolean Optimization

1. はじめに

競プロ界隈では「燃やす埋める」や「Project Selection Problem」などと呼ばれる,最小 $s$-$t$ カットに帰着して解く問題がしばしば登場します.
これらの問題では,最小 $s$-$t$ カットに帰着できることが分かっていても,「どの頂点からどの頂点へ辺を張るか」「どの変数をフリップするか」といったグラフの構成で混乱しがちです.

そこで,問題から最小 $s$-$t$ カットのグラフを直接構成するのではなく,これらの問題を Quadratic Pseudo-Boolean Optimization(以下 QPBO)として統一的に扱うことを考えます.
QPBO として問題を定式化しておけば,そこからグラフを構成する処理をライブラリに任せることができます.また,変数のフリップによって解ける問題についても,フリップを明示的に考える必要がなくなります.
これにより,問題を解く際には「問題をどのように QPBO として定式化するか」に注力できます.

本記事では,まず「単純な関数の場合」として,単純な形の QPBO がそのまま最小 $s$-$t$ カットに帰着できることを確認します.
次に「劣モジュラ関数の場合」として,すべての二変数項が劣モジュラである QPBO は再パラメータ化という処理によってこの単純な形に帰着できることを示します.
さらに「変数のフリップで解ける場合」として,一見すると劣モジュラでない場合でも,変数をフリップすることで劣モジュラにできる場合を考えます.
最後に「一般の関数の場合」で,一般の QPBO に対する QPBO 法を説明します.
実装は QPBO.hpp に記載します.

2. Quadratic Pseudo-Boolean Optimization

QPBO では,$x_p \in \lbrace 0, 1 \rbrace$ として,次の形の関数 $E(\bold x)$ の最小化を目的とします.

$$ \begin{equation} E(\bold x) = \theta_{\mathrm{const}} + \sum_p \theta_p(x_p) + \sum_{p < q} \theta_{pq}(x_p, x_q). \end{equation} $$

ここで,$\theta_{\mathrm{const}}$ は定数項,$\theta_p$ は一つの変数に依存するコスト,$\theta_{pq}$ は二つの変数に依存するコストです.

例えば,変数が $x_1,x_2$ の 2 つの場合を考えます.
各関数値は次の通りだとします.

$\theta_{\mathrm{const}} = 2$
$\theta_{1}(0) = 3,\quad \theta_{1}(1) = 1$
$\theta_{2}(0) = 0,\quad \theta_{2}(1) = 4$

$\begin{array}{c|cc} \theta_{12}(x_1,x_2) & x_2=0 & x_2=1\\ \hline x_1=0 & 0 & 2 \\ x_1=1 & 5 & 1 \end{array}$

このとき,目的関数は

$$ E(x_1,x_2) = 2 + \theta_1(x_1) + \theta_2(x_2) + \theta_{12}(x_1,x_2) $$

です.例えば $(x_1,x_2)=(0,1)$ のとき,

$$ E(0,1)=2+3+4+2=11 $$

となります.

変数が $2$ つなので,全部で $2^2 = 4$ 個の解が考えられます.

$$ \begin{array}{c|c} \ (x_1, x_2) & E(\bold{x}) \\ \hline (0, 0) & 2 + \theta_1(0) + \theta_2(0) + \theta_{12}(0, 0) = 5 \\ (0, 1) & 2 + \theta_1(0) + \theta_2(1) + \theta_{12}(0, 1) = 11 \\ (1, 0) & 2 + \theta_1(1) + \theta_2(0) + \theta_{12}(1, 0) = 8 \\ (1, 1) & 2 + \theta_1(1) + \theta_2(1) + \theta_{12}(1, 1) = 8 \end{array} $$

今回の場合,$x_1 = 0, x_2 = 0$ が最適解で,その目的関数値は $5$ です.
このように QPBO では,各変数の値だけで決まるコスト $\theta_p(x_p)$ と,2 つの変数の組み合わせで決まるコスト $\theta_{pq}(x_p,x_q)$ の和として目的関数を表します.

競プロの問題を解く時は,この形式に帰着できないかを考えます.
例えば,$x_p,x_q\in\lbrace0,1\rbrace$ について,「$x_p$ と $x_q$ の値が異なる場合にコスト $c$ がかかる」という条件は,

$$ \theta_{pq}(0,0)=0,\qquad \theta_{pq}(0,1)=c,\qquad \theta_{pq}(1,0)=c,\qquad \theta_{pq}(1,1)=0 $$

と定義することで表せます.
このように,QPBO を用いると「どのような組み合わせにどれだけのコストを与えるか」を一変数項 $\theta_p$ と二変数項 $\theta_{pq}$ として記述できます.

3. $s$-$t$ カット

QPBO を解くには $s$-$t$ カットを求める必要があるため,$s$-$t$ カットを導入します.
頂点集合 $V$ と有向辺集合 $E$ からなる有向グラフ $G = (V, E)$ が与えられます.辺 $(i, j)$ には容量 $c_{ij} \ge 0$ が定まっているものとします.
頂点集合 $V$ を 2 つの部分集合 $S$ と $T = V \backslash S$ に分割します.2 つの頂点 $s$ と $t$ について $s \in S$,$t \in T$ となるような分割を $s$-$t$ カットと呼びます.
$S$ から出て $T$ に入るような辺の容量の総和を $s$-$t$ カットの容量と呼びます.その値は次の式で定義されます.すべての $s$-$t$ カットのうち最小のものを最小 $s$-$t$ カットと呼びます.

$$ \begin{aligned} c(S) = \sum_{(i, j) \in (S, T)} c_{ij} \end{aligned} $$

下のグラフの $s$-$t$ カットの例をいくつか見ていきます 1.

s-t カットの数値例に用いる有向グラフ


3.1. $s$-$t$ カットの例 1

頂点の部分集合として,$S = \lbrace s, 0, 1 \rbrace$ を選んだとします.
$S$ に属する頂点を赤,$T = V \backslash S$ に属する頂点を青で示します. $S$ から出て $T$ に入るような辺は辺 (0, 2) と辺 (1, 3) です.よって,この $s$-$t$ カットの容量は 3 + 2 = 5 となります.
すべての $s$-$t$ カットの中でこのカットより容量の小さい $s$-$t$ カットは存在しないのでこれは最小 $s$-$t$ カットです.

頂点 0 と 1 を S 側に置いた容量 5 の最小 s-t カット

3.2. $s$-$t$ カットの例 2

頂点の部分集合として $S = \lbrace s, 0, 1, 2, 3 \rbrace$ を選んだとします.
この $s$-$t$ カットの容量は 2 + 3 = 5 となります.
このカットも最小 $s$-$t$ カットです.このように最小 $s$-$t$ カットは複数存在することがあります.

頂点 0、1、2、3 を S 側に置いた容量 5 の最小 s-t カット

3.3. $s$-$t$ カットの例 3

頂点の部分集合として $S = \lbrace s, 1, 2 \rbrace$ を選んだとします.
この $s$-$t$ カットの容量は,3 + 2 + 4 + 2 = 11 となります.
辺 (0, 1) や辺 (0, 2) は $T$ から $S$ に入る辺なので含まれません.

頂点 1 と 2 を S 側に置いた容量 11 の s-t カット


最小 $s$-$t$ カットは最大流問題を解き,残余ネットワーク上で頂点 $s$ から到達できる頂点集合を $S$ とすることで求められます.詳しくは 最大フロー最小カット定理 などを参照してください.
次節から $s$-$t$ カットを使って $E(\bold x)$ を最小化する方法を見ていきます.

4. 単純な関数の場合

ここから,関数 $E(\bold x)$ を最小化する方法を考えていきます.
$x_{p} \in \lbrace 0, 1 \rbrace$ なので,変数の個数が $n$ 個のとき解は $2^n$ 個存在します.
この $2^n$ 個の解のなかから目的関数値を最小にする 0/1 割り当てを見つけることが目標です.

$$ \begin{equation} E(\bold x) = \theta_{\mathrm{const}} + \sum_{p} \theta_{p}(x_p) + \sum_{p \lt q} \theta_{pq}(x_p, x_q) \tag {1} \end{equation} $$

単純な関数の場合を考えたいので,$\theta_{pq}(0,0) = \theta_{pq}(1,1) = 0$ であり,すべての一変数項 $\theta_p$ と二変数項 $\theta_{pq}$ が $0$ 以上の値をとると仮定します.
実は,この仮定を満たす関数の場合は $E(\bold x)$ の解と $s$-$t$ カットの解が 1 対 1 対応するグラフを作成できます.よって,グラフの最小 $s$-$t$ カットがわかれば $E(\bold x)$ の最適解を求めることができます. グラフは各変数を頂点とし,これに特別な頂点 $s$ と $t$ を加えた $n + 2$ 個の頂点から構成されます.辺は下記のルールにしたがって張ります.

関数辺容量
$\theta_{p}(0)$$p \rightarrow t$$\theta_{p}(0)$
$\theta_{p}(1)$$s \rightarrow p$$\theta_{p}(1)$
$\theta_{pq}(0, 1)$$p \rightarrow q$$\theta_{pq}(0, 1)$
$\theta_{pq}(1, 0)$$q \rightarrow p$$\theta_{pq}(1, 0)$

具体例として,変数が $a$ と $b$ の 2 つだけの場合を見てみます.
各変数の値に対応する $E(\bold x) = \theta_a(a) + \theta_b(b) + \theta_{ab}(a, b)$ の値は以下のように定まります.
ただし,$\theta_{ab}(0, 0)$ と $\theta_{ab}(1, 1)$ の値は $0$ であり,$\theta_{\mathrm{const}}$ は定数のため省略しています.

ab$E(\bold x)$
00$\theta_{a}(0) + \theta_{b}(0)$
01$\theta_{a}(0) + \theta_{b}(1) + \theta_{ab}(0, 1)$
10$\theta_{a}(1) + \theta_{b}(0) + \theta_{ab}(1, 0)$
11$\theta_{a}(1) + \theta_{b}(1)$

ルールに従うと下のグラフが構築されます.
このグラフの $s$-$t$ カットをいくつか見ていきます.

単純な関数から構成する 2 変数のグラフ

4.1. $s$-$t$ カットの例 1

$S = \lbrace s, a, b \rbrace$ とします.この $s$-$t$ カットの容量は $\theta_{a}(0) + \theta_{b}(0)$ です.
また,$a = 0$,$b = 0$ としたとき $E(\bold x)$ の値は $\theta_{a}(0) + \theta_{b}(0)$ です.
よって,$S = \lbrace s, a, b \rbrace$ としたときの $s$-$t$ カットの容量と,$a = 0$,$b = 0$ としたときの $E(\bold x)$ の関数値は一致しています.

a と b がともに 0 となる s-t カット

4.2. $s$-$t$ カットの例 2

$S = \lbrace s, a \rbrace$ とします.この $s$-$t$ カットの容量は $\theta_{a}(0) + \theta_{b}(1) + \theta_{ab}(0, 1)$ です.
また,$a = 0$,$b = 1$ としたとき $E(\bold x)$ の値は $\theta_{a}(0) + \theta_{b}(1) + \theta_{ab}(0, 1)$ です.
よって,$S = \lbrace s, a \rbrace$ としたときの $s$-$t$ カットの容量と,$a = 0$,$b = 1$ としたときの $E(\bold x)$ の関数値は一致しています.

a が 0、b が 1 となる s-t カット

4.3. $s$-$t$ カットの全パターン

変数が 2 つの場合は $s$-$t$ カットは $2^2$ 通りあります.すべてのパターンは以下の通りです.

2 変数のグラフにおける s-t カットの全 4 パターン

このように $s$-$t$ カットの構成と各変数への 0/1 の割り当てが 1 対 1 対応するため,最小 $s$-$t$ カットを計算することで $E(\bold x)$ の最適解を求めることができます.
最小 $s$-$t$ カットを計算し,$S$ に属する頂点に対応する変数の値を $0$,$T$ に属する頂点に対応する変数の値を $1$ と設定することで最適な $\bold x$ を構成できます.

5. 劣モジュラ関数の場合

以降は表記を簡潔にするため,$i, j \in \lbrace 0, 1 \rbrace$ に対して,関数値 $\theta_p(i)$ と $\theta_{pq}(i, j)$ をそれぞれ $\theta_{p;i}$ と $\theta_{pq;ij}$ と表記する場合があります.
対応関係は以下の通りです.

関数値以降の表記
$\theta_{p}(0)$$\theta_{p;0}$
$\theta_{p}(1)$$\theta_{p;1}$
$\theta_{pq}(0,0)$$\theta_{pq;00}$
$\theta_{pq}(0,1)$$\theta_{pq;01}$
$\theta_{pq}(1,0)$$\theta_{pq;10}$
$\theta_{pq}(1,1)$$\theta_{pq;11}$

「単純な関数の場合」では,$\theta_{pq;00} = \theta_{pq;11} = 0$ とし,どの関数も $0$ 以上の値を返すことを仮定していました.
この節ではこの仮定を外し,すべての二変数項 $\theta_{pq}$ が劣モジュラであることのみを仮定します.今回は $2$ 値変数を考えているので,$\theta_{pq;01} + \theta_{pq;10} \ge \theta_{pq;00} + \theta_{pq;11}$ を満たすことになります. また,逆に,$\theta_{pq;01} + \theta_{pq;10} \le \theta_{pq;00} + \theta_{pq;11}$ を満たす場合を優モジュラと呼びます.

$\theta_{pq;00}$ や $\theta_{pq;11}$ が $0$ 以外の値をとったり関数値が負の値をとる場合があるので,今回は節 4 で説明したルール通りにグラフを作ることはできません.

5.1. 再パラメータ化

再パラメータ化という操作を行うことでこの問題に対処します.再パラメータ化とは,$E(\bold x)$ の関数値を保ちつつ $\theta_{pq;10}$ などの各関数値を変化させる操作です.

再パラメータ化を行うと標準形とよばれる以下の条件を満たす形になります2.標準形では,定数項を除く各関数が $0$ 以上の値をとります.

  • $\min \lbrace \theta_{pq;00}, \theta_{pq;10} \rbrace = 0$
  • $\min \lbrace \theta_{pq;01}, \theta_{pq;11} \rbrace = 0$
  • $\min \lbrace \theta_{pq;00}, \theta_{pq;01} \rbrace = 0$
  • $\min \lbrace \theta_{pq;10}, \theta_{pq;11} \rbrace = 0$
  • $\min \lbrace \theta_{p;0}, \theta_{p;1} \rbrace = 0$

再パラメータ化をすると,関数 $\theta_{pq}(x_{p}, x_{q})$ が劣モジュラの場合は $\theta_{pq;00} = \theta_{pq;11} = 0$ に,優モジュラの場合は $\theta_{pq;01} = \theta_{pq;10} = 0$ になります.
よって,すべての $\theta_{pq}(x_{p}, x_{q})$ が劣モジュラ関数の場合は,再パラメータ化をすることで「単純な関数の場合」に帰着できます.
再パラメータ化の手続きは以下の通りです.

$\theta_{pq}(x_p,x_q)$ は,次のように $x_p$ を行,$x_q$ を列とする表として考えます.

$$ \begin{array}{c|cc} \theta_{pq}(x_p,x_q) & x_q=0 & x_q=1\\ \hline x_p=0 & \theta_{pq;00} & \theta_{pq;01}\\ x_p=1 & \theta_{pq;10} & \theta_{pq;11} \end{array} $$
  • step 1

    • すべての関数 $\theta_{pq}$ について以下を行う
      • 各列 $j \in \lbrace 0,1 \rbrace$ について
        • $\delta = \min \lbrace \theta_{pq;0j}, \theta_{pq;1j} \rbrace$
        • $\theta_{pq;0j}, \theta_{pq;1j}$ から $\delta$ を引く
        • $\theta_{q;j}$ に $\delta$ を加える
      • 各行 $i \in \lbrace 0,1 \rbrace$ について
        • $\delta = \min \lbrace \theta_{pq;i0}, \theta_{pq;i1} \rbrace$
        • $\theta_{pq;i0}, \theta_{pq;i1}$ から $\delta$ を引く
        • $\theta_{p;i}$ に $\delta$ を加える
  • step 2

    • すべての関数 $\theta_{p}$ について
      • $\delta = \min \lbrace \theta_{p;0}, \theta_{p;1} \rbrace$
      • $\theta_{p;0}, \theta_{p;1}$ から $\delta$ を引く
      • $\theta_{\mathrm{const}}$ に $\delta$ を加える

この再パラメータ化によって目的関数値は変化しません.
step 1 では,例えば列 $j$ から $\delta$ を引くと,$x_q=j$ のときに $\theta_{pq}$ の値が $\delta$ 小さくなりますが,同時に $\theta_{q;j}$ に $\delta$ を加えているので打ち消し合います. 行についても同様です.
step 2 では,$\theta_p$ から引いた $\delta$ を $\theta_{\mathrm{const}}$ に加えているので,やはり目的関数値は変化しません.

5.2. 再パラメータ化の例

具体例として,次の関数を再パラメータ化してみます.

$\theta_{\mathrm{const}} = 0$
$\theta_{p;0}=5,\quad \theta_{p;1}=2$
$\theta_{q;0}=1,\quad \theta_{q;1}=3$

$\begin{array}{c|cc} \theta_{pq}(x_p,x_q) & x_q=0 & x_q=1\\ \hline x_p=0 & 5 & 4\\ x_p=1 & 3 & 1 \end{array}$

この関数は

$$ \theta_{pq;01}+\theta_{pq;10} = 4+3 \ge 5+1 = \theta_{pq;00}+\theta_{pq;11} $$

を満たすので劣モジュラ関数です.
再パラメータ化することにより,目的関数値を保ったまま $\theta_{pq;00}=\theta_{pq;11}=0$ となることを確認します.

5.2.1. step 1 の列 $j=0$

まず,$j=0$ の列について

$$ \delta = \min \lbrace \theta_{pq;00},\theta_{pq;10} \rbrace = \min \lbrace 5,3 \rbrace = 3 $$

とします.

$\theta_{pq;00},\theta_{pq;10}$ から $3$ を引き,$\theta_{q;0}$ に $3$ を加えます.
すると,

$\theta_{\mathrm{const}} = 0$
$\theta_{p;0}=5,\quad \theta_{p;1}=2$
$\theta_{q;0}=4,\quad \theta_{q;1}=3$

$\begin{array}{c|cc} \theta_{pq}(x_p,x_q) & x_q=0 & x_q=1\\ \hline x_p=0 & 2 & 4\\ x_p=1 & 0 & 1 \end{array}$

となります.

5.2.2. step 1 の列 $j=1$

次に,$j=1$ の列について

$$ \delta = \min \lbrace \theta_{pq;01},\theta_{pq;11} \rbrace = \min \lbrace 4,1 \rbrace = 1 $$

とします.

$\theta_{pq;01},\theta_{pq;11}$ から $1$ を引き,$\theta_{q;1}$ に $1$ を加えます.
すると,

$\theta_{\mathrm{const}} = 0$
$\theta_{p;0}=5,\quad \theta_{p;1}=2$
$\theta_{q;0}=4,\quad \theta_{q;1}=4$

$\begin{array}{c|cc} \theta_{pq}(x_p,x_q) & x_q=0 & x_q=1\\ \hline x_p=0 & 2 & 3\\ x_p=1 & 0 & 0 \end{array}$

となります.

続いて,各行について同じ操作を行います.

5.2.3. step 1 の行 $i=0$

$i=0$ の行について

$$ \delta = \min \lbrace \theta_{pq;00},\theta_{pq;01} \rbrace = \min \lbrace 2,3 \rbrace = 2 $$

とします.

$\theta_{pq;00},\theta_{pq;01}$ から $2$ を引き,$\theta_{p;0}$ に $2$ を加えます.
すると,

$\theta_{\mathrm{const}} = 0$
$\theta_{p;0}=7,\quad \theta_{p;1}=2$
$\theta_{q;0}=4,\quad \theta_{q;1}=4$

$\begin{array}{c|cc} \theta_{pq}(x_p,x_q) & x_q=0 & x_q=1\\ \hline x_p=0 & 0 & 1\\ x_p=1 & 0 & 0 \end{array}$

となります.

$i=1$ の行については,$\delta=0$ なので何も変化しません.

5.2.4. step 2

次に一変数項を再パラメータ化します.

まず $\theta_p$ について,

$$ \delta = \min \lbrace \theta_{p;0},\theta_{p;1} \rbrace = \min \lbrace 7,2 \rbrace = 2 $$

なので,$\theta_{p;0},\theta_{p;1}$ から $2$ を引き,$\theta_{\mathrm{const}}$ に $2$ を加えます.

次に $\theta_q$ について,

$$ \delta = \min \lbrace \theta_{q;0},\theta_{q;1} \rbrace = \min \lbrace 4,4 \rbrace = 4 $$

なので,$\theta_{q;0},\theta_{q;1}$ から $4$ を引き,$\theta_{\mathrm{const}}$ に $4$ を加えます.

最終的に再パラメータ化によって以下のようになりました.
$\theta_{pq;00} = \theta_{pq;11} = 0$ となっており,「単純な関数の場合」の形に帰着できていることがわかります.

再パラメータ化前:

$\theta_{\mathrm{const}} = 0$
$\theta_{p;0}=5,\quad \theta_{p;1}=2$
$\theta_{q;0}=1,\quad \theta_{q;1}=3$

$\begin{array}{c|cc} \theta_{pq}(x_p,x_q) & x_q=0 & x_q=1\\ \hline x_p=0 & 5 & 4\\ x_p=1 & 3 & 1 \end{array}$

再パラメータ化後:

$\theta^{\prime}_{\mathrm{const}} = 6$
$\theta^{\prime}_{p;0}=5,\quad \theta^{\prime}_{p;1}=0$
$\theta^{\prime}_{q;0}=0,\quad \theta^{\prime}_{q;1}=0$

$\begin{array}{c|cc} \theta^{\prime}_{pq}(x_p,x_q) & x_q=0 & x_q=1\\ \hline x_p=0 & 0 & 1\\ x_p=1 & 0 & 0 \end{array}$

最後に,再パラメータ化によって目的関数値が変化していないことを確認します.

再パラメータ化前の目的関数は

$$ E(\bold x) = 0 + \theta_p(x_p) + \theta_q(x_q) + \theta_{pq}(x_p,x_q) $$

であり,再パラメータ化後は

$$ E^{\prime}(\bold x) = 6 + \theta_p^{\prime}(x_p) + \theta_q^{\prime}(x_q) + \theta^{\prime}_{pq}(x_p,x_q) $$

です.各 $(x_p,x_q)$ について値を計算すると次のようになります.

$(x_p,x_q)$再パラメータ化前再パラメータ化後
$(0,0)$$0+5+1+5=11$$6+5+0+0=11$
$(0,1)$$0+5+3+4=12$$6+5+0+1=12$
$(1,0)$$0+2+1+3=6$$6+0+0+0=6$
$(1,1)$$0+2+3+1=6$$6+0+0+0=6$

すべての $(x_p,x_q)$ について目的関数値が一致しているので,再パラメータ化によって目的関数値が変化していないことが確認できます.

6. 変数のフリップで解ける場合

前節では,すべての二変数項 $\theta_{pq}$ が劣モジュラである場合について考えました.
この場合,再パラメータ化によって「単純な関数の場合」に帰着できるため,最小 $s$-$t$ カットを用いて $E(\bold x)$ を最小化できます.

この節では,劣モジュラでない二変数項が含まれる場合について考えます.

6.1. 再パラメータ化だけでは解けない場合

二変数項 $\theta_{pq}$ について,

$$ D_{pq} = \theta_{pq;01} + \theta_{pq;10} - \theta_{pq;00} - \theta_{pq;11} $$

を定義します.

$D_{pq} \ge 0$ のとき $\theta_{pq}$ は劣モジュラであり,$D_{pq} \le 0$ のとき優モジュラです. $D_{pq}=0$ の場合,$\theta_{pq}$ は劣モジュラかつ優モジュラであり,変数をフリップしても $D_{pq}$ の符号は変化しません.
また,この二変数項はフリップの選び方に制約を与えません.したがって,以下では $D_{pq} \neq 0$ の二変数項のみを考えます.

前節では再パラメータ化によって劣モジュラ関数を単純な形に変形しましたが,再パラメータ化によって $\theta_{pq}$ が劣モジュラか優モジュラかを変化させることはできません.

例えば,$\theta_{pq}$ の行 $i$ から $\delta$ を引く再パラメータ化を考えます. $i=0$ の場合,

$$ \begin{aligned} D_{pq}^{\prime} &= (\theta_{pq;01}-\delta)+\theta_{pq;10} -(\theta_{pq;00}-\delta)-\theta_{pq;11}\\ &= \theta_{pq;01}+\theta_{pq;10} -\theta_{pq;00}-\theta_{pq;11}\\ &= D_{pq} \end{aligned} $$

となります. 列から値を引く場合についても同様なので,再パラメータ化によって $D_{pq}$ の値は変化しません.
したがって,優モジュラな二変数項を再パラメータ化だけで劣モジュラにすることはできません.

6.2. 変数のフリップ

そこで,変数の $0$ と $1$ の意味を入れ替える操作を考えます. 変数 $x_p$ に対して

$$ x_p^{\prime}=1-x_p $$

とする操作を,$x_p$ のフリップと呼ぶことにします.

二変数項 $\theta_{pq}$ を

$$ \begin{array}{c|cc} \theta_{pq}(x_p,x_q) & x_q=0 & x_q=1\\ \hline x_p=0 & \theta_{pq;00} & \theta_{pq;01}\\ x_p=1 & \theta_{pq;10} & \theta_{pq;11} \end{array} $$

とします.

$x_p$ だけをフリップすると,$x_p=0$ と $x_p=1$ が入れ替わるので,

$$ \begin{array}{c|cc} \theta_{pq}^{\prime}(x_p^{\prime},x_q) & x_q=0 & x_q=1\\ \hline x_p^{\prime}=0 & \theta_{pq;10} & \theta_{pq;11}\\ x_p^{\prime}=1 & \theta_{pq;00} & \theta_{pq;01} \end{array} $$

となります.

このとき,

$$ \begin{aligned} D_{pq}^{\prime} &= \theta_{pq;11}+\theta_{pq;00} -\theta_{pq;10}-\theta_{pq;01}\\ &= -D_{pq} \end{aligned} $$

となります.

したがって,二変数項に含まれる 2 つの変数のうち片方だけをフリップすると,劣モジュラと優モジュラが入れ替わります.

同様に考えると,

  • $x_p, x_q$ のどちらもフリップしない場合,劣モジュラ性は変化しない
  • $x_p, x_q$ の片方だけをフリップする場合,劣モジュラと優モジュラが入れ替わる
  • $x_p, x_q$ の両方をフリップする場合,劣モジュラ性は変化しない

ことがわかります.

よって,優モジュラな二変数項について片方の変数だけをフリップできれば,その二変数項を劣モジュラにできます.

6.3. すべての二変数項を劣モジュラにできる条件

一つの二変数項だけであれば,優モジュラな場合にどちらか一方の変数をフリップすれば劣モジュラにできます. しかし,一つの変数が複数の二変数項に含まれている場合,それぞれの二変数項について独立にフリップすることはできません.

どのようなときにすべての二変数項を同時に劣モジュラにできるかを考えます.

変数 $x_p$ をフリップするかどうかを $f_p \in \lbrace 0,1\rbrace$ で表し,

$$ f_p= \begin{cases} 0 & x_p\text{ をフリップしない}\\ 1 & x_p\text{ をフリップする} \end{cases} $$

とします.

$D_{pq} \gt 0$ の場合,フリップ後も劣モジュラにするには $x_p$ と $x_q$ のフリップの有無を同じにする必要があります. したがって,

$$ f_p \oplus f_q = 0 $$

を満たす必要があります.

一方,$D_{pq} < 0$ の場合,$\theta_{pq}$ は優モジュラなので,片方だけをフリップする必要があります.
したがって,

$$ f_p \oplus f_q = 1 $$

を満たす必要があります.
したがって,すべての二変数項を劣モジュラにできるかどうかは,

$$ f_p \oplus f_q= \begin{cases} 0 & D_{pq} \gt 0 \\ 1 & D_{pq} \lt 0 \end{cases} $$

という制約をすべて同時に満たす $\bold f$ が存在するかどうかに帰着できます.

これはグラフの問題として考えることができます.
各変数を頂点とし,二変数項 $\theta_{pq}$ について頂点 $p, q$ の間に辺を張ります.
$D_{pq} \gt 0$ の二変数項に対応する辺では両端の $f$ を同じ値にし,$D_{pq} \lt 0$ の二変数項に対応する辺では異なる値にします.

連結成分ごとに一つの頂点の $f_p$ を決めれば,DFS や BFS によって隣接する頂点の値を順に決めることができます.
途中で既に決まっている値と矛盾すれば,すべての二変数項を劣モジュラにするフリップは存在しません.

例えば,3 変数 $x_1, x_2, x_3$ について,$\theta_{12}, \theta_{23}, \theta_{13}$ がすべて優モジュラである場合を考えます. このとき,

$$ f_1\oplus f_2=1,\qquad f_2\oplus f_3=1,\qquad f_3\oplus f_1=1 $$

を同時に満たす必要があります. 最初の 2 式から $f_1 = f_3$ となりますが,3 式目では $f_1 \neq f_3$ が必要なので矛盾します.

したがって,この場合はすべての二変数項を同時に劣モジュラにするようなフリップは存在しません.

適切なフリップが存在する場合は,ここまでの方法でそれを求めてから前節の方法で解くことができます.
次節ではこのフリップを明示的に求める必要がなくなる QPBO 法について説明します.

7. 一般の関数の場合

前節では,変数を適切にフリップすることで,すべての二変数項を劣モジュラにできる場合について考えました.
このようなフリップが存在する場合は,フリップ後の問題を「劣モジュラ関数の場合」に帰着することで,最小 $s$-$t$ カットを用いて最適解を求めることができます.

一方,一般の QPBO では,すべての二変数項を劣モジュラにするようなフリップが存在するとは限りません. このような一般の QPBO の最小化は NP-hard であり,常に最小 $s$-$t$ カットに帰着して最適解を求められるわけではありません.

そこで,問題を緩和して劣モジュラな問題に帰着する QPBO 法を考えます. QPBO 法では,必ずしもすべての変数の値を決定できるわけではありませんが,最適解の一部を求めることができます3.

QPBO 法は解として,各変数 $x_p$ に $0,1,\emptyset$ のいずれかを与えます.
$\emptyset$ は $x_p$ の値を決定できなかったことを表します. $x_p$ に $0$ または $1$ が与えられたとき,$x_p$ はラベル付けされたといいます.

QPBO 法には次の性質があります.

  1. アルゴリズムの出力を $\bold x$ とする.完全にラベル付けされた任意の解 $\bold y$ に対して,
$$ z_p = \begin{cases} x_p & x_p \in \lbrace 0,1\rbrace\\ y_p & x_p=\emptyset \end{cases} $$

と定めると,常に $E(\bold z)\le E(\bold y)$ を満たす.

  1. すべての二変数項が劣モジュラの場合は,すべての変数がラベル付けされ,最適解を求めることができる.
  2. QPBO 法は,変数のフリップに対して不変である.
    すなわち,変数をフリップしてから QPBO 法を適用した場合でも,そのフリップに合わせてラベルの $0/1$ を読み替えれば,元の問題に対応する結果が得られる.

性質 1 で $\bold y$ として最適解を選ぶと,QPBO 法によって $0$ または $1$ が与えられた変数については,その値を保つ最適解が存在することがわかります.
また,性質 2 と 3 から,前節で考えた「変数を適切にフリップすることですべての二変数項を劣モジュラにできる場合」には,実際にそのフリップを求めなくても QPBO 法によってすべての変数をラベル付けできます.
したがって,QPBO 法を利用すれば,問題を解く側でどの変数をフリップするかを考える必要がありません.

7.1. QPBO 法

QPBO 法を説明します.
あらかじめ $E(\bold x)$ を標準形に再パラメータ化します.以降の $\theta$ は再パラメータ化後の係数を表すものとします.

まず $x_{\bar{p}} = 1 - x_{p}$ を導入します.$x_{\bar{p}}$ は $x_{p}$ が決まれば一意に定まります.
$x_p$ と $x_{\bar p}$ を使って各項を 2 通りに表し,その平均をとる形に変形します.
この変形の目的は,各項を「単純な関数の場合」で扱った,2 つの変数の値が異なるときだけコストが発生する形に分解することです.
$E(\bold x) = \theta_{\mathrm{const}} + \sum \theta_{p}(x_p) + \sum \theta_{pq}(x_p, x_q)$ は次のように変形できます.

$$ \begin{alignedat}{2} E(\bold x) &= \theta_{\mathrm{const}} + \sum \theta_{p}(x_p) &&+ \sum \theta_{pq}(x_p, x_q) \\ &= \theta_{\mathrm{const}} + \sum \big( \theta_{p;1} x_{p} + \theta_{p;0}(1 - x_{p}) \big) \\ &\quad &&+ \sum \big( \theta_{pq;00} (1 - x_{p})(1 - x_{q}) \\ &\quad &&\quad + \theta_{pq;01} (1 - x_{p}) x_{q} \\ &\quad &&\quad + \theta_{pq;10} x_{p}(1 - x_{q}) \\ &\quad &&\quad + \theta_{pq;11} x_{p} x_{q} \big) \\ &= \theta_{\mathrm{const}} + \sum \bigg( \frac{\theta_{p;1}}{2}(x_p + (1 - x_{\bar{p}})) && + \frac{\theta_{p;0}}{2}(x_{\bar{p}} + (1 - x_p)) \bigg) \\ &\quad &&+ \sum \bigg( \frac{\theta_{pq;00}}{2} \big(x_{\bar{p}} (1 - x_q) + (1 - x_p) x_{\bar{q}} \big) \\ &\quad &&\quad + \frac{\theta_{pq;01}}{2} \big((1 - x_p) x_q + x_{\bar{p}} (1 - x_{\bar{q}}) \big) \\ &\quad &&\quad + \frac{\theta_{pq;11}}{2} \big(x_p (1 - x_{\bar{q}}) + (1 - x_{\bar{p}}) x_q \big) \\ &\quad &&\quad + \frac{\theta_{pq;10}}{2} \big(x_p (1 - x_q) + (1 - x_{\bar{p}}) x_{\bar{q}} \big) \bigg) \end{alignedat} $$

ここで,$x_{\bar{p}} = 1 - x_p$ という制約を緩和し,$x_p$ と $x_{\bar{p}}$ が独立に値をとれる緩和問題を考えます.
同じ変数の組に依存する項をまとめると下の関数に分割できることがわかります.
標準形では定数項を除く各係数が $0$ 以上なので,分割して得られる各関数も非負の劣モジュラ関数になります.したがって「劣モジュラ関数の場合」に帰着できます.

$x_p$
0$\frac{1}{2}\theta_{p;0}$
1$\frac{1}{2}\theta_{p;1}$
$x_{\bar{p}}$
0$\frac{1}{2}\theta_{p;1}$
1$\frac{1}{2}\theta_{p;0}$
$x_p$$x_q$
000
01$\frac{1}{2}\theta_{pq;01}$
10$\frac{1}{2}\theta_{pq;10}$
110
$x_p$$x_{\bar{q}}$
000
01$\frac{1}{2}\theta_{pq;00}$
10$\frac{1}{2}\theta_{pq;11}$
110
$x_{\bar{p}}$$x_q$
000
01$\frac{1}{2}\theta_{pq;11}$
10$\frac{1}{2}\theta_{pq;00}$
110
$x_{\bar{p}}$$x_{\bar{q}}$
000
01$\frac{1}{2}\theta_{pq;10}$
10$\frac{1}{2}\theta_{pq;01}$
110

上記関数について,ルールに従ってグラフを構築します. 整理すると以下のルールに従ってグラフを構築すればいいことがわかります.

$\theta$辺容量4
$\theta_{p;0}$$(p \rightarrow t), (s \rightarrow \bar p)$$\frac{1}{2} \theta_{p;0}$
$\theta_{p;1}$$(s \rightarrow p), (\bar p \rightarrow t)$$\frac{1}{2} \theta_{p;1}$
$\theta_{pq;01}$$(p \rightarrow q), (\bar q \rightarrow \bar p)$$\frac{1}{2} \theta_{pq;01}$
$\theta_{pq;10}$$(q \rightarrow p), (\bar p \rightarrow \bar q)$$\frac{1}{2} \theta_{pq;10}$
$\theta_{pq;00}$$(p \rightarrow \bar q), (q \rightarrow \bar p)$$\frac{1}{2} \theta_{pq;00}$
$\theta_{pq;11}$$(\bar q \rightarrow p), (\bar p \rightarrow q)$$\frac{1}{2} \theta_{pq;11}$

このグラフの最小 $s$-$t$ カットを計算します.
元の制約 $x_{\bar p}=1-x_p$ は緩和しているので,頂点 $p$ と $\bar p$ が異なる側に分かれた場合だけ $x_p$ の値を決定します.
最小 $s$-$t$ カットが複数ある場合は,できるだけ多くの変数の値が定まるものを選びます.この選び方により,すべての二変数項が劣モジュラの場合は全変数の値が定まります.
選び方については参考文献の “Choosing a minimum cut” を参照してください.
よって,$\bold x$ は次のように構成されます.

$$ x_{p} = \left\{ \begin{array}{ll} 0 & \text{if} \space p \in S, \bar p \in T \\ 1 & \text{if} \space p \in T, \bar p \in S \\ \emptyset & \text{otherwise} \end{array} \right. $$

7.2. $s$-$t$ カットとラベルの対応例

具体例として,変数が $a$ と $b$ の 2 つだけの場合を考えます.
各変数の値に対応する $E(\bold x) = \theta_a(a) + \theta_b(b) + \theta_{ab}(a, b)$ の値は以下のように定まります.

ab$E(\bold x)$
00$\theta_{a}(0) + \theta_{b}(0) + \theta_{ab}(0, 0)$
01$\theta_{a}(0) + \theta_{b}(1) + \theta_{ab}(0, 1)$
10$\theta_{a}(1) + \theta_{b}(0) + \theta_{ab}(1, 0)$
11$\theta_{a}(1) + \theta_{b}(1) + \theta_{ab}(1, 1)$

対応するグラフは以下のようになります.表記が煩雑になるので図では $\frac{1}{2}$ を除外しています.
このグラフの(最小とは限らない) $s$-$t$ カットの例をいくつか見ていきます.

QPBO 法で構成する 2 変数のグラフ

$S = \lbrace s, a, b \rbrace$ とします.この $s$-$t$ カットの容量は $\frac{1}{2} (\theta_{a}(0) + \theta_{a}(0) + \theta_{b}(0) + \theta_{b}(0) + \theta_{ab}(0, 0) + \theta_{ab}(0, 0))$ です.
この値は $a = 0$,$b = 0$ としたときの $E(\bold x)$ の目的関数値と一致します.

QPBO 法で a と b がともに 0 とラベル付けされるカット

$S = \lbrace s, a, \bar b \rbrace$ とします.この $s$-$t$ カットの容量は $\frac{1}{2} (\theta_{a}(0) + \theta_{a}(0) + \theta_{b}(1) + \theta_{b}(1) + \theta_{ab}(0, 1) + \theta_{ab}(0, 1))$ です.
この値は $a = 0$,$b = 1$ としたときの $E(\bold x)$ の目的関数値と一致します.

QPBO 法で a が 0、b が 1 とラベル付けされるカット

$S = \lbrace s, a, b, \bar{b} \rbrace$ とします.この場合,$a = 0$,$b = \emptyset$ とし,$b$ のラベルは未定となります.

QPBO 法で a が 0、b が未定となるカット

8. 問題

ここからは,問題から最小 $s$-$t$ カットのグラフを直接構成するのではなく,QPBO として定式化することで競プロの問題を解いていきます.

8.1. ARC085 E - MUL

宝石が $N$ 個あり,それぞれ $1,2,\cdots,N$ と数が書かれています。
あなたは,以下の操作を好きなだけ行うことが出来ます(一度も行わなくてもよいです)。

  • 正整数 $x$ を選ぶ。$x$ の倍数が書かれた宝石を全て叩き割る。

そして,$i$ が書かれていた宝石が割られずに残っていた場合,$a_i$ 円貰います。 ただし,この $a_i$ は負の場合もあり,その場合はお金を払わなくてはいけません。
うまく操作を行った時,あなたは最大で何円お金を貰えるでしょうか?

まず変数を定義します.
宝石 $i$ が残っているかどうかを $x_i$ で表します.宝石を割る場合 $1$ を,残す場合は $0$ をとります.

次に関数を定義します.
QPBO は目的関数値の最小化を目指すのでコストがいくらかかるかで表します.
宝石 $i$ が残っている場合 $a_i$ 円貰えます.これは $-a_i$ 円のコストを払うということなので,次のように定義できます.

  • $\theta_{i}(0) = -a_i$
  • $\theta_{i}(1) = 0$

また,宝石 $i$ を割るにもかかわらず $i$ で割り切れる値が書かれた宝石 $j$ を残すことは許されないので,この場合は無限のコストがかかるとします.よって,次のように定義できます.

  • $\theta_{ij}(0, 0) = 0$
  • $\theta_{ij}(0, 1) = 0$
  • $\theta_{ij}(1, 0) = \infty$
  • $\theta_{ij}(1, 1) = 0$

この関数は $\theta_{ij}(0, 1) + \theta_{ij}(1, 0) \ge \theta_{ij}(0, 0) + \theta_{ij}(1, 1)$ を満たしているので劣モジュラ関数です.
あとは,すべての $i$ と $i$ で割り切れる $j$ について上記関数を定義すれば問題を解けます.
求める答えは $-\min E(\bold x)$ です.

提出コード

8.2. ABC193 F - Zebraness

縦 $N$ マス、横 $N$ マスのマス目があります。上から $i$ 行目、左から $j$ 列目のマスをマス $(i,j)$ と表すことにします。 マス $(i,j)$ の色の情報が文字 $c_{i,j}$ により与えられます。
$B$ はマスが黒で塗られていることを、 $W$ はマスが白で塗られていることを、 $?$ はマスにまだ色が塗られていないことを表します。
高橋くんは、まだ色が塗られていないマスをそれぞれ黒または白で塗り、白黒のマス目を作ります。マス目のしまうま度を、辺で接する黒マスと白マスの組の個数と定義します。高橋くんが達成できるしまうま度の最大値を求めてください。

まず変数を定義します.
添字を簡潔にするため,マス $(i, j)$ を $p = i \times N + j$ で表します. マス $p$ の色を変数 $x_{p}$ で表します.白の場合 $0$ をとり,黒の場合 $1$ をとります.

次に関数を定義します.
与えられている色の変更はできないので白から黒や黒から白に変更すると無限のコストがかかるとします.次のように定義できます.

  • マス $p$ の色が黒の場合

    • $\theta_{p}(0) = \infty$
    • $\theta_{p}(1) = 0$
  • マス $p$ の色が白の場合

    • $\theta_{p}(0) = 0$
    • $\theta_{p}(1) = \infty$
  • マス $p$ が未着色の場合

    • $\theta_p(0) = 0$
    • $\theta_p(1) = 0$

マス $p$ と辺で接するマス $q$ が異なる色だと -1 のコストがかかります.

  • $\theta_{pq}(0, 0) = 0$
  • $\theta_{pq}(0, 1) = -1$
  • $\theta_{pq}(1, 0) = -1$
  • $\theta_{pq}(1, 1) = 0$

この二変数項は優モジュラですが,$(i + j)$ が奇数のマスだけ変数をフリップすると,隣接する 2 マスのうち必ず片方だけがフリップされるため,すべての二変数項を劣モジュラにできます. したがって節 6 で説明した条件を満たしています.
また,QPBO 法は変数のフリップに対して不変なので実際にこのフリップを行う必要はなく,上記の関数をそのまま定義すれば問題を解けます.
求める答えは $-\min E(\bold x)$ です.

提出コード

8.3. その他の問題

9. 参考


  1. この数値例は最大フロー最小カット定理から引用しています ↩︎

  2. 標準形は一意に定まるとは限りません ↩︎

  3. Minimizing non-submodular functions with graph cuts – a review ↩︎

  4. 実装では容量に $\frac{1}{2}$ をかけるのではなく,容量を 2 倍したグラフを構築しています.全変数の値が定まった場合,最小 $s$-$t$ カット値を $C$ とすると目的関数値は $\theta_{\mathrm{const}} + \frac{C}{2}$ です. ↩︎

Hugo で構築されています。
テーマ Stack は Jimmy によって設計されています。