From: Spencer Graves <spencer.graves_at_pdf.com>

Date: Sat 11 Feb 2006 - 15:05:31 EST

iteration 2

iteration 3

iteration 4

iteration 5

iteration 6

Linear mixed-effects model fit by maximum likelihood

> (denDF.week <- anova(fit)$denDF[1])

[1] 166

> (denDF.week <- anova(fit)$denDF[1])

[1] 166

> (par.week <- fixef(fit)[1:5])

factor(week)0 factor(week)2 factor(week)4 factor(week)6 factor(week)11

* > CM <- array(0, dim=c(5*4/2, 5))
*

> i1 <- 0

> for(i in 1:4)for(j in (i+1):5){

[8,] 0 0 -1 1 0

[9,] 0 0 -1 0 1

[10,] 0 0 0 -1 1

> library(multcomp)

> csimint(par.week, df=denDF.week, covm=vc.week,cmatrix=CM)

[10,] 0.001 -1.206 1.209

[10,] 0 0 0 -1 1

[10,] 0.997

> sessionInfo()

R version 2.2.1, 2005-12-20, i386-pc-mingw32

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Sat Feb 11 15:08:26 2006

Date: Sat 11 Feb 2006 - 15:05:31 EST

The following appears to be an answer to your question, though I'd be pleased to receive critiques from others. Since your example is NOT self contained, I modified an example in the "glmmPQL" help file:

(fit <- glmmPQL(y ~ factor(week)-1+trt, random = ~ 1 | ID,

+ family = binomial, data = bacteria))iteration 1

iteration 2

iteration 3

iteration 4

iteration 5

iteration 6

Linear mixed-effects model fit by maximum likelihood

Data: bacteria

Log-likelihood: -551.1184

Fixed: y ~ factor(week) - 1 + trt

factor(week)0 factor(week)2 factor(week)4 factor(week)6 factor(week)11

3.3459650 3.5262521 1.9102037 1.7645881 1.7660845 trtdrug trtdrug+ -1.2527642 -0.7570441

Random effects:

Formula: ~1 | ID

(Intercept) Residual

StdDev: 1.426534 0.7747477

Variance function:

Structure: fixed weights

Formula: ~invwt

Number of Observations: 220

Number of Groups: 50

> anova(fit)

numDF denDF F-value p-value factor(week) 5 166 10.821682 <.0001 trt 2 48 1.889473 0.1622

> (denDF.week <- anova(fit)$denDF[1])

[1] 166

> (denDF.week <- anova(fit)$denDF[1])

[1] 166

> (par.week <- fixef(fit)[1:5])

factor(week)0 factor(week)2 factor(week)4 factor(week)6 factor(week)11

3.345965 3.526252 1.910204 1.764588 1.766085

> (vc.week <- vcov(fit)[1:5, 1:5])

factor(week)0 factor(week)2 factor(week)4 factor(week)6 factor(week)0 0.3351649 0.1799365 0.1705898 0.1694884 factor(week)2 0.1799365 0.3709887 0.1683038 0.1684096 factor(week)4 0.1705898 0.1683038 0.2655072 0.1655673 factor(week)6 0.1694884 0.1684096 0.1655673 0.2674647 factor(week)11 0.1668450 0.1665177 0.1616748 0.1638169 factor(week)11 factor(week)0 0.1668450 factor(week)2 0.1665177 factor(week)4 0.1616748 factor(week)6 0.1638169 factor(week)11 0.2525962

> i1 <- 0

> for(i in 1:4)for(j in (i+1):5){

+ i1 <- i1+1 + CM[i1, c(i, j)] <- c(-1, 1) + }

> CM

[,1] [,2] [,3] [,4] [,5] [1,] -1 1 0 0 0 [2,] -1 0 1 0 0 [3,] -1 0 0 1 0 [4,] -1 0 0 0 1 [5,] 0 -1 1 0 0 [6,] 0 -1 0 1 0 [7,] 0 -1 0 0 1

[8,] 0 0 -1 1 0

[9,] 0 0 -1 0 1

[10,] 0 0 0 -1 1

> library(multcomp)

> csimint(par.week, df=denDF.week, covm=vc.week,cmatrix=CM)

Simultaneous confidence intervals: user-defined contrasts

95 % confidence intervals

Estimate 2.5 % 97.5 % [1,] 0.180 -1.439 1.800 [2,] -1.436 -2.838 -0.034 [3,] -1.581 -2.995 -0.168 [4,] -1.580 -2.967 -0.193 [5,] -1.616 -3.123 -0.109 [6,] -1.762 -3.273 -0.250 [7,] -1.760 -3.244 -0.277 [8,] -0.146 -1.382 1.091 [9,] -0.144 -1.359 1.070

[10,] 0.001 -1.206 1.209

> csimtest(par.week, df=denDF.week, covm=vc.week,cmatrix=CM)

Simultaneous tests: user-defined contrasts

Contrast matrix:

[,1] [,2] [,3] [,4] [,5]

[1,] -1 1 0 0 0 [2,] -1 0 1 0 0 [3,] -1 0 0 1 0 [4,] -1 0 0 0 1 [5,] 0 -1 1 0 0 [6,] 0 -1 0 1 0 [7,] 0 -1 0 0 1 [8,] 0 0 -1 1 0 [9,] 0 0 -1 0 1

[10,] 0 0 0 -1 1

Adjusted P-Values

p adj

[1,] 0.011 [2,] 0.013 [3,] 0.014 [4,] 0.015 [5,] 0.020 [6,] 0.024 [7,] 0.985 [8,] 0.985 [9,] 0.985

[10,] 0.997

> sessionInfo()

R version 2.2.1, 2005-12-20, i386-pc-mingw32

attached base packages:

[1] "methods" "stats" "graphics" "grDevices" "utils" "datasets"
[7] "base"

other attached packages:

multcomp mvtnorm MASS statmod nlme "0.4-8" "0.7-2" "7.2-24" "1.2.4" "3.1-68.1"

If this does NOT answer your question (or even if it does), PLEASE do read the posting guide! "www.R-project.org/posting-guide.html". I'd prefer not to have to guess whether you would think the example I chose was relevant.

hope this helps, spencer graves

Michaël Coeurdassier wrote:

> Dear R community, > > I performed a generalized linear mixed model using glmmPQL (MASS > library) to analyse my data i.e : y is the response with a poisson > distribution, t and Trait are the independent variables which are > continuous and categorical (3 categories C, M and F) respectively, ind > is the random variable. > > mydata<-glmmPQL(y~t+Trait,random=~1|ind,family=poisson,data=tab) > Do you think it is OK? > > Trait is significant (p < 0.0001) and I would like to perform post-hoc > comparisons to check where the difference among Trait categories but > I did not find a solution in R help list or others. > > Thank you in advance for your help > > Michael > > ______________________________________________ > R-help@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 ______________________________________________R-help@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 Received on Sat Feb 11 15:08:26 2006

*
This archive was generated by hypermail 2.1.8
: Fri 03 Mar 2006 - 03:42:27 EST
*