2.2. Moving Least Squares Method (MLSM)
In the traditional least squares method (LSM), the response surface function is formed as:
where
is a vector of
polynomial basis functions, and
is an evaluated vector of the coefficient. It is noted that these evaluated coefficients are constant with global regression, while
changes. That is to say, the whole sample points are considered and weighted equally. However, in practical application, a precise approximated model requires greater weight on the sample points closer to the prediction point. Consequently, a weighted regression method called the moving least squares method (MLSM) is proposed to improve the coefficients
. This implies that the coefficients
become functions of
. And the coefficients
can be calculated by minimizing the residual error between the output value and the response surface value.
To construct this function,
sample points of variables in the design space must be chosen based on optimal Latin hypercube design (OLHD). The output value can be calculated from these OLHD sample points using computer simulation, such as the finite element method (FEM). Once the required sample data was obtained, the polynomial constructed for all sample points can be shown as:
where
can be defined as:
After that, the experimental residual error between the output value and the response surface value can be expressed as:
where can
be defined as:
In Equation (6), the weight
can be presented as:
where
is the Euclid distance between the prediction point and the sample point, and
is the influence range whose appropriate size should be selected among a sufficient number of neighboring sample points to avoid singularity in the solution. It is useful to remember that
should contain at least
sample points, and the weight function should vanish outside the
influence range. And
is the coefficient of the exponential and the Gaussian weight function. Among all the weight functions, the Gaussian function is considered to be an effective weight by most scholars for MLSM [
30,
31].
Figure 2 shows the shape of different weight functions using Equation (8). From
Figure 2, the exponential weight function and the Gaussian weight function may have a better weight effect both in the domain range close to the sample point and away from the sample point. Nevertheless, using the exponential weight function, the weight value is not zero near the boundary of the domain area. In this regard, the Gaussian weight function may have a more suitable performance for the MLSM compared to the exponential weight function.
Figure 3 demonstrates the shape of the Gaussian weight function under different coefficients
. Specifically, the value of this coefficient suggested
has a nice distribution characteristic from
Figure 3 [
30].
The minimum of the experimental residual error
can be gained by the partial derivatives with respect to coefficients, as below:
Then, the evaluated coefficient can be achieved from Equation (9):
Therefore, the contribution of the sample points will be decreased as the distance between the prediction point and is increased. This implies the MLSM can denote both a local and global approximation over the entire region.
2.4. Progressive Trigonometric Mixed Response Surface Method
Usually, most RSM models in the literature take second-order polynomials as their basic form. TMRSM, though, has successfully improved the accuracy of the surrogate model, and the determination of the higher-order polynomial should be tested in advance for the given sample data. To circumvent this deficiency, a progressive trigonometric mixed response surface method (PTMRSM) is proposed in this paper to approximate these sample points and test the necessity of the higher-order terms using a t-test.
According to Equations (4) and (10), the experimental residual error can be written as:
If
, then the residual error can be abbreviated to:
Hence, the mathematical expectation and the variance can be expressed as:
where
is the expectation of the evaluated coefficient
. For the MLSM, the expectation of the estimated coefficient is the actual coefficient for the sample points. If
, the variance of
is defined as
. In which the variance of
is indicated as
.
Furthermore, the expectation of the sum of the squared residual error is:
If the trace is simplified as
, there is the variance of
:
Therefore, the unbiased estimation of
is:
The squared residual error can comply with
, and then it can be written as
. Hence, this squared residual error obeys the Chi-square distribution:
According to Equations (18) and (19), the summation of the square residual error should abide by
:
From Equation (10), the mathematical expectation of this evaluated coefficient can be expressed as:
So, the variance of the evaluated coefficient can be indicated as:
Make
, then we have
. And a standard normal distribution can be achieved:
Combining Equations (20) and (23), the t distribution can be represented as:
where is a
t-distribution of the
degree of freedom.
Determining the significance of the higher-order terms means measuring the importance of these evaluated coefficients. This indicates that the insignificant coefficients stand for the unimportant higher-order terms. Equation (10) shows the vector of the coefficients, and it is easy to find the coefficient values
of the highest order. Then, the focus of this work becomes to check whether the highest coefficients have significant effects on the output response. It is equivalent to whether the highest order coefficients
are equal to the values of 0 (0 shows the insignificant impact). Assume a null hypothesis
, and the test can be performed by calculating the
t-test value.
where
indicates t distribution of
degree of the freedom (see Equation (24));
is the standard deviation of the coefficients;
is the number of the terms of the model, which can be selected from the first number of the highest-order terms to the last number of the highest-order terms.
For this
t-test, a two-sided test and 90% confidence intervals are usually chosen in this paper [
34]. Then, the boundary value
can be calculated by
t-distribution functions. If
, the null hypothesis
can be rejected. This implies these coefficients
and the relevant high-order terms may have a significant influence. Otherwise, the null hypothesis
cannot be dismissed. That is to say, the coefficients
and these high-order terms may not have a substantial effect.
In addition to the
t-test criterion, other criteria should exist to ensure the relationships in which the odd or even order terms have a significant effect on the response results. Only testing whether the highest order has an important effect may easily lead to inaccurate approximations of these special functional relationships. To avoid this deficiency, the determination of the coefficient
and the mean relative error
should be introduced for the
sample points, respectively.
where
represents the sum of squares for error, which is caused by experimental error. And
expresses the sum of squares for the total, representing the total variation of values of
sample observations. These two values can be shown as:
where
is the mean value of output response
, which can be expressed as:
The determination of the coefficient is a statistical indicator that reflects the reliability of the response surface models to account for changes in terms of dependent variables. The closer this value is to 1, the better the approximated effect will be. If the determination of the coefficient
for the m sample points, it means that the approximation degree is sufficient [
38].
The mean relative error is another indicator that can illustrate the relative error between the actual response and the approximate response. The lesser the value
is, the smaller the error is. For most practical engineering problems, it can be viewed as an appropriate limit when
[
38].
Eventually, another stopping criterion is the highest order limit. Given the possibility of over-fitting for the approximate functions, this model takes a stopping criterion with respect to the sixth highest order. That is due to the fact that the highest terms are really rare to exceed the fifth order for most engineering problems [
34].
The proposed PTMRSM can be started with the order
, and the program can determine whether the
t-test values of the third-order terms are sufficient for the 90% confidence intervals. If not, the determination of the coefficient and the mean relative error will be utilized to check the reliability of the fitting performance and the relative error. If these criteria are satisfied, stop this process, and output the response. If not, continue this process until the highest order reaches the sixth order. Then, the response of PTMRSM will be outputted, displaying that “Warning: the highest order has reached the maximum value (6th order)”.
Figure 4 shows the step-by-step procedure of this algorithm:
Perform the OLHD sample points, and obtain the response with respect to these sample points based on FEM or specific functions;
Transform all the sample points and the input parameters into the trigonometric functions and . Then, choose order as the beginning, and construct this third-order TMRSM;
Calculate the values of the t-test criterion . Meanwhile, the determination of the coefficient and the mean relative error can also be calculated from Equations (28) and (29);
Check the t-test criterion and decide whether this operation should be stopped. If the coefficients , it is concluded that the highest-order terms are insignificant. Then, stop the process, and construct the model with (k − 1)-th order TMRSM;
If the t-test criterion is not satisfied, check the criteria for the determination of the coefficient and the mean relative error . If these criteria are fulfilled, it implies the performance and the accuracy of the approximated degree are sufficient. Then, stop this operation and build the model with k-th order TMRSM;
Otherwise, select (k + 1)-th order as the next step order, and check whether the highest order is satisfied. If not, repeat steps 3 to step 6 until one of these criteria is satisfied;
When the highest order is fulfilled, establish the model with the sixth order TMRSM, “Warning: the highest order number has reached the maximum value (6th order)”.
The pseudocode for the PTMRSM process is listed in the Algorithm 1.
Algorithm 1: The pseudocode for the PTMRSM process: |
1: | Read input variable matrix and output response matrix: |
2: | load Sample.mat; |
3: | Set the confidence level of t-test criterion, the determination coefficient criterion and |
4: | the mean relative error criterion. |
5: | Transform input variables and sample matrix into trigonometric functions sin θ and cos θ. |
6: | for sn = 1:Number of response column |
7: | for k = 3:6 |
8: | Construct the X matrix with MLSM using k order |
9: | Evaluated coefficient: b = (P.′ ∗ W ∗ P)\(P.′ ∗ W ∗ Y); |
10: | Calculate errors between actual response values and predicted values |
11: | Calculate the value of t-test criterion: |
12: | H= P *(P.′ ∗ W ∗ P)\(P.′ ∗ W); |
13: | s = trace((eye(m)-H) ∗ (eye(m)-H).’); |
14: | SSE = e.’ ∗ e; %e is the error between actual value and evaluated value |
15: | C_jj = diag((P.′ ∗ W ∗ P)\(P.′ ∗ W ∗ W ∗ P)/(P.′ ∗ W ∗ P)); |
16: | s_d = (SSE./s). ∗ C_jj; |
17: | tj = b./(s_d).^0.5; |
18: | Calculate the determination coefficient criterion value: R2 = 1-SSE./SST; |
19: | Calculate the mean relative error criterion value: MSE = mean(abs(e./Y)); |
20: | if (t < =tinv(confidence level,s) ∗ ones(highest order number,1)) |
21: | k = k − 1; |
22: | break |
23: | end |
24: | if (R2 > = 0.99)&(MSE < = 0.02) |
25: | k = k; |
26: | break |
27: | end |
28: | end |
29: | Adopt the new order and reconstruct the TMRSM model using updated order |
30: | Output the response value: |
31: | f (sn) = p ∗ b; |
32: | h_order(sn) = k; |
33: | sn = sn + 1; |
34: | end |