雙因子變異數分析(Two-way ANOVA)

最近接觸到滿多次雙因子變異數分析，想說就記錄一下了解這個統計方法的過程吧

還是先說明什麼是變異數分析(analysis of variance, ANOVA)，它是一種分析類別自變數對數值應變數影響的統計方法。
舉例說給雞添加兩種營養物質(A, B)，並分析這兩種處裡各自對體增重的影響。
添加營養物質就是類別自變數; 體增重為數值應變數。

跳回雙因子變異數分析吧，就是指兩項類別自變數對數值應變數的影響。
從剛才例子中衍生，原本只有一項類別自變數(營養物質)，又多加一項環境溫度(20, 25, 30℃)
也就是一隻雞同時接受兩種因子影響(如餵飼A物質且處於25℃)，組合總共有2X3= 6組處理

接下用R進行演示一邊說明如何使用two-way anova
對R有興趣可以搭配程式碼，一般直接看統計結果即可。
資料呈現如下表，每組有10隻雞的體增重(由於只是演示，請忽略體增重的合理性)

	20°C	25°C	30°C
A	4.2	16.5	23.6
	11.5	16.5	18.5
	7.3	15.2	33.9
	5.8	17.3	25.5
	6.4	22.5	32.5
	10	17.3	32.5
	11.2	13.6	26.7
	11.2	14.5	21.5
	5.2	18.8	23.3
	7	15.5	29.5
B	15.2	19.7	25.5
	21.5	23.3	26.4
	17.6	23.6	22.4
	9.7	26.4	24.5
	14.5	20	24.8
	10	25.2	30.9
	8.2	25.8	26.4
	9.4	21.2	27.3
	16.5	14.5	29.4
	9.7	27.3	23

一、輸入資料

setwd("C:/R/biostatistics")
data <- read.csv("data_twoway_anova.csv", header=T)
data$temp <- as.factor(data$temp)
```  

二、描述性統計(盒鬚圖、平均值及標準差)  
```R
#boxplot
boxplot(gain ~ supp + temp, col=c("red", "orange"), ylab="Body weight gain", xlab="supplement and temperature", data=data)
libraray(ggplot2)
ggplot(data, aes(x=supp, y=gain, col=temp)) + geom_boxplot()
 ```  
boxplot可使用內建函式(boxplot)或用ggplot繪圖
![](boxplot.png)
 
 ```R
# class by supp and temp, and then count means and standard error
library(reshape); library(agricolae)
mystats <- function(x) {c(n=length(x), mean=mean(x, na.rm=T), sd=sd(x))}
dfm <- melt(data, measure.vars="gain", id.vars=c("supp", "temp"))
cast(dfm, supp + temp + variable ~ ., mystates)

六個組別各自的平均值(mean)與標準差(sd)

	supp	temp	variable	n	mean	sd
1	A	20	gain	10	7.98	2.746634
2	A	25	gain	10	16.77	2.515309	
3	A	30	gain	10	26.14	4.797731	
4	B	20	gain	10	13.23	4.459709
5	B	25	gain	10	22.70	3.910953
6	B	30	gain	10	26.06	2.655058

三、主效果變異數分析
先不考慮交互作用，只分析營養物質(supp)和環境溫度(temp)個別對體增重(gain)的影響效果

# main effect ANOVA with no interaction
model_aov <- aov(data=data, gain ~ supp + temp)
summary(model_aov)
model_aov$cofficients  
# if we don't care supp and temp, mean of gain is 12.455 
# if suppOJ = 0. and then suppVC = -3.7