例：血液型問題

血液型は、両親から遺伝した2つの遺伝子のペアから決定します。例えば、遺伝子がAA,AO,OAならば血液型はA型になります。遺伝子A,B,Oが出現する確率をそれぞれ$p_a, p_b, p_o$とします。（つまり$p_a+p_b+p_o = 1$)。そして、両親がランダムな血液型を持っているとしたときの子の血液型がA,B,O,ABになる確率をそれぞれ$P_A, P_B, P_O, P_{AB}$とします。つまり、

\begin{aligned}
P_A &= p_a^2 + 2p_ap_o \\
P_B &= p_b^2 + 2p_bp_o \\
P_O &= p_o^2 \\
P_{AB} &= 2p_ap_b 
\end{aligned}

次の二つの問題を考えてみましょう

(問題1)$p_a = p_b = p_c$の時、$P_A, P_B, P_O, P_{AB}$は何か？
(問題2)A型が4人、B型が2人、O型が3人、AB型が1人だったというデータがあった時、$p_a, p_b, p_o$は何である可能性が高いか？

一つ目は式に当てはめたら終わりそうですが、二つ目は方程式を解かないと行けなさそうですね。そこで、どちらも計算できる確率論理プログラミングを触ってみました。（初心者です）

Problog

インストール方法

OS Xだとうまく行かなかったのでUbuntuで動かしました。 残念です

pip install problog

問題1

コード

blood.pl

0.333::parent(a, X); 0.333::parent(b, X); 0.334::parent(o, X).

genotype(X, Y) :- parent(X, 1), parent(Y, 2).

bloodtype(P) :- genotype(X, Y), (
  X = Y, P = X;
  X = o, P = Y;
  Y = o, P = X;
  X = a, Y = b, P = ab;
  X = b, Y = a, P = ab
).

query(bloodtype(X)).

実行方法

problog blood.pl

実行結果

 bloodtype(a):	0.333333  
bloodtype(ab):	0.221778  
 bloodtype(b):	0.333333  
 bloodtype(o):	0.111556

0.333::parent(a, X)はparent(a, X)が任意のXで確率0.333で真であるという意味です。また、;でつなげることで、それらのアトムが排他的であるということ(つまり$p_a + p_b + p_c = 1$であること)を示しています。query/1は確率を計算したいアトムの対象です

問題2

bloodln.pl

t(_)::parent(a, X); t(_)::parent(b, X); t(_)::parent(o, X).

genotype(X, Y) :- parent(X, 1), parent(Y, 2).

bloodtype(P) :- genotype(X, Y), (
  X = Y, P = X;
  X = o, P = Y;
  Y = o, P = X;
  X = a, Y = b, P = ab;
  X = b, Y = a, P = ab
).

query(bloodtype(X)).

evidence.pl

evidence(bloodtype(a), true).
-----
evidence(bloodtype(a), true).
-----
evidence(bloodtype(a), true).
-----
evidence(bloodtype(a), true).
-----
evidence(bloodtype(b), true).
-----
evidence(bloodtype(b), true).
-----
evidence(bloodtype(o), true).
-----
evidence(bloodtype(o), true).
-----
evidence(bloodtype(o), true).
-----
evidence(bloodtype(ab), true).

t(_)::parent(a, X)がparent(a, X)が学習対象であるという意味です。また、evidence.plは例のリストで、-----で区切られています。

実行方法

problog lfi blood.pl evidence.pl

実行結果

-12.79944937914457 [0.29231361, 0.16301681, 0.54466957] [t(_)::parent(a,X), t(_)::parent(b,X), t(_)::parent(o,X)] 21

$p_a = 0.292, p_b= 0.163, p_o = 0.545$という推論結果が出ました。

PRISM

インストール方法

wget https://rjida.meijo-u.ac.jp/prism/download/prism23_macx.tar.gz   
tar zxvf prism23_macx.tar.gz

問題1

コード

PRISMをダウンロードしたフォルダに次のblood.psmを作成してください。

blood.psm

values(gene, [a, b, o]).

genotype(X, Y) :- msw(gene, X), msw(gene, Y).

bloodtype(P) :- genotype(X, Y), (
  X = Y -> P = X;
  X = o -> P = Y;
  Y = o -> P = X;
  P = ab
).

述語values(gene, [a, b, o])はリスト[a, b, o]の要素が選ばれる確率の和が1となる空間geneを宣言しています。述語msw(gene, X)は、その空間geneからその要素Xをその確率に沿って選ぶ述語です（つまり、msw(gene, a)となる確率は$p_a$になります)。また、その確率はデフォルトでは一様です。

実行方法

ダウンロードしたフォルダで、./prism/bin/prismを実行してください。（あるいはパスを通してprismを実行してください。）

次を実行して、blood.psmを読み込みます。

| ?- prism(blood).
loading::blood.psm.out

yes

あるアトムの確率を求めるにはprob/1を利用します。前述の通り、まだデータをいれていないので、$p_a = p_b = p_o$であるという前提で計算されます。

| ?- prob(bloodtype(a)).
Probability of bloodtype(a) is: 0.333333333333333

yes
| ?- prob(bloodtype(b)).
Probability of bloodtype(b) is: 0.333333333333333

yes
| ?- prob(bloodtype(o)).
Probability of bloodtype(o) is: 0.111111111111111

yes
| ?- prob(bloodtype(ab)).
Probability of bloodtype(ab) is: 0.222222222222222

yes

問題2

なんと、全く同じファイルblood.psmで動きます。

サンプルを入れるにはクエリで、learn/1とcount/2を利用します。

| ?- learn([count(bloodtype(a), 4), count(bloodtype(b), 2), 
            count(bloodtype(o), 3), count(bloodtype(ab), 1)]). 
#goals: 0(4)
Exporting switch information to the EM routine ... done
#em-iters: 0(6) (Converged: -128.004797160)
Statistics on learning:
        Graph size: 27
        Number of switches: 1
        Number of switch instances: 3
        Number of iterations: 6
        Final log likelihood: -128.004797160
        Total learning time: 0.000 seconds
        Explanation search time: 0.000 seconds
        Total table space used: 4736 bytes
Type show_sw to show the probability distributions.

show_sw/0で推論結果を表示します。

| ?- show_sw.
Switch gene: unfixed_p: a (p: 0.292330954) b (p: 0.163020599) o (p: 0.544648448)

$p_a = 0.292, p_b= 0.163, p_o = 0.545$という推論結果が出ました。

こいつら何が違うの？

Problogの公式サイト曰く、PRISMは

特定の問題形式に限定することで、計算量を少なくしている
確率的なメモ化をしていない
否定に特別な扱いが必要である

うーん私にはよくわからないですし、

私の感想では、problogは複雑なベイジアンネットワークを構成するには簡単だと感じましたが、今回取り上げた問題のように例から元の確率を推論する問題に関してはPRISMの方が書きやすかったように感じます。

参考

Prolog公式サイト
 PRISM公式サイト