bbduk for mirna mapping

[danova]

Hi,
I have run a miseq run for microrna detection, filtered the data (removed adaptors and ncRNAs (others than microRNAs) and i´m now going to identify reads mapping to microRNAs in Human.

As for the strategy not sure if i have to map to the human genome or i can directly map to the mature.fa database from mirbase. Any suggestion ? The reason i have choosen not to map to the human genome is that not really interested in denovo microRNAS but on the contrary i´m not sure if i have first to remove reads that can map to other parts of the genome.

So far i have decided to align to mature mirbase database.

I´m using the following command:

bbduk.sh ref=mature.fa.gz in=1.fa.gz stats=statsmature.txt out=mirna.fa hdist=0


1/ What is the best approach to identify full match microRNAs (no mismatch) with bbduk.
2/ How to identify isomirs ? Do i need to cluster reads first with vsearch with a 99% identity and then map to mirbase with 100% full match.

Thanks for your comments

[Brian Bushnell]

By default BBDuk will allow one mismatch of the middle base in a kmer (even with "hdist=0"). You need to add the flag "mm=f" to require exact matches. mm stands for "maskmiddle".

If you want to require every kmer in a sequence to also be present in the reference, you can also add "mkf=1" (which stands for min kmer fraction). That will make the process completely intolerant of any errors or remaining adapter sequence, though.

Be sure to specify a kmer length. If you are looking for RNAs as short as 17bp, then set "k=17"; the default is much longer.

So, the full command might look like:

bbduk.sh ref=mature.fa.gz in=1.fa.gz stats=statsmature.txt out=mirna.fa hdist=0 [B]mm=f mkf=1 k=17[/B]

...but you should determine for yourself what value you want for K.