Entering edit mode
2.5 years ago
MolGeek
▴
80
Hi all.
I want to perform read count from several bams (Alignment was performed using mm10) using the htseq-count.
I am using the Mus_musculus.GRCm38.101.ccds.gtf:
head Mus_musculus.GRCm38.101.ccds.gtf
1 ensembl_havana exon 3214482 3216968 . - . gene_id "ENSMUSG00000051951"; gene_version "5"; transcript_id "ENSMUST00000070533"; transcript_version "4"; exon_number "3"; gene_name "Xkr4"; gene_source "ensembl_havana"; gene_biotype "protein_coding"; havana_gene "OTTMUSG00000026353"; havana_gene_version "2"; transcript_name "Xkr4-201"; transcript_source "ensembl_havana"; transcript_biotype "protein_coding"; tag "CCDS"; ccds_id "CCDS14803"; havana_transcript "OTTMUST00000065166"; havana_transcript_version "1"; exon_id "ENSMUSE00000448840"; exon_version "2"; tag "basic"; transcript_support_level "1":
The command is:
htseq-count -f bam -t exon File1.bam Mus_musculus.GRCm38.101.ccds.gtf > File1.txt
The output is:
ENSMUSG00000118396 0
ENSMUSG00000118409 0
ENSMUSG00000118459 0
ENSMUSG00000118506 0
ENSMUSG00000118560 0
__no_feature 30945311
__ambiguous 0
__too_low_aQual 0
__not_aligned 6999135
__alignment_not_unique 7287645
The code runs without errors but feature counts are 0.
How could i fix it?
Along with -t option you also need to provide -i option where you can mention the attribute for which you need the count (e:g exon_id).
I already tried but it doesn't work.