# Split Hierarchical Variational Compression

Tom Ryder Chen Zhang∗† Ning Kang† Shifeng Zhang†
†Huawei Noah’s Ark Lab
{chenzhang10, kang.ning2<EMAIL_ADDRESS>co-first author. The work
of Tom Ryder is conducted during his employment at Huawei Technologies R&D UK.

###### Abstract

Variational autoencoders (VAEs) have witnessed great success in performing the
compression of image datasets. This success, made possible by the bits-back
coding framework, has produced competitive compression performance across many
benchmarks. However, despite this, VAE architectures are currently limited by
a combination of coding practicalities and compression ratios. That is, not
only do state-of-the-art methods, such as normalizing flows, often demonstrate
out-performance, but the initial bits required in coding makes single and
parallel image compression challenging. To remedy this, we introduce Split
Hierarchical Variational Compression (SHVC). SHVC introduces two novelties.
Firstly, we propose an efficient autoregressive prior, the autoregressive sub-
pixel convolution, that allows a generalisation between per-pixel
autoregressions and fully factorised probability models. Secondly, we define
our coding framework, the autoregressive initial bits, that flexibly supports
parallel coding and avoids – for the first time – many of the practicalities
commonly associated with bits-back coding. In our experiments, we demonstrate
SHVC is able to achieve state-of-the-art compression performance across full-
resolution lossless image compression tasks, with up to 100x fewer model
parameters than competing VAE approaches.

##  1 Introduction

The volume of data, measured in terms of IP traffic, is currently witnessing
an exponential year-on-year growth [9]. Fuelled by the demand for high-
resolution media content, it is estimated that 80% of this data is in the form
of images and video [9]. Data service providers, such as cloud and streaming
platforms, have consequently seen costs associated with transmission and
storage become prohibitively expensive. For example, an increased demand for
streaming services forced major providers to throttle the maximum resolution
of video content to 720p during the coronavirus pandemic. As such, these
challenges have renewed the need for the development of high-performance data
compression codecs.

One solution to this problem has been the development of approaches using
likelihood-based generative models capable of discrete density estimation [22,
35, 14, 3, 13, 36, 18, 24, 6, 43, 42]. Such methods operate by learning a deep
probabilistic model of the data distribution, which, in combination with
entropy coders, can be used to compress data. Here, according to Shannon’s
source coding theorem [21], the minimal required average codelength is bounded
by the expected negative log-likelihood of the data distribution.

From this family of generative models, there have emerged three dominant modes
for data compression: normalizing flows [14, 3, 43, 42], variational
autoencoders [36, 18, 24] and autoregressive models [31, 37, 15] 111Most
recently, score-based generative models have been adapted for data compression
[16], but current approaches require an impractical number of operations at
inference time.. In fact, each of these approaches can be thought of as a
traversal on the Pareto frontier of inference speed and compression
performance. With broad generality, autoregressive models can often be the
most powerful but the slowest; variational autoencoders are often the weakest
but the fastest; and normalizing flows – depending on the variant – sit
somewhere in between.

In this paper, we consider data compression with VAEs, and focus on extending
the efficient frontier; obtaining solutions faster than popular VAEs that
achieve state-of-the-art compression ratios. Use of VAEs, however, poses two
outstanding challenges. Firstly, we should achieve competitive coding ratios
without greatly sacrificing time complexity. For example, best iterates
currently require one of two ingredients to improve performance: building
either a deep hierarchy of latent variables [8] or use of autoregressive
priors [29, 11]. The latter idea, especially popular in the codecs of the
lossy compression community [25], posits a model that flexibly learns both
local (via autoregression) and global (via hierarchical latent representation)
data modalities (e.g. low-frequency information). Whilst these approaches,
such as MS-PixelCNN [29] and PixelVAE [11], have had some success in achieving
more efficient trade-offs, generation of even moderately sized images is still
to the order of minutes [22].

Secondly, there should exist a practical means by which to efficiently perform
single-image compression. Single-image compression then permits parallel
coding, which is highly desirable. However, translating a VAE into a lossless
codec is currently achieved using the bits-back coding framework
(predominantly, bits-back ANS), which requires a large number of initial bits
[36, 39, 12] (see Section 3.1). Whilst this is a trivial number of bits on
large image datasets (where we can amortize this cost), it renders bits-back
an impractical approach for single-image compression. Furthermore, even large
datasets are often coded such that images are interlinked. Access to a single
image in the middle of a sequence would therefore require all prior images in
the bitstream to be additionally decompressed.

To that end, we propose two novelties for use in VAE-based compression
designed to address these challenges. The first, our autoregressive sub-pixel
convolution, introduces a simple autoregressive factorisation – not dissimilar
from the transformations used in normalizing flows [3, 43, 42] – designed to
present an efficient interpolation between fully-factorised probability
distributions and the impractical per-pixel autoregressions. Built from a
modified space-to-depth convolution operator, we losslessly downsample data
variables before performing a computationally efficient autoregression along
the channel dimension. Our autoregressive operator is then advantaged by a
number of network evaluations invariant to data dimensions, with each
autoregression crucially performed on a downsampled version on the input
tensor. More broadly, we view this framework as a generalisation of many
popular autoregressive “context” models used in data compression [31, 26, 11,
41].

Our second contribution, autoregressive initial bits, presents a general
framework for avoiding the impracticalities of bits-back ANS, allowing for
eminently parallelizable coding. This technique, highly compatible with our
autoregressive model, partitions the data variable into two splits such that
the second partition is conditionally independent of the latent variable(s),
given the first. In this way, we illustrate how we can use the entropy coding
of the conditionally independent partition to both supply and remove the
initial bits necessitated by bits-back ANS. We demonstrate that this approach
reduces the bit overhead on a per-image basis by close to 20x.

Finally, we combine the above contributions to present our codec, Split
Hierarchical Variational Compression (SHVC). SHVC posits a hierarchical VAE of
general-form autoregressive priors that permits parallel coding. Using our
framework, we outperform all other VAE-based compression approaches with fewer
latent variables and a comparable number of neural network evaluations. We
further illustrate the effectiveness of our architecture by training a small
model which is able to outperform similar VAE approach Bit-Swap [18] – but
with 100x fewer model parameters.

##  2 Related Work

Compression with VAEs can be separated into those with [18, 36] and those
without [22] stochastic posterior sampling (the latter uses a discrete
distribution of one symbol, assumed to have a probability of one). Whilst
obtaining theoretically superior compression ratios, approaches adopting
stochastic posteriors, such as HiLLoC [36] and Bit-Swap [18], must entropy
code using derivatives of the bits-back argument [12]. These approaches, under
the umbrella of bits-back ANS (bb-ANS), require access to an initial
bitstream. Although it is possible to amortize the cost of the initial
bitstream across a large dataset, single and parallel data compression has
presently proven challenging (see Section 3.1). In HiLLoC [36] the authors
propose use of a conventional codec to compress and send parts of a dataset,
which is then used as the basis for the initial bitstream. Whilst this
partially solves some of the coding challenges, it requires the implementation
of sub-optimal traditional techniques and still does not permit practical
coding of a single image. In contrast, our approach avoids all of these
challenges, requiring fewer latent variables, for a negligible additional bit
cost. Similarly, approaches such as L3C [22] that leverage deterministic
posteriors can also avoid the challenges associated with bb-ANS with use of
arithmetic (or adaptive arithmetic) coding (AC) [40]. Closely related approach
RC [24] – which can be loosely thought of as a VAE – uses a lossy compressed
image as a de facto latent variable to condition the data distribution.
However these techniques pay for their comparative practicality with a penalty
in compression ratio as they must explicitly code the joint distribution of
data and latent variables (see Section 3.2). Since this work focuses on VAEs
based codecs, we refer readers to [14, 3, 43, 42] and references therein for
compression with alternative deep generative models.

Figure 1: Left: The bb-ANS coding framework as discussed in Section 3.1.
Right: The coding framework of ArIB as proposed in Section 4.2. Here we nest a
bb-ANS coder inside of a block-based autoregressive structure that removes the
need for initial bits.

Autoregressive Models are a popular means to extend the independence
assumptions of fully factorised models to high-dimensional multivariate
densities. They are popular as both stand-alone models [38, 31] or in
application with VAEs [25, 29, 11]. In their most computationally expensive
forms, such as PixelCNN++ (and variants) [31], they posit a pixel-by-pixel
autoregression, which then codes in raster-scan order. Whilst of broad
academic interest, the $\mathcal{O}(n^{2})$ time complexity makes them
prohibitive for application in data compression. One proposed solution to this
problem has been to parametrise the priors of VAEs with autoregressive
densities. Here probability estimation proceeds by combining hierarchical
latent representations with decoded autoregressive context. Supplementing
autoregressive components with auxiliary latent features permits causality
restrictions that reduce time complexity, without greatly diminishing
performance. These restrictions include channel-wise autoregressions [26];
independent, block-based models [29]; “checkerboard” context [41]; and small
neural networks [11], amongst others. In fact, these restrictions are likely
of dual purpose: combining powerful autoregressive models with VAEs will
likely expedite posterior collapse (see Section 5.3) [5, 28, 11, 20]. Like
these techniques, our approach combines a VAE with a restricted autoregressive
model. Our method can be thought as most similar to [26] and [41]. Like the
former, we perform a channel-wise autoregression, but do so after our
autoregressive operator downsamples the data tensor. As such, our causality
more closely resembles that of [41]. However, in contrast to the authors of
[41], who enforce their causality with a binary mask, we do so using our sub-
pixel convolution. This precipitates a greater degree of parallelism and
presents the flexibility to efficiently recover a number of causal dependency
schemes, such as PixelCNN++.

##  3 Background

Suppose access to a dataset of size $n$, $\\{x_{1},x_{2},\dots,x_{n}\\}$,
drawn from some intractable $p_{\text{data}}({x})$ that we wish to compress.
In order to achieve this, we introduce a discrete probability distribution
$p({x})$ that, in combination with entropy coding, requires a codelength of
$\sum_{i=1}^{n}-\log_{2}p(x_{i})$ bits to represent. Ideally, $p({x})$ should
closely resemble $p_{\text{data}}({x})$. In such a case, the average
codelength in the limit of $n\longrightarrow\infty$ is given by
$\mathbb{E}_{p_{\text{data}}}\left[-\log_{2}p(x)\right]\longrightarrow
H({x})$, where $H({x})$ is the entropy of the data. Here, the compression
scheme is said to be optimal under Shannon’s source coding theorem [32].

###  3.1 Bits-back ANS for VAEs

As discussed in Section 1, variational autoencoders (VAEs) are one popular
approach to estimating $p_{\text{data}}({x})$, which defines a latent variable
model such that

$\displaystyle p({x})=\int p({x},{z})d{z}=\int p({x}|{z})p({z})d{z},$ (1)

where $p({z})$ is the prior distribution over latent variable ${z}$. As
$p({x})$ is normally intractable, VAEs introduce an approximate posterior
$q(z|x)$, which is optimised to maximise a lower bound on the marginal
evidence, the Evidence Lower Bound (ELBO)

$\displaystyle\log p(x)\geq$
$\displaystyle\mathcal{L}=\mathbb{E}_{q}\bigg{[}-\log q(z|x)+\log
p(x|z)p(z)\bigg{]},$ (2)

where low-variance estimates of the expectation in (2) are via Monte-Carlo
integration and the reparametrization trick [17].

Entropy coding requires explicit probabilities of data symbols. However, in
VAEs, the model is factorised into a prior and a likelihood, and therefore
does not allow direct coding of the data. To remedy this, several authors have
proposed coding variants of bb-ANS [35, 18, 36]. This process is outlined as
follows, which, without loss of generality, we describe for a model with a
single latent variable. During compression, one decodes $z$ from some
auxiliary initial bits with $q(z|x)$; encodes $x$ with $p(x|z)$; and encodes
$z$ with $p(z)$ to obtain the complete bitstream. In the decompression stage,
one decodes $z$ from the bitstream with $p(z)$; decodes $x$ with $p(x|z)$;
encodes $z$ with $q(z|x)$ and thus returns the initial bits (hence bits-back
coding). This technique is visualised in Fg. 1 Left.

The first decoding step common in existing bb-ANS codecs requires access to an
initial bitstream. This requirement leads to several disadvantages when
compared to AC. Firstly, while the initial bits are returned after
decompression, the same bits are occupied and not readily readable beforehand.
Secondly, although we can amortize the cost of the initial bits across a large
dataset by chaining the compressed data, access to any given data point
requires decompression of all data points posterior to the one of our target
in the original data sequence. As such, compression of a single image carries
substantial overhead – and, by extension, so too do parallel coding
implementations.

###  3.2 Deterministic vs Stochastic Posterior Sampling

Within VAEs based lossless compression, the impracticalities associated with
bits-back coding are not the only choice. Indeed, approaches that use
“deterministic” posterior sampling [22] may eschew bb-ANS for AC. This
approach is almost ubiquitous in the lossless codecs of lossy approaches [2,
25, 7, 26], where eminently parallelizable, low-latency codecs are especially
preferable (e.g. streaming media). We note that when using deterministic
posterior sampling, the likelihoods associated with the prior in Eq. (2)
trivially cancel to zero, such that the objective resolves to maximum
likelihood estimation of the joint distribution (see e.g. [2]).

Whilst gaining notable practical coding advantages, in sacrificing a
stochastic posterior one also sacrifices the ability to minimize the cost of
sending latent variables. To offset this limitation, models will repeatedly
downsample (RDS) the number of symbols available to latent representations
between each layer. Whilst this limits posterior expressiveness, there is
little experimental evidence to support how much this matters in practice. In
addition, models with stochastic posteriors require large $L$ to excel (where
$L$ is the number of latent variables), which hinders run-time.

Figure 2: Average BPD savings against the number of latent variables
($L=1,...,5$) used in three VAE lossless compression models. The baseline is
the deterministic posterior with $L=1$ where the BPD saving is zero. Results
are presented as an average for the three approaches trained across CIFAR10,
ImageNet32 and ImageNet64.

To that end, we display the results of a simple experiment in Fg. 2 designed
to investigate this difference further. Here we train three VAEs across
CIFAR10, ImageNet32 and ImageNet64: two with stochastic posteriors (one with
and one without RDS) and a deterministic posterior (with RDS). The
architectures in each model are identical (with the exception of downsampling
operations) and we quantify compression performance in bits per dimension
(BPD). Further experimental details can be found in the Appendix. Here we
observe that, even with $L\leq 3$, both stochastic posteriors outperform by
$\sim 5\%$. This difference grows as $L$ increases – but it does not matter if
RDS is used (at least, for $L\leq 5$). This result is of important
consequence: the best current approach to avoid the impracticalities of bb-ANS
(ie. using a deterministic posterior) carries a 5% BPD penalty. For single-
image compression with stochastic posteriors, the initial bits required would
typically be much larger than this. Likewise, unless extending to a deep
hierarchy of latent variables, RDS seems like a compute-efficient choice that
does not limit performance.

##  4 Method

Our method posits a hierarchical VAE where we parameterize the priors using an
autoregressive factorisation. We begin by defining a lossless downsampling
convolution operator, before describing its application to density estimation
using both weak and strong autoregressive models. We then describe how this
autoregressive structure can be leveraged to avoid many of the challenges
associated with bb-ANS without sacrificing the performance of stochastic
posteriors. Finally, we describe how these contributions can be combined to
form our SHVC codec.

Figure 3: Left: a $3\times 4\times 4$ input RGB image. Centre: the image $x$
downsampled using the convolution operator of (15) with $k=2$. Right: the
image $x$ downsampled using the convolution operator of (8) with $k=2$.

###  4.1 Autoregressive Sub-Pixel Convolutions

The space-to-depth and depth-to-space transformations are popular operations
across image analysis, from generative modelling [14, 3] to super-resolution
[33]. They define adjacent operations for efficient up and downsampling
transformations by folding spatial dimensions into channel dimensions – and
vice versa. Unlike learned operations, they greatly reduce computational
complexity, allowing for greater parallelism by losslessly moving computation
(and data) into the channels. Indeed, these operations have become an
essential component in papers seeking real-time execution (e.g. [33, 30, 10,
19]). Specifically, given a tensor of $C$ channels, $H$ height and $W$ width,
we define the space-to-depth and depth-to-space transformations, $f$ and
$f^{-1}$, such that

$\displaystyle f$ $\displaystyle:\mathbb{R}^{C\times H\times
W}\longrightarrow\mathbb{R}^{Ck^{2}\times\frac{H}{k}\times\frac{W}{k}},$ (3)
$\displaystyle f^{-1}$ $\displaystyle:\mathbb{R}^{C\times H\times
W}\longrightarrow\mathbb{R}^{\frac{C}{k^{2}}\times Hk\times Wk},$ (4)

where $k$ is the scale factor.

As described in [33], these operations can be efficiently performed using sub-
pixel convolutions, which are referred to as pixel unshuffle and pixel
shuffle. In particular, their space-to-depth transformation, pixel unshuffle,
is performed using a $k$-stride depthwise convolution where the $n^{th}$
element of $Ck^{2}$ $k\times k$ filters has one non-zero element such that

$K_{h,w}^{(n)}=\begin{cases}1&\text{if }h=\lfloor
n\mathbin{/}k\rfloor\,\text{mod}\,k,\,w=n\,\text{mod}\,k\\\
0&\text{else}\end{cases},$ (5)

where $h,w$ are the indices over spatial dimensions. The result of this
operation is visualised in Fg. 3 Centre.

Defining a channel-wise autoregression over the resulting tensor would posit a
checkerboard autoregressive structure over each of the channels in the
original tensor, sequentially. However, as identified in PixelCNN++ [31], sub-
pixels in adjacent channels, sharing the same spatial location in the original
tensor, have high correlation and therefore do not require complex models to
describe the dependency structure. As such, the authors of PixelCNN++ use a
linear model predicted by a single network evaluation, conditioned on decoded
context, to define the joint distribution across channels. In this way, they
obfuscate the need for separate RGB network evaluations. (We note that in our
setting, context refers to previously decoded pixels in either the current or
previous hierarchical latent variable.) From henceforth what we refer to as a
weak autoregression, is then defined similarly to [31] according to

$\displaystyle
p\left(x_{0:C,h,w}|D\right)=p\left(x_{0,h,w}|D\right)\prod_{c=1}^{C}p\left(x_{c,h,w}|x_{<c,h,w},D\right)$
(6)

where $D$ is the decoded context and $p$ is some parametric probability mass
function (pmf), obtained via integrating a probability density function (pdf)
over discretization bins, with mean at channel $c$ location $h,w$ given by

$\mu_{c,h,w}=\alpha_{c,h,w}+\sum_{i=0}^{c-1}\beta_{c,h,w}^{(i)}x_{i,h,w}.$ (7)

Here $\alpha$ and $\beta$ are scalars predicted for all channels and spatial
locations by a single network evaluated on decoded context, and $i$ is the
index over channels in decoded context such that $\beta_{c,h,w}^{(i)}$ is the
scalar for prediction of the mean associated with pixel at channel $i$,
spatial location $h,w$.

Inspired by this, we introduce a new space-to-depth convolution such that the
resulting autoregression is alternatively re-ordered into $k^{2}$ sub-blocks
of $C$ channels each. Crucially, the resulting channels in each sub-block
share the same spatial index allowing application of the autoregression
detailed in (6) and (7). We note that should $k=H=W$ we return an equivalence
to the per-pixel autoregression of PixelCNN++ but perform an autoregression
exclusively in the channel dimension. Likewise should $k<H$, we define a
block-based context model in raster scan order where, unlike MS-PixelCNN,
adjacent blocks are dependent.

To achieve our desired downsampling operation, which we denote by $g(\cdot)$,
we expand the depthwise convolutions of (15) into regular three-dimensional
kernels where the $n^{th}$ of $Ck^{2}$ $C\times k\times k$ filters has one
non-zero element such that

$K_{c,h,w}^{(n)}=\begin{cases}1&\text{if }c=n\,\,\text{mod}\,\,C,\,h=\lfloor
n\mathbin{/}Ck\rfloor\,\text{mod}\,k,\\\ &\,\,\,\,w=\lfloor
n\mathbin{/}C\rfloor\,\text{mod}\,k\\\ 0&\text{else}\end{cases}.$ (8)

We further visualise this operation in Fg. 3 Right. The resulting density of
the downsampled tensor for spatial location $h,w$ is then given by

$\displaystyle
p\left(g\left(x;k\right)_{0:Ck^{2},h,w}|D\right)=\prod_{i=0}^{k^{2}}\bigg{[}p\left(g\left(x;k\right)_{iC,h,w}|D\right)\times$
$\displaystyle\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\times\prod_{j=iC+1}^{(i+1)C-1}p\left(g\left(x;k\right)_{j,h,w}|g\left(x;k\right)_{<j,h,w},D\right)\bigg{]},$
(9)

where $i$ is the index over sub-blocks (i.e. a strong autoregression evaluated
using neural networks).

#### Masked 3D Convolutions

Whilst we are restricted to $k^{2}$ evaluations per latent variable at
inference time, the same does not have to be true during training. One
efficient parallel training scheme is to use 3-dimensional convolutions
applied to the downsampled tensors by expanding them into a $d\times
Ck^{2}\times H\times W$ volume [27, 23], where $d$ is some auxiliary
dimension. Here we can apply zero-masking along the channel dimension of the
kernels to enforce the causality condition, along with k-stride channel
convolutions on the input. Full details are available in the Appendix.

#### Choice of Distribution

For our choices of $p$ and $q$, we use a discretized mixture of logistic
distributions for $x$ and a discretized univariate logistic distribution for
all $z^{(l)}$ [31]. That is, given some mean $\mu$, scale $s$ and uniform
discretization bin-width $b$, one can obtain the univariate pmf by integrating
the logistic pdf over the discretiztion bin. For $x$, we typically use a
mixture of 5 discrete logistic distributions as defined above.

###  4.2 Autoregressive Initial Bits

As discussed in Section 3.1, bb-ANS is able to achieve efficient codelengths,
but can lead to several shortcomings. Fortunately, our proposed autoregressive
model naturally accommodates the possibility to bypass the auxiliary bits
needed in other bb-ANS methods. We achieve this by exploiting the block-based
autoregressive structure on the data variable. We outline this process below,
which we refer to as autoregressive initial bits (ArIB). Different from the
models considered in existing VAE-based codecs, we remove the direct causality
between the latent variable $z$ and some partition of the data variable $x$.
In practice, we simply remove $z$ from $D$ in Eq. (9) for the final $n$ sub-
blocks in $p(x|z)$, along with the partition from $x$ in
$q\left(z^{(1)}|x\right)$. As a result, we factorise the likelihood as
$p(x|z)=p(x_{s+1:k^{2}}|x_{1:s})p(x_{1:s}|z)$ with the approximate posterior
as $q(z|x)=q(z|x_{1:s})$, where $s$ is our ‘split’ index. Instead of
conducting the first step by decoding $z$ from $q(z|x)$, one can encode
$x_{s+1:k^{2}}$ with $p(x_{s+1:k^{2}}|x_{1:s})$ and thus obtain the bitstream
from which to decode $z$. Then one decodes $z$ with $q(z|x_{1:s})$, encodes
$x_{1:s}$ with $p(x_{1:s}|z)$ and encodes $z$ with $p(z)$. At the
decompression stage, one decodes $z$ with $p(z)$, decodes $x_{1:s}$ with
$p(x_{1:s}|z)$, encodes $z$ with $q(z|x_{1:s})$ and decodes $x_{s+1:k^{2}}$
with $p(x_{s+1:k^{2}}|x_{1:s})$. We illustrate this technique in Fg. 1 Right.

For this approach to be valid, we require the satisfaction of two criteria:

  1. 1.

There exists some $s$, $k$ and $z$ such that imposing
$\left(x_{s+1:k^{2}}\perp z|x_{1:s}\right)$ does not greatly hinder
performance.

  2. 2.

The entropy of $p\left(x_{s+1:k^{2}}|x_{1:s}\right)$ and
$q\left(\hat{z}|x_{1:s}\right)$, where $\hat{z}$ is the discretized analogue
of $z$, should be such that
$\mathcal{H}_{p(x_{s+1:k^{2}}|x_{1:s})}\geq\mathcal{H}_{q(\hat{z}|x_{1:s})}$.

In our experiments, we demonstrate that the performance costs associated with
criteria one are negligible. Crucially, we demonstrate that it is both orders
of magnitude less that initial bits required of vanilla bb-ANS and a
parameterisation of our approach using deterministic posteriors.

For criteria two, we formulate the optimization of (2) as a constrained
problem subject to
$\mathcal{H}_{p(x_{s+1:k^{2}}|x_{1:s})}\geq\mathcal{H}_{q(\hat{z}|x_{1:s})}$,
where we estimate the respective expectations during training using Monte-
Carlo integration. Whilst a variety of techniques from optimization theory may
be applied, we found it sufficient to simply penalise (2) according to

$\mathcal{L}_{pen}=\mathcal{L}+\lambda\max\left(0,\mathcal{H}_{q(\hat{z}|x_{1:s})}-\mathcal{H}_{p(x_{s+1:k^{2}}|x_{1:s})}\right),$
(10)

where $\lambda$ is some Lagrange multiplier. We find that this further
presents flexibility when choosing $s$, with a variety of choices yielding the
same result.

| Compression Model | CIFAR10 | ImageNet32 | ImageNet64 | CLIC.mobile | CLIC.pro | DIV2K
---|---|---|---|---|---|---|---
Generic | PNG [4] | 5.71 | 5.87 | 6.39 | 3.90 | 4.00 | 3.09
| FLIF [34] | 4.19 | 4.19 | 4.52 | 2.49 | 2.78 | 2.91
| JPEG-XL [1] | 5.74 | 5.89 | 6.39 | 2.36 | 2.63 | 2.79
VAE-Based | L3C [22] | - | 4.76 | 4.42 | 2.64 | 2.94 | 3.09
| Bit-Swap [18] | 3.82 | 4.50 | - | - | - | -
| HiLLoC [36] | 3.56‡ | 4.20‡ | 3.90‡ | - | - | -
| SHVC | 3.16/3.41‡ | 3.98 | 3.68/3.71‡ | 1.96∗ | 2.02∗ | 2.57∗
| SHVC Lite | 3.76 | 4.49 | 4.16 | - | - | -
Flow-Based | IDF [14] | 3.34/3.60‡ | 4.18 | 3.90/3.94‡ | - | - | -
| IDF++ [3] | 3.26 | 4.12 | 3.81 | - | - | -
| LBB [13] | 3.12 | 3.88 | 3.70 | - | - | -
| iVPF [43] | 3.20/ 3.49‡ | 4.03 | 3.75/3.79‡ | 2.39∗ | 2.54∗ | 2.68∗
| iFlow [42] | 3.12/3.36‡ | 3.88 | 3.70/3.65‡ | 2.26∗ | 2.44∗ | 2.57∗
Table 1: Compression results in BPD for SHVC and other popular codecs across
three low-resolution datasets and three full-resolution datasets. Lower is
better. Here ‡ and ∗ denote models trained on ImageNet32 and ImageNet64,
respectivley.

###  4.3 Split Hierarchical Variational Compression

SHVC formulates a hierarchical VAE built from the components described above.
Here we partition the latent variable into a simple disjoint hierarchy of $L$
layers, such that $z=\\{z^{(1)},...,z^{(L)}\\}$. We define the prior and
posterior according to

$\displaystyle p\left(x,z^{(1:L)}\right)$
$\displaystyle=p\left(x|z^{(1)}\right)p\left(z^{(L)}\right)\prod_{i=1}^{L-1}p\left(z^{(i)}|z^{(i+1)}\right),$
(11) $\displaystyle q\left(z^{(1:L)}|x\right)$
$\displaystyle=q\left(z^{(1)}|x\right)\prod_{i=1}^{L-1}q\left(z^{(i+1)}|z^{(i)}\right),$
(12)

where we parameterise every conditional density in Eq. (11) as per Eq. (9).
While this factorisation naturally fits the coding scheme proposed in Bit-Swap
[18], we additionally introduce a local reverse encoding to accommodate the
autoregressive structure for factors in Eq. (11). In more detail, for the
encoding of $z^{(i)}=[z^{(i)}_{1},...,z^{(i)}_{k^{2}}]$ with
$p(z^{(i)}|z^{(i+1)})$, one needs to encode in the reserved order of
$z^{(i)}_{k^{2}},...,z^{(i)}_{1}$, to accommodate the first-in-last-out nature
of ANS based codecs.

For purposes of experimentation, we define two versions of our model: one with
and one without the dependency structure permitting ArIB. From henceforth we
shall refer to these models as SHVC and SHVC-ArIB, respectivley. For SHVC, one
can encode $x$ with $p(x|z^{(1)})$, along with other variables in Eq. (11) as
discussed above. For SHVC-ArIB, one performs encoding and decoding for $x$ as
discussed in Section 4.2; and applies local reverse encoding for slices in
$x_{s+1:k^{2}}$ and $x_{1:s}$, respectively. We note that the only difference
in SHVC-ArIB is that, whilst $p(x_{s+1:k^{2}}|x_{1:s})$ and
$p(x_{1:s}|z^{(1)})$ are both modelled using (9), the former evidently omits
$z^{(1)}$ from $D$. In addition, we restrict the posterior such that
$q(z^{(1)}|x_{1:s})$. We visualise the overall architecture along with the
coding scheme for SHVC-ArIB in the Appendix.

##  5 Experiments

In this Section, we perform a series of experiments to evaluate the
effectiveness of compression with SHVC. We begin by discussing the
architecture and training details in Section 5.1. In Section 5.2, we evaluate
compression performance in terms on both low and full-resolution images. Here
we additionally evaluate the effect of the ArIB constraints on compression
performance and assess inference speed. Finally, in 5.3 we perform a series of
ablation studies.

###  5.1 Architecture and Training Details

As discussed in Section 3, we observe an efficient trade-off between inference
speed and performance where we repeatedly downsample and limit the number of
latent variable layers. Inspired by this, we apply repeated downsampling,
using three and four latent variables when training $32\times 32$ and
$64\times 64$ images, respectively. With these settings, we require 3 and 16
network evaluations for posterior and prior inference on CIFAR10. In contrast,
leading VAE approach HiLLoC [36] requires 24 prior and posterior evaluations.

For all of our experiments, we set $k=2$ in our autoregressive model. Whilst
this choice is discussed more extensively in Section 5.3, we found the $k=2$
presented an effective compromise between inference speed and representational
power. Indeed, increasing $k$ can often result in worse performance. We
observe that increasing $k$ increases the possibility of posterior collapse
[5, 28, 11, 20], making training a hierarchy of powerful autoregressive models
challenging. We further detail comprehensive architectures for all of our
experiments in the Appendix.

###  5.2 Compression Performance

#### Low-Resolution Images

We begin by testing our method on three toy datasets: CIFAR10, Imagenet32 and
Imagenet64. We compare our method against leading approaches from traditional
codecs, normalizing flows and VAEs. (Given the time complexity associated with
per-pixel autoregressive factorisations, we follow the broader compression
community and eschew them from our comparisons). As an ad hoc test of
generalizability, we additionally follow the authors of [36, 43, 14, 42] by
training a model on Imagenet32 and testing it across all other datasets. We
present our results in terms of bits per dimension (BPD), which we display in
Table 1.

#### Full-Resolution Images

We further compare our approach against full-resolution algorithms, i.e., L3C
[22]. Here we follow [43] and adopt our Imagenet64 model using a patch-based
evaluation protocol in which images are cropped to $64\times 64$. From Table
1, we demonstrate reliable out-performance of every other considered codec.
Further, when comparing the performance of SHVC to e.g. iFlow across small and
full-resolution images we interestingly note that difference in out-
performance becomes greater. Indeed, we hypothesise that autoregressive
”context” becomes an increasingly important inductive bias as the resolution
of the image increases – and likewise the extent of the spatial redundancy.

#### ArIB

ArIB poses non-trivial constraints on the data variable, such that a
sufficiently large partition should be: 1. conditionally independent of latent
variables; 2. encoded with an entropy larger than that with which the first
latent variable is decoded with. Using the objective of Eq. (10), we train
models capable of single-image compression: SHVC, SHVC-ArIB and a
parameterisation of SHVC with deterministic posteriors using AC (henceforth
known as Deterministic SHVC). We display results across CIFAR10, Imagenet32
and Imagenet64 in Fg. 4, and quantity our results as additional bits on a per-
image basis against SHVC. In more detail, the overhead in SHVC comes from
initial bits, whilst other models incur a performance cost. Here we see that
our ArIB adds minimal additional bits – less than the number in SHVC by a
factor of $\sim$20\. Crucially, our approach also outperforms Deterministic
SHVC, which we hope may serve as motivation to adapt the codec to lossy
compressors.

| SHVC | Bit-Swap | HiLLoC | IDF
---|---|---|---|---
BPD | 3.18 | 3.82 | 3.32 | 3.34
Time (s) | 4.63 | 5.86 | 10.20 | 20.58
Table 2: Inference time and BPD across CIFAR10. Inference time measures
evaluation of 10,000 test images with a batch size of 100.

#### Inference Speed

To better evaluate our approach in the context of VAEs, we compare SHVC
against popular publicly available compression models, HiLLoC 222In contrast
to SHVC and Bit-Swap, HiLLoC is implemented in TensorFlow. As such, HiLLoC is
likely advantaged in terms of run-time environment. [36], Bit-Swap [18] and
Integer Discrete Flow (IDF) [14] 333We note follow-up paper IDF++ [3] has no
publicly available code but uses largely the same architecture as IDF. As such
IDF++ inference times can be loosely inferred from IDF speeds.. We report
achieved BPDs in Table 2 and measure inference time in seconds (s) to evaluate
the 10,000 CIFAR10 test images with a batch size of 100. Here we observe that
our model is faster than HiLLoC for a lower BPD and achieves lower BPD for the
same speed as Bit-Swap. We further note that, unlike both considered
approaches, we are able to easily achieve parallel coding with minimal
overhead.

Figure 4: Additional bits per single image compressed across SHVC, SHVC-ArIB
and Deterministic SHVC. The overhead in SHVC comes from initial bits, whilst
other models incur a performance cost.

###  5.3 Ablation Studies

Figure 5: The theoretically minimal BPD against iteration, comparing the
common space-to-depth operation vs ours.

#### Space-to-depth operation

One alternative in our approach is to replace our downsampling operator of (8)
with the usual space-to-depth transformation described in Eq. (15). However,
as discussed, this would instead define a weak autoregressive property over
spatially adjacent pixels, channel-by-channel. In Fg. 5 we demonstrate the
differences resulting from the choice of spatial downsampling transformation
by training two models on CIFAR10. Here we observe that our convolutional
operator provides non-trivial benefits over the vanilla space-to-depth
transformation.

#### Choice of $k$

As discussed, one important hyper-parameter choice is that of $k$. As $k$
grows, the prior becomes more powerful but the time complexity grows. In Fg. 6
we visualise the effect of increasing $k$ on the compression ratio, which is
displayed as an average across three models trained on CIFAR10, ImageNet32 and
ImageNet64. As discussed, we note that performance of the model peaks at
$k=4$, before becoming worse as $k$ increases. This non-intuitive behaviour
can be better explained further in Fg. 6, where we see evidence of a posterior
collapse common in hierarchical VAEs – especially those with hierarchical
autoregressive priors [5, 28, 11, 20].

Figure 6: Top: (left) BPD savings against the choice of $k$ (the baseline is
$k=1$), and (right) the reconstruction loss as a percent of the total BPD.
Bottom: marginal densities of latent variables for different $k$.

##  6 Conclusion

We have proposed and evaluated a new VAE model for data compression. SHVC is
able to outperform existing VAE approaches in terms of both speed and
compression ratios. Additionally, unlike competing approaches using bb-ANS,
SHVC-ArIB is able to support parallel coding with minimal overhead. As such,
we believe ArIB represents the most promising means to convert a VAE into a
lossless codec. Motivated by this, one promising area of future work could
consider the application of SHVC to lossy compression, which has traditionally
ignored bits-back schemes.

Whilst our approach does not introduce any new negative societal consequences,
we remain acutely aware of the issues precipitated by this area of research.
These concerns are two-fold. Firstly, improving the capacity to store data
could greatly increase the amount of non-essential personal data held by
third-parties. This may be viewed as a challenge on broader digital liberties.
Secondly, any class of generative model trained on sensitive data will learn
to closely approximate the distribution of that data. As such, the model
itself might extend to malignant use-cases beyond its intended purpose., e.g.
classifiers.

## References

  * [1] Jyrki Alakuijala, Ruud van Asseldonk, Sami Boukortt, Martin Bruse, Iulia-Maria Comșa, Moritz Firsching, Thomas Fischbacher, Evgenii Kliuchnikov, Sebastian Gomez, Robert Obryk, et al.  Jpeg xl next-generation image compression architecture and coding tools.  In Applications of Digital Image Processing XLII, volume 11137, page 111370K. International Society for Optics and Photonics, 2019.
  * [2] Johannes Ballé, David Minnen, Saurabh Singh, Sung Jin Hwang, and Nick Johnston.  Variational image compression with a scale hyperprior.  In International Conference on Learning Representations, 2018.
  * [3] Rianne van den Berg, Alexey A Gritsenko, Mostafa Dehghani, Casper Kaae Sønderby, and Tim Salimans.  IDF++: Analyzing and improving integer discrete flows for lossless compression.  In International Conference on Learning Representations, 2021.
  * [4] Thomas Boutell and T Lane.  Png (portable network graphics) specification version 1.0.  Network Working Group, pages 1–102, 1997.
  * [5] Samuel R Bowman, Luke Vilnis, Oriol Vinyals, Andrew M Dai, Rafal Jozefowicz, and Samy Bengio.  Generating sentences from a continuous space.  2016\.
  * [6] Sheng Cao, Chao-Yuan Wu, and Philipp Krähenbühl.  Lossless image compression through super-resolution.  arXiv preprint arXiv:2004.02872, 2020.
  * [7] Zhengxue Cheng, Heming Sun, Masaru Takeuchi, and Jiro Katto.  Learned image compression with discretized gaussian mixture likelihoods and attention modules.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition (CVPR), 2020.
  * [8] Rewon Child.  Very deep VAEs generalize autoregressive models and can outperform them on images.  In International Conference on Learning Representations, 2021.
  * [9] Cisco.  Cisco visual networking index: global mobile data traffic forecast update, 2018–2023.  Update, 2018:2023, 2020.
  * [10] Tiago Cortinhal, George Tzelepis, and Eren Erdal Aksoy.  Salsanext: fast, uncertainty-aware semantic segmentation of lidar point clouds for autonomous driving.  arXiv preprint arXiv:2003.03653, 2020.
  * [11] Ishaan Gulrajani, Kundan Kumar, Faruk Ahmed, Adrien Ali Taiga, Francesco Visin, David Vazquez, and Aaron Courville.  Pixelvae: A latent variable model for natural images.  International Conference on Learning Representations, 2017.
  * [12] GE Hinton and Drew van Camp.  Keeping neural networks simple by minimising the description length of weights. 1993.  In Proceedings of COLT-93, pages 5–13.
  * [13] Jonathan Ho, Evan Lohn, and Pieter Abbeel.  Compression with flows via local bits-back coding.  In Advances in Neural Information Processing Systems, pages 3879–3888, 2019.
  * [14] Emiel Hoogeboom, Jorn Peters, Rianne van den Berg, and Max Welling.  Integer discrete flows and lossless compression.  In Advances in Neural Information Processing Systems, pages 12134–12144, 2019.
  * [15] Heewoo Jun, Rewon Child, Mark Chen, John Schulman, Aditya Ramesh, Alec Radford, and Ilya Sutskever.  Distribution augmentation for generative modeling.  In International Conference on Machine Learning, pages 5006–5019. PMLR, 2020.
  * [16] Diederik P Kingma, Tim Salimans, Ben Poole, and Jonathan Ho.  Variational diffusion models.  Advances in Neural Information Processing Systems, 2021.
  * [17] Diederik P Kingma and Max Welling.  Auto-encoding variational bayes.  International Conference on Learning Representations, 2014.
  * [18] Friso H Kingma, Pieter Abbeel, and Jonathan Ho.  Bit-swap: Recursive bits-back coding for lossless compression with hierarchical latent variables.  International Conference on Machine Learning, 2019.
  * [19] Haojie Liu, Tong Chen, Qiu Shen, Tao Yue, and Zhan Ma.  Deep image compression via end-to-end learning.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition Workshops, 2018.
  * [20] James Lucas, George Tucker, Roger B Grosse, and Mohammad Norouzi.  Don’t blame the elbo! a linear vae perspective on posterior collapse.  Advances in Neural Information Processing Systems, 32:9408–9418, 2019.
  * [21] David JC MacKay and David JC Mac Kay.  Information theory, inference and learning algorithms.  Cambridge university press, 2003.
  * [22] Fabian Mentzer, Eirikur Agustsson, Michael Tschannen, Radu Timofte, and Luc Van Gool.  Practical full resolution learned lossless image compression.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pages 10629–10638, 2019.
  * [23] Fabian Mentzer, Eirikur Agustsson, Michael Tschannen, Radu Timofte, and Luc Van Gool.  Conditional probability models for deep image compression.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition (CVPR), 2018.
  * [24] Fabian Mentzer, Luc Van Gool, and Michael Tschannen.  Learning better lossless compression using lossy compression.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 6638–6647, 2020.
  * [25] David Minnen, Johannes Ballé, and George D Toderici.  Joint autoregressive and hierarchical priors for learned image compression.  In S. Bengio, H. Wallach, H. Larochelle, K. Grauman, N. Cesa-Bianchi, and R. Garnett, editors, Advances in Neural Information Processing Systems, 2018.
  * [26] David Minnen and Saurabh Singh.  Channel-wise autoregressive entropy models for learned image compression.  In IEEE International Conference on Image Processing, ICIP, 2020\.
  * [27] Aaron Van Oord, Nal Kalchbrenner, and Koray Kavukcuoglu.  Pixel recurrent neural networks.  In Proceedings of The 33rd International Conference on Machine Learning, 2016.
  * [28] Ali Razavi, Aaron van den Oord, Ben Poole, and Oriol Vinyals.  Preventing posterior collapse with delta-VAEs.  In International Conference on Learning Representations, 2019.
  * [29] Scott Reed, Aäron van den Oord, Nal Kalchbrenner, Sergio Gómez Colmenarejo, Ziyu Wang, Yutian Chen, Dan Belov, and Nando de Freitas.  Parallel multiscale autoregressive density estimation.  In Proceedings of the 34th International Conference on Machine Learning, 2017.
  * [30] Oren Rippel, Alexander G. Anderson, Kedar Tatwawadi, Sanjay Nair, Craig Lytle, and Lubomir Bourdev.  Elf-vc: Efficient learned flexible-rate video coding.  2021 IEEE/CVF International Conference on Computer Vision (ICCV), 2021.
  * [31] Tim Salimans, Andrej Karpathy, Xi Chen, and Diederik P Kingma.  Pixelcnn++: Improving the pixelcnn with discretized logistic mixture likelihood and other modifications.  International Conference on Learning Representations, 2017.
  * [32] C. E. Shannon.  A mathematical theory of communication.  SIGMOBILE Mob. Comput. Commun. Rev., 2001.
  * [33] Wenzhe Shi, Jose Caballero, Ferenc Huszár, Johannes Totz, Andrew P Aitken, Rob Bishop, Daniel Rueckert, and Zehan Wang.  Real-time single image and video super-resolution using an efficient sub-pixel convolutional neural network.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 1874–1883, 2016.
  * [34] Jon Sneyers and Pieter Wuille.  Flif: Free lossless image format based on maniac compression.  In 2016 IEEE International Conference on Image Processing (ICIP), pages 66–70. IEEE, 2016.
  * [35] James Townsend, Tom Bird, and David Barber.  Practical lossless compression with latent variables using bits back coding.  International Conference on Learning Representations, 2019.
  * [36] James Townsend, Thomas Bird, Julius Kunze, and David Barber.  Hilloc: Lossless image compression with hierarchical latent variable models.  International Conference on Learning Representations, 2020.
  * [37] Aaron Van den Oord, Nal Kalchbrenner, Lasse Espeholt, Oriol Vinyals, Alex Graves, et al.  Conditional image generation with pixelcnn decoders.  In Advances in Neural Information Processing Systems, pages 4790–4798, 2016.
  * [38] Aaron Van Oord, Nal Kalchbrenner, and Koray Kavukcuoglu.  Pixel recurrent neural networks.  In International Conference on Machine Learning, 2016.
  * [39] Chris S Wallace.  Classification by minimum-message-length inference.  In International Conference on Computing and Information, pages 72–81. Springer, 1990.
  * [40] Ian H. Witten, Radford M. Neal, and John G. Cleary.  Arithmetic coding for data compression.  Commun. ACM, 1987.
  * [41] Honglei Zhang, Francesco Cricri, Hamed R. Tavakoli, Nannan Zou, Emre Aksu, and Miska M. Hannuksela.  Lossless image compression using a multi-scale progressive statistical model.  In Proceedings of the Asian Conference on Computer Vision (ACCV), 2020.
  * [42] Shifeng Zhang, Ning Kang, Tom Ryder, and Zhenguo Li.  iflow: Numerically invertible flows for efficient lossless compression via a uniform coder.  Advances in Neural Information Processing Systems, 2021.
  * [43] Shifeng Zhang, Chen Zhang, Ning Kang, and Zhenguo Li.  ivpf: Numerical invertible volume preserving flow for efficient lossless compression.  Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, 2021.

## A.1 Experimental Details for Section 3.2

For both the encoders and decoders, we use four fully-connected two-
dimensional convolutional layers with 128 channels and a 3x3 kernel. We
additionally use weight-normalization at each layer and PReLU activation
units. Where relevant, we downsample with 2-stride convolutions at the third
convolution, and upsample using a transposed convolution at the third
convolution.

For training, we use the Adam optimizer with default learning settings and an
initial learning rate of $5\times 10^{-4}$. We exponentially anneal this
learning rate to $1\times 10^{-5}$ during training.

Unlike stochastic posterior sampling, where we can train with continuous
latent variables because discretization schemes cancel across distributions,
deterministic posterior sampling requires discretization during training.
Because the discrteization operation (i.e. rounding) is not differentiable, we
adopt the popular technique of adding uniform noise during training, such that
our discretized latent variable is defined by

$\hat{z}^{(l)}=z^{(l)}+\epsilon,\quad\epsilon\sim\mathcal{U}\left(-\frac{\delta}{2},\frac{\delta}{2}\right),$
(13)

where $\delta$ is the uniform discretization bin and $\mathcal{U}$ represents
a uniform distribution. In practice, we take $\delta=1$.

## A.2 Masked 3D Convolutions

For large $k$ it becomes impractical to train using two-dimensional
convolutions. Doing so typically necessitates a serial scheme across data
partitions at a given latent variable. One approach to train our models in
parallel is to use masked three-dimensional convolutions. We achieve this by
expanding our downsampled data tensors into a $d\times Ck^{2}\times H\times W$
volume, where $d$ is some auxiliary dimension.

In order to retain the causality constraints, we build our approach of two
steps:

  1. 1.

We use an off-center convolution of stride $C$ to enforce the autoregressive
structure within sub-blocks. We define this operation as one that concatenates
a zero-tensor of dimension $d\times C\times H\times W$ to the data variable
along the auxiliary dimension and then applies the convolution as described.
The result of applying this convolution is a $f\times k^{2}\times H\times W$
tensor, where $f$ are the output channels of the convolution.

  2. 2.

We then apply repeated masked three-dimensional convolutions to the output of
the the off-center convolution. To enforce the causality constraint between
sub-blocks we apply a point-wise mask to the kernels prior to convolution. We
define two types of masks: type ‘A’ and ‘B’. We use B-type masks at all
locations apart from the input, where use an A-type mask. We describe these
masks in more detail below, and visualise them for a $3\times 3\times 3$
kernel in Figure 7.

For a three-dimensional kernel of depth $d$, height $h$ and width $w$,
consider the following masks that we apply as a point-wise multiplication to
the kernel.

#### A-Type Mask

$M_{d,h,w}=\begin{cases}1&\text{if }d\leq\lfloor k/2\rfloor\\\
0&\text{else}\end{cases}.$ (14)

#### B-Type Mask

$M_{d,h,w}=\begin{cases}1&\text{if }d\leq\lceil k/2\rceil\\\
0&\text{else}\end{cases}.$ (15)

Figure 7: A $3\times 3\times 3$ kernel masking strategies. Mask type ‘B’ left;
type ‘A’ right. Dark grey elements indicates zeros; light elements indicate
ones. Masks are applied as point-wise multiplications to the kernel.

## A.3 Visualisation and Coding Scheme for SHVC-ArIB

#### Variable Dependencies

In Figure 8 and Figure 9, we illustrate the differences in the dependency
structures in the priors and posteriors of SHVC and SHVC-ArIB. For ease of
presentation, we do so using one latent variable (i.e. $L=1$) and assume
$k=2$. We further assume $z^{(1)}$ is two-times smaller in spatial resolution
that $x$ but has the same number of channels, i.e. $C=3$.

#### Coding Schemes

Here we use the above model specification as an example to illustrate the
encoding and decoding processes of SHVC and SHVC-ArIB. Encoding and decoding
algorithms for SHVC can be found in Algorithms 1 and 2. Encoding and decoding
algorithms for SHVC-ArIB can be found in Algorithms 3 and 4. At the global
level, the coding algorithm is consistent with that of Bit-Swap, and at the
local level, the encoding of slices in latent and the data is conducted in the
reverse order. Since the above model only involves one latent variable, the
global level Bit-Swap degenerates to the original bb-ANS.

Algorithm 1 SHVC Encoding

Input: data to compress $x$

Step 0: Get auxiliary initial bits $c_{0}$

Step 1: Decode $z^{(1)}$ with $q(z^{(1)}|x)$

Step 2: Encode $x$ with $p(x|z^{(1)})$

Encode $x_{12}$ with $p(x_{12}|x_{1:11},z^{(1)})$

Encode $x_{11}$ with $p(x_{11}|x_{1:10},z^{(1)})$

Encode $x_{10}$ with $p(x_{10}|x_{1:9},z^{(1)})$

…

Encode $x_{2}$ with $p(x_{2}|x_{1},z^{(1)})$

Encode $x_{1}$ with $p(x_{1}|z^{(1)})$

Step 3: Encode $z^{(1)}$ with $p(z^{(1)})$

Encode $z^{(1)}_{12}$ with $p(z^{(1)}_{12}|z^{(1)}_{1:11})$

Encode $z^{(1)}_{11}$ with $p(z^{(1)}_{11}|z^{(1)}_{1:10})$

Encode $z^{(1)}_{10}$ with $p(z^{(1)}_{10}|z^{(1)}_{1:9})$

…

Encode $z^{(1)}_{2}$ with $p(z^{(1)}_{2}|z^{(1)}_{1})$

Encode $z^{(1)}_{1}$ with $p(z^{(1)}_{1})$

Output: final bit stream $c$

Algorithm 2 SHVC Decoding

Input: bit stream $c$

Step 1: Decode $z^{(1)}$ with $p(z^{(1)})$

Decode $z^{(1)}_{1}$ with $p(z^{(1)}_{1})$

Decode $z^{(1)}_{2}$ with $p(z^{(1)}_{2}|z^{(1)}_{1})$

…

Decode $z^{(1)}_{10}$ with $p(z^{(1)}_{10}|z^{(1)}_{1:9})$

Decode $z^{(1)}_{11}$ with $p(z^{(1)}_{11}|z^{(1)}_{1:10})$

Decode $z^{(1)}_{12}$ with $p(z^{(1)}_{12}|z^{(1)}_{1:11})$

Step 2: Decode $x$ with $p(x|z^{(1)})$

Decode $x_{1}$ with $p(x_{1}|z^{(1)})$

Decode $x_{2}$ with $p(x_{2}|x_{1},z^{(1)})$

…

Decode $x_{10}$ with $p(x_{10}|x_{1:9},z^{(1)})$

Decode $x_{11}$ with $p(x_{11}|x_{1:10},z^{(1)})$

Decode $x_{12}$ with $p(x_{12}|x_{1:11},z^{(1)})$

Step 3: Encode $z^{(1)}$ with $q(z^{(1)}|x)$

Output: auxiliary initial bit stream $c_{0}$, data to decompress $x$

Algorithm 3 SHVC-ArIB Encoding

Input: data to compress $x$

Step 0: Get autoregressive initial bits by encoding $x_{7:12}$

Encode $x_{12}$ with $p(x_{12}|x_{1:11})$

…

Encode $x_{7}$ with $p(x_{7}|x_{1:6})$

Step 1: Decode $z^{(1)}$ with $q(z^{(1)}|x_{1:6})$

Step 2: Encode $x_{1:6}$ with $p(x_{1:6}|z^{(1)})$

Encode $x_{6}$ with $p(x_{6}|x_{1:5},z^{(1)})$

…

Encode $x_{1}$ with $p(x_{1}|z^{(1)})$

Step 3: Encode $z^{(1)}$ with $p(z^{(1)})$

Encode $z^{(1)}_{12}$ with $p(z^{(1)}_{12}|z^{(1)}_{1:11})$

Encode $z^{(1)}_{11}$ with $p(z^{(1)}_{11}|z^{(1)}_{1:10})$

Encode $z^{(1)}_{10}$ with $p(z^{(1)}_{10}|z^{(1)}_{1:9})$

…

Encode $z^{(1)}_{2}$ with $p(z^{(1)}_{2}|z^{(1)}_{1})$

Encode $z^{(1)}_{1}$ with $p(z^{(1)}_{1})$

Output: final bit stream $c$

Algorithm 4 SHVC-ArIB Decoding

Input: bit stream $c$

Step 1: Decode $z^{(1)}$ with $p(z^{(1)})$

Decode $z^{(1)}_{1}$ with $p(z^{(1)}_{1})$

Decode $z^{(1)}_{2}$ with $p(z^{(1)}_{2}|z^{(1)}_{1})$

…

Decode $z^{(1)}_{10}$ with $p(z^{(1)}_{10}|z^{(1)}_{1:9})$

Decode $z^{(1)}_{11}$ with $p(z^{(1)}_{11}|z^{(1)}_{1:10})$

Decode $z^{(1)}_{12}$ with $p(z^{(1)}_{12}|z^{(1)}_{1:11})$

Step 2: Decode $x_{1:6}$ with $p(x_{1:6}|z^{(1)})$

Decode $x_{1}$ with $p(x_{1}|z^{(1)})$

…

Decode $x_{6}$ with $p(x_{6}|x_{1:5},z^{(1)})$

Step 3: Encode $z^{(1)}$ with $q(z^{(1)}|x_{1:6})$

Step 4: Decode $x_{7:12}$ with $p(x_{7:12}|x_{1:6})$

Decode $x_{7}$ with $p(x_{7}|x_{1:6})$

…

Decode $x_{12}$ with $p(x_{12}|x_{1:11})$

Output: data to decompress $x$

Figure 8: A comparison of the factorisation used in the priors of SHVC (left)
and SHVC-ArIB (right). Variable groupings are represented by coloured blocks.
Arrows indicate explicit dependencies. In SHVC-ArIB, there is no direct link
between $z^{(1)}$ (blue) and $x_{7:12}$ (green). Figure 9: A comparison of the
factorisation used in the posteriors of SHVC (left) and SHVC-ArIB (right).
Variable groupings are represented by coloured blocks. Arrows indicate
explicit dependencies. In SHVC-ArIB, there is no direct link between $z^{(1)}$
(blue) and $x_{7:12}$ (green).

## A.4 SHVC Architecture and Experimental Details

For both our encoder and decoder architectures, we use 8-layer Residual
networks with PReLU activation units and weight-normalization. For CIFAR10 we
additionally use dropout layers between residual connections to prevent
overfitting.

To highlight the effectiveness of our approach, we additionally train a small
”Lite” model, which uses four fully-connected convolutional layers. Here we
reduce the number of channels as we downsample the latent variables across
layers. These are detailed as follows:

  * •

$p(x|z^{(1)})$ uses 32 channels.

  * •

$p(z^{(1)}|z^{(2)})$ uses 24 channels.

  * •

$p(z^{(2)}|z^{(3)})$ uses 16 channels.

  * •

$p(z^{(3)}|z^{(4)})$ uses 8 channels, if it exists.

For training, we use the Adam optimizer with default learning settings and an
initial learning rate of $5\times 10^{-4}$. We exponentially anneal this
learning rate to $1\times 10^{-5}$ during training. We further use gradient-
clipping to control for numerical stability.

We run all of our experiments on a single NVIDIA Tesla V100.


[a]Matteo Di Carlo [c,d]Nils Hermansson-Truedsson [a]Antonin Portelli

# Isospin-breaking corrections to light leptonic decays in lattice QCD+QED at
the physical point

Peter Boyle  Felix Erben  Vera Gülpers  Maxwell T. Hansen  Tim Harris  Raoul
Hodgson  Andreas Jüttner  Fionn Ó hÓgáin  James Richings  Andrew Z. N. Yong

###### Abstract

We report on the physical-point RBC/UKQCD calculation of the leading isospin-
breaking corrections to light-meson leptonic decays. This is highly relevant
for future precision tests in the flavour physics sector, in particular the
first-row unitarity of the Cabibbo-Kobayashi-Maskawa matrix containing the
elements $V_{us}$ and $V_{ud}$. The simulations were performed using Domain-
Wall fermions for $2+1$ flavours, and with isospin-breaking effects included
perturbatively in the path integral through order $\alpha$ and
$(m_{u}-m_{d})/\Lambda_{\mathrm{QCD}}$. We use QEDL for the inclusion of
electromagnetism, and discuss here the non-locality of this prescription which
has significant impact on the infinite-volume extrapolation.

##  1 Introduction

Leptonic decays of pions and kaons are important for flavour physics precision
tests of the Standard Model (SM), since they give access to the ratio
$|V_{us}|/|V_{ud}|$. The Cabibbo-Kobayashi-Maskawa (CKM) matrix elements
$V_{us}$ and $V_{ud}$ are fundamental parameters in the SM expected to satisfy
a unitarity relation with also $V_{ub}$, namely
$|V_{ud}|^{2}+|V_{us}|^{2}+|V_{ub}|^{2}=1$ [1]. With a precision goal at the
(sub-)per cent level, isospin-breaking effects in the strong and
electromagnetic sector have to be included in theoretical predictions from
e.g. lattice field theory [2].

In the following we consider leptonic decays
$P^{+}\rightarrow\ell^{+}\nu_{\ell}$, where $P^{+}$ is a pion or kaon, and the
lepton $\ell^{+}$ a muon, including isospin-breaking effects in the strong
sector ($m_{u}\neq m_{d}$) and from quantum electrodynamics (QED). To handle
the problematic photon zero-momentum modes prohibiting QED to be
straightforwardly defined in a finite volume with periodic boundary
conditions, we employ the QEDL prescription [3] where the spatial zero-
momentum modes are dropped. Due to the long-range nature of QED, only the
inclusive decay rate
$\Gamma\left[P^{+}\rightarrow\ell^{+}\nu_{\ell}(\gamma)\right]=\Gamma_{0}\left[P^{+}\rightarrow\ell^{+}\nu_{\ell}\right]+\Gamma_{1}\left[P^{+}\rightarrow\ell^{+}\nu_{\ell}\gamma\right]$
yields an infrared (IR) finite result [4]. As originally laid out in Ref. [5],
one may separately study the virtual and real radiative corrected decay rates,
$\Gamma_{0}$ and $\Gamma_{1}$, respectively, by adding and subtracting a
perturbatively calculated quantity with the same IR divergence. For instance,
it is therefore possible to evaluate the virtual decay rate on the lattice in
a spacetime of volume $T\times L^{3}$, and the real radiative decay
perturbatively with a photon of mass $m_{\gamma}$ and energy less than some
cut-off $\omega_{\gamma}$, so that

$\displaystyle\Gamma\left[P^{+}\rightarrow\ell^{+}\nu_{\ell}(\gamma)\right]=\lim_{L\rightarrow\infty}\left\\{\Gamma_{0}(L)-\Gamma_{0}^{\textrm{uni}}(L)\right\\}+\lim_{m_{\gamma}\rightarrow
0}\left\\{\Gamma_{0}^{\textrm{uni}}(m_{\gamma})+\Gamma_{1}(\omega_{\gamma},m_{\gamma})\right\\}\,.$
(1)

Here each term in brackets is IR-finite, with $L$ and $m_{\gamma}$ being used
as IR regulators, and $\Gamma_{0}^{\textrm{uni}}$ is a perturbatively
calculated virtual decay rate independent of the internal meson structure
(hence the superscript universal) [5, 6]. The quantity
$\Gamma_{0}^{\textrm{uni}}(L)$ was calculated for the first time in Refs. [6,
7]. As discussed in Section 3, $\Gamma_{0}^{\textrm{uni}}$ subtracts the
logarithmic IR-divergence and finite-volume effects through order $1/L$. In
Ref. [8], an extension to a function $\Gamma_{0}^{(n)}(L)$ subtracting
structure-dependent higher-order finite-volume effects was proposed. Here we
use $\Gamma_{0}^{(2)}(L)$ and therefore study the quantity

$\displaystyle\Gamma\left[P^{+}\rightarrow\ell^{+}\nu_{\ell}(\gamma)\right]$
$\displaystyle=\lim_{L\rightarrow\infty}\left\\{\Gamma_{0}(L)-\Gamma_{0}^{(2)}(L)\right\\}+\lim_{m_{\gamma}\rightarrow
0}\left\\{\Gamma_{0}^{\textrm{uni}}(m_{\gamma})+\Gamma_{1}(\omega_{\gamma},m_{\gamma})\right\\}$
$\displaystyle\equiv\lim_{L\rightarrow\infty}\Gamma_{P}^{\mathrm{latt}}(L)+\lim_{m_{\gamma}\rightarrow
0}\Gamma_{P}^{\mathrm{pert}}(m_{\gamma})\,,$ (2)

with $\Gamma_{0}(L)$ evaluated non-perturbatively and at leading order in
$\alpha$ and $m_{u}-m_{d}$ on the lattice.

The leading isospin-breaking corrections to the decay rates are defined
through
$\Gamma_{P}^{\mathrm{latt}}(L)=\Gamma_{\mathrm{P}}^{\textrm{tree}}\left[1+\delta
R_{P}^{\textrm{latt}}(L)-\delta R_{P}^{\textrm{(2)}}(L)\right]$ and
$\Gamma_{P}^{\mathrm{pert}}(m_{\gamma})=\Gamma_{\mathrm{P}}^{\textrm{tree}}\,\delta
R_{P}^{\textrm{pert}}(\omega_{\gamma},m_{\gamma})$ where the tree level decay
rate is given by [1]

$\displaystyle\Gamma^{\textrm{tree}}_{P}=\frac{G_{F}^{2}}{8\pi}\,|V_{q_{1}q_{2}}|^{2}\,m_{\ell}^{2}\left(1-\frac{m_{\ell}^{2}}{m_{P}^{2}}\right)^{2}\,m_{P}\,{f_{P}}^{2}\,.$
(3)

Here, $G_{F}$ is the Fermi constant, $V_{q_{1}q_{2}}$ is the CKM element
associated to the decaying meson $P$ comprised of valence quarks $q_{1,2}$,
and ${f_{P}}$ is the isospin-symmetric decay constant of $P$. Introducing the
combination $\delta R_{P}=\delta R_{P}^{\textrm{latt}}(L)-\delta
R_{P}^{\textrm{(2)}}(L)+\delta
R_{P}^{\textrm{pert}}(\omega_{\gamma},m_{\gamma})$, the ratio
$|V_{us}|^{2}/|V_{ud}|^{2}$ can be obtained from the ratio of inclusive kaon
and pion decay rates according to

$\displaystyle\frac{|V_{us}|^{2}}{|V_{ud}|^{2}}=\frac{\Gamma\left[K^{+}\rightarrow\mu^{+}\nu_{\mu}(\gamma)\right]}{\Gamma\left[\pi^{+}\rightarrow\mu^{+}\nu_{\mu}(\gamma)\right]}\left(\frac{m_{\pi}^{2}-m_{\mu}^{2}}{m_{K}^{2}-m_{\mu}^{2}}\right)^{2}\frac{m_{K}^{3}}{m_{\pi}^{3}}\,\left(\frac{{f_{\pi}}}{{f_{K}}}\right)^{2}\frac{1}{1+\delta
R_{K}-\delta R_{\pi}}\,.$ (4)

The ratio of decay constants as well as $\delta R_{K\pi}=\delta R_{K}-\delta
R_{\pi}$ can be predicted from the lattice [2], with the remaining parts taken
from experiments. In the following we will focus on our calculation of $\delta
R_{K\pi}$. For further details, we refer the reader to our paper Ref. [9]. So
far, there is only one other lattice prediction for this quantity, by the
RM123/Southampton (RM123S) collaboration [10, 11].

##  2 Calculating $\delta R_{K\pi}$ on the lattice

We consider the leptonic decay $P^{+}\rightarrow\ell^{+}\nu_{\ell}$, with
associated Euclidean 4-momenta in the rest frame of $P^{+}$,
$p=(im_{P},\mathbf{0})$, $p_{\ell}=(i\omega_{\ell},\mathbf{p}_{\ell})$ and
$p_{\nu}=(i\omega_{\nu},-\mathbf{p}_{\ell})$. We further denote polarisations
of the final state leptons by $r$ and $s$, which have to be summed over to get
the non-polarised decay rate. The finite-volume decay rate in the full QCD+QED
theory can thus be written
$\Gamma_{0}(L)=\mathcal{K}\,\left|\mathcal{M}\right|^{2}$ where the
kinematical factor $\mathcal{K}$ and the QCD+QED squared matrix element
$\left|\mathcal{M}\right|^{2}=\sum_{r,s}\left|\mathcal{M}^{rs}\right|^{2}$ are
given by

$\displaystyle\mathcal{K}=\frac{G_{F}^{2}}{16\pi}\,\left|V_{q_{1}q_{2}}\right|^{2}\,\frac{1}{2m_{P}}\left(1-\frac{m_{\ell}^{2}}{m_{P}^{2}}\right)\,,\quad\mathcal{M}^{rs}=\mathcal{Z}\,\left\langle\ell^{+},\,r,\,\mathbf{p}_{\ell};\nu_{\ell},\,s,\,\mathbf{p}_{\nu}\right|O_{W}\left|P^{+},\,\mathbf{0}\right\rangle\,.$
(5)

Here we introduced the four-fermion operator $O_{W}$ entering the effective
Hamiltonian responsible for the decay, and its renormalisation factor
$\mathcal{Z}$. We do not enter into the details regarding the renormalisation
procedure here, and refer the reader to Ref. [11].

Next define
$\left|\mathcal{M}\right|^{2}=\mathcal{Z}^{2}\,4\,m_{\ell}^{2}\left(1-m_{\ell}^{2}/m_{P}^{2}\right)\,\left|\mathcal{A}_{P}\right|^{2}$.
Plugging this into $\Gamma_{0}(L)=\mathcal{K}\,\left|\mathcal{M}\right|^{2}$
and expanding to leading order in isospin breaking effects yields

$\displaystyle\delta
R_{P}^{\textrm{latt}}=2\left(\frac{\delta\mathcal{A}_{P}}{{\mathcal{A}}_{P}^{(0)}}+\frac{\delta\mathcal{Z}}{\mathcal{Z}_{0}}-\frac{\delta
m_{P}}{{m_{P}^{(0)}}}\right)\,,$ (6)

where $\delta\mathcal{A}_{P}$ are the isospin-breaking corrections in the
decay amplitude $\mathcal{A}_{P}$, and similar definitions hold for $\delta
m_{P}$ and $\delta Z_{P}$. The isospin-symmetric quantities are
$\mathcal{Z}_{0}\,{\mathcal{A}}^{(0)}_{P}=\mathcal{Z}_{0}\,{\left\langle
0\right|A^{0}\left|P\right\rangle}^{(0)}=i\,{m_{P}}^{(0)}\,{f_{P}}$, where
$A^{0}=\bar{q}_{2}\gamma_{0}\gamma_{5}q_{1}$ is the zeroth component of the
axial current. From the mass-independence of $\mathcal{Z}$, it follows that
the difference $\delta R_{K\pi}$ is given by

$\displaystyle\delta R_{K\pi}$
$\displaystyle=2\left(\frac{\delta\mathcal{A}_{K}}{{\mathcal{A}^{(0)}}_{K}}-\frac{\delta
m_{K}}{m_{K}^{(0)}}\right)-2\left(\frac{\delta\mathcal{A}_{\pi}}{{\mathcal{A}^{(0)}}_{\pi}}-\frac{\delta
m_{\pi}}{m_{\pi}^{(0)}}\right)$ $\displaystyle-\left[\delta
R_{K}^{\textrm{(2)}}(L)-\delta R_{\pi}^{\textrm{(2)}}(L)\right]+\left[\delta
R_{K}^{\textrm{pert}}(\omega_{\gamma},m_{\gamma})-\delta
R_{\pi}^{\textrm{pert}}(\omega_{\gamma},m_{\gamma})\right]\,.$ (7)

We can get an understanding of the various isospin-breaking contributions
$\delta\mathcal{A}_{P}$ by expanding the squared matrix element around the
iso-symmetric point according to
$\left|\mathcal{M}\right|^{2}=\left|{\mathcal{M}}^{(0)}\right|^{2}+\delta_{\textrm{f}}\left|\mathcal{M}\right|^{2}+\delta_{\textrm{nf}}\left|\mathcal{M}\right|^{2}$.
The term $\delta_{\textrm{f}}\left|\mathcal{M}\right|^{2}$ corresponds to the
leading isospin-breaking corrections not involving the final-state lepton in
the decay, which are diagrammatically depicted in diagrams (a)–(e) in Fig. 1.
As can be seen, the leptonic parts are factorised. Conversely,
$\delta_{\textrm{nf}}\left|\mathcal{M}\right|^{2}$ contains the non-
factorisable diagrams (f)–(g). The remaining diagram, (h), does not need to be
considered as it has purely perturbative nature and thus cancels in the
difference in (1). Note that we employ an electroquenched approximation where
sea quarks are neutral. These neglected disconnected diagrams are expected to
be small [11].

Figure 1: The leading connected isospin-breaking corrections to the decay
$P^{+}\rightarrow\ell^{+}\nu_{\ell}$, where $P^{+}$ has valence quarks $q_{1}$
and $q_{2}$. The filled circles represent the meson interpolating operator,
the squares the weak current and the diamonds scalar insertions.

###  2.1 Extracting $\delta R_{K\pi}$ from correlators

The non-perturbative quantities needed for $\delta R_{P}^{\textrm{latt}}$ in
(6) can be extracted from Euclidean correlation functions on the lattice. We
denote the finite time-extent by $T$, and here choose the meson to be created
at time $-t$, the weak current to be at the origin and the lepton a time-
separation $t_{\ell}$ away from the weak current. The factorisable
contributions, which do not couple to the leptons, can be obtained by studying
the two correlators

$\displaystyle C_{\textrm{PP}}(t)$ $\displaystyle=\int
d^{3}\mathbf{x}\,\left\langle
0\right|\mathbb{T}\left\\{A^{0}(0,\mathbf{0})\,\phi^{\dagger}(-t,\mathbf{x})\right\\}\left|0\right\rangle\,,$
(8) $\displaystyle C_{\textrm{PA}}(t)$ $\displaystyle=\int
d^{3}\mathbf{x}\,\left\langle
0\right|\mathbb{T}\left\\{\phi(0,\mathbf{0})\,\phi^{\dagger}(-t,\mathbf{x})\right\\}\left|0\right\rangle\,,$
(9)

where $\phi$ is a pseudoscalar meson interpolator with the appropriate flavour
structure. The above correlators are defined in the full theory, i.e. QCD+QED,
but we will expand them around the iso-symmetric point through leading order
in isospin breaking. Performing spectral decompositions of the above
correlators in the iso-symmetric limit, and retaining only the leading forward
and backward propagating exponentials, yields

$\displaystyle C_{\textrm{PP}}^{(0)}(t)$
$\displaystyle=\frac{|Z_{P}|^{2}}{2\,m_{P}^{(0)}}\,\left[e^{-m_{P}^{(0)}t}+e^{-m_{P}^{(0)}\left(T-t\right)}\right]+\ldots\,,$
(10) $\displaystyle C_{\textrm{PA}}^{(0)}(t)$
$\displaystyle=\frac{Z_{P}\,{\mathcal{A}^{(0)}}_{P}}{2\,m_{P}^{(0)}}\,\left[e^{-m_{P}^{(0)}t}-e^{-m_{P}^{(0)}\left(T-t\right)}\right]+\ldots\,,$
(11)

where $Z_{P}=\langle P,\mathbf{0}|\phi^{\dagger}(0)|0\rangle$ and the ellipses
contain sub-leading exponentials. It is from the above spectral decompositions
that the iso-symmetric quantities in (6) can be fitted. Defining $\delta
C_{\textrm{PA}}(t)$ as the leading isospin-breaking contribution to
$C_{\textrm{PA}}(t)$ gives the factorisable correction
$\delta_{\textrm{f}}\mathcal{A}_{P}$ from the leading finite-$T$ behaviour of
the ratio

$\displaystyle\frac{\delta C_{\textrm{PA}}(t)}{C_{\textrm{PA}}^{(0)}(t)}$
$\displaystyle=\frac{\delta_{\textrm{f}}\mathcal{A}_{P}}{{\mathcal{A}^{(0)}}_{P}}+\frac{\delta
Z_{P}}{Z_{P}}-\frac{\delta m_{P}}{m_{P}^{(0)}}\,f_{PA}(t)+\ldots\,,$ (12)
$\displaystyle f_{PA}(t)$
$\displaystyle=1+m_{P}^{(0)}\left\\{\frac{T}{2}-\left(t-\frac{T}{2}\right)\,\coth\left[m_{P}^{(0)}\left(t-\frac{T}{2}\right)\right]\right\\}\,.$
(13)

A similar procedure can be used for $C_{\textrm{PP}}(t)$.

What remains is to determine the non-factorisable contribution
$\delta_{\textrm{nf}}\mathcal{A}_{P}$. In a fashion similar to above, we
define a correlation function $C_{P\ell}(t)$ containing the meson, the weak
current as well as the leptonic part. This can be appropriately traced with
spinors to yield a scalar object. However, the lepton propagator between the
weak current at the origin and the external lepton at $t_{\ell}$ affects the
finite-$T$ behaviour severely. By constructing projectors onto either the
forward or backward propagating signal, the numerical analysis can be
simplified. For brevity, we refrain from defining this correlator in detail,
and simply show the leading finite-$T$ behaviour when projecting onto the
forward propagating lepton, namely

$\displaystyle\frac{\delta
C_{P\ell}(t)}{\overline{C}_{P\ell}(t)}=\frac{\delta_{\textrm{nf}}\mathcal{A}_{P}}{{\mathcal{A}^{(0)}}_{P}}\,f_{P\ell}(t)+\ldots\,,\quad
f_{P\ell}(t)=\frac{1}{2}\left\\{1+\kappa_{\ell}-(1-\kappa_{\ell})\,\coth\left[m_{P}^{(0)}\left(t-\frac{T}{2}\right)\right]\right\\}\,,$
(14)

where $\kappa_{\ell}$ parametrises the backward propagating signal.

We thus see that all the ingredients needed for $\delta R_{P}^{\textrm{latt}}$
in (6) can be obtained from Euclidean correlation functions on the lattice. We
do simultaneous correlated fits of the factorisable and non-factorisable
correlators, respectively. We choose optimal fit ranges through a genetic
algorithm for the factorisable case, and construct $\delta R_{K\pi}$ from
these using AIC-based model averaging as in e.g. Refs. [12, 13]. This approach
results in a distribution of possible $\delta R_{K\pi}$, from which we choose
our prediction as the median and estimate the associated statistical and
systematic uncertainties.

###  2.2 Simulation details

For the simulations needed to generate the correlators discussed above we make
use of Grid [14] and Hadrons [15]. We have temporal and spatial extents
$T/a=96$ and $L/a=48$, respectively, and use domain-wall fermions close to the
physical point for $L_{s}/a=24$ and $aM_{5}=1.8$ [16]. We use 60 statistically
independent QCD gauge configurations generated using the Iwasaki action [17]
by the RBC/UKQCD collaboration. Our valence- and sea-quark masses are
$\hat{m}_{u}=\hat{m}_{d}=0.00078$ and $\hat{m}_{s}=0.0362$, where hat denotes
lattice units. The ensemble pion mass is $m_{\pi}=139.15(36)$ MeV, and the
inverse lattice spacing $a^{-1}=1.7295(38)$ GeV. We thus have a single lattice
spacing and volume.

The correlators are created using sequential propagators, with 96 Coulomb
gauge-fixed wall sources per configuration. The muon momentum
$\mathbf{p}_{\ell}$ fixed in the direction $(1,1,1)$ is injected via twisted
boundary conditions, and the muon propagator is evaluated for 8 different
source-sink separations. With the weak current at time $0$, these are
$t_{\ell}\in\\{12,16,20,24,28,32,36,40\\}$. The photons are implemented in
$\textrm{QED}_{\textrm{L}}$, with the photon fields sampled from a Gaussian
distribution, and we use a renormalised local vector current.

###  2.3 Defining the iso-symmetric theory

The bare parameters of the full QCD+QED theory can be unambiguously determined
by requiring that a complete set of experimental hadronic masses be
reproduced. Once determined, any physical quantity can be predicted.
Determining the isospin-breaking corrections to an observable, on the other
hand, requires the definition of an iso-symmetric point. Consider a physical
observable $X^{\phi}$. This can be separated into three terms according to
$X^{\phi}={X}^{(0)}+\delta X^{\gamma}+\delta X^{\textrm{SIB}}$, where
${X}^{(0)}$ is the iso-symmetric quantity, $\delta X^{\gamma}$ contains the
electromagnetic corrections where $\alpha\neq 0$, and $\delta
X^{\textrm{SIB}}$ includes the strong isospin-breaking effects from $m_{u}\neq
m_{d}$ for $\alpha=0$. Such a separation, however, requires additional
conditions and constitutes the choice of a scheme. Therefore, when comparing
predictions of isospin-breaking corrections, scheme ambiguities have to be
taken into account. This is particularly important since the ambiguities in
principle can be of the same size as the predictions themselves. Relations
between different schemes can be obtained due to the small numerical size of
isospin-breaking effects, which implies that unphysical theories to a good
approximation are within a linear deviation from the physical point.

In our calculation we have 3 flavours. Denoting lattice units by a hat, we
tune our parameters in the full theory to reproduce the ratios

$\displaystyle\frac{\hat{m}_{P}}{\hat{M}_{\Omega^{-}}}=\frac{m^{\textrm{PDG}}_{P}}{M^{\textrm{PDG}}_{\Omega^{-}}}\,,$
(15)

where $P=\pi^{+},K^{+},K^{0}$ and PDG denotes the Particle Data Group values
from Ref. [1]. In addition, we can define the lattice spacing as
$a=(\hat{M}_{\Omega^{-}})/M^{\textrm{PDG}}_{\Omega^{-}}$. To define the iso-
symmetric point we employ the same scheme as the BMW collaboration in Ref.
[18]. In short, this is done by considering non-physical, neutral, purely
connected mesons with quark content $\bar{q}q$ and masses $M_{\bar{q}q}$,
which are combined into the three variables (in the full theory)

$\displaystyle
M^{2}_{ud}=\frac{1}{2}\left(M^{2}_{\bar{u}u}+M^{2}_{\bar{d}d}\right)\,,\quad\Delta
M^{2}=M^{2}_{\bar{u}u}-M^{2}_{\bar{d}d}\,,\quad
2\,M^{2}_{K\chi}=M^{2}_{K^{+}}+M_{K^{0}}^{2}-M_{\pi^{+}}^{2}\,.$ (16)

The scheme is then defined by requiring that the meson masses above coincide
in the full theory and QCD, as well as imposing the conditions
$(\hat{M}^{(0)}_{ud}/\hat{M}_{\Omega^{-}})^{2}=(\hat{M}_{ud}/\hat{M}_{\Omega^{-}})^{2}$
,
$(\hat{M}_{K\chi}^{(0)}/\hat{M}_{\Omega^{-}})^{2}=(\hat{M}_{K\chi}/\hat{M}_{\Omega^{-}})^{2}$
and $(\Delta\hat{M}^{(0)}/\hat{M}_{\Omega^{-}})^{2}=0$ between the full and
iso-symmetric theories.

##  3 Finite-volume effects

The finite-volume effects in the decay rate $\Gamma_{0}^{(n)}(L)$ are
contained in the correction $\delta R_{P}^{(n)}(L)$ defined through

$\displaystyle\Gamma_{0}^{(n)}(L)$
$\displaystyle=\Gamma_{P}^{\textrm{tree}}\left[1+\delta
R_{P}^{(n)}(L)\right]\,,$ (17) $\displaystyle\delta R_{P}^{(n)}(L)$
$\displaystyle=2\,\frac{\alpha}{4\pi}\left[\widetilde{Y}_{P}(L)+\sum_{j=0}^{n}\frac{1}{(m_{P}L)^{j}}\,Y_{P,\,j}\right]\,.$
(18)

The universal $\widetilde{Y}_{P}(L)$ contains the infinite-volume contribution
as well as the logarithmic divergence in $L$. The numerical coefficients
$Y_{P,\,j}$ in general depend on the structure of the decaying mesons.
However, it was proven in Ref. [6] that through order $1/L$ the coefficients
are universal, i.e. structure-independent, and they were determined within
point-like scalar QED.

The structure-dependent contributions at higher order can be determined using
the approach in Ref. [8]. At present, the coefficients are known fully through
order $1/L^{2}$ and can be found in Ref. [8]. The structure dependence in
$Y_{P,\,2}$ is encoded in the on-shell axial-vector form factor $F_{A}^{P}$
from the real radiative decay $P\rightarrow\ell\nu_{\ell}\gamma$, which is
known experimentally [1], from chiral perturbation theory [19] and the lattice
[20]. However, it was shown in Ref. [8] that due to the smallness of
$F_{A}^{P}$, the structure-dependence at order $1/L^{2}$ is only a per cent
level contribution to $Y_{P,\,2}$, implying that $Y_{P,\,2}\approx
Y_{P,\,2}^{\textrm{pt}}$ is a good approximation.

At order $1/L^{3}$, only the point-like approximation
$Y_{P,\,3}^{\textrm{pt}}$ of the full coefficient
$Y_{P,\,3}=Y_{P,\,3}^{\textrm{pt}}+Y_{P,\,3}^{\textrm{sd}}$ is currently
known. The structure-dependent piece $Y_{P,\,3}^{\textrm{sd}}$ has yet to be
determined, which is complicated due to non-local effects arising at order
$1/L^{3}$ in $\textrm{QED}_{\textrm{L}}$ [8]. The point-like result is

$\displaystyle
Y_{P,\,3}^{\textrm{pt}}=\frac{32\pi^{2}c_{0}\,[2+(m_{\ell}/m_{P})^{2}]}{[1+(m_{\ell}/m_{P})^{2}]^{3}}\,.$
(19)

Here $c_{0}=-1$ is a finite-volume coefficient. As a final remark on the
finite-volume expansion, all the point-like contributions at order $1/L^{4}$
and higher vanish as is easily shown using the approach in Ref. [8].

In our analysis we subtract $\delta R_{P}^{(2)}(L)$ from $\delta
R_{P}^{\textrm{latt}}(L)$ to yield $\delta R_{K\pi}$ as in (2). At our
simulated volume $L/a=48$ we then use the point-like contribution $\delta
R_{P}^{(3),\,\textrm{pt}}(L)-\delta R_{P}^{(2)}(L)$ as a systematic error.
Denoting $\delta R_{K\pi}^{(n)}(L)=\delta R_{K}^{(n)}(L)-\delta
R_{\pi}^{(n)}(L)$ we find at the simulation point for $L_{48}=48a$

$\displaystyle\delta R^{(1)}_{K\pi}(L_{48})\approx-0.00468\,,\quad\delta
R^{(2)}_{K\pi}(L_{48})\approx-0.00730\,,\quad\delta
R^{(3),\,\textrm{pt}}_{K\pi}(L_{48})\approx-0.00337\,.$ (20)

The relative shift in going from $1/L^{2}$ to $1/L^{3}$ is roughly $-54\%$,
which means that we have a large systematic error from not knowing the full
coefficient $Y_{P,\,3}$. We are currently studying possible ways to determine
the unknown structure-dependent $Y_{P,\,3}^{\textrm{sd}}$.

Finally, we present the importance of analytical knowledge of the finite-
volume dependence for the infinite-volume extrapolation. For this we use the
published $\delta R_{K\pi}$ data of Ref. [11], which was kindly provided to us
by the authors. The data was produced for non-physical pions and kaons of
masses $m_{\pi}\approx 320$ MeV and $m_{K}\approx 580$ MeV, respectively, and
for several volumes. We here investigate the impact on the infinite-volume
extrapolation of this data from the recent observation that the structure-
dependence in $Y_{P,\,2}$ is negligible, which was not known at the time of
Ref. [11].

In Figs. 2(a)–(b) the volume-dependence of the pion and kaon data for $\delta
R_{P}$ is shown, with different orders of $\delta R_{P}^{(n)}(L)$ subtracted.
The circular points in the figures correspond to the results of Ref. [11]
where the universal $\delta R_{P}^{(1)}(L)$ are subtracted (labelled $1/L$
subtracted). As can be seen, a linear extrapolation to the infinite-volume
limit (dashed line) describes the data well. In Ref. [11] the point-like
approximation of $\delta R_{P}^{(2)}(L)$ was then subtracted, which here
corresponds to the square points, and due to the residual slope the authors of
Ref. [11] concluded that there is significant structure-dependence at order
$1/L^{2}$. However, from Ref. [8] we now know that the structure-dependence in
$Y_{P,\,2}$ is negligible, so that the slope of the $1/L^{2}$ subtracted data
in Fig. 2 must be from the $1/L^{3}$ coefficient. We therefore also subtract
the point-like $\delta R_{P}^{\textrm{pt},\,(3)}(L)$ to yield the diamond
points, and perform a fit of the form $a+b/L^{3}$ to those points (solid line
in the figure). The $1/L^{3}$ ansatz for the volume-dependence describes the
data well and estimates the remaining $1/L^{3}$-dependence. Using the fitted
curve and adding back the analytically known coefficients yields the remaining
solid lines, which also give a good description of the data. This means that a
sizeable $1/L^{3}$ dependence mimics a $1/L^{2}$ behaviour in the range of
volumes in the figure.

Although the infinite-volume extrapolations for the dashed and the solid
curves do not agree in the figure, we stress that this difference can be
washed out in predictions of $\delta R_{K\pi}$. For instance, the numerical
analysis to obtain the final value of $\delta R_{K\pi}$ in Ref. [11] involves
several steps beyond the infinite-volume extrapolation, such as a simultaneous
chiral extrapolation to physical meson masses, and the difference observed
here might be well within the associated uncertainties. In conclusion,
analytical knowledge of the finite-volume dependence is crucial when
extrapolating to infinite volume, and there is clear need to determine the
full coefficient $Y_{P,\,3}$.

(a)

(b)

Figure 2: The volume dependence of $\delta R_{P}$ for (a) pions and (b) kaons,
subtracted with $\delta R_{P}^{(n)}(L)$ for different $n$. The lattice data
points are from Ref. [11].

##  4 A prediction for $\delta R_{K\pi}$

We now present our prediction for the leading isospin-breaking to ratio of
light leptonic decays,

$\displaystyle\delta
R_{K\pi}=-0.0086(3)_{\textrm{stat.}}\left({}^{+11}_{-4}\right)_{\textrm{fit}}(5)_{\textrm{disc.}}(5)_{\textrm{quench.}}(39)_{\textrm{vol.}}\,.$
(21)

As can be seen the dominant contribution comes from the finite-volume
uncertainty. Further details on the error budget is provided in our paper [9].
This value for $\delta R_{K\pi}$ can be compared to those from RM123S [11] and
chiral perturbation theory ($\chi$PT) [21], which respectively are

$\displaystyle\delta R_{K\pi}^{\textrm{RM123S}}=-0.0126(14)\,,\quad\delta
R_{K\pi}^{\chi\textrm{PT}}=-0.0112(21)\,.$ (22)

All three values are in agreement, although our error bars are very large due
to the finite-volume systematics. It is clear that the finite-volume issue has
to be resolved in the future.

##  5 Conclusions

In these proceedings we have presented the RBC/UKQCD prediction for the
leading isospin-breaking corrections to leptonic decays of pions and kaons,
encoded in $\delta R_{K\pi}$. Further details on our calculation are presented
in Ref. [9]. Knowledge of this quantity is essential for future flavour
physics precision tests of the Standard Model, particularly as it is required
for testing CKM-matrix unitarity at (sub-)per cent precision.

## Acknowledgements

We thank the RM123S collaboration for kindly sharing their lattice data for
our study of the finite-volume dependence. This work used the DiRAC Extreme
Scaling service at the University of Edinburgh, operated by the Edinburgh
Parallel Computing Centre on behalf of the STFC DiRAC HPC Facility
(www.dirac.ac.uk). This equipment was funded by BEIS capital funding via STFC
capital grant ST/R00238X/1 and STFC DiRAC Operations grant ST/R001006/1. DiRAC
is part of the National e-Infrastructure. P.B. has been supported in part by
the U.S. Department of Energy, Office of Science, Office of Nuclear Physics
under the Contract No. DE-SC-0012704 (BNL). M.D.C., F.E., T.H., V.G., M.T.H.,
and A.P. are supported in part by UK STFC grant ST/P000630/1. Additionally
M.T.H. is supported by UKRI Future Leader Fellowship MR/T019956/1. F.E., V.G.,
R.H., F.Óh., A.P. and A.Z.N.Y. received funding from the European Research
Council (ERC) under the European Union’s Horizon 2020 research and innovation
programme under grant agreement No 757646 and A.P. additionally under grant
agreement No 813942. N. H.-T. is funded in part by the Albert Einstein Center
for Fundamental Physics at the University of Bern, and in part by the Swedish
Research Council, project number 2021-06638. A.J. and J.R have been supported
in part by UK STFC grant ST/P000711/1 and ST/T000775/1. J.R. is also supported
in part by UK STFC DiRAC operational grants ST/S003762/1 and ST/W002701/1.

## References

  * [1] Particle Data Group collaboration, _Review of Particle Physics_ , _PTEP_ 2022 (2022) 083C01.
  * [2] Flavour Lattice Averaging Group (FLAG) collaboration, _FLAG Review 2021_ , _Eur. Phys. J. C_ 82 (2022) 869 [2111.09849].
  * [3] M. Hayakawa and S. Uno, _QED in finite volume and finite size scaling effect on electromagnetic properties of hadrons_ , _Prog. Theor. Phys._ 120 (2008) 413 [0804.2044].
  * [4] F. Bloch and A. Nordsieck, _Note on the Radiation Field of the electron_ , _Phys. Rev._ 52 (1937) 54.
  * [5] N. Carrasco, V. Lubicz, G. Martinelli, C.T. Sachrajda, N. Tantalo, C. Tarantino et al., _QED Corrections to Hadronic Processes in Lattice QCD_ , _Phys. Rev. D_ 91 (2015) 074506 [1502.00257].
  * [6] V. Lubicz, G. Martinelli, C.T. Sachrajda, F. Sanfilippo, S. Simula and N. Tantalo, _Finite-Volume QED Corrections to Decay Amplitudes in Lattice QCD_ , _Phys. Rev. D_ 95 (2017) 034504 [1611.08497].
  * [7] N. Tantalo, V. Lubicz, G. Martinelli, C.T. Sachrajda, F. Sanfilippo and S. Simula, _Electromagnetic corrections to leptonic decay rates of charged pseudoscalar mesons: finite-volume effects_ , 1612.00199.
  * [8] M. Di Carlo, M.T. Hansen, N. Hermansson-Truedsson and A. Portelli, _Relativistic, model-independent determination of electromagnetic finite-size effects beyond the pointlike approximation_ , _Phys. Rev. D_ 105 (2022) 074509 [2109.05002].
  * [9] P. Boyle et al., _Isospin-breaking corrections to light-meson leptonic decays from lattice simulations at physical quark masses_ , 2211.12865.
  * [10] D. Giusti, V. Lubicz, G. Martinelli, C.T. Sachrajda, F. Sanfilippo, S. Simula et al., _First lattice calculation of the QED corrections to leptonic decay rates_ , _Phys. Rev. Lett._ 120 (2018) 072001 [1711.06537].
  * [11] M. Di Carlo, D. Giusti, V. Lubicz, G. Martinelli, C.T. Sachrajda, F. Sanfilippo et al., _Light-meson leptonic decay rates in lattice QCD+QED_ , _Phys. Rev. D_ 100 (2019) 034514 [1904.08731].
  * [12] S. Borsanyi et al., _Ab initio calculation of the neutron-proton mass difference_ , _Science_ 347 (2015) 1452 [1406.4088].
  * [13] W.I. Jay and E.T. Neil, _Bayesian model averaging for analysis of lattice field theory results_ , _Phys. Rev. D_ 103 (2021) 114502 [2008.01069].
  * [14] P.A. Boyle, G. Cossu, A. Yamaguchi and A. Portelli, _Grid: A next generation data parallel C++ QCD library_ , _PoS_ LATTICE2015 (2016) 023.
  * [15] A. Portelli, R. Abott, N. Asmussen, A. Barone, P.A. Boyle, F. Erben et al., _aportelli/hadrons: Hadrons v1.3_ , Mar., 2022.  10.5281/zenodo.6382460.
  * [16] RBC, UKQCD collaboration, _Domain wall QCD with physical quark masses_ , _Phys. Rev. D_ 93 (2016) 074505 [1411.7017].
  * [17] Y. Iwasaki, _Renormalization group analysis of lattice theories and improved lattice action: Two-dimensional non-linear O(N) sigma model_ , _Nuclear Physics B_ 258 (1985) 141.
  * [18] Budapest-Marseille-Wuppertal Collaboration collaboration, _Isospin Splittings in the Light-Baryon Octet from Lattice QCD and QED_ , _Phys. Rev. Lett._ 111 (2013) 252001.
  * [19] V. Cirigliano, G. Ecker, H. Neufeld, A. Pich and J. Portolés, _Kaon decays in the standard model_ , _Rev. Mod. Phys._ 84 (2012) 399.
  * [20] A. Desiderio et al., _First lattice calculation of radiative leptonic decay rates of pseudoscalar mesons_ , _Phys. Rev. D_ 103 (2021) 014502 [2006.05358].
  * [21] V. Cirigliano and H. Neufeld, _A note on isospin violation in $P\ell 2(\gamma)$ decays_, _Phys. Lett. B_ 700 (2011) 7 [1102.0563].