2016-05-14 4 views
0

目的関数でx1x2を表現する方法について質問があります。 ここはインターネット上の例です。x1 * x2をquadprogで表す方法R

例から
## 
## min x1^2 +2x2^2 + 4x3^2 - x1 - x2 + 5x3 
## x1 + x3 <= 1 
## x1 >= 5 
## x2 <= 0 
## 

P = 2*diag (c (1, 2, 4)); 
d = c (-1, -1, 5); 
A = matrix (0, nrow=3, ncol=3); 
A[1,] = c(-1, 0, -1); 
A[2,] = c(1, 0, 0); 
A[3,] = c(0, -1, 0); 
b = c(-1, 5, 0); 

、目的finctionは、それは、しかしP = 2*diag (c (1, 2, 4)); d = c (-1, -1, 5);

であるRでx1^2 +2x2^2 + 4x3^2 - x1 - x2 + 5x3

である私は、どのようにコマンドをキーにX1X2またはX1^2 のような目的関数を持っている場合R?

ありがとうございます。

+0

BTW、私はなぜP =​​ 2 * diag(c(1,2,4))この関数 "2" diagを掛ける必要がある別の質問を思いますか? – GGININDERRRRRRRRRRRRRRRRRRRRRR

+0

目的関数の非線形部分に1/2を掛けるため、 "2"が必要です。詳細はsolve.QPのドキュメントと私の答えを参照してください。 – digEmAll

答えて

0

solve.QP関数は、次式を最小化:

min 1/2 * (t(x) %*% D %*% x) - t(d) %*% x 
xは n変数(例えば、X1、X2 ...)、 n x n係数の D対称行列と n係数の dベクトルのベクトルである

。これらの行列積は、私が入力与えられたquadprogの目標を拡張するには、次の簡単な関数を作成しました、目的関数式に変換する方法を理解するために

# Given matrix D and vector d that will be passed to solve.QP, 
# it returns the objective function expression as string 
expandQPObj <- function(D, d, var.name = 'x') { 
    # helper function that expands a matrix product. NB it does not compute the result, 
    # it just returns the expressions of each element of the resulting matrix 
    expandMatrixProd <- function(m1, m2) { 
    m1 <- as.matrix(m1) 
    m2 <- as.matrix(m2) 
    if (ncol(m1) != nrow(m2)) { stop("incompatible dimensions") } 
    n <- nrow(m1) 
    m <- ncol(m2) 
    res <- matrix('', nrow = n, ncol = m) 
    for (i in 1:n) { 
     for (j in 1:m) { 
     a <- m1[i, ] 
     b <- m2[, j] 
     a <- ifelse(grepl('[+*]', a), paste0('(', a, ')'), a) 
     b <- ifelse(grepl('[+*]', b), paste0('(', b, ')'), b) 
     res[i, j] <- gsub('+-','-',paste(a, b, sep = '*', collapse = '+'),fixed = TRUE) 
     } 
    } 
    return(res) 
    } 
    D <- as.matrix(D) 
    d <- as.vector(d) 
    n <- length(d) 
    if (!all(dim(D) == n) || n == 0) { 
    stop('Dimensions problem: D should be an nxn matrix and d a vector of length n (n>0)') 
    } 
    xvec <- paste(var.name, 1:n, sep = '') 
    quadComp <- as.vector(Reduce(list(t(xvec), D, xvec),f=expandMatrixProd)) 
    linearComp <- as.vector(strMatrixMult(t(d), xvec)) 
    return(paste0('1/2*(', quadComp, ') - (', linearComp, ')')) 
} 

リテラル係数でそれを呼び出すことによって、我々得るために、そう

Dliteral <- matrix(paste0('D',1:9),nrow=3,byrow = T) 
#  [,1] [,2] [,3] 
#[1,] "D1" "D2" "D3" 
#[2,] "D4" "D5" "D6" 
#[3,] "D7" "D8" "D9" 

dliteral <- paste0('d',1:3) 
#[1] "d1" "d2" "d3" 

expandQPObj(Dliteral,dliteral) 
#"1/2*((x1*D1+x2*D4+x3*D7)*x1+(x1*D2+x2*D5+x3*D8)*x2+(x1*D3+x2*D6+x3*D9)*x3) - d1*x1+d2*x2+d3*x3" 

:簡単に私たちの目的の結果を得るために、これらの係数を設定する方法を理解することができます

min x1*x2 
(マトリックス Dは対称でなければならないためだけでなく、 D2=2!)、従って他のすべての係数= 0、我々は D2=1 and D4=1を設定する必要が

:得るために、同様

D = rbind(c(0,1,0),c(1,0,0),c(0,0,0)) 
#  [,1] [,2] [,3] 
#[1,] 0 1 0 
#[2,] 1 0 0 
#[3,] 0 0 0 
d <- rep(0,3) 
#[1] 0 0 0 

expandQPObj(D,d) 
[1] "1/2*((x1*0+x2*1+x3*0)*x1+(x1*1+x2*0+x3*0)*x2+(x1*0+x2*0+x3*0)*x3) - 0*x1+0*x2+0*x3" 

min x1^2 (== min x1*x1) 

D1=2と他のすべての係数= 0を設定する必要があります。

D = rbind(c(2,0,0),c(0,0,0),c(0,0,0)) 
#  [,1] [,2] [,3] 
#[1,] 2 0 0 
#[2,] 0 0 0 
#[3,] 0 0 0 

d <- rep(0,3) 
#[1] 0 0 0 

expandQPObj(D,d) 
# "1/2*((x1*2+x2*0+x3*0)*x1+(x1*0+x2*0+x3*0)*x2+(x1*0+x2*0+x3*0)*x3) - 0*x1+0*x2+0*x3" 
関連する問題