Consensus sequence from mapping
Here’s a set of commands to generate consensus sequence from mapping reads onto reference. This is implementation of the answer to my question “Software for gene recontruction from WGS data” on BioStar.
|
1 2 3 4 5 6 7 |
bwa index database.fasta bwa aln -t $num_threads database.fasta reads.fastq > reads.sai bwa samse database.fasta reads.sai reads.fastq > aln.sam samtools view -b -o aln.bam aln.sam samtools sort aln.bam aln.sorted samtools index aln.sorted.bam samtools mpileup -uf database.fasta aln.sorted.bam |bcftools view -cg - |vcfutils.pl vcf2fq >aln.fastq |
The above code is also available as gist.
