call-somatic-dna-vars

Call somatic DNA variants in a long-read WGS BAM file.

Usage

Template:

exacto call-somatic-dna-vars \
    --bam-file <bam_file> \
    --bai-file <bai_file> \
    --control-bam-files <control_bam_files> \
    --control-bai-files <control_bai_files> \
    --fasta-file <fasta_file> \
    --output-tsv-file <output_tsv_file> \
    [--preset PRESET] \
    [--num-threads NUM_THREADS] \
    [--regions REGIONS [REGIONS ...]] \
    [--min-reads MIN_READS] \
    [--min-mapping-quality MIN_MAPPING_QUALITY] \
    [--min-base-quality MIN_BASE_QUALITY] \
    [--min-total-depth MIN_TOTAL_DEPTH] \
    [--min-alt-allele-fraction MIN_ALT_ALLELE_FRACTION] \
    [--min-size-proportion MIN_SIZE_PROPORTION] \
    [--max-ins-norm-edit-distance MAX_INS_NORM_EDIT_DISTANCE] \
    [--max-intrachromosomal-distance MAX_INTRACHROMOSOMAL_DISTANCE] \
    [--max-slippage-repeat-length MAX_SLIPPAGE_REPEAT_LENGTH] \
    [--apply-infinite-sites-assumption APPLY_INFINITE_SITES_ASSUMPTION] \
    [--chunk-size CHUNK_SIZE] \
    [--max-records MAX_RECORDS] \
    [--expected-variant-allele-fraction EXPECTED_VARIANT_ALLELE_FRACTION] \
    [--expected-mutation-rate EXPECTED_MUTATION_RATE] \
    [--expected-sequencing-error EXPECTED_SEQUENCING_ERROR] \
    [--expected-slippage-probability EXPECTED_SLIPPAGE_PROBABILITY] \
    [--max-fpr MAX_FPR] \
    [--temp-dir TEMP_DIR] \
    [--min-terminal-softclip-length MIN_TERMINAL_SOFTCLIP_LENGTH] \
    [--bkpt-rescue BKPT_RESCUE] \
    [--bkpt-rescue-min-ins-len BKPT_RESCUE_MIN_INS_LEN] \
    [--bkpt-rescue-max-ins-len BKPT_RESCUE_MAX_INS_LEN] \
    [--bkpt-rescue-search-distance BKPT_RESCUE_SEARCH_DISTANCE] \
    [--bkpt-rescue-realignment-gap-open-score BKPT_RESCUE_REALIGNMENT_GAP_OPEN_SCORE] \
    [--bkpt-rescue-realignment-gap-extend-score BKPT_RESCUE_REALIGNMENT_GAP_EXTEND_SCORE] \
    [--bkpt-rescue-realignment-k BKPT_RESCUE_REALIGNMENT_K] \
    [--bkpt-rescue-realignment-band-width BKPT_RESCUE_REALIGNMENT_BAND_WIDTH] \
    [--bkpt-rescue-realignment-min-score-fraction BKPT_RESCUE_REALIGNMENT_MIN_SCORE_FRACTION] \
    [--bkpt-rescue-realignment-min-query-coverage BKPT_RESCUE_REALIGNMENT_MIN_QUERY_COVERAGE] \
    [--bkpt-rescue-realignment-min-span-proportion BKPT_RESCUE_REALIGNMENT_MIN_SPAN_PROPORTION] \
    [--poa-match-score POA_MATCH_SCORE] \
    [--poa-mismatch-score POA_MISMATCH_SCORE] \
    [--poa-gap-open-score POA_GAP_OPEN_SCORE] \
    [--poa-gap-extend-score POA_GAP_EXTEND_SCORE] \
    [--read-depth-max-merge-distance READ_DEPTH_MAX_MERGE_DISTANCE] \
    [--min-homopolymer-len MIN_HOMOPOLYMER_LEN] \
    [--min-dinucleotide-context-len MIN_DINUCLEOTIDE_CONTEXT_LEN] \
    [--max-control-reads MAX_CONTROL_READS]

Example:

exacto call-somatic-dna-vars \
    --bam-file tumor_dna.sorted.bam \
    --bai-file tumor_dna.sorted.bam.bai \
    --control-bam-files normal_dna.sorted.bam \
    --control-bai-files normal_dna.sorted.bam.bai \
    --fasta-file reference_genome.fasta \
    --preset pb \
    --output-tsv-file tumor_somatic_dna_variants.tsv

Description

Call somatic DNA variants in a long-read WGS BAM file.

NoteAt a glance

Inputs: *.bam, *.bam.bai (case), one or more control *.bam / *.bam.bai, *.fasta (reference genome)

Outputs: *.tsv (one row per somatic DNA variant call)

Typical next step: integrate-vars, translate-transcripts

Required arguments

Flag Type Description
--bam-file str Input BAM file.
--bai-file str Input BAM.BAI file.
--control-bam-files str Input control BAM file(s) (e.g. –control-bam-files BAM_FILE_1 BAM_FILE_2).
--control-bai-files str Input control BAM.BAI file(s) (e.g. –control-bam-bai-files BAM_BAI_FILE_1 BAM_BAI_FILE_2).
--fasta-file str Input reference genome FASTA file.
--output-tsv-file str Output TSV file.

Optional arguments

Flag Type Default Description
--preset str (pb|ont) Sequencing-platform preset that fills platform-typical defaults for –expected-sequencing-error, –expected-slippage-probability, –min-reads, and –min-total-depth. Choices: ‘pb’ (PacBio HiFi) or ‘ont’ (Oxford Nanopore). Any explicit parameter wins over the preset.
--num-threads int 4 Number of threads.
--regions str Genomic regions in which to identify variants (e.g. –regions chr1 chr2 or –regions chr1:1-1000000 chr2:1-1000000). If unspecified, Exacto identifies variants in all contigs found in the BAM file (–bam-file BAM_FILE).
--min-reads int Minimum number of supporting reads.
--min-mapping-quality int 4 Minimum mapping quality.
--min-base-quality int 20 Minimum base quality.
--min-total-depth int Minimum total depth.
--min-alt-allele-fraction float 0.05 Minimum alternate allele fraction.
--min-size-proportion float 0.5 Minimum size proportion between two variants. Size proportion = smaller variant size / longer variant size.
--max-ins-norm-edit-distance float 0.5 Maximum insertion normalized edit (Levenshtein) distance. Normalized edit distance = edit distance / longer insertion size.
--max-intrachromosomal-distance int 1000 Maximum distance for clustering intrachromomsomal variants.
--max-slippage-repeat-length int 30 Maximum slippage repeat length.
--apply-infinite-sites-assumption str2bool yes If ‘yes’, apply infinite sites assumption to the variant calling. That is, if a variant in the BAM file shares breakpoint with any of the variant inany of the control BAM files, filter it out.
--chunk-size int 100000 Chunk size for variant calling.
--max-records int 7 Maximum number of records. Read names having more than this value will be excluded.
--expected-variant-allele-fraction float 0.25 Expected variant allele fraction.
--expected-mutation-rate float 1e-06 Expected mutation rate.
--expected-sequencing-error float Expected sequencing error rate.
--expected-slippage-probability float Expected slippage probability. Suggested: 2x the expected sequencing error rate.
--max-fpr float 1e-06 Maximum false positive rate.
--temp-dir str Temp directory.
--min-terminal-softclip-length int 4 Terminal soft-clip runs shorter than this many bases emit no variant record; 1-2 bp terminal clips are aligner anchor jitter, not junction or insertion evidence. Set 0 to disable. Currently has no effect: the somatic pipeline always uses 0 and leaves clip evidence to breakend clustering.
--bkpt-rescue str2bool yes Whether to rescue breakpoints from insertions by re-aligning long insertion sequences.
--bkpt-rescue-min-ins-len int 100 Minimum insertion length to be considered for re-alignment (to rescue breakpoints).
--bkpt-rescue-max-ins-len int 500 Maximum insertion length to be considered for re-alignment (to rescue breakpoints).
--bkpt-rescue-search-distance int 1000 Search distance in bases around an insertion for re-alignment (to rescue breakpoints).
--bkpt-rescue-realignment-gap-open-score int -5 Gap open score for insertion re-alignment (to rescue breakpoints).
--bkpt-rescue-realignment-gap-extend-score int -1 Gap extend score for insertion re-alignment (to rescue breakpoints).
--bkpt-rescue-realignment-k int 19 K-mer size for insertion re-alignment (to rescue breakpoints).
--bkpt-rescue-realignment-band-width int 100 Band width for insertion re-alignment (to rescue breakpoints).
--bkpt-rescue-realignment-min-score-fraction float 0.7 Minimum score fraction for insertion re-alignment (to rescue breakpoints).
--bkpt-rescue-realignment-min-query-coverage float 0.5 Minimum query coverage for insertion re-alignment (to rescue breakpoints).
--bkpt-rescue-realignment-min-span-proportion float 0.9 Minimum span proportion for insertion re-alignment (to rescue breakpoints): the shorter of the insertion length and the reference span its re-aligned sequence covers, over the longer.
--poa-match-score int 0 Partial order alignment match score (to determine consensus insertion sequence).
--poa-mismatch-score int 4 Partial order alignment mismatch score (to determine consensus insertion sequence).
--poa-gap-open-score int 6 Partial order alignment gap open score (to determine consensus insertion sequence).
--poa-gap-extend-score int 2 Partial order alignment gap extend score (to determine consensus insertion sequence).
--read-depth-max-merge-distance int 1000 Positions less than this many bases apart have their read depths counted from one query of the BAM file. Only the number of queries depends on it, not the depths.
--min-homopolymer-len int 4 Minimum homopolymer length for a variant to be treated as homopolymer-context.
--min-dinucleotide-context-len int 4 Minimum dinucleotide repeat context length for a variant to be treated as dinucleotide-context.
--max-control-reads int 0 Maximum number of control reads a variant may have and still be called somatic.