`virgil.likelihood`

Likelihoods of source models given data, and numpyro models.

Likelihoods of source models given interferometric data.

A model is given either as a template SourceModel, whose parameters at dot-separated zodiax paths (e.g. "comp.flux") are replaced, or as a class/callable called with the parameters as keyword arguments (see build_model).

Every likelihood, grid, limit and fit goes through one residual vector, whitened_residuals: the residuals divided by their uncertainties, with unprojected phases measured as a chord, 2 sin(Δ/2), so that the likelihood is smooth where phases wrap at ±π. Closure phases from four or more telescopes are correlated: their independent combinations are whitened together, using the sine sin Δ (not the chord) of each residual, with a periodic penalty 2 sin²(Δ/2) / σ added per closure phase. Both terms repeat every 2π and are smooth, so the likelihood is continuous everywhere.

`whitened_residuals(model, data, **noise)`

Residuals of a model divided by the data uncertainties.

The one residual vector behind every likelihood in virgil: model_loglike is -0.5 * sum(whitened_residuals**2) plus the normalization, and least-squares fits minimize its sum of squares. Visibility and projected-phase residuals are (model - data) / σ; an unprojected phase residual Δ is the chord 2 sin(Δ/2) / σ; correlated closure phases are replaced by whitened independent combinations plus periodic penalty rows. See conventions and closure phases with four or more telescopes.

Parameters:

Name	Type	Description	Default
`model`	`SourceModel`	Model to evaluate.	required
`data`	`OIData`	Data to compare with.	required
`**noise`		Error terms (see common arguments).	`{}`

Returns:

Type	Description
`array - like`	One dimensionless residual per independent observable, in the order of `flatten_data`; the length is `n_residuals`.

`build_model(model, params, values)`

Build a model from parameter names and values.

model is either a class/callable, called as model(**dict(zip(params, values))), or a SourceModel instance used as a template whose leaves at the (dot-separated) paths params are replaced by values.

`model_loglike(model, data, *, reject_unphysical=False, **noise)`

Evaluate the log likelihood for an instantiated model object.

This is -0.5 * sum(r**2) for the residuals r of whitened_residuals, plus each density's normalization (Gaussian, or von Mises for unprojected phases; see likelihood normalization).

Parameters:

Name	Type	Description	Default
`model`	`SourceModel`	Model to evaluate.	required
`data`	`OIData`	Data to compare with.	required
`reject_unphysical`	`bool`	If True, return `-inf` when `is_physical` is false, e.g. for a negative flux. Works inside `jax.jit`, as a hard prior boundary.	`False`
`**noise`		Error terms (see common arguments). The normalization uses the inflated errors; for data with gains it includes the log-determinant of the marginalized covariance.	`{}`

`inflated_errors(data, prediction, vis_error_rel=None, phi_error=None, vis_scale=None, phi_scale=None, vis_error=None, *, where='model', combine='quadrature')`

The data uncertainties, scaled and with extra terms in quadrature.

The visibility errors become hypot(vis_scale σ, vis_error, vis_error_rel V), for the model visibility observable V, and the phase errors hypot(phi_scale σ, phi_error). Terms left as None are not applied. Extra observables (OIData.extras) are left unchanged; give them floors with OIData.with_error_floor.

Parameters:

Name	Type	Description	Default
`data`	`OIData`	Data whose uncertainties are inflated.	required
`prediction`	`array - like`	Model vector, e.g. from `OIData.model`.	required
`vis_error_rel`	`float`	Extra visibility error, as a fraction of the model visibility observable (e.g. of the model V² for squared visibilities): a calibration error.	`None`
`phi_error`	`float`	Extra phase error in radians.	`None`
`vis_scale`	`float`	Factors multiplying the visibility and phase uncertainties.	`None`
`phi_scale`	`float`	Factors multiplying the visibility and phase uncertainties.	`None`
`vis_error`	`float`	Extra absolute visibility error, in the units of the observable.	`None`
`where`	`(model, data)`	What `vis_error_rel` is relative to: the model (default, the right choice for a fitted term, as it does not reward the model for low data points) or the data (as error floors are).	`"model"`
`combine`	`(quadrature, max)`	Add the terms in quadrature (default), or take the largest (a floor, as `with_error_floor` does; both use `_utils.inflate_errors`).	`"quadrature"`

Returns:

Type	Description
`array - like`	Uncertainties matching `flatten_data`.

`flux_scale_posterior(model, data)`

The grey scales of the extra spectra, given a model.

The scales (and polynomial coefficients) of OI_FLUX spectra and correlated fluxes are marginalized in the likelihood (see FluxSpectrum); this is their Gaussian posterior conditional on model, for reporting.

Parameters:

Name	Type	Description	Default
`model`	`SourceModel`	E.g. the best fit.	required
`data`	`OIData`	Data with extra spectra (`extras=("flux",)` and the like).	required

Returns:

Type	Description
`dict`	Per kind (`"flux"`, `"nflux"`, `"corrflux"`): `mean` `(n_group, p)` and `cov` `(n_group, p, p)` of the weights, whose first is the scale k multiplying the model's template (normalized to a mean of 1 per group); and `scale`, the scale of `total_spectrum` itself (k over the template's normalization), for `"flux"` and `"corrflux"`: data ≈ scale × model. `groups` gives each sample's group.

`noise_sites(noise, n_datasets)`

Expand a noise specification into named sites.

noise maps error-inflation terms (NOISE_TERMS) to priors and applies to every dataset, giving sites "noise.<term>"; a list of such dicts, one per dataset, gives sites "noise[i].<term>". An entry may also be tied: a function of the dict of sampled parameter values (keyed like priors) that returns the term's value, which then has no prior of its own (see is_tied), unless it has a log_prior(values) method (see tied_log_prior).

Returns:

Type	Description
`dict`	`{site: (prior, datasets, term)}`, `datasets` being the indices of the datasets the term applies to; `prior` is the function for a tied term.

`is_tied(spec)`

Whether a noise entry is tied: a function of the sampled parameters (rather than a prior), e.g. one epoch's error scale drawn from a population (see hierarchical_scales).

`tied_log_prior(sites, values)`

The summed log_prior(values) of the distinct tied noise terms that have one, such as the population density of a centred hierarchical_scales member (a term used twice counts once).

`noise_for(sites, values, index)`

The error-inflation terms of dataset index, from site values.

`loglike(values, params, model, data, **options)`

Gaussian log-likelihood of a model with the given parameter values, assuming Gaussian errors.

Parameters:

Name	Type	Description	Default
`values`	`array - like`	Values of the model parameters.	required
`params`	`list`	List of parameter names.	required
`model`	`SourceModel or callable`	Template model whose parameters at the dot-separated paths `params` are replaced by `values`, or a class/callable called as `model(**dict(zip(params, values)))` (see `build_model`).	required
`data`	`OIData`	Object containing the data to be fitted.	required
`**options`		Error terms and `reject_unphysical`, passed to `model_loglike`.	`{}`

Returns:

Type	Description
`float`	Log-likelihood value.

`joint_prediction(params, model, data)`

Concatenate predictions for a parameter pytree and multiple data.

model(params, index) defines which parameters are shared and which are specific to each observation.

`joint_data(data)`

Concatenate observed vectors in the same order as joint_prediction.

`joint_errors(data)`

Concatenate uncertainty vectors in the same order as joint_prediction.

`joint_loglike(params, model, data, **options)`

Sum independent Gaussian log likelihoods over multiple data.

options (error terms and reject_unphysical) are passed to model_loglike.

`numpyro_model(model, priors, data, regularizers=(), noise=None, likelihoods=(), *, flat_coordinates=False, **options)`

Return a numpyro model sampling the parameters in priors.

Arguments are shared with fit; see common arguments.

Parameters:

Name	Type	Description	Default
`model`	`SourceModel or callable`	A template whose leaves at the paths in `priors` are sampled, or a function of the sampled values returning a `SourceModel` or a list of them, one per dataset (see `build_model`).	required
`priors`	`dict[str, Distribution]`	Parameter path (or function argument name) to prior; each key is also the sample-site name. An `AngleVector` prior samples a 2-D vector at `"<path>_vec"` and records the angle at `"<path>"`.	required
`data`	`OIData or sequence of OIData`	Gaussian log likelihood, added with `numpyro.factor`; may be `()` when `likelihoods` holds all the data.	required
`regularizers`	`sequence`	Log-prior terms, e.g. `Centroid`.	`()`
`noise`	`dict or list of dict`	Priors on error and calibration terms, sites `"noise.<term>"` (`"noise[i].<term>"` for a list); see common arguments.	`None`
`likelihoods`	`sequence`	Extra data terms, sites `"likelihood_<i>"`; see common arguments.	`()`
`flat_coordinates`	`bool`	Sample each prior in its flat coordinate (default `False`; see common arguments). Pass the same value to `gauss_newton_mass`.	`False`
`**options`		Fixed error terms and `reject_unphysical`, passed to `model_loglike`.	`{}`

Returns:

Type	Description
`callable`	Zero-argument numpyro model, e.g. for `numpyro.infer.NUTS`.

Examples:

>>> posterior = numpyro_model(model, priors, data)
>>> mcmc = MCMC(NUTS(posterior), num_warmup=500, num_samples=500)

`chain_init_params(model, starts, key=None, mass=None)`

Initial parameters for numpyro's MCMC, one chain per start.

Converts one dict of values per chain (e.g. the distinct fits of OrbitStart.chain_values) into numpyro's init_params: unconstrained, with a leading axis over chains.

Parameters:

Name	Type	Description	Default
`model`	`callable`	The numpyro model, e.g. from `numpyro_model`.	required
`starts`	`sequence of dict`	One dict of values per chain, keyed by sample site as `init_to_value` takes them (e.g. `FitResult.values`). Missing sites start uniformly at random.	required
`key`	`Array`	Random key for those sites (default `PRNGKey(0)`).	`None`
`mass`	`dict`	The result of `gauss_newton_mass`; a `ValueError` is raised if it was computed with another `flat_coordinates` than `model`.	`None`

Returns:

Type	Description
`dict`	Unconstrained parameters by site, with a leading axis of `len(starts)` (none for a single start), for `mcmc.run(key, init_params=...)` with `num_chains=len(starts)`.

Examples:

>>> posterior = numpyro_model(model, priors, data)
>>> init = chain_init_params(posterior, start.chain_values(4))
>>> mcmc = MCMC(NUTS(posterior), num_warmup=500, num_samples=500,
...             num_chains=4, chain_method="vectorized")
>>> mcmc.run(jax.random.PRNGKey(1), init_params=init)

`posterior_predictive_summary(samples, model, data, params=None)`

Mean and spread of the model observables over posterior samples.

Parameters:

Name	Type	Description	Default
`samples`	`dict[str, array - like]`	Posterior samples, as equal-length 1D arrays keyed by parameter name or path (e.g. from `mcmc.get_samples()`).	required
`model`	`SourceModel or callable`	Template model or class, as for :func:`loglike`.	required
`data`	`OIData`	Data defining the observables.	required
`params`	`list[str]`	Which keys of `samples` to use (default: all of them).	`None`

Returns:

Type	Description
`dict`	`vis_mean`, `vis_std`, `phi_mean` and `phi_std`: the mean and standard deviation over the samples of each visibility and phase observable.

`group_outlier_mixture(loglikes, outlier_loglikes)`

Log likelihood of groups that may each be outliers, the outlier rate marginalized exactly.

Group k has likelihood (1 - ε) L_k + ε B_k, with ε shared by all groups and Jeffreys prior Beta(½, ½). The integral over ε is a polynomial in cos θ, so Gauss–Chebyshev quadrature with K + 1 nodes is exact. Differentiable.

Provisional

May change in any 0.4.x release without a deprecation cycle.

Parameters:

Name	Type	Description	Default
`loglikes`	`(array, shape(K))`	`ln L_k`, each group's log likelihood under the model.	required
`outlier_loglikes`	`(array, shape(K) or scalar)`	`ln B_k`, each group's log likelihood as an outlier.	required

Returns:

Type	Description
`Array`	`ln ∫ Π_k [(1 - ε) L_k + ε B_k] π(ε) dε`.

Notes

Source: Gauss–Chebyshev quadrature, Abramowitz & Stegun (1964) 25.4.38 (references); ε prior, Jeffreys (1946), references.

`group_outlier_probabilities(loglikes, outlier_loglikes)`

Each group's posterior probability of being an outlier, with ε marginalized as in group_outlier_mixture.

Provisional

May change in any 0.4.x release without a deprecation cycle.

Returns:

Type	Description
`(Array, shape(K))`	`q_k = Σ_j W_j ε_j B_k / [(1 - ε_j) L_k + ε_j B_k]`, with `W_j` the posterior weights of the ε nodes.