How to run the rpkm function of edgeR
1
0
Entering edit mode
Gary ▴ 20
@gary-7967
Last seen 5.7 years ago

Hi,

I would like to use edgeR to estimate the RPKM values.  I have (1) read counts files estimated by HTSeq-count, and (2) a transcript length file. I know how to estimate CPM in edgeR, using below command lines. However, I don’t know how to estimate RPKM values based on the files I have. Could you show me how to run the command lines in edgeR? Many thanks.

Best,

Gary

 

Transcript count files

http://68.181.92.180/~Gary/temporal/T8_Early_R_T.txt

http://68.181.92.180/~Gary/temporal/T8_Late_R_T.txt

http://68.181.92.180/~Gary/temporal/T9_Early_R_T.txt

http://68.181.92.180/~Gary/temporal/T9_Late_R_T.txt

Transcript length file

http://68.181.92.180/~Gary/temporal/hg19refGene20160206_TranscriptLength.txt

 

> targets=readTargets()

> targets

             files group  description

1 T8_Early_R_T.txt Early T8_Early_R_T

2  T8_Late_R_T.txt  Late  T8_Late_R_T

3 T9_Early_R_T.txt Early T9_Early_R_T

4  T9_Late_R_T.txt  Late  T9_Late_R_T

> d=readDGE(targets,skip=0,comment.char="!")

> keep=rowSums(cpm(d)>1)>=2

> d=d[keep,,keep.lib.sizes=FALSE]

> d=calcNormFactors(d,method="TMM")

> dim(d)

[1] 8202    4

> detags=rownames(d)

> detagsfile=cpm(d)[detags,]

> write.table(detagsfile,file="CPM_MSC_8202ene.txt")

> rpkm<-rpkm(d)

Error in rpkm.DGEList(d) : Gene lengths not found

> rpkm<-rpkm(d,gene.lengths)

Error in rpkm.DGEList(d, gene.lengths) : object 'gene.lengths' not found

edger rpkm • 5.6k views
ADD COMMENT
4
Entering edit mode
Aaron Lun ★ 28k
@alun
Last seen 21 hours ago
The city by the bay

Well, it would be a good start to load the gene lengths into memory.

gene.data <- read.table("hg19refGene20160206_TranscriptLength.txt", header=TRUE)

You then need to match the gene names in d to the ordering of gene.data:

m <- match(rownames(d), gene.data$Transcript)
gene.lengths <- gene.data$TranscriptLength[m]

... and then you should be able to run rpkm on that.

ADD COMMENT
0
Entering edit mode

Hi Aaron,

Thank you so much.

ADD REPLY

Login before adding your answer.

Traffic: 680 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6