Bayesian Hierarchical Analysis

• Bayesian hierarchical analysis is a statistical framework that uses nested priors and hyperparameters to share information across different model levels.
• It integrates low- and high-fidelity simulations in surrogate modeling, achieving calibrated uncertainty and accurate predictions.
• The approach reduces computational cost through closed-form posteriors and recursive covariance inversion, benefiting complex engineering applications.

Bayesian hierarchical analysis is a statistical framework in which parameters are assigned probability distributions that themselves depend on further (“hyper”) parameters equipped with their own priors. This multi-level structure formalizes the notion of sharing information across units or levels of a system, supports principled quantification of uncertainty, and enables flexible modelling of complex systems, such as multi-fidelity computer simulations and stochastic physical processes. In the context of surrogate modeling for expensive computer codes, Bayesian hierarchical analysis provides a tractable and scalable methodology that incorporates both fast, low-fidelity approximations and accurate, high-fidelity models in a unified, uncertainty-aware predictive scheme (Gratiet, 2011). The hierarchical paradigm underpins closed-form inference, computational efficiency, and calibrated uncertainty for real-world engineering and scientific applications.

1. Gaussian Process–Based Hierarchical Modeling for Multi-Fidelity Codes

When modeling a complex physical system by computer codes of different fidelity, the outputs at successive fidelity levels are coupled using an autoregressive Gaussian process (GP) structure. In the basic two-level case, the high-fidelity code response $Z_2(x)$ is modeled as a linear transformation (“scale factor” or adjustment) of the lower-fidelity code response $Z_1(x)$ plus an independent GP discrepancy term: $Z_2(x) = \rho\,Z_1(x) + \delta(x)$ where $\rho$ is the scale factor and $\delta(x)$ is a zero-mean, independent GP with its own covariance structure. Both $Z_1(x)$ and $\delta(x)$ are assigned GP priors with specified correlation kernels (e.g., Matern 5/2). This model generalizes to $s$ levels recursively; each level $t$ is modeled as: $Z_t(x) = \rho_{t-1}Z_{t-1}(x) + \delta_t(x)$ with independent GPs $\delta_t(x)$. The hierarchical structure allows uncertainty from all code levels to be propagated explicitly.

2. Bayesian Parameter Inference and Closed-Form Solutions

The hierarchical approach treats all unknowns—including regression effects $\beta$, scale parameters $\{\rho\}$, variances $\{\sigma_t^2\}$, and correlation hyperparameters $\{\theta_t\}$—as random variables. Bayesian inference is then applied, with Gaussian, inverse-gamma, or non-informative (e.g., Jeffreys) priors chosen depending on prior knowledge and computational tractability.

Key results include the derivation of closed-form expressions for joint posteriors. For example, the joint posterior of the scale factor and regression bias for the high-fidelity response given the covariates and observed outputs is: $$[(\rho, \beta_2)\mid z_1, z_2, \sigma_2^2, \theta_2] \sim \mathcal{N}_{p_2+1}\left( \left[F^T R_2^{-1} F\right]^{-1} F^T R_2^{-1} z_2, \,\sigma_2^2 \left[F^T R_2^{-1} F\right]^{-1} \right)$$ where $F = [z_1(D_2) \;\; F_2(D_2)]$ and $R_2$ is the correlation matrix at level 2. Separate inverse-gamma posteriors for $\sigma^2$ parameters are also available. This enables analytic, non-MCMC estimation of adjustment parameters and credible intervals, quantifying uncertainty without intensive computation.

3. Computational Complexity and Efficient Covariance Inversion

Hierarchical multi-fidelity co-kriging poses scalability challenges due to the necessity of inverting large, block-structured covariance matrices. The paper leverages nested experimental designs (where $D_t \subset D_{t-1}$ for all $t$) and the autoregressive model hierarchy to construct a block inverse recursion: $$V_s^{-1} = \begin{pmatrix} V_{s-1}^{-1} + C & -B^T \ -B & D \end{pmatrix}$$ with explicit formulas for each block. The base case is $V_1^{-1} = \sigma_1^{-2}R_1^{-1}$. This approach reduces the inversion cost from $O((\sum_t n_t)^3)$ (the full stacked covariance) to $O(\sum_t n_t^3)$ (sum of individual level costs). Computational savings are dramatic when employing cross-validation or when the number of code levels is large.

4. Joint Posterior Representation and Prediction

A distinctive feature of this hierarchical Bayesian approach is the explicit joint factorization of the posterior: $$p(\beta_1, \beta_2, \rho, \sigma_1^2, \sigma_2^2\mid z_1, z_2) = p(\beta_1|\sigma_1^2, z_1)\; p(\beta_2, \rho | z_1, z_2, \sigma_2^2)\; p(\sigma_1^2 | z_1)\; p(\sigma_2^2 | z_1, z_2)$$ Each term is a normal or inverse-gamma density with closed-form solutions. Conditioning on fixed hyperparameters (typically estimated via restricted likelihood), full predictive distributions are obtained by integrating over all uncertainties, which can be executed by low-dimensional quadrature. This explicit representation removes the need for Markov Chain Monte Carlo to approximate the posterior, resulting in considerable gains in speed and reproducibility.

5. Practical Application: Thermodynamic Surrogate Modeling

A thermodynamic case paper involving a fluidized-bed process demonstrates the Bayesian hierarchical co-kriging framework. Four response types (experimental $T_\text{exp}$, high-fidelity $T_3$, medium-fidelity $T_2$, low-fidelity $T_1$) are considered, with the aim of predicting $T_\text{exp}$ given limited data. The three-level hierarchical model—using Matern 5/2 kernels, level-specific hyperparameters, and nested experimental designs—outperforms the two-level version:

• The three-level model attains significantly lower RMSE and improved Q2 compared to two-level or single-level models.
• Prediction intervals from the model are well calibrated, reflecting integrated parameter uncertainty.
• Hyperparameter estimates reveal that lower-fidelity simulations can be much smoother, depending strongly on input dimension reductions.

These results confirm the substantial advantages of Bayesian hierarchical analysis in accurately propagating the uncertainty and integrating heterogeneous computational data for surrogate modeling.

6. Methodological Impact and Generalizability

The analytic formulation for the multi-fidelity hierarchical Bayesian model permits rigorous uncertainty quantification in complex simulation environments, especially when observational data are very limited and computational codes are expensive to run at full fidelity. The closed-form posteriors enable calibration and prediction without the need for high-dimensional sampling, and the efficient covariance inversion scales to settings with many approximation levels.

Key methodological contributions include:

• Analytical joint posteriors for key adjustment parameters (scale, regression, variance) with uncertainty quantification.
• Computational reductions via recursive covariance matrix inversion benefiting nested designs.
• Generalizable framework accommodating any number of simulation fidelity levels.

This paradigm is directly extensible to surrogate modeling in engineering, uncertainty quantification in scientific computation, and design of computer experiments where multi-level code hierarchies are present.

---

Summary Table: Bayesian Hierarchical Co-kriging Framework

Component	Mathematical object / method	Role in framework
Model structure	Autoregressive GP hierarchy	Functional coupling of code fidelities
Parameter inference	Closed-form Bayesian posteriors	Uncertainty quantification for all parameters
Covariance inversion	Recursive block partitioning	Reduces cubic cost; enables scalability
Predictive uncertainty	Explicit joint posterior factorization	Direct integration; well-calibrated intervals, no need for MCMC
Real-world utility	Thermodynamic surrogate modeling example	Demonstrates improved accuracy, uncertainty quantification in sparse-data

The Bayesian hierarchical analysis of multi-fidelity models establishes an analytically tractable and computationally efficient approach that is particularly advantageous for engineering and physical sciences applications requiring accurate, uncertainty-aware surrogate modeling with limited high-fidelity data (Gratiet, 2011).