From 546d26a6759eddd5ef8723cbd4fec4627765a9ed Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 10:10:57 +0800 Subject: [PATCH 01/25] Add certified maximin cutting-plane solver --- src/bayesian_ach/design_certificate.py | 356 +++++++++++++++++++++++++ 1 file changed, 356 insertions(+) create mode 100644 src/bayesian_ach/design_certificate.py diff --git a/src/bayesian_ach/design_certificate.py b/src/bayesian_ach/design_certificate.py new file mode 100644 index 0000000..ec3c46c --- /dev/null +++ b/src/bayesian_ach/design_certificate.py @@ -0,0 +1,356 @@ +"""Certified cutting-plane maximin allocation over a finite design grid.""" + +from __future__ import annotations + +import math +from dataclasses import dataclass +from typing import Any + +import numpy as np +from numpy.typing import NDArray +from scipy.optimize import Bounds, LinearConstraint, milp + +from bayesian_ach.design_geometry import pairwise_residual_matrix +from bayesian_ach.design_optimizer import optimize_maximin_design + + +@dataclass(frozen=True, slots=True) +class CertifiedDesignConfig: + """Configuration for continuous or exact-integer maximin certification.""" + + budget: int = 60 + max_point_fraction: float = 0.15 + integer: bool = True + absolute_gap_tolerance: float = 1.0e-8 + relative_gap_tolerance: float = 1.0e-7 + cut_violation_tolerance: float = 1.0e-9 + max_iterations: int = 100 + master_time_limit_s: float = 120.0 + master_mip_relative_gap: float = 1.0e-9 + + def validate(self, point_count: int, candidate_count: int) -> None: + if self.budget < candidate_count + 1: + raise ValueError("budget must exceed the number of candidate signals") + if point_count < 1: + raise ValueError("point_count must be positive") + if not 0.0 < self.max_point_fraction <= 1.0: + raise ValueError("max_point_fraction must lie in (0, 1]") + if self.absolute_gap_tolerance < 0.0: + raise ValueError("absolute_gap_tolerance must be nonnegative") + if self.relative_gap_tolerance < 0.0: + raise ValueError("relative_gap_tolerance must be nonnegative") + if self.cut_violation_tolerance < 0.0: + raise ValueError("cut_violation_tolerance must be nonnegative") + if self.max_iterations < 1: + raise ValueError("max_iterations must be positive") + if self.master_time_limit_s <= 0.0: + raise ValueError("master_time_limit_s must be positive") + if self.master_mip_relative_gap < 0.0: + raise ValueError("master_mip_relative_gap must be nonnegative") + maximum_count = max(1, int(math.ceil(self.max_point_fraction * self.budget))) + if point_count * maximum_count < self.budget: + raise ValueError("allocation cap makes the requested budget infeasible") + + +@dataclass(frozen=True, slots=True) +class MaximinCertificate: + """Certified objective bounds and the best feasible allocation found.""" + + allocation: NDArray[np.float64] + integer: bool + certified: bool + lower_bound: float + upper_bound: float + absolute_gap: float + relative_gap: float + heuristic_lower_bound: float + iterations: int + cut_count: int + last_master_status: int + last_master_message: str + trace: tuple[dict[str, Any], ...] + + +@dataclass(frozen=True, slots=True) +class _PairCut: + generator: int + alternative: int + intercept: float + slope: float + losses: NDArray[np.float64] + residual: float + + +def _pairwise_oracle( + signals: NDArray[np.float64], + allocation: NDArray[np.float64], +) -> tuple[float, tuple[_PairCut, ...]]: + total = float(np.sum(allocation)) + if not np.isfinite(total) or total <= 0.0: + raise ValueError("allocation must have positive finite mass") + weights = allocation / total + root_weight = np.sqrt(weights) + cuts: list[_PairCut] = [] + minimum = math.inf + for generator in range(signals.shape[1]): + response = signals[:, generator] + for alternative in range(signals.shape[1]): + if generator == alternative: + continue + design = np.column_stack((np.ones(signals.shape[0]), signals[:, alternative])) + weighted_design = design * root_weight[:, None] + weighted_response = response * root_weight + coefficients, _, _, _ = np.linalg.lstsq( + weighted_design, + weighted_response, + rcond=None, + ) + residual_vector = response - design @ coefficients + losses = np.asarray(residual_vector**2, dtype=float) + residual = float(weights @ losses) + minimum = min(minimum, residual) + cuts.append( + _PairCut( + generator=generator, + alternative=alternative, + intercept=float(coefficients[0]), + slope=float(coefficients[1]), + losses=losses, + residual=residual, + ) + ) + return float(minimum), tuple(cuts) + + +def _validate_initial_allocation( + allocation: NDArray[np.float64], + *, + point_count: int, + config: CertifiedDesignConfig, + maximum_count: int, +) -> NDArray[np.float64]: + result = np.asarray(allocation, dtype=float) + if result.shape != (point_count,) or not np.all(np.isfinite(result)): + raise ValueError("initial_allocation must be a finite vector matching the grid") + if np.any(result < -1.0e-9) or np.any(result > maximum_count + 1.0e-9): + raise ValueError("initial_allocation violates its bounds") + if not math.isclose(float(np.sum(result)), config.budget, abs_tol=1.0e-7): + raise ValueError("initial_allocation must sum to budget") + if config.integer and not np.allclose(result, np.rint(result), atol=1.0e-9): + raise ValueError("integer certification requires integer initial_allocation") + return np.rint(result) if config.integer else result + + +def _master_upper_bound(result: Any, *, integer: bool, incumbent: float) -> float: + if bool(result.success): + return incumbent + if integer: + dual = getattr(result, "mip_dual_bound", None) + if dual is not None and np.isfinite(float(dual)): + return float(-float(dual)) + return math.inf + + +def certify_maximin_design( + standardized_signals: NDArray[np.float64], + config: CertifiedDesignConfig | None = None, + *, + initial_allocation: NDArray[np.float64] | None = None, +) -> MaximinCertificate: + """Certify continuous or integer maximin residual geometry with OLS cuts. + + Each ordered residual is the infimum, over an intercept and slope, of a + loss linear in the allocation weights. A finite set of such losses defines + a master LP/MILP upper bound. Weighted least squares at each master + solution supplies a separating cut and a feasible lower bound. + """ + + signals = np.asarray(standardized_signals, dtype=float) + if signals.ndim != 2 or signals.shape[1] < 2 or not np.all(np.isfinite(signals)): + raise ValueError("standardized_signals must be a finite two-dimensional matrix") + config = CertifiedDesignConfig() if config is None else config + config.validate(signals.shape[0], signals.shape[1]) + maximum_count = max(1, int(math.ceil(config.max_point_fraction * config.budget))) + + if initial_allocation is None: + heuristic = optimize_maximin_design( + signals, + config.budget, + max_point_fraction=config.max_point_fraction, + ) + initial = np.asarray(heuristic.counts, dtype=float) + else: + initial = np.asarray(initial_allocation, dtype=float) + initial = _validate_initial_allocation( + initial, + point_count=signals.shape[0], + config=config, + maximum_count=maximum_count, + ) + best_lower, initial_cuts = _pairwise_oracle(signals, initial) + best_allocation = initial.copy() + heuristic_lower = best_lower + + cuts: list[_PairCut] = [] + cut_keys: set[tuple[int, int, float, float]] = set() + + def add_cut(cut: _PairCut) -> bool: + key = ( + cut.generator, + cut.alternative, + round(cut.intercept, 12), + round(cut.slope, 12), + ) + if key in cut_keys: + return False + cut_keys.add(key) + cuts.append(cut) + return True + + for cut in initial_cuts: + add_cut(cut) + + variable_count = signals.shape[0] + 1 + lower_bounds = np.zeros(variable_count, dtype=float) + upper_bounds = np.full(variable_count, float(maximum_count), dtype=float) + upper_bounds[-1] = float(np.max(signals**2)) + objective = np.zeros(variable_count, dtype=float) + objective[-1] = -1.0 + integrality = np.zeros(variable_count, dtype=np.int32) + if config.integer: + integrality[:-1] = 1 + + global_upper = float(upper_bounds[-1]) + trace: list[dict[str, Any]] = [] + last_status = -1 + last_message = "master was not run" + certified = False + + for iteration in range(1, config.max_iterations + 1): + rows = np.zeros((len(cuts) + 1, variable_count), dtype=float) + lower = np.full(len(cuts) + 1, -np.inf, dtype=float) + upper = np.zeros(len(cuts) + 1, dtype=float) + rows[0, :-1] = 1.0 + lower[0] = float(config.budget) + upper[0] = float(config.budget) + for row_index, cut in enumerate(cuts, start=1): + rows[row_index, :-1] = -cut.losses + rows[row_index, -1] = float(config.budget) + + options: dict[str, float | bool] = { + "presolve": True, + "time_limit": float(config.master_time_limit_s), + } + if config.integer: + options["mip_rel_gap"] = float(config.master_mip_relative_gap) + result = milp( + objective, + integrality=integrality, + bounds=Bounds(lower_bounds, upper_bounds), + constraints=LinearConstraint(rows, lower, upper), + options=options, + ) + last_status = int(result.status) + last_message = str(result.message) + if result.x is None: + trace.append( + { + "iteration": iteration, + "status": last_status, + "message": last_message, + "cut_count": len(cuts), + "master_incumbent": None, + "master_upper_bound": global_upper, + "best_lower_bound": best_lower, + "absolute_gap": max(0.0, global_upper - best_lower), + "violated_pair_count": None, + "new_cut_count": 0, + } + ) + break + + master_incumbent = float(result.x[-1]) + master_upper = _master_upper_bound( + result, + integer=config.integer, + incumbent=master_incumbent, + ) + global_upper = min(global_upper, master_upper) + candidate = np.asarray(result.x[:-1], dtype=float) + if config.integer: + candidate = np.rint(candidate) + candidate_lower, oracle_cuts = _pairwise_oracle(signals, candidate) + if candidate_lower > best_lower: + best_lower = candidate_lower + best_allocation = candidate.copy() + + absolute_gap = max(0.0, global_upper - best_lower) + relative_gap = absolute_gap / max(1.0, abs(best_lower)) + violated = [ + cut + for cut in oracle_cuts + if master_incumbent > cut.residual + config.cut_violation_tolerance + ] + new_cut_count = sum(int(add_cut(cut)) for cut in violated) + trace.append( + { + "iteration": iteration, + "status": last_status, + "message": last_message, + "cut_count": len(cuts), + "master_incumbent": master_incumbent, + "master_upper_bound": master_upper, + "global_upper_bound": global_upper, + "candidate_lower_bound": candidate_lower, + "best_lower_bound": best_lower, + "absolute_gap": absolute_gap, + "relative_gap": relative_gap, + "violated_pair_count": len(violated), + "new_cut_count": new_cut_count, + } + ) + + tolerance = config.absolute_gap_tolerance + ( + config.relative_gap_tolerance * max(1.0, abs(best_lower)) + ) + if np.isfinite(global_upper) and absolute_gap <= tolerance: + certified = True + break + if not violated or new_cut_count == 0: + break + + absolute_gap = max(0.0, global_upper - best_lower) + relative_gap = absolute_gap / max(1.0, abs(best_lower)) + return MaximinCertificate( + allocation=np.asarray(best_allocation, dtype=float), + integer=config.integer, + certified=certified, + lower_bound=float(best_lower), + upper_bound=float(global_upper), + absolute_gap=float(absolute_gap), + relative_gap=float(relative_gap), + heuristic_lower_bound=float(heuristic_lower), + iterations=len(trace), + cut_count=len(cuts), + last_master_status=last_status, + last_master_message=last_message, + trace=tuple(trace), + ) + + +def certificate_matches_geometry( + signals: NDArray[np.float64], + certificate: MaximinCertificate, + *, + tolerance: float = 1.0e-9, +) -> bool: + """Return whether the frozen lower bound matches direct residual geometry.""" + + direct = pairwise_residual_matrix(signals, certificate.allocation) + off_diagonal = ~np.eye(direct.shape[0], dtype=bool) + return math.isclose( + certificate.lower_bound, + float(np.min(direct[off_diagonal])), + rel_tol=tolerance, + abs_tol=tolerance, + ) From e7634870889fbd9106bb9305060d67cfa689ee8d Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 10:11:47 +0800 Subject: [PATCH 02/25] Add certified design artifact CLI --- src/bayesian_ach/design_certificate_cli.py | 156 +++++++++++++++++++++ 1 file changed, 156 insertions(+) create mode 100644 src/bayesian_ach/design_certificate_cli.py diff --git a/src/bayesian_ach/design_certificate_cli.py b/src/bayesian_ach/design_certificate_cli.py new file mode 100644 index 0000000..9861977 --- /dev/null +++ b/src/bayesian_ach/design_certificate_cli.py @@ -0,0 +1,156 @@ +"""Command-line export for certified finite-grid maximin allocations.""" + +from __future__ import annotations + +import argparse +import csv +import hashlib +import json +from dataclasses import asdict +from pathlib import Path +from typing import Any, Sequence + +import numpy as np +import scipy + +from bayesian_ach.design_certificate import ( + CertifiedDesignConfig, + certificate_matches_geometry, + certify_maximin_design, +) +from bayesian_ach.design_grid import DESIGN_CANDIDATE_NAMES, generate_transition_design_grid + + +def _write_json(path: Path, value: dict[str, Any]) -> None: + path.write_text( + json.dumps(value, indent=2, sort_keys=True, allow_nan=False) + "\n", + encoding="utf-8", + ) + + +def _write_csv(path: Path, rows: list[dict[str, Any]]) -> None: + if not rows: + raise ValueError(f"cannot write empty table: {path}") + with path.open("w", newline="", encoding="utf-8") as handle: + writer = csv.DictWriter(handle, fieldnames=list(rows[0])) + writer.writeheader() + writer.writerows(rows) + + +def _manifest(output: Path, filenames: Sequence[str]) -> None: + rows = [] + for filename in sorted(filenames): + path = output / filename + rows.append( + { + "file": filename, + "bytes": path.stat().st_size, + "sha256": hashlib.sha256(path.read_bytes()).hexdigest(), + } + ) + _write_csv(output / "SHA256SUMS.csv", rows) + + +def _parser() -> argparse.ArgumentParser: + parser = argparse.ArgumentParser() + parser.add_argument("--output", type=Path, required=True) + parser.add_argument("--code-sha", required=True) + parser.add_argument("--budget", type=int, default=60) + parser.add_argument("--mode", choices=("integer", "continuous"), default="integer") + parser.add_argument("--max-point-fraction", type=float, default=0.15) + parser.add_argument("--absolute-gap-tolerance", type=float, default=1.0e-8) + parser.add_argument("--relative-gap-tolerance", type=float, default=1.0e-7) + parser.add_argument("--cut-violation-tolerance", type=float, default=1.0e-9) + parser.add_argument("--max-iterations", type=int, default=100) + parser.add_argument("--master-time-limit", type=float, default=120.0) + parser.add_argument("--master-mip-relative-gap", type=float, default=1.0e-9) + parser.add_argument("--require-certificate", action="store_true") + return parser + + +def main(argv: Sequence[str] | None = None) -> int: + """Run the certificate and write a checksum-locked result package.""" + + args = _parser().parse_args(argv) + config = CertifiedDesignConfig( + budget=args.budget, + max_point_fraction=args.max_point_fraction, + integer=args.mode == "integer", + absolute_gap_tolerance=args.absolute_gap_tolerance, + relative_gap_tolerance=args.relative_gap_tolerance, + cut_violation_tolerance=args.cut_violation_tolerance, + max_iterations=args.max_iterations, + master_time_limit_s=args.master_time_limit, + master_mip_relative_gap=args.master_mip_relative_gap, + ) + rows, _, standardized = generate_transition_design_grid() + result = certify_maximin_design(standardized, config) + output = args.output.resolve() + output.mkdir(parents=True, exist_ok=False) + + allocation_rows: list[dict[str, Any]] = [] + for point_index, mass in enumerate(result.allocation): + if mass <= 1.0e-12: + continue + allocation_rows.append( + { + **rows[point_index], + "allocation": int(round(float(mass))) + if result.integer + else float(mass), + } + ) + trace_rows = [dict(row) for row in result.trace] + summary = { + "schema_version": 1, + "experiment": "certified_finite_grid_maximin_allocation", + "code_repository": "IPS-Stuttgart/Bayesian-ACh", + "code_sha": args.code_sha, + "scipy_version": scipy.__version__, + "solver": "scipy.optimize.milp (HiGHS) with OLS cutting planes", + "candidate_names": list(DESIGN_CANDIDATE_NAMES), + "grid_point_count": len(rows), + "config": asdict(config), + "certified": result.certified, + "lower_bound": result.lower_bound, + "upper_bound": result.upper_bound, + "absolute_gap": result.absolute_gap, + "relative_gap": result.relative_gap, + "heuristic_lower_bound": result.heuristic_lower_bound, + "certified_improvement_over_heuristic": ( + result.lower_bound - result.heuristic_lower_bound + ), + "iterations": result.iterations, + "cut_count": result.cut_count, + "last_master_status": result.last_master_status, + "last_master_message": result.last_master_message, + "direct_geometry_matches_lower_bound": certificate_matches_geometry( + standardized, + result, + ), + "scope": ( + "This certificate is conditional on the independently instantiated finite " + "grid, global candidate standardization, budget, and per-cell cap. The " + "allocation is a count target, not a validated sequential protocol; no " + "history constructor, washout, reset, or carry-over feasibility is implied." + ), + } + _write_json(output / "certificate_summary.json", summary) + _write_csv(output / "certified_allocation.csv", allocation_rows) + _write_csv(output / "cut_trace.csv", trace_rows) + _manifest( + output, + ( + "certificate_summary.json", + "certified_allocation.csv", + "cut_trace.csv", + ), + ) + print(json.dumps(summary, indent=2, sort_keys=True)) + if args.require_certificate and not result.certified: + return 2 + return 0 + + +if __name__ == "__main__": + raise SystemExit(main()) From 5fa1ba12c8da5a963344068c2bf6b125da2abf08 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 10:12:20 +0800 Subject: [PATCH 03/25] Preserve formal master upper bounds --- src/bayesian_ach/design_certificate.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/bayesian_ach/design_certificate.py b/src/bayesian_ach/design_certificate.py index ec3c46c..7362b61 100644 --- a/src/bayesian_ach/design_certificate.py +++ b/src/bayesian_ach/design_certificate.py @@ -192,7 +192,7 @@ def certify_maximin_design( heuristic_lower = best_lower cuts: list[_PairCut] = [] - cut_keys: set[tuple[int, int, float, float]] = set() + cut_keys: set[tuple[int, int, str, str]] = set() def add_cut(cut: _PairCut) -> bool: key = ( From 6882d0f4ce126d41ae3701fb9bb29a75510354b6 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 10:12:55 +0800 Subject: [PATCH 04/25] Test certified maximin solver and artifacts --- tests/test_design_certificate.py | 189 +++++++++++++++++++++++++++++++ 1 file changed, 189 insertions(+) create mode 100644 tests/test_design_certificate.py diff --git a/tests/test_design_certificate.py b/tests/test_design_certificate.py new file mode 100644 index 0000000..b765e7e --- /dev/null +++ b/tests/test_design_certificate.py @@ -0,0 +1,189 @@ +from __future__ import annotations + +import csv +import itertools +import json +from pathlib import Path +from types import SimpleNamespace + +import numpy as np +import pytest + +from bayesian_ach.design_certificate import ( + CertifiedDesignConfig, + MaximinCertificate, + _master_upper_bound, + certificate_matches_geometry, + certify_maximin_design, +) +from bayesian_ach.design_certificate_cli import main +from bayesian_ach.design_geometry import pairwise_residual_matrix + + +def _objective(signals: np.ndarray, counts: np.ndarray) -> float: + residuals = pairwise_residual_matrix(signals, counts) + off_diagonal = ~np.eye(residuals.shape[0], dtype=bool) + return float(np.min(residuals[off_diagonal])) + + +def test_integer_cutting_plane_matches_brute_force() -> None: + signals = np.asarray( + [ + [-1.2, -0.7], + [-0.4, 0.9], + [0.1, -0.2], + [0.8, 0.3], + [1.3, 1.1], + ], + dtype=float, + ) + config = CertifiedDesignConfig( + budget=4, + max_point_fraction=0.5, + integer=True, + absolute_gap_tolerance=1.0e-9, + relative_gap_tolerance=1.0e-9, + max_iterations=100, + master_time_limit_s=30.0, + master_mip_relative_gap=1.0e-10, + ) + certificate = certify_maximin_design(signals, config) + + feasible = ( + np.asarray(values, dtype=np.int64) + for values in itertools.product(range(3), repeat=signals.shape[0]) + if sum(values) == config.budget + ) + brute_value = max(_objective(signals, counts) for counts in feasible) + + assert certificate.certified + assert certificate.lower_bound == pytest.approx(brute_value, abs=1.0e-8) + assert certificate.upper_bound - certificate.lower_bound <= 1.0e-8 + assert np.allclose(certificate.allocation, np.rint(certificate.allocation)) + assert int(np.sum(certificate.allocation)) == config.budget + assert certificate_matches_geometry(signals, certificate) + + +def test_continuous_certificate_bounds_integer_optimum() -> None: + signals = np.asarray( + [ + [-1.0, -0.2], + [-0.5, 1.0], + [0.0, -0.8], + [0.5, 0.7], + [1.0, 0.1], + ], + dtype=float, + ) + integer = certify_maximin_design( + signals, + CertifiedDesignConfig( + budget=4, + max_point_fraction=0.5, + integer=True, + max_iterations=100, + master_time_limit_s=30.0, + ), + ) + continuous = certify_maximin_design( + signals, + CertifiedDesignConfig( + budget=4, + max_point_fraction=0.5, + integer=False, + max_iterations=100, + master_time_limit_s=30.0, + ), + ) + + assert integer.certified + assert continuous.certified + assert continuous.lower_bound + 1.0e-8 >= integer.lower_bound + assert continuous.upper_bound - continuous.lower_bound <= 2.0e-7 + + +def test_mip_dual_bound_is_used_even_for_success_status() -> None: + result = SimpleNamespace(success=True, mip_dual_bound=-0.31) + assert _master_upper_bound(result, integer=True, incumbent=0.30) == pytest.approx( + 0.31 + ) + + +def test_invalid_initial_integer_allocation_is_rejected() -> None: + signals = np.asarray([[-1.0, 0.0], [0.0, 1.0], [1.0, -1.0]]) + with pytest.raises(ValueError, match="integer initial_allocation"): + certify_maximin_design( + signals, + CertifiedDesignConfig(budget=3, max_point_fraction=1.0), + initial_allocation=np.asarray([0.5, 1.0, 1.5]), + ) + + +def test_certificate_cli_writes_hash_bound_artifacts( + tmp_path: Path, + monkeypatch: pytest.MonkeyPatch, +) -> None: + import bayesian_ach.design_certificate_cli as cli + + rows = ( + {"point_id": 0, "condition": "a"}, + {"point_id": 1, "condition": "b"}, + ) + signals = np.asarray([[-1.0, 0.5], [1.0, -0.5]]) + certificate = MaximinCertificate( + allocation=np.asarray([1.0, 1.0]), + integer=True, + certified=True, + lower_bound=0.2, + upper_bound=0.2, + absolute_gap=0.0, + relative_gap=0.0, + heuristic_lower_bound=0.1, + iterations=1, + cut_count=2, + last_master_status=0, + last_master_message="optimal", + trace=( + { + "iteration": 1, + "status": 0, + "message": "optimal", + "cut_count": 2, + }, + ), + ) + monkeypatch.setattr( + cli, + "generate_transition_design_grid", + lambda: (rows, signals, signals), + ) + monkeypatch.setattr(cli, "certify_maximin_design", lambda *_args, **_kwargs: certificate) + monkeypatch.setattr(cli, "certificate_matches_geometry", lambda *_args: True) + + output = tmp_path / "certificate" + assert ( + main( + [ + "--output", + str(output), + "--code-sha", + "deadbeef", + "--budget", + "2", + "--max-point-fraction", + "1.0", + "--require-certificate", + ] + ) + == 0 + ) + summary = json.loads((output / "certificate_summary.json").read_text()) + assert summary["certified"] + assert summary["code_sha"] == "deadbeef" + with (output / "SHA256SUMS.csv").open(newline="", encoding="utf-8") as handle: + manifest = list(csv.DictReader(handle)) + assert {row["file"] for row in manifest} == { + "certificate_summary.json", + "certified_allocation.csv", + "cut_trace.csv", + } From ecfa33e82749f6165fabb1b09a8668f024c59c07 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 10:13:23 +0800 Subject: [PATCH 05/25] Expose certified design API --- src/bayesian_ach/design.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/bayesian_ach/design.py b/src/bayesian_ach/design.py index 95af0ad..54d4d92 100644 --- a/src/bayesian_ach/design.py +++ b/src/bayesian_ach/design.py @@ -1,6 +1,6 @@ """Public prospective experimental-design API.""" -from bayesian_ach.design_geometry import ( +from bayesian_ach.design_certificate import (\n CertifiedDesignConfig,\n MaximinCertificate,\n certificate_matches_geometry,\n certify_maximin_design,\n)\nfrom bayesian_ach.design_geometry import ( DesignDiagnostics, design_diagnostics, pairwise_residual_matrix, From 7c8b2a6c3fde2ea38cbee4f5b87b8b0f681e3530 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 10:13:38 +0800 Subject: [PATCH 06/25] Register certified design CLI From b9d0bc46157c320bdef318232318005943251c10 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 10:14:22 +0800 Subject: [PATCH 07/25] Fix certified design API exports --- src/bayesian_ach/design.py | 12 +++++++++++- 1 file changed, 11 insertions(+), 1 deletion(-) diff --git a/src/bayesian_ach/design.py b/src/bayesian_ach/design.py index 54d4d92..bb1c1e0 100644 --- a/src/bayesian_ach/design.py +++ b/src/bayesian_ach/design.py @@ -1,6 +1,12 @@ """Public prospective experimental-design API.""" -from bayesian_ach.design_certificate import (\n CertifiedDesignConfig,\n MaximinCertificate,\n certificate_matches_geometry,\n certify_maximin_design,\n)\nfrom bayesian_ach.design_geometry import ( +from bayesian_ach.design_certificate import ( + CertifiedDesignConfig, + MaximinCertificate, + certificate_matches_geometry, + certify_maximin_design, +) +from bayesian_ach.design_geometry import ( DesignDiagnostics, design_diagnostics, pairwise_residual_matrix, @@ -16,10 +22,14 @@ from bayesian_ach.design_optimizer import OptimizedDesign, optimize_maximin_design __all__ = [ + "CertifiedDesignConfig", "DESIGN_CANDIDATE_NAMES", "DesignDiagnostics", + "MaximinCertificate", "OptimizedDesign", "TransitionDesignGridConfig", + "certificate_matches_geometry", + "certify_maximin_design", "coupled_novelty_design", "design_diagnostics", "generate_transition_design_grid", From cb4c5c375f31b0b43b87a76e22fa01d2e75b8ef0 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 10:14:34 +0800 Subject: [PATCH 08/25] Register certified design CLI --- pyproject.toml | 1 + 1 file changed, 1 insertion(+) diff --git a/pyproject.toml b/pyproject.toml index a0c6107..477ef4d 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -70,6 +70,7 @@ plot = ["matplotlib>=3.8"] [project.scripts] bayesian-ach = "bayesian_ach.cli_ext:main" bayesian-ach-design = "bayesian_ach.design_cli:main" +bayesian-ach-design-certify = "bayesian_ach.design_certificate_cli:main" bayesian-ach-replay = "bayesian_ach.replay_cli:main" [tool.setuptools.packages.find] From 61d9863c740239fd617eeecdd86544148dcc6a08 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 10:14:55 +0800 Subject: [PATCH 09/25] Use rigorous master bounds and exact cut keys --- src/bayesian_ach/design_certificate.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/bayesian_ach/design_certificate.py b/src/bayesian_ach/design_certificate.py index 7362b61..982a811 100644 --- a/src/bayesian_ach/design_certificate.py +++ b/src/bayesian_ach/design_certificate.py @@ -142,12 +142,12 @@ def _validate_initial_allocation( def _master_upper_bound(result: Any, *, integer: bool, incumbent: float) -> float: - if bool(result.success): - return incumbent if integer: dual = getattr(result, "mip_dual_bound", None) if dual is not None and np.isfinite(float(dual)): - return float(-float(dual)) + return max(incumbent, float(-float(dual))) + if bool(result.success): + return incumbent return math.inf @@ -198,8 +198,8 @@ def add_cut(cut: _PairCut) -> bool: key = ( cut.generator, cut.alternative, - round(cut.intercept, 12), - round(cut.slope, 12), + cut.intercept.hex(), + cut.slope.hex(), ) if key in cut_keys: return False From 017c0e8956ea4eba5747abbb7eb55b4c550d101c Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 10:15:23 +0800 Subject: [PATCH 10/25] Document certified design scope and bounds --- docs/optimal_design.md | 42 ++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 42 insertions(+) diff --git a/docs/optimal_design.md b/docs/optimal_design.md index b243b5b..6b3fb57 100644 --- a/docs/optimal_design.md +++ b/docs/optimal_design.md @@ -202,6 +202,48 @@ closed-loop timing can be incorporated by replacing the candidate matrix with the corresponding forward-model predictions while retaining the same auditable finite-allocation structure. +## Certified continuous and integer follow-up + +The original `maximin_optimized` allocation is a deterministic greedy +construction with at most three one-for-one exchanges. Its trace is auditable, +but it is not a certificate of global optimality. + +`certify_maximin_design` supplies a separate global certificate. For each +ordered pair, + +```math +R_{k\mid l}(w) += +\inf_{b,c}\sum_d w_d[x_k(d)-b-cx_l(d)]^2. +``` + +For any fixed `(b,c)`, the loss is linear in the allocation. The certificate +iteratively solves a HiGHS LP or MILP master problem, evaluates exact weighted +least squares for all ordered pairs, and adds every violated loss cut. The +finite master gives an upper bound; the directly evaluated allocation gives a +feasible lower bound. Results report both bounds and stop as certified only +when their declared tolerance is met. + +The integer mode certifies the stated budget and per-cell count cap. The +continuous mode certifies the capped-simplex relaxation and therefore also +provides an upper bound for every integer allocation. Neither mode constructs a +sequential behavioral history. The 240 cells are independently instantiated +belief conditions; reset, washout, carry-over, and physical ordering constraints +must be encoded before calling any count vector an executable protocol. + +```bash +bayesian-ach-design-certify \ + --output results/certified-design-n60 \ + --code-sha \ + --budget 60 \ + --mode integer \ + --require-certificate +``` + +Certificate artifacts are versioned separately from the previously frozen +greedy/exchange evidence. A changed certified allocation must trigger new +recovery evidence; it must never silently replace the earlier artifact. + ## Use ```bash From 5a5acf2db26a2200a83a4c03bc10b0d87e36f687 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 10:17:21 +0800 Subject: [PATCH 11/25] Validate certificate with shared OLS oracle --- src/bayesian_ach/design_certificate.py | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/src/bayesian_ach/design_certificate.py b/src/bayesian_ach/design_certificate.py index 982a811..66101d1 100644 --- a/src/bayesian_ach/design_certificate.py +++ b/src/bayesian_ach/design_certificate.py @@ -10,7 +10,6 @@ from numpy.typing import NDArray from scipy.optimize import Bounds, LinearConstraint, milp -from bayesian_ach.design_geometry import pairwise_residual_matrix from bayesian_ach.design_optimizer import optimize_maximin_design @@ -346,11 +345,13 @@ def certificate_matches_geometry( ) -> bool: """Return whether the frozen lower bound matches direct residual geometry.""" - direct = pairwise_residual_matrix(signals, certificate.allocation) - off_diagonal = ~np.eye(direct.shape[0], dtype=bool) + direct, _ = _pairwise_oracle( + np.asarray(signals, dtype=float), + certificate.allocation, + ) return math.isclose( certificate.lower_bound, - float(np.min(direct[off_diagonal])), + direct, rel_tol=tolerance, abs_tol=tolerance, ) From e6e1f2270ba699f39651f866aea62f6562ac2792 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 10:19:57 +0800 Subject: [PATCH 12/25] Fix certified CLI lint imports --- src/bayesian_ach/design_certificate_cli.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/bayesian_ach/design_certificate_cli.py b/src/bayesian_ach/design_certificate_cli.py index 9861977..6a9f7c6 100644 --- a/src/bayesian_ach/design_certificate_cli.py +++ b/src/bayesian_ach/design_certificate_cli.py @@ -3,14 +3,14 @@ from __future__ import annotations import argparse +from collections.abc import Sequence import csv import hashlib import json from dataclasses import asdict from pathlib import Path -from typing import Any, Sequence +from typing import Any -import numpy as np import scipy from bayesian_ach.design_certificate import ( From 521a752460f680435802d9843478b97bf27dc3b4 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 10:21:51 +0800 Subject: [PATCH 13/25] Sort certified CLI imports --- src/bayesian_ach/design_certificate_cli.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/bayesian_ach/design_certificate_cli.py b/src/bayesian_ach/design_certificate_cli.py index 6a9f7c6..6d1fc41 100644 --- a/src/bayesian_ach/design_certificate_cli.py +++ b/src/bayesian_ach/design_certificate_cli.py @@ -3,10 +3,10 @@ from __future__ import annotations import argparse -from collections.abc import Sequence import csv import hashlib import json +from collections.abc import Sequence from dataclasses import asdict from pathlib import Path from typing import Any From 362305cbcbf48754130116a2389b79a5f0075179 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 11:11:53 +0800 Subject: [PATCH 14/25] Derive rigorous N_eff brackets from certificate bounds --- src/bayesian_ach/design_certificate.py | 62 ++++++++++++++++++++++++++ 1 file changed, 62 insertions(+) diff --git a/src/bayesian_ach/design_certificate.py b/src/bayesian_ach/design_certificate.py index 66101d1..9b797be 100644 --- a/src/bayesian_ach/design_certificate.py +++ b/src/bayesian_ach/design_certificate.py @@ -10,6 +10,7 @@ from numpy.typing import NDArray from scipy.optimize import Bounds, LinearConstraint, milp +from bayesian_ach.design_geometry import profiled_gaussian_log_score_gap from bayesian_ach.design_optimizer import optimize_maximin_design @@ -337,6 +338,67 @@ def add_cut(cut: _PairCut) -> bool: ) + +def population_n_eff_index( + residual_variance: float, + *, + effect_size: float = 1.0, + noise_std: float = 1.0, + target_log_score_gap: float = 5.0, +) -> int: + """Map positive residual geometry to a population observation-equivalent index.""" + + values = ( + residual_variance, + effect_size, + noise_std, + target_log_score_gap, + ) + if not all(np.isfinite(value) for value in values): + raise ValueError("N_eff inputs must be finite") + if residual_variance <= 0.0: + raise ValueError("residual_variance must be positive") + if effect_size <= 0.0 or noise_std <= 0.0 or target_log_score_gap <= 0.0: + raise ValueError("effect_size, noise_std, and target_log_score_gap must be positive") + per_observation = profiled_gaussian_log_score_gap( + residual_variance, + effect_size=effect_size, + noise_std=noise_std, + ) + return int(math.ceil(target_log_score_gap / per_observation)) + + +def population_n_eff_bracket( + residual_lower_bound: float, + residual_upper_bound: float, + *, + effect_size: float = 1.0, + noise_std: float = 1.0, + target_log_score_gap: float = 5.0, +) -> tuple[int, int]: + """Return rigorous lower/upper indices induced by residual objective bounds.""" + + if ( + not np.isfinite(residual_lower_bound) + or not np.isfinite(residual_upper_bound) + or residual_lower_bound <= 0.0 + or residual_upper_bound < residual_lower_bound + ): + raise ValueError("residual bounds must be finite, positive, and ordered") + optimistic = population_n_eff_index( + residual_upper_bound, + effect_size=effect_size, + noise_std=noise_std, + target_log_score_gap=target_log_score_gap, + ) + conservative = population_n_eff_index( + residual_lower_bound, + effect_size=effect_size, + noise_std=noise_std, + target_log_score_gap=target_log_score_gap, + ) + return optimistic, conservative + def certificate_matches_geometry( signals: NDArray[np.float64], certificate: MaximinCertificate, From 3e1fa12f42cefa0149d7b5713a7d13f670bd966e Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 11:12:10 +0800 Subject: [PATCH 15/25] Wrap N_eff validation diagnostics --- src/bayesian_ach/design_certificate.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/src/bayesian_ach/design_certificate.py b/src/bayesian_ach/design_certificate.py index 9b797be..67c9fed 100644 --- a/src/bayesian_ach/design_certificate.py +++ b/src/bayesian_ach/design_certificate.py @@ -359,7 +359,9 @@ def population_n_eff_index( if residual_variance <= 0.0: raise ValueError("residual_variance must be positive") if effect_size <= 0.0 or noise_std <= 0.0 or target_log_score_gap <= 0.0: - raise ValueError("effect_size, noise_std, and target_log_score_gap must be positive") + raise ValueError( + "effect_size, noise_std, and target_log_score_gap must be positive" + ) per_observation = profiled_gaussian_log_score_gap( residual_variance, effect_size=effect_size, From 019177f1af2af568f007837558aa30b068a015f9 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 11:12:34 +0800 Subject: [PATCH 16/25] Record the certified population N_eff bracket --- src/bayesian_ach/design_certificate_cli.py | 37 ++++++++++++++++++++++ 1 file changed, 37 insertions(+) diff --git a/src/bayesian_ach/design_certificate_cli.py b/src/bayesian_ach/design_certificate_cli.py index 6d1fc41..702d8f1 100644 --- a/src/bayesian_ach/design_certificate_cli.py +++ b/src/bayesian_ach/design_certificate_cli.py @@ -6,6 +6,7 @@ import csv import hashlib import json +import math from collections.abc import Sequence from dataclasses import asdict from pathlib import Path @@ -17,6 +18,7 @@ CertifiedDesignConfig, certificate_matches_geometry, certify_maximin_design, + population_n_eff_bracket, ) from bayesian_ach.design_grid import DESIGN_CANDIDATE_NAMES, generate_transition_design_grid @@ -64,6 +66,9 @@ def _parser() -> argparse.ArgumentParser: parser.add_argument("--max-iterations", type=int, default=100) parser.add_argument("--master-time-limit", type=float, default=120.0) parser.add_argument("--master-mip-relative-gap", type=float, default=1.0e-9) + parser.add_argument("--effect-size", type=float, default=1.0) + parser.add_argument("--noise-std", type=float, default=1.0) + parser.add_argument("--target-log-score-gap", type=float, default=5.0) parser.add_argument("--require-certificate", action="store_true") return parser @@ -83,8 +88,27 @@ def main(argv: Sequence[str] | None = None) -> int: master_time_limit_s=args.master_time_limit, master_mip_relative_gap=args.master_mip_relative_gap, ) + population_target = ( + float(args.effect_size), + float(args.noise_std), + float(args.target_log_score_gap), + ) + if ( + not all(math.isfinite(value) for value in population_target) + or any(value <= 0.0 for value in population_target) + ): + raise ValueError( + "effect-size, noise-std, and target-log-score-gap must be finite and positive" + ) rows, _, standardized = generate_transition_design_grid() result = certify_maximin_design(standardized, config) + n_eff_lower, n_eff_upper = population_n_eff_bracket( + result.lower_bound, + result.upper_bound, + effect_size=args.effect_size, + noise_std=args.noise_std, + target_log_score_gap=args.target_log_score_gap, + ) output = args.output.resolve() output.mkdir(parents=True, exist_ok=False) @@ -128,6 +152,19 @@ def main(argv: Sequence[str] | None = None) -> int: standardized, result, ), + "population_n_eff": { + "effect_size": args.effect_size, + "noise_std": args.noise_std, + "target_log_score_gap": args.target_log_score_gap, + "lower_index_from_residual_upper_bound": n_eff_lower, + "upper_index_from_residual_lower_bound": n_eff_upper, + "index_certified": n_eff_lower == n_eff_upper, + "interpretation": ( + "Population observation-equivalent index under effectively " + "independent Gaussian observations; not a physical trial, " + "time-bin, session, or animal count." + ), + }, "scope": ( "This certificate is conditional on the independently instantiated finite " "grid, global candidate standardization, budget, and per-cell cap. The " From ac2f42c19024ca2d5f11312f91e8f56c6cd5ab86 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 11:13:25 +0800 Subject: [PATCH 17/25] Add independent certificate package verification --- src/bayesian_ach/design_certificate_verify.py | 206 ++++++++++++++++++ 1 file changed, 206 insertions(+) create mode 100644 src/bayesian_ach/design_certificate_verify.py diff --git a/src/bayesian_ach/design_certificate_verify.py b/src/bayesian_ach/design_certificate_verify.py new file mode 100644 index 0000000..10cfe7a --- /dev/null +++ b/src/bayesian_ach/design_certificate_verify.py @@ -0,0 +1,206 @@ +"""Independent verification for frozen integer maximin certificate packages.""" + +from __future__ import annotations + +import argparse +import csv +import hashlib +import json +import math +import re +from collections.abc import Sequence +from pathlib import Path +from typing import Any + +import numpy as np + +from bayesian_ach.design_certificate import population_n_eff_bracket +from bayesian_ach.design_geometry import pairwise_residual_matrix +from bayesian_ach.design_grid import generate_transition_design_grid + +_EXPECTED_PAYLOADS = { + "certificate_summary.json", + "certified_allocation.csv", + "cut_trace.csv", +} +_SHA_PATTERN = re.compile(r"^[0-9a-f]{40}$") + + +def _sha256(path: Path) -> str: + digest = hashlib.sha256() + with path.open("rb") as handle: + for chunk in iter(lambda: handle.read(1024 * 1024), b""): + digest.update(chunk) + return digest.hexdigest() + + +def _locked_payloads(directory: Path) -> dict[str, dict[str, str]]: + checksum_path = directory / "SHA256SUMS.csv" + if not checksum_path.is_file(): + raise ValueError(f"missing checksum table: {checksum_path}") + with checksum_path.open(newline="", encoding="utf-8") as handle: + rows = list(csv.DictReader(handle)) + locked: dict[str, dict[str, str]] = {} + for row in rows: + name = str(row["file"]) + if name in locked: + raise ValueError(f"duplicate checksum row: {name}") + locked[name] = {key: str(value) for key, value in row.items()} + if set(locked) != _EXPECTED_PAYLOADS: + raise ValueError("checksum table does not bind the exact certificate payload set") + for name, row in locked.items(): + path = directory / name + if not path.is_file(): + raise ValueError(f"missing locked payload: {name}") + if int(row["bytes"]) != path.stat().st_size: + raise ValueError(f"locked byte count mismatch: {name}") + if row["sha256"] != _sha256(path): + raise ValueError(f"locked SHA-256 mismatch: {name}") + return locked + + +def _allocation_counts( + path: Path, + *, + point_count: int, +) -> np.ndarray: + counts = np.zeros(point_count, dtype=np.int64) + with path.open(newline="", encoding="utf-8") as handle: + rows = list(csv.DictReader(handle)) + if not rows: + raise ValueError("certified allocation is empty") + seen: set[int] = set() + for row in rows: + point = int(row["point_id"]) + value = float(row["allocation"]) + rounded = int(round(value)) + if point < 0 or point >= point_count or point in seen: + raise ValueError("certified allocation has an invalid or duplicate point") + if rounded <= 0 or not math.isclose(value, rounded, abs_tol=1.0e-9): + raise ValueError("certified integer allocation contains a noninteger count") + counts[point] = rounded + seen.add(point) + return counts + + +def verify_certificate_package(directory: Path) -> dict[str, Any]: + """Recompute hashes, finite-grid geometry, gap, and N_eff rounding.""" + + directory = directory.resolve() + _locked_payloads(directory) + summary = json.loads( + (directory / "certificate_summary.json").read_text(encoding="utf-8") + ) + if ( + summary.get("schema_version") != 1 + or summary.get("experiment") != "certified_finite_grid_maximin_allocation" + or summary.get("certified") is not True + or summary.get("direct_geometry_matches_lower_bound") is not True + ): + raise ValueError("package is not a certified finite-grid allocation") + if _SHA_PATTERN.fullmatch(str(summary.get("code_sha", ""))) is None: + raise ValueError("certificate code_sha is not a full lowercase commit SHA") + + config = summary["config"] + if config.get("integer") is not True: + raise ValueError("claim-bearing package must certify an integer allocation") + budget = int(config["budget"]) + maximum_count = max( + 1, + int(math.ceil(float(config["max_point_fraction"]) * budget)), + ) + rows, _, signals = generate_transition_design_grid() + if int(summary["grid_point_count"]) != len(rows): + raise ValueError("certificate grid-point count does not match the canonical grid") + counts = _allocation_counts( + directory / "certified_allocation.csv", + point_count=len(rows), + ) + if int(np.sum(counts)) != budget or np.any(counts > maximum_count): + raise ValueError("certified allocation violates its budget or per-cell cap") + + geometry = pairwise_residual_matrix(signals, counts) + off_diagonal = ~np.eye(geometry.shape[0], dtype=bool) + direct_lower = float(np.min(geometry[off_diagonal])) + lower = float(summary["lower_bound"]) + upper = float(summary["upper_bound"]) + absolute_gap = float(summary["absolute_gap"]) + relative_gap = float(summary["relative_gap"]) + if ( + not np.isfinite([lower, upper, absolute_gap, relative_gap]).all() + or lower <= 0.0 + or upper < lower + ): + raise ValueError("certificate bounds are not finite, positive, and ordered") + if not math.isclose(direct_lower, lower, rel_tol=1.0e-9, abs_tol=1.0e-9): + raise ValueError("direct finite-grid geometry does not reproduce the lower bound") + recomputed_gap = max(0.0, upper - lower) + if not math.isclose(absolute_gap, recomputed_gap, rel_tol=1.0e-9, abs_tol=1.0e-12): + raise ValueError("reported absolute certificate gap is inconsistent") + expected_relative = recomputed_gap / max(1.0, abs(lower)) + if not math.isclose(relative_gap, expected_relative, rel_tol=1.0e-9, abs_tol=1.0e-12): + raise ValueError("reported relative certificate gap is inconsistent") + tolerance = float(config["absolute_gap_tolerance"]) + float( + config["relative_gap_tolerance"] + ) * max(1.0, abs(lower)) + if absolute_gap > tolerance: + raise ValueError("reported certified gap exceeds the frozen tolerance") + + target = summary["population_n_eff"] + n_eff = population_n_eff_bracket( + lower, + upper, + effect_size=float(target["effect_size"]), + noise_std=float(target["noise_std"]), + target_log_score_gap=float(target["target_log_score_gap"]), + ) + reported_n_eff = ( + int(target["lower_index_from_residual_upper_bound"]), + int(target["upper_index_from_residual_lower_bound"]), + ) + if n_eff != reported_n_eff or bool(target["index_certified"]) != (n_eff[0] == n_eff[1]): + raise ValueError("reported population N_eff bracket is inconsistent") + + with (directory / "cut_trace.csv").open(newline="", encoding="utf-8") as handle: + trace = list(csv.DictReader(handle)) + if len(trace) != int(summary["iterations"]): + raise ValueError("cut trace length does not match the reported iteration count") + if int(trace[-1]["cut_count"]) != int(summary["cut_count"]): + raise ValueError("cut trace does not end at the reported cut count") + + return { + "package": directory.name, + "verified": True, + "producer_commit": summary["code_sha"], + "budget": budget, + "support_size": int(np.count_nonzero(counts)), + "maximum_cell_count": int(np.max(counts)), + "lower_bound": lower, + "upper_bound": upper, + "absolute_gap": absolute_gap, + "population_n_eff_lower": n_eff[0], + "population_n_eff_upper": n_eff[1], + "sha256sums_sha256": _sha256(directory / "SHA256SUMS.csv"), + } + + +def _parser() -> argparse.ArgumentParser: + parser = argparse.ArgumentParser( + prog="bayesian-ach-design-certificate-verify", + description="Verify frozen integer maximin certificate packages.", + ) + parser.add_argument("packages", nargs="+", type=Path) + return parser + + +def main(argv: Sequence[str] | None = None) -> int: + reports = [ + verify_certificate_package(path) + for path in _parser().parse_args(argv).packages + ] + print(json.dumps(reports, indent=2, sort_keys=True, allow_nan=False)) + return 0 + + +if __name__ == "__main__": # pragma: no cover + raise SystemExit(main()) From 5c76df81b1050e4826f9263ed58afb926e945360 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 11:14:05 +0800 Subject: [PATCH 18/25] Expose certificate package verification --- pyproject.toml | 1 + 1 file changed, 1 insertion(+) diff --git a/pyproject.toml b/pyproject.toml index 477ef4d..acf6f9c 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -70,6 +70,7 @@ plot = ["matplotlib>=3.8"] [project.scripts] bayesian-ach = "bayesian_ach.cli_ext:main" bayesian-ach-design = "bayesian_ach.design_cli:main" +bayesian-ach-design-certificate-verify = "bayesian_ach.design_certificate_verify:main" bayesian-ach-design-certify = "bayesian_ach.design_certificate_cli:main" bayesian-ach-replay = "bayesian_ach.replay_cli:main" From a72861b2878b5092b3e27b462e13e32ed47442cb Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 11:14:37 +0800 Subject: [PATCH 19/25] Test N_eff rounding and certificate package verification --- tests/test_design_certificate.py | 39 ++++++++++++++++++++++++++++++-- 1 file changed, 37 insertions(+), 2 deletions(-) diff --git a/tests/test_design_certificate.py b/tests/test_design_certificate.py index b765e7e..c0b507b 100644 --- a/tests/test_design_certificate.py +++ b/tests/test_design_certificate.py @@ -15,6 +15,8 @@ _master_upper_bound, certificate_matches_geometry, certify_maximin_design, + population_n_eff_bracket, + population_n_eff_index, ) from bayesian_ach.design_certificate_cli import main from bayesian_ach.design_geometry import pairwise_residual_matrix @@ -102,6 +104,15 @@ def test_continuous_certificate_bounds_integer_optimum() -> None: assert continuous.upper_bound - continuous.lower_bound <= 2.0e-7 +def test_population_n_eff_bracket_certifies_rounding_without_objective_tolerance() -> None: + lower = 0.253106863558372 + upper = 0.2532561178966192 + assert population_n_eff_index(lower) == 45 + assert population_n_eff_bracket(lower, upper) == (45, 45) + with pytest.raises(ValueError, match="ordered"): + population_n_eff_bracket(upper, lower) + + def test_mip_dual_bound_is_used_even_for_success_status() -> None: result = SimpleNamespace(success=True, mip_dual_bound=-0.31) assert _master_upper_bound(result, integer=True, incumbent=0.30) == pytest.approx( @@ -167,7 +178,7 @@ def test_certificate_cli_writes_hash_bound_artifacts( "--output", str(output), "--code-sha", - "deadbeef", + "d" * 40, "--budget", "2", "--max-point-fraction", @@ -179,7 +190,8 @@ def test_certificate_cli_writes_hash_bound_artifacts( ) summary = json.loads((output / "certificate_summary.json").read_text()) assert summary["certified"] - assert summary["code_sha"] == "deadbeef" + assert summary["code_sha"] == "d" * 40 + assert summary["population_n_eff"]["index_certified"] is True with (output / "SHA256SUMS.csv").open(newline="", encoding="utf-8") as handle: manifest = list(csv.DictReader(handle)) assert {row["file"] for row in manifest} == { @@ -187,3 +199,26 @@ def test_certificate_cli_writes_hash_bound_artifacts( "certified_allocation.csv", "cut_trace.csv", } + + import bayesian_ach.design_certificate_verify as verifier + + monkeypatch.setattr( + verifier, + "generate_transition_design_grid", + lambda: (rows, signals, signals), + ) + monkeypatch.setattr( + verifier, + "pairwise_residual_matrix", + lambda *_args: np.asarray([[0.0, 0.2], [0.2, 0.0]]), + ) + report = verifier.verify_certificate_package(output) + assert report["verified"] is True + assert report["population_n_eff_lower"] == report["population_n_eff_upper"] + + (output / "certified_allocation.csv").write_text( + "point_id,condition,allocation\n0,a,2\n", + encoding="utf-8", + ) + with pytest.raises(ValueError, match="SHA-256 mismatch"): + verifier.verify_certificate_package(output) From dba80ced740659d0721d70276f112f9763de31ad Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 11:14:56 +0800 Subject: [PATCH 20/25] Normalize certificate verifier test imports --- tests/test_design_certificate.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/tests/test_design_certificate.py b/tests/test_design_certificate.py index c0b507b..d7b6da6 100644 --- a/tests/test_design_certificate.py +++ b/tests/test_design_certificate.py @@ -9,6 +9,7 @@ import numpy as np import pytest +import bayesian_ach.design_certificate_verify as verifier from bayesian_ach.design_certificate import ( CertifiedDesignConfig, MaximinCertificate, @@ -200,8 +201,6 @@ def test_certificate_cli_writes_hash_bound_artifacts( "cut_trace.csv", } - import bayesian_ach.design_certificate_verify as verifier - monkeypatch.setattr( verifier, "generate_transition_design_grid", From 4c952e6c599534f8717362c71c7eda6efca95e88 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 11:15:36 +0800 Subject: [PATCH 21/25] Document objective and N_eff certificate boundaries --- docs/optimal_design.md | 25 +++++++++++++++++++++++++ 1 file changed, 25 insertions(+) diff --git a/docs/optimal_design.md b/docs/optimal_design.md index 6b3fb57..1b347eb 100644 --- a/docs/optimal_design.md +++ b/docs/optimal_design.md @@ -244,6 +244,31 @@ Certificate artifacts are versioned separately from the previously frozen greedy/exchange evidence. A changed certified allocation must trigger new recovery evidence; it must never silently replace the earlier artifact. +Objective certification and planning-index certification are reported +separately. If rigorous residual bounds are `[R_L,R_U]`, monotonicity of the +profiled Gaussian gap gives + +```math +\left\lceil\frac{B}{G(R_U)}\right\rceil +\le N_{\mathrm{eff}} \le +\left\lceil\frac{B}{G(R_L)}\right\rceil. +``` + +The rounded `N_eff` index is certified when these endpoints agree, even if a +continuous objective run has not met a much tighter numerical residual-gap +tolerance. This does not relabel that continuous objective as certified. The +claim-bearing finite schedules use exact integer certificates. + +Frozen integer packages can be independently checked against their SHA-256 +table, allocation budget and cap, canonical 240-cell geometry, objective gap, +cut trace, and `N_eff` rounding: + +```bash +bayesian-ach-design-certificate-verify \ + results/certified-maximin-design/n60 \ + results/certified-maximin-design/n45 +``` + ## Use ```bash From e82f48d335fb3f17f843d16bf2442693b6c83eef Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 11:17:54 +0800 Subject: [PATCH 22/25] Accept either tamper-detection boundary in verifier test --- tests/test_design_certificate.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_design_certificate.py b/tests/test_design_certificate.py index d7b6da6..1fa158d 100644 --- a/tests/test_design_certificate.py +++ b/tests/test_design_certificate.py @@ -219,5 +219,5 @@ def test_certificate_cli_writes_hash_bound_artifacts( "point_id,condition,allocation\n0,a,2\n", encoding="utf-8", ) - with pytest.raises(ValueError, match="SHA-256 mismatch"): + with pytest.raises(ValueError, match="byte count mismatch|SHA-256 mismatch"): verifier.verify_certificate_package(output) From 80565d7ac6d3ebea509fcbea133929faf65185fe Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 05:21:58 +0200 Subject: [PATCH 23/25] Freeze exact integer maximin certificates Add independently verifiable N=60 and N=45 finite-grid integer certificate packages generated at e82f48d335fb3f17f843d16bf2442693b6c83eef. Both retain the population N_eff bracket [45, 45] and remain count targets rather than executable protocols. --- .../n45/SHA256SUMS.csv | 4 ++ .../n45/certificate_summary.json | 50 +++++++++++++++++++ .../n45/certified_allocation.csv | 16 ++++++ .../n45/cut_trace.csv | 8 +++ .../n60/SHA256SUMS.csv | 4 ++ .../n60/certificate_summary.json | 50 +++++++++++++++++++ .../n60/certified_allocation.csv | 13 +++++ .../n60/cut_trace.csv | 6 +++ 8 files changed, 151 insertions(+) create mode 100644 results/certified-maximin-design/n45/SHA256SUMS.csv create mode 100644 results/certified-maximin-design/n45/certificate_summary.json create mode 100644 results/certified-maximin-design/n45/certified_allocation.csv create mode 100644 results/certified-maximin-design/n45/cut_trace.csv create mode 100644 results/certified-maximin-design/n60/SHA256SUMS.csv create mode 100644 results/certified-maximin-design/n60/certificate_summary.json create mode 100644 results/certified-maximin-design/n60/certified_allocation.csv create mode 100644 results/certified-maximin-design/n60/cut_trace.csv diff --git a/results/certified-maximin-design/n45/SHA256SUMS.csv b/results/certified-maximin-design/n45/SHA256SUMS.csv new file mode 100644 index 0000000..14902cc --- /dev/null +++ b/results/certified-maximin-design/n45/SHA256SUMS.csv @@ -0,0 +1,4 @@ +file,bytes,sha256 +certificate_summary.json,1970,b527c4373533650c9259e92846ad7680e7cff7d7b7cb169885aeba3663dae2f7 +certified_allocation.csv,2382,bbd0796cc7c767c3388cd6c3b9522cebaed78993586aa5dfd3b353e856abe856 +cut_trace.csv,1702,4d365c5f8635fbc17c2d39390326ceeace7ddc1db4a2eb9ca3529252859289ba diff --git a/results/certified-maximin-design/n45/certificate_summary.json b/results/certified-maximin-design/n45/certificate_summary.json new file mode 100644 index 0000000..7c05255 --- /dev/null +++ b/results/certified-maximin-design/n45/certificate_summary.json @@ -0,0 +1,50 @@ +{ + "absolute_gap": 2.7755575615628914e-16, + "candidate_names": [ + "innovation_l2", + "surprise", + "gain", + "update_l2", + "information_gain", + "change_probability" + ], + "certified": true, + "certified_improvement_over_heuristic": 0.0025656042479573204, + "code_repository": "IPS-Stuttgart/Bayesian-ACh", + "code_sha": "e82f48d335fb3f17f843d16bf2442693b6c83eef", + "config": { + "absolute_gap_tolerance": 1e-08, + "budget": 45, + "cut_violation_tolerance": 1e-09, + "integer": true, + "master_mip_relative_gap": 1e-09, + "master_time_limit_s": 120.0, + "max_iterations": 100, + "max_point_fraction": 0.15, + "relative_gap_tolerance": 1e-07 + }, + "cut_count": 53, + "direct_geometry_matches_lower_bound": true, + "experiment": "certified_finite_grid_maximin_allocation", + "grid_point_count": 240, + "heuristic_lower_bound": 0.24983023368504692, + "iterations": 7, + "last_master_message": "Optimization terminated successfully. (HiGHS Status 7: Optimal)", + "last_master_status": 0, + "lower_bound": 0.25239583793300424, + "population_n_eff": { + "effect_size": 1.0, + "index_certified": true, + "interpretation": "Population observation-equivalent index under effectively independent Gaussian observations; not a physical trial, time-bin, session, or animal count.", + "lower_index_from_residual_upper_bound": 45, + "noise_std": 1.0, + "target_log_score_gap": 5.0, + "upper_index_from_residual_lower_bound": 45 + }, + "relative_gap": 2.7755575615628914e-16, + "schema_version": 1, + "scipy_version": "1.15.3", + "scope": "This certificate is conditional on the independently instantiated finite grid, global candidate standardization, budget, and per-cell cap. The allocation is a count target, not a validated sequential protocol; no history constructor, washout, reset, or carry-over feasibility is implied.", + "solver": "scipy.optimize.milp (HiGHS) with OLS cutting planes", + "upper_bound": 0.2523958379330045 +} diff --git a/results/certified-maximin-design/n45/certified_allocation.csv b/results/certified-maximin-design/n45/certified_allocation.csv new file mode 100644 index 0000000..85e75ed --- /dev/null +++ b/results/certified-maximin-design/n45/certified_allocation.csv @@ -0,0 +1,16 @@ +point_id,observed_probability,residual_shape,concentration,reset_observed_probability,hazard,innovation_l2,surprise,gain,update_l2,information_gain,change_probability,allocation +5,0.05,0.5,2.0,0.95,0.15,1.1635076278220096,2.995732273553991,0.3333333333333333,0.3878358759406699,1.6491929980444475,0.7702702702702702,1 +10,0.05,0.5,8.0,0.95,0.01,1.1635076278220096,2.995732273553991,0.1111111111111111,0.1292786253135566,0.793706251013266,0.16101694915254236,4 +13,0.05,0.5,32.0,0.05,0.15,1.1635076278220096,2.995732273553991,0.030303030303030304,0.03525780690369727,0.2655001957924723,0.15,1 +16,0.05,0.5,32.0,0.95,0.01,1.1635076278220096,2.995732273553991,0.030303030303030304,0.03525780690369727,0.2655001957924723,0.16101694915254236,4 +120,0.35,0.9,2.0,0.05,0.01,0.8768979416100827,1.0498221244986778,0.3333333333333333,0.29229931387002756,0.3355856642737046,0.0014409221902017292,1 +123,0.35,0.9,2.0,0.5,0.15,0.8768979416100827,1.0498221244986778,0.3333333333333333,0.29229931387002756,0.3355856642737046,0.20134228187919462,3 +125,0.35,0.9,2.0,0.95,0.15,0.8768979416100827,1.0498221244986778,0.3333333333333333,0.29229931387002756,0.3355856642737046,0.32386363636363635,2 +128,0.35,0.9,8.0,0.5,0.01,0.8768979416100827,1.0498221244986778,0.1111111111111111,0.09743310462334258,0.10687049253534275,0.014224751066856332,1 +129,0.35,0.9,8.0,0.5,0.15,0.8768979416100827,1.0498221244986778,0.1111111111111111,0.09743310462334258,0.10687049253534275,0.20134228187919462,7 +130,0.35,0.9,8.0,0.95,0.01,0.8768979416100827,1.0498221244986778,0.1111111111111111,0.09743310462334258,0.10687049253534275,0.026685393258426966,5 +143,0.35,0.9,128.0,0.95,0.15,0.8768979416100827,1.0498221244986778,0.007751937984496124,0.0067976584620936574,0.007218032071614644,0.32386363636363635,2 +197,0.9,0.5,2.0,0.95,0.15,0.12247448713915887,0.10536051565782628,0.3333333333333333,0.0408248290463862,0.023123169408776692,0.15702479338842976,7 +201,0.9,0.5,8.0,0.5,0.15,0.12247448713915887,0.10536051565782628,0.1111111111111111,0.013608276348795313,0.006640070861282682,0.08928571428571429,4 +207,0.9,0.5,32.0,0.5,0.15,0.12247448713915887,0.10536051565782628,0.030303030303030304,0.0037113480951260457,0.0017170260888290656,0.08928571428571429,2 +211,0.9,0.5,128.0,0.05,0.15,0.12247448713915887,0.10536051565782628,0.007751937984496124,0.000949414628985708,0.00043283472005928303,0.00970873786407767,1 diff --git a/results/certified-maximin-design/n45/cut_trace.csv b/results/certified-maximin-design/n45/cut_trace.csv new file mode 100644 index 0000000..b0d9984 --- /dev/null +++ b/results/certified-maximin-design/n45/cut_trace.csv @@ -0,0 +1,8 @@ +iteration,status,message,cut_count,master_incumbent,master_upper_bound,global_upper_bound,candidate_lower_bound,best_lower_bound,absolute_gap,relative_gap,violated_pair_count,new_cut_count +1,0,Optimization terminated successfully. (HiGHS Status 7: Optimal),40,0.29969132172795643,0.29969132172795643,0.29969132172795643,0.018338177433384985,0.24983023368504692,0.04986108804290951,0.04986108804290951,10,10 +2,0,Optimization terminated successfully. (HiGHS Status 7: Optimal),43,0.2736078369096001,0.2736078369096001,0.2736078369096001,0.20997088490389842,0.24983023368504692,0.02377760322455319,0.02377760322455319,3,3 +3,0,Optimization terminated successfully. (HiGHS Status 7: Optimal),46,0.2601656411587341,0.2601656411587341,0.2601656411587341,0.2472594709475595,0.24983023368504692,0.010335407473687164,0.010335407473687164,3,3 +4,0,Optimization terminated successfully. (HiGHS Status 7: Optimal),50,0.25432101764237947,0.25432101764237947,0.25432101764237947,0.16308942112856467,0.24983023368504692,0.004490783957332545,0.004490783957332545,4,4 +5,0,Optimization terminated successfully. (HiGHS Status 7: Optimal),51,0.2523958379330045,0.2523958379330045,0.2523958379330045,0.2216371103803767,0.24983023368504692,0.002565604247957598,0.002565604247957598,1,1 +6,0,Optimization terminated successfully. (HiGHS Status 7: Optimal),53,0.2523958379330047,0.2523958379330047,0.2523958379330045,0.23336287252923216,0.24983023368504692,0.002565604247957598,0.002565604247957598,2,2 +7,0,Optimization terminated successfully. (HiGHS Status 7: Optimal),53,0.2523958379330047,0.2523958379330047,0.2523958379330045,0.25239583793300424,0.25239583793300424,2.7755575615628914e-16,2.7755575615628914e-16,0,0 diff --git a/results/certified-maximin-design/n60/SHA256SUMS.csv b/results/certified-maximin-design/n60/SHA256SUMS.csv new file mode 100644 index 0000000..368f037 --- /dev/null +++ b/results/certified-maximin-design/n60/SHA256SUMS.csv @@ -0,0 +1,4 @@ +file,bytes,sha256 +certificate_summary.json,1968,7c89b45552568665fd22103feaa301d9d69030bc1a5101b20c3ce6d438b3743c +certified_allocation.csv,1925,694a84107c1ba94f39146a05675def1b5813621fab486d7567f6711510a6fc08 +cut_trace.csv,1264,34d2b9598f96960ca483bfc4054f5fc975c4aecbdf0fade0a02d28f734ac759f diff --git a/results/certified-maximin-design/n60/certificate_summary.json b/results/certified-maximin-design/n60/certificate_summary.json new file mode 100644 index 0000000..0f74aa3 --- /dev/null +++ b/results/certified-maximin-design/n60/certificate_summary.json @@ -0,0 +1,50 @@ +{ + "absolute_gap": 5.551115123125783e-17, + "candidate_names": [ + "innovation_l2", + "surprise", + "gain", + "update_l2", + "information_gain", + "change_probability" + ], + "certified": true, + "certified_improvement_over_heuristic": 0.0028311171816116643, + "code_repository": "IPS-Stuttgart/Bayesian-ACh", + "code_sha": "e82f48d335fb3f17f843d16bf2442693b6c83eef", + "config": { + "absolute_gap_tolerance": 1e-08, + "budget": 60, + "cut_violation_tolerance": 1e-09, + "integer": true, + "master_mip_relative_gap": 1e-09, + "master_time_limit_s": 120.0, + "max_iterations": 100, + "max_point_fraction": 0.15, + "relative_gap_tolerance": 1e-07 + }, + "cut_count": 52, + "direct_geometry_matches_lower_bound": true, + "experiment": "certified_finite_grid_maximin_allocation", + "grid_point_count": 240, + "heuristic_lower_bound": 0.25003654826964467, + "iterations": 5, + "last_master_message": "Optimization terminated successfully. (HiGHS Status 7: Optimal)", + "last_master_status": 0, + "lower_bound": 0.25286766545125633, + "population_n_eff": { + "effect_size": 1.0, + "index_certified": true, + "interpretation": "Population observation-equivalent index under effectively independent Gaussian observations; not a physical trial, time-bin, session, or animal count.", + "lower_index_from_residual_upper_bound": 45, + "noise_std": 1.0, + "target_log_score_gap": 5.0, + "upper_index_from_residual_lower_bound": 45 + }, + "relative_gap": 5.551115123125783e-17, + "schema_version": 1, + "scipy_version": "1.15.3", + "scope": "This certificate is conditional on the independently instantiated finite grid, global candidate standardization, budget, and per-cell cap. The allocation is a count target, not a validated sequential protocol; no history constructor, washout, reset, or carry-over feasibility is implied.", + "solver": "scipy.optimize.milp (HiGHS) with OLS cutting planes", + "upper_bound": 0.2528676654512564 +} diff --git a/results/certified-maximin-design/n60/certified_allocation.csv b/results/certified-maximin-design/n60/certified_allocation.csv new file mode 100644 index 0000000..2fbb7fd --- /dev/null +++ b/results/certified-maximin-design/n60/certified_allocation.csv @@ -0,0 +1,13 @@ +point_id,observed_probability,residual_shape,concentration,reset_observed_probability,hazard,innovation_l2,surprise,gain,update_l2,information_gain,change_probability,allocation +0,0.05,0.5,2.0,0.05,0.01,1.1635076278220096,2.995732273553991,0.3333333333333333,0.3878358759406699,1.6491929980444475,0.01,1 +3,0.05,0.5,2.0,0.5,0.15,1.1635076278220096,2.995732273553991,0.3333333333333333,0.3878358759406699,1.6491929980444475,0.6382978723404256,1 +7,0.05,0.5,8.0,0.05,0.15,1.1635076278220096,2.995732273553991,0.1111111111111111,0.1292786253135566,0.793706251013266,0.15,2 +11,0.05,0.5,8.0,0.95,0.15,1.1635076278220096,2.995732273553991,0.1111111111111111,0.1292786253135566,0.793706251013266,0.7702702702702702,6 +16,0.05,0.5,32.0,0.95,0.01,1.1635076278220096,2.995732273553991,0.030303030303030304,0.03525780690369727,0.2655001957924723,0.16101694915254236,3 +121,0.35,0.9,2.0,0.05,0.15,0.8768979416100827,1.0498221244986778,0.3333333333333333,0.29229931387002756,0.3355856642737046,0.02459016393442623,5 +134,0.35,0.9,32.0,0.5,0.01,0.8768979416100827,1.0498221244986778,0.030303030303030304,0.02657266489727523,0.02843542878084193,0.014224751066856332,6 +139,0.35,0.9,128.0,0.05,0.15,0.8768979416100827,1.0498221244986778,0.007751937984496124,0.0067976584620936574,0.007218032071614644,0.02459016393442623,9 +140,0.35,0.9,128.0,0.5,0.01,0.8768979416100827,1.0498221244986778,0.007751937984496124,0.0067976584620936574,0.007218032071614644,0.014224751066856332,9 +197,0.9,0.5,2.0,0.95,0.15,0.12247448713915887,0.10536051565782628,0.3333333333333333,0.0408248290463862,0.023123169408776692,0.15702479338842976,4 +205,0.9,0.5,32.0,0.05,0.15,0.12247448713915887,0.10536051565782628,0.030303030303030304,0.0037113480951260457,0.0017170260888290656,0.00970873786407767,5 +210,0.9,0.5,128.0,0.05,0.01,0.12247448713915887,0.10536051565782628,0.007751937984496124,0.000949414628985708,0.00043283472005928303,0.0005608524957936063,9 diff --git a/results/certified-maximin-design/n60/cut_trace.csv b/results/certified-maximin-design/n60/cut_trace.csv new file mode 100644 index 0000000..8c832ab --- /dev/null +++ b/results/certified-maximin-design/n60/cut_trace.csv @@ -0,0 +1,6 @@ +iteration,status,message,cut_count,master_incumbent,master_upper_bound,global_upper_bound,candidate_lower_bound,best_lower_bound,absolute_gap,relative_gap,violated_pair_count,new_cut_count +1,0,Optimization terminated successfully. (HiGHS Status 7: Optimal),36,0.3003019448952873,0.3003019448952873,0.3003019448952873,1.8899792520920073e-32,0.25003654826964467,0.05026539662564261,0.05026539662564261,6,6 +2,0,Optimization terminated successfully. (HiGHS Status 7: Optimal),40,0.2746209089988777,0.2746209089988777,0.2746209089988777,0.1998940141548362,0.25003654826964467,0.024584360729233035,0.024584360729233035,4,4 +3,0,Optimization terminated successfully. (HiGHS Status 7: Optimal),43,0.2625384062229467,0.2625384062229467,0.2625384062229467,0.1845235932375922,0.25003654826964467,0.012501857953302054,0.012501857953302054,3,3 +4,0,Optimization terminated successfully. (HiGHS Status 7: Optimal),45,0.2556738948533919,0.2556738948533919,0.2556738948533919,0.25286766545125633,0.25286766545125633,0.002806229402135596,0.002806229402135596,2,2 +5,0,Optimization terminated successfully. (HiGHS Status 7: Optimal),52,0.2528676654512564,0.2528676654512564,0.2528676654512564,0.166514444222626,0.25286766545125633,5.551115123125783e-17,5.551115123125783e-17,7,7 From c6dea54d1e8a1b702c212c83535af8a5de3ab3ea Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 06:04:29 +0200 Subject: [PATCH 24/25] Add paired maximin recovery freeze contract --- docs/certified_maximin_paired_recovery.md | 20 + pyproject.toml | 2 + src/bayesian_ach/design_paired_recovery.py | 396 ++++++++++++++++++ .../design_paired_recovery_verify.py | 241 +++++++++++ tests/test_design_paired_recovery.py | 96 +++++ 5 files changed, 755 insertions(+) create mode 100644 docs/certified_maximin_paired_recovery.md create mode 100644 src/bayesian_ach/design_paired_recovery.py create mode 100644 src/bayesian_ach/design_paired_recovery_verify.py create mode 100644 tests/test_design_paired_recovery.py diff --git a/docs/certified_maximin_paired_recovery.md b/docs/certified_maximin_paired_recovery.md new file mode 100644 index 0000000..eb6701f --- /dev/null +++ b/docs/certified_maximin_paired_recovery.md @@ -0,0 +1,20 @@ +# Paired heuristic-versus-certified recovery diagnostic + +The exact integer certificate maximizes the declared finite-grid worst ordered-pair +residual. It does not follow that the certified allocation must outperform the earlier +heuristic allocation in a finite held-out recovery simulation. + +`bayesian-ach-design-paired-recovery` therefore compares the chronologically locked +heuristic N=60 allocation with the exact certified N=60 allocation. The comparison uses +the same five seeds (`7, 11, 19, 23, 31`), 200 replicates per generating candidate and +seed, a 0.35 held-out fraction, unit effect and noise scales, and an identically reset RNG +stream for both allocations. Both source allocations are checksum-bound; the heuristic +constructor and integer certificate are independently reverified before simulation. + +`bayesian-ach-design-paired-recovery-verify` rechecks every output checksum, both source +allocations, certificate geometry, all 60 recovery rows, five seed summaries, and every +headline allocation/recovery diagnostic by rerunning the deterministic simulation. + +This is a diagnostic sensitivity analysis. It does not replace the chronologically locked +paper allocation, select a post-result schedule, establish empirical superiority, guarantee +recovery power, or specify a physical trial, time-bin, or animal protocol. diff --git a/pyproject.toml b/pyproject.toml index acf6f9c..4861caa 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -72,6 +72,8 @@ bayesian-ach = "bayesian_ach.cli_ext:main" bayesian-ach-design = "bayesian_ach.design_cli:main" bayesian-ach-design-certificate-verify = "bayesian_ach.design_certificate_verify:main" bayesian-ach-design-certify = "bayesian_ach.design_certificate_cli:main" +bayesian-ach-design-paired-recovery = "bayesian_ach.design_paired_recovery:main" +bayesian-ach-design-paired-recovery-verify = "bayesian_ach.design_paired_recovery_verify:main" bayesian-ach-replay = "bayesian_ach.replay_cli:main" [tool.setuptools.packages.find] diff --git a/src/bayesian_ach/design_paired_recovery.py b/src/bayesian_ach/design_paired_recovery.py new file mode 100644 index 0000000..a767a2a --- /dev/null +++ b/src/bayesian_ach/design_paired_recovery.py @@ -0,0 +1,396 @@ +"""Checksum-bound paired recovery for heuristic and certified maximin designs.""" + +from __future__ import annotations + +import argparse +import csv +import hashlib +import json +import math +import re +import shutil +import subprocess +from collections.abc import Mapping, Sequence +from dataclasses import asdict, dataclass +from pathlib import Path +from typing import Any + +import numpy as np +from numpy.typing import NDArray + +from bayesian_ach.design_certificate_verify import verify_certificate_package +from bayesian_ach.design_grid import generate_transition_design_grid +from bayesian_ach.design_optimizer import optimize_maximin_design +from bayesian_ach.design_recovery import recover_design + +_SHA40 = re.compile(r"^[0-9a-f]{40}$") +_SHA64 = re.compile(r"^[0-9a-f]{64}$") +_BUDGET = 60 +_MAX_POINT_FRACTION = 0.15 +_VARIANTS = ("heuristic_maximin", "certified_integer_maximin") + + +@dataclass(frozen=True, slots=True) +class PairedRecoveryConfig: + """Frozen paired-recovery settings.""" + + seeds: tuple[int, ...] = (7, 11, 19, 23, 31) + replicates: int = 200 + test_fraction: float = 0.35 + effect_size: float = 1.0 + noise_std: float = 1.0 + recovery_seed_offset: int = 100_003 + + def validate(self) -> None: + if len(self.seeds) != 5 or len(set(self.seeds)) != len(self.seeds): + raise ValueError("paired recovery requires five unique seeds") + if self.replicates < 20: + raise ValueError("replicates must be at least 20") + if not 0.0 < self.test_fraction < 1.0: + raise ValueError("test_fraction must lie in (0, 1)") + if ( + not math.isfinite(self.effect_size) + or not math.isfinite(self.noise_std) + or self.effect_size <= 0.0 + or self.noise_std <= 0.0 + ): + raise ValueError("effect_size and noise_std must be finite and positive") + + +def sha256(path: Path) -> str: + """Return a streaming SHA-256 digest.""" + + digest = hashlib.sha256() + with path.open("rb") as handle: + for chunk in iter(lambda: handle.read(1024 * 1024), b""): + digest.update(chunk) + return digest.hexdigest() + + +def _git_provenance(repo_root: Path, expected_sha: str) -> None: + if _SHA40.fullmatch(expected_sha) is None: + raise ValueError("--code-sha must be a lowercase 40-character commit SHA") + head = subprocess.run( + ["git", "-C", str(repo_root), "rev-parse", "HEAD"], + check=True, + capture_output=True, + text=True, + ).stdout.strip() + status = subprocess.run( + ["git", "-C", str(repo_root), "status", "--porcelain", "--untracked-files=all"], + check=True, + capture_output=True, + text=True, + ).stdout + if head != expected_sha or status: + raise RuntimeError("paired recovery requires the clean exact --code-sha checkout") + + +def _integer_count(value: str) -> int: + number = float(value) + rounded = int(round(number)) + if rounded <= 0 or not math.isclose(number, rounded, abs_tol=1.0e-9): + raise ValueError("allocation contains a nonpositive or noninteger count") + return rounded + + +def load_heuristic_allocation( + path: Path, + *, + expected_sha256: str, + source_code_sha: str, +) -> tuple[NDArray[np.int64], dict[str, Any]]: + """Load and reconstruct the chronologically locked heuristic N=60 allocation.""" + + path = path.resolve() + if _SHA64.fullmatch(expected_sha256) is None or sha256(path) != expected_sha256: + raise ValueError("heuristic allocation SHA-256 mismatch") + if _SHA40.fullmatch(source_code_sha) is None: + raise ValueError("heuristic source code SHA must contain 40 lowercase hex characters") + _, _, signals = generate_transition_design_grid() + counts = np.zeros(signals.shape[0], dtype=np.int64) + seen: set[int] = set() + with path.open(newline="", encoding="utf-8") as handle: + rows = list(csv.DictReader(handle)) + if not rows: + raise ValueError("heuristic allocation table is empty") + for row in rows: + if row.get("design") != "maximin_optimized": + continue + point = int(row["point_id"]) + if point < 0 or point >= counts.size or point in seen: + raise ValueError("heuristic allocation has an invalid or duplicate point") + counts[point] = _integer_count(str(row["count"])) + seen.add(point) + if int(np.sum(counts)) != _BUDGET: + raise ValueError("heuristic maximin allocation does not sum to N=60") + reconstructed = optimize_maximin_design( + signals, + _BUDGET, + max_point_fraction=_MAX_POINT_FRACTION, + effect_size=1.0, + noise_std=1.0, + target_log_score_gap=5.0, + ).counts + if not np.array_equal(counts, reconstructed): + raise ValueError("heuristic allocation does not match the frozen constructor") + return counts, { + "kind": "chronologically_locked_heuristic_maximin_n60", + "source_repository": "IPS-Stuttgart/Bayesian-ACh", + "source_code_sha": source_code_sha, + "source_file": path.name, + "source_bytes": path.stat().st_size, + "source_sha256": expected_sha256, + "budget": _BUDGET, + "max_point_fraction": _MAX_POINT_FRACTION, + "constructor_reproduced": True, + } + + +def _certified_counts(path: Path, point_count: int) -> NDArray[np.int64]: + counts = np.zeros(point_count, dtype=np.int64) + seen: set[int] = set() + with path.open(newline="", encoding="utf-8") as handle: + rows = list(csv.DictReader(handle)) + for row in rows: + point = int(row["point_id"]) + if point < 0 or point >= point_count or point in seen: + raise ValueError("certified allocation has an invalid or duplicate point") + counts[point] = _integer_count(str(row["allocation"])) + seen.add(point) + if int(np.sum(counts)) != _BUDGET: + raise ValueError("certified allocation does not sum to N=60") + return counts + + +def load_certified_allocation( + package: Path, +) -> tuple[NDArray[np.int64], dict[str, Any]]: + """Verify and load an exact integer N=60 certificate package.""" + + package = package.resolve() + report = verify_certificate_package(package) + if int(report["budget"]) != _BUDGET: + raise ValueError("paired recovery requires an N=60 certificate package") + summary_path = package / "certificate_summary.json" + allocation_path = package / "certified_allocation.csv" + summary = json.loads(summary_path.read_text(encoding="utf-8")) + _, _, signals = generate_transition_design_grid() + counts = _certified_counts(allocation_path, signals.shape[0]) + return counts, { + "kind": "certified_integer_maximin_n60", + "source_repository": "IPS-Stuttgart/Bayesian-ACh", + "certificate_code_sha": summary["code_sha"], + "certificate_summary_sha256": sha256(summary_path), + "certificate_allocation_sha256": sha256(allocation_path), + "certificate_sha256sums_sha256": sha256(package / "SHA256SUMS.csv"), + "budget": _BUDGET, + "lower_bound": summary["lower_bound"], + "upper_bound": summary["upper_bound"], + "certified": True, + } + + +def paired_recovery_rows( + heuristic: NDArray[np.int64], + certified: NDArray[np.int64], + config: PairedRecoveryConfig, +) -> tuple[list[dict[str, Any]], list[dict[str, Any]], dict[str, Any]]: + """Run both schedules with identical random streams for every frozen seed.""" + + config.validate() + _, _, signals = generate_transition_design_grid() + variants = { + "heuristic_maximin": heuristic, + "certified_integer_maximin": certified, + } + recovery_rows: list[dict[str, Any]] = [] + seed_rows: list[dict[str, Any]] = [] + for seed in config.seeds: + by_variant: dict[str, list[float]] = {} + for name, counts in variants.items(): + records = recover_design( + name, + signals, + counts, + replicates=config.replicates, + test_fraction=config.test_fraction, + effect_size=config.effect_size, + noise_std=config.noise_std, + seed=seed + config.recovery_seed_offset, + ) + rates = [float(record.recovery_rate) for record in records] + by_variant[name] = rates + recovery_rows.extend({"seed": seed, **record.as_dict()} for record in records) + heuristic_rates = by_variant["heuristic_maximin"] + certified_rates = by_variant["certified_integer_maximin"] + seed_rows.append( + { + "seed": seed, + "heuristic_minimum_recovery": min(heuristic_rates), + "certified_minimum_recovery": min(certified_rates), + "minimum_recovery_delta_certified_minus_heuristic": ( + min(certified_rates) - min(heuristic_rates) + ), + "heuristic_mean_recovery": float(np.mean(heuristic_rates)), + "certified_mean_recovery": float(np.mean(certified_rates)), + "mean_recovery_delta_certified_minus_heuristic": ( + float(np.mean(certified_rates)) - float(np.mean(heuristic_rates)) + ), + } + ) + difference = certified - heuristic + overlap = int(np.sum((certified > 0) & (heuristic > 0))) + diagnostics = { + "heuristic_support_count": int(np.count_nonzero(heuristic)), + "certified_support_count": int(np.count_nonzero(certified)), + "support_overlap_count": overlap, + "support_union_count": int(np.sum((certified > 0) | (heuristic > 0))), + "allocation_l1_distance": int(np.sum(np.abs(difference))), + "maximum_absolute_cell_change": int(np.max(np.abs(difference))), + "minimum_recovery_across_seeds": { + "heuristic_maximin": min( + float(row["heuristic_minimum_recovery"]) for row in seed_rows + ), + "certified_integer_maximin": min( + float(row["certified_minimum_recovery"]) for row in seed_rows + ), + }, + "all_seed_minima_certified_minus_heuristic_negative": all( + float(row["minimum_recovery_delta_certified_minus_heuristic"]) < 0.0 + for row in seed_rows + ), + } + return recovery_rows, seed_rows, diagnostics + + +def _write_csv(path: Path, rows: Sequence[Mapping[str, Any]]) -> None: + if not rows: + raise ValueError(f"cannot write empty table: {path}") + with path.open("w", newline="", encoding="utf-8") as handle: + writer = csv.DictWriter(handle, fieldnames=list(rows[0]), lineterminator="\n") + writer.writeheader() + writer.writerows(rows) + + +def freeze_paired_recovery( + *, + output: Path, + repo_root: Path, + code_sha: str, + heuristic_path: Path, + heuristic_sha256: str, + heuristic_source_code_sha: str, + certificate_package: Path, + config: PairedRecoveryConfig, +) -> dict[str, Any]: + """Generate a self-checksummed paired-recovery evidence package.""" + + _git_provenance(repo_root.resolve(), code_sha) + heuristic, heuristic_provenance = load_heuristic_allocation( + heuristic_path, + expected_sha256=heuristic_sha256, + source_code_sha=heuristic_source_code_sha, + ) + certified, certificate_provenance = load_certified_allocation(certificate_package) + recovery, seeds, diagnostics = paired_recovery_rows(heuristic, certified, config) + summary = { + "schema_version": 1, + "experiment": "paired_heuristic_vs_certified_n60_recovery_diagnostic", + "producer_commit": code_sha, + "producer_clean_worktree": True, + "config": asdict(config), + "paired_rng": True, + "recovery_seed_rule": ( + "seed + recovery_seed_offset; RNG reset identically for each allocation" + ), + "input_provenance": [heuristic_provenance, certificate_provenance], + **diagnostics, + "interpretation": ( + "The exact integer certificate optimizes the frozen asymptotic worst-residual " + "objective. This paired finite-sample diagnostic does not select a replacement " + "schedule, establish empirical superiority, or define a physical trial/animal protocol." + ), + } + output = output.resolve() + output.mkdir(parents=True, exist_ok=False) + shutil.copyfile(heuristic_path, output / "heuristic_allocation_source.csv") + _write_csv(output / "paired_recovery.csv", recovery) + _write_csv(output / "seed_summary.csv", seeds) + (output / "summary.json").write_text( + json.dumps(summary, indent=2, sort_keys=True, allow_nan=False) + "\n", + encoding="utf-8", + ) + payload_names = ( + "heuristic_allocation_source.csv", + "paired_recovery.csv", + "seed_summary.csv", + "summary.json", + ) + manifest = { + "schema_version": 1, + "artifact": summary["experiment"], + "producer_commit": code_sha, + "producer_git_dirty": False, + "inputs": summary["input_provenance"], + "files": [ + { + "path": name, + "bytes": (output / name).stat().st_size, + "sha256": sha256(output / name), + } + for name in payload_names + ], + } + (output / "artifact_manifest.json").write_text( + json.dumps(manifest, indent=2, sort_keys=True, allow_nan=False) + "\n", + encoding="utf-8", + ) + checksum_names = (*payload_names, "artifact_manifest.json") + _write_csv( + output / "SHA256SUMS.csv", + [ + { + "file": name, + "bytes": (output / name).stat().st_size, + "sha256": sha256(output / name), + } + for name in checksum_names + ], + ) + return summary + + +def _parser() -> argparse.ArgumentParser: + parser = argparse.ArgumentParser( + prog="bayesian-ach-design-paired-recovery", + description="Freeze paired heuristic-versus-certified N=60 recovery diagnostics.", + ) + parser.add_argument("--output", required=True, type=Path) + parser.add_argument("--repo-root", default=Path.cwd(), type=Path) + parser.add_argument("--code-sha", required=True) + parser.add_argument("--heuristic-allocation", required=True, type=Path) + parser.add_argument("--heuristic-allocation-sha256", required=True) + parser.add_argument("--heuristic-source-code-sha", required=True) + parser.add_argument("--certificate-package", required=True, type=Path) + return parser + + +def main(argv: Sequence[str] | None = None) -> int: + args = _parser().parse_args(argv) + summary = freeze_paired_recovery( + output=args.output, + repo_root=args.repo_root, + code_sha=args.code_sha, + heuristic_path=args.heuristic_allocation, + heuristic_sha256=args.heuristic_allocation_sha256, + heuristic_source_code_sha=args.heuristic_source_code_sha, + certificate_package=args.certificate_package, + config=PairedRecoveryConfig(), + ) + print(json.dumps(summary, indent=2, sort_keys=True, allow_nan=False)) + return 0 + + +if __name__ == "__main__": # pragma: no cover + raise SystemExit(main()) diff --git a/src/bayesian_ach/design_paired_recovery_verify.py b/src/bayesian_ach/design_paired_recovery_verify.py new file mode 100644 index 0000000..bbbb295 --- /dev/null +++ b/src/bayesian_ach/design_paired_recovery_verify.py @@ -0,0 +1,241 @@ +"""Independent verification for frozen paired maximin recovery evidence.""" + +from __future__ import annotations + +import argparse +import csv +import json +import math +import re +from collections.abc import Mapping, Sequence +from pathlib import Path +from typing import Any + +from bayesian_ach.design_paired_recovery import ( + PairedRecoveryConfig, + load_certified_allocation, + load_heuristic_allocation, + paired_recovery_rows, + sha256, +) + +_PAYLOADS = { + "heuristic_allocation_source.csv", + "paired_recovery.csv", + "seed_summary.csv", + "summary.json", + "artifact_manifest.json", +} +_SHA40 = re.compile(r"^[0-9a-f]{40}$") + + +def _locked_payloads(directory: Path) -> dict[str, dict[str, str]]: + checksum = directory / "SHA256SUMS.csv" + with checksum.open(newline="", encoding="utf-8") as handle: + rows = list(csv.DictReader(handle)) + locked: dict[str, dict[str, str]] = {} + for row in rows: + name = str(row["file"]) + if name in locked: + raise ValueError(f"duplicate checksum row: {name}") + locked[name] = {key: str(value) for key, value in row.items()} + if set(locked) != _PAYLOADS: + raise ValueError("checksum table does not bind the exact paired-recovery payload set") + for name, row in locked.items(): + path = directory / name + if not path.is_file(): + raise ValueError(f"missing locked payload: {name}") + if int(row["bytes"]) != path.stat().st_size or row["sha256"] != sha256(path): + raise ValueError(f"locked paired-recovery payload mismatch: {name}") + return locked + + +def _provenance_item(summary: Mapping[str, Any], kind: str) -> Mapping[str, Any]: + items = [item for item in summary["input_provenance"] if item.get("kind") == kind] + if len(items) != 1: + raise ValueError(f"expected one provenance item of kind {kind}") + return items[0] + + +def _compare_rows( + path: Path, + expected: Sequence[Mapping[str, Any]], +) -> None: + with path.open(newline="", encoding="utf-8") as handle: + rows = list(csv.DictReader(handle)) + if len(rows) != len(expected): + raise ValueError(f"row count mismatch: {path.name}") + for observed, target in zip(rows, expected, strict=True): + if set(observed) != set(target): + raise ValueError(f"column mismatch: {path.name}") + for key, value in target.items(): + raw = str(observed[key]) + if isinstance(value, bool): + equal = raw == str(value) + elif isinstance(value, int): + equal = int(raw) == value + elif isinstance(value, float): + equal = math.isclose( + float(raw), + value, + rel_tol=1.0e-12, + abs_tol=1.0e-12, + ) + else: + equal = raw == str(value) + if not equal: + raise ValueError(f"value mismatch in {path.name}: {key}") + + +def _compare_summary_value(observed: Any, expected: Any, name: str) -> None: + if isinstance(expected, bool): + equal = observed is expected + elif isinstance(expected, int): + equal = int(observed) == expected + elif isinstance(expected, float): + equal = math.isclose( + float(observed), + expected, + rel_tol=1.0e-12, + abs_tol=1.0e-12, + ) + elif isinstance(expected, Mapping): + equal = set(observed) == set(expected) + if equal: + for key, value in expected.items(): + _compare_summary_value(observed[key], value, f"{name}.{key}") + return + else: + equal = observed == expected + if not equal: + raise ValueError(f"summary mismatch: {name}") + + +def verify_paired_recovery_package( + directory: Path, + *, + certificate_package: Path, +) -> dict[str, Any]: + """Recompute checksums, inputs, simulations, and all headline diagnostics.""" + + directory = directory.resolve() + locked = _locked_payloads(directory) + summary = json.loads((directory / "summary.json").read_text(encoding="utf-8")) + manifest = json.loads( + (directory / "artifact_manifest.json").read_text(encoding="utf-8") + ) + producer = str(summary.get("producer_commit", "")) + if ( + summary.get("schema_version") != 1 + or summary.get("experiment") + != "paired_heuristic_vs_certified_n60_recovery_diagnostic" + or summary.get("producer_clean_worktree") is not True + or _SHA40.fullmatch(producer) is None + ): + raise ValueError("invalid paired-recovery summary contract") + if ( + manifest.get("schema_version") != 1 + or manifest.get("producer_commit") != producer + or manifest.get("producer_git_dirty") is not False + or set(item["path"] for item in manifest["files"]) + != _PAYLOADS - {"artifact_manifest.json"} + ): + raise ValueError("invalid paired-recovery manifest contract") + for item in manifest["files"]: + row = locked[item["path"]] + if int(item["bytes"]) != int(row["bytes"]) or item["sha256"] != row["sha256"]: + raise ValueError("manifest and checksum table disagree") + + config_values = dict(summary["config"]) + config_values["seeds"] = tuple(config_values["seeds"]) + config = PairedRecoveryConfig(**config_values) + config.validate() + if config != PairedRecoveryConfig(): + raise ValueError("paired-recovery config differs from the frozen default") + heuristic_input = _provenance_item( + summary, + "chronologically_locked_heuristic_maximin_n60", + ) + heuristic, heuristic_report = load_heuristic_allocation( + directory / "heuristic_allocation_source.csv", + expected_sha256=str(heuristic_input["source_sha256"]), + source_code_sha=str(heuristic_input["source_code_sha"]), + ) + certified, certificate_report = load_certified_allocation(certificate_package) + for key in ( + "source_code_sha", + "source_sha256", + "budget", + "max_point_fraction", + "constructor_reproduced", + ): + _compare_summary_value(heuristic_input[key], heuristic_report[key], f"heuristic.{key}") + certificate_input = _provenance_item(summary, "certified_integer_maximin_n60") + for key in ( + "certificate_code_sha", + "certificate_summary_sha256", + "certificate_allocation_sha256", + "certificate_sha256sums_sha256", + "budget", + "lower_bound", + "upper_bound", + "certified", + ): + _compare_summary_value( + certificate_input[key], + certificate_report[key], + f"certificate.{key}", + ) + + recovery_rows, seed_rows, diagnostics = paired_recovery_rows( + heuristic, + certified, + config, + ) + _compare_rows(directory / "paired_recovery.csv", recovery_rows) + _compare_rows(directory / "seed_summary.csv", seed_rows) + for key, value in diagnostics.items(): + _compare_summary_value(summary[key], value, key) + expected_count = len(config.seeds) * 2 * 6 + if len(recovery_rows) != expected_count: + raise ValueError("paired-recovery table does not cover every seed/design/generator") + return { + "package": directory.name, + "verified": True, + "producer_commit": producer, + "seed_count": len(config.seeds), + "replicates_per_generator_per_seed": config.replicates, + "recovery_row_count": expected_count, + "allocation_l1_distance": diagnostics["allocation_l1_distance"], + "heuristic_minimum_recovery": diagnostics["minimum_recovery_across_seeds"][ + "heuristic_maximin" + ], + "certified_minimum_recovery": diagnostics["minimum_recovery_across_seeds"][ + "certified_integer_maximin" + ], + "sha256sums_sha256": sha256(directory / "SHA256SUMS.csv"), + } + + +def _parser() -> argparse.ArgumentParser: + parser = argparse.ArgumentParser( + prog="bayesian-ach-design-paired-recovery-verify", + description="Verify frozen paired maximin recovery evidence.", + ) + parser.add_argument("package", type=Path) + parser.add_argument("--certificate-package", required=True, type=Path) + return parser + + +def main(argv: Sequence[str] | None = None) -> int: + args = _parser().parse_args(argv) + report = verify_paired_recovery_package( + args.package, + certificate_package=args.certificate_package, + ) + print(json.dumps(report, indent=2, sort_keys=True, allow_nan=False)) + return 0 + + +if __name__ == "__main__": # pragma: no cover + raise SystemExit(main()) diff --git a/tests/test_design_paired_recovery.py b/tests/test_design_paired_recovery.py new file mode 100644 index 0000000..1af06a9 --- /dev/null +++ b/tests/test_design_paired_recovery.py @@ -0,0 +1,96 @@ +from __future__ import annotations + +import csv +from pathlib import Path + +import numpy as np +import pytest + +from bayesian_ach import design_paired_recovery as paired +from bayesian_ach.design_grid import generate_transition_design_grid +from bayesian_ach.design_optimizer import optimize_maximin_design +from bayesian_ach.design_paired_recovery_verify import verify_paired_recovery_package + + +def _heuristic_source(path: Path) -> Path: + _, _, signals = generate_transition_design_grid() + counts = optimize_maximin_design( + signals, + 60, + max_point_fraction=0.15, + effect_size=1.0, + noise_std=1.0, + target_log_score_gap=5.0, + ).counts + with path.open("w", newline="", encoding="utf-8") as handle: + writer = csv.DictWriter(handle, fieldnames=["design", "count", "point_id"]) + writer.writeheader() + for point in np.flatnonzero(counts): + writer.writerow( + { + "design": "maximin_optimized", + "count": int(counts[point]), + "point_id": int(point), + } + ) + return path + + +def _certificate_package() -> Path: + return Path(__file__).resolve().parents[1] / "results/certified-maximin-design/n60" + + +def test_locked_heuristic_is_hash_bound_and_reconstructed(tmp_path: Path) -> None: + path = _heuristic_source(tmp_path / "heuristic.csv") + digest = paired.sha256(path) + counts, report = paired.load_heuristic_allocation( + path, + expected_sha256=digest, + source_code_sha="1" * 40, + ) + assert int(np.sum(counts)) == 60 + assert report["constructor_reproduced"] is True + path.write_text(path.read_text(encoding="utf-8") + "\n", encoding="utf-8") + with pytest.raises(ValueError, match="SHA-256"): + paired.load_heuristic_allocation( + path, + expected_sha256=digest, + source_code_sha="1" * 40, + ) + + +def test_paired_package_recomputes_and_detects_tamper( + tmp_path: Path, + monkeypatch: pytest.MonkeyPatch, +) -> None: + heuristic = _heuristic_source(tmp_path / "heuristic.csv") + monkeypatch.setattr(paired, "_git_provenance", lambda *_args, **_kwargs: None) + output = tmp_path / "evidence" + summary = paired.freeze_paired_recovery( + output=output, + repo_root=tmp_path, + code_sha="2" * 40, + heuristic_path=heuristic, + heuristic_sha256=paired.sha256(heuristic), + heuristic_source_code_sha="1" * 40, + certificate_package=_certificate_package(), + config=paired.PairedRecoveryConfig(), + ) + assert summary["paired_rng"] is True + report = verify_paired_recovery_package( + output, + certificate_package=_certificate_package(), + ) + assert report["verified"] is True + assert report["recovery_row_count"] == 60 + + seed_summary = output / "seed_summary.csv" + seed_summary.write_text( + seed_summary.read_text(encoding="utf-8").replace("0.83", "0.84", 1), + encoding="utf-8", + ) + with pytest.raises(ValueError, match="payload mismatch"): + verify_paired_recovery_package( + output, + certificate_package=_certificate_package(), + ) From b7e95c0d6d2fb6d2f4c3500c9f36095f8b5eb521 Mon Sep 17 00:00:00 2001 From: Florian Pfaff <6773539+FlorianPfaff@users.noreply.github.com> Date: Tue, 25 Aug 2026 06:11:51 +0200 Subject: [PATCH 25/25] Freeze paired maximin recovery diagnostic --- .gitattributes | 1 + ...ertified_maximin_paired_recovery_result.md | 22 +++++ .../paired-recovery/SHA256SUMS.csv | 6 ++ .../paired-recovery/artifact_manifest.json | 53 +++++++++++ .../heuristic_allocation_source.csv | 87 +++++++++++++++++++ .../paired-recovery/paired_recovery.csv | 61 +++++++++++++ .../paired-recovery/seed_summary.csv | 6 ++ .../paired-recovery/summary.json | 59 +++++++++++++ 8 files changed, 295 insertions(+) create mode 100644 .gitattributes create mode 100644 docs/certified_maximin_paired_recovery_result.md create mode 100644 results/certified-maximin-design/paired-recovery/SHA256SUMS.csv create mode 100644 results/certified-maximin-design/paired-recovery/artifact_manifest.json create mode 100644 results/certified-maximin-design/paired-recovery/heuristic_allocation_source.csv create mode 100644 results/certified-maximin-design/paired-recovery/paired_recovery.csv create mode 100644 results/certified-maximin-design/paired-recovery/seed_summary.csv create mode 100644 results/certified-maximin-design/paired-recovery/summary.json diff --git a/.gitattributes b/.gitattributes new file mode 100644 index 0000000..53c89e8 --- /dev/null +++ b/.gitattributes @@ -0,0 +1 @@ +results/certified-maximin-design/paired-recovery/heuristic_allocation_source.csv -text -diff diff --git a/docs/certified_maximin_paired_recovery_result.md b/docs/certified_maximin_paired_recovery_result.md new file mode 100644 index 0000000..a70c871 --- /dev/null +++ b/docs/certified_maximin_paired_recovery_result.md @@ -0,0 +1,22 @@ +# Frozen paired N=60 recovery result + +The checksum-bound paired diagnostic was produced from clean commit +`c6dea54d1e8a1b702c212c83535af8a5de3ab3ea`. It compares the chronologically locked +heuristic allocation (SHA-256 `a823be49faf6c6cbebf60b11d4b5ca895cf7734d6e9c577ee98f97a5907b69b2`) +with the exact certified allocation (SHA-256 +`694a84107c1ba94f39146a05675def1b5813621fab486d7567f6711510a6fc08`). + +The schedules are materially different: allocation L1 distance 100, support 15 versus 12, +support overlap 4, and maximum absolute cell change 9. With identical random streams for +five seeds and 200 replicates per generating candidate, the minimum recovery rate changed +from 0.750 for the heuristic schedule to 0.665 for the certified schedule. The certified +schedule's minimum was lower for every seed. + +This result is a negative sensitivity boundary. The certificate remains an exact statement +about the declared asymptotic worst-residual objective; it is not evidence that the certified +schedule is superior under finite recovery. The diagnostic does not replace the locked +heuristic schedule, select a post-result design, guarantee power, or define a physical +trial, time-bin, or animal protocol. + +The frozen evidence is in `results/certified-maximin-design/paired-recovery/`. Its verifier +recomputes both inputs, every recovery row, all headline diagnostics, and every checksum. diff --git a/results/certified-maximin-design/paired-recovery/SHA256SUMS.csv b/results/certified-maximin-design/paired-recovery/SHA256SUMS.csv new file mode 100644 index 0000000..c8abc30 --- /dev/null +++ b/results/certified-maximin-design/paired-recovery/SHA256SUMS.csv @@ -0,0 +1,6 @@ +file,bytes,sha256 +heuristic_allocation_source.csv,14450,a823be49faf6c6cbebf60b11d4b5ca895cf7734d6e9c577ee98f97a5907b69b2 +paired_recovery.csv,5515,c04675c240b34a110be3c6a2376775ced4381396026240e96a756bf4e1c94c3c +seed_summary.csv,644,0893d3da81d6d26cd94121dfefc6b19c406edf068184019e1d51a363adff44c8 +summary.json,2309,2abde4359fcf9e6ada150931e9dbb6b74e70cd06ea989998dfe25e42ee7191ed +artifact_manifest.json,1954,b20585a4e29c1a0776a0b327f25bd38c493800f54fbcc0f81b64a5a5feaa4296 diff --git a/results/certified-maximin-design/paired-recovery/artifact_manifest.json b/results/certified-maximin-design/paired-recovery/artifact_manifest.json new file mode 100644 index 0000000..6836462 --- /dev/null +++ b/results/certified-maximin-design/paired-recovery/artifact_manifest.json @@ -0,0 +1,53 @@ +{ + "artifact": "paired_heuristic_vs_certified_n60_recovery_diagnostic", + "files": [ + { + "bytes": 14450, + "path": "heuristic_allocation_source.csv", + "sha256": "a823be49faf6c6cbebf60b11d4b5ca895cf7734d6e9c577ee98f97a5907b69b2" + }, + { + "bytes": 5515, + "path": "paired_recovery.csv", + "sha256": "c04675c240b34a110be3c6a2376775ced4381396026240e96a756bf4e1c94c3c" + }, + { + "bytes": 644, + "path": "seed_summary.csv", + "sha256": "0893d3da81d6d26cd94121dfefc6b19c406edf068184019e1d51a363adff44c8" + }, + { + "bytes": 2309, + "path": "summary.json", + "sha256": "2abde4359fcf9e6ada150931e9dbb6b74e70cd06ea989998dfe25e42ee7191ed" + } + ], + "inputs": [ + { + "budget": 60, + "constructor_reproduced": true, + "kind": "chronologically_locked_heuristic_maximin_n60", + "max_point_fraction": 0.15, + "source_bytes": 14450, + "source_code_sha": "1b2028929ac6ebc1cce0882f0c22af9918044342", + "source_file": "optimal_design_allocation_seed7.csv", + "source_repository": "IPS-Stuttgart/Bayesian-ACh", + "source_sha256": "a823be49faf6c6cbebf60b11d4b5ca895cf7734d6e9c577ee98f97a5907b69b2" + }, + { + "budget": 60, + "certificate_allocation_sha256": "694a84107c1ba94f39146a05675def1b5813621fab486d7567f6711510a6fc08", + "certificate_code_sha": "e82f48d335fb3f17f843d16bf2442693b6c83eef", + "certificate_sha256sums_sha256": "ae337d44c4b5b0ae3bbaaaf1a4981a29bffe7cbb5458104f05d19de13c313938", + "certificate_summary_sha256": "7c89b45552568665fd22103feaa301d9d69030bc1a5101b20c3ce6d438b3743c", + "certified": true, + "kind": "certified_integer_maximin_n60", + "lower_bound": 0.25286766545125633, + "source_repository": "IPS-Stuttgart/Bayesian-ACh", + "upper_bound": 0.2528676654512564 + } + ], + "producer_commit": "c6dea54d1e8a1b702c212c83535af8a5de3ab3ea", + "producer_git_dirty": false, + "schema_version": 1 +} diff --git a/results/certified-maximin-design/paired-recovery/heuristic_allocation_source.csv b/results/certified-maximin-design/paired-recovery/heuristic_allocation_source.csv new file mode 100644 index 0000000..af3261a --- /dev/null +++ b/results/certified-maximin-design/paired-recovery/heuristic_allocation_source.csv @@ -0,0 +1,87 @@ +design,count,point_id,observed_probability,residual_shape,concentration,reset_observed_probability,hazard,innovation_l2,surprise,gain,update_l2,information_gain,change_probability +coupled_novelty,4,29,0.05,0.9,2.0,0.95,0.15,1.281620068507044,2.995732273553991,0.3333333333333333,0.427206689502348,1.6491929980444475,0.7702702702702702 +coupled_novelty,5,77,0.15,0.9,2.0,0.95,0.15,1.1467126928747235,1.8971199848858813,0.3333333333333333,0.3822375642915745,0.8051447609206146,0.5277777777777777 +coupled_novelty,5,83,0.15,0.9,8.0,0.95,0.15,1.1467126928747235,1.8971199848858813,0.1111111111111111,0.12741252143052484,0.30077194367141624,0.5277777777777777 +coupled_novelty,12,129,0.35,0.9,8.0,0.5,0.15,0.8768979416100827,1.0498221244986778,0.1111111111111111,0.09743310462334258,0.10687049253534275,0.20134228187919462 +coupled_novelty,1,131,0.35,0.9,8.0,0.95,0.15,0.8768979416100827,1.0498221244986778,0.1111111111111111,0.09743310462334258,0.10687049253534275,0.32386363636363635 +coupled_novelty,3,134,0.35,0.9,32.0,0.5,0.01,0.8768979416100827,1.0498221244986778,0.030303030303030304,0.02657266489727523,0.02843542878084193,0.014224751066856332 +coupled_novelty,2,135,0.35,0.9,32.0,0.5,0.15,0.8768979416100827,1.0498221244986778,0.030303030303030304,0.02657266489727523,0.02843542878084193,0.20134228187919462 +coupled_novelty,12,158,0.65,0.5,32.0,0.5,0.01,0.4286607049870561,0.4307829160924542,0.030303030303030304,0.012989718332941111,0.008302262397641869,0.007710100231303007 +coupled_novelty,6,162,0.65,0.5,128.0,0.05,0.01,0.4286607049870561,0.4307829160924542,0.007751937984496124,0.0033229512014500266,0.0020964132981635686,0.0007763975155279505 +coupled_novelty,1,182,0.65,0.9,32.0,0.5,0.01,0.4721758147131214,0.4307829160924542,0.030303030303030304,0.014308358021609735,0.008302262397641869,0.007710100231303007 +coupled_novelty,9,210,0.9,0.5,128.0,0.05,0.01,0.12247448713915887,0.10536051565782628,0.007751937984496124,0.000949414628985708,0.00043283472005928303,0.0005608524957936063 +uniform_factorial,1,0,0.05,0.5,2.0,0.05,0.01,1.1635076278220096,2.995732273553991,0.3333333333333333,0.3878358759406699,1.6491929980444475,0.01 +uniform_factorial,1,2,0.05,0.5,2.0,0.5,0.01,1.1635076278220096,2.995732273553991,0.3333333333333333,0.3878358759406699,1.6491929980444475,0.09174311926605505 +uniform_factorial,1,4,0.05,0.5,2.0,0.95,0.01,1.1635076278220096,2.995732273553991,0.3333333333333333,0.3878358759406699,1.6491929980444475,0.16101694915254236 +uniform_factorial,1,5,0.05,0.5,2.0,0.95,0.15,1.1635076278220096,2.995732273553991,0.3333333333333333,0.3878358759406699,1.6491929980444475,0.7702702702702702 +uniform_factorial,1,7,0.05,0.5,8.0,0.05,0.15,1.1635076278220096,2.995732273553991,0.1111111111111111,0.1292786253135566,0.793706251013266,0.15 +uniform_factorial,1,10,0.05,0.5,8.0,0.95,0.01,1.1635076278220096,2.995732273553991,0.1111111111111111,0.1292786253135566,0.793706251013266,0.16101694915254236 +uniform_factorial,1,12,0.05,0.5,32.0,0.05,0.01,1.1635076278220096,2.995732273553991,0.030303030303030304,0.03525780690369727,0.2655001957924723,0.01 +uniform_factorial,1,14,0.05,0.5,32.0,0.5,0.01,1.1635076278220096,2.995732273553991,0.030303030303030304,0.03525780690369727,0.2655001957924723,0.09174311926605505 +uniform_factorial,1,21,0.05,0.5,128.0,0.5,0.15,1.1635076278220096,2.995732273553991,0.007751937984496124,0.009019438975364405,0.07219424175040157,0.6382978723404256 +uniform_factorial,1,28,0.05,0.9,2.0,0.95,0.01,1.281620068507044,2.995732273553991,0.3333333333333333,0.427206689502348,1.6491929980444475,0.16101694915254236 +uniform_factorial,1,35,0.05,0.9,8.0,0.95,0.15,1.281620068507044,2.995732273553991,0.1111111111111111,0.142402229834116,0.7937062510132655,0.7702702702702702 +uniform_factorial,1,41,0.05,0.9,32.0,0.95,0.15,1.281620068507044,2.995732273553991,0.030303030303030304,0.03883697177294071,0.2655001957924723,0.7702702702702702 +uniform_factorial,1,46,0.05,0.9,128.0,0.95,0.01,1.281620068507044,2.995732273553991,0.007751937984496124,0.009935039290752292,0.07219424175045841,0.16101694915254236 +uniform_factorial,1,49,0.15,0.5,2.0,0.05,0.15,1.0410331406828506,1.8971199848858813,0.3333333333333333,0.3470110468942835,0.8051447609206146,0.05555555555555555 +uniform_factorial,1,54,0.15,0.5,8.0,0.05,0.01,1.0410331406828506,1.8971199848858813,0.1111111111111111,0.11567034896476121,0.3007719436714167,0.0033557046979865775 +uniform_factorial,1,55,0.15,0.5,8.0,0.05,0.15,1.0410331406828506,1.8971199848858813,0.1111111111111111,0.11567034896476121,0.3007719436714167,0.05555555555555555 +uniform_factorial,1,57,0.15,0.5,8.0,0.5,0.15,1.0410331406828506,1.8971199848858813,0.1111111111111111,0.11567034896476121,0.3007719436714167,0.37037037037037035 +uniform_factorial,1,58,0.15,0.5,8.0,0.95,0.01,1.0410331406828506,1.8971199848858813,0.1111111111111111,0.11567034896476121,0.3007719436714167,0.060126582278481014 +uniform_factorial,1,62,0.15,0.5,32.0,0.5,0.01,1.0410331406828506,1.8971199848858813,0.030303030303030304,0.03154645880857123,0.08502152841333821,0.03257328990228013 +uniform_factorial,1,70,0.15,0.5,128.0,0.95,0.01,1.0410331406828506,1.8971199848858813,0.007751937984496124,0.008070024346378674,0.02191450800715744,0.060126582278481014 +uniform_factorial,1,71,0.15,0.5,128.0,0.95,0.15,1.0410331406828506,1.8971199848858813,0.007751937984496124,0.008070024346378674,0.02191450800715744,0.5277777777777777 +uniform_factorial,1,74,0.15,0.9,2.0,0.5,0.01,1.1467126928747235,1.8971199848858813,0.3333333333333333,0.3822375642915745,0.8051447609206146,0.03257328990228013 +uniform_factorial,1,90,0.15,0.9,128.0,0.05,0.01,1.1467126928747235,1.8971199848858813,0.007751937984496124,0.008889245681199386,0.02191450800715744,0.0033557046979865775 +uniform_factorial,1,103,0.35,0.5,8.0,0.05,0.15,0.796084166404533,1.0498221244986778,0.1111111111111111,0.08845379626717033,0.10687049253534298,0.02459016393442623 +uniform_factorial,1,108,0.35,0.5,32.0,0.05,0.01,0.796084166404533,1.0498221244986778,0.030303030303030304,0.024123762618319177,0.028435428780849037,0.0014409221902017292 +uniform_factorial,1,111,0.35,0.5,32.0,0.5,0.15,0.796084166404533,1.0498221244986778,0.030303030303030304,0.024123762618319177,0.028435428780849037,0.20134228187919462 +uniform_factorial,1,115,0.35,0.5,128.0,0.05,0.15,0.796084166404533,1.0498221244986778,0.007751937984496124,0.006171195088407225,0.007218032071614644,0.02459016393442623 +uniform_factorial,1,117,0.35,0.5,128.0,0.5,0.15,0.796084166404533,1.0498221244986778,0.007751937984496124,0.006171195088407225,0.007218032071614644,0.20134228187919462 +uniform_factorial,1,132,0.35,0.9,32.0,0.05,0.01,0.8768979416100827,1.0498221244986778,0.030303030303030304,0.02657266489727523,0.02843542878084193,0.0014409221902017292 +uniform_factorial,1,140,0.35,0.9,128.0,0.5,0.01,0.8768979416100827,1.0498221244986778,0.007751937984496124,0.0067976584620936574,0.007218032071614644,0.014224751066856332 +uniform_factorial,1,145,0.65,0.5,2.0,0.05,0.15,0.4286607049870561,0.4307829160924542,0.3333333333333333,0.142886901662352,0.10803846135795697,0.013392857142857144 +uniform_factorial,1,157,0.65,0.5,32.0,0.05,0.15,0.4286607049870561,0.4307829160924542,0.030303030303030304,0.012989718332941111,0.008302262397641869,0.013392857142857144 +uniform_factorial,1,160,0.65,0.5,32.0,0.95,0.01,0.4286607049870561,0.4307829160924542,0.030303030303030304,0.012989718332941111,0.008302262397641869,0.014548238897396632 +uniform_factorial,1,163,0.65,0.5,128.0,0.05,0.15,0.4286607049870561,0.4307829160924542,0.007751937984496124,0.0033229512014500266,0.0020964132981635686,0.013392857142857144 +uniform_factorial,1,164,0.65,0.5,128.0,0.5,0.01,0.4286607049870561,0.4307829160924542,0.007751937984496124,0.0033229512014500266,0.0020964132981635686,0.007710100231303007 +uniform_factorial,1,170,0.65,0.9,2.0,0.5,0.01,0.4721758147131214,0.4307829160924542,0.3333333333333333,0.15739193823770714,0.10803846135795675,0.007710100231303007 +uniform_factorial,1,175,0.65,0.9,8.0,0.05,0.15,0.4721758147131214,0.4307829160924542,0.1111111111111111,0.052463979412569015,0.03188325982637452,0.013392857142857144 +uniform_factorial,1,176,0.65,0.9,8.0,0.5,0.01,0.4721758147131214,0.4307829160924542,0.1111111111111111,0.052463979412569015,0.03188325982637452,0.007710100231303007 +uniform_factorial,1,177,0.65,0.9,8.0,0.5,0.15,0.4721758147131214,0.4307829160924542,0.1111111111111111,0.052463979412569015,0.03188325982637452,0.11952191235059761 +uniform_factorial,1,178,0.65,0.9,8.0,0.95,0.01,0.4721758147131214,0.4307829160924542,0.1111111111111111,0.052463979412569015,0.03188325982637452,0.014548238897396632 +uniform_factorial,1,179,0.65,0.9,8.0,0.95,0.15,0.4721758147131214,0.4307829160924542,0.1111111111111111,0.052463979412569015,0.03188325982637452,0.20503597122302158 +uniform_factorial,1,180,0.65,0.9,32.0,0.05,0.01,0.4721758147131214,0.4307829160924542,0.030303030303030304,0.014308358021609735,0.008302262397641869,0.0007763975155279505 +uniform_factorial,1,181,0.65,0.9,32.0,0.05,0.15,0.4721758147131214,0.4307829160924542,0.030303030303030304,0.014308358021609735,0.008302262397641869,0.013392857142857144 +uniform_factorial,1,186,0.65,0.9,128.0,0.05,0.01,0.4721758147131214,0.4307829160924542,0.007751937984496124,0.0036602776334350487,0.0020964132981635686,0.0007763975155279505 +uniform_factorial,1,195,0.9,0.5,2.0,0.5,0.15,0.12247448713915887,0.10536051565782628,0.3333333333333333,0.0408248290463862,0.023123169408776692,0.08928571428571429 +uniform_factorial,1,197,0.9,0.5,2.0,0.95,0.15,0.12247448713915887,0.10536051565782628,0.3333333333333333,0.0408248290463862,0.023123169408776692,0.15702479338842976 +uniform_factorial,1,200,0.9,0.5,8.0,0.5,0.01,0.12247448713915887,0.10536051565782628,0.1111111111111111,0.013608276348795313,0.006640070861282682,0.005580357142857143 +uniform_factorial,1,202,0.9,0.5,8.0,0.95,0.01,0.12247448713915887,0.10536051565782628,0.1111111111111111,0.013608276348795313,0.006640070861282682,0.010549694614103275 +uniform_factorial,1,204,0.9,0.5,32.0,0.05,0.01,0.12247448713915887,0.10536051565782628,0.030303030303030304,0.0037113480951260457,0.0017170260888290656,0.0005608524957936063 +uniform_factorial,1,205,0.9,0.5,32.0,0.05,0.15,0.12247448713915887,0.10536051565782628,0.030303030303030304,0.0037113480951260457,0.0017170260888290656,0.00970873786407767 +uniform_factorial,1,209,0.9,0.5,32.0,0.95,0.15,0.12247448713915887,0.10536051565782628,0.030303030303030304,0.0037113480951260457,0.0017170260888290656,0.15702479338842976 +uniform_factorial,1,211,0.9,0.5,128.0,0.05,0.15,0.12247448713915887,0.10536051565782628,0.007751937984496124,0.000949414628985708,0.00043283472005928303,0.00970873786407767 +uniform_factorial,1,217,0.9,0.9,2.0,0.05,0.15,0.1349073756323204,0.10536051565782628,0.3333333333333333,0.04496912521077338,0.023123169408776692,0.00970873786407767 +uniform_factorial,1,219,0.9,0.9,2.0,0.5,0.15,0.1349073756323204,0.10536051565782628,0.3333333333333333,0.04496912521077338,0.023123169408776692,0.08928571428571429 +uniform_factorial,1,227,0.9,0.9,8.0,0.95,0.15,0.1349073756323204,0.10536051565782628,0.1111111111111111,0.014989708403591043,0.006640070861284458,0.15702479338842976 +uniform_factorial,1,230,0.9,0.9,32.0,0.5,0.01,0.1349073756323204,0.10536051565782628,0.030303030303030304,0.004088102291888515,0.0017170260888290656,0.005580357142857143 +uniform_factorial,1,231,0.9,0.9,32.0,0.5,0.15,0.1349073756323204,0.10536051565782628,0.030303030303030304,0.004088102291888515,0.0017170260888290656,0.08928571428571429 +uniform_factorial,1,233,0.9,0.9,32.0,0.95,0.15,0.1349073756323204,0.10536051565782628,0.030303030303030304,0.004088102291888515,0.0017170260888290656,0.15702479338842976 +uniform_factorial,1,234,0.9,0.9,128.0,0.05,0.01,0.1349073756323204,0.10536051565782628,0.007751937984496124,0.0010457936095528526,0.00043283472005928303,0.0005608524957936063 +uniform_factorial,1,237,0.9,0.9,128.0,0.5,0.15,0.1349073756323204,0.10536051565782628,0.007751937984496124,0.0010457936095528526,0.00043283472005928303,0.08928571428571429 +maximin_optimized,3,0,0.05,0.5,2.0,0.05,0.01,1.1635076278220096,2.995732273553991,0.3333333333333333,0.3878358759406699,1.6491929980444475,0.01 +maximin_optimized,2,5,0.05,0.5,2.0,0.95,0.15,1.1635076278220096,2.995732273553991,0.3333333333333333,0.3878358759406699,1.6491929980444475,0.7702702702702702 +maximin_optimized,1,11,0.05,0.5,8.0,0.95,0.15,1.1635076278220096,2.995732273553991,0.1111111111111111,0.1292786253135566,0.793706251013266,0.7702702702702702 +maximin_optimized,3,18,0.05,0.5,128.0,0.05,0.01,1.1635076278220096,2.995732273553991,0.007751937984496124,0.009019438975364405,0.07219424175040157,0.01 +maximin_optimized,3,23,0.05,0.5,128.0,0.95,0.15,1.1635076278220096,2.995732273553991,0.007751937984496124,0.009019438975364405,0.07219424175040157,0.7702702702702702 +maximin_optimized,8,120,0.35,0.9,2.0,0.05,0.01,0.8768979416100827,1.0498221244986778,0.3333333333333333,0.29229931387002756,0.3355856642737046,0.0014409221902017292 +maximin_optimized,8,125,0.35,0.9,2.0,0.95,0.15,0.8768979416100827,1.0498221244986778,0.3333333333333333,0.29229931387002756,0.3355856642737046,0.32386363636363635 +maximin_optimized,1,131,0.35,0.9,8.0,0.95,0.15,0.8768979416100827,1.0498221244986778,0.1111111111111111,0.09743310462334258,0.10687049253534275,0.32386363636363635 +maximin_optimized,8,138,0.35,0.9,128.0,0.05,0.01,0.8768979416100827,1.0498221244986778,0.007751937984496124,0.0067976584620936574,0.007218032071614644,0.0014409221902017292 +maximin_optimized,7,143,0.35,0.9,128.0,0.95,0.15,0.8768979416100827,1.0498221244986778,0.007751937984496124,0.0067976584620936574,0.007218032071614644,0.32386363636363635 +maximin_optimized,3,192,0.9,0.5,2.0,0.05,0.01,0.12247448713915887,0.10536051565782628,0.3333333333333333,0.0408248290463862,0.023123169408776692,0.0005608524957936063 +maximin_optimized,4,197,0.9,0.5,2.0,0.95,0.15,0.12247448713915887,0.10536051565782628,0.3333333333333333,0.0408248290463862,0.023123169408776692,0.15702479338842976 +maximin_optimized,4,210,0.9,0.5,128.0,0.05,0.01,0.12247448713915887,0.10536051565782628,0.007751937984496124,0.000949414628985708,0.00043283472005928303,0.0005608524957936063 +maximin_optimized,4,215,0.9,0.5,128.0,0.95,0.15,0.12247448713915887,0.10536051565782628,0.007751937984496124,0.000949414628985708,0.00043283472005928303,0.15702479338842976 +maximin_optimized,1,219,0.9,0.9,2.0,0.5,0.15,0.1349073756323204,0.10536051565782628,0.3333333333333333,0.04496912521077338,0.023123169408776692,0.08928571428571429 diff --git a/results/certified-maximin-design/paired-recovery/paired_recovery.csv b/results/certified-maximin-design/paired-recovery/paired_recovery.csv new file mode 100644 index 0000000..a683515 --- /dev/null +++ b/results/certified-maximin-design/paired-recovery/paired_recovery.csv @@ -0,0 +1,61 @@ +seed,design,generator,replicate_count,correct_count,recovery_rate,median_log_evidence_margin,minimum_log_evidence_margin +7,heuristic_maximin,innovation_l2,200,168,0.84,2.370863092071815,0.015989971664971847 +7,heuristic_maximin,surprise,200,166,0.83,2.231980360659616,0.057823111828959384 +7,heuristic_maximin,gain,200,184,0.92,4.481044949432173,4.146743496136196e-05 +7,heuristic_maximin,update_l2,200,183,0.915,3.6377256267286757,0.07554977876299063 +7,heuristic_maximin,information_gain,200,178,0.89,3.7894046971621016,0.002462782765924487 +7,heuristic_maximin,change_probability,200,195,0.975,8.48263787995816,0.014963841319399762 +7,certified_integer_maximin,innovation_l2,200,170,0.85,2.544359373404669,0.0012180282734988168 +7,certified_integer_maximin,surprise,200,153,0.765,1.763398782921188,0.0010887618280790434 +7,certified_integer_maximin,gain,200,163,0.815,2.9126658322986447,0.029952119646040387 +7,certified_integer_maximin,update_l2,200,151,0.755,1.6357895610814186,0.0005664630917436853 +7,certified_integer_maximin,information_gain,200,146,0.73,1.8053502038969018,0.01787382673903437 +7,certified_integer_maximin,change_probability,200,191,0.955,5.063948113851465,0.15624436037130351 +11,heuristic_maximin,innovation_l2,200,166,0.83,2.3024986642055367,0.004581505820588916 +11,heuristic_maximin,surprise,200,155,0.775,2.1856774921034177,0.018913756171507146 +11,heuristic_maximin,gain,200,189,0.945,4.454267869790888,0.02624997808282359 +11,heuristic_maximin,update_l2,200,186,0.93,3.7035164034001653,0.015589964133553735 +11,heuristic_maximin,information_gain,200,170,0.85,4.022788756418635,0.12414149665198693 +11,heuristic_maximin,change_probability,200,195,0.975,7.392091023319658,0.047698147319572115 +11,certified_integer_maximin,innovation_l2,200,169,0.845,2.413170996967807,0.004635748437902265 +11,certified_integer_maximin,surprise,200,143,0.715,1.7625061433579798,0.02673033480230913 +11,certified_integer_maximin,gain,200,171,0.855,2.714811081653588,0.02662005762579156 +11,certified_integer_maximin,update_l2,200,153,0.765,1.5946771740664474,0.031209347540261945 +11,certified_integer_maximin,information_gain,200,143,0.715,2.0066179951054384,0.0008556368456247299 +11,certified_integer_maximin,change_probability,200,184,0.92,4.880238456331114,0.14273230654612945 +19,heuristic_maximin,innovation_l2,200,168,0.84,2.2663990487777426,0.01718431694705913 +19,heuristic_maximin,surprise,200,150,0.75,2.1654303717821985,0.024464177346505522 +19,heuristic_maximin,gain,200,184,0.92,4.818207899063886,0.011962138771430375 +19,heuristic_maximin,update_l2,200,167,0.835,3.3799030416968012,0.016508853014620684 +19,heuristic_maximin,information_gain,200,184,0.92,4.625042042377174,0.03273894211902473 +19,heuristic_maximin,change_probability,200,191,0.955,8.18246635934211,0.00313147462063057 +19,certified_integer_maximin,innovation_l2,200,166,0.83,2.4163963897329737,0.0008561319843636284 +19,certified_integer_maximin,surprise,200,151,0.755,1.80185092580707,0.026437310936017866 +19,certified_integer_maximin,gain,200,172,0.86,3.1242290725561155,0.0026623369757743376 +19,certified_integer_maximin,update_l2,200,133,0.665,1.6571953314403594,0.03789957994277415 +19,certified_integer_maximin,information_gain,200,146,0.73,2.030181611815136,0.005589918690425577 +19,certified_integer_maximin,change_probability,200,182,0.91,5.668587772335426,0.025401509147279455 +23,heuristic_maximin,innovation_l2,200,164,0.82,2.5163316176315895,0.002974129428505279 +23,heuristic_maximin,surprise,200,159,0.795,2.205472551459371,0.005146518921758769 +23,heuristic_maximin,gain,200,173,0.865,4.739307560921329,0.006971053873847666 +23,heuristic_maximin,update_l2,200,183,0.915,3.4667460840671858,0.021056531448742533 +23,heuristic_maximin,information_gain,200,175,0.875,4.258551239967007,0.0012091898314423588 +23,heuristic_maximin,change_probability,200,192,0.96,7.720564828864887,0.018597892423983353 +23,certified_integer_maximin,innovation_l2,200,159,0.795,2.782291634205338,0.0009879176030267445 +23,certified_integer_maximin,surprise,200,150,0.75,1.732900030949148,0.006299553614422138 +23,certified_integer_maximin,gain,200,173,0.865,2.5669721898872826,0.0032445805976735187 +23,certified_integer_maximin,update_l2,200,148,0.74,1.752343444429858,0.03837867263441197 +23,certified_integer_maximin,information_gain,200,148,0.74,2.1045538468034444,0.012118246894512907 +23,certified_integer_maximin,change_probability,200,185,0.925,4.8889904825337585,0.08761845525354772 +31,heuristic_maximin,innovation_l2,200,169,0.845,2.216623342313243,0.020811543625612217 +31,heuristic_maximin,surprise,200,157,0.785,2.2485447247125663,0.01372465957428659 +31,heuristic_maximin,gain,200,188,0.94,5.278517797038381,0.10002786879880432 +31,heuristic_maximin,update_l2,200,172,0.86,3.3189306226121946,0.029335133508567424 +31,heuristic_maximin,information_gain,200,170,0.85,3.5590486269567965,0.017401882997731377 +31,heuristic_maximin,change_probability,200,194,0.97,7.746689205126872,0.08641809500532815 +31,certified_integer_maximin,innovation_l2,200,167,0.835,2.449945971344926,0.060487250483234334 +31,certified_integer_maximin,surprise,200,153,0.765,1.7568475230602179,0.018640464802992085 +31,certified_integer_maximin,gain,200,171,0.855,2.922935267347187,0.049261356807939194 +31,certified_integer_maximin,update_l2,200,142,0.71,1.4796783686643238,0.01668079180099724 +31,certified_integer_maximin,information_gain,200,148,0.74,1.8146322100636105,0.01842186360967446 +31,certified_integer_maximin,change_probability,200,185,0.925,5.532824461795972,0.00958530454727935 diff --git a/results/certified-maximin-design/paired-recovery/seed_summary.csv b/results/certified-maximin-design/paired-recovery/seed_summary.csv new file mode 100644 index 0000000..b160f54 --- /dev/null +++ b/results/certified-maximin-design/paired-recovery/seed_summary.csv @@ -0,0 +1,6 @@ +seed,heuristic_minimum_recovery,certified_minimum_recovery,minimum_recovery_delta_certified_minus_heuristic,heuristic_mean_recovery,certified_mean_recovery,mean_recovery_delta_certified_minus_heuristic +7,0.83,0.73,-0.09999999999999998,0.8949999999999999,0.8116666666666665,-0.08333333333333337 +11,0.775,0.715,-0.06000000000000005,0.8841666666666667,0.8025000000000001,-0.08166666666666655 +19,0.75,0.665,-0.08499999999999996,0.87,0.7916666666666666,-0.07833333333333337 +23,0.795,0.74,-0.05500000000000005,0.8716666666666666,0.8025000000000001,-0.06916666666666649 +31,0.785,0.71,-0.07500000000000007,0.8749999999999999,0.805,-0.06999999999999984 diff --git a/results/certified-maximin-design/paired-recovery/summary.json b/results/certified-maximin-design/paired-recovery/summary.json new file mode 100644 index 0000000..1e52bba --- /dev/null +++ b/results/certified-maximin-design/paired-recovery/summary.json @@ -0,0 +1,59 @@ +{ + "all_seed_minima_certified_minus_heuristic_negative": true, + "allocation_l1_distance": 100, + "certified_support_count": 12, + "config": { + "effect_size": 1.0, + "noise_std": 1.0, + "recovery_seed_offset": 100003, + "replicates": 200, + "seeds": [ + 7, + 11, + 19, + 23, + 31 + ], + "test_fraction": 0.35 + }, + "experiment": "paired_heuristic_vs_certified_n60_recovery_diagnostic", + "heuristic_support_count": 15, + "input_provenance": [ + { + "budget": 60, + "constructor_reproduced": true, + "kind": "chronologically_locked_heuristic_maximin_n60", + "max_point_fraction": 0.15, + "source_bytes": 14450, + "source_code_sha": "1b2028929ac6ebc1cce0882f0c22af9918044342", + "source_file": "optimal_design_allocation_seed7.csv", + "source_repository": "IPS-Stuttgart/Bayesian-ACh", + "source_sha256": "a823be49faf6c6cbebf60b11d4b5ca895cf7734d6e9c577ee98f97a5907b69b2" + }, + { + "budget": 60, + "certificate_allocation_sha256": "694a84107c1ba94f39146a05675def1b5813621fab486d7567f6711510a6fc08", + "certificate_code_sha": "e82f48d335fb3f17f843d16bf2442693b6c83eef", + "certificate_sha256sums_sha256": "ae337d44c4b5b0ae3bbaaaf1a4981a29bffe7cbb5458104f05d19de13c313938", + "certificate_summary_sha256": "7c89b45552568665fd22103feaa301d9d69030bc1a5101b20c3ce6d438b3743c", + "certified": true, + "kind": "certified_integer_maximin_n60", + "lower_bound": 0.25286766545125633, + "source_repository": "IPS-Stuttgart/Bayesian-ACh", + "upper_bound": 0.2528676654512564 + } + ], + "interpretation": "The exact integer certificate optimizes the frozen asymptotic worst-residual objective. This paired finite-sample diagnostic does not select a replacement schedule, establish empirical superiority, or define a physical trial/animal protocol.", + "maximum_absolute_cell_change": 9, + "minimum_recovery_across_seeds": { + "certified_integer_maximin": 0.665, + "heuristic_maximin": 0.75 + }, + "paired_rng": true, + "producer_clean_worktree": true, + "producer_commit": "c6dea54d1e8a1b702c212c83535af8a5de3ab3ea", + "recovery_seed_rule": "seed + recovery_seed_offset; RNG reset identically for each allocation", + "schema_version": 1, + "support_overlap_count": 4, + "support_union_count": 23 +}