2

我一直在尝试使用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？

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

4

1 回答 1

2

包的作者 Terry Therneaucoxme通过电子邮件回答——他让我把这个贴在这里。

下面是我对您的示例的重写，删除了 Surv 间接并在 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 参数）进行拟合来直接检查方差是否为零。如果随着方差达到 1e-6 或更小的值，可能性是恒定的或增加的，那么我假设 MLE 为零并从模型中删除该随机项。特里*

于 2018-09-06T08:00:30.800 回答