TSS plot with htseq 0.6.1 gives segmentation fault

Question

[jellevandewege]

I am trying to run htseq with a bam and gtf file to create a TSS plot. After a while got a segmentation fault (core dumped). Im running the code shown in the example. I tired running it on a smaller file and this gave me no error...

bamfile = HTSeq.BAM_Reader( "s_8_1_I12-1251-03_DP-H3K27Ac.chipseq.hg19.srt.bam" )
gtffile = HTSeq.GFF_Reader( "hg19_IlluminaAnnotation_genes.gtf" )
coverage = HTSeq.GenomicArray( "auto", stranded=False, typecode="i" )

for almnt in bamfile:
    if almnt.aligned:
        coverage[ almnt.iv ] += 1

tsspos = set()
for feature in gtffile:
    if feature.type == "start_codon" and feature.:
        tsspos.add( feature.iv.start_d_as_pos )

p = HTSeq.GenomicPosition( "1", 145439814, "+" )

print p in tsspos

halfwinwidth = 3000
window = HTSeq.GenomicInterval( p.chrom, p.pos - halfwinwidth, p.pos + halfwinwidth, "." )

list( coverage[window] )  
wincvg = numpy.fromiter( coverage[window], dtype='i', count=2*halfwinwidth )
print wincvg

pyplot.plot( wincvg )
pyplot.show()

I have no idea what is error means. It is a memory problem? or is there might be something wrong with the file?

[Devon Ryan]

Did it create a core dump when it crashed? Also, is it possible that you ran out of memory (unlikely, but it's possible). Finally, are you able to run the full version line by line to determine at what point it breaks?

Answer 1

[Ming Tommy Tang]

I have a post using htseq to generate similar figures http://crazyhottommy.blogspot.com/2013/04/how-to-make-tss-plot-using-rna-seq-and.html

or you can use CGAT tool

https://github.com/CGATOxford/cgat

or ngsplot