Recipe View

This recipe compares alignment before and after QC using two different tools, bwa and bowtie2

0 results • updated 19 months ago by Bioinformatics Admin • needs authorization

This recipe compares alignment before and after QC using two different tools, bwa and bowtie2

Requirements the bioinfo environment and:

pip bio install --upgrade

Download the recipe and run it at command line:

 bash Compare_alignment_rates.sh CP050134 SRR653447

The first parameter is the SRR accession number the second number is the GenBank id for the reference genome.

Interesting cases

 # Ebola virus, bwa better by default
  bash Compare_alignment_rates.sh AF086833  SRR1972739

 # Zika virus Bowtie works better by default
 bash Compare_alignment_rates.sh KY075939 SRR5122835

 # Tobacco etch virus genome
 bash Compare_alignment_rates.sh NC_001555 SRR3466599

 # The alignment is worse after trimming! Whoops!
 bash Compare_alignment_rates.sh VWSY01000001  SRR16249304

If you want to know more about each sample run:

bio search AF086833  SRR1972739

Copy recipe

You need write access to the project to edit.

# Short read aligner - SRA accession, reference GenBank accession
#
# Usage example:
#
#     bash align.sh CP050134 SRR653447
#

# Stop on error
set -uex

# Reference genome accession number.
ACC=$1

# SRA accession number.
SRR=$2

# Reference FASTA.
REF=ref/$ACC.fa

# Create the directory.
mkdir -p ref reads

# Get the reference accession.
bio fetch $ACC -format fasta > $REF

# Get FASTQ data from your query.
fastq-dump -X 10000 --split-files --outdir reads $SRR 1>> log.txt 

# Assign variable names to your FASTQ data.
R1=reads/${SRR}_1.fastq
R2=reads/${SRR}_2.fastq

# The trimmed reads.
TP1=reads/${SRR}_trimmed_1.fastq
TU1=reads/${SRR}_unpaired_1.fastq

TP2=reads/${SRR}_trimmed_2.fastq
TU2=reads/${SRR}_unpaired_2.fastq

# Index the two genomes
bwa index $REF > log.txt 2> log.txt
bowtie2-build $REF $REF > log.txt 2>>log.txt

# Align before doing any quality control on your FASTQ data.
bwa mem $REF $R1 $R2 > bwa.sam 2>> log.txt
bowtie2 -x $REF -1 $R1 -2 $R2 > bowtie.sam 2>> log.txt

# Define your adaptors.
echo ">illumina" > adapter.fa
echo "AGATCGGAAGAG" >> adapter.fa

# Trim off adaptors and low quality reads.
trimmomatic PE ${R1} ${R2} $TP1 $TU1 $TP2 $TU2 ILLUMINACLIP:adapter.fa:2:30:5 SLIDINGWINDOW:4:30 >> log.txt 2>> log.txt

# Align after doing any quality control on your FASTQ data.
bwa mem $REF $TP1 $TP2 > bwa_qc.sam 2>> log.txt
bowtie2 -x $REF -1 $TP1 -2 $TP2 > bowtie_qc.sam 2>> log.txt

# Display certain lines from flagstat output
ls -1 *.sam | parallel -n 1 '(echo {} && samtools flagstat {} | grep mapped)'

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

Name

Recipe display name

Identifier

Unique identifier for the recipe.

Image :

Optional image for the recipe ( 500px Maximum ).

Rank:

Used to order recipes (optional).

Insert Image

From the web

From your computer

Cancel

Back

Compare alignment rates

Insert Image