Visualizing the Gilbreath expectation sequence

One byproduct of learning how to use coding agents to create visualization apps is that it now becomes straightforward to convert any figure in one’s papers that had already been generated by code (e.g., in Python) into a more interactive, animated applet.

I can illustrate this with Figure 1 from my recent paper on the Gilbreath conjecture with Chase and Hunter, reproduced below:

This plot displays both exact and numerically simulated values of a certain poorly understood sequence relating to the Gilbreath conjecture, which I will call the “Gilbreath expectation sequence” here for lack of a better name. The definition of the sequence is as follows. Consider a “Gilbreath array” which is an inverted pyramid, where the top entries are independent exponential random variables of mean 1, and all the other entries are the absolute values of the differences of the two entries immediately above it. Thanks to the visualizer app, I can quickly give an example (with ):

The left diagonal entries are then random variables; the sequence are defined to be the expectation of these values. (The process is stationary, so in fact any entry on the row will have expectation .)

If one starts with the first normalized prime gaps (which have expectation about , and are conjecturally distributed asymptotically according to a geometric distribution), then standard conjectures (e.g., the prime tuples conjecture) predict that the row entries should decay like , at least for small . So the Gilbreath conjecture appears to be tied to how fast the sequence decays with .

One can in principle work out each value of as an explicit rational number by performing a certain complicated multivariate integral, but in the paper we only did this for (the orange line in the above figure); for the remaining we performed a Monte Carlo simulation with Gilbreath arrays to obtain a numerical approximation (in blue), which (as per the law of large numbers) agreed well with the theoretical values. A later calculation of Michael Ross extended the theoretical values to , maintaining the good fit:

The asymptotic behavior of the sequence remains mysterious. Clearly, it is not monotonic; in fact we cannot even prove it is bounded. The best we could do in our paper was establish an inequality which, roughly speaking, showed that cannot decay faster than .

In a recent preprint of Ross, these numerics were extended, and a rough empirical prediction

was proposed for some constants 0" class="latex" /> and 1" class="latex" /> (empirically ), where is the number of 1’s in the binary expansion of ; in particular, it is the fluctuation in this quantity that is intended to explain much of the non-monotonic behavior of . These are now all displayed in the following companion applet, which was a routine matter to generate in about an hour by the coding agent (which by this point has extensive experience with creating such apps, encoded via a “skill” markdown file that it maintains):

The appearance of the quantity may initially appear mysterious, but it is related to Lucas’s theorem, Kummer’s theorem, and the Sierpinski gasket. Consider for instance a Gilbreath array where all the entries are zero except for a single “spike”. Then the following Sierpinski pattern emerges:

Here is what an version of this picture looks like (with the spike positioned at the 32th entry):

The number of 1s in the row is then (if we index the rows starting from zero), which is at least of the same shape as the empirical prediction, albeit with different constants.

Numerically, we seem to observe fragments of Sierpinski gaskets being generated before decaying (often due to “collisions” with other gaskets):

However, it is not clear to me at all what the asymptotic probabilistic model should be, even heuristically; it does not resemble any random shape model that I am familiar with. But perhaps there are readers more expert in probability theory or statistical physics who may be able to suggest such an asymptotic limit?