Probabilistic Graphical Models – Semana 8

Este post apresenta as notas de aula da semana 8 do curso online Probabilistic Graphical Models.

Aprendizado estruturado

Na aula passada vimos maneiras de aprender as distribuições de probabilidade a partir de dados, onde a estrutura da rede já era dada. Agora, estamos interessados em aprender a estrutura da rede.

Os método que discutiremos a seguir são baseados em pontuação (score). Definimos uma função de pontuação que avalia quão bem a estrutura representa com os dados. O objetivo é encontrar uma estrutura que maximiza essa função.

Pontuação usando verossimilhança

Podemos definir uma função de pontuação denotada por \mbox{score}_L({\cal G}: {\cal D}), como a log-verossimilhança,

\mbox{score}_L({\cal G}: {\cal D}) = \ell((\hat \theta, {\cal G}): {\cal D})

Onde \hat \theta representa a MLE para \ell({\cal G}: {\cal D})

Definimos a distribuição de probabilidade empírica, denotada por \hat P(x, y) como a frequência de um dado resultado no conjunto de dados {\cal D} = \{x[m], y[m]\}_{m=1}^M, como

\hat P(x, y) = \dfrac{M[x,y]}{M}

Onde M[x,y] é o número de instâncias em {\cal D} tal que x[m]=x, y[m]=y. Para o caso da distribuição multinomial, temos que $\hat P$ é a MLE.

Definimos agora a informação mútua de uma distribuição \hat P, denotada por \hat I_{\hat P}(X; Y), como

\hat I_{\hat P}(X; Y) = \sum_{x, y} \hat P(x, y) \log \dfrac{\hat P(x,y)}{\hat P(x) \hat P(y)}

Temos que \hat I_{\hat P}(X; Y) \ge 0, sendo 0 se e somente se X e Y são independentes.

Definimos a entropia de uma distribuição \hat P, denotada por \hat H_P(X), como

H_P(X) = - \sum_{x} P(x) \log P(x)

Podemos então escrever a função de pontuação em função dessas novas definições

\mbox{score}_L({\cal G} : {\cal D}) = M\sum_{i=1}^n I_{\hat P}(X_i; Pa_{X_i}^{\cal G}) - M\sum_{i}H_{\hat P}(X_i)

Onde Pa_{X_i}^{\cal G} são os pais de X_i na rede \cal G.

Uma limitação dessa função é que para a maioria dos casos, a rede que a maximiza é o grafo completo. Para evitar problemas com overfitting, costuma-se inserir restrições como limitar o número de pais que um vértice pode ter ou o número de parâmetros da rede.

Há também a alternativas mais flexíveis como penalizar redes complexas.

Pontuação BIC e consistência assintótica

Uma maneira para penalizar a pontuação usando log-verossimilhança é através da pontuação BIC (Bayesian Information Criterion), que é dada por

\mbox{score}_{BIC}({\cal G}: {\cal D}) = \ell(\hat \theta_{\cal G}: {\cal D}) - \dfrac{\log M}{2}\mbox{Dim}[{\cal G}]

O último termo é o fator de penalização. \mbox{Dim}[{\cal G}] é equivalente ao número de parâmetros independentes de \cal G.

Essa penalização é um tradeoff entre o fitting dos dados e a complexidade do modelo.

Comportamento assintótico. Detalhando os termos da pontuação BIC teremos

= M\sum_{i=1}^n I_{\hat P}(X_i; Pa_{X_i}^{\cal G}) - M\sum_{i}H_{\hat P}(X_i) - \dfrac{\log M}{2}\mbox{Dim}[{\cal G}]

Note que a entropia \hat H não depende de G. Além disso, a informação mútua cresce linearmente com o número de amostras, enquanto a penalidade cresce logaritmicamente. Conforme M aumenta, maior peso é dado ao fitting e menos à simplicidade do modelo.

Consistência. Seja G^* a rede que representa verdadeiramente as amostras. Pode se mostrar que se M tende a infinito, $G^*$ maximiza essa função de pontuação.

Intuitivamente, arestas sobressalentes não contribuem com a informação mútua, mas aumentam a complexidade do modelo e portanto são penalizadas. Por outro lado, arestas necessárias são adicionadas por causa do crescimento linear comparado com a penalização logarítmica.

Pontuação Bayesiana

Essa função de pontuação é derivada do teorema de Bayes, que no nosso contexto é dado por

P({\cal G} | {\cal D}) = \dfrac{P({\cal D}|{\cal G}) P({\cal G})}{P({\cal D})}

Notando que a distribuição marginal dos dados P({\cal D}) é independente de \cal G, podemos escrever a pontuação Bayesiana como o logarítmo do numerador acima, resultando em

(1) \mbox{score}_B({\cal G} : {\cal D}) = \log P({\cal D}|{\cal G}) + \log(P({\cal G}))

Vamos analisar cada um dos termos com mais detalhes. Primeiramente, definimos a verossimilhança marginal P({\cal D}|{\cal G}) como

(2) P({\cal D}|{\cal G}) = \int P({\cal D}|{\cal G},\theta_{\cal G}) P(\theta_{\cal G}|{\cal G})d\theta_{\cal G}

O primeiro termo dessa integral é a verossimilhança, enquanto o segundo termo é a distribuição a priori dos parâmetros.

Para o caso de redes Bayesianas multinomiais, a expressão acima pode ser reescrita de uma maneira mais fácil de se computar. Para simplificar, vamos definir

\displaystyle A_i^{\cal G} = \prod_{u_i \in Val(Pa_{X_i}^{\cal G})} \dfrac{\Gamma(\alpha_{X_i|u_i})}{\Gamma(\alpha_{X_i|u_i} + M[u_i])}

onde Val(Pa_{X_i}^{\cal G}) é o conjunto de todas as atribuições possíveis dos pais de X_i, \alpha_{X_i|u_i} são os parâmetros da distribuição de Dirichlet, M[u_i] é a contagem das instâncias com o valor u_i e \Gamma é a função gamma, que é uma generalização da função de Fibonacci para o domínio contínuo.

Também definimos

\displaystyle B_i^{\cal G} = \prod_{x_i^j \in Val(X_i)} \dfrac{\Gamma(\alpha_{x_i^j|u_i} + M[x_i^j, u_i])}{\Gamma(\alpha_{x_i^j|u_i})}

Onde Val(X_i) é o conjunto de todos os possíveis valores de X_i. Com isso, podemos escrever a verossimilhança marginal como

P({\cal D}|{\cal G}) = \prod_i (A_i^{\cal G} B_i^{\cal G})

Aplicando o logaritmo na expressão acima, podemos decompô-la em famílias de pontuação,

\log P({\cal D}|{\cal G}) = \sum_i \mbox{FamScore}_B(X_i|Pa_{X_i}^{\cal G} : {\cal D})

Analisada a verossimilhança marginal, vamos prosseguir com o segundo termo de (1), P({\cal G}), denominada estrutura a priori.

Podemos escolher a distribuição inicial de \cal G como sendo uniforme ou penalizando a complexidade da rede como por exemplo o número de arestas ou o número de parâmetros.

A constante normalizadora da estrutura a priori é a mesma para qualquer rede. Como aplicamos o logaritmo, essa constante é somada na função de pontuação independentemente da rede e portanto podemos ignorá-la.

Parâmetros a priori. Voltando a falar sobre a verossimilhança marginal, temos a distribuição a priori dos parâmetros, P(\theta_{\cal G} | {\cal G}). Da mesma forma que na estimativa de parâmetros, vamos utilizar a distribuição de Dirichlet (BDe prior).

Lembrando, temos o tamanho do conjunto equivalente, denotado por \alpha, que representa o número de amostras que foram hipoteticamente observadas. Temos também a rede B_0, contendo as distribuições a priori.

O número de amostras hipoteticamente observadas para cada combinação de valores é dada por \alpha[x_i, Pa_{x_i}^{\cal G}] = \alpha P(x_i, Pa_{x_i}^{\cal G} | B_0). Note que não necessariamente Pa_{x_i}^{\cal G} = Pa_{x_i}^{B_0}

A vantagem da rede Bayesiana a priori B_0 é que podemos usá-la para prover as distribuições a priori para qualquer rede \cal G.

O uso dessa distribuição para a distribuição a priori dos parâmetros tem a propriedade de que rede I-equivalentes possuem a mesma pontuação Bayesiana.

Além disso, conforme o número de instâncias amostradas M tende a infinito, uma rede com distribuição a priori de Dirichlet satisfaz a pontuação de BIC mais uma constante

\mbox{score}_{BIC}({\cal G} : {\cal D}) + c

Aprendizado de árvores

Agora que definimos algumas possibilidades de se avaliar a qualidade de uma rede através de funções de pontuação, vamos focar no problema de encontrar a rede que maximiza uma dada função de pontuação, que é essencialmente um problema de otimização.

Inicialmente, vamos nos restringir a redes que são árvores. Para esse caso, temos que cada vértice tem no máximo um pai. Seja p(i) > 0 o pai de X_i ou p(i) = 0 se X_i não possui pai. Assim, podemos re-escrever uma dada função de pontuação como

\mbox{score}({\cal G} : {\cal D}) = \sum_{i} \mbox{score}(X_i|Pa_{X_i}^{\cal G} : {\cal D})

\displaystyle = \sum_{i:p(i)>0} \mbox{score}(X_i | X_{p(i)} : {\cal D}) + \sum_{i:p(i)=0} \mbox{score}(X_i : {\cal D})

Somando \mbox{score}(X_i : {\cal D}) para i : p(i) > 0 no segundo somatório e subtraindo do primeiro temos

\displaystyle = \sum_{i:p(i)>0} (\mbox{score}(X_i | X_{p(i)} : {\cal D}) - \mbox{score}(X_i : {\cal D})) +

\displaystyle \sum_{i} \mbox{score}(X_i : {\cal D})

Agora o segundo somatório é independente da estrutura da árvore. Vamos definir pesos para as arestas

w(i \rightarrow j) = \mbox{score}(X_i | X_{p(i)} : {\cal D}) - \mbox{score}(X_i : {\cal D})

O objetivo é encontrar o conjunto de arestas que maximiza a soma desses pesos.

Se a função de pontuação for a verossimilhança, temos que w(i \rightarrow j)  = M I_{\hat P} (X_i ; X_j) e todos os pesos são não negativos e portanto a solução ótima é uma árvore (ou seja, um grafo conexo).

Para as pontuações BIC e BDe, a penalização pode fazer com que os pesos nas arestas sejam negativos. Dessa forma, a solução ótima é uma floresta (ou seja, não necessariamente conexo).

Dizemos que uma função de pontuação satisfaz a equivalência de pontuação se estruturas I-equivalentes possuem a mesma pontuação. Esse é o caso das pontuações de verossimilhança, BIC e BDe. Nesse caso, podemos mostrar que w(i \rightarrow j) = w(j \rightarrow i). Portanto, para fins de computar a pontuação de uma árvore, podemos ignorar o sentido das arestas.

Com isso, podemos reduzir esse problema para o da árvore geradora mínima, que pode ser resolvido através de algoritmos como o Kruskal e o Prim.

A entrada é um grafo completo com custo dado por w(i \rightarrow j). Uma vez computada a árvore geradora mínima desse grafo, removemos as arestas de custo negativo dessa árvore, obtendo uma floresta (conexa ou não). Essa floresta representa a estrutura que maximiza a função de pontuação.

Aprendizado em grafos gerais

Vimos que é possível determinar a estrutura que maximiza as funções de pontuação baseadas em verossimilhança, BIC e BDe de maneira exata se nos restringirmos a florestas. Entretanto, para grafos gerais, o problema é NP-difícil, como enuncia o tereoma baixo.

Teorema. Encontrar uma estrutura de pontuação máxima onde os vértices têm no máximo k pais é NP-difícil para k > 1.

Podemos então usar heurísticas tais como a busca local. Alguns tipos de vizinhança comumente utilizados são remoção de aresta, adição de aresta e inversão do sentido da aresta.

Como solução inicial podemos usar uma rede vazia (sem arestas), a melhor árvore, uma rede aleatória ou alguma outra estrutura construída a partir de conhecimento a priori.

Aqui também temos os conhecidos problemas de busca local como ótimos locais e platôs. Aliás, para pontuações que respeitam a equivalência de pontuação o problema de platôs é comum porque estruturas I-equivalentes são geralmente vizinhas com relação à reversão de arestas.

Variantes que visam contornar esses efeitos são a busca tabu e o GRASP, por exemplo.

Complexidade

Vamos analisar a complexidade de um passo da busca local. Há n(n-1) arestas que podem sofrer remoção ou inversão de sentido ou podem ser adicionadas, levando a um número O(n^2) de vizinhos. A complexidade total é dada por O(n^3 M).

Podemos usar a propriedade da “decomponibilidade” da função de pontuação para melhorar a complexidade da busca local. Lembrando, temos que,

\mbox{score}({\cal G} : {\cal D}) = \sum_i \mbox{score}(X_i | Pa_{X_i}^{\cal G} : {\cal D})

É possível mostrar que os únicos componentes do somatório que são modificados entre duas soluções vizinhas são correspondes aos vértices incidentes à aresta alterada. Portanto, a ideia é manter um conjunto de deltas representando a diferença entre o custo de uma solução vizinha e a solução atual. Cada vértice tem O(n) deltas que correspondem às possíveis arestas incidentes a ele. Portanto temos que manter O(n^2) deltas.

Ao realizarmos um movimento, no máximo O(n) deltas devem ser atualizados, correspondendo aos vértices incidentes à aresta alterada. Para atualizar um delta, precisamos re-avaliar todas as amostras, ficando com a complexidade O(M).

O conjunto de deltas pode ser mantido em uma fila de prioridade, de modo que podemos escolher o próximo movimento (o de maior delta) em O(1). Os O(n) deltas recomputados devem ser atualizados na fila, o que pode ser feito em O(n\log n).

A complexidade total de um movimento explorando essas ideias é O(nM + n\log n)

Leitura adicional

[1] Wikipedia – Bayesian Information Criterion
[2] Wikipedia – Minimum Description Length
[3] Wolfram – Gamma Function

Os comentários estão fechados.

%d bloggers like this: