Rstudio——绘制生存曲线

# 设置PDF绘制路径并加载
pdfFile=paste(resultPath, "PINS_", dataset, ".pdf" ,sep="")
pdf(pdfFile)
# 加载绘制数据
file=paste(dataPath,dataset, "/", dataset, "_ProcessedData.RData" ,sep="")
load(file)
# 提取聚类结果
groups = result$groups
groups2=result$groups2
# 绘图

plot(result$dataTypeResult[[1]]$Discrepancy$AUC, ylab= "Area under the curve", xlab="Cluster number", main=paste("AUC of gene expression for ", dataset, " data", sep=""))
lines(1:Kmax, result$dataTypeResult[[1]]$Discrepancy$AUC)
points(result$dataTypeResult[[1]]$k, result$dataTypeResult[[1]]$Discrepancy$AUC[result$dataTypeResult[[1]]$k],col="red")

plot(result$dataTypeResult[[2]]$Discrepancy$AUC, ylab= "Area under the curve", xlab="Cluster number", main=paste("AUC of methylation for ", dataset, " data", sep=""))
lines(1:Kmax, result$dataTypeResult[[2]]$Discrepancy$AUC)
points(result$dataTypeResult[[2]]$k, result$dataTypeResult[[2]]$Discrepancy$AUC[result$dataTypeResult[[2]]$k],col="red")

plot(result$dataTypeResult[[3]]$Discrepancy$AUC, ylab= "Area under the curve", xlab="Cluster number", main=paste("AUC of microRNA for ", dataset, " data", sep=""))
lines(1:Kmax, result$dataTypeResult[[3]]$Discrepancy$AUC)
points(result$dataTypeResult[[3]]$k, result$dataTypeResult[[3]]$Discrepancy$AUC[result$dataTypeResult[[3]]$k],col="red")

coxFit <- coxph(Surv(time = Survival, event = Death) ~ as.factor(result$dataTypeResult[[1]]$groups), data = survival, ties="exact")
mfit <- survfit(Surv(Survival, Death == 1) ~ as.factor(result$dataTypeResult[[1]]$groups), data = survival)
plot(mfit, col=unique(result$dataTypeResult[[1]]$groups), main = paste("Survival curves for gene expression of ",dataset, " (PertCluster)", sep=""), xlab = "Days", ylab="Survival", lwd=2)
legend("topright", legend = paste("Cox p-value:", round(summary(coxFit)$sctest[3],digits = 5), sep=""))

coxFit <- coxph(Surv(time = Survival, event = Death) ~ as.factor(result$dataTypeResult[[2]]$groups), data = survival, ties="exact")
mfit <- survfit(Surv(Survival, Death == 1) ~ as.factor(result$dataTypeResult[[2]]$groups), data = survival)
plot(mfit, col=unique(result$dataTypeResult[[2]]$groups), main = paste("Survival curves for methylation of ", dataset, " (PertCluster)", sep=""), xlab = "Days", ylab="Survival", lwd=2)
legend("topright", legend = paste("Cox p-value:", round(summary(coxFit)$sctest[3],digits = 5), sep=""))

coxFit <- coxph(Surv(time = Survival, event = Death) ~ as.factor(result$dataTypeResult[[3]]$groups), data = survival, ties="exact")
mfit <- survfit(Surv(Survival, Death == 1) ~ as.factor(result$dataTypeResult[[3]]$groups), data = survival)
plot(mfit, col=unique(result$dataTypeResult[[3]]$groups), main = paste("Survival curves for microRNA of ", dataset, " (PertCluster)", sep=""), xlab = "Days", ylab="Survival", lwd=2)
legend("topright", legend = paste("Cox p-value:", round(summary(coxFit)$sctest[3],digits = 5), sep=""))

ageCol <- abs(as.numeric(clinical$"birth_days_to"))/365
names(ageCol) <- rownames(clinical)
age <- list()
for (j in levels(factor(groups))) {
age[[j]] <- ageCol[names(groups[groups==j])]
}
boxplot(age, main=paste("Age distribution, ", dataset, sep=""), xlab="Groups", ylab="Age")

coxFit <- coxph(Surv(time = Survival, event = Death) ~ as.factor(groups), data = survival[names(groups),], ties="exact")
mfit <- survfit(Surv(Survival, Death == 1) ~ as.factor(groups), data = survival[names(groups),])
plot(mfit, col=levels(factor(groups)), main = paste("Survival curves for ", dataset, ", level 1 (PertCluster)", sep=""), xlab = "Days", ylab="Survival", lwd=2)
legend("top", legend = paste("Cox p-value:", round(summary(coxFit)$sctest[3],digits = 5), sep=""))
legend("topright", fill=levels(factor(groups)), legend=paste("Group ",levels(factor(groups)), ": ", table(groups)[levels(factor(groups))], sep=""))

age <- list()
for (j in levels(factor(groups2))) {
age[[j]] <- ageCol[names(groups2[groups2==j])]
}
boxplot(age, main=paste("Age distribution, ", dataset, sep=""), xlab="Groups", ylab="Age")

coxFit <- coxph(Surv(time = Survival, event = Death) ~ as.factor(groups2), data = survival, ties="exact")
mfit <- survfit(Surv(Survival, Death == 1) ~ as.factor(groups2), data = survival)
a <-intersect(unique(groups2), unique(groups));names(a) <- intersect(unique(groups2), unique(groups)); a[setdiff(unique(groups2), unique(groups))] <- seq(setdiff(unique(groups2), unique(groups)))+max(groups)
colors <- a[levels(factor(groups2))]
plot(mfit, col=colors, main = paste("Survival curves for ", dataset, ", level 2 (PertCluster)", sep=""), xlab = "Days", ylab="Survival", lwd=2)
legend("top", legend = paste("Cox p-value:", round(summary(coxFit)$sctest[3],digits = 5), sep=""))
legend("topright", fill=colors, legend=paste("Group ",levels(factor(groups2)), ": ", table(groups2)[levels(factor(groups2))], sep=""))

# 关闭 pdf

dev.off()