feat(analysis): PhILR, SBP and balance-dendrogram ILR bases with Agda proofs (#20)

.github/workflows/ci.yml

@@ -162,6 +162,40 @@ jobs:
           fi
           echo "commit convention ok: $count non-merge commit(s) graded in $range"
 
+  # Agda proof gate for the compositional transforms (issues #20, #21); scope,
+  # layering and residue are in docs/formal/verification-plan.md. The toolchain
+  # is the estate's pin (epistemic-types, residual-evidence-types): Debian 13
+  # packages agda-bin 2.6.4.3-1+b2 and agda-stdlib 2.1-4, installed from the
+  # signed Debian archive rather than a downloaded binary. The same script runs
+  # locally with AGDA=... AGDA_STDLIB_LIB=... overrides.
+  proofs:
+    name: Proofs (Agda)
+    runs-on: ubuntu-24.04
+    timeout-minutes: 30
+    container: debian:13-slim@sha256:d7e12182ce18b85b93007c1dedf31f2d29e01ccf3182cc4017c709b6259bc132
+    permissions:
+      contents: read
+    env:
+      LANG: C.UTF-8
+      LC_ALL: C.UTF-8
+    steps:
+      - name: Install the pinned Agda toolchain (Debian archive)
+        run: |
+          apt-get update
+          apt-get install --no-install-recommends -y ca-certificates git \
+            agda-bin=2.6.4.3-1+b2 agda-stdlib=2.1-4
+
+      - uses: actions/checkout@3d3c42e5aac5ba805825da76410c181273ba90b1  # v7.0.1
+        with:
+          persist-credentials: false
+
+      - name: Guard, type-check, negative controls
+        run: |
+          test "$(agda --version)" = 'Agda version 2.6.4.3'
+          AGDA_STDLIB_LIB=$(dpkg -L agda-stdlib | grep 'standard-library\.agda-lib$' | head -1)
+          export AGDA_STDLIB_LIB
+          bash scripts/check-proofs.sh
+
   test:
     # The name is a fixed string, it carries no version, and the job has NO matrix.
     # Both are required, and the first without the second is a trap that has already
@@ -616,6 +650,11 @@ jobs:
           # over the same shadowing. Running it means the next AnalysisConfig API
           # change breaks this step immediately instead of at test time.
           julia --project=. bench/analysis_config/benchmark.jl
+          # ILR bases (issue #20): CLR, default ILR and the three new bases at 100,
+          # 1 000 and 10 000 taxa. Emits ::warning:: for any workload over 5 minutes
+          # or over 1 GiB (allocation or peak-RSS growth); never fails the job — the
+          # hard CLR/ILR gate is the same-runner base-vs-head step below.
+          julia --project=. bench/ilr_bases/benchmark.jl
           julia --project=. bench/comprehensive_benchmark.jl
 
       - name: Summarise Julia benchmark deltas (informational)
@@ -632,6 +671,34 @@ jobs:
           if [ -f bench/results/comprehensive_results.json ]; then
             echo "Comprehensive results: $(cat bench/results/comprehensive_results.json | head -c 500)"
           fi
+          if [ -f bench/results/ilr_bases_results.json ]; then
+            echo "ILR bases results: $(head -c 500 bench/results/ilr_bases_results.json)"
+          fi
+
+      # Issue #20: "fail CI if CLR/ILR performance regresses more than 10%". A committed
+      # baseline cannot carry that (cross-host noise, as above), so this measures the PR's
+      # base commit and its head on this runner, interleaved base/head/base/head, and fails
+      # if head allocates >10% more or its minimum time is >10% slower. Rationale and the
+      # statistic are documented in bench/ilr_bases/regression_gate.jl. Pull requests only:
+      # a push to main has no base to compare against.
+      - name: CLR/ILR regression gate (same runner, base vs head)
+        if: github.event_name == 'pull_request'
+        env:
+          BASE_SHA: ${{ github.event.pull_request.base.sha }}
+          R_LIBS_SITE: /usr/local/lib/R/site-library:/usr/lib/R/site-library:/usr/lib/R/library
+        run: |
+          set -euo pipefail
+          git fetch --no-tags --depth=1 origin "$BASE_SHA"
+          git worktree add --detach "$RUNNER_TEMP/base" "$BASE_SHA"
+          julia --project="$RUNNER_TEMP/base" -e 'using Pkg; Pkg.instantiate()'
+          gate=bench/ilr_bases/regression_gate.jl
+          out="$RUNNER_TEMP/ilr_gate"
+          mkdir -p "$out"
+          for round in 1 2; do
+            julia --project="$RUNNER_TEMP/base" "$gate" measure base "$out/base_$round.json"
+            julia --project=. "$gate" measure head "$out/head_$round.json"
+          done
+          julia --project=. "$gate" compare --base "$out"/base_*.json --head "$out"/head_*.json
 
       - name: Upload Julia benchmark artifacts
         if: always()
@@ -641,6 +708,7 @@ jobs:
           path: |
             bench/*/baseline.json
             bench/results/comprehensive_results.json
+            bench/results/ilr_bases_results.json
             bench/**/baseline.json
           if-no-files-found: warn
 

.gitignore

@@ -266,6 +266,8 @@ docs/*
 !docs/owner-review-2026-09-25.md
 !docs/integration/
 !docs/integration/**
+!docs/formal/
+!docs/formal/**
 
 # Test + benchmark run artifacts (not the committed fixtures/baselines)
 frontend/tests/results/
@@ -358,3 +360,8 @@ archive/
 # Track the Stipple/Vue migration recon and implementation plan.
 !docs/migration/
 !docs/migration/**
+
+### Agda ###
+# Interface files (regenerated by scripts/check-proofs.sh)
+proofs/agda/_build/
+*.agdai

CHANGELOG.md

@@ -12,6 +12,39 @@ types, tests, infrastructure, and alignment.
 
 ## [Unreleased]
 
+### Added — the three deferred ILR bases are implemented, with proofs (#20, 2026-09-26)
+
+- **`phylogenetic` (PhILR), `sequential_binary_partition` and `balance_dendrogram`** are
+  computed, no longer refused (`DEFERRED_ILR_BASIS` is now empty; it stays, so a future
+  basis can be deferred the same way). `src/analysis/ilr_basis.jl` is one engine for all
+  three: each basis is a rooted binary tree, and balances are clade sums in one post-order
+  pass, `O(D)` per sample, with no dense basis matrix. Pure Julia; no new dependency.
+- **Inputs and contracts.** `advanced.ilr_phylo_tree_path`, `ilr_sbp_matrix_path`,
+  `ilr_balance_dendrogram_method` (plus philr's part/balance weights and the SBP history)
+  are identical in the Julia validator, the Nickel contract, the JSON schema, DEED and the
+  frontend: each basis requires its input and refuses the others'. Configurations that
+  use none of them hash exactly as before; the default Helmert balances are byte-identical.
+- **Validity is refused, not repaired.** Trees must be rooted, bifurcating and cover the
+  retained taxa (tips outside them are pruned and recorded); an SBP must be the SBP of a
+  binary tree (Egozcue & Pawlowsky-Glahn 2005), over exactly the retained taxa.
+- **Provenance and guards.** Each run records the tree or SBP SHA-256, the dendrogram
+  method, the weights and the balance-id rule. More than 3 distinct SBPs tried in a project
+  raises a DANGER (p-hacking guard). BH stays mandatory.
+- **Evidence.** Known answers; an independent Julia reference
+  (`test/fixtures/ilr/ilr_reference.jl`) that recomputes all 17 committed expectations
+  (philr, `compositions::ilr`, R `hclust`) every run; R cross-checks where philr,
+  compositions and robCompositions are installed; negative controls. The balance algebra
+  (contrast sums, orthonormality, injectivity with its positive-weight hypothesis, scale
+  and perturbation invariance, comb = Helmert, SBP validity) is proved in Agda
+  (`proofs/agda/`, new `proofs` CI job) and mapped to the tests in
+  `docs/formal/verification-plan.md`. The conditions document is
+  `docs/statistics/method-conditions/ilr-bases.md`.
+- **Benchmarks.** `bench/ilr_bases/benchmark.jl` (100 / 1 000 / 10 000 taxa; warns over
+  5 minutes or 1 GiB) and a pull-request gate that fails when CLR or default ILR allocate
+  or run more than 10 % worse than the base commit on the same runner.
+- **Removed** the Python fixture generator: Python is not permitted by the estate language
+  policy; the Julia reference replaces it (`docs/compliance/standards-alignment.md`).
+
 ### Fixed — the NB test fixture is data a negative binomial describes (2026-09-26)
 
 - The estimation tests' synthetic table was **under-dispersed** (variance below the mean,

Justfile

@@ -274,6 +274,11 @@ lint:
 hygiene: spdx format lint
     @echo "hygiene: OK"
 
+# Agda proof gate (docs/formal/verification-plan.md): guard, type-check,
+# negative controls. Honours AGDA=... and AGDA_STDLIB_LIB=... overrides.
+proofs:
+    ./scripts/check-proofs.sh
+
 # Lint a commit message against the canonical format (default: HEAD).
 commit-check msg="":
     #!/usr/bin/env bash
@@ -364,6 +369,25 @@ bench-julia data="":
     fi
     $JULIA_CMD --project=. -t4 bench/layer1_mock_recovery/runner.jl {{data}}
 
+# ILR-basis scaling benchmark (issue #20); taxa="100,1000" for a quick run.
+bench-ilr taxa="100,1000,10000":
+    ILR_BENCH_TAXA={{taxa}} $JULIA_CMD --project=. bench/ilr_bases/benchmark.jl
+
+# The CI CLR/ILR regression gate, locally: this checkout vs `base` (a git ref),
+# interleaved base/head/base/head on this machine; fails on >10% (time or allocation).
+bench-ilr-gate base="origin/main":
+    #!/usr/bin/env bash
+    set -euo pipefail
+    tmp=$(mktemp -d)
+    trap 'git worktree remove --force "$tmp/base" >/dev/null 2>&1 || true; rm -rf "$tmp"' EXIT
+    git worktree add --detach "$tmp/base" "{{base}}"
+    $JULIA_CMD --project="$tmp/base" -e 'using Pkg; Pkg.instantiate()'
+    for round in 1 2; do
+        $JULIA_CMD --project="$tmp/base" bench/ilr_bases/regression_gate.jl measure base "$tmp/base_$round.json"
+        $JULIA_CMD --project=. bench/ilr_bases/regression_gate.jl measure head "$tmp/head_$round.json"
+    done
+    $JULIA_CMD --project=. bench/ilr_bases/regression_gate.jl compare --base "$tmp"/base_*.json --head "$tmp"/head_*.json
+
 # ----------------------------------------------------------------------- #
 # Composites
 # ----------------------------------------------------------------------- #

README.md

@@ -553,6 +553,7 @@ clean checkout (`frontend/`):
 | Licence headers | `scripts/check-spdx.sh` | gated |
 | Formatting | `scripts/check-format.sh` | gated |
 | Lint (tsc semantics + shell) | `scripts/check-lint.sh` | gated |
+| Agda proofs (ILR bases, compositions) | `scripts/check-proofs.sh` | gated (`proofs` job) |
 
 CI runs the same gates (see `.github/workflows/ci.yml`: repo-hygiene job,
 then the pinned Julia and frontend jobs). Contributor setup, commit and

bench/ilr_bases/benchmark.jl

@@ -0,0 +1,141 @@
+# SPDX-License-Identifier: AGPL-3.0-only
+# SPDX-FileCopyrightText: 2026 Jonathan D.A. Jewell (hyperpolymath) <j.d.a.jewell@open.ac.uk>
+"""
+ILR-basis scaling benchmark (issue #20) — 100, 1 000 and 10 000 taxa.
+
+Measures, per size, the CLR transform and the default (Helmert) ILR through
+`Execution.prepare_analysis_table`, and each new basis through the engine
+(`ILRBasis.ilr_transform`, file read and SHA-256 included):
+
+  * phylogenetic — a balanced rooted bifurcating tree;
+  * sequential_binary_partition — the SBP CSV of the same tree (D x (D-1) cells: the
+    format the issue specifies is dense, 10^8 cells at 10 000 taxa; the reader streams it);
+  * balance_dendrogram — ward, complete and average on the variation matrix.
+
+Reported per workload: wall time, bytes allocated, and growth of the process's peak
+resident set. A workload over 5 minutes, or over 1 GiB of allocation or of peak-RSS growth,
+is WARNED about (a GitHub `::warning::` annotation in CI), as the issue asks. Like every
+other Julia benchmark in this repository the result is informational (cross-host timings
+are noise; see the comment above the frontend benchmark step in ci.yml): the hard >10%
+CLR/ILR regression gate is the same-runner base-vs-head comparison in
+bench/ilr_bases/regression_gate.jl.
+
+Run:  julia --project=. bench/ilr_bases/benchmark.jl
+      ILR_BENCH_TAXA=100,1000 julia --project=. bench/ilr_bases/benchmark.jl
+Writes bench/results/ilr_bases_results.json.
+"""
+
+using MetaManifold
+using MetaManifold: AnalysisConfig, Execution
+using JSON3
+using Logging
+using OrderedCollections
+
+const ILR = MetaManifold.ILRBasis
+const WARN_SECONDS = 300.0
+const WARN_BYTES = 1 << 30
+const SAMPLES = 20
+
+sizes() = parse.(Int, split(get(ENV, "ILR_BENCH_TAXA", "100,1000,10000"), ','))
+
+# Deterministic strictly positive taxa x samples table (no RNG: identical on every host).
+synth(D, n) = [1.5 + mod(i * 7919 + j * 104729, 997) + 0.25 * mod(i * j, 7) for i in 1:D, j in 1:n]
+
+function balanced_newick(names::Vector{String}, lo::Int, hi::Int)
+    lo == hi && return names[lo]
+    mid = (lo + hi) >>> 1
+    return "(" * balanced_newick(names, lo, mid) * ":0.1," * balanced_newick(names, mid + 1, hi) * ":0.1)"
+end
+
+function write_sbp(path, taxa, W)
+    open(path, "w") do io
+        println(io, "taxon,", join(("b$k" for k in 1:size(W, 2)), ","))
+        for i in eachindex(taxa)
+            print(io, taxa[i])
+            for k in 1:size(W, 2)
+                print(io, ',', W[i, k])
+            end
+            println(io)
+        end
+    end
+end
+
+function config(norm::String, method::String)
+    AnalysisConfig.AnalysisConfig(
+        method = method, formula = "~ group", metadata_columns = ["group"],
+        normalization = AnalysisConfig.NormalizationConfig(method = norm, pseudocount = 0.5),
+        advanced = AnalysisConfig.AdvancedConfig(min_prevalence = 0.0, min_abundance = 0.0),
+        created_by = "bench/ilr_bases")
+end
+
+function measure(f)
+    GC.gc()
+    rss0 = Sys.maxrss()
+    stats = @timed f()
+    return (seconds = stats.time, allocated_bytes = stats.bytes, gc_seconds = stats.gctime,
+            peak_rss_growth_bytes = max(0, Int(Sys.maxrss()) - Int(rss0)))
+end
+
+function workloads(D::Int, dir::String)
+    X = synth(D, SAMPLES)
+    taxa = ["t$i" for i in 1:D]
+    samples = ["s$j" for j in 1:SAMPLES]
+    tree = joinpath(dir, "tree_$D.nwk")
+    write(tree, balanced_newick(taxa, 1, D) * ";")
+    sbp = joinpath(dir, "sbp_$D.csv")
+    write_sbp(sbp, taxa, ILR.sbp_matrix(first(ILR.phylo_balance_tree(ILR.parse_newick(read(tree, String)), taxa))))
+    prep(cfg) = () -> Execution.prepare_analysis_table(cfg, X; sample_ids = samples, taxa_ids = taxa, drop_policy = "drop")
+    clr, ilr = config("clr", "clr_lm"), config("ilr", "ilr_lm")
+    w = OrderedDict{String,Any}(
+        "clr (prepare_analysis_table)" => prep(clr),
+        "ilr default Helmert (prepare_analysis_table)" => prep(ilr),
+        "phylogenetic" => () -> ILR.ilr_transform(X, taxa; basis = "phylogenetic", tree_path = tree),
+        "sequential_binary_partition" => () -> ILR.ilr_transform(X, taxa; basis = "sequential_binary_partition", sbp_path = sbp),
+    )
+    for m in ILR.VALID_DENDROGRAM_METHODS
+        w["balance_dendrogram $m"] = () -> ILR.ilr_transform(X, taxa; basis = "balance_dendrogram", dendrogram_method = m)
+    end
+    return w
+end
+
+function main()
+    in_ci = haskey(ENV, "GITHUB_ACTIONS")
+    results = OrderedDict{String,Any}("julia" => string(VERSION), "samples" => SAMPLES,
+                                      "warn_seconds" => WARN_SECONDS, "warn_bytes" => WARN_BYTES,
+                                      "runs" => Any[])
+    warnings = String[]
+    mktempdir() do dir
+        with_logger(NullLogger()) do
+            # compile everything once on a small problem so no size pays for compilation
+            foreach(f -> f(), values(workloads(20, dir)))
+            for D in sizes()
+                for (name, f) in workloads(D, dir)
+                    m = measure(f)
+                    push!(results["runs"], OrderedDict{String,Any}("taxa" => D, "workload" => name, pairs(m)...))
+                    flags = String[]
+                    m.seconds > WARN_SECONDS && push!(flags, "took $(round(m.seconds; digits = 1)) s (> 5 min)")
+                    m.allocated_bytes > WARN_BYTES && push!(flags, "allocated $(round(m.allocated_bytes / 2^30; digits = 2)) GiB (> 1 GiB)")
+                    m.peak_rss_growth_bytes > WARN_BYTES && push!(flags, "grew peak RSS by $(round(m.peak_rss_growth_bytes / 2^30; digits = 2)) GiB (> 1 GiB)")
+                    isempty(flags) || push!(warnings, "ILR benchmark, $D taxa, $name: " * join(flags, "; "))
+                end
+            end
+        end
+    end
+    println("=== ILR bases: scaling (", SAMPLES, " samples) ===")
+    println(rpad("taxa", 7), rpad("workload", 48), lpad("time s", 10), lpad("alloc MiB", 12), lpad("ΔpeakRSS MiB", 14))
+    for r in results["runs"]
+        println(rpad(string(r["taxa"]), 7), rpad(r["workload"], 48), lpad(string(round(r["seconds"]; digits = 3)), 10),
+                lpad(string(round(r["allocated_bytes"] / 2^20; digits = 1)), 12),
+                lpad(string(round(r["peak_rss_growth_bytes"] / 2^20; digits = 1)), 14))
+    end
+    results["warnings"] = warnings
+    for w in warnings
+        println(in_ci ? "::warning::$w" : "WARNING: $w")
+    end
+    out = joinpath(@__DIR__, "..", "results", "ilr_bases_results.json")
+    mkpath(dirname(out))
+    open(io -> JSON3.pretty(io, results), out, "w")
+    println("results: ", normpath(out))
+end
+
+main()