HTSeq count : outsam

[barak]

Barak,

I think this is a bug in HTSeq when trying to use BAM files directly as input; it has happened to me as well. The only way around it is to fall back to using SAM input. You can at least pipe the SAM stream to htseq-count so you don't have to save the uncompressed SAM on disk.

samtools view -h accepted_hits.bam | htseq-count -f sam --stranded yes -a 10 -t exon -i gene_id -m union -o samout3.sam - genes.3utr.gtf > counts3.txt

Thanks, it works!
While working on this I realized that I could use the method: get_sam_line() and change the code in HTSeq-count. I'm still trying to figure out a reason NOT to do so, but it seems to hold the same info as with the original sam line.
B

[kmcarr]

Hi
I'm running HTSeq-count version 0.6.1p1. Trying out outsam option: I do get the sam file with the XF tags but the reads are blanks.
example:
XF:Z:__alignment_not_unique
XF:Z:__no_feature
XF:Z:__no_feature
...

Also I tried with the python framework, imported the lib, read the bam file, everything looks OK, but the attribute: original_sam_line is indeed a blank string.

Any ideas?

The input bam file is from tophat (v2.0.10)


This is my command:
htseq-count -f bam --stranded yes -a 10 -t exon -i gene_id -m union -o samout3.sam accepted_hits.bam genes.3utr.gtf > counts3.txt

thanks,
B

Barak,

I think this is a bug in HTSeq when trying to use BAM files directly as input; it has happened to me as well. The only way around it is to fall back to using SAM input. You can at least pipe the SAM stream to htseq-count so you don't have to save the uncompressed SAM on disk.

samtools view -h accepted_hits.bam | htseq-count -f sam --stranded yes -a 10 -t exon -i gene_id -m union -o samout3.sam - genes.3utr.gtf > counts3.txt

[barak]

Hi
I'm running HTSeq-count version 0.6.1p1. Trying out outsam option: I do get the sam file with the XF tags but the reads are blanks.
example:
XF:Z:__alignment_not_unique
XF:Z:__no_feature
XF:Z:__no_feature
...

Also I tried with the python framework, imported the lib, read the bam file, everything looks OK, but the attribute: original_sam_line is indeed a blank string.

Any ideas?

The input bam file is from tophat (v2.0.10)


This is my command:
htseq-count -f bam --stranded yes -a 10 -t exon -i gene_id -m union -o samout3.sam accepted_hits.bam genes.3utr.gtf > counts3.txt

thanks,
B