Recipe View

Bacterial variant caller

2 results • updated 4.0 years ago by Istvan Albert

Bacterial variant caller

This recipe contains the code presented by Torsten Seemann in the blog post titled A Unix one-liner to call bacterial variants

The code from the above website has been generalized a bit with additional utility.

This recipe provides code that:

Downloads a reference genome (default value: NZ_CP008918 Pasteurella multocida strain ATCC 43137, complete genome)
Downloads a sequencing run from SRA (default value: SRR4124989)
Generates statistics on the sequencing data
Runs minimap2 to align the sequences and create a SAM output
Runs bcftools mpileup to generate the genotype likelihoods of each base
Runs bcftools call to filter for multiallelic variants only
Runs bcftools norm to normalize each variant to a standard form
Runs bcftools filter to remove variants with low quality (QUAL) or low coverage (DP)

Copy recipe

You need write access to the project to edit.

# Stop on errors.
set -uexo pipefail

# The reference genome.
ACC=NZ_CP008918

# Select sequencing run.
SRR=SRR4124989

# The directory that stores the reference sequence.
mkdir -p ref

# The file that stores the reference genome.
REF=ref/${ACC}.fa

# Download the reference accession number.
efetch -db nuccore -id $ACC -format fasta > $REF

# The directory that stores the sequencing reads.
mkdir -p reads

# Download the sequencing data.
fastq-dump --origfmt  --split-files -O reads $SRR

# Optional step. Generate a report on the data size
echo "*** Sequence statistics"
seqkit stat reads/${SRR}_?.fastq
echo "*** "
#
# The one-liner code comes from:
#
# https://thegenomefactory.blogspot.com/2018/10/a-unix-one-liner-to-call-bacterial.html
#

# The number of CPUs to use.
CPUS=4

# The names for the read pairs.
R1=reads/${SRR}_1.fastq
R2=reads/${SRR}_2.fastq

# Unix trick: A single command line that is too long may be broken
# into several lines as long as each ends with a \ symbol.
#
# Call the variants.
minimap2 -a -x sr -t $CPUS $REF $R1 $R2 \
 | samtools sort -l 0 --threads $CPUS \
 | bcftools mpileup -Ou -B --min-MQ 60 -f $REF - \
 | bcftools call -Ou -v -m - \
 | bcftools norm -Ou -f $REF -d all - \
 | bcftools filter -Ov -e 'QUAL<40 || DP<10 || GT!="1/1"' \
 > variants.vcf

# Generate a statistics on variants.
bcftools stats variants.vcf > variants-stats.txt

# Plot the variant statistics.
plot-vcfstats -P -p plots variants-stats.txt

You need write access to the original recipe to edit.

Click the buttons on the right to create new fields.

Add text field Add float field Add data field Add checkbox Add dropdown Add upload field Add integer field Add radio button

Edit the content of each interface element.

[acc]
label = "Genome Accession Number"
display = "TEXTBOX"
value = "NZ_CP008918"
regex = "^\\w{1,20}$"
help = "An NCBI accession number"

[srr]
label = "SRA Run Number"
display = "TEXTBOX"
value = "SRR4124989"
help = "An SRR Run number"

[settings]
name = "Bacterial Variant Calling"
template = "Bacterial_Variant_Calling_cookbook_6.sh"
image = "Bacterial_Variant_Calling_cookbook_6.png"
id = 6
recipe_uid = "bacterialsnp"
uid = "bacterialsnp"
help = "Bacterial variant caller\n\nThis recipe contains the code presented by Torsten Seemann in the blog post titled [A Unix one-liner to call bacterial variants][code]\n\nThe code from the above website has been generalized a bit with additional utility.\n\nThis recipe provides code that:\n\n1. Downloads a reference genome (default value: [NZ_CP008918][NZ_CP008918] *Pasteurella multocida strain ATCC 43137, complete genome*)\n1. Downloads a sequencing run from SRA (default value: [SRR4124989][SRR4124989])\n1. Generates statistics on the sequencing data\n1. Runs `minimap2` to align the sequences and create a SAM output\n1. Runs `bcftools mpileup` to generate the genotype likelihoods of each base\n1. Runs `bcftools call` to filter for multiallelic variants only\n1. Runs `bcftools norm` to normalize each variant to a standard form\n1. Runs `bcftools filter` to remove variants with low quality (QUAL) or low coverage (DP)\n\n[NZ_CP008918]: https://www.ncbi.nlm.nih.gov/nuccore/NZ_CP008918.1\n[SRR4124989]: https://www.ncbi.nlm.nih.gov/sra/?term=SRR4124989\n\n[code]: https://thegenomefactory.blogspot.com/2018/10/a-unix-one-liner-to-call-bacterial.html"
url = "http://localhost8000"

You need write access to the original recipe to edit.

Name

Recipe display name

Identifier

Unique identifier for the recipe.

Image :

Optional image for the recipe ( 500px Maximum ).

Rank:

Used to order recipes (optional).

Bacterial variant caller

This recipe contains the code presented by Torsten Seemann in the blog post titled [A Unix one-liner to call bacterial variants][code]

The code from the above website has been generalized a bit with additional utility.

This recipe provides code that:

1. Downloads a reference genome (default value: [NZ_CP008918][NZ_CP008918] *Pasteurella multocida strain ATCC 43137, complete genome*)
1. Downloads a sequencing run from SRA (default value: [SRR4124989][SRR4124989])
1. Generates statistics on the sequencing data
1. Runs `minimap2` to align the sequences and create a SAM output
1. Runs `bcftools mpileup` to generate the genotype likelihoods of each base
1. Runs `bcftools call` to filter for multiallelic variants only
1. Runs `bcftools norm` to normalize each variant to a standard form
1. Runs `bcftools filter` to remove variants with low quality (QUAL) or low coverage (DP)

[NZ_CP008918]: https://www.ncbi.nlm.nih.gov/nuccore/NZ_CP008918.1
[SRR4124989]: https://www.ncbi.nlm.nih.gov/sra/?term=SRR4124989

[code]: https://thegenomefactory.blogspot.com/2018/10/a-unix-one-liner-to-call-bacterial.html

Insert Image

From the web

From your computer

Cancel

Back