Probabilistic Graphical Models – Semana 9

Maio 21, 2012

Este post apresenta as notas de aula da semana 9, a última do curso online Probabilistic Graphical Models.

Aprendizado com dados incompletos

Visão geral

Na prática as amostras com as quais trabalhamos muitas vezes são incompletas, podendo ser oriundas de falhas na aquisição dos dados ou mesmo a presença de variáveis latentes (ocultas) no modelo.

Relembrando, utilizamos variáveis latentes no modelo para simplificá-lo. Considere o exemplo abaixo:

Figura 1: Uso de variáveis latentes simplifica o modelo

No exemplo podemos ver que a introdução da variável H reduz o número de arestas no modelo e consequentemente o número de parâmetros.

Lidando com dados incompletos

Podemos considerar dois casos de falta de informação. Se ela ocorre de maneira independente do valor observado, então basta ignoramos as amostras com dados incompletos. Por outro lado, se existe uma dependência entre a falta de informação e o valor observado, então não podemos simplesmente ignorar amostras defeituosas. Nesse caso precisamos considerar o mecanismo que modela a perda desses dados.

Perda ao acaso ou MAR (Missing At Random) é quando a probabilidade da ausência de um dado é independente de seu valor. Um exemplo simples é quando obtemos amostras de lançamentos de uma moeda. Se em ocasiões aleatórias deixamos de registrar o resultado de um lançamento, então temos um MAR.

Por outro lado, se por algum motivo só deixamos de anotar lançamentos que resultaram em cara, então a probabilidade de ausência do dado é influenciada pelo valor da amostra e nesse caso não temos um MAR.

Com a presença de dados restantes, a função de verossimilhança possui múltiplos ótimos locais. Uma intuição da razão disso é que substituindo o valor de uma amostra está faltando, o valor da verossimilhança continua o mesmo, embora o ponto no espaço que gera esse valor seja diferente. Essa característica é denominada identificabilidade.

Métodos para a otimização de verossimilhança

Uma vez que a função de verossimilhança pode ter vários ótimos locais, podemos usar métodos baseados em gradientes para estimar parâmetros que resultam em bons valores de verossimilhança.

Subida de gradiente

Podemos obter o gradiente usando o resultado do seguinte teorema:

\displaystyle \dfrac{\partial \log P(D| \Theta)}{\partial \theta_{x_i, u_i}} = \dfrac{1}{\theta_{x_i|u_i}} \sum_m P(x_i, u_i | d[m], \Theta)

Onde \theta_{x_i, u_i} representa o parâmetro para a atribuição X_i = x_i e Y_i = y_i no CPD de X_i. Temos que P(x_i, u_i | D, \Theta) representa a probabilidade da atribuição X_i = x_i e Y_i = y_i dadas as amostras e os parâmetros \Theta.

Para cada cálculo do gradiente, é preciso computar essa probabilidade para cada atribuição e para cada amostra d[m]. Usando inferência em uma árvore-clique apenas uma vez para um dado d[m], conseguimos obter P(x_i, u_i | D, \Theta).

As vantagens desse método é que ele pode ser generalizado para CPDs que não são na forma tabular. As desvantagens é que é preciso garantir que os parâmetros encontrados são de fato CPDs e para acelerar a convergência pode ser preciso combinar com métodos mais avançados.

Maximização de Esperança

Como alternativa ao método usando gradientes, existe um método específico para maximizar a verossimilhança que é chamado de maximização de esperança ou EM (Expectation Maximization).

A ideia por trás do método é que estimar parâmetros com dados completos é fácil. Por outro lado, computar a distribuição dos dados incompletos é fácil (via inferência) dados os parâmetros.

O algoritmo pode ser descrito pelos seguintes passos:

  • Escolha um ponto inicial para os parâmetros
  • Itere
    • Passo-E (Esperança): completar os dados restantes fazendo amostragem com os parâmetros atuais
    • Passo-M (Maximização): estimar novos parâmetros usando os dados completos

Vamos detalhar mais os Passo-E e Passo-M:

Esperança (Passo-E). Para cada conjunto de dados d[m] calcule P(X,U | d[m], \theta^t) para todas as atribuições de (X|U). Compute as estatísticas suficientes esperadas ou ESS (Expected Sufficient Statistics) para cada atribuição x, u, dadas por

\displaystyle \bar M_{\theta^t} [x, u] = \sum_{m = 1}^M P(x, u| d[m], \theta^t)

Seja \bar x[m] (\bar u[m]) o conjunto de variáveis de x (u) que estão faltando em d[m]. Podemos concluir que

\displaystyle P(x, u| d[m], \theta^t) = P(\bar x[m], \bar u[m] | d[m], \theta)

Maximização (Passo-M). Calcule a MLE usando as estatísticas suficientes esperadas como se fossem as estatísticas suficientes reais, ou seja,

\theta_{x|u}^{t+1} = \dfrac{\bar M_{\theta^t}[x,u]}{\bar M_{\theta^t}[u]}

A vantagem deste método é que podemos aproveitar algoritmos para calcular a MLE e na prática este método obtém bons resultados rapidamente, nas primeiras iterações. A desvantagem é que nas iterações finais a convergência é lenta.

É possível mostrar que a verossimilhança L(\theta: D) melhora a cada iteração nesse método.

Aspectos práticos da maximização da esperança

Algumas outras considerações práticas do algoritmo EM:

  1. Convergência da verossimilhança não necessariamente implica na convergência dos parâmetros
  2. Executar muitas iterações do algoritmo sobre um conjunto de treinamento pode levar a overfitting
  3. O chute inicial de parâmetros pode afetar bastante a qualidade das soluções obtidas

Referências

[1] Wikipedia – Expectation Maximization

Comentários Desativados em Probabilistic Graphical Models – Semana 9 | inteligência artificial, stanford | Permalink
Publicado por kunigami

Probabilistic Graphical Models – Semana 8

Maio 13, 2012

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

Comentários Desativados em Probabilistic Graphical Models – Semana 8 | inteligência artificial, stanford | Permalink
Publicado por kunigami

Probabilistic Graphical Models – Semana 7

Maio 6, 2012

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

Estimativa de parâmetros em Redes Bayesianas

Estimativa de verossimilhança máxima

Dada uma amostra {\cal D} = \{x[1], \cdots, x[M] \} obtida de maneira independente de uma distribuição desconhecida, queremos estimar \theta que represente \cal D.

Uma medida da qualidade dessa distribuição é a verossimilhança, denotada por L(\theta:{\cal D}) e definida como a probabilidade de se obter \cal D dada a distribuição \theta, ou seja

L(\theta:{\cal D}) = P({\cal D}|\theta) = \prod_{i=1}^{M} P(x[i]|\theta)

Definimos estatística suficiente como uma estatística que provê a mesma informação de uma amostra, sendo suficiente para calcularmos a verossimilhança. Vamos defini-la em termos de uma função que recebe uma amostra x e retorne uma tupla \mathbb{R}^k.

Um exemplo é a estatística suficiente para uma distribuição multinomial de k valores. Podemos definir a função s(x_i) = (0, \cdots, 0, 1, 0, \cdots, 0) que recebe um valor x_i e retorna uma tupla apenas com a i-ésima posição contendo 1 e o resto 0.

Agora calculamos a verossimilhança definindo \bar M = (M_1, \cdots, M_k) = \sum_{x[i] \in {\cal D}} s(x[i]). Note que M_k representa o número de amostras com o valor k. Logo, a verossimilhança é dada por

L(\bar{\theta}:{\cal D}) = \prod_{i=1}^{k} \theta_i^{M_i}

Outro exemplo é a distribuição Gaussiana com média \mu e variância \sigma^2. A distribuição pode ser dada por

p(X) = \dfrac{1}{\sqrt{2\pi \sigma^2}} \exp \left(-x^2 \dfrac{1}{2\sigma^2} + x \dfrac{\mu}{\sigma^2} - \dfrac{\mu^2}{\sigma^2}\right)

Podemos definir a estatística suficiente como s(x) = (1, x, x^2). Seja

\bar M = (M_1, M_2, M_3) = \sum_{x[i] \in {\cal D}} s(x[i])

podemos escrever a função de verossimilhança como

L(\theta: {\cal D}) = \left(\dfrac{1}{\sqrt{2\pi \sigma^2}}\right)^k \exp \left(-M_3 \dfrac{1}{2\sigma^2} + M_2 \dfrac{\mu}{\sigma^2} - M_1 \dfrac{\mu^2}{\sigma^2}\right)

O problema da estimativa de verossimilhança máxima (Maximum Likelihood Estimation, MLE) consiste em, dado um conjunto de amostras \cal D, encontrar a distribuição \theta que maximiza a verossimilhança L(\theta: D).

No caso da distribuição multinomial, podemos mostrar que o MLE é

\theta^*_i = \frac{M_i}{\sum_{i=1}^m M_i}

No caso da distribuição Gaussiana, temos

\mu^* = \dfrac{1}{M}

\sigma^* = \sqrt{\dfrac{1}{M} \sum_m (x[m] - \mu^*)^2}

Estimativa de verossimilhança máxima em Redes Bayesianas

Podemos generalizar a estimativa de verossimilhança para as distribuições em Redes Bayesianas. Suponha que tenhamos uma rede Bayesiana e queremos estimar as distribuições de probabilidades, nesse caso, os CPDs.

Novamente, dado um conjunto de amostras, \cal D, definimos a verossimilhança

L(\Theta: D) = \prod_{m} P(x[m]:\Theta)

Podemos fatorar essa distribuição em um produto dos CPDs de cada nó i usando a regra da cadeia.

= \prod_{m} \prod_i P(x_i[m] | \pmb{u}_i[m] : \theta_i)

Neste caso \pmb{u}_i representa os pais da variável i na rede. Re-arranjando os fatores, podemos trocar a ordem dos produtórios e obter,

= \prod_{i} \prod_m P(x_i[m] | \pmb{u}_i[m] : \theta_i)

Definindo a verossimilhança local como

L_i({\cal D}: \theta_i) = \prod_m P(x_i[m] | \pmb{u}_i[m] : \theta_i)

Temos finalmente que

L(\Theta: D) = \prod_{i} L_i(D:\theta_i)

Supondo inicialmente que os CPDs não são compartilhados por nenhuma variável, temos um \theta_i para cada variável. Vamos mostrar como determinar esses CPDs de maneira independente.

L_i({\cal D}: \theta_i) = \prod_{m} P(x[m]|\pmb{u}_[m] : \theta) = \prod_{m} P(x[m]|\pmb{u}_[m] : \theta_{X|\pmb{U}})

Re-arranjando os produtos para agrupá-los por valores de x e \pmb{u} (buckets), temos

\displaystyle = \prod_{x, \pmb{u}} \prod_{\substack{m:x[m] = x\\\pmb{u}[m] = \pmb{u}}} P(x[m]|\pmb{u}_[m] : \theta_{X|\pmb{U}})

Observamos que P(x|\pmb{u}) corresponde a uma entrada do CPD dado por \theta_{X|\pmb{U}}, que denotamos por \theta_{x|\pmb{u}}, portanto,

\displaystyle = \prod_{x, \pmb{u}} \prod_{\substack{m:x[m] = x\\\pmb{u}[m] = \pmb{u}}} \theta_{x|\pmb{u}}

Note que o termo \theta_{x|\pmb{u}} será multiplicado pelo número de amostras tais que x[m] = x e \pmb{u}[m] = \pmb{u}, que denotamos por M[x,\pmb{u}],

\displaystyle = \prod_{x, \pmb{u}} \theta_{x|\pmb{u}}^{M[x,\pmb{u}]}

Finalmente, vemos que a distribuição P(X|\pmb{U} = \pmb{u}) corresponde à distribuição multinomial, então podemos calcular o MLE \theta^*_{x|\pmb{u}} como

\theta^*_{x|\pmb{u}} = \dfrac{M[x,\pmb{u}]}{\sum_{x'} M[x',\pmb{u}]} = \dfrac{M[x,\pmb{u}]}{M[\pmb{u}]}

CPDs compartilhados. Note que a fórmula que obtivemos depende apenas dos valores de x e \pmb{u} e do CPD \theta_{x|\pmb{u}}. Assim, se tivermos CPDs compartilhados entre variáveis, como no modelo temporal, a única diferença é que podemos contar as ocorrências em M[x,\pmb{u}] múltiplas vezes para uma mesma amostra.

Observamos finalmente que o número de buckets cresce exponencialmente com o número de pais da variável. Isso faz com que os valores de \frac{M[x,\pmb{u}]}{M[\pmb{u}]} fiquem muito pequenos, levando à ocorrência da fragmentação e requerendo um maior número de amostras.

Estimativa Bayesiana

Uma maneira alternativa de se fazer estimativas das distribuições dos CPDs é usando uma outra rede Bayesiana. Neste caso supomos que cada amostra de \cal D representa um nó da rede e elas dependem de uma variável contínua que representa a distribuição \theta.

Neste caso a variável \theta tem um CPD correspondente a uma distribuição inicial P(\theta). A ideia é usar a rede para determinar a distribuição a posteriori de P(\theta) uma vez que observamos \cal D. Usando a regra de Bayes temos que,

(1) P(\theta|x[1], \cdots, x[M]) = \dfrac{P(x[1], \cdots, x[M]|\theta) P(\theta)}{P(x[1], \cdots, x[M])}

Veja que o primeiro termo do numerador corresponde à verossimilhança. P(\theta) é a distribuição a priori e P(\theta|x[1], \cdots, x[M]) a distribuição a posteriori. Uma vez que temos o conjunto de amostragem, o denominador é uma constante normalizadora.

Portanto, podemos reescrever (1) como

(2) P(\theta|{\cal D}) \propto P({\cal D}|\theta) P(\theta)

Distribuição multinomial

Suponha que \theta é uma distribuição multinomial sobre k valores. Como um valor a priori da distribuição de \theta, em geral utiliza-se a distribuição de Dirichlet, denotada por Dir(\alpha_1, \cdots, \alpha_k), onde \alpha_i é chamado de hiper-parâmetro. Basicamente, a distribuição é proporcional a um produtório,

Dir(\alpha_1, \cdots, \alpha_k) \propto \prod_{i=1}^{k} {\theta_i}^{\alpha_i - 1}

Intuitivamente os hiper-parâmetros correspondem às frequências de variáveis amostradas discriminadas por valor.

Usando a distribuição de Dirichlet como a distribuição a priori e lembrando que para uma distribuição multinomial a verossimilhança é dada por

P({\cal D}|\theta) = \prod_{i=1}^{k} \theta_i^{M_i}

Por (2), concluímos que a distribuição a posteriori é dada por,

P(\theta|{\cal D}) = \prod_{i=1}^{k} \theta_i^{\alpha_i + M_i - 1}

Note que P(\theta|{\cal D}) = Dir(\alpha_1 + M_1, \cdots, \alpha_k + M_k), ou seja, a distribuição a posteriori continuou sendo uma distribuição de Dirichlet. Devido a essa característica, dizemos que a distribuição de Dirichlet é um conjugado a priori para a distribuição multinomial.

Previsão Bayesiana

Considere uma rede bayesiana simples como a da Figura 1. Suponha que a variável \theta tem distribuição Dir(\alpha_1,\cdots,\alpha_k) = \prod_{j=1}^{k} \theta_j^{\alpha_j - 1}.

Figura 1: Variável condicionada por uma distribuição de Dirichlet

Podemos calcular a distribuição em X multiplicando os CPD’s e marginalizando \theta, ou seja,

(3) P(X) = \int_{\theta} P(X|\theta) P(\theta) d\theta

Para um valor específico de X = x_i, temos que P(X = x_i|\theta) = \theta_i. Substituindo também P(\theta) pela definição da eq. de Dirichlet teremos

P(X = x_i) = \frac{1}{Z} \int_{\theta} \theta_i \left(\prod_{j=1}^{k} \theta_j^{\alpha_j - 1}\right)d\theta

Que segundo as aulas é igual a

(4) = \dfrac{\alpha_i}{\sum_{j=1}^{k} \alpha_j}

Agora, suponha que tenhamos um conjunto de M+1 variáveis condicionadas por \theta, que novamente supomos possuir a distribuição Dir(\alpha_1,\cdots,\alpha_k), conforme a Figura 2.

Figura 2: Rede Bayesiana com M+1 variáveis condicionadas a uma distribuição de Dirichlet

Imagine que observamos as M primeiras variáveis (conjunto \cal D) e queremos prever a distribuição da M+1-ésima, temos que

P(X[M+1] | {\cal D}) = \int_{\theta} P(X[M+1]|{\cal D},\theta) P(\theta|{\cal D}) d\theta

Dado \theta, X[M+1] é independente de \cal D, logo

P(X[M+1] | {\cal D}) = \int_{\theta} P(X[M+1]|\theta) P(\theta|{\cal D}) d\theta

Vimos anteriormente que P(\theta|{\cal D}) = Dir(\alpha_1 + M_1, \cdots, \alpha_k + M_k). Note então que o lado direito da equação acima é idêntico ao lado direito da equação (3). Portanto, podemos calcular $P(X[M+1] = x_i |{\cal D})$, mas usando \alpha_i + M_i, levando a

(5) P(X[M+1] = x_i |{\cal D}) = \dfrac{\alpha_i + M_i}{\sum_{j=1}^{k} (\alpha_j + M_j)}

Definimos como o tamanho da amostra equivalente a soma \alpha = \alpha_1 + \cdots + \alpha_k. Observe que se \alpha_j = 0, teríamos o MLE.

Usando redes Bayesianas, damos um “chute” inicial na distribuição, que eventualmente converge para o mesmo valor que o MLE conforme o número de iterações aumenta. A vantagem é que temos uma maior robustez nas iterações iniciais, quando temos poucas amostras. O gráfico da Figura 3 mostra bem esse fenômeno.

Figura 3: Gráfico de estimativa de uma distribuição de Bernoulli, usando diferentes métodos.

Note que a magnitude de \alpha define o quanto a distribuição inicial terá efeito sobre as distribuições obtidas a priori. Veja como a Dir(.5,.5) (verde-claro) quase não tem diferença sobre a MLE, enquanto Dir(10,10) (preto) muda bastante.

Estimativas Bayesianas para Redes Bayesianas

Podemos generalizar as ideias apresentadas anteriormente para estimar as distribuições dos CPDs de toda uma rede Bayesiana. Consideremos uma rede Bayesiana com uma variável X e uma variável Y que depende condicionalmente de X.

Para estimar os CPDs dessas duas variáveis usando um conjunto de M amostras, podemos construir outra rede Bayesiana conforme a Figura 4.

Figura 4: Estimativa de parâmetros para uma rede Bayesiana”

Note que \theta_X corresponde ao CPD de X e \theta_{Y|X} corresponde ao CPD de Y. Se partirmos de distribuições iniciais para \theta_X e \theta_{Y|X}, temos que as amostras (X[m], Y[m]) são independentes.

Dado o conjunto de amostras {\cal D}, temos que \theta_X e \theta_{Y|X} são independentes. Assim, a distribuição a posteriori de \theta pode ser dada por

P(\theta_X, \theta_{Y|X} | {\cal D}) = P(\theta_X|{\cal D}) P(\theta_{Y|X}|{\cal D})

Podemos refinar o modelo e considerar os possíveis valores de X, criando um CPD de Y para cada valor. Por exemplo, se X for uma variável binária teremos a rede Bayesiana da Figura 5

Figura 5: Estimativa refinada de parâmetros para uma rede Bayesiana

Neste caso, dado o conjunto \cal D, não fica evidente pela rede que \theta_{Y|x^0} é independente de \theta_{Y|x^1}, mas aqui temos independência específica por contexto, pois quando o valor de X é conhecido, uma das duas arestas incidentes a cada Y[m] some.

De qualquer maneira, podemos computar as probabilidades a posteriori de \theta independentemente. Se a distribuição a priori de \theta_{X|\pmb{u}} for Dir(\alpha_{x^1|\pmb{u}}, \cdots, \alpha_{x^k|\pmb{u}}), a distribuição a posteriori será Dir(\alpha_{x^1|\pmb{u}} + M[x^1|\pmb{u}], \cdots, \alpha_{x^k|\pmb{u}} + M[x^k|\pmb{u}])

Uma escolha inicial dos hiper-parâmetros das distribuições de Dirichlet, consiste em escolher um tamanho de amostra equivalente \alpha e para cada CPD distribuí-lo pelo número de entradas. Por exemplo, se o CPD for sobre uma única variável binária, teremos Dir(\alpha/2, \alpha/2). Se tivermos uma variável binária dependendo de outra, fazemos a divisão de \alpha por 4.

Estimativa de parâmetros em Redes de Markov

Vamos apresentar métodos para a estimativa de parâmetros em redes de Markov ou MRF (Markov Random Fields) e para redes de Markov condicionais ou CRF (Conditional Random Fields). Ao invés de um produto de fatores, vamos considerar a representação mais genéricas através de modelos log-lineares.

Verossimilhança máxima para modelos log-lineares

Considere uma distribuição conjunta de um modelo log-linear

P(\pmb{X} : \theta) = \dfrac{1}{Z(\theta)} \tilde P(\pmb{X} : \theta) = \dfrac{1}{Z(\theta)} \exp \left\{\sum_{i=1}^{k} \theta_i f_i(\pmb{x}) \right\}

Note que para simplificar a notação, estamos supondo que f_i ignora todas as entradas de \pmb{x} que não pertencem a seu domínimo D_i. Aqui Z(\theta) é a constante normalizadora dada pela soma de \tilde P(\pmb{X}=\pmb{x}) para todos as possíveis atribuições \pmb{x}, ou seja

Z(\theta) = \sum_{\pmb{x}} \exp (\sum_i \theta_i f_i (\pmb{x}))

Podemos definir a log-verossimilhança de um conjunto de amostras \cal D, denotado por \ell(\theta:{\cal D}) como

\ell(\theta:{\cal D}) = \ln (L(\theta:{\cal D})) = \sum_m(\ln P(\pmb{x}[m] : \theta))

Substituindo P(\pmb{x}[m] : \theta) pela distribuição log-linear temos

= \sum_m (\sum_i \theta_i f_i(\pmb{x}[m]) - \ln Z(\theta))

podemos re-arrajar os somatórios e obter:

= (\sum_i \theta_i (\sum_m f_i(\pmb{x}[m]))) - M\ln Z(\theta)

Diferentemente do caso de redes Bayesianas, aqui não podemos decompor em verossimilhanças locais porque o termo Z(\theta) acopla os termos \theta_1, \cdots, \theta_k.

Porém, podemos mostrar que essa função é côncava, podendo ser resolvida numericamente por um método de descida de gradiente, por exemplo. Na prática, costuma-se utilizar o L-BFGS (Low Memory Broyden Fletcher Goldfarb Shanno).

Para computar os gradientes, vamos considerar a derivada da log-verossimilhança.

\dfrac{\partial}{\partial\theta_i} \dfrac{1}{M} \ell(\theta:{\cal D}) = \dfrac{\partial}{\partial\theta_i} \left((\sum_i \theta_i \dfrac{1}{M} (\sum_m f_i(\pmb{x}[m]))) - \ln Z(\theta)\right)

Note que dividimos a log-verossimilhança pela constante M antes de derivar. Isso é para simplificar o resultado. Como estamos derivando com relacão a \theta_i, o único termo que restará do primeiro somatório é

\dfrac{1}{M} \sum_m f_i(x[m])

que corresponde à média dos valores de f_i considerando todas as amostras. Definimos esse termo como esperança empírica, denotada por E_{\cal D}[f_i(\pmb{X})]. Veja que este valor independe da distribuição \theta.

Já para o segundo termo, aplicamos inicialmente a regra da cadeia para obter

\dfrac{\partial}{\partial\theta_i} \ln Z(\theta)
= \dfrac{1}{Z(\theta)} \sum_{\pmb{x}} \left( \dfrac{\partial}{\partial \theta_i} \exp \{\sum_j \theta_j f_j(\pmb{x}) \} \right)

Podemos aplicar novamente a regra da cadeia para \exp, obtendo

= \dfrac{1}{Z(\theta)} \sum_{\pmb{x}} \left( f_i(\pmb{x}) \exp \{\sum_j \theta_j f_j(\pmb{x}) \} \right)

Rearranjando os termos ficamos com

= \sum_{\pmb{x}} \left( \dfrac{1}{Z(\theta)} \exp \{\sum_j \theta_j f_j(\pmb{x}) \} \right) f_i(\pmb{x})

Veja que o priimeiro termo do somatório é exatamente nossa distribuição original, P(\pmb{X} : \theta) com \pmb{X} = \pmb{x}. Reescrevendo temos:

\dfrac{\partial}{\partial\theta_i} \ln Z(\theta) = \sum_{\pmb{x}} P(\pmb{X} = \pmb{x}) f_i(\pmb{x})

Esse termo representa a esperança de f_i considerando a distribuição \theta e será denotado por E_\theta[f_i]. Veja que para calcular a probabilidade de cada possível atribuição \pmb{x}, podemos usar os métodos de inferência vistos em aulas passadas.

Resumindo, temos o seguinte resultado:

\dfrac{\partial}{\partial\theta_i} \dfrac{1}{M} \ell(\theta:{\cal D}) = E_{\cal D}[f_i(\pmb{X})] - E_\theta[f_i]

Usando essa diferença para construir o gradiente, obtemos finalmente o seguinte teorema:

Teorema. A distribuição \hat \theta resulta na MLE se e somente se E_{\cal D}[f_i(X)] = E_\theta[f_i] para todo i=1,\cdots,k

MLE para CRFs

No caso de CRFs (Conditional Random Fields), temos a distribuição conjunta condicional

P_\theta(\pmb{Y}|\pmb{x}) = \dfrac{1}{Z_{\pmb{x}}(\theta)}\tilde P_\theta(\pmb{x}, \pmb{Y}),

onde Z_{\pmb{x}}(\theta) = \sum_{\pmb{Y}} \tilde P_\theta(\pmb{x}, \pmb{Y})

Considere o conjunto de amostras {\cal D} = \{(x[m], y[m])\}_{m=1}^M. Podemos novamente definir a log-verossimilhança:

\ell_{Y|x}(\theta:{\cal D}) = \sum_{m=1}^M \ln P_\theta(y[m], x[m] : \theta)

= \sum_m (\sum_i \theta_i f_i(\pmb{x}[m]) - \ln Z_{x[m]}(\theta))

Note que neste caso não podemos extrair o termo Z_{x[m]}(\theta) para fora do somatório porque ele depende de x[m].

Para determinar o gradiente, podemos novamente calcular a derivada em relação a \theta_i:

\dfrac{\partial}{\partial \theta_i} \dfrac{1}{M} \ell_{Y|x}(\theta:{\cal D}) = \dfrac{\partial}{\partial \theta_i} \dfrac{1}{M} \sum_m (\sum_i \theta_i f_i(\pmb{x}[m]) - \ln Z_{x[m]}(\theta))

Usando as mesmas ideias do caso de MRF’s, concluiremos que

\dfrac{\partial}{\partial \theta_i} \dfrac{1}{M} \ell_{Y|X}(\theta:{\cal D}) = \dfrac{1}{M} \sum_{m=1}^{M} (f_i(x[m],y[m]) - E_\theta[f_i(x[m], Y)])

Onde E_\theta[f_i(x[m], Y)] = \sum_Y P_\theta(Y|x[m]) f_j(Y, x[m])

A diferença aqui é que a cada passo do algoritmo, devemos recalcular P_\theta(Y|x[m])] para cada amostra x[m]. Porém, em geral o conjunto Y será menor para CRF’s.

Estimativa MAP para MRFs e CRFs

Da mesma forma que para redes Bayesianas, podemos “chutar” uma distribuição a priori e obter uma distribuição a posteriori. Entretanto, devido ao acoplamento presente na verossimilhança de redes de Markov, não há uma fórmula fechada para a distribuição, como a de Dirichlet para redes Bayesianas.

O que fazemos é manter um conjunto de parâmetros que maximizem a dstribuição a posteriori, ou seja, uma MAP (Maximum A Posteriori).

Uma possível distribuição a pirori que pode ser usada é a distribuição gaussiana com métia 0.

P(\theta : \sigma^2) = \prod_{i=1}^{k} \dfrac{1}{\sqrt{2\pi}\sigma} \exp \left\{ - \dfrac{\theta_i^2}{2\sigma^2} \right\}

Intuitivamente, quanto maior o valor de \sigma, maior a confiança de que o valor de \theta_i esteja próximo de 0.

Outra alternativa é usar a distribuição de Laplace,

P(\theta:\beta) = \prod_{i=1}^{k} \dfrac{1}{2\beta} \exp \left\{-\dfrac{|\theta_i|}{\beta} \right\}

Queremos encontrar a distribuição \theta que maximiza a distribuição conjunta dado o conjunto de amostras, ou seja,

\mbox{argmax}_\theta P({\cal D}, \theta) = \mbox{argmax}_\theta P({\cal D} | \theta) P(\theta)

Neste caso P({\cal D} | \theta) é a verossimilhança e P(\theta) é a distribuição a priori. Como o logaritmo é monotomicamente crescente com sua entrada, podemos aplicar o log ao termo dentro de argmax:

= \mbox{argmax}_\theta(\ell(\theta:{\cal D}) + \log P(\theta))

O plot de \log P(\theta) para a distribuição Gaussiana, é uma função quadrática. Esta curva corresponde a uma penalização da log-verossimilhança e é conhecida como regularização-L2.

De maneira análoga, o plot para a distribuição de Laplace é uma função linear e é equivalente à regularização-L1.

Comentários Desativados em Probabilistic Graphical Models – Semana 7 | inteligência artificial, stanford | Permalink
Publicado por kunigami