An Efficient Imperialist Competitive Algorithm for Solving the QFD Decision Problem

It is an important QFD decision problem to determine the engineering characteristics and their corresponding actual fulfillment levels. With the increasing complexity of actual engineering problems, the corresponding QFD matrixes become much huger, and the time spent on analyzing thesematrixes andmaking decisions will be unacceptable. In this paper, a solution for efficiently solving the QFD decision problem is proposed. The QFD decision problem is reformulated as a mixed integer nonlinear programming (MINLP) model, which aims to maximize overall customer satisfaction with the consideration of the enterprises’ capability, cost, and resource constraints. And then an improved algorithm G-ICA, a combination of Imperialist Competitive Algorithm (ICA) and genetic algorithm (GA), is proposed to tackle this model. The G-ICA is compared with other mature algorithms by solving 7 numerical MINLP problems and 4 adapted QFD decision problems with different scales. The results verify a satisfied global optimization performance and time performance of the G-ICA. Meanwhile, the proposed algorithm’s better capabilities to guarantee decision-making accuracy and efficiency are also proved.


Introduction
Quality Function Deployment (QFD), as a well-known planning methodology originated in the late 1960s in Japan by Akao and Mazur [1], focuses on translating customer requirements (CRs) into a series of relevant design requirements. During this planning process, a lot of decisions should be made. One of them is to determine the engineering characteristics (ECs) to be implemented among numerous possible selections and their corresponding actual fulfillment levels, which not only directly influences the deployment of other design requirements and the final quality of the product but also influences overall customer satisfaction and the distribution of resources.
In practice, it is a complex process to analyze and organize a lot of information related to the QFD decision problem, including multifacets of customer needs, market driven factors, price/revenue streams of enterprises, and various regulation compliance. Multiple variables and complex relationships bring many difficulties for decision makers to make decisions. Mathematical programming, an optimization technique, such as integer programming (IP), linear programming (LP), nonlinear programming (NLP), and mixed integer nonlinear programming (MINLP), has been widely used in helping people solve these difficulties. Firstly, in order to clarify those complex relationships, the problem is abstracted as a mathematical model. Then the model is calculated by an algorithm to offer a result for decision makers as a valuable reference [2][3][4][5][6][7][8][9][10][11][12][13].
However, a rapidly changing decision environment makes people not only satisfied with obtaining a decision outcome but also want to seek for a higher decision-making efficiency, in other words, a shorter time spent on making decisions. Especially in such a society with highly developed economy, technology, and culture, customer requirements tend to be more and more diverse, product systems become gradually complex, and science technologies are continuously developing. All of these may lead to a huge and complex QFD matrix. As an important step in tackling the above QFD decision problem, the selection of algorithms influences  not only the accuracy of decisions but also total decision time. And if the calculating time of the algorithm is too long, the decision efficiency and time cannot be guaranteed.
The QFD decision problem in this paper is reformulated as a MINLP model aiming to maximize overall customer satisfaction with the consideration of the enterprises' capability, cost, and resource constraints. Actually, many decision problems in various application fields, such as chemical engineering [14,15], gas/water distribution networks optimal design [16,17], and engineering design [18], can be formulated as MINLP models. MINLP is a NP-hard problem and the calculating algorithms to MINLPs have been developed since 1980s, which can be broadly divided into two main classes: deterministic algorithms and stochastic algorithms. Deterministic algorithms, such as branch and bound method, outer approximation method, and extended cutting plane method, can provide accurate solutions for problems in a continuous space but require an additional effort of analysis (such as the gradient analysis of objective function and constraints), which may considerably add to the work [19,20]. Considering the existence of integer variables in QFD decision problems and easy-operation of engineers, stochastic algorithms, which are more capable of handling NP-hard problems and need not much professional mathematical knowledge, will be focused on in this paper. A lot of stochastic algorithms have been proposed. An algorithm called M-SIMPSA suitable for the optimization of MINLP is presented by Cardoso et al. [19]. Mohan and Nguyen [21] developed an improved computational algorithm called RST2ANU, which was based on the original RST2U [22]. Yiqing et al. [23] improved the particle swarm optimization (PSO) algorithm for solving nonconvex NLP/MINLP problem with equality and/or inequality constraints. Schlüter et al. [24] introduced two novel extensions for ant colony optimization (ACO) framework. Deep et al. [25] proposed a real coded genetic algorithm named MI-LXPM.
In this paper, Imperialist Competitive Algorithm (ICA), one of the stochastic algorithms proposed by Atashpaz-Gargari and Lucas [26], which simulates the social-political process of imperialism and imperialist competition, is applied to solve MINLP models. The key features of ICA are its fast convergent rate to reach global optimum, which has been proved in dealing with various optimization problems, such as structure design [27,28], distribution network optimization [29,30], and process planning and scheduling [31]. The advantages of ICA are beneficial to improvement of decision efficiency. In this QFD decision problem, considering the existence of integer variables in this problem, an improvement by combining genetic algorithm (GA) is carried out because the original ICA is only applied to real variables. The improved ICA, called G-IGA in the following, is compared with other mature algorithms to MINLPs, including RST2U and RST2ANU mentioned above, by calculating 7 simple numerical MINLP problems and 4 different scale QFD decision problems, and the results verify the advantages of G-ICA.
The remainder of the paper is organized as follows: the QFD decision problem is defined in Section 2. The G-ICA to resolve the problem is introduced in Section 3. Some numerical experiments and case studies are shown in Section 4, and Section 5 summarizes the paper.

Problem Definition
The QFD matrix is just as in Figure 1. Before the decision problem based on QFD is defined, some assumptions should be made in advance: (i) The QFD matrix for a certain product, including CRs, possible ECs, weight vectors, relationship, and correlationship matrixes, is already established. cost for every EC), and resource (maximum total resources and required resource to implement every EC) is offered by the enterprise.
(iii) The fulfillment level is a percentage from 0 to 100%. For every EC, actual fulfillment level is only influenced by the corresponding investment cost, whose value is from 0 to 10 with no units. And the enterprise can determine the influence relationship based on long-term practical experience all by itself.
(iv) There are no natural units or maximum values for customer satisfaction. It is the result of synthesizing all the actual fulfillment levels of engineering characteristics to be implemented.
Then the mathematical model is expressed as follows: = ( ) in which the decision variables and parameters are defined as follows: CS: the overall customer satisfaction : the required resource corresponding to engineering characteristic : the maximum total budget : the maximum total resources.

Objective Function and Decision
Variables. The objective function in this model is to maximize the overall customer satisfaction CS, just as in formula (1). Whether the engineering characteristics are necessary can be viewed as 0 or 1 integer variables. If engineering characteristic is selected, is taken as 1. Otherwise, is taken as 0. Another decision variable , the actual fulfillment level for engineering characteristic , is directly influenced by , the corresponding investment cost. The detailed relationship expressions will be introduced in next section. Then these fulfillment levels of selected engineering characteristics will be synthesized according to their weights, to get the overall customer satisfaction. A traditional calculation method to get the weight of engineering characteristic is as follows: where is the initial importance weight of the engineering characteristic , is the relative importance weight of the requirements , and is the relationship between the customer requirement and the engineering characteristic , ∈ {±1, ±3, ±5}. Moreover, the correlationship between different engineering characteristics is taken into consideration as in the following: with denoting the correlation between customer requirement and , ∈ {±1, ±3, ±5}.

Relationship between Investment Cost and Actual Fulfillment Level.
In reality, in formula (7), the relationship between investment cost and actual fulfillment levels is too complex to describe quantitatively. In order to simplify the problem, the relationship was generally viewed as linear in some past research [6,32]. Here, the relationship is updated by a new expression, which can represent linear, concave, convex, or mixed relationship type. Considering the current known conditions firstly (i) the independent variable ranges from 0 to 10; (ii) the dependent variable ranges from 0 to 100%; (iii) when is taken as 0, is taken as 0. When is taken as 10, is taken as 100%.
A preliminary expression is suggested with the help of an auxiliary function as follows: Here, a tangent function is chosen as the auxiliary function . Moreover, in order to express more relationship types in the domain of definition with only one formula, a shape parameter is introduced. The final expression is as follows: The shape parameter varies with different engineering characteristics. Actually, can be understood as a special investment cost value. At this point, the fulfillment level will get the fastest growing; that is, the slope of ( ) is maximum.
When goes to +∞ or −∞, ( ) is approximately linear, just as in Figure 2(a), which is an ideal situation.

When
≥ 10, ( ) is a concave function, just as in Figure 2(b), which represents the case that the actual fulfillment level grows faster and faster with the increase of investment cost.
When ≤ 0, ( ) is a convex function, just as in Figure 2(c), which represents the case that the actual fulfillment level grows slower and slower with the increase of investment cost.
When 0 < < 10, ( ) is concave at the beginning and convex afterward, just as in Figure 2(d). The investment cost value at the inflection point is equal to . The fulfillment level will increase sharply around this point.  corresponding to the lowest and the highest levels achieved by the enterprise. The actual fulfillment level should be within this range, just as in formula (2). Formula (3) and formula (4) are the cost constraints. The investment cost of engineering characteristic should not surpass the upper bound , and the total investment cost should not exceed the total budget limit .

Capability, Cost, and Resource
Resource for every engineering characteristic is assumed to be known. The total resource should not exceed the resource limit for product development as in formula (5).

Solving Algorithm
The detailed G-ICA algorithm steps are shown in Figure 3.
Step 1. A hybrid real-binary concatenated coding method is applied for the encoding and decoding of a solution.
Step 2. A penalty function method is used as the constraints handling method.
Step 3 (empires construction). The initial countries are generated and their corresponding objective function values are calculated. According to the power value of each country, which is a normalized expression of the objective function value, the most powerful countries are chosen as imperialists and the remained other countries are chosen as colonies. After the colonies are assigned to every imperialist randomly, the empires are constructed.
Step 4 (assimilation). Assimilation occurs inside every empire. According to the basic ICA, imperialists try to assimilate their colonies and make them similar to themselves, which aims to make solutions represented by colonies tending to solutions represented by imperialists. Here, the assimilation of real variables is the same as the original ICA algorithm.
But for integer variables, a certain modification has to be conducted. A multipoint crossover operator and a multipoint mutation operator are introduced, which refers to GA.
The crossover probability is cr , whose value is generally larger than 0.6. The mutation probability is mu , which is generally taken less than 0.1.
Step 5 (revolution). In each decade, revolution operators are performed on some of the colonies in each empire. There are many ways to realize revolution, such as swap, insertion, reversion, and perturbation. In this paper, some colonies are replaced by newly generated colonies. The population percentage that revolts is shown by re .
Step 6 (imperialist update). After a series of operators inside every empire, the cost of every imperialist should be updated if there exists any colony whose cost is lower than the imperialist's. In this case, the colony will become the imperialist in the current empire and vice versa.
Step 7. Calculate the total cost of each empire.
Step 8 (empires competition). The weakest colony of all empires will be selected and assigned to some empire. The lucky empire is chosen randomly; that is to say, the strongest empire does not always gain the initiative.
Step 9 (elimination of the powerless empire). When an empire loses all of the colonies, the empire will be eliminated from the population.
Step 10 (stop criterion). The stopping criterion is considered when there is only one empire that governs all of the colonies or the decades reached the limitations.  Table 1.   Index: -the percentage of success, t-average time in seconds used by the algorithm in achieving the optimal solution in the case of the successful runs, and FEs-the average number of function evaluations in the case of successful runs.

Case Studies and Discussions
Performance of G-ICA is compared with GA, ( + )-ES, ( , )-ES, and standard simulated annealing algorithm (SSA). Another two improved algorithms for MINLPs, RST2U, and RST2ANU algorithm, mentioned in Section 1, also take part in the comparison. All the algorithms were programmed by MATLAB R2014b, and these test problems are done on a personal computer equipped with Intel Core i7-5500U CUP, 2.4 GHz processor, 4 G RAM under Windows 8.1 operating system. The results obtained are given in Table 2.
Each problem is executed 100 runs with all the seven algorithms. Each run is initiated using a different set of initial population. A run is considered a success if achieved value of the objective function is within 1% of the known optimal value (in case the optimal value of the objective is zero, a run is considered success if the achieved absolute value of the objective function is less than 0.01). For each problem, the percentage of success (obtained as the ratio of the number of successful runs to total number of runs), the average number of function evaluations (FEs) in the case of successful runs, and the average time in seconds used by the algorithm in achieving the optimal solution in the case of the successful runs are also listed.
From the results listed in Table 2, it can be seen that, for most problems, the G-ICA can get an optimal solution in a much shorter time. For Problem 1, most algorithms can get 100% success rate, but G-ICA's running time is the shortest and FEs is the smallest. For Problems 4, 5, and 7, on the premise of similar FEs, G-ICA can get a higher success rate in a much shorter time. For Problem 2, when the maximum number of function evaluations (MAXFEs) is 1000, G-ICA gets 90% success rate, which is lower than RST2U and RST2ANU. But when the MAXFEs reaches 1500, G-ICA gets 100% success rate, whose running time is still shorter than RST2U and RST2ANU. The same situation occurs in solving Problem 3.
An interesting phenomenon is found in solving Problem 6. The global optimal solution of Problem 6 given in the Appendix is ( 1 , 2 , 3 , 1 , 2 , 3 , 4 , ) = (0.2,  1.280624, 1.954483, 1, 0, 0, 1, 3.557463). And the 37 successful runs of G-ICA are all around the global optimal solution. But RST2U and RST2ANU find another local optimal solution in searching process, ( 1 , 2 , 3 , 1 , 2 , 3 , 4 , ) = (0. 2, 0.8, 1.9079, 1, 1, 0, 1, 3.5796). In the 89 successful runs of RST2ANU, there are about 70% points around the local optimal solution, and the following are around the global optimal solution. In our calculating method of percentage of success, all these points are viewed as successful. So the data in Table 2 shows that G-ICA is far behind in the percentage of success index. This may be an evidence to prove that G-ICA can guarantee global optimal.
In fact, the number of initial populations has a great influence on the performances of algorithms, which has been proposed by Engelbrecht [33]. In this paper, each run for each algorithm is initiated using a different set of initial population, and every index's average value is taken as reference, in order to compare these algorithms fairly. And for G-ICA, both the number of initial countries pop and number of imperialists imp , that is, the number of groups which the initial countries are divided into, are set differently. Here, an experiment is done for Problem 2 to research how these parameters influence the algorithm's performances. The iteration stop criteria of G-ICA is to reach the maximum FEs. The results are shown in Figures 4 and 5.
For the percentage of success, it can be seen from Figures 4(a) and 5(a) that more FEs provide higher percentage of success. But the parameters pop and imp 's influences to this index are hard to describe, because the relationship between pop ( imp ) and the percentage of success is wavy.
From the view of index FEs, the performance has become deteriorative with the increase of pop , just as in Figure 4(b), because the average FEs increase. But the average FEs will decrease when imp becomes larger, just as in Figure 5(b).
The parameters' influences to the average running time are shown clearly in the figures. On the premise of the same maximum FEs, the larger pop is, the fewer iteration times are, so the average time is shorter, just as in Figure 4(c). But the larger imp is, the more times of empires competitions are, so the average time is longer, just as Figure 5(c). In conclusion, if you want to pursue a shorter running time, you may set a larger pop and a smaller imp .

Case of the Washing Machine and Discussions.
The data of this case is adapted from a washing machine case [13,32,34]. The QFD of a washing machine is as in Figure 6.
The importance weight matrix of engineering characteristics is as follows: The relationship matrix between requirements and engineering characteristics is as follows: ] . ( The correlationship matrix between engineering characteristics is as follows: ] . The initial importance weight of engineering characteristic is calculated by formula (8). And the result is as follows:
The final importance weight of engineering characteristic is calculated by formula (9). And the result is as follows: Here, the corresponding shape parameters are −600000, 200000, −1.3, and −0.5, 2.
The object function of the optimization problem is as the following formula:  The model is solved using G-ICA and RST2ANU in the same computer and every algorithm runs 10 times. The results are as the first row in Table 3. On the premise of getting similar solutions, the average FEs of RST2ANU are much smaller than that of G-ICA, but G-ICA's average time is still shorter. All of these two algorithms' solutions are clustered together. And one solution is taken as the decision reference: 99, 2.31, 5.54, 4.53, 2.80, 1, 0, 1, 1, 1) .

Mathematical Problems in Engineering
The solution means that, considering the limited resource in this design phase, the engineering characteristic "noise level" can be ignored for the time being. The investment costs of the other four engineering characteristics are 1.99, 5.54, 4.53, and 2.80, respectively. And the overall customer satisfaction is 1375.4.
Beyond that, another three cases adapted from the above case are also implemented. 10, 15, and 25 possible engineering characteristics are assumed in cases 2, 3, and 4. The other problem parameters are generated randomly. Cases 2 and 3 are implemented 10 times, but case 4 is only implemented 1   time because of its relatively long running time. The results are listed in Table 3, and the time advantage of G-ICA is very obvious.

Conclusions and Future Works
Determining the engineering characteristics to be implemented among so many possible selections and their corresponding actual fulfillment levels is an important QFD decision problem. But actual engineering problems' increasing complexity may lead to a huge QFD matrix, which is timeconsuming to analyze and tackle it. In order to improve QFD decision efficiency, an improved algorithm combining ICA and GA, called G-ICA, is proposed and applied to solve a MINLP model abstracted from the practical decision problem.
In the process of defining problems, the relationship between investment cost and fulfillment levels of engineering characteristics is focused. This relationship is generally very complicated in practice and is often simplified to linear. Here, a new relationship expression is proposed to express linear, concave, convex, and mixed type with different values of shape parameter . This shape parameter is assumed to be determined by the enterprises according to their experience, but with no more discussions in this paper. In the future work, more influence factors in practice will be taken into account while establishing QFD models, in order to accord with the actual product development situations better.
The G-ICA algorithm is compared with GA, ( + )-ES, ( , )-ES, SSA, RST2U, and RST2ANU by calculating 7 MINLP problems in the same personal computer. Three indexes, percentage of success , average time of successful runs , and average FEs in the case of successful runs, are used to evaluate these algorithms. The results show that, comparing with other algorithms, the G-ICA algorithm has a satisfactory global optimization and runtime performance Mathematical Problems in Engineering 11 in solving MINLP problems. Finally, 4 cases of washing machine with different scales are solved by G-ICA and RST2ANU. For every case, both algorithms' optimal results are clustered around the same point and can be as references for designers, to determine engineering characteristics and actual fulfillment levels. Meanwhile, the results show the G-ICA's better runtime performance than RST2ANU again.