The absorption-shape test — does a ratio spectrum have a gas's shape?

Write $C_n = F_0 \cos\theta_s\, T_{\text{atm}}^{m} / \pi$ for the factor that pair $n$’s two pixels share. They share it when the out-pixels lie near the in-pixels, so both pixels see the same sun, elevation and water vapour; different pairs need not share it. Leaving out the noise,

with $E_n = E(\lambda; \Delta X_n)$. This is an identity, not an approximation, and it says three things.

$\eta$ is a ratio of brightness-weighted mean reflectances. Pairing on bands that methane does not absorb makes $\rho_n \approx \rho'_n$ pair by pair outside methane’s band. Pairing cannot see a mismatch that exists only inside methane’s band, such as a carbonate’s dip at 2.3 μm, so that mismatch stays in $\eta$, and it is the fit and the score that must reject it.

$\bar E$ is a mean of transmittances, not the transmittance of the mean column. $E$ is decreasing and convex in $\Delta X$, so Jensen’s inequality gives $\bar E \ge E(\overline{\Delta X})$ with $\overline{\Delta X} = \sum_n w_n \Delta X_n / \sum_n w_n$. The column that reproduces $\bar E$ is therefore at most $\overline{\Delta X}$: an enhancement whose column varies across the pairs is fitted with a smaller column than its weighted mean.

The ratio is taken after averaging. For independent noisy quantities $a$ and $b$ with means $\mu_a, \mu_b$ and relative standard deviation $s_b$ of $b$, a second-order expansion of $1/b$ gives $\mathbb{E}[a/b] \approx (\mu_a / \mu_b)(1 + s_b^2)$. Averaging $N$ pairs first divides $s_b^2$ by $N$; dividing each pair and then averaging leaves the bias at the size of a single pixel’s $s_b^2$. The bias is largest for dark pixels, whose $L_{\text{out}}$ has the lowest signal-to-noise ratio.

Band $b$ averages the wavenumbers $\nu$ of its response with a weight $w_b(\nu) \ge 0$, normalised so that $\int w_b\, d\nu = 1$, that combines the response, the solar irradiance and the ambient transmittance. With $\sigma(\nu)$ methane’s cross-section in m²/mol at the enhancement’s layer,

The exponential is convex, so Jensen’s inequality gives

with equality only if $\sigma$ is constant across the band. Physically, the line cores saturate first: once $e^{-m\sigma\Delta X} \approx 0$ at a line’s centre, more methane removes no more light there, and the band darkens only through the line wings. A model $e^{-\varepsilon_b \Delta X}$ predicts more absorption than there is for every $\Delta X > 0$, so it reaches the observed transmittance with a smaller column, and the deficit grows with the enhancement’s strength.

Taking logarithms does not remove the effect. $\ln E_b(\Delta X)$ is the cumulant generating function of $m\sigma$ under the weight $w_b$, evaluated at $-\Delta X$. That function is convex, and it is linear only when $\sigma$ is constant across the band. A model that is linear in $\Delta X$ in log-radiance, such as the matched filter on $\ln L$, has the same low bias for strong enhancements.

Pairing uses only bands that methane does not absorb. If methane’s bands took part, pairing would favour out-pixels that happen to be dark there and cancel part of the absorption. The Euclidean distance on radiance matches brightness as well as shape, so it makes $\eta$ close to one and not only flat.

One to one. If one out-pixel served many in-pixels, its noise and any surface mismatch would enter many pairs at once. The injective assignment (7) is a linear assignment problem, solved exactly by the Hungarian algorithm (scipy.optimize.linear_sum_assignment).

The best-matched half. Dropping the worst-matched pairs trades a factor of $\sqrt{2}$ in noise for a smaller and smoother surface mismatch in $\eta$. The reference code does not give its reason for choosing a half.

Out-pixels free of the gas. An out-pixel that carries the enhancement too cancels the signal being tested. Choosing out-pixels with a small matched-filter score, below the window’s mf_background_mol_m2, keeps them out; that choice is the caller’s.

Let $V \in \mathbb{R}^{B \times K}$ hold the Legendre polynomials at the window’s bands, $V_{\lambda k} = P_k(x_\lambda)$ with $x_\lambda = 2(\lambda - \lambda_{\min}) / (\lambda_{\max} - \lambda_{\min}) - 1$. For fixed $\Delta X$, the model in (8) is $A(\Delta X)\, c$ with $A(\Delta X) = \operatorname{diag}(E(\Delta X))\, V$, so

with $\Pi_A = A (A^\top A)^{-1} A^\top$ the orthogonal projector onto the columns of $A$. The grid point $k$ of least $\mathrm{RSS}$ needs no starting point: it is the global minimum over the grid. $\widehat{\Delta X}$ then refines it between the nodes: it is the vertex of the parabola through $\mathrm{RSS}$ at three neighbouring nodes around $k$ (the first or last three at a grid end), clipped to their span; $E(\widehat{\Delta X})$ interpolates $\ln E$ linearly between the two rows around it, and $\hat c$ is solved again there. The refined fit replaces the node’s only where the parabola is convex and the refit lowers the residual, so a fit exact on a node stays there. If $\mathrm{RSS} \approx a u^2 + b u^3$ near its minimum, with $u$ the distance to it, and the nodes are $h$ apart, the vertex is off by $(b / 2a)(3\delta^2 - h^2)$, with $\delta$ the middle node’s distance to the minimum: at most $|b| h^2 / a$, against $h / 2$ for the node alone. That matters for weak enhancements: on a 0.005 mol/m² grid, a noise-free 62 ppm·m (0.0025 mol/m²) ratio would otherwise sit on a node and score $D_{\text{norm}} \approx 0.5$; refined, it scores below 10^-3. The grid starts at 0, which imposes $\Delta X \ge 0$: a ratio whose best fit wants no methane, or less than none, sits on that boundary at $\widehat{\Delta X} = 0$.

Least squares is the maximum-likelihood fit when the errors in $R$ are independent across bands and have equal variance. Neither holds exactly. The noise in $R$ is larger in bands with less light, and the leftover surface mismatch in $\eta$ is smooth in $\lambda$, so it is correlated across bands.

Linearise for a weak enhancement in log space, where the continuum is additive: $\ln R \approx V a - \varepsilon\, \Delta X + n$, treating $\ln P$ as lying in the span of $V$ and taking $n \sim \mathcal{N}(0, \sigma^2 I)$. Profiling out $a$ (the Frisch–Waugh–Lovell theorem) leaves

Only the part of methane’s band shape that a polynomial of degree $K - 1$ cannot reproduce informs $\widehat{\Delta X}$. Each added degree removes another direction from $\varepsilon_\perp$ and inflates the variance; with $K = B$ the continuum fits anything, and $\Delta X$ is not identifiable. This is the sense in which eleven continuum terms on 46 bands leave little of methane’s shape to score (7. Where radtrax departs from the reference code).

The basis does not change this. Legendre polynomials and monomials of the same degree span the same space, so they give the same $\varepsilon_\perp$ and, in exact arithmetic, the same fit. What the Legendre basis changes is conditioning: monomials of $\lambda$ near 2300 nm are nearly collinear, and the least-squares solve loses digits to that.

The same algebra says what the test can reject. A confounder that adds $d$ to $\ln R$ shifts $\widehat{\Delta X}$ by $-\varepsilon_\perp^\top d / \lVert \varepsilon_\perp \rVert^2$ and leaves the rest of $(I - \Pi_V) d$ in the residual. The test rejects a confounder only to the extent that its shape, once the continuum is removed, differs from methane’s.

Write $q = E(\widehat{\Delta X}) + r$, with $r$ the residual. The numerator of (9) is $\operatorname{mean}_\lambda |r|$, the part of $q$ the fit leaves unexplained. The denominator is $q$’s spread about its own mean, which is the absorption the fit explains plus the residual. $D_{\text{norm}}$ is an L1 analogue of the fraction of variance unexplained, $1 - R^2$.

It is not purely a test of shape. Let the residual have size $\sigma_q$ per band and the fitted transmittance vary across the window by $\delta$. When $\delta \gg \sigma_q$, $D_{\text{norm}}$ is of order $\sigma_q / \delta$; when $\delta \ll \sigma_q$, numerator and denominator are both about $\sigma_q$ and $D_{\text{norm}}$ tends to 1. A real but weak enhancement therefore scores like none. For $N$ pairs of pixels with signal-to-noise ratio SNR, $\sigma_q \approx \sqrt{2} / (\mathrm{SNR}\sqrt{N})$, so a fixed threshold such as 0.3 also depends on the sensor’s noise, on the number of pairs and on the window, besides the continuum’s degree.

It can exceed 1, when the fitted transmittance has structure that $q$ lacks.

The fit and the score optimise different things. The fit minimises the L2 norm of $R - P E$; the score is an L1 norm on $q = R / P$. $\widehat{\Delta X}$ therefore does not minimise $D_{\text{norm}}$'s numerator, though the difference is small when $P$ is close to 1.