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.
______________________________________________
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.
