FlexSurv中的自定义生存函数

3ks5zfa0  于 2023-04-09  发布在  其他
关注(0)|答案(1)|浏览(145)

我想在flexsurv包中实现一个自定义生存函数,类似于第4节(第12页):
https://cran.r-project.org/web/packages/flexsurv/vignettes/flexsurv.pdf
假设我有一个 Dataframe df,其中status确定一个项目是否失败。1=failed0=non-failed

> library(flexsurv)
> library(survival)
> df
equip_id      age status
       1 20.33812      0
       2 28.44830      1
       3 22.47019      0
       4 47.21489      1
       5 30.42093      1

现在假设我有一个密度函数:

# must be named 'dcustom.weibull' to denote density function
dcustom.weibull <- function(shape, scale, t)
{
  f_t <- (shape/scale) * ((t/scale)^(shape - 1)) * exp(-(t/scale)^shape)
  return(f_t)
}

我尝试创建一个自定义函数,用于flexsurv包:

custom.weibull <- list(name = "custom.weibull", pars = c("shape", "scale"),
                       location = "scale",
                       transforms = c(log, log),
                       inv.transforms = c(exp, exp),
                       inits = function(t){ c(1, median(t)) })

我尝试用flexsurv包调用它:

flexsurvreg(Surv(age, status)~1, data = df, dist = custom.weibull)

我得到以下错误:

Forming integrated rmst function...
Forming integrated mean function...
Error in (function (shape, scale, t)  : 
  unused arguments (x = c(28.4482952329068, 47.2148908312751, 30.4209324295688), log = TRUE)

现在我意识到有很多软件包可以为2参数Weibull函数(例如fitdistplus)做这件事。我试图了解如何使用已知的密度函数设置它,以便最终实现一个不同的,不太传统的密度函数。

djmepvbi

djmepvbi1#

一般来说,我们必须为一个新的分布建立一个合适的函数族。我只看到了一个d-分布函数。此外,密度函数的“x”参数通常放在第一位,所以对我来说(R经验)眼睛,将其视为“t”并且处于第三位置是相当出乎意料的。所以我构建了一个自定义dweibull并添加了一个log参数,因为我的第一个版本抛出了一个错误,抱怨一个未使用的log参数。

library(flexsurv)
# it's poor practice to name a data object with the same token as an existing density fun
df1 <- read.table(text=" equip_id      age status
 1 20.33812      0
 2 28.44830      1
 3 22.47019      0
 4 47.21489      1
 5 30.42093      1", header=TRUE)
dcustom.weibull <- function(x, shape, scale, log=FALSE)
{
    f_x <- (shape/scale) * ((x/scale)^(shape - 1)) * exp(-(x/scale)^shape)
    return(f_x)
}
flexsurvreg(Surv(age, status)~1, data = df1, dist = custom.weibull)
#-------------------------
Forming cumulative distribution function...
Forming integrated rmst function...
Forming integrated mean function...
Call:
flexsurvreg(formula = Surv(age, status) ~ 1, data = df1, dist = custom.weibull)

Estimates: 
       est     L95%    U95%    se    
shape    4.07      NA      NA      NA
scale  196.43      NA      NA      NA

N = 5,  Events: 3,  Censored: 2
Total time at risk: 148.8924
Log-likelihood = 0.0001369537, df = 2
AIC = 3.999726

Warning message:
In flexsurvreg(Surv(age, status) ~ 1, data = df1, dist = custom.weibull) :
  Optimisation has probably not converged to the maximum likelihood - Hessian is not positive definite. 

So the density function was apparently sufficient without either a p- or q- function to complete the custom family, so the rewrite and addition of a log parameter was sufficient to cure the issue.

相关问题