tensiometer.chains_convergence¶
This file contains some functions to study convergence of the chains and to compare the two posteriors.

tensiometer.chains_convergence.
GR_test
(chains, param_names=None)[source]¶ Function performing the Gelman Rubin (GR) test (described in Gelman and Rubin 92 and Brooks and Gelman 98) on a list of
MCSamples
or on a singleMCSamples
with different subchains. This test compares the variation of the mean across a pool of chains with the expected variation of the mean under the pdf that is being sampled. If we define the covariance of the mean as:\[C_{ij} \equiv {\rm Cov}_c({\rm Mean}_s(\theta))_{ij}\]and the mean covariance as:
\[M_{ij} \equiv {\rm Mean}_c[{\rm Cov}_s(\theta)_{ij}]\]then we seek to maximize:
\[R1 = {\rm max_{\theta}}\frac{C_{ij} \theta^i \theta^j} {M_{ij}\theta^i \theta^j}\]where the subscript \(c\) means that the statistics is computed across chains while the subscrit \(s\) indicates that it is computed across samples. In this case the maximization is solved by finding the maximum eigenvalue of \(C M^{1}\).
Parameters:  chains – single or list of
MCSamples
 param_names – names of the parameters involved in the test. By default uses all nonderived parameters.
Returns: value of the GR test and corresponding parameter combination
 chains – single or list of

tensiometer.chains_convergence.
GR_test_from_samples
(samples, weights)[source]¶ Lower level function to perform the Gelman Rubin (GR) test. This works on a list of samples from different chains and corresponding weights. Refer to
tensiometer.chains_convergence.GR_test()
for more details of what this function is doing.Parameters:  samples – list of samples from different chains
 weights – weights of the samples for each chain
Returns: value of the GR test and corresponding parameter combination

tensiometer.chains_convergence.
GRn_test
(chains, n, theta0=None, param_names=None, feedback=0, optimizer='ParticleSwarm', **kwargs)[source]¶ Multi dimensional higher order moments convergence test. Compares the variation of a given moment among the population of chains with the expected variation of that quantity from the samples pdf.
We first build the \(k\) order tensor of parameter differences around a point \(\tilde{\theta}\):
\[Q^{(k)} \equiv Q_{i_1, \dots, i_k} \equiv (\theta_{i_1} \tilde{\theta}_{i_1}) \cdots (\theta_{i_k} \tilde{\theta}_{i_k})\]then we build the tensor encoding its covariance across chains
\[V_M = {\rm Var}_c (E_s [Q^{(k)}])\]which is a \(2k\) and then build the second tensor encoding the mean in chain moment:
\[M_V = {\rm Mean}_c (E_s[(Q^{(k)}E_s[Q^{(k)}]) \otimes(Q^{(k)}E_s[Q^{(k)}])])\]where we have suppressed all indexes to not crowd the notation.
Then we maximize over parameters:
\[R_n 1 \equiv {\rm max}_\theta \frac{V_M \theta^{2k}}{M_V \theta^{2k}}\]where \(\theta^{2k}\) is the tensor product of \(\theta\) for \(2k\) times.
Differently from the 2D case this problem has no solution in terms of eigenvalues of tensors so far and the solution is obtained by numerical minimization with the pymanopt library.
Parameters:  chains – single or list of
MCSamples
 n – order of the moment
 theta0 – center of the moments. By default equal to the mean
 param_names – names of the parameters involved in the test. By default uses all nonderived parameters.
 feedback – level of feedback. 0=no feedback, >0 increasingly chatty
 optimizer – choice of optimization algorithm for pymanopt. Default is ParticleSwarm, other possibility is TrustRegions.
 kwargs – keyword arguments for the optimizer.
Returns: value of the GR moment test and corresponding parameter combination
 chains – single or list of

tensiometer.chains_convergence.
GRn_test_1D
(chains, n, param_name, theta0=None)[source]¶ One dimensional higher moments test. Compares the variation of a given moment among the population of chains with the expected variation of that quantity from the samples pdf.
This test is defined by:
\[R_n(\theta_0)1 = \frac{{\rm Var}_c ({\rm Mean}_s(\theta\theta_0)^n)}{{\rm Mean}_c ({\rm Var}_s(\theta\theta_0)^n) }\]where the subscript \(c\) means that the statistics is computed across chains while the subscrit \(s\) indicates that it is computed across samples.
Parameters:  chains – single or list of
MCSamples
 n – order of the moment
 param_name – names of the parameter involved in the test.
 theta0 – center of the moments. By default equal to the mean.
Returns: value of the GR moment test and corresponding parameter combination (an array with one since this works in 1D)
 chains – single or list of

tensiometer.chains_convergence.
GRn_test_1D_samples
(samples, weights, n, theta0=None)[source]¶ Lower level function to compute the one dimensional higher moments test. This works on a list of samples from different chains and corresponding weights. Refer to
tensiometer.chains_convergence.GRn_test_1D()
for more details of what this function is doing.Parameters:  samples – list of samples from different chains
 weights – weights of the samples for each chain
 n – order of the moment
 theta0 – center of the moments. By default equal to the mean.
Returns: value of the GR moment test and corresponding parameter combination (an array with one since this works in 1D)

tensiometer.chains_convergence.
GRn_test_from_samples
(samples, weights, n, theta0=None, feedback=0, optimizer='ParticleSwarm', **kwargs)[source]¶ Lower level function to compute the multi dimensional higher moments test. This works on a list of samples from different chains and corresponding weights. Refer to
tensiometer.chains_convergence.GRn_test()
for more details of what this function is doing.Parameters:  samples – list of samples from different chains
 weights – weights of the samples for each chain
 n – order of the moment
 theta0 – center of the moments. By default equal to the mean
 feedback – level of feedback. 0=no feedback, >0 increasingly chatty
 optimizer – choice of optimization algorithm for pymanopt. Default is ParticleSwarm, other possibility is TrustRegions.
 kwargs – keyword arguments for the optimizer.
Returns: value of the GR moment test and corresponding parameter combination