前面给大家介绍过
在上面这篇文章里面我们直接根据样本类型来分成了两类,即tumor和normal。
但更多的时候,我们其实想研究这些基因的表达跟临床特征的关系,如组织病理分期stage I, stage II, stage III, stage IV,TNM分期,T0, T1, T2, T3, T4,N0, N1, N2, N3(☞肿瘤TNM分期)等等。这个就跟前面讲的不太一样了,这些临床特征都大于两类,那么今天我们就来看看如何根据临床特征对样本进行分组绘制boxplot,并显示p值。
我们首先来捋一捋思路
1.从TCGA数据库下载想要研究肿瘤的表达谱数据和临床数据,这里以CHOL(胆管癌为例)
具体可以参考
☞如何从TCGA数据库下载RNAseq数据以及临床信息(一)
2.合并得到表达谱矩阵,并进行归一化
具体可以参考
3.搜集m6A甲基化相关的基因,根据需要做基因的ID转换
这里我已经帮大家从一些相关的研究中搜集了m6A甲基化相关的16个基因,当然大家也可以根据自己的需要增,删相关基因。我们讲的主要是方法,代码都是一样的跑,多几个基因少几个基因都无所谓。甚至你还可以换成其他的任何基因列表都是可以的。
"METTL3","METTL14","METTL16","WTAP","VIRMA","RBM15","ZC3H13","YTHDF1","YTHDF2","YTHDF3","YTHDC1","YTHDC2","HNRNPA2B1","HNRNPC","FTO","ALKBH5"
这里需要注意的是TCGA中gene的ID一般是ensembl gene ID,而我们从paper中搜集的gene一般是gene名字,所以这里就需要做一个gene ID转换。这个我在前面也讲过
4.提取m6A基因相关的表达矩阵
这个其实就是R里面最常见的数据框取下标操作了。这里所有分析都是基于R来做的,所以了解一些R的基础知识还是很有必要的。
具体可以参考
m6A基因相关的表达矩阵如下
5.添加样本的类型信息,如肿瘤,癌旁正常对照
前面也用视频给大家详细介绍过如何从TCGA数据下载临床信息
☞如何从TCGA数据库下载RNAseq数据以及临床信息(一)
6.下载并处理临床信息
前面我们也给大家讲解过
也分享过如果处理多分组的临床信息
这里我们将组织病理分期分成三组之后跟m6a基因的表达值合并起来
7.对数据格式进行转换
8.根据分组绘制boxplot
9.计算p值
在m6a甲基化相关基因boxplot并显示p值一文中,我们是根据样本类型来分类的,只有两类,所以我们可以做t检验来获取p值。本文中我们将样本根据组织病理分期分成了三类(当然四类也可以),stage I,stage II和stage III/IV。所以这里就不能再使用t检验来获取p值了,而应该使用aov方差检验来获取p值。
10.将计算得到的p值转换成对应的*,添加到boxplot上
可以参考☞【R语言】P值转换成***将p值转换成相应的*
对TCGA数据库还不了解的小伙伴,我们也有相关的视频课程详细介绍TCGA数据库挖掘
参考资料:
☞如何从TCGA数据库下载RNAseq数据以及临床信息(一)
完整代码参考