TCGA的原始数据下载后,只有一串字母混合的代号,没有样本号,所有的信息都存在metadata中,但是他的形式是json
现在我们来读取json,需要安装jsonlite包。
metadata <- jsonlite::fromJSON("metadata.cart.2017-11-15T09_56_59.722935.json")
library(tidyverse)
metadata_id <- metadata %>%
select(c(file_name,associated_entities))
最终读到1208个观测,15个变量,我需要的是file_name和样本名称,样本名称藏在了associated_entities 列表中
里面包括了entity_id,case_id,entity_submitter_id,entity_type这四个项目,查看第一个了解一下
metadata$associated_entities[1]
[[1]]
entity_id case_id
1 52033f64-1e6f-4657-a4fb-7cfeffc61951 39de7761-e762-4811-b95c-8216b79ae06b
entity_submitter_id entity_type
1 TCGA-AN-A0XW-01A-11R-A109-07 aliquot
现在的想法是我把ilename和associated_entities中的entity_submitter_id提取出来,放入列表中,形成1208个小列表
每个小列表中包含两个元素,即lename和entity_submitter_id
name_id <- list()
for (i in 1:1208){
name_id[[i]] <- substr(metadata_id$file_name[i],1,nchar(metadata_id$file_name[i])-3)
name_id[[i]][2] <- metadata_id$associated_entities[i][[1]]$entity_submitter_id
}
不是很方便尝试直接把list转化成datafram,感受一下do.call的另类用法
naid_df <- as.data.frame(do.call(rbind, name_id))
得到的naid_df是一个2列1208行的数据框,实际上一开始避开list也可以完成,只是我选择了list,正好使用一下do.call
现在把1208个小文件读入一个矩阵文件,并且给每一个文件加上filename和entity_submitter_id
之前有个帖子介绍过
在R语言中将多个同样的行列式文件合并起来
但是当时不知道,TCGA的单个文件是没有列名的,导致无法合并,本次要复杂一点
#读入所有解压的文件 1208个,参考[R语言中选取多个文件夹中的文件合并到新的文件夹](http://guoshipeng.com/2017/11/13/15-multiple-files-into-one-by-R/)
nameList <- list.files("data_unzip/")
location <- which(naid_df==nameList[1],arr.ind = TRUE) ##which函数有一个已知value返回坐标的功能
TCGA_id <- as.character(naid_df[location[1],2]) ##通过坐标,获取TCGA_id
expr_df<- read.table(paste0("data_unzip/",nameList[1]),stringsAsFactors = F, header = F) #读入第一个文件,保存为data.frame
names(expr_df) <- c("gene_id",TCGA_id) #给刚才数据库命名
这边开始批量作业
for (i in 2:length(nameList)){
location <- which(naid_df==nameList[i],arr.ind = TRUE)
TCGA_id <- as.character(naid_df[location[1],2])
dfnew <- read.table(paste0("data_unzip/",nameList[i]),stringsAsFactors = F,header = F)
names(dfnew) <- c("gene_id",TCGA_id)
expr_df <- inner_join(expr_df,dfnew,by="gene_id")
}
晚上走的时候没运行完,早上来的时候已经完毕,限速环节应该是read.table,早上再来尝试运行一次总是说内存不够
我想到之前阅读GTF时,fread的速度很快,就试了一下
library(data.table)
nameList <- list.files("data_unzip/")
location <- which(naid_df==nameList[1],arr.ind = TRUE)
TCGA_id <- as.character(naid_df[location[1],2])
expr_df<- fread(paste0("data_unzip/",nameList[1]))
names(expr_df) <- c("gene_id",TCGA_id)
for (i in 2:length(nameList)){
location <- which(naid_df==nameList[i],arr.ind = TRUE)
TCGA_id <- as.character(naid_df[location[1],2])
dfnew <- fread(paste0("data_unzip/",nameList[i]))
names(dfnew) <- c("gene_id",TCGA_id)
expr_df <- inner_join(expr_df,dfnew,by="gene_id")
}
大概花了15分钟左右
save(expr_df,file = "brca_expr_df.Rda")
如果下次使用
load(file = "brca_expr_df.Rda")
更新:无意间看一个帖子说list的获取可以使用$符号,但是在储存的时候要有特定的格式,相当于把第一个变量当做名称,直接¥调用
如果成功的话,之前汇总文件的过程就会简单一点,应该是这种格式的循环
name_id <- list()
for (i in 1:1208){
list_name <- substr(metadata_id$file_name[i],1,nchar(metadata_id$file_name[i])-3)
name_id$list_name <- metadata_id$associated_entities[i][[1]]$entity_submitter_id
}
使用第一个文件试一下
list_name <- substr(metadata_id$file_name[1],1,nchar(metadata_id$file_name[1])-3)
name_id$list_name <- metadata_id$associated_entities[1][[1]]$entity_submitter_id
测试 一下
name_id
$list_name
[1] "TCGA-AN-A0XW-01A-11R-A109-07"
不对,直接把list_name的表达式放进去又不识别,换成下面的表达式还是不行
name_id <- list(substr(metadata_id$file_name[1],1,nchar(metadata_id$file_name[1])-3) = metadata_id$associated_entities[1][[1]]$entity_submitter_id)
看到$符号,我想到我可以把这些数据弄进数据框,只有两行,行名是filename,把TCGAid弄成一行
试一下第一个数据
name_id <- data.frame()
list_name <- substr(metadata_id$file_name[1],1,nchar(metadata_id$file_name[1])-3)
name_id$list_name <- metadata_id$associated_entities[1][[1]]$entity_submitter_id
不能运行,换一种表达方式
name_id <- data.frame()
name_id[1,1] <- metadata_id$associated_entities[1][[1]]$entity_submitter_id
names(name_id)[1] <- substr(metadata_id$file_name[1],1,nchar(metadata_id$file_name[1])-3)
成功,批量运行
name_id <- data.frame()
for (i in 1:1208){
name_id[1,i] <- metadata_id$associated_entities[i][[1]]$entity_submitter_id
names(name_id)[i] <- substr(metadata_id$file_name[i],1,nchar(metadata_id$file_name[i])-3)
}
运行成功,速度也很快,下面提取每个数据
nameList <- list.files("data_unzip/")
TCGA_id <- name_id$nameList[1]#这一步失败,最后发现是字符串中有-号,暂时没能解决
总结:
fread速度很快
R语言的data.frame如果列名包含“-”,无法使用$获取元素,在list中也是一样的情况
等以后看能不能解决这个问题吧
网友评论