R语言constrOptim约束优化报初始值不在可行域内部求解咨询

• constrOptim 仅支持严格不等式约束，要求初始值必须满足 ui %*% theta > ci（严格大于，不能取等），可行域必须存在非空内部。你把m2元素和为1这个等式约束拆成sum(m2)>=1和-sum(m2)>=-1两个不等式，所有满足约束的点都刚好取等，完全位于约束边界上，不存在符合要求的内点，必然触发报错。
• 目标函数参数顺序错误：constrOptim 会把待优化参数作为第一个入参传给目标函数，你写的const.fun第一个形参是固定矩阵m1，第二个才是待优化的m2，传参逻辑完全错位。

利用等式约束做变量替换，把3维优化问题降为2维，消除等式约束，让可行域存在明确内部点：

1. 由m2[1]+m2[2]+m2[3] = 1，直接令m2[3] = 1 - m2[1] - m2[2]，待优化参数变为2维向量theta = c(m2[1], m2[2])
2. 剩余约束全部为存在严格内点的不等式形式：m2[1]>=0、m2[2]>=0、1 - m2[1] - m2[2] >=0，可行域是二维平面上的三角形，存在大量内点可选做初始值。

对应实现代码：

# 重写目标函数，第一个参数固定为待优化的theta
const.fun_fix <- function(theta, m1){
  m2 <- c(theta[1], theta[2], 1 - theta[1] - theta[2])
  as.numeric(t(m2) %*% m1 %*% m2)
}
# 重写约束矩阵，对应ui%*%theta >= ci
ui_fix <- rbind(
  c(1, 0),   # theta[1] >= 0
  c(0, 1),   # theta[2] >= 0
  c(-1, -1)  #  -theta[1] - theta[2] >= -1 即 m2[3] >=0
)
ci_fix <- c(0, 0, -1)
# 选严格内点作为初始值，比如c(0.3,0.3)，对应m2=c(0.3,0.3,0.4)，所有元素为正、和为1，严格满足所有不等式
init <- c(0.3, 0.3)
# 调用优化，求最小值直接运行，求最大值需在control参数里设置fnscale=-1
res <- constrOptim(
  theta = init,
  f = const.fun_fix,
  grad = NULL,
  ui = ui_fix,
  ci = ci_fix,
  m1 = m1 # 传入固定的3*3矩阵m1
)
# 还原最优的m2
m2_opt <- c(res$par[1], res$par[2], 1 - res$par[1] - res$par[2])

你的目标函数是标准二次型t(m2) %*% m1 %*% m2，属于典型的二次规划问题，用quadprog包可以直接同时处理等式和不等式约束，不需要做变量替换，求解精度和效率都更高。
对应实现代码：

# 首次使用先运行 install.packages("quadprog")
library(quadprog)
# 组装二次规划标准形式参数：min 0.5*t(x)%*%Dmat%*%x - t(dvec)%*%x
Dmat <- 2 * m1 # 匹配原二次型x'm1x，对应0.5x'Dx的形式D=2m1
dvec <- rep(0, 3) # 目标无一次项，d取零向量
# 约束矩阵，前1个为等式约束，后3个为非负不等式约束
Amat <- t(rbind(
  c(1,1,1), # sum(m2) = 1
  c(1,0,0), # m2[1] >=0
  c(0,1,0), # m2[2] >=0
  c(0,0,1)  # m2[3] >=0
))
bvec <- c(1, 0, 0, 0) # 约束右侧值，满足t(Amat)%*%x >= bvec
# meq=1指定前1个约束为等式
res_qp <- solve.QP(Dmat = Dmat, dvec = dvec, Amat = Amat, bvec = bvec, meq = 1)
# 直接取最优m2
m2_opt_qp <- res_qp$solution

注意：如果你的优化目标是最大化二次型，只需要把目标函数取负转为最小化问题即可，也就是将代码里的Dmat设为-2*m1。

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