Robust parameter design of mixed multiple responses based on a latent variable Gaussian process model

Traditional robust parameter design methods mainly focus on the optimization of quantitative quality characteristics. However, computer experiments involving qualitative and quantitative mixed input and output occur frequently in the manufacturing industry, which prompted the authors to develop an effective meta-modeling and optimization technique for such experiments. This article combines the latent variable Gaussian process (LVGP) model and fuzzy set theory to create a mixed multi-response LVGP (MMR-LVGP) model involving qualitative and quantitative mixed input and output. Then, the optimization scheme is established by comprehensively weighing the location and dispersion effects of each quality characteristic to find the joint optimal solution of qualitative and quantitative factors. Numerical and industrial cases are used to illustrate the validity of the proposed method in the modeling and optimization of experimental data with qualitative and quantitative mixed input and output or spatio-temporal structure. The comparison results indicate that the proposed method is preferred over existing methods.


Introduction
Robust parameter design (RPD) is an important method for improving product quality proposed by Taguchi (1986) based on experimental design and signal-to-noise ratio.Its purpose is to minimize the variability of the quality characteristics while keeping them close to their target values by setting the levels of controllable variables.RPD is mainly used in the quality design phase of a product, which can reduce the uncertainty related to the product or process performance at the source and has produced significant economic benefits (Li, Xiao, and Gao 2019a;Li et al. 2019bLi et al. , 2020)).
It is supposed that all input factors and output responses are quantitative in the standard statistical framework of optimization techniques (Taguchi 1985;Johnson and Wichern 2007;Myers, Montgomery, and Anderson-Cook 2016;Ouyang et al. 2021).However, owing to the intrinsic features of quality characteristics or the limitations of measuring tools, the coexistence of qualitative and quantitative is common in industrial practice (Liu and Zhou 2006;Han et al. 2009;Li, Tsung, and Zou 2014;Huang et al. 2016;Ozdemir and Cho 2017).The following examples describe the three situations of qualitative and quantitative mixing.
(1) Qualitative and quantitative mixed inputs.In the simulation experiment of the solid end milling process, the radius of the anterior angle, posterior angle and helix angle are continuous-numeric-type quantitative factors.The number of chip removal grooves of the milling cutters (possibly 2, 3 or 4) is a quantitative factor that takes discrete values.Different from the previous two numeric-type factors, the kind of material used to make the workpiece, such as titanium alloy or high-speed steel, is a typical qualitative factor of the nominal type.If the situation of the tool is graded as 'bad', 'normal' or 'good', this will be another ordered qualitative factor (Joseph, Gul, and Ba 2020).
(2) Qualitative and quantitative mixed outputs.In the semiconductor manufacturing process (Spanos and Chen 1997), the printed pattern is 'fuzzy' or 'sharp', the surface is 'smooth' or 'rough' and the adhesion between the layers and the substrate is 'good' or 'bad'.These qualitative observations are usually captured together with observations of quantitative characteristics such as the thickness of each layer, the threshold voltage and the line width.
(3) Qualitative and quantitative mixed inputs and mixed outputs.Take the ion implantation process experiment given by Hsieh and Tong (2001) as an example.Six control factors are considered in this experiment.Among them, factor A is a two-level qualitative factor that takes the value 'Type 1' or 'Type 2', and B ∼ F are three-level quantitative factors of continuous-numeric type.
The ion implantation amount is a quantitative response with a target value of 1000.The uniformity of the ion implantation is subjectively divided into five categories, expressed by I ∼ V, which is a typical qualitative response.In addition to the above three cases, the spatio-temporal response or signal response is another data type that implies qualitative and quantitative mixed variables (Xiong, Li, and Wu 2020).For example, in the alternator experiment by Nair, Taam, and Ye (2002), if the rotational speed of the time factor is regarded as a qualitative factor, it can be recombined with the quantitative spatial factors, which can be classified as the first situation.
Imagine if the spatio-temporal data also contains qualitative and quantitative mixed response variables.Such spatio-temporal data is the type of data that will be studied with qualitative and quantitative mixed input and output.In many fields, such as manufacturing and engineering design, the input factors and output responses may be continuous or discrete numerical variables, or may be nominal or ordered categorical variables.Therefore, the meta-modeling and RPD of computer experiments with multiple types of factors and responses have received more and more attention in continuous quality improvement.
Since the acquisition of qualitative information usually requires expensive analytical techniques and the traditional global optimization algorithms for dual-response and triple-response systems cannot be applied directly to the optimization of qualitative and quantitative mixed multi-responses (Del Castillo, Fan, and Semple 1999;Fan, Fan, and Huang 2013;Del Castillo, Fan, and Semple 2018).Therefore, in manufacturing, qualitative observation is usually informally collected or subjectively regarded as tacit knowledge.Ordered categorical or linguistic descriptions are usually used to characterize the performance of qualitative features to capture important information in qualitative observations.For the ordered categorical response, Jeng and Guo (1995) incorporated both location effects and dispersion effects into the mean square deviation (MSD) and constructed a weighted probability scoring scheme (WPSS) to acquire the optimal level combination of control factors by minimizing the MSD.Hsieh (2007) took the fuzzy quality loss function (FQLF) as the criterion for calculating the optimal parameter settings and developed a new optimization method for ordered categorical responses by using fuzzy sets.In addition, Yang et al. (2021) proposed an hierarchical Bayesian model and used the satisfaction function as the optimization criterion to obtain the optimal solution for the parameters.However, the above methods only focus on a single qualitative response, which does not apply to the optimization of a qualitative and quantitative mixed multi-response.For the qualitative responses described by linguistics, the expert scoring method is usually used to quantify them.Since the conventional binary set cannot define the uncertainty involved in the linguistic description (Zadeh 1973), the expert scoring method often results in information loss due to ignoring the uncertainties of qualitative characteristics.These facts motivated the present authors to develop a new technology to optimize qualitative and quantitative mixed quality characteristics.
The Gaussian process (GP) can use finite samples to approximate the expensive calculationintensive function, which has been widely used to fit the nonlinear relationship between inputs and outputs (Sacks et al. 1989;Williams and Rasmussen 2006;Rougier 2008;Kleijnen 2010;Alshraideh and Del Castillo 2014;Gu and Berger 2016;Feng et a. 2020;Gu and Xu 2020).However, there is no specific distance measure between the different levels of qualitative factors (Deng et al. 2017;Kang and Deng 2020).There are three key challenges to incorporating qualitative and quantitative mixed input factors into GP modeling.
(1) How to build an appropriate covariance function for qualitative factors.
(2) How to establish a flexible relationship between the correlation functions of qualitative factors and quantitative factors.
(3) When the number of classification levels is large, how to accurately estimate the covariance structure between different levels of qualitative factors.
Some covariance structures have been proposed for qualitative and quantitative mixed input variables (Qian, Wu, and Wu 2008;Han et al. 2009;Zhou, Qian, and Zhou 2011;Deng et al. 2017;Zhang et al. 2020), and these methods will be reviewed in Section 2. However, the existing GP models assume that the experimental data only contains qualitative and quantitative mixed input factors.Therefore, establishing the model and optimization framework for experimental data with qualitative and quantitative mixed input and output is an important but unresolved problem.
This work aims to establish a model and optimization framework for computer experiments with qualitative and quantitative mixed input and output.The specific contributions are summarized as follows.
(1) The latent variable Gaussian process (LVGP) model proposed by Zhang et al. (2020) is combined with fuzzy set theory to develop a data-driven model for computer experiments with qualitative and quantitative mixed input and output.(2) An optimization scheme is established by comprehensively balancing the location and dispersion effects to find the joint optimal solutions for qualitative factors and quantitative factors.
(3) The proposed method can effectively deal with the parameter optimization problem of highdimensional spatio-temporal responses or signal responses, and obtains more robust optimization results while reducing the computational cost.
The remainder of this article is organized as follows.The purpose of Section 2 is to provide a brief review of existing Gaussian process models and the basics of fuzzy set theory.Section 3 introduces the modeling and optimization methods proposed in this article and their realization process.Then, the effectiveness and applicability of the proposed method are illustrated through a simulation experiment and two industrial cases in Section 4. Section 5 summarizes this article, including the limitations of the proposed method and possible future research directions.

Notation and background
This section briefly reviews the existing GP models, providing symbols and background for expanding the subsequent model.The article considers a computer experiment with input variables v = (x T , t T ) T ∈ R p+q and output response y(v), where x = (x 1 , . . ., x p ) T and t = (t 1 , . . ., t q ) T denote p quantitative factors and q qualitative factors, respectively.t j ∈ {1, 2, . . ., m j } represents the jth qualitative factor with m j classification levels, j = 1, . . ., q.

Gaussian process model with only quantitative inputs
When the input factors are all quantitative, the relationship between response y(x) and factor x is usually expressed as where μ is a constant mean , Z(x) is a GP with mean zero and covariance function φ(•, Here, σ 2 is the prior variance, R(•, • | θ) represents the correlation function for which the correlation parameter vector is θ = (θ 1 , . . ., θ p ).The Gaussian correlation function is commonly used for computer experiments with only quantitative input factors: which quantifies the correlation between Z(x i ) and Z(x j ) under any two input level combinations x i and x j .The parameters θ, together with μ and σ 2 , need to be estimated via the maximum likelihood estimation (MLE) method.

Gaussian process models with qualitative and quantitative mixed inputs
Similar to Equation (1), the relationship between qualitative and quantitative mixed input factors v and output response y(v) can be expressed as where μ is a constant mean and Z(v) is a GP with mean value zero and variance σ 2 .The following is a brief review of existing covariance structures for qualitative and quantitative mixed input factors.Qian, Wu, and Wu (2008) and Zhou, Qian, and Zhou (2011) used the product of the covariance matrix of qualitative and quantitative mixed input factors to establish a covariance model.Specifically, for any two inputs v 1 = (x 1 , t 1 ) and v 2 = (x 2 , t 2 ), they assumed that the covariance structure of y(v 1 ) and y(v 2 ) was

Unrestrictive covariance (UC)
The 2), and τ (j) t 1j ,t 2j denotes the correlation between the two levels of t 1j and t 2j for the jth qualitative input variable.There is a total of q j=1 m j (m j − 1)/2 parameters to be estimated in Equation ( 4).The UC model assumes that all quantitative input variables use a fixed covariance structure for qualitative variables, which is not flexible enough.

Multiplicative covariance (MC)
The function MC is simplified based on the function UC, and its assumption is that, for all t 1 = t 2 , τ (j) where φ j l > 0 is a parameter related to the lth level of the qualitative factor t j .It is easy to see from Equation ( 5) that each level of each qualitative factor corresponds to one parameter, so a total of q j=1 m j parameters need to be estimated during the MLE optimization process.

Additive Gaussian process (AGP)
To overcome the limitations of the multiplicative covariance structure in both the UC and MC methods, a more flexible additive covariance structure was presented by Deng et al. (2017): where the Z j 's are independent GPs with mean value zero and covariance function φ j .R(x 1 , x 2 | θ (j) ) is a Gaussian correlation function, as defined in Equation ( 2).θ (j) and σ 2 j are the correlation parameter vector and the prior variance of the qualitative factor t j , respectively.The definition of τ (j) t 1j ,t 2j is the same as Equation ( 4).
Although the AGP is more flexible than UC when multiple qualitative factors are present, it also has some limitations.For example, the total number of (1 + p) × q + q j=1 m j (m j − 1)/2 parameters to be estimated via MLE will increase rapidly as the dimensionality of the problem increases.

Latent variable Gaussian process (LVGP)
The LVGP method maps the m j levels of each qualitative factor t j to m j points in the jth latent space.The distance between any two points in a latent space can be defined.That is, any standard correlation function is applicable for points in the latent space, such as the following Gaussian correlation function: where T is the point corresponding to the lth level of the qualitative factor t j in the latent space, d z represents the dimension of the latent space and • 2 represents the Euclidean 2-norm.

Fuzzy set theory
The attribute values of the qualitative and quantitative mixed quality characteristics may include various types of expression such as an exact number, a fuzzy number and a fuzzy linguistic.Considering the uncertainty of the experimenter in expressing the experiment, introducing Zadeh's fuzzy set theory into robust parameter design can provide a flexible and appropriate form of information expression for multiple types of response (Zadeh 1965;Lee et al. 2016).
then A is said to be a triangular fuzzy number, and 0 < a L ≤ a M ≤ a U .If the triangular fuzzy number ) are both triangular fuzzy numbers, the distance between A and B is defined as 3. The proposed model and optimization method

The process of building the proposed model
This article extends the LVGP model to create a mixed multi-response LVGP (MMR-LVGP) model for computer experiments with qualitative and quantitative mixed input and output.Compared with non-spatio-temporal data, when using the proposed method to model spatio-temporal data, the structure of the spatio-temporal data needs to be converted first.While adding a signal factor that takes the same value for all levels, such as s = (1, 1, . . ., 1), as a time factor, the non-spatio-temporal data can also be regarded as spatio-temporal data.Therefore, taking spatio-temporal data as an example, the construction process of the proposed model can be explained more comprehensively.In the following, taking the spatio-temporal data containing qualitative and quantitative mixed input and output as an example, the construction process of the MMR-LVGP model is introduced through the following four steps.
Step 1: Quantify the qualitative response.In order to maintain the uncertainty of the qualitative response, it is common to use triangular fuzzy numbers to quantify the qualitative response described by linguistics.Refer to Table 1 for the specific conversion relationship between a linguistic fuzzy number and a triangular fuzzy number.
The quantitative and quantified qualitative responses form a mixed response matrix A = (a ij ) m×n , where a ij represents the observation at the ith spatial factor combination x i and the jth time factor t j .
To incorporate qualitative and quantitative mixed characteristics into MMR-LVGP modeling, the triangular fuzzy numbers in the mixed multi-response matrix A can be converted into an exact value using Equation ( 10): , and μ a ij (x) is the MF of the triangular fuzzy number a ij .
Step 2: Convert the data structure.Let Y = (y ij ) m×n be the response matrix after defuzzification, where y ij = y(x i , t j ) denotes the observation at the spatio-temporal location (x i , t j ).
is an (mn) × 1 column vector, and vec(•) is an operation that connects the columns of a matrix into a vector in left-to-right order (Alshraideh and Del Castillo 2014).Equations ( 10) and ( 11) cleverly convert qualitative and quantitative mixed spatio-temporal response data into the form of quantitative single-response data.This conversion overcomes the shortcoming of information loss caused by ignoring the influence of the time factor on output responses in traditional modeling methods.Simultaneously, it reduces the difficulty and computational cost of spatio-temporal data modeling.
Step 3: Latent variable representation of qualitative factors.The LVGP method offers a straightforward way of processing time factors by mapping their levels onto a two-dimensional (2D) numerical latent space.It has robust physical principles since the impact of any qualitative factor on a quantitative response is ultimately caused by a group of underlying numerical variables {h 1 , h 2 , . . .= h 1 (t), h 2 (t), . ..}.As shown in Figure 1, the three levels of the time factor correspond to points in the high-dimensional underlying space of {h 1 (t), h 2 (t), . ..}.Although the dimensionality of these potential numerical variables can be very high, their integral contribution to the response can be approximated by a 2D combination z(t) = g(h 1 (t), h 2 (t), . ..).
Step 4: Mixed multi-response latent variable Gaussian process modeling.After the latent representation in Step 3, the distance between any two points in the latent space can easily be defined.Under the model ( 7) given in Section 2.2.4, the log-likelihood function is where ln(•) is the natural logarithm, mn is the sample capacity, 1 is an mn-by-1 vector with elements all one, y is the mn-by-1 vector of observations, U q j=1 {Z j (1), . . ., Z j (m j )} is the set of latent mapping values for each level of qualitative factors, and R = R(Z, θ ) is the mn-by-mn correlation matrix calculated by inserting the sample value pairs of (x, t) into (7) for the elements.The MLE estimation of μ and σ 2 in ( 12) can be expressed by the correlation matrix: After substituting ( 13) and ( 14) into ( 12), the estimated values of Ẑ and θ can be calculated by maximizing the log-likelihood function: After obtaining the MLE estimates of Ẑ and θ, the MMR-LVGP model prediction at any new input where r(v * ) = (R(v * , v (1) ), . . ., R(v * , v (mn) )) is a vector composed of the pairwise correlation between v * = (x * , t * ) and each point v (ij) = (x i , t j ), i = 1, 2, . . ., m; j = 1, 2, . . ., n.
In addition, in order to quantify the uncertainty of prediction, the prediction error variance of Equation ( 16) is When there is uncertainty in the actual model, in order to consider the noise of the response, a random parameter λ can be added to all diagonal elements in the correlation matrix R, and the parameter λ is estimated together with Z and θ.

The proposed optimization method
In the mixed response matrix A = (a ij ) m×n , there is incommensurability between the responses (Qin 2003).On the one hand, each response has different units, orders of magnitude, and forms of taking values.On the other hand, according to the optimization objectives, the response can be divided into three different types: the larger-the-better (LTB) type, the nominal-the-best (NTB) type, and the smaller-the-better (STB) type.Therefore, directly calculating the contribution of each response to the variation is complicated and easily causes unreasonable parameter optimization results.To eliminate the influence of different physical dimensions, data types and optimization objectives on the optimization results, the range transformation method is used to normalize various responses and obtain a normalized response matrix B = (b ij ) m×n in this article. (

1) Normalization of the real number response
If the response value a ij is a real number, its normalized value is (2) Normalization of the fuzzy number response where M = 1, 2, . . ., m is the subscript set of different combinations of spatial factor levels, J 1 , J 2 and J 3 are the subscript sets of LTB-type, NTB-type and STB-type of the quality characteristics, and a * j is the target value of the jth quality characteristic.
The more the difference in the values of a quality characteristic across all combinations of input levels, the more unstable it is.Therefore, the quality loss weight of the jth quality characteristic can be calculated using the deviation maximization method: where where Q 1 and Q 2 are the subscript sets for quantitative and qualitative quality characteristics, respectively.
The deviation maximization method only uses the degree of difference in the values taken for the quality characteristic to determine its loss weight.It has clear principles and simple calculations, which can reduce the interference of experts' subjective factors to obtain a more objective loss weight result.
To improve the applicability of the mean squared error (MSE) method presented by Lin and Tu (1995) and comprehensively weigh the influence of the location and dispersion effects on the optimization results, a quadratic quality loss function is established for the three types of quality characteristics, LTB-type, NTB-type and STB-type: where y j (x, t), T j and σ 2 j are the predicted value, target value, and predicted variance of the jth quality characteristic, respectively.J 1 , J 2 and J 3 are the subscript sets for the quality characteristics of LTBtype, NTB-type and STB-type.
Under the modeling framework of MMR-LVGP, combining the loss weights and the quadratic quality loss function of each quality characteristic, the following optimization model is established: where and are the feasible regions of the quantitative factor x and qualitative factor t, respectively.Then, the global optimal solution of model( 24) can be found by genetic algorithm.
The basic flow of the proposed method is presented in Figure 2, and the specific steps are summarized as follows.
Step 1: Develop an experimental design plan and collect the experimental data.
Step 2: Refer to Table 1, converting the linguistic fuzzy numbers in the mixed response matrix into triangular fuzzy numbers.Transform the mixed spatio-temporal response into a quantitative single response format by Equations ( 10) and (11).
Step 3: Map the qualitative factors into a two-dimensional numerical latent space.
Step 4: Fit the MMR-LVGP model using the experimental data processed in Steps 2 and 3, and create a quality loss function based on the output responses and variances.
Step 5: Use the range transformation method to normalize the mixed response matrix.Then, each quality characteristic's loss weight is calculated using the deviation maximization method.
Step 6: The optimization scheme is established by comprehensively weighing the quality loss of each quality characteristic, using the genetic algorithm to find the jointly optimal solution for qualitative and quantitative factors.

Illustrative example
In this part, a simulation experiment and two industrial cases are used to illustrate the effectiveness and applicability of the proposed method.Univariate GP, multivariate GP (MGP) and the existing models for qualitative and quantitative mixed input given in Section 2.2 are selected as comparison methods, and the root mean square error (RMSE) given by Equation ( 25) is used to compare the prediction performance of different models: where ŷj (v * i ) and y j (v * i ) are the predicted value and true value of the jth quality characteristic at the ith input location v * i , respectively.

Numerical example
The specific definition of the simulation function is as follows: where and is an indicator variable with a value of zero or one.x i ∈ [0, 1], i = 1, 2, are quantitative spatial factors, and t ∈ {−1, −0.5, 0, 0.5, 1} is a qualitative time factor.
To test the dependence of the proposed model on the sample distribution, i.e. the model's generalization ability, the test samples were generated in two steps.In the first stage, the value of the quantitative factor x i is restricted to [0.1, 0.9].The maximum projection (MaxPro) design method is used to generate 14 design points for the qualitative and quantitative mixed input factors (Joseph, Gul, and Ba 2020).The output data is obtained by Equation ( 26) to generate 14 training samples.Then, from the 10,000 quantitative factor candidate points in the range [0, 1], an additional 16 design points were selected to augment the first design using the MaxPro criteria, generating a total of 30 test samples.
The distribution of the design points and augmented points of the quantitative factors are shown in Figure 3, and the generated simulation data are shown in Table A1 of the online supplemental data, which can be accessed at https://doi.org/10.1080/0305215X.2022.2124982.As can be seen in Figure 3, the first generated 14 design points for training samples fall within the square region [0.1, 0.9] × [0.1, 0.9] enclosed by the dotted line.The 30 design points used as test samples are generated by adding 16 design points in the region [0, 1] × [0, 1] based on the first 14 design points.
The values of the function corresponding to t = 0 and t = 0 are qualitative and quantitative, respectively.It is a typical qualitative and quantitative mixed multi-response spatio-temporal data set.The surface graphs and scatter plots of C 1 (x), C 2 (x) and C 3 (x) are drawn by using the transformed simulation data.It is not difficult to see from Figure 4 that the test function (26) has a highly nonlinear input-output relationship.
The simulation data in Table A1 of the online supplemental data are modeled and optimized using the meta-modeling and parameter optimization methods proposed in Section 3. The optimization results are shown in Table 2. Comparing the RMSE and quality loss results under the UC, MC, AGP and MMR-LVGP models in Table 2, the proposed method can obtain the minimum RMSE and quality loss, which are 0.4977 and 0.1515, respectively.Therefore, the proposed method has higher prediction    accuracy.To verify further the effectiveness of the extrapolation of the prediction error result of the proposed method, the above simulation is repeated 30 times.The 30 RMSEs under the method UC, MC, AGP and MMR-LVGP are shown in Figure 5, and their average values are 0.6835, 0.7121, 0.7697 and 0.5857.The RMSE box plots and their average values under different methods show that the proposed model has a better predictive performance.
In addition to high prediction accuracy, the MMR-LVGP model can also visualize the influence of each level of the qualitative factor on the response.It can be seen from the definition of the simulation function and the simulation data in Table A1 that, when the qualitative factor t takes the third level, the corresponding response is qualitative, and the responses corresponding to other levels are quantitative.As shown in Figure 6, the third level of qualitative factor t is far from the other levels.Therefore, the ranking and relative distance of each level of qualitative factors are very close to their actual levels, which is convenient for an in-depth understanding of their impact on response.

Ion implantation process experiment
In this section, the ion implantation process described in the introduction (Section 1) is used as an example to illustrate the effectiveness of the proposed method in dealing with the modeling and RPD problems of computer experiments with qualitative and quantitative mixed input and output.The specific settings and levels of control factors in this experiment are listed in Table 3.According to  engineering knowledge, the defect grade of the sensitive area after ion implantation often affects the uniformity of the implanted ion.In order to simplify the analysis, the defects of the sensitive area are divided into five grades I, II, III, IV and V.During the ion implantation process, 36 sensitive areas were assigned to each wafer.Among the 36 sensitive areas, the grade that obtains the largest count is used as the final grade of ion uniformity.The detailed classification definition of ion uniformity is given in Table 4.The experimental data of the Taguchi L 18 orthogonal array are shown in Table B1 of the online supplemental data.The experimental data are modeled and optimized according to the steps given in Section 3.2.The optimization results under different methods are shown in Table 5.The joint optimal solution of the qualitative and quantitative factors obtained by the proposed method is (A,B,C,D,E,F) opt = (Type2, 6.2025, 100.0897, 12.9796, 4.2002, 73.4114), and the prediction value corresponding to the optimal solution is (ŷ 1 , ŷ2 ) opt = (II, 1000).The predicted value of uniformity of the optimal solution is 'good', and the predicted ion implantation amount is just equal to its target value of 1000.In addition, the proposed method obtains a minimum quality loss value of 0.8666 and a minimum RMSE value of 38.1396.Therefore, the proposed method can obtain more robust optimization results compared with the existing methods.

Electric alternator experiment
This section uses the alternator experiment given by Nair, Taam, and Ye (2002) to show that the proposed method can solve the modeling and optimization problem of spatio-temporal data.The specific settings of time and spatial factors are shown in Table 6.The interesting response is the current produced by the alternator at seven different rotational speeds.The upper and lower specification limits (U and L) and the target value (T) of the electric current at different rotational speeds are given in Table 7.The experimenter used the Taguchi L 18 orthogonal array design and repeated the experiment six times.Here, the mean of the six observations at each combination of the factor levels is used as the final experimental result.The experimental data are presented in Table C1 of the online supplemental data.The loss weights of electric current at seven different rotational speeds are w = (0.2654, 0.2280, 0.1594, 0.1237, 0.0860, 0.0706, 0.0669), and the optimization results under different methods are given in Table 8.
As shown in Table 8, under the proposed method the joint optimal solution of the spatial factor and the time factor is (x opt , s opt ) = (−0.20,0.34, 0.32, 0.10, 0.52, −0.43, 0.58, 0.88, 5000).The RMSE, quality loss and deviation between the predicted value and the target value of the optimal solution obtained by the proposed method are 2.2334, 128.9871, and zero, respectively.The proposed model has the smallest values under these comparison criteria compared with the comparison models.Neither the GP nor the MGP methods can obtain the optimal solution of the time factor, and the prediction deviation and quality loss of the two methods are significantly larger than those of a 'PD' represents the deviation between the predicted value and the target value of the optimal solution.b 'QL' represents the quality loss of the predicted response value corresponding to the optimal setting of input factors.c '×' indicates that the corresponding method cannot optimize the time factor s.
the other methods.The prediction value of the optimal solution under the UC method exceeded its lower specification limit, and its quality loss is as high as 488.7830.The predicted value of the optimal solution of the MC method is the same as the target value, but the quality loss is 4.72 times that of the proposed method.Although the quality loss of the AGP method is almost the same as that of the proposed method, it has a higher prediction deviation.In summary, the proposed method can effectively deal with the modeling and optimization problems of spatio-temporal data and obtain more robust optimization results.

Conclusions
All input factors and output responses are assumed to be quantitative in the standard statistical framework of GP modeling in RPD .However, the coexistence of qualitative and quantitative mixed inputs and outputs frequently occurs in the manufacturing industry.Therefore, it is necessary to develop effective meta-modeling and optimization techniques for computer experiments with qualitative and quantitative mixed input and output.This article combines the LVGP model and fuzzy set theory to develop a data-driven model involving qualitative and quantitative mixed input and output.Then, an optimization scheme is established by comprehensively balancing the location and dispersion effects to find the jointly optimal solution of qualitative factors and quantitative factors.The effectiveness of the proposed method in dealing with the modeling and optimization of experimental data with qualitative and quantitative mixed input and output or spatio-temporal structure is verified by a numerical example and two industrial cases.A comparison of the results indicates that the proposed method can obtain more robust optimization results while reducing the computational cost.
Although the proposed method is useful for computer experiments with qualitative and quantitative mixed input and output, it has some drawbacks that need to be improved.For example, owing to the limitations of experimental conditions and statistical techniques, there is a certain degree of deviation between the true response values from physical experiments and the predicted response value based on the model under the same input parameter setting.As a result, there is some bias in comparing the merits of different meta-models in robust parameter design based on the quality loss calculated from the predicted response of the model.In addition, the proposed method mainly focuses on experimental data involving exact numbers, fuzzy numbers or fuzzy linguistics.However, progress in sensing technology has led to most sensing systems equipping manufacturing processes being able to obtain rich data streams of multiple types, having high speeds, being high-dimensional and having complex structures.The focus of future research will be how to design efficient meta-modeling and RPD methods for such high-dimensional complex data.

Disclosure statement
No potential conflict of interest was reported by the authors.

Figure 2 .
Figure 2. Flowchart of the proposed method.

Figure 3 .
Figure 3. Distribution of design points and augmented points of the quantitative factors.

Figure 4 .
Figure 4.The relationship between C 1 , C 2 , C 3 and the quantitative input variables x 1 , x 2 .(a) Surface plot of C 1 (x).(b) Surface plot of C 2 (x).(c) Surface plot of C 3 (x).

Figure 5 .
Figure 5.Comparison of RMSEs for different methods.

Figure 6 .
Figure 6.Latent space of qualitative factor t.

Table 1 .
The conversion relations between linguistic fuzzy numbers and triangular fuzzy numbers.Mapping from the high-dimensional underlying variables to the 2D latent variables.

Table 2 .
Optimization results of different methods in the simulation experiment.
b 'QL' represents the quality loss of the true response value corresponding to the optimal setting of the input factors.

Table 3 .
Control factors and their levels.

Table 4 .
Categories definition of ion uniformity.

Table 5 .
Optimization results of different methods in the ion implantation process experiment.'QL' represents the quality loss of the predicted response value corresponding to the optimal setting of the input factors. a

Table 6 .
List of spatial factors and time factors s.

Table 7 .
Specification limits and target values at different rotational speeds.

Table 8 .
Optimization results of different research methods in the alternator experiment.