반응형

6.1 회귀분석

) 페인트의 불순도는 페인트를 얼마나 빨리 저어주느냐에 따라 달라진다. 아래표는 휘젓는 장치의 회전율과 불순도를 측정한 데이터이다.

x(회전율)

20

22

24

26

28

30

32

34

36

38

40

42

y(불순도)

8.4

9.5

11.8

10.4

13.3

14.8

13.2

14.7

16.4

16.5

18.9

18.5

> setwd('c:/Rwork')

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

> paint

x    y

1  20  8.4

2  22  9.5

3  24 11.8

4  26 10.4

5  28 13.3

6  30 14.8

7  32 13.2

8  34 14.7

9  36 16.4

10 38 16.5

11 40 18.9

12 42 18.5

> out=lm(y~x,data=paint)

> plot(y~x,data=paint)

> abline(out)


> out

 

Call:

  lm(formula = y ~ x, data = paint)

 

Coefficients:

  (Intercept)            x 

-0.2893       0.4566 

 

> summary(out)

 

Call:

  lm(formula = y ~ x, data = paint)

 

Residuals:

  Min      1Q  Median      3Q     Max

-1.1834 -0.5432 -0.3233  0.8333  1.3900

 

Coefficients:

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

(Intercept) -0.28928    1.22079  -0.237    0.817   

x            0.45664    0.03844  11.880 3.21e-07 ***

  ---

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

 

Residual standard error: 0.9193 on 10 degrees of freedom

Multiple R-squared:  0.9338,           Adjusted R-squared:  0.9272

F-statistic: 141.1 on 1 and 10 DF,  p-value: 3.211e-07

 

회귀식 Y = -0.28928 + 0.45664X

기울기 0.45664 p-value 3.21e-07로 유의.

결정계수R2 0.9338, 총 변동 중 93.38%가 회귀모형에 의해 설명되고 있다. 이는 두 변수 사이의 피어슨 상관계수의 제곱이다.


> with(paint,cor(x,y))^2

[1] 0.933832

 

 F(df1=1,df2=10) 141.1 p-value 3.21e-07로 유의.


cor.test()로 구한 Pearson 상관계수의 t(df=10)=11.88, p-value 3.211e-07은 위의 p-value값과 일치한다. T 값인 11.88의 제곱  11.88^2 값은 141.1344로 모형의 적합성을 나타내는 F 141.1과 일치하고 당연히 p-value3.21e-07로 동일하다.

> with(paint,cor.test(x,y))

 

Pearson's product-moment correlation

 

data:  x and y

t = 11.88, df = 10, p-value = 3.211e-07

alternative hypothesis: true correlation is not equal to 0

95 percent confidence interval:

0.8810937 0.9907768

sample estimates:

cor

0.9663498

 

out으로 저장된 결과를 plot()을 이용하여 회귀진단에 필요한 그래프를 얻을 수 있다.

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

> plot(out)

 

해석:표준화된 잔차의 Normal Q-Q plot이 직선에 가깝고 다른 그래프들도 별다른 추세를 보이지 않는다.

 

resid()을 이용하면 회귀분석 결과물로부터 잔차를 어을 수 있다. 이렇게 받은 잔차에 아래의 명령어를 이용하여 normal Q-Q Plot을 그릴 수 있다.

> qqnorm(resid(out))

> qqline(resid(out))


 

정확한 p-value를 알고 싶다면 shapiro.test()를 이용하여 잔차가 정규분포를 따르는지 검정한다.

> shapiro.test(resid(out))

 

Shapiro-Wilk normality test

 

data:  resid(out)

W = 0.9196, p-value = 0.2826

H0 정규분포를 따른다

H1 정규분포를 따르지 않는다

 

p-value 0.2826로 잔차는 정규분포 가정을 만족시킨다분산분석을 통하여도 회귀의 유의성을 검정할 수 있다.

> anova(out)

Analysis of Variance Table

 

Response: y

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

x          1 119.275 119.275  141.13 3.211e-07 ***

Residuals  10   8.451   0.845                     

---

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

 


6.2 공분산분석

실험이 잘 통제되어 종속변수Y의 변동을 설명하는 데 그룹변수 이외에 다른 변인이 없다면 Two-sample t-test나 일원분산분석을 하면 된다. 그러나 현실적으로 통제는 쉽지 않기 때문에 이 경우에는 여러가지 다른 변수들을 통제해주어 조사하고자 하는 변수만의 효과를 조사해야 한다. 공분산분석은 분산분석에 연속형 변수를 추가한 것이다. 궁극적인 목적인 각 그룹 간 평균들의 차이가 있는 검정하는 것으로 분산분석과 동일하나, 통제가 안 되는 연속형 변수(covariate)를 추가하여 오차를 줄이고 검정력을 높이는 것이 차이점이다.

기계1

기계2

기계3

y 강도

x 섬유 두께

y 강도

x 섬유 두께

y 강도

x 섬유 두께

36

20

40

22

35

21

41

25

48

28

37

23

39

24

39

22

42

26

42

25

45

30

34

21

49

32

44

28

32

15

 

lm(y~공변량변수 + 그룹변수)

> machine<-read.csv('machine.csv')

> head(machine)

machine  y  x

1      m1 36 20

2      m1 41 25

3      m1 39 24

4      m1 42 25

5      m1 49 32

6      m2 40 22

> levels(machine$machine)

[1] "m1" "m2" "m3"

해석: m1기계1이 대조군이고 m2 m3은 비교군이다.

 

> out=lm(y~x+machine,data=machine)

> anova(out)

Analysis of Variance Table

 

Response: y

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

x          1 305.130 305.130 119.9330 2.96e-07 ***

machine    2  13.284   6.642   2.6106   0.1181   

Residuals 11  27.986   2.544                     

---

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

 해석: p-value 0.1181로 기계 간의 차이가 있다고 결론을 내릴 수 없다



summary()로 자세한 결과물을 보자.

> summary(out)

 

Call:

  lm(formula = y ~ x + machine, data = machine)

 

Residuals:

  Min      1Q  Median      3Q     Max

-2.0160 -0.9586 -0.3841  0.9518  2.8920

 

Coefficients:

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

(Intercept)   17.360      2.961   5.862 0.000109 ***

x              0.954      0.114   8.365 4.26e-06 ***

machinem2      1.037      1.013   1.024 0.328012   

machinem3     -1.584      1.107  -1.431 0.180292   

---

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

 

Residual standard error: 1.595 on 11 degrees of freedom

Multiple R-squared:  0.9192,        Adjusted R-squared:  0.8972

F-statistic: 41.72 on 3 and 11 DF,  p-value: 2.665e-06

해석: m2 reference m1와 비교한 machinem2 p-value값이 0.328012이고 m3 m1과 비교한 machinem3 p-value 0.180292이다. 그러나 이 p-value들을 조정없이 사용하여 두 machine 모두 m1와 유의한 차이가 없다고 결론을 내리지는 않는다.

 

1개의 대조군(m1) 2개의 비교군과 비교하므로 Dunnett의 방법으로 p-value를 조정한다.

> library(multcomp)

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

> summary(dunnett)

 

Simultaneous Tests for General Linear Hypotheses

 

Multiple Comparisons of Means: Dunnett Contrasts

 

 

Fit: lm(formula = y ~ x + machine, data = machine)

 

Linear Hypotheses:

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

m2 - m1 == 0    1.037      1.013   1.024    0.518

m3 - m1 == 0   -1.584      1.107  -1.431    0.304

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

> plot(dunnett)

 


 해석:  p-value 모두가 원래 p-value보다 커졌으며, 여전히 유의한 차이가 있다고 결론을 지을 수 없다. Plot 그래프로 확인할 수도 있다. 신뢰구간이 0을 포함하고 있으므로 유의하다고 볼 수 없다.

 

summary(out)

dunnett

m2

0.328012

0.518

m3

0.180292

0.304

 

> with(data=machine,tapply(y,machine,mean))

m1   m2   m3

41.4 43.2 36.0

> with(data=machine,tapply(y,machine,sd))

m1       m2       m3

4.827007 3.701351 3.807887

해석: 각 기계에서 생산된 섬유 제품의 강도 평균은 41.4, 43.2, 36.0으로 얻어져 혹시 차이가 있다면 기계3 강도가 떨어지는 하다.


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

반응형
Posted by 마르띤
,
반응형

이원산분석은 그룹변수가 하나인 일원산분석의 확장으로 두 개의 그룹변수를 가진다. 이 분석의 특징은 두 그룹변수들의 효과뿐만 아니라 두 그룹 변수들이 서로 어떤 영향을 미치는지 교호작용을 볼 수 있다는 것이다

 

1) 교호 작용이 있는 경우: lm(y~변수a*변수b, data=)

2) 교호 작용이 없는 경우: summary(lm(y~변수a+변수b, data = )

 

예시) 각기 다른 사료 a1~a4와 돼지품종 b1~3이 체중 증가에 미치는 영향을 조사

 

b1

b2

b3

a1

64

66

70

72

81

64

74

51

65

a2

65

63

58

57

43

52

47

58

67

a3

59

68

65

66

71

59

58

39

42

a4

58

41

46

57

61

53

53

59

38

 

> setwd('c:/Rwork')

> pigs<-read.csv('pigs.csv')

> head(pigs)

result breed feed

1     64    b1   a1

2     66    b1   a1

3     70    b1   a1

4     65    b1   a2

5     63    b1   a2

6     58    b1   a2

> with(pigs,tapply(result,breed,mean))

b1       b2       b3

60.25000 61.33333 54.25000

> with(pigs,tapply(result,feed,mean))

a1       a2       a3       a4

67.44444 56.66667 58.55556 51.77778

> with(pigs,tapply(result,list(breed,feed),mean))

a1       a2       a3       a4

b1 66.66667 62.00000 64.00000 48.33333

b2 72.33333 50.66667 65.33333 57.00000

b3 63.33333 57.33333 46.33333 50.00000

> boxplot(result~breed+feed,col='blue',data=pigs)



 

기술통계를 위해 tapply 함수를 활용하여 품종, 사료, 품종x사료 조합의 평균을 구하였고, boxplot 함수를 통해 변수 조합의 result 분포 차이를 나타냈다.


 

교호 작용 검사

> with(pigs,interaction.plot(x.factor=feed,trace.factor=breed,response=result,fun=mean,type='b',main='Interaction Plot',pch=c(1,2,16)))

 


 

#교호 작용이 있는 경우

> out=lm(result~breed*feed,data=pigs)

> out

 

 

 

> anova(out)

Analysis of Variance Table

 

Response: result

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

breed       2  349.39  174.69  2.7926 0.081214 .

feed        3 1156.56  385.52  6.1628 0.002939 **

breed:feed   6  771.28  128.55  2.0549 0.097118 .

Residuals  24 1501.33   62.56                   

---

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

품종(breed)는 유의하지 않지만 p-value 0.05에 근접하고, 사료(feed)는 유의하다. 품종에 따라 사료가 주는 효과의 크기가 조금씩 다를 수 있지만, Interaction breed:feed p-value 0.09로 유의하지 않으므로 교호작용이 있다는 결론을 내릴 수 없다. 유의성 검정 결과를 보면 사료의 종류에 따라 체중증가에 유의한 차이가 있으나, 돼지 품종에 따라서는 체중증가에 차이가 있다는 결론은 내릴 수 없다. 품종에 대한 유의확률은 0.081로 약간의 영향력을 행사하는 듯하다.



#교호 작용이 없는 경우

 to be updated..



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


반응형
Posted by 마르띤
,
반응형

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

 

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

제품

반복

A1

1.93

2.38

2.2

2.25

A2

2.55

2.72

2.75

2.7

A3

2.4

2.68

2.32

2.28

A4

2.33

2.38

2.28

2.25

> 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-a295% 신뢰구간이 0을 포함하지 않으므로 유의한 차이가 있음을 알 수 있다.

 

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

 

 

Dunnett

Tukey

a2 vs a1

0.001

0.00187

a3 vs a1

0.0962

0.15096

a4 vs a1

0.5085

0.63632

a3 vs a2

 

0.09215

a4 vs a2

 

0.01375

a4 vs a3

 

0.69444

 

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

반응형
Posted by 마르띤
,
반응형

2.3 짝지어진 비교 Paired T-Test

Paired T-Test는 쌍을 이룬 두 변수의 차이를 보는 검정이다. 가장 흔한 예는 한 집단을 대상으로 어떤 개입의 효과를 보기 위해 개입 전-후 값을 비교하여 개입의 효과를 측정하는 것이다.

Paired t-test는 다음 순서를 따른다.


 1) 정규분포를 따르는지 검정 with(데이터, shapiro.test(사후-사전)

 2) 정규분포를 따르면 사후-사전 차이의 평균이 0인지 검정하는 one sample t-test를 이용

   with(데이터, t.test(사후-사전))

 3) 정규분포를 따르지 않으면 비모수적 방법을 적용 with(데이터, wilcox.test(사후-사전)

 

) 운동화 밑창에 사용되는 B 재질과 A 재질 간 마모도에 차이가 있는지 알아보고자 한다. 이를 위하여 10명의 아이들을 대상으로 임의로 한 쪽 발에는 재질 A의 밑창을 단 운동화를 한 쪽 발에는 재질 B의 밑창을 단 운동화를 신겨 일정기간 사용하도록 한 다음 밑창이 얼마나 마모되었는지를 측정한 결과 다음 자료를 얻었다.

아이

재질A

재질B

차이B-A

1

13.2

14.0

0.8

2

8.2

8.8

0.6

3

10.9

11.2

0.3

4

14.3

14.2

-0.1

5

10.7

11.8

1.1

6

6.6

6.4

-0.2

7

9.5

9.8

0.3

8

10.8

11.3

0.5

9

8.8

9.3

0.5

10

13.3

13.6

0.3

 

 

 

평균차이=0.41

 

> shoes<-read.csv('shoes.csv',header=T)

> head(shoes)

A    B

1 13.2 14.0

2  8.2  8.8

3 10.9 11.2

4 14.3 14.2

5 10.7 11.8

6  6.6  6.4

#paired t-test 1. 정규분포를 따르는지 검정

> with(shoes,shapiro.test(B-A))

 

Shapiro-Wilk normality test

 

data:  B - A

W = 0.96132, p-value = 0.8009

 

귀무가설 H0: B-A의 차이가 정규분포를 따른다.

대립가설 H1: B-A의 차이가 정규분포를 따르지 않는다.

검정통계량 W = 0.96132

p-value = 0.8009

결정 : 귀무가설 기각할 수 없다. B-A의 차이가 정규 분포를 따른다.

 

#paired t-test 2. 정규분포를 따르면 one sample t-test 이용

> with(shoes,t.test(B-A))

 

One Sample t-test

 

data:  B - A

t = 3.3489, df = 9, p-value = 0.008539

alternative hypothesis: true mean is not equal to 0

95 percent confidence interval:

  0.1330461 0.6869539

sample estimates:

  mean of x

0.41

귀무가설 H0: μ = 0

대립가설 H1: μ ≠ 0

검정통계량 t(df=9) = 3.3489

p-value = 0.008539

결정 : 귀무가설 기각한다. 두 신발 재질의 마모도에 유의한 차이가 있다.

 

#paired t-test 3. 정규분포를 따르지 않는 경우 wilcox.test 사용한다.

> with(shoes,wilcox.test(B-A)

귀무가설 H0: 전후 차이의 median 0이다

대립가설 H1: 전후 차이의 median 0이 아니다

검정통계량 V =

p-value =

결정: if p-value > 0.05 then 귀무가설 기각할 수 없다, 전후 차이가 없다.


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

반응형
Posted by 마르띤
,
반응형

2.2 독립표본을 이용한 두 모평균 차이에 대한 추론. Two-sample T-test

독립표본을 바탕으로 두 개의 모집단의 평균을 비교. 가장 흔한 실험 연구는 실험군과 대조군에 서로 다른 개입(intervention)을 적용시킨 후 두 집단의 평균이 같은지를 비교하여 개입 효과의 차이를 평가하는 것이다. 이 경우 two-sample t-test를 사용하는데, 서로 독립적인 두 변수 간에 차이의 평균이 0인지를 검정한다.

Two-sample t-test는 다음 순서를 따른다.

 1) 두 집단의 분산이 같은지 검정한다. var.test(y~그룹변수)

 2) 분산이 다르면 Welch t-test를 적용한다. t.test(y~그룹변수)

 3) 분산이 같으면 pooled variance를 이용한 t-test를 적용한다. t.test(y~그룹변수, var.equal=TRUE)

 

예제) 제약회사에서 어떤 약을 오래 보관해도 약효가 지속되는지를 검사하려고 한다. 표본1 2를 랜덤추출한 결과가 아래와 같다.

표본1

10.2

10.5

10.3

10.8

9.8

10.6

10.7

10.2

10.0

10.1

표본2

9.8

9.6

10.1

10.2

10.1

9.7

9.5

9.6

9.8

9.9

 

> medical<-read.csv('medical.csv',header=T)

> head(medical,3)

sample result

1 sample1   10.2

2 sample1   10.5

3 sample1   10.3

> tail(medical,3)

sample result

18 sample2    9.6

19 sample2    9.8

20 sample2    9.9

> boxplot(result~sample,data=medical)


-> 해석: Sample1의 분산이 sample2보다 더 큼을 알 수 있다.

 

#two sample test 1. 등분산 검정

> var.test(result~sample,data=medical)

F test to compare two variances

 

data:  result by sample

F = 1.7965, num df = 9, denom df = 9, p-value = 0.3959

alternative hypothesis: true ratio of variances is not equal to 1

95 percent confidence interval:

  0.4462364 7.2328801

sample estimates:

  ratio of variances

1.796545  

귀무가설 H0: σ21 = σ22

대립가설 H1: σ21 ≠ σ22

검정통계량 F(df1=9, df2=9) = 1.7965

p-value = 0.3959

결정 : 귀무가설 기각할 수 없다. 등분산을 가정한다.

 

> 1/1.7965  #F 1.7965의 역수는 0.55 sample2의 분산이 sample1대비 0.55배임을 알 수 있다.

[1] 0.5566379

 

#two sample test 2. 분산이 같은 경우, pooled variance사용

> t.test(result~sample,var.equal=TRUE,data=medical)

 

Two Sample t-test

 

data:  result by sample

t = 3.8511, df = 18, p-value = 0.00117

alternative hypothesis: true difference in means is not equal to 0

95 percent confidence interval:

  0.222688 0.757312

sample estimates:

  mean in group sample1 mean in group sample2

10.32                  9.83

 

귀무가설 H0: μ1 = μ2

대립가설 H1: μ1 ≠ μ2

검정통계량 t(df=18) = 3.8511

p-value = 0.00117

결정 : 귀무가설 기각한다. 두 모집단의 평균이 다르다.

 

#two sample test 3. 분산이 다른 경우, Welch t-test 한다.

> t.test(result~sample,data=medical)

 

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

반응형
Posted by 마르띤
,