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.

Gaussian variational approximations are widely used for summarizing posterior distributions in Bayesian models, especially in high-dimensional settings. However, a drawback of such approximations is the inability to capture skewness or more complex features of the posterior. Recent work suggests applying skewness corrections to existing Gaussian or other symmetric approximations to address this limitation. We propose to incorporate the skewness correction into the definition of an approximating variational family. We consider approximating the posterior for hierarchical models, in which there are “global” and “local” parameters. A baseline variational approximation is defined as the product of a Gaussian marginal posterior for global parameters and a Gaussian conditional posterior for local parameters given the global ones. Skewness corrections are then considered. The adjustment of the conditional posterior term for local variables is adaptive to the global parameter value. Optimization of baseline variational parameters is performed jointly with the skewness correction. Our approach allows the location, scale and skewness to be captured separately, without using additional parameters for skewness adjustments. The proposed method substantially improves accuracy for only a modest increase in computational cost compared to state-of-the-art Gaussian approximations. Good performance is demonstrated in generalized linear mixed models and multinomial logit discrete choice models.

Modular Bayesian methods perform inference in models that are specified through a collection of coupled sub-models, known as modules. These modules often arise from modelling different data sources or from combining domain knowledge from different disciplines. “Cutting feedback” is a Bayesian inference method that ensures misspecification of one module does not affect inferences for parameters in other modules, and produces what is known as the cut posterior. However, choosing between the cut posterior and the standard Bayesian posterior is challenging. When misspecification is not severe, cutting feedback can greatly increase posterior uncertainty without a large reduction of estimation bias, leading to a bias-variance trade-off. This trade-off motivates semi-modular posteriors, which interpolate between standard and cut posteriors based on a tuning parameter. In this work, we provide the first precise formulation of the bias-variance trade-off that is present in cutting feedback, and we propose a new semi-modular posterior that takes advantage of it. Under general regularity conditions, we prove that this semi-modular posterior is more accurate than the cut posterior according to a notion of posterior risk. An important implication of this result is that point inferences made under the cut posterior are inadmissable. The new method is demonstrated in a number of examples.

We propose a new distributional regression model for a multivariate response vector based on a copula process over the covariate space. It uses the implicit copula of a Gaussian multivariate regression, which we call a "regression copula process." To allow for large covariate vectors the regression coefficients are regularized using a novel multivariate extension of the horseshoe prior. Bayesian inference and distributional predictions are evaluated using efficient variational inference methods, allowing application to large datasets. An advantage of the approach is that the marginal distributions of the response vector can be estimated separately and accurately, resulting in predictive distributions that are marginally calibrated. Two applications of the methodology illustrate its efficacy. The first is the econometric modeling and prediction of half-hourly regional Australian electricity prices. The second is the evaluation of multivariate posteriors in likelihood-free inference (LFI) for a model of tree species abundance. This extends a previous univariate regression copula LFI method for evaluating univariate posteriors. In both applications, we demonstrate that our new approach exhibits a desirable marginal calibration property. Supplementary materials for this article are available online.

Likelihood-free inference (LFI) methods, such as approximate Bayesian computation, have become commonplace for conducting inference in complex models. Many approaches are based on summary statistics or discrepancies derived from synthetic data. However, determining which summary statistics or discrepancies to use for constructing the posterior remains a challenging question, both practically and theoretically. Instead of relying on a single vector of summaries for inference, we propose a new pooled posterior that optimally combines inferences from multiple LFI posteriors. This pooled approach eliminates the need to select a single vector of summaries or even a specific LFI algorithm. Our approach is straightforward to implement and avoids performing a high-dimensional LFI analysis involving all summary statistics. We give theoretical guarantees for the improved performance of the pooled posterior mean in terms of asymptotic frequentist risk and demonstrate the effectiveness of the approach in a number of benchmark examples.

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.

In copula models the marginal distributions and copula function are specified separately. We treat these as two modules in a modular Bayesian inference framework, and propose conducting modified Bayesian inference by “cutting feedback”. Cutting feedback limits the influence of potentially misspecified modules in posterior inference. We consider two types of cuts. The first limits the influence of a misspecified copula on inference for the marginals, which is a Bayesian analogue of the popular Inference for Margins (IFM) estimator. The second limits the influence of misspecified marginals on inference for the copula parameters by using a rank likelihood to define the cut model. We establish that if only one of the modules is misspecified, then the appropriate cut posterior gives accurate uncertainty quantification asymptotically for the parameters in the other module. Computation of the cut posteriors is difficult, and new variational inference methods to do so are proposed. The efficacy of the new methodology is demonstrated using both simulated data and a substantive multivariate time series application from macroeconomic forecasting. In the latter, cutting feedback from misspecified marginals to a 1096 dimension copula improves posterior inference and predictive accuracy greatly, compared to conventional Bayesian inference.

Simulation-based Bayesian inference (SBI) methods are widely used for parameter estimation in complex models where evaluating the likelihood is challenging but generating simulations is relatively straightforward. However, these methods commonly assume that the simulation model accurately reflects the true data-generating process, an assumption that is frequently violated in realistic scenarios. In this paper, we focus on the challenges faced by SBI methods under model misspecification. We consolidate recent research aimed at mitigating the effects of misspecification, highlighting three key strategies: i) robust summary statistics, ii) generalised Bayesian inference, and iii) error modelling and adjustment parameters. To illustrate both the vulnerabilities of popular SBI methods and the effectiveness of misspecification-robust alternatives, we present empirical results on an illustrative example.

Statistical models are often defined by a generative process for simulating synthetic data, but this can lead to intractable likelihoods. Likelihood free inference (LFI) methods enable Bayesian inference to be performed in this case. Extending a popular approach to simulation-efficient LFI for single-source data, we propose Multi-objective Bayesian Optimization for Likelihood Free Inference (MOBOLFI) to perform LFI using multi-source data. MOBOLFI models a multi-dimensional discrepancy between observed and simulated data, using a separate discrepancy for each data source. The use of a multivariate discrepancy allows for approximations to individual data source likelihoods in addition to the joint likelihood, enabling detection of conflicting information and deeper understanding of the importance of different data sources in estimating individual parameters. The adaptive choice of simulation parameters using multi-objective Bayesian optimization ensures simulation efficient approximation of likelihood components for all data sources. We illustrate our approach in sequential sampling models (SSMs), which are widely used in psychology and consumer-behavior modeling. SSMs are often fitted using multi-source data, such as choice and response time. The advantages of our approach are illustrated in comparison with a single discrepancy for an SSM fitted to data assessing preferences of ride-hailing drivers in Singapore to rent electric vehicles.

Bayesian simulation-based inference (SBI) methods are used in statistical models where simulation is feasible but the likelihood is intractable. Standard SBI methods can perform poorly in cases of model misspecification, and there has been much recent work on modified SBI approaches which are robust to misspecified likelihoods. However, less attention has been given to the issue of inappropriate prior specification, which is the focus of this work. In conventional Bayesian modelling, there will often be a wide range of prior distributions consistent with limited prior knowledge expressed by an expert. Choosing a single prior can lead to an inappropriate choice, possibly conflicting with the likelihood information. Robust Bayesian methods, where a class of priors is considered instead of a single prior, can address this issue. For each density in the prior class, a posterior can be computed, and the range of the resulting inferences is informative about posterior sensitivity to the prior imprecision. We consider density ratio classes for the prior and implement robust Bayesian SBI using amortized neural methods developed recently in the literature. We also discuss methods for checking for conflict between a density ratio class of priors and the likelihood, and sequential updating methods for examining conflict between different groups of summary statistics. The methods are illustrated for several simulated and real examples.

The Bayesian Synthetic Likelihood (BSL) method is a widely-used tool for likelihood-free Bayesian inference. This method assumes that some summary statistics are normally distributed, which can be incorrect in many applications. We propose a transformation, called the Wasserstein Gaussianization transformation, that uses a Wasserstein gradient flow to approximately transform the distribution of the summary statistics into a Gaussian distribution. BSL also implicitly requires compatibility between simulated summary statistics under the working model and the observed summary statistics. A robust BSL variant which achieves this has been developed in the recent literature. We combine the Wasserstein Gaussianization transformation with robust BSL, and an efficient Variational Bayes procedure for posterior approximation, to develop a highly efficient and reliable approximate Bayesian inference method for likelihood-free problems.

Scientists continue to develop increasingly complex mechanistic models to reflect their knowledge more realistically. Statistical inference using these models can be challenging since the corresponding likelihood function is often intractable and model simulation may be computationally burdensome. Fortunately, in many of these situations it is possible to adopt a surrogate model or approximate likelihood function. It may be convenient to conduct Bayesian inference directly with a surrogate, but this can result in a posterior with poor uncertainty quantification. In this paper, we propose a new method for adjusting approximate posterior samples to reduce bias and improve posterior coverage properties. We do this by optimizing a transformation of the approximate posterior, the result of which maximizes a scoring rule. Our approach requires only a (fixed) small number of complex model simulations and is numerically stable. We develop supporting theory for our method and demonstrate beneficial corrections to approximate posteriors across several examples of increasing complexity.

We consider the problem of estimating complex statistical latent variable models using variational Bayes methods. These methods are used when exact posterior inference is either infeasible or computationally expensive, and they approximate the posterior density with a family of tractable distributions. The parameters of the approximating distribution are estimated using optimisation methods. This article develops a flexible Gaussian mixture variational approximation, where we impose sparsity in the precision matrix of each Gaussian component to reflect the appropriate conditional independence structure in the model. By introducing sparsity in the precision matrix and parameterising it using the Cholesky factor, each Gaussian mixture component becomes parsimonious (with a reduced number of non-zero parameters), while still capturing the dependence in the posterior distribution. Fast estimation methods based on global and local variational boosting moves combined with natural gradients and variance reduction methods are developed. The local boosting moves adjust an existing mixture component, and optimisation is only carried out on a subset of the variational parameters of a new component. The subset is chosen to target improvement of the current approximation in aspects where it is poor. The local boosting moves are fast because only a small number of variational parameters need to be optimised. The efficacy of the approach is illustrated by using simulated and real datasets to estimate generalised linear mixed models and state space models.

Mixtures of linear mixed models are widely used for modelling longitudinal data for which observation times differ between subjects. In typical applications, temporal trends are described using a basis expansion, with basis coefficients treated as random effects varying by subject. Additional random effects can describe variation between mixture components, or other known sources of variation in complex experimental designs. A key advantage of these models is that they provide a natural mechanism for clustering, which can be helpful for interpretation in many applications. Current versions of mixtures of linear mixed models are not specifically designed for the case where there are many observations per subject and a complex temporal trend, which requires a large number of basis functions to capture. In this case, the subject-specific basis coefficients are a high-dimensional random effects vector, for which the covariance matrix is hard to specify and estimate, especially if it varies between mixture components. To address this issue, we consider the use of recently-developed deep mixture of factor analyzers models as the prior for the random effects. The resulting deep mixture of linear mixed models is well-suited to high-dimensional settings, and we describe an efficient variational inference approach to posterior computation. The efficacy of the method is demonstrated on both real and simulated data.

Modern statistics deals with complex models from which the joint model used for inference is built by coupling submodels, called modules. We consider modular inference where the modules may depend on parametric and nonparametric components. In such cases, a joint Bayesian inference is highly susceptible to misspecification across any module, and inappropriate priors for nonparametric components may deliver subpar inferences for parametric components, and vice versa. We propose a novel “optimization-centric” approach to cutting feedback for semiparametric modular inference, which can address misspecification and prior-data conflicts. The proposed generalized cut posteriors are defined through a variational optimization problem for generalized posteriors where regularization is based on Rényi divergence, rather than Kullback-Leibler divergence (KLD), and variational computational methods are developed. We show empirically that using Rényi divergence to define the cut posterior delivers more robust inferences than KLD. We derive novel posterior concentration results that accommodate the Rényi divergence and allow for semiparametric components, greatly extending existing results for cut posteriors that were derived for parametric models and KLD. We demonstrate these new methods in a benchmark toy example and two real examples: Gaussian process adjustments for confounding in causal inference and misspecified copula models with nonparametric marginals.

Despite the increasing attention given to frequentist model averaging methods, few studies address the theoretical aspects of semiparametric model averaging. In this paper, we introduce an exponentially weighted model averaging (EWA) method for the varying-coefficient partially linear model (VCPLM). Our weighting method is straightforward to implement and has a natural Bayesian interpretation. We establish oracle inequalities for the finite sample error bound of the empirical estimators and thoroughly investigate the asymptotic optimality of the proposed weighting criterion, along with other statistical properties. Some numerical studies are provided which support the theoretical findings.

Bayesian inference has many advantages for complex models, but standard Monte Carlo methods for summarizing the posterior can be computationally demanding, and it is attractive to consider optimization-based variational methods. Our work considers Gaussian approximations with sparse precision matrices which are tractable to optimize in high-dimensions. The optimal Gaussian approximation is usually defined as being closest to the posterior in Kullback-Leibler divergence, but it is useful to consider other divergences when the Gaussian assumption is crude, to capture important posterior features for given applications. Our work studies the weighted Fisher divergence, which focuses on gradient differences between the target posterior and its approximation, with the Fisher and score-based divergences as special cases. We make three main contributions. First, we compare approximations for weighted Fisher divergences under mean-field assumptions for Gaussian and non-Gaussian targets with Kullback-Leibler approximations. Second, we go beyond mean-field and consider approximations with sparse precision matrices reflecting posterior conditional independence structure for hierarchical models. Using stochastic gradient descent to enforce sparsity, we develop two approaches to minimize the Fisher and score-based divergences, based on the reparametrization trick and a batch approximation of the objective. Finally, we study the performances of our methods using logistic regression, generalized linear mixed models and stochastic volatility models.

Bayesian inference is a powerful tool for combining information in complex settings, a task of increasing importance in modern applications. However, Bayesian inference with a flawed model can produce unreliable conclusions. This review discusses approaches to performing Bayesian inference when the model is misspecified, where, by misspecified, we mean that the analyst is unwilling to act as if the model is correct. Much has been written about this topic, and in most cases we do not believe that a conventional Bayesian analysis is meaningful when there is serious model misspecification. Nevertheless, in some cases it is possible to use a well-specified model to give meaning to a Bayesian analysis of a misspecified model, and we focus on such cases. Three main classes of methods are discussed: restricted likelihood methods, which use a model based on an insufficient summary of the original data; modular inference methods, which use a model constructed from coupled submodels, with some of the submodels correctly specified; and the use of a reference model to construct a projected posterior or predictive distribution for a simplified model considered to be useful for prediction or interpretation.

Variational Bayes methods approximate the posterior density by a family of tractable distributions whose parameters are estimated by optimisation. Variational approximation is useful when exact inference is intractable or very costly. Our article develops a flexible variational approximation based on a copula of a mixture, which is implemented by combining boosting, natural gradient, and a variance reduction method. The efficacy of the approach is illustrated by using simulated and real datasets to approximate multimodal, skewed and heavy-tailed posterior distributions, including an application to Bayesian deep feedforward neural network regression models.

We analyse the behaviour of the synthetic likelihood (SL) method when the model generating the simulated data differs from the actual data generating process. One of the most common methods to obtain SL-based inferences is via the Bayesian posterior distribution, with this method often referred to as Bayesian synthetic likelihood (BSL). We demonstrate that when the model is misspecified, the BSL posterior can be poorly behaved, placing significant posterior mass on values of the model parameters that do not represent the true features observed in the data. Theoretical results demonstrate that in misspecified models the BSL posterior can display a wide range of behaviours depending on the level of model misspecification, including being asymptotically non-Gaussian. Our results suggest that a recently proposed robust BSL approach can ameliorate this behavior and deliver reliable posterior inference under model misspecification. We document all theoretical results using a simple running example.