FLASH not pairing obvious reads

[JamieHeather]

I've been trying to use FLASH to combine reads that should overlap in read 1 and 2 in a recent MiSeq 2x251 run, but it seems to be missing obvious overlaps.

My data is a pool of four amplicons of different length (roughly 320, 340, 420 and 430 bp long), so I should be getting overlaps from ~70 to 180bp, depending on the amplicon.

When I run FLASH however, it only combines about 30% of the reads. If I look in the files of reads it doesn't combine, I can see reads that should easily be detected - e.g. reads that have 200 bases of 100% identical sequence.

Here's the command I'm using:

flash -m 30 -M 220 -r 251 -f 376 -s 56 file1.fq file2.fq

So -m 30 -M 220 is a minimum overlap size of 30 and maximum overlap of 220 bases (I'm being generous in the size allowance, as my amplicons have a variable length section in the middle). -r is then read length, f is average fragment length, with the standard deviation of s.

If I BLAST some of the paired reads that don't get overlaps found, I can easily find ones with huge overlaps (as well as finding some that did get combined that shouldn't have!).

Any ideas what I might be doing wrong?

[JackieBadger]

I'm guessing (from your previous post) that these are MHC amplicons?
I would say that your different length amplicons should be conserved enough where you do not need to specify the over lap.
Have you not tried running flash with no preset options?
I have done this on MHC and it correctly joins say 99% of PE reads.

[JamieHeather]

They're actually TCR amplicons (hence the MHC interest!), so there is a fair bit of length deviation even within an amplicon.

I did try just default options, which produces a similar 30% success rate (a % or so less than the one above) - which I find odd, as the default max overlap is 70bp, which should only really just allow some of my smaller amplicons.

[JackieBadger]

I would try to demultiplex your amplicon products by size as much as possible just so you can try to figure out the issue.
First off, are all of your PE1 and PE2 sequences in corresponding order between the files?
Can you de-multiplex these based on primer sequence i.e. does each amplicon of a particular length have its own primer set?

[JamieHeather]

Yep, the reads are in order in each file. Thanks, I'll try demultiplexing on the primers and seeing how they do!

[JamieHeather]

Ah-ha, now this is interesting - once demultiplexed two of the amplicons give combine at >90%, while the other two don't even manage a % between them! Demultiplexing was great advice.

It's the shorter two (~330 base amplicons, therefore ~170 bp overlap) that fail. Still having no luck tweaking these with the different parameters either.

[GenoMax]

There are a few additional tools posted that you could give a try: http://thegenomefactory.blogspot.com...aired-end.html

[JackieBadger]

jamie ping me a mail...I'd like to discuss what you're up to [email protected]

[JamieHeather]

Cheers Genomax - I had tried a few of the others from that, particularly Cope and Pandaseq. Sadly Cope did worse than FLASH, and I'm not convinced Pandaseq installed correctly - it kept throwing up odd errors that the read IDs weren't identical, when they clearly were (even if I made the read names exactly the same, so both :1, or removing that whole block).

I'm just now trying to see how Cope works on the different demultiplexed files.

[JamieHeather]

In case anyone has a similar problem and finds this thread in the future, here's how it got solved.

Basically, as the reads are so long and the amplicons are so short, the end of each read goes past the start of the other read, so the 3' end of read has no equivalent in the other file. As this is what FLASH (and presumably the others) looks for, it doesn't find the huge overlap in the 5'.

It's easily solved - I just used FASTX-Toolkit to trim the reads down to a size where they will only overlap in their 3's, before running FLASH, e.g.:

fastx_trimmer -Q33 -f 1 -l 175 -i R1.fq -o shortR1.fq
fastx_trimmer -Q33 -f 1 -l 175 -i R2.fq -o shortR2.fq
flash shortR1.fq shortR2.fq

Works great!

[JackieBadger]

Now the 'fun' really begins for you!