##### One-way ANOVA
> PlantGrowth = read.csv("PlantGrowth.csv")
> with(PlantGrowth, tapply(weight, group, mean))
ctrl trt1 trt2
5.032 4.661 5.526
> with(PlantGrowth, tapply(weight,group,sd))
ctrl trt1 trt2
0.5830914 0.7936757 0.4425733
# boxplot
> boxplot(weight ~ group, col="red", data=PlantGrowth)
![](https://t1.daumcdn.net/cfile/cafe/1324ED1D4B2903B778)
# fit one-way anova
> out = lm(weight ~ group, data=PlantGrowth)
> anova(out)
Analysis of Variance Table
Response: weight
Df Sum Sq Mean Sq F value Pr(>F)
group 2 3.7663 1.8832 4.8461 0.01591 *
Residuals 27 10.4921 0.3886
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# 회귀진단
> par(mfrow=c(2,2))
> plot(out)
> shapiro.test(resid(out))
# Install "multcomp" library and call the library
install.packages("multcomp")
library(multcomp)
### Dunnett
> out = lm(weight ~ group, data=PlantGrowth)
> dunnett = glht(out,linfct=mcp(group="Dunnett"))
> summary(dunnett)
Simultaneous Tests for General Linear Hypotheses
Multiple Comparisons of Means: Dunnett Contrasts
Fit: lm(formula = weight ~ group, data = PlantGrowth)
Linear Hypotheses:
Estimate Std. Error t value Pr(>|t|)
trt1 - ctrl == 0 -0.3710 0.2788 -1.331 0.323
trt2 - ctrl == 0 0.4940 0.2788 1.772 0.153
(Adjusted p values reported -- single-step method)
> plot(dunnett)
![](https://t1.daumcdn.net/cfile/cafe/111A4A1A4B2905691E)
### Tukey
> tukey = glht(out,linfct=mcp(group="Tukey"))
> summary(tukey)
Simultaneous Tests for General Linear Hypotheses
Multiple Comparisons of Means: Tukey Contrasts
Fit: lm(formula = weight ~ group, data = PlantGrowth)
Linear Hypotheses:
Estimate Std. Error t value Pr(>|t|)
trt1 - ctrl == 0 -0.3710 0.2788 -1.331 0.3908
trt2 - ctrl == 0 0.4940 0.2788 1.772 0.1979
trt2 - trt1 == 0 0.8650 0.2788 3.103 0.0119 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)
![](https://t1.daumcdn.net/cfile/cafe/121A4A1A4B2905691F)
> out = lm(weight ~ group, data=PlantGrowth)
> weight = rnorm(n=nrow(PlantGrowth),mean=fitted(out),sd=summary(out)$sigma)
##### Kurskal-Wallis Test
> kruskal.test(weight~group, data=PlantGrowth)
Kruskal-Wallis rank sum test
data: weight by group
Kruskal-Wallis chi-squared = 7.9882, df = 2, p-value = 0.01842
첫댓글 음.. 이해가 잘 안돼요. 복사해서 붙여봐도 에러메세지가...
> dunnett=glht(out,linfct=mcp(group='Dunnett'))
Error in mcp2matrix(model, linfct = linfct) :
Variable(s) ‘group’ have been specified in ‘linfct’ but cannot be found i$
저는 잘 되는데요. 이걸 한꺼번에 돌려도 안되나요?
install.packages("multcomp")
library(multcomp)
out = lm(weight ~ group, data=PlantGrowth)
dunnett = glht(out,linfct=mcp(group="Dunnett"))
install.packages("multcomp") 부분은 RStudio나 R 버전에 따라서 실행이 안되기도 하네요
저도 책에서 install.packages("multcomp") 보고 입력했더니
Warning message:
package ‘mulcomp’ is not available (for R version 3.1.0)
라고 뜨네요.
그래서 바로 위에 써주신
install.packages("multcomp")
library(multcomp)
out = lm(weight ~ group, data=PlantGrowth)
dunnett = glht(out,linfct=mcp(group="Dunnett"))
복사하고 엔터 치니깐 되네요. 왜그런거죠?
이상하네요. 저는 3.1.0에서 해봤는데 잘 되는데요.
RStudio 버전은 상관없을텐데, R은 최신버전에서 해야할겁니다.
안선생님 질문있습니다!
13장에 Kruskal-Wallis test 가 one-way ANOVA 의 비모수 방법이라고 하셨는데, 그럼 one-way ANOVA의 회귀진단으로부터 잔차의 가정이 만족되지 않는경우에 Kruskal-Wallis test를 쓰는게 아닌가요? 12장에서 PlantGrowth 데이터를 이용해서 weight~group로 lm하고 plot(out)과 shapiro.test(resid(out))으로 잔차의 가정을 체크했을 때 모두 만족되었었는데, 13장에서 왜 weight~group 으로 Kruskal-Wallis test를 하는 건가 해서요...답변 주시면 감사하겠습니다!
Kruskal-Wallis test에 적합한 데이터가 없어서 그냥 one-way anova에 썼던 데이터를 사용했습니다^^
@안재형 아하. 그랬군요. ^^ 그럼 원래대로라면 one-way ANOVA의 회귀진단으로 부터 가정이 만족되지 않았을 때 Kruskal-wallis test를 쓰는게 맞는거군요. ㅎㅎ 답변 감사합니다! 이제 책의 끝을 향해 달려가고 있습니다. 이론적으로만 배웠던 것들을 직접 응용하면서 해보니 정리가 잘 되는거 같아요!! XD
@NiKe 책에 설명을 좀 했었어야했는데요^^ 2판이 나가게 되면 적절한 데이터를 찾던지, 그런 얘기를 추가하겠습니다. 감사합니다~
안녕하세요,
12.5에서
dunnett = glht(out,linfct=mcp(group="Dunnett"))
명령어를 실행하면
Error in glht(out, linfct = mcp(group = "Dunnett")) :
could not find function "glht"
Error in mcp2matrix(model, linfct = linfct) :
Variable(s) ‘group’ of class ‘character’ is/are not contained as a factor in ‘model’.
라고 오류가 나오는데요, 왜그런지 알 수 있을까요?
해결법도 궁금합니다.
install.packages("multcomp")
설치하셨나요?