Annealing Harmony Search Algorithm to Solve the Nurse Rostering Problem

: A real-life problem is the rostering of nurses at hospitals. It is a famous nondeterministic, polynomial time (NP) -hard combinatorial optimization problem. Handling the real-world nurse rostering problem (NRP) constraints in distributing workload equally between available nurses is still a difficult task to achieve. The international shortage of nurses, in addition to the spread of COVID-19, has made it more difficult to provide convenient rosters for nurses. Based on the literature, heuristic-based methods are the most commonly used methods to solve the NRP due to its computational complexity, especially for large rosters. Heuristic-based algorithms in general have problems striking the balance between diversification and intensification. Therefore, this paper aims to introduce a novel metaheuristic hybridization that combines the enhanced harmony search algorithm (EHSA) with the simulated annealing (SA) algorithm called the annealing harmony search algorithm (AHSA). The AHSA is used to solve NRP from a Malaysian hospital. The AHSA performance is compared to the EHSA, climbing harmony search algorithm (CHSA), deluge harmony search algorithm (DHSA), and harmony annealing search algorithm (HAS). The results show that the AHSA performs better than the other compared algorithms for all the tested instances where the best ever results reported for the UKMMC

because of the wide range of hard and soft constraints in addition to the variety of individual requests that change from one rostering period to another depending on the needs of the staff nurses.
Nurse rostering is a popular research areas due to several factors that are intensifying the NRP, including an international shortage of nurses due to high turnover, inconvenient working environment, and the high risks that nurses face due to direct interaction with COVID-19 [2]. According to [3], "Nine million shortage estimated in 2014 is predicted to decrease by two million by 2030". An effective nurse roster algorithm can help achieve better working rosters, which might decrease high turnover due to an inconvenient working environment. The aim of this research is to introduce a better nurse rostering algorithm by investigating and exploring metaheuristic algorithms and their hybridization to solve the NRP.
Inspired by the relationship with the music spontaneous creation measure, the harmony search algorithm (HSA) is an evolutionary algorithm that was presented in 2001 by Geem et al. [4] to improve instrument pitches to achieve better harmonies. The original HSA is classified as a population-based search method that has few local-based search aspects. The HSA is similar to other populationbased search algorithms because it satisfactorily distinguishes the superior solution search space when compared to existing local-based search algorithms [1]. A challenge for the HSA is to balance the processes of diversification and intensification. Several attempts have been made in the literature to hybridize the original and enhanced HSA with other local search algorithms to improve its performance [5]. The HSA has been hybridized with other heuristic-based algorithms, such as hillclimbing [6,7], great deluge [8], simulated annealing [9,10], and particle swarm optimization [11]. For the process of diversification and intensification, the hybridization of population-based methods with local search-based methods has been considered by many researchers to strike the balance between these processes [12][13][14]. More details about the algorithms that have been hybridized within the framework of the HSA and the hybridization of the HSA with other metaheuristics can be found in [14,15].
To overcome some of the drawbacks of the original HSA, the research in [16] presented a modification to the HSA called the enhanced harmony search algorithm (EHSA). In [16], the modifications applied to the EHSA were the use of the SCSPA mechanism [17] in building up the HM and the use of a dynamic updating mechanism for the parameters of the harmony memory considering rate (HMCR) and the pitch adjustment rate (PAR).
The main contribution in this research is the novel hybridization introduced between the enhanced harmony search algorithm (EHSA) with a simulated annealing (SA) called the annealing harmony search algorithm (AHSA). AHSA then used to solve a real-world NRP belonging to the medical centre of National University Malaysia (UKMMC). To show the strength of the introduced algorithm, the AHSA is evaluated and compared against previous variants of the HSA in [6,8,16] in addition to another algorithm from the literature that hybridizes the HSA with the SA and is called the HAS [9]. This research is an extension of previous works aimed at improving the EHSA [16] by hybridizing it with other local-search algorithms, such as the hill-climbing [6] and great deluge [8] algorithms. These algorithms were selected because of their diverse acceptance criteria and the different mechanisms to the use to avoid being trapped in local optima. enthused by the annealing process in metallurgy. Annealing is the practice of cooling solid items gradually until they are frozen at the configuration of minimum energy. SA is considered a local search-based algorithm. Researchers in [18] introduced SA, which was inspired by the annealing physical process in statistical mechanics to solve complex optimization problems. To employ SA to solve optimization problems, a suitable disorder mechanism, cost function, solution search space, and cooling plan are needed to escape local optima and find high-quality solutions. Isken et al. [19] discussed the implementation of SA as a solution for solving the NRP in terms of the time of day and urgency of demand in a post anaesthetic care unit. Thompson (1996) [20] applied SA using noncontinuously available employees for shift scheduling. One hundred forty-four test problems were used to examine the proposed method, and another twenty test problems were used to evaluate the SA performance with an integer programming (IP) model. Bailey et al. (1997) [21] employed SA and the genetic algorithm (GA) to generate the rosters of 27 nurses over a period of 6 weeks with 3 daily shifts. Each nurse who is categorized into one out of three levels of skills is assigned to one of four teams. In general, they found that SA was capable of yielding a better-quality solution in less computational time than GA. In [22], stepwise adaptation of weights (SAWing) and noising with SA is proposed to overcome the problem of weights adjusting on constraint binding in solving the NRP. SAWing was used as an intelligent tool to escape local optima by changing the weights to offer precedence to the constraints that were confirmed to be more difficult than others. The method of noising used with SA was random adjusting based on the weights of the evaluation function. The proposed method was tested on a dataset from several Welsh hospitals that included 11 instances. Comparing the results of using SA alone, SA with SAWing, SA with the noising method, and solutions produced manually showed that SA with the noising method performs better than all the other methods. In [23], two solutions to the NRP were introduced: the exact approach and the metaheuristic approach. For metaheuristics, SA was used to further enhance the initial solutions and to solve the small instances for optimality in a short time. The proposed methods were tested on real-world and benchmark NRP datasets.
HSA applied successfully to solve NRP, recently as in [24,25]. Due to the great success of the HSA in approaching different challenging optimization problems, it has been hybridized with other metaheuristic algorithms to overcome its weakness. Some surveys of hybridizing the HSA with other metaheuristic algorithms and its applications can be found in [26][27][28]. According to [29], there are two main ways for hybridization: (i) including some components of one metaheuristic to another and (2) combining components or techniques from other fields into metaheuristic algorithms as operations research and other artificial intelligence methods.
Researchers in [30] hybridized the HSA with particle swarm optimization (PSO) by substituting the random function with the global-best function of PSO. To enhance the HSA diversification ability, the HSA was hybridized with a memetic algorithm approach to solve the problem of unit commitment for electrical power systems. To reach the minimum cost of the systems life cycle, the HSA was hybridized with SA in [31] to optimize the source of residential electrical load through mixed energy systems. In [32], the HSA was also hybridized with SA to all HSAs to accept the inferior harmonies using the probability of the SA temperature. The proposed method was tested on the IEEE CEC 2014 benchmark in addition to a computer vision real-life problem. eight instances. It also has seven hard constraints and five soft constraints. For further details about the mathematical model, constraints, coverage demand, types of nurses, and penalty values of the UKMMC dataset, please refer to [33,34]. 3 W is the set of weights; P is the set of possible shifts of night patterns;

Mathematical Model of UKMMC
R ds is the coverage demand for day d on shift s; d p is the starting date of night pattern p.

Decision Variables
Eqs. (1) and (2) show the decision variables of every nurse regarding night shift patterns.
In this mathematical model, we have four working shifts: morning (M), evening (E), night (N) and a day off (o). These 4 shifts are represented by the set S = {1, 2, 3, 4}. To simplify the modelling, the day off is considered the fourth shift.

UKMMC Constraints
As with other scheduling problems, constraints can be classified as hard and soft. Hard constraints cannot be violated under any condition while soft constraints can be violated, but we add a penalty for each constraint violated based on how many times the violation occurs. The UKMMC dataset has the following hard and soft constraints.

Hard Constraints
1). The coverage demand for nurses must be fulfilled.
2). One-shift a day is the highest number of shifts for each nurse each day.
3). It is obligatory to have at least one senior-nurse for each working-shift.  4). An isolated working shift is not allowed where the nurse has off before and after a single working shift. 5). For a roster period of 14 days, 12 days is the maximum number of working days whilst 10 days is the minimum. 6). A nurse is not allowed to work more than four working days consecutively. 7). If the nurse works four consecutive night shifts, then s/he is entitled to get two days off.

Soft Constraints
1). Equal allocation for working days and days off for all working nurses. 2). It is highly preferable to assign a day off during weekends (Saturday or Sunday). 3). Four consecutive morning shifts followed by a day off is highly preferable. 4). Four evening shifts consecutively followed by a day off is highly preferable. 5). Four night shifts followed by evening shift or day off after the 2 days off that follows the patterns of night shift (NNNNxxE) or (NNNNxxx). (1) and (2) are the most desired constraints to be satisfied according to the extracted information from the UKMMC Nursing Department. For soft constraints 1, 2, 3, 4, and 5, the associated penalty values (PVs) are 100, 100, 10, 10, and 1, respectively. The constraint that has a higher PV specifies the significance of that constraint. For example, it is more important that the constraint with PV 100 is satisfied than that with PV 10 or 1.

Coverage Demand
For the coverage demand required to cover selected instances of the UKMMC dataset, the number of nurses, senior nurses, and required nurses for the morning (M), evening (E), and night (N) shifts on weekdays and weekends are presented in Tab. 1.

Objective Function
The objective function for the UKMMC dataset is to reduce all penalty values that occur because of violating the soft constraints. Eq. (3) presents the objective function used in this research.

Hybridizing the EHSA with SA for the NRP (AHSA)
The EHSA is an enhanced variant of the basic HSA that uses the semicyclic shift patterns approach (SCSPA) introduced in [18] to generate feasible shift patterns instead of a random mechanism to fill the harmony memory (HM). The SCSPA is proposed to address the night shift since it is the most problematic shift to be assigned due to the extra constraints that it has. SCSPA first assigns predesigned night-shift patterns in a cyclic manner and then assigns morning and evening shift patterns in a noncyclic manner. Fig. 1 shows the main steps of SCSPA. After the allocation, the process is done and then the repair mechanism is used to repair the infeasible initial solution produced by the SCSPA. Regarding the patterns containing night shifts, we select a shift pattern with the same night shifts and different M and E shifts than the allocated pattern to replace one of the allocated patterns that contain night shifts. One rule that governs the replacement and swapping process is that they should be performed within the same duty roster. If the shift pattern to be replaced or swapped is from the duty roster of the first week, it must be swapped with another shift pattern from the first week N shift patterns. The same rule applies when replacing or swapping a shift pattern in the duty roster of the second week. An example of applying the moves of the repair mechanism of the SCSPA is presented in Fig. 2.
Rather than using fixed parameter values for the HMCR and the PAR, it is recognized that the performance of any population-based method is affected by the values of its parameters, which are problem dependent. Thus, to overcome these shortcomings, several methods have been proposed to control the problem of parameter tuning. To overcome the above-addressed weaknesses in the EHSA in this research, the HMCR and PAR were updated dynamically instead of through the fixed methods used in the basic HSA [35]. The initialization mechanism of the HM is replaced with a semicyclic mechanism in generating the initial population. The employment of a dynamic method to adjust the values of the HMCR and PAR parameters throughout the instance solving, instead of the fixed parameter values of the BHSA. This modification is based on the idea of adaptive harmony search (AHS) Hasancebi et al. which is inspired by the dynamic modification of the PAR and the BW in [36]. At the beginning of the search process, the parameter values of the HMCR should be lower in the range of HMCR min and HMCR max . This facilitates more exploration. Towards the end of the search process, the value of the HMCR should be increased to facilitate more exploration. Furthermore, the value of the PAR should be as high as possible in the range of PAR min and PAR max at the beginning of the search process. This facilitates making higher exploration. Towards the end, these values of the PAR are decreased to facilitate exploring the candidate solutions. Fig. 3 shows the pseudo-code of the improvisation step only with the proposed modifications. In this modifying step, the parameter values of the HMCR and the PAR keep changing through the execution based on Eqs. (4) and (5). These equations are based on the work in [37] to solve RNA secondary structure prediction. The main idea in this research is to hybridize the EHSA in [17] with standard SA [38]. In this section, the hybridization of the EHSA with SA is discussed in detail. The SA procedure used in this research is the same as the standard SA procedure used in [38]. The basic pseudo-code of SA is presented in Fig. 4.  The harmony annealing search (HAS) [9] is implemented and compared to the AHSA. The HAS is an HSA variant that applies the idea of the dynamic HMCR and PAR along with SA. The AHSA is compared with the HAS because both algorithms use the same idea to update the parameter settings of the HMCR and the PAR. Furthermore, the AHSA and HAS use SA in their solution framework. However, the only two differences are: 1) In the initial stage of filling the HM, the HAS uses the random mechanism, whereas the AHSA uses the SCSPA mechanism.
2) The HAS applied SA to improve the solutions that come from the random selection based on the probability of the HMCR, while the AHSA uses SA to improve the best solution after the EHSA reaches one of the stopping criteria.

Results
For each instance, 20 runs were performed, and the obtained results are reported. The same experimental design is used for both the AHSA and HAS to fairly compare them. The stopping criteria are as follows: 1) 1000 iterations (based on a preliminary test), 2) PV = 0, and 3) 600 s of execution time. A 32-bit Intel laptop with 1.73 GHz, and 2 GB RAM were used to perform the experiments. The parameter values of HMCR min and HMCR max are 0.1 and 0.95, respectively [39]. The parameter values used for PAR min and PAR max are 0.01 and 0.99, respectively, [40,41]. For the BW, we set it experimentally to be −5 and 5. Section 5.1 presents the results obtained by the AHSA, 5.2 presents the results obtained by the HAS, and 5.3 presents the comparison between the AHSA, HAS, EHSA, CHSA, and DHSA.

Results for the UKMMC Dataset Using the AHSA
The results obtained by applying the AHSA to solve the UKMMC dataset are presented in Tab. 2, where the best, average, worst PVs, standard deviation (Stddv.), desirable patterns (DPs) and execution time are recorded.
Based on the obtained results, AHSA was able to produce an optimal solution for small instances as produced roster didn't break any of the soft constraints. For medium and large instances, very promising results were reported. Soft constraints 1 and 2 that have the higher value didn't break.

Results for the UKMMC Dataset Using the HAS
The computational results obtained by applying the HAS algorithm to solve the UKMMC dataset are presented in Tab. 3, where the best PVs, average PV, worst PVs, standard deviation (Stddv.), desirable patterns (DPs) and execution time are recorded.    6 depicts the average PVs for AHSA and the other four compared algorithms. It is apparently that AHSA was able to generate the lowest PVs for all tested instances. Even though HAS could obtain feasible solutions for the tested instances; the quality of the obtained results was very low.       Based on the results in Tab. 4, the AHSA succeeds in producing better rosters than the EHSA. The AHSA remarkably reduced the PVs compared to the EHSA PVs. The AHSA needs almost twice the time that is required by the EHSA. Due to the acceptance criteria of SA, it escapes the local optima, and as a result, it further decreases the PVs for all instances. For the small instances (CIUC, SGY5, and MD1), the AHSA shows that it is very promising in reaching optimal solutions for NRP as it reaches zero PVs in each run. For medium instances (NICU, N50, ED) and large instances (SGY and CICU), AHSA generates very high-quality rosters with very low PVs. To further evaluate the performance of the AHSA, it was compared to previously published work of other HSA hybridized algorithms the CHSA [6] and DHSA [8]). The results in Tab. 4 show that the AHSA outperformed the CHSA and DHSA in all of the evaluation criteria except for execution time, which increased slightly. This is because the AHSA obtained better PVs and a higher number of DPs for all instances.
The obtained results in Tab. 4 show that the AHSA performed better than the HAS for all instances. This is because the AHSA was able to generate better-quality rosters with less execution time than HAS. Doing the selection randomly by SA out of HM caused two shortcomings: (i) premature convergence because of selecting new-solution vectors out of HM that has been improved by SA, and (ii) the additional consumed time by random selection to generate an entire new harmony and then improve it by SA.

Conclusion
In this paper, a novel hybridization between EHSA and SA is introduced. The AHSA is applied to solve a real-world dataset collected from UKMMC Malaysia. AHSA was employed to construct ideal duty rosters with minimum PVs and the maximum number of DPs, taking into consideration all of the associated constraints. The AHSA is then compared to the EHSA alone, the CHSA, and the DHSA algorithms to demonstrate its performance. The results showed that the AHSA yields the best results. It is observed that the success of the AHSA is because of the ability to maintain the balance between diversification and intensification in the solution search pool, in addition to the strong mechanism of SA in escaping local optima. Moreover, the comparison between the AHSA and the HAS (both use SA to enhance the constructed solution but in different stages) showed that the AHSA performed extremely better than the HAS in all UKMMC instances for all comparison criteria. This is due to the premature convergence that occurred for the population of the solutions in the HM caused by running the SA several times in the HAS.
For future work, it is highly recommended to hybridize the EHSA with enhanced variants of SA, which might further enhance the mechanism of the HSA and help construct high-quality rosters for nurses.