To appear in the second edition of the MCMC handbook, S. P. Brooks, A. Gelman, G. Jones and X.-L. Meng (eds), Chapman & Hall.

Symbolic data analysis (SDA) aggregates large individual-level datasets into a small number of distributional summaries, such as random rectangles or random histograms. The inference is carried out using these summaries in place of the original dataset, resulting in computational gains at the loss of some information. In likelihood-based SDA, the likelihood function is characterised by an integral with a large exponent, which limits the method's utility as for typical models the integral is unavailable in closed form. In addition, the likelihood function is known to produce biased parameter estimates in some circumstances. Our article develops a Bayesian framework for SDA methods in these settings that resolves the issues resulting from integral intractability and biased parameter estimation using pseudo-marginal Markov chain Monte Carlo methods. We develop an exact but computationally expensive method based on path sampling and the Poisson estimator, and a much faster, but approximate, method based on a Taylor expansion. Through simulation and real-data examples we demonstrate the performance of the developed methods, showing large reductions in computation time compared to the full-data analysis, with only a small loss of information.

Calibration ensures that predicted uncertainties align with observed uncertainties. While there is an extensive literature on recalibration methods for univariate probabilistic forecasts, work on calibration for multivariate forecasts is much more limited. This paper introduces a novel post-hoc recalibration approach that addresses multivariate calibration for potentially misspecified models. Our method involves constructing local mappings between vectors of marginal probability integral transform values and the space of observations, providing a flexible and model free solution applicable to continuous, discrete, and mixed responses. We present two versions of our approach: one uses K-nearest neighbors, and the other uses normalizing flows. Each method has its own strengths in different situations. We demonstrate the effectiveness of our approach on two real data applications: recalibrating a deep neural network's currency exchange rate forecast and improving a regression model for childhood malnutrition in India for which the multivariate response has both discrete and continuous components.

Stochastic partition processes divide a multi-dimensional space into a number of regions, such that the data within each region exhibit some form of homogeneity. Due to the nature of their partition strategies, partition processes can often create many unnecessary divisions in sparse regions when trying to describe data in dense regions. To avoid this problem we introduce a parsimonious partition model - the Rectangular Bounding Process (RBP) - to efficiently partition multi-dimensional spaces, by employing a bounding strategy to enclose data points within rectangular bounding boxes. The RBP is self-consistent and as such can be directly extended from a finite hypercube to an infinite (unbounded) space. We extend the RBP to establish a data-dependent RBP (data-RBP) to generate bounding boxes only over existing data points in a sequential manner, which can effectively reduce model complexity and enable online learning. To achieve this, we design an alternative way to generate bounding boxes and prove the distributional equivalence between the data-RBP and the RBP when empty boxes are removed. We demonstrate application of the RBP and the data-RBP in three scenarios: regression trees, relational modelling, and random feature construction for online learning. Extensive experimental results validate the performance of the RBP and the data-RBP for both accuracy and efficiency.

Positional Encoder Graph Neural Networks (PE-GNNs) are a leading approach for modeling continuous spatial data. However, they often fail to produce calibrated predictive distributions, limiting their effectiveness for uncertainty quantification. We introduce the Positional Encoder Graph Quantile Neural Network (PE-GQNN), a novel method that integrates PE-GNNs, Quantile Neural Networks, and recalibration techniques in a fully nonparametric framework, requiring minimal assumptions about the predictive distributions. We propose a new network architecture that, when combined with a quantile-based loss function, yields accurate and reliable probabilistic models without increasing computational complexity. Our approach provides a flexible, robust framework for conditional density estimation, applicable beyond spatial data contexts. We further introduce a structured method for incorporating a KNN predictor into the model while avoiding data leakage through the GNN layer operation. Experiments on benchmark datasets demonstrate that PE-GQNN significantly outperforms existing state-of-the-art methods in both predictive accuracy and uncertainty quantification.

The analysis of unperturbed tumor growth kinetics, particularly the estimation of parameters for S-shaped equations used to describe growth, requires an appropriate likelihood function that accounts for the increasing error in solid tumor measurements as tumor size grows over time. This study aims to propose suitable likelihood functions for parameter estimation in S-shaped models of unperturbed tumor growth. Five different likelihood functions are evaluated and compared using three Bayesian criteria (the Bayesian Information Criterion, Deviance Information Criterion, and Bayes Factor) along with hypothesis tests on residuals. These functions are applied to fit data from unperturbed Ehrlich, fibrosarcoma Sa-37, and F3II tumors using the Gompertz equation, though they are generalizable to other S-shaped growth models for solid tumors or analogous systems (e.g., microorganisms, viruses). Results indicate that error models with tumor volume-dependent dispersion outperform standard constant-variance models in capturing the variability of tumor measurements, particularly the Thres model, which provides interpretable parameters for tumor growth. Additionally, constant-variance models, such as those assuming a normal error distribution, remain valuable as complementary benchmarks in analysis. It is concluded that models incorporating volume-dependent dispersion are preferred for accurate and clinically meaningful tumor growth modeling, whereas constant-dispersion models serve as useful complements for consistency and historical comparability.

Doubly intractable models are encountered in a number of fields, e.g. social networks, ecology and epidemiology. Inference for such models requires the evaluation of a likelihood function, whose normalising function depends on the model parameters and is typically computationally intractable. The normalising constant of the posterior distribution and the additional normalising function of the likelihood function result in a so-called doubly intractable posterior, for which it is difficult to directly apply Markov chain Monte Carlo (MCMC) methods. We propose a signed pseudo-marginal Metropolis-Hastings (PMMH) algorithm with an unbiased block-Poisson estimator to sample from the posterior distribution of doubly intractable models. As the estimator can be negative, the algorithm targets the absolute value of the estimated posterior and uses an importance sampling correction to ensure simulation consistent estimates of the posterior mean of any function. The advantages of our estimator over previous approaches are that its form is ideal for correlated pseudo-marginal methods which are well known to dramatically increase sampling efficiency. Moreover, we develop analytically derived heuristic guidelines for optimally tuning the hyperparameters of the estimator. We demonstrate the algorithm on the Ising model and a Kent distribution model for spherical data.

Machine learning is increasingly being applied in polymer chemistry to link chemical structures to macroscopic properties of polymers and to identify chemical patterns in the polymer structures that help improve specific properties. To facilitate this, a chemical dataset needs to be translated into machine readable descriptors. However, limited and inadequately curated datasets, broad molecular weight distributions, and irregular polymer configurations pose significant challenges. Most off the shelf mathematical models often need refinement for specific applications. Addressing these challenges demand a close collaboration between chemists and mathematicians as chemists must formulate research questions in mathematical terms while mathematicians are required to refine models for specific applications. This review unites both disciplines to address dataset curation hurdles and highlight advances in polymer synthesis and modeling that enhance data availability. It then surveys ML approaches used to predict solid-state properties, solution behavior, composite performance, and emerging applications such as drug delivery and the polymer-biology interface. A perspective of the field is concluded and the importance of FAIR (findability, accessibility, interoperability, and reusability) data and the integration of polymer theory and data are discussed, and the thoughts on the machine-human interface are shared.

This study aimed to identify solvent characteristics that enhance drug loading in polymeric micelles. Polyethylene glycol-block-polystyrene (PEG-b-PS) and curcumin were used as model compounds to investigate the impact of 40 different solvent mixtures on drug loading during flow-based assembly. We tested five algorithms: Random Forest (RF), Gradient Boosting (GP), XGBoost, Support Vector Regression (SVR), and Multilayer Perceptron (MLP), with the MLP model proving to be the most effective among them. To explain the model's predictions, we utilized SHapley Additive exPlanations (SHAP) values to identify solvent properties that contribute to high drug loading. Of the nine descriptors examined-curcumin solubility, polarity, Hildebrand solubility parameters, dipole moment, dielectric constants, viscosity, and Hansen solubility parameters (δD, δP, and δH)-solubility emerged as the most critical factor. Therefore, to achieve optimal drug loading, researchers should prioritize solvents with the highest solubility.

The expressiveness of flow-based models combined with stochastic variational inference (SVI) has expanded the application of optimization-based Bayesian inference to highly complex problems. However, despite the importance of multi-model Bayesian inference, defined over a transdimensional joint model and parameter space, flow-based SVI has been limited to problems defined over a fixed-dimensional parameter space. We introduce CoSMIC normalizing flows (COntextually-Specified Masking for Identity-mapped Components), an extension to neural autoregressive conditional normalizing flow architectures that enables use of a single amortized variational density for inference over a transdimensional (multi-model) conditional target distribution. We propose a combined stochastic variational transdimensional inference (VTI) approach to training CoSMIC flows using ideas from Bayesian optimization and Monte Carlo gradient estimation. Numerical experiments show the performance of VTI on challenging problems that scale to high-cardinality model spaces.

Sampling from multivariate normal distributions, subjected to a variety of restrictions, is a problem that is recurrent in statistics and computing. In the present work, we demonstrate a general framework to efficiently sample a multivariate normal distribution subject to any set of linear inequality constraints and/or linear equality constraints simultaneously. In the approach we detail, sampling a multivariate random variable from the domain formed by the intersection of linear constraints proceeds via a combination of elliptical slice sampling to address the inequality constraints, and linear mapping to address the equality constraints. We also detail a linear programming method for finding an initial sample on the linearly constrained domain; such a method is critical for sampling problems where the domain has small probability. We demonstrate the validity of our methods on an arbitrarily chosen four-dimensional multivariate normal distribution subject to five inequality constraints and/or two equality constraints. Our approach compares favourably to direct sampling and/or accept-reject sampling methods; the latter methods vary widely in their efficiency, whereas the methods in the present work are rejection-free. Where practical we compare predictions of probability density functions between our sampling methods and analytical computation. For all simulations we demonstrate that our methods yield accurate computation of the mean and covariance of the multivariate normal distributions restricted by the imposed linear constraints. MATLAB codes to implement our methods are readily available at https://dx.doi.org/10.6084/m9.figshare.29956304 .

Quantum computers promise to surpass the most powerful classical supercomputers when it comes to solving many critically important practical problems, such as pharmaceutical and fertilizer design, supply chain and traffic optimization, or optimization for machine learning tasks. Because quantum computers function fundamentally differently from classical computers, the emergence of quantum computing technology will lead to a new evolutionary branch of statistical and data analytics methodologies. This review provides an introduction to quantum computing designed to be accessible to statisticians and data scientists, aiming to equip them with an overarching framework of quantum computing, the basic language and building blocks of quantum algorithms, and an overview of existing quantum applications in statistics and data analysis. Our goal is to enable statisticians and data scientists to follow quantum computing literature relevant to their fields, to collaborate with quantum algorithm designers, and, ultimately, to bring forth the next generation of statistical and data analytics tools.

Max-stable processes serve as the fundamental distributional family in extreme value theory. However, likelihood-based inference methods for max-stable processes still heavily rely on composite likelihoods, rendering them intractable in high dimensions due to their intractable densities. In this paper, we introduce a fast and efficient inference method for max-stable processes based on their angular densities for a class of max-stable processes whose angular densities do not put mass on the boundary space of the simplex. This class can also be used to construct r-Pareto processes. We demonstrate the efficiency of the proposed method through two new max-stable processes: the truncated extremal-t process and the skewed Brown-Resnick process. The skewed Brown-Resnick process contains the popular Brown-Resnick model as a special case and possesses nonstationary extremal dependence structures. The proposed method is shown to be computationally efficient and can be applied to large datasets. We showcase the new max-stable processes on simulated and real data.

Abstract The Gompertz model, a mainstay in tumor growth kinetics analysis, requires an accurate likelihood function for its parameterestimation, applicable in both classical and Bayesian methodologies. This study compares five distinct error models, eachrepresenting a different likelihood function. Our comparative analysis employs the Bayesian Information Criterion (BIC), theDeviance Information Criterion (DIC), the Bayes Factor (BF), and hypothesis tests on residuals. Applying these criteria to fitthe Gompertz model to Ehrlich and fibrosarcoma Sa-37 tumor data, we find that error models with tumor volume-dependentdispersion consistently outperform others in quantitative evaluations. However, the conventional Normal error model withconstant variance remains a vital tool, offering significant clinical insights. This study underscores the complexity of likelihoodmodel selection in tumor growth kinetics and highlights the need for a multifaceted approach in such analysis.

Statistical modelling of spatial extreme events has gained increasing attention over the last few decades with max-stable processes, and more recently r-Pareto processes, becoming the reference tools for the statistical analysis of asymptotically dependent data. Although inference for r-Pareto processes is easier than for max-stable processes, there remain major hurdles for their application to high dimensional datasets within a reasonable timeframe. In addition, both approaches have almost exclusively focused on the Brown-Resnick model, for its Gaussian foundations, and for the continuity of its exponent measure. In this paper, we derive a class of models for which this continuity property holds and present the skewed Brown-Resnick model, an extension of the Brown-Resnick that allows for non-stationarity in the dependence structure, and the truncated extremal-t model, a refinement of the well-known extremal-t model. We use an inference methodology based on the intensity function of the process which is derived from the exponent measure, and demonstrate the statistical and computational efficiency of this approach. Applications to two real-world problems illustrate valuable gains in modelling flexibility as well as appealing computational gains over reference methodologies.

Artificial neural networks (ANNs) are highly flexible predictive models. However, reliably quantifying uncertainty for their predictions is a continuing challenge. There has been much recent work on "recalibration" of predictive distributions for ANNs, so that forecast probabilities for events of interest are consistent with certain frequency evaluations of them. Uncalibrated probabilistic forecasts are of limited use for many important decision-making tasks. To address this issue, we propose a localized recalibration of ANN predictive distributions using the dimension-reduced representation of the input provided by the ANN hidden layers. Our novel method draws inspiration from recalibration techniques used in the literature on approximate Bayesian computation and likelihood-free inference methods. Most existing calibration methods for ANNs can be thought of as calibrating either on the input layer, which is difficult when the input is high-dimensional, or the output layer, which may not be sufficiently flexible. Through a simulation study, we demonstrate that our method has good performance compared to alternative approaches, and explore the benefits that can be achieved by localizing the calibration based on different layers of the network. Finally, we apply our proposed method to a diamond price prediction problem, demonstrating the potential of our approach to improve prediction and uncertainty quantification in real-world applications.

Stein variational gradient descent (SVGD) is a particle based approximate inference algorithm with largely well understood theoretical properties. In recent years, many variants of SVGD have been proposed and shown to share those properties. A preliminary test of the hybrid kernel variant (h-SVGD) has demonstrated promising results on image classification with deep neural network ensembles. However, the theoretical properties of h-SVGD have not yet been established, and its practical advantages have not been fully explored. In this paper, we define a hybrid kernelised Stein discrepancy (h-KSD) and prove that the h-SVGD update direction is optimal within an appropriate reproducing kernel Hilbert space. We also prove a descent lemma that guarantees a decrease in the KL divergence at each step along with other limit results. Numerical results demonstrate that h-SVGD mitigates the variance collapse behaviour of SVGD at no additional computational cost whilst remaining competitive at inference tasks.

Gaussian process state-space models (GPSSMs) provide a principled and flexible approach to modeling the dynamics of a latent state, which is observed at discrete-time points via a likelihood model. However, inference in GPSSMs is computationally and statistically challenging due to the large number of latent variables in the model and the strong temporal dependencies between them. In this paper, we propose a new method for inference in Bayesian GPSSMs, which overcomes the drawbacks of previous approaches, namely over-simplified assumptions, and high computational requirements. Our method is based on free-form variational inference via stochastic gradient Hamiltonian Monte Carlo within the inducing-variable formalism. Furthermore, by exploiting our proposed variational distribution, we provide a collapsed extension of our method where the inducing variables are marginalized analytically. We also showcase results when combining our framework with particle MCMC methods. We show that, on six real-world datasets, our approach can learn transition dynamics and latent states more accurately than competing methods.

There has been much recent interest in modifying Bayesian inference for misspecified models so that it is useful for specific purposes. One popular modified Bayesian inference method is "cutting feedback" which can be used when the model consists of a number of coupled modules, with only some of the modules being misspecified. Cutting feedback methods represent the full posterior distribution in terms of conditional and sequential components, and then modify some terms in such a representation based on the modular structure for specification or computation of a modified posterior distribution. The main goal of this is to avoid contamination of inferences for parameters of interest by misspecified modules. Computation for cut posterior distributions is challenging, and here we consider cutting feedback for likelihood-free inference based on Gaussian mixture approximations to the joint distribution of parameters and data summary statistics. We exploit the fact that marginal and conditional distributions of a Gaussian mixture are Gaussian mixtures to give explicit approximations to marginal or conditional posterior distributions so that we can easily approximate cut posterior analyses. The mixture approach allows repeated approximation of posterior distributions for different data based on a single mixture fit. This is important for model checks which aid in the decision of whether to "cut". A semi-modular approach to likelihood-free inference where feedback is partially cut is also developed. The benefits of the method are illustrated on two challenging examples, a collective cell spreading model and a continuous time model for asset returns with jumps.