Dear all:for the folowing data, a two-period, two treatment (A=1 vs. B=2) cross-over is fitted using the folowing SAS code.? data one; input? Sbj? Seq ?Per? Trt? PEF; cards; 1???? ?1?? 1?? 1?? 310 1???? ?1?? 2?? 2?? 270 4???? ?1?? 1?? 1?? 310 4???? ?1?? 2?? 2?? 260 6???? ?1?? 1?? 1?? 370 6???? ?1?? 2?? 2?? 300 7????? 1?? 1?? 1?? 410 7????? 1?? 2?? 2?? 390 10??? 1?? 1?? 1?? 250 10??? 1?? 2?? 2?? 210 11??? 1?? 1?? 1?? 380 11??? 1?? 2?? 2?? 350 14??? 1?? 1?? 1?? 330 14??? 1?? 2?? 2?? 365 2????? 2?? 1?? 2?? 370 2????? 2?? 2?? 1?? 385 3????? 2?? 1?? 2?? 310 3???? ?2?? 2?? 1?? 400 5???? ?2?? 1?? 2?? 380 5????? 2?? 2?? 1?? 410 9????? 2?? 1?? 2?? 290 9????? 2?? 2?? 1?? 320 12??? 2?? 1?? 2?? 260 12??? 2?? 2?? 1?? 340 13??? 2?? 1?? 2??? 90 13??? 2?? 2?? 1?? 220 ; run; proc mixed data=one method=reml; class Sbj Per Trt; ?? model PEF = Per Trt /ddfm=kr; ?? repeated Trt / sub=Sbj type=un r; ?? lsmeans Trt / cl alpha=0.05; ?? estimate 'B vs. A' Trt -1? 1 / alpha=0.1 cl; run; (where kr option is for Kenward-Roger method).I need to use R to reproduce the results similar to what the above SAS code generates. I have used several R functions including lme, lmer with no success so far.Any advice will be greatly appreciated,Sincerely, Keramat [[alternative HTML version deleted]]
knouri <nouri4 <at> yahoo.com> writes:> > Dear all:for the folowing data, a two-period, two treatment (A=1 vs. B=2) > cross-over is fitted > using the folowing SAS code.? > data one;[snip]> run; > proc mixed data=one method=reml; > class Sbj Per Trt; > ?? model PEF = Per Trt /ddfm=kr; > ?? repeated Trt / sub=Sbj type=un r; > ?? lsmeans Trt / cl alpha=0.05; > ?? estimate 'B vs. A' Trt -1? 1 / alpha=0.1 cl; > run;> (where kr option is for Kenward-Roger method).I need to use R to > reproduce the results similar to what the above SAS code generates. > I have used several R functions including lme, lmer with no success > so far.Any advice will be greatly appreciated,Sincerely, KeramatThis is more appropriate for r-sig-mixed-models at r-project.org. Please post followups there. The lmerTest and lsmeans packages will probably be useful. As a statistical point, I don't understand why you can't just do a paired t-test on these data?? dat <- read.table(header=TRUE,text"Sbj Seq Per Trt PEF 1 1 1 1 310 1 1 2 2 270 4 1 1 1 310 4 1 2 2 260 6 1 1 1 370 6 1 2 2 300 7 1 1 1 410 7 1 2 2 390 10 1 1 1 250 10 1 2 2 210 11 1 1 1 380 11 1 2 2 350 14 1 1 1 330 14 1 2 2 365 2 2 1 2 370 2 2 2 1 385 3 2 1 2 310 3 2 2 1 400 5 2 1 2 380 5 2 2 1 410 9 2 1 2 290 9 2 2 1 320 12 2 1 2 260 12 2 2 1 340 13 2 1 2 90 13 2 2 1 220") library(lmerTest) library(ggplot2); theme_set(theme_bw()) ggplot(dat,aes(x=Per,y=PEF,colour=factor(Trt)))+geom_point()+ geom_line(colour="gray",aes(group=Sbj))+ scale_x_continuous(breaks=c(1,2)) m1 <- lmer(PEF~Per+Trt +(Trt|Sbj), data=dat) ## warning about unidentifiability
On 08 Jun 2015, at 22:59 , Ben Bolker <bbolker at gmail.com> wrote:> As a statistical point, I don't understand why you can't just > do a paired t-test on these data??It's a cross-over trial, so you need to account for the period effect. However, a standard way of treating it is by to calculate period differences within subjects and use a t-test (welch or classical) to see if the differ according to sequence (trt1-trt2 or v.v.). -- Peter Dalgaard, Professor, Center for Statistics, Copenhagen Business School Solbjerg Plads 3, 2000 Frederiksberg, Denmark Phone: (+45)38153501 Office: A 4.23 Email: pd.mes at cbs.dk Priv: PDalgd at gmail.com