A good curve and a bad verdict

A column-pivoted QR produces one object, and two entirely different readers take it away. The first reads a single number off it — the last diagonal entry, thresholded, turned into an integer — and calls that integer the rank. The second reads the whole diagonal as a shape: how fast the magnitudes fall, where they cross a tolerance, what would be lost by stopping at column twenty rather than column forty.

The first reader has been dealt with. The cheap rank and what it cannot see shows a triangular matrix on which that last entry sits eight orders of magnitude above the smallest singular value, with no interchange made anywhere and no threshold that would have helped, and shows that the discrepancy grows without limit in the size. The verdict has no floor.

Nothing has asked the second question. It is not the same question, and the answer turns out not to be the same answer: over the whole leading block of the same factorisations, |r_kk| sits within a factor of about two of the σ_k it stands in for, and it does so on the counterexample as well as on the matrices the method was designed for. The number that settles it is the width of the band that the ratio |r_kk| / σ_k occupies, and it is measured below on sixty-four diagonal entries of planted-rank matrices and on the first thirty-nine entries of the matrix the whole subject is built around.

The two curves do not merely share a cliff — they share a level #

The figure above is normally read for its cliff. Both curves drop eleven orders of magnitude at the same index, so any threshold dropped into the gap gives the same integer from either, which is the whole of what “rank-revealing” is usually taken to mean. That reading throws away almost everything the picture contains, because it uses two points out of eighty.

Read instead every entry of the leading block. On the matrix above the twelve ratios |r_kk| / σ_k run 0.576, 0.541, 0.551, 0.536, 0.544, 0.560, 0.560, 0.567, 0.524, 0.518, 0.506 and 0.431. The diagonal is not tracking the singular values loosely; it is tracking them at a nearly constant fraction of about a half, and the fraction wanders by a factor of 1.34 across the entire block.

That is a different kind of statement from the one about the cliff, and it holds at every planted rank the generator draws. At planted rank 4 the band is 0.407 to 0.515; at 8, 0.465 to 0.533; at 12, 0.431 to 0.576; at 16, 0.553 to 0.639; at 24, 0.489 to 0.775. Sixty-four diagonal entries in total, across five matrices with ranks spanning a factor of six, and every one of them lies in [0.407, 0.775]. The worst single entry understates its singular value by a factor of 2.456 and the best by 1.29.

A factor of 2.46 across sixty-four entries is not a tight bound in the sense that a numerical analyst uses the word — there is no theorem here, and the classical worst case for pivoted QR is exponential in the index. It is, however, two significant figures of the spectrum obtained from a factorisation that costs a fraction of the one that would produce the spectrum exactly. That is the claim this essay is about, and it is a claim about a curve rather than about a gap.

The band drifts upward with the rank, and does not widen #

Planted rank 4 is the worst case in the sample and planted rank 24 is the largest block, and the comparison between them says something about how the approximation behaves as more of the matrix is being described.

The ratios in the larger block sit higher — around 0.7 rather than around 0.5 — and the spread inside the block is 1.58 rather than 1.26. Neither number moves much, and the direction of the drift is worth naming: a longer leading block puts the diagonal closer to the singular values, not further from them. The greedy pivot rule picks the largest remaining column norm at every step, and a column norm is an upper bound on nothing and a lower bound on the singular value it is competing with; the more genuine directions there are left to choose from, the better that choice is.

The one place the ratio falls away inside a block is its final entry: 0.431 at index 12 of the twelve, 0.489 at index 24 of the twenty-four, 0.515 at index 4 of the four. The last column before the cliff is the one the rule has the least freedom about, because everything it might have swapped in has already been taken. That is the first appearance of a pattern the rest of the essay is about — the diagonal’s disagreement with the spectrum concentrates at the ends of the block, and the middle is where it is reliable.

The counterexample, read along its whole length #

Kahan’s matrix is the construction that removes the word guarantee from all of this. Its trailing column norms are exactly equal at every step, so the greedy rule has nothing to be greedy about, it makes no interchange at any size, and the diagonal it produces is the matrix’s own diagonal unchanged. At n = 40 and c = 0.5 that diagonal ends at 0.00366 while σ_min is 9.07·10⁻¹⁰: a ratio of 4.04·10⁶, and the essay that priced it shows the ratio multiplying by 57.7 for every ten unknowns added.

Now measure the ratio at every index rather than at the last one. Over the first thirty-nine indices the band is [0.173, 0.816]. Discard index 1, which the next section shows is a separate and bounded effect, and the band over indices 2 through 39 is [0.669, 0.816] — a factor of 1.221 across thirty-eight consecutive entries, tighter than anything measured on the planted matrices.

That is the finding. On the matrix that is the standing counterexample to the whole method, the pivoted diagonal reproduces the singular value curve to within twenty-two per cent over ninety-five per cent of its length, and is wrong by six orders of magnitude at exactly one index. The failure is not diffuse and it is not a degradation. It is a single entry, and it is the last one, and the last one is the only place a rank verdict looks.

The reason is visible in the picture and is a property of the matrix rather than of the rule. The diagonal is a pure geometric sequence — sⁿ⁻¹ down the whole length, with s = 0.866 — and the singular values follow that same geometric decay until the very end, where the accumulated off-diagonal mass of the construction collapses the smallest one far below the trend. The diagonal has no way to report a departure that happens after it has stopped, and the departure happens in the last singular value.

Only one end of the curve is unbounded, and the other is a theorem #

Index 1 is the other place the ratio leaves the band, and it leaves it in a way that is bounded in advance.

|r_11| is the largest column norm of the matrix, because that is what the first pivot selects. A largest column norm is at most ‖A‖₂, since every column norm is, and it is at least ‖A‖_F / √n, since the Frobenius norm is the root of the sum of the n squared column norms and the largest is at least the root mean square. So σ_1 / √n ≤ |r_11| ≤ σ_1 for every matrix, with no assumption about pivoting, conditioning or the arithmetic.

Kahan’s matrix nearly attains the lower half of that. At n = 40, c = 0.5 the measured ratio at index 1 is 0.1730 against a bound of 1/√40 = 0.1581 — within nine per cent of the worst any matrix could do. At n = 60 it is 0.1369 against 0.1291, within six per cent. At n = 40 with c = 0.8 it is 0.1619 against 0.1581, within two and a half per cent. The construction that defeats the rank verdict also drives the first entry of the curve as low as the algebra allows, and the algebra allows only √n.

Moving c moves the last entry by four orders of magnitude and moves the band by a fifth. Moving the size does the same thing in the other direction:

Across the whole family measured — sizes 20, 40 and 60, parameters 0.2, 0.5 and 0.8 — the ratio over the first n − 1 indices stays inside [0.137, 0.913], while the ratio at index n runs from 1.2·10³ to 1.3·10¹⁰. Two quantities computed from the same factorisation, one of which moves by a factor of seven across the family and one of which moves by seven orders of magnitude.

What reads the curve is a truncation, and it is not asking for an integer #

The reason any of this matters is that the rank verdict is the minority consumer. Rank is a decision settles what a rank means on a machine and its answer is that there is no rank, only a spectrum and a threshold somebody chose. Most of the code that factorises a matrix and looks at the diagonal is not choosing that threshold at all. It is asking how much of the matrix it can throw away for a stated error, which is a question about the shape of the decay and not about where it ends.

So the honest test of the curve is a truncation error rather than an integer. Stop the pivoted factorisation after k columns, form Q₁R₁ from what has been computed, and measure how far that is from the permuted matrix in the 2-norm. The best approximation there is fixes what the answer would be if the singular values had been computed instead: the error of the best rank-k approximation is not bounded by the next singular value, it is exactly equal to it.

Measured against that standard, the pivoted truncation is close to optimal almost everywhere. On the planted-rank-12 matrix the ratio of the pivoted error to σₖ₊₁ runs between 1.041 and 1.192 for every k from 1 to 11. On a 20 × 20 Hilbert matrix — a smooth decay with no gap anywhere and fifteen interchanges — it runs between 1.102 and 2.472 for k from 1 to 13. And on Kahan’s matrix at n = 40, c = 0.5, the matrix on which the rank verdict is out by four million, it runs between 1.000 and 3.808 for every k from 1 to 38, and reaches 4.04·10⁶ at k = 39.

At n = 60 the same measurement gives 1.000 to 4.813 across the first fifty-eight cuts and 1.34·10¹⁰ at the fifty-ninth; at c = 0.2, 1.000 to 3.089 across thirty-eight and 640 at the thirty-ninth. The blow-up is at the last cut in every case, and it is the same number as the rank verdict’s error because at k = n − 1 the truncation error is |r_nn| and the optimum is σ_min.

The practical reading is blunt. A caller who truncates a Kahan matrix at any rank they would plausibly want gets an approximation within a factor of four of the best possible one, on the matrix that the literature holds up as the reason not to trust this factorisation. A caller who asks the same factorisation whether the matrix is singular gets an answer that is wrong by four million. Same object, same code, same run.

A rank chosen from a tolerance moves by two indices, not by six orders #

Between those two consumers sits a third, and it is the common one in practice: choose k as the number of diagonal entries above a tolerance. That is a rank, so it looks like the first consumer, but the tolerance is usually an accuracy budget rather than a singularity test, so it behaves like the second.

Sweeping forty-five tolerances from 10⁻¹⁵ to 0.5 and comparing the count from the diagonal against the count from the singular values makes the distinction quantitative. On the Hilbert matrix the two counts differ at four of the forty-five tolerances and never by more than one index. That is what a factor-of-two error in a curve buys when the curve is steep: nothing, almost always.

On Kahan’s matrix at n = 40, c = 0.5 the two counts differ at twenty-seven of the forty-five, and the disagreements split cleanly into the two regimes this essay has been separating. At every tolerance from 10⁻⁹ up to 2·10⁻³ — twenty of them — the diagonal says 40 and the singular values say 39. That is the rank verdict, off by exactly one index, and the one index is the last. At every tolerance coarser than that the diagonal says fewer than the singular values do: 37 against 39, 33 against 35, 28 against 30, 21 against 24, 17 against 19, 12 against 14, and 5 against 8.

The second regime is the factor of two showing up as an index shift, and its size is predictable from the decay rate. The diagonal falls by 0.866 per index, so a magnitude ratio r costs ln r / ln(1/0.866) indices; the median ratio over the first thirty-nine entries is 0.688, which is a factor of 1.454, which is 2.60 indices. The measured shifts are two and three. A curve that understates the spectrum by a constant fraction crosses a horizontal tolerance a fixed number of indices early, and on a curve this steep that number is small.

The two errors are therefore of completely different kinds and it is a mistake to quote one for the other. At the last index the error is a factor of four million and no tolerance repairs it. In the body of the curve the error is a factor of one and a half and it costs two or three columns of a truncation whose error budget is a factor of a hundred wide.

Where the under-report stops being a rule #

Two limits, because the band measured above is a measurement on three families and not a theorem about matrices.

The under-report is not a property of the curve. On planted matrices, on Hilbert and on Kahan’s, |r_kk| sits below σ_k almost everywhere — but on 120 Gaussian 12 × 9 matrices the ratio spans [0.627, 1.987] and exceeds one at 465 of the 1,080 diagonal entries. The diagonal of a pivoted R routinely overstates an interior singular value. It cannot overstate the last one, and that is the only place the one-sidedness is a theorem: eᵀ_n R⁻¹ e_n = 1/r_nn gives |r_nn| ≥ σ_min for every triangular factor whatever the pivoting did. Across the same 120 Gaussian matrices the last entry sits above σ_min every time, with the smallest ratio in the sample at 1.087 and the median at 1.454, and the refusal a one-sided rank test claimed to err in both directions is fed the claim that some of them go the other way and is required to fail.

And below the arithmetic’s own floor neither curve means anything. The Hilbert band above stops at k = 13 for a reason. That matrix has ‖A‖₂ = 1.91, so a singular value below about 4.2·10⁻¹⁶ is below the rounding of the largest entry; σ₁₄ is 1.6·10⁻¹⁷ and σ₂₀ is 3.9·10⁻¹⁹. The measured truncation error stops falling at 8.0·10⁻¹⁶ and stays there, so the ratio of measured error to σₖ₊₁ climbs to 2.06·10³ at k = 19 — not because the factorisation degraded, but because the denominator is a number the format cannot represent a difference at. That is the same boundary deciding that a zero has arrived draws around detection thresholds, arriving here as a limit on what a ratio can be computed about.

What follows for anything that reads a decay #

A hierarchical block’s rank is chosen from exactly this curve. A block nobody can call sparse is the shape: ninety-six singular values, five that matter, and a factor of forty a column. The number wanted there is where the decay crosses an accuracy budget, which the measurements above put within two or three columns of the SVD’s own answer, and the cost of being three columns high is three columns of storage rather than a wrong answer.

A rank stated as a number of digits is the curve and nothing else. A rank that is a number of digits measures 0.55 columns a decade for a kernel block; a factor of 1.5 in the magnitudes is 0.1 of a decade, which is one twentieth of a column. The error the verdict cannot survive is invisible on that axis.

A null result that depends on the curve keeps working. The size the rank does not notice turns on a rank being the same number across a factor of eight in sampling, and a constant multiplicative offset in the diagonal shifts every one of those counts identically or not at all. The same holds for the growth in the kernel with nothing to compress, where the quantity of interest is a rate rather than a level.

A truncation performed inside an arithmetic is unaffected. The rounding that was not the problem does ninety-eight truncations inside one formatted Cholesky and finds their cost below the representation error. Each of those truncations is a cut on a decay curve, and a two- or three-column shift in where the cut falls changes the storage and not the residual.

A singularity test is where none of this transfers. An estimate that can be fooled is the other cheap certificate whose error is always in the flattering direction, and it fails the same way: a quantity that is usually right and one-sidedly wrong is safe when it is priced as an approximation and unsafe when it is thresholded into a yes.

A verdict that already depends on how the sums were ordered is not made safer by a good curve. A rank that depends on the thread count gets seven different numerical ranks from one matrix and one threshold, and the gap decides the eigenvector shows the same threshold sensitivity in a neighbouring quantity. Both are about the last index. Neither is helped by the diagonal’s accuracy sixteen entries earlier.

And the distinction is not the one the word “rank-revealing” makes. Gu and Eisenstat’s factorisation earns the name by checking after the greedy pass, which is what buys a polynomial bound rather than an exponential one. What the measurements here separate is different and cheaper: the same greedy pass, unchecked, is a two-figure picture of a spectrum and a six-order-of-magnitude error about a rank. The word has been attached to the object, and the object does two jobs.