如何模拟审查数据

机器算法验证 r 生存 模拟 随机生成
2022-02-27 05:00:21

我想知道如何模拟包含 I 型右删失观察的 n 个 Weibull 分布寿命的样本。例如,让 n = 3,shape = 3,scale = 1,审查率 = .15,审查时间 = .88。我知道如何生成 Weibull 样本,但我不知道如何生成 R 中类型 I 右删失的删失数据。

T = rweibull(3, shape=.5, scale=1)
2个回答

(作为 R 编码风格问题,最好不要 T 用作变量名,因为它是 的别名 TRUE,那样做难免会出问题。)


你的问题有点模棱两可;有几种解释它的方法。让我们来看看它们:

  1. 您规定要模拟类型 1 审查这通常被认为意味着实验运行了一段时间,并且到那时任何没有发生过事件的研究单位都会被审查。如果这就是您的意思,那么(必然)不可能同时规定形状和比例参数,以及审查时间和速率。已经规定了任何三个,最后一个必然是固定的。

    (尝试)求解形状参数:
    失败;无论形状参数是什么,似乎不可能在 0.88 的审查时间与韦布尔分布的尺度参数保持为 1 的情况下获得 15% 的审查率。

    optim(.5, fn=function(shp){(pweibull(.88, shape=shp, scale=1, lower.tail=F)-.15)^2})
    # $par
    # [1] 4.768372e-08
    # ...
    # There were 46 warnings (use warnings() to see them)
    pweibull(.88, shape=4.768372e-08, scale=1, lower.tail=F)
    # [1] 0.3678794
    
    optim(.5, fn=function(shp){(pweibull(.88, shape=shp, scale=1, lower.tail=F)-.15)^2},
          control=list(reltol=1e-16))
    # $par
    # [1] 9.769963e-16
    # ...
    # There were 50 or more warnings (use warnings() to see the first 50)
    pweibull(.88, shape=9.769963e-16, scale=1, lower.tail=F)
    # [1] 0.3678794
    

    求解尺度参数:

    optim(1, fn=function(scl){(pweibull(.88, shape=.5, scale=scl, lower.tail=F)-.15)^2})
    # $par
    # [1] 0.2445312
    # ...
    pweibull(.88, shape=.5, scale=0.2445312, lower.tail=F)
    # [1] 0.1500135
    

    求解审查时间:

    qweibull(.15, shape=.5, scale=1, lower.tail=F)
    # [1] 3.599064
    

    求解审查率:

    pweibull(.88, shape=.5, scale=1, lower.tail=F)
    # [1] 0.3913773
    
  2. 另一方面,我们可以将审查视为在整个研究过程中随机(通常是独立地)发生的,例如由于辍学。在这种情况下,程序是模拟两组 Weibull 变量。然后您只需注意哪个先出现:您使用较小的值作为端点,如果较小的值是审查时间,则将该单位称为审查单位。例如:

    set.seed(0775)  
    t    = rweibull(3, shape=.5, scale=1)
    t      # [1] 0.7433678 1.1325749 0.2784812
    c    = rweibull(3, shape=.5, scale=1.5)
    c      # [1] 3.3242417 2.8866217 0.9779436
    time = pmin(t, c)
    time   # [1] 0.7433678 1.1325749 0.2784812
    cens = ifelse(c<t, 1, 0)
    cens   # [1] 0 0 0
    

只是为了确保我们谈论的是同一件事,I 型审查是什么时候

...一个实验有一定数量的受试者或项目,并在预定时间停止实验,此时剩余的任何受试者都将被右删失。

要使用censoring time = 0.88生成正确的删失数据,您只需使用以下min函数:

T <- rweibull(3, shape=.5, scale=1)
censoring_time <- 0.88
T_censored <- min(censoring_time, T)

但是,我不完全确定您所说的“审查率 = 0.15 ”是什么意思……您的意思是说您的科目中有 15% 是正确审查的吗?这些关于审查的说明似乎表明I 型审查所需的唯一参数是审查时间,所以我不确定这个比率是如何影响的。