-
-
Notifications
You must be signed in to change notification settings - Fork 2
/
Copy pathanovaSPFpq.R
78 lines (63 loc) · 3.19 KB
/
anovaSPFpq.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
## ------------------------------------------------------------------------
wants <- c("car", "DescTools", "multcomp")
has <- wants %in% rownames(installed.packages())
if(any(!has)) install.packages(wants[!has])
## ------------------------------------------------------------------------
set.seed(123)
Nj <- 10
P <- 3
Q <- 3
muJK <- c(rep(c(1,-1,-2), Nj), rep(c(2,1,-1), Nj), rep(c(3,3,0), Nj))
dfSPFpqL <- data.frame(id=factor(rep(1:(P*Nj), times=Q)),
IVbtw=factor(rep(LETTERS[1:P], times=Q*Nj)),
IVwth=factor(rep(1:Q, each=P*Nj)),
DV=rnorm(Nj*P*Q, muJK, 3))
## ------------------------------------------------------------------------
aovSPFpq <- aov(DV ~ IVbtw*IVwth + Error(id/IVwth), data=dfSPFpqL)
summary(aovSPFpq)
## ------------------------------------------------------------------------
dfSPFpqW <- reshape(dfSPFpqL, v.names="DV", timevar="IVwth",
idvar=c("id", "IVbtw"), direction="wide")
## ------------------------------------------------------------------------
library(car)
fitSPFpq <- lm(cbind(DV.1, DV.2, DV.3) ~ IVbtw, data=dfSPFpqW)
inSPFpq <- data.frame(IVwth=gl(Q, 1))
AnovaSPFpq <- Anova(fitSPFpq, idata=inSPFpq, idesign=~IVwth)
summary(AnovaSPFpq, multivariate=FALSE, univariate=TRUE)
## ------------------------------------------------------------------------
anova(fitSPFpq, M=~1, X=~0, idata=inSPFpq, test="Spherical")
anova(fitSPFpq, M=~IVwth, X=~1, idata=inSPFpq, test="Spherical")
## ------------------------------------------------------------------------
mauchly.test(fitSPFpq, M=~IVwth, X=~1, idata=inSPFpq)
## ------------------------------------------------------------------------
library(DescTools)
EtaSq(aovSPFpq, type=1)
## ------------------------------------------------------------------------
summary(aov(DV ~ IVbtw, data=dfSPFpqL, subset=(IVwth==1)))
summary(aov(DV ~ IVbtw, data=dfSPFpqL, subset=(IVwth==2)))
summary(aov(DV ~ IVbtw, data=dfSPFpqL, subset=(IVwth==3)))
## ------------------------------------------------------------------------
summary(aov(DV ~ IVwth + Error(id/IVwth), data=dfSPFpqL,
subset=(IVbtw=="A")))
summary(aov(DV ~ IVwth + Error(id/IVwth), data=dfSPFpqL,
subset=(IVbtw=="B")))
summary(aov(DV ~ IVwth + Error(id/IVwth), data=dfSPFpqL,
subset=(IVbtw=="C")))
## ------------------------------------------------------------------------
mDf <- aggregate(DV ~ id + IVbtw, data=dfSPFpqL, FUN=mean)
aovRes <- aov(DV ~ IVbtw, data=mDf)
## ------------------------------------------------------------------------
cMat <- rbind("-0.5*(A+B)+C"=c(-1/2, -1/2, 1),
"A-C"=c(-1, 0, 1))
## ------------------------------------------------------------------------
library(multcomp)
summary(glht(aovRes, linfct=mcp(IVbtw=cMat), alternative="greater"),
test=adjusted("none"))
## ------------------------------------------------------------------------
try(detach(package:DescTools))
try(detach(package:multcomp))
try(detach(package:mvtnorm))
try(detach(package:car))
try(detach(package:survival))
try(detach(package:splines))
try(detach(package:TH.data))