Bayesian Separation Logic

# Bayesian Hoare triple: {P} c {Q}
# P, Q are predicates on distributions (not just values)

import random

def sample_bernoulli(p):
    return 1 if random.random() < p else 0

# Program: sample a coin bias, then observe heads
# {prior: bias ~ Uniform} coin; observe(heads) {posterior: bias shifted up}

# Simplified: discrete approximation
prior = {0.2: 0.2, 0.4: 0.2, 0.5: 0.2, 0.7: 0.2, 0.9: 0.2}

# Observe heads: weight by likelihood
def observe_heads(prior):
    weighted = {bias: prob * bias for bias, prob in prior.items()}
    total = sum(weighted.values())
    return {bias: p / total for bias, p in weighted.items()}

posterior = observe_heads(prior)
print("prior:    ", {k: round(v, 3) for k, v in prior.items()})
print("posterior:", {k: round(v, 3) for k, v in posterior.items()})
print()
print("After observing heads, higher biases are more probable.")
print("That's Bayesian conditioning — the program's semantics.")

# s-finite kernels: allow unnormalized measures
# After conditioning, weights don't sum to 1 until normalized

def s_finite_kernel(prior, likelihood, evidence):
    """Unnormalized posterior — an s-finite kernel"""
    unnorm = {h: ph * likelihood(h, evidence)
              for h, ph in prior.items()}
    return unnorm  # NOT normalized — that's the point

def normalize(kernel):
    total = sum(kernel.values())
    if total == 0:
        return None
    return {k: v / total for k, v in kernel.items()}

prior = {"fair": 0.5, "biased": 0.5}
def lk(h, e):
    if h == "fair": return 0.5
    return 0.8 if e == "heads" else 0.2

unnorm = s_finite_kernel(prior, lk, "heads")
norm = normalize(unnorm)
print(f"unnormalized: {unnorm}")
print(f"total weight: {sum(unnorm.values())}")
print(f"normalized:   {norm}")
print()
print("s-finite = work with unnormalized weights,")
print("normalize only when you need probabilities.")

# Frame rule: independent variables stay independent
# {P} c {Q} implies {P * R} c {Q * R}
# where * means "on independent variables"

import random

# Two independent coins: x ~ Bernoulli(0.5), y ~ Bernoulli(0.7)
# Program c: observe x = heads (only touches x)

# Before: x and y are independent
prior_x = {0: 0.5, 1: 0.5}
prior_y = {0: 0.3, 1: 0.7}

# Joint = product (independence)
joint_before = {(xv, yv): px * py
                for xv, px in prior_x.items()
                for yv, py in prior_y.items()}

# Observe x = 1 (only touches x)
def observe_x1(joint):
    filtered = {k: v for k, v in joint.items() if k[0] == 1}
    total = sum(filtered.values())
    return {k: v / total for k, v in filtered.items()}

joint_after = observe_x1(joint_before)

# Check: y's marginal is unchanged
marginal_y_before = {0: 0.3, 1: 0.7}
marginal_y_after = {}
for (x, y), p in joint_after.items():
    marginal_y_after[y] = marginal_y_after.get(y, 0) + p

print("y marginal before:", marginal_y_before)
print("y marginal after: ", {k: round(v, 4) for k, v in marginal_y_after.items()})
print()
print("y is unchanged — the frame rule in action.")
print("Independence (separation) is preserved.")

# COMP for Bayesian programs
# {P} c1 {R} and {R} c2 {Q} implies {P} c1;c2 {Q}

# c1: sample bias ~ Uniform{0.3, 0.5, 0.7}
# c2: observe 2 heads given bias

prior = {0.3: 1/3, 0.5: 1/3, 0.7: 1/3}

# c1: just the prior (identity)
after_c1 = dict(prior)

# c2: observe 2 heads — likelihood = bias^2
def observe_2_heads(dist):
    weighted = {b: p * b**2 for b, p in dist.items()}
    total = sum(weighted.values())
    return {b: p / total for b, p in weighted.items()}

after_c2 = observe_2_heads(after_c1)

print("prior:          ", {k: round(v, 3) for k, v in prior.items()})
print("after 2 heads:  ", {k: round(v, 3) for k, v in after_c2.items()})
print()
print("Higher biases now dominate — evidence compounds.")
print("COMP chains evidence accumulation.")