mixle.stats.directional.watson module

Watson distribution – a rotationally symmetric distribution for axial data on the sphere.

Axial data are unit vectors identified with their antipodes (x and -x are the same), e.g. fibre or crystal orientations, where the von Mises-Fisher (which distinguishes x from -x) does not apply. The Watson distribution on S^{p-1} concentrates around an axis mu with shape kappa:

f(x; mu, kappa) = M(1/2, p/2, kappa)^{-1} / omega_p * exp(kappa (mu^T x)^2),

where M is Kummer’s confluent hypergeometric function and omega_p = 2 pi^{p/2} / Gamma(p/2) is the sphere’s surface area. kappa > 0 is bipolar (mass near the axis +/-mu), kappa < 0 girdle (mass on the equator orthogonal to mu); both are antipodally symmetric. It is fit by maximum likelihood: mu is the leading (kappa>0) or trailing (kappa<0) eigenvector of the scatter matrix, and kappa solves E[(mu^T x)^2] = mu^T S mu (a monotone 1-D equation in the Kummer ratio).

Reference: Mardia & Jupp, Directional Statistics (Wiley, 2000).

class WatsonDistribution(mu, kappa, name=None, keys=None)[source]

Bases: SequenceEncodableProbabilityDistribution

Watson distribution on the unit sphere S^{p-1} with axis mu and concentration kappa.

Parameters:
density(x)[source]

Return the density at a single unit vector x.

Parameters:

x (ndarray)

Return type:

float

log_density(x)[source]

Return the log-density at a single unit vector x.

Parameters:

x (ndarray)

Return type:

float

seq_log_density(x)[source]

Vectorized log-density for a stack of unit vectors, shape (N, p).

Parameters:

x (ndarray)

Return type:

ndarray

classmethod compute_capabilities()[source]

Describe backend support for generated Watson scoring kernels.

backend_seq_log_density(x, engine)[source]

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

Parameters:
Return type:

Any

sampler(seed=None)[source]

Return a sampler for drawing unit vectors from this distribution.

Parameters:

seed (int | None)

Return type:

WatsonSampler

estimator(pseudo_count=None)[source]

Return a maximum-likelihood estimator (scatter eigenvector + Kummer-ratio kappa solve).

Parameters:

pseudo_count (float | None)

Return type:

WatsonEstimator

dist_to_encoder()[source]

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

Return type:

WatsonDataEncoder

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

Bases: DistributionSampler

Draw axes by sampling s = |mu^T x| via a numerical inverse-CDF, then a random orthogonal tangent.

Parameters:

  • dist (WatsonDistribution)

  • seed (int | None)

sample(size=None)[source]

Draw one unit vector or a stack of iid unit vectors.

Parameters:

size (int | None)

Return type:

ndarray

class WatsonAccumulator(dim, name=None, keys=None)[source]

Bases: SequenceEncodableStatisticAccumulator

Accumulate the weighted scatter matrix S = sum_i w_i x_i x_i^T and total weight.

Parameters:

  • dim (int)

  • name (str | None)

  • keys (str | None)

update(x, weight, estimate)[source]

Accumulate one weighted outer product into the scatter matrix.

Parameters:

  • x (ndarray)

  • weight (float)

  • estimate (WatsonDistribution | None)

Return type:

None

initialize(x, weight, rng)[source]

Initialize statistics from one unit vector.

Parameters:
Return type:

None

seq_update(x, weights, estimate)[source]

Accumulate weighted scatter statistics from encoded unit vectors.

Parameters:
Return type:

None

seq_initialize(x, weights, rng)[source]

Initialize statistics from encoded unit vectors.

Parameters:
Return type:

None

combine(suff_stat)[source]

Merge another Watson sufficient-statistic tuple.

Parameters:

suff_stat (tuple[ndarray, float])

Return type:

WatsonAccumulator

value()[source]

Return the scatter matrix and total weight.

Return type:

tuple[ndarray, float]

from_value(x)[source]

Replace accumulator contents from scatter statistics.

Parameters:

x (tuple[ndarray, float])

Return type:

WatsonAccumulator

key_merge(stats_dict)[source]

Merge keyed statistics into stats_dict when keys are configured.

Parameters:

stats_dict (dict[str, Any])

Return type:

None

key_replace(stats_dict)[source]

Replace this accumulator from keyed statistics when available.

Parameters:

stats_dict (dict[str, Any])

Return type:

None

acc_to_encoder()[source]

Return the encoder used by this accumulator.

Return type:

WatsonDataEncoder

class WatsonAccumulatorFactory(dim, name=None, keys=None)[source]

Bases: StatisticAccumulatorFactory

Factory for WatsonAccumulator.

Parameters:

  • dim (int)

  • name (str | None)

  • keys (str | None)

make()[source]

Create a fresh Watson accumulator.

Return type:

WatsonAccumulator

class WatsonEstimator(dim, name=None, keys=None)[source]

Bases: ParameterEstimator

Maximum-likelihood estimator: scatter eigenvector for the axis, Kummer-ratio solve for kappa.

Parameters:

  • dim (int)

  • name (str | None)

  • keys (str | None)

accumulator_factory()[source]

Return an accumulator factory for Watson scatter statistics.

Return type:

WatsonAccumulatorFactory

estimate(nobs, suff_stat)[source]

Estimate the Watson axis and concentration from weighted scatter.

Parameters:
Return type:

WatsonDistribution

class WatsonDataEncoder[source]

Bases: DataSequenceEncoder

Encode a sequence of unit vectors as an (N, p) float array.

seq_encode(x)[source]

Encode unit vectors as an (N, p) floating-point array.

Parameters:

x (Sequence[ndarray])

Return type:

ndarray