# Optimal probes and error-correction schemes in multi-parameter quantum
metrology

Wojciech Górecki  Faculty of Physics, University of Warsaw, Pasteura 5, 02-093
Warsaw, Poland  Sisi Zhou  Departments of Applied Physics and Physics, Yale
University, New Haven, Connecticut 06511, USA  Yale Quantum Institute, Yale
University, New Haven, Connecticut 06511, USA  Pritzker School of Molecular
Engineering, University of Chicago, Chicago, IL 60637, USA  Liang Jiang
Departments of Applied Physics and Physics, Yale University, New Haven,
Connecticut 06511, USA  Yale Quantum Institute, Yale University, New Haven,
Connecticut 06511, USA  Pritzker School of Molecular Engineering, University
of Chicago, Chicago, IL 60637, USA  Rafał Demkowicz-Dobrzański  Faculty of
Physics, University of Warsaw, Pasteura 5, 02-093 Warsaw, Poland

###### Abstract

We derive a necessary and sufficient condition for the possibility of
achieving the Heisenberg scaling in general adaptive multi-parameter
estimation schemes in presence of Markovian noise. In situations where the
Heisenberg scaling is achievable, we provide a semidefinite program to
identify the optimal quantum error correcting (QEC) protocol that yields the
best estimation precision. We overcome the technical challenges associated
with potential incompatibility of the measurement optimally extracting
information on different parameters by utilizing the Holevo Cramér-Rao (HCR)
bound for pure states. We provide examples of significant advantages offered
by our joint-QEC protocols, that sense all the parameters utilizing a single
error-corrected subspace, over separate-QEC protocols where each parameter is
effectively sensed in a separate subspace.

††footnotetext: *These two authors provided key and equal contributions to the
project.

##  1 Introduction

Figure 1: General adaptive mutli-parameter quantum metrological scheme, where
$P$ parameters $\bm{\omega}=[\omega_{i}]_{i=1}^{P}$ are to be estimated. Total
probe system evolution time $T$ is divided into a number $m$ of $t$-long steps
of probe evolution $\mathcal{E}_{t}^{\bm{\omega}}$ interleaved with general
unitary controls $U_{i}$. In the end a general POVM $\\{M_{\ell}\\}$ is
performed yielding estimated value of all parameters
$\tilde{\bm{\omega}}(\ell)$ with probability
$p(\ell)=\textrm{Tr}(\rho_{\bm{\omega}}M_{\ell})$.

Quantum metrology aims at exploiting all possible features of quantum systems,
such as coherence or entanglement, in order to boost the precision of
measurements beyond that achievable by metrological schemes that operate
within classical or semi-classical paradigms [1, 2, 3, 4, 5, 6, 7, 8, 9]. The
most persuasive promise of quantum metrology is the possibility of obtaining
the so-called Heisenberg scaling (HS), which manifests itself in the
quadratically improved scaling of precision as a function of number of
elementary probe systems involved in the experiment [10, 11, 12, 13, 14, 15,
16, 17, 18, 19] or the total interrogation time of a probe system [20]. In
either of these cases, the presence of decoherence typically restricts the
quadratic improvement to a small particle number or a short-time regime,
whereas in the asymptotic regime the quantum-enhancement amounts to constant
factor improvements [21, 22, 23, 24, 25, 26] even under the most general
adaptive schemes [27]. Still, there are specific models where even in the
presence of decoherence the asymptotic HS is achievable via application of
appropriate quantum error correction (QEC) protocols [28, 29, 30, 31, 32, 33,
34, 35, 36, 37, 38, 39, 40, 41, 42].

Recently, a general theory providing a necessary and sufficient condition, the
HNLS condition (an acronym for “Hamiltonian-Not-in-Lindblad-Span”), for
achieving the HS in a finite-dimensional system in the most general adaptive
quantum metrological protocols under Markovian noise, has been developed [35,
36]. The theory allows for a quick identification of the most promising
quantum metrological models and provides a clear recipe for designing the
optimal adaptive schemes based on appropriately tailored QEC protocols.
However, HNLS is restricted to the single-parameter estimation case, while a
lot of relevant metrological problems, like vector field sensing (e.g.
magnetic field) [43], imaging [44], multiple-arm interferometry [45, 46] or
waveform estimation [47, 48] are inherently multi-parameter estimation
problems. Multi-parameter estimation problems drew a lot of attention in
recent years [49, 50, 51, 52, 53, 54, 55, 56], yet the fundamental questions
regarding the achievability of the HS as well as the form of the optimal
metrological protocols in multiple-parameter estimation in presence of noise
have not been answered so far. The aim of this paper is to fill this gap.

The methods developed for the single parameter estimation case, in particular
the semidefinite program (SDP) that allows to identify the optimal QEC
protocol [36], are not applicable in the multi-parameter regime. The reasons
are threefold.

First, the widely used quantum Cramér-Rao (CR) bound for multiple parameters
is not in general saturable, due to the incompatibility of the optimal
measurements for different parameters [51, 52, 49]. Therefore, unlike in the
single-parameter case, the quantum Fisher information (QFI) does not provide
the full insight into the problem [51, 57, 58, 59, 60]. On the other hand,
stronger bounds, such as the HCR bound [59, 61, 62, 63, 64, 60], have no
closed-form expressions (except for specific cases [65]). Moreover, the HCR
bound, although solvable via an SDP [66], does not shed light on the
corresponding optimal measurements saturating it. In general, the HCR bound is
only saturable when collective measurements on all copies of quantum states
are allowed, which is a demanding condition in practice [63, 64]. As a result,
the optimal measurements on the output quantum states are hard to identify.

Second, there is no general recipe to find the optimal input states in multi-
parameter estimation even in the noiseless case [52, 53], unlike in the
single-parameter estimation case where the optimal input state is simply the
equally weighted superposition between the eigenstates corresponding to the
maximum and minimum eigenvalues of the Hamiltonian.

Finally, in the single-parameter case [36], all valid two-dimensional QEC
codes were mapped into a traceless Hermitian matrix representing the
difference between logical zero and one codes. In the multi-parameter case,
however, it is not clear whether valid QEC codes could be mapped into a convex
set when the code dimension is large, which is inevitable in multi-parameter
estimation.

In this paper, we generalize the HNLS condition to multi-parameter scenarios,
and provide an SDP to find the best possible QEC protocol in the situations
when the HNLS condition is satisfied (including all noiseless cases). The
solution yields an explicit form of the optimal input state, QEC codes and
measurements. No collective measurements are required on the output states.
Our protocol goes beyond the typically used QFI-based formalism and overcomes
all the challenges related with the multi-parameter aspect of the problem
mentioned above. Our work reveals the advantage of QEC protocols in multi-
parameter estimation and we expect that the SDP formulation of our problem
will also be an inspiration for other research areas in quantum error
correction and quantum metrology.

##  2 Formulation of the model

We assume the dynamics of a $d$-dimensional probe system $\mathcal{H}_{S}$ is
given by a general quantum master equation [67, 68, 69]:

$\frac{d\rho}{dt}=-i[H,\rho]+\sum_{k=1}^{r}(L_{k}\rho
L_{k}^{\dagger}-\frac{1}{2}\\{L_{k}^{\dagger}L_{k},\rho\\}),$ (1)

where the parameters to be estimated
$\bm{\omega}=[\omega_{1},\dots,\omega_{P}]$ enter linearly into the
Hamiltonian of the evolution via Hermitian generators
$\bm{G}=[G_{1},\dots,G_{P}]^{T}$ (where T denotes transpose) so that
$H=\bm{\omega}\cdot\bm{G}\equiv\sum_{k=1}^{P}\omega_{k}G_{k}$, and $L_{k}$ are
operators representing a general Markovian noise. Similar to the previous
investigations [35, 36] we consider the most general adaptive scheme (see Fig.
1) [27] with an unlimited number of ancillae (denoted jointly as
$\mathcal{H}_{A}$), instantaneous perfect intermediate unitary operations
$U_{i}$ and a general POVM on the final state $\rho_{\bm{\omega}}$.
$\mathcal{E}_{t}^{\bm{\omega}}$ represents the probe system dynamics
integrated over time $t$, whereas the total probe interrogation time is $T$.
Such schemes are the most general schemes of probing quantum dynamics,
assuming the total interrogation time is $T$, and encompass in particular all
QEC procedures.

In single-parameter estimation the optimal protocol is the one that yields the
minimum estimation variance. In multi-parameter case the estimator covariance
matrix is the key object capturing estimation precision, defined as [58, 59]:

$\Sigma_{ij}=\textstyle{\sum}_{\ell}\,\mathrm{Tr}(\rho_{\bm{\omega}}M_{\ell})(\tilde{\omega}_{i}(\ell)-\omega_{i})(\tilde{\omega}_{j}(\ell)-\omega_{j})$
(2)

for $i,j=1,\ldots,P$, where $M_{\ell}\geq 0$,
$\sum_{\ell}M_{\ell}=\mathbbm{1}$, are measurement operators (“$\geq 0$” for
matrices means positive semidefinite) and
$\tilde{\text{\boldmath$\omega$}}(\ell)$ is an estimator function mapping a
measurement result $\ell$ to the parameter space.

Diagonal entries of $\Sigma$ represent variances of estimators of respective
parameters while off-diagonal terms represent covariance between different
parameters. As a figure of merit one may simply choose $\mathrm{Tr}(\Sigma)$
which will be the sum of all individual parameter variance, or more generally
$\mathrm{Tr}(W\Sigma)$, where $W$ is a real positive cost matrix that
determines the weight we associate with each parameter in the effective scalar
cost function

${\Delta_{W}^{2}\tilde{\bm{\omega}}}\equiv\mathrm{Tr}(W\Sigma).$ (3)

Note that we require strict positivity of $W$ which is equivalent to saying
that it is an estimation problem of all $P$ parameters, and not a problem
where effectively only a smaller number of parameters are relevant. We assume
the measurement-estimation strategy to be locally unbiased at some fixed
parameter point ${\bm{\omega}}$, i.e.

$\sum_{\ell}\tilde{\omega}_{j}(\ell)\mathrm{Tr}(\rho_{\text{\boldmath$\omega$}}M_{\ell})=\omega_{j},\quad\sum_{\ell}\tilde{\omega}_{j}(\ell)\mathrm{Tr}\big{(}\partial_{i}\rho_{\text{\boldmath$\omega$}}M_{\ell}\big{)}=\delta_{ij},$
(4)

where
$\partial_{i}\rho_{\text{\boldmath$\omega$}}=\frac{\partial\rho_{\text{\boldmath$\omega$}}}{\partial\omega_{i}}$,
which is a standard assumption necessary to obtain meaningful precision bounds
within the frequentist estimation framework [70, 71].

##  3 The necessary and sufficient condition for the HS

We say that the HS in a multi-parameter estimation problem is achieved when
there exists an adaptive protocol such that for every $W>0$,
${\Delta_{W}^{2}\tilde{\bm{\omega}}}\propto 1/T^{2}$ in the limit
$T\rightarrow\infty$. This is equivalent to a requirement that all parameters
(and any combination of parameters) are estimated with precision that scales
like the HS. The following theorem generalize the HNLS condition [35, 36] to
multi-parameter scenarios.

######  Theorem 1 (Multi-parameter HNLS).

The HS can be achieved in a multi-parameter estimation problem if and only if
$\\{(G_{i})_{\perp},i=1,\dots,P\\}$ are linearly independent operators. Here
$(G_{i})_{\perp}$ are orthogonal projections of $G_{i}$ onto space
$\mathcal{S}^{\perp}$ which is the orthogonal complement of the Lindblad span

$\mathcal{S}=\mathrm{span}_{\mathbb{R}}\\{\mathbbm{1},L_{k}^{{\mathrm{H}}},iL_{k}^{{\rm
AH}},(L_{k}^{\dagger}L_{j})^{{\mathrm{H}}},i(L_{k}^{\dagger}L_{j})^{{\rm
AH}},\,\forall j,k\\},$ (5)

in the Hilbert space of Hermitian matrices under the standard Hilbert-Schmidt
scalar product, whereas the superscripts H, AH denote the Hermitian and anti-
Hermitian part of an operator respectively.

###### Proof.

Let us start with a brief review of the HNLS condition in the single-parameter
case, where $H=\omega G$ involves only a single generator $G$. As shown in
[35, 36], the necessary and sufficient condition to achieve the HS is that
$G\notin\mathcal{S}$, or in other words that $G_{\perp}\neq 0$. In particular,
following [36] (see the section named “QEC code for HL scaling when HNLS
holds”), an explicit construction of the optimal QEC code was provided, where
the code space
$\mathcal{H}_{\mathcal{C}}\subseteq\mathcal{H}_{S}\otimes\mathcal{H}_{A}$ is
defined on the Hilbert space of the probe system $\mathcal{H}_{S}$ extended by
an ancillary space $\mathcal{H}_{A}\cong\mathcal{H}_{S}$. The code space
satisfies the QEC condition [36, 72]:

$\Pi_{\mathcal{H}_{\mathcal{C}}}(S\otimes\mathbbm{1})\Pi_{\mathcal{H}_{\mathcal{C}}}\propto\Pi_{\mathcal{H}_{\mathcal{C}}},\forall
S\in\mathcal{S},$ (6)

where the operator $S$ acting on $\mathcal{H}_{S}$ was tensored with identity
on $\mathcal{H}_{A}$ and $\Pi_{\mathcal{H}_{\mathcal{C}}}$ denotes the
projection onto $\mathcal{H}_{\mathcal{C}}$. Metrological sensitivity is
guaranteed by the fact that $G$ acts non-trivially on
$\mathcal{H}_{\mathcal{C}}$:

$G^{\mathcal{C}}=\Pi_{\mathcal{H}_{\mathcal{C}}}(G\otimes\mathbbm{1})\Pi_{\mathcal{H}_{\mathcal{C}}}\not\propto\Pi_{\mathcal{H}_{\mathcal{C}}}.$
(7)

As a result we obtain a noiseless unitary evolution generated by
$G^{\mathcal{C}}$ leading to the HS in the estimation precision of $\omega$.

(Necessity) Suppose $(G_{i})_{\perp}$’s are linearly dependent. Then there
exists a linear (invertible) transformation on the parameter space
$A\in\mathbb{R}^{P\times P}$: ${\bm{\omega}}^{\prime}={\bm{\omega}}A^{-1}$,
(where we also modify accordingly the generators $\bm{G}^{\prime}=A\bm{G}$ and
the cost matrix $W^{\prime}=AWA^{T}$, so that $H$ and
${\Delta_{W}^{2}\tilde{\bm{\omega}}}$ remain unchanged), such that
$(G^{\prime}_{i})_{\perp}=0$ for some $i$. Then, from the single-parameter
theorem, $\omega^{\prime}_{i}$ cannot be estimated with precision better than
$\Delta^{2}\tilde{\omega}^{\prime}_{i}\sim 1/T$ which contradicts the HS
requirements.

(Sufficiency) Suppose $(G_{i})_{\perp}$’s are linearly independent. We assume
the ancillary space to be a direct sum of $P$ subspaces $\mathcal{H}_{A_{i}}$
so that the whole Hilbert space is
$\mathcal{H}_{S}\otimes(\mathcal{H}_{A_{1}}\oplus\cdots\oplus\mathcal{H}_{A_{P}})$
(see Fig. 2a). We may construct separate code spaces for each parameter using
orthogonal ancillary subspace
$\mathcal{H}_{\mathcal{C}_{i}}\subseteq\mathcal{H}_{S}\otimes\mathcal{H}_{A_{i}}$
so that the QEC conditions Eq. (6) are satisfied within each code space
$\mathcal{H}_{\mathcal{C}_{i}}$ separately. While constructing the code space
for the $i$-th parameter, we include all the remaining generators $G_{j}$
($j\neq i$) in the Lindblad span, so effectively treating them as noise i.e.
$\mathcal{S}_{i}=\mathrm{span}_{\mathbb{R}}\\{\\{\mathbbm{1},L_{k}^{{\mathrm{H}}},iL_{k}^{{\rm
AH}},(L_{k}^{\dagger}L_{j})^{{\mathrm{H}}},i(L_{k}^{\dagger}L_{j})^{{\rm
AH}}\\}_{j,k}\cup\\{G_{j}\\}_{j\neq i}\\}$. As a result thanks to the QEC
condition it follows that $\forall_{i\neq
j}\Pi_{\mathcal{H}_{\mathcal{C}_{i}}}\left(G_{j}\otimes\mathbbm{1}\right)\Pi_{\mathcal{H}_{\mathcal{C}_{i}}}\propto\Pi_{\mathcal{H}_{\mathcal{C}_{i}}}$
and hence within a given subspace only one parameter is being sensed via the
effective generator
$G_{i}^{\mathcal{C}_{i}}=\Pi_{\mathcal{H}_{\mathcal{C}_{i}}}(G_{i}\otimes\mathbbm{1})\Pi_{\mathcal{H}_{\mathcal{C}_{i}}}$,
while all other generators act trivially. If
$\ket{\psi_{i}}\in\mathcal{H}_{\mathcal{C}_{i}}$ is the optimal state for
measuring $\omega_{i}$, the state to be used in order to obtain HS for all
parameters which is not affected by noise reads
$\rho_{\text{in}}=\frac{1}{P}\sum_{i=1}^{P}\ket{\psi_{i}}\bra{\psi_{i}}\in\mathcal{H}_{S}\otimes\left(\bigoplus_{i=1}^{P}\mathcal{H}_{A_{i}}\right)$—there
is no measurement incompatibility issue because different parameters are
encoded on orthogonal subspaces. ∎

Similar to the single-parameter case, it must be admitted that in realistic
situations with generic noise, HNLS is often violated [24, 35]. Therefore, a
more pragmatic approach is required taking into account the fact that in a
real experiment the total time of evolution $T$ is always finite. Let us
consider a situation where $H\in\mathcal{S}$, but where some noise components
are weak [36]. Specifically, we divide Lindblad operators in Eq. (1) into two
sets—strong noise generators $\\{L_{k}\\}$ and weak noise generators
$\\{J_{m}\\}$ satisfying $\epsilon:=\|\sum_{m}J^{\dagger}_{m}J_{m}\|$ where
$\|\cdot\|$ denotes operator norm. If the HNLS condition is satisfied for the
strong noise part, we could choose the code space $\mathcal{H}_{\mathcal{C}}$
which allows to completely erase the strong noise $\\{L_{k}\\}$ and the
resulting effective noise rate would be upper bounded by $\epsilon$ [36]. This
means that the distance between state of the error-corrected probe and the
state evolving under ideal noiseless evolution will be of the order
$\Theta(\epsilon T)$. Therefore, for sufficiently short evolution,
$T=o(1/\epsilon)$, the precision of estimation will still scale quadratically
with the total time
${\Delta_{W}^{2}\tilde{\bm{\omega}}}\propto\frac{1}{T^{2}}$ whereas for larger
$T$, it will gradually approach the standard $1/T$ scaling.

Figure 2: Schematic diagrams of relations between Hilbert spaces
$\mathcal{H}_{S},\mathcal{H}_{A},\mathcal{H}_{\mathcal{C}},\mathcal{H}_{A_{i}},\mathcal{H}_{\mathcal{C}_{i}},\mathcal{H}_{\mathcal{R}}$.
(a) In SEP-QEC, we use $P$ mutually orthogonal ancillary subspaces
$\mathcal{H}_{A_{i}}$ to sense each parameter $\omega_{i}$.
$\mathcal{H}_{A}=\bigoplus_{i=1}^{P}\mathcal{H}_{A_{i}}$ and
$\mathcal{H}_{\mathcal{C}}=\bigoplus_{i=1}^{P}\mathcal{H}_{\mathcal{C}_{i}}$.
$\dim(\mathcal{H}_{A_{i}})=\dim(\mathcal{H}_{S})=d$ and
$\dim(\mathcal{H}_{\mathcal{C}_{i}})=2$. (b) In JNT-QEC, we use a single code
space
$\mathcal{H}_{\mathcal{C}}\subseteq\mathcal{H}_{S}\otimes\mathcal{H}_{A}$ to
estimate all parameters jointly. $\dim(\mathcal{H}_{A})=(P+1)d$ and
$\dim(\mathcal{H}_{\mathcal{C}})=P+1$. (c) We use $\mathcal{H}_{\mathcal{R}}$
to represent the logical space ${\rm
span}\\{\ket{0},\ket{1},\ldots,\ket{P}\\}$ which is encoded into the physical
space $\mathcal{H}_{\mathcal{C}}={\rm
span}\\{\ket{c_{0}},\ket{c_{1}},\ldots,\ket{c_{P}}\\}$.

##  4 Optimal probes, error correction schemes and measurements

In Sec. 3, we provided a QEC code where each parameter is sensed separately in
different error-corrected subspaces (see Fig. 2a). Such protocols will be
referred as separate-QEC protocols (SEP-QEC). In contrast to this
construction, we will now consider QEC strategies which allow for simultaneous
estimation of all the parameters in a single coherent protocol by utilizing
states within a single protected code space, which we will call the joint-
parameter QEC scheme (JNT-QEC). In this section we provide a general method to
identify the optimal JNT-QEC, while its potential advantages over the optimal
SEP-QEC will be discussed in Sec. 5, as well as in Sec. 6 when studying
concrete estimation models.

From now on, we assume the multi-parameter HNLS condition is satisfied.
Without loss of generality, we assume the generators
$\\{G_{i}\\}_{i=1}^{P}\subset\mathcal{S}^{\perp}$ are orthonormal, since the
components in $\mathcal{S}$ do not contribute and there is always a linear
transformation $A$ on parameters leading to orthonormality. The following
theorem provides a recipe to find the optimal JNT-QEC.

######  Theorem 2 (Optimal JNT-QEC).

Given a cost matrix $W$. If the multi-parameter HNLS condition is satisfied
with generators $\\{G_{i}\\}^{P}_{i=1}\subset S^{\perp}$, the minimum cost
${\Delta_{W}^{2}\tilde{\bm{\omega}}}$ that can be achieved in a JNT-QEC reads

$\begin{split}&{\Delta_{W}^{2}\tilde{\bm{\omega}}}=\frac{P}{4T^{2}}\min_{G^{\mathcal{R}}_{i},B_{i},\nu_{i},K,w}w,\\\
\text{\leavevmode\nobreak\ subject to\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
}&\mathbbm{1}_{\scriptscriptstyle
P+1}\otimes\frac{\mathbbm{1}_{d}}{d}+\sum_{i=1}^{P}(G^{\mathcal{R}}_{i})^{T}\otimes
G_{i}+\sum_{i=P+1}^{P^{\prime}}\nu_{i}\mathbbm{1}_{\scriptscriptstyle
P+1}\otimes S_{i}+\sum_{i=P^{\prime}+1}^{d^{2}-1}B_{i}\otimes R_{i}\geq 0,\\\
&\Gamma_{ij}=\mathrm{Im}[G^{\mathcal{R}}_{j}]_{i0},\quad\begin{pmatrix}w\mathbbm{1}_{\scriptscriptstyle
P}&K\\\ K&\mathbbm{1}_{\scriptscriptstyle P}\end{pmatrix}\geq
0,\quad\begin{pmatrix}K&\mathbbm{1}_{\scriptscriptstyle P}\\\
\mathbbm{1}_{\scriptscriptstyle P}&\Gamma\sqrt{W^{-1}}\end{pmatrix}\geq
0,\end{split}$ (8)

where $\mathbbm{1}_{d}/\sqrt{d}$, $\\{G_{i}\\}_{i=1}^{P}$,
$\\{S_{i}\\}_{i=P+1}^{P^{\prime}}$, $\\{R_{i}\\}_{i=P^{\prime}+1}^{d^{2}-1}$
form an orthonormal basis of Hermitian operators acting on $\mathcal{H}_{S}$
such that
$\mathcal{S}=\mathrm{span}_{\mathbb{R}}\\{\mathbbm{1}_{d},{(S_{i})}_{i=P+1}^{P^{\prime}}\\}$.
Moreover, $G^{\mathcal{R}}_{i}$, $B_{i}$ are Hermitian $(P+1)\times(P+1)$
matrices (with matrix indices taking values from $0$ to $P$), and
$\nu_{i}\in\mathbb{R}$. $\Gamma$ and $K$ are real $P\times P$ matrices (with
matrix indices from $1$ to $P$).

As a semidefinite program (SDP) it could be easily solved numerically, for
example, using the Matlab-based package CVX [73]. Before giving a formal
proof, let us briefly review some existing bounds in multi-parameter metrology
and discuss their saturability.

###  4.1 General bounds in multi-parameter metrology

In this section, we discuss bounds for general parameter estimation problems
on fixed quantum states. We use $\mathcal{H}$ to denote a general finite-
dimensional Hilbert space.

Most commonly, quantum multi-parameter estimation problems are analyzed
utilizing the standard quantum CR bound [58, 57, 59]:

${\Delta_{W}^{2}\tilde{\bm{\omega}}}\geq\mathrm{Tr}(WF^{-1}),\quad
F_{ij}=(\mathrm{Tr}(\rho_{\bm{\omega}}\Lambda_{i}\Lambda_{j})),$ (9)

where $F$ is a $P\times P$ QFI matrix and $\Lambda_{i}$ (symmetric logarithmic
derivatives) satisfy
$\partial_{i}\rho_{\bm{\omega}}=\frac{1}{2}(\Lambda_{i}\rho_{\bm{\omega}}+\rho_{\bm{\omega}}\Lambda_{i})$.
This bound is not saturable in general, due to potential non-compatibility of
the optimal measurements, unless
$\mathrm{Im}[\mathrm{Tr}(\rho_{\bm{\omega}}\Lambda_{i}\Lambda_{j})]=0$ [51].
Moreover, a direct minimization of the CR bound over all JNT-QEC, with the
saturability constraint imposed, does not necessarily guarantee the
identification of the optimal protocol—there is a possibility that the optimal
protocol does not meet the saturability constraint for the CR bound.

Therefore, in order to identify truly optimal protocols we need to resort to a
stronger HCR bound [59, 61, 62]:

$\begin{split}&{\Delta_{W}^{2}\tilde{\bm{\omega}}}\geq\min_{\\{X_{i}\\}}(\mathrm{Tr}\left(W\cdot{\rm
Re}V)+\mathrm{Tr}\left[{\rm
abs}(W\cdot\mathrm{Im}V)\right]\right),\quad\text{where}\leavevmode\nobreak\
V_{ij}=\mathrm{Tr}(X_{i}X_{j}\rho_{{\bm{\omega}}}),\\\ &\text{for
Hermitian}\leavevmode\nobreak\ X_{i}\in{\cal
L}(\mathcal{H}),\quad\text{subject\leavevmode\nobreak\ to}\leavevmode\nobreak\
\mathrm{Tr}(X_{i}\partial_{j}\rho_{{\bm{\omega}}})=\delta_{ij},\end{split}$
(10)

where ${\cal L}(\circ)$ denotes the set of all linear operators acting on
$\circ$, ${\rm Re}$ and ${\rm Im}$ denote the real and imaginary part of a
matrix (not to be confused with the Hermitian and anti-Hermitian part of a
matrix), and $\mathrm{Tr}[{\rm abs}(\cdot)]$ is the sum of absolute values of
the eigenvalues of a matrix. When the second term is dropped the HCR bound
reduces to the standard CR bound [59, 51]. Unlike the CR bound this bound is
saturable in general using collective measurements on many copies [63, 64]. On
the other hand, the HCR is defined via an optimization problem, making it
usually more difficult to deal with than the standard quantum CR bound with a
closed-form expression.

In the case of pure states
$\rho_{{\bm{\omega}}}=\ket{\psi_{\bm{\omega}}}\bra{\psi_{\bm{\omega}}}$,
however, we note that the HCR is exactly equivalent to [49]:

$\begin{split}&{\Delta_{W}^{2}\tilde{\bm{\omega}}}\geq\min_{\\{\ket{x_{i}}\\}}\mathrm{Tr}(WV),\quad{\textrm{where}}\leavevmode\nobreak\
\leavevmode\nobreak\ V_{ij}=\braket{x_{i}}{x_{j}},\\\
&\text{for}\leavevmode\nobreak\
\ket{x_{i}}\in\mathrm{span}\\{\ket{\psi_{\bm{\omega}}},\partial_{1}\ket{\psi_{\bm{\omega}}},...,\partial_{P}\ket{\psi_{\bm{\omega}}}\\}\oplus\mathbb{C}^{P},\\\
&\text{subject\leavevmode\nobreak\ to}\leavevmode\nobreak\
2(\braket{x_{i}}{\partial_{j}}{\psi_{\bm{\omega}}})=\delta_{ij},\,\braket{x_{i}}{\psi_{\bm{\omega}}}=0,\;\mathrm{Im}(V)=0,\end{split}$
(11)

which we will call the Matsumoto bound. It was also shown in [49] that the
bound is always saturable using individual measurements (which relaxes the
requirement of collective measurements in the pure state case). However, there
are no known efficient numerical algorithms to find the solutions
$\\{\ket{x_{i}}\\}$ of the Mastumoto bound and it is not clear yet whether or
not the Mastumoto bound will be useful in finding the optimal QEC protocols in
our situation. In the following, we will start the proof of Theorem 2 by first
sketching the proof of the Matsumoto bound and then reformulating it in such a
way that the final optimization problem becomes an SDP, which will eventually
be incorporated into the QEC protocol optimization procedure.

###  4.2 Proof of Theorem 2

####  4.2.1 Proof and reformulation of the Matsumoto bound (Eq. (11))

According to the Naimark’s theorem [59], for any general measurement
$\\{M_{\ell}\\}$ on $\mathcal{H}$ there exists a projective measurement
$\\{E_{\ell}\\}$ on an extended space $\mathcal{H}_{M}$ (where
$\mathcal{H}\subseteq\mathcal{H}_{M}$) satisfying
$E_{\ell}E_{\ell^{\prime}}=\delta_{\ell\ell^{\prime}}E_{\ell}$ and
$M_{\ell}=\Pi_{\mathcal{H}}E_{\ell}\Pi_{\mathcal{H}}$. We now define a set of
vectors $\ket{x_{i}}\in\mathcal{H}_{M}$:

$\ket{x_{i}}=\sum_{\ell}(\tilde{\omega}_{i}(\ell)-\omega_{i})E_{\ell}\ket{\psi_{\bm{\omega}}}.$
(12)

One may see that, thanks to the projective nature of measurements
$\\{E_{\ell}\\}$, scalar products of vectors $\ket{x_{i}}$ yield the
covariance matrix of the estimator:

$V_{ij}=\braket{x_{i}}{x_{j}}=\sum_{\ell,\ell^{\prime}}\bra{\psi_{\bm{\omega}}}(\tilde{\omega}_{i}(\ell)-\omega_{i})E_{\ell}E_{\ell^{\prime}}(\tilde{\omega}_{j}(\ell^{\prime})-\omega_{j})\ket{\psi_{\bm{\omega}}}=\Sigma_{ij}.$
(13)

Now, instead of minimizing over the measurement $\\{M_{\ell}\\}$ on
$\mathcal{H}$, we can perform the minimization directly over the vectors
$\ket{x_{i}}\in\mathcal{H}_{M}$, imposing the following constraints:

$\mathrm{Im}(\braket{x_{i}}{x_{j}})=0,\quad\braket{x_{i}}{\psi_{\bm{\omega}}}=0,\quad
2(\braket{x_{i}}{\partial_{j}}{\psi_{\bm{\omega}}})=\delta_{ij}.$ (14)

These constraints correspond respectively to the projective nature of the
measurement $\\{E_{\ell}\\}$ and the unbiasedness conditions as given in Eq.
(4). At this point one may wonder how big the space $\mathcal{H}_{M}$ should
be (as for a general measurement it might be arbitrary large). However, we can
always map
$\mathrm{span}\\{\ket{\psi_{\bm{\omega}}},\\{\partial_{i}\ket{\psi_{\bm{\omega}}},\ket{x_{i}}\\}_{i=1}^{P}\\}\subseteq\mathcal{H}_{M}$
isometrically to a $(2P+1)$-dimensional space. Therefore when looking for the
bound, under the constraint Eq. (14), it is enough to perform the minimization
over
$\ket{x_{i}}\in\mathrm{span}\\{\ket{\psi_{\bm{\omega}}},\partial_{1}\ket{\psi_{\bm{\omega}}},...,\partial_{P}\ket{\psi_{\bm{\omega}}}\\}\oplus\mathbb{C}^{P}$,
which results in equation Eq. (11).

Finally, we show that indeed for any set of $\ket{x_{i}}$ satisfying Eq. (14)
there exists a proper projective measurement on
$\mathcal{H}\oplus\mathbb{C}^{P}$ and a locally unbiased estimator satisfying
Eq. (12), and consequently there exists a corresponding general measurement on
$\mathcal{H}$. To see this, notice that since
$\forall_{i}\braket{\psi_{\bm{\omega}}}{x_{i}}=0$ and
$\forall_{i,j}\braket{x_{i}}{x_{j}}\in\mathbb{R}$ one may choose an
orthonormal basis $\\{\ket{b_{i}}\\}$ of
$\textrm{span}\\{\ket{\psi_{\bm{\omega}}},\ket{x_{1}},\ldots,\ket{x_{P}}\\}$
satisfying:
$\forall_{i}\braket{\psi_{\bm{\omega}}}{b_{i}}\in\mathbb{R}\backslash\\{0\\}$
and $\forall_{i,j}\braket{x_{i}}{b_{j}}\in\mathbb{R}$. Then one can define a
projective measurement:

$E_{\ell}=\ket{b_{\ell}}\bra{b_{\ell}}\;(\ell=1,\ldots,P+1),\quad
E_{0}=\mathbbm{1}_{\textrm{dim}(\mathcal{H}_{M})}-\textstyle{\sum}_{\ell=1}^{P+1}\ket{b_{\ell}}\bra{b_{\ell}},$
(15)

with the corresponding estimator:

$\tilde{\omega}_{i}(\ell)=\frac{\braket{b_{\ell}}{x_{i}}}{\braket{b_{\ell}}{\psi_{\bm{\omega}}}}+\omega_{i},\;\ell\geq
1,\quad\tilde{\omega}_{i}(0)=0,$ (16)

which is locally unbiased and satisfies

$\ket{x_{i}}=\sum_{\ell=0}^{P+1}(\tilde{\omega}_{i}(\ell)-\omega_{i})E_{\ell}\ket{\psi_{\bm{\omega}}}.$
(17)

This proves Eq. (11).

Specifically, if $\dim(\mathcal{H})\geq 2P+1$ we may choose
$\mathrm{span}\\{\ket{\psi_{\bm{\omega}}},\partial_{1}\ket{\psi_{\bm{\omega}}},...,\partial_{P}\ket{\psi_{\bm{\omega}}}\\}\oplus\mathbb{C}^{P}$
as a subspace of $\mathcal{H}$ and optimize over $\ket{x_{i}}\in\mathcal{H}$.
In this case we may also reformulate the Matsumoto bound in a slightly
different form. First, note that any vectors $\\{\ket{x_{i}}\\}$ satisfying
Eq. (14) need to be linearly independent. Let $\\{\ket{c_{i}}\\}_{i=1}^{P}$ be
an orthonormal basis of $\mathrm{span}\\{\ket{x_{1}},...,\ket{x_{P}}\\}$,
satisfying $\forall_{i,j}\mathrm{Im}\braket{x_{i}}{c_{j}}=0$ (such a set may
be generated using the Gram-Schmidt orthonormalization procedure). The locally
unbiased conditions may now be rewritten as:

$2(\braket{x_{i}}{\partial_{j}}{\psi_{\bm{\omega}}})=\sum_{k=1}^{P}2(\braket{x_{i}}{c_{k}}\braket{c_{k}}{\partial_{j}}{\psi_{\bm{\omega}}})=\sum_{k=1}^{P}\braket{x_{i}}{c_{k}}2(\braket{c_{k}}{\partial_{j}}{\psi_{\bm{\omega}}})=\delta_{ij},$
(18)

which (after introducing matrices $\mathcal{X}_{ki}=\braket{c_{k}}{x_{i}}$,
$D_{kj}=2(\braket{c_{k}}{\partial_{j}}{\psi_{\bm{\omega}}})$ is equivalent to
the matrix equality $\mathcal{X}^{T}D=\mathbbm{1}_{\scriptscriptstyle P}$.
From $\mathcal{X}^{T}D=\mathbbm{1}_{\scriptscriptstyle P}$ we have
$\mathcal{X}^{T}=D^{-1}\Rightarrow\mathrm{Tr}(W\cdot
V)=\mathrm{Tr}(W\cdot\mathcal{X}^{T}\mathcal{X})=\mathrm{Tr}(W\cdot(D^{T}D)^{-1})$,
which gives

$\begin{split}&\min_{\ket{c_{1}},...,\ket{c_{P}}\in\mathcal{H}}\mathrm{Tr}(W\cdot(D^{T}D)^{-1})),\\\
&\text{where}\leavevmode\nobreak\
D_{ij}=2\braket{c_{i}}{\partial_{j}}{\psi_{\bm{\omega}}},\quad\text{subject
to}\leavevmode\nobreak\ \braket{c_{i}}{c_{j}}=\delta_{ij}.\end{split}$ (19)

This formulation will be more convenient to use when we will formulate the QEC
protocol optimization problem as an SDP.

####  4.2.2 Optimizing the error-correction codes

Now we apply the reformulated Matsumoto bound to our task of identification of
the optimal JNT-QEC. Consider a given input state $\ket{\psi_{\rm in}}$. Let
$\mathcal{H}_{\mathcal{C}}$ be any code subspace of
$\mathcal{H}_{S}\otimes\mathcal{H}_{A}$ containing $\ket{\psi_{\rm in}}$ and
satisfying the QEC conditions Eq. (6)—in order to be in accordance with the
reformulated Matsumoto bound, this space may be required to be at least $2P+1$
dimensional, but as we show in the following it will always be possible to
reduce its dimensionality to $P+1$ effectively. Using QEC, our goal is to
preserve an effective unitary evolution in the encoded space and coherently
acquire the sensing signal. Therefore, we are effectively dealing with pure
state $\ket{\psi_{{\bm{\omega}}}}$, which allows us to utilize Eq.
(LABEL:eq:orthomat) as a formula for the minimal cost of sensing multiple
parameters.

The effective evolution after implementing QEC is given by

$\ket{\psi_{\bm{\omega}}}=\exp\left(-iT\sum_{j=1}^{P}\omega_{j}\Pi_{\mathcal{H}_{\mathcal{C}}}\left(G_{j}\otimes\mathbbm{1}_{\textrm{dim}(\mathcal{H}_{A})}\right)\Pi_{\mathcal{H}_{\mathcal{C}}}\right)\ket{\psi_{\rm
in}}.$ (20)

We focus on the estimation around point ${\bm{\omega}}=[0,\ldots,0]$ (which
can always be achieved by applying inverse Hamiltonian dynamics [52]) and
denote $\ket{c_{0}}=\ket{\psi_{{\bm{\omega}}=0}}$ for notational simplicity.
Then for any $\ket{c_{i}}\in\mathcal{H}_{\mathcal{C}}$ we have
$2\braket{c_{i}}{\partial_{j}}{\psi_{{\bm{\omega}}=0}}=2T\mathrm{Im}\braket{c_{i}}{(G_{j}\otimes\mathbbm{1}_{\textrm{dim}(\mathcal{H}_{A})})}{c_{0}}$,
and according to Eq. (LABEL:eq:orthomat) the minimum achievable cost for a
fixed code space $\mathcal{H}_{\mathcal{C}}$ is given by:

$\begin{split}&\min_{\ket{c_{1}},...,\ket{c_{P}}\in\mathcal{H}_{\mathcal{C}}}\mathrm{Tr}(W\cdot(D^{T}D)^{-1})),\\\
&\text{where}\leavevmode\nobreak\
D_{ij}=2T\mathrm{Im}[\braket{c_{i}}{(G_{j}\otimes\mathbbm{1}_{\textrm{dim}(\mathcal{H}_{A})})}{\psi_{\bm{\omega}}}],\quad\text{subject
to}\leavevmode\nobreak\ \braket{c_{i}}{c_{j}}=\delta_{ij}.\end{split}$ (21)

From the above formulation it is clear that we may always reduce the code
space $\mathcal{H}_{\mathcal{C}}$ to
$\mathrm{span}\\{\ket{c_{k}}\\}_{k=0}^{P}$ without increasing the cost. Hence,
the problem of optimization over both probes and error-correction protocols is
now equivalent to identification of the set $\\{\ket{c_{k}}\\}_{k=0}^{P}$ that
minimizes the cost with the constraint that
$\mathcal{H}_{\mathcal{C}}=\mathrm{span}\\{\ket{c_{k}}\\}_{k=0}^{P}$ satisfies
the QEC conditions.

To solve this problem, it will be convenient to formally extend the Hilbert
space $\mathcal{H}_{S}\otimes\mathcal{H}_{A}$ by tensoring it with a
$(P+1)$-dimensional reference space
${\mathcal{H}_{\mathcal{R}}}=\textrm{span}\\{\ket{0},\dots,\ket{P}\\}$ (see
Fig. 2c). This reference space will be representing the effective evolution of
the probe state that happens within the code space and it will allow us to
encode QEC conditions in a compact and numerically friendly way.

First, we introduce a matrix $Q\in{\cal
L}({\mathcal{H}_{\mathcal{R}}}\otimes\mathcal{H}_{S})$ that represents a code

$Q=\mathrm{Tr}_{\mathcal{H}_{A}}\begin{pmatrix}\begin{bmatrix}\ket{c_{0}}\\\
\vdots\\\
\ket{c_{P}}\end{bmatrix}\begin{bmatrix}\bra{c_{0}}&\cdots&\bra{c_{P}}\end{bmatrix}\end{pmatrix},$
(22)

which, more formally, will be written as
$Q=\mathrm{Tr}_{\mathcal{H}_{A}}(\sum_{k,l=0}^{P}\ket{k}\bra{l}_{\mathcal{H}_{\mathcal{R}}}\otimes\ket{c_{k}}\bra{c_{l}}_{\mathcal{H}_{S}\otimes\mathcal{H}_{A}})$.
This matrix is proportional to the reduced density matrix of the maximum
entangled state between ${\mathcal{H}_{\mathcal{R}}}$ and
$\mathcal{H}_{\mathcal{C}}$. By its construction $Q\geq 0$ and contains all
relevant information on the code states in $\mathcal{H}_{\mathcal{C}}$.

Next, we introduce effective generators $G^{\mathcal{R}}_{i}$ acting on
${\mathcal{H}_{\mathcal{R}}}$ so that they represent properly the action of
the physical generators on the code space
$[G^{\mathcal{R}}_{i}]_{kl}=\braket{c_{k}}{G_{i}\otimes\mathbbm{1}_{\textrm{dim}(\mathcal{H}_{A})}}{c_{l}}$.
The effective evolution generators are related with the $Q$ matrix via:

$(G^{\mathcal{R}}_{i})^{T}=\mathrm{Tr}_{\mathcal{H}_{S}}\left[Q(\mathbbm{1}_{\scriptscriptstyle
P+1}\otimes G_{i})\right]\quad i=1,\ldots,P.$ (23)

Note that the identity operator here acts on the reference space
${\mathcal{H}_{\mathcal{R}}}$, and _not_ on the ancillary space
$\mathcal{H}_{A}$. Taking into account the orthonormality of $\ket{c_{k}}$ and
the QEC condition Eq. (6), we obtain the following constraints on $Q$

$\mathrm{Tr}_{\mathcal{H}_{S}}(Q)=\mathbbm{1}_{\scriptscriptstyle
P+1},\;\forall_{S_{i}\in\mathcal{S}}\leavevmode\nobreak\
\mathrm{Tr}_{\mathcal{H}_{S}}\left[Q(\mathbbm{1}_{\scriptscriptstyle
P+1}\otimes S_{i})\right]\propto\mathbbm{1}_{\scriptscriptstyle P+1}.$ (24)

Let $\mathbbm{1}_{d}/\sqrt{d}$, $\\{G_{i}\\}_{i=1}^{P}$,
$\\{S_{i}\\}_{i=P+1}^{P^{\prime}}$, $\\{R_{i}\\}_{i=P^{\prime}+1}^{d^{2}-1}$
form an orthonormal basis of Hermitian operators in ${\cal
L}(\mathcal{H}_{S})$ such that
$\mathcal{S}=\mathrm{span}_{\mathbb{R}}\\{\mathbbm{1}_{d},{(S_{i})}_{i=P+1}^{P^{\prime}}\\}$.
Any non-negative $Q$ satisfying Eqs. (23)-(24) has the following form:

$Q=\mathbbm{1}_{\scriptscriptstyle
P+1}\otimes\frac{\mathbbm{1}_{d}}{d}+\sum_{i=1}^{P}(G^{\mathcal{R}}_{i})^{T}\otimes
G_{i}+\sum_{i=P+1}^{P^{\prime}}\nu_{i}\mathbbm{1}_{\scriptscriptstyle
P+1}\otimes S_{i}+\sum_{i=P^{\prime}+1}^{d^{2}-1}B_{i}\otimes R_{i}\geq 0,$
(25)

where $\nu_{i}\in\mathbb{R}$ and $B_{i}$ are Hermitian. Conversely, for any
nonnegative defined $Q\geq 0$, we can consider its purification
$\ket{Q}\in{\mathcal{H}_{\mathcal{R}}}\otimes\mathcal{H}_{S}\otimes\mathcal{H}_{A}$,
which when written as
$\ket{Q}=\sum_{k=0}^{P}\ket{k}_{{\mathcal{H}_{\mathcal{R}}}}\otimes\ket{c_{k}}_{\mathcal{H}_{S}\otimes\mathcal{H}_{A}}$
yields the code states $\ket{c_{k}}$. Note that it implies that the rank of
$Q$ corresponding to the dimension of the ancillary space. It is always
sufficient to assume the dimension of the ancillary space to be
$\dim\mathcal{H}_{A}=(P+1)d$. Therefore $\\{G^{\mathcal{R}}_{i}\\}$ is an
achievable set of effective generators in ${\cal
L}({\mathcal{H}_{\mathcal{R}}})$ (satisfying the QEC condition) if and only if
there exist such $\nu_{i}\in\mathbb{R}$ and $B_{i}$, for which $Q\geq 0$.

Finally, in order to have an explicit dependence of the cost on the total time
parameter $T$, we introduce a matrix $\Gamma=\frac{1}{2T}D$, i.e.
$\Gamma_{ij}=\mathrm{Im}[\braket{c_{i}}{G_{j}\otimes\mathbbm{1}_{\textrm{dim}(\mathcal{H}_{A})}}{c_{0}}]=\mathrm{Im}[G^{\mathcal{R}}_{j}]_{i0}$,
and we end up with:

$\begin{split}&\frac{1}{4T^{2}}\min_{G^{\mathcal{R}}_{i},B_{i},\nu_{i}}\mathrm{Tr}\left(W(\Gamma^{T}\Gamma)^{-1}\right),\quad\text{where\leavevmode\nobreak\
}\Gamma_{ij}=\mathrm{Im}[G^{\mathcal{R}}_{j}]_{i0},\\\ &\text{subject
to}\leavevmode\nobreak\ \mathbbm{1}_{\scriptscriptstyle
P+1}\otimes\frac{\mathbbm{1}_{d}}{d}+\sum_{i=1}^{P}(G^{\mathcal{R}}_{i})^{T}\otimes
G_{i}+\sum_{i=P+1}^{P^{\prime}}\nu_{i}\mathbbm{1}_{\scriptscriptstyle
P+1}\otimes S_{i}+\sum_{i=P^{\prime}+1}^{d^{2}-1}B_{i}\otimes R_{i}\geq
0.\end{split}$ (26)

####  4.2.3 Reduction to an SDP

In order to reformulate Eq. (LABEL:almostend) as an SDP, we first show that we
may assume without loss of generality that $\Gamma\sqrt{W^{-1}}\geq 0$. Note
that for any full rank matrix $\Gamma$, the polar decomposition theorem
implies that there always exists an orthonormal matrix $O$ such that
$O\Gamma\sqrt{W^{-1}}\geq 0$. Next, as
$\Gamma_{ij}=\mathrm{Im}[\braket{i}{G^{\mathcal{R}}_{j}}{0}]$, multiplication
$\Gamma$ by $O$ is equivalent to rotating the base in the reference space
${\mathcal{H}_{\mathcal{R}}}$. Since, according to Eq. (22) such a rotation
cannot change the non-negativity of $Q$ and at the same time it does not
affect the figure of merit $\mathrm{Tr}\left(W(\Gamma^{T}\Gamma)^{-1}\right)$,
the statement is proven. To put Eq. (LABEL:almostend) in the form of an SDP,
we introduce a positive matrix $K\in{\mathbb{R}}^{P\times P}$ and a positive
real number $w$. Now, using the following two relations,

$\begin{bmatrix}K&\mathbbm{1}_{\scriptscriptstyle P}\\\
\mathbbm{1}_{\scriptscriptstyle P}&\Gamma\sqrt{W^{-1}}\\\ \end{bmatrix}\geq
0\leavevmode\nobreak\ \leavevmode\nobreak\ \Leftrightarrow\leavevmode\nobreak\
\leavevmode\nobreak\ K\geq(\Gamma\sqrt{W^{-1}})^{-1},$ (27)
$\begin{bmatrix}w\mathbbm{1}_{\scriptscriptstyle P}&K\\\
K&\mathbbm{1}_{\scriptscriptstyle P}\\\ \end{bmatrix}\geq
0\leavevmode\nobreak\ \leavevmode\nobreak\ \Leftrightarrow\leavevmode\nobreak\
\leavevmode\nobreak\ w\mathbbm{1}_{\scriptscriptstyle P}\geq K^{2},$ (28)

we see that $P\min
w=\min\mathrm{Tr}(K^{2})=\min\mathrm{Tr}\left(W(\Gamma^{T}\Gamma)^{-1}\right)$
in Eq. (LABEL:eq:thm2), making it equivalent to Eq. (LABEL:almostend). Hence
the problem takes the form of an SDP.

###  4.3 Discussion

It should be remarked that JNT-QEC do not contain SEP-QEC as a subclass. In
SEP-QEC, unlike in JNT-QEC, the noises are not fully corrected in the entire
space, and the decoherence is only avoided by choosing a properly mixed state
input. In general one could combine both these approaches in a unified
framework by dividing the set of all parameters into smaller subsets and then
applying JNT-QEC for each of these subset separately—in this approach the SEP-
QEC case would correspond to the situation where JNT-QEC optimization is
applied to single parameter subsets. Such an optimization is in principle
doable, but will involve much large numerical effort and it is not clear that
it will lead to better protocols.

It is also worth noting that, apart from the improved metrological performance
provided by QEC protocols when dealing with noisy systems, the above algorithm
is also applicable in the noiseless scenario when
$\mathcal{S}=\textrm{span}_{\mathbb{R}}\\{\mathbbm{1}_{d}\\}$. In such a
situation _no QEC is required_ (for simplicity we may still use
$\mathcal{H}_{\mathcal{C}}$ for $\mathrm{span}\\{\ket{c_{k}}\\}_{k=0}^{P}$,
but no recovery operation or projection $\Pi_{\mathcal{H}_{\mathcal{C}}}$ is
needed during evolution), but the condition ${\bm{\omega}}=[0,\ldots,0]$
(which is achievable by applying inverse Hamiltonian dynamics [52]) is still
required, as otherwise the derivatives of the state may not scale linearly
with $T$. In such situations, the solution of JNT-QEC yields an optimally
ancilla-assisted sensing protcol under arbitrary system dynamics
(Hamiltonitans) that resolves the potential incompatibility issues between
sensing of different parameters. It should be stressed that our approach is
universal and unlike existing approaches [43, 52, 53] does not assume any
specific structures of the Hamiltonians.

##  5 Advantages of JNT-QEC over SEP-QEC

In this section we investigate the potential advantages JNT-QEC provide
compared with SEP-QEC. The characteristic feature of the SEP-QEC is the
division resources so that each part of the resources is used to measure a
given parameter independently of the others. This is reflected in the form of
the input probe state
$\rho_{\textrm{in}}=\frac{1}{P}\sum_{i=1}^{P}\ket{\psi_{i}}\bra{\psi_{i}}$. As
a consequence we effectively measure each parameter only once in every $P$
repetitions of an experiment (corresponding to the $1/P$ factor in the
$\rho_{\textrm{in}}$). Therefore for a fixed total number of measurements, the
uncertainty of estimating a given parameter will grow proportionally to $P$.
Contrastingly, in JNT-QEC there is a chance to avoid the division of resources
and use a single pure state and a single measurement to estimate all of them
simultaneously. If there existed a code space, a state and a measurement which
were all simultaneously optimal for all parameters, then the decrease the cost
by a factor $P$ would be achievable. This would be the largest possible
advantage offered by JNT-QEC. In general we can write (see Appx. A for more
formal derivation):

${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm
JNT}\geq\tfrac{1}{P}{\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm SEP}.$ (29)

This inequality may be saturated only in special examples.

It is important to note that Theorem 2 gives us a recipe for identification of
the optimal JNT-QEC, while the explicit construction of SEP-QEC discussed in
the proof of Theorem 1 was only aimed at demonstrating the possibility of the
HS and hence the protocol was not optimized. Therefore, to make a fair
comparison between the two approaches, we need to compare the performance of
the best JNT-QEC with the best SEP-QEC.

Below we present a way to find the optimal SEP-QEC, i.e. the protocol for
which all parameters are measured independently on mutually orthogonal
subspaces. We will also provide a useful lower bound for the minimal cost
achievable in such protocols.

###  5.1 Optimization of SEP-QEC

Adapting the results from [36] (see the section named “Geometrical picture”),
we can infer that for a given set of $\\{G_{i}\\}$, the minimal variance
achievable in estimation of a single parameter $\omega_{i}$ in SEP-QEC is
given by:

$\Delta^{2}\tilde{\omega}_{i}=\frac{1}{F_{i}},\quad
F_{i}=4T^{2}\min_{\widetilde{G}_{i\parallel}\in\mathcal{S}\oplus\mathrm{span}_{\mathbb{R}}\\{G_{j}\\}_{j\neq
i}}\|G_{i}-\widetilde{G}_{i\parallel}\|,$ (30)

where $\|\cdot\|$ denotes operator norm. Let $\ket{\psi_{i}}$ be the optimal
state for measuring $\omega_{i}$. Therefore, using $\rho_{\rm
in}=\sum^{P}_{i=1}p_{i}\ket{\psi_{i}}\bra{\psi_{i}}$ (where in the naive
approach we would set all $p_{i}=\frac{1}{P}$) leads to the total cost
${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{{\rm
SEP}}=\sum_{i=1}^{P}\frac{1}{p_{i}}\frac{W_{ii}}{F_{i}}$. After optimization
over $p_{i}$ is performed (keeping in mind that $\sum_{i=1}^{P}p_{i}=1$,
$p_{i}\geq 0$) we get ${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{{\rm
SEP}}=\left(\sum_{i=1}^{P}\sqrt{\frac{W_{ii}}{F_{i}}}\right)^{2}$.

Next, we may still improve the peromance of SEP-QEC by choosing a different
QEC code based on a new set of generators obtained via a linear transformation
on the parameters $A\in{\mathbb{R}}^{P\times P}$:

${\bm{\omega}}^{\prime}={\bm{\omega}}A^{-1},\quad\bm{G}^{\prime}=A\bm{G},\quad
W^{\prime}=AWA^{T},$ (31)

so that the cost function remains unchanged.

Note also, that rescaling any generator by constant factor has no impact on
the results, so we may restrict to linear transformations satisfying
$(AA^{T})_{ii}=1$. Therefore the cost for the optimal SEP-QEC is given by:

$\begin{split}{\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm
SEP}&=\min_{A:\forall_{i}(AA^{T})_{ii}=1}\left(\sum_{i=1}^{P}\sqrt{\frac{(AWA^{T})_{ii}}{F_{i}(A)}}\right)^{2},\\\
&{\rm where}\quad
F_{i}(A)=4T^{2}\min_{\widetilde{G}_{i\parallel}\in\mathcal{S}\oplus\mathrm{span}_{\mathbb{R}}\\{(A\bm{G})_{j}\\}_{j\neq
i}}\|(A\bm{G})_{i}-\widetilde{G}_{i\parallel}\|.\end{split}$ (32)

We want to stress that introducing this (relatively complicated) procedure of
SEP-QEC optimization is necessary in order to distniguish the cases where the
true advantage is offered by the joint multi-parameter approach compared with
the situation where the advatage is only apparent and results from a
suboptimal choice of separate protocols.

Formula Eq. (32) is rather complicated and hard to compute in general.
However, in case $W=\mathbbm{1}$ the following lower bound is valid:

${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm SEP}\geq\min_{a_{i}:\sum
a_{i}^{2}=1}\frac{P^{2}}{T^{2}(\lambda^{*}_{+}-\lambda^{*}_{-})^{2}},\quad
G^{*}=\sum_{i=1}^{P}a_{i}G_{i},$ (33)

where $\lambda^{*}_{\pm}$ are extreme eigenvalues of $G^{*}$.

Since the reasoning leading to Eq. (29) holds for any set of generators, not
necessary the ones optimal for SEP-QEC, therefore, the optimal joint-
estimation cost, in case $W=\mathbbm{1}$, may be bounded as (see Appx. A for a
formal derivation):

${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{{\rm
JNT}}\geq\min_{i}\frac{P}{T^{2}(\lambda_{i+}-\lambda_{i-})^{2}},\quad{\rm
for\,any\,set\,of}\,G_{i}$ (34)

where $\lambda_{i\pm}$ are extreme eigenvalues of $G_{i}$.

In order to provide the reader with some intuition on the concept presented
above, we discuss below two extreme cases illustrating the apparent and
maximal advantage of joint-estimation protocols over separate ones. These
examples deal with noiseless scenarios where JNT-QEC is always optimal and are
aimed to prepare the reader for more physical noisy examples discussed in Sec.
6. Still, in order not to introduce additional abbreviations we will still use
the acronyms JNT-QEC, SEP-QEC to describe separate and joint estimation
schemes, even though the role of QEC in these examples is trivial.

###  5.2 Example: Apparent advantage of JNT-QEC

Consider a noiseless physical system comprising $P=2^{r}$ qubits (for
technical reasons we require the number of qubits to be the integer power of
$2$). Each qubit, regarded as a spin 1/2, senses the $z$-component of a local
magnetic field which is assumed to be independent for different qubit
locations. The corresponding sensing Hamiltonian for this system reads:

$H=\sum_{i=1}^{P}\omega_{i}G_{i},\quad
G_{i}=\sigma_{z}^{(i)}=\mathbbm{1}^{\otimes(i-1)}\otimes\sigma_{z}\otimes\mathbbm{1}^{\otimes(P-i)}$
(35)

and the cost matrix is assumed to be $W=\mathbbm{1}$. The minimum variance of
estimating each parameter independently is lower bounded by
$\Delta^{2}\tilde{\omega}_{i}\geq\frac{1}{T^{2}(\lambda_{+}-\lambda_{-})^{2}}$.
As the maximum and the minimum eigenvalues of each $G_{i}$ are
$\lambda_{i\pm}=\pm 1$, the minimal variance reads
$\Delta^{2}\tilde{\omega}_{i}=\frac{1}{4T^{2}}$. According to Eq. (34)) this
implies that the optimal JNT-QEC cost is lower bounded by $\frac{P}{4T^{2}}$.
Moreover, since the state
$\ket{\psi}=\tfrac{1}{\sqrt{2^{P}}}\left(\ket{+}+\ket{-}\right)^{\otimes P}$
is simultaneously optimal for all the parameters (and the tensor structure
guarantees no measurement incompatibility issue), this bound is saturable and
hence

${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{{\rm JNT}}=\frac{P}{4T^{2}}.$ (36)

On the other hand, the estimation strategy where each parameter $\omega_{i}$
is estimated separately leads to the cost equal $\frac{P^{2}}{4T^{2}}$.
Therefore, one may naively think that it is an example of superiority JNT-QEC
over SEP-QEC.

However, the separate protocol may be significantly improved here, by
estimating different combinations of the parameters. According to Eq. (33), we
should be looking for a proper linear combination of $G_{i}$ with the biggest
difference of extreme eigenvalues. The most obvious choice is
$G^{*}=\frac{1}{\sqrt{P}}\sum_{i=1}^{P}G_{i}$, for which
$\lambda^{*}_{\pm}=\sqrt{P}$, and therefore
${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm SEP}\geq\frac{P}{4T^{2}}$, which is
exactly equal to ${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{{\rm JNT}}$! Below we
show, that this bound may be saturated.

Let the matrix $A$ defining the transformation to the new set of parameters be
proportional to a Hadamard matrix of size $P\times P$ (i.e. square matrix
whose entries are either $+1$ or $-1$ and whose rows are mutually orthogonal)
defined as follows:

$A_{ij}=\frac{1}{\sqrt{P}}\prod_{k=0}^{r-1}(-1)^{i_{k}j_{k}},$ (37)

where the indices of the matrix are written using the binary representation,
$i=\sum_{k=0}^{r-1}i_{k}2^{k}$, where $i_{k}$ represent the binary digit of
$i$ at position $k$. Then for $\bm{G}^{\prime}=A\bm{G}$ we have:

$\lambda^{\prime}_{i\pm}=\pm\sqrt{P},\quad\ket{\lambda^{\prime}_{i\pm}}=\bigotimes_{j=1}^{P}\ket{\pm
A_{ij}},\quad\forall_{j\neq
i}\braket{\lambda^{\prime}_{i\pm}}{G^{\prime}_{j}}{\lambda^{\prime}_{i\pm}}=\braket{\lambda^{\prime}_{i\pm}}{G^{\prime}_{j}}{\lambda^{\prime}_{i\mp}}=0.$
(38)

As a result every $\omega_{i}^{\prime}$ may be measured separately on subspace
$\mathrm{span}\\{\ket{\lambda^{\prime}_{i+}},\ket{\lambda^{\prime}_{i-}}\\}$
with precision $\frac{1}{4PT^{2}}$. Therefore, for the initial sate $\rho_{\rm
in}=\frac{1}{P}\sum_{i=1}^{P}\ket{\psi_{i}}\bra{\psi_{i}}$ with
$\ket{\psi_{i}}=\tfrac{1}{\sqrt{2}}(\ket{\lambda_{i+}}+\ket{\lambda_{i-}})$ we
get

${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{{\rm SEP}}=\frac{P}{4T^{2}},$ (39)

which is exactly the same as the optimal
${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{{\rm JNT}}$. We see that the apparent
advantage of a JNT-QEC disappears once a proper combinations of parameters are
estimated in the SEP-QEC.

###  5.3 Example: Maximal advantage of JNT-QEC

As a contrasting example, here we discuss a situation where the advantage of
JNT-QEC is genuine and is the maximal possible. This example also provides an
intuition, what relation between the generators $G_{i}$ is responsible for
this.

Consider a noiseless system
$\mathcal{H}_{S}=\mathrm{span}\\{\ket{i}\\}_{i=0}^{P}$ with Hamiltonian:

$H=\sum_{i=1}^{P}\omega_{i}G_{i},\quad
G_{i}=\frac{1}{\sqrt{2}}(\ket{0}\bra{i}+\ket{i}\bra{0}).$ (40)

We focus on the estimation around point ${\bm{\omega}}=[0,\ldots,0]$ with the
cost matrix $W=\mathbbm{1}$. As all generators are orthonormal
$\mathrm{Tr}(G_{i}G_{j})=\delta_{ij}$, the difference between extreme
eigenvalues of any normalized combination of them, such as $G^{*}$ in Eq.
(33), is at most equal to $\sqrt{2}$. Therefore, the optimal SEP-QEC cost is
bounded by

${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm SEP}\geq\frac{P^{2}}{2T^{2}}.$ (41)

Importantly, the state $\ket{\psi_{\rm in}}=\ket{0}$ is simultaneously optimal
for measuring all the parameters. There is also no measurement incompatibility
problem as the optimal measurement for all the parameters is the measurement
in the basis $\ket{i}$. Therefore,

${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm JNT}=\frac{P}{2T^{2}}$ (42)

is achievable, and hence the there is a factor $P$ decrease in the cost in
case of JNT-QEC compared with the optimal SEP-QEC protocols.

##  6 Examples

Here we provide representative examples of a large class of multi-parameter
estimation models, where there are unavoidable tradeoffs in determining the
optimal states and measurements that arise due to the multi-parameter nature
of the problem. Our methods are useful in identifying optimal strategies in
all such models, provided the structure of noise admits the achievability of
the HS via application of the most general QEC schemes. By the construction of
our algorithm (Theorem 2), we have the guarantee that the solutions found are
the optimal ones. For an interested reader, a broader discussion and
generalizations of the results presented in this section, including proofs and
analytical constructions of the codes are provided in Appx. B, Appx. C, Appx.
D.

###  6.1 Single qubit case

Consider first the simplest single-qubit case with $d=2$. The HS is achievable
via QEC only in the case of single-rank Pauli noise (specified by a single
Hermitian Lindbladian $L$) [34]. Without loss of generality we can set
$L=\sigma_{z}$ (the Pauli-Z matrix). Since
$\mathcal{S}=\textrm{span}\\{\mathbbm{1},\sigma_{z}\\}$, at most two
parameters may be estimated in a qubit system with the HS (as
$\dim(\mathcal{S}^{\perp})=2$). However, it turns out that when the multi-
parameter HNLS condition is met there is no benefit in performing the more
sophisticated JNT-QEC compared to SEP-QEC, which is shown analytically in
Appx. B.

###  6.2 Two qubits in a magnetic field

In order to appreciate the superiority of JNT-QEC over SEP-QEC, let us
consider a two-qubit model which is a multi-parameter generalization of the
one from [37]. Consider two localized qubits, coupled to a magnetic field,
which is constant in both time and space, apart from some small fluctuations
in the $z$ direction. These fluctuations are assumed to be uncorrelated in
time, but maximally anticorrelated in space (for the two qubits they have
always opposite signs). Such a system may be effectively described by Eq. (1)
with $H=\frac{1}{2}\sum_{i=1}^{2}{\bm{\omega}}\cdot{\bm{\sigma}}^{(i)}$ (where
${\bm{\sigma}}^{(i)}=[\sigma_{x}^{(i)},\sigma_{y}^{(i)},\sigma_{z}^{(i)}]$
acts on the $i^{\text{th}}$ atom) and a single Lindblad operator
$L=\sqrt{2\gamma}(\sigma_{z}^{(1)}-\sigma_{z}^{(2)})$. It can be shown that
the minimal cost when each parameter is estimated with the optimal individual
parameter strategy is $\Delta^{2}\tilde{\omega}_{x,y,z}=\frac{1}{4T^{2}}$. At
the same time, in accordance with the discussion presented in Sec. 5, the best
precision achievable using SEP-QEC for the standard cost matrix
$W=\mathbbm{1}$, is
${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\textrm{SEP}}=\frac{P^{2}}{4T^{2}}=\frac{9}{4T^{2}}$—see
Appx. C for the formal derivation of the above formulas.

Below we present the result of numerical optimization of the JNT-QEC approach
(found by the algorithm presented in Theorem 2 and reconstructed to its
analytical form). We will use the standard Bell states notation:

$\begin{split}\ket{\Phi_{+}}=\frac{1}{\sqrt{2}}(\ket{\uparrow\uparrow}+\ket{\downarrow\downarrow}),\quad\ket{\Phi_{-}}=\frac{1}{\sqrt{2}}(\ket{\uparrow\uparrow}-\ket{\downarrow\downarrow}),\\\
\ket{\Psi_{+}}=\frac{1}{\sqrt{2}}(\ket{\uparrow\downarrow}+\ket{\downarrow\uparrow}),\quad\ket{\Psi_{-}}=\frac{1}{\sqrt{2}}(\ket{\uparrow\downarrow}-\ket{\downarrow\uparrow}).\end{split}$
(43)

Entanglement with ancilla will be abbreviated in the subscript
$\ket{\psi}\otimes\ket{i}_{A}\equiv\ket{\psi}_{i}$. Using the numerical
algorithm we have found out, that the optimal code space has the form

$\begin{split}\ket{c_{0}}&=-\cos(\varphi)\ket{\Phi_{+}}_{1}+\frac{i}{\sqrt{2}}\sin(\varphi)(\ket{\Phi_{+}}_{2}+\ket{\Phi_{-}}_{3}),\\\
\ket{c_{1}}&=-i\sin(\varphi)\ket{\Phi_{+}}_{1}-\cos(\varphi)\ket{\Psi_{+}}_{2},\\\
\ket{c_{2}}&=-\sin(\varphi)\ket{\Phi_{-}}_{1}-i\cos(\varphi)\ket{\Psi_{+}}_{3},\\\
\ket{c_{3}}&=-\frac{1}{\sqrt{2}}\sin(\varphi)(\ket{\Phi_{-}}_{2}+\ket{\Phi_{+}}_{3})+\cos(\varphi)\ket{\Psi_{+}}_{4},\\\
\end{split}$ (44)

where the input state is $\ket{\psi_{{\bm{\omega}}=\bm{0}}}=\ket{c_{0}}$. Note
that the presence of the last term in $\ket{c_{3}}$ (entangled with
$\ket{4}_{A}$) is necessary to satisfy the QEC conditions. The value of
$\varphi$ can be found analytically and the minimal total cost of estimation
${\Delta_{W}^{2}\tilde{\bm{\omega}}}$ is achieved for:

$\cos(\varphi)=\sqrt{\frac{\sqrt{7+4\sqrt{2}}-3}{4\sqrt{2}-2}}\approx 0.39,$
(45)

while the corresponding optimal cost is:

${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm JNT}\approx\frac{5.31}{4T^{2}}.$
(46)

As we see a significant improvement has been achieved here compared to SEP-
QEC.

Figure 3:  Numerical results for the optimal JNT-QEC strategy for estimating
the $SU(d)$ generators under the noise $J_{z}$ (black points joined by dashed
line), in contrast with the optimal precision asymptotically achievable by
SEP-QEC ${\Delta_{W}^{2}\tilde{\bm{\omega}}}=\frac{P^{2}}{2T^{2}}$ (black
solid line) and the lower bound
${\Delta_{W}^{2}\tilde{\bm{\omega}}}=\frac{P^{3/2}}{4T^{2}}$ asymptotically
achievable only in the noiseless case (solid gray line).

###  6.3 $SU(d)$ generators’ estimation

Finally let us consider an example of estimating parameters associated with
$SU(d)$ generators which shows an asymptotic advantage (with the number of
parameters) of the JNT-QEC over the SEP-QEC protocol. Let $\mathcal{H}_{S}$ be
$d$-dimensional Hilbert space, which for a more intuitive notation may be
regarded as the one associated with a spin-$j$ particle (where $d=2j+1$).

First, let us recall the noiseless case, where the Hamiltonian
$H=\sum_{i=1}^{d^{2}-1}\omega_{i}G_{i}$ is composed of all $P=d^{2}-1$ $SU(d)$
generators. The generators together with the identity
$\\{\frac{1}{\sqrt{d}}\mathbbm{1},G_{1},\ldots,G_{d^{2}-1}\\}$ form an
orthonormal basis of Hermitian operators on $\mathcal{H}_{S}$. We focus on the
estimation around point ${\bm{\omega}}=[0,\ldots,0]$ with the cost matrix
$W=\mathbbm{1}$. Since all $G_{i}$ are orthonormal the maximum spread between
their maximum and minimal eigenvalues as well as any normalized combination of
them is $\sqrt{2}$. Using the same line of reasoning as presented in Sec. 5.3
we conclude that

${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm
SEP}\geq\frac{P^{2}}{2T^{2}}=\frac{(d^{2}-1)^{2}}{2T^{2}}.$ (47)

The bound for ${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm JNT}$ may be derived
analytically. To achieve that, we use the following chain of inequalities
involving the QFI matrix:

$\sum_{i=1}^{d^{2}-1}\Delta^{2}\tilde{\omega}_{i}\geq\sum_{i=1}^{d^{2}-1}(F^{-1})_{ii}\geq\sum_{i=1}^{d^{2}-1}\frac{1}{F_{ii}}\geq\frac{(d^{2}-1)^{2}}{\sum_{i=1}^{d^{2}-1}F_{ii}},$
(48)

where the first one is the CR inequality and the rest are general algebraic
properties of positive semidefinite matrices. What remains to be done is to
derive a proper bound for the trace of the QFI matrix. For any input state
$\ket{\psi_{\rm in}}\in\mathcal{H}_{S}\otimes\mathcal{H}_{A}$ we have:

$F_{ii}=4T^{2}\left(\braket{\psi_{\bm{\omega}}}{G_{i}^{2}\otimes\mathbbm{1}}{\psi_{\bm{\omega}}}-\braket{\psi_{\bm{\omega}}}{G_{i}\otimes\mathbbm{1}}{\psi_{\bm{\omega}}}^{2}\right)\leq
4T^{2}\braket{\psi_{\bm{\omega}}}{G_{i}^{2}\otimes\mathbbm{1}}{\psi_{\bm{\omega}}}.$
(49)

Taking into account the normalization of $G_{i}$ and noting that
$\sum_{i=1}^{d^{2}-1}G_{i}^{2}$ is the Casimir operator of the $SU(d)$
algebra, and as such is proportional to the identity, we get that
$\sum_{i=1}^{d^{2}-1}G_{i}^{2}=\frac{d^{2}-1}{d}\mathbbm{1}$. Therefore:

$\sum_{i=1}^{d^{2}-1}F_{ii}\leq
4T^{2}\braket{\psi_{\bm{\omega}}}{\sum_{i=1}^{d^{2}-1}G_{i}^{2}\otimes\mathbbm{1}}{\psi_{\bm{\omega}}}=4T^{2}\frac{d^{2}-1}{d}.$
(50)

After substituting the above to Eq. (48) we get

${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm
JNT}=\sum_{i=1}^{d^{2}-1}\Delta^{2}\tilde{\omega}_{i}\geq\frac{d(d^{2}-1)}{4T^{2}}\approx\frac{P^{3/2}}{4T^{2}}.$
(51)

The example of a state which saturates the above bound is $\ket{\psi_{\rm
in}}=\frac{1}{\sqrt{d}}\sum_{k=1}^{d}\ket{k}_{S}\otimes\ket{k}_{A}$[53, 52].
For such a state the QFI matrix is given by
$F_{ij}=\delta_{ij}\frac{4T^{2}}{d}$, so the second and the third inequalities
in Eq. (48) become equalities. As
$\mathrm{Im}(\braket{\psi_{\bm{\omega}}}{\Lambda_{i}\Lambda_{j}}{\psi_{\bm{\omega}}})\propto\braket{\psi_{\bm{\omega}}}{[G_{i}\otimes\mathbbm{1},G_{j}\otimes\mathbbm{1}]}{\psi_{\bm{\omega}}}=0$,
the first one (the CR bound) is saturable as well. Note that the role of
ancillae here is to make optimal measurements with respect to different
parameters compatible. As a result, we see that optimal
${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm JNT}$ is approximately $2\sqrt{P}$
times smaller than ${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm SEP}$.

Let us now consider a noisy version with a single Lindblad operator
$J_{z}=\sum_{k=-j}^{j}k\ket{k}\bra{k}$. From Theorem 1 we know that only
parameters associated with generators
$G_{i}\notin\mathrm{span}_{\mathbb{R}}\\{\mathbbm{1},J_{z},J_{z}^{2}\\}$ may
be estimated with the HS. Therefore we consider the Hamiltonian
$H=\sum_{i=1}^{P}\omega_{i}G_{i}$ composed of $P=d^{2}-3$ of $SU(d)$
generators orthogonal to
$\mathrm{span}_{\mathbb{R}}\\{\mathbbm{1},J_{z},J_{z}^{2}\\}$), with the
standard cost matrix $W=\mathbbm{1}$ (such cost makes the problem independent
on choosing peculiar set of $\\{G_{i}\\}$). In Fig. 3, we present numerical
results for such a problem, and we observe a significant advantage over the
SEP-QEC protocol as well as strong indication of the asymptotic $P^{3/2}$
scaling identical to the noiseless case. Even though the optimal JNT-QEC code
cannot be written down analytically in a concise way, in Appx. D we provide an
analytical suboptimal construction achieving the $P^{3/2}$ scaling, supporting
the numerical results. The scaling advantage we prove here is not trivial
because there are no decoherence-free subspaces in the system.

##  7 Conclusions and outlook

We have generalized previous results on single-parameter error-corrected
metrology to the multi-parameter scheme, obtaining a necessary and sufficient
condition (HNLS) for the achievability of the HS. In case of scenarios when
HNLS is satisfied, we developed an efficient numerical algorithm (formulated
as an SDP) to find the optimal QEC protocol, including the optimal input
state, QEC codes and measurements. Our algorithm is applicable to arbitrary
system dynamics as long as the HS is achievable (including noiseless cases),
which contrasts previous works where special dynamics of quantum system is
assumed [43, 52, 53] or the estimation is performed on fixed quantum states
[66]. However, it still remains open in which situations the requirement of
noiseless ancillae could be removed [37, 38] and whether QEC is still helpful
in the case when HNLS is violated, as in the single-parameter case [74, 41].

We also remark that our way of formulating the Knill-Laflamme conditions [72]
as a positive semidefinite constraint Eq. (25) is novel and may have
applications beyond error-corrected quantum metrology.

Finally, we note that in this paper we have followed the frequentist
estimation approach, and in principle more stringent HS bounds might be
derived when following the Bayesian approach as was demonstrated recently in
the single parameter case [75].

## Acknowledgements

We thank Mengzhen Zhang, John Preskill and Francesco Albarelli for inspiring
discussions and Joseph Renes for his input in formulating the final form of
the SDP. WG and RDD acknowledge support from the National Science Center
(Poland) grant No. 2016/22/E/ST2/00559. SZ and LJ acknowledge support from the
ARL-CDQI (W911NF-15-2-0067, W911NF-18-2-0237), ARO (W911NF-18-1-0020,
W911NF-18-1-0212), ARO MURI (W911NF-16-1-0349), AFOSR MURI (FA9550-14-1-0052,
FA9550-15-1-0015), DOE (DE-SC0019406), NSF (EFMA-1640959), and the Packard
Foundation (2013-39273).

##  Appendix A The proof of ${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm
JNT}\geq\frac{1}{P}{\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm SEP}$

In the reasoning presented in the beginning of Sec. 5 we have assumed that,
when estimating a single parameter, the minimal variance of the estimator
$\Delta^{2}\tilde{\omega}_{i}$ achievable in SEP-QEC is no bigger than the one
achievable in JNT-QEC (even if we focus only on this particular parameter).
This statement is not so obvious as both protocols impose different
constraints on the code space:

$\begin{split}{\textrm{JNT-
QEC}}&:\,\forall_{j}\Pi_{\mathcal{H}_{C}}(G_{j}\otimes\mathbbm{1})\Pi_{\mathcal{H}_{C}}\not\propto\Pi_{\mathcal{H}_{C}},\\\
{\textrm{SEP-QEC}}&:\,\forall_{j\neq
i}\Pi_{\mathcal{H}_{C_{i}}}(G_{j}\otimes\mathbbm{1})\Pi_{\mathcal{H}_{C_{i}}}\propto\Pi_{\mathcal{H}_{C_{i}}},\quad\Pi_{\mathcal{H}_{C_{i}}}(G_{i}\otimes\mathbbm{1})\Pi_{\mathcal{H}_{C_{i}}}\not\propto\Pi_{\mathcal{H}_{C_{i}}},\end{split}$
(52)

where none of them is weaker or stronger than the other one.

Therefore, here we show directly that having JNT-QEC with accuracy
${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm JNT}$, one may always construct SEP-
QEC giving ${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm
SEP}=P\cdot{\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm JNT}$. Let
$\mathcal{H}_{\mathcal{C}}$ be the code space used in JNT-QEC. Without loss of
generality let us assume that $W$ is diagonal (otherwise we apply a
transformation of parameters that diagonalizes $W$) and
$\dim(\mathcal{H}_{\mathcal{C}})\geq 2P+1$ (otherwise we trivially extend it
using an additional ancilla). From the Matsumoto bound (see Eq. (11)) the
optimal cost is given by:

$\begin{split}&{\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm
JNT}=\min_{\\{\ket{x_{i}}\\}}\mathrm{Tr}(WV),\quad{\textrm{where}}\leavevmode\nobreak\
\leavevmode\nobreak\ V_{ij}=\braket{x_{i}}{x_{j}},\\\
&\text{for}\leavevmode\nobreak\ \ket{x_{i}}\in\mathcal{H}_{\mathcal{C}},\\\
&\text{subject\leavevmode\nobreak\ to}\leavevmode\nobreak\
2(\braket{x_{i}}{\partial_{j}}{\psi_{\bm{\omega}}})=\delta_{ij},\,\braket{x_{i}}{\psi_{\bm{\omega}}}=0,\;\mathrm{Im}(V)=0,\end{split}$
(53)

where for $\omega=0$ we have
$\partial_{j}\ket{\psi_{\bm{\omega}}}\big{|}_{\omega=0}=-iTG_{j}^{\mathcal{C}}\ket{\psi_{{\bm{\omega}}=\bm{0}}}$.
We define
$\mathcal{H}_{\mathcal{C}_{i}}=\mathrm{span}_{\mathbb{R}}\\{\ket{c_{0}^{i}},\ket{c_{1}^{i}}\\}$
with:

$\ket{c_{0}^{i}}=\tfrac{1}{\sqrt{2}}\left(\ket{\psi_{{\bm{\omega}}=\bm{0}}}\ket{0}_{A_{i}}+\tfrac{1}{\sqrt{\braket{x_{i}}{x_{i}}}}\ket{x_{i}}\ket{1}_{A_{i}}\right),\quad\ket{c_{1}^{i}}=\tfrac{1}{\sqrt{2}}\left(\tfrac{1}{\sqrt{\braket{x_{i}}{x_{i}}}}\ket{x_{i}}\ket{0}_{A_{i}}+\ket{\psi_{{\bm{\omega}}=\bm{0}}}\ket{1}_{A_{i}}\right),$
(54)

where $\ket{x_{i}}$ are the result of optimization Eq. (53). As
$\ket{\psi_{{\bm{\omega}}=\bm{0}}},\ket{x_{i}}\in\mathcal{H}_{\mathcal{C}}$,
obviously the noise acts trivially on $\mathcal{H}_{\mathcal{C}_{i}}$.
Moreover:

$\begin{split}\bra{c_{0}^{i}}G_{j}\otimes\mathbbm{1}\ket{c_{1}^{i}}&=\tfrac{1}{2\sqrt{\braket{x_{i}}{x_{i}}}}\left(\bra{\psi_{{\bm{\omega}}=\bm{0}}}G_{j}^{\mathcal{C}}\ket{x_{i}}+\bra{x_{i}}G_{j}^{\mathcal{C}}\ket{\psi_{{\bm{\omega}}=\bm{0}}}\right)=\tfrac{1}{T\sqrt{\braket{x_{i}}{x_{i}}}}\delta_{ij}\\\
\bra{c_{0}^{i}}G_{j}\otimes\mathbbm{1}\ket{c_{0}^{i}}&=\tfrac{1}{2}\left(\bra{\psi_{{\bm{\omega}}=\bm{0}}}G_{j}^{\mathcal{C}}\ket{\psi_{{\bm{\omega}}=\bm{0}}}+\tfrac{1}{\braket{x_{i}}{x_{i}}}\bra{x_{i}}G_{j}^{\mathcal{C}}\ket{x_{i}}\right)=\bra{c_{1}^{i}}G_{j}\otimes\mathbbm{1}\ket{c_{1}^{i}},\end{split}$
(55)

so also condition $\forall_{j\neq
i}\Pi_{\mathcal{H}_{\mathcal{C}_{i}}}G_{j}\Pi_{\mathcal{H}_{\mathcal{C}_{i}}}\propto\Pi_{\mathcal{H}_{\mathcal{C}_{i}}}$
is satisfied.

Therefore SEP-QEC with initial state $\rho_{\rm
in}=\frac{1}{P}\sum_{i=1}^{P}\ket{c_{0}^{i}}\bra{c_{0}^{i}}$ and
$\mathcal{H}_{\mathcal{C}_{i}}$ defined above leads to the total cost
${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm SEP}=\sum_{i=1}^{P}P\cdot
W_{ii}\braket{x_{i}}{x_{i}}=P\cdot{\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm
JNT}$. From this reasoning we see clearly that the largest possible advantage
of JNT-QEC over SEP-QEC is to decreasing the total cost by a factor $P$.

##  Appendix B Optimality of SEP-QEC in the single-qubit model

Let us consider the most general case of a two-parameter quantum estimation
problem under Markovian noise in a $2$-dimensional Hilbert space when the HS
achievability condition is satisfied. Without loss of generality we assume
that the Lindblad operator is $L=\sigma_{z}$, the Hamiltonian is
$H=\omega_{x}\sigma_{x}+\omega_{y}\sigma_{y}$ and the cost matrix is diagonal:

$\begin{bmatrix}W_{xx}&0\\\ 0&W_{yy}\end{bmatrix}$ (56)

(otherwise, one can always apply a proper transformation in the parameter
space ${\bm{\omega}}^{\prime}={\bm{\omega}}A^{-1}$ which diagonalizes the cost
matrix, without changing orthonormality of the generators). Below we show that
for such a problem there is no advantage of JNT-QEC over SEP-QEC.

First let us consider the optimal SEP-QEC. Each generator has eigenvalues
$\ket{\lambda_{x/y\pm}}=\pm 1$ and each parameter may be estimated with
precision $\Delta^{2}\tilde{\omega}_{x/y}=\frac{1}{4T^{2}}$ [34]. Moreover, as
the cost is diagonal, there is no point in applying an additional
transformation in the optimization procedure given in Eq. (32)—indeed, any
matrix $A$ satisfying $\forall_{i}(AA^{T})_{ii}=1$ preserves the eigenvalues
of the generators as well as $\mathrm{Tr}(AWA^{T})$. Therefore the optimal
SEP-QEC cost is simply given by ${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm
SEP}=\frac{(\sqrt{W_{xx}}+\sqrt{W_{yy}})^{2}}{4T^{2}}$. We will show below
that this cannot be outperformed by the optimal JNT-QECs.

For the sake of notation simplicity, when tensoring an operator acting on
$\mathcal{H}_{S}$ with the identity on $\mathcal{H}_{A}$, the part
$\otimes\mathbbm{1}$ will be omitted and we will denote
$\sigma_{i}\otimes\mathbbm{1}$ simply as $\sigma_{i}$ (unless this leads to
ambiguity).

First, we note that the diagonal elements of the QFI matrix for the state
$\ket{\psi_{\bm{\omega}}}$ are

$\begin{split}F_{ii}&=4T^{2}(\braket{\psi_{\bm{\omega}}}{\sigma_{i}\Pi_{\mathcal{H}_{\mathcal{C}}}\sigma_{i}}{\psi_{\bm{\omega}}}-\left|{\braket{\psi_{\bm{\omega}}}{\Pi_{\mathcal{H}_{\mathcal{C}}}\sigma_{i}}{\psi_{\bm{\omega}}}}\right|^{2})\\\
&\leq
4T^{2}\braket{\psi_{\bm{\omega}}}{\sigma_{i}\Pi_{\mathcal{H}_{\mathcal{C}}}\sigma_{i}}{\psi_{\bm{\omega}}}\quad(i=x,y).\end{split}$
(57)

Moreover

$\begin{split}{\Delta_{W}^{2}\tilde{\bm{\omega}}}&=W_{xx}\Delta^{2}\tilde{\omega}_{x}+W_{yy}\Delta^{2}\tilde{\omega}_{y}\geq\frac{W_{xx}}{F_{xx}}+\frac{W_{yy}}{F_{yy}}=\frac{1}{F_{xx}+F_{yy}}\left(\frac{W_{xx}}{\frac{F_{xx}}{F_{xx}+F_{yy}}}+\frac{W_{yy}}{\frac{F_{yy}}{F_{xx}+F_{yy}}}\right)\\\
&\geq\frac{1}{F_{xx}+F_{yy}}\min_{p\in[0,1]}\left(\frac{W_{xx}}{p}+\frac{W_{yy}}{1-p}\right)=\frac{(\sqrt{W_{xx}}+\sqrt{W_{yy}})^{2}}{F_{xx}+F_{yy}}.\end{split}$
(58)

Therefore we may focus on an upper bound for
$\sum_{i=x,y}\braket{\psi_{\bm{\omega}}}{\sigma_{i}\Pi_{\mathcal{H}_{\mathcal{C}}}\sigma_{i}}{\psi_{\bm{\omega}}}$.

Let $\\{\ket{c_{0}},\ket{c_{1}},\ket{c_{2}}\\}$ be an orthonormal basis of
$\mathcal{H}_{\mathcal{C}}\subseteq\mathcal{H}_{S}\otimes\mathcal{H}_{A}$.
These vectors can be written down as

$\ket{c_{i}}=\cos(\varphi^{i})\ket{0}\ket{A^{i}_{0}}+\sin(\varphi^{i})\ket{1}\ket{A^{i}_{1}},$
(59)

where $\ket{A^{i}_{0/1}}$ are normalized states in $\mathcal{H}_{A}$ and
$\varphi^{i}\in[0,\frac{\pi}{2}]$ (a potential complex phase is incorporated
in the definition of$\ket{A^{i}_{0/1}}$). The QEC condition requires
$\forall_{i,j}\braket{c_{i}}{\sigma_{z}}{c_{j}}=\lambda\delta_{ij}$, which
leads to the following two constraints: (i)
$\forall_{i}\cos^{2}(\varphi^{i})-\sin^{2}(\varphi^{i})=\lambda$ means all
$\varphi^{i}$ are equal (therefore superscript $i$ will be omitted); (ii)
$\forall_{i\neq
j}\cos^{2}(\varphi)\braket{A^{i}_{0}}{A^{j}_{0}}-\sin^{2}(\varphi)\braket{A^{i}_{1}}{A^{j}_{1}}=0$.
Together with the orthonormality condition $\forall_{i\neq
j}\cos^{2}(\varphi)\braket{A^{i}_{0}}{A^{j}_{0}}+\sin^{2}(\varphi)\braket{A^{i}_{1}}{A^{j}_{1}}=0$,
we have

$\begin{split}&\ket{c_{i}}=\cos(\varphi)\ket{0}\ket{A^{i}_{0}}+\sin(\varphi)\ket{1}\ket{A^{i}_{1}},\\\
&\quad\forall_{i,j}\braket{A_{0}^{i}}{A_{0}^{j}}=\delta_{ij},\quad\braket{A_{1}^{i}}{A_{1}^{j}}=\delta_{ij}.\end{split}$
(60)

Note that there is no fixed relationship between sets
$\\{\ket{A_{0}^{i}}\\}_{i=0}^{i=2}$ and $\\{\ket{A_{1}^{i}}\\}_{i=0}^{i=2}$—in
particular it may happen that
$\mathrm{span}\\{\ket{A_{0}^{i}}\\}\neq\mathrm{span}\\{\ket{A_{1}^{i}}\\}$.
Effective generators in the chosen basis are given as:

$\displaystyle(G^{\mathcal{R}}_{x})_{ji}$
$\displaystyle=\braket{c_{j}}{\sigma_{x}}{c_{i}}=\frac{\sin(2\varphi)}{2}(\braket{A_{0}^{j}}{A_{1}^{i}}+\braket{A_{1}^{j}}{A_{0}^{i}}),$
(61) $\displaystyle(G^{\mathcal{R}}_{y})_{ji}$
$\displaystyle=\braket{c_{j}}{\sigma_{y}}{c_{i}}=i\frac{\sin(2\varphi)}{2}(\braket{A_{1}^{j}}{A_{0}^{i}}-\braket{A_{0}^{j}}{A_{1}^{i}}).$
(62)

We focus on estimation around point ${\bm{\omega}}=[0,...,0]$ for which
$\ket{\psi_{\bm{\omega}}}=\ket{\psi_{\rm in}}=\ket{c_{0}}$. Then

$\braket{\psi_{\bm{\omega}}}{\sigma_{{x}/{y}}\Pi_{\mathcal{H}_{\mathcal{C}}}\sigma_{{x}/{y}}}{\psi_{\bm{\omega}}}=\sum_{i=0}^{2}\braket{c_{0}}{\sigma_{{x}/{y}}}{c_{i}}\braket{c_{i}}{\sigma_{{x}/{y}}}{c_{0}}=\\\
\frac{\sin^{2}(2\varphi)}{4}\big{(}\sum_{i=0}^{2}\left|{\braket{A_{0}^{0}}{A_{1}^{i}}}\right|^{2}+\left|{\braket{A_{1}^{0}}{A_{0}^{i}}}\right|^{2}\pm
2(\braket{A_{1}^{0}}{A_{0}^{i}}\braket{A_{1}^{i}}{A_{0}^{0}})\big{)}.$ (63)

Since for each $k=0/1$, states $\\{\ket{A_{k}^{i}}\\}_{i=0}^{2}$ are mutually
orthonormal,

$\sum_{i=x,y}\braket{\psi}{\sigma_{i}\Pi_{\mathcal{H}_{\mathcal{C}}}\sigma_{i}}{\psi}=\frac{\sin^{2}(2\varphi)}{2}\cdot\sum_{i=0}^{2}(\left|{\braket{A_{0}^{0}}{A_{1}^{i}}}\right|^{2}+\left|{\braket{A_{1}^{0}}{A_{0}^{i}}}\right|^{2})\leq\sin^{2}(2\varphi)\leq
1,$ (64)

where the first inequality is saturated if and only if both
$\ket{A^{0}_{0}}\in\mathrm{span}\\{\ket{A^{i}_{1}}\\}_{i=0}^{2}$ and
$\ket{A^{0}_{1}}\in\mathrm{span}\\{\ket{A^{i}_{0}}\\}_{i=0}^{2}$. Using Eq.
(57) and Eq. (58) we get
${\Delta_{W}^{2}\tilde{\bm{\omega}}}\geq\frac{(\sqrt{W_{xx}}+\sqrt{W_{yy}})^{2}}{4T^{2}}$.
This implies that the JNT-QEC stategy cannot outperform the best SEP-QEC
strategy.

##  Appendix C The optimal SEP protocol for sensing all magnetic field
components in presence of correlated dephasing noise in the two-qubit model

Here we prove formally, that the optimal precision achievable in the second
example from the main text is ${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm
SEP}=\frac{9}{4T^{2}}$.

Let us briefly review the problem. We consider a two-atom system with
Hamiltonian $H=\frac{1}{2}\sum_{i=1}^{2}{\bm{\omega}}\cdot{\bm{\sigma}}^{(i)}$
(where
${\bm{\sigma}}^{(i)}=[\sigma_{x}^{(i)},\sigma_{y}^{(i)},\sigma_{z}^{(i)}]$
acts on the $i^{\text{th}}$ atom) and a single Lindblad operator
$L=\sqrt{2\gamma}(\sigma_{z}^{(1)}-\sigma_{z}^{(2)})$, therefore

$\mathcal{S}=\mathrm{span}_{\mathbb{R}}\\{\mathbbm{1},L,L^{2}\\}=\mathrm{span}_{\mathbb{R}}\\{\mathbbm{1},\sigma_{z}^{(1)}-\sigma_{z}^{(2)},\sigma_{z}^{(1)}\sigma_{z}^{(2)}\\}.$
(65)

As for all three generators the minimal and the maximal eigenvalues are
respectively $-1$ and $+1$, we immediately see that
$\Delta^{2}\tilde{\omega}_{x,y,z}\geq\frac{1}{4T^{2}}$ and from that
${\Delta_{W}^{2}\tilde{\bm{\omega}}}_{\rm
SEP}\geq\frac{P^{2}}{4T^{2}}=\frac{9}{4T^{2}}$ (which cannot be improved by
applying transformation $A$). Below we show that such a precision is indeed
achievable.

First, $\omega_{z}$ can be estimated using a decoherence-free subspace [76]
$\mathrm{span}\\{\ket{\psi_{z}},\frac{1}{2}(\sigma_{z}^{(1)}+\sigma_{z}^{(2)})\ket{\psi_{z}}\\}$,
where
$\ket{\psi_{z}}=\frac{1}{\sqrt{2}}(\ket{\uparrow\uparrow}+\ket{\downarrow\downarrow})$,
which leads to the precision $\Delta^{2}\omega_{z}=\frac{1}{4T^{2}}$.

In case of $\omega_{x}$ the situation is slightly more complicated, as using
the analogue approach the subspace
$\mathrm{span}\\{\ket{\psi_{x}},\frac{1}{2}(\sigma_{x}^{(1)}+\sigma_{x}^{(2)})\ket{\psi_{x}}\\}$
would not satisfy the QEC conditions. To get the desired estimation precision
we need to find the state $\ket{\psi_{x}}$ which is optimal for measuring
$\omega_{x}$ (from the point of view of noiseless, single-parameter
estimation) and for which

$\mathcal{H}_{\mathcal{C}_{x}}=\mathrm{span}\\{\ket{\psi_{x}},\frac{1}{2}(\sigma_{x}^{(1)}+\sigma_{x}^{(2)})\ket{\psi_{x}}\\}$
(66)

satisfies the QEC condition:

$\begin{split}\Pi_{\mathcal{H}_{\mathcal{C}_{x}}}(\sigma_{z}^{(1)}-\sigma_{z}^{(2)})\Pi_{\mathcal{H}_{\mathcal{C}_{x}}}\propto\Pi_{\mathcal{H}_{\mathcal{C}_{x}}},\\\
\Pi_{\mathcal{H}_{\mathcal{C}_{x}}}\sigma_{z}^{(1)}\sigma_{z}^{(2)}\Pi_{\mathcal{H}_{\mathcal{C}_{x}}}\propto\Pi_{\mathcal{H}_{\mathcal{C}_{x}}},\end{split}$
(67)

and, moreover, other generators act trivially inside
$\mathcal{H}_{\mathcal{C}_{x}}$:

$\begin{split}\Pi_{\mathcal{H}_{\mathcal{C}_{y}}}(\sigma_{y}^{(1)}+\sigma_{y}^{(2)})\Pi_{\mathcal{H}_{\mathcal{C}_{x}}}\propto\Pi_{\mathcal{H}_{\mathcal{C}_{x}}},\\\
\Pi_{\mathcal{H}_{\mathcal{C}_{z}}}(\sigma_{z}^{(1)}+\sigma_{z}^{(2)})\Pi_{\mathcal{H}_{\mathcal{C}_{x}}}\propto\Pi_{\mathcal{H}_{\mathcal{C}_{x}}}.\end{split}$
(68)

It is known that for single-parameter frequency estimation the optimal state
corresponds to an equally weighted superposition of states with minimal and
maximal eigenvalues of the generator. For
$\frac{1}{2}(\sigma_{x}^{(1)}+\sigma_{x}^{(2)})$ it will be:

$\ket{\psi_{x}^{\varphi}}=\frac{1}{2\sqrt{2}}\cdot\left((\ket{\uparrow}+\ket{\downarrow})(\ket{\uparrow}+\ket{\downarrow})+e^{i\varphi}(\ket{\uparrow}-\ket{\downarrow})(\ket{\uparrow}-\ket{\downarrow})\right)$
(69)

for any $\varphi\in\mathbb{R}$. Note, that any superposition of
$\ket{\psi_{x}^{\varphi}}$ (with different $\varphi$) entangled with separated
ancillae is still optimal for sensing $\omega_{x}$. Therefore, we can take
$\ket{\psi_{x}}=\frac{1}{\sqrt{2}}(\ket{\psi_{x}^{0}}\ket{1}_{A}+\ket{\psi_{x}^{\pi}}\ket{2}_{A})$.
We have:

$\ket{\psi_{x}}=\frac{1}{2}\big{(}(\ket{\uparrow\uparrow}+\ket{\downarrow\downarrow})\ket{1}_{A}+(\ket{\uparrow\downarrow}+\ket{\downarrow\uparrow})\ket{2}_{A}\big{)},$
(70)
$\frac{1}{2}(\sigma_{x}^{(1)}+\sigma_{x}^{(2)})\ket{\psi_{x}}=\frac{1}{2}\big{(}(\ket{\uparrow\downarrow}+\ket{\downarrow\uparrow})\ket{1}_{A}+(\ket{\uparrow\uparrow}+\ket{\downarrow\downarrow})\ket{2}_{A}\big{)}.$
(71)

and then the code space Eq. (66) satisfies Eq. (67) and Eq. (68). It gives as
$\Delta^{2}\omega_{x}=\frac{1}{4T^{2}}$. Analogous reasoning could be provided
for $\omega_{y}$.

Therefore for $\rho_{\rm
in}=\frac{1}{3}\sum_{i=x,y,z}\ket{\psi_{i}}\bra{\psi_{i}}$ we have
${\Delta_{W}^{2}\tilde{\bm{\omega}}}=\frac{9}{4T^{2}}$ in line with general
considerations on the performance of the SEP-QEC codes as given in Sec. 5.

##  Appendix D Estimating the $SU(d)$ generators

Below we present an example of a JNT-QEC protocol allowing one to achieve the
total cost ${\Delta_{W}^{2}\tilde{\bm{\omega}}}=\Theta(P^{\frac{3}{2}})$ for
the last example in the main text. For clarification, we treat $d$-dimensional
Hilbert space as a single spin-$j$ particle ($d=2j+1$) and we use the notation
where $\\{\ket{k}\\}_{k=-j}^{j}$ is the eigenbasis of the $J_{z}$ operator.

We consider a problem where the noise generator $J_{z}$ and the unitary
evolution $H$ read:

$J_{z}=\sum_{k=-j}^{j}k\ket{k}\bra{k},\qquad H=\sum_{i=1}^{P}\omega_{i}G_{i},$
(72)

where $G_{i}$ is an orthonormal basis of $\mathcal{S}^{\perp}$—the orthogonal
complement of $\mathcal{S}=\mathrm{span}\\{\mathbbm{1},J_{z},J_{z}^{2}\\}$
(therefore $P=d^{2}-3$). For technical reasons we distinguish three groups of
operators that form the basis $\\{G_{i}\\}$:

  * •

Real off-diagonal:
$G^{R}_{kl}=\frac{1}{\sqrt{2}}(\ket{k}\bra{l}+\ket{l}\bra{k})$

  * •

Imaginary off-diagonal:
$G^{I}_{kl}=\frac{i}{\sqrt{2}}(\ket{k}\bra{l}-\ket{l}\bra{k})$

  * •

Diagonal: $G^{D}_{i}=\sum_{k=-j}^{j}g_{i}^{k}\ket{k}\bra{k}$

and in what follows we prove the scaling $\Theta(P^{\frac{3}{2}})$ for each
group. For simplicity, we assume, that $j$ is an integer (for half-integer $j$
the proof remains almost the same) and in this section we focus on the
estimation around point ${\bm{\omega}}=[0,\ldots,0]$ and set $T=1$.

Real off-diagonal generators. We take
$\dim(\mathcal{H}_{A})=\dim(\mathcal{H}_{S})$ and the state
$\ket{\psi_{\bm{\omega}}}=\ket{\psi^{R}}=\frac{1}{\sqrt{2j+1}}\sum_{k=-j}^{j}\ket{k}\ket{k}_{A}\in\mathcal{H}_{S}\otimes\mathcal{H}_{A}$,
we have

$\braket{\psi^{R}}{J_{z}}{\psi^{R}}=0,\quad\braket{\psi^{R}}{J_{z}^{2}}{\psi^{R}}=\frac{j(j+1)}{3}.$
(73)

We construct the code space in the following way. First, we act on
$\ket{\psi}$ with generators:
$G_{kl}^{R}\ket{\psi}=\frac{1}{\sqrt{2}}(\ket{k}\ket{l}_{A}+\ket{l}\ket{k}_{A})$
and then we “fix it” to satisfy the QEC condition by extending ancilla
$\mathcal{H}_{A}\rightarrow\mathcal{H}_{A}\oplus\mathcal{H}_{B}$ and adding
more terms:

$\ket{c^{R}_{kl}}=\frac{p}{\sqrt{2}}(\ket{k}\ket{l}_{A}+\ket{l}\ket{k}_{A})++q\ket{j}\ket{klj}_{B}+r\ket{-j}\ket{kl(-j)}_{B}+s\ket{0}\ket{kl0}_{B},$
(74)

where
$\braket{klm}{k^{\prime}l^{\prime}m^{\prime}}_{B}=\delta_{(klm)(k^{\prime}l^{\prime}m^{\prime})}$.
Then the QEC condition is equivalent to:

$\begin{split}&\braket{c^{R}_{kl}}{c^{R}_{kl}}=p^{2}+q^{2}+r^{2}+s^{2}=1,\\\
&\braket{c^{R}_{kl}}{J_{z}}{c^{R}_{kl}}=\frac{p^{2}}{2}(k+l)+(q^{2}-r^{2})j=0,\\\
&\braket{c^{R}_{kl}}{J_{z}^{2}}{c^{R}_{kl}}=\frac{p^{2}}{2}(k^{2}+l^{2})+(q^{2}+r^{2})j^{2}=\frac{j(j+1)}{3}.\\\
\end{split}$ (75)

The off-diagonal terms are automatically zero, no matter what $p,q,r,s$ are.
We can write down $q^{2}$, $r^{2}$ and $s^{2}$ as linear functions of $p^{2}$:

$\begin{split}q^{2}&=\frac{1}{2j^{2}}\Big{(}\frac{j(j+1)}{3}-\frac{p^{2}}{2}(k^{2}+l^{2}+j(k+l))\Big{)},\\\
r^{2}&=\frac{1}{2j^{2}}\Big{(}\frac{j(j+1)}{3}-\frac{p^{2}}{2}(k^{2}+l^{2}-j(k+l))\Big{)},\\\
s^{2}&=1-p^{2}-\frac{1}{j^{2}}\Big{(}\frac{j(j+1)}{3}-\frac{p^{2}}{2}(k^{2}+l^{2})\Big{)}.\end{split}$
(76)

Note that $p$ is a valid coefficient if the above set of equations has a
solution (i.e. if the right-hand sides are positive). As $-2j\leq k+l\leq 2j$,
$k^{2}+l^{2}\leq 2j^{2}$, this always holds provided $p^{2}=\frac{1}{6}$. For
the code space ${\mathcal{H}_{\mathcal{C}_{R}}}$ spanned by vectors
constructed in such a way, we have

$(G_{kl}^{R})^{{\mathcal{H}_{\mathcal{C}_{R}}}}\ket{\psi^{R}}=\frac{p}{\sqrt{2j+1}}\ket{c^{R}_{kl}}.$
(77)

The QFIs are

$F^{R}_{(kl)(k^{\prime}l^{\prime})}=4(\braket{\psi^{R}}{(G_{kl}^{R})^{{\mathcal{H}_{\mathcal{C}_{R}}}}(G_{k^{\prime}l^{\prime}}^{R})^{{\mathcal{H}_{\mathcal{C}_{R}}}}}{\psi^{R}}-\braket{\psi^{R}}{(G_{kl}^{R})^{{\mathcal{H}_{\mathcal{C}_{R}}}}}{\psi^{R}}\braket{\psi^{R}}{(G_{k^{\prime}l^{\prime}}^{R})^{{\mathcal{H}_{\mathcal{C}_{R}}}}}{\psi^{R}})$
(78)

which in our case simplifies to:

$F^{R}_{(kl)(k^{\prime}l^{\prime})}=\delta_{(kl)(k^{\prime}l^{\prime})}4\braket{\psi^{R}}{((G_{kl}^{R})^{{\mathcal{H}_{\mathcal{C}_{R}}}})^{2}}{\psi^{R}}=\frac{4p^{2}}{2j+1}.$
(79)

As
$\braket{\psi^{R}}{[(G_{kl}^{R})^{{\mathcal{H}_{\mathcal{C}_{R}}}},(G_{k^{\prime}l^{\prime}}^{R})^{{\mathcal{H}_{\mathcal{C}_{R}}}}]}{\psi^{R}}=0$,
the CR bound is saturable and the total cost is

$\sum_{k>l}\Delta^{2}\omega^{R}_{kl}=\frac{2j+1}{4}\sum_{k>l}\frac{1}{p^{2}}=\frac{3j(2j+1)^{2}}{2}=\Theta(P^{\frac{3}{2}}).$
(80)

Imaginary off-diagonal generators. The reasoning is analogous to the previous
case. Note that using different ancillary spaces for real and imaginary
generators is needed. Even though
$\braket{\psi^{R}}{G^{R}_{kl}J_{z}G^{R}_{k^{\prime}l^{\prime}}}{\psi^{R}}\propto\delta_{kl,k^{\prime}l^{\prime}}$
and
$\braket{\psi^{R}}{G^{I}_{kl}J_{z}G^{I}_{k^{\prime}l^{\prime}}}{\psi^{R}}\propto\delta_{kl,k^{\prime}l^{\prime}}$
are satisfied automatically,
$\braket{\psi^{R}}{G^{R}_{kl}J_{z}G^{I}_{k^{\prime}l^{\prime}}}{\psi^{R}}\propto\delta_{kl,k^{\prime}l^{\prime}}$
may not be true.

Diagonal generators. The number of diagonal generators scales like $\Theta(j)$
(whereas for off-diagonal elements the scaling is $\Theta(j^{2})$), implying
that the estimation with respect to diagonal generators does not contribute
significantly to the overall scaling. Therefore we could simply use the SEP-
QEC approach. Following [36], any traceless generator may by written down as:

$G^{D}_{i}=\frac{1}{2}\mathrm{Tr}(|G^{D}_{i}|)(\rho_{i+}-\rho_{i-}).$ (81)

We define states $\ket{c_{i+}},\ket{c_{i-}}$ as purifications of these density
matrices by using mutually orthogonal ancillary subspaces
$\mathcal{H}_{Ai+},\mathcal{H}_{Ai-}$:

$\rho_{i+/-}=\mathrm{Tr}_{Ai+/-}(\ket{c_{i+/-}}\bra{c_{i+/-}}).$ (82)

Therefore

$\begin{split}\bra{c_{i+}}G^{D}_{i}\ket{c_{i+}}&=\frac{1}{2}\mathrm{Tr}(|G^{D}_{i}|)\geq\frac{1}{2\sqrt{2j+1}},\\\
\bra{c_{i-}}G^{D}_{i}\ket{c_{i-}}&=-\frac{1}{2}\mathrm{Tr}(|G^{D}_{i}|)\leq-\frac{1}{2\sqrt{2j+1}},\end{split}$
(83)

and from that for code space $\mathrm{span}\\{\ket{c_{i+}},\ket{c_{i-}}\\}$
and input state $\frac{1}{\sqrt{2}}(\ket{c_{i+}}+\ket{c_{i-}})$ we have
$F_{\omega_{i}}\geq\frac{1}{2j+1}$. For a single-parameter problem the CR
bound is always saturable so using the SEP-QEC approach for the input state

$\rho_{\rm
in}^{D}=\frac{1}{2j-2}\sum_{i=1}^{2j-2}\frac{1}{2}\left(\ket{c_{i+}}+\ket{c_{i-}}\right)\left(\bra{c_{i+}}+\bra{c_{i-}}\right)$
(84)

we have

$\sum_{i=1}^{2j-2}\Delta^{2}\omega^{D}_{i}=(2j-2)\sum_{i=1}^{2j-2}\frac{1}{F_{\omega_{i}}}\leq(2j+1)(2j-2)^{2}=\Theta(P^{\frac{3}{2}}).$
(85)

Results. Finally, combining all above, for the input state:

$\rho_{\rm
in}=\frac{1}{3}(\ket{\psi^{R}}\bra{\psi^{R}}+\ket{\psi^{I}}\bra{\psi^{I}}+\rho^{D})$
(86)

(with properly applied QEC protocol) we get

$\sum_{i=1}^{P}\Delta^{2}\omega_{i}=\Theta(P^{\frac{3}{2}}).$ (87)

## References

  * Giovannetti  _et al._ [2006] V. Giovannetti, S. Lloyd, and L. Maccone,  _Quantum metrology,_Phys. Rev. Lett. 96, 010401 (2006).
  * Paris [2009] M. G. A. Paris,  _Quantum estimation for quantum technologies,_Int. J. Quantum Inf. 07, 125 (2009).
  * Giovannetti  _et al._ [2011] V. Giovannetti, S. Lloyd, and L. Maccone,  _Advances in quantum metrology,_Nat. Photonics 5, 222 (2011).
  * Toth and Apellaniz [2014] G. Toth and I. Apellaniz,  _Quantum metrology from a quantum information science perspective,_J. Phys. A: Math. Theor. 47, 424006 (2014).
  * Demkowicz-Dobrzański  _et al._ [2015] R. Demkowicz-Dobrzański, M. Jarzyna, and J. Kołodyński, in  _Prog. Optics_, Vol. 60, edited by E. Wolf (Elsevier, 2015) pp. 345–435.
  * Schnabel [2017] R. Schnabel,  _Squeezed states of light and their applications in laser interferometers,_Phys. Rep. 684, 1 (2017).
  * Degen  _et al._ [2017] C. L. Degen, F. Reinhard, and P. Cappellaro,  _Quantum sensing,_Rev. Mod. Phys. 89, 035002 (2017).
  * Pezzè  _et al._ [2018] L. Pezzè, A. Smerzi, M. K. Oberthaler, R. Schmied, and P. Treutlein,  _Quantum metrology with nonclassical states of atomic ensembles,_Rev. Mod. Phys. 90, 035005 (2018).
  * Pirandola  _et al._ [2018] S. Pirandola, B. R. Bardhan, T. Gehring, C. Weedbrook, and S. Lloyd,  _Advances in photonic quantum sensing,_Nat. Photonics 12, 724 (2018).
  * Caves [1981] C. M. Caves,  _Quantum-mechanical noise in an interferometer,_Phys. Rev. D 23, 1693 (1981).
  * Holland and Burnett [1993] M. Holland and K. Burnett,  _Interferometric detection of optical phase shifts at the heisenberg limit,_Phys. Rev. Lett. 71, 1355 (1993).
  * Lee  _et al._ [2002] H. Lee, P. Kok, and J. P. Dowling,  _A quantum rosetta stone for interferometry,_J. Mod. Optic. 49, 2325 (2002).
  * Wineland  _et al._ [1992] D. Wineland, J. Bollinger, W. Itano, F. Moore, and D. Heinzen,  _Spin squeezing and reduced quantum noise in spectroscopy,_Phys. Rev. A 46, R6797 (1992).
  * McKenzie  _et al._ [2002] K. McKenzie, D. A. Shaddock, D. E. McClelland, B. C. Buchler, and P. K. Lam,  _Experimental demonstration of a squeezing-enhanced power-recycled michelson interferometer for gravitational wave detection,_Phys. Rev. Lett. 88, 231102 (2002).
  * Bollinger  _et al._ [1996] J. Bollinger, W. M. Itano, D. Wineland, and D. Heinzen,  _Optimal frequency measurements with maximally correlated states,_Phys. Rev. A 54, R4649 (1996).
  * Leibfried  _et al._ [2004] D. Leibfried, M. Barrett, T. Schaetz, J. Britton, J. Chiaverini, W. Itano, J. Jost, C. Langer, and D. Wineland,  _Toward heisenberg-limited spectroscopy with multiparticle entangled states,_Science 304, 1476 (2004).
  * Giovannetti  _et al._ [2004] V. Giovannetti, S. Lloyd, and L. Maccone,  _Quantum-enhanced measurements: beating the standard quantum limit,_Science 306, 1330 (2004).
  * Huelga  _et al._ [1997] S. F. Huelga, C. Macchiavello, T. Pellizzari, A. K. Ekert, M. B. Plenio, and J. I. Cirac,  _Improvement of frequency standards with quantum entanglement,_Phys. Rev. Lett. 79, 3865 (1997).
  * Berry and Wiseman [2000] D. W. Berry and H. M. Wiseman,  _Optimal states and almost optimal adaptive measurements for quantum interferometry,_Phys. Rev. Lett. 85, 5098 (2000).
  * de Burgh and Bartlett [2005] M. de Burgh and S. D. Bartlett,  _Quantum methods for clock synchronization: Beating the standard quantum limit without entanglement,_Phys. Rev. A 72, 042301 (2005).
  * Fujiwara and Imai [2008] A. Fujiwara and H. Imai,  _A fibre bundle over manifolds of quantum channels and its application to quantum statistics,_J. Phys. A: Math. Theor. 41, 255304 (2008).
  * Demkowicz-Dobrzański  _et al._ [2009] R. Demkowicz-Dobrzański, U. Dorner, B. Smith, J. Lundeen, W. Wasilewski, K. Banaszek, and I. Walmsley,  _Quantum phase estimation with lossy interferometers,_Phys. Rev. A 80, 013825 (2009).
  * Escher  _et al._ [2011] B. Escher, R. de Matos Filho, and L. Davidovich,  _General framework for estimating the ultimate precision limit in noisy quantum-enhanced metrology,_Nat. Phys. 7, 406 (2011).
  * Demkowicz-Dobrzański  _et al._ [2012] R. Demkowicz-Dobrzański, J. Kołodyński, and M. Guţă,  _The elusive heisenberg limit in quantum-enhanced metrology,_Nat. Commun. 3, 1063 (2012).
  * Kołodyński and Demkowicz-Dobrzański [2013] J. Kołodyński and R. Demkowicz-Dobrzański,  _Efficient tools for quantum metrology with uncorrelated noise,_New J. Phys. 15, 073043 (2013).
  * Knysh  _et al._ [2014] S. I. Knysh, E. H. Chen, and G. A. Durkin,  _True limits to precision via unique quantum probe,_arXiv:1402.0495 (2014).
  * Demkowicz-Dobrzański and Maccone [2014] R. Demkowicz-Dobrzański and L. Maccone,  _Using entanglement against noise in quantum metrology,_Phys. Rev. Lett. 113, 250801 (2014).
  * Kessler  _et al._ [2014] E. M. Kessler, I. Lovchinsky, A. O. Sushkov, and M. D. Lukin,  _Quantum error correction for metrology,_Phys. Rev. Lett. 112, 150802 (2014).
  * Dür  _et al._ [2014] W. Dür, M. Skotiniotis, F. Froewis, and B. Kraus,  _Improved quantum metrology using quantum error correction,_Phys. Rev. Lett. 112, 080801 (2014).
  * Ozeri [2013] R. Ozeri,  _Heisenberg limited metrology using quantum error-correction codes._arXiv:1310.3432 (2013).
  * Arrad  _et al._ [2014] G. Arrad, Y. Vinkler, D. Aharonov, and A. Retzker,  _Increasing sensing resolution with error correction,_Phys. Rev. Lett. 112, 150801 (2014).
  * Unden  _et al._ [2016] T. Unden, P. Balasubramanian, D. Louzon, Y. Vinkler, M. B. Plenio, M. Markham, D. Twitchen, A. Stacey, I. Lovchinsky, A. O. Sushkov, _et al._ ,  _Quantum metrology enhanced by repetitive quantum error correction,_Phys. Rev. Lett. 116, 230502 (2016).
  * Reiter  _et al._ [2017] F. Reiter, A. S. Sørensen, P. Zoller, and C. A. Muschik,  _Dissipative quantum error correction and application to quantum sensing with trapped ions,_Nat. Commun. 8, 1822 (2017).
  * Sekatski  _et al._ [2017] P. Sekatski, M. Skotiniotis, J. Kołodyński, and W. Dür,  _Quantum metrology with full and fast quantum control,_Quantum 1, 27 (2017).
  * Demkowicz-Dobrzański  _et al._ [2017] R. Demkowicz-Dobrzański, J. Czajkowski, and P. Sekatski,  _Adaptive quantum metrology under general markovian noise,_Phys. Rev. X 7, 041009 (2017).
  * Zhou  _et al._ [2018] S. Zhou, M. Zhang, J. Preskill, and L. Jiang,  _Achieving the heisenberg limit in quantum metrology using quantum error correction,_Nat. Commun. 9, 78 (2018).
  * Layden and Cappellaro [2018] D. Layden and P. Cappellaro,  _Spatial noise filtering through error correction for quantum sensing,_npj Quantum Inf. 4, 30 (2018).
  * Layden  _et al._ [2019] D. Layden, S. Zhou, P. Cappellaro, and L. Jiang,  _Ancilla-free quantum error correction codes for quantum metrology,_Phys. Rev. Lett. 122, 040502 (2019).
  * Kapourniotis and Datta [2019] T. Kapourniotis and A. Datta,  _Fault-tolerant quantum metrology,_Phys. Rev. A 100, 022335 (2019).
  * Tan  _et al._ [2019] K. C. Tan, S. Omkar, and H. Jeong,  _Quantum-error-correction-assisted quantum metrology without entanglement,_Phys. Rev. A 100, 022312 (2019).
  * Zhou and Jiang [2020a] S. Zhou and L. Jiang,  _The theory of entanglement-assisted metrology for quantum channels,_arXiv:2003.10559 (2020a).
  * Chen  _et al._ [2020] Y. Chen, H. Chen, J. Liu, Z. Miao, and H. Yuan,  _Fluctuation-enhanced quantum metrology,_arXiv:2003.13010 (2020).
  * Baumgratz and Datta [2016] T. Baumgratz and A. Datta,  _Quantum enhanced estimation of a multidimensional field,_Phys. Rev. Lett. 116, 030801 (2016).
  * Tsang  _et al._ [2016] M. Tsang, R. Nair, and X.-M. Lu,  _Quantum theory of superresolution for two incoherent optical point sources,_Phys. Rev. X 6, 031033 (2016).
  * Humphreys  _et al._ [2013] P. C. Humphreys, M. Barbieri, A. Datta, and I. A. Walmsley,  _Quantum enhanced multiple phase estimation,_Phys. Rev. Lett. 111, 070403 (2013).
  * Gessner  _et al._ [2018] M. Gessner, L. Pezzè, and A. Smerzi,  _Sensitivity bounds for multiparameter quantum metrology,_Phys. Rev. Lett. 121, 130503 (2018).
  * Tsang  _et al._ [2011] M. Tsang, H. M. Wiseman, and C. M. Caves,  _Fundamental quantum limit to waveform estimation,_Phys. Rev. Lett. 106, 090401 (2011).
  * Berry  _et al._ [2013] D. W. Berry, M. J. W. Hall, and H. M. Wiseman,  _Stochastic heisenberg limit: Optimal estimation of a fluctuating phase,_Phys. Rev. Lett. 111, 113601 (2013).
  * Matsumoto [2002] K. Matsumoto,  _A new approach to the cramér-rao-type bound of the pure-state model,_J. Phys. A.: Math. Theor. 35, 3111 (2002).
  * Genoni  _et al._ [2013] M. G. Genoni, M. G. A. Paris, G. Adesso, H. Nha, P. L. Knight, and M. S. Kim,  _Optimal estimation of joint parameters in phase space,_Phys. Rev. A 87, 012107 (2013).
  * Ragy  _et al._ [2016] S. Ragy, M. Jarzyna, and R. Demkowicz-Dobrzański,  _Compatibility in multiparameter quantum metrology,_Phys. Rev. A 94, 052108 (2016).
  * Yuan [2016] H. Yuan,  _Sequential feedback scheme outperforms the parallel scheme for hamiltonian parameter estimation,_Phys. Rev. Lett. 117, 160801 (2016).
  * Kura and Ueda [2018] N. Kura and M. Ueda,  _Finite-error metrological bounds on multiparameter hamiltonian estimation,_Phys. Rev. A 97, 012101 (2018).
  * Liu and Yuan [2017] J. Liu and H. Yuan,  _Control-enhanced multiparameter quantum estimation,_Phys. Rev. A 96, 042114 (2017).
  * Nichols  _et al._ [2018] R. Nichols, P. Liuzzo-Scorpo, P. A. Knott, and G. Adesso,  _Multiparameter gaussian quantum metrology,_Phys. Rev. A 98, 012114 (2018).
  * Ge  _et al._ [2018] W. Ge, K. Jacobs, Z. Eldredge, A. V. Gorshkov, and M. Foss-Feig,  _Distributed quantum metrology with linear networks and separable inputs,_Phys. Rev. Lett. 121, 043604 (2018).
  * Braunstein and Caves [1994] S. L. Braunstein and C. M. Caves,  _Statistical distance and the geometry of quantum states,_Phys. Rev. Lett. 72, 3439 (1994).
  * Helstrom [1976] C. W. Helstrom,  _Quantum detection and estimation theory_ (Academic press, 1976).
  * Holevo [1982] A. S. Holevo,  _Probabilistic and Statistical Aspects of Quantum Theory_ (North Holland, Amsterdam, 1982).
  * Demkowicz-Dobrzanski  _et al._ [2020] R. Demkowicz-Dobrzanski, W. Gorecki, and M. Guta,  _Multi-parameter estimation beyond quantum fisher information,_Journal of Physics A: Mathematical and Theoretical (2020).
  * Nagaoka and Hayashi [2005] H. Nagaoka and M. Hayashi,  _Asymptotic Theory of Quantum Statistical Inference_ (World Scientific Singapore, 2005) Chap. 8.
  * Suzuki [2016] J. Suzuki,  _Explicit formula for the holevo bound for two-parameter qubit-state estimation problem,_J. Math. Phys. 57, 042201 (2016).
  * Guţă and Jenčová [2007] M. Guţă and A. Jenčová,  _Local asymptotic normality in quantum statistics,_Comm. Math. Phys. 276, 341 (2007).
  * Yamagata  _et al._ [2013] K. Yamagata, A. Fujiwara, R. D. Gill, _et al._ ,  _Quantum local asymptotic normality based on a new quantum likelihood ratio,_Ann. Statist. 41, 2197 (2013).
  * Fujiwara [1994] A. Fujiwara,  _Multi-parameter pure state estimation based on the right logarithmic derivative,_ Math. Eng. Tech. Rep 94, 94 (1994).
  * Albarelli  _et al._ [2019] F. Albarelli, J. F. Friel, and A. Datta,  _Evaluating the holevo cramér-rao bound for multiparameter quantum metrology,_Phys. Rev. Lett. 123, 200503 (2019).
  * Lindblad [1976] G. Lindblad,  _On the generators of quantum dynamical semigroups,_Comm. Math. Phys. 48, 119 (1976).
  * Gorini  _et al._ [1976] V. Gorini, A. Kossakowski, and E. C. G. Sudarshan,  _Completely positive dynamical semigroups of n-level systems,_J. Math. Phys. 17, 821 (1976).
  * Breuer  _et al._ [2002] H.-P. Breuer, F. Petruccione, _et al._ ,  _The theory of open quantum systems_ (Oxford University Press on Demand, 2002).
  * Kay [1993] S. M. Kay,  _Fundamentals of statistical signal processing: estimation theory_ (Prentice Hall, 1993).
  * Gill and Massar [2000] R. Gill and S. Massar,  _State estimation for large ensembles,_Phys. Rev. A 61, 042312 (2000).
  * Knill and Laflamme [1997] E. Knill and R. Laflamme,  _Theory of quantum error-correcting codes,_Phys. Rev. A 55, 900 (1997).
  * [73] M. Grant and S. Boyd, Cvx: Matlab software for disciplined convex programming, .
  * Zhou and Jiang [2020b] S. Zhou and L. Jiang,  _Optimal approximate quantum error correction for quantum metrology,_Phys. Rev. Research 2, 013235 (2020b).
  * Górecki  _et al._ [2020] W. Górecki, R. Demkowicz-Dobrzański, H. M. Wiseman, and D. W. Berry,  _$\pi$ -corrected heisenberg limit, _Phys. Rev. Lett. 124, 030501 (2020).
  * Lidar  _et al._ [1998] D. A. Lidar, I. L. Chuang, and K. B. Whaley,  _Decoherence-free subspaces for quantum computation,_Phys. Rev. Lett. 81, 2594 (1998).