# A Large-Scale Study of Probabilistic Calibration
in Neural Network Regression

Victor Dheur  Souhaib Ben Taieb

###### Abstract

Accurate probabilistic predictions are essential for optimal decision making.
While neural network miscalibration has been studied primarily in
classification, we investigate this in the less-explored domain of regression.
We conduct the largest empirical study to date to assess the probabilistic
calibration of neural networks. We also analyze the performance of
recalibration, conformal, and regularization methods to enhance probabilistic
calibration. Additionally, we introduce novel differentiable recalibration and
regularization methods, uncovering new insights into their effectiveness. Our
findings reveal that regularization methods offer a favorable tradeoff between
calibration and sharpness. Post-hoc methods exhibit superior probabilistic
calibration, which we attribute to the finite-sample coverage guarantee of
conformal prediction. Furthermore, we demonstrate that quantile recalibration
can be considered as a specific case of conformal prediction. Our study is
fully reproducible and implemented in a common code base for fair comparisons.

Calibration, Regression, Study, Regularization, Post-hoc, Neural Networks,
Deep Learning, Machine Learning, ICML

##  1 Introduction

Neural network predictions affect critical decisions in many applications,
including medical diagnostics and autonomous driving (Gulshan2016-kw, 24, 23).
However, effective decision making often requires accurate probabilistic
predictions (Gawlikowski2021-ty, 16, 1). For example, consider a probabilistic
regression model that produces 90% prediction intervals. An important property
would be that 90% of these prediction intervals contain the realizations.

For models that output a predictive distribution, probabilistic calibration is
an important property that states that all quantiles must be calibrated, i.e.,
the frequency of realizations below these quantiles must match the
corresponding quantile level. Additionally, predictive distributions should be
sufficiently sharp (i.e., concentrated around the realizations) and leverage
the information in the inputs.

In the classification setting, Guo2017-ow (25) found that common neural
architectures trained on image and text data were miscalibrated, sparking
increased interest in neural network calibration. In a follow-up study,
Minderer2021-xw (35) showed that more recent neural architectures demonstrate
improved calibration. However, there has been less research on calibration for
neural probabilistic regression models compared to classification. Therefore,
it remains uncertain whether the same results apply to the regression setting.
This paper addresses this gap by conducting a comprehensive study on
probabilistic calibration for regression using tabular data. We explore
various calibration methods, including quantile recalibration
(Kuleshov2018-tb, 33) and conformalized quantile regression (Romano2019-kp,
39). We also consider regularization methods, which have been shown to perform
well in the classification setting (Karandikar2021-ml, 31, 38, 47).

We make the following main contributions:

  1. 1.

We conduct the largest empirical study to date on probabilistic calibration of
neural regression models using 57 tabular datasets (Sections 4 and 6). We
consider multiple state-of-the-art calibration methods (Section 5), including
post-hoc recalibration, conformal prediction, and regularization methods, with
various scoring rules and predictive models.

  2. 2.

Building on quantile recalibration, we propose a new differentiable
calibration map using kernel density estimation, which provides improved
negative log-likelihood compared to baselines. We also introduce two new
regularization objectives based on the probabilistic calibration error
(Section 5).

  3. 3.

We show that quantile recalibration is a special case of conformal prediction,
providing an explanation for their superior performance in terms of
probabilistic calibration (Section 6).

##  2 Background

We consider a univariate regression problem where the target variable
$Y\in\mathcal{Y}$ depends on an input variable $X\in\mathcal{X}$, with
$\mathcal{Y}=\mathbb{R}$ representing the target space and $\mathcal{X}$
representing the input space. Our objective is to approximate the conditional
distribution $P_{Y\mid X}$ using training data
$\mathcal{D}=\Set{(X_{i},Y_{i})}_{i=1}^{N}$ where
$(X_{i},Y_{i})\overset{\text{i.i.d.}}{\sim}P\equiv P_{X}\times P_{Y\mid X}$.

A probabilistic predictor $F_{\theta}:\mathcal{X}\to\mathcal{F}$ is a function
parametrized by $\theta\in\Theta$ that maps an input $x\in\mathcal{X}$ to a
predictive cumulative distribution function (CDF) $F_{\theta}(\cdot\mid x)$ in
the space $\mathcal{F}$ of distributions over $\mathbb{R}$. Additionally,
given $x\in\mathcal{X}$, we denote the predictive quantile function (QF) by
$Q_{\theta}(\cdot\mid x)$, and probability density function (PDF) by
$f_{\theta}(\cdot\mid x)$. Similarly, the marginal CDF, QF, or PDF of a random
variable $R$ is denoted by $F_{R}$, $Q_{R}$, or $f_{R}$, respectively.

##### Probabilistic calibration.

Given an input $x\in\mathcal{X}$, the model $F_{\theta}$ is ideal if it
precisely matches the conditional distribution $P_{Y\mid X}$. However,
learning the ideal model based on finite data is not possible without
additional (strong) assumptions (Foygel_Barber2021-ig, 13). To avoid
additional assumptions, we can instead enforce certain desirable properties
that are attainable in practice and that a good or ideal forecaster should
exhibit. One such property is probabilistic calibration.

Let $Z=F_{\theta}(Y\mid X)\in[0,1]$ denote the probability integral transform
(PIT) of $Y$ conditional on $X$. The model $F_{\theta}$ is _probabilistically
calibrated_ (also known as PIT-calibrated) if $\forall\alpha\in[0,1]$,

$F_{Z}(\alpha)\doteq\text{Pr}(Z\leq\alpha)=\alpha.$ (1)

Let $U\in[0,1]$ be a uniform random variable independent of $Z$. The left and
right hand sides of 1 can be interpreted as the CDF of $Z$ and $U$,
respectively, as a function of $\alpha$. This shows that the uniformity of the
PIT is equivalent to probabilistic calibration (Dawid1984-gn, 8).

Since the ideal forecaster is probabilistically calibrated, we can require
this property from any competent forecaster. However, probabilistic
calibration, though necessary, is not sufficient for making accurate
probabilistic predictions. Additionally, as discussed by Gneiting2021-vu (18),
probabilistic calibration primarily addresses unconditional aspects of
predictive performance and is implied by more robust conditional notions of
calibration, such as auto-calibration.

##### Probabilistic calibration error.

The most common approach for evaluating probabilistic calibration is to
consider distances of the form
$\int_{0}^{1}|F_{Z}(\alpha)-F_{U}(\alpha)|^{p}d\alpha$ where $p>0$. The
particular cases of $p=1$ and $p=2$ are known as the 1-Wasserstein distance
and Cramér-von Mises distance, respectively. We denote the empirical CDF of
the PIT as
$\hat{F}_{Z}(\alpha)=\frac{1}{N}\sum_{i=1}^{N}\mathbbm{1}(Z_{i}\leq\alpha)$
where $Z_{i}=F_{\theta}(Y_{i}\mid X_{i})$ are PIT realizations. A common
approach to assess probabilistic calibration using Monte Carlo estimation is
to evaluate it at equidistant values $\alpha_{1}<\dots<\alpha_{M}$ as follows:

$\text{PCE}_{p}(F_{\theta},\mathcal{D})=\frac{1}{M}\sum_{j=1}^{M}\left|\alpha_{j}-\hat{F}_{Z}(\alpha_{j})\right|^{p}.$
(2)

This metric has been previously employed in literature such as Zhao2020-ze
(48, 49) with $p=1$, and Kuleshov2018-tb (33, 43) with $p=2$. It is important
to note that, unlike the classical definition of the $p$-norm, we do not
exponentiate $\frac{1}{p}$ in 2 to maintain consistency with prior literature.
In the subsequent sections, we focus our analysis on PCE1 and use the
abbreviation PCE for brevity.

One limitation of scalar metrics like PCE is their inability to provide
detailed information regarding calibration errors at individual quantile
levels, $\alpha_{1},\dots,\alpha_{M}$. Instead, PIT reliability diagrams offer
a visual assessment of probabilistic calibration across all quantile levels by
plotting the empirical CDF of the PIT $Z$. These diagrams display the right
side of 1 against its left side, with a perfectly calibrated model represented
by a diagonal line (in the asymptotic case). Figure 2 provides examples of
such reliability diagrams, which have been employed in studies by
Pinson2012-ba (37) and Kuleshov2018-tb (33).

##  3 Related Work

Figure 1: Multiple regression benchmark datasets with references. Datasets
inside parentheses have not been considered in this study. Full dataset names
are available in Table 1.

Post-hoc calibration approaches involve adjusting the predictions of a trained
model using a mapping learned from a separate calibration dataset. In the
context of classification, temperature scaling (Guo2017-ow, 25) is a simple
and effective method that adjusts predictive confidence while maintaining
accuracy. For regression tasks, quantile recalibration (Kuleshov2018-tb, 33)
aims to achieve probabilistic calibration. Conformal prediction (Vovk2020-pg,
45) is a general approach that provides prediction sets with a finite-sample
coverage guarantee. Notable methods applied with deep learning include
Conformal Quantile Regression (Romano2019-kp, 39) and Distributional Conformal
Prediction (Izbicki2020-ed, 29, 5). Furthermore, post-hoc approaches have also
been proposed for conditional notions of calibration (Song2019-bk, 40, 32).

Regularization approaches aim to improve calibration during training by
incorporating regularization techniques. Some methods, proposed by Zhao2020-ze
(48) and Feldman2021-ew (12), utilize regularization to target different
conditional notions of calibration based on the inputs. Zhou2021-hx (49)
introduced an alternative loss function involving the simultaneous training of
two neural networks, while Pearce2018-lo (36), Chung2021-rh (6), and
Thiagarajan2020-gm (42) proposed objectives that allow control over the
tradeoff between coverage and sharpness of prediction intervals. To our
knowledge, the only regularization objective specifically targeting
probabilistic calibration is quantile regularization (Utpala2020-nw, 43).
Other types of uncertainty quantification methods include ensembling
(Lakshminarayanan2017-zg, 34) and Bayesian methods (Jospin2022-zi, 30)

##  4 Are Neural Regression Models Probabilistically Calibrated?

We conduct an extensive empirical study to evaluate the probabilistic
calibration of neural regression models. To this end, we calculate the
probabilistic calibration error defined in (2) for various state-of-the-art
models across multiple benchmark datasets.

Benchmark datasets. We analyze a total of 57 datasets, including 27 from the
OpenML curated benchmark (Grinsztajn2022-nu, 21), 18 from the AutoML
Repository (Gijsbers2019-xk, 17), and 12 from the UCI Machine Learning
Repository (Dua2017-ut, 10).

These datasets are widely used in the evaluation of deep probabilistic models
and uncertainty quantification, as evidenced by previous studies such as
Fakoor2021-os (11, 6, 49, 43, 15).

Figure 1 provides an overview of the utilization of these datasets in previous
studies. To the best of our knowledge, our study represents the most
comprehensive assessment of probabilistic calibration for neural regression
models published to date.

Neural probabilistic regression models. We consider three state-of-the-art
neural probabilistic regression models. The first model predicts a parametric
distribution, where the parameters are obtained as outputs of a hypernetwork.
Previous studies have often focused on the Gaussian distribution
(Lakshminarayanan2017-zg, 34, 43, 48). To introduce more flexibility, we
consider a mixture of $K$ Gaussian distributions. Given an input
$x\in\mathcal{X}$, the hypernetwork parametrizes the means $\mu_{k}(x)$,
standard deviations $\sigma_{k}(x)$, and weights $w_{k}(x)$ for each component
$k=1,...,K$. To ensure positive standard deviations and that the mixture
weights form a discrete probability distribution, we use the Softplus and
Softmax activations, respectively. We have two variants of this model
depending on the scoring rule used for training: the negative log-likelihood
(NLL) or the continuous ranked probability score (CRPS). These models are
denoted as MIX-NLL and MIX-CRPS, respectively. It is worth noting that the
CRPS of a mixture of Gaussians has a closed-form expression (Grimit2006-dg,
20).

The second model predicts quantiles of the distribution (Tagasovska2019-ts,
41, 6, 12). Specifically, given an input $x\in\mathcal{X}$ and a quantile
level $\alpha\in[0,1]$, the model outputs a quantile $Q_{\theta}(\alpha\mid
x)$. The full quantile function can be obtained by evaluating the model at
multiple quantile levels. The model is trained by minimizing the quantile
score at multiple levels, which is asymptotically equivalent to minimizing the
CRPS (Bracher2021-ut, 4). We denote this model as SQR-CRPS, where SQR stands
for simultaneous quantile regression (Tagasovska2019-ts, 41).

Figure 2: The top row shows the PCE for different datasets with one standard
error (error bar). The bottom row gives examples of PIT reliability diagrams
for five datasets.

Experimental setup. We adopt the large-sized regime introduced by
Grinsztajn2022-nu (21), which involves truncating the datasets to a maximum of
50,000 examples. Among the 57 datasets, the number of examples ranges from 135
to 50,000, and the number of features ranges from 3 to 3,611111Please refer to
Appendix D, specifically Table 1, for a detailed summary of each dataset..
Each of the 57 datasets is divided into four sets: training (65%), validation
(10%), calibration (15%), and test (10%). We normalize the input $X$ and
target $Y$ using the mean and standard deviation from the training split. The
final predictions are then transformed back to the original scale. For our
neural network models, we employ the same fully-connected architecture as
previous studies conducted by Kuleshov2018-tb (33), Chung2021-rh (6), and
Fakoor2021-os (11). Further details regarding the model hyperparameters can be
found in Appendix C.

Results. In Figure 2, the first row displays the PCE (averaged over five
random train-validation-test splits) for MIX-NLL in blue on each of the 57
datasets. For comparison, the PCE of a perfectly calibrated model, i.e. with
uniformly distributed PITs, computed using $5\times 10^{4}$ simulated values
is shown in orange. The second row presents reliability diagrams for five
datasets, with 90% consistency bands as in Gneiting2023-at (19). Similar
information is provided for MIX-CRPS and SQR-CRPS in Figures 12 and 13 in
Section B.4, respectively. Additionally, reliability diagrams for all datasets
can be found in Figure 15 in Section B.6.

The analysis reveals that the (average) PCE is generally high across many
datasets, although there are significant variations between datasets. To test
the statistical significance of these results, $10^{4}$ samples were generated
from the sampling distribution of the average PCE under the null hypothesis of
probabilistic calibration. The resulting sampling distribution for all
datasets is presented in Section B.5.

By computing the p-value associated with a one-sided test in the upper tail of
the distribution (as illustrated in Section B.5), it was observed that most
datasets have a p-value of zero. This indicates that the average PCE obtained
for the considered model is higher than all the simulated average PCEs of the
probabilistically calibrated model. Applying a threshold of 0.01 and a Holm
correction for the 57 hypothesis tests, the null hypothesis is rejected for 11
datasets out of the 57.

Overall, the results indicate that the neural models considered in this study
are generally probabilistically miscalibrated on a significant number of
benchmark tabular datasets. In Section 6, we will further explore how
calibration methods can substantially improve the PCE of neural models.

##  5 Calibration Methods

We begin by discussing the three main approaches to calibration: quantile
recalibration, conformal prediction, and regularization-based calibration.
Following that, we introduce two novel variants of regularization-based
calibration.

Quantile recalibration and conformal prediction are post-hoc methods, meaning
they are applied after model training. These approaches utilize a separate
calibration dataset
$\mathcal{D}^{\prime}=\Set{(X^{\prime}_{i},Y^{\prime}_{i})}_{i=1}^{N^{\prime}}$,
where $(X^{\prime}_{i},Y^{\prime}_{i})\overset{\text{i.i.d.}}{\sim}P_{X,Y}$.
On the other hand, regularization-based calibration operates directly during
training and relies solely on the training data $\mathcal{D}$.

###  5.1 Quantile Recalibration

Quantile recalibration aims to transform a potentially miscalibrated CDF
$F_{\theta}$ into a probabilistically calibrated CDF
$F^{\prime}_{\theta}=F_{Z}\circ F_{\theta}$, using the calibration map $F_{Z}$
which represents the CDF of the PITs for $F_{\theta}$. For a given quantile
level $\alpha\in[0,1]$, the recalibrated CDF $F^{\prime}_{\theta}$ satisfies:

$\displaystyle\text{Pr}(F^{\prime}_{\theta}(Y\mid X)\leq\alpha)$
$\displaystyle=\text{Pr}(F_{\theta}(Y\mid X)\leq Q_{Z}(\alpha))$ (3)
$\displaystyle=F_{Z}(Q_{Z}(\alpha))$ (4) $\displaystyle=\alpha.$ (5)

In practice, $F_{Z}$ is not directly available and needs to be estimated from
data. Kuleshov2018-tb (33) proposed estimating it using isotonic regression,
while Utpala2020-nw (43) showed that computing the empirical CDF is an
equivalent and simpler method. Specifically, given a set of PIT values
$Z^{\prime}_{i}=F_{\theta}(Y^{\prime}_{i}\mid
X^{\prime}_{i}),i=1,\dots,N^{\prime}$, the calibration map $\phi^{\text{EMP}}$
is computed as:

$\phi^{\text{EMP}}(\alpha;\Set{Z^{\prime}_{i}}_{i=1}^{N^{\prime}})=\frac{1}{N^{\prime}}\sum_{i=1}^{N^{\prime}}\mathbbm{1}(Z^{\prime}_{i}\leq\alpha),$
(6)

where $\alpha\in[0,1]$.

Similarly to Utpala2020-nw (43), we also consider a linear calibration map
$\phi^{\text{LIN}}$, which is continuous, and corresponds to a linear
interpolation between the points
$\Set{(0,0),(Z^{\prime}_{(1)},\nicefrac{{1}}{{N^{\prime}+1}}),\dots,(Z^{\prime}_{(N^{\prime})},\nicefrac{{N^{\prime}}}{{N^{\prime}+1}}),(1,1)}$,
where $Z^{\prime}_{(k)}$ is the $k$th order statistic of
$Z^{\prime}_{1},\dots,Z^{\prime}_{N^{\prime}}$.

In addition, we propose a calibration map based on kernel density estimation
(KDE), denoted as $\phi^{\text{KDE}}$. This calibration map offers the
advantage of being differentiable and can lead to improved NLL performance.
The key idea is to use a relaxed approximation of the indicator function,
which allows us to make the PIT CDF 6 differentiable. Specifically, we compute

$\mathbbm{1}_{\tau}(a\leq b)=\sigma(\tau(b-a))\approx\mathbbm{1}(a\leq b),$

where $\tau>0$ is an hyperparameter and $\sigma(x)=\frac{1}{1+e^{-x}}$ denotes
the sigmoid function. The resulting smoothed empirical CDF is given by

$\phi^{\text{KDE}}(\alpha;\Set{Z^{\prime}_{i}}_{i=1}^{N^{\prime}})=\frac{1}{N^{\prime}}\sum_{i=1}^{N^{\prime}}\mathbbm{1}_{\tau}(Z^{\prime}_{i}\leq\alpha).$
(7)

This corresponds to estimating the CDF $F_{Z}$ using KDE based on $N^{\prime}$
realizations of $Z$ ($\Set{Z^{\prime}_{i}}_{i=1}^{N^{\prime}}$). Since
$\sigma$ is the CDF of the logistic distribution, we use the PDF of the
logistic distribution as the kernel in the KDE. Algorithm 1 summarises this
method.

Algorithm 1 Quantile recalibration

Input: Predictive CDF $F_{\theta}$ and
$\mathcal{D}^{\prime}=\Set{(X^{\prime}_{i},Y^{\prime}_{i})}_{i=1}^{N^{\prime}}$.

Compute $Z^{\prime}_{i}=F_{\theta}(Y^{\prime}_{i}\mid
X^{\prime}_{i})\quad(i=1,\dots,N^{\prime})$

Compute the calibration map $\phi$, either $\phi^{\text{EMP}}$,
$\phi^{\text{LIN}}$, or $\phi^{\text{KDE}}$

Return: Recalibrated CDF $F^{\prime}_{\theta}=\phi\circ F_{\theta}$.

###  5.2 Conformal Prediction

Let us assume the realizations of our calibration dataset
$\mathcal{D}^{\prime}$ are drawn exchangeably from $P_{X,Y}$222This is implied
by the common i.i.d. assumption.. Given a predictive model $M_{\theta}$ and a
coverage level $\alpha\in[0,1]$, (inductive) conformal prediction allows us to
construct a prediction set $C_{\alpha}(X)\subseteq\mathcal{Y}$ for any input
$X$, satisfying the property:

$\displaystyle\text{Pr}(Y\in C_{\alpha}(X))$
$\displaystyle=\frac{\lceil(N^{\prime}+1)\alpha\rceil}{N^{\prime}+1}$ (8)
$\displaystyle\approx\alpha.$ (9)

Conformal prediction achieves this by utilizing a conformity score
$s_{\theta}(Y\mid X)$, which intuitively quantifies the similarity between new
samples and previously observed samples. When the conformity score increases
with $Y$, an interval $C_{\alpha}(X)=(-\infty,s_{\theta}^{-1}(\alpha\mid X)]$
can be constructed, ensuring the conformal guarantee (8) at level $\alpha$.

Let $Q^{\prime}_{\theta}(\alpha\mid X)=s_{\theta}^{-1}(\alpha\mid X)$
represent the (revised) model obtained through conformal prediction from
$Q_{\theta}(\alpha\mid X)$. Under the assumption that $Q^{\prime}_{\theta}$ is
continuous and strictly increasing, the conformal guarantee implies that
$\text{Pr}(Y\leq Q^{\prime}_{\theta}(\alpha\mid X))\approx\alpha$, which
indicates approximate probabilistic calibration at level $\alpha$.

Conformalized Quantile Regression (Romano2019-kp, 39) is an example of a
conformal procedure, where the conformity score is defined as
$s_{\theta}(Y\mid X)=Y-Q_{\theta}(\alpha\mid X)$, representing the quantile
residual. Another example is Distributional Conformal Prediction
(Izbicki2020-ed, 29, 5), which employs the conformity score $s_{\theta}(Y\mid
X)=F_{\theta}(Y\mid X)$, referring to the PIT. Algorithm 2 provides a summary
of how to compute calibrated quantiles using inductive conformal prediction.

Algorithm 2 Calibrated quantiles with conformal prediction

Input: Trained model $M_{\theta}$,
$\mathcal{D}^{\prime}=\Set{(X^{\prime}_{i},Y^{\prime}_{i})}_{i=1}^{N^{\prime}}$,
strictly increasing conformity score $s$, quantile level $\alpha\in[0,1]$,
input $X$.

Compute $S_{i}=s_{\theta}(Y^{\prime}_{i}\mid
X^{\prime}_{i})\quad(i=1,\dots,N^{\prime})$

Compute $\hat{q}=S_{(\lceil(N^{\prime}+1)\alpha\rceil)}$ where $S_{(k)}$
denote the $k$th smallest value among
$\Set{S_{1},\dots,S_{N^{\prime}},+\infty}$

Return: Calibrated quantile $Q^{\prime}_{\theta}(\alpha\mid
X)=s_{\theta}^{-1}(\hat{q}\mid X)$

###  5.3 Regularization-based Calibration

Regularization-based calibration methods aim to enhance model calibration by
incorporating a regularization term into the training objective. Although
widely used in classification, there are relatively fewer methods specifically
designed for regression problems. In this section, we discuss two approaches:
quantile regularization (Utpala2020-nw, 43) and the truncation method
(Chung2021-rh, 6). The main steps of regularization-based calibration are
summarized in Algorithm 3.

Algorithm 3 Regularization-based calibration

Input: Model $M_{\theta}$, calibration regularizer $\mathcal{R}(\theta)$ and
tuning parameter $\lambda\geq 0$.

Compute
$\theta^{*}=\operatorname*{arg\,min}_{\theta\in\Theta}\mathcal{L}^{{}^{\prime}}(\theta;\mathcal{D})$
where

$\mathcal{L}^{{}^{\prime}}(\theta;\mathcal{D})=\nicefrac{{1}}{{N}}\sum_{i=1}^{N}\mathcal{L}(M_{\theta}(\cdot\mid
X_{i}),Y_{i})+\lambda\mathcal{R}(\theta;\mathcal{D})$

Return: Regularized model $M_{\theta^{*}}$

####  5.3.1 Quantile Regularization

The regularizer proposed by Utpala2020-nw (43) aims to measure the deviation
of the PIT variable $Z$ from a uniform distribution, which is characteristic
of a probabilistically calibrated model. This regularization penalty
encourages the selection of calibrated models during training.

The authors observed that the KL divergence between $Z$ and a uniform random
variable is equivalent to the negative differential entropy of $Z$, denoted as
$H(Z)$. To approximate $H(Z)$, they employed sample-spacing entropy estimation
(Vasicek1976-fa, 44), resulting in the following regularizer:

$\displaystyle\mathcal{R}_{\text{QR}}(\theta;\mathcal{D})$ (10)
$\displaystyle=\frac{1}{N-k}\sum_{i=1}^{N-k}\log\left[\frac{N+1}{k}(Z_{(i+k)}-Z_{(i)})\right]$
(11) $\displaystyle\approx H(Z),$ (12)

where $k$ is a hyperparameter satisfying $1\leq k\leq N$, and $Z_{(i)}$
represents the $i$th order statistic of $Z$.

To ensure differentiability during optimization, the authors employed a
differentiable relaxation technique called NeuralSort (Grover2019-iw, 22), as
sorting is a non-differentiable operation.

####  5.3.2 Truncation-based Calibration

The regularization approach introduced by Chung2021-rh (6), that we denote
Trunc, involves truncating the predictive distribution based on the current
level of calibration.

Given a quantile model $Q_{\theta}$, let
$\hat{F}_{Z}(\alpha)=\frac{1}{N}\sum_{i=1}^{N}\mathbbm{1}(Y_{i}\leq
Q_{\theta}(\alpha\mid X_{i}))$ be the estimated PIT CDF evaluated at $\alpha$
and $\rho(x,y)=(y-x)\mathbbm{1}(x<y)$. The regularization objective for level
$\alpha$ is defined as follows:

$\displaystyle\mathcal{R}_{\text{Trunc}}(\theta;\mathcal{D},\alpha)$ (13)
$\displaystyle=\begin{cases}\frac{1}{N}\sum_{i=1}^{N}\rho(Q_{\theta}(\alpha\mid
X_{i}),Y_{i})&\text{if $\hat{F}_{Z}(\alpha)<\alpha$}\\\
\frac{1}{N}\sum_{i=1}^{N}\rho(Y_{i},Q_{\theta}(\alpha\mid
X_{i}))&\text{otherwise}\end{cases}$ (14)

This regularization objective adjusts $\hat{F}_{Z}(\alpha)$ to match $\alpha$
by increasing it when $\hat{F}_{Z}(\alpha)<\alpha$, and vice versa. The final
regularization objective is computed by averaging
$\mathcal{R}_{\text{Trunc}}(\theta;\mathcal{D},\alpha)$ over multiple quantile
levels $\Set{\alpha_{j}}_{j=1}^{M}$:

$\mathcal{R}_{\text{Trunc}}(\theta;\mathcal{D})=\frac{1}{M}\sum_{j=1}^{M}\mathcal{R}_{\text{trunc}}(\theta;\mathcal{D},\alpha_{j}).$
(15)

It is worth noting that Chung2021-rh (6) combine the previous regularization
objective with a sharpness objective that penalizes the width between the
quantile predictions, given by
$\frac{1}{M}\sum_{j=1}^{M}\frac{1}{N}\sum_{i=1}^{N}|Q_{\theta}(\alpha_{j}\mid
X_{i})-Q_{\theta}(1-\alpha_{j}\mid X_{i})|$. Instead, we combine it with a
strictly proper scoring rule.

###  5.4 New Regularization-based Calibration Methods

Building upon the quantile calibration method discussed in Section 5.3.1, we
propose two new regularization objectives which compute a differentiable
$\text{PCE}_{p}$ using alternative statistical distances.

The first approach, named PCE-KDE, leverages the differentiable calibration
map $\phi^{\text{KDE}}$ 7 based on KDE. Given a set of quantile levels
$\Set{\alpha_{j}}_{j=1}^{M}$, the regularization objective is given by

$\displaystyle\mathcal{R}_{\text{PCE-
KDE}}(\theta;\mathcal{D})=\frac{1}{M}\sum_{j=1}^{M}\left|\alpha_{j}-\phi^{\text{KDE}}(\alpha_{j};\Set{Z_{i}}_{i=1}^{N})\right|^{p},$
(16)

where $p>0$. Note that $\mathcal{R}_{\text{PCE-KDE}}$ reduces to
$\text{PCE}_{p}$ in (2) when $\tau$ in (7) goes to $\infty$.

The second approach considers distances of the form
$\int_{0}^{1}|Q_{Z}(\alpha)-Q_{U}(\alpha)|^{p}d\alpha$, where $Q_{Z}$ and
$Q_{U}$ denote the quantile functions of the true and uniform distributions,
respectively. When $p=1$, this distance reduces to the 1-Wasserstein distance,
equivalent to $\int_{0}^{1}|F_{Z}(\alpha)-F_{U}(\alpha)|d\alpha$, which aligns
with PCE (see Proposition 1 in Section A.1).

By exploiting the fact that
$\mathbb{E}_{\,\\!\\!}\left[F_{Z}(Z_{(i)})\right]=\nicefrac{{i}}{{N+1}}$, we
approximate $Q_{Z}(\nicefrac{{i}}{{N+1}})$ using the $i$-th order statistic
$Z_{(i)}$. The regularization objective is given by

$\mathcal{R}_{\text{PCE-
Sort}}(\theta;\mathcal{D})=\frac{1}{N}\sum_{i=1}^{N}\left\lvert
Z_{(i)}-\frac{i}{N+1}\right\rvert^{p},$ (17)

where $p>0$. Differentiable relaxations to sorting, such as those proposed by
Blondel2020-sf (3) and Cuturi2019-dv (7), can be employed to obtain the order
statistics.

##  6 A Comparative Study of Probabilistic Calibration Methods

In continuation of the empirical study described in Section 4, we now proceed
to evaluate the performance of the probabilistic calibration methods outlined
in the previous section. Specifically, we apply eight distinct calibration
methods to the three neural regression models introduced in Section 4. These
methods are evenly divided into two categories: post-hoc methods and
regularization-based methods.

To assess the effectiveness of these calibration methods, we employ four
different evaluation metrics. The evaluation is conducted on a set of 57
datasets, utilizing the same experimental setup detailed in Section 4. To
ensure a fair and consistent comparison, all the methods have been implemented
within a unified codebase333The code can be accessed at the following GitHub
repository: https://github.com/Vekteur/probabilistic-calibration-study.

###  6.1 Experimental Setup

##### Base probabilistic models and calibration methods.

We consider the three probabilistic models presented in Section 4, namely MIX-
NLL, MIX-CRPS, and SQR-CRPS. For the MIX models, when applying quantile
recalibration, we transform the CDF using the empirical CDF estimator (Rec-
EMP), the linear estimator (Rec-LIN), or the KDE estimator (Rec-KDE). For SQR-
CRPS, we transform multiple quantiles using conformalized quantile regression
(CQR). For the three models, we consider the four regularization objectives
presented in Sections 5.3 and 5.4 (with $p=1$), namely
$\mathcal{R}_{\text{PCE-KDE}}$ (PCE-KDE), $\mathcal{R}_{\text{PCE-Sort}}$
(PCE-Sort), $\mathcal{R}_{\text{QR}}$ (QR), and $\mathcal{R}_{\text{Trunc}}$
(Trunc). PCE-Sort is only shown in the Appendix because it performs similarly
to PCE-KDE.

##### Metrics.

We measure the accuracy of the probabilistic predictions using NLL and CRPS.
For the SQR model, we estimate CRPS by averaging the quantile score at 64
equidistant quantile levels. Probabilistic calibration is measured using PCE,
defined in (1). Finally, we measure sharpness using the mean standard
deviation of the predictive distributions, denoted by STD.

##### Hyperparameters.

In our experiments, MIX-NLL and MIX-CRPS output a mixture of 3 Gaussians, and
SQR-CRPS outputs 64 quantiles. We justify the choice of these hyperparameters
in Appendix C. The hyperparameter $\tau$ of Rec-KDE and PCE-KDE is fixed at
100, which was found to perform well empirically. For regularization methods,
an important hyperparameter is the regularization factor $\lambda$. As
previously observed in classification (Karandikar2021-ml, 31), we found that
higher values of $\lambda$ tend to improve calibration but worsen NLL, CRPS,
and STD. Karandikar2021-ml (31) proposed to limit the loss in accuracy by a
maximum of 1%. We adopt a similar strategy by selecting $\lambda$ which
minimizes PCE with a maximum increase in CRPS of 10% in the validation set.
For each dataset, we select $\lambda$ in the set $\\{0,0.01,0.05,0.2,1,5\\}$,
which corresponds to various degrees of calibration regularization.

##### Comparison of multiple models over many datasets.

As in Karandikar2021-ml (31), and since NLL, CRPS and STD have different
scales across datasets, we report Cohen’s d, which is a standardized effect
size comparing the mean of one method (over 5 runs, in our case) against a
baseline. Values of $-0.8$ and $-2$ are considered large and huge,
respectively. Due to the heterogeneity of the datasets that we consider, the
performance of our models can vary widely across datasets. To visualize the
results, we show the distribution of Cohen’s d using letter-value plots
(Hofmann2011-rg, 26), which indicate the quantiles at levels
$\nicefrac{{1}}{{8}}$, $\nicefrac{{1}}{{4}}$, $\nicefrac{{1}}{{2}}$,
$\nicefrac{{3}}{{4}}$ and $\nicefrac{{7}}{{8}}$, as well as outliers. A median
value below zero indicates that the model improved the metric on more than
half the datasets.

In order to assess whether significant differences exist between different
methods, we follow the recommendations of Ismail_Fawaz2019-od (28), which are
based on (Demsar2006-ed, 9). First, we test for a significant difference among
model performances using the Friedman test (Friedman1940-vi, 14). Then, we use
the pairwise post-hoc analysis recommended by Benavoli2016-mh (2) using a
Wilcoxon signed-rank test (Wilcoxon1945-cp, 46) with Holm’s alpha correction
(Holm1979-rk, 27). The results of this procedure are shown using a critical
difference diagram (Demsar2006-ed, 9). The lower the rank (further to the
right), the better performance of a model. A thick horizontal line shows a
group of models whose performance is not significantly different, with a
significance level of 0.05.

###  6.2 Results

Figure 3 shows the letter-values plots for the Cohen’s d of PCE (top panel) as
well as the associated critical diagram (bottom panel), for all methods and
datasets. The reference model is MIX-NLL. The results with other models as
reference are available in Section B.1. Blue, green, and red colors are used
for the post-hoc methods, the regularization-based methods, and the base
models, respectively. The same information is given in Figures 4 and 5 for the
CRPS and the NLL, respectively.

##### Comparison of PCE.

As expected, Figure 3 shows that the PCE of calibration methods is improved
compared to the base models. Furthermore, independently of the base model, we
can see that post-hoc methods achieve significantly better PCE than
regularization methods. When comparing PCE-KDE with QR, we can see that there
is a significantly larger decrease in PCE with the MIX-CRPS base model
compared to MIX-NLL. Finally, both PCE-KDE and Trunc decrease PCE for SQR-
CRPS, without a significant difference between them.

((a)) Cohen’s d of PCE with respect to the MIX-NLL model  ((b)) Critical
difference diagram

Figure 3: Comparison of PCE with multiple base losses and calibration methods.

##### Comparison of CRPS.

While post-hoc methods outperform regularization methods in terms of PCE,
Figure 4 shows they have a higher CRPS (except for the SQR base model). This
can be explained by the fact that regularization methods prevent the CRPS from
increasing exceedingly due to the selection criterion for $\lambda$.

##### Comparison of NLL.

Figure 5 shows the importance of the calibration map. In fact, quantile
recalibration with a linear map significantly increases the NLL, while smooth
interpolation decreases PCE without a large increase in NLL. Note that we only
consider MIX models since we cannot compute the NLL for SQR.

((a)) Cohen’s d of CRPS with respect to the MIX-NLL model  ((b)) Critical
difference diagram

Figure 4: Comparison of CRPS with multiple base losses and calibration
methods.

##### On the choice of a calibration method.

If probabilistic calibration is critical to the application, our experiments
suggest that post-hoc methods such as quantile recalibration and conformal
prediction should be preferred. However, when we also want to control the CRPS
or the NLL, regularization methods can offer a better trade-off in terms of
calibration and sharpness. In fact, as shown in Figure 6 in Section B.1, when
the base model is MIX-NLL, all regularization methods provide a significant
improvement in probabilistic calibration without deteriorating the CRPS, NLL
or STD. For the MIX-CRPS model, Figure 7 shows that QR has limited impact on
CRPS and NLL, while providing better calibration. For the SQR-CRPS base model,
Figure 8 shows that the SQR-CRPS + CQR conformal method significantly
outperforms the Trunc and PCE-KDE regularization methods both in terms of PCE
and CRPS. Overall, Section B.1 suggests that MIX-NLL + PCE-KDE, MIX-CRPS + QR
and SQR-CRPS + CQR are good choices for practitioners aiming to improve PCE
without significantly impacting other aspects of the conditional distribution.
Finally, since both regularization and post-hoc methods are able to improve
calibration, we investigate whether a combination of these two methods can
lead to better performance. Figure 9 in Section B.2 shows that such a
combination does not significantly improve probabilistic calibration, with an
increase in CRPS and NLL. This indicates that practitioners should exercise
caution when applying regularization to a model that is already well-
calibrated.

((a)) Cohen’s d of NLL with respect to the MIX-NLL model  ((b)) Critical
difference diagram

Figure 5: Comparison of NLL with multiple base losses and calibration methods.

###  6.3 Link between Quantile Recalibration and Conformal Prediction

Conformal prediction methods are well-known for their finite-sample coverage
guarantee. Interestingly, a specific implementation of quantile recalibration
can be considered a special case of conformal prediction. This implies that
quantile recalibration can also provide a finite-sample coverage guarantee.
This observation could potentially explain why both methods, conformal
prediction and quantile recalibration, are effective in improving
probabilistic calibration.

######  Theorem 1.

Quantile recalibration is equivalent to Distributional Conformal Prediction
(DCP) of left intervals at each coverage level $\alpha\in[0,1]$. The
equivalence is obtained when the estimator of the calibration map is defined
by a slightly different estimator than the conventional one in (6), namely
$\phi_{\text{DCP}}(\alpha)=\frac{1}{N^{\prime}+1}\sum_{i=1}^{N^{\prime}}\mathbbm{1}(Z^{\prime}_{i}\leq\alpha)$.

###### Proof.

Given a predictive distribution $F_{\theta}$ learned from a training dataset
$\mathcal{D}=\Set{(X_{i},Y_{i})}_{i=1}^{N}$ where
$(X_{i},Y_{i})\overset{\text{i.i.d.}}{\sim}P_{X,Y}$, let
$Z^{\prime}_{i}=F_{\theta}(Y^{\prime}_{i}\mid X^{\prime}_{i})$ represent the
PIT values computed on a separate calibration dataset
$\mathcal{D}^{\prime}=\Set{(X^{\prime}_{i},Y^{\prime}_{i})}_{i=1}^{N^{\prime}}$,
where $(X^{\prime}_{i},Y^{\prime}_{i})\overset{\text{i.i.d.}}{\sim}P_{X,Y}$.

In the DCP approach, as outlined in Algorithm 2, the conformal scores are
given by the PIT values $Z^{\prime}_{i}$. DCP first computes the $\alpha$
empirical quantile of the scores as
$\hat{q}=Z^{\prime}_{(\lceil(N^{\prime}+1)\alpha\rceil)}$, where
$Z^{\prime}_{(k)}$ represents the $k$th smallest value among
$\Set{Z^{\prime}_{1},\dots,Z^{\prime}_{N^{\prime}},+\infty}$. Then, the
conformalized quantile is computed as $Q^{\prime}_{\theta}(\alpha\mid
X)=Q_{\theta}(\hat{q}\mid X)$, which corresponds to conformal prediction with
coverage $\alpha$ for the left interval
$(-\infty,Q^{\prime}_{\theta}(\alpha\mid X)]$.

Let us consider quantile recalibration with the calibration map $\phi$ in
Algorithm 1 given by
$\phi_{\text{DCP}}(\alpha)=\frac{1}{N^{\prime}+1}\sum_{i=1}^{N^{\prime}}\mathbbm{1}(Z^{\prime}_{i}\leq\alpha)$.
It computes a recalibrated CDF $F^{\prime}_{\theta}$ by composing the original
CDF $F_{\theta}$ with $\phi_{\text{DCP}}$, yielding $F^{\prime}_{\theta}(y\mid
X)=\phi_{\text{DCP}}(F_{\theta}(y\mid X))$.

We observe that $\phi_{\text{DCP}}$ is the CDF of a discrete random variable,
with
$\phi_{\text{DCP}}^{-1}(\alpha)=Z^{\prime}_{(\lceil(N^{\prime}+1)\alpha\rceil)}$
representing its empirical quantile function. Furthermore, the composition
$\phi_{\text{DCP}}\circ F_{\theta}(\cdot\mid X)$ acts as the inverse function
of $Q_{\theta}(\cdot\mid X)\circ\phi_{\text{DCP}}^{-1}$. As a result, both the
DCP approach and quantile recalibration yield QFs and CDFs that correspond to
the same underlying distribution.

∎

Quantile recalibration with other recalibration maps (e.g.,
$\phi^{\text{EMP}}$, $\phi^{\text{LIN}}$, or $\phi^{\text{KDE}}$) would
correspond to DCP where the empirical quantile $\hat{q}$ is selected using
other strategies which does not provide the exact conformal guarantee 8.

##  7 Conclusion

The observation that neural network classifiers tend to be miscalibrated
(Guo2017-ow, 25) has prompted the development of various approaches for
calibrating these models. In this paper, we present the largest empirical
study conducted to date on the probabilistic calibration of neural regression
models. Our study provides valuable insights into their performance and the
selection of calibration methods. Notably, we introduce a novel differentiable
calibration map based on kernel density estimation for quantile recalibration,
as well as two novel regularization objectives derived from the PCE.

Our study reveals that regularization methods can provide a favorable tradeoff
between calibration and sharpness. However, post-hoc methods demonstrate
superior performance in terms of PCE. We attribute this finding to the finite-
sample coverage guarantee offered by conformal prediction and demonstrate that
quantile recalibration can be viewed as a specific case of conformal
prediction.

Future investigations may extend the study of probabilistic calibration to
other models, such as tree-based models, and explore alternative notions of
calibration (Gneiting2021-vu, 18). Notably, distribution calibration
represents a promising direction, as it has inspired the development of
calibration methods (Song2019-bk, 40, 32).

## Acknowledgement

This work was supported by the Fonds de la Recherche Scientifique - FNRS under
Grants T.0011.21 and J.0011.20.

## References

  * (1) Moloud Abdar et al.  “A review of uncertainty quantification in deep learning: Techniques, applications and challenges”  In _An international journal on information fusion_ 76, 2021, pp. 243–297  DOI: 10.1016/j.inffus.2021.05.008
  * (2) Alessio Benavoli, Giorgio Corani and Francesca Mangili  “Should We Really Use Post-Hoc Tests Based on Mean-Ranks?”  In _Journal of machine learning research: JMLR_ 17.5, 2016, pp. 1–10  URL: https://jmlr.org/papers/v17/benavoli16a.html
  * (3) Mathieu Blondel, Olivier Teboul, Quentin Berthet and Josip Djolonga  “Fast Differentiable Sorting and Ranking”  In _Proceedings of the 37th International Conference on Machine Learning_ 119, Proceedings of Machine Learning Research  PMLR, 2020, pp. 950–959  URL: https://proceedings.mlr.press/v119/blondel20a.html
  * (4) Johannes Bracher, Evan L Ray, Tilmann Gneiting and Nicholas G Reich  “Evaluating epidemic forecasts in an interval format”  In _PLoS computational biology_ 17.2 journals.plos.org, 2021, pp. e1008618  DOI: 10.1371/journal.pcbi.1008618
  * (5) Victor Chernozhukov, Kaspar Wüthrich and Yinchu Zhu  “Distributional conformal prediction”  In _Proceedings of the National Academy of Sciences of the United States of America_ 118.48 National Acad Sciences, 2021  DOI: 10.1073/pnas.2107794118
  * (6) Youngseog Chung, Willie Neiswanger, Ian Char and Jeff Schneider  “Beyond Pinball Loss: Quantile Methods for Calibrated Uncertainty Quantification”  In _Advances in Neural Information Processing Systems_ 34 Curran Associates, Inc., 2021, pp. 10971–10984  URL: https://proceedings.neurips.cc/paper/2021/file/5b168fdba5ee5ea262cc2d4c0b457697-Paper.pdf
  * (7) Marco Cuturi, Olivier Teboul and Jean-Philippe Vert  “Differentiable Ranks and Sorting using Optimal Transport”, 2019  arXiv: https://proceedings.neurips.cc/paper/2019/file/d8c24ca8f23c562a5600876ca2a550ce-Paper.pdf
  * (8) A P Dawid  “Present Position and Potential Developments: Some Personal Views: Statistical Theory: The Prequential Approach”  In _Journal of the Royal Statistical Society. Series A_ 147.2 [Royal Statistical Society, Wiley], 1984, pp. 278–292  DOI: 10.2307/2981683
  * (9) Janez Demšar  “Statistical Comparisons of Classifiers over Multiple Data Sets”  In _Journal of machine learning research: JMLR_ 7 JMLR.org, 2006, pp. 1–30  URL: https://dl.acm.org/doi/10.5555/1248547.1248548
  * (10) Dheeru Dua and Casey Graff  “UCI Machine Learning Repository”, 2017  URL: http://archive.ics.uci.edu/ml
  * (11) Rasool Fakoor et al.  “Flexible Model Aggregation for Quantile Regression”, 2021  arXiv: http://arxiv.org/abs/2103.00083
  * (12) Feldman, Bates and Romano  “Improving conditional coverage via orthogonal quantile regression”  In _Advances in neural information processing systems_ proceedings.neurips.cc, 2021  URL: https://proceedings.neurips.cc/paper/2021/hash/1006ff12c465532f8c574aeaa4461b16-Abstract.html
  * (13) Rina Foygel Barber, Emmanuel J Candès, Aaditya Ramdas and Ryan J Tibshirani  “The limits of distribution-free conditional predictive inference”  In _Information and Inference: A Journal of the IMA_ 10.2 Oxford Academic, 2021, pp. 455–482  DOI: 10.1093/imaiai/iaaa017
  * (14) Milton Friedman  “A Comparison of Alternative Tests of Significance for the Problem of $m$ Rankings”  In _The Annals of Mathematical Statistics_ 11.1 Institute of Mathematical Statistics, 1940, pp. 86–92  DOI: 10.1214/aoms/1177731944
  * (15) Yarin Gal and Zoubin Ghahramani  “Dropout as a Bayesian Approximation: Representing Model Uncertainty in Deep Learning”  In _Proceedings of The 33rd International Conference on Machine Learning_ 48, Proceedings of Machine Learning Research  New York, USA: PMLR, 2016, pp. 1050–1059  URL: https://proceedings.mlr.press/v48/gal16.html
  * (16) Jakob Gawlikowski et al.  “A Survey of Uncertainty in Deep Neural Networks”, 2021  arXiv: http://arxiv.org/abs/2107.03342
  * (17) Pieter Gijsbers et al.  “An Open Source AutoML Benchmark”, 2019  arXiv: http://arxiv.org/abs/1907.00909
  * (18) Tilmann Gneiting and Johannes Resin  “Regression Diagnostics meets Forecast Evaluation: Conditional Calibration, Reliability Diagrams, and Coefficient of Determination”, 2021  arXiv: http://arxiv.org/abs/2108.03210
  * (19) Tilmann Gneiting et al.  “Model Diagnostics and Forecast Evaluation for Quantiles”  In _Annual Review of Statistics and Its Application_ Annual Reviews, 2023  DOI: 10.1146/annurev-statistics-032921-020240
  * (20) E P Grimit, T Gneiting, V J Berrocal and N A Johnson  “The continuous ranked probability score for circular variables and its application to mesoscale forecast ensemble verification”  In _Quarterly Journal of the Royal Meteorological Society_ 132.621C Wiley, 2006, pp. 2925–2942  DOI: 10.1256/qj.05.235
  * (21) Léo Grinsztajn, Edouard Oyallon and Gaël Varoquaux  “Why do tree-based models still outperform deep learning on tabular data?”, 2022  arXiv: http://arxiv.org/abs/2207.08815
  * (22) Aditya Grover, Eric Wang, Aaron Zweig and Stefano Ermon  “Stochastic Optimization of Sorting Networks via Continuous Relaxations”  In _International Conference on Learning Representations_ , 2019  URL: https://openreview.net/forum?id=H1eSS3CcKX
  * (23) Vitor Guizilini et al.  “3D Packing for Self-Supervised Monocular Depth Estimation”  In _IEEE Conference on Computer Vision and Pattern Recognition (CVPR)_ , 2020
  * (24) Varun Gulshan et al.  “Development and Validation of a Deep Learning Algorithm for Detection of Diabetic Retinopathy in Retinal Fundus Photographs”  In _JAMA: the journal of the American Medical Association_ 316.22, 2016, pp. 2402–2410  DOI: 10.1001/jama.2016.17216
  * (25) Chuan Guo, Geoff Pleiss, Yu Sun and Kilian Q Weinberger  “On Calibration of Modern Neural Networks”  In _Proceedings of the 34th International Conference on Machine Learning_ 70, Proceedings of Machine Learning Research  PMLR, 2017, pp. 1321–1330  URL: https://proceedings.mlr.press/v70/guo17a.html
  * (26) Heike Hofmann, Karen Kafadar and Hadley Wickham  “Letter-value plots: Boxplots for large data”, 2011
  * (27) Sture Holm  “A Simple Sequentially Rejective Multiple Test Procedure”  In _Scandinavian journal of statistics, theory and applications_ 6.2 [Board of the Foundation of the Scandinavian Journal of Statistics, Wiley], 1979, pp. 65–70  URL: http://www.jstor.org/stable/4615733
  * (28) Hassan Ismail Fawaz et al.  “Deep learning for time series classification: a review”  In _Data mining and knowledge discovery_ 33.4, 2019, pp. 917–963  DOI: 10.1007/s10618-019-00619-1
  * (29) Rafael Izbicki, Gilson Shimizu and Rafael Stern  “Flexible distribution-free conditional predictive bands using density estimators”  In _Proceedings of the Twenty Third International Conference on Artificial Intelligence and Statistics_ 108, Proceedings of Machine Learning Research  PMLR, 2020, pp. 3068–3077  URL: https://proceedings.mlr.press/v108/izbicki20a.html
  * (30) Laurent Valentin Jospin et al.  “Hands-On Bayesian Neural Networks—A Tutorial for Deep Learning Users”  In _IEEE Computational Intelligence Magazine_ 17.2 ieeexplore.ieee.org, 2022, pp. 29–48  DOI: 10.1109/MCI.2022.3155327
  * (31) Archit Karandikar et al.  “Soft Calibration Objectives for Neural Networks”, 2021  arXiv: http://arxiv.org/abs/2108.00106
  * (32) Volodymyr Kuleshov and Shachi Deshpande  “Calibrated and Sharp Uncertainties in Deep Learning via Density Estimation”  In _Proceedings of the 39th International Conference on Machine Learning_ 162, Proceedings of Machine Learning Research  PMLR, 2022, pp. 11683–11693  URL: https://proceedings.mlr.press/v162/kuleshov22a.html
  * (33) Volodymyr Kuleshov, Nathan Fenner and Stefano Ermon  “Accurate Uncertainties for Deep Learning Using Calibrated Regression”  In _Proceedings of the 35th International Conference on Machine Learning_ 80, Proceedings of Machine Learning Research  PMLR, 2018, pp. 2796–2804  URL: https://proceedings.mlr.press/v80/kuleshov18a.html
  * (34) Lakshminarayanan and Pritzel  “Simple and scalable predictive uncertainty estimation using deep ensembles”  In _Advances in neural information processing systems_ proceedings.neurips.cc, 2017  URL: https://proceedings.neurips.cc/paper/2017/hash/9ef2ed4b7fd2c810847ffa5fa85bce38-Abstract.html
  * (35) Matthias Minderer et al.  “Revisiting the Calibration of Modern Neural Networks”  In _Advances in Neural Information Processing Systems_ 34 Curran Associates, Inc., 2021, pp. 15682–15694  URL: https://proceedings.neurips.cc/paper/2021/file/8420d359404024567b5aefda1231af24-Paper.pdf
  * (36) Tim Pearce, Alexandra Brintrup, Mohamed Zaki and Andy Neely  “High-Quality Prediction Intervals for Deep Learning: A Distribution-Free, Ensembled Approach”  In _Proceedings of the 35th International Conference on Machine Learning_ 80, Proceedings of Machine Learning Research  PMLR, 2018, pp. 4075–4084  URL: https://proceedings.mlr.press/v80/pearce18a.html
  * (37) Pierre Pinson and Renate Hagedorn  “Verification of the ECMWF ensemble forecasts of wind speed against analyses and observations”  In _Meteorological Applications_ 19.4 Wiley, 2012, pp. 484–500  DOI: 10.1002/met.283
  * (38) Teodora Popordanoska, Raphael Sayer and Matthew B Blaschko  “A Consistent and Differentiable Lp Canonical Calibration Error Estimator”, 2022  URL: https://openreview.net/pdf?id=HMs5pxZq1If
  * (39) Romano and Patterson  “Conformalized quantile regression”  In _Advances in neural information processing systems_ proceedings.neurips.cc, 2019  URL: https://proceedings.neurips.cc/paper/2019/hash/5103c3584b063c431bd1268e9b5e76fb-Abstract.html
  * (40) Hao Song, Tom Diethe, Meelis Kull and Peter Flach  “Distribution calibration for regression”  In _Proceedings of the 36th International Conference on Machine Learning_ 97, Proceedings of Machine Learning Research  PMLR, 2019, pp. 5897–5906  URL: https://proceedings.mlr.press/v97/song19a.html
  * (41) Tagasovska and Lopez-Paz  “Single-model uncertainties for deep learning”  In _Advances in neural information processing systems_ proceedings.neurips.cc, 2019  URL: https://proceedings.neurips.cc/paper/2019/hash/73c03186765e199c116224b68adc5fa0-Abstract.html
  * (42) Jayaraman J Thiagarajan et al.  “Designing accurate emulators for scientific processes using calibration-driven deep models”  In _Nature communications_ 11.1, 2020, pp. 5622  DOI: 10.1038/s41467-020-19448-8
  * (43) Saiteja Utpala and Piyush Rai  “Quantile Regularization: Towards Implicit Calibration of Regression Models”, 2020  arXiv: http://arxiv.org/abs/2002.12860
  * (44) Oldrich Vasicek  “A Test for Normality Based on Sample Entropy”  In _Journal of the Royal Statistical Society. Series B, Statistical methodology_ 38.1 [Royal Statistical Society, Wiley], 1976, pp. 54–59  URL: http://www.jstor.org/stable/2984828
  * (45) Vladimir Vovk et al.  “Conformal calibrators”  In _Proceedings of the Ninth Symposium on Conformal and Probabilistic Prediction and Applications_ 128, Proceedings of Machine Learning Research  PMLR, 2020, pp. 84–99  URL: https://proceedings.mlr.press/v128/vovk20a.html
  * (46) Frank Wilcoxon  “Individual Comparisons by Ranking Methods”  In _Biometrics Bulletin_ 1.6 [International Biometric Society, Wiley], 1945, pp. 80–83  DOI: 10.2307/3001968
  * (47) Hee Suk Yoon et al.  “ESD: Expected Squared Difference as a Tuning-Free Trainable Calibration Measure”  In _The Eleventh International Conference on Learning Representations_ , 2023  URL: https://openreview.net/forum?id=bHW9njOSON
  * (48) Shengjia Zhao, Tengyu Ma and Stefano Ermon  “Individual Calibration with Randomized Forecasting”  In _Proceedings of the 37th International Conference on Machine Learning_ 119, Proceedings of Machine Learning Research  PMLR, 2020, pp. 11387–11397  URL: https://proceedings.mlr.press/v119/zhao20e.html
  * (49) Tianhui Zhou, Yitong Li, Yuan Wu and David Carlson  “Estimating Uncertainty Intervals from Collaborating Networks”  In _Journal of machine learning research: JMLR_ 22, 2021  URL: https://www.ncbi.nlm.nih.gov/pubmed/35754923

##  Appendix A Proofs

###  A.1 Integral of the Absolute Difference between CDFs or QFs

######  Proposition 1.

Let $F_{A},F_{B}:[0,1]\to[0,1]$ denote two strictly increasing CDFs of random
variables defined on $[0,1]$ with corresponding QFs $Q_{A}$ and $Q_{B}$. Then,

$\int_{0}^{1}|F_{A}(q)-F_{B}(q)|\,dq=\int_{0}^{1}|Q_{A}(p)-Q_{B}(p)|\,dp.$
(18)

###### Proof.

We define two functions $r,s:[0,1]\times[0,1]\to\Set{0,1}$ where

$\displaystyle r(q,p)$ $\displaystyle=\begin{cases}1&\text{if }F_{A}(q)\leq
p\leq F_{B}(q)\text{ or }F_{B}(q)\leq p\leq F_{A}(q)\\\
0&\text{otherwise},\end{cases}$ (19) $\displaystyle\text{and }s(q,p)$
$\displaystyle=\begin{cases}1&\text{if }Q_{A}(p)\leq q\leq Q_{B}(p)\text{ or
}Q_{B}(p)\leq q\leq Q_{A}(p)\\\ 0&\text{otherwise}.\end{cases}$ (20)

Let us show that $r$ and $s$ are equal. Considering $q\in[0,1]$ and
$p\in[0,1]$, we can write

$\displaystyle F_{A}(q)\leq p\leq F_{B}(q)$ (21) $\displaystyle\iff$
$\displaystyle(F_{A}(q)\leq p)\land(p\leq F_{B}(q))$ (22) $\displaystyle\iff$
$\displaystyle(q\leq Q_{A}(p))\land(Q_{B}(p)\leq q)$ (23) $\displaystyle\iff$
$\displaystyle Q_{B}(p)\leq q\leq Q_{A}(p),$ (24)

where 23 holds since both $F_{A}$ and $F_{B}$ are strictly increasing.

Similarly, $F_{B}(q)\leq p\leq F_{A}(q)\iff Q_{A}(p)\leq q\leq Q_{B}(p)$.
Hence $r(q,p)=1\iff s(q,p)=1$ and $r$ and $s$ are equal.

By Fubini’s theorem, we have

$\int_{0}^{1}\int_{0}^{1}r(q,p)\,dp\,dq=\int_{0}^{1}\int_{0}^{1}s(q,p)\,dq\,dp.$
(25)

Furthermore, upon evaluating the inner integrals, we obtain

$\displaystyle\int_{0}^{1}r(q,p)\,dp$
$\displaystyle=\begin{cases}\displaystyle\int_{F_{A}(q)}^{F_{B}(q)}1\,dp&\text{if
}F_{A}(q)\leq F_{B}(q)\\\
\displaystyle\int_{F_{B}(q)}^{F_{A}(q)}1\,dp&\text{otherwise}\\\ \end{cases}$
(26) $\displaystyle=|F_{A}(q)-F_{B}(q)|.$ (27)

Similarly, we have $\displaystyle\int_{0}^{1}s(q,p)\,dq=|Q_{A}(p)-Q_{B}(p)|$.
Finally, by substituting these results in 25, we prove 18.

∎

##  Appendix B Detailed Results

This section presents additional experimental results.

###  B.1 Comparison between Recalibration, Conformal Prediction and
Regularization Approaches per Base Model

First, we present the results of our experiments comparing recalibration,
conformal prediction, and regularization approaches. Our objective is to
determine which metrics are improved by these methods compared to a vanilla
model. We divide our comparisons based on the three base models considered:
MIX-NLL (Figure 6), MIX-CRPS (Figure 7) and SQR-CRPS (Figure 8).

Since NLL, CRPS, and standard deviation cannot be directly compared across
different datasets, we utilize Cohen’s d as an effect size measure, with the
baseline being a vanilla model of the same base model. For instance, the
baseline for MIX-CRPS + Rec-EMP is MIX-CRPS. Additionally, we provide critical
difference diagrams to assess the significance of differences.

Overall, recalibration and conformal prediction demonstrate significantly
improved PCE compared to the baseline, although there is a trade-off with
other metrics. For both MIX-NLL and MIX-CRPS, Rec-EMP yields infinite NLL,
Rec-LIN substantially increases NLL, while Rec-KDE has a lesser impact on NLL.
However, Rec-KDE results in a significant degradation of CRPS compared to
other recalibration methods when the base model is MIX-CRPS. In the case of
quantile predictions, CQR significantly improves PCE.

While regularization methods generally lead to improved PCE, they are still
outperformed by recalibration and conformal prediction in this regard.
However, we observe that with the MIX-NLL base model, regularization methods
(PCE-KDE, PCE-Sort and QR) have minimal impact on CRPS, NLL, and STD compared
to recalibration methods. With the MIX-CRPS base model, the difference in CRPS
between recalibration and regularization is less pronounced. Nevertheless, it
is evident that regularization methods PCE-KDE and PCE-Sort, which rely on
PCE, result in less sharp predictions compared to recalibration methods, which
produce sharper predictions.

Regarding quantile predictions, the case is reversed: conformal prediction
(SQR-CRPS + CQR) yields less sharp predictions, while regularization with SQR-
CRPS + Trunc leads to sharper predictions.

((a)) Boxplots of Cohen’s d of different metrics on all datasets, with respect
to MIX-NLL.  ((b)) PCE  ((c)) CRPS  ((d)) NLL  ((e)) STD

Figure 6: Comparison of different metrics where the base model is MIX-NLL.

((a)) Boxplots of Cohen’s d of different metrics on all datasets, with respect
to MIX-CRPS.  ((b)) PCE  ((c)) CRPS  ((d)) NLL  ((e)) STD

Figure 7: Comparison of different metrics where the base model is MIX-CRPS.

((a)) Boxplots of Cohen’s d of different metrics on all datasets, with respect
to SQR-CRPS.  ((b)) PCE

((c)) CRPS  ((d)) STD

Figure 8: Comparison of different metrics where the base model is SQR-CRPS.

###  B.2 Combining Regularization and Post-hoc Methods

In this paper, we have established that post-hoc methods are generally more
favorable than regularization methods when the primary objective is to enhance
probabilistic calibration. Since regularization methods operate during
training and do not alter the form of predictions (e.g., Gaussian mixture
predictions), they can be easily combined with post-hoc methods. In this
section, we address the question: "Which metrics do regularization methods
improve when combined with a post-hoc method compared to the same model
without regularization?"

To ensure clarity, we focus our presentation on a selection of paired
regularization and post-hoc methods. Figure 9 illustrates the impact of
regularization on various metrics for these pairs. In Figure 9(a), the
baseline corresponds to the same post-hoc method without regularization,
enabling a direct measurement of the effect of adding regularization to a
post-hoc method. It is important to note that the boxplots in this figure
cannot be directly compared due to the different baselines.

The critical difference diagrams provide a comparison of all methods, with and
without regularization. Overall, when combined with post-hoc methods,
regularization has a negative impact: no regularization method significantly
improves probabilistic calibration, and they tend to negatively affect CRPS,
NLL, and STD metrics.

((a)) Boxplots of Cohen’s d of different metrics on all datasets, with respect
to the same model except that regularization is not applied.  ((b)) PCE  ((c))
CRPS  ((d)) NLL  ((e)) STD

Figure 9: Comparison of different metrics showing the effect of regularization
when combined with a post-hoc method, compared to the same model without
regularization.

###  B.3 Post-hoc Calibration based on the Training Dataset

In this paper, the calibration map or conformity scores have been computed on
a separate calibration dataset, following common practice in the literature.
However, holding out data for post-hoc calibration reduces the quantity of
training data. For the sake of clarity, we focus our analysis on the MIX-NLL
and SQR-CRPS base losses

In this section, we compare post-hoc calibration based on the training dataset
to post-hoc calibration based on the calibration dataset. We aim to answer the
question: "Can it be beneficial to use post-hoc calibration based on the
training dataset, and should it be preferred over regularization methods when
there is no calibration dataset available?" One advantage of regularization
methods and post-hoc calibration methods based on the training dataset is that
the base model can be trained on more data (80% in our experiments, compared
to 65% when holding out the calibration dataset).

Figure 10 presents a comparison of different methods, with post-hoc methods
trained on the calibration dataset indicated by (calib) and those trained on
the training dataset indicated by (train). We observe that post-hoc methods
based on the calibration dataset tend to significantly outperform their
counterparts based on the training dataset in terms of probabilistic
calibration. Specifically, MIX-NLL + Rec-LIN and MIX-NLL + Rec-KDE achieve
significantly better calibration when the calibration map is learned on the
calibration dataset. Similarly, SQR-CRPS + CQR tends to improve calibration
when conformal prediction is based on the calibration dataset. It is worth
noting that even without a calibration dataset, post-hoc methods tend to be
better calibrated than regularization methods.

Finally, we observe that post-hoc methods based on the training dataset tend
to achieve better CRPS and NLL scores, although not significantly.
Additionally, they are also significantly sharper. This may be attributed to
the larger training dataset available to the base model when there is no held-
out dataset.

((a)) Boxplots of Cohen’s d of different metrics on all datasets, with respect
to MIX-NLL.  ((b)) PCE  ((c)) CRPS  ((d)) NLL  ((e)) STD

Figure 10: Comparison of different metrics.

###  B.4 Calibration of Vanilla Models

Figure 12 and Figure 13 provide additional results from our empirical study in
Section 4, specifically focusing on the PCE obtained with MIX-CRPS and SQR-
CRPS. The datasets are ordered in the same manner as shown in Figure 2 for
comparison. We observe that SQR-CRPS is less calibrated compared to MIX-NLL.

Figure 11: PCE obtained on different datasets, with examples of reliability
diagrams. The height of each bar is the mean PCE of 5 runs with different
dataset splits while the error bar represents the standard error of the mean.
For 5 datasets, the PIT reliability diagrams of 5 runs are displayed in the
bottom row.  Figure 12: PCE of SQR-CRPS, on all datasets. Figure 13: PCE of
MIX-CRPS, on all datasets.

###  B.5 Distribution of the Test Statistic

Figure 14 shows the distribution of the test statistic, as described in
Section 4. We observe that, in a lot of cases, the average PCE of the compared
models is larger than all the $10^{4}$ samples of the average PCE from a
probabilistically calibrated model. Among the different calibration methods,
post-hoc calibration with MIX-NLL + Rec-EMP achieves the highest level of
calibration performance in the majority of cases.

Figure 14: Distribution of the test statistic on all datasets for different
models.

###  B.6 Reliability Diagrams

Figure 15 and Figure 16 compare reliability diagrams obtained on models with
and without post-hoc calibration, respectively. With only a few exceptions,
the post-hoc calibrated models exhibit a visual proximity to the diagonal
line.

Figure 15: Reliability diagrams on all datasets for different models. Figure
16: Reliability diagrams on all datasets for different models with post-hoc
calibration.

##  Appendix C Hyperparameters

In our experiments, we adopt a specific architecture consisting of 3 hidden
layers with 100 units per layer, ReLU non-linearities, and a dropout rate of
0.2 on the last hidden layer. Early stopping with a patience of 30 is applied
to select the epoch with the lowest base loss on the validation dataset.

In this section, we delve into the performance of different model parameters,
including the number of components in Gaussian mixture predictions, the number
of quantiles in quantile predictions, and the number of hidden layers in the
underlying models.

Figure 17 compares models that predict mixtures with varying numbers of
components compared to the reference of 3 components. Notably, when there is
only 1 component (yielding a single Gaussian prediction), the model’s
performance significantly deteriorates in terms of CRPS, NLL, and sharpness.
However, as the number of components increases beyond 3, the differences
become less pronounced.

Figure 18 compares models with different numbers of quantiles compared to a
reference of 64 quantiles. The results reveal a consistent pattern: predicting
more quantiles consistently enhances performance in terms of probabilistic
calibration, CRPS, and sharpness.

Figure 19 compares models with different numbers of layers relative to a
3-layer model. It highlights that models with 2, 3, or 5 layers tend to yield
superior performance in terms of CRPS and NLL.

Figure 17: Comparison of models whose predictions are Gaussian mixtures with
different numbers of components. All models are trained with NLL loss, without
regularization or post-hoc method. The box plots show Cohen’s d of different
metrics on all datasets. Cohen’s d is computed with respect to a model whose
predictions are Gaussian mixtures with 3 components. Figure 18: Comparison of
models whose predictions are different numbers of quantiles. All models are
trained with CRPS loss, without regularization or post-hoc method. The box
plots show Cohen’s d of different metrics on all datasets. Cohen’s d is
computed with respect to a model whose predictions are 64 quantiles. Figure
19: Comparison of models with different number of layers. All models predict
Gaussian mixtures and are trained with NLL loss, without regularization or
post-hoc method. The box plots show Cohen’s d of different metrics on all
datasets. Cohen’s d is computed with respect to a model with 3 hidden layers.

##  Appendix D Tabular Regression Datasets

Table 1 presents the datasets considered in our experiments. To ensure
consistency, when datasets are available from multiple sources, we select one
specific source per dataset, as indicated in Figure 1. Our selection
prioritizes the suites 297, 299, and 269 of OpenML, followed by UCI datasets.

In the OpenML suite 297, we discovered that the datasets houses and california
are identical, and thus, we only included the california dataset in our
analysis. Moreover, the UCI archive for the dataset wine_quality contains two
separate datasets for red and white wine. As there was no indication regarding
the specific dataset(s) used in previous studies, we followed the approach of
Grinsztajn2022-nu (21) and solely considered the dataset related to white
wine. In Figure 1, other studies may have employed the alternative dataset or
a combination of both datasets.

Table 1: Datasets |  |  |  Nb of training instances |  Nb of features
---|---|---|---|---
Group | Dataset | Abbrev. |  |
UCI | CPU | CP1 | 135 | 7
Yacht | YAC | 200 | 6
MPG | MPG | 254 | 7
Energy | ENE | 499 | 9
Crime | CRI | 531 | 104
Fish | FIS | 590 | 6
Concrete | CON | 669 | 8
Airfoil | AI1 | 976 | 5
Kin8nm | KIN | 5324 | 8
Power | POW | 6219 | 4
Naval | NAV | 7757 | 17
Protein | PRO | 29724 | 9
OpenML 297 | wine_quality | WIN | 4223 | 11
isolet | ISO | 5068 | 613
cpu_act | CP2 | 5324 | 21
sulfur | SUL | 6552 | 6
Brazilian_houses | BRA | 6949 | 8
Ailerons | AIL | 8937 | 33
MiamiHousing2016 | MIA | 9055 | 13
pol | POL | 9750 | 26
elevators | ELE | 10789 | 16
Bike_Sharing_Demand | BIK | 11296 | 6
fifa | FIF | 11740 | 5
california | CAL | 13416 | 8
superconduct | SUP | 13820 | 79
house_sales | HO3 | 14048 | 15
house_16H | HO1 | 14809 | 16
diamonds | DIA | 35061 | 6
medical_charges | MED | 50000 | 3
year | YEA | 50000 | 90
nyc-taxi-green-dec-2016 | NYC | 50000 | 9
OpenML 299 | analcatdata_supreme | ANA | 2633 | 12
Mercedes_Benz _Greener_Manufacturing | MER | 2735 | 735
visualizing_soil | VIS | 5616 | 5
yprop_4_1 | YPR | 5775 | 82
OnlineNewsPopularity | ONL | 25768 | 73
black_friday | BLA | 50000 | 23
SGEMM_GPU _kernel_performance | SGE | 50000 | 15
particulate-matter -ukair-2017 | PAR | 50000 | 26
OpenML 269 | tecator | TEC | 156 | 124
boston | BOS | 328 | 22
MIP-2016-regression | MIP | 708 | 111
socmob | SOC | 751 | 39
Moneyball | MON | 800 | 18
house_prices_nominal | HO2 | 711 | 234
us_crime | US_ | 1295 | 101
quake | QUA | 1415 | 3
space_ga | SPA | 2019 | 6
abalone | ABA | 2715 | 10
SAT11-HAND- runtime-regression | SAT | 2886 | 118
Santander_transaction _value | SAN | 2898 | 3611
colleges | COL | 4351 | 34
topo_2_1 | TOP | 5775 | 252
Allstate_Claims_Severity | ALL | 50000 | 477
Yolanda | YOL | 50000 | 100
Buzzinsocialmedia_Twitter | BUZ | 50000 | 70
Airlines_DepDelay_10M | AI2 | 50000 | 5