cluster-rna-reads

Cluster RNA reads (full-length long RNA reads).

Usage

Template:

exacto cluster-rna-reads \
    --bam-file <bam_file> \
    --bai-file <bai_file> \
    --reference-genome-fasta-file <reference_genome_fasta_file> \
    --reference-gene-annotation-file <reference_gene_annotation_file> \
    --reference-gene-annotation-source <reference_gene_annotation_source> \
    --reference-gene-annotation-assembly <reference_gene_annotation_assembly> \
    --reference-gene-annotation-version <reference_gene_annotation_version> \
    --output-dir <output_dir> \
    --output-prefix <output_prefix> \
    [--preset PRESET] \
    [--analyte-type ANALYTE_TYPE] \
    [--remove-unspliced-rnas REMOVE_UNSPLICED_RNAS] \
    [--max-batch-reads MAX_BATCH_READS] \
    [--max-locus-gap MAX_LOCUS_GAP] \
    [--unspliced-bin-size UNSPLICED_BIN_SIZE] \
    [--min-mapping-quality MIN_MAPPING_QUALITY] \
    [--min-terminal-softclip-length MIN_TERMINAL_SOFTCLIP_LENGTH] \
    [--min-reads-per-cluster MIN_READS_PER_CLUSTER] \
    [--max-records MAX_RECORDS] \
    [--expected-sequencing-error EXPECTED_SEQUENCING_ERROR] \
    [--expected-slippage-probability EXPECTED_SLIPPAGE_PROBABILITY] \
    [--max-fpr MAX_FPR] \
    [--max-slippage-repeat-length MAX_SLIPPAGE_REPEAT_LENGTH] \
    [--min-size-proportion MIN_SIZE_PROPORTION] \
    [--max-ins-norm-edit-distance MAX_INS_NORM_EDIT_DISTANCE] \
    [--max-intrachromosomal-distance MAX_INTRACHROMOSOMAL_DISTANCE] \
    [--dna-variants-tsv-file DNA_VARIANTS_TSV_FILE] \
    [--allowed-variants-tsv-file ALLOWED_VARIANTS_TSV_FILE] \
    [--allowed-variant-max-distance ALLOWED_VARIANT_MAX_DISTANCE] \
    [--mec-max-k MEC_MAX_K] \
    [--mec-num-restarts MEC_NUM_RESTARTS] \
    [--mec-max-iter MEC_MAX_ITER] \
    [--mec-seed MEC_SEED] \
    [--ts-min-homology TS_MIN_HOMOLOGY] \
    [--ts-soft-min-homology TS_SOFT_MIN_HOMOLOGY] \
    [--ts-max-breakpoint-dispersion TS_MAX_BREAKPOINT_DISPERSION] \
    [--num-threads NUM_THREADS] \
    [--temp-dir TEMP_DIR] \
    [--max-unspliced-locus-gap MAX_UNSPLICED_LOCUS_GAP] \
    [--soft-clip-removal SOFT_CLIP_REMOVAL] \
    [--soft-clip-max-boundary-distance SOFT_CLIP_MAX_BOUNDARY_DISTANCE] \
    [--soft-clip-bases-per-edit SOFT_CLIP_BASES_PER_EDIT] \
    [--soft-clip-min-partner-bases SOFT_CLIP_MIN_PARTNER_BASES] \
    [--unplaced-tail-removal UNPLACED_TAIL_REMOVAL] \
    [--unplaced-tail-min-ins-len UNPLACED_TAIL_MIN_INS_LEN] \
    [--chunk-size CHUNK_SIZE] \
    [--bkpt-rescue BKPT_RESCUE] \
    [--bkpt-rescue-min-ins-len BKPT_RESCUE_MIN_INS_LEN] \
    [--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-placed-fraction BKPT_RESCUE_REALIGNMENT_MIN_PLACED_FRACTION] \
    [--bkpt-rescue-realignment-max-pieces BKPT_RESCUE_REALIGNMENT_MAX_PIECES] \
    [--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] \
    [--min-reads MIN_READS] \
    [--min-total-depth MIN_TOTAL_DEPTH] \
    [--min-homopolymer-len MIN_HOMOPOLYMER_LEN] \
    [--min-dinucleotide-context-len MIN_DINUCLEOTIDE_CONTEXT_LEN] \
    [--template-switch-flank TEMPLATE_SWITCH_FLANK] \
    [--template-switch-foldback-min-stem TEMPLATE_SWITCH_FOLDBACK_MIN_STEM] \
    [--template-switch-foldback-max-loop-len TEMPLATE_SWITCH_FOLDBACK_MAX_LOOP_LEN] \
    [--template-switch-foldback-max-distance TEMPLATE_SWITCH_FOLDBACK_MAX_DISTANCE] \
    [--template-switch-foldback-slack TEMPLATE_SWITCH_FOLDBACK_SLACK]

Example:

exacto cluster-rna-reads \
    --bam-file tumor_rna.sorted.bam \
    --bai-file tumor_rna.sorted.bam.bai \
    --reference-genome-fasta-file reference_genome.fasta \
    --reference-gene-annotation-file gencode.gtf.gz \
    --reference-gene-annotation-source gencode \
    --reference-gene-annotation-assembly hg38 \
    --reference-gene-annotation-version v45 \
    --output-dir pass1/ \
    --output-prefix tumor_pass1 \
    --preset pb

Description

Cluster RNA reads (full-length long RNA reads).

NoteAt a glance

Inputs: *.bam, *.bam.bai (spliced long-read RNA alignment carrying cs tags), *.fasta (reference genome), *.gtf.gz (gene annotation)

Outputs: {prefix}_exacto_rna_clusters.tsv (read-to-cluster assignments) and {prefix}_exacto_rna_clusters_{summary,reference_transcripts,splice_junctions,template_switch,variants_passed,variants_failed}.tsv under --output-dir

Typical next step: correct-rna-reads after pass 1 on the raw reads (--preset pb or ont); determine-rna-consensus after pass 2 on the realigned corrected reads (--preset corrected --allowed-variants-tsv-file the pass-1 _variants_passed.tsv)

Required arguments

Flag	Type	Description
--bam-file	str	Input BAM file.
--bai-file	str	Input BAM.BAI file.
--reference-genome-fasta-file	str	Reference genome FASTA file.
--reference-gene-annotation-file	str	Reference gene annotation file.
--reference-gene-annotation-source	str	Reference gene annotation source.
--reference-gene-annotation-assembly	str	Reference gene annotation assembly (e.g. ‘hg38’).
--reference-gene-annotation-version	str	Reference gene annotation version (e.g. ‘v41’).
--output-dir	str	Output directory.
--output-prefix	str	Output prefix.

Optional arguments

Flag	Type	Default	Description
--preset	str (pb|ont|corrected)		Error-profile preset that fills typical defaults for –expected-sequencing-error, –expected-slippage-probability, –min-reads-per-cluster, –min-reads, and –min-total-depth. Choices: ‘pb’ (PacBio HiFi), ‘ont’ (Oxford Nanopore), or ‘corrected’ (exacto correct-rna-reads output, for pass-2 re-clustering). Any explicit parameter wins over the preset.
--analyte-type	str (cdna|rna)	cdna	Analyte sequenced: ‘cdna’ (made by reverse transcription) or ‘rna’ (direct RNA sequencing, which turns reverse transcriptase template-switch filtering off).
--remove-unspliced-rnas	str2bool	yes	Remove unspliced RNAs, as remove-unspliced-rnas does, before clustering; pass ‘no’ when the input BAM is already its output.
--max-batch-reads	int	100000	Maximum number of reads to process in a batch.
--max-locus-gap	int	1000000	Maximum locus gap for initial read clustering.
--unspliced-bin-size	int	100000	Unspliced read bin size for initial read clustering.
--min-mapping-quality	int	0	Minimum mapping quality.
--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.
--min-reads-per-cluster	int		Minimum number of reads per splice-junction cluster. This is the junction-cluster read floor only; variant calls are gated separately by –min-reads and –min-total-depth. A flat floor, and the only read-support test a splice-junction cluster has to pass.
--max-records	int	7	Maximum number of records. Read names having more than this value will be excluded.
--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 of a variant call within a cluster. A call needs enough supporting reads that P(X >= reads | site depth, –expected-sequencing-error) falls below this; inside a repeat the error rate is the slippage probability (–expected-slippage-probability, –max-slippage-repeat-length). Retained introns and novel splice junctions, including unannotated consecutive-junction pairs, use the same beta-binomial error test at local read depth, with mean and dispersion both set to –expected-sequencing-error. The separate –min-reads-per-cluster floor also applies.
--max-slippage-repeat-length	int	30	Maximum slippage repeat length.
--min-size-proportion	float	0.5	Minimum size proportion between two insertions for them to cluster. Size proportion = smaller variant size / longer variant size.
--max-ins-norm-edit-distance	float	0.5	Maximum normalized edit distance between two insertion sequences for them to cluster.
--max-intrachromosomal-distance	int	1000	Maximum distance between two insertions for them to cluster.
--dna-variants-tsv-file	str		DNA variants TSV file (an output of exacto call-germline-dna-vars or call-somatic-dna-vars). When supplied, an RNA variant call that matches a listed DNA variant exactly (chromosomes, positions, operations and sequence) passes without the depth, read-support and template-switch tests. Cluster membership and –min-reads-per-cluster are unaffected.
--allowed-variants-tsv-file	str		Allowed variants TSV file: the _exacto_rna_clusters_variants_passed.tsv of a pass-1 exacto cluster-rna-reads run, for pass 2 of the two-pass RNA pipeline. When supplied, a variant call that matches no listed variant (same type and chromosomes, both positions within –allowed-variant-max-distance, an insertion’s sequence within –max-ins-norm-edit-distance) fails as ‘not_in_allowed_list’ and is not phased on. A listed call still faces every other test.
--allowed-variant-max-distance	int	10	Largest distance in bases between a call’s positions and those of the allowed variant it matches; correction and realignment can move a variant between passes.
--mec-max-k	int	7	Maximum number of clusters to identify at the single isoform level for minimum error correction.
--mec-num-restarts	int	1000	Number of restarts for the minimum error correction algorithm.
--mec-max-iter	int	10000	Maximum number of iterations the minimum error correction algorithm.
--mec-seed	int	42	Seed for the minimum error correction algorithm.
--ts-min-homology	int	8	Junction-spanning flank homology (bases) that flags a junction on its own, unless it carries the GT-AG spliceosomal signature.
--ts-soft-min-homology	int	5	Homology at or above this flags only with corroboration: breakpoint dispersion, fold-back geometry, or hub promiscuity.
--ts-max-breakpoint-dispersion	int	2	Pooled member breakpoints spreading more than this many bases beyond the homology interval count as dispersion corroboration.
--num-threads	int	4	Number of threads.
--temp-dir	str		Temp directory.
--max-unspliced-locus-gap	int	0	Maximum allowed distance between two unspliced RNA reads for them to share a locus.
--soft-clip-removal	str2bool	yes	Whether to remove the insertions made of a terminal soft clip that spells the reference across an intron. Such a clip is the end of a degraded read reaching into the neighbouring exon, too short for the aligner to splice. Also resolve a clip that spells a fusion, breakpoint or translocation junction other reads align, as those reads spell it past the junction (see –soft-clip-min-partner-bases), as that junction, and move its read to the one cluster carrying it.
--soft-clip-max-boundary-distance	int	3	Maximum distance in bases between the end of the alignment and the boundary of the intron or junction a terminal soft clip is compared across (for –soft-clip-removal).
--soft-clip-bases-per-edit	int	8	Number of bases of a terminal soft clip for every edit allowed between the clip and the reference across the intron, or the supporting reads across the junction (for –soft-clip-removal). Set 0 to allow no edit.
--soft-clip-min-partner-bases	int	20	Minimum number of bases of a terminal soft clip, past the untemplated bases of a fusion, breakpoint or translocation junction, that must spell the junction’s partner arm for the clip to be resolved as that junction (for –soft-clip-removal). A shorter stretch can match an unrelated arm by chance.
--unplaced-tail-removal	str2bool	yes	Whether to remove the stretch of a read the aligner left unplaced around a short aligned block: a large insertion, a splice junction, the block and a terminal soft clip, where the block is shorter than the insertion and the clip. The junction into the block and the variant records of the stretch are removed.
--unplaced-tail-min-ins-len	int	100	Minimum insertion length for a stretch to be removed as an unplaced tail (for –unplaced-tail-removal).
--chunk-size	int	10000	Number of reads to process in parallel per chunk.
--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-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	11	K-mer size for insertion re-alignment (to rescue breakpoints).
--bkpt-rescue-realignment-band-width	int	50	Band width for insertion re-alignment (to rescue breakpoints).
--bkpt-rescue-realignment-min-score-fraction	float	0.8	Minimum score fraction for insertion re-alignment (to rescue breakpoints).
--bkpt-rescue-realignment-min-query-coverage	float	0.1	Minimum query coverage for insertion re-alignment (to rescue breakpoints).
--bkpt-rescue-realignment-min-placed-fraction	float	0.9	Minimum placed fraction for insertion re-alignment (to rescue breakpoints).
--bkpt-rescue-realignment-max-pieces	int	3	Maximum number of pieces for insertion re-alignment (to rescue breakpoints).
--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).
--min-reads	int		Minimum number of reads supporting a variant call.
--min-total-depth	int		Minimum total depth at a variant site.
--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.
--template-switch-flank	int	20	Flank length in bases examined on each side of a junction for reverse transcriptase template-switch homology.
--template-switch-foldback-min-stem	int	5	Minimum stem length for the reverse transcriptase template-switch fold-back signature.
--template-switch-foldback-max-loop-len	int	20	Maximum loop length for the reverse transcriptase template-switch fold-back signature.
--template-switch-foldback-max-distance	int	1000	Maximum distance (stem + loop + stem) for the reverse transcriptase template-switch fold-back signature.
--template-switch-foldback-slack	int	20	Slack in bases for the reverse transcriptase template-switch fold-back signature.