2

次の方程式を最小化したい:

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

1 に答える 1

4

ここでの最適化問題と前の質問の類似点を考えると、その質問に対する私の回答から直接言葉を借ります。ただし、それらはかなり異なるため (前の問題は線形計画問題で、これは二次計画問題であり、制約が異なります)、重複していません。

得られた最適化目標を展開しますQuw*ruw^2 - 2*Quw*ruw*yuw + Quw*yuw^2。これは決定変数 の 2 次関数であることがわかりますyuw。したがって、パッケージのsolve.QPメソッドをquadProg使用して最適化問題を解くことができます。

これを少し抽象化するために、入力行列の次元を想定R=20して説明しましょう。C=10次に、R*C決定変数があり、それらに order を割り当ててy11, y21, ... yR1, y12, y22, ... yR2, ..., y1C, y2C, ..., yRC、変数の行列の列を読み下げることができます。

から?solve.QP、目的は-d'b + 0.5b'Db決定変数 の形式をとることが読み取れますbd決定変数に対応するの要素yuwは値2*Quw*ruwを持ち、決定変数に対応する要素が値を取るD対角行列です。この関数では行列が正定値である必要があることに注意してください。したがって、すべてのペアに必要です。yuw2*Quwsolve.QPDQuw > 0u, w

最初のR*(C-1)制約は制約yuw >= yu,w+1に対応し、次の(R-1)*C制約は制約に対応しyuw >= yu-1,wます。次の2*R制約は制約 (およびyuC = 0として入力される) に対応し、最後の制約は(論理的には と同等) です。yuC >= 0-yuC >= 0-yR1 >= -100yR0 = 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
于 2015-06-21T22:33:03.897 に答える