使用R中solnp函数带约束优化双参数函数时遇错误求助

先来说说你遇到的两个核心问题：线性化无可行解和海森矩阵求逆失败，以及结果不稳定的原因：

1. 「线性化问题无可行解」的根源：初始值不满足等式约束 #

你设置的初始参数是c(1000000, 750000)，但等式约束要求x[1]+x[2]=17000000，而初始值的和只有1750000，离约束值差了一个数量级！solnp在迭代初期会尝试线性化问题寻找可行方向，但初始点完全不在可行域内，导致算法一开始就找不到满足所有约束的可行解，直接抛出这个错误。

另外，你的上下界和等式约束还有个隐含的逻辑：根据x[2] = 17000000 - x[1]和x[2] >=750000，可以推导出x[1] <= 16250000，你给x[1]设的上界16000000本身没问题，但初始值不满足等式约束才是触发错误的核心。

2. 「海森矩阵求逆失败」的原因：目标函数非光滑（不可导） #

你的目标函数是用多层ifelse构建的分段线性函数，在每个分段点（比如x[1]=1000000、x[2]=750000这类节点），函数的斜率会突然变化，也就是函数在这些点不可导。

solnp是依赖梯度的优化算法，需要目标函数具备光滑性（连续可导）才能正常计算海森矩阵（二阶导数矩阵），进而寻找最优方向。当算法碰到不可导的分段点时，海森矩阵会变得奇异（不可逆），就会出现「求逆失败」的错误。

3. 改变初始值导致结果变化的原因 #

因为目标函数是分段线性的，整个可行域内的最优解大概率出现在分段点或者约束边界上。不同的初始值会让算法落在不同的分段区间里，加上如果初始值不满足约束，算法的搜索方向会混乱，最终可能收敛到局部的边界点，而不是全局最优，自然会出现结果不稳定的情况。

解决建议：简化问题，降维优化 #

既然有严格的等式约束x[1]+x[2]=17000000，我们可以直接把问题降维成单变量优化，避开多变量约束的麻烦：

1. 用x[2] = 17000000 - x[1]替换目标函数里的x[2]
2. 推导x[1]的可行范围：
   • 由x[1] >=1000000
   • 由x[2] =17000000 -x[1] >=750000 → x[1] <=16250000
   • 结合你设置的上界x[1] <=16000000，最终x[1]的范围是[1000000, 16000000]
3. 用单变量优化函数（比如基础的optimise）来寻找最大值，因为分段线性函数的最大值必然在区间端点或者分段点上，optimise可以很好处理这类问题。

修改后的目标函数示例：

my_fn_single <- function(x1) {
  x2 <- 17000000 - x1
  # 计算x1对应的收益
  term1 <- ifelse(x1 < 1000000, 0, 
                  ifelse(x1 < 19500000, x1 * 0.0175, 
                         ifelse(x1 < 20500000, x1 * 0.0190, 
                                ifelse(x1 < 21500000, x1 * 0.0220, 
                                       ifelse(x1 < 22500000, x1 * 0.0270, x1 * 0.0340)))))
  # 计算x2对应的收益
  term2 <- ifelse(x2 < 750000, 0, 
                  ifelse(x2 < 15750000, x2 * 0.0100, 
                         ifelse(x2 < 16500000, x2 * 0.0150, 
                                ifelse(x2 < 17250000, x2 * 0.0275, 
                                       ifelse(x2 < 18000000, x2 * 0.0375, 
                                              ifelse(x2 < 18750000, x2 * 0.0450, x2 * 0.0550))))))
  return(term1 + term2) # 直接最大化，不用取负
}

# 调用优化
optimise(my_fn_single, interval = c(1000000, 16000000), maximum = TRUE)

这样既避开了多变量约束的问题，也不用依赖需要光滑函数的solnp，结果会更稳定可靠。

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