數值的Hessian使用梯度函數R

Question

我需要計算使用我的梯度函數（由式，不是數字編程）我的函數數值的的Hessian。像numDeriv或rootSolve這樣的軟件包使用數字計算黑森州不符合我的需求的梯度。我需要執行的內容內部（我可以調用沒有單獨的方法）在optim包，但唯一的方法來處理我的優化任務良好實施在nlopt包和傳遞其最佳值爲optim爲了得到粗麻布太昂貴我的程序。

所以我需要一些功能，使用非數字梯度（見例如這些公式https://neos-guide.org/content/difference-approximations）計算黑森州。我不能做出這樣的功能，因爲我不明白如何選擇參數h（增量），我的功能非常敏感。我可以在R中找到這樣的函數嗎？或者從optim包中以某種方式找回它？或者可能有人至少解釋如何選擇最小化h的錯誤值，然後我會自己發佈這個函數？

Answer 1

如果我理解正確，您應該使用numDeriv::jacobian()，它採用矢量值函數並計算矩陣（每個元素相對於每個輸入的導數），並將其應用於您的解析梯度函數。 jacobian()確實使用數值近似（理查德森外推，準確而言），但我沒有看到任何其他方式可以從黑盒梯度函數轉換爲Hessian？

你需要指定（或使用默認值）數值「三角洲」功能（1E-4默認）。在另一方面，該optim()使用計算的Hessian內部代碼還使用有限差：見here和here ...

下面我定義一個函數，它的梯度的功能，其Hessian矩陣;此代碼顯示jacobian(grad(x))與Hessian（針對特定測試用例）相同。

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

測試，我沒有搞砸的梯度和黑森州功能：

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

我黑森州的數值等價和梯度的雅可比：

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) 
}