I am working on a pooled bacterial sequence reads. My final goal is to create a protein database. The reads have already been mapped to the reference genome. I can see the aligned, indexed and sorted reads and the respective reference genome with the annotation(genes etc) on IGV. I just don't know how to retrieve the annotation information (not just visualize it)of the mapped reads(transfer the annotation of the Reference genome to my reads) for further analysis. Is there an already existing method or tool for achieving this. All the annotation tools I have come across annotate assembled genomes not reads. It seems an easy task but since I am new to the field, I haven't been successful in my search for a solution. If I need to write a script to do it, I can; I just thought I should check if there is an already existing method/tool before I do that. I would appreciate any help. Thank you.
•
written
by
eva •
It's not entirely clear what you mean by "annotating reads", but could you use the reference annotation to identify the regions corresponding to the features you are interested in (i.e. genes), and then use 'samtools view' (which takes an optional argument corresponding to region) to extract the set of reads mapping to each such region?
I agree its not clear what you are trying to do, but bedtools intersect might be what you are looking for:
http://bedtools.readthedocs.org/en/latest/content/tools/intersect.html
Also check this:
http://gettinggeneticsdone.blogspot.com/2014/03/visualize-coverage-exome-targeted-ngs-bedtools.html
Log in to answer this question.
