joh.tex

% TEXINPUTS=.:/home/hei2/Documents/LaTeX/myTexStyles:
% created for Journal of Heuristics 
% Time-stamp: "2011-07-18 14:17 hei2"
% if all fonts computer modern use -G1
% dvips -Ppdf -G0 <filename>
% http://www.springer.de/comp/lncs/authors.html

%!TEX TS-program = xelatex
%!TEX encoding = UTF-8 Unicode

%! xelatex joh.tex

\RequirePackage{fix-cm} 

\documentclass[twocolumn]{svjour3} 

\makeatletter
\def\cl@chapter{\@elt {theorem}}
\makeatother

\usepackage{mathptmx} % use Times fonts if available on your TeX system

\journalname{Journal of Scheduling}

\title{Discovering dispatching rules from data using imitation learning}
\subtitle{Case study for the job-shop problem}

\author{Helga Ingimundardottir \and Thomas Philip Runarsson}

\institute{H. Ingimundardottir \at
Dunhaga 5, IS-107 Reykjavik, Iceland \\
Tel.: +354-525-4704\\
Fax: +354-525-4632\\
\email{hei2@hi.is}\\
\and
T.P. Runarsson \at
Hjardarhagi 2-6, IS-107 Reykjavik, Iceland \\
Tel.: +354-525-4733\\
Fax: +354-525-4632\\
\email{tpr@hi.is}\\
}
\date{Received: November 6, 2015 / Accepted: date}
% The correct dates will be entered by the editor

\input{shorthand} % put your own shorthand declarations in this document

\begin{document}
\sloppy % fixes many bad boxes
\maketitle

\selectlanguage{english}

\begin{abstract}
    Dispatching rules can be automatically generated from scheduling data.
  This paper will demonstrate that the key to learning an effective
  dispatching rule is through the careful construction of the training data,  
  $\{\vec{x}_i(k),y_i(k)\}_{k=1}^K\in\mathcal{D}$, where
  \begin{enumerate*}
    \item features of partially constructed schedules $\vec{x}_i$ should 
    necessarily reflect the induced 
    data distribution $\mathcal{D}$ for when the rule is applied. This is 
    achieved by updating the learned model in 
    an active imitation learning fashion
    \item $y_i$ is labelled optimally using a MIP solver
    \item data needs to be balanced, as the set is unbalanced with respect to 
    the dispatching step $k$
  \end{enumerate*}
  Using the guidelines set by our framework the design of custom \dr s, for  
  a particular scheduling application, will become more effective. In the 
  study presented three different distributions of the \jsp\ will be 
  considered. The machine learning approach considered is based on preference 
  learning, i.e., which dispatch (post-decision state) is preferable to another.
 
  \keywords{Scheduling \and Composite \dr s \and Performance Analysis \and 
   Imitation Learning \and DAgger \and Preference Learning}
\end{abstract}

\section{Introduction}\label{sec:introduction}

Hand crafting heuristics for scheduling is an ad hoc approach to finding 
approximate solutions to problems. The practice is time-consuming and its 
performance can even vary dramatically between different problem instances. The 
aim of this work is to increase our understanding of this process. In 
particular, the learning of new problem specific priority \dr s (DR) will be 
addressed for a subclass of scheduling problems known as the \jsp\ scheduling 
problem (\JSP). 

A recent editorial the state-of-the-art approaches \cite{Chen13} in advanced 
\dr s for large-scale manufacturing systems reminds us that:
\lq\lq ... most traditional \dr s are based on historical data. 
With the emergence of data mining and on-line analytic processing, dispatching 
rules can now take predictive information into account.\rq\rq~The importance of 
automated discovery of \dr s was also emphasised by \cite{Monch13}. 
Data for learning can also be generated using a known heuristic on a set of 
problem instances.
Such an approach is taken in \cite{Siggi05} for a single machine where
a decision tree is learned from the data to have similar logic to the guiding 
\dr. 
However, the learned method cannot outperform the original \dr\ used for the 
data generation. 
This drawback is confronted in \cite{Malik08,Russell09,Siggi10} by using an 
optimal scheduler or policy, computed off-line. The 
resulting \dr s, as decision trees, gave significantly better schedules than 
using popular heuristics in that field, and a lower worst-case factor from 
optimality. 
Although using optimal policies for creating training data gives vital 
information on how to learn good scheduling rules an experimental study will 
show that this is not sufficient. 
Once these rules make a suboptimal dispatch then they are in uncharted 
territory and the effects are relatively unknown. 
This work will illustrate the sensitivity of learned \dr's performance in the 
way the training data is sampled.
For this purpose, \JSP\ is used as a case study to illustrate a methodology for 
generating meaningful training data, which can be successfully 
learned using preference-based imitation learning (IL).

The competing alternative to learning dispatching rules from data would be to 
search the \dr\ space directly. The prevalent approach in this case would be 
using an evolutionary algorithm, such as genetic programming (GP). 
The predominant approach in hyper-heuristics is a framework for creating 
new heuristics from a set of predefined heuristics via genetic algorithm 
optimisation \cite{Burke10}. 
Adopting a two-stage hyper-heuristic approach to generate a set of 
machine-specific DRs for dynamic \jsp, \cite{Pickardt2013} used genetic 
programming (GP) to evolve CDRs from basic features, along with an evolutionary 
algorithm to assign a CDR to a specific machine. 
The problem space consists of \jsp s in semiconductor manufacturing, with 
additional shop constraints, as machines are grouped into similar work centres, 
which can have different set-up times, workloads, etc. 
In fact, the GP emphasised efficient dispatch at the work centres with 
set-up requirements and batching capabilities, which are rules that are 
non-trivial to determine manually.

With meta-heuristics existing DRs can be used, and for example 
portfolio-based algorithm selection \cite{Rice76,Gomes01,Xu07}, either based on 
a single instance or class of instances, to determine which DR to choose from. 
Implementing ant colony optimisation to select the best DR 
from a selection of nine DRs for \JSP, experiments from \cite{Korytkowski13} 
showed that the choice of DR does affect the results and for all performance 
measures considered. They showed that it was better to have all the DRs to 
choose from rather than just a single DR at a time.
A simpler and more straightforward way to automate selection of \cdr s (CDR), 
\cite{InRu14}, translated \dr s into measurable features which describe the 
partial schedule and optimise directly what their contribution should be via 
evolutionary search. 

Using case based reasoning for timetable scheduling, training data in 
\cite{Burke06} is guided by the two best heuristics in the literature.
They point out that in order for their framework to be successful, problem 
features need to be sufficiently explanatory and training data needs to be 
selected carefully so they can suggest the appropriate solution for a specific 
range of new cases. 
When learning new \dr s there are several important 
factors to consider. First and foremost the context in which the training data 
is constructed 
will influence the quality of the learned \dr\ \cite{Burke06}. 
Since the training data consists of a collection of features, the quality of 
training data is interchangeable with the predictability of features. 
The training data is necessarily also problem instance specific. 
In addition to addressing these aspects, the paper will show that during 
the scheduling process, the most critical moment to make the `right' dispatch 
will vary. Furthermore, depending on the distribution of problem 
instances, these critical moments can vary greatly. 
Moreover, a supervised learning algorithm will optimize classification 
accuracy, while it is the actual end-performance of the dispatching rule 
learned that will determine the success of the learning method. 

The outline of the paper is as follows, with \cref{sec:problemdef} giving the 
mathematical formalities of the scheduling problem and 
\cref{sec:DR} describing the main features for \jsp\ 
and illustrating how schedules are created with \dr s. 
\Cref{sec:learnOPT} sets up the framework for learning from optimal schedules,  
in particular, the probability of choosing optimal decisions and the effects of 
making a suboptimal decision. Furthermore, the optimality of common \sdr s is 
investigated.
With these guidelines presented, \cref{ch:expr:CDR} goes into detail on how to 
create 
meaningful \cdr s using preference learning, focusing on how to 
compare operations and collect training data with the importance of  the
sampling strategy applied. 
\Cref{sec:il:active,sec:il:passive} explain the trajectories for 
sampling meaningful schedules used in preference learning, either 
using passive or active imitation learning. 
Experimental results are jointly presented in \cref{sec:il:expr} with 
comparison for a randomly generated problem space. 
Furthermore, some general adjustments for performance boost are also considered.
Following a discussion in \cref{sec:dis} the main conclusion and future work are presented in \cref{sec:con}.

\section{\Jsp~Scheduling}\label{sec:problemdef}
\JSP\ involves the scheduling of jobs for a set of 
machines. Each job consists of a number of operations which are then processed 
on the machines in a predetermined order. An optimal solution to the problem 
will depend on the specific objective. 

This study will consider the $n\times m$ \JSP, where $n$ jobs, 
$\mathcal{J}=\{J_j\}_{j=1}^n$, are scheduled on a finite set, 
$\mathcal{M}=\{M_a\}_{a=1}^m$, of $m$ machines. The index $j$ refers to a job 
$J_j\in\mathcal{J}$ while the index $a$ refers to a machine 
$M_a\in\mathcal{M}$. 
Each job requires a number of processing steps or operations, and the pair 
$(j,a)$ refers to the operation, i.e., processing the task of job $J_j$ on 
machine $M_a$. 

Each job $J_j$ has an indivisible operation time (or cost) on machine $M_a$, 
$p_{ja}$, which is assumed to be integral and finite. 
Starting time of job $J_j$ on machine $M_a$ is denoted $x_s(j,a)$ and its 
end time is denoted $x_e(j,a)$. 
Each job $J_j$ has a specified processing order through the machines. It is a 
permutation vector, $\vsigma_j$, of $\{1,\ldots,m\}$, representing a job $J_j$ 
which can be processed on $M_{\vsigma_j(a)}$ only after it has been completely 
processed on $M_{\vsigma_j(a-1)}$, namely:
\begin{equation}\quad\label{eq:permutation}
x_s(j,\vsigma_j(a)) \geq x_e(j,\vsigma_j(a-1)) 
\end{equation}
for all $J_j\in\mathcal{J}$ and $a\in\{2,..,m\}$. 
Note, that each job can have its own distinctive flow pattern through the 
machines which is independent of the other jobs. 
However, in the case that all jobs share the same fixed permutation 
route, it is referred to as \fsp~(\FSP). 

The objective function is to minimise the schedule's maximum completion times 
for all tasks, commonly referred to as the makespan, $C_{\max}$. 
This family of scheduling problems is denoted by $J||C_{\max}$ 
\cite{Pinedo08}.
Additional constraints commonly considered are job release-dates and due-dates 
or sequence dependent set-up times; however, these will not be considered here. 

In order to find an optimal (or near optimal) solution for scheduling problems 
it is possible to use either exact methods or heuristics methods. Exact methods 
guarantee an optimal solution. However, \jsp\ scheduling is strongly NP-hard 
\cite{Garey76:NPhard}. Any exact algorithm generally suffers from the curse of 
dimensionality, which impedes the application in finding the global optimum in 
a reasonable amount of time. 
Using state-of-the-art software for solving scheduling problems, such as 
LiSA 
(A Library of Scheduling Algorithms) \cite{LiSA}, which includes a specialised 
version of branch and bound that manages to find optimums for \jsp\ problems of 
up to $14\times14$ \cite{Ru12}. However, problems that are of greater size 
become intractable. 
Heuristics are generally more time efficient but 
do not necessarily attain the global optimum. Therefore, \jsp\ has the 
reputation of being notoriously difficult to solve. 
As a result, it's been widely studied in deterministic scheduling theory and 
its class of problems has been tested on a plethora of different solution 
methodologies from various research fields \cite{Meeran12}, all from simple and 
straight forward \dr s to highly sophisticated frameworks.

\section{Priority Dispatching Rules} \label{sec:DR}

Priority \dr s determine, from a list of incomplete jobs, 
$\mathcal{L}$, which job should be dispatched next. This process, where an 
example of 
a temporal partial schedule of six jobs scheduled on five machines, is 
illustrated in \cref{fig:jssp:example}.
The numbers in the boxes represent the job identification $j$. 
The width of the box illustrates the processing times for a given job for a 
particular machine $M_a$ (on the vertical axis). 
The dashed boxes represent the resulting partial schedule for when a particular 
job is scheduled next. 
Moreover, the current $C_{\max}$ is denoted by a dotted vertical line. 
The object is to keep this value as small as possible once all operations are 
complete. As shown in the example there are $15$ operations already scheduled. 
The \textit{sequence}, $\vchi$, of dispatches used to create this partial 
schedule is:
\begin{equation}\quad
\vchi=\left(J_3,J_3,J_3,J_3,J_4,J_4,J_5,J_1,J_1,J_2,J_4,J_6,J_4,J_5,J_3\right)
\end{equation}
This refers to the sequential ordering of job dispatches to machines, i.e., 
$(j,a)$; 
the collective set of allocated jobs to machines is interpreted by its 
sequence, which is referred to as a schedule.
A scheduling policy will pertain to the manner in which 
the sequence is determined from the available jobs to be scheduled. 
In our example, the available jobs are given by the job-list
$\mathcal{L}^{(k)}=\{J_1,J_2,J_4,J_5,J_6\}$ with the five potential jobs 
to be dispatched at step $k=16$ (note that $J_3$ is completed).

However, deciding which job to dispatch is not sufficient as one must also know 
where to place it. In order to build tight schedules, it is sensible to place a 
job as soon as it becomes available and such that the machine idle time is 
minimal, i.e., schedules are non-delay. 
There may also be a number of different options for such a placement. 
In \cref{fig:jssp:example} one observes that $J_2$, to be scheduled on $M_3$, 
could be placed immediately in a slot between $J_3$ and $J_4$, or after $J_4$ 
on this machine. 
If $J_6$ had been placed earlier, a slot would have been created between it and 
$J_4$, thus creating a third alternative, namely scheduling $J_2$ after $J_6$. 
The time in which machine $M_a$ is idle between consecutive jobs $J_j$ and 
$J_{j'}$ is called idle time or slack:
\begin{equation}\quad 
s(a,j):=x_s(j,a)-x_e(j',a) \label{eq:slack}
\end{equation}
where $J_j$ is the immediate successor of $J_{j'}$ on $M_a$. 

Construction heuristics are designed in such a way that it limits the search 
space in a logical manner while not excluding the optimum. Here, the 
construction heuristic, $\Upsilon$, is to schedule the dispatches as closely 
together as possible, i.e., minimise the schedule's idle time. 
More specifically, once an operation $(j,a)$ has been chosen from the job-list 
$\mathcal{L}$ by some dispatching rule, it can then be placed immediately after 
(but not prior) to $x_e(j,\vsigma_j(a-1))$ on machine $M_a$ due to constraint 
\cref{eq:permutation}, this could be considered the release time for $J_j$. 
However, to guarantee the disjunctive condition (i.e. no overlapping of jobs on 
machines) is not violated, idle times $M_a$ from \cref{eq:slack} are inspected 
as they create a slot which $J_j$ can occupy, but only considering those that 
occur after $J_j$ has been released from its previous machine, namely:
\begin{equation}\quad
\tilde{s}(a,j'):= x_s(j'',a)-\max\{x_e(j',a),x_e(j,\vsigma_j(a-1))\} 
\end{equation}
for all already dispatched jobs, $J_{j'},J_{j''}\in \mathcal{J}_a$ where 
$J_{j''}$ is $J_{j'}$ successor on $M_a$. Moreover, since preemption is not 
allowed, the only applicable slots for a new dispatch, $\tilde{S}_{ja}$, are 
those idle times that can process the entire operation, namely:
\begin{equation}\quad
\tilde{S}_{ja} := 
\left\{J_{j'}\in \mathcal{J}_a \;:\; \tilde{s}(a,j')\geq p_{ja} \right\}
\label{eq:slots}. \end{equation} 
The placement rule applied to decide which slot from \cref{eq:slots} is used to 
place $J_j$ is intrinsic to the construction heuristic, which is chosen 
independently of the priority dispatching rule that is applied. 
Different placement rules could be considered for selecting a slot, e.g., if 
the main concern were to utilise the slot space, then choosing the slot with 
the smallest idle time would yield a closer-fitted schedule and leave greater 
idle times undiminished for subsequent dispatches on $M_a$.
In our experiments, cases were discovered where such a placement could rule out 
the possibility of constructing optimal solutions.
However, this problem did not occur when jobs are simply placed as early as 
possible, which is beneficial for subsequent dispatches for $J_j$. 
For this reason, it will be the placement rule applied here.

\begin{figure}[t!]\centering
\includegraphics[width=\columnwidth]{figures/jssp_example}
\caption[Gantt chart of a partial \JSP\ schedule]{Gantt chart of a
partial \JSP\ schedule after 15 dispatches: Solid and dashed boxes
represent $\vchi$ and $\mathcal{L}^{(16)}$, respectively. Current
$C_{\max}$ denoted as dotted line.}
\label{fig:jssp:example}
\end{figure}

Priority \dr s will use features of operations, such as processing time, 
in order to determine the job with the highest priority. 
Consider again \cref{fig:jssp:example}; if the job with the shortest processing 
time (SPT) were to be scheduled next, then $J_2$ would be dispatched. 
Similarly, for the longest processing time (LPT) heuristic, $J_5$ would have 
the highest priority. 
Dispatching can also be based on features related to the partial schedule. 
Examples of these are dispatching the job with the most work remaining (MWR) or 
alternatively the least work remaining (LWR). A survey of more than $100$ of 
such rules are presented in \cite{Panwalkar77}. 
However, the reader is referred to an in-depth survey for simple or 
\emph{\sdr} (SDR) by \cite{Haupt89}. 
The SDRs assign an index to each job in the job-list and is generally only 
based on a few features and simple mathematical operations.

\begin{table}[t!] \centering
\caption[Feature space $\mathcal{F}$ for \JSP]{Feature space 
$\mathcal{F}$ for \JSP\ where job $J_j$ on machine $M_a$ given the 
resulting temporal schedule after operation $(j,a)$.
}
\label{tbl:jssp:feat}
{\setlength{\tabcolsep}{3pt} \input{tables/features-description}}
\end{table}

Designing priority \dr s requires recognising the important features of the 
partial schedules needed to create a reasonable scheduling rule. 
These features attempt to grasp key attributes of the schedule being 
constructed. Which features are most important will necessarily depend on the
objectives of the scheduling problem. 
The features used in this study applied for each possible operation encountered 
are given in \cref{tbl:jssp:feat}, where the set of machines already dispatched 
for $J_j$ is $\mathcal{M}_j\subset\mathcal{M}$, and similarly, $M_a$ has 
already had the jobs $\mathcal{J}_a\subset\mathcal{J}$ previously dispatched.
The features of particular interest were obtained by inspecting the 
aforementioned SDRs. Features \phiJobRelated\ and \phiMacRelated\ are 
job-related and machine-related, respectively.
In fact, \cite{Pickardt2013} note that in the current literature, there is a 
lack of global perspective in the feature space, as omitting them won't 
address the possible negative impact an operation $(j,a)$ might have on other 
machines at a later time. It is for this reason features such as 
\phiSlackRelated\ are considered, since they are slack related and are a means 
of indicating the current quality of the schedule.

Priority \dr s are attractive since they are relatively easy to 
implement, perform fast, and find reasonable schedules. In addition, they are 
relatively easy to interpret, which makes them desirable for the end-user.
However, they can also fail unpredictably. 
A careful combination of \dr s has been shown to perform significantly better 
\cite{Jayamohan04}. These are referred to as \emph{\cdr s} 
(CDR), where the priority ranking is an expression of several \dr s. 
CDRs deal with a greater number of more complicated functions and are 
constructed from the schedule features. In short, a CDR is a combination 
of several DRs. 
For instance, let $\pi$ be a CDR comprised of $d$ DRs, then the index $I$ for 
$J_j\in\mathcal{L}^{(k)}$ using $\pi$ is:
\begin{equation}\quad I_j^{\pi} = \sum_{i=1}^d w_i \pi_i(\vchi^j) 
\label{eq:CDR}
\end{equation}
where $w_i>0$ and $\sum_{i=0}^d w_i = 1$ with $w_i$ giving the weight of the 
influence of $\pi_i$ (which could be an SDR or another CDR) to $\pi$. Note: 
each $\pi_i$ is a function of $J_j$'s features from the current sequence 
$\vchi$, where $\vchi^j$ implies that $J_j$ was the latest dispatch, i.e., the 
partial schedule given $\chi_k=J_j$.

At each step $k$, an operation is dispatched which has the highest 
priority. If there is a tie, some other priority measure is used. Generally 
the \dr s are static during the entire scheduling process. However, ties could 
also be broken randomly (RND). 

While investigating 11 SDRs for JSP, \cite{Lu13} a pool of 33 CDRs was created. 
This pool strongly outperformed the original CDRs by using multi-contextual 
functions based on either job waiting time or machine idle time (similar 
to \phiwait\ and \phimacSlack\ in \cref{tbl:jssp:feat}), i.e., the CDRs are a 
combination of either one or both of these key features and then the SDRs.
However, no combinations of the basic SDRs were explored, only those two 
features. 
Similarly, using priority rules to combine 12 existing DRs from the literature, 
\cite{Yu13} had 48 CDR combinations which yielded 48 different models 
to implement and test. 
It is intuitive to get a boost in performance by introducing new CDRs, since 
where one DR might be failing, another could be excelling, so combining them 
together should yield a better CDR. However, these approaches introduce fairly 
ad hoc solutions and there is no guarantee the optimal combination of 
\dr s have been found.

The \cdr\ presented in \cref{eq:CDR} can be considered as a special case of 
the following general linear value function:
\begin{equation}\quad\label{eq:jssp:linweights}
\pi(\vchi^j)=\sum_{i=1}^d w_i \phi_i(\vchi^j).
\end{equation}
when $\pi_i=\phi_i$, i.e., a composite function of the features 
from \cref{tbl:jssp:feat}. Finally, the job to be dispatched, $J_{j^*}$, 
corresponds to the one with the highest value, namely:
\begin{equation}\quad\label{eq:jstar}
J_{j^*}=\argmax_{J_j\in \mathcal{L}}\; \pi(\vchi^j)
\end{equation}
Similarly, \sdr s may be described by this linear model. For instance, let all 
$w_i=0$, but with following exceptions: $w_1=-1$ for SPT, $w_1=+1$ for LPT, 
$w_7=-1$ for LWR and $w_7=+1$ for MWR. Generally, the weights $\vec{w}$ are 
chosen by the designer or the 
rule apriori. A more attractive approach would be to learn these weights from 
problem examples directly. The following \lcnamecref{sec:learnOPT} will  
investigate how this may be accomplished.

\section{Performance Analysis of Priority Dispatching Rules}\label{sec:learnOPT}

In order to create successful \dr s, a good starting point is to 
investigate the properties of optimal solutions and hopefully be able to learn 
how to mimic the construction of such solutions. For this, optimal 
solutions (obtained by using a commercial software package \cite{gurobi}) are 
followed and the probability of SDRs being optimal is inspected. 
This serves as an indicator of how hard it is to put our objective up as a 
machine learning problem. 
However, the end goal, which is minimising \namerho, must also be taken into 
consideration because its relationship to stepwise optimality is not fully 
understood.

In this \lcnamecref{sec:learnOPT} the concerns of learning new priority \dr s 
will be addressed. At the same the time experimental set-up used in the study 
is described. 

\subsection{Problem Instances}\label{sec:data:sim}

The class of problem instances used in our studies was the \jsp\ scheduling 
problem described in \cref{sec:problemdef}. Each instance will have 
different processing times and machine ordering. Each instance will 
therefore create different challenges for a priority dispatching rule. 
Dispatching rules learned will be customised for the problems used for their 
training. For real world application using historical data would be most 
appropriate. The aim would be to learn a dispatching rule that works well on 
average for a given distribution of problem instances. 
To illustrate the performance difference of priority \dr s on different 
problem distributions within the same class of problems, 
consider the following three cases.
Problem instances for JSP are generated stochastically by fixing the number of 
jobs and machines to ten. A discrete processing time is sampled independently 
from a discrete uniform distribution from the interval $I=[u_1,u_2]$, i.e., 
$\vec{p}\sim \mathcal{U}(u_1,u_2)$. 
The machine order was a random permutation of all of the machines in the 
\jsp. Two different processing times distributions were explored, namely 
\jrnd{n}{m} where $I=[1,99]$ and \jrndn{n}{m} where $I=[45,55]$. These 
instances are referred to as random and random-narrow, respectively. 
In addition, the case where the machine order is fixed and the same for 
all jobs, i.e. $\sigma_j(a)=a$ for all $J_j\in\mathcal{J}$ and where 
$\vec{p}\sim\mathcal{U}(1,99)$, was also considered. 
These jobs are denoted by \frnd{n}{m} and are analogous to \jrnd{n}{m}.
The problem spaces are summarised in \cref{tbl:data:sim}.

The goal is to minimise the makespan, $C_{\max}$. The optimum 
makespan is denoted $C_{\max}^{\pi_\star}$ (using the expert policy 
$\pi_\star$), and the makespan obtained from the 
scheduling policy $\pi$ under inspection by $C_{\max}^{\pi}$. Since the optimal 
makespan varies between problem instances the performance measure is the 
following:
\begin{equation}\quad\label{eq:rho}
\rho=\frac{C_{\max}^{\pi}-C_{\max}^{\pi_\star}}{C_{\max}^{\pi_\star}}\cdot
100\%
\end{equation}
which indicates the percentage relative deviation from optimality. 
Note: \cref{eq:rho} measures the discrepancy between predicted value and true 
outcome, and is commonly referred to as a loss function, which should be 
minimised for policy $\pi$.

\Cref{fig:boxplot:SDR} depicts the box plot for \cref{eq:rho} when using the 
SDRs from \cref{sec:DR} for all of the problem spaces from \cref{tbl:data:sim}.
These box plots show the difference in performance of the various SDRs. The 
rule MWR performs on average the best on the \jrnd{n}{m} and \jrndn{n}{m} 
problem instances, whereas for \frnd{n}{m} it is LWR that performs best. It is 
also interesting to observe that all but MWR perform statistically worse than a 
random job dispatching on the \jrnd{n}{m} and \jrndn{n}{m} problem instances.

\begin{table}\centering
\caption[Problem space distributions used in experimental studies.]{Problem 
space distributions used in experimental studies. Note, problem instances 
are synthetic and each problem space is i.i.d. 
}\label{tbl:data:sim}
\input{tables/data-sim}
\end{table} 

\begin{figure}
\centering
\includegraphics[width=\columnwidth]{figures/{boxplotRho_SDR_10x10}.pdf}
\caption{Box plot for deviation from optimality, $\rho$, (\%) for SDRs}
\label{fig:boxplot:SDR}
\end{figure}

\subsection{Reconstructing optimal solutions}\label{sec:opt:sdr}

When building a complete schedule, $K=n\cdot m$ dispatches must be made 
sequentially. A job is placed at the earliest available time slot for its next 
machine, whilst still responding to the fact that each machine can handle at 
most one job at each time, and jobs need to have been finished by their 
previous machines according 
to their machine order. Unfinished jobs from the job list $\mathcal{L}$ are 
dispatched one at a time according to a deterministic 
scheduling policy (or heuristic). This process is given as a pseudo-code in 
\cref{pseudo:constructJSP}. 
After each dispatch\footnote{Dispatch and time step are used interchangeably.} 
the schedule's current features are updated based on 
the half-finished schedule, $\vchi$. 
For each possible post-decision state, $J_j\in\mathcal{L}^{(k)}$, the temporal 
features, $\vphi_j$, are collected (cf. 
\cref{pseudo:constructJSP:phi}) forming the feature set, $\Phi$, based on the 
union of all 
$N_{\text{train}}$ problem instances available, namely:
\begin{equation}\quad \label{eq:Phi}
\Phi := \bigcup_{\{\vec{x}_i\}_{i=1}^{N_{\text{train}}}} 
\left\{\vphi^j \;:\; J_j\in\mathcal{L}^{(k)}\right\}_{k=1}^K
\subset\mathcal{F}
\end{equation}
where the feature space $\mathcal{F}$ is described in \cref{tbl:jssp:feat}, and 
are based on job and machine features which are widespread in practice.

\input{pseudocode/constructJSP}

It is easy to see that the sequence of task assignments is by no means unique. 
Inspecting a partial schedule further along in the dispatching process such as 
in \cref{fig:jssp:example}, then let's say $J_1$ would be dispatched next, and 
in the next iteration $J_2$. Now this sequence would yield the same schedule as 
if $J_2$ had been dispatched first and then $J_1$ in the next iteration, 
i.e., these are jobs with non-conflicting machines. 
In this particular scenario, one cannot infer that choosing $J_1$ is better 
and $J_2$ is worse (or vice versa) since they can both yield the same solution.
Furthermore, 
there may be multiple optimal solutions to the same 
problem instance. Hence not only is the sequence representation `flawed' in the 
sense that slight permutations on the sequence are in fact equivalent in terms 
of the end-result, but very varying permutations on the dispatching sequence 
(although given the same partial initial sequence) can result in very different 
complete schedules yet can still achieve the same makespan. 

The redundancy in building optimal solutions using dispatching rules means that 
many different dispatches may yield an optimal solution to the problem instance.
Let's formalise the probability of optimality (or stepwise 
classification accuracy) for a given policy $\pi$, yet following the expert 
policy $\pi_\star$, as:
\begin{equation}\quad \label{eq:tracc:opt}
\xi^\star_{\pi} := \mathbb{E}_{\pi_\star}\left\{\pi_{\star} = \pi \right\}
\end{equation}
that is to say, the mean likelihood of our policy $\pi$ being equivalent to the 
expert policy $\pi_\star$.
The probability that a job chosen by a SDR yields an optimal makespan on a 
step-by-step basis, i.e., $\xi^\star_{\langle \text{SDR} \rangle}$, is depicted 
in \cref{fig:opt:SDR:xistar}. These probabilities
vary quite a bit between the different problem instance distributions studied. 
From \cref{fig:opt:SDR:xistar} it is observed that $\xi^\star_{\text{MWR}}$ has 
a higher probability than random guessing, in choosing a dispatch which may 
result in an optimal schedule. This is especially true towards the end of the 
schedule building process. 
Similarly, $\xi^\star_{\text{LWR}}$ selects dispatches 
resulting in optimal schedules with a higher probability. This would appear to 
support the idea that the higher the probability of dispatching job that 
may lead to an optimal schedule, the better the SDRs performance, as 
illustrated by \cref{fig:boxplot:SDR}. However, there is a counter example; 
$\xi^\star_{\text{SPT}}$ has a higher probability than random dispatching of 
selecting a jobs that may lead to an optimal solution. Nevertheless, the random 
dispatching performs better than SPT on problem instances \jrnd{10}{10} and 
\jrndn{10}{10}. 

\begin{figure}[t]
\centering
\includegraphics[width=\columnwidth]{figures/{trdat_prob_moveIsOptimal_10x10_SDR_xistar}.pdf}
\caption{Probability of SDR being optimal, 
$\xi^\star_{\langle\text{SDR}\rangle}$}
\label{fig:opt:SDR:xistar}
\end{figure}

As shown in \cref{fig:opt:SDR:xistar}, \jrnd{10}{10} has a relatively high 
probability ($70\%$ and above) of selecting an optimal job at random. 
However, it is imperative to keep making optimal decisions, because the 
consequences of making suboptimal dispatches are unknown. 
To demonstrate this \cref{fig:case} depicts mean worst and best case scenario 
of the resulting \namerho, once off the optimal track, defined as follows:
\begin{subequations}\label{eq:bwc:opt}
\begin{eqnarray} 
\quad \zeta_{\min}^{\star}(k) &:=& \mathbb{E}_{\pi_\star}\left\{
\min_{J_j\in\mathcal{L}^{(k)}}(\rho) \;:\;
\forall C_{\max}^{\vchi^j} \gneq C_{\max}^{\pi_\star} \right\} \\
\quad \zeta_{\max}^{\star}(k) &:=& \mathbb{E}_{\pi_\star}\left\{
\max_{J_j\in\mathcal{L}^{(k)}}(\rho) \;:\;
\forall C_{\max}^{\vchi^j} \gneq C_{\max}^{\pi_\star} \right\}
\end{eqnarray}
\end{subequations}
Note, that it is a given that there is only one non-optimal dispatch. 
Generally, there will be more, and then the compound effects of making 
suboptimal decisions are cumulative. 

It is interesting to observe that for \jrnd{10}{10} and \jrndn{10}{10} making 
suboptimal decisions later impacts on the resulting makespan more than making a 
mistake early. 
The opposite seems to be the case for \frnd{10}{10}. 
In this case it is imperative to make good decisions right from the start. This 
is due to the major structural differences between \JSP\ and \FSP, namely the 
latter having a homogeneous machine ordering and therefore constricting the 
solution immensely. 

\begin{figure}[t]
\centering
\includegraphics[width=\columnwidth]{figures/{stepwise_10x10_OPT_casescenario}.pdf}
\caption{Mean deviation from optimality, $\rho$, (\%), for best and worst 
case scenarios of making one suboptimal dispatch (i.e. 
$\zeta^{\star}_{\min}$ and $\zeta^{\star}_{\max}$), depicted as lower and 
upper bounds, respectively, for \jrnd{10}{10}, \jrndn{10}{10} and 
\frnd{10}{10}. The mean suboptimal move is given as a dashed line.}
\label{fig:case}
\end{figure}

\subsection{Blended \dr s}\label{sec:opt:bdr}
A naive approach to creating a simple blended \dr~(BDR) would be to 
switch between SDRs at a predetermined time.
Observing again \cref{fig:opt:SDR:xistar}, a presumably good BDR for 
\jrnd{10}{10} would be to start with $\xi^\star_{\text{SPT}}$ and then switch 
over to $\xi^\star_{\text{MWR}}$ at around time step $k=40$, where the SDRs 
change places in outperforming one another. 
A box plot for $\rho$ for the BDR compared with MWR and SPT is depicted in 
\cref{fig:boxplot:BDR} and its main statistics are reported in 
\cref{tbl:BDR:stats}. 
This simple swap between SDRs does outperform the SPT heuristic, yet 
doesn't manage to gain the performance edge of MWR. Using SPT downgrades the 
performance of MWR.
A reason for this lack of performance of our proposed BDR is perhaps that by 
starting out with SPT in the beginning, it sets up the schedules in such a way 
that it's quite greedy and only takes into consideration jobs with the shortest 
immediate processing times. Now, even though it is possible to find optimal 
schedules from this scenario, as \cref{fig:opt:SDR:xistar} shows, the inherent 
structure that's already taking place might make it hard to come across by 
simple methods. Therefore, it is by no means guaranteed that by simply swapping 
over to MWR will handle the situation that applying SPT has already created. 
\Cref{fig:boxplot:BDR} does show, however, that by applying MWR instead of SPT 
in the latter stages does help the schedule to be more compact in terms of SPT. 
However, the fact remains that the schedules have diverged too far from what 
MWR would have been able to achieve on its own. 

\begin{figure}[t]
\centering
\includegraphics[width=\columnwidth]{figures/j_rnd/{boxplotRho_BDR_10x10}.pdf}
\caption{Box plot for deviation from optimality, $\rho$, (\%) for BDR where 
SPT is applied for the first 10\%, 15\%, 20\%, 30\% or 40\% of the 
dispatches, followed by MWR}
\label{fig:boxplot:BDR}
\end{figure}
\input{tables/stats.BDR.10x10.tex}

In \cref{fig:opt:SDR:xistar} the stepwise optimality was inspected, given that 
all committed dispatches were based on the optimal trajectory. 
As mistakes are bound to be made at some points, it is interesting to see how 
the stepwise optimality evolves for its intended trajectory, $\pi$, thereby 
updating \cref{eq:tracc:opt} to:
\begin{equation}\quad \label{eq:tracc:track}
\xi_{\pi} := \mathbb{E}_{\pi}\left\{\pi_{\star} = \pi \right\}
\end{equation}
\Cref{fig:opt:SDR:xi} shows the log likelihood for $\xi_{\langle 
\text{SDR} \rangle}$ using \jrnd{10}{10}.
There one can see that even though $\xi_{\text{SPT}}$ is generally more 
likely to find optimal dispatches in the initial steps, then shortly after 
$k=15$, $\xi_{\text{MWR}}$ becomes a 
contender again. This could explain why our BDR switch at $k=40$ from 
\cref{fig:boxplot:BDR} was unsuccessful. However, changing to MWR at $k\leq20$
is not statistically significant from MWR (the boost in mean $\rho$ is at most 
-0.5\%). 
But as pointed out for \cref{fig:case}, it's not so fatal to 
make bad moves in the very first dispatches for \jrnd{10}{10}, hence there is 
little gain with improved classification accuracy in that region. 
However, after $k>20$ then the BDR performance starts diverging from that of 
MWR. 

\begin{figure}[t]
\centering
\includegraphics[width=\columnwidth]{figures/j_rnd/{trdat_prob_moveIsOptimal_10x10_SDR_xi}.pdf}
\caption{Log likelihood of SDR being optimal for \jrnd{10}{10}\ when 
following its corresponding SDR trajectory, i.e.,  
$\log\left(\xi_{\langle\text{SDR}\rangle}\right)$}
\label{fig:opt:SDR:xi}
\end{figure}

\section{Preference Learning}\label{ch:expr:CDR}
\Cref{sec:opt:bdr} demonstrated there is something to be gained by
trying out different combinations of DRs; however, it is non-trivial. This 
section presents one approach to learning how such combinations can work.
Learning models considered in this study are based on ordinal regression in 
which the learning task is formulated as learning preferences, and in the case 
of scheduling, learning which operations are preferred to others. Ordinal 
regression has been previously presented in \cite{Ru06:PPSN} and in 
\cite{InRu11a} for \JSP, and given here for completeness. 

The optimum makespan is known for each problem instance. At each time step $k$, 
a number of feature pairs are created. 
Let $\vphi^{o}\in\mathcal{F}$ denote the post-decision state when dispatching 
$J_o\in\mathcal{O}^{(k)}$ corresponds to an optimal schedule being built. 
All post-decision states corresponding to suboptimal dispatches, 
$J_s\in\mathcal{S}^{(k)}$, are denoted by $\vphi^{s}\in\mathcal{F}$.

The approach taken here is to verify analytically, at each time step, by fixing 
the current temporal schedule as an initial state, whether it is possible to 
somehow yield an optimal schedule by manipulating the remainder of the 
sequence. This also takes care of the scenario that having dispatched a job 
resulting in a different temporal makespan would have resulted in the same 
final makespan if another optimal dispatching sequence had been chosen. 
That is to say, the training data generation takes into consideration when 
there are multiple optimal solutions\footnote{
There can be several optimal solutions available for each problem instance. 
However, it is deemed sufficient to inspect only one optimal trajectory per 
problem instance as there are $N_{\text{train}}=300$ independent instances, 
which gives the training data variety.} 
to the same problem instance. 

Let's compare features from \cref{eq:Phi}, and define the difference of optimal 
from suboptimal ones as, 
\mbox{$\vpsi^{o}=\vphi^{o}-\vphi^{s}$}, and vice versa, 
\mbox{$\vpsi^{s}=\vphi^{s}-\vphi^{o}$}, and label the differences by $y_o=+1$ 
and $y_s=-1$ respectively. 
Then, the preference learning problem is specified by a set of preference pairs:
\begin{equation}
\begin{split}
\quad \Psi & = \left\{\left(\vpsi^o,y_o\right),\left(\vpsi^s,y_s\right)
\;:\; 
\forall \left(J_o,J_s\right) \in \mathcal{O}^{(k)} \times 
\mathcal{S}^{(k)}\right\}_{k=1}^{K} \\
& \subset  \Phi\times Y
\end{split} \label{eq:prefset}
\end{equation}
where $\Phi\subset \mathbb{R}^d$ is the training set of $d=\NrFeatLocal$ 
features (cf. \cref{tbl:jssp:feat}), $Y=\{+1,-1\}$ is the outcome space from 
job pairs $J_o\in\mathcal{O}^{(k)}$ and $J_s\in\mathcal{S}^{(k)}$, for all 
dispatch steps $k$.

To summarise, each job is compared against another job of the job list, 
$\mathcal{L}^{(k)}$, and if the makespan differs, i.e., 
$C_{\max}^{\pi_\star(\vchi^s)} \gneq C_{\max}^{\pi_\star(\vchi^o)}$ an 
optimal/suboptimal pair is created. 
However, if the makespans are identical the pair is omitted since they give the 
same optimal makespan. 
In this way, only features from a dispatch resulting in a suboptimal solution 
is labelled undesirable.

Now let's consider mappings from solutions to ranks, such 
function $\pi$ induces an ordering on the solutions by the following rule:
\begin{equation}\quad
\vchi^i \succ \vchi^j \quad \iff \quad \pi(\vchi^i) > 
\pi(\vchi^j)
\end{equation}
where the symbol $\succ$ denotes \emph{is preferred to} with respect to the 
solutions' ranks. Referring to the ranks in a discrete  manner as follows:
\begin{equation}\quad\label{eq:rankorder}
r_1 \succ r_2 \succ \cdots \succ r_{n'} \quad (n' \leq n)
\end{equation}
implies $r_1$ is preferable to $r_2$, and $r_2$ is preferable to $r_3$, etc. 

The function used to bring about preference is a linear function in the feature 
space, previously defined in \cref{eq:jssp:linweights}.
Logistic regression learns the optimal parameters $\vec{w}^*\in\mathbb{R}^d$. 
For this study, L2-regularised logistic regression from the \textsc{liblinear} 
package \cite{liblinear} without bias is used to learn the preference set 
$\Psi$, defined by \cref{eq:prefset}.
Hence the job chosen to be dispatched, $J_{j^*}$, is the one corresponding to 
the highest preference estimate, i.e., \cref{eq:jstar} where $\pi$ is 
the classification model obtained by the preference set.

Preliminary experiments for creating step-by-step model was done in 
\cite{InRu11a} resulting in a local linear model for each dispatch for a total 
of $K$ linear models for solving $n\times m$ JSP. 
However, the experiments showed that by fixing the weights to each mean 
value throughout the dispatching sequence the results remained satisfactory.
A more sophisticated way would be to create a new linear model, where 
the preference set, $\Psi$, is the union of the preference pairs across the 
$K$ dispatches, as described in \cref{eq:prefset}. 
This would amount to a substantial preference set, and for $\Psi$ to be 
computationally feasible to learn, $\Psi$ has to be reduced. For this several 
ranking strategies were explored in \cite{InRu15a}, and the results showed 
that it is sufficient to use partial subsequent rankings with combinations 
of $r_i$ and $r_{i+1}$ for $i\in\{1,\ldots,n'\}$, added to the preference 
set in such a manner that there is more than one 
operation with the same ranking and only one from that rank is needed to be 
compared to the subsequent rank. 
Moreover, for this study, which dealt with $10\times 10$ problem instances 
instead of $6\times5$, the partial subsequent ranking became necessary, as 
full ranking was computationally infeasible due to its size. 
Defining the size of the preference set as $l=\abs{\Psi}$, then if $l$ is too 
large, re-sampling to size $l_{\max}$ may be needed in order for the 
ordinal regression to be computationally feasible. 

The training data from \cite{InRu11a} was created from optimal solutions of 
randomly generated problem instances, i.e., traditional \emph{passive
imitation learning} (PIL). 
As \JSP\ is a sequential decision making process, errors are bound to emerge. 
Due to the compound effect of making suboptimal dispatches, the model leads the 
schedule astray from learned feature-space, resulting in the new input being 
foreign to the learned model. 
Alternatively, training data could be generated using suboptimal solution 
trajectories as well, as was done in \cite{InRu15a}, where the training data 
was also incorporated following the trajectories obtained by applying 
successful SDRs from the literature. 
The reasoning behind it was that they would be beneficial for learning, 
as they might help the model to escape from local minima once off the coveted 
optimal path. 
Simply aggregating training data obtained by following the trajectories of 
well-known SDRs yielded better models with lower \namerho. 

Inspired by the work of \cite{RossB10,RossGB11}, the methodology of generating 
training data will now be such that it will iteratively improve upon the model, 
such that the feature-space learned will be representative of the feature-space 
the eventual model would likely encounter, known as DAgger for \emph{active 
imitation learning} (AIL) and thereby, eliminating the ad hoc nature of 
choosing trajectories to learn, by instead letting the model lead its own way 
in a self-perpetuating manner until it converges.

Furthermore, in order to boost training accuracy, two strategies were explored 
\begin{enumerate}[after={{}}, leftmargin=*,
label={\textbf{Boost.\arabic*}}, ref={{Boost.\arabic*}}]
\item \label{expr:boost:varylmax} increasing number of preferences used 
in training (i.e. varying \mbox{$l_{\max} \leq \abs{\Psi}$}),
\item \label{expr:boost:newdata} introducing more problem instances (denoted 
EXT in experimental setting).
\end{enumerate}
Note, the following experimental studies will address 
\ref{expr:boost:newdata}, whereas preliminary experiments for 
\ref{expr:boost:varylmax} showed no statistical significance in boost of 
performance. Hence, the default set-up will be $l_{\max}=5 \cdot 10^5$ which is 
roughly the amount of features encountered from one pass of sampling a 
\mbox{$K$-stepped} trajectory using a fixed policy $\pi$ for the default 
$N_{\text{train}}=300$. 

Another way to adjust training accuracy is to give different weights to various 
time steps. To address this problem, two different stepwise sampling biases (or 
data balancing techniques) are here considered
\begin{enumerate}[after={{}}, leftmargin=*,
label={\textbf{Bias.\arabic*}}, ref={{Bias.\arabic*}}]
\item \label{bias:equal} \textbf{(equal)} where each time step has equal 
probability. This was used in \cite{InRu14,InRu15a} and serves as a baseline.
\item \label{bias:adjdbl2nd} \textbf{(adjdbl2nd)} where each time step is 
adjusted to the number of preference pairs for that particular step (i.e. 
each step now has equal probability irrespective of quantity of encountered 
features). This is done with re-sampling.
In addition, there is superimposed twice as much likelihood of 
choosing pairs from the latter half of the dispatching process. 
\end{enumerate}

\noindent Remark: as the following \lcnamecref{sec:il:passive}s require 
repeated collection of training data and since labelling is a very 
time intensive task, the remainder of the paper will solely focus 
on \jrnd{10}{10}.
The experiments in \cite{InRu15a} generally took several hours to collect for 
$N_{\text{train}}^{6\times5}=500$. 
Going to higher dimension and using $N_{\text{train}}^{10\times10}=300$ 
computational intensity becomes an issue, as \jrnd{10}{10} needs a few days and 
\jrndn{10}{10} or \frnd{10}{10} require several weeks using a personal computer.

\section{Passive Imitation Learning}\label{sec:il:passive}
Using the terms from game theory used in \cite{CesaBianchi06}, 
then our problem is a basic version of the sequential prediction problem where 
the predictor (or forecaster), $\pi$, observes each element of a sequence 
$\vchi$ of jobs, where at each time step $k \in \{1,...,K\}$, before the 
$k$-th job of the sequence is revealed, the predictor guesses its value 
$\chi_k$ on the basis of the previous $k-1$ observations. 

\subsection{Prediction with Expert Advice}\label{sec:expertPolicy}
Let us assume the expert policy $\pi^\star$ is known, which can query what 
is the optimal choice of $\chi_k={j^*}$ at any given time step $k$. 
Now let's use \cref{eq:jstar} to back-propagate the relationship between 
post-decision states and $\hat{\pi}$ with preference learning via our collected 
feature set, denoted $\Phi^\text{OPT}$, i.e., collecting the feature set 
corresponding to the following optimal tasks $J_{j^*}$ from $\pi^\star$ in 
\cref{pseudo:constructJSP}.
This baseline sampling trajectory originally introduced in \cite{InRu11a} for 
adding features to the feature set is a pure strategy where at each dispatch an 
optimal task is dispatched.

By querying the expert policy, $\pi_\star$, the ranking of the job-list, 
$\mathcal{L}$, satisfies \cref{eq:rankorder}. 
In this study, then it's known that the best rank $r_1 \propto 
C_{\max}^{\pi_\star}$, hence the optimal job-list is the following:
\begin{equation}\quad
\mathcal{O}=\left\{J_j \;:\; r_1 \propto \min_{J_j \in \mathcal{L}} 
C_{\max}^{\pi_\star(\vchi^j)}\right\}
\end{equation}
found by solving the current partial schedule to optimality using a MIP solver.

When $\abs{\mathcal{O}^{(k)}}>1$, there can be several trajectories worth 
exploring. However, only one is chosen at random. This is deemed sufficient as 
the number of problem instances, $N_{\text{train}}$, is relatively large.

\subsection{Follow the Perturbed Leader}\label{sec:perturbedLeader}
By allowing a predictor to randomise it's possible to achieve improved 
performance \cite{CesaBianchi06,Hannan57}. This is the inspiration for our next 
strategy called Follow the Perturbed Leader, denoted OPT$\epsilon$. 
Its pseudo code is given in \cref{pseudo:perturbedLeader} and describes how the 
expert policy (i.e. optimal trajectory) from \cref{sec:expertPolicy} is subtly
``perturbed'' with $\epsilon=10\%$ likelihood, by choosing a job corresponding 
to the second best $C_{\max}$ instead of an optimal one with some small 
probability. 

\input{pseudocode/perturbedLeader}

\subsection{Experimental study}\label{sec:pil:expr}

\begin{figure}
\includegraphics[width=\columnwidth]{figures/{j_rnd}/{boxplot_passive_10x10}.pdf}
\caption{Box plot for \jrnd{10}{10} \namerho, using either expert policy and 
following perturbed leader.}\label{fig:passive:boxplot} 
\end{figure}

Results for \cref{sec:expertPolicy,sec:perturbedLeader} using the \jrnd{10}{10} 
box plot of \namerho, is given in \cref{fig:passive:boxplot} and the main 
statistics are reported in \cref{tbl:IL:stats}. 
To address \ref{expr:boost:newdata}, the extended training set was simply 
obtained by iterating over more examples, namely $N^{\text{OPT}}_{\text{train, 
EXT}}=1000$. However, one can see that the increased number of varied features 
dissuades the preference models from achieving a good performance in terms of 
$\rho$. 
It's preferable to use the default $N^{\text{OPT}}_{\text{train}}=300$ 
allowing slight perturbations of the optimal trajectory, as was done for 
$\Phi^{\text{OPT}\epsilon}$. Unfortunately, all this overhead has not managed 
to surpass MWR in performance, except for $\Phi^{\text{OPT}\epsilon}$ using 
\ref{bias:adjdbl2nd} with a $\Delta\rho\approx-4.24\%$ boost in mean 
performance. Otherwise, for \ref{bias:equal}, there is a loss of 
$\Delta\rho\approx+6.23\%$ in mean performance. 
This is likely due to the fact that if equal probability is used for stepwise 
sampling, then there is hardly any emphasis given to the final dispatches as 
there are relatively few (compared to previous steps) preference pairs 
belonging to those final stages.
Revisiting \cref{fig:case}, then the band for 
$\{\zeta^{\star}_{\min},\zeta^{\star}_{\max}\}$ is quite tight as the problem 
is immensely constricted and few operations to choose from. However, the 
empirical evidence from using \ref{bias:adjdbl2nd} shows that it is 
imperative to make the right decisions at the very end.

Based on the results from \cite{InRu11a} the expert policy is a promising 
starting point.
However, that was for $6\times5$ dimensionality (i.e. $K=30$), which is a much 
simpler problem space. Notice that in \cref{fig:opt:SDR:xi} there was 
virtually no chance for $\xi_\pi(k)$ of choosing a job resulting in optimal 
makespan after step $k=28$.
Since \jsp\ is a sequential prediction problem, all future observations are 
dependent on previous operations. 
Therefore, learning sampled features that correspond only to optimal or 
near-optimal schedules isn't of much use when the preference model has 
diverged too far. \Cref{sec:opt:bdr} showed that good classification 
accuracy based on $\xi^\star_\pi$ does not necessarily mean a low mean \namerho.
This is due to the learner's predictions affects future input observations 
during its execution, which violates the crucial i.i.d. assumptions of the 
learning approach, and ignoring this interaction leads to poor performance.
In fact, \cite{RossB10} proves that assuming the preference model has a 
training error of $\epsilon$, then the total compound error (for all $K$ 
dispatches), the classifier induces itself and grows quadratically, 
$\bigOh{\epsilon K^2}$, for the entire schedule, rather than having linear 
loss, $\bigOh{\epsilon K}$, as if it were i.i.d.

\section{Active Imitation Learning}\label{sec:il:active}

To amend performance from $\Phi^{\text{OPT}}$-based models, suboptimal 
partial schedules were explored in \cite{InRu15a} by inspecting the features 
from successful SDRs, $\Phi^{\langle\text{SDR}\rangle}$, by passively observing 
a full execution of following the task chosen by the corresponding SDR. 
This required some trial-and-error as the experiments showed that features 
obtained by SDR trajectories were not equally useful for learning.

To automate this process, inspiration from AIL presented in \cite{RossGB11} was 
sought, called \emph{Dataset Aggregation} (DAgger) method, which addresses a 
no-regret algorithm in an on-line learning setting. 
The novel meta-algorithm for IL learns a deterministic policy guaranteed to 
perform well under its induced distribution of states. 
The method is closely related to Follow-the-leader (cf. 
\cref{sec:il:passive}), only with a more sophisticated leverage to the 
expert policy. 
In short, it entails the model $\pi_i$ that queries an expert policy as 
in \cref{sec:expertPolicy}, $\pi_\star$, that it's trying to mimic, 
but also ensuring the learned model updates itself in an iterative fashion, 
until it converges. 
The benefit of this approach is that the feature-states that are likely to 
occur in practice are also investigated and as such used to dissuade the model 
from making poor choices. In fact, the method queries the expert about the 
desired action at individual post-decision states which are both based on past 
queries, and the learner's interaction with the current environment.

DAgger has been proven successful in a variety of benchmarks 
\cite{RossGB11,Ross13}, such as the video games Super Tux Kart and Super 
Mario Bros., handwriting recognition and autonomous navigation for large 
unmanned aerial vehicles. In all cases it greatly improved on traditional 
supervised IL approaches.

\subsection{DAgger}
The policy of AIL at iteration $i>0$ is a mixed strategy given as follows:
\begin{equation}\quad\label{eq:il}
\pi_i = \beta_i\pi_\star + (1-\beta_i)\hat{\pi}_{i-1}
\end{equation}
where $\pi_\star$ is the expert policy and $\hat{\pi}_{i-1}$ is the learned 
model from the previous iteration. 
Note, for the initial iteration, $i=0$, a pure strategy of $\pi_\star$ is 
followed. Hence, $\hat{\pi}_0$ corresponds to the preference model from 
\cref{sec:expertPolicy} (i.e. $\Phi^{\text{IL}0}=\Phi^{\text{OPT}}$). 

\Cref{eq:il} shows that $\beta_i$ controls the probability distribution of 
querying the expert policy $\pi_\star$ instead of the previous imitation model, 
$\hat{\pi}_{i-1}$. 
The only requirement for $\{\beta_i\}_i^\infty$ according to \cite{RossGB11} is 
that $\limit{\frac{1}{T}\sum_{i=0}^T\beta_i}{T\to\infty}{0}$ to guarantee 
finding a policy $\hat{\pi}_i$ that achieves $\epsilon$ surrogate loss under 
its own state distribution limit.

\Cref{pseudo:activeIL} explains the pseudo code for how to collect 
partial training set, $\Phi^{\text{IL}i}$ for $i$-th iteration of AIL.
Subsequently, the resulting preference model, $\hat{\pi}_i$, learns on the 
aggregated datasets from all previous iterations, its update procedure is 
detailed in \cref{pseudo:DAgger}.

\input{pseudocode/imitationLearning}
\input{pseudocode/DAgger}

\subsection{Results}\label{sec:ail:expr}
Due to time constraints, only $T=3$ iterations will be inspected here.
In addition, preliminary experiments using DAgger for \JSP\ favoured a 
simple parameter-free version of $\beta_i$ in \cref{eq:il}. 
Namely, the mixed strategy for $\{\beta_i\}_{i=0}^T$ was unsupervised 
with $\beta_i=I(i=0)$, where $I$ is the indicator 
function.\footnote{$\beta_0=1$ and $\beta_i=0,\forall i>0$.}

Regarding \ref{expr:boost:newdata} strategy, \cref{sec:il:passive} showed that 
adding new problem instances did not boost performance for the expert policy 
(which is equivalent for the initial iteration of DAgger). 
Hence, for active IL, the extended set now consists of each iteration 
encountering $N_{\text{train}}$ new problem instances.
This way, the extended training data is used sparingly, as labelling for each 
problem instance is computationally intensive. As a result, the computational 
budget for DAgger is the same regardless of whether or not new problem 
instances are used, i.e., 
$\abs{\Phi^{\text{DA}i}}\approx\abs{\Phi^{\text{DA}i}_{\text{EXT}}}$.

The results for \jrnd{10}{10} box plot of \namerho, are given in 
\cref{fig:active:boxplot} and the main statistics are reported in 
\cref{tbl:IL:stats}. As can be seen, DAgger is not fruitful when the same 
problem instances are continually used. This is due to the fact that there is 
not enough variance between $\Phi^{\text{IL}i}$ and $\Phi^{\text{IL}(i-1)}$, 
hence the aggregated feature set $\Phi^{\text{DA}i}$ is only slightly perturbed 
with each iteration. \Cref{sec:pil:expr} has showed this was not a 
very successful modification for the expert policy; even though it's noted that 
by introducing suboptimal feature-space the preference model is not as 
drastically bad as the extended optimal policy, even though 
$\abs{\Phi^{\text{DA}i}}\approx\abs{\Phi^{\text{OPT}}_{\text{EXT}}}$.
However, when using new problem instances at each iteration, the feature set 
becomes varied enough that situations arise that can be learned to achieve a 
better represented classification problem which yields a lower mean \namerho.

\begin{figure}
\includegraphics[width=\columnwidth]{figures/{j_rnd}/{boxplot_active_10x10}.pdf}
\caption{Box plot for \jrnd{10}{10} \namerho, using DAgger for \JSP}
\label{fig:active:boxplot} 
\end{figure}

\section{Summary of Imitation Learning}\label{sec:il:expr}

A summary of \jrnd{10}{10} best PIL and AIL models in terms of \namerho, from 
\cref{sec:pil:expr,sec:ail:expr}, respectively, are illustrated in 
\cref{fig:all:boxplot}, and the main statistics are given in 
\cref{tbl:IL:stats}. To summarise, the following trajectories were used
\begin{enumerate*}
\item expert policy, trained on $\Phi^{\text{OPT}}$
\item perturbed leader, trained on $\Phi^{\text{OPT}\epsilon}$
\item imitation learning, trained on $\Phi^{\text{DA}i}_{\text{EXT}}$ for 
iterations $i=\{1,\ldots,3\}$ using an extended training set
\end{enumerate*}
As a reference, the \sdr\ MWR is shown at the edges of 
\cref{fig:all:boxplot}.

At first one can see that the perturbed leader ever-so-slightly improves the 
mean for $\rho$, rather than using the baseline expert policy. 
However, AIL is by far the best improvement. With each iteration of DAgger, the 
models improve upon the previous iteration
\begin{enumerate*} 
    \item for \ref{bias:equal} with \ref{expr:boost:newdata} then $i=1$ starts 
    with increasing $\Delta\rho\approx+1.39\%$. However, after that first  
    iteration there is a performance boost of $\Delta\rho\approx-15.11\%$ after 
    $i=2$ and $\Delta\rho\approx-0.19\%$ for the final iteration $i=3$
    \item on the other hand when using \ref{bias:adjdbl2nd} with 
    \ref{expr:boost:newdata}, only one iteration is needed, as 
    $\Delta\rho\approx-11.68$ for $i=1$, and after that it stagnates with 
    $\Delta\rho\approx+0.55\%$ for $i=2$ and for $i=3$ it is significantly 
    worse than the previous iteration by $\Delta\rho\approx+0.75\%$
\end{enumerate*}
In both cases, DAgger outperforms MWR
\begin{enumerate*}
    \item after $i=3$ iterations by $\Delta\rho\approx-5.31\%$ for 
    \ref{bias:equal} with \ref{expr:boost:newdata}
    \item after $i=1$ iteration by $\Delta\rho\approx-9.31\%$ for 
    \ref{bias:adjdbl2nd} with \ref{expr:boost:newdata}
\end{enumerate*}
Note, for \ref{bias:equal} without \ref{expr:boost:newdata},  DAgger is 
unsuccessful, and the aggregated data set downgrades the performance of the 
previous iterations, making it best to learn solely on the initial expert 
policy for that model configuration.

Regarding \ref{expr:boost:newdata}, it is not then successful for the expert 
policy as $\rho$ increased approximately 10\%. This could most likely be 
counteracted by increasing $l_{\max}$ to reflect the 700 additional examples. 
What is interesting, though, is that \ref{expr:boost:newdata} is well suited 
for AIL, using the same $l_{\max}$ as before. 
Note, the number of problems used for $N^{\text{OPT}}_{\text{train,EXT}}$ is 
equivalent to $T=2\tfrac{1}{3}$ iterations of extended DAgger. 
The new varied data give the aggregated feature set more information 
on what is important to learn in subsequent iterations, as those new 
feature-states are more likely to be encountered `in practice.' 
Not only does the AIL converge faster, it also consistently improves with each 
iteration.

\begin{figure}
\includegraphics[width=\columnwidth]{figures/{j_rnd}/{boxplot_summary_10x10}.pdf}
\caption{Box plot for \jrnd{10}{10} \namerho, using either expert policy, 
DAgger or following perturbed leader strategies. 
MWR shown for reference.}\label{fig:all:boxplot} 
\end{figure}

\input{tables/stats.IL.10x10.tex}

\section{Discussion}\label{sec:dis}

Experiments in \cref{sec:pil:expr} show that following the optimal
policy is not without its faults. There are many obstacles to consider in order 
to improve model configurations. 
When training the learning model, there is a trade-off between making the 
over-all best decisions (in terms of highest mean validation accuracy) versus 
making the right decision on crucial time points in the scheduling process, as 
\cref{fig:case} clearly illustrates. 
Moreover, before training the learning model, the preference set $\Psi$ 
needs to be re-sampled to size $l_{\max}$. 
As the effects of making suboptimal choices varies as a function of time, the 
stepwise bias should rather take into account the disproportional amount of 
features during the dispatching process. 
As the experimental studies in \cref{sec:pil:expr,sec:ail:expr,sec:il:expr} 
showed, instead of equal probability (i.e. \ref{bias:equal}) it was much more 
fruitful to adjust the set to its number of preference and doubling the 
emphasis on the second half (i.e. \ref{bias:adjdbl2nd}).
However, there are many other stepwise sampling strategies based on our 
analysis that could have been chosen instead, as here only a 
simplification of the trend from \cref{fig:case} was chosen. 
This also opens up the question of how validation accuracy should be measured, 
taking into account that the model is based on learning preferences, both based 
on optimal versus suboptimal, and then varying degrees of sub-optimality. 
Since ranks are only looked at in a black and white fashion, such that the 
makespans need to be strictly greater to belong to a higher rank, then it can 
be argued that some ranks should be grouped together if their makespans are 
sufficiently close. 
This would simplify the training set, making it (presumably) have fewer
contradictions and be more appropriate for linear learning. Or simply the 
validation accuracy could be weighted in terms of the difference in makespan.
During the dispatching process, there are some significant time points which 
need to be especially taken care of. 
\Cref{fig:case} shows how making suboptimal decisions is especially critical 
during the later stages for job-shop, whereas for flow-shop the earlier stages 
of dispatches are more critical. 

Despite the information gathered by following an optimal trajectory, the 
knowledge obtained is not enough by itself. Since the learning model isn't 
perfect, it is bound to make a suboptimal dispatch eventually. 
When it does, the model is in uncharted territory as there is no certainty the 
samples already collected are able to explain the current situation. 
For this we propose investigating partial schedules from suboptimal 
trajectories as well, since future observations depend on previous 
predictions. 
A straight forward approach would be to inspect the trajectories of promising 
SDRs or CDRs. However, more information is gained when applying AIL inspired by 
the work of \cite{RossB10,RossGB11}, such that the 
learned policy following an optimal trajectory is used to collect training 
data, and the learned model is iteratively updated. 
This can be done over several iterations, with the benefit being that the 
scheduling features that are likely to occur in practice are investigated and 
as such used to dissuade the model from making poor choices in the future.

The main drawback of DAgger is that it quite aggressively queries the expert, 
making it impractical for some problems, especially if it involves human 
experts. A way to confront that, \cite{Kim13,Judah12} propose frameworks to 
minimise the expert's labelling effort.
Or even circumvent the expert policy altogether by using a `poorer' reference 
policy instead (i.e. $\pi_\star$ in \cref{eq:il} is suboptimal) 
\cite{ChangKADL15}.


\section{Conclusion}\label{sec:con}

The \sdr s remain a popular approach to scheduling, as they are simple to implement and quite efficient. By inspecting optimal schedules, and investigating the probability that an optimal dispatch could be chosen by chance, and by looking at the impact of choosing sub-optimal dispatches,  some light is shed on how SDRs vary in performance. Furthermore, the problem instance space was varied, giving an even better understanding of the behaviour of the SDRs. This analysis, however, also revealed that creating effective \dr s from data is by no means trivial. 

Our analysis has shown that training data must be generated not only from optimal trajectories but also from the actual dispatching rule undergoing the training. The learning of effective dispatching rules is an iterative process with many calls to the expert (MIP solver). Techniques for reducing computationally expensive calls to the expert need to be addressed in the future. This study has been structured around the \jsp\ scheduling problem. However, it can easily be extended to other types of deterministic optimisation problems that involve sequential decision making. Applying the technique to other types of
scheduling problems is of interest as future work.


\bibliographystyle{spmpsci} 
\bibliography{references} 

\end{document}