It looks like your example can be found here:
http://books.google.com/books?id=pm-DsOVKY78C&pg=PT166

This page is also very useful:
http://www.ats.ucla.edu/stat/r/seminars/repeated_measures/repeated_measures.htm

Note that "pat" should be a factor.

For the compound symmetry model, the PROC MIXED code has no "random" term,
so instead of using "lme", try using "gls".

In the code below, first I plot the data, then create a groupedData object,
then finally fit the Comp Symm model.

# Plot
alz2 <- unialz
alz2$pat <- factor(alz2$pat)
alz2$month <- factor(alz2$month, levels=c("ADAS02", "ADAS04", "ADAS06",
"ADAS08", "ADAS10", "ADAS12"))
alz2 <- alz2[ order(alz2$month), ]
xyplot(adascog ~ month|trt , data=alz2, groups=pat, type='l',layout=c(3,1),
       scales=list(x=list(rot=90, cex=.75)))

alz2 <- groupedData(adascog ~ as.numeric(month)|pat, data=alz2)
plot(alz2, inner=~trt)

# CS - nlme needs to use "groupedData" with gls
cs2 <- gls(adascog ~ trt * month, data=alz2,
           correlation=corCompSymm(form = ~ 1 | pat),
           na.action=na.omit)
logLik(cs2)*2 # match to SAS
summary(cs2$modelStruct$corStruct)
corMatrix(cs2$modelStruct$corStruct)[1] # match to SAS


Kevin Wright



On Wed, Aug 10, 2011 at 9:23 AM, Paul Miller <pjmiller_57@yahoo.com> wrote:

> Hello All,
>
> I sent the message below about two weeks ago and got no reply. I'm told
> that many people on the list were away at JSM in Miami around that time. My
> hope is that some people who otherwise would have replied were away and that
> my message simply fell through the cracks. So I'm going to try posting it
> again.
>
> Hope no one will mind.
>
> Thanks,
>
> Paul
>
>
> Hello Everyone,
>
> I'm A SAS user who is in the process of learning R. I've read about a dozen
> R-related books and have been testing my understanding of how to perform
> various analyses by replicating the results from one of my SAS books. Most
> of the time, I can replicate the results from SAS using R. And when I can't,
> it's ususally clear that there are differences in the way SAS and R handle
> the analysis.
>
> I've replicated almost all of the analyses in the of 22 chapters in my
> book, but am currently stuck on a couple of things in the repeated measures
> ANOVA chapter. Specifically, I'm having trouble with repeated measures ANOVA
> using mixed models and GEE.
>
> Below is some R code and output along with corresponding code and output
> from SAS. I tried using both lmer and nlme and found I was able to come
> closest to replicating the mixed model results using the latter. The results
> in R and SAS look very different and this makes me suspect that I'm using R
> incorrectly.
>
> Can anyone help me understand what I'm doing wrong?
>
> Thanks,
>
> Paul
>
>
> ###############################################################
> #### Example 8.3: Disease Progression in Alzheimer's Trial ####
> ###############################################################
>
> #### Data in Wide Format ####
>
> connection <- textConnection("
> TREAT PAT ADAS02 ADAS04 ADAS06 ADAS08 ADAS10 ADAS12
> L  1  22  30   .  33  28  30
> L  5  34  35  46  37  31  35
> L  8  40  41  41  46  52  48
> L 12  24   .  21  28  30  27
> L 13  29  26  29  26   .  36
> L 15  31  36  41  46  52  57
> L 19  22  27  28  24  27  28
> L 21  43  49  42  48  48  46
> L 24  18  28  29   .  25  28
> L 28  25  24  27  18  21  22
> L 31 37 35 35 38  42   .
> L 34  24  27  28  24  27  25
> L 37  45  50  58  59  60  58
> L 40  33  32  35  30  31  35
> L 44  34  37  43  44  39  38
> L 47  25 27 29 28 31   .
> L 51  30   .  36  32  34  38
> L 54  23   .  33  28  32  32
> L 57  35  37  39  38  41  43
> L 59  44  48  48  45  50  52
> L 63  28  30  32  31  35  32
> L 67  24  22  23  24  27  30
> L 68   .  49  51  48  55  54
> L 73  26  28  30  27  30  33
> L 76  30  32  35  35  36  38
> L 78  40  42  44  43  45  46
> H  2  31  36  35  31  31  31
> H  6  24  27  28  21  27  26
> H  7  31  31  39  37  41   .
> H 14  45  48  46  52  48  42
> H 17  24  28  26  23  24  29
> H 20  21  32  39  36  33  30
> H 22  32  34  45  42  37  32
> H 25  18  22  26  26  27  24
> H 27  51  47   .  43  43  43
> H 33  20  22  29  24  29  30
> H 38  41  34  37  29  35  33
> H 42  24  35  39  32  24   .
> H 45  23   .  33  36  33  30
> H 50  25  28  25  28  28  30
> H 52  31  34   .  33  34  35
> H 56  27  31  26  33  33  34
> H 60  37  43  39  42  43  36
> H 62  41  42  51  45  46  51
> H 66  35  33  34  35  36  41
> H 69  30  31  27  34  33  36
> H 72  54  60  55  58   .  65
> H 75  35  37  39  41  39  44
> H 79  18  21  19  19  20  27
> H 80  40  35  33  39  38  41
> P  3  31  36  37  41  39  44
> P  4  20  26  32  35  25  29
> P  9  33  33  29  33  39  41
> P 10  35  39  40  38  40  38
> P 11  26  24  31  42  50   .
> P 16  44  48  44  37  36  47
> P 18  25  31  21  27  41  32
> P 23  28  34  26  26  36  35
> P 26  27   .  28  35  40   .
> P 29  20  30  30  27  33  29
> P 30  49   .  43  48  44  53
> P 32  26  29  31  30  35  38
> P 35  30  33  41  .   41  44
> P 36  31  34  44  44  50  56
> P 39  42  46  36  43  48  48
> P 41  31  30  31   .  41  38
> P 43  27  22  36  45  54  60
> P 46  24  37  41  31  36  44
> P 48  33  31  38  41  31   .
> P 49  27  30  36  36  32  33
> P 53  35  34  45  44  38  40
> P 55  39  40  38  44  43  44
> P 58  32  34  40  45  36  38
> P 61  45  50   .  54  50  53
> P 64  21  23  31  34  27  27
> P 65  26  30  37  37  30  32
> P 70  53  50  55  57   .   .
> P 71  32  34  27  30  36  35
> P 74   .  50  52  56  52  54
> P 77  24  32  31  37  35  30
> ")
>
> alzhmrs <- read.table(connection, header=TRUE, sep="", na.strings=".")
> alzhmrs$TREAT <- factor(alzhmrs$TREAT, labels=c("High", "Low", "Placebo"))
> head(alzhmrs)
>
> #### Reshape Data to Long Format ####
>
> library(reshape)
> unialz <- melt(alzhmrs, id=c("TREAT","PAT"))
> names(unialz) <- c("TREAT", "PAT", "MONTH", "ADASCOG")
> unialz$TREAT <- relevel(unialz$TREAT,"Placebo")
> unialz$MONTH <- relevel(unialz$MONTH,"ADAS12")
> unialz <- sort_df(unialz,vars=c("PAT","MONTH"))
> head(unialz)
>
> #############################################
> #### Mixed Model Approach (Unstructured) ####
> #############################################
>
> #### R code and output ####
>
> library("nlme")
> model.un <- lme(ADASCOG ~ TREAT * MONTH, random = ~ TREAT | PAT,
>            correlation=corSymm(form = ~ 1 | PAT),
>            data=unialz, na.action=na.omit)
> summary(model.un)
> anova(model.un)
>
> > anova(model.un)
>            numDF denDF   F-value p-value
> (Intercept)     1   359 1560.5256  <.0001
> TREAT           2    77    0.8724  0.4220
> MONTH           5   359   24.6385  <.0001
> TREAT:MONTH    10   359    1.9795  0.0346
>
> #### SAS code and output ####
>
> ODS SELECT Mixed.RCorr Mixed.Dimensions Mixed.FitStatistics Mixed.Tests3;
>
> PROC MIXED DATA = UNIALZ;
>    CLASS TREAT MONTH PAT;
>    MODEL ADASCOG = TREAT MONTH TREAT*MONTH;
>    REPEATED / TYPE=UN SUBJECT=PAT(TREAT) RCORR;
>    TITLE4 'PROC MIXED Using Unstructured Covariance (UN)';
> RUN;
>
>   PROC MIXED Using Unstructured Covariance (UN)
>
>               The Mixed Procedure
>
>         Type 3 Tests of Fixed Effects
>
>                Num     Den
> Effect           DF      DF    F Value    Pr > F
>
> TREAT             2      77       0.94    0.3967
> MONTH             5      77      21.55    <.0001
> TREAT*MONTH      10      77       2.08    0.0357
>
> model.cs <- update(model.un, correlation=corCompSymm(form = ~ 1 | PAT))
> summary(model.cs)
> anova(model.cs)
>
> ####################################
> #### Mixed Model Approach (AR1) ####
> ####################################
>
> #### R code and output ####
>
> model.ar1 <- update(model.un, correlation=corAR1(form = ~ 1 | PAT))
> summary(model.ar1)
> anova(model.ar1)
>
> > anova(model.ar1)
>            numDF denDF   F-value p-value
> (Intercept)     1   359 1597.7569  <.0001
> TREAT           2    77    1.1192  0.3318
> MONTH           5   359   30.9931  <.0001
> TREAT:MONTH    10   359    2.3784  0.0097
>
> #### SAS code and output ####
>
> ODS SELECT Mixed.RCorr Mixed.Dimensions Mixed.FitStatistics Mixed.Tests3;
>
> PROC MIXED DATA = UNIALZ;
>    CLASS TREAT MONTH PAT;
>    MODEL ADASCOG = TREAT MONTH TREAT*MONTH;
>    REPEATED / TYPE=AR(1) SUBJECT=PAT(TREAT);
>    TITLE4 'PROC MIXED Using 1st-Order Auto-Regressive Covariance (AR(1))';
> RUN;
>
>         Type 3 Tests of Fixed Effects
>
>                Num     Den
> Effect           DF      DF    F Value    Pr > F
>
> TREAT             2      77       0.94    0.3952
> MONTH             5     359      11.21    <.0001
> TREAT*MONTH      10     359       1.43    0.1671
>
> ######################
> #### GEE Approach ####
> ######################
>
> #### R code and output ####
>
> library("geepack")
> model <- geeglm(ADASCOG ~ TREAT * MONTH, id = PAT,
>        data=unialz, family = gaussian, corstr = "ar1")
> summary(model)
> anova(model)
>
> > anova(model)
> Analysis of 'Wald statistic' Table
> Model: gaussian, link: identity
> Response: ADASCOG
> Terms added sequentially (first to last)
>
>            Df   X2 P(>|Chi|)
> TREAT        2  3.7   0.15595
> MONTH        5 95.0   < 2e-16 ***
> TREAT:MONTH 10 30.7   0.00065 ***
> ---
> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
> #### SAS code and output ####
>
> ODS SELECT Genmod.Type3;
> PROC GENMOD DATA = UNIALZ;
>    CLASS TREAT MONTH PAT;
>    MODEL ADASCOG = TREAT MONTH TREAT*MONTH / DIST=NORMAL TYPE3;
>    REPEATED SUBJECT=PAT / TYPE=AR(1);
>    TITLE4 'GEE Analysis Using PROC GENMOD';
>    TITLE5 'Autoregressive Correlation (AR(1)) Working Correlation';
> RUN;
>
>            GEE Analysis Using PROC GENMOD
> Autoregressive Correlation (AR(1)) Working Correlation
>
>                 The GENMOD Procedure
>
>       Score Statistics For Type 3 GEE Analysis
>
>                                Chi-
>    Source             DF     Square    Pr > ChiSq
>
>    TREAT               2       2.02        0.3639
>    MONTH               5      46.14        <.0001
>    TREAT*MONTH        10      21.76        0.0164
>
> _______________________________________________
> R-sig-mixed-models@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>

	[[alternative HTML version deleted]]


