Profile an observable's additive offset — and the coupled (scale, offset) pair — out of the search, completing the ADR-0066/ADR-0108 family (6 of 23 benchmark slugs, up to 42% of the search): method proposal + evaluation requested before we commit

[wshlavacek]

ADR-0066/ADR-0099 profile a declared column's optimal multiplicative scale out analytically, and ADR-0108 (#562) profiles out an estimated noise scale. The third member of that family is missing: an additive offset, and — the case that is not just "one more scalar" — the coupled (scale, offset) pair, which cannot be profiled one parameter at a time because the two are jointly determined.

This issue proposes a method and asks for an evaluation before we commit to it. See "The decision this issue is actually asking for" at the bottom; there is a real counter-hypothesis, and I do not think the answer is obvious.

What is in the box today

A PEtab observableParameters scale/offset declared as a fit parameter is an ordinary free parameter searched in the box. MeasurementModel.prediction_sensitivity (measurement/base.py:147) resolves it as "a directly-named free parameter ... a scale/offset estimated as a fit parameter", contributing its ∂f/∂param column and nothing more. Nothing profiles it.

Scope across the Grein et al. 2026 subset-I corpus, counting only observable-layer parameters that enter the prediction linearly:

slug	linear observable params in the box	form	k	share
Schwen_PONE2015	scale, offset, offset_nExpID1–4, scaleElisa_nExpID1–4	joint, nested + per-condition	24	42%
Weber_BMC2015	scale_yCERTpRN24, scale_yPI4K3BpRN24, scale_yPKDpN0/24/25	pure scale	36	14%
Elowitz_Nature2000	scale, background	joint	21	10%
Borghans_BiophysChem1997	scale, offset	joint	23	9%
Brannmark_JBC2010	k_IRSiP_1Step, k_IRSiP_2Step	pure scale, per-condition	22	9%
Laske_PLOSComputBiol2019	vRNA_offset	pure offset	13	8%

6 of 23 slugs, 22 parameters. Three of the six carry a genuine coupled pair.

Two details worth pulling out of that table:

• Brannmark's two look like rate constants by name, but neither appears in model_Brannmark_JBC2010.xml — they are pure observable scales bound per condition through the measparams table.
• Laske's vRNA_offset is searched over loguniform 1e-10 1e10. Twenty decades of box for a quantity with a closed form.

And one slug must be refused. Fiedler_BMCSystBiol2016 binds the same eight tokens (s_pErk_*, s_pMek_*) to both observableParameter1_* and noiseParameter1_* — each parameter is simultaneously the per-gel scale and the per-gel σ. It is not a free linear coefficient, and profiling it as one would silently change the fit. Any implementation needs to detect this double binding and refuse it by name, the way noise_profiling already refuses formula/data-column noise sources.

The two slugs this lands on are the two slugs #562 already singled out

That is the part I did not expect. ADR-0108's motivating evidence was Borghans (converges to the no-dynamics attractor) and Schwen (σ runs into its box bound). Both of those also carry an unprofiled (scale, offset) pair — Borghans is literally observable: Ca, formula: Z_state*scale + offset, and Schwen carries ten linear observable parameters, 42% of its search.

So on Borghans the proposal here would take the search from 23 dimensions to 20 (σ via #562, plus scale and offset), and on Schwen from 24 to 12.

Proposed method

For one series with prediction y_i(θ) = a·s_i(θ) + b over raw simulated column s, data d, and weights W = diag(1/σ_i²), stack

Φ(θ) = [ s(θ)  1 ] ∈ R^{n×2}
c*(θ) = (Φᵀ W Φ)⁻¹ Φᵀ W d          c = (a, b)
r(θ)  = W^{1/2} ( d − Φ(θ) c*(θ) ) = W^{1/2} (I − P_Φ) d

which is exactly variable projection (Golub & Pereyra 1973): score the data against its residual after orthogonal projection onto span(Φ). Pure-offset is the same construction with Φ = [1], pure-scale with Φ = [s] — so the existing ADR-0066 scale falls out as the 1-column special case rather than being a parallel code path.

The Golub–Pereyra derivative is

∂r/∂θ = − W^{1/2} [ (I − P)(∂Φ/∂θ)Φ⁺ + ((I − P)(∂Φ/∂θ)Φ⁺)ᵀ ] d

with Kaufman's variant dropping the second term (cheaper per iteration, near-identical iteration counts in practice).

This needs no new sensitivity machinery. ∂Φ/∂θ has exactly one nonzero column — ∂s/∂θ, the forward sensitivity already assembled — because the intercept column is constant in θ. The added cost is one 2×2 weighted solve per series per evaluation.

The weighting question, which is the interesting one

ADR-0099 had to carry the full ∂c*/∂θ product rule, and explained why: PyBNF's existing profiling criterion is σ-unweighted and location-agnostic, so it is not the objective's own conditional minimizer over c, so ∂Obj/∂c* ≠ 0 and the envelope theorem does not apply.

If the new joint solve is W-weighted — i.e. it genuinely is the conditional minimizer of the Gaussian objective over (a, b) — then the envelope theorem does apply, ∂Obj/∂c* = 0, and the reduced gradient is just the partial. That is strictly simpler than what ADR-0099 already ships, not harder.

The cost is that it makes the new joint profiling inconsistent with the existing normalization = ..., scale, which stays σ-unweighted. Options: (a) weighted joint solve, accept two different criteria in the codebase; (b) unweighted, stay consistent with ADR-0066 and pay the product rule; (c) weight both and migrate ADR-0066, which changes existing scale results. This is a decision, not a detail, and it should be made deliberately in an ADR rather than fallen into.

Other implementation notes

• Reparametrization. Schwen declares scale*(IR1 + IR1in + offset) — not linear in the declared (scale, offset), but spanning the same 2-column space, since scale*x + scale*offset is linear in (scale, scale*offset). VarPro solves for span coefficients; mapping back to the user's declared names for reporting is a separate step, and is ill-posed when scale → 0. Needs a stated convention.
• Partial separability is real and should be supported. Schwen's observable_Insulin is observableParameter1 + observableParameter2 * g(·; observableParameter3, observableParameter4) — two linear parameters and two nonlinear ones in one observable. Classifying which declared parameters are linear (rather than assuming the whole observable is) is the part that makes this a general facility instead of a special case.
• Ordering against noise_profiling. Profiled σ and profiled linear parameters are nested — σ's MLE depends on residuals that depend on the linear solve. Linear-first, then σ in closed form (Loos et al. 2018's order). Pin it explicitly.
• information_criteria.txt must keep counting profiled linear parameters in k, exactly as Profile an estimated noise scale out analytically, as ADR-0066 already does for a multiplicative scale (13 of 23 benchmark slugs, up to 14% of the search) #562 required for σ. They remain estimated quantities; only the search drops them.
• Bayesian refusal carries over. config.py:2766 already refuses noise_profiling for Bayesian samplers because profiling collapses a posterior dimension. Same argument, same refusal.
• profile_likelihood interaction. A profiled-out parameter has no profile curve. Decide whether it is reportable (the closed form makes a conditional CI cheap) or simply absent from the summary.
• Rank deficiency. Φ = [s, 1] is singular when s is constant over the series — which happens at exactly the degenerate θ a global sampler will visit. Needs the _warn_degenerate_profile treatment (objective.py:1329) and a pseudo-inverse fallback, not a raw solve.

The decision this issue is actually asking for

I want this evaluated before it is built, because there is a counter-hypothesis that cuts the other way and I cannot tell from first principles which wins.

The case for. Fewer search dimensions (up to 42% on Schwen), better conditioning, every draw a global sampler ranks is linear-optimal, no box artifacts on quantities that have a closed form, and no more 20-decade boxes on an offset.

The case against, and it is not weak. #562's argument was that σ is the cheapest descent direction in the box, so the sampler ranks draws by how wrong their σ is. An offset is a second cheap descent direction of exactly that kind — and worse, Borghans's known attractor (a flat line at the best constant, scoring -51.204092) is reachable purely through the linear parameters: it is a = 0, b = mean(d). Profiling makes that solution available for free at every θ. By the envelope theorem the joint optimum does not move, but the reduced surface is now everywhere at-or-below the flat-line score, which could compress the ranking a global method sees rather than sharpening it. Profiling σ removed a distractor; profiling the offset might instead hand every candidate the same floor.

I genuinely do not know which of those dominates. So:

Requested evaluation, with kill criteria.

1. Offline, no fitting: on the six slugs above, compute the reduced (linear-profiled) objective and compare its landscape against the searched version along the box_probe.py sections. Does the reduced surface separate the reference optimum from the flat-line floor better or worse?
2. On Borghans specifically: does the -51.204092 attractor remain a converged endpoint when scale, offset, and σ are all profiled? This is the single most informative measurement.
3. Head-to-head fits on all six at matched budget, profiled vs searched, reporting best objective and iterations-to-converge.
4. Confirm the Fiedler double-binding refusal triggers and that no other slug in the wider corpus has the same pattern undetected.

Kill criterion. If (1) and (2) show the reduced landscape is no better — or worse — for global search, and (3) shows no reliable improvement, then close this as a tidiness item and do not implement. Fewer dimensions is not on its own a reason to ship a change to what the objective means.

Non-goals

General variable projection over the dynamics parameters. Rate constants enter the ODE right-hand side and a trajectory is a nonlinear functional of them; there is no Φ(α)c factorization to exploit. The separable structure in this domain lives entirely in the observation model, and this issue is scoped to it. Prior art: hierarchical optimization for ODE models (Loos et al. 2018), and pyPESTO/AMICI, which profile scale, offset, and σ analytically by default — the same family ADR-0066 and ADR-0108 already borrowed the other two thirds from.