Optimization in R - Job Assignment for highest productivity

Introduction

Business units often must decide on job assignments under various constrains, to minimize time taken to complete a set of job. Mathematics can help with such decisions, as shown below using a simplified example.

Problem Definition

Suppose you have 5 contractors (C1 to C5) and 5 buildings (B1 to B5). Each contractor provides a quote in terms of cost and time to complete each building. The goal is to select a set of contractors such that the total cost and time are minimized.

Data

Create matrices for costs and times for each contractor for each building. Also, define the budget and time limits.

library(lpSolve)
library(ggplot2)

# Cost matrix (rows: contractors, columns: buildings)
costs <- matrix(c(
  100, 120, 110, 130, 140,
  90,  100, 95,  110, 105,
  80,  85,  90,  95,  100,
  150, 160, 155, 165, 170,
  110, 115, 120, 125, 130
), nrow = 5, byrow = TRUE)

# Time matrix (rows: contractors, columns: buildings)
times <- matrix(c(
  10, 12, 11, 13, 11,
  9,  10, 9.5, 11, 10.5,
  8,  8.5, 9,  9, 10,
  17, 16, 15.5, 16.5, 15,
  11, 11.5, 12, 12.5, 13
), nrow = 5, byrow = TRUE)

# Adjusted constraints
budget_limit <- 600  # Increase budget limit
time_limit <- 60     # Increase time limit

# Number of contractors and buildings
num_contractors <- 5
num_buildings <- 5

# Convert matrices to vectors
costs_vec <- as.vector(t(costs))
times_vec <- as.vector(t(times))

Linear Programing formula

Define decision variables, where \[ x[i, j] \] is 1 if contractor i is assigned to building j, and 0 otherwise. The objective is to minimize the total cost and time.

Objective Function:

Minimize \[ sum(costs[i, j] * x[i, j]) + sum(times[i, j] * x[i, j]) \]

Constraints:

Use at least 2 contractors: \[ sum(y[i]) >= 2 \] where \[y[i]\] is 1 if contractor i is selected, 0 otherwise. Total cost: \[ sum(costs[i, j] * x[i, j]) <= budget_limit \] Total time: \[ sum(times[i, j] * x[i, j]) <= time_limit \] Each building must be assigned to exactly one contractor: \[ sum(x[i, j]) == 1 \] for all j Logical constraint: \[ x[i, j] <= y[i] \] to ensure if a contractor is not selected for any building, they are not used.

Solution

# Define the constraints matrix
num_vars <- num_contractors * num_buildings
constraints <- matrix(0, nrow = num_contractors + num_buildings + 2, ncol = num_vars + num_contractors)

# Contractor selection constraints
for (i in 1:num_contractors) {
  constraints[i, ((i-1)*num_buildings+1):(i*num_buildings)] <- 1
  constraints[i, num_vars + i] <- -1
}

# Building assignment constraints
for (j in 1:num_buildings) {
  constraints[num_contractors + j, seq(j, num_vars, by = num_buildings)] <- 1
}

# Total cost and time constraints
constraints[num_contractors + num_buildings + 1, 1:num_vars] <- costs_vec
constraints[num_contractors + num_buildings + 2, 1:num_vars] <- times_vec

# Right-hand side of constraints
rhs <- c(rep(0, num_contractors), rep(1, num_buildings), budget_limit, time_limit)

# Directions of the constraints
directions <- c(rep("<=", num_contractors), rep("=", num_buildings), "<=", "<=")

# Objective coefficients
objective <- c(costs_vec + times_vec, rep(0, num_contractors))

# Solve the ILP
solution <- lp("min", objective, constraints, directions, rhs, binary.vec = 1:(num_vars + num_contractors))

# Check if a solution is found

if (solution$status == 0) {
  selected_vars <- solution$solution[1:num_vars]
  selected_contractors <- solution$solution[(num_vars+1):(num_vars+num_contractors)]

  selected_vars_matrix <- matrix(selected_vars, nrow = num_contractors, byrow = TRUE)

  cat("Selected contractor-building assignments:\n")
  print(selected_vars_matrix)
  cat("Selected contractors:\n")
  print(selected_contractors)
  
  # Create a data frame for visualization
  assignments_df <- data.frame(
    Contractor = rep(1:num_contractors, each = num_buildings),
    Building = rep(1:num_buildings, times = num_contractors),
    Assigned = as.factor(selected_vars)
  )
  
  # Visualization of the assignments
  ggplot(assignments_df, aes(x = Building, y = Contractor, fill = Assigned)) +
    geom_tile() +
    scale_fill_manual(values = c("0" = "white", "1" = "blue"), labels = c("Not Assigned", "Assigned")) +
    labs(title = "Contractor-Building Assignments", x = "Building", y = "Contractor", fill = "Assignment") +
    theme(axis.text.x = element_text(angle = 45, hjust = 1))
  
} else {
  cat("No feasible solution found with the current constraints.\n")
}

Selected contractor-building assignments:
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    0    0    0    0
[2,]    0    0    0    0    1
[3,]    0    0    0    1    0
[4,]    0    0    1    0    0
[5,]    0    1    0    0    0
Selected contractors:
[1] 1 1 1 1 1

sensitivity_analysis <- function(budget_limits, time_limits) {
  results <- data.frame(Budget = numeric(), TimeLimit = numeric(), TotalCost = numeric(), TotalTime = numeric(), ContractorsUsed = numeric())

  for (budget in budget_limits) {
    for (time in time_limits) {
      rhs <- c(rep(0, num_contractors), rep(1, num_buildings), budget, time)

      solution <- lp("min", objective, constraints, directions, rhs, binary.vec = 1:(num_vars + num_contractors))
      
      if (solution$status == 0) {
        total_cost <- sum(costs_vec * solution$solution[1:num_vars])
        total_time <- sum(times_vec * solution$solution[1:num_vars])
        contractors_used <- sum(solution$solution[(num_vars+1):(num_vars+num_contractors)])
        
        results <- rbind(results, data.frame(Budget = budget, TimeLimit = time, TotalCost = total_cost, TotalTime = total_time, ContractorsUsed = contractors_used))
      }
    }
  }
  
  return(results)
}

budget_limits <- seq(500, 700, by = 50)
time_limits <- seq(50, 70, by = 5)
sensitivity_results <- sensitivity_analysis(budget_limits, time_limits)

# Print sensitivity results
print(sensitivity_results)

   Budget TimeLimit TotalCost TotalTime ContractorsUsed
1     600        55       575      55.0               5
2     600        60       570      56.5               5
3     600        65       570      56.5               5
4     600        70       570      56.5               5
5     650        55       575      55.0               5
6     650        60       570      56.5               5
7     650        65       570      56.5               5
8     650        70       570      56.5               5
9     700        55       575      55.0               5
10    700        60       570      56.5               5
11    700        65       570      56.5               5
12    700        70       570      56.5               5

If time limit is increased by 5 units, cost can also be reduced by 5 units (actual increase in time is 1.5 units)!