Esercizio - Epidemia di colera di Broad Street (laboratorio 5)

Soluzione

Query: $\mathbf P(Y\mid\neg x)$ (presa a monte: $\neg x$). Le variabili nascoste sono $Z$ e $W$. Come nell'inferenza per enumerazione (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 →): $$\mathbf P(Y\mid\neg x)=\alpha\,\mathbf P(Y,\neg x)=\alpha\sum_{z,w}\mathbf P(Y,w,\neg x,z)=\alpha\sum_z\sum_w\mathbf P(Y\mid w,z)P(w\mid\neg x,z)P(\neg x)P(z)=\alpha\,P(\neg x)\sum_zP(z)\sum_w\mathbf P(Y\mid w,z)P(w\mid\neg x,z).$$ $P(\neg x)$ esce dalle due somme e $P(z)$ dalla somma su $w$.

Passo 1: somma su $W$, per ciascun valore di $Z$. Il vettore ha le componenti $\langle y,\neg y\rangle$.

• $z$: $P(w\mid\neg x,z)=0{,}10$, $P(\neg w\mid\neg x,z)=0{,}90$: $$\langle0{,}80;\,0{,}20\rangle\cdot0{,}10+\langle0{,}15;\,0{,}85\rangle\cdot0{,}90=\langle0{,}08;\,0{,}02\rangle+\langle0{,}135;\,0{,}765\rangle=\langle0{,}215;\ 0{,}785\rangle.$$
• $\neg z$: $P(w\mid\neg x,\neg z)=0{,}02$, $P(\neg w\mid\neg x,\neg z)=0{,}98$: $$\langle0{,}75;\,0{,}25\rangle\cdot0{,}02+\langle0{,}05;\,0{,}95\rangle\cdot0{,}98=\langle0{,}015;\,0{,}005\rangle+\langle0{,}049;\,0{,}931\rangle=\langle0{,}064;\ 0{,}936\rangle.$$

Passo 2: somma su $Z$, pesando con $P(z)=0{,}25$ e $P(\neg z)=0{,}75$: $$0{,}25\cdot\langle0{,}215;\,0{,}785\rangle+0{,}75\cdot\langle0{,}064;\,0{,}936\rangle=\langle0{,}05375;\,0{,}19625\rangle+\langle0{,}048;\,0{,}702\rangle=\langle0{,}10175;\ 0{,}89825\rangle.$$

Passo 3: moltiplicare per $P(\neg x)=0{,}5$ e normalizzare. $0{,}5\cdot\langle0{,}10175;\,0{,}89825\rangle=\langle0{,}050875;\,0{,}449125\rangle$, somma $0{,}5$ (è già $P(\neg x)$): dividendo si ottiene $$\mathbf P(Y\mid\neg x)=\langle0{,}10175;\ 0{,}89825\rangle.$$ Con presa a monte la probabilità di colera è circa $10{,}2\%$. Per confronto, con presa a valle ($x$) lo stesso calcolo dà $\mathbf P(Y\mid x)=\langle0{,}6675;\ 0{,}3325\rangle$: il $66{,}8\%$ (con probabilità finte, ma il confronto mostra l'effetto della posizione della presa attraverso l'acqua contaminata).

Codice

import numpy as np

t, f = 0, 1                                        # indici di 'vero' e 'falso'
P_x = np.array([0.5, 0.5])                         # P(X)
P_z = np.array([0.25, 0.75])                       # P(Z)
# P_w_xz[x][z] = distribuzione di W: <P(w|x,z), P(non w|x,z)>
P_w_xz = np.array([[[0.90, 0.10], [0.85, 0.15]],
                   [[0.10, 0.90], [0.02, 0.98]]])
# P_y_wz[w][z] = distribuzione di Y: <P(y|w,z), P(non y|w,z)>
P_y_wz = np.array([[[0.80, 0.20], [0.75, 0.25]],
                   [[0.15, 0.85], [0.05, 0.95]]])

def P_Y_dato_x(x):
    """P(Y | X=x) = alfa P(x) sum_z P(z) sum_w P(Y|w,z) P(w|x,z)."""
    totale = np.zeros(2)
    for z in (t, f):
        interno = np.zeros(2)
        for w in (t, f):
            interno += P_y_wz[w][z] * P_w_xz[x][z][w]       # somma su W
        totale += P_z[z] * interno                          # somma su Z
    totale *= P_x[x]
    return totale / totale.sum()                            # normalizzazione

print(P_Y_dato_x(f))      # presa a monte (non x)
print(P_Y_dato_x(t))      # presa a valle (x)

Output: [0.10175 0.89825] e [0.6675 0.3325]. Verificato sommando la congiunta completa di tutte le $16$ combinazioni di $X,Z,W,Y$ (somma totale $1$). I valori intermedi [0.215 0.785] ($z$) e [0.064 0.936] ($\neg z$) coincidono con il calcolo a mano. Il codice del corso fa gli stessi passi, con indicizzazioni diverse (P_y_wz[t][:][f]), e ottiene lo stesso [0.10175 0.89825].