如何在R的ROI包中修改二次优化目标函数以最大化x·V/(x^T Q x)

要在ROI中最大化目标函数 x·V/(x^T Q x)，由于ROI默认的二次目标无法直接表达这个分式形式，你可以利用支持非线性优化的solver="alabama"，通过自定义非线性目标函数来实现。具体步骤如下：

核心思路 #

原目标是最大化 f(x) = sum(x*V) / (x^T Q x)，等价于最小化 -f(x)（最大化一个函数等于最小化它的相反数）。我们可以用ROI的nonlinear_objective定义这个自定义目标，同时保留原有的线性约束和变量边界。

修改后的可运行代码 #

library(ROI)
library(Matrix)  # 用于生成对称矩阵
library(ROI.plugin.alabama)  # 加载alabama非线性求解插件

# 生成协方差矩阵和V向量（设置随机种子保证结果可复现）
set.seed(123)
C1 <- matrix(runif(100, min=10, max=100), nrow=10)
C <- as.matrix(forceSymmetric(C1))
V <- sqrt(diag(C))

# 定义全额投资约束：sum(x) = 1
full_invest <- L_constraint(rep(1, 10), "==", 1)

# 定义非线性目标函数：将最大化原目标转为最小化其相反数
nonlin_obj <- nonlinear_objective(
  fun = function(x) {
    numerator <- sum(x * V)
    denominator <- as.numeric(t(x) %*% C %*% x)
    -numerator / denominator
  },
  # 可选：提供解析梯度，提升求解稳定性和速度（省略则自动用数值梯度）
  grad = function(x) {
    numerator <- sum(x * V)
    denominator <- as.numeric(t(x) %*% C %*% x)
    grad_num <- V
    grad_den <- 2 * C %*% x
    # 应用商数求导法则计算梯度
    -(grad_num * denominator - numerator * grad_den) / (denominator^2)
  }
)

# 构建优化问题：补充投资组合权重的非负约束
qcqp <- OP(
  objective = nonlin_obj,
  constraints = full_invest,
  bounds = V_bound(ui = seq_len(10), lb = rep(0, 10), ub = rep(0.25, 10)),
  max = FALSE  # 目标是最小化-f(x)，对应原问题的最大化
)

# 求解优化问题
sol1 <- ROI_solve(qcqp, solver = "alabama", start = rep(1/10, 10))

# 查看结果
print(sol1)
x_opt <- solution(sol1)
# 计算原目标函数的最优值
cat("原目标函数最优值：", sum(x_opt * V) / as.numeric(t(x_opt) %*% C %*% x_opt), "\n")

关键说明 #

1. 自定义目标函数：通过nonlinear_objective直接定义分式目标的相反数，将原最大化问题转为ROI默认支持的最小化问题。
2. 解析梯度：提供解析梯度可以让alabama求解器更高效稳定，若不想编写梯度逻辑，可直接省略grad参数，求解器会自动用数值方法计算梯度。
3. 变量边界：补充了lb=rep(0,10)的非负约束（符合投资组合权重的实际意义），原代码仅设置了上限，建议保留该约束。
4. 初始点：选择均匀权重rep(1/10,10)作为初始点，是投资组合优化的常规初始值选择。

替代分式规划变形方案 #

如果你偏好将问题转化为带约束的形式，可通过引入额外变量t，将原问题等价转化为：

最大化 t
约束：
V^T x ≥ t * (x^T C x)
sum(x) = 1
0 ≤ x_i ≤ 0.25
t ≥ 0

对应的代码实现如下：

# 定义非线性约束：V^T x - t*(x^T C x) ≥ 0
nonlin_constr <- NL_constraint(
  fun = function(x_t) {
    x <- x_t[1:10]
    t_val <- x_t[11]
    sum(x*V) - t_val * as.numeric(t(x)%*%C%*%x)
  },
  dir = ">=",
  rhs = 0
)

# 构建优化问题（变量为c(x, t)）
qcqp2 <- OP(
  objective = L_objective(c(rep(0,10), 1)),  # 最大化t
  constraints = c(full_invest, nonlin_constr),
  bounds = V_bound(ui = 1:11, lb = c(rep(0,10), 0), ub = c(rep(0.25,10), Inf)),
  max = TRUE
)

sol2 <- ROI_solve(qcqp2, solver = "alabama", start = c(rep(1/10,10), 0.1))

不过这种方法需要额外引入变量t，不如直接自定义目标函数直观。

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