simultaneously optimizing the effect of epidemic control and effect of common disease treatment. There are three main challenges to solving this problem. The first is to evaluate the effects of epidemic control and common disease treatment in a relatively accurate manner. The second is to meet the budget constraint, which is often strict in the pandemic. The third is to approximate the Pareto front of the problem in an efficient manner. For the first challenge, we develop a procedure to simulate the arrival and treatment of cases of infection and cases of common diseases according to a general principle of disease treatment and medical supplies usage. However, this also makes the objective functions expensive. Moreover, a major hospital often involves tens of thousands of medical supplies, which makes the dimension of the solution space too high. The combination of these reasons makes the problem very difficult to solve.
Decomposition is a general approach to solving large complex problems that are beyond the reach of standard techniques. Decomposition in optimization appears in early work on large-scale linear programming problems from the 1960s [1] . Many problems with separable objective functions are trivial to solve by mathematical methods. For additively decomposed functions that are not able to optimize by standard genetic algorithms, Mühlenbein and Mahnig [2] proposed the factorized distribution algorithm that factors the distribution into conditional and marginal distributions based on function structures. Strasser et al. [3] proposed factored evolutionary algorithms, which factors an optimization problem by creating overlapping subpopulations that optimize over subsets of variables. Many combinatorial optimization problems can be solved by using efficient methods to solve subproblems and combining the results to obtain solutions to the original problems. Zheng and Xue [4] utilized these characteristics to automatically derive efficient problem-solving algorithms, including evolutionary algorithms (EAs) that are mainly used for NP-hard problems. Unfortunately, the medical supplies procurement problem considered in this paper does not satisfy the basic conditions of decomposition, because it is quite common that one supply can be used in multiple diseases, and the treatment of a disease can involve many supplies.
For complex problems whose subcomponents interact and affect each other, Potter and De Jong [5] proposed a cooperative coevolution architecture that decomposes a problem into subcomponents, which are then evolved as a collection of cooperating species. Different decomposition methods have different performance in cooperative coevolution for different problems. Yang et al. [6] proposed a cooperative coevolution framework that uses random grouping and adaptive weighting in problem decomposition and coevolution for optimizing large nonseparable problems. Chandra [7] presented a competitive cooperative coevolution method for training recurrent neural networks, where different decomposition methods compete with different features they have in terms of diversity and degree of non-separability. To capture the interdependency among variables, Hu et al. [8] proposed a fast interdependency identification algorithm, which avoids obtaining the full interdependency information of nonseparable variables for cooperative coevolution. Yang et al. [9] proposed a data-driven approach, which exploits historical data to mine the evolution consistency among variables for dynamic variable grouping. Gomes et al. [10] extended cooperative coevolution with operators that put the number of coevolving subpopulations under evolutionary control. Wang et al. [11] proposed a formulabased grouping strategy for grouping variables directly based on the separable and nonseparable operations in the formula of an objective function. Mahdavi et al. [12] proposed an incremental cooperative coevolution algorithm in which subcomponents are dynamically added to an integrated subcomponent being evolved. For multimodal optimization problems, Peng et al. [13] proposed a method that concurrently searches for multiple optima as informative representatives to be exchanged among subcomponents for compensation in coevolution. For dynamic optimization problems, Peng et al. [14] used multiple populations in cooperative coevolution to compensate information in dynamic landscapes. Zhang et al. [15] proposed a dynamic cooperative coevolution framework, which allocates computational resources to elitist subcomponents with superior variables. Although cooperative coevolution has shown its efficiency in solving many engineering optimization problems [16] [17] [18] [19] [20] , for most nonseparable problems, decomposition causes the loss of a great deal of information, and the algorithms easily gravitate towards sub-optima represented by Nash equilibria rather than global optima [21] .
A multiobjective optimization problem is much more difficult than its single-objective counterpart, and decomposition is also a basic strategy in multiobjective optimization. Zhang and Li [22] proposed a multiobjective evolutionary algorithms based on decomposition (MOEA/D), which decomposes a multiobjective optimization problem into a set of singleobjective optimization subproblems using decomposition approaches such as weighted sum, weighted Tchebycheff, and penalty-based boundary interaction. However, for some problems, these decomposition approaches may not be suitable for balancing the diversity and convergence. Wang et al. [23] revolved this difficulty by imposing a constraint to an unconstrained subproblem, where the improvement region of each subproblem is determined by an adaptive control parameter. MOEA/D makes an assumption that two neighboring subproblems have similar optimal solutions, but some problems do not satisfy this assumption. Mei et al. [24] proposed a decomposition-based memetic algorithm with neighborhood search for multiobjective capacitated arc routing problem, which combines decomposition-based and domination-based techniques for solution selection. An EA proposed by Cai et al. [25] also combined domination-based sorting and decomposition, the former for evolving an internal population and the latter for maintaining an external archive. Jan and Zhang [26] introduced a penalty function to MOEA/D for multiobjective constrained optimization. Konstantinidis and Yang [27] adapted MOEA/D to solve a K-connected deployment and power assignment problem by introducing a problemspecific repair heuristic for infeasible solutions. Zhang et al. [28] extended MOEA/D for big optimization problems by embedding a gradient-based local search. Chen et al. [29] extended MOEA/D by assigning each subproblem with an upper bound vector based on -constraint for constrained optimization. Qiao et al. equipped MOEA/D with an anglebased adaptive penalty scheme. Fang et al. [30] proposed coordinate transformation in the objective space to accelerate the convergence process of multiobjective EA. Jiang et al. [31] studied the effect of scalarizing functions and presented two new functions for improving decomposition in MOEA/D. There have been many other extensions and applications of MOEA/D in recent years [32] . The cooperative coevolution architecture that decomposes a problem in the decision space has also been extended for multiobjective optimization [33] [34] [35] [36] . To our knowledge, there is only one study, by He et el. [37] , combining decomposition in both the decision space and objective space, which is similar to our work. The key difference is that our method utilizes problem-specific knowledge to ensure the separability of the transformed problem. Unfortunately, we found that, although using decompositionbased strategies in MOEA/D and other similar algorithms [38] [39] [40] [41] [42] can reduce the complexity to a certain degree, the performance of those algorithms is still far from satisfactory in solving the medical supplies procurement problem in practice.
In this study, we present a transform-and-divide approach to efficiently solve the problem. First, we transform the original problem of determining the purchase quantity of each supply to a new problem of distributing the budget to epidemic control and all common diseases. In our case studies, the dimension of the transformed problem is only one to two percent of that of the original problem. However, evaluating each solution to the transformed problem is itself a nontrivial optimization problem. Second, we divide the evaluation problem into a set of low-dimensional, single-objective optimization subproblems. We propose a hybrid evolutionary optimization approach, which employs a multiobjective EA to evolve a population of main solutions to the transformed problem and uses a tabu search algorithm to solve the subproblems. During the peak of COVID-19, we applied the proposed approach to six hospitals in Zhejiang Province, China. Results demonstrated that the transform-and-divide evolutionary optimization approach exhibits significantly better performance than that of directly using multiobjective EAs to solve the original problem. The main contributions of this paper are twofold:
• We propose a transform-and-divide evolutionary optimization approach to medical supplies procurement and demonstrate its practicability and efficiency during COVID-19.
• We show that using problem-specific knowledge to transform and divide a complex optimization problem can lead to competitive EAs for the problem. This approach can be extended to many other problems and enlarge the application field of EAs. The remainder of this paper is organized as follows. Section II presents the medical supplies procurement problem. Section III simply describes how to directly use basic multiobjective EAs to solve the original problem. Section IV proposes the transform-and-divide evolutionary optimization approach. section V presents the computational results. Section VI concludes with a discussion.
We consider a medical supplies procurement problem formulated as follows. In a pandemic, a hospital plans to procure medical supplies, including a set S = {S 1 , S 2 , . . . , S n } of n supplies for epidemic control, and a set S = {S n+1 , S n+2 , . . . , S n+n } of n supplies for normal disease treatment. For each supply S k , the current inventory is a k , the unit price is c k , and the unit volume is v k . The problem is to determine the purchase quantity x k of each supply (1 ≤ k ≤ n + n ), such that the effects of epidemic control and normal disease treatment are simultaneously optimized.
The supplies for epidemic control can be divided into two classes. The first class consists of supplies such as latex gloves and normal saline that must be used in the treatment of a suspected case of infection; we use Ψ 0 to denote the set of these supplies, and use q 0,k to denote the quantity of each S k ∈ Ψ 0 required to treat a case. The second class consists of supplies that are alternative in some treatment items. Table I presents six treatment items and their alternative supplies used for COVID-19 control in this study. The six sets of supplies are denoted by Ψ 1 , Ψ 2 , . . . , Ψ 6 , respectively, and the quantity of each alternative S k ∈ Ψ j required to treat a case is denoted by q j,k . Different alternatives have different treatment effects. The treatment effect of using each alternative S k ∈ Ψ j is estimated as e j,k . For example, the effects of peroxide, impermeable gown, and normal gown in "body protection" item are estimated as 1, 0.9 and 0.7, respectively. If we choose the S kj ∈ Ψ j for the j-th treatment item (1 ≤ j ≤ 6), the corresponding epidemic control effect on the case is empirically estimated as: E(k 1 , ...k 6 ) = (0.4e 1,k1 +0.6e 2,k2 )e 3,k3 (0.2e 4,k4 +0.8e 5,k5 )e 6,k6 (1) The hospital is capable of treating a set D = {D 1 , D 2 , . . . , D m } of m diseases. Similarly, for each disease D i , the set of supplies that must be used is denoted by Ψ i,0 , and the set of supplies that are alternative in J i treatment items are denoted by Ψ i,1 , Ψ i,2 , . . . , Ψ i,Ji , respectively. The quantity of each S k ∈ Ψ i,0 required to treat a case is q i,0,k , and the quantity of each alternative S k ∈ Ψ i,j required to treat a case is q i,j,k . Different alternatives have different treatment effects. The treatment effect of using S k ∈ Ψ i,j is estimated as e i,j,k . If we choose the S kj ∈ Ψ i,j for the j-th treatment item (1 ≤ j ≤ J i ), the corresponding treatment effect on the case is empirically estimated by a therapeutic effect function E i (k 1 , k 2 , ..., k Ji ). Like Eq. (1), the typical expression of E i is a weighted sum or product of e i,j,k [43] .
Let T be the procurement decision cycle. In our case studies, the hospital procures medical supplies every 15 days. The supply quantities are determined based on the estimation of the number of hospital visits in the next decision cycle. For the number of cases of each disease D i , we estimate three values: the expected value r i , lower limit (optimistic value) r i , and upper limit (pessimistic value) r i . The values can be obtained based on historical morbidity data and environmental influence factors [44] [45] [46] [47] .
The number of suspected cases of epidemic infection is estimated based on the number of hospital visits of different diseases. For each disease D i , we estimate a probability p i that a patient of D i is a suspected case of COVID-19. In general, a disease having more similar symptoms with the epidemic has a higher p i . For example, an acute respiratory infectious disease has a high p i . For a disease (such as fracture) that is unrelated to the epidemic, we set p i to the current incidence p e of infection (including suspected infection) in the local region. We also estimate an average number r i of accompanying persons of a patient of D i ; in general, a critical disease has a large r i . The probability that an accompanying person of a patient of D i is a suspected case of COVID-19 is p i , which is set to 0.5p i if D i has similar symptoms with the epidemic and p e otherwise. The total number of suspected cases of infection in the next decision cycle is estimated as follows (we use r i as we take a serious or pessimistic view of epidemic control):
A solution to the medical supplies procurement problem can be represented by a (n + n )-dimensional vector x = {x 1 , . . . , x n , x n+1 , . . . , x n+n }. The fitness of x is evaluated by two objective functions: (1) the epidemic control effect Υ(x), which is the sum of treatment effects of all suspected cases of infection; (2) the common disease treatment effect Υ (x), which is the weighted sum of treatment effects of all common cases, where the weight of each D i is w i . It is assumed that the arrival time of cases follows a uniform distribution. That is, for each disease D i , as the expected number of cases in a decision cycle of 15 days is r i , then there is a case arriving every (15× h w )/r i hours, while h w is the daily working hours (24 for emergency diseases and 8 for non-emergency diseases in our study). Moreover, there is one suspected case of infection in every 1/(p i + p i r i ) cases of D i .
A general principle of disease treatment is "focusing on the current patient", i.e., whenever a new case arrives, the physician always chooses the most effective supply from the available alternatives, as he does not know how many 
in nonincreasing order of e j,k or e i,j,k , and simulate the arrival and treatment of all cases according to the procedure shown in Algorithm 1 to calculate the values of Υ(x) and Υ (x). In Algorithm 1, the boolean variable tr denotes whether the remaining supplies are capable of treating a suspected case of epidemic infection, and tr i denotes whether the remaining supplies are capable of treating a case of D i .
Algorithm 1: Procedure for evaluating the effects of epidemic control and disease treatment for the original problem. if Ψi,j = ∅ then tri ← false;
if the case is a suspected infected case and tr = true then 
A procurement solution x must satisfy problem constraints. First, the total procurement cost cannot exceed the budget C:
The hospital should perform its normal functions. In this study, it is required that the hospital is able to treat r i (the lower number) cases of each disease D i . These constraints can be tested by simulating the arrival and treatment of the lower numbers of cases in Algorithm 1: if any new case cannot be treated, i.e., whenever the condition tr i = false (Line 6 of Algorithm 1) is triggered while the number of simulated cases of D i is less than r i , the constraint is violated.
It is also required that the hospital is able to treat r 0 suspected cases of infection. Whenever the condition tr = false (Line 18) in Algorithm 1 is triggered, the constraint is violated.
For the above (n+n )-dimensional, constrained bi-objective optimization problem, we can use evolutionary constrained multiobjective algorithms to search for the Pareto optimal solutions. The search range of each dimension k is [x k , x k ]. The lower limit x k is set to the total quantity of S k required in non-alternative treatment items for all cases:
n+1 ≤ k ≤ n+n (4) The upper limit x k can be set to the total required quantity of S k under the assumption that S k is always chosen whenever S k is an alternative. That is, if 1 ≤ k ≤ n, we set
Otherwise, we set
We adopt the following five well-known evolutionary constrained multiobjective algorithms to solve the problem:
• The nondominated sorting genetic algorithm II (NSGA-II) with the constrained-domination principle [48] . • The constrained multiobjective evolutionary algorithm (CMOEA) based on an adaptive penalty function and a distance measure [49] . • The multiobjective evolutionary algorithm based on decomposition (MOEA/D) [22] with a penalty function for constrain handling [26] .
• The differential evolution with self-adaptation and local search for constrained multiobjective optimization (DECMOSA) [50] , which combines constraineddomination and penalty function for constrain handling. • The multi-objective particle swarm optimizer based on dominance with decomposition (D 2 MOPSO) [38] . CMOEA, MOEA/D, and DECMOSA employ penalty functions for constrain handling. Violation of constraint (3) is calculated as max(0, n+n k=1 c k x k − C). For constraints that all suspected cases of infection and the lower number of cases of each common disease must be treated, we set the violation of each constraint equal to the budget C, i.e., the violation is C times the number of false tr i and tr in Algorithm 1.
Nevertheless, the performance of all the above algorithms is not satisfying, mainly because the dimension (n+n ) is very high (approximately 10,000∼40,000 in a major hospital in our case studies) and the evaluation of a solution using Algorithm 1 is computationally expensive.
In this section, we propose a new transform-and-divide approach to efficiently solve the problem. First, we transform the original high-dimensional, constrained bi-objective optimization problem to a low-dimensional, unconstrained biobjective optimization problem, which can be solved using evolutionary (unconstrained) multiobjective algorithms. The evaluation of each solution to the transformed problem can be divided into a set of low-dimensional, single-objective optimization subproblems, which can be solved using a tabu search algorithm.
We transform the original problem of determining the purchase quantity of each supply to a problem of determining the purchase budget for epidemic control and the purchase budget for each disease. First of all, we calculate the cost for purchasing the supplies that must be used in the non-alternative treatment items and, therefore, obtain the remaining budget as:
Consequently, the transformed problem is to distribute C to m+1 components, denoted by {y 0 , y 1 , . . . , y m }, where y 0 is the budget for purchasing alternative supplies for epidemic control, and y i is the budget for purchasing alternative supplies for treating disease D i (1 ≤ i ≤ m). The dimension of the transformed problem is m (approximately 200∼600 in a major hospital in our case studies), which is significantly smaller than the dimension n+n of the original problem.
Moreover, the search range of each dimension of the transformed problem is also much smaller than that of the original problem. For epidemic control, the lower limit y 0 of budget y 0 can be obtained using the following steps: 1) Use supplies in storage to treat as many suspected cases of infection as possible;
2) If there is no remaining case, set y 0 = 0; 3) Else, for each remaining case, always purchase the cheapest supply among the alternatives, and set y 0 to the total purchase cost. And the upper limit y 0 of budget y 0 can be obtained using the following steps:
1) Treat each suspected case in the most effective way, i.e, always select the supply with the maximum treatment effect e j,k among the alternatives, and calculate the total required quantity of each supply; 2) Calculate the purchase quantity of each supply, and set y 0 to the total purchase cost. Therefore, the search range of y 0 is limited to [y 0 , y 0 ]. We can obtain the search range [y i , y i ] of each y i in a similar manner. In our case studies, the average value of (y i −y i ) is approximately 95 (in unit of 100 RMB), while that of (x k −x k ) is approximately 1100 (in minimum order quantity).
The solution space of the transformed problem is significantly smaller than that of the original problem. But how to evaluate a solution y = {y 0 , y 1 , . . . , y m } to the transformed problem? The task can be divided into m+1 optimization subproblems. The first subproblem is to determine the purchase quantities under the budget y 0 so as to maximize the epidemic control effect. Each of the remaining m subproblems is to determine the purchase quantities under the budget y i so as to maximize the treatment effect of disease D i (1 ≤ i ≤ m).
However, the division leads to difficulty in allocating supplies in storage to different diseases. We overcome this difficulty by employing a procedure similar to Algorithm 1 to simulate the arrival and treatment of all cases. But the procedure has two differences from Algorithm 1:
• Initially, we only consider supplies in storage, i.e., Line 3 of Algorithm 1 is not executed. • If there is no supply in storage that can be used for a treatment item (i.e., the condition in Line 14 or Line 26 is satisfied), we temporarily purchase "in advance" the cheapest alternative supply for the item. The procedure also produces the "cheapest" solution to each subproblem, which can be evolved to an optimal or nearoptimal solution, as described in the next subsection.
The proposed method employs an evolutionary multiobjective algorithm to evolve a population of main solutions to the transformed problem, and employs a tabu search algorithm to solve the subproblems for evaluating each main solution.
For the first subproblem, each solution z can be represented by six vectors as follows (the vector lengths do not need to be the same): where z j,k denotes the number of cases that use the k-th alternative supply for the j-th treatment item, and each vector satisfies |Ψj | k=1 z j,k = r 0 . The procedure described in Sec. IV-B produces the cheapest solution to the subproblem, denoted by z † . First, we continually use the following steps to improve z † by replacing an alternative supply to a more effective alternative for a randomly selected case until z † cannot be further improved: 1) Randomly selecting two components z j,k and z j,k in a vector satisfying z j,k > 0; 2) Set z j,k = z j,k −1 and z j,k = z j,k +1 if doing so would not violate the budget constraint.
Starting from the improved z † , the tabu search algorithm continually uses the following steps to search around and improve z † until the stopping condition is satisfied: 1) Generate k N neighboring solutions of the current z † , each being obtained by randomly selecting two components z j,k and z j ,k satisfying k < |Ψ j |, k < |Ψ j |, z j,k > 0, and z j ,k +1 > 0, and setting z j,k = z j,k −1, z j,k+1 = z j,k+1 +1, z j ,k = z j ,k +1, and z j ,k +1 = z j ,k +1 −1, if doing so would not violate the budget constraint; 2) Select the best neighbor that is not tabued or is better than the current z † , make z † move to this neighbor, and add this move to the tabu list.
The remaining m subproblems can be solved by tabu search in a similar way. As demonstrated by the experiments, the tabu search algorithm can quickly obtain optimal solutions for most subproblem instances, given that the dimensions of the subproblems are relatively small. For example, as we can observe from Table I , the dimension of the first subproblem is 18 (note that the last dimension of each vector can be determined by other dimensions of the vector, and the actual dimension in the solution space is only 12). Therefore, the tabu search algorithm is very suitable for the subproblems, as it will be invoked many times to evaluate main solutions.
For the main transformed problem, we adopt the following evolutionary multiobjective algorithms to evolve main solutions and invoke the tabu search algorithm:
• NSGA-II [48] . • MOEA/D [22] . • A differential evolution for multiobjective optimization with self-adaptation (DEMOwSA) [51] . • A multiobjective particle swarm optimization (MOPSO) algorithm [52] which extends comprehensive learning [53] for multiobjective optimization.
In this subsection, we theoretically compare the complexities of the original problem and the transformed problem. For notational simplicity, we regard suspected infection as a disease with subscript i = 0. The number of all possible solutions to the original problem is
And the algorithm 1 for evaluating each solution to the original problem has a time complexity
After transformation, the number of all possible solutions to the new problem is
To evaluate each solution to the transformed problem, we need to solve m + 1 subproblems. The number of possible solutions to the i-th subproblem is
where z j,k and z j,k denote the upper and lower limits of decision variable z j,k . And the time complexity of evaluating each solution to the i-th subproblem is
Consequently, the total time complexity of the transformed problem is N m i=0 N i O(f i ) , while that of the original problem is N · O(f ), and the complexity reduction ratio of transformation is 
In our case studies, the average values are x ≈ 1100, y ≈ 95, z ≈ 66, |Φ| ≈ 37, n + n ≈ 25000, and m ≈ 450. Therefore, the average complexity reduction ratio on the instances is approximately 58.
We use the proposed method for medical supplies procurement in Zhejiang Hospital of Traditional Chinese Medicine (ZJHTCM) from 15 Feb to 15 Apr, 2020, the peak of COVID-19 in Zhejiang Province, China. Since 15 Mar, we also extend the application to other five hospitals (denoted by H1-H5) .
Therefore, there are 14 real-world instances of the medical supplies procurement problem. Table II summarizes the main characteristics of the instances, where i r i denotes the total expected number of cases of all common diseases, J denotes the average treatment items per disease, |Φ| denotes the average number of alternatives per treatment item, and the budget C is in RMB. The instances are solved on a workstation with an i7-6500 2.5GH CPU, 8GB DDR4 RAM, and an NVIDIA Quadro M500M card.
Before testing the algorithms for solving the main problem, we first test the performance of the tabu search algorithm for subproblems. From the above real-world main problem instances, we select 16 subproblem instances, the dimensions D of which range from 12 to 72. For the algorithm, we set the neighborhood size k N to 2D, tabu length to 12, and the maximum number of iterations to 50D. On each instance, we run the algorithm 50 times to test whether and how long it can obtain the exact optimal solution (validated by an exact branch-and-bound algorithm [54] ). Fig. 1 presents the convergency curves (averaged over the 50 runs) of the tabu search algorithm on the subproblem instances. The algorithm reaches the optima within 100 iterations (10 ms in our computing environment) when the problem dimension is smaller than 24, within 200 iterations (30 ms) when the dimension is smaller than 40, and within 400 iterations (120 ms) on all instances. In our case studies, the average dimension of the instances is approximately 37, which can be solved using approximately 160 iterations (25 ms); the dimension of the largest instance is 72, which can be solved using 369 iterations (116 ms). Using multithreading and GPU acceleration, the average CPU time for evaluating a main solution to a problem of 400 diseases is approximately 600 ms.
For each main problem instance, we use five evolutionary constrained multiobjective algorithms, including NSGA-II with constraint handling (denoted by NSGA-II-C) [48] , CMOEA [49] , MOEA/D with constraint handling (denoted by MOEA/D-C) [26] , DECMOSA [50] , and D 2 MOPSO [38] , to solve the original problem; we also use four evolutionary multiobjective algorithms, including NSGA-II [48] , MOEA/D [22] , DEMOwSA [51] , and MOPSO [52] , all combined with tabu search, to solve the transformed problem. The control parameters of all algorithms are tuned on the whole set of instances. For a fair comparison, all the algorithms use the same stopping criterion that the CPU time does not exceed 90 minutes, which is also applied in our practice. On each instance, each algorithm is run 30 times. Fig. 2 compares the hyperarea (the area under the Paretoapproximated front in objective space, also known as the hypervolume) [55] , [56] obtained by each algorithm on each main problem instance. It is clear that the last four algorithms using transform-and-divide exhibit significant performance advantages over the first five algorithms. On large-size instances of ZJHTCM and H5, the median hyperareas of the four transform-and-divide EAs are approximately eight to ten times of those of the basic EAs; on the other instances, the median hyperareas of the transform-and-divide EAs are approximately six to eight times of those of the basic EAs. In general, the performance of the transform-and-divide EAs is mainly affected by m (the number of diseases) and C (the total budget), while that of the basic EAs is mainly affected by n+n (the number of supplies). This is why all algorithms obtain relatively high hyperareas on the instances of H1 and H2, where m, n+n , and C are relatively small. Nevertheless, the performance advantages of the transform-and-divide EAs over the basic EAs are very significant on these relatively smallsize instances. On the last two instances of H5, the values of these parameters are large, and the performance advantages of the transform-and-divide EAs over the basic EAs are not so significant. On each instance, the minimum hyperareas of the transform-and-divide EAs are still significantly larger (about four to five times) than the maximum hyperareas of the basic EAs. The results demonstrate that the proposed transform-anddivide method can greatly reduce the difficulty of solving the complex original problem.
We also make a pairwise comparison between each transform-and-divide EA and its constrained version for the original problem, e.g., NSGA-II vs. NSGA-II-C, NSGA-II vs. CMOEA, MOEA/D vs. MOEA/D-C, DEMOwSA vs DECMOSA, and MOPSO vs D 2 MOPSO. That is, in each pair, the first is an unconstrained multiobjective EA used in transform-and-divide, and the second is a multiobjective EA adding constraint handling mechanisms (such as constraineddomination or penalty functions) to the first one. The comparison is made based on the ratio of the hyperarea obtained by the first algorithm to the maximum hyperarea obtained by the second algorithm. Fig. 3 presents the changes of hyperarea ratios over the CPU time on each problem instance. In general, it takes about five minutes for a transform-anddivide EA to reach the maximum hyperarea obtained by its counterpart after 90 minutes; on large-size instances of ZJHTCM and H5, the ratio is between 3 and 4 at ten minutes and exceeds 5 before twenty minutes; on the other instances, the ratio is approximately 2 at ten minutes and between 3 and 4 at twenty minutes. The high hyperarea ratios demonstrate that the transform-and-divide EAs approximate the Parato front of the problem significantly more efficiently than the basic EAs, mainly because the transform-and-divide strategy significantly reduces the solution space, and also because the unconstrained EAs for the transformed problem can explore the solution space more effectively than their constrained versions. Comparatively, the ratios of NSGA-II vs. its counterparts (NSGA-II-C using constrained-domination and CMOEA using a penalty) grow more slowly than those of DEMOwSA and MOPSO vs. their counterparts in early stages of evolution. As there is no significant difference among the five basic EAs for the original problem, the results indicate that DEMOwSA and MOPSO are more efficient in exploring the solution space to approximate the Parato front of the transformed problem.
Besides the hyperarea metric, we also compare the al- divide EAs, while solutions obtained by the basic EAs can hardly provide reference. In practice, decision-makers choose final solutions for implementation as follows:
• The best solutions of MOEA/D on the fourth instance of ZJHTCM and the second instance of H4. • The best solutions of DEMOwSA on the second instance of ZJHTCM, the second instance of H1, the second instance of H2, the second instance of H3, and the second instance of H5. • The best solutions of MOPSO on the remaining seven instances.
This paper presents a transform-and-divide evolutionary optimization approach to medical supplies procurement under the background of COVID-19. Our approach first transforms the original high-dimensional, constrained multiobjective optimization problem to a low-dimensional, unconstrained multiobjective optimization problem, and then evaluates each solution to the transformed problem by solving a set of simple single-objective optimization subproblems, such that the problem can be efficiently solved by existing evolutionary multiobjective algorithms. We applied the transform-anddivide evolutionary optimization approach to six hospitals in Zhejiang Province, China, during the peak of COVID-19. Results showed that our approach exhibits significantly better performance than that of directly solving the original problem. Decision-makers of the hospitals always choose the best solutions produced by the transform-and-divide method for implementation and achieve promising results in balancing epidemic control and common disease treatment in practice.
The proposed transform-and-divide evolutionary optimization based on problem-specific knowledge can be an efficient solution approach to many other complex problems. For example, considering a problem of personalizing healthcare solutions for a large number of residents. As the number of candidate healthcare medicines and treatment items are large, the dimension of the problem is high. However, by clustering the residents based on their health status and limit the medicines and treatment items to each cluster, the problem dimension can be significantly reduced, and we can solve the subproblem of personalizing healthcare solutions for each cluster much more efficiently. Another example is to distribute m types of disaster relief supplies from n 1 suppliers to n 2 demanders, we need to determine the quantity of each supply from each supplier to each demander. The dimension of the problem is mn 1 n 2 . However, we can establish a virtual "intermediary" and transform the problem to a new problem of determining the quantity of each supply from each supplier to the intermediary and that from the intermediary to each demander. Therefore, the dimension of the transformed problem is m(n 1 +n 2 ), but we need to solve additional subproblems of determining at the intermediary which parts of supplies should be sent to each demander. In many cases, the transform-anddivide strategy can greatly reduce the difficulty of problemsolving, but it often requires effective discovery and utilization of problem-specific knowledge for problem transformation and division. Although EAs are regarded as robust problemindependent search heuristics for a large variety of optimization problems [57] , [58] , we argue that proper exploitation of problem-dependent knowledge can significantly improve the efficiency of EAs in solving highly complex problems and, therefore, enlarge the application field of EAs.