Abstract
A cubic spline approximation-Bayesian composite quantile regression algorithm is proposed to estimate parameters and structure of the Wiener model with internal noise. Firstly, an ARX model with a high order is taken to represent the linear block; meanwhile, the nonlinear block (reversibility) is approximated by a cubic spline function. Then, parameters are estimated by using the Bayesian composite quantile regression algorithm. In order to reduce the computational burden, the Markov Chain Monte Carlo algorithm is introduced to calculate the expectation of parameters’ posterior distribution. To determine the structure order, the Final Output Error and the Akaike Information Criterion are used in the nonlinear block and the linear block, respectively. Finally, a numerical simulation and an industrial case verify the effectiveness of the proposed algorithm.
1. Introduction
The Wiener model is a nonlinear system, which is composed of a dynamic linear block and a static nonlinear block. Due to its simple structure and strong adaptability, the Wiener model is widely used in industrial processes, such as chemical engineering processes [1], oil recovery [2], biological plants [3], and fluid control units [4]. However, how to identify the parameters and structure of the Wiener model is still a challenging issue.
Most of the existing works are related to identify the Wiener model with a known structure. Xu et al. introduced a sigmoid function to improve the mutation rate F in the Adaptive Differential Evolution algorithm [5]. Furthermore, a blind identification method for the Wiener model was also investigated [6]. The input signal of the system was given by a cyclostationary signal instead of a Gaussian random signal and the internal variables were recovered only based on the output. Then, the order and parameters of the Wiener model were estimated by using the support vector machine regression.
However, the structure and order of the Wiener model are usually unknown in the industrial process. It is difficult to identify the Wiener model because of its nonlinear characteristic. Additionally, the existing research studies are mostly related to the system without noise or only with the output noise [7–9]. Lamia et al. described the Wiener model using the polynomial nonlinear state space (PNLSS) model and developed an output error identification method for the nonlinear block [10]. Riccardo et al. presented a kernel-based identification to estimate parameters of the Wiener system [11]. The impulse response of the linear block and the static nonlinearity were modelled with a Gaussian process and combination of basis functions, respectively. Then, an iterative algorithm using the expectation-maximization method was developed to estimate parameters of the Wiener model. Al-Dhaifallah applied twin support vector regressions to identify the nonlinear Wiener system, including a linear dynamic block [12]. The linear block was expanded in terms of basis functions, while the nonlinear part is determined by twin support vector machine regressions.
Although mentioned algorithms work well, the industrial process with the internal noise is less involved. Lindsten et al. handled the internal white noise in a state space systematic manner and used the nonparametric Gaussian process model for the static nonlinearity [13]. Jing et al. proposed a variable knot spline approximation recursive Bayesian algorithm to reduce the influence of the internal white noise [14]. Li et al. presented a neurofuzzy-based single-input-single-output (SISO) Wiener model identification method for colour noises [15]. Zhang and Mao proposed a robust recursive least squares algorithm with a dead zone weighted factor based on the inverse of the nonlinear function block, which took process noises and measurement noises into consideration [16]. Most of related articles are on the assumption that the internal process noise satisfies the Gaussian distribution or approximately satisfies the symmetrical distribution. Due to the internal noise change with the gain variation of the system, it is difficult to eliminate affections.(a)Variables between the linear block and the nonlinear block cannot be measured, and the existing methods cannot estimate parameters accurately [9].(b)Owing to the nonlinear characteristic of the Wiener model, the internal noise will be amplified with the gain of the nonlinear block [14].(c)The inverse nonlinear block is commonly approximated by using the polynomial function. However, there are some oscillation phenomena in the polynomial with a high order [17]. As a result, the existing algorithms cannot achieve a satisfied estimation when the effect of external interpolation is not good.(d)The best unbiased estimation of parameters can be obtained only if the modelling error satisfies the distribution with the mean zero and the same variance. Because the error cannot be guaranteed to satisfy the Gaussian distribution, the least squares algorithm no longer shows its robustness [18].
To get an unbiased estimation of the Wiener model with the internal noise, a cubic spline approximation-Bayesian composite quantile regression (CSA-BCQR) is presented in this work. The contributions are as follows:(a)The internal process noise with the unknown asymmetric distribution form in the Wiener model is addressed(b)Repeated sampling by using the Markov Chain Monte Carlo (MCMC) method is implemented to achieve a faster convergence of parameters(c)Overcome the high-order oscillation caused by using the polynomial approximation, and a structure detection framework is also considered
The structure is organized as follows. Section 2 gives the problem description. Section 3 states the principle of quantile regression and Bayesian composite quantile regression in detail. Section 4 describes the order determination of the nonlinear block and linear block of the Wiener model. Section 5 presents a numerical simulation and an industrial case to evaluate the proposed algorithm, respectively. Finally, conclusions are summarized in Section 6.
2. Problem Description
As shown in Figure 1, the input signal firstly passes through the linear block of the Wiener model and generates . Under the disturbance of the internal noise , is further developed into and becomes the input of the nonlinear block . The overall output is finally generated.

Assume that(a)The nonlinear function is continuous, monotonic, and reversible(b) is considered as a Gaussian white noise and independent of the input signal (c)Unknown polynomial order can be defined by the Final Output Error (FOE) and the Akaike Information Criterion (AIC)
Here, an ARX model with a high order is employed to approximate the linear block of the Wiener model:
Considering the nonlinear part is reversible, equation (1) is rewritten as follows:
The cost function is defined as follows:where is the size of samples. Because the specific form of the expression is unknown, a stable linear transfer function is approximated by the finite impulse response (FIR) model.where is the order of the linear block.
For the nonlinear block of the Wiener model, a cubic spline approximation (CSA) function is applied to fit the inverse function of the nonlinear clock:where is the order of spline function, are the estimated parameters, are the internal gathering points in the knot set and satisfy , , is the size of samples, and .
Because is unmeasurable, an arbitrary size gain is allocated between the linear module and the nonlinear module. Set and combine equations (2) and (5):
Then, the Wiener model is parameterized as follows:with
3. Bayesian Composite Quantile Regression
3.1. Quantile Regression
According to the theory of the least squares method, the estimate is unbiased when the modelling error . If modelling error and , the estimate is unbiased with the minimum variance. When the data has a sharp/thick tail distribution and a significant heteroscedastic, the least squares algorithm no longer has the abovementioned superiority. To make up those deficiencies, a Quantile Regression (QR) has been proposed by Bang et al. [19].
The nature of the QR is to describe the linear relationship between the independent variable and the dependent variable in the quantile formation. When the quantile changes from 0 to 1, the position and steering of the corresponding regression plane will be adjusted. Therefore, the dependent variable under different quantile conditions can be obtained. QR not only describes the range of variation of the dependent variable but also measures the influence of the conditional distribution shape of the regression variable. Unlike the least squares regression, QR has the best linear unbiased estimates and modelling residuals satisfy the normal distribution.
If the conditional quantile of is a linear function of , then . The quantile estimate of is shown as follows:where , is the linear loss function, and is the indicator function.
Equivalently, . Equation (9) can also be expressed as follows:
When takes different values, the estimated parameter is also different, as well as the expression of the regression equation.
However, the objective function in QR is nondifferentiable or nonconvex, which makes the optimal solution process complicate [20]. Both equations (9) and (10) are piecewise linear functions, so the solution cannot be solved directly from the objective function. It is necessary to find another way to obtain the posterior density function of each QR’s parameter.
Definition 1. A random variable has an Asymmetric Laplace Distribution (ALD) which is noted as with and [21]. Its probability density function (pdf) is given bywhere is a location parameter, is a scale parameter, and is an asymmetry parameter. When , the distribution simplifies to the Laplace distribution.
Assuming that satisfies the ALD, the corresponding pdf is written as follows [22]:and the likelihood function of is defined as follows:Then, minimization of the cost function equation (10) is equivalent to the maximization of the likelihood function equation (13) based on ALD [23].
3.2. Bayesian Composite Quantile Regression
Compared with least squares algorithm, QR has much more flexibility in assessing the effect of predictors on different locations of the response distribution [22]. It also solves the problem that the least squares regression can only describe the local influence of the dependent variable to the dependent variable. Zou and Yuan pointed out that QR can lead to arbitrarily small relative efficiency compared with the least square [24]. In addition, the median regression, a special case of QR, may not be the best choice for some abnormal errors. Therefore, the Composite Quantile Regression (CQR) can aggregate multiple quantile information together to produce a robust parameter estimation. Meanwhile, CQR can achieve an estimated efficiency gain based on a single quantile regression [25].
With , , the CQR of iswhich is also rewritten by
Obviously, the estimated can be calculated by the QR of a combination of . However, the computational burdens of CQR and QR are high. It is necessary to find a good solution.
The posterior distribution of , is given bywhere is the prior distribution of , is the likelihood function:
Lemma 1. If the likelihood function is , then the posterior distribution of , , will have a proper distribution [25], i.e.,
If , , still exists and the expectation of the posterior distribution is the Bayesian estimate of .
In order to improve the computational efficiency, the Markov Chain Monte Carlo (MCMC) method, which is convenient for researchers to sample from complex distributions, has been developed [26]. There are two main MCMC sampling methods: Metropolis-Hastings (M-H) sampling and Gibbs sampling. Here, the M-H algorithm is used. The procedure of M-H sampling is shown in Algorithm 1.
The Bayesian CQR (BCQR) algorithm used in this paper does not depend on the actual distribution of data, but on the likelihood function formed by the ALD [27]. The essence of BCQR is that the estimated parameter is regarded as a random variable, and the sampling distribution of parameter can be obtained by repeated sampling. When a stable distribution is obtained, the mean and standard deviation of the parameters at each quantile can be determined.
|
4. Order Selection of the Wiener Model
The FOE and AIC are commonly used to determine the order of the nonlinear block and linear block, respectively.
4.1. Order Determination of the Nonlinear Block
As shown in equation (5), the order of the inverse nonlinear function is and the corresponding FOE is as follows:where is the size of samples and is the number of parameters. is the factor which overcomes the overfitting issue. However, equation (19) requires the estimated value . Therefore, an equivalent criterion is introduced:where .
Assume
Since the nonlinear function is monotonic, then .
4.2. Order Determination of the Linear Block
FOE criterion needs to know the input of the system and the output of the linear block, but the intermediate variables of the system are unmeasurable. So, an indirect method is adopted to estimate the intermediate variables firstly and the order is determined lately.
According to the estimated parameter , the linear module iswith
Equation (22) can also be rewritten by
The likelihood function of is . The maximum likelihood estimation of is given by
Using AIC criterion,
Finally, the order is determined with .
5. Case Studies
5.1. Numerical Simulation
The discrete Wiener model with an internal noise is given bywhere , , satisfies the ALD, and .
5.1.1. Estimation of
The order of the Cubic Spline Function is estimated by using the FOE criterion under different variances. The variation of the FOE of the cubic spline function with different quantiles was shown in Figure 2. When , FOE is minimum.

5.1.2. Estimation of
As shown in Figure 3, AIC is minimum when , which is completely consistent with the true order of the Wiener model.

5.1.3. Parameter Estimation Using CSA-BCQR
To estimate the error of the linear part and nonlinear part, absolute relative error () and mean square error () are defined as follows:and the comprehensive error () is
The estimated results were shown in Table 1 and Figure 4. To show the superiority of CAS-BCQR algorithm, CSA-Bayesian quantile regression (CSA-BQR) and CSA-recursive least squares (CSA-RLS) algorithms were taken for comparison.

When using CSA-BCQR algorithm, converge at the iteration of 195, 79, and 83, respectively. However, in CSA-BQR algorithm, does not converge when the iteration time is 300 and converge at 91 and 249, respectively. In CSA-RLS algorithm, also fails to converge at 300 and converges at 190 and 298, respectively. Obviously, the convergence speed of CSA-BCQR algorithm is faster than those of the CSA-BQR and CSA-RLS algorithms.
5.2. Industrial Case
A water tank model was taken for an example (shown in Figure 5). is the water inlet flow, is the theoretical liquid level, is the actual liquid level, and is the water outlet flow [28].

The relationship between the variables of the model is as follows:where is a persistent excitation signal sequence, . Here, and and parameters estimated by CSA-BCQR algorithm were shown in Table 2.
As shown in Table 2, the estimation errors of the CSA-BCQR, CSA-BQR, and CSA-RLS algorithms decrease gradually with increasing. Meanwhile the estimated value was getting closer to the true value. When , the ARE of CSA-BCQR algorithm was 0.0282 and the ARE of CSA-BQR and CSA-RLS were 0.1136 and 0.1642, respectively.
As shown in Figure 6, the MSE of CSA-BCQR algorithm has a significant reduction and converges after 100 iterations. Using CSA-BQR and CSA-RLS, the accuracy of the nonlinear block parameter identification is less than that of CSA-BCQR algorithm. It can be seen that the CSA-BCQR algorithm has higher accuracy and faster convergence.

6. Conclusion
A Cubic Spline Approximation-Bayesian Composite Quantile Regression algorithm is presented to estimate the structure and parameter of the Wiener system with internal noise. Using a cubic spline function to approximate the nonlinear block, overcome the high-order oscillation caused by the polynomial approximation. Then, using Bayesian Composite Quantile Regression algorithm to consider different information of quantiles can effectively improve the accuracy of parameter estimation. A numerical simulation and an industrial case show that the CSA-BCQR algorithm has faster convergence speed and higher parameter identification accuracy compared with CSA-BQR and CSA-RLS algorithms. Furthermore, the CSA-BCQR algorithm may also be applied to other block-oriented models, such as a Hammerstein system or a Hammerstein-Wiener system.
Data Availability
The detailed mechanism model and model parameters of the Wiener model are given in the manuscript. The results are computed on the MATLAB software with the model and given parameters, while the relevant results are also given in the manuscript.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
This work was supported by the National Nature Science Foundation of China (Grant no. 61873113), the 333 High-level Talent Training Program of Jiangsu Province, the Key R&D Program of Jiangsu (Grant no. BE2018370), and the Postgraduate Research and Practice Innovation Program of Jiangsu Province (Grant no. KYCX17_1785).