Reads on contigs

[George Willian Condomitti]

Hi guys,

I assembled with SSQ=1, pe-sam mp-sam options, and then I got the .SAM output file, which gives me detailed information regarding reads on contigs mapping.

I used samtools to generate BAM file from the generated "-3.sam" file, and then indexed it and finally called "samtools idxstats", which (if I did not misunderstand) gives me the how many reads mapped to each contig.

I'm a bit confused now because some contigs have 0 as mapped reads. Shouldn't it have at least 1 read mapped back to it?

Thanks a lot.


Cheers,
Condomitti.

[Tony Raymond]

Hi Condomitti,

The SAM file generated by abyss-pe will only contain reads mapping to different contigs. This would be why you'll see some contigs with 0 mapped reads. 

The tool that filters the reads like this is abyss-fixmate. If you want to get all of the alignments you'll need to run abyss-map separately. 

Cheers,

Tony

--
You received this message because you are subscribed to the Google Groups "ABySS" group.
To unsubscribe from this group and stop receiving emails from it, send an email to abyss-users...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Tony Raymond]

Hi Condomitti,

It was bugging me that there wasn't a way to make abyss-pe behave like you wanted. So I've made a change to abyss-fixmate that will let you do what you were expecting abyss-pe to do.

If you felt adventurous, I've attached the patch that you could apply to the abyss source and then recompile. It adds an option to abyss-fixmate (--all) to print all of the alignments to the SAM files. To apply the patch run the following in the source root directory:

patch -p1 <all_opt.patch

Then to assemble, getting all of the alignments, you'd specify all the same options, but add `FIXMATE_OPTIONS=--all`. 

Cheers,

Tony

[George Willian Condomitti]

Tony, thanks a lot mate! I really appreciate your help!


I'm gonna ask the administrator to apply the patch and recompile or give me privileges to do so.
I'll get back to you soon with the new results.

Thanks again mate!

Cheers,
Condomitti.

[Tony Raymond]

Hi Condomitti,

I sent out that patch a little hastily. There seems to be some issues that I didn't deal with downstream of abyss-fixmate. I'll send you another patch shortly.

Tony

 

[George Willian Condomitti]

Hi Tony,

Ok, no problem. Thanks a lot!


Cheers,
Condomitti.

[Tony Raymond]

Hi Condomitti,

Attached is the patch that I was mentioning earlier as well as the original one. It makes it so that abyss-fixmate handles mateless reads correctly. You'll want to apply both patches one after another and recompile. i.e.: 

patch -p1 <all_opt.patch && patch -p1 <mateless.patch

Let me know how it goes.

Cheers,

Tony

[George Willian Condomitti]

Hello Tony,

Thanks a lot for the patches! I have applied them and recompiled Abyss.

I started a new assembly three days ago... I'll let you know the results as soon as the assembly finishes.

Have a great week.

Cheers,
Condomitti.

[George Willian Condomitti]

Hello Tony,

I just got this error message from abyss-fixmate:


abyss-fixmate: error: All reads are mateless. This can happen when first and second read IDs do not match.
make: *** [mp1-6.sam.gz] Error 1
make: *** Deleting file `mp1-6.sam.gz'


The weird thing here is that I haven't changed anything of the input reads. They are exactly the same I was using initially and getting to the end of the assembly with no error messages.
Any clue?



Thanks.

Cheers,
George Condomitti.

[George Willian Condomitti]

From the stdout I can see the call to abyss-fixmate:

abyss-map-ssq -v -j70 -l32  MP_READS_CORRECTED_J1_R1.fq MP_READS_CORRECTED_J1_R2.fq snake-6.fa \
                |abyss-fixmate-ssq -v -l32 --all -h mp1-6.hist \
                |sort -snk3 -k4 \
                |gzip >mp1-6.sam.gz


Those mate-pair reads are the ones I was using before.



Condomitti.

[Tony Raymond]

Hi Condomitti,

That's weird! Please send the first few lines of each file (i.e. `head MP_READS_CORRECTED_J1_R1.fq MP_READS_CORRECTED_J1_R2.fq).

Thanks,

Tony

[George Willian Condomitti]

Hi Tony,

It seems there was a misspelling in the file names. Sorry about that!

Once fixed, will abyss-pe be able to continue from the point it stopped (once intermediary files are already there) or should I manually call subsequent scripts (-fixmate etc)?

Thanks.

Cheers,
Condomitti.

[Tony Raymond]

Hi Condomitti,

I think you'll want to delete the *-6.hist and *-6.sam.gz files if present. Other than that, abyss-pe should pick up where it left off.

Cheers,

Tony

[George Willian Condomitti]

Thanks Tony,

Everything worked like a charm thanks to your great help.

I'm processing the .SAM file now to obtain all the necessary information.=)

Cheers,
Condomitti.

[George Willian Condomitti]

Hello Tony,

I've been working and analyzing the SAM file generated after your patch and everything seems ok. Thanks!

I'm just a bit confused because there are still loads of entries with no reads mapped to it. I'm using samtools idxstats to get to the mapping stage, and then I'm looking at the third column which shows the numbers of mapped reads. This is the one I'm referring to when I say there are some with 0 on it. Is this the expected behavior? If yes, was that something wrong I've set to abyss so that it didn't map all of the contigs? I'm using the option you told me to use: FIXMATE_OPTIONS=--all




Cheers,
Condomitti.

[Ka Ming Nip]

Hi Condomitti,

abyss-map is a very specialized aligner for ABySS assemblies.
It only produces a single local alignment for each read and it soft-clips at mismatch positions.
If a read maps to multiple locations, then abyss-map only report the first (best) alignment if I remember correctly.
In addition, abyss-map will not report alignments for the soft-clipped regions of the primary alignment by default.
You can use the '--multi' option to turn on this feature, but abyss-fixmate wouldn't work with these alignments as it expects one alignment per read...

Regards,
Ka Ming

[Tony Raymond]

Hi Condomitti,

I don't think you are doing anything wrong :)

As Ka Ming says, when a read multi maps with abyss-map, only one alignment is reported. The first/best matching sequence is used, and the mapping quality is set to 0. Every contig really should contain some unique sequence so hopefully this isn't the source of conitigs having 0 reads mapped.

Another possibility for a contig to have 0 reads mapped would be when you have a small contig (say, k-mer length) with longer contigs on either side. In this case, all of the reads that could align to the small contig would align to the other flanking contigs better. If you only considered contigs that are longer than the read length for your analysis, do you still see contigs with 0 reads mapped? If so, what proportion of them are there? I can imagine there are cases where the unique portion of a contig is also very small and your reads map someplace else better, but this should be quite rare.

Cheers,

Tony

[George Willian Condomitti]

Thanks a lot, Ka Ming and Tony for your replies!

Now it's clear about the "0" mapped =).

Tony I haven't done an analysis of the contigs that are longer than the read length but I'll do it now - good point, thanks!
I'll eliminate those very small contigs and check the numbers.


Have a good weekend mates!


Cheers,
Condomitti.