思思热思思,亚州无码第9

經(jīng)過悉心的研究以及很多老師，朋友的幫助，順利的進行了自己預想的分析，但仍有寫小細節(jié)仍有缺陷。特別感謝genefu 包的作者Benjamin Haibe-Kains教授，謹以本文記錄自己的探索過程，如可對后來者有幫助也算是功德一件了。
####step18-與基于芯片的模型比較(TCGA cohort，AJCC stage)
####比較AUC值，以ROC曲線形式展示
##載入已準備好的數(shù)據(jù)
## code used to compute the risk score
library(genefu)
library(survival)
library(survivalROC)
library(survminer)

##先準備注釋文件
ddata<-eset##準備好col為symbol的表達矩陣
##
s<-colnames(ddata)
##id 轉(zhuǎn)換
library(org.Hs.eg.db)
s2g=toTable(org.Hs.egSYMBOL)
g=s2g[match(s,s2g$symbol),1];
head(g)
##構(gòu)造注釋文件
#  probe Gene.symbol Gene.ID
dannot=data.frame(probe=s,
                  "Gene.Symbol" =s, 
                  "EntrezGene.ID"=g)
dannot[1:3,1:3]
##保留已注釋的基因及注釋文件
ddata=ddata[,!is.na(dannot$EntrezGene.ID)]
dannot=dannot[!is.na(dannot$EntrezGene.ID),] 
head(dannot)

##如已準備好數(shù)據(jù)
##my data 
ddata[1:5,1:5]
dim(ddata)
dannot[1:5,1:3]
dim(dannot)
all(colnames(ddata)==dannot$probe)
##
rownames(dannot) <- as.character(dannot[,'probe'])
data(sig.oncotypedx)
sig.oncotypedx
rs.onco <- oncotypedx(data=ddata,annot=dannot, do.mapping=TRUE)
table(rs.onco$risk)##

##計算OncotypedROC曲線
test_expr <- exprSet
new_exprSet <- test_expr[,c('sampleID','OS_time','OS_event')]
new_exprSet$multi_model <- rs.onco$score

##
cutoff3 =60##5年
##oncotypedx風險分數(shù)
Mayo.fit1= survivalROC(Stime=new_exprSet$OS_time,
                     status=new_exprSet$OS_event,
                     marker =new_exprSet$multi_model,
                     predict.time = cutoff3,method="KM")
pdf("TCGA_oncotypedx_ROC_OS.pdf")
plot(Mayo.fit1$FP,Mayo.fit1$TP,type="l",
     xlim=c(0,1),ylim=c(0,1),
     xlab="1-Specificity",ylab="Sensitivity",main="TCGA cohort ",col="#1c61b6",lwd=2)
abline(0,1,col="grey",lwd=2)

#my gene signature
##按自己的模型計算
model_predcit_values <- predict(Multi_cox3,test_expr[,multi_sign_gene],type="risk")
new_exprSet <- test_expr[,c('sampleID','OS_time','OS_event')]
new_exprSet$multi_model <- model_predcit_values
Mayo.fit3= survivalROC(Stime=new_exprSet$OS_time,
                       status=new_exprSet$OS_event,
                       marker =new_exprSet$multi_model,
                       predict.time = cutoff3,method="KM")
lines(Mayo.fit3$FP,Mayo.fit3$TP,type="l",
      xlim=c(0,1),ylim=c(0,1),col="#840000",lwd=2)

##添加圖注
legend(0.25,0.25,
       legend=c(paste("Gene21 AUC",round(Mayo.fit1$AUC,3)),
                 paste("Gene5 AUC",round(Mayo.fit3$AUC,3))),
       title.adj=1,##圖例標題字體位置1為右
       bty='n',##不畫出圖例框
       lwd=2,col=c("#1c61b6","#840000"))
dev.off()

####################################################################3
####RFS模型比較
##計算OncotypedROC曲線
test_expr <- exprSet
new_exprSet <- test_expr[,c('sampleID','RFS_time','RFS_event')]
new_exprSet$multi_model <- rs.onco$score

##
cutoff3 =60##5年
##oncotypedx風險分數(shù)
Mayo.fit1= survivalROC(Stime=new_exprSet$RFS_time,
                       status=new_exprSet$RFS_event,
                       marker =new_exprSet$multi_model,
                       predict.time = cutoff3,method="KM")
pdf("TCGA_oncotypedx_ROC_RFS.pdf")
plot(Mayo.fit1$FP,Mayo.fit1$TP,type="l",
     xlim=c(0,1),ylim=c(0,1),
     xlab="1-Specificity",ylab="Sensitivity",main="TCGA cohort RFS ",col="#1c61b6",lwd=2)
abline(0,1,col="grey",lwd=2)

#my gene signature
##按自己的模型計算
model_predcit_values <- predict(Multi_cox3,test_expr[,multi_sign_gene],type="risk")
new_exprSet <- test_expr[,c('sampleID','RFS_time','RFS_event')]
new_exprSet$multi_model <- model_predcit_values
Mayo.fit3= survivalROC(Stime=new_exprSet$RFS_time,
                       status=new_exprSet$RFS_event,
                       marker =new_exprSet$multi_model,
                       predict.time = cutoff3,method="KM")
lines(Mayo.fit3$FP,Mayo.fit3$TP,type="l",
      xlim=c(0,1),ylim=c(0,1),col="#840000",lwd=2)

##添加圖注
legend(0.25,0.25,
       legend=c(paste("Gene21 AUC",round(Mayo.fit1$AUC,3)),
                paste("Gene5 AUC",round(Mayo.fit3$AUC,3))),
       title.adj=1,##圖例標題字體位置1為右
       bty='n',##不畫出圖例框
       lwd=2,col=c("#1c61b6","#840000"))
##
dev.off()
#####################################################################
##early stage(I/II) ROC comparision
##提取早期組數(shù)據(jù)
test_expr <- exprSet
dim(test_expr)
test_expr[1:5,1:5]
table(test_expr$stage)##I/II期589sample
table(test_expr$Stage)##

##提取stage為I/II期的sample
early_index<-which(test_expr$stage=='0')
test_expr<-test_expr[early_index,]
dim(test_expr)
new_exprSet <- test_expr[,c('sampleID','OS_time','OS_event')]

##oncotypedx風險分數(shù)
ddata[1:5,1:5]
data(sig.oncotypedx)
sig.oncotypedx
##注意sample矩陣的改變
ddata<-ddata[early_index,]
dim(ddata)
rs.onco <- oncotypedx(data=ddata,annot=dannot, do.mapping=TRUE)
table(rs.onco$risk)##
new_exprSet$multi_model <- rs.onco$score
##
cutoff3 =60##5年

Mayo.fit1= survivalROC(Stime=new_exprSet$OS_time,
                       status=new_exprSet$OS_event,
                       marker =new_exprSet$multi_model,
                       predict.time = cutoff3,method="KM")
pdf("Early_stage_oncotypedx_ROC_OS.pdf",width = 5,height = 5)
plot(Mayo.fit1$FP,Mayo.fit1$TP,type="l",
     xlim=c(0,1),ylim=c(0,1),
     xlab="1-Specificity",ylab="Sensitivity",main="Early stage cohort ",col="#1c61b6",lwd=2)
abline(0,1,col="grey",lwd=2)

#my gene signature
##按自己的模型計算
model_predcit_values <- predict(Multi_cox3,test_expr[,multi_sign_gene],type="risk")
new_exprSet <- test_expr[,c('sampleID','OS_time','OS_event')]
new_exprSet$multi_model <- model_predcit_values
Mayo.fit3= survivalROC(Stime=new_exprSet$OS_time,
                       status=new_exprSet$OS_event,
                       marker =new_exprSet$multi_model,
                       predict.time = cutoff3,method="KM")
lines(Mayo.fit3$FP,Mayo.fit3$TP,type="l",
      xlim=c(0,1),ylim=c(0,1),col="#840000",lwd=2)

##添加圖注
legend(0.25,0.25,
       legend=c(paste("Gene21 AUC",round(Mayo.fit1$AUC,3)),
                paste("Gene5 AUC",round(Mayo.fit3$AUC,3))),
       title.adj=1,##圖例標題字體位置1為右
       bty='n',##不畫出圖例框
       lwd=2,col=c("#1c61b6","#840000"))
dev.off()

###########################################################
####Earlt stage RFS比較模型
new_exprSet <- test_expr[,c('sampleID','RFS_time','RFS_event')]
new_exprSet$multi_model <- rs.onco$score

##
cutoff3 =60##5年
##oncotypedx風險分數(shù)
Mayo.fit1= survivalROC(Stime=new_exprSet$RFS_time,
                       status=new_exprSet$RFS_event,
                       marker =new_exprSet$multi_model,
                       predict.time = cutoff3,method="KM")
pdf("Early_oncotypedx_ROC_RFS.pdf")
plot(Mayo.fit1$FP,Mayo.fit1$TP,type="l",
     xlim=c(0,1),ylim=c(0,1),
     xlab="1-Specificity",ylab="Sensitivity",main="Early stage cohort RFS ",col="#1c61b6",lwd=2)
abline(0,1,col="grey",lwd=2)

#my gene signature
##按自己的模型計算
model_predcit_values <- predict(Multi_cox3,test_expr[,multi_sign_gene],type="risk")
new_exprSet <- test_expr[,c('sampleID','RFS_time','RFS_event')]
new_exprSet$multi_model <- model_predcit_values
Mayo.fit3= survivalROC(Stime=new_exprSet$RFS_time,
                       status=new_exprSet$RFS_event,
                       marker =new_exprSet$multi_model,
                       predict.time = cutoff3,method="KM")
lines(Mayo.fit3$FP,Mayo.fit3$TP,type="l",
      xlim=c(0,1),ylim=c(0,1),col="#840000",lwd=2)

##添加圖注
legend(0.25,0.25,
       legend=c(paste("Gene21 AUC",round(Mayo.fit1$AUC,3)),
                paste("Gene5 AUC",round(Mayo.fit3$AUC,3))),
       title.adj=1,##圖例標題字體位置1為右
       bty='n',##不畫出圖例框
       lwd=2,col=c("#1c61b6","#840000"))
##
dev.off()

###############################################################
##Advanced stage ROC comparision
test_expr <- exprSet
dim(test_expr)
test_expr[1:5,1:5]
colnames(test_expr)[1:60]
table(test_expr$stage)
table(test_expr$Stage)
##late 168

##提取III IV期sample 168例
early_index<-which(test_expr$stage=='1')
test_expr<-test_expr[early_index,]
dim(test_expr)

##提取生存數(shù)據(jù)
new_exprSet <- test_expr[,c('sampleID','OS_time','OS_event')]
##oncotypedx風險分數(shù)
ddata[1:5,1:5]
data(sig.oncotypedx)
sig.oncotypedx
##注意sample矩陣的改變
ddata<-ddata[early_index,]
dim(ddata)
rs.onco <- oncotypedx(data=ddata,annot=dannot, do.mapping=TRUE)
table(rs.onco$risk)##
new_exprSet$multi_model <- rs.onco$score
##
cutoff3 =60##5年

Mayo.fit1= survivalROC(Stime=new_exprSet$OS_time,
                       status=new_exprSet$OS_event,
                       marker =new_exprSet$multi_model,
                       predict.time = cutoff3,method="KM")
pdf("Advanced_stage_oncotypedx_ROC_OS.pdf",width = 5,height = 5)
plot(Mayo.fit1$FP,Mayo.fit1$TP,type="l",
     xlim=c(0,1),ylim=c(0,1),
     xlab="1-Specificity",ylab="Sensitivity",main="Advanced stage cohort ",col="#1c61b6",lwd=2)
abline(0,1,col="grey",lwd=2)

#my gene signature
##按自己的模型計算
model_predcit_values <- predict(Multi_cox3,test_expr[,multi_sign_gene],type="risk")
new_exprSet <- test_expr[,c('sampleID','OS_time','OS_event')]
new_exprSet$multi_model <- model_predcit_values
Mayo.fit3= survivalROC(Stime=new_exprSet$OS_time,
                       status=new_exprSet$OS_event,
                       marker =new_exprSet$multi_model,
                       predict.time = cutoff3,method="KM")
lines(Mayo.fit3$FP,Mayo.fit3$TP,type="l",
      xlim=c(0,1),ylim=c(0,1),col="#840000",lwd=2)

##添加圖注
legend(0.25,0.25,
       legend=c(paste("Gene21 AUC",round(Mayo.fit1$AUC,3)),
                paste("Gene5 AUC",round(Mayo.fit3$AUC,3))),
       title.adj=1,##圖例標題字體位置1為右
       bty='n',##不畫出圖例框
       lwd=2,col=c("#1c61b6","#840000"))
dev.off()

###########################################################
####Earlt stage RFS比較模型
new_exprSet <- test_expr[,c('sampleID','RFS_time','RFS_event')]
new_exprSet$multi_model <- rs.onco$score

##
cutoff3 =60##5年
##oncotypedx風險分數(shù)
Mayo.fit1= survivalROC(Stime=new_exprSet$RFS_time,
                       status=new_exprSet$RFS_event,
                       marker =new_exprSet$multi_model,
                       predict.time = cutoff3,method="KM")
pdf("Advanced_oncotypedx_ROC_RFS.pdf",width = 5,height = 5)
plot(Mayo.fit1$FP,Mayo.fit1$TP,type="l",
     xlim=c(0,1),ylim=c(0,1),
     xlab="1-Specificity",ylab="Sensitivity",main="Advanced stage cohort RFS ",col="#1c61b6",lwd=2)
abline(0,1,col="grey",lwd=2)

#my gene signature
##按自己的模型計算
model_predcit_values <- predict(Multi_cox3,test_expr[,multi_sign_gene],type="risk")
new_exprSet <- test_expr[,c('sampleID','RFS_time','RFS_event')]
new_exprSet$multi_model <- model_predcit_values
Mayo.fit3= survivalROC(Stime=new_exprSet$RFS_time,
                       status=new_exprSet$RFS_event,
                       marker =new_exprSet$multi_model,
                       predict.time = cutoff3,method="KM")
lines(Mayo.fit3$FP,Mayo.fit3$TP,type="l",
      xlim=c(0,1),ylim=c(0,1),col="#840000",lwd=2)

##添加圖注
legend(0.25,0.25,
       legend=c(paste("Gene21 AUC",round(Mayo.fit1$AUC,3)),
                paste("Gene5 AUC",round(Mayo.fit3$AUC,3))),
       title.adj=1,##圖例標題字體位置1為右
       bty='n',##不畫出圖例框
       lwd=2,col=c("#1c61b6","#840000"))
##
dev.off()
色偷偷精品伊人,欧洲久久精品,欧美综合婷婷骚逼,国产AV主播,国产最新探花在线,九色在线视频一区,伊人大交九欧美,1769亚洲,黄色成人av

genefu package多基因預后模型的比較

genefu package多基因預后模型的比較

相關閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容

色偷偷精品伊人,欧洲久久精品,欧美综合婷婷骚逼,国产AV主播,国产最新探花在线,九色在线视频一区,伊人大交九 欧美,1769亚洲,黄色成人av

genefu package多基因預后模型的比較

相關閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容

色偷偷精品伊人,欧洲久久精品,欧美综合婷婷骚逼,国产AV主播,国产最新探花在线,九色在线视频一区,伊人大交九欧美,1769亚洲,黄色成人av
