Optimal Control for a Linear System Subject to a General ARIMA Disturbance

Abstract

A novel run-to-run (R2R) control algorithm based on Kalman filter approach is proposed to deal with a linear system with a general ARIMA( p,d,q ) process, in the presence of measurement error and adjustment error together with a random initial bias. We mathematically prove its optimality. The performance of the newly proposed controller and the proportional-integral ( PI ) controller is evaluated and compared under multiple scenarios through Monte Carlo simulations. Almost all the results reflect the new controller's superiority.

Hongyan Xie 1 and Fangyi He 2

Academic Editor:Alain Vande Wouwer

1, School of Economics, Southwestern University of Finance and Economics, Chengdu, Sichuan 611130, China
2, School of Finance, Southwestern University of Finance and Economics, Chengdu, Sichuan 611130, China

Received 28 November 2014; Revised 24 January 2015; Accepted 16 February 2015; 25 February 2015

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction and Motivation

Optimizing a sequence of actions to obtain some future goal is the general topic of control theory (see [1]). Many optimal control methods have been widely used in actuarial science (e.g., [2, 3]), production management (e.g., [4-6]), and quality management (e.g., [7-9]). Specifically, disturbance rejection is one of the major concerns in quality control problems, such as machine setup adjustment problems and R2R process control problems in semiconductor manufacturing (see [10, 11]). Grubbs [12] originally studies the machine setup adjustment problem and develops an optimal adjustment rule for a linear system with normal disturbances, which is referred to by Trietsch [13] as the "harmonic rule." Vander Wiel et al. [14] give an optimal control algorithm for a linear system with ARMA( [figure omitted; refer to PDF] ) disturbances. For a linear system with IMA( [figure omitted; refer to PDF] ) disturbances, Ingolfsson and Sachs [15] and Box et al. [16] introduce the exponential weighted moving average (EWMA) control algorithm and prove its optimality; He et al. [17] recently develop an optimal control algorithm, the ARMA controller, for a linear system with a general ARMA( [figure omitted; refer to PDF] ) disturbance.

In this paper, we will extend the results in [17] and develop an optimal control algorithm for a linear system with a general ARIMA( [figure omitted; refer to PDF] ) disturbance, where [figure omitted; refer to PDF] . The ARIMA disturbance has been widely used to describe process dynamics (see [18-21]). Since ARIMA( [figure omitted; refer to PDF] ) processes can be used to model a large class of nonstationary disturbances, many machine setup problems and R2R control problems can be solved under this new framework. Similar to [17], the newly proposed controller will consider a more realistic case such that both measurement error and adjustment error exist, and the initial bias of the process is a random variable. We present the problem as follows.

Supposing a process to be controlled [figure omitted; refer to PDF] can be expressed as [figure omitted; refer to PDF] where [figure omitted; refer to PDF] is the state of [figure omitted; refer to PDF] at time [figure omitted; refer to PDF] ; [figure omitted; refer to PDF] is an ARIMA( [figure omitted; refer to PDF] ) process satisfying [figure omitted; refer to PDF] where [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , and [figure omitted; refer to PDF] is a white noise process with mean 0 and variance [figure omitted; refer to PDF] , that is, [figure omitted; refer to PDF] ; [figure omitted; refer to PDF] is the measurement error for [figure omitted; refer to PDF] and [figure omitted; refer to PDF] . [figure omitted; refer to PDF] is the backward shift operator defined by [figure omitted; refer to PDF] Neither [figure omitted; refer to PDF] nor [figure omitted; refer to PDF] can be measured or observed directly. At time [figure omitted; refer to PDF] , suppose we need to make an adjustment of magnitude [figure omitted; refer to PDF] to [figure omitted; refer to PDF] to bring the process output to target [figure omitted; refer to PDF] in the next run. That is, [figure omitted; refer to PDF] where [figure omitted; refer to PDF] is the adjustment error and [figure omitted; refer to PDF] . In practice, the process adjustment [figure omitted; refer to PDF] is assumed to be done by adjusting one controllable factor, [figure omitted; refer to PDF] , via the following model: [figure omitted; refer to PDF] where [figure omitted; refer to PDF] is called the process gain and [figure omitted; refer to PDF] . Assume [figure omitted; refer to PDF] , [figure omitted; refer to PDF] . The initial value [figure omitted; refer to PDF] is assumed to be a random variable with [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , and [figure omitted; refer to PDF] . Without loss of generality, the target [figure omitted; refer to PDF] is assumed to be 0 in the rest of this paper. For a positive time index [figure omitted; refer to PDF] , we hope to determine the optimal [figure omitted; refer to PDF] , [figure omitted; refer to PDF] that satisfy [figure omitted; refer to PDF] where [figure omitted; refer to PDF] is the expectation operator conditional on all the information until time [figure omitted; refer to PDF] . Problem (6) is a finite-horizon problem. When [figure omitted; refer to PDF] , we get the infinite-horizon problem that searches the optimal [figure omitted; refer to PDF] satisfying [figure omitted; refer to PDF] The purpose of this paper is to find the solutions to both problems (6) and (7).

The rest of this paper is organized as follows. In Section 2 we derive the state-space representation of the system and get the recursive estimation formulae of the system states using Kalman filter. In Section 3, we develop a control algorithm for the system and prove its optimality without normal distribution assumptions. We further give the implementation steps of the controller in practice. Simulation studies are done in Section 4 under multiple scenarios. Section 5 gives an illustrative example for the application of the control algorithm. Concluding remarks are included in Section 6.

2. The State-Space Representation

Let [figure omitted; refer to PDF] and [figure omitted; refer to PDF] be the coefficients of [figure omitted; refer to PDF] in the power series expansion of [figure omitted; refer to PDF] , [figure omitted; refer to PDF] . Brockwell and Davis [22] gave a state-space representation for a general ARIMA( [figure omitted; refer to PDF] ) process. By extending their results we can get the following Theorem 1 easily.

Theorem 1.

The linear system (1)-(5) has a state-space representation as follows: [figure omitted; refer to PDF] where [figure omitted; refer to PDF] and [figure omitted; refer to PDF] are the [figure omitted; refer to PDF] and [figure omitted; refer to PDF] matrices defined by [figure omitted; refer to PDF] if [figure omitted; refer to PDF] and [figure omitted; refer to PDF] if [figure omitted; refer to PDF] ; [figure omitted; refer to PDF] and [figure omitted; refer to PDF] are the [figure omitted; refer to PDF] and [figure omitted; refer to PDF] matrices defined by [figure omitted; refer to PDF] , [figure omitted; refer to PDF] if [figure omitted; refer to PDF] and [figure omitted; refer to PDF] if [figure omitted; refer to PDF] ; [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , [figure omitted; refer to PDF] .

A proof of Theorem 1 is presented in Appendix A. For simplicity, in the rest of this paper we denote [figure omitted; refer to PDF] where [figure omitted; refer to PDF] is defined to be a [figure omitted; refer to PDF] vector of ones.

That the problems of optimal control and state estimation can be decoupled in certain cases is one of the most fundamental principles in feedback control theory. This is known as the separation principle (see [23]). The Kalman filter produces the statistically optimal estimate of the state of system (8)-(9). Note the fact that the disturbances of (8) and (9) are correlated to each other, so we should use the Kalman filter formulae for correlated measurement and process noise. At any time [figure omitted; refer to PDF] , let us define [figure omitted; refer to PDF] and [figure omitted; refer to PDF] and assume [figure omitted; refer to PDF] and [figure omitted; refer to PDF] . Using the results on page 123 in Lewis' book [24], we can directly get the following Lemma 2.

Lemma 2.

For the systems (8)-(9), we have the following recursive formulae: [figure omitted; refer to PDF] where [figure omitted; refer to PDF] is called the modified Kalman Gain that [figure omitted; refer to PDF]

3. The Optimal Control Algorithm

In this section, we will derive an optimal control algorithm without the normal distribution assumptions on [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , and [figure omitted; refer to PDF] , [figure omitted; refer to PDF] . This optimal control algorithm can be applied to both the finite-horizon and infinite-horizon problems.

Theorem 3.

The optimal control algorithm that solves both [figure omitted; refer to PDF] for the system (1)-(5) is [figure omitted; refer to PDF] where [figure omitted; refer to PDF] is updated according to [figure omitted; refer to PDF] and (14), (15) are with initial values [figure omitted; refer to PDF] and [figure omitted; refer to PDF] .

Proof.

For a given positive time index [figure omitted; refer to PDF] , let [figure omitted; refer to PDF] and [figure omitted; refer to PDF] . We can get the Bellman equation as [figure omitted; refer to PDF] From (8), it is evident that [figure omitted; refer to PDF] . Using the fact that [figure omitted; refer to PDF] and the property of trace operator, we obtain that [figure omitted; refer to PDF] From (13) and the fact that [figure omitted; refer to PDF] , we have [figure omitted; refer to PDF] Let [figure omitted; refer to PDF] . It is evident that [figure omitted; refer to PDF] is irrelated to [figure omitted; refer to PDF] since [figure omitted; refer to PDF] is updated based on (14) and (15). Now (19) changes to [figure omitted; refer to PDF] When [figure omitted; refer to PDF] , we get [figure omitted; refer to PDF] Solving the first order condition, we get that [figure omitted; refer to PDF] Putting (25) into (24), we get [figure omitted; refer to PDF] . When [figure omitted; refer to PDF] , we have [figure omitted; refer to PDF] Solving the first order condition, we get [figure omitted; refer to PDF] and [figure omitted; refer to PDF] . Repeating the procedures above, we can prove that (17) is the optimal adjustment strategy at time [figure omitted; refer to PDF] and get [figure omitted; refer to PDF] At last, we can derive (18) by putting (17) into (13). Since (17) is irrelated to the given time index [figure omitted; refer to PDF] , the control strategy also solves the infinite-horizon problem. Then we finish the proof.

In practice, the steps to implement the proposed control algorithm are as follows.

Step 0. Setup the algorithm's parameters and initial values, such as [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , based on experience or history data.

Step 1. Compute [figure omitted; refer to PDF] based on (15).

Step 2. Collect the new observation [figure omitted; refer to PDF] and compute the adjustment [figure omitted; refer to PDF] based on (17).

Step 3. Update [figure omitted; refer to PDF] based on (18).

Step 4. Update [figure omitted; refer to PDF] based on (14); let [figure omitted; refer to PDF] and go back to Step 1.

Note that Step 0 is an off-line procedure and Step 1 to Step 4 are on-line procedures. As the newly proposed control algorithm is specially designed for adjusting any ARIMA disturbances, we call it ARIMA controller in the rest of this paper.

4. Simulation Study

In this section, we will study the performance of the ARIMA controller under multiple scenarios through Monte Carlo simulations. For simplicity, we only focus on the system (1)-(5) with [figure omitted; refer to PDF] being an ARIMA( [figure omitted; refer to PDF] ) disturbance, although the ARIMA controller can be applied to any general ARIMA( [figure omitted; refer to PDF] ) disturbance. Without loss of generality, we set the parameters in the ARIMA( [figure omitted; refer to PDF] ) disturbance as [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , and [figure omitted; refer to PDF] . For comparison, a [figure omitted; refer to PDF] controller's performance is also evaluated. The [figure omitted; refer to PDF] controller is widely used in feedback control, and it involves two separate constant parameters: the proportional and the integral values, denoted by [figure omitted; refer to PDF] and [figure omitted; refer to PDF] . The steps to implement a [figure omitted; refer to PDF] controller can be found in [25].

Tuning a [figure omitted; refer to PDF] controller could be a challenging task. There are many ways of tuning, such as Ziegler-Nichols tuning, lambda tuning, robust loop shaping, optimization methods, and others (see [26]). Among them, optimization methods are powerful and direct ways. In order to do fair comparisons with the ARIMA controller, we need to choose [figure omitted; refer to PDF] and [figure omitted; refer to PDF] to make the [figure omitted; refer to PDF] controller have the best performance for controlling the specified ARIMA( [figure omitted; refer to PDF] ) disturbance. To the best of our knowledge, there are no closed-form expressions for [figure omitted; refer to PDF] and [figure omitted; refer to PDF] to optimally control a ARIMA( [figure omitted; refer to PDF] ) process; then experimental tuning of the [figure omitted; refer to PDF] parameters has to be performed. The procedure of optimally choosing [figure omitted; refer to PDF] and [figure omitted; refer to PDF] in this paper is described in Appendix B. For the ARIMA( [figure omitted; refer to PDF] ) disturbance with [figure omitted; refer to PDF] and [figure omitted; refer to PDF] , we get [figure omitted; refer to PDF] and [figure omitted; refer to PDF] according to the optimization procedure.

Four scenarios are examined in the following simulations. In Scenario 1, the effect of the measurement error and the adjustment error on the ARIMA controller and the [figure omitted; refer to PDF] controller is investigated; in Scenario 2, the effect of the process initial bias [figure omitted; refer to PDF] on the two controllers is explored; in Scenario 3, we study how the estimate of the process gain [figure omitted; refer to PDF] affects the two controllers; in Scenario 4, we investigate both controllers' performance when the disturbance parameters are not estimated accurately.

We do 1000 replications for each case and run 100 steps for [figure omitted; refer to PDF] in each replication. The first 100-run mean square error (MSE) of [figure omitted; refer to PDF] is computed. The Average MSE (AMSE) of the 1000 replications is reported in Tables 1-4. Also the standard error in the AMSE (SEAMSE) is computed, which is [figure omitted; refer to PDF] where SDMSE is the standard deviation of mean square errors, and the number of replicates here is 1000. AMSE measures the performance of the controllers, and SEAMSE reflects the variability of the AMSE. The smaller the AMSE is, the better performance the controller has.

Table 1: AMSE when both types of error exist. [figure omitted; refer to PDF] is the ARIMA(1, 1, 1) disturbance with [figure omitted; refer to PDF] and [figure omitted; refer to PDF] .

ARIMA controller [figure omitted; refer to PDF] [figure omitted; refer to PDF]
[figure omitted; refer to PDF]	[figure omitted; refer to PDF]	0.5	1.0	1.5	2.0
0	1.014 (.005)	1.298 (.006)	2.104 (.009)	3.392 (.015)	5.155 (.023)
0.5	1.657 (.008)	1.922 (.009)	2.695 (.012)	3.959 (.018)	5.709 (.026)
1.0	3.165 (.014)	3.431 (.016)	4.205 (.019)	5.471 (.025)	7.225 (.032)
1.5	5.286 (.024)	5.569 (.025)	6.381 (.029)	7.698 (.035)	9.503 (.043)
2.0	7.966 (.037)	8.269 (.038)	9.136 (.041)	10.530 (.047)	12.420 (.056)
[figure omitted; refer to PDF] controller [figure omitted; refer to PDF] [figure omitted; refer to PDF]
[figure omitted; refer to PDF]	[figure omitted; refer to PDF]	0.5	1.0	1.5	2.0
0	1.084 (.005)	1.379 (.006)	2.265 (.010)	3.742 (.017)	5.809 (.027)
0.5	1.872 (.009)	2.172 (.010)	3.063 (.015)	4.544 (.022)	6.616 (.032)
1.0	4.247 (.023)	4.552 (.025)	5.448 (.029)	6.934 (.036)	9.010 (.047)
1.5	8.211 (.048)	8.520 (.049)	9.420 (.053)	10.911 (.060)	12.992 (.071)
2.0	13.761 (.082)	14.075 (.083)	14.980 (.088)	16.475 (.095)	18.561 (.105)

The corresponding SEAMSE's are enclosed in the parentheses.

Table 2: AMSE when the initial bias' mean and standard deviation vary. [figure omitted; refer to PDF] is the ARIMA(1, 1, 1) disturbance with [figure omitted; refer to PDF] and [figure omitted; refer to PDF] .

ARIMA controller [figure omitted; refer to PDF] [figure omitted; refer to PDF]
[figure omitted; refer to PDF]	[figure omitted; refer to PDF]	1	2	3	4
-4	1.176 (.006)	1.180 (.007)	1.203 (.008)	1.244 (.010)	1.304 (.013)
-3	1.106 (.006)	1.110 (.006)	1.134 (.007)	1.176 (.009)	1.237 (.011)
-2	1.055 (.005)	1.061 (.005)	1.085 (.006)	1.128 (.008)	1.190 (.010)
-1	1.025 (.005)	1.031 (.005)	1.057 (.005)	1.101 (.007)	1.163 (.009)
0	1.014 (.005)	1.021 (.005)	1.048 (.005)	1.093 (.006)	1.156 (.008)
1	1.024 (.005)	1.032 (.005)	1.059 (.005)	1.105 (.007)	1.169 (.009)
2	1.053 (.005)	1.062 (.005)	1.091 (.006)	1.137 (.008)	1.202 (.010)
3	1.103 (.006)	1.113 (.006)	1.142 (.007)	1.189 (.009)	1.255 (.011)
4	1.172 (.007)	1.183 (.007)	1.213 (.008)	1.261 (.011)	1.328 (.013)
[figure omitted; refer to PDF] controller [figure omitted; refer to PDF] [figure omitted; refer to PDF]
[figure omitted; refer to PDF]	[figure omitted; refer to PDF]	1	2	3	4
-4	1.279 (.007)	1.284 (.008)	1.311 (.010)	1.360 (.012)	1.432 (.016)
-3	1.194 (.006)	1.200 (.007)	1.228 (.008)	1.279 (.010)	1.351 (.013)
-2	1.133 (.006)	1.140 (.006)	1.170 (.007)	1.221 (.009)	1.294 (.011)
-1	1.096 (.005)	1.105 (.005)	1.135 (.006)	1.187 (.008)	1.261 (.010)
0	1.084 (.005)	1.093 (.005)	1.124 (.006)	1.177 (.007)	1.253 (.010)
1	1.095 (.005)	1.105 (.005)	1.137 (.006)	1.192 (.008)	1.268 (.010)
2	1.130 (.006)	1.141 (.006)	1.174 (.007)	1.230 (.009)	1.307 (.011)
3	1.189 (.007)	1.201 (.007)	1.235 (.008)	1.292 (.010)	1.370 (.013)
4	1.272 (.007)	1.285 (.008)	1.320 (.009)	1.378 (.012)	1.458 (.016)

The corresponding SEAMSE's are enclosed in the parentheses. Neither measurement errors nor adjustment errors are considered in the simulations.

Table 3: AMSE with different estimates of the process gain. The true process gain [figure omitted; refer to PDF] .

[figure omitted; refer to PDF] is ARIMA(1, 1, 1) [figure omitted; refer to PDF] [figure omitted; refer to PDF]
[figure omitted; refer to PDF]	ARIMA controller	[figure omitted; refer to PDF] controller	[figure omitted; refer to PDF]	ARIMA controller	[figure omitted; refer to PDF] controller
( [figure omitted; refer to PDF] [figure omitted; refer to PDF] )	( [figure omitted; refer to PDF] [figure omitted; refer to PDF] )	( [figure omitted; refer to PDF] [figure omitted; refer to PDF] )	( [figure omitted; refer to PDF] [figure omitted; refer to PDF] )
AMSE (SEAMSE)	AMSE (SEAMSE)	AMSE (SEAMSE)	AMSE (SEAMSE)
1.7	2.207 (.020)	5.113 (.084)	2.6	1.017 (.005)	1.088 (.005)
1.8	1.583 (.011)	2.213 (.021)	2.7	1.024 (.005)	1.099 (.005)
1.9	1.320 (.008)	1.597 (.011)	2.8	1.036 (.005)	1.115 (.006)
2.0	1.183 (.006)	1.343 (.008)	2.9	1.051 (.005)	1.135 (.006)
2.1	1.104 (.006)	1.215 (.006)	3.0	1.068 (.005)	1.159 (.006)
2.2	1.058 (.005)	1.145 (.006)	3.1	1.087 (.005)	1.186 (.006)
2.3	1.032 (.005)	1.107 (.005)	3.2	1.108 (.005)	1.215 (.007)
2.4	1.018 (.005)	1.089 (.005)	3.3	1.131 (.006)	1.246 (.007)
2.5	1.014 (.005)	1.084 (.005)	3.4	1.155 (.006)	1.279 (.007)

Table 4: AMSE with uncertainties in the ARIMA parameters.

ARIMA controller [figure omitted; refer to PDF] [figure omitted; refer to PDF]
[figure omitted; refer to PDF]	[figure omitted; refer to PDF]	-0.4	-0.3	-0.2	-0.1
0.5	1.177 (.006)	1.109 (.005)	1.061 (.005)	1.034 (.005)	1.027 (.005)
0.6	1.106 (.005)	1.056 (.005)	1.027 (.005)	1.018 (.005)	1.031 (.005)
0.7	1.052 (.005)	1.022 (.005)	1.014 (.005)	1.028 (.005)	1.064 (.005)
0.8	1.024 (.005)	1.021 (.005)	1.044 (.005)	1.091 (.006)	1.162 (.006)
0.9	1.073 (.006)	1.129 (.007)	1.218 (.008)	1.340 (.010)	1.495 (.012)
0.95	1.248 (.010)	1.398 (.013)	1.596 (.018)	1.843 (.023)	2.138 (.029)
[figure omitted; refer to PDF] controller [figure omitted; refer to PDF] [figure omitted; refer to PDF]
[figure omitted; refer to PDF]	[figure omitted; refer to PDF]	-0.4	-0.3	-0.2	-0.1
0.5	1.199 (.006)	1.137 (.006)	1.097 (.005)	1.078 (.005)	1.081 (.005)
0.6	1.136 (.006)	1.095 (.005)	1.077 (.005)	1.081 (.005)	1.109 (.005)
0.7	1.092 (.005)	1.076 (.005)	1.084 (.005)	1.116 (.006)	1.172 (.006)
0.8	1.074 (.005)	1.087 (.005)	1.128 (.006)	1.197 (.007)	1.293 (.008)
0.9	1.094 (.005)	1.151 (.006)	1.242 (.007)	1.366 (.009)	1.523 (.011)
0.95	1.129 (.006)	1.220 (.007)	1.348 (.009)	1.515 (.011)	1.720 (.014)

The corresponding SEAMSE's are included in the parentheses.

4.1. Effects of Measurement Error and Adjustment Error

In order to focus on the effects of the measurement error [figure omitted; refer to PDF] and the adjustment error [figure omitted; refer to PDF] on the controllers, we set the initial process bias [figure omitted; refer to PDF] , that is, [figure omitted; refer to PDF] , and set the estimate of the process gain [figure omitted; refer to PDF] to the true value [figure omitted; refer to PDF] . We also assume all the ARIMA( [figure omitted; refer to PDF] ) parameters can be accurately estimated, that is, [figure omitted; refer to PDF] and [figure omitted; refer to PDF] .

Table 1 and Figure 1 show the performance of the two controllers affected by the measurement errors and adjustment errors. We set the values of [figure omitted; refer to PDF] and [figure omitted; refer to PDF] to be 0, 0.5, 1.0, 1.5 and 2.0, respectively, so there are 25 pairs of ( [figure omitted; refer to PDF] ). For each pair of ( [figure omitted; refer to PDF] ), we repeat the simulations. We can observe that for both controllers the AMSE of [figure omitted; refer to PDF] increases when [figure omitted; refer to PDF] or [figure omitted; refer to PDF] increases. However, the AMSE of [figure omitted; refer to PDF] using the ARIMA controller is always smaller than that using the [figure omitted; refer to PDF] controller at a significant level for all the pairs of [figure omitted; refer to PDF] .

Figure 1: The performance of the ARIMA controller and the [figure omitted; refer to PDF] controller when both measurement error ( [figure omitted; refer to PDF] ) and adjustment error ( [figure omitted; refer to PDF] ) exist and vary.

[figure omitted; refer to PDF]

4.2. Effects of Initial Bias

For varying sizes of [figure omitted; refer to PDF] and [figure omitted; refer to PDF] in the prior distribution of the process initial bias, Table 2 and Figure 2 show the results of the AMSE of [figure omitted; refer to PDF] using the ARIMA controller or the [figure omitted; refer to PDF] controller for the ARIMA( [figure omitted; refer to PDF] ) disturbance. In order to focus on the effects of the initial bias, we set all the [figure omitted; refer to PDF] and [figure omitted; refer to PDF] to 0, [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , and [figure omitted; refer to PDF] .

Figure 2: The performance of the ARIMA controller and the [figure omitted; refer to PDF] controller when the process initial bias varies.

[figure omitted; refer to PDF]

It can be seen that for both controllers, the AMSE increases with [figure omitted; refer to PDF] when [figure omitted; refer to PDF] is fixed; and the AMSE increases with [figure omitted; refer to PDF] when [figure omitted; refer to PDF] is fixed. For the same pair of [figure omitted; refer to PDF] , the AMSE of [figure omitted; refer to PDF] using the ARIMA controller is always smaller than that using the [figure omitted; refer to PDF] controller at a significant level.

4.3. Effect of Estimation Uncertainties in the Process Gain

In order to focus on the effect of the estimate of the process gain brought to the controller's performance, we assume both measurement error and adjustment error are absent, that is, [figure omitted; refer to PDF] ; the initial bias [figure omitted; refer to PDF] is 0, that is, [figure omitted; refer to PDF] ; and all the ARIMA parameters are accurately estimated, that is, [figure omitted; refer to PDF] and [figure omitted; refer to PDF] .

We set the process gain [figure omitted; refer to PDF] to 2.5 and assume [figure omitted; refer to PDF] . Table 3 and Figure 3 present the performance of the two controllers when [figure omitted; refer to PDF] varies. It is shown that the ARIMA controller always outperforms the [figure omitted; refer to PDF] controller for all the [figure omitted; refer to PDF] s we examined at a significant level. Additionally, underestimated [figure omitted; refer to PDF] (i.e., [figure omitted; refer to PDF] ) will hurt both controllers' performance more than overestimated [figure omitted; refer to PDF] (i.e., [figure omitted; refer to PDF] ). It is shown that if the process gain [figure omitted; refer to PDF] is underestimated to some degree, both controllers would completely fail.

Figure 3: The performance of the ARIMA controller and the [figure omitted; refer to PDF] controller when the process gain is not accurately estimated.

[figure omitted; refer to PDF]

4.4. Effect of Estimation Uncertainties in the Disturbance

Again, in order to focus on the effect of the estimation uncertainties in the disturbance, we assume [figure omitted; refer to PDF] and [figure omitted; refer to PDF] . Suppose the disturbance parameters are estimated as [figure omitted; refer to PDF] and [figure omitted; refer to PDF] . However, the disturbance's true parameters may not be equal to the estimates, that is, [figure omitted; refer to PDF] and/or [figure omitted; refer to PDF] . We choose 5 values of [figure omitted; refer to PDF] and 5 values of [figure omitted; refer to PDF] as the true parameters for the ARIMA( [figure omitted; refer to PDF] ) disturbance, so there are a total of 25 disturbance ARIMA( [figure omitted; refer to PDF] ) models that are studied to investigate both controllers' robustness.

Table 4 reports the results of the AMSE of [figure omitted; refer to PDF] using the ARIMA controller and the [figure omitted; refer to PDF] controller for different ARIMA( [figure omitted; refer to PDF] ) disturbances with varying [figure omitted; refer to PDF] and [figure omitted; refer to PDF] . It is shown that the ARIMA controller's performance is better than the [figure omitted; refer to PDF] controller's in most of the cases at a significant level. Only when [figure omitted; refer to PDF] is the [figure omitted; refer to PDF] controller's performance better than the ARIMA controller's. This is because a ARIMA( [figure omitted; refer to PDF] ) process will become a ARIMA( [figure omitted; refer to PDF] ) process when [figure omitted; refer to PDF] . As this ARIMA controller is purposely designed to control a ARIMA( [figure omitted; refer to PDF] ) disturbance, its performance will be contaminated when [figure omitted; refer to PDF] is close to 1. The results also imply that if a ARIMA( [figure omitted; refer to PDF] ) disturbance is misidentified to be a ARIMA( [figure omitted; refer to PDF] ) disturbance, the ARIMA controller would have a bad performance.

For illustration, Figure 4 further presents the contour plots for the AMSE as a function of [figure omitted; refer to PDF] and [figure omitted; refer to PDF] by using the ARIMA controller and the [figure omitted; refer to PDF] controller, respectively. It can be generally observed that both controllers that are designed based on underestimated disturbance parameters (i.e., [figure omitted; refer to PDF] , [figure omitted; refer to PDF] ) will have worse performance than the ones designed based on overestimated disturbance parameters (i.e., [figure omitted; refer to PDF] , [figure omitted; refer to PDF] ). Specially, as can be seen from Figure 4(a), when [figure omitted; refer to PDF] is not greater than 0.9, the AMSE of [figure omitted; refer to PDF] by using the ARIMA controller changes very slowly with varying [figure omitted; refer to PDF] and [figure omitted; refer to PDF] . However, when [figure omitted; refer to PDF] is greater than 0.9, the AMSE increases very quickly with [figure omitted; refer to PDF] increasing. The parameter [figure omitted; refer to PDF] has less effects than [figure omitted; refer to PDF] on the performance of the ARIMA controller. In comparison with the ARIMA controller, the [figure omitted; refer to PDF] controller's performance is affected by both [figure omitted; refer to PDF] and [figure omitted; refer to PDF] in almost the same degree as shown in Figure 4(b). In all the change area of [figure omitted; refer to PDF] and [figure omitted; refer to PDF] except the area [figure omitted; refer to PDF] , the [figure omitted; refer to PDF] controller's performance is worse than the ARIMA controller's. However, the maximum AMSE (the worst case) used by the [figure omitted; refer to PDF] controller in the whole change area of [figure omitted; refer to PDF] and [figure omitted; refer to PDF] is smaller than that used by the ARIMA controller. In this sense, we can say that the [figure omitted; refer to PDF] controller's performance is more robust than the ARIMA controller's when uncertainties exist in the ARIMA parameters.

Contour plots of AMSE for the ARIMA controller and the [figure omitted; refer to PDF] controller when the disturbance parameters [figure omitted; refer to PDF] vary.

(a) ARIMA controller

[figure omitted; refer to PDF]

(b) [figure omitted; refer to PDF] controller

[figure omitted; refer to PDF]

5. An Illustrative Example

The studies in the last section are carried out based on an ARIMA( [figure omitted; refer to PDF] ) disturbance. As aforementioned, the proposed ARIMA controller can be applied to any general ARIMA( [figure omitted; refer to PDF] ) disturbance. For illustration purpose, we assume a process as follows: [figure omitted; refer to PDF] and [figure omitted; refer to PDF] , [figure omitted; refer to PDF] ; [figure omitted; refer to PDF] is the deviation from target. The disturbance model [figure omitted; refer to PDF] is the same as one of the ARIMA( [figure omitted; refer to PDF] ) disturbances presented in [21]. For the sake of simplicity, we further set [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , and assume that all the disturbance parameters are known. That is, we ignore measurement errors, adjustment errors, initial bias uncertainties, and disturbance parameter uncertainties and only show the ARIMA controller's superiority in controlling higher order ARIMA disturbance models.

For the above ARIMA( [figure omitted; refer to PDF] ) disturbance, we can get [figure omitted; refer to PDF] and [figure omitted; refer to PDF] So, in the offline procedure, the ARIMA controller's parameters are given as follows: [figure omitted; refer to PDF] The initial values [figure omitted; refer to PDF] . Using the same optimization tuning techniques described in Appendix B, we get [figure omitted; refer to PDF] and [figure omitted; refer to PDF] .

We randomly draw one simulation and show the paths of the process output [figure omitted; refer to PDF] in Figure 5 when the two control algorithms are implemented, respectively. The MSE of [figure omitted; refer to PDF] implemented by the [figure omitted; refer to PDF] controller is 0.819 while that implemented by the ARIMA controller is 0.646. The paths suggest that the ARIMA controller maintains the process output [figure omitted; refer to PDF] closer to the target than the [figure omitted; refer to PDF] controller in much more places in the first 100 simulated runs.

Figure 5: One simulated path by using the ARIMA and [figure omitted; refer to PDF] controllers, respectively, for the ARIMA(1,1,3) disturbance that [figure omitted; refer to PDF] with [figure omitted; refer to PDF] .

[figure omitted; refer to PDF]

6. Concluding Remarks

In this paper, we have developed an optimal control algorithm, the ARIMA controller, for a linear system with a general ARIMA( [figure omitted; refer to PDF] ) disturbance, in the presence of measurement errors and adjustment errors together with a random initial bias. We theoretically prove that the ARIMA controller is optimal to both the finite-horizon and infinite-horizon problems. The performances of the ARIMA and [figure omitted; refer to PDF] controllers have been evaluated and compared via Monte Carlo simulations under multiple scenarios. In almost all the analyzed scenarios, the ARIMA controller outperforms the [figure omitted; refer to PDF] controller. Only when uncertainties exist in the ARIMA parameters does the [figure omitted; refer to PDF] controller have more robust performance than the ARIMA controller. Our simulation studies also show that when designing the two controllers, underestimating the process parameters, including the process gain and the ARIMA parameters, will have more negative impact on the controllers' performance than overestimating the parameters. Although the steps to implement the [figure omitted; refer to PDF] controller are simple, a tuning process is always needed to determine the [figure omitted; refer to PDF] parameters [figure omitted; refer to PDF] and [figure omitted; refer to PDF] for each given disturbance. It is usually time-consuming. Conversely, the ARIMA controller can automatically determine its algorithm parameters for each given disturbance.

The ARIMA controller is a complement to the ARMA controller which is designed to optimally control any general ARMA( [figure omitted; refer to PDF] ) disturbance. The ARMA controller can be used to control a large class of weakly stationary disturbances while the ARIMA controller can handle a large class of nonstationary or periodic disturbances. Such disturbances are common in rapid thermal processing, reactive ion etching, I-line lithography, and lapping processes in the semiconductor manufacturing.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (Grant no. 71102145) and the Fundamental Research Funds for the Central Universities of China (Grant no. JBK150144). The authors are grateful to the editors and the anonymous reviewers for their valuable comments and suggestions that have greatly improved the quality of this paper.

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

References

[1] R. F. Stengel Optimal Control and Estimation , Dover Publications, New York, NY, USA, 1993.

[2] A. Cadenillas, T. Choulli, M. Taksar, L. Zhang, "Classical and impulse stochastic control for the optimization of the dividend and risk policies of an insurance firm," Mathematical Finance , vol. 16, no. 1, pp. 181-202, 2006.

[3] N. Kulenko, H. Schmidli, "Optimal dividend strategies in a Cramer-Lundberg model with capital injections," Insurance: Mathematics and Economics , vol. 43, no. 2, pp. 270-278, 2008.

[4] K. Nakashima, H. Arimitsu, T. Nose, S. Kuriyama, "Optimal control of a remanufacturing system," International Journal of Production Research , vol. 42, no. 17, pp. 3619-3625, 2004.

[5] H. Liu, Z. B. Zabinsky, W. Kohn, "Rule-based forecasting and production control system design utilizing a feedback control architecture," IIE Transactions , vol. 43, no. 2, pp. 143-152, 2011.

[6] G. Singer, E. Khmelnitsky, "A finite-horizon, stochastic optimal control policy for a production-inventory system with backlog-dependent lost sales," IIE Transactions , vol. 42, no. 12, pp. 855-864, 2010.

[7] F. He, K. Wang, W. Jiang, "A general harmonic rule controller for run-to-run process control," IEEE Transactions on Semiconductor Manufacturing , vol. 22, no. 2, pp. 232-244, 2009.

[8] Z. K. Nagy, R. D. Braatz, "Open-loop and closed-loop robust optimal control of batch processes using distributional and worst-case analysis," Journal of Process Control , vol. 14, no. 4, pp. 411-422, 2004.

[9] J. Shi, S. Zhou, "Quality control and improvement for multistage systems: a survey," IIE Transactions (Institute of Industrial Engineers) , vol. 41, no. 9, pp. 744-753, 2009.

[10] E. del Castillo, "Statistical process adjustment: a brief retrospective, current status, and some opportunities for further work," Statistica Neerlandica , vol. 60, no. 3, pp. 309-326, 2006.

[11] A.-J. Su, J.-C. Jeng, H.-P. Huang, C.-C. Yu, S.-Y. Hung, C.-K. Chao, "Control relevant issues in semiconductor manufacturing: overview with some new results," Control Engineering Practice , vol. 15, no. 10, pp. 1268-1279, 2007.

[12] F. E. Grubbs, "An optimum procedure for setting machines or adjusting processes," Journal of Quality Technology , vol. 15, pp. 186-189, 1983.

[13] D. Trietsch, "The harmonic rule for process setup adjustment with quadratic loss," Journal of Quality Technology , vol. 30, no. 1, pp. 75-84, 1998.

[14] S. A. Vander Wiel, W. T. Tucker, F. W. Faltin, N. Doganaksoy, "Algorithmic statistical process control. Concepts and an application," Technometrics , vol. 34, no. 3, pp. 286-291, 1992.

[15] A. Ingolfsson, E. Sachs, "Stability and sensitivity of an EWMA controller," Journal of Quality Technology , vol. 25, pp. 271-287, 1993.

[16] G. E. P. Box, G. M. Jenkins, G. C. Reinsel Time Series Analysis: Forecasting and Control , Prentice Hall, Upper Saddle River, NJ, USA, 1994., 3rd.

[17] F. He, H. Xie, K. Wang, "Optimal setup adjustment and control of a process under ARMA disturbances," IIE Transactions , vol. 47, no. 3, pp. 230-244, 2015.

[18] J. F. MacGregor, T. J. Harris, J. D. Wright, "Duality between the control of processes subject to randomly occurring deterministic disturbances and ARIMA stochastic disturbances," Technometrics , vol. 26, no. 4, pp. 389-397, 1984.

[19] F. Tsung, H. Wu, V. N. Nair, "On the efficiency and robustness of discrete proportional-integral control schemes," Technometrics , vol. 40, no. 3, pp. 214-222, 1998.

[20] D. W. Apley, "A cautious minimum variance controller with ARIMA disturbances," IIE Transactions , vol. 36, no. 5, pp. 417-432, 2004.

[21] M.-D. Ma, C.-C. Chang, S.-S. Jang, D. S.-H. Wong, "Mixed product run-to-run process control-an ANOVA model with ARIMA disturbance approach," Journal of Process Control , vol. 19, no. 4, pp. 604-614, 2009.

[22] P. J. Brockwell, R. A. Davis Time Series: Theory and Methods , Springer, New York, NY, USA, 2009.

[23] K. J. Åström Introduction To Stochastic Control Theory , of Mathematics in Science and Engineering, Vol. 70, Academic Press, New York, NY, USA, 1970.

[24] F. L. Lewis Optimal Estimation With An Introduction To Stochastic Control Theory , John Wiley & Sons, New York, NY, USA, 1986.

[25] E. del Castillo Statistical Process Adjustment for Quality Control , of Wiley Series in Probability and Statistics, John Wiley & Sons, New York, NY, USA, 2002.

[26] K. J. Åströmand, T. Hagglund Advanced PID Control , ISA, Research Triangle Park, NC, USA, 2006.

[27] R Core Team http://www.R-project.org/ R: A Language and Environment for Statistical Computing , R Foundation for Statistical Computing, Vienna, Austria, 2014.

[28] J. A. Nelder, R. Mead, "A simplex algorithm for function minimization," Computer Journal , vol. 7, pp. 308-313, 1965.

[29] C. G. Broyden, "The convergence of a class of double-rank minimization algorithms," Journal of the Institute of Mathematics and Its Applications , vol. 6, pp. 76-90, 1970.

[30] R. Fletcher, "New approach to variable metric algorithms," Computer Journal , vol. 13, no. 3, pp. 317-322, 1970.

[31] D. Goldfarb, "A family of variable-metric methods derived by variational means," Mathematics of Computation , vol. 24, pp. 23-26, 1970.

[32] D. F. Shanno, "Conditioning of quasi-Newton methods for function minimization," Mathematics of Computation , vol. 24, pp. 647-656, 1970.

Appendix

Appendices

A. Proof of Theorem 1

In the following proof, all the notations are the same as defined in Theorem 1. For the general ARIMA( [figure omitted; refer to PDF] ) process [figure omitted; refer to PDF] in (2), Brockwell and Davis [22] already provided its state-space representation. Referring to Example [figure omitted; refer to PDF] on page 471 in the book [22], we can get [figure omitted; refer to PDF] 's state-space representation as follows: [figure omitted; refer to PDF]

Now, we will give the state-space representation of the linear system (1)-(5). Based on (1) and (A.1), we can get [figure omitted; refer to PDF] So, we proved (8) in Theorem 1. Then, based on (4) and (5), we can get [figure omitted; refer to PDF] Combining (A.4) with (A.2), we can obtain the state equation: [figure omitted; refer to PDF] Thus, we proved (9) in Theorem 1. Then, we finish the proof of Theorem 1.

B. Tuning Procedure for [figure omitted; refer to PDF] Parameters

Given each pair of [figure omitted; refer to PDF] and [figure omitted; refer to PDF] , the MSE of the first 100 outputs [figure omitted; refer to PDF] is computed. We do 1000 replicates of this procedure and use the AMSE as the performance of the [figure omitted; refer to PDF] controller. Obviously, AMSE is a function of [figure omitted; refer to PDF] and [figure omitted; refer to PDF] . Our way to determine [figure omitted; refer to PDF] and [figure omitted; refer to PDF] is to solve the following optimization problem: [figure omitted; refer to PDF] To solve (B.1), we use the function optim in the base package of R software (see [27]). Two types of numerical optimization techniques, the Nelder-Mead method (see [28]) and the Broyden-Fletcher-Goldfarb-Shanno algorithm (see [29-32]), are used, respectively. For the two different numerical optimization techniques, we get the same optimal solution that [figure omitted; refer to PDF] and [figure omitted; refer to PDF] for the ARIMA( [figure omitted; refer to PDF] ) disturbance with [figure omitted; refer to PDF] and [figure omitted; refer to PDF] investigated in Section 4; and the same optimal solution that [figure omitted; refer to PDF] and [figure omitted; refer to PDF] for the ARIMA( [figure omitted; refer to PDF] ) disturbance with [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , [figure omitted; refer to PDF] , and [figure omitted; refer to PDF] is investigated in Section 5.