|
| 1 | +# Search-agnostic prior-support enforcement: a Clipper class |
| 2 | + |
| 3 | +Type: feature |
| 4 | +Target: PyAutoFit |
| 5 | +Repos: |
| 6 | +- PyAutoFit |
| 7 | +Difficulty: medium |
| 8 | +Autonomy: supervised |
| 9 | +Priority: high |
| 10 | +Status: formalised |
| 11 | + |
| 12 | +## Why |
| 13 | + |
| 14 | +`@PyAutoFit/autofit/non_linear/search/mle/multi_start_gradient/search.py` builds |
| 15 | +its objective as |
| 16 | + |
| 17 | +``` |
| 18 | +fom = -2 * (log_likelihood + sum(log_prior_list)) |
| 19 | +``` |
| 20 | + |
| 21 | +A `UniformPrior` returns `log_prior = -inf` outside its box, and the search steps |
| 22 | +in **physical** parameter space with nothing constraining it to that box. A lane |
| 23 | +that oversteps a hard prior edge reads as non-finite, is marked dead, and with |
| 24 | +`resurrect=False` is never redrawn. |
| 25 | + |
| 26 | +Measured on the real `imaging/mge` profiling cell (16 starts x 150 steps, cloud |
| 27 | +CPU) — full investigation and evidence in autolens_profiling#128: |
| 28 | + |
| 29 | +| arm | value-NaN lane-steps | lanes dead | alive at end | |
| 30 | +|---|---:|---:|---:| |
| 31 | +| baseline | 1446 (60.25%) | 14/16 | 2 | |
| 32 | +| shear box widened to ±1 | 1422 (59.25%) | 15/16 | 1 | |
| 33 | +| prior term neutered (diagnostic) | 215 (8.96%) | 3/16 | 13 | |
| 34 | +| **clip to prior box (prototype)** | **425 (17.71%)** | **5/16** | **11** | |
| 35 | + |
| 36 | +**The likelihood never went non-finite** in ~7200 lane-steps. This is entirely a |
| 37 | +prior-support problem. |
| 38 | + |
| 39 | +The behaviour is worse than "frozen": the overshoot is tiny (median 3% of box |
| 40 | +width, min 0.21%), and because `log_prior = -inf` is *constant* outside the box |
| 41 | +its derivative is zero, so the total gradient is the finite **likelihood** |
| 42 | +gradient. `optax.apply_if_finite` therefore never fires and the dead lane **keeps |
| 43 | +stepping forever** — full likelihood-and-gradient cost every step, output |
| 44 | +discarded, wandering far (one parameter went `0.30 -> -1.76`). 0/16 lanes ever |
| 45 | +revive. |
| 46 | + |
| 47 | +## The exposure is not MultiStart-only |
| 48 | + |
| 49 | +This is why the fix should not live inside one search: |
| 50 | + |
| 51 | +- **`MultiStartGradient`** (`MultiStartAdam` / `MultiStartADABelief` / |
| 52 | + `MultiStartLion` / `MultiStartProdigy` all share one `_fit`) — measured above. |
| 53 | +- **`@PyAutoFit/autofit/non_linear/search/mle/bfgs/search.py`** — same |
| 54 | + `Fitness(fom_is_log_likelihood=False, resample_figure_of_merit=-np.inf, |
| 55 | + convert_to_chi_squared=True)`, steps in physical space, and calls |
| 56 | + `optimize.minimize(fun=..., x0=..., method=self.method, options=..., tol=...)` |
| 57 | + with **no `bounds=` argument**. `L-BFGS-B` supports box bounds natively; they |
| 58 | + are simply not passed. Being single-start, this presents as a failed fit rather |
| 59 | + than a population collapse, so it is easier to misattribute. |
| 60 | +- **NUTS** (`@PyAutoFit/autofit/non_linear/search/mcmc/blackjax/nuts/search.py`) |
| 61 | + also targets the log posterior from a physical `initial_position`. HMC entering |
| 62 | + a `-inf` region *diverges* rather than freezing. **Out of scope here** — different |
| 63 | + mechanism, needs its own investigation. See "Deliberately out of scope". |
| 64 | + |
| 65 | +Not exposed, and correctly so: the nested samplers already work in unit-cube |
| 66 | +coordinates, and the MCMC samplers reject `-inf` proposals so the walker stays |
| 67 | +put. **Rejection is the restoring mechanism that gradient methods lack.** |
| 68 | + |
| 69 | +## The design |
| 70 | + |
| 71 | +A `Clipper`, modelled on `@PyAutoFit/autofit/non_linear/initializer.py` — a |
| 72 | +pluggable, per-search strategy object with a config-resolved default. |
| 73 | + |
| 74 | +**One place the `Initializer` analogy does not carry.** `Initializer` has a single |
| 75 | +consumption pattern (`samples_from_model`). `Clipper` has **two structurally |
| 76 | +different consumers** and must serve both from one source of truth: |
| 77 | + |
| 78 | +- `MultiStartGradient` enforces the constraint itself, every step → wants an |
| 79 | + imperative `project(...)`. |
| 80 | +- `LBFGS` hands bounds to scipy and lets *scipy* enforce → wants a declarative |
| 81 | + `bounds`. |
| 82 | + |
| 83 | +Proposed contract: |
| 84 | + |
| 85 | +```python |
| 86 | +class AbstractClipper(ABC): |
| 87 | + @abstractmethod |
| 88 | + def bounds_from_model(self, model) -> tuple[np.ndarray, np.ndarray]: |
| 89 | + """(lower, upper) in PHYSICAL parameter order. Unbounded -> -inf/+inf.""" |
| 90 | + |
| 91 | + @abstractmethod |
| 92 | + def project(self, vector, model, xp=np): |
| 93 | + """Return (projected_vector, clipped_mask). Identity where unbounded.""" |
| 94 | + |
| 95 | + |
| 96 | +class ClipperNone(AbstractClipper): |
| 97 | + """No-op. Bounds are ±inf, project is the identity. THE DEFAULT (see below).""" |
| 98 | + |
| 99 | + |
| 100 | +class ClipperPriorBox(AbstractClipper): |
| 101 | + """Hard projection onto the prior support, inset by a margin.""" |
| 102 | +``` |
| 103 | + |
| 104 | +`project` **must return which coordinates it clipped**, not just the new vector. |
| 105 | +That mask is what lets a caller zero the optimiser momentum along clipped |
| 106 | +directions. It is needed: the prototype left 5 of 16 lanes pinned to a bound at |
| 107 | +the end of the run because the parameters were projected while Prodigy's |
| 108 | +accumulated state kept pushing outward. The `Clipper` cannot fix that itself — it |
| 109 | +does not own `opt_state` — so it must expose enough for the search to. |
| 110 | + |
| 111 | +Later strategies (`ClipperReflect`, a soft-wall variant) drop in without touching |
| 112 | +callers. A soft wall must be a *Clipper* (search-local), **never** a change to the |
| 113 | +`Prior` classes — that would silently alter the objective for the nested samplers, |
| 114 | +where the hard box currently works correctly. |
| 115 | + |
| 116 | +## Scope — PR 1 (this task) |
| 117 | + |
| 118 | +1. `AbstractClipper` + `ClipperNone` + `ClipperPriorBox` in a new |
| 119 | + `@PyAutoFit/autofit/non_linear/clipper.py`. |
| 120 | +2. Bounds extraction covering **every** prior type. Confirmed present in the |
| 121 | + reference model: `UniformPrior` (finite both sides), `TruncatedGaussianPrior` |
| 122 | + (finite both sides, e.g. `(-1, 1)` for `ell_comps`), `GaussianPrior` |
| 123 | + (`±inf` — must pass through untouched). Audit the rest (`LogUniformPrior`, |
| 124 | + `LogGaussianPrior`, any `Constant`/deterministic entries). |
| 125 | +3. Wire into `AbstractMultiStartGradient._fit`, applied after |
| 126 | + `optax.apply_updates`, **opt-in**. |
| 127 | +4. Wire into `LBFGS`, passing `bounds=` through to `optimize.minimize`, |
| 128 | + **opt-in**. Only valid for bound-supporting methods (`L-BFGS-B`, `TNC`, |
| 129 | + `SLSQP`) — guard or warn for plain `BFGS`. |
| 130 | +5. `clipper: AbstractClipper = None` constructor arg on the searches, resolved |
| 131 | + like `initializer`. |
| 132 | + |
| 133 | +**Default is `ClipperNone`, and PR 1 must be bit-identical with it.** Follow the |
| 134 | +precedent set by PyAutoFit#1475, whose models declaring no constraint |
| 135 | +short-circuit to bit-identical behaviour. Flipping the default is a real |
| 136 | +behaviour change that shifts every stored multi-start benchmark, which is exactly |
| 137 | +the comparability argument PyAutoFit#1472 made when it deferred its own policy |
| 138 | +change. |
| 139 | + |
| 140 | +## Scope — PR 2 (separate prompt, file after PR 1 lands) |
| 141 | + |
| 142 | +Flip `MultiStartGradient`'s default to `ClipperPriorBox`, **with** the benchmark |
| 143 | +re-baseline, plus the momentum-reset-on-clip decision informed by how bad the |
| 144 | +pinning actually is at production budget. |
| 145 | + |
| 146 | +## Traps, measured |
| 147 | + |
| 148 | +- **Parameter ordering is load-bearing and silent if wrong.** |
| 149 | + `model.priors_ordered_by_id` was used for the prototype and lined up correctly |
| 150 | + with `model.instance_from_vector`, but a mismatch would clip the *wrong |
| 151 | + parameter* with no error. Assert the correspondence in a test rather than |
| 152 | + trusting it. |
| 153 | +- **Boundary semantics.** Decide and document whether `log_prior` at *exactly* the |
| 154 | + limit is finite. The prototype inset by `1e-6` of the box width to stay strictly |
| 155 | + inside; that margin is a guess and should be a justified constant. |
| 156 | +- **Pinning is correct behaviour, not a bug.** Where the likelihood genuinely |
| 157 | + prefers a value outside the prior, a clipped lane sitting on the bound is the |
| 158 | + correct MAP answer under the declared prior. It is worth surfacing (it says the |
| 159 | + prior is fighting the data) rather than hiding. In the reference cell the shear |
| 160 | + escapes were mixed-sign (`+0.353`, `-0.341`, `+0.301`, `-0.312`), which reads |
| 161 | + more like a poorly-constrained parameter diffusing out than a true value sitting |
| 162 | + outside. |
| 163 | +- **Clipping does not fix every death.** 5/16 lanes still died in the prototype; |
| 164 | + those are the NaN-gradient population (likelihood NaN in the *jitted* path, |
| 165 | + which the `Fitness` guard maps to `-inf` and whose `where` makes the gradient |
| 166 | + NaN). Separate mechanism, do not expect this task to remove it. |
| 167 | + |
| 168 | +## Two incidental bugs found while investigating — do not lose these |
| 169 | + |
| 170 | +Both surfaced only because clipping let lanes *survive*, i.e. on a code path this |
| 171 | +cell had apparently never taken: |
| 172 | + |
| 173 | +1. **`float32` is not JSON serializable in result output.** |
| 174 | + `@PyAutoFit/autofit/non_linear/paths/directory.py:80` `save_json` raises |
| 175 | + `TypeError: Object of type float32 is not JSON serializable` at the end of a |
| 176 | + successful clipped run. Did not fire on the baseline runs, where 14/16 lanes |
| 177 | + were dead. File separately if confirmed. |
| 178 | +2. **A crashed run poisons the next run of the same name.** The half-written |
| 179 | + output left by (1) caused the next search with the same `name` to fail with |
| 180 | + `JSONDecodeError` while trying to resume — a 4-second no-op run that *looked |
| 181 | + like* a clean result (zero deaths, because zero steps). This is a new form of |
| 182 | + the cached-result hazard already recorded in |
| 183 | + `complete/2026/08/multistart-nan-step-diagnostics.md`. |
| 184 | + |
| 185 | +## Deliberately out of scope |
| 186 | + |
| 187 | +- **NUTS.** Divergence, not lane death; may need a transform or a soft wall rather |
| 188 | + than projection. Its own task. |
| 189 | +- **Unit-cube stepping.** The more principled long-term fix — PyAutoFit's prior |
| 190 | + machinery is already unit-cube and the nested samplers work that way, and it |
| 191 | + would also normalise parameter scales (`einstein_radius ∈ [0,8]` alongside |
| 192 | + `ell_comps ∈ [-1,1]`). Rejected *for now* on three grounds: a logit |
| 193 | + reparameterisation sends the optimum to infinity when it genuinely sits on a |
| 194 | + boundary, which this cell demonstrably has; the inverse-CDF transform for |
| 195 | + non-uniform priors has `∂θ/∂u -> ∞` at the cube faces, trading one numerical |
| 196 | + hazard for another; and it invalidates every stored benchmark. If pursued, note |
| 197 | + that reparameterising the *search path* does not move the optimum **provided the |
| 198 | + objective is still the physical-space posterior evaluated at `θ(u)`** — optimise |
| 199 | + the density *of u* instead and the Jacobian makes the MAP non-invariant, which |
| 200 | + fails silently. |
| 201 | +- **Changing `resurrect` defaults.** Not the fix: a redrawn lane walks out again. |
| 202 | + |
| 203 | +## Testing |
| 204 | + |
| 205 | +- Bounds extraction per prior type, including `±inf` passthrough for `GaussianPrior`. |
| 206 | +- Ordering assertion (see traps). |
| 207 | +- `ClipperNone` is bit-identical: same seed, same final parameters, on both |
| 208 | + `MultiStartGradient` and `LBFGS`. |
| 209 | +- A lane deliberately stepped across a boundary is projected back inside, and the |
| 210 | + returned mask names exactly the crossed coordinates. |
| 211 | +- `LBFGS` passes bounds through and rejects/warns for non-bound-supporting methods. |
| 212 | +- Regression: with `ClipperPriorBox` on a model with a tight `UniformPrior`, the |
| 213 | + value-NaN rate falls substantially. The reference numbers above are CPU/float32, |
| 214 | + single seed — assert a direction and a large margin, not an exact figure. |
0 commit comments