Online Knowledge Extraction and Preference Guided Multi-Objective Optimization in Manufacturing

The integration of simulation-based optimization and data mining is an emerging approach to support decision-making in the design and improvement of manufacturing systems. In such an approach, knowledge extracted from the optimal solutions generated by the simulation-based optimization process can provide important information to decision makers, such as the importance of the decision variables and their influence on the design objectives, which cannot easily be obtained by other means. However, can the extracted knowledge be directly used during the optimization process to further enhance the quality of the solutions? This paper proposes such an online knowledge extraction approach that is used together with a preference-guided multi-objective optimization algorithm on simulation models of manufacturing systems. Specifically, it introduces a combination of the multi-objective evolutionary optimization algorithm, NSGA-II, and a customized data mining algorithm, called Flexible Pattern Mining (FPM), which can extract knowledge in the form of rules in an online and automatic manner, in order to guide the optimization to converge towards a decision maker’s preferred region in the objective space. Through a set of application problems, this paper demonstrates how the proposed FPM-NSGA-II can be used to support higher quality decision-making in manufacturing.


I. INTRODUCTION
It is widely acknowledged that manufacturing is the engine of the modern economy, due to its leading contribution to overall productivity and its multiple effects on growth in the rest of the economy [1]. For manufacturing companies to stay competitive, efficient manufacturing operations are required, both in their design phase and during continuous improvements that necessitate changes and investments to improve their efficiency. Central to the decision-making process in designing or improving manufacturing operations is the concept of knowledge -the more a decision maker understands the relationship between how to design and run a system and how it performs, the better are the decisions that can be made [2]. However, knowing what actions to The associate editor coordinating the review of this manuscript and approving it for publication was Yongming Li . perform to maintain and improve efficiency can be hard, due to the complex nature of manufacturing processes. One of the tools for analyzing manufacturing processes is simulation, where a model is built to depict a real-world process or system. Simulation has historically been used for prediction and performance analysis of different aspects of manufacturing, for example, workforce planning, supply chain management and capacity planning [3], as well as improving resource efficiency [4]. Some recent comprehensive reviews of the use of simulation in manufacturing can be found in [5] and [6]. Although simulation is often used as an evaluative tool, it cannot generate any optimal solution when used on its own. Consequently, researchers have shown a significant interest in merging optimization and simulation in order to find the optimal or near-optimal design and/or operating policies for manufacturing systems [5]. In such a simulation-based optimization (SBO) procedure, an optimization algorithm is used  together with a simulation model for evaluating candidate solutions in an iterative manner as shown in Figure 1. It has huge potential to generate optimal or near-optimal solutions that can be used to assist decision-making in designing and improving manufacturing systems [7], [8].
Optimization problems within manufacturing as well as many other domains often involve multiple conflicting objectives [9]. A multi-objective optimization problem can be defined as: Minimize F(x) = {f 1 (x), . . . , f M (x)} subject to g j (x) ≥ 0 ∀j = 1, 2, . . . , J h k (x) = 0 ∀k = 1, 2, . . . , K where f i : R n → R are two or more objectives (M ) that are conflicting and need to be minimized. The vector of variables x = [x 1 , x 2 , . . . , x n ] T belongs to the non-empty feasible region formed by the inequality constraints g j (x), equality constraints h k (x) and variable bounds. The feasible region is defined by the constraints of the problem and bounds on the variables. An important concept within multi-objective optimization is that of dominance which can be used to compare two solutions. A variable vector x 1 is said to dominate another variable vector x 2 (x 1 x 2 ) if and only if: , ∀i ∈ {1, 2, . . . , M }, 2) ∃j ∈ {1, 2, . . . , M } such that f j (x 1 ) < f j (x 2 ). If neither x 1 x 2 nor x 2 x 1 , then they are said to be non-dominated with respect to each other and denoted as x 1 || x 2 . A feasible solution x * is said to be Pareto-optimal if there does not exist any other feasible solution dominating x * . The manifold containing the Pareto-optimal solutions is called the Pareto-optimal front.
The n-dimensional space formed by the variables is called the decision space. Each solution in the decision space maps to a point in the M -dimensional space formed by the objectives, which is called the objective space. Figure 2 shows the mapping between the decision space and objective space in a problem with two variables and two conflicting minimization objectives. The Pareto-optimal solutions are shown in dark.
Simulation-based multi-objective optimization problems are difficult to solve using classical optimization methods due to various reasons. Firstly, as explained above, multiple objectives lead to multiple solutions. Most classical optimization methods can only generate one solution at a time, therefore, in order to find the Pareto-optimal set, several optimization runs have to be carried out. Secondly, by definition, a simulation model has non-analytical objective functions. Classical optimization techniques often rely on objective and constraint function gradients which are either unavailable or difficult to calculate for simulated outputs. Thirdly, due to the stochastic nature of objective and constraint function, the resulting multi-objective optimization problem can have several locally optimal solutions and classical optimization techniques are prone to such solutions. The three main issues described above can be overcome using evolutionary algorithms. Evolutionary algorithms aim to mimic principles of natural evolution observed in the nature by using a population of individuals (solutions) that evolve over several generations (iterations) to generate the optimal solutions [10]. Evolutionary algorithms do not require gradient information and, due to their stochastic nature, they also have the ability to recover from local optima. Evolutionary algorithms are especially popular in solving multi-objective optimization problems due to their population-based approach which allows obtaining multiple optimal solutions simultaneously. Evolutionary algorithms are also better suited for parallelization.
A decision maker's expert knowledge can be advantageous when solving complex multi-objective optimization problems, if his/her preferences can somehow be embedded into an evolutionary multi-objective (EMO) algorithm [11]- [13]. There are many different ways of incorporating decision maker's preferences. For example, by using different weights for each objective, the decision maker can define the relative importance of one objective over another, and then minimize the weighted sum of the objectives, essentially converting the multi-objective optimization problem into a single objective optimization problem. However, such scalarization methods have some drawbacks and do not always work as expected. For example, the simple weighted sum approach described here is only applicable to problems where the Pareto-optimal front is convex [14].
Combining SBO with evolutionary algorithms and data mining can be advantageous as a tool for extracting knowledge about a manufacturing system. The manual use of these techniques can be a complex undertaking and, therefore, they should be integrated into a decision support system. However, to our best knowledge, there are hitherto no decision support systems that combine these techniques [15].
The aim of this article is therefore to study how knowledge extracted from combining data mining methods with SBO can be used to better understand the implications of optimal solutions in manufacturing system problems for the purpose of decision-making. The study also explores how online knowledge extraction can be used to guide the SBO process based on the decision maker's preference. The paper is structured as follows: Section II reviews previous related research studies. Section III describes a technique that uses data mining to extract knowledge during and after an optimization run. Section IV proposes an extension of a multi-objective evolutionary algorithm that uses knowledge extracted through data mining to automatically guide an optimization towards the preferences indicated by a decision maker. Then, in Section V, the proposed methods are demonstrated with applications from the manufacturing industry. The obtained results are discussed in Section VI, and finally, the paper concludes in Section VII.

II. RELATED WORK
In this section, we present a review of previous works related to the three methodological aspects of this work. Section II-A discusses data mining methods that can generate knowledge in a form suitable for integration with an optimization algorithm. Section II-B covers existing preference-based evolutionary algorithms and their corresponding preference articulation methods. And finally, Section II-C explores existing approaches to combining data mining methods with evolutionary algorithms.

A. DATA MINING METHODS FOR KNOWLEDGE DISCOVERY
Several data mining approaches have been used in the literature for knowledge discovery from the solutions obtained through multi-objective optimization. A recent survey of such methods [16] classifies them based on the representation of the obtained knowledge: (i) descriptive statistics, which are simple numerical measures that summarize data, leads to explicit knowledge, (ii) visual data mining techniques represent knowledge through various graphical means such as clusters or heatmaps, thus they have an implicit representation, (iii) machine learning methods can represent knowledge in both implicit and explicit forms. Implicit knowledge is knowledge without a formal notation and therefore it is often presented visually, which not only leads to subjective interpretations, but also makes it difficult for an algorithm to use the knowledge automatically. Explicit knowledge, on the other hand, uses mathematical notation which is well-defined and complete so that no information is lost [17]. Explicit forms of knowledge have a fixed mathematical notation and are therefore easy to store, transfer and process automatically within a computer program. Hence, here we only focus on data mining methods that can discover explicit knowledge.
The most common type of explicit knowledge takes the form of if A then B rules involving the variables in the antecedent A and objectives in the consequent B. The earliest application of rule extraction in multi-objective optimization can be found in [18], where each objective is discretized into different levels and the rules involving the variables are obtained using classification trees. A similar method involving regression trees is used in [19] and [20] to generate rules describing a preferred region in the objective space. Both works demonstrate how the obtained rules can provide a better understanding of complex SBO problems to decision makers. Rough set theory and association rule mining have also been used to obtain if-then rules. The former has been applied to multi-objective design optimization of a centrifugal impeller in [18], [21], while the latter has been shown to be effective in [22] on the same problem. Another approach called Multi-Objective Rule Extraction, specifically designed for obtaining rules for SBO problems in production systems, is proposed in [23]. The method extracts rules by finding critical values for each variable that give significant jumps in objective function values. This gives the decision maker information about which configuration is best suited for production system given the desired level of investment. The authors show the application of the method in real-world production system optimization.
Another common type of explicit knowledge takes the form of mathematical relationships involving variables, objective functions and constraint functions. The simplest form of mathematical relationships are linear correlations and these can easily be obtained through methods like Principal Components Analysis (PCA) and Linear Discriminant Analysis (LDA). These methods have traditionally been used for graphical representation of solutions [24], [25], but can also be used for knowledge extraction. For example, a closely related method called Proper Orthogonal Decomposition has been employed in [26] and [27] to extract design knowledge from the Pareto-optimal solutions of an airfoil shape optimization problem.
While linear correlations are indeed insightful when present in optimization data, often the relationships between any given set of variables are nonlinear. Moreover, the nature of the relationships may be different in different regions of the search space. Innovization or (innovation through optimization) proposed in [28] can handle such nonlinearities through a manual approach of analyzing scatter plots of different variable combinations and performing regression with appropriate models on the correlated parts of the data. Later, an automated innovization approach was developed [29], [30] which uses unsupervised machine learning to extract multiple nonlinear relationships simultaneously without human intervention [31]. This approach was further extended in [32] to obtain mathematical relationships of arbitrary kinds using genetic programming. Automated innovization has been used for knowledge discovery in several real-world multi-objective optimization problems [33].

B. PREFERENCE-BASED EVOLUTIONARY ALGORITHMS
Preference-based evolutionary algorithms use information about the decision maker's preferences to obtain a preferred set of solutions on the Pareto-optimal front. These algorithms can be classified into three categories based on when the user input is needed, a priori, a posteriori and interactive methods. A priori methods take the decision maker's preferences into account right from the start of an optimization run. It is assumed that the decision maker is able to elicit the preference information accurately at the outset. On the other hand, a posteriori methods require preferences only after obtaining a representative set of trade-off solutions. While this helps the decision maker to realistically set his/her expectations, generating the entire trade-off front in the first place can itself be a challenging task, particularly when many objectives are involved. Interactive methods try to combine aspects of a priori and a posteriori methods by involving the decision maker during the optimization run in an iterative manner to guide search towards his/her preferences [34]. While this generally reduces the computational effort required to obtain preferred solutions, it drastically increases the number of times the decision maker's inputs are needed, thus risking cognitive overload.
Preference-based evolutionary algorithms can also be classified based on how the user specifies his/her preferences. One of the most common ways of articulating preferences is to use one or more reference points. A reference point is a point in the objective space which specifies the aspirations of the decision maker for each objective function. Depending on its location, the reference point can either be attainable (achievable) or unattainable (unachievable). An attainable point is one that lies in the feasible space and it is therefore possible for an optimization algorithm to reach the point. An unattainable point is one that is specified in the infeasible region, in which case the optimization algorithm can approach it only as close as the feasible space allows. Preferences can also be obtained through pairwise comparisons between representative solutions. However, in practice, this approach can be time-consuming and too demanding of the decision maker. A more recent method of preference articulation involves specifying a region of the objective space through the action of brushing using mouse input as presented in [35]. This is a visual and interactive method and as in the case of reference points, the decision maker does not need to have prior knowledge of the range of values in the objective space. Other common ways of articulating preferences is through reference directions, light beams and value functions [36]. Here, we specifically review preference-based evolutionary algorithms that use reference point(s), since the algorithm proposed in this paper is of a similar nature.
The reference point-based NSGA-II or R-NSGA-II introduced in [37], uses the weighted Euclidean distance of each solution from the reference point as a metric to promote diversity close to the reference point. It also employs a clearing strategy to avoid overcrowding of solutions at a particular reference point by controlling the extent and distribution of obtained solutions. While R-NSGA-II uses a priori elicitation of preferences, an interactive version of the algorithm proposed in [38] relies on pairwise comparisons between selected solutions at regular intervals. Based on the inputs from a decision maker, the algorithm creates a utility function which in turn guides search.
The Preference-Based Evolutionary Algorithm (PBEA), introduced in [39], incorporates preferences into an evolutionary algorithm using modified quality indicators. The decision maker's preference is expressed as reference points which are used in the selection operator of the algorithm together with an achievement scalarizing function to guide the search towards the region of interest. A weighted form of achievement scalarizing function (WASF-GA) [40] has also been developed where the weights are chosen so as to project the reference point onto a region of interest on the Paretooptimal front.
The concept of g-dominance proposed in [41] works by redefining Pareto-dominance such that solutions that dominate or are dominated by the reference point are considered to be preferable to solutions that non-dominated with the reference point. The proposed g-NSGA-II algorithm requires minimal changes to the base algorithm and can be utilized in both a priori and interactive manner. However, the reference point needs to be fairly close to the Pareto-optimal front in order to avoid obtaining a large region of interest.
The concept of r-dominance introduced in [42] also redefines Pareto-dominance such that a strict partial order can be obtained among non-dominated solutions. This allows solutions closer to the reference point to be preferred over other solutions of the same rank. The authors demonstrate that the approach is effective in both a priori and interactive forms.
The recently proposed RMEAD2 [43] and MOEA/D-NUMS [44] integrate user-preferences into the decomposition-based evolutionary algorithm MOEA/D by biasing the weight vectors towards the given reference point.
It is noteworthy that most of the preference-based evolutionary algorithms in the literature use NSGA-II as the underlying multi-objective optimizer. Since this is also true for the algorithm proposed in this paper, for completeness we briefly describe NSGA-II later in Section III-A.

C. COMBINING EVOLUTIONARY ALGORITHMS WITH DATA MINING
The main purpose of combining evolutionary algorithms with data mining methods is to improve the performance of the algorithm through the knowledge extracted during VOLUME 9, 2021 the optimization process. Such algorithms may be referred to as online knowledge-driven optimization algorithms as proposed in [16]. Since the knowledge needs to be processed automatically by an algorithm, only explicit representations can be used.
One of the early works in this area is the Learnable Evolution Model (LEM) proposed in [45]. It uses machine learning to obtain rules that differentiate between high-performing and low-performing solutions. These rules are used in turn to generate new solutions. While LEM was only applied to single-objective optimization problems, it showed the potential of embedding extracted knowledge into the optimization algorithm. Later, the method was extended in [46] to handle multi-objective optimization problems.
Linear correlations obtained by PCA can also be utilized to enhance search. For example, the approach proposed in [47] detects correlations in the decision space and builds a probability model that is used to generate solutions in subsequent iterations. Identification of linear and nonlinear correlations in the objective space can also help the search process when redundant objectives are removed online as demonstrated in [48]. This approach is especially useful for solving many-objective optimization problems.
The interactive method presented in [49] uses Dominancebased Rough Set approach to build a preference model based on the decision maker's ranking of certain representative solutions. The preference model consists of if-then-else rules which are in turn added as constraints to the original optimization problem, thus guiding the algorithm towards preferred regions of the search space.

III. ONLINE KNOWLEDGE EXTRACTION
In this section, we first describe the NSGA-II algorithm, a well-known multi-objective evolutionary algorithm, to facilitate our discussion later on how it can be extended with online knowledge extraction based on the decision maker's preference. Next, we present the Flexible Pattern Mining (FPM) approach that is used in this work for extracting knowledge. Thereafter, we first describe how FPM is used in a decision support system to extract knowledge in an offline manner, i.e. after the optimization run is complete. And finally, we present the extension of FPM-based decision support system to facilitate online knowledge extraction, i.e. during the optimization run.
A. THE NSGA-II ALGORITHM Genetic algorithms are inspired by the process of evolution seen in nature. A typical genetic algorithm starts with initialization of a set of random solutions called the parent population. Based on their fitness with respect to the objective function, a subset of these solutions is selected to form the so called mating pool. A child population is generated from the mating pool by applying the genetic operators, crossover and mutation on the solutions. Based on a crossover probability, each crossover operation recombines two parent solutions to form two recombinant solutions. Similarly, based on the mutation probability, the mutation operator mutates each recombinant to form a child solution. The parent and the child population are then merged and the new parent population is formed by selecting the best fitness solutions while maintaining the population size. This process is repeated for a fixed number of iterations or generations.
There are many multi-objective optimization techniques based on genetic algorithms. Non-dominated Sorting Genetic Algorithm II (NSGA-II) [50] is once such technique that has been used as a template for many other multi-objective optimization algorithms in several different application areas [9]. NSGA-II has been chosen in this work because it is well-known and has relatively good performance on multi-objective test problems as well as real world optimization problems [50]. Like genetic algorithms, NSGA-II is also a population based technique, which allows it to converge to multiple Pareto-optimal solutions in a single algorithmic run. The solutions in the population are evolved with crossover and mutation and each new set of individuals is called a generation. The individuals in a generation is evaluated with, in our case, a test function for a test problem or a simulation model. To create a set of individuals for a new generation the current generation is sorted based on the concept of dominance and crowding distance. Based on this sorting a new population is selected and the process is repeated until a stopping criterion is met, for example, maximum number of generations or running time. This process of selecting a new generation is shown in Figure 3 [50]. The child population C t is created from the population P t . These two populations is combined (R t ) and are sorted with non-dominated sorting into fronts (F 1 , F 2 , F 3 ,. . . , F l ) where the individuals from the best fronts are used to create the new generation. If a front does not fit into the new generation it is sorted with crowding distance and the individuals with the highest crowding distance is kept for the new generation (P t+1 ). The algorithm is described in detail in [50] and pseudo-code for the algorithm can be seen in Algorithm 1.  [51] to the popular sequential pattern mining (SPM) algorithm proposed in [52]. The original SPM algorithm was designed to extract frequent customer sequences (also known as sequential patterns) from market basket data. Given a dataset of transactions of different customers, called the transaction database, the SPM algorithm first generates a database of customer sequences. A customer sequence is essentially a chronologically ordered list of itemsets, each itemset itself being a group of items bought together by the customer. Frequent customer sequences are those that meet a minimum support threshold. FPM uses a similar approach, by considering the solutions generated during optimization as transactions and the variables as items. Thus, each transaction will have the same number of items, which is not the case in SPM. When the variable values are continuous or ordinal in nature, FPM can also generate sequential rules instead of sequential patterns by transforming the optimization data into a truth matrix. For each possible value d ij of a variable x i present in the solution set, FPM adds three binary columns in the truth matrix corresponding to the logical expressions x i < d ij , x i == d ij and x i > d ij . Frequent patterns can then be found within the truth matrix. For example, if x 1 > d 12 is true for 100% of the solutions and x 2 < d 24 is true for 80% of the solutions, then the sequential rule x 1 > d 12 ∧ x 2 < d 24 has a support of 80%.

C. DECISION SUPPORT SYSTEM FOR KNOWLEDGE DISCOVERY
To support decision makers with the task of SBO in manufacturing applications, the authors have previously proposed a Decision Support System (DSS) in [15]. This DSS, implemented in C++, has the ability to perform multi-objective optimization using NSGA-II on discrete event simulation models of real-world production systems and thereafter explore the obtained Pareto-optimal solutions. Here, we describe how FPM was integrated into this DSS to enable both offline and online knowledge discovery.  The FPM algorithm is added as the third component of the DSS as shown in Figure 4. It uses the traditional SBO process to feed optimization data to the FPM algorithm in real-time through an SQL database containing an archive of solutions from all generations. The archive also contains the rules generated by the FPM algorithm at the end of each generation of NSGA-II in Algorithm 1. The solutions and rules from any generation can be visualized by the decision maker through another interface called the Experiment Browser.
The Experiment Browser is simply a piece of software, implemented in C#, that connects to the internal storage of the DSS and provides the user with different tools to visualize and analyze the results. The Experiment Browser has scatter plots for two and three dimensional optimization data, parallel coordinates, as well as a table that displays the raw data from experiments. Figure 5 shows a screen-shot of the Experiment Browser with three open diagrams: parallel coordinates, scatter plot, and FPM rule plot. The FPM rule plot (shown as 'Online Analysis') is used to visualize the rules obtained by the FPM algorithm at the last generation. Each vertical bar represents the bounds of a variable in the optimization problem. Since FPM rules are of the form x i < d ij or x i > d ij , they can be interpreted as subsets of the original variable bounds and shown in color on the vertical bars.

IV. FPM-NSGA-II: A PREFERENCE GUIDED EVOLUTIONARY ALGORITHM
The FPM algorithm described in Section III-B is capable of generating knowledge in the form of explicit rules from VOLUME 9, 2021 intermediate generations of an evolutionary multi-objective optimization algorithm. As mentioned before, a unique advantage of explicit knowledge is that it can be processed automatically by an algorithm. Our proposed preference guided evolutionary algorithm is called of FPM-NSGA-II due to the integration of FPM into NSGA-II. Like most other a priori and interactive methods, FPM-NSGA-II uses a reference point for articulating user preferences. After every T generations, where T is set by the user, FPM-NSGA-II calls a function called RunFPM() shown in Algorithm 2. This function first retrieves the current non-dominated from, currentNDF from the archive A of evaluated solutions. After normalizing the objectives to [0, 1], the k nearest Euclidean neighbors of the reference point, refPoint, are selected from the currentNDF as shown in Figure 6 to form the set selectedSols. This set of solutions is assumed to represent the current preferences of the decision maker. Any significant rules extracted from this subset of non-dominated solutions will thus capture the decision maker's preferences in the decision space. The solutions that do not represent the decision maker's preferences, unselectedSols, are obtained by taking the set different between currentNDF and selectedSols. Next, FPM is used to obtain rules corresponding to selectedSols. The support for each rule in both selectedSols and unselectedSols is also calculated.

Algorithm 2 Pseudo-Code for RunFPM(A, k)
currentNDF ← GetNonDominated(A) selectedSols ← kNN(currentNDF, refPoint) unselectedSols ← currentNDF \ selectedSols rules ← FPM(selectedSols, unselectedSols) In order to guide the optimization towards the preferred regions (in this case, towards refPoint), the rules obtained by FPM are added as probabilistic constraints to the current generation. In subsequent generations, the probabilistic constraints penalize solutions that do not conform to the decision maker's preference, while promoting solutions close to refPoint. Probabilistic constraints act like regular inequality constraints, with the only difference that not all solutions are affected by the same constraints. Each rule from FPM that is added as a probabilistic constraint has an associated probability equal to its support. Thus, a support value of 95% means that the rule acts as a constraint with a probability of 0.95 on each solution in the current generation. A solution that adheres to the rule is not affected by this constraint. However, on the average 95% of the solutions that do not follow this rule will be penalized. Such probabilistic constraints help in maintaining some diversity among the solutions, that avoiding premature convergence.
The crossover and mutation operators are also modified using the rules to promote creation of children that adhere to the rules obtained from FPM. Simulated binary crossover (SBX) [50] and polynomial mutation operators are used in this work. Both the operators are modified by altering the bounds of the variable being recombined or mutated based on the FPM rule corresponding to that variable. For example, if a variable x i with original bounds x is associated with the FPM rule x i > p i , then the bounds are modified as p i , x (U ) i in SBX and polynomial mutation. The support of the rule is again used as the probability with which this modification is applied. The net effect is a change in the probability distributions of the offspring in the current and subsequent generations.
Together, the probabilistic constraints and modified SBX and polynomial mutation operators guide the solutions towards the preferred region as specified by the reference point. The complete FPM-NSGA-II algorithm is shown in Figure 7.
In this work, we set the number of generations between FPM calls to T = 10 generations. After every T generations, existing probabilistic constraints are purged and all modifications to the crossover and mutation operators are reversed. Thus, theoretically the algorithm has the ability to recover from premature convergence towards refPoint.  Note that since the online knowledge extraction procedure is a module attached to the regular NSGA-II algorithm, the latter can therefore easily be replaced with any other state-ofthe-art multi-objective optimization algorithm.
FPM-NSGA-II is fundamentally different from existing preference-based evolutionary algorithms in two main ways. Firstly, as discussed in Section II-B, algorithms such R-NSGA-II, PBEA and g-NSGA-II work via the selection mechanism by modifying the dominance operator to emphasize solutions closer to refPoint. On the other hand, FPM-NSGA-II works via the variation mechanisms by modifying crossover and mutation operators as described above. It also affects the selection process indirectly through probabilistic constraints that change every T generations. A more important distinction however is the fact that unlike existing preference-based evolutionary algorithms, FPM-NSGA-II also generates a set of rules at the end that capture decision maker's preferences in the variable space.

A. COMPUTATIONAL COMPLEXITY OF FPM-NSGA-II
The time complexity of NSGA-II algorithm is determined by the non-dominated sorting procedure depicted in Figure 3, which is known to be O(MN 2 ) [50], where M is the number of objectives and N is the population size. The computational runtime of FPM is determined by the underlying Apriori algorithm [53] and depends on four factors [54], namely (i) support threshold, (ii) number of transactions (N ), (iii) number of unique items (d) and (iv) average transaction width (w).
As mentioned before, in FPM all transactions (solutions) have the same number of unique items (variables). Therefore, the average transaction width is w = d = n, n being the number of variables. The worst case complexity of O(2 n ) occurs when support values for all possible variable combinations need to be calculated. However, this upper bound is practically only possible when support threshold is zero. In FPM-NSGA-II, only unary rules of the forms x i < d ij , x i == d ij and x i > d ij are sought, which only requires O(Nn) operations. Thus, the overall computational complexity of FPM-NSGA-II is either O(MN 2 ) or O(Nn), whichever is greater. From this we can conclude that the FPM part of FPM-NSGA-II scales linearly with the number of solutions and number of variables. In the application problems used in this paper, we have N > n, so the runtime of FPM-NSGA-II is limited by the non-dominated sorting procedure itself.

V. APPLICATION STUDIES
In this section, we describe two different production system design problems that are used in this work to demonstrate (i) online knowledge extraction through the developed DSS, and (ii) the proposed FPM-NSGA-II algorithm.
The first problem is based on a machining cell model extended from the one presented in [55] and is called 'PSU2'. Although the model is hypothetical, it is realistically designed based on the authors' experience in working with real-world complex models and decision-making situations in industry. This model makes use of three so-called selection objects in the simulation software FACTS Analyzer (or FACTS hereafter) [56] that can enable several design alternatives to co-exist in a model. By changing the values of the corresponding decision variables, the optimization module in FACTS can then alter the selection objects into various combinations. The analysis of such a model can then be conducted efficiently, by letting the optimization algorithm seek the optimal combination of design alternatives with their optimal settings, together with the optimal values of other decision variables (e.g. buffer spaces), in a single optimization run. The model is shown in Figure 8. The optimization problem formulation can be found in [15] with a complete description of the objectives and variables. Here it suffices to say that the model involves three objectives, minimize payoff time for investments (calculated as the ratio of total investment to profit per year), minimize work-in-process (WIP) and maximize throughput. The variables, their possible values and investments costs associated with them are shown in Table 1. The second problem under consideration is based on a scalable simulation model depicting simple flow lines that manufacture a single product type. Although the simulation model is simple, it has the interesting property of being scalable to any desired number of stations (operations). This allows the model to be used for benchmarking SBO algorithms and for understanding the behavior of other complex models under scaling. In this paper, we consider simulation models with four different number of stations: 5, 10, 15 and 20, hereafter referred to as 5S, 10S, 15S and 20S respectively. Between every two stations, there is one buffer whose capacity can be adjusted. Thus, for s stations, the total number of buffers is s − 1. Figure 9 illustrates the 5S scalable model which starts and ends with a source and a sink, and with the operations and buffers alternating in between.
An optimization problem is set-up for each of the scaled models with the same type of configuration. The variables for each station are the processing time, availability and mean time to repair. The buffer capacities are also included as variables. Table 2 shows the notation and possible values for each variable. The subscript i in the notation refers to the operation or buffer number, starting from i = 1 immediately after the source. The total number of variables is n = 3s + (s − 1). Thus, 5S, 10S, 15S and 20S respectively involve 19, 39, 59 and 79 variables to be optimized. The optimization problem has three objectives: (i) maximizing the throughput of the line, (ii) minimizing the total number of improvements in processing time, availability and mean time to repair from the base values of T i = 80, A i = 90 and R i = 300 ∀ i, and (iii) minimizing the total number of buffer spaces required. An optimization of this problem should achieve the optimal trade-off between the throughput, number of improvement actions, and total buffer capacity.

VI. RESULTS AND DISCUSSION
Each experiment is composed of a simulation model, the experiment type and the parameter settings. To account for the stochastic nature of evolutionary algorithms, each experiment is replicated 10 times. The evaluations and the FPM rules from each generation of all optimization runs are stored in the database, which can be retrieved and plotted either during or after the optimization run. Both the original NSGA-II algorithm and the proposed FPM-NSGA-II algorithm use the parameter settings shown in Table 3. For the mutation probability, the 1/n rule-of-thumb, where n is the number of variables, is used.
In addition to the above mentioned NSGA-II parameters, FPM-NSGA-II involves two other parameters, namely, (i) the number of nearest neighbors k from the given reference point refPoint to be selected from the current non-dominated front currentNDF, and (ii) the support threshold ST for FPM rule generation. The value of k determines how well the decision maker's preferences are captured in the objective space. When too many solutions are selected, the preferences may not be adequately emphasized to guide search towards the reference point. On the other hand, selecting too few solutions may mean that no rules are found that meet the support threshold. Ideally, the two parameters should be adaptively controlled during runtime, with a high k and low ST at the beginning when relatively fewer solutions are non-dominated, and a low k and high ST towards the end when a majority of the population is expected to be non-dominated. This way the emphasis on the preferences is gradually increased with each generation. However, the control algorithm for an optimal balance between the two parameters is a matter for future research. In this paper, we keep both k and ST static thoruhgout the runtime. In our experiments, we found that using a low k and moderate ST allows preferences to be emphasized right from the start while still allowing rules satisfying the support threshold to be found. The following results were generated with k = 20% the size of currentNDF and ST = 50%. Small changes to k ∈ [20%, 30%] and ST ∈ [50%, 60%] have shown to have negligible effect on the region and rate of convergence for all problems considered. The DSS developed in this work has an architecture that enables parallelization to distribute the evaluations of solutions on several computers, which reduces the total runtime. The implementation uses 100 computation nodes to parallely evaluate all solutions in a given generation, thus achieving a speed-up of nearly 100x. The generation of rules with FPM cannot be efficiently parallelized, and is hence performed on the local machine.
In Section VI-A, we demonstrate the online knowledge extraction process through the developed DSS from the solutions obtained by solving the above problems with the original NSGA-II algorithm. In this first scenario, the obtained rules apply to the whole non-dominated front and therefore can be interpreted as IF 'non-dominated', THEN 'rule'. Section VI-B presents how the FPM-NSGA-II algorithm uses the same DSS for online knowledge extraction and incorporates that knowledge to converge towards a specified reference point, thus meeting the decision maker's preferences. Thus, the second scenario allows finding rules that apply to a preferred region of the Pareto-optimal front, without having to generate a representation for the whole front. These rules can be interpreted as IF 'preferred', THEN 'rule'.

A. VISUALIZATION OF FPM RULES
When the experiment type is set to 'Online knowledge extraction', the non-dominated fronts from all generations of a complete optimization are analyzed. The final non-dominated front contains the trade-off solutions obtained by the original NSGA-II algorithm. By applying FPM to these solutions, knowledge about what makes a solution perform well can be discovered.
In order to visualize the FPM rules obtained from multiple optimization runs, a visualization method called the rule map has been developed. Figure 10 shows the rule map for the final non-dominated front from all ten replications of the 5S model. The diagram shows, for each input parameter, the range of variable values described by the FPM rules. The figure is divided into four parts representing the four groups of input parameters: processing time, availability, mean time to repair (MTTR), and buffer capacity. The y-axis represents the station/buffer number. For each rule in a replication, a translucent rectangle representing the value range of the rule and parameter is added to the diagram. Due to the stochastic nature of optimization, the rules from the last nondominated front may vary from one replication to another. However, when rules from multiple replications are added on top of each other, the overlapping areas become darker indicating the most reliable values for each variable.
In Figure 10, it can be seen that MTTR variables follow several rules which makes their map darker. For process times, availability, and buffer capacities, there are very few rules, which means that FPM could not find many significant rules for these variables when considering the whole trade-off front. Figure 11 shows the same analysis for the 10S variant of the scalable models. The behavior of the MTTR parameters is similar to the 5S model, but it has more significant rules for the other parameters. Rules obtained from the 15S model, presented in Figure 12, shows even more rules for process times, availability, and buffer capacities. Most processing time values are centered around 70 seconds. For the 20 station model, Figure 13 shows that the rules are approximately the same as for the 15S model. Such rule maps can help decision makers easily visualize the most frequent optimal values when dealing with a large number of variables. For example, the rule map of PSU2 in Figure 14 shows that B 1 and B 2 usually take values between 1 and 10, while B 3 and B 5 mostly take values between 10 and 30.

B. RESULTS FROM FPM-NSGA-II
The proposed preference guided evolutionary algorithm, FPM-NSGA-II, is now tested on the five application problems. The experimental setup is the same as described before, except that FPM-NSGA-II also requires a reference point for each simulation model. These reference points, which are all unattainable, are shown in Table 5. Our DSS can currently only handle one reference point as this is the most common use-case. The algorithm itself is not restricted with respect to the number of reference points.     Table 6. The table also shows the number of times the rules appear in 10 replications of FPM-NSGA-II. Note that due to the specification of a  preferred region (through reference points), the obtained rules are very different from those in Table 4.
Similar to the previous sub-section, the rule maps can be generated for all five models. These are shown in Figures 15-18. Note that there are more distinct and darker regions in these figures than in Figures 10-14. This is a graphical manifestation of the difference in rules observed in Table 6. The visual nature of rule maps helps decision makers to easily make sense of multiple variables at a single glance. For example, Figure 17 shows that processing times        can also generate rules that describe these solutions. In the classical sense, decision making is almost always performed in the objective space. FPM-NSGA-II allows decision makers to simultaneously understand the implications of their decisions in the decision space. This is not possible with any of the reference-point based evolutionary algorithms that exist in the literature.

VII. CONCLUSIONS
Although a combination of simulation and multi-objective optimization has previously been shown to be effective for supporting decision making in manufacturing, higher quality decisions can be assured when the optimal solutions are supplemented with the knowledge about what makes them optimal in terms of the variable values. Such knowledge can be automatically extracted using data mining techniques.
In this paper, we first extended an existing DSS, that combines a discrete event simulation engine and an optimization module, with a data mining module that implements the recently proposed Flexible Pattern Mining algorithm. The DSS is capable of replicating optimization runs on a parallel architecture and online knowledge extraction. Also part of the DSS is an Experiment Browser for visualizing the solutions and the obtained knowledge. In order to visualize the rules generated by FPM, we developed the concept of rule map, which allows FPM rules from multiple replications of the experiments to be shown on a single plot. The decision maker can obtain an overview of knowledge pertaining to the trade-off solutions at a glance. The paper also proposes a new preference guided multiobjective optimization algorithm called FPM-NSGA-II. The new algorithm incorporates the decision maker's preferences by (i) selecting solutions that are close to the specified reference point in intermediate generations, (ii) using the DSS to obtain FPM rules specific to the selected solutions, and (iii) adding probabilistic constraints and modifying the crossover/mutation operators of NSGA-II to promote the generation of solutions that satisfy the rules. In this way, the algorithm is able to converge towards a preferred region on the Pareto-optimal front. While there are several preferencebased evolutionary algorithms in the literature, the advantage of FPM-NSGA-II is that in addition to the set of preferred solutions, the algorithm also generates knowledge that describe those solutions in the form of explicit and easily interpretable rules.
The DSS for online knowledge extraction and FPM-NSGA-II, are applied to SBO problems involving discrete event simulation models of realistic production systems. We observed that knowledge extracted from the complete trade-off front can be different from that pertaining to a specific region of the trade-off front. In either case, the rule maps were found to be an effective way of communicating the extracted knowledge to the decision maker. The results shows that FPM-NSGA-II is able to extract relevant knowledge in an automatic and iterative manner to guide the evolutionary algorithm towards the decision maker's preferences. In each application study, we found that the rules generated by FPM-NSGA-II throw new light on understanding the implications of choosing a solution from the preferred optimal set. Such an approach to decision making is believed to provide the decision maker with additional information that cannot be obtained by any existing industrial practice today.
INGEMAR KARLSSON received the B.Sc. degree in computer science, the M.Sc. degree in automation engineering, and the Ph.D. degree from the University of Skövde, in 2006, 2011, and 2018, respectively. He is a Senior Software Developer with the University of Skövde. His research interests include simulation-based optimization, data mining, augmented reality, and cloud technologies.
SUNITH BANDARU (Senior Member, IEEE) received the bachelor's degree in mechanical engineering from Jawaharlal Nehru Technological University, India, in 2006, and the master's and Ph.D. degrees from IIT Kanpur, India, in 2008 and 2013, respectively. He is an Associate Professor of production and automation engineering with the University of Skövde, Sweden. He is also the Chair of the Committee for Research Education in Informatics. His research interests include the intersection of evolutionary computation, machine learning, and data mining, where he develops methods for knowledge discovery in multi-objective optimization.
AMOS H. C. NG received the Ph.D. degree in computing sciences and engineering from De Montfort University, Leicester, U.K. He is a Professor of production and automation engineering with the University of Skövde, Sweden. He is also a Visiting Professor with the Department of Civil and Industrial Engineering, Uppsala University, Sweden, and the CEO of Evoma AB. His main research interests include applying simulation, multi-objective optimization, and prescriptive analytics for manufacturing/service/health-care systems design, analysis, and improvement.