Esercizio - Eliminazione di variabili nella rete burglary (slide del corso)

Soluzione

Fattori. Ogni CPT è una tabella (Reti bayesiane - semantica, indipendenze e inferenza esattaUna rete bayesiana è un grafo diretto aciclico (DAG) con una variabile per nodo e una tabella delle probabilità condizionate (CPT) $P(X_i\mid\text{Parents}(X_i))$ per nodo; rappresenta la congiunta come prodotto $P(x_1,\dots,x_n)=\prod_i P(x_i\mid\text{parents}(X_i))$, con molti meno numeri della tabella completa (rete dell'allarme: 10 contro 31). Tre strutture elementari: catena $X\to Y\to Z$ e fork $Y\leftarrow X\to Z$ (estremi indipendenti dato il centro), collider $X\to Z\leftarrow Y$ (indipendenti, ma dipendenti se si osserva $Z$ o un suo discendente). Modelli causali (genitori = cause) danno reti più piccole. Inferenza esatta per enumerazione: $\mathbf P(X\mid\mathbf e)=\alpha\sum_{\mathbf y}\mathbf P(X,\mathbf e,\mathbf y)$ calcolata come somma di prodotti, in profondità e da sinistra a destra; per l'allarme $\mathbf P(B\mid j,m)=\langle0{,}284;\ 0{,}716\rangle$. Molti prodotti si ripetono: l'eliminazione di variabili li salva come fattori e li riusa.Reti bayesiane - semantica, indipendenze e inferenza esatta →): $f_B(B)$, $f_E(E)$, $f_A(A,B,E)=P(A\mid B,E)$, $f_J(J,A)=P(J\mid A)$, $f_M(M,A)=P(M\mid A)$. L'evidenza $m$ fissa $M$ al valore vero: $f_M$ diventa il fattore $f_m(A)=\langle P(m\mid a),\,P(m\mid\neg a)\rangle=\langle0{,}7;\ 0{,}01\rangle$.

Query 1: $\mathbf P(J\mid m)$. Variabili nascoste: $B,E,A$. Si sceglie l'ordine di eliminazione $B,E,A$.

Passo 1: eliminare $B$. I fattori che contengono $B$ sono $f_B$ e $f_A$. Si moltiplicano e si somma su $B$: $g_1(A,E)=\sum_bf_B(b)\,f_A(A,b,E)$.

• $g_1(a,e)=0{,}001\cdot0{,}95+0{,}999\cdot0{,}29=0{,}00095+0{,}28971=0{,}29066$;
• $g_1(a,\neg e)=0{,}001\cdot0{,}94+0{,}999\cdot0{,}001=0{,}00094+0{,}000999=0{,}001939$;
• $g_1(\neg a,e)=1-0{,}29066=0{,}70934$ e $g_1(\neg a,\neg e)=1-0{,}001939=0{,}998061$ (le colonne sommano a $1$ perché è ancora una probabilità di $A$ dato $E$).

Passo 2: eliminare $E$. Fattori con $E$: $f_E$ e $g_1$. $g_2(A)=\sum_ef_E(e)\,g_1(A,e)$: $$g_2(a)=0{,}002\cdot0{,}29066+0{,}998\cdot0{,}001939=0{,}000581+0{,}001935=0{,}0025164,\qquad g_2(\neg a)=0{,}9974836.$$ Il fattore $g_2$ è la distribuzione marginale dell'allarme: $P(a)=0{,}0025164$.

Passo 3: eliminare $A$. Fattori con $A$: $g_2$, $f_m$ e $f_J$. Il prodotto è una tabella su $J$ e $A$ e la somma su $A$ la riduce a un vettore su $J$: $h(J)=\sum_ag_2(a)\,f_m(a)\,f_J(J,a)$.

• $h(j)=0{,}0025164\cdot0{,}7\cdot0{,}9+0{,}9974836\cdot0{,}01\cdot0{,}05=0{,}0015853+0{,}0004987=0{,}0020841$;
• $h(\neg j)=0{,}0025164\cdot0{,}7\cdot0{,}1+0{,}9974836\cdot0{,}01\cdot0{,}95=0{,}0001761+0{,}0094761=0{,}0096522$.

Passo 4: normalizzare. $h(j)+h(\neg j)=0{,}0117363=P(m)$ e $$\mathbf P(J\mid m)=\langle0{,}0020841;\ 0{,}0096522\rangle/0{,}0117363=\langle0{,}1776;\ 0{,}8224\rangle.$$ Interpretazione: $P(j)=0{,}0521<P(j\mid m)=0{,}1776$: se Maria chiama, è più probabile che l'allarme sia suonato e che chiami anche Giovanni (fork $J\leftarrow A\rightarrow M$: $J$ e $M$ sono dipendenti se non si osserva $A$).

Query 2: quali prodotti si riusano. Nell'enumerazione di $P(b,j,m)=P(b)\sum_eP(e)\sum_aP(a\mid b,e)P(j\mid a)P(m\mid a)$ compaiono, sotto $e$ e sotto $\neg e$ (e di nuovo per $\neg b$), gli stessi prodotti $P(j\mid a)P(m\mid a)=0{,}9\cdot0{,}7=0{,}63$ e $P(j\mid\neg a)P(m\mid\neg a)=0{,}05\cdot0{,}01=0{,}0005$: nell'albero delle slide sono i rami cerchiati in rosa e in verde. L'eliminazione li calcola una sola volta, come fattore $f_J(A)f_M(A)=\langle0{,}63;\ 0{,}0005\rangle$, poi forma la tabella $g(B,E)=\sum_af_A\,f_Jf_M$ (otto numeri: $0{,}598525$, $0{,}59223$, $0{,}183055$, $0{,}0011295$ per $(b,e),(b,\neg e),(\neg b,e),(\neg b,\neg e)$) e infine $h(B)=\sum_ef_E(e)g(B,e)=\langle0{,}5922426;\ 0{,}0014934\rangle$. Moltiplicando per $f_B$: $\langle0{,}0005922;\ 0{,}0014919\rangle\to\langle0{,}284;\ 0{,}716\rangle$, uguale al risultato dell'enumerazione.

Ordine di eliminazione. L'ordine $E,B,A$ dà lo stesso risultato: eliminando prima $E$ si forma $g_1'(A,B)=\sum_ef_E(e)f_A(A,B,e)$ con $g_1'(a,b)=0{,}002\cdot0{,}95+0{,}998\cdot0{,}94=0{,}94002$ e $g_1'(a,\neg b)=0{,}002\cdot0{,}29+0{,}998\cdot0{,}001=0{,}001578$; poi si elimina $B$ e infine $A$, e si ottiene ancora $\langle0{,}1776;\ 0{,}8224\rangle$ (verificato). Cambia il costo, perché in reti più grandi alcuni ordini creano fattori più grandi: il costo dipende dalla dimensione del fattore più grande che si forma, non dal risultato.

Codice

import numpy as np

# fattori della rete burglary; ogni asse e' una variabile, indice 0 = vero, 1 = falso
f_B = np.array([0.001, 0.999])                              # f_B(B)
f_E = np.array([0.002, 0.998])                              # f_E(E)
f_A = np.zeros((2, 2, 2))                                   # f_A(A, B, E) = P(A | B, E)
for (b, e), p in {(0, 0): .95, (0, 1): .94, (1, 0): .29, (1, 1): .001}.items():
    f_A[0, b, e], f_A[1, b, e] = p, 1 - p
f_J = np.array([[0.90, 0.05], [0.10, 0.95]])                # f_J(J, A) = P(J | A)
f_M = np.array([[0.70, 0.01], [0.30, 0.99]])                # f_M(M, A) = P(M | A)

# Query: P(J | m).  Evidenza m: si fissa M = vero, f_M diventa un fattore su A soltanto
f_m = f_M[0, :]                                             # <P(m|a), P(m|~a)> = <0.7, 0.01>

# 1) eliminare B: g1(A, E) = sum_b f_B(b) f_A(A, b, E)
g1 = np.einsum('b,abe->ae', f_B, f_A)
# 2) eliminare E: g2(A) = sum_e f_E(e) g1(A, e)   (e' la marginale P(A))
g2 = np.einsum('e,ae->a', f_E, g1)
# 3) eliminare A: h(J) = sum_a g2(a) f_m(a) f_J(J, a)
h = np.einsum('a,a,ja->j', g2, f_m, f_J)
print("g1\n", g1); print("g2", g2); print("h", h, h / h.sum())

Output: g1 = [[0.29066 0.001939] [0.70934 0.998061]], g2 = [0.0025164 0.9974836], h = [0.0020841 0.0096522] e distribuzione [0.1775766 0.8224234]. np.einsum esegue in un colpo prodotto puntuale e somma: 'b,abe->ae' moltiplica $f_B(b)$ per $f_A(a,b,e)$ e somma sull'indice $b$ che non compare nel risultato. Verificato con l'enumerazione sulla congiunta completa: $\mathbf P(J\mid m)=\langle0{,}1775766;\ 0{,}8224234\rangle$.