# Large eddy simulation of a muffler with the high-order spectral difference
method

Matteo Parsani Current Institution: Computational Aerosciences Branch, NASA
Langley Research Center, Hampton, VA 23681, USA
(matteo.parsani@nasa.gov)Division of Computer, Electrical and Mathematical
Sciences & Engineering, King Abdullah University of Science and Technology,
Thuwal, 23955-6900, KSA Michael Bilka University of Notre Dame, Department of
Aerospace and Mechanical Engineering, 365 Fitzpatrick Hall, Notre Dame, IN
46556-5637, USA<EMAIL_ADDRESS>Chris Lacor Department of Mechanical
Engineering, Vrije Universiteit Brussel, Pleinlaan 2, 1050 Brussels, Belgium
<EMAIL_ADDRESS>
###### Abstract

The combination of the high-order accurate spectral difference discretization
on unstructured grids with subgrid-scale modelling is investigated for large
eddy simulation of a muffler at $Re=4.64\cdot 10^{4}$ and low Mach number. The
subgrid-scale stress tensor is modelled by the wall-adapting local eddy-
viscosity model with a cut-off length which is a decreasing function of the
order of accuracy of the scheme. Numerical results indicate that although the
high-order solver without subgrid-scale modelling is already able to capture
well the features of the flow, the coupling with the wall-adapting local eddy-
viscosity model improves the quality of the solution.

##  1 Introduction

Throughout the past two decades, the development of high-order accurate
spatial discretization has been one of the major fields of research in
computational fluid dynamics (CFD), computational aeroacoustics (CAA),
computational electromagnetism (CEM) and in general computational physics
characterized by linear and nonlinear wave propagation phenomena. High-order
accurate discretizations have the potential to improve the computational
efficiency required to achieve a desired error level. In fact, compared with
low order schemes, high order methods offer better wave propagation properties
and increased accuracy for a comparable number of degrees of freedom (DOFs).
Therefore, it may be advantageous to use such schemes for problems that
require very low numerical dissipation and small error levels [1]. Moreover,
since computational science is increasingly used as an industrial design and
analysis tool, high accuracy must be achieved on unstructured grids which are
required for efficient meshing. These needs have been the driving force for
the development of a variety of higher order schemes for unstructured meshes
such as the Discontinuous Galerkin (DG) method [2, 3], the Spectral Volume
(SV) method [4], the Spectral Difference (SD) method [5, 6], the Energy Stable
Flux Reconstruction [7] and many others.

In this study we focus on a SD solver for unstructured hexahedral grids
(tensorial cells). The SD method has been proposed as an alternative high
order collocation-based method using local interpolation of the strong form of
the equations. Therefore, the SD scheme has an important advantage over
classical DG and SV methods, that no integrals have to be evaluated to compute
the residuals, thus avoiding the need for costly high-order accurate
quadrature formulas.

Although the formulation of high-order accurate spatial discretization is now
fairly mature, their application for the simulation of general turbulent flows
implies that particular attention has still to be paid to subgrid-scale (SGS)
models. So far, the combination of the SD method with SGS models for LES has
not been widely investigated. In 2010, Parsani et al. [8] reported the first
implementation in study of a two-dimensional (2D) third-order accurate SD
solver coupled with the Wall-Adapting Local Eddy-viscosity (WALE) model [14]
and a cut-off length which is a decreasing function of the order of accuracy.
A successful extension of that approach to a three-dimensional (3D) second-
order accurate SD solver has been reported in [12]. Very recently, Lodato and
Jameson [13] have presented an alternative technique to model the unresolved
scales in the flow field: A structural SGS approach with the WALE Similarity
Mixed model (WSM), where constrained explicit filtering represents a key
element to approximate subgrid-scale interactions. The performance of such an
algorithm has been also satisfactory.

In this study, we couple for the first time the approach proposed in [8] with
a 3D fourth-order accurate SD solver, for the simulation of the turbulent flow
in an industrial-type muffler at $Re=4.64\cdot 10^{4}$. The goal is to
investigate if the coupling of a high-order SD scheme with a sub-grid closure
model improves the quality of the results when the grid resolution is
relatively low. The latter requirement is often desirable when a high-order
accurate spatial discretization is used.

##  2 Physical model and numerical algorithm

In this study the system of the Navier-Stokes equations for a compressible
flow are discretized in space using the SD method and the subgrid-scale stress
tensor is modelled by the WALE approach.

###  2.1 Filtered Navier-Stokes equations

The three physical conservation laws for a general Newtonian fluid, i.e., the
continuity, the momentum and energy equations, are introduced using the
following notation: $\rho$ for the mass density, $\vec{u}\in\mathbb{R}^{dim}$
for the velocity vector in a physical space with $dim$ dimensions, $P$ for the
static pressure and $E$ for the specific total energy which is related to the
pressure and the velocity vector field by
$E=\frac{1}{\gamma-1}\frac{P}{\rho}+\frac{|\vec{u}|^{2}}{2}$, where $\gamma$
is the constant ratio of specific heats and it is $1.4$ for air in standard
conditions.

The system, written in divergence form and equipped with suitable initial-
boundary conditions, is

$\frac{\partial\textbf{w}}{\partial
t}+\vec{\nabla}\cdot\left(\vec{\textbf{f}}_{C}\left(\textbf{w}\right)-\vec{\textbf{f}}_{D}\left(\textbf{w},\vec{\nabla}\textbf{w}\right)\right)=\frac{\partial\textbf{w}}{\partial
t}+\vec{\nabla}\cdot\vec{\textbf{f}}=0,$ (1)

where
$\textbf{w}=\left(\overline{\rho},\overline{\rho}\tilde{\vec{u}},\overline{\rho}\tilde{E}\right)^{T}$
is the vector of the filtered conservative variables and
$\vec{\textbf{f}}_{C}=\vec{\textbf{f}}_{C}\left(\textbf{w}\right)$ and
$\vec{\textbf{f}}_{D}=\vec{\textbf{f}}_{D}\left(\textbf{w},\vec{\nabla}\textbf{w}\right)$
represent the convective and the diffusive fluxes, respectively. Here the
symbols $(\overline{\cdot})$ and $(\tilde{\cdot})$ represent the spatially
filtered field and the Favre filtered field defined as
$\tilde{\vec{u}}=\overline{\rho\vec{u}}/\overline{\rho}$.

In a general 3D ($dim=3$) Cartesian space,
$\vec{x}=\left[x_{1},x_{2},x_{3}\right]^{T}$, the components of the flux
vector
$\vec{\textbf{f}}\left(\textbf{w},\vec{\nabla}\textbf{w}\right)=\left[\textbf{f}_{1},\textbf{f}_{2},\textbf{f}_{3}\right]^{T}$
are given by

$\textbf{f}_{1}=\left(\begin{array}[]{c}\overline{\rho}\tilde{u}_{1}\\\
\overline{\rho}\tilde{u}_{1}^{2}+\overline{P}-\tilde{\sigma}_{11}+\tau_{11}^{sgs}\\\
\overline{\rho}\tilde{u}_{1}\tilde{u}_{2}-\tilde{\sigma}_{21}+\tau_{21}^{sgs}\\\
\overline{\rho}\tilde{u}_{1}\tilde{u}_{3}-\tilde{\sigma}_{31}+\tau_{31}^{sgs}\\\
\tilde{u}_{1}\left(\overline{\rho}\tilde{E}+\overline{P}\right)-\tilde{u}_{1}\left(\tilde{\sigma}_{11}-\tau_{11}^{sgs}\right)-\tilde{u}_{2}\left(\tilde{\sigma}_{21}-\tau_{21}^{sgs}\right)-\tilde{u}_{3}\left(\tilde{\sigma}_{31}-\tau_{31}^{sgs}\right)-c_{P}\frac{\mu}{Pr}\frac{\partial\tilde{T}}{\partial
x_{1}}+q_{1}^{sgs}\end{array}\right),$
$\textbf{f}_{2}=\left(\begin{array}[]{c}\overline{\rho}\tilde{u}_{2}\\\
\overline{\rho}\tilde{u}_{1}\tilde{u}_{2}-\tilde{\sigma}_{12}+\tau_{12}^{sgs}\\\
\overline{\rho}\tilde{u}_{2}^{2}+\overline{P}-\tilde{\sigma}_{22}+\tau_{22}^{sgs}\\\
\overline{\rho}\tilde{u}_{2}\tilde{u}_{3}-\tilde{\sigma}_{32}+\tau_{32}^{sgs}\\\
\tilde{u}_{2}\left(\overline{\rho}\tilde{E}+\overline{P}\right)-\tilde{u}_{1}\left(\tilde{\sigma}_{12}-\tau_{12}^{sgs}\right)-\tilde{u}_{2}\left(\tilde{\sigma}_{22}-\tau_{22}^{sgs}\right)-\tilde{u}_{3}\left(\tilde{\sigma}_{32}-\tau_{32}^{sgs}\right)-c_{P}\frac{\mu}{Pr}\frac{\partial\tilde{T}}{\partial
x_{2}}+q_{2}^{sgs}\end{array}\right),$
$\textbf{f}_{3}=\left(\begin{array}[]{c}\overline{\rho}\tilde{u}_{3}\\\
\overline{\rho}\tilde{u}_{1}\tilde{u}_{3}-\tilde{\sigma}_{13}+\tau_{13}^{sgs}\\\
\overline{\rho}\tilde{u}_{2}\tilde{u}_{3}-\tilde{\sigma}_{23}+\tau_{23}^{sgs}\\\
\overline{\rho}\tilde{u}_{3}^{2}+\overline{P}-\tilde{\sigma}_{33}+\tau_{33}^{sgs}\\\
\tilde{u}_{3}\left(\overline{\rho}\tilde{E}+\overline{P}\right)-\tilde{u}_{1}\left(\tilde{\sigma}_{13}-\tau_{13}^{sgs}\right)-\tilde{u}_{2}\left(\tilde{\sigma}_{23}-\tau_{23}^{sgs}\right)-\tilde{u}_{3}\left(\tilde{\sigma}_{33}-\tau_{33}^{sgs}\right)-c_{P}\frac{\mu}{Pr}\frac{\partial\tilde{T}}{\partial
x_{3}}+q_{3}^{sgs}\end{array}\right),$

where $c_{P}$, $\mu$, $Pr$ and $T$ represent respectively the specific heat
capacity at constant pressure, the dynamic viscosity, the Prandtl number and
the temperature of the fluid. Moreover, $\sigma_{ij}$ represents the
$ij-$component of the resolved viscous stress tensor [15].

Both momentum and energy equations differ from the classical fluid dynamic
equations only for two terms which take into account the contributions from
the unresolved scales. These contributions, represented by the specific
subgrid-scale stress tensor $\tau_{ij}^{sgs}$ and by the subgrid heat flux
vector defined $q_{i}^{sgs}$, appear when the spatial filter is applied to the
convective terms [15]. The interactions of $\tau_{ij}^{sgs}$ and $q_{i}^{sgs}$
with the resolved scales have to be modeled through a subgrid-scale closure
model because they cannot be determined using only the resolved flow field w.

####  2.1.1 The wall-adapted local eddy-viscosity closure model

The smallest scales present in the flow field and their interaction with the
resolved scales have to be modeled through the subgrid-scale term
$\tau_{ij}^{sgs}$. The most common approach to model such a tensor is based on
the eddy-viscosity concept in which one assumes that the residual stress is
proportional to a measure of the filtered local strain rate [15], which is
defined as follows:

$\tau_{ij}^{sgs}-\tau_{kk}^{sgs}\delta_{ij}=-2\,\overline{\rho}\,\nu_{t}\left(\tilde{S}_{ij}-\frac{\delta_{ij}}{3}\tilde{S}_{kk}\right).$
(2)

In the WALE model, it is assumed that the eddy-viscosity $\nu_{t}$ is
proportional to the square of the length scale of the cut-off length (or width
of the grid filter) and the filtered local rate of strain. Although the model
was originally developed for incompressible flows, it can also be used for
variable density flows by giving the formulation as follows

$\nu_{t}=\left(C\Delta\right)^{2}\left|\tilde{S}\right|.$ (3)

Here $\left|\tilde{S}\right|$ is defined as

$\left|\tilde{S}\right|=\frac{\left[\tilde{S}_{ij}^{d}\,\tilde{S}_{ij}^{d}\right]^{3/2}}{\left[\tilde{S}_{ij}\,\tilde{S}_{ij}\right]^{5/2}+\left[\tilde{S}_{ij}^{d}\,\tilde{S}_{ij}^{d}\right]^{5/4}},$
(4)

where $\tilde{S}_{ij}^{d}$ is the traceless symmetric part of the square of
the resolved velocity gradient tensor
$\tilde{g}_{ij}=\frac{\partial\tilde{u}_{i}}{\partial x_{j}}$. Note that in
Equation (3) $\Delta$, i.e., the cut-off length, is an unknown function. Often
the cut-off length is taken proportional to the smallest resolvable length
scale of the discretization. In the present work, the definition of the grid
filter function is given in Section 2.2, where the SD method is discussed.

###  2.2 Spectral difference method

Consider a problem governed by a general system of conservation laws given by
Equation (1) and valid on a domain $\Omega\subset\mathbb{R}^{dim}$ with
boundary $\partial\Omega$ and completed with consistent initial and boundary
conditions. The domain is divided into $N$ non-overlapping cells, with cell
index $i$.

In order to achieve an efficient implementation of the SD method, all
hexahedral cells in the physical domain are mapped into cubic elements using
local coordinates $\vec{\xi}=\left[\xi_{1},\xi_{2},\xi_{3}\right]^{T}$. Such a
transformation is characterized by the Jacobian matrix
$\vec{\vec{\left.\mathrm{J}\right.}}_{i}$ with determinant
$det(\vec{\vec{\left.\mathrm{J}\right.}}_{i})$. Therefore, system (1) can be
written in the mapped coordinate system as

$\frac{\partial\mathbf{w}_{i}^{\hskip 1.13791pt\vec{\xi}}}{\partial
t}=-\frac{\partial\mathbf{f}_{1,i}^{\hskip
1.13791pt\vec{\xi}}}{\partial\xi_{1}}-\frac{\partial\mathbf{f}_{2,i}^{\hskip
1.13791pt\vec{\xi}}}{\partial\xi_{2}}-\frac{\partial\mathbf{f}_{3,i}^{\hskip
1.13791pt\vec{\xi}}}{\partial\xi_{3}}=-\vec{\nabla}^{\hskip
1.13791pt\vec{\xi}}\cdot\vec{\textbf{f}}^{\hskip 1.13791pt\vec{\xi}}_{i},$ (5)

where $\mathbf{w}^{\hskip 1.13791pt\vec{\xi}}_{i}\equiv
det(\vec{\vec{\left.\mathrm{J}\right.}}_{i})\,\mathbf{w}$ and
$\vec{\nabla}^{\hskip 1.13791pt\vec{\xi}}$ are the conserved variables and the
generalized divergence differential operator in the mapped coordinate system,
respectively.

For a $\left(p+1\right)$-th-order accurate $dim$-dimensional scheme, $N^{s}$
_solution collocation points_ with index $j$ are introduced at positions
$\vec{\xi}_{j}^{s}$ in each cell $i$, with $N^{s}$ given by
$N^{s}=\left(p+1\right)^{dim}$. Given the values at these points, a polynomial
approximation of degree $p$ of the solution in cell $i$ can be constructed.
This polynomial is called the _solution polynomial_ and is usually composed of
a set of Lagrangian basis polynomial $L_{j}^{s}\left(\vec{\xi}\right)$ of
degree $p$:

$\mathbf{W}_{i}\left(\vec{\xi}\right)=\sum_{j=1}^{N^{s}}\mathbf{W}_{i,j}\,L_{j}^{s}\left(\vec{\xi}\right).$
(6)

Therefore, the unknowns of the SD method are the interpolation coefficients
$\mathbf{W}_{i,j}=\mathbf{W}_{i}\left(\vec{\xi}_{j}^{s}\right)$ which are the
approximated values of the conserved variables $\mathbf{w}_{i}$ at the
solution points.

The divergence of the mapped fluxes $\vec{\nabla}^{\hskip
1.13791pt\vec{\xi}}\cdot\vec{\textbf{f}}^{\hskip 1.13791pt\vec{\xi}}$ at the
solution points is computed by introducing a set of $N^{f}$ flux collocation
points with index $l$ and at positions $\vec{\xi}_{l}^{f}$, supporting a
polynomial of degree $p+1$. The evolution of the mapped flux vector
$\vec{\textbf{f}}^{\hskip 1.13791pt\vec{\xi}}$ in cell $i$ is then
approximated by a flux polynomial $\vec{\textbf{F}}^{\hskip
1.13791pt\vec{\xi}}_{i}$, which is obtained by reconstructing the solution
variables at the flux points and evaluating the fluxes
$\vec{\textbf{F}}^{\hskip 1.13791pt\vec{\xi}}_{i,l}$ at these points. The flux
is also represented by a Lagrange polynomial:

$\vec{\textbf{F}}^{\hskip
1.13791pt\vec{\xi}}_{i}\left(\vec{\xi}\right)=\sum_{l=1}^{N^{f}}\vec{\textbf{F}}^{\hskip
1.13791pt\vec{\xi}}_{i,l}\,L_{l}^{f}\left(\vec{\xi}\right),$ (7)

where the coefficients of the flux interpolation are defined as

$\vec{\textbf{F}}^{\hskip
1.13791pt\vec{\xi}}_{i,l}=\begin{cases}\vec{\textbf{F}}^{\hskip
1.13791pt\vec{\xi}}_{i}\left(\vec{\xi}_{l}^{f}\right),&\quad\vec{\xi}_{l}^{f}\in\Omega_{i},\\\
\vec{\textbf{F}}_{\mathrm{num}}^{\hskip
1.13791pt\vec{\xi}}\left(\vec{\xi}_{l}^{f}\right),&\quad\vec{\xi}_{l}^{f}\in\partial\Omega_{i}.\end{cases}$
(8)

Here $\vec{\textbf{F}}_{\mathrm{num}}^{\hskip 1.13791pt\vec{\xi}}$ is the
numerical flux vector at the cell interface. In fact, the solution at a face
is in general not continuous and requires the solution of a Riemann problem to
maintain conservation at a cell level (i.e., the flux component normal to a
face $\vec{\textbf{F}}_{\mathrm{num}}^{\hskip
1.13791pt\vec{\xi}}\cdot\vec{n}^{\hskip 1.13791pt\vec{\xi}}$ must be
continuous between two neighboring cells). Approximate Riemann solvers are
typically used (e.g. Rusanov Riemann solver). The tangential component of
$\vec{\textbf{F}}_{\mathrm{num}}^{\hskip 1.13791pt\vec{\xi}}$ is usually taken
from the interior cell.

Taking the divergence of the flux polynomial $\vec{\nabla}^{\hskip
1.13791pt\vec{\xi}}\cdot\vec{\textbf{F}}^{\hskip 1.13791pt\vec{\xi}}_{i}$ in
the solution points results in the following modified form of (5), describing
the evolution of the conservative variables at the solution points:

$\frac{d\mathbf{W}_{i,j}}{dt}=-\left.\vec{\nabla}\cdot\vec{\textbf{F}}_{i}\right|_{j}=-\frac{1}{J_{i,j}}\left.\vec{\nabla}^{\hskip
1.13791pt\vec{\xi}}\cdot\vec{\textbf{F}}^{\hskip
1.13791pt\vec{\xi}}_{i}\right|_{j}=\mathbf{R}_{i,j},$ (9)

where $\vec{\textbf{F}}_{i}$ is the flux polynomial vector in the physical
space whereas $\mathbf{R}_{i,j}$ is the SD residual associated with
$\mathbf{W}_{i,j}$. This is a system of ODEs, in time, for the unknowns
$\mathbf{W}_{i,j}$. In this work, the optimized explicit eighteen-stages
fourth-order Runge-Kutta schemes presented in [16] is used to solve such a
system at each time step.

####  2.2.1 Solution and flux points distributions

In 2007, Huynh [9] showed that for quadrilateral and hexahedral cells, tensor
product flux point distributions based on a one-dimensional flux point
distribution consisting of the end points and the Legendre-Gauss quadrature
points lead to stable schemes for arbitrary order of accuracy. In 2008, Van
den Abeele et al. [10] showed an interesting property of the SD method, namely
that it is independent of the positions of its solution points in most general
circumstances. This property implies an important improvement in efficiency,
since the solution points can be placed at flux point positions and thus a
significant number of solution reconstructions can be avoided. Recently, this
property has been proved by Jameson [11].

####  2.2.2 Cut-off length $\Delta$

In Section 2.1.1 we have seen that in the WALE model the cut-off length
$\Delta$ is used to compute the turbulent eddy-viscosity $\nu_{t}$, i.e.,
Equation (3). Following the approach presented in [8], for each cell with
index $i$ and each flux points with index $l$ and positions
$\boldsymbol{\xi}_{l}^{f}$, we use the following definition of filter width

$\Delta_{i,l}=\left[\frac{1}{N^{s}}det\left(\left.\vec{\vec{\left.\mathrm{J}\right.}}_{i}\right|_{\boldsymbol{\xi}_{l}^{f}}\right)\right]^{1/dim}=\left(\frac{det(J_{i,l})}{N^{s}}\right)^{1/dim}.$
(10)

Notice that the cell filter width is not constant in one cell, but it varies
because the Jacobian matrix is a function of the positions of the flux points.
Moreover, for a given mesh, the number of solution points depends on the order
of the SD scheme, so that the grid filter width decreases by increasing the
polynomial order of the approximation.

##  3 Numerical results

The main purpose of this section is to evaluate the accuracy and the
reliability of the fourth-order SD-LES solver for simulating a 3D turbulent
flow in an industrial-type muffler. The results are compared with the particle
image velocimetry (PIV) measurement performed at the Department of
Environmental and Applied Fluid Dynamics of the von Karman Institute for Fluid
Dynamics [17]. In Figure 1, the geometry of the muffler and its characteristic
dimensions are illustrated, where the flow is from left to right.

$5\,d$$7.5\,d$$5\,d$$4.75\,d$

$d$

$0.625\,d$

$0.625\,d$

$\chi_{3}$$\chi_{2}$

Figure 1: Configuration of the 3D muffler test case.

At the inlet, mass density and velocity profiles are imposed. The inlet
velocity profile in the $x_{3}$ direction is given by

$u_{3}=u_{max}\left\\{\frac{1}{2}-\frac{1}{2}\tanh{\left[2.2\left(\frac{r}{d/2}-\frac{d/2}{r}\right)\right]}\right\\}.$

At the outlet only the pressure is prescribed. In accordance to the
experiments, the inlet Mach number and the Reynolds number, based on maximum
velocity at the inlet $u_{max}$ and the diameter of the inlet/outlet $d$
($d=4\,cm$), are set respectively to $M_{inlet}=0.05$ and $Re=4.64\cdot
10^{4}$.

The flow is computed using fourth-order ($p=3$) SD scheme on a grid with
$36,612$ hexahedral elements which was generated with the open source software
Gmsh [18]. Second-order boundary elements are used to approximate the curved
geometry. The total number of DOFs is approximately $2.3\cdot 10^{6}$ (i.e.,
$36,612\cdot(p+1)^{3}$). The maximum CFL number used for the computations
started from $0.1$ and increased up to $0.65$. After the flow field was fully
developed, time averaging is performed for a period corresponding to about 25
flow-through times.

The computation is validated on the center plane of the expansion coinciding
with the center planes of the inlet and outlet pipes using the PIV results
from [17]. All of the measurements are taken on the symmetrical center plane
of the muffler. The reference cross section corresponds to the entrance of the
expansion chamber. It should be noted that the circular nature of the geometry
acts as a lens causing a change in magnification in the radial direction
($x_{2}$) which prevents from capturing images close to the wall. It is found
that outside $~{}1\,cm$ from the wall the magnification effect is negligible
and as the mean stream-wise direction is in the direction of constant
magnification and has only little effect on the particle correlations no
corrections are deemed necessary.

In Figure 2, the non-dimensional mean velocity profile in the axial direction
$\langle\tilde{u}_{3}\rangle/u_{max}$ is shown for four different cross
sections in the expansion chamber, where the PIV measurements were done. In
this figure, the PIV data are also plotted for comparison. Figure 3 shows the
non-dimensional Reynolds stress $\langle
u_{2}^{\prime}u_{3}^{\prime}\rangle/{u_{max}^{2}}$ at the same cross sections.
Although the high-order implicit LES is already able to capture well the
features of the flow field, the use of the WALE model improves the results. In
particular, when the SGS model is active, the local extrema of the time-
averaged velocity profiles and the second-order statistical moment (which get
fairly oscillatory by moving far away from the inlet pipe) are better
captured.

##  4 Conclusions

The fourth-order SD method in combination with the WALE model and the variable
filter width performed well. The numerical results confirm that the model is
correctly accounting for the unresolved shear stress computed from the
resolved field, for the present internal flow. However, it should be noted
that the SD discretization without subgrid-scale modelling also worked rather
well, at least for the grid resolution used in this study.

Work is currently under way to test both approaches for different orders of
accuracy, grid resolutions and other realistic turbulent flows. We believe
that the flexibility of the high-order SD scheme on unstructured grids
together with the development of robust sub-grid closure models for highly
separated flows and efficient grid generators for high-order accurate schemes
will allow to perform LES of industrial-type flows in the near future.

(a) $1d$ downstream. (b) $4d$ downstream. (c) $6d$ downstream. (d) $7d$
downstream.

Figure 2: Time-averaged velocity profile in the axial direction
$\langle\tilde{u}_{3}\rangle/u_{max}$ at four cross sections in the expansion
chamber, obtained with fourth-order ($p=3$) SD-LES method. Comparison with
experimental measurements (PIV) [17].

(a) $1d$ downstream. (b) $4d$ downstream. (c) $6d$ downstream. (d) $7d$
downstream.

Figure 3: Reynolds stress $\langle
u_{2}^{\prime}u_{3}^{\prime}\rangle/u_{max}^{2}$ in the axial direction at
four cross sections in the expansion chamber, obtained with fourth-order
($p=3$) SD-LES method. Comparison with experimental measurements (PIV) [17].

## Acknowledgement

The authors would like to thank Professor David I. Ketcheson and Professor
Mark H. Carpenter for their support. This research used partially the
resources of the KAUST Supercomputing Laboratory and was supported in part by
an appointment to the NASA Postdoctoral Program at Langley Research Center,
administered by Oak Ridge Associates Universities. These supports are
gratefully acknowledged.

## References

  * [1] Wang, Z. J.: Adaptive High-order Methods in Computational Fluid Dynamics (Advances in Computational Fluid Dynamics), World Scientific Publishing Company (2011).
  * [2] Busch, K., König, M., and Niegemann, J.: Discontinuous Galerkin methods in nanophotonics. Laser Photonics Rev. 5(6), 773–809 (2011).
  * [3] Hesthaven, Jan S., and Warburton, Tim: Nodal Discontinuous Galerkin Methods: Algorithms, Analysis, and Applications. Texts in Applied Mathematics, 54. Computational Science & Engineering, Springer Publishing Company, Incorporated (2007).
  * [4] Sun, Y., Wang, Z. J., and Liu Y.: Spectral (finite) volume method for conservation laws on unstructured grids VI: extension to viscous flow. J. of Comput. Phys., 215(1), 41–58 (2006).
  * [5] May, G., and Jameson, A.: A spectral difference method for the Euler and Navier-Stokes equations on unstructured meshes. _AIAA paper_ 2006-304. 44th AIAA Aerospace Sciences Meeting, Reno, Nevada, USA, January 9-12, 2006.
  * [6] Sun, Y., Wang, Z. J., and Liu, Y.: High-order multidomain spectral difference method for the Navier-Stokes equations on unstructured hexahedral grids. Commun. Comput. Phys, 2(2), 310–333 (2007).
  * [7] Castonguay, P., Vincent P., and Jameson, A.: Application of High-Order Energy Stable Flux Reconstruction Schemes to the Euler Equations. _AIAA paper_ 2011-686. 49th AIAA Aerospace Sciences Meeting including the New Horizons Forum and Aerospace Exposition, Orlando, Florida, USA, January 4-7, 2011.
  * [8] Parsani, M., Ghorbaniasl, G., Lacor C., and Turkel, E.: An implicit high-order spectral difference approach for large eddy simulation. J. Comput. Phys., 229(14), 5373–5393 (2010).
  * [9] Huynh, H. T.: A flux reconstruction approach to high-order schemes including discontinuous Galerkin methods. _AIAA paper_ 2007-4079. 18th AIAA Computational Fluid Dynamics Conference, Miami, Florida, USA, June 25-28, 2007.
  * [10] Van den Abeele, K., Lacor, C., and Wang, Z. J.: On the stability and accuracy of the spectral difference method, J. Sci. Comput., 37(2), 162–188 (2008).
  * [11] Jameson, A.: A proof of the stability of the spectral difference method for all orders of accuracy, J. Sci. Comput., 45(1-3), 348–358 (2010).
  * [12] Parsani, M., Ghorbaniasl, G., and Lacor C.: Validation and application of an high-order spectral difference method for flow induced noise simulation. J. Comput. Acoust., 19(3), 241–268 (2011).
  * [13] Lodato, G., and Jameson, A.: LES modeling with high-order flux reconstruction and spectral difference schemes. _ICCFD paper_ 2201\. 7th ICCFD Conference, Big Island, Hawaii, July 9-13, 2012.
  * [14] Nicoud, F., and Ducros, F.: Subgrid-scale stress modelling based on the square of the velocity gradient tensor. Flow Turbul. Combust., 62(3), 183–200 (1999).
  * [15] Pope, Stephen B.: Turbulent flows. Cambridge University Press (2003).
  * [16] Parsani, M., Ketcheson, David I., and Deconinck, W.: Optimized explicit Runge-Kutta schemes for the spectral difference method applied to wave propagation problems, SIAM J. Sci. Comput, 35(2):A957-A986 (2013).
  * [17] Bilka, M., and Anthoine, J.: Experimental investigation of flow noise in a circular expansion using PIV and acoustic measurements. _AIAA paper_ 2008-2952. 14th AIAA/CEAS Aeroacoustics Conference, Vancouver, British Columbia, Canada, May 5-6, 2008.
  * [18] Geuzaine, C., and Remacle J.-F.: Gmsh: a three-dimensional finite element mesh generator with built-in pre- and post-processing facilities. Int. J. Numer. Meth. Eng., 79(11), 1309–1331 (2009).


# Adapting Conformal Prediction to Distribution Shifts Without Labels

Kevin Kasa Work done while visiting Harvard University. School of Engineering,
University of Guelph  Vector Institute for Artificial Intelligence  Zhiyu
Zhang  Harvard University  Heng Yang  Harvard University  Graham W. Taylor
School of Engineering, University of Guelph  Vector Institute for Artificial
Intelligence

###### Abstract

Conformal prediction (CP) enables machine learning models to output prediction
sets with guaranteed coverage rate, assuming exchangeable data. Unfortunately,
the exchangeability assumption is frequently violated due to distribution
shifts in practice, and the challenge is often compounded by the lack of
ground truth labels at test time. Focusing on classification in this paper,
our goal is to improve the quality of CP-generated prediction sets using only
unlabeled data from the test domain. This is achieved by two new methods
called ECP and EACP, that adjust the score function in CP according to the
base model’s uncertainty on the unlabeled test data. Through extensive
experiments on a number of large-scale datasets and neural network
architectures, we show that our methods provide consistent improvement over
existing baselines and nearly match the performance of supervised algorithms.

##  1 Introduction

Advances in deep learning are fundamentally changing the autonomous decision
making pipeline. While most works have focused on accurate point predictions,
calibrating the _uncertainty_ of the model is arguably as important. Taking
autonomous driving as an example: if a detection model predicts the existence
of an obstacle, it would be reasonable to take different maneuvering
strategies depending on the _confidence_ of the prediction. But is that
reliable? In a possible failure mode, the model could report $60\%$ (resp.
$99\%$) confidence, but the _probability_ of an obstacle actually showing up
is $99\%$ (resp. $60\%$). Such discrepancy between the model’s own uncertainty
evaluation and the ground truth probability (or post-hoc frequency) is
actually very common (Guo et al., 2017; Liang et al., 2023), and can
significantly compromise the safety in downstream decision making.

_Conformal Prediction_ (CP; Vovk et al., 2005) is a simple post-hoc and
distribution-free method that addresses this problem. Given a fixed black-box
machine learning model (called the _base model_) and a covariate
$x_{\mathrm{test}}$, the goal of CP is to generate a small prediction set
$\mathcal{C}_{\mathrm{test}}$ that contains (or, _covers_) the unknown ground
truth label $y_{\mathrm{test}}$ with a pre-specified probability. Crucially,
CP relies on the assumption that the distribution of the data stream is
_exchangeable_ (a weaker variant of i.i.d.), which allows the fairly
straightforward inference of $y_{\mathrm{test}}$ from $x_{\mathrm{test}}$ and
the base model’s performance on a pre-collected _calibration dataset_. Note
that the target label $y_{\mathrm{test}}$ does not need to be revealed after
the set prediction is made: exchangeability together with a large enough
calibration dataset is sufficient to ensure the desirable coverage
probability. This is particularly important for autonomous decision making,
where real-time data annotation is expensive or even infeasible.

However, real-world data streams are usually corrupted by all sorts of
_distribution shifts_ , violating the exchangeability assumption. Even when
the data stream itself is exchangeable, we often want to continually update
the base model rather than keeping it fixed, and this can be effectively
understood as a distribution shift in the CP context. In such cases, simply
applying exchangeability-based CP methods could fail miserably (Tibshirani et
al., 2019; Bhatnagar et al., 2023a; Kasa and Taylor, 2023). Therefore, making
CP compatible with distribution shifts has become a focal point of recent
works.

A number of solutions have been proposed, but the key challenge still remains.
For example, Gibbs and Candes (2021) formulated the connection between CP and
_Online Convex Optimization_ (OCO; Zinkevich, 2003), and the latter is able to
handle arbitrarily distribution-shifted environments. The weakness is that
ground truth labels are now required at test time, as opposed to the standard
CP procedure. On the other direction, there are CP methods that combat
distribution shifts without test time labels (Tibshirani et al., 2019; Roth,
2022; Barber et al., 2023; Cauchois et al., 2024), but they assume the
distribution shifts are “easy”, such that even without labels, we can still
rigorously infer the test distribution to a certain extent using the labeled
calibration dataset. Overall, it appears that handling both difficulties –
distribution shifts and the lack of test time labels – is a formidable but
important challenge remaining in the CP literature.

#### Contributions

In this paper, we focus on classification, and address the above problem by
exploiting the uncertainty evaluation of the base model itself. Although such
a quantity is not necessarily calibrated in a strict sense, it has been
empirically shown to strongly correlate with the magnitude of distribution
shifts (Wang et al., 2021), thus providing a valuable way to probe the test
distribution without label access. Under this high level idea, we make the
following contributions:

  * •

First, we propose a new method named ECP (Entropy scaled Conformal
Prediction), which scales the _score function_ in CP by an “entropy
quantile”111Enforced to be larger than 1 to ensure that the prediction sets do
not become smaller on in-distribution data. of the base model on the test
distribution.

In the CP jargon, the score function is determined by the fixed base machine
learning model, and assigns each candidate label a “propensity score” after
observing the covariate $x_{\mathrm{test}}$. Then, the CP prediction set
$\mathcal{C}_{\mathrm{test}}$ simply includes all the labels whose score is
above a certain threshold.222Since we mainly focus on classification, the
score functions are assumed to be _positively oriented_ (Sadinle et al.,
2019). By scaling up the score function, ECP makes the prediction sets larger,
which naturally corresponds to the intuition that the uncertainty of the base
model should be inflated under distribution shifts. Moreover, the amount of
such inflation is correlated with the magnitude of the distribution shift,
through the use of the entropy quantile.

  * •

Second, we refine the above using techniques from unsupervised _Test Time
Adaptation_ (TTA) (Niu et al., 2022), and the resulting method is named as
EACP (Entropy base-Adapted Conformal Prediction). The key idea is that while
ECP keeps the base model fixed during testing, we may concurrently update it
using _entropy minimization_ (Grandvalet and Bengio, 2004; Wang et al., 2021)
– a widely adopted idea in unsupervised TTA, alongside the aforementioned
entropy scaling. This “adaptively” reduces the scaling effect that ECP applies
on the score function, thus making the prediction sets of ECP smaller.

  * •

Finally, despite the loss of statistical guarantees, we evaluate the proposed
methods on a wide range of large-scale datasets under distribution shifts, as
well as different neural network architectures. We find that exchangeability-
based CP (with and without TTA on the base model) consistently leads to lower-
than-specified coverage frequency, but our methods can effectively mitigate
this under-coverage issue while keeping the sizes of the prediction sets
moderate. Furthermore, our methods also significantly improve the prediction
sets generated by the base model itself (without CP). It shows that by
bridging the CP procedure (which is statistically sound) and the base model’s
own uncertainty evaluation (which is often informative), our methods enjoy the
practical benefit from both worlds.

##  2 Related works

#### CP in decision making

Our interest in the considered setting – distribution shifts without test time
labels – is mainly motivated by the growing applications of CP in autonomous
decision making. A very much incomplete list: see (Lekeufack et al., 2023) for
a generic treatment; (Lindemann et al., 2023) for trajectory optimization in
robotics; (Yang and Pavone, 2023; Gao et al., 2024) for 3D vision; (Kumar et
al., 2023; Quach et al., 2024) for large language models (LLM); and (Ren et
al., 2023) for LLM-powered robotics.

#### CP under distribution shifts

Considerable efforts have been devoted to developing CP methods robust to
distribution shifts. We now discuss two possible directions and their
respective limitations.

  * •

The first direction does not require test time labels, but the distribution
shift is assumed to be simple in some sense. For example, Tibshirani et al.
(2019) studied CP under _covariate shifts_ , where the distribution of the
label $Y$ conditioned on the covariate $X$ remains unchanged. Here, it
suffices to use the classical _likelihood ratio reweighting_ on the
calibration dataset, but due to the difficulty of likelihood ratio estimation,
one could also take a robust optimization perspective (Cauchois et al., 2024)
and apply _multicalibration_ on the possible test domain, as sketched in
(Roth, 2022, Chapter 8). Barber et al. (2023) generalized the reweighting
approach to handle mild general distribution shifts, but the choice of weights
is generally unclear in practice.

  * •

The second direction is connecting CP to adversarial online learning. A line
of works (Gibbs and Candes, 2021; Angelopoulos et al., 2023; Gibbs and Candès,
2023; Bhatnagar et al., 2023b; Zhang et al., 2024) applied regret minimization
algorithms in OCO to select the score threshold in CP, and Bastani et al.
(2022) achieved this task using multicalibration. By relaxing the CP objective
from the _coverage probability_ to the _post-hoc coverage frequency_ , these
methods can handle arbitrary _continual_ distribution shifts. However, they
require the true label to be provided after every prediction, which is a
limiting requirement for many use cases. Our experiments show that it is
possible to achieve comparable performance in these settings without this
limitation, i.e., “label free”.

#### Test Time Adaptation (TTA)

Our EACP method is inspired by core ideas in TTA. The goal of TTA is to update
a trained machine learning model at test time, using unlabeled data from
shifted distributions. To achieve this, one could update the batch-norm
statistics on the test data (Nado et al., 2020; Schneider et al., 2020;
Khurana et al., 2021), or minimize the test-time prediction entropy (Wang et
al., 2021; Zhang et al., 2022; Niu et al., 2022; Song et al., 2023; Press et
al., 2024). Notably, these methods can be applied on any probabilistic and
differentiable model (such as modern neural networks), which is naturally
congruent with the strength of CP. However, to date this line of work has been
disconnected from the conformal prediction literature, and to our knowledge
there are no studies investigating the benefit of TTA to the distribution-
shifted CP setting.

##  3 Preliminaries of CP

We begin by introducing the standard background of CP. For clarity, we assume
i.i.d. data in our exposition, rather than the slightly weaker notion of
exchangeability. Also see (Roth, 2022; Angelopoulos and Bates, 2023;
Tibshirani, 2023).

Let $\mathcal{D}$ be an unknown distribution on the space
$\mathcal{X}\times\mathcal{Y}$ of covariate-label pairs, and let
$\alpha\in(0,1)$ be the _error rate_ we aim for. Given a dataset $D$
consisting of $n$ i.i.d. samples $(x_{i},y_{i})_{i\in[n]}\sim\mathcal{D}^{n}$,
the goal of CP is to generate a set-valued function
$\mathcal{C}_{D}:\mathcal{X}\rightarrow 2^{\mathcal{Y}}$, such that

$\mathbb{P}_{(x,y)\sim\mathcal{D},D\sim\mathcal{D}^{n}}\left[y\in\mathcal{C}_{D}(x)\right]\geq
1-\alpha.$ (1)

That is, for a fresh sample $(x,y)\sim\mathcal{D}$, our prediction set
$\mathcal{C}_{D}(x)$ covers the label $y$ with guaranteed high
probability.333In the literature this is called a _marginal_ (w.r.t. the
covariate $x$), _in-expectation_ (w.r.t. the dataset $D$) guarantee – see
(Roth, 2022, Chapter 7) for various ways to strengthen it. Notice that Eq.(1)
alone is a trivial objective, since it suffices to predict the entire label
space $\mathcal{C}_{D}(x)=\mathcal{Y}$ for all $x$. Therefore, CP is
essentially a bi-objective problem: as long as Eq.(1) is satisfied, we want
the prediction set $\mathcal{C}_{D}(x)$ to be small.

The main difficulty of this set-valued prediction problem is that the range of
output $2^{\mathcal{Y}}$ is too large. In this regard, the key idea of CP is
reducing the problem to 1D prediction via a trained machine learning model,
such as a neural network. Specifically, we assume access to a (positively
oriented; i.e., larger is better) _score function_
$s:\mathcal{X}\times\mathcal{Y}\rightarrow\mathbb{R}_{+}$ given by the base
model, such that for each covariate $x\in\mathcal{X}$ and _candidate label_
$y^{\prime}\in\mathcal{Y}$, $s(x,y^{\prime})$ measures how likely the model
believes that $y^{\prime}$ is the true label $y$. Then, all there is left for
CP is to pick a threshold $\tau_{D}\in\mathbb{R}$ that depends on the dataset
$D$, and predict the label set (if the score function is negatively oriented,
then $\geq$ is replaced by $\leq$)

$\mathcal{C}_{D}\left(x;\tau_{D}\right):=\left\\{y^{\prime}\in\mathcal{Y}:s\left(x,y^{\prime}\right)\geq\tau_{D}\right\\}.$
(2)

Under the i.i.d. assumption, the coverage objective Eq.(1) is satisfied by
picking $\tau_{D}$ as the $\alpha(1+n^{-1})$-quantile of the _empirical
scores_ $\\{s(x_{i},y_{i})\\}_{i\in[n]}$. Since the training data of the base
model is split from the calibration dataset used to determine $\tau_{D}$, this
approach is commonly known as _split conformal prediction_ , which we also
refer to as SplitCP. Notably, $\tau_{D}$ is determined by the calibration
dataset; once the latter is fixed, there is no need to access the ground truth
labels at test time.

#### Example

In classification (which we focus on), the standard practice is to define the
score function $s(x,y)=\pi_{\theta}(x)_{y}$ (Sadinle et al., 2019), where
$\pi_{\theta}$ is a trained neural network parameterized by $\theta$, and
$\pi_{\theta}(x)_{y}\in[0,1]$ is the softmax score corresponding to one of the
$k$-classes $y\in[k]$. Such a score function is positively oriented (i.e.,
larger means more likely to be included in the prediction set).

#### Distribution shift

For the rest of this paper, we study a deviation of the above problem where
the test sample $(x,y)$ is _not_ drawn from $\mathcal{D}$, but rather some new
unknown distribution $\widetilde{\mathcal{D}}_{\mathrm{test}}$. In general, it
is impossible to guarantee Eq.(1) without restrictions on the shift, but we
will show that with more help from the base model, the above CP procedure can
be modified to work well in practice.

##  4 Our methods

In this section, we first propose a method called ECP (Entropy scaled
Conformal Prediction), which improves the coverage rate of CP by enlarging its
prediction sets using the uncertainty evaluation of the base model itself.
Crucially, this notion of uncertainty can be directly minimized and refined
through unsupervised TTA, leading to an improved method called EACP (Entropy
base-Adapted Conformal Prediction). It is able to both i) recover the desired
error rate on many challenging distribution-shifted datasets, while ii)
significantly reduce inflated set sizes under increased uncertainty.

###  4.1 Scaling conformal scores by uncertainty

Let us start by motivating our idea. A common criticism of standard SplitCP
(although typically studied in the regression setting rather than
classification) is that the size of the CP prediction set is invariant to the
covariate $x$. To provide greater “local adaptivity” (w.r.t. $x$),
uncertainty-aware techniques have been proposed to adjust (or scale) the score
function $s(x,y)$ based on some notion of uncertainty (or difficulty) that the
base model decides at each $x$ (Papadopoulos et al., 2008; Johansson et al.,
2015; Lei et al., 2018; Izbicki et al., 2019; Romano et al., 2019; Seedat et
al., 2023; Rossellini et al., 2024). This has the effect of inflating the
prediction set on the base model’s uncertain regions, and has been shown to
improve the _conditional coverage rate_ of CP (Angelopoulos and Bates, 2023;
Tibshirani, 2023).

Figure 1: Entropy vs. the softmax score of true label, averaged on each
dataset.

#### Key idea

Building on these results, our key idea is to apply an analogous uncertainty
scaling on the score function, in order to improve the performance of CP under
distribution shifts. However, instead of using the uncertainty of the base
model at each covariate $x$, we crucially evaluate its uncertainty on the
whole distribution-shifted _test dataset_ , effectively aggregating its
“localized” uncertainty at each $x$. In other words, instead of aiming for
“local adaptivity” (w.r.t. $x$) as in prior works, we aim for the adaptivity
w.r.t. the unknown distribution shift.

#### Prediction entropy

More concretely, which uncertainty measure do we use on the base model? As
discussed above, the ideal dataset-specific uncertainty measure would follow
from a “localized” uncertainty measure at each $x$, and a particularly useful
one for the latter is the _entropy_ of the model’s predictions,

$h(x)=-\sum_{y\in[k]}\pi_{\theta}(x)_{y}\text{log}\pi_{\theta}(x)_{y}.$

Previous works have established the relation between such an entropy notion
(of the base model) and the magnitude of the distribution shift, showing that
larger shifts are strongly correlated with higher entropy (thus uncertainty in
the base model) (Wang et al., 2021). We provide a consistent but unique
observation in Figure 1, which plots the relation between the entropy
(averaged over all $x$ values in the dataset) and the softmax score of the
true label (also averaged over $x$), across a range of different
datasets.444Darker shades represent greater severity levels of ImageNet-C
corruptions; see Section 5 for more details on the datasets. For the true
label to be included in the CP prediction set, which is eventually what we aim
for, its softmax value should be greater than the CP threshold $\tau_{D}$.
From Figure 1, one could see that an increase in entropy is associated with a
decrease in the softmax score of the true label, which crucially means that
the prediction sets should be larger when the entropy is higher.

Figure 2: The targeted coverage rate $1-\alpha$ vs. the empirical coverage
rate, induced by ECP with different $\beta$ values. It shows that simply
setting $\beta=1-\alpha$ in ECP consistently works well for all but the most
severe distribution shifts (e.g., ImageNet-R and ImageNet-C Contrast Severity
5). Such an observation also holds across various $\alpha$ values, suggesting
the general effectiveness of this hyperparameter choice.

Now consider going from the “localized” uncertainty measure $h(x)$ to an
uncertainty measure on the test dataset, denoted as $u_{D_{\mathrm{test}}}$.
We introduce another design choice: defining $u_{D_{\mathrm{test}}}$ as a
specified, $\beta$-quantile of $h(x)$ on the test dataset. Surprisingly, we
find that simply setting $\beta$ to the desired coverage rate $1-\alpha$ is a
fairly reliable choice in practice (see Figure 2), which gives a robust
(over)-estimate of typical $h(x)$ values on the test dataset. We perform all
the experiments with this direct relationship to avoid excessive
hyperparameter tuning, but it can be further refined if desired.

#### Method: ECP

Now we are ready to use $u_{D_{\mathrm{test}}}$ above to scale the score
functions on the test dataset, without label access. The resulting method is
named as ECP (Entropy scaled Conformal Prediction).

Formally, let us denote $q_{\beta}(\cdot)$ as the $\beta$-th quantile of its
argument, and
$h_{D_{\mathrm{test}}}=\\{h(x_{1}),...,h(x)_{N_{\mathrm{test}}}\\}$ as the
entropy values on the test dataset. Then, we compute the “entropy quantile”

$u_{D_{\mathrm{test}}}=q_{1-\alpha}\left(h_{D_{\mathrm{test}}}\right),$ (3)

and on any test sample, modified from Eq.(2), we scale the score function by
$\max(1,u_{D_{\mathrm{test}}})$ to form the CP prediction set

$\mathcal{C}_{D_{\mathrm{test}}}\left(x;\tau_{D}\right):=\left\\{y^{\prime}\in[k]:s(x,y^{\prime})\cdot\max(1,u_{D_{\mathrm{test}}})\geq\tau_{D}\right\\}.$
(4)

Here, we take a maximum with $1$ to ensure that the prediction sets of ECP
cannot be smaller than those of standard SplitCP.

The intuition is that a larger distribution shift will result in larger
entropy, and thus lead to a correspondingly larger up-scaling of the score
function. In this way, more labels have score larger than the fixed CP
threshold $\tau_{D}$, and the prediction set grows. This can help mitigate the
under-coverage issue of standard SplitCP under distribution shifts (without
accessing test labels), and further details are provided in our experiments
(Section 5).

###  4.2 Optimizing uncertainty using TTA

While our ECP method already improves the coverage of SplitCP on several
datasets, it inevitably leads to larger set sizes and, like all post-hoc CP
methods, still relies on a fixed base model. To remedy this, we refine ECP
using techniques from unsupervised TTA, which are able to update the base
model’s parameters on the unlabeled test dataset. While these techniques have
been often investigated in the context of top-1 accuracy, we take a different
perspective and study their ability to improve set-valued classifiers like
conformal prediction.

#### Key idea

Concretely, consider the following loss function (of the base model’s
parameter $\theta$),

$\mathcal{L}(x;\theta):=h(x)=-\sum_{y\in[k]}\pi_{\theta}(x)_{y}\text{log}\pi_{\theta}(x)_{y}.$
(5)

Our insight is that by directly minimizing this loss (through updating
$\theta$), we implicitly minimize the entropy quantile $u_{D_{\mathrm{test}}}$
used to scale the score functions. Consequently, the sizes of ECP prediction
sets can be reduced, and the coverage rate can be increased.

#### Method: EACP

A number of specific TTA methods have been developed to minimize Eq.(5), while
ensuring certain notions of stability. In this work, we leverage the recent
ETA (Efficient Test-time Adaptation) method (Niu et al., 2022) as it is a
simple approach shown to be effective even in continuously shifting
distributions (Press et al., 2023). Adopting this with ECP results in a new CP
method, which we call EACP (Entropy base-Adapted Conformal Prediction).

In practice, one could simply call ETA as a subroutine, so we only present its
high level idea for the sake of completeness. On a batch of new test data
(i.e., $x_{i}$ with a collection of indices $i$), ETA filters the base model’s
outputs (i.e., softmax scores) $s(x_{i},\cdot)$ by excluding the outputs
similar to those already seen. Then, it reweights the remaining outputs based
the associated entropy $h(x_{i})$, with lower entropy (less uncertain) indices
receiving higher weights. This leads to a weighted variant of the loss
function Eq.(5), which is then minimized by performing a single gradient
update. Subsequently, the updated base model is applied to ECP to form the
prediction sets of EACP, according to Eq.(4).

In Section 5, we demonstrate that EACP can further improve the performance of
ECP, by increasing the coverage rate while maintaining informative set sizes.

###  4.3 Uncertainty scaling function

Previously, in Eq.(4) we scale the score functions linearly by the entropy
quantile $u_{D_{\mathrm{test}}}$ of the base model. However, in principle this
can be adjusted by any (potentially non-linear) function $f(\cdot)$ of the
entropy quantile, $f(u_{D_{\mathrm{test}}})$. The best choice of $f(\cdot)$
depends on the unknown relation between $u_{D_{\mathrm{test}}}$ and the
$(1-\alpha)$-quantile of the true labels’ conformal scores, denoted as

$\tau_{D_{\mathrm{test}}}:=q_{1-\alpha}\left[s(x,y);(x,y)\in
D_{\mathrm{test}}\right].$

Specifically, we aim to find an $f(\cdot)$ satisfying
$f(u_{D_{\mathrm{test}}})=\tau_{D}/\tau_{D_{\mathrm{test}}}$.

Figure 3: $u_{D_{\mathrm{test}}}$ versus $\tau_{D}/\tau_{D_{\mathrm{test}}}$
on a log-log scale. On mild and moderate distribution shifts, the linear fit
on the log-log plot has slope between 1 and 2. This suggests the effectiveness
of using a linear or quadratic function as $f(\cdot)$, which acts on the
entropy quantile. However, we also observe that a higher-order polynomial is
required on more difficult shifts, such as ImageNet-R.

While finding this optimal $f(\cdot)$ is obviously infeasible without
observing the labels in $D_{\mathrm{test}}$, in Figure 3 we empirically
evaluate the ideal choice in a post-hoc manner, across different datasets.
Recall that both $u_{D_{\mathrm{test}}}$ and $\tau_{D_{\mathrm{test}}}$ depend
on the desired error rate $\alpha$. Therefore, for each dataset, we vary
$\alpha$ and plot the resulting $u_{D_{\mathrm{test}}}$ versus
$\tau_{D}/\tau_{D_{\mathrm{test}}}$ on a log-log scale. Furthermore, if we
mildly restrict $f(\cdot)$ to polynomials, then its optimal order can be
approximated by the slope of a linear fit on the log-log plot. We only use
not-too-small $\alpha$ values for the linear fit, since it is closer to the
practice of CP and less prone to noise.

Figure 3 shows that the optimal polynomial order generally increases with the
severity of the distribution shift. A larger polynomial order would lead to
larger prediction sets. While end-users can refine $f(\cdot)$ based on a
preference towards ensuring coverage or small set sizes, we empirically
validate that both linear ($\texttt{ECP}_{1}$ / $\texttt{E{\small A}CP}_{1}$)
and quadratic ($\texttt{ECP}_{2}$ / $\texttt{E{\small A}CP}_{2}$) scaling are
sufficient for coverage and moderate set sizes, on a diverse range of
distribution shifts.

##  5 Experiments

We conduct experiments across a number of large-scale datasets and neural
network architectures. Our setup builds on the standard SplitCP procedure
introduced in Section 3, which relies on a held-out, in-distribution,
“development set” for calibrating the CP threshold. On ImageNet variants, we
split the original ImageNet development set, into a calibration and validation
set, resulting in 25,000 calibration points. we refer to Appendix A.1 for more
details on the calibration splits and sizes for all datasets.

The conformal threshold is found on the calibration set, and used for all
subsequent distribution-shifted datasets. _After the conformal threshold is
estimated in-distribution, all subsequent steps are unsupervised._ We show
results on both stationary and continuously shifting distributions.

#### Baselines

We compare our proposed ECP and EACP methods to a number of possible
baselines. This includes using a pre-calibrated threshold on new distribution-
shifted data as in the standard SplitCP approach, using this calibrated
threshold while updating the base model using ETA, and an approach we call
NAIVE which forms prediction sets by including classes until their cumulative
probability is greater or equal to the desired $1-\alpha$ coverage level,
without any prior calibration. In all experiments, the target coverage is
0.90. We also analyze linear and quadratic scaling of the uncertainty
adjustment, as described in Section 4.3.

#### Datasets

We investigate a number of ImageNet (Deng et al., 2009) variants including:
ImageNet-V2 (Recht et al., 2019), ImageNet-R (Hendrycks et al., 2021a), and
ImageNet-A (Hendrycks et al., 2021b). We also test our approach on datasets
from the WILDS Benchmark (Koh et al., 2021) which represent in-the-wild
distribution shifts across many real world applications, including iWildCam
(animal trap images), RXRX1 (cellular images), and FMOW (satellite images).

While the previous datasets present a single distribution shift, the
ImageNet-C (Hendrycks and Dietterich, 2019) dataset allows us to investigate
shifts across many types and severities. Specifically, ImageNet-C applies 19
visual corruptions to the ImageNet validation set across four corruption
categories — noise, blur, weather, and digital, with five severity levels for
each corruption. Please refer to Appendix A.1 for more details on the studied
datasets.

###  5.1 Stationary shifts

Table 1: ECP and EACP can achieve the desired coverage on a number of
distribution-shifted datasets, across a variety of imaging domains
(ecological, cellular, satellite, etc.). The coverage gap is significantly
reduced even on extremely challenging datasets. All results are from ResNet-50
models except FMOW, which uses a DenseNet-121 (Huang et al., 2016). Quadratic
uncertainty scaling provides better coverage rates, however, linear scaling
results in smaller set sizes.

| Coverage | Set Size
---|---|---
Dataset | SplitCP | NAIVE | ETA | $\texttt{ECP}_{1}$ | $\texttt{ECP}_{2}$ | $\texttt{E{\small A}CP}_{1}$ | $\texttt{E{\small A}CP}_{2}$ | SplitCP | NAIVE | ETA | $\texttt{ECP}_{1}$ | $\texttt{ECP}_{2}$ | $\texttt{E{\small A}CP}_{1}$ | $\texttt{E{\small A}CP}_{2}$
ImageNet-V2 | 0.81 | 0.88 | 0.81 | 0.86 | 0.91 | 0.86 | 0.91 | 2.5 | 11.7 | 2.5 | 4.2 | 7.6 | 4.5 | 8.7
ImageNet-R | 0.50 | 0.69 | 0.62 | 0.61 | 0.72 | 0.71 | 0.80 | 3.4 | 20.9 | 3.0 | 9.1 | 23.3 | 6.8 | 16.1
ImageNet-A | 0.03 | 0.14 | 0.05 | 0.10 | 0.27 | 0.14 | 0.30 | 3.4 | 12.7 | 3.6 | 7.4 | 15.1 | 8.7 | 10.1
iWildCam | 0.84 | 0.76 | 0.84 | 0.84 | 0.88 | 0.84 | 0.89 | 3.9 | 2.5 | 3.8 | 3.8 | 5.5 | 3.7 | 5.6
RXRX1 | 0.84 | 0.48 | 0.87 | 0.87 | 0.90 | 0.90 | 0.93 | 81.8 | 6.4 | 100 | 105 | 137 | 133 | 177
FMOW | 0.87 | 0.83 | 0.87 | 0.93 | 0.96 | 0.93 | 0.94 | 6.2 | 5.8 | 6.5 | 10.3 | 15.3 | 11.1 | 16.4

Table 1 summarizes our performance on various natural distribution-shifted
datasets. We observe that ETA alone is unable to fully recover the desired
coverage, however, ECP closes the gap significantly. Coverage is further
improved by adaptation through EACP, which helps to reduce set sizes on some
datasets – it recovers the desired coverage on all but the most challenging
datasets. Here, we also see the trade-off between linear and quadratic
uncertainty scaling. We observe $\texttt{E{\small A}CP}_{2}$ consistently
achieves higher coverage rates, however this also leads to “over-coverage” on
some datasets and thus larger sets. In contrast, $\texttt{E{\small A}CP}_{1}$
leads to lower coverage but also smaller set sizes. This trade-off can be
selected by end-users depending on their preference for more accurate or more
efficient prediction sets. In subsequent experiments, we will focus on
demonstrating the benefit of $\texttt{E{\small A}CP}_{2}$ on coverage, while
noting that the observed set sizes are nonetheless practically useful and far
from trivial.

Table 2: Coverage on four different corruption types representing each
ImageNet-C category. By explicitly leveraging uncertainty under distribution
shift, our $\texttt{ECP}_{2}$ method closes the coverage gap on most severity
levels. $\texttt{E{\small A}CP}_{2}$ improves on this by refining the
uncertainty at test-time, resulting in a recovered coverage rate of 0.90 on
nearly all corruption types and severities.

Method | Contrast | Brightness | Gaussian Noise | Motion Blur
---|---|---|---|---
1 | 2 | 3 | 4 | 5 | 1 | 2 | 3 | 4 | 5 | 1 | 2 | 3 | 4 | 5 | 1 | 2 | 3 | 4 | 5
NAIVE | 0.91 | 0.89 | 0.87 | 0.83 | 0.76 | 0.92 | 0.92 | 0.91 | 0.91 | 0.90 | 0.88 | 0.85 | 0.79 | 0.69 | 0.79 | 0.91 | 0.90 | 0.85 | 0.77 | 0.71
SplitCP (Sadinle et al., 2019) | 0.83 | 0.78 | 0.66 | 0.36 | 0.09 | 0.88 | 0.87 | 0.86 | 0.83 | 0.78 | 0.79 | 0.69 | 0.50 | 0.26 | 0.07 | 0.83 | 0.74 | 0.57 | 0.37 | 0.27
ETA (Niu et al., 2022) | 0.87 | 0.86 | 0.84 | 0.79 | 0.63 | 0.88 | 0.88 | 0.87 | 0.86 | 0.84 | 0.86 | 0.82 | 0.76 | 0.69 | 0.54 | 0.86 | 0.84 | 0.80 | 0.73 | 0.68
$\texttt{ECP}_{2}$ (ours) | 0.93 | 0.92 | 0.89 | 0.79 | 0.60 | 0.94 | 0.94 | 0.94 | 0.93 | 0.92 | 0.92 | 0.88 | 0.80 | 0.86 | 0.38 | 0.94 | 0.92 | 0.86 | 0.75 | 0.68
$\texttt{E{\small A}CP}_{2}$ (ours) | 0.93 | 0.93 | 0.93 | 0.92 | 0.87 | 0.93 | 0.93 | 0.93 | 0.93 | 0.93 | 0.93 | 0.93 | 0.92 | 0.90 | 0.84 | 0.93 | 0.92 | 0.92 | 0.91 | 0.89

In Table 2, we show fine-grained results of our method’s performance on one
corruption type for each ImageNet-C category and across each severity level.
Here, we can see the benefit of leveraging an uncertainty notion that can be
directly minimized and refined on new test samples. By both adjusting the
conformal scores by an entropy-based notion of uncertainty while
simultaneously refining this uncertainty, EACP is able to recover the target
coverage rate on most corruption types and severities.

Figure 4: Our $\texttt{E{\small A}CP}_{2}$ method is able to improve coverage
using various neural network models and architectures, as well as a large and
diverse range of distribution shifts. It consistently “hugs” the desired
coverage rate, while maintaining practical set sizes. Results are averaged
across five severity levels for each corruption type in the ImageNet-C
dataset. We zoom in on the right to clearly see the benefit of adapting at
test time on coverage and set sizes.

We report performance using a number of neural networks architectures and
sizes on all 19 corruption types averaged across severities for each
corruption in Figure 4. Here, the methodology and hyper-parameters are exactly
the same for all models, including the vision transformers. We can see that
even with no hyper-parameter tuning (including on the TTA algorithm), our
method considerably reduces the coverage gap across a variety of neural
networks. An immediate observation is that while ETA does not recover
coverage, it also consistently leads to smaller set sizes. While minimizing
set sizes is indeed desirable, this may be of secondary importance to the
primary goal of achieving the desired coverage level. Further, we note that
the set sizes without ECP are in fact consistent and with little variance,
regardless of the empirical coverage. We argue that this is actually an
undesirable quality. As an uncertainty quantification framework, conformal
prediction is considered to convey uncertainty through its set sizes; greater
uncertainty is supposed to be reflected in larger sets, and vice versa. Our
results suggest that explicitly incorporating uncertainty into the conformal
predictor seems necessary to accurately convey the underlying uncertainty
under distribution shift. While even large mis-coverage still results in small
set sizes using ETA, larger sets under EACP consistently correspond to lower
coverage — an important and practical quality for end-users.

###  5.2 Continuous shifts

Table 3: We evaluate performance on the challenging setting of continuously
shifting distributions. The “label free” column denotes whether a method
relies on labels at test-time from the target data. We recall that SplitCP
does not adapt to new data. In addition to the average coverage ($\uparrow$)
and average size ($\downarrow$), we also measure the worst local corruption
error LCE128 ($\downarrow$) and worst local set size LSS128, ($\downarrow$) on
a sliding window of 128 test points.

|  | Gradual shift | Sudden shift
---|---|---|---
Label Free | Method | Avg. Cov | Avg. Size |  LCE128 |  LSS128 | Avg. Cov | Avg. Size |  LCE128 |  LSS128
- |  SplitCP (Sadinle et al., 2019) | 0.59 | 3.1 | 0.70 | 3.6 | 0.59 | 2.8 | 0.71 | 3.5
✗ |  SAOCP (Seedat et al., 2023) | 0.79 | 145 | 0.24 | 353 | 0.78 | 139 | 0.28 | 349
✗ |  FACI (Gibbs and Candes, 2021) | 0.90 | 101 | 0.07 | 455 | 0.90 | 142 | 0.09 | 450
✗ |  MAGL-D (Zhang et al., 2024) | 0.90 | 403 | 0.05 | 856 | 0.90 | 355 | 0.05 | 844
✗ |  MAGL (Zhang et al., 2024) | 0.90 | 117 | 0.06 | 573 | 0.90 | 168 | 0.3 | 704
✗ |  MAGDIS (Zhang et al., 2024) | 0.90 | 417 | 0.06 | 841 | 0.90 | 372 | 0.07 | 852
✓ |  ETA (Niu et al., 2022) | 0.69 | 2.9 | 0.52 | 3.4 | 0.67 | 2.7 | 0.54 | 3.5
✓ |  $\texttt{ECP}_{2}$ (ours) | 0.84 | 36.6 | 0.35 | 90.4 | 0.82 | 37.5 | 0.38 | 88.5
✓ |  $\texttt{E{\small A}CP}_{2}$ (ours) | 0.88 | 22.4 | 0.20 | 47.8 | 0.86 | 23.1 | 0.28 | 55.7

We further investigate the setting of continuous distribution shift in Table
3. This setting is similar to that studied in previous Online Conformal
Prediction (OCP) works (Zhang et al., 2024; Bhatnagar et al., 2023b; Gibbs and
Candès, 2023; Gibbs and Candes, 2021). This involves shifting between
ImageNet-C severity level 1 to level 5 (either suddenly or gradually), while
sampling random corruptions at each corresponding severity. We emphasize that
this is a particularly challenging setting, as it presents a continuous shift
in both the magnitude as well as type of corruption. See Appendix A.2 for more
details on this set-up. In addition to unsupervised methods, we also compare
with supervised OCP methods that rely on the correct label being revealed
after every prediction. We notice that on gradual shifts in severity,
$\texttt{E{\small A}CP}_{2}$ achieves nearly the same coverage results as
these strong baselines, with greatly reduced set-sizes. We believe this
considerable reduction in set size is due to two factors. First,
$\texttt{E{\small A}CP}_{2}$ is able to update the underlying neural network
on new data. Second, while OCP is able to learn the (potentially expansive)
necessary threshold, $\texttt{E{\small A}CP}_{2}$ instead directly updates the
conformal scores themselves. In addition to average coverage and set size, we
also measure the local coverage error LCE128 across the worst sliding window
of 128 samples, and similarly the worst local set size, LSS128. While the
supervised methods unsurprisingly result in better local coverage, they also
lead to local set sizes that are an order of magnitude larger than those
produced by $\texttt{E{\small A}CP}_{2}$.

##  6 Conclusion

This paper studies how to adapt CP methods to distribution-shifted data,
without relying on labels from the target dataset. This is an important
challenge in many real world settings, where exchangeability assumptions are
violated and labels may be difficult to attain. We propose an uncertainty-
aware method based on the prediction entropy (ECP), and leverage unsupervised
test time adaptation to update the base model and refine its uncertainty
(EACP). We demonstrate that the proposed methods are able to recover the
desired error rate on a wide range of distribution shifts, while maintaining
efficient set sizes. Furthermore, they are even competitive with supervised
approaches on challenging and continuously shifting distributions. We hope
this inspires future works continuing to tackle this important challenge.

## References

  * Angelopoulos et al. (2023) Anastasios Angelopoulos, Emmanuel Candes, and Ryan J Tibshirani.  Conformal pid control for time series prediction.  _Advances in Neural Information Processing Systems_ , 36, 2023.
  * Angelopoulos and Bates (2023) Anastasios N. Angelopoulos and Stephen Bates.  Conformal prediction: A gentle introduction.  _Foundations and Trends® in Machine Learning_ , 16(4):494–591, 2023.  ISSN 1935-8237.  doi: 10.1561/2200000101.  URL http://dx.doi.org/10.1561/2200000101.
  * Angelopoulos et al. (2021) Anastasios Nikolas Angelopoulos, Stephen Bates, Michael Jordan, and Jitendra Malik.  Uncertainty sets for image classifiers using conformal prediction.  In _International Conference on Learning Representations_ , 2021.  URL https://openreview.net/forum?id=eNdiU_DbM9.
  * Barber et al. (2023) Rina Foygel Barber, Emmanuel J Candes, Aaditya Ramdas, and Ryan J Tibshirani.  Conformal prediction beyond exchangeability.  _The Annals of Statistics_ , 51(2):816–845, 2023.
  * Bastani et al. (2022) Osbert Bastani, Varun Gupta, Christopher Jung, Georgy Noarov, Ramya Ramalingam, and Aaron Roth.  Practical adversarial multivalid conformal prediction.  _Advances in Neural Information Processing Systems_ , 35:29362–29373, 2022.
  * Beery et al. (2020) Sara Beery, Elijah Cole, and Arvi Gjoka.  The iwildcam 2020 competition dataset.  _arXiv preprint arXiv:2004.10340_ , 2020.
  * Bhatnagar et al. (2023a) Aadyot Bhatnagar, Huan Wang, Caiming Xiong, and Yu Bai.  Improved online conformal prediction via strongly adaptive online learning.  In _International Conference on Machine Learning_ , pages 2337–2363. PMLR, 2023a.
  * Bhatnagar et al. (2023b) Aadyot Bhatnagar, Huan Wang, Caiming Xiong, and Yu Bai.  Improved online conformal prediction via strongly adaptive online learning.  In _Proceedings of the 40th International Conference on Machine Learning_ , ICML’23. JMLR.org, 2023b.
  * Cauchois et al. (2024) Maxime Cauchois, Suyash Gupta, Alnur Ali, and John C Duchi.  Robust validation: Confident predictions even when distributions shift.  _Journal of the American Statistical Association_ , pages 1–66, 2024.
  * Christie et al. (2018) Gordon Christie, Neil Fendley, James Wilson, and Ryan Mukherjee.  Functional map of the world.  In _Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition_ , 2018.
  * Deng et al. (2009) Jia Deng, Wei Dong, Richard Socher, Li-Jia Li, Kai Li, and Li Fei-Fei.  Imagenet: A large-scale hierarchical image database.  In _2009 IEEE conference on computer vision and pattern recognition_ , pages 248–255. Ieee, 2009.
  * Gao et al. (2024) Yihuai Gao, Yukai Tang, Han Qi, and Heng Yang.  Closure: Fast quantification of pose uncertainty sets.  _arXiv preprint arXiv:2403.09990_ , 2024.
  * Gibbs and Candes (2021) Isaac Gibbs and Emmanuel Candes.  Adaptive conformal inference under distribution shift.  _Advances in Neural Information Processing Systems_ , 34:1660–1672, 2021.
  * Gibbs and Candès (2023) Isaac Gibbs and Emmanuel Candès.  Conformal inference for online prediction with arbitrary distribution shifts, 2023.
  * Grandvalet and Bengio (2004) Yves Grandvalet and Yoshua Bengio.  Semi-supervised learning by entropy minimization.  _Advances in Neural Information Processing Systems_ , 17, 2004.
  * Guo et al. (2017) Chuan Guo, Geoff Pleiss, Yu Sun, and Kilian Q Weinberger.  On calibration of modern neural networks.  In _International conference on machine learning_ , pages 1321–1330. PMLR, 2017.
  * Hendrycks and Dietterich (2019) Dan Hendrycks and Thomas Dietterich.  Benchmarking neural network robustness to common corruptions and perturbations.  _Proceedings of the International Conference on Learning Representations_ , 2019.
  * Hendrycks et al. (2021a) Dan Hendrycks, Steven Basart, Norman Mu, Saurav Kadavath, Frank Wang, Evan Dorundo, Rahul Desai, Tyler Zhu, Samyak Parajuli, Mike Guo, Dawn Song, Jacob Steinhardt, and Justin Gilmer.  The many faces of robustness: A critical analysis of out-of-distribution generalization.  _ICCV_ , 2021a.
  * Hendrycks et al. (2021b) Dan Hendrycks, Kevin Zhao, Steven Basart, Jacob Steinhardt, and Dawn Song.  Natural adversarial examples.  _CVPR_ , 2021b.
  * Huang et al. (2016) Gao Huang, Zhuang Liu, and Kilian Q. Weinberger.  Densely connected convolutional networks.  _2017 IEEE Conference on Computer Vision and Pattern Recognition (CVPR)_ , pages 2261–2269, 2016.  URL https://api.semanticscholar.org/CorpusID:9433631.
  * Izbicki et al. (2019) Rafael Izbicki, Gilson T Shimizu, and Rafael B Stern.  Flexible distribution-free conditional predictive bands using density estimators.  _arXiv preprint arXiv:1910.05575_ , 2019.
  * Johansson et al. (2015) Ulf Johansson, Cecilia Sönströd, and Henrik Linusson.  Efficient conformal regressors using bagged neural nets.  In _2015 International Joint Conference on Neural Networks (IJCNN)_ , pages 1–8, 2015.  doi: 10.1109/IJCNN.2015.7280763.
  * Kasa and Taylor (2023) Kevin Kasa and Graham W. Taylor.  Empirically validating conformal prediction on modern vision architectures under distribution shift and long-tailed data.  In _ICML 2023 Workshop on Structured Probabilistic Inference & Generative Modeling_, 2023.  URL https://openreview.net/forum?id=TMxpy0aluL.
  * Khurana et al. (2021) Ansh Khurana, Sujoy Paul, Piyush Rai, Soma Biswas, and Gaurav Aggarwal.  Sita: Single image test-time adaptation.  _arXiv preprint arXiv:2112.02355_ , 2021.
  * Koh et al. (2021) Pang Wei Koh, Shiori Sagawa, Henrik Marklund, Sang Michael Xie, Marvin Zhang, Akshay Balsubramani, Weihua Hu, Michihiro Yasunaga, Richard Lanas Phillips, Irena Gao, Tony Lee, Etienne David, Ian Stavness, Wei Guo, Berton A. Earnshaw, Imran S. Haque, Sara Beery, Jure Leskovec, Anshul Kundaje, Emma Pierson, Sergey Levine, Chelsea Finn, and Percy Liang.  WILDS: A benchmark of in-the-wild distribution shifts.  In _International Conference on Machine Learning (ICML)_ , 2021.
  * Kumar et al. (2023) Bhawesh Kumar, Charlie Lu, Gauri Gupta, Anil Palepu, David Bellamy, Ramesh Raskar, and Andrew Beam.  Conformal prediction with large language models for multi-choice question answering.  _arXiv preprint arXiv:2305.18404_ , 2023.
  * Lei et al. (2018) Jing Lei, Max G’Sell, Alessandro Rinaldo, Ryan J Tibshirani, and Larry Wasserman.  Distribution-free predictive inference for regression.  _Journal of the American Statistical Association_ , 113(523):1094–1111, 2018.
  * Lekeufack et al. (2023) Jordan Lekeufack, Anastasios A Angelopoulos, Andrea Bajcsy, Michael I Jordan, and Jitendra Malik.  Conformal decision theory: Safe autonomous decisions from imperfect predictions.  _arXiv preprint arXiv:2310.05921_ , 2023.
  * Liang et al. (2023) Percy Liang, Rishi Bommasani, Tony Lee, Dimitris Tsipras, Dilara Soylu, Michihiro Yasunaga, Yian Zhang, Deepak Narayanan, Yuhuai Wu, Ananya Kumar, et al.  Holistic evaluation of language models.  _Transactions on Machine Learning Research_ , 2023.
  * Lindemann et al. (2023) Lars Lindemann, Matthew Cleaveland, Gihyun Shim, and George J Pappas.  Safe planning in dynamic environments using conformal prediction.  _IEEE Robotics and Automation Letters_ , 2023.
  * Minderer et al. (2021) Matthias Minderer, Josip Djolonga, Rob Romijnders, Frances Hubis, Xiaohua Zhai, Neil Houlsby, Dustin Tran, and Mario Lucic.  Revisiting the calibration of modern neural networks.  _Advances in Neural Information Processing Systems_ , 34:15682–15694, 2021.
  * Nado et al. (2020) Zachary Nado, Shreyas Padhy, D Sculley, Alexander D’Amour, Balaji Lakshminarayanan, and Jasper Snoek.  Evaluating prediction-time batch normalization for robustness under covariate shift.  _arXiv preprint arXiv:2006.10963_ , 2020.
  * Niu et al. (2022) Shuaicheng Niu, Jiaxiang Wu, Yifan Zhang, Yaofo Chen, Shijian Zheng, Peilin Zhao, and Mingkui Tan.  Efficient test-time model adaptation without forgetting.  In _The Internetional Conference on Machine Learning_ , 2022.
  * Papadopoulos et al. (2008) Harris Papadopoulos, Alex Gammerman, and Volodya Vovk.  Normalized nonconformity measures for regression conformal prediction.  In _Proceedings of the 26th IASTED International Conference on Artificial Intelligence and Applications_ , AIA ’08, page 64–69, USA, 2008. ACTA Press.  ISBN 9780889867109.
  * Press et al. (2023) Ori Press, Steffen Schneider, Matthias Kümmerer, and Matthias Bethge.  Rdumb: A simple approach that questions our progress in continual test-time adaptation.  _Advances in Neural Information Processing Systems_ , 36, 2023.
  * Press et al. (2024) Ori Press, Ravid Shwartz-Ziv, Yann LeCun, and Matthias Bethge.  The entropy enigma: Success and failure of entropy minimization.  In _International Conference on Machine Learning_. PMLR, 2024.
  * Quach et al. (2024) Victor Quach, Adam Fisch, Tal Schuster, Adam Yala, Jae Ho Sohn, Tommi S. Jaakkola, and Regina Barzilay.  Conformal language modeling.  In _The Twelfth International Conference on Learning Representations_ , 2024.  URL https://openreview.net/forum?id=pzUhfQ74c5.
  * Recht et al. (2019) Benjamin Recht, Rebecca Roelofs, Ludwig Schmidt, and Vaishaal Shankar.  Do imagenet classifiers generalize to imagenet?  In _International Conference on Machine Learning_ , 2019.
  * Ren et al. (2023) Allen Z. Ren, Anushri Dixit, Alexandra Bodrova, Sumeet Singh, Stephen Tu, Noah Brown, Peng Xu, Leila Takayama, Fei Xia, Jake Varley, Zhenjia Xu, Dorsa Sadigh, Andy Zeng, and Anirudha Majumdar.  Robots that ask for help: Uncertainty alignment for large language model planners.  In _Proceedings of the Conference on Robot Learning (CoRL)_ , 2023.
  * Romano et al. (2019) Yaniv Romano, Evan Patterson, and Emmanuel Candes.  Conformalized quantile regression.  _Advances in Neural Information Processing Systems_ , 32, 2019.
  * Rossellini et al. (2024) Raphael Rossellini, Rina Foygel Barber, and Rebecca Willett.  Integrating uncertainty awareness into conformalized quantile regression.  In Sanjoy Dasgupta, Stephan Mandt, and Yingzhen Li, editors, _Proceedings of The 27th International Conference on Artificial Intelligence and Statistics_ , volume 238 of _Proceedings of Machine Learning Research_ , pages 1540–1548. PMLR, 02–04 May 2024.  URL https://proceedings.mlr.press/v238/rossellini24a.html.
  * Roth (2022) Aaron Roth.  Uncertain: Modern topics in uncertainty estimation.  _Unpublished Lecture Notes_ , 2022.  URL https://www.cis.upenn.edu/~aaroth/uncertainty-notes.pdf.
  * Sadinle et al. (2019) Mauricio Sadinle, Jing Lei, and Larry Wasserman.  Least Ambiguous Set-Valued Classifiers with Bounded Error Levels.  _Journal of the American Statistical Association_ , 114(525):223–234, January 2019.  ISSN 0162-1459, 1537-274X.  doi: 10.1080/01621459.2017.1395341.  URL http://arxiv.org/abs/1609.00451.  arXiv:1609.00451 [cs, stat].
  * Schneider et al. (2020) Steffen Schneider, Evgenia Rusak, Luisa Eck, Oliver Bringmann, Wieland Brendel, and Matthias Bethge.  Improving robustness against common corruptions by covariate shift adaptation.  _Advances in Neural Information Processing Systems_ , 33:11539–11551, 2020.
  * Seedat et al. (2023) Nabeel Seedat, Alan Jeffares, Fergus Imrie, and Mihaela van der Schaar.  Improving adaptive conformal prediction using self-supervised learning.  In _International Conference on Artificial Intelligence and Statistics_ , pages 10160–10177. PMLR, 2023.
  * Song et al. (2023) Junha Song, Jungsoo Lee, In So Kweon, and Sungha Choi.  Ecotta: Memory-efficient continual test-time adaptation via self-distilled regularization.  In _Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition_ , pages 11920–11929, 2023.
  * Taylor et al. (2019) J. Taylor, B. Earnshaw, B. Mabey, M. Victors, and J. Yosinski.  Rxrx1: An image set for cellular morphological variation across many experimental batches.  In _International Conference on Learning Representations (ICLR)_ , 2019.
  * Tibshirani (2023) Ryan Tibshirani.  Advanced topics in statistical learning: Conformal prediction.  https://www.stat.berkeley.edu/~ryantibs/statlearn-s23/lectures/conformal.pdf, 2023.
  * Tibshirani et al. (2019) Ryan J Tibshirani, Rina Foygel Barber, Emmanuel Candes, and Aaditya Ramdas.  Conformal prediction under covariate shift.  _Advances in Neural Information Processing Systems_ , 32, 2019.
  * Vovk et al. (2005) Vladimir Vovk, Alexander Gammerman, and Glenn Shafer.  _Algorithmic learning in a random world_ , volume 29.  Springer, 2005.
  * Wang et al. (2021) Dequan Wang, Evan Shelhamer, Shaoteng Liu, Bruno Olshausen, and Trevor Darrell.  Tent: Fully test-time adaptation by entropy minimization.  In _International Conference on Learning Representations_ , 2021.  URL https://openreview.net/forum?id=uXl3bZLkr3c.
  * Yang and Pavone (2023) Heng Yang and Marco Pavone.  Object pose estimation with statistical guarantees: Conformal keypoint detection and geometric uncertainty propagation.  In _Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition_ , pages 8947–8958, 2023.
  * Zhang et al. (2022) Marvin Zhang, Sergey Levine, and Chelsea Finn.  Memo: Test time robustness via adaptation and augmentation.  _Advances in Neural Information Processing Systems_ , 35:38629–38642, 2022.
  * Zhang et al. (2024) Zhiyu Zhang, David Bombara, and Heng Yang.  Discounted adaptive online prediction.  _International Conference on Machine Learning_ , 2024.
  * Zinkevich (2003) Martin Zinkevich.  Online convex programming and generalized infinitesimal gradient ascent.  In _Machine Learning, Proceedings of the Twentieth International Conference (ICML 2003), August 21-24, 2003, Washington, DC, USA_ , pages 928–936. AAAI Press, 2003.  URL http://www.aaai.org/Library/ICML/2003/icml03-120.php.

## Appendix

Appendix A contains detailed information on our studied datasets and
experimental protocols, including TTA hyper-parameters, CP procedure, and
continuous shift set-up. Appendix B discusses other possible uncertainty
measures and their deficiencies. Appendix C presents results with another base
TTA method, Appendix D has results of more neural network architectures on
natural distribution shifts, Appendix E has graphs performance on continuous-
shifts at each time-step, and Appendix F has full ImageNet-C results on each
severity level.

##  Appendix A Experimental details

###  A.1 Dataset details

We perform experiments on a number of large-scale datasets that are frequently
used to evaluate deep learning performance under distribution shift [Koh et
al., 2021, Wang et al., 2021, Minderer et al., 2021, Niu et al., 2022, Zhang
et al., 2022, Bhatnagar et al., 2023a, Zhang et al., 2024]:

  * •

ImageNet-V2 [Recht et al., 2019] is a new ImageNet test-set that contains
10,000 images that were collected by closely following the original ImageNet
data collection process.

  * •

ImageNet-R [Hendrycks et al., 2021a] includes renditions (e.g., paintings,
sculptures, drawings, etc.) of 200 ImageNet classes, resulting in a test set
of 30,000 images.

  * •

ImageNet-A [Hendrycks et al., 2021b] consists of 7,500 real-world, unmodified,
and naturally occurring adversarial images which a ResNet-50 model failed to
correctly classify.

  * •

ImageNet-C [Hendrycks and Dietterich, 2019] applies 19 visual corruptions
across four categories and at five severity levels to the original ImageNet
validation set.

  * •

iWildCam [Koh et al., 2021, Beery et al., 2020] contains camera-trap images
from different areas of the world, representing geographic distribution-shift.
It includes a validation set of 7,314 images from the same camera traps the
model was trained on, which is used as our calibration data, as well as 42,791
images from different camera traps that is used as our test set. The images
contain one of the 182 possible animal species.

  * •

RXRX1 [Koh et al., 2021, Taylor et al., 2019] consists of high resolution
fluorescent microscopy images of human cells which have been given one of
1,139 genetic treatments, with the goal of generalizing across experimental
batches. It is split into a 40,612 in-distribution validation set and 34,432
test set.

  * •

FMOW [Koh et al., 2021, Christie et al., 2018] is a satellite imaging dataset
with the goal of classifying images into one of 62 different land use or
building types. It consists of 11,483 validation images from the years from
2002–2013, and 22,108 test images from the years from 2016–2018.

###  A.2 Experimental protocols

#### Conformal prediction

Our split conformal prediction set-up follow previous works [Angelopoulos et
al., 2021, Angelopoulos and Bates, 2023], which divides a held-out dataset
into a calibration and test set. On ImageNet variants, we split the original
validation set in half to produce 25,000 calibration points and 25,000 in-
distribution test points. The calibrated scores and / or threshold are then
used for subsequent distribution-shifted data. On the WILDS datasets, we
similarly split the in-distribution validation sets.

#### Adaptation procedure

Our ImageNet-based experiments are conducted on pre-trained ResNets provided
by the torchvision library555https://github.com/pytorch/vision, and ViTs
provided by the timm library 666https://github.com/huggingface/pytorch-image-
models. Experiments on WILDS datasets are conducted using pre-trained models
provided by the authors of that study 777https://github.com/p-lambda/wilds.
For EACP and ETA, we closely follow the optimization hyperparameters from the
original paper [Niu et al., 2022]: we use SGD optimizer with a momentum of 0.9
and learning rate of 0.00025. We use a batch size of 64 for all ImageNet
experiments, 128 for RXRX1 and FMOW, and 42 for iWildCam. Our experiments are
conducted using a single NVIDIA A40 GPU.

#### Continuous shift

We adopt a slightly modified version of the experimental design for continuous
distribution shift presented in previous works [Bhatnagar et al., 2023a, Zhang
et al., 2024]. This involves sampling random corruptions from the ImageNet-C
dataset under two regimes: gradual shifts where the severity level first
increases in order from $\\{1,...,5\\}$ then decreases from $\\{5,...,1\\}$,
and sudden shifts where the severity level alternates between 1 and 5. In
addition to sampling random corruptions, we also consider in Figure 7 results
on the “easier” setting of shifting severity’s on a single corruption type.

##  Appendix B What is the right measure of uncertainty?

(a) Softmax variance

(b) 1 - Max softmax

Figure 5: Similarly to Figure 1, we present the relation between different
uncertainty measures and the average score of the true label. We see that
softmax variance (left) has an inverse relation with distribution shift, and
$1-$ maximum softmax is a bounded metric that may provide an insufficient
adjustment.

Although in Section 4.3 we propose adjusting the conformal scores by the
prediction entropy of the base model, it is worth asking if there exist other
notions of uncertainty that may instead be used. Here, we consider two
additional uncertainty measures and their relation with the softmax value of
the true label (which is ultimately what we would like to include in our
prediction set), and show they are ill-suited for our task. Firstly, in Figure
5(a) we consider the variance of the softmax scores. Perhaps surprisingly, we
see that distribution shift most often leads to a _smaller_ variance, thus
conveying that the base model is _less_ uncertain. This suggests that softmax
variance is a deficient uncertainty measure as it fails to capture the actual
underlying uncertainty on distribution-shifted data.

We also consider $1-$ maximum softmax score as another possible uncertainty
measure, and see in Figure 5(b) that distribution shift is associated with
smaller maximum softmax values. Unlike softmax variance, this does appear to
better capture the uncertainty, as we would expect the base model to be less
confident on distribution-shifted data. However, this uncertainty measure can
still only take a maximum value of one and thus may not provide necessary
adjustment magnitude, and it is unknown if it can reliably update the base
model label-free.

While there may exist better uncertainty measures that future works can
explore, these results suggest that the prediction entropy is a simple and
reliable measure for conformal adjustments that can effectively capture the
underlying uncertainty.

##  Appendix C Other TTA methods

We investigate our methods performance with another base TTA method in Table
4. Here, we use the Tent update [Wang et al., 2021], which is a simpler
version of ETA with no re-weighing of the entropy loss. While our proposed
methods are also compatible with Tent, we notice that the more powerful ETA
leads to better coverage and set sizes as seen in Table 1. We can expect that
additional improvements in TTA will similarly lead to improvements in our EACP
method.

Table 4: Our proposed EACP performs well with other TTA methods, as seen here
using Tent Wang et al. [2021] as the TTA update.

Dataset | SplitCP | Tent | $\texttt{ECP}_{2}$ | $\texttt{E{\small A}CP}_{2}$
---|---|---|---|---
(coverage / set size) | (coverage / set size) | (coverage / set size) | (coverage / set size)
ImageNet-V2 | 0.81 / 2.5 | 0.81 / 2.6 | 0.91 / 8.0 | 0.92 / 9.6
ImageNet-R | 0.50 / 3.2 | 0.58 / 3.3 | 0.73 / 23 | 0.77 / 17
ImageNet-A | 0.07 / 1.5 | 0.21 / 3.1 | 0.58 / 204 | 0.40 / 24
iWildCam | 0.83 / 3.5 | 0.81 / 2.6 | 0.89 / 5.7 | 0.85 / 3.4
RXRX1 | 0.85 / 83 | 0.87 / 101 | 0.90 / 136.7 | 0.92 / 176
FMOW | 0.87 / 6.3 | 0.85 / 5.7 | 0.96 / 15.6 | 0.94 / 13.4

##  Appendix D More architecture comparisons

In Table 5, we further demonstrate our methods improvements to coverage loss
on natural distribution shifts using diverse neural network architectures. As
expected, larger and more accurate neural networks result in better coverage
and smaller set sizes using ECP and EACP Ṫhis is encouraging as it
demonstrates our methods can scale along with the underlying model.

Table 5: On natural distribution shifts, the performance of our methods scale
well with the performance of the base classifier. This is encouraging as it
suggests compatibility

Dataset | Model | SplitCP | ETA | $\texttt{ECP}_{2}$ | $\texttt{E{\small A}CP}_{2}$
---|---|---|---|---|---
(coverage / set size) | (coverage / set size) | (coverage / set size) | (coverage / set size)
ImageNet-V2 | Resnet50 | 0.81 / 2.5 | 0.81 / 2.5 | 0.91 / 7.6 | 0.91 / 8.7
Resnet152 | 0.81 / 2.0 | 0.81 / 2.1 | 0.89 / 4.6 | 0.91 / 6.3
Vit-S | 0.80 / 1.5 | 0.80 / 1.5 | 0.90 / 3.4 | 0.90 / 3.4
ViT-B | 0.80 / 1.2 | 0.80 / 1.2 | 0.90 / 2.4 | 0.90 / 2.4
ImageNet-R | Resnet50 | 0.50 / 3.4 | 0.62 / 3.0 | 0.72 / 23.3 | 0.80 / 16.1
Resnet152 | 0.53 / 2.7 | 0.60 / 2.6 | 0.71 / 15.3 | 0.79 / 17.3
Vit-S | 0.52 / 1.3 | 0.53 / 1.3 | 0.74 / 12.3 | 0.75 / 11.8
ViT-B | 0.58 / 0.9 | 0.59 / 0.9 | 0.78 / 8.3 | 0.79 / 8.0
ImageNet-A | Resnet50 | 0.03 / 3.4 | 0.05 / 3.6 | 0.27 / 15.1 | 0.30 / 19.1
Resnet152 | 0.18 / 3.0 | 0.17 / 3.3 | 0.43 / 11.8 | 0.50 / 19.6
Vit-S | 0.37 / 1.7 | 0.37 / 1.7 | 0.65 / 8.4 | 0.66 / 8.3
ViT-B | 0.47 / 1.2 | 0.47 / 1.2 | 0.76 / 6.5 | 0.76 / 6.4

##  Appendix E Continuous shifts

In Figure 7, we visualize the coverage and set-sizes of our unsupervised
methods and a number of supervised baselines on the previously described
continuous distribution shifts. We show results on random corruption types as
well as fixed corruption types. Our proposed methods perform well across all
these settings; they closely maintain coverage even on sudden and severe
shifts, while leading to substantially smaller set sizes than the baselines.

(a) Shifting (random) corruptions

(b) Contrast corruption

(c) Brightness corruption

(a) Gaussian noise corruption

(b) Motion blur corruption

Figure 7: Our unsupervised methods ECP and EACP are able to provide nearly the
same empirical coverage, and considerably smaller set sizes, that supervised
methods on continuously shifting distributions. Dashed lines denote methods
that rely on a ground truth label being revealed at test time.

##  Appendix F ImageNet-C all severity levels

In Figure 8, we present full results across all ImageNet-C severity levels. We
see that our method is effective in recovering coverage even under many highly
severe distribution shifts, and nearly always recovers the desired coverage on
less severe shifts.

(a) Severity level 1

(b) Severity level 2

(c) Severity level 3

(d) Severity level 4

(e) Severity level 5

Figure 8: Performance on 19 ImageNet-C corruptions on each severity level.
$\texttt{E{\small A}CP}_{2}$ hugs the desired coverage line on nearly all
severity levels.