[BioC] Error in estimateCommonDisp

Gordon K Smyth smyth at wehi.EDU.AU
Fri Oct 18 03:23:17 CEST 2013


Dear Xiaowan Lu,

The estimateCommonDisp() function in edgeR is used a great deal by us and 
others without giving an error, so it is likely that there is some special 
feature of your data that is causing a problem.  Although you have given 
some code in your email, it doesn't provide us any way to reproduce the 
error you have or to track down the cause.

First, the error message that you give could not have arisen from the code 
given in your email.  The error message has arisen from the 
edgeRFindChangedGenes() function in the GROseq package, which is commented 
out in the code you give.

It appears that you are not using edgeR directly but are instead using the 
GROseq package.  GROseq is not available from either Bioconductor or CRAN. 
If you get errors from GROseq functions, you should write directly to the 
authors of that package.

If you want to pursue this further you could email your data and edgeR 
code offline to me or to one of the other edgeR authors.  However you need 
to provide us with a data example that we can run ourselves using edgeR 
commands.  We cannot debug functions in the GROseq package.

Best wishes
Gordon

---------------------------------------------
Professor Gordon K Smyth,
Bioinformatics Division,
Walter and Eliza Hall Institute of Medical Research,
1G Royal Parade, Parkville, Vic 3052, Australia.
http://www.statsci.org/smyth


> Date: Fri, 11 Oct 2013 11:28:52 -0400
> From: xiaowan.lu at uhnresearch.ca
> To: <bioconductor at r-project.org>
> Subject: [BioC] Error in estimateCommonDisp
> Message-ID: <FNVshYDb.1381505332.8512780.xiaowan.lu at uhnresearch.ca>
>
> To whom it may concern,
>
> I am trying to use the function estimateCommonDisp(object) in edgeR package to run some analysis, but I keep getting the error:
>
> Error in if (any(input.mean < 0)) stop("input.mean must be non-negative") :
>  missing value where TRUE/FALSE needed
> Calls: edgeRFindChangedGenes ... estimateCommonDisp -> equalizeLibSizes -> q2qnbinom
>
> I double checked the object I used for the function fun and it looked 
> fine for me. I looked at input.mean, it returns whole bunch of NAs, I am 
> not sure what went wrong, please help.
>
> Here I attached the code I used:
>
> require(GROseq)
> require(edgeR)
> source("FitModel.R")
> Gbed <- read.table("FinalTranscripts.bed")
> tID <- paste(Gbed[[1]],Gbed[[2]],Gbed[[6]], sep="")
> Gbed <- data.frame(Chr=as.character(Gbed[[1]]), Start=as.integer(Gbed[[2]]), End=as.integer(Gbed[[3]]), Str=as.character(Gbed[[6]]), RefSeqID=as.character(tID), MGI=as.character(tID))
>
> G <- LimitToXkb(Gbed, SIZE=13000)
>
>
> load("SOAP.RData")
> nrNH00 <- NROW(NH00)
> nrNH10 <- NROW(NH10)
> nrNH40 <- NROW(NH40)
> nrNH160<- NROW(NH160)
> nrLC00 <- NROW(LC00)
> nrLC10 <- NROW(LC10)
> nrLC40 <- NROW(LC40)
> nrLC160<- NROW(LC160)
>
> RefSeqNH00 <- CountReadsInInterval(f= G, p= NH00[,c(1:3,6)])
> RefSeqLC00 <- CountReadsInInterval(f= G, p= LC00[,c(1:3,6)])
> RefSeqNH10 <- CountReadsInInterval(f= G, p= NH10[,c(1:3,6)])
> RefSeqLC10 <- CountReadsInInterval(f= G, p= LC10[,c(1:3,6)])
>
> #ChangedGenes <- edgeRFindChangedGenes(RefSeqNH00, nrNH00, 
> RefSeqLC00,nrLC00, RefSeqNH10,nrNH10, RefSeqLC10,nrLC10, 
> FILENAME="T/T.ModelFit-10m.png", G=Gbed)
>
> nZI <- !(RefSeqNH00 == 0 & RefSeqLC00 == 0 & RefSeqNH10 == 0 & RefSeqLC10 == 0)
> cat("nZI=", nZI, "/n")
>
> d <- list()
> d$counts <- as.matrix(data.frame(NH00= RefSeqNH00[nZI], LC00= RefSeqLC00[nZI], NHE2= RefSeqNH10[nZI], LCE2= RefSeqLC10[nZI]))
> dim(d$counts) <- c(NROW(d$counts), NCOL(d$counts))
>
> colnames(d$counts) <- c("NH00", "LC00", "NHE2", "LCE2")
>
> d$samples$files <- c("NH00", "LC00", "NHE2", "LCE2")
> d$samples$group <- as.factor(c("VEH", "VEH", "E2", "E2"))
> d$samples$description <- c("VEH", "VEH", "E2", "E2")
> d$samples$lib.size <- c(nrNH00, nrLC00, nrNH10, nrLC10)
> d <- new("DGEList",d)
> d <- estimateCommonDisp(d)
>
>
> Xiaowan Lu
> Bioinformatics Research Assistant,
> CO-OP Student
> University Health Network

______________________________________________________________________
The information in this email is confidential and intend...{{dropped:4}}



More information about the Bioconductor mailing list