Hi all,
Yesterday I asked a question about IGV but I realized that I asked from the wrong side. So I am asking again by revising my question.
I just only want to find the count of mapped alignments shown in the red box for the egfr gene.
And to do that; I created a bed file that only span EGFR gene then I run
bedtools multicov bedtools multicov -S -bams HCC827-1_Aligned_sorted.out.bam -bed ~/refgen/chr7EGFR_downloaded.bed
command above but my output was like this;
NC_000007.13 55086709 55279321 EGFR 0 + 0
NC_000007.13 55247442 55256642 EGFR-AS1 0 - 0
I know there are mapped alignments to that region but the code does not give anything. If someone could help me I would be appreciated. Thanks
•
written
by
biology_inform •
Are you sure about the -S option ? we can not appreciate strandness from the igv snapshot.
Are the chrosome names exactly the same between the .bed and the .bam files ? You can quickly verify that with the commands below:
samtools idxstats HCC827-1_Aligned_sorted.out.bam | cut -f 1 # Gives chrom names in .bam files
cut -f 1 ~/refgen/chr7EGFR_downloaded.bed | sort | uniq # Gives chrom names in .bed file
Also note that bedtools multicov filter out low quality reads and multimaping reads. You can suppress this behaviour by using the -F ad -D options, respectively. It is not recommended to do so, but it might explain the discrepancy between igv and bedtools.
Finally, I suppose that the repetition of bedtools multicov (you wrote "bedtools multicov bedtools multicov") in your command is a typo.
It is a copy paste mistake :) Okay I am going to check with the commands. But if they are not same how can I do that ?
@Carlo Yague they are not the same :) in bed file NC_000007.13 like this but bam file it is chr7. But I downloaded my bed file from Genome data viewer NCBI (https://www.ncbi.nlm.nih.gov/genome/gdv/browser/gene) should I choose different track?
You should always try to download the bed annotation from the same source as the genome used for the alignment.
Log in to answer this question.