參考：公眾號(hào) 果子學(xué)生信 https://mp.weixin.qq.com/s/rdFnq6jCMIjmrWI8A8fS5g
GDC網(wǎng)站
點(diǎn)擊repository
點(diǎn)擊case
確定數(shù)據(jù)下載的組織
點(diǎn)擊file選擇下載數(shù)據(jù)的類型
下載數(shù)據(jù)
此時(shí)點(diǎn)擊購物車出吹，就會(huì)進(jìn)入下載頁面，因?yàn)楫?dāng)前數(shù)據(jù)挺大的励负，不建議直接下載而是采用他推薦的GDC Data Transfer Tool來進(jìn)行而使用這個(gè)工具温亲，需要兩個(gè)文件棚壁，一個(gè)是Manifest，一個(gè)是metadata
R下載數(shù)據(jù)
### 1.創(chuàng)建數(shù)據(jù)下載文件夾
if(!file.exists("rawdata")) dir.create("rawdata")
### 2.生成下載代碼
### 這個(gè)就是之前下載的manifest文件魂务，每個(gè)人的不一樣在刺，注意修改
manifest <- "gdc_manifest_20220929_023209.txt"
### 這個(gè)是下載的文件夾
rawdata <- "rawdata"
command <- sprintf("./gdc-client download -m %s -d %s",
manifest,
rawdata)
### 3.執(zhí)行下載操作
system(command = command)
整理數(shù)據(jù)
這時(shí)候所有數(shù)據(jù)都在rawdata這個(gè)文件夾中，總共有1226個(gè)文件夾相艇，每一個(gè)里面點(diǎn)開后有一個(gè)文件
其中第一列 gene_id我們是需要的颖杏，第四列unstranded 就是傳統(tǒng)的counts 數(shù)據(jù)，以前的TCGA數(shù)據(jù)只提供兩列坛芽。
其中第二列是轉(zhuǎn)換過后的基因名稱留储，第三列是基因類型，包括編碼和非編碼信息咙轩，這個(gè)可以保留一份用來做基因注釋
我們只需要第1列和第4列
1.把所有的tsv文件获讳，拷貝到新的文件夾data_in_one
#合并文件
if(!file.exists("data_in_one")) dir.create("data_in_one")
### 使用for循環(huán)來批量做，回顧for循環(huán)的要點(diǎn)
for (dirname in dir("rawdata/")){
## 要查看的單個(gè)文件夾的絕對(duì)路徑
mydir <- paste0(getwd(),"/rawdata/",dirname)
##找到對(duì)應(yīng)文件夾中的文件并提取名稱活喊，pattern表示模式丐膝，可以是正則表達(dá)式
file <- list.files(mydir,pattern = "_counts")
## 當(dāng)前文件的絕對(duì)路徑是
myfile <- paste0(mydir,"/",file)
#復(fù)制這個(gè)文件到目的文件夾
file.copy(myfile,"data_in_one")
}
2.提取文件名稱和TCGAID 的對(duì)應(yīng)關(guān)系【metadta文件中】
#提取文件名稱和TCGAID 的對(duì)應(yīng)關(guān)系
metadata <- jsonlite::fromJSON("metadata.cart.2022-09-29.json")
library(dplyr)
metadata_id <- metadata %>%
dplyr::select(c(file_name,associated_entities))
metadata
metadataid
## 提取對(duì)應(yīng)的文件
naid_df <- data.frame()
for (i in 1:nrow(metadata)){
naid_df[i,1] <- metadata_id$file_name[i]
naid_df[i,2] <- metadata_id$associated_entities[i][[1]]$entity_submitter_id
}
colnames(naid_df) <- c("filename","TCGA_id")
3.批量讀取數(shù)據(jù)
##第三步钾菊，批量讀取數(shù)據(jù)
myfread <- function(files){
data = data.table::fread(paste0("data_in_one/",files))
data = data[-seq(1,4),]
data = data$unstranded
}
files <- dir("data_in_one")
system.time(test <- myfread(files[1]))測試讀取一個(gè)文件的所需時(shí)間
system.time(f <- lapply(files,myfread)) #lapply批量讀取
expr_df <- as.data.frame(do.call(cbind,f)) #完成后是個(gè)列表帅矗，需要轉(zhuǎn)換為數(shù)據(jù)框
當(dāng)前的數(shù)據(jù)是裸的滞详，沒有行凛俱，行應(yīng)該是基因，但是沒有料饥，列應(yīng)該是TCGA ID也沒有蒲犬，接下來添加一下。
先添加列名, 對(duì)應(yīng)關(guān)系之前已經(jīng)提取到了稀火，只是要限定以下行的順序暖哨，需要跟讀入的樣本順序一致
rownames(naid_df) <- naid_df[,1]
naid_df <- naid_df[files,]
colnames(expr_df) <- naid_df$TCGA_id
test=expr_df[1:20,1:20]
最后添加基因名稱，因?yàn)檫@些基因名稱和類型凰狞，在每一個(gè)樣本中都有篇裁，所以我們隨便讀取一個(gè)數(shù)據(jù)，然后添加即可赡若。
gene_id <- data.table::fread(paste0("data_in_one/",files[1]))$gene_id
gene_id <- gene_id[-seq(1,4)]
expr_df <- cbind(gene_id=gene_id,expr_df)
順利完成达布，該數(shù)據(jù)格式，第一列是基因逾冬，后面都是樣本黍聂，轉(zhuǎn)錄組下游分析躺苦，無論是公共數(shù)據(jù)，還是自測數(shù)據(jù)产还，都應(yīng)該獲取這樣的數(shù)據(jù)
總代碼·
### 1.從GDC下載counts數(shù)據(jù)
## https://portal.gdc.cancer.gov/repository
### 1.創(chuàng)建數(shù)據(jù)下載文件夾
if(!file.exists("rawdata")) dir.create("rawdata")
### 2.生成下載代碼
### 這個(gè)就是之前下載的manifest文件匹厘，每個(gè)人的不一樣，注意修改
manifest <- "gdc_manifest_20220929_023209.txt"
### 這個(gè)是下載的文件夾
rawdata <- "rawdata"
command <- sprintf("./gdc-client download -m %s -d %s",
manifest,
rawdata)
### 3.執(zhí)行下載操作
system(command = command)
#合并文件
if(!file.exists("data_in_one")) dir.create("data_in_one")
### 使用for循環(huán)來批量做脐区，回顧for循環(huán)的要點(diǎn)
for (dirname in dir("rawdata/")){
## 要查看的單個(gè)文件夾的絕對(duì)路徑
mydir <- paste0(getwd(),"/rawdata/",dirname)
##找到對(duì)應(yīng)文件夾中的文件并提取名稱愈诚，pattern表示模式，可以是正則表達(dá)式
file <- list.files(mydir,pattern = "_counts")
## 當(dāng)前文件的絕對(duì)路徑是
myfile <- paste0(mydir,"/",file)
#復(fù)制這個(gè)文件到目的文件夾
file.copy(myfile,"data_in_one")
}
#提取文件名稱和TCGAID 的對(duì)應(yīng)關(guān)系
metadata <- jsonlite::fromJSON("metadata.cart.2022-09-29.json")
library(dplyr)
metadata_id <- metadata %>%
dplyr::select(c(file_name,associated_entities))
## 提取對(duì)應(yīng)的文件
naid_df <- data.frame()
for (i in 1:nrow(metadata)){
naid_df[i,1] <- metadata_id$file_name[i]
naid_df[i,2] <- metadata_id$associated_entities[i][[1]]$entity_submitter_id
}
colnames(naid_df) <- c("filename","TCGA_id")
##第三步牛隅，批量讀取數(shù)據(jù)
myfread <- function(files){
data = data.table::fread(paste0("data_in_one/",files))
data = data[-seq(1,4),]
data = data$unstranded
}
files <- dir("data_in_one")
system.time(test <- myfread(files[1]))
system.time(f <- lapply(files,myfread))
expr_df <- as.data.frame(do.call(cbind,f))
test=expr_df[1:20,1:20]
rownames(naid_df) <- naid_df[,1]
naid_df <- naid_df[files,]
colnames(expr_df) <- naid_df$TCGA_id
test=expr_df[1:20,1:20]
gene_id <- data.table::fread(paste0("data_in_one/",files[1]))$gene_id
gene_id <- gene_id[-seq(1,4)]
expr_df <- cbind(gene_id=gene_id,expr_df)
save(expr_df,naid_df,file="exp.Rdata")
