AHC 016 - Kiri8128の日記

本稿について

AHC016 Graphorean の考察（主に理論面）について書いています。最尤推定（ ※ 追記：あとで調べたら事前分布を仮定するものはベイズ推定と呼ぶのが正しいようです *1）に類する方法で復元をします。正規分布近似によるグラフを構築した際の誤判定率の計算方法についても記載しています。

残念ながらこれをフルに活用した解法を提出するところまでいけませんでしたが、理論面の参考になればと思います。

方針

本問は

構築フェイズ
復元フェイズ

に分けることができます。

復元フェイズについては最尤推定（後述）が理論的に最善であることが証明できます。

本稿では、厳密な最尤推定の代わりに、特徴量をいくつか用いて連続関数の対数尤度の最大化に帰着する方法を紹介します。

構築フェイズは、（復元フェイズで最尤推定を用いる前提で）期待得点の最大化をすれば良いです。グラフを選択した段階で理論的には期待得点が決まりますが、これを現実的に計算することはできるでしょうか？

実は、特徴量を表す確率変数を独立正規分布で近似すると、構築したふたつのグラフについて誤判定を起こす確率が計算できます *2 。ここから期待得点を簡易的に予測することができます。あとは焼きなまし等で期待得点が最大となるグラフのセットを構築すれば良さそうです。

実はこの方法には問題（後述）があるのですが、途中までは理論として面白いと思うので紹介します。

最尤推定

最尤推定とは

最尤推定とは、本問の例では反転・並び替え後の結果（ $H_k$ ）を見て、 $1$ から $M$ のグラフのうち、それを選んだ場合に $H_k$ になる条件付き確率 *3 が最大のものを選ぶ方法です *4 。

なお復元フェイズでは最尤推定による選択が理論的に最善となることが証明できます *5 。これは、（通常の最尤推定では事前分布が不明であることが多いが）本問では事前分布として $1$ から $M$ までのグラフが等確率で選ばれることが与えられているので、各グラフが正しい確率は尤度に厳密に比例することから従います。

特徴量

実際に、各グラフから条件付き確率を計算することは、グラフの頂点数 $n$ *6 がよほど小さくないと現実的でないです。このためいくつか特徴量を使って判定しました。特徴量としては以下を採用しました。

全体の辺の本数。
各頂点 $i$ を固定したときの、 $i$ と隣接する辺の本数。
各頂点 $i$ を固定したときの、 $(j, k)$ の組であって $i$ から $j$ 、 $i$ から $k$ 、 $j$ から $k$ のいずれにも辺が存在するものの個数
各頂点 $i$ を固定したときの、 $(j, k)$ の組であって $i$ から $j$ 、 $i$ から $k$ 、 $j$ から $k$ のいずれにも辺が存在しないものの個数

最初のものはスカラー、それ以外はいずれも $n$ 元ベクトルになります。

なお、グラフは $n \times n$ 行列で表されているとして考えていた *7 *8 ので、この考え方だと例えば 2. は行ごとの 1 のマスの個数と言い換えられます。

3. と 4. は、 $i$ を含む 3 点からなるクリーク・反クリークの個数とも言い換えられます。

判定の方法

行ごとの 1 のマスの個数の遷移行列は、愚直に DP で prob(from, to) を全パターン計算できます *9。実装上はこれの対数を持っておくとよいです。

3 点クリークは実際に計算するのは難しいので、期待値と分散のみを計算することにしました。あとは中心極限定理を信じて正規分布を仮定すれば良いです。期待値は 3 点を固定して独立に計算してから合計すればよい（期待値の線形性）ので簡単です。分散は少し工夫が必要ですが、「共分散」の概念を知っていれば計算できると思います。

頂点の対応

ひとつ問題なのは、シャッフル前後でどの頂点とどの頂点が対応するかが分からないところです *10。

簡単のため、頂点ごとの特徴量をソートして、大小順にぴったり対応すると考えました。 *11 。

期待値・分散・共分散について（復習）

確率変数 $X$ および $Y$ について、 $E(X)$ を $X$ の期待値、 $V(X)$ を $X$ の分散、 $\mathrm{Cov}(X, Y)$ を $X$ と $Y$ の共分散とします。

共分散は

$\mathrm{Cov}(X, Y) = E(XY) - E(X)E(Y)$

で定義されます。

このとき、次の性質が成立します。

$\mathrm{Cov}(X, X) = V(X)$

$X$ と $Y$ が独立のとき $\mathrm{Cov}(X, Y) = 0$

$V(X+Y) = V(X) + 2\mathrm{Cov}(X, Y) + V(Y)$

$\mathrm{Cov}(X + Y, Z) = \mathrm{Cov}(X, Z) + \mathrm{Cov}(Y, Z)$

$V(aX) = a^2V(X)$

これを繰り返し使うと、多数の確率変数の和である確率変数について、それぞれの分散と任意の 2 つの共分散が分かっていれば、その合計の期待値・分散を求めることができます *12。

各特徴量の期待値・分散の計算

反転率（問題文中の eps）を $p$ 、反転しない確率を $q=1-p$ とします。

1. 全体の本数

全体の辺の本数については、最初の辺の本数を $c$ とすると

平均 $p * c + q * (\displaystyle\frac{n(n-1)}{2} -c)$

分散 $(\displaystyle\frac{n(n-1)}{2})pq$

の二項分布に従います。計算上は正規分布に近似します *13 。

2. 各頂点から出ている辺の本数

各頂点から出ている辺の本数は、最初の該当の本数を $c$ とすると

平均 $pn + q(n - 1 - c)$

分散 $(n-1)pq$

です。

これも二項分布に従うので正規分布で近似できそうですが、少し気になる点があります。それは各頂点から出ている辺の本数を全部合計すると、辺を 2 回ずつカウントしていることになるため、独立性が担保されなくなります。

全く同じものを 2 回合計すると平均は 2 倍、分散は 4 倍になるので、分散を最大 2 倍ぐらい過小評価している可能性があります *14 。

3. 各頂点を含む 3 点クリークの個数

各頂点（ $i$ とする）に隣接する 3 点クリークの個数は、3 辺の組ごとに遷移確率が次の表で定まります。期待値はここから直接計算できます。

from\to	0	1	2	3
0	q^3	3p*q^2	3p^2*q	p^3
1	p*q^2	q^3 + 2p^2*q	2p*q^2 + p^3	p^2*q
2	p^2*q	2p*q^2 + p^3	q^3 + 2p^2*q	p*q^2
3	p^3	3p^2*q	3p*q^2	q^3

分散はちょっとややこしいですが、共分散の考え方を使うと計算できます。つまり、 3 頂点の組すべてについて期待値と分散が求まっていて、 2 つの「3 頂点の組」に対して共分散が求まっていれば全体の期待値と分散が求まります。

分散については、 $i,j,k$ の組が 3 点クリークになる場合に 1 、それ以外の場合に 0 を取る確率変数を $X$ とすると、 $V(X) = E(X)(1 - E(X))$ なので上の表から求まります。

頂点 $i$ を固定するとき、 $i,j,k$ の組と $i,l,m$ の組は相関がないので共分散はゼロです *15 。 $i,j,k$ と $i,j,l$ の組については $i,j$ 間の辺が共通なので相関があります。

上の公式に当てはめて頑張って計算すると共分散を求めることができます。

例えば、 $i,j,k,l$ の 4 点について、任意の 2 点間に最初辺があるとして、

$i,j,k$ の組が 3 点クリークになる場合に 1 、それ以外の場合に 0 を取る確率変数 $X$ と $i,j,l$ の組が 3 点クリークになる場合に 1 、それ以外の場合に 0 を取る確率変数 $Y$ の共分散 $\mathrm{Cov}(X, Y)$ を求めてみます。上の定義を思い出すと

$\mathrm{Cov}(X, Y) = E(XY) - E(X)E(Y) = q^5 - q^3q^3=q^5(1-q)=pq^5$ と計算できます。

なお、上記の $X,\ Y$ および

$i,k,l$ の組が 3 点クリークになる場合に 1 、それ以外の場合に 0 を取る確率変数 $Z$

について、 $\mathrm{Cov}(X, Y)+\mathrm{Cov}(X, Z)+\mathrm{Cov}(Y, Z)$ をまとめて考えると考察・実装がラクになります。

結果を言ってしまうと、この結果は $i,j,k,l$ 間の 6 個の候補にある辺の本数と場所によって、次のように決まります。

6 本の辺がある場合： $3pq^5$

5 本の辺がある場合： $2p^2q^4+pq^5$

4 本の辺がある場合で、それらがサイズ 4 のループをなす場合： $pq^5+2p^3q^3$

4 本の辺がある場合で、それらのうち 3 つがループをなす場合： $p^3q^3+2p^2q^4$

3 本の辺がある場合で、それらが直線状に並んでいる場合： $p^4q^2+p^3q^3+p^2q^4$

3 本の辺がある場合で、それらが長さ 3 のループをなすまたは 1 点から 3 本出ている場合： $3p^3q^3$

2 本の辺がある場合で、それらがつながっている場合： $p^3q^3 + 2p^4q^2$

2 本の辺がある場合で、それらが離れている場合： $p^5q + 2p^3q^3$

1 本の辺がある場合： $2p^4q^2+p^5q$

0 本の辺がある場合： $3p^5q$

反クリークについても同様です。

計算量が気になるかもしれませんが、グラフ盤面を予め固定しておけばパラメトリックな方法で $O(1)$ で計算できる場合もあります。受け取ったグラフのカウントは bit 演算でグラフごとに $O(N^3/w)$ などでできるので十分高速です。

盤面候補の選出

いくつか候補を出して、なるべくお互いが「近くならない」（誤判定する確率が低い）ように選びたいです。

下記の 4 パターンを使いました。

Type 0（左上）：一部のノードのみ全結合 *16

Type 1（左下）：Type 0 の逆 *17

Type 2（右上）：差が一定以内のものを結合 *18

Type 3（右下）：3 段タイプ *19

誤判定率

最尤推定による判定では、各初期グラフ $g=1,2,\cdots, M$ について各特徴量 $i=1,2,\cdots, 3n+1$ による尤度の積を計算して、その合計が最大となる $g$ を選ぶという方法でした。

ここではすべての分布を正規分布で近似して、対数尤度の和の最大化と考えます。正規分布の確率密度関数は $f(x)=\displaystyle\frac{1}{2\pi\sigma ^2}e^{-(x - \mu)^2/2}$ で表されることを思い出しましょう。

平均 $\mu$ 、分散 $\sigma ^2$ の正規分布を $N(\mu, \sigma ^2)$ で表します。確率変数 $X$ がこの分布に従うことを $X \sim N(\mu, \sigma ^2)$ のように書きます。

各特徴量 $i$ の分布を $X_i \sim N(\mu_i, \sigma_i ^2)$ とすると、 $X_i=x$ における対数尤度は $-\displaystyle\frac{(x-\mu_i)^2}{2\sigma_i ^2} - \displaystyle\frac{1}{2}\log 2\pi \sigma_i ^2$ で表されます。復元フェイズではこの和が最大の $i$ を選ぶのが最善です。

次にこれが誤判定を起こす確率を計算 *20 します。

誤判定をするというのは $g$ 番目のグラフについて、 $g$ 番目のグラフの特徴量に基づく対数尤度の合計よりもあるほかのグラフ（ $g'$ とします）の特徴量に基づく対数尤度の合計の方が大きくなることです。まず $g$ と $g'$ を固定して、 $g$ を $g'$ であると誤判定する確率を計算します。

グラフ $g$ の $i$ 番目の特徴量の確率変数を $X_{g,i} \sim N(\mu_{g,i}, \sigma_{g,i} ^2)$ とすると、次の式が成立する確率を求めれば良いです。

$\displaystyle\sum_i \left( -\displaystyle\frac{(X_{g,i}-\mu_{g,i})^2}{2\sigma_{g,i} ^2} - \displaystyle\frac{1}{2}\log 2\pi \sigma_{g,i} ^2 \right) \lt \displaystyle\sum_i \left( -\displaystyle\frac{(X_{g,i}-\mu_{g',i})^2}{2\sigma_{g',i} ^2} - \displaystyle\frac{1}{2}\log 2\pi \sigma_{g',i} ^2 \right)$

$\displaystyle\sum_i \left(\displaystyle\frac{(X_{g,i}-\mu_{g,i})^2}{\sigma_{g,i} ^2} + \log \sigma_{g,i} ^2 \right) \gt \displaystyle\sum_i \left(\displaystyle\frac{(X_{g,i}-\mu_{g',i})^2}{\sigma_{g',i} ^2} + \log \sigma_{g',i} ^2 \right)$

右辺はパラメータが $g$ ではなく $g'$ に基づいていますが、確率変数は $g$ に基づく $X_{g,i}$ であることに注意してください。

ここで $Z_{g,i}= \displaystyle\frac{X_{g,i}-\mu_{g,i}}{\sigma_{g,i}}$ とおくと、これは標準正規分布に従う（ $Z_{g,i} \sim N(0, 1)$ ）ので、 $E(Z_{g,i})=0,\ E(Z_{g,i}^2)=V(Z_{g,i})=1$ が成立します。

これを使って上式をさらに変形すると

$\displaystyle\sum_i \left(\displaystyle\frac{(X_{g,i}-\mu_{g,i})^2}{\sigma_{g,i} ^2} + \log \sigma_{g,i} ^2 \right) \gt \displaystyle\sum_i \left(\Big(\displaystyle\frac{X_{g,i}-\mu_{g,i}}{\sigma_{g,i}}\displaystyle\frac{\sigma_{g,i}}{\sigma_{g',i}} + \displaystyle\frac{\mu_{g,i} - \mu_{g',i}}{\sigma_{g',i}}\Big)^2 + \log \sigma_{g',i} ^2 \right)$

$\displaystyle\sum_i \left(Z_{g,i}^2 + \log \sigma_{g,i} ^2 \right) \gt \displaystyle\sum_i \left(\Big(\displaystyle\frac{\sigma_{g,i}}{\sigma_{g',i}}Z + \displaystyle\frac{\mu_{g,i} - \mu_{g',i}}{\sigma_{g',i}}\Big)^2 + \log \sigma_{g',i} ^2 \right)$

$\displaystyle\sum_i \left(\Big(\displaystyle\frac{\sigma_{g,i} ^2}{\sigma_{g',i} ^2}-1\Big) Z_{g,i}^2 + \displaystyle\frac{2\sigma_{g,i}(\mu_{g,i} - \mu_{g',i})}{\sigma_{g',i} ^2} Z_{g,i} + \Big(\displaystyle\frac{\mu_{g,i}-\mu_{g',i}}{\sigma_{g',i}}\Big)^2+ \log \displaystyle\frac{\sigma_{g',i} ^2}{\sigma_{g,i} ^2} \right) \lt 0$

よってこの左辺を $A$ とおくと、この期待値 $e$ と分散 $\sigma^2$ は

$e=E(A) = \displaystyle\sum_i \left(\Big(\displaystyle\frac{\sigma_{g,i} ^2}{\sigma_{g',i} ^2}-1\Big) + \Big(\displaystyle\frac{\mu_{g,i}-\mu_{g',i}}{\sigma_{g',i}}\Big)^2+ \log \displaystyle\frac{\sigma_{g',i} ^2}{\sigma_{g,i} ^2} \right)$

$\sigma ^2=V(A) = \displaystyle\sum_{i} \left(2\Big(\displaystyle\frac{\sigma_{g,i} ^2}{\sigma_{g',i} ^2}-1\Big)^2 + \displaystyle\frac{4\sigma_{g,i}^2(\mu_{g,i} - \mu_{g',i})^2}{\sigma_{g',i} ^4} \right)$

と計算できます。これも正規分布だと思うと、求める誤判定率は $P(A \lt 0) = 1 - F(\displaystyle\frac{e}{\sigma})$ （ただし $F$ は正規分布の累積分布関数）となります。

問題点

紹介した方法にはいくつか問題点があります。

二項分布などに従う特徴量を正規分布で近似している
誤判定率の計算において、対数尤度の確率変数を正規分布で近似している
行ごとの特徴量をソートして使っている

実験してみると 1. と 2. の問題はあまり大きくないことが分かります。一方、条件によっては 3. の問題がかなり大きいようです。実際、場合によっては本稿の方法で計算した誤判定率がかなり小さい *21 場合でも、正規分布を仮定して実験するとほとんどの場合で誤判定をする *22 こともありました。順序統計量（ソート後の確率変数）が満たす性質に近い分布を選びがちになる場合があるということかなと思います。

ソートして使う以外の方法としては、「〇〇以下である特徴量の個数」のような確率変数を新たに特徴量とすれば、これはソートによる影響を受けないのでもう少し良い判定ができそうです *23 。