R - 优化 X，使得所有 n 的 min(sum(X_i))，其中 X_n + X_(n-1) + X_(n-2) >= Y_n，其中 Y 对于所有 n 都是已知的

Question

我有一个长度为 50 的未知向量 X 和一个已知的常量向量 Y（长度为 50）。

我希望找到 X，使得对于 X_i>=0，sum(X_i) 最小化，约束条件为：

X_n + X_{n-1} >= Y_n

我不知道从哪里开始使用 R。

Answer 1

我想你可以尝试CVXR解决优化问题。

首先，让我们定义一个矩阵M，如下所示

M <- matrix(0,nrow = 10,ncol = 11)
for (i in 1:nrow(M)) {
  for (j in 1:ncol(M)) {
    if (j %in% (i+(0:1))) M[i,j] <- 1
  }
}

看起来像

> M
      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11]
 [1,]    1    1    0    0    0    0    0    0    0     0     0
 [2,]    0    1    1    0    0    0    0    0    0     0     0
 [3,]    0    0    1    1    0    0    0    0    0     0     0
 [4,]    0    0    0    1    1    0    0    0    0     0     0
 [5,]    0    0    0    0    1    1    0    0    0     0     0
 [6,]    0    0    0    0    0    1    1    0    0     0     0
 [7,]    0    0    0    0    0    0    1    1    0     0     0
 [8,]    0    0    0    0    0    0    0    1    1     0     0
 [9,]    0    0    0    0    0    0    0    0    1     1     0
[10,]    0    0    0    0    0    0    0    0    0     1     1

然后，我们构建objective函数以及约束

library(CVXR)
X <- Variable(11)
objective <- Minimize(sum(X))
constraints <- list( X>=0, M%*%X >= Y)
problem <- Problem(objective,constraints)
res <- solve(problem)

最后，我们可以通过res$getValue(X)

X

例子给定 Y 如下

set.seed(1)
Y <- runif(10)

我们可以得到

Xopt <- res$getValue(X)

> Xopt
              [,1]
 [1,] 1.667850e-07
 [2,] 3.072356e-01
 [3,] 6.488860e-02
 [4,] 6.214644e-01
 [5,] 2.867441e-01
 [6,] 1.486883e-02
 [7,] 8.835218e-01
 [8,] 7.476340e-02
 [9,] 5.860353e-01
[10,] 5.264372e-02
[11,] 9.142897e-03

另一个可能的选项可能是 pracma::fmincon，例如，

pracma:: fmincon(rep(0, 11),
  function(x) sum(x),
  A = -M,
  b = -Y,
  lb = 0,
)

Answer 2

这可以表示为线性规划

min 1'x
s.t. Ax >= y
     x >= 0

我们可以使用lpSolve包来解决这个问题。令 n 为 x 的长度，k 为每个约束中 x 的数量。问题正文中定义的约束对应于 k=2，但问题的主题对约束的定义不同，对应于 k=3。将有 n-k+1 个约束条件。

embed(1:n,k) 创建一个 n-k+1 x k 矩阵，其中每一行是 A 的对应行中的列索引。例如，对于 k=2 第一行embed(...) 是 1:2 因为 A 第一行的元素 1 和 2 是 1，而 A 行的其余部分是零。 embed 输出的第二行是 2:3，第三行是 3:4，等等。然后我们在行上应用替换以替换零的 n-vector，numeric(n)，在那些位置。 apply 的工作方式是它给出了我们想要的转置，所以我们将它转回得到 A.

最后我们运行线性程序。我们可以使用 str(out) 检查返回的输出组件。特别是我们将解决方案显示为 out$solution.

library(lpSolve)

# inputs - replace these three with your inputs
n <- 5 # no of x variables
k <- 2 # no of x's in each constraint
Y <- seq_len(n - k + 1)

A <- t(apply(embed(1:n, k), 1, replace, x = numeric(n), values = 1))
out <- lp("min", rep(1, n), A, ">=", Y)
out
## Success: the objective function is 6 

out$solution
## [1] 0 2 0 4 0

R - 优化 X，使得所有 n 的 min(sum(X_i))，其中 X_n + X_(n-1) + X_(n-2) >= Y_n，其中 Y 对于所有 n 都是已知的

R - Optimize X, such that min(sum(X_i)) for all n, where X_n + X_(n-1) + X_(n-2) >= Y_n, where Y is known for all n

r

optimization

vector

lag

inequality