Search This Blog

Saturday, March 22, 2014

Binary Linear Programming using R

In the following series we will perform Linear optimizations using R. The first thing to do is to install the necessary packages in R. We will be using a package called lpSolveAPI for this exercise. We will start by downloading and installing the necessary packages followed by selecting a problem for project selection, and solving it within R

Installing the necessary R packages

First start R as an administrator. In Ubuntu, it will be sudo R.


Type command to install lpSolveAPI package as shown below.

install.packages(“lpSolveAPI”)



R will prompt us to choose an appropriate mirror. Choose the mirror closest to you.



Next install ggplot2



Next install reshape



Next install gridExtra


We can quit running R as an administrator, and start it again as a normal user, as shown below.


Now we are ready to formulate the linear programming problem and solve it within R.

Project Portfolio Selection Using R

The next thing we want to do is a pick up a relevant problem. For purpose of illustration, we are going to use an example of project portfolio selection example from a relevant paper “Strategic Selection of the Most Feasible Projects Using Linear Programming Models” by Yara Hamdan. The purpose is to pick up projects that are going to have the maximum benefit, and minimize costs along multiple dimensions.



In the paper, the author has cited 20 hypothetical projects that compete with each other for cost, labour and other similar constraints. The objective is to select as many projects as possible within the given constraints.

The overall data for all the projects is given in the following table, that is identical to the data given in the paper.



Project
Cost
Labor
Priority
Quality
Risk
Historical Info
New Skills
Commercial Complex
Tech. Advance
1
A
100000
300000
0
70
80
1200
720
1400
230
2
B
2600000
550000
0
130
130
2300
830
1950
410
3
C
500000
90000
0
420
40
550
560
670
1230
4
D
3200000
700000
0
360
140
3500
1000
3200
1200
5
E
2700000
500000
0
400
100
2500
870
2300
1000
6
F
4500000
850000
0
200
175
3800
1200
3600
750
7
G
700000
120000
0
150
65
850
660
800
560
8
H
250000
75000
0
90
30
400
400
470
320
9
I
750000
200000
0
320
70
900
680
890
900
10
J
2300000
450000
0
120
120
2000
880
2200
290
11
K
3500000
560000
0
300
155
3700
1100
3400
1000
12
L
1700000
350000
0
90
100
1500
750
1800
320
13
M
900000
270000
1
360
70
1000
700
1300
980
14
N
3600000
750000
0
400
160
2100
950
2400
1500
15
O
1200000
320000
0
250
110
1300
730
1000
900
16
P
2200000
480000
0
220
125
2100
860
2400
720
17
Q
600000
100000
0
160
60
780
700
800
520
18
R
800000
160000
0
100
75
950
790
900
250
19
S
3000000
500000
0
340
135
2900
920
3200
900
20
T
542663
70000
0
210
60
620
560
590
640

The paper provides the formulation in LINGO language. We shall solve the same problem using R programming language.

The R script for the above is as follows

# Load the lpSolveAPI library
library(lpSolveAPI)
maxrow <-9
maxcolumn <-20

# We create a new empty model with 9 rows and 20 columns
lpmodel <-make.lp(maxrow,maxcolumn)

# Setting the appropriate data for each project, that translates to a column in the lpSolve equation
set.column(lpmodel, 1, c(100000,300000,0,70,80,1200,720,1400,230))
set.column(lpmodel, 2, c(2600000,550000,0,130,130,2300,830,1950,410))
set.column(lpmodel, 3, c(500000,90000,0,420,40,550,560,670,1230))
set.column(lpmodel, 4, c(3200000,700000,0,360,140,3500,1000,3200,1200))
set.column(lpmodel, 5, c(2700000,500000,0,400,100,2500,870,2300,1000))
set.column(lpmodel, 6, c(4500000,850000,0,200,175,3800,1200,3600,750))
set.column(lpmodel, 7, c(700000,120000,0,150,65,850,660,800,560))
set.column(lpmodel, 8, c(250000,75000,0,90,30,400,400,470,320))
set.column(lpmodel, 9, c(750000,200000,0,320,70,900,680,890,900))
set.column(lpmodel, 10, c(2300000,450000,0,120,120,2000,880,2200,290))
set.column(lpmodel, 11, c(3500000,560000,0,300,155,3700,1100,3400,1000))
set.column(lpmodel, 12, c(1700000,350000,0,90,100,1500,750,1800,320))
set.column(lpmodel, 13, c(900000,270000,1,360,70,1000,700,1300,980))
set.column(lpmodel, 14, c(3600000,750000,0,400,160,2100,950,2400,1500))
set.column(lpmodel, 15, c(1200000,320000,0,250,110,1300,730,1000,900))
set.column(lpmodel, 16, c(2200000,480000,0,220,125,2100,860,2400,720))
set.column(lpmodel, 17, c(600000,100000,0,160,60,780,700,800,520))
set.column(lpmodel, 18, c(800000,160000,0,100,75,950,790,900,250))
set.column(lpmodel, 19, c(3000000,500000,0,340,135,2900,920,3200,900))
set.column(lpmodel, 20, c(542663,70000,0,210,60,620,560,590,640))

# Declaring the rhs values
cost <-10000000
hours <-2000000
priority  <- 1
quality <-1000
risk <-600
histinfo <-10000
newexp <-5000
comm <-10000
tech <-3000

# To start, each project has uniform weightage
set.objfn(lpmodel,rep(1,maxcolumn))

# Each constraint provides the upper bound of rhs
set.constr.type(lpmodel,rep("<=", maxrow))

# Set the RHS values
set.rhs(lpmodel, c(cost,hours,priority,quality,risk,histinfo,newexp,comm,tech))

# Setting all variable types to be binary
set.type(lpmodel,1,"binary")
set.type(lpmodel,2,"binary")
set.type(lpmodel,3,"binary")
set.type(lpmodel,4,"binary")
set.type(lpmodel,5,"binary")
set.type(lpmodel,6,"binary")
set.type(lpmodel,7,"binary")
set.type(lpmodel,8,"binary")
set.type(lpmodel,9,"binary")
set.type(lpmodel,10,"binary")
set.type(lpmodel,11,"binary")
set.type(lpmodel,12,"binary")
set.type(lpmodel,13,"binary")
set.type(lpmodel,14,"binary")
set.type(lpmodel,15,"binary")
set.type(lpmodel,16,"binary")
set.type(lpmodel,17,"binary")
set.type(lpmodel,18,"binary")
set.type(lpmodel,19,"binary")
set.type(lpmodel,20,"binary")

#Setting the objective function direction
lp.control(lpmodel,sense='max')

#Finally lets name the model variables
Rownames <- c("Cost","Labor","Priority","Quality","Risk","Historical Info","New Skills","Commercial Complex","Tech. Advance")

Colnames <- c("Project A","Project B","Project C","Project D","Project E","Project F","Project G","Project H","Project I","Project J","Project K","Project L","Project M","Project N","Project O","Project P","Project Q","Project R","Project S","Project T")
dimnames(lpmodel) <- list(Rownames, Colnames)

# Solving the model
solve(lpmodel)

get.objective(lpmodel)
get.variables(lpmodel)
get.constraints(lpmodel)

Solving the model, we can see that the same 7 projects got selected using R lpSolveAPI

> solve(lpmodel)
[1] 0
> get.objective(lpmodel)
[1] 7
> get.variables(lpmodel)
 [1] 1 1 0 0 0 0 1 1 0 0 0 1 0 0 0 0 1 1 0 0
> get.constraints(lpmodel)
[1] 6750000 1655000       0     790     540    7980    4850    8120    2610

Here is the screenshot


That's it folks.... We solved a Binary LP problem using R.