# Improved Hybrid Response Surface Method Based on Double Weighted Regression and Vector Projection

Authors: Yu Xia; Yeming Wang
Journal: Mathematical Problems in Engineering
 (2022)

---

## Abstract

In order to increase the accuracy and stability of the classical response surface method and relevant method, a new improved response surface method based on the idea of double weighting factors and vector projection method is proposed. The response surface is fitted by the weighted regression technique, which allows the sampling points to be weighted by their distance from the true failure surface and that from the estimated design point. It uses the vector of the gradient projection method to get new sampling points in the process of iteration, in order to make the sampling points closer to the design point, and the value of deviation coefficient is constantly adjusted. To some extent, these strategies increase the accuracy and stability of the response surface method, while the calculation time is decreased. At last, the rationality and efficiency of the proposed method are demonstrated through five examples. Besides, as revealed from this investigation, compared with other conventional algorithms, this method has a few obvious advantages; this algorithm not only has high precision and efficiency, but also has solid stability.


---

## Body

## 1. Introduction

In the research of engineering structural reliability, common methods include first-order reliability method, Monte Carlo method, and response surface method. The first-order reliability method has many applications in engineering practice because it is easy to understand and put in use [1–3]. For the first-order reliability method, the problem of robustness and efficiency has always been the core problem of the method. In order to improve its robustness and efficiency, researchers have put forward their own ideas [4–8]. However, in general, the drawbacks of these methods have not been effectively resolved. In other words, they are very suitable in some cases of high nonlinearity, but in other cases, they may not even converge. And Monte Carlo method is recognized as a general method for structural reliability, but its efficiency is very low. Therefore, the research on Monte Carlo method is mainly on how to reduce its calculation amount [9, 10].In summary, the accuracy of the first order reliability method is low, and the stability is poor; the Monte Carlo method is too computationally expensive. The response surface method can effectively overcome the disadvantages of the above two methods. Because the response surface method is simple in principle, efficient in calculation, and easy to operate, scholars in the field of reliability have paid extensive attention to the response surface method [11–13]. At present, many scholars have improved the response surface method in the two aspects of weighted regression and vector projection.As early as 2005, Kaymaz [14] et al. proposed the idea of weighted regression and applied it to the response surface method, which triggered a large number of scholars to enable further exploration of their research results. Subsequently, Chen et al. [15] constructed an improved response surface method based on weighted regression. Xiong et al. [16] proposed a weighted stochastic response surface method that considered the sample probabilistic weights in regression. Dong et al. [17] proposed an improved weighted regression distributed collaborative surrogate model method. Lu et al. [18] proposed a weighted regression-based extremum response surface method. Recently, Li et al. [19] introduced the idea of weighted regression in the quasisparse response surface method, and the problem of overfitting in the process of atom selection was reduced. The above research results are all improvements to the weighted response surface method that contains only a single weighting factor.Inspired by the idea of weights, Nguyen et al. [20] added a weighting factor (the distance from the fitting points to the design point) based on a single weighting factor to improve the single weighted response surface method proposed by Kaymaz and then proposed a double weighted response surface method. Xiong et al. [21] proposed a new double weighted stochastic response surface method. Guimarães et al. [22] used a cross-validation method and double weighted regression technique to improve the accuracy and efficiency of response surface method. These studies can be regarded as an extension of the single weighted response surface method. And the single weighted response surface improvement method series and the double weighted response surface improvement method series belong to the use of weighted regression ideas to improve the response surface method.In the field of vector projection, Kim [23] is one of the first scholars to apply the vector projection method to the response surface method, and the vector projection response surface method is applied to practical engineering applications by many scholars due to its high efficiency and accuracy [24, 25]. The vector projection response surface method proposed by Basaga et al. [26] searches for points close to the true limit state function around the design point and then uses the Second-Order Reliability Method to obtain the reliability index and probability of failure. Li et al. [27] proposed a new class of hybrid perturbation-Galerkin method based on the combination of perturbation technique and Bubnov-Galerkin projection. The core of the above research is the improvement and application of the vector projection response surface method.All the above research results are based on the weighted regression or the vector projection method to improve the accuracy and efficiency of the response surface method, but there are few researches that combine the two to study together. The most representative of these is the improved response surface method developed by Zhang [28]. However, in the weights setting of the sampling points, Zhang only considered the distance between the sampling points and the true failure surface and did not take into account that the distance between the sampling points, and the design point will also affect the model fitting accuracy of the final response surface. In addition, in its iterative process, although the linear interpolation method can achieve higher accuracy, it will increase the amount of calculation to a certain extent.To address the above issues, on the basis of Zhang’s research [28], this paper further applies vector projection method and weighted regression to the classical response surface method and then adds a weight factor considering the distance between the sampling points and the design point. The overall weight takes the average weight of the two weights. Besides, this paper no longer simply uses the small step size method of linear interpolation to search for the design point because this way will cause a lot of duplication of work. According to the characteristics of the design point, this paper makes multiple iterations. At the beginning, the first iteration was performed in the case of a linear response surface combined with the negative gradient direction. Starting from the second iteration, the quadratic polynomial function without cross terms is combined with the negative gradient direction. The effect of this is that it can quickly lean to the true failure surface. Then, iterate alternately with the linear interpolation method of small step size and the negative gradient direction of large step size, so that the calculation efficiency and accuracy are improved. This method solves the problem that the response surface series method has always been affected by the deviation coefficient change, which makes low stability, and the method in this paper further improves the calculation accuracy of the response surface method while ensuring the efficiency. Finally, the rationality of the method in this paper is verified by the relatively accurate solution of the Monte Carlo method and five examples.


## 2. Improved Weighted Response Surface Method

### 2.1. Choice of Function Form

There are many choices for the form of response surface functions, such as neural networks, splines, and polynomials. However, polynomials are the most commonly used, especially simple linear polynomials and quadratic polynomials without cross terms in the weighted response surface method.The expression of the linear response surface is as follows:(1)Zr=g¯X=a0+∑i=1naixi.Among them,a0 and ai represent unknown coefficients, and xi represents random variables.The expression of the quadratic polynomial response surface function without cross terms is as follows:(2)Zr=g¯X=a0+∑i=1naixi+∑i=1nbixi2.Among them,a0, ai and bi all represent unknown coefficients.


### 2.2. Determination of Function Coefficients

The coefficients of the classical response surface method are obtained by the following formula:(3)A=B′B−1B′y^,where B represents the design matrix at sampling points, and y^ represents the structural response obtained at sampling points. From the above formula, it can be found that all sampling points are equalized, and the weights are exactly the same, which directly leads to the low degree of the fitting of the response surface in our classical response surface method. And it further explains the reason for the lower accuracy of the classical response surface.In order to construct a response surface closer to the true failure surface, Kaymaz proposed a weighted regression idea. The coefficient matrix of the response surface function proposed by him is constructed as follows:(4)A=B′CB−1B′Cy^,where C is the weight matrix.


### 2.3. Weights of Sampling Points

#### 2.3.1. The First Weight

According to Kaymaz’s research [14], each sampling point has a different effect on the response surface function in fitting the true failure surface. The closer the point is to the true failure surface, the greater the effect on the fitting is. Therefore, the value of the limit function at sampling point describes the distance from sampling point to the true failure surface in theory. It can be used as the first weight. And the smaller the absolute value of the limit function at sampling point, the greater the effect of this point. Based on this, Kaymaz [14] gives a form of weight function as follows:(5)Cj=exp−Gxj−GminGmin,where j represents the j-th sampling point, and Gmin represents the response value of the sampling point closest to the true failure surface.Although this weighting form greatly improves the accuracy and efficiency of the classical response surface method, when the minimum value ofG is close to 0, the ill-conditioned problem of the regression matrix often occurs, which limits its application.In order to solve this problem, Wu [29] proposed a better form of weight function:(6)Cj=exp−θGxj−Gx0Gx0,where Gx0 represents the response value of the limit state function at the origin, and θ is a non-zero constant used to adjust the smoothness of the fitted curve.This formula solves the ill-conditioned problem of the regression matrix very well, but the randomness of the valueθ is its obvious shortcoming. Based on the ideas of the above two scholars, this paper constructs a more stable weight form that will not cause matrix ill-conditioned problems. When Gxj=0, Cj=1.2 . This is because the closer the test point is to the limit state surface, the greater its weight is.The weight of the test point atGmin is 1. In theory, except for the sampling points on the true failure surface, it is the point closest to the true failure surface. Therefore, the weight of the sampling points on the true failure surface should be greater than that at Gmin, and the weight value of the sampling points on the limit state surface is 1.2 here.(7)Cj=GminGxj,Gxj≠0,1.2,Gxj=0.


#### 2.3.2. The Second Weight

Considering that the closer the sampling point is to the design point, the more obvious its fitting effect will be. Following the functional form in [22], the expression for constructing the second weight is as follows:(8)dj=x∗−xj,dmax=maxdj,fi=djdmax,C2j=1−fj23,fj≤1,0,fj>1,where x∗ represents the design point, and xj represents the sampling points.


## 2.1. Choice of Function Form

There are many choices for the form of response surface functions, such as neural networks, splines, and polynomials. However, polynomials are the most commonly used, especially simple linear polynomials and quadratic polynomials without cross terms in the weighted response surface method.The expression of the linear response surface is as follows:(1)Zr=g¯X=a0+∑i=1naixi.Among them,a0 and ai represent unknown coefficients, and xi represents random variables.The expression of the quadratic polynomial response surface function without cross terms is as follows:(2)Zr=g¯X=a0+∑i=1naixi+∑i=1nbixi2.Among them,a0, ai and bi all represent unknown coefficients.


## 2.2. Determination of Function Coefficients

The coefficients of the classical response surface method are obtained by the following formula:(3)A=B′B−1B′y^,where B represents the design matrix at sampling points, and y^ represents the structural response obtained at sampling points. From the above formula, it can be found that all sampling points are equalized, and the weights are exactly the same, which directly leads to the low degree of the fitting of the response surface in our classical response surface method. And it further explains the reason for the lower accuracy of the classical response surface.In order to construct a response surface closer to the true failure surface, Kaymaz proposed a weighted regression idea. The coefficient matrix of the response surface function proposed by him is constructed as follows:(4)A=B′CB−1B′Cy^,where C is the weight matrix.


## 2.3. Weights of Sampling Points

### 2.3.1. The First Weight

According to Kaymaz’s research [14], each sampling point has a different effect on the response surface function in fitting the true failure surface. The closer the point is to the true failure surface, the greater the effect on the fitting is. Therefore, the value of the limit function at sampling point describes the distance from sampling point to the true failure surface in theory. It can be used as the first weight. And the smaller the absolute value of the limit function at sampling point, the greater the effect of this point. Based on this, Kaymaz [14] gives a form of weight function as follows:(5)Cj=exp−Gxj−GminGmin,where j represents the j-th sampling point, and Gmin represents the response value of the sampling point closest to the true failure surface.Although this weighting form greatly improves the accuracy and efficiency of the classical response surface method, when the minimum value ofG is close to 0, the ill-conditioned problem of the regression matrix often occurs, which limits its application.In order to solve this problem, Wu [29] proposed a better form of weight function:(6)Cj=exp−θGxj−Gx0Gx0,where Gx0 represents the response value of the limit state function at the origin, and θ is a non-zero constant used to adjust the smoothness of the fitted curve.This formula solves the ill-conditioned problem of the regression matrix very well, but the randomness of the valueθ is its obvious shortcoming. Based on the ideas of the above two scholars, this paper constructs a more stable weight form that will not cause matrix ill-conditioned problems. When Gxj=0, Cj=1.2 . This is because the closer the test point is to the limit state surface, the greater its weight is.The weight of the test point atGmin is 1. In theory, except for the sampling points on the true failure surface, it is the point closest to the true failure surface. Therefore, the weight of the sampling points on the true failure surface should be greater than that at Gmin, and the weight value of the sampling points on the limit state surface is 1.2 here.(7)Cj=GminGxj,Gxj≠0,1.2,Gxj=0.


### 2.3.2. The Second Weight

Considering that the closer the sampling point is to the design point, the more obvious its fitting effect will be. Following the functional form in [22], the expression for constructing the second weight is as follows:(8)dj=x∗−xj,dmax=maxdj,fi=djdmax,C2j=1−fj23,fj≤1,0,fj>1,where x∗ represents the design point, and xj represents the sampling points.


## 2.3.1. The First Weight

According to Kaymaz’s research [14], each sampling point has a different effect on the response surface function in fitting the true failure surface. The closer the point is to the true failure surface, the greater the effect on the fitting is. Therefore, the value of the limit function at sampling point describes the distance from sampling point to the true failure surface in theory. It can be used as the first weight. And the smaller the absolute value of the limit function at sampling point, the greater the effect of this point. Based on this, Kaymaz [14] gives a form of weight function as follows:(5)Cj=exp−Gxj−GminGmin,where j represents the j-th sampling point, and Gmin represents the response value of the sampling point closest to the true failure surface.Although this weighting form greatly improves the accuracy and efficiency of the classical response surface method, when the minimum value ofG is close to 0, the ill-conditioned problem of the regression matrix often occurs, which limits its application.In order to solve this problem, Wu [29] proposed a better form of weight function:(6)Cj=exp−θGxj−Gx0Gx0,where Gx0 represents the response value of the limit state function at the origin, and θ is a non-zero constant used to adjust the smoothness of the fitted curve.This formula solves the ill-conditioned problem of the regression matrix very well, but the randomness of the valueθ is its obvious shortcoming. Based on the ideas of the above two scholars, this paper constructs a more stable weight form that will not cause matrix ill-conditioned problems. When Gxj=0, Cj=1.2 . This is because the closer the test point is to the limit state surface, the greater its weight is.The weight of the test point atGmin is 1. In theory, except for the sampling points on the true failure surface, it is the point closest to the true failure surface. Therefore, the weight of the sampling points on the true failure surface should be greater than that at Gmin, and the weight value of the sampling points on the limit state surface is 1.2 here.(7)Cj=GminGxj,Gxj≠0,1.2,Gxj=0.


## 2.3.2. The Second Weight

Considering that the closer the sampling point is to the design point, the more obvious its fitting effect will be. Following the functional form in [22], the expression for constructing the second weight is as follows:(8)dj=x∗−xj,dmax=maxdj,fi=djdmax,C2j=1−fj23,fj≤1,0,fj>1,where x∗ represents the design point, and xj represents the sampling points.


## 3. Vector Projection Method

Assuming that the design point obtained last time isx∗, the expression of the response surface function is (1). Calculate the unit vector αr of the response surface function at the design point:(9)αr=−∇g¯x∗∇g¯x∗.Thei-th random variable is projected onto the response surface, and the expression of the projection vector is as follows:(10)Ti=ei−αrTeiαr,i=1,2,…,n,ti=TiTi,Ti≠0,0,Ti=0,where ei represents the unit basis vector along the coordinate axis Xi; that is, the i-th component of ei is 1, and the others are 0.In order to project the sampling points as close to the response surface as possible and make the sampling points close to the true failure surface, the expression for the projection direction is shown as follows:(11)Qi=ti+εαr,qi=QiQi,where ε represents a small constant.The coordinates of the sampling points on theXi axis finally obtained become(12)xi=xi∗±fσXiqiTei,where f represents the deviation coefficient.Figure1 shows a schematic diagram of sampling points for selection using vector projection technology in a two-dimensional situation.Figure 1
Sampling points selected by vector projection.


## 4. The Method Proposed of This Paper

Based on the classic response surface method, this paper proposes the following five improvement measures: ① A more reasonable weight form is proposed, so that the method of this paper will not cause the ill-conditioned problem of the matrix and is more stable. In particular, considering that the weight coefficient should be larger when the sampling point is on the failure surface, Therefore, the weight of the sampling points on the true failure surface should be greater than that atGmin, and their weight value on the limit state surface is 1.2 here. ② This paper proposes a strategy for searching design points in steps. The specific operation is as follows: according to the characteristics of the design point, the first iteration was performed in the case of a linear response surface combined with the negative gradient direction. Starting from the second iteration, the quadratic polynomial function without cross terms is combined with the negative gradient direction. The effect of this is that it can quickly lean to the true failure surface. Then, iterate alternately with the linear interpolation method of small step size and the negative gradient direction of large step size, so that the calculation efficiency and accuracy are improved. ③ In this paper, the vector projection method is used to optimize the selection of sampling points, so that the sampling points are more likely to fall near the failure surface, so as to improve the accuracy of the classical response surface method; ④ Introduce a new weight factor, so that the sampling points close to the design point have greater weight, thereby achieving the effect of improving the fitting degree of the response surface; ⑤ In the iterative process, the value of the deviation coefficient can be adjusted adaptively, so that the distribution range of the sampling points is first large and then small, and finally, it can stably converge near the design point, enhancing the robustness of the classical response surface method. Finally, the algorithm steps and flowchart of this method are as follows.


### 4.1. Steps of Algorithm

(i)
Step 1: Assuming that the initial design point is the mean point, construct the remaining 2n axial sampling points using the traditional experimental design method, and then use the traditional linear response surface method to obtain the next design point, where the form of the linear response surface function is shown in (1).(ii)
Step 2: The design point obtained in the previous step is the center point, and the corresponding sampling points is obtained using the central composite design method with only the axis point, the corresponding coefficient matrixA is constructed, and then the finite element method or numerical calculation method is used to get estimated value of the performance function at sampling points.(iii)
Step 3: According to the two weighting factors of the distance between the sampling point and the true failure surface and the distance between the sampling point and the design point, construct the corresponding double weighted matrix. The weight is the average weight of the two weights. Then, use the response surface function without cross terms to perform a weighted regression analysis to determine the undetermined coefficients of the response surface function.(iv)
Step 4: Use the First-Order Reliability Method to calculate the reliability index, and perform the iteration of the design point in the negative gradient direction of the large step size. The iteration form is as follows:(13)xi∗=μXi−∇gx.β.σi2∑i=1n∇gx2⋅σXi2,i=1,2,…,n.(v)
Step 5: Use the vector projection method to construct the sampling points for the next weighted least squares regression analysis, and, at the same time, raise the sixth power of the deviation coefficientf.(vi)
Step 6: Repeat steps 2 to 5 until the distance between the two design points meets a certain accuracy, and then proceed to the next step after the iteration converges.(vii)
Step 7: Perform linear interpolation between the starting pointx1,gx1 and the design point xi∗,gxi∗ obtained in the previous iteration, and then get the next design point xi+1∗.(14)xi+1∗=x1+gx1gx1−gxi∗⋅xi∗−x1,(viii)
Step 8: Construct a weighting matrix according to the idea of double weighting factors, perform regression analysis, and then determine the unknown coefficients of the response surface.(ix)
Step 9: Use the First Order Reliability Method to calculate the reliability index and get the next design point.(x)
Step 10: Perform linear interpolation from the mean pointμx,gμx and the design point xi∗,gxi∗ obtained in the previous iteration.(15)xi+1∗=μx+gμxgμx−gxi∗⋅xi∗−μx,(xi)
Step 11: Use the vector projection method to construct new sampling points.(xii)
Step 12: Repeat steps 8 to 11 until the reliability index meets a certain accuracy.The calculation flow chart of the method in this paper is shown in Figure2.Figure 2
Flow chart of this proposed method.


## 4.1. Steps of Algorithm

(i)
Step 1: Assuming that the initial design point is the mean point, construct the remaining 2n axial sampling points using the traditional experimental design method, and then use the traditional linear response surface method to obtain the next design point, where the form of the linear response surface function is shown in (1).(ii)
Step 2: The design point obtained in the previous step is the center point, and the corresponding sampling points is obtained using the central composite design method with only the axis point, the corresponding coefficient matrixA is constructed, and then the finite element method or numerical calculation method is used to get estimated value of the performance function at sampling points.(iii)
Step 3: According to the two weighting factors of the distance between the sampling point and the true failure surface and the distance between the sampling point and the design point, construct the corresponding double weighted matrix. The weight is the average weight of the two weights. Then, use the response surface function without cross terms to perform a weighted regression analysis to determine the undetermined coefficients of the response surface function.(iv)
Step 4: Use the First-Order Reliability Method to calculate the reliability index, and perform the iteration of the design point in the negative gradient direction of the large step size. The iteration form is as follows:(13)xi∗=μXi−∇gx.β.σi2∑i=1n∇gx2⋅σXi2,i=1,2,…,n.(v)
Step 5: Use the vector projection method to construct the sampling points for the next weighted least squares regression analysis, and, at the same time, raise the sixth power of the deviation coefficientf.(vi)
Step 6: Repeat steps 2 to 5 until the distance between the two design points meets a certain accuracy, and then proceed to the next step after the iteration converges.(vii)
Step 7: Perform linear interpolation between the starting pointx1,gx1 and the design point xi∗,gxi∗ obtained in the previous iteration, and then get the next design point xi+1∗.(14)xi+1∗=x1+gx1gx1−gxi∗⋅xi∗−x1,(viii)
Step 8: Construct a weighting matrix according to the idea of double weighting factors, perform regression analysis, and then determine the unknown coefficients of the response surface.(ix)
Step 9: Use the First Order Reliability Method to calculate the reliability index and get the next design point.(x)
Step 10: Perform linear interpolation from the mean pointμx,gμx and the design point xi∗,gxi∗ obtained in the previous iteration.(15)xi+1∗=μx+gμxgμx−gxi∗⋅xi∗−μx,(xi)
Step 11: Use the vector projection method to construct new sampling points.(xii)
Step 12: Repeat steps 8 to 11 until the reliability index meets a certain accuracy.The calculation flow chart of the method in this paper is shown in Figure2.Figure 2
Flow chart of this proposed method.


## 5. Example Analysis

Five different types of examples are given and calculated by the proposed method (PM) of this paper. The calculation results of this method are compared with those of classical response surface method (RSM), linear weighted response surface method (WRSM), and vector projection response surface method (WPRSM). Besides, the results of Monte Carlo method (MCS) are also used as relatively accurate solutions to verify the effectiveness and applicability of the method in this paper.Among them, the absolute value of the relative error difference between different interpolation coefficients is used as the stability index of the response surface function. The smaller the difference, the higher the stability of the corresponding method.


### 5.1. Example 1: Simple Low Nonlinear Function

The performance function is shown below, wherex and y both obey the standard normal distribution. The calculation results of the above methods are listed in Table 1, and the iterative process is shown in Figure 3.(16)gx,y=exp0.2x+6.2−exp0.47y+5.0.Table 1
Comparison of calculation results of Example 1.
MethodDeviation coefficientfReliability indexβRelative error (%)StabilityIteration times (s)MCS2.3497RSM32.35750.33201.28955.1680e − 0452.38781.62156.3970e − 04WRSM32.35170.08510.48529.1380e − 0452.36310.57031.1853e − 03VPRSM32.34870.04260.05546.3904e − 0452.34740.09801.3155e − 03PM32.34930.017006.1710e − 0452.34930.01706.1850e − 04Figure 3
Iterative process of Example 1.It can be seen from Table1 that the accurate result of the reliable index obtained by MCS is 2.3497. Compared with MCS, PM gives the most accurate results with an error of only 0.0170%. It is worth noting that WRSM and VPRSM not only have a larger relative error than PM, but also require longer iteration time. Although the iteration time of RSM is similar to that of PM, the relative error of its reliability index is the largest. When the deviation coefficient is equal to 5, the iterative process diagram of each method is shown in Figure 1. It can be seen from the figure that the convergence solution of PM is closest to the exact solution of MCS, and PM only needs 7 iterations to reach the convergence solution, with the least number of iterations. This phenomenon is mutually corroborated by the results of the reliability index and relative error shown in Table 1. According to Table 1, the stability index value of the method in this paper is 0, which shows that the reliability index obtained by the method in this paper is almost not affected by the change of the deviation coefficient. At the same time, it can be obtained that the order of the stability of each method from high to low is PM, VPRSM, WRSM, RSM.In summary, in the case of low nonlinear functions, PM is not only superior in accuracy to other method, but also the iteration convergence speed of the method in this paper is faster and higher in efficiency. In addition, PM has high stability and strong robustness.The reason for the above phenomenon is that the PM method uses the strategy of searching the design point in steps. In the iterative process, according to the characteristics of the design point, it is first iterated in the negative gradient direction with a large step size, so that it can approach the design point faster, and then it is alternately iterated in the linear interpolation with a small step size and the negative gradient direction with a large step size afterwards, which makes the computational efficiency increase and the iteration time decrease.


### 5.2. Example 2: Functional Function with Cross terms

Example 2 is to study the influence of the cross term on the algorithm. The calculation results are listed in Table2, and the iteration process is shown in Figure 4.(17)gx,y=570−xy,x∼N1000,33,y∼N2,0.1.Table 2
Comparison of calculation results of Example 2.
MethodDeviation coefficientfReliability indexβRelative error (%)StabilityIteration times (s)MCS2.1257RSM32.12580.00050.05123.8610e − 0452.12680.05173.8490e − 04WRSM32.12560.00050.03713.9210e − 0452.12490.03763.9780ee − 04VPRSM32.08292.01353.53763.4970e − 0452.00775.55113.5350e − 04PM32.12560.000509.4610e − 0452.12560.00059.5260e − 04Figure 4
Iterative process of Example 2.It can be seen from Table2 that although the iteration time of PM is the longest, the reliability index obtained by PM is closest to the solution of MCS, so the accuracy of PM is the highest. Secondly, as the deviation coefficient changes, the stability index of PM is always 0, indicating that the stability of PM is very high. It can be found intuitively from Figure 4 that the reliability index obtained by PM, WRSM, and RSM is all close to the exact solution, when the relative error of the reliability of VPRSM is the largest.Therefore, it can be concluded that when the performance function contains cross terms, the efficiency of this method is slightly lower than other methods, but in terms of accuracy and stability, this method is still superior to other methods.


### 5.3. Example 3: Nonlinear Function of Concave toward Failure Region

The performance function of this example is shown below, and bothx and y obey the standard normal distribution. The calculation results are shown in Table 3, and the iterative process and the fitting effect are shown in Figure 5 and Figure 6, respectively. The fitting effect of the response surface function obtained in the last iteration of each method to the true failure surface at the most probable point is shown in Figure 6. Among them, LSF represents the true failure surface, and MPP represents the most probable point, that is, the design point.(18)gx,y=exp0.2x+1.4−y.Table 3
Comparison of calculation results of Example 3.
MethodDeviation coefficientfReliability indexβRelative error (%)StabilityIteration times (s)MCS3.3811RSM33.35060.91030.18575.1300e − 0453.35660.72465.0850e − 04WRSM33.34990.92280.03068.5460e − 0453.35120.89228.6010e − 04VPRSM33.48403.04346.95935.3020e − 0453.719310.00275.9950e − 04PM33.34970.928705.1780e − 0453.34970.92875.3230e − 04Figure 5
Iterative process of example 3.Figure 6
Fitting effect of each method.From the research of Nguyen, it can be known that the key of the response surface method is that the response surface function fits the sampling points well, especially near the design point. It can be seen from Figure6 that, except for VPRSM, other methods have a relatively good fitting effect on the true failure surface at the design point, while the vector projection response surface is still a certain distance from the design point. The fitting effect at the design point is not very good, so it can be concluded that the accuracy of PM in this example is much higher than that of VPRSM. This conclusion is also verified in column relative error of Table 3. Figure 5 can also draw conclusions similar to Figure 6. In Figure 5, the reliability index of PM, RSM, and WRSM is slightly smaller than the exact solution, while the reliability of VPRSM is much larger than the exact solution, which also shows that the high accuracy of PM.Table3 accurately shows the characteristics of each method from a quantitative perspective. The reliability index of PM, RSM, and WRSM has similar results, but according to the iteration time, PM and RSM are better than WRSM. Although the iteration time of VPRSM is also very small, the relative error of its reliability is very large, even as high as 10.0027%. Judging from the stability index, we can draw conclusions almost the same as those in the first and second examples. The method in this paper has high stability and strong robustness.


### 5.4. Example 4: Nonlinear Function of Convex toward Failure Region

The type of example 4 is just the opposite of that of example 3. It is a convex function, and its performance function is as follows. Bothx and y obey the standard normal distribution. The calculation results are shown in Table 4, and the iteration path diagram of each method is shown in Figure 7.(19)gx,y=exp0.4x+2+6.2−exp0.3y+5.0−200.Table 4
Comparison of calculation results of example 4.
MethodDeviation coefficientfFailure probabilityRelative error (%)StabilityIteration times (s)MCS2.6855RSM32.71130.96070.01495.5980e − 0452.71170.97566.1820e − 04WRSM32.71030.92350.06336.9080e − 0452.71200.98686.9330e − 04VPRSM32.63192.03660.29895.3510e − 0452.51916.06558.5430e − 04PM32.70990.908605.8180e − 0452.70990.90865.9790e − 04Figure 7
Iterative process of each method.It can be seen from Table4 that the relative errors of the PM, RSM, and WRSM methods are very low, and the iteration time is similar. Although the iteration time of VPRSM is similar, its relative error is the largest. At the same time, it can be concluded that the stability of the method in this paper is still the best (0 stability).It can be seen from the figure above that, except for the vector projection response surface method, other methods all converge to the true failure surface. This indicates that the PM is more accurate than the VPRSM. From the iterative path, it can be seen that the method in this paper can quickly approach the convergence point in only 4 steps, which is significantly faster than other methods in efficiency. This is because the method in this paper adopts a strategy of iterating in the negative gradient direction with a large step size first, and then iterating with a small step size.


### 5.5. Example 5: Ten-Bar Truss Structure

In this example, a ten-bar truss structure（Figure8）is selected. This structure is fixed with hinge support at node 1 and node 2, where E is the elastic modulus E=107psi1psi=6.895kPa, L is the length of the horizontal and vertical rods L=360in1in=2.54cm, and P is the vertical load P=100000lb1lb=0.4536kg. And the section areas A1, A2, and A3 are all independent normal random variables A1=N13,1.3A2∼N2,0.2A3=N9,0.9. The failure state with respect to the serviceability of the structure is reached when the horizontal displacement dreal of the node 5 exceeded a threshold value 3.5 in. Obviously, this is a case of an implicit performance function. And the results are shown in Table 5.(20)gA1,A2,A3=dallow−dreal.Figure 8
Ten-bar truss structure.Table 5
Comparison of calculation results of Example 5.
MethodDeviation coefficientfFailure probabilityRelative error (%)StabilityIteration times (s)MCS0.836895RSM30.8199432.02560.00050.015650.8199472.02510.0155WRSM30.8199462.02520.00050.026250.8199502.02470.0248VPRSM30.8634803.17667.23180.019150.92400210.40840.0192PM30.8199422.025700.010850.8199422.02570.0105In this example, the calculation error obtained by VPRSM is too large and is no longer applicable to this example. Regardless of whether it is in the relative error of the failure probability or in the iteration time, the results displayed by the PM are satisfactory.


## 5.1. Example 1: Simple Low Nonlinear Function

The performance function is shown below, wherex and y both obey the standard normal distribution. The calculation results of the above methods are listed in Table 1, and the iterative process is shown in Figure 3.(16)gx,y=exp0.2x+6.2−exp0.47y+5.0.Table 1
Comparison of calculation results of Example 1.
MethodDeviation coefficientfReliability indexβRelative error (%)StabilityIteration times (s)MCS2.3497RSM32.35750.33201.28955.1680e − 0452.38781.62156.3970e − 04WRSM32.35170.08510.48529.1380e − 0452.36310.57031.1853e − 03VPRSM32.34870.04260.05546.3904e − 0452.34740.09801.3155e − 03PM32.34930.017006.1710e − 0452.34930.01706.1850e − 04Figure 3
Iterative process of Example 1.It can be seen from Table1 that the accurate result of the reliable index obtained by MCS is 2.3497. Compared with MCS, PM gives the most accurate results with an error of only 0.0170%. It is worth noting that WRSM and VPRSM not only have a larger relative error than PM, but also require longer iteration time. Although the iteration time of RSM is similar to that of PM, the relative error of its reliability index is the largest. When the deviation coefficient is equal to 5, the iterative process diagram of each method is shown in Figure 1. It can be seen from the figure that the convergence solution of PM is closest to the exact solution of MCS, and PM only needs 7 iterations to reach the convergence solution, with the least number of iterations. This phenomenon is mutually corroborated by the results of the reliability index and relative error shown in Table 1. According to Table 1, the stability index value of the method in this paper is 0, which shows that the reliability index obtained by the method in this paper is almost not affected by the change of the deviation coefficient. At the same time, it can be obtained that the order of the stability of each method from high to low is PM, VPRSM, WRSM, RSM.In summary, in the case of low nonlinear functions, PM is not only superior in accuracy to other method, but also the iteration convergence speed of the method in this paper is faster and higher in efficiency. In addition, PM has high stability and strong robustness.The reason for the above phenomenon is that the PM method uses the strategy of searching the design point in steps. In the iterative process, according to the characteristics of the design point, it is first iterated in the negative gradient direction with a large step size, so that it can approach the design point faster, and then it is alternately iterated in the linear interpolation with a small step size and the negative gradient direction with a large step size afterwards, which makes the computational efficiency increase and the iteration time decrease.


## 5.2. Example 2: Functional Function with Cross terms

Example 2 is to study the influence of the cross term on the algorithm. The calculation results are listed in Table2, and the iteration process is shown in Figure 4.(17)gx,y=570−xy,x∼N1000,33,y∼N2,0.1.Table 2
Comparison of calculation results of Example 2.
MethodDeviation coefficientfReliability indexβRelative error (%)StabilityIteration times (s)MCS2.1257RSM32.12580.00050.05123.8610e − 0452.12680.05173.8490e − 04WRSM32.12560.00050.03713.9210e − 0452.12490.03763.9780ee − 04VPRSM32.08292.01353.53763.4970e − 0452.00775.55113.5350e − 04PM32.12560.000509.4610e − 0452.12560.00059.5260e − 04Figure 4
Iterative process of Example 2.It can be seen from Table2 that although the iteration time of PM is the longest, the reliability index obtained by PM is closest to the solution of MCS, so the accuracy of PM is the highest. Secondly, as the deviation coefficient changes, the stability index of PM is always 0, indicating that the stability of PM is very high. It can be found intuitively from Figure 4 that the reliability index obtained by PM, WRSM, and RSM is all close to the exact solution, when the relative error of the reliability of VPRSM is the largest.Therefore, it can be concluded that when the performance function contains cross terms, the efficiency of this method is slightly lower than other methods, but in terms of accuracy and stability, this method is still superior to other methods.


## 5.3. Example 3: Nonlinear Function of Concave toward Failure Region

The performance function of this example is shown below, and bothx and y obey the standard normal distribution. The calculation results are shown in Table 3, and the iterative process and the fitting effect are shown in Figure 5 and Figure 6, respectively. The fitting effect of the response surface function obtained in the last iteration of each method to the true failure surface at the most probable point is shown in Figure 6. Among them, LSF represents the true failure surface, and MPP represents the most probable point, that is, the design point.(18)gx,y=exp0.2x+1.4−y.Table 3
Comparison of calculation results of Example 3.
MethodDeviation coefficientfReliability indexβRelative error (%)StabilityIteration times (s)MCS3.3811RSM33.35060.91030.18575.1300e − 0453.35660.72465.0850e − 04WRSM33.34990.92280.03068.5460e − 0453.35120.89228.6010e − 04VPRSM33.48403.04346.95935.3020e − 0453.719310.00275.9950e − 04PM33.34970.928705.1780e − 0453.34970.92875.3230e − 04Figure 5
Iterative process of example 3.Figure 6
Fitting effect of each method.From the research of Nguyen, it can be known that the key of the response surface method is that the response surface function fits the sampling points well, especially near the design point. It can be seen from Figure6 that, except for VPRSM, other methods have a relatively good fitting effect on the true failure surface at the design point, while the vector projection response surface is still a certain distance from the design point. The fitting effect at the design point is not very good, so it can be concluded that the accuracy of PM in this example is much higher than that of VPRSM. This conclusion is also verified in column relative error of Table 3. Figure 5 can also draw conclusions similar to Figure 6. In Figure 5, the reliability index of PM, RSM, and WRSM is slightly smaller than the exact solution, while the reliability of VPRSM is much larger than the exact solution, which also shows that the high accuracy of PM.Table3 accurately shows the characteristics of each method from a quantitative perspective. The reliability index of PM, RSM, and WRSM has similar results, but according to the iteration time, PM and RSM are better than WRSM. Although the iteration time of VPRSM is also very small, the relative error of its reliability is very large, even as high as 10.0027%. Judging from the stability index, we can draw conclusions almost the same as those in the first and second examples. The method in this paper has high stability and strong robustness.


## 5.4. Example 4: Nonlinear Function of Convex toward Failure Region

The type of example 4 is just the opposite of that of example 3. It is a convex function, and its performance function is as follows. Bothx and y obey the standard normal distribution. The calculation results are shown in Table 4, and the iteration path diagram of each method is shown in Figure 7.(19)gx,y=exp0.4x+2+6.2−exp0.3y+5.0−200.Table 4
Comparison of calculation results of example 4.
MethodDeviation coefficientfFailure probabilityRelative error (%)StabilityIteration times (s)MCS2.6855RSM32.71130.96070.01495.5980e − 0452.71170.97566.1820e − 04WRSM32.71030.92350.06336.9080e − 0452.71200.98686.9330e − 04VPRSM32.63192.03660.29895.3510e − 0452.51916.06558.5430e − 04PM32.70990.908605.8180e − 0452.70990.90865.9790e − 04Figure 7
Iterative process of each method.It can be seen from Table4 that the relative errors of the PM, RSM, and WRSM methods are very low, and the iteration time is similar. Although the iteration time of VPRSM is similar, its relative error is the largest. At the same time, it can be concluded that the stability of the method in this paper is still the best (0 stability).It can be seen from the figure above that, except for the vector projection response surface method, other methods all converge to the true failure surface. This indicates that the PM is more accurate than the VPRSM. From the iterative path, it can be seen that the method in this paper can quickly approach the convergence point in only 4 steps, which is significantly faster than other methods in efficiency. This is because the method in this paper adopts a strategy of iterating in the negative gradient direction with a large step size first, and then iterating with a small step size.


## 5.5. Example 5: Ten-Bar Truss Structure

In this example, a ten-bar truss structure（Figure8）is selected. This structure is fixed with hinge support at node 1 and node 2, where E is the elastic modulus E=107psi1psi=6.895kPa, L is the length of the horizontal and vertical rods L=360in1in=2.54cm, and P is the vertical load P=100000lb1lb=0.4536kg. And the section areas A1, A2, and A3 are all independent normal random variables A1=N13,1.3A2∼N2,0.2A3=N9,0.9. The failure state with respect to the serviceability of the structure is reached when the horizontal displacement dreal of the node 5 exceeded a threshold value 3.5 in. Obviously, this is a case of an implicit performance function. And the results are shown in Table 5.(20)gA1,A2,A3=dallow−dreal.Figure 8
Ten-bar truss structure.Table 5
Comparison of calculation results of Example 5.
MethodDeviation coefficientfFailure probabilityRelative error (%)StabilityIteration times (s)MCS0.836895RSM30.8199432.02560.00050.015650.8199472.02510.0155WRSM30.8199462.02520.00050.026250.8199502.02470.0248VPRSM30.8634803.17667.23180.019150.92400210.40840.0192PM30.8199422.025700.010850.8199422.02570.0105In this example, the calculation error obtained by VPRSM is too large and is no longer applicable to this example. Regardless of whether it is in the relative error of the failure probability or in the iteration time, the results displayed by the PM are satisfactory.


## 6. Conclusions

Aiming at the problem of the classic response surface method and its relevant method with low calculation accuracy in some cases, this paper proposes a new method based on double weighting factors idea and vector projection technology to improve the accuracy and stability.Regarding accuracy and stability, the accuracy and stability of the method in this paper are better than the classical response surface method, linear weighted response surface, and vector projection response surface method in the case of low nonlinear performance function, performance function with cross term, concave function, or ten-bar truss structure. In the case of concave functions, the accuracy of the method in this paper is slightly better than that of the vector projection response surface method, which is slightly inferior to other methods. But its stability and robustness are still higher than other methods.Regarding efficiency, although the steps of the method in this paper are more complicated, the method in this paper reduces the iteration of invalid points, which makes the program running time less, and the calculation speed is improved. Compared with other methods, the method in this paper requires less calculation, and this method is more efficient.In summary, the method in this paper improves the accuracy and stability of the response surface method while ensuring the calculation efficiency, and it has certain research significance and application value.



---

*Source: 5104027-2022-05-03.xml*