De novo genome assembly for non-model organisms face the issue of bacterial contamination. For assembled contigs with mostly bacterial-like sequences (based on BLAST search), the entire contig can be dropped off. But for a long contig with small fragments of bacterial-like sequences, one way to solve it is to break the contigs into smaller pieces and drop the bacterial-like pieces. It would work better if the raw reads could be removed BEFORE assembly. This was commonly done for Illumina reads, but for PacBio long reads (non-HiFi reads) with a higher error rate, is it possible to filter contaminated bacterial sequences based on sequence similarity, like BLAST, Kraken? Of course, the potentially contaminated region should be long enough to drop the entire sequence, for PacBio long reads. Maybe 30%? Any suggestion of the processing details?
2 Answers
I don't think that in principle the issue is any different for PacBio and Illumina; it's just that the process might be noisier with higher error rates. Note that PacBio errors aren't randomly distributed substitutions and probably many (most?) short-ish k-mers (<30nt) will likely be intact. I would suggest doing k-mer filtering of your reads as suggested here or here or here.
Alternately, use your favorite alignment-based technique, if you have a good idea of what your contaminant is. In principle PB reads will be better filtered by this method than Illumina as local alignment will be able to figure out the errors quite well with longer queries.
If you go the route of trying to do alignment to some set of contaminant references, the obvious answer would be to use minimap2 as your aligner as this is exactly the task that it is designed for. You can then use samtools to filter out anything that maps well.
Here is a related question/answer that also addresses the contaminant issue and suggests using Blobology, which I'm not familiar with but may be helpful (granted it also seems to be Illumina-centric). Genome-size estimation using k-mers on PB reads appears to be pretty accurate too, so I think that it's likely ok to work with raw PB k-mers.
-
1$\begingroup$ Right, the PacBio CLR reads errors are randomly distributed. Using smaller kmer is a great idea. I find the BBsketch is the fastest comparison using kmer analysis, although the documentation is embedded deep in a readme file. BBSketch could build a local NCBI nr sketch and compare reads against it in less than an hour. I am comparing the accuracy of BBSketch and Kraken (micro-organisms) now. It might be hard to track down the sources for de novo assembly, within unknown bacteria, parasites, and even plants. Blobology is now Blobtools2, a very nice tool I have been exploring! Thanks! $\endgroup$ Jan 7, 2021 at 19:12
You can use the kraken database. It will classify your reads based on what organism it matches with. Then can write a custom script to extract the reads you are interested in for downstream analysis.
