##  Chapter 0 Estimation Theory

Estimation theory is the mathematical language of metrology. Statistical error
in classical estimation theory is ultimately constrained by the central limit
theorem. Quantum metrology overcomes this limitation thanks to quantum
entanglement. With the vast number of applications and straightforward proof
of principle, it is unsurprising that quantum metrology is witnessing a boon
of theoretical and experimental developments [giovannetti2011, degen2017,
pirandola2018].

This chapter is divided into three sections. The first section summarizes
important concepts from classical estimation theory [kay1993, cox2006,
sheskin2003, rice2006, poor2013]. The second section is devoted to the
analogous concepts of quantum estimation theory formalized by Helstrom
[helstrom1967, helstrom1968, helstrom1969] and Holevo [holevo1973,
holevo1982]. The final section examines example applications of quantum
metrology (phase estimation and amplitude estimation) to put into perspective
the mathematical tools and concepts introduced throughout the first two
sections. For a quantum information perspective on quantum metrology see
[toth2014]. For a more mathematical rigorous review of quantum metrology and
quantum estimation theory see [simon2017]. For more information on estimation
theory and statistical inference see [kay1993, cox2006].

> An experiment is a question which science poses to Nature and a measurement
> is the recording of Nature’s answer.
>
> -Max Planck

###  1 Classical Estimation Theory

In an abstract sense, the scientific and mathematical knowledge of humankind
is reflected in the mathematical models used to describe the contents of the
universe: planetary orbits, bacterial growth in a petri dish, even social
constructs like financial trends. These models, are not fabricated
haphazardly, instead they are a manifestation of a multitude of observations
and tested by making predictions. As our efficiency of gathering and
interpreting data increases, so do the mathematical models, and in turn our
understanding of the universe. For example, the theory of gravity has evolved
along with the capabilities of telescopes; from Galilean and Newtonian gravity
to Einstein’s theory of general relativity to the (currently unconfirmed)
theory of dark matter and dark energy.

Estimation theory is a branch of statistics at the heart of mathematical
modelling. It addresses the question: ‘What is the most efficient way of
extracting information from a set of data?’. This seemingly simple question is
difficult to answer. Typically, the variables used to describe a mathematical
model can be partitioned in two categories

  1. 1.

observables - an attribute which can be inherently measured (e.g. position and
speed).

  2. 2.

latent parameters - an attribute which cannot be inherently measured, (e.g.
strength of an electromagnetic field).

The parameter estimation problem is concerned with the extent at which
collected data (observables) can be used to estimate the unknown latent
parameters [sheskin2003]. With respect to the listed examples, one could
observe the dynamics of a charged particle to estimate the strength of an
electromagnetic field.

Formally, observed data $\mathbf{x}=\\{x_{1},\ldots,x_{N}\\}$ is treated as a
realisation of $N$ independent and identically distributed (iid) random
variables $X$. A probably density function $p(X|\theta)$ dictates the
distribution of observed data, where $\theta$ is a latent parameter. The goal
of the parameter estimation problem is to construct an estimator
$\hat{\theta}(\mathbf{x})$, which should be interpreted as a function whose
input is the collected data $\mathbf{x}$ and outputs an estimate of $\theta$.
The explicit dependence on $\mathbf{x}$ is sometimes dropped for clarity,
$\hat{\theta}(\mathbf{x})\rightarrow\hat{\theta}$. Estimators are subjected to
two conditions. The first condition is that the expected estimate is the true
value of the parameter, this is known as having an unbiased estimator

$\expval*{\hat{\theta}}=\int
p(\mathbf{x}|\theta)\hat{\theta}(\mathbf{x})\mathrm{d}\mathbf{x}=\theta.$ (1)

The integral equation is used for observed data which can take on a continuum
of values, it is interchangeable with a sum in the discrete case. The second
condition is that an estimator tends towards the correct value as the amount
of data increases, this is known as being consistent

$\lim_{N\rightarrow\infty}\hat{\theta}=\theta.$ (2)

An estimator is a manifestation of random variables, and is thus also a random
variable, hence, statistical moments such as mean and variance are well-
defined.

The statistical inference process adopted to the parameter estimation problem
is dependent on the nature of the latent variable: deterministic or
stochastic. Usually, a frequentist inference approach is taken for
deterministic parameters and a Bayesian inference approach is taken for
stochastic parameters [li2018]. Mathematically, these two approaches vary
greatly, the primary differences are listed in Tab. (1), but they are not
mutually exclusive. The subsequent chapters of this thesis employ the
frequentist approach, and therefore the frequentist approach is summarized in
greater detail in this chapter. That being said, the Bayesian approach has
been adapted to the realm of quantum information [holevo1982, tsang2011], and
has been gaining traction in the community [berry2009, gammelmark2013,
jarzyna2015, wiebe2016, rubio2020b]. Specifically, to circumvent problems of
the frequentist approach: i) lack of a priori knowledge [kolodynski2010,
demkowicz2011] and ii) inaccuracies with limited resources [rubio2020b]. Even
though it is not applied to the research presented in this thesis, for the
sake of completeness, a brief summary of the Bayesian approach used in
classical parameter estimation problems and its adaptation to quantum
parameter estimation problems is included in this chapter.

| Frequentist Approach | Bayesian Approach
---|---|---
Parameter(s) | Deterministic | Stochastic
Figure of Merit | Mean squared error | Cost function
Optimization | Local | Global
Table 1: The main differences between the frequentist approach and Bayesian
approach for statistical inference. This is a broad perspective and the
statistical inference approaches are not restricted by this table.

####  1 The Frequentist Approach

The frequentist approach is typically used when $\mathbf{\theta}$ is
deterministic (sometimes called static). As $N\rightarrow\infty$ the frequency
of collected data tends to reflect the probability density function, hence the
etymology. Therefore with a sufficient amount of collected data, the unknown
parameter can be estimated to any desired precision. The figure of merit used
by the frequentist approach is the mean-squared error (MSE)

$\Delta^{2}\hat{\theta}=\expval*{(\hat{\theta}-\theta)^{2}}=\int
p(\mathbf{x}|\theta)\big{(}\hat{\theta}(\mathbf{x})-\theta\big{)}^{2}\mathrm{d}\mathbf{x},$
(3)

in which the aim is to find an estimator which minimizes the above equation.
Because the estimator is assumed to be unbiased, the MSE is equal to the
variance, which is often a more significant statistical quantity.

The first controversy of the the frequentist approach arises due to the fact
that an optimal estimator (one where Eq. eq:MSE is minimized) is potentially
dependent on $\theta$. Some estimators may be optimal for specific values of
$\theta$ (local), whereas an estimator which is optimal for all values of
$\theta$ (global) can only be worse than ones which are locally optimized. At
first glance, this appears counter intuitive because a locally optimized
estimator requires exact knowledge of $\theta$, which defeats the purpose of
parameter estimation. However, it is reasonable to assume that a priori
approximate knowledge $\theta\approx\theta_{0}$ is often known because of
theory or previous estimates. In the absence of a priori knowledge, one can
construct a locally efficient estimator by increasing $N$. To do so, a
fraction of the results are first used to obtain a local approximation
$\theta_{0}$, and the remaining are used within the locally optimized
estimator. Unfortunately, the frequentist approach does not provide a method
on bounding $N$ such that the local regime can be assured; thus the saturation
of an optimal estimator may not be possible without the ability to infinitely
increase $N$.

####  2 Cramér-Rao Bound and Fisher Information

The Cramér-Rao Bound (CRB) is an inequality which assigns a lower bound to the
MSE of unbiased estimators [cramer1946], the derivation of which is
straightforward. The unbiased condition, Eq. eq:unbiased, can be re-written as

$\int
p(\mathbf{x}|\theta)\big{(}\hat{\theta}(\mathbf{x})-\theta\big{)}\mathrm{d}\mathbf{x}=0,$
(4)

from which it follows that

$\begin{split}0&=\frac{\partial}{\partial\theta}\int
p(\mathbf{x}|\theta)\big{(}\hat{\theta}(\mathbf{x})-\theta\big{)}\mathrm{d}\mathbf{x}\\\
&=\int\frac{\partial
p(\mathbf{x}|\theta)}{\partial\theta}\big{(}\hat{\theta}(\mathbf{x})-\theta\big{)}\mathrm{d}\mathbf{x}-\int
p(\mathbf{x}|\theta)\mathrm{d}\mathbf{x}\\\ &=\int
p(\mathbf{x}|\theta)\frac{\partial\ln
p(\mathbf{x}|\theta)}{\partial\theta}\big{(}\hat{\theta}(\mathbf{x})-\theta\big{)}\mathrm{d}\mathbf{x}-1.\end{split}$
(5)

Using the Cauchy–Schwarz inequality

$\bigg{|}\int f(x)g(x)\mathrm{d}x\bigg{|}^{2}\leq\bigg{(}\int
f(x)^{2}\mathrm{d}x\bigg{)}\cdot\bigg{(}\int g(x)^{2}\mathrm{d}x\bigg{)},$ (6)

with $x\rightarrow\mathbf{x}$,
$f(x)\rightarrow\sqrt{p(\mathbf{x}|\theta)}\frac{\partial\ln
p(\mathbf{x}|\theta)}{\partial\theta}$ and
$g(x)\rightarrow\sqrt{p(\mathbf{x}|\theta)}\big{(}\hat{\theta}(\mathbf{x})-\theta\big{)}$,
Eq. eq:crbproof is transformed into the inequality

$1\leq\bigg{(}\int
p(\mathbf{x}|\theta)\big{(}\hat{\theta}(\mathbf{x})-\theta\big{)}^{2}\mathrm{d}\mathbf{x}\bigg{)}\cdot\bigg{(}\int
p(\mathbf{x}|\theta)\Big{(}\frac{\partial\ln
p(\mathbf{x}|\theta)}{\partial\theta}\Big{)}^{2}\mathrm{d}\mathbf{x}\bigg{)}.$
(7)

The above can be manipulated to obtain the CRB

$\Delta^{2}\hat{\theta}\geq\frac{1}{\mathcal{I}\big{(}p(\mathbf{x}|\theta)\big{)}},$
(8)

where

$\split\mathcal{I}\big{(}p(\mathbf{x}|\theta)\big{)}&=\int
p(\mathbf{x}|\theta)\Big{(}\frac{\partial\ln
p(\mathbf{x}|\theta)}{\partial\theta}\Big{)}^{2}\mathrm{d}\mathbf{x}\\\
=\int\frac{1}{p(\mathbf{x}|\theta)}\Big{(}\frac{\partial
p(\mathbf{x}|\theta)}{\partial\theta}\Big{)}^{2}\mathrm{d}\mathbf{x}\\\ =-\int
p(\mathbf{x}|\theta)\frac{\partial^{2}\ln
p(\mathbf{x}|\theta)}{\partial\theta^{2}}\mathrm{d}\mathbf{x}$ (9)

is the Fisher Information (FI), where three equivalent (assuming that $p$ is
twice differentiable) expressions given. The FI is a non-negative and additive
quantity. Because $\mathbf{x}$ is $N$ independent realisations of the random
variable $X$, the CRB can be equivalently expressed as

$\Delta^{2}\hat{\theta}\geq\frac{1}{N\mathcal{I}\big{(}p(X|\theta)\big{)}}.$
(10)

The above form of the CRB reflects the limitations of central limit theorem:
as $N\rightarrow\infty$ the sample average will take on a normal distribution
with a variance of $\mathcal{O}\big{(}N^{-1}\big{)}$.

The FI is often interpreted as a measure of how much information about an
unknown parameter can be extracted from a probability density function
[fisher1925]. In particular, $\theta$ can be learned perfectly when
$\mathcal{I}\rightarrow\infty$, and conversely no information can be learned
about $\theta$ when $\mathcal{I}=0$. In fact, when viewing probability density
functions as points on a manifold (parameterized by $\theta$), the FI is a
Riemannian metric between neighbouring probability density functions
$p(X|\theta)$ and $p(X|\theta+\delta\theta)$ [nielsen2013]. Similarly, the
statistical angle111This is the classical version of the Bures angle
[wootters1981]. between probability density functions

$D\big{(}p_{1}(x),p_{2}(x)\big{)}=\arccos\int\sqrt{p_{1}(x)p_{2}(x)}\mathrm{d}x,$
(11)

can be expressed as [barndorff1986]

$D\big{(}p(X|\theta),p(X|\theta+\delta\theta)\big{)}=\frac{1}{2}\sqrt{\mathcal{I}\big{(}p(X|\theta)\big{)}}\delta\theta+\mathcal{O}\big{(}\delta\theta^{2}\big{)}.$
(12)

Hence, a probability density function with a high FI will deviate more upon
small perturbations $\delta\theta$ than the opposing case of a probability
density function with a small FI.

Figure 1: The FI for a normal distribution
$p(X|\mu,\sigma)=\frac{1}{\sqrt{2\pi}\sigma}e^{-\frac{1}{2}(\frac{x-\mu}{\sigma})^{2}}$.
Regardless of whether $\mu$ or $\sigma$ is the latent parameter, the FI is
exclusively dependent on the standard deviation of the normal distribution.
This is logical because it is more difficult to interpret data with a larger
standard deviation; which is in accordance with the interpretation of the FI
being a measure of extractable information. One could equally consider the
scenario in which both $\mu$ and $\sigma$ are unknown parameters. Here
multiparameter parameter estimation techniques are needed - which are
discussed in a latter part of this chapter.

The Cauchy-Schwarz inequality, Eq. eq:CauchySchwarz, is saturated if

$\frac{|f(x)|}{|g(x)|}=\frac{\int f(x)^{2}\mathrm{d}x}{\int
g(x)^{2}\mathrm{d}x}$ (13)

Therefore an estimator which saturates the CRB for all $\theta$ (global)
satisfies

$\frac{\partial\ln
p(\mathbf{x}|\theta)}{\partial\theta}=\mathcal{I}\big{(}\mathbf{x}|\theta\big{)}\big{(}\hat{\theta}(\mathbf{x})-\theta\big{)}.$
(14)

An estimator which saturated the CRB is said to be efficient. The above
expression can be equivalently written as

$\split
p(\mathbf{x}|\theta)&=\exp\bigg{(}\int\mathcal{I}\big{(}\mathbf{x}|\theta\big{)}\big{(}\hat{\theta}(\mathbf{x})-\theta\big{)}\mathrm{d}\theta\bigg{)}\\\
=\exp\bigg{(}\frac{\partial\mathcal{J}\big{(}\mathbf{x}|\theta)}{\partial\theta}\big{(}\hat{\theta}(\mathbf{x})-\theta\big{)}+\mathcal{J}\big{(}\mathbf{x}|\theta)+c(\mathbf{x})\bigg{)},$
(15)

where $\mathcal{J}\big{(}\mathbf{x}|\theta\big{)}$ is a function which
satisfies
$\frac{\partial^{2}\mathcal{J}\big{(}\mathbf{x}|\theta\big{)}}{\partial\theta^{2}}=\mathcal{I}\big{(}\mathbf{x}|\theta\big{)}$
and $c(\mathbf{x})$ is an arbitrary function independent of $\theta$, both of
which are chosen such that the unbiased condition, Eq. eq:unbiased, is
satisfied. This general expression for a probability density function can
correspond to a multitude of well-known distributions in statistics with
exponential tendencies: Gaussian, Bernoulli, Poisson, et cetera. It should be
stressed that an efficient global estimator does not necessarily exist,
further it may encounter the earlier stated problem of having a dependence on
$\theta$. A locally (approximately) efficient estimator can be constructed
with prior knowledge that $\theta\approx\theta_{0}$ by re-arranging Eq.
eq:CRBcond

$\hat{\theta}_{\text{Local}}=\theta_{0}+\frac{1}{\mathcal{I}\big{(}\mathbf{x}|\theta_{0}\big{)}}\frac{\partial\ln
p(\mathbf{x}|\theta)}{\partial\theta}\bigg{|}_{\theta\rightarrow\theta_{0}}.$
(16)

Unfortunately, the locally approximate estimator is ultimately constrained by
ones prior knowledge, as shifting
$\theta_{0}\rightarrow\theta_{0}+\delta\theta_{0}$ will similarly shift Eq.
eq:localestimator by $\mathcal{O}(\delta\theta_{0})$. Furthermore, the locally
approximate estimator may be ill defined on certain domains, for example one
of circular symmetry (such as the problem of phase estimation which is
discussed in a later section of this chapter).

####  3 Maximum Likelihood Estimation

The likelihood function $L(\theta|\mathbf{x})=p(\mathbf{x}|\theta)$ is a
goodness of fit between a model and the sampled data. It should be understood
that the likelihood function is not a probability density function; the
observed data $\mathbf{x}$ is held fixed and the latent parameter $\theta$ is
considered a variable. The intuition is simplistic: if
$L(\theta_{1}|\mathbf{x})>L(\theta_{2}|\mathbf{x})$, then it is more likely
that the true value of $\theta$ is $\theta_{1}$ rather than $\theta_{2}$. This
is the principal idea of maximum likelihood estimation [berger1988]. The
maximum likelihood estimator, $\hat{\theta}_{\text{ML}}$, outputs the value of
$\theta$ which maximizes $L(\theta|\mathbf{x})$

$\hat{\theta}_{\text{ML}}(\mathbf{x})=\underset{\theta}{\text{argmax}}L(\theta|\mathbf{x})=\underset{\theta}{\text{argmax}}\ln
L(\theta|\mathbf{x}).$ (17)

Because $p(\mathbf{x}|\theta)$ is a joint probability density function of $N$
independent probability density functions,
$p(\mathbf{x}|\theta)=\prod_{j=1}^{N}p(x_{j}|\theta)$, it is often simpler to
maximize the log of the likelihood function, $\ln L(\theta|\mathbf{x})$,
sometimes shortened to the log-likelihood.

.495 .495

Figure 2: $N=4$ Figure 3: $N=100$ Figure 4: Likelihood function for a biased
coin where $p_{H}$ is an unknown probability of the coin toss resulting in
heads. The maximum likelihood estimation strategy outputs
$\hat{p}_{H}=\frac{\text{\\#Heads}}{N}$. The estimate is sensitive to small
fluctuations in the observed data for small $N$ (a), but becomes more robust
to fluctuations as $N$ increases (b).

One controversy with maximum likelihood estimation is that
$\hat{\theta}_{\text}{ML}$ does not generally satisfy the unbiased condition,
Eq. eq:unbiased. Specifically, for small $N$, where the estimator is much more
susceptive to statistical outliers within the collected data. However, as
$N\rightarrow\infty$ the estimator becomes more unbiased,
$\expval*{\hat{\theta}_{\text{ML}}}\rightarrow\theta$. The sensitivity of the
maximum likelihood estimator to small fluctuations in $\mathbf{x}$ is
illustrated in Fig. (4). Additionally, the MSE of the maximum likelihood
estimator tends to saturate the CRB as it becomes more unbiased [kay1993,
van2000]. It is important to remark that there is no general formula to
determine an appropriate value of $N$. However, within the framework of
quantum metrology, unknown parameters are encoded into quantum resources;
because of the abundance of these resources the issue of small $N$ is often
ignored.

####  4 Example: Biased Coin

Consider a biased coin, which when flipped results in heads with an unknown
probability $p_{H}$ and tails with probability $1-p_{H}$. For the sake of
creating a locally optimized estimator, previous coin tosses suggest that the
bias is $p_{H}\approx p_{H,0}$. The FI of a single flip is easy to compute

$\mathcal{I}_{\text}{coin}=\frac{1}{p_{H}}\Big{(}\frac{\partial
p_{H}}{\partial
p_{H}}\Big{)}^{2}+\frac{1}{1-p_{H}}\Big{(}\frac{\partial(1-p_{H})}{\partial
p_{H}}\Big{)}^{2}=\frac{1}{p_{H}(1-p_{H})},$ (18)

thus the CRB imposes that the MSE of an unbiased estimator using $N$ outcomes
is bounded by

$\Delta^{2}\hat{p}_{H}\geq\frac{p_{H}(1-p_{H})}{N}.$ (19)

To remain somewhat general, the data collected is from $N$ coin tosses, $h$ of
which resulted in heads and $N-h$ of which resulted in tails, which occurs
with probability $\binom{N}{h}p_{H}^{h}(1-p_{H})^{N-h}$. Using the locally
optimized estimation strategy, Eq. eq:localestimator, the estimator is

$\hat{p}_{H}^{\text{Local}}=p_{H,0}+\frac{1}{N\mathcal{I}_{\text}{coin}}\frac{\partial\ln
p_{H}^{h}(1-p_{H})^{N-h}}{\partial p_{H}}\bigg{|}_{p_{H}\rightarrow
p_{H,0}}=\frac{h}{N},$ (20)

which is unbiased because

$\expval*{\hat{p}_{H}^{\text{Local}}}=\sum_{h=0}^{N}\binom{N}{h}p_{H}^{h}(1-p_{H})^{N-h}\frac{h}{N}=p_{H}.$
(21)

Furthermore, the estimator is efficient because it saturates the CRB

$\Delta^{2}\hat{p}_{H}^{\text{Local}}=\sum_{h=0}^{N}\binom{N}{h}p_{H}^{h}(1-p_{H})^{N-h}\big{(}\frac{h}{N}-p_{H})^{2}=\frac{p_{H}(1-p_{H})}{N}=\frac{1}{N\mathcal{I}_{\text}{coin}}.$
(22)

Despite the fact that the estimator was initially constructed using a local
approximation, the estimator is independent of $p_{H,0}$, and is thus globally
optimized. In addition, the same estimator is realized using the maximum
likelihood estimation strategy, see Fig. (4). The biased coin exemplifies the
underlying nature of the frequentist approach: as $N$ increases, the quantity
$\frac{h}{N}$ converges to the quantity $p_{H}$, from which the (albeit
simple) probability density function can be reverse engineered.

####  5 The Bayesian Approach

The Bayesian approach is typically used to estimate unknown parameters which
are stochastic. In other words, the latent parameters are themselves a random
variable and have an intrinsic probability distribution $p(\theta)$ \- which
should to be confused with $p(\mathbf{x}|\theta)$. Therefore, the observed
data $\mathbf{x}$ is dependent on specific realisations of $\theta$.
Consequently, a well-constructed estimator within the Bayesian approach aims
to minimize the MSE for all values (global) of $\theta$, and not subjected to
local values like the frequentist approach. To achieve this, the Bayesian
approach minimizes the average of a cost function $C(\hat{\theta},\theta)$
[kay1993, trees2007]

$\expval*{C(\hat{\theta},\theta)}=\int p(\theta)\bigg{(}\int
p(\mathbf{x}|\theta)\big{(}C(\hat{\theta}(\mathbf{x}),\theta)\mathrm{d}\mathbf{x}\bigg{)}\mathrm{d}\theta.$
(23)

In principle, a cost function is a generalisation of the MSE for the
frequentist approach. It is a function which decreases as $\hat{\theta}$
approaches $\theta$. The MSE is an example of a cost function, so too is the
absolute error $C=|\hat{\theta}-\theta|$. Different cost functions are
tailored to specific probability density functions to take advantage of
specific symmetries or properties.

By merging the two probability distributions, the average cost can be
interpreted as an average over the simultaneous realisations of $X$ and
$\theta$. According to Bayes’ theorem (hence the name of this approach), the
joint probability distribution can be interpreted in two ways

$p(\mathbf{x},\theta)=p(\mathbf{x}|\theta)p(\theta)=p(\theta|\mathbf{x})p(\mathbf{x}),$
(24)

thus the average cost can be written as

$\expval*{C(\hat{\theta},\theta)}=\int
p(\mathbf{x})\bigg{(}\int\mathrm{d}\theta
p(\theta|\mathbf{x})C(\hat{\theta}(\mathbf{x}),\theta)\mathrm{d}\theta\bigg{)}\mathrm{d}\mathbf{x}.$
(25)

The average cost can then be minimized through standard optimization
techniques, i.e by solving the equation

$\frac{\partial}{\partial\hat{\theta}}\int
p(\theta|\mathbf{x})C(\hat{\theta}(\mathbf{x}),\theta)\mathrm{d}\theta=0,$
(26)

where the quantity $p(\theta|\mathbf{x})$ can be computed using Bayes’ theorem

$p(\theta|\mathbf{x})=\frac{p(\mathbf{x}|\theta)p(\theta)}{p(\mathbf{x})}=\frac{p(\mathbf{x}|\theta)p(\theta)}{\int
p(\mathbf{x}|\theta)p(\theta)\mathrm{d}\theta}.$ (27)

A priori knowledge of $p(\theta)$ is needed to evaluate Eq.
eq:bayesianestimator, which is why the Bayesian approach is often used in
tandem with adaptive techniques. The estimator continually outputs a new
probability density function $p(\theta)$ based on the previous density
function and collected data, and as the number of repetitions increases it
will converge towards the correct value. There are precision bounds similar to
the CRB within the Bayesian framework, but they are dependent on the cost
function [bobrovsky1987]. More information about Bayesian inference can be
found in [trees2007].

###  2 Quantum Estimation Theory

In the quantum setting, the foundations of the parameter estimation problem
remains mostly unchanged from the classical setting [helstrom1969,
holevo1982]. An unknown parameter $\theta$ governs an $n$ qubit quantum state
$\rho_{\theta}$, the individual qubits can be measured with respect to a PVM
$M$, and the measurement outcomes $m_{1},\ldots,m_{n}$ are used to construct
an estimate $\hat{\theta}$222The assumptions that the qubits are acted on
independently and identically (both the encoding and the measurement) are
unnecessary and impose a limit on the most general framework of a quantum
parameter estimation scheme, see Fig. (7). These assumptions are introduced to
provide a natural extension from a classical framework to a quantum
framework.. The main difference from the classical setting is that the
measurement outcomes (analogous to $\mathbf{x}$) are not necessarily
independent from each other because of entanglement. As a result, estimates
can be made with a super-classical precision known as the Heisenberg limit.

.49  .49

Figure 5: Phase estimation. Figure 6: General framework. Figure 7: Diagrams of
the prepare, encode and measure segment of a quantum parameter estimation
problem. (a) The canonical example of quantum metrology is phase estimation
[giovannetti2006], in which a phase $\theta$ is independently and identically
encoded into each of the $n$ qubits of $\rho$ through a unitary $U_{\theta}$.
Each of the qubits are individually measured in accordance with a PVM $M$. (b)
A more general framework involves $\theta$ being encoded through a general
CPTP map $\Lambda_{\theta}$, which does not necessarily act identically on the
$n$ qubits. Additionally, the PVM $M$ is replaced by a POVM $\mathcal{M}$. The
generalized setting depicted in (b) is allowable in the realm of quantum
mechanics, but unlike the problem of phase estimation, it is difficult to
compare to classical setting. Further, highly entangling operations and
measurements are not feasible for current quantum technologies
[calsamiglia2001], so it is often more practical to consider the simplistic
setting of phase estimation as a benchmark for quantum metrology.

A quantum parameter estimation problem can be viewed as a two step process.
The first is the ‘prepare, encode and measure’ step, which is inherently
quantum by construction and depicted in Fig. (7). The second is the
statistical inference step, which is uniquely classical, thus the techniques
discussed in the the previous section can be applied. Therefore, using a
frequentist approach with an unbiased estimator, if the quantum portion is
repeated $\nu$ times, the MSE is bounded by a quantum version of the CRB,
otherwise known as the quantum Cramér-Rao bound (QCRB)

$\Delta^{2}\hat{\theta}\geq\frac{1}{\nu\mathcal{I}(\rho_{\theta},\mathcal{M})}\geq\frac{1}{\nu\mathcal{Q(\rho_{\theta})}},$
(28)

where $\mathcal{Q}$ is the quantum Fisher information (QFI), which is the FI
maximized over all POVM’s $\mathcal{M}$ [braunstein1994]. Evidently, the goal
of finding an optimal estimator $\hat{\theta}$ naturally divides into a
classical goal and a quantum goal. The classical goal is to devise an optimal
estimation technique, e.g. a locally optimized estimator or the maximum
likelihood estimator, whilst the quantum goal is to find an optimal
combination of initialized states $\rho$ and POVM $\mathcal{M}$. For the task
of phase estimation, this the QCRB can be saturated using highly entangled
states, such as the GHZ state or NOON states, and a local measurement strategy
[giovannetti2004]. In general, the QFI is a highly non-linear equation, and
there is no universal optimization strategy which is applicable to an
arbitrary encoding $\Lambda_{\theta}$. There are different mathematical
techniques to approximately solve this optimization problem [gill2013,
koczor2020, meyer2021].

The quantum metrology schematics in Fig. (7) are idealized settings. In
reality, it is much more complicated: environmental decoherence occurs in
simultaneity with the parameter encoding, resulting in noisy measurement
statistics and added uncertainty [escher2011a, escher2011b, demkowicz2012].
More so, quantum technologies are not perfect, and an error may be introduced
in either the quantum state preparation step or quantum measurement step. This
more realistic noisy scenario is explored in greater detail in Chapter 5.

####  1 Inferring an Estimate from an Observable

A simple frequentist estimation strategy used in quantum metrology is to
construct an estimator for the expectation value of an observable $O$ and
infer the value of the latent parameter from this estimate [toth2014].
Assuming that $O$ is chosen appropriately, the expectation value
$\expval{O}=\Tr(O\rho_{\theta})$ will be a function of $\theta$, denoted by
$f(\theta)=\Tr(O\rho_{\theta})$. An estimate of $f(\theta)$, $\hat{f}$, can be
inverted to obtain $\hat{\theta}=f^{-1}(\hat{f})$. The estimator $\hat{f}$ is
designed using the frequentist philosophy: with sufficient data $\nu\gg 1$,
the frequency of the measurement results will mimic the true probability
density function. Denote the eigenvalues of $O$ as $\\{\lambda_{j}\\}$ with
corresponding eigenvectors $\\{\ket{\phi_{j}}\\}$

$O=\sum_{j}\lambda_{j}\dyad{\phi_{j}}.$ (29)

The state $\rho_{\theta}$ is measured with respect to the eigenbasis of $O$.
The results are recorded as $m_{1},\ldots,m_{\nu}$: if the $k$th measurement
results in $\ket{\phi_{j}}$, then $m_{k}=\lambda_{j}$ and the maximum
likelihood estimate can be written as

$\hat{f}=\frac{1}{\nu}\sum_{k=1}^{\nu}m_{k}.$ (30)

This is an unbiased estimate because

$\mathbb{E}(\hat{f})=\frac{1}{\nu}\sum_{k=1}^{\nu}\mathbb{E}(m_{k})=\frac{1}{\nu}\sum_{k=1}^{\nu}\sum_{j}\lambda_{j}\Tr(\rho_{\theta}\dyad{\phi_{j}})=\frac{1}{\nu}\sum_{k=1}^{\nu}\expval{O}=\expval{O},$
(31)

and the MSE is proportional to the variance of $O$

$\Delta^{2}\hat{f}=\frac{\Delta^{2}O}{\nu}=\frac{\Tr(O^{2}\rho_{\theta})-\Tr(O\rho_{\theta})^{2}}{\nu}.$
(32)

An issue with this estimation technique is that $f(\theta)$ is not necessarily
an invertible function, and thus $\hat{f}=f^{-1}(\hat{f})$ may be ambiguous.
That is of course, unless one has a priori knowledge of
$\theta\approx\theta_{0}$ such that one can properly define a local inverse in
the region surrounding $f(\theta_{0})$. Assuming this is true and that the MSE
is small, $\Delta^{2}\hat{f}\ll 1$, then by the central limit theorem
$\hat{f}$ fluctuates close to $\expval{O}$, validating the first order Taylor
approximation

$\hat{\theta}=f^{-1}(\hat{f})\approx f^{-1}(\expval{O})+\frac{\partial
f^{-1}(\hat{f})}{\partial\hat{f}}\Big{|}_{\hat{f}\rightarrow\expval{O}}(\hat{f}-\expval{O})=\theta+\frac{1}{\frac{\partial\expval{O}}{\partial\theta}}(\hat{f}-\expval{O}).$
(33)

It follows from the above approximation that the estimator $\hat{\theta}$ is
unbiased and has MSE

$\Delta^{2}\hat{\theta}=\frac{\Delta^{2}\hat{f}}{|\frac{\partial\expval{O}}{\partial\theta}|^{2}}=\frac{\Delta^{2}O}{\nu|\frac{\partial\expval{O}}{\partial\theta}|^{2}}.$
(34)

Figure 8: Graphical calculation of the MSE $\Delta^{2}\hat{\theta}$ using the
error propagation formula. The solid red curve depicts
$\hat{\theta}=f^{-1}(\hat{f})$, which at the point
$\hat{f}\rightarrow\expval{O}$ has a tangent with angle $\alpha$, therefore,
$\frac{\Delta\hat{\theta}}{\Delta\hat{f}}=|\tan\alpha|=|\frac{\partial\expval{O}}{\partial\theta}|^{-1}$.

Eq. eq:ErrorPropagation is the error propagation formula, which quantifies the
amount that $\hat{\theta}$ fluctuates around $\theta$ in terms of the
fluctuations of $\hat{f}$ around $\expval{O}$ [ku1966]. A geometric intuition
of the formula is depicted in Fig. (8). The term in the denominator
$|\frac{\partial\expval{O}}{\partial\theta}|^{2}$ encapsulates the difficulty
of inverting a function when there is uncertainty. The effects of uncertainty
are amplified near a local maxima or minima, but diminish as
$|\frac{\partial\expval{O}}{\partial\theta}|\rightarrow\infty$.

####  2 Quantum Fisher Information

Using the semantics of quantum information theory, the explicit expression for
the FI with respect to a POVM $\mathcal{M}$ with outcomes $\\{E_{m}\\}$ can be
written as

$\mathcal{I}(\rho_{\theta},\mathcal{M})=\int\frac{\big{(}\Tr(E_{m}\dot{\rho}_{\theta})\big{)}^{2}}{\Tr(E_{m}\rho_{\theta})}\mathrm{d}m,$
(35)

where the notation $\dot{\square}=\frac{\partial\square}{\partial\theta}$ is
used for conciseness. Just as the FI is interpreted as an information measure,
so too is the QFI [barndorff2000]. Eq. eq:AngularDistance suggests that the
POVM which maximizes the distinguishability between the probability density
functions associated to $\rho_{\theta}$ and $\rho_{\theta+\delta\theta}$ will
similarly maximize the FI. This is a principle idea behind the derivation of
the closed form expression of the QFI [braunstein1994].

The derivation begins by defining the superoperator

$\mathcal{R}_{\rho_{\theta}}(O)=\frac{1}{2}(\rho_{\theta}O+O\rho_{\theta}),$
(36)

whose inverse333The inverse $\mathcal{R}_{\rho_{\theta}}^{-1}(O)$ is not
always well defined for all $O$, however the quantity used in the derivation
of the QFI, $\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})$, always
converges to a well-defined Hermitian operator. is

$\mathcal{R}_{\rho_{\theta}}^{-1}(O)=\sum_{j,k}\frac{2}{\lambda_{j}+\lambda_{k}}O_{jk}\dyad{j}{k},$
(37)

where $\rho_{\theta}=\sum_{j}\lambda_{j}\dyad{j}$ is the orthonormal expansion
of $\rho_{\theta}$ and $O_{jk}=\mel{j}{O}{k}$. A property of $\mathcal{R}$ is
that for any Hermitian $A$ and $B$,
$\Tr\big{(}AB\big{)}=\text{Re}\big{[}\Tr\big{(}\rho_{\theta}A\mathcal{R}_{\rho_{\theta}}^{-1}(B)\big{)}\big{]}$,
from which it follows that the FI can be written as

$\split\mathcal{I}(\rho_{\theta},\mathcal{M})&=\int\frac{\text{Re}\big{[}\Tr\big{(}\rho_{\theta}E_{m}\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})\big{)}\big{]}^{2}}{\Tr(E_{m}\rho_{\theta})}\mathrm{d}m\\\
\leq\int\frac{\big{|}\Tr\big{(}\rho_{\theta}E_{m}\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})\big{)}\big{|}^{2}}{\Tr(E_{m}\rho_{\theta})}\mathrm{d}m\\\
=\int\frac{\big{|}\Tr\big{(}\sqrt{\rho_{\theta}}\sqrt{E_{m}}\sqrt{E_{m}}\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})\sqrt{\rho_{\theta}}\big{)}\big{|}^{2}}{\Tr(E_{m}\rho_{\theta})}\mathrm{d}m.$
(38)

The final step in the derivation uses the Cauchy-Schwarz inequality
$|\Tr(A^{\dagger}B)|^{2}\leq\Tr(AA^{\dagger})\Tr(BB^{\dagger})$ with
$A=\sqrt{E_{m}}\sqrt{\rho_{\theta}}$ and
$B=\sqrt{E_{m}}\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})\sqrt{\rho_{\theta}}$,

$\split\mathcal{I}(\rho_{\theta},\mathcal{M})&\leq\int\Tr\big{(}E_{m}\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})\rho_{\theta}\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})\big{)}\mathrm{d}m\\\
=\Tr\big{(}\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})\rho_{\theta}\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})\big{)}\\\
=\mathcal{Q}(\rho_{\theta}).$ (39)

The Hermitian operator $\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})$
is the symmetric logarithmic derivative. The QCRB can be saturated by setting
$\mathcal{M}$ to be the measurement in the eigenbasis of
$\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})$ [braunstein1994,
luo2000, matsumoto2002]. Unfortunately, but not surprisingly, such a
measurement is encumbered by the usual quandary of the frequentist approach:
the measurement basis is dependent on $\theta$444Similar to how a locally
optimized estimator, Eq. eq:localestimator, approximately saturates the CRB,
measuring in the eigenbasis of
$\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})|_{\theta\rightarrow\theta_{0}}$
approximately saturates the QCRB.. Furthermore, this measurement strategy is
very sophisticated and out of reach for current technologies
[calsamiglia2001]. Fortunately, this is not the unique measurement strategy
which saturates the QCRB [giovannetti2004]. As mentioned, the quantum goal of
parameter estimation problems is to determine feasible measurement schemes
which best saturate the QCRB.

A closed form expression for the QFI can be derived using the definition of
$\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})$, Eq. eq:Rinverse,

$\mathcal{Q}(\rho_{\theta})=\sum_{j}\frac{\dot{\lambda}_{j}^{2}}{\lambda_{j}}+2\sum_{j,k}\frac{(\lambda_{j}-\lambda_{k})^{2}}{\lambda_{j}+\lambda_{k}}|\braket{\dot{j}}{k}|^{2}.$
(40)

The first sum is reminiscent of the classical FI and quantifies the amount of
extractable information from the eigenvalues $\\{\lambda_{j}\\}$. Whilst the
second sum accounts for quantum effects such as superposition and entanglement
and quantifies the amount of extractable information from the quantum states
$\\{\ket{j}\\}$. To a certain extent, the classical term is limited to
‘amplitudes’, while the quantum term has access to ‘amplitudes’ and ‘phases’.
As such, the quantum term is significantly more influential than the classical
term, this is reinforced by the convexity property of the QFI [alipour2015]

$\mathcal{Q}\big{(}p\rho_{1}+(1-p)\rho_{2}\big{)}\leq
p\mathcal{Q}(\rho_{1})+(1-p)\mathcal{Q}(\rho_{2}).$ (41)

For the special case of pure states $\rho_{\theta}=\dyad{\psi_{\theta}}$, the
expression is much more aesthetically pleasing. It follows from
$\rho_{\theta}^{2}=\rho_{\theta}$ that
$\dot{\rho}_{\theta}=\rho_{\theta}\dot{\rho}_{\theta}+\dot{\rho}_{\theta}\rho_{\theta}$
and thus
$\mathcal{R}_{\rho_{\theta}}^{-1}(\dot{\rho}_{\theta})=2\dot{\rho}_{\theta}$.
The QFI simplifies greatly to

$\mathcal{Q}(\ket{\psi_{\theta}})=4\Tr(\rho_{\theta}\dot{\rho}_{\theta}^{2})=4\big{(}\braket*{\dot{\psi}_{\theta}}-|\braket*{\dot{\psi}_{\theta}}{\psi_{\theta}}|^{2}\big{)}.$
(42)

In fact, it was shown that a similar expression holds for arbitrary mixed
states [escher2011b]

$\mathcal{Q}(\rho_{\theta})=\min_{\ket{\Psi_{\theta}}}4\big{(}\braket*{\dot{\Psi}_{\theta}}-|\braket*{\dot{\Psi}_{\theta}}{\Psi_{\theta}}|^{2}\big{)},$
(43)

where the minimization is taken over all possible purifications, Eq.
eq:purification, of $\rho_{\theta}$.

####  3 Geometric Perspectives of the QFI

The representation introduced in Chapter 2 is that a quantum state
$\ket{\psi}$ can be thought of as a vector which is an element of a Hilbert
space $\mathcal{H}$. An alternative to this is a geometric representation,
where $n$ qubit quantum states are thought to be elements of the complex
projective space $\mathbb{CP}^{n}$ [wootters1981, petz1996, grabowski2005].
Pure states reside on the surface of this Riemannian manifold and mixed states
in the interior, the $n=1$ case is the well-known Bloch sphere portrayed in
Fig. (LABEL:fig:blochsphere). $\mathbb{CP}^{n}$ is equipped with an
infinitesimal metric called the Fubini-Study metric $ds^{2}$, which is called
the Bures metric [bures1969, sommers2003] when it is extended to include the
interior. Such a metric allows one to compare neighbouring quantum states
$\rho_{\theta}$ and $\rho_{\theta+\delta\theta}$, analogous to the FI metric
for (classical) statistical manifolds, it can be shown that
$ds^{2}=\frac{1}{4}\mathcal{Q(\rho_{\theta})}\delta\theta^{2}$ [facchi2010,
sidhu2020].

Figure 9: The length of the chord between $\ket{\psi_{\theta}}$ and
$\ket{\psi_{\theta+\delta\theta}}$ is $2\sin\frac{d_{A}}{2}$. For small
$d_{A}$, the length of the chord is approximately equal to the length of the
geodesic, $ds\approx 2\sin\frac{d_{A}}{2}\approx d_{A}$. This idea generalizes
to higher dimensional abstract surfaces $\mathbb{CP}^{n}$ as well as their
interior points (mixed states).

The Bures angle $d_{A}$ is the angle between the rays of $\rho_{1}$ and
$\rho_{2}$, explicitly [amari2016, bengtsson2017]

$d_{A}(\rho_{1},\rho_{2})=\arccos\sqrt{\mathscr{F}(\rho_{1},\rho_{2})},$ (44)

where $\mathscr{F}$ is the fidelity, Eq. eq:fidelity. For neighbouring quantum
states, the Bures angle can be approximated two different ways. The first way
is by using a first order Taylor expansion

$d_{A}(\rho_{\theta},\rho_{\theta+\delta\theta})=\sqrt{2-2\sqrt{\mathscr{F}(\rho_{\theta},\rho_{\theta+\delta\theta})}}+\mathcal{O}(\delta\theta^{2}).$
(45)

The second is a geometric approximation using the Bures metric (and by
extension the QFI), the intuition of which is given in Fig. (9)

$d_{A}(\rho_{\theta},\rho_{\theta+\delta\theta})=ds+\mathcal{O}(\delta\theta^{2})=\frac{1}{2}\sqrt{\mathcal{Q}(\rho_{\theta},\rho_{\theta+\delta\theta})}\delta\theta+\mathcal{O}(\delta\theta^{2}).$
(46)

A new expression for the QFI is obtained by merging the two equations
[sidhu2020]

$\mathcal{Q}(\rho_{\theta})=\lim_{\delta\theta\rightarrow
0}8\frac{1-\sqrt{\mathscr{F}(\rho_{\theta},\rho_{\theta+\delta\theta})}}{\delta\theta^{2}},$
(47)

which can be useful to derive analytic bounds for the QFI and other
information theoretic quantities [suzuki2019, tsang2020]. A corollary of Eq.
eq:QFIlimit is the concavity of the QFI under CPTP maps $\mathcal{E}$

$\mathcal{Q}\big{(}\mathcal{E}(\rho_{\theta})\big{)}\leq\mathcal{Q}\big{(}\rho_{\theta}\big{)},$
(48)

which follows from the monotonicity of the fidelity
$\mathscr{F}\big{(}\mathcal{E}(\rho_{1}),\mathcal{E}(\rho_{2})\big{)}\geq\mathscr{F}\big{(}\rho_{1},\rho_{2}\big{)}$.
If $\mathcal{E}$ is thought of as an interaction with an environment (Chapter
5) or a malicious adversary (Chapter 6), then the concavity of the QFI can be
understood as information about $\theta$ being lost to these outside sources.

####  4 Ultimate Precision: The Heisenberg Limit

To recapitulate: the CRB is a bound on the MSE by optimizing over estimation
strategies, and the QCRB extends the bound by optimizing over measurement
strategies. The next natural extension is to optimize over initialized quantum
states, to find the true limit of precision attainable through quantum
mechanics. The upper bound for which is referred to as the Heisenberg Limit
(HL) [yurke1986, holland1993].

Originally, the HL was derived within the framework of phase estimation. In
the phase estimation problem, a phase $\theta$ is encoded into each qubit of
an $n$ qubit pure state $\ket{\psi}$ by a unitary $U_{\theta}=e^{-i\theta H}$,
where the Hamiltonian $H$ acts independently and identically on all $n$
qubits. The QFI can be calculated to be

$\mathcal{Q}=4\big{(}\mel*{\psi}{H^{2}}{\psi}-|\mel*{\psi}{H}{\psi}|^{2}\big{)}=4\Delta^{2}H.$
(49)

The etymology of the term ‘Heisenberg limit’ stems from the fact that the QCRB
(with $\nu=1$) can be manipulated to mimic the the Heisenberg uncertainty
principle

$\Delta^{2}\hat{\theta}\Delta^{2}H\geq\frac{1}{4}.$ (50)

The QFI for phase estimation can be maximized by setting $\ket{\psi}$ to be a
highly entangled state, such as the GHZ state for qubit systems or the NOON
state for photonic systems, which results in $\mathcal{Q}=n^{2}$. Hence, the
ultimate allowable precision by quantum mechanics (the HL) is

$\Delta^{2}\hat{\theta}_{\text}{HL}=\frac{1}{\nu n^{2}}.$ (51)

The HL offers a quadratic improvement compared to the standard quantum limit
(SQL), where the $\ket{\psi}$ is limited to separable states

$\Delta^{2}\hat{\theta}_{\text}{SQL}=\frac{1}{\nu n}.$ (52)

The SQL is also referred to as the classical limit or the shot-noise limit
[xiao1987].

For qubit (and qudit) systems555For CV systems a quantum advantage can be
achieved with squeezing [yurke1986, ono2010]., entanglement is a crucial
resource for quantum metrology [pezze2009, pezze2018]. In fact, the quadratic
tendencies of the QFI of a quantum state for phase estimation can be bounded
with respect to the geometric measure of entanglement $G$666The geometric
measure of entanglement for a pure state $\ket{\psi}$ is
$G(\ket{\psi})=1-\max_{\ket{\phi}}|\braket{\phi}{\psi}|^{2}$, where
$\ket{\phi}$ is maximized over all fully separable states. The definition is
extended to mixed states by finding the convex roof of the geometric measure
of entanglement over all possible statistical ensembles [wei2003].
[augusiak2016]

$\mathcal{Q}(\rho_{\theta})\leq n+8n^{2}\sqrt{G(\rho_{\theta})}.$ (53)

It is worth stressing that entanglement may be a necessary condition to
surpass the SQL but it is not a sufficient condition [hyllus2010,
oszmaniec2016]. Additionally, the bounds in Eq. eq:HL and Eq. eq:SQL are
exclusive to the problem of phase estimation with an iid encoding. The QFI can
surpass $n^{2}$ for non-linear $H$ [luis2004, boixo2007, choi2008, braun2018],
and scenarios can be devised in which entanglement is not a necessary resource
[tilma2010].

####  5 Bayesian Approach to Quantum Metrology

In the quantum version of the frequentist approach, the MSE is minimized by
optimizing over all possible POVM’s and input quantum states. The quantum
version of the Bayesian approach [holevo1982, jarzyna2015, rubio2018] is
enhanced in an analogous fashion. As the estimator is updated adaptively, so
too can the initialized quantum state as well as choice of POVM.

For parameter estimation problems which exhibit periodicity, such as phase
estimation, the circular cost function

$C_{\circ}(\hat{\theta},\theta)=4\sin^{2}\Big{(}\frac{\hat{\theta}-\theta}{2}\Big{)}$
(54)

is a natural choice as a figure of merit [demkowicz2011, demkowicz2015], and
converges to the MSE as $\hat{\theta}$ approaches $\theta$. If the initial
choice of input quantum state and POVM are
$\rho_{\theta}=U_{\theta}\rho_{0}U_{\theta}^{\dagger}$ and $\int
E_{m}\mathrm{d}m$ respectively, then the average cost is

$\expval{C_{\circ}}=\int
p(\theta)\bigg{(}\int\Tr(\rho_{\theta}E_{m})C_{\circ}(\hat{\theta}(m),\theta)\mathrm{d}m\bigg{)}\mathrm{d}\theta,$
(55)

which is invariant when replacing the POVM $\\{E_{m}\\}$ with a covariant POVM
$\\{E_{\hat{\theta}}\\}$ [holevo1982, derka1998, chiribella2004,
chiribella2005]

$E_{\hat{\theta}}=U_{\hat{\theta}}\Sigma U_{\hat{\theta}}^{\dagger}$ (56)

and $\Sigma$ is the positive-semi definite operator defined for a specific
$\hat{\theta}$

$\Sigma=\int
U_{\hat{\theta}(m)}^{\dagger}E_{m}U_{\hat{\theta}(m)}\mathrm{d}m.$ (57)

This re-parametrization allows the average cost to be expressed as

$\split\expval{C_{\circ}}&=\int\int
p(\theta)\Tr(\rho_{\theta}E_{\hat{\theta}})C_{\circ}(\hat{\theta},\theta)\mathrm{d}\hat{\theta}\mathrm{d}\theta\\\
=\int
p(\theta)\Tr(\rho_{\theta}\Sigma)4\sin^{2}\frac{\theta}{2}\mathrm{d}\theta.$
(58)

By optimizing the above expression, the initialized quantum state $\rho_{0}$
and POVM characterized by $\Sigma$ can be updated adaptively [derka1998,
chiribella2005].

As mentioned, the Bayesian statistical inference approach addresses the issues
inherent to the frequentist approach: lack of a priori knowledge
[kolodynski2010, demkowicz2011] and limited resources [rubio2020b]. The work
presented in the subsequent chapters exclusively focus on the frequentist
approach, as such, an interesting future perspective would be to generalize
some of the findings to the Bayesian approach. Specifically in Chapter 6,
where the estimation process is adapted in some capacity to account for the
cryptographic framework.

####  6 Multiple Parameters

In the interest of simplicity, this chapter introduced the problem of
parameter estimation with a single unknown parameter. The problem naturally
generalizes to include multiple latent parameters
$\theta\rightarrow\boldsymbol{\theta}=\\{\theta_{1},\ldots,\theta_{m}\\}$,
where the goal extends to devising estimators for each parameter
$\hat{\theta}(\mathbf{x})\rightarrow\hat{\boldsymbol{\theta}}(\mathbf{x})=\\{\hat{\theta}_{1}(\mathbf{x}),\ldots,\hat{\theta}_{m}(\mathbf{x})\\}$.
There are two major quandaries which arise in the multiparameter setting.
First, the parameters may be statistically dependent on one another, which
adds ambiguity when trying to interpret the observed data $\mathbf{x}$.
Second, it is not always possible to simultaneously construct an efficient
estimator for each unknown parameter.

Within the frequentist inference framework, assuming that each estimator
satisfies the unbiased estimator constraint,
$\mathbb{E}(\hat{\boldsymbol{\theta}})=\boldsymbol{\theta}$, the
generalization of the QCRB is the matrix equation [yuen1973, helstrom1974,
tsang2011]

$\textbf{Cov}(\boldsymbol{\theta})\geq\frac{1}{\nu}\boldsymbol{\mathcal{I}}^{-1}(\boldsymbol{\theta})\geq\frac{1}{\nu}\boldsymbol{\mathcal{Q}}^{-1}(\boldsymbol{\theta}),$
(59)

where $\textbf{Cov}(\boldsymbol{\theta})$ is the covariance matrix with
entries
$\textbf{Cov}(\boldsymbol{\theta})_{i,j}=\expval*{(\hat{\theta}_{i}-\theta_{i})(\hat{\theta}_{j}-\theta_{j})}$,
$\boldsymbol{\mathcal{I}}(\boldsymbol{\theta})$ is the FI matrix

$\boldsymbol{\mathcal{I}}(\boldsymbol{\theta})_{i,j}=\int
p(\mathbf{x}|\boldsymbol{\theta})\Big{(}\frac{\partial\ln
p(\mathbf{x}|\boldsymbol{\theta})}{\partial\theta_{i}}\Big{)}\Big{(}\frac{\partial\ln
p(\mathbf{x}|\boldsymbol{\theta})}{\partial\theta_{j}}\Big{)}\mathrm{d}\mathbf{x},$
(60)

and $\boldsymbol{\mathcal{Q}}(\boldsymbol{\theta})$ is the QFI matrix

$\boldsymbol{\mathcal{Q}}(\boldsymbol{\theta})_{i,j}=\Tr\Big{(}\mathcal{R}^{-1}_{\rho_{\boldsymbol{\theta}}}\big{(}\frac{\partial\rho_{\boldsymbol{\theta}}}{\partial\theta_{i}}\big{)}\rho_{\boldsymbol{\theta}}\mathcal{R}^{-1}_{\rho_{\boldsymbol{\theta}}}\big{(}\frac{\partial\rho_{\boldsymbol{\theta}}}{\partial\theta_{j}}\big{)}\Big{)}.$
(61)

The diagonal elements of an invertible positive semi-definite matrix $M$
satisfy $M_{i,i}\geq 1/M_{i,i}$, and equality holds for each $i$ when $M$ is
diagonal. Hence, when the $m$ parameters are statistically independent from
each other, Eq. eq:multiQCRB reduces to the QCRB for each individual parameter
$\theta_{i}$. Multiparameter estimation sustains additional complications in
the quantum context with respect to optimizing quantum states and
measurements. The superoperator
$\mathcal{R}^{-1}_{\rho_{\boldsymbol{\theta}}}\big{(}\frac{\partial\rho_{\boldsymbol{\theta}}}{\partial\theta_{i}}\big{)}$
is the symmetric logarithmic derivative with respect to $\theta_{i}$, and
measuring in all of these bases will saturate the QCRB. However, these
measurements may not be compatible and thus cannot be realized in simultaneity
[vidrighin2014, crowley2014, ragy2016].

Despite the additional complexity of simultaneously estimating multiple
parameters, multiparameter quantum metrology is an active research topic
[szczykulska2016, nichols2018, albarelli2019, rubio2020b, meyer2021]. With
respect to phase estimation, when the phase encoding unitaries do not commute,
it is more efficient to estimate them in simultaneity rather than
independently [baumgratz2016]; when the encoding unitaries do commute, one can
devise a simultaneous estimation strategy which is at least as efficient as
estimating the phases independently [humphreys2013]. Multiparameter quantum
metrology is a natural framework for eigenvalue estimation of higher
dimensional unitaries [fujiwara2001, ballester2004, berry2015, baumgratz2016]
and for spatially distributed estimation problems [eldredge2018, ge2018,
proctor2018, zhuang2018, rubio2020a, guo2020]. The subsequent research
chapters focus on the single parameter setting. Nonetheless, the mathematical
techniques and derivations can easily be adapted to the multiple parameter
setting. Formally addressing these generalizations is a future perspective of
the works presented.

###  3 Example Applications of Quantum Metrology

The final section of this chapter explores well-known applications of quantum
metrology where the concepts and tools that were introduced are put into
practise. The examples chosen, phase estimation and amplitude estimation, have
a simple mathematical formalism and highlight the novelty of a quantum
parameter estimation problem. Specifically, phase estimation, which is
indisputably the canonical usage of quantum metrology [caves1981,
giovannetti2004], clearly showcases the advantages a quantum system can
provide. Additionally, phase estimation is the core problem of Chapter 4 and
Chapter 5. The example of amplitude estimation, although not present in the
subsequent research chapters, is included to showcase a simple usage which is
not phase estimation.

####  1 Phase Estimation (Photonic Interferometry)

Phase estimation is a benchmarking problem for quantum metrology
[holland1993]. It encapsulates a variety of applications for quantum sensors,
notably magnetometry [taylor2008, wasilewski2010, sewell2012, brask2015,
razzoli2019], frequency estimation for optomechanical sensors [zheng2016,
djorwe2019, tsang2013], spectroscopy [meyer2001, leibfried2004, kira2011,
dorfman2016, shaniv2018] and unitary tomography [svore2013, o2019]. The
premise of phase estimation is simple, yet the results are an elegant display
of a quantum advantage. A relative phase is encoded into a quantum system via
a physical interaction, and using highly entangled states and a simple
measurement strategy, the QCRB can be saturated to attain the HL
[giovannetti2006].

Figure 10: Schematic of a Mach-Zehnder interferometer. A photonic source is
fired onto a 50:50 beam splitter, which reflects a photon with 50% probability
(path $a$), otherwise the photon is transmitted (path $b$). A phase shift of
$\theta$ is introduced uniquely on path $a$, after the paths interfere via a
second beam splitter. The quantum state is measured with photon-counters,
whose outcomes $j$ are dependent on the relative phase $\theta$.

For photonic sources, a relative phase can be introduced using a Mach-Zehnder
interferometer [zetie2000], depicted in Fig. 10. A photonic quantum state
passes through a 50:50 beam splitter, in which photons are either transmitted
or reflected, both with probability of one half, where reflection induces a
phase shift of $\pi/2$ [dowling2008]. After passing through the first beam
splitter, a single photon will be in a superposition of the two possible
modes, labelled by the respective paths $\ket{a}$ and $\ket{b}$. Since photons
are indistinguishable particles, it is customary to write an $n$ photon
quantum state as

$\ket{\psi}=\sum_{k=0}^{n}\alpha_{k}\ket{k,n-k},$ (62)

where the notation $\ket{k,n-k}$ denotes the quantum state with $k$ photons in
mode $\ket{a}$ and $n-k$ photons in mode $\ket{b}$. The path-dependent phase
shift is represented by the unitary

$U_{\theta}=\exp(-i\frac{\theta}{2}\big{(}\dyad{a}-\dyad{b}\big{)})^{\otimes
n},$ (63)

thus

$U_{\theta}\ket{k,n-k}=e^{-i\frac{2k-n}{2}\theta}\ket{k,n-k}.$ (64)

After the phase shifts, the two paths interfere by passing through a second
beam splitter. This causes the detection probabilities to be dependent on
$\theta$, from which the statistics can be used to generate an estimate
$\hat{\theta}$. The precision of the estimate is ultimately bounded by the QFI

$\mathcal{Q}(\ket{\psi})=\sum_{k=0}^{n}|\alpha_{k}|^{2}(2k-n)^{2}-\Big{(}\sum_{k=0}^{n}|\alpha_{k}|^{2}(2k-n)\Big{)}^{2}=4\sum_{k=0}^{n}|\alpha_{k}|^{2}k^{2}-4\Big{(}\sum_{k=0}^{n}|\alpha_{k}|^{2}k\Big{)}^{2}.$
(65)

Suppose that the source of photons is uncorrelated such that the quantum state
after passing through the first beam splitter is

$\ket{\psi}_{\text}{sep}=\frac{1}{\sqrt{2^{n}}}\Big{(}\ket{a}+\ket{b}\Big{)}^{\otimes
n}=\sum_{k=0}^{n}\sqrt{2^{-n}\binom{n}{k}}\ket{k,n-k}.$ (66)

Using Eq. eq:qfiEx1, the QFI can be computed to be equal to the SQL,
$\mathcal{Q}_{\text}{sep}=n$. In contrast, if the quantum state is initialized
in the NOON state [dowling2008, matthews2016, zhang2018]

$\ket{\psi}_{\text}{ent}=\frac{1}{\sqrt{2}}\Big{(}\ket{n,0}+\ket{0,n}\Big{)},$
(67)

then the HL is attained777This has been experimentally accomplished using
$n=4$ photons [nagata2007], $\mathcal{Q}_{\text}{ent}=n^{2}$. The QCRB can be
saturated using the previously described method of inferring an estimate from
an observable. The chosen observable is the parity of the detected photons in
detector $a$, labelled $D_{a}$ in, Fig. 10,

$O_{\text}{parity}=\sum_{k=0}^{n}(-1)^{k}\dyad{k,n-k}.$ (68)

Assuming no dark counts, if $D_{a}$ detects $j$ photons, $D_{b}$ is fixed to
$n-j$ photons, hence the parity of $D_{b}$ could have equally been considered.
After intercepting the second beam splitter, the quantum state (up to a global
phase) is

$\split\ket{\psi_{\theta}}_{\text}{ent}&=\frac{1}{\sqrt{2^{n}}}\sum_{k=0}^{n}\sqrt{\binom{n}{k}}\Big{(}e^{-in\theta/2}i^{n-k}+e^{in\theta/2}i^{k}\Big{)}\ket{k,n-k}\\\
=\frac{1}{\sqrt{2^{n-1}}}\sum_{k=0}^{n}\sqrt{\binom{n}{k}}\cos\Big{(}\frac{n\theta}{2}+\frac{\pi}{4}(2k-n)\Big{)}\ket{k,n-k},$
(69)

from which it can be computed that

$\expval*{O}_{\text}{ent}=\cos\big{(}n\theta-\frac{\pi n}{2}\big{)}.$ (70)

If the prepare and measure protocol is repeated $\nu$ times, the MSE of the
estimator is

$\Delta^{2}\hat{\theta}_{\text}{ent}=\frac{\Delta^{2}O_{\text}{ent}}{\nu\big{|}\frac{\partial\expval*{O}_{\text}{ent}}{\partial\theta}\big{|}^{2}}=\frac{1}{\nu
n^{2}},$ (71)

thus the QCRB is saturated with a simple measurement strategy. This example
highlights the achievability of a quantum advantage, but simultaneously the
locality of the frequentist approach. The periodicity of
$\expval*{O}_{\text}{ent}$ implies that a priori knowledge of $\theta$ is
required to an order of $2\pi/n$. If $\theta$ is completely unknown, a
frequentist approach can still be employed using varying $n=2^{m}$, where
successive rounds of estimation are used to estimate the $m$th (binary) digit
of $\theta$ [kitaev1995].

In this example, a relative phase is encoded using a photonic source by means
of a Mach-Zehnder interferometer. Notably though, analogous results are
obtained using spin systems [toth2014]. In fact, the NOON state, Eq. eq:N00N,
can be interpreted as a GHZ state

$\ket{\text{GHZ}}=\frac{1}{\sqrt{2}}(\ket{0}^{\otimes n}+\ket{1}^{\otimes
n}),$ (72)

which when used for phase estimation, similarly saturates the HL.

####  2 Amplitude Estimation (Thermometry)

After phase estimation, the next obvious example of a quantum metrology
problem is amplitude estimation. A well-known example of amplitude estimation
is quantum thermometry [correa2015, de2018, mehboudi2019], where the unknown
parameter in question is temperature. Temperature is a seemingly intuitive
notion ever present in our daily lives, and measuring this quantity may appear
trivial. For every day objects, an infra-red thermometer converts infra-red
radiation into a voltage which is converted into a temperature. This is done
extremely quickly and accurately. Even so, in the ‘very cold’ regime near zero
Kelvin, measuring temperature is a complicated task, but very necessary for
modern technologies such as superconductors.

Experimental implementations of quantum thermometry differ greatly
[neumann2013, toyli2013, kucsko2013], yet the underlying principle is
straightforward [mehboudi2019]. An $N$ level system interacts with an external
source at temperature $T$. Eventually, the collective ensemble of the system
and the bath will reach thermal equilibrium, and the state of the system is
given by the Gibbs ensemble

$\rho_{T}=\frac{1}{Z}e^{-\frac{H}{k_{B}T}},$ (73)

where $k_{B}$ is the Boltzmann constant, $H$ is the system Hamiltonian

$H=\sum_{k=1}^{N}\epsilon_{k}\dyad{\epsilon_{k}}$ (74)

with energy eigenvalues $\epsilon_{k}$ and eigenstates $\ket{\epsilon_{k}}$,
and $Z=\Tr e^{-\frac{H}{k_{B}T}}$ is the partition function. Within the
standard convention of setting $k_{B}=1$, the derivative of the quantum system
with respect to temperature is

$\dot{\rho}_{T}=\frac{H}{T^{2}}\rho_{T}-\frac{\expval*{H}}{T^{2}}\rho_{T}=\frac{1}{2T^{2}}\big{(}(H-\expval*{H})\rho_{T}+\rho_{t}(H-\expval*{H})\big{)}.$
(75)

From which it is clear that the symmetric logarithmic derivative is

$\mathcal{R}_{\rho_{T}}^{-1}(\dot{\rho}_{T})=\frac{1}{T^{2}}(H-\expval*{H}),$
(76)

therefore

$\mathcal{Q}(\rho_{T})=\Tr\big{(}\mathcal{R}_{\rho_{T}}^{-1}(\dot{\rho}_{T})\rho_{T}\mathcal{R}_{\rho_{T}}^{-1}(\dot{\rho}_{T})\big{)}=\frac{1}{T^{4}}\Tr\big{(}(H-\expval{H})^{2}\rho_{T}\big{)}=\frac{\Delta^{2}H}{T^{4}}.$
(77)

In fact, the heat capacity

$\frac{\partial\expval{H}}{\partial T}=\frac{\Delta^{2}H}{T^{2}}$ (78)

is directly proportional to the QFI [phillies1984], hence the QCRB can be re-
written as

$\Delta^{2}\hat{T}\geq\frac{1}{\nu\mathcal{Q}}=\frac{\Delta^{2}H}{\nu|\frac{\partial\expval{H}}{\partial
T}|^{2}}.$ (79)

Which suggests that the QCRB can be saturated by using the previously
described estimation strategy of inferring the temperature from an observable,
Eq. eq:ErrorPropagation. In this instance the observable is the energy of the
system [jahnke2011], $H$, which may not be surprising because of how
intertwined energy and temperature are as quantities within the realm of
statistical mechanics888On a macroscopic scale, temperature is an average
quantity of a system composed of many many particles. This definition is
somewhat ambiguous on a microscopic scale. In Eq. eq:thermalequilibrium,
temperature can be interpreted as a variable which governs the probability of
the quantum system occupying a specific energy eigenstate..

Analogous to a GHZ state or NOON state being the optimal probe for phase
estimation, Eq. eq:QFIthermometry can be maximized to find the optimal probe
for thermometry. The solution is an effective two level system with a single
eigenstate having an energy of $\epsilon_{-}$ and $N-1$ eigenstates having a
degenerate energy of $\epsilon_{+}$, the relative error of such a probe is
$\Delta^{2}\hat{T}/T^{2}=\mathcal{O}(1/\log N)$ [correa2015, mehboudi2019].
Note that the information presented in this example holds only for fully
thermalized systems, which may be a time consuming process. The analysis is
significantly more complex for partially thermalized systems [correa2015].