bayesmix/algorithms

The Algorithm class handles other class objects and performs the MCMC simulation. There are two types of Algorithm: marginal and conditional, each of which can only be used with the matching type of Mixing.

Algorithms¶

class BaseAlgorithm¶

Abstract template class representing a Gibbs sampler for mixture models. Gibbs samplers are a particular class of Markov chain Monte Carlo (MCMC) algorithms that are used to perform posterior inference in Bayesian models. For the specific class of models that can be fit using these algorithm, see MarginalAlgorithm and ConditionalAlgorithm.

This class receives a number of initialization parameters related to the MCMC algorithm itself (such as the number of iterations) in a Protobuf message. It also utilizes several class objects related to the underlying model: a Mixing object, and a vector of Hierarchy objects. It includes methods for running the MCMC simulation and for estimating the posterior density on a given grid of points. In particular, density estimation is performed by averaging the estimate associated to each single MCMC iteration. The specific method depends on whether the algorithm is marginal or conditional (see derived classes MarginalAlgorithm and ConditionalAlgorithm). Results for this class’ methods are saved either to Collector objects, or to text files.

Subclassed by ConditionalAlgorithm, MarginalAlgorithm

Public Functions

virtual bool is_conditional() const = 0¶: Returns whether the algorithm is conditional or marginal.

inline virtual bool requires_conjugate_hierarchy() const¶: Returns whether it is restricted to conjugate Hierarchy objects.

inline void run(BaseCollector *collector)¶: Runs the algorithm and saves the MCMC chain to the given Collector

virtual Eigen::MatrixXd eval_lpdf(BaseCollector *const collector, const Eigen::MatrixXd &grid, const Eigen::RowVectorXd &hier_covariate = Eigen::RowVectorXd(0), const Eigen::RowVectorXd &mix_covariate = Eigen::RowVectorXd(0))¶

Evaluates the estimated log-pdf on all iterations over a grid of points

Parameters:

collector – Collector object from which the MCMC is read
grid – Grid of row points on which the density is to be evaluated
hier_covariate – Optional covariates related to the Hierarchy
mix_covariate – Optional covariates related to the Mixing

Returns:

The estimation over all iterations (rows) and points (cols)

virtual bayesmix::AlgorithmId get_id() const = 0¶: Returns the Protobuf ID associated to this class.

virtual void read_params_from_proto(const bayesmix::AlgorithmParams &params)¶: Reads and sets algorithm parameters from an appropriate Protobuf message.

virtual Eigen::VectorXd lpdf_from_state(const Eigen::MatrixXd &grid, const Eigen::RowVectorXd &hier_covariate, const Eigen::RowVectorXd &mix_covariate) = 0¶

Evaluates the estimated log-lpdf on the state contained in curr_state

Parameters:

grid – Grid of row points on which the density is to be evaluated
hier_covariate – Optional covariates related to the Hierarchy
mix_covariate – Optional covariates related to the Mixing

Returns:

The estimation on the iteration of curr_state over all points

class MarginalAlgorithm : public BaseAlgorithm ¶

Template class for a marginal sampler deriving from BaseAlgorithm.

This template class implements a generic Gibbs sampling marginal algorithm as the child of the BaseAlgorithm class. A mixture model sampled from a Marginal Algorithm can be expressed as

$x_i \mid c_i, \theta_1, \dots, \theta_k &\sim f(x_i \mid \theta_{c_i}) \\ \theta_1, \dots, \theta_k &\sim G_0 \\ c_1, \dots, c_n &\sim EPPF(c_1, \dots, c_n)$

where $f(x \mid \theta_j)$ is a density for each value of $\theta_j$ and $c_i$ take values in ${1, \dots, k}$ . Depending on the actual implementation, the algorithm might require the kernel/likelihood $f(x \mid \theta)$ and $G_0(\phi)$ to be conjugagte or not. In the former case, a conjugate hierarchy must be specified. In this library, each $\theta_j$ is represented as an Hierarchy object (which inherits from AbstractHierarchy), that also holds the information related to the base measure $G_0$ is (see AbstractHierarchy). The $EPPF$ is instead represented as a Mixing object, which inherits from AbstractMixing.

The state of a marginal algorithm only consists of allocations and unique values. In this class of algorithms, the local lpdf estimate for a single iteration is a weighted average of likelihood values corresponding to each component (i.e. cluster), with the weights being based on its cardinality, and of the marginal component, which depends on the specific algorithm.

Subclassed by Neal2Algorithm, SplitAndMergeAlgorithm

Public Functions

inline virtual bool is_conditional() const override¶: Returns whether the algorithm is conditional or marginal.

virtual Eigen::VectorXd lpdf_from_state(const Eigen::MatrixXd &grid, const Eigen::RowVectorXd &hier_covariate, const Eigen::RowVectorXd &mix_covariate) override¶

Evaluates the estimated log-lpdf on the state contained in curr_state

Parameters:

grid – Grid of row points on which the density is to be evaluated
hier_covariate – Optional covariates related to the Hierarchy
mix_covariate – Optional covariates related to the Mixing

Returns:

The estimation on the iteration of curr_state over all points

class Neal2Algorithm : public MarginalAlgorithm ¶

Template class for Neal’s algorithm 2 for conjugate hierarchies.

This class implements Neal’s Gibbs sampling algorithm 2 from Neal (2000) that generates a Markov chain on the clustering of the provided data.

This algorithm requires the use of a ConjugateHierarchy object.

Subclassed by Neal3Algorithm, Neal8Algorithm

Public Functions

inline virtual bool requires_conjugate_hierarchy() const override¶: Returns whether it is restricted to conjugate Hierarchy objects.

inline virtual bayesmix::AlgorithmId get_id() const override¶: Returns the Protobuf ID associated to this class.

class Neal3Algorithm : public Neal2Algorithm ¶

Template class for Neal’s algorithm 3 for conjugate hierarchies.

This class implements Neal’s Gibbs sampling algorithm 3 from Neal (2000) that generates a Markov chain on the clustering of the provided data.

This algorithm requires the use of a ConjugateHierarchy object. Algorithm 3 from Neal (2000) is almost identical to Algorithm 2, except that conjugacy is further exploied by marginalizing the unique values from the state when updating the cluster allocations, which leads to improved efficiency in terms of effective sample size, but might require longer runtimes. For more information, please refer to the Neal2Algorithm class, as well as BaseAlgorithm and MarginalAlgorithm on which it is based.

Public Functions

inline virtual bayesmix::AlgorithmId get_id() const override¶: Returns the Protobuf ID associated to this class.

class Neal8Algorithm : public Neal2Algorithm ¶

Template class for Neal’s algorithm 8 for conjugate hierarchies.

This class implements Neal’s Gibbs sampling algorithm 8 from Neal (2000) that generates a Markov chain on the clustering of the provided data.

It extends Neal’s algorithm 2 to also deal with cases when the kernel/likelihood f(x | phi) is not conjugate to G, thanks to the introduction of additional, auxiliary unique values.

Public Functions

inline virtual bool requires_conjugate_hierarchy() const override¶: Returns whether it is restricted to conjugate Hierarchy objects.

inline unsigned int get_n_aux() const¶: Returns number of auxiliary blocks.

inline void set_n_aux(const unsigned int n_aux_)¶: Sets number of auxiliary blocks.

inline virtual bayesmix::AlgorithmId get_id() const override¶: Returns the Protobuf ID associated to this class.

virtual void read_params_from_proto(const bayesmix::AlgorithmParams &params) override¶: Reads and sets algorithm parameters from an appropriate Protobuf message.

class SplitAndMergeAlgorithm : public MarginalAlgorithm ¶

Template class for Split and Merge for conjugate hierarchies.

This class implements Split and Merge algorithm from Jain and Neal (2004) that generates a Markov chain on the clustering of the provided data.

This algorithm requires the use of a ConjugateHierarchy object. After picking at random two points, the Split and Merge algorithm takes all the points that are in the same clusters as the first two and assigns to each of them, at random, one of two temporary clusters. Then, it updates these assignments by doing a restricted Gibbs sampling (GS), restricted in the sense that it considers only the points cited above and only the two temporary clusters. The result of this operation is used to compute the acceptance probability for the split or merge proposal. These Metropolis-Hastings (MH) steps are alternated with full Gibbs sampling steps, similar to the ones that are computed in Neal3 algorithm.

Public Functions

inline virtual bool requires_conjugate_hierarchy() const override¶: Returns whether it is restricted to conjugate Hierarchy objects.

inline virtual bayesmix::AlgorithmId get_id() const override¶: Returns the Protobuf ID associated to this class.

virtual void read_params_from_proto(const bayesmix::AlgorithmParams &params) override¶: Reads and sets algorithm parameters from an appropriate Protobuf message.

class ConditionalAlgorithm : public BaseAlgorithm ¶

Template class for a conditional sampler deriving from BaseAlgorithm.

This template class implements a generic Gibbs sampling conditional algorithm as the child of the BaseAlgorithm class. A mixture model sampled from a conditional algorithm can be expressed as

$x_i \mid c_i, \theta_1, \dots, \theta_k &\sim f(x_i \mid \theta_{c_i}) \\ \theta_1, \dots, \theta_k &\sim G_0 \\ c_1, \dots, c_n \mid w_1, \dots, w_k &\sim \text{Cat}(w_1, \dots, w_k) \\ w_1, \dots, w_k &\sim p(w_1, \dots, w_k)$

where $f(x \mid \theta_j)$ is a density for each value of $\theta_j$ , $c_i$ take values in $\{1, \dots, k\}$ and $w_1, \dots, w_k$ are nonnegative weights whose sum is a.s. 1, i.e. $p(w_1, \dots, w_k)$ is a probability distribution on the k-1 dimensional unit simplex). In this library, each $\theta_j$ is represented as an Hierarchy object (which inherits from AbstractHierarchy), that also holds the information related to the base measure $G$ is (see AbstractHierarchy). The weights $(w_1, \dots, w_k)$ are represented as a Mixing object, which inherits from AbstractMixing.

The state of a conditional algorithm consists of the unique values, the cluster allocations and the mixture weights. The former two are stored in this class, while the weights are stored in the Mixing object.

Subclassed by BlockedGibbsAlgorithm, SliceSampler

Public Functions

inline virtual bool is_conditional() const override¶: Returns whether the algorithm is conditional or marginal.

virtual Eigen::VectorXd lpdf_from_state(const Eigen::MatrixXd &grid, const Eigen::RowVectorXd &hier_covariate, const Eigen::RowVectorXd &mix_covariate) override¶

Evaluates the estimated log-lpdf on the state contained in curr_state

Parameters:

grid – Grid of row points on which the density is to be evaluated
hier_covariate – Optional covariates related to the Hierarchy
mix_covariate – Optional covariates related to the Mixing

Returns:

The estimation on the iteration of curr_state over all points

class BlockedGibbsAlgorithm : public ConditionalAlgorithm ¶

Template class for the blocked Gibbs sampling algorithm.

This class implement the blocked Gibbs sampling procedure from [1].

[1] Ishwaran, H., & James, L. F. (2001). Gibbs sampling methods for stick-breaking priors. Journal of the American Statistical Association, 96(453), 161-173.

Public Functions

inline virtual bayesmix::AlgorithmId get_id() const override¶: Returns the Protobuf ID associated to this class.