獨斷論

SAS 초급7: GLM을 이용한 이원분산분석 본문

과학과 기술/SAS

SAS 초급7: GLM을 이용한 이원분산분석

부르칸 2021. 7. 10. 01:28

SAS GLM을 이용하여 이원분산분석(two-way ANOVA)를 수행해보도록 한다.

1. 데이터파일:

drugeffect.csv
0.00MB

쉼표로 구분된 csv파일이고 첫번째 행에 변수가 포함되어 있다.

data hdlchol;
	infile 'd:\work\statistics\book\sas\tutorial\drugeffect.csv' 
		dlm="," firstobs=2;
	input treat $ age $ response;
run;

HDL콜레스테롤을 증가시키는 약을 개발하여 t1 group에는 가짜약을 넣고 t2에는 5mg을 투약하고 t3에는 10mg을 투약한 후에 나이대별로 HDL콜레스테롤이 어떻게 변화하는지 알아본 데이터이다. 데이터파일의 첫번째 열은 t1, t2, t3 그룹을 나타내고 두번째 열은 연령별 그룹이고 세번째 열은 콜레스테롤 증가량이다

 

2. 이원분산분석

proc glm data = hdlchol;
  class treat age;
  model response = treat age treat*age;
  means treat age treat*age / bon;
run; quit;

 

AMOVA 모델은 독립변수로 treat과 age, 그리고 interaction effect로  treat*age를 설정하고

Bonferoni t test로 multiple 그룹비교를 수행하였다.

SAS 코드에서 bon 대신에 scheffe 또는 tukey를 사용할수 있다.

 

treat의 main effect는 통계적으로 의미를 갖는다(F=25.68, p-value < 0.0001) 그러나 age는 p-value=0.086으로 통계적으로 유의미하다고 볼만한 충분한 증거를 가지 못하였다. interaction effect는 보이지 않는다(F = 1.08, p-value=0.37).

위 그림은  Bonferroni 테스트 결과로 오른쪽의 색깔이 서로 다른 treat 그룹을 통계적으로 유의미한 차이를 갖는다. 그러나 아래 age 그룹을 보면 모두 같은 색깔을 가지며 앞서 본 anova결과와 동일하다.

 

3. Residual plot

Residual plot은 glm procedure에 직접 그릴수 없으므로 결과 데이터를 따로 불러내 아래와 같이 수행한다.

proc glm data = hdlchol;
  class treat age; 
  model response = treat age treat*age;
  output out=anovaresults p=yhat r=resid;
run;

symbol color=blue value=dot;
proc gplot data = anovaresults;
  plot resid * yhat;
run;quit;

 

4. Interaction plot

interaction plot을 그리기 위해서 proc means를 이용하여 각 그룹의 level별 평균을 따로 구하고, proc gplot을 이용하여 앞에서 구한 평균값의 그래프를 그린다.

proc means nway data = hdlchol;
	class treat age;
	var response;
	output out=resultsMean mean=HDLmean;
run;

proc gplot data = resultsMean;
	plot HDLmean * treat = age;
run;

plot 명령어에서 HDLmean이 세로축, treat가 가로축을 나타내고 이 둘을  * 로 분리하였다. 그리고 age에 따라서 또 다른 그래프를 중첩하여 그리도록 하였는데 이들을 =로 분리하였다.

symbol을 십자가에서 원으로 바꾸고 각 점을 선으로 연결하려면 다음과 같이 한다.

symbol1 value=square color=black interpol=join;
symbol2 value=circle color=red interpol=join;
proc gplot data = resultsMean;
	plot HDLmean * treat = age;
run;

symbol1은 age의 첫번째 레벨이고 symbol2는 age의 두번째 레벨에 대한 그래프 옵션이다.

 

Comments