Assigning students to rank-ordered group projects with linear optimization

It is fairly common in a classroom setting to have students who must be assiged to projects, either to work alone or in groups. If the students have each ranked the projects from most to least desirable, one typically wants to assign all of the students to their most desired projects.

This is a perfect scenario for optimization, in that we can come up with a "satisfaction" metric based on the rankings and let the solver try to find the best configuration. In this post, we will set up the problem using the proprietary Gurobi solver via their Python bindings, but it can be just as easily set up using PuLP or OR-Tools.

Here, we will consider the problem of assigning students to groups of at most 4 for projects that they have ranked from 1 (most desirable) up to the number projects, and we will try to make everyone as happy as possible in the resulting solution. The data is anonymized but is real data from a recent class (though slightly modified to make the problem harder).

Input the data

Here's some data that was collected using a web form asking each student to rate each of the nine projects from 1 (most desired) to 9 (least desired). In our case, the form required a complete ordering from each student, but if that was not available we could easily create a constraint that prevents students from being assigned to a project they wanted to little that they refused to rate it.

In [1]:
import random
import pandas as pd
import gurobipy as gp
from gurobipy import GRB

# number of students
I = 16
students = [f"Student {i:02d}" for i in range(I)]

# number of projects
J = 9
projects = [f"Project {j:02d}" for j in range(J)]

# students rank each project
rank_df = pd.DataFrame([
        [6, 4, 1, 9, 8, 2, 7, 5, 3],
        [1, 2, 5, 3, 7, 6, 9, 4, 8],
        [4, 5, 6, 8, 9, 3, 2, 1, 7],
        [2, 4, 3, 7, 9, 1, 5, 8, 6],
        [3, 2, 5, 1, 9, 8, 4, 6, 7],
        [2, 1, 7, 4, 9, 8, 6, 3, 5],
        [5, 2, 3, 1, 9, 4, 8, 7, 6],
        [4, 1, 7, 2, 9, 3, 6, 5, 8],
        [8, 6, 3, 5, 2, 4, 7, 1, 9],
        [2, 6, 7, 3, 5, 1, 4, 9, 8],
        [1, 5, 3, 4, 9, 7, 2, 8, 6],
        [4, 5, 8, 7, 9, 3, 2, 1, 6],
        [7, 5, 8, 1, 6, 9, 4, 3, 2],
        [4, 5, 2, 9, 8, 1, 3, 7, 6],
        [3, 5, 2, 4, 8, 1, 9, 6, 7],
        [1, 3, 4, 2, 8, 5, 7, 6, 9]
    ],
    index=students,
    columns=projects,
)
rank_df
Out[1]:
Project 00 Project 01 Project 02 Project 03 Project 04 Project 05 Project 06 Project 07 Project 08
Student 00 6 4 1 9 8 2 7 5 3
Student 01 1 2 5 3 7 6 9 4 8
Student 02 4 5 6 8 9 3 2 1 7
Student 03 2 4 3 7 9 1 5 8 6
Student 04 3 2 5 1 9 8 4 6 7
Student 05 2 1 7 4 9 8 6 3 5
Student 06 5 2 3 1 9 4 8 7 6
Student 07 4 1 7 2 9 3 6 5 8
Student 08 8 6 3 5 2 4 7 1 9
Student 09 2 6 7 3 5 1 4 9 8
Student 10 1 5 3 4 9 7 2 8 6
Student 11 4 5 8 7 9 3 2 1 6
Student 12 7 5 8 1 6 9 4 3 2
Student 13 4 5 2 9 8 1 3 7 6
Student 14 3 5 2 4 8 1 9 6 7
Student 15 1 3 4 2 8 5 7 6 9

And here are the students who have indicated that they would specifically like to be together with or apart from another students.

In [2]:
together = [
    ("Student 01", "Student 13"),
    ("Student 05", "Student 06")
]
apart = [
    ("Student 00", "Student 05")
]

Setting up the optimization

This setup naturally leads to a classic assignment ILP (integer linear program), where each decision variable is a one (assign a certain student to a certain project) or a zero (not assigned). Here's how the ILP breaks down:

Data

  • List of students indexed by $i \in I$
  • List of projects indexed by $j \in J$
  • Student preferences $\mathrm{rank}(i,j) \in {1, \ldots, J}$ available $\forall i \in I,\, j \in J$ or else assumed to be infeasible (set equal to 100)

Decision variables

  • Assignment matrix (binary):
$$ X_{ij} = \begin{cases} 1 & \mathrm{if\ student\ } i \mathrm{\ assigned\ to\ project\ } j \\ 0 & \mathrm{otherwise} \end{cases} $$

Constraints

  • Each student is assigned to one and only one project
$$ \sum_j X_{ij} = 1 \qquad \forall \, i \in I$$
  • Each project has less than or equal to four students
$$ \sum_i X_{ij} \leq 4 \qquad \forall \, j \in J$$
  • If two students mutually request on another, they are put in the same group.

  • If any student requests not to work with another, they are not put together.

  • The number of projects is less than or equal to some value $P$ (see below for discussion).

Objective function

Place students in their most desired project, if possible. Higher (worse) preference rank results in a squared loss (making lower places much worse).

$$ \mathrm{minimize}\,Z = \sum_i \sum_j X_{ij} \cdot \left( \mathrm{rank}(i,j) - 1\right)^2 $$

Example:

  • student $i$ gets 1st choice $\implies Z_i = 0^2 = 0$
  • student $i$ gets 2nd choice $\implies Z_i = 1^2 = 1$
  • student $i$ gets 3rd choice $\implies Z_i = 2^2 = 4$
  • ...

Implementation

Given the problem definition above, we can start to set this up with code. Gurobi has a multidict object which allows for convenient summing and other relevant operations by indices:

In [3]:
permutations, ratings = gp.multidict({
    (i, j): rank_df.loc[i, j]
    for i in students
    for j in projects
})

The first object returned is a list of all of indices of decision variables, basically the $i$ and $j$ for each $X_{ij}$ but with the Python-native string keys instead of the actual indexes:

In [4]:
permutations[:15]
Out[4]:
[('Student 00', 'Project 00'),
 ('Student 00', 'Project 01'),
 ('Student 00', 'Project 02'),
 ('Student 00', 'Project 03'),
 ('Student 00', 'Project 04'),
 ('Student 00', 'Project 05'),
 ('Student 00', 'Project 06'),
 ('Student 00', 'Project 07'),
 ('Student 00', 'Project 08'),
 ('Student 01', 'Project 00'),
 ('Student 01', 'Project 01'),
 ('Student 01', 'Project 02'),
 ('Student 01', 'Project 03'),
 ('Student 01', 'Project 04'),
 ('Student 01', 'Project 05')]

The ratings part is the same dictionary comprehension defined above that has all of these permutations as dict keys and then the currently assigned values as the dict values, except as a multidict it now has some extra functions that come in handy for defining the ILP:

In [5]:
ratings.sum("*", "Project 00")
Out[5]:
<gurobi.LinExpr: 57.0>
In [6]:
ratings.sum("*", "Project 01")
Out[6]:
<gurobi.LinExpr: 61.0>

Now we'll use the Gurobi built-ins to set this up. Here's the twist — above we set a maximum number of projects $P$, except we're going to put the problem solving inside a function and let that be the one variable we specifcy dynamically.

The motivation for this is that you can imagine having a soft constraint on the number of projects while still having a preference for fewer, so it could be helpful to see a few scenarios of number of projects and how that affects the objective function.

This could be set as a variable in the optimization, but in many assignment problems there is value to putting a "human in the loop" so we can make informed decisions about whether it's worth having some groups of three and one extra project with only two if that means everyone is much happier about their assignments.

In [7]:
MAX_STUDENTS_PER_PROJECT = 4
MIN_STUDENTS_PER_PROJECT = 2

def solve(max_projects):
    # initialize the model object
    m = gp.Model(f"project_assignment_{max_projects}")

    assign = m.addVars(permutations, vtype=GRB.BINARY, name="assign")
    use_project = m.addVars(projects, vtype=GRB.BINARY, name="use_project")

    # each student has one and only one project group
    m.addConstrs(
        (assign.sum(student, "*") == 1 for student in students),
        name="EachStudentAssignedToOneProject"
    )

    # projects can't exceed the maximum number of students
    m.addConstrs(
        (assign.sum("*", project) <= MAX_STUDENTS_PER_PROJECT for project in projects),
        name="LimitGroupSize"
    )

    # projects must be considered 'in use' if any students are assigned
    m.addConstrs(
        (use_project[project] >= assign[(student, project)] for student in students for project in projects),
        name="ProjectMustBeInUseIfAnyStudentsAssigned"
    )

    # don't exceed max number of projects
    m.addConstr(use_project.sum() <= max_projects, name="MaxProjects")

    # if any students are assigned to a project, the project must have at least 2 students
    m.addConstrs(
        (assign.sum("*", project) >= use_project[project] * MIN_STUDENTS_PER_PROJECT for project in projects),
        name="ProjectsInUseMustHaveAtLeastTwoStudents"
    )

    # put students together who both indicated the other
    for student1, student2 in together:
        m.addConstrs(
            (assign[(student1, project)] == assign[(student2, project)] for project in projects),
            name=f"PairStudents[{(student1, student1)}]"
        )
    
    # keep students apart who contraindicated another
    for student1, student2 in apart:
        m.addConstrs(
            (
                (assign[(student1, project)] + assign[(student2, project)]
            ) <= 1 for project in projects),
            name=f"ApartStudents[{(student1, student1)}]"
        )

    # set the objective function to be minimized
    m.setObjective(
        (ratings.prod(assign) - 1) ** 2,
        sense=GRB.MINIMIZE,
    )

    m.optimize()
    return m, assign

A first pass

Let's solve this for a maximum of 4 projects:

In [8]:
m, assign = solve(max_projects=4)
Academic license - for non-commercial use only - expires 2021-03-24
Using license file /home/isaac/gurobi.lic
Gurobi Optimizer version 9.1.1 build v9.1.1rc0 (linux64)
Thread count: 4 physical cores, 8 logical processors, using up to 8 threads
Optimize a model with 206 rows, 153 columns and 792 nonzeros
Model fingerprint: 0x9eddfe84
Model has 10440 quadratic objective terms
Variable types: 0 continuous, 153 integer (153 binary)
Coefficient statistics:
  Matrix range     [1e+00, 2e+00]
  Objective range  [2e+00, 2e+01]
  QObjective range [2e+00, 3e+02]
  Bounds range     [1e+00, 1e+00]
  RHS range        [1e+00, 4e+00]
Found heuristic solution: objective 8100.0000000
Presolve removed 56 rows and 18 columns
Presolve time: 0.01s
Presolved: 150 rows, 135 columns, 639 nonzeros
Presolved model has 8001 quadratic objective terms
Variable types: 0 continuous, 135 integer (135 binary)

Root relaxation: objective 5.290000e+02, 86 iterations, 0.00 seconds

    Nodes    |    Current Node    |     Objective Bounds      |     Work
 Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time

     0     0  529.00000    0   10 8100.00000  529.00000  93.5%     -    0s
H    0     0                     676.0000000  529.00000  21.7%     -    0s
*    0     0               0     529.0000000  529.00000  0.00%     -    0s

Cutting planes:
  Gomory: 1
  Cover: 1
  Zero half: 1

Explored 1 nodes (93 simplex iterations) in 0.04 seconds
Thread count was 8 (of 8 available processors)

Solution count 3: 529 676 8100 

Optimal solution found (tolerance 1.00e-04)
Best objective 5.290000000000e+02, best bound 5.290000000000e+02, gap 0.0000%

This solved almost immediately, and in fairness this problem is small enough that a person could get pretty close to optimal just going by hand. But you can probably imagine how hard this would be for a person to do once the class size got to be in the dozens of students. Meanwhile, the solver will still be very fast for almost any reasonably sized classroom.

The solve function gives us back a model object and an assignment dictionary, but we'll set up a function to turn this into some helpful dataframes:

In [9]:
def get_results(assign):
    """ Take the dict of results and turn it into useful DataFrames """
    
    # create df with impossible placeholder
    assign_df = pd.DataFrame(-1, index=students, columns=projects)

    # fill in the decision variable results
    for (i, j), x_ij in assign.items():
        assign_df.loc[i, j] = int(x_ij.X)

    # sanity check that none were missed
    assert ((assign_df == 0) | (assign_df == 1)).all().all()
    
    # count how many students got their nth choice
    choices = (assign_df * rank_df).values.ravel()
    choices = choices[choices > 0]
    n_ranks = pd.Series(choices).value_counts().rename(index=lambda x: f"choice {x}")

    # count up how big the group sizes are
    group_sizes = assign_df.sum(axis="rows").sort_values(ascending=False).rename("n").sort_values(ascending=False)
    
    return assign_df, n_ranks, group_sizes

assign_df, n_ranks, group_sizes = get_results(assign)
assign_df
Out[9]:
Project 00 Project 01 Project 02 Project 03 Project 04 Project 05 Project 06 Project 07 Project 08
Student 00 0 0 0 0 0 1 0 0 0
Student 01 1 0 0 0 0 0 0 0 0
Student 02 0 0 0 0 0 0 0 1 0
Student 03 0 0 0 0 0 1 0 0 0
Student 04 0 1 0 0 0 0 0 0 0
Student 05 0 1 0 0 0 0 0 0 0
Student 06 0 1 0 0 0 0 0 0 0
Student 07 0 1 0 0 0 0 0 0 0
Student 08 0 0 0 0 0 0 0 1 0
Student 09 0 0 0 0 0 1 0 0 0
Student 10 1 0 0 0 0 0 0 0 0
Student 11 0 0 0 0 0 0 0 1 0
Student 12 0 0 0 0 0 0 0 1 0
Student 13 1 0 0 0 0 0 0 0 0
Student 14 0 0 0 0 0 1 0 0 0
Student 15 1 0 0 0 0 0 0 0 0
In [10]:
n_ranks
Out[10]:
choice 1    11
choice 2     3
choice 3     1
choice 4     1
dtype: int64

So most students got their first choice, a few got their second, and then a couple students were lower than that.

In [11]:
group_sizes
Out[11]:
Project 00    4
Project 01    4
Project 05    4
Project 07    4
Project 02    0
Project 03    0
Project 04    0
Project 06    0
Project 08    0
Name: n, dtype: int64

Meanwhile, as specified this time around we only used allowed four of the projects to be used, so with 16 students each project ended up with the maximum number of students.

Could we get better assignments if we allow another project instead of dense packing students into the maximum group size of 4 for each of 4 projects?

Solving again with a bit more slack

In [12]:
m, assign = solve(max_projects=5)
assign_df, n_ranks, group_sizes = get_results(assign)
Gurobi Optimizer version 9.1.1 build v9.1.1rc0 (linux64)
Thread count: 4 physical cores, 8 logical processors, using up to 8 threads
Optimize a model with 206 rows, 153 columns and 792 nonzeros
Model fingerprint: 0x403ba5fa
Model has 10440 quadratic objective terms
Variable types: 0 continuous, 153 integer (153 binary)
Coefficient statistics:
  Matrix range     [1e+00, 2e+00]
  Objective range  [2e+00, 2e+01]
  QObjective range [2e+00, 3e+02]
  Bounds range     [1e+00, 1e+00]
  RHS range        [1e+00, 5e+00]
Found heuristic solution: objective 8836.0000000
Presolve removed 56 rows and 18 columns
Presolve time: 0.01s
Presolved: 150 rows, 135 columns, 639 nonzeros
Presolved model has 8001 quadratic objective terms
Variable types: 0 continuous, 135 integer (135 binary)

Root relaxation: objective 4.000000e+02, 54 iterations, 0.00 seconds

    Nodes    |    Current Node    |     Objective Bounds      |     Work
 Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time

*    0     0               0     400.0000000  400.00000  0.00%     -    0s

Explored 0 nodes (54 simplex iterations) in 0.02 seconds
Thread count was 8 (of 8 available processors)

Solution count 2: 400 8836 

Optimal solution found (tolerance 1.00e-04)
Best objective 4.000000000000e+02, best bound 4.000000000000e+02, gap 0.0000%
In [13]:
n_ranks
Out[13]:
choice 1    13
choice 2     2
choice 4     1
dtype: int64
In [14]:
group_sizes
Out[14]:
Project 00    4
Project 05    4
Project 01    3
Project 07    3
Project 03    2
Project 02    0
Project 04    0
Project 06    0
Project 08    0
Name: n, dtype: int64

This way we have demonstrated that going from 4 to 5 projects spreads the students out without causing much improvement in satisfaction. At this point, a human instructor could decide whether it is worthwhile or not to take on supervising another project.

Conclusion

This is the type of problem that fits neatly into an optimization approach. It is extremely difficult to do by hand once the problem size grows past a handful of decision variables — try looking at the rank dataframe above and thinking about how you would do this manually if necessary. Meanwhile it is trivial to represent and solve as an ILP.

The approach can easily be extended to a situation where only one student can work on a given "thing" they have ranked and there are at least as many of the assignable things as students. For example, instead of groups they might have ranked physical resources or project topics they want to work on.

In the real world, assignment problems like this often need a "fudge factor." After all, you can optimize a problem like nurse scheduling or project assignment to 4 decimal places in the objective function, but people do messy things like calling in sick or dropping classes, and in practice our objective function isn't really one-to-one with actual happiness. Tiny changes rarely matter.

A better strategy than coming up with one perfect solution is to jump past all the hard work of figuring this out on paper but then put a human in the loop by showing a nunber of feasible and nearly-optimal solutions and let a person decide based on the totality of the circumstances.

Previous: Solving a logic puzzle with the Z3 theorem prover

Any comments or suggestions? Let me know.