You are viewing an old version of this page. View the current version.

Compare with Current View Page History

« Previous Version 29 Next »

Samtools

Setup output directory.

mkdir bowtie_samtools

If you do not have an alignment file in the SAM format you may want to start with Introduction to mapping.

cp bowtie/REL606.5.sam bowtie_samtools/ 
cp bowtie/REL606.5.fasta bowtie_samtools/

Index the reference file.

samtools faidx bowtie_samtools/REL606.5.fasta

Convert from SAM to BAM format.

samtools view -bS -o bowtie_samtools/REL606.5.bam bowtie/REL606.5.sam |borderStyle=solid}

Sort the BAM file.

samtools sort bowtie_samtools/REL606.5.bam bowtie_samtools/sorted_REL606.5

Output VCF file.

samtools mpileup -uf bowtie_samtools/REL606.5.fasta bowtie_samtools/sorted_REL606.5.bam \|bcftools view -vcg - \> bowtie_samtools/output.vcf 

Produces output.vcf from Bowtie and output.vcf from BWA.
Move all 3(question) bam files and all 3(question) vcf files to lonestar
introduce bedtools.

Developing a pipeline.

Entering all these commands to first create an alignment and then using Samtools to create a VCF format quickly becomes tedious... Let's create a simple bashscript to streamline things.

Mapping with Bowtie - bowtie_samtools.sh
#!/usr/bin/env bash

#Bowtie.
mkdir bowtie
bp_seqconvert --from genbank --to fasta <REL606.5.gbk > bowtie/REL606.5.fasta
bowtie-build bowtie/REL606.5.fasta bowtie/REL606.5
bowtie -t --sam bowtie/REL606.5 -1 SRR030257_1.fastq -2 SRR030257_2.fastq bowtie/REL606.5.sam

#Samtools
mkdir bowtie_samtools
cp bowtie/REL606.5.sam bowtie_samtools/ 
cp bowtie/REL606.5.fasta bowtie_samtools/
samtools faidx bowtie_samtools/REL606.5.fasta
samtools view -bS -o bowtie_samtools/REL606.5.bam bowtie/REL606.5.sam
samtools sort bowtie_samtools/REL606.5.bam bowtie_samtools/sorted_REL606.5
samtools mpileup -uf bowtie_samtools/REL606.5.fasta bowtie_samtools/sorted_REL606.5.bam \|bcftools view -vcg - \> bowtie_samtools/output.vcf 
Mapping with BWA - bwa_samtools.sh
#!/usr/bin/env bash

#BWA
mkdir bwa 
bp_seqconvert --from genbank --to fasta <REL606.5.gbk > bwa/REL606.5.fasta
bwa index bwa/REL606.5.fasta
bwa aln -f bwa/SRR030257_1.sai bwa/REL606.5.fasta SRR030257_1.fastq
bwa aln -f bwa/SRR030257_2.sai bwa/REL606.5.fasta SRR030257_2.fastq
bwa sampe -f bwa/REL606.5.sam bwa/REL606.5.fasta bwa/SRR030257_1.sai bwa/SRR030257_2.sai SRR030257_1.fastq SRR030257_2.fastq

#Samtools
mkdir bwa_samtools
cp bwa/REL606.5.sam bwa_samtools/ 
cp bwa/REL606.5.fasta bwa_samtools/
samtools faidx bwa_samtools/REL606.5.fasta
samtools view -bS -o bwa_samtools/REL606.5.bam bwa/REL606.5.sam 
samtools sort bwa_samtools/REL606.5.bam bwa_samtools/sorted_REL606.5
samtools mpileup -uf bwa_samtools/REL606.5.fasta bwa_samtools/sorted_REL606.5.bam \|bcftools view -vcg - \> bwa_samtools/output.vcf 

Both are available at

/corral-repl/utexas/BioITeam/ngs_course/local/data/mapping_example/scripts
  • No labels