# Coupled Least Squares Identification Algorithms for Multivariate Output-Error Systems
## The Multivariate Auxiliary Model Coupled Identification Algorithm
### The Coupled Auxiliary Model Based Recursive Least Squares Algorithm  
In order to avoid the redundant parameter estimates, we use the coupling identification concept to derive a coupled AM-RLS algorithm based on the C-S-AM-RLS algorithm. Referring to the partially coupled SG identification method [24], and with the help of the Jacobi or Gauss-Seidel iterative algorithm, replacing \(\hat{\boldsymbol{\vartheta}}_{m}(t-1)\) with \(\hat{\boldsymbol{\vartheta}}_{1}(t-1)\) for \(i=1\), replacing \(\hat{\boldsymbol{\vartheta}}_{i}(t-1)\) with \(\hat{\boldsymbol{\vartheta}}_{i-1}(t)\) for \(i=2,3,\cdots,m\), we can obtain the following coupled auxiliary model based recursive least squares (C-AM-RLS) identification algorithm:  
\[\hat{\boldsymbol{\vartheta}}_{1}(t) =\hat{\boldsymbol{\vartheta}}_{m}(t-1)+\boldsymbol{L}_{1}(t)[y_{1 }(t)-\hat{\boldsymbol{\vartheta}}_{1}^{\top}(t)\hat{\boldsymbol{\vartheta}}_{ m}(t-1)], \tag{26}\] \[\boldsymbol{L}_{1}(t) =\boldsymbol{P}_{m}(t-1)\hat{\boldsymbol{\vartheta}}_{1}(t)[1+ \hat{\boldsymbol{\varphi}}_{1}^{\top}(t)\boldsymbol{P}_{m}(t-1)\hat{ \boldsymbol{\vartheta}}_{1}(t)]^{-1},\] (27) \[\boldsymbol{P}_{1}(t) =[\boldsymbol{I}_{n_{0}}-\boldsymbol{L}_{1}(t)\hat{\boldsymbol{ \vartheta}}_{1}^{\top}(t)]\boldsymbol{P}_{m}(t-1),\] (28) \[\hat{\boldsymbol{\vartheta}}_{i}(t) =\hat{\boldsymbol{\vartheta}}_{i-1}(t)+\boldsymbol{L}_{i}(t)[y_{ i}(t)-\hat{\boldsymbol{\vartheta}}_{i}^{\top}(t)\hat{\boldsymbol{\vartheta}}_{i-1}(t)],\] (29) \[\boldsymbol{L}_{i}(t) =\boldsymbol{P}_{i-1}(t)\hat{\boldsymbol{\vartheta}}_{i}(t)[1+ \hat{\boldsymbol{\varphi}}_{i}^{\top}(t)\boldsymbol{P}_{i-1}(t)\hat{ \boldsymbol{\vartheta}}_{i}(t)]^{-1},\] (30) \[\boldsymbol{P}_{i}(t) =[\boldsymbol{I}_{n_{0}}-\boldsymbol{L}_{i}(t)\hat{\boldsymbol{ \vartheta}}_{i}^{\top}(t)]\boldsymbol{P}_{i-1}(t),\quad i=2,3,\cdots,m,\] (31) \[\hat{\boldsymbol{\Phi}}_{x}(t) =[-\mathbf{x}_{\text{a}}(t-1),-\mathbf{x}_{\text{a}}(t-2),\cdots,- \mathbf{x}_{\text{a}}(t-n_{a})],\] (32) \[\hat{\boldsymbol{\Phi}}(t) =[\hat{\boldsymbol{\Phi}}_{x}(t),\boldsymbol{\Phi}_{s}(t)]\] (33) \[=[\hat{\boldsymbol{\varphi}}_{1}(t),\hat{\boldsymbol{\varphi}}_{ 2}(t),\cdots,\hat{\boldsymbol{\varphi}}_{m}(t)]^{\top},\] (34) \[\mathbf{x}_{\text{a}}(t) =\hat{\boldsymbol{\Phi}}(t)\hat{\boldsymbol{\vartheta}}(t). \tag{35}\]In the above algorithm in (26)-(35), \(\hat{\mathbf{\vartheta}}_{i}(t)\in\mathbb{R}^{n_{0}}\) is the parameter estimation vector of the \(i\)th subsystem at time \(t\), \(\mathbf{L}_{i}(t)\in\mathbb{R}^{n_{0}}\) is the gain vector of the \(i\)th subsystem at time \(t\), \(\mathbf{P}_{i}(t)\in\mathbb{R}^{n_{0}\times n_{0}}\) is the covariance matrix of the \(i\)th subsystem at time \(t\). \(\hat{\mathbf{\vartheta}}_{i-1}(t)\) and \(\mathbf{P}_{i-1}(t)\) are the parameter estimation vector and the covariance matrix of the \((i-1)\)th subsystem at time \(t\), respectively; \(\hat{\mathbf{\vartheta}}_{m}(t-1)\) and \(\mathbf{P}_{m}(t-1)\) are the parameter estimation vector and the covariance matrix of the \(m\)th subsystem at time \(t-1\), respectively, the system parameter estimation vector is defined by the parameter estimation vector of the \(m\)th subsystem at time \(t\): \(\hat{\mathbf{\vartheta}}(t)=\hat{\mathbf{\vartheta}}_{m}(t)\). The procedure of computing the parameter estimation vector \(\hat{\mathbf{\vartheta}}_{m}(t)\) in (26)-(35) is as follows. 1. Set the initial values: \(t=1\), \(\hat{\mathbf{\vartheta}}_{m}(0)=\mathbf{1}_{n_{0}}/p_{0}\), \(\mathbf{P}_{m}(0)=p_{0}\mathbf{I}_{n_{0}}\), \(\mathbf{x}_{\text{a}}(t-j)=\mathbf{1}_{m}/p_{0}\), \(j=1\), \(2\), \(\cdots\), \(n_{a}\), \(p_{0}=10^{6}\). Set the data length \(L\). 2. Collect the observation data \(y(t)\) and \(\mathbf{\Phi}_{s}(t)\), and construct \(\hat{\mathbf{\Phi}}_{x}(t)\) and \(\hat{\mathbf{\Phi}}(t)\) by (32) and (33). 3. Read \(\hat{\mathbf{\vartheta}}_{i}(t)\) from \(\hat{\mathbf{\Phi}}(t)\) in (34), compute \(\mathbf{L}_{1}(t)\) and \(\mathbf{P}_{1}(t)\) by (27) and (28), and update the parameter estimation vector \(\hat{\mathbf{\vartheta}}_{1}(t)\) by (26). 4.