Application of job shop scheduling approach in green patient flow optimization using a hybrid swarm intelligence

With the increasing demand for hospital services amidst the COVID-19 pandemic, allocation of limited public resources and management of healthcare services are of paramount importance. In the field of patient flow scheduling, previous research primarily focused on classical-based objective functions, while ignoring environmental-based objective functions. This study presents a flexible job shop scheduling problem to optimize patient flow and, thereby, minimize the total carbon footprint, as the sustainability-based objective function. Since flexible job shop scheduling is an NP-hard problem, a metaheuristic optimization algorithm, called Chaotic Salp Swarm Algorithm Enhanced with Opposition-Based Learning and Sine Cosine (CSSAOS), was developed. The proposed algorithm integrates the Salp Swarm Algorithm (SSA) with chaotic maps to update the position of followers, the sine cosine algorithm to update the leader position, and opposition-based learning for a better exploration of the search space. generating more accurate solutions. The proposed method was successfully applied in a real-world case study and demonstrated better performance than other well-known metaheuristic algorithms, including differential evolution, genetic algorithm, grasshopper optimization algorithm, SSA based on opposition-based learning, quantum evolutionary SSA, and whale optimization algorithm. In addition, it was found that the proposed method is scalable to different sizes and complexities.


Introduction
Optimal flow, in terms of patient flow, is critical in providing quality care in healthcare environments, particularly hospitals. Enhancing patient flow is not only beneficial to healthcare providers, it also provides a way to refine health services and improve patient safety, outcomes, and satisfaction (Bacelar-Silva, Cox III, & Rodrigues, 2022;Leviner, 2020;Modi, 2007).
The wide use of complex equipment and technologies in various treatments, especially in hospitals, consumes a large amount of electricity and, thus, can increase CO 2 emissions (Brown, Buettner, & Canyon, 2012;MacNeill, Lillywhite, & Brown, 2017). CO 2 emissions from healthcare in the world's largest economies account for about 4 % of their national carbon footprints. Hospitals consume more energy than other nonresidential buildings per square meter of floor space, in part because of their continuous operation (Gaglia et al., 2007). According to Lancet Commission on Health and Climate Change, greenhouse gas (GHG) emissions of healthcare systems must be included as an indicator in assessments of health and climate (Watts et al., 2017). However, few studies have examined the emissions caused by the healthcare sector as well as potential mitigation strategies (McMichael, Neira, Bertollini, Campbell-Lendrum, & Hales, 2009;Watts et al., 2015).
Scheduling is used to allocate machines for various industrial processes and to determine processing sequences of products. Considering the emphasis on optimizing scheduling to reduce carbon emissions (Fang, Uhan, Zhao, & Sutherland, 2011;Liu & Huang, 2014;Yi, Li, Tang, & Wang, 2012;Zheng, Wang, & Wang, 2014), job shop scheduling problems (JSPs) are among the hardest combinatorial optimization problems even in a deterministic environment, in which all data are assumed to be fixed and precisely known in advance (Jain & Meeran, 1999). Since job shop scheduling optimization problems are NP-hard problems, intelligent algorithms (Jarosław, Czesław, & Dominik, 2013;Pan, 2012;Verma & Kaushal, 2017;Xiao & Konak, 2017) are commonly used. Pollard et. al. (Pollard, Paddle, Taylor, & Tillyard, 2014) proposed a bottom-up modeling framework to help in the decision-making for both cost and carbon in healthcare, using data from a case study in Cornwall, UK. Research findings confirm that a bottom-up model is an efficient tool in the process of estimating and modeling the carbon footprint (CFP) of healthcare.
Most previous studies on patient flow scheduling primarily focused on classical objective functions, such as waiting time, length of stay, and patient throughput (Pham & Klinkert, 2008;Tai & Williams, 2012;White, Froehle, & Klassen, 2011;Wojtys, Schley, Overgaard, & Agbabian, 2009), but ignored environmental concerns. Therefore, this work applied a flexible job shop scheduling approach to model the green patient flow problem (GPFP). Considering that the job shop scheduling problem is NP-hard and complex, an improved intelligent algorithm is developed to prepare an efficient model.
The scientific contributions of this research are delivered through three different novelties in both modeling and solution. Firstly, to the best of our knowledge, it is the first paper that considers carbon emissions when patients are both receiving and waiting for treatments. It is important to note that while patients are kept waiting for the next treatment step, some electrical equipment is normally attached to them, resulting in electricity usage and CO2 production. Not only minimizing the sum of carbon emissions during the two periods is novel, but it is also critical from an application standpoint. Secondly, this study addresses carbon footprint in the patient flow, using an approach analogous to the FJSP. Thirdly, to solve the optimization problem we propose an improved evolutionary algorithm.
The remaining contents of this research are organized as follows. In Section 2, a literature review on the current studies is addressed. In Section 3, the case study is introduced. In Section 4, a bi-criterion green patient flow flexible job shop scheduling problem is described, and its mixed-integer programming model is constructed. In Section 5, the proposed algorithm, i.e. Chaotic Salp Swarm Algorithm Enhanced with Opposition-Based Learning and Sine Cosine (CSSAOS), is described. In Section 6, the datasets, parameter settings, and computational results are described, and the sensitivity analysis in terms of the number of patients in each category of ESI is presented. Finally, Section 7 discusses the findings of the research and draws conclusions based on the outputs of the research.

Related works
The delivery of healthcare services produces a surprising amount of greenhouse gas (GHG) emissions, which play an unequivocal factor in climate change and worldwide warming. Previous studies have addressed the environmental impacts of GHGs caused by the healthcare sector, primarily contributed to the large energy consumption of treatment to procedures. To illustrate, Gilliam et al. (Gilliam, Davidson, & Guest, 2008) estimated direct CO 2 emissions from laparoscopic surgeries, while Ryan and Nielsen (Ryan & Nielsen, 2010) determined the 20-year global warming potentials of three common anesthetic gases, including sevoflurane, isoflurane, and desflurane, through clinical scenarios to estimate the impacts on the environment. The consumption and generation of energy are associated with significant damages to the climate, environment, and, consequently, economy. In spite of the fact that power system dispatch incorporates a significate part in GHGs, especially CO 2 emissions, carbon capture power plants as a basic critical low-carbon generation option will have a vital effect on power system operation and dispatch. Ji et al. (Ji et al., 2013) introduced a model for low-carbon power system dispatch incorporating carbon capture power plants, demonstrating its effectiveness and validity using numerical examples based on an IEEE 118-bus tested system.
Based on a United Nations Framework Convention on Climate Change (UNFCCC) report, electricity production is responsible for 22 % of GHG emissions, 3 % of which is due to electric consumption by hospitals (Eckelman, Sherman, & MacNeill, 2018). In the United States, the healthcare system produces about 8-10 % of global GHG emissions (Chung & Meltzer, 2009;UNEP, 2012) while Canada is responsible for around 25 % (UNEP, 2012) climate change is one of the most important problems facing public health, healthcare services will continue to significantly generate GHGs. Therefore, reducing the environmental impacts and, of course, GHGs caused by the healthcare sector is one of the key responsibilities of the health sector in preventing global warming (Rossati, 2017). A previous study (Becker, 2012) estimated that American hospitals contributed 5.5 % of the total delivered energy to the commercial sector.
Concerning reports on the usage and amount of energy consumed in the healthcare sector (Becker, 2012;CBECS, 2012), optimization of this energy utilization can have a large impact on economic profit and, more importantly, on reducing greenhouse gases. Since the healthcare system has a high demand for electricity, such as for lighting, heating, ventilation, equipment, and air conditioning (Bi & Hansen, 2018;Bujak, 2010;Chirarattananon, Chaiwiwatworakul, Hien, Rakkwamsuk, & Kubaha, 2010;Chung & Meltzer, 2009;Renedo, Ortiz, Mañana, Silio, & Perez, 2006), medical centers have been recognized as one of the largest energy consumers and, thus, emitter of GHGs in the world. Particularly, hospitals consume more energy than other nonresidential buildings per square meter of floor space, as a result of their continuous operation (Gaglia et al., 2007). Based on a study by Yale School of Medicine (New Haven, CT) and Northeastern University (Boston, MA), the US healthcare system is a top producer of GHGs due to its energy consumption for electricity and heating (Review, 2016). To comply with the international (Kyoto and Paris) and national conventions to prevent global warming by reducing and managing emissions in such places is more probable. Ulli Weisz et al. (Weisz et al., 2020) proposed numerous untapped possibilities for reducing GHGs in healthcare services and determined six concrete steps toward sustainable healthcare that apply to most industrial countries. Moreover, life cycle assessments (LCA), as the most established approach to estimate GHG footprints from various treatments, individual products, locations, and industries, have been performed for associated industries (Belboom, Renzoni, Verjans, Léonard, & Germain, 2011;Eckelman, Mosher, Gonzalez, & Sherman, 2012;McGain & Naylor, 2014;Thiel et al., 2015;Usubharatana & Phungrassami, 2018), specific pharmaceuticals (McAlister et al., 2016;Parvatker et al., 2019;Wernet, Conradt, Isenring, Jiménez-González, & Hungerbühler, 2010), and medical procedures (Campion et al., 2012;Connor, Lillywhite, & Cooke, 2010;Danesh-Meyer, 2011;MacNeill et al., 2017).
Scheduling methods are based on operations research techniques, including optimization, mathematical modeling, forecasting, stochastic processes, and queue model. These techniques are used in scheduling staff, managing the patient flow, planning surgeries, and setting appointments.
Patient flow scheduling continuously remains one of the foremost vital issues within the healthcare framework (Hall, 2012). Liang et al. (Liang, Turkcan, Ceyhan, & Stuart, 2015) developed a discrete event simulation and mathematical programming model to evaluate the operational performance in an oncology clinic, which showed to reduce the total working times of clinics and patient waiting times while balancing resource utilization. Gupta and Denton (Gupta & Denton, 2008) introduced a state-of-the-art appointment scheduling system to manage access to service providers and demonstrated its potential for novel applications of Industrial Engineering and Operations Research (IE/OR) techniques. Burdett and Kozan (Burdett & Kozan, 2018) proposed a flexible job shop scheduling (FJSS) model to effectively utilize hospital beds, operating rooms (OR), and other treatment spaces, considering patients, beds, hospital wards, and healthcare activities as jobs, single machines, parallel machines, and operations, respectively. The latter researchers developed hybrid and constructive metaheuristic algorithms to solve the FJSS problem and verified the potential of the scheduling model via integration in actual hospital information systems. Erhard et al. (Erhard, Schoenfelder, Fügener, & Brunner, 2018) reviewed the quantitative methods for physician scheduling in hospitals, including physician scheduling problems (e.g. staffing, rostering, re-planning, and personnel qualification) as well as shift types. The review, including 68 publications in operations research and management science fields, revealed the gaps that require future research activities.
As a result of increasing patient care and satisfaction, while reducing costs, Wright and Mahar (Wright & Mahar, 2013) examined the effect of centralizing scheduling decisions over departments in a clinic. The performance result of the centralized model showed an improvement in nurse management by 34 %, reduced overtime by 80 %, and minimized costs by just under 11 %. Other works (Piroozfard, Wong, & Wong, 2018;Zhang, Gu, & Jiang, 2015) on manufacturing emphasized GHGs reduction, especially CO 2 (environmental influence). For instance, Zhang et al. (Zhang et al., 2015) displayed a low-carbon scheduling model for the flexible job shop problem). The model was solved using a hybrid non-dominated sorting genetic algorithm II, and model efficiency was considering production components (i.e. makespan and machine workload) and environmental impact (i.e. carbon emission proven with some well-known benchmark instances and an actual case. In another research, Wang, et al. (Wang, Ding, Qiu, & Dong, 2011) displayed a lowcarbon production scheduling system to minimize total CO 2 during the whole planning horizon, which was verified by computational experiments considering renewable energy.
Consequently, the main point of this investigation was to develop a mixed-integer programming model as a flexible job shop scheduling problem that can optimize an environmental-based objective function (total CFP) all over the scheduling time.

Case study
Our case study was conducted in Bushehr Heart Hospital (BHH) in southern Iran. The BHH has ten care units, including triage, cardiopulmonary resuscitation (CPR), inpatient emergency department (IED), coronary care unit I and II (CCUI and CCU II), post coronary care unit (PCCU), intensive care unit I and II (ICU I and ICU II), Catheterization Laboratory (Cath lab), and operating rooms (ORs). Also, it has two administration units, namely reception and discharge units, which have been conceptualized as workstations.
Based on the Emergency Severity Index (ESI), patients were categorized into resuscitation (ESI1), emergent (ESI2), urgent (ESI3), less urgent (ESI4), and non-urgent (ESI5) cases (Tanabe, Gimbel, Yarnold, & Adams, 2004). Using the obtained ESI of all patients, the routes that patients must follow from arrival time to discharge from the hospital were determined and are shown in Table 11 (in Appendix). As an example, a patient who is categorized in ESI1 category is immediately moved to CPR. If treatment is successful, the patient is transferred to the inpatient unit and then to the Cath lab. A patient with acute cardiovascular disease (ACS) and experiencing severe chest pain is categorized into ESI2 and is taken to Cath lab for angiography, then if necessary, receives Percutaneous Coronary Intervention (PCI). A patient who is categorized in ESI3 is transferred to the inpatient unit, but if needed angiography, he/she is taken to the Cath lab. A patient categorized into ESI4 is admitted to the emergency department for up to a 6-h stay and will be discharged if the treatment is successful, else he/she is transferred to the inpatient unit. A patient who is in the ESI5 category is referred to an outpatient clinic. The whole process is conceptually considered as patient flow. The average time of using the medical equipment in each care unit according to one patient and the amount of electricity consumed by the equipment (in KW/h) are given in Table 10 and Table 11 (in Appendix), respectively.

Mathematical formulation of green patient flow scheduling
In this research, an application of flexible job shop scheduling in green patient flow management (GPFM) was studied. A flexible job shop scheduling problem (FJSP) could be a generalized form of classical job shop scheduling (Pezzella, Morganti, & Ciaschetti, 2008), where each operation can be executed on a qualified machine. To comply with the methodology, we considered patients as jobs, treatment operations as operations, and medical equipment as machines.
To model the problem, we have used the analogy of the flexible job shop problem. It means that we see the problem as an FJSP instance, where patients resemble jobs, and treatment operations were seen as operations of a specific job. However, we did not intend to use all modeling aspects of the FJSP. In fact, the FJSP analogy helped us to understand the system and also helped in defining the solution procedure.
The FJSP is known to be an NP-hard combinatorial optimization problem (Garey, Johnson, & Sethi, 1976). FJSP consists of two subproblems, the assignment and the scheduling problems. Using this FJSP analogy, we have considered that patient flow problem comprises two sub-problems, which are patient routing and sequencing. The patient routing problem, with the pseudocode represented in Algorithm 4, is concerned with assigning proper medical equipment to treatment operations from a set of eligible medical equipment to execute the treatment operations; the second sub-problem (sequencing) including algorithms 2 and 5 of section 5.6 involves ordering all treatment operations on all medical equipment. The total amount of CO 2 delivered by medical equipment is calculated from both processes.

Assumptions of the proposed method
The green patient flow as a flexible job shop problem is defined as follows. Consider a set of n given patients J (1). Simultaneously executing more than one treatment operation on one piece of medical equipment is not permitted, i.e. to avoid overlap in treatment operations. (2). A piece of medical equipment cannot be used in more than one treatment operation. Thus, more than one patient may be under the same treatment operation only if the required equipment is available for each of them, but a single piece of equipment cannot be shared between two or more operations or patients. (3). Treatment operations are non-preemptable, i.e. once a treatment operation is commenced on medical equipment, it cannot be hindered or delayed. (4). Medical equipment or patients are independent of each other, i.e.
there is no relation among distinctive medical equipment or patients. In any case, priority connections and mechanical arrangements must be considered between the treatment operations of the same patients. (5). A boundless buffer size between medical equipment is accepted. (6). Medical equipment and patients are accepted to be accessible from the start. (7). Medical equipment breakdown and preventive maintenance are not considered, i.e. types of medical equipment are ceaselessly accessible. (8). Emitted CFP per kilowatt-hour is accepted to be constant. (9). The length of stay in each care unit is equal to the length of bed occupancy in that care unit.

Mathematical model
To find the optimal solution for the green patient flow management scheduling (GPFMS) at BHH, mixed-integer linear programming (MIP) model was formulated. The GPFMS is similar to the flexible job shop scheduling problem (FJSP), which was used as a general framework for modeling patient flow.
The following lists include the indices, parameters, decision variables, and the mathematical model considered in this work. In addition, the constraints, notation, and parameters of the mathematical model are provided. Indices: Medical equipment related to treatment operation of bed q (1,2,…,mq) Care and administration unit index (1, 2,…, K)

Parameters:
Aj Arrival time of patient j m Total number of medical equipment n Total number of patients sj Total number of treatment operation s for the j th patient MC l,j A set of capable medical equipment for the l th treatment operation of patient j P i l,j Processing time for the l th the treatment operation of patient j that should be executed on medical equipment i Pmi Power consumption of medical equipment i in the processing condition EF Quantity of emitted CFP per kilowatt-hour B It is assumed as a big number C k The capacity of the care unit (Number of beds in the treatment unit(

Decision variables:
The mathematical model consists of the objective function and constraints to capture the reality of the system. The modeling requirements of FJSP are considered in constraints defined in Eqs. (4, 6, 7, 11, 12). Constraints (8, 9, 13) represent the assignment sub-problem while the sequencing constraints are defined in constraints (5, 10, 14, 15).
MC l,j ⊂M, ∀j = 1, 2, ⋯, n; ∀l = 1, 2, ⋯, s j (11) The objective function, defined in Eq. (3), accounts for the minimization of the total emitted CFP obtained by calculating the emitted CFP of medical equipment for each patient and the carbon produced due to the patient's waiting time in the care units. While it is possible to assign different weights to the two components of the objective function, we have used here equal weights to reflect the hospital authority's preference not to prioritize either component. Eq.
(4) guarantees that the beginning time of the first patient's treatment operation will be longer than the arrival time. Eq. (5) ensures that the start time of a preceding treatment operation (BO i l,j ) plus its execution time (P i l,j ) is less than or equal to the start time of the beginning time of the consequent treatment operation of the same patient to fulfill the priority relationship between diverse treatment operations of the same understanding whereas guaranteeing no relation between treatment operations of distinctive patients.
Eq. (6) guarantees that medical equipment cannot commence preparing diverse treatment operations at the same time, i.e. to prevent overlap in medical equipment. Eq. (7) guarantees that the medical equipment for processing the l th the treatment operation of patient j (O l,j ) is empty and its previous treatment operation O l− 1,j is already processed, i.e. to prevent overlap in treatment operations. In Eq. (8), the appropriate type of medical equipment is determined for each treatment operation of patients. Eq. (9) assigns the treatment operations of the patients to their medical equipment, and subsequently, all treatment operations are arranged on all medical equipment while taking into account the priority of treatment operations on medical equipment. Eq. (10) limits each treatment operation to be executed on a single capable medical equipment with one priority, i.e. each assigned medical equipment for processing a treatment operation has one priority.
Eq. (11) ensures that the eligible set of medical equipment for executing O l,j is from the given set of medical equipment MC l,j . Eq. (12) specifies that once a treatment operation is commenced, it should be performed without interruption, i.e. treatment operations are nonpreemptable. Eq. (13) indicates that the number of occupied beds in each hospitalized care unit shall be less than or equal to the capacity of that hospitalized care unit. Eq. (14) shows the waiting time of patient j on an inpatient bed q. Eq. (15) indicates that the start time of O l,j on medical equipment i and start time of medical equipment i in priority p is bigger than or equal to zero. Eq. (16) express that the finish time of O l,j on medical equipment i should be bigger than or equal to zero. Eqs. (17)-(18) specify that the decision variables including X l,j,i,p and h l,j,i are binary.
As stated in Section 4, the objective function minimizes carbon emitted both during the treatment process and during the time patients are waiting for the next treatment step. As the objective function has two components, it could be modeled with a bi-objective optimization model and then solved using multi-objective metaheuristic algorithms. However both components are of environmental concern and are independent of each other, i.e. an increase in one component does not result in a decrease in the other. In fact, a multi-objective solution algorithm could be used if the optimal decisions need to be taken in the presence of tradeoffs between two or more conflicting objectives.
Although in an FJSP, machine idle time is of concern, in this research beds and equipment idle time have not been concerned. The reason is that there is no CO 2 emission as long as they are not used in the treatment process. The usage time of equipment is based on the patient's ESI, as shown in Table 11.

Proposed CSSAOS algorithm
To solve the problem at hand, an evolutionary algorithm, called Chaotic Salp Swarm Algorithm Enhanced with Opposition-Based Learning and Sine Cosine (CSSAOS), is proposed to solve the presented objective flexible job shop scheduling problem (FJSP). The SSA, recently introduced by Seyedali Mirjalili et al. , is prone to problems, such as a slow convergence rate and local optimal solution, like most metaheuristic algorithms. To solve these problems, SSA is integrated with the chaotic maps, sine cosine algorithm (SCA), and opposition-based learning (OBL).
In this area, the proposed strategy is presented, at that point, the encoding and decoding of a solution are completely clarified as a vital portion of the algorithm. Hence, a stepwise outline of the calculation is given. To depict the proposed algorithm, an outline of SSA, SCA, OBL, and chaotic maps are brought as takes after.

Salp swarm algorithm
SSA mirrors the swarm behavior of salps, alter their position design utilizing quick agreeable changes to rummage around for nourishment (Anderson & Bone, 1980;Sutherland & Weihs, 2017). The populace of salps is partitioned into two bunches: leader and followers. The leader flies at the front of the chain to assign another movement for the rest of the salps, or devotees, to mimic. The position of salps is characterized in an n-dimensional search space, where n is the number of variables of a given problem and food source F is the target. Thus, the position of all salps is stored in a two-dimensional matrix called x. Based on the taking after equation, the position of the leader is overhauled: where x 1 j is the position of the primary salp (called leader) in the j th dimension; F j is the position of the food source within the j th dimension; ub j is the upper bound of the j th dimension; lb j is the lower bound of the j th dimension and c 1 , c 2 , and c 3 are random numbers. The coefficient c 1 is the most important parameter in SSA because it balances exploration and exploitation, which is characterized as follows: where l is the current iteration; and L is the maximum number of iterations. Newton's law of motion is used to update the position of the follower as: where i ≥ 2, x i j is the position of the i th follower salp within the j th dimension; t is time; v 0 is the initial speed; and a = v final v0 where v = x− x0 t . Since the time for optimization is an iteration, the discrepancy between iterations is equal to 1. Considering v 0 = 0, the position of the i th follower salp in the j th dimension is expressed as follows: Logistic map Sinusoidal map Tent map Bernoulli map

Sine cosine calculation for upgrading the leader's position
The SCA could be a population-based optimization method and a modern metaheuristic calculation (Mirjalili, 2016), the solutions are updated based on the sine or cosine function as in equations (23): where P i is the destination solution, X i is the current solution; || demonstrates the absolute value; and r 1 , r 2 , r 3 and r 4 are random variables. The parameter r 1 is upgraded utilizing Eq. (24) to adjust exploration and exploitation (Mirjalili, 2016), as follows: where T is the maximum number of iterations; a is a constant; and t is the current iteration. The r 2 is a random variable utilized to discover the direction of the movement of the next solution (i.e. if it moves towards or away from P i ). Moreover, r 3 is a random variable that gives random weights for P i to stochastically emphasize (r 3 > 1) or deemphasize (r 3 < 1) the impact of desalination in characterizing the distance. Moreover, r 4 is utilized to switch between the sine and cosine functions as in Eq. (23).

Chaotic maps
In recent studies, chaotic generators have been chosen over random number generators (RNGs) as RNGs are not completely random (Caponetto, Fortuna, Fazzino, & Xibilia, 2003). Therefore, this work employed various chaotic maps to update the followers' position and form a new solution, as listed in Table 1.
The chaotic maps are used to update the followers' position. For each follower i, its next position is calculated using Eq. (25).
where N(i +1) is calculated using a chaotic map taken from Table 1; and x i j is the position of the i th follower salp within the j th dimension. Other than that, the salp population X is built by the chaotic maps utilizing Eq. (26).

Opposition-based learning
OBL describes a contrary solution to the current solution, and after that assesses the fitness function to accept or reject the new solution (Tizhoosh, 2005), as follows: where x i j is the opposite vector from the real vector x i j ; and x i j is characterized as a real number over the interval x i j ∈ [L j , U j ].

Chaotic salp swarm algorithm enhanced with opposition-based learning and sine cosine
Within the proposed strategy, SCA is utilized for upgrading the position of the leader, the chaotic maps to update the position of followers, and OBL for a better exploration of the search space, generating more accurate solutions. Based on the 10 chaotic maps provided in Table 1, different CSSAOS algorithms, CSSAOS1 to CSSAOS10, are introduced.
The pseudo-code of the CSSAOS algorithm is outlined in Algorithm 1.

Chromosomes representation
The FJSP is an NP-hard problem consisting of two sub-problems, which are the assignment and the scheduling problems. The FJSP analogy is utilized to see the system in which a different number of jobs (patients) is be processed (treated) on different numbers of machines (treatment units) at the same time.
In this research, two chromosomes are defined to represent a solution. First, the routing chromosome captures the path that will be followed by patients based on their ESIs. Second, the ranking chromosome is used to prioritize bed allocation to patients in all units.
The whole process can be summarized in three steps as follows: (1). Define the patient path chromosome and priority chromosomes.
(2). Schedule bed allocation to the patient based on priority chromosomes. (3). Determine the amount of carbon consumed by each patient, which is the total carbon produced by the services and waiting time.
In the first step, the algorithm sets a uniformly distributed random number r ∈ U(0, 1) to each gene of the chromosomes (random numbers are different for the routing and ranking chromosomes). The chromosomes are divided into n (total number of patients) segments, each representing one patient. The number of genes in each segment is equal to the number of treatment units that each patient will go through in the treatment operation.
In the second step, the beds are assigned to the patients by decoding the first chromosome using Eq. (28): where ⌊.⌋ is the floor function; and q l and q h indicate the first and last bed index of the corresponding care unit, respectively. The amount of carbon dioxide produced is calculated in the third step via Algorithm 6. The most rationale of the proposed calculation is portrayed in Algorithm 2, which consists of six steps, each with a specific set of operations. Both routing and ranking chromosomes are initialized with random numbers, which are defined as Grand. The rank chromosome is valued using Algorithm 3.
Algorithm 3. Pseudo-code of setGenesOrder() END FOR 6. //assigns a rank value to each gene based on its randomly-set value 7.

END DO
Bed assignment is based on patients' priorities. A patient with higher priority should be allocated to a bed before other patients. Algorithm 4 is used to determine the priority of patients. This algorithm defines the ordered set of patients based on their priority. The priority then will be used for sequencing patients and allocating beds.  In allocating beds to the patient, it is critical to consider the plausibility of patient boarding, due to the unavailability of required beds in the destination care unit. Based on the assignment sub-problem of the FJSP, Algorithm 5 is used to allocate beds to the patient and determine if the patient needs to wait for an unoccupied bed. This algorithm ensures that simultaneous execution of more than one treatment operation using a specific medical device is not allowed.

DO 2.
FOR each patient p in patientsOrder 3.
FOR each unit required for patient p 4.

5.
Allocate bed b to patient p 6.
SET patient p in boarding 9.

END DO
The timing of the whole care process needs to be scheduled. Using the concept of the scheduling sub-problem of the FJSP, Algorithm 6 determines the total time that a bed is occupied by a patient and the boarding time of a patient. It shows that the length of stay of a patient in each treatment unit is equal to the duration of bed occupancy in that unit. Also, the priority relationship between different treatment operations shows that there is no connection between different treatment operations. If the next ward does not have an unoccupied bed, the patient will wait in the previous unit until a bed becomes available in the destination unit.

DO 2.
FOR each patient p 3.

END DO
The amount of carbon consumed due to the use of electrical equipment during services and boarding time is calculated using Algorithm 7.
FOR each bed assigned to each patient j 3.
The produced carbon is calculated by using Eq. (29) 4.

5.
The total carbon produced per bed is calculated by using Eq. (30) 6. End DO Eq. (29) shows the amount of carbon dioxide produced during the treatment and boarding of each patient for each bed. Eq. (30) determines the total CO 2 produced from the treatment of a patient.

A non-trivial example
Suppose we have 5 patients with 4 treatment units and 15 beds, and 5 treatment pathways are defined based on the ESI of patients as represented in Table 2.
The bed in each treatment unit is marked with a numerical index. Beds with index numbers 1 to 4 belong to treatment unit 1; beds with

Treatment pathways ESI
treatment unit1 → treatment unit 3 →treatment unit2→treatment unit 4 1 treatment unit 2 →treatment unit 3→treatment unit 1 2 treatment unit 1 →treatment unit 3→treatment unit 4 3 treatment unit 2→treatment unit 4 4 treatment unit 1→treatment unit 2 5 index numbers 5 to 6 belong to treatment unit 2; beds with index numbers 7 to 10 belong to treatment unit 3; and beds with index numbers 10 to 15 belong to treatment unit 4. Also, assume that patients 1 and 5 are categorized in ESI1, patient 2 in ESI4, patient 3 in ESI3, and patient 4 in ESI2 categories. The arrival times of patients 1-5 are 67, 1034, 108, 97, and 19 min, respectively. The chromosome for the above five patients, generated using setRandomGene(), is illustrated in Fig. 1. As seen in Fig. 1, the chromosome is divided into five segments, each segment represents one patient. Therefore, the chromosome represents five patients. The number of genes in each segment is equal to the number of treatment units that the patient will go through in the treatment process.
In the next stage of the procedure, the beds are assigned to patients by decoding the following matrix: In Fig. 2, the first row of the matrix shows the treatment pathways for the five patients. The second row represents the index number of the first bed in the corresponding treatment unit (B lm ), while the third row represents the index of the last bed of the treatment unit (B ln ). The chromosome is represented in the fourth row of the matrix. Equation (31) is used to convert the matrix to the decoded chromosome: where ⌊ ⌋ is the correct component function; B lm is the index number of the first bed; and B ln is the index number of the last bed of the corresponding treatment unit l. Fig. 3 presents the allocation of beds for each patient.
To consider the priority of each patient in assigning beds and performing the required treatment, the ranking chromosome is defined in Fig. 4.
As shown in Fig. 5, using rank-based rankGenes(), the priority of the patient is determined.
For example, patient 5 with rank 1 in locus 3 is the first to be allocated on bed number 3, as seen in Fig. 5. If the target gene of the corresponding segment of the decoded allocation chromosome is in locus 1, the bed will be allocated directly. Otherwise, if the target gene of the corresponding segment is in locus 2 or greater, all the genes before the target should be allocated as well to ensure that the patient is allocated all the required beds in the corresponding treatment units. Using the same procedures, all beds are allocated to all patients.
In Table 3, cell [82, 1-3] shows that a patient in ESI1 (row 1) goes to unit 1, and the total time of the care process will be 82 min. If available, beds 1, 2, or 3 can be allocated to the patient. The other cells can be read     accordingly. As seen in Table 4, patient 5 arrives at time 19 min in unit 2 and is allocated bed 2. He then goes to unit 4 with no waiting time and then back to unit 1 at time 4421. The care process is finished in unit 1 at time 4503, then the patient goes to unit 4 and is allocated bed 14. Since there is no vacant bed, the patient must wait 131 min (his boarding time) until the bed is available at 4634. The whole care process is finally completed at time 8964. In order to calculate the total carbon footprint for the patient, Eq. (30) is used.

Parameter tuning
To detect the optimum level for parameters of algorithms, the Taguchi method (Li & Kwan, 2004) was used. Also, the signal-to-noise ratio (Fattahi, Hajipour, & Nobari, 2015), as stated in Eq. (32), was where N and Y indicate the number of orthogonal arrays and the response, respectively. Table 5 gives a general view of the level of parameters for the Taguchi method. Orthogonal arrays are used in the Taguchi method with the aim of studying all factors concurrently. The L9 design is used for all CSSAOS1, CSSAOS2, …, and CSSAOS10 algorithms. As illustrated in Fig. 6, the best Table 7 Aggregate results of the Wilcoxon Signed-Rank Test.  According to the results of the Taguchi method, the best level for MaxIt is the median level (100 repetitions) and the best level for the number of agents is the lower level, i.e. 30.

Results
The standard deviation, elapsed time, and mean value of the best solutions of 30 independent runs were employed as the performance metrics. Specifically, standard deviation indicates the stability of CSSAOS during all the runs, and the mean value indicates the expected optimal value between all the independent runs. Elapsed time refers to    the average total time (in seconds) that each algorithm needs to run in order to determine the computation cost of each algorithm. For providing a fair comparison, the main controlling parameters of these algorithms, i.e. maximum iteration and the number of search agents, were set equal to 100 and 30, respectively. The average number of patients admitted to the emergency department of BHH between August 2018 and August 2019 was 6996, of which 10 % were categorized into ESI1 (3 % for route 11, 4 % for route 12, and 3 % for route 13) and another 15 % into ESI2 (2 % for route 21, 5 % for route 22, and 8 % for route 23). For the others, 20 % went to ESI3 (15 % for route 31 and 5 % for route 32), 25 % to ESI4, and 30 % to ESI5.
To examine the applicability of the method, nine different test problems grouped as small, medium, and large-size problems were considered. Each test problem exemplifies a different level of complexity, different planning horizon, different number of patients, and different percentages of patients with different ESIs. The number of beds and electrical equipment and their electricity consumption were also considered, which is provided in Table 10 and Table 11 (in  Appendix).
The comparison results in Table 6 verify that CSSAOS performs slightly better than the other well-known metaheuristic algorithms. The * and ** symbols represent the smallest and the second smallest value in each column, respectively.
The Wilcoxon Signed-Rank test is used in order to perform a pairwise comparison and the comparison are made between CSSAOS1, CSSAOS2, …, and CSSAOS10 algorithms using rank-sum values. The result of the pairwise comparison is given in Table 7. In summary, from the results of 9 test problems, it could be easily concluded that CSSAOS8 is very efficient and produces very competitive results based on quality and reliability criteria. In fact, CSSAOS8 strongly competes with the current state-of-the-art algorithms.
The convergence curves of different CSSAOS algorithms on different test problems TP1,…, TP9 (small, medium, and big sizes) are shown in Fig. 7. The behavior of the average fitness of all agents demonstrates the capability of CSSAOS in achieving very good results during the solution process.
The hypothesis test is performed to show whether the mean values of the superior algorithms are significantly different. The Wilcoxon Signed-Rank Test (James & Li, 2015) reveals the performance of the proposed algorithm with the other well-known algorithms. To this aim, the null hypothesis (H 0 ) and the alternative hypothesis (H 1 ) can be used to determine the significant level of rejecting the null hypothesis, which is 0.01: { H 0 : μ 1 = μ 2 H 1 : μ 1 ∕ = μ 2 According to Table 8, H 0 is rejected at the 99 % significance level, while the acceptance of H 1 implies that the obtained optimal values of our proposed algorithm are distinctive from those of the other algorithms.

Sensitivity analysis
Considering the growing trend of cardiovascular diseases worldwide (such as arrhythmias, aorta disease, congenital heart disease, coronary artery, heart attack, heart failure, and cardiomyopathy), it is pertinent to study the effect of the increasing number of high-risk patients (ESI1, ESI2) on the amount of CO 2 emissions in the hospital. For this, two scenarios are defined. In the first scenario, the percentage of patients in ESI1, ESI2, ESI3, ESI4, and ESI5 are 12 %, 17 %, 18 %, 23 %, and 30 %, respectively. In the second scenario, the percentages are 15 %, 20 %, 15 %, 20 %, and 30 %, accordingly. The results of 30 independent runs are reported in Table 9.
As can be seen in Table 9, for each test problem and each scenario, there is an algorithm that performs superior to the other algorithms. Considering the increase in high-risk patients due to the COVID-19 pandemic, two scenarios have been designed. In the second scenario, the percentage of high-risk patients is higher than in the first scenario. According to Fig. 8, as the number of patients in a serious ESI category increases, the amount of CO 2 emitted grows exponentially. Although, the rate of growth in scenario 2 is sharper than in scenario 1 and the normal situation.

Managerial insights
This study provided the hospital with managerial insights for better patient flow and the environmental impact of the treatment process. The issue with the current operating conditions and the flow of the patient system we identified during the course of the study led to a meeting with hospital authorities to highlight the hospital's responsibility in complying with low-carbon healthcare. Our optimization results confirmed that patient waiting times would be reduced if the proposed  algorithm was used. This would assist the hospital to achieve a better service level with less waiting time and length of stay. Furthermore, hospital administrators agreed that the research results could be applied, as they addressed both operational and environmental concerns. They highlighted the importance of implementing an integrated web-based scheduling system across the hospital to help the hospital offer timely services to patients while respecting its environmental commitment.

Conclusion
The complexity and cost of healthcare systems, especially in hospitals, are becoming increasingly concerning due to the emergence of complex equipment and technologies. In addition, numerous treatments consume large amounts of electricity and, therefore, indirectly increase CO 2 emissions (MacNeill et al., 2017). Therefore, this research proposes an NP-hard carbon-efficient flexible job shop scheduling problem and a metaheuristic optimization algorithm called CSSAOS. The algorithm integrates SSA with chaotic maps to update the position of followers, the sine cosine algorithm to update the leader position, and oppositionbased learning for a better exploration of the search space, generating more accurate solutions. Based on a real-world case study, the results indicate that the proposed CSSAOS algorithm exhibited better performance to solve problems with a complex search space compared to DE, GA, GOA, OSSA, QSSA, SSA, and WOA.
We have applied the CSSAOS algorithm to solve a green patient flow problem. This research could be extended in two directions. From the problem viewpoint, more details, such as resource management and staff rostering, could be added to the patient flow. Also, energy management and policies, such as replacing equipment with different technology, could be considered. From a solution viewpoint, one can compare the proposed algorithm and its variants with other multi-objective metaheuristics such as NSGA-II and MOPSO. One could also consider extending the proposed CSSAOS algorithm to a multi-objective algorithm.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Data availability
Data will be made available on request.  Fig. 10. Patient flow routes based on his/her ESI.