# How to set limits using constrOptim in R?

I am using constrOptim to minimize a log likelihood function for maximum likelihood estimation of parameters.

I wish to set the bounds on my parameters, but to not understand the constrOptim definition of the feasibility region.

The feasible region is defined by ui %*% theta - ci >= 0

I have a set of parameters with bounds [lower, upper]

a[0,5]  ie 0<a<5
b[0,Inf]
c[0,Inf]
e[0,1]

theta (starting values) = c(1, 1, 0.01,0.1)


What are the ui (constraint matrix (k x p)) and ci (constraint vector of length k) for these parameter bounds?

Is there a straightforward way to get from a list of upper and lower bounds to a ui and ci value?

Here's an example that we can use to illustrate ui and ci, with some extraneous output removed for brevity. It's maximizing the log likelihood of a normal distribution. In the first part, we use the optim function with box constraints, and in the second part, we use the constrOptim function with its version of the same box constraints.

# function to be optimized
> foo.unconstr <- function(par, x) -sum(dnorm(x, par, par, log=TRUE))

> x <- rnorm(100,1,1)

> optim(c(1,1), foo.unconstr, lower=c(0,0), upper=c(5,5), method="L-BFGS-B", x=x)
$par  1.147652 1.077654$value
 149.3724

>
> # constrOptim example
>
> ui <- cbind(c(1,-1,0,0),c(0,0,1,-1))
> ui
[,1] [,2]
[1,]    1    0
[2,]   -1    0
[3,]    0    1
[4,]    0   -1
> ci <- c(0, -5, 0, -5)
>
> constrOptim(c(1,1), foo.unconstr, grad=NULL, ui=u1, ci=c1, x=x)
$par  1.147690 1.077712$value
 149.3724

... blah blah blah ...

outer.iterations
 2

$barrier.value  -0.001079475 >  If you look at the ui matrix and imagine multiplying by the parameter vector to be optimized, call it$\theta$, you'll see that the result has four rows, the first of which is$\theta_1$, the second$-\theta_1$, the third$\theta_2$, and the fourth$-\theta_2$. Subtracting off the ci vector and enforcing the$\ge 0$constraint on each row results in$\theta_1 \ge 0$,$-\theta_1 + 5 \ge 0$,$\theta_2 \ge 0$and$-\theta_2 + 5 \ge 0$. Obviously, multiplying the second and fourth constraints by -1 and moving the constant to the right hand side gets you to$\theta_1 \le 5$and$\theta_2 \le 5$, the upper bound constraints. Just substitute your own values into the ci vector and add appropriate columns (if any) to the ui vector to get the box constraint set you want. • I am still slightly confused about how to construct ui and ci, but getting there. Thank you very much Apr 25, 2012 at 21:22 • This is a useful link to understand how to set up your ui and ci for constrOptim(): youtube.com/watch?v=MCvz-c6UUkw. The key takeaway is that all your constraints have to be set up in the Ax >= b format. Aug 7, 2015 at 8:11 Your constraints are of two types, either$\theta_i \geq a_i$, or$\theta_i \leq b_i$. The first ones are already in the right form (and the matrix ui is just the identity matrix), while the others can be written as$-\theta_i \geq - b_i$: ui is then$-I_n$and ci is$-b\$.

# Constraints
bounds <- matrix(c(
0,5,
0,Inf,
0,Inf,
0,1
), nc=2, byrow=TRUE)
colnames(bounds) <- c("lower", "upper")

# Convert the constraints to the ui and ci matrices
n <- nrow(bounds)
ui <- rbind( diag(n), -diag(n) )
ci <- c( bounds[,1], - bounds[,2] )

# Remove the infinite values
i <- as.vector(is.finite(bounds))
ui <- ui[i,]
ci <- ci[i]

# Constrained minimization
f <- function(u) sum((u+1)^2)


We can check how the constraint matrices ci and ui are interpreted:

# Print the constraints
k <- length(ci)
n <- dim(ui)
for(i in seq_len(k)) {
j <- which( ui[i,] != 0 )
cat(paste( ui[i,j], " * ", "x[", (1:n)[j], "]", sep="", collapse=" + " ))
cat(" >= " )
cat( ci[i], "\n" )
}
# 1 * x >= 0
# 1 * x >= 0
# 1 * x >= 0
# 1 * x >= 0
# -1 * x >= -5
# -1 * x >= -1


Some of the algorithms in optim allow you to specify the lower and upper bounds directly: that is probably easier to use.

• Does't removing the infinities affect which parameter is being bounded by which value? Apr 25, 2012 at 21:20
• The parameters correspond to the columns of ui, but what I remove I rows, i.e., constraints of the form x[i] <= Inf. Apr 25, 2012 at 22:43