exp03_ledger.py
#!/usr/bin/env python3
"""Experiment 0.3: ledger calibration.
Eight systems solve goal-reaching on a ring with barriers. For each: competency
in bits (episode-level and behaviour-level), substrate bias, and three effort
terms (code, parameters, training signal). Four accounting rules combine them
into a residual. plan.md fixes the pass criterion and the expected outcome.
Usage: exp03_ledger.py [--episodes 4000] [--train 100000] [--out results/exp03]
exp03_ledger.py --selftest
"""
import argparse
import csv
import inspect
import lzma
import math
import os
import random
import sys
import zlib
from collections import Counter
from itertools import combinations
L = 12 # ring length
T = 10 # step budget
ACTIONS = (-1, 0, 1)
MAX_BARRIERS = 2
SEED_BITS = 32
# --- environment ---
def states():
"""Every (pos, goal, barriers) with a path of length <= T."""
cells = range(L)
for k in range(MAX_BARRIERS + 1):
for bars in combinations(cells, k):
bset = frozenset(bars)
for pos in cells:
if pos in bset:
continue
for goal in cells:
if goal in bset or dist(pos, goal, bset) > T:
continue
yield (pos, goal, bset)
def step(pos, a, bars):
nxt = (pos + a) % L
return pos if nxt in bars else nxt
def dist(pos, goal, bars):
"""Shortest path length on the ring avoiding bars; L if unreachable."""
if pos == goal:
return 0
seen, frontier, d = {pos}, [pos], 0
while frontier:
d += 1
nxt = []
for p in frontier:
for a in (-1, 1):
q = (p + a) % L
if q in bars or q in seen:
continue
if q == goal:
return d
seen.add(q)
nxt.append(q)
frontier = nxt
return L
def optimal(pos, goal, bars):
"""Canonical optimal action: 0 at goal, else the direction of the shortest
path, +1 on ties."""
if pos == goal:
return 0
best, best_d = 0, L + 1
for a in (1, -1):
q = step(pos, a, bars)
if q == pos:
continue
d = dist(q, goal, bars)
if d < best_d:
best, best_d = a, d
return best
# --- systems: policy(pos, goal, bars) -> action ---
def pol_solver(pos, goal, bars):
return optimal(pos, goal, bars)
def make_table():
return {s: optimal(*s) for s in states()}
def pol_greedy(pos, goal, bars):
if pos == goal:
return 0
cw = (goal - pos) % L
a = 1 if cw <= L - cw else -1
if step(pos, a, bars) == pos:
a = -a
return a
def pol_greedy_buggy(pos, goal, bars):
if pos == goal:
return 0
cw = (goal - pos) % L
return 1 if cw <= L - cw else -1
def q_train(episodes, rng, alpha=0.2, gamma=0.95, eps=0.2):
q = {}
for _ in range(episodes):
pos, goal, bars = random_state(rng)
for _ in range(T):
if pos == goal:
break
key = (pos, goal, bars)
row = q.setdefault(key, [0.0, 0.0, 0.0])
i = rng.randrange(3) if rng.random() < eps else row.index(max(row))
nxt = step(pos, ACTIONS[i], bars)
r = 1.0 if nxt == goal else -0.01
nrow = q.get((nxt, goal, bars), [0.0, 0.0, 0.0])
row[i] += alpha * (r + gamma * max(nrow) - row[i])
pos = nxt
return q
def pol_from_q(q, rng):
def pol(pos, goal, bars):
row = q.get((pos, goal, bars), [0.0, 0.0, 0.0])
m = max(row)
idx = [i for i, v in enumerate(row) if v == m]
return ACTIONS[rng.choice(idx)]
return pol
def make_reservoir(rng):
return {s: rng.choice(ACTIONS) for s in states()}
def random_state(rng):
while True:
k = rng.randrange(MAX_BARRIERS + 1)
bars = frozenset(rng.sample(range(L), k))
free = [c for c in range(L) if c not in bars]
pos, goal = rng.choice(free), rng.choice(free)
if dist(pos, goal, bars) <= T:
return pos, goal, bars
# --- measurement ---
def success_rate(pol, episodes, rng):
ok = 0
for _ in range(episodes):
pos, goal, bars = random_state(rng)
for _ in range(T):
if pos == goal:
break
pos = step(pos, pol(pos, goal, bars), bars)
ok += pos == goal
return ok / episodes
def mutual_info(pairs):
n = len(pairs)
jx, jy, jxy = Counter(), Counter(), Counter()
for x, y in pairs:
jx[x] += 1
jy[y] += 1
jxy[(x, y)] += 1
return sum(c / n * math.log2(c * n / (jx[x] * jy[y])) for (x, y), c in jxy.items())
def behaviour_bits(pol, state_list):
pairs = [(pol(*s), optimal(*s)) for s in state_list]
return len(state_list) * mutual_info(pairs)
def bits(obj):
"""Compressed size in bits; the smaller of zlib and lzma."""
raw = repr(obj).encode() if not isinstance(obj, bytes) else obj
return 8 * min(len(zlib.compress(raw, 9)), len(lzma.compress(raw)))
def policy_bytes(pol, state_list):
"""Canonical parameter encoding: one action byte per enumerated state."""
return bytes(ACTIONS.index(pol(*s)) for s in state_list)
def q_bytes(q, state_list):
"""Q-values quantised to 8 bits each, in canonical state order."""
vals = [v for s in state_list for v in q.get(s, [0.0, 0.0, 0.0])]
lo, hi = min(vals), max(vals)
scale = 255 / (hi - lo) if hi > lo else 0.0
return bytes(int((v - lo) * scale) for v in vals)
def code_bits(*fns):
return bits("".join(inspect.getsource(f) for f in fns).encode())
def run(episodes, train, out_dir):
rng = random.Random(7)
state_list = list(states())
p_rand = success_rate(lambda p, g, b: rng.choice(ACTIONS), episodes, random.Random(1))
env_code = (step, dist, optimal)
q = q_train(train, random.Random(2))
feedback_bits = 1 + math.log2(T) # success flag plus step count per episode
table = make_table()
reservoir = make_reservoir(random.Random(3))
systems = [
# name, policy, class, E_code fns, params, E_train
("random", lambda p, g, b: rng.choice(ACTIONS), "random", (), None, 0),
("solver_program", pol_solver, "solver", (pol_solver, *env_code), None, 0),
("solver_table", lambda p, g, b: table[(p, g, b)], "solver", (step,), "policy", 0),
("greedy_correct", pol_greedy, "greedy", (pol_greedy, step), None, 0),
("greedy_buggy", pol_greedy_buggy, "greedy", (pol_greedy_buggy,), None, 0),
("q_trained", pol_from_q(q, random.Random(4)), "q", (q_train, pol_from_q, step, random_state), "q", train * feedback_bits),
("q_untrained", pol_from_q({}, random.Random(5)), "q", (pol_from_q,), "policy", 0),
("reservoir_untrained", lambda p, g, b: reservoir[(p, g, b)], "reservoir", (make_reservoir,), "policy", 0),
]
rows = []
e_params_full = 0
for name, pol, cls, fns, params, e_train in systems:
p = success_rate(pol, episodes, random.Random(11))
c_ep = math.log2(p / p_rand) if p > 0 else float("-inf")
c_beh = behaviour_bits(pol, state_list)
if params == "policy":
e_params = bits(policy_bytes(pol, state_list))
elif params == "q":
e_params = bits(policy_bytes(pol, state_list))
e_params_full = bits(q_bytes(q, state_list))
else:
e_params = 0
rows.append({
"system": name,
"class": cls,
"p_success": p,
"C_ep": c_ep,
"C_beh": c_beh,
"E_code": code_bits(*fns) if fns else 0,
"E_params": e_params,
"E_params_q8": e_params_full if params == "q" else "",
"E_seed": SEED_BITS if params is not None or cls == "random" else 0,
"E_train": e_train,
})
# substrate bias: competency of the untrained instance of the class
base = {"q": "q_untrained", "reservoir": "reservoir_untrained"}
by_name = {r["system"]: r for r in rows}
for r in rows:
r["B"] = by_name[base[r["class"]]]["C_beh"] if r["class"] in base else 0.0
r["R1"] = r["C_beh"] - r["B"] - (r["E_code"] + r["E_params"])
r["R2"] = r["C_beh"] - r["B"] - (r["E_code"] + r["E_seed"] + r["E_train"])
r["R3"] = r["C_beh"] - r["B"] - (r["E_code"] + r["E_params"] + r["E_train"])
# R4: leave-one-out linear fit of C_beh on R3's effort
pts = [(r["E_code"] + r["E_params"] + r["E_train"], r["C_beh"]) for r in rows]
for i, r in enumerate(rows):
xs = [x for j, (x, _) in enumerate(pts) if j != i]
ys = [y for j, (_, y) in enumerate(pts) if j != i]
mx, my = sum(xs) / len(xs), sum(ys) / len(ys)
sxx = sum((x - mx) ** 2 for x in xs)
slope = sum((x - mx) * (y - my) for x, y in zip(xs, ys)) / sxx if sxx else 0.0
r["R4"] = r["C_beh"] - (my + slope * (pts[i][0] - mx))
os.makedirs(out_dir, exist_ok=True)
with open(os.path.join(out_dir, "ledger.csv"), "w", newline="") as f:
w = csv.DictWriter(f, fieldnames=list(rows[0].keys()))
w.writeheader()
w.writerows(rows)
print(f"states={len(state_list)} p_random={p_rand:.3f} episodes={episodes} train={train}")
hdr = f"{'system':20} {'p':>5} {'C_ep':>6} {'C_beh':>8} {'B':>7} {'E_code':>7} {'E_par':>7} {'E_train':>8} {'R1':>9} {'R2':>9} {'R3':>9} {'R4':>9}"
print(hdr)
for r in rows:
print(
f"{r['system']:20} {r['p_success']:5.3f} {r['C_ep']:6.2f} {r['C_beh']:8.0f} {r['B']:7.0f} "
f"{r['E_code']:7.0f} {r['E_params']:7.0f} {r['E_train']:8.0f} {r['R1']:9.0f} {r['R2']:9.0f} {r['R3']:9.0f} {r['R4']:9.0f}"
)
return rows
def selftest():
assert dist(0, 6, frozenset()) == 6
assert dist(0, 1, frozenset({1})) == L # unreachable
assert dist(0, 2, frozenset({1})) == L - 2 # forced the long way
assert optimal(0, 2, frozenset({1})) == -1
assert optimal(3, 3, frozenset()) == 0
assert pol_greedy_buggy(0, 2, frozenset({1})) == 1 # walks into the wall
assert pol_greedy(0, 2, frozenset({1})) == -1
n = sum(1 for _ in states())
assert n > 5000, n
assert abs(mutual_info([(0, 0), (1, 1)] * 10) - 1.0) < 1e-9
assert mutual_info([(0, 0), (0, 1)] * 10) == 0.0
rng = random.Random(0)
q = q_train(3000, rng)
p = success_rate(pol_from_q(q, rng), 300, rng)
assert p > 0.5, p
print("selftest ok")
def main():
ap = argparse.ArgumentParser()
ap.add_argument("--episodes", type=int, default=4000)
ap.add_argument("--train", type=int, default=100000)
ap.add_argument("--out", default=os.path.join(os.path.dirname(__file__), "results", "exp03"))
ap.add_argument("--selftest", action="store_true")
a = ap.parse_args()
if a.selftest:
selftest()
return
run(a.episodes, a.train, a.out)
if __name__ == "__main__":
main()