API Reference

API Reference#

Core Module#

Bayesian Matrix-Based Kinematic Deconvolution#

This module infers the intrinsic Line-of-Sight Velocity Distribution (LOSVD) from discrete, heteroscedastic stellar observations using a pre-computed linear design matrix and a hierarchical smoothness prior.

class veldist.veldist.KinematicSolver[source]#

Bases: object

A high-level interface for performing Bayesian kinematic deconvolution.

This class manages the full inference workflow: 1. Defining the velocity grid (setup_grid). 2. Ingesting data and building the design matrix (add_data). 3. Running the MCMC sampler (run). 4. Visualizing the results (plot_result).

matrix#

The pre-computed design matrix of shape (N_stars, N_bins).

Type:: jnp.ndarray or None

grid#

Metadata defining the velocity grid (centers, edges, width).

Type:: dict

n_stars#

Number of stars loaded via add_data. Used as the bin_flux analog when writing Dynamite output (see write_dynamite_kinematics).

Type:: int or None

samples#

Posterior samples from the MCMC run.

Type:: dict or None

clipped_samples#

Per-bin summary statistics (median LOSVD and clipped uncertainties) populated by clip_uncertainties.

Type:: dict or None

add_data(vel, err)[source]#

Load observations and pre-compute the design matrix.

Parameters:

vel (array-like) – Observed velocities of stars.
err (array-like) – Measurement errors associated with the velocities.

Returns:

Sets self.matrix.

Return type:

None

clip_uncertainties(floor_fraction=0.01, abs_floor=1e-10)[source]#

Apply uncertainty floors and store LOSVD summary statistics.

This is a post-processing step that does not modify the raw posterior samples in self.samples. It summarises the posterior as per-bin marginal medians and half-CI-widths in probability-mass space, then raises the uncertainties to a floor so that no bin carries a zero into the Dynamite output writer.

Output format matches the Dynamite BayesLOSVD ECSV convention (see context/dynamite_format_spec.md):

losvd_median stores the per-bin marginal median of the posterior probability mass. Because the joint posterior is a simplex but marginals are taken independently, the median values typically sum to 0.85–0.95, not 1. This is expected and correct.
losvd_uncertainty stores the half-width of the 68% credible interval: (p84 − p16) / 2. Used as symmetric ±error bars.

Both quantities are dimensionless probability mass per bin — they are not divided by the bin width.

Motivation#

Zero uncertainties in LOSVD bins propagate into Dynamite’s internal NNLS projection matrices and produce econ zeros that cause weight- solving failures in large orbit-library runs. The relative floor (floor_fraction * max_uncertainty) is the primary safeguard; the absolute floor is a numerical backstop for channels where the posterior is pathologically tight across the board.

param floor_fraction:

Relative floor as a fraction of the maximum per-bin half-CI-width across all bins. Default 0.01 (1%).

type floor_fraction:

float

param abs_floor:

Absolute floor applied after the relative floor. Default 1e-10.

type abs_floor:

float

returns:

Sets self.clipped_samples as a dict with keys:

'losvd_median' — per-bin marginal median, probability mass (dimensionless); shape (n_bins,).
'losvd_uncertainty' — clipped half-width of 68% CI, probability mass; shape (n_bins,).

rtype:

None

plot_result(ax=None, true_intrinsic=None)[source]#

Visualize the inferred LOSVD.

Parameters:

ax (matplotlib.axes.Axes, optional) – Axes to plot on. If None, creates a new figure.
true_intrinsic (array-like, optional) – True intrinsic velocities for comparison (if available).

Returns:

ax – The plot axes.

Return type:

matplotlib.axes.Axes

run(num_warmup=500, num_samples=1000, gpu=True, seed=5567)[source]#

Run the NUTS sampler.

Parameters:

num_warmup (int) – Number of warmup (burn-in) steps.
num_samples (int) – Number of MCMC samples to draw.
gpu (bool) – Whether to use GPU acceleration (if available).
seed (int) – RNG seed for the NUTS sampler. Default 5567 (kept for backwards compatibility). When running many bins in a batch, pass distinct seeds per bin to avoid any correlation in the sampling chains; a simple convention is seed + bin_index (see fit_all_bins).

Returns:

samples – Posterior samples (e.g., “intrinsic_pdf”, “smoothness_sigma”).

Return type:

dict

setup_grid(center, width, n_bins)[source]#

Define the velocity grid (histogram bins).

Parameters:

center (float) – Center of the velocity grid.
width (float) – Total width of the velocity grid.
n_bins (int) – Number of bins in the grid.

Returns:

Sets self.grid.

Return type:

None

truncate_losvd(n_sigma=3.0, abs_floor=1e-10)[source]#

Suppress LOSVD bins beyond n_sigma dispersions from the bulk mean.

Note

This method should only be called if tail contamination is observed to be causing problems (e.g. non-negligible posterior mass in channels well outside the range of the input data, or econ zeros in Dynamite that persist after clip_uncertainties). It is not part of the standard pipeline and is not called by fit_all_bins(). When in doubt, omit it — modifying the posterior summary without a clear diagnostic reason introduces unnecessary bias.

Velocity channels far from the data carry unphysical posterior weight due to the random-walk prior leaking into unconstrained tails. This method zeros the mean density in those channels and sets their uncertainty to abs_floor, preventing Dynamite from trying to fit orbits in unphysical regimes.

Approach follows the velocity-range truncation strategy discussed in Falcón-Barroso & Martig (2021, A&A).

This is a post-processing step and does not modify self.samples. It operates on self.clipped_samples, creating it first (with default floors) if it is not already present.

Parameters:

n_sigma (float) – Number of velocity dispersions beyond which to truncate. Default 3.0.
abs_floor (float) – Uncertainty value assigned to truncated bins (must be > 0 to avoid econ zeros in Dynamite). Default 1e-10.

Returns:

Updates self.clipped_samples in-place.

Return type:

None

veldist.veldist.fit_all_bins(bin_data_list, grid_kwargs, run_kwargs=None, min_stars=10)[source]#

Run the full inference pipeline for a list of Voronoi bins.

For each bin, this executes the setup_grid → add_data → run → clip_uncertainties pipeline and returns a list of KinematicSolver instances ready for the Dynamite output writer. Bins with too few stars are skipped (returning None at that position) so the writer can mask them.

truncate_losvd() is deliberately not called here. It is an optional diagnostic repair step; call it manually on individual solvers only if tail contamination is observed to be causing problems.

The same velocity grid is used for every bin (grid_kwargs is shared). Each bin receives a unique RNG seed derived as base_seed + bin_index to avoid correlations between sampling chains.

Parameters:

bin_data_list (list of dict) –
One dict per Voronoi bin. Required keys:
- 'vel' — array of observed stellar velocities.
- 'err' — array of per-star measurement errors.
Any additional keys (e.g. spatial metadata) are ignored here and can be passed separately to the output writer.
grid_kwargs (dict) – Keyword arguments forwarded to KinematicSolver.setup_grid() (center, width, n_bins). Shared across all bins.
run_kwargs (dict, optional) – Keyword arguments forwarded to KinematicSolver.run() (e.g. num_warmup, num_samples, gpu). The seed key, if present, is used as the base seed; each bin then receives seed + bin_index. Defaults to {} (all run defaults apply).
min_stars (int) – Minimum number of stars required to attempt inference. Bins with fewer stars are skipped with a warning. Default 10.

Returns:

solvers – One entry per input bin. Entries are either a fully solved KinematicSolver (with samples and clipped_samples populated) or None for skipped bins.

Return type:

list

Examples

>>> solvers = fit_all_bins(
...     bin_data_list,
...     grid_kwargs={"center": 0.0, "width": 600.0, "n_bins": 60},
...     run_kwargs={"num_warmup": 500, "num_samples": 1000, "gpu": False},
... )
>>> # Pass to the output writer (Task 2):
>>> solver.write_dynamite_kinematics(output_dir, voronoi_bin_metadata)

veldist.veldist.generate_smooth_curve(N_bins, smoothness_sigma)[source]#

Generates a 1D smooth curve using a Gaussian Random Walk.

Model:: step[i] ~ Normal(0, sigma) curve[i] = curve[i-1] + step[i]

Parameters:

N_bins (int) – Number of velocity bins.
smoothness_sigma (float) – Controls the flexibility of the curve. - Low sigma: Stiff, very smooth. - High sigma: Flexible, jagged.

Returns:

curve – Latent log-density curve.

Return type:

jnp.ndarray (N_bins,)

veldist.veldist.model(matrix, n_bins)[source]#

The Model definition for NumPyro.

Parameters:

matrix (jnp.ndarray (N_stars, N_bins)) – The pre-computed Design Matrix M.
n_bins (int) – Number of velocity bins.

Returns:

This function defines the probabilistic graph and has no return value.

Return type:

None

veldist.veldist.precompute_design_matrix(obs_val, obs_err, bin_centers, bin_width=None)[source]#

Computes the Probability Design Matrix (M) for 1D LOSVD inference.

This function bakes the observations and uncertainties into a static matrix, converting the deconvolution problem into a single matrix multiplication during inference.

Parameters:

obs_val (array-like (N,)) – Observed velocities of individual stars.
obs_err (array-like (N,)) – Standard deviation (measurement error) for each star.
bin_centers (array-like (K,)) – The velocity grid centers (the intrinsic histogram bins).
bin_width (float, optional) – Width of bins. If None, inferred from centers.

Returns:

M – Probability matrix. M[i, j] = P(Star i | True Velocity is in Bin j)

Return type:

jnp.ndarray (N, K)

veldist.veldist.write_dynamite_kinematics(solvers, output_dir, voronoi_bin_metadata, kin_filename='bayes_losvd_kins.ecsv', aperture_filename='aperture.dat', bins_filename='bins.dat', bin_flux_mode='nstars')[source]#

Write Dynamite-compatible BayesLOSVD input files from a list of solved bins.

Produces three files that Dynamite expects for its BayesLOSVD kinematics representation (see context/dynamite_format_spec.md for full format details):

{kin_filename} — Astropy ECSV, one row per solved Voronoi bin, containing the per-bin marginal median LOSVD and ±half-CI uncertainties.
{aperture_filename} — pixel grid geometry.
{bins_filename} — pixel-to-bin mapping.

Any None entries in solvers (bins skipped by fit_all_bins()) are automatically masked: their pixels are written as 0 in the bins file and they are omitted from the kinematics table. The remaining bins are re-numbered sequentially (1-indexed) as required by Dynamite.

clip_uncertainties() is called automatically on any solver that has not already had its clipped_samples populated.

Parameters:

solvers (list) – Solved KinematicSolver instances (or None for skipped bins), as returned by fit_all_bins(). All non-None entries must share the same velocity grid.
output_dir (str or path-like) – Directory in which to write the three output files. Created if it does not exist.

voronoi_bin_metadata (dict) –

Spatial and observational metadata. Required structure:

{
    'bins': [
        {
            'xbin': float,       # arcsec, bin centre x
            'ybin': float,       # arcsec, bin centre y
            # 'bin_flux' only read when bin_flux_mode='custom'
            'bin_flux': float,
        },
        ...  # one entry per entry in solvers
    ],
    'aperture': {
        'x_start':   float,   # arcsec, lower-left corner x
        'y_start':   float,   # arcsec, lower-left corner y
        'x_size':    float,   # arcsec, total x extent
        'y_size':    float,   # arcsec, total y extent
        'angle_deg': float,   # degrees (= 90 - position_angle)
        'nx':        int,     # pixels along x
        'ny':        int,     # pixels along y
    },
    # 1-indexed bin IDs (bin 1 = solvers[0], bin 2 = solvers[1], …)
    # 0 = masked.  Skipped bins are re-mapped to 0 automatically.
    'pixel_bin_ids': array-like,   # shape (nx*ny,) or (ny, nx)
    'psf': {
        'sigma':  [float, ...],   # Gaussian sigma(s) in arcsec
        'weight': [float, ...],   # must sum to 1
    },
}

kin_filename (str) – File name for the kinematics ECSV. Default 'bayes_losvd_kins.ecsv'.
aperture_filename (str) – File name for the aperture file. Default 'aperture.dat'.
bins_filename (str) – File name for the bins file. Default 'bins.dat'.
bin_flux_mode ({'nstars', 'uniform', 'custom'}) –
Controls what is written to the bin_flux column.

'nstars' (default)
Use the number of stars in each bin (solver.n_stars, set by add_data()). This is the physically meaningful analog of IFU surface brightness for discrete stellar kinematic data. bin_flux is used by Dynamite only for flux-weighted systemic velocity centering (center_v_systemic); it does not enter the NNLS chi² (confirmed from NNLS.construct_nnls_matrix_and_rhs and BayesLOSVD.get_observed_values_and_uncertainties in Dynamite’s source). N_stars is the physically appropriate quantity to pass for discrete stellar data.

'uniform'
Write 1.0 for every bin. Use this if you want equal weighting of all bins in any flux-weighted calculation, or if you are uncertain about the right quantity and want a neutral default.

'custom'
Read bin_flux from voronoi_bin_metadata['bins'][i]['bin_flux'] for each bin. Useful if you have an external flux estimate (e.g. sum of photometric counts in each Voronoi cell).

Raises:

ValueError – If no solved bins are found, or if solvers share inconsistent grids.
AssertionError – If any LOSVD uncertainty is ≤ 0 after clipping (would cause econ zeros in Dynamite).

Return type:

None

Analysis Module#

Statistical analysis utilities for inferred velocity distributions.

veldist.analysis.bimodality_score(pdf_samples)[source]#

Count the number of peaks in the smoothed posterior-mean LOSVD.

This is a diagnostic integer flag, not a posterior quantity. It is computed from the posterior-mean LOSVD (not from individual samples), so no credible interval is available. Bins with a score of 2 or more indicate a potentially multimodal distribution where GH-analogue summary statistics (mean, sigma, skewness, kurtosis) may be misleading and the full histogram shape should be inspected.

Peak detection uses a 3-point boxcar to smooth over single-bin noise, and requires each candidate peak to exceed 1% of the global maximum to suppress spurious detections in low-probability tails.

Parameters:

pdf_samples (array-like, shape (n_samples, n_bins)) – MCMC samples of the probability mass function. Each row sums to 1.

Returns:

Number of local maxima detected. Typical values:

1: Unimodal distribution (normal case).
2: Bimodal — could indicate counter-rotation, two kinematic components, or a contaminating population.
≥ 3: Highly irregular; inspect visually before interpreting any scalar summaries.

Return type:

int

Notes

A more principled treatment would count peaks on every individual MCMC sample and return a distribution over the count. This is deferred as a future improvement; the posterior-mean approach is sufficient for identifying bins that require visual inspection.

Examples

>>> score = bimodality_score(solver.samples["intrinsic_pdf"])
>>> if score >= 2:
...     print("Multimodal — inspect full histogram before trusting moments")

veldist.analysis.cdf_percentile(pdf_samples, grid_centers, p)[source]#

Interpolate velocity values at cumulative probability level(s) p.

Builds the empirical CDF for each MCMC sample by cumulative summation of the probability mass, then uses linear interpolation to find the velocity at which the CDF equals p. This is the building block for compute_summary()’s median and IQR calculations, but is also useful on its own for constructing credible bands on cumulative profiles.

Parameters:

pdf_samples (array-like, shape (n_samples, n_bins)) – MCMC samples of the probability mass function. Each row must sum to 1 and grid_centers must be monotonically increasing.
grid_centers (array-like, shape (n_bins,)) – Bin-centre velocities in ascending order.
p (float or array-like) – Cumulative probability level(s) in the interval [0, 1]. Scalar input returns a 1-D result; array input returns a 2-D result.

Returns:

Shape (n_samples,) when p is a scalar, or (n_samples, len(p)) when p is array-like.

Return type:

ndarray

Examples

Compute the posterior median velocity for each MCMC sample:

>>> v_median_samples = cdf_percentile(pdf_samples, grid_centers, 0.5)

Compute Q25 and Q75 simultaneously:

>>> q25_q75 = cdf_percentile(pdf_samples, grid_centers, [0.25, 0.75])
>>> iqr_samples = q25_q75[:, 1] - q25_q75[:, 0]

veldist.analysis.compute_moments(pdf_samples, grid_centers)[source]#

Compute statistical moments from posterior LOSVD samples.

Deprecated since version This: function is retained for backward compatibility. For new code, prefer compute_summary(), which returns the same quantities plus robust alternatives (median, IQR, tail weight, bimodality score) and uses the posterior half-68% CI as its uncertainty convention instead of the posterior standard deviation.

Parameters:

pdf_samples (array-like, shape (n_samples, n_bins)) – MCMC samples of the probability mass function. Each row must sum to 1.
grid_centers (array-like, shape (n_bins,)) – Centres of the velocity bins (km/s or consistent velocity unit).

Returns:

Keys and values:

'mean'(float, float): (posterior_mean, posterior_std) of the flux-weighted mean velocity across MCMC samples.
'std'(float, float): (posterior_mean, posterior_std) of the velocity dispersion.
'skewness'(float, float): (posterior_mean, posterior_std) of the normalised third central moment.
'kurtosis'(float, float): (posterior_mean, posterior_std) of the excess kurtosis (fourth central moment / σ⁴ − 3).

Return type:

dict

veldist.analysis.compute_summary(pdf_samples, grid_centers)[source]#

Compute spatially-mappable scalar summaries from posterior LOSVD samples.

This is the primary analysis function for extracting kinematic maps from KinematicSolver output. All metrics except bimodality_score are evaluated independently on every MCMC sample, so the full posterior uncertainty — measurement noise, finite star count, and prior regularisation — is propagated automatically with no separate bootstrap or error-propagation step.

Each metric is summarised as (posterior_median, half_68ci) using the same convention as the LOSVD itself: the reported uncertainty is (p84 − p16) / 2.

Parameters:

pdf_samples (array-like, shape (n_samples, n_bins)) – MCMC samples of the probability mass function. Each row must sum to 1.
grid_centers (array-like, shape (n_bins,)) – Centres of the velocity bins (km/s or consistent velocity unit).

Returns:

Each key maps to a (median, half_68ci) tuple of floats — both in the same units as grid_centers for velocity quantities, and dimensionless for shape metrics — except 'bimodality_score', which is a plain int.

Location

'v_mean': Flux-weighted mean velocity. GH analogue: V. Sensitive to tail contamination; compare with v_median as a cross-check.
'v_median': Median velocity (CDF = 0.5). Robust against edge-bin contamination and heavy tails.
'v_asymmetry': Mean minus median. Near zero for symmetric LOSVDs; the sign mirrors that of the low-velocity tail (positive = mean pulled toward higher velocities by a trailing tail). Closely related to h₃, but does not require computing higher-order moments.

Dispersion

'sigma': Standard deviation of the LOSVD. GH analogue: σ.
'iqr': Interquartile range Q75 − Q25. Robust dispersion estimate insensitive to tail contamination.
'sigma_iqr': IQR / 1.3490 — the Gaussian-equivalent dispersion derived from the IQR. For a Gaussian, sigma_iqr ≈ sigma. sigma_iqr < sigma → heavy tails (radial anisotropy); sigma_iqr > sigma → flat top (tangential anisotropy).

Shape

'skewness': Normalised third central moment γ₁. Zero for symmetric distributions. GH analogue: h₃ ≈ −*γ₁* / √6. Note sign: γ₁ > 0 (right-skewed, trailing tail) → h₃ < 0, the expected signal on the receding side of a rotating system.
'kurtosis': Excess kurtosis κ = fourth central moment / σ⁴ − 3. Zero for a Gaussian. GH analogue: h₄ ≈ κ / √24. Positive (leptokurtic) → radially anisotropic; negative (platykurtic) → tangentially anisotropic / flat-topped.
'tail_weight': Fraction of probability mass outside ±1*σ* of the mean. Gaussian reference: 0.3173. A more direct anisotropy diagnostic than h₄: no expansion assumption, always interpretable even for non-Gaussian shapes. See also tail_weight().

Diagnostic

'bimodality_score': Integer number of peaks in the smoothed posterior-mean LOSVD (see bimodality_score()). Score 1 = unimodal; ≥ 2 = inspect visually. No uncertainty is returned for this metric.

Return type:

dict

Notes

The approximate Gauss-Hermite conversions (valid for |h₃|, |h₄| ≲ 0.2):

h3 ≈ -skewness / sqrt(6)
h4 ≈  kurtosis / sqrt(24)

These allow cross-validation against GH-based models and literature maps.

For spatial bins where bimodality_score ≥ 2, the mean, sigma, skewness, and kurtosis should be treated with caution: the mean lands between two peaks, sigma is inflated by their separation, and skewness reflects which peak is taller rather than any genuine asymmetry.

Examples

>>> solver = KinematicSolver(v_grid)
>>> solver.add_data(velocities, uncertainties)
>>> solver.run()
>>> summary = compute_summary(solver.samples["intrinsic_pdf"],
...                           solver.grid["centers"])
>>> v, dv = summary["v_mean"]
>>> s, ds = summary["sigma"]
>>> print(f"V = {v:.1f} ± {dv:.1f}  σ = {s:.1f} ± {ds:.1f}  km/s")

veldist.analysis.compute_summary_maps(solvers)[source]#

Compute summary statistics for all solved bins from fit_all_bins().

Iterates over a list of KinematicSolver instances, calls compute_summary() on each one, and assembles the results into arrays shaped (n_bins,) — one entry per spatial bin — ready to pass directly to a spatial map plotting function.

Bins that were skipped during inference (None entries in solvers, as returned for bins below min_stars by fit_all_bins()) produce NaN in all output arrays so that the spatial indexing is preserved.

Parameters:

solvers (list of KinematicSolver or None) – As returned by fit_all_bins(). None entries (skipped bins) are silently mapped to NaN. The list must contain at least one non-None entry.

Returns:

One key per metric name from compute_summary(). Each value is a sub-dict with two keys:

'median'ndarray, shape (n_bins,): Posterior median of the metric. NaN for skipped bins. For bimodality_score, the integer score is cast to float.
'uncertainty'ndarray, shape (n_bins,): Half-width of the 68% credible interval. NaN for skipped bins and for bimodality_score (which has no posterior CI).

Return type:

dict

Raises:

ValueError – If all entries in solvers are None.

Examples

Build velocity and dispersion maps and plot them side by side:

>>> maps = compute_summary_maps(solvers)
>>> import matplotlib.pyplot as plt
>>> fig, axes = plt.subplots(1, 2, figsize=(10, 4))
>>> for ax, key, label in zip(axes,
...                           ["v_mean", "sigma"],
...                           ["Mean velocity (km/s)", "Dispersion (km/s)"]):
...     sc = ax.scatter(xbin, ybin, c=maps[key]["median"], cmap="RdBu_r")
...     plt.colorbar(sc, ax=ax, label=label)

Flag multimodal bins before plotting moments:

>>> bscore = maps["bimodality_score"]["median"]
>>> multimodal = bscore >= 2
>>> print(f"{multimodal.sum():.0f} bins have bimodality_score ≥ 2")

Notes

The bimodality_score sub-dict has 'median' filled with the integer score cast to float and 'uncertainty' filled with NaN. This keeps the return structure uniform so callers can iterate over all keys without special-casing the diagnostic metric.

veldist.analysis.half_68ci(samples)[source]#

Half-width of the 68% posterior credible interval.

Computes (p84 − p16) / 2 — the symmetric ± uncertainty reported throughout veldist for both LOSVD values and scalar summary metrics. This matches the convention used by BayesLOSVD and is the natural non-parametric analogue of a 1-σ Gaussian error bar.

Parameters:: samples (array-like, shape (n_samples,)) – Posterior samples of a scalar quantity (e.g., one metric evaluated over all MCMC draws).
Returns:: (p84 − p16) / 2, in the same units as samples.
Return type:: float

Notes

For a Gaussian posterior this equals the posterior standard deviation. For skewed or heavy-tailed posteriors it can differ substantially, but it always has the interpretation “the true value lies within ±half_68ci of the median with approximately 68% posterior probability.”

Examples

>>> v_mean_samples = pdf_samples @ grid_centers
>>> uncertainty = half_68ci(v_mean_samples)
>>> median = float(np.median(v_mean_samples))
>>> print(f"v_mean = {median:.1f} ± {uncertainty:.1f} km/s")

veldist.analysis.tail_weight(pdf_samples, grid_centers, means, stds)[source]#

Fraction of probability mass outside ±1σ of the mean, per MCMC sample.

This is a direct, model-free measure of tail heaviness that can be interpreted without any expansion assumption. It is the non-parametric analogue of the Gauss-Hermite h₄ coefficient, but more intuitive: for a Gaussian the value is exactly 1 − erf(1/√2) ≈ 0.3173.

Parameters:

pdf_samples (array-like, shape (n_samples, n_bins)) – MCMC samples of the probability mass function. Each row sums to 1.
grid_centers (array-like, shape (n_bins,)) – Bin-centre velocities.
means (array-like, shape (n_samples,)) – Per-sample flux-weighted mean velocities. Compute as pdf_samples @ grid_centers if not already available.
stds (array-like, shape (n_samples,)) – Per-sample velocity dispersions (standard deviations of the LOSVD).

Returns:

Tail-weight value for each MCMC sample. Values greater than 0.317 indicate heavier tails than Gaussian (associated with radial velocity anisotropy); values below 0.317 indicate lighter tails / flat top (tangential anisotropy).

Return type:

ndarray, shape (n_samples,)

Notes

The Gaussian reference value can be computed as:

from math import erf, sqrt
gaussian_tail_weight = 1 - erf(1 / sqrt(2))   # ≈ 0.3173

Examples

>>> means = pdf_samples @ grid_centers
>>> stds  = np.sqrt(np.einsum('ij,ij->i', pdf_samples,
...                           (grid_centers - means[:, None])**2))
>>> tw_samples = tail_weight(pdf_samples, grid_centers, means, stds)
>>> print(f"tail weight = {np.median(tw_samples):.4f}")

API Reference

Contents

API Reference#

Core Module#

Bayesian Matrix-Based Kinematic Deconvolution#

Motivation#

Analysis Module#