Quantum Multi-Level Estimation of Functionals of Discrete Distributions

In this paper, to estimate functionals of an unknown discrete distribution ${(p_i)}_{i=1}^n$, we employ the purified quantum query access model, where one is given a unitary oracle, ${𝒪}_p$, for the distribution $p$ such that

with each $|{\varphi }_i⟩$ an arbitrary normalized pure state. This quantum query model is the now standard input model for quantum property testing, commonly used in quantum computational complexity [45] and quantum algorithms [19].

To establish our quantum estimators, we provide a quantum multi-level estimation framework for general functionals of the form ${\sum }_{i=1}^{n}f(p_i)$. While building on the partitioning strategy of [36], our framework significantly reduces the circuit implementation cost. We introduce the function $g(x)≔f(x)/x$ to rewrite the target functional as ${\sum }_{i=1}^n{p}_{i}g(p_i)$.

As a direct application of Theorem 1, we obtain a quantum estimator for Shannon entropy with query complexity $O(\sqrt{n}{\log }^{4.5}(n/\epsilon )/\epsilon )$, reproducing the result of [36] with explicit polylogarithmic factors (see Corollary 14).

Our specific main applications of Theorem 1 are a series of quantum estimators for the $q$-Tsallis entropy [39]

for all $q\in (0,1)\cup (1,+\mathrm{\infty })$, where ${\mathrm{F}}_q\left(p\right)={\sum }_{i=1}^n{p}_i^q$ is the key component to estimate. Classical estimators for the $q$-Tsallis entropy have been investigated in [5, 24, 25]. We introduce our quantum estimators below and compare them with prior results in Table 1.

To the best of our knowledge, Theorem 2 and Theorem 3 present the first near-optimal quantum estimators for the (parameterized) $q$-entropy of probability distributions for non-integer $q$.³³3Near-optimal quantum estimators for the $q$-entropy of probability distributions are known for integer $q$ in certain cases, e.g., the $2$-Rényi entropy ($q=2$) [26], the $q$-Tsallis entropy for integer $q\ge 2$ [40], and the Shannon entropy ($q=1$) [36]. It is known that the sample complexity of estimating the $q$-Tsallis entropy is $\mathrm{\Theta }(1/{\epsilon }^2)$ for $q\ge 1.5$ [25] and $\tilde{\mathrm{\Theta }}(1/{\epsilon }^{\frac{1}{q-1}})$ for [24]. In comparison, Theorem 2 and Theorem 3 achieve a strict quantum advantage over the corresponding classical estimators. Moreover, our results also improve the previous best query complexity $O(1/{\epsilon }^{1+\frac{1}{q-1}})$ for $q>1.5$ due to [28].⁴⁴4In [28], they considered the estimation of the $q$-Tsallis entropy of quantum states, which is a more general task. Even though, their query complexity remains unchanged when the quantum states degenerate to probability distributions.

Theorem 4 demonstrates a quantum speedup over the near-optimal classical estimator [24] for $q$-Tsallis entropy with sample complexity $\tilde{\mathrm{\Theta }}\left(\frac{n^{1/q}}{{\epsilon }^{1/q}}+\frac{n^{2-2q}}{{\epsilon }^2}\right)$ for . We also present a quantum query lower bound $\tilde{\mathrm{\Omega }}(n^{1/2q})$ (see Theorem 19 for details), suggesting that there is only room for sub-quadratic speedup on the distribution size $n$. We note that when $q=1-o(1)$, Theorem 4 achieves optimal dependence on $n$ up to a quasi-polynomial factor.

Moreover, given that an estimate of $q$-Tsallis entropy is always an estimate of $q$-Rényi entropy to the same precision when , Theorem 4 also implies an improved quantum query complexity of estimating the $\alpha$-Rényi entropy in certain regimes, combined with the results in [44, Theorem 1].

Corollary 5 improves the quantum query complexity $\tilde{O}\left(\frac{n^{\frac{1}{2\alpha }}}{{\epsilon }^{\frac{1}{2\alpha }+1}}\right)$ in [44] when $\epsilon \le O(n^{\frac{1-\alpha }{1-2\alpha }})$ and .

For more details of the theorems, propositions and their proofs, we refer to the full version [14].

Our paper leaves several natural open questions for future investigation:

• $\mathrm{■}$

  The bounds we obtain for Tsallis entropy estimation in the regime are not optimal. Is it possible to derive improved upper or lower bounds on the corresponding quantum query complexity? For example, can we show better lower bounds using other lower bound techniques such as the polynomial method [6], adversary method [2], and compressed oracle method [46]? More broadly, can the idea of multi-level estimation be applied to design optimal quantum algorithms for ${\mathrm{\ell }}_{\alpha }$ closeness testing or for estimating other notions of relative entropy between discrete distributions?

• $\mathrm{■}$

  There are results on the complexity of entropy estimation of quantum states [1, 19, 21, 37, 44, 41, 43, 28, 15]. The technique of [19] can be used to reduce the estimation of quantum entropies to the estimation of entropies of discrete probability distributions [44]. It remains an open question whether our quantum multi-level estimation technique can be applied in this setting, or adapted appropriately, to obtain optimal algorithms for estimating quantum entropies. A potential difficulty is that the eigenbasis of the quantum state is unknown.

• $\mathrm{■}$

  Our multi-level amplitude estimation shares several conceptual similarities with variable-time amplitude amplification and estimation algorithms [3, 12, 4, 29]. Are there deeper or more formal connections between these approaches? Moreover, given the broad range of applications of variable-time amplitude amplification, it is natural to ask whether our multi-level amplitude estimation framework can be extended to, or applied in, more general settings.

Unless otherwise stated, all operator approximations are measured in the spectral norm $\parallel \cdot \parallel$, and state approximations are measured in the Euclidean norm $\parallel \cdot {\parallel }_2$. We use $|\mathrm{𝟎}⟩$ to denote the all-zero state ${|0⟩}^{\otimes k}$ on a quantum register, where the number of qubits $k$ is clear from the context.

In this subsection, we review the basic quantum subroutines that serve as building blocks for our algorithms.

To estimate $a$ to within a given additive error $\epsilon$, we propose the following two-stage algorithm that first decides whether $a\le \epsilon$ and then determines $t$ accordingly.

A projected unitary encoding $(\tilde{\mathrm{\Pi }},\mathrm{\Pi },U)$ [20] of a matrix $A\in {ℂ}^{d\times d}$ is a unitary $U$ and two orthogonal projections $\mathrm{\Pi },\tilde{\mathrm{\Pi }}$ such that $A=\tilde{\mathrm{\Pi }}U\mathrm{\Pi }$. Let

be the singular value decomposition (SVD) of $A$. Given a scalar function $f$, the singular value transformation of $A$ is defined by applying $f$ to its singular values.

Given a projected unitary encoding of $A$, we can implement the singular value transformation of $A$ corresponding to a bounded polynomial with definite parity.

We recall the following lemmas about polynomial approximations of power functions.

Let $U$ be a projected unitary encoding of a matrix $A$ with SVD $A={\sum }_{i=1}^d{\sigma }_i|{\tilde{\psi }}_i⟩⟨{\psi }_i|$. A standard method for singular value discrimination used in [44, 36] utilizes QSVT with a polynomial $P$ that approximates a threshold function. Specifically, to distinguish whether a singular value ${\sigma }_i$ lies below $\gamma /\rho$ or above $\gamma$ for a threshold $\gamma \in (0,1)$, consider an even polynomial $P$ satisfying $P(x)\approx 0$ on $[0,\gamma /\rho ]$ and $P(x)\approx 1$ on $[\gamma ,1]$. Applying the associated QSVT unitary to the singular vector $|{\psi }_i⟩$ produces the state

where $|{\psi }_i^{⟂}⟩$ denotes a garbage state orthogonal to the subspace of interest. To partition ${\sigma }_i$ into finer intervals, multiple discrimination steps are necessary. Note that if , the produced state contains a non-negligible component of the garbage state $|{\psi }_i^{⟂}⟩$. The main drawback is that this garbage state cannot be reused for subsequent tests, since QSVT leaves the orthogonal subspace undefined. This necessitates preparing a fresh copy of $|{\psi }_i⟩$ in an additional ancilla register for every new discrimination step, incurring a substantial space overhead. To overcome this and discriminate singular values using constant ancilla qubits, we apply the gapped phase estimation with branch marking technique proposed in [29]. Since the algorithm only applies to Hermitian matrices, we convert the singular value problem into an eigenvalue problem using a Hermitian dilation. Consider the unitary $U_H=|0⟩⟨1|\otimes U+|1⟩⟨0|\otimes U^{†}$ acting on an extended system. With respect to the projector ${\mathrm{\Pi }}_H=|0⟩⟨0|\otimes \tilde{\mathrm{\Pi }}+|1⟩⟨1|\otimes \mathrm{\Pi }$, the operator $U_H$ encodes the Hermitian matrix

satisfying ${\mathrm{\Pi }}_H{U}_H{\mathrm{\Pi }}_H=H$. The spectral decomposition of $H$ takes the form $H={\sum }_{i,\nu }\nu {\sigma }_i|{\varphi }_{i,\nu }⟩⟨{\varphi }_{i,\nu }|$, where the eigenvectors are given by $|{\varphi }_{i,\nu }⟩=\frac{1}{\sqrt{2}}(|0⟩|{\psi }_i⟩+\nu |1⟩|{\tilde{\psi }}_i⟩)$ for $\nu \in \{+1,-1\}$.

The pseudocode of the algorithm is presented in Algorithm 2. For each $j\in [m]$, we isolate singular values in the $j$-th loop via two non-destructive singular value discriminators $D_{j-1}$ and $D_j$.

In the following theorem, we show that this process effectively partitions the $i$-th probability of the distribution into $m$ groups. Specifically, each $\sqrt{p_i}$ is assigned a probability weight ${|{\beta }_{i,j}|}^2$ for the $j$-th group, satisfying two key properties. Suppose $\sqrt{p_i}\in {\mathrm{\Delta }}_{j_i}$, and then we have

1. 1.

   (Localization of ${\beta }_{i,j}$) ${|{\beta }_{i,j}|}^2\approx 0$ for any $j\notin \{j_i-1,j_i\}$,

2. 2.

   (Completeness of ${\beta }_{i,j}$) ${|{\beta }_{i,j_i-1}|}^2+{|{\beta }_{i,j_i}|}^2\approx 1$.

We then apply a QSVT with polynomial $P_j$ to the isolated singular values, followed by amplitude estimation to estimate the squared norm of the output state. Taking the sum over all groups yields an estimate of:

Due to the localization of ${\beta }_{i,j}$, the polynomial $P_j^2(\sqrt{p_i}/2)$ only needs to approximate $g(p_i)$ for $\sqrt{p_i}\in {\mathrm{\Delta }}_j\cup {\mathrm{\Delta }}_{j+1}$. The completeness property then ensures that Equation (3) approximates the target functional ${\sum }_{i=1}^n{p}_{i}g(p_i)$.

A few remarks on the technical conditions are in order to clarify their meaning:

• $\mathrm{■}$

  Condition 1 (QSVT validity): This ensures that the polynomials are bounded and valid for implementation via QSVT.

• $\mathrm{■}$

  Condition 2 (Local approximation): This guarantees that the scaled polynomial approximates the target function with high accuracy within the $j$-th interval and its transition boundaries.

• $\mathrm{■}$

  Condition 3 and Condition 4 (Well-bounded tail): These conditions handle the residual interval $[0,{\phi }_m]$ closest to the singularity at the origin. In this region, we do not require an accurate polynomial approximation. Instead, Condition 4 ensures the true functional contribution from these small probabilities is sufficiently small, while Condition 3 restricts the magnitude of the polynomial estimate. Together, they guarantee via the triangle inequality that the total estimation error from this tail region is well-bounded.