||
转录组数据结构
2、NCBI的 GEO数据集
(1)转录组数据主要有2种,一种是 Expression profiling by array,可下载的数据格式是 CEL,在 Series Matrix File(s)中下载的数据,可以用R读取:
gse112389.data <- read.table("GSE112389_series_matrix.txt",header = T,sep = "\t",comment.char = "!")
注释部分是以“!”开始的,可以用文本阅读器读取。
读出的数据中,是各探针的id检测到的表达值,探针ID(如244906_at)与基因ID(AT2G07783)之间的对应关系在GPL文件中,该表格可以下载,
读取gpl文件,需要Bioconductor包
library(Biobase)
library(GEOquery)
如果没有,用BiocManager::install("GEOquery", version = "3.8")来安装。
从网上直接获取GPL:
gpl198 <- getGEO("GPL198" , destdir = ".")
再次使用时,可以从本地读取:gpl198 <- getGEO(filename = "GPL198.soft")
提取探针id与基因id的关系部分,这里分了两步,
gpl198.dataTable <- gpl198@dataTable
gpl198.table <- gpl198.dataTable@table
可以一步完成:gpl198.table <- gpl198@dataTable@table
把表达数据和探针信息合并
两个表的探针id顺序一致时,可以简单的合并:
gse112389.express <- cbind(gse112389.data,gpl198.table)
或者,先统一两个数据集共同列的名称
colnames(gpl198.table)[1] <- "ID_REF"
然后融合
gse112389.express <- merge(gse112389.data,gpl198.table,by = "ID_REF")
(2)另一种转录组数据是 Expression profiling by high throughput sequencing,这种直接下载raw文件,解包后自己组合所需要的数据,拟南芥中直接就是基因id,不需要其它处理。
d1 <- read.table("GSM3373992_D_1.txt", sep = "\t")
d2 <- read.table("GSM3373993_D_2.txt", sep = "\t")
d3 <- read.table("GSM3373994_D_3.txt", sep = "\t")
基因名称列表顺序相同时,可以直接合并
d0 <- cbind(d1,d2,d3)
或者:d0 <- merge(merge(d1,d2, by = "V1"),d3,by = "V1"),
更多文件时,放在一个目录下,用dir()获取文件名后,用for循环处理。
expfiles = dir()
expdata = read.table(expfiles[1],sep = "\t")
colnames(expdata) = c("ID",substr(expfiles[1],1,14) )
for ( i in 2:length(expfiles)){
expdata0 = read.table(expfiles[i],sep = "\t")
colnames(expdata0) = c("ID",substr(expfiles[i],1,14) )
expdata = merge(expdata,expdata0,by = "ID")
}
最后可以输出成TXT文件
write.table(expdata,"expr.txt")
Archiver|手机版|科学网 ( 京ICP备07017567号-12 )
GMT+8, 2024-9-1 17:21
Powered by ScienceNet.cn
Copyright © 2007- 中国科学报社