定義生存分佈:: survreg（）

繼?survreg.distributions我定義我喜歡這個自定義分佈：

gamma <- list(name = 'gamma', 
      parms = c(2,2), 
      init = function(x, weights, ...){ 
      c(median(x), mad(x)) 
      }, 
      density = function(x, parms){ 
      shape <- parms[1] 
      scale <- parms[2] 
      cbind(pgamma(x, shape=shape, scale=scale), 
        1-pgamma(x, shape=shape, scale=scale), 
        dgamma(x, shape=shape, scale=scale), 
        (shape-1)/x - 1/scale, 
        (shape-1)*(shape-2)/x^2 - 2*(shape-1)/(x*scale) + 1/scale^2) 
      }, 
      quantile = function(p, parms) { 
      qgamma(p, shape=parms[1], scale=parms[2]) 
      }, 
      deviance = function(...) stop('deviance residuals not defined') 
)

但是我不能讓它運行：

require(survival) 
survreg(Surv(log(time), status) ~ ph.ecog + sex, lung, dist=gamma) 
#Error in coxph.wtest(t(x) %*% (wt * x), c((wt * eta + weights * deriv$dg) %*% : 
# NA/NaN/Inf in foreign function call (arg 3)

這個錯誤來自一些C代碼，但我認爲它生成得早得多...

任何提示/建議/替代survreg？

來源

2013-04-08 EDi

我看到兩種可能性。您可能會在停止呼叫時拋出錯誤，或者如果您將負數傳遞給log（），您會希望NaN的人知道，除非您提供數據？ – 2013-04-09 03:44:09

@DWin：謝謝，我會嘗試調試C代碼並檢查輸入。我的例子是可重複的，肺數據與生存包一起提供。日誌（時間）很好，都是正面的。 – EDi 2013-04-09 08:52:48

在Markmail中進行一些搜索我發現，survreg旨在用於具有位置參數的分佈，並且該gamma不在該族中。 http://markmail.org/search/?q=list%3Aorg.r-project.r-help++survreg+gamma+therneau#query:list%3Aorg.r-project.r-help%20%20survreg% 20gamma％20therneau +頁面：3 +中：t2pcdpq6p2pimpvl +狀態：結果 – 2013-04-09 12:51:32

我發現了flexsurv包實現了廣義伽馬分佈。

威布爾分佈從survreg和flexsurvreg估計都差不多（但要注意不同的參數化：

require(survival) 
summary(survreg(Surv(log(time), status) ~ ph.ecog + sex, data = lung, dist='weibull')) 

Call: 
survreg(formula = Surv(log(time), status) ~ ph.ecog + sex, data = lung, 
    dist = "weibull") 
       Value Std. Error  z   p 
(Intercept) 1.7504  0.0364 48.13 0.00e+00 
ph.ecog  -0.0660  0.0158 -4.17 3.10e-05 
sex   0.0763  0.0237 3.22 1.27e-03 
Log(scale) -1.9670  0.0639 -30.77 6.36e-208 

Scale= 0.14 

Weibull distribution 
Loglik(model)= -270.5 Loglik(intercept only)= -284.3 
    Chisq= 27.62 on 2 degrees of freedom, p= 1e-06 
Number of Newton-Raphson Iterations: 6 
n=227 (1 observation deleted due to missingness) 


require(flexsurv) 
flexsurvreg(Surv(log(time), status) ~ ph.ecog + sex, data = lung, dist='weibull') 

Call: 
flexsurvreg(formula = Surv(log(time), status) ~ ph.ecog + sex,  data = lung, dist = "weibull") 

Maximum likelihood estimates: 
      est L95% U95% 
shape 7.1500 6.3100 8.1000 
scale 5.7600 5.3600 6.1800 
ph.ecog -0.0660 -0.0970 -0.0349 
sex  0.0763 0.0299 0.1230 

N = 227, Events: 164, Censored: 63 
Total time at risk: 1232.1 
Log-likelihood = -270.5, df = 4 
AIC = 549

隨着flexsurvreg我們可以適應通用Gamma分佈這樣的數據：

flexsurvreg(Surv(log(time), status) ~ ph.ecog + sex, data = lung, dist='gengamma') 

Call: 
flexsurvreg(formula = Surv(log(time), status) ~ ph.ecog + sex,  data = lung, dist = "gengamma") 

Maximum likelihood estimates: 
      est L95% U95% 
mu  1.7800 1.7100 1.8600 
sigma 0.1180 0.0971 0.1440 
Q  1.4600 1.0200 1.9100 
ph.ecog -0.0559 -0.0853 -0.0266 
sex  0.0621 0.0178 0.1060 

N = 227, Events: 164, Censored: 63 
Total time at risk: 1232.1 
Log-likelihood = -267.57, df = 5 
AIC = 545.15

的loglogistic分佈是（與survreg不同），但可以很容易地定製（參見flexsurvreg的示例）

我還沒有測試太多，但flexsurv似乎是survival的一個很好的選擇。

來源

2013-04-09 15:37:50 EDi

感謝您發佈一個'倖存'問題的答案。 – 2014-06-26 05:33:02

定義生存分佈:: survreg（）

回答

相關問題