# Spherical analysis with saved posteriors

After the [step-by-step tutorials](index.md), use this worked example to compare
the complete NumPy/emcee and JAX/NumPyro storage workflows on the same data.
Unlike the minimal Quickstart, these scripts include checkpoint restart and
additional diagnostics.

Build a **Plummer + generalized NFW + constant anisotropy** model, generate
mock velocities, run MCMC, resume from disk, and inspect the saved posterior.
The generalized NFW halo uses `ZhaoModel` with fixed `a=1`, `b=3` and a sampled
inner slope `gamma`. The small dataset lets the complete example run on a CPU.

Install the package using the [installation guide](../installation.md). The emcee
example needs the `plotting` extra; the NumPyro example also needs
`numpyro_cpu`. Run the scripts below from the source checkout. Use a new output
directory for each example.

## Model and mock observations

Both backends use the same 32 mock stars. Lengths are in pc, velocities in km/s,
and halo densities in solar masses per cubic parsec. The Plummer tracer weights
the kinematics and contributes no gravitational mass. The tracer scale, Zhao
exponents `alpha=1, beta=3` and untruncated halo cutoff are fixed. The mock has `gamma=1`
(standard NFW), while inference allows the inner slope to vary from 0 to 2.

```{literalinclude} ../../../examples/docs_quickstart_data.py
:language: python
:start-after: mock-start
:end-before: mock-end
```

```{figure} ../_static/quickstart/observations.png
:alt: Synthetic line-of-sight velocities versus projected radius with measurement error bars.
:width: 640px

Mock observations generated by the code above, with random seed 123.
```

We sample five parameters with the same uniform priors in both backends:

| Sampling coordinate | Lower | Upper | Physical meaning |
| --- | ---: | ---: | --- |
| `log10_rs_pc` | 1.5 | 4.0 | Logarithm of the halo scale radius |
| `log10_rhos_Msunpc3` | -3.0 | 1.0 | Logarithm of the halo density scale |
| `gamma` | 0.0 | 2.0 | Inner halo density slope |
| `beta_ani` | -1.0 | 0.75 | Constant spherical anisotropy |
| `vmem_kms` | -50 | 50 | Systemic velocity in km/s |

The priors on the logarithms are log-uniform in the physical scales. These
bounds illustrate this mock analysis; choose priors for the scientific problem
being studied. With five parameters and only 32 stars, degeneracies and prior
sensitivity are expected; see the [MCMC notebook](inference.ipynb). The standalone
scripts share one mock dataset generated by the NumPy solver at `n=256`,
`n_kernel=64`; the notebooks generate their own mocks with their displayed JAX
settings, so their numerical outputs need not match these scripts exactly.

## Sample and inspect the posterior

`````{tab-set}
````{tab-item} NumPy / SciPy / emcee
:sync: emcee

Run the complete example:

```bash
python examples/docs_quickstart_emcee.py --output-dir /tmp/jeanspy-emcee-example
```

**Construct the model.** NumPy/SciPy components store their physical parameters.

```{literalinclude} ../../../examples/docs_quickstart_emcee.py
:language: python
:start-after: imports-start
:end-before: imports-end
```

```{literalinclude} ../../../examples/docs_quickstart_emcee.py
:language: python
:start-after: model-start
:end-before: model-end
```

**Define the posterior, sample and resume.** The Gaussian velocity likelihood
adds the measurement-error variance to the Jeans prediction. The HDF5 backend
stores the ensemble and its random state. Reopening that backend and passing
`None` continues from the saved state. Keep the model, priors and observations
identical when using this emcee interface; the JeansPy
[`Sampler`](../api/generated/jeanspy.sampler.Sampler.rst) adds analysis-identity
checks for its supported estimation-model classes.

```{literalinclude} ../../../examples/docs_quickstart_emcee.py
:language: python
:start-after: inference-start
:end-before: inference-end
```

**Recorded execution result.** The first run stores 256 steps; the resumed run adds 128.
The first 128 steps are discarded when plotting and summarizing this example.

```{literalinclude} ../_static/quickstart/emcee.txt
:language: text
```

```{figure} ../_static/quickstart/emcee-trace.png
:alt: Posterior histograms and stored emcee walker traces for five sampled parameters.
:width: 100%

Posterior histograms and walker traces. Dashed red lines mark the generating values.
```

```{figure} ../_static/quickstart/emcee-autocorrelation.png
:alt: Autocorrelation versus lag for each of the five emcee sampling coordinates.
:width: 640px

Mean autocorrelation across walkers. A short-chain estimate of its integrated
time is exploratory; interacting walkers are not independent chains for R-hat.
```

```{figure} ../_static/quickstart/emcee-posterior.png
:alt: Marginal and pairwise emcee posterior distributions, with generating values and central 68 percent intervals.
:width: 100%

The corner plot shows 68% and 95% enclosed-probability contours in each
2-D marginal. The diagonal dashed lines mark 16%, 50% and 84% quantiles;
black lines mark the generating values. Contours are estimated from smoothed
histograms. A mock-data posterior need not peak at the generating value.
```
````
````{tab-item} JAX / NumPyro
:sync: numpyro

Run the complete example in double precision:

```bash
JEANSPY_JAX_PLATFORM=cpu JEANSPY_JAX_ENABLE_X64=true \
python examples/docs_quickstart_numpyro.py --output-dir /tmp/jeanspy-numpyro-example
```

**Construct the model.** JAX components receive physical parameters explicitly.
Import the components from `jeanspy.model_jax` and the likelihood and
sampler from `jeanspy.sampler_numpyro`.

```{literalinclude} ../../../examples/docs_quickstart_numpyro.py
:language: python
:start-after: imports-start
:end-before: imports-end
```

```{literalinclude} ../../../examples/docs_quickstart_numpyro.py
:language: python
:start-after: model-start
:end-before: model-end
```

**Define priors and run NUTS.** `ParameterSpec` associates sampling coordinates
with their physical parameters. `NumPyroSampler` saves a checkpoint and a
backend-native ArviZ chunk. A fresh sampler instance resumes without repeating
warmup, after checking that the analysis and environment match the stored state.

```{literalinclude} ../../../examples/docs_quickstart_numpyro.py
:language: python
:start-after: inference-start
:end-before: inference-end
```

**Recorded execution result.** Two chains run sequentially, each with 200 warmup steps
and 256 posterior draws. Resuming adds another 256 draws per chain. The ArviZ
summary and divergences below are calculated from the combined saved chunks.

```{literalinclude} ../_static/quickstart/numpyro.txt
:language: text
```

```{figure} ../_static/quickstart/numpyro-trace.png
:alt: Posterior histograms and two stored NumPyro chain traces for five sampled parameters.
:width: 100%

Posterior histograms and the two independent chain traces. Dashed red lines mark
the generating values.
```

```{figure} ../_static/quickstart/numpyro-autocorrelation.png
:alt: Autocorrelation versus lag for each of the five NumPyro sampling coordinates.
:width: 640px

Mean autocorrelation across the two chains.
```

```{figure} ../_static/quickstart/numpyro-posterior.png
:alt: Marginal and pairwise NumPyro posterior distributions, with generating values and central 68 percent intervals.
:width: 100%

Marginal and pairwise posterior distributions, with the same conventions as the
emcee panels.
```
````
`````

For rank-normalized R-hat, bulk/tail ESS, MCSE, BFMI, tree depth and locating
NUTS divergences in a corner plot, follow the [MCMC diagnostics tutorial](inference.ipynb).

The displayed outputs are retained from the execution identified by the
[metadata](../_static/quickstart/execution.json); a documentation rebuild does
not rerun them. These are short tutorial runs. Inspect chain mixing, effective sample sizes,
R-hat for independent chains and divergences, and extend the analysis as needed.
Good sampler diagnostics do not establish the physical adequacy of a Jeans
model or the calibration of its intervals.

## Plotting and saved files

The figures use the posterior draws loaded from disk. The shared plotting
function also makes the dimensional convention explicit:

````{dropdown} Plotting code
```{literalinclude} ../../../examples/docs_quickstart_plots.py
:language: python
```
````

Each output directory contains `observations.csv`, the three posterior figures
and the mock-observation figure. The emcee directory contains `chain.h5`;
the NumPyro directory contains `metadata.json`, `last_state.pkl` and two
`chunks/*.nc` stores. See [saving and resuming](storage.ipynb) for longer
analyses and [axisymmetric models](../guides/axisymmetric.md) for flattened systems.

The [execution metadata](../_static/quickstart/execution.json) records package
versions, source hashes and output hashes. To regenerate the displayed results
from a checkout, run `python scripts/run_quickstart.py` in the locked docs
environment. Ordinary documentation builds reuse these recorded outputs; the
metadata identifies the source and environment used to produce them. Both
workflows run again on a published GitHub release or a manual Documentation
workflow run with `run_mcmc=true`.

<a id="classical-estimation-model-wrapper"></a>

(numpy-scipy-estimation-model-wrapper)=
## NumPy/SciPy estimation-model wrapper

The analysis above supplies its own log posterior to emcee. When using
JeansPy estimation models, {class}`~jeanspy.sampler.Sampler` adds HDF5
storage and restart identity checks. Its `FlatPriorModel` uses an ordered
table with finite `lower` and `upper` bounds. Spherical convenience models
also use an explicit photometric prior on `log10_re_pc`; a systemic-velocity
prior is not inferred from the velocities unless the caller explicitly
requests the data-derived option.

This separate example writes and resumes six steps in a temporary directory.
Its autocorrelation estimate may be undefined because the chain is deliberately
short. The wrapper still computes that diagnostic when early stopping is
disabled; the example only exercises the storage workflow.

```{literalinclude} ../../../examples/docs_inference.py
:language: python
:end-before: numpy-inference-end
```

Choose the same likelihood and prior measure when comparing samplers. Keep
raw draws and diagnostics when a stopping criterion fails; convergence and
calibration need their own checks as discussed in the [MCMC tutorial](inference.ipynb).