Matching algorithms in R (bipartite matching, Hung

Question

I wonder how to set up some example some fundamental matching procedures in R. There are many examples in various programming languages, but I have not yet found a good example for R.

Let’s say I want to match students to projects and I would consider 3 alternative approaches which I came across when googling on this issue:

1) Bipartite matching case: I ask each student to name 3 projects to work on (without stating any preference ranking among those 3).

ID  T.1 T.2 T.3 T.4 T.5 T.6 T.7
1   1   1   1   0   0   0   0
2   0   0   0   0   1   1   1
3   0   1   1   1   0   0   0
4   0   0   0   1   1   1   0
5   1   0   1   0   1   0   0
6   0   1   0   0   0   1   1
7   0   1   1   0   1   0   0

--

d.1 <- structure(list(Student.ID = 1:7, Project.1 = c(1L, 0L, 0L, 0L, 
1L, 0L, 0L), Project.2 = c(1L, 0L, 1L, 0L, 0L, 1L, 1L), Project.3 = c(1L, 
0L, 1L, 0L, 1L, 0L, 1L), Project.4 = c(0L, 0L, 1L, 1L, 0L, 0L, 
0L), Project.5 = c(0L, 1L, 0L, 1L, 1L, 0L, 1L), Project.6 = c(0L, 
1L, 0L, 1L, 0L, 1L, 0L), Project.7 = c(0L, 1L, 0L, 0L, 0L, 1L, 
0L)), .Names = c("Student.ID", "Project.1", "Project.2", "Project.3", 
"Project.4", "Project.5", "Project.6", "Project.7"), class = "data.frame", row.names = c(NA, 
-7L))

2) Hungarian algorithm: I ask each student name 3 projects to work on WITH stating a preference ranking among those 3. As far as I understood the reasoning when applying the algorithm in this case would be something like: the better the rank the lower the “costs” to the student.

ID  T.1 T.2 T.3 T.4 T.5 T.6 T.7
1   3   2   1   na  na  na  na
2   na  na  na  na  1   2   3
3   na  1   3   2   na  na  na
4   na  na  na  1   2   3   na
5   2   na  3   na  1   na  na
6   na  3   na  na  na  2   1
7   na  1   2   na  3   na  na

--

d.2 <- structure(list(Student.ID = 1:7, Project.1 = structure(c(2L, 3L, 
3L, 3L, 1L, 3L, 3L), .Label = c("2", "3", "na"), class = "factor"), 
    Project.2 = structure(c(2L, 4L, 1L, 4L, 4L, 3L, 1L), .Label = c("1", 
    "2", "3", "na"), class = "factor"), Project.3 = structure(c(1L, 
    4L, 3L, 4L, 3L, 4L, 2L), .Label = c("1", "2", "3", "na"), class = "factor"), 
    Project.4 = structure(c(3L, 3L, 2L, 1L, 3L, 3L, 3L), .Label = c("1", 
    "2", "na"), class = "factor"), Project.5 = structure(c(4L, 
    1L, 4L, 2L, 1L, 4L, 3L), .Label = c("1", "2", "3", "na"), class = "factor"), 
    Project.6 = structure(c(3L, 1L, 3L, 2L, 3L, 1L, 3L), .Label = c("2", 
    "3", "na"), class = "factor"), Project.7 = structure(c(3L, 
    2L, 3L, 3L, 3L, 1L, 3L), .Label = c("1", "3", "na"), class = "factor")), .Names = c("Student.ID", 
"Project.1", "Project.2", "Project.3", "Project.4", "Project.5", 
"Project.6", "Project.7"), class = "data.frame", row.names = c(NA, 
-7L))

3) ??? approach: This should be pretty much related to (2). However, I think it is probably a better/ fairer approach (at least in the setting of the example). The students cannot pick projects, they even don’t know about the projects, but they have rate their qualifications (1 “not existent” to 10 “professional level”) with regards to a certain skillset. Further, the lecturer has rated the required skillset for every project. In addition to (2), a first step would be to calculate a similarity matrix and then to run the optimization routine from above.

PS: Programming Skills
SK: Statistical Knowledge
IE: Industry Experience

ID  PS  SK  IE
1   10  9   8
2   1   2   10
3   10  2   5
4   2   5   3
5   10  2   10
6   1   10  1
7   5   5   5

--

d.3a <- structure(list(Student.ID = 1:7, Programming.Skills = c(10L, 1L, 
10L, 2L, 10L, 1L, 5L), Statistical.knowlegde = c(9L, 2L, 2L, 
5L, 2L, 10L, 5L), Industry.Expertise = c(8L, 10L, 5L, 3L, 10L, 
1L, 5L)), .Names = c("Student.ID", "Programming.Skills", "Statistical.knowlegde", 
"Industry.Expertise"), class = "data.frame", row.names = c(NA, 
-7L))

--

T: Topic ID
PS: Programming Skills
SK: Statistical Knowledge
IE: Industry Experience

T  PS   SK  IE
1   10  5   1
2   1   1   5
3   10  10  10
4   2   8   3
5   4   3   2
6   1   1   1
7   5   7   2

--

d.3b <- structure(list(Project.ID = 1:7, Programming.Skills = c(10L, 
1L, 10L, 2L, 4L, 1L, 5L), Statistical.Knowlegde = c(5L, 1L, 10L, 
8L, 3L, 1L, 7L), Industry.Expertise = c(1L, 5L, 10L, 3L, 2L, 
1L, 2L)), .Names = c("Project.ID", "Programming.Skills", "Statistical.Knowlegde", 
"Industry.Expertise"), class = "data.frame", row.names = c(NA, 
-7L))

I would appreciate any help in implementing those 3 approaches in R. Thank you.

UPDATE: The following questions seem to be related, but none show how to solve it in R: https://math.stackexchange.com/questions/132829/group-membership-assignment-by-preferences-optimization-problem https://superuser.com/questions/467577/using-optimization-to-assign-by-preference

Answer 1

Here are possible solutions using bipartite matching and the Hungarian algorithm.

My proposed solution using bipartite matching might not be what you have in mind. All the code below does is sample randomly for a specified number of iterations, after which at least one solution hopefully will have been identified. This might require a large number of iterations and a long time with large problems. The code below found three solutions to your example problem within 200 iterations.

matrix1 <- matrix(c( 1,   1,   1,  NA,  NA,  NA,  NA,
                    NA,  NA,  NA,  NA,   1,   1,   1,
                    NA,   1,   1,   1,  NA,  NA,  NA,
                    NA,  NA,  NA,   1,   1,   1,  NA,
                     1,  NA,   1,  NA,   1,  NA,  NA,
                    NA,   1,  NA,  NA,  NA,   1,   1,
                    NA,   1,   1,  NA,   1,  NA,  NA), nrow=7, byrow=TRUE)

set.seed(1234)

iters <- 200

my.match <- matrix(NA, nrow=iters, ncol=ncol(matrix1))

for(i in 1:iters) {

     for(j in 1:nrow(matrix1)) {

          my.match[i,j] <- sample(which(matrix1[j,] == 1), 1)

     }
}

n.unique <- apply(my.match, 1, function(x) length(unique(x)))

my.match[n.unique==ncol(matrix1),]

#      [,1] [,2] [,3] [,4] [,5] [,6] [,7]
# [1,]    3    7    4    6    1    2    5
# [2,]    1    7    4    5    3    6    2
# [3,]    3    5    4    6    1    7    2

Here is code for the Hungarian algorithm using package clue and solve_LSAP() as @jackStinger suggested. For this to work I had to replace the missing observations and I arbitrarily replaced them with 4. Person 5 did not get their first choice and Person 7 did not get any of their three choices.

library(clue)

matrix1 <- matrix(c( 3,   2,   1,   4,   4,   4,   4,
                     4,   4,   4,   4,   1,   2,   3,
                     4,   1,   3,   2,   4,   4,   4,
                     4,   4,   4,   1,   2,   3,   4,
                     2,   4,   3,   4,   1,   4,   4,
                     4,   3,   4,   4,   4,   2,   1,
                     4,   1,   2,   4,   3,   4,   4), nrow=7, byrow=TRUE)

matrix1

solve_LSAP(matrix1, maximum = FALSE)

# Optimal assignment:
# 1 => 3, 2 => 5, 3 => 2, 4 => 4, 5 => 1, 6 => 7, 7 => 6

Here is a link to a site showing how the Hungarian algorithm works: http://www.wikihow.com/Use-the-Hungarian-Algorithm

EDIT: June 5, 2014

Here is my first stab at optimizing the third scenario. I randomly assign each student to a project, then calculate the cost for that set of assignments. Cost is calculated by finding the difference between a student's skill set and the project's required skills. The absolute values of those differences are summed to give a total cost for the seven assignments.

Below I repeat the process 10,000 times and identify which of those 10,000 assignments results in the lowest cost.

An alternative approach would be to do an exhaustive search of all possible student-project assignments.

Neither the random search nor the exhaustive search is likely what you had in mind. However, the former will give an approximate optimal solution and the latter would give an exact optimal solution.

I might return to this problem later.

students <- matrix(c(10,   9,   8,
                      1,   2,  10,
                     10,   2,   5,
                      2,   5,   3,
                     10,   2,  10,
                      1,  10,   1,
                      5,   5,   5), nrow=7, ncol=3, byrow=TRUE)

projects <- matrix(c(10,   5,    1,
                      1,   1,    5,
                     10,  10,   10,
                      2,   8,    3,
                      4,   3,    2,
                      1,   1,    1,
                      5,   7,    2), nrow=7, ncol=3, byrow=TRUE)

iters <- 10000

# col = student, cell = project
assignments <- matrix(NA, nrow=iters, ncol=nrow(students))

for(i in 1:iters) {
      assignments[i,1:7] <- sample(7,7,replace=FALSE)
}

cost <- matrix(NA, nrow=iters, ncol=nrow(students))

for(i in 1:iters) {

     for(j in 1:nrow(students)) {

          student <- j
          project <- assignments[i,student]

          student.cost <- rep(NA,3)

          for(k in 1:3) {     

               student.cost[k] <- abs(students[student,k] - projects[project,k])

          } 

          cost[i,j] <- sum(student.cost)

     }

}


total.costs <- rowSums(cost)

assignment.costs <- cbind(assignments, total.costs)
head(assignment.costs)

assignment.costs[assignment.costs[,8]==min(assignment.costs[,8]),]

#                    total.costs
# [1,] 3 2 1 4 5 6 7          48
# [2,] 3 2 1 6 5 4 7          48
# [3,] 3 2 1 6 5 4 7          48

# student 1, project 3, cost = 3
# student 2, project 2, cost = 6
# student 3, project 1, cost = 7
# student 4, project 4, cost = 3
# student 5, project 5, cost = 15
# student 6, project 6, cost = 9
# student 7, project 7, cost = 5

3+6+7+3+15+9+5

# [1] 48

EDIT: June 6, 2014

Here is the exhaustive search. There are only 5040 possible ways to assign projects to the seven students. This search returns four optimal solutions:

students <- matrix(c(10,   9,   8,
                      1,   2,  10,
                     10,   2,   5,
                      2,   5,   3,
                     10,   2,  10,
                      1,  10,   1,
                      5,   5,   5), nrow=7, ncol=3, byrow=TRUE)

projects <- matrix(c(10,   5,    1,
                      1,   1,    5,
                     10,  10,   10,
                      2,   8,    3,
                      4,   3,    2,
                      1,   1,    1,
                      5,   7,    2), nrow=7, ncol=3, byrow=TRUE)

library(combinat)

n <- nrow(students)

assignments <- permn(1:n)
assignments <- do.call(rbind, assignments)
dim(assignments)

# column of assignments = student
# row of assignments = iteration
# cell of assignments = project

cost <- matrix(NA, nrow=nrow(assignments), ncol=n)

for(i in 1:(nrow(assignments))) {
     for(student in 1:n) {

          project      <- assignments[i,student]
          student.cost <- rep(NA,3)

          for(k in 1:3) {     
               student.cost[k] <- abs(students[student,k] - projects[project,k])
          } 

          cost[i,student] <- sum(student.cost)
     }
}


total.costs <- rowSums(cost)

assignment.costs <- cbind(assignments, total.costs)
head(assignment.costs)

assignment.costs[assignment.costs[,(n+1)]==min(assignment.costs[,(n+1)]),]

                   total.costs
[1,] 3 2 5 4 1 6 7          48
[2,] 3 2 5 6 1 4 7          48
[3,] 3 2 1 6 5 4 7          48
[4,] 3 2 1 4 5 6 7          48