使用梯度函数 R 的数值 Hessian

Question

我需要使用梯度函数（由公式编程，而不是数字）以数值方式计算函数的Hessian矩阵。像或使用数值梯度计算粗麻布的包不能满足我的需求。我需要在包中内部执行的事情（没有可以调用的单独方法），但是处理我的优化任务的唯一方法在包中很好地实现并传递其最佳值以获得粗麻布对于我的程序来说太昂贵了。numDerivrootSolveoptimnloptoptim

所以我需要一些不使用数字梯度来计算 Hessian 的函数（例如，请参阅这些公式https://neos-guide.org/content/difference-approximations）。我无法制作这样的函数，因为我不明白如何选择h我的函数非常敏感的参数（增量）。我可以在 R 中找到这样的函数或以某种方式从optim包中检索它吗？或者有人至少可以解释如何选择误差最小化值，h然后我自己发布这个函数？

Answer 1

如果我理解正确，您应该使用numDeriv::jacobian()，它采用向量值函数并计算矩阵（每个元素相对于每个输入的导数），并将其应用于您的分析梯度函数。jacobian()确实使用数值近似（准确地说是理查森外推法），但我没有看到任何其他方法可以从黑盒梯度函数到 Hessian 函数？

您确实需要指定（或使用默认值）数字“delta”函数（默认为 1e-4）。另一方面，optim()用于计算 Hessian 矩阵的内部代码也使用有限差分：请参见此处和此处...

下面我定义了一个函数，它的梯度函数，以及它的Hessian；此代码显示它jacobian(grad(x))与 Hessian 相同（对于特定测试用例）。

library(numDeriv)
x1 <- c(1,1,1)

测试一下我没有搞砸梯度和 Hessian 函数：

all.equal(grad(f,x1),g(x1)) ## TRUE
all.equal(h(x1),hessian(f,x1)) ## TRUE

我的 Hessian 矩阵和梯度的 Jacobian 矩阵的数值等价：

all.equal(h(x1),jacobian(g,x1)) ## TRUE

测试功能：

f <- function(x) {
  sin(x[1])*exp(x[2])*x[3]^2
}

g <- function(x) {
  c(cos(x[1])*exp(x[2])*x[3]^2,
    sin(x[1])*exp(x[2])*x[3]^2,
    sin(x[1])*exp(x[2])*2*x[3])
}  

h <- function(x) {
    m <- matrix(NA,3,3)
    m[lower.tri(m,diag=TRUE)] <-
        c(-sin(x[1])*exp(x[2])*x[3]^2,
          cos(x[1])*exp(x[2])*x[3]^2,
          cos(x[1])*exp(x[2])*2*x[3],
    # col 2
           sin(x[1])*exp(x[2])*x[3]^2,
           sin(x[1])*exp(x[3])*2*x[3],
    # col 3
           sin(x[1])*exp(x[2])*2)
    m <- Matrix::forceSymmetric(m,"L")
    m <- unname(as.matrix(m))
    return(m)
}