20 millions reads per sample for SNP

Question

[Omics data mining]

Dear all,

I have RNASEQ sequencing data (tumour) in which each sample has approximately 20 million reads. My target is to identify SNPs and predict the isoform expression variability using RNAseq data. Before doing all analysis , I wanted to check coverage of sequence data. Just to check coverage , I used formaula N*L/G where N no of reads per sample, L read length (150 as paired end data 75+75), and human genome size. For sample, with 26130832 reads count, I got coverage 1.3. Does this 1.3 indicates sequence coverage is 1X ? Am I right ? My target is to identity SNPs. Can I proceed with this much amount of data per sample?

Thank you in advance

Answer 1

[swbarnes2]

For what it's worth, in my lab, we shoot for 30M reads just to assess expression differences; you'd want more for variants. For variant calling, the average coverage is not helpful. You need a finer breakdown...you need to know what percentage of the transcriptome is adequately covered (say, 25x). You will need something like BEDTools to assess that. My guess is, you might have enough depth to examine the most highly expressed genes, but not much else.

Answer 2

[Titus]

Hi ,

You should use G = EXOM , you count non coding part in your formula :)

Best

[Omics data mining]

Dear Titus,

Yes. You are right. G=EXOM part will represent the actual targeted genome size . Thanks :)

Based on this, I think my data is good enough to proceed for analysis.

[WouterDeCoster]

RNA-seq has a very uneven coverage pattern. Some highly abundant genes will have high coverage (>1000 reads) while you'll also have genes with some expression (~40-80 reads), low abundant genes with just a few reads (<5) and many genes without expression in your tissue of interest (0 reads). As such, simply taking an average of the read count is pointless and meaningless.

I'm not saying that you can't, but RNA-seq is far from optimal for variant calling.