The cheap rank and what it cannot see

Rank is a decision settled what rank means on a machine. There is no rank; there is a list of singular values, a gap somewhere in it with luck, and a threshold somebody chose. The essay’s figure is the spectrum with the threshold drawn across it, and its argument is that the threshold is the content.

What it did not say is that almost nobody computes those singular values.

An SVD — the best approximation there is — costs several times what a factorisation costs — a bidiagonalisation and then an iterative phase, against one sweep of reflections. The standard substitute is QR with column pivoting: eliminate greedily, always taking the column with the largest norm in the part still to be eliminated, and read the rank off |rₖₖ| — a decision, as always. LAPACK’s xGEQP3, MATLAB’s qr(A, 0) with a permutation, scipy.linalg.qr(pivoting=True) — all the same algorithm, and all described as rank-revealing.

What the diagonal is worth, in one line #

There is a guarantee, it is one line of algebra, and it holds for every triangular factor — pivoted, unpivoted, from any algorithm at all.

eₙᵀ R⁻¹ eₙ = 1/rₙₙ, so ‖R⁻¹‖₂ ≥ 1/|rₙₙ|, so σₘᵢₙ(R) = 1/‖R⁻¹‖₂ ≤ |rₙₙ|.

|rₙₙ| is an upper bound on σ_min — the one quantity on this page that a gap could be read from with no arithmetic caveat. The cheap test can report a matrix as further from singular than it is, and never as closer. Its rank verdict can be too high and never too low.

That is the same asymmetry as an estimate that can be fooled, reached by a completely different mechanism, and the two together are the phase’s clearest instance of a pattern: a cheap certificate of ill-conditioning is always a witness, and a witness gives a one-sided bound.

assertTheDiagonalIsAnUpperBoundOnSigmaMin checks it across forty Gaussian matrices, three planted ranks and three Kahan matrices, with the worst overestimate in the sample at 2.7·10⁴.

And it works #

The essay would be dishonest without this half. On a matrix with a genuine gap, the diagonal of R finds it.

assertTheCheapTestFindsAPlantedGap runs planted ranks 2 through n − 2 at three shapes, and the two verdicts agree every time, through a gap in the diagonal of at least 10⁵. The generator is greedy, the greedy choice is the right one on a matrix with a real gap, and the whole thing works.

What the counterexample removes is the word guarantee.

Kahan’s matrix #

Let s² + c² = 1 and build the upper triangular matrix with sⁱ down the diagonal and −c·sⁱ everywhere above it:

K = diag(1, s, s², …) · [ 1  −c  −c  ⋯ ]
                        [    1  −c  ⋯ ]
                        [        1  ⋯ ]

Now compute the norm of column j restricted to rows k and below, which is exactly what the pivot rule looks at at step k:

‖K(k:, j)‖²  =  s²ʲ + c²(s²ᵏ + s²ᵏ⁺² + ⋯ + s²ʲ⁻²)
             =  s²ʲ + s²ᵏ − s²ʲ
             =  s²ᵏ

Every trailing column norm is exactly s²ᵏ, for every j ≥ k, at every step. They are all equal. The greedy rule is choosing between numbers that are the same number.

So it makes no interchange, R is K unchanged, and the diagonal ends at sⁿ⁻¹. At n = 50 and c = 0.5 that is 8.7·10⁻⁴, and σ_min is 3.7·10⁻¹²: a factor of 2.3·10⁸.

Not a bad choice — no choice #

Three things are asserted about that, because each rules out an explanation a reader will reach for.

The rule is not choosing badly. assertKahanDefeatsColumnPivoting checks that the interchange count is zero. There was no better column at any step, because there was no different column at any step.

The factorisation is not inaccurate. ‖AP − QR‖/‖A‖ is at rounding. This is not a case of cancellation destroying the diagonal; the diagonal is exactly what the matrix put there.

And there is no threshold that would have worked. Every consecutive ratio in the diagonal is the same s, to within 10⁻⁴ across the whole of it, so the diagonal has no gap in it anywhere for a threshold to sit in. The refusal a threshold offered as the repair for a diagonal with no gap in it is fed exactly that claim and fails.

And it has no floor #

A single counterexample is a curiosity. The statement worth making is that the failure is unbounded.

At c = 0.5 the ratio runs 1.2·10³, 7.0·10⁴, 4.0·10⁶ and 2.3·10⁸ at n = 20, 30, 40 and 50. Along the other axis, at n = 40, it runs 6.4·10² , 6.5·10⁴ and 4.0·10⁶ as c goes 0.2, 0.35, 0.5. Both directions are asserted to be monotone.

The verdict a caller receives #

Ratios become decisions when a threshold is applied, and the decision is where the argument lands.

Put the threshold in the middle of the gap between σ_min and σₙ₋₁ — the only place it could go if the singular values were known — and the singular values report rank 49 out of 50 while the diagonal of R reports 50. The two routines disagree about whether the matrix is singular, and neither has made a mistake.

assertTheTwoRoutinesDisagreeAboutTheRank states it in exactly that form: the threshold is checked to be inside the gap, and the two counts are checked to differ.

The second finding: a strict comparison of two equal numbers #

The trailing column norms are equal in exact arithmetic. Once computed, they differ by an ulp.

So a pivot rule written with a strict > — which is the obvious way to write it, and the way it is written in every textbook — is comparing two rounding errors and calling the answer a pivot choice.

The strict rule does better. At n = 50 it returns a ratio of 51.6 instead of 2.3·10⁸.

It is not a repair, and the count is the evidence. Fourteen interchanges at n = 40 and six at n = 50 — not monotone in the size, because it is not a function of the size. It is a function of which ulp fell where, so it does not survive a different summation order, a different blocking, or a different precision. LAPACK, which downdates its column norms rather than recomputing them, is a different summation order.

assertTheStrictRuleReordersOnRounding asserts the non-monotonicity directly, and the refusal a rounding-driven interchange count read as a pivot strategy is fed the claim that the count grows with the size.

This is the site’s absolute-tolerance rule turned around, and it is the same asymmetry the swap that is not optional prices for pivoting: a comparison written without a tolerance is a claim about the scale of what is being compared. The steady-06 phase found that an absolute epsilon inside a tick loop is a claim about the scale of the data. Here a zero tolerance inside a comparison is a claim that two computed quantities can be meaningfully ordered, and on this matrix they cannot.

Three summations, one comparison, eight orders #

The paragraph above says the strict rule’s improvement will not survive a different summation order and names LAPACK’s downdating as an example. That is a prediction, so it was run. Three ways of computing the same trailing column norms — sum the squares from row k downwards, keep a running squared norm and subtract the eliminated entry (which is what LAPACK does), or sum from the bottom row upwards — with one strict comparison, on Kahan’s matrix at n = 50, c = 0.5:

Every one of those is a correct factorisation — the residual is at rounding in all three — and the rank test they support spans eight orders of magnitude. The 51.6 quoted above is the middle of the three. Summing the same squares from the bottom up instead of the top down restores the whole of the failure the strict rule was said to repair; downdating them lands within a factor of three of the truth, which is better than any of this essay’s other routes and is not a property of the algorithm at all.

Two further readings come off the same table.

The strict ratios do not move with the size. 51.58 at n = 30, 40 and 50 for the forward recompute; 3.019 at all three for the downdate. Meanwhile the tie-tolerant ratio grows like sⁿ, from 7.0·10⁴ to 2.3·10⁸. So the rounding-driven interchanges are not a small perturbation of the failure — they put the routine in a different regime, one whose answer happens not to degrade with n, which is exactly the property somebody would report as rank-revealing if they had measured only that column.

And the tie tolerance is not summation-independent either. Under downdating it still makes three interchanges at n = 40 and thirteen at n = 50, because the downdated norms drift apart by more than the 10⁻¹² the tie test allows. The repair for a rounding-driven comparison is itself decided by how the numbers were accumulated, which is the second-order version of the same defect and is the one a library would ship.

There is a practical instruction in it that the essay’s earlier sections cannot give. A caller who needs a rank verdict they can defend has three options and the table orders them: take the SVD, which costs more and answers the question asked; take |rₙₙ| and read it as the upper bound it provably is, which never claims a matrix is worse than it is; or read a pivoted diagonal as a rank and accept that the answer is a function of the library’s blocking. The third is what every default does, and nothing in any of the three interfaces distinguishes them.

What survives all of it is the theorem. |rₙₙ| is an upper bound on σ_min whatever the pivot rule did, and it is a bound for every one of the six factorisations above. The bound is the only thing on this page that does not depend on which ulp fell where. assertTheStrictRulesAdvantageIsAnArtefactOfTheSummation holds the residuals, the eight-order spread, the zero interchanges under reversal, the constancy in n of both strict ratios, and the tie tolerance’s failure under downdating.

Where the rank verdict is actually used #

It would be easy to read all of this as being about a diagnostic printed and forgotten. The verdict is usually an input to something.

A least-squares solve with a rank-deficient design. xGELSY runs column-pivoted QR, decides the rank from the diagonal, and solves the reduced problem. Overstate the rank and the solve inverts a direction that is not there, amplifying noise by the reciprocal of a singular value that should have been discarded. The routine’s own documentation offers xGELSD, which uses the SVD, as the alternative, and the choice between them is a choice between speed and this essay.

A rank-deficient constraint set. An optimisation code that detects redundant constraints does it by factorising the constraint Jacobian and reading off a rank. Overstating it keeps a redundant constraint, which makes the KKT system singular in a way the next factorisation has to deal with — and that is when symmetry is not enough inheriting a decision made here.

Model-order reduction and subset selection. Column-pivoted QR is the standard way to choose k columns of a matrix to keep, and the pivot order is the answer, not a diagnostic. On a matrix where the rule has nothing to choose between, the columns it keeps are the first k, which was not a choice.

In all three, the failure mode is the same shape: the rank comes out too high, so nothing is discarded, so a direction that carries no information is kept and inverted. The consequence appears later, as an answer with a large component along a direction the data never determined — which is where the answer stops being in the data arriving by a route that essay does not mention.

A ratio at one size is a number; a ratio across sizes is a rate, and the rate is what decides whether this family is a curiosity or a hazard:

Ten more again, to see whether 57.7 is a coincidence of that particular pair of sizes — and this is the reason to draw the sweep rather than quote its endpoints, because two points cannot tell an exponential from anything else.

Three identical multipliers is a law rather than a trend. The last frame runs it to the largest size the sweep affords:

Fifty-seven point seven, three times. The ratio multiplies by exactly 57.7 for every ten unknowns added — 1,214 to 7·10⁴ to 4.04·10⁶ to 2.33·10⁸ — and the last step, fourteen unknowns rather than ten, gives 292, which is 1.486¹⁴ to three figures. So the underestimate is exponential in the size with a rate of 1.486 per unknown, and it is exponential to the digit rather than approximately.

And the interchange count is zero at every size. That is the sharper half of the essay’s title. The cheap rank estimate is not making a bad choice among the columns; it is declining to choose, because at every step the remaining column norms are equal to working precision and the pivoting rule has nothing to break the tie with. A test that never fires is not a test that is sometimes wrong — it is a test whose failure mode is invisible from inside it.

Which is why the size matters so much here. At twenty unknowns the estimate is out by three orders of magnitude, which a careful caller might notice; at sixty-four it is out by eleven, which no tolerance would survive. The same code, the same rule, the same zero interchanges, and an error that grows by a factor of a hundred billion across a range of sizes any real problem passes through.

What it costs to be sure #

The honest accounting, because “just use the SVD” is the obvious response and it is not free.

A Householder QR of an m×n matrix costs about 2mn² − 2n³/3. Column pivoting adds the norm bookkeeping, which is O(mn) if the norms are downdated and O(mn²) if they are recomputed — the implementation here recomputes, deliberately, so that what is measured is the pivot rule rather than one implementation’s cancellation in the downdate.

An SVD of the same matrix costs the bidiagonalisation, about 4mn² − 4n³/3, plus an iterative phase. Call it two to three times the pivoted QR in practice, and more if singular vectors are wanted.

So the substitution buys a factor of two or three. That is worth having on a large problem and it is not worth having on a small one, and the sentence a reader can use is: if the matrix is small enough that the difference would not be noticed, compute the singular values. The cases where the cheap route matters are exactly the cases where checking it is expensive, which is an uncomfortable place to leave it and is where it is.

What is actually rank-revealing #

Strong rank-revealing factorisations exist — Gu and Eisenstat’s is the standard one — and they work by doing what column pivoting does and then checking: after the greedy pass, look for a column swap that would increase |det R₁₁| by more than a factor f, and perform it if there is one. Iterate until none exists.

That gives a genuine bound, σₖ / (|rₖₖ|) ≤ p(n, k) with p polynomial, and it costs the greedy pass plus however many corrective swaps are found. On most matrices there are none and the cost is the greedy pass; on Kahan’s matrix there are many.

The distinction is exactly the one this essay is about. Column pivoting is a heuristic that usually finds the gap; Gu–Eisenstat is an algorithm that verifies it found one. The word “rank-revealing” is attached to both in common usage, and only the second has a theorem.

What LAPACK’s downdate changes #

One implementation detail is worth separating from the argument, because it is the most common objection to this essay’s second finding.

xGEQP3 does not recompute the trailing column norms at every step. It downdates them: having eliminated one entry from a column, it subtracts that entry’s square from the stored norm rather than re-adding the rest. That is O(n) a step instead of O(mn), and it is why the routine is affordable.

It is also a subtraction of two nearly equal quantities, which on a matrix whose columns shrink steadily loses relative accuracy exactly where the accuracy matters. LAPACK guards against it by recomputing whenever the downdated norm falls below a fraction of the original — the standard trick, and one more absolute-versus-relative tolerance of the kind this site has now met four times.

So the ulps this essay’s second finding is about are, in the shipped routine, a different set of ulps: they come from the downdate rather than from the sum. That changes which interchanges get made and does not change the conclusion, which is that the comparison being made is between numbers that are equal, and its outcome is therefore a property of the arithmetic rather than of the matrix.

qrColumnPivoted here recomputes, deliberately, so that what is measured is the pivot rule and not one implementation’s cancellation in its bookkeeping.

The other direction of the bound #

|rₙₙ| ≥ σₘᵢₙ is the half this essay is about. There is a bound the other way, and knowing what it costs is what makes the first one worth stating.

For column-pivoted QR the classical result is σₖ(A) ≤ |rₖₖ|·2^k at worst — an exponential factor, and Kahan’s matrix is the construction showing it is not merely an artefact of the proof. So the two-sided statement is: |rₙₙ| is above σ_min and within an exponential factor of it, which brackets nothing useful.

Gu and Eisenstat’s strong rank-revealing factorisation replaces the exponential with a polynomial, at the cost of the corrective swaps described above. That is the whole difference between a heuristic and an algorithm here: not the greedy pass, which both perform, but whether anything checks the result.

What is worth carrying #

The diagonal of a pivoted R bounds σ_min from above, always, and the bound can be off by any factor. The direction is structural and the size is not bounded by anything.

Kahan’s matrix does not defeat the rule by being adversarially pivoted. It defeats it by making every choice identical, so the greedy rule has nothing to be greedy about, and no gap appears in the diagonal for a threshold to find.

A strict comparison of two quantities that are equal in exact arithmetic is a comparison of their rounding errors, and an improvement obtained that way is not reproducible across implementations even when it is real.

And “rank-revealing” names two different things. One of them has a theorem.

The cheapest verdict of all, and what it is worth #

This page’s diagonal is a cheap certificate that errs in the flattering direction. There is a cheaper one still — the determinant — and it does not err in a direction at all.