Count the bytes, not the FLOPs

TL;DR: tsbootstrap’s fast path used to build every bootstrap replicate into one giant batched tensor: on a routine job, a 160 MB tensor pushed out to main memory and hauled back to produce one summary number per replicate. We rebuilt the hot path around two rules: never materialize an array whose only consumer is a reduction, and never carry state you can derive. Every head-to-head benchmark cell now favors the fused path; the settled benchmark grid sustains 3.1x to 20x on the longer series. The last stubborn cell, where the runtime was mostly Python building random-number seeds, dropped from 13.1 ms to 0.55 ms in the redesign’s validation run once we deleted the seed objects. The full replicate tensor is still available when you want it, and on that materializing path our default backend remains slower than the incumbent. That row is printed below with the rest. Early this summer, in the middle of a performance push I thought was going well, we ran tsbootstrap’s block-bootstrap engine head-to-head against arch, a mature econometrics library whose bootstrap loop has had years of polish, for the first time. We lost. At a realistic workload, a few thousand replicates of a few-thousand-point series, the whole call came back several times slower end-to-end, on the machine we developed on, before the compiled engine this article describes existed. Every release we had ever shipped carried that loss; users hitting exactly that workload had been paying it the whole time, and nothing in our tooling would ever have told us.

I had not seen it coming: our benchmark suite had only ever compared the library against its own history. A regression suite tells you when you get slower than yesterday; it cannot tell you that you ship slower than the other library, at the exact workload where users would notice. The day a head-to-head existed, the loss appeared. It had been there all along. The profiler’s first answer was the ordinary one. Most of the measured deficit was removable interpreter tax: a per-replicate Python index loop, thousands of seed objects spawned per call, stacked object wrapping. But beneath it sat a second, structural problem, the reason the reflexive fix could not win either. The reflexive fix for Python overhead is to batch: replace the loop with a single vectorized NumPy operation.

Batching removes the interpreter from the inner loop and swaps it for B-fold traffic through main memory (B being the number of bootstrap replicates, the count every cost in this story scales with) because the batched design materializes every replicate at once. The incumbent’s replicate loop is ordinary Python too. It wins by keeping its working set (the data its loop actively touches) small enough to live in cache, consuming each resample the moment it exists. We were not going to beat that by batching harder. Here is the machine’s side of the argument. A CPU core computes only on what it holds in its registers. Think of that as your hands. Behind them sit the caches: L1, a desk you can reach without looking up; L2, the shelf behind you; L3, shared with the other cores, a cabinet down the hall. Behind all of that is DRAM: main memory, the warehouse across the street, bigger and slower at every step, orders of magnitude in latency top to bottom. The hardware hides none of that from your runtime; it hides it only from your source code. Now put a bootstrap on that machine. A bootstrap, one time for anyone who has not run one: it re-draws your data thousands of times to see how much a statistic would wobble across plausible alternate samples, and each redrawn copy is a replicate. Take a modest, entirely realistic job: a series of n = 10,000 observations, B = 2,000 bootstrap replicates, float64 (8 bytes per number). Build all the replicates at once, the way a batched NumPy design wants to, and you have asked for a B x n tensor: 10,000 x 2,000 x 8 bytes = 160 MB. That is far beyond the tens of megabytes of L3 a core complex sees. Every byte of that tensor is written out to main memory, then hauled back in so the statistic can read it once, and if the statistic is a mean, the answer you keep is 2,000 numbers. Sixteen kilobytes, bought with 160 MB of round trips.

Every replicate is a copy of your data. The tensor is B times larger than the data. It lives in DRAM. Gather and reduce in the same pass.

Keep only the answer. 20 MB where 1.94 GB stood, B = 50,000. measured 2026-07-04 on v0.4.0; the streaming path is unchanged since The batched design's tensor: assembled, flooding past the cache boundary, and collapsed into the answer. The boundary is schematic; the byte counts are arithmetic on the worked example. At a larger job (fifty thousand replicates of a shorter series) the materializing path peaked at 1.94 GB where the streaming path held 20 MB, roughly 97x (measured 2026-07-04 on version 0.4.0; the streaming path is unchanged since). data table quantityvaluereplicate tensor, worked example (n = 10,000, B = 2,000, float64)160 MBanswer kept (one mean per replicate)16 kBpeak memory at B = 50,000 on an n = 2,000 series: materialize vs streaming1.94 GB vs 20 MB (~97x), measured 2026-07-04 on v0.4.0peak resident memory at the worked example’s shape, stationary-block resample (measured 2026-07-11): materialize vs streaming617 MB vs 234 MB, where 234 MB is the bare process baseline; the streaming run is indistinguishable from a run that computes nothing. The materialize delta exceeds the 160 MB tensor because the path also holds a returned copy and the index matrix The tsbootstrap that ships today never builds that tensor on its fast path. The story of how it stopped: The first wall is memory traffic: the bytes themselves, moving through the hierarchy. The second took longer to see, because it hides inside the language runtime and the library’s own design, not the hardware: the cost of state (allocated objects and seeded generators) that exists only to make randomness happen in the right order. We stopped materializing arrays, and we stopped carrying state. Wall one: memory traffic¶ Start with what the batched design actually does.

To bootstrap a series with blocks (contiguous runs resampled whole, so the short-range dependence that makes a series a series survives inside each block; the first piece in this series is about what happens when it doesn’t), you draw random block starts, expand them into resample indices, gather the data through those indices, and compute your statistic on each replicate. The batched version does each of those as one array operation over all B replicates: build a (B, n) index matrix, gather into a (B, n, d) values tensor (d is just the number of columns in your series), then reduce along the middle axis. Count the work in each cell of that tensor. Eight bytes written on the way out. Eight bytes read back for the reduction. Roughly one floating-point addition, one addition per sixteen bytes moved, where a compute-bound kernel like a matrix multiply does hundreds of operations for every byte it pulls. That ratio, the arithmetic intensity (how much computing you do per byte you move), is the number that decides whether your workload is limited by the processor or by the memory system. Ours is about as low as arithmetic intensity gets. The cores were never the constraint; they spent the benchmark waiting on memory. The gather is slow not because the code is interpreted but because its arithmetic intensity is far too low to keep the processor busy: it is a bytes problem being graded as a FLOPs problem. The incumbent’s loop, seen this way, stops looking old-fashioned and starts looking correct. One replicate at a time means one (n, d) resample at a time: a working set that fits in cache, is consumed by the statistic immediately, and is gone before the next one arrives. That working set never grows with B, so DRAM barely enters its story. The batched design’s traffic grows linearly with B, and once the tensor outgrows the last cache level, every additional replicate’s bytes take the round trip through main memory. NumPy’s programming model quietly couples two decisions that are really separate: batching the dispatch and batching the bytes. The loop was winning because it never built the thing we were building. Here is the whole argument as accounting, on the worked example: n = 10,000, B = 2,000, float64, one mean per replicate.

Redo it on a napkin for your own workload; the ratio is the point.

quantitymaterialize pathfused, streaming pathintermediate allocated160 MB, the full replicate tensorone (n,) index row + one (n, d) scratch per thread (~80 kB each)bytes through DRAMat least twice the tensor: written out, read backconsumed in cache, per work itemuseful output kept16 kB, one mean per replicate16 kB, the same answerbytes moved per byte of answer~20,000 : 1~1 : 1

None of this is a new insight: compilers have fused loops to avoid intermediate arrays for decades, and every numerical programmer eventually meets the memory wall. What is worth knowing is that the same wall now dominates transformer training and serving. FlashAttention, the kernel that reshaped both, describes its own core move in the caption of its first figure: it uses “tiling to prevent materialization of the large N × N attention matrix”, and that one refusal, plus fusing the surrounding operations into a single kernel, yielded a 7.6x speedup on GPT-2’s attention. The hierarchy it plays against is the GPU’s version of ours: an A100’s main memory streams at up to 2.0 TB/s, while the tiny on-chip SRAM (192 KB per streaming multiprocessor) runs roughly ten times faster. Same shape, same conclusion: keep the working set in fast memory; never write the big intermediate at all. FlashAttention is the proof that the old compiler principle still binds at datacenter scale. One pass, and nothing left behind¶ The fix is a single compiled kernel that does everything the batched pipeline did, per replicate and in one pass, without ever building the B-wide tensor. Per-worker scratch still exists, but it is n-sized and transient, not B times the data. The standard name for this is kernel fusion: instead of separate passes that hand each other arrays, one loop body builds, gathers, and reduces for one replicate and emits only the statistic. Each replicate is one unit of parallel work.