Thinking in Waves Episode 51 / Part XVIII — Putting Error Bars On It

Episode 49's handbook was only ever tested on noiseless functions

Error Bars on the Decade Count —
One periodogram ordinate always scatters by 5.57 dB This is the first line of Episode 50's "what to take up next" — saying "a power law over how many decades" from noisy real data.
★ One periodogram ordinate scatters by 5.57 dB. Independent of signal strength, frequency and instrument.
★★ And what we found — Episode 49's handbook cannot be used on real data as written. Even averaging 10,000 ordinates, a 6-decade power law was reported as only 2.6 decades.

Tools needed: Episodes 9, 22, 27, 42, 47, 48, 49Computation: kenshou/gosa.js

Since Episode 27 we have spoken in "how many decades."
★ It is only right to put an error bar on the decades themselves — so we start on Episode 50 §05's first line.

01One periodogram ordinate always scatters the same

Each ordinate is "the true spectrum × \(\chi^2_2/2\)" $$\ln\!\left(\frac{\chi^2_2}{2}\right):\quad \text{mean}=-\gamma=-0.5772,\qquad \text{variance}=\frac{\pi^2}{6},\quad \text{sd}=\frac{\pi}{\sqrt6}=1.2825$$
QuantityMonte CarloTheory
mean \(E[\ln(\chi^2_2/2)]\)−0.575227−0.577216
standard deviation1.2790371.282550
★ the same thing in dB5.5548 dB★ 5.5700 dB

The point of §01

One periodogram ordinate scatters by 5.57 dB.
Independent of signal strength, frequency and instrument.

★ This is not "because the measurement is poor" — one ordinate carries only two degrees of freedom, a real and an imaginary part.
★ The famous \(-0.5772\) bias (\(-2.51\) dB) is there too, but it affects only the intercept, never the slope.

02The slope's error is \(1.2825/\sqrt M\)

For frequencies uniform on \(0\sim F\), \(\mathrm{Var}(\ln f)=1\) exactly $$\ln S=c-\beta\ln f+\varepsilon\quad\Longrightarrow\quad \mathrm{sd}(\beta)=\frac{1.2825}{\sqrt M}$$
\(M\) (number of ordinates)\(\mathrm{sd}(\beta)\) (measured)\(1.2825/\sqrt M\)ratio
2000.0895320.0906900.9872
10000.0394020.0405580.9715
200000.0092590.0090691.0209

★ \(\mathrm{Var}(\ln f)=1\) is a property of linearly spaced frequencies — widening the decade span does not change it. Only the count \(M\) matters.

03★ How much data does Episode 42's and 47's shift need?

Since \(\text{dB/oct}=3\beta\), \(\mathrm{sd}(\text{dB/oct})=3.85/\sqrt M\).

Difference to be measured\(M\) for \(1\sigma\)\(M\) for \(3\sigma\)record length \(N\)
★ turbulent intermittency (Ep. 42) 0.08122.3e+32.0e+44.0e+4
★ the universe's \(n_s\) (Ep. 47) 0.10531.3e+31.2e+42.4e+4
one whole order (3 dB/oct)1.61530

The point of §03

★ To claim 0.081 dB/oct at \(3\sigma\) takes 20,000 ordinates — a 40,000-point record. Surprisingly modest: statistics alone is within reach.
★★★ But that is not where the difficulty lay.

04★★ Systematics — first, notice that they cancel

A real spectrum is not an infinite power law; it bends at both ends (call the inertial range \(R\) decades wide). We place a window of width \(D\) at its centre and fit a line.

(a) symmetric model\(D=0.5\)\(D=1\)\(D=2\)\(D=3\)
\(R=3\) decades1.3e−143.3e−152.4e−148.0e−15
\(R=6\) decades6.7e−161.8e−142.1e−141.1e−14

§04 (a) — ★★★ a result we did not expect

Essentially exactly zero — the two end curvatures pull the slope in opposite directions and cancel.

★ Simply centring the window in the inertial range removes the leading systematic.
★★ This is the single most useful practical note here. We had not noticed it before writing this episode.

But it depends on the symmetry. Real turbulence cuts off exponentially at high frequency (the dissipation range) and bends gently at low frequency — make it asymmetric and the cancellation fails.

(b) asymmetric model\(D=0.5\)\(D=1\)\(D=2\)\(D=3\)
\(R=3\) decades0.09120.09820.12340.1380
\(R=4\) decades0.03030.03320.04590.0692
\(R=6\) decades0.00310.00340.00490.0083
(c) shifting the window off-centre (symmetric model, \(R=4,D=2\))★ systematic [dB/oct]
shift 0 decades−0.00000
0.250.00937
0.50★ 0.03201
1.000.26611

The point of §04

★★★ Half a decade off-centre already gives 0.03 dB/oct — the cancellation breaks quickly.
★ With only 3 decades of inertial range and an asymmetric spectrum, the systematic is 0.09–0.14 dB/oct — as large as or larger than Episode 42's 0.081.
★★ Statistical error shrinks with \(M\); systematic error does not.

05★★★ There is an optimal window width

Statistics falls with a wider window; systematics rises $$\text{stat}=\frac{3\times1.2825}{\sqrt{nD}\cdot\frac{\ln10}{\sqrt{12}}\cdot D}\ \propto D^{-3/2},\qquad \text{syst}=\text{the table of §04}\ \nearrow$$
inertial range \(R\)★ optimal \(D\)statsyst★ total [dB/oct]
3 decades2.950.11420.1382★ 0.1793
4 decades3.950.07370.08880.1155
6 decades5.350.04680.0355★ 0.0587

The point of §05

★★★ At the optimum, statistical and systematic errors are about equal. Shrinking only one does not reduce the total — so "just measure longer" plateaus.
★★ With only 3 decades of inertial range, the best achievable is 0.18 dB/oct — short of Episode 42's 0.081. With 6 decades it falls to 0.059 and finally reaches.
★ So "choose a system with a wide inertial range" beats "measure longer" — as experiment design, that comes first.

06★ Episode 49's handbook cannot be used as written

We run Episode 49's procedure (the longest run where the local order stays within \(\pm0.1\)) with noise. The true signal is 1/f over 6 decades.

averaged ordinates \(K\)sd per point [dB]★ decades reportedtruth
15.5700.049 ± 0.0516.0
1000.5570.205 ± 0.0766.0
100000.056★ 2.592 ± 1.0396.0

The point of §06 — ★★★ the core of this episode

Even averaging 10,000 ordinates down to 0.056 dB per point,
a 6-decade power law is reported as only 2.6 decades.

★ The cause is building the local order from differences between neighbours. A difference multiplies the noise by \(\sqrt2\), and dividing by the spacing multiplies it again by \(1/\Delta\).
★★ Episode 49's handbook was only tested on noiseless functions — applied to real data it fails.
★ This is a defect this series found in itself: we did what Episode 50 said to do next, and out it came.

07★ The fix — fit over a window instead of differencing

window \(w\) [decades]sd of the order★ decades reportedceiling \(6-w\)
0.100.698410.051 ± 0.0175.9
0.300.127950.302 ± 0.0875.7
0.600.044663.409 ± 1.3285.4
1.000.02005★ 5.0175.0
1.500.011284.5354.5

The point of §07

★ Build the local order by least squares over a window of width \(w\) instead of by differencing, and it works. The noise falls as the square root of the number of points in the window.
★★★ But a window of width \(w\) consumes \(w/2\) at each end, so the reportable decades are capped at \(6-w\). Indeed \(w=1.0\) gives 5.017 (ceiling 5.0) and \(w=1.5\) gives 4.535 (ceiling 4.5) — it is the edge, not the noise, that limits.
★★ Too narrow and noise fragments the run; too wide and you lose the ends — the window width itself is Episodes 9 and 22's "allocation of resolution."

08Amendments to Episode 49's handbook

StepNote
③-a take the per-point error as 5.57 dBfor a periodogram it is always this
③-b averaging \(K\) gives \(5.57/\sqrt K\) dBaverage across frequency or across trials, either way
★ ③-c build the local order by window fitting★ never by differencing neighbours (§06)
★ ③-d use \(w=0.5\)–1 decade★ the reportable span is capped at \(D_{\rm total}-w\)
★ ③-e centre the window in the inertial range★ the end curvatures then cancel (§04)
③-f stop at the \(M\) where stat = systbeyond that, more data buys nothing

09What wave language can and cannot say

ItemVerdictBasis
★ Scatter of one ordinate◎ numerics★ 5.57 dB, universally
\(\mathrm{sd}(\beta)=1.2825/\sqrt M\)◎ numerics★ set by count, not decade span
Record length for 0.081 dB/oct◎ numerics40,000 points at \(3\sigma\) (statistics only)
★ Systematics cancel in a centred window◎ numerics★ centring alone removes it
Half a decade off-centre breaks it◎ numerics0.032 dB/oct
★★ There is an optimal window width◎ numerics★ stat = syst there
★ The naive decade count fails◎ numerics★ 2.6 decades even after 10,000 averages
Window fitting repairs it◎ numerics5.0 decades at \(w=1\) (the ceiling)
★ Validation on real data× not done★ synthetic data only, here too

Exercises

  1. Why is §01's 5.57 dB independent of the instrument?
    Show answer
    Because one periodogram ordinate has only two degrees of freedom, a real and an imaginary part. However good the instrument, one estimate at one frequency carries the information of a single complex number.
    ★ Its squared magnitude follows \(\chi^2_2\), and the variance of its logarithm is the universal constant \(\pi^2/6\).
    ★★ There is exactly one way to reduce it — average more independent estimates (\(5.57/\sqrt K\) dB for \(K\) of them). Binning in frequency or repeating trials works equally. ★ The same structure as Episode 9's "the allocation is choosable but the total is not" — here the total is the number of degrees of freedom.
  2. §04's systematics cancelled. What must you watch for in using that on real data?
    Show answer
    Three things.
    ★ ① Are both ends the same shape? As §04 (b) shows, an exponential cutoff on one side kills the cancellation. Turbulence's dissipation range is exactly that.
    ★ ② Is the window really centred? Half a decade off gives 0.03 dB/oct (§04 (c)). You need independent knowledge of where the inertial range ends.
    ★ ③ Do not rely on it. The cancellation is luck, not a guarantee, so always check how much the answer moves when you shift the window (sensitivity analysis).
    ★★ In this article's stance, the cancellation is "a bonus if it happens," never "a premise you may assume."
  3. §06 found a defect in the handbook. Does that make Episode 49 wrong?
    Show answer
    Not wrong, but under-declared. Episode 49's own text says "only noiseless functions were treated; with noisy real data the decade estimates acquire statistical error."
    ★ In practice, though, what they acquire is not "statistical error" but a large systematic underestimate — the answer does not gain an error bar; it drops from 6 decades to 2.6.
    ★★ So Episode 49's caveat was too soft. "Acquires error" and "cannot be used" are different, and the truth was the latter.
    ★ Having decided in Episode 15 never to delete a calculation we lost, this is the way to write it — and since the fix (§07) came with it, the net is still forward.

Summary — the error bar came in two parts

One periodogram ordinate scatters by 5.57 dB (\(\pi/\sqrt6=1.2825\) nats) — independent of signal, frequency and instrument. The \(-2.51\) dB mean bias affects only the intercept, never the slope.

The slope's error is \(\mathrm{sd}(\beta)=1.2825/\sqrt M\) (confirmed numerically) — set by the count, not by the decade span. On statistics alone, 0.081 dB/oct at \(3\sigma\) needs 40,000 points — within reach.

★★★ But — centre the window in the inertial range and
the two end curvatures pull opposite ways and cancel (to \(10^{-14}\)).
★ Half a decade off-centre and 0.032 dB/oct is back.

When the ends are asymmetric (an exponential high-frequency cutoff — real turbulence) the cancellation fails, giving 0.09–0.14 dB/oct over a 3-decade inertial range. Combining statistics and systematics gives an optimal window width, at which the two are about equal — best 0.18 dB/oct at \(R=3\), 0.059 at \(R=6\). ★ "Choose a system with a wide inertial range" beats "measure longer."

★★★ And the core — Episode 49's handbook cannot be used on real data as written.
Even averaged down to 0.056 dB per point, a 6-decade power law reads as 2.6.

★ The cause was building the order from differences between neighbours — differencing multiplies noise by \(\sqrt2/\Delta\). The fix is least squares over a window of width \(w\), which recovered 5.0 decades at \(w=1\).
★★ But a width-\(w\) window eats \(w/2\) at each end, capping the reportable span at \(D_{\rm total}-w\) — too narrow and noise fragments it, too wide and you lose the ends. The window width itself is Episodes 9 and 22's "allocation of resolution."

★ Honestly: this article too used only synthetic data; nothing was applied to a real measurement. ★ And Episode 49's caveat was too soft — not "acquires error" but "cannot be used."

This document is Episode 51 of the "Thinking in Waves" series. Computations are in kenshou/gosa.js. On the numbers here: the Monte Carlo of log-periodogram statistics, the slope standard deviations, the required record lengths, the systematic-error calculations, the decade-count simulations, and the improvement from window fitting were all computed here. That the periodogram is \(\chi^2_2\) distributed, that the log-periodogram has mean \(-\gamma\) and variance \(\pi^2/6\) (the basis of the Geweke–Porter-Hudak method), and the resulting slope-estimation error are all well-known results and none are discoveries of this article — what this article did was translate them into this series' language of "dB/oct" and "decades," apply them to Episode 49's handbook, find where it breaks, and fix it. §04's "the systematics cancel" holds only when both ends have the same functional form and the window is exactly centred; it does not hold in general — (b) and (c) are the counterexamples. The optimization in §05 idealizes "\(n=100\) independent estimates per decade"; in practice log-binning leaves fewer estimates at low frequency. The record lengths in §03 assume a near-white process and are a guide only; real data adds windowing, leakage and pre-processing effects. This article too used only synthetic data and was not applied to any real measurement. — To print, use your browser's Print and "Save as PDF."

Printing / PDF: ⌘+P (Ctrl+P on Windows). "Show answer" opens each solution.