On Jun 7, 2013, at 10:36 AM, Daofeng Li <lid...@gmail.com> wrote:
> Thank you Sarah and Marc for your fast and nice response.
> Apology for didn't include all information.
>
> I have a input file like following:
>
> gene1 18 15 13 13 16 9 20 24
> gene2 15 8 8 7 0 12 18 4
> gene3 10 9 8 12 9 11 12 12
> gene4 4 0 4 3 0 5 0 0
> gene5 0 1 0 0 1 5 1 0
> gene6 3 3 3 3 4 4 4 4
> gene7 0 4 0 2 2 2 0 0
> gene8 4 4 7 3 0 6 6 12
> gene9 11 6 13 10 13 7 12 9
> gene10 6 3 6 3 4 7 6 3
>
> I am using following R code to compare the difference between the 1st 4
> numbers against 2nd 4 numbers:
>
> library(MASS)
> library(coin)
> suppressPackageStartupMessages(suppressWarnings(library(tcltk)))
> library(qvalue)
> library(plyr)
> dat =
> read.table("test",col.names=c("gene","b1","b2","b3","b4","c1","c2","c3","c4"))
> index=(apply(dat[,-1],1,sum)>0)
> data = dat[index,]
> group=c(1,1,1,1,0,0,0,0)
> n=nrow(data)
> result=NULL
> for (i in 1:n){
> print(i)
> y=as.numeric(data[i,-1])
> if (all((y-mean(y))==0))
> result=rbind(result,c(0,0,0,1))
> else {
> fit=glm.nb(y~group)
> result=rbind(result,summary(fit)$coef[2,])
> }
> }
> qval = qvalue(result[,4])
> fdr=0.1
> index=(qval$qvalues<fdr)
> dat.result = data[index,]
> write.table(dat.result,file="test_result",row.names=F,quote=F)
>
> If you use this input file and code, would reproduce the same error:
>
> Loading required package: methods
> Loading required package: survival
> Loading required package: splines
> Loading required package: mvtnorm
> Loading required package: modeltools
> Loading required package: stats4
>
> Attaching package: ‘plyr’
>
> The following object is masked from ‘package:modeltools’:
>
> empty
>
> [1] 1
> [1] 2
> [1] 3
> [1] 4
> [1] 5
> [1] 6
> Error in while ((it <- it + 1) < limit && abs(del) > eps) { :
> missing value where TRUE/FALSE needed
> Calls: glm.nb -> as.vector -> theta.ml
> In addition: Warning messages:
> 1: In theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace = control$trace
> > :
> iteration limit reached
> 2: In theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace = control$trace
> > :
> iteration limit reached
> Execution halted
>
> So might be the error was in 6th line, not the line I pasted before (5th
> line)? Sorry about that.
>
> Thanks.
>
> Daofeng
Your 'y' at that point in the loop is:
> y
[1] 3 3 3 3 4 4 4 4
and 'group' is:
> group
[1] 1 1 1 1 0 0 0 0
> glm.nb(y ~ group)
Error in while ((it <- it + 1) < limit && abs(del) > eps) { :
missing value where TRUE/FALSE needed
Think about it... :-)
Regards,
Marc Schwartz
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
