Re-calling¶
The recallers re-run the likelihood-ratio passes on a BAM that already has HMM calls, without running the HMM again. Use them to:
- add TF calls to the output of
fiberhmm-apply; - re-call with a refit or custom emission table;
- change recall settings (
--min-llr, the nucleosome split thresholds, the CpG policy) on an existing call set.
For DddA output of fiberhmm-call, use fiberhmm-call again instead (see
DddA call output); that includes re-calling after
fiberhmm-tag-m5c adds CpG-island states
(DAF-seq).
| Command | Does |
|---|---|
fiberhmm-recall-tfs |
TF recall over the input's MSPs and short footprints; nucleosomes pass through (nuc.Q) |
fiberhmm-recall-nucs |
nucleosome recall, MSP re-derivation, then TF recall and promotion (same as fiberhmm-recall-tfs --recall-nucs) |
Both rebuild each read's observations from its own MM/ML (or DAF
evidence) and sequence, and take the HMM footprints from the input's
ns/nl/as/al, or from its MA when the legacy tags are absent
(--no-legacy-tags output, fibertools Ma). On fiberhmm-apply output,
fiberhmm-recall-nucs gives the same footprint tags as fiberhmm-call with
matching options.
After fiberhmm-apply¶
fiberhmm-apply -i demo/hia5_pacbio.bam --enzyme hia5 --seq pacbio -o out/apply -c 2
fiberhmm-recall-tfs -i out/apply/hia5_pacbio_footprints.bam -o out/apply.tfs.bam \
--enzyme hia5 --seq pacbio -c 2
fiberhmm-recall-nucs -i out/apply/hia5_pacbio_footprints.bam -o out/apply.recalled.bam \
--enzyme hia5 --seq pacbio -c 2
[recall_tfs] enzyme=hia5 mode=pacbio-fiber k=3 min_llr=5.00 uplift=1.00 tf_decoder=multi_interval_v1 cpg_mask=off unify_threshold=90 cores=2 numba=on
[recall_tfs] ML threshold for MM/ML calls: 125
[recall_tfs] processed 300 reads; 300 carried v2 tags; ... TF calls emitted; ... v2 short nucs demoted to tf+
The input needs ns/nl that list every HMM footprint, which is the
fiberhmm-apply default. Output of fiberhmm-apply --harmonize-legacy-tags
holds only the MA-view nucleosomes in ns/nl and marks each such read
with hl:i:1; both recallers stop at the first read with that tag:
error: read 'pacbio_0275' carries hl: the input is fiberhmm-apply --harmonize-legacy-tags output, whose ns/nl no longer list every HMM footprint that recall needs; re-run fiberhmm-apply without --harmonize-legacy-tags (or use fiberhmm-call)
A writer that replaces ns/nl (the recallers, fiberhmm-call) removes
hl; tools that pass the tags through keep it.
fiberhmm-apply output declares its chemistry like fiberhmm-call output,
so a refit table given with -m alone inherits it (see below). The recall
output declares its chemistry and MA-TYPES, keeps the input order, and is
not indexed; run samtools index if the input was sorted.
With a refit or custom table¶
On a BAM that fiberhmm-call or fiberhmm-apply produced, give only the model:
NOTE: fiberhmm-recall-tfs: custom model on an input declaring [assay=daf enzyme=dddb platform=nanopore mode=daf]; using the defaults of --enzyme dddb --seq nanopore with the given model file(s).
The custom model inherits the declared enzyme and platform because its
observation mode (daf) matches, and the run uses that enzyme's defaults
(ML threshold, CpG masking, keep-one, TF presets) with
your table. The output keeps the input's chemistry declaration.
If the table's mode does not match, or you pass a conflicting
--enzyme/--seq, the run stops:
error: fiberhmm-recall-tfs: the input BAM declares chemistry [assay=daf enzyme=dddb platform=nanopore mode=daf] but this run would declare [assay=fiber-seq enzyme=hia5 platform=nanopore mode=nanopore-fiber]. Pass --enzyme/--seq matching the input (with --model for a custom table), or --replace-chemistry to re-declare the output deliberately.
--replace-chemistry drops the input's declaration and declares this run's
own chemistry instead; for a custom table that is
enzyme=custom;platform=unknown, and no enzyme defaults are inherited. Use it
only for a deliberate change of chemistry.
A custom DddA table given without an inherited DddA declaration does not
turn on CpG-aware recall automatically; add --use-m5c.
What is kept and what is regenerated¶
- Regenerated:
nuc(fromrecall-nucs;recall-tfspasses nucleosomes through),msp,tf, and the legacy tags (unless--no-legacy-tags). - Kept unchanged, with their quality bytes: any other
MAgroup, for exampleddda_mcg/ddda_ucgand thedeam+/deam-coverage of a merged duplex. - Short HMM footprints (
nl<--unify-threshold) that a recaller call overlaps move totf; unmatched ones stay innucwithnq= 0. - Records the recaller skips lose this run's call tags rather than keeping stale ones.
The recallers re-derive deaminations from each read. They do not re-apply a
SNP mask or a --reference used at calling time, and say so:
WARNING: the input was called by fiberhmm-call with SNP mask (daf_snp_mask=on/0sites). Recall re-derives deaminations from each read and does not re-apply the SNP mask or reference, so masked sites count as hits again. Re-run fiberhmm-call with the same options for SNP-masked calls.
DddA call output¶
Recalling a DddA BAM from fiberhmm-call does not reproduce the call. With
DddA, fiberhmm-call chooses the TF scan space from the HMM's own footprints
and drops TFs that only the radial nucleosome refinement exposed. Its output
keeps only the refined footprints, so a recall scans different space: on the
demo data, recall-tfs with identical settings changes the TF calls on 370 of
377 reads. The recallers print a warning on such input. Re-run
fiberhmm-call on the BAM instead: it re-runs the HMM, gives the same calls
on its own output, and keeps the ddda_ucg/ddda_mcg island calls, so it
also re-calls DddA after fiberhmm-tag-m5c. Hia5 and DddB recalls of
fiberhmm-call output reproduce the call at the call's ML threshold: for Hia5
PacBio pass --prob-threshold 128, since the recallers default to 125 where
fiberhmm-call uses 128 (see Options that differ).
Options that differ from fiberhmm-call¶
| Option | Default | Notes |
|---|---|---|
--prob-threshold |
125; 248 for Hia5 Nanopore | fiberhmm-call uses 128 for non-Nanopore data |
--input-frame |
auto |
frame of the input's legacy tags. auto applies the shared frame rule: the last footprint writer on each @PG PP chain decides (fibertools: molecular; FiberHMM call/apply/recall: molecular when it declares coord=molecular, else query); with no writer, a coord=molecular declaration means molecular, otherwise query (FiberHMM 2.12 output). When the header has no evidence at all (no footprint writer, no declaration, for example a stripped @PG), the tags are taken as query frame with a warning, except on reads whose MA has a . nuc/msp/tf group, which are molecular. Merged histories that disagree with nothing declared stop the run at the first read with these tags. fibertools Ma tags are always molecular. |
-c/--cores |
1 | 0 = all CPUs |
--chunk-size |
1024 | reads per worker chunk |
--recall-nucs |
off (recall-tfs), on (recall-nucs) |
Input and output may be - (stdin/stdout) for piping. fiberhmm-recall-nucs
handles linear reads only; circular molecules need fiberhmm-call -r.
Output modes¶
- Spec mode (default):
MA/AQcarrynuc,mspandtfwith their scores; the legacyns/nlare refreshed and hold nucleosomes only. - Downstream-compat mode (
--downstream-compat): TF calls go intons/nlnext to the nucleosomes and noMA/AQis written (an existing one is removed). Per-TFtq/el/erare lost. Use it only for tools that read legacy tags alone.
Every option: fiberhmm-recall-tfs,
fiberhmm-recall-nucs.