MachineSex/paper/pnas/make_figs.py
Giorgio Gilestro 96902e87f0 figures: publication-ready — unified, lettered, codename-free
paper/pnas/make_figs.py re-plots every panel directly from the committed
results artifacts into six single-file figures (figs/fig1..fig6.pdf):
no experiment codenames or suptitles (interpretation moved to captions),
bold panel letters, plain-language axis labels and legend entries, one
consistent style (8pt, no top/right spines). Panels: fig1 A-B (grounding
equilibrium + MNIST montage with its baked-in title cropped), fig2 A-B
(blending cancellation + Fisher-Muller), fig3 A-D (outbreeding, directed
recombination, mating breadth champion + diversity), fig4 A-C (society
ablation trajectories), fig5 A-F (speciation: analytic curve + cliff,
MLP decomposition + conflict sweep, LLM coherence + duration null),
fig6 A-D (seed-replicated merging, 7B-hard routing vs averaging,
predictive-test scatter, predictor comparison). build.py now places the
single PDFs; captions rewritten per lettered panel; in-text panel refs
updated (5B->5C-D, 5C->5E-F); stale stacked copies removed. Document
20pp -> 18pp.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01BkRLcc18rwT2Lysu6PbG7v
2026-09-07 09:09:46 +01:00

381 lines
17 KiB
Python
Raw Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

"""Publication figures for the PNAS draft — unified, lettered, codename-free.
Re-plots every panel directly from the committed results artifacts into six single-file figures
(figs/fig1.pdf .. fig6.pdf): no experiment codenames, no suptitles, no per-panel headline titles
(interpretation lives in the captions), bold panel letters, one consistent style. The per-experiment
figures under results/ remain the exploratory versions; these are the manuscript's.
Usage: python paper/pnas/make_figs.py
"""
from __future__ import annotations
import sys
from pathlib import Path
import matplotlib.pyplot as plt
import numpy as np
ROOT = Path(__file__).resolve().parents[2]
sys.path.insert(0, str(ROOT / "figures"))
sys.path.insert(0, str(ROOT / "src"))
import os
os.chdir(ROOT) # load_bundle uses repo-relative paths
from _figlib import load_bundle, mean_ci # noqa: E402
OUT = ROOT / "paper" / "pnas" / "figs"
plt.rcParams.update({
"font.size": 8, "axes.labelsize": 8.5, "legend.fontsize": 7, "legend.frameon": False,
"lines.markersize": 3.6, "axes.spines.top": False, "axes.spines.right": False,
})
def letter(ax, s, x=-0.14):
ax.text(x, 1.02, s, transform=ax.transAxes, fontsize=12, fontweight="bold", va="bottom")
def save(fig, name):
OUT.mkdir(exist_ok=True)
fig.savefig(OUT / f"{name}.pdf", bbox_inches="tight")
plt.close(fig)
print("wrote", OUT / f"{name}.pdf")
# ---------------------------------------------------------------- fig 1: grounding + MNIST
def fig1():
from knowledge.analysis import critical_grounding, reduce_to_stationary
from knowledge.metrics import heterozygosity
from knowledge.truth import make_true_distribution
df, cfg = load_bundle("results/E2")
n = cfg["dynamics"]["n"]
td = make_true_distribution(cfg["truth"]["K"], 1, "zipf", cfg["truth"]["tail_frac"],
cfg["truth"]["zipf_s"], 0, tail_threshold=cfg["truth"]["tail_threshold"])
H_star = heterozygosity(td.p_star)
last = int(cfg["generations"] * 0.8)
stat = df[df["generation"] >= last]
fig, axes = plt.subplots(1, 2, figsize=(10.6, 3.5), gridspec_kw={"width_ratios": [1, 1.35]})
ax = axes[0]
st = reduce_to_stationary(stat, value_col="heterozygosity", replicate_col="replicate", last_frac=1.0)
gg, Hm, Hci = mean_ci(stat, "g", "heterozygosity")
nz = gg > 0
ax.errorbar(gg[nz], Hm[nz], yerr=Hci[nz], fmt="o", color="#1f77b4", capsize=2, label="simulation")
ax.plot(gg[~nz], Hm[~nz], "o", mfc="white", mec="#1f77b4")
m_of_g = stat.groupby("g")["m"].first().to_numpy()
m_grid = np.linspace(0, m_of_g.max(), 400)
def H_eq(m):
m = np.asarray(m, float)
return np.where(m <= 0, 0.0, H_star * m * (2 * n + m - 1) / (n + 2 * n * m + m * m))
ax.plot(m_grid / (n + m_grid), H_eq(m_grid), "k--", lw=1, label="exact equilibrium")
ax.axhline(H_star, ls=":", color="gray", lw=1, label="source diversity $H^*$")
r = critical_grounding(st, H_star=H_star, frac=0.95, seed=7)
ax.axvspan(r["ci_low"], r["ci_high"], color="#d62728", alpha=0.15)
ax.axvline(r["g_star"], color="#d62728", lw=1.1,
label=f"95%-retention threshold $g\\approx{r['g_star']:.3f}$")
ax.set(xlabel="grounding fraction $g = m/(n+m)$", ylabel="stationary diversity $H$")
ax.legend()
letter(ax, "A")
ax = axes[1]
from PIL import Image
im = np.asarray(Image.open("results/mnist_collapse/mnist_montage.png"))
crop = int(im.shape[0] * 0.085) # remove the baked-in title band
ax.imshow(im[crop:], interpolation="bilinear")
ax.set_axis_off()
letter(ax, "B", x=-0.02)
save(fig, "fig1")
# ---------------------------------------------------------------- fig 2: blending vs union + FisherMuller
def fig2():
fig, axes = plt.subplots(1, 2, figsize=(10.6, 3.5))
df, _ = load_bundle("results/E4")
r0 = df[(df["g"] == 0.0) & (df["rho"] == 0.0)]
mx = r0.groupby("K_T")["surviving_max"].agg(["mean", "sem"])
mn = r0.groupby("K_T")["surviving_mean"].agg(["mean", "sem"])
ax = axes[0]
ax.errorbar(mx.index, mx["mean"], yerr=1.96 * mx["sem"], fmt="-o", color="#1f77b4",
capsize=2, label="union operator (strongest source)")
ax.errorbar(mn.index, mn["mean"], yerr=1.96 * mn["sem"], fmt="--s", color="#d62728",
capsize=2, label="output-mean (blending)")
ax.set(xlabel="number of parents", ylabel="rare capabilities surviving in the child",
xticks=sorted(r0["K_T"].unique()))
ax.legend()
letter(ax, "A")
df8, cfg8 = load_bundle("results/E8")
L = cfg8["society"]["L"]
d0 = df8[df8["rho"] == 0.0]
ax = axes[1]
for col, c, lab in [("best_parent", "#7f7f7f", "best single parent"),
("average", "#1f77b4", "blended average"),
("sexual", "#d62728", "recombined offspring")]:
k, m, ci = mean_ci(d0, "K_T", col)
ax.errorbar(k, m, yerr=ci, fmt="-o", color=c, capsize=2, label=lab)
ax.axhline(L, ls=":", color="green", lw=1, label="optimum")
ax.set(xlabel="number of parents", ylabel="offspring capability")
ax.legend()
letter(ax, "B")
save(fig, "fig2")
# ---------------------------------------------------------------- fig 3: rugged landscapes
def fig3():
fig, axes = plt.subplots(2, 2, figsize=(10.6, 6.8))
df9, _ = load_bundle("results/E9")
Ks = sorted(df9["K"].unique())
colors = plt.cm.viridis(np.linspace(0, 0.85, len(Ks)))
bp = df9.groupby("K")["best_parent"].mean()
ax = axes[0, 0]
for K, c in zip(Ks, colors):
s = df9[df9["K"] == K].groupby("rate")["mean_offspring"].mean() - bp[K]
ax.plot(s.index, s.values, "-o", color=c, label=f"$K$={K}")
ax.axhline(0, ls=":", color="gray", lw=1)
ax.set(xlabel="recombination rate", ylabel="mean offspring best parent")
ax.legend(title="ruggedness", ncol=2)
letter(ax, "A")
df10, _ = load_bundle("results/E10")
ax = axes[0, 1]
for col, c, lab in [("global_opt", "green", "global optimum"),
("directed_sex", "#d62728", "screened recombination (directed)"),
("best_parent", "#7f7f7f", "best single parent"),
("random_sex", "#1f77b4", "blind recombination")]:
k, m, ci = mean_ci(df10, "K", col)
if col == "global_opt":
ax.plot(k, m, ":", color=c, label=lab)
else:
ax.errorbar(k, m, yerr=ci, fmt="-o", color=c, capsize=2, label=lab)
ax.set(xlabel="landscape ruggedness $K$", ylabel="offspring capability")
ax.legend()
letter(ax, "B")
df14, _ = load_bundle("results/E14")
last = df14[df14["generation"] == df14["generation"].max()].copy()
last["best_n"] = last["best_fitness"] / last["global_opt"]
K14 = sorted(last["K"].unique())
cmap = plt.get_cmap("viridis")
c14 = {K: cmap(i / max(1, len(K14) - 1)) for i, K in enumerate(K14)}
for ax, col, ylab, L in [(axes[1, 0], "best_n", "best fitness / optimum", "C"),
(axes[1, 1], "diversity", "population diversity", "D")]:
for K in K14:
g = last[last["K"] == K].groupby("breadth")[col].agg(["mean", "sem"]).reset_index()
ax.errorbar(g["breadth"], g["mean"], yerr=1.96 * g["sem"].fillna(0), fmt="-o",
color=c14[K], capsize=2, label=f"$K$={K}")
ax.set_xscale("log")
ax.set(xlabel="mate-pool breadth (monogamous → panmictic)", ylabel=ylab)
ax.legend(title="ruggedness")
letter(ax, L)
save(fig, "fig3")
# ---------------------------------------------------------------- fig 4: the society
def fig4():
df, _ = load_bundle("results/E11")
arms = [("full", "#2ca02c", "full system"),
("no_sex", "#ff7f0e", "no recombination"),
("no_diversity", "#9467bd", "no diversity preservation"),
("no_grounding", "#d62728", "no grounded evaluation")]
arms = [a for a in arms if a[0] in set(df["arm"].unique())]
g_opt = df["global_opt"].mean()
fig, axes = plt.subplots(1, 3, figsize=(11.4, 3.2))
panels = [("best_fitness", "best real fitness", "A"),
("diversity", "population diversity", "B"),
("conformity_true_gap", "conformity true fitness", "C")]
for ax, (col, ylab, L) in zip(axes, panels):
for name, c, lab in arms:
sub = df[df["arm"] == name]
g, m, ci = mean_ci(sub, "generation", col)
ax.plot(g, m, "-", color=c, lw=1.6, label=lab)
ax.fill_between(g, m - ci, m + ci, color=c, alpha=0.15)
if col == "best_fitness":
ax.axhline(g_opt, ls=":", color="gray", lw=1, label="global optimum")
ax.legend()
ax.set(xlabel="generation", ylabel=ylab)
letter(ax, L)
save(fig, "fig4")
# ---------------------------------------------------------------- fig 5: speciation, three tiers
def fig5():
fig, axes = plt.subplots(2, 3, figsize=(11.4, 6.6))
bdm, _ = load_bundle("results/E12")
rhos = sorted(bdm["rho"].unique())
colors = plt.cm.viridis(np.linspace(0.15, 0.85, len(rhos)))
def agg(df, keys, value):
g = df.groupby(keys)[value].agg(["mean", "std", "count"]).reset_index()
g["se"] = g["std"] / np.sqrt(g["count"].clip(lower=1))
return g
ax = axes[0, 0]
par = agg(bdm, "divergence", "parent_fitness")
ax.plot(par["divergence"], par["mean"], "k--", lw=1.3, label="parents")
for rho, c in zip(rhos, colors):
g = agg(bdm[bdm["rho"] == rho], "divergence", "offspring_fitness")
ax.plot(g["divergence"], g["mean"], "-o", color=c, label=f"hybrid, density {rho:g}")
ax.fill_between(g["divergence"], g["mean"] - g["se"], g["mean"] + g["se"], color=c, alpha=0.15)
ax.axhline(0, color="#999", lw=0.7, ls=":")
ax.set(xlabel="parental divergence (substitutions)", ylabel="fitness")
ax.legend()
letter(ax, "A")
ax = axes[0, 1]
for rho, c in zip(rhos, colors):
g = agg(bdm[bdm["rho"] == rho], "divergence", "isolation")
ax.plot(g["divergence"], g["mean"], "-o", color=c, label=f"{rho:g}")
ax.set(xlabel="parental divergence (substitutions)", ylabel="P(hybrid inviable)", ylim=(-0.02, 1.02))
ax.legend(title="incompatibility density")
letter(ax, "B")
dec, _ = load_bundle("results/speciation_real")
order = [c for c in ["shared", "independent", "conflict"] if c in set(dec["condition"])]
g = dec.groupby("condition").agg(naive=("barrier_naive", "mean"),
res=("residual_scale", "mean")).reindex(order)
ax = axes[0, 2]
x = np.arange(len(order)); w = 0.38
ax.bar(x - w / 2, g["naive"], w, color="#9ecae1", label="before alignment")
ax.bar(x + w / 2, g["res"], w, color="#d62728", label="after alignment (residual)")
ax.set_xticks(x)
ax.set_xticklabels(["same task,\nshared start", "same task,\ndifferent start", "conflicting\ntasks"])
ax.set(ylabel="merge error barrier")
ax.legend()
letter(ax, "C")
cliff, _ = load_bundle("results/speciation_real_cliff")
cg = cliff.groupby("conflict_frac").agg(res=("residual_scale", "mean"),
hyb=("acc_merge_scale", "mean")).reset_index()
ax = axes[1, 0]
ax.plot(cg["conflict_frac"], cg["res"], "-o", color="#d62728", label="residual barrier")
ax2 = ax.twinx()
ax2.plot(cg["conflict_frac"], cg["hyb"], "-s", color="#2c7fb8", label="merged-model accuracy")
ax2.set_ylabel("merged accuracy", color="#2c7fb8")
ax2.tick_params(axis="y", labelcolor="#2c7fb8")
ax2.spines["right"].set_visible(True)
ax.set(xlabel="fraction of classes in conflict", ylabel="residual barrier")
l1, la1 = ax.get_legend_handles_labels(); l2, la2 = ax2.get_legend_handles_labels()
ax.legend(l1 + l2, la1 + la2, loc="center left")
letter(ax, "D")
rep, _ = load_bundle("results/llm_speciation")
def series(df, mode, model, metric):
sub = df[(df["mode"] == mode) & (df["model"] == model) & (df["metric"] == metric)]
g = sub.groupby("x")["accuracy"].mean().reset_index()
return g["x"], g["accuracy"]
ax = axes[1, 1]
x_, y_ = series(rep, "conflict", "parent_a", "ambig_asc")
ax.plot(x_, y_, "--o", color="#9ecae1", label="parent A, own convention")
x_, y_ = series(rep, "conflict", "parent_b", "ambig_desc")
ax.plot(x_, y_, "--o", color="#a1d99b", label="parent B, own convention")
x_, y_ = series(rep, "conflict", "merge_soup", "coherence")
ax.plot(x_, y_, "-s", color="#d62728", label="merge, best convention")
ax.set(xlabel="fraction of training in conflict", ylabel="accuracy, shared prompts")
ax.legend()
letter(ax, "E")
ax = axes[1, 2]
x_, y_ = series(rep, "duration", "merge_soup", "mean_private")
ax.plot(x_, y_, "-o", color="#d62728", label="merged model")
x_, y_ = series(rep, "duration", "parent_a", "strings")
ax.plot(x_, y_, "--o", color="#9ecae1", label="parent A, own task")
x_, y_ = series(rep, "duration", "parent_b", "arith")
ax.plot(x_, y_, "--o", color="#a1d99b", label="parent B, own task")
ax.set(xlabel="specialist training (epochs)", ylabel="accuracy", ylim=(0, 1.02))
ax.legend()
letter(ax, "F")
save(fig, "fig5")
# ---------------------------------------------------------------- fig 6: the language-model tier
def fig6():
import pandas as pd
from scipy.stats import spearmanr
fig, axes = plt.subplots(2, 2, figsize=(10.6, 6.6))
dfm, _ = load_bundle("results/llm_merge_seeds")
specs = sorted(m for m in dfm["model"].unique() if m.startswith("spec_"))
rows = []
for s, sub in dfm.groupby("seed"):
ov = {m: sub[(sub["model"] == m) & (sub["metric"] == "overall")]["accuracy"].mean() for m in specs}
b = sub[sub["model"] == max(ov, key=ov.get)].copy(); b["model"] = "best_specialist"
rows.append(b)
dfm = pd.concat([dfm] + rows, ignore_index=True)
models = ["base", "best_specialist", "merge_soup", "merge_ties"]
labels = ["base", "best\nspecialist", "merged\n(average)", "merged\n(interference-aware)"]
ax = axes[0, 0]
x = np.arange(len(models))
for off, metric, c, lab in ((-0.19, "overall", "#2c7fb8", "overall"),
(0.19, "worst_family", "#d62728", "worst task family")):
vals, errs = [], []
for m in models:
v = dfm[(dfm["model"] == m) & (dfm["metric"] == metric)].groupby("seed")["accuracy"].mean()
vals.append(v.mean()); errs.append(1.96 * v.std(ddof=1) / max(1, np.sqrt(len(v))))
ax.bar(x + off, vals, 0.36, yerr=errs, capsize=2, color=c, label=lab)
ax.set_xticks(x); ax.set_xticklabels(labels, fontsize=7)
ax.set(ylabel="verifier accuracy")
ax.legend()
letter(ax, "A")
df7, _ = load_bundle("results/llm_moe_hard_hpc")
def acc7(model, metric):
r = df7[(df7["model"] == model) & (df7["metric"] == metric)]["accuracy"]
return float(r.iloc[0]) if len(r) else np.nan
specs7 = sorted(m for m in df7["model"].unique() if m.startswith("spec_"))
best7 = max(specs7, key=lambda m: acc7(m, "overall"))
models7 = [best7, "merge_soup", "merge_ties", "moe_oracle"]
labels7 = ["best\nspecialist", "merged\n(average)", "merged\n(interference-aware)", "routed\n(kept separate)"]
ax = axes[0, 1]
x = np.arange(len(models7))
for off, metric, c, lab in ((-0.19, "overall", "#2c7fb8", "overall"),
(0.19, "worst_family", "#d62728", "worst task family")):
ax.bar(x + off, [acc7(m, metric) for m in models7], 0.36, color=c, label=lab)
ax.set_xticks(x); ax.set_xticklabels(labels7, fontsize=7)
ax.set(ylabel="verifier accuracy")
ax.legend()
letter(ax, "B")
a = pd.read_parquet("results/llm_epistasis/results.parquet")
b = pd.read_parquet("results/llm_epistasis_compat/results.parquet")
dfe = pd.concat([a, b], ignore_index=True)
ax = axes[1, 0]
for mode, c, mk, lab in (("conflict", "#d62728", "o", "conflicting conventions"),
("duration", "#2c7fb8", "s", "divergence only"),
("compat", "#41ab5d", "^", "overlap, no conflict")):
sub = dfe[dfe["mode"] == mode]
ax.scatter(sub["epi_conf"], sub["merge_penalty"], c=c, marker=mk, s=26, alpha=0.75, label=lab)
ax.axhline(0, color="#999", lw=0.6)
ax.set(xlabel="pre-merge functional conflict (confidence-weighted)",
ylabel="merge penalty")
ax.legend()
letter(ax, "C")
preds = [("dis_raw", "raw\ndisagreement"), ("epi_conf", "conf-weighted\nconflict"),
("grad_cos", "gradient\nalignment"), ("delta_cos", "weight\ncosine"),
("delta_l2", "weight\ndistance"), ("cross_perf", "cross-task\naccuracy")]
ax = axes[1, 1]
rhos_ = [abs(spearmanr(dfe[c], dfe["merge_penalty"])[0]) for c, _ in preds]
cols = ["#fc9272", "#d62728", "#9ecae1", "#9ecae1", "#9ecae1", "#9ecae1"]
ax.bar(np.arange(len(preds)), rhos_, 0.6, color=cols)
ax.set_xticks(np.arange(len(preds)))
ax.set_xticklabels([l for _, l in preds], fontsize=6.5)
ax.set(ylabel="|Spearman ρ| vs merge penalty", ylim=(0, 0.8))
letter(ax, "D")
save(fig, "fig6")
if __name__ == "__main__":
for f in (fig1, fig2, fig3, fig4, fig5, fig6):
f()