데이터마이너를 꿈꾸며

제3장 일원배치 완전확률화 계획법,

예제) 직물의 긁힘에 대한 저항력을 측정하기 위해 원단 A1,2,3,4의 제품 마모도를 비교하조가 한다.

> setwd('c:/Rwork')

> fabric = read.csv('fabric.csv')

> head(fabric)

result group

1   1.93    a1

2   2.38    a1

3   2.20    a1

4   2.25    a1

5   2.55    a2

6   2.72    a2

> boxplot(result~group,data=fabric)

> with(fabric,tapply(result,group,mean))

a1   a2   a3   a4

2.19 2.68 2.42 2.31

> round(with(fabric,tapply(result,group,mean)))

a1 a2 a3 a4

2  3  2  2

> with(fabric,tapply(result,group,sd))

a1         a2         a3         a4

0.18920888 0.08906926 0.18036999 0.05715476

> out=lm(result~group,data=fabric)

> out

Call:

  lm(formula = result ~ group, data = fabric)

Coefficients:

  (Intercept)      groupa2      groupa3      groupa4 

2.19         0.49         0.23         0.12 

Lm의 결과물을 보면 intercept 추정치가 a1의 평균과 같으므로 a1이 reference로 사용되었음을알 수 있다. A1와 a2 평균이 각각 2.19, 2.68으로 그 차이가 groupa2 임을 알 수 있다.

> summary(out)

Call:

  lm(formula = result ~ group, data = fabric)

Residuals:

  Min      1Q  Median      3Q     Max

-0.2600 -0.0700  0.0150  0.0625  0.2600

Coefficients:

  Estimate Std. Error t value Pr(>|t|)   

(Intercept)  2.19000    0.07050  31.062 7.79e-13 ***

  groupa2      0.49000    0.09971   4.914 0.000357 ***

  groupa3      0.23000    0.09971   2.307 0.039710 * 

  groupa4      0.12000    0.09971   1.204 0.251982   

---

  Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.141 on 12 degrees of freedom

Multiple R-squared:  0.6871,           Adjusted R-squared:  0.6089

F-statistic: 8.785 on 3 and 12 DF,  p-value: 0.002353

> anova(out)

Analysis of Variance Table

Response: result

Df Sum Sq  Mean Sq F value   Pr(>F)  

group     3 0.5240 0.174667  8.7846 0.002353 **

Residuals 12 0.2386 0.019883                   

---

  Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

분산분석표에서 group의 F(df1=3,df2=12)=8.7846 (p-value=0.002353)으로 세 그룹의 평균들이 모두 같다는 귀무가설을 기각한다.

귀무가설 H0: μa1 = μa2 = μa3 = μa4

대립가설 H1: not H0

검정통계량 F(df1=3,df2=12)=8.7846

p-value = 0.002353 < 0.05

결정 : 모든 그룹의 평균이 같다는 귀무가설을 기각한다. 4개의 그룹 중 어느 쌍에 차이가 있는지 다중비교로 좀 더 진단할 수 있다.

> par(mfrow=c(2,2))

> plot(out)

> abline(out,col='blue')

Normal Q-Q Plot이 직선에 가깝고 다른 그림들도 별다른 추세를 보이지 않으므로 문제가 없다는 결론을 내린다.

> shapiro.test(resid(out))

Shapiro-Wilk normality test

data:  resid(out)

W = 0.97318, p-value = 0.8872

Shapiro.test()를 이용하여 잔차의 정규성검정 결과는 p-value=0.8872로 정규분포 가정을 만족시킨다.

모든 평균이 같다는 귀무가설이 기각되었다는 말은 그룹 중 최소한 하나는 0이 아니다라는 말이다. 어느 쌍의 차이로 귀무가설이 기각되었는지 조사하기 위해 다중비교를 한다. 분산분석에서 많이 쓰이는 다중비교 방법은 Dunnett와 Tukey이다. Tukey는 가능한 모든 조합의 쌍을, Dunnett는 하나의 대조군(reference)을 나머지 비교군(treatment)들과 비교한다.

> library(multcomp)

> dunnett=glht(out,linfct=mcp(group='Dunnett'))

> dunnett

General Linear Hypotheses

Multiple Comparisons of Means: Dunnett Contrasts

Linear Hypotheses:

  Estimate

a2 - a1 == 0     0.49

a3 - a1 == 0     0.23

a4 - a1 == 0     0.12

> summary(dunnett)

Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Dunnett Contrasts

Fit: lm(formula = result ~ group, data = fabric)

Linear Hypotheses:

  Estimate Std. Error t value Pr(>|t|)   

a2 - a1 == 0  0.49000    0.09971   4.914   <0.001 ***

a3 - a1 == 0  0.23000    0.09971   2.307   0.0962 . 

a4 - a1 == 0  0.12000    0.09971   1.204   0.5085   

---

  Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Adjusted p values reported -- single-step method)

> plot(dunnett)

a2-a1의 p-value만 0.001로 유의

> tukey=glht(out,linfct=mcp(group='Tukey'))

> tukey

General Linear Hypotheses

Multiple Comparisons of Means: Tukey Contrasts

Linear Hypotheses:

  Estimate

a2 - a1 == 0     0.49

a3 - a1 == 0     0.23

a4 - a1 == 0     0.12

a3 - a2 == 0    -0.26

a4 - a2 == 0    -0.37

a4 - a3 == 0    -0.11

> summary(tukey)

Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Tukey Contrasts

Fit: lm(formula = result ~ group, data = fabric)

Linear Hypotheses:

  Estimate Std. Error t value Pr(>|t|)  

a2 - a1 == 0  0.49000    0.09971   4.914  0.00187 **

a3 - a1 == 0  0.23000    0.09971   2.307  0.15096  

a4 - a1 == 0  0.12000    0.09971   1.204  0.63632  

a3 - a2 == 0 -0.26000    0.09971  -2.608  0.09215 .

a4 - a2 == 0 -0.37000    0.09971  -3.711  0.01375 *

a4 - a3 == 0 -0.11000    0.09971  -1.103  0.69444  

---

  Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Adjusted p values reported -- single-step method)

> plot(tukey)

a2-a1, a4-a2간 95% 신뢰구간이 0을 포함하지 않으므로 유의한 차이가 있음을 알 수 있다.

분산분석은 모든 그룹의 평균이 같은지를 유의수준 0.05로 한 번에 검정한다. 그러나 여러 쌍의 평균이 같은지 다시 검정하면 비교하는 쌍이 하나 이상으로 늘어나 실제로 유의하지 않은데 유의하게 나올 가능성이 높아지므로 유의수준을 0.05보다 낮추거나, p-value를 높게 조종해야 한다. 아래 표처럼 실제 Tukey의 p-value가 Dunnett의 p-value보다 높게 조정되었음을 알 수 있다.

출처: 실험 계획과 응용, R로 하는 통계 분석