Compensated Summation for GLM Objectives

Why Summation Matters in GLM Computation

A GLM objective is assembled from contributions: observation-level log-likelihood terms are added, score contributions are summed into a vector, and information contributions are summed into a matrix. In exact arithmetic, the order of addition does not matter. In floating-point arithmetic, each addition may round, so changing the order can change the result. The effect is especially noticeable when many terms are accumulated or when large positive and negative contributions nearly cancel.

A small summation error can affect more than the displayed objective. An inaccurate score can make a computed update appear too small or point in a slightly different direction; inaccurate information entries can alter later calculations. Compensated summation reduces the error caused by accumulation. It does not change the GLM or correct errors already present in the individual contributions.

A Bound for Sequential Summation

Use the standard floating-point model for an addition: if the exact result is finite and normal, then its computed value is \(\operatorname{fl}(a+b)=(a+b)(1+\delta)\), where \(|\delta|\leq u\) and \(u\) is the unit roundoff. Assume here that no overflow or underflow occurs. Define \(\gamma_k=ku/(1-ku)\) for \(ku<1\), and let \(\gamma_0=0\).

Proof. Write \(s_i=(s_{i-1}+x_i)(1+\delta_i)\), with \(|\delta_i|\leq u\). Since \(x_1\) is itself a floating-point number, the first addition is exact: \(\operatorname{fl}(0+x_1)=x_1\). Thus \(\delta_1=0\). Expanding the recurrence gives \[ s_n=x_1\prod_{i=2}^n(1+\delta_i)+\sum_{j=2}^n x_j\prod_{i=j}^n(1+\delta_i). \] The coefficient of \(x_1\) therefore has at most \(n-1\) nontrivial factors; the coefficient of \(x_j\), for \(j\geq 2\), has \(n-j+1\leq n-1\) factors. This factor count uses the exact first addition; without \(\delta_1=0\), the coefficient of \(x_1\) would include another factor.

For any product of \(k\) factors \(1+\delta_i\), with \(ku<1\), its distance from 1 is at most \(\gamma_k\). Indeed, the product is at most \((1+u)^k\leq 1/(1-ku)\), and at least \((1-u)^k\geq 1-ku\). These bounds imply that the absolute difference from 1 is at most \(ku/(1-ku)\). Since \(\gamma_k\leq\gamma_{n-1}\), subtracting the exact sum from the expanded expression and applying the triangle inequality gives the stated bound. \(\square\)

Compensated Summation

Compensated summation keeps track of low-order information that ordinary addition would discard. One practical variant is Neumaier summation. It maintains a running sum \(s\) and a correction \(c\). For each input \(x\), it computes \(t=\operatorname{fl}(s+x)\), estimates the part lost in forming \(t\), and adds that residual to \(c\). A final result is formed from \(s+c\).

The branch chooses an ordering that helps estimate the rounding residual of the addition; \(c\) accumulates those residual estimates. In common floating-point settings this can substantially improve accuracy, particularly when small terms are added to a much larger running sum. It is not a promise of an exact answer for every input: the residual and the correction are themselves computed in floating-point arithmetic.

What the Condition Number Tells You

For inputs with exact sum \(S=\sum_i x_i\neq 0\), define the summation condition number by \[ \kappa_{\mathrm{sum}}=\frac{\sum_i |x_i|}{|S|}. \] It is at least 1 by the triangle inequality. It is close to 1 when the terms mostly reinforce one another and large when they largely cancel.

Proof. Divide the theorem’s absolute error bound by \(|S|\), which is nonzero. The resulting right-hand side is \(\gamma_{n-1}\sum_i|x_i|/|S|=\gamma_{n-1}\kappa_{\mathrm{sum}}\). \(\square\)

This result separates two sources of difficulty. The factor \(\gamma_{n-1}\) reflects the accumulation process and the number of additions. The condition number reflects the inputs themselves. Compensation can reduce accumulation error, but a sum whose terms nearly cancel remains sensitive to small perturbations in the contributions.

Using Compensation in Practice

Compensation is most useful when many contributions are accumulated, their magnitudes differ substantially, or cancellation is expected. For reproducible computation, also keep the summation order fixed: compensation improves accuracy but does not make floating-point addition associative. In vector and matrix calculations, compensation may be applied entry by entry, using an independent running sum and correction for each entry.

A useful implementation practice is to compare ordinary and compensated totals on the same finite contributions. A discrepancy can indicate that ordinary accumulation discarded meaningful low-order terms. It is also useful to inspect the scale of \(\sum_i|x_i|\) relative to the scale of the final sum. These checks help identify a fragile accumulation without claiming that a particular result is wrong merely because the two totals differ.

Common Mistakes

• Assuming the first addition contributes a rounding factor. In the stated recurrence, the inputs are floating-point numbers, so \(\operatorname{fl}(0+x_1)=x_1\). The first addition is exact; this is why the coefficient of \(x_1\) has at most \(n-1\), not \(n\), nontrivial factors.
• Assuming addition order does not matter. Exact addition is associative, but floating-point addition generally is not. Reordering contributions can change the computed total.
• Treating compensation as a cure for every numerical problem. It addresses loss from summation. It cannot recover information lost when a contribution was computed inaccurately, nor can it prevent overflow in that computation.
• Reading an absolute error bound as a useful relative bound. When the exact sum is near zero, the condition number can be large, so a small absolute error may still be a large relative error.