meaning of every column of accepted_hit.bam from tophat

[chenyao]

Dear all,

I got the results accepted_hit.bam from tophat:

SRR039211RIP-seq2.1012865 0 chr1 3040608 255 36M * 0 0 TCTTCTTCTTCTTCTCTATCTCCCTGTAACCCCCTC IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII NM:i:0 NH:i:1
SRR039211RIP-seq2.1976722 0 chr1 3043067 255 36M * 0 0 GCCCCCCATTTGGTTATTAGGAGGAGGTCCCCCCAG IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII NM:i:2 NH:i:1

----
what's the meaning of every colum

[iansealy]

http://samtools.sourceforge.net/SAM1.pdf should help.

Cheers,
Ian

[chenyao]

http://samtools.sourceforge.net/SAM1.pdf should help.

Cheers,
Ian

Thank you. Another question, I don't find the count of read number in this file.

[thurisaz]

Each line represents an aligned read. Samtools can provide a summary of the file, but if you've allowed multiple alignments in tophat that will be a count of the number of alignments, rather than the number of unique reads aligned. To get the latter, you can use unix tools like awk, sort, uniq, wc...

[DZhang]

samtools view provides a read count function and you can specify MAPQ value. I believe the options are -c (count the reads) and -q INT (MAPQ value, choose 1 for "unique" reads.

[chenyao]

samtools view provides a read count function and you can specify MAPQ value. I believe the options are -c (count the reads) and -q INT (MAPQ value, choose 1 for "unique" reads.

Dear Zhang,

I think the samtools view -c is cout the reads, for example, if an alignment may be given multiple times if it is overlapping several regions. But I want to now the reads number or you can call it the copy number of reads, and I don't see this in the output from Tophat.

[DZhang]

Hi Chenyao,

Do you want to know how many times a particular read sequence is present in the data? If that's the case, you may use FastQC or any data QC program to get the information. I believe FastQC does not sample every reads but a small portion of them to give you a sense of duplication. If you know perl, it should take a few minutes to write up a script to count that.

[chenyao]

Hi Chenyao,

Do you want to know how many times a particular read sequence is present in the data? If that's the case, you may use FastQC or any data QC program to get the information. I believe FastQC does not sample every reads but a small portion of them to give you a sense of duplication. If you know perl, it should take a few minutes to write up a script to count that.

Yes, that's what I mean. I am confused, if the output of tophat didn't contain this information, how can we use the accept_hit.bam to calculate the RPKM or FPKM value.
I may get wrong about the RNA-seq. I thought I didn't want to find duplication. I want to know the number of some transcript. So if some transcript were transcripted many times, which I think is quite often, then a particular read from this transcript should have large numbers. Suppose read A was only mapped to gene B and gene B only has this only read. If there are 200 read A, so roughly, not normalized, the expression of gene B is 200. That's what we get expression value from RNA-seq, isn't it?