SNP calling with paired and single read data

[mht]

Hi all,

I'm working with a set of Illumina paired-end data and my reference genome is a non-model organism made up of only scaffolds and contigs at the moment.
My overall aim is to use the reference as a guide for assembly as well as SNP calling.

My first question is, should I do assembly first and use it for SNP calling, or are these two processes separate and independent?

For now, I've decided to work on the SNP calling part first, by mapping my sequence reads (using BWA) to the set of publicly available reference scaffold and contigs. But before that, I did a quality trim step which left me with paired and unpaired reads (single).

My second question is, how do I deal with the unpaired reads? Can I also use them for SNP calling and later merge my SNP results with the ones I obtain with paired reads?

Any help/suggestion is much appreciated.
Thanks!

[TiborNagy]

1. You can use assembly or BWA mapping for SNP calling. So the two process are separate.
2. I usually does not use unpaired reads, because BWA can make higher quality mappings with paired reads. But you have only scaffolds, so I think you could use those reads.

[mht]

Much thanks TiborNagy for your reply.

I'm wondering what the difference will be if I use my assembly for SNP calling versus BWA mapping. I'm thinking of pros and cons of doing SNP calling based on assembly:

pro: if my sample is not so similar to the reference, then calling SNPs after assembly would be better.
con: the SNP results I get will be dependent on the quality of my assembly.

Any advice on this?

[TiborNagy]

Well, if there is no consensus, I think your assembly maybe has much more benefits than a simple BWA mapping. Do not forget that the quality of the contigs also limit the SNP results.

[mht]

Hi TiborNagy,

I was wondering if you know of any good open source reference assembly programs that I can use for references which contain scaffolds/contigs and not full assemblies?

[TiborNagy]

Hi mht,

My personal favourite is MIRA.

[mht]

Hi TiborNagy,

I've used MIRA and am now wondering what I should do next and how to judge the quality of my assembly. I can calculate the N50 and total size of the original reference, and compare it with my assembly. Should the N50 and total size increase as compared to the original reference? I notice that MIRA did mostly mapping back to the backbone, but also did some extensions to each backbone contig.

Also, since it is already a reference assembly, is it necessary to try to join the contigs to see if any overlap (since some contigs have been extended). Are there any programs you can recommend?

Thanks so much for your help.

[TiborNagy]

It is hard to say anything about the comparison of your assembly and the original one. It depends on many factors. For your second question: I can not recommend other programs.