GTEx數(shù)據(jù)庫獲取表達矩陣.tpm
一、下載數(shù)據(jù)
共要下載三個數(shù)據(jù),分別為表達矩陣、樣本信息、注釋信息
進入網(wǎng)站:UCSC Xena
點擊“Launch Xena”,選擇“DATA SETs”
點擊“GTEX(11 datasets)”
下載框中的兩個數(shù)據(jù),上面一個是表達矩陣,下面一個是樣本信息。還差一個注釋信息,下載地址:https://toil.xenahubs.net/download/probeMap/gencode.v23.annotation.gene.probemap?
需要注意的是:
表達矩陣中數(shù)據(jù)格式為log2(tpm+0.001)
下載完成后,三個文件的文件名分別為:
- gtex_RSEM_gene_tpm.gz
- GTEX_phenotype.gz
- gencode.v23.annotation.gene.probemap
二、載入數(shù)據(jù)
library(data.table) #載入數(shù)據(jù)用
#表達矩陣
exp_gtex.tpm=fread("gtex_RSEM_gene_tpm.gz",header = T, sep = '\t',data.table = F)
rownames(exp_gtex.tpm)=exp_gtex.tpm[,1]
exp_gtex.tpm=exp_gtex.tpm[,-1]
#樣本信息
data_cl=fread("GTEX_phenotype.gz",header = T, sep = '\t',data.table = F)
data_cl=data_cl[,c(1,3)]
names(data_cl)=c('Barcode','Tissue')
data_cl=data_cl[data_cl$Tissue == 'Prostate',] #篩選出Prostate的數(shù)據(jù)
#注釋信息
annotat=fread("gencode.v23.annotation.gene.probemap",header = T, sep = '\t',data.table = F)
annotat=annotat[,c(1,2)]
rownames(annotat)=annotat[,1] #這里沒有選擇刪去id這一列
?View(exp_gtex.tpm)
?View(data_cl)
?樣本信息中有122個barcode來自Prostate組織
View(annotat)
?三、處理數(shù)據(jù)
1 篩選出exp_gtex.tpm中的Prostate組織數(shù)據(jù),并還原為TPM
#篩選,篩選之后還剩100個barcode
exp_gtex.tpm=exp_gtex.tpm[,colnames(exp_gtex.tpm) %in% data_cl$Barcode]
#還原為TPM
exp_gtex.tpm=2^exp_gtex.tpm-0.001
2 基因注釋,去重復(fù)基因名,讀出表達矩陣
#基因注釋
exp_gtex.tpm=as.matrix(exp_gtex.tpm)
t_index=intersect(rownames(exp_gtex.tpm),rownames(annotat)) #行名取交集,t_index中是能夠進行注釋的probe_id
exp_gtex.tpm=exp_gtex.tpm[t_index,]
annotat=annotat[t_index,]
rownames(exp_gtex.tpm)=annotat$gene
#去除重復(fù)基因名
t_index1=order(rowMeans(exp_gtex.tpm),decreasing = T)
t_data_order=exp_gtex.tpm[t_index1,]
keep=!duplicated(rownames(t_data_order))#對于有重復(fù)的基因,保留第一次出現(xiàn)的那個,即行平均值大的那個
exp_gtex.tpm=t_data_order[keep,]#得到最后處理之后的表達譜矩陣
#讀出
write.csv(exp_gtex.tpm,file = "exp_gtex.tpm.csv",quote = FALSE)
?View(exp_gtex.tpm)
TCGA數(shù)據(jù)庫獲取表達矩陣.tpm
?TCGA_改版后STAR-count處理方法_老實人謝耳朵的博客-CSDN博客
result <- fromJSON(file = "E:/R/PRAD Data Mining/PRAD_data_mining/TCGA/Results/DESeq2差異分析/TP vs NT/GDCdata_star_count_TP&NT/metadata.cart.2022-05-01.json")
metadata <- data.frame(t(sapply(result,function(x){
id <- x$associated_entities[[1]]$entity_submitter_id
file_name <- x$file_name
all <- cbind(id,file_name)
})))
rownames(metadata) <- metadata[,2]
#獲取raw
t_dir <- 'E:/R/PRAD Data Mining/PRAD_data_mining/TCGA/Results/DESeq2差異分析/TP vs NT/GDCdata_star_count_TP&NT/all/'
t_samples=list.files(t_dir)
sampledir <- paste0(t_dir,t_samples) #各個文件路徑
example <- data.table::fread('E:/R/PRAD Data Mining/PRAD_data_mining/TCGA/Results/DESeq2差異分析/TP vs NT/GDCdata_star_count_TP&NT/all/005d2b9e-722c-40bd-aa5c-bd4e8842cb04.rna_seq.augmented_star_gene_counts.tsv',data.table = F)#讀入一個tsv文件,查看需要的列數(shù),“unstranded”
raw <- do.call(cbind,lapply(sampledir, function(x){
rt <- data.table::fread(x,data.table = F) #data.table::fread函數(shù)
rownames(rt) <- rt[,1]
rt <- rt[,7]###第7列為“tpm_unstranded”
}))
#替換行名、列名
colnames(raw)=sapply(strsplit(sampledir,'/'),'[',11)###列名,11為文件名005d2b9e-722c-40bd-aa5c-bd4e8842cb04.rna_seq.augmented_star_gene_counts.tsv
rownames(raw) <- example$gene_id ##行名
raw_t <- t(raw)
t_same <- intersect(row.names(metadata),row.names(raw_t))
dataPrep2 <- cbind(metadata[t_same,],raw_t[t_same,])
rownames(dataPrep2) <- dataPrep2[,1]
dataPrep2 <- t(dataPrep2)
dataPrep2 <-dataPrep2[-c(1:6),] #dataPrep2為未注釋count矩陣
#dataPrep2中數(shù)據(jù)類型為“character”,需要轉(zhuǎn)為“numeric”
puried_data=apply(dataPrep2,2,as.numeric)
#基因注釋
rownames(puried_data)=example[5:nrow(example),'gene_name']
#去除重復(fù)基因名
t_index=order(rowMeans(puried_data),decreasing = T)#計算所有行平均值,按降序排列
t_data_order=puried_data[t_index,]#調(diào)整表達譜的基因順序
keep=!duplicated(rownames(t_data_order))#對于有重復(fù)的基因,保留第一次出現(xiàn)的那個,即行平均值大的那個
exp_tcga.tpm=t_data_order[keep,]#得到最后處理之后的表達譜矩陣
write.csv(exp_tcga.tpm,file = "exp_tcga.tpm.csv",quote = FALSE)
?View(exp_tcga.tpm)文章來源:http://www.zghlxwxcb.cn/news/detail-416269.html
文章來源地址http://www.zghlxwxcb.cn/news/detail-416269.html
到了這里,關(guān)于TCGA_聯(lián)合GTEx分析1_得到表達矩陣.tpm的文章就介紹完了。如果您還想了解更多內(nèi)容,請在右上角搜索TOY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關(guān)文章,希望大家以后多多支持TOY模板網(wǎng)!