下面给它改成非均衡设计：
weight <- c(0.82,0.65,0.51,0.73,0.54,0.23,0.43,0.34,0.28,0.41,0.21,0.31,0.68,0.43,0.24)block <- c(rep(c("1","2","3","4","5"),each=3))# 每组样本量不一样group <- c(rep(c("A"),2),rep("B",9),rep("C",4))data4_4 <- data.frame(weight,block,group)
下面再看一下3种类型方差分析的结果：
# 1型fit1 <- aov(weight ~ block + group, data = http://www.kingceram.com/post/data4_4)summary(fit)##Df Sum Sq Mean Sq F valuePr(>F)## block4 0.2284 0.057095.978 0.01579 * ## group2 0.2280 0.1140011.937 0.00397 **## Residuals8 0.0764 0.00955## ---## Signif. codes:0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# 2型car::Anova(fit1, type = 2)## Anova Table (Type II tests)## ## Response: weight##Sum Sq Df F value Pr(>F)## block0.07978940.5896 0.6798## group0.03375020.4988 0.6250## Residuals 0.2706508# 3型car::Anova(fit1, type = 3)## Anova Table (Type III tests)## ## Response: weight##Sum Sq Df F valuePr(>F)## (Intercept) 1.080451 31.9364 0.0004807 ***## block0.0797940.5896 0.6798021## group0.0337520.4988 0.6249619## Residuals0.270658## ---## Signif. codes:0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
可以看到结果完全不一样了哦！
协方差分析
就用一个简单的完全随机设计资料的协方差分析进行演示 ， 示例数据来自课本例13-1 。
df13_1 <- data.frame(x1=c(10.8,11.6,10.6,9.0,11.2,9.9,10.6,10.4,9.6,10.5,10.6,9.9,9.5,9.7,10.7,9.2,10.5,11.0,10.1,10.7,8.5,10.0,10.4,9.7,9.4,9.2,10.5,11.2,9.6,8.0),y1=c(9.4,9.7,8.7,7.2,10.0,8.5,8.3,8.1,8.5,9.1,9.2,8.4,7.6,7.9,8.8,7.4,8.6,9.2,8.0,8.5,7.3,8.3,8.6,8.7,7.6,8.0,8.8,9.5,8.2,7.2),x2=c(10.4,9.7,9.9,9.8,11.1,8.2,8.8,10.0,9.0,9.4,8.9,10.3,9.3,9.2,10.9,9.2,9.2,10.4,11.2,11.1,11.0,8.6,9.3,10.3,10.3,9.8,10.5,10.7,10.4,9.4),y2=c(9.2,9.1,8.9,8.6,9.9,7.1,7.8,7.9,8.0,9.0,7.9,8.9,8.9,8.1,10.2,8.5,9.0,8.9,9.8,10.1,8.5,8.1,8.6,8.9,9.6,8.1,9.9,9.3,8.7,8.7),x3=c(9.8,11.2,10.7,9.6,10.1,9.8,10.1,10.3,11.0,10.5,9.2,10.1,10.4,10.0,8.4,10.1,9.3,10.5,11.1,10.5,9.7,9.2,9.3,10.4,10.0,10.3,9.9,9.4,8.3,9.2),y3=c(7.6,7.9,9.0,7.8,8.5,7.5,8.3,8.2,8.4,8.1,7.0,7.7,8.0,6.6,6.1,8.1,7.8,8.4,8.2,8.0,7.6,6.9,6.7,8.1,7.4,8.2,7.6,7.8,6.6,7.2))suppressPackageStartupMessages(library(tidyverse))df13_11 <- df13_1 %>% pivot_longer(cols = everything(), # 变长names_to = c(".value","group"),names_pattern = "(.)(.)") %>% mutate(group = as.factor(group)) # 组别变为因子型glimpse(df13_11)## Rows: 90## Columns: 3## $ group 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1…## $ x 10.8, 10.4, 9.8, 11.6, 9.7, 11.2, 10.6, 9.9, 10.7, 9.0, 9.8, 9.6…## $ y 9.4, 9.2, 7.6, 9.7, 9.1, 7.9, 8.7, 8.9, 9.0, 7.2, 8.6, 7.8, 10.0…
group是分组因素 ， x是协变量 ， y是因变量 。
下面进行3种类型的协方差分析：
# 1型fit <- aov(y ~ x + group, data = http://www.kingceram.com/post/df13_11) summary(fit)##Df Sum Sq Mean Sq F value Pr(>F)## x129.0629.057171.20 <2e-16 ***## group219.859.92558.48 <2e-16 ***## Residuals8614.600.170## ---## Signif. codes:0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# 2型car::Anova(fit, type = 2)## Anova Table (Type II tests)## ## Response: y##Sum Sq Df F valuePr(>F)## x30.1831177.83 < 2.2e-16 ***## group19.851258.48 < 2.2e-16 ***## Residuals 14.596 86## ---## Signif. codes:0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# 3型car::Anova(fit, type = 3)## Anova Table (Type III tests)## ## Response: y##Sum Sq DfF value Pr(>F)## (Intercept)0.381812.2493 0.1373## x30.18301 177.8346 <2e-16 ***## group19.8510258.4798 <2e-16 ***## Residuals14.5963 86## ---## Signif. codes:0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1