Update indexing constraints and CLI options
This commit introduces new constraints for kmer size and minimizer selection, defines super-kmers, and adds extensive new command-line options for filtering, conversion, merging, and packing indices. Documentation across the codebase has also been updated.
1 parent
c1e139c597
commit
5313788f7c
10 files changed
+18
-25
No files matched your search
@@ -26,7 +26,6 @@ Two verification modes are available, selected at build time (`index --approx`)
|
|||||||
|
|
||||||
## On-disk layout
|
## On-disk layout
|
||||||
|
|
||||||
```
|
|
||||||
<index_root>/
|
<index_root>/
|
||||||
index.meta global configuration (k, minimizer size, partition count,
|
index.meta global configuration (k, minimizer size, partition count,
|
||||||
evidence mode, whether counts are stored) and genome list/metadata
|
evidence mode, whether counts are stored) and genome list/metadata
|
||||||
@@ -45,7 +44,6 @@ Two verification modes are available, selected at build time (`index --approx`)
|
|||||||
counts/ per-genome kmer counts (if counts were requested)
|
counts/ per-genome kmer counts (if counts were requested)
|
||||||
presence/ per-genome presence/absence bits
|
presence/ per-genome presence/absence bits
|
||||||
layer_1/, layer_2/, ... added by later merges, same internal structure
|
layer_1/, layer_2/, ... added by later merges, same internal structure
|
||||||
```
|
|
||||||
|
|
||||||
`unitigs.bin` is the only file from which the indexed kmer content can be fully recovered; it is always retained. Every other file (MPHF, evidence, counts) is derived from it.
|
`unitigs.bin` is the only file from which the indexed kmer content can be fully recovered; it is always retained. Every other file (MPHF, evidence, counts) is derived from it.
|
||||||
|
|
||||||
|
|||||||
+1
-1
@@ -5,7 +5,7 @@
|
|||||||
Every nucleotide is encoded on 2 bits, most-significant-bit first within each word:
|
Every nucleotide is encoded on 2 bits, most-significant-bit first within each word:
|
||||||
|
|
||||||
| Base | Encoding |
|
| Base | Encoding |
|
||||||
|------|----------|
|
| --- | --- |
|
||||||
| A | `00` |
|
| A | `00` |
|
||||||
| C | `01` |
|
| C | `01` |
|
||||||
| G | `10` |
|
| G | `10` |
|
||||||
|
|||||||
+1
-1
@@ -22,7 +22,7 @@ A value near 0 indicates low complexity (e.g. a homopolymer run); near 1 indicat
|
|||||||
|
|
||||||
## Final score
|
## Final score
|
||||||
|
|
||||||
The filter evaluates $\hat{H}(ws)$ for every word size from 1 to ws_max and keeps the minimum:
|
The filter evaluates $\hat{H}(ws)$ for every word size from 1 to ws\_max and keeps the minimum:
|
||||||
|
|
||||||
$$\text{entropy}(kmer) = \min_{ws=1}^{ws_{\max}} \hat{H}(ws)$$
|
$$\text{entropy}(kmer) = \min_{ws=1}^{ws_{\max}} \hat{H}(ws)$$
|
||||||
|
|
||||||
|
|||||||
@@ -6,9 +6,7 @@ An index is split into a fixed number of **partitions**, each handling an indepe
|
|||||||
|
|
||||||
The canonical minimizer of a super-kmer (see [Minimizer selection](theory-minimizer_selection)) is hashed to produce a $p$-bit routing value that selects the destination partition:
|
The canonical minimizer of a super-kmer (see [Minimizer selection](theory-minimizer_selection)) is hashed to produce a $p$-bit routing value that selects the destination partition:
|
||||||
|
|
||||||
```
|
|
||||||
canonical minimizer → hash(minimizer) → p-bit value → partition index
|
canonical minimizer → hash(minimizer) → p-bit value → partition index
|
||||||
```
|
|
||||||
|
|
||||||
Within a partition, kmers are indexed as plain values via a minimal perfect hash function (see [On-disk storage](formats-index_layout)); the minimizer plays no further role once a super-kmer has reached its partition.
|
Within a partition, kmers are indexed as plain values via a minimal perfect hash function (see [On-disk storage](formats-index_layout)); the minimizer plays no further role once a super-kmer has reached its partition.
|
||||||
|
|
||||||
@@ -17,7 +15,7 @@ Within a partition, kmers are indexed as plain values via a minimal perfect hash
|
|||||||
Even though $H$ already makes minimizer values well-distributed (see [Minimizer selection](theory-minimizer_selection)), choosing $p$ well below the number of bits available in the minimizer ($2m$) leaves a comfortable entropy margin, provided the number of distinct minimizers actually observed is much larger than the number of partitions.
|
Even though $H$ already makes minimizer values well-distributed (see [Minimizer selection](theory-minimizer_selection)), choosing $p$ well below the number of bits available in the minimizer ($2m$) leaves a comfortable entropy margin, provided the number of distinct minimizers actually observed is much larger than the number of partitions.
|
||||||
|
|
||||||
| Minimizer size $m$ | Minimizer bits ($2m$) | Typical partition-index bits $p$ | Partitions |
|
| Minimizer size $m$ | Minimizer bits ($2m$) | Typical partition-index bits $p$ | Partitions |
|
||||||
|----|-----------|-----------|------------|
|
| --- | --- | --- | --- |
|
||||||
| 9 | 18 | 6–8 | 64–256 |
|
| 9 | 18 | 6–8 | 64–256 |
|
||||||
| 11 | 22 | 8–10 | 256–1 024 |
|
| 11 | 22 | 8–10 | 256–1 024 |
|
||||||
| 13 | 26 | 10–12 | 1 024–4 096 |
|
| 13 | 26 | 10–12 | 1 024–4 096 |
|
||||||
|
|||||||
@@ -11,7 +11,7 @@ A **kmer** is a DNA subsequence of fixed length $k$. Two constraints apply to $k
|
|||||||
|
|
||||||
A **super-kmer** is a maximal run of consecutive, overlapping kmers from a read that share the same canonical minimizer (see [Minimizer selection](theory-minimizer_selection)). Each kmer in the run overlaps the next by $k-1$ nucleotides. A super-kmer is capped at 256 nucleotides; a longer run is split at that boundary.
|
A **super-kmer** is a maximal run of consecutive, overlapping kmers from a read that share the same canonical minimizer (see [Minimizer selection](theory-minimizer_selection)). Each kmer in the run overlaps the next by $k-1$ nucleotides. A super-kmer is capped at 256 nucleotides; a longer run is split at that boundary.
|
||||||
|
|
||||||
For a random minimizer of length $m$ over kmers of length $k$, the expected length of a super-kmer is approximately ([Golan & Shur 2025](#ref-Golan2025-xf); [Zheng *et al.* 2020](#ref-Zheng2020-ji)):
|
For a random minimizer of length $m$ over kmers of length $k$, the expected length of a super-kmer is approximately ([Golan & Shur 2025](#ref-Golan2025-xf); [Zheng et al. 2020](#ref-Zheng2020-ji)):
|
||||||
|
|
||||||
$$L_{\text{nt}} \approx \frac{k-m+2}{2} + k - 1$$
|
$$L_{\text{nt}} \approx \frac{k-m+2}{2} + k - 1$$
|
||||||
|
|
||||||
@@ -25,17 +25,17 @@ Super-kmers are the unit of work used throughout construction and querying: sequ
|
|||||||
|
|
||||||
## Bibliography
|
## Bibliography
|
||||||
|
|
||||||
<div id="refs" class="references csl-bib-body hanging-indent" data-entry-spacing="0">
|
<div id="refs" class="references csl-bib-body hanging-indent">
|
||||||
|
|
||||||
<div id="ref-Golan2025-xf" class="csl-entry">
|
<div id="ref-Golan2025-xf" class="csl-entry">
|
||||||
|
|
||||||
Golan, S. & Shur, A.M. (2025). [Expected density of random minimizers](https://doi.org/10.1007/978-3-031-82670-2\_25). In: *Lecture notes in computer science*, Lecture notes in computer science. Springer Nature Switzerland, Cham, pp. 347–360.
|
Golan, S. & Shur, A.M. (2025). <a href="https://doi.org/10.1007/978-3-031-82670-2\_25">Expected density of random minimizers</a>. In: <span style="font-style: italic;">Lecture Notes in Computer Science</span>, Lecture Notes in Computer Science. Springer Nature Switzerland, Cham, pp. 347–360.
|
||||||
|
|
||||||
</div>
|
</div>
|
||||||
|
|
||||||
<div id="ref-Zheng2020-ji" class="csl-entry">
|
<div id="ref-Zheng2020-ji" class="csl-entry">
|
||||||
|
|
||||||
Zheng, H., Kingsford, C. & Marçais, G. (2020). [Improved design and analysis of practical minimizers](https://doi.org/10.1093/bioinformatics/btaa472). *Bioinformatics (Oxford, England)*, 36, i119–i127.
|
Zheng, H., Kingsford, C. & Marçais, G. (2020). <a href="https://doi.org/10.1093/bioinformatics/btaa472">Improved design and analysis of practical minimizers</a>. <span style="font-style: italic;">Bioinformatics (Oxford, England)</span>, 36, i119–i127.
|
||||||
|
|
||||||
</div>
|
</div>
|
||||||
|
|
||||||
|
|||||||
@@ -8,7 +8,7 @@ The minimizer partitions a sequence into super-kmers: maximal runs of overlappin
|
|||||||
|
|
||||||
## Hash-based ("random") minimizer
|
## Hash-based ("random") minimizer
|
||||||
|
|
||||||
`obikmer` selects minimizers by hash order rather than plain lexicographic order. Ordering m-mers lexicographically on their 2-bit encoding systematically favors AT-rich m-mers (an all-A m-mer always encodes to 0), which causes low-complexity regions to dominate as minimizers and produces unbalanced partitions ([Golan & Shur 2025](#ref-Golan2025-xf); [Kille *et al.* 2023](#ref-Kille2023-px); [Pan & Reinert 2024](#ref-Pan2024-hb); [Zheng *et al.* 2020](#ref-Zheng2020-ji), [2021](#ref-Zheng2021-cc)).
|
`obikmer` selects minimizers by hash order rather than plain lexicographic order. Ordering m-mers lexicographically on their 2-bit encoding systematically favors AT-rich m-mers (an all-A m-mer always encodes to 0), which causes low-complexity regions to dominate as minimizers and produces unbalanced partitions ([Golan & Shur 2025](#ref-Golan2025-xf); [Kille et al. 2023](#ref-Kille2023-px); [Pan & Reinert 2024](#ref-Pan2024-hb); [Zheng et al. 2020](#ref-Zheng2020-ji); [2021](#ref-Zheng2021-cc)).
|
||||||
|
|
||||||
Instead, a well-distributed hash function $H$ is applied to the canonical (lexicographically minimal) form of each m-mer, and the m-mer with the smallest $H$ value wins. Because $H$ is a bijection with good avalanche properties, every distinct m-mer in a window has an equal chance of holding the minimum hash value, independent of its nucleotide composition.
|
Instead, a well-distributed hash function $H$ is applied to the canonical (lexicographically minimal) form of each m-mer, and the m-mer with the smallest $H$ value wins. Because $H$ is a bijection with good avalanche properties, every distinct m-mer in a window has an equal chance of holding the minimum hash value, independent of its nucleotide composition.
|
||||||
|
|
||||||
@@ -37,7 +37,7 @@ H(x):
|
|||||||
The choice of $s$ is not arbitrary. Low-complexity m-mers (homopolymers, short tandem repeats) are disproportionately abundant in real genomes; if one of them happened to be the perpetual argmin of $H$ — as the all-A m-mer is when $s = 0$, since $\text{mix64}(0) = 0$ is a fixed point — it would win far more windows than the composition-uniform behavior established above predicts, not because $H$ favors it, but because that pathological input keeps recurring in real sequence data. Exhaustive checks confirm that, with this seed, the argmin is never a homopolymer or any periodic repeat, for every tested $m$:
|
The choice of $s$ is not arbitrary. Low-complexity m-mers (homopolymers, short tandem repeats) are disproportionately abundant in real genomes; if one of them happened to be the perpetual argmin of $H$ — as the all-A m-mer is when $s = 0$, since $\text{mix64}(0) = 0$ is a fixed point — it would win far more windows than the composition-uniform behavior established above predicts, not because $H$ favors it, but because that pathological input keeps recurring in real sequence data. Exhaustive checks confirm that, with this seed, the argmin is never a homopolymer or any periodic repeat, for every tested $m$:
|
||||||
|
|
||||||
| $m$ | argmin (canonical) | decoded sequence | minimal period |
|
| $m$ | argmin (canonical) | decoded sequence | minimal period |
|
||||||
|-----|--------------------|-------------------|----------------|
|
| --- | --- | --- | --- |
|
||||||
| 3 | 16 | `CAA` | 3 |
|
| 3 | 16 | `CAA` | 3 |
|
||||||
| 5 | 78 | `ACATG` | 5 |
|
| 5 | 78 | `ACATG` | 5 |
|
||||||
| 7 | 5512 | `CCCGAGA` | 7 |
|
| 7 | 5512 | `CCCGAGA` | 7 |
|
||||||
@@ -58,35 +58,35 @@ See [Partitioning and indexing architecture](theory-indexing_architecture) for m
|
|||||||
|
|
||||||
## Bibliography
|
## Bibliography
|
||||||
|
|
||||||
<div id="refs" class="references csl-bib-body hanging-indent" data-entry-spacing="0">
|
<div id="refs" class="references csl-bib-body hanging-indent">
|
||||||
|
|
||||||
<div id="ref-Golan2025-xf" class="csl-entry">
|
<div id="ref-Golan2025-xf" class="csl-entry">
|
||||||
|
|
||||||
Golan, S. & Shur, A.M. (2025). [Expected density of random minimizers](https://doi.org/10.1007/978-3-031-82670-2\_25). In: *Lecture notes in computer science*, Lecture notes in computer science. Springer Nature Switzerland, Cham, pp. 347–360.
|
Golan, S. & Shur, A.M. (2025). <a href="https://doi.org/10.1007/978-3-031-82670-2\_25">Expected density of random minimizers</a>. In: <span style="font-style: italic;">Lecture Notes in Computer Science</span>, Lecture Notes in Computer Science. Springer Nature Switzerland, Cham, pp. 347–360.
|
||||||
|
|
||||||
</div>
|
</div>
|
||||||
|
|
||||||
<div id="ref-Kille2023-px" class="csl-entry">
|
<div id="ref-Kille2023-px" class="csl-entry">
|
||||||
|
|
||||||
Kille, B., Garrison, E., Treangen, T.J. & Phillippy, A.M. (2023). [Minmers are a generalization of minimizers that enable unbiased local jaccard estimation](https://doi.org/10.1093/bioinformatics/btad512). *Bioinformatics (Oxford, England)*, 39.
|
Kille, B., Garrison, E., Treangen, T.J. & Phillippy, A.M. (2023). <a href="https://doi.org/10.1093/bioinformatics/btad512">Minmers are a generalization of minimizers that enable unbiased local Jaccard estimation</a>. <span style="font-style: italic;">Bioinformatics (Oxford, England)</span>, 39.
|
||||||
|
|
||||||
</div>
|
</div>
|
||||||
|
|
||||||
<div id="ref-Pan2024-hb" class="csl-entry">
|
<div id="ref-Pan2024-hb" class="csl-entry">
|
||||||
|
|
||||||
Pan, C. & Reinert, K. (2024). [A simple refined DNA minimizer operator enables 2-fold faster computation](https://doi.org/10.1093/bioinformatics/btae045). *Bioinformatics (Oxford, England)*, 40.
|
Pan, C. & Reinert, K. (2024). <a href="https://doi.org/10.1093/bioinformatics/btae045">A simple refined DNA minimizer operator enables 2-fold faster computation</a>. <span style="font-style: italic;">Bioinformatics (Oxford, England)</span>, 40.
|
||||||
|
|
||||||
</div>
|
</div>
|
||||||
|
|
||||||
<div id="ref-Zheng2020-ji" class="csl-entry">
|
<div id="ref-Zheng2020-ji" class="csl-entry">
|
||||||
|
|
||||||
Zheng, H., Kingsford, C. & Marçais, G. (2020). [Improved design and analysis of practical minimizers](https://doi.org/10.1093/bioinformatics/btaa472). *Bioinformatics (Oxford, England)*, 36, i119–i127.
|
Zheng, H., Kingsford, C. & Marçais, G. (2020). <a href="https://doi.org/10.1093/bioinformatics/btaa472">Improved design and analysis of practical minimizers</a>. <span style="font-style: italic;">Bioinformatics (Oxford, England)</span>, 36, i119–i127.
|
||||||
|
|
||||||
</div>
|
</div>
|
||||||
|
|
||||||
<div id="ref-Zheng2021-cc" class="csl-entry">
|
<div id="ref-Zheng2021-cc" class="csl-entry">
|
||||||
|
|
||||||
Zheng, H., Kingsford, C. & Marçais, G. (2021). [Sequence-specific minimizers via polar sets](https://doi.org/10.1093/bioinformatics/btab313). *Bioinformatics (Oxford, England)*, 37, i187–i195.
|
Zheng, H., Kingsford, C. & Marçais, G. (2021). <a href="https://doi.org/10.1093/bioinformatics/btab313">Sequence-specific minimizers via polar sets</a>. <span style="font-style: italic;">Bioinformatics (Oxford, England)</span>, 37, i187–i195.
|
||||||
|
|
||||||
</div>
|
</div>
|
||||||
|
|
||||||
|
|||||||
+1
-1
@@ -20,7 +20,7 @@ obikmer index -o OUTPUT [OPTIONS] [INPUTS...]
|
|||||||
| `--force` | off | Overwrite an existing output directory |
|
| `--force` | off | Overwrite an existing output directory |
|
||||||
| `--label` | input file name without extension | Genome label stored in the index |
|
| `--label` | input file name without extension | Genome label stored in the index |
|
||||||
| `--meta KEY=VALUE` | none | Attach a categorical metadata field to the genome (repeatable) |
|
| `--meta KEY=VALUE` | none | Attach a categorical metadata field to the genome (repeatable) |
|
||||||
| `-k, --kmer-size` | `31` | Kmer size (odd, in [11, 31]) |
|
| `-k, --kmer-size` | `31` | Kmer size (odd, in \[11, 31\]) |
|
||||||
| `-m, --minimizer-size` | `11` | Minimizer size (odd, in $[3, k-1]$) |
|
| `-m, --minimizer-size` | `11` | Minimizer size (odd, in $[3, k-1]$) |
|
||||||
| `--theta` | `0.7` | Entropy threshold for the low-complexity filter |
|
| `--theta` | `0.7` | Entropy threshold for the low-complexity filter |
|
||||||
| `--level-max` | `6` | Maximum sub-word size for the entropy score |
|
| `--level-max` | `6` | Maximum sub-word size for the entropy score |
|
||||||
|
|||||||
+1
-2
@@ -295,9 +295,8 @@ CSV matrix layout (`_dist.csv`, `_shared.csv`, `_family_overlap.csv`): header `g
|
|||||||
**`_iqtree.model`** (`--iqtree`) — lower-triangular exchangeability matrix $R(a,b) = e^{-\text{cost}(a,b)}$ (one row of increasing length per state, whitespace-separated, PAML order), followed by one line of empirical state frequencies. Only states actually occurring in the alignment are kept, compactly renumbered `0..k-1`.
|
**`_iqtree.model`** (`--iqtree`) — lower-triangular exchangeability matrix $R(a,b) = e^{-\text{cost}(a,b)}$ (one row of increasing length per state, whitespace-separated, PAML order), followed by one line of empirical state frequencies. Only states actually occurring in the alignment are kept, compactly renumbered `0..k-1`.
|
||||||
|
|
||||||
**`_iqtree.fasta`** (`--iqtree`) — alignment recoded to that same compact `0..k-1` alphabet (symbols `0-9A-F`). Under `--free-loss`, non-detection becomes `?` and columns left non-informative once missing calls are ignored are dropped first (required for `+ASC`); with `--iqtree-min-freq` also set (the default), any state rarer than that threshold is folded into the same `?` treatment, and non-informative columns are re-checked and dropped again. Run with:
|
**`_iqtree.fasta`** (`--iqtree`) — alignment recoded to that same compact `0..k-1` alphabet (symbols `0-9A-F`). Under `--free-loss`, non-detection becomes `?` and columns left non-informative once missing calls are ignored are dropped first (required for `+ASC`); with `--iqtree-min-freq` also set (the default), any state rarer than that threshold is folded into the same `?` treatment, and non-informative columns are re-checked and dropped again. Run with:
|
||||||
```
|
|
||||||
iqtree3 -s <prefix>_iqtree.fasta --seqtype MORPH -m <prefix>_iqtree.model+ASC --prefix <prefix>_iqtree -T AUTO
|
iqtree3 -s <prefix>_iqtree.fasta --seqtype MORPH -m <prefix>_iqtree.model+ASC --prefix <prefix>_iqtree -T AUTO
|
||||||
```
|
|
||||||
|
|
||||||
**`_iqtree_states.csv`** (`--iqtree`) — one row per state actually kept in `_iqtree.model`/`_iqtree.fasta` (header `iqtree_symbol,canonical_symbol,frequency`): `iqtree_symbol` is the compact `0-9A-F` symbol as written in those two files, `canonical_symbol` is the matching `_sankoff_matrix.csv` state, `frequency` is that state's empirical frequency at full precision. Under `--free-loss`, absent (`0`/`?`) is never a kept state, so it never appears here — nor does any state `--iqtree-min-freq` folded away for being too rare.
|
**`_iqtree_states.csv`** (`--iqtree`) — one row per state actually kept in `_iqtree.model`/`_iqtree.fasta` (header `iqtree_symbol,canonical_symbol,frequency`): `iqtree_symbol` is the compact `0-9A-F` symbol as written in those two files, `canonical_symbol` is the matching `_sankoff_matrix.csv` state, `frequency` is that state's empirical frequency at full precision. Under `--free-loss`, absent (`0`/`?`) is never a kept state, so it never appears here — nor does any state `--iqtree-min-freq` folded away for being too rare.
|
||||||
|
|
||||||
|
|||||||
@@ -20,9 +20,7 @@ Multiple `--ingroup` predicates are combined with AND; multiple `--outgroup` pre
|
|||||||
|
|
||||||
A metadata value is treated as a taxonomy path when it starts with the literal prefix `taxonomy:/`; any other value is treated as a plain string and only supports `=`/`!=`.
|
A metadata value is treated as a taxonomy path when it starts with the literal prefix `taxonomy:/`; any other value is treated as a plain string and only supports `=`/`!=`.
|
||||||
|
|
||||||
```
|
|
||||||
taxonomy:/segment1@rank1/segment2@rank2/...
|
taxonomy:/segment1@rank1/segment2@rank2/...
|
||||||
```
|
|
||||||
|
|
||||||
Each segment is a name, optionally annotated with a rank (e.g. `@family`, `@genus`, `@species`); ranks are optional and can be mixed within a path. The `@` character is reserved inside taxonomy paths and cannot appear in segment names or rank labels.
|
Each segment is a name, optionally annotated with a rank (e.g. `@family`, `@genus`, `@species`); ranks are optional and can be mixed within a path. The `@` character is reserved inside taxonomy paths and cannot appear in segment names or rank labels.
|
||||||
|
|
||||||
|
|||||||
+1
-1
@@ -16,7 +16,7 @@ obikmer superkmer [OPTIONS] [INPUTS...]
|
|||||||
|
|
||||||
| Option | Default | Description |
|
| Option | Default | Description |
|
||||||
| --- | --- | --- |
|
| --- | --- | --- |
|
||||||
| `-k, --kmer-size` | `31` | Kmer size (must be odd, in [11, 31]) |
|
| `-k, --kmer-size` | `31` | Kmer size (must be odd, in \[11, 31\]) |
|
||||||
| `-m, --minimizer-size` | `11` | Minimizer size (must be odd, in $[3, k-1]$) |
|
| `-m, --minimizer-size` | `11` | Minimizer size (must be odd, in $[3, k-1]$) |
|
||||||
| `--theta` | `0.7` | Entropy threshold; kmers with a normalized entropy at or below this value are excluded |
|
| `--theta` | `0.7` | Entropy threshold; kmers with a normalized entropy at or below this value are excluded |
|
||||||
| `--level-max` | `6` | Maximum sub-word size used for the entropy score |
|
| `--level-max` | `6` | Maximum sub-word size used for the entropy score |
|
||||||
|
|||||||
Reference in new issue
Block a user