User cookbook

In order to process BAM files and generate file index, the software samtools and tabix need to be installed in the system.

The package sequenza-utils includes several programs and it should support the generation of seqz files using commonly available input files, such as fasta, BAM and vcf files.

To write your own program using the sequenza-utils library, lease refer to the API library interface

Generate GC reference file

The GC content source required to generate seqz files must to be in the wiggle track format (WIG). In order to generate the wig file from any fasta file use the gc_wiggle program.

sequenza-utils gc_wiggle --fasta genome.fa.gz -w 50 -o genome_gc50.wig.gz

From BAM files

Normal and tumor BAM files

sequenza-utils bam2seqz --normal normal_sample.bam --tumor tumor_sample.bam \
    --fasta genome.fa.gz -gc genome_gc50.wig.gz --output sample.seqz.gz

Normal and tumor pileup files

sequenza-utils bam2seqz --normal normal_sample.pielup.gz \
    --tumor tumor_sample.pielup.gz --fasta genome.fa.gz \
    -gc genome_gc50.wig.gz --output sample.seqz.gz --pileup

Without normal, workaround

sequenza-utils bam2seqz --normal tumor_sample.bam --tumor tumor_sample.bam \
    --normal2 non_matching_normal_sample.bam --fasta genome.fa.gz \
    -gc genome_gc50.wig.gz --output sample.seqz.gz

Binning seqz, reduce memory

sequenza-utils seqz_binning --seqz sample.seqz.gz --window 50 \
    -o sample_bin50.seqz.gz

Seqz from VCF files

VCF files with DP and AD tags

sequenza-utils snp2seqz --vcf sample_calls.vcf.gz -gc genome_gc50.wig.gz \
    --output samples.seqz.gz

Mutect/Caveman/Strelka2 preset

sequenza-utils snp2seqz --vcf sample_calls.vcf.gz -gc genome_gc50.wig.gz \
    --preset mutect --output samples.seqz.gz
sequenza-utils snp2seqz --vcf sample_calls.vcf.gz -gc genome_gc50.wig.gz \
    --preset caveman --output samples.seqz.gz
sequenza-utils snp2seqz --vcf sample_calls.vcf.gz -gc genome_gc50.wig.gz \
    --preset strelka2_som --output samples.seqz.gz

Merge seqz files

Non overlapping calls (eg different chromosomes)

gzcat sample_chr1.seqz.gz sample_chr1.seqz.gz | \
    gawk '{if (NR!=1 && $1 != "chromosome") {print $0}}' | bgzip > \
tabix -f -s 1 -b 2 -e 2 -S 1 sample.seqz.gz

Overlapping sample_calls

sequenza-utils seqz_merge --seqz1 sample_somatic.seqz.gz \
    --seqz2 sample_snps.seqz.gz --output samples.seqz.gz