A block nobody can call sparse
Worth reading first: Rank is a decision · The factor is not sparse · An operator with no entries.
This collection has two answers to the sentence the matrix is too large to store, and both of them are things the problem hands over rather than choices anybody makes.
The first is sparsity. Most of the entries are zero, so keep the ones that are not, and the whole
sparsity field follows from what that costs at elimination. The
second is not forming it at all: an operator can be a function, y = A(x), with
no entries anywhere, and the essay on that ends where
this one begins — an operator that can only be applied is an operator that cannot be factorised, which is why
every solver in it is an iteration.
There is a third, and it is the one this field is about. It applies to matrices that are not sparse and not small, it is the only one of the three whose cost is chosen rather than handed over, and the choosing is the whole subject.
The matrix
Take 96 points spread evenly over the interval [0, 1] and another 96 over [2, 3], and build the matrix whose (i, j) entry is
1 / |xᵢ − yⱼ| .
Every entry is nonzero. The smallest is about 1/3 and the largest about 1, so there is not a single one that could be dropped and no threshold at which the picture becomes sparse. It is not a Toeplitz matrix either — the points are the same spacing, so in this particular arrangement it happens to be one, but move the points off a uniform grid and it stops being, while everything below survives unchanged. The structure field’s whole subject is a matrix described by fewer numbers than it contains, and this is not that.
So neither existing answer applies. There are 9,216 entries and no reason yet to think fewer than 9,216 numbers are needed.
The measurement
Take its singular values and divide by the largest:
| k | σₖ ⁄ σ₁ |
|---|---|
| 1 | 1 |
| 2 | 2.40·10⁻² |
| 3 | 4.59·10⁻⁴ |
| 4 | 8.47·10⁻⁶ |
| 5 | 1.54·10⁻⁷ |
| 6 | 2.79·10⁻⁹ |
| 7 | 5.01·10⁻¹¹ |
| 8 | 8.97·10⁻¹³ |
| 9 | 1.60·10⁻¹⁴ |
| 10 | 2.84·10⁻¹⁶ |
A factor of about fifty a column, held for nine columns, until it hits the unit roundoff and stops. Fitted over the first five it is 1.71 decades a column.
Five of the ninety-six are above 10⁻⁸ of the first. By the theorem this collection already has an essay about — the best rank-k approximation is the one the singular value decomposition writes down, and its error is the next singular value — a rank-five approximation of this block is wrong in the eighth digit. Two factors of 96 × 5 are 960 numbers, against 9,216 entries.
The control, which is the argument
A curve that falls is not evidence of anything. Every matrix has decaying singular values in the sense that they are sorted, and a reader shown one falling line has been shown a convention.
So the figure carries a second block, of exactly the same size, on exactly the same points, whose entries are independent draws from a normal distribution instead of a function of the distance between two points. Its singular values fall too — from 1 to about 10⁻² over ninety-six columns, which is what a random matrix does — and 96 of its 96 are above 10⁻⁸.
That is the whole comparison. Same size, same density, no zero entries in either, and one of them needs five columns while the other needs all ninety-six. Nothing about the shape of the matrix distinguishes them. What distinguishes them is that one of them is a function of the geometry of two sets of points and the other is not.
Why, and the answer is one line of algebra
The kernel 1/(x − y) can be expanded about the two intervals’ centres. Write x = cx + u and y = cy + v, with d = cx − cy the distance between the centres. Then
1/(x − y) = 1/(d + u − v) = (1/d) · Σₖ (−1)^k ((u − v)/d)^k ,
which converges whenever
q = (rx + ry) / d < 1 ,
with rx and ry the half-widths of the two intervals. Truncate the sum at order p and expand (u − v)^k by the binomial theorem, and every term separates into a function of u times a function of v. That is a rank p + 1 matrix, written down without any linear algebra at all — two tables of powers of a shifted coordinate — and its relative error is bounded by q^(p+1)/(1 − q).
For the two intervals here, rx = ry = 1/2 and d = 2, so q = 1/2 and the expansion converges.
Two things follow and they are worth separating, because the rest of this field lives in the gap between them.
- The block is approximately low rank, and the reason is a property of the geometry rather than of the entries. That is the finding.
- The rank the expansion asks for is not the rank the matrix needs. The expansion loses a factor of two per column and the matrix loses a factor of fifty. Both are geometric; they are not the same geometric.
Two routes, which is why the middle curve is on the figure
The expansion above is a second route to a number the decomposition already gives, and this collection’s standing habit is that a number arrived at by one route has been wrong every time it has been published.
So it is evaluated rather than described. The two tables of powers are built, multiplied together, and subtracted from the block, and what comes out is a relative error of 0.224 at order 0, falling to 8.2·10⁻⁶ at order 12 — a slope of 0.365 decades a column, against the log₁₀(1/q) = 0.301 the bound promises. The expansion behaves as its own analysis says it does, and it never beats the decomposition, which it cannot: one of them is optimal by construction.
The three curves are therefore an ordering and a disagreement. The ordering — bound above expansion above decomposition — holds at every rank measured, and it has to. The disagreement is that the constant is wrong by a factor in the exponent, and that is a measurement rather than an error.
The factor of five is one point on an axis
The two slopes above — 0.365 decades a column for the written expansion, 1.71 for the decomposition — divide to 4.7, and the sentence they were used to write said the bound is loose by a factor of five in rank. That sentence has a hidden argument. Both numbers were measured on one pair of intervals, [0, 1] against [2, 3], where q is exactly 1/2, and a factor is only a summary of a relationship if it holds along the axis the relationship lives on.
The axis here is the separation, and it is not an incidental parameter: the whole business of deciding which pairs of clusters are allowed to be compressed is a rule about where a pair sits on it. So the sweep is the thing to do. Hold the block at 128 square and the accuracy at eight digits, slide the second interval from touching to far away, and compare the rank the four-number bound asks for against the rank the block turns out to need.
| gap | q | rank needed | rank the bound asks for | factor |
|---|---|---|---|---|
| 0.05 | 0.952 | 10 | 440 | 44.0 |
| 0.1 | 0.909 | 9 | 219 | 24.3 |
| 0.25 | 0.800 | 7 | 90 | 12.9 |
| 0.5 | 0.667 | 6 | 49 | 8.2 |
| 1 | 0.500 | 5 | 28 | 5.6 |
| 2 | 0.333 | 4 | 18 | 4.5 |
| 4 | 0.200 | 4 | 12 | 3.0 |
The bound remains a bound at every row, which is the half of it that is a theorem and is worth no less for being unsurprising. What it is not is loose by a constant. Across pairs that are all admissible — every q below 1, every expansion convergent — its looseness runs from 3x to 44x, and it runs monotonically. The factor of five quoted above is neither a typical value nor an upper one. It is the value at the single geometry the figures happen to be drawn on.
The mechanism is visible in the two rates rather than in the two ranks. As the intervals close, q rises to 1 and the bound’s decay rate log₁₀(1/q) goes to zero with it — at q = 0.952 the bound is promising a fifth of a decade every ten columns. The block’s own decay does not go to zero at all: it is still 0.88 decades a column there, and the matrix still compresses to ten columns at eight digits. A barely-admissible block is not a barely-compressible one. The bound thinks it is because the quantity the bound is built out of, the ratio of radii to distance, is the wrong quantity for a spectrum near the touching limit.
Those seven rows are seven spectra, and drawn rather than tabulated they say something the rank column cannot. The same 96 × 96 block, the same four kernels, the same eight-digit cut — only the second interval moves.
The ranks are 10 and 7 against the sweep’s 10 and 7, measured here at 96 square where the table was taken at 128, and the whole column agrees row for row. What is not in the table is the rate. σ₂/σ₁ runs 0.192, 0.143, 0.0846, 0.0489, 0.024, 0.00981 and 0.0034 across the seven separations — the first column of the cliff steepens from a factor of 5.2 to a factor of 294, which is a factor of 56 in the decay rate against a factor of 2.5 in the rank.
That ratio is the answer to why a barely-admissible block is not a barely-compressible one, and it is arithmetic rather than luck: a rank is a number of digits divided by a rate, and a rate is already a logarithm. Eight digits against the 0.88 decades a column fitted at q = 0.952 is nine or ten columns; the same eight digits against the 2.47 decades of the first step at q = 0.2 is three or four. So a factor of fifty-six in the decay buys a factor of two and a half in the rank — and the bound, which is built out of q directly rather than out of its logarithm, moves by a factor of forty-four over the same range.
The control does not move at any of them. The block of independent draws keeps 96 of its 96 singular values above the cut at every separation on the sweep, because there is no geometry in it to separate — which is what makes the four falling curves a measurement of the two intervals rather than of the sorting. A reader who suspects the cliff of being an artefact of the particular pair [0, 1] against [2, 3] can watch the pair move by a factor of eighty in q and watch the control not notice.
That is the direction that matters for a partition. An admissibility rule accepts a pair the moment q drops below its threshold, so the pairs it accepts pile up near the threshold rather than far past it — the comfortably separated blocks are the rare ones. The blocks whose rank a bound is actually asked to predict are the blocks it is worst at predicting, and the essay on what the test costs against what it saves is where that becomes a number rather than an inconvenience.
assertTheBoundsLoosenessIsNotAConstant holds the whole sweep: that the bound never goes under, that
its looseness varies by more than a factor of five across the range, that it tightens monotonically
with the gap, and that at the tightest pair the rank is under twelve while the prediction is over two
hundred. The habit it enforces is the one this collection keeps arriving at from different
directions — a
quantity measured at one setting and quoted as a property
is a quantity nobody has swept. The kernel that
has nothing to compress at all is the other end of the
same warning, and the next essay’s straight line is only straight because somebody moved the knob.
What the third answer actually is
Sparsity is a property of the entries. An entry is zero or it is not, and no decision anybody makes changes how many are.
Matrix-free is a property of the code. Either a routine that applies the operator exists or the code have the entries; the problem decided which before a numerical analyst arrived.
Numerical rank is neither. It is a property of a matrix and an accuracy together, and the accuracy is a number somebody types. Ask for two digits and this block costs two columns; ask for fourteen and it costs nine. That is the sentence this field is organised around, and the next essay is about the shape of that trade — which turns out to be a straight line, so a digit costs the same handful of columns wherever it is bought.
What it does not say
Three things this measurement is not, since each of them is a reading somebody will take.
It is not a statement about the matrix’s rank. The block has rank 96, exactly, in the algebra. The singular values below the tenth are not zero, they are 10⁻¹⁶ — and this collection has an essay about why calling those zero is a decision with two ways to be wrong. Nothing here changes that. What is being claimed is that a rank-five matrix is within 10⁻⁷ of it, and the two sentences are different.
It is not a statement about conditioning. A block whose singular values span fourteen decades is not a well-conditioned object, and a solve with it would be in trouble. Nothing here proposes to. The blocks a hierarchical format compresses are the off-diagonal ones; the diagonal blocks — where the points of the two clusters are the same points — are kept dense precisely because none of this applies to them.
And it is not free. Computing those singular values costs a decomposition of the block, which is cubic. The next four essays in this field are about what to do with the fact; the one after them is about a case where none of it works; and the one on building the representation from products alone is about the fact that the cheapest available route does not compute the decomposition at all.
Where else this shape has appeared here
The site has met a matrix with hidden structure three times before and named the structure differently each time.
The matrix that is one row is a circulant: n numbers, no zero entries, diagonalised exactly by a transform. A limit the matrix never reaches is Toeplitz: 2n − 1 numbers, a symbol, and a spectrum that has a limit it never attains. Both of those are matrices described by fewer numbers than they contain, which is a statement that can be checked by counting.
This is not that, and the difference matters. A kernel matrix is described by exactly n² numbers, each of them the value of a function at a pair of points, and there is no shorter description. What it has instead is that large sub-blocks of it are close to matrices with short descriptions, and close is a word with a parameter in it.
The refusal
The claim under test is the one this essay opens by contradicting: that a matrix with no zero entries has to be stored entry by entry.
It is refuted by the measurement, and it would be refuted by any falling curve, which is why the
refusal is fed the other block. The assertion rank ≤ 12 at 10⁻⁸ is handed the 96 × 96 block of
independent draws on the same points, and it fails — 96 of 96 survive the cut — which is what it must
do. Without that, the essay’s first figure is a picture of a sorted list.
The refusal has a second half worth stating. Compressibility is not a property that a matrix either has or lacks; it is a property a matrix has with respect to an ordering of its indices, and the essay that numbers the same matrix twice takes the same 256 × 256 matrix, renumbers it with one permutation, and finds that it stores 118,208 numbers where the original stored 27,008. Nothing a norm can see moves. The refusal above rules out one wrong reading; that essay rules out the other.
What it costs to find out
There is an awkwardness in everything above and it is worth naming before the field is four essays further on, because a reader who notices it and is not told will assume it was missed.
The five columns were found by taking a singular value decomposition of the block. That costs O(n³) operations and reads every one of the 9,216 entries — so the result is a 960-number object and the route to it touched every entry of the thing it was replacing. On the problems this format exists for, where the matrix is a boundary integral operator and n is in the millions, forming n² entries is the one thing nobody can afford, and a method whose first step is to form them has answered no question at all.
There are three ways out and this field takes all of them.
The first is that most blocks are never formed. A partition splits the matrix so that the compressed blocks are large and few, and the dense blocks are small and near the diagonal; the total number of entries ever visited is the thing that has to be counted, not n².
The second is the expansion above. It writes down a rank p + 1 approximation from four numbers, so a code that knows its kernel analytically can build the factors and never look at the block. That it is loose by a factor of five in rank is the price of not looking, and for a great many kernels it is worth paying.
The third is the one that took the longest to find and is the most general: sample the block. Apply the operator to a handful of random vectors supported on the block’s columns, read the block’s rows, and what comes back spans the block’s range. The whole representation of a 512-square can be built from 256 applications of an operator that is never assembled — a measurement that belongs to a later essay, and the reason this one is careful not to claim that a decomposition is the only route.
The number to keep
Five columns, 96 × 96, eight digits. Nine hundred and sixty numbers instead of nine thousand two hundred and sixteen, for a matrix in which every entry is nonzero and no entry is small.
The rest of this field is what follows from that, and each of the next steps is a question the measurement above leaves open. What does an extra digit cost? Does the five grow when the block does? Which pairs of clusters is it true of, and who decides? And — the question that turns a storage format into a numerical method — if the whole matrix is now an approximation, what exactly has been approximated?
What links here
Computed from the collection, not written here: the essays that point at this one.
Reads more easily once this is understood
Essays that name this one as worth reading first.
- The kernel with nothing to compress
- The size the rank does not notice
- The test that costs what it saves
- A model that is a rational function
- A knob calibrated in residuals
- The cliff behind the count
- The offset that moved the slope
- A tolerance is priced by the problem
- Bracketing an error nobody can measure
- A decomposition made only of SVDs
- A factorisation that is unique for once
- A nearest point that is not there
- A rank that is a number of digits
- A rate that belongs to the matrix
- An accuracy that is a backward error
- An index that is a pair
- Built from products alone
- The fill that is not independent
- The knob that moved two things
- Which pairs are allowed to be small
- Why a Gramian can be truncated at all
- A quarter of the leaf
- The smaller cluster sets the rank
- One build tells the leaf
- An offset that bends twice
- Each halving reads its own width
Shares its objects with
Essays that name at least two of the same things, and that neither author linked.
- A rank that is a number of digits — both name kernel matrix, low-rank approximation, numerical rank, off-diagonal rank, spectral decay, truncated svd
- A nearest point that is not there — both name eckart–young, low-rank approximation, numerical rank, truncated svd
- The size the rank does not notice — both name kernel matrix, low-rank approximation, numerical rank, off-diagonal rank
- A bound that holds with probability — both name eckart–young, low-rank approximation, singular values
- A decomposition made only of SVDs — both name eckart–young, low-rank approximation, truncated svd
- A knob calibrated in residuals — both name low-rank approximation, off-diagonal rank, truncated svd
Named objects
A flat tag is an object no other essay names yet.
Eckart–YoungKernel matrixLow-rank approximationMatrix-freeNumerical rankOff-diagonal rankSingular valuesSparsitySpectral decayTruncated SVD