Roman Hahn
Roman Hahn

Reputation: 83

Constraint Optimization with one parameter included in the constraint of the other

I want to calculate the following

Optimization Problem

So I want to find Theta and Sigma that maximizes the function. The constraints are:

> Theta>-Sigma
> -1<Sigma<1

So one of my problem is that I dont know how to deal with the fact that one parameter is included in the constraint of the other Parameter, that I want to optimize over.

I tried with optim(), constrOptim and dfoptim!

Using optim():

k=8
i=1:(k-1)
x=c(5,0.2)
n=24
nj=c(3,4,8,1,1,4,2,1)


EPPF <- function(x,n,nj) {
  y=(x[1]+1):(x[1]+1+(n-1)-1)
  z=-(prod(x[1]+i*x[2])/(prod(y))*prod(sapply(nj, hfun)))
  return(z)}

hfun <- function(p){
  h=(1-x[2]):((1-x[2])+p-1)
  hfun=prod(h)
  return(hfun)
}

> optim(c(6,0.3), fn=EPPF,method = "L-BFGS-B", n=n,nj=nj, lower = c(-x[1],-1), upper = c(Inf,1))
$par
[1] 6.0 0.3

$value
[1] -1.258458e-15

$counts
function gradient 
       2        2 

$convergence
[1] 0

$message
[1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH"

I also tried using a constructor function:

make.EPPF <- function(n,nj,fixed=c(FALSE,FALSE)){
  params <-fixed
  function(p) {
    hfun <- function(y){
      h=(1-sigma):((1-sigma)+y-1)
      hfun=prod(h)
      return(hfun)
    }

    params[!fixed] <- p
    theta <- params[1]
    sigma <- params[2]
    y=(theta+1):(theta+1+(n-1)-1)
    z=(prod(theta+i*sigma)/(prod(y))*prod(sapply(nj, hfun)))
    z
  }
}

EPPF <-make.EPPF (n,nj)

> optim(c(theta=6, sigma=0.5), fn=EPPF,method = "L-BFGS-B",lower = c(-sigma,-1), upper = c(Inf,1))
Error in optim(c(theta = 6, sigma = 0.5), fn = EPPF, method = "L-BFGS-B",  : 
  object 'sigma' not found

Using constrOptim():

> A <- matrix(c(1,1,0,1,0,-1),3,2,byrow=T)
> b <- c(0,-1,-1)
> 
> constrOptim(c(3,0.3),EPPF,NULL,A,b, control=list(fnscale=-1))
$par
[1] 3.0 0.3

$value
[1] 9.712117e-16

$counts
[1] 0

$convergence
[1] 0

$message
NULL

$outer.iterations
[1] 1

$barrier.value
[1] 7.313452e-05

Using Package dfoptim:

> library(dfoptim)
> nmkb(x=c(6,0.3), EPPF, lower=c(-x[2],-1), upper=c(Inf, 1 ))
Error in par < lower : 
  comparison (3) is possible only for atomic and list types

Either there is for some reasons no convergence or some other Errors.

I am relativ new to programming and R and would really appreciate if someone could help me. Thanks!

Upvotes: 0

Views: 1193

Answers (1)

Ott Toomet
Ott Toomet

Reputation: 1956

These are 3 linear inequality constraints:

sigma + theta > 0
sigma + 1 > 0
-sigma + 1 > 0

You can do this in maxLik. But note that maxLik maximizes the function, hence remove the '-' in front of z. Here is the code that works for me (using Rscript):

k=8
i=1:(k-1)
x=c(5,0.2)
n=24
nj=c(3,4,8,1,1,4,2,1)

EPPF <- function(x,n,nj) {
   theta <- x[1]
   sigma <- x[2]
  y=(x[1]+1):(x[1]+1+(n-1)-1)
   z <- prod(x[1]+i*x[2])/(prod(y))*prod(sapply(nj, hfun))
   z <- log(z)
   return(z)
}

hfun <- function(p){
  h=(1-x[2]):((1-x[2])+p-1)
  hfun=prod(h)
  return(hfun)
}

library(maxLik)
constraints <- list(ineqA=matrix(c(1,0,0,1,1,-1),3,2), ineqB=c(0,1,1))
m <- maxBFGS(EPPF, start=c(6,0.3), constraints=constraints, n=n, nj=nj)
print(summary(m))

I also took logarithm of the result as this leads to more "human" numbers. Otherwise you have to re-tune the stopping conditions. The answer seems to be -1, 1.

Upvotes: 1

Related Questions