Standard Quantum Limit and Heisenberg Limit in Function Estimation
Naoto Kura Masahito Ueda






 Unlike well-established parameter estimation, function estimation faces conceptual and mathematical difficulties despite its enormous potential utility. 
 We establish the fundamental error bounds on function estimation in quantum metrology for a spatially varying phase operator, where various degrees of smooth functions are considered.
 
 The error bounds are identified in both cases of absence and presence of interparticle entanglement, which correspond to the standard quantum limit and the Heisenberg limit, respectively.
 Notably, these error bounds can be reached by either position-localized states or wavenumber-localized ones.
 In fact, we show that these error bounds are theoretically optimal for any type of probe states, indicating that quantum metrology on functions is also subject to the Nyquist-Shannon sampling theorem, even if classical detection is replaced by quantum measurement.  

 
    Masahito Ueda
    December 30, 2023
=====================




 Accurate estimation of signals with a limited amount of resource is a fundamental problem in physics.
 Quantum metrology has made a profound contribution to this problem by demonstrating a classically unattainable scaling of the estimation error 
 <cit.>.
 This non-classical accuracy, called the Heisenberg limit, can be reached by various forms of quantum features including entanglement <cit.>, quantum circuit <cit.>, bosonic <cit.> and spin <cit.> squeezing, and quantum statistics <cit.>.
Applications of the Heisenberg-limited measurement range from detection of fundamental signals such as atomic clocks <cit.>, gravitational waves <cit.>, a scalable cat state <cit.> and polariton condensates <cit.>.


 The Heisenberg limit is significantly better than the standard quantum limit (SQL) δ = O(N^-1/2), which sets the accuracy bound arising from uncorrelated noises, where N is the size of resource.
 While the Heisenberg limit δ = O(N^-1) and the SQL δ = O(N^-1/2) apply to both scalar estimation and vector estimation <cit.>, it is generally considered that the continuity of the signal may alter the scaling law.
 Such a problem can be categorized as function estimation, which has attracted growing attention.
 For example, atomic clocks <cit.> and graviational waves <cit.> involves time-varying singals, which offer richer information when treated as functions.
 More generally, exploration for new phenomena involves observing structures in a continuous in space and/or time, which can be represented as functions.
 This indicates that functional structures play crucial roles in generic continuous systems including measurements on magnetometry <cit.>, nanostructured materials <cit.>, live cells <cit.>, and event horizons <cit.>.
 Hence, it is not only a fundamental question but also relevant to a wide range of applications to ask how quantum metrology can contribute to the detection of functions with ultimate accuracy.

 The quantum version of function estimation has been investigated in terms of the signal detection theory in Refs. <cit.>.  In fact, weaker scaling laws are implied when the target parameter can change continuously in time, such as a Gaussian signal.
 The demonstration of such unconventional limits has recently become within the experimental reach due to the realization of, e.g., high-N00N states <cit.> and optical phase tracking <cit.>.
 Although the detection theory is applicable to stochastic noises, it does not support the case where the relevant parameter is not inherently stochastic, which is often the case with quantum imaging and quantum signal processing <cit.>.

 In this Letter, we present a fundamental framework of quantum metrology on functions.
 Unlike parametric estimation, function estimation involves infinite degrees of freedom and inevitably requires further assumptions on the target function.
 Assuming only the smoothness of the function, we find the SQL of O(N^-q/(2q+1)) and the Heisenberg limit of O(N^-q/(q+1)), where q indicates the degree of smoothness of the function.
 Our framework allows analysis of estimation errors of data series under given smoothness, such as a bound on the amplitude of derivatives.
 This includes the previous results on Gaussian processes through computation of their smoothness <cit.>, as demonstrated later.
 The data series requires neither to have a prior distribution nor even to be continuous, allowing, for example, a sample with a finite number of discontinuous points <cit.>.


 Moreover, we have found that the error limit can equally be saturated by states which are localized in position or wavenumber.
 This result implies the equivalence between space discretization and momentum cutoff in quantum information processing, reminiscent of the Nyquist–Shannon sampling theorem in classical statistics.

 

  
Setup.
 
 We consider the estimation of an unknown function φ(x) defined over the interval 0≤ x≤ L by using the position-dependent phase-shift gate Û_φ shown in Fig. <ref>.
 For simplicity, we also assume the periodic boundary condition.

 To estimate φ, we prepare an N-particle state as a probe, which evolves according to the unitary operator Û_φ and then is measured.
 These particles are distributed in the interval [0,L] and have two internal states: one state |+⟩ interacts with the phase-shift gate Û_φ, and the other state |-⟩ does not.

 We work in the first-quantization formalism and denote by |x;±⟩, which is the position eigenstate at x with internal state ±.
 Let the phase-shift gate act as
 
  Û_φ |x; +⟩ = e^iφ(x) |x; +⟩
 
 and
 
  Û_φ |x; -⟩ = |x; -⟩.
 
 The unknown function φ can be estimated from measurement on this output.
 When the probe is composed of N separable particles, the error of the function estimation is bounded by the SQL.
 A probe with appropriately entangled N particles, on the other hand, leads to the Heisenberg limit.

 We recall that the estimation error of a scalar parameter θ is computed as δ^2 = 3E[θ̃-θ^2], where θ̃ is the estimator depending on the stochastic nature of measurement outcome.
 Similarly, we consider a stochastic estimator φ̃ for the function, and compute the mean-square periodic error (MSPE)  <cit.> as
 
    δ^2 = 3E[ ∫_0^L dx/L [φ̃(x)-φ(x)]_2π^2 ].

 In other words, the estimation error is averaged over x and the modulus is replaced by [φ̃(x)-φ(x)]_2π = min_n∈3Z[]φ̃(x)-φ(x)+2π n, i.e. the minimal modulus regarding the 2π periodicity of a phase.

 The main difficulty in function estimation lies in the fact that the problem involves an the infinite degrees of freedom.
 In particular, the lower bound on δ^2 cannot be established for an arbitrary function, since we cannot exclude any rapidly fluctuating functions from a finite number of measurements.
 Hence, we impose the following constraint on the target function φ:
 
    ∫_0^L dx/L[]φ'(x)^2 ≤M^2/L^2

 for some positive number M>0.
 With this constraint, we can establish a suitable lower bound on sufficiently smooth and slowly varying functions φ.

 The condition (<ref>) can be applied only when the target function is differentiable.
 In our framework, we consider more general functions without differentiability: the Hölder continuity φ(x+ϵ)-φ(x) = O(ϵ^q) for fixed 0<q≤ 1 <cit.>.
 More rigorously, we impose a general constraint:
 
    sup_0<ϵ<a∫_0^L dx/L[]φ(x+ϵ)-φ(x)/ϵ^q^2
      ≤M^2/L^2q,

 where a > 0 is a constant that does not affect the estimation error in the limit of large N.
 The special case with q=1 reduces to Eq. (<ref>).
 
 

 

  
Estimation methods.
 Given the target function φ under the constraint (<ref>), there exist estimation methods that ensure a finite estimation error δ defined in Eq. (<ref>).
 We here compare the following two different methods:

 Position-state (PS) method – We estimate the individual phases φ(x_j) at several positions x_j, and then computationally reconstruct the entire function.
             
             

 Wavenumber-state (WS) method – We prepare a sufficiently large number of wavefunctions ψ(x) ∝ e^iφ(x), and estimate the function φ by reconstructing the quantum state |ψ⟩ by the quantum tomography.

 We find that the numbers of particles N required for these two methods are the same up to a constant factor.
 
 

  
Position-state method.
 The PS method can be used when the target function is relatively small, say, φ(x)≤π/3 for all x.
 In this case, we can circumvent the phase wrapping problem and employ a method analogous to the kernel density estimation <cit.>.


 In the first step, we sample n_1 positions x_1, x_2, …, x_n_1 in the interval 0≤ x≤ L with equal spacing. 
 Then, the phase φ_j = φ(x_j) at each position x_j is measured by using n_2 particles localized at x_j.
 The estimation error δ_1ind of the individual phase φ_j is known <cit.> as it is the quantum metrology on a scalar parameter; the SQL δ_1ind = O(n_2^-1/2) is established by the probe 1/√(2)(|x_j;-⟩+|x_j;+⟩)^⊗ n_2 and the Heisenberg limit δ_1ind = O(n_2^-1) by 1/√(2)(|x_j;-⟩^⊗ n_2+|x_j;+⟩^⊗ n_2).
 
 Finally, the function estimator φ̃(x) is computed from the individual estimators φ̃_j by local linear smoothing <cit.>:
 
    φ̃(x) = ∑_j=1^n_1φ̃_jf(x-x_j),

 where f is a smoothing function.
 In the present case, we may just set f(x) = 1 for x ≤L/2n_1 and f(x)=0 otherwise.
 This corresponds to the approximation by the value at the nearest site, i.e., we set φ̃(x) = φ̃_j where x_j is the point nearest to x.

 The estimation error can be decomposed into two parts: the statistical error δ_1stat caused by the measurement, and the deterministic error δ_1det due to smoothing.
 The balance between these errors can be tuned by the width l of smoothing.
 The estimated value φ̃_j(x) is of the same order as δ_1ind: δ_1stat = O(n_2^-1/2) for the SQL and δ_1stat = O(n_2^-1) for the Heisenberg limit.




 On the other hand, the deterministic error δ_1det is the variation of φ(x) within the width L/2n_1, which turns out to be
 
  δ_1det = O(n_1^-qM)
 
 by virtue of the constraint (<ref>).
 
 For a given number of particles N = n_1n_2, the optimal accuracy is determined by the trade-off between δ_1stat and δ_1det.
 As a consequence of Young's inequality, we obtain
 
    δ≥ O(n_1^-q M) + O(n_2^-1/2)
       ≥ O( (M^1/qN^-1)^q/(2q+1))

 for the SQL and
 
    δ≥ O(n_1^-q M) + O(n_2^-1)
       ≥ O( (M^1/qN^-1)^q/(q+1))

 for the Heisenberg limit.








 Therefore, the overall estimation error δ is significantly larger than the traditional quantum limit, which is an expected feature of the function estimation.
 We note that entanglement of particles in different positions is not necessary to achieve the Heisenberg limit; such intersite entanglement does not enhance the estimation of linear parameters, as suggested in studies of quantum network sensors <cit.>.
 

  
Wavenumber-state method.
 In the WS method, we begin with the wavenumber eigenstate with zero eigenvalue: 
 
  ∫_0^L dx/√(2L)[|x;-⟩^⊗n_1p + |x;+⟩^⊗n_1p].
 
 We use the one-particle state (n_1p=1) for the SQL and a multipartite EPR state (n_1p>1) for the Heisenberg limit.
 
 By the phase-shift gate Û_φ, one obtains the output probe state
 
    |S_φ⟩ =
        ∫_0^L dx/√(2L)[|x;-⟩^⊗ n_1p + e^in_1pφ(x)|x;+⟩^⊗ n_1p].

 The estimation is conducted by reconstructing |S_φ⟩ as accurate as possible by measuring n_1c copies of the probe state.
 For this purpose, we consider the projection P_K onto the subspace of wavenumbers k such that k≤ 2π K/L.
 Since the postselected state |S^*_φ⟩∝ P_K|S_φ⟩ belongs to a (2K+1)-dimensional Hilbert space, it can be identified by the quantum tomography.


 The error of the state reconstruction can be quantified by the infidelity
 
  1 - ⟨ S_φ | S_φ̃ ⟩,
  
 
 where |S_φ̃⟩ denotes the reconstructed state.
 In fact, we show in the Supplemental Material 
 [See Supplemental Material for the rigorous derivation of 
 the biased Cramér-Rao bound, details of the explicit estimation procedures, and an application of Kitaev's method on function estimation.]
 that the MSPE has can be bounded by the expected infidelity as
 
    δ^2    ≤π^2/n_1p^23E[
        1 - ⟨ S_φ| S_φ̃⟩
      ] 
       ≤π^2/n_1p^2(1 - ⟨ S_φ| S^*_φ⟩^2) +
       π^2/n_1p^23E[
        1 - ⟨ S_φ^* | S_φ̃⟩^2
      ] 
       = δ^2_1PS + δ^2_1QT.

 
 Here, the error is divided into the postselection part δ^2_1PS and the quantum-tomography part δ^2_1QT.





















 The postselection error can be bounded by the constraint (<ref>) as δ_1PS≤ O(K^-qM) <cit.>,
 while the results of the finite-dimensional tomography imply n_1pδ_1QT≤ O(K^1/2 n_1c^-1/2) <cit.>.

 When n_1p = 1 and N = n_1c, we have the trade-off relation between δ_1PS and δ_1QT for the SQL:
 
    δ_1PS = O(K^-qM),  δ_1QT = O(K^1/2N^-1/2).

 By setting n_1=K and n_2 = N/K, the errors
 δ_1PS and δ_1QT can be mapped to the errors
 δ_1det and δ_1stat in the PS method, respectively.
 Therefore, the SQL in the WS method reduces to that in Eq. (<ref>) obtained by the PS method.

 The error bound can be lowered for n_1p > 1 and N = n_1p n_1c,
 while K ≤ O(n_1c) must be maintained in order to robustly conduct the quantum tomography.
 Therefore, the optimal trade-off relation for the Heisenberg limit is
 
    δ_1PS = O(n_1c^-qM),  δ_1QT = O(n_1p^-1).

 By setting n_1=n_1c and n_2=n_1p, this tradeoff relation corresponds exactly to the PS method, and we obtain the same Heisenberg limit as Eq. (<ref>).
 However, there is a caveat that, in the output state (<ref>) the phase ambiguity of φ(x) modulo 2π/n_1p needs to be removed.
 We show in the Supplemental Material <cit.> that this removal can be handled by analogy with the Kitaev's method <cit.>.
 



 
 

  
Optimality of the SQL.
 We have preposed the SQL (<ref>) and the Heisenberg limits (<ref>) that can be achieved by both the PS and the WS methods.
 We show that these limits are in fact optimal; 
 any theoretical method is subject to the same bounds on the estimation error.

 We first derive the theoretical lower bound on the SQL.
 We consider the Fourier transform of the function φ(x):
 
    φ_k = ∫_0^L dx/L e^-2π ikx/Lφ(x).

 
 On the wavenumber basis, the constraint (<ref>) corresponds to the suppression of high-wavenumber components: φ_k = o(k^-q).
 In particular, the special constraint (<ref>) is equivalent to
 
  ∑_k=1^∞ k^2φ_k^2 ≤M^2/8π^2,
  
 
 which can be seen from Perseval's equality.
 
 A generalization of this argument leads to a sufficient condition on the constraint (<ref>):







 
    ∑_k=1^∞ k^2qφ_k^2 ≤M^2/2c_0^2,

 where c_0 = 2π^q max_0≤ x≤π [x^-qsin x] is a q-dependent constant.

 To utilize the known results in the discrete parameter estimation <cit.>, we focus on the functions with only some low-wavenumber components.
 Using a K-dimensional vector 5u=(u_1,…,u_K), we parametrize the function φ as
 
    φ_5u(x) = ∑_k=1^K√(2)u_ksin (2π ikx/L).

 Such function φ_5u meets the constraint if 5u≤ρ is satisfied for ρ = c_0^-1MK^-q.
 
 
 Since ρ=o(1) and √(2)sin (2π ikx/L) forms an orthonormal basis, the MSPE δ^2 of the function φ can be bounded by MSE δ_5u^2 of the vector 5u. 

 Hence, instead of the function estimation, we may consider the vector estimation in which case the error can be evaluated by the quantum Cramér-Rao bound (QCRB) <cit.>.
 For an unbiased estimation, the QCRB is given as
 
    δ_5u≥δ_1UUB := K( [J(5u)])^-1/2,

 where [J(5u)] is the Fisher information matrix defined for the output probe state |ψ_5u⟩ as
 
    [J(5u)]_jk = 4⟨[|]∂∂ u_jψ_5u[1-|ψ_5u⟩⟨ψ_5u|] |[⟩]∂∂ u_kψ_5u.

 Since the Fisher information is bounded from above by the SQL <cit.>:
 [J(5u)]_jj≤ 8N for each 1≤ j≤ K, we obtain a uniform, unbiased bound
 
  δ_1UUB = (K/8N)^1/2.
 

 Although such uniform bound is not applicable to the biased estimation,
 there exists the worst-case biased bound δ_1WBB <cit.>.
 Within the region 5u≤ρ, one can find a vector 5u satisfying
 
    δ_5u≥δ_1WBB,   δ_1WBB^-1 := δ_1UUB^-1 + ρ^-1.

 Since this biased version of QCRB holds for any integer K≥ 1, we choose K that gives the maximal bound δ_1WBB.
 This is satisfied when ρ and δ_1UUB are comparable to each other — this case holds when we set K=O((M^2 N)^1/(2q+1)). 
 Hence the SQL is given as
 
    δ≥ c_1 (M^1/qN^-1)^q/(2q+1),

 where the constant factor c_1 is lower-bounded by q/2(2q+1)^2 [c_0^-1q]^2q/(2q+1).
 

  
Optimality of the Heisenberg limit.
 We consider the case in which entanglement between at most n_1p (≥ 1) particles is allowed.
 It is known that the quantum information of a probe state is maximal when their wavefunction is completely symmetric <cit.>.
 With completely symmetric probe states, the problem becomes equivalent to the estimation of an effective phase n_1pφ(x) with n_1c=n_1p^-1N separate particles;
 the probe states in (<ref>) serve as an example for the WS method.














 Since the function of interest φ is replaced by its effective one n_1pφ, the MSPE δ^2 and the normalization constant M are replaced by n_1p^2δ^2 and n_1p M, respectively.  This argument leads to a generalized limit:
 
    n_1pδ   ≥ c_1 [(n_1p M)^1/q(n_1p^-1N)^-1]^q/(2q+1)
       = c_1 (M n_1p^q+1 N^-q)^1/(2q+1).

 To restore the original function φ(x) from the estimate of n_1pφ(x), we need to resolve the phase ambiguity by 2π/n_1p.
 For this purpose, the left-hand side of (<ref>) should not exceed π, giving
 
    n_1p≤ [(π c_1^-1)^(2q+1)/2M^-1N^q]^1/(q+1).

 With the maximal n_1p substituted in (<ref>), we obtain the Heisenberg limit
 
    δ≥ c_2 (M^1/qN^-1)^q/(q+1),

 where the constant c_2 is at least (π^-1c_1)^(2q+1)/(q+1).
 

  
Extension to smoother functions.
 
 The degree of smoothness can further be extended into q>1, where the target function is known to be more than just differentiable.
 For an integer m and 0<σ≤ 1 satisfying q=m+σ, the constraint for smoother functions is given as
 
    sup_0<ϵ<a∫_0^L dx/L[]φ^(m)(x+ϵ)-φ^(m)(x)/ϵ^σ^2
       ≤M^2/L^2q.


 Our results in the PS method and the optimality are also valid for q>1, thus leaving the quantum limits (<ref>) and (<ref>) unchanged.
 On the other hand, the straightforward extension of the WS method into q>1 does not work.
 Therefore, the asymptotic equivalence between the PS method and the WS method can be obtained only for 0<q≤ 1.
 See the Supplemental Material for more details <cit.>.


 

  
Comparison with Gaussian signal estimation.
 The error bounds we have obtained here are related to that of the Gaussian signal estimation <cit.>, in which the time-dependent phase φ_t is subject to a Gaussian process with the power spectrum I(ω)∼ω^-p.
 The estimation error of an instantaneous phase φ_t=0 is O(2N^-(p-1)/2p) for a coherent state and O(2N^-(p-1)/(p+1)) for a squeezed state, where 2N is the photon flux <cit.>.
 This can exactly be mapped into the SQL O(N^-q/(2q+1)) and the Heisenberg limit O(N^-q/(q+1)) in our study by setting p=2q+1.
 In fact, almost all sample functions of the Gaussian process φ_t satisfy Eq. (<ref>) by taking the large time span <cit.>.






 
 Therefore, the estimation error of the Gaussian process is subject to the quantum metrology of function estimation.
 We suspect that this fact is related to the minmax theorem <cit.>, which explains the consistency between the Bayesian and non-Bayesian estimation methods, though no clear connection is established yet.
 
 

  
Conclusion and outlook.
 In this Letter, we have established the fundamental limits on function estimation subject to a bounded qth-order differentiability.
 The estimation error is bounded from below by O( N^-q/(2q+1) ) in the standard quantum limit and O( N^-q/(q+1) ) in the Heisenberg limit.
 These results can be reduced to the previous studies on the signal estimation in quantum optics <cit.> when the target function is an infinitely extended stochastic process.
 We have presented two theoretical methods of the functional quantum metrology, both of which saturate the fundamental limits for 0<q≤ 1.
 
 
 This is a fundamental result for the effecient detection of functional structures — continuous signals and images, for example — which is a common target of estimation today.
 In fact, our results set theoretical bounds on various types of analysis relying on the function structure, such as model prediction or feature extraction <cit.>.
 On one hand, these bounds indicate the turning point where quantum methods outperform classical methods on functional data, with the scaling laws different from those obtained from parameter estimation.
 On the other hand, our result shows the optimal strategies for the quantum estimation of functions, such as an appropriate choice of temporal/spatial resolution or the size of entanglement.
 We note that choice of resolution is crucial in the real application <cit.>, and what is more in the quantum case, we have seen that larger entanglement does not necessarily mean better accuracy.

 The framework presented here also enables further quantum information-theoretic analysis on functions, such as a the quantum version of the Nyquist–Shannon sampling theorem which concerns the exact equivalence between the position- and wavenumber-states in the signal detection, including the O(1) prefactor that has remained undetermined. 


 
 We gratefully acknowledge D. W. Berry for helpful discussions on the quantum signal estimation.  Special thanks are due to Y. Ashida and R. Hamazaki for critical advice. 
 This work was supported by a Grand-in-Aid for Scientific Research on innovative Areas “Topological Materials Science” (KAKENHI Grant No. JP15H05855) from MEXT of Japan.
 N. K. was supported by the Leading Graduate Schools “ALPS.”
  

