Recursive Algorithms for Multivariable Output-Error-Like ARMA Systems

Abstract

This paper studies the parameter identification problems for multivariable output-error-like systems with colored noises. Based on the hierarchical identification principle, the original system is decomposed into several subsystems. However, each subsystem contains the same parameter vector, which leads to redundant computation. By taking the average of the parameter estimation vectors of each subsystem, a partially-coupled subsystem recursive generalized extended least squares (PC-S-RGELS) algorithm is presented to cut down the redundant parameter estimates. Furthermore, a partially-coupled recursive generalized extended least squares (PC-RGELS) algorithm is presented to further reduce the computational cost and the redundant estimates by using the coupling identification concept. Finally, an example indicates the effectiveness of the derived algorithms.

1. Introduction

System identification is an important branch in the field of modern control and is an important method to establish systematic mathematical models from the combination of observation data and prior knowledge [1,2,3,4,5,6,7,8], and has been applied in many fields for decades, such as controller design [9,10,11,12,13,14,15] and system analysis [16,17,18,19,20]. Parameter identification is an important part of system identification and is to estimate the parameters by using the measurable data [21,22,23,24,25,26,27]. Parameter estimation methods can be applied to many areas [28,29,30,31]. Recently, in the literature of parameter identification, Wan et al. studied the problem of parameter estimation for multivariable equation-error systems with colored noises and proposed a hierarchical gradient-based iterative identification algorithm by using the hierarchical identification principle [32]. Chen et al. transformed the time-delay rational model into an augmented model based on the redundant rule and proposed a biased compensation recursive least squares-based threshold algorithm [33]. Other identification methods can be found in [34,35,36,37,38,39,40].

Multivariable systems widely consist in practical industrial processes on account of multi-input multi-output systems can describe modern industrial process more accurately [41,42,43,44,45]. Parameter estimation of multivariable systems has attracted extensive research attention over the past decades, and many different identification approaches have been proposed to solve the parameter identification problems of multivariable systems, such as the hierarchical identification principle and the coupling identification concept [46,47]. The core idea of the hierarchical identification principle is to decompose the original model into several submodels, and to combine other approaches to estimate the parameters of the submodels [48,49]. The coupled identification methods have been derived to identify the parameters of multivariable systems and were first presented in [50]. The basic idea of the coupling identification concept is to decompose the original system into several subsystems, and to estimate the parameters based on the coupled parameter relationships between these subsystems [51,52,53].

In the field of system modeling and control, the recursive identification and iterative identification methods are basic [54,55,56,57,58,59]. The recursive least squares methods are the commonly used parameter estimation approaches among many different parameter estimation techniques [60,61,62]. Recently, the recursive least squares (RLS) algorithm is always combined other methods to identify the complex systems. For instance, Zhang et al. proposed a filtering-based two-stage RLS algorithm for a bilinear system which is described by the state space form based on the filtering technique [63]. Liu et al. sudied the parameter estimation problems of multivariate output-error autoregressive systems and derived a filtering-based auxiliary model recursive generalized least squares algorithm based on the data filtering technique and the auxiliary model identification idea [64].

Multivariable output-error-like systems are special type of multivariable systems, which contain not only multiple inputs and multiple outputs, but also more complex parameter forms, i.e., the parameter vector and the parameter matrix [65,66]. Hence, the multivariable output-error like systems can describe modern industrial process more accurately than other types of multivariable systems. Recently, for multivariate output-error systems, Wang et al. proposed a decomposition based recursive least squares identification algorithm by using the auxiliary model, and analyzed its convergence through the stochastic process theory [67]. Ding proposed a hierarchical iterative identification algorithm for multi-input output-error systems with autoregressive noise [68]. Different from the methods in [67,68], this paper studies the parameter identification problems of multivariable output-error-like (M-OE-like) systems with colored noises which is described by the autoregressive moving average (ARMA) model by means of the decomposition technique and the coupling identification concept [69,70]. The main contributions of this paper lie in the following aspects.

• Based on the hierarchical identification principle, this paper decomposes the original system into m subsystems.
• Based on the coupled relationships between subsystems, a partially-coupled recursive generalized extended least squares (PC-RGELS) algorithm is proposed to identify the parameters of M-OEARMA-like systems.
• The derived PC-RGELS algorithm has higher computation efficiency and higher estimation accuracy than the recursive generalized extended least squares (RGELS) algorithm.

This paper is organized as follows. Section 2 describes the identification model. A RGELS algorithm is proposed to give some comparisons with the proposed algorithms in Section 3. Section 4 proposes a partially-coupled recursive generalized extended least squares algorithm. Section 5 provides the numerical simulation results to illustrate the performance of the proposed algorithms. Finally, Section 6 gives some conclusions.

2. The System Description

Let us introduce some symbols. “$B=:Y$” or “$Y:=B$” stands for “B is defined as Y”; the superscript T stands for the vector/matrix transpose; the symbol ${\mathit{I}}_n$ denotes an identity matrix of size $n\times n$; ${\mathbf{1}}_n$ stands for an n-dimensional column vector whose elements are 1; the symbol ⊗ represents the Kronecker product, for example, $\mathit{C}=\left[c_{ij}\right]\in {\mathbb{R}}^{m\times n}$, $\mathit{D}=\left[d_{ij}\right]{\mathbb{R}}^{p\times q}$, $\mathit{C}\otimes \mathit{D}=\left[c_{ij}\mathit{D}\right]\in {\mathbb{R}}^{\left(mp\right)\times \left(nq\right)}$, in general $\mathit{C}\otimes \mathit{D}\ne \mathit{D}\otimes \mathit{C}$; $\mathrm{col}\left[\mathit{Y}\right]$ is defined as the vector formed by all columns of matrix $\mathit{Y}$ arranged in order, for example, $\mathit{Y}=[{\mathit{y}}_1,{\mathit{y}}_2,\cdots ,{\mathit{y}}_n]{\mathbb{R}}^{m\times n}$, $\mathrm{col}\left[\mathit{Y}\right]={[{\mathit{y}}_1^{\mathrm{T}},{\mathit{y}}_2^{\mathrm{T}},\cdots ,{\mathit{y}}_n^{\mathrm{T}}]}^{\mathrm{T}}\in {\mathbb{R}}^{mn}$; $\widehat{\mathit{Y}}\left(s\right)$ denotes the estimate of $\mathit{Y}$ at time s; the norm of a matrix (or a column vector) $\mathit{Y}$ is defined by ${\parallel \mathit{Y}\parallel }^2:=\mathrm{tr}\left[\mathit{Y}{\mathit{Y}}^{\mathrm{T}}\right]$.

Recently, for the multivariable output-error system:

$$\mathit{y}\left(s\right)=\frac{{\mathit{\Phi }}_s\left(s\right)\mathit{\theta }}{A\left(z\right)}+\mathit{v}\left(s\right),$$

(1)

where $\mathit{y}\left(s\right):={[y_1\left(s\right),y_2\left(s\right),\cdots ,y_m\left(s\right)]}^{\mathrm{T}}\in {\mathbb{R}}^m$ refers to the output vector of the system, $\mathit{v}\left(s\right)\in {\mathbb{R}}^m$ is the white noise vector with zero mean, ${\mathit{\Phi }}_s\left(s\right){\mathbb{R}}^{m\times n}$ is the information matrix consisting of the input-output data, $\mathit{\theta }\in {\mathbb{R}}^n$ is the parameter vector, $z^{-1}$ is a unit delay operator: $z^{-1}y\left(s\right)=y(s-1)$, $A\left(z\right)$ is a monic polynomial in $z^{-1}$, and

$$A\left(z\right):=1+a_1{z}^{-1}+a_2{z}^{-2}+\cdots +a_{n_a}{z}^{n_a},a_i\in \mathbb{R},$$

a coupled stochastic gradient identification algorithm has been proposed to estimate the parameters of this system [71].

Different from the system in [71], this paper considers the multivariable output-error-like system with autoregressive moving average noise:

$$\mathit{y}\left(s\right)=\frac{\mathit{Q}\left(z\right)}{\alpha \left(z\right)}\mathit{u}\left(s\right)+\frac{D\left(z\right)}{C\left(z\right)}\mathit{v}\left(s\right),$$

(2)

where the definitions of $\mathit{y}\left(s\right)$ and $\mathit{v}\left(s\right)$ are same to above, $\mathit{u}\left(s\right):={[u_1\left(s\right),u_2\left(s\right),\cdots ,u_r\left(s\right)]}^{\mathrm{T}}\in {\mathbb{R}}^r$ is the system input vector, $\alpha \left(z\right)$, $C\left(z\right)$ and $D\left(z\right)$ are monic polynomials in $z^{-1}$ and $\mathit{Q}\left(z\right)$ is a matrix-coefficient polynomial in $z^{-1}$, and defined by

$$\begin{array}{lll}\alpha \left(z\right)& :=& 1+a_1{z}^{-1}+a_2{z}^{-2}+\cdots +a_n{z}^{-n},a_i\in \mathbb{R},\\ \mathit{Q}\left(z\right)& :=& {\mathit{Q}}_1{z}^{-1}+{\mathit{Q}}_2{z}^{-2}+\cdots +{\mathit{Q}}_n{z}^{-n},{\mathit{Q}}_i\in {\mathbb{R}}^{m\times r},\\ C\left(z\right)& :=& 1+c_1{z}^{-1}+c_2{z}^{-2}+\cdots +c_{n_c}{z}^{n_c},c_i\in \mathbb{R},\\ D\left(z\right)& :=& 1+d_1{z}^{-1}+d_2{z}^{-2}+\cdots +d_{n_d}{z}^{n_d},d_i\in \mathbb{R}.\end{array}$$

In order to focus on the essence of the parameter estimation, we assume that the orders m, r, n, $n_c$ and $n_d$ are known, and $\mathit{y}\left(s\right)=\mathbf{0}$, $\mathit{u}\left(s\right)=\mathbf{0}$ and $\mathit{v}\left(s\right)=\mathbf{0}$ for $s⩽0$. Define the actual output vector of the system,

$$\mathit{x}\left(s\right):=\frac{\mathit{Q}\left(z\right)}{\alpha \left(z\right)}\mathit{u}\left(s\right)\in {\mathbb{R}}^m.$$

(3)

Define the parameter vector $\mathit{\alpha }$ and the parameter matrix $\mathit{\theta }$, and the information vector $\mathit{\phi }\left(s\right)$ and the information matrix ${\mathit{\psi }}_s\left(s\right)$ as

$$\begin{array}{ccc}\mathit{\alpha }& :=& {[a_1,a_2,\cdots ,a_n]}^{\mathrm{T}}\in {\mathbb{R}}^n,\hfill \\ {\mathit{\theta }}^{\mathrm{T}}& :=& [{\mathit{Q}}_1,{\mathit{Q}}_2,\cdots ,{\mathit{Q}}_n]\in {\mathbb{R}}^{m\times \left(rn\right)},\hfill \\ \mathit{\phi }\left(s\right)& :=& {[{\mathit{u}}^{\mathrm{T}}(s-1),{\mathit{u}}^{\mathrm{T}}(s-2),\cdots ,{\mathit{u}}^{\mathrm{T}}(s-n)]}^{\mathrm{T}}\in {\mathbb{R}}^{rn},\hfill \\ {\mathit{\psi }}_s\left(s\right)& :=& [-\mathit{x}(s-1),-\mathit{x}(s-2),\cdots ,-\mathit{x}(s-n)]\in {\mathbb{R}}^{m\times n}.\hfill \end{array}$$

Equation (3) can be expressed as

$$\begin{array}{lll}\mathit{x}\left(s\right)& =& [1-\alpha (z\left)\right]\mathit{x}\left(s\right)+\mathit{Q}\left(z\right)\mathit{u}\left(s\right),\\ & =& {\mathit{\psi }}_s\left(s\right)\mathit{\alpha }+{\mathit{\theta }}^{\mathrm{T}}\mathit{\phi }\left(s\right).\end{array}$$

(4)

Define an internal vector of the system,

$$\mathit{w}\left(s\right):=\frac{D\left(z\right)}{C\left(z\right)}\mathit{v}\left(s\right)\in {\mathbb{R}}^m.$$

(5)

Let $n_1:=n_c+n_d$, define the parameter vector $\mathit{\rho }$ and the information matrix ${\mathit{\psi }}_n\left(s\right)$ of the system as

$$\begin{array}{ccc}\mathit{\rho }& :=& {[c_1,c_2,\cdots ,c_{n_c},d_1,d_2,\cdots ,d_{n_d}]}^{\mathrm{T}}\in {\mathbb{R}}^{n_1},\hfill \\ {\mathit{\psi }}_n\left(s\right)& :=& [-\mathit{w}(s-1),-\mathit{w}(s-2),\cdots ,-\mathit{w}(s-n_c),\mathit{v}(s-1),\mathit{v}(s-2),\cdots ,\mathit{v}(s-n_d)]\in {\mathbb{R}}^{m\times n_1}.\end{array}$$

Equation (5) can be rewritten as

$$\begin{array}{lll}\mathit{w}\left(s\right)& =& [1-C(z\left)\right]\mathit{w}\left(s\right)+D\left(z\right)\mathit{v}\left(s\right),\end{array}$$

$$\begin{array}{lll}& =& \mathit{y}\left(s\right)-{\mathit{\psi }}_s\left(s\right)\mathit{\alpha }-{\mathit{\theta }}^{\mathrm{T}}\mathit{\phi }\left(s\right),\end{array}$$

(6)

$$\begin{array}{lll}& =& {\mathit{\psi }}_n\left(s\right)\mathit{\rho }+\mathit{v}\left(s\right).\end{array}$$

(7)

Using (4) and (7), Equation (2) can be equivalently written as

$$\begin{array}{lll}\mathit{y}\left(s\right)& =& \mathit{x}\left(s\right)+\mathit{w}\left(s\right),\end{array}$$

$$\begin{array}{lll}& =& {\mathit{\psi }}_s\left(s\right)\mathit{\alpha }+{\mathit{\theta }}^{\mathrm{T}}\mathit{\phi }\left(s\right)+{\mathit{\psi }}_n\left(s\right)\mathit{\rho }+\mathit{v}\left(s\right),\end{array}$$

(8)

$$\begin{array}{lll}& =& [{\mathit{\psi }}_s\left(s\right),{\mathit{\psi }}_n\left(s\right)]\left[\begin{array}{c}\mathit{\alpha }\\ \mathit{\rho }\end{array}\right]+{\mathit{\theta }}^{\mathrm{T}}\mathit{\phi }\left(s\right)+\mathit{v}\left(s\right).\end{array}$$

(9)

Let $n_0:=n_1+n$, define the information matrix $\mathit{\psi }\left(s\right)$ and the parameter vector $\mathit{\beta }$ as

$$\begin{array}{ccc}\mathit{\psi }\left(s\right)& :=& [{\mathit{\psi }}_s\left(s\right),{\mathit{\psi }}_n\left(s\right)]\in {\mathbb{R}}^{m\times n_0},\end{array}$$

(10)

$$\begin{array}{ccc}\mathit{\beta }& :=& {[{\mathit{\alpha }}^{\mathrm{T}},{\mathit{\rho }}^{\mathrm{T}}]}^{\mathrm{T}}\in {\mathbb{R}}^{n_0}.\end{array}$$

(11)

Substituting (10) and (11) into (9) gives the hierarchical identification model

$$\mathit{y}\left(s\right)=\mathit{\psi }\left(s\right)\mathit{\beta }+{\mathit{\theta }}^{\mathrm{T}}\mathit{\phi }\left(s\right)+\mathit{v}\left(s\right).$$

(12)

For convenience, we define an information matrix $\mathit{\Psi }\left(s\right)$ by making use of the Kronecker product of the information matrix $\mathit{\psi }\left(s\right)$ and the information vector $\mathit{\phi }\left(s\right)$ as

$$\mathit{\Psi }\left(s\right):=[\mathit{\psi }\left(s\right),{\mathit{\phi }}^{\mathrm{T}}\left(s\right)\otimes {\mathit{I}}_m]\in {\mathbb{R}}^{m\times n_2},n_2:=n_0+mrn.$$

Hence, a parameter vector $\mathit{\vartheta }$ is defined by using the parameter vector $\mathit{\beta }$ and the parameter matrix $\mathit{\theta }$,

$$\mathit{\vartheta }:=\left[\begin{array}{c}\mathit{\beta }\\ \mathrm{col}\left[{\mathit{\theta }}^{\mathrm{T}}\right]\end{array}\right]\in {\mathbb{R}}^{n_2}.$$

Then Equation (12) can be equivalently expressed as

$$\mathit{y}\left(s\right)=\mathit{\Psi }\left(s\right)\mathit{\vartheta }+\mathit{v}\left(s\right).$$

(13)

Therefore, we get the identification model (13) of the M-OEARMA-like system in (2), where $\mathit{\vartheta }$ is the parameter vector to be identified and contains all the parameters of the system (2)

3. The RGELS Algorithm

Based on Equation (13), define a criterion function,

$$J_1\left(\mathit{\vartheta }\right):={\displaystyle \sum _{j=1}^s}{\parallel \mathit{y}\left(j\right)-\mathit{\Psi }\left(j\right)\mathit{\vartheta }\parallel }^2.$$

(14)

Let $\widehat{\mathit{\vartheta }}\left(s\right)$ be the estimate of $\mathit{\vartheta }$ at time s. Minimizing $J_1\left(\mathit{\vartheta }\right)$ gives

$$\frac{\partial J_1\left(\mathit{\vartheta }\right)}{\partial \mathit{\vartheta }}{|}_{\mathit{\vartheta }=\widehat{\mathit{\vartheta }}\left(s\right)}=\frac{\partial {({\mathit{Y}}_s-{\mathit{H}}_s\mathit{\vartheta })}^{\mathrm{T}}({\mathit{Y}}_s-{\mathit{H}}_s\mathit{\vartheta })}{\partial \mathit{\vartheta }}{|}_{\mathit{\vartheta }=\widehat{\mathit{\vartheta }}\left(s\right)}=\mathbf{0},$$

(15)

where

$$\begin{array}{lll}{\mathit{Y}}_s& :=& {[{\mathit{y}}^{\mathrm{T}}\left(1\right),{\mathit{y}}^{\mathrm{T}}\left(2\right),{\mathit{y}}^{\mathrm{T}}\left(3\right),\cdots ,{\mathit{y}}^{\mathrm{T}}\left(s\right)]}^{\mathrm{T}}\in {\mathbb{R}}^{ms},\\ {\mathit{H}}_s& :=& {[{\mathit{\Psi }}^{\mathrm{T}}\left(1\right),{\mathit{\Psi }}^{\mathrm{T}}\left(2\right),{\mathit{\Psi }}^{\mathrm{T}}\left(3\right),\cdots ,{\mathit{\Psi }}^{\mathrm{T}}\left(s\right)]}^{\mathrm{T}}\in {\mathbb{R}}^{\left(ms\right)\times n_2},\end{array}$$

Parameter estimate $\widehat{\mathit{\vartheta }}\left(s\right)$ of $\mathit{\vartheta }$ can be obtained from (15) as

$$\begin{array}{lll}\widehat{\mathit{\vartheta }}\left(s\right)& =& {\left({\mathit{H}}_s^{\mathrm{T}}{\mathit{H}}_s\right)}^{-1}{\mathit{H}}_s^{\mathrm{T}}{\mathit{Y}}_s,\end{array}$$

(16)

$$\begin{array}{lll}& =& {\left[{\displaystyle \sum _{j=1}^s}{\mathit{\Psi }}^{\mathrm{T}}\left(j\right)\mathit{\Psi }\left(j\right)\right]}^{-1}\left[{\displaystyle \sum _{j=1}^s}{\mathit{\Psi }}^{\mathrm{T}}\left(j\right)\mathit{y}\left(j\right)\right].\end{array}$$

(17)

Define the covariance matrix

$$\begin{array}{lll}{\mathit{P}}^{-1}\left(s\right)& :=& {\displaystyle \sum _{j=1}^s}{\mathit{\Psi }}^{\mathrm{T}}\left(j\right)\mathit{\Psi }\left(j\right)\in {\mathbb{R}}^{n_2\times n_2},\\ & =& {\mathit{P}}^{-1}(s-1)+{\mathit{\Psi }}^{\mathrm{T}}\left(s\right)\mathit{\Psi }\left(s\right).\end{array}$$

(18)

Let $\mathit{L}\left(s\right):=\mathit{P}\left(s\right){\mathit{\Psi }}^{\mathrm{T}}\left(s\right)\in {\mathbb{R}}^{n_2\times m}$ be the gain matrix. Based on the derivation of the RLS algorithm in [72,73], we can easily get the RLS relations:

$$\begin{array}{ccc}\widehat{\mathit{\vartheta }}\left(s\right)& =& \widehat{\mathit{\vartheta }}(s-1)+\mathit{L}\left(s\right)[\mathit{y}\left(s\right)-\mathit{\Psi }\left(s\right)\widehat{\mathit{\vartheta }}(s-1)],\end{array}$$

(19)

$$\begin{array}{ccc}\mathit{L}\left(s\right)& =& \mathit{P}(s-1){\mathit{\Psi }}^{\mathrm{T}}\left(s\right){[{\mathit{I}}_m+\mathit{\Psi }\left(s\right)\mathit{P}(s-1){\mathit{\Psi }}^{\mathrm{T}}\left(s\right)]}^{-1},\end{array}$$

(20)

$$\begin{array}{ccc}\mathit{P}\left(s\right)& =& \mathit{P}(s-1)-\mathit{L}\left(s\right){\left[\mathit{P}(s-1){\mathit{\Psi }}^{\mathrm{T}}\left(s\right)\right]}^{\mathrm{T}}.\end{array}$$

(21)

However, Equations (19)–(21) cannot figure out the parameter estimate $\widehat{\mathit{\vartheta }}\left(s\right)$ because of the information matrix $\mathit{\Psi }\left(s\right)$ contains the unknown vectors $\mathit{x}(s-i)$, $\mathit{w}(s-i)$ and $\mathit{v}(s-i)$. The solution is to replace these unknown vectors in $\mathit{\Psi }\left(s\right)$ with their corresponding estimates $\widehat{\mathit{x}}(s-i)$, $\widehat{\mathit{w}}(s-i)$ and $\widehat{\mathit{v}}(s-i)$ by using the auxiliary model. Define the estimates of $\mathit{\Psi }\left(s\right)$, $\mathit{\psi }\left(s\right)$, ${\mathit{\psi }}_s\left(s\right)$ and ${\mathit{\psi }}_n\left(s\right)$ as

$$\begin{array}{ccc}\widehat{\mathit{\Psi }}\left(s\right)& :=& [\widehat{\mathit{\psi }}\left(s\right),{\mathit{\phi }}^{\mathrm{T}}\left(s\right)\otimes {\mathit{I}}_m]\in {\mathbb{R}}^{m\times n_2},\end{array}$$

(22)

$$\begin{array}{ccc}\widehat{\mathit{\psi }}\left(s\right)& :=& [{\widehat{\mathit{\psi }}}_s\left(s\right),{\widehat{\mathit{\psi }}}_n\left(s\right)]\in {\mathbb{R}}^{m\times n_0},\end{array}$$

(23)

$$\begin{array}{ccc}{\widehat{\mathit{\psi }}}_s\left(s\right)& :=& [-\widehat{\mathit{x}}(s-1),-\widehat{\mathit{x}}(s-2),\cdots ,-\widehat{\mathit{x}}(s-n)]\in {\mathbb{R}}^{m\times n},\end{array}$$

(24)

$$\begin{array}{ccc}{\widehat{\mathit{\psi }}}_n\left(s\right)& :=& [-\widehat{\mathit{w}}(s-1),\cdots ,-\widehat{\mathit{w}}(s-n_c),\widehat{\mathit{v}}(s-1),\cdots ,\widehat{\mathit{v}}(s-n_d)]\in {\mathbb{R}}^{m\times n_1}.\end{array}$$

(25)

Replacing ${\mathit{\psi }}_s\left(s\right)$, $\mathit{\alpha }$ and $\mathit{\theta }$ in (4) and (6) with their estimates ${\widehat{\mathit{\psi }}}_s\left(s\right)$, $\widehat{\mathit{\alpha }}\left(s\right)$ and $\widehat{\mathit{\theta }}\left(s\right)$, the estimates $\widehat{\mathit{x}}\left(s\right)$ and $\widehat{\mathit{w}}\left(s\right)$ can be calculated by two auxiliary models:

$$\begin{array}{ccc}\widehat{\mathit{x}}\left(s\right)& :=& {\widehat{\mathit{\psi }}}_s\left(s\right)\widehat{\mathit{\alpha }}\left(s\right)+{\widehat{\mathit{\theta }}}^{\mathrm{T}}\left(s\right)\mathit{\phi }\left(s\right),\end{array}$$

(26)

$$\begin{array}{ccc}\widehat{\mathit{w}}\left(s\right)& :=& \mathit{y}\left(s\right)-\widehat{\mathit{x}}\left(s\right)=\mathit{y}\left(s\right)-{\widehat{\mathit{\psi }}}_s\left(s\right)\widehat{\mathit{\alpha }}\left(s\right)-{\widehat{\mathit{\theta }}}^{\mathrm{T}}\left(s\right)\mathit{\phi }\left(s\right).\end{array}$$

(27)

From (13), use the estimates $\widehat{\mathit{\Psi }}\left(s\right)$ and $\widehat{\mathit{\vartheta }}\left(s\right)$ of $\mathit{\Psi }\left(s\right)$ and $\mathit{\vartheta }$ to define the estimate of $\mathit{v}\left(s\right)$ as

$$\widehat{\mathit{v}}\left(s\right):=\mathit{y}\left(s\right)-\widehat{\mathit{\Psi }}\left(s\right)\widehat{\mathit{\vartheta }}\left(s\right).$$

(28)

Combining (22)–(28) and replacing $\mathit{\Psi }\left(s\right)$ in (19)–(21) with $\widehat{\mathit{\Psi }}\left(s\right)$ yield the following recursive generalized extended least squares (RGELS) algorithm:

$$\begin{array}{ccc}\widehat{\mathit{\vartheta }}\left(s\right)& =& \widehat{\mathit{\vartheta }}(s-1)+\mathit{L}\left(s\right)[\mathit{y}\left(s\right)-\widehat{\mathit{\Psi }}\left(s\right)\widehat{\mathit{\vartheta }}(s-1)],\end{array}$$

(29)

$$\begin{array}{ccc}\mathit{L}\left(s\right)& =& \mathit{P}(s-1){\widehat{\mathit{\Psi }}}^{\mathrm{T}}\left(s\right){[{\mathit{I}}_m+\widehat{\mathit{\Psi }}\left(s\right)\mathit{P}(s-1){\widehat{\mathit{\Psi }}}^{\mathrm{T}}\left(s\right)]}^{-1},\end{array}$$

(30)

$$\begin{array}{ccc}\mathit{P}\left(s\right)& =& \mathit{P}(s-1)-\mathit{L}\left(s\right){\left[\mathit{P}(s-1){\widehat{\mathit{\Psi }}}^{\mathrm{T}}\left(s\right)\right]}^{\mathrm{T}},\end{array}$$

(31)

$$\begin{array}{ccc}\widehat{\mathit{\Psi }}\left(s\right)& =& [\widehat{\mathit{\psi }}\left(s\right),{\mathit{\phi }}^{\mathrm{T}}\left(s\right)\otimes {\mathit{I}}_m],\end{array}$$

(32)

$$\begin{array}{ccc}\mathit{\phi }\left(s\right)& =& {[{\mathit{u}}^{\mathrm{T}}(s-1),{\mathit{u}}^{\mathrm{T}}(s-2),\cdots ,{\mathit{u}}^{\mathrm{T}}(s-n)]}^{\mathrm{T}},\end{array}$$

(33)

$$\begin{array}{ccc}\widehat{\mathit{\psi }}\left(s\right)& =& [{\widehat{\mathit{\psi }}}_s\left(s\right),{\widehat{\mathit{\psi }}}_n\left(s\right)],\end{array}$$

(34)

$$\begin{array}{ccc}{\widehat{\mathit{\psi }}}_s\left(s\right)& =& [-\widehat{\mathit{x}}(s-1),-\widehat{\mathit{x}}(s-2),\cdots ,-\widehat{\mathit{x}}(s-n)],\end{array}$$

(35)

$$\begin{array}{ccc}{\widehat{\mathit{\psi }}}_n\left(s\right)& =& [-\widehat{\mathit{w}}(s-1),-\widehat{\mathit{w}}(s-2),\cdots ,-\widehat{\mathit{w}}(s-n_c),\widehat{\mathit{v}}(s-1),\widehat{\mathit{v}}(s-2),\cdots ,\widehat{\mathit{v}}(s-n_d)],\end{array}$$

(36)

$$\begin{array}{ccc}\widehat{\mathit{x}}\left(s\right)& =& {\widehat{\mathit{\psi }}}_s\left(s\right)\widehat{\mathit{\alpha }}\left(s\right)+{\widehat{\mathit{\theta }}}^{\mathrm{T}}\left(s\right)\mathit{\phi }\left(s\right),\end{array}$$

(37)

$$\begin{array}{ccc}\widehat{\mathit{w}}\left(s\right)& =& \mathit{y}\left(s\right)-{\widehat{\mathit{\psi }}}_s\left(s\right)\widehat{\mathit{\alpha }}\left(s\right)-{\widehat{\mathit{\theta }}}^{\mathrm{T}}\left(s\right)\mathit{\phi }\left(s\right),\end{array}$$

(38)

$$\begin{array}{ccc}\widehat{\mathit{v}}\left(s\right)& =& \mathit{y}\left(s\right)-\widehat{\mathit{\Psi }}\left(s\right)\widehat{\mathit{\vartheta }}\left(s\right).\end{array}$$

(39)

The procedure contained in the RGELS algorithm in (29)–(39) as follows.

• For $s⩽0$, all variables are set to zero. Set the data length L. Let $s=1$, set the initial values $\widehat{\mathit{\vartheta }}\left(0\right)={\mathbf{1}}_{n_2}/p_0$, $\widehat{\mathit{x}}\left(0\right)={\mathbf{1}}_m/p_0$, $\widehat{\mathit{w}}\left(0\right)={\mathbf{1}}_m/p_0$, $\widehat{\mathit{v}}\left(0\right)={\mathbf{1}}_m/p_0$, $\mathit{P}\left(0\right)=p_0{\mathit{I}}_{n_3}$, $p_0={10}^6$.
• Collect the input-output data $\mathit{u}\left(s\right)$ and $\mathit{y}\left(s\right)$, and construct $\mathit{\phi }\left(s\right)$ using (33).
• Form ${\widehat{\mathit{\psi }}}_s\left(s\right)$, ${\widehat{\mathit{\psi }}}_n\left(s\right)$ and $\widehat{\mathit{\psi }}\left(s\right)$ using (35)–(36) and (34), and form $\widehat{\mathit{\Psi }}\left(s\right)$ using (32).
• Calculate the covariance matrix $\mathit{P}\left(s\right)$ and the gain matrix $\mathit{L}\left(s\right)$ using (31) and (30), and update the estimate $\widehat{\mathit{\vartheta }}\left(s\right)$ using (29).
• Figure the estimates $\widehat{\mathit{x}}\left(s\right)$, $\widehat{\mathit{w}}\left(s\right)$ and $\widehat{\mathit{v}}\left(s\right)$ using (37)–(39).
• Compare s with L: if $s⩽L$, increase s by 1 and go to Step 2; otherwise obtain the parameter estimate $\widehat{\mathit{\vartheta }}\left(L\right)$ of $\mathit{\vartheta }$ and break up the program.

The flowchart of computing $\widehat{\mathit{\vartheta }}\left(L\right)$ in the RGELS algorithm is shown in Figure 1. The RGELS algorithm is basic in system identification, and can be extended to study the parameter estimation problems of different systems such as signal modeling and communication networked systems [74,75,76,77,78,79].

4. The PC-RGELS Algorithm

In this part, a partially-coupled recursive generalized extended least squares (PC-RGELS) identification algorithm is studied to cut down the redundant estimates and improve the computational efficiency of the RGELS algorithm based on the decomposition technique and the coupling identification concept.

The identification model in (12) of system (2) is rewritten as follows:

$$\mathit{y}\left(s\right)=\mathit{\psi }\left(s\right)\mathit{\beta }+{\mathit{\theta }}^{\mathrm{T}}\mathit{\phi }\left(s\right)+\mathit{v}\left(s\right).$$

(40)

Referring to the decomposition methods in [51,52,71], let ${\mathit{\psi }}_i^{\mathrm{T}}\left(s\right)\in {\mathbb{R}}^{n_0}$ be the ith row of the information matrix $\mathit{\psi }\left(s\right)$, that is

$$\mathit{\psi }\left(s\right):={[{\mathit{\psi }}_1\left(s\right),{\mathit{\psi }}_2\left(s\right),\cdots ,{\mathit{\psi }}_m\left(s\right)]}^{\mathrm{T}}\in {\mathbb{R}}^{m\times n_0}.$$

Similarly, let ${\mathit{\theta }}_i\in {\mathbb{R}}^{rn}$ be the ith column of the parameter matrix $\mathit{\theta }$:

$$\mathit{\theta }:=[{\mathit{\theta }}_1,{\mathit{\theta }}_2,\cdots ,{\mathit{\theta }}_m]\in {\mathbb{R}}^{\left(rn\right)\times m}.$$

Then Equation (40) can be decomposed into m subsystem identification models:

$$\begin{array}{lll}{y}_i\left(s\right)& =& {\mathit{\psi }}_i^{\mathrm{T}}\left(s\right)\mathit{\beta }+{\mathit{\theta }}_i^{\mathrm{T}}\mathit{\phi }\left(s\right)+v_i\left(s\right),\\ & =& {\mathit{\psi }}_i^{\mathrm{T}}\left(s\right)\mathit{\beta }+{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\mathit{\theta }}_i+v_i\left(s\right),i=1,2,\cdots ,m.\end{array}$$

(41)

According to the identification model in (41), define a gradient criterion function,

$$\begin{array}{c}{J}_2(\mathit{\beta },{\mathit{\theta }}_i):={\displaystyle \sum _{j=1}^s}{[y_i\left(j\right)-{\mathit{\psi }}_i^{\mathrm{T}}\left(j\right)\mathit{\beta }-{\mathit{\phi }}^{\mathrm{T}}\left(j\right){\mathit{\theta }}_i]}^2.\end{array}$$

(42)

Let $\widehat{\mathit{\beta }}\left(s\right)$ and ${\widehat{\mathit{\theta }}}_i\left(s\right)$ be the estimates of $\mathit{\beta }$ and ${\mathit{\theta }}_i$ at time s. Minimizing $J_2(\mathit{\beta },{\mathit{\theta }}_i)$ gives

$$\begin{array}{ccc}\frac{\partial J_2(\mathit{\beta },{\mathit{\theta }}_i)}{\partial \mathit{\beta }}{|}_{\mathit{\beta }=\widehat{\mathit{\beta }}\left(s\right)}& =& \frac{\partial \parallel {\mathit{Y}}_{i,s}-{\mathit{H}}_{i,s}\mathit{\beta }-{\mathit{Z}}_{i,s}{\mathit{\theta }}_i{\parallel }^2}{\partial \mathit{\beta }}{|}_{\mathit{\beta }=\widehat{\mathit{\beta }}\left(s\right)}=\mathbf{0},\end{array}$$

(43)

$$\begin{array}{ccc}\frac{\partial J_2(\mathit{\beta },{\mathit{\theta }}_i)}{\partial {\mathit{\theta }}_i}{|}_{{\mathit{\theta }}_i={\widehat{\mathit{\theta }}}_i\left(s\right)}& =& \frac{\partial \parallel {\mathit{Y}}_{i,s}-{\mathit{H}}_{i,s}\mathit{\beta }-{\mathit{Z}}_{i,s}{\mathit{\theta }}_i{\parallel }^2}{\partial {\mathit{\theta }}_i}{|}_{{\mathit{\theta }}_i={\widehat{\mathit{\theta }}}_i\left(s\right)}=\mathbf{0},\end{array}$$

(44)

where

$$\begin{array}{lll}{\mathit{Y}}_{i,s}& :=& {[y_i\left(1\right),y_i\left(2\right),y_i\left(3\right),\cdots ,y_i\left(s\right)]}^{\mathrm{T}}\in {\mathbb{R}}^s,\\ {\mathit{H}}_{i,s}& :=& {[{\mathit{\psi }}_i\left(1\right),{\mathit{\psi }}_i\left(2\right),{\mathit{\psi }}_i\left(3\right),\cdots ,{\mathit{\psi }}_i\left(s\right)]}^{\mathrm{T}}\in {\mathbb{R}}^{s\times n_0},\\ {\mathit{Z}}_{i,s}& :=& {[\mathit{\phi }\left(1\right),\mathit{\phi }\left(2\right),\mathit{\phi }\left(3\right),\cdots ,\mathit{\phi }\left(s\right)]}^{\mathrm{T}}\in {\mathbb{R}}^{s\times \left(rn\right)},\end{array}$$

From (43) and (44), we can get the least squares estimates $\widehat{\mathit{\beta }}\left(s\right)$ and ${\widehat{\mathit{\theta }}}_i\left(s\right)$ of $\mathit{\beta }$ and ${\mathit{\theta }}_i$:

$$\begin{array}{lll}\widehat{\mathit{\beta }}\left(s\right)& =& {\left({\mathit{H}}_{i,s}^{\mathrm{T}}{\mathit{H}}_{i,s}\right)}^{-1}({\mathit{H}}_{i,s}^{\mathrm{T}}{\mathit{Y}}_{i,s}-{\mathit{H}}_{i,s}^{\mathrm{T}}{\mathit{Z}}_{i,s}{\mathit{\theta }}_i),\end{array}$$

(45)

$$\begin{array}{lll}& =& {\left[{\displaystyle \sum _{j=1}^s}{\mathit{\psi }}_i\left(j\right){\mathit{\psi }}_i^{\mathrm{T}}\left(j\right)\right]}^{-1}\left[{\displaystyle \sum _{j=1}^s}{\mathit{\psi }}_i\left(j\right)y_i\left(j\right)-{\displaystyle \sum _{j=1}^s}{\mathit{\psi }}_i\left(j\right){\mathit{\phi }}^{\mathrm{T}}\left(j\right){\mathit{\theta }}_i\right],\end{array}$$

(46)

$$\begin{array}{lll}{\widehat{\mathit{\theta }}}_i\left(s\right)& =& {\left({\mathit{Z}}_{i,s}^{\mathrm{T}}{\mathit{Z}}_{i,s}\right)}^{-1}({\mathit{Z}}_{i,s}^{\mathrm{T}}{\mathit{Y}}_{i,s}-{\mathit{Z}}_{i,s}^{\mathrm{T}}{\mathit{H}}_{i,s}\mathit{\beta }),\end{array}$$

(47)

$$\begin{array}{lll}& =& {\left[{\displaystyle \sum _{j=1}^s}\mathit{\phi }\left(j\right){\mathit{\phi }}^{\mathrm{T}}\left(j\right)\right]}^{-1}\left[{\displaystyle \sum _{j=1}^s}\mathit{\phi }\left(j\right)y_i\left(j\right)-{\displaystyle \sum _{j=1}^s}\mathit{\phi }\left(j\right){\mathit{\psi }}_i^{\mathrm{T}}\left(j\right)\mathit{\beta }\right].\end{array}$$

(48)

Define the covariance matrixes ${\mathit{P}}_{\beta ,i}\left(s\right)$ and ${\mathit{P}}_{\theta ,i}\left(s\right)$, and the gain matrixes ${\mathit{L}}_{\beta ,i}\left(s\right)$ and ${\mathit{L}}_{\theta ,i}\left(s\right)$ as

$$\begin{array}{ccc}{\mathit{P}}_{\beta ,i}^{-1}\left(s\right)& :=& {\displaystyle \sum _{j=1}^s}{\mathit{\psi }}_i\left(j\right){\mathit{\psi }}_i^{\mathrm{T}}\left(j\right)\in {\mathbb{R}}^{n_0\times n_0},\end{array}$$

(49)

$$\begin{array}{ccc}& =& {\mathit{P}}_{\beta ,i}^{-1}(s-1)+{\mathit{\psi }}_i\left(s\right){\mathit{\psi }}_i^{\mathrm{T}}\left(s\right),\end{array}$$

(50)

$$\begin{array}{ccc}{\mathit{P}}_{\theta ,i}^{-1}\left(s\right)& :=& {\displaystyle \sum _{j=1}^s}\mathit{\phi }\left(j\right){\mathit{\phi }}^{\mathrm{T}}\left(j\right)\in {\mathbb{R}}^{\left(rn\right)\times \left(rn\right)},\end{array}$$

(51)

$$\begin{array}{ccc}& =& {\mathit{P}}_{\theta ,i}^{-1}(s-1)+\mathit{\phi }\left(s\right){\mathit{\phi }}^{\mathrm{T}}\left(s\right),\end{array}$$

(52)

$$\begin{array}{ccc}{\mathit{L}}_{\beta ,i}\left(s\right)& :=& {\mathit{P}}_{\beta ,i}\left(s\right){\mathit{\psi }}_i\left(s\right)\in {\mathbb{R}}^{n_0},\end{array}$$

(53)

$$\begin{array}{ccc}{\mathit{L}}_{\theta ,i}\left(s\right)& :=& {\mathit{P}}_{\theta ,i}\left(s\right)\mathit{\phi }\left(s\right)\in {\mathbb{R}}^{rn}.\end{array}$$

(54)

Based on the derivation of the RLS algorithm in [72,73], we can summarize the RLS estimates $\widehat{\mathit{\beta }}\left(s\right)$ and ${\widehat{\mathit{\theta }}}_i\left(s\right)$ of ${\mathit{\beta }}_i$ and ${\mathit{\theta }}_i$:

$$\begin{array}{ccc}\widehat{\mathit{\beta }}\left(s\right)& =& \widehat{\mathit{\beta }}(s-1)+{\mathit{L}}_{\beta ,i}\left(s\right)[y_i\left(s\right)-{\mathit{\psi }}_i^{\mathrm{T}}\left(s\right)\widehat{\mathit{\beta }}(s-1)-{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\widehat{\mathit{\theta }}}_i(s-1)],\end{array}$$

(55)

$$\begin{array}{ccc}{\mathit{L}}_{\beta ,i}\left(s\right)& =& \frac{{\mathit{P}}_{\beta ,i}(s-1){\mathit{\psi }}_i\left(s\right)}{1+{\mathit{\psi }}_i^{\mathrm{T}}\left(s\right){\mathit{P}}_{\beta ,i}(s-1){\mathit{\psi }}_i\left(s\right)},\end{array}$$

(56)

$$\begin{array}{ccc}{\mathit{P}}_{\beta ,i}\left(s\right)& =& {\mathit{P}}_{\beta ,i}(s-1)-{\mathit{L}}_{\beta ,i}\left(s\right){\left[{\mathit{P}}_{\beta ,i}(s-1){\mathit{\psi }}_i\left(s\right)\right]}^{\mathrm{T}},\end{array}$$

(57)

$$\begin{array}{ccc}{\widehat{\mathit{\theta }}}_i\left(s\right)& =& {\widehat{\mathit{\theta }}}_i(s-1)+{\mathit{L}}_{\theta ,i}\left(s\right)[y_i\left(s\right)-{\mathit{\psi }}_i^{\mathrm{T}}\left(s\right)\widehat{\mathit{\beta }}(s-1)-{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\widehat{\mathit{\theta }}}_i(s-1)],\end{array}$$

(58)

$$\begin{array}{ccc}{\mathit{L}}_{\theta ,i}\left(s\right)& =& \frac{{\mathit{P}}_{\theta ,i}(s-1)\mathit{\phi }\left(s\right)}{1+{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\mathit{P}}_{\theta ,i}(s-1)\mathit{\phi }\left(s\right)},\end{array}$$

(59)

$$\begin{array}{ccc}{\mathit{P}}_{\theta ,i}\left(s\right)& =& {\mathit{P}}_{\theta ,i}(s-1)-{\mathit{L}}_{\theta ,i}\left(s\right){\left[{\mathit{P}}_{\theta ,i}(s-1)\mathit{\phi }\left(s\right)\right]}^{\mathrm{T}}.\end{array}$$

(60)

However, Equations (55)–(60) cannot figure out the estimates $\widehat{\mathit{\beta }}\left(s\right)$ and ${\widehat{\mathit{\theta }}}_i\left(s\right)$. Because the subsystem information vector ${\mathit{\psi }}_i\left(s\right)$ includes the unknown vectors $\mathit{x}(s-i)$, $\mathit{w}(s-i)$ and $\mathit{v}(s-i)$. In order to deal with this problem, we replace these unknown vectors in ${\mathit{\psi }}_i\left(s\right)$ with their corresponding estimates $\widehat{\mathit{x}}(s-i)$, $\widehat{\mathit{w}}(s-i)$ and $\widehat{\mathit{v}}(s-i)$ by making use of the auxiliary model identification idea. Then the estimates of $\mathit{\psi }\left(s\right)$, ${\mathit{\psi }}_s\left(s\right)$ and ${\mathit{\psi }}_n\left(s\right)$ can be constructed by

$$\begin{array}{ccc}\widehat{\mathit{\psi }}\left(s\right)& :=& [{\widehat{\mathit{\psi }}}_s\left(s\right),{\widehat{\mathit{\psi }}}_n\left(s\right)],\end{array}$$

(61)

$$\begin{array}{cc}=& {[{\widehat{\mathit{\psi }}}_1\left(s\right),{\widehat{\mathit{\psi }}}_2\left(s\right),\cdots ,{\widehat{\mathit{\psi }}}_m\left(s\right)]}^{\mathrm{T}}\in {\mathbb{R}}^{m\times n_0},\end{array}$$

(62)

$$\begin{array}{ccc}{\widehat{\mathit{\psi }}}_s\left(s\right)& :=& [-\widehat{\mathit{x}}(s-1),-\widehat{\mathit{x}}(s-2),\cdots ,-\widehat{\mathit{x}}(s-n)]\in {\mathbb{R}}^{m\times n},\end{array}$$

(63)

$$\begin{array}{ccc}{\widehat{\mathit{\psi }}}_n\left(s\right)& :=& [-\widehat{\mathit{w}}(s-1),\cdots ,-\widehat{\mathit{w}}(s-n_c),\widehat{\mathit{v}}(s-1),\cdots ,\widehat{\mathit{v}}(s-n_d)]\in {\mathbb{R}}^{m\times n_1}.\end{array}$$

(64)

Based on (4), (6) and (41), the estimates $\widehat{\mathit{x}}\left(s\right)$, $\widehat{\mathit{w}}\left(s\right)$ and ${\widehat{v}}_i\left(s\right)$ at time s can be calculated through three auxiliary models:

$$\begin{array}{ccc}\widehat{\mathit{x}}\left(s\right)& :=& {\widehat{\mathit{\psi }}}_s\left(s\right)\widehat{\mathit{\alpha }}\left(s\right)+{\widehat{\mathit{\theta }}}^{\mathrm{T}}\left(s\right)\mathit{\phi }\left(s\right),\end{array}$$

(65)

$$\begin{array}{ccc}\widehat{\mathit{w}}\left(s\right)& :=& \mathit{y}\left(s\right)-{\widehat{\mathit{\psi }}}_s\left(s\right)\widehat{\mathit{\alpha }}\left(s\right)-{\widehat{\mathit{\theta }}}^{\mathrm{T}}\left(s\right)\mathit{\phi }\left(s\right),\end{array}$$

(66)

$$\begin{array}{ccc}{\widehat{v}}_i\left(s\right)& :=& y_i\left(s\right)-{\widehat{\mathit{\psi }}}_i^{\mathrm{T}}\left(s\right)\widehat{\mathit{\beta }}\left(s\right)-{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\widehat{\mathit{\theta }}}_i\left(s\right).\end{array}$$

(67)

For clarity, let ${\widehat{\mathit{\beta }}}_i\left(s\right)$ represent the estimate of $\mathit{\beta }$ in Subsystem i at time s. Combining (61)–(67) and replacing ${\mathit{\psi }}_i\left(s\right)$ in (55)–(60) with its estimate ${\widehat{\mathit{\psi }}}_i\left(s\right)$ and the same parameter vector $\widehat{\mathit{\beta }}\left(s\right)$ in (55)–(60) with ${\widehat{\mathit{\beta }}}_i\left(s\right)$ give the subsystem recursive generalized extended least squares (S-RGELS) algorithm:

$$\begin{array}{ccc}{\widehat{\mathit{\beta }}}_i\left(s\right)& =& {\widehat{\mathit{\beta }}}_i(s-1)+{\mathit{L}}_{\beta ,i}\left(s\right)[y_i\left(s\right)-{\widehat{\mathit{\psi }}}_i^{\mathrm{T}}\left(s\right){\widehat{\mathit{\beta }}}_i(s-1)-{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\widehat{\mathit{\theta }}}_i(s-1)],i=1,2,\cdots ,m,\end{array}$$

(68)

$$\begin{array}{ccc}{\mathit{L}}_{\beta ,i}\left(s\right)& =& {\mathit{P}}_{\beta ,i}(s-1){\widehat{\mathit{\psi }}}_i\left(s\right){[1+{\widehat{\mathit{\psi }}}_i^{\mathrm{T}}\left(s\right){\mathit{P}}_{\beta ,i}(s-1){\widehat{\mathit{\psi }}}_i\left(s\right)]}^{-1},\end{array}$$

(69)

$$\begin{array}{ccc}{\mathit{P}}_{\beta ,i}\left(s\right)& =& {\mathit{P}}_{\beta ,i}(s-1)-{\mathit{L}}_{\beta ,i}\left(s\right){\left[{\mathit{P}}_{\beta ,i}(s-1){\widehat{\mathit{\psi }}}_i\left(s\right)\right]}^{\mathrm{T}},\end{array}$$

(70)

$$\begin{array}{ccc}{\widehat{\mathit{\theta }}}_i\left(s\right)& =& {\widehat{\mathit{\theta }}}_i(s-1)+{\mathit{L}}_{\theta ,i}\left(s\right)[y_i\left(s\right)-{\widehat{\mathit{\psi }}}_i^{\mathrm{T}}\left(s\right){\widehat{\mathit{\beta }}}_i(s-1)-{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\widehat{\mathit{\theta }}}_i(s-1)],\end{array}$$

(71)

$$\begin{array}{ccc}{\mathit{L}}_{\theta ,i}\left(s\right)& =& {\mathit{P}}_{\theta ,i}(s-1)\mathit{\phi }\left(s\right){[1+{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\mathit{P}}_{\theta ,i}(s-1)\mathit{\phi }\left(s\right)]}^{-1},\end{array}$$

(72)

$$\begin{array}{ccc}{\mathit{P}}_{\theta ,i}\left(s\right)& =& {\mathit{P}}_{\theta ,i}(s-1)-{\mathit{L}}_{\theta ,i}\left(s\right){\left[{\mathit{P}}_{\theta ,i}(s-1)\mathit{\phi }\left(s\right)\right]}^{\mathrm{T}},\end{array}$$

(73)

$$\begin{array}{ccc}\mathit{\phi }\left(s\right)& =& {[{\mathit{u}}^{\mathrm{T}}(s-1),{\mathit{u}}^{\mathrm{T}}(s-2),\cdots ,{\mathit{u}}^{\mathrm{T}}(s-n)]}^{\mathrm{T}},\end{array}$$

(74)

$$\begin{array}{ccc}\widehat{\mathit{\psi }}\left(s\right)& =& [{\widehat{\mathit{\psi }}}_s\left(s\right),{\widehat{\mathit{\psi }}}_n\left(s\right)],\end{array}$$

(75)

$$\begin{array}{ccc}& =& {[{\widehat{\mathit{\psi }}}_1\left(s\right),{\widehat{\mathit{\psi }}}_2\left(s\right),\cdots ,{\widehat{\mathit{\psi }}}_m\left(s\right)]}^{\mathrm{T}},\end{array}$$

(76)

$$\begin{array}{ccc}{\widehat{\mathit{\psi }}}_s\left(s\right)& =& [-\widehat{\mathit{x}}(s-1),-\widehat{\mathit{x}}(s-2),\cdots ,-\widehat{\mathit{x}}(s-n)],\end{array}$$

(77)

$$\begin{array}{ccc}{\widehat{\mathit{\psi }}}_n\left(s\right)& =& [-\widehat{\mathit{w}}(s-1),-\widehat{\mathit{w}}(s-2),\cdots ,-\widehat{\mathit{w}}(s-n_c),\widehat{\mathit{v}}(s-1),\widehat{\mathit{v}}(s-2),\cdots ,\widehat{\mathit{v}}(s-n_d)],\end{array}$$

(78)

$$\begin{array}{ccc}\widehat{\mathit{x}}\left(s\right)& =& {\widehat{\mathit{\psi }}}_s\left(s\right)\widehat{\mathit{\alpha }}\left(s\right)+{\widehat{\mathit{\theta }}}^{\mathrm{T}}\left(s\right)\mathit{\phi }\left(s\right),\end{array}$$

(79)

$$\begin{array}{ccc}\widehat{\mathit{w}}\left(s\right)& =& \mathit{y}\left(s\right)-{\widehat{\mathit{\psi }}}_s\left(s\right)\widehat{\mathit{\alpha }}\left(s\right)-{\widehat{\mathit{\theta }}}^{\mathrm{T}}\left(s\right)\mathit{\phi }\left(s\right),\end{array}$$

(80)

$$\begin{array}{ccc}{\widehat{v}}_i\left(s\right)& =& y_i\left(s\right)-{\widehat{\mathit{\psi }}}_i^{\mathrm{T}}\left(s\right){\widehat{\mathit{\beta }}}_i\left(s\right)-{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\widehat{\mathit{\theta }}}_i\left(s\right),\end{array}$$

(81)

$$\begin{array}{ccc}\widehat{\mathit{v}}\left(s\right)& =& {[{\widehat{v}}_1\left(s\right),{\widehat{v}}_2\left(s\right),\cdots ,{\widehat{v}}_m\left(s\right)]}^{\mathrm{T}},\end{array}$$

(82)

$$\begin{array}{ccc}\widehat{\mathit{\theta }}\left(s\right)& =& [{\widehat{\mathit{\theta }}}_1\left(s\right),{\widehat{\mathit{\theta }}}_2\left(s\right),\cdots ,{\widehat{\mathit{\theta }}}_m\left(s\right)].\end{array}$$

(83)

From the S-RGELS algorithm in (68)–(83), we can acquire m estimates ${\widehat{\mathit{\beta }}}_i\left(s\right)$ of $\mathit{\beta }$ from (68)–(83), it leads to a lot of redundant computation. However, we only need one parameter estimate of $\mathit{\beta }$. In order to cut down the redundant parameter estimates and improve the parameter estimation accuracy, the first way is to take their average value as the estimate of $\mathit{\beta }$:

$$\widehat{\mathit{\beta }}\left(s\right)=\frac{{\widehat{\mathit{\beta }}}_1\left(s\right)+{\widehat{\mathit{\beta }}}_2\left(s\right)+\cdots +{\widehat{\mathit{\beta }}}_m\left(s\right)}{m}\in {\mathbb{R}}^{n_0}.$$

(84)

Replacing ${\widehat{\mathit{\beta }}}_i(s-1)$ in (68)–(83) with $\widehat{\mathit{\beta }}(s-1)$ gives the partially-coupled subsystem recursive generalized extended least squares (PC-S-RGELS) algorithm:

$$\begin{array}{ccc}{\widehat{\mathit{\beta }}}_i\left(s\right)& =& \widehat{\mathit{\beta }}(s-1)+{\mathit{L}}_{\beta ,i}\left(s\right)[y_i\left(s\right)-{\widehat{\mathit{\psi }}}_i^{\mathrm{T}}\left(s\right)\widehat{\mathit{\beta }}(s-1)-{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\widehat{\mathit{\theta }}}_i(s-1)],i=1,2,\cdots ,m,\end{array}$$

(85)

$$\begin{array}{ccc}{\mathit{L}}_{\beta ,i}\left(s\right)& =& {\mathit{P}}_{\beta ,i}(s-1){\widehat{\mathit{\psi }}}_i\left(s\right){[1+{\widehat{\mathit{\psi }}}_i^{\mathrm{T}}\left(s\right){\mathit{P}}_{\beta ,i}(s-1){\widehat{\mathit{\psi }}}_i\left(s\right)]}^{-1},\end{array}$$

(86)

$$\begin{array}{ccc}{\mathit{P}}_{\beta ,i}\left(s\right)& =& {\mathit{P}}_{\beta ,i}(s-1)-{\mathit{L}}_{\beta ,i}\left(s\right){\left[{\mathit{P}}_{\beta ,i}(s-1){\widehat{\mathit{\psi }}}_i\left(s\right)\right]}^{\mathrm{T}},\end{array}$$

(87)

$$\begin{array}{ccc}{\widehat{\mathit{\theta }}}_i\left(s\right)& =& {\widehat{\mathit{\theta }}}_i(s-1)+{\mathit{L}}_{\theta ,i}\left(s\right)[y_i\left(s\right)-{\widehat{\mathit{\psi }}}_i^{\mathrm{T}}\left(s\right)\widehat{\mathit{\beta }}(s-1)-{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\widehat{\mathit{\theta }}}_i(s-1)],\end{array}$$

(88)

$$\begin{array}{ccc}{\mathit{L}}_{\theta ,i}\left(s\right)& =& {\mathit{P}}_{\theta ,i}(s-1)\mathit{\phi }\left(s\right){[1+{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\mathit{P}}_{\theta ,i}(s-1)\mathit{\phi }\left(s\right)]}^{-1},\end{array}$$

(89)

$$\begin{array}{ccc}{\mathit{P}}_{\theta ,i}\left(s\right)& =& {\mathit{P}}_{\theta ,i}(s-1)-{\mathit{L}}_{\theta ,i}\left(s\right){\left[{\mathit{P}}_{\theta ,i}(s-1)\mathit{\phi }\left(s\right)\right]}^{\mathrm{T}},\end{array}$$

(90)

$$\begin{array}{ccc}\mathit{\phi }\left(s\right)& =& {[{\mathit{u}}^{\mathrm{T}}(s-1),{\mathit{u}}^{\mathrm{T}}(s-2),\cdots ,{\mathit{u}}^{\mathrm{T}}(s-n)]}^{\mathrm{T}},\end{array}$$

(91)

$$\begin{array}{ccc}\widehat{\mathit{\psi }}\left(s\right)& =& [{\widehat{\mathit{\psi }}}_s\left(s\right),{\widehat{\mathit{\psi }}}_n\left(s\right)],\end{array}$$

(92)

$$\begin{array}{ccc}& =& {[{\widehat{\mathit{\psi }}}_1\left(s\right),{\widehat{\mathit{\psi }}}_2\left(s\right),\cdots ,{\widehat{\mathit{\psi }}}_m\left(s\right)]}^{\mathrm{T}},\end{array}$$

(93)

$$\begin{array}{ccc}{\widehat{\mathit{\psi }}}_s\left(s\right)& =& [-\widehat{\mathit{x}}(s-1),-\widehat{\mathit{x}}(s-2),\cdots ,-\widehat{\mathit{x}}(s-n)],\end{array}$$

(94)

$$\begin{array}{ccc}{\widehat{\mathit{\psi }}}_n\left(s\right)& =& [-\widehat{\mathit{w}}(s-1),-\widehat{\mathit{w}}(s-2),\cdots ,-\widehat{\mathit{w}}(s-n_c),\widehat{\mathit{v}}(s-1),\widehat{\mathit{v}}(s-2),\cdots ,\widehat{\mathit{v}}(s-n_d)],\end{array}$$

(95)

$$\begin{array}{ccc}\widehat{\mathit{x}}\left(s\right)& =& {\widehat{\mathit{\psi }}}_s\left(s\right)\widehat{\mathit{\alpha }}\left(s\right)+{\widehat{\mathit{\theta }}}^{\mathrm{T}}\left(s\right)\mathit{\phi }\left(s\right),\end{array}$$

(96)

$$\begin{array}{ccc}\widehat{\mathit{w}}\left(s\right)& =& \mathit{y}\left(s\right)-{\widehat{\mathit{\psi }}}_s\left(s\right)\widehat{\mathit{\alpha }}\left(s\right)-{\widehat{\mathit{\theta }}}^{\mathrm{T}}\left(s\right)\mathit{\phi }\left(s\right),\end{array}$$

(97)

$$\begin{array}{ccc}\widehat{\mathit{v}}\left(s\right)& =& \mathit{y}\left(s\right)-\widehat{\mathit{\psi }}\left(s\right)\widehat{\mathit{\beta }}\left(s\right)-{\widehat{\mathit{\theta }}}^{\mathrm{T}}\left(s\right)\mathit{\phi }\left(s\right),\end{array}$$

(98)

$$\begin{array}{ccc}\widehat{\mathit{\beta }}\left(s\right)& =& \frac{{\widehat{\mathit{\beta }}}_1\left(s\right)+{\widehat{\mathit{\beta }}}_2\left(s\right)+\cdots +{\widehat{\mathit{\beta }}}_m\left(s\right)}{m},\end{array}$$

(99)

$$\begin{array}{ccc}\widehat{\mathit{\theta }}\left(s\right)& =& [{\widehat{\mathit{\theta }}}_1\left(s\right),{\widehat{\mathit{\theta }}}_2\left(s\right),\cdots ,{\widehat{\mathit{\theta }}}_m\left(s\right)].\end{array}$$

(100)

Generally, for the recursive algorithms, it is considered that the parameter estimates are close to their true parameters as the time s increasing. Therefore, the estimate ${\widehat{\mathit{\beta }}}_{i-1}\left(s\right)$ is closer to the true parameter than the estimate ${\widehat{\mathit{\beta }}}_i(s-1)$. Combining (91)–(100) and replacing ${\widehat{\mathit{\beta }}}_i(s-1)$ in (68) and (71) when $i=1$ with ${\widehat{\mathit{\beta }}}_m(s-1)$ and ${\widehat{\mathit{\beta }}}_i(s-1)$ in (68) and (71) when $i=2,3,\cdots ,m$ with ${\widehat{\mathit{\beta }}}_{i-1}\left(s\right)$ give the partially-coupled recursive generalized extended least squares (PC-RGELS) algorithm:

$$\begin{array}{ccc}{\widehat{\mathit{\beta }}}_1\left(s\right)& =& {\widehat{\mathit{\beta }}}_m(s-1)+{\mathit{L}}_{\beta ,1}\left(s\right)[y_1\left(s\right)-{\widehat{\mathit{\psi }}}_1^{\mathrm{T}}\left(s\right){\widehat{\mathit{\beta }}}_m(s-1)-{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\widehat{\mathit{\theta }}}_1(s-1)],\end{array}$$

(101)

$$\begin{array}{ccc}{\mathit{L}}_{\beta ,1}\left(s\right)& =& {\mathit{P}}_{\beta ,1}(s-1){\widehat{\mathit{\psi }}}_i\left(s\right){[1+{\widehat{\mathit{\psi }}}_1^{\mathrm{T}}\left(s\right){\mathit{P}}_{\beta ,1}(s-1){\widehat{\mathit{\psi }}}_1\left(s\right)]}^{-1},\end{array}$$

(102)

$$\begin{array}{ccc}{\mathit{P}}_{\beta ,1}\left(s\right)& =& {\mathit{P}}_{\beta ,1}(s-1)-{\mathit{L}}_{\beta ,1}\left(s\right){\left[{\mathit{P}}_{\beta ,1}(s-1){\widehat{\mathit{\psi }}}_1\left(s\right)\right]}^{\mathrm{T}},\end{array}$$

(103)

$$\begin{array}{ccc}{\widehat{\mathit{\theta }}}_1\left(s\right)& =& {\widehat{\mathit{\theta }}}_1(s-1)+{\mathit{L}}_{\theta ,1}\left(s\right)[y_1\left(s\right)-{\widehat{\mathit{\psi }}}_1^{\mathrm{T}}\left(s\right){\widehat{\mathit{\beta }}}_m(s-1)-{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\widehat{\mathit{\theta }}}_1(s-1)],\end{array}$$

(104)

$$\begin{array}{ccc}{\mathit{L}}_{\theta ,1}\left(s\right)& =& {\mathit{P}}_{\theta ,1}(s-1)\mathit{\phi }\left(s\right){[1+{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\mathit{P}}_{\theta ,1}(s-1)\mathit{\phi }\left(s\right)]}^{-1},\end{array}$$

(105)

$$\begin{array}{ccc}{\mathit{P}}_{\theta ,1}\left(s\right)& =& {\mathit{P}}_{\theta ,1}(s-1)-{\mathit{L}}_{\theta ,i}\left(s\right){\left[{\mathit{P}}_{\theta ,1}(s-1)\mathit{\phi }\left(s\right)\right]}^{\mathrm{T}},\end{array}$$

(106)

$$\begin{array}{ccc}{\widehat{\mathit{\beta }}}_i\left(s\right)& =& {\widehat{\mathit{\beta }}}_{i-1}\left(s\right)+{\mathit{L}}_{\beta ,i}\left(s\right)[y_i\left(s\right)-{\widehat{\mathit{\psi }}}_i^{\mathrm{T}}\left(s\right){\widehat{\mathit{\beta }}}_{i-1}\left(s\right)-{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\widehat{\mathit{\theta }}}_i(s-1)],i=2,3,\cdots ,m,\end{array}$$

(107)

$$\begin{array}{ccc}{\mathit{L}}_{\beta ,i}\left(s\right)& =& {\mathit{P}}_{\beta ,i}(s-1){\widehat{\mathit{\psi }}}_i\left(s\right){[1+{\widehat{\mathit{\psi }}}_i^{\mathrm{T}}\left(s\right){\mathit{P}}_{\beta ,i}(s-1){\widehat{\mathit{\psi }}}_i\left(s\right)]}^{-1},\end{array}$$

(108)

$$\begin{array}{ccc}{\mathit{P}}_{\beta ,i}\left(s\right)& =& {\mathit{P}}_{\beta ,i}(s-1)-{\mathit{L}}_{\beta ,i}\left(s\right){\left[{\mathit{P}}_{\beta ,i}(s-1){\widehat{\mathit{\psi }}}_i\left(s\right)\right]}^{\mathrm{T}},\end{array}$$

(109)

$$\begin{array}{ccc}{\widehat{\mathit{\theta }}}_i\left(s\right)& =& {\widehat{\mathit{\theta }}}_i(s-1)+{\mathit{L}}_{\theta ,i}\left(s\right)[y_i\left(s\right)-{\widehat{\mathit{\psi }}}_i^{\mathrm{T}}\left(s\right){\widehat{\mathit{\beta }}}_{i-1}\left(s\right)-{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\widehat{\mathit{\theta }}}_i(s-1)],\end{array}$$

(110)

$$\begin{array}{ccc}{\mathit{L}}_{\theta ,i}\left(s\right)& =& {\mathit{P}}_{\theta ,i}(s-1)\mathit{\phi }\left(s\right){[1+{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\mathit{P}}_{\theta ,i}(s-1)\mathit{\phi }\left(s\right)]}^{-1},\end{array}$$

(111)

$$\begin{array}{ccc}{\mathit{P}}_{\theta ,i}\left(s\right)& =& {\mathit{P}}_{\theta ,i}(s-1)-{\mathit{L}}_{\theta ,i}\left(s\right){\left[{\mathit{P}}_{\theta ,i}(s-1)\mathit{\phi }\left(s\right)\right]}^{\mathrm{T}},\end{array}$$

(112)

$$\begin{array}{ccc}\mathit{\phi }\left(s\right)& =& {[{\mathit{u}}^{\mathrm{T}}(s-1),{\mathit{u}}^{\mathrm{T}}(s-2),\cdots ,{\mathit{u}}^{\mathrm{T}}(s-n)]}^{\mathrm{T}},\end{array}$$

(113)

$$\begin{array}{ccc}\widehat{\mathit{\psi }}\left(s\right)& =& [{\widehat{\mathit{\psi }}}_s\left(s\right),{\widehat{\mathit{\psi }}}_n\left(s\right)],\end{array}$$

(114)

$$\begin{array}{ccc}& =& {[{\widehat{\mathit{\psi }}}_1\left(s\right),{\widehat{\mathit{\psi }}}_2\left(s\right),\cdots ,{\widehat{\mathit{\psi }}}_m\left(s\right)]}^{\mathrm{T}},\end{array}$$

(115)

$$\begin{array}{ccc}{\widehat{\mathit{\psi }}}_s\left(s\right)& =& [-\widehat{\mathit{x}}(s-1),-\widehat{\mathit{x}}(s-2),\cdots ,-\widehat{\mathit{x}}(s-n)],\end{array}$$

(116)

$$\begin{array}{ccc}{\widehat{\mathit{\psi }}}_n\left(s\right)& =& [-\widehat{\mathit{w}}(s-1),-\widehat{\mathit{w}}(s-2),\cdots ,-\widehat{\mathit{w}}(s-n_c),\widehat{\mathit{v}}(s-1),\widehat{\mathit{v}}(s-2),\cdots ,\widehat{\mathit{v}}(s-n_d)],\end{array}$$

(117)

$$\begin{array}{ccc}\widehat{\mathit{x}}\left(s\right)& =& {\widehat{\mathit{\psi }}}_s\left(s\right)\widehat{\mathit{\alpha }}\left(s\right)+{\widehat{\mathit{\theta }}}^{\mathrm{T}}\left(s\right)\mathit{\phi }\left(s\right),\end{array}$$

(118)

$$\begin{array}{ccc}\widehat{\mathit{w}}\left(s\right)& =& \mathit{y}\left(s\right)-{\widehat{\mathit{\psi }}}_s\left(s\right)\widehat{\mathit{\alpha }}\left(s\right)-{\widehat{\mathit{\theta }}}^{\mathrm{T}}\left(s\right)\mathit{\phi }\left(s\right),\end{array}$$

(119)

$$\begin{array}{ccc}{\widehat{v}}_i\left(s\right)& =& y_i\left(s\right)-{\widehat{\mathit{\psi }}}_i^{\mathrm{T}}\left(s\right){\widehat{\mathit{\beta }}}_m\left(s\right)-{\mathit{\phi }}^{\mathrm{T}}\left(s\right){\widehat{\mathit{\theta }}}_i\left(s\right),\end{array}$$

(120)

$$\begin{array}{ccc}\widehat{\mathit{v}}\left(s\right)& =& {[{\widehat{v}}_1\left(s\right),{\widehat{v}}_2\left(s\right),\cdots ,{\widehat{v}}_m\left(s\right)]}^{\mathrm{T}},\end{array}$$

(121)

$$\begin{array}{ccc}\widehat{\mathit{\theta }}\left(s\right)& =& [{\widehat{\mathit{\theta }}}_1\left(s\right),{\widehat{\mathit{\theta }}}_2\left(s\right),\cdots ,{\widehat{\mathit{\theta }}}_m\left(s\right)].\end{array}$$

(122)

The procedures for achieving the PC-RGELS algorithm in (101)–(122) are as follows.

• For $s⩽0$, all variables are set to zero. Set the data length L. Let $s=1$, set the initial values ${\widehat{\mathit{\beta }}}_m\left(0\right)={\mathbf{1}}_{n_0}/p_0$, ${\widehat{\mathit{\theta }}}_i\left(0\right)={\mathbf{1}}_{rn}/p_0$, $\widehat{\mathit{x}}\left(0\right)={\mathbf{1}}_m/p_0$, $\widehat{\mathit{w}}\left(0\right)={\mathbf{1}}_m/p_0$, $\widehat{\mathit{v}}\left(0\right)={\mathbf{1}}_m/p_0$, ${\mathit{P}}_{\beta ,i}\left(0\right)=p_0{\mathit{I}}_{n_0}$, ${\mathit{P}}_{\theta ,i}\left(0\right)=p_0{\mathit{I}}_{rn}$, $p_0={10}^6$.
• Collect the input-output data $\mathit{u}\left(s\right)={[u_1\left(s\right),u_2\left(s\right),\cdots ,u_r\left(s\right)]}^{\mathrm{T}}$ and $\mathit{y}\left(s\right)={[y_1\left(s\right),y_2\left(s\right),\cdots ,y_m\left(s\right)]}^{\mathrm{T}}$, and construct $\mathit{\phi }\left(s\right)$ using (113).
• Form ${\widehat{\mathit{\psi }}}_s\left(s\right)$ and ${\widehat{\mathit{\psi }}}_n\left(s\right)$ using (116)–(117) and construct $\widehat{\mathit{\psi }}\left(s\right)$ using (114), and read ${\widehat{\mathit{\psi }}}_i\left(s\right)$ from $\widehat{\mathit{\psi }}\left(s\right)$ in (115), $i=1,2,3,\cdots ,m$.
• Compute ${\mathit{P}}_{\beta ,1}\left(s\right)$ and ${\mathit{P}}_{\theta ,1}\left(s\right)$ using (103) and (106), and compute ${\mathit{L}}_{\mathit{\beta },1}\left(s\right)$ and ${\mathit{L}}_{\theta ,1}\left(s\right)$ using (102) and (105), and update the estimates ${\widehat{\mathit{\beta }}}_1\left(s\right)$ and ${\widehat{\mathit{\theta }}}_1\left(s\right)$ using (101) and (104).
• For $i=2,3,\cdots ,m$, calculate ${\mathit{P}}_{\beta ,i}\left(s\right)$ and ${\mathit{P}}_{\theta ,i}\left(s\right)$ using (109) and (112), and compute ${\mathit{L}}_{\mathit{\beta },i}\left(s\right)$ and ${\mathit{L}}_{\theta ,i}\left(s\right)$ using (108) and (111), and refresh the estimates ${\widehat{\mathit{\beta }}}_i\left(s\right)$ and ${\widehat{\mathit{\theta }}}_i\left(s\right)$ using (107) and (1110).
• Construct $\widehat{\mathit{\theta }}\left(s\right)$ by (112), calculate the estimates $\widehat{\mathit{x}}\left(s\right)$, $\widehat{\mathit{w}}\left(s\right)$ and $\widehat{\mathit{v}}\left(s\right)$ using (118)–(120).
• Compare s with L: if $s⩽L$, increase s by 1 and go to Step 2; otherwise obtain the estimate $\widehat{\mathit{\beta }}\left(L\right)$ and $\widehat{\mathit{\theta }}\left(L\right)$ and terminate this procedure.

The flowchart of computing $\widehat{\mathit{\beta }}\left(L\right)$ and $\widehat{\mathit{\theta }}\left(L\right)$ in the PC-RGELS algorithm is shown in Figure 2.

5. Example

The numerical simulation is based on the following M-OEARMA-like system:

$$\begin{array}{lll}\mathit{y}\left(s\right)& =& \frac{\mathit{Q}\left(z\right)}{\alpha \left(z\right)}\mathit{u}\left(s\right)+\frac{D\left(z\right)}{C\left(z\right)}\mathit{v}\left(s\right),\\ \alpha \left(z\right)& :=& 1+a_1{z}^{-1}=1+0.26z^{-1},\\ \mathit{Q}\left(z\right)& :=& {\mathit{A}}_1{z}^{-1}=\left[\begin{array}{cc}{a}_2& a_3\\ a_4& a_5\end{array}\right]z^{-1}=\left[\begin{array}{cc}-0.57z^{-1}& 0.93z^{-1}\\ 0.87z^{-1}& {0.65}^{-1}\end{array}\right],\\ C\left(z\right)& :=& 1+c_1{z}^{-1}=1-0.56z^{-1},\\ D\left(z\right)& :=& 1+d_1{z}^{-1}=1+0.32z^{-1},\end{array}$$

the parameter vectors of the system are

$$\begin{array}{lll}\mathit{\alpha }& :=& {[a_1,a_2,a_3,a_4,a_5]}^{\mathrm{T}}={[0.26,-0.57,0.93,0.87,0.65]}^{\mathrm{T}},\\ \mathit{\beta }& :=& {[c_1,d_1]}^{\mathrm{T}}={[-0.56,0.32]}^{\mathrm{T}},\\ \mathit{\vartheta }& :=& \left[\begin{array}{c}\mathit{\alpha }\\ \mathit{\beta }\end{array}\right]={[a_1,a_2,a_3,a_4,a_5,c_1,d_1]}^{\mathrm{T}},\\ & =& {[0.26,-0.57,0.93,0.87,0.65,-0.56,0.32]}^{\mathrm{T}}.\end{array}$$

Here, the inputs $\left\{u_1\left(s\right)\right\}$ and $\left\{u_2\left(s\right)\right\}$ are taken as two independent persistent excitation signal sequences with zero mean and unit variances, $\left\{v_1\left(s\right)\right\}$ and $\left\{v_2\left(s\right)\right\}$ are taken as two white noise sequences with zero mean and variances ${\sigma }_1^2$ for $v_1\left(s\right)$ and ${\sigma }_2^2$ for $v_2\left(s\right)$. Set ${\sigma }_1^2={\sigma }_2^2={0.30}^2$. Then we can use them to generate the output vector $\mathit{y}\left(s\right)={[y_1\left(s\right),y_2\left(s\right)]}^{\mathrm{T}}$. Set the data length $L=3000$. Applying the RGELS, PC-S-RGELS and PC-RGELS algorithms to identy the parameters of the given system. The parameter estimates and their estimation errors are shown in

Table 1. The RGELS estimates and their errors.

s	${\mathit{a}}_1$	${\mathit{a}}_2$	${\mathit{a}}_3$	${\mathit{a}}_4$	${\mathit{a}}_5$	${\mathit{c}}_1$	${\mathit{d}}_1$	$\mathit{\delta }(\%)$
100	0.21713	−0.55471	0.90627	0.84332	0.71040	−0.48095	0.18013	10.72413
200	0.27110	−0.56221	0.94490	0.86938	0.64018	−0.47391	0.23325	7.35692
500	0.24836	−0.54410	0.93994	0.86873	0.66404	−0.55063	0.06758	15.08382
1000	0.26144	−0.55626	0.94061	0.86873	0.64676	−0.59024	0.03508	16.99606
2000	0.25755	−0.54862	0.94853	0.86843	0.65790	−0.57538	0.03620	16.91802
3000	0.25925	−0.55582	0.94269	0.86686	0.66788	−0.56765	0.05916	15.52793
True values	0.26000	−0.57000	0.93000	0.87000	0.65000	−0.56000	0.32000

,

Table 2. The PC-S-RGELS estimates and their errors.

s	${\mathit{a}}_1$	${\mathit{a}}_2$	${\mathit{a}}_3$	${\mathit{a}}_4$	${\mathit{a}}_5$	${\mathit{c}}_1$	${\mathit{d}}_1$	$\mathit{\delta }(\%)$
100	0.20140	−0.56387	0.92631	0.70446	0.79768	−0.35762	0.68920	28.39191
200	0.24432	−0.56531	0.93373	0.79687	0.69082	−0.38330	0.64577	22.51734
500	0.25446	−0.56183	0.94070	0.83956	0.67391	−0.47448	0.47151	10.58749
1000	0.26117	−0.56373	0.93438	0.85347	0.65752	−0.53154	0.40254	5.29978
2000	0.25977	−0.55982	0.93912	0.86161	0.65870	−0.54046	0.37127	3.42342
3000	0.26036	−0.56350	0.93576	0.86339	0.66214	−0.54830	0.37027	3.20492
True values	0.26000	−0.57000	0.93000	0.87000	0.65000	−0.56000	0.32000

and

Table 3. The PC-RGELS estimates and their errors.

s	${\mathit{a}}_1$	${\mathit{a}}_2$	${\mathit{a}}_3$	${\mathit{a}}_4$	${\mathit{a}}_5$	${\mathit{c}}_1$	${\mathit{d}}_1$	$\mathit{\delta }(\%)$
100	0.24978	−0.60667	0.87031	0.78115	0.56724	−0.39507	0.30167	12.87537
200	0.27882	−0.58559	0.91259	0.82014	0.61386	−0.42910	0.41254	10.32160
500	0.25835	−0.56532	0.93012	0.85414	0.64496	−0.58010	0.25299	4.26796
1000	0.26982	−0.56619	0.93045	0.85993	0.64143	−0.60694	0.26914	4.21835
2000	0.26104	−0.56177	0.93780	0.86480	0.65036	−0.57656	0.30440	1.53744
3000	0.26119	−0.56487	0.93488	0.86510	0.65647	−0.56786	0.34146	1.49741
True values	0.26000	−0.57000	0.93000	0.87000	0.65000	−0.56000	0.32000

. The parameter estimation errors $\delta :=\parallel \widehat{\mathit{\theta }}\left(s\right)-\mathit{\theta }\parallel /\parallel \mathit{\theta }\parallel$ versus s is shown in Figure 3.

Taking the simulation conditions with ${\sigma }_1^2={\sigma }_2^2={0.20}^2$ and ${\sigma }_1^2={\sigma }_2^2={0.60}^2$, using the PC-RGELS algorithm to identify the parameters of the given systems, respectively, the parameter estimates and their estimation errors are shown in

Table 4. The PC-RGELS estimates and their errors with ${\sigma }^2={0.20}^2$.

s	${\mathit{a}}_1$	${\mathit{a}}_2$	${\mathit{a}}_3$	${\mathit{a}}_4$	${\mathit{a}}_5$	${\mathit{c}}_1$	${\mathit{d}}_1$	$\mathit{\delta }(\%)$
100	0.24988	−0.56806	0.90461	0.89323	0.60879	−0.59192	0.23442	6.30236
200	0.26570	−0.55866	0.90290	0.89140	0.62628	−0.56227	0.30809	2.68948
500	0.26807	−0.57330	0.92228	0.86927	0.62782	−0.55656	0.31763	1.50427
1000	0.26418	−0.57144	0.92234	0.87348	0.63193	−0.56993	0.32145	1.34724
2000	0.26205	−0.57095	0.92456	0.87546	0.64561	−0.54789	0.32853	1.03059
3000	0.26050	−0.56927	0.92491	0.87164	0.64979	−0.54752	0.31805	0.81408
True values	0.26000	−0.57000	0.93000	0.87000	0.65000	−0.56000	0.32000

and

Table 5. The PC-RGELS estimates and their errors with ${\sigma }^2={0.60}^2$.

s	${\mathit{a}}_1$	${\mathit{a}}_2$	${\mathit{a}}_3$	${\mathit{a}}_4$	${\mathit{a}}_5$	${\mathit{c}}_1$	${\mathit{d}}_1$	$\mathit{\delta }(\%)$
100	−0.07666	−1.10323	1.30348	0.68217	0.72605	−0.49016	−0.09463	51.44481
200	0.02363	−0.82959	1.02696	0.84406	0.70880	−0.49353	−0.03312	30.52969
500	0.16820	−0.71220	1.00389	0.83283	0.64901	−0.50435	0.06928	18.85790
1000	0.20791	−0.63623	0.95860	0.86602	0.63011	−0.54443	0.15480	11.21196
2000	0.23053	−0.60174	0.94024	0.87889	0.65490	−0.53307	0.25093	5.15448
3000	0.24053	−0.58686	0.93174	0.87014	0.66124	−0.53546	0.27576	3.42763
True values	0.26000	−0.57000	0.93000	0.87000	0.65000	−0.56000	0.32000

, and the estimation errors $\delta$ versus s is shown in Figure 4.

From Table 1, Table 2, Table 3, Table 4 and Table 5 and Figure 3 and Figure 4, we can draw the following conclusions.

• The parameter estimation errors given by the RGELS, PC-S-RGELS and PC-RGELS algorithms become smaller as s increasing. Thus the proposed algorithms for multivariable OEARMA-like system are effective.
• Under the same simulation conditions, the PC-S-RGELS and PC-RGELS algorithms can give more accurate parameter estimates compared with the RGELS algorithm.
• A lower noise level leads to a higher parameter estimation accuracy by the PC-RGELS algorithm under the same data length.

6. Conclusions

In this paper, we have dealt with the parameter identification problems of the M-OEARMA-like systems. A partially-coupled recursive generalized extended least squares algorithm is presented based on the hierarchical identification principle and the coupling identification concept. At the last of this paper, we discussed the comparison between the RGELS algorithm and the derived algorithms, and found that the PC-S-RGELS and PC-RGELS algorithms have higher computational efficiency than the RGELS algorithm on account of the proposed algorithms avoid calculating the inverse of the covariance matrix. The analysis and numerical example indicate that the PC-RGELS algorithm can give more accurate parameter estimates than the RGELS algorithm. Furthermore, the proposed methods can be extended to other fields by means of some other mathematical tools and approaches [80,81,82,83,84,85,86] to model industrial processes [87,88,89,90,91,92,93,94,95,96,97,98].