# How to derive the expression for the probability in quantum phase estimation? ((5.27) Nielsen & Chuang)

I'm trying to understand the QPE algorithm that is presented in the Nielsen and Chuang textbook. More precisely, I do not understand Equation $$(5.27)$$.

Context: In the following, let $$b$$ be a natural number such that $$\frac{b}{2^t}$$ denotes the best $$t$$ bit approximation of $$\varphi$$, which is less than $$\varphi$$. Next, define $$\delta := \varphi - \frac{b}{2^t}$$. At the end of the QPE, we obtain the following state:

$$\sum_{l = 0}^{2^t-1} a_l \left|(b+l) (\text{mod } 2^t)\right>$$ where $$a_l := \frac{1}{2^t}(\frac{1 - e^{i2\pi(2^t\delta - l)}}{1 - e^{i2\pi(\delta - l/2^t)}})$$. Now, if a measurement yields $$m$$, our goal is to bound the probability of $$\left|m - b\right| > e$$, where $$e$$ is a natural number. The textbook now proceeds as follows: The probability of observing such an $$m$$ is given by

$$p(\left|m - b\right| > e) = \sum_{-2^{t-1} < l \le -(e+1)} \left|a_l\right|^2 + \sum_{e + 1 \le l \le 2^{t-1}} \left|a_l\right|^2$$ (Equation $$(5.27)$$ in the Nielsen and Chuang textbook).

To be honest, I really fail to see why this holds in general. Naively, I would say, the probability is given by $$p(\left|m - b\right| > e) = \sum_{l = -b}^{-(e+1)} \left|a_l\right|^2 + \sum_{l = e + 1}^{2^t -1 - b} \left|a_l\right|^2$$ (that is, the sum over all offets $$l$$ such that $$b - l \ge e + 1$$ or $$b + l \ge e + 1$$).

For certain values of $$t, b, e$$ both equations yield different results:
For instance, consider $$t := 3$$, $$b := 6$$ and $$e := 2$$. Then,

$$p(\left|m - b\right| > e) = \sum_{-2^{3-1} < l \le -(2+1)} \left|a_l\right|^2 + \sum_{2 + 1 \le l \le 2^{3-1}} \left|a_l\right|^2 = \sum_{-4 < l \le -3} \left|a_l\right|^2 + \sum_{3 \le l \le 4} \left|a_l\right|^2 = \left|a_{-3}\right|^2 + \left|a_3\right|^2 + \left|a_4\right|^2$$ VS. $$p(\left|m - b\right| > e) = \sum_{l = -6}^{-(2+1)} \left|a_l\right|^2 + \sum_{l = 2 + 1}^{2^3 - 1 - 6} \left|a_l\right|^2 = \sum_{l = -6}^{-3} \left|a_l\right|^2 + \sum_{l = 3}^{1} \left|a_l\right|^2 = \left|a_{-6}\right|^2 + \left|a_{-5}\right|^2 + \left|a_{-4}\right|^2 + \left|a_{-3}\right|^2 + 0$$

Since $$\left|a_{4}\right|^2 = \left|a_{4 - 2^3}\right|^2 = \left|a_{-4}\right|^2$$ and $$\left|a_{3}\right|^2 = \left|a_{3-2^3}\right|^2 = \left|a_{-5}\right|^2$$, only the $$\left|a_{-6}\right|^2$$ summand is missing.

Can you give me a hint why Equation (5.27) is correct and my naive approach is wrong?

You observe a value $$m=l+b\text{ mod }2^t$$ (I'm not sure if you've directly equated $$m$$ with $$l$$, which is why you still have a $$b$$ floating around). So, $$m-b\equiv l\text{ mod }2^t.$$ Hence, we're interested in $$|l|>e\text{ mod }2^t.$$ It's probably easiest to think about $$l$$ being in the range $$-2^{t-1}+1$$ to $$2^{t-1}$$, which of course is equivalent under the modulus. The range that we're excluding is then $$l=-e$$ to $$e$$. Thus, we're keeping 2 ranges: $$-2^{t-1}+1$$ to $$-e-1$$ and $$e+1$$ to $$2^{t-1}$$.

• Thanks a lot for your answer! I still have to think about it (btw, I just assume that you meant $2^{t}$, $2^{t-1}$ instead of $2^l$, $2^{l-1}$). What I don't quite understand is $\left|l\right| > e \mod 2^t$ In my example with $t := 3$ and $e := 2$, we have for $l \in \{3, 4, 5, 6, 7\}$ that $\left|l\right| > 2 \mod 2^3$ and, hence, we would end up with $\left|a_3\right|^2 + \left|a_4\right|^2 + \left|a_5\right|^2 + \left|a_6\right|^2 + \left|a_7\right|^2$ (which is not what Equation (5.27) would give us). Note that at this point I still assume that $l \in \{0, \ldots, 2^t - 1\}$. Aug 23, 2023 at 20:20
• The range is important here. The notation I've used only makes proper mathematical sense for the range I've used. Otherwise, we have two conditions: $l>e\text{ mod }2^t$ and $-l>e\text{ mod }2^t$. That second condition failing is what cuts out the $l=6,7$ terms. Aug 24, 2023 at 1:42
• Ah, I got the distance $\text{mod } 2^t$ wrong. For $l \in \{0, \ldots, 2^t - 1\}$, define $\left|l\right|_{\text{mod } 2^t} := \min \{ l \text{ mod } 2^t, -l \text{ mod } 2^t\}$. That is, we are interested in all $m = b + l$ such that $\left|m - b\right|_{\text{mod } 2^t} = \left|b + l - b\right|_{\text{mod } 2^t} = \left| l \right|_{\text{mod } 2^t} > e$. I'll come up with a small write-up that puts everything together within the next few days. Thanks a lot for your help! Aug 24, 2023 at 18:06

I think it's sometimes useful to see the same thing derived with a (slightly) different notation. The more direct answer to the original question is in the last paragraph here: one has to realise that $$\ell$$ is a way to reparametrise the qubit states, and that its range is $$-N/2\le \ell\le N/2$$, $$N\equiv 2^n$$. This simplifies considering indices $$k$$ that are "effectively" further away from the max by some quantity $$e$$, which is otherwise nontrivial due to the periodicity of $$a_k$$ with respect to $$k$$.

Brief overview of QPE — As a brief reminder, doing QPE for a unitary $$U$$ with eigenvector $$|\psi\rangle$$ such that $$U|\psi\rangle=e^{2\pi i\varphi}|\psi\rangle$$ amounts to a scheme to apply the unitary $$\sum_x \mathbb{P}_x\otimes U^x$$, $$\mathbb{P}_x\equiv|x\rangle\!\langle x|$$, to the input state $$|+\rangle\otimes|\psi\rangle$$ with $$|+\rangle\equiv N^{-1/2}\sum_{x=0}^{N-1}|x\rangle$$, where $$N\equiv 2^n$$ and $$n$$ is the number of qubits (in the first register). This gives you on the first register the state $$N^{-1/2}\sum_x e^{2\pi i x \varphi}|x\rangle$$, and you then apply the (inverse) QFT to it to get $$N^{-1}\sum_{k=0}^{N-1}\underbrace{\left(\sum_{x=0}^{N-1}e^{2\pi ix(\varphi-k/N)}\right)}_{\equiv a_k}|k\rangle.$$ Simple algebra now leads to the expression $$a_k = \frac{1- e^{2\pi i(\varphi-k/N)N}}{1-e^{2\pi i(\varphi-k/N)}} \implies |a_k| = \left\lvert\frac{\sin[\pi(N\varphi-k)]}{\sin[\pi(N\varphi-k)/N]}\right\rvert.$$ Remember that here $$k,N\in\mathbb{N}$$ and $$0\le k\le N-1$$. On the other hand $$\varphi\in\mathbb{R}$$ with $$\varphi\in[0,1]$$. Also, by definition of $$a_k$$ we must have $$\sum_k |a_k|^2 = N^2$$. Note that the periodic nature of $$\sin$$ is such that if $$N\varphi$$ is an integer, then $$a_{N\varphi}=N$$ and all other coefficients vanish.

Finding the most likely outcome — The first question is: what's the max of $$a_k$$? Ideally this would correspond to $$k=N\varphi$$, for which we'd have $$|a_k|=N$$ and thus $$a_j=0$$ for all $$j\neq k$$ (technically this expression for $$a_k$$ is not defined for $$k=N\varphi$$ but you can easily see this as a removable singularity via Taylor expansion). However, in general, $$\varphi$$ is not even rational and therefore $$k=N\varphi$$ is not possible. A more general observation is that the function $$\sin(\pi y)/\sin(\pi y/N)$$ always achieves its maximum at $$y=0$$ in the domain $$y\in[-N,N]$$, so the max $$|a_k|$$ will correspond to the $$k$$ that minimises the absolute value of $$N\varphi-k$$.

Suppose then that $$\varphi=0.\varphi_1\varphi_2\cdots$$, which is concise notation to say $$\varphi=\sum_{k\ge1} 2^{-k}\varphi_k$$. It follows that $$N\varphi=\varphi_1\cdots\varphi_n.\varphi_{n+1}\cdots$$, that is, the integer part of $$N\varphi$$ is $$\varphi_1\cdots\varphi_n\equiv\sum_{k=1}^n 2^{n-k}\varphi_k$$, where $$\varphi_k\in\{0,1\}$$. Given the constraint $$k\in\mathbb{N}$$, it's clear that the min of $$|N\varphi-k|$$ corresponds to either $$k=\lfloor N\varphi\rfloor\equiv\varphi_1\cdots\varphi_n$$, where I denoted with $$\lfloor N\varphi\rfloor$$ the integer part of $$N\varphi$$, or $$k=\lfloor N\varphi\rfloor+1$$, with the caveat of remembering that $$k$$ should be understood mod $$N$$, it being a label for an $$n$$-qubit state. For example, if $$\varphi=0.1111=1/2+1/4+1/8+1/16=15/16$$ and $$N=2^2$$, then $$N\varphi=11.11=3+3/4$$, $$\lfloor N\varphi\rfloor=3$$, but the best integer approximation is $$N\varphi \simeq 4\equiv 0\pmod 4$$. So more generally, let's denote with $$\tilde\varphi$$ the best $$n$$-bit integer approximation to $$N\varphi$$, so that $$|\tilde\varphi-N\varphi|\le 1/2$$.

To gain a more concrete idea of what we're talking about, let's consider as an explicit example the case $$N=2^4$$, $$\varphi=0.0100101\equiv 37/128$$. Then $$N \varphi=100.101=4.624$$ and thus $$\tilde\varphi=2^2+1=5$$. The thus expect $$|a_k|$$ to be maximal at $$k=5$$, and hopefully decay fast far from this value. In fact, this is exactly what we find plotting $$|a_k|^2/N^2$$, as shown in the following:

Here we see clearly how $$|a_k|^2/N^2$$ is largest at $$k=5$$ and decays fast moving away from it in any direction. The purpose of most of the following will be to show analytically that the probability of getting values far from the max is sufficiently low.

Upper bound the error probability — Here's where things get tricky, in part due to the periodic nature of the quantities we're dealing with. For starters, note that $$a_k$$ is compatible with the periodic nature of $$k$$, as $$a_{k+jN}=a_k$$ for any $$j\in\mathbb{Z}$$. Furthermore, observe the following general upper bound to the amplitudes: $$|a_k| \le \frac{N}{2|N\varphi-k|},$$ where I used $$\lvert \sin(x)\rvert\le1$$ and $$\lvert\sin(x)\rvert\ge 2|x|/\pi$$ for $$|x|\le\pi/2$$, and this assumption translates to $$|N\varphi-k|\le N/2$$ in our case. Applying this bound to our case is however trickier than it seems, because $$k$$ is defined mod $$N$$, and therefore rather than $$N\varphi-k$$ we would actually want there $$N\varphi-\tilde k$$ where $$\tilde k\equiv k\pmod N$$ is the integer such that $$|N\varphi-\tilde k|$$ is closest to $$0$$. As an example, imagine $$N=2^2$$, $$\varphi=0.0001$$, $$N\varphi=0.01=1/4$$, $$\tilde\varphi=0$$, and $$k=3$$. Then $$|N\varphi-k|=|3-1/4|>N/2$$, which means our approximation would seem to not work. However, due to the periodicity of the $$\sin$$ we're trying to approximate, we can arbitrarily translate $$k$$ by a multiple of $$N$$ without changing its value. In this example, this means we can replace $$k=3$$ with $$\tilde k=k-4=-1$$. Now $$|N\varphi-\tilde k|, and we thus have the upper bound $$\frac{1}{\lvert\sin(\pi(N\varphi-\tilde k)/N)\rvert }\le \frac{N}{2\lvert N\varphi-\tilde k\vert}.$$

Reparametrise indices to simplify matters — The above observation makes things messy: we'd have to take care when doing the upper bound to use $$\tilde k$$ instead of $$k$$ every time, but then it's not trivial to write concisely the expression for $$\tilde k$$. The trick employed in N&C to avoid this is to reparametrise indices, effectively "centering" our "coordinate system" around $$k=\tilde\varphi$$. To this end, define $$-N/2 \le \ell\le N/2$$ such that $$k=(\tilde\varphi+\ell)\bmod N$$. Now $$\ell$$ tells us how far we are from the true value, correctly taking into account the periodicity of the sine. In particular, we have $$\lvert a_k\rvert \le \frac{1}{\lvert\sin[\pi(N\varphi - \tilde\varphi-\ell)/N]\rvert} \le \frac{N}{2\lvert N\varphi - \tilde\varphi-\ell \rvert },$$ where we didn't need to keep writing the mod $$N$$ part in the sin, because the sin is intrinsically periodic. However, we can safely say that $$(N\varphi - \tilde\varphi-\ell) \le N/2$$, and therefore the upper bound on $$|a_k|$$ works without issues.

Actually, there is another caveat in the above argument that makes it slightly incorrect. You might notice that $$\ell=-N/2$$ and $$\ell=N/2$$ actually correspond to the same value of $$k$$, and therefore the range of $$\ell$$ doesn't include one of the two extremes. Which extreme is to be excluded, however, depends on whether $$\tilde\varphi\ge N\varphi$$ or $$\tilde\varphi\le N\varphi$$, as only one of these will give $$(N\varphi - \tilde\varphi-\ell) \le N/2$$. As always, an explicit example might help in understanding things. Say we have $$N=2^2$$, $$N\varphi=0.01\equiv 1/4$$, and thus $$\tilde\varphi=0$$. This gives the following correspondence between $$k$$ and $$\ell$$: $$k=0\simeq \ell=0, \qquad k=1 \simeq \ell=1, \qquad k=2 \simeq \ell=\pm2, \qquad k=3 \simeq \ell=-1,$$ where clearly we can chose either $$\ell=2$$ or $$\ell=-2$$ to correspond to $$k=2$$. But only one of these values works, because $$|N\varphi-\tilde\varphi - 2| = |1/4-2| < 2, \\ |N\varphi-\tilde\varphi - (-2)| = |1/4+2| > 2.$$ In conclusion, the range of $$\ell$$ is either $$-N/2 < \ell \le N/2$$ or $$-N/2\le \ell < N/2$$, depending on whether $$\tilde\varphi\le N\varphi$$ or $$\tilde\varphi\ge N\varphi$$.

Get the final upper bound — Having worked out the above boring technicalities, let's try and derive the final upper bound. We want to upper bound the probability of getting an outcome $$k$$ that has distance $$e\in\mathbb{N}$$ from the best value $$\tilde\varphi$$, that is, $$p(|k-\tilde\varphi|>e)$$. This is by definition equal to the sum of $$|a_k|^2/N^2$$ for all $$k$$ such that $$|k-\tilde\varphi|>e$$, except, again, we don't really care about $$|k-\tilde\varphi|$$, but rather about $$|\tilde k-\tilde\varphi|$$ with $$\tilde k$$ the integer with $$\tilde k\equiv k \pmod N$$ that is closest to $$\tilde\varphi$$ taking into account the mod $$N$$ periodicity. This statement is again made much simpler using our reparametrisation with $$\ell$$, because it amounts to just looking at the indices corresponding to $$|\ell|>e$$. We thus have $$p(|k-\tilde\varphi|>e) = \sum_{k: \,\, |\ell|>e} \frac{\lvert a_{k} \rvert^2}{N^2} \le \frac{1}{4}\sum_{\ell: \,\, |\ell|>e} \frac{1}{\lvert N\varphi-\tilde\varphi-\ell \rvert^2} \\ \le \frac{1}{4}\left[\sum_{\ell=-N/2}^{-e-1} \frac{1}{\lvert N\varphi-\tilde\varphi-\ell \rvert^2} + \sum_{\ell=e+1}^{N/2} \frac{1}{\lvert N\varphi-\tilde\varphi-\ell \rvert^2}\right],$$ where in the last step I willfully neglected the issue with the range of $$\ell$$ not including one of $$-N/2$$ or $$N/2$$, as we can still get a decent upper bound including both.

Because by definition of $$\tilde\varphi$$ we have $$|N\varphi-\tilde\varphi|\le1/2$$, then for $$\ell<0$$ we can write $$|N\varphi-\tilde\varphi-\ell|\ge |-\ell-1/2|$$, while for $$\ell>0$$ we can write $$|N\varphi-\tilde\varphi-\ell|\ge |-\ell+1/2|$$. Therefore $$p(|k-\tilde\varphi|>e) \le \frac{1}{4}\left[\sum_{\ell=-N/2}^{-e-1} \frac{1}{( \ell+1/2 )^2} + \sum_{\ell=e+1}^{N/2} \frac{1}{( \ell-1/2 )^2}\right] \\ \le \frac14\left[ \sum_{\ell=-N/2}^{-e-1}\frac{1}{(\ell+1)^2} + \sum_{\ell=e+1}^{N/2}\frac{1}{(\ell-1)^2} \right] = \frac12 \sum_{\ell=e}^{N/2-1}\frac{1}{\ell^2} \le \frac{1}{2(e-1)}.$$ In the last step I used the same way to approximate sums with integral employed in N&C (interesting, although I used a slightly different procedure, using $$\tilde\varphi$$ rather than just approximating $$N\varphi$$ from below, and ignoring the associated issue with the range of $$\ell$$, I ended up with the same identical final result).