Beyond the swap test: optimal estimation of quantum state overlap

TL;DR

This paper develops optimal collective measurement strategies for estimating the overlap between two unknown quantum states, outperforming traditional swap tests, especially for small overlaps and in noisy conditions.

Contribution

It introduces the optimal measurement framework for quantum state overlap estimation, surpassing swap tests in accuracy and robustness, with detailed analysis of performance and noise resilience.

Findings

Optimal measurements reduce mean square error compared to swap tests.

Performance advantage increases with system dimension and for small overlaps.

Optimal measurement is less invasive and more noise-robust than swap tests.

Abstract

We study the estimation of the overlap between two unknown pure quantum states of a finite dimensional system, given $M$ and $N$ copies of each type. This is a fundamental primitive in quantum information processing that is commonly accomplished from the outcomes of $N$ swap-tests, a joint measurement on one copy of each type whose outcome probability is a linear function of the squared overlap. We show that a more precise estimate can be obtained by allowing for general collective measurements on all copies. We derive the statistics of the optimal measurement and compute the optimal mean square error in the asymptotic pointwise and finite Bayesian estimation settings. Besides, we consider two strategies relying on the estimation of one or both the states, and show that, although they are suboptimal, they outperform the swap test. In particular, the swap test is extremely inefficient for small values of the overlap, which become exponentially more likely as the dimension increases. Finally, we show that the optimal measurement is less invasive than the swap test and study the robustness to depolarizing noise for qubit states.

Figures

Tables

Table 1. Table 1: Local MSE and global MSE attainable via the optimal, EP and EE strategies in two asymptotic limits. In all the cases the global MSEs coincide with the corresponding average local MSE values, apart from asymptotically vanishing corrections.

Local est.	${𝐯}_{\mathrm{op}}\left(𝐜\right)$	${𝐯}_{\mathrm{ep}}\left(𝐜\right)$	${𝐯}_{\mathrm{ee}}\left(𝐜\right)$
$M=N\to \infty$	$\frac{4c\left(1-c\right)}{M+N}$	$\frac{3}{2}{v}_{\mathrm{op}}\left(c\right)$	$2v_{\mathrm{op}}\left(c\right)$
$M\to \infty$	$\frac{c\left(1-c\right)}{N}$	$v_{\mathrm{op}}\left(c\right)$	$2v_{\mathrm{op}}\left(c\right)$
Bayesian est.	${𝐯}_{\mathrm{op}}$	${𝐯}_{\mathrm{ep}}$	${𝐯}_{\mathrm{ee}}$
$M=N\to \infty$	$\frac{4\left(d-1\right)}{d\left(d+1\right)\left(M+N\right)}$	$\frac{3}{2}{v}_{\mathrm{op}}$	$2v_{\mathrm{op}}$
$M\to \infty$	$\frac{\left(d-1\right)}{d\left(d+1\right)\left(d+N\right)}$	$v_{\mathrm{op}}$	$\frac{d+2N}{2+N}{v}_{\mathrm{op}}$

Equations

$$v=\sum _k\int dUdV(c(k)-c{)}^2\mathrm{Tr}\left[E_k|\Psi ⟩⟨\Psi |\right].$$

$$\rho (c)=\int \mathrm{d}U{U}^{\otimes (N+M)}|{\Psi }_0⟩⟨{\Psi }_0|U^{†\otimes (N+M)}.$$

$$v=\sum _k\int dcp(c)(c(k)-c{)}^2\mathrm{Tr}\left[E_k\rho (c)\right]$$

$$p(c)=\int dU\delta (c-|⟨0|U|0⟩{|}^2)=(d-1)(1-c{)}^{d-2}.$$

$$v(c)=\sum _k(c(k)-c{)}^2\mathrm{Tr}\left[E_k\rho (c)\right],$$

$$\rho(c)=\sum_{J=J_{\min}}^{J_{\max}}p(J|c)\frac{{\mbox{$1{\bf l}$}}_{J}}{\chi_{J}}\otimes\left|\sigma\right\rangle\left\langle\sigma\right|_{J},$$

$$p(J|c)=\frac{(2J+1)N!M!(1-c{)}^M{P}_{J+J_{\min }}^{(0,-2J_{\min })}\left(\frac{1+c}{1-c}\right)}{(J_{\max }-J)!(J_{\max }+1+J)!},$$

$$c_{\mathrm{op}}^{\mathrm{loc}}(J)={\left(\frac{2J}{M+N}\right)}^2,v_{\mathrm{op}}(c)=\frac{4c(1-c)}{M+N},$$

$$c_{\mathrm{op}}^{\mathrm{bay}}(J)$$

$$v_{\mathrm{op}}$$

$$F(c)={\int }_{U\in \mathrm{SU}(d)}dU\sum _k{\left|⟨{\Psi }_U|E_k|{\Psi }_U⟩\right|}^2,$$

$$F_{\mathrm{op}}(c)=\sum _{J=J_{\min }}^{J_{\max }}p(J|c{)}^2,F_{\mathrm{sw}}(c)={\left(\frac{1+c^2}{2}\right)}^N,$$

$$u_{N+M}(U)|i_1⟩\otimes |i_2⟩\otimes ...\otimes |i_{N+M}⟩$$

$$=U|i_1⟩\otimes U|i_2⟩\otimes ...\otimes U|i_{M+N}⟩,\forall U\in \mathrm{SU}(d)$$

$$s_{N+M}(\sigma )|i_1⟩\otimes |i_3⟩\otimes ...\otimes |i_{N+M}⟩$$

$${\mathcal{H}}_d^{\otimes (N+M)}\cong \underset{Y}{⨁}{\mathcal{U}}^{(Y)}(\mathrm{SU}(d))\otimes {\mathcal{U}}^{(Y)}(S_{N+M}),$$

$$u_{N+M}(U)$$

$$s(\sigma )$$

$${\mathcal{G}}_{\mathrm{SU}(d)}[\rho ]=\int \mathrm{d}U{U}^{\otimes (N+M)}\rho U^{†\otimes (N+M)},$$

$${\mathcal{G}}_{S_{N+M}}[\rho ]=\frac{1}{(N+M)!}\sum _{\sigma }{s}_{N+M}(\sigma )\rho s_{N+M}^{†}(\sigma ).$$

$${\mathcal{G}}_{\mathrm{SU}(d)}[\rho ]$$

$${\mathcal{G}}_{S_{N+M}}[\rho ]$$

$${\mathrm{Sym}}_{N,M}^{(d)}:={\mathrm{Sym}}_N^{(d)}\otimes {\mathrm{Sym}}_M^{(d)}\subset {\mathcal{H}}_d^{\otimes N+M}.$$

$${\mathrm{Sym}}_{N,M}^{(d)}={\oplus }_J\left({\mathcal{U}}^{(J)}(\mathrm{SU}(d))\otimes {\mathcal{K}}^J\right),$$

$${\mathrm{Sym}}_{N,M}^{(2)}:={\oplus }_J\left({\mathcal{U}}^{(J)}(\mathrm{SU}(2))\otimes {\mathcal{K}}_J^{(2)}\right)\subseteq {\mathrm{Sym}}_{N,M}^{(d)}.$$

$${\mathrm{Sym}}_{N,M}^{\mathcal{E}}:={\oplus }_J\left({\mathcal{U}}_{\mathcal{E}}^{(J)}(\mathrm{SU}(2))\otimes {\mathcal{K}}_{\mathcal{E}}^J\right)\subseteq {\mathrm{Sym}}_{N,M}^{(d)}.$$

$${\mathcal{U}}_{\mathcal{E}}^{(J)}(\mathrm{SU}(2))\otimes {\mathcal{K}}_{\mathcal{E}}^J\subseteq {\mathcal{U}}^{(J)}(\mathrm{SU}(d))\otimes {\mathcal{K}}^J$$

$$|\varphi ⟩=V(\theta )|\psi ⟩\equiv \sqrt{c}|\psi ⟩+\sqrt{1-c}|{\psi }_{\perp }⟩$$

$${\mathcal{G}}_{\mathrm{SU}(d)}\left[{\left(|\psi ⟩⟨\psi |\right)}^{\otimes N}\otimes {\left(|\varphi ⟩⟨\varphi |\right)}^{\otimes M}\right]\equiv {\int }_{\mathrm{SU}(d)}\mathrm{d}U{\left(U|\psi ⟩⟨\psi |U^{†}\right)}^{\otimes N}\otimes {\left(U|\varphi ⟩⟨\varphi |U^{†}\right)}^{\otimes M}$$

$${\mathcal{G}}_{\mathrm{SU}(d)}\left[{\left(|\psi ⟩⟨\psi |\right)}^{\otimes N}\otimes {\left(|\varphi ⟩⟨\varphi |\right)}^{\otimes M}\right]={\int }_{\mathrm{SU}(d)}\mathrm{d}U{U}^{\otimes (N+M)}\left[{\left(|\psi ⟩⟨\psi |\right)}^{\otimes N}\otimes {\left(|\varphi ⟩⟨\varphi |\right)}^{\otimes M}\right]{U^{†}}^{\otimes (N+M)}$$

$$={\int }_{\mathrm{SU}(d)}\mathrm{d}U{U}^{\otimes (N+M)}\left({\int }_{\mathrm{SU}(\mathcal{E})}\mathrm{d}V{V}^{\otimes (N+M)}\left[{\left(|\psi ⟩⟨\psi |\right)}^{\otimes N}\otimes {\left(|\varphi ⟩⟨\varphi |\right)}^{\otimes M}\right]{V^{†}}^{\otimes (N+M)}\right){U^{†}}^{\otimes (N+M)}$$

$$={\mathcal{G}}_{\mathrm{SU}(d)}\left[{\mathcal{G}}_{\mathrm{SU}(\mathcal{E})}\left[{\left(|\psi ⟩⟨\psi |\right)}^{\otimes N}\otimes {\left(|\varphi ⟩⟨\varphi |\right)}^{\otimes M}\right]\right]$$

$${|\psi ⟩}^{\otimes N}\otimes {|\varphi ⟩}^{\otimes M}$$

$$=\sum _{J=J_{\min }}^{J_{\max }}\sum _{k=-\frac{M}{2}}^{\frac{M}{2}}{C}_{\frac{N}{2},\frac{N}{2};\frac{M}{2},k}^{J,\frac{N}{2}+k}{d}_{k,\frac{M}{2}}^{\left(\frac{M}{2}\right)}(\theta )|J,\frac{N}{2}+k⟩,$$

Full text

Introduction.— The overlap between two unknown quantum states is an archetypical instance of quantum relative information OverlapLSB ; Unspeakable ; ReviewFramesInfo ; ResourceAsymmetry ; AsymmetryBasic ; AsymmetryModes ; genSchurEst ; ChiribellaMo and the estimation of the overlap is a basic primitive in quantum information processing, with applications ranging from quantum fingerprinting SwapTest ; QFinger1 ; QFinger2 , entanglement estimation FuncOverlap ; EntDet1 ; EntDet2 ; EntDet3 , communication without a shared reference frame OverlapBRS ; OverlapBIM ; OverlapLSB ; OverlapGI to quantum machine learning HHL ; LloydAlgo ; QuantumSVM ; ProgrammableProjective ; LearningSwap ; GaussianProc ; NNAlgorithm ; IBMSupervised ; Nana ; variationalnear . Recently, with the advent of quantum machine learning protocols HHL ; LloydAlgo ; QuantumSVM ; GaussianProc ; NNAlgorithm ; IBMSupervised ; Nana , overlap estimation (OvE) has attracted renewed interest as a fundamental primitive and its efficient implementation and generalization on near-term quantum computers are subjects of current research ProgrammableProjective ; LearningSwap . In most applications, OvE is carried out through the swap test (SWT) SwapTest ; LearningSwap ; ProgrammableProjective : given two systems in the state $\left|\psi\right\rangle\left|\phi\right\rangle$, the probability of projecting it on its symmetric or antisymmetric part is determined by the overlap between $\left|\psi\right\rangle$ and $\left|\phi\right\rangle$. By repeating this measurement on several pairs of copies one can obtain a good estimate of this probability, and hence of the overlap. It is then natural to ask whether, for the same number of copies, one could reach a larger accuracy via a collective strategy that extracts the relevant information using a joint and less-destructive measurement. In this article we answer in the positive, evaluating the ultimate precision attainable in the OvE of two pure quantum states, given a number of copies of each and assuming no prior knowledge about them.

The task we consider is as follows: given $N$ and $M\geq N$ copies of two unknown pure states $\left|\psi\right\rangle,\left|\phi\right\rangle$ of a $d$-dimensional quantum system we are requested to provide an estimate of their (squared) overlap $\lvert\left\langle\psi\lvert\phi\right\rangle\rvert^{2}$ which is fixed, but unknown to us. The task is carried out by a machine that performs a measurement on the state $\left|\Psi\right\rangle=\left|\psi\right\rangle^{\otimes N}\otimes\left|\phi\right\rangle^{\otimes M}$ of $M+N$ qudits and produces an estimate with maximum precision, as quantified by the mean square error (MSE). Furthermore we consider the case of unlabeled states, i.e., when the machine receives $U_{\sigma}\left|\Psi\right\rangle$, with $U_{\sigma}$ an unknown permutation of the qudits. Note that in this case OvE constitutes in itself an instance of unsupervised quantum-classical learning problem, in a setting similar to SentisUnsuper .

The measurements optimizing the average information gain OverlapBRS and the average error OverlapBIM ; OverlapLSB have been derived for the case of qubits, with only numerical solutions OverlapGI for higher dimensions.

Here we tackle OvE in full generality, characterizing the optimal estimation within both local (pointwise) and global (Bayesian) approaches. For local estimation, we provide an asymptotically achievable lower bound using the quantum Fisher information (QFI) HelstromBOOK ; FisherParis , whereas for Bayesian estimation HelstromBOOK ; Personick we provide an exact solution, generalizing the results of OverlapGI . We find that the optimal local strategy is also Bayesian-optimal asymptotically and that it performs identically in the labeled and unlabeled scenarios. We compare our results with the SWT and with two LOCC strategies based on estimating either one or both $\left|\psi\right\rangle$ and $\left|\phi\right\rangle$, see Fig. 1. Such strategies are useful in distributed scenarios where copies of $\left|\psi\right\rangle$ and $\left|\phi\right\rangle$ are produced in different and distant laboratories. We show that in the limit of large $M+N$ and $|M-N|$ constant the optimal strategy displays a finite asymptotic gap with respect to all the others. Moreover, we show that the optimal measurement is less invasive than the SWT and robust against single-qubit noise.

Assessing the machine’s performance.— The states $\left|\psi\right\rangle=U\left|0\right\rangle$ and $\left|\phi\right\rangle=V\left|0\right\rangle$ are drawn uniformly at random, i.e., with Haar-distributed $U,V\in SU(d)$. Upon performing a measurement $\{E_{k}\}$ on $\left|\Psi\right\rangle$ with outcome $k$, the machine outputs an estimate $c(k)$ of the overlap $c=\lvert\left\langle\psi\lvert\phi\right\rangle\rvert^{2}$, with squared error $(c(k)-c)^{2}$.

In the global approach, the machine’s performance is quantified by averaging the squared error over all possible states and outcomes:

[TABLE]

We refer to $v$ as global MSE. Writing $V=UW$ and using the Haar-measure invariance $dV=dW$ and $c=|\left\langle 0\right|W\left|0\right\rangle|^{2}$ the average mean square error can be written as: $v=\sum_{k}\int dW\,(c(k)-c)^{2}\mathrm{Tr}\left[E_{k}\rho(c)\right]$, where we have defined the effective state

[TABLE]

where \left|\Psi_{0}\right\rangle={\mbox{1{\bf l}}}\otimes W^{\otimes M}\left|0\right\rangle^{\otimes(N+M)}. In addition, we can write the above integral over $W$ as an integral over the overlap such that

[TABLE]

where the distribution over overlaps is given by OverlapStatistics ; SentisUnsuper

[TABLE]

From the above discussion we see that the average over both types of states, i.e. over $U$ and $V$, is equivalent to an average over overlaps with weight $p(c)$ and over different orientations, $U^{\otimes(N+M)}$. This is a direct consequence of the fact that if the states are completely unknown, then all pairs of states with equal overlap are equally probable and are related by a rigid unitary: $\lvert\left\langle\psi\lvert\phi\right\rangle\rvert^{2}=\lvert\left\langle\psi^{\prime}\lvert\phi\right\rangle^{\prime}\rvert^{2}$ if and only if it exists $U$ such that $\left|\psi^{\prime}\right\rangle=U\left|\psi\right\rangle$ and $\left|\phi^{\prime}\right\rangle=U\left|\phi\right\rangle$.

At variance with the global approach, where the overlap is a random variable, in the local approach the overlap is considered to be fixed. We then assess the performance of the machine by computing the average of the square error over all states with fixed overlap $c$ and over all outcomes:

[TABLE]

also in terms of the average state for a fixed overlap $\rho(c)$. We refer to $v(c)$ as local MSE. As shown in the supplemental material (SM), the integral in Eq. (2) can be performed using $SU(d)$ representation theory and $SU(2)$ Clebsch-Gordan coefficients, obtaining the block-diagonal form

[TABLE]

with $J_{\min}=\frac{\left\lvert M-N\right\rvert}{2}$, $J_{\max}=\frac{M+N}{2}$ and

[TABLE]

where $P_{n}^{(\alpha,\beta)}(x)$ is the nth-degree Jacobi polynomial. In the previous equations, for $d=2$, $J$ is the familiar total-angular-momentum label and {\mbox{1{\bf l}}}_{J}=\sum_{M=-J}^{J}\left|J,M\right\rangle\!\!\left\langle J,M\right| is the projector on the subspace of total angular momentum $J$, of dimension $\chi_{J}=2J+1$. In general, for $d>2$, {\mbox{1{\bf l}}}_{J} are projectors over the subspaces of dimension $\chi_{J}(d)$ hosting irreducible representations (irreps) of $\mathrm{SU}(d)$ arising from the tensor product of two completely symmetric representations of $M$ and $N$ qudits; these irreps are still indexed by an (half)-integer $J\in[J_{\min},J_{\max}]$. Finally, $\left|\sigma\right\rangle_{J}$ is a state representing the known labeling of the states and it belongs to the irrep-space of the permutation group, also labeled by $J$. Note that, in the unlabeled scenario, the average over qudit permutations acts only on $\left|\sigma\right\rangle_{J}$ for each $J$, depolarizing it to a projector on the whole irrep space. Importantly, note that all the information about the overlap is contained in the $J$-statistics $p(J|c)$, which is independent of dimension and labeling. In particular, the optimal measurement is given by the projectors $\Pi_{J}$ on the subspaces labeled by $J$, and it can be implemented via weak Schur sampling Schur3 . Indeed, for any POVM $\{E_{k}\}_{k}$, we can get the same outcome probability distribution if we use the POVM $\{\Pi_{J}E_{k}\Pi_{J}\}_{k,J}$ and then ignore the $J$ label. When this POVM is applied to $\rho(c)$ the outcome probabilities are $p_{k,J}:=p(k|J)p(J|c)$. The same outcome probabilities $p_{k,J}$ can be generated by applying directly the POVM $\{\Pi_{J}\}_{J}$ followed by classical post-processing. The latter can only increase the variance of the estimator, by convexity of the figure of merit: $\sum_{k}p(k|J)(c_{k,J}-c)^{2}\geq(\sum_{k}p(k|J)c_{k,J}-c)^{2}$, which follows from the Cauchy-Schwartz inequality. Therefore $\{\Pi_{J}\}_{J}$ is optimal both for local and global estimation and one can replace $\mathrm{Tr}\left[E_{J}\rho(c)\right]$ with $p(J|c)$ in Eqs. (1,5), effectively reducing our problem to one of classical estimation, i.e., optimizing the function $c(J)$.

Local estimation.—

The classical Cramer-Rao bound HolevoBook places a lower bound on the MSE of all local unbiased estimators $c(J)$ as $v(c)\geq H(c)^{-1}$, where $H(c)=\sum_{J}(\partial_{c}p(J|c))^{2}/p(J|c)$ is the Fisher information of the measurement statistics. In the limit $M+N\to\infty$ and $M-N\ll(M+N)\sqrt{c}$, we can use an approximation of the Jacobi polynomial given in Szego1939 to obtain the following asymptotically-unbiased estimator and its associated MSE:

[TABLE]

In the SM we show that $v_{\rm op}(c)$ coincides with $H(c)^{-1}$ to leading order in $\frac{1}{M+N}$ and hence the Cramer-Rao bound is achievable in this limit. If instead $M\rightarrow\infty$ and $N$ is finite, it is clear that $\left|\phi\right\rangle$ can be estimated perfectly and hence the optimal strategy is to project the copies of $\left|\psi\right\rangle$ in this known direction, with resulting $v_{\rm op}(c)=\frac{c(1-c)}{N}$.

Bayesian estimation.—

The optimal classical Bayesian (global) estimator is given HelstromBOOK by $c_{\rm op}^{\rm bay}(J)=\frac{\int dc\,c\,p(c)p(J|c)}{\int dc\,p(c)p(J|c)}$. Using graphical calculus techniques for the recoupling theory of Clebsch-Gordan coefficients Varshalovich1988 , as explained in the SM, we obtain the following optimal global estimator and corresponding MSE:

[TABLE]

We pause to highlight the following facts: i) when $d$ is fixed and the number of copies is large, the prior distribution of the states is little informative with respect to the information that can be obtained by the actual measurement; indeed we can see that when $M+N\rightarrow\infty$, $M-N$ constant, $c_{\rm op}^{\rm bay}(J)\approx c_{\rm op}^{\rm loc}(J)$ implying that the local optimal estimator is also a good Bayesian estimator and viceversa; ii) contrarily to the local estimation results, the global MSE of Eq. (10) is exact for all $M$, $N$ and depends on $d$ due to the prior, Eq. (4); iii) in particular, $v_{\rm op}$ decays as $d^{-2}$ if one of either $M$ or $N$ is kept finite, whereas it decays only as $d^{-1}$ when $M,N\gg 1$.

1-LOCC strategies.—

We now consider a family of intermediate strategies that employ 1-LOCC on $\left|\psi\right\rangle^{\otimes N}$ and $\left|\phi\right\rangle^{\otimes M}$. The estimate-and-project (EP) strategy consists in estimating $\left|\phi\right\rangle$ from its $M$ copies, then projecting each copy of $\left|\psi\right\rangle$ on this estimate and counting the fraction of successful projections. When $\left|\phi\right\rangle$ is known, projecting $\left|\psi\right\rangle$ on $\left|\phi\right\rangle$ is optimal HolevoBook . However, EP is not necessarily the optimal 1-LOCC strategy. The corresponding POVM elements can be written as $E_{V,k}^{\rm(ep)}=dVE_{V}^{(M)}\otimes V^{\otimes N}\Pi_{k}^{(N)}V^{{\dagger}\otimes N}$, where $E_{V}^{(M)}=\chi_{\frac{M}{2}}(d)(V\left|0\right\rangle\left\langle 0\right|V^{{\dagger}})^{\otimes M}$ is the optimal covariant measurement to estimate $\left|\phi\right\rangle$ Hayashi1997 ; HayashiGroupRep and $\Pi_{k}^{(N)}$ represents $k$ successful projections of the copies of $\left|\psi\right\rangle$ on the estimate of $\left|\phi\right\rangle$. The estimator is $c_{\rm ep}^{\rm loc}(k)=\frac{k}{N}$. The estimate-and-estimate (EE) strategy instead consists in estimating both $\left|\psi\right\rangle$ and $\left|\phi\right\rangle$ separately, then computing the overlap between the estimated states. The corresponding POVM elements can be written as $E_{V,W}^{\rm(ee)}=dVdWE_{V}^{(M)}\otimes E_{W}^{(N)}$, i.e., a product of two covariant measurements to estimate $\left|\phi\right\rangle$ and $\left|\psi\right\rangle$. We take as local estimator $c_{\rm ee}^{\rm loc}(V,W)=\left\lvert\left\langle 0\right|V^{{\dagger}}W\left|0\right\rangle\right\rvert^{2}$. In the SM we provide exact results for local and Bayesian estimation using EP and EE. Table (1) compares the performance of these strategies with the optimal one in two asymptotic limits. We find that, for both local and Bayesian estimation, the EE strategy is always worse than the optimal by a factor of $2$, whereas the EP strategy attains a MSE equal to the optimal in the limit $M\rightarrow\infty$, $N$ finite.

Performances comparison.— We now compare the strategies discussed so far with the traditional SWT SwapTest . Note that all these strategies except the optimal one require labeling of the states. The latter consists in projecting the state $\left|\psi\right\rangle\otimes\left|\phi\right\rangle$ on its triplet/singlet components, hence it coincides with the optimal measurement for $M=N=1$. As the SWT acts on couples of states, we restrict to the case $M=N$. The probability of a triplet projection $p(c)=\frac{1}{2}\left(1+c\right)$ and the ensuing statistics of $k$ successful projections out of $N$ trials is given by the binomial distribution. The optimal local MSE attainable by this test is well-known, $v_{\rm sw}(c)=\frac{1-c^{2}}{N},$ while for the optimal global MSE $v_{\rm sw}$ one can derive an exact expression for each value of $k$, then compute the sum numerically, as detailed in the SM. In the asymptotic limit of $M=N\gg d$ a good approximation is provided by averaging the optimal local MSE: $v_{\rm sw}\simeq\int dc\;p(c)v_{\rm sw}(c)=(d+2)(d-1)/(d(d+1)N)$, which is $\sim d$ times larger than $v_{\rm op}$.

In the same limit, we can compare the local MSE of all the strategies, see Fig. 2. First, we observe a gap between the optimal strategy, that attains the QFI, and all the other strategies. This means that, even with a large number of copies, the collective measurement on $\left|\psi\right\rangle^{\otimes N}\otimes\left|\phi\right\rangle^{\otimes M}$ has a clear advantage over a non-collective one. Second, we observe that the relative error $\frac{\sqrt{v(c)}}{c}$ for small $c$ scales as $\frac{1}{c\sqrt{N}}$ for the SWT and as $\frac{1}{\sqrt{cN}}$ for the other strategies, implying a quadratic improvement in $\frac{1}{\sqrt{c}}$ in the number copies needed to reach a fixed relative error, while the optimal measurement is still computationally efficient (see next section). This is particularly relevant since for large $d$ small overlaps are exponentially more likely, see Eq. 4. This phenomenon is also at the source of the so-called “barren plateau” problem googleplateau ; glassyplateau for quantum variational circuits, and other types of strategies have been proposed to address this issue plateaulearning ; quantumplateau ; unitaryoverlap .

We notice similar features for the global MSE, plotted in Fig. 3 as a function of $N$ for $M$ fixed and increasing $d$ (inset). We observe that the SWT is comparable with EE for $M\sim N$ and $d=2$, but with a small increase in dimension this feature disappears. Moreover, there is in general a gap between the EP and EE strategies, the former being closer to the optimal one.

Gate complexity.— The advantage in the precision of the optimal estimation comes with the tradeoff that the optimal measurement requires entangling operations over the whole system of $N+M$ qudits. The Schur transform  Schur1 ; Schur2 ; Schur3 is a way to perform the optimal measurement, and requires $O(\mathrm{poly}(N+M,\log d,\log\frac{1}{\epsilon}))$ qudit gates for precision $\epsilon$. The resulting algorithm is efficient, but still unfeasible without error correction. The SWT instead requires $N$ independent circuits of fixed depth, and may still be convenient for large overlaps or very noisy gates.

A mid-term solution is to divide input data in $R$ groups of $S$ copies of $\left|\phi\right\rangle$ and $\left|\psi\right\rangle$, such that $S$ is the largest integer for which the given architecture can perform the optimal measurement with high fidelity, repeat the measurement $R$ times and do classical post-processing. The performances of these intermediate protocols are between SWTs and optimal measurement. See the SM for a more detailed discussion of these issues.

Measurement invasiveness.— Another relevant figure of merit for applications is the fidelity between the post-measurement state and the initial one, averaged over the measurement outcomes. Both the optimal measurement and the SWT are projective measurements. We assume that the post-measurement states are given by the result of such projections and hence the average post-measurement fidelity can be written as

[TABLE]

with \{E_{k}\equiv{\mbox{1{\bf l}}}_{J}\} for the optimal measurement and \{E_{k}=\mathcal{G}_{{\cal S}_{N}}({\mbox{1{\bf l}}}_{2}^{S\,\otimes k}\otimes{\mbox{1{\bf l}}}_{2}^{A\,\otimes N-k})\} for the SWT, where {\mbox{1{\bf l}}}^{S/A}_{2} are the projectors on the singlet/triplet components of $\mathcal{H}_{2}^{\otimes 2}$. Then Eq. (11) is simply given by

[TABLE]

as shown in the SM. In Fig. 3 we plot these two quantities as a function of $c$, showing that the optimal measurement is less invasive than the SWT, especially for small overlap values.

Noise-robustness.— Finally, we consider how the optimal strategy changes when the states, which are expected to be pure, are affected by depolarizing noise acting independently on each qudit before reaching the measurement stage. Note that if the noisy channel is of a different kind, one can at least reach the optimal MSE for the depolarizing channel by performing a twirling operation, realizable by pre- and post-processing with random unitaries on each qudit plus classical forward communication. This operation is $\int dUU^{{\dagger}}{\cal N}(U\rho U^{{\dagger}})U=\Delta_{r}(\rho)$ for some $r$, where we have defined the depolarizing channel as \Delta_{r}=r\mathcal{I}+(1-r)\frac{{\mbox{1{\bf l}}}}{d}\mathrm{Tr} and $\mathcal{I}$ is the identity channel. After this operation the overall state of the system can now be written as $\Delta_{r_{0}}(\psi)^{\otimes N}\otimes\Delta_{r_{1}}(\phi)^{\otimes M}$.

In the SM we compute the optimal MSE in this case, restricting to $d=2$ for simplicity. In the limit $M,N\rightarrow\infty$ with $\frac{M}{N}$ finite, the global MSE at leading order is $v_{\rm op,mix}=\frac{1}{6Mr_{0}^{2}}+\frac{1}{6Nr_{1}^{2}}$, which agrees with the previously found limit of Eq. (10) for zero-noise, $r_{i}=1$. Hence the net effect of white noise is to rescale the MSE by a factor $r_{i}^{-2}$ for each state.

Conclusions.— In this article we have computed the ultimate precision attainable in estimating the overlap of two arbitrary pure quantum states, as a function of the dimension of their Hilbert space and their number of copies. We showed that the commonly-used SWT is highly inefficient for small values of the overlap and also on average over Haar-distributed random states. The optimal strategy is a collective measurement on all the copies and can be implemented efficiently using the Schur transform, although it remains experimentally challenging. A practical alternative is to do Schur sampling on subsets of the dataset, followed by classical post-processing. In addition, we proposed two intuitive strategies that estimate separately one or both states and showed that they also outperform the SWT. Finally, we showed that the optimal measurement is less invasive than the SWT and robust to white noise. The strategies we introduced provide several clear advantages over the SWT, and they could become a standard tool for various quantum technologies, as well as providing improvements in the runtime of quantum algorithms.

Note added after publication.— We thank one of the referees of TQC 2020 for pointing out BOW to us, where the authors find a minimum variance unbiased estimator of the Hilbert-Schmidt distance of two unknown mixed states, and compute its variance. Our results, valid for more general estimators, complement their analysis in the case of pure states.

Acknowledgments.— M. F. and V. G. acknowledge support from PRIN 2017 “Taming complexity with quantum strategies”. M. R., M. S. and J. C. acknowledge support from the Spanish MINECO, project FIS2016-80681-P with the support of AEI/FEDER funds; the Generalitat de Catalunya, project CIRIT 2017-SGR-1127. M. R. also acknowledges partial financial support by the Baidu-UAB collaborative project “Learning of Quantum Hidden Markov Models”. M. S. also acknowledges support from the Spanish MINECO project IJCI-2015-24643. M. F. thanks Matthias Christandl for helpful comments. M. S. and J. C. acknowledge useful discussion with Nana Liu.

Appendix A Schur-Weyl duality and the irreducible representations of $\mathrm{SU(d)}$ and $S_{N+M}$

We begin by first reviewing some key concepts and techniques in group representation theory that we will frequently make use of throughout this Supplemental material. The most important ingredient is Schur-Weyl duality Weyl . Consider the state space of $N+M$, $d$-dimensional systems, $\mathcal{H}_{d}^{\otimes(N+M)}$. This space carries the action of two different groups; the special unitary group of $d\times d$ complex matrices, $\mathrm{SU}(d)$, and the permutation group of $N+M$ objects, $S_{N+M}$. Specifically, the groups $\mathrm{SU}(d)$ and $S_{N+M}$ act on a basis $\{\left|i_{1}\right\rangle\otimes\left|i_{2}\right\rangle\otimes...\otimes\left|i_{N+M}\right\rangle\}_{i_{1},i_{2}....,i_{N+M}}$ of $\mathcal{H}_{d}^{\otimes(N+M)}$ via unitary representations $u_{N+M}:\mathrm{SU}(d)\to\mathrm{U}(\mathcal{H}_{d}^{\otimes(N+M)})$, and $s_{N+M}:S_{N+M}\to\mathrm{U}(\mathcal{H}_{d}^{\otimes(N+M)})$ as follows

[TABLE]

Observe that $[U^{\otimes(N+M)},s(\sigma)]=0,\;\forall U\in\mathrm{SU}(d),\,\mathrm{and}\,\forall\sigma\in S_{N+M}$. Schur-Weyl duality states that the total state space $\mathcal{H}_{d}^{\otimes(N+M)}$ can be decomposed as

[TABLE]

where $\mathcal{U}^{(Y)}(\mathrm{SU}(d))$ is the space of dimension $\chi_{Y}$ upon which the unitary irreducible representation (irrep) $u^{(Y)}$ of $\mathrm{SU}(d)$ acts and $\mathcal{U}^{(Y)}(S_{N+M})$ is the space of dimension $\omega_{Y}$ on which the irrep $s^{(Y)}$ of $S_{N+M}$ acts. The use of the congruence sign in Eq. (15) indicates that this block decomposition is accomplished by a unitary transformation; in the case considered here this unitary is the Schur transform Schur1 ; Schur2 .

The label $Y$ indexes the integer partitions of $N+M$ in at most $d$ parts written in descending order, pictorially represented by Young diagrams, where $N+M$ boxes are arranged into at most $d$ rows. For $d=2$ the block decomposition of Eq. (15) is the familiar decomposition of $N+M$ spin-$\nicefrac{{1}}{{2}}$ systems into the total angular momentum label $J$. The latter is related to Young diagrams of at most two rows where the number of boxes in any of the two rows is given by $(\frac{N+M}{2}+J,\frac{N+M}{2}-J)$.

With the help of Eq. (15) we can express the action of $u_{N+M}$ and $s_{M+N}$ as

[TABLE]

where we have made implicit the fact that $u^{(Y)}(U)$ acts on $\mathcal{U}^{(Y)}(\mathrm{SU}(d))$ and $s^{(Y)}(\sigma)$ acts on $\mathcal{U}^{(Y)}(S_{N+M})$ respectively. Given a state $\rho\in\mathcal{B}(\mathcal{H}_{d}^{\otimes(N+M)})$, and applying $u_{N+M}(U)$ with $U$ extracted from the Haar measure $\mathrm{d}U$ of $\mathrm{SU}(d)$ gives rise to the following completely positive, trace-preserving (CPTP) map

[TABLE]

Similarly, applying $s(\sigma)$ chosen uniformly at random gives rise to the following CPTP map

[TABLE]

By decomposing the representations $u_{N+M}$, and $s_{N+M}$ into their irreps as in Eq. (16), and making use of Schur’s lemmas Sternberg Eqs (17, 18) can be conveniently written as ReviewFramesInfo

[TABLE]

where $\mathcal{I}$ is the identity map, i.e., $\mathcal{I}_{\mathcal{H}}[A]=A$, $\mathcal{D}$ is the completely depolarizing map, i.e., \mathcal{D}_{\mathcal{H}}[A]=\frac{\mathrm{Tr}\left[A\right]}{\mathrm{dim}(\mathcal{H})}{\mbox{1{\bf l}}}_{\mathcal{H}}, and $\Pi_{Y}$ is the projector onto the block $\mathcal{U}^{(Y)}(\mathrm{SU}(d))\otimes\mathcal{U}^{(Y)}(S_{N+M})$.

In the following we will need some more observations. Define $\mathrm{Sym}_{N}^{(d)}$ as the completely symmetric subspace of $\mathcal{H}_{d}^{\otimes N}$; in particular, if taken $Y_{max}$ to be the label for a Young diagram of just one row, one has $\mathrm{Sym}_{N}^{(d)}=\mathcal{U}^{(Y_{max})}(\mathrm{SU}(d))\otimes\mathcal{U}^{(Y_{max})}(S_{N})$. Moreover, the tensor product space

[TABLE]

admits a decomposition

[TABLE]

where $J$ are the labels of Young diagram with two rows and $\mathcal{K}^{J}:=\mathrm{span}\{\left|J\right\rangle\}$, $\left|J\right\rangle\in\mathcal{U}^{(J)}(S_{M+N})$, $\mathcal{K}^{J}$ being a one dimensional space because the multiplicity of each irrep $\mathcal{U}^{(J)}(\mathrm{SU}(d))$ is one. Notice that the label $J$ can be indexed by an half-integer $\frac{M-N}{2}\leq J\leq\frac{M+N}{2}$.

The states of the form $\left(\left|\psi\right\rangle\!\!\left\langle\psi\right|\right)^{\otimes N}\otimes\left(\left|\phi\right\rangle\!\!\left\langle\phi\right|\right)^{\otimes M}$ are supported on $\mathrm{Sym}_{N,M}^{(d)}$. Define also $\mathcal{E}:=\mathrm{span}\{\left|\psi\right\rangle,\left|\psi_{\perp}\right\rangle\}\subseteq\mathcal{H}_{d}$. Since $\mathcal{E}$ has dimension 2, it is isomorphic to $\mathcal{H}_{2}$. Using this isomorphism we can define $\mathrm{Sym}^{\mathcal{E}}_{N}\cong\mathrm{Sym}^{(2)}_{N}$. Notice that $\left(\left|\psi\right\rangle\!\!\left\langle\psi\right|\right)^{\otimes N}\otimes\left(\left|\phi\right\rangle\!\!\left\langle\phi\right|\right)^{\otimes M}$ is also supported on the product of the two completely symmetric subspaces of $\mathcal{E}^{\otimes N}$ and $\mathcal{E}^{\otimes M}$, $\mathrm{Sym}^{\mathcal{E}}_{N,M}:=\mathrm{Sym}^{\mathcal{E}}_{N}\otimes\mathrm{Sym}^{\mathcal{E}}_{M}$, which is isomorphic to

[TABLE]

where now $K^{(2)}_{J}=\mathrm{span}\{\left|J\right\rangle_{(2)}\}$, $\left|J\right\rangle_{(2)}\in\mathcal{U}^{(J)}(S_{M+N})$. Using this isomorphism, we can write the decomposition

[TABLE]

and now $K_{\mathcal{E}}^{J}=\mathrm{span}\{\left|J\right\rangle_{\mathcal{E}}\}$, $\left|J\right\rangle_{\mathcal{E}}\in\mathcal{U}^{(J)}(S_{M+N})$. What is more, $\mathcal{E}^{\otimes N}\otimes\mathcal{E}^{\otimes M}$ is an invariant subspace of $\mathcal{H}_{d}^{\otimes M+N}$ under the action of the symmetric group, therefore the irreps of the symmetric group supported on $\mathcal{E}^{\otimes N}\otimes\mathcal{E}^{\otimes M}$ can also be taken as irreps in the decomposition of the representation of $S_{M+N}$ on $\mathcal{H}_{d}^{\otimes M+N}$, and it follows that $\mathcal{U}^{(J)}_{\mathcal{E}}(\mathrm{SU}(2))\otimes\mathcal{K}_{\mathcal{E}}^{J}\subseteq\mathcal{U}^{(J)}(\mathrm{SU}(d))\otimes\mathcal{U}^{(J)}(\mathcal{S}_{N+M})$. Finally, since $\mathrm{Sym}_{N,M}^{(2)}\subseteq\mathrm{Sym}_{N,M}^{(d)}$ and each of the $\mathcal{U}^{(J)}_{\mathcal{E}}(\mathrm{SU}(2))\otimes\mathcal{K}^{J}_{\mathcal{E}}$ and $\mathcal{U}^{(J)}(\mathrm{SU}(d))\otimes\mathcal{K}^{J}$ are included in the same subspace $\mathcal{U}^{(J)}(\mathrm{SU}(d))\otimes\mathcal{U}^{(J)}(\mathcal{S}_{N+M})$, it follows that

[TABLE]

and $\left|J\right\rangle=\left|J\right\rangle_{\mathcal{E}}$ for each $J$.

Appendix B Average state at fixed overlap

In this section we derive the form of the average state of $N$ copies of a pure state $\left|\psi\right\rangle\in\mathcal{H}_{d}$ and $M\geq N$ copies of a pure state $\left|\phi\right\rangle\in\mathcal{H}_{d}$ with $\left|\psi\right\rangle$ and $\left|\phi\right\rangle$ having fixed overlap. Using group representation theory techniques we find a basis in which the average state is diagonal and compute the eigenvalues. We generalize the result to the average state in the unlabeled scenario.

Without loss of generality we may write the latter of the two states as

[TABLE]

where $c=\cos^{2}\frac{\theta}{2}\in(0,1),\,\theta\in(0,\pi),\,\left\langle\psi_{\perp}\lvert\psi\right\rangle=0$.

As both $\left|\psi\right\rangle,\,\left|\phi\right\rangle\in\mathcal{H}_{d}$ are randomly chosen the global state describing the $N+M$ qudits is given by

[TABLE]

where $dU$ is the Haar measure of $SU(d)$. We can decompose the integral over $\mathrm{SU}(d)$ in Eq. (26) as follows. First we consider the two-dimensional subspace $\mathcal{E}\equiv\mathrm{span}\{\left|\psi\right\rangle,\left|\psi_{\perp}\right\rangle\}\subseteq\mathcal{H}_{d}$. Using the invariance of the Haar measure, we first perform the group average over the $\mathrm{SU}(2)$ subgroup of $\mathrm{SU}(d)$ which acts non-trivially only on $\mathcal{E}$, which we denote as $\mathrm{SU}(\mathcal{E})$. Afterwards, we can average over all $\mathrm{SU}(d)$. This implies that Eq. (26) can be written as

[TABLE]

First of all we compute the action of $\mathcal{G}_{\mathrm{SU}(2)}$ on $\left(\left|\psi\right\rangle\!\!\left\langle\psi\right|\right)^{\otimes N}\otimes\left(\left|\phi\right\rangle\!\!\left\langle\phi\right|\right)^{\otimes M}$, when $d=2$. We use the addition rules for angular momentum on $\mathcal{H}_{2}^{\otimes N+M}$ to write

[TABLE]

where $C^{J,\frac{N}{2}+k}_{\frac{N}{2},\frac{N}{2};\frac{M}{2},k}=\left\langle J,\frac{N}{2}+k\lvert\frac{N}{2},\frac{N}{2};\frac{M}{2},k\right\rangle$ are the Clebsch-Gordan coefficients. Using Eq. (17) and the decomposition in Eq. (23) ($\left(\left|\psi\right\rangle\!\!\left\langle\psi\right|\right)^{\otimes N}\otimes\left(\left|\phi\right\rangle\!\!\left\langle\phi\right|\right)^{\otimes M}$ is supported on $\mathrm{Sym}_{N,M}^{\mathcal{E}}$) the first Haar-measure group average of Eq. (27) reads

[TABLE]

where we have used the total angular momentum $J$, instead of the Young frame label $Y$, $\chi_{J}=2J+1$, and $\left|J\right\rangle\!\!\left\langle J\right|_{\mathcal{U}^{(J)}(S_{N+M})}$ is a pure state. The coefficients of Eq. (29) are:

[TABLE]

and we have made use of the following expression for the Wigner $D$ matrix in going from the second to the third line in Eq. (30)

[TABLE]

with $P_{n}^{(\alpha,\beta)}(x)$ the Jacobi polynomials, defined in general as

[TABLE]

Let’s consider now $d>2$. Notice that $\mathcal{E}^{\otimes M+N}\cong\oplus_{J}\mathcal{U}^{(J)}_{\mathcal{E}}(\mathrm{SU}(2))\otimes\mathcal{U}_{\mathcal{E}}^{(J)}(S_{N+M})$. Performing the average over $\mathrm{SU}(d)$ on the state given by Eq. (29) results in

[TABLE]

where in the last equality we used the fact that the support of {\mbox{1{\bf l}}}_{\mathcal{U}^{(J)}_{\mathcal{E}}(\mathrm{SU}(2))}\otimes\left|J\right\rangle\!\!\left\langle J\right|_{\mathcal{U}^{(J)}_{\mathcal{E}}(S_{N+M})} is a subspace of the support of {\mbox{1{\bf l}}}_{\mathcal{U}^{(J)}(\mathrm{SU}(d))}\otimes\left|J\right\rangle\!\!\left\langle J\right|_{\mathcal{U}^{(J)}(S_{N+M})}, as stated in Eq. 24.

Hitherto in the computation we tacitly assumed that the copies of states are labelled, i.e., that the first $N$ states are all $\left|\psi\right\rangle\in\mathcal{H}_{d}$ and the remaining $M$ states are all $\left|\phi\right\rangle\in\mathcal{H}_{d}$. We now lift this assumption and derive the form of the average state in the case where the copies of the states are unlabelled. This is equivalent to averaging over all possible permutations of the $N+M$ copies, i.e., by applying the map $\mathcal{G}_{S_{N+M}}[\rho(c)]$ of Eq. (19). The final state can be easily shown to be

[TABLE]

by observing that the maps $\mathcal{G}_{\mathrm{SU}(d)},\,\mathcal{G}_{S_{N+M}}$ both project on the same irrep label $J$, but depolarize the states on different irrep subspaces $\mathcal{U}^{(J)}(\mathrm{SU}(d))$ and $\mathcal{U}^{(J)}(S_{N+M})$ respectively. As all the information concerning the overlap between the two states is extracted by the projective measurement $\{\Pi_{J}\}_{J}$, and the action of any permutation–or global unitary rotation–does not alter the statistics $p(J|c)$. Consequently, all our results apply equally well to both labelled and unlabelled scenarios.

Appendix C Fisher Information of $p(J|c)$

In this section we derive the Fisher information of the probability distribution given by Eq. (30), and analyze its asymptotics. We can obtain the first non trivial order by using an asymptotic approximation for the Jacobi polynomials and large deviation evaluations. Using the same techniques, we identify an asymptotically unbiased estimator that saturates the Cramer-Rao bound.

The result we obtain is summarized by the following theorem:

Theorem**.**

Given $N$ and $M$ copies of two Haar-random states, respectively $\left|\psi\right\rangle$ and $\left|\phi\right\rangle$, define $\left|\Psi\right\rangle=\left|\psi\right\rangle^{\otimes N}\otimes\left|\phi\right\rangle^{\otimes M}$. The mean square error on the estimation of the overlap $c=|\left\langle\psi\lvert\phi\right\rangle|^{2}$ attained by the weak Schur sampling projectors $\{E_{J}\}_{J}$ with the estimator $c_{J}=\left(\frac{J}{M+N}\right)^{2}$ is

[TABLE]

This measurement is optimal in the sense that saturates the quantum Cramer-Rao bound for the problem of estimating of $c$ from the family of average states at fixed overlap

[TABLE]

This implies that $c$ can be estimated at a relative error $\frac{\sqrt{v(c)}}{c}=\epsilon$ with $N=M=\mathcal{O}(\frac{1}{\epsilon^{2}c})$ samples.

We begin by recalling the definition of the Fisher information

[TABLE]

Using the identity

[TABLE]

it follows that

[TABLE]

For $x\geq 1$ the following asymptotic expansion for the Jacobi polynomials holds Szego1939 ; Darboux1878 . Defining the function

[TABLE]

one has

[TABLE]

where the convergence is uniform on any half-line $x\in[1+\delta,+\infty)$, $\delta>0$. We want to use the leading order as an approximation to the Jacobi polynomial. Notice that the remainder in the expansion does not depend on $J_{\max}$. Take $J_{0}$ such that for each $J>J_{0}$ and a certain $C>0$

[TABLE]

We can then split the sum in Eq. 37 in two parts. The first is

[TABLE]

Writting $p(J|c)$ as

[TABLE]

and noting that

[TABLE]

it follows that for $J\leq J_{0}$

[TABLE]

On the other hand, isolating in $p(J|c)$ the terms depending on $J_{\max}$, and defining $w(J,J_{\min},c)=(2J+1)\sum_{s=0}^{J-J_{\min}}\binom{J+J_{\min}}{s}\binom{J-J_{\min}}{s}c^{s}$ which is bounded for $J_{\min}\leq J\leq J_{0}$, we have

[TABLE]

By introducing the following binomial distribution

[TABLE]

with $p=\frac{1-\sqrt{c}}{2}$, whose mean and variance are given by

[TABLE]

Eq. (47) can be written as

[TABLE]

with $w^{\prime}(J,J_{\min},c)$ still bounded.

Now by Sanov’s theorem Sanov1957 , the binomial distribution $q(J)$ is exponentially suppressed as $J_{\max}\rightarrow\infty$ and $J<J_{0}$. In particular, if $J_{\max}>\frac{1}{\sqrt{c}}J_{0}$

[TABLE]

where $D(\frac{J_{\max}-J_{0}}{2J_{\max}+1}||\frac{1-\sqrt{c}}{2})$ is the relative entropy between the two Bernoulli distributions with probabilities $p_{1}=\frac{J_{\max}-J_{0}}{2J_{\max+1}}$ and $p_{2}=\frac{1-\sqrt{c}}{2}$. For $J_{\max}\rightarrow\infty$, $p_{1}\rightarrow\frac{1}{2},\,D(\frac{J_{\max}-J_{0}}{2J_{\max}+1}||\frac{1-\sqrt{c}}{2})\rightarrow D(\frac{1}{2}||\frac{1-\sqrt{c}}{2})>0$ unless $c=0$. As the rest of the terms in both $p(J|c)$ and $\left|\left(\frac{\frac{\mathrm{d}p(J|c)}{\mathrm{d}c}}{p(J|c)}\right)\right|$ tend to a power law in $J_{max}$, $H(c)_{1}$ is exponentially suppressed.

As an indication for which ratio $\frac{J_{\max}}{J_{\min}}$ is expected to give exponential suppression when $c$ is small, we can look at the first order of the Taylor expansion of $D(p_{1}||p_{2})$ in $\frac{1}{J_{\max}}$ and $\sqrt{c}$ around zero: $D(p_{1}||p_{2})=-\frac{1}{2}\log(1-c)-\frac{1}{4}\frac{2J_{0}+1}{J_{max}}\log\frac{1+\sqrt{c}}{1-\sqrt{c}}$, which for small $c$ requires $J_{\max}c>>1$ from the zeroth order and $J_{\max}>>\frac{1}{\sqrt{c}}J_{0}$ from the first order, implying also $J_{\max}>>\frac{1}{\sqrt{c}}J_{\min}$.

It remains to evaluate

[TABLE]

Using Eq. (41) one obtains, after some algebra

[TABLE]

where the same binomial distribution, $q(J)$ of Eq. (48), appears again. Writing

[TABLE]

where

[TABLE]

we have

[TABLE]

In the limit $J_{\max}\to\infty$

[TABLE]

whilst $g(J,J_{\max},J_{\min},c)=\mathcal{O}(J_{\max}^{\frac{3}{2}})$ by inspection.

Now consider the first $\sqrt{d}J_{\max}-(J_{0}+1)$ terms in Eq. (52), where $d<c$,

[TABLE]

Sanov’s theorem tells us that these are exponentially suppressed,

[TABLE]

and as a result we are only left with having to evaluate

[TABLE]

As we are interested in computing the $QFI$ in the limit $J_{\max}\to\infty$ we will perform a Taylor series expansion of $H(c)_{2b}$, around the mean value $\langle J\rangle=J_{\max}-\mu=\sqrt{c}\left(J_{\max}+\frac{1}{2}\right)-\frac{1}{2}$ of $q(J)$. For any analytical function $h(J)$, using the Lagrange remainder, an expansion about $\langle J\rangle$ is given by

[TABLE]

with $J\leq\xi(J)\leq\langle J\rangle$ (so that in any case $\sqrt{d}J_{\max}\leq\xi(J)\leq J_{\max}$). Using Eq. (61) and the fact that $g(J,J_{\max},J_{\min},c)$ is analytic for $J>J_{\min}$ by inspection, we have

[TABLE]

where the change in the order of the summation is allowed since one of the sums is finite. Applying Sanov’s theorem, Eq. (51), gives

[TABLE]

It follows that we can change the summation in $J$ of Eq. (62) to start from $-J_{\max}-1$ just by adding a term exponentially suppressed in $J_{\max}$. For the Lagrange remainder, before extending the summation, we notice that for $\sqrt{d}J_{\max}\leq\xi(J)\leq J_{\max}$, $\left.\frac{d^{s}g(J^{\prime},M,N,c)}{dJ^{\prime s}}\right\rvert_{J^{\prime}=\xi(J)}=\mathcal{O}(J_{\max}^{\frac{3}{2}-n})$ (by inspection).

Finally, using Eq. (49) and some algebra one obtains

[TABLE]

Multiplying Eqs. ((57), (64)) gives the final result

[TABLE]

The neglected term $\mathcal{O}\left(\frac{1}{J+J_{\min}}\right)$ in Eq. (56) contributes to $H(c)$ with a next to leading order term, by simple power counting, adapting the steps above.

Using similar techniques it can be shown that $H(c)$ is achievable by the estimator of Eq. (7) of the main text. First note that the estimator is indeed asymptotically unbiased since

[TABLE]

apart from an exponentially suppressed term, and we have again used the asymptotic expansion of Eq. (41) and identified the binomial distribution of Eq. (48). In the limit $J_{\max}\to\infty$ one can proceed as above, separating the $J$-dependence and applying Eq. (61), to obtain

[TABLE]

It follows that

[TABLE]

A similar calculation gives for the MSE of the estimator

[TABLE]

Appendix D Optimal global mean squared error

In this section we derive the optimal estimator and corresponding global mean squared error (glMSE) for the case where the overlap $c$ is a random variable with a distribution induced by the Haar-uniform measure of $\mathrm{SU}(d)$.

The result we obtain is summarized by the following theorem:

Theorem**.**

Given $N$ and $M$ copies of two Haar-random states, respectively $\left|\psi\right\rangle$ and the state $\left|\phi\right\rangle$, define $\left|\Psi\right\rangle=\left|\phi\right\rangle^{\otimes N}\otimes\left|\psi\right\rangle^{\otimes M}$. The global mean square error on the estimation of the overlap $c=|\left\langle\psi\lvert\phi\right\rangle|^{2}$ attained by a POVM $\{E_{k}\}_{k}$ and an estimator $c_{k}$ is

[TABLE]

The minimum of $v$ is attainable by a projective measurement and it reads

[TABLE]

We follow the standard treatment in Personick to compute the estimator that minimizes glMSE, and use representation theory to perform the integrations.

The probability distribution of the overlap, induced by the Haar measure, is  (OverlapStatistics, , Eq. (13))

[TABLE]

Following Personick the optimal estimator, $S$, satisfies

[TABLE]

where

[TABLE]

and is explicitly given by

[TABLE]

Plugging Eq. (29) into Eqs. (74,75) gives

[TABLE]

where $p(J,c)=p(c)p(J|c)$, $\int p(J,c)\,\mathrm{d}c=p(J)$ and we have used Bayes’ theorem in going from the second to the third line of Eq. (76). Again, the optimal measurement corresponds to measuring the total irrep label $J$. Upon a given outcome the estimator that minimizes the glMSE is $c_{opt}^{bay}(j)=\langle c\rangle_{J}$ where the expectation value is taken with respect to the conditional probability distribution $p(c|J)$.

The operators $\Gamma$ and $\eta$ are given by

[TABLE]

where we have made use of the fact that $W_{V}V=T(c)$ for a unitary $W_{V}$ such that $W_{V}\left|\psi\right\rangle\left\langle\psi\right|W_{V}^{\dagger}=\left|\psi\right\rangle\left\langle\psi\right|$ in the third equality, the invariance of the Haar measure for the fourth equality, and used the addition rules for $\mathrm{SU}(d)$ representations for the last equality.

To compute $\eta$ we make use of

[TABLE]

with $D^{j}_{m,n}(g)$ being Wigner matrices, so that

[TABLE]

It is now trivial to compute the optimal Bayesian estimator for a given measurement outcome $J$:

[TABLE]

$\chi_{J}$ can be calculated with the Hook formula HayashiGroupRep . For Young diagrams with two rows one has

[TABLE]

To simplify the denominator we employ the graphical calculus techniques in  Varshalovich1988 :

[TABLE]

where the term in curly brackets is the Wigner 6-j symbol. Plugging everything together the optimal glMSE estimator for a given measurement outcome $J$ is given by

[TABLE]

with its corresponding glMSE

[TABLE]

Appendix E 1-LOCC strategies

Here we derive the pointwise MSE attainable by the estimate-and-project (EP) and estimate-and-estimate (EE) strategies, and study the optimal Bayesian estimator for both.

E.0.1 Estimate-and-project: local estimation

The EP strategy consists in first estimating at best one of the states using an optimal collective measurement on all of its copies, and then projecting each copy of the other state on the estimate of the first one. The optimal estimation of a random state $\left|\phi\right\rangle$, given $N$ copies of it, is provided by the covariant measurement of Ref. Hayashi1997 , $\{M_{dV}\}$, that produces an estimate $\left|\phi_{V}\right\rangle=V\left|0\right\rangle$, $V\in\mathrm{SU}(d)$, with probability density

[TABLE]

Finally, we perform the projective measurement \{\left|\phi_{V}\right\rangle\left\langle\phi_{V}\right|,{\mbox{1{\bf l}}}-\left|\phi_{V}\right\rangle\left\langle\phi_{V}\right|\} on each copy of $\left|\psi\right\rangle$. This succeeds with probability

[TABLE]

The overall measurement operator is then $E_{V,k}^{(ep)}=\mathrm{d}VE_{V}^{(M)}\otimes V^{\otimes N}\Pi_{k}^{(N)}V^{{\dagger}\otimes N}$, as defined in the main text. Its outcome statistics for a fixed value of the overlap can be written as

[TABLE]

where we have set without loss of generality $\left|\psi\right\rangle=\left|0\right\rangle$, $\left|\phi\right\rangle=U_{c}\left|0\right\rangle$ thanks to the presence of the reference-frame average. We now consider both local and Bayesian estimation with this strategy.

For local estimation with EP, we employ as an estimator $c_{\rm ep}^{\rm loc}(k)=\frac{k}{N}$, i.e., the fraction of successful projections. This estimator is in general biased, except in the limit of large $N$, and it is not necessarily optimal but provides a natural guess for the overlap given our strategy. Its MSE is

[TABLE]

where in the first equality we have introduced the binomial distribution $\mathrm{Bin}(k,N,p)$ of $m$ successes out of $N$ trials, with a single-trial success probability $p$, in the second equality we expanded the square and used the mean and variance of $\mathrm{Bin}(k,N,p)$, and in the third equality we defined the integrals

[TABLE]

These can be computed by expanding the scalar products and writing the states in a collective-spin basis, obtaining the expectation value of the operator in Eq. (26), with the substitution $N\mapsto i$:

[TABLE]

where we have eliminated the integral over $V$ by including this rotation into the group average. The result is then given by Eqs. (30,33) with the same substitution. Plugging in the expression of \chi_{J}=\left(\begin{array}[]{c}2J+d-1\\ d-1\end{array}\right), Ref. Hayashi1997 , we have

[TABLE]

which can be computed explicitly for $i=1,2$. Inserting these expressions in Eq. (86) we finally obtain

[TABLE]

In the limit $M\rightarrow\infty$, $N$ constant we have

[TABLE]

which coincides with the optimal strategy, corresponding to a projection on the known direction of $\left|\phi\right\rangle$. In the limit $M+N\rightarrow\infty$, $M-N$ fixed we have instead

[TABLE]

which is $3/2$ times larger than the optimal strategy.

E.0.2 Estimate-and-project: Bayesian estimation

For Bayesian estimation with EP, as in the previous section, the optimal classical estimator is given by

[TABLE]

We start by computing the probability distribution of the outcomes, using Eq. (85):

[TABLE]

where we performed the average over $c\in\mathrm{SU}(d)$ in the third equality, defining {\mbox{1{\bf l}}}_{M}^{sym}={\mbox{1{\bf l}}}_{\mathcal{U}^{\left(\frac{M}{2}\right)}(\mathrm{SU}(d))} as the projector on the completely symmetric subspace of $M$ qudits, and employed its invariance under $U^{\otimes M}$. We are then left to compute the overlap of $\Pi_{k}^{(N)}$ with {\mbox{1{\bf l}}}^{sym}_{N}. In order to do so, we recall that the latter can be written as the average of all permutation operators $V(\sigma),\,\sigma\in S_{N}$ of $N$ $d$-level systems. Then the symmetrization of $\Pi_{k}^{(N)}$ gives a trivial binomial factor and we can write

[TABLE]

where $\beta_{j}(\underline{i})$ is the number of times the integer $j$ appears in the sequence $\underline{i}$. The sum over $\underline{i}$ can then be broken up into the sum over all partitions of $N-k$ systems in $d-1$ sectors, i.e., the sum over all possible vectors $\underline{\beta}$ of $d-1$ components that add up to $N-k$, times the sum of all permutations of $N-k$ systems which are equal in groups of size $\beta_{j}$. The latter can be carried out immediately since the summand is invariant under permutation of the $j$’s:

[TABLE]

where in the last equality we summed over $\underline{\beta}$. For the first moment of the distribution we get, similarly,

[TABLE]

where we have introduced $U^{{\dagger}}U$ and its conjugate on the additional subsystem, then employed the invariance of {\mbox{1{\bf l}}}_{M+1}^{sym} under $U^{\otimes(M+1)}$. In order to proceed, we first compute the value of the following operator:

[TABLE]

where in the first equality we have taken the partial trace over $M$ subsystems, while in the second one we have written {\mbox{1{\bf l}}}_{M+1}^{sym} as an average of permutations, like before, and written the explicit basis representation of the last subsystem. The third equality follows by evaluating the only non-zero elements in the sums: the first term contains all the permutations of $M$ subsystems times the identity on the remaining subsystem; the second term considers the additional permutations in the case $i=j=0$, where the last subsystem can be exchanged with any of the other $M$ subsystems. By substituting this expression in Eq. (E.0.2) we obtain

[TABLE]

where we have used Eq. (94). We then just need to compute the second term in the sum above, which is very similar to Eq. (95) with the change $\left|0\right\rangle^{\otimes k}\mapsto\left|0\right\rangle^{\otimes(k+1)}$. We finally obtain

[TABLE]

and the optimal Bayesian estimator for each $k$ is

[TABLE]

The corresponding glMSE is given by

[TABLE]

In the limit $M\rightarrow\infty$, $N$ constant we have

[TABLE]

which again coincides with the optimal Bayesian strategy. In the limit $M+N\rightarrow\infty$, $M-N$ fixed we have instead

[TABLE]

which again is $3/2$ times larger than the optimal Bayesian strategy.

E.0.3 Estimate-and-estimate: Bayesian estimation

The EE strategy instead consists in estimating both states with a covariant measurement, hence it is described by overall POVM operators $E_{V,k}^{(ep)}=\mathrm{d}V\mathrm{d}WE_{V}^{(M)}\otimes E_{W}^{(N)}$, as mentioned in the main text. Its success probability can be written as

[TABLE]

where again we could include one of the outcomes into the unitary average and, since by redefining $W\mapsto V^{{\dagger}}W$ the dependence on $V$ is a constant, we averaged over $V$ without loss of generality.

In the Bayesian case for EE we proceed as before and compute first

[TABLE]

then

[TABLE]

The second equality above comes from averaging over $U_{c}$, the third one from Eq. (98) and introducing $U^{{\dagger}}U$ and its conjugate on the additional subsystem, the fourth one from redefining $U\mapsto W^{{\dagger}}U$, switching the operators acting on the $N$ subsystems and averaging over $U$, while the fifth one from applying Eq. (98) again and defining $w=\left\lvert\left\langle 0\right|W\left|0\right\rangle\right\rvert^{2}$.

Then the optimal EE Bayesian estimator for each $W$ is simply $c(W)=\tilde{c}(W)/dW$ and the minimum glMSE attained by it is

[TABLE]

where we have carried out the group averages in the usual way:

[TABLE]

In the limit $M\rightarrow\infty$, $N$ constant we have

[TABLE]

which is $(d+2N)/(d+N)$ times larger than the optimal Bayesian strategy. In the limit $M+N\rightarrow\infty$, $M-N$ fixed we have instead

[TABLE]

which is $2$ times larger than the optimal Bayesian strategy.

E.0.4 Estimate-and-estimate: local estimation

Finally, for the local EE estimation the estimator $\tilde{c}_{W}=w=\left\lvert\left\langle 0\right|W\left|0\right\rangle\right\rvert^{2}$ is a natural guess. Its variance can be computed in terms of its first and second moments according to the distribution $p(W|c)$:

[TABLE]

where

[TABLE]

Then the first moment is straightforward to compute by inserting Eq. (98):

[TABLE]

For the second moment we first need to evaluate the following operator:

[TABLE]

As before, the third equality follows by evaluating the only non-zero elements in the sums: the first term contains all the permutations of $N$ subsystems times the identity and the swap on the remaining subsystems, which add up to the projector on the completely symmetric subspace of the two subsystems,

[TABLE]

the second term considers the additional permutations in the case $i_{1}=j_{1}=0$ and $i_{2}=j_{2}=0$, where one of the remaining subsystems can be swapped or not with the other, then permuted with any of the other $N$ subsystems; analogously, the third term considers the additional permutations in the case $\underline{i}=\underline{j}=\underline{0}$, where each remaining subsystem can be permuted respectively with $N$ and $N-1$ of the others. Hence the second moment of $w$ can be written as

[TABLE]

By plugging the expressions of Eqs. (114,117) into Eq. (112) we finally get

[TABLE]

In the limit $M\rightarrow\infty$, $N$ finite we have

[TABLE]

which is twice as large as the optimal strategy in the leading order of $N$. In the limit $M+N\rightarrow\infty$, $M-N$ fixed we have instead

[TABLE]

which is $2$ times larger than the optimal strategy.

Appendix F Bayesian estimation using the swap test

Here we derive the formulas employed for the plot of the swap-test performance in Fig. 3 of the main text. The measurement statistics is given by a binomial distribution $\mathrm{Bin}(k,N,p(c))$, of $k$ events out of $N$, with single-event probability $p(c)=\frac{1+c}{2}$. The corresponding optimal classical Bayesian estimator is given by Eq. (93), with $p(k|c)\rightarrow\mathrm{Bin}(k,N,p(c))$. We have

[TABLE]

where ${}_{2}F_{1}$ and ${}_{2}\tilde{F}_{1}$ are hypergeometric and regularized hypergeometric functions respectively. Following the derivations of the previous section, the minimum glMSE attainable with the swap test can be written as

[TABLE]

Appendix G Comparison between Schur transform and swap test with imperfect implementations

In this section we sketch an evaluation of the effect of imperfect gates on the accuracy of the estimate of the overlap. First of all we model the error of each iteration of the swap test as white noise for each iteration: $\mathcal{N}_{sw}(p(c))=(1-\epsilon_{sw})p(c)+\epsilon_{sw}\frac{1}{2}$, the Fisher information becomes

[TABLE]

For $N$ repetitions, one gets a resulting MSE

[TABLE]

We model the noise on the Schur transform measurement outcomes also as mixing with a probability distribution $q(c)$: $\mathcal{N}_{sw}(p(J|c))=(1-\epsilon_{Sch})p(J|c)+\epsilon_{Sch}q(c)$, with a probability of mixing that scales exponentially in the number of gates, $1-\epsilon_{Sch}\approx(1-\epsilon)^{g}$, where $\epsilon$ is the error per gate, and $g$ is the total number of gates. We recall the joint convexity property of the Fisher Information, coming from its monotonicity:

[TABLE]

If we assume $q(c)$ to be overlap independent, we obtain the bound

[TABLE]

so that

[TABLE]

This is a very conservative estimate, as we are assuming we acquiring useful information with exponentially small probability. Hence the Swap test outperforms our optimal strategy, based on the Schur transform, when the respective implementation errors satisfy the following relation

[TABLE]

One can express $\epsilon_{Sch}$ and $\epsilon_{Sw}$ in terms of the error per gate, $\epsilon$, raised to gate complexity of their respective circuits. An intermediate strategy could be to divide the $N$ copies of both $\left|\phi\right\rangle$ and $\left|\psi\right\rangle$, into $R$ groups of $S$ copies, and perform the optimal measurement on each group, followed by classical post-processing. If $N=M=RS$ and $F(J|c,S)$ is the optimal Fisher information for the case with $M=N=S$ copies, the Cramer-Rao bound reads

[TABLE]

The best option would be to choose $S$ as the highest number of copies such that the architecture can perform the optimal measurement in a sufficiently precise way. On the other hand, if one requires to be in the asymptotic regime of the approximation for $c>c_{0}$, one can just find the minimum $S$ for which the approximation works, and perform the optimal measurement with $S$ copies $R$ times. The classical post processing will have the optimal asymptotic performance for $c>c_{0}$. In any case the bound (130) is asymptotically achieved by a maximum likelihood estimator when $R\rightarrow\infty$.

Appendix H Average post-measurement fidelity

Here we compute the average post-measurement fidelity for the optimal strategy and the swap test. For the former we have measurement operators $\{\Pi_{J}\}_{J}$, so that Eq. (8) of the main text reads

[TABLE]

where we have used the $\mathrm{SU}(d)$-invariance of $\Pi_{J}$. For the swap test we restrict to $M=N$ as usual and we consider that the measurement is separable and identical on each couple of copies. Moreover, the measurement on a single pair of copies is a triplet/singlet projection, which is again $\mathrm{SU}(d)$-invariant, and succeeds/fails with probability $(1\pm c)/2$. Hence we have

[TABLE]

Appendix I Estimating the overlap between two arbitrary mixed qubits

In this appendix we derive the optimal estimator and corresponding mean squared error for the case where we are given $N$ and $M$ copies of mixed states. We shall restrict our attention to qubit mixed states and for ease of notation we shall revert to the standard angular momentum notation for irrep labels.

The mixed states whose overlap we wish to estimate are

[TABLE]

where $r_{0(1)}$ denotes the corresponding purity of the states. Following EstimationGRC the states in Eq. (133) can be written in the total angular momentum basis, after tracing out multiplicities, as

[TABLE]

where

[TABLE]

with $R_{i}=\frac{1+r_{i}}{1-r_{i}},\,Z^{(i)}_{J_{i}}=\frac{R_{i}^{J_{i}+1}-R_{i}^{-J_{i}}}{R_{i}-1}$, and just as for the case of pure states, we have chosen $\vec{n}_{0}=\vec{z}$ without loss of generality. Moreover,

[TABLE]

and similarly for $p_{J_{1}}$. Using Eq. (17) we obtain

[TABLE]

where $\sigma^{(J_{0},J_{1})}\in\mathcal{B}(\mathcal{U}^{(J)}(S_{N+M}))$ and they are orthogonal for different couples $(J_{0},J_{1})$. To calculate the glMSE we need to compute the operators $\Gamma,\,\eta$ of Eq. (74). A similar calculation as in Eq. (77), (78) gives

[TABLE]

For $\eta$ one obtains

[TABLE]

We finally obtain

[TABLE]

For a given $J,J_{0},\,J_{1}$ and overlap $c$ the estimator is given by

[TABLE]

and the glMSE reads

[TABLE]

The sums in $L,L^{\prime},m,k$ can be done exactly. The sums in $J$ at the leading order in $R_{0},R_{1}$ can be done by keeping track of the non-exponentially decaying (in $J$) contributions. The final sum in $J_{0}$ and $J_{1}$ can be done using the fact that $p_{J_{0}}\frac{R_{(0)}^{J_{0}}}{Z_{J_{0}}}$ can be written as

[TABLE]

and in the limit $M=\alpha Z,\,N=\beta Z,\,Z\to\infty$ one can approximate the glMSE expanding in moments around the mean of the binomial distribution. The final result reads

[TABLE]

in agreement with the pure state case for $d=2,\,r_{0}=r_{1}=1$.