We propose a Bayesian nonparametric model for grouped data called the Shared Atoms Model (SAM). SAM generates dependent random discrete distributions through “atom skipping,” where a mixture component may have zero weight in some groups. This allows the model to represent unique, shared, and common clusters, and to estimate the extent of overlap across groups. As a result, interpretable posterior inference is feasible, such as reporting the posterior probability of a unique cluster exclusive to a single group or a shared cluster belonging to some but not all groups. We discuss the theoretical properties of the proposed and related models. Minor extensions of the proposed model for multivariate or count data are presented. Simulation studies and applications using real-world datasets illustrate the performance of the new models in comparison with existing models. The code is publicly available at https://github.com/edwardbi/SharedAtomsModel.
Bayesian nonparametrics (BNP) is a major success story in statistics. The field has grown from a small, niche area of a handful of active researchers in the 1970s and 1980s to the current large, vibrant community with a major presence at all of the big statistical conferences and beyond. While there are many reasons for this growth, this paper offers a modern perspective on one specific contribution: a class of BNP models for, in the language of regression, Y∣x, or more formally, the collection of conditional distributions that we refer to as Fx. Designed with minimal technical detail, this paper is intended to accompany MacEachern (2026). Within it, we emphasize a framework for modelling and inference that aligns both with this class of models and with traditional modelling strategies. Finally, we highlight several promising but underdeveloped applications of these models.
Dependent Dirichlet processes provide a means of modelling a collection of random distributions as related to, but not identical to, each other. A key point is that the realizations of these processes are dependent. This article develops the processes through definition and results. It also describes the basic tools for computations with the processes and suggests a number of applications in which they will be useful. The extension of the Dirichlet process to the dependent Dirichlet process can be paralleled for essentially any other nonparametric process which, for a single random distribution, relies on a countable collection of random variables. A parallel development for two such processes is indicated.
Preferential sampling is a common feature in geostatistics and occurs when the locations are sampled based on information about the phenomena under study. In this case, point pattern models are commonly used as the probability law for the distribution of the locations. However, analytic intractability of the point process likelihood prevents its direct calculation. Many Bayesian (and non-Bayesian) approaches in non-parametric model specifications handle this difficulty with approximations, both to the model and to the computations required for drawing inference. Procedures to approximate the model lead to errors that are sometimes difficult to quantify and can lead to biased inference. This paper presents an approach for performing exact Bayesian inference for this setting without the need for model approximation. A qualitatively minor change on the traditional model is proposed to circumvent the likelihood intractability. This change enables the use of an augmented model strategy. Recent work on Bayesian inference for point pattern models can be adapted to the geostatistics setting and renders computational tractability for exact inference for the proposed methodology. Estimation of model parameters and prediction of the response at unsampled locations can then be obtained from the joint posterior distribution. Simulated studies showed good quality of the proposed model for estimation and prediction in a variety of scenarios. The performance of our approach is illustrated in the analysis of simulated and real datasets and also compares favourably against approximation-based approaches. The paper is concluded with comments regarding extensions and improvements to the proposed methodology.
We consider the problem of boundary detection for areal data, focusing on situations where for each areal unit multiple observations are available. We propose a Bayesian nonparametric mixture model for the area-specific population densities, with spatially dependent weights and a random number of components. Contrary to previously proposed methods for boundary detection, which consider one observation per areal unit, ours does not require external information such as area-specific covariates or dissimilarity metrics. Instead, by exploiting information from multiple samples per area, it is able to identify boundaries between areas that exhibit different densities. Crucially, the number of mixture components needs to be learned from data to obtain meaningful boundary detection, due to the non-identifiability of overfitted mixtures. Therefore, we assume it random by placing a prior on it. The motivating application is the analysis of economic inequality in the greater Los Angeles region, which typically yields social inequality and unrest. Efficient posterior computation is facilitated by a transdimensional Markov Chain Monte Carlo sampler which exploits the recently introduced optimal auxiliary priors to improve the mixing. The methodology is validated via extensive simulations and applied to the income data in the greater Los Angeles region. We identify several boundaries in the income distributions, which can be explained ex-post in terms of the percentage of the population without health insurance, though not in terms of the total number of crimes, showing the usefulness of such an analysis to policymakers.
We introduce an extension of a recently developed semi-parametric generalized linear model by adding a nonparametric Bayesian prior on the baseline distribution of the GLM. The construction is similar to a varying-weight dependent Dirichlet process. We show that the resulting model takes the form of an inhomogeneous completely random measure that arises from exponential tilting of a normalized completely random measure. Building on familiar posterior sampling methods for mixtures with respect to normalized random measures, we introduce posterior simulation in the resulting model. We validate the proposed methodology through extensive simulation studies and illustrate its application using data from a speech intelligibility study.
Posterior samples of factor loadings in Bayesian exploratory factor analysis are not directly comparable across MCMC iterations due to rotational indeterminacy. As a result, their posterior means are typically close to zero because they cancel across rotationally equivalent orientations, yielding uninterpretable factor loading estimates. Alignment-based post-processing is therefore required, but existing approaches face a trade-off: exact alignment methods do not scale well beyond low-dimensional settings, whereas scalable alternatives rely on approximations that can reduce accuracy. We introduce an efficient version of the Rotation-Sign-Permutation (RSP) algorithm of Papastamoulis and Ntzoufras (2022) that overcomes these limitations by making exact alignment scalable to a large number of factors. Simulations and empirical examples demonstrate that this algorithm achieves higher alignment accuracy in high-dimensional settings with negligible computational overhead. An optimized C++ implementation is available in the open-source R package BayesianEFA.
We propose a fully Bayesian approach for causal inference with multivariate categorical data based on staged tree models, a class of probabilistic graphical models capable of representing asymmetric and context-specific dependencies. To account for uncertainty in both structure and parameters, we introduce a flexible family of prior distributions over staged trees. These include product partition models to encourage parsimony, a novel distance-based prior to promote interpretable dependence patterns, and an extension that incorporates continuous covariates into the learning process. Posterior inference is achieved via a tailored Markov Chain Monte Carlo algorithm with split-and-merge moves, yielding posterior samples of staged trees from which average treatment effects and uncertainty measures are derived. Posterior summaries and uncertainty measures are obtained via techniques from the Bayesian nonparametrics literature. Two case studies on electronic fetal monitoring and cesarean delivery and on anthracycline therapy and cardiac dysfunction in breast cancer illustrate the methods.
We develop a Bayesian tree ensemble model to estimate heterogeneous treatment effects in censored survival data with high-dimensional covariates. Instead of imposing sparsity through the tree structure, we place a horseshoe prior directly on the step heights to achieve adaptive global-local shrinkage. This strategy allows flexible regularisation and reduces noise. We develop a reversible jump Gibbs sampler to accommodate the non-conjugate horseshoe prior within the tree ensemble framework. We show through extensive simulations that the method accurately estimates treatment effects in high-dimensional covariate spaces, at various sparsity levels, and under non-linear treatment effect functions. We further illustrate the practical utility of the proposed approach by a re-analysis of pancreatic ductal adenocarcinoma (PDAC) survival data from The Cancer Genome Atlas.
The hierarchical Dirichlet process is the cornerstone of Bayesian nonparametric multilevel models. Its generative model can be described through a set of latent variables, commonly referred to as tables within the popular restaurant franchise metaphor. The latent tables simplify the expression of the posterior and allow for the implementation of Gibbs sampling algorithms to approximately draw posterior samples. However, managing their assignments can become computationally expensive, especially as the size of the dataset and the number of levels increase. In this work, we identify a prior for the concentration parameter of the hierarchical Dirichlet process that (i) induces a quasi-conjugate posterior distribution, and (ii) removes the need for tables, leading to more interpretable expressions for the posterior, with both a scalable and an exact algorithm to sample from it. Remarkably, this construction extends beyond the Dirichlet process, leading to a new framework for defining normalized hierarchical random measures and a new class of algorithms to sample from their posteriors. The key analytical tool is the independence of multivariate increments, that is, their representation as completely random vectors.
In many research fields, researchers aim to identify significant associations between a set of explanatory variables and a response while controlling the false discovery rate (FDR). The Knockoff filter has been recently proposed in the frequentist paradigm to introduce controlled noise in a model by cleverly constructing copies of the predictors as auxiliary variables. In this paper, we develop a fully Bayesian generalization of the classical model-X knockoff filter for normally distributed covariates. In our approach, we consider a joint model for the covariates and the response, where the conditional independence structure of the covariates is captured through a Gaussian graphical model and used to define a latent knockoff layer through a parameter-expanded representation of the response model. Estimating the covariate graph informs the knockoff construction and improves inference on the covariate effects. We use a modified spike-and-slab prior on the regression coefficients, which avoids the increase of the model dimension typical of the classical knockoff filter. We also address extensions to settings with non-Gaussian responses. Our model performs variable selection using an upper bound on the posterior probability of non-inclusion. We show that the induced latent knockoff layer defines valid Gaussian model-X knockoffs under the proposed construction and that the resulting procedure controls the Bayesian FDR at an arbitrary level, in finite samples, if the distribution of the covariates is fully known; under an estimated graphical structure, it satisfies an asymptotic FDR guarantee. We use simulated data to demonstrate that our proposal increases the stability of the selection with respect to classical knockoff methods, as it relies on the entire posterior distribution of the latent knockoff variables instead of a single sample. With respect to Bayesian variable selection methods, we show that our selection procedure achieves comparable or better performance, while maintaining control over the FDR. Finally, we show the usefulness of the proposed model with an application to real data.
The use of hierarchical mixture priors with shared atoms has recently flourished in the Bayesian literature for partially exchangeable data. Leveraging on nested levels of mixtures, these models allow the estimation of a two-layered data partition: across groups and across observations. This paper discusses and compares the properties of such modeling strategies when the mixing weights are assigned either a finite-dimensional Dirichlet distribution or a Dirichlet process prior. Based on these considerations, we introduce a novel hierarchical nonparametric prior based on a finite set of shared atoms, a specification that enhances the flexibility of the induced random measures and the availability of fast posterior inference. To support these findings, we analytically derive the induced prior correlation structure and partially exchangeable partition probability function. Additionally, we develop a novel mean-field variational algorithm for posterior inference to boost the applicability of our nested model to large multivariate data. We then assess and compare the performance of the different shared-atom specifications via simulation. We also show that our variational proposal is highly scalable and that the accuracy of the posterior density estimate and the estimated partition is comparable with state-of-the-art Gibbs sampler algorithms. Finally, we apply our model to a real dataset of Spotify's song features, simultaneously segmenting artists and songs with similar characteristics.
Analysing non-Gaussian spatial-temporal data requires introducing spatial as well as temporal dependence in generalised linear models through the link function of an exponential family distribution. Unlike in Gaussian likelihoods, inference is considerably encumbered by the inability to analytically integrate out the random effects and reduce the dimension of the parameter space. Iterative estimation algorithms struggle to converge due to the presence of weakly identified parameters. We devise Bayesian inference using predictive stacking that assimilates inference from analytically tractable conditional posterior distributions. We achieve this by expanding upon the Diaconis-Ylvisaker family of conjugate priors and exploiting generalised conjugate multivariate (GCM) distribution theory for exponential families, which enables exact sampling from analytically available posterior distributions conditional upon some process parameters. Subsequently, we assimilate inference over a range of values of these parameters using Bayesian predictive stacking. We evaluate inferential performance on simulated data, compare with full Bayesian inference using Markov chain Monte Carlo (MCMC) and apply our method to analyse spatially-temporally referenced avian count data from the North American Breeding Bird Survey database.
We develop a Bayesian inference method for discretely-observed stochastic differential equations (SDEs). Inference is challenging for most SDEs, due to the analytical intractability of the likelihood function. Nevertheless, forward simulation via numerical methods is straightforward, motivating the use of approximate Bayesian computation (ABC). We propose a conditional simulation scheme for SDEs that is based on lookahead strategies for sequential Monte Carlo (SMC) and particle smoothing using backward simulation. This leads to the simulation of trajectories that are consistent with the observed trajectory, thereby increasing the ABC acceptance rate. We additionally employ an invariant neural network, previously developed for Markov processes, to learn the summary statistics function required in ABC. The neural network is incrementally retrained by exploiting an ABC-SMC sampler, which provides new training data at each round. Since the SDEs simulation scheme differs from standard forward simulation, we propose a suitable importance sampling correction, which has the added advantage of guiding the parameters towards regions of high posterior density, especially in the first ABC-SMC round. Our approach achieves accurate inference and is about three times faster than standard (forward-only) ABC-SMC. We illustrate our method in five simulation studies, including three examples from the Chan-Karaolyi-Longstaff-Sanders SDE family, a stochastic bi-stable model (Schl{\"o}gl) that is notoriously challenging for ABC methods, and a two dimensional biochemical reaction network.
Likelihood-based inference in stochastic non-linear dynamical systems, such as those found in chemical reaction networks and biological clock systems, is inherently complex and has largely been limited to small and unrealistically simple systems. Recent advances in analytically tractable approximations to the underlying conditional probability distributions enable long-term dynamics to be accurately modelled, and make the large number of model evaluations required for exact Bayesian inference much more feasible. We propose a new methodology for inference in stochastic non-linear dynamical systems exhibiting oscillatory behaviour and show the parameters in these models can be realistically estimated from simulated data. Preliminary analyses based on the Fisher Information Matrix of the model can guide the implementation of Bayesian inference. We show that this parameter sensitivity analysis can predict which parameters are practically identifiable. Several Markov chain Monte Carlo algorithms are compared, with our results suggesting a parallel tempering algorithm consistently gives the best approach for these systems, which are shown to frequently exhibit multi-modal posterior distributions.
Meta-analysis aims to generalize results from multiple related statistical analyses through a combined analysis. While the natural outcome of a Bayesian study is a posterior distribution, traditional Bayesian meta-analyses proceed by combining summary statistics (i.e., point-valued estimates) computed from data. In this paper, we develop a framework for combining posterior distributions from multiple related Bayesian studies into a meta-analysis. Importantly, the method is capable of reusing pre-computed posteriors from computationally costly analyses, without needing the implementation details from each study. Besides providing a consensus across studies, the method enables updating the local posteriors post-hoc and therefore refining them by sharing statistical strength between the studies, without rerunning the original analyses. We illustrate the wide applicability of the framework by combining results from likelihood-free Bayesian analyses, which would be difficult to carry out using standard methodology.
We consider estimation of a normal mean matrix under the Frobenius loss. Motivated by the Efron--Morris estimator, a generalization of Stein's prior has been recently developed, which is superharmonic and shrinks the singular values towards zero. The generalized Bayes estimator with respect to this prior is minimax and dominates the maximum likelihood estimator. However, here we show that it is inadmissible by using Brown's condition. Then, we develop two types of priors that provide improved generalized Bayes estimators and examine their performance numerically. The proposed priors attain risk reduction by adding scalar shrinkage or column-wise shrinkage to singular value shrinkage. Parallel results for Bayesian predictive densities are also given.
Gaussian processes (GPs) are a highly flexible, nonparametric statistical model that are commonly used to fit nonlinear relationships or account for correlation between observations. However, the computational load of fitting a Gaussian process is $\mathcal{O}(n^3)$ making them infeasible for use on large datasets. To make GPs more feasible for large datasets, this research focuses on the use of minibatching to estimate GP parameters. Specifically, we outline both approximate and exact minibatch Markov chain Monte Carlo algorithms that substantially reduce the computation of fitting a GP by only considering small subsets of the data at a time. We demonstrate and compare this methodology using various simulations and real datasets.
Bayesian hypothesis testing leverages posterior probabilities, Bayes factors, or credible intervals to assess characteristics that summarize data. We propose a framework for power curve approximation with such hypothesis tests that assumes data are generated using statistical models with fixed parameters for the purposes of sample size determination. We present a fast approach to explore the sampling distribution of posterior probabilities when the conditions for the Bernstein-von Mises theorem are satisfied. We extend that approach to facilitate targeted sampling from the approximate sampling distribution of posterior probabilities for each sample size explored. These sampling distributions are used to construct power curves for various types of posterior analyses. Our resulting method for power curve approximation is orders of magnitude faster than conventional power curve estimation for Bayesian hypothesis tests. We also prove the consistency of the corresponding power estimates and sample size recommendations under certain conditions.
