# Adaptive measurement filter: efficient strategy for optimal estimation of
quantum Markov chains

Alfred Godley and Mădălin Guţă  School of Mathematical Sciences, University of
Nottingham, University Park, NG82SD Nottingham, United Kingdom

###### Abstract

Continuous-time measurements are instrumental for a multitude of tasks in
quantum engineering and quantum control, including the estimation of dynamical
parameters of open quantum systems monitored through the environment. However,
such measurements do not extract the maximum amount of information available
in the output state, so finding alternative optimal measurement strategies is
a major open problem.

In this paper we solve this problem in the setting of discrete-time input-
output quantum Markov chains. We present an efficient algorithm for optimal
estimation of one-dimensional dynamical parameters which consists of an
iterative procedure for updating a ‘measurement filter’ operator and
determining successive measurement bases for the output units. A key
ingredient of the scheme is the use of a coherent quantum absorber as a way to
post-process the output after the interaction with the system. This is
designed adaptively such that the joint system and absorber stationary state
is pure at a reference parameter value. The scheme offers an exciting prospect
for optimal continuous-time adaptive measurements, but more work is needed to
find realistic practical implementations.

##  1 Introduction

The quantum input-output (I-O) formalism is an effective framework for
describing the evolution, monitoring and control of Markovian quantum open
systems [21, 58, 12]. In this setting, the interaction with the environment is
modelled by coupling the system of interest with a quantum transmission line
(channel) represented by a Gaussian bosonic field. The output field carries
information about the open system’s dynamics which can be accessed by
performing continuous-time measurements, and the corresponding conditional
system evolution is described in terms of stochastic Schrödinger or filtering
equations [3, 7, 9, 13, 57].

While these theories are key to quantum engineering applications, they rely on
the precise knowledge of the system’s dynamical parameters (e.g. Hamiltonian
of field coupling), which are often uncertain, or completely unknown, and
therefore need to be estimated from measurement data. The I-O formalism is
ideally suited for this statistical inference task, and more generally for
implementing online system identification methods [40]. Unlike direct
measurement techniques which require repeated system re-preparations and fast
control operations [32, 17, 15, 4, 50, 14, 60], the parameters can be
estimated continuously from the output measurement trajectory, even if the
system is not directly accessible or it is involved in an information
processing task. The first investigation in parameter estimation for
continuously-observed quantum systems considered the estimation of the Rabi
frequency of an atom in a cavity mode, while a photon counting measurement is
performed on the cavity output [41]. Subsequent works have addressed a variety
of related problems including the dependence on measurement choice [18],
adaptive estimation [5, 46] filtering with uncertain parameters [47] particle
filters for estimation [11, 49] achieving Heisenberg scaling [42, 1, 2],
sensing with error correction [45], Bayesian estimation [19, 43, 33, 48, 59],
quantum smoothing [52, 55, 56, 27], waveform estimation [54, 53] estimation of
linear systems [26, 22, 39, 38], classical and quantum Fisher informations of
the output channel [28, 19, 10, 34, 20, 22, 29, 25].

An upshot of these studies is that standard measurements such as counting,
homodyne or heterodyne generally do no achieve the ultimate limit given by the
quantum Cramér–Rao bound [30, 31, 8], while the optimal measurement prescribed
by the symmetric logarithmic derivative requires collective operations on the
output state. In this work we make a first step towards addressing the key
issue of devising _realistic and statistically effective_ measurement
strategies within the the framework of the I-O theory. By realistic we mean
procedures which involve sequential continuous-time measurements (as opposed
to general non-separable measurements on the output state), possibly combined
with more advanced but theoretically well understood operations such as series
connections and feedback [24].

For conceptual clarity we focus primarily on discrete-time dynamics, but we
will indicate how the techniques may be extended to continuous-time. In the
discrete-time setting, the I-O dynamics consists of a $d$-dimensional system
of interest interacting sequentially with a chain of $k$-dimensional ‘noise’
input units, which are identically and independently prepared in a state
$|\chi\rangle$, cf. Figure 1.

Figure 1: Quantum input-output discrete-time dynamics with $\theta$ dependent
unitary interaction $U_{\theta}$

We assume that the interaction unitary $U_{\theta}$ acting on
$\mathbb{C}^{d}\otimes\mathbb{C}^{k}$ depends on a parameter
$\theta\in\mathbb{R}$, which we would like to estimate by measuring the output
after $n$ evolution steps. In principle, this can be done by applying the
adaptive, separable measurement scheme developed in [61] to the joint pure
system-output state; indeed this has been shown to attain the quantum Cramér-
Rao bound. However, while theoretically applicable, the algorithm involves
manipulating multi-partite operators, making it unsuitable for processing
output states with a large number $n$ of noise units. In addition, it is not
clear how the algorithm can be applied to continuous-time dynamics.

Our main contribution is to eliminate these drawbacks by devising a scheme
which exploits the intrinsic Markovian structure of the problem. Concretely,
we propose an algorithm which finds optimal measurement bases for each of the
output units by only performing computations on the space of a doubled-up
system and a noise unit, i.e.
$\mathbb{C}^{d}\otimes\mathbb{C}^{d}\otimes\mathbb{C}^{k}$. Our algorithm has
a similar structure to that of the quantum state filter describing the
system’s conditional evolution, and can be run ‘in real-time’ without having
to specify the time length $n$ in advance. This eliminates the need for
storing the quantum output while the computation is performed. While our
general algorithm requires measurements on both the output and system in order
to achieve optimality, in Proposition 6.1 we prove that by measuring only the
output we incur a loss of Fisher information which is bounded by a constant,
independent of the time $n$. Since the quantum Fisher information scales
linearly in time, this implies that output measurement is optimal in the
leading contribution to the QFI.

We now describe our scheme in more detail. In the first stage of the protocol
we use a small proportion of the output units (of samplesize $\tilde{n}\approx
n^{1-\epsilon}$ with small $\epsilon$) in order to compute a preliminary
‘rough estimator’ $\theta_{0}$ of the true parameter $\theta$, by performing a
standard sequential measurement. This step is necessary in _any_ quantum
estimation problem in which the optimal measurement depends on the unknown
parameter. In particular, this means that strictly speaking one can only
attain optimality in the limit of large sample sizes, as $\theta-\theta_{0}$
decays as $n^{-(1-\epsilon)/2}$ [23].

Figure 2:

In the second (main) stage of the protocol we use $\theta_{0}$ to design an
sequential measurement which achieves the output QFI at $\theta=\theta_{0}$,
and in conjunction with the first stage, at any parameter value $\theta$, for
large $n$. This is illustrated in Figure 2: each output unit undergoes a
‘quantum post-processing’ transformation followed by an adaptive projective
measurement whose basis is computed according to the ‘measurement filter’
algorithm described below. More specifically, after interacting with the
system, a unitary $V_{\theta_{0}}$ is applied to the noise unit together with
an additional ancilla of the same size $d$ as the system. The system and
ancilla can be regarded as a single open system (denotes ‘s+a’) of dimension
$D=d^{2}$ coupled to the noise units via the unitary
$W_{\theta}=V_{\theta_{0}}U_{\theta}$. The unitary $V_{\theta_{0}}$ is chosen
such that s+a has a _pure stationary state_ $|\psi\rangle\in\mathbb{C}^{D}$ at
$\theta=\theta_{0}$, and the output state is identical to that of the input.
This is a discrete-time analogue of the notion of _coherent quantum absorber_
introduced in [51], and it insures that the ‘reference’ output state at
$\theta=\theta_{0}$ is the same as the product input state (the ‘vacuum’),
while deviations from $\theta_{0}$ produce ’excitations’ in the output. After
the interaction with the ancilla (absorber), the noise unit is measured in a
basis determined by an iterative algorithm which we briefly outline here ( see
section 5 for more details).

Initialisation step ($j=1$). The first measurement basis
$\left\\{|e^{[1]}_{i}\rangle\right\\}_{i=1}^{k}$ in $\mathbb{C}^{k}$ is
determined by the conditions

$\left\langle
e^{[1]}_{i}\right|B_{1}\left|e^{[1]}_{i}\right\rangle=0,\qquad\left|\left\langle\left.e^{[1]}_{i}\right|\chi\right\rangle\right|^{2}=\frac{1}{k},\qquad{\rm
for~{}all~{}}i=1,\dots k.$

where $B_{1}={\rm Tr}_{s+a}(A_{1})\in M(\mathbb{C}^{k})$ and
$A_{1}=\dot{W}P-P\dot{W}^{*}\in M(\mathbb{C}^{D}\otimes\mathbb{C}^{k})$ with
$\dot{W}=\left.\frac{dW_{\theta}}{d\theta}\right|_{\theta_{0}}$ and
$P=|\psi\otimes\chi\rangle\langle\psi\otimes\chi|.$ By measuring the first
output unit in the basis $\left\\{|e^{[1]}_{i}\rangle\right\\}_{i=1}^{k}$ we
obtain the outcome $X_{1}=i_{1}$, and then compute the ‘measurement filter’
operator $\Pi_{1}=\left\langle
e^{[1]}_{i_{1}}\right|A_{1}\left|e^{[1]}_{i_{1}}\right\rangle$ on
$\mathbb{C}^{D}$.

Iterative step (with $j=2,\dots n$). Given $\Pi_{j-1}$ from the previous step,
we compute the operators

$A_{j}=W_{\theta_{0}}\left(\Pi_{j-1}\otimes
P\right)W_{\theta_{0}}^{*}+\frac{1}{D^{j-1}}A_{1},\qquad B_{j}={\rm
Tr}_{s+a}A_{j}$

and find the $j$-th measurement basis
$\left\\{|e^{[j]}_{i}\rangle\right\\}_{i=1}^{k}$ by applying the conditions

$\left\langle
e^{[j]}_{i}\right|B_{j}\left|e^{[j]}_{i}\right\rangle=0,\qquad{\rm
and}\qquad\left|\braket{e^{[j]}_{i}}{\chi}\right|^{2}=\frac{1}{k},\qquad{\rm
for~{}all~{}}i=1,\dots k.$

By measuring the $j$-th output unit in the basis
$\left\\{|e^{[j]}_{i}\rangle\right\\}_{i=1}^{k}$ we obtain the outcome
$X_{j}=i_{j}$, and compute $\Pi_{j}=\left\langle
e^{[j]}_{i_{j}}\right|A_{j}\left|e^{[j]}_{i_{j}}\right\rangle$.

Final measurement on s+a (optional). In order to achieve the full system-
output QFI, a final measurement can be applied to s+a, whose basis
$\\{|e^{[s+a]}_{i}\rangle\\}_{i=1}^{D}$ is determined by the conditions

$\bra{e^{[s+a]}_{i}}\Pi_{n}\ket{e^{[s+a]}_{i}}=0,\qquad\left|\left\langle
e^{[s+a]}_{i}|\psi\right\rangle\right|^{2}=\frac{1}{d}.$

As stated above, the final measurement improves the estimation error for small
$n$, but is only marginally beneficial for large $n$.

In section 7 we describe the results of a simulation study for a simplified
model where the system plus absorber are represented by a two-dimensional
system with a pure stationary state. While this sidesteps the preliminary
estimation stage of the protocol, it allows us to specifically test the key
features of the adaptive measurement algorithm with a reasonably large
trajectory length and a high number of repetitions. For this model, we can
explicitly compute the system-output QFI (cf. Lemma 7.1), while the classical
Fisher information of any output measurement strategy can be estimated by
sampling techniques. The results confirm that the adaptive measurement attains
the QFI when the system is measured at the end, while the output-only strategy
is only worse by a constant independent of trajectory length. On the other
hand, simple measurements (same fixed basis for each unit) perform strictly
worse even when the measurement basis is optimised. While the improvement here
is not dramatic, our preliminary investigations indicate that the gap
increases significantly with the system dimension, depending on the chosen
model. We further test the performance of the maximum likelihood estimator and
find that its mean square error approaches the inverse of the classical Fisher
information in the long-time limit, which agrees broadly with the Cramér-Rao
bound.

Finally, we note that our scheme can be extended to continuous-time dynamics
by using standard time-discretisation techniques [44, 37]. Although we do not
treat this in detail here, we comment on this extension at the end of the
paper.

The paper is organised as follows. In section 2 we briefly review the adaptive
algorithm for optimal separable measurements developed in [61]. Section 3
introduces the Markov dynamics setting and reviews a key result on the
asymptotic QFI of the output. Section 4 explains how the use of ‘post-
processing’ by quantum absorber reduces the general estimation problem to one
concerning a system with a pure stationary state. This is then used in section
5, which details the adaptive measurement procedure including the key
‘measurement filter’ algorithm. In section 6.2 we show that the proposed
adaptive output measurement achieves the optimal QFI rate even if the system
is not measured. We also devise a scheme to estimate the classical Fisher
information of the measurement process by sampling over trajectories. Section
7 presents simulation results using an elaboration of an amplitude decay qubit
model.

##  2 Optimal separable measurements

In this section we review a result by Zhou, Zou and Jiang [61] concerning
optimal parameter estimation for multipartite pure states, using separable
measurements (local measurements and classical communication). Their method
will then be applied to the problem of estimating parameters of discrete time
quantum input-output systems. By exploiting the Markovian nature of the
dynamics, we will show that the algorithm can be recast in a simpler procedure
akin to that of a quantum state filter.

Consider a one parameter quantum statistical model
$\\{\rho_{\theta}:\theta\in\mathbb{R}\\}$ where $\rho_{\theta}$ is a state on
a Hilbert space $\mathcal{H}$ which depends smoothly on the unknown parameter
$\theta$. To estimate $\theta$ we perform a measurement on the state
$\rho_{\theta}$ and compute an estimator $\hat{\theta}$ based on the
measurement outcome. According to the quantum Cramér-Rao bound (QCRB) [30, 31,
8], the variance of any unbiased estimator $\hat{\theta}$ is lower bounded as

${\rm Var}(\hat{\theta})=\mathbb{E}\left[(\hat{\theta}-\theta)^{2}\right]\geq
F_{\theta}^{-1}$

where $F_{\theta}$ is the quantum Fisher information (QFI) defined as
$F(\theta)={\rm Tr}(\rho_{\theta}\mathcal{L}_{\theta}^{2})$, and
$\mathcal{L}_{\theta}$ is the symmetric logarithmic derivative (SLD)
satisfying
$\frac{d}{d\theta}\rho_{\theta}=\mathcal{L}_{\theta}\circ\rho_{\theta}$. In
general, for any given parameter value $\theta_{0}$, the QCRB is saturated111
This achievability argument can be made rigorous in an asymptotic setting
where the experimenter has an ensemble of $n$ independent, identically
prepared systems and employs an adaptive procedure for ‘locating’ the
parameter [23, 16].  by measuring the SLD $\mathcal{L}_{\theta_{0}}$ and
constructing a locally unbiased estimator
$\hat{\theta}=\theta_{0}+X/F(\theta_{0})$ where $X$ is the measurement
outcome.

While for full rank states the optimal measurement is essentially unique, for
rank deficient states this is not the case and a necessary and sufficient
condition for a measurement to saturate the QCRB has been derived in [8]. This
has practical relevance for multipartite systems where the measurement of the
SLD may not be easy to implement. Motivated by this limitation, the
saturability condition has been further investigated in [61] where it is shown
that the QCRB for pure states of multipartite systems is achievable using
separate measurements constructed in an adaptive fashion which we now proceed
to describe.

Consider the pure state model
$\rho_{\theta}=\ket{\psi_{\theta}}\bra{\psi_{\theta}}$ with
$|\psi_{\theta}\rangle\in\mathcal{H}$. We denote
$|\dot{\psi}_{\theta}\rangle=\frac{d}{d\theta}|\psi_{\theta}\rangle$ and
assume that $\langle\psi_{\theta}|\dot{\psi}_{\theta}\rangle=0.$ This can
generally be arranged by choosing the (unphysical) phase of
$|\psi_{\theta}\rangle$ to have an appropriate dependence on $\theta$. In
particular, in this case we have
$\ket{\psi_{\theta}^{\perp}}:=(1-\ket{\psi_{\theta}}\bra{\psi_{\theta}})\ket{\dot{\psi_{\theta}}}=\ket{\dot{\psi_{\theta}}}$.
Under this assumption the QFI is given by

$F_{\theta}=4\|\dot{\psi}_{\theta}\|^{2}.$ (1)

Further, we define the operator $M$ which will play a key role in the analysis

$M=\ket{\psi_{\theta}}\bra{\dot{\psi_{\theta}}}-\ket{\dot{\psi_{\theta}}}\bra{\psi_{\theta}}.$
(2)

The authors of [61] note that if a projective rank-one measurement
$\\{E_{i}=|e_{i}\rangle\langle e_{i}|\\}$ satisfies the conditions

$\langle e_{i}|M|e_{i}\rangle=0,\qquad{\rm and}\qquad p_{\theta}(i)=|\langle
e_{i}|\psi_{\theta}\rangle|^{2}=1/D,\qquad D={\rm dim}(\mathcal{H})$ (3)

then it fulfils the general criteria of [8] and therefore it achieves the
QCRB. Note that the second condition can be relaxed to $p(i)\neq 0$ for all
$i$, but we will stick to the chosen expression for concreteness.

We now assume that we deal with a multipartite system such that
$\mathcal{H}=\mathcal{H}_{1}\otimes\mathcal{H}_{2}\otimes\dots\otimes\mathcal{H}_{n}$,
with ${\rm dim}(\mathcal{H}_{i})=d_{i}$, and follow [61] to show that a
separable measurement satisfying the above conditions can be constructed by
using the algorithm outlined below. For any set of indices
$A\subset\\{1,\dots,n\\}$ we denote its complement by $A^{c}$ and by
$\rho_{A}={\rm Tr}_{A^{c}}(\rho_{\theta})$ the partial state of the sub-
systems with indices in $A$. For $k>1$ we denote by $\underline{k}$ the set
$\\{1,\dots,k\\}$. Similarly, we denote $M_{A}={\rm Tr}_{A^{c}}(M)$ and for
single sub-systems ($A=\\{i\\}$) we use the notation $\rho_{i}$ and $M_{i}$.

We measure the sub-systems sequentially, such that each individual measurement
basis depends on the outcomes of the previous measurements, as follows. In the
first step, the measurement basis $\\{|e_{i}^{[1]}\rangle\\}_{i=1}^{d_{1}}$ of
system $\mathcal{H}_{1}$ is chosen such that (see [61] for how to construct
such a basis)

$\left\langle
e_{i}^{[1]}\right|M_{1}\left|e_{i}^{[1]}\right\rangle=0\qquad{\rm and}\qquad
p_{1}(i)=\left\langle
e_{i}^{[1]}\right|\rho_{1}\left|e_{i}^{[1]}\right\rangle=\frac{1}{d_{1}}.$

After this, the following procedure is applied sequentially to determine the
measurement basis for system $k+1$ with $k=1,\dots,n-1$: given the outcomes
$i_{\underline{k}}:=\\{i_{1},\dots,i_{k}\\}$ of the first $k$ measurements, we
choose the basis
$\\{|e^{[k+1]}_{i}(i_{\underline{k}})\rangle\\}_{i=1}^{d_{k+1}}$ in
$\mathcal{H}_{k+1}$ such that

$\left\langle
e_{i}^{[k+1]}(i_{\underline{k}})\right|M_{k+1}(i_{\underline{k}})\left|e_{i}^{[k+1]}(i_{\underline{k}})\right\rangle=0\quad{\rm
and}\quad p_{k+1}(i_{\underline{k}},i)=\left\langle
e_{i}^{[k+1]}(i_{\underline{k}})\right|\rho_{k+1}(i_{\underline{k}})\left|e_{i}^{[k+1]}(i_{\underline{k}})\right\rangle=\frac{1}{d_{1}\cdot
d_{2}\dots\cdot d_{k}}$

for all $i=1,\dots,d_{k+1}$, where

$M_{k+1}(i_{\underline{k}})=\left\langle
e^{[\underline{k}]}_{i_{\underline{k}}}\right|M_{\underline{k+1}}\left|e^{[\underline{k}]}_{i_{\underline{k}}}\right\rangle,\quad\rho_{k+1}(i_{\underline{k}})=\left\langle
e^{[\underline{k}]}_{i_{\underline{k}}}\right|\rho_{\underline{k+1}}\left|e^{[\underline{k}]}_{i_{\underline{k}}}\right\rangle,\quad{\rm
and}\quad\left|e^{[\underline{k}]}_{i_{\underline{k}}}\right\rangle=\left|e^{[1]}_{i_{1}}\right\rangle\otimes\dots\otimes\left|e^{[k]}_{i_{k}}\right\rangle$

Note that the second condition means that for each $k$ the outcome $i_{k}$ is
independent of the others and has equal probabilities $1/d_{k}$.

After $n$ steps we have defined in an adpative fashion a product measurement
basis

$\left|e^{[\underline{n}]}_{i_{\underline{n}}}\right\rangle=\left|e^{[1]}_{i_{1}}\right\rangle\otimes\dots\otimes\left|e^{[n]}_{i_{n}}\right\rangle$

and one can verify that such a measurement satisfies the general condition
(3).

##  3 Discrete quantum Markov chains and the output QFI

In the input-output formalism the dynamics of a discrete-time quantum open
system $\mathcal{H}_{s}\cong\mathbb{C}^{d}$ is modeled by successive unitary
interactions with independent ‘noise units’, identically prepared in a state
$|\chi\rangle\in\mathcal{H}_{u}\cong\mathbb{C}^{k}$. This can be pictured as a
conveyor belt where the incoming ‘noise units’ constitute the input, while the
outgoing ‘noise units’ make up the output of the process, cf. Figure 1. If
$|\phi\rangle\in\mathcal{H}_{s}$ is the initial state of the system, and $U$
is the unitary on $\mathcal{H}_{s}\otimes\mathcal{H}_{u}$ describing the
interaction between system and a noise unit, then the state of the system and
output after $n$ time units is

$|\Psi(n)\rangle=U(n)|\phi\otimes\chi^{\otimes
n}\rangle=U^{(n)}\cdot\dots\cdot U^{(2)}\cdot U^{(1)}|\phi\otimes\chi^{\otimes
n}\rangle\in\mathcal{H}_{s}\otimes\mathcal{H}_{u}^{\otimes n}$ (4)

where $U^{(i)}$ is the unitary acting on the system and the $i$-th noise unit.
From equation (4) we find that the reduced state of the system at time $n$ is
given by

$\rho(n):={\rm Tr}_{\rm out}(|\Psi(n)\rangle\langle\Psi(n)|)=T^{n}(\rho_{\rm
in}),\qquad\rho_{\rm in}=|\phi\rangle\langle\phi|,$

where the partial trace is taken over the output (noise units), and
$T:\mathcal{T}_{1}(\mathcal{H}_{s})\to\mathcal{T}_{1}(\mathcal{H}_{s})$ is the
Markov transition operator

$T:\rho\mapsto{\rm
Tr}_{u}(U(\rho\otimes\tau)U^{*}),\qquad\tau:=|\chi\rangle\langle\chi|.$

Fixing an orthonormal basis $\\{|1\rangle,\dots,|k\rangle\\}$ in
$\mathcal{H}_{u}$, we can express the system-output state as a matrix product
state

$|\Psi(n)\rangle=\sum_{i_{1},\dots,i_{n}=1}^{k}K_{i_{n}}\dots
K_{i_{1}}|\phi\rangle\otimes|i_{1}\rangle\otimes\dots\otimes|i_{n}\rangle$ (5)

where $K_{i}=\langle i|U|\chi\rangle$ are the Kraus operators of $T$, so that
$T(\rho)=\sum_{i}K_{i}\rho K_{i}^{*}$.

Now, let us assume that the dynamics depends on a parameter
$\theta\in\mathbb{R}$ which we would like to estimate, so that $U=U_{\theta}$
and $|\Psi(n)\rangle=|\Psi_{\theta}(n)\rangle$. In the input-output formalism
it is usually assumed that the experimenter can measure the output (noise
units after the interaction) but may not have access to the system. In this
case the relevant quantum statistical model is that of the (mixed) output
state given by

$\rho^{\rm out}_{\theta}(n)={\rm
Tr}_{s}(|\Psi_{\theta}(n)\rangle\langle\Psi_{\theta}(n)|).$

The problem of estimating $\theta$ in this formulation has been investigate in
both the discrete time [28, 29] and the continuous time [10, 20, 25] settings.
For our purposes, we summarise here the relevant results of [29]. We will
assume that the Markov chain is _primitive_ , i.e. the transition operator $T$
has a unique full-rank steady state $\rho_{\rm ss}$ (i.e. $T(\rho_{\rm
ss})=\rho_{\rm ss}$) , and is aperiodic (i.e. the only eigenvalue of $T$ on
the unit circle is $1$). In particular, for any initial state $\rho_{\rm in}$,
the system converges to the stationary state $T^{n}(\rho_{\rm in})\to\rho_{\rm
ss}$ in the large $n$ limit. Therefore, for the asymptotic analysis we can
assume that the dynamics is in the _stationary regime_ and focus on the large
time properties of the stationary output state. The following Theorem shows
that the output QFI scales linearly with time and provides an explicit
expression of the rate.

######  Theorem 3.1.

Consider a primitive discrete time Markov chain as described above, whose
unitary depends smoothly on a one-dimensional parameter $\theta$, so that
$U=U_{\theta}$. The quantum Fisher information of the output $F_{\theta}(n)$
state $\rho^{\rm out}_{\theta}(n)$ scales linearly with $n$ and its rate is
equal to

$\lim_{n\to\infty}\frac{1}{n}F_{\theta}(n)=f_{\theta}=4\sum_{i=1}^{k}\left[{\rm
Tr}\left[\rho_{ss}\dot{K}_{i}^{*}\dot{K}_{i}\right]+2{\rm Tr}\left[{\rm
Im}(K_{i}\rho_{ss}\dot{K}_{i}^{*})\cdot\mathcal{R}({\rm
Im}\sum_{j}\dot{K}_{j}^{*}K_{j})\right]\right]$ (6)

where $\mathcal{R}$ is the Moore-Penrose inverse of ${\rm Id}-T_{\theta}$.

Following standard quantum Cramér-Rao theory [30, 31, 8], the theorem implies
that the variance of any (unbiased) output-based estimator is bounded from
below by $n^{-1}/f(\theta)$ for large $n$. A more in depth analysis [29] shows
that the output model satisfies the property of _local asymptotic normality_
which pertains to a certain quantum Gaussian approximation of the output state
and implies that there exists an estimator $\hat{\theta}_{n}$ which achieves
the CR bound asymptotically and has normally distributed errors:

$\sqrt{n}\left(\hat{\theta}_{n}-\theta\right)\longrightarrow
N(0,f(\theta)^{-1})$

where the convergence is in distribution to a normal variable with variance
$f(\theta)^{-1}$. Below, we will make use of an extension of Theorem 3.1 which
shows that the same result holds for rank-deficient stationary states of
ergodic chains, and in particular for pure states [35].

Having identified the output QFI rate, we would like to investigate
measurement schemes which can provide good accuracy for estimating the
parameter $\theta$. As noted before, the QCRB can be achieved by measuring the
SLD of the statistical model. However, the SLD of the output state is
generally a complicated operator whose measurement requires collective
operations on the noise units. On the other hand, one can consider separate
measurements of the same observable on the different noise units and system.
The average statistic may provide an efficient estimator and its (asymptotic)
Fisher information can be computed explicitly [10]. However, such measurements
are in general not optimal. Here would would like to ask the more fundamental
question: is it possible to achieve the QCRB using simpler ‘local’
manipulation of the output units which involve operations on single, rather
than multiple units.

##  4 Output post-processing using quantum coherent absorber

In this section we introduce a key tool used in our estimation scheme, which
allows us to recast the estimation problem for a general primitive Markov
chain with unitary $U$, into one concerning a Markov chain with a ‘doubled-up’
system having a pure stationary state.

Consider the input-output system of section 3 characterised by a unitary $U$
on $\mathcal{H}_{s}\otimes\mathcal{H}_{u}$. We now modify the setup as
illustrated in Figure 2 by inserting an additional $d$-dimensional system
$\mathcal{H}_{a}$ called a _coherent quantum absorber_ which interacts with
each of the noise units via a fixed unitary $V$ on
$\mathcal{H}_{a}\otimes\mathcal{H}_{u}$ applied immediately after $U$. This
can be seen as a _quantum post-processing_ of the output prior to the
measurement. The original system and the absorber can be considered a single
open system with space $\mathcal{H}_{s}\otimes\mathcal{H}_{a}$ which interacts
with the same conveyor belt of noise units via the unitary $W_{\theta}=V\cdot
U$ where $V,U$ are now understood as the ampliations of the unitaries to the
tensor product $\mathcal{H}_{s}\otimes\mathcal{H}_{a}\otimes\mathcal{H}_{u}$.
The following lemma shows that for certain choices of $V$ the auxiliary system
forms a pure stationary state together with the original one, and the noise
units pass unperturbed from input to output. This explains the ‘absorber’
terminology, which was originally introduced in the context of continuous-time
input-output dynamics [51].

######  Lemma 4.1.

Any given primitive Markov with unitary $U$ can be extended to a quantum
Markov chain including an absorber with unitary $V$, such that the doubled-up
system has a pure stationary state
$|\tilde{\psi}\rangle\in\mathcal{H}_{s}\otimes\mathcal{H}_{a}$ and

$W:|\tilde{\psi}\rangle\otimes|\chi\rangle\mapsto|\tilde{\psi}\rangle\otimes|\chi\rangle,\qquad
W=VU.$

In particular, if the initial state of the doubled-up system is
$|\tilde{\psi}\rangle$, then the $n$-steps output state of the doubled-up
system is identical to the input state $|\chi\rangle^{\otimes n}$.

_Proof._ Let $\rho_{\rm ss}=\sum_{i}\lambda_{i}|f_{i}\rangle\langle f_{i}|$ be
the spectral decomposition of the stationary state of the original system with
unitary $U$. We construct the purification

$|\tilde{\psi}\rangle=\sum\sqrt{\lambda_{i}}|f_{i}\rangle\otimes|f_{i}\rangle\in\mathcal{H}_{s}\otimes\mathcal{H}_{a}$

which will play the role of stationary state of the extended system. Let
$|\phi\rangle:=U|\tilde{\psi}\otimes\chi\rangle\in\mathcal{H}_{s}\otimes\mathcal{H}_{a}\otimes\mathcal{H}_{u}$
be the state after applying $U$. We therefore look for unitary $V$ on
$\mathcal{H}_{a}\otimes\mathcal{H}_{u}$ (ampliated by identity on
$\mathcal{H}_{s}$) such that $V$ reverts the action of $U$

$V:|\phi\rangle\mapsto|\tilde{\psi}\otimes\chi\rangle.$

Since ${\rm Tr}_{a}(|\tilde{\psi}\rangle\langle\tilde{\psi}|)=\rho_{\rm ss}$
this means that the reduced state of the system after applying $U$ is still
the stationary state, so that

$|\phi\rangle=\sum_{i}\sqrt{\lambda_{i}}|f_{i}\rangle\otimes|g_{i}\rangle$

where $|g_{i}\rangle$ are mutually orthogonal unit vectors in
$\mathcal{H}_{a}\otimes\mathcal{H}_{u}$. We now choose a unitary $V$ such that
$V|g_{i}\rangle=|f_{i}\otimes\chi\rangle$, for all $i$, which is always
possible due to orthogonality.

∎

##  5 Adaptive measurement algorithm

In this section we describe our adaptive output measurement protocol for
estimating an unknown one-dimensional dynamical parameter $\theta$ of a
discrete time quantum Markov chain with unitary $U_{\theta}$, as described in
section 3. The protocol has two stages. In the first stage we use a small
proportion (e.g. $\tilde{n}=n^{1-\epsilon}$, with $0<\epsilon\ll 1$) of the
output units in order to compute a preliminary ‘rough estimator’ $\theta_{0}$
of the true parameter $\theta$ by performing a standard sequential
measurement. This step is necessary in any quantum estimation problem in which
the optimal measurement depends on the unknown parameter [23, 16], and will
inform the second stage of the protocol. In the second stage we use
$\theta_{0}$ to design a sequential measurement which achieves the output QFI
at $\theta=\theta_{0}$. Since
$\delta\theta=\theta-\theta_{0}=O(n^{-1/2+\epsilon})$ [10], this implies that
the procedure is optimal for any parameter value $\theta$ for large n. This
stage has two key ingredients (cf. Figure 2): a ‘quantum post-processing’
operation on output units immediately after interacting with the system,
followed by an adaptive projective measurement whose basis is computed
according to a ‘measurement filter’ algorithm inspired by [61]. We now
describe the two steps in detail.

_Quantum post-processing._ After the interaction $U_{\theta}$ with the system,
the output units interact sequentially with a d-dimensional ancillary system,
cf. Figure 2. The interaction unitary $V_{\theta_{0}}$ is chosen such that the
ancillary system becomes a coherent quantum absorber for $\theta=\theta_{0}$,
cf. section 4. This means that the system plus absorber (s+a) can be regarded
as a single open system with a unitary $W_{\theta}=V_{\theta_{0}}U_{\theta}$
and pure stationary state $|\psi\rangle\in\mathbb{C}^{D}$, with $D=d^{2}$,
whose output state is identical to the input. The general estimation problem
for $U_{\theta}$ has been reduced to a special one for a doubled-up system
with unitary $W_{\theta}$ which features a pure stationary state at
$\theta=\theta_{0}$.

_Adaptive measurement algorithm._ We will assume for simplicity that the
inital state of $s+a$ is the stationary state $|\psi\rangle$ such that the
full (s+a)-output state is
$|\Psi_{\theta}(n)\rangle=W_{\theta}(n)|\psi\otimes\chi^{\otimes n}\rangle$ as
defined in (4). In principle we could now apply the algorithm described in
section 2 to our model to achieve the system-output QFI (noting that
postprocessing does not change the system-output QFI). However, this procedure
has some drawbacks. Indeed, in order to compute the measurement bases one
needs to work with large dimensional spaces which becomes unfeasible in an
asymptotic setting. Secondly, it is not clear a priori whether the output
units can be measured immediately after the interaction with the system, and
whether the measurements depend on the length of the output (sample size). In
addition, the procedure requires a final measurement on the system, which may
be impractical in the context of input-output dynamics. We will show that all
these issues can be addressed by taking into account the Markovian structure
of our model, and exploiting the pure stationary state property. Let us denote
$W=W_{\theta_{0}}$, $\dot{W}=\left.\frac{dW}{d\theta}\right|_{\theta_{0}}$ and

$A_{1}=M^{(1)}=\dot{W}|\psi\otimes\chi\rangle\langle\psi\otimes\chi|W^{*}-W|\psi\otimes\chi\rangle\langle\psi\otimes\chi|\dot{W}^{*}=\dot{W}P_{\psi\otimes\chi}-P_{\psi\otimes\chi}\dot{W}^{*},$

and

$B_{1}={\rm Tr}_{s+a}A_{1}=KP_{\chi}-P_{\chi}K^{*},\qquad{\rm with}\qquad
K=\langle\psi|\dot{W}|\psi\rangle.$

In Appendix 9 we show that the adaptive measurement of [61] reduces to the
following iterative algorithm which is conceptually similar to quantum state
filtering [3, 7], and involves individual measurements on the output units
immediately after the interaction with the system, and computations with
operators on $\mathbb{C}^{D}\otimes\mathbb{C}^{k}$ at each step.

Initialisation step (j=1). The first measurement basis
$\left\\{\ket{e^{[1]}_{i}}\right\\}$ in $\mathbb{C}^{k}$ is chosen such that
the following conditions are fulfilled:

$\left\langle
e^{[1]}_{i}\right|B_{1}\left|e^{[1]}_{i}\right\rangle=0,\qquad{\rm
and}\qquad\left|\braket{e^{[1]}_{i}}{\chi}\right|^{2}=\frac{1}{k},\qquad{\rm
for~{}all~{}}i=1\dots,k.$

The first noise unit is measured in this basis and the outcome $X_{1}=i_{1}$
is obtained. The filter at time $j=1$ is defined as the (trace zero) s+a
operator

$\Pi_{1}=\left\langle
e^{[1]}_{i_{1}}\right|A_{1}\left|e^{[1]}_{i_{1}}\right\rangle.$

Iterative step. The following step is iterated for $j=2,\dots,n$. Given the
filter operator $\Pi_{j-1}$ of the previous step, we define

$A_{j}=\frac{1}{D^{j-1}}A_{1}+W\left(\Pi_{j-1}\otimes
P_{\chi}\right)W^{*},\qquad B_{j}={\rm Tr}_{s+a}A_{j}$

The $j$-th measurement basis $\left\\{\ket{e^{[j]}_{i}}\right\\}$ is chosen to
fulfill the conditions

$\left\langle
e^{[j]}_{i}\right|B_{j}\left|e^{[j]}_{i}\right\rangle=0,\qquad{\rm
and}\qquad\left|\braket{e^{[j]}_{i}}{\chi}\right|^{2}=\frac{1}{k},\qquad{\rm
for~{}all~{}}i=1,\dots k.$ (7)

We measure the $j$-th noise unit in the basis
$\left\\{|e^{[j]}_{i}\rangle\right\\}$ and obtain the result $X_{j}=i_{j}$.
The filter at time $j$ is updated to

$\Pi_{j}=\left\langle
e^{[j]}_{i_{j}}\right|A_{j}\left|e^{[j]}_{i_{j}}\right\rangle.$

Final s+a measurement. This is an optional step which involves a final joint
measurement on system and absorber. The basis
$\left\\{|e^{[s+a]}_{i}\rangle\right\\}_{i=1}^{D}$ is determined by the
following conditions

$\left\langle
e^{[s+a]}_{i}\right|\Pi_{n}\left|e^{[s+a]}_{i}\right\rangle=0,\qquad
p(i|i_{1},\dots,i_{n})=\left|\langle
e^{[s+a]}_{i}|\psi\rangle\right|^{2}=1/D.$

The system and absorber is measured in this basis and the outcome $X=i_{0}$ is
obtained.

The output measurement record $\\{i_{1},\dots,i_{n},i_{0}\\}$ is collected and
used for estimating the parameter $\theta$. The likelihood function is given
by

$p_{\theta}(i_{1},\dots i_{n},i_{0})=\left|\left\langle
e^{[s+a]}_{i_{0}}\otimes e^{[1]}_{i_{1}}\otimes\dots\otimes
e^{[n]}_{i_{n}}\left.\right|\Psi_{\theta}(n)\right\rangle\right|^{2}.$ (8)

For later use, we denote by $p_{\theta}(i_{1},\dots i_{n})$ the marginal
distribution of the output measurement record only.

##  6 Fisher informations considerations

In this section we investigate the relationship between the classical Fisher
information (CFI) of the output adaptive measurement process and the system-
output QFI. We prove that both scale with the same rate and the latter may be
larger than the former by at most a constant, independent of time.

We also provide an expression of the CFI of sequential (adaptive or standard)
output measurements, which is amenable to estimation by sampling. This tool
will be used to confirm the optimality of our adaptive algorithm in numerical
simulations.

###  6.1 Achievability of the QFI with adaptive output measurements

The adaptive measurement scheme described in section 5 insures the CFI of the
full measurement (output and s+a) is equal to the QFI of the full pure state
model. However, in certain physical implementations the system may not be
accessible for measurements, so the more interesting scenario is that in which
only the output state is measured. In this case the CFI will generally be
strictly smaller that the QFI, and the question is whether by measuring only
the output we incur a significant loss of information. In proposition 6.1 we
show that this is not the case: the difference between the QFI and the output
CFI is bounded by a constant, so for large times the loss of information is
negligible compared to both QFI and output CFI, which scale linearly with
time. In section 7 we will illustrate the result on a specific model.

######  Proposition 6.1.

Consider the setup described in section 5, and let $F(n)$ be the system-
absober-output QFI, and $I^{\rm(o)}(n)$ be the output CFI for the optimal
adaptive measurement, at $\theta=\theta_{0}$. Then $F(n)-I^{\rm(o)}(n)<c$ for
all $n$ where $c$ is a constant depending only on the model $U_{\theta}$.
Consequently,

$\lim_{n\to\infty}\frac{1}{n}I^{\rm(o)}(n)=\lim_{n\to\infty}\frac{1}{n}F(n)=f>0$

where $f=f_{\theta_{0}}$ is the QFI rate (6).

###  6.2 Computing the classical Fisher information of the output

The _classical Fisher information_ of the output measurement process at
$\theta_{0}$ is

$I^{{\rm(o)}}(n)=\mathbb{E}_{\theta_{0}}\left(\frac{d\log
p_{\theta}}{d\theta}\right)^{2}=\sum_{i_{1},\dots,i_{n}}p_{\theta_{0}}(i_{1},\dots
i_{n})^{-1}\left(\left.\frac{dp_{\theta}(i_{1},\dots,i_{n})}{d\theta}\right|_{\theta_{0}}\right)^{2}$

where the sum runs over indices such that $p_{\theta}(i_{1},\dots,i_{n})>0$.
In general

$I^{{\rm(o)}}(n)\leq F^{{\rm(o)}}(n)\leq F^{({\rm s+o})}(n)$

where the successive upper bounds are output and system-output QFIs
respectively.

In our simulation study we will be interested to study to what extent these
bounds are saturated in the adaptive and non-adaptive scenarios, and in
particular, to verify the prediction of Proposition 6.1. Since the classical
Fisher information is difficult to compute for long trajectories, we will
recast it as an expectation which can be estimated by sampling measurement
trajectories. In Lemma 6.2 below, we will use the fact that at
$\theta=\theta_{0}$ the vector $|\psi\rangle$ is the stationary state, and
therefore $K^{[j]}_{i}|\psi\rangle=c^{[j]}_{i}|\psi\rangle$ for _any_ Kraus
decomposition $K_{i}^{[j]}=\langle e_{j}^{[i]}|W|\chi\rangle$ (for simplicity
we use the same notation for s+a Kraus operators as in section 3). The proof
of Lemma 6.2 can be found in Appendix 11.

######  Lemma 6.2.

Consider the setup described in section 5. The output CFI at
$\theta=\theta_{0}$ is given by

$I^{{\rm(o)}}(n)=\mathbb{E}_{\theta_{0}}(f^{2})=\sum_{i_{1},\dots,i_{n}}p_{\theta_{0}}(i_{1},\dots,i_{n})f^{2}(i_{1},\dots,i_{n})$
(9)

where $f$ is the function

$f(i_{1},\dots,i_{n})=2{\rm
Re}\sum_{j=1}^{n}\frac{\langle\psi|K^{[n]}_{i_{n}}\dots
K^{[j+1]}_{i_{j+1}}\dot{K}^{[j]}_{i_{j}}|\psi\rangle}{c^{[j]}_{i_{j}}\dots
c^{[n]}_{i_{n}}}$

In particular, $I^{{\rm(o)}}(n)$ can be estimated by computing the empirical
average of $f^{2}$ over sampled trajectories.

We now consider the case where a (projective) measurement $\\{P^{({\rm
s+a})}_{i}\\}$ is performed on s+a, after obtaining the output measurement
record $(i_{1},\dots,i_{n})$. We denote the additional outcome by $i_{0}$. The
CFI of the full process is

$I^{({\rm s+a+o})}(n)=\mathbb{E}_{\theta_{0}}\left(\left.\frac{d\log
p_{\theta}}{d\theta}\right|_{\theta_{0}}\right)^{2}$

where

$p_{\theta}(i_{1},\dots,i_{n},i_{0})=\left\|P^{({\rm
s+a})}_{i_{0}}K^{[n]}_{i_{n}}\dots K^{[1]}_{i_{1}}\psi\right\|^{2}$

is the likelihood function of a trajectory augmented by the system measurement
outcome $i_{0}$. The relevant upper bound in this case is

$I^{({\rm s+a+o})}(n)\leq F^{({\rm s+a+o})}(n)=F^{({\rm s+o})}(n).$ (10)

A similar computation to that of Lemma 6.2 gives the system-output classical
Fisher information

$I^{({\rm s+a+o})}(n)=\mathbb{E}_{\theta_{0}}(\tilde{f}^{2})$

where $\tilde{f}$ is the function

$\displaystyle\tilde{f}(i_{1},\dots,i_{n},i_{0})$ $\displaystyle=$
$\displaystyle 2D\,{\rm Re}\sum_{j=1}^{n}\frac{\langle\psi|P^{({\rm
s})}_{i_{0}}K^{[n]}_{i_{n}}\dots
K^{[j+1]}_{i_{j+1}}\dot{K}^{[j]}_{i_{j}}|\psi\rangle}{c^{[j]}_{i_{j}}\dots
c^{[n]}_{i_{n}}}$

For fixed (non-adaptive) measurements, the bound (10) is generally not
saturated except for special models (e.g. if the state coefficients in the
measurement basis are real for all $\theta$). In contrast, the system-
absorber-output classical Fisher information for adaptive measurements is
equal to the QFI

$I^{({\rm s+a+o})}_{\theta_{0}}(n)=F^{({\rm s+o})}_{\theta_{0}}(n)$

thanks to the optimality of the adaptive measurement procedure. This is
confirmed by our simulation study which also investigates the performance of
the fixed measurement scenario.

##  7 Numerical simulations

We now test the key properties of the adaptive measurement scheme developed in
section 5. For reasons outlined below we do this with a simplified Markov
model and leave a comprehensive numerical study for a future work. The first
simplification is that we will bypass stage-one of the scheme (computing a
rough estimator $\theta_{0}$) and simulate data at $\theta=\theta_{0}$. This
allows us to directly study the performance of the algorithm itself (stage
two), rather than that of the combination of the two stages. The second
simplification is that we choose a system which has a pure stationary state at
$\theta_{0}$; this means that no absorber is required, so the system can be
seen as a surogate for the system+absorber in the general scheme. The reason
for this is mainly practical, as it allows us to use a two-dimensional system
while a system-absorber would have dimension at least four.

The model consists of a two-dimensional system coupled via a unitary
$U_{\theta}$ to two dimensional noise units in state $|\chi\rangle=|0\rangle$.
The input state and the unitary are designed such that the stationary state at
$\theta_{0}=0$ is $|\psi\rangle=|0\rangle$. Since the input is prepared in a
fixed state, we only need to define the action of $U_{\theta}$ on the basis
vectors $|0\rangle\otimes|0\rangle$ and $|1\rangle\otimes|0\rangle$. The
following choice has unknown parameter $\theta$ and two known parameters
$\lambda$ and $\phi$

$\displaystyle U_{\theta}:$
$\displaystyle\ket{00}\longrightarrow\cos(\theta)\sqrt{1-\theta^{2}}\ket{00}+i\sin(\theta)\sqrt{1-\theta^{2}}\ket{10}+\theta\ket{11},$
$\displaystyle U_{\theta}:$ $\displaystyle\ket{10}\longrightarrow
i\sin(\theta)\sqrt{1-\lambda}\ket{00}+\cos(\theta)\sqrt{1-\lambda}\ket{10}+\sqrt{\lambda}e^{i\phi}\ket{01}.$
(11)

A non-zero value of the phase parameter $\phi$ insures that the system-output
state does not have real coefficients in the standard basis, in which case the
standard basis measurement would be optimal. Note that at $\theta_{0}$ the
dynamics provides a simple model for ‘photon-decay’ with decay parameter
$\lambda$.

We compare two measurements scenarios. In the first, non-adaptive scenario,
the noise units are measured in a fixed orthonormal basis
$\\{|f_{0}\rangle=(|0\rangle+|1\rangle)/\sqrt{2},|f_{1}\rangle=(|0\rangle-|1\rangle)/\sqrt{2}\\}$
while in the second scenario they are measured adaptively following the
algorithm described in section 5. In both cases, given the measurement record
$\\{i_{1},\dots,i_{n}\\}\in\\{0,1\\}^{n}$, the conditional state of the system
is

$|\psi_{n}(i_{1},\dots i_{n})\rangle=\frac{K^{[n]}_{i_{n}}\dots
K^{[1]}_{i_{1}}|\psi\rangle}{\|K^{[n]}_{i_{n}}\dots K^{[1]}_{i_{1}}\psi\|}$
(12)

where the Kraus operators are $\\{K_{0}=\langle
f_{0}|U|\chi\rangle,K_{1}=\langle f_{1}|U|\chi\rangle\\}$ in the first
scenario, and $K^{[j]}_{i}=\langle e^{[j]}_{i}|U|\chi\rangle$ in the second
one. The likelihood of the output measurement trajectory is

$p_{\theta}(i_{1},\dots i_{n})=\|K^{[n]}_{i_{n}}\dots
K^{[1]}_{i_{1}}\psi\|^{2}.$ (13)

For two dimensional noise units, the computation of the measurement basis at
each step reduces to the following scheme. We express the traceless, anti-
Hermitian matrix $B_{j}$ defined in section 5 as $B_{j}=i\vec{r}^{\
[j]}\cdot\vec{\sigma}$ where $\vec{r}^{\ [j]}=(r^{[j]}_{x},r^{[j]}_{y},0)$ is
its Bloch vector, and similarly, we let $\pm\vec{s}^{\ [j]}$ be the Bloch
vectors of the basis vectors $\\{|e_{1}^{[j]}\rangle,|e_{2}^{[j]}\rangle\\}$
satisfying the conditions (7). Then one finds that $\vec{s}^{\ [j]}$ can be
taken to be $\vec{s}^{\ [j]}=(r^{[j]}_{y},-r^{[j]}_{x},0)$.

Finally, we study the scenario in which the system of interest is measured
after obtaining the output trajectory. In this case the algorithm guarantees
that the classical Fisher information of the full measurement record is equal
to the QFI of the system-output state. This claim is verified numerically by
comparing the estimated classical Fisher information with the explicit
expression of the QFI computed in section 7.1

###  7.1 The system-output quantum Fisher information

To verify the optimality of our adaptive measurement strategy we will compare
the _classical_ Fisher information computed in section 6.2 with the _quantum_
Fisher information of the system-output state. As shown in Theorem 3.1, the
system-output QFI scales linearly with time, i.e.

$F^{({\rm s+o})}_{\theta}(n)=nf_{\theta}+o(n)$

where the QFI _rate_ is given by the equation (6) and the $o(n)$ term depends
on the specific system parameters. Applying this to our model we obtain

$f_{\theta_{0}}=\frac{8}{1-\sqrt{1-\lambda}}.$ (14)

However, it turns out that for this specific model the QFI can be computed
explicitly for _any_ fixed $n$. The proofs of (14) and of the following lemma
can be found in Appendix 12.

######  Lemma 7.1.

The system-output QFI at $\theta=\theta_{0}$ is given by the formula

$\displaystyle F^{({\rm s+o})}_{\theta}(n)$ $\displaystyle=$
$\displaystyle=\frac{8n}{1-a}+4\left[\frac{2(a^{2}-a^{n})}{(1-a)^{2}}-2\frac{b^{2}(a-a^{n})}{(1-a)^{3}}+\frac{2(a^{2}-a^{2n})}{b^{2}}\right]$
(15)

where $a=\sqrt{1-\lambda}$ and $b=\sqrt{\lambda}$. In particular, the leading
term in $n$ is given by (14) while the remaining terms are bounded.

###  7.2 Simulation studies

The first simulation study focuses on the comparison between the different
notions of Fisher information: the system-output QFI, the CFI of the output
trajectory in the non-adapted and adapted scenarios, and the CFI of the
system-output measurement process in the adapted measurement scenario. The QFI
is computed using the formula in Lemma 7.1 while the CFIs are estimated by
sampling using the expression in Lemma 6.2.

The results are illustrated in Figure 3 where the different informations are
plotted as a function of time (trajectory length) $n$, for
$\lambda=0.8,\phi=\pi/4$. The simulation confirms the fact that the adaptive
algorithm achieves the QFI when the system is measured together with the
output, while the CFI of the output trajectory (without measuring the system)
provides a close approximation which only differs by a constant factor. In
contrast, the CFI of the standard measurement has a smaller rate of increase;
additional numerical work shows that the CFI rate can be improved by
optimising the basis of the standard measurement but it does not achieve the
QFI.

Figure 3: Fisher informations as function of output lenght $n$: quantum Fisher
information (blue), classical Fisher information for the adaptive measurement
with/without system measurement (orange/green), classical Fisher information
for a regular (non-adaptive) measurement.

The second simulation study focuses on the ‘measurement trajectory’, i.e. the
sequence of measurement settings produced in the adaptive measurement
scenario. Since all measurement bases consist of vectors in the equatorial
plane on the Bloch sphere, one can parametrise each basis by the polar
coordinate $\varphi$ of the basis vector $(|0\rangle\pm
e^{i\varphi}|1\rangle)/\sqrt{2}$ for which $\varphi$ belongs to a specified
interval of length $\pi$. This is illustrated in the left panel of Figure 4.
This parametrisation has the disadvantage that it does not reproduce the
topology of the space of measurements which is that of a circle, leading to
some jumps in measurement angles appearing to be larger than the actual
‘distance’ between measurement bases. To remedy this, in the right panel of
Figure 4 we plot $2\varphi$ on circles of radius increasing linearly with
time.

Figure 4: Adaptive measurement trajectory with basis angle $\phi$ on left and
$2\phi$ plotted on the circle on the right.

We note that the initial steps of the trajectory show large variations which
then ‘stabilise’ around a certain range of values. This has to do with the
fact that the initial angle can be chosen arbitrarily in this model as
$B_{1}=0$, cf. section 5. Understanding the nature of this stochastic process
remains an interesting topic of future research.

Figure 5: QFI (blue), CFI of adapted measurement (green) CFI of the standard
measurement (purple), inverse MLE error for adapted measurement (orange) and
for standard measurement (red), as function of trajectory lenght at
$\lambda=0.8,\varphi=\pi/4.$

The third simulation study concerns the performance of the maximum likelihood
estimator (MLE) in an adapted and non-adapted output measurement scenarios.
The MLE is defined by

$\hat{\theta}_{n}:=\underset{\tau}{\arg\max}\,\,p_{\tau}(i_{1},\dots i_{n}).$

where the likelihood $p_{\tau}(i_{1},\dots i_{n})$ is computed as in equation
(13). In numerics, we maximise the log-likelihood function which can be
computed as the sum

$\log p_{\tau}(i_{i},\dots
i_{n})=\sum_{j=1}^{n}\log\|K^{[j]}_{i_{j}}\psi(i_{1},\dots i_{j-1})\|^{2}$

where $\psi(i_{1},\dots i_{j-1})$ is the system’s state conditional on the
output trajectory (filter), cf. equation (12). The MLE accuracy is quantified
by the mean square error (MSE)
$E(n)=\mathbb{E}_{\theta}(\hat{\theta}_{n}-\theta)^{2}$, which is estimated
empirically by averaging over a number of simulations runs to obtain
$\hat{E}(n)$. To verify that the MSE scales as $n^{-1}$ we plot the inverse
empirical error $\hat{E}^{-1}(n)$ as a function of time. Figure 5 shows the
inverse error for the adaptive and non-adaptive measurements together with the
QFI (which is equal to the CFI of the adaptive measurement) and the CFI of the
non-adaptive measurement. We note that for small values of $n$ the inverse
error is significantly lower that the corresponding Fisher information, but it
approaches the latter for larger values of $n$. This suggests that the MLE
achieves the Cramér-Rao bound asymptotically, which is not surprising since
the MLE is known to be asymptotically optimal for independent samples as well
as for certain classes of hidden Markov chains [36, 6]. However, proving its
optimality for the adaptive measurement process remains an open problem. In
addition to the mean square error, we looked at the distribution of the MLE.
Figure 6 shows a histogram of the MLE based on $N=10000$ simulations with
$n=200$, and it indicates that the MLE is approximately normally distributed.
Based on this, it is reasonable to conjecture that the MLE is an _efficient
estimator_ [36, 6] (i.e. has asymptotically normal distribution with variance
equal to the inverse of the Fisher information).

Figure 6: Histogram of the MLE distribution in the adaptive measurement
scenario, with $n=200,\theta=0$ and $10000$ samples, at
$\lambda=0.8,\varphi=\pi/4.$

##  8 Conclusions and Outlook

In this paper we developed an efficient iterative algorithm for optimal
estimation of dynamical parameters of a discrete quantum Markov chain, using
adaptive sequential measurements on the output. The algorithm builds on the
general measurement scheme of [61] which achieves the quantum Fisher
information for pure state models of multipartite systems with one dimensional
unknown parameters. However, unlike the scheme [61] which requires
manipulations involving the full multipartite state, the proposed algorithm
only involves computations on $d^{2}\cdot k$-dimensional systems where $d$ and
$k$ are the dimensions of the open system and noise unit respectively.
Therefore, the method can be readily applied to Markov parameter estimation
for large output sizes. The algorithm exploits the Markovian structure of the
dynamics to sequentially compute optimal measurement bases in terms of a
single time-dependent ‘measurement filter’ operator, which is updated in a way
that is reminiscent of a state filter. One of the key ingredients of the
proposed scheme is the use of a coherent quantum absorber [51] which reduces
the estimation problem to one concerning a system with a pure stationary
state. We considered both output and output-system measurements scenarios and
we showed that while the former achieve the full quantum Fisher information,
the latter is short of this by just a fixed constant, and in particular both
have the same scaling with time. Our theoretical results are confirmed by
numerical simulations using a simplified model related to the amplitude decay
channel.

Our discrete-time procedure raises interesting questions about the the
possibility to design realistic optimal sequential measurements in continuous-
time dynamics. In principle the scheme can be applied to continuous time by
using time-discretisation techniques [44, 37]. Although we did not treat this
in detail here, we can readily make several observations in the special case
of a single Bosonic I-O channel. For small enough time intervals $\delta t$,
the field can be approximated by a noise unit $\mathbb{C}^{2}$ with basis
vectors representing the vacuum and a one-photon state respectively. The
proposed measurements consist of projections whose corresponding Bloch vectors
are in the equatorial plane; loosely speaking, this corresponds to an adaptive
homodyne measurement with a time-dependent angle. Our preliminary
investigations indicate that the behavior of the time-dependent angle ranges
from deterministic evolution to a noisy stochastic process which does not
appear to be of diffusive type. This raises the question whether the remaining
freedom in choosing the measurement bases can be used to improve the adaptive
algorithm and produce a more regular measurement processes. Indeed, the second
defining condition for the measurement vectors can be relaxed, allowing for
more general classes of optimal measurements, which may be more suitable for
continuous-time direction. This will be the topic of a future investigation.

From a theoretical viewpoint, it is important to understand the mathematical
properties of the stochastic processes introduced here, the adaptive
measurement and the measurement trajectory. Finally, it is intriguing to
consider to what extent the proposed method can adapted to multi-parameter
estimation and to general (time-dependent) matrix product states as opposed to
stationary output states of Markov processes.

Acknowledgements: MG was supported by the EPSRC grant EP/T022140/1.

## References

  * [1] F. Albarelli, M. A. C. Rossi, M. G. A. Paris, and M. G. Genoni.  Ultimate limits for quantum magnetometry via time-continuous measurements.  New J. Phys., 19:123011, 2017.
  * [2] F. Albarelli, M. A. C. Rossi, D. Tamascelli, and M. G. Genoni.  Restoring heisenberg scaling in noisy quantum metrology by monitoring the environment.  Quantum, 2:110, 2018.
  * [3] V. P. Belavkin.  Nondemolition principle of quantum measurement theory.  Found. Phys., 24:685–714, 1994.
  * [4] F. Benatti, S. Alipour, and A. T. Rezakhani.  Dissipative quantum metrology in manybody systems of identical particles.  New J. Phys., 16:015023, 2014.
  * [5] D. W. Berry and H. M. Wiseman.  Adaptive quantum measurements of a continuously varying phase.  Phys. Rev. A, 65:043803, 2002.
  * [6] P. J. Bickel, Y. Ritov, and T. Rydén.  Asymptotic normality of the maximum likelihood estimator for general hidden markov models.  Ann. Statist., 26:1614–1635, 1998.
  * [7] L. Bouten, R. van Handel, and M. R. James.  An introduction to quantum filtering.  SIAM J. Control Optim., 46:2199–2241, 2007.
  * [8] Samuel L. Braunstein and Carlton M. Caves.  Statistical distance and the geometry of quantum states.  Phys. Rev. Lett., 72(22):3439–3443, May 1994.
  * [9] H. J. Carmichael.  H.J. Carmichael, An Open Systems Approach to Quantum Optics.  (Springer-Verlag, Berlin, 1993.
  * [10] C. Catana, L. Bouten, and M. Guta.  Fisher informations and local asymptotic normality for continuous-time quantum markov processes.  Journal of Physics A: Mathematical and Theoretical, 48:365301, 2015\.
  * [11] B. Chase and J. M. Geremia.  Single-shot parameter estimation via continuous quantum measurement.  Phys. Rev. A, 79:022314, 2009.
  * [12] J. Combes, J. Kerckhoff, and M. Sarovar.  The SLH framework for modeling quantum input-output networks.  Adv. Phys. X, 2(784), 2017.
  * [13] Jean Dalibard, Yvan Castin, and Klaus Mølmer.  Wave-function approach to dissipative processes in quantum optics.  Phys. Rev. Lett., 68:580–583, 1992.
  * [14] R. Demkowicz-Dobrzanski, J. Czajkowski, and P. Sekatski.  Adaptive quantum metrology under general markovian noise.  Phys. Rev. X, 7:041009, 2017.
  * [15] R. Demkowicz-Dobrzański and M. Guta.  The elusive heisenberg limit in quantum-enhanced metrology.  Nature Communications, 3:1063, 2012.
  * [16] Rafał Demkowicz-Dobrzański, Wojciech Górecki, and Mădălin Guţă.  Multi-parameter estimation beyond quantum fisher information.  Journal of Physics A: Mathematical and Theoretical, 53:363001, 2020\.
  * [17] B. Escher, R. de Matos Filho, and L. Davidovich.  General framework for estimating the ultimate precision limit in noisy quantum-enhanced metrology.  Nature Phys., 7:406–411, 2011.
  * [18] J. Gambetta and H. M. Wiseman.  State and dynamical parameter estimation for open quantum systems.  Phys. Rev. A, 64:042105, 2001.
  * [19] S. Gammelmark and K. Mölmer.  Bayesian parameter inference from continuously monitored quantum systems.  Phys. Rev. A, 87:032115, 2013.
  * [20] S. Gammelmark and K. Mølmer.  Fisher information and the quantum cramér-rao sensitivity limit of continuous measurements.  Phys. Rev. Lett., 112:170401, 2014.
  * [21] C. W. Gardiner and P. Zoller.  Quantum Noise: A Handbook of Markovian and Non-Markovian Quantum Stochastic Methods with Applications to Quantum Optics.  Springer-Verlag, Berlin and New York, 3 edition, 2004.
  * [22] Marco G. Genoni.  Cramér-rao bound for time-continuous measurements in linear gaussian quantum systems.  Phys. Rev. A, 95:012116, 2017.
  * [23] Richard D. Gill and Serge Massar.  State estimation for large ensembles.  Phys. Rev. A, 61:042312, Mar 2000.
  * [24] J. Gough and M. R. James.  The series product and its application to quantum feedforward and feedback networks.  IEEE Trans. Automat. Control, 54(11):2530–2544, 2009.
  * [25] M. Guţǎ and J. Kiukas.  Information geometry and local asymptotic normality for multi-parameter estimation of quantum markov dynamics.  J. Math. Phys., 58(052201):052201, 2017.
  * [26] M. Guţǎ and N. Yamamoto.  System identification for passive linear quantum systems.  IEEE Trans. Automat. Contr., 61(4):921–936, 2016.
  * [27] Ivonne Guevara and Howard Wiseman.  Quantum state smoothing.  Phys. Rev. Lett., 115:180407, 2015.
  * [28] M. Guta.  Fisher information and asymptotic normality in system identification for quantum markov chains.  Physical Review A, 83, 2011.
  * [29] M. Guta and J. Kiukas.  Equivalence classes and local asymptotic normality in system identification for quantum markov chains.  Communications in Mathematical Physics, 335(3):1397–1428, 2014\.
  * [30] C. W. Helstrom.  Quantum detection and estimation theory.  Academic press, 1976.
  * [31] A. S. Holevo.  Probabilistic and Statistical Aspects of Quantum Theory.  North Holland, Amsterdam, 1982.
  * [32] 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–3868, 1997.
  * [33] A. H. Kiilerich and K. Mølmer.  Bayesian parameter estimation by continuous homodyne detection.  Phys. Rev. A, 94:032103, 2016.
  * [34] Alexander Holm Kiilerich and Klaus Mølmer.  Estimation of atomic interaction parameters by photon counting.  Phys. Rev. A, 89:052110, 2014.
  * [35] J. Kiukas and M. Guta.  in preparation.  2022\.
  * [36] E. L. Lehmann and G. Casella.  Theory of point estimation.  Springer, 1998.
  * [37] I. Lesanovsky, M. van Horssen, M. Guta, and J. P. Garrahan.  Characterization of dynamical phase transitions in quantum jump trajectories beyond the properties of the stationary state.  Phys. Rev. Lett., 110:150401, 2013.
  * [38] M. Levitt, M. Guţǎ, and H. I. Nurdin.  Power spectrum identification for quantum linear systems.  Automatica, 90:255–262, 2018.
  * [39] Matthew Levitt and M. Guta.  Identification of single-input–single-output quantum linear systems.  Phys. Rev. A, 95:033825, 2017.
  * [40] L. Ljung.  System Identification: Theory for the User.  Prentice-Hall, 2 edition, 1999.
  * [41] H. Mabuchi.  Dynamical identification of open quantum systems.  Quantum Semiclass. Opt., 8:1103, 1996.
  * [42] Katarzyna Macieszczak, Mădălin Guţă, Igor Lesanovsky, and Juan P. Garrahan.  Dynamical phase transitions as a resource for quantum enhanced metrology.  Physical Review A, 93(2), Feb 2016.
  * [43] A. Negretti and K. Molmer.  Estimation of classical parameters via continuous probing of complementary quantum observables.  New J. Phys., 15:125002, 2013.
  * [44] Tobias J. Osborne, Jens Eisert, and Frank Verstraete.  Holographic quantum states.  Phys. Rev. Lett., 105:260401, 2010.
  * [45] Martin B. Plenio and Susana F. Huelga.  Sensing in the presence of an observed environment.  Phys. Rev. A, 93:032123, 2016.
  * [46] D. T. Pope, H. M. Wiseman, and N. K. Langford.  Adaptive phase estimation is more accurate than nonadaptive phase estimation for continuous beams of light.  Phys. Rev. A, 70:043812, 2004.
  * [47] Jason F. Ralph, Kurt Jacobs, and Charles D. Hill.  Frequency tracking and parameter estimation for robust quantum state estimation.  Phys. Rev. A, 84:052119, 2011.
  * [48] Jason F. Ralph, Simon Maskell, and Kurt Jacobs.  Multiparameter estimation along quantum trajectories with sequential monte carlo methods.  Phys. Rev. A, 96:052306, 2017.
  * [49] P. Six, Ph. Campagne-Ibarcq, L. Bretheau, B. Huard, and P. Rouchon.  Parameter estimation from measurements along quantum trajectories.  In 2015 54th IEEE Conference on Decision and Control (CDC), pages 7742–7748, 2015.
  * [50] Andrea Smirne, J. Kolodynski, S. F. Huelga, and R. Demkowicz-Dobrzanski.  Ultimate precision limits for noisy frequency estimation.  Phys. Rev. Lett., 116:120801, 2016.
  * [51] K. Stannigel, P. Rabl, and P. Zoller.  Driven-dissipative preparation of entangled states in cascaded quantum-optical networks.  New J. Phys., 12:063014, 2012.
  * [52] M. Tsang.  Time-symmetric quantum theory of smoothing.  Phys. Rev. Lett., 102:250403, 2009.
  * [53] M. Tsang and R. Nair.  Fundamental quantum limits to waveform detection.  Phys. Rev. A, 86:042115, 2012.
  * [54] M. Tsang, H. M. Wiseman, and C. Caves.  Fundamental quantum limit to waveform estimation.  Phys. Rev. Lett., 106:090401, 2011.
  * [55] Mankei Tsang.  Optimal waveform estimation for classical and quantum systems via time-symmetric smoothing.  Phys. Rev. A, 80:033840, 2009.
  * [56] Mankei Tsang.  Optimal waveform estimation for classical and quantum systems via time-symmetric smoothing. ii. applications to atomic magnetometry and hardy’s paradox.  Phys. Rev. A, 81:013824, 2010.
  * [57] H. M. Wiseman and G. J. Milburn.  Quantum theory of field-quadrature measurements.  Phys. Rev. A, 47:642–662, 1993.
  * [58] H. M. Wiseman and G. J. Milburn.  Quantum Measurement and Control.  Cambridge University Press, 2010.
  * [59] Cheng Zhang, Kai Zhou, Wei Feng, and Xin-Qi Li.  Estimation of parameters in circuit qed by continuous quantum measurement.  Phys. Rev. A, 99:022114, 2019.
  * [60] S. Zhou, M. Zhang, J. Preskill, and L. Jiang.  Achieving the heisenberg limit in quantum metrology using quantum error correction.  Nature Communications, 78:78, 2018.
  * [61] S. Zhou, C.-L. Zou, and L. Jiang.  Saturating the quantum cramér–rao bound using locc.  Quantum Sci. Technol., 5:025005, 2020.

##  9 Appendix: Derivation of the measurement filter

We start by applying the algorithm [61] to our estimation problem. The natural
setup is to measure the noise units in the order in which they emerge in the
output, followed by a final measurement on the system+absorber. For simplicity
we refer to the latter as the system. We start by defining

$M^{(n)}\coloneqq\ket{\Psi(n)}\bra{\dot{\Psi}(n)}-\ket{\dot{\Psi}(n)}\bra{\Psi(n)}.$
(16)

where the index $n$ keeps track of the output length. In this section we will
use the label $0$ for the system, and $1,\dots,n$ for the noise units of the
output. The algorithm prescribes measurement bases
$\left\\{|e_{i}^{[l]}\rangle\right\\}_{i=1}^{k}$ for each of the output units
$l=1,\dots,n$ in an adaptive, sequential fashion. The first basis satisfies
the equations

$\left\langle
e^{[1]}_{i}\right|M^{(n)}_{1}\left|e^{[1]}_{i}\right\rangle=0\qquad\text{and}\qquad\left|\left\langle\left.e^{[1]}_{i}\right|\chi\right\rangle\right|^{2}=\frac{1}{d},\qquad{\rm
for~{}all~{}}i.$

where $M^{(n)}_{1}={\rm Tr}_{0,2,\dots,n}M^{(n)}$. The second equality follows
from the fact that at $\theta_{0}$ we have
$|\Psi(n)\rangle=|\psi\otimes\chi^{\otimes n}\rangle$. The next basis depends
on the outcome $i_{1}$ of the first measurement and satisfies the constraints

$\left\langle
e^{[2]}_{i}\right|M^{(n)}_{2}(i_{1})\left|e^{[2]}_{i}\right\rangle=0\qquad\text{and}\qquad\left|\left\langle\left.e^{[2]}_{i}\right|\chi\right\rangle\right|^{2}=\frac{1}{d},\qquad{\rm
for~{}all~{}}i$

where

$M^{(n)}_{2}(i_{1})={\rm Tr}_{0,3,\dots,n}\left\langle
e^{[1]}_{i_{1}}\right|M^{(n)}\left|e^{[1]}_{i_{1}}\right\rangle.$

Assuming the first $j<n$ units have been measured and a measurement record
$i_{\underline{j}}:=\\{i_{1},\dots,i_{j}\\}$ has been obtained, we denote

$M^{(n)}({i}_{\underline{j}})\coloneqq\left\langle e^{[1]}_{i_{1}}\otimes
e^{[2]}_{i_{2}}\otimes\dots\otimes
e^{[j]}_{i_{j}}\right|M^{(n)}\left|e^{[1]}_{i_{1}}\otimes
e^{[2]}_{i_{2}}\otimes\dots\otimes e^{[j]}_{i_{j}}\right\rangle,\qquad
M^{(n)}_{j+1}({i}_{\underline{j}}):={\rm
Tr}_{0,j+2,\dots,n}M^{(n)}({i}_{\underline{j}}).$

The measurement basis
$\left\\{\left|e^{[j+1]}_{i}({i}_{\underline{j}})\right\rangle\right\\}$ for
unit $j+1$ is then obtained by solving the constraints

$\left\langle
e^{[j+1]}_{i}({i}_{\underline{j}})\right|M^{(n)}_{j+1}({i}_{\underline{j}})\left|e^{[j+1]}_{i}({i}_{\underline{j}})\right\rangle=0\qquad\text{and}\qquad\left|\left\langle\left.e^{[j+1]}_{i}({i}_{\underline{j}})\right|\chi\right\rangle\right|^{2}=\frac{1}{D},\qquad{\rm
for~{}all~{}}i.$

The last step consists of measuring the system using the same procedure as for
the output units.

A priori, the procedure depends on the size $n$ of the output. The following
lemma shows that the optimal bases obtained for different output sizes
coincide.

######  Lemma 9.1.

Let $j,n$ be two output lengths with $j<n$ and consider applying the above
procedure to the corresponding states $|\Psi^{(j)}\rangle$ and respectively
$|\Psi^{(n)}\rangle$. The optimal measurement bases
$\\{|e^{[k]}_{i}\rangle\\}$ for the units $k=1,\dots,j$ satisfy the same
constraints and can be chosen to be the same. In addition we have

$M^{(n)}_{j}({i}_{\underline{j-1}})=M^{(j)}_{j}({i}_{\underline{j-1}}).$

_Proof._ Recall that $|\Psi(n)\rangle=W^{(n)}\dots
W^{(1)}|\psi\otimes\chi^{\otimes n}\rangle$, and let us denote

$P(n):=\ket{\psi\otimes\chi^{\otimes n}}\bra{\psi\otimes\chi^{\otimes
n}}=|\Psi_{\theta_{0}}(n)\rangle\langle\Psi_{\theta_{0}}(n)|.$

Then

$|\dot{\Psi}(n)\rangle=\sum_{i=1}^{n}W^{(n)}\dots\dot{W}^{(i)}\dots
W^{(1)}|\psi\otimes\chi^{\otimes
n}\rangle=\sum_{i=1}^{n}W^{(n)}\dots\dot{W}^{(i)}|\psi\otimes\chi^{\otimes
n}\rangle$

where we used the fact that $W$ leaves $|\psi\otimes\chi\rangle$ invariant.

We first show that the matrix $M_{1}^{(n)}$ does not depend on $n$, and
therefore the first measurement basis does not depend on the length of the
output.

$\displaystyle M^{(n)}_{1}={\rm Tr}_{0,2,\dots,n}\left[M^{(n)}\right]$
$\displaystyle=\sum_{i}\text{Tr}_{0,2,\dots,n}\left[\ W^{(n)}\cdots
W^{(i)}\cdots W^{(1)}P(n)W^{(1)*}\cdots\dot{W}^{(i)*}\cdots W^{(n)*}-c.c\
\right]$ $\displaystyle=\sum_{i}\text{Tr}_{0,2,\dots,n}\left[\ W^{(n)}\cdots
W^{(i)}P(n)\dot{W}^{(i)*}\cdots W^{(n)*}-c.c\ \right]$
$\displaystyle=\sum_{i}\text{Tr}_{0,2,\dots,n}\left[\
W^{(i)}P(n)\dot{W}^{(i)*}-c.c\ \right]$ $\displaystyle={\rm
Tr}_{0}\left[W^{(1)}P(1)\dot{W}^{(1)*}-c.c\right]+\sum_{i=2}\text{Tr}_{0,2,\dots,n}\left[W{i}P(n)\dot{W}^{(i)*}-c.c\right]$
$\displaystyle={\rm Tr}_{0}\left[M^{(1)}\right]=M^{(1)}_{1},$ (17)

where $c.c,$ denotes the adjoint. In the third equality we used

${\rm Tr}_{0,2}[W_{0,2}A_{0,1,2}W^{*}_{0,2}]={\rm Tr}_{0,2}[A_{0,1,2}]$

where $W_{0,2}$ is a unitary acting on subsystems $0,2$ of a tripartite
system, and $A_{0,1}$ acts on subsystems $0,1$. In the last equality we used
our assumption $\bra{\psi\otimes\chi}W^{*}\dot{W}\ket{\psi\otimes\chi}=0$.

Let $\left\\{\left|e^{[1]}_{i}\right\rangle\right\\}$ be the measurement basis
determined from $M_{1}=M_{1}^{(n)}=M_{1}^{(j)}$. We now show that, conditional
on the outcome $i_{1}$ of this measurement, the second basis
$\left\\{\left|e^{[2]}_{i}\right\rangle\right\\}$ does not depend on the
length of the output, which can be take to be equal to 2.

$\displaystyle M_{2}^{(n)}(i_{1})$ $\displaystyle={\rm
Tr}_{0,3,\dots,n}\left[M^{(n)}(i_{1})\right]$ $\displaystyle=\sum_{i}{\rm
Tr}_{0,3,\dots,n}\left[\left\langle e^{[1]}_{i_{1}}\right|W^{(n)}\cdots
W^{(i)}\cdots W^{(1)}P(n)W^{(1)*}\cdots\dot{W}^{(i)*}\cdots
W^{(n)*}\left|e_{i_{1}}^{[1]}\right\rangle-c.c\right]$ $\displaystyle={\rm
Tr}_{0,3,\dots,n}\left[\left\langle e^{[1]}_{i_{1}}\right|W^{(n)}\cdots
W^{(1)}P(n)\dot{W}^{(1)*}\cdots
W^{(n)*}\left|e_{i_{1}}^{[1]}\right\rangle-c.c\right]$
$\displaystyle\quad+{\rm Tr}_{0,3,\dots,n}\left[\left\langle
e^{[1]}_{i_{1}}\right|W^{(n)}\cdots
W^{(2)}W^{(1)}P(n)W^{(1)*}\dot{W}^{(2)*}\cdots
W^{(n)*}\left|e_{i_{1}}^{[1]}\right\rangle-c.c\right]$
$\displaystyle\qquad+\sum_{i=3}^{n}{\rm Tr}_{0,3,\dots,n}\left[\left\langle
e^{[1]}_{i_{1}}\right|W^{(n)}\cdots W^{(i)}\cdots
W^{(1)}P(n)W^{(1)*}\cdots\dot{W}^{(i)*}\cdots
W^{(n)*}\left|e_{i_{1}}^{[1]}\right\rangle-c.c\right]$ $\displaystyle={\rm
Tr}_{0}\left[W^{(2)}\left\langle
e^{[1]}_{i_{1}}\right|W^{(1)}P(2)\dot{W}^{(1)*}\left|e_{i_{1}}^{[1]}\right\rangle
W^{(2)*}-c.c\right]$ $\displaystyle\quad+\frac{1}{D}{\rm
Tr}_{0}\left[W^{(1)}P(1)\dot{W}^{(1)*}-c.c\right]$ $\displaystyle={\rm
Tr}_{0}\left[W^{(2)}\left(\left\langle
e^{[1]}_{i_{1}}\right|M^{(1)}\left|e_{i_{1}}^{[1]}\right\rangle\otimes
P_{\chi}\right)W^{(2)*}\right]+\frac{1}{D}{\rm
Tr}_{0}\left[M^{(1)}\right]=M_{2}^{(2)}(i_{1}).$

The equalities follow in the same way as in (17), and in addition we used
$|\braket{e^{[1]}_{i_{1}}}{\chi}|^{2}=\frac{1}{D}$ in the third equality.

Using the same techniques as above we obtain the general statement

$\displaystyle M_{j}^{(n)}(i_{\underline{j-1}})={\rm
Tr}_{0,j+1,\dots,n}\left[M^{(n)}({i}_{\underline{j-1}})\right]$
$\displaystyle=\frac{1}{D^{j-1}}{\rm Tr}_{0}\left[M^{(1)}\right]$
$\displaystyle+{\rm
Tr}_{0}\left[W^{(j)}\left(\bra{e^{[\underline{j-1}]}_{{i}_{\underline{j-1}}}}M^{(j-1)}({i}_{\underline{j-2}})\ket{e^{[\underline{j-1}]}_{{i}_{\underline{j-1}}}}\otimes
P_{\chi}\right)W^{(j)*}\right]$ $\displaystyle={\rm
Tr}_{0}\left[M^{(j)}({i}_{\underline{j-1}})\right]=M_{j}^{(j)}(i_{\underline{j-1}}).$
(18)

∎

Next, we show that we can express $M^{(j)}({i}_{\underline{j-1}})$ in terms of
$M^{(1)}$, $M^{(j-1)}({i}_{\underline{j-2}})$ and $\ket{e^{[j-1]}_{j-1}}$.
Indeed by writing $W(j)=W^{(j)}W(j-1)$ we have

$\displaystyle M^{(j)}({i}_{\underline{j-1}})$ $\displaystyle:=\left\langle
e^{[\underline{j-1}]}_{{i}_{\underline{j-1}}}\right|M^{(j)}\left|e^{[\underline{j-1}]}_{{i}_{\underline{j-1}}}\right\rangle$
$\displaystyle=\left\langle
e^{[\underline{j-1}]}_{{i}_{\underline{j-1}}}\right|W^{(j)}W(j-1)P(j)W^{*}(j-1)\dot{W}^{(j)*}\left|e^{[\underline{j-1}]}_{{i}_{\underline{j-1}}}\right\rangle-c.c$
$\displaystyle\qquad+\left\langle
e^{[\underline{j-1}]}_{{i}_{\underline{j-1}}}\right|W^{(j)}W(j-1)P(j)\dot{W}(j-1)^{*}W^{(j)*}\left|e^{[\underline{j-1}]}_{{i}_{\underline{j-1}}}\right\rangle-c.c$
$\displaystyle=\frac{1}{D^{j-1}}M^{(1)}+W^{(j)}\left(\left\langle
e^{[j-1]}_{i_{j-1}}\right|M^{(j-1)}(i_{\underline{j-2}})\left|e^{[j-1]}_{i_{j-1}}\right\rangle\otimes
P_{\chi}\right)W^{(j)*}$ (19)

This can then be used to determine the next measurement basis, producing the
iterative procedure described in section 5, which consists in updating the
‘filter’ that determines the optimal basis at each time step using the last
measurement outcome $|e^{[j]}_{i_{j}}\rangle$.

Let $A_{1}=M^{(1)}$ and $A_{j}:=M^{(j)}(i_{\underline{j-1}})$ for $j>1$, and
denote

$\Pi_{j}:=\left\langle
e^{[j]}_{i_{j}}\right|A_{j}\left|e^{[j]}_{i_{j}}\right\rangle=\left\langle
e^{[j]}_{i_{j}}\right|M^{(j)}({i}_{\underline{j-1}})\left|e^{[j]}_{i_{j}}\right\rangle$

Then equation (19) can be written as

$A_{j}=\frac{1}{D^{j-1}}A_{1}+U^{(j)}\left(\Pi_{j-1}\otimes
P_{\chi}\right)U^{(j)*}$

The optimal measurement is obtained by applying the conditions to the operator
$B_{j}:={\rm Tr}_{0}A_{j}$.

##  10 Appendix: Proof of Proposition 6.1

We start by using a general Fisher information identity for bipartite systems.
Consider a generic pure state model
$|\psi_{\theta}\rangle\in\mathcal{H}_{s}\otimes\mathcal{H}_{o}$ and let
$\\{|e^{o}_{i}\rangle\\}$ and $\\{|e^{s}_{j}\rangle\\}$ be optimal bases in
the ‘output’ and ‘system’ subsystems, for estimating $\theta$ at a particular
value $\theta_{0}$, as prescribed by [61]. Assume we perform the ‘output’
measurement and let $X$ denote the outcome whose distribution is

$\mathbb{P}_{\theta}(X=i)=\langle\psi_{\theta}|\mathbf{1}\otimes
P_{i}|\psi_{\theta}\rangle,\qquad P_{i}=|e^{o}_{i}\rangle\langle e^{o}_{i}|.$

The conditional state of the ‘system’ given $X=i$ is

$|\psi_{\theta}(i)\rangle=\frac{\mathbf{1}\otimes
P_{i}|\psi_{\theta}\rangle}{\|\mathbf{1}\otimes P_{i}\psi_{\theta}\|}$

and this state contains the ‘remaining’ information about $\theta$. The
following inequality bounds the total available information as

$I_{\theta}(X)+\mathbb{E}_{X}F(\psi_{\theta}(X))\leq F_{\theta}$

where the first term on the left side is the classical Fisher information of
the outcome distribution $\mathbb{P}_{\theta}$ and the second is the expected
QFI of the conditional state $|\phi^{i}_{\theta}\rangle$. Consider now the
second measurement $|e^{s}_{j}\rangle$ on the system and let $Y$ be its
outcome. Then

$F_{\theta}=I_{\theta}(X)+\mathbb{E}_{X}I_{\theta}(Y|X)\leq
I_{\theta}(X)+\mathbb{E}_{X}F(\psi_{\theta}(X))\leq F_{\theta}$

where the first equality is due to measurement optimality, while the second is
the inequality between classical and quantum information. This implies that

$F_{\theta}-I_{\theta}(X)=\mathbb{E}_{X}F_{\theta}(\psi_{\theta}(X)).$ (20)

We now consider the Markov setup in which the system+absorber play the role of
‘system’ while the $n$ noise units are the ‘output’. We assume that the the
output and the system+absorber are measured according to the optimal scheme
presented in section 5. The joint state is given by

$|\Psi_{\theta}(n)\rangle=W_{\theta}^{(n)}\cdot\dots\cdot
W_{\theta}^{(1)}|\psi\otimes\chi^{\otimes n}\rangle$

and the conditional states are

$|\psi_{\theta}(i_{1},\dots,i_{n})\rangle=\frac{K^{[n]}_{\theta,i_{n}}\dots
K^{[1]}_{\theta,i_{1}}|\psi\rangle}{\sqrt{p_{\theta}(i_{1},\dots i_{n})}}.$

We will show that at $\theta=\theta_{0}$ the left side of (20) is bounded by a
constant which does not depend on $n$. For simplicity, whenever possible we
will use the compact notations such as

$K_{\bf i}:=K^{[n]}_{\theta,i_{n}}\dots
K^{[1]}_{\theta,i_{1}},\qquad\mathrm{and}\qquad|e_{\bf
i}\rangle=|e^{[n]}_{i_{n}}\otimes\dots\otimes e^{[1]}_{i_{1}}\rangle.$

Recall that by design the following condition holds
$W_{\theta_{0}}|\psi\otimes\chi\rangle=|\psi\otimes\chi\rangle$, which implies
$|\Psi_{\theta_{0}}(n)\rangle=|\psi\otimes\chi^{\otimes n}\rangle$ and also
$K^{[j]}_{i_{j}}|\psi\rangle=c^{[j]}_{i_{j}}|\psi\rangle$ for some constants
$c^{[j]}_{i_{j}}$.

Recall that for a pure state model $|\psi_{\theta}\rangle$ the QFI is given by

$F_{\theta}=4\left(\|\dot{\psi}_{\theta}\|^{2}-|\langle\dot{\psi}_{\theta}|\psi_{\theta}\rangle|^{2}\right)=4\|\psi_{\theta}^{\perp}\|^{2},\qquad|\psi_{\theta}^{\perp}\rangle=|\dot{\psi}_{\theta}\rangle-
P_{\psi_{\theta}}|\dot{\psi}_{\theta}\rangle.$

Therefore, expected system QFI on the left side of (20) is given by

$F_{s}(\theta_{0})=4\sum_{\bf i}p({\bf i})\|\psi^{\perp}_{\theta_{0}}({\bf
i})\|^{2}$ (21)

where

$|\psi^{\perp}_{\theta_{0}}({\bf i})\rangle=|\dot{\psi}_{\theta_{0}}({\bf
i})\rangle-P_{\psi}|\dot{\psi}_{\theta_{0}}({\bf i})\rangle$

For simplicity we now drop the subscript $\theta_{0}$ and we have

$|\dot{\psi}({\bf i})\rangle=\frac{\dot{K}_{\bf i}|\psi\rangle}{\sqrt{p({\bf
i})}}-\frac{1}{2}\frac{K_{\bf i}|\psi\rangle}{p^{3/2}({\bf i})}\dot{p}({\bf
i})=\frac{\dot{K}_{\bf i}|\psi\rangle}{\sqrt{p({\bf
i})}}-\frac{1}{2}\frac{c({\bf i})\dot{p}({\bf i})}{p^{3/2}({\bf
i})}|\psi\rangle$

and therefore

$|\psi^{\perp}({\bf i})\rangle=\frac{1}{\sqrt{p({\bf
i})}}\left(I-P_{\psi}\right)\dot{K}_{\bf i}|\psi\rangle.$

Equation (21) becomes

$F_{s}(\theta_{0})=4\sum_{\bf i}\left\|P^{\perp}_{\psi}\dot{K}_{\bf
i}\psi\right\|^{2}=4\left\|(P^{\perp}_{\psi}\otimes
I)\dot{\Psi}(n)\right\|^{2}.$

Now

$|\dot{\Psi}(n)\rangle=\sum_{j=1}^{n}W^{(n)}\dots
W^{(j+1)}\dot{W}^{(j)}W^{(j-1)}\dots W^{(1)}|\psi\otimes\chi^{\otimes
n}\rangle=\sum_{j=1}^{n}W^{(n)}\dots
W^{(j+1)}\dot{W}^{(j)}|\psi\otimes\chi^{\otimes n}\rangle$

and by triangle inequality

$F_{s}(\theta_{0})\leq 4\left(\sum_{j=1}^{n}\left\|(P^{\perp}_{\psi}\otimes
I)W^{(n)}\dots W^{(j+1)}\dot{W}^{(j)}\psi\otimes\chi^{\otimes
n}\right\|\right)^{2}$

Let $|\dot{\Psi}(1)\rangle:=\dot{W}|\psi\otimes\chi\rangle$ and let
$\tau:={\rm Tr}_{1}(\dot{\Psi}(1)\rangle\langle\dot{\Psi}(1)|)$. Then

$\displaystyle\left\|(P^{\perp}_{\psi}\otimes I)W^{(n)}\dots
W^{(j+1)}\dot{W}^{(j)}\psi\otimes\chi^{\otimes n}\right\|^{2}$
$\displaystyle=$ $\displaystyle{\rm Tr}_{0,1,\dots
n-j}\left((P^{\perp}_{\psi}\otimes I)W^{(n-j)}\dots W^{(1)}\tau W^{(1)*}\dots
W^{(n-j)*}\right)$ $\displaystyle=$ $\displaystyle{\rm
Tr}_{0}(P^{\perp}_{\psi}T^{n-j}(\tau))$

where in the last equality we have uses the definition of the transition
operator $T$ of the system+absorber. Assuming that $T$ is ergodic we have

$T^{n}(\tau)\to P_{\psi}$

exponentially fast with $n$ so that

${\rm Tr}_{0}(P^{\perp}_{\psi}T^{n-j}(\tau))\leq a^{2(n-j)}$

for some $a<1$. Therefore

$F_{s}(\theta_{0})\leq 4\left(\sum_{j=1}^{n}a^{n-j}\right)^{2}\leq
4\frac{1}{(1-a)^{2}}$

Note that if the spectral gap of $T$ becomes small then the convergence to
stationarity is slower and the upper bound increases.

∎

##  11 Appendix: Proof of Lemma 6.2

The CFI of any output measurement can be computed explicitly by writing

$\displaystyle\left.\frac{d}{d\theta}p_{\theta}(i_{1},\dots
i_{n})\right|_{\theta_{0}}$ $\displaystyle=$
$\displaystyle\frac{d}{d\theta}\left.\left\|K^{[n]}_{i_{n}}\cdots
K^{[1]}_{i_{1}}\psi\right\|^{2}\right|_{\theta_{0}}$ $\displaystyle=$
$\displaystyle 2{\rm Re}\sum_{j=1}^{n}\langle\psi|K^{[1]*}_{i_{1}}\cdots
K^{[n]*}_{i_{n}}K^{[n]}_{i_{n}}\cdots\dot{K}^{[j]}_{i_{j}}\cdots
K^{[1]}_{i_{1}}|\psi\rangle$ $\displaystyle=$
$\displaystyle\left|c^{[1]}_{i_{1}}\cdots c^{[n]}_{i_{n}}\right|^{2}\cdot
2{\rm Re}\sum_{j=1}^{n}\frac{\langle\psi|K^{[n]}_{i_{n}}\dots
K^{[j+1]}_{i_{j+1}}\dot{K}^{[j]}_{i_{j}}|\psi\rangle}{c^{[j]}_{i_{j}}\dots
c^{[n]}_{i_{n}}}$

Therefore

$\displaystyle I^{({\rm out})}_{\theta_{0}}(n)$ $\displaystyle=$
$\displaystyle\sum_{i_{1},\dots,i_{n}}\left|c^{[1]}_{i_{1}}\cdots
c^{[n]}_{i_{n}}\right|^{2}\left(2{\rm
Re}\sum_{j=1}^{n}\frac{\langle\psi|K^{[n]}_{i_{n}}\dots
K^{[j+1]}_{i_{j+1}}\dot{K}^{[j]}_{i_{j}}|\psi\rangle}{c^{[j]}_{i_{j}}\dots
c^{[n]}_{i_{n}}}\right)^{2}$ (22) $\displaystyle=$
$\displaystyle\sum_{i_{1},\dots,i_{n}}p_{\theta_{0}}(i_{1},\dots,i_{n})f^{2}(i_{1},\dots,i_{n})=\mathbb{E}_{\theta_{0}}(f^{2})$

where $f$ is the function

$f(i_{1},\dots,i_{n})=2{\rm
Re}\sum_{j=1}^{n}\frac{\langle\psi|K^{[n]}_{i_{n}}\dots
K^{[j+1]}_{i_{j+1}}\dot{K}^{[j]}_{i_{j}}|\psi\rangle}{c^{[j]}_{i_{j}}\dots
c^{[n]}_{i_{n}}}$

∎

##  12 Appendix: Computation of finite time system-output QFI

From (5) and (1) we have

$\displaystyle F^{({\rm s+o})}_{\theta}(n)$ $\displaystyle=$ $\displaystyle
4\|\dot{\Psi}_{\theta}(n)\|^{2}$ (23) $\displaystyle=$ $\displaystyle
4\sum_{i_{1},\dots
i_{n}}\left\|\sum_{j=1}^{n}K_{i_{n}}\dots\dot{K}_{i_{j}}\dots
K_{i_{1}}\psi\right\|^{2}$

where $K_{i}$ are the (fixed) Kraus operators with respect to the standard
basis, and $|psi\rangle=|0\rangle$. Our specific model ha the feature that
both $K_{i}$ and $\dot{K}_{i}$ map the basis vectors into each other:

$\displaystyle K_{0}|0\rangle$ $\displaystyle=|0\rangle$
$\displaystyle\dot{K}_{0}|0\rangle$ $\displaystyle=i|1\rangle$ $\displaystyle
K_{0}|1\rangle$ $\displaystyle=\sqrt{1-\lambda}|1\rangle$
$\displaystyle\dot{K}_{0}|1\rangle$ $\displaystyle=i\sqrt{1-\lambda}|0\rangle$
$\displaystyle K_{1}|0\rangle$ $\displaystyle=0$
$\displaystyle\dot{K}_{1}|0\rangle$ $\displaystyle=|1\rangle$ $\displaystyle
K_{1}|1\rangle$ $\displaystyle=\sqrt{\lambda}e^{i\phi}|0\rangle$
$\displaystyle\dot{K}_{1}|1\rangle$ $\displaystyle=0$

This allows to compute the QFI explicitly by noting the terms in the sum (23)
with more that two indices equal to $1$ have zero contribution. The remaining
terms can be computed as follows. The term with only zero indices is

$F^{(0)}=4\left\|\sum_{j=1}^{n}K_{0}\dots\dot{K}_{0}\dots
K_{0}|0\rangle\right\|^{2}=4\left|\sum_{j=1}^{n}ia^{n-j}\right|^{2}=4\left(\frac{1-a^{n}}{1-a}\right)^{2}$
(24)

where $a=\sqrt{1-\lambda}$.

Consider now a sequence $(0,\dots,,0,1,0,\dots 0)$ with a single one on
position $l$. Since $K_{0}\dots\dot{K}_{0}\dots K_{1}\dots K_{0}|0\rangle=0$
the only contributing terms will be those with derivative on the first $(l-1)$
$K_{0}$s or on $K_{1}$. This gives

$\displaystyle F^{(1)}$ $\displaystyle=$ $\displaystyle
4\sum_{l=2}^{n}\left\|i\sum_{r=1}^{l-1}e^{i\phi}a^{l-1-r}b|0\rangle+a^{n-l}|1\rangle\right\|^{2}+4\|a^{n-1}|1\rangle\|^{2}$
(25) $\displaystyle=$ $\displaystyle
4\sum_{l=2}^{n}\left(b^{2}\frac{1-a^{l-1}}{1-a}+a^{2(n-l)}\right)+4a^{2(n-1)}$
$\displaystyle=$ $\displaystyle
4(n-1)\frac{b^{2}}{(1-a)^{2}}+4\frac{a^{2}-a^{2n}}{(1-a)^{2}}-8\frac{b^{2}(a-a^{n})}{(1-a)^{3}}+4\frac{1-a^{2(n-1)}}{b^{2}}+4a^{2(n-1)}$

where $b=\sqrt{\lambda}$.

Finally, consider the sequences of the type $(0,\dots 0,1,0,\dots 0,1,0\dots
0)$ with 1s on positions $1\leq i<k\leq n$. In this case the nonzero
contributions come from terms where the derivative is on positions $j=i$. The
Fisher contribution is

$\displaystyle F^{(2)}$ $\displaystyle=$ $\displaystyle 4\sum_{1\leq i<k\leq
n}\|e^{i\phi}\sqrt{\lambda}(1-\lambda)^{(k-i-1)/2}|0\rangle\|^{2}$ (26)
$\displaystyle=$ $\displaystyle 4\sum_{1\leq i<k\leq
n}b^{2}a^{2(k-i-1)}=4(n-1)-4\frac{a^{2}-a^{2n}}{b^{2}}.$

Adding together the contributions (24), (25) and (26) we obtain the total QFI

$\displaystyle F_{\theta}^{({\rm s+o})}(n)$ $\displaystyle=$ $\displaystyle
F^{(0)}+F^{(1)}+F^{(2)}$ (27) $\displaystyle=\frac{8n}{1-a}$
$\displaystyle+4\left[\left(\frac{1-a^{n}}{1-a}\right)^{2}-\frac{b^{2}}{(1-a)^{2}}+\frac{a^{2}-a^{2n}}{(1-a)^{2}}-2\frac{b^{2}(a-a^{n})}{(1-a)^{3}}+\frac{1-a^{2(n-1)}}{b^{2}}+a^{2(n-1)}\right]$
$\displaystyle-4\left[1+\frac{a^{2}-a^{2n}}{b^{2}}\right]$

where the leading term is consistent with the QFI rate formula (14).