{
  "cells": [
    {
      "cell_type": "markdown",
      "id": "markdown-000",
      "metadata": {},
      "source": [
        "# Measurement and probability\n",
        "\n",
        "[Read the lesson](https://tirta.co/a-byte-of-physics/measurement-and-probability). This notebook includes the lesson text and its Python examples in reading order. Run the setup cell, then run the remaining cells from top to bottom. In Colab, use File → Save a copy in Drive to keep your edits. You can run and edit this view without changing the course source. Examples use illustrative data unless stated otherwise."
      ]
    },
    {
      "cell_type": "code",
      "id": "code-001",
      "metadata": {},
      "source": [
        "import importlib.metadata\n",
        "import subprocess\n",
        "import sys\n",
        "\n",
        "requirements = {\"numpy\": (2, 0), \"sympy\": (1, 13), \"matplotlib\": (3, 8)}\n",
        "to_install = []\n",
        "for package, minimum in requirements.items():\n",
        "    try:\n",
        "        installed = importlib.metadata.version(package)\n",
        "        version = tuple(int(part) for part in installed.split(\".\")[:2])\n",
        "    except importlib.metadata.PackageNotFoundError:\n",
        "        version = (0, 0)\n",
        "    if version < minimum:\n",
        "        to_install.append(f\"{package}>={minimum[0]}.{minimum[1]}\")\n",
        "if to_install:\n",
        "    subprocess.check_call([sys.executable, \"-m\", \"pip\", \"install\", *to_install])\n",
        "\n",
        "import numpy as np\n",
        "import sympy as sp\n",
        "import matplotlib\n",
        "import matplotlib.pyplot as plt\n",
        "\n",
        "print(\"Python\", sys.version.split()[0])\n",
        "print(\"NumPy\", np.__version__)\n",
        "print(\"SymPy\", sp.__version__)\n",
        "print(\"Matplotlib\", matplotlib.__version__)\n"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "id": "markdown-002",
      "metadata": {},
      "source": [
        "Quantum mechanics predicts probabilities for possible measurement outcomes. A\n",
        "single measurement returns one outcome, not a percentage. Probabilities become\n",
        "testable by preparing many comparable systems, applying the same measurement,\n",
        "and comparing observed frequencies with the predicted distribution.\n",
        "\n",
        "The probability rule uses complex amplitudes. Amplitudes can add or cancel\n",
        "before their squared magnitudes become probabilities, which is why relative\n",
        "phase affects interference even though probabilities are non-negative.\n",
        "\n",
        "## The Born rule for discrete outcomes\n",
        "\n",
        "Suppose an observable\n",
        "$\\hat A$ has orthonormal,\n",
        "nondegenerate eigenstates\n",
        "$|a_n\\rangle$:\n",
        "\n",
        "\n",
        "$$\n",
        "\\hat A|a_n\\rangle=a_n|a_n\\rangle\n",
        "$$\n",
        "\n",
        "\n",
        "A normalized state can be expanded in this basis:\n",
        "\n",
        "\n",
        "$$\n",
        "|\\psi\\rangle=\\sum_n c_n|a_n\\rangle,\\qquad c_n=\\langle a_n|\\psi\\rangle\n",
        "$$\n",
        "\n",
        "\n",
        "The Born rule says that an ideal measurement returns eigenvalue\n",
        "$a_n$ with probability\n",
        "\n",
        "\n",
        "$$\n",
        "P(a_n)=|c_n|^2=|\\langle a_n|\\psi\\rangle|^2\n",
        "$$\n",
        "\n",
        "\n",
        "Completeness and normalization imply\n",
        "\n",
        "\n",
        "$$\n",
        "\\sum_nP(a_n)=1\n",
        "$$\n",
        "\n",
        "\n",
        "If input coefficients are not normalized, normalize their squared magnitudes:\n",
        "\n",
        "\n",
        "$$\n",
        "P(a_n)=\\frac{|c_n|^2}{\\sum_j|c_j|^2}\n",
        "$$\n",
        "\n",
        "\n",
        "This last formula is a convenience for a nonzero state representation. The\n",
        "physical state itself should normally be normalized.\n",
        "\n",
        "## Complex magnitude is not ordinary squaring\n",
        "\n",
        "If $c_n=x_n+iy_n$,\n",
        "\n",
        "\n",
        "$$\n",
        "|c_n|^2=c_n^*c_n=x_n^2+y_n^2\n",
        "$$\n",
        "\n",
        "\n",
        "Squaring $c_n$ without\n",
        "conjugation can produce a complex number and is not a probability. For\n",
        "coefficients with magnitudes one and two, the weights are one and four, so the\n",
        "probabilities are\n",
        "$1/5$ and\n",
        "$4/5$. They are not\n",
        "one third and two thirds."
      ]
    },
    {
      "cell_type": "code",
      "id": "code-003",
      "metadata": {},
      "source": [
        "import numpy as np\n",
        "\n",
        "\n",
        "def born_probabilities(state, basis):\n",
        "    if not np.allclose(basis.conj().T @ basis, np.eye(len(state))):\n",
        "        raise ValueError(\"Basis columns must be orthonormal.\")\n",
        "    amplitudes = basis.conj().T @ state\n",
        "    weights = abs(amplitudes) ** 2\n",
        "    if weights.sum() <= 0:\n",
        "        raise ValueError(\"State must have nonzero norm.\")\n",
        "    return weights / weights.sum()\n",
        "\n",
        "\n",
        "print(born_probabilities(np.array([1.0, 1j]), np.eye(2)))\n"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "id": "markdown-004",
      "metadata": {},
      "source": [
        "The returned values should all be non-negative and sum to one within a scaled\n",
        "floating-point tolerance. A large normalization discrepancy indicates an\n",
        "incomplete or nonorthonormal basis, a bad inner product, or arithmetic error.\n",
        "\n",
        "## The measurement basis matters\n",
        "\n",
        "Consider a two-state system with basis states\n",
        "$|0\\rangle$ and\n",
        "$|1\\rangle$. Define\n",
        "\n",
        "\n",
        "$$\n",
        "|+\\rangle=\\frac{|0\\rangle+|1\\rangle}{\\sqrt2},\\qquad |-\\rangle=\\frac{|0\\rangle-|1\\rangle}{\\sqrt2}\n",
        "$$\n",
        "\n",
        "\n",
        "If the state is $|+\\rangle$,\n",
        "a measurement in the plus/minus basis returns plus with certainty. A measurement\n",
        "in the zero/one basis returns either result with probability one half.\n",
        "\n",
        "The state has not become vague; the question changed. Quantum probabilities are\n",
        "always probabilities for a specified measurement. Listing amplitudes without\n",
        "the basis is incomplete.\n",
        "\n",
        "## Relative phase changes other-basis probabilities\n",
        "\n",
        "The states\n",
        "\n",
        "\n",
        "$$\n",
        "|\\psi_1\\rangle=\\frac{|0\\rangle+|1\\rangle}{\\sqrt2},\\qquad |\\psi_2\\rangle=\\frac{|0\\rangle-|1\\rangle}{\\sqrt2}\n",
        "$$\n",
        "\n",
        "\n",
        "have identical zero/one probabilities. Yet they are orthogonal: the first is\n",
        "$|+\\rangle$ and the\n",
        "second is $|-\\rangle$.\n",
        "A plus/minus measurement distinguishes them perfectly.\n",
        "\n",
        "By contrast, multiplying the whole state by a global phase\n",
        "\n",
        "\n",
        "$$\n",
        "|\\psi\\rangle\\longrightarrow e^{i\\chi}|\\psi\\rangle\n",
        "$$\n",
        "\n",
        "\n",
        "does not change any projection's squared magnitude. Global phase has no\n",
        "observable effect; relative phase does."
      ]
    },
    {
      "cell_type": "code",
      "id": "code-005",
      "metadata": {},
      "source": [
        "plus = np.array([1.0, 1.0]) / np.sqrt(2)\n",
        "minus = np.array([1.0, -1.0]) / np.sqrt(2)\n",
        "plus_minus_basis = np.column_stack([plus, minus])\n",
        "for state in [plus, minus]:\n",
        "    print(\"Zero/one basis:\", born_probabilities(state, np.eye(2)))\n",
        "    print(\"Plus/minus basis:\", born_probabilities(state, plus_minus_basis))\n"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "id": "markdown-006",
      "metadata": {},
      "source": [
        "## Degenerate outcomes and projectors\n",
        "\n",
        "If several eigenstates share the same eigenvalue, a measurement of that\n",
        "eigenvalue does not identify one arbitrary vector inside the subspace. Let\n",
        "$\\hat P_a$ project\n",
        "onto the entire eigenspace belonging to result\n",
        "$a$. Then\n",
        "\n",
        "\n",
        "$$\n",
        "P(a)=\\langle\\psi|\\hat P_a|\\psi\\rangle\n",
        "$$\n",
        "\n",
        "\n",
        "for a normalized state. The projectors for a complete ideal measurement obey\n",
        "\n",
        "\n",
        "$$\n",
        "\\hat P_a^\\dagger=\\hat P_a,\\qquad \\hat P_a^2=\\hat P_a,\\qquad \\sum_a\\hat P_a=I\n",
        "$$\n",
        "\n",
        "\n",
        "This form handles degeneracy and makes clear that probabilities belong to\n",
        "physical alternatives, not necessarily individual basis-vector labels.\n",
        "\n",
        "## Conditional state after an ideal measurement\n",
        "\n",
        "In the standard projective measurement model, if outcome\n",
        "$a$ is obtained, the\n",
        "conditional post-measurement state is\n",
        "\n",
        "\n",
        "$$\n",
        "|\\psi_a\\rangle=\\frac{\\hat P_a|\\psi\\rangle}{\\sqrt{\\langle\\psi|\\hat P_a|\\psi\\rangle}}\n",
        "$$\n",
        "\n",
        "\n",
        "provided the denominator is nonzero. An immediate repetition of the same ideal\n",
        "measurement then returns the same result with certainty. This update rule is not\n",
        "the unitary Schrödinger evolution of the measured subsystem.\n",
        "\n",
        "Real instruments can be destructive, noisy, weak, or only partially resolving.\n",
        "General measurements require a broader operator framework. The projective rule\n",
        "is an ideal model whose applicability must be stated.\n",
        "\n",
        "## Continuous position probability\n",
        "\n",
        "For a normalized position-space wavefunction,\n",
        "\n",
        "\n",
        "$$\n",
        "\\int_{-\\infty}^{\\infty}|\\psi(x)|^2\\,dx=1\n",
        "$$\n",
        "\n",
        "\n",
        "The probability of finding a particle in the finite interval from\n",
        "$a$ to\n",
        "$b$ is\n",
        "\n",
        "\n",
        "$$\n",
        "P(a\\leq x\\leq b)=\\int_a^b|\\psi(x)|^2\\,dx\n",
        "$$\n",
        "\n",
        "\n",
        "For a continuous variable, probability at one exact point is generally zero;\n",
        "$|\\psi(x)|^2$\n",
        "is a probability **density** with units of inverse length in one dimension.\n",
        "Only its integral over an interval is dimensionless probability.\n",
        "\n",
        "On a grid with cell widths\n",
        "$w_j$, use\n",
        "\n",
        "\n",
        "$$\n",
        "P_j\\approx|\\psi_j|^2w_j,\\qquad \\sum_jP_j\\approx1\n",
        "$$\n",
        "\n",
        "\n",
        "Omitting the cell width confuses density with probability and fails on a\n",
        "nonuniform grid.\n",
        "\n",
        "## Expectation and spread\n",
        "\n",
        "For discrete outcomes, the expectation is\n",
        "\n",
        "\n",
        "$$\n",
        "\\langle A\\rangle=\\sum_a aP(a)\n",
        "$$\n",
        "\n",
        "\n",
        "and the variance is\n",
        "\n",
        "\n",
        "$$\n",
        "(\\Delta A)^2=\\sum_a(a-\\langle A\\rangle)^2P(a)=\\langle A^2\\rangle-\\langle A\\rangle^2\n",
        "$$\n",
        "\n",
        "\n",
        "Expectation is the average approached across many comparable trials; it need\n",
        "not be a possible individual outcome. For a two-outcome observable with values\n",
        "positive one and negative one, an equal distribution has expectation zero even\n",
        "though zero is never measured.\n",
        "\n",
        "Variance describes intrinsic spread in the predicted outcome distribution.\n",
        "Experimental scatter can additionally include imperfect preparation, detector\n",
        "noise, calibration error, and environmental disturbance.\n",
        "\n",
        "## Comparing predictions with repeated trials\n",
        "\n",
        "Suppose outcome $a$ has\n",
        "predicted probability $p_a$\n",
        "and is observed $n_a$ times\n",
        "in $N$ trials. The observed\n",
        "frequency is\n",
        "\n",
        "\n",
        "$$\n",
        "f_a=\\frac{n_a}{N}\n",
        "$$\n",
        "\n",
        "\n",
        "For independent, identically prepared Bernoulli trials, the standard deviation\n",
        "of the frequency is\n",
        "\n",
        "\n",
        "$$\n",
        "\\sigma_{f_a}=\\sqrt{\\frac{p_a(1-p_a)}{N}}\n",
        "$$\n",
        "\n",
        "\n",
        "The typical absolute fluctuation shrinks like\n",
        "$1/\\sqrt N$.\n",
        "A finite experiment need not reproduce the predicted percentage exactly.\n",
        "Agreement should be judged with a statistical model, not by demanding exact\n",
        "counts."
      ]
    },
    {
      "cell_type": "code",
      "id": "code-007",
      "metadata": {},
      "source": [
        "rng = np.random.default_rng(42)\n",
        "probabilities = np.array([0.25, 0.75])\n",
        "counts = rng.multinomial(1000, probabilities)\n",
        "frequencies = counts / counts.sum()\n",
        "standard_errors = np.sqrt(probabilities * (1 - probabilities) / counts.sum())\n",
        "print(\"Observed frequencies:\", frequencies)\n",
        "print(\"Expected standard errors:\", standard_errors)\n"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "id": "markdown-008",
      "metadata": {},
      "source": [
        "Independence and identical preparation are assumptions. Drift, correlated\n",
        "trials, post-selection, background counts, and detector inefficiency require a\n",
        "more detailed likelihood.\n",
        "\n",
        "## Measurement uncertainty is not one idea\n",
        "\n",
        "The position–momentum relation\n",
        "\n",
        "\n",
        "$$\n",
        "\\Delta x\\,\\Delta p_x\\geq\\frac{\\hbar}{2}\n",
        "$$\n",
        "\n",
        "\n",
        "concerns the standard deviations of outcome distributions in one quantum state.\n",
        "It is not the same as a device's calibration uncertainty or the statistical\n",
        "uncertainty of estimating a mean from finite trials. All three can matter in an\n",
        "experiment, but they enter the analysis differently.\n",
        "\n",
        "Measurement disturbance is also a separate question. An ideal projective\n",
        "measurement changes the conditional state, while the Robertson uncertainty\n",
        "relation follows from noncommuting observables even before a particular\n",
        "measurement apparatus is specified.\n",
        "\n",
        "## Validation checklist for a probability calculation\n",
        "\n",
        "A reliable calculation demonstrates that:\n",
        "\n",
        "- the state and measurement basis have matching dimensions;\n",
        "- the basis is orthonormal or a valid projector description is used;\n",
        "- every computed probability is real and non-negative;\n",
        "- probabilities for a complete outcome set sum to one;\n",
        "- global phase leaves all probabilities unchanged;\n",
        "- a known eigenstate gives its eigenvalue with probability one;\n",
        "- continuous densities include the integration measure;\n",
        "- observed counts are not mistaken for exact theoretical probabilities.\n",
        "\n",
        "This lesson uses pure states and ideal projective measurements. Mixed states use\n",
        "density operators, and realistic generalized measurements use positive\n",
        "operator-valued measures. Decoherence, detector response, and interpretations\n",
        "of what constitutes an outcome are deeper topics. The Born rule and its\n",
        "experimentally testable distributions remain common to the standard\n",
        "calculations presented here.\n",
        "\n",
        "## A worked Python notebook\n",
        "\n",
        "Normalise a two-state vector with real amplitudes cosine and sine of a mixing angle. The question for this calculation is precise: how does **probability of outcome one**\n",
        "change when we change **amplitude mixing angle**, while keeping the other choices fixed?\n",
        "The earlier sections explain the physical relationship. Here we turn one part of that\n",
        "relationship into a calculation that you can run, inspect, and change in Google Colab.\n",
        "\n",
        "Probabilities use squared complex magnitudes, while a different measurement basis can reveal phase information. Keep that sentence beside the results. A formula can return a number outside\n",
        "its useful range, and a plotting library will usually draw that number without asking\n",
        "whether the physical assumptions still hold. The selected interval is from\n",
        "0 to 1.5708 rad. It gives us a specific domain to discuss\n",
        "instead of leaving the graph's range to a software default.\n",
        "\n",
        "Use the notebook link above the lesson, run the setup cell, and then run the cells in\n",
        "order. The code below is also included there. NumPy handles the arrays of numerical\n",
        "values, SymPy keeps the equation symbolic, and Matplotlib draws the results. If these\n",
        "libraries or notebook cells are unfamiliar, start with\n",
        "[Your first physics notebook](https://tirta.co/a-byte-of-physics/your-first-physics-notebook).\n",
        "\n",
        "### State the relationship before choosing a graph\n",
        "\n",
        "In this small example, `q` stands for amplitude mixing angle, expressed in rad. The output is\n",
        "probability of outcome one, expressed in dimensionless. The equation below uses the numerical values in those\n",
        "units. Any fixed coefficients therefore carry the units needed to make the relation\n",
        "consistent. A value supplied in another unit must be converted before it enters this\n",
        "expression.\n",
        "\n",
        "\n",
        "$$\n",
        "f(q) = \\sin^{2}{\\left(q \\right)}\n",
        "$$\n",
        "\n",
        "\n",
        "This is a function of one changing input. It is not a statement that the physical\n",
        "system has only one variable. We deliberately hold the other inputs fixed so that one\n",
        "cause of variation is visible. Changing several inputs at once would make it harder\n",
        "to tell which change produced the output difference. Write down the held-fixed choices\n",
        "before trying a different experiment.\n",
        "\n",
        "At 0.523599 rad, the expression gives 0.25 dimensionless. This is our\n",
        "reference calculation. It is a prediction of the stated example, not a new measurement.\n",
        "Its displayed digits help us compare implementations. They should not be read as a\n",
        "claim that physical inputs are known to six significant figures."
      ]
    },
    {
      "cell_type": "code",
      "id": "code-009",
      "metadata": {},
      "source": [
        "import numpy as np\n",
        "import sympy as sp\n",
        "import matplotlib.pyplot as plt\n",
        "\n",
        "q = sp.Symbol(\"q\", real=True)\n",
        "expression = sp.sympify(\"sin(q)**2\", locals={\"q\": q})\n",
        "reference_input = 0.5235987755982989\n",
        "reference_output = float(expression.subs(q, reference_input))\n",
        "print(\"Input:\", reference_input, \"rad\")\n",
        "print(\"Predicted probability of outcome one:\", reference_output, \"dimensionless\")\n"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "id": "markdown-010",
      "metadata": {},
      "source": [
        "`sp.Symbol` declares a mathematical symbol. It does not assign a numerical measurement\n",
        "to that symbol. `subs` then replaces the symbol with our chosen input, and `float`\n",
        "converts the resulting symbolic number into a standard Python numerical value.\n",
        "Keeping these steps separate lets us inspect the equation before numerical substitution.\n",
        "It also prevents a later plot from becoming the only surviving record of the model.\n",
        "\n",
        "### Compare a few cases you can explain\n",
        "\n",
        "Before drawing hundreds of points, calculate three. The lower endpoint gives\n",
        "0 dimensionless. The midpoint, at 0.785398 rad, gives\n",
        "0.5 dimensionless. The upper endpoint gives 1 dimensionless.\n",
        "Read these together with the assumptions: they describe different choices of amplitude mixing angle\n",
        "within the same example, with the other inputs held fixed.\n",
        "\n",
        "These three values are a useful first check, but they do not establish the shape between\n",
        "them. Two endpoints can be equal even when the curve bends or oscillates in between.\n",
        "Likewise, a large change between endpoints does not tell us whether the relationship is\n",
        "linear. We need the equation and intermediate samples to answer those questions."
      ]
    },
    {
      "cell_type": "code",
      "id": "code-011",
      "metadata": {},
      "source": [
        "evaluate = sp.lambdify(q, expression, modules=\"numpy\")\n",
        "inputs = np.array([0.0, 0.7853981633974483, 1.5707963267948966])\n",
        "outputs = np.asarray(evaluate(inputs), dtype=float)\n",
        "outputs = np.broadcast_to(outputs, inputs.shape)\n",
        "for input_value, output_value in zip(inputs, outputs):\n",
        "    print(f\"{input_value:.6g} rad -> {output_value:.6g} dimensionless\")\n",
        "assert np.all(np.isfinite(outputs))\n"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "id": "markdown-012",
      "metadata": {},
      "source": [
        "`lambdify` turns the symbolic expression into a numerical function. We explicitly select\n",
        "NumPy as its numerical backend so it can accept an array. The same expression now has\n",
        "two useful forms: one for algebra and one for calculating many cases. The assertion\n",
        "checks that the chosen inputs produce finite outputs. It would catch an infinity or an\n",
        "undefined numerical result, but a finite value can still be physically inappropriate.\n",
        "\n",
        "### Make a graph that answers the question\n",
        "\n",
        "The horizontal axis shows amplitude mixing angle; the vertical axis shows probability of outcome one. Every plotted\n",
        "point is a prediction of the same expression. We mark the reference case separately\n",
        "so that a reader can connect the earlier printed value to its place on the curve.\n",
        "The axis labels include units because a shape without a scale cannot tell us the size\n",
        "of a physical effect."
      ]
    },
    {
      "cell_type": "code",
      "id": "code-013",
      "metadata": {},
      "source": [
        "input_grid = np.linspace(0.0, 1.5707963267948966, 301)\n",
        "output_grid = np.asarray(evaluate(input_grid), dtype=float)\n",
        "output_grid = np.broadcast_to(output_grid, input_grid.shape)\n",
        "fig, ax = plt.subplots(figsize=(7, 4))\n",
        "ax.plot(input_grid, output_grid, label=\"Stated model\")\n",
        "ax.scatter([reference_input], [reference_output], label=\"Reference case\", zorder=3)\n",
        "ax.set_xlabel(\"Amplitude mixing angle (rad)\")\n",
        "ax.set_ylabel(\"Probability of outcome one (dimensionless)\")\n",
        "ax.legend()\n",
        "ax.grid(alpha=0.25)\n",
        "fig.tight_layout()\n",
        "plt.show()\n"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "id": "markdown-014",
      "metadata": {},
      "source": [
        "Look for a sign change, a turning point, or a region where the curve becomes steep.\n",
        "Connect each feature to the expression rather than to the colour of the line. If the\n",
        "curve has a minimum or maximum, changing the input in opposite directions may have\n",
        "similar effects on the output near that point. If it is nearly straight over a short\n",
        "interval, a local linear approximation may be useful even when the complete expression\n",
        "is nonlinear.\n",
        "\n",
        "Increasing the number of samples makes the line smoother and may reveal features that\n",
        "a coarse grid missed. It does not improve the physical inputs or enlarge the model's\n",
        "validity range. For quantities that can take only discrete values, use the line as a\n",
        "guide to the relationship and interpret only the allowed input values. A plotting\n",
        "function does not enforce such restrictions for you.\n",
        "\n",
        "### Measure the local sensitivity\n",
        "\n",
        "A derivative asks a more specific question than the whole graph: how quickly does\n",
        "probability of outcome one change near the reference input? Its units are the output unit divided by\n",
        "the input unit. For this calculation, the derivative at 0.523599 rad is\n",
        "0.866025 in those ratio units. Locally, the output increases as the\n",
        "input increases. This statement is local; the sign or size of the slope can be different\n",
        "elsewhere on the curve.\n",
        "\n",
        "\n",
        "$$\n",
        "f'(q) = 2 \\sin{\\left(q \\right)} \\cos{\\left(q \\right)}\n",
        "$$\n",
        "\n",
        "\n",
        "For a small input change, multiply the derivative by that change to estimate the output\n",
        "change. Then compare that estimate with a fresh evaluation of the full expression.\n",
        "The distinction between these two calculations matters: the derivative gives a tangent\n",
        "approximation, while substitution evaluates the chosen expression at the new point."
      ]
    },
    {
      "cell_type": "code",
      "id": "code-015",
      "metadata": {},
      "source": [
        "slope_expression = sp.diff(expression, q)\n",
        "slope_at_reference = float(slope_expression.subs(q, reference_input))\n",
        "changed_input = 0.5314527572322734\n",
        "linear_prediction = reference_output + slope_at_reference * (\n",
        "    changed_input - reference_input\n",
        ")\n",
        "full_prediction = float(evaluate(changed_input))\n",
        "print(\"Symbolic derivative:\", slope_expression)\n",
        "print(\"Tangent prediction:\", linear_prediction, \"dimensionless\")\n",
        "print(\"Full expression:\", full_prediction, \"dimensionless\")\n",
        "print(\"Approximation difference:\", full_prediction - linear_prediction, \"dimensionless\")\n"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "id": "markdown-016",
      "metadata": {},
      "source": [
        "Here the changed input is 0.531453 rad. The tangent estimate is\n",
        "0.256802 dimensionless, while the full expression gives\n",
        "0.256832 dimensionless. If the two values are close, it means the chosen\n",
        "input change is small enough for this local approximation at this point. It does not\n",
        "mean the entire model is linear. If the slope is zero, a quadratic or higher-order term\n",
        "may supply the leading change instead.\n",
        "\n",
        "This is also the beginning of uncertainty propagation. If the input has a small stated\n",
        "standard uncertainty, the magnitude of the slope times that uncertainty estimates one\n",
        "contribution to output uncertainty. Do not automatically use the same rule for a wide\n",
        "interval, a discontinuity, a turning point with zero first derivative, or several\n",
        "correlated uncertain inputs. Those situations need more of the model than one slope.\n",
        "\n",
        "### Check the derivative with numerical differences\n",
        "\n",
        "We can estimate the same slope from values on either side of the reference point.\n",
        "Subtract the lower output from the upper output and divide by twice the input step.\n",
        "Parentheses are essential: the entire output difference must be divided by the entire\n",
        "input interval. Without them, Python can execute a different mathematical expression\n",
        "while still producing a plausible number.\n",
        "\n",
        "The calculation below uses three decreasing step sizes. Each estimate is compared with\n",
        "the symbolic derivative, so the reference is independent of the numerical difference\n",
        "formula. Near an endpoint, it is acceptable here to evaluate a small extension only\n",
        "because the displayed mathematical expression is smooth there; such a check does not\n",
        "assert that extrapolated values describe additional physical cases."
      ]
    },
    {
      "cell_type": "code",
      "id": "code-017",
      "metadata": {},
      "source": [
        "input_span = 1.5707963267948966\n",
        "for fraction in [1e-2, 1e-3, 1e-4]:\n",
        "    step = input_span * fraction\n",
        "    upper = float(evaluate(reference_input + step))\n",
        "    lower = float(evaluate(reference_input - step))\n",
        "    numerical_slope = (upper - lower) / (2 * step)\n",
        "    difference = numerical_slope - slope_at_reference\n",
        "    print(f\"step={step:.6g}: slope={numerical_slope:.6g}, difference={difference:.6g}\")\n"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "id": "markdown-018",
      "metadata": {},
      "source": [
        "For a smooth function, reducing the step often reduces the approximation error at first.\n",
        "Very small steps can eventually make subtraction of nearly equal floating-point values\n",
        "less reliable. Therefore “smaller” is not an unlimited guarantee of “better.” Compare\n",
        "the observed behaviour over a useful range of steps instead of choosing the smallest\n",
        "number the computer can represent.\n",
        "\n",
        "An absolute slope difference is useful when the expected slope is near zero. A relative\n",
        "difference would divide by a very small reference and could look enormous even when the\n",
        "absolute discrepancy is tiny. Away from zero, a relative difference can help compare\n",
        "quantities with very different scales. The tolerance should reflect the question and\n",
        "the numerical method, rather than a universal decimal threshold copied into every test.\n",
        "\n",
        "### Change the example and preserve the reasoning\n",
        "\n",
        "Make one copy of the reference calculation before editing it. Choose a second value of\n",
        "amplitude mixing angle inside the stated interval, predict whether probability of outcome one will rise or fall,\n",
        "and then evaluate the expression. Place the prediction and the result in a text cell.\n",
        "If they disagree, identify whether the mistake was in the sign, the unit conversion,\n",
        "the algebra, or the assumption about how the curve behaves.\n",
        "\n",
        "Next choose a shorter plotting interval around the reference point. Compare the local\n",
        "curve with its tangent approximation. Then widen the interval again. The difference\n",
        "shows why a useful local estimate need not describe the whole range. Keep the axis\n",
        "labels and numerical scale visible when comparing the plots; automatic axis limits can\n",
        "make two different effects appear equally large.\n",
        "\n",
        "Finally, return to the physical limit stated at the start: Probabilities use squared complex magnitudes, while a different measurement basis can reveal phase information. Describe one change\n",
        "to the setup that would violate that limit. State what new input, equation, or observation\n",
        "would be needed before trusting an answer in the changed situation. This last step\n",
        "connects the Python experiment to the physics instead of treating successful execution\n",
        "as the only criterion for a useful result.\n",
        "\n",
        "Before saving, restart the runtime and run all cells in order. A result that depends on\n",
        "an earlier, now-deleted input cell is not reproducible from the saved notebook. Keep\n",
        "unrounded values for later calculations, round only the reported output, and label any\n",
        "new data as measured, assumed, or generated. The graph and printed values should remain\n",
        "traceable to the same expression and the same choices of inputs.\n",
        "\n",
        "### Recover the output change by integrating the slope\n",
        "\n",
        "There is another way to check that the expression and its derivative describe the same\n",
        "relationship. Add up the small output changes predicted by the slope across the chosen\n",
        "input interval. In the continuous limit, this integral equals the output at the upper\n",
        "endpoint minus the output at the lower endpoint. The result has units of dimensionless:\n",
        "the input unit in the integration cancels the input unit in the denominator of the slope.\n",
        "\n",
        "This statement is about the mathematics of this particular relationship. It does not\n",
        "assign a new physical meaning to an arbitrary area under a graph. For example, an area\n",
        "under force versus distance represents work, while an area under a different pair of\n",
        "quantities may represent something else or have no standard name. Always carry the\n",
        "axis units through the calculation before naming an integral."
      ]
    },
    {
      "cell_type": "code",
      "id": "code-019",
      "metadata": {},
      "source": [
        "slope_function = sp.lambdify(q, slope_expression, modules=\"numpy\")\n",
        "endpoint_change = float(evaluate(1.5707963267948966) - evaluate(0.0))\n",
        "for points in [101, 201, 401]:\n",
        "    grid = np.linspace(0.0, 1.5707963267948966, points)\n",
        "    slopes = np.asarray(slope_function(grid), dtype=float)\n",
        "    slopes = np.broadcast_to(slopes, grid.shape)\n",
        "    integrated_change = np.trapezoid(slopes, grid)\n",
        "    print(\n",
        "        points,\n",
        "        \"integrated change:\",\n",
        "        integrated_change,\n",
        "        \"endpoint change:\",\n",
        "        endpoint_change,\n",
        "    )\n"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "id": "markdown-020",
      "metadata": {},
      "source": [
        "The trapezoidal rule joins neighbouring slope samples with straight segments and adds\n",
        "their signed areas. A negative contribution remains negative; replacing every slope\n",
        "by its magnitude would calculate a different quantity. Increasing the sample count\n",
        "provides a practical refinement check. If the values do not settle toward the endpoint\n",
        "difference, first inspect the interval, singularities, and array ordering.\n",
        "\n",
        "Do not conclude that two agreeing calculations prove the model is physically complete.\n",
        "Both use the same selected relationship, so an omitted effect would be omitted from\n",
        "both. Their agreement is evidence that the numerical tools implement that relationship\n",
        "consistently. Comparing with an independent experiment or a more appropriate physical\n",
        "model is a separate test.\n",
        "\n",
        "Keep this comparison in the notebook when you change the expression. It is particularly\n",
        "useful after adding a coefficient, changing a sign convention, or converting an input\n",
        "unit. A code cell that produces a smooth graph can conceal these mistakes; a second\n",
        "calculation of the same output change gives you another place to notice them.\n",
        "\n",
        "## Sources\n",
        "\n",
        "- [OpenStax, *University Physics Volume 3*, §7.1, “Wave Functions”](https://openstax.org/books/university-physics-volume-3/pages/7-1-wave-functions).\n",
        "- [OpenStax, *University Physics Volume 3*, §7.2, “The Heisenberg Uncertainty Principle”](https://openstax.org/books/university-physics-volume-3/pages/7-2-the-heisenberg-uncertainty-principle).\n",
        "- [OpenStax, *University Physics Volume 3*, §7.3, “The Schrödinger Equation”](https://openstax.org/books/university-physics-volume-3/pages/7-3-the-schrodinger-equation).\n",
        "- [Tipler and Llewellyn, *Modern Physics*, quantum probability and measurement](https://web.pdx.edu/~pmoeck/books/Tipler_Llewellyn.pdf).\n",
        "- [Benjamin Crowell, *Light and Matter*, quantum physics](https://open.umn.edu/opentextbooks/textbooks/light-and-matter).\n",
        "- [Physics LibreTexts, Quantum Mechanics bookshelf](https://phys.libretexts.org/Bookshelves/Quantum_Mechanics)."
      ]
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "name": "python",
      "version": "3.12"
    },
    "colab": {
      "name": "measurement-and-probability.ipynb",
      "provenance": []
    },
    "tirta": {
      "lesson": "https://tirta.co/a-byte-of-physics/measurement-and-probability",
      "source_sha256": "3a8fa34a363c66a1b641135bd5a2131e46677d46ece5bf9d9031640b9dcd9952"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 5
}
