1.问题起源
predict的参数type默认值就是“lp”，代表线性模型，理论上和大量文章里描述的公式应该是结果相同的。比如这个：
但是我突然发现他俩对不上号(⊙﹏⊙)
2.探索原因
没去搜索，直接发了个朋友圈，五分钟后收到了这个链接：https://mp.weixin.qq.com/s/Fg-Nefyz8BIzN17oJ23wxQ
哈哈！我神奇的几千人朋友圈，真的是藏龙卧虎，这也是我解决困惑的基本操作，近水楼台先得月。而且在这个链接里找到了我的名字~美滋滋。
简单来说，就是加权乘积比predict的预测结果多加了一个常数。
3.比较二者的结果
上文种已经对公式做了很好的推导，我用R语言完成一下咯：
rm(list = ls())
library(survival)
load("surv_model.Rdata")
names(model$coefficients)
## [1] "HFE" "SHISA5" "BRCA1" "EPM2A" "ERLIN1" "GPX7" "SLN"
## [8] "DNAJB11" "MMP2" "NOL3" "CP" "ATP2A2" "GLA" "MAPK3"
## [15] "SREBF2" "CASP2" "SNCA" "DDIT3" "BAG1" "ANXA5"
model$coefficients
## HFE SHISA5 BRCA1 EPM2A ERLIN1 GPX7 SLN
## 0.3024894 0.6928951 0.7909278 -0.5709679 -0.7810678 0.4341572 -0.1456389
## DNAJB11 MMP2 NOL3 CP ATP2A2 GLA MAPK3
## -0.8620475 0.1276739 0.3153988 0.1535990 0.6495272 0.3665937 0.8236297
## SREBF2 CASP2 SNCA DDIT3 BAG1 ANXA5
## -0.4178062 0.3980908 0.1977293 0.1403852 -0.2884173 0.1684487
a = predict(model,surv1)
b = apply(surv1[,names(model$coefficients)], 1,function(k){
sum(model$coefficients * k)
})
dat = data.frame(a,b,k = a-b)
head(dat)
## a b k
## TCGA-06-0878-01A 2.232431 12.03249 -9.800062
## TCGA-26-5135-01A 2.736399 12.53646 -9.800062
## TCGA-06-5859-01A 1.992062 11.79212 -9.800062
## TCGA-06-2563-01A 1.153412 10.95347 -9.800062
## TCGA-41-2571-01A 2.982385 12.78245 -9.800062
## TCGA-28-5207-01A 1.966414 11.76648 -9.800062
plot(a,b)
a是predict计算的结果，b是加权乘积计算的结果。如果是lasso模型，这两个数值是一致的。但cox模型，确实不一致。
4.不擅长公式推导，勉强试试
需要认识一下R语言里的函数exp和log
exp(1)
## [1] 2.718282
exp(10)
## [1] 22026.47
log(exp(1))
## [1] 1
log(exp(10))
## [1] 10
exp是自然常数e的指数运算，exp(n)就是e的n次方，log则是数学公式里的ln，取自然对数的意思，二者刚好可逆
所以公式可以写作
h(t,x) = h0(t)*exp(b)
log(h(t,x))-log(h0(t)) = b
所以predict计算出的预测值a就是公式里的log(h(t,x))，而上面的a-b的差-9.800062这个数字，就是公式里的常数log(h0(t))咯。
