Differentially expressed genes from RPKM

Question

[Dataminer]

Hi!

So, I have 3 control samples and 13 treatment samples. I have only RPKMs of these samples and not the tag counts or saw access to alignment file. The only information available to me is about the RPKMs of the genes in control and treatment samples.

Is their a way to get a list of genes that are upregulated and down regulated between control and treatment.

I know, if I had information on tag/reads counts application of software like DE-Seq, edgeR etc would have easily solved my problem. However, I did apply DE-Seq on my samples but the results were far from ideal (as I was using RPKMs instead of readcount).

Moreover, visualizing the samples using heatmap along with the application of Kmeans produced some clusters of DE genes. Somehow I am not convinced if this is better way or not.

Any suggestions or guidance will be of great help.

Thank you

Answer 1

[Dataminer]

So, I googled up the same terms with different combinations and found this package in R , sRAP.

It exactly does what I was looking for here is the link http://www.bioconductor.org/packages/release/bioc/html/sRAP.html

It also has a manual explaining everything one would like to know.

Thank you all for your comments and answer.

Answer 2

[Devon Ryan]

You might be able to just log2 transform the values and then use limma. Alternatively, log2 transform them and then just use a T-test, though this will have lower power.

[Dataminer]

Could you please elaborate a bit on your answer with some example or so. I have 3 control samples and 13 treated samples. Even a link to discussion or a small tutorial will be more than fine as how to proceed. I have transformed my data to log2.

[Devon Ryan]

I doubt there will be any example anywhere since this is pretty straight forward. If you have these values in a matrix, then pretty much any limma example will apply (it has a very extensive user's guide with numerous example experiments).

The T-test method would just be using the apply() function on the matrix.

[poisonAlien]

Yeah limma would be perfect. Or if you just want to use t.test, these commands would be enough to begin with.

Assuming x is the data frame (log2 transformed vales) with rownames corresponding to gene names.

x$p.values=apply(x,1,function(x){t.test(x[1:3],x[4:16])$p.value}) #this will add a column named p.value with calculated p.values
x$fdr=p.adjust(x$p.values) #adjust p values for FDR