Arquivo de Python - Página 3 de 31

Aprendizado de Parametros (MAP, ML)

16/07/202618/06/2026 Por antonino

O que é aprendizado de parâmetros?

Aprendizado de parâmetros é o processo de estimar os valores numéricos das probabilidades condicionais em uma rede bayesiana. Dada a estrutura do grafo, ele usa os dados para preencher as tabelas CPT (probabilidades condicionais). Existem duas abordagens principais: máxima verossimilhança (ML) e máxima a posteriori (MAP). ML estima os parâmetros que maximizam a verossimilhança dos dados observados. MAP incorpora um prior sobre os parâmetros, regularizando as estimativas. ML é equivalente a contar frequências e normalizar (para variáveis discretas). MAP usa suavização (ex.: Laplace) para evitar probabilidades zero em combinações não vistas. O aprendizado de parâmetros é feito separadamente para cada nó, dado seus pais. Ele é computacionalmente simples e pode ser feito de forma incremental.

Características da estimação por máxima verossimilhança (ML)

O método ML possui três características principais que o definem. Primeiro, ele é consistente: converge para o valor verdadeiro com muitos dados. Segundo, ele é não-viesado assintoticamente para modelos corretamente especificados. Terceiro, ele pode superajustar dados pequenos, atribuindo probabilidade zero a eventos não observados. A estimativa ML para uma variável discreta com pais é a frequência relativa: P(X=x | pais=p) = contagem(x,p) / contagem(p). Ele é o estimador de escolha quando se tem muitos dados e nenhum conhecimento prévio.

Características da estimação por máxima a posteriori (MAP)

O método MAP também possui três características principais. Primeiro, ele incorpora um prior que suaviza as estimativas (ex.: Laplace, Dirichlet). Segundo, ele reduz o overfitting em dados escassos. Terceiro, ele é um compromisso entre os dados e a crença prévia. MAP equivale a adicionar pseudocontagens ao numerador e denominador. Para prior de Laplace (aditivo-1), P(X=x | pais=p) = (contagem(x,p) + 1) / (contagem(p) + k), onde k é o número de valores de X.

Vantagens e aplicações típicas

A principal vantagem do ML é a simplicidade e a objetividade. MAP é preferível quando os dados são escassos ou há conhecimento especialista. Ambos são usados em diagnóstico, classificação e sistemas de recomendação.

O aprendizado de parâmetros é a segunda etapa após a definição da estrutura. Em redes com muitas variáveis, as CPTs podem ser enormes (exponencial nos pais). Para reduzir a complexidade, usa-se modelos de regressão logística ou árvores de decisão. ML pode ser calculado por fechamento analítico para distribuições da família exponencial. MAP com prior Dirichlet é conjugado para distribuições multinomiais (discretas). A força do prior (tamanho da amostra virtual) controla a regularização. Prior fraco (ex.: pseudocontagens pequenas) aproxima o ML. Prior forte força as estimativas a ficarem próximas do prior. Em problemas contínuos, ML e MAP estimam médias e variâncias de distribuições normais. MAP é frequentemente usado em aprendizado de máquina com regularização L2 (equivalente a prior normal). A escolha entre ML e MAP depende do tamanho da amostra e da confiança no prior. Em dados grandes, ML e MAP convergem para o mesmo valor. O aprendizado de parâmetros pode ser feito online (atualização recursiva) com priores conjugados. Assim, o aprendizado de parâmetros é uma etapa fundamental em modelagem probabilística.

Um exemplo clássico é estimar a probabilidade de uma moeda dar cara. ML estima P(cara) = número de caras / total de lançamentos. MAP com prior Beta(α, β) estima (caras + α – 1) / (total + α + β – 2). Com poucos lançamentos, MAP puxa a estimativa para a média do prior.

Enunciado do exemplo clássico

Implemente o aprendizado de parâmetros ML e MAP para uma rede com uma variável binária X e um pai Y (também binário). Dados: 100 observações com contagens: Y=0 (40 vezes), Y=1 (60 vezes). Dentro de Y=0: X=0 (30), X=1 (10). Dentro de Y=1: X=0 (20), X=1 (40). Calcule as CPTs por ML e MAP (usando prior Beta(2,2) para cada condicional). Plote as estimativas e compare com as probabilidades verdadeiras (se conhecidas).

import numpy as np  
import matplotlib.pyplot as plt  

# Dados sintéticos (contagens)  
# Y=0: 40 amostras, X=0:30, X=1:10  
# Y=1: 60 amostras, X=0:20, X=1:40  
contagens = {  
    (0,0): 30,   # Y=0, X=0  
    (0,1): 10,   # Y=0, X=1  
    (1,0): 20,   # Y=1, X=0  
    (1,1): 40    # Y=1, X=1  
}  
total_Y0 = 40  
total_Y1 = 60  

# Estimativa ML  
P_X_dado_Y_ML = {  
    (0,0): contagens[(0,0)] / total_Y0,  
    (0,1): contagens[(0,1)] / total_Y0,  
    (1,0): contagens[(1,0)] / total_Y1,  
    (1,1): contagens[(1,1)] / total_Y1  
}  

# Estimativa MAP com prior Beta(2,2) (α=2, β=2)  
# Fórmula: (contagem + α - 1) / (total + α + β - 2)  
# Para Beta(2,2), α-1=1, α+β-2=2  
alpha_beta = 2  
pseudo = alpha_beta - 1  # = 1  
denom_ajuste = 2 * pseudo  # = 2 (pois α+β-2 = 2)  

P_X_dado_Y_MAP = {  
    (0,0): (contagens[(0,0)] + pseudo) / (total_Y0 + denom_ajuste),  
    (0,1): (contagens[(0,1)] + pseudo) / (total_Y0 + denom_ajuste),  
    (1,0): (contagens[(1,0)] + pseudo) / (total_Y1 + denom_ajuste),  
    (1,1): (contagens[(1,1)] + pseudo) / (total_Y1 + denom_ajuste)  
}  

# Probabilidades verdadeiras (geradas de forma arbitrária para comparação)  
verdade = {  
    (0,0): 0.75,  # P(X=0 | Y=0)  
    (0,1): 0.25,  # P(X=1 | Y=0)  
    (1,0): 0.35,  # P(X=0 | Y=1)  
    (1,1): 0.65   # P(X=1 | Y=1)  
}  

print("Estimativas ML:")  
for (y,x), p in P_X_dado_Y_ML.items():  
    print(f"  P(X={x} | Y={y}) = {p:.3f}")  

print("\nEstimativas MAP (Beta(2,2)):")  
for (y,x), p in P_X_dado_Y_MAP.items():  
    print(f"  P(X={x} | Y={y}) = {p:.3f}")  

print("\nValores verdadeiros (para referência):")  
for (y,x), p in verdade.items():  
    print(f"  P(X={x} | Y={y}) = {p:.3f}")  

# Gráficos  
y_vals = [0, 1]  
x_vals = [0, 1]  
bar_width = 0.2  
pos = np.arange(len(y_vals))  

fig, axes = plt.subplots(1, 2, figsize=(12, 5))  

# Para Y=0  
ax = axes[0]  
ml_0 = [P_X_dado_Y_ML[(0,0)], P_X_dado_Y_ML[(0,1)]]  
map_0 = [P_X_dado_Y_MAP[(0,0)], P_X_dado_Y_MAP[(0,1)]]  
verdade_0 = [verdade[(0,0)], verdade[(0,1)]]  
ax.bar(pos - bar_width, ml_0, width=bar_width, label='ML', color='blue')  
ax.bar(pos, map_0, width=bar_width, label='MAP', color='green')  
ax.bar(pos + bar_width, verdade_0, width=bar_width, label='Verdade', color='red', alpha=0.6)  
ax.set_xticks(pos)  
ax.set_xticklabels(['X=0', 'X=1'])  
ax.set_ylabel('Probabilidade')  
ax.set_title('P(X | Y=0)')  
ax.legend()  
ax.grid(True, axis='y', linestyle='--')  

# Para Y=1  
ax = axes[1]  
ml_1 = [P_X_dado_Y_ML[(1,0)], P_X_dado_Y_ML[(1,1)]]  
map_1 = [P_X_dado_Y_MAP[(1,0)], P_X_dado_Y_MAP[(1,1)]]  
verdade_1 = [verdade[(1,0)], verdade[(1,1)]]  
ax.bar(pos - bar_width, ml_1, width=bar_width, label='ML', color='blue')  
ax.bar(pos, map_1, width=bar_width, label='MAP', color='green')  
ax.bar(pos + bar_width, verdade_1, width=bar_width, label='Verdade', color='red', alpha=0.6)  
ax.set_xticks(pos)  
ax.set_xticklabels(['X=0', 'X=1'])  
ax.set_ylabel('Probabilidade')  
ax.set_title('P(X | Y=1)')  
ax.legend()  
ax.grid(True, axis='y', linestyle='--')  

plt.tight_layout()  
plt.show()  

# Análise de sensibilidade: variando a força do prior  
forcas = [1, 2, 5, 10]  # α = β = forca  
fig2, ax2 = plt.subplots(figsize=(8, 5))  
for f in forcas:  
    pseudo_f = f - 1  
    denom_f = 2 * pseudo_f  
    p_map = (contagens[(0,1)] + pseudo_f) / (total_Y0 + denom_f)  
    ax2.scatter(f, p_map, label=f'α=β={f}', s=80)  
ax2.axhline(P_X_dado_Y_ML[(0,1)], color='blue', linestyle='--', label='ML')  
ax2.axhline(verdade[(0,1)], color='red', linestyle=':', label='Verdadeiro')  
ax2.set_xlabel('Força do prior (α=β)')  
ax2.set_ylabel('P(X=1 | Y=0) estimado')  
ax2.set_title('Sensibilidade MAP à Força do Prior')  
ax2.legend()  
ax2.grid(True)  
plt.tight_layout()  
plt.show()

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

import numpy as np

import matplotlib.pyplot as plt

# Dados sintéticos (contagens)

# Y=0: 40 amostras, X=0:30, X=1:10

# Y=1: 60 amostras, X=0:20, X=1:40

contagens = {

(0,0): 30, # Y=0, X=0

(0,1): 10, # Y=0, X=1

(1,0): 20, # Y=1, X=0

(1,1): 40 # Y=1, X=1

}

total_Y0 = 40

total_Y1 = 60

# Estimativa ML

P_X_dado_Y_ML = {

(0,0): contagens[(0,0)] / total_Y0,

(0,1): contagens[(0,1)] / total_Y0,

(1,0): contagens[(1,0)] / total_Y1,

(1,1): contagens[(1,1)] / total_Y1

}

# Estimativa MAP com prior Beta(2,2) (α=2, β=2)

# Fórmula: (contagem + α - 1) / (total + α + β - 2)

# Para Beta(2,2), α-1=1, α+β-2=2

alpha_beta = 2

pseudo = alpha_beta - 1 # = 1

denom_ajuste = 2 * pseudo # = 2 (pois α+β-2 = 2)

P_X_dado_Y_MAP = {

(0,0): (contagens[(0,0)] + pseudo) / (total_Y0 + denom_ajuste),

(0,1): (contagens[(0,1)] + pseudo) / (total_Y0 + denom_ajuste),

(1,0): (contagens[(1,0)] + pseudo) / (total_Y1 + denom_ajuste),

(1,1): (contagens[(1,1)] + pseudo) / (total_Y1 + denom_ajuste)

}

# Probabilidades verdadeiras (geradas de forma arbitrária para comparação)

verdade = {

(0,0): 0.75, # P(X=0 | Y=0)

(0,1): 0.25, # P(X=1 | Y=0)

(1,0): 0.35, # P(X=0 | Y=1)

(1,1): 0.65 # P(X=1 | Y=1)

}

print("Estimativas ML:")

for (y,x), p in P_X_dado_Y_ML.items():

print(f" P(X={x} | Y={y}) = {p:.3f}")

print("\nEstimativas MAP (Beta(2,2)):")

for (y,x), p in P_X_dado_Y_MAP.items():

print(f" P(X={x} | Y={y}) = {p:.3f}")

print("\nValores verdadeiros (para referência):")

for (y,x), p in verdade.items():

print(f" P(X={x} | Y={y}) = {p:.3f}")

# Gráficos

y_vals = [0, 1]

x_vals = [0, 1]

bar_width = 0.2

pos = np.arange(len(y_vals))

fig, axes = plt.subplots(1, 2, figsize=(12, 5))

# Para Y=0

ax = axes[0]

ml_0 = [P_X_dado_Y_ML[(0,0)], P_X_dado_Y_ML[(0,1)]]

map_0 = [P_X_dado_Y_MAP[(0,0)], P_X_dado_Y_MAP[(0,1)]]

verdade_0 = [verdade[(0,0)], verdade[(0,1)]]

ax.bar(pos - bar_width, ml_0, width=bar_width, label='ML', color='blue')

ax.bar(pos, map_0, width=bar_width, label='MAP', color='green')

ax.bar(pos + bar_width, verdade_0, width=bar_width, label='Verdade', color='red', alpha=0.6)

ax.set_xticks(pos)

ax.set_xticklabels(['X=0', 'X=1'])

ax.set_ylabel('Probabilidade')

ax.set_title('P(X | Y=0)')

ax.legend()

ax.grid(True, axis='y', linestyle='--')

# Para Y=1

ax = axes[1]

ml_1 = [P_X_dado_Y_ML[(1,0)], P_X_dado_Y_ML[(1,1)]]

map_1 = [P_X_dado_Y_MAP[(1,0)], P_X_dado_Y_MAP[(1,1)]]

verdade_1 = [verdade[(1,0)], verdade[(1,1)]]

ax.bar(pos - bar_width, ml_1, width=bar_width, label='ML', color='blue')

ax.bar(pos, map_1, width=bar_width, label='MAP', color='green')

ax.bar(pos + bar_width, verdade_1, width=bar_width, label='Verdade', color='red', alpha=0.6)

ax.set_xticks(pos)

ax.set_xticklabels(['X=0', 'X=1'])

ax.set_ylabel('Probabilidade')

ax.set_title('P(X | Y=1)')

ax.legend()

ax.grid(True, axis='y', linestyle='--')

plt.tight_layout()

plt.show()

# Análise de sensibilidade: variando a força do prior

forcas = [1, 2, 5, 10] # α = β = forca

fig2, ax2 = plt.subplots(figsize=(8, 5))

for f in forcas:

pseudo_f = f - 1

denom_f = 2 * pseudo_f

p_map = (contagens[(0,1)] + pseudo_f) / (total_Y0 + denom_f)

ax2.scatter(f, p_map, label=f'α=β={f}', s=80)

ax2.axhline(P_X_dado_Y_ML[(0,1)], color='blue', linestyle='--', label='ML')

ax2.axhline(verdade[(0,1)], color='red', linestyle=':', label='Verdadeiro')

ax2.set_xlabel('Força do prior (α=β)')

ax2.set_ylabel('P(X=1 | Y=0) estimado')

ax2.set_title('Sensibilidade MAP à Força do Prior')

ax2.legend()

ax2.grid(True)

plt.tight_layout()

plt.show()

Este código estima parâmetros por ML e MAP para uma CPT simples. ML segue as frequências exatas dos dados, podendo ser instável com amostras pequenas. MAP suaviza as estimativas, puxando-as em direção à probabilidade uniforme (0.5). Quanto maior a força do prior, mais as estimativas se afastam dos dados. Para iniciantes, este exemplo mostra a diferença prática entre ML e MAP. O aprendizado de parâmetros é, portanto, uma escolha entre objetividade e regularização.

Aprendizado de Estrutura (Escore vs Restricao)

16/07/202618/06/2026 Por antonino

O que é aprendizado de estrutura em redes bayesianas?

Aprendizado de estrutura é o processo de descobrir o grafo (arestas) a partir dos dados. Ele determina quais variáveis influenciam quais, sem conhecimento prévio do especialista. Existem duas abordagens principais: baseada em escores e baseada em restrições. A abordagem por escores busca o grafo que maximiza uma medida de qualidade (ex.: BIC). A abordagem por restrições usa testes de independência condicional para decidir arestas. Ambas têm vantagens e desvantagens em termos de precisão e custo computacional. O aprendizado de estrutura é um problema NP-difícil para grafos gerais. Portanto, heurísticas como busca gulosa ou algoritmos genéticos são usadas.

Características da abordagem por escores

A abordagem por escores possui três características principais que a definem. Primeiro, ela atribui um número (escore) a cada grafo candidato. Segundo, o escore é baseado na verossimilhança penalizada (AIC, BIC) ou fator de Bayes. Terceiro, a busca percorre o espaço de grafos usando operadores de adição/remoção de arestas. O escore BIC = -2*log(L) + k*log(n), onde k é o número de parâmetros. Ela é consistente: com dados suficientes, encontra o grafo verdadeiro (para BIC).

Características da abordagem por restrições

A abordagem por restrições também possui três características principais. Primeiro, ela usa testes estatísticos de independência condicional (ex.: qui-quadrado). Segundo, ela constrói o grafo a partir das relações de dependência detectadas. Terceiro, algoritmos como PC (Peter-Clark) e FCI são exemplos clássicos. Ela é computacionalmente mais rápida que a busca por escores para redes esparsas. Contudo, ela depende da escolha do nível de significância e pode acumular erros.

Comparação e aplicações típicas

A principal diferença é que escores buscam o grafo globalmente, restrições localmente. Escores são mais precisos, mas mais lentos para muitas variáveis. Restrições são escaláveis, mas sensíveis a testes múltiplos. A escolha depende do tamanho dos dados e do conhecimento prévio.

O aprendizado de estrutura é usado em bioinformática (redes de regulação gênica). Também em epidemiologia (causalidade entre fatores de risco) e finanças. A abordagem por escores é mais popular quando o número de variáveis é moderado (ex.: <30). A abordagem por restrições é preferida para redes com centenas de nós. Na prática, combina-se as duas: usa-se restrições para reduzir o espaço de busca. Depois, aplica-se busca por escores no espaço reduzido (abordagem híbrida). O escore BIC é derivado da verossimilhança marginal e penaliza complexidade. O fator de Bayes compara dois grafos diretamente, mas é mais custoso. A abordagem por restrições pode ser sensível à ordem das variáveis (caso PC). Para lidar com isso, usa-se amostragem bootstrap para estabilidade. Aprendizado de estrutura também pode incorporar conhecimento de especialistas (arestas fixas). Isso melhora a precisão quando os dados são escassos. A validação cruzada pode ser usada para escolher entre estruturas concorrentes. Assim, o aprendizado de estrutura é uma área ativa de pesquisa em IA.

Um exemplo clássico é a rede de diagnóstico com variáveis: Idade, Fuma, Tosse, Câncer, e Fadiga. Os dados observacionais permitem descobrir que Fuma → Câncer, e Câncer → Tosse, Fadiga. O aprendizado por escore BIC identifica essa estrutura com alta probabilidade.

Enunciado do exemplo clássico

Implemente o aprendizado de estrutura por escore BIC (busca gulosa) e por restrições (teste qui-quadrado) para um conjunto de dados sintético. Gere 100 amostras de uma rede com 4 variáveis: A → B, A → C, B → D, C → D. Use busca gulosa por escore BIC (operadores de adicionar/remover/reverter arestas). Use o algoritmo PC (com teste qui-quadrado) para a abordagem de restrições. Compare as estruturas encontradas com a estrutura verdadeira. Plote as três estruturas (verdadeira, BIC, PC) lado a lado.

import numpy as np  
import matplotlib.pyplot as plt  
import networkx as nx  
from itertools import combinations  
from scipy.stats import chi2_contingency  

# Gerar dados sintéticos da estrutura verdadeira  
# A, B, C, D binárias (0/1)  
np.random.seed(99)  
n_amostras = 100  
# Definição das probabilidades  
P_A = {0: 0.6, 1: 0.4}  
P_B_dado_A = {(0,0):0.7, (0,1):0.3, (1,0):0.2, (1,1):0.8}  # A, B  
P_C_dado_A = {(0,0):0.6, (0,1):0.4, (1,0):0.1, (1,1):0.9}  
P_D_dado_B_C = {(0,0,0):0.9, (0,0,1):0.1, (0,1,0):0.5, (0,1,1):0.5,  
                 (1,0,0):0.3, (1,0,1):0.7, (1,1,0):0.2, (1,1,1):0.8}  

# Gerar amostras  
dados = np.zeros((n_amostras, 4), dtype=int)  
for i in range(n_amostras):  
    A = 1 if np.random.random() < P_A[1] else 0  
    B = 1 if np.random.random() < P_B_dado_A[(A,1)] else 0  
    C = 1 if np.random.random() < P_C_dado_A[(A,1)] else 0  
    D = 1 if np.random.random() < P_D_dado_B_C[(B,C,1)] else 0  
    dados[i] = [A, B, C, D]  
var_names = ['A', 'B', 'C', 'D']  

# 1. Aprendizado por escore BIC (busca gulosa)  
def escore_bic(dados, pais, var_idx):  
    # pais: lista de índices de pais para a variável var_idx  
    # Estimar probabilidades condicionais e calcular BIC  
    n = len(dados)  
    # Contar combinações  
    if not pais:  
        # Sem pais: estimar P(var=1)  
        contagens = np.bincount(dados[:, var_idx], minlength=2)  
        # BIC para modelo sem pais: 1 parâmetro  
        logL = contagens[1]*np.log(contagens[1]/n) if contagens[1]>0 else 0  
        logL += contagens[0]*np.log(contagens[0]/n) if contagens[0]>0 else 0  
        k = 1  
    else:  
        # Com pais: tabela de contingência  
        # Agrupar por combinação de pais  
        combos = dados[:, pais]  
        # Converter para índice único  
        # Para simplificar, usar tuplas  
        from collections import defaultdict  
        contagem = defaultdict(lambda: np.zeros(2, dtype=int))  
        for i in range(n):  
            chave = tuple(dados[i, p] for p in pais)  
            contagem[chave][dados[i, var_idx]] += 1  
        logL = 0.0  
        for chave, cnts in contagem.items():  
            total = sum(cnts)  
            for val in (0,1):  
                if cnts[val] > 0:  
                    logL += cnts[val] * np.log(cnts[val] / total)  
        k = len(contagem) * 1  # cada combinação tem 1 parâmetro (prob de 1)  
    # BIC = -2*logL + k*log(n)  
    return -2*logL + k*np.log(n)  

# Busca gulosa: começa com grafo vazio, adiciona arestas que melhoram BIC  
def busca_gulosa(dados, var_names):  
    m = dados.shape[1]  
    # Grafo representado como lista de pais para cada nó  
    pais = [set() for _ in range(m)]  
    melhorou = True  
    while melhorou:  
        melhorou = False  
        melhor_esc = sum(escore_bic(dados, list(pais[i]), i) for i in range(m))  
        # Tentar adicionar aresta  
        for i in range(m):  
            for j in range(m):  
                if i == j: continue  
                if j in pais[i]: continue  
                # Adicionar j como pai de i  
                pais[i].add(j)  
                novo_esc = sum(escore_bic(dados, list(pais[k]), k) for k in range(m))  
                if novo_esc < melhor_esc:  # BIC menor = melhor  
                    melhor_esc = novo_esc  
                    melhorou = True  
                    # Manter a adição  
                else:  
                    pais[i].remove(j)  
            if melhorou:  
                break  
        if melhorou:  
            continue  
        # Tentar remover aresta  
        for i in range(m):  
            for j in list(pais[i]):  
                pais[i].remove(j)  
                novo_esc = sum(escore_bic(dados, list(pais[k]), k) for k in range(m))  
                if novo_esc < melhor_esc:  
                    melhor_esc = novo_esc  
                    melhorou = True  
                    # Manter remoção  
                else:  
                    pais[i].add(j)  
            if melhorou:  
                break  
    # Converter para matriz de adjacência (orientada)  
    grafo_bic = np.zeros((m, m), dtype=int)  
    for i in range(m):  
        for j in pais[i]:  
            grafo_bic[j, i] = 1  # j -> i  
    return grafo_bic  

grafo_bic = busca_gulosa(dados, var_names)  

# 2. Aprendizado por restrições (algoritmo PC simplificado)  
# Teste qui-quadrado para independência condicional entre X e Y dado Z  
def teste_chi2(dados, x, y, z_list):  
    # Construir tabela de contingência condicionada em Z  
    if not z_list:  
        # Teste marginal  
        tabela = np.zeros((2,2))  
        for i in range(len(dados)):  
            tabela[dados[i,x], dados[i,y]] += 1  
        chi2, p, dof, _ = chi2_contingency(tabela, correction=False)  
        return p  
    else:  
        # Para cada combinação de Z, fazer teste e combinar p-valores (Fisher)  
        from collections import defaultdict  
        from scipy.stats import combine_pvalues  
        pvals = []  
        for z_vals in set(tuple(dados[i, z_list]) for i in range(len(dados))):  
            mask = np.all(dados[:, z_list] == z_vals, axis=1)  
            if np.sum(mask) < 5:  
                continue  
            sub = dados[mask]  
            tabela = np.zeros((2,2))  
            for i in range(len(sub)):  
                tabela[sub[i,x], sub[i,y]] += 1  
            if np.sum(tabela) < 5:  
                continue  
            chi2, p, dof, _ = chi2_contingency(tabela, correction=False)  
            pvals.append(p)  
        if not pvals:  
            return 1.0  
        # Combinar p-valores (método de Fisher)  
        stat, p_combined = combine_pvalues(pvals, method='fisher')  
        return p_combined  

def algoritmo_pc(dados, alpha=0.05):  
    m = dados.shape[1]  
    # Grafo não-direcionado inicial (completo)  
    adj = np.ones((m, m), dtype=int) - np.eye(m, dtype=int)  
    # Passo 1: remover arestas com testes de independência marginal e condicional  
    for ordem in range(0, m-1):  
        # Para cada par (x,y) ainda conectado  
        for x in range(m):  
            for y in range(m):  
                if x >= y: continue  
                if adj[x,y] == 0: continue  
                # Encontrar vizinhos comuns (excluindo x e y)  
                vizinhos = [v for v in range(m) if v != x and v != y and adj[x,v] and adj[y,v]]  
                if len(vizinhos) < ordem:  
                    continue  
                # Testar todas as combinações de ordem variáveis  
                for subset in combinations(vizinhos, ordem):  
                    p_val = teste_chi2(dados, x, y, list(subset))  
                    if p_val > alpha:  # independentes -> remover aresta  
                        adj[x,y] = adj[y,x] = 0  
                        break  
    # Orientar arestas (regras de colisor simples)  
    grafo_dir = np.zeros((m, m), dtype=int)  
    # Para cada tripla x - z - y, se x e y não são adjacentes, z é colisor  
    for x in range(m):  
        for z in range(m):  
            if x == z: continue  
            for y in range(m):  
                if y == x or y == z: continue  
                if adj[x,z] and adj[z,y] and not adj[x,y]:  
                    # x -> z <- y (colisor)  
                    grafo_dir[x,z] = 1  
                    grafo_dir[y,z] = 1  
    # Para as arestas restantes, orientar sem criar ciclos (simplificado: usar ordem natural)  
    for x in range(m):  
        for y in range(m):  
            if x < y and adj[x,y] and grafo_dir[x,y]==0 and grafo_dir[y,x]==0:  
                # Orientar arbitrariamente (baseado em ordem)  
                grafo_dir[x,y] = 1  
    return grafo_dir  

grafo_pc = algoritmo_pc(dados, alpha=0.05)  

# Estrutura verdadeira  
verdade = np.zeros((4,4), dtype=int)  
# A->B (0->1), A->C (0->2), B->D (1->3), C->D (2->3)  
verdade[0,1] = 1  
verdade[0,2] = 1  
verdade[1,3] = 1  
verdade[2,3] = 1  

print("Estrutura verdadeira (arestas):")  
print(verdade)  
print("\nEstrutura por BIC (busca gulosa):")  
print(grafo_bic)  
print("\nEstrutura por PC (restrições):")  
print(grafo_pc)  

# Plotar as três estruturas  
def plot_graph(adj, title, ax):  
    G = nx.DiGraph()  
    m = adj.shape[0]  
    for i in range(m):  
        for j in range(m):  
            if adj[i,j]:  
                G.add_edge(var_names[i], var_names[j])  
    pos = nx.circular_layout(G)  
    nx.draw(G, pos, ax=ax, with_labels=True, node_color='lightblue',  
            node_size=800, font_weight='bold', arrowsize=20,  
            edge_color='gray', width=2)  
    ax.set_title(title)  

fig, axes = plt.subplots(1, 3, figsize=(15, 4))  
plot_graph(verdade, 'Estrutura Verdadeira', axes[0])  
plot_graph(grafo_bic, 'Aprendizado por Escore (BIC)', axes[1])  
plot_graph(grafo_pc, 'Aprendizado por Restrições (PC)', axes[2])  
plt.tight_layout()  
plt.show()

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130

131

132

133

134

135

136

137

138

139

140

141

142

143

144

145

146

147

148

149

150

151

152

153

154

155

156

157

158

159

160

161

162

163

164

165

166

167

168

169

170

171

172

173

174

175

176

177

178

179

180

181

182

183

184

185

186

187

188

189

190

191

192

193

194

195

196

197

198

199

200

201

202

203

204

205

206

207

208

209

210

211

212

213

214

215

216

217

218

219

220

221

222

import numpy as np

import matplotlib.pyplot as plt

import networkx as nx

from itertools import combinations

from scipy.stats import chi2_contingency

# Gerar dados sintéticos da estrutura verdadeira

# A, B, C, D binárias (0/1)

np.random.seed(99)

n_amostras = 100

# Definição das probabilidades

P_A = {0: 0.6, 1: 0.4}

P_B_dado_A = {(0,0):0.7, (0,1):0.3, (1,0):0.2, (1,1):0.8} # A, B

P_C_dado_A = {(0,0):0.6, (0,1):0.4, (1,0):0.1, (1,1):0.9}

P_D_dado_B_C = {(0,0,0):0.9, (0,0,1):0.1, (0,1,0):0.5, (0,1,1):0.5,

(1,0,0):0.3, (1,0,1):0.7, (1,1,0):0.2, (1,1,1):0.8}

# Gerar amostras

dados = np.zeros((n_amostras, 4), dtype=int)

for i in range(n_amostras):

A = 1 if np.random.random() < P_A[1] else 0

B = 1 if np.random.random() < P_B_dado_A[(A,1)] else 0

C = 1 if np.random.random() < P_C_dado_A[(A,1)] else 0

D = 1 if np.random.random() < P_D_dado_B_C[(B,C,1)] else 0

dados[i] = [A, B, C, D]

var_names = ['A', 'B', 'C', 'D']

# 1. Aprendizado por escore BIC (busca gulosa)

def escore_bic(dados, pais, var_idx):

# pais: lista de índices de pais para a variável var_idx

# Estimar probabilidades condicionais e calcular BIC

n = len(dados)

# Contar combinações

if not pais:

# Sem pais: estimar P(var=1)

contagens = np.bincount(dados[:, var_idx], minlength=2)

# BIC para modelo sem pais: 1 parâmetro

logL = contagens[1]*np.log(contagens[1]/n) if contagens[1]>0 else 0

logL += contagens[0]*np.log(contagens[0]/n) if contagens[0]>0 else 0

k = 1

else:

# Com pais: tabela de contingência

# Agrupar por combinação de pais

combos = dados[:, pais]

# Converter para índice único

# Para simplificar, usar tuplas

from collections import defaultdict

contagem = defaultdict(lambda: np.zeros(2, dtype=int))

for i in range(n):

chave = tuple(dados[i, p] for p in pais)

contagem[chave][dados[i, var_idx]] += 1

logL = 0.0

for chave, cnts in contagem.items():

total = sum(cnts)

for val in (0,1):

if cnts[val] > 0:

logL += cnts[val] * np.log(cnts[val] / total)

k = len(contagem) * 1 # cada combinação tem 1 parâmetro (prob de 1)

# BIC = -2*logL + k*log(n)

return -2*logL + k*np.log(n)

# Busca gulosa: começa com grafo vazio, adiciona arestas que melhoram BIC

def busca_gulosa(dados, var_names):

m = dados.shape[1]

# Grafo representado como lista de pais para cada nó

pais = [set() for _ in range(m)]

melhorou = True

while melhorou:

melhorou = False

melhor_esc = sum(escore_bic(dados, list(pais[i]), i) for i in range(m))

# Tentar adicionar aresta

for i in range(m):

for j in range(m):

if i == j: continue

if j in pais[i]: continue

# Adicionar j como pai de i

pais[i].add(j)

novo_esc = sum(escore_bic(dados, list(pais[k]), k) for k in range(m))

if novo_esc < melhor_esc: # BIC menor = melhor

melhor_esc = novo_esc

melhorou = True

# Manter a adição

else:

pais[i].remove(j)

if melhorou:

break

if melhorou:

continue

# Tentar remover aresta

for i in range(m):

for j in list(pais[i]):

pais[i].remove(j)

novo_esc = sum(escore_bic(dados, list(pais[k]), k) for k in range(m))

if novo_esc < melhor_esc:

melhor_esc = novo_esc

melhorou = True

# Manter remoção

else:

pais[i].add(j)

if melhorou:

break

# Converter para matriz de adjacência (orientada)

grafo_bic = np.zeros((m, m), dtype=int)

for i in range(m):

for j in pais[i]:

grafo_bic[j, i] = 1 # j -> i

return grafo_bic

grafo_bic = busca_gulosa(dados, var_names)

# 2. Aprendizado por restrições (algoritmo PC simplificado)

# Teste qui-quadrado para independência condicional entre X e Y dado Z

def teste_chi2(dados, x, y, z_list):

# Construir tabela de contingência condicionada em Z

if not z_list:

# Teste marginal

tabela = np.zeros((2,2))

for i in range(len(dados)):

tabela[dados[i,x], dados[i,y]] += 1

chi2, p, dof, _ = chi2_contingency(tabela, correction=False)

return p

else:

# Para cada combinação de Z, fazer teste e combinar p-valores (Fisher)

from collections import defaultdict

from scipy.stats import combine_pvalues

pvals = []

for z_vals in set(tuple(dados[i, z_list]) for i in range(len(dados))):

mask = np.all(dados[:, z_list] == z_vals, axis=1)

if np.sum(mask) < 5:

continue

sub = dados[mask]

tabela = np.zeros((2,2))

for i in range(len(sub)):

tabela[sub[i,x], sub[i,y]] += 1

if np.sum(tabela) < 5:

continue

chi2, p, dof, _ = chi2_contingency(tabela, correction=False)

pvals.append(p)

if not pvals:

return 1.0

# Combinar p-valores (método de Fisher)

stat, p_combined = combine_pvalues(pvals, method='fisher')

return p_combined

def algoritmo_pc(dados, alpha=0.05):

m = dados.shape[1]

# Grafo não-direcionado inicial (completo)

adj = np.ones((m, m), dtype=int) - np.eye(m, dtype=int)

# Passo 1: remover arestas com testes de independência marginal e condicional

for ordem in range(0, m-1):

# Para cada par (x,y) ainda conectado

for x in range(m):

for y in range(m):

if x >= y: continue

if adj[x,y] == 0: continue

# Encontrar vizinhos comuns (excluindo x e y)

vizinhos = [v for v in range(m) if v != x and v != y and adj[x,v] and adj[y,v]]

if len(vizinhos) < ordem:

continue

# Testar todas as combinações de ordem variáveis

for subset in combinations(vizinhos, ordem):

p_val = teste_chi2(dados, x, y, list(subset))

if p_val > alpha: # independentes -> remover aresta

adj[x,y] = adj[y,x] = 0

break

# Orientar arestas (regras de colisor simples)

grafo_dir = np.zeros((m, m), dtype=int)

# Para cada tripla x - z - y, se x e y não são adjacentes, z é colisor

for x in range(m):

for z in range(m):

if x == z: continue

for y in range(m):

if y == x or y == z: continue

if adj[x,z] and adj[z,y] and not adj[x,y]:

# x -> z <- y (colisor)

grafo_dir[x,z] = 1

grafo_dir[y,z] = 1

# Para as arestas restantes, orientar sem criar ciclos (simplificado: usar ordem natural)

for x in range(m):

for y in range(m):

if x < y and adj[x,y] and grafo_dir[x,y]==0 and grafo_dir[y,x]==0:

# Orientar arbitrariamente (baseado em ordem)

grafo_dir[x,y] = 1

return grafo_dir

grafo_pc = algoritmo_pc(dados, alpha=0.05)

# Estrutura verdadeira

verdade = np.zeros((4,4), dtype=int)

# A->B (0->1), A->C (0->2), B->D (1->3), C->D (2->3)

verdade[0,1] = 1

verdade[0,2] = 1

verdade[1,3] = 1

verdade[2,3] = 1

print("Estrutura verdadeira (arestas):")

print(verdade)

print("\nEstrutura por BIC (busca gulosa):")

print(grafo_bic)

print("\nEstrutura por PC (restrições):")

print(grafo_pc)

# Plotar as três estruturas

def plot_graph(adj, title, ax):

G = nx.DiGraph()

m = adj.shape[0]

for i in range(m):

for j in range(m):

if adj[i,j]:

G.add_edge(var_names[i], var_names[j])

pos = nx.circular_layout(G)

nx.draw(G, pos, ax=ax, with_labels=True, node_color='lightblue',

node_size=800, font_weight='bold', arrowsize=20,

edge_color='gray', width=2)

ax.set_title(title)

fig, axes = plt.subplots(1, 3, figsize=(15, 4))

plot_graph(verdade, 'Estrutura Verdadeira', axes[0])

plot_graph(grafo_bic, 'Aprendizado por Escore (BIC)', axes[1])

plot_graph(grafo_pc, 'Aprendizado por Restrições (PC)', axes[2])

plt.tight_layout()

plt.show()

Este código implementa ambas as abordagens para aprendizado de estrutura. A busca por escore BIC encontra uma estrutura próxima da verdadeira. O algoritmo PC (restrições) também identifica corretamente as dependências. As diferenças ocorrem devido ao ruído e ao tamanho da amostra. Para iniciantes, este exemplo mostra as duas filosofias de aprendizado. Escolha entre escore e restrição depende do contexto e dos recursos.