mixle.stats.directional.bingham module

Bingham distribution on the 2-sphere (antipodally symmetric axial data).

The Bingham distribution models axes (unit vectors identified with their antipodes, x ~ -x) on S^2. Its density is

f(x) = c(Z)^{-1} exp( sum_i z_i * (m_i . x)^2 ),

where M = [m_1, m_2, m_3] is a 3 x 3 orthonormal orientation and z = (z_1, z_2, z_3) are the concentrations (defined only up to a common additive shift, since sum_i (m_i . x)^2 = 1; by convention the largest is taken to be 0). It is the antipodally symmetric analogue of the von Mises-Fisher / Kent laws and the standard model for undirected orientation data (crystallography, geology, principal axes). Equal concentrations give the uniform distribution on the sphere; a large gap concentrates mass on a great circle (girdle) or a pair of poles (bipolar).

Normalizer (reduced to a stable 1-D integral, verified to <1e-7 against arbitrary-precision integration over S^2):

c(Z) = 2 pi integral_0^pi exp(z_3 cos^2 t + (z_1+z_2)/2 sin^2 t) I_0((z_1-z_2)/2 sin^2 t) sin t dt.

Sampling is by exact angular-central-Gaussian-envelope rejection (Kent, Ganeiber & Mardia 2013); the orientation and concentrations are fit by maximum likelihood (the orientation is the eigenbasis of the scatter matrix; the concentrations solve the concave moment-matching problem d log c / d z_i = E[(m_i . x)^2]).

Reference: Bingham, “An antipodally symmetric distribution on the sphere”, Annals of Statistics 2 (1974); Kent, Ganeiber & Mardia, “A new method to simulate the Bingham and related distributions”, (2013); Mardia & Jupp, Directional Statistics (2000).

class BinghamDistribution(m, z, name=None, keys=None)[source]

Bases: SequenceEncodableProbabilityDistribution

Bingham distribution on S^2 with orientation m (3x3) and concentrations z (length 3).

Parameters:
density(x)[source]

Return the probability density or mass at a single observation.

Concrete default: exponentiate log_density (the abstract method subclasses must provide). Leaves with a cheaper closed form may override this.

Parameters:

x (Any)

Return type:

float

log_density(x)[source]

Return the log-density at a unit 3-vector x (the same value at x and -x).

Parameters:

x (Any)

Return type:

float

seq_log_density(x)[source]

Return vectorized log-density for a sequence-encoded (n, 3) array of unit vectors.

Parameters:

x (ndarray)

Return type:

ndarray

classmethod compute_capabilities()[source]
backend_seq_log_density(x, engine)[source]

Engine-neutral vectorized log-density for (N, 3) unit vectors.

Parameters:
Return type:

Any

sampler(seed=None)[source]

Return a sampler for drawing observations from this distribution.

Parameters:

seed (int | None)

Return type:

BinghamSampler

estimator(pseudo_count=None)[source]

Return an estimator for fitting this distribution from data.

Parameters:

pseudo_count (float | None)

Return type:

BinghamEstimator

dist_to_encoder()[source]

Return the data encoder used by this distribution for vectorized methods.

Return type:

BinghamDataEncoder

class BinghamSampler(dist, seed=None)[source]

Bases: DistributionSampler

Sample by angular-central-Gaussian-envelope rejection (Kent, Ganeiber & Mardia 2013).

Parameters:
  • dist (BinghamDistribution)

  • seed (int | None)

sample(size=None)[source]

Draw observations.

Combinator samplers (mixture/sequence/…) accept batched. With batched=True (the default) each child stream is drawn in one vectorized call instead of a per-draw Python loop – far faster. Because every child sampler owns an independent RandomState, batching consumes each stream in the same order as the loop, so the draws are identical to the legacy path. batched=False forces that legacy per-draw loop as a guaranteed- stable reference. Leaf samplers are already vectorized and ignore the flag.

Parameters:

size (int | None)

Return type:

Any

class BinghamAccumulator(name=None, keys=None)[source]

Bases: SequenceEncodableStatisticAccumulator

Accumulate (count, sum_xx (3,3)) – the scatter matrix is sufficient for the Bingham fit.

Parameters:
  • name (str | None)

  • keys (str | None)

update(x, weight, estimate)[source]
Parameters:
  • x (Any)

  • weight (float)

  • estimate (BinghamDistribution | None)

Return type:

None

initialize(x, weight, rng)[source]
Parameters:
Return type:

None

seq_update(x, weights, estimate)[source]
Parameters:
Return type:

None

seq_initialize(x, weights, rng)[source]
Parameters:
Return type:

None

combine(suff_stat)[source]
Parameters:

suff_stat (tuple[float, ndarray])

Return type:

BinghamAccumulator

value()[source]
Return type:

tuple[float, ndarray]

from_value(x)[source]
Parameters:

x (tuple[float, ndarray])

Return type:

BinghamAccumulator

key_merge(stats_dict)[source]

Pool this accumulator’s statistics into stats_dict under its merge key.

The structural default implements the common single-key pattern: store the accumulator under self.keys the first time the key is seen, else combine into the one already there. Accumulators with several named keys (e.g. an HMM’s init/trans/state keys) or a non-accumulator stats payload override this. A keys of None (the default) is a no-op.

Parameters:

stats_dict (dict[str, Any])

Return type:

None

key_replace(stats_dict)[source]

Replace this accumulator’s statistics from the pooled stats_dict entry (see key_merge).

Parameters:

stats_dict (dict[str, Any])

Return type:

None

acc_to_encoder()[source]
Return type:

BinghamDataEncoder

class BinghamAccumulatorFactory(name=None, keys=None)[source]

Bases: StatisticAccumulatorFactory

Factory for BinghamAccumulator.

Parameters:
  • name (str | None)

  • keys (str | None)

make()[source]
Return type:

BinghamAccumulator

class BinghamEstimator(name=None, keys=None)[source]

Bases: ParameterEstimator

Maximum-likelihood Bingham fit: scatter eigenbasis + concave concentration moment-matching.

Parameters:
  • name (str | None)

  • keys (str | None)

estimate(nobs, suff_stat)[source]
Parameters:
Return type:

BinghamDistribution

accumulator_factory()[source]
Return type:

BinghamAccumulatorFactory

class BinghamDataEncoder[source]

Bases: DataSequenceEncoder

Encode axial data as a normalized (n, 3) float array (sign is irrelevant: x ~ -x).

seq_encode(x)[source]

Encode the iid observation sequence x for vectorized evaluation.

Parameters:

x (Sequence[Any])

Return type:

ndarray