Skip to content

What every module shares

How long the per-site work takes

Every module here reads the same thing: one comma-separated list of read counts per pool per site, in the depth tables. On a genome running to tens of millions of called sites, splitting those strings is most of the time a module spends, and for basicstats and association it is essentially all of it. So each of the derivations that does per-site work has two implementations, and the compiled one is the default. Both produce the same numbers; the choice is only ever about time.

Measured on an Intel i7-7700HQ at 2.8 GHz with R 4.6.1 and GCC 16.2, over a corpus of 3.2 million sites — 92% biallelic, 7% triallelic, 1% tetrallelic, which is what a called cohort looks like. CPU time, scaled to what 100 million called sites would cost in one pass:

Derivation Reads Vectorized R Compiled Ratio
site_diversity one pool 135 s 11 s 12×
allele_frequencies one pool 104 s 16 s
allele_frequencies six pools 675 s 69 s 10×
nei_distance six pools 710 s 22 s 32×

The last row is a different kind of work and its ratio should not be averaged with the others. The three parsers split strings, and what limits them is memory traffic; nei_distance reads what a parser already produced and does arithmetic on it. Its vectorized form makes one pass per pair of pools — fifteen at six pools — so what compiling removes is interpreted call overhead that grows with the square of the pool count, not with the site count. Expect its ratio to climb with more pools where the parsers' will not. It is also the reason mds is the one module whose statistic costs about as much as reading the table for it.

For the three parsers, take the ratio as about ten, not as a precise figure. It is not constant, and it is not constant in a direction worth knowing about: on a corpus small enough to sit in cache the compiled path looks two to three times better than this, and below about ten thousand sites the whole call is shorter than the clock can resolve, so a figure extrapolated from a small test corpus will flatter it by a wide margin. At three million sites one pool's cells are already 48 MB — past any current L3 — and both implementations spend their time waiting on memory rather than on arithmetic. That is the regime a genome is in. A machine with more memory bandwidth than a 2017 laptop does better than this table on both columns.

What the table is really for is the decision it supports. For association, reading a hundred million sites across six pools is eleven minutes of splitting strings against one. For mds, which pays for the parse and the distance both, it is twenty-three minutes against a minute and a half. That is why the compiled path is the default rather than something to ask for, and why a module offering it stops rather than quietly falling back when it cannot build one.

dev/scripts/bench-compiled-paths.R is what produced the table, and re-running it on your own machine is how you find out what these numbers are where you work.

The analysis environment already has a compiler. Conda's r-base depends on one — GCC on Linux, clang on macOS — because R needs a toolchain to build packages from source, so an environment built by PoolSeqFlow analysis install can compile on every platform this ships to.

Each module publishes its compiled source into the results folder whether or not the run used it, and the header of the module's own script beside it names the path that produced the numbers.

When the compiled path will not build

Having a compiler and being able to use it are two things. Rcpp::sourceCpp writes a source file, builds it, and then executes the result out of a temporary directory — so a machine that forbids executing from temporary storage fails here even with the toolchain installed. The usual causes:

What you see Why
The tools required to build C++ code for R were not found the environment is not activated, so its compiler is not on PATH. Run the module through PoolSeqFlow analysis, which activates it, rather than calling Rscript yourself
cannot open shared object file, or a build that succeeds and then fails to load /tmp mounted noexec, which is common on hardened clusters. Point TMPDIR at a filesystem you may execute from
a permission error while compiling a read-only home, or a read-only conda environment shared across users

The run stops and says so; it never quietly drops to the other path. A run that silently took a different implementation is a run whose timings mean nothing, and you would have no way to know which produced your numbers. Add nocpp after the module name and it will finish — slower, and with identical output.