I am trying to align short reads DNA sequences on a reference database comprising all the prokaryotic genomes reported at the ncbi. Unfortunately I noticed that both BOWTIE and BWA have a maximum limit for the database size and my database cannot be indexed. Is it possible to split the database in more parts and then perform separate indexing? If this is possible, does the alignment program consider all the separate indexing obtained?
It is possible, you literally split the references into multiple sections and then align to each separately, then merge the BAM files. You would however have to recalculate the mapq scores in the merged output from the AS and XS scores. The latest versions of BWA will support larger databases.
One has to ask why you're trying to align to so many genomes? Even if it were possible to index the genomes your searches would probably take forever and need huge amounts of RAM. Plus, how would you deal with reads that match to tens of genomes...?
If this is a metagenomics study, it's not how the data are typically analysed. See this paper for example:
http://www.ploscompbiol.org/article/info%3Adoi%2F10.1371%2Fjournal.pcbi.1000667
In a trade-off between RAM used and hard drive space used, you could split up all prokaryotic genomes into separate, indexed files. Then construct a shell script to call bowtie/bwa on your sequences to each genome serially.
You could then search for alignments of your specific sequence of interest to see which species has a hit. If all species' genomes are concatenated and indexed together, you may filter out hits that you would consider valid because you needed to increase the match filter to get your aligner to run faster through your massive database file.
Sometimes one has to answer a question, not to judge the appropriateness of what one wants to do. First, because in science sometimes who does things in non-classical ways often discover new things, second because maybe you simply don't know or even imagine why one wants to do something. Now, for instance this may be the case. Imagine that you want to get the genome sequence of an organism that you can't grow in the lab and that therefore you also sequence thousands of bacteria that live with it in the environment. If you don't have a close enough reference genome, how do you filter out the sequences that are not likely coming from the organism of interest?
Third because for instance the bowtie1 database size limit was set to 3.6 billions characters (just enough for human?) but it can be a problem if you use e.g. spruce genome which is some 20 G. Also, if you want to index the latest Rfam release (for instance to filter out non-coding RNAs in a small RNA sequencing experiment), you have to split the file and then perform separate alignments!
Therefore I think the question is a good one. Solution with bowtie, is as it has been said, splitting, align separately and cat the outputs. Clearly this can be interesting when using a cluster, such that you can run separate alignment runs on different machines. Anyway, a downvote for bowtie because it should have been really easy for them to add a wrapper to allow to do the splitting and separate alignments and concatenation without asking the user to do it.
- Badges
- Science topic
- Similar topics
- Computer Science and Engineering
- Bioinformatics
- Bioinformatic Tools