# Random Variables & Expectation > Every number you ever read off a quantum computer is an *estimate* of an expectation value from a > finite number of repetitions. This lesson builds the machinery — expectation, variance, > concentration — that tells you exactly how good that estimate is, and how many shots you must pay > for. ## Learning Objectives By the end of this lesson you will be able to: 1. **Define** discrete and continuous random variables and their PMF/PDF/CDF, and **compute** expectation and variance from each. 2. **Use** linearity of expectation and the variance of a sum to analyze estimators. 3. **Identify** and **parameterize** the Bernoulli, Binomial, Geometric, and Gaussian distributions and state their means and variances. 4. **State** the Law of Large Numbers and the Central Limit Theorem and **explain** the $1/\sqrt{M}$ standard-error scaling of a sample mean. 5. **Apply** Markov, Chebyshev, and Hoeffding inequalities to bound deviation probabilities and **derive** a shot count for a target confidence. 6. **Simulate** sample-mean convergence in NumPy and **verify** the $1/\sqrt{M}$ error scaling empirically. --- ## Intuition A **random variable** is a number attached to each outcome of an experiment — the *measurement* of the random world. Roll two dice; the sum is a random variable. Run a circuit and assign $+1$ to outcome $0$ and $-1$ to outcome $1$; that $\pm1$ reading is a random variable, and its **expectation value** is precisely $\langle Z\rangle$, the quantity variational algorithms minimize. Two questions dominate everything downstream: - **Where is it centered?** The **expectation** $\mathbb E[X]$ — the long-run average. - **How much does it wobble?** The **variance** $\mathrm{Var}(X)$ — the spread around the center. The punchline of the whole lesson is this. You cannot read $\mathbb E[X]$ off a quantum computer directly; you can only *sample* $X$ a finite number $M$ of times ("shots") and average. That sample average is itself random. The Law of Large Numbers says it converges to $\mathbb E[X]$; the Central Limit Theorem says its fluctuations shrink like $1/\sqrt M$; and concentration inequalities (Hoeffding) turn that into a hard promise: *"with probability $\ge 1-\delta$, my estimate is within $\varepsilon$ of the truth, provided $M \ge \dots$"*. That inequality is, quite literally, the cost model for running on a QPU. --- ## Theory ### Random variables, PMF, PDF, CDF Given a probability space $(\Omega, \mathcal F, \mathbb P)$, a **random variable** is a (measurable) function $X : \Omega \to \mathbb R$. Measurability means $\{\omega : X(\omega) \le x\} \in \mathcal F$ for every $x$, so we may assign it a probability. The **cumulative distribution function (CDF)** is $$ F_X(x) := \mathbb P(X \le x), $$ non-decreasing, right-continuous, with $F_X(-\infty)=0$, $F_X(+\infty)=1$. A **discrete** RV takes values in a countable set; its **probability mass function (PMF)** is $p_X(x) = \mathbb P(X = x)$, with $p_X(x) \ge 0$ and $\sum_x p_X(x) = 1$. A **continuous** RV has a **probability density function (PDF)** $f_X \ge 0$ with $\int_{-\infty}^{\infty} f_X(x)\,dx = 1$ and $\mathbb P(X \in [a,b]) = \int_a^b f_X(x)\,dx$. Note $f_X(x)$ is a *density*, not a probability — it can exceed $1$; only its integral over a set is a probability. The CDF and PDF are linked by $f_X = F_X'$ (where the derivative exists). ### Expectation The **expectation** (mean) of $X$ is its probability-weighted average: $$ \mathbb E[X] = \sum_x x\,p_X(x) \quad\text{(discrete)}, \qquad \mathbb E[X] = \int_{-\infty}^{\infty} x\,f_X(x)\,dx \quad\text{(continuous)}, $$ when the sum/integral converges absolutely. For any (measurable) function $g$, the **law of the unconscious statistician** gives $\mathbb E[g(X)] = \sum_x g(x)p_X(x)$ — you do *not* need the distribution of $g(X)$, just that of $X$. > **Quantum tie-in.** For a Pauli-$Z$ measurement with outcomes $\pm1$, $\langle Z\rangle = \mathbb > E[\text{reading}] = (+1)p_0 + (-1)p_1 = p_0 - p_1$. Estimating $\langle Z\rangle$ on hardware *is* > estimating an expectation from samples. We meet $\langle A\rangle = \langle\psi|A|\psi\rangle$ > formally in [Term 1.3](../../term-1-quantum-mechanics/c03-measurement/02-expectation-uncertainty.md). ### Linearity of expectation For *any* random variables $X, Y$ (independent or not) and constants $a, b$, $$ \mathbb E[aX + bY] = a\,\mathbb E[X] + b\,\mathbb E[Y]. $$ **Derivation (discrete, finite $\Omega$).** Expectation can be written as a sum over outcomes, $\mathbb E[X] = \sum_{\omega} X(\omega)\,\mathbb P(\{\omega\})$. Then $$ \mathbb E[aX+bY] = \sum_\omega \big(aX(\omega)+bY(\omega)\big)\mathbb P(\{\omega\}) = a\sum_\omega X(\omega)\mathbb P(\{\omega\}) + b\sum_\omega Y(\omega)\mathbb P(\{\omega\}) = a\mathbb E[X] + b\mathbb E[Y]. $$ The power of linearity is that it needs **no independence** — it is the single most useful identity in applied probability. ### Variance and standard deviation The **variance** measures spread about the mean: $$ \mathrm{Var}(X) := \mathbb E\big[(X - \mathbb E[X])^2\big] = \mathbb E[X^2] - (\mathbb E[X])^2, $$ where the second form follows by expanding the square and using linearity: $\mathbb E[(X-\mu)^2] = \mathbb E[X^2] - 2\mu\mathbb E[X] + \mu^2 = \mathbb E[X^2] - \mu^2$ with $\mu = \mathbb E[X]$. The **standard deviation** is $\sigma_X = \sqrt{\mathrm{Var}(X)}$, in the same units as $X$. Scaling: $\mathrm{Var}(aX+b) = a^2\,\mathrm{Var}(X)$ (shifts don't change spread; scales square). For the sum of random variables, $$ \mathrm{Var}\Big(\textstyle\sum_i X_i\Big) = \sum_i \mathrm{Var}(X_i) + \sum_{i\ne j}\mathrm{Cov}(X_i, X_j), \qquad \mathrm{Cov}(X,Y) := \mathbb E[(X-\mathbb E X)(Y-\mathbb E Y)]. $$ When the $X_i$ are **pairwise uncorrelated** (in particular, independent), the covariance terms vanish and **variances add**: $\mathrm{Var}(\sum_i X_i) = \sum_i \mathrm{Var}(X_i)$. This single fact drives the $1/\sqrt M$ scaling below. ### Key distributions | Distribution | PMF / PDF | Mean | Variance | Models | |---|---|---|---|---| | **Bernoulli$(p)$** | $p_X(1)=p,\ p_X(0)=1-p$ | $p$ | $p(1-p)$ | one shot, one bit | | **Binomial$(n,p)$** | $\binom{n}{k}p^k(1-p)^{n-k}$ | $np$ | $np(1-p)$ | # of $1$s in $n$ i.i.d. shots | | **Geometric$(p)$** | $(1-p)^{k-1}p,\ k\ge 1$ | $1/p$ | $(1-p)/p^2$ | # shots until first success | | **Gaussian$(\mu,\sigma^2)$** | $\frac{1}{\sqrt{2\pi}\sigma}e^{-(x-\mu)^2/2\sigma^2}$ | $\mu$ | $\sigma^2$ | CLT limit; measurement noise | **Bernoulli variance derivation.** $X^2 = X$ for a $0/1$ variable, so $\mathbb E[X^2] = \mathbb E[X] = p$ and $\mathrm{Var}(X) = p - p^2 = p(1-p)$, maximized at $p=\tfrac12$ (a fair coin / an equal superposition is the *noisiest* single bit). A Binomial$(n,p)$ is a sum of $n$ independent Bernoulli$(p)$, so by linearity and additivity of variance its mean is $np$ and variance $np(1-p)$ — no binomial-coefficient sums required. ### Law of Large Numbers (LLN) Let $X_1, X_2, \dots$ be i.i.d. with mean $\mu = \mathbb E[X_i]$ and finite variance $\sigma^2$. Define the **sample mean** $\bar X_M = \frac1M\sum_{i=1}^M X_i$. The **(weak) law of large numbers** states that for every $\varepsilon > 0$, $$ \mathbb P\big(|\bar X_M - \mu| \ge \varepsilon\big) \xrightarrow[M\to\infty]{} 0 . $$ **Derivation via Chebyshev.** By linearity, $\mathbb E[\bar X_M] = \mu$. By additivity of variance for independent terms, $$ \mathrm{Var}(\bar X_M) = \frac{1}{M^2}\sum_{i=1}^M \mathrm{Var}(X_i) = \frac{\sigma^2}{M}. $$ Chebyshev's inequality (below) gives $\mathbb P(|\bar X_M - \mu| \ge \varepsilon) \le \sigma^2/(M\varepsilon^2) \to 0$. $\blacksquare$ This is *why* averaging more shots works. ### The $1/\sqrt M$ standard error — derived The line $\mathrm{Var}(\bar X_M) = \sigma^2/M$ is the most important formula in this whole course for practical work. The **standard error** of the sample mean — the standard deviation of the estimator — is its square root: $$ \boxed{\ \mathrm{SE}(\bar X_M) = \sqrt{\mathrm{Var}(\bar X_M)} = \frac{\sigma}{\sqrt M}\ } $$ So to **halve** your statistical error you must **quadruple** your shots. To estimate a Pauli expectation $\langle Z\rangle$ with outcomes $\pm1$, the per-shot variance is $\sigma^2 = 1 - \langle Z\rangle^2 \le 1$, so $\mathrm{SE} \le 1/\sqrt M$. This $1/\sqrt M$ wall is fundamental: it governs runtime and dollar cost of *every* expectation-value estimation on a QPU, and it is the reason variational algorithms are shot-hungry. ### Central Limit Theorem (CLT) The LLN says $\bar X_M \to \mu$; the CLT says *how* the fluctuations are shaped. For i.i.d. $X_i$ with mean $\mu$ and finite variance $\sigma^2$, the standardized sample mean converges *in distribution* to a standard Gaussian: $$ \frac{\bar X_M - \mu}{\sigma/\sqrt M} \;\xrightarrow[M\to\infty]{d}\; \mathcal N(0,1). $$ Equivalently, $\bar X_M \approx \mathcal N(\mu, \sigma^2/M)$ for large $M$. Two consequences we use constantly: (i) the error bars on a sampled expectation are Gaussian, so a $\pm 1.96\,\sigma/\sqrt M$ interval is an approximate $95\%$ confidence interval; (ii) the *width* of that interval is $\propto 1/\sqrt M$, consistent with the standard error above. ### Concentration inequalities These bound the probability of large deviations and let us choose $M$ *before* running. **Markov's inequality.** For a non-negative RV $Y \ge 0$ and $a > 0$, $$ \mathbb P(Y \ge a) \le \frac{\mathbb E[Y]}{a}. $$ *Derivation.* $a\,\mathbf 1\{Y\ge a\} \le Y$ pointwise (check both cases $Y\ge a$ and $Y = p0 - p1 = 0.4 per_shot_var = 1 - true_expectation**2 # = 1 - 0.16 = 0.84 def sample_readings(M, rng, p0): """Return M samples of the +/-1 reading for a qubit with P(0)=p0.""" bits = rng.random(M) < p0 # True (outcome 0) with prob p0 return np.where(bits, +1.0, -1.0) # map 0 -> +1, 1 -> -1 # Sample-mean convergence: estimate at increasing shot counts. print(f"true = {true_expectation:.4f}\n") print(f"{'M':>9} {'estimate':>10} {'|error|':>10} {'SE = sqrt(var/M)':>18}") for M in (10, 100, 1_000, 10_000, 100_000, 1_000_000): est = sample_readings(M, rng, p0).mean() se = np.sqrt(per_shot_var / M) # predicted standard error print(f"{M:>9} {est:>10.4f} {abs(est-true_expectation):>10.4f} {se:>18.5f}") ``` ```python # Verify the 1/sqrt(M) scaling: estimate the EMPIRICAL standard error by repeating # the whole experiment many times at each M and taking the std of the estimates. # A clean 1/sqrt(M) law means error * sqrt(M) is (approximately) constant. Ms = np.array([16, 64, 256, 1024, 4096, 16384]) trials = 2000 emp_se = np.empty(len(Ms)) for i, M in enumerate(Ms): estimates = np.array([sample_readings(M, rng, p0).mean() for _ in range(trials)]) emp_se[i] = estimates.std(ddof=1) # spread of the estimator theory_se = np.sqrt(per_shot_var / Ms) print(f"\n{'M':>7} {'empirical SE':>14} {'theory SE':>12} {'SE*sqrt(M)':>12}") for M, e, t in zip(Ms, emp_se, theory_se): print(f"{M:>7} {e:>14.5f} {t:>12.5f} {e*np.sqrt(M):>12.5f}") # The last column ~ sqrt(per_shot_var) = sqrt(0.84) ~ 0.9165 for every M, # confirming SE ∝ 1/sqrt(M). print(f"\nsqrt(per-shot var) = {np.sqrt(per_shot_var):.5f} (the 1/sqrt(M) constant)") ``` **Plot description.** Plot `emp_se` against `Ms` on **log–log axes**: the points fall on a straight line of slope $-\tfrac12$, the signature of $1/\sqrt M$ decay. Overlaying `theory_se` gives a line the empirical points sit on. (To render: `plt.loglog(Ms, emp_se, 'o'); plt.loglog(Ms, theory_se, '-'); plt.xlabel("shots M"); plt.ylabel("standard error")`.) Equivalently, the rightmost column of the table is flat at $\sqrt{0.84}\approx 0.9165$ — error times $\sqrt M$ is constant, which is the $1/\sqrt M$ law restated. ```python # Sanity-check the Hoeffding shot count from Example 2. def hoeffding_shots(eps, delta, R=2.0): return (R**2) / (2 * eps**2) * np.log(2 / delta) print(f"\nShots for eps=0.02, delta=0.05 (95% conf): {hoeffding_shots(0.02, 0.05):.0f}") print(f"Shots for eps=0.01, delta=0.01 (99% conf): {hoeffding_shots(0.01, 0.01):.0f}") ``` Expected output (abridged; reproducible under the seed): ```text true = 0.4000 M estimate |error| SE = sqrt(var/M) 10 0.0000 0.4000 0.28983 100 0.4600 0.0600 0.09165 1000 0.3900 0.0100 0.02898 10000 0.4080 0.0080 0.00917 100000 0.3977 0.0023 0.00290 1000000 0.4005 0.0005 0.00092 M empirical SE theory SE SE*sqrt(M) 16 0.22689 0.22913 0.90756 64 0.11403 0.11456 0.91227 256 0.05776 0.05728 0.92421 1024 0.02923 0.02864 0.93545 4096 0.01445 0.01432 0.92453 16384 0.00692 0.00716 0.88535 sqrt(per-shot var) = 0.91652 (the 1/sqrt(M) constant) Shots for eps=0.02, delta=0.05 (95% conf): 18444 Shots for eps=0.01, delta=0.01 (99% conf): 105966 ``` This is the simulated version of running a circuit with `shots=M` on Braket and computing an expectation value from the counts — the convergence you see here is the convergence you pay for there. --- ## Exercises **1. (Easy) Mean and variance of a die.** Let $X$ be the result of a fair six-sided die. Compute $\mathbb E[X]$ and $\mathrm{Var}(X)$.
Solution $\mathbb E[X] = \frac16(1+2+3+4+5+6) = \frac{21}{6} = 3.5$. $\mathbb E[X^2] = \frac16(1+4+9+16+25+36) = \frac{91}{6} \approx 15.1\overline{6}$. $\mathrm{Var}(X) = \frac{91}{6} - 3.5^2 = 15.1\overline{6} - 12.25 = \frac{35}{12} \approx 2.9167$.
**2. (Easy) Binomial via linearity.** A circuit is run $n=1000$ times; each shot reads $1$ with probability $p=0.3$ independently. Let $K$ be the number of $1$s. Give $\mathbb E[K]$ and $\mathrm{Var}(K)$ without summing any binomial coefficients.
Solution $K = \sum_{i=1}^{n} X_i$ with $X_i \sim \text{Bernoulli}(p)$ independent. By linearity $\mathbb E[K] = np = 300$. By additivity of variance (independence) $\mathrm{Var}(K) = np(1-p) = 1000(0.3)(0.7) = 210$, so $\sigma_K = \sqrt{210}\approx 14.5$.
**3. (Medium) Geometric mean.** Shots are repeated until the first $1$ appears; each shot is $1$ with probability $p$. Let $N$ be the shot index of the first $1$. Show $\mathbb E[N] = 1/p$.
Solution $N\sim$ Geometric$(p)$ with $\mathbb P(N=k) = (1-p)^{k-1}p$. Let $q=1-p$. Then $$ \mathbb E[N] = \sum_{k=1}^{\infty} k\,q^{k-1}p = p\sum_{k=1}^{\infty} k q^{k-1} = p\cdot\frac{1}{(1-q)^2} = p\cdot\frac{1}{p^2} = \frac1p, $$ using $\sum_{k\ge1} k q^{k-1} = \frac{d}{dq}\sum_{k\ge0}q^k = \frac{d}{dq}\frac{1}{1-q} = \frac{1}{(1-q)^2}$ for $|q|<1$. So if a "success" has probability $p$, you wait $1/p$ shots on average — e.g. amplitude $\sqrt p$ needs $\sim 1/p$ repetitions, the intuition behind why Grover's $O(1/\sqrt N)$ is a real speedup over the classical $O(1/p)=O(N)$.
**4. (Medium) Chebyshev vs. Hoeffding shot counts.** You estimate the mean of a $[0,1]$-bounded reading to within $\varepsilon = 0.05$ with failure probability $\delta = 0.05$. (a) Use Chebyshev with the worst-case variance for a $[0,1]$ variable to bound $M$. (b) Use Hoeffding ($R=1$). Compare.
Solution (a) Worst-case variance of a $[0,1]$ variable is $\tfrac14$ (Bernoulli at $p=\tfrac12$). For the sample mean, $\mathrm{Var}(\bar X_M) = \sigma^2/M \le 1/(4M)$. Chebyshev: $\mathbb P(|\bar X_M-\mu|\ge\varepsilon) \le \frac{1}{4M\varepsilon^2} \le \delta$ gives $M \ge \frac{1}{4\delta\varepsilon^2} = \frac{1}{4(0.05)(0.0025)} = 2000$. (b) Hoeffding with $R=1$: $M \ge \frac{R^2}{2\varepsilon^2}\ln\frac{2}{\delta} = \frac{1}{2(0.0025)}\ln(40) = 200\times3.689 \approx 738$. Hoeffding needs $\sim 738$ vs. Chebyshev's $2000$ — for *bounded* variables the exponential bound is much tighter (here $\sim 2.7\times$ fewer shots). Chebyshev wins only when you *can't* bound the range.
**5. (Hard) Variance reduction by negation/antithetic pairing.** Suppose you can sample a reading $X\in[-1,1]$ and also its "mirror" $X'$ with the *same* mean but $\mathrm{Cov}(X, X') < 0$. Show the average $\frac{X+X'}{2}$ has variance strictly less than $\frac12\mathrm{Var}(X)$, and explain why negatively correlated samples beat independent ones for fixed shot count.
Solution With $\mathrm{Var}(X)=\mathrm{Var}(X')=\sigma^2$ and $c = \mathrm{Cov}(X,X')$, $$ \mathrm{Var}\!\Big(\frac{X+X'}{2}\Big) = \frac14\big(\mathrm{Var}(X)+\mathrm{Var}(X')+2\,\mathrm{Cov}(X,X')\big) = \frac14(2\sigma^2 + 2c) = \frac{\sigma^2 + c}{2}. $$ If $c < 0$ this is $< \sigma^2/2$, which is what *independent* pairs would give ($c=0$). Negatively correlated samples cancel each other's fluctuations, so you reach a target precision in fewer total samples. This is the principle behind **antithetic variates** and, in quantum estimation, behind clever measurement groupings that lower the effective shot count for a fixed precision.
**6. (Hard) CLT confidence interval coverage.** Simulate (in NumPy) the estimation of $\langle Z\rangle$ for $p_0 = 0.6$ using $M = 500$ shots, repeated $10{,}000$ times. For each repetition form the $95\%$ CLT interval $\bar X \pm 1.96\,\hat\sigma/\sqrt M$ (using the *sample* std $\hat\sigma$). What fraction of intervals contain the true $\langle Z\rangle = 0.2$? Should it be near $0.95$?
Solution ```python import numpy as np rng = np.random.default_rng(7) p0, M, reps = 0.6, 500, 10_000 true_Z = 2*p0 - 1 # = 0.2 covered = 0 for _ in range(reps): x = np.where(rng.random(M) < p0, 1.0, -1.0) xbar, s = x.mean(), x.std(ddof=1) half = 1.96 * s / np.sqrt(M) if xbar - half <= true_Z <= xbar + half: covered += 1 print(covered / reps) # ~0.95 ``` The coverage is close to $0.95$ because, by the CLT, $\bar X$ is approximately $\mathcal N(\mu, \sigma^2/M)$ at $M=500$, and the $\pm1.96\sigma$ band of a Gaussian captures $95\%$ of its mass. Slight under-coverage can appear because we plug in the *estimated* $\hat\sigma$; for small $M$ a Student-$t$ critical value would be more honest, but at $M=500$ the difference is negligible. This is exactly how you'd report error bars on a measured expectation value.
--- ## Checkpoint 1. State linearity of expectation. Does it require independence? What about additivity of variance? 2. Write the standard error of a sample mean of $M$ i.i.d. samples, and explain in one sentence why halving the error costs $4\times$ the shots. 3. For a $\pm1$ Pauli reading, what is the per-shot variance in terms of $\langle Z\rangle$, and which state maximizes it? 4. Contrast what the LLN and the CLT each tell you about $\bar X_M$. 5. When is Hoeffding's bound preferable to Chebyshev's, and why is it so much tighter there? 6. Give the Hoeffding shot count $M$ to reach precision $\varepsilon$ with confidence $1-\delta$ for a $\pm1$ observable.
Answers 1. $\mathbb E[aX+bY]=a\mathbb E[X]+b\mathbb E[Y]$, for *any* $X,Y$ — **no** independence needed. Additivity of variance, $\mathrm{Var}(\sum X_i)=\sum\mathrm{Var}(X_i)$, *does* require the terms to be uncorrelated (independence suffices); otherwise covariance terms appear. 2. $\mathrm{SE} = \sigma/\sqrt M$. Error $\propto M^{-1/2}$, so to cut it in half you need $(2)^2 = 4$ times as many samples. 3. $\mathrm{Var}(z) = 1 - \langle Z\rangle^2$, maximized ($=1$) at $\langle Z\rangle = 0$, i.e. the equal-superposition state $p_0=p_1=\tfrac12$. 4. LLN: $\bar X_M \to \mu$ (the estimate converges to the truth). CLT: the *fluctuations* of $\bar X_M$ around $\mu$ are approximately Gaussian with std $\sigma/\sqrt M$ — it shapes the error bars. 5. When the variables are *bounded*. Hoeffding then gives an *exponential* (in $M\varepsilon^2$) tail bound vs. Chebyshev's merely polynomial $1/(M\varepsilon^2)$, so it needs far fewer samples for a target $\delta$. Chebyshev only needs finite variance, so it wins for unbounded variables. 6. $M \ge \dfrac{R^2}{2\varepsilon^2}\ln\dfrac{2}{\delta}$ with $R=2$, i.e. $M \ge \dfrac{2}{\varepsilon^2}\ln\dfrac{2}{\delta}$.
--- ## Further Reading - **[CT, §2.1, 3.1]** Cover & Thomas — expectation, AEP, and the LLN as a warm-up to entropy. - **[NC, §2.2.5–2.2.6, Box 2.1]** Nielsen & Chuang — expectation values of observables and measurement statistics. - **[Pre, Ch. 10]** Preskill — sampling, estimation, and the role of shot noise in quantum experiments. - Hoeffding's original: W. Hoeffding, "Probability inequalities for sums of bounded random variables", *J. Amer. Statist. Assoc.* **58** (1963) — the source of the concentration bound used here. --- ← Prev: [Probability Spaces](01-probability-spaces.md) · Up: [Term 0](../README.md) · Next: [Classical Information Theory](03-classical-information-theory.md) →