Elisabeth,
Thank you for the code. I will give it a shot and see how it goes.
I have also seen that glmmPQL {MASS} can use the quasipoisson family
(which is usually fairly close to the negative binomial (?)). So this
looks like an "easier" way to accomplish what I want. [I would like
to thank Ben Bolker for his suggestion to use the quasipoisson. As a
matter of fact, I have used it for glm models in the past (along with
negative binomial) and gave similar results. Thank you also for the
useful reference.]
Thank you also to Prof. Ripley for his response. In the past I had
used glm.nb, but I have to admit I had not realized that there is
also a negative.binomial family (to use as an argument for "family")!
Once again, thank you for the help and the wonderful book (it adds
some weight to my backpack, but it is worth it!).
In a somewhat related comment, I have stopped receiving emails from R-
help this year, i.e. after Dec 31, 2005 (for some unknown reason... I
subscribed with a different address, so I am getting the messages
again today, but I cannot reply to the two emails that Ben Bolker and
Prof. Ripley sent directly. I am replying above, as I only saw them
in the archives).
Best regards,
Costas
On Jan 3, 2006, at 11:54 PM, Elizabeth Lawson wrote:
> Have you tried nlme?
>
> I tried to something similar.
>
> Here is the code that I used for a negative binomial random effects
> model
>
>
> library(nlme)
> mydata<-read.table("C:\\Plx\\plx.all\\plxall.txt",header=TRUE)
> loglike = function
> (PLX_NRX,PD4_42D,GAT_34D,VIS_42D,MSL_42D,SPE_ROL,XM2_DUM,THX_DUM,b0,b1
> ,b2,b3,b4,b5,b6,b7,alpha){
> lambda = exp(b0 + b1*GAT_34D+b2*VIS_42D+b3*MSL_42D+b4*PD4_42D
> +b5*SPE_ROL+b6*XM2_DUM+b7*THX_DUM)
> y=round(PLX_NRX)
> y <- table(y)
> freq <- as.vector(y)
> count <- as.numeric(names(y))
> count <- count[!(freq < 1)]
> freq <- freq[!(freq < 1)]
> n <- length(count)
> df <- n - 1
> &nbs! p;
>
>
>
> df <- df - 2
> xbar <- weighted.mean(count, freq)
> s2 <- var(rep(count, freq))
> p <- xbar/s2
> alpha <- xbar^2/(s2 - xbar)
>
> ( dnbinom(y,alpha,(alpha/(alpha+lambda))) )
> }
>
>
> plx.nlme<-nlme(PLX_NRX~loglike
> (PLX_NRX,PD4_42D,GAT_34D,VIS_42D,MSL_42D,SPE_ROL,XM2_DUM,THX_DUM,b0,b1
> ,b2,b3,b4,b5,b6,b7,alpha),
> data=mydata,
> fixed=list(b0 + b1+b2+b3+b4+b5+b6+b7+alpha~1),
> random=b0~1|menum,
>
> start=c(b0=0,b1=0,b2=0,b3=0,b4=0,b5=0,b6=0,b7=0,alpha=5)
> )
>
> I am not sure that this is what you are looking for , but I hope
> this helps!
>
> Elizabeth Lawson
> Constantinos Antoniou <[EMAIL PROTECTED]> wrote:
> Hello all,
>
> I would like to fit a mixed effects model, but my re! sponse is of the
> negative binomial (or overdispersed poisson) family. The only (?)
> package that looks like it can do this is glmm.ADMB (but it cannot
> run on Mac OS X - please correct me if I am wrong!) [1]
>
> I think that glmmML {glmmML}, lmer {Matrix}, and glmmPQL {MASS} do
> not provide this "family" (i.e. nbinom, or overdispersed poisson). Is
> there any other package that offers this functionality?
>
> Thanking you in advance,
>
> Costas
>
>
> [1] Yes, I know I can use this on another OS. But it is kind of a
> nuisance, as I have my whole workflow setup on a mac, including emacs
> +ess, the data etc etc. It will be non-trivial to start moving/
> syncing files between >1 computers, in order to use this package...
>
> --
> Constantinos Antoniou, Ph.D.
> Department of Transportation Planning and Engineering
> National Technical University of Athens
> 5, Iroon Polytechniou str. GR-15773, Athens, Greece
>
> ______________________________________________
> [email protected] mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-
> guide.html
>
>
> Yahoo! for Good - Make a difference this year.
--
Constantinos Antoniou, Ph.D.
Department of Transportation Planning and Engineering
National Technical University of Athens
5, Iroon Polytechniou str. GR-15773, Athens, Greece
[[alternative HTML version deleted]]
______________________________________________
[email protected] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
