migec 2.5.1

The composition table says which half of the barcode a position is in

checkout.umi_composition.tsv gains a segment column (cell / umi). The table has spanned the whole barcode since 2.5.0 keyed the counters on cell + UMI; without this it is one anonymous run of 26 positions. The halves are read differently and that is the point: a cell position is drawn from a whitelist so its skew is the whitelist's, a UMI position is a free synthesiser mix so its skew is a defect you are paying for.

Three figures were lying, and one run was reported as broken

set format in a panel script needs a single brace and a single percent. Panel scripts are substituted with str.replace, not str.format — only _PREAMBLE is an f-string — so "10^{{%%T}}" reached gnuplot verbatim and every y tick of the committed assets/consensus_error.svg read a literal 10%T. That axis spans 10^0 to 10^-250 and no tick said so. Three panels affected: consensus_error, quality_calibration, coverage.

An empty table is skipped, not failed. A purely positional pattern has no constant bases, so checkout.quality_calibration.tsv is written with a header and no rows — the calibration is measured against the pattern's own constant bases and 10x has none. migec plot handed that to gnuplot, got x range is invalid, and reported a failure. Empty and missing are now reported apart, because the advice differs: a missing table means run the stage, an empty one means the stage ran and had nothing to put in it.

coverage is now log2 on x, since checkout.coverage.tsv is written in MIGEC's own doubling bins and a base-10 axis puts the ticks somewhere other than the data. umi_pwm autoscales y instead of pinning [0:1]: a well-made barcode sits at 1/4 everywhere, so a fixed unit axis spent three quarters of the panel on emptiness.

Two README panels, on a real library

assets/coverage.svg and assets/umi_pwm.PBMC.svg, drawn off checkout's own tables on sc5p_v2_hs_PBMC_1k VDJ-T lane 1: 3,155,166 read pairs, 100% assigned, 311,421 molecules, effective length 25.85 of 26 nt. Experimental, not simulated, and SOURCES.md says so next to the simulated panels it sits beside.

The PWM is the case for the segment column: the 16 cell positions swing between 0.198 and 0.305 while the 10 UMI positions sit tight around 1/4. The MIG size panel shows 195,251 molecules seen once against 1.80 M reads piled onto 5,118 molecules in the 256 bin.

Also in this line, from 2.5.0's follow-up audit

• A barcode too long for the key is refused. 20 nt of cell and 16 of UMI each clear their own 32-base bound; their 36 nt sum does not fit the key, and two molecules sharing their first 32 bases merged silently — measured, 1 distinct barcode where there were 2. Refused at pattern-compile time, and across both mates for dual-end.
• >> 64 is undefined, not zero. x86 masks the count to 6 bits, so it became >> 0 and the UMI landed back on top of the cell barcode. Confirmed under UBSan that the old expression is flagged and returns a different value.
• place_reads indexes each read once, not once per pair: 40,786 -> 150,152 reads/s, output MD5-identical.

Full detail in CHANGELOG.md.

🤖 Generated with Claude Code

migec 2.5.0

A molecule is sample + cell + UMI, and checkout's counters were keyed on the UMI alone

w.umis.emplace_back(s, umi_key) merged every cell's copy of one UMI into a single entry. On sc5p_v2_hs_PBMC_1k that read 221,026 pooled UMIs where there are 311,962 molecules, put the depth 1.41x high, and fired the saturation warning at an apparent 21% occupancy against a true 3e-6. The key is now cell then UMI.

The follow-through: effective_length is a property of the counter, so it spans 26 nt on a 10x run. Weighed against the 10 nt UMI, the composition-skew warning read "25.4 of 10 nt usable" — which cannot be true, so it was dead on every single-cell library. It is weighed against barcode_length now. checkout.summary.tsv carries all three lengths as columns, and checkout.barcode_space.tsv's second column is renamed umi_length -> barcode_length, which is what its value always was.

Cell Ranger, on all three axes

scripts/compare_cellranger{,_chains,_contigs}.py, docs/single_cell.rst. Run against both 5.0.0 (10x's published calls) and 10.1.0 (run here) — the two agree at Jaccard 0.9938, which is the control saying the migec gap is not version drift.

• Cells: migec 888 against 479, sharing 470, and it loses 1.9 points of reads-in-cells. The extra barcodes are nearly empty, so a cell count on its own is not an accuracy figure.
• Per-cell chains, migec + arda against a per-cell assembler: TRA 426/426, TRB 468/469 recall, junction agreement 0.9507 / 0.9915, 22 s over 47,584 consensuses.
• Reference-free contigs (arda cells): 933/943 CDR3s recovered verbatim, k-mer coverage 0.9759, chain recall 0.9777.

The knee is Kneedle at its global maximum, and it refuses when there is no knee

Kneedle's published rule — walk the local maxima, stop at the first that falls by a sensitivity step — is degenerate without the paper's smoothing spline: on 136,032 barcodes the unsmoothed difference curve has 287 local maxima and the walk stops at fifteen cells. The global maximum needs no smoothing parameter, and is now guarded at 10x the mean molecules per observed barcode, so an ambient-only library reports no knee instead of a meaningless rank.

place_reads asks union-find before it scans: 640.9 s -> 32.6 s

Exact, not a heuristic: join already returns immediately when the two roots match, so the pairs skipped are the ones whose answer it discarded. Consensus FASTQ and mig.tsv are MD5-identical over 47,584 molecules.

Full detail in CHANGELOG.md.

🤖 Generated with Claude Code

migec 2.2.1 — the ctDNA ground truth was found, not built

The wheel is unchanged from 2.2.0. No command, flag, output file or number that migec itself
produces is different. Everything in this release is in the repository: a new documentation page,
the SLURM templates, the Nextflow downstream layer, two more notebooks, and the two scripts behind
them.

The ctDNA ground truth was found, not built

ROADMAP.md item 7 said a ctDNA ground truth had to be constructed, because the runs used by the
published UMI-aware caller benchmark carry no UMI. That is true of those two runs and not of the
archive. Screening SRA read structure instead turned up two studies that kept a real 12 nt inline
UMI — one that migec suggest recovers from base composition alone, with no knowledge of the
protocol:

study	what it is
PRJNA788522	72 runs, cell-free DNA reference material at certified 0 / 0.125 / 0.25 / 1% VAF, crossed with 5/20/80 ng input and 3.3/10/30x depth, three replicates
PRJNA507366	28 runs, six polymerases on the same material, plus 0.031% and 0.0625% VAF

Run over all 100 (scripts/ctdna_titration.py): at 0.125% VAF the input mass decides the
outcome, not the caller.

input	depth	molecules at the site	variant molecules	P(at least 3 copies)
5 ng	3.3x	2,699	3.4	0.65
5 ng	10x	4,530	5.7	0.89
20 ng	3.3x	7,162	9.0	0.99
20 ng	10x	11,961	15.0	1.00

At 5 ng a 0.125% variant is a coin flip and no caller fixes it — a third of replicates simply do
not contain three copies of it. The WT arms are 33 runs and 6.1 M molecules in which the truth is
zero, which is a false-positive floor measured on real chemistry rather than on a simulation.

The barcode error rate tracks polymerase fidelity

PRJNA507366 varies the polymerase while holding template and protocol fixed, which makes it an
independent check on refine's barcode-error estimate. It separates the high-fidelity enzymes from
the standard ones by 3 Phred, a factor of two — Accuprime, Platinum and Phusion at Q30.8–31.4
against Platinum HiFi at Q34.0 and Platinum SuperFi at Q34.8–35.1 — from the barcodes alone, with
no reference, no alignment and no knowledge of which enzyme it was reading.

Also in this release

• docs/variants.rst — which variant caller, and the molecule count that decides whether any
  of them can see it. Standard callers (Mutect2, LoFreq, FreeBayes, VarDict) compose with migec;
  UMI-aware ones (UMI-VarCal, UMIErrorCorrect) replace a stage. Never set a family-size filter
  after assemble: every family has size 1 by construction, so --min-family-size 3 discards the
  entire library and reports zero variants without an error.
• integrations/slurm/ — two sbatch templates and a sample sheet. Both run as ordinary bash
  without SLURM, which is how they are tested and how a layout should be checked before a cohort is
  queued.
• integrations/nextflow/ — a runnable entry point (--mode consensus | ctdna | airr) and a
  downstream layer (align, call variants, arda). The align module exits non-zero if MI:Z: does
  not reach the BAM, because that failure is otherwise silent.
• scripts/sra_fetch.py — fetch benchmark runs on demand instead of mirroring them. NCBI S3 is
  33x faster than ENA here (6.7 MB/s on 8 connections against 200 kB/s), and some runs list their
  original submitted file ahead of the .sra, so fasterq-dump gets handed a gzipped FASTQ
  unless you select on the file type.
• notebooks/ — exome_capture.py, airr_repertoire.py and a README. The PEP 723 headers
  were under-declared: altair renders a polars frame through pyarrow and then pandas, so a clean
  machine failed three imports deep with an error naming neither.

Corrected

"Sequencing deeper buys reads, not molecules" was wrong, and the data says so. The same 20 ng
library gives 6,310 / 10,299 / 16,809 molecules per amplicon at 3.3 / 10 / 30x — depth buys
molecules up to the input ceiling, past which it only raises reads per molecule.

Verified

bash setup.sh from a clean venv, ruff, 217 pytest, ctest and sphinx-build -W all clean; the
SLURM templates were run end to end. Nextflow was not run — it is not installed on the machine
this was written on — so those modules are reviewed against the nf-core module spec rather than
verified by a pipeline run.

Install: pip install migec or uv pip install migec. Wheels for CPython 3.10–3.13 on Linux
x86-64 and macOS arm64.

v2.2.0 — the counter handing out the work, and the error rate at every depth

Throughput

A sampling profile of assemble put 21% of all CPU samples, across every thread, on one instruction: the atomic fetch_add in parallel_for that claimed one item at a time. When an item is one read's tag scan, sixteen cores serialise on one cache line. Items are claimed in batches now, sized so each worker takes ~8 turns, collapsing to 1 when there are few items — the uneven case (one bucket per item) the counter exists for.

Two serial blocks went with it, both read-only 3L-binary-search scans of the barcode table that tally integers: the distance-1 census in estimate_umi_error, which is checkout's per-sample statistics tail, and refine's residual-FDR scan, 0.53 s of a 2.17 s run on one core after everything around it had been parallelised.

stage	1 thread	16 threads	was, at 16
checkout	213,880	1,548,835	1,056,472
refine	617,802	1,554,156	1,012,368
assemble	554,106	2,470,928	2,051,937

reads/s. checkout's gap between end-to-end and matching throughput fell from 20% to 9%.

The barcode error rate, measured at every depth

New table <sample>.umi_errors.tsv, one row per exact parent depth. A parent carrying c reads offered c·L barcode bases to be miscalled, so the same eps falls out two ways:

distinct children  u(c) = 3L (1 - exp(-c eps / 3))    saturates at 3L
reads in children  r(c) = c L eps                     no ceiling

Where they part company is where the library's barcode neighbourhood filled up, measured rather than predicted. On 1e-3 injected: distance-1 excess 9.73e-04, children 9.89e-04 = Q30, all depths 9.98e-04.

New fields error_at_depth, error_phred, error_from_children, error_depth, and two migec plot panels.

Neither estimator is saturation-free — both are bounded by the merges correction actually made. As a fraction of an injected truth:

occupancy	0.2%	2.3%	9.8%	33%	100%
distance-1 excess	0.97	0.96	0.76	0.45	0.001
from the children	0.99	0.95	0.88	0.62	0.00

At 100% both are zero, because correct_umis refuses to merge on a full barcode space and is right to. The saturated flag is what says the answer is a floor.

A sparse spectrum is points, never a line

mig_size_spectrum drew reads in them with lines over a table at exact sizes, where past the head almost every size holds one molecule — so reads == size and the line drew the y = x diagonal as the figure's most prominent feature, a tautology that reads as a second mode.

Compatibility

No breaking change: every 2.1.0 command, flag and output file behaves as it did. The additions are one table, four JSON fields and two plot panels.

Full notes in CHANGELOG.md; the method is written up in docs/umi_errors.rst.

migec 2.1.0 — the four figures everyone already knows how to read

First non-alpha of the rewrite. The three stages, the eight commands and the on-disk formats have
been stable across a3 and a4; what this release adds is the QC layer catching up with them.

Four figures you already know how to read. They come off new tables rather than out of a new
computation, so each is still redrawable from a committed TSV:

• The barcode rank plot, on Cell Ranger's
  axes, because it is the figure every user of a droplet protocol has seen: barcodes sorted by
  content, log-log, knee where cells stop. Never: the y axis is unique UMIs, never reads — one
  over-amplified molecule would otherwise put an empty droplet high on the curve, which is the
  exact artefact the plot exists to show. The call is drawn on the curve, not described in a
  caption. New table <sample>.cell_rank.tsv.
• The MIG size spectrum, molecules and the reads they account for, against log(1 + size).
  Both series, because they peak in different places the moment a library is over-sequenced: most
  molecules are shallow, most reads are in the deep ones, and a figure with only one of them says
  the opposite of a figure with only the other. log1p, so a molecule seen once has a place on
  the axis.
• The rank/Zipf curve, molecule size against rank on log-log. Never: this is why the new
  <sample>.sizes.tsv is written at exact sizes and not power-of-two bins — four bins make
  four steps, and a straight line cannot be told from a bent one. It costs one row per distinct
  depth, a few thousand on a real library, not one row per molecule.
• Unique UMIs and reads per sample barcode, off checkout.summary.tsv. The multiplexed
  analogue of the same question.

Consensus quality is a box, not a thinned scatter. Emitted quality is discrete and capped at
the RT floor, so at any real depth every molecule sits on one or two integers: a cloud of dots
draws that as a flat line whether the bin holds ten molecules or ten million, and the every 17
that kept the SVG small threw away the tails that were the only thing the cloud could have shown.
assemble now accumulates the exact joint distribution of (depth bin, rounded Phred) — both are
small integers, so it is 61 counters per bin per bucket — and assemble.quality_by_depth.tsv
carries real order statistics over every molecule.

Publication defaults on every panel. Transparent background, so one SVG serves a light README,
a dark README and print. One ink colour (#808080), which reads on both. The legend is inside
the plot box: a key in the margin makes every figure wider than its data and is the first thing a
journal asks you to move. Frame is 760x520 rather than 900x560 wide.

The pipeline figure is page-shaped. rankdir = TB and three rank = same groups pinning each
side tool level with the data it reads; it was a 5:1 strip that filled the README column with air.
Transparent, and minibwa is in the downstream box.

Twenty panels now, from sixteen. docs/formats.rst documents every QC table and its columns.

migec 2.0.0a4 — the partition threads, the docs navigate, and the mapping order is weighed

The partition threads, and the memory estimate behind it was wrong. With the consensus already
parallel, assemble's partition was 2.07 s of a 2.69 s run — 77% of the stage, on one thread.
gzip -dc on the same file takes 0.23 s, so five sixths of it was not the inflate: it was the tag
scan, the barcode packing, the record serialisation and the level-1 deflate of each bucket block.
All four run on the workers now, by ownership rather than locking — worker w owns every
bucket with bucket % threads == w for the whole run, so a bucket file has exactly one writer and
no bucket state is shared.

4 M reads, -t 16	before	after
wall clock	2.70 s	1.95 s
partition	2.06 s	1.45 s
reads/s	1,481,946	2,051,937
peak RSS	1,479 MB	789 MB

Half of that was the reader rather than the threading: assign into the chunk rather than clearing
it, because clear() destroys four std::string per record and the reader ends up spending its
time in malloc instead of inflate. One chunk is held rather than one per worker, so it costs ~2 MB
at any -t.

Note: a bigger chunk is 22% faster again — 2,510,241 reads/s at 64 k reads, because
parallel_for starts and joins its threads per call and 4 M reads at 8 k a chunk pays ~15,000
thread creations. It is not taken: 16 MB of resident chunk makes the partition the memory peak on
a finely partitioned shallow library, which breaks the property that a finer partition costs less
rather than more, and test_shallow_memory_is_still_bounded_by_the_bucket catches it. The upgrade
path is a persistent worker pool, not a bigger chunk.

Never: an estimate that nothing checks will be wrong. The constant deciding how finely to cut
the input said a gzipped FASTQ goes resident at 8x its on-disk size. Measured, it is 19x —
a resident record is two heap std::string with their allocator headers and rounded-up buckets,
plus three 8-byte keys, not the 180 bytes of payload. Guessing low is the expensive direction,
because it picks too few buckets and pass 2 holds sixteen of them at once. That single wrong
number is where the 1,479 MB came from.

subsample says what it did. The report now carries the median and the deepest reads per kept
barcode next to the mean, and five kept barcodes with their depths. Never: in key order, not
first-seen order — a barcode with 100 reads appears early about 100x more often than a singleton,
so the head of a file is a sample of the deep MIGs and of nothing else, which is the same trap
subsample exists to avoid, one level down.

minibwa is in the downstream contract, run and counted like the rest: 600/600 records keep
RX/CB/MI through a sorted BAM. Note: the comment flag is -y on minibwa map (the
minimap2 spelling) and -C on the legacy minibwa mem (bwa's), and each rejects the other's
flag with a non-zero exit rather than dropping the tags quietly.

Map first, or collapse first? docs/downstream.rst now works the question through: what the
chromosomal position buys as extra key bits (and when — a 5 nt TSO500 UMI is 1,024 barcodes), what
it costs (N x the alignment, a mismapping becoming a grouping error, needing a reference at all),
and why assemble's linkage sub-clustering recovers most of it from the payload without an
aligner. With a table separating tools that transport a UMI from tools that deduplicate on one:
the first compose with migec, the second replace a stage of it.

The docs navigate. Twenty pages of flat toctree put every long page title in the header; they
are grouped into seven sections now — Installation, Examples, Layouts, Commands, Downstream,
Method, Reference — with landing pages that say what each page answers. Every command page carries
a subtitle (assemble -- one consensus per molecule), and docs/nextflow.rst is new.

Nextflow. Never: containsKey, not ?:. Groovy's elvis treats false as absent, so a
per-sample contig: false against a params default of true silently meant its opposite — the one
direction a per-sample override exists to make possible. --limit-read/--limit-umi now reach
assemble, and the container tag is one param rather than four literals. Note: nextflow is not
installed on the machine this was measured on, so the modules are reviewed against the nf-core
spec, not verified by a run, and the docs say so.

migec 2.0.0a3 — every stage threads, named RT floor, counting mode, sixteen QC panels

Every stage threads now, and the first fix was not a thread. refine and assemble were the
pipeline's bottleneck at ~200 k reads/s against checkout's 1.06 M. The measurement said why: zlib
at its default level 6 spent 1.78 s of refine's 2.14 s run compressing an intermediate that the
next stage decompresses immediately. Level 1 costs 21% more bytes and gave 3x on its own. Then the
parallelism:

stage	before	after (16 threads)
checkout	1,056,472	1,056,472
refine	222,017	1,012,368
assemble	202,977	1,434,573

refine --threads splits the neighbourhood scan, which is a pure function of the barcode table --
it reads no union-find state -- so it parallelises exactly, and the merges it finds are applied
serially afterwards in the original smallest-first order. The result is identical, not merely
equivalent. assemble --threads gives each worker its own bucket; the buckets are independent by
construction because the partition is on the barcode.

Never: -t still changes nothing but the wall clock, on all three stages. That is why the
bucket count is a constant floor of 16 rather than a function of --threads -- if -t chose how
finely the input was cut, it would choose the gzip member boundaries too, and two runs would
produce byte-different files holding identical records. Asserted three ways: in C++ per stage
(tests/cpp/test_parallel_stages.cpp), at the CLI over a full three-stage chain
(tests/synthetic/test_thread_invariance.py), and under the thread sanitizer, which reports
no data race across 104 test cases and 224,116 assertions -- with the instrumentation proven to
fire on a deliberate race in the same helper.

--limit-read N and --limit-umi N on every stage: stop after N reads, or after N distinct
barcodes. For getting an answer out of a 400 GB run in a minute. Never: a limited run is not a
sample, and says so in its own report -- the first N reads of a FASTQ are one corner of one
flowcell. subsample remains the sampler.

The nextflow integration is rebuilt as three modules and a subworkflow
(integrations/nextflow/), nf-core layout, with meta.yml and a stub: for each. One process
meant one resource label for three stages with different shapes and no resume between them: a
failed assemble re-ran the whole demultiplex. refine now carries process_high_memory because
its memory is set by distinct barcodes and by nothing else, and every stage passes task.cpus,
which is safe precisely because the output does not depend on it.

The pre-amplification error floor is named, not guessed. --rt-error takes a fidelity class:

migec assemble ... --rt-error rt        # 1e-4, caps at Q40 -- anything with an RT step (default)
migec assemble ... --rt-error medium    # 1e-5, caps at Q50 -- no RT, an ordinary polymerase
migec assemble ... --rt-error high      # 1e-6, caps at Q60 -- no RT, a proofreading polymerase
migec assemble ... --rt-error 7.37e-5   # or the rate, e.g. TruSight Oncology 500 v2

Which class applies is a property of the protocol, so it is declared rather than fitted, and every
value is cited in SOURCES.md. The default is 1e-4 because that is 10x's stated figure for the
V(D)J RT reaction and what X2 measured here independently (1.54e-4 on SRR1763769). Never: it is
the one-molecule floor and every record migec emits is one molecule. 10x assign Q40 to bases
covered by a single UMI and Q60 only to bases covered by two or more — an RT error is common-mode
within a molecule and independent between them — and combining molecules is arda's job. A
per-molecule record claiming Q50 is claiming two-UMI confidence on one-UMI evidence.

migec assemble --fast: counting mode. The group's most frequent exact sequence, with each
base carrying the best quality any read of that sequence reported for it. No column model, so no
per-base error correction and no sub-clustering — and the RT floor still caps what it claims. For
expression and clonotype abundance, where the deliverable is a molecule count. Measured against
the full path at 8 reads and 5e-3 per base: the column posterior clears essentially every
sequencing error and the majority string keeps what it carried. Refused with --contig, whose
tiling reads share no exact sequence to take a majority over. The per-molecule table gains a
support column: how many of the molecule's reads carried what was emitted.

Coverage into the consensus is capped at 10,000 reads per barcode, which is 10x's rule and
their reasoning — past that the column posterior has long since saturated while the group still
costs time and memory. Never: the cap applies to the reads that are consensed, never to the reads
that are counted. cD and the table's reads column stay the molecule's true depth, because
capping a count would flatten the abundance of exactly the most-amplified molecules.

migec plot: sixteen QC panels, drawn with gnuplot from the tables the stages already wrote.
UMI PWM and information content, barcode quality and its calibration, coverage, trimming, barcode
space, the suggest cycle trace, overrepresented k-mers, the cell rank curve, consensus quality,
error and layout, and thread scaling. It reads no reads and computes nothing, so a figure can be
redrawn from the table beside it long after the FASTQ is gone and can never disagree with the
number in the report. gnuplot is not a Python dependency: without it the .gp scripts are still
written. migec plot joins info and sheet outside the five-pipeline-command rule.

migec suggest now reports overrepresented k-mers, and stitches them back into the sequence
they came from. Run it on a stage's output to find what the trim left behind: an 8-mer occurs by
chance every ~65 kb, so a surviving primer appears as a run of k-mers each shifted one base from
the last. Counted exactly in a flat 4^8 array — no hash map, nothing to size — and measured
against the reads' own base composition, never a flat 1/4, so the table is not just a description
of GC content.

checkout writes checkout.trimming.tsv, the payload length distribution after trimming. A
pattern matched one base off still matches and still trims; it just leaves every payload one base
short, which no counter of matched reads can show. The mean payload length is now a column in the
report.

Fixes

• checkout's report divided by zero on an empty input file, so a run with no reads ended in a
  traceback rather than a summary.
• The benchmark corpus assigned every read to one sample of four (i % 4 with four reads per
  molecule is always 0). The matcher still scored all four patterns, so the throughput figures
  stand, but the per-sample counters and the memory figure were a single sample's. Fixed in
  tests/benchmark/ and in the new scripts/benchmark_threads.py, and the published numbers are
  re-measured: 1,056,472 reads/s end to end and 1,684,654 matching at 16 threads, 217 MB.
• A sample that received no reads reported an effective barcode length of inf. Infinity is not a
  length; it is the absence of one, and it was being written into a numeric TSV column.
• The consensus log-likelihood tables were rebuilt per molecule — 122 transcendentals per group
  against the ~3 per emitted base the posterior actually needs.
• polars was a runtime dependency of the wheel and nothing in the package imported it; it is now
  in the notebooks extra with marimo and matplotlib. The pipeline's only runtime dependency
  is typer.
• Stale documentation corrected: docs/roadmap.rst still called assemble and refine planned;
  docs/checkout.rst said dual-end barcodes were unimplemented; docs/fragmented.rst documented
  --mode {amplicon,fragmented}, which shipped as --contig; the README quoted the retracted
  1.04x position-independence excess (it is 1.0103x) and linked the docs badge at a host the docs
  are not published on.
• _not_yet() and its claim in the skill that subsample exits 2 are gone; it has worked since
  M4.

migec 2.0.0a2 — positional layouts, presets, and a measured downstream contract

Install: pip install --pre migec

Positional is the primary mode

Most libraries fix the barcode at an offset in one read. Saying so no longer takes a sample sheet
or a flag:

migec checkout reads.fq.gz --bc-pattern '^NNNNNNNN' -o out/      # a caret anchors it
migec checkout reads.fq.gz --bc-pattern '0:8'       -o out/      # or a half-open slice
migec checkout R1.fq.gz R2.fq.gz --bc-pattern 'cell:0:16,16:26' -o out/
migec checkout R1.fq.gz R2.fq.gz --preset 10x-v2 -o out/

Slices are half-open and 0-based like Python's, each a UMI slice unless prefixed cell:; gaps
between them become skipped bases, which is what a spacer is.

--max-offset is now automatic and should not be passed. A caret, a slice list, a read
structure and a pattern with nothing to score all anchor at the first base. This was the sharpest
edge in the interface: 10x and TSO500 needed --max-offset 0 typed by hand, and without it every
read was refused — correct, but only after the run had failed.

Presets

umi, migec, primerid, duplex, 10x, 10x-v2, tso500, smarter-umi.
migec sheet --presets prints each with its source. Every one is a published chemistry with a
citable source and a test that compiles it.

Two carry warnings that are part of the preset:

• duplex emits single-strand consensuses. Duplex pairing is not implemented; do not quote a
  duplex error rate from it.
• TSO500's UMI is 5 nt, on R1 only (5M5S+T +T). 1,024 barcodes does not identify a molecule
  on a real ctDNA panel, and TSO500's own pipeline does not claim it does — it groups on the UMI
  and the mapping position, after alignment. migec groups on the barcode and will report the
  space as saturated; on this chemistry that warning is the correct answer.

The downstream contract is measured, not asserted

New docs/downstream.rst, run against real assemble output:

tool	what arrives	result
minimap2 -ax sr -y	name + tags	600/600 records carry RX, CB, MI; valid sorted BAM
bwa mem -C	name + tags	600/600
arda amplicon	name only	AIRR sequence_id is the molecule id
salmon / kallisto	sequence only	counts are molecule counts

Never run alevin, bustools or STARsolo on a consensus FASTQ — they deduplicate from a raw
barcode read that no longer exists. STAR is unverified here: the Homebrew arm64 build reports zero
input reads for any FASTQ, so the failure is the build's.

Also

• A run that matches nothing says so, and says nothing else. Zero assigned reads used to print
  three warnings computed from reads that never arrived. Found by the smarter-umi preset, whose
  first draft scored a three-base anchor — 6.0 bits against an anchored bar of 6.64, so it refused
  every read and reported that as a library property.
• SRP150352 cannot be reprocessed from SRA: its pipeline moves the UMI into the FASTQ header
  and SRA rewrites headers. Confirmed on three runs; migec suggest reports it unprompted.
• Translating from zUMIs (so, from NASC-seq2): ranges there are 1-based and inclusive, slices
  here are 0-based and half-open, so UMI(12-19) is 11:19.
• README badges, docs/layouts.rst, and a Nextflow module that takes migec_preset.

Full notes in CHANGELOG.md.

migec 2.0.0a1

First published build of the C++20 rewrite of MIGEC and MAGERI.

pip install --pre migec

What works

All three stages, validated against real data on four layouts:

layout	data	result
bulk amplicon with an anchor	HIV-1 Primer ID, SRR1763769	pattern recovered by migec suggest unaided; consensus places at HXB2 2,328-2,595
10x droplet, cell barcodes	sc5p_v2_hs_PBMC_1k VDJ-T	100% of 3,155,166 reads assigned, 221,024 barcodes, 813 cells called
dual-end barcodes	MAGERI's design	24 nt UMI from twelve bases at each end
TSO500 read structures	5M5S+T	fgbio/Picard syntax, UMI on both mates

• checkout - barcode extraction with a pattern table, an inline --bc-pattern, or an fgbio --read-structure. 1.18 M reads/s at 16 threads, output byte-identical at any thread count.
• refine - barcode error correction using the evidence that survives at one read per molecule: the barcode's own base quality and payload agreement, not only the count ratio. Whitelists with a measured background prior. OrdMag cell calling.
• assemble - one consensus per molecule, with emitted quality capped by a measured RT/first-cycle-PCR floor of 1.5e-4 rather than the 1e-6 that gets assumed.
• suggest, subsample, sheet.

Numbers that are measured rather than assumed

• The RT/PCR quality floor: 1.54e-4 [1.36e-4, 1.74e-4], so no emitted quality above ~Q38.
• The MIG-split threshold: 8.68, bootstrap 95% CI [8.42, 9.14] - the nominal p < 0.01 over-calls by 19x because reads are not exchangeable.
• What the reported Phred is worth: slope 1.04 over 46.3 M constant pattern bases.
• At ~1 read per UMI, 80% of barcode errors have no observable parent and cannot be corrected by any method. migec says so instead of reporting a corrected count it has no evidence for.

Why an alpha

The roadmap's remaining items are real: index hopping from the i7 x i5 table, .mig bucket output from checkout (which is what would bound the UMI counters at NovaSeq scale), and the published benchmark comparisons against MIGEC v1, MAGERI, UMI-tools, Calib and fgbio.

Java users of the Groovy MIGEC 1.2.9 want the jars on the 1.2.9 release; that code is archived on branch legacy-v1 and at tag v1-final.

1.2.9

• Update logging for assemble routing
• Docs update
• Fix issues with sequence chunk allocation in CdrBlast