Hello all,
I need some help interpreting what I'm seeing in my assembly data.
I'm doing some variant calling in a small number of closely related bacterial genomes: I assembled each strain de novo (there is no finished reference) using the SPAdes assembler, chose one as the reference and aligned the other's reads against it using bwa mem. I then used freebayes (--ploidy 1 --standard-filters -F 0.95 -C 10) to search for SNV's. As a kind of positive control, I also included the bam alignment of the reference genome itself, i.e. reads from the reference strain mapped to itself.
In a perfect world I suppose the resultant VCF would show total congruence between the reference allele at a given site and itself, since both should be referring to the exact same nucleotide. But in fact there were often sites where the two differed. Puzzled, I eyeballed some of these sites in Tablet, looking at the reference's reads mapped to itself, and sure enough there are sites where the nucleotide in the reference sequence and the bases in the reads mapped underneath it differ - sometimes completely (e.g. see screengrab in attachments).*
I understand that neither the assembly nor the mapping processes are going to be perfect, and that there may be errors in both across the genome. But my question is, why would the assembler show an "A" (for example) when there is an abundance of evidence, from the mapper, that in fact it should be a different base?*
Are these kinds of patterns what we might expect when looking at repeated regions in the genome (TE's etc) - i.e. where bwa has incorrectly placed reads from one repeat underneath another, resulting in what appears to be a SNP in the bam alignment? Or is it evidence of a more fundamental cock-up somewhere in the assembly/mapping process...?
Finally, given all of the above, would it be a sensible idea to discard any putative SNPs that are variable within the reference alignment, on the basis that these sites are prone to misalignment and cannot be trusted in any of the samples?
I hope at least some of that makes sense,
Cheers!
I need some help interpreting what I'm seeing in my assembly data.
I'm doing some variant calling in a small number of closely related bacterial genomes: I assembled each strain de novo (there is no finished reference) using the SPAdes assembler, chose one as the reference and aligned the other's reads against it using bwa mem. I then used freebayes (--ploidy 1 --standard-filters -F 0.95 -C 10) to search for SNV's. As a kind of positive control, I also included the bam alignment of the reference genome itself, i.e. reads from the reference strain mapped to itself.
In a perfect world I suppose the resultant VCF would show total congruence between the reference allele at a given site and itself, since both should be referring to the exact same nucleotide. But in fact there were often sites where the two differed. Puzzled, I eyeballed some of these sites in Tablet, looking at the reference's reads mapped to itself, and sure enough there are sites where the nucleotide in the reference sequence and the bases in the reads mapped underneath it differ - sometimes completely (e.g. see screengrab in attachments).*
I understand that neither the assembly nor the mapping processes are going to be perfect, and that there may be errors in both across the genome. But my question is, why would the assembler show an "A" (for example) when there is an abundance of evidence, from the mapper, that in fact it should be a different base?*
Are these kinds of patterns what we might expect when looking at repeated regions in the genome (TE's etc) - i.e. where bwa has incorrectly placed reads from one repeat underneath another, resulting in what appears to be a SNP in the bam alignment? Or is it evidence of a more fundamental cock-up somewhere in the assembly/mapping process...?
Finally, given all of the above, would it be a sensible idea to discard any putative SNPs that are variable within the reference alignment, on the basis that these sites are prone to misalignment and cannot be trusted in any of the samples?
I hope at least some of that makes sense,
Cheers!
Comment