Apologies if this question is previously asked. I have a
BAM file alignment to
hg19 and miRNA gff3 file from
mirbase for annotation of BAM file. I want to extract information in the following format (e.g.):
miRNA chr chr_start chr_end expression_count sequence
hsa-mir-155 chr2 ... .... 1020 ATGC....
I know that I can use htseq-count for BAM annotation and for expression counts but how can I get their chromosomal location. The reason for doing this specific way is that I want to show them in a CIRCOS plot.
Just counting is not the goal, their chromosomal location also needs to be extracted. So just subread/htseq-count is not sufficient for extracting genomic location information.
Previously, I had the same idea about processing the sam records. If I want to annotate chromosomal information with miRNA using .gff file, I have to check for each line in .sam records whether their chr_start and chr_end lie in the bandwidth of start and end of any miRNA in .gff file. This process can be very heavy and slow because .sam records are of more than 5 GB. That's why I want to know if there is any extremely fast and generalized way to do it.
I want to avoid processing sam records in python because my experience with extremely heavy data like fastq or sam in python is not good. I would rather prefer awk/bash or any low-level language.