Bayesian Predictive Synthesis (BPS)

• Bayesian Predictive Synthesis (BPS) is a Bayesian framework for combining predictive distributions from multiple models by using latent state calibration to correct biases and learn dependencies.
• BPS extends traditional Bayesian model averaging and linear pooling by dynamically adjusting forecast weights and modeling interdependencies, leading to improved calibration and bias correction.
• BPS supports dynamic, multivariate, and sequential forecasting through scalable computational methods, delivering robust performance even during structural changes.

Bayesian Predictive Synthesis (BPS) is a general Bayesian framework for combining predictive distributions from multiple models, forecasters, or information sources in a way that goes beyond standard Bayesian model averaging (BMA) and beyond ordinary linear opinion pools. Its central idea is that a decision maker does not simply average forecasts; instead, the forecasts are treated as information about latent model states, and those latent states are then synthesized into a predictive distribution for the target outcome. In this sense, BPS expands standard Bayesian model uncertainty analysis into a broader theory of calibration, combination, dependence learning, and model-set incompleteness (1803.01984).

1. Historical emergence and conceptual basis

An important precursor is Bayesian Synthesis, in which several data analysts work independently on portions of a data set, elicit separate Bayesian models, update them on data not used in model building, and then combine them using weights derived from geometric means of pairwise Bayes factors. In that formulation, the resulting weighted mixture is a “hyper-model,” and the procedure was explicitly framed as a generalization of BMA in which the competing units are not just models but analysts plus their models. The same paper also introduced Convex Synthesis, a fixed-weight analogue. In modern BPS language, this earlier framework is conceptually very close: each analyst contributes a predictive component, and an overview layer learns how much trust to place in each component (Yu et al., 2011).

A later technical development clarified the formal status of the now-standard predictive synthesis formula. The 2024 report showed that modern BPS applications had been using the identity

$$p(y\mid \mathcal H)=\int \alpha(y\mid \mathbf x)\prod_{k=1}^K h_k(x_k)\,d\mathbf x$$

without an exact proof in the multi-agent, general-density setting, and then supplied two proofs. The report traced the formula back to discrete-event coherence results associated with Genest–Schervish and West–Crosse, and showed that the continuous-density identity arises as the limit of those discrete synthesis constructions. This also clarified a recurrent misconception: the product form is part of the predictive representation and is not, by itself, a substantive assumption that the agents are independent in reality (Masuda et al., 2024).

2. Formal representation and coherence

In its canonical form, BPS begins with predictive densities

$h_j(y), \qquad j=1,\dots,J,$

supplied by models $M_j$. BPS posits a latent state vector

$\mathbf{x}=(x_1,\dots,x_J)',$

where each $x_j$ represents the latent forecast contribution of model $M_j$. If

$h(\mathbf{x})=\prod_{j=1}^J h_j(x_j),$

then the BPS posterior predictive density is

$$p(y\mid \mathcal H)=\int \alpha(y\mid \mathbf{x})\, h(\mathbf{x})\, d\mathbf{x},$$

where $\mathcal H=\{h_1,\ldots,h_J\}$ and $\alpha(y\mid \mathbf{x})$ is the synthesis function. This is the central BPS representation: the models are not just weighted directly; rather, latent predictive draws are mapped through $h_j(y), \qquad j=1,\dots,J,$0 to produce the combined forecast (1803.01984).

A consistency condition links the decision maker’s prior predictive to the synthesis model: $h_j(y), \qquad j=1,\dots,J,$1 Updating from $h_j(y), \qquad j=1,\dots,J,$2 to the realized $h_j(y), \qquad j=1,\dots,J,$3 is a Jeffrey-style Bayesian update, not a standard likelihood-based Bayes update. This makes BPS an overview operator for predictive information rather than a conventional posterior obtained from a single fully specified joint likelihood. The synthesis function $h_j(y), \qquad j=1,\dots,J,$4 is therefore the critical modeling object: it expresses how the decision maker calibrates, corrects, and combines latent agent states (1803.01984).

3. Pooling, calibration, dependence, and model incompleteness

BPS contains BMA and linear pooling as special cases. If the synthesis function takes the form

$h_j(y), \qquad j=1,\dots,J,$5

with $h_j(y), \qquad j=1,\dots,J,$6 and $h_j(y), \qquad j=1,\dots,J,$7, then

$h_j(y), \qquad j=1,\dots,J,$8

This is a linear pool, and when $h_j(y), \qquad j=1,\dots,J,$9 it reduces to the standard linear pool. BMA arises as a further special case when the weights are interpreted as posterior model probabilities under a standard model-space assumption. BPS nonetheless generalizes these constructions by allowing a baseline model $M_j$0, allowing the combination to depend on latent forecast values, allowing flexible calibration and bias correction, allowing dependence across model forecasts, and allowing dynamic extensions (1803.01984).

The generalized linear pool

$M_j$1

is especially important. In BPS, outcome-dependent weights receive a Bayesian generative justification rather than functioning as heuristic devices. For the synthesis function

$M_j$2

one obtains

$M_j$3

with

$M_j$4

Thus BPS does not merely pool forecasts; it recalibrates them into transformed densities $M_j$5. The framework explicitly supports bias correction and dispersion adjustment, and it includes specific constructions such as Gaussian weights, Gaussian well weights, and cross-model weights based on

$M_j$6

so that synthesis can reward agreement, penalize herding, discount outliers, or favor diversity, depending on the substantive forecasting goal (1803.01984).

Dependence is handled through the full latent forecast vector. In the more general form

$M_j$7

the weights may depend on all of $M_j$8, not only on a single $M_j$9. This allows BPS to down-weight redundant or herding models, recognize consensus among models, discount outlier forecasts, and adapt weights according to how a forecast compares with the rest of the set. At the same time, the factorization $\mathbf{x}=(x_1,\dots,x_J)',$0 does not mean the forecasts are substantively independent; the prior structure $\mathbf{x}=(x_1,\dots,x_J)',$1 can encode dependence arbitrarily. The baseline density $\mathbf{x}=(x_1,\dots,x_J)',$2 addresses the $\mathbf{x}=(x_1,\dots,x_J)',$3-open problem by providing a principled “safe haven” when the available models are unreliable, too extreme, or simply incomplete (1803.01984).

A particularly important parametric instance is linear regression synthesis,

$\mathbf{x}=(x_1,\dots,x_J)',$4

in which $\mathbf{x}=(x_1,\dots,x_J)',$5 corrects overall bias, $\mathbf{x}=(x_1,\dots,x_J)',$6 captures cross-model dependence and relative influence, and $\mathbf{x}=(x_1,\dots,x_J)',$7 captures residual uncertainty after synthesis. This representation is often the starting point for dynamic BPS models (1803.01984).

4. Dynamic, multivariate, and sequential formulations

Dynamic BPS extends the static synthesis identity to sequential forecasting. In the univariate time-series formulation,

$\mathbf{x}=(x_1,\dots,x_J)',$8

and the main implementation uses a dynamic linear model synthesis function,

$\mathbf{x}=(x_1,\dots,x_J)',$9

The state evolves as

$x_j$0

with discount-factor evolution for both coefficients and volatility. In this construction, the intercept captures time-varying bias, the coefficients capture time-varying calibration and interdependence, and the model naturally supports horizon-specific calibration through BPS$x_j$1, where synthesis is trained directly on $x_j$2-step-ahead targets rather than extrapolated from 1-step performance (McAlinn et al., 2016).

The multivariate extension specializes the synthesis function to a dynamic latent factor or seemingly unrelated regression structure,

$x_j$3

Each outcome series has its own intercept and agent coefficients, while the equations are coupled through the full residual covariance matrix $x_j$4. This permits series-specific bias and calibration effects: one agent can be good for inflation and poor for wages, and BPS can accommodate that heterogeneity. The multivariate framework also uses posterior correlations of latent agent states and a Kullback–Leibler divergence between prior and posterior agent forecasts to summarize time-varying recalibration, dependence, and redundancy (McAlinn et al., 2017).

Sequential BPS has also been reformulated for online computation. In the linear-Gaussian case, the synthesis function is written as

$x_j$5

and the final synthesized prediction averages over posterior uncertainty in the calibration state $x_j$6. In this formulation, the intercept is crucial for bias correction, while the regression-like coefficients on the agent forecasts are not constrained to be nonnegative or to sum to one; they can be negative, can exceed one, and can adaptively correct location and uncertainty. A custom Rao-Blackwellized particle filter, supplemented by periodic MCMC interventions, was developed to make sequential DBPS fast enough for online use and to improve robustness under sudden distribution shifts. The same work introduced discount-factor averaging based on power-discounted likelihoods to improve adaptation during structural change (Masuda et al., 2023).

5. Computation and empirical behavior

BPS computation has been developed through several simulation-based and filtering-based schemes. For dynamic mixture-based BPS, the posterior is not analytically tractable, and implementations use Monte Carlo sampling, a Gibbs sampler with latent component indicators, and variational Bayes to refit normal-inverse-Wishart and Dirichlet approximations after each update. Dynamic linear BPS implementations reduce to state-space models conditional on the latent agent states, so forward filtering backward sampling is available for the synthesis states, while the latent states can be sampled directly under Gaussian agent densities or through scale-mixture augmentation when the agent densities are Student-$x_j$7. In the multivariate case, three-block Gibbs sampling is used for dynamic coefficients, volatilities, and latent agent states; in the sequential online case, Rao–Blackwellization collapses the particle weights to a predictive density term, and effective sample size diagnostics trigger periodic MCMC refreshes when particle degeneracy becomes severe (1803.01984, McAlinn et al., 2016, McAlinn et al., 2017, Masuda et al., 2023).

Empirical studies emphasize both predictive gains and learned structure. In the euro/USD foreign-exchange example, using daily euro log-prices against the US dollar in the second half of 2016 and combining forecasts from TVAR(1), TVAR(2), TVAR(5), and a locally linear dynamic linear model at a 5-trading-day horizon, BPS learned time-varying biases, learned cross-model dependence patterns, adapted around the U.S. presidential election when volatility rose, and used the baseline density when model forecasts became too extreme; it achieved the best RMSE and log score among the methods considered. In quarterly U.S. inflation forecasting, dynamic BPS improved forecast accuracy over all individual agents and pooling benchmarks, while horizon-specific BPS(4) outperformed direct extrapolation for 4-quarter-ahead forecasts. In monthly multivariate U.S. macroeconomic forecasting, BPS improved 1-step point forecasting for 5 of the 6 series, and BPS$x_j$8 and BPS$x_j$9 outperformed all the agent models for almost all series, with gains increasing with horizon. In sequential U.S. inflation forecasting through the 2020–2022 inflation burst, fixed-discount DBPS adapted slowly, whereas discount-factor averaging selected more adaptive discount factors during the burst period, yielding wider predictive intervals and better predictive scores in 2020–2022 (1803.01984, McAlinn et al., 2016, McAlinn et al., 2017, Masuda et al., 2023).

6. Specialized extensions and adjacent frameworks

Several specialized variants retain the BPS synthesis logic while altering the synthesis function or the inferential target. Bayesian Spatial Predictive Synthesis (BSPS) makes the synthesis function spatially varying: $M_j$0 with Gaussian process priors on the coefficient surfaces. Under the paper’s assumptions, this form is the best approximation of the data generating process given the predictive variables, and the resulting predictive distribution is exact minimax in Kullback–Leibler risk. In applied work, BSPS learned heterogeneous model performance across space, delivered improved predictive accuracy and uncertainty quantification, and remained interpretable through posterior distributions over spatially varying coefficients. A different extension replaced the usual parametric synthesis function with a regression-tree-based nonparametric version, BPS-RT, in which tree or BART structures define prior means for static and dynamic weights as functions of weight modifiers. In euro-area GDP growth and U.S. inflation applications, BPS-RT improved density forecast accuracy, especially in CRPS and tail-weighted CRPS, while retaining an interpretable weight structure (Cabel et al., 2022, Chernis et al., 2023).

Another line of work shifts the focus from prediction to decisions. Bayesian Predictive Decision Synthesis (BPDS) extends BPS by defining decision-dependent synthesis weights of the form

$M_j$1

where the score function reflects model-specific decision outcomes. This construction is justified through entropic tilting and generalized Bayesian updating, and it was developed in applications such as optimal design for regression prediction and sequential portfolio decisions. Predictive Decision Synthesis for portfolios used a closely related relaxed entropic tilting framework, with model-specific portfolio scores and a baseline model to represent the possibility that all models are wrong. General Bayesian Predictive Synthesis (GBPS) takes a different route by defining a loss-based generalized posterior

$M_j$2

so that ensemble weights are learned through loss minimization rather than parameter estimation in the usual predictive sense (Tallman et al., 2022, Tallman et al., 2024, Kato, 2024).

Quantile-focused synthesis extends BPS beyond mean or full-density combination. Dynamic Bayesian regression quantile synthesis (DRQS) uses an asymmetric Laplace synthesis function,

$M_j$3

so that the conditional quantile is directly

$M_j$4

Its multivariate generalization, factor DRQS (FDRQS), introduces a time-varying latent factor structure for the synthesis weights, enabling the model to leverage cross-sectional dependencies and shared information across multiple time series simultaneously. In empirical applications to U.S. inflation and global GDP growth, these models improved quantile forecasting, and FDRQS showed superior resilience during periods of extreme economic stress such as the COVID-19 pandemic by adaptively rebalancing agent contributions and capturing emergent global dependencies (Kobayashi et al., 12 Mar 2026).

7. Graphon and network generalizations

BPS has also been extended from scalar or vector outcomes to random networks. At the graphon level, agent models induce candidate graphons $M_j$5, and synthesis is defined through the linear combination

$M_j$6

In this setting, BPS corresponds to the $M_j$7 projection of the true graphon onto the agent span, and least-squares graphon-BPS achieves the minimax $M_j$8 rate over that span. The same theory establishes an oracle inequality that separates approximation error from estimation error, shows that single-agent selection is uniformly inconsistent on a nontrivial subset of the convex hull of the agents, and formalizes a “combination beats components” phenomenon. It also proves Lipschitz-type bounds that transfer graphon error into error for edge density, degree distributions, subgraph densities, clustering coefficients, and giant component phase transitions, and shows that entropic tilting keeps exponential random graph model agents within their family under log-linear tilting with Kullback–Leibler optimal moment calibration (Papamichalis et al., 21 Dec 2025).

Dynamic network BPS treats structural mechanisms such as communities, latent geometry, degree heterogeneity, and local link heuristics as agents that forecast dyad-level edge probabilities. These are combined through an affine logit pooling layer with a sparse-safe parametrization separating density from structure. A major theoretical result is that a single graph snapshot can identify and estimate the synthesis weights, with rate

$M_j$9

because the effective information scale is the number of dyads rather than the length of the time series. The same theory gives a sharp distinguishability threshold

$h(\mathbf{x})=\prod_{j=1}^J h_j(x_j),$0

shows that changes in the active mechanism can be tracked at an optimal per-switch cost, and proves that static link prediction is the $h(\mathbf{x})=\prod_{j=1}^J h_j(x_j),$1 special case. Under misspecification, BMA collapses to the single KL-closest agent, whereas BPS converges to the log-odds projection onto the mechanism class and can strictly dominate BMA whenever that affine-logit projection improves on every individual agent. Across these network formulations, the defining BPS characteristic remains unchanged: predictive combination is treated as synthesis over latent agent structures rather than as direct algebraic averaging (Papamichalis et al., 18 Jun 2026).