Как минимизировать функцию в R с двумя ограничениями? - PullRequest
0 голосов
/ 17 апреля 2019

Я хочу решить задачу минимизации с двумя ограничениями в R.

Задача: мин x ^ T H x st

(1) e ^ T *x = 1,

(2) 0

H - матрица 20x20, а e - вектор единицс длиной 20

Как я могу это сделать?Я посмотрел на optimize, optim и optimix, но почему-то не знаю, с чего начать.

Я начал создавать свою целевую функцию:

f<- function(x) {t(x)%*%H%*%x}

и ограничение (1)

g<- function(x) {t(e)*x=1}

Но я не знаю, как сформулировать ограничение (2),

Я также не знаю, какая функция оптимизации подходит для этой проблемы.Я благодарен за любой совет.

В конце я хочу получить вектор x с 20 значениями.

Ответы [ 2 ]

0 голосов
/ 18 апреля 2019

Вот как это сделать с пакетом CVXR.

Обратите внимание, что (1) совпадает с sum(x_i) = 1, а (2) можно упростить до x_i > 0, поскольку sum(x_i) = 1 и x_i > 0 => x_i <= 1.

H <- rWishart(1, df = 30, diag(20))[,,1] # a 20x20 symmetric positive matrix

library(CVXR)

# the variable
x <- Variable(20)

# objective 
objective <- Minimize(quad_form(x, H))

# define problem
constraint1 <- sum(x) == 1
constraint2 <- x > 0
problem <- Problem(objective, constraints = list(constraint1, constraint2))

# solve problem
result <- solve(problem)

# results
result$getValue(x) # optimal x
result$value # value of objective at optimal x
0 голосов
/ 17 апреля 2019

проверьте пакет quadprog.Он имеет функцию solve.QP для решения следующего квадратичного программирования:

min (-d ^ T b + 1/2 b ^ TD b) с ограничениями A ^ T b> = b_0.

В вашем случае D = H * 2, d = 0.

Вот как вы можете построить матрицу A и вектор b0 для вашей задачи:

Amat <- t(rbind(rep(1, n), diag(1, 20), -diag(1, 20))
b0 <- c(1, rep(0, 20), rep(-1, 20))

Тогда вы можетеrun

library(quadprog)
solve.QP(Dmat = 2*H, dvec = 0, Amat = Amat, bvec = b0, meq = 1)

Обратите внимание, что meq = 1 указывает, что первое неравенство должно быть равенством.

Если вам нужны только оптимальные значения x, попробуйте следующее:

solve.QP(Dmat = 2*H, dvec = 0, Amat = Amat, bvec = b0, meq = 1)$solution

Надеюсь, это поможет.

...