syndrome and consider the MMSE estimator
$\tilde{\bm{x}}_{d}\left(\bm{s}\right)=\bm{f}_{\rm MMSE}(\bm{s})$. Likewise,
let
${\bm{s}}^{\prime}=R_{\sqrt{2\pi}}(\bm{N}\bm{M}^{\top}\bm{\Omega}\bm{x}_{a})$
be the error syndrome in the new basis and the corresponding MMSE estimator
$\tilde{\bm{x}}_{d}^{\prime}\left(\bm{s}^{\prime}\right)=\bm{f}_{\rm
MMSE}(\bm{s}^{\prime})$. Then,
$\tilde{\bm{x}}_{d}^{\prime}\left(\bm{s}^{\prime}\right)=\tilde{\bm{x}}_{d}\left(\bm{s}\right)$.

###### Proof.

From Theorem 10, we have that
$P\left(\bm{x}_{d},\bm{s}\right)=P^{\prime}(\bm{x}_{d},\bm{s}^{\prime})$ and
thus $P\left(\bm{x}_{d}|\bm{s}\right)=P^{\prime}(\bm{x}_{d}|\bm{s}^{\prime})$
for the conditional distribution, such that

$\displaystyle\tilde{\bm{x}}_{d}^{\prime}(\bm{s}^{\prime})$
$\displaystyle=\int_{\mathbb{R}^{2N}}\differential{\bm{x}_{d}}\;\bm{x}_{d}P^{\prime}(\bm{x}_{d}|\bm{s}^{\prime})$
$\displaystyle=\int_{\mathbb{R}^{2N}}\differential{\bm{x}_{d}}\;\bm{x}_{d}P\left(\bm{x}_{d}|R_{\sqrt{2\pi}}\bm{s}\right)$
$\displaystyle=\tilde{\bm{x}}_{d}\left(\bm{s}\right),$

which was to be proved. ∎

######  Corollary 12.1.

The output covariance matrix $\bm{V}_{\rm out}$ describing the noise on the
data modes (i.e., the second moments of data error) is invariant under a basis
transformation on the ancillary lattice.

This corollary follows immediately from Theorem 12, but we explicitly derive
such for clarity.

###### Proof.

Consider the output covariance matrix in the new basis $\bm{V}_{\rm
out}^{\prime}$ [defined in a similar fashion as Eq. (97)] and the following
set of equalities,

$\displaystyle[\bm{V}_{\rm out}^{\prime}]_{ij}$
$\displaystyle=\int_{\mathbb{R}^{2N}}\differential{\bm{x}_{d}}\int_{\mathfrak{I}^{2M}}\differential{\bm{s}^{\prime}}\;\left[\bm{x}_{d}-\tilde{\bm{x}}_{d}^{\prime}\left(\bm{s}^{\prime}\right)\right]_{i\times
j}P^{\prime}(\bm{x}_{d},\bm{s}^{\prime})$
$\displaystyle=\int_{\mathbb{R}^{2N}}\differential{\bm{x}_{d}}\int_{\mathfrak{I}^{2M}}\differential{\bm{s}^{\prime}}\;\left[\bm{x}_{d}-\tilde{\bm{x}}_{d}\left(\bm{s}\right)\right]_{i\times
j}P(\bm{x}_{d},\bm{s})$
$\displaystyle=\int_{\mathbb{R}^{2N}}\differential{\bm{x}_{d}}\int_{\mathfrak{I}^{2M}}\differential{\left(R_{\sqrt{2\pi}}(\bm{N}\bm{s})\right)}\;\left(\bm{x}_{d}-\tilde{\bm{x}}_{d}\left(\bm{s}\right)\right)_{i\times
j}P(\bm{x}_{d},\bm{s})$
$\displaystyle=\int_{\mathbb{R}^{2N}}\differential{\bm{x}_{d}}\int_{\mathfrak{I}^{2M}}\differential{\bm{s}}\;\left(\bm{x}_{d}-\tilde{\bm{x}}_{d}\left(\bm{s}\right)\right)_{i\times
j}P(\bm{x}_{d},\bm{s})$ $\displaystyle=[\bm{V}_{\rm out}]_{ij}.$

We have used Theorems 10 and 12 to go from the first equality to the second
equality and used the relation
$\bm{s}^{\prime}=R_{\sqrt{2\pi}}\left(\bm{N}{\bm{s}}\right)$ and Lemma 11 to
go from the third equality to the fourth equality. The final equality follows
by definition of the covariance matrix. ∎

We now discuss linear estimation and show that such is not invariant under a
lattice basis transformation.

######  Theorem 13.

Linear estimation is not invariant under a lattice basis transformation.

###### Proof.

Considering only the $\bm{n}=\bm{0}$ element in the sum of Eq. (96), we get
the linear estimator

$\bm{f}_{\rm{Linear}}(\bm{s})=-\bm{V}_{d}^{-1}\bm{V}_{da}\bm{s},$ (116)

which is obviously not invariant under a basis transformation. Indeed,
consider a basis transformation such that
${\bm{s}}^{\prime}=R_{\sqrt{2\pi}}(\bm{N}\bm{s})$. The estimator in Eq. (116)
then changes to

$\displaystyle\bm{f}_{\rm{Linear}}(\bm{s}^{\prime})$
$\displaystyle=-\bm{V}_{d}^{-1}\bm{V}_{da}\bm{s}^{\prime}$
$\displaystyle=-\bm{V}_{d}^{-1}\bm{V}_{da}R_{\sqrt{2\pi}}(\bm{N}\bm{s}).$
(117)

and $\bm{f}_{\rm{Linear}}(\bm{s})\neq\bm{f}_{\rm{Linear}}(\bm{s}^{\prime})$. ∎

As a consequence, the output covariance matrix for linear estimation also
changes under a change of lattice basis. To show this explicitly, let us first
write the output covariance matrix $\bm{V}_{\rm out}$ in full,

$\displaystyle[\bm{V}_{\rm out}]_{ij}$
$\displaystyle=\int_{\mathbb{R}^{2N}}\differential{\bm{x}_{d}}\int_{\mathfrak{I}^{2M}}\differential{\bm{s}}\;{(\bm{x}_{d}-\bm{f}_{\rm{Linear}}(\bm{s}))}_{i\times
j}P(\bm{x}_{d},\bm{s})$
$\displaystyle=\sum_{\bm{k}}\int_{\mathfrak{I}^{2M}}\differential{\bm{s}}\int_{\mathbb{R}^{2N}}\differential{\bm{x}_{d}}\;{(\bm{x}_{d}-\bm{f}_{\rm{Linear}}(\bm{s}))}_{i\times
j}g(\bm{V}_{d}^{-1},\bm{x}_{d}+\bm{V}^{-1}_{d}\bm{V}_{da}(\bm{s}-\bm{k}\sqrt{2\pi}))g(\bm{V}_{d|a}^{-1},\bm{s}-\bm{k}\sqrt{2\pi})$
$\displaystyle=\sum_{\bm{k}}\int_{\mathfrak{I}^{2M}}\differential{\bm{s}}\;{[\bm{V}_{d}^{-1}]}_{ij}g(\bm{V}_{d|a}^{-1},\bm{s}-\bm{k}\sqrt{2\pi})+\left[-\bm{V}^{-1}_{d}\bm{V}_{da}\bm{s}+\bm{V}^{-1}_{d}\bm{V}_{da}(\bm{s}-\bm{k}\sqrt{2\pi})\right]_{i\times
j}g(\bm{V}_{d|a}^{-1},\bm{s}-\bm{k}\sqrt{2\pi})$
$\displaystyle=[\bm{V}_{d}^{-1}]_{ij}+\sum_{\bm{k}}2\pi\left[\bm{V}_{d}^{-1}\bm{V}_{da}\bm{k}\right]_{i\times
j}\int_{\mathfrak{I}^{2M}}\mathop{}\\!\mathrm{d}{\bm{s}}\;g(\bm{V}_{d|a}^{-1},\bm{s}-\bm{k}\sqrt{2\pi}),$
(118)

where we have used the expression of $P(\bm{x}_{d},\bm{s})$ in (93) and (98)
in the first and second steps, respectively. Under change of basis, the new
elements of the covariance matrix $\bm{V}_{\rm out}^{\prime}$ are

$\displaystyle[\bm{V}_{\rm out}^{\prime}]_{ij}$
$\displaystyle=\int_{\mathbb{R}^{2N}}\differential{\bm{x}_{d}}\int_{\mathfrak{I}^{2M}}\differential{\bm{s}^{\prime}}\;\left[\bm{x}_{d}-\bm{f}_{\rm{Linear}}(\bm{s}^{\prime})\right]_{i\times
j}P^{\prime}(\bm{x}_{d},\bm{s}^{\prime})$ (119)
$\displaystyle=\int_{\mathbb{R}^{2N}}\differential{\bm{x}_{d}}\int_{\mathfrak{I}^{2M}}\differential{\left(R_{\sqrt{2\pi}}\left(\bm{N}\bm{s}\right)\right)}\;[\bm{x}_{d}+\bm{V}_{d}^{-1}\bm{V}_{da}R_{\sqrt{2\pi}}(\bm{N}\bm{s})]_{i\times
j}P\left(\bm{x}_{d},\bm{s}\right)$ (120)
$\displaystyle=\int_{\mathbb{R}^{2N}}\differential{\bm{x}_{d}}\int_{\mathfrak{I}^{2M}}\differential{\bm{s}}\;[\bm{x}_{d}+\bm{V}_{d}^{-1}\bm{V}_{da}R_{\sqrt{2\pi}}(\bm{N}\bm{s})]_{i\times
j}\;P\left(\bm{x}_{d},\bm{s}\right)$ (121)
$\displaystyle=\sum_{\bm{k}}\int_{\mathbb{R}^{2N}}\differential{\bm{x}_{d}}\int_{\mathfrak{I}^{2M}}\differential{\bm{s}}\;[\bm{x}_{d}+\bm{V}_{d}^{-1}\bm{V}_{da}R_{\sqrt{2\pi}}(\bm{N}\bm{s})]_{i\times
j}\;g(\bm{V}_{d}^{-1},\bm{x}_{d}+\bm{V}^{-1}_{d}\bm{V}_{da}(\bm{s}-\bm{k}\sqrt{2\pi}))g(\bm{V}_{d|a}^{-1},\bm{s}-\bm{k}\sqrt{2\pi})$
(122)
$\displaystyle=\sum_{\bm{k}}\int_{\mathfrak{I}^{2M}}\differential{\bm{s}}\;g(\bm{V}_{d|a}^{-1},\bm{s}-\bm{k}\sqrt{2\pi})\left\\{[\bm{V}_{d}^{-1}]_{ij}+\left[\bm{V}_{d}^{-1}\bm{V}_{da}R_{\sqrt{2\pi}}(\bm{N}\bm{s})-\bm{V}^{-1}_{d}\bm{V}_{da}(\bm{s}-\bm{k}\sqrt{2\pi})\right]_{i\times
j}\right\\}$ (123)
$\displaystyle=[\bm{V}_{d}^{-1}]_{ij}+\sum_{\bm{k}}\int_{\mathfrak{I}^{2M}}\differential{\bm{s}}\;g(\bm{V}_{d|a}^{-1},\bm{s}-\bm{k}\sqrt{2\pi})\left[\bm{V}_{d}^{-1}\bm{V}_{da}R_{\sqrt{2\pi}}(\bm{N}\bm{s})-\bm{V}^{-1}_{d}\bm{V}_{da}(\bm{s}-\bm{k}\sqrt{2\pi})\right]_{i\times
j}$ (124)
$\displaystyle=[\bm{V}_{d}^{-1}]_{ij}+\sum_{\bm{k}}\int_{\mathfrak{I}^{2M}}\differential{\bm{s}}\;g(\bm{V}_{d|a}^{-1},\bm{s}-\bm{k}\sqrt{2\pi})\left[\bm{V}^{-1}_{d}\bm{V}_{da}\bm{k}\sqrt{2\pi}+\bm{V}_{d}^{-1}\bm{V}_{da}\left(R_{\sqrt{2\pi}}(\bm{N}\bm{s})-\bm{s}\right)\right]_{i\times
j}$ (125)
$\displaystyle=\underbrace{[\bm{V}_{d}^{-1}]_{ij}+\sum_{\bm{k}}2\pi\left[\bm{V}_{d}^{-1}\bm{V}_{da}\bm{k}\right]_{i\times
j}\int_{\mathfrak{I}^{2M}}\mathop{}\\!\mathrm{d}{\bm{s}}\;g(\bm{V}_{d|a}^{-1},\bm{s}-\bm{k}\sqrt{2\pi})}_{=[\bm{V}_{\rm
out}]_{ij}}$
$\displaystyle\qquad\qquad\;\;+\sum_{\bm{k}}\int_{\mathfrak{I}^{2M}}\mathop{}\\!\mathrm{d}{\bm{s}}\;g(\bm{V}_{d|a}^{-1},\bm{s}-\bm{k}\sqrt{2\pi})\left[\bm{V}_{d}^{-1}\bm{V}_{da}\left(R_{\sqrt{2\pi}}(\bm{N}\bm{s})-\bm{s}\right)\right]_{i\times
j}$
$\displaystyle\qquad\qquad\;\;+\sum_{\bm{k}}\sqrt{2\pi}\left[\bm{V}_{d}^{-1}\bm{V}_{da}\bm{k}\right]_{i}\int_{\mathfrak{I}^{2M}}\mathop{}\\!\mathrm{d}{\bm{s}}\;g(\bm{V}_{d|a}^{-1},\bm{s}-\bm{k}\sqrt{2\pi})\left[\bm{V}_{d}^{-1}\bm{V}_{da}\left(R_{\sqrt{2\pi}}(\bm{N}\bm{s})-\bm{s}\right)\right]_{j}$
$\displaystyle\qquad\qquad\;\;+\sum_{\bm{k}}\sqrt{2\pi}\left[\bm{V}_{d}^{-1}\bm{V}_{da}\bm{k}\right]_{j}\int_{\mathfrak{I}^{2M}}\mathop{}\\!\mathrm{d}{\bm{s}}\;g(\bm{V}_{d|a}^{-1},\bm{s}-\bm{k}\sqrt{2\pi})\left[\bm{V}_{d}^{-1}\bm{V}_{da}\left(R_{\sqrt{2\pi}}(\bm{N}\bm{s})-\bm{s}\right)\right]_{i}.$
(126)

From (119) to (120), we use
$\bm{s}^{\prime}=R_{\sqrt{2\pi}}\left(\bm{N}\bm{s}\right)$,
$P^{\prime}(\bm{x}_{d},\bm{s}^{\prime})=P(\bm{x}_{d},\bm{s})$ and the
expression in (117). From (120) to (121), we apply Lemma 11. To obtain (122),
we input the explicit equation of (93). The integral of $\bm{x}_{d}$ is
simplified with (98) from (122) to (123). From (123) to (125), we first
recognize that the $\bm{V}_{d}^{-1}$ part does not depend on $\bm{s}$ and then
reorganized the integrand. To obtain (125), we make use of the definition
$\bm{x}_{i\times j}=\bm{x}_{i}\bm{x}_{j}$, where $\bm{x}_{i}$ is the $i$th
component of $\bm{x}$, so that $(\bm{a}+\bm{b})_{i\times
j}=(\bm{a}+\bm{b})_{i}(\bm{a}+\bm{b})_{j}=\bm{a}_{i\times j}+\bm{b}_{i\times
j}+\bm{a}_{i}\bm{b}_{j}+\bm{a}_{j}\bm{b}_{i}$. Comparing equations (118) and
(126), we see that the linear estimation is not invariant under a lattice
basis transformation.

##  Appendix H Lower bound on output variance for multimode codes

We find a lower bound on the output variance for multimode GKP codes (with $N$
data modes and $M\geq N$ GKP ancilla modes), in a similar vein to the single-
mode lower bound found in Ref. [11]. We first provide an effective reduction
from multimode additive noise to independent single-mode additive noises which
we utilize to prove the result.

Given a $N$ mode input data state $\rho$, the output data state
$\rho^{\prime}$ after a round of error correction with a multimode GKP
stabilizer code (which protects against AGN) is related to the input by an
additive non-Gaussian noise channel $\widetilde{\mathcal{N}}$ via

$\rho^{\prime}=\widetilde{\mathcal{N}}\left(\rho\right)=\int_{\mathbb{R}^{2N}}\differential{\bm{e}_{d}}P(\bm{e}_{d})D_{\bm{e}_{d}}\rho
D_{\bm{e}_{d}}^{\dagger},$ (127)

where $P(\bm{e}_{d})$ is a multivariate non-Gaussian pdf that describes the
residual data noise $\bm{e}_{d}$—e.g., $\expectationvalue{\bm{e}_{d}}=0$ and
$\expectationvalue{\anticommutator{(\bm{e}_{d})_{i}}{(\bm{e}_{d})_{j}}}=(\bm{V}_{\rm
out})_{ij}$; for brevity in what follows, we let $\bm{V}\equiv\bm{V}_{\rm
out}$. Consider the Gaussian approximations $\rho_{G}$ and $\rho^{\prime}_{G}$
for the input and output states $\rho$ and $\rho^{\prime}$, which agree with
$\rho$ and $\rho^{\prime}$ at the level of first and second moments. We have
that

$\rho^{\prime}_{G}=\mathcal{N}_{\bm{V}}\left(\rho_{G}\right)=\int_{\mathbb{R}^{2N}}\differential{\bm{e}_{d}}g(\bm{V},\bm{e}_{d})D_{\bm{e}_{d}}\rho_{G}D_{\bm{e}_{d}}^{\dagger},$
(128)

where $g(\bm{V},\bm{e}_{d})$ is a Gaussian distribution according to Eq. (83)
with the same first and second moments as the non-Gaussian pdf $P(\bm{e}_{d})$
and $\mathcal{N}_{\bm{V}}$ is a multimode AGN channel with noise matrix
$\bm{V}$.

Let $\bm{S}$ be the symplectic transformation that diagonalizes $\bm{V}$ such
that $\bm{S}\bm{V}\bm{S}^{\top}=\bigoplus_{i=1}^{N}\nu_{i}\bm{I}_{2}$, where
$\nu_{i}\geq 0$ are the symplectic eigenvalues (single-mode variances) of
$\bm{V}$. We append Gaussian unitaries $U_{\bm{S}}^{\dagger}$ and $U_{\bm{S}}$
before and after error correction, respectively, and define the resulting
channel
$\widetilde{\mathcal{N}}^{(\bm{S})}\equiv\mathcal{U}_{\bm{S}}\circ\widetilde{\mathcal{N}}\circ\mathcal{U}_{\bm{S}}^{-1}$,
where $\mathcal{U}_{\bm{S}}$ is the unitary channel of $U_{\bm{S}}$. Then,

$\displaystyle\rho^{\prime}(\bm{S})$
$\displaystyle=\widetilde{\mathcal{N}}^{(\bm{S})}(\rho)$
$\displaystyle=\int_{\mathbb{R}^{2N}}\differential{\bm{e}_{d}}P(\bm{e}_{d})\left(U_{\bm{S}}D_{\bm{e}_{d}}U_{\bm{S}}^{\dagger}\right)\rho\left(U_{\bm{S}}D_{\bm{e}_{d}}^{\dagger}U_{\bm{S}}^{\dagger}\right)$
$\displaystyle=\int_{\mathbb{R}^{2N}}\differential{\bm{e}_{d}}P(\bm{e}_{d})D_{\bm{S}\bm{e}_{d}}\rho
D_{\bm{S}\bm{e}_{d}}^{\dagger}$
$\displaystyle=\int_{\mathbb{R}^{2N}}\differential{\bm{e}_{d}}P(\bm{S}^{-1}\bm{e}_{d})D_{\bm{e}_{d}}\rho
D_{\bm{e}_{d}}^{\dagger},$ (129)

where Eq. (26) was used in the second line and a change of variables was used
in the third line, along with the fact that
$\differential(\bm{S}^{-1}\bm{e}_{d})=\differential{\bm{e}_{d}}$. The Gaussian
approximation for $\rho^{\prime}(\bm{S})$ is then simply,

$\rho^{\prime}_{G}(\bm{S})=\mathcal{N}_{\bm{S}\bm{V}\bm{S}^{\top}}\left(\rho_{G}\right)=\bigotimes_{i=1}^{N}\mathcal{N}_{\nu_{i}}(\rho_{G}).$
(130)

This is due to the fact that displacements are always local,
$D_{\bm{e}_{d}}=\bigotimes_{i}D_{\bm{e}_{d_{i}}}$, and that
$P(\bm{S}^{-1}\bm{e}_{d})$ has a corresponding Gaussian distribution
satisfying

$\displaystyle g(\bm{V},\bm{S}^{-1}\bm{e}_{d})$
$\displaystyle=g(\bm{S}\bm{V}\bm{S}^{\top},\bm{e}_{d})$
$\displaystyle=g\Big{(}{\bigoplus}_{i}\nu_{i}\bm{I}_{2},\bm{e}_{d}\Big{)}$
$\displaystyle=\prod_{i}g(\nu_{i}\bm{I}_{2},\bm{e}_{d_{i}}),$ (131)

where $\bm{e}_{d_{i}}\in\mathbb{R}^{2}$. In other words, at the level of first
and second moments, we have decorrelated the residual noises on the data modes
via $\bm{S}$ and reduced the channel to a direct product of independent AGN
channels. We will use this simplification in what follows to prove a lower
bound on the quantum capacity of a general non-Gaussian additive noise channel
$\widetilde{\mathcal{N}}$, which includes the single-mode case proven in Ref.
[11].

######  Lemma 14.

Consider a $N$ mode additive non-Gaussian noise channel
$\widetilde{\mathcal{N}}$ [see Eq. (127)]. Let the geometric mean error of the
channel be defined as ${\bar{\sigma}_{\rm GM}^{2}=\sqrt[2N]{\det\bm{V}}}$,
where $\bm{V}$ is the $2N\times 2N$ covariance matrix of
$\widetilde{\mathcal{N}}$. Then the quantum capacity $C_{\mathcal{Q}}$ of
$\widetilde{\mathcal{N}}$ has the following lower bound,

$C_{\mathcal{Q}}(\widetilde{\mathcal{N}})\geq\max\left[0,N\log_{2}\left(\frac{1}{e\bar{\sigma}^{2}_{\rm
GM}}\right)\right].$ (132)

###### Proof.

To prove this result, we first consider the channel
$\widetilde{\mathcal{N}}^{(\bm{S})}$ of Eq. (129). Since
$\widetilde{\mathcal{N}}^{(\bm{S})}$ is related to $\widetilde{\mathcal{N}}$
by unitary pre- and post-processing, the quantum capacities of the two
channels are equivalent,
$C_{\mathcal{Q}}(\widetilde{\mathcal{N}})=C_{\mathcal{Q}}(\widetilde{\mathcal{N}}^{(\bm{S})})$.
We focus on the channel $\widetilde{\mathcal{N}}^{(\bm{S})}$ from hereon.

Observe that the single-shot coherent information places a lower bound on the
quantum capacity [41, 42, 43]

$C_{\mathcal{Q}}(\widetilde{\mathcal{N}}^{(\bm{S})})\geq\max_{\rho}I_{c}(\rho,\widetilde{\mathcal{N}}^{(\bm{S})})$
(133)

where
$I_{c}(\rho,\widetilde{\mathcal{N}}^{(\bm{S})})=S(\widetilde{\mathcal{N}}^{(\bm{S})}(\rho))-S(\widetilde{\mathcal{N}}^{(\bm{S})\,c}(\rho))$
is the coherent information, $S$ is the von Neumann entropy and
$\widetilde{\mathcal{N}}^{(\bm{S})\,c}$ is the complementary channel of
$\widetilde{\mathcal{N}}^{(\bm{S})}$. By Gaussian extremality [44], one can
show that $I_{c}(\rho,\widetilde{\mathcal{N}}^{(\bm{S})})\geq
S(\widetilde{\mathcal{N}}^{(\bm{S})}(\rho)_{G})-S(\widetilde{\mathcal{N}}^{(\bm{S})\,c}(\rho)_{G})$
where, e.g., $\widetilde{\mathcal{N}}^{(\bm{S})}(\rho)_{G}$ is a Gaussian
state with the same first and second moments as
$\widetilde{\mathcal{N}}^{(\bm{S})}(\rho)$. In fact,
$\widetilde{\mathcal{N}}^{(\bm{S})}(\rho)_{G}=\bigotimes_{i=1}^{N}\mathcal{N}_{\nu_{i}}(\rho_{G})$
from Eq. (130), which is just a direct product of AGN channels acting on an
input Gaussian state $\rho_{G}$; similarly,
$\widetilde{\mathcal{N}}^{(\bm{S})\,c}(\rho)_{G}=\bigotimes_{i=1}^{N}\mathcal{N}_{\nu_{i}}^{c}(\rho_{G})$.
Therefore, $I_{c}(\rho,\widetilde{\mathcal{N}}^{(\bm{S})})\geq
I_{c}(\rho_{G},\bigotimes_{i=1}^{N}\mathcal{N}_{\nu_{i}})$. In other words, we
have reduced the problem to analyzing a set of independent AGN channels with
Gaussian input states; this has been done in previous works for a single mode
with an input thermal state $\Theta_{\bar{n}}$ [45, 25]. In the energy
unconstrained setting ($\bar{n}\rightarrow\infty$), the single mode result is
$\max_{\rho_{G}}I_{c}(\rho_{G},\mathcal{N}_{\nu_{i}})=\log(1/e\nu_{i})$.
Hence, for the multimode channel $\bigotimes\mathcal{N}_{\nu_{i}}$,

$\displaystyle\max_{\rho_{G}}I_{c}\left(\rho_{G},\bigotimes_{i=1}^{N}\mathcal{N}_{\nu_{i}}\right)$
$\displaystyle\geq\sum_{i=1}^{N}\log(\frac{1}{e\nu_{i}})$ (134)
$\displaystyle=\log(\frac{1}{e^{N}\prod_{i=1}^{N}\nu_{i}}).$ (135)

Since $\nu_{i}$ are the symplectic eigenvalues of the channel covariance
matrix $\bm{V}$, $\det\bm{V}=\prod_{i=1}^{N}\nu_{i}^{2}$, and the result is
proved.  ∎

Now we consider the upper bound of the quantum capacity. We begin by quoting a
lemma in Refs. [25, 20].

######  Lemma 15.

Given an $N+M$ mode AGN channel, $\mathcal{N}_{\bm{Y}}$, with noise covariance
matrix $\bm{Y}\geq 0$, the quantum capacity $C_{\mathcal{Q}}$ then has the
following upper bound [25, 20],

$C_{\mathcal{Q}}\left(\mathcal{N}_{\bm{Y}}\right)\leq\sum_{i=1}^{N+M}\log_{2}\left(\frac{1-\sigma_{i}^{2}}{\sigma_{i}^{2}}\right),$
(136)

where $\sigma^{2}_{i}$ are the symplectic eigenvalues of $\bm{Y}$ (i.e., the
single-mode variances) .

From the above results of lower and upper bounds, we can have the following
theorem.

######  Theorem 16 (Lower bound for output noise of multimode codes).

Consider a multimode GKP stabilizer code which consumes $M\geq N$ ancillary
GKP modes to protect $N$ data modes against independent AGN. The joint noise
channel across all modes is given by
$\bigotimes_{i=1}^{N+M}\mathcal{N}_{\sigma_{i}^{2}}$ with variances
$\sigma_{i}^{2}$. Let the $2N\times 2N$ output noise matrix of the data be
$\bm{V}_{\rm out}$, with RMS and GM errors given as $\bar{\sigma}^{2}_{\rm
RMS}=\Tr{\bm{V}_{\rm out}}/2N$ and $\bar{\sigma}_{\rm
GM}^{2}=\sqrt[2N]{\det\bm{V}_{\rm out}}$, respectively. Then,

$\bar{\sigma}_{\rm RMS}\geq\bar{\sigma}_{\rm
GM}\geq\frac{1}{\sqrt{e}}\sqrt[2N]{\left(\prod_{i=1}^{N+M}\frac{\sigma_{i}^{2}}{1-\sigma_{i}^{2}}\right)}.$
(137)

###### Proof.

The proof for the geometric mean error follows by concurrently considering the
lower bound and upper bound of Lemmas 14 and 15, respectively. In detail, we
can use a multimode GKP code—and thus the resultant additive non-Gaussian
channel $\widetilde{\mathcal{N}}$—to transmit information across the AGN
channels $\bigotimes_{i=1}^{N+M}\mathcal{N}_{\sigma_{i}^{2}}$. The
transmission rate is bounded from below by the quantum capacity
$C_{\mathcal{Q}}(\widetilde{\mathcal{N}})$ since
$C_{\mathcal{Q}}(\widetilde{\mathcal{N}})$ is achievable, which in turn is
bounded from below by Lemma 14. Likewise, the transmission rate is bounded
from above by the quantum capacity of the original channel,
$C_{\mathcal{Q}}(\bigotimes_{i=1}^{N+M}\mathcal{N}_{\sigma_{i}^{2}})$, which
in turn is bounded from above by Lemma 15. This establishes the second
inequality in Eq. (137). The bound on the RMS error [the first inequality in
Eq. (137)] follows from the fact that the arithmetic mean is an upper bound on
the geometric mean. ∎

######  Corollary 16.1.

For iid AGN, such that $\sigma_{i}=\sigma\,\forall\,i\in\\{1,\dots,M+N\\}$,
$\bar{\sigma}_{\rm RMS}\geq\bar{\sigma}_{\rm
GM}\geq\frac{1}{\sqrt{e}}\left(\frac{\sigma^{2}}{1-\sigma^{2}}\right)^{\frac{N+M}{2N}}$.

Thus, for general multimode codes, if the number of ancilla modes is equal to
the number of data modes ($M=N$), then error suppression is at most quadratic
in the initial noise $\sigma$. The performance can be further enhanced with
concatenated codes ($M>N$), i.e. $\bar{\sigma}_{\rm
GM}\sim\sigma^{1+\frac{M}{N}}$. We find an upper bound on the break-even
point, $\sigma^{\star}\leq 1/\sqrt{2}$; see Fig. 5 of the main text. For
$\sigma>1/\sqrt{2}$, there is no gain to be had from QEC. This agrees with the
fact that the upper bound for the quantum capacity of the AGN channel [Eq.
(136)] vanishes as $\sigma\rightarrow 1/\sqrt{2}$.

##  Appendix I Classical entropies for the iid Gaussian channel

The differential entropy of multivariate normal distribution is [23]

$S(\bm{x})=\frac{1}{2}\ln\left[(2\pi e)^{K}\det\bm{V}\right],$ (138)

where $K$ is the dimension of the random variable $\bm{x}$. In the case of $N$
data modes and $M$ ancilla modes [$K=2(M+N)$], the correlated data and ancilla
noises (just prior to corrective displacements on the data) has a covariance
matrix

$\displaystyle\begin{split}\bm{V}_{\bm{x}}&=\bm{S}_{\rm{enc}}^{-1}\bm{V}_{\xi}\bm{S}_{\rm{enc}}^{-\top}\\\
&=(\bm{\Lambda}_{d}^{-1}\oplus\bm{\Lambda}_{a}^{-1})\left[\sigma^{2}\left(\bigoplus_{i=1}^{N}\bm{S}_{G_{i}}\bm{S}_{G_{i}}^{\top}\right)\oplus\bm{I}_{2(M-N)}\right](\bm{\Lambda}_{d}^{-\top}\oplus\bm{\Lambda}_{a}^{-\top}),\end{split}$
(139)

where $\bm{\Lambda}_{d}$ and $\bm{\Lambda}_{a}$ are local symplectic
transformations on the data and ancilla, respectively. The general form in
terms of two-mode squeezing blocks, for an arbitrary Gaussian encoding
$\bm{S}_{\rm enc}$, is a consquence of Lemma 1.

Since the local symplectic transforms do not change the determinants of the
subsystems, we can easily calculate all entropies of interest,

$\displaystyle S(\bm{x}_{d})=\sum_{i=1}^{N}\ln{[2\pi e(2G_{i}-1)\sigma^{2}]},$
(140) $\displaystyle S(\bm{x}_{a})=\sum_{i=1}^{N}\ln{[2\pi
e(2G_{i}-1)\sigma^{2}]}+(M-N)\ln(2\pi e\sigma^{2}),$ (141) $\displaystyle
S(\bm{x}_{d},\bm{x}_{a})=(N+M)\ln(2\pi e\sigma^{2}).$ (142)

The mutual information of the joint distribution is then

$\displaystyle I(\bm{x}_{d};\bm{x}_{a})$
$\displaystyle=S(\bm{x}_{d})+S(\bm{x}_{a})-S(\bm{x}_{d},\bm{x}_{a})$
$\displaystyle=2\sum_{i=1}^{N}\ln{(2G_{i}-1)},$ (143)

and the conditional differential entroy is

$\displaystyle S(\bm{x}_{d}|\bm{x}_{a})$
$\displaystyle=S(\bm{x}_{d},\bm{x}_{a})-S(\bm{x}_{a})$
$\displaystyle=\sum_{i=1}^{N}\ln{\left(\frac{2\pi
e\sigma^{2}}{2G_{i}-1}\right)}.$ (144)

In obtaining these results, we have utilized the fact that
$\bm{\Lambda}_{d}\oplus\bm{\Lambda}_{a}$ does not change any of the above
quantities. So in the calculations, we can actually consider the random
variables $\bm{\Lambda}_{d}\bm{x}_{d}$ and $\bm{\Lambda}_{a}\bm{x}_{a}$ for
simplicity. In terms of these variables, the data modes (and ancilla modes)
are uncorrelated amongst themselves, and there are only two-mode correlations
for each data-ancilla mode pair. Moreover, due to the structure of two-mode
squeezing, there are only $qq$ and $pp$ correlations and no $qp$ correlations.
Therefore, the total entropy, mutual information, and conditional entropy
originate from: (1) two-mode pairs (one data and one ancilla) due to two-mode
squeezing and (2) individual quadratures ($qq$ or $pp$ correlations but no
$qp$ correlations). As an example,
$S(\bm{x}_{d}|\bm{x}_{a})=S(\bm{\Lambda}_{d}\bm{x}_{d}|\bm{\Lambda}_{a}\bm{x}_{a})=\sum_{i=1}^{2N}S(\bm{\Lambda}_{d}\bm{x}_{d}|\bm{\Lambda}_{a}\bm{x}_{a})_{i}$,
where
$S(\bm{\Lambda}_{d}\bm{x}_{d}|\bm{\Lambda}_{a}\bm{x}_{a})_{i}=\ln{\left(\frac{2\pi
e\sigma^{2}}{2G_{i}-1}\right)}/2$ is the conditional entropy for the $qq$
($pp$) quadrature correlations of the $i$th data-ancilla mode pair. Everything
considered, we see that we can treat each data quadrature independently and
sum their individual contributions to calculate global properties.

##  Appendix J Error from approximated states

Given the noise model for finite GKP squeezing discussed in Section 9, it is
not hard to show that the error syndrome changes to

$\displaystyle\bm{s}$
$\displaystyle=R_{\sqrt{2\pi}}\left(\bm{M}^{\top}\bm{\Omega}(\bm{\xi}^{(1)}_{\rm{GKP}}+\bm{x}_{a})+\bm{\xi}^{(2)}_{\rm{GKP}}\right)$
(145)
$\displaystyle=R_{\sqrt{2\pi}}\left(\bm{M}^{\top}\bm{\Omega}\bm{x}_{a}+\bm{\delta}\right),$
(146)

where $\bm{\xi}^{(1)}_{\rm{GKP}}\sim\mathcal{N}(0,\sigma_{\rm GKP}^{2})$ is
the GKP noise from the ancilla modes,
$\bm{\xi}^{(2)}_{\rm{GKP}}\sim\mathcal{N}(0,\sigma_{\rm GKP}^{2})$ is the GKP
noise from the noisy stabilizer measurements, and
$\bm{\delta}=\bm{M}^{\top}\bm{\Omega}\bm{\xi}^{(1)}_{\rm{GKP}}+\bm{\xi}^{(2)}_{\rm{GKP}}$.
The PDF of $\bm{\delta}$ is

$P(\bm{\delta})=g\left(\sigma_{\rm{GKP}}^{2}(\bm{M}^{\top}\bm{M}+\bm{I}_{2M}),\bm{\delta}\right).$
(147)

We can show that the joint distribution of $\bm{x}_{d}$ and the error syndrome
$\bm{s}$ is

$\displaystyle P(\bm{x}_{d},\bm{s})$
$\displaystyle=\int_{\mathbb{R}^{4M}}\differential\left(\bm{M}^{\top}\bm{\Omega}\bm{x}_{a}\right)\differential\bm{\delta}\;\;P(\bm{x}_{d},\bm{M}^{\top}\bm{\Omega}\bm{x}_{a})P(\bm{\delta})\sum_{\bm{k}}\delta\left(\bm{s}-\bm{M}^{\top}\bm{\Omega}\bm{x}_{a}-\bm{\delta}-\bm{k}\sqrt{2\pi}\right)$
(148)
$\displaystyle=\sum_{\bm{k}}\int_{\mathbb{R}^{2M}}\differential\bm{\delta}\;\;g\left((\bm{I}_{2N}\oplus{\bm{M}^{\top}\bm{\Omega}})\bm{V}_{\bm{x}}(\bm{I}_{2N}\oplus{(\bm{M}^{\top}\bm{\Omega})}^{\top}),(\bm{x}_{d},\bm{s}-\bm{\delta}-\bm{k}\sqrt{2\pi})\right)g\left(\sigma_{\rm{GKP}}^{2}(\bm{M}^{\top}\bm{M}+\bm{I}_{2M}),\bm{\delta}\right)$
(149)
$\displaystyle=\sum_{\bm{k}}g\left((\bm{I}_{2N}\oplus{\bm{M}^{\top}\bm{\Omega}})\bm{V}_{\bm{x}}(\bm{I}_{2N}\oplus{(\bm{M}^{\top}\bm{\Omega})}^{\top}+(\bm{0}_{2N}\oplus\sigma_{\rm{GKP}}^{2}(\bm{M}^{\top}\bm{M}+\bm{I}_{2M}))),(\bm{x}_{d},\bm{s}-\bm{k}\sqrt{2\pi})\right)$
(150)
$\displaystyle\equiv\sum_{\bm{k}}g(\bm{V}_{d}^{-1},\bm{x}_{d}+\bm{V}^{-1}_{d}\bm{V}_{da}(\bm{s}-\bm{k}\sqrt{2\pi}))g(\bm{V}_{d|a}^{-1},\bm{s}-\bm{k}\sqrt{2\pi}),$
(151)

where we have redefined the notation $\bm{V}_{d},\bm{V}_{da},\bm{V}_{d|a}$ to
incorporate finite-squeezing GKP noise.

We use a moment-generating function to prove the identities (149) and (150).
Let $\bm{X}$ and $\bm{Y}$ be the vectors of independent random Gaussian
variables, then

$\displaystyle
P\left(\bm{X}=(\bm{x}_{d},\bm{s}-\bm{k}\sqrt{2\pi})\right)=g\left((\bm{I}_{2N}\oplus{\bm{M}^{\top}\bm{\Omega}})\bm{V}_{\bm{x}}(\bm{I}_{2N}\oplus{(\bm{M}^{\top}\bm{\Omega})}^{\top}),(\bm{x}_{d},\bm{s}-\bm{k}\sqrt{2\pi})\right),$
(152) $\displaystyle
P\left(\bm{Y}=(\bm{0}_{2N},\bm{\delta})\right)=g\left((\bm{0}_{2N}\oplus\sigma_{\rm{GKP}}^{2}(\bm{M}^{\top}\bm{M}+\bm{I}_{2M})),(\bm{0}_{2N},\bm{\delta})\right).$
(153)

The corresponding moment generating functions are

$\displaystyle
M_{\bm{X}}(\bm{t})=\exp{-\frac{1}{2}\bm{t}^{\top}(\bm{I}_{2N}\oplus{\bm{M}^{\top}\bm{\Omega}})\bm{V}_{\bm{x}}(\bm{I}_{2N}\oplus{(\bm{M}^{\top}\bm{\Omega})}^{\top})\bm{t}},$
(154) $\displaystyle
M_{\bm{Y}}(\bm{t})=\exp{-\frac{1}{2}\bm{t}^{\top}(\bm{0}_{2N}\oplus\sigma_{\rm{GKP}}^{2}(\bm{M}^{\top}\bm{M}+\bm{I}_{2M}))\bm{t}}.$
(155)

Let $\bm{Z}=\bm{X}+\bm{Y}$. Since $\bm{X}$ and $\bm{Y}$ are independent, then

$\displaystyle M_{\bm{Z}}(\bm{t})$
$\displaystyle=M_{\bm{X}}(\bm{t})M_{\bm{Y}}(\bm{t})$ (156)
$\displaystyle=\exp{-\frac{1}{2}\bm{t}^{\top}\left[(\bm{I}_{2N}\oplus{\bm{M}^{\top}\bm{\Omega}})\bm{V}_{\bm{x}}(\bm{I}_{2N}\oplus{(\bm{M}^{\top}\bm{\Omega})}^{\top})+(\bm{0}_{2N}\oplus\sigma_{\rm{GKP}}^{2}(\bm{M}^{\top}\bm{M}+\bm{I}_{2M}))\right]\bm{t}},$
(157)

and we have

$\displaystyle\eqref{eq:xd-s-joint-distribution-
step2}=P(\bm{Z}=\bm{X}+\bm{Y})=\eqref{eq:xd-s-joint-distribution-step3_appx}.$
(158)

This complete the proof.