mixle.stats.matrix.lkj module¶

LKJ distribution over correlation matrices.

The Lewandowski-Kurowicka-Joe (LKJ) law places density f(R) = c_d(eta) * det(R)^(eta-1) on d x d correlation matrices (symmetric, unit diagonal, positive definite). The concentration eta > 0 tilts mass toward the identity: eta = 1 is uniform over correlation matrices, eta > 1 favours weak correlations (R near I), and eta < 1 favours strong ones. It is the standard prior on correlation matrices in hierarchical Bayesian models (Stan’s default), separating a covariance into scales times a correlation.

Normalizer (C-vine derivation, verified to high precision against arbitrary-precision integration over the correlation elliptope for d = 2, 3):

Z_d(eta) = prod_{k=1}^{d-1} B(eta + (d-1-k)/2, 1/2)^(d-k), c_d(eta) = 1 / Z_d(eta).

It samples exactly by the onion method (Lewandowski et al. 2009, sec. 3.2): each off-diagonal entry then has the exact marginal (r + 1)/2 ~ Beta(eta + (d-2)/2, eta + (d-2)/2). Because the density depends on R only through det(R), observations are encoded as log det(R), and eta is fit by maximum likelihood (a 1-D root find: the mean log-determinant equals sum_k (d-k)[psi(eta+e_k) - psi(eta+e_k+1/2)] with e_k = (d-1-k)/2).

Reference: Lewandowski, Kurowicka & Joe, “Generating random correlation matrices based on vines and extended onion method”, J. Multivariate Analysis 100 (2009).

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

Bases: SequenceEncodableProbabilityDistribution

LKJ distribution over dim x dim correlation matrices with concentration eta > 0.

Parameters:

dim (int)
eta (float)
name (str | None)
keys (str | None)

density(x)[source]

Return the probability density at a correlation matrix x.

Parameters:: x (Any)
Return type:: float

log_density(x)[source]

Return the log-density at a dim x dim correlation matrix (-inf if not positive definite).

Parameters:: x (Any)
Return type:: float

seq_log_density(x)[source]

Return vectorized log-density for a sequence-encoded array of log det(R) values.

Parameters:: x (ndarray)
Return type:: ndarray

sampler(seed=None)[source]

Return an onion-method sampler for correlation matrices.

Parameters:: seed (int | None)
Return type:: LKJSampler

estimator(pseudo_count=None)[source]

Return a maximum-likelihood estimator for the concentration eta (dim fixed).

Parameters:: pseudo_count (float | None)
Return type:: LKJEstimator

dist_to_encoder()[source]

Return the data encoder (a correlation matrix is encoded as its log-determinant).

Return type:: LKJDataEncoder

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

Bases: DistributionSampler

Sample correlation matrices by the onion method (Lewandowski-Kurowicka-Joe 2009).

Parameters:

dist (LKJDistribution)
seed (int | None)

sample(size=None)[source]

Draw one correlation matrix or a list of independent correlation matrices.

Parameters:: size (int | None)
Return type:: Any

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

Bases: SequenceEncodableStatisticAccumulator

Accumulate (count, sum of log det(R)) – the sufficient statistics for the eta-MLE.

Parameters:

name (str | None)
keys (str | None)

update(x, weight, estimate)[source]

Accumulate the weighted log determinant for one correlation matrix.

Parameters:

x (Any)
weight (float)
estimate (LKJDistribution | None)

Return type:

None

initialize(x, weight, rng)[source]

Initialize the sufficient statistics with one weighted matrix.

Parameters:

x (Any)
weight (float)
rng (RandomState | None)

Return type:

None

seq_update(x, weights, estimate)[source]

Accumulate weighted log determinants from encoded matrices.

Parameters:

x (ndarray)
weights (ndarray)
estimate (Any)

Return type:

None

seq_initialize(x, weights, rng)[source]

Initialize the sufficient statistics from encoded log determinants.

Parameters:

x (ndarray)
weights (ndarray)
rng (RandomState | None)

Return type:

None

combine(suff_stat)[source]

Merge serialized LKJ sufficient statistics into this accumulator.

Parameters:: suff_stat (tuple[float, float])
Return type:: LKJAccumulator

value()[source]

Return the total weight and weighted sum of log determinants.

Return type:: tuple[float, float]

from_value(x)[source]

Restore the accumulator from serialized LKJ sufficient statistics.

Parameters:: x (tuple[float, float])
Return type:: LKJAccumulator

key_merge(stats_dict)[source]

Merge this accumulator into a keyed statistics dictionary.

Parameters:: stats_dict (dict[str, Any])
Return type:: None

key_replace(stats_dict)[source]

Replace this accumulator from a keyed statistics dictionary.

Parameters:: stats_dict (dict[str, Any])
Return type:: None

acc_to_encoder()[source]

Return an encoder that reduces correlation matrices to log determinants.

Return type:: LKJDataEncoder

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

Bases: StatisticAccumulatorFactory

Factory for LKJAccumulator.

Parameters:

name (str | None)
keys (str | None)

make()[source]

Create an empty LKJ accumulator.

Return type:: LKJAccumulator

class LKJEstimator(dim, eta_bounds=(0.05, 1.0e4), name=None, keys=None)[source]

Bases: ParameterEstimator

Maximum-likelihood estimator for the concentration eta at fixed dimension dim.

Parameters:

dim (int)
eta_bounds (tuple[float, float])
name (str | None)
keys (str | None)

accumulator_factory()[source]

Return a factory for LKJ sufficient-statistic accumulators.

Return type:: LKJAccumulatorFactory

estimate(nobs, suff_stat)[source]

Estimate the LKJ concentration from the mean log determinant.

Parameters:

nobs (float | None)
suff_stat (tuple[float, float])

Return type:

LKJDistribution

class LKJDataEncoder[source]

Bases: DataSequenceEncoder

Encode each correlation matrix as its log-determinant (the only data the density depends on).

seq_encode(x)[source]

Encode correlation matrices as their log determinants.

Parameters:: x (Sequence[Any])
Return type:: ndarray