有些GEO平台的探针转换比较麻烦
这个帖子是为了回答一个学员的提问,
GPL4381的 平台的探针和gene symbol转换怎么做?
这个平台的信息在这个页面:
https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GPL4381
找不到ID和genesymbol的转换,准确地说是没有直接的对应关系
但是我们知道gene symbol藏在那串长长的名字中,此刻,即使不知道如何提取,也要想到R语言是万能的。
我们下载这个平台的信息
然后读入R语言
GPL4381 <- data.table::fread("GPL4381-4306.txt")
大概是这个样子的
我们先观察一下,提取一个名称该用什么逻辑。举个例子
Homo sapiens myelin expression factor 2 (MYEF2), mRNA
首先,他必须在两个括号之内,括号里面第一个元素是大写字母,后面是大写字母和数字的组合,那么stringr这个包中的str_extract可以实现。
str_extract("Homo sapiens myelin expression factor 2 (MYEF2), mRNA","\\([A-Z]{1}[A-Z,0-9]{2,}\\)")
[1] "(MYEF2)"
我这里编写的不好,能够实现功能,这种情况下如果碰到这样的元素,有两个括号
Homo sapiens intraflagellar transport 80 homolog (Chlamydomonas) (IFT80), mRNA
这个表达式也是能够顺利实现的
str_extract("Homo sapiens intraflagellar transport 80 homolog (Chlamydomonas) (IFT80), mRNA","\\([A-Z]{1}[A-Z,0-9]{2,}\\)")
(IFT80)
既然我们明确地定义了做一件事情的步骤,那么我们就可以批量做。
library(dplyr)
library(tidyr)
library(stringr)
probeGenesymbol <- GPL4381 %>%
select(ID,GB_DEFINITION) %>%
filter(GB_DEFINITION != "") %>%
mutate(gene_symbol = str_extract(.$GB_DEFINITION, "\\([A-Z]{1}[A-Z,0-9]{2,}\\)")) %>%
filter(gene_symbol !="") %>%
separate(gene_symbol,into=("gene_symbol"),sep = "\\)")%>%
separate(gene_symbol,into = c("drop","gene_symbol"),sep = "\\(")%>%
select(ID,gene_symbol,GB_DEFINITION)
最终把行数从14112变成了10715
我们来看看得到的probeGenesymbol
我们觉得效果还可以,这时候就可以用merge,inner_join等方法跟表达矩阵合并了
参看这个帖子:
skr!GEO芯片数据的探针ID转换
我们刚才发现,行数从14112变成了10715,理论上我们得到的是mRNA的信息,肯定是损失了非编码RNA的信息,这个信息如何获得呢?
回看平台信息
红框内是每一个探针的序列,如果把这个序列比对到基因组,可以对每一行的信息进行注释,
方法可以参考这个:
芯片探针序列的基因组注释
也可以在微店中检索小丫画图,其第58个作品就是如何提取lncRNA。
但是,对于初学者来说,都相对困难,需要有比较好的R语言基础才能驾驭。
好了总结一下今天帖子的三个注意点,
1 今天读入数据的函数是data.table中的fread,原因有两个,1是速度块,2是智能,大多时候不需要设置参数。我十分喜欢。
2 今天出现的字符串截取函数,str_match 来自stringr这个包,我们的读者应该是第一次看到,这里面涉及了一点正则表达式的知识,曾经有一位大神说过这样的话,“正则表达式是我们认识这个世界的哲学”,他的名字叫洲更
3 点号的特殊使用,在R语言的管道符号%>% 中,我们不需要设置中间变量,所有的中间变量就像水一样顺着管道往下流淌,那么如何在中间某一步使用上一步函数产生的数据呢?就用点号“.”代替。上过课的学员们应该知道,我们在GEO的数据处理时,探针去重那一步用到了这个方法。
而这个用法也隐藏在了这个帖子里面
skr!GEO芯片数据的探针ID转换
想想我出去培训学员,所有的中间材料和录屏全部赠送给学员,
资料:
视频:
这本质上是我让渡了讲师的一部分权益,但为了学员们能够掌握技能我觉得无所谓,而且我根本不怕有人录屏传播和抄袭。
原因有3个:
第一,很少有人会连续观看17个小时的教学视频,除非那是电视剧。
第二,大部分群众购买教程只是为了满足收藏癖,收藏就是雪藏,不会打开第二次。
第三,果子是活的,定期解答学员的问题。
下次见。