A data-driven hierarchical MILP approach for scheduling clinical pathways: a real-world case study from a German university hospital

Facing economic pressure and case-based compensation systems, hospitals strive for effectively planning patient hospitalization and making efficient use of their resources. To support this endeavor, this paper proposes a flexible hierarchical mixed-integer linear programming (MILP)-based approach for the day-level scheduling of clinical pathways (CP). CP form sequences of ward stays and treatments to be performed during a patient’s hospitalization under consideration of all relevant resources such as beds, operating rooms and clinical staff. Since in most hospitals CP-related information needed for planning is not readily available, we propose a data-driven approach in which the structure of the CP to be scheduled including all CP-related constraints is automatically extracted from standardized hospital billing data available in every German hospital. The approach uses a flexible multi-criteria objective function considering several patient- and hospital-related aspects which makes our approach applicable in various scenarios. Furthermore, in contrast to other approaches, it considers several practically relevant aspects ensuring the implementability of the scheduling results such as multiple ward stays per hospitalization and gender-separated room assignments. Regarding the treatment resources such as operation rooms and clinical staff, it considers the eligibility of resources for treatments based on information such as special equipment or qualification and represents complex resources individually to avoid disaggregation problems. To allow solving the resulting complex and large-scale scheduling problem for realistically dimensioned problem instances, we propose a hierarchical two-stage MILP approach involving carefully designed anticipation components in the first-stage model. We evaluate our approach in a case study with real-world data from a German university hospital showing that our approach is able to solve instances with a planning horizon of 1 month exhibiting 1088 treatments and 302 ward stays of 286 patients. In addition to comparing our approach to a monolithic MILP approach, we provide a detailed discussion of the scheduling results for two practically motivated scenarios.


Introduction
In Germany, hospitals cause around 25% of the total costs of the health care system (Destatis 2014). To reduce these costs of about 76bn Euro per year, the German Diagnosis-Related Groups (DRG) system was established in 2003. In contrast to the earlier length-of-stay-based reimbursement system, hospitals now receive a diagnosis-dependent lump sum per case. As a consequence, hospitals have a strong incentive to reduce their patients' lengths of stay and to effectively manage clinical processes and resources.
At the time of writing, however, every third hospital is still unprofitable (Augurzky 2015) which, together with demographic change, yields a great need for increasing the efficiency of hospitals. According to Villa et al. (2009), clinical pathways (CP) are among the most promising instruments for achieving an efficient utilization of both human and physical hospital resources. CP are specific sets of time-constrained treatments and ward stays to be performed between a patient's admission and discharge to cure a certain disease. Besides increasing the transparency and the standardization of medical processes, CP are instrumental in hospitalization scheduling and controlling (Jacobs 2007). Despite this often-stated potential, however, CP and patient flow scheduling are still waiting to form an integral part of day-to-day practice in many hospitals.
According to Roeder et al. (2003), Kirschner et al. (2007), Küttner and Roeder (2007), and Salfeld et al. (2009), an important obstacle for a widespread utilization of CP is the fact that manually creating CP forms a tedious endeavor. In addition, since CP are mostly used in a descriptive way, only few CP formalisms discussed in the literature enable both the automatic extraction of CP structures from hospital data and their employment for scheduling. Furthermore, many important parameters affecting pathway schedules such as emergency patient arrivals, (exact) treatment durations, resource availability or complications arising in surgeries are subject to uncertainty. Moreover, given that the pathway-scheduling problem has to consider complex constraints as well as all relevant resources, it constitutes a complex largescale optimization problem. As a consequence, many pathway-scheduling approaches rely on simplifications and operate on a high level of aggregation which makes putting the resulting schedules into practice difficult: The approach from Saadani et al. (2014), for example, completely ignores bed resources. Finally, while hospital scheduling typically considers multiple objectives related to economics, patients and staff members, many approaches focus on one-dimensional performance indicators, typically based on cost minimization or profit maximization (Hulshof et al. 2012).
Motivated from a case study with a German university hospital, this paper addresses most of the issues discussed above and proposes a tactical pathwayscheduling approach intended to support hospital case managers responsible for the coordination and scheduling of treatments and ward stays. The main contributions of this paper can be summarized as follows.
First, we propose a data-driven approach to pathway scheduling. Our approach is based on scheduling relevant pathway information automatically extracted from standardized hospital billing data available for every hospital in Germany using the pathway-mining approach introduced by Helbig et al. (2015). As a result, one of the most time-consuming data preparation steps necessary for pathway schedulingdefining pathways and their constraints-is automated to a large extent. To the best of our knowledge, the resulting approach is the first to combine process mining techniques with optimization-based decision support for scheduling clinical pathways.
Second, we propose a day-level pathway-scheduling approach considering all pathway constraints and hospital resources on an adequate level of detail which is flexible enough to consider multiple scheduling criteria: it considers multiple ward stays during hospitalization and gender-separated room assignments. Moreover, it explicitly accounts for resource eligibility constraints for treatments, e.g., based on specific feature and skill requirements. In addition, it allows considering important resources individually while aggregating identical resources if appropriate. Finally, to support scheduling pathways for a broad range of possible objective structures, the objective function in our approach involves a weighted combination of multiple hospital-and patient-related criteria including fairness and workload-balancing criteria.
Third, to tackle the resulting complex large-scale scheduling problem for realistic problem instances involving up to 300 patients with a planning horizon of a full month, we propose a two-stage mixed-integer linear programming (MILP) approach based on a hierarchical decomposition: The first-stage model determines the admission dates and schedules complex treatments with a high duration. Given instructions from the first model, the second model schedules the ward stays, assigns all patients to gender-separated rooms and schedules the remaining treatments. To align the first-stage instructions with the second-stage requirements, the first-stage model involves anticipation components to account for the gender separation requirement and to consider pathway constraints involving the full set of treatments.
Fourth, we evaluate the overall approach including the automatic CP extraction in a set of experiments conducted with real-world data from a German university hospital involving 286 elective patients and 1088 treatments. On the one hand, we evaluate our hierarchical approach in comparison with a monolithic mixed integer linear programming (MILP) approach: we compare solution times as well as solution quality of both approaches; furthermore, we investigate the effect of employing different degrees of anticipation within the hierarchical approach. On the other hand, we provide a detailed evaluation of two scheduling scenarios: in the first scenario, our approach is used to establish balanced resource utilization throughout the planning horizon; in the second scenario, the main goal is to obtain a high resource utilization at the beginning of the planning horizon.
Note that in this paper, we do not explicitly address the stochastic character of the pathway-scheduling problem: on the one hand, our approach only deals with scheduling elective patients for which most planning-relevant information is known ahead. Nonetheless, one of the goals of the planning approach is to establish balanced resource utilization and resource buffers to ensure that emergency cases can be handled adequately. On the other hand, our approach employs conservative point estimates for uncertain treatment times to obtain schedules which are robust with regard to treatment time variations. Like the other authors dealing with clinical pathway planning (see e.g., Gartner andKolisch 2014 andBurdett andKozan 2018), we assume that additional uncertainties such as absence of medical staff, resource breakdowns or stochastic patient arrivals are addressed by embedding our pathwayscheduling approach in a rolling planning process.
The remainder of this paper is structured as follows. Section 2 provides an overview of the rich literature dealing with offline operational scheduling in hospitals and identifies possible objectives, efficient modeling techniques and aspects to be considered in clinical pathway scheduling. Section 3 introduces the constraint-based CP concept according to which the scheduling is carried out and sketches the automated CP mining approach mentioned above. In Sect. 4, we introduce our hierarchical MILP approach consisting of an aggregated first-stage model for scheduling complex treatments and determining the admission day and a detailed second-stage model for scheduling the remaining treatment and assigning patients to rooms. In Sect. 5, we explain the real-world dataset from a university hospital used in our computational experiments. In addition, we illustrate the results from mining CP for the hospital's largest department. Section 6 presents the results from using our MILP approach for scheduling the extracted CP for all elective patients of the department of urology. To demonstrate the effects of varying the weights in the multi-criteria objective function, we discuss two scenarios: in the first scenario, the resource allocation is balanced across the planning horizon; in the second scenario, resources are preferably allocated at the beginning of the horizon. Section 7 concludes the paper and points to further research opportunities.

Related work
We employ the taxonomy of planning decisions in health care introduced by Hulshof et al. (2012) to structure the following literature review. According to this taxonomy, the problem considered in this paper can be characterized as an offline operational planning problem: complete CP for elective patients are scheduled within a planning horizon of about 1 month considering activities at the level of individual patients and resources. Since our approach takes all kinds of hospital resources into account, the following literature review does not only regard CPrelated literature but also considers more general hospital resource scheduling approaches to identify aspects relevant for our research. Since the assignment of staff to shifts and inventory management issues are not in the scope of this paper, articles dealing with these topics are not included in the following review. Furthermore, the review is restricted to articles involving mathematical programming as solution technique published after 2005. Table 1 shows the related work classified by the combinations of considered resources and treatments.
Since surgeries form the greatest source of income for most hospitals (Denton et al. 2007) and at the same time operating rooms (OR) constitute the most expensive type of resource using more than 10% of a typical hospital budget (Jebali et al. 2006;Chaabane et al. 2008), there is a rich body of literature dealing with surgery scheduling. The main goals of surgery scheduling are to reduce costs (overtime costs, penalties for idle time or fixed costs for opening an OR) and to increase the utilization of the available OR. This is achieved by either allocating surgeries to a given set of OR (Lamiri et al. 2007;Lamiri et al. 2008a, b;Fei et al. 2009;Denton et al. 2010;Riise and Burke 2011;Meskens et al. 2013), sequencing a set of given surgeries (Denton et al. 2007;Cardoen et al. 2009a;Cardoen and Demeulemeester 2011) or both (Jebali et al. 2006;Testi and Tànfani 2009;Roland et al. 2010;Batun et al. 2011;Marques et al. 2012;Clavel et al. 2017). In most cases, mathematical programming, heuristics, column generation or combinations of these methods are used to find a good or even optimal solution with respect to a typically multi-criteria objective function. With regard to the scope of the present article, the main insights from the surgery scheduling literature can be summarized as follows: Overtime should be avoided and the well-being of the medical staff should be considered (Roland et al. 2010;Meskens et al. 2013), a flexible pooling strategy of OR has significant benefits (Batun et al. 2011), the eligibility of  Treatments Vlah Jerić andFigueira (2010, 2012) and Schimmelpfeng et al. (2012) Operating room and beds Pham and Klinkert (2008), Cardoen et al. (2009b), Augusto et al. (2010), Fei et al. (2010), Chow et al. (2011), Banditori et al. (2013), Sun et al. (2013), Vancroonenburg et al. (2013), Ceschia and Schaerf (2014) and Li et al. (2015) Clinical pathways Vissers (2005), Conforti et al. (2011), Helbig (2011, Gartner and Kolisch (2014) and Saadani et al. (2014) Business Research (2019) 12:597-636 601 operating rooms for surgeries may depend on special operating room equipment (Denton et al. 2007) and the patients' welfare should be taken into account (Testi and Tànfani 2009). Regarding the solution methods, carefully decomposing the problem into solving multiple small problems often yields a good performance in terms of solution time without having a significantly negative impact on solution quality (Jebali et al. 2006;Lamiri et al. 2007;Lamiri et al. 2008a, b;Fei et al. 2009;Denton et al. 2010;Batun et al. 2011;Marques et al. 2012). For a recent survey on OR scheduling, see Samudra et al. (2016). As stated by Helm and Van Oyen (2014), a proper bed management is necessary to avoid having to dismiss patients because of blocked beds, surgical cancelations and operational chaos. In the literature, two kinds of bed management problems are discussed both of which typically consider medical needs as well as patient preferences: elective admission scheduling and the patient-to-room assignment. Elective admission scheduling problems are solved, e.g., by local search heuristics (Ceschia and Schaerf 2011), MILP (Helm and Van Oyen 2014) or random forest models (Schmidt et al. 2013). Demeester et al. (2010) solve the patient-to-room assignment problem using a tabu search heuristic. Vancroonenburg et al. (2014) develop an integrated model to solve both problems at once to increase both planning flexibility and operational efficiency. Gender-separated room assignment is typically considered as a crucial issue, see Demeester et al. (2010, Ceschia andSchaerf (2011), Schmidt et al. (2013), and Vancroonenburg et al. (2014).
According to Schimmelpfeng et al. (2012), scheduling treatments (encompassing not only surgeries but all types of medical activities) manually often leads to an inefficient resource allocation and can have negative effects on quality of care and patient satisfaction. To avoid this, Vlah Jerić and Figueira (2010) develop a decision support system (DSS) to support the construction of a daily schedule of medical treatments considering available resources and other criteria. They use a scatter search heuristic to construct a set of Pareto-optimal solutions among which the user can interactively choose. Another approach by the same authors presented in Vlah Jerić and Figueira (2012) determines daily treatment schedules considering the availability of medical equipment and physicians. Their approach is based on a multi-objective binary programming formulation for which different solution approaches such as variable neighborhood search, scatter search and non-dominated sorting genetic algorithms are compared. Schimmelpfeng et al. (2012) outline a conceptual framework for a DSS to support the scheduling process of treatments in a rehabilitation hospital. To handle realistic problem dimensions, a hierarchical approach is developed using a daily model for aggregated long-term scheduling and an intra-day model (either time-or resource-oriented) for short-term scheduling.
As argued by Chow et al. (2011), trying to achieve a high OR utilization without considering further surgery-related resources often results in issues such as staff overtime, surgical cancelations and long surgical waiting times. In particular, scheduling surgeries without taking the availability of recovery beds into account often leads to problems. To avoid this issue, several authors propose to simultaneously schedule both the surgery and the following stay in a recovery bed or on a ward. For instance, Pham and Klinkert (2008) propose to allocate hospital resources to individual surgical cases divided into preoperative, perioperative (= surgery), postoperative intensive care unit (ICU) modes and formulate this problem as a MILP multi-mode blocking job shop model. Cardoen et al. (2009b) solve a surgical case scheduling problem in a day-care facility. They formulate a multi-objective model to minimize peaks in recovery bed usage, the occurrence of recovery overtime and violations of staff and patients' preferences. Their approach is based on a combination of column generation and dynamic programming. To analyze the impact of recovery in OR when no recovery bed is available, Augusto et al. (2010) use a Lagrangian relaxation-based approach. Their results show a high benefit of this improved flexibility in resource usage, in particular in case of high demands for recovery beds. Fei et al. (2010) construct weekly surgery schedules for OR considering recovery beds to maximize resource utilization while minimizing overtime and idle time between surgeries. They use a hierarchical two-stage solution approach which first computes the date of the surgeries based on a set partitioning formulation solved by a column generation technique. In the second step, taking into account the availability of recovery beds the daily surgery scheduling problem in which the sequence of surgeries is determined is formulated as a flow-shop problem and solved by a hybrid genetic heuristic. The authors show that their overall approach allows reducing idle-and overtime increasing resource utilization. Chow et al. (2011) use Monte Carlo Simulation to predict bed requirements from historical data. Based on this, they formulate a MILP to reduce peaks in bed occupancy by scheduling surgeon blocks and patient types. Another simulation optimization approach for the surgery scheduling problem aiming at a maximal patient throughput is developed by Banditori et al. (2013). The authors first use a MILP model to determine the number of cases to be treated by surgery groups in each time slot of a month. The model is then fine-tuned with a simulation approach to obtain both robust and easy-toimplement schedules. Sun et al. (2013) propose a weekly scheduling approach to achieve a high OR utilization. The authors formulate a MILP taking into account limited key resources like ICU-beds and workload of surgeons. Vancroonenburg et al. (2013) incorporate OR and gender separation constraints in their approach for the elective patient admission scheduling problem. Their goal is to determine how scheduling surgical and non-surgical admissions impacts room assignments at wards. Ceschia and Schaerf (2014) address a similar patient admission problem considering OR utilization constraints, gender-separated room assignment, a flexible planning horizon and the notion of patient delay using a local search heuristic. Li et al. (2015) develop a MILP-based lexicographic goal programming approach of scheduling OR and beds considering competing resources such as surgical and nursing staff, anesthesiologists and recovery beds.
The most promising way to avoid bottlenecks, idle time and overtime during hospital patient flow is to schedule not only surgeries, possibly augmented with induced bed resources, but each patient's full CP. Vissers (2005) proposes a master surgery scheduling approach for CP involving cardiothoracic surgeries based on computing the patient mix with a MILP model. He divides the care process into four successive steps (medium care unit (MCU), surgery, ICU, MCU) and computes the master schedule for a planning cycle of 4 weeks taking OR, MCU-beds, ICU-beds and nursing staff for ICU-beds into account. An approach to schedule full CP of elective patients for a week hospital 1 is introduced by Conforti et al. (2011). The authors show their results from scheduling clinical services (diagnostics and surgeries) as well as patient admissions to maximize the patient flow using a MILP model in a case study with 20 patients. Helbig (2011) presents a MILP model for the operational hour-based scheduling of interdisciplinary CP containing treatments, surgeries, order constraints and ward stays. The aim of the model is to minimize waiting time and to finish each CP as early as possible within the planning horizon under consideration of practical constraints like gender separation in rooms and a fair distribution of waiting time among all patients. The paper reports results from experiments with artificial small-scale problem instances with up to four patients. A patient flow planning approach to maximize the length of stay (LOS)-based contribution margin for DRG-based reimbursement policies is developed by Gartner and Kolisch (2014). The authors formulate two MILP models to schedule CP, one with fixed and one with variable patient admission dates. They experiment with these models in a static and in a dynamic setting in which the MILP is embedded in a rolling horizon approach. The results dealing with a planning horizon of 28 days and 150 patients show the potential to achieve a higher contribution margin and a significantly reduced time between admission and surgery compared to given manual solutions. Saadani et al. (2014) propose a MILP-based CP scheduling approach in which treatments using different resources are scheduled to minimize the patients' length of stay. Their CP contain multiple treatments but no bed requirements. The authors present results for an instance of 20 patients and a planning horizon of 14 days. Recently, Burdett and Kozan (2018) proposed an approach using constructive heuristics and metaheuristics for scheduling CP operating on a very high level of detail: They schedule not only the day, but also the time of day of the treatments. Their approach takes into account all major treatment resources under consideration of resource eligibility constraints; in addition, it considers multiple ward stays and bed resources on wards.
Summarizing the results of the literature review, the most promising way to improve the patient flow in hospitals is to consider the full CP of each patient. As shown above, scheduling CP involves various complex aspects to be considered by a CP scheduling approach to achieve practically useful and implementable schedules: • Given the multitude of relevant objectives, the approach should allow to include multiple objective criteria related to both patients and resources: -Important patient-related objectives encompass the consideration of LOS (the main driver of schedule-related costs in a DRG system), admission day preferences and hospitalization without delays. -Among the resource-related objectives, minimization of overtime and idle time, establishing a balanced workload and achieving a high level of resource utilization are the most important.
• To obtain implementable pathway schedules, patients should be assigned to concrete rooms taking gender separation into account. • For an adequate representation of resource requirements, scheduling should consider the eligibility of resources for treatments which depends on resource properties such as special qualifications or equipment. Furthermore, to avoid disaggregation problems as sketched in the introduction, important treatment resources such as OR should be considered individually instead of in an aggregated way. • It is unlikely to find a monolithic approach for exactly solving a CP scheduling problem considering all of the aspects described above. However, decomposing the problem into multiple stages, if carried out carefully, can be a promising approach to achieve good results within acceptable solution times also for largescale instances involving more than 200 patients, a 30-day planning horizon and 1000 treatments.
We are not aware of any mathematical programming-based CP scheduling approach taking into account all the aspects mentioned above. As a consequence, this paper proposes a novel hierarchical MILP approach to schedule whole CP with hospital-wide ward stays able to handle the mentioned aspects for real-world problem instances.
An additional important feature of our approach is its data-driven nature: by employing the pathway-mining approach and the constraint-based representation of CP introduced by Helbig et al. (2015), the structure of the CP to be scheduled can be automatically extracted from standardized data available in every German hospital. It becomes clear from the recent survey articles Yang and Su (2014) and Rojas et al. (2016) that process mining for clinical pathways and in healthcare in general forms a very active field of research. As described in Rojas et al. (2016), however, the results of the process mining are typically used for analyzing and improving healthcare processes using techniques from business process management. Our approach in contrast is, to the best of our knowledge, the first one tightly coupling process mining with mathematical optimization for scheduling clinical pathways.

Constraint-based CP representation
Even though the general concept of CP is well-known and widely accepted, there is no standard formalism for representing CP (Vanhaecht et al. 2006). While the formal representation of a CP is of secondary importance if it serves merely descriptive purposes, it is crucial when CP form the basis of scheduling patient flow. Conforti et al. (2011), for example, consider a CP as a set of treatments to be scheduled without imposing inter-treatment dependencies. Gartner and Kolisch (2014) represent the possible schedules for each CP on a directed graph in which the arcs form minimum time lags between the clinical activities. The CP scheduling approach proposed in this work relies on the constraint-based CP representation introduced in (Helbig et al. 2015). We give a brief overview and example of this representation in this section; for details regarding the representation as well as regarding the pathway-mining approach, we refer the reader to (Helbig et al. 2015).
As illustrated in Fig. 1, the time scale of the pathway representation is days and a CP is considered on two levels: the level of ward stays and the level of treatments. On the level of ward stays, a CP forms a chronological sequence of ward stays S. Note that while in many cases, the set of ward stays S is a singleton, there are cases in which ward changes are part of the planned pathway of elective patients: as an example, in the department of urology considered in our case study, there are several cases with complex surgeries involving both a stay in an intensive care unit (ICU) and on a regular ward; in total, 6% of the cases in the dataset involve multiple ward stays.
In the CP representation, a ward stay s 2 S is associated with a set of eligible wards W s ; this means that for stay s, the patient may be assigned to one of the wards in w s . In Fig. 1, for example, the ward stay s 2 may either be assigned to ward w 2 or to ward w 3 . A ward stay is characterized by an admission event v ad s and a discharge event v dis s . The admission event v ad s 1 of the first stay s 1 corresponds to the hospital admission and always takes place at the first day of the CP. All other ward stay events have feasible time intervals indicating when they can be scheduled (e.g., the discharge v dis s 1 can be scheduled from day 3 to day 6, see Fig. 1). The number of days between admission and discharge on a ward forms the length of stay (LOS) on that ward. For every ward stay, the LOS is constrained by a minimum and a maximum number of days. During hospitalization, each ward stay requires a bed resource and it is assumed that there are no gap days between two adjacent ward stays.
On the second level of the CP representation, a treatment t corresponds to a surgery or any other clinical service to be performed to cure a disease. A CP Fig. 1 Constraint-based CP concept contains a set T of different treatments. Each of these treatments lasts at most 1 day and can only be performed if the required amount of all necessary resource types (e.g., surgeons, special rooms, theater nurses or the patient) is available. To facilitate scheduling, treatments are clustered into treatment groups; the treatments within the same group g have to be scheduled at the same day (see e.g., t 2 and t 3 in Fig. 1). Some treatment types may occur in multiple groups which means that they have to be performed once per group (t 2 , for example, occurs in two groups). To ensure medically correct schedules, each treatment group g has a feasible time interval d g ; d g Â Ã relative to the first admission event. In Fig. 1, this interval corresponds to the length of the bar associated with the treatment group (e.g., g 2 has the feasible interval [2;3]). In addition, a medically correct order can be enforced by precedence constraints imposing a minimum and a maximum time lag c g;g 0 ; c g;g 0 h i between two related treatment groups g and g 0 (e.g., in Fig. 1, g 5 needs to fall into the interval [1;3] relative to the scheduled day of g 4 ). Moreover, a group can be associated with a certain ward stay s which means that the group can only be scheduled during that ward stay (e.g., g 1 and g 2 from Fig. 1 can only be scheduled during the first ward stay). Note that by grouping treatments as described, many CPrelated constraints are aggregated and in some cases even considered implicitly. On the one hand, this contributes to a tractable size of the MILP models presented in the next section, on the other hand, the modeling itself is facilitated by the grouping.

A hierarchical MILP approach for scheduling clinical pathways
The pathway-scheduling problem considered in this paper consists of determining the patient's admission and discharge days, the start and end days of the ward stays and the days on which treatments are performed. The scheduling accounts for precedence relations between treatments and considers all critical hospital resources such as beds, medical staff as well as treatment and operating rooms. Bed resources are considered on the level of rooms, accounting for the requirement of gender separation. Treatments are associated with resource requirements. It is assumed that for each requirement there is a set of eligible resources and a given amount of resource time needed for the treatment. Each resource, in turn, is associated with a time capacity for each given day. Identical resources may be considered in an aggregated way; however, to avoid disaggregation problems, important complex resources such as operating rooms are considered individually. With respect to the level of detail considered, the pathway-scheduling problem addressed in this paper is situated between the problems considered in Gartner and Kolisch (2014) and in Burdett and Kozan (2018): both the problems considered in Gartner and Kolisch (2014) and in this paper schedule treatments on the timeaggregated level of days, but in contrast to the problem addressed here, Gartner and Kolisch (2014) do not consider multiple ward stays, bed assignment to genderseparated rooms and resource eligibility constraints. Burdett and Kozan (2018), in contrast, consider multiple ward stays and resource eligibility constraints such that the level of detail for representing resources is similar to our work. However, the problem addressed in Burdett and Kozan (2018) deals with the exact timing of treatments instead of only determining the day on which a treatment is performed. Gartner and Kolisch (2014) point out that the pathway-scheduling problem they address is related to the resource-constrained project scheduling problem (RCPSP), with the important difference that treatments are assumed to take a fraction of a the time unit (in this case: a day) while in project scheduling, activities are typically assumed to take multiple units of time. Using the terminology from Artigues et al. (2015), the model presented by Gartner and Kolisch (2014) can be characterized as a time-indexed formulation, that is, the main decision variables for scheduling activities are binary variables with a day index. The authors assume that each resource requirement is associated with a single set of aggregated resources. As a result, no treatment-to-resource assignment is needed and the resource capacity constraint can be formulated using the time-indexed treatment scheduling variables. The resulting model instances are fairly compact and can be solved to optimality by standard MILP solvers even for realistically sized instances.
Using the analogy to RCPSP, extending the problem addressed in Gartner and Kolisch (2014) by the consideration of treatment-specific resource eligibility is very similar to the extension of the RCPSP to a project staffing and scheduling problem in which resources are considered as flexible or multi-skilled and each task is associated with a set of eligible resources. As discussed in Correia and Saldanha-da-Gama (2015), this setting adds a whole new dimension to the RCPSP since tasks are not only scheduled but also explicitly assigned to resources. Incorporating resource eligibility in MILP formulations for the RCPSP typically results in additional variables and constraints, see, e.g., the formulation proposed in Correia and Saldanha-da-Gama (2015) and, as also discussed therein, makes the problems much harder to solve. Note that the modeling framework for project staffing and scheduling problems proposed in Correia and Saldanha-da-Gama (2015) is based on natural-date variables, that is, on variables representing start times of tasks. Our formulation presented below employs both natural-date variables (e.g., for the start day of ward stays) and time-indexed variables (e.g., for resource assignments). Note that Burdett and Kozan (2018) do not provide a MILP formulation of their very detailed CP scheduling problem-their solution approach uses constructive heuristics and metaheuristics.
Instead of resorting to metaheuristic approaches to deal with the computational complexity of the problem considered in this paper, similar to Schimmelpfeng et al. (2012), we employ a two-stage MILP approach which can be viewed as a constructional hierarchical decomposition in the terminology of hierarchical planning, see, e.g., (Schneeweiss 1998). The first stage determines the ''cornerstones'' of the clinical pathway, that is, the patient admission and discharge dates as well as the dates for the complex key treatments of each clinical pathway. Considering the first-stage decisions as fixed instructions, the second stage solves the full problem sketched above.
To obtain a high-quality solution for the full problem, the first-stage model anticipates major aspects of the full planning problem using anticipation model components obtained by the means of aggregation and relaxation. While the first stage does not determine the ward stays, it anticipates ward stay scheduling under consideration of the gender separation requirement on the level of aggregated room types. Furthermore, while it only schedules complex treatments (under consideration of all required treatment resources), to make this scheduling consistent with the pathway constraints involving all treatments, the first-stage model also involves all remaining treatment groups (neglecting the respective resource requirements). Before presenting the MILP formulations for the first-stage and the second-stage model, we give a short overview of the notational conventions used in the presentation. A full list of all symbols used in the models is given in the Appendix; in addition, all symbols are explained in the description of the models.

Notation conventions and variable domains
We use calligraphic symbols for representing sets, in some cases refined with descriptive superscripts. As an example, P denotes the set of patients and P m denotes the set of male patients. Descriptive superscripts are also used for variables and parameters; subscripts are reserved for indices. Indices and parameters are in lowercase. Greek letters are used for path-relative time information (see Sect. 3).
Decision variables are capitalized and we utilize the following conventions to enhance the readability of the models: the letter X is used for variables for which the value corresponds to a day in the planning period, that is, the domain of these variables is the set D ¼ 1; . . .; D j j f g . Superscripts are used to refer to specific groups of variables, e.g., X ws p;s denotes the start day of ward stay s of patient p. Dayindexed variables are mostly denoted with Y. For example, the binary variable Y ward p;w;d indicates whether patient p stays at ward w at day d or not. Bookkeeping variables representing a duration, a slack or a violation of a soft constraint are denoted with U; again specified with a descriptive superscript and possibly with an additional superscript signifying the direction of a violation. As an example, the variable U adþ p represents a positive violation from the desired admission date of patient p. If not specified otherwise, the bookkeeping variables are continuous non-negative variables. The letter Z is used for integer variables mostly representing assignment decisions-for example, the indicator variable Z room p;s;q takes the value of 1 if patient p stays in room q throughout her ward stay s.

First-stage model
In the following exposition, we present the first-stage model used to determine the admission and discharge dates of the patients as well as the dates of the (complex) key treatments in the clinical pathways. Note that all other decisions considered in the following model such as the start and end days of the ward stays, the aggregated gender-separated room assignments and the scheduled dates of the non-complex treatment groups are not transferred to the second-stage model as instructions but merely serve the purpose of anticipation.

Scheduling admission, ward stays and discharge
For each ward stay s from the ordered set S p ¼ 1. . . S p È É of ward stays of patient p, the decision variable X ws p;s represents the start day of s. The set of constraints (1) both compute the length U los p;s of each ward stay and establish the correct sequencing of each patient's ward stays. Note that in the constraint sets associated with the last ward stay S p of patient p, the index S p þ 1 arises as a subscript of the ward stay start variable. This variable can be interpreted as the first day after the end of hospitalization. The set of constraints (2) enforces the hard lower bound b los p;s as well as the soft upper bound b los p;s along with the upper bound violation variable U losþ p;s for the length of stay of ward stay s of patient p: Every patient p has a desired admission day d desAd p . In the set of constraints (3), the deviation of the scheduled admission day X ws p;1 from this day is registered by the variables U adÀ p and U adþ p both of which are bound by the maximum allowed deviation b devAd enforced by the constraint set (4):

Aggregated consideration of bed resources and gender separation
For anticipation purposes, the first-level model considers bed resources as well as the gender separation requirement on the level of aggregated room types. To model these requirements, the day-indexed variables Y ward p;s;w;d indicate whether on day d; patient p occupies a bed on ward w during ward stay s on one of the eligible wards in set W p;s . The three sets of constraints (5-7) link the bed occupation variables Y ward p;s;w;d to the scheduling of the ward stays: Constraint (5) ensures that the number of bed occupation days equals the length of the ward stays; constraints (6) and (7) To ensure that a patient is assigned to a single ward during the full ward stay, we introduce the binary ward assignment variable Z ward p;s;w indicating whether the ward stay s of patient p is assigned to ward w. The set of bundle constraints (8) ensures that this assignment is unique; constraints (9) ensure that the bed occupation variables only take a value of one for the assigned ward: X While the first-stage model does not involve an assignment of patients to a concrete room, it enforces that for each day, there exists an assignment of patients to rooms respecting the gender separation requirement. To achieve this, the set of rooms Q w on ward w is partitioned into subsets Q w;b according to the number of beds b in the rooms; the set B contains the possible number of beds in rooms. For each day d 2 D, the number of rooms on ward w with b beds used as female (male) room is decided by the integer decision variable Z fRoom w;b;d (Z mRoom w;b;d ). The constraints (10) ensure that on each day d on each ward w and for room capacity b, the total number of male and female rooms does not exceed the number Q w;b of rooms with the corresponding capacity. The constraints (11) and (12)

Scheduling treatments considering pathway constraints
As noted in Sect. 3, our approach aggregates treatments occurring on the same day form a treatment group. In the mathematical model, the integer decision variable X treat p;g represents the day on which the treatment group g from the set G p of all treatment groups of a patient p is scheduled. To anticipate the pathway constraints in the first-stage model, the set G p encompasses all treatment groups (not only those containing key treatments to be scheduled in the first stage).
The scheduling of treatments is linked to the scheduling of ward stays in two ways: first, as modelled in constraints (13), for each treatment group g, there is a time window d p;g ; d p;g Â Ã relative to the admission date of patient p associated with the variable X ws p;1 ; violations of the time window are registered by the variables and U twþ p;g . Second, as modelled by constraints (14), a subset G p;s of the treatment groups G p of a patient p may require to be performed during the ward stay s of p: X ws p;s X treat p;g X ws p;sþ1 À 1 8p 2 P; s 2 S p ; g 2 G p;s : In addition to the relations enforced by (13) and (14), the constraint-based pathway model described in Sect. 4 considers sequence constraints providing an interval c g;g 0 ; c g;g 0 h i between a treatment group g and a successive treatment group g 0 . The set of successors of a treatment group g in the pathway of patient p subject to this type of constraint is denoted with G succ p;g ; the respective constraints in the model can be formulated as follows: X treat p;g þ c p;g;g 0 X treat p;g 0 X treat p;g þ c p;g;g 0 8p 2 P; g 2 G p ; g 0 2 G succ

Modelling treatment resource requirements
To perform the treatments in the aggregated treatment groups, certain resources such as operating rooms, physicians and nurses are needed. As mentioned above, the first-stage model only considers the resource requirement for a subset of treatment groups. In the following exposition, this subset is denoted with G res p . When it comes to the representation of these resource requirements, most publications dealing with pathway scheduling rely on a simple representation of resources and requirements: resources have a unique resource type (or a unique identifier) and the resource requirements are formulated in terms of time needed per resource type (or per specific resource). In reality, however, the eligibility of resources for treatment requirements depends on a set of resource properties such as qualifications (e.g., the capability of conducting a certain surgery) or features (e.g., special equipment only available in certain operating rooms). In this work, this more complex type of representing resource eligibility is modelled as follows: We represent the resource requirements on the level of treatment groups. For each group, the requirements of all treatments in the group are aggregated. The set of resource requirements for treatment group g of patient p is denoted with K p;g ; the amount of resource time needed for the requirement k 2 K p;g is given by the parameter a p;g;k . A resource requirement k is further associated with a set R k of eligible resources.
For each day d and for each eligible resource k 2 K p;g for group g of patient p, there is a binary assignment variable Z res p;g;k;r;d indicating whether resource r is used to cover the resource requirement k of treatment group g of patient p scheduled at day d. The set of bundle constraints (16) ensures that exactly one resource is chosen for each requirement; constraints (17) Note that our formulation of resource requirements and capacities given by the constraints (16-19) is valid both for individual and aggregated resources. Resources can be aggregated if they have the same properties, e.g., the same qualifications or features. For example, in our case study, the OR-nurses can be considered as an aggregated resource. In case of an aggregated resource r 2 R, the resource capacity b res r;d per day corresponds to the sum of the capacity of all individual resources on that day.
Note, however, that in addition to adequately modelling resources with different features and skills, there is another reason for individually considering (complex) resources: the requirements for these resources often comprise multiple hours and it is not possible to switch the assigned resource during a treatment. Aggregating these resources thus easily leads to a situation where there is no feasible disaggregated allocation to individual resources for an aggregated solution.

Objective function
The objective function of the first-stage model involves multiple patient-and resource-related objectives combined in a single function. Depending on the preferences of the decision maker, the objective function coefficients can be chosen accordingly. An important part of the objective function deals with patient delays. The constraints sets (20) and (21) compute the total delay U del p of a patient: constraint set (20) ensured that the total delay of a patient is at least as large as the number of days exceeding the medically required minimum length of stay b los p . In addition to that, the constraints (21) force the total delay to be greater than the sum of the ward stay delays U losþ p;s and the treatment delays U twþ p;t . To establish fairness among the patients and to avoid huge delays of single patients, a min-max objective involving the maximum total delay U maxDel of all patients is used; U maxDel is forced to the maximum of all patient delays by constraints (22). Similarly, to be able to enforce a smooth resource allocation along the full planning horizon, for each resource r 2 R, the constraints (23) The first patient-related term in the first-stage minimization objective function (24) is the maximum total delay U maxDel weighted with the penalty c maxDel . The other patient-related terms involve the individual delays U del p associated with the penalty c del as well as the deviations from the desired admission days U adÀ p and U adþ p of each patient p weighted with the coefficient c ad . The remaining elements of the objective function deal with resource-related aspects: a penalty c bedþ is imposed for each female and male extra bed needed per day. Regarding the treatment resources, overtime and idle time of each resource r 2 R are considered in the objective function two ways: first, for each resource, the maximum overtime U maxOt r and idle time U maxIt r over the full planning period are considered using the penalties c maxOt and c maxIt . Second, for each day d 2 D, the overtime U ot r;d and idle time U it r;d , weighted with the penalties c ot d and c it d , contribute to the objective function:

Second-stage model
Considering the instructions from the first-stage model, that is, the patient admission and discharge dates as well as the dates of the complex treatments, the second-stage model schedules the remaining treatments and performs the detailed allocation of patients to rooms on the wards respecting gender separation constraints. Besides the fact that in the second-stage model, certain decisions are considered as fixed, the main differences between the first-stage and the second-stage model are the level of detail on which the bed resources are modeled and the fact that in the second-stage model, all treatment groups are scheduled under consideration of their resource requirements. As a result, assuming that in the second-stage model, the sets of treatment groups G res p contains all relevant treatment groups instead of only the subset of complex key treatments, the constraint sets (1-4) dealing with scheduling ward stays (implying admission and discharge days) and (13-23) mainly dealing with treatment scheduling and the allocation treatment resources also appear in the second-stage model. In the following, we present the structurally different parts of the second-stage model, that is, the constraints for fixing the first-level decisions, the model component for representing bed resources on a high level of detail and the objective function in which the bed-related terms make use of the variables introduced in the second-stage model.

Enforcing the instructions from the first-stage model
To enforce the instructions from the first-stage model in the second-stage model, we impose the following constraints fixing the respective variables: constraints set (25) fixes the admission day d ad p and the discharge day d dis p of each patient. Assuming that for each patient p, the set of treatment groups scheduled by the first-stage model is denoted with G fixed p and that d treat p;g corresponds to the day on which treatment group g 2 G fixed p is scheduled, constraints (26) fix the respective scheduling decisions in the second-stage model:

Bed resources, room assignment and gender separation based on individual rooms
In contrast to the first-stage model in which on each ward, the rooms with a common bed capacity are aggregated, the second-stage model considers each room individually. For each day d, the binary day-indexed room stay variables Y room p;s;q;d determine whether a patient p stays in room q in ward stay s; Q p;s denotes the set of rooms which can be used by patient p for ward stay s. In analogy to the constraints (5-7) in the first-stage model, constraints (27-29) link the Y room p;s;q;d variables to the ward stay start variables X ws p;s ; enforcing the allocation of a bed in a room for each day of each ward stay: X In the second-stage model, it is assumed that a patient resides in the same room during the full ward stay. The assignment of a room q to a patient p at her ward stay s is represented in the binary decision variable Z room p;s;q ; the constraints set (30) ensures that each patient is assigned exactly one room during each ward stay. Constraints (31) enforce that the values of the day-indexed room assignment variables Y room p;s;q;d are in line with the same-room requirement.
The gender separation requirement means that on each day in each room, all occupants share the same gender. In the model, the binary decision variable Y mRoom q;d is introduced; Y mRoom q;d equals 1 if room q is a ''male room'' and 0 if room q is a ''female room'' on day d. The constraints (32) and (33) enforce that for each day d and in each room q, the gender of the room occupants is consistent with the value of the variable Y mRoom q;d and at the same time ensure the room capacity is respected. It is assumed that the room capacity can be extended by extra beds registered by the variable U mBedþ q;d for beds occupied by males and U mBedþ q;d occupied by females. The set of constraints (34) ensures that for each room and for each day, the upper bound b add q for the number of extra beds is respected and that the gender of the extra bed occupants matches the gender assignment of the rooms on each day:

Objective function
The objective function of the second-stage model is given by (35). While the scheduling objectives of both levels are basically identical, the objective function (35) reflects the main differences between the models on the two levels: first, since the admission decisions are fixed in the second stage, (35) does not consider admission day preferences. Second, since the level of detail of the bed allocation in the second-stage model is higher than in the first-stage model, the corresponding term is adapted: 5 Real-world case study from a German university hospital The research presented in this paper is based on a partnership with a German University Hospital providing us both with a real-world dataset and with an opportunity to discuss our assumptions and results with domain experts such as the case manager of the department of urology considered in our investigation. In 2011, the case study hospital had around 36,000 admissions divided into around 25,000 elective and 11,000 emergency inpatients. Figure 2 shows the amount of elective admissions on each ward in 2011. The department of urology (HA2200) has the highest number of cases; Fig. 3 shows the distribution of these cases grouped by month of admission. We agreed with the partner hospital to base our analysis and our experiments on the peak month March 2011 exhibiting the highest number of cases. In close collaboration with the case manager, we prepared and validated the input data required for our approach including two sets objective function weights representing two different planning scenarios. The remainder of this section describes and illustrates these input data: Sect. 5.1 deals with the cases and their CP information automatically extracted by the pathway-mining approach sketched in Sect. 3. Section 5.2 describes how the data regarding resource requirements and resource capacities were obtained. Finally, Sect. 5.3 explains the planning scenarios used for the computational experiments.

Case data and automated pathway extraction
As explained above, we conduct our investigation with the real-world data for cases admitted in March 2011 to the department of urology of a German university hospital. For these cases, we applied the CP mining approach proposed in Helbig et al. (2015) and shortly sketched in Sect. 3 of this paper. Since this mining approach operates with standardized hospital billing data according to the §21KHEntG 2 which has to be reported by each hospital in Germany, it is readily applicable for each German hospital. The pathway-mining approach identifies similar treatment profiles for each primary diagnosis and creates corresponding homogenous case groups. Based on these groups, feasible intervals for treatments and ward stays as well as order constraints are determined. Note that if only a single case with a certain treatment profile occurs in the database, the mining approach returns treatment chains without scheduling flexibility for treatments and ward stays.
In total, the data in the considered month encompass 286 cases (217 male and 69 female patients) with 90 different primary diagnoses and 69 different DRGs.
2 More details about structure and containing data as well as an example dataset can be found at: http:// www.g-drg.de/cms/Datenveroeffentlichung_gem._21_KHEntgG. Applying the mining approach to the described dataset resulted in 229 different CP among which 136 are based on unique treatment profiles and 93 represent profiles occurring more than once. In total, the cases involve 1088 treatment groups and 302 ward stays; around 6% of all cases involve more than one ward stay. Figures 4, 5, and 6 show different CP for the primary diagnosis benign neoplasm: rectum based on a single treatment profile (Fig. 4), prostatic hyperplasia, a profile occurring three times (Fig. 5) and the most frequent profile (Fig. 6) for the primary diagnosis of malignant tumors of the prostate. The CP depicted in Fig. 4 exhibits two ward stays; the first on the department of urology and the second on an ICU (after the ICU the patient is relocated to a rehabilitation). Figures 5 and 6 exhibit a single ward stay on the department of urology. Each CP contains a specific set of treatments based on the OPS-Codes (operations and procedures) documented in the billing data according to §21KHEntG. Figure 4 shows a unique CP without scheduling flexibility for the ward stays and treatment groups, Figs. 5 and 6 depict profiles occurring more than once exhibiting scheduling flexibility for the treatment groups.

Resource requirements and availability
In addition to the CP structure discussed so far, information regarding required and available resources is needed for our approach. To obtain a realistic set of parameters, we discussed all 259 different OPS-Codes that had occurred during our planning horizon (March 2011) with the case manager of the department for urology who provided us with information regarding the required resources and resource times for each code. Table 2 shows an excerpt of the required time on resource types needed by a treatment.
Since OPS-Codes are very detailed and surgeries usually involve multiple codes, considering each individual set of codes would typically result in a unique type of pathway for each surgery. To avoid this, all surgery codes (starting with a ''5'') appearing at the same day of a CP are aggregated to a single artificial code. To obtain a robust approximation the resource requirements for such a surgery code, we consider the maximum time a single resource unit was required for any of the OPS-Codes aggregated to the surgery code.
This approximation is based on the assumption provided by the case manager that all required resources need to be available for the full surgery duration. This assumption is reflected by all of her resource requirement estimates for individual surgery OPS-Codes. Using the maximum duration over all involved OPS-Codes, however, tends to overestimate the actual surgery duration: Taking all cases from March 2011 considered in our study, on average the duration of a surgery is overestimated by 4.6 min. In accordance with the case manager from the involved department, we chose this slightly pessimistic estimation procedure for the surgery duration and resource consumption instead of a more optimistic estimate to make  the scheduling results more robust against uncertainty with respect to surgery durations. Table 3 shows an excerpt of required resource types and their required time for some example surgeries; each surgery belongs to a primary diagnosis and a CP variant. For example, the surgery with the primary diagnosis C61 in variant 1 is based on the OPS-Codes 557,840 and 560,401. Using our assumptions, the required capacity for physicians and anesthetists is 8 h (2 physicians working 4 h) and 4 h.
The case study department has a capacity of 88 patients during workdays and 73 at weekends. It has 7 single-, 5 two-, 9 three-and 11 four-bed rooms. We assume that small rooms are not used for scheduling elective patients because of their flexibility in handling emergencies. Based on this assumption, we reduced the amount of available beds to 53 by keeping only the 11 four-bed rooms and 3 threebed rooms. With regard to the CP to be scheduled, there are four additional departments in which beds are needed for ward stays. For these departments, we assume a bed capacity based on the observed number patient cases available in form of two-bed rooms.
Regarding the treatment resources, we consider scarce and special resources such as operating rooms and the urography individually. On the one hand, this accurately This is a special kind of operating room located in the department of urology models resource eligibility for treatments requiring special equipment: as an example, the central OR is eligible for 139 treatments whereas the OR of the urology department is only eligible for 66 of these treatments. A more complex example arises in the Urography: There are ten treatments for which both the Uro-1 and the Uro-2 resource is eligible; however, there are 19 treatments for which only Uro-1 and 13 treatments for which only Uro-2 is eligible. On the other hand, considering complex resources individually avoids disaggregation issues such as scheduling three 4-h surgeries on a day while there are only two OR available for 6 h on the same day. In addition, complex resources often exhibit different daily capacities since they may be shared among different departments. Other resources such as human resources with identical skills (e.g., OR-nurses) are aggregated.
In the considered planning month (March 2011), approximately 60% of all patients were elective. Based on this observation, we assume that 60% of the total available time of all treatment resources owned by the urology department can be scheduled for elective patients. Furthermore, the per-day availability of shared hospital resources such as central OR, MRT or computer tomography for the urology department was determined based on interviews with the department case manager. Table 4 shows an excerpt of the capacity available for scheduling elective patients per resource type for each weekday.

Scheduling scenarios considered in the computational experiments
The computational experiments conducted with the data described above are based on two main scenarios differing with respect to the resource allocation objectives. In the first scenario referred to as ''smooth allocation'', the hospital aims at establishing a smooth resource allocation throughout the planning horizon. As noted by the case manager at our partner hospital department, a smooth resource allocation is favored by the staff members and establishes a similar level of flexibility for reacting to emergency cases on each day of the planning horizon. To achieve a smooth allocation, the objective function coefficients associated with the min-max objectives as well as the penalty associated with daily resource overtime are both ten while the other objectives are smaller (see Table 5). Since extra beds on rooms are considered as worse than idle time or patient delays, we set extra bed penalty to five and the delay and idle time penalties to two. In the second scenario referred to as ''early allocation'', the hospital aims at establishing a high level of resource utilization at the beginning of the planning horizon. The scenario can be motivated by the observation from the case manager that the level of uncertainty affecting scheduling-relevant information regarding patient cases and resource availability increases towards the end of the planning horizon. This scenario is particularly useful in a setting where our approach is embedded in a rolling horizon approach. Establishing a high level of resource utilization at the beginning of the planning period aims at providing more flexibility for coping with the higher amount of uncertainty towards the end of the planning period. In addition, freeing future resource capacity can be instrumental for admitting a higher number of elective patients. The set of objective function coefficients for this scenario is depicted in the second row of Table 5; note that in particular, this scenario involves a quadratic decreasing day-dependent penalty for the idle time of resources and zero penalties for the maximum idle time and for the deviation from admission day penalties.

Experimental results
The real-world data described in the previous section forms the basis of a series of computational experiments with our hierarchical approach. The first set of experiments, discussed in Sect. 6.1, aims at assessing the hierarchical approach proposed in this paper. The second sets of experiments aim at discussing the results of the two scenarios in detail.
For the experiments, the MILP models presented in Sect. 4 were implemented AMPL; all treatments exhibiting a duration of more than 30 min (about 50% of the treatments) were considered to form the set G res p of treatment groups to be scheduled in the first-stage model. The model instances were solved with Gurobi 6.5 on an Intel i7-3770 CPU with 3.40 GHz and 8.00 GB RAM using 6 threads. For solving the aggregated first-stage model, the relative MIP gap tolerance was set to 1%, the Early allocation 10 2 0 5 10 0 10 5 þ 31Àd ð Þ 2 100 time limit was set to 2 h and the MIP search strategy was set to focus on improving bounds. The second-stage model was solved using default Gurobi parameter settings.

Evaluating the hierarchical approach
Compared to approaches presented in the literature, our approach addresses a very detailed pathway-scheduling problem considering multiple ward stays, genderseparated room assignment and treatment-specific resource eligibility. In this section, we first compare the hierarchical two-stage approach to a monolithic MILP formulation. Next, we evaluate the effect of ignoring the gender separation requirement on the quality of the solution, and finally, we present experimental results showing the importance of carefully anticipating the second-stage problem in the first stage of our hierarchical solution approach.

Monolithic model vs hierarchical approach
The main reason for introducing a hierarchical approach for the pathway-scheduling problem considered in this paper is that solving a monolithic MILP formulation of the problem with a standard solver is not possible in a reasonable amount of time for realistically sized problem instances. Table 6 depicts the results from experiments with both a monolithic model and the hierarchical approach proposed in Sect. 4 for the dataset described in Sect. 5. For the computation with the monolithic model, we allowed for a maximal solution time of 12 h and varied the maximum allowed deviation b devAd of the desired admission from 3 up to 10 days. The monolithic model instances exhibit 19,500 rows, 112,000 columns and 700,000 nonzeros; the instances of the first-stage model in the hierarchical approach contain 19,000 rows, 29,000 columns and 530,000 nonzeros; the instances of the second-stage model contain 2,200 rows, 6,200 columns and 25,000 nonzeros.
The results in Table 6 show that within the desired solution time of 2 h, the MILP solver was unable to find a feasible solution to any of the instances of the monolithic model. For instances permitting a higher deviation from the desired admission day, this even holds after 12 h of solution time: it was only possible to find integer feasible solutions of allowed admission day deviations of up to 5 days in the smooth allocation scenario and for up to 3 days in the early allocation scenario. Furthermore, only for a single instance, a solution meeting the predefined 1% optimality gap is found within 12 h. For the hierarchical approach, the solution time is less than 10 min for all but one instance for which the solution time is less than 73 min. The last column in Table 6 aims at evaluating the price at which this considerable decrease in solving time comes: For those instances for which a feasible solution could be found for the monolithic model, it depicts the percentage objective function gap between the solution obtained with the hierarchical approach and the best integer solution found for the monolithic model. It turns out that for this gap is between 2.5 and 6.2%; the highest gap is observed in the smooth allocation scenario with the smallest flexibility regarding the deviation from the desired admission day. The objective function gaps between the monolithic and the hierarchical approach can mainly be attributed to the fact that the first stage model determining the cornerstones of the CP schedule is based on an aggregate representation of room assignments and only considers the resource requirements for complex treatments. However, as will become clear in the next section, in contrast to the monolithic model, the hierarchical approach allows obtaining high-quality feasible solutions with a fully flexible admission date within a time limit of 2 h.

The effect of ignoring gender separation
Compared to other pathway-scheduling approaches, our model operates on a higher level of detail. With regard to treatments, our model considers the eligibility of resources based on resource properties (qualifications, available equipment) instead of characterizing resources using a single resource type (nurse, operating room) and stating requirements in terms of these resource types. While the advantage of this more realistic representation for matching treatments and resources cannot be evaluated experimentally in a straightforward way, in the following, we discuss results from comparing different levels of detail of representing the allocation of bed resources. To evaluate the relevance of considering gender-separated rooms, we compare the solution obtained with our approach (a priori consideration of gender separation) to a solution obtained by first scheduling pathways ignoring the gender separation requirement and afterwards determining a gender-separated room assignment is performed based on these schedules (a posteriori consideration of gender separation). Table 7 shows the results from using this comparison for two instances, both based on the ''early bed allocation'' scenario, with an increased penalty of 50 for extra beds: in the first instance, the number of beds and rooms in the department of urology is as described in Sect. 5 (11 four-bed rooms and 3 three-bed rooms) while in the second instance, only 10 four-bed rooms are considered. It turns out that in both instances, the overall solution quality is better and the number of needed extra beds is smaller when considering gender separation a priori. While the objective function differences are comparably small, note that a significant amount of extra beds typically represents solutions which are not implementable in practice.

The role of anticipation in the hierarchical approach
In our hierarchical approach, we carefully designed the anticipation of the secondstage problem in the first stage. As an example, we propose to consider all treatment groups (not only those to be scheduled in the first-stage model) to anticipate the full set of precedence constraints in the first-stage model. In our experiments, this anticipation turned out to be crucial: without this anticipation, almost all the secondstage model instances were infeasible. As depicted in Table 8 showing results from based on the ''early allocation'' scenario with the reduced room capacity setting described above, it is also important to anticipate the gender separation requirement for patient rooms: if the gender separation anticipation is dropped from the firststage model, the final solution obtained from the second-stage model involves considerably more extra beds.

Smooth resource allocation scenario: detailed results
In the first scenario to be discussed in detail, the objective structure favors a smooth allocation of resource utilization over the planning horizon. Table 9 shows the solution times and objective function values for different settings of the parameter b devAd limiting the deviation from the preferred admission day from 0 to 30. The results show that the solution time of the first-stage model is very sensitive to the flexibility regarding the admission day: The model with a flexibility of 5 days is  solved in less than 10 min whereas the model with 10 days needs an eight times higher solution time. Surprisingly, the solution time of the last instance (30 days) is shorter than for the 10-day instance. We suspect that that this is due the fact that the increased admission day flexibility avoids resource bottlenecks. The second-stage model could easily be solved in less than 20 s in all cases but the first. Table 10 provides more details regarding the effects of increasing admission day flexibility: Even a small amount of flexibility considerably improves resource-and patient-related performance indicators. The total amount of overtime for OR, ORnurses and urography is reduced by around 92% in case of allowing the admission day to deviate from the desired day by 3 days and up to 96% in the 30-day case. Considering the OR values in detail (two central OR and one ward-owned OR) shows that the total overtime in the 30-day case is allocated to both central OR with a maximum overtime of 0.7 h (the maximum overtime in the zero-days-case was 6 h). The ward-owned OR which cannot perform all kinds of surgeries needs no additional overtime at all in the 30-day case anymore.
With an increasing flexibility, the idle times of the OR are highly reduced. This contrasts with the comparably small change of idle time for OR-nurses and urography-we suspect that these resources are not scarce. Taking a closer look, we can confirm this observation for OR-nurses and for one of the urography resources. The other urography, however, forms a scarce resource since it is eligible for a wider range of treatments. With regard to patient-related indicators, the maximum delay is reduced from 3 to 0 days and the total amount of all delays decreases from 46 to 0 days. On the other hand, the distance to the preferred admission day increases. Figure 7 depicts the frequencies of the deviations from the desired admission days for all patients for the last instance: from all 286 patients, 231 can be admitted at their preferred day; 7% of all patients are admitted before and 11% are admitted after their preferred admission day. Figure 8 visualizes the allocation profiles of the resource type physician, operating room and urography as given by empirical ex-post data (left hand side) and after the scheduling with an admission day flexibility of 30 days (right hand side). The blue background depicts the available capacity of treatment time and the visible area corresponds to the remaining capacity. The black line on top is the goal ratio of 80% capacity allocation. The red and the green areas are the allocated resource times at each day. For each resource, the allocation obtained by our approach is very smooth compared to the empirical data. The resource allocation of the physicians shows an even workload distribution. Given the large blue area, considering only elective patients, physicians and OR-nurses obviously do not form scarce resources. Regarding the combined allocation profiles of all operating rooms, the results obtained by our approach exhibit a better fit to the available capacity profile with less overtime occurring during the month compared to the empirical data. Since there is no idle time for the central OR, the visible available capacity can be attributed exclusively to the ward-owned OR.
To quantify the smoothing effect obtained by our scheduling approach, let us consider the standard deviation of the daily allocation of each resource: For physicians, the standard deviation change from 12.7 to 7.1, for operating rooms from 7 to 4, for OR-nurses from 8.4 to 5.2 and for the urography from 4.6 to 2.2. Given these values, it can be summarized that the resource allocation is successfully smoothed in this scenario. Furthermore, even under consideration of gender separation, no extra beds are required to hospitalize all 286 elective patients. While the scenario presented in the previous section aimed at obtaining a smooth allocation of resources, the present section focuses on a scenario in which scheduling favors a high resource utilization at the beginning of the planning horizon. As in the first scenario, we vary the maximum allowed deviation from the preferred admission day by changing the parameter b devAd from 0 to 30 days for the first-stage model. Table 11 shows the objective function values, solving times and gaps for each considered value of b devAd . It turns out that even the largest possible instance with completely flexible patient admission days was solved in less than 2 h. We suspect that this considerable decrease in solving time compared to the first scenario results from the fact that in the second scenario, the objective function is less symmetric. For all variants, the second model is solved to optimality in less than 10 min. Table 12 gives an overview on how selected indicators are affected by increasing admission date flexibility from 0 up to 30 days. The patient-related indicators (maximum and individual delay) can be reduced to zero. The resource-related indicators, however, exhibit a more complex pattern: while the total over-and idle time is reduced from the 0-day to the 10-day case, these values increase in the 30-day case. The reason for this behavior is patients shifted from the end of the scheduling horizon to the beginning to avoid most expensive idle times on all resources. Figure 9 shows the absolute frequencies of the deviations between scheduled and the preferred admission days for the case of a fully flexible admission day. The scheduled admission days deviate up to 28 days from the preferred admission day. 65% (186) of the patients are admitted earlier then preferred; 30% (80 patients) are admitted more than 10 days earlier. On average, each patient is 4.5 days earlier admitted.
The effect of the early resource allocation is visualized in Fig. 10 showing the time-dependent resource allocation profile of physicians, OR, OR-nurses and urography based on historical ex-post data (left hand side) as well as based on scheduling with full admission day flexibility (right hand side). The results show that the goal of establishing a high resource utilization at the beginning of the planning horizon could be achieved.   Since at the beginning of the planning horizon, the penalty for idle time substantially outweighs the penalties for overtime, the scarce OR resources tend to exhibit overtime while other resources are highly utilized. Note that an even higher utilization is prevented by the scarcity of bed resources at the admission ward, see Fig. 11 for a visualization of the bed allocation over time. At the first weekend starting at day 5, the OR exhibits a high amount of idle time: in this case, physicians and OR-nurses form the scarce resources, and reducing OR idle time would induce overtime for several resources at once. Interestingly, the scheduling results on the first weekend also reflect the fact that in contrast to OR, the urography resources do not require the presence of OR-nurses: their utilization does not decrease as heavily as the OR utilization; furthermore, given the reduced physician availability on Sundays, the urography utilization is lowest on Sunday.  Summarizing the results of this scenario, the desired maximization of the resource allocation at the beginning of the planning horizon can be achieved. However, putting the described results into practice may turn out to be risky since in case of a rate of arriving emergency patients higher than the expected ratio of 40%, capacity bottlenecks will arise at least during the first week. Nevertheless, our scheduling approach can be used to identify such bottlenecks and to take measures to avoid them or to mitigate their consequences. In the current case, operating rooms and OR-nurses seem to be the scarcest resource types. This becomes particularly obvious after scheduling, showing noticeable available capacity only during the last days of the planning horizon.

Conclusion
This paper presents a hierarchical MILP approach for CP scheduling aiming at practical applicability in several ways. The CP information including the pathway constraints can be automatically extracted from standardized hospital billing data. The data-driven character of our approach substantially reduces the amount of work to be performed to obtain the CP information needed for scheduling. Moreover, the CP scheduling approach considers all relevant resources on an adequate level of detail; in particular, it models bed assignment on the level of rooms and allows considering treatment resources individually to avoid disaggregation problems of complex resources such as operating rooms. Furthermore, it considers several practically relevant aspects such as multiple ward stays per hospitalization, genderseparated rooms and eligibility of resources for treatments (e.g., special equipment in operating rooms needed for a certain surgery or special qualifications of nurses).
To support a broad range of application scenarios and hospital-specific objective structures, our approach employs a multi-criteria objective function accounting for both patient-and hospital-related objectives. Since the detailed representation of practically relevant aspects yields complex large-scale model instances, we propose a carefully designed hierarchical two-stage approach involving anticipation components aligning the first-stage solution with aspects relevant to the secondstage problem.
The experimental results based a case study with real-world data involving a planning horizon of 30 days, 1088 treatments and 302 ward stays of 286 patients show that our approach is capable of obtaining high-quality solutions within a reasonable amount of time: even the most complex instances in which the admission days of the patients can be freely chosen can be solved within less than 90 min. Furthermore, the detailed evaluation shows that our model is capable of considering different scenarios such as achieving a smooth allocation of resources over the planning horizon and establishing a high resource allocation in the beginning of the horizon.
The results presented in this paper offer multiple opportunities for future research: The robustness of the schedules obtained with our approach could be assessed by a simulation study accounting for sources of uncertainty such as such as emergency patients, stochastic treatment times and short-term resource unavailability. Moreover, evaluating our approach in additional settings would be instrumental for understanding in how far the results presented in this paper can be generalized. From a practical perspective, the main goal of our efforts is to embed our model into the information systems of our partner university hospital to provide decision support for hospital-wide scheduling of clinical pathways. While we are aware that this goal is presently out of reach mainly due to the present organizational structure in hospitals and due to a lack of integrated information systems, we believe that providing an approach that is data-driven and considers many practically relevant aspects is a step in the right direction.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.