Hello Dr. Bates and group, I understand, the attached data file did not accompany my original message. I have listed below the code used to create that file.
data.1 <- data.frame(subject = factor(rep(c("one", "two", "three", "four", "five", "six", "seven", "eight"), each = 4), levels = c("one", "two", "three", "four", "five", "six", "seven", "eight")), day = factor(rep(c("one", "two", "three", "four"), times = 8), levels = c("one", "two", "three", "four")), expt = rep(c("control", "treatment"), each = 16), response = c(58, 63, 57, 54, 63, 59, 61, 53, 52, 62, 46, 55, 59, 63, 58, 59, 62, 59, 64, 53, 63, 75, 62, 64, 53, 58, 62, 53, 64, 72, 65, 74)) mtrx.1 <- matrix(apply(data.1[, -4], 2, function(x) rep(x, 100 - data.1$response)), ncol = 3, byrow = F) mtrx.2 <- matrix(apply(data.1[, -4], 2, function(x) rep(x, data.1$response)), ncol = 3, byrow = F) data.2 <- data.frame(subject = factor(c(mtrx.1[,1], mtrx.2[,1]), levels = c("one", "two", "three", "four", "five", "six", "seven", "eight")), day = factor(c(mtrx.1[,2], mtrx.2[,2]), levels = c("one", "two", "three", "four")), expt = factor(c(mtrx.1[,3], mtrx.2[,3]), levels = c("control", "treatment")), response = factor(c(rep("yes", nrow(mtrx.1)), rep("no", nrow(mtrx.2))), levels = c("yes", "no"))) #-------------------------------------------------------------------------------# Douglas Bates wrote: >On 9/4/05, Andrew R. Criswell <[EMAIL PROTECTED]> wrote: > >>Hello All, >> >>I have a question regarding how glmmPQL should be specified. Which of >>these two is correct? >> >>summary(fm.3 <- glmmPQL(cbind(response, 100 - response) ~ expt, >> data = data.1, random = ~ 1 | subject, >> family = binomial)) >> >>summary(fm.4 <- glmmPQL(response ~ expt, data = data.2, >> random = ~ 1 | subject, family = binomial)) >> >>One might think it makes no difference, but it does. >> >>I have an experiment in which 8 individuals were subjected to two types >>of treatment, 100 times per day for 4 consecutive days. The response >>given is binary--yes or no--for each treatment. >> >>I constructed two types of data sets. On Rfile-01.Rdata (attached here) >>are data frames, data.1 and data.2. The information is identical but the >>data are arranged differently between these two data frames. Data frame, >>data.1, groups frequencies by subject, day and treatment. Data frame, >>data.2, is ungrouped. >> > >I don't think your attached .Rdata file made it through the various >filters on the lists or on my receiving email. Could you send me a >copy in a separate email message? > > >>Consistency of these data frames is substantiated by computing these >>tables: >> >>ftable(xtabs(response ~ expt + subject + day, >> data = data.1)) >>ftable(xtabs(as.numeric(response) - 1 ~ expt + subject + day, >> data = data.2)) >> >>If I ignore the repeated measurement aspect of the data, I get, using >>glm, identical results (but for deviance and df). >> >>summary(fm.1 <- glm(cbind(response, 100 - response) ~ expt, >> data = data.1, family = binomial)) >> >>summary(fm.2 <- glm(response ~ expt, data = data.2, >> family = binomial)) >> >>However, if I estimate these two equations as a mixed model using >>glmPQL, I get completely different results between the two >>specifications, fm.3 and fm.4. Which one is right? The example which >>accompanies help(glmmPQL) suggests fm.4. >> >>summary(fm.3 <- glmmPQL(cbind(response, 100 - response) ~ expt, >> data = data.1, random = ~ 1 | subject, >> family = binomial)) >> >>summary(fm.4 <- glmmPQL(response ~ expt, data = data.2, >> random = ~ 1 | subject, family = binomial)) >> >>Thank you, >>Andrew >> >> >> >> >>______________________________________________ >>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