Likelihood-ratio-based confidence intervals for neural networks
This paper introduces a first implementation of a novel
likelihood-ratio-based approach for constructing confidence intervals for
neural networks. Our method, called DeepLR, offers several qualitative
advantages: most notably, the ability to construct asymmetric intervals that
expand in regions with a limited amount of data, and the inherent incorporation
of factors such as the amount of training time, network architecture, and
regularization techniques. While acknowledging that the current implementation
of the method is prohibitively expensive for many deep-learning applications,
the high cost may already be justified in specific fields like medical
predictions or astrophysics, where a reliable uncertainty estimate for a single
prediction is essential. This work highlights the significant potential of a
likelihood-ratio-based uncertainty estimate and establishes a promising avenue
for future research.
Laurens Sluijterman, Eric Cator, Tom Heskes
# Likelihood-ratio-based confidence intervals

###### Abstract

This paper introduces a first implementation of a novel likelihood-ratio-based approach for constructing confidence intervals for neural networks. Our method, called DeepLR, offers several qualitative advantages: most notably, the ability to construct asymmetric intervals that expand in regions with a limited amount of data, and the inherent incorporation of factors such as the amount of training time, network architecture, and regularization techniques. While acknowledging that the current implementation of the method is prohibitively expensive for many deep-learning applications, the high cost may already be justified in specific fields like medical predictions or astrophysics, where a reliable uncertainty estimate for a single prediction is essential. This work highlights the significant potential of a likelihood-ratio-based uncertainty estimate and establishes a promising avenue for future research.

**Keywords:** Likelihood ratio, Uncertainty estimation, Neural network,

Regression, Classification

## 1 Introduction

Over the past two decades, neural networks have seen an enormous rise in popularity and are currently being used in almost every area of science and
industry. In light of this widespread usage, it has become increasingly clear that trustworthy uncertainty estimates are essential (Gal, 2016).

Many uncertainty estimation methods have been developed using Bayesian techniques (Neal et al., 2011; MacKay, 1992; Gal, 2016), ensembling techniques (Heskes, 1997; Lakshminarayanan et al., 2017), or applications of frequentist techniques such as the delta method (Kallus and McInerney, 2022).

Many of the resulting confidence intervals (for the frequentist methods) and credible regions (for the Bayesian methods) have two common issues. Firstly, most methods result in symmetric intervals around the prediction which can be overly restrictive and can lead to very low coverage in biased regions (Sluijterman et al., 2022). Secondly, most methods rely heavily on asymptotic theorems (such as the central limit theorem or the Bernstein-von-Mises theorem) and can therefore only be trusted in the asymptotic regime where we have many more data points than model parameters, the exact opposite scenario of where we typically find ourselves within machine learning.

#### Contribution

In this paper, we demonstrate how the likelihood-ratio test can be leveraged to combat the two previously mentioned issues. We provide a first implementation of a likelihood-ratio-based approach, called DeepLR, that has the ability to produce asymmetric intervals that are more appropriately justified in the scenarios where we have more parameters than data points. Furthermore, these intervals exhibit desirable behavior in regions far removed from the data, as evidenced in Figure 1.

Figure 1: A comparison of the confidence intervals of our likelihood-ratio approach (DeepLR), an ensembling approach, and MC-Dropout on the two-moon data set. The colorbar represents the width of 95% confidence intervals, where yellow indicates greater uncertainty. The orange and the blue circles indicate the location of the data points of the two different classes. Crucially, DeepLR presents high levels of uncertainty in regions far away from the data, generating confidence intervals of [0.00, 1.00], unlike the other methods that display extreme certainty in those regions.
#### Organisation

This paper is structured into four sections, with this introduction being the first. Section 2 explains our method in detail and also contains the related work section which is simultaneously used to highlight the advantages and disadvantages of our method. Section 3 presents experimental results that illustrate the desirable properties of a likelihood-ratio-based approach. Finally, Section 4 summarizes and discusses the results and outlines possible directions for future work.

## 2 DeepLR: Deep Likelihood-Ratio-based confidence intervals

In this section, we present our method, named DeepLR, for constructing confidence intervals for neural networks using the likelihood-ratio test. We first formalize the problem that we are considering in Subsection 2.1. Subsection 2.2 explains the general idea behind constructing a confidence interval via the likelihood-ratio test. Subsequently, in Subsection 2.3, we outline the high level idea for translating this general procedure to neural networks. The details regarding the distribution and the calculation of the test statistic are provided in Subsections 2.4 and 2.5. Finally, in Subsection 2.6, we compare our method to related work while simultaneously highlighting its strengths and limitations.

### Problem formulation

We consider a data set \(\mathcal{D}=\{(X_{1},Y_{1}),\ldots,(X_{n},Y_{n})\}\), consisting of \(n\) independent observations of the random variable pair \((X,Y)\). We consider networks that provide an estimate for the conditional density of \(Y\mid X\). This is achieved by assuming a distribution and having the network output the parameter(s) of that distribution.

Three well-known types of networks that fall in this class are: (1) A regression setting where the network outputs a mean estimate and is trained using a mean-squared error loss. This is equivalent to assuming a normal distribution with homoscedastic noise. (2) Alternatively, the network could output both a mean and a variance estimate and be optimized by minimizing the negative loglikelihood assuming a normal distribution. (3) In a classification setting, the network could output logits that are transformed to class probabilities while assuming a categorical distribution.

The network is parametrized by \(\theta\in\mathbb{R}^{p}\), where \(p\) is typically much larger than \(n\). With \(\Theta\), we denote the set containing all the \(\theta\) that are reachable for a network with a specific training process. This includes choices such as training time, batch size, optimizer, and regularization techniques. With \(p_{\theta}\), we denote the predicted conditional density. Additionally, we assume that the true conditional density is given by \(p_{\theta_{0}}\) for some \(\theta_{0}\) in \(\Theta\). In other words, we assume that our model is well specified.
The objective of our method is to construct a confidence interval for one of the output nodes of the network for a specific input of interest \(X_{0}\). We denote this output of interest with \(f_{\theta_{0}}(X_{0})\) for the remainder of the paper. In the context of a regression setting, this output of interest is the true regression function value at \(X_{0}\) and in a classification setting it is the true class probability for input \(X_{0}\).

We define a \((1-\alpha)\cdot 100\%\) confidence interval for \(f_{\theta_{0}}(X_{0})\) as an interval, \(\mathrm{CI}(f_{\theta_{0}}(X_{0}))\), which is random since it depends on the random realization of the data, such that the probability (taken with respect to the random data set) that \(\mathrm{CI}(f_{\theta_{0}}(X_{0}))\) contains the true value \(f_{\theta_{0}}(X_{0})\) is \((1-\alpha)\cdot 100\%\).

### Confidence interval based on the likelihood ratio

We explain the general idea behind constructing a confidence interval with the likelihood-ratio by working through a well-known example. We consider \(n\) observations \(Y_{i}\) that are assumed to be normally distributed with unknown mean \(\mu\) and unknown variance \(\sigma^{2}\). Our goal is to create a confidence interval for \(\mu\) by using the likelihood-ratio test.

The duality between a confidence interval and hypothesis testing states that we can create a \((1-\alpha)\cdot 100\%\) confidence interval for \(\mu\) by including all the values \(c\) for which the hypothesis \(\mu=c\) cannot be rejected at a \((1-\alpha)\cdot 100\%\) confidence level. We must therefore test for what values \(c\) we can accept the hypothesis \(\mu=c\).

The general approach to test a hypothesis is to create a test statistic of which we know the distribution under the null hypothesis and to reject this hypothesis if the probability of finding the observed test statistic or an extremer value is smaller than \(\alpha\).

As our test statistic, we take two times the log of the likelihood ratio:

\[T(c):=2\left(\sup_{\Theta}\left(\sum_{i=1}^{n}\log(L(Y_{i};\theta))\right)- \sup_{\Theta_{0}}\left(\sum_{i=1}^{n}\log(L(Y_{i};\theta))\right)\right),\]

where \(L\) denotes the likelihood function \(\theta\mapsto L(Y_{i};\theta)\), \(\Theta\) is the full parameter space and \(\Theta_{0}\) the restricted parameter space. In our example, we have

\[\Theta=\{(\mu,\sigma^{2})\mid\mu\in\mathbb{R},\sigma^{2}\in\mathbb{R}_{>0}\},\]

and

\[\Theta_{0}=\{(c,\sigma^{2})\mid\sigma^{2}\in\mathbb{R}_{>0}\}.\]

Wilks (1938) proved that \(T(c)\) weakly converges to a \(\chi^{2}(1)\) distribution under the null hypothesis that \(\mu=c\). We therefore reject if \(T(c)>\chi^{2}_{1-\alpha}(1)\) and our confidence interval for \(\mu\) becomes the set \(\{c\mid T(c)\leq\chi^{2}_{1-\alpha}(1)\}\), where \(\chi^{2}_{1-\alpha}(1)\) is the \((1-\alpha)\)-quantile of a \(\chi^{2}(1)\) distribution.

In our example, this results in the well-known interval
\[\bar{Y}\pm z_{1-\alpha/2}\sqrt{\frac{1}{n}\frac{1}{n-1}\sum_{i=1}^{n}(Y_{i}-\bar{Y} )^{2}},\]

where \(z_{1-\alpha/2}\) is the \((1-\alpha/2)\)-quantile of a standard-normal distribution.

### High-level idea of DeepLR

Our goal is to apply the likelihood-ratio testing procedure outlined in the previous subsection to construct a confidence interval for \(f_{\theta_{0}}(X_{0})\), the value of one of the output nodes given input \(X_{0}\). We create this confidence interval by including all the values \(c\) for which the hypothesis \(f_{\theta_{0}}(X_{0})=c\) cannot be rejected. The testing of the hypothesis is done with the likelihood-ratio test. Specifically, we use two times the log likelihood ratio as our test statistic:

\[T(c) :=2\bigg{(}\sup_{\Theta}\big{(}\sum_{i=1}^{n}\log(L(X_{i},Y_{i}; \theta))\big{)}-\sup_{\Theta_{0}(c)}\big{(}\sum_{i=1}^{n}\log(L(X_{i},Y_{i}; \theta))\big{)}\bigg{)}\] \[=2\bigg{(}\sup_{\Theta}\big{(}\sum_{i=1}^{n}\log(p_{\theta}(Y_{i} \mid X_{i}))\big{)}-\sup_{\Theta_{0}(c)}\big{(}\sum_{i=1}^{n}\log(p_{\theta}(Y _{i}\mid X_{i}))\big{)}\bigg{)}, \tag{1}\]

and we construct a confidence interval for \(f_{\theta_{0}}(X_{0})\) by including all values \(c\) for which the test statistic is not larger than \(\chi^{2}_{1-\alpha}(1)\):

\[\text{CI}(f_{\theta_{0}}(X_{0}))=\{c\mid T(c)\leq\chi^{2}_{1-\alpha}(1)\}. \tag{2}\]

Here, we consider \(\Theta\subset\mathbb{R}^{p}\) to be the set containing all reachable parameters, and \(\Theta_{0}(c)=\{\theta\in\Theta\mid f_{\theta_{0}}(X_{0})=c\}\). The set \(\Theta\) is explicitly not equal to all parameter combinations. Due to explicit (e.g., early stopping) and implicit regularization (e.g., lazy training (Chizat et al., 2019)) not all parameter combinations can be reached. The set \(\Theta\) should therefore be seen as the set containing the parameters of all neural networks that can be found given the optimizer, training time, regularization techniques, and network architecture.

Intuitively, this approach answers the question: _What values could the network have reached at location \(X_{0}\) while still explaining the data well?_ This is a sensible question for a highly flexible and typically overparameterized machine-learning approach. After training, the model ends up with a certain prediction at location \(X_{0}\). However, since the model is typically very complex, it is likely that the model could just as well have made other predictions at that location while still explaining the data well. Therefore, all those other function values should also be considered as possibilities. Inherently, all modeling choices are taken into account by asking this question. A more flexible model, for instance, is likely able to reach more values without affecting the likelihood of the training data, leading to a larger confidence interval.

The construction of the confidence interval in Equation (2) assumes that the test statistic, \(T(c)\), has a \(\chi^{2}(1)\) distribution. We discuss this assumption in
the following subsection. The subsection thereafter describes how to calculate the test statistic.

### Distribution of the test statistic

In the classical setting, Wilks (1938) proved that the likelihood-ratio test statistic asymptotically has a \(\chi^{2}(1)\) distribution when the submodel has one degree of freedom less than the full model. We are, however, not in this classical regime. We have many more parameters than data points and therefore need a similar result for this setting.

It has been shown that the likelihood-ratio test statistic converges to a \(\chi^{2}\) distribution for a wide range of settings, which is referred to as the Wilks-phenomenon by later authors (Fan et al., 2001; Boucheron and Massart, 2011). For a semi-parametric model, which more closely resembles our situation, it has been shown that the test statistic also converges in distribution to a \(\chi^{2}\) distribution under appropriate regularity conditions (Murphy and van der Vaart, 1997).

We prove a similar result for our setting in the appendix. The theorem states that, under suitable assumptions, our test statistic has a \(\chi^{2}(1)\) distribution. Intuitively, this results from the fact that we added a single constraint, namely that \(f_{\theta}(X_{0})=c\). We emphasize that even if the test statistic does not exactly follow a \(\chi^{2}(1)\) distribution, the qualitative characteristics of the confidence intervals will remain evident, albeit with inaccurate coverage levels.

### Calculating the test statistic

Calculating the two terms in equation (1) presents certain challenges. The first term, \(\sup_{\Theta}\big{(}\sum_{i=1}^{n}\log(p_{\theta}(Y_{i}\mid X_{i}))\big{)}\), is relatively straightforward. We train a network that maximizes the likelihood, which gives the conditional densities \(p_{\hat{\theta}}(Y_{i}\mid X_{i})\).

The second term is substantially more complex. Ideally, we would optimize over the set \(\Theta_{0}(c)=\{\theta\in\Theta\mid f_{\theta}(X_{0})=c\}\). This is problematic for two reasons. Firstly, it is unclear how we can add this constraint to the network and secondly, this would necessitate training our network for every distinct value \(c\) that we wish to test. Even when employing an efficient bisection algorithm, this could easily result in needing to train upwards of 10 additional networks.

We address this problem as follows. We first create a network that is perturbed in the direction of a relatively large value (\(c_{\max}\)) at \(X_{0}\) and a network that is perturbed in the direction a relatively small value (\(c_{\min}\)) at \(X_{0}\) while maximizing the likelihood of the data. We denote the resulting network parameters with \(\hat{\theta}_{\pm}\):

\[\hat{\theta}_{+}=\underset{\begin{subarray}{c}\theta\in\Theta,\\ f_{\theta}(X_{0})\approx c_{\min}\end{subarray}}{\arg\max}L(\mathcal{D};\theta),\quad\text{and}\quad\hat{\theta}_{-}=\underset{\begin{subarray}{c}\theta \in\Theta,\\ f_{\theta}(X_{0})\approx c_{\min}\end{subarray}}{\arg\max}L(\mathcal{D};\theta).\]
Subsequently, the network that maximizes the likelihood under the constraint \(f_{\theta}(X_{0})=c\) is approximated using a linear combination. Specifically, suppose we want the network that maximizes the likelihood of the training data while passing through \(c\) at \(X_{0}\). In the case that \(c>f_{\hat{\theta}}(X_{0})\), we approximate this network by taking a linear combination of the outputs such that

\[(1-\lambda)f_{\hat{\theta}}(X_{0})+\lambda f_{\hat{\theta}_{+}}(X_{0})=c,\]

and we define \(p_{c}\) as the density that we get by using the same linear combinations for the distributional parameters that are predicted by the networks parametrized by \(\hat{\theta}\) and \(\hat{\theta}_{+}\). We then approximate the second term in equation (1) as follows:

\[\sup_{\Theta_{0}(c)}\big{(}\sum_{i=1}^{n}\log(p_{\theta}(Y_{i}\mid X_{i})) \big{)}\approx\sum_{i=1}^{n}\log(p_{c}(Y_{i}\mid X_{i})). \tag{3}\]

This procedure is visualized in Figure 2 for a regression setting. Details on the second step, finding the perturbed networks, are provided below both for a regression and binary-classification setting.

RegressionSuppose we completed step 1 and we have a network that maximizes the data's likelihood. For step 2, our objective is to adjust this network such that it goes through a relatively large value at \(X_{0}\) while continuing to maximize the data's likelihood. Moreover, we aim to achieve this in as stable a manner as possible, given that minor differences can significantly influence the test statistic, particularly for large data sets.

We accomplish this by copying the original network and training it on the objective to maintain the original predictions - as those maximized the likelihood - while predicting \(f_{\hat{\theta}}(X_{0})+1\) for input \(X_{0}\). For the perturbation in the negative direction we use \(-1\). These values \(\pm 1\) are chosen assuming that the targets are normalized prior to training.

We obtain this training objective by using modified training set, \(\tilde{\mathcal{D}}\), that is constructed by replacing the targets \(Y_{i}\) in the original training set with the predictions \(f_{\hat{\theta}}(X_{i})\) of the original network and adding the data point \((X_{0},f_{\hat{\theta}}(X_{0})+1)\).

The resulting problem is very imbalanced. We want the network to change the prediction at location \(X_{0}\), which is only present in the data once. This makes the training very unstable since, especially when training for a small number of epochs, it is very influential which specific batch contains the new point. To remedy this, we use a combination of upsampling and downweighting of the new data point \((X_{0},f_{\hat{\theta}}(X_{0})+1)\). Merely upsampling the new data point - i.e., adding it multiple times - is undesirable as this can introduce significant biases (Van den Goorbergh et al., 2022). Hence, we add many copies of \((X_{0},f_{\hat{\theta}}(X_{0})+1)\) but we reweigh the loss contributions of these added data
points such that they have a total contribution to the loss that is equivalent to that of a single data point.

We propose to add \(2n/\)batch size extra data points such that each batch is expected to have 2 new data points. The same training procedure is used as for the original network. We found slightly larger or smaller number of added

Figure 2: Illustration of the steps of our method for the positive direction in a regression setting. Steps 2, 3, and 4 are also carried out in the negative direction.
data points to perform very similarly. This setting worked for a wide variety of data sets and architectures.

Binary classificationConsider a data set where the targets are either \(1\) or \(0\). Our network outputs logits, denoted with \(f_{\hat{\theta}}(X)\), that are transformed to probabilities via a sigmoid function. The procedure is nearly identical for this binary classification setting: We create a positively perturbed network, parametrized by \(\hat{\theta}_{+}\), and a negatively perturbed network, parametrized by \(\hat{\theta}_{-}\).

The only difference is in the construction of the augmented data sets. We again replace the targets \(Y_{i}\) by the predictions of the original network, \(f_{\hat{\theta}}(X_{i})\), but now we add multiple copies of the data point \((X_{0},1)\) for the positive direction, and multiple copies of the data point \((X_{0},0)\) for the negative direction.

The entire method is summarized in Algorithm 1. In summary, we want to test what values the network can reach while still explaining the data well. We do this by perturbing the network in a positive direction and a negative direction and subsequently testing which linear combinations would still explain the data reasonably well, i.e., linear combinations with a test statistic smaller than \(\chi^{2}_{1-\alpha}(1)\).

### Related work

In this subsection, we place our work in context of existing work while simultaneously highlighting the strengths and limitations of our approach. Our aim is not to give a complete overview of all uncertainty quantification methods, for which we refer to the various reviews and surveys on the subject (Abdar et al., 2021; Gawlikowski et al., 2022; He and Jiang, 2023). Instead, we discuss several broad groups in which most methods can be categorized: Ensembling methods, Bayesian methods, frequentist methods, and distance-aware methods.

Ensembling methodstrain multiple models and use the variance of the predictions as an estimate for model uncertainty (Heskes, 1997; Lakshminarayanan et al., 2017; Zhang et al., 2017; Wenzel et al., 2020; Jain et al., 2020; Dwaracherla et al., 2022). While being extremely easy to implement, they can be computationally expensive due to the need to train multiple networks. Moreover, the resulting confidence intervals can behave poorly in regions with a limited amount of data where the predictor is likely biased. Additionally, ensemble members may interpolate in a very similar manner, potentially leading to unreasonably narrow confidence intervals.

Bayesian approachesplace a prior distribution on the model parameters and aim to simulate from the resulting posterior distribution given the observed data (MacKay, 1992; Neal, 2012; Hernandez-Lobato and Adams, 2015). Since this posterior is generally intractable, it is often approximated, with MC-Dropout being a notable example (Gal and Ghahramani, 2016; Gal et al., 2017).
```
0:\(\mathcal{D}=\{(X_{1},Y_{1}),\ldots,(X_{n},Y_{n})\}\), \(\hat{\theta},X_{0}\), \(\alpha\);
1:\(n_{\text{extra}}=2n/\text{(batch size)}\);
2: For binary classification \(c_{\text{max}}=1\) and \(c_{\text{min}}=0\), and for regression \(c_{\text{max}}=f_{\hat{\theta}}(X_{0})+\delta\) and \(c_{\text{min}}=f_{\hat{\theta}}(X_{0})-\delta\);
3:\(\tilde{\mathcal{D}}_{+}=\{(X_{1},f_{\hat{\theta}}(X_{1})),\ldots,(X_{n},f_{ \hat{\theta}}(X_{n})),\overbrace{(X_{0},c_{\text{max}}),\ldots,(X_{0},c_{ \text{max}})}^{n_{\text{extra}}\text{times}}\}\);
4:\(\tilde{\mathcal{D}}_{-}=\{(X_{1},f_{\hat{\theta}}(X_{1})),\ldots,(X_{n},f_{ \hat{\theta}}(X_{n})),\overbrace{(X_{0},c_{\text{min}}),\ldots,(X_{0},c_{ \text{min}})}^{n_{\text{extra}}\text{times}}\}\);
5: Make two copies of network parametrized by \(\hat{\theta}\) and train them on \(\tilde{\mathcal{D}}_{+}\) and \(\tilde{\mathcal{D}}_{-}\) using the original training procedure. During the training, the loss contribution of the added data points is divided by \(n_{\text{extra}}\). Denote the parameters of the resulting networks with \(\hat{\theta}_{+}\) and \(\hat{\theta}_{-}\);
6:for\(c\in\mathbb{R}\)do\(\triangleright\) Since testing all possible \(c\) is impossible, we propose to use some variation of a bisection search algorithm.
7:if\(c>f_{\hat{\theta}}(X_{0})\)then
8: Pick \(\lambda\) such that \((1-\lambda)f_{\hat{\theta}}(X_{0})+\lambda f_{\hat{\theta}_{+}}(X_{0})=c\) with corresponding density \(p_{c}\); \(\triangleright\) The density \(p_{c}\) is obtained by taking the same linear combination of the predicted distribution parameters outputted by the networks parametrized by \(\hat{\theta}\) and \(\hat{\theta}_{+}\).
9:endif
10:if\(c<f_{\hat{\theta}}(X_{0})\)then
11: Pick \(\lambda\) such that \((1-\lambda)f_{\hat{\theta}}(X_{0})+\lambda f_{\hat{\theta}_{-}}(X_{0})=c\) with corresponding density \(p_{c}\);
12:endif
13:if\(2\bigg{(}\sum_{i=1}^{n}\log(p_{\hat{\theta}}(Y_{i}\mid X_{i}))-\sum_{i=1}^{ n}\log(p_{c}(Y_{i}\mid X_{i}))\bigg{)}\leq\chi_{1-\alpha}^{2}(1)\)then
14: Include \(c\) in \(\text{CI}(f_{\theta_{0}}(X_{0}))\);
15:endif
16:endfor
17:return\(\text{CI}(f_{\theta_{0}}(X_{0}))\)
```

**Algorithm 1** Pseudocode for the construction of \(\text{CI}(f_{\theta_{0}}(X_{0}))\) using DeepLR

A downside is that these methods can be challenging to train, potentially resulting in a lower accuracy. Our proposed approach does not change the optimization procedure and therefore has no accuracy loss. Another downside is that while asymptotically Bayesian credible sets become confidence sets by the Bernstein-von Mises theorem (see Van der Vaart (2000, Chapter 10))? However, this theorem generally does not apply for a neural network where the dimension of the parameter space typically exceeds the number of data points. Moreover, the prior distribution is often chosen out of computational convenience instead of being motivated by domain knowledge.

_Distance-aware methods_ have a more pragmatic nature. They use the dissimilarity of a new input compared to the training data as a metric for model uncertainty (Lee et al., 2018; Van Amersfoort et al., 2020; Ren et al., 2021).
As we will see in the next section, our method also exhibits this distance-aware property, albeit for a different reason: The further away from the training data, the easier it becomes for the network to change the predictions without negatively affecting the likelihood of the training data.

_Frequentist methods_ use classical parametric statistics to obtain model uncertainty estimates. The typical approach - more elaborately explained in textbooks on parametric statistics, e.g. Seber and Wild (2003) - involves obtaining (an estimate of) the variance of the model parameters using asymptotic theory and then converting this variance to the variance of the model predictions using the delta method. This approach has been used by various authors to create confidence intervals for neural networks (Kallus and McInerney, 2022; Nilsen et al., 2022; Deng et al., 2023; Khosravi et al., 2011).

Confidence intervals of this type are often referred to as Wald-type intervals. These intervals are necessarily symmetric. Various authors have noted that, for classical models, Wald-type intervals often behave worse than likelihood-ratio type intervals in the low-data regime (Hall and La Scala, 1990; Andersen et al., 2012; Murphy, 1995; Murphy and van der Vaart, 1997). Specifically, when the loglikelihood cannot be effectively approximated with a quadratic function, Wald-type intervals may behave very poorly (Pawitan, 2001, Chapter 2). The significant advantage of Wald-type intervals in the classical setting is the easier computation. However, while only a single model needs to be fitted, the necessary inversion of a high-dimensional \(p\times p\) matrix and the quadratic approximation of the likelihood strongly rely on being in the asymptotic regime.

Conversely, the construction of the DeepLR confidence interval does not rely on a quadratic approximation of the likelihood, which is in general only valid asymptotically. While we still utilize asymptotic theory to determine the distribution of our test statistic (see our proof of Theorem 1 in Appendix A), we do impose the extremely strong requirement that the second derivative of the loglikelihood must converge and be invertible. We only require this second derivative to behave nicely in a single direction, which is a much weaker requirement.

Another benefit of likelihood-ratio-based confidence intervals is that they are transformation invariant (Pawitan, 2000). In other words, a different parametrization of the distribution does not alter the resulting confidence intervals. This is not the case for most Wald-type intervals or Bayesian approaches.

The main limitation of DeepLR in its current form is the computational cost. We compare the computational costs of the various approaches in Table 1. Our method comes at no extra training cost but requires two additional networks to be trained for every confidence interval. Ensembling methods also need to train multiple networks, typically from five to ten, but these networks can be reused for different confidence intervals. Frequentist methods typically come at no additional training cost but require the inversion of a \(p\times p\) matrix
to construct the confidence interval. The cost of Bayesian methods varies drastically from method to method. The training process can be substantially more involved and the construction of a credible set requires multiple samples from the (approximate) posterior.

In summary, a likelihood-ratio-based method is distance aware, transformation invariant, has no accuracy loss, and is capable of creating asymmetric confidence intervals. However, it comes with the downside of being computationally expensive. Producing a confidence interval for a single input requires the training of two additional networks.

## 3 Experimental results

In this section, we present the results of various experiments that showcase the desirable properties of DeepLR, such as its distance-aware nature and capability to create asymmetric confidence intervals. The high computational cost of our method prohibits any large-scale experiments. Nevertheless, the following experiments demonstrate the effectiveness of a likelihood-ratio-based approach.

All code used in the following experiments can be found at [https://github.com/LaurensSluyterman/Likelihood_ratio_intervals](https://github.com/LaurensSluyterman/Likelihood_ratio_intervals)

### Toy examples

To start, we present two one-dimensional toy examples that effectively illustrate the behavior of our method. Specifically, these examples illustrate the capability of our method to produce asymmetric confidence intervals that expand in regions with a limited amount of data points, both during interpolation and extrapolation.

RegressionThe data set consists of 80 realisations of the random variable pair \((X,Y)\). Half of the \(x\)-values are sampled uniformly from the interval \([-1,-0.2]\), while the

\begin{table}
\begin{tabular}{|p{113.8pt}|p{113.8pt}|p{113.8pt}|} \hline Approach & Additional training cost & Additional inference cost \\ \hline DeepLR & No additional costs & Two additional networks per new input and multiple additional forward passes. A forward pass through each of the ensemble members. \\ Frequentist Bayesian & None & Inversion of a \(p\times p\) matrix. Varies from method to method. Typically, a large number of forward passes have to be made for every new input. \\ \hline \end{tabular}
\end{table}
Table 1: A comparison of the computational costs of different types of methods that create confidence intervals or credible regions.
remaining half are sampled uniformly from the interval \([0.2,1]\). The \(y\)-values are subsequently sampled using

\[Y\mid X=x\sim\mathcal{N}\left(2x^{2},0.1^{2}\right).\]

On this training set, we train a mean-variance estimation network (MVE) (Nix and Weigend, 1994). This particular type of network provides both mean and variance estimates and is trained by minimizing the negative loglikelihood under the assumption of a normal distribution. The network is trained for 400 epochs, using a default Adam optimizer (Kingma and Ba, 2014) and a batch size of 32. The network consists of 3 hidden layers for the mean estimation network with 40, 30, and 20 hidden units respectively and 2 hidden units for the variance using 5 and 2 hidden units respectively. All layers have elu activation functions (Clevert et al., 2015) with \(l_{2}\)-regularization applied in each dense layer with a constant value of 1e-4.

Figure 3 gives the 95% confidence intervals for mean predictions. Notably, in the biased region around 0, the intervals become highly asymmetric. In contrast, most other methods produce symmetric interval around the original network. In a region with a bias, this can easily lead to intervals with very poor coverage. This is exemplified by the biased symmetric intervals generated by an ensemble consisting of 10 ensemble members (we used the ensembling strategy employed by Lakshminarayanan et al. (2017)).

**Fig. 3**: This figure illustrates DeepLR for a regression problem. The blue dots indicate the locations of the training points, the dotted orange line represents the true function, and the solid green line represents the predicted regression function of the network. The shaded blue region gives the 95% CI of the regression function. DeepLR exhibits two desirable properties when compared to an ensemble approach. Firstly, the intervals expand in regions where data is sparse. Secondly, the intervals can be asymmetric, allowing for the compensation of potential bias.
### Binary classification

The data set consists of 60 realisations of the random variable pair \((X,Y)\), where half of the \(x\)-values are sampled uniformly from the interval \([0,0.2]\) and the other half are sampled uniformly from the interval \([0.8,1]\). The \(y\)-values are subsequently simulated using

\[Y\mid X=x\sim\text{Ber}(p(x)),\quad\text{with }p(x)=0.5+0.4\cos(6x).\]

On this training set, we train a fully connected network with three hidden layers consisting of 30 hidden units with elu activations functions. The final layer outputs a logit that is transformed using a sigmoid to yield a class probability. The network is trained for 300 epochs using a binary-crossy-entropy loss function and the Adam optimizer with a batch size of 32.

The resulting 95% confidence intervals of the predicted probability of class 1 are given in Figure 4. We carried out the experiment with two amounts of regularization to illustrate how this affects the result. For comparison, we also implemented an ensembling approach and MC-Dropout (see Lakshminarayanan et al. (2017) and Gal and Ghahramani (2016) for details). We used ten ensemble members and a standard dropout rate of 0.2. All networks were trained using the same training procedure.

We observe the same desirable properties as in the regression example. The intervals get much larger in regions with a limited amount of data, also when interpolating, and can become asymmetrical. Additionally, the intervals get smaller when we increase the amount of regularization. The model class becomes smaller (fewer parameters can be reached) making it more difficult for the model to change the predictions without affecting the likelihood. This, in turn, leads to smaller confidence intervals. If the model is overly regularized, it will become miss specified (\(\theta_{0}\notin\Theta\)) and the resulting confidence intervals will not be correct.

The other approaches do not share the same qualitative properties. The ensembling approach results in confidence intervals that are far too narrow. All ensemble members behave more or less the same, especially when interpolating. The MC-Dropout credible regions do not expand in the regions in between the data and also only moderately expand when extrapolating.

### Two-moon example

The data set consists of 80 data points, generated using the make_moons function from the scikit-learn package, which creates a binary classification problem with two interleaving half circles (Pedregosa et al., 2011).

We utilize the same network architecture as in the toy classification example. The network is trained for 500 epochs using the Adam optimizer with default learning rate and batch size of 32, while also applying \(l_{2}\)-regularization in each layer with a constant value of 1e-3.
Figure 1 presents the 95% confidence intervals for the predicted class probabilities. The results illustrate that DeepLR becomes extremely uncertain in regions farther from the data (i.e., the confidence interval for the class probability spans the full range of \([0,1]\)).

In contrast, both an ensembling approach and MC-Dropout report excessively high certainty in the upper left and lower right regions. The ensemble's behavior can be attributed to all ensemble members extrapolating in the

Figure 4: This figure illustrates DeepLR for a binary classification problem. The blue dots represent the training data, the dotted orange line represents the true probability of class 1 and the solid green line represents the predicted probability of class 1. The shaded blue region provides the 95% CI for the predicted probability of class 1. The top two figures (a) and (b) demonstrate the behavior of DeepLR for varying amounts of regularization. The more regularization, the smaller the class of admissible functions becomes, which naturally results in smaller intervals. Additionally, the top left figure demonstrates that the intervals expand when interpolating, a feature not shared by the dropout and ensembling approach.
same direction causing all ensemble members to report more or less the same class probability. For MC-dropout, a saturated sigmoid is causing the narrow credible intervals.

This comparison underscores the unique capability of DeepLR to provide more accurate uncertainty estimates in regions less well represented by data - a crucial capability in practical applications.

### MNIST binary example

For a more difficult task, we train a small convolutional network on the first two classes of the MNIST data set, consisting of handwritten digits. In this binary classification task, the 0's labeled as class 0 and the 1's as class 1.

The CNN architecture consists of two pairs of convolutional layers (with 28 filters and 3x3 kernels) and max-pooling layers (2x2 kernel), followed by a densely connected network with two hidden layers with 30 hidden units each and elu activation functions.

The network is trained for 10 epochs, using the SGD optimizer with a batch size of 32, default learning rate, and \(l_{2}\) regularization with a constant value of 1e-5, and binary cross-entropy loss function. The amounts of training time and regularization were determined from a manual grid search using an 80/20 split of training data. For the actual experiment, the entire training set was utilized.

Figure 5 presents 95% CIs for a number of different training points, test points, and OoD points. As shown, the OoD points have wider confidence intervals than the training and test points, reflecting greater uncertainty.

In an additional experiment, we rotated one of the test-points and created confidence intervals for the rotated images. As Figure 6 illustrates, increasing the rotation angle results in larger confidence intervals. This behavior is also seen with the ensembling and MC-dropout, but to a significantly lesser extent.

These larger intervals can be explained as follows: Our approach essentially asks the intuitive question how much the network can change the prediction for this new input without overly affecting the likelihood of the training data. A heavily rotated image of a number 1 deviates greatly from the typical input, utilizing pixels that are almost never used by the training inputs. It is therefore relatively straightforward for the network to change the prediction of this rotated input without dramatically lowering the likelihood of the training data. This, in turn, results in very large confidence intervals, accurately indicating that it is a very unfamiliar input.
**Fig. 5**: 95% CIs for different training points (first column), test points (second column), and OoD points (third column). Each subcaption displays the 95% CI for the probability of class 1 made with DeepLR (LR), MC-dropout (DR), and ensembling (EN). For in-distribution points (zeros and ones), all three methods provide extremely narrow intervals. For the OoD points, our method provides much larger CIs, a property also present in the other methods, albeit to a lesser extent.
Figure 6: This figure provides 95% CIs for different amounts of rotation. For rotations of 0 and 30 degrees, all methods produce very narrow confidence intervals, implying high certainty. At 60 degrees of rotation, DeepLR outputs high uncertainty whereas the ensembling approach and MC-dropout remain fairly certain. At 90 degrees of rotation, DeepLR outputs very high uncertainty, a confidence interval of [0.05, 1.00], a behavior also observed to a lesser extent in both ensembling and MC-dropout.
### CIFAR binary example

We extend the previous experiment to the CIFAR10 data set, using the first two classes - planes and cars - as a binary classification problem.

We use a CNN consisting of two pairs of convolutional layers (with 32 filters and 3x3 kernels) and max-pooling layers (2x2 kernel), followed by a densely connected network with three hidden layers with 30 hidden units each and elu activation functions.

The CNN is trained for 15 epochs using the SGD optimizer with default learning rate, a batch size of 32, and \(l_{2}\)-regularization with a constant value of 1e-5. The training time and regularization are determined the same way as for the MNIST experiment, using an 80/20 split of the training data.

Figure 7 presents 95% confidence intervals for several training, test, and OoD points. Our method is uncertain for out-of-distribution inputs. However, contrary to the MNIST example, we also see that the model is uncertain for various in-distribution points. Figures 7(e) and 7(g) provide examples of such uncertain predictions. The exact reason why the model is uncertain for those inputs remains speculation. Possible explanations might be the open hood of the car in (e) or the large amount of blue sky in (g). An interesting avenue for future work would be to investigate what specific features cause DeepLR to output greater uncertainty.

### Adversarial example

In addition to the previous experiments, we briefly tested how the method deals with adversarial examples (Goodfellow et al., 2014). Adversarial examples are modified inputs that are specifically designed to mislead the model, typically by adding small perturbations to the input. While virtually imperceptible to humans, these perturbations can dramatically alter the model's prediction.

We constructed adversarial versions of the two most confident test-inputs in Figure 7 using the FGSM method (Goodfellow et al., 2014), which works by using the gradients of the networks' loss function with respect to the input to create an input that maximizes the loss.

As illustrated in Figure 8, DeepLR demonstrates a higher uncertainty for the two adversarial inputs. This effect can be explained as follows. Although very similar to a human observer, an adversarial example is significantly different to a neural network. The difference was so large that both adversarial examples were wrongly classified by the network. Since these adversarial examples differ significantly from the training data, the network can more easily change the prediction at this location without changing the other predictions and thus the likelihood of the training data too much. This results in much larger confidence intervals.

These results provide an encouraging sign that our method could also be able to handle adversarial examples, offering the unique capability to create confidence intervals, detect OoD examples, and offer robustness against
adversarial examples, all with a single method. Unfortunately, the high computational cost of the method prevents more large-scale comparisons with other methods to more firmly establish this potential.

Figure 7: This figures gives 95% CIs for various training points (first column), test points (second column), and OoD points (third column). Each sub-caption provides the CIs for the probability of class 1 (cars) made with DeepLR (LR), MC-dropout (DR), and ensembling (EN). All methods demonstrate greater uncertainty than in the MNIST example, which is sensible as the CIFAR data set is significantly harder. Our method is very uncertain for all OoD points, which is not always the case for MC-dropout (c and i) and ensembling (i). We also observe relatively uncertain predictions by all methods for some in-distribution points, notably (e) and (g).
## 4 Discussion and conclusion

In this paper, we demonstrated the potential of a likelihood-ratio-based uncertainty estimate for neural networks. This approach is capable of producing asymmetric confidence intervals that are better motivated in cases where we have fewer data points than parameters, i.e. most deep learning applications.

The experimental results verify the theoretical advantages of a likelihood-ratio-based approach. The intervals are larger in regions with fewer data points, can get asymmetric in biased regions, and get larger for OoD and adversarial outputs. However, not being specifically designed for OoD detection or robustness against adversarial attacks, we do not claim it to be competitive in this regard against tailor-made alternatives.

While we made an effort to reduce it, our method still has some variance and can produce slightly different intervals upon repetition due to the randomness of the optimization procedure. This effect is greater for larger data sets where small differences can have a large effect on the test statistic.

Figure 8: Adversarial inputs generated by the FGSM method result in more uncertain predictions. The same network was used as for the CIFAR example illustrated in Figure 7. The intuition behind the larger CIs is as follows. While very similar to humans, the adversarial examples are significantly different to a neural network. This allows the network to change the prediction for the adversarial example without changing the predictions of the training data, resulting in larger confidence intervals.
Furthermore, it is essential that the model is well specified. This is the case for any model and not specific to our method. If the true density \(p_{\theta_{0}}\) cannot be reached, the resulting confidence intervals will surely be wrong. A model can be miss specified if it is overly regularized or if incorrect distributional assumptions are made (e.g. incorrectly assuming Gaussian noise).

In its current form, the high computational cost makes DeepLR unsuitable for many deep learning applications. A self-driving car that is approaching a cross-section cannot stop and wait for an hour until it has an uncertainty estimate. Nevertheless, a trustworthy uncertainty estimate may be critical in certain situations, or only a limited number of confidence intervals may be required. For instance, for medical applications, the extra computational time may be worthwhile. Alternatively, for some applications within astrophysics, only very few confidence intervals may be needed. If only a single interval is needed, for instance, our method is cheaper than an ensemble.

### Future work

Overall, our findings highlight the potential of a likelihood-ratio-based approach as a new branch of uncertainty estimation methods. We hope that our work will inspire further research in this direction. Several areas for improvement include:

* although not all
- applications and we hope that the proof of concept in this paper motivates further research in this direction that may result in reduced computational cost.
* Improved approximation of the test statistic: The calculation of the test statistic uses an approximation for the second term in equation (1). Further research could focus on finding better and possibly cheaper approximations. It may also be worthwhile to investigate the use of a Bartlett correction.
* Application to other machine-learning models: We applied this approach to neural networks. However, the methodology should also be applicable to other models, for example random forests. Especially for models that are relatively cheap to train, this approach could be very promising.
* Development of the theory on the distribution of the test statistic: It would be interesting to develop the theory surrounding the distribution of the test statistic in greater generality, possibly also when explicitly considering a regularization term. We constructed a reparametrization that showed that, under some assumptions, the test statistic has a \(\chi^{2}(1)\) distribution. It would be interesting to study these assumptions further.
* A better understanding of what causes DeepLR to become uncertain: We saw, for example, that various planes and cars had rather large accompanying confidence intervals. It would be interesting to study what causes certain input to become more uncertain than others.