下面，我将R函数的结果与自己的代码进行比较。该算法仅包含最大化多个参数的函数（此处为19）。我的代码定义了功能，并使用nlm进行了优化。幸运的是，两者都返回相同的结果。但是，R功能非常快。因此，我怀疑我可以比使用nlm（或R中类似的优化例程）做得更好。任何想法？
这是一些survival data可以与Cox model配合的。为此，需要最大化部分对数似然性（Wikipedia链接中的第三个方程）。
在R中，这可以通过coxph()（survival包的一部分）完成：
> library(survival)
> fmla <- as.formula(paste("Surv(time, event) ~ ",
+ paste(names(data)[-(1:3)], collapse=" +")))
> mod <- coxph(formula=fmla, data=data)
> round(mod$coef, 3)
x1 x2 x3 x4 x5 x6 x7 x8 x9 x10 x11 x12 x13 x14 x15
-0.246 -0.760 0.089 -0.033 -0.138 -0.051 -0.484 -0.537 -0.620 -0.446 -0.204 -0.112 -0.089 -0.451 0.043
x16 x17 x18 x19
0.106 -0.015 -0.245 -0.653
可以通过显式编写部分对数可能性并使用一些数值优化例程来检查。这是完成此工作的一些原始代码。
该代码已根据我收到的评论进行了编辑
> #------ minus partial log-lik ------
> Mpll <- function(beta, data)
+ #!!!data must be ordered by increasing time!!!
+ #--> data <- data[order(data$time), ]
+ {
+ #preparation
+ N <- nrow(data)
+ linpred <- as.matrix(data[, -(1:3)]) %*% beta
+
+ #pll
+ pll <- sum(sapply(X=which(data$event == 1), FUN=function(j)
+ linpred[j] - log(sum(exp(linpred[j:N])))))
+
+ #output
+ return(- pll)
+ }
> #-----------------------------------
>
> data <- data[order(data$time), ]
> round(nlm(f=Mpll, p=rep(0, 19), data=data)$estimate, 3)
[1] -0.246 -0.760 0.089 -0.033 -0.138 -0.051 -0.484 -0.537 -0.620 -0.446 -0.204 -0.112 -0.089 -0.451
[15] 0.043 0.106 -0.015 -0.245 -0.653
好的，它可以工作...但是要慢得多！
是否有人对
coxph()中的处理方法有一个想法，以使其快速完成？ 最佳答案
这是代码的向量化版本。
Mpll2 <- function(beta, data) {
X <- as.matrix(data[, -(1:3)])
a <- X %*% beta
b <- log(rev(cumsum(rev(exp(a)))))
-sum((a - b)[data$event==1])
}
这是运行时间的简单测试。
data <- data[order(data$time), ] # No reason to order every time
# Yours
system.time(round(nlm(f=Mpll, p=rep(0, 19), data=data)$estimate, 3))
# user system elapsed
# 2.77 0.01 2.79
# Vectorized
system.time(round(nlm(f=Mpll2, p=rep(0, 19), data=data)$estimate, 3))
# user system elapsed
# 0.28 0.00 0.28
# Optimized C code
fmla <- as.formula(paste("Surv(time, event) ~ ",
paste(names(data)[-(1:3)], collapse=" +")))
system.time(round(coxph(formula=fmla, data=data)$coef,3))
# user system elapsed
# 0.02 0.00 0.03
因此，每种类型之间大约有一个数量级的差异。 C的速度非常快，并且您永远都不会接近R中的那些速度。但是C很难编写。
关于r - 加快手工制作Cox模型的拟合速度(相对于`survival::coxph`)，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/14153393/