[R] Fitting GLMM models with glmer

Ben Bolker bbolker at gmail.com
Sat Sep 25 20:51:35 CEST 2010


Clodio Almeida <clodioalm <at> gmail.com> writes:

> 
> Hi everybody:
> 
[snip]

> 
> I obtained the same results for parameters estimates and
> standarderrors, by using:
> 
> glmer(cens ~ treat + heart + offset(log(SURVT) + (1 | INST),
> data=liver, family=poisson)

  It's nice that this check works.  I would like to see a lot
more NLMIXED/glmer comparison ...
 
> After that, the authors present the same model, but now bi = log(gi)
> where gi has the following gamma distribution:
> 
> f(gi | theta1) = gi ^ (1/theta1 - 1) exp(-gi/theta1)/
> GAMMA(1/theta1)(theta1^(1/theta1)
> 
> They used a set of transformations proposed in the paper, and with the
> following rotine fitted the model achieving
> estimates for the same fixed parameters and for theta1:
> 
> proc nlmixed data=liver;
> parms theta1 = 1 b0 = -2:8 btrt = -0.54 bhrt =0:18;
> pi = CDF('NORMAL',ai);
> if(pi > 0:999999) then pi =0:999999;
> gi2 = quantile('GAMMA',pi, 1/theta1);
> gi = theta1 * gi2;
> xb= log_s + b0 + btrt * treat + bhrt * heart+ log(gi);
> lambda=exp(xb);
> model cens ~ poisson(lambda);
> random ai ~ normal(0,1) subject=INST;
> run;
> 
> This time I' m lost. Could anyone please give me a hint?
> 
> The data set is:
> liver <- as.data.frame(matrix(c(1,1,1,1,2,2,2,3,3,3,3,4,4,4,4,
> + 5,5,6,6,7,7,23.286,6.429,26.857,11.143,3.857,9.000,8.714,
> + 1.143,23.143,2.571,2.857,76.429,35.857,25.857,52.286,25.143,
> + 25.143,29.286,28.857,1.857,14.286,1,0,1,0,0,1,0,0,1,1,0,1,
> + 1,0,0,1,0,1,0,0,1,1,0,0,1,1,0,0,1,0,0,1,0,1,0,0,1,0,1,0,1,0),
> + ncol=4,dimnames=list(NULL,c("INST","SURVT","treat","heart"))))

  It's not immediately obvious that you can do this in glmer;
in fact, I suspect you can't. I don't know of an R package that
can fit non-normal random effects 'out of the box'; the only solutions
I know of other than SAS PROC NLMIXED are alternative programs
like WinBUGS and AD Model Builder (which do have R interfaces,
but have their own fairly significant learning curves).

  You might try asking this question on r-sig-mixed-models
instead.



More information about the R-help mailing list