This is a test version of Biostars. For the public version, visit https://www.biostars.org.
read coverage and transcripts

Can somebody suggest me how to obtain read coverage of a particular transcript in a transcript assembly done by TopHat/Cufflinks?

rna-seq sequencing

a transcript assembly done by TopHat/Cufflinks?

Correct me if I'm wrong, but are you sure this is an assembly? TopHat is an aligner, not an assembler.

@WouterDeCoster, TopHat is an aligner and Cufflinks is used for assembly.

Right, but it doesn't modify the alignment, does it? You want to assess the coverage in the alignment. An assembly doesn't have a coverage.

Yes I want to assess the coverage in alignment. Actually in my question, I wanted to make it clear that what software I used

Alright very well. Excuse me for my pedantic nitpicking here, but correct terminology is quite important for quickly getting the right answers.

3 answers

Here, you need to find the number of reads that are mapped to the given transcripts. For this purpose you need mapped BAM file, the transcript chromosome number and it's start-end co-ordinates. If you have this information, you can find read coverage using samtools as follows

samtools view BAM_file Chromosome:start-end

It will give you all reads that mapped to given transcript region. To count the read coverage, you can pipe wc -l command

samtools view BAM_file Chromosome:start-end | wc -l

Renesh, what if a gene has more than one transcript and each transcript has same co-ordinate as in case of exon skipping or intron retention. How one would find the coverage for each different transcript?

Alongside samtools, bedtools offers ways of reading Tophat/cufflinks output too http://bedtools.readthedocs.io/en/latest/index.html

In addition to the tools mentioned earlier, commonly used tools for counting reads in a genomic interval (gene, exon,...) are featureCounts and htseq-counts. I would recommend featureCounts, because it's very fast, has convenient options, but htseq-counts also works fine. Both are nicely documented.

Log in to answer this question.