ryuyan⽅差分析_R语⾔中实现⽅差分析及其可视化
⽅差分析(Analysis of Variance，简称ANOVA)，⼜称“变异数分析”，是R.A.Fisher发明的，⽤于两个及两个以上样本均数差别的显著性检验。 由于各种因素的影响，研究所得的数据呈现波动状。 造成波动的原因可分成两类，⼀是不可控的随机因素，另⼀是研究中施加的对结果形成影响的可控因素。
⽅差分析分类：
⽅差分析常⽤的为单因素⽅差分析和双因素⽅差分析，⽽多因素⽅差分析⽐较复杂⽤的不多，这⾥主要介绍前两种⽅差分析。
单因素⽅差分析
单因素⽅差分析⽐较的是分类因⼦定义的两个或多个组别中的因变量均值。与t检验区别：t检验适⽤于两列数据的均值⽐较。单因素⽅差分析适⽤于两列或更多列数据的均值⽐较。但对于两列数据的均值⽐较，单因素⽅差分析=等⽅差假设的双尾t检验。
以R中内置数据PlantGrowth集为例。不同地块( ‘ctrl’, ‘trt1’, and ‘trt2’)种植的植物的产量(weight)数据；解决的问题如下
1.不同地块产量是否存在差异
<1与trt2是哪块差异更⼤
##正态性检验
my_data
#结果表明，p=0.8915 >0.05，接受原假设，说明数据在多个⽔平下都是正态分布的。
##⽅差齐性检验
#结果表明，p=0.2371 >0.05，接受原假设，说明数据在不同⽔平下是等⽅差的。
##单因素⽅差分析
fit
summary(fit)
#            Df Sum Sq Mean Sq F value Pr(>F)
#group        2  3.766  1.8832  4.846 0.0159 *
#Residuals  27 10.492  0.3886
#---
#Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#结果表明，不同地块见产量差异⾮常显著。
#如果⽤⾮参数检验法 Kruskal-Wallis test
##多重⽐较(多重T检验)，哪块差异更明显
#R语⾔中TukeyHSD()函数提供了对各组均值差异的成对检验。
TukeyHSD(fit)
#如果⽤KW检验多重⽐较，
st(PlantGrowth$weight, PlantGrowth$group,
hod = "BH")
更好的⽅法+绘图展⽰：
利⽤ggpubr包做⽅差分析，并可视化分析结果：
library("ggpubr")
compare_means(weight ~ group, data = PlantGrowth,method="anova")
compare_means(weight ~ group, data = PlantGrowth,method="st")
my_comparisons
c("trt1", "trt2"))
ggboxplot(my_data, x = "group", y = "weight",
color = "group", palette = "jco",
order = c("ctrl", "trt1", "trt2"),
ylab = "Weight", xlab = "Treatment")+
stat_compare_means(label.y = 11,method = "anova")+
stat_compare_means(comparisons = my_comparisons,label.y = c(7,9,8))
或者⼩提琴图
ggviolin(my_data, x = "group", y = "weight",
fill = "group", palette = "jco",
order = c("ctrl", "trt1", "trt2"),
ylab = "Weight", xlab = "Treatment",add = "boxplot", add.params = list(fill="white"))+
stat_compare_means(label.y = 11,method = "anova")+    # Global p-value
stat_compare_means(comparisons = my_comparisons,label.y = c(7,9,8))
或者柱状图
ggbarplot(PlantGrowth, x = "group", y = "weight", add = "mean_se",
color = "group",fill="group", palette = "jco")+
stat_compare_means(label.y = 11,method = "anova") +    # Global p-value
stat_compare_means(comparisons = my_comparisons,label.y = c(7,9, 8),method="t.test")+ scale_y_continuous(expand=c(0,0))
comparisons或者线图
ggline(my_data, x = "group", y = "weight",
add = c("mean_se", "jitter"),
order = c("ctrl", "trt1", "trt2"),
ylab = "Weight", xlab = "Treatment")+
stat_compare_means(label.y = 11,method = "anova")+  # Global p-value
stat_compare_means(comparisons = my_comparisons,label.y = c(7,9,8))
延伸阅读
更多⽣物信息课程：
1. ⽂章越来越难发？是你没发现新思路，基因家族分析发2-4分⽂章简单快速，学习链接：基因家族分析实操课程、基因家族⽂献思路解读
2. 转录组数据理解不深⼊？图表看不懂？点击链接学习深⼊解读数据结果⽂件，学习链接：转录组(有参)结果解读；转录组(⽆参)结果解读
3. 转录组数据深⼊挖掘技能-WGCNA，提升你的⽂章档次，学习链接：WGCNA-加权基因共表达⽹络分析

版权声明：本站内容均来自互联网，仅供演示用，请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198，我们将在24小时内删除。