【问题标题】:inflection point for 3-parameter weibull in R?R中3参数weibull的拐点?
【发布时间】:2021-09-29 19:52:30
【问题描述】:

R 中的 drc 包包含自启动函数,用于将 3 参数 Weibull 模型拟合到数据。该软件包提供了 3 参数 weibull (https://cran.r-project.org/web/packages/drc/drc.pdf) 的 2 个参数化:

f(x) = 0 + (d − 0) exp(− exp(b(log(x) − log(e))))
f(x) = 0 + (d − 0)(1 − exp(− exp(b(log(x) − log(e)))))

包的文档表明d 代表上渐近线,b 是速率参数,并且曲线在剂量e 处具有拐点。但是,绘制这些函数并使用 R 来求解何时二阶导数 = 0 给出的值与拐点的 e 不同。例如:

dd=1
bb = -5
ee = 30
curve(dd*exp(-exp(bb * (log(x) - log(ee))) ) , xlim=c(0,100))
abline(v=ee)

# 1st deriv
g <- function(x) {}
body(g) <- D( expression(dd*exp(-exp(bb * (log(x) - log(ee))) )), "x")
curve(g, xlim=c(0, 100))
abline(v=ee)

# 2nd deriv
g <- function(x) {}
body(g) <- D(D( expression(dd*exp(-exp(bb * (log(x) - log(ee))) )), "x"), "x")
curve(g, xlim=c(0, 100))
abline(v=ee)

uniroot(g, c(20,50))  # should be 30, but is not?

文档不正确吗?或者这是否与 weibull 的特定记录参数化有关(将其限制为正数)?

【问题讨论】:

    标签: r nonlinear-functions weibull drc


    【解决方案1】:

    我认为文档有误。

    exp(-exp(x)) 在 0 处有一个拐点:

    DD <- function(expr, name, order = 1) {
      if(order < 1) stop("'order' must be >= 1")
      if(order == 1) D(expr, name)
      else DD(D(expr, name), name, order - 1)
    }
    
    DD( expression(exp(-exp(x))) , "x" , 2)
    #-(exp(-exp(x)) * exp(x) - exp(-exp(x)) * exp(x) * exp(x))
    

    exp(0) 等于 1。因此,我们有: -(1/e - 1*e) = 0

    我可以看到当内部指数的参数变为零时,如何假设exp(-exp(log(x) - log(c))) 也会有一个拐点。但事实并非如此,因为log 改变了曲率和二阶导数:

    DD( expression(exp(-exp(log(x) - log(c)))) , "x" , 2)
    #-(exp(-exp(log(x) - log(c))) * (exp(log(x) - log(c)) * (1/x) * 
    #    (1/x) - exp(log(x) - log(c)) * (1/x^2)) - exp(-exp(log(x) - 
    #    log(c))) * (exp(log(x) - log(c)) * (1/x)) * (exp(log(x) - 
    #    log(c)) * (1/x)))
    

    如果我们把 x = c,我们有:

    -(1/e * (1 * 1/c * 1/c - 1 * 1/c^2) - 1/e * 1 * 1/c * 1 * 1/c)

    = -(0 - 1/e * 1/c^2)

    = 1/e * 1/c^2

    那不是 0。

    只是为了检查我的数学:

    g <- function(x) {}
    body(g) <- DD( expression(exp(-exp(log(x) - log(c)))) , "x" , 2)
    c <- 2
    g(2)
    #[1] 0.09196986
    1/exp(1) * 1/2^2
    #[1] 0.09196986
    

    (为了获得额外的学分,您应该手动进行二阶导数。我没有那么多时间。)

    【讨论】:

    • 这是证实我怀疑的好方法。您是否清楚c 参数代表什么?它确实清楚地沿 x 轴移动了曲线。或者这可能是对 Weibull 函数的错误描述?
    猜你喜欢
    • 2019-04-30
    • 2014-06-23
    • 2012-08-02
    • 2011-09-12
    • 1970-01-01
    • 2017-12-25
    • 1970-01-01
    • 2020-10-29
    • 2013-06-23
    相关资源
    最近更新 更多