TCGA2.GDC数据整理

测试开头

今天是生信星球陪你的第509天

大神一句话，菜鸟跑半年。我不是大神，但我可以缩短你走弯路的半年~

就像歌儿唱的那样，如果你不知道该往哪儿走，就留在这学点生信好不好~

这里有豆豆和花花的学习历程，从新手到进阶，生信路上有你有我！

TCGA2.GDC数据整理

花花写于2020.1.2 今天去看牙医了，痛苦啊。但是不耽误学习。

本系列是我的TCGA学习记录，跟着生信技能树B站课程学的，已获得授权（嗯，真的^_^）。课程链接：https://www.bilibili.com/video/av49363776
目录：TCGA-1.数据下载

1.xml文件探索

使用R包XML：

library(XML)
result <- xmlParse("./clinical/00e01b05-d939-49e6-b808-e1bab0d6773b/nationwidechildrens.org_clinical.TCGA-J2-8192.xml")
rootnode <- xmlRoot(result)
rootsize <- xmlSize(rootnode)
print(rootnode[2])

探索发现，xml共有两个节点，第二个节点中储存着病人的信息。

xmldataframe <- xmlToDataFrame(rootnode[2])
t(xmlToDataFrame(rootnode[2]))
                   [,1]                                                
additional_studies ""                                                  
tumor_tissue_site  "Lung"                                              
histological_type  "Lung Adenocarcinoma- Not Otherwise Specified (NOS)"
other_dx           "Yes, History of Prior Malignancy"                  
gender             "FEMALE"                                            
vital_status       "Alive"

2.写循环，临床信息整理为数据框

xmls = dir("clinical/",pattern = "*.xml$",recursive = T)

td = function(x){
  result <- xmlParse(file.path("clinical/",x))
  rootnode <- xmlRoot(result)
  xmldataframe <- xmlToDataFrame(rootnode[2])
  return(t(xmldataframe))
}

cl = lapply(xmls,td)
cl_df <- t(do.call(cbind,cl))
 cl_df[1:3,1:3]
# additional_studies tumor_tissue_site histological_type                                   
# [1,] ""                 "Lung"            "Lung Adenocarcinoma- Not Otherwise Specified (NOS)"
# [2,] ""                 "Lung"            "Lung Adenocarcinoma- Not Otherwise Specified (NOS)"
# [3,] ""                 "Lung"            "Lung Adenocarcinoma- Not Otherwise Specified (NOS)"

3.探索miRNA信息

x = read.table(file = "mirna/00c46e8b-f303-4d20-bd6d-d650c36895f5/af29b644-b3a8-455b-9f07-b956d41f6ec4.mirbase21.mirnas.quantification.txt",header = T,sep = "t") 
head(x)
# miRNA_ID read_count reads_per_million_miRNA_mapped cross.mapped
# 1 hsa-let-7a-1      75848                     21277.9626            N
# 2 hsa-let-7a-2      75451                     21166.5905            Y
# 3 hsa-let-7a-3      75769                     21255.8004            N
# 4   hsa-let-7b      54045                     15161.4741            N
# 5   hsa-let-7c      23828                      6684.5704            Y
# 6   hsa-let-7d       1794                       503.2785            N

第一列和第二列有用。探索一下会发现，每个病人的第一列miRNA_ID内容和顺序都是一样的。

4.写循环合成表达矩阵

mis = dir("mirna/",pattern = "*mirnas.quantification.txt$",recursive = T)

ex = function(x){
  result <- read.table(file.path("mirna/",x),sep = "t",header = T)[,1:2]
  return(result)
}

mi = lapply(mis,ex)
mi_df <- t(do.call(cbind,mi))
mi_df[1:4,1:4]
# [,1]           [,2]           [,3]           [,4]        
# miRNA_ID   "hsa-let-7a-1" "hsa-let-7a-2" "hsa-let-7a-3" "hsa-let-7b"
# read_count "  43563"      "  44020"      "  43680"      "  65476"   
# miRNA_ID   "hsa-let-7a-1" "hsa-let-7a-2" "hsa-let-7a-3" "hsa-let-7b"
# read_count "  75848"      "  75451"      "  75769"      "  54045"   
colnames(mi_df) <- mi_df[1,]
#奇数列是多余的，只保留偶数列
mi_df <- mi_df[seq(2,nrow(mi_df),2),]
dim(mi_df)
#[1]  567 1881
mi_df[1:4,1:4]
# hsa-let-7a-1 hsa-let-7a-2 hsa-let-7a-3 hsa-let-7b
# read_count "  43563"    "  44020"    "  43680"    "  65476" 
# read_count "  75848"    "  75451"    "  75769"    "  54045" 
# read_count " 56189"     " 56527"     " 56437"     " 44829"  
# read_count "  7956"     "  7951"     "  7834"     " 15317" 
#转为数值型
mi_df <- apply(mi_df, 2, as.numeric)
mi_df[1:4,1:4]
# hsa-let-7a-1 hsa-let-7a-2 hsa-let-7a-3 hsa-let-7b
# [1,]        43563        44020        43680      65476
# [2,]        75848        75451        75769      54045
# [3,]        56189        56527        56437      44829
# [4,]         7956         7951         7834      15317