做泛癌时,20多个差异基因集需要取交集怎么办
R学习往期回顾:
R学习 从Tidyverse学起,入门R语言 dplyr合并数据
R学习 从Tidyverse学起,入门R语言(tidyr和stringr)
R学习 从Tidyverse学起,入门R语言(tibble,readr和dplyr)
今日正文
在做泛癌研究的时候,TCGA有30多个肿瘤,其中20多个是有正常组织的,可以用这些肿瘤做差异分析。当我们做了差异分析以后,需要选择共同上调或者下调的基因,而且20多个基因集做出的交集极有可能啥也没有,这要怎么办呢?
我们可以设定一个标准,比如在25个肿瘤中至少在9个肿瘤中是上调的,这样的基因我们认为是普遍上调的基因,那么问题来了,如何实现这样的计算?
我们先构建模拟数据,考虑到电脑的感受,就构建12个集合,在12个集合中选择至少在4个集合中存在的数据。
rm(list = ls())
set.seed(202000)
a<-sample(1:1500,100)
b<-sample(1:1500,100)
c<-sample(1:1500,100)
d<-sample(1:1500,100)
e<-sample(1:1500,100)
f<-sample(1:1500,100)
g<-sample(1:1500,100)
h<-sample(1:1500,100)
i<-sample(1:1500,100)
j<-sample(1:1500,100)
k<-sample(1:1500,100)
l<-sample(1:1500,100)
我们可以尝试做两两的交集
intersect(a,b)
intersect函数只能做两两的交集,那么如何做12个基因集的交集呢?有人说可以写个循环,现在还没这个必要,下面这段代码可以实现(不要问我咋知道的,我搜索的)
Reduce(intersect, list(a,b,c,d,e,f,g,h,i,j,k,l))
看到没有,交集啥都没有
下面就尝试一下4个数据集的交集,这会可以复习一下排列组合了,12选4有多少种组合呢?R语言种有个choose函数可以实现(不要问我咋知道的,我搜索的)
choose(12,4)
495种组合,还好,25选9可是有200多万种啊。
下面我们要知道所有的组合情况,我们先建立一个变量
labels<-c(a,b,c,d,e,f,g,h,i,j,k,l)
很显然这样是不对的,我们需要知道的是a到l的排列组合
labels<-c('a','b','c','d','e','f','g','h','i','j','k','l')
这就对了
接下来这段代码,可以知道所有的组合方式!(不要问我咋知道的,我搜索的)
res = t(combn(labels, 4))
res <-as.data.frame(res)
head(res)
每一行代表一种组合方式,总共有495行,跟上面算的一样。
接下来要干啥,当然是写个循环计算咯,我们可以把行当作变量,列当作常量,第一种组合是这样的
res[1,1]
res[1,2]
res[1,3]
res[1,4]
把1换成变量i,然后取每一行的这4个字符所对应的向量的交集。(不知道啥叫向量?看看前面的基础或者看本书)
这里值得注意的是,用res[1,1]取出来的是“factor”,而不是变量名,R语言中也是有一个函数可以将字符串转变成相对应的变量名(不要问我咋知道的,我搜索的)
get(as.character(res[1,1]))
解决了这个关键问题,咱就可以写个循环了计算一下了
output<-vector()
for (i in 1:495) {
value = Reduce(intersect, list(get(as.character(res[i,1])),get(as.character(res[i,2])),get(as.character(res[i,3])),
get(as.character(res[i,4])) ))
output<-c(value,output)
}
output
结果就是有4个数字在这12个数据集中的至少4个数据集中出现。
稍微验证一下
1453%in%a
1453%in%b
1453%in%c
1453%in%d
1453%in%e
1453%in%f
1453%in%g
1453%in%h
1453%in%i
1453%in%j
1453%in%k
1453%in%l
四个TRUE,结果是没错的,可怜了我的电脑,现在还在运行那200多万次的交集运算和200多万次的output,我也不知道要多久,主机也不是很热,反正拿个风扇对着吹。应该还有其他方法可以更快的运算,尽情搜索吧。
今日感触:只要关键词找的对,总会有答案的。
公众号“生信小课堂”