Information-Based Calibration of Uncertainty Quantification in Product-of-Experts Gaussian Process Models

Yean Hoon Ong; Paolo Barucca; Wei Pan; Jun Wang

Figure 1. Illustration of posterior predictions produced by a single global GP model (GP-glo) and a product-of-experts GP model (GP-pro) for the one-dimensional Ackley function. The white dots represent the training data points. The plot on the left shows the posterior prediction of the GP-glo model, with the shaded region indicating the standard deviation at each input value

x

. The plot on the right shows a GP-pro model comprising two local experts. The balltree assignment method is used to allocate training data points to the local GP models. The green and red curves, together with their corresponding shaded regions, represent the posterior predictions and standard deviations of the two local GP models. The GP-pro model employs an entropy-weighted aggregation scheme to combine the predictions of the local GP models. The blue curve and shaded region represent the aggregated posterior prediction and its associated standard deviation. Compared with the GP-glo model, the GP-pro model exhibits overestimated posterior variance, as evidenced by the wider blue uncertainty intervals. Two line graphs (left and right) show the Ackley function values, with the y-axis ranging from −4 to 2 and the x-axis representing input values from −10 to 10 in increments of 5. The left graph presents the function values and their corresponding confidence intervals predicted by a single global Gaussian process (GP) model. The right graph presents the function values and their corresponding confidence intervals predicted by a product of two local GP models. The confidence intervals produced by the product of local GP models are noticeably wider than those of the single global GP model, illustrating the issue of variance overestimation that can arise in product-of-Gaussian-process models.

Dataset	Size	Dim	GP-pro-bt-ent			GP-pro-full			GP-full
			RMSE $\downarrow$	NLL $\downarrow$	ENCE $\downarrow$	RMSE $\downarrow$	NLL $\downarrow$	ENCE $\downarrow$	RMSE $\downarrow$	NLL $\downarrow$	ENCE $\downarrow$
Ackley	500	1	0.300	0.267	0.296	0.290	0.204	0.252	0.290	0.204	0.252
Ackley	5000	10	0.487	0.775	0.245	0.468	0.663	0.057	0.468	0.663	0.057
Levy	5000	10	0.681	1.021	0.160	0.653	0.982	0.137	0.653	0.982	0.137
Rastrigin	5000	10	0.745	1.138	0.127	0.713	1.081	0.039	0.713	1.081	0.039
Average $\downarrow$			0.553	0.800	0.207	0.531	0.732	0.121	0.531	0.732	0.121

Table 1. Average RMSE, NLL, and ENCE results for the GP-pro, GP-pro-full, and GP-full models on four synthetic functions. GP-pro-full and GP-full yield identical results, as expected, since GP-pro-full trains each expert on the entire training dataset, resulting in identical posteriors whose product recovers the global GP solution. GP-pro exhibits higher NLL and ENCE than GP-pro-full and GP-full, indicating variance overestimation in the GP-pro model.

Figure 2. Reliability diagrams for the GP-pro, GP-pro-full, and GP-full models on four synthetic functions. Well-calibrated models yield curves closer to the ideal diagonal line. The curves of GP-pro-full overlap with those of GP-full, as expected. For the

10d

Levy function, GP-full also exhibits noticeable deviation from the diagonal line, which can be attributed to kernel misspecification and the high-dimensional, nonstationary nature of the function. In contrast, GP-pro shows consistently larger deviations across all functions, indicating systematic overestimation of posterior variances in the GP-pro models. Four line graphs, labelled (a), (b), (c), and (d), show reliability diagrams for three types of Gaussian process models evaluated on the one-dimensional Ackley function, the 10-dimensional Ackley function, the 10-dimensional Levy function, and the 10-dimensional Rastrigin function, respectively. In each graph, the horizontal axis represents the Root Mean Variance (RMV), ranging from 0 to 1 in increments of 0.2, and the vertical axis represents the Root Mean Square Error (RMSE), also ranging from 0 to 1 in increments of 0.2. Each graph contains four lines: an ideal dashed diagonal line extending from the lower-left corner to the upper-right corner, and three model-specific lines corresponding to (i) a product-of-experts Gaussian process model, (ii) a product-of-experts Gaussian process model in which each local model was trained on the entire training dataset, and (iii) a single global Gaussian process model. Across all four graphs, the single global Gaussian process model produces a line that is closest to the ideal diagonal, indicating better calibration. In contrast, the product-of-experts Gaussian process model produces a line that deviates noticeably from the diagonal, reflecting poorer calibration.

Dataset	Size	Dim	GP-pro	GP-pro-full	GP-full
			training, test time	training, test time	training, test time
Ackley	500	1	0.666, 0.018	0.677, 0.021	0.569, 0.013
Ackley	5000	10	3.101, 1.425	963.947, 13.745	47.905, 0.885
Levy	5000	10	3.380, 1.450	982.690, 14.161	46.986, 0.892
Rastrigin	5000	10	3.119, 1.452	972.641, 14.850	50.469, 1.053
Average $\downarrow$			2.566, 1.086	729.988, 10.694	36.482, 0.711

Table 2. Average training and test times (in seconds) across four synthetic functions for the GP-pro, GP-pro-full, and GP-full models.

Figure 3. Illustration of posterior predictions before and after applying the calibration ratio to GP-pro-c models comprising two local experts for the one-dimensional Ackley function. The GP-pro-c model on the left uses the random assignment method, whereas the model on the right uses balltree assignment. Both GP-pro-c models employ an entropy-weighted aggregation scheme. The posterior predictions of the two local GP models are represented by the green and red curves, with the corresponding light green and light red shaded regions indicating the standard deviation at each input value

x

. The blue curve and shaded region represent the aggregated posterior prediction. The darker green, red, and blue shaded regions indicate the reduced standard deviations obtained after applying the calibration ratio. Two line graphs, shown side by side, illustrate posterior predictions for the one-dimensional Ackley function. In both graphs, the y-axis shows function values ranging from -4 to 2, and the x-axis shows input values from -10 to 10 in increments of 2. Each graph contains three prediction curves and their corresponding uncertainty intervals. The first and second curves represent the posterior mean predictions of two local Gaussian process models, together with their uncertainty intervals before and after applying the calibration ratio. The third curve represents the aggregated posterior prediction obtained by combining the outputs of the two local Gaussian process models, along with its corresponding uncertainty interval. The darker shaded regions indicate the calibrated uncertainties after applying the calibration ratio, while the lighter shaded regions indicate the original uncertainties before calibration.

Method	Training	Test-time	Test-time	Notes
	complexity	mean	variance
Global GP	$\mathcal{O}(n^{3})$	$\mathcal{O}(n)$	$\mathcal{O}(n^{2})$	Full dataset, accurate variance
Vanilla GP-pro ( $M$ experts, $n_{m}$ points each)	$\mathcal{O}(Mn_{m}^{3})$	$\mathcal{O}(Mn_{m})$	$\mathcal{O}(Mn_{m}^{2})$	Distributed experts, standard PoE
GP-pro-c (proposed)	$\mathcal{O}(Mn_{m}^{3})$	$\mathcal{O}(Mn_{m})$	$\mathcal{O}(Mn_{m}^{2})$	Adds calibrated variance via $\alpha_{m}({\bm{x}})$ with negligible cost

Table 3. Computational complexity comparison.

Model	Data assignment method	Weighting scheme	Uncertainty calibration
bt-ent	balltree	entropy	no
bt-var	balltree	variance	no
bt-uni	balltree	uniform	no
kx-ent	$k$ -means	entropy	no
kx-var	$k$ -means	variance	no
kx-uni	$k$ -means	uniform	no

Table 4. Summary of the models used in the experiments for the comparison and evaluation of the GP-pro-c and GP-pro models.

Table 5. Average RMSE values across four synthetic functions and six benchmark datasets using balltree,

k

-means, and random assignment methods, as well as entropy-, variance-, and uniform-weighted schemes for the GP-pro-c and GP-pro models. GP-pro-c models with information-based calibration achieve RMSE values comparable to those of GP-pro models, indicating that information-based calibration improves uncertainty quantification without affecting predictive accuracy.

Table 6. Average NLL values across four synthetic functions and six benchmark datasets using balltree,

k

-means, and random assignment methods, as well as entropy-, variance-, and uniform-weighted schemes for the GP-pro-c and GP-pro models. GP-pro-c models with information-based calibration achieve lower NLL values than GP-pro models, demonstrating that information-based calibration improves uncertainty quantification without affecting predictive accuracy.

Table 7. Average ENCE values across four synthetic functions and six benchmark datasets using balltree,

k

-means, and random assignment methods, as well as entropy-, variance-, and uniform-weighted schemes for the GP-pro-c and GP-pro models. GP-pro-c models with information-based calibration achieve lower ENCE values than GP-pro models, demonstrating improved uncertainty quantification.

Figure 4. Reliability diagrams comparing GP-pro-c and GP-pro models using balltree,

k

-means, and random assignment methods across four synthetic functions and six regression datasets. GP-pro-c consistently yields curves closer to the ideal diagonal than GP-pro, indicating improved uncertainty calibration. The degree of deviation from the diagonal varies across datasets: smoother or lower-dimensional problems (e.g., Concrete, Airfoil, and Power) show near-diagonal curves, whereas noisier or higher-dimensional problems (e.g., Levy, Rosenbrock, Parkinsons, and Kin40K) show larger residual deviations. Ten line graphs show reliability diagrams for six variants of product-of-experts Gaussian process models evaluated on the 10-dimensional Ackley function, 10-dimensional Levy function, 10-dimensional Rastrigin function, 50-dimensional Rosenbrock function, the 8-dimensional Concrete dataset, the 5-dimensional Airfoil dataset, the 20-dimensional Parkinsons dataset, the 4-dimensional Power dataset, the 8-dimensional Kin40K dataset, and the 9-dimensional Protein dataset. \parIn each graph, the horizontal axis represents the Root Mean Variance (RMV), ranging from 0 to 1 in increments of 0.2, and the vertical axis represents the Root Mean Square Error (RMSE), also ranging from 0 to 1 in increments of 0.2. \parEach graph contains seven lines: an ideal dashed diagonal line extending from the lower-left corner to the upper-right corner, and six model-specific lines corresponding to (i) a product-of-experts Gaussian process model using the balltree assignment method, an entropy-weighted aggregation scheme, and uncertainty calibration; (ii) the same model without uncertainty calibration; (iii) a product-of-experts Gaussian-process model using the

k

-means assignment method, an entropy-weighted aggregation scheme, and uncertainty calibration; (iv) the same model without uncertainty calibration; (v) a product-of-experts Gaussian process model using random assignment, an entropy-weighted aggregation scheme, and uncertainty calibration; and (vi) the same model without uncertainty calibration. \parAcross all graphs, the calibrated product-of-experts Gaussian process model using the balltree assignment method generally produces a curve that lies closest to the ideal diagonal, indicating better uncertainty calibration than the other model variants. Curves corresponding to uncalibrated models tend to deviate more substantially from the diagonal.

Table 8. Average training and test times (in seconds) across four synthetic functions and six regression datasets using balltree,

k

-means, and random assignment methods, as well as entropy- and variance-weighted schemes (the uniform-weighted scheme is omitted for brevity) for the GP-pro-c and GP-pro models. Variants of GP-pro-c have training and test times similar to those of GP-pro.

Figure 5. Average NLL, RMSE, ENCE, and training time for GP-pro-c models using balltree,

k

-means, and random assignment methods with 100-600 data points per local GP, evaluated on four synthetic functions and four benchmark datasets. Random assignment yields lower ENCE but higher RMSE and NLL due to conservative, weakly correlated local GPs. Balltree and

k

-means preserve locality, yielding lower RMSE and NLL but higher ENCE.

k

-means performs poorly with small cluster sizes because small, high-dimensional clusters can be heterogeneous and poorly conditioned; performance improves sharply as cluster size increases.

k

-means can also incur slightly higher computational overhead due to uneven cluster sizes, whereas balltree and random assignments are more balanced. Overall, balltree assignment with 200-400 points per local GP provides a good trade-off between accuracy, calibration, and computational cost. Four line graphs labelled (a), (b), (c), and (d). \parIn all graphs, the horizontal x-axis shows the number of data points per local Gaussian process, ranging from 100 to 600 in increments of 100. Each graph contains three lines representing product-of-experts Gaussian process models constructed using the balltree,

k

-means, and random assignment methods, respectively. All models use the entropy-weighted scheme with uncertainty calibration. \parIn (a), the vertical y-axis shows the average Root Mean Square Error (RMSE), ranging from 0.40 to 0.60 in increments of 0.05. All three lines decrease consistently as the number of data points increases. The balltree line remains the lowest throughout. The

k

-means line is initially the highest at 100 data points, but from 200 data points onward it becomes the second lowest, overtaking the random assignment line. \parIn (b), the vertical y-axis shows the average Negative Log-Likelihood (NLL), ranging from 0.2 to 1.2 in increments of 0.2. All three lines decrease consistently. The balltree line remains the lowest throughout. The

k

-means line is initially the highest at 100 data points, but from 200 data points onward it becomes the second lowest, overtaking the random assignment line. \parIn (c), the vertical y-axis shows the average Expected Normalised Calibration Error (ENCE), ranging from 0.1 to 0.5 in increments of 0.1. The balltree and random assignment lines remain the lowest, at approximately 0.2 ENCE, throughout the graph. The

k

-means line begins at 0.45 ENCE with 100 data points, then decreases sharply to approximately 0.25 ENCE at 200 data points and remains the highest thereafter. \parIn (d), the vertical y-axis shows the average training time in seconds, ranging from 10 to 35 in increments of 5. All three lines increase exponentially as the number of data points grows. The balltree and random assignment lines remain the lowest throughout, while the

k

-means line remains the highest.

Dataset	Size	Dim	GP-pro-bt-ent-c			GP-pro-bt-ent			GP-pro-sm			GP-full
			RMSE $\downarrow$	NLL $\downarrow$	ENCE $\downarrow$	RMSE $\downarrow$	NLL $\downarrow$	ENCE $\downarrow$	RMSE $\downarrow$	NLL $\downarrow$	ENCE $\downarrow$	RMSE $\downarrow$	NLL $\downarrow$	ENCE $\downarrow$
Ackley	5000	10	0.487	0.742	0.190	0.487	0.775	0.245	0.527	0.801	0.115	0.468	0.663	0.057
Levy	5000	10	0.681	1.012	0.103	0.681	1.021	0.160	0.731	1.097	0.115	0.653	0.982	0.137
Rastrigin	5000	10	0.745	1.124	0.060	0.745	1.138	0.127	0.801	1.211	0.116	0.713	1.081	0.039
Rosenbrock	50000	50	0.610	0.918	0.116	0.610	0.929	0.167	0.711	1.106	0.162	N/A	N/A	N/A
Concrete	1030	8	0.336	0.287	0.162	0.336	0.287	0.166	0.341	0.298	0.186	0.335	0.282	0.190

Table 9. Average RMSE, NLL, and ENCE values across four synthetic functions and six regression benchmark datasets, comparing GP-pro-c models with uncertainty calibration, GP-pro models without calibration, GP-pro-sm models using softmax-weight calibration, and a single global GP (GP-full). GP-pro-c models achieve lower NLL and ENCE scores, indicating that the proposed information-based calibration maintains predictive accuracy while improving uncertainty calibration.

Figure 6. Reliability diagrams for the GP-pro-c, GP-pro, GP-pro-sm, and GP-full models on four synthetic benchmark functions and six regression datasets. The GP-full models consistently exhibit curves closest to the diagonal. The GP-pro-c models demonstrate improved calibration compared to GP-pro across most datasets. On the Kin40K dataset, the GP-pro-sm model exhibits a curve closest to the diagonal. However, on the Protein dataset, GP-pro-sm exhibits noticeable deviation from the diagonal, whereas GP-pro-c maintains more stable calibration behaviour across datasets. Ten line graphs showing reliability diagrams for six types of Gaussian process models on the 10-dimensional Ackley function, 10-dimensional Levy function, 10-dimensional Rastrigin function, 50-dimensional Rosenbrock function, 8-dimensional Concrete dataset, 5-dimensional Airfoil dataset, 20-dimensional Parkinsons dataset, 4-dimensional Power dataset, 8-dimensional Kin40K dataset, and 9-dimensional Protein dataset. \parOn each line graph, the horizontal x-axis shows Root Mean Variance ranging from 0 to 1 in increments of 0.2, and the vertical y-axis shows Root Mean Square Error ranging from 0 to 1 in increments of 0.2. \parEach plot contains five lines: a dashed diagonal line from the lower-left to the upper-right corner representing perfect calibration, and four curves corresponding to (i) a product-of-experts Gaussian process model using the balltree assignment method with entropy-weighted scheme and uncertainty calibration, (ii) the same model without uncertainty calibration, (iii) the same model with softmax weight calibration, and (iv) a single global Gaussian process model. \parIn the Ackley, Levy, Rastrigin, Concrete, Airfoil, and Parkinsons datasets, the single global Gaussian process model produces curves closest to the ideal diagonal line. \parAcross all datasets except Kin40K, the balltree-based product-of-experts Gaussian process model with entropy weighting and uncertainty calibration yields curves closer to the ideal diagonal line than the corresponding model with softmax weight calibration.

Dataset	Size	Dim	GP-pro-bt-ent-c	GP-pro-bt-ent	GP-pro-sm	GP-full
			training, test time	training, test time	training, test time	training, test time
Ackley	5000	10	3.101, 1.427	3.101, 1.425	3.101, 1.425	47.905, 0.885
Levy	5000	10	3.380, 1.452	3.380, 1.450	3.380, 1.450	46.986, 0.892
Rastrigin	5000	10	3.119, 1.454	3.119, 1.452	3.119, 1.453	50.469, 1.053
Rosenbrock	50000	50	32.522, 85.483	32.522, 85.375	32.522, 85.695	N/A
Concrete	1030	8	1.263, 0.037	1.263, 0.036	1.263, 0.037	1.717, 0.033

Table 10. Average training and test times (in seconds) across four synthetic functions and six regression benchmark datasets for GP-pro-c models with uncertainty calibration, GP-pro models without calibration, GP-pro-sm models using softmax-weight calibration, and a single global GP (GP-full).

Information-Based Calibration of Uncertainty Quantification in Product-of-Experts Gaussian Process Models

Abstract

Figures & tables

Appendix

Explore similar work

Bayesian Optimisation Using Product-of-Experts Gaussian Process Models with Uncertainty Calibration

Predictively Oriented Gaussian Process Posteriors

Wasserstein-type Gaussian Process Regressions for Input Measurement Uncertainty