[R] Glmm for multiple outcomes

Spencer Graves spencer.graves at pdf.com
Sun Oct 9 18:05:03 CEST 2005


Does the following help:

n.subjects <- 3
J <- 4
K <- 5
n.ijk <- rep(2, each=n.subjects*J*K)

x <- rep(1:K, n.subjects, each=J)

subj <- factor(rep(1:n.subjects, each=K*J))

sa.subject <- 1
sb.subject <- 1

set.seed(2)
a.subj <- rep(sa.subject*rnorm(n.subjects), each=K*J)
b.subj <- rep(sb.subject*rnorm(n.subjects), each=K*J)

Z <- a.subj+b.subj*x

library(boot)
Y <- (rbinom(n.subjects*K*J, n.ijk, inv.logit(Z))
          /n.ijk)

Dat <- data.frame(subj=subj, x=x, y=Y)

library(lme4)

fit <- lmer(y~x+(x|subj), Dat)
Linear mixed-effects model fit by REML
Formula: y ~ x + (x | subj)
    Data: Dat
       AIC      BIC    logLik MLdeviance REMLdeviance
  51.63172 64.19779 -19.81586    33.1066     39.63172
Random effects:
  Groups   Name        Variance  Std.Dev. Corr
  subj     (Intercept) 0.0446346 0.211269
           x           0.0032613 0.057108 1.000
  Residual             0.0879438 0.296553
# of obs: 60, groups: subj, 3

Fixed effects:
              Estimate Std. Error DF t value Pr(>|t|)
(Intercept)  0.350000   0.151459 58  2.3109  0.02442 *
x            0.033333   0.042661 58  0.7814  0.43777
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

	  spencer graves

Abderrahim Oulhaj wrote:

> Dear All,
> 
> I wonder if there is an efficient way to fit the  generalized linear mixed model  for multivariate outcomes.
> 
> More specifically, Suppose that for a given subject i and at a  given time j we observe a multivariate  outcome Yij = (Y_ij1, Y_ij2, ..., Y_ijK). 
>  where Y_ijk is a binomial(n_ijk, p_ijk). 
> 
> One way to jointly model  the data is to use the following specification:
> 
> g(p_ijk) = beta_0k + b_0ik + (beta_1k + b_1ik)*x_ijk  with k = 1,2 ...., K , g is a specified link function and (b_0ik,b_1ik) k=1,...K are random effects ...
> 
>  I my case, the glmmPQL converges only  and give good results when k is less than 3 (i.e. for a small number of random effects). I also used the gee (generalized estimating equations) to estimate the fixed effects and the same probleme ariseed with k.
> 
> Is there any help?
> 
> Thank you in advance, 
>  
> Abderrahim Oulhaj
> 	[[alternative HTML version deleted]]
> 
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

-- 
Spencer Graves, PhD
Senior Development Engineer
PDF Solutions, Inc.
333 West San Carlos Street Suite 700
San Jose, CA 95110, USA

spencer.graves at pdf.com
www.pdf.com <http://www.pdf.com>
Tel:  408-938-4420
Fax: 408-280-7915




More information about the R-help mailing list