R语言：生成Hessian矩阵函数后，传入环境中表达式求值失败求助

问题分析 #

你在生成Hessian矩阵表达式并指定点求值时，遇到两个核心问题：

1. 最初的矩阵存储方案中，return语句放在内层循环里，导致函数提前终止，只返回第一行第一列结果；且未正确处理表达式的求值逻辑。
2. 使用环境存储表达式后，求值时错误使用do.call（该函数第一个参数应为函数/调用对象，而非环境变量名的字符串），同时变量赋值的方式也不符合R的环境操作规则。

正确实现方案 #

我们可以分两步完成：先生成存储二阶偏导表达式的矩阵，再编写函数将指定点的值代入计算数值Hessian矩阵。

1. 生成Hessian表达式矩阵

修正循环逻辑，直接将二阶偏导的表达式存入矩阵：

hessian_expr_matrix <- function(gx, respect_to) {
  n <- length(respect_to)
  # 创建空的表达式矩阵
  out_mat <- matrix(expression(), nrow = n, ncol = n)
  rownames(out_mat) <- respect_to
  colnames(out_mat) <- respect_to
  
  for(i in 1:n) {
    # 先对第i个变量求一阶偏导
    first_deriv <- D(gx, respect_to[i])
    for(j in 1:n) {
      # 再对第j个变量求二阶偏导
      second_deriv <- D(first_deriv, respect_to[j])
      out_mat[i, j] <- second_deriv
    }
  }
  return(out_mat)
}

# 测试生成表达式矩阵
g <- expression(x^3 - 2*x*y - y^6)
hess_expr <- hessian_expr_matrix(g, respect_to = c('x', 'y'))

此时hess_expr是存储二阶偏导表达式的矩阵，输出如下：

x           y          
x    6 * x       -2         
y    -2          -30 * y^5  

2. 在指定点求值Hessian矩阵

编写求值函数，将变量值代入表达式矩阵计算数值结果：

evaluate_hessian <- function(hess_expr, at_val) {
  # 创建临时环境，存入指定的变量值
  eval_env <- new.env()
  list2env(at_val, envir = eval_env)
  
  n <- nrow(hess_expr)
  num_mat <- matrix(0, nrow = n, ncol = n)
  rownames(num_mat) <- rownames(hess_expr)
  colnames(num_mat) <- colnames(hess_expr)
  
  # 遍历表达式矩阵，逐个求值
  for(i in 1:n) {
    for(j in 1:n) {
      num_mat[i, j] <- eval(hess_expr[i, j], envir = eval_env)
    }
  }
  return(num_mat)
}

# 测试在(x=1, y=2)处求值
evaluate_hessian(hess_expr, at_val = list(x=1, y=2))

输出结果：

x    y
x    6   -2
y   -2 -960

针对你原有环境方案的修正 #

如果坚持使用环境存储表达式，修正后的代码如下：

# 修正后的Hessian环境生成函数
hessian_env <- function(gx, respect_to) {
  out_env <- new.env()
  n <- length(respect_to)
  for(i in 1:n) {
    first_deriv <- D(gx, respect_to[i])
    for(j in 1:n) {
      second_deriv <- D(first_deriv, respect_to[j])
      assign(paste0(i,j), second_deriv, out_env)
    }
  }
  # 存入变量名，方便后续映射矩阵位置
  assign("var_names", respect_to, out_env)
  return(out_env)
}

# 修正后的求值函数
fisher_observed <- function(h_env, at_val) {
  # 创建包含变量值的环境，继承原环境的表达式
  eval_env <- new.env(parent = h_env)
  list2env(at_val, envir = eval_env)
  
  expr_names <- setdiff(ls(h_env), "var_names")
  n <- length(h_env$var_names)
  num_mat <- matrix(0, nrow = n, ncol = n)
  rownames(num_mat) <- h_env$var_names
  colnames(num_mat) <- h_env$var_names
  
  # 遍历每个表达式求值并填充矩阵
  for(name in expr_names) {
    expr <- get(name, envir = h_env)
    i <- as.integer(substr(name, 1, 1))
    j <- as.integer(substr(name, 2, 2))
    num_mat[i, j] <- eval(expr, envir = eval_env)
  }
  return(num_mat)
}

# 测试
h_g <- hessian_env(g, respect_to = c('x', 'y'))
fisher_observed(h_g, at_val = list(x=1, y=2))

同样会得到正确的数值Hessian矩阵。

关键知识点总结 #

• D()函数返回的是表达式，可直接存入矩阵或环境，使用eval()时需提供包含变量值的环境。
• list2env()能快速将列表形式的变量值存入环境，避免手动循环赋值。
• 循环中的return语句必须放在所有循环完成之后，否则会提前终止函数执行。

内容的提问来源于stack exchange，提问作者Casey Jayne