CLI Reference

This document contains the help content for the twitcher command-line program.

Command Overview:

twitcher↴
twitcher vcf↴
twitcher reads↴
twitcher viz↴

`twitcher`

Template sWITCH alignER – Reannotate variant calls and try to explain mutation clusters with short-range template switches.

Usage: twitcher <COMMAND>

Examples: Annotate a VCF, writing the result to a new VCF: twitcher vcf input.vcf ref.fa -o annotated.vcf

Restrict to a few regions and resolve ambiguous phasing from a BAM: twitcher vcf input.vcf ref.fa -r chr1:1-1000000,chr7 –phasing-bam reads.bam -o annotated.vcf

Scan reads for template switches, then render the hits as SVGs: twitcher reads reads.bam ref.fa -o hits.csv twitcher viz hits.csv -r ref.fa -o svgs/

Subcommands:

vcf — Annotate VCF files with template switch mutations
reads — Find and extract template switch mutations in reads
viz — Visualize template switch mutations

`twitcher vcf`

Annotate VCF files with template switch mutations

Usage: twitcher vcf [OPTIONS] <REFERENCE|--reference <FILE>> <INPUT>

Arguments:

<INPUT> — The input vcf file, or - to read from stdin
<REFERENCE> — The FASTA file containing the reference as the second positional argument. It must be indexed (.fai)

Options:

--reference <FILE> — The reference may also be provided as an argument with –reference at any position
--soft-mask — When enabled, lowercase letters in the reference sequence will be excluded from any alignments
-r, --regions <LIST> — A list of regions to process, comma-separated or by repeating the flag. Each region shall be given in the format chr1:100-200 or chr1. The reader will use the index of the input file to jump to the specified regions. Therefore, the input file must be indexed
-R, --regions-file <FILE> — Same as -r, but the list of regions is specified in a file. Each line corresponds to a region; each region may be either be a BED-like tab-separated definition, or a “position” in the form of chr1:100-200
-t, --targets <LIST> — A list of targets, comma-separated or by repeating the flag. Format is the same as for -r, but the reader will stream the whole input and skip non-matching records, instead of jumping to the specified regions. This may be more performant if the list of targets is long, or if the input has no index. A record counts as targeted when it overlaps with a target region
-T, --targets-file <FILE> — Same as -t, but the list of targets is specified in a file. See also -R
-o, --output <FILE> — The output VCF file. If none is given, or it is -, stdout is used
--only-realigned — When enabled, only the records of clusters in which a template switch was found are written to the output. Every cluster is realigned; this filters on the outcome, not on whether a realignment happened
--realigned-and-context <BASES> — Like –only-realigned, but also write any record within the given number of bases of a record that describes a template switch
--no-version — Do not append version and command line to the header
--output-regions <FILE> — Optionally, a file where to write the regions where there was a template switch. Can be - to output on stdout, but then the output file (-o/–output) must be specified as a file. By default, the regions will be output as a tab-delimited three-column file: (CHROM, POS, END), 1-based and closed. If the specified filename ends with “.bed”, the output will be 0-based and half-open, like usual with the BED file standard
--csv <FILE> — Optionally, a file where to write one CSV row per template switch that was found, in the same format as the output of the reads subcommand. This is the format that the viz subcommand reads
--phasing-bam <FILE> — Optional BAM file (with index) used for local read-based phasing of unphased heterozygous records
--phasing-min-reads <NUM> — Minimum number of reads that must support the winning haplotype assignment. Assignments backed by fewer reads are rejected, guarding against decisions made on too little coverage; the record then stays unphased

Default value: 3
--phasing-confidence <FRACTION> — Minimum fraction of informative reads that must agree on the assignment [0.5–1.0]. “Informative” means the read covers the position and carries a recognisable allele. Higher values reject conflicting or noisy evidence; lower values phase more aggressively. Records below this threshold remain unphased

Default value: 0.8
--phasing-min-mapq <NUM> — Minimum read mapping quality (MAPQ) to use a read as phasing evidence. Reads below this are skipped entirely, so mis-mapped reads cannot contribute spurious allele calls

Default value: 20
--phasing-min-base-qual <NUM> — Minimum base quality (Phred) at the variant position for a read’s allele call to count. A read whose base(s) at the locus fall below this are treated as uninformative there, so likely sequencing errors do not become votes. Deletions (no read base) are unaffected

Default value: 20
-f, --require-flags <FLAG> — Only use as evidence reads that have ALL of the FLAGs present

Default value: 0
-F, --exclude-flags <FLAG> — Only use as evidence reads that have NONE of the FLAGs present

Default value: 3852
--cluster-strategy <STRATEGY> — Clustering strategy. legacy (default) uses log-scaled complexity with a reference-only span and a two-record minimum; edit-mass uses linear edit mass with a symmetric span, a one-record minimum, and the --cluster-min-mass gate

Default value: legacy

Possible values:
- legacy: Original behaviour: log-scaled variant complexity, reference-only span, a minimum of two records, and no edit-mass gate. This is the default
- edit-mass: Linear edit mass, symmetric span (max(reference extent, query extent)), a minimum of one record, and an additional --cluster-min-mass gate. Lets a single large indel or MNV form a cluster on its own
--cluster-min-records <NUM> — Minimum number of mutation events required to form a cluster. When unset, defaults to 2 under the legacy strategy and 1 under edit-mass (where a single sufficiently large event can form a cluster on its own)
--cluster-min-mass <NUM> — Minimum edit mass required for a cluster to be reported (edit-mass strategy only). Edit mass is the sum of edited bases across all events (see [edit_len]): a SNP contributes 1, an n-base MNV/indel contributes n. Defaults to 2.0, which excludes a lone SNP or 1-base indel while admitting any larger single edit or multi-event cluster
--cluster-max-gap <BASES> — Maximum gap in bases between consecutive mutations that still allows them to be in the same cluster. A gap larger than this value will break the cluster

Default value: 20
--threads <NUM> — Number of threads to use for the alignments. Defaults to the number of available cores
-m, --memory-limit <EXPR> — Set an overall memory limit, for example 4gb or 500mb. It is split evenly across the threads. Defaults to the memory currently available on the machine, so results may differ between machines and runs unless it is set explicitly
--aligner-timeout <EXPR> — Set a timeout for each individual run of the aligner, as a number with a unit: ms, s, m or h (e.g. 30s). If an alignment takes longer than that, it will fail and no template switch is found. Unlimited by default
-p, --padding <BASES> — How much context to add on each side of a mutation cluster. This bounds where a template switch may “jump” to (e.g. 100 for +-100 bases around the cluster), but not where it may start and end – that is what –range-extension bounds

Default value: 200
--range-extension <BASES> — How many extra bases on each side of the cluster the alignment range covers. A template switch must start and end inside that range, so this bounds where one may be placed, while –padding bounds where it may jump to. Raising –padding alone does not help an event whose entrance or exit falls outside the range

Default value: 5
-c, --costs <FILE> — A config file for the costs (in custom .tsa format), this is used by the aligner
--no-ts — Force the aligner to NOT align with template switches, so no template switch is ever reported. Useful for performance benchmarks, useless otherwise
--fpa — Use the Four-Point-Aligner, based on dynamic programming. While it has limited capabilities, it is generally faster. It cannot, however, find alignments that have indels in the inner part of the template switch or that have more than one template switch. A best effort is made to keep the costs consistent with the --costs argument, with some exceptions. The most notable one is that since this aligner is not gap-affine, there is no gap open cost. Only the cost for extending a gap is taken into account
--min-length-strategy <ML_STRATEGY> — Ignored with –fpa and –no-ts: Select a strategy to determine lower bounds on the length of a template switch. This is an implementation detail and not part of the cost function. To enforce a minimum (or maximum) length of the inner part of a template switch, use the cost file

Default value: preprocess-filter

Possible values:
- none: Do not compute lower bounds. Turns A* into Dijkstra
- lookahead: Compute a lower bound by looking ahead from each candidate jump
- preprocess-filter: Check for perfect matches of minimum length and filter out imperfect matches. Note that this is only implemented for reverse template switches; with forward ones enabled, prefer lookahead
- preprocess-price: Check for perfect matches of minimum length and compute a lower bound for imperfect matches
- preprocess-lookahead: Check for perfect matches of minimum length and cache the costs of imperfect matches
--uncertainty-range-mode <UR_MODE> — Ignored with –fpa and –no-ts: Select how to compute the uncertainty ranges around template switches. See the aligner repo for details

Default value: equal-cost

Possible values:
- none: Do not compute uncertainty ranges; they are reported as invalid ([-])
- equal-cost: Extend TS entrance and exit as far as possible without increasing the alignment cost
- equal-cost-ignore-geometry: Like equal-cost, but ignore cost increases that are caused by changing the template switch geometry. For example, with this it is possible to report the minimum and maximum length of a perfect match, even if the minimum (maximum) length is smaller than what is specified via --costs
--allow-mixed-descendants — Ignored with –fpa and –no-ts: Allow multiple template switches in the same solution to have different descendants. By default, all template switches in the solution must have the same descendant, but it can be either of the sequences (ref / alt)
--output-failed-alignments <FILE> — If specified, write the region of every failed alignment to this file, one per line, in samtools region format (chr:start-end, 1-based inclusive). The file can be passed directly to samtools/bcftools via -R or -T to retry or inspect failed regions
--output-alignment-stats <FILE> — If specified, write one CSV row per alignment to this file, covering runtime, peak memory, where the result came from (a fresh computation, the in-memory cache or the database) and the statistics reported by the aligner. Clusters that map to the same alignment yield one row each; filter on source == computed to see only the alignments that were actually run
--db <FILE> — Path to a database file (sqlite) used as a persistent cache of alignment results. Results are keyed by the aligner configuration and the aligned sequences, so a rerun with the same settings reuses them instead of recomputing. The file is created if it does not exist
--db-mode <MODE> — How to use the cache file given with –db

Default value: read-write

Possible values:
- read-only: Reuse cached alignments, but do not store newly computed ones
- write-only: Store newly computed alignments, but do not reuse cached ones (recomputes everything)
- read-write: Reuse cached alignments and store newly computed ones
--output-unresolvable <FILE> — Write a TSV report of clusters that could not be resolved into haplotypes.

One line per unresolvable cluster, with three tab-separated columns:
1. the cluster region as contig:start-end (1-based, inclusive)
2. the phasing status of the cluster: trivial (no phasing needed), no reads (no --phasing-bam given), error(...) (phasing failed, the parenthesised text names the cause), or phased[D,D,...] with one decision per record that phasing attempted — OK_DIRECT, OK_FLIP, NO_EVIDENCE, MIN_READS or MIN_CONFIDENCE
3. why the cluster is unresolvable: MoreThanOneUnphasedHet or MissingAllele
The records of these clusters are still passed through to the output VCF unchanged.
--buffer <NUM> — The capacity of the in-memory buffer channels connecting the reader, aligner and writer tasks, in messages. Larger values let the tasks run further ahead of each other at the cost of memory

Default value: 8192
--sample <NAME>
--stat-output <FILE> — Write the counters that are otherwise logged at the end of the run to this file instead, one name: value pair per line
-v, --verbose — Increase logging verbosity
-q, --quiet — Decrease logging verbosity

`twitcher reads`

Find and extract template switch mutations in reads

Usage: twitcher reads [OPTIONS] <REFERENCE|--reference <FILE>> <BAM_FILE>

Arguments:

<BAM_FILE> — The input BAM file
<REFERENCE> — The FASTA file containing the reference as the second positional argument. It must be indexed (.fai)

Options:

-n, --number-of-reads <NUM> — Only analyze the first n reads. Respects -r/-R but not -t/-T. Intended only for debug/testing purposes
--reference <FILE> — The reference may also be provided as an argument with –reference at any position
--soft-mask — When enabled, lowercase letters in the reference sequence will be excluded from any alignments
-r, --regions <LIST> — A list of regions to process, comma-separated or by repeating the flag. Each region shall be given in the format chr1:100-200 or chr1. The reader will use the index of the input file to jump to the specified regions. Therefore, the input file must be indexed
-R, --regions-file <FILE> — Same as -r, but the list of regions is specified in a file. Each line corresponds to a region; each region may be either be a BED-like tab-separated definition, or a “position” in the form of chr1:100-200
-t, --targets <LIST> — A list of targets, comma-separated or by repeating the flag. Format is the same as for -r, but the reader will stream the whole input and skip non-matching records, instead of jumping to the specified regions. This may be more performant if the list of targets is long, or if the input has no index. A record counts as targeted when it overlaps with a target region
-T, --targets-file <FILE> — Same as -t, but the list of targets is specified in a file. See also -R
--cluster-strategy <STRATEGY> — Clustering strategy. legacy (default) uses log-scaled complexity with a reference-only span and a two-record minimum; edit-mass uses linear edit mass with a symmetric span, a one-record minimum, and the --cluster-min-mass gate

Default value: legacy

Possible values:
- legacy: Original behaviour: log-scaled variant complexity, reference-only span, a minimum of two records, and no edit-mass gate. This is the default
- edit-mass: Linear edit mass, symmetric span (max(reference extent, query extent)), a minimum of one record, and an additional --cluster-min-mass gate. Lets a single large indel or MNV form a cluster on its own
--cluster-min-records <NUM> — Minimum number of mutation events required to form a cluster. When unset, defaults to 2 under the legacy strategy and 1 under edit-mass (where a single sufficiently large event can form a cluster on its own)
--cluster-min-mass <NUM> — Minimum edit mass required for a cluster to be reported (edit-mass strategy only). Edit mass is the sum of edited bases across all events (see [edit_len]): a SNP contributes 1, an n-base MNV/indel contributes n. Defaults to 2.0, which excludes a lone SNP or 1-base indel while admitting any larger single edit or multi-event cluster
--cluster-max-gap <BASES> — Maximum gap in bases between consecutive mutations that still allows them to be in the same cluster. A gap larger than this value will break the cluster

Default value: 20
-o, --output <FILE> — Output CSV file; If not specified or “-”, csv will be printed to stdout
-a, --all-clusters — If enabled, write all clusters to the output. Otherwise (the default) only write clusters that are realigned with a template switch
--threads <NUM> — Number of threads to use for the alignments. Defaults to the number of available cores
-m, --memory-limit <EXPR> — Set an overall memory limit, for example 4gb or 500mb. It is split evenly across the threads. Defaults to the memory currently available on the machine, so results may differ between machines and runs unless it is set explicitly
--aligner-timeout <EXPR> — Set a timeout for each individual run of the aligner, as a number with a unit: ms, s, m or h (e.g. 30s). If an alignment takes longer than that, it will fail and no template switch is found. Unlimited by default
-p, --padding <BASES> — How much context to add on each side of a mutation cluster. This bounds where a template switch may “jump” to (e.g. 100 for +-100 bases around the cluster), but not where it may start and end – that is what –range-extension bounds

Default value: 200
--range-extension <BASES> — How many extra bases on each side of the cluster the alignment range covers. A template switch must start and end inside that range, so this bounds where one may be placed, while –padding bounds where it may jump to. Raising –padding alone does not help an event whose entrance or exit falls outside the range

Default value: 5
-c, --costs <FILE> — A config file for the costs (in custom .tsa format), this is used by the aligner
--no-ts — Force the aligner to NOT align with template switches, so no template switch is ever reported. Useful for performance benchmarks, useless otherwise
--fpa — Use the Four-Point-Aligner, based on dynamic programming. While it has limited capabilities, it is generally faster. It cannot, however, find alignments that have indels in the inner part of the template switch or that have more than one template switch. A best effort is made to keep the costs consistent with the --costs argument, with some exceptions. The most notable one is that since this aligner is not gap-affine, there is no gap open cost. Only the cost for extending a gap is taken into account
--min-length-strategy <ML_STRATEGY> — Ignored with –fpa and –no-ts: Select a strategy to determine lower bounds on the length of a template switch. This is an implementation detail and not part of the cost function. To enforce a minimum (or maximum) length of the inner part of a template switch, use the cost file

Default value: preprocess-filter

Possible values:
- none: Do not compute lower bounds. Turns A* into Dijkstra
- lookahead: Compute a lower bound by looking ahead from each candidate jump
- preprocess-filter: Check for perfect matches of minimum length and filter out imperfect matches. Note that this is only implemented for reverse template switches; with forward ones enabled, prefer lookahead
- preprocess-price: Check for perfect matches of minimum length and compute a lower bound for imperfect matches
- preprocess-lookahead: Check for perfect matches of minimum length and cache the costs of imperfect matches
--uncertainty-range-mode <UR_MODE> — Ignored with –fpa and –no-ts: Select how to compute the uncertainty ranges around template switches. See the aligner repo for details

Default value: equal-cost

Possible values:
- none: Do not compute uncertainty ranges; they are reported as invalid ([-])
- equal-cost: Extend TS entrance and exit as far as possible without increasing the alignment cost
- equal-cost-ignore-geometry: Like equal-cost, but ignore cost increases that are caused by changing the template switch geometry. For example, with this it is possible to report the minimum and maximum length of a perfect match, even if the minimum (maximum) length is smaller than what is specified via --costs
--allow-mixed-descendants — Ignored with –fpa and –no-ts: Allow multiple template switches in the same solution to have different descendants. By default, all template switches in the solution must have the same descendant, but it can be either of the sequences (ref / alt)
--output-failed-alignments <FILE> — If specified, write the region of every failed alignment to this file, one per line, in samtools region format (chr:start-end, 1-based inclusive). The file can be passed directly to samtools/bcftools via -R or -T to retry or inspect failed regions
--output-alignment-stats <FILE> — If specified, write one CSV row per alignment to this file, covering runtime, peak memory, where the result came from (a fresh computation, the in-memory cache or the database) and the statistics reported by the aligner. Clusters that map to the same alignment yield one row each; filter on source == computed to see only the alignments that were actually run
--db <FILE> — Path to a database file (sqlite) used as a persistent cache of alignment results. Results are keyed by the aligner configuration and the aligned sequences, so a rerun with the same settings reuses them instead of recomputing. The file is created if it does not exist
--db-mode <MODE> — How to use the cache file given with –db

Default value: read-write

Possible values:
- read-only: Reuse cached alignments, but do not store newly computed ones
- write-only: Store newly computed alignments, but do not reuse cached ones (recomputes everything)
- read-write: Reuse cached alignments and store newly computed ones
--buffer <NUM> — The capacity of the in-memory buffer channels connecting the reader, aligner and writer tasks, in messages. Larger values let the tasks run further ahead of each other at the cost of memory

Default value: 8192
--no-dedup — When specified, no effort is made to deduplicate clusters before aligning them

Default value: false
--min-quality <NUM> — Specify a minimum base quality (Phred-scaled) that all bases in a cluster must have
--stat-output <FILE> — Write the counters that are otherwise logged at the end of the run to this file instead, one name: value pair per line
-v, --verbose — Increase logging verbosity
-q, --quiet — Decrease logging verbosity

`twitcher viz`

Visualize template switch mutations

Usage: twitcher viz [OPTIONS] --output-dir <DIR> --reference <FILE> [INPUT]

Arguments:

<INPUT> — The input csv file. If omitted, read from stdin

Options:

-o, --output-dir <DIR> — The output directory, where one svg file per csv line will be created
-r, --reference <FILE> — The FASTA file containing the reference. It must be indexed (.fai)
-a, --arrows — Render arrows
-c, --more-complement — Render more complement
-z, --zoom <BASES> — “Zoom” in to show only the specified number of bases as context. Default: no zoom
-e, --uncertainty-range-mode <MODE> [alias: equal-cost-range-mode] — Visualize the uncertainty ranges

Default value: none

Possible values:
- none: Do not render TSM uncertainty ranges
- inner-only: Render the TSM uncertainty ranges only on the inner sequences
- full: Render the TSM uncertainty ranges on all sequences. This may be messy if there are overlaps
-f, --overwrite — Overwrite existing output files
-i, --ids <LIST> — IDs to process. Accepts individual IDs or a comma-separated list. See also -I
-I, --id-file <FILE> — File with IDs to process. One ID per line
--stat-output <FILE> — Write the counters that are otherwise logged at the end of the run to this file instead, one name: value pair per line
-v, --verbose — Increase logging verbosity
-q, --quiet — Decrease logging verbosity