Reproducibility

Reproducibility is the ability to take the same raw data, run the same analysis, and get the same answer. Nothing about the biology changes. The only thing being tested is whether the computation and the measurement are stable enough to be trusted.

Four words get used loosely and mean different things:

• Reproducibility: same data, same code, same result. A computational property.
• Replicability: new data, same procedure, same conclusion. A scientific property.
• Repeatability: same sample, same instrument, same operator, same day. A measurement property, usually quoted as a coefficient of variation.
• Accuracy: how close the result is to truth. Independent of all three above.

Accurate and reproducible are orthogonal. A pipeline that assigns the wrong reference allele at a position will do so identically every run: perfectly reproducible, wrong. A perfectly accurate assay with 25% intra-assay CV will not reproduce its own number next Tuesday. You need both, and you diagnose them differently.

How it works

Computational reproducibility fails for a small, enumerable set of reasons, and each has a fix.

Software versions. bwa-mem2 2.2.1 and bwa 0.7.17 do not produce identical BAMs. GATK HaplotypeCaller output changes across minor versions. The fix is containers plus a package manager with a lockfile: Docker or Singularity images built from pinned Conda/Bioconda environments, which is the pattern the life sciences have converged on¹. Workflow managers (Nextflow, Snakemake, WDL) bind each process to a specific container digest, not a tag, because :latest is not a version.

Reference data. GRCh37 versus GRCh38 versus T2T-CHM13 changes coordinates. Within GRCh38, the analysis set with alt contigs behaves differently from the no-alt one. Annotation is worse: VEP output depends on the Ensembl cache release and on which transcript set you flag as canonical, so --everything with cache 110 and cache 112 will give you different consequence calls for the same variant². Record the cache release in your notes the way you record the aligner version.

Nondeterminism. Multithreaded aligners can reorder equal-scoring alignments. Anything sampling (bootstrap in salmon/kallisto, permutation in gene set enrichment, UMAP) needs a fixed seed. Enrichment results in particular vary substantially by algorithm choice on the same input gene list, so the method name and its parameters are part of the result, not a footnote³.

Provenance capture. The strongest version records the full DAG that produced a file inside the artifact itself, so anyone can replay how it was made rather than trusting a README. QIIME 2’s Provenance Replay does exactly this, reconstructing an executable script from the output archive⁴. Most pipelines do not, so you have to do it yourself by keeping the workflow run directory.

In your own data

Concrete things to check when you hold a Personal Molecular Profile.

Your VCF header is the receipt. Read the ##reference=, ##source=, and ##GATKCommandLine lines. If you have two VCFs from two vendors, diff those headers first. Confirm contig names match the style you expect (chr1 versus 1), and confirm ##contig=<ID=chr1,length=248956422> for GRCh38. Mixing conventions silently drops variants in downstream bedtools intersect.

Re-annotate rather than trust a shipped annotation. Run VEP yourself with an explicit cache version and --assembly GRCh38, and keep the command line in the same directory as the output. When a variant’s ClinVar classification changes in six months, you want to know whether the change came from ClinVar or from your own tooling².

RNA-seq counts. Store the transcriptome FASTA name and version (GENCODE v44 is not v39) alongside the count matrix. Gene symbols are unstable: symbols get renamed between releases, which breaks joins across time. Key on Ensembl gene IDs and use HGNC-approved symbols only for display, which is the point of standardized nomenclature⁵. Strip version suffixes (ENSG00000141510.17 to ENSG00000141510) consistently or not at all, but decide once.

Assay reproducibility. For blood biomarkers and proteomics, the number you care about is CV, not the assay’s accuracy claim. Ask your vendor for intra- and inter-assay CV per analyte. For Olink NPX or SomaScan RFU values, a 15% inter-assay CV means a 10% change between two draws is noise. Continuous glucose data has its own version: sensor-to-sensor bias can shift your mean by 5-10 mg/dL, so compare within a sensor session before comparing across them.

Common mistakes: running an analysis in a Jupyter kernel with hand-installed packages and no environment.yml, keeping intermediate files without the command that made them, and re-running “the same” pipeline after a conda update. Use a cloud-optimized, versioned workflow set if you do not want to maintain your own⁶, or containerize per-step as CoBRA does for ChIP/ATAC⁷.

Limitations

Reproducibility does not get you truth. A result that reproduces exactly can still be an artifact of one reference, one annotation release, or one batch. Replication on independent data is the test that catches that, and at biobank scale, systematically checking which genotype-phenotype associations replicate is how spurious ones get filtered out⁸. Consortium-scale work has repeatedly found that apparent signal depends on assay and analysis choices, which is why cross-platform agreement matters more than internal consistency⁹.

It also does not tell you what a finding means for you. Variant interpretation and biomarker changes need a clinician and, for anything germline and actionable, a genetic counselor. A reproducible pipeline gets the number right in the sense that it is stable. Someone qualified still has to say what it implies.

Woolf Software builds longitudinal molecular profiles of individuals: whole-genome sequencing, RNA sequencing, proteomics, blood biomarkers, and continuous glucose data, integrated into one model of you. Build your profile.