{"id":1285,"date":"2016-01-06T07:43:12","date_gmt":"2016-01-06T12:43:12","guid":{"rendered":"https:\/\/blogs.mathworks.com\/loren\/?p=1285"},"modified":"2016-06-16T10:35:05","modified_gmt":"2016-06-16T15:35:05","slug":"generating-an-optimal-employee-work-schedule-using-integer-linear-programming","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/loren\/2016\/01\/06\/generating-an-optimal-employee-work-schedule-using-integer-linear-programming\/","title":{"rendered":"Generating an Optimal Employee Work Schedule Using Integer Linear Programming"},"content":{"rendered":"<div class=\"content\"><!--introduction--><p>Today's guest blogger is Teja Muppirala, who is a member of our Consulting Services group. Teja has been with MathWorks for 6 years and is based in our Tokyo office.<\/p><!--\/introduction--><h3>Contents<\/h3><div><ul><li><a href=\"#70d024fd-92d8-47cf-af03-8a6d398c131b\">Mixed-Integer Linear Programming and The Nurse Scheduling Problem<\/a><\/li><li><a href=\"#264d84bd-9e77-4595-a91e-8e1a972b6d5d\">Problem Statement<\/a><\/li><li><a href=\"#4ddc1667-693d-4808-9ab3-89706de94079\">1. Read in the requirements and employee data from the Excel sheet<\/a><\/li><li><a href=\"#04a71d52-4161-40f1-a92c-60ba31b460f4\">2. Generate the f, A, and b matrices based on the the constraints and objectives<\/a><\/li><li><a href=\"#d4520395-e909-4164-8d45-ccdb3ca6b8e6\">3. Call intlinprog with every variable as an integer 0 or 1<\/a><\/li><li><a href=\"#8821edca-c7b6-4b3b-8882-25610f5f8653\">4. Gather and display the results<\/a><\/li><li><a href=\"#d44f2e37-6892-40b6-84a6-065d703446f6\">5. Conclusion<\/a><\/li><\/ul><\/div><h4>Mixed-Integer Linear Programming and The Nurse Scheduling Problem<a name=\"70d024fd-92d8-47cf-af03-8a6d398c131b\"><\/a><\/h4><p>Since it's introduction in release R2014a, we've had <a href=\"https:\/\/blogs.mathworks.com\/loren\/?s=intlinprog\">several blog posts now<\/a> showing some applications of <a href=\"https:\/\/www.mathworks.com\/help\/optim\/ug\/intlinprog.html\"><tt>intlinprog<\/tt><\/a> , the mixed-integer linear programming (MILP) solver found in the <a href=\"https:\/\/www.mathworks.com\/help\/optim\/index.html\">Optimization Toolbox<\/a>. This had been one of our most requested features, as MILP has trememdous application in a variety of areas such as scheduling and resource optimization.<\/p><p>There are a number of classic optimization problems that can be framed using MILP, such as the <a href=\"https:\/\/www.mathworks.com\/help\/optim\/ug\/travelling-salesman-problem.html\">Travelling Salesman Problem<\/a> and <a href=\"https:\/\/blogs.mathworks.com\/pick\/2014\/03\/28\/solve-it-and-solve-it-right\/\">Knapsack Problem<\/a>.<\/p><p>The work in this blog post is based loosely on a discussion I recently had with a customer, who wanted to make an optimal shift schedule for his employees while satisfying certain availability and staffing constraints. This is a variation of what is known as the <a href=\"https:\/\/en.wikipedia.org\/wiki\/Nurse_scheduling_problem\">Nurse Scheduling Problem.<\/a><\/p><h4>Problem Statement<a name=\"264d84bd-9e77-4595-a91e-8e1a972b6d5d\"><\/a><\/h4><p>Obviously the \"Nurse Scheduling Problem\" is not limited to \"nurses\" as an occupation, so I will just use the generic term \"employee\" here. The customer's problem statement was as follows.<\/p><p><b>Given:<\/b><\/p><div><ul><li>A list of employees with their available work hours, and hourly   salaries<\/li><li>A prescibed minimum number of staff needed to be at work at a given hour   (fewer staff are needed at night, more staff are needed during   peak hours)<\/li><\/ul><\/div><p><b>Find a schedule which MINIMIZES:<\/b><\/p><div><ul><li>The total daily wages the employer must pay out.<\/li><\/ul><\/div><p><b>While meeting the hard CONSTRAINTS:<\/b><\/p><div><ul><li>At any given hour, you must meet the minimum staffing requirement<\/li><\/ul><\/div><div><ul><li>An employee can only work one shift a day<\/li><\/ul><\/div><div><ul><li>An employee must work within his available hours<\/li><\/ul><\/div><div><ul><li>If an employee is called for duty, they must work at least a specified minimum number of hours, and no more than a specified maximum number of hours<\/li><\/ul><\/div><p>This is probably a bit abstract, but let's load in some concrete data so you can get a clearer idea of what exactly we are trying to do here.<\/p><h4>1. Read in the requirements and employee data from the Excel sheet<a name=\"4ddc1667-693d-4808-9ab3-89706de94079\"><\/a><\/h4><p>The staff information and scheduling requirements are in an Excel file, so we will first need to import the data. This is a fictional dataset and the Excel file is available along with the MATLAB code for this blog post so you can feel free to try it out on your own.<\/p><p>The first sheet in the Excel file (available <a href=\"https:\/\/blogs.mathworks.com\/images\/loren\/2016\/scheduling.xlsx\">here<\/a>) contains the staff information, and it is in a tabular format suitable to be imported directly as a MATLAB <a href=\"https:\/\/www.mathworks.com\/help\/matlab\/ref\/table.html\"><tt>table<\/tt><\/a> using the <tt>readtable<\/tt> function. Let's take a look at what's inside.<\/p><pre class=\"codeinput\">filename = <span class=\"string\">'scheduling.xlsx'<\/span>;\r\nstaffTable = readtable(filename)\r\n<\/pre><pre class=\"codeoutput\">staffTable = \r\n    EmployeeName    MinHours    MaxHours    HourlyWage    Availability\r\n    ____________    ________    ________    __________    ____________\r\n    'SMITH'         6           8           30            '6-20'      \r\n    'JOHNSON'       6           8           50            ''          \r\n    'WILLIAMS'      6           8           30            ''          \r\n    'JONES'         6           8           30            ''          \r\n    'BROWN'         6           8           40            ''          \r\n    'DAVIS'         6           8           50            ''          \r\n    'MILLER'        6           8           45            '6-18'      \r\n    'WILSON'        6           8           30            ''          \r\n    'MOORE'         6           8           35            ''          \r\n    'TAYLOR'        6           8           40            ''          \r\n    'ANDERSON'      2           3           60            '0-6, 18-24'\r\n    'THOMAS'        2           4           40            ''          \r\n    'JACKSON'       2           4           60            '8-16'      \r\n    'WHITE'         2           6           55            ''          \r\n    'HARRIS'        2           6           45            ''          \r\n    'MARTIN'        2           3           40            ''          \r\n    'THOMPSON'      2           5           50            '12-24'     \r\n    'GARCIA'        2           4           50            ''          \r\n    'MARTINEZ'      2           4           40            ''          \r\n    'ROBINSON'      2           5           50            ''          \r\n<\/pre><p>The variable <tt>staffTable<\/tt> has a list of each employee, along with the minimum hours they must work (if they are called in), maximum hours they may work, their hourly wage, and any limits on availability if there are any.<\/p><p>For example, the first employee SMITH, if called in, must work at least 6 hours, no more than 8 hours, earns $30\/hour, and is only available between 6am and 8pm.<\/p><p>We can also see the required staffing requirements for each hour of the 24-hour day. I don't need this to be table, I just want to import it as a numeric array, so <tt>xlsread<\/tt> will work just fine.<\/p><pre class=\"codeinput\">requirements = xlsread(filename,2); <span class=\"comment\">% Second sheet has minimum staff requirements<\/span>\r\ndisp(<span class=\"string\">'Required staffing (hourly from 0:00 - 24:00): '<\/span>);\r\ndisp(mat2str(requirements(2,:)));\r\n<\/pre><pre class=\"codeoutput\">Required staffing (hourly from 0:00 - 24:00): \r\n[1 1 2 3 6 6 7 8 9 8 8 8 7 6 6 5 5 4 4 3 2 2 2 2]\r\n<\/pre><p>Late at night, only 1 or 2 employees are needed, while during peak hours in the morning to afternoon, we may need as many as 9 employees on duty.<\/p><h4>2. Generate the f, A, and b matrices based on the the constraints and objectives<a name=\"04a71d52-4161-40f1-a92c-60ba31b460f4\"><\/a><\/h4><p>Generating a MILP formulation of a particular problem involves expressing the minimization objective and constraints using linear equations, and these are typically written using matrix notation. The specifics of this are covered thoroughly in the <a href=\"https:\/\/www.mathworks.com\/help\/optim\/ug\/intlinprog.html\">documentation<\/a>.<\/p><p>There is a bit of machinery involved in generating the appropriate matrices in this particular case, and it would be a bit dense to go through all the details in this blog post. My aim is more to show an application of MILP, rather than discussing the intricate details of problem formulation. If you are interested in that though, I definitely encourage you to <a href=\"https:\/\/blogs.mathworks.com\/images\/loren\/2016\/makeMILPMatrices.m\">download<\/a> the MATLAB files for this example and check it out for yourself.<\/p><p>However, I do want to discuss briefly the structure of the constraint matrices and decision variables for this particular problem. Let's take a look at the constraint matrix <b>A<\/b>. Here I used the <tt>spy<\/tt> function to view the sparsity pattern.<\/p><pre class=\"codeinput\"><span class=\"comment\">% f,A,b are inputs for the MILP Problem<\/span>\r\n<span class=\"comment\">% staffNumberVector is needed later when plotting the results<\/span>\r\n[f,A,b,staffNumberVector] = makeMILPMatrices(staffTable,requirements);\r\n\r\nfigure,\r\nspy(A,1); axis <span class=\"string\">fill<\/span>;\r\ntitle({<span class=\"string\">'Sparsity pattern of the constraint matrix A:'<\/span>;\r\n       <span class=\"string\">'The top portion is from the minimum hour requirement,'<\/span>;\r\n       <span class=\"string\">'The bottom portion is from the one-shift-a-day requirement'<\/span>});\r\n<\/pre><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/images\/loren\/2016\/shiftScheduleDemoBlog_01.png\" alt=\"\"> <p>There are two parts to the inequality constraint matrix <b>A<\/b>, which are clearly visible. The top portion has 24 rows (because there are 24 hours in a day) and each row represents the constraint that at a particular hour, you need a minimum number of staff.<\/p><p>The bottom section of <b>A<\/b>, in this case, has 20 rows (because there are 20 available employees). The k-th row of this section is associated with the constraint that the k-th employee may only work one shift per day.<\/p><p><b>A<\/b> will in general be much wider than tall, and the columns of <b>A<\/b> represent every possible shift for every possible employee. Every column of <b>A<\/b> is also associated with a decision variable that we constrain to be either 0 or 1 (a binary variable).<\/p><p>$\\textbf{x} = [x_1,x_2,...x_k]$<\/p><p>${x_n} \\in \\{0,1\\}$ for all n<\/p><p>A value of 1 for $x_n$ means that that n-th shift possibility is implemented, and a value of 0 indicates it is not. For this particular problem, <b>A<\/b> has 1295 columns, and so <b>x<\/b> is a decision vector of length k = 1295. It is this <b>x<\/b> that <tt>intlinprog<\/tt> will solve for.<\/p><p>As the number of employees and possible shifts increases, <b>A<\/b> may consist of many thousands of columns. Although I didn't use sparse matrices in this demo, <tt>intlinprog<\/tt> is capable of working with sparse matrices to make more efficient use of memory and allow for larger problem sizes.<\/p><h4>3. Call intlinprog with every variable as an integer 0 or 1<a name=\"d4520395-e909-4164-8d45-ccdb3ca6b8e6\"><\/a><\/h4><p>Once the relevant matrices are created, all that's left to do is actually call <tt>intlinprog<\/tt>. I specify that all the decision variables are integers, and have a value between 0 and 1.<\/p><pre class=\"codeinput\">nVars = numel(f);  <span class=\"comment\">% Number of variables in the problem<\/span>\r\nlb = zeros(nVars,1);    <span class=\"comment\">% The upper and lower bounds are 0 and 1<\/span>\r\nub = ones(nVars,1);\r\n[x, cost] = intlinprog(f,1:nVars,A,b,[],[],lb,ub);\r\n<\/pre><pre class=\"codeoutput\">LP:                Optimal objective value is 4670.000000.                                          \r\n\r\nCut Generation:    Applied 1 zero-half cut.                                                         \r\n                   Lower bound is 4670.000000.                                                      \r\n\r\nBranch and Bound:\r\n\r\n   nodes     total   num int        integer       relative                                          \r\nexplored  time (s)  solution           fval        gap (%)                                         \r\n      49      0.07         1   4.670000e+03   0.000000e+00                                          \r\n\r\nOptimal solution found.\r\n\r\nIntlinprog stopped because the objective value is within a gap tolerance of the\r\noptimal value, options.TolGapAbs = 0 (the default value). The intcon variables\r\nare integer within tolerance, options.TolInteger = 1e-05 (the default value).\r\n\r\n<\/pre><p>I am returning two arguments, the first is the decision vector itself, <tt>x<\/tt>, and the second argument, <tt>cost<\/tt>, is the minimized daily cost associated with this solution. You could also find <tt>cost<\/tt> as<\/p><p><tt>cost = f*x<\/tt>.<\/p><p><tt>intlinprog<\/tt> by default outputs certain diagnostic information to the command window (though there are options to hide this if preferred). For example, the miminized cost is shown in the command window to be 4670.0 (in real-world terms, $4760 dollars in daily wages). Also, one can see the total time that it took to solve this problem.<\/p><p><tt>intlinprog<\/tt> is <i>very fast!<\/i><\/p><h4>4. Gather and display the results<a name=\"8821edca-c7b6-4b3b-8882-25610f5f8653\"><\/a><\/h4><p><tt>intlinprog<\/tt> will give me my optimal <tt>x<\/tt> vector, but this is simply a long vector of of ones and zeros saying which shifts are implemented. I will go back and find out which index of <tt>x<\/tt> corresponds with which employee and shift, and then display the data in a more human-friendly form.<\/p><pre class=\"codeinput\">numStaff = size(staffTable,1); <span class=\"comment\">% Total number of staff available<\/span>\r\n\r\n<span class=\"comment\">% Convert from indices in x to employee and shift information<\/span>\r\nselected = find(x);\r\nhoursMatrix = zeros(numStaff,24);\r\n\r\n<span class=\"keyword\">for<\/span> n = 1:numel(selected);\r\n    thisEntry =selected(n);\r\n    thisStaff = staffNumberVector(thisEntry);\r\n    hoursOnDuty = -A(1:24,thisEntry);\r\n    hoursMatrix(thisStaff,:) = hoursOnDuty;\r\n<span class=\"keyword\">end<\/span>\r\n\r\n<span class=\"comment\">% Make a figure in a suitable location<\/span>\r\nhf = figure(<span class=\"string\">'visible'<\/span>,<span class=\"string\">'off'<\/span>,<span class=\"string\">'units'<\/span>,<span class=\"string\">'pixel'<\/span>,<span class=\"string\">'position'<\/span>,[100 100 560 700]);\r\nmovegui(hf,<span class=\"string\">'center'<\/span>);\r\nset(hf,<span class=\"string\">'visible'<\/span>,<span class=\"string\">'on'<\/span>);\r\n\r\n<span class=\"comment\">% Draw the employee work schedule matrix<\/span>\r\nsubplot(211);\r\nimagesc(0.5:23.5,1:20,hoursMatrix)\r\nset(gca,<span class=\"string\">'xtick'<\/span>,0:24);\r\nset(gca,<span class=\"string\">'ytick'<\/span>,1:numStaff,<span class=\"string\">'yticklabel'<\/span>,staffTable.EmployeeName);\r\naxis <span class=\"string\">tight<\/span>;\r\nset(gca,<span class=\"string\">'xlimmode'<\/span>,<span class=\"string\">'manual'<\/span>,<span class=\"string\">'ylimmode'<\/span>,<span class=\"string\">'manual'<\/span>,<span class=\"string\">'layer'<\/span>,<span class=\"string\">'bottom'<\/span>);\r\nhold <span class=\"string\">all<\/span>;\r\nylabel(<span class=\"string\">'Employee Name'<\/span>,<span class=\"string\">'FontSize'<\/span>,12);\r\n\r\n<span class=\"comment\">% Let's make the grid lines manually...<\/span>\r\n<span class=\"keyword\">for<\/span> n = 0.5:numStaff+0.5\r\n    plot(xlim,[n n],<span class=\"string\">'k'<\/span>);\r\n<span class=\"keyword\">end<\/span>\r\n<span class=\"keyword\">for<\/span> n = 0:24\r\n    plot([n n], ylim,<span class=\"string\">'k'<\/span>);\r\n<span class=\"keyword\">end<\/span>\r\n\r\ntitle([<span class=\"string\">'Employee shift schedule'<\/span> 10 <span class=\"keyword\">...<\/span>\r\n       <span class=\"string\">'Total wages over 24 hours: $'<\/span> num2str(cost)],<span class=\"string\">'FontSize'<\/span>,16);\r\ncolormap([1 1 1; 0 0.5 1]);\r\n\r\n<span class=\"comment\">% Draw the required vs. actual staff<\/span>\r\nsubplot(212);\r\nplot(0.5:23.5,requirements(2,:),<span class=\"string\">'b.-'<\/span>,<span class=\"string\">'linewidth'<\/span>,3,<span class=\"string\">'markersize'<\/span>,30);\r\nactualHours = -A(1:24,:)*x;\r\nhold <span class=\"string\">on<\/span>;\r\nplot(0.5:23.5,actualHours,<span class=\"string\">'r:.'<\/span>,<span class=\"string\">'linewidth'<\/span>,2,<span class=\"string\">'markersize'<\/span>,16);\r\nxlim([0 24]);\r\nset(gca,<span class=\"string\">'xtick'<\/span>,0:24,<span class=\"string\">'ytick'<\/span>,0:10);\r\nylim([0 10]);\r\ngrid <span class=\"string\">on<\/span>;\r\ntitle([<span class=\"string\">'Hourly Staff Count'<\/span>],<span class=\"string\">'FontSize'<\/span>,16);\r\nlegend({<span class=\"string\">'Employees Required'<\/span>, <span class=\"string\">'Employees Scheduled'<\/span>});\r\nxlabel(<span class=\"string\">'Time of day'<\/span>,<span class=\"string\">'fontsize'<\/span>,12);\r\nylabel(<span class=\"string\">'Employee count'<\/span>,<span class=\"string\">'fontsize'<\/span>,12);\r\n<\/pre><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/images\/loren\/2016\/shiftScheduleDemoBlog_02.png\" alt=\"\"> <p>We have two plots here. The first shows the actual work schedule for each employee (blue indicates when on shift). Note that there are a couple of employees (ANDERSON and JACKSON) that do not get called in. Also, while it seems that JOHNSON has two shifts, this is because we are wrapping around and the shift actually goes from 10pm to 6am the next morning.<\/p><p>The second plot shows the number of staff on duty as compared to our minimum staffing requirements. We see that we meet the requirements exactly, without overstaffing or understaffing at any point. This makes sense, because if the goal is to minimize paid wages, then you'd try to stay as close to the minimum staff as possible.<\/p><h4>5. Conclusion<a name=\"d44f2e37-6892-40b6-84a6-065d703446f6\"><\/a><\/h4><p>The problem outlined in the example is not a trivial one to solve. To be sure, it takes a bit of study and experience to be able to know how to easily convert a real world problem into its equivalent MILP formulation.<\/p><p>But for problems like this, especially at larger problem sizes, attempting to solve it using other methods such trial-and-error or genetic algorithms will usually take longer, and generally not yield as good results.<\/p><p>Have you tried using MATLAB's mixed-integer linear programming solver for any of your own work? If you want to share your thoughts, or if you have any questions regarding this example, please feel free to let us know in <a href=\"https:\/\/blogs.mathworks.com\/loren\/?p=1285#respond\">here<\/a>.<\/p><script language=\"JavaScript\"> <!-- \r\n    function grabCode_929c965db0ed419baab512796124b3e0() {\r\n        \/\/ Remember the title so we can use it in the new page\r\n        title = document.title;\r\n\r\n        \/\/ Break up these strings so that their presence\r\n        \/\/ in the Javascript doesn't mess up the search for\r\n        \/\/ the MATLAB code.\r\n        t1='929c965db0ed419baab512796124b3e0 ' + '##### ' + 'SOURCE BEGIN' + ' #####';\r\n        t2='##### ' + 'SOURCE END' + ' #####' + ' 929c965db0ed419baab512796124b3e0';\r\n    \r\n        b=document.getElementsByTagName('body')[0];\r\n        i1=b.innerHTML.indexOf(t1)+t1.length;\r\n        i2=b.innerHTML.indexOf(t2);\r\n \r\n        code_string = b.innerHTML.substring(i1, i2);\r\n        code_string = code_string.replace(\/REPLACE_WITH_DASH_DASH\/g,'--');\r\n\r\n        \/\/ Use \/x3C\/g instead of the less-than character to avoid errors \r\n        \/\/ in the XML parser.\r\n        \/\/ Use '\\x26#60;' instead of '<' so that the XML parser\r\n        \/\/ doesn't go ahead and substitute the less-than character. \r\n        code_string = code_string.replace(\/\\x3C\/g, '\\x26#60;');\r\n\r\n        copyright = 'Copyright 2016 The MathWorks, Inc.';\r\n\r\n        w = window.open();\r\n        d = w.document;\r\n        d.write('<pre>\\n');\r\n        d.write(code_string);\r\n\r\n        \/\/ Add copyright line at the bottom if specified.\r\n        if (copyright.length > 0) {\r\n            d.writeln('');\r\n            d.writeln('%%');\r\n            if (copyright.length > 0) {\r\n                d.writeln('% _' + copyright + '_');\r\n            }\r\n        }\r\n\r\n        d.write('<\/pre>\\n');\r\n\r\n        d.title = title + ' (MATLAB code)';\r\n        d.close();\r\n    }   \r\n     --> <\/script><p style=\"text-align: right; font-size: xx-small; font-weight:lighter;   font-style: italic; color: gray\"><br><a href=\"javascript:grabCode_929c965db0ed419baab512796124b3e0()\"><span style=\"font-size: x-small;        font-style: italic;\">Get \r\n      the MATLAB code <noscript>(requires JavaScript)<\/noscript><\/span><\/a><br><br>\r\n      Published with MATLAB&reg; R2015b<br><\/p><\/div><!--\r\n929c965db0ed419baab512796124b3e0 ##### SOURCE BEGIN #####\r\n%% Generating an Optimal Employee Work Schedule Using Integer Linear Programming\r\n% Today's guest blogger is Teja Muppirala, who is a member of our\r\n% Consulting Services group. Teja has been with MathWorks for 6 years and\r\n% is based in our Tokyo office.\r\n\r\n%% Mixed-Integer Linear Programming and The Nurse Scheduling Problem\r\n% Since its introduction in release R2014a, we've had\r\n% <https:\/\/blogs.mathworks.com\/loren\/?s=intlinprog several blog posts now>\r\n% showing some applications of\r\n% <https:\/\/www.mathworks.com\/help\/optim\/ug\/intlinprog.html |intlinprog|> ,\r\n% the mixed-integer linear programming (MILP) solver found in the\r\n% <https:\/\/www.mathworks.com\/help\/optim\/index.html Optimization Toolbox>.\r\n% This had been one of our most requested features, as MILP has trememdous\r\n% application in a variety of areas such as scheduling and resource\r\n% optimization.\r\n%\r\n% There are a number of classic optimization problems that can be framed\r\n% using MILP, such as the\r\n% <https:\/\/www.mathworks.com\/help\/optim\/ug\/travelling-salesman-problem.html\r\n% Travelling Salesman Problem> and\r\n% <https:\/\/blogs.mathworks.com\/pick\/2014\/03\/28\/solve-it-and-solve-it-right\/\r\n% Knapsack Problem>.\r\n%\r\n% The work in this blog post is based loosely on a discussion I recently\r\n% had with a customer, who wanted to make an optimal shift schedule for his\r\n% employees while satisfying certain availability and staffing constraints.\r\n% This is a variation of what is known as the\r\n% <https:\/\/en.wikipedia.org\/wiki\/Nurse_scheduling_problem Nurse Scheduling\r\n% Problem.>\r\n%\r\n\r\n%% Problem Statement\r\n% Obviously the \"Nurse Scheduling Problem\" is not limited to \"nurses\" as an\r\n% occupation, so I will just use the generic term \"employee\" here.\r\n% The customer's problem statement was as follows.\r\n\r\n%%\r\n% *Given:*\r\n%\r\n% * A list of employees with their available work hours, and hourly\r\n%   salaries\r\n% * A prescibed minimum number of staff needed to be at work at a given hour \r\n%   (fewer staff are needed at night, more staff are needed during\r\n%   peak hours)\r\n% \r\n\r\n%% \r\n% *Find a schedule which MINIMIZES:*\r\n%\r\n% * The total daily wages the employer must pay out.\r\n%\r\n\r\n%%\r\n% *While meeting the hard CONSTRAINTS:*\r\n%\r\n% * At any given hour, you must meet the minimum staffing requirement\r\n%   \r\n% * An employee can only work one shift a day\r\n%\r\n% * An employee must work within his available hours\r\n%\r\n% * If an employee is called for duty, they must work at least a specified\r\n% minimum number of hours, and no more than a specified maximum number of\r\n% hours\r\n% \r\n% This is probably a bit abstract, but let's load in some concrete data so\r\n% you can get a clearer idea of what exactly we are trying to do here.\r\n\r\n%% 1. Read in the requirements and employee data from the Excel sheet\r\n% The staff information and scheduling requirements are in an Excel file,\r\n% so we will first need to import the data. This is a fictional dataset\r\n% and the Excel file is available along with the MATLAB code for this blog \r\n% post so you can feel free to try it out on your own.\r\n%%\r\n% The first sheet in the Excel file (available\r\n% <https:\/\/blogs.mathworks.com\/images\/loren\/2016\/scheduling.xlsx here>)\r\n% contains the staff information, and it is in a tabular format suitable to\r\n% be imported directly as a MATLAB\r\n% <https:\/\/www.mathworks.com\/help\/matlab\/ref\/table.html |table|> using the\r\n% |readtable| function. Let's take a look at what's inside.\r\n\r\nfilename = 'scheduling.xlsx';\r\nstaffTable = readtable(filename)\r\n\r\n%%\r\n% The variable |staffTable| has a list of each employee, along with the\r\n% minimum hours they must work (if they are called in), maximum hours they\r\n% may work, their hourly wage, and any limits on availability if there are\r\n% any.\r\n%\r\n% For example, the first employee SMITH, if called in, must work at least 6\r\n% hours, no more than 8 hours, earns $30\/hour, and is only available\r\n% between 6am and 8pm.\r\n\r\n%%\r\n% We can also see the required staffing requirements for each hour of the\r\n% 24-hour day. I don't need this to be table, I just want to import it as a\r\n% numeric array, so |xlsread| will work just fine.\r\n\r\nrequirements = xlsread(filename,2); % Second sheet has minimum staff requirements\r\ndisp('Required staffing (hourly from 0:00 - 24:00): ');\r\ndisp(mat2str(requirements(2,:)));\r\n\r\n\r\n%%\r\n% Late at night, only 1 or 2 employees are needed, while during peak hours\r\n% in the morning to afternoon, we may need as many as 9 employees on duty.\r\n\r\n%% 2. Generate the f, A, and b matrices based on the the constraints and objectives\r\n%\r\n% Generating a MILP formulation of a particular problem involves expressing\r\n% the minimization objective and constraints using linear equations, and\r\n% these are typically written using matrix notation. The specifics of this\r\n% are covered thoroughly in the <https:\/\/www.mathworks.com\/help\/optim\/ug\/intlinprog.html\r\n% documentation>. \r\n\r\n%%\r\n% There is a bit of machinery involved in generating the appropriate\r\n% matrices in this particular case, and it would be a bit dense to go\r\n% through all the details in this blog post. My aim is more to show an\r\n% application of MILP, rather than discussing the intricate details of\r\n% problem formulation. If you are interested in that though, I definitely\r\n% encourage you to\r\n% <https:\/\/blogs.mathworks.com\/images\/loren\/2016\/makeMILPMatrices.m\r\n% download> the MATLAB files for this example and check it out for\r\n% yourself.\r\n\r\n%%\r\n% However, I do want to discuss briefly the structure of the constraint \r\n% matrices and decision variables for this particular problem. Let's take a\r\n% look at the constraint matrix *A*. Here I used the |spy| function to view \r\n% the sparsity pattern.\r\n\r\n% f,A,b are inputs for the MILP Problem\r\n% staffNumberVector is needed later when plotting the results\r\n[f,A,b,staffNumberVector] = makeMILPMatrices(staffTable,requirements);\r\n\r\nfigure,\r\nspy(A,1); axis fill;\r\ntitle({'Sparsity pattern of the constraint matrix A:';\r\n       'The top portion is from the minimum hour requirement,';\r\n       'The bottom portion is from the one-shift-a-day requirement'});\r\n   \r\n%%\r\n% There are two parts to the inequality constraint matrix *A*, \r\n% which are clearly visible. The top portion has 24 rows (because there are\r\n% 24 hours in a day) and each row represents the constraint that at a\r\n% particular hour, you need a minimum number of staff.\r\n%\r\n% The bottom section of *A*, in this case, has 20 rows (because there are \r\n% 20 available employees). The k-th row of this section is associated with the \r\n% constraint that the k-th employee may only work one shift per day.\r\n%\r\n% *A* will in general be much wider than tall, and the columns of *A* \r\n% represent every possible shift for every possible employee. Every column \r\n% of *A* is also associated with a decision variable that we constrain to\r\n% be either 0 or 1 (a binary variable).\r\n%\r\n% $\\textbf{x} = [x_1,x_2,...x_k]$\r\n%\r\n% ${x_n} \\in \\{0,1\\}$ for all n\r\n%\r\n% A value of 1 for $x_n$ means that that n-th shift possibility is \r\n% implemented, and a value of 0 indicates it is not. For this particular \r\n% problem, *A* has 1295 columns, and so *x* is a decision vector of length \r\n% k = 1295. It is this *x* that |intlinprog| will solve for.\r\n%\r\n% As the number of employees and possible shifts increases, *A* may consist\r\n% of many thousands of columns. Although I didn't use sparse matrices in\r\n% this demo, |intlinprog| is capable of working with sparse matrices to \r\n% make more efficient use of memory and allow for larger problem sizes.\r\n\r\n%% 3. Call intlinprog with every variable as an integer 0 or 1\r\n% Once the relevant matrices are created, all that's left to do is actually\r\n% call |intlinprog|. I specify that all the decision variables are\r\n% integers, and have a value between 0 and 1.\r\n\r\nnVars = numel(f);  % Number of variables in the problem\r\nlb = zeros(nVars,1);    % The upper and lower bounds are 0 and 1\r\nub = ones(nVars,1);\r\n[x, cost] = intlinprog(f,1:nVars,A,b,[],[],lb,ub);\r\n\r\n%%\r\n% I am returning two arguments, the first is the decision vector itself,\r\n% |x|, and the second argument, |cost|, is the minimized daily cost \r\n% associated with this solution. You could also find |cost| as \r\n%\r\n% |cost = f*x|.\r\n% \r\n% |intlinprog| by default outputs certain diagnostic information to the command\r\n% window (though there are options to hide this if preferred). For \r\n% example, the miminized cost is shown in the command window to be 4670.0\r\n% (in real-world terms, $4760 dollars in daily wages). Also, one can see \r\n% the total time that it took to solve this problem.\r\n%\r\n% |intlinprog| is _very fast!_\r\n\r\n%% 4. Gather and display the results\r\n% |intlinprog| will give me my optimal |x| vector, but this is simply a\r\n% long vector of of ones and zeros saying which shifts are implemented. I\r\n% will go back and find out which index of |x| corresponds with which\r\n% employee and shift, and then display the data in a more human-friendly\r\n% form.\r\n\r\n\r\nnumStaff = size(staffTable,1); % Total number of staff available\r\n\r\n% Convert from indices in x to employee and shift information\r\nselected = find(x);\r\nhoursMatrix = zeros(numStaff,24);\r\n\r\nfor n = 1:numel(selected);\r\n    thisEntry =selected(n);\r\n    thisStaff = staffNumberVector(thisEntry);\r\n    hoursOnDuty = -A(1:24,thisEntry);\r\n    hoursMatrix(thisStaff,:) = hoursOnDuty;\r\nend\r\n\r\n% Make a figure in a suitable location\r\nhf = figure('visible','off','units','pixel','position',[100 100 560 700]);\r\nmovegui(hf,'center');\r\nset(hf,'visible','on');\r\n\r\n% Draw the employee work schedule matrix\r\nsubplot(211);\r\nimagesc(0.5:23.5,1:20,hoursMatrix)\r\nset(gca,'xtick',0:24);\r\nset(gca,'ytick',1:numStaff,'yticklabel',staffTable.EmployeeName);\r\naxis tight; \r\nset(gca,'xlimmode','manual','ylimmode','manual','layer','bottom');\r\nhold all;\r\nylabel('Employee Name','FontSize',12);\r\n\r\n% Let's make the grid lines manually...\r\nfor n = 0.5:numStaff+0.5 \r\n    plot(xlim,[n n],'k');\r\nend\r\nfor n = 0:24\r\n    plot([n n], ylim,'k');\r\nend\r\n\r\ntitle(['Employee shift schedule' 10 ...\r\n       'Total wages over 24 hours: $' num2str(cost)],'FontSize',16);\r\ncolormap([1 1 1; 0 0.5 1]);\r\n\r\n% Draw the required vs. actual staff\r\nsubplot(212);\r\nplot(0.5:23.5,requirements(2,:),'b.-','linewidth',3,'markersize',30);\r\nactualHours = -A(1:24,:)*x;\r\nhold on;\r\nplot(0.5:23.5,actualHours,'r:.','linewidth',2,'markersize',16);\r\nxlim([0 24]);\r\nset(gca,'xtick',0:24,'ytick',0:10);\r\nylim([0 10]);\r\ngrid on;\r\ntitle(['Hourly Staff Count'],'FontSize',16);\r\nlegend({'Employees Required', 'Employees Scheduled'});\r\nxlabel('Time of day','fontsize',12);\r\nylabel('Employee count','fontsize',12);\r\n\r\n%%\r\n% We have two plots here. The first shows the actual work schedule for each\r\n% employee (blue indicates when on shift). Note that there are a \r\n% couple of employees (ANDERSON and JACKSON) that do not get called in.\r\n% Also, while it seems that JOHNSON has two shifts, this is because we\r\n% are wrapping around and the shift actually goes from 10pm to 6am the next\r\n% morning.\r\n%\r\n% The second plot shows the number of staff on duty as compared to our\r\n% minimum staffing requirements.\r\n% We see that we meet the requirements exactly, without overstaffing or\r\n% understaffing at any point. This makes sense, because if the goal is to \r\n% minimize paid wages, then you'd try to stay as close to the minimum staff \r\n% as possible.\r\n\r\n%% 5. Conclusion\r\n% \r\n% The problem outlined in the example is not a trivial one to solve. To be \r\n% sure, it takes a bit of study and experience to be able to know how to\r\n% easily convert a real world problem into its equivalent MILP formulation.\r\n% \r\n% But for problems like this, especially at larger problem sizes, \r\n% attempting to solve it using other methods such trial-and-error or \r\n% genetic algorithms will usually take longer, and generally not yield as \r\n% good results. \r\n%\r\n% Have you tried using MATLAB's mixed-integer linear programming solver for\r\n% any of your own work? If you want to share your thoughts, or if you have \r\n% any questions regarding this example, please feel free to let us know in\r\n% <https:\/\/blogs.mathworks.com\/loren\/?p=1285#respond here>.\r\n%\r\n\r\n##### SOURCE END ##### 929c965db0ed419baab512796124b3e0\r\n-->","protected":false},"excerpt":{"rendered":"<div class=\"overview-image\"><img decoding=\"async\"  class=\"img-responsive\" src=\"https:\/\/blogs.mathworks.com\/images\/loren\/2016\/shiftScheduleDemoBlog_02.png\" onError=\"this.style.display ='none';\" \/><\/div><!--introduction--><p>Today's guest blogger is Teja Muppirala, who is a member of our Consulting Services group. Teja has been with MathWorks for 6 years and is based in our Tokyo office.... <a class=\"read-more\" href=\"https:\/\/blogs.mathworks.com\/loren\/2016\/01\/06\/generating-an-optimal-employee-work-schedule-using-integer-linear-programming\/\">read more >><\/a><\/p>","protected":false},"author":39,"featured_media":0,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[39,60,1],"tags":[],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/posts\/1285"}],"collection":[{"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/users\/39"}],"replies":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/comments?post=1285"}],"version-history":[{"count":6,"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/posts\/1285\/revisions"}],"predecessor-version":[{"id":1302,"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/posts\/1285\/revisions\/1302"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/media?parent=1285"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/categories?post=1285"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/tags?post=1285"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}