2

我希望尽量减少以下等式:R:二次规划/等渗回归

F=SUM{u 1:20}sum{w 1:10} Quw(ruw-yuw)^2 

有以下限制:

yuw >= yu,w+1 
yuw >= yu-1,w 
y20,0 = 100 
y0,10 = 0 
yu,10 = 0 

我有一个20×10 RUW和20 * 10 QUW矩阵,我现在需要生成一个遵循约束的yuw矩阵。我在R编码,熟悉lpsolve,optimx和quadprog软件包,但不知道如何将它们用于这个特定的问题。我知道我必须使用quadprog软件包,因为这是一个二次编程问题。我并不是在寻找完整的答案,我想要了解如何构建约束矩阵以及解决问题的最佳方法。

回答

4

鉴于此处优化问题与your previous question的相似之处,我将直接从我对该问题的回答中借用一些语言。然而它们有很大的不同(前面的问题是线性规划问题,这是一个二次规划问题,约束不同),所以它们不是重复的。

扩大优化目标我们得到Quw*ruw^2 - 2*Quw*ruw*yuw + Quw*yuw^2。我们看到这是决策变量yuw的二次函数,因此quadProg包的solve.QP方法可用于解决优化问题。

为了抽象出一点,我们假设R=20C=10描述了输入矩阵的维数。然后有R*C决策变量,我们可以将它们分配的顺序为y11, y21, ... yR1, y12, y22, ... yR2, ..., y1C, y2C, ..., yRC,读取变量矩阵的列。

?solve.QP,我们看到目标的形式为-d'b + 0.5b'Db决策变量b。对应于决策变量yuwd的元素具有值2*Quw*ruwD是对角矩阵,其中元素对应于取值2*Quw的决策变量yuw。请注意0​​函数要求D矩阵是正定的,所以我们需要Quw > 0对每个u, w对。

第一个R*(C-1)约束对应于yuw >= yu,w+1约束,并且下一个(R-1)*C约束对应于约束yuw >= yu-1,w。下一个约束条件对应于yuC = 0约束条件(作为yuC >= 0-yuC >= 0输入),最后一个约束条件为-yR1 >= -100(逻辑上等于yR0 = 100)。

我们可以输入该模型到quadProg包与下述R命令,使用随机输入数据:

# Sample data 
set.seed(144) 
Quw <- matrix(runif(200), nrow=20) 
ruw <- matrix(runif(200), nrow=20) 
R <- nrow(Quw) 
C <- ncol(Quw) 

# Build constraint matrix 
part1 <- matrix(0, nrow=R*(C-1), ncol=R*C) 
part1[cbind(1:(R*C-R), 1:(R*C-R))] <- 1 
part1[cbind(1:(R*C-R), (R+1):(R*C))] <- -1 
part2 <- matrix(0, nrow=(R-1)*C, ncol=R*C) 
pos2 <- as.vector(sapply(2:R, function(r) r+seq(0, R*(C-1), by=R))) 
part2[cbind(1:nrow(part2), pos2)] <- 1 
part2[cbind(1:nrow(part2), pos2-1)] <- -1 
part3 <- matrix(0, nrow=2*R, ncol=R*C) 
part3[cbind(1:R, (R*C-R+1):(R*C))] <- 1 
part3[cbind((R+1):(2*R), (R*C-R+1):(R*C))] <- -1 
part4 <- rep(0, R*C) 
part4[R] <- -1 
const.mat <- rbind(part1, part2, part3, part4) 

library(quadProg) 
mod <- solve.QP(Dmat = 2*diag(as.vector(Quw)), 
       dvec = 2*as.vector(ruw)*as.vector(Quw), 
       Amat = t(const.mat), 
       bvec = c(rep(0, nrow(const.mat)-1), -100)) 

我们现在可以访问模型的解决方案:

# Objective (including the constant term): 
mod$value + sum(Quw*ruw^2) 
# [1] 9.14478 
matrix(mod$solution, nrow=R) 
#   [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9]  [,10] 
# [1,] 0.4346995 0.4346995 0.4346995 0.4346995 0.4346995 0.3215992 0.1818095 0.1818095 0.1818095 0.000000e+00 
# [2,] 0.4346995 0.4346995 0.4346995 0.4346995 0.4346995 0.4346995 0.4346995 0.2882339 0.2882339 0.000000e+00 
# [3,] 0.4346995 0.4346995 0.4346995 0.4346995 0.4346995 0.4346995 0.4346995 0.2882339 0.2882339 2.775558e-17 
# [4,] 0.5728478 0.4346995 0.4346995 0.4346995 0.4346995 0.4346995 0.4346995 0.2882339 0.2882339 0.000000e+00 
# [5,] 0.5728478 0.5111456 0.5111456 0.4699046 0.4346995 0.4346995 0.4346995 0.2882339 0.2882339 0.000000e+00 
# [6,] 0.5728478 0.5111456 0.5111456 0.4699046 0.4346995 0.4346995 0.4346995 0.2882339 0.2882339 0.000000e+00 
# [7,] 0.5728478 0.5111456 0.5111456 0.4699046 0.4346995 0.4346995 0.4346995 0.2882339 0.2882339 0.000000e+00 
# [8,] 0.5728478 0.5111456 0.5111456 0.5111456 0.4346995 0.4346995 0.4346995 0.2882339 0.2882339 0.000000e+00 
# [9,] 0.5728478 0.5111456 0.5111456 0.5111456 0.4346995 0.4346995 0.4346995 0.4346995 0.2882339 1.110223e-16 
# [10,] 0.5728478 0.5111456 0.5111456 0.5111456 0.4346995 0.4346995 0.4346995 0.4346995 0.4346995 0.000000e+00 
# [11,] 0.6298100 0.5111456 0.5111456 0.5111456 0.4518205 0.4346995 0.4346995 0.4346995 0.4346995 0.000000e+00 
# [12,] 0.6298100 0.5111456 0.5111456 0.5111456 0.4518205 0.4346995 0.4346995 0.4346995 0.4346995 0.000000e+00 
# [13,] 0.6298100 0.5111456 0.5111456 0.5111456 0.4518205 0.4346995 0.4346995 0.4346995 0.4346995 0.000000e+00 
# [14,] 0.6298100 0.5111456 0.5111456 0.5111456 0.4518205 0.4346995 0.4346995 0.4346995 0.4346995 0.000000e+00 
# [15,] 0.6298100 0.6009718 0.5111456 0.5111456 0.4518205 0.4346995 0.4346995 0.4346995 0.4346995 0.000000e+00 
# [16,] 0.6298100 0.6009718 0.6009718 0.6009718 0.4518205 0.4346995 0.4346995 0.4346995 0.4346995 0.000000e+00 
# [17,] 0.6298100 0.6009718 0.6009718 0.6009718 0.6009718 0.4346995 0.4346995 0.4346995 0.4346995 0.000000e+00 
# [18,] 0.6298100 0.6009718 0.6009718 0.6009718 0.6009718 0.6009718 0.4346995 0.4346995 0.4346995 0.000000e+00 
# [19,] 0.6298100 0.6009718 0.6009718 0.6009718 0.6009718 0.6009718 0.4346995 0.4346995 0.4346995 0.000000e+00 
# [20,] 0.6298100 0.6009718 0.6009718 0.6009718 0.6009718 0.6009718 0.5643033 0.5643033 0.5643033 0.000000e+00 
+0

您好Josilber,一些yuw值与沿着行的值相同,也与沿列向下的值相同。这是否遵守yuw> = yu,w + 1和yu,w> = yu-1,w约束? – Ankit

+0

另外,第一个值(20,1)应该是100,并且y0,10 = 0 yu,10 = 0约束不成立,请问为什么? – Ankit

+0

@Ankit元素(1,1)位于左上方,所以'yu,w> = yu,w + 1'表示数字在您右移至左移时不减少,'yu,w> = yu -1,w'表示数字在您从上到下移动时不减少。我不会显示'y20,0',因为它没有出现在客观值中。从'y20,0 = 100'和'yu,w> = yu,w + 1'我得出结论:'y20,1 <= 100',这是我添加的约束。正如你所看到的,'y20,1'的最终选择值是0.6298100。 – josliber