Skip to content

New ADVI API - #635

Merged
zaxtax merged 8 commits into
pymc-devs:mainfrom
jessegrabowski:advi-refactor
Aug 27, 2026
Merged

New ADVI API #635
zaxtax merged 8 commits into
pymc-devs:mainfrom
jessegrabowski:advi-refactor

Conversation

@jessegrabowski

@jessegrabowski jessegrabowski commented Feb 2, 2026

Copy link
Copy Markdown
Member

This PR moves the work from pymc-devs/pymc#7799 over here to extras. The key idea is to copy the numpyro guide model API, but with our own PyMC flair.

I also added:

  • stick the landing estimator
  • forward sampling helper
  • LLM sketch for a training API patterned after pytorch-lightning

I updated the notebook to use the proposed training API. It obviously needs a lot of work (being that it's llm trash) but hopefully it can get some ideas flowing.

@review-notebook-app

Copy link
Copy Markdown

Check out this pull request on  ReviewNB

See visual diffs & provide feedback on Jupyter Notebooks.


Powered by ReviewNB

@codecov-commenter

codecov-commenter commented Feb 2, 2026

Copy link
Copy Markdown

Codecov Report

❌ Patch coverage is 93.42857% with 23 lines in your changes missing coverage. Please review.
✅ Project coverage is 92.99%. Comparing base (86fac3c) to head (380120a).
⚠️ Report is 53 commits behind head on main.

Files with missing lines Patch % Lines
pymc_extras/inference/advi/training.py 93.37% 11 Missing ⚠️
pymc_extras/inference/advi/optimizers.py 92.36% 10 Missing ⚠️
pymc_extras/inference/advi/compile.py 96.55% 1 Missing ⚠️
pymc_extras/inference/advi/fit.py 94.44% 1 Missing ⚠️
Additional details and impacted files

Impacted file tree graph

@@             Coverage Diff             @@
##             main     #635       +/-   ##
===========================================
+ Coverage   51.60%   92.99%   +41.38%     
===========================================
  Files          73      107       +34     
  Lines        8003     9749     +1746     
===========================================
+ Hits         4130     9066     +4936     
+ Misses       3873      683     -3190     
Files with missing lines Coverage Δ
pymc_extras/inference/advi/__init__.py 100.00% <100.00%> (ø)
pymc_extras/inference/advi/autoguide.py 98.11% <100.00%> (+98.11%) ⬆️
pymc_extras/inference/advi/compile.py 97.77% <96.55%> (ø)
pymc_extras/inference/advi/fit.py 94.44% <94.44%> (ø)
pymc_extras/inference/advi/optimizers.py 92.36% <92.36%> (ø)
pymc_extras/inference/advi/training.py 93.41% <93.37%> (+93.41%) ⬆️

... and 44 files with indirect coverage changes

🚀 New features to boost your workflow:
  • ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

@zaxtax zaxtax left a comment

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think this is pretty neat! After some rebasing and small API changes we should definitely try to get this merged!

Comment thread pymc_extras/inference/advi/objective.py Outdated
The probabilistic model.
guide : AutoGuideModel
The variational guide.
stick_the_landing : bool, optional

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I feel like we need a better name than stick_the_landing . Also the function is supposed to return the logp and logq terms but the STL estimator is about returning only the path derivative component of the gradient.

Comment thread pymc_extras/inference/advi/training.py Outdated
Number of MC draws per step for gradient estimation, by default 10.
model : Model
The PyMC model to fit. If None, the model is inferred from context.
state : SVIState, optional

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'd prefer we copy what PyTorch Lightning did here and have the SVI State live in the Trainer object and be something we can pass upon initalisation.

Comment thread pymc_extras/inference/advi/training.py Outdated
return state

def sample_posterior(
self, draws: int, state: SVIState, model: Model | None = None

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

drop SVIState and move into the Trainer object

@ricardoV94
ricardoV94 force-pushed the advi-refactor branch 3 times, most recently from 6525424 to a661c84 Compare June 12, 2026 16:29
Comment thread pymc_extras/inference/advi/fit.py Outdated
model: Model | None = None,
*,
n_steps: int = 10_000,
draws_per_step: int = 8,

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
draws_per_step: int = 8,
draws_per_step: int = 1,

@ricardoV94
ricardoV94 force-pushed the advi-refactor branch 3 times, most recently from 68e83fe to eefecfb Compare June 16, 2026 14:23
@fonnesbeck

Copy link
Copy Markdown
Member

Is the plan to deprecate the OPVI fit function when this merged? (or matriculates to the main project?)

Deterministic(
rv.name,
loc + pt.softplus(scale) * z,
loc + pt.exp(scale) * z,

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We had it with softplus since that's what numpyro does? Why change?

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

2x faster, logdet becomes -scale, instead of log(softplus)

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It's not a sure choice, just something I'm exploring

@ricardoV94

Copy link
Copy Markdown
Member

Is the plan to deprecate the OPVI fit function when this merged? (or matriculates to the main project?)

There's no big plan on my end, still very much WIP. Maybe if it proves itself. The VI part of this is def smaller in scope (but also complexity)

@jessegrabowski

jessegrabowski commented Jun 18, 2026

Copy link
Copy Markdown
Member Author

Is the plan to deprecate the OPVI fit function when this merged? (or matriculates to the main project?)

From my part, this refactor was an effort to explore a different design. Could maybe lead to a full replacement of pm.fit, but the scope so far is just ADVI. I'm not sure the particle methods would be in scope for what we have in this PR.

But I think the guide model setup is really interesting, and works really well with the pymc Model API + pytensor graph mutations. So its appealing in and of itself.

@ricardoV94

Copy link
Copy Markdown
Member

@zaxtax let's leave the minibatch out of this. I already had quite some discussion with @jessegrabowski I don't think PyMC models have the right API for us to eagerly go and guess what's batched data and what's not, and assume it's over the 0-th axis and so forth. Check the over-enginereered #700 I started for what I think is the minimum work a library built on top of PyMC would have to do to "auto minibatch". I would not put any minibatch stuff on this PR. For now if users want to minibatch they'd use the lover-or-hate it pm.Minibatch path

@zaxtax

zaxtax commented Jul 5, 2026 via email

Copy link
Copy Markdown
Contributor

@ricardoV94

ricardoV94 commented Jul 5, 2026

Copy link
Copy Markdown
Member

But the proposed API here offers a high degree of control, so I'm not sure what's ambiguous here

That's exactly the sort of discussion I didn't want to have here... Anyway I see no place to control which X variables you stream and how. Your proposal is a bit more low-level and less debatable than an "auto minibatch", but I still would prefer to leave it to its own PR. The scaling thing is also a bit naive (in terms of calling model.logp() per var which can be costly). I didn't check how it handles univariate vs multivariates. Also why merging observed transform into the fitting routine?

@zaxtax

zaxtax commented Jul 5, 2026 via email

Copy link
Copy Markdown
Contributor

@zaxtax

zaxtax commented Jul 17, 2026 via email

Copy link
Copy Markdown
Contributor

@ricardoV94

Copy link
Copy Markdown
Member

@ricardoV94 can I close this PR as superseded by the other ones?

no, but this needs rebasing of things that went elsewhere

@ricardoV94

Copy link
Copy Markdown
Member

Did my Trainer class commit end up in another PR?

I reverted the minibatch commit, the idea was to discuss in another PR yes. I didn't open yet

@zaxtax

zaxtax commented Jul 17, 2026 via email

Copy link
Copy Markdown
Contributor

@ricardoV94

ricardoV94 commented Jul 17, 2026

Copy link
Copy Markdown
Member

I tried to be careful to only revert, not delete stuff that wasn't mine. But may have failed

@zaxtax

zaxtax commented Jul 17, 2026

Copy link
Copy Markdown
Contributor

I mean this one c186229

@zaxtax

zaxtax commented Jul 17, 2026 via email

Copy link
Copy Markdown
Contributor

@ricardoV94

Copy link
Copy Markdown
Member

I'll take a look, but I suspect we may have 3 chefs who want 3 slight different things here.

Personally, I'd distinguish fit machinery state and the fit result. For the fit machinery state we should try to reuse the pymc step stampler state api/machinery (not the step samplers, those are different things). https://github.com/pymc-devs/pymc/blob/df1a369d5fbdc0cc49082b0c0ef7cbd511cf76e9/pymc/step_methods/state.py#L54

We can reject if we find a good reason why it doesn't fit the mold. It doesn't need to happen in this PR either, but I'd like to plan on checking it.

Or data structures / classes, I don't want to incuur overhead of storing each step result in a fancy dataclass or promise a state that was stale (i.e. read the shared variables content on the first iteration and pretend those are up to date in the pure pytensor path). I don't know if I pushed the latest state of the api notebook, but most my decisions were based on bringing inference speed down, and they checked out there (in the last example brought from @juanitorduz blogpost). When the inner step function gets cheap I don't want to be losing to jax/numpyro.

My one criteria is I don't want @juanitorduz to have an excuse to use numpyro for SVI. Or at least not a good one.

On the control over the fit loop, and what we expose vs not I don't have much opinion, but both you @jessegrabowski and @zaxtax need to check if you're on the same page. I didn't look into that at all but may have accidentally reverted or erased something (sorry if erased, it wasn't supposed to if it happened)

@zaxtax

zaxtax commented Jul 17, 2026

Copy link
Copy Markdown
Contributor

I'll take a look, but I suspect we may have 3 chefs who want 3 slight different things here.

Personally, I'd distinguish fit machinery state and the fit result. For the fit machinery state we should try to reuse the pymc step stampler state api/machinery (not the step samplers, those are different things). https://github.com/pymc-devs/pymc/blob/df1a369d5fbdc0cc49082b0c0ef7cbd511cf76e9/pymc/step_methods/state.py#L54

We can reject if we find a good reason why it doesn't fit the mold. It doesn't need to happen in this PR either, but I'd like to plan on checking it.

Or data structures / classes, I don't want to incuur overhead of storing each step result in a fancy dataclass or promise a state that was stale (i.e. read the shared variables content on the first iteration and pretend those are up to date in the pure pytensor path). I don't know if I pushed the latest state of the api notebook, but most my decisions were based on bringing inference speed down, and they checked out there (in the last example brought from @juanitorduz blogpost). When the inner step function gets cheap I don't want to be losing to jax/numpyro.

My one criteria is I don't want @juanitorduz to have an excuse to use numpyro for SVI. Or at least not a good one.

On the control over the fit loop, and what we expose vs not I don't have much opinion, but both you @jessegrabowski and @zaxtax need to check if you're on the same page. I didn't look into that at all but may have accidentally reverted or erased something (sorry if erased, it wasn't supposed to if it happened)

I don't have a strong opinion here but I worry sampler state is a poor fit for the different state representations ADVI methods can have. As for the fit loop, all I ask is @jessegrabowski have a peek about what I made in the existing commit (possibly the reverted one) and we confirm if we are on the same page. If we are not, I'm happy for us to just pick something and go with it

@juanitorduz

juanitorduz commented Jul 17, 2026

Copy link
Copy Markdown

I will be happy to test this against a real-life example https://github.com/pymc-labs/pymc-forecast/blob/main/docs/examples/retail_stockouts.ipynb :D.

In the meantime,I wanna try to run this (as in, the current api, not this one) on GPU + Modal (I made the numpyro version fit this model in less than 10 min with SVI and 50K time series), https://github.com/pymc-labs/pymc-forecast/blob/main/docs/examples/retail_stockouts.ipynb is the same example, and we tried to port the batching machinery pymc-labs/pymc-forecast#48

@zaxtax

zaxtax commented Jul 17, 2026

Copy link
Copy Markdown
Contributor

And yes I want @juanitorduz to love whatever we make

@juanitorduz

Copy link
Copy Markdown

And yes I want @juanitorduz to love whatever we make

I AM loving it ... that is why we are pushing with @twiecki a testing field : https://github.com/pymc-labs/pymc-forecast :D

@ricardoV94

Copy link
Copy Markdown
Member

I will be happy to test this against a real-life example https://github.com/pymc-labs/pymc-forecast/blob/main/docs/examples/retail_stockouts.ipynb :D.

In the meantime,I wanna try to run this (as in, the current api, not this one) on GPU + Modal (I made the numpyro version fit this model in less than 10 min with SVI and 50K time series), https://github.com/pymc-labs/pymc-forecast/blob/main/docs/examples/retail_stockouts.ipynb is the same example, and we tried to port the batching machinery pymc-labs/pymc-forecast#48

what batching machinery? minibatch or vectorization -adjacent?

@juanitorduz

juanitorduz commented Jul 17, 2026

Copy link
Copy Markdown

I will be happy to test this against a real-life example https://github.com/pymc-labs/pymc-forecast/blob/main/docs/examples/retail_stockouts.ipynb :D.
In the meantime,I wanna try to run this (as in, the current api, not this one) on GPU + Modal (I made the numpyro version fit this model in less than 10 min with SVI and 50K time series), https://github.com/pymc-labs/pymc-forecast/blob/main/docs/examples/retail_stockouts.ipynb is the same example, and we tried to port the batching machinery pymc-labs/pymc-forecast#48

what batching machinery? minibatch or vectorization -adjacent?

For posterior predictive sampling. In the NumPyro function Predictive() this was killing the memory, and the "trick" that worked was to move this computation to the CPU and do it by batches.

@zaxtax

zaxtax commented Jul 17, 2026 via email

Copy link
Copy Markdown
Contributor

self.convergence_window = convergence_window
self.relative_tolerance = relative_tolerance

self.model = model

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should we make Trainer object modelcontext aware?

Comment thread pymc_extras/inference/advi/training.py Outdated
return self._guide

@property
def clip_norm(self) -> float | None:

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do we need this?

Comment thread tests/inference/advi/test_schedules.py Outdated

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do we need learning rate schedulers in this PR?

Comment thread tests/inference/advi/test_fit.py Outdated
assert not np.array_equal(draws_a, draws_c)


def test_fit_advi_early_stopping(conjugate_model):

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do we need early stopping support in this PR?

@zaxtax

zaxtax commented Aug 20, 2026

Copy link
Copy Markdown
Contributor

@ricardoV94 what's left to do on this one other than maybe remove learning rate and clip norm from the Trainer object?

ricardoV94 and others added 7 commits August 26, 2026 20:53
fit_advi entry point with an optax-like optimizer API and numpyro-style
defaults, a compiled SVI step that bakes in draws and runs optimizer updates
(Adam bias-correction kept in floatX), the training loop, and the objective.
Applies to the mean-field, full-rank and low-rank guides, so the scale
parameterization is the same across all three.
Replace the SVIModule/ADVIModule/SVITrainer split with one Trainer that
owns the training loop, following the design of pymc-devs/pymc#8333 and
PyTorch Lightning: all configuration (guide, optimizer, learning rate,
convergence-based early stopping, model, backend) lives at construction,
there are no user-facing hooks or callbacks, and fit(n) just runs.

The duplicate fit/fit_jitted loops are merged into a single fit that
internally picks the compiled fast path (clipped Adam baked into the
step function; the default) or the Python-side update loop (when an
optax-like optimizer is passed). The fast path now also supports
resuming parameters from a passed SVIState. The trainer keeps the last
state, so sample_posterior() works without arguments. fit_advi becomes
a thin wrapper over Trainer, and draws_per_step is renamed to
n_particles and moved to the constructor.
fit(n) now always continues from the current parameters and Adam moments,
replacing a state=None argument that meant "start fresh" on the python path
and "silently resume" on the compiled one. reset(), load_state() and the
state property cover those cases explicitly, and SVIState carries the Adam
moments so resuming in another trainer matches continuing in place.

State is read out of the shared variables once at the end of a fit rather
than rebuilt per step, so the loop no longer allocates a state object it
cannot fill in honestly. Learning rates are resolved up front for the whole
call, lifting the schedule (and its np.interp) out of the hot loop.

Configuration splits the same way: what gets compiled into the step function
is read-only, what is per-run policy is read afresh by each fit call, with a
learning_rate override for a single call.

Drops the python-side optimizer loop, which cost 61us/step against 3.6us for
the compiled path, was used by no test, notebook or default, and was where
both divergences came from. optimizers.py loses the unreachable optax-like
transformations and becomes schedules.py.

fit_advi is unchanged bit-for-bit.
Cut it from 76 cells to 26, dropping the SVML/llvmlite probes, the timeit
scratch cells, the dead SGD and RMSProp optimizer classes, and the duplicate
benchmark runs. What remains is fit_advi and Trainer on three models:
a linear regression against NUTS, the Trainer API itself (fit continuing,
state snapshots, load_state, reset, per-call learning rate), radon with an
LKJ covariance, and the 116k-parameter forecasting benchmark.

The radon comparison is a two-panel scatter over all 85 counties rather than
an 85-row forest plot: same information, and 3.4MB less of it.

The benchmark runs 10_000 steps instead of 3_000. At 3_000 the one-cycle rate
has fully annealed while sigma is still descending, which read as fit_advi
losing to pm.fit on quality when it was only losing on step count. All three
now land on sigma 0.313 and seasonal correlation 0.876.

Also drops the claim that numba loses to XLA here. On pytensor
catch_up_with_jax at 4be21562a the ADVI step is 6.7ms on numba against 7.7ms
on jax; what is left is numba's cold compile, which the notes now say instead.
fit(n) now runs n steps. The window check interacted badly with the default
schedule: the one-cycle is sized to the requested n, so stopping early exits
wherever that schedule had reached, which for an early stop is near the peak.
On radon it fired at 2800 of 10000 steps at lr 0.00706 against a 0.008 peak,
i.e. at close to maximum step size. Spending the same 2800 steps with the
schedule sized to 2800 lands fully annealed and gets 2.2x lower mean error
against a nutpie reference (0.0122 vs 0.0273), so the heuristic was worse than
simply asking for fewer steps.

That leaves n as a pacing parameter and not just a cap, which the docstrings
now say. Deciding when to stop by hand is better served by what the trainer
already does: fit in chunks, look at state.loss_history, and call fit again to
continue from the same parameters and Adam moments. KeyboardInterrupt also
keeps the state, so an interrupted run is resumable.

Schedules are now functions of the trainer's global step rather than of the
offset within a call, and the default one-cycle is sized to start_step + n. A
follow-up fit anneals the rest of one cycle instead of ramping a second time,
which was the last thing about resuming that did not carry over. Resuming from
a snapshot still matches continuing in place, since load_state restores the
step count the schedule reads.

Removes convergence_window and relative_tolerance from Trainer and fit_advi,
and reruns the notebook.
…imizers, add sgd

Port the optimizer refactor from advi-minibatch into advi-refactor.

- Add optimizers.py with GradientTransformation (init/update/pytensor), adam,
  clipped_adam, sgd, rmsprop, chain, clip_by_global_norm, and
  linear_onecycle_schedule; remove schedules.py.
- compile_svi_step_fn now takes an optimizer and returns
  (step_fn, shared_params, shared_optimizer_state); step_fn takes no inputs.
- Trainer takes an optimizer instead of learning_rate/clip_norm and defaults
  to clipped_adam(); fit() no longer accepts a per-call learning_rate.
- fit_advi takes an optimizer instead of learning_rate/clip_norm.
- Preserve the optimizer-state snapshot/resume feature: SVIState.optimizer_state
  is read from the shared variables, so state/load_state/reset keep working
  (empty for stateless optimizers like sgd).
@zaxtax
zaxtax merged commit d48872a into pymc-devs:main Aug 27, 2026
29 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

enhancements New feature or request help wanted Extra attention is needed inference

Projects

None yet

Development

Successfully merging this pull request may close these issues.

6 participants