Validation against R

Reference. The EcoTyper R pipeline (R 4.4.3, NMF 0.28, HiClimR 2.2.1; docker/Dockerfile.r) on x86-64/glibc without FMA contraction. The port ran on the same machine with ECOTYPER_NO_FMA=1, with the same seed and the same inputs. CIBERSORTx outputs are taken as inputs, not recomputed.

Method. Every output file of every tutorial is compared with R’s (validation/tutorials.py, ecotyper.validate):

  • byte-identical: the files are identical byte for byte;

  • discrete-equal: all non-floating-point cells (labels, ranks, assignments, gene and state sets, row and column names) are identical and the shapes match;

  • max rel. diff: the largest relative difference over all numeric cells, compared as parsed doubles.

Result. 2,608 files compared: all discrete-equal, none missing, numeric differences ≤ 8e-11 relative (target 1e-8). The remaining differences are last-digit noise from the order of floating-point operations in a few reductions (e.g. BLAS vs R’s loops). None of them changes a state, ecotype, rank or gene set. The per-tutorial breakdown and a chart are on the Benchmarks page.

What was needed for exactness on x86

  • R on x86-64 accumulates sum, mean, colMeans, rowMeans, cov/cor/var/scale in 80-bit long double. src/ecotyper/rcompat/_ld.c reimplements these kernels (transcribed from R 4.4.3 summary.c, array.c, cov.c); used automatically when long double is wider than double (ECOTYPER_LD=0/1 overrides). validation/ld_check.sh: 12/12 kernels bit-identical to R.

  • R’s write.table number formatting (scientific() in format.c) uses long double too; ld_scientific reproduces it. validation/fmt_check.sh: 462,063 doubles, 0 differences. This matters because R re-reads its own files between steps, so a last-digit difference in text propagates into ties and ranks.

  • Brunet NMF, RNG (Mersenne Twister + R’s sampling), hclust, etc. are bit-exact (unit tests).

Not exact / caveats

  • On arm64 (Apple Silicon) long double == double, so the numba paths are used; results are discrete-equal to the x86 R reference but not byte-identical, because libm log2 differs from glibc by 1 ulp for a few inputs and FMA contraction differs. The reference platform for exactness is x86-64/glibc.

  • Mac-R vs x86-R outputs differ from each other for the same reason (compiler/libm); the x86 R run is the primary reference.

  • ECOTYPER_NO_FMA=1 is needed on CPUs where numba would otherwise emit FMA, to match R’s non-FMA conda build.

  • CIBERSORTx is not reimplemented and was not tested.

  • Exact-mode speed-up on t5 is small (1.18x): most of its wall time is serial parent steps (per-cell-type cross-correlation, DE) that are kept in R’s exact reduction order; they were not parallelised because that would risk bit differences. t5 peak RSS summed over processes (380 GB) double-counts shared pages (PSS 15 GB).