cox by R

1. 数据格式

   
   image.png
2. 基线表

clinTumor2 <- read.table("clinTumor2.txt",header = T,stringsAsFactors = F,sep = "\t")
## 频数和频率
#clinTumor2是表名
Char <- function(x){
  Count <- table(clinTumor2[,x])
  Prop <- round(prop.table(Count) * 100,1)
  Num <- paste0(Count,' (',Prop,')')
  
  Unit <- data.frame(Characteristics = c(x,names(Count)),
                     Number_of_Cases = c(NA,Num))
  return(Unit)
}
Char('vitalstatus')

CharNames <- colnames(clinTumor2)[-c(1,3,11,12)]  ##选择需要计算基线表信息的变量
Characteristics <- lapply(CharNames, Char)
Characteristics <- ldply(Characteristics)
write.csv(Characteristics,'Characteristics.csv',row.names = F,na = '')

3. 单因素cox回归

library(survival)
library(survminer)
library(plyr)

Sur2 <- Surv(time = clinTumor2$daystolastfollowup,event = clinTumor2$vitalstatus)   ##提取生存时间、生存状态
##将因子转换为数值
clinTumor2[,c(4:10)] <- lapply(clinTumor2[,c(4:10)], as.numeric)
##建立函数，用于lapply
UniCoxReg <- function(x){
  fml <- as.formula(paste0('Sur2~',x))
  uniMod <- coxph(fml,data = clinTumor2)
  SuniMod <- summary(uniMod)
  HR <- round(SuniMod$coefficients[,2],2)
  P_Value <- round(SuniMod$coefficients[,5],3)
  LCI <- round(SuniMod$conf.int[,3],2)
  UCI <- round(SuniMod$conf.int[,4],2)
  CI95 <- paste(LCI,'-',UCI)
  UniCox <- data.frame(Characteristics = x,
                       Hazard_Ratio = HR,
                       CI95 = CI95,
                       P_Value = P_Value)
  return(UniCox)
  
}
##提取需要算单因素cox的变量
UniNames <- colnames(clinTumor2)[c(1,4:7,9:12)]
Uni <- lapply(UniNames, UniCoxReg)
Uni <- ldply(Uni)

4. 多因素回归

##multicox
##所有变量都入 
MultiNames <- colnames(clinTumor2)[c(1,4:7,9:12)]  ##c里面是需要输入的变量

FML <- as.formula(paste0('Sur2~',  paste(MultiNames,collapse = '+')))
MultiMod <- coxph(FML,data = clinTumor2)
SMultiMod <- summary(MultiMod)
HR <- round(SMultiMod$coefficients[,2],2)
P_Value <- round(SMultiMod$coefficients[,5],3)
LCI <- round(SMultiMod$conf.int[,3],2)
UCI <- round(SMultiMod$conf.int[,4],2)
CI95 <- paste(LCI,'-',UCI)
Multi <- data.frame(Characteristics = MultiNames,
                    Hazard_Ratio = HR,
                    CI95 = CI95,
                    P_Value = P_Value)

##合并表格与单因素cox 
CoxReg <- merge(Uni,Multi,by = 'Characteristics',all = T)
write.csv(CoxReg,'CoxReg.csv',row.names = F,na = '')

5. km curve

library(survival)
clinTumor<-read.csv("eg.csv")
## 整理临床信息
clinTumor$daystolastfollowup <- as.numeric(clinTumor$daystolastfollowup)
clinTumor$vitalstatus <- as.numeric(clinTumor$vitalstatus)
median(clinTumor$gene)
clinTumor$group <- ifelse(clinTumor$gene < median(clinTumor$gene),'< Median','>= Median')

## 拟合模型
Sur <- Surv(time = clinTumor$daystolastfollowup,event = clinTumor$vitalstatus)
Mod <- survfit(Sur~group,data = clinTumor)
survdiff(Sur~group,data = clinTumor)

## 作图
plot(Mod,col = c('red','blue'),lty = 1,lwd = 2,
     xlab = 'Time',ylab = 'Survival',main = 'mRNAArray - CISD2')
legend('bottomleft',legend = c('< Median','>= Median','p = 0.4'),
       lty = 1,col = c('red','blue',0),cex = 0.7)

##加载数据
library(survival)

library(ggplot2)
library(survminer)
library(ggsci)
##构建生存模型
fit <- survfit(Surv(survival, status) ~ lpa_stage,
data = subset)
##绘图
#方式一
p<-ggsurvplot(fit, data = subset, risk.table = TRUE,tables.height = 0.2,tables.theme = theme_cleantable(), palette = c("lancet"), ggtheme = theme_bw())
#方式2
p<-ggsurvplot(fit, data = subset, surv.median.line = "hv",pval = TRUE,conf.int = TRUE,risk.table = TRUE,tables.height = 0.2,tables.theme = theme_cleantable(), palette = c("lancet"), ggtheme = theme_bw())
p

image.png