What a jump costs evenly spread points

An error bar for points that are not random ended with a measurement and a warning. The measurement: Sobol’ points with their binary digits scrambled at random integrate a smooth function with an error that falls like $N^{-3/2}$ in the number of points, against $N^{-1/2}$ for independent random points — at four thousand points, a thousand times more accurate. The warning: that rate belongs to the smooth integrand it was measured on. For an integrand with a jump inside the square, scrambled points still improve on random ones, but by a smaller power.

This essay measures the smaller power, and finds where it comes from. The answer is short enough to draw. Evenly spread points earn their accuracy by putting exactly one point in every small box of the square. A box in which the integrand is smooth is then integrated very well. A box cut in two by a jump is not: its one point lands on one side or the other, and the box’s contribution is right or wrong by a whole box’s worth of value, at random. Everything else follows from counting how many boxes a jump cuts.

A box the edge crosses #

The integrand in the figure is $e^{x+y}$ below the line $y = 0.15 + 0.7x$ and nought above it. That is the simplest kind of jump: a smooth function, cut off sharply along a straight edge that is tilted against the sides of the square. Such integrands are common. The probability that a random quantity lies below a threshold is the integral of a function that is one on one side of a surface and nought on the other; Buffon’s needle estimates $\pi$ from exactly such a function, one when the needle crosses a line and nought when it does not.

The scrambled Sobol’ points in the figure have the property the previous essay established. Divide the square into $N$ equal boxes in any of the binary ways — four by four, two by eight, sixteen by one — and every box holds exactly one point. With $N = 256$ the square-shaped boxes are sixteen by sixteen, and each holds one point, placed at random inside it by the scrambling.

For a box the line misses, one point at a random place in a smooth piece of the function gives an estimate of that box’s contribution whose error is very small, and the errors of different boxes partly cancel. That is the source of the $N^{-3/2}$ rate. For a box the line crosses, the single point lands above the line with a probability equal to the share of the box above it, and then contributes nothing; or below, and contributes as though the whole box were below. The estimate of that box is right on average and wrong, in any one run, by an amount comparable to the box’s whole value. Those boxes behave like a sample of independent coin tosses, one for each crossed box, and the coin tosses are not averaged away by anything.

A line crosses a number of boxes proportional to the number of boxes along one side: $6$, $12$, $26$ in the three panels, between one and two times $\sqrt N$. Each contributes an error of about one box’s area, $1/N$, independently, so the variance is about $\sqrt N \times (1/N)^2 = N^{-3/2}$, and the error — the square root of the variance — is about $N^{-3/4}$.

Half the exponent, and no worse #

The prediction can be checked directly, and it can be compared with what a jump along an axis does.

The smooth integrand improves with slope $-1.47$, close to the $-\tfrac32$ of Owen’s theorem for scrambled nets. The tilted jump has slope $-0.76$, and the circle, a curved jump, $-0.75$: half the smooth exponent, as the box count predicts. Random points have slope $-0.49$, the usual $-\tfrac12$. So the jump costs evenly spread points half their exponent and leaves them still well ahead: at sixteen thousand points the scrambled estimate of the tilted integrand is more than ten times as accurate as the random one.

The jump along the vertical line $x = 0.37$ behaves differently, with slope $-1.00$. The box picture explains that too. The binary boxes of a Sobol’ set include vertical strips of width $1/N$, each of which holds exactly one point. A vertical jump lies inside exactly one of those strips. Every other strip lies wholly on one side of it, so only one strip is left to chance, and its contribution is wrong by about $1/N$ — an error of order $1/N$, not $N^{-3/4}$. A jump aligned with a coordinate is almost free, because the point set was built to be even in exactly those directions.

For the simplest evenly spread design the rate is a theorem. Put one independent random point in each of $N$ equal squares — jittered sampling — and the box count is exact: the variance of the estimate of an area bounded by a smooth curve is of order $N^{-3/2}$, an error of $N^{-3/4}$. Scrambled Sobol’ points are harder to analyse, because their boxes come in many shapes at once and the points in different boxes are not independent. Zhijian He and Xiaoqun Wang proved in 2015 that for a smooth function cut off by a jump the mean-square error of a scrambled net is at most of order $N^{-4/3}$ in the plane, up to an arbitrarily small loss in the exponent — an error of $N^{-2/3}$. The measured slope is steeper than their bound and matches the jittered rate, which suggests that the box count, not the bound, is the truth here; a proof of that for scrambled nets is not available.

The bound that has nothing to say #

The classical theory of evenly spread points cannot reach this. Points too even to be random explained the Koksma–Hlawka inequality: the error of any deterministic point set is at most its discrepancy, a measure of unevenness, times the variation of the integrand in the sense of Hardy and Krause. For a smooth function that variation is finite and the bound gives a rate of about $1/N$.

For a jump along a tilted line the Hardy–Krause variation is infinite. The variation is built from rectangles aligned with the axes: it adds up how much the function changes across the corners of every such rectangle, at every scale. A tilted edge passes through the corners of arbitrarily small rectangles in a way that makes each one contribute about the same amount, and there are infinitely many, so the total diverges. A jump along an axis does not, because it never cuts a rectangle diagonally. The bound therefore says nothing at all about the integrand in the figure — not a weak rate, but no rate — while the measured error is falling steadily like $N^{-3/4}$.

This is the same distinction that which curves have a length at all draws for curves: some perfectly reasonable objects have infinite variation, and a theorem that needs finite variation is silent about them. Here the silence is an artefact of the bound, not of the method. Scrambling turns the question from one about the worst case into one about variance, and the variance can be counted box by box even when the variation is infinite.

Only the axes are safe #

The two extremes in the figure are a jump at an angle of about $35°$ and one at $90°$. Turning the edge from one to the other shows how sharply the saving depends on direction.

Measured against random points with the same budget, the scrambled points make about a fifth of the error at 256 points, a seventh at 1,024 and a tenth at 4,096 for every edge well away from the axes. That share shrinks only like $N^{-1/4}$, the ratio of $N^{-3/4}$ to $N^{-1/2}$: it takes sixteen times the points to halve it. At the two ends of the range, where the edge is horizontal or vertical, the share falls to a few per cent at 256 points and to about one per cent at 4,096, and it keeps falling like $N^{-1/2}$.

The striking feature is how narrow the safe directions are. A tilt of three degrees already gives up most of what the axis saved: the share there is 11, 8 and 5 per cent at the three sizes, far nearer the tilted values than the axis ones. A line at three degrees from horizontal crosses only a few of the horizontal strips of width $1/N$ — but it crosses a great many of the square boxes of side $1/\sqrt N$, and once it crosses more square boxes than strips, it behaves like any other tilted line. The exact angle at which the cost appears moves towards the axis as $N$ grows; in the limit, any tilt at all gives $N^{-3/4}$.

Lattice rules have no safe direction at all. They work by being blind to a thin set of frequencies that a smooth periodic function hardly uses, and a jump in any direction spreads its weight over frequencies that decay only like one over their size, many of them in the lattice’s blind set. Sobol’ points at least have their axes, the directions in which they were built to be even.

This matters because the coordinates in which a problem is posed are often arbitrary. An integrand whose jump is aligned with a coordinate in one formulation can be tilted in another, and the scrambled points cannot tell which formulation they were given, only how the edge lies against their boxes.

The length of the edge, not the size of the region #

The box count says that the error depends on how many boxes the jump crosses, and that is set by the jump’s length. Random points see something else. Their error in estimating the area of a region is the binomial error, the square root of $p(1 - p)/N$ where $p$ is the area: it depends on how big the region is and not on its shape. The two can be separated by estimating the area of circles of different sizes.

The scrambled error grows like $\sqrt r$, with a measured slope of $0.49$: like the square root of the circumference, because the number of crossed boxes is proportional to the circumference. Divided by $\sqrt r$, the error is the same to within seven per cent from the smallest circle to the largest. Random points’ error grows like $r$ itself for small circles, slope $0.97$, because for a small region $p(1 - p)$ is nearly $p$, proportional to the area $\pi r^2$, and its square root to $r$.

So the gain from even spreading depends on the shape. On the smallest circle, of radius 0.04, the scrambled points beat random ones by a factor of 3.6; at radius 0.35, by 7.8. A region with a lot of edge for its area — small, thin or ragged — is the case in which evenly spread points gain least. A long thin strip at an angle across the square is close to the worst case: little area, much edge.

The same arithmetic governs a much older question. The dots a circle catches counts the points of the square grid inside a circle of radius $R$ and compares the count with the area $\pi R^2$. The difference comes entirely from the grid squares the circle cuts — about $R$ of them — and if each were an independent coin toss the difference would be about $\sqrt R$. Gauss’s circle problem asks for the true size of that difference, which is conjectured to be close to $\sqrt R$ and proved only to be well below $R^{2/3}$. A fixed grid is not a scrambled one, and its errors along the edge do not behave like independent tosses, which is exactly why the problem is hard. Scrambling makes the coin-toss picture true by construction, and so makes the answer easy.

More dimensions, less advantage #

In $d$ dimensions, $N$ evenly spread points fill $N$ boxes of side $N^{-1/d}$, and a jump is a surface of dimension $d - 1$. The number of boxes a surface cuts is about its area divided by the area of a box face, $N^{(d-1)/d}$, which is a share $N^{-1/d}$ of all the boxes.

The exact counts follow the prediction: the share of boxes cut falls like $N^{-1/2}$ in the plane, $N^{-1/3}$ in three dimensions, $N^{-1/4}$ in four and close to $N^{-1/8}$ in eight. Each cut box contributes a variance of about $1/N^2$, so the error is about the square root of $N^{(d-1)/d}/N^2$, which is $N^{-1/2 - 1/(2d)}$ — the jittered-sampling rate, and what the measurements in the plane suggest for scrambled points too. The exponents are $-0.75$, $-0.667$, $-0.625$ and $-0.5625$.

The extra exponent beyond random sampling’s $-\tfrac12$ is $1/(2d)$, and it shrinks as the dimension grows. In eight dimensions a jump leaves scrambled points with an exponent only $\tfrac1{16}$ better than random points: sixteen times the points to halve the remaining advantage. In fifty dimensions, the size of many problems in finance, the rate is indistinguishable from random sampling’s.

That is not the end of the matter, as the error that does not care how many dimensions and its successors found for smooth integrands: what matters in high dimension is often how many coordinates the integrand really depends on, and how the jump lies against those. A jump that depends only on the first two coordinates of a fifty-dimensional problem behaves like a two-dimensional jump. But the box count says plainly what a jump genuinely spread across many coordinates costs.

Integrating across the jump first #

The damage comes from one thing: within a cut box, a single point decides which side the box is on. If the side could be settled exactly instead of sampled, the damage would go. For the tilted integrand it can. For each $x$, the integral of $e^{x+y}$ from $y = 0$ up to the jump at $y = 0.15 + 0.7x$ can be done by hand:

$$\int_0^{0.15 + 0.7x} e^{x+y}\, dy = e^x \left(e^{0.15 + 0.7x} - 1\right).$$

What is left is a smooth function of $x$ alone, with the same integral as the original, and no jump anywhere.

The scrambled points on the smoothed integrand have slope $-1.52$, back to the smooth rate, and at sixteen thousand points their error is about two thousand times smaller than on the original. Nothing was approximated. The jump was integrated over exactly, where it lies, and only the smooth remainder was sampled.

The method is old in Monte Carlo, where it is called conditioning: replace a random quantity by its average over some of the randomness, done exactly, and the variance can only go down. For random points the gain is a constant factor. For evenly spread points it is a change of rate, because it removes the one feature the points cannot handle. Recent work in high dimensions, under the name preintegration, does this for integrands from finance and from partial differential equations with random coefficients: choose one coordinate along which the jump can be located, find the jump along that coordinate at each sampled point by solving an equation, and integrate across it exactly. Griewank, Kuo, Leövey and Sloan gave conditions in 2018 under which the result is smooth enough for the full rates to return.

The price is that the location of the jump must be found, at every point, along one direction in which it crosses only once. When the jump is the boundary of a complicated region that may not be possible. It is the same trade as in sampling where the answer lives: knowledge of where the integrand changes buys accuracy, and the method is only as good as that knowledge.

Where jumps come from #

Jumps are not exotic. Any probability is the integral of an indicator, which is all jump. In computational finance, an option that pays a fixed amount if a price ends above a level — a digital option — has a discontinuous payoff, and so does a barrier option that is cancelled if the price ever crosses a level; their prices are integrals over many dimensions of functions with jumps, and the slow convergence of quasi-Monte Carlo on them was noticed long before the rate was proved. One response in that literature is to rotate the coordinates so that the jump lies as nearly as possible along one of them, which turns a tilted jump towards the cheap axis-parallel case the angle figure shows.

In computer graphics, the colour of a pixel is the average of the scene over the pixel’s area, and the edge of an object crossing the pixel is a jump. Renderers sample each pixel at several points, and stratified or quasi-random placement is standard. Don Mitchell observed in 1996 that stratified sampling of a pixel with an edge has error falling like $N^{-3/4}$, against $N^{-1}$ or better in smooth regions — the same count of cut cells, met by people trying to draw a clean line. Smooth regions of an image converge quickly and edges are where the noise remains.

Still open: the rate for rougher jumps and what to do in high dimension #

The rates above need the jump’s boundary to be reasonably regular — a smooth surface, or one made of finitely many smooth pieces. For boundaries that are fractal or have infinitely many wiggles, the number of boxes cut grows faster than $N^{(d-1)/d}$, and the variance of scrambled nets depends on the boundary’s box-counting dimension in a way that has been worked out only in special cases. How much of the advantage survives for the boundary of a set like the Mandelbrot set, or the level set of a rough random function, is not generally known, and for scrambled nets even the smooth case rests on the gap between He and Wang’s proved $N^{-2/3}$ and the observed $N^{-3/4}$.

In high dimension the open questions are practical. Preintegration works when one coordinate can be found along which the jump is crossed once and can be located cheaply; when no such coordinate exists, or when there are many jumps in many directions, there is no general method that restores the smooth rate, and it is not known whether one can exist. For deterministic point sets the situation is worse. Since the Hardy–Krause variation of a tilted jump is infinite, the classical theory guarantees no rate at all, and the natural replacement is discrepancy measured against half-planes or discs instead of boxes. Its best possible order is known — József Beck and Ralph Alexander proved lower bounds between 1987 and 1990 and Jiří Matoušek a matching upper bound in 1995, giving $N^{-1/2-1/(2d)}$, the same exponent as the box count — but the point sets that achieve it come from combinatorial constructions rather than from any of the sequences used in practice, and how close Sobol’ or lattice points come to it, deterministically, is not known.