爾云間? 一個(gè)專門做科研的團(tuán)隊(duì)

原創(chuàng)?小果?生信果

經(jīng)常做癌癥數(shù)據(jù)預(yù)后分析的小伙伴肯定經(jīng)常接觸像單因數(shù)cox分析，多因數(shù)cox分析，LASSO-cox回歸分析這些分析方法，其中在單因數(shù)cox回歸篩選之后，往往還需要進(jìn)一步進(jìn)行篩選，而這個(gè)進(jìn)一步篩選的方法一般就是LASSO或者多因數(shù)逐步cox分析，根據(jù)小云這段時(shí)間分析數(shù)據(jù)的經(jīng)驗(yàn)，多因素逐步cox分析往往是起決定作用的。

?

下來讓我們一起來看看多因素逐步cox分析的具體代碼吧

library(survival)
library(survminer) #載入R包
rt=read.table("逐步.txt",header=T,sep="\t",check.names=F,row.names=1) #讀入數(shù)據(jù)，前三列分別為樣本名，時(shí)間，狀態(tài)，后面為各基因表達(dá)量。
multiCox=coxph(Surv(time, status) ~ ., data = rt) #計(jì)算與時(shí)間狀態(tài)的預(yù)后關(guān)系
multiCox=step(multiCox,direction = "both") #進(jìn)行逐步分析
#以下內(nèi)容為逐步分析顯示內(nèi)容
Version:1.0 StartHTML:0000000107 EndHTML:0000002932 StartFragment:0000000127 EndFragment:0000002914
Start:? AIC=837.78
Surv(time, status) ~ NXPE4 + TAGLN2 + PSAT1 + MAGEA1 + LAMP5 +
??? CD44 + CTAG2 + RCBTB2 + SLC40A1

????????? Df??? AIC
- NXPE4??? 1 835.78
- CTAG2??? 1 835.93
- MAGEA1?? 1 836.20
- PSAT1??? 1 836.34
- LAMP5??? 1 837.15
- RCBTB2?? 1 837.76
<none>?????? 837.78
- CD44???? 1 838.12
- SLC40A1? 1 840.83
- TAGLN2?? 1 849.15

Step:? AIC=835.78
Surv(time, status) ~ TAGLN2 + PSAT1 + MAGEA1 + LAMP5 + CD44 +
??? CTAG2 + RCBTB2 + SLC40A1

????????? Df??? AIC
- CTAG2??? 1 833.94
- MAGEA1?? 1 834.22
- PSAT1??? 1 834.35
- LAMP5??? 1 835.15
- RCBTB2?? 1 835.76
<none>?????? 835.78
- CD44???? 1 836.12
+ NXPE4??? 1 837.78
- SLC40A1? 1 838.86
- TAGLN2?? 1 848.53

Step:? AIC=833.94
Surv(time, status) ~ TAGLN2 + PSAT1 + MAGEA1 + LAMP5 + CD44 +
??? RCBTB2 + SLC40A1

????????? Df??? AIC
- PSAT1??? 1 832.53
- MAGEA1?? 1 832.84
- LAMP5??? 1 833.32
<none>?????? 833.94
- RCBTB2?? 1 833.99
- CD44???? 1 834.32
+ CTAG2??? 1 835.78
+ NXPE4??? 1 835.93
- SLC40A1? 1 837.30
- TAGLN2?? 1 846.71

Step:? AIC=832.53
Surv(time, status) ~ TAGLN2 + MAGEA1 + LAMP5 + CD44 + RCBTB2 +
??? SLC40A1

????????? Df??? AIC
- MAGEA1?? 1 831.37
- LAMP5??? 1 831.73
<none>?????? 832.53
- RCBTB2?? 1 832.62
- CD44???? 1 833.00
+ PSAT1??? 1 833.94
+ CTAG2??? 1 834.35
+ NXPE4??? 1 834.50
- SLC40A1? 1 835.86
- TAGLN2?? 1 844.72

Step:? AIC=831.37
Surv(time, status) ~ TAGLN2 + LAMP5 + CD44 + RCBTB2 + SLC40A1

????????? Df??? AIC
- LAMP5??? 1 830.80
<none>?????? 831.37
- RCBTB2?? 1 831.60
- CD44???? 1 832.16
+ MAGEA1?? 1 832.53
+ CTAG2??? 1 832.74
+ PSAT1??? 1 832.84
+ NXPE4??? 1 833.26
- SLC40A1? 1 834.76
- TAGLN2?? 1 844.23

Step:? AIC=830.8
Surv(time, status) ~ TAGLN2 + CD44 + RCBTB2 + SLC40A1

????????? Df??? AIC
<none>?????? 830.80
- RCBTB2?? 1 831.28
+ LAMP5??? 1 831.37
+ MAGEA1?? 1 831.73
- CD44???? 1 832.01
+ CTAG2??? 1 832.05
+ PSAT1??? 1 832.44
+ NXPE4??? 1 832.71
- SLC40A1? 1 836.28
- TAGLN2?? 1 843.0

multiCoxSum=summary(multiCox)

outTab=data.frame()
outTab=cbind(
? coef=multiCoxSum$coefficients[,"coef"],
? HR=multiCoxSum$conf.int[,"exp(coef)"],
? HR.95L=multiCoxSum$conf.int[,"lower .95"],
? HR.95H=multiCoxSum$conf.int[,"upper .95"],
? pvalue=multiCoxSum$coefficients[,"Pr(>|z|)"])
outTab=cbind(id=row.names(outTab),outTab) #整理結(jié)果
write.table(outTab,file="multiCox.xls",sep="\t",row.names=F,quote=F) #寫入數(shù)據(jù)

好了，今天的主要內(nèi)容就是這樣的，歡迎小伙伴們來和小云分享討論啊。