生存分析只需要tumor数据，不要normal，将其去掉，新表达矩阵数据命名为exprSet；GTEX可得到正常样本做补充（xena将GTEX和TCGA数据对应起来了）
clinical信息需要进一步整理，成为生存分析需要的格式，新临床信息数据命名为meta。
由于不同癌症的临床信息表格组织形式不同，这里的代码需要根据实际情况修改。
rm(list=ls())
options(stringsAsFactors = F)
load("TCGA-CHOL_gdc.Rdata")
library(stringr)
1.简化临床信息，选出要用的列
tmp = data.frame(colnames(clinical))#将clinical的列名组成一个数据框，以便于在阅览视图中找到它并复制
clinical = clinical[,c(
'bcr_patient_barcode',
'vital_status',
'days_to_death',
'days_to_last_followup',
'race_list',
'days_to_birth',
'gender' ,
'stage_event'
)]#只选关心的列
#其实days_to_last_followup应该是找age_at_initial_pathologic_diagnosis，这表格里没有，用days_to_birth计算一下年龄，暂且替代。
dim(clinical)
#> [1] 48 8
#rownames(clinical) <- NULL
rownames(clinical) <- clinical$bcr_patient_barcode#将clinical行名改为病人ID
clinical[1:4,1:4]
#> bcr_patient_barcode vital_status days_to_death
#> TCGA-W5-AA2Q TCGA-W5-AA2Q Alive
#> TCGA-W5-AA2O TCGA-W5-AA2O Dead 640
#> TCGA-W5-AA39 TCGA-W5-AA39 Dead 170
#> TCGA-W6-AA0S TCGA-W6-AA0S Alive
#> days_to_last_followup
#> TCGA-W5-AA2Q 50
#> TCGA-W5-AA2O
#> TCGA-W5-AA39
#> TCGA-W6-AA0S 352
meta = clinical
exprSet=exp[,Group=='tumor']#生存分析不需要正常样本
exprSet = log2(edgeR::cpm(exprSet)+1)#取cpm，后面用cpm来做生存分析
#简化meta的列名
colnames(meta)=c('ID','event','death','last_followup','race','age','gender','stage')
#空着的值改为NA
meta[meta==""]=NA
2.实现表达矩阵与临床信息的匹配
有的病人会有两个或两个以上的肿瘤样本，就有重复。两种可行的办法：
（1）以病人为中心，对表达矩阵的列按照病人ID去重复，每个病人只保留一个样本。（本文档）
（2）以样本为中心，把meta里的病人ID替换成样本ID，这样同一个病人的两个样本就会有两行完全一致的临床信息。（在zz.R中）
# 以病人为中心，表达矩阵按病人ID去重复
k = !duplicated(str_sub(colnames(exprSet),1,12));table(k)
#> k
#> TRUE
#> 36
exprSet = exprSet[,k]#按逻辑值不重复的取一下
#调整meta的ID顺序与exprSet列名一致
meta=meta[match(str_sub(colnames(exprSet),1,12),meta$ID),]
identical(meta$ID,str_sub(colnames(exprSet),1,12))
#> [1] TRUE
3.整理生存分析的输入数据
#1.1由随访时间和死亡时间计算生存时间(月)(自己考虑生存月份特别少的数据是否要去掉)
meta$time = ifelse(meta$event=="Alive",
meta$last_followup,
meta$death)
meta$time = as.numeric(meta$time)/30
#1.2 根据生死定义event，活着是0，死的是1
meta$event=ifelse(meta$event=='Alive',
0,
1)#注意Alive的大小写
table(meta$event)
#>
#> 0 1
#> 20 16
#1.3 年龄和年龄分组
meta$age=ceiling(abs(as.numeric(meta$age))/365)#ceiling向上取整
meta$age_group=ifelse(meta$age>median(meta$age,na.rm = T),'older','younger')#看是否大于中位数判断是older还是younger
table(meta$age_group)
#>
#> older younger
#> 18 18
#1.4 stage
library(stringr)
head(meta$stage)#版本不同，严谨的提取分期时应该按照TNM（T3N0M1部分）一个个的找
#> [1] "6thStage IVT3N0M1" "7thStage IIIT3N0M0" "7thStage IIT2bNXMX"
#> [4] "7thStage IT1N0M0" "7thStage IIT2N0M0" "7thStage IT1N0M0"
a = str_extract_all(meta$stage,"I|V");head(a)
#> [[1]]
#> [1] "I" "V"
#>
#> [[2]]
#> [1] "I" "I" "I"
#>
#> [[3]]
#> [1] "I" "I"
#>
#> [[4]]
#> [1] "I"
#>
#> [[5]]
#> [1] "I" "I"
#>
#> [[6]]
#> [1] "I"
b = sapply(a,paste,collapse = "");head(b)
#> [1] "IV" "III" "II" "I" "II" "I"
meta$stage = b
table(meta$stage)
#1.5 race
table(meta$race)
#>
#> ASIAN BLACK OR AFRICAN AMERICAN WHITE
#> 3 2 31
save(meta,exprSet,cancer_type,file = paste0(cancer_type,"_sur_model.Rdata"))
*生信技能树课程笔记
