1.Kaplan-Meier
library(survival)
library(survminer)
res.cut <- surv_cutpoint(svdata, time = "time", event = "event", variables = names(svdata)[3:ncol(svdata)], minprop = 0.3)#修改此处可以更改分组:minprop = 0.3
res.cat <- surv_categorize(res.cut)
my.surv <- Surv(res.cat$time, res.cat$event)
pl<-list()          
for (i in colnames(res.cat)[3:ncol(svdata)]) {
  group <- res.cat[,i] 
  survival_dat <- data.frame(group = group)
  fit <- survfit(my.surv ~ group)
  group <- factor(group, levels = c("low", "high"))
  data.survdiff <- survdiff(my.surv ~ group)
  p.val = 1 - pchisq(data.survdiff$chisq, length(data.survdiff$n) - 1)
  HR = (data.survdiff$obs[2]/data.survdiff$exp[2])/(data.survdiff$obs[1]/data.survdiff$exp[1])
  up95 = exp(log(HR) + qnorm(0.975)*sqrt(1/data.survdiff$exp[2]+1/data.survdiff$exp[1]))
  low95 = exp(log(HR) - qnorm(0.975)*sqrt(1/data.survdiff$exp[2]+1/data.survdiff$exp[1]))
  HR <- paste("Hazard Ratio = ", round(HR,2), sep = "")
  CI <- paste("95% CI: ", paste(round(low95,2), round(up95,2), sep = " - "), sep = "")
  svsort <- svdata[order(svdata[,i]),]
  pl[[i]]<-ggsurvplot(fit, data = survival_dat ,conf.int = F, censor = F, palette = c("#E7B800","#00AFBB"), legend.title = i,font.legend = 11,legend.labs=c(paste0(">",round(svsort[fit$n[2],i],2),"(",fit$n[1],")"),paste0("<",round(svsort[fit$n[2],i],2),"(",fit$n[2],")")),pval = paste(pval = ifelse(p.val < 0.001, "p < 0.001", paste("p = ",round(p.val,3), sep = "")),HR, CI, sep = "\n"))}
res <- arrange_ggsurvplots(pl, print = T,ncol = 3, nrow = 4)
ggsave(".pdf",res,width = 12,height = 16)

2. Forest
rm(list=ls())
inputFile="input.txt"        
outFile="forest.pdf"         

rt=read.table(inputFile,header=T,sep="\t",row.names=1,check.names=F)
gene=rownames(rt)
hr=sprintf("%.3f",rt$"HR")
hrLow=sprintf("%.3f",rt$"HR.95L")
hrHigh=sprintf("%.3f",rt$"HR.95H")
Hazard.ratio=paste0(hr,"(",hrLow,"-",hrHigh,")")
pVal=ifelse(rt$pvalue<0.001, "<0.001", sprintf("%.3f", rt$pvalue))

pdf(file=outFile, width = 6, height =4.5)
n=nrow(rt)
nRow=n+1
ylim=c(1,nRow)
layout(matrix(c(1,2),nc=2),width=c(3,2))

xlim = c(0,3)
par(mar=c(4,2,1.5,1.5))
plot(1,xlim=xlim,ylim=ylim,type="n",axes=F,xlab="",ylab="")
text.cex=0.8
text(0,n:1,gene,adj=0,cex=text.cex)
text(1.5-0.5*0.2,n:1,pVal,adj=1,cex=text.cex);text(1.5-0.5*0.2,n+1,'pvalue',cex=text.cex,font=2,adj=1)
text(3,n:1,Hazard.ratio,adj=1,cex=text.cex);text(3,n+1,'Hazard ratio',cex=text.cex,font=2,adj=1,)

par(mar=c(4,1,1.5,1),mgp=c(2,0.5,0))
xlim = c(0,max(as.numeric(hrLow),as.numeric(hrHigh)))
plot(1,xlim=xlim,ylim=ylim,type="n",axes=F,ylab="",xaxs="i",xlab="Hazard ratio")
arrows(as.numeric(hrLow),n:1,as.numeric(hrHigh),n:1,angle=90,code=3,length=0.03,col="darkblue",lwd=2.5)
abline(v=1,col="black",lty=2,lwd=2)
boxcolor = ifelse(as.numeric(hr) > 1, 'red', 'blue')
points(as.numeric(hr), n:1, pch = 15, col = boxcolor, cex=1.3)
axis(1)
dev.off()

3.ROC curve
#install.packages("pROC")
rm(list = ls())
library(pROC)             
inputFile="input.txt"      
outFile="ROC.pdf"         
x="RiskScore"
setwd("")   

rt=read.table(inputFile,header=T,sep="\t",check.names=F,row.names=1)
y=colnames(rt)[1]

rocobj1=roc(rt[,y], as.vector(rt[,x]))
pdf(file=outFile,width=5,height=5)
plot(rocobj1, print.auc=TRUE, col="red")
dev.off()

4.violin plot
rm(list = ls())
library(ggpubr)
setwd("")
data <- read.csv("Gide.csv")  
head(data)
my_comparisons <- list( c("LCK(R)", "LCK(NR)")) 
pdf("Gide.pdf", width = 5, height = 6)
ggviolin(data, x="svdata", y="TPM", fill = "svdata", 
         palette = c("#00AFBB", "#E7B800"), 
         add = "boxplot", add.params = list(fill="white"))+ 
  stat_compare_means(comparisons = my_comparisons)
dev.off()

5. Unsupervised Clustering
library(ConsensusClusterPlus)
tempdir <- function() ""
title = tempdir()
setwd("")
d <- read.table("input.txt", sep = "\t", header = T, row.names = 1)
d <- data.matrix(d)
dim(d)
### sweep
d = sweep(d, 1, apply(d, 1, median,na.rm=T))
maxK <- 3 
results = ConsensusClusterPlus(d, maxK = maxK, reps = 1000, pItem = 0.8, pFeature = 1,
                               title = title, clusterAlg = "hc", distance = "pearson", seed = 1262118388.71279, plot = "pdf")
icl = calcICL(results, title=title, plot = "pdf")