What is the relation between number of measurements and additive error, when estimating the expectation value of an observable?

Question

Given an unknown state $\rho$ I would like to measure it using an observable $O$. I would perform measurements using the orthogonal basis of $O$ on multiple copies of $\rho$ and my answer would approach $\text{tr}[O\rho]$.

I imagine that if I have $n$ copies of the state $\rho$ then my error, $\epsilon = |\text{tr}[O\rho] - \sum_i \mu_i|$ would depend on $n$. I would like to know the formula that considers $O$ and $n$.

I am not measuring in computational basis but the basis of $O$. Further I am not trying to estimate the entire state-vector making it much easier from the standard tomography problem.

Answer 1

Let us denote $\left\{\lambda_i\right\}_i$ the eigenvalues of $O$. We can write $O$ as $UDU^\dagger$ with $U$ being unitary and $D$ being diagonal. We then have: $$\DeclareMathOperator{\tr}{tr}\tr\left(O\rho\right)=\tr\left(UDU^\dagger\rho\right)=\tr\left(DU^\dagger\rho U\right)$$ So, what you want to do is to apply the $U^\dagger$ unitary on $\rho$ and then measure in the "computational basis" (not exactly since each measurement outcome is associated with $\lambda_i$, not $\pm1$). Thus, it is equivalent to answer your question for an arbitrary diagonal observable $D$. By linearity of the trace and of the expectation value, we can also afford to assume that $\rho$ is pure: $$\rho=|\psi\rangle\!\langle\psi|$$ We thus want to evaluate $\tr\left(D|\psi\rangle\langle\psi|\right)=\langle\psi|D|\psi\rangle$. Let us denote: $$|\psi\rangle=\sum_i\psi_i|i\rangle$$ We thus have: $$\langle\psi|D|\psi\rangle=\sum_i\lambda_i\left|\psi_i\right|^2$$ Note that when measuring $|\psi\rangle$ in the $\left\{|i\rangle\right\}_i$ basis, we obtain $\lambda_i$ with probability $\left|\psi_i\right|^2$. Thus, we can see this as sampling a random variable $X$ whose expectation value is given by $\langle\psi|D|\psi\rangle$.

Now we only have to deal with probability theory.

---

First result: Using Bienaymé–Chebyshev inequality

Given $n$ i.i.d. measurements $X_1,\cdots,X_n$, we have, using the Bienaymé–Chebyshev inequality: $$\mathbb{P}\left[\left|\frac1n\sum_{i=1}^nX_i-\mathbb{E}\left[X_1\right]\right|\geqslant\varepsilon\right]\leqslant\frac{1}{\varepsilon^2}\mathbb{V}\left[\frac1n\sum_{i=1}^nX_i\right]=\frac{\mathbb{V}\left[X_1\right]}{n\varepsilon^2}$$

---

Second result: Upper-bounding the mean error

Given $n$ i.i.d. measurements $X_1,\cdots,X_n$, we want to estimate $$\mathbb{E}\left[\left|\frac1n\sum_{i=1}^nX_i-\mathbb{E}\left[X_1\right]\right|\right]=\mathbb{E}\left[\sqrt{\left(\frac1n\sum_{i=1}^nX_i-\mathbb{E}\left[X_1\right]\right)^2}\right]$$ Using Jensen's inequality, we thus have: $$\mathbb{E}\left[\left|\frac1n\sum_{i=1}^nX_i-\mathbb{E}\left[X_1\right]\right|\right]\leqslant\sqrt{\mathbb{E}\left[\left(\frac1n\sum_{i=1}^nX_i-\mathbb{E}\left[X_1\right]\right)^2\right]}=\sqrt{\mathbb{V}\left[\frac1n\sum_{i=1}^nX_i\right]}=\sqrt{\frac{\mathbb{V}\left[X_1\right]}{n}}$$

---

In both cases, we would like to upper-bound $\mathbb{V}\left[X_1\right]$ independently of $|\psi\rangle$. We can use Popoviciu's inequality for this purpose, giving us: $$\mathbb{V}\left[X_1\right]\leqslant\frac14\left(\lambda_{\text{max}}-\lambda_{\text{min}}\right)^2$$ Where $\lambda_{\text{max}}$ and $\lambda_{\text{min}}$ are respectively the largest and smallest eigenvalues of $O$ (which are real as a recall).

All in all, both results intuitively say that in order to achieve an (expected) error of $\varepsilon$, you need $n=\mathcal{O}\left(\frac{1}{\varepsilon^2}\right)$ samples.

Answer 2

You can frame this question as follows: define an estimator which gives an estimate for the expectation value you seek for each possible measurement outcome. We are assuming to be measuring in the eigenbasis of the operator, therefore a natural choice of estimator is $$\hat o(b) = \langle b|O|b\rangle=\lambda_b,$$ where the observable's eigendecomposition reads $O=\sum_b\lambda_b |b\rangle\!\langle b|$. Note that this does what you want because on average $\hat o$ is exactly $\operatorname{tr}(O\rho)$, when the input state is $\rho$. More precisely you see this observing that $$\mathbb{E}[\hat o|\rho] = \sum_b \langle b|\rho|b\rangle \hat o(b)= \operatorname{tr}(O\rho). $$ To figure out the additive errors obtained with a given number of measurements, the standard approach is to compute the variance of this estimator, which directly translates into the mean squared error (MSE) for unbiased estimators. The variance reads $$\mathbb{E}[\hat o^2]-\mathbb{E}[\hat o]^2 = \sum_b \langle b|\rho|b\rangle \hat o(b)^2-\operatorname{tr}(O\rho)^2 = \langle O^2\rangle_\rho - \langle O\rangle_\rho^2 \equiv \sigma^2_O.$$ In other words, for this simple choice of measurement and estimator, the variance is precisely the "intrinsic variance" of the observable (note that this is not the case for the estimators you obtain for more general POVMs; you can see e.g. https://arxiv.org/abs/2301.13229 and references therein for how the analysis is performed in those instances).

Once you have the variance, it's only a matter of figuring out which bound you can apply to get your estimates. A standard approach is to compute a lower bound on the number of measurements required to obtain an additive error lower than $\epsilon$ with "high probability" (ie probability higher than $1-\delta$). A standard approach is to use Chebyshev's inequality, which tells you that $$\operatorname{Prob}(|\overline o_N-\mathbb{E}[\hat o]|\ge \epsilon) \le \frac{\sigma_O^2}{N \epsilon^2}.$$ Or in other words, this tells you that to have $\operatorname{Prob}(|\overline o_N-\mathbb{E}[\hat o]|\ge\epsilon)\le\delta$ you need $$N \ge \frac{\sigma_O^2}{\epsilon^2\delta}.$$

Another approach is to use bounds that exploit the sub-Gaussian nature of the kinds of distributions you usually get. In other words, these are bounds that work when the distribution does not have "heavy tails", i.e. "large oscillations". For example, using Hoeffding's inequality, denoting with $\overline o_N$ the standard mean of the estimates obtained from $N$ measurement rounds, $\overline o_N\equiv \frac{1}{N}\sum_{k=1}^N \hat o(b_k)$, you have $$\operatorname{Prob}(|\overline o_N-\mathbb{E}[\hat o]|\ge \epsilon) \le 2\exp\left(-\frac{2N \epsilon^2}{(\lambda_M-\lambda_m)^2}\right) \le 2\exp\left(-\frac{N \epsilon^2}{\Delta_O^2}\right),$$ where $\lambda_m,\lambda_M$ are lower and higher eigenvalues of $O$, and we used the notation $\Delta_O\equiv |\lambda_M-\lambda_m|$. Note that $\Delta_O\le 2\|O\|_\infty$. This bound ensures that you have $\operatorname{Prob}(|\overline o_N-\operatorname{tr}(O\rho)|>\epsilon) < \delta,$ in words, the error probability is larger than (some small) $\epsilon$ with (some small) probability $\delta$, if the number $N$ of measurements satisfies $$N \ge \frac{\Delta_O^2}{\epsilon^2}\log(2/\delta).$$ Note how this bound doesn't involve the variance, but rather how large the range of possible values of the estimator is, and that you get a better scaling with respect to the error probability $\delta$.

Both this and the other inequality can be used, which one is more useful depends on the context. If for example you know that the variance is very small compared to the range (eg you know the measured state is close to an eigenstate) then Chebyshev's might give you better estimates for the required $N$. You can also use more sophisticated bounds of "Hoeffding type" that involve the variance, e.g. Bernstein's. For a nice more thorough overview of some of these ideas see e.g. PRXQuantum.2.010201, and in particular Section 2C.