专栏首页生信技能树十步搞定单因素方差分析

十步搞定单因素方差分析

看过我TCGA肿瘤数据库知识图谱的小伙伴都只是如何在任意癌症查询指定感兴趣基因的表达量,并且对样本进行分组比较,网站是:https://xenabrowser.net/heatmap/ 悄咪咪的上线了TCGA知识图谱视频教程(B站和YouTube直达)

image-20190506132847002

根据视频教程拿到数据,很容易可视化如下:

image-20190506132955029

可以得出结论,我们感兴趣的基因(这里是CUL5)在乳腺癌的正常组织及癌症组织(原位和转移)表达量,使用单因素方差分析,得到了统计学显著的结果。

定义

单因素方差分析是两个样本平均数比较的引伸,它是用来检验多个平均数之间的差异,从而确定因素对试验结果有无显著性影响的一种统计方法。

  • 因素:影响研究对象的某一指标、变量。
  • 水平:因素变化的各种状态或因素变化所分的等级或组别。
  • 单因素试验:考虑的因素只有一个的试验叫单因素试验。

了解数据

数据文件可以下载,然后读入R里面进行可视化,代码如下:

rm(list = ls())
options(stringsAsFactors = F)
library(ggpubr)
df=read.table('CUL5-BRCA-type.tsv',header = T,sep = '\t')
colnames(df)

p = ggboxplot(df, "sample_type", "CUL5",
          color = "sample_type", palette =c("#00AFBB", "#E7B800", "#FC4E07"),
          add = "jitter", shape = "sample_type")
p

这个文件 CUL5-BRCA-type.tsv 如果你没有看我的TCGA肿瘤数据库知识图谱可能不知道如何下载,可以发邮件给我找我申请这个测试数据 ( 邮箱: jmzeng1314@163.com )

可以看到,比网页工具出图要好看:

image-20190506133842005

但是,这个时候还没有进行统计分析,可以添加的统计学检验包括:

image-20190506133655190

代码也很简单:

p+stat_compare_means(method = "anova", label.y = 10)+      # Add global p-value
  stat_compare_means(label = "p.signif", method = "t.test",
                     ref.group = ".all.")                  # Pairwise comparison against all

image-20190506134047019

可以看到,跟网页工具结果一模一样,而且出图更漂亮,下面我们就手把手带领大家完成这个分析,把这个ggpubr一步就完成的工作拆解开来。

第1步:计算各组内样本均值

lapply(split(df,df$sample_type),function(x) mean(x$CUL5))

第2步:计算所有样本均值

mean(df$CUL5 )

第3步:计算各组内部误差平方和

tmp=lapply(split(df,df$sample_type),function(x) sum((x$CUL5-mean(x$CUL5))^2) )
sse = sum(unlist(tmp))

第4步:计算各组间误差平方和

tmp=lapply(split(df,df$sample_type),function(x) nrow(x)*(mean(x$CUL5) - mean(df$CUL5 ))^2 )
ssb = sum(unlist(tmp))

第5步:计算各组内部均方误

mse=sse/(nrow(df)-length(unique(df$sample_type)))

第6步:计算组间均方误

msd=ssb/length(unique(df$sample_type))-1

第7步:计算F比率

f= msb/mse
f

第8步:查找F临界值

df1=(length(unique(df$sample_type))-1)
df2=(nrow(df)-length(unique(df$sample_type)))
qf(0.05,2,1215)

差别可知这里的F值是0.05,远小于我们真实情况,所以非常显著了。

第9步:判断是否显著

1-pf(f,2,1215)

现在我们已经知道了,在选定的显著水平为0.05时候,这个F统计是显著的,但是仍然是不知道哪组之间不一样, 所以可以选择tukey检验

第10步:进行tukey检验,多重比较

J·W·图凯(Tukey)于1953年提出一种能将所有各对平均值同时比较的方法,这种方法现在已被广泛采用,一般称之为“HSD检验法”,或称“W法”。 Tukey (John Wilder Tukey) for multiple comparisons 主要应用于3组或以上的多重比较。比如说一共有4组数据,两两比较产生6个统计值,Tukey test用于生成一个critical value来控制总体误差(Familywise error rate,FER);与Tukey test相类似的是Dunnett test,它是控制多对一比较(即3组同时和一个参照组比较)的FER。

这个多重比较算法还蛮多的,参考:https://zhuanlan.zhihu.com/p/44880434 这里超纲了,我就不具体介绍了!

本文分享自微信公众号 - 生信技能树(biotrainee),作者:生信技能树

原文出处及转载信息见文内详细说明,如有侵权,请联系 yunjia_community@tencent.com 删除。

原始发表时间:2019-05-07

本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。

我来说两句

0 条评论
登录 后参与评论

相关文章

  • NPC的突变特性(逆向收费读文献2019-13)

    本周尝试一种新的文献解读方式,就是一次性解读3篇文献,总结他们的异同点,整体把握该领域发展情况。

    生信技能树
  • 【直播】我的基因组74:快速给测序reads比对到物种

    其实这一讲只是把未比对到人类基因组的序列快速比对到细菌基因组,并得到各个种类的菌的占比。 在这之前我们讲的是对几亿条reads定位到指定参考基因组的具体某个坐标...

    生信技能树
  • FDA致病菌检测挑战赛结果发布,厦门极元科技总分排名第二

    近期,美国食品药品监督局FDA通过在线平台precisionFDA (https://precision.fda.gov) 举办了CFSAN Pathogen ...

    生信技能树
  • [数据分析工具] Pandas 功能介绍(二)

    条件过滤 我们需要看第一季度的数据是怎样的,就需要使用条件过滤 ? 体感的舒适适湿度是40-70,我们试着过滤出体感舒适湿度的数据 ? 最后整合上面两种条件,在...

    数据分析
  • Pandas从小白到大师

    说明:本文内容翻译、节选自外文From Pandas-wan to Pandas-master[1],原作者Rudolf Höhn小哥,实验数据来自kaggle...

    zhangqibot
  • 基础知识 | 每日一面(16)

    小林:为了确保分配连续的结构数组时正确对齐, 结构可能有这种尾部填充。即使结构不是数组的成员, 填充也会保持, 以便 sizeof 能够总是返回一致的大小。

    闫小林
  • [译]Python正则表达式拾珠

    Python标准库中有很多非常恶心的模块,但Python的re模块不是其中之一。虽然它已经很老了而且多年未更新,它仍是我认为的众多动态语言中最好的(正则表达式模...

    岂不美哉Frost
  • Python+django网页设计入门(14):使用中间件对抗爬虫

    Python+django网页设计入门(12):使用Bootstrap和jQuery

    Python小屋屋主
  • 微课|玩转Python轻松过二级(2.4节):常用内置函数用法精要2

    Python小屋屋主
  • 从根上理解 Docker 是个啥

    早都想写 Docker 了,只是一直太忙!直到最近越来越多的网友再咨询我,Docker 为什么有优势?Docker 和虚拟机有什么本质上的区别?。。。所以,整体...

    业余草

扫码关注云+社区

领取腾讯云代金券