随机斜率Cox比例风险

Question

我一直试图使用coxme为我的模型中的每个协变量提取随机斜率。

library (coxme)
Start <- runif(5000, 1985, 2015)
Stop <- Start + runif(5000, 2, 10)
S <- data.frame (
  X1 <- runif(5000, 5.0, 7.5),
  X2 <- runif(5000, 5.0, 7.5),
  D <- rbinom(5000, 1, 0.8), 
  Letters <- sample(LETTERS, 5000, replace = TRUE), 
  Start <- Start,
  Stop <- Stop
  )

S_ind1 <- Surv (time = S$Start, time2 = S$Stop, event = S$D)
a <- coxme (S_ind1 ~ X1 + X2 + (X1 + X2|Letters), data = S)

我所得到的就是：

Error in gchol(kfun(theta, varlist, vparm, ntheta, ncoef)) : 
NA/NaN/Inf in foreign function call (arg 5)
In addition: Warning messages:
1: In sqrt(xvar * zvar) : NaNs produced
2: In sqrt(xvar * zvar) : NaNs produced

当我使用我自己的数据时，我经常会得到：

Error in coxme.fit(X, Y, strats, offset, init, control, weights = weights,  : 
  'Calloc' could not allocate memory (56076596 of 8 bytes)

是否可以使用coxme包括随机斜率？

如果没有，是否还有其他选择使用其他软件包？

Answer 1

coxme包裹的作者Therneau通过电子邮件回答--他让我把这个贴在这里。

下面是我对您的示例的重写，删除data.frame调用中的“=”，并在该上下文中使用“=”(我对<-在该上下文中工作有点惊讶)，并添加set.seed，以便该示例可以复制。

library (coxme)

set.seed(1953)
time1 <- runif(5000, 1985, 2015)
time2 <- time1 + runif(5000, 2, 10)
test <- data.frame (
  x1 = runif(5000, 5.0, 7.5),
  x2 = runif(5000, 5.0, 7.5),
  death = rbinom(5000, 1, 0.8), 
  letters = sample(LETTERS, 5000, replace = TRUE), 
  time1 = time1,
  time2 = time2)

fit1 <- coxme(Surv(time1, time2, death) ~ x1 + x2 + (1|letters), data=test)
fit2 <- coxme(Surv(time1, time2, death) ~ x1 + x2 + (1+x1 | letters), test)
fit3 <- coxme(Surv(time1, time2, death) ~ x1 + x2 + (1+x2 | letters), test)
fit4 <- coxme(Surv(time1, time2, death) ~ x1 + x2 + (1+ x1 + x2 | letters),
              data=test, vinit= c(1e-6, 1e-8, 1e-8))

*1.所有模型直到fit4为止。

1. 我觉得你的模型令人担忧，因为它有一个随机的斜率，但没有随机的拦截，就像所有通过原点的回归都让我担心:我很难去交换结果。虽然lme默认设置拦截条款，但coxme没有。
2. 我希望fit4能起作用，也许有更好的开始估计它会起作用。coxme的底层代码是我在所有生存工作中遇到的最困难的最大化问题，在某种意义上，最大化器很容易迷路，永远也找不到它的方向。这是一个功能，有时需要手握，通过有限的迭代计数和/或开始评估。我希望它不是这样，我有一些长期的计划来改善这一点，增加一个备用的基于MCMC的最大化器，这在理论上不会丢失，但牺牲了更长的计算时间。
3. 如果任何方差太接近于零，那么sqrt()消息就会作为舍入错误的函数出现。当然，在您的测试用例中，实际MLE的方差为0。当这种情况发生时，我通常会通过与固定变量序列(vfixed参数)相匹配来直接检查零方差。如果似然为常数或随着方差的增加而增加，则假设MLE为零，并将该随机项从模型中移除。特里·T.*