打开APP
userphoto
未登录

开通VIP,畅享免费电子书等14项超值服

开通VIP
R读取gmt文件

    前面我们简单介绍过什么是gmt文件,基因矩阵转置文件格式(* .gmt)。今天我们就用R来去读gmt文件。

    首先我们从GESA(https://www.gsea-msigdb.org/gsea/downloads.jsp)的官网上,下载一个gmt文件。这里以KEGG的gmt文件为例,其他gmt文件的读取方法一样。

c2.cp.kegg.v7.0.symbols.gmt这个文件里面保存的是基因的名字,

而c2.cp.kegg.v7.0.entrez.gmt这个文件里面保存的是基因的Entrez gene id,

Entrez gene id一般是以一串数字代表一个基因,这串数字直接贴到NCBI里面就可以查找到对应的基因名字了。

下面我们会用两种不同的方法来将KEGG symbol的gmt文件读到R里,并转换成列表。由于gmt文件的每一行都是不一样长的,所以传统的read.table在这里是毫无用武之地的。

方法一:

x <- readLines("c2.cp.kegg.v7.0.symbols.gmt") res <- strsplit(x, "\t") names(res) <- vapply(res, function(y) y[1], character(1))    res <- lapply(res, "[", -c(1:2))

该方法会将KEGG通路的名字作为列表中每个元素的名字,然后将前两列删掉,剩下的基因名字作为列表的元素

方法二:

dat = readLines("c2.cp.kegg.v7.0.symbols.gmt") n = length(dat) res = list(genesets = vector(mode = "list", length = n), geneset.names = vector(mode = "character", length = n), geneset.descriptions = vector(mode = "character", length = n)) for (i in 1:n) { s = strsplit(dat[i], "\t")[[1]] res$genesets[[i]] = s[-c(1:2)] res$geneset.names[i] = s[1] res$geneset.descriptions[i] = s[2] } names(res$genesets) = res$geneset.names    res

该方法,会保留gmt文件中的所有信息,结果会生成一个复杂的数据结构,列表里面嵌套列表。res为列表,长度为3,分别保存genesets,KEGG通路名字和数据来源,而geneset也是一个列表,里面保存186条KEGG通路上的所有基因名字。

本站仅提供存储服务,所有内容均由用户发布,如发现有害或侵权内容,请点击举报
打开APP,阅读全文并永久保存 查看更多类似文章
猜你喜欢
类似文章
【热】打开小程序,算一算2024你的财运
对象何必到处乱找,自己创造即可
GO富集可视化:clusterProfiler GOplot
TCGA转录组差异分析后多种基因功能富集分析:从GO/KEGG到GSEA和GSVA/ssGSEA(含基因ID转换)
0055
wikipathway : 代谢通路专用数据库
一个affymetrix表达芯片实战
更多类似文章 >>
生活服务
热点新闻
分享 收藏 导长图 关注 下载文章
绑定账号成功
后续可登录账号畅享VIP特权!
如果VIP功能使用有故障,
可点击这里联系客服!

联系客服