count reads or count base-pairs

yuelics · 07-26-2011, 06:48 PM

Hi all,

I am new to the field. Here is a silly question:

For quantification of defined genomic regions (let it be genes, windows, etc) based on sequencing data, why do people rarely count base-pairs (a convenient function for doing that is "coverage" in R) but rather count reads?

Thanks in advance,
Yue

westerman · 07-28-2011, 07:52 AM

By "quantification" I presume you mean something like a transcriptome project where you want to find out the expression level of genes. For such a project every biological transcript should equal to exactly one read. E.g., if a cell is producing 50 mRNAs for protein 'A' and 10 mRNAs for protein 'B' then there will be a total of 60 reads generated, 50 of which will map to gene 'A' and 10 to gene 'B'. Depending on the sequencing technology, the reads may vary in length. Therefore the number of reads is a better count of the transcriptomes than the base pair coverage.

For other projects you may want (or at least get by) with a base-pair count. For example when discovering SNPs via whole-genome-mapping or exon-mapping then, since they are (duh!) single-base the number of reads will equal the number of base-pairs. So either number could suffice. However for SNP quality control, it is useful to know the number of *unique* reads covering the SNP thus one might as well use the read count instead of the base-pair coverage count. (BTW: the number of unique reads should roughly equal the number of reads otherwise you potentially have problems.)

For other projects, such as general mapping or de-novo assembly then the coverage via base-pairs is an important statistic.

Hope that the above helps.

yuelics · 07-29-2011, 05:50 AM

Thanks westerman. Based on your analysis, if we know the transcripts along the genome, then it is better to use read coverage.

In de novo approach, where a sliding window of fixed size is used, should the base coverage be equivalent to window size of 1?

I have never seen anyone set this window size to 1 (25 bp at least in some literature, but usually 1kb). But using window size of 1 eliminates the concern for reads spanning across two windows. But the complication may be that we introduce highly dependent event for bases within a read, which might make it hard to model with simple Markov chain of limited order (e.g., HMM of order 1).

If using read coverage, what would be an effective window width, and should the adjacent windows of large size (e.g., 1kb) be overlapping each other to some extent (e.g., 500 bp)?

Thanks for the helps!

westerman · 07-29-2011, 06:41 AM

If you are talking about de-novo transcriptomes (which I do a lot of, mainly from the 454 platform) here too I would use read counts. Where appropriate I bundle all of the reads together to create the de-novo transcripts. I then figure how how many of each sample's reads map to the transcripts. This is the count that we use. In many ways this work is no different than when you actually know the transcripts ahead of time ... except for de-novo you have no idea (without extra work) of what the transcripts actually represent or if they actually exist in nature.

I can not see using a sliding window on transcripts. Except, perhaps, to determine if the transcripts are likely to exist. In other words transcripts that have a large number of bases supporting them are much more likely to be correct than those without a number of bases. The sliding window count of bases would even out disparities in read coverage and could allow one to differentiate between false transcripts versus those just expressed at a very low level.

But then you may be doing a more in-depth analysis (or an entirely different analysis) than what I normally do. Basically what we give our customers is a list of their generated transcripts (plus blast/GO terms) with the associated read counts for each of their samples. It is up to them to figure out if those transcripts and counts signify anything biologically important or how to model them further.

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
Quantification: count reads or count base pairs?	yuelics	Bioinformatics	0	07-27-2011 05:48 AM
count reads across exon junctions	suninsky	Bioinformatics	4	02-23-2011 08:34 AM
looking for solution to count the short reads in rice	dingkai0564	RNA Sequencing	0	11-04-2010 06:19 PM
cufflinks exonic reads count	repinementer	Bioinformatics	0	08-13-2010 08:25 AM

07-26-2011, 06:48 PM	#1
yuelics Member Location: Toronto, Canada Join Date: Apr 2011 Posts: 13	count reads or count base-pairs Hi all, I am new to the field. Here is a silly question: For quantification of defined genomic regions (let it be genes, windows, etc) based on sequencing data, why do people rarely count base-pairs (a convenient function for doing that is "coverage" in R) but rather count reads? Thanks in advance, Yue

07-29-2011, 05:50 AM	#3
yuelics Member Location: Toronto, Canada Join Date: Apr 2011 Posts: 13	Thanks westerman. Based on your analysis, if we know the transcripts along the genome, then it is better to use read coverage. In de novo approach, where a sliding window of fixed size is used, should the base coverage be equivalent to window size of 1? I have never seen anyone set this window size to 1 (25 bp at least in some literature, but usually 1kb). But using window size of 1 eliminates the concern for reads spanning across two windows. But the complication may be that we introduce highly dependent event for bases within a read, which might make it hard to model with simple Markov chain of limited order (e.g., HMM of order 1). If using read coverage, what would be an effective window width, and should the adjacent windows of large size (e.g., 1kb) be overlapping each other to some extent (e.g., 500 bp)? Thanks for the helps! Last edited by yuelics; 07-29-2011 at 05:58 AM.

07-28-2011, 07:52 AM	#2
westerman Rick Westerman Location: Purdue University, Indiana, USA Join Date: Jun 2008 Posts: 1,104	By "quantification" I presume you mean something like a transcriptome project where you want to find out the expression level of genes. For such a project every biological transcript should equal to exactly one read. E.g., if a cell is producing 50 mRNAs for protein 'A' and 10 mRNAs for protein 'B' then there will be a total of 60 reads generated, 50 of which will map to gene 'A' and 10 to gene 'B'. Depending on the sequencing technology, the reads may vary in length. Therefore the number of reads is a better count of the transcriptomes than the base pair coverage. For other projects you may want (or at least get by) with a base-pair count. For example when discovering SNPs via whole-genome-mapping or exon-mapping then, since they are (duh!) single-base the number of reads will equal the number of base-pairs. So either number could suffice. However for SNP quality control, it is useful to know the number of unique reads covering the SNP thus one might as well use the read count instead of the base-pair coverage count. (BTW: the number of unique reads should roughly equal the number of reads otherwise you potentially have problems.) For other projects, such as general mapping or de-novo assembly then the coverage via base-pairs is an important statistic. Hope that the above helps.

07-29-2011, 06:41 AM	#4
westerman Rick Westerman Location: Purdue University, Indiana, USA Join Date: Jun 2008 Posts: 1,104	If you are talking about de-novo transcriptomes (which I do a lot of, mainly from the 454 platform) here too I would use read counts. Where appropriate I bundle all of the reads together to create the de-novo transcripts. I then figure how how many of each sample's reads map to the transcripts. This is the count that we use. In many ways this work is no different than when you actually know the transcripts ahead of time ... except for de-novo you have no idea (without extra work) of what the transcripts actually represent or if they actually exist in nature. I can not see using a sliding window on transcripts. Except, perhaps, to determine if the transcripts are likely to exist. In other words transcripts that have a large number of bases supporting them are much more likely to be correct than those without a number of bases. The sliding window count of bases would even out disparities in read coverage and could allow one to differentiate between false transcripts versus those just expressed at a very low level. But then you may be doing a more in-depth analysis (or an entirely different analysis) than what I normally do. Basically what we give our customers is a list of their generated transcripts (plus blast/GO terms) with the associated read counts for each of their samples. It is up to them to figure out if those transcripts and counts signify anything biologically important or how to model them further.