## 3 Identification Algorithm
Minimizing the cost function:  
\[\varepsilon(i)=\min\lVert\boldsymbol{r}_{k-1}-\boldsymbol{\phi}_{i}\boldsymbol {\theta}_{i}\rVert^{2},\ i=1,2,\cdots,N\]0  
leads to the \(k\)-th parameter estimate:  
\[\varepsilon(i)=\min\lVert\boldsymbol{r}_{k-1}-\boldsymbol{\phi}_{i}\boldsymbol {\theta}_{i}\rVert^{2},\ i=1,2,\cdots,N\]1  
Then, the \(k\)-th residual can be computed by:  
\[\varepsilon(i)=\min\lVert\boldsymbol{r}_{k-1}-\boldsymbol{\phi}_{i}\boldsymbol {\theta}_{i}\rVert^{2},\ i=1,2,\cdots,N\]2  
If the system is noise free, then the residual will be zero after \(K\) iterations, and the non-zero parameters, as well as their positions can be exactly determined. This is the so-called orthogonal matching pursuit (OMP) algorithm. However, for systems disturbed by noises, the OMP algorithm cannot give the sparsest solution, because the residual \(\mathbf{r}_{K}\) is non-zero after \(K\) iterations, and the iteration continues. In order to solve this problem, an appropriate small threshold \(\epsilon\) can be set to filter the parameter estimate \(\hat{\mathbf{\theta}}_{\Lambda_{k}}\). If \(|\hat{\Theta}_{\Lambda_{\epsilon}}|<\epsilon\), let \(\hat{\Theta}_{\Lambda_{\epsilon}}=0\), where \(\hat{\Theta}_{\Lambda_{\epsilon}}\) is the \(i\)-th element of \(\hat{\mathbf{\theta}}_{\Lambda_{k}}\). The choice of \(\epsilon\) will be discussed in the simulation part. Denote \(\hat{\mathbf{\theta}}_{\Lambda_{k\epsilon}}\) as the parameter estimate after the filtering. Then, the residual can be computed by:  
\[\varepsilon(i)=\min\lVert\boldsymbol{r}_{k-1}-\boldsymbol{\phi}_{i}\boldsymbol {\theta}_{i}\rVert^{2},\ i=1,2,\cdots,N\]3if \(\|\mathbf{r}_{k}\|<\varepsilon\), then stop the iteration; otherwise, go to apply another iteration. If the iteration stops at \(k\), we can recover the parameter estimate \(\hat{\boldsymbol{\theta}}\in\mathbb{R}^{N}\) from \(\hat{\boldsymbol{\theta}}_{\Lambda_{kc}}\in\mathbb{R}^{k}\) based on the support set \(\Lambda_{k}\). It is worth noting that taking the derivative of \(J(\boldsymbol{\theta}_{\Lambda_{k}})\) with respect to \(\hat{\boldsymbol{\theta}}_{\Lambda_{k}}\) gives:  
\[\varepsilon(i)=\min\lVert\boldsymbol{r}_{k-1}-\boldsymbol{\phi}_{i}\boldsymbol {\theta}_{i}\rVert^{2},\ i=1,2,\cdots,N\]4 \[\varepsilon(i)=\min\lVert\boldsymbol{r}_{k-1}-\boldsymbol{\phi}_{i}\boldsymbol {\theta}_{i}\rVert^{2},\ i=1,2,\cdots,N\]5  
The above equation shows that the residual \(\mathbf{r}_{k}\) is orthogonal to the sub-information matrix \(\boldsymbol{\Phi}_{\Lambda_{k}}\). Therefore, it does not require computing the inner products of \(\mathbf{r}_{k}\) and the columns in \(\boldsymbol{\Phi}_{\Lambda_{k}}\) during the \((k+1)\)-th iteration. Modify Equation (7) as:  
\[\varepsilon(i)=\min\lVert\boldsymbol{r}_{k-1}-\boldsymbol{\phi}_{i}\boldsymbol {\theta}_{i}\rVert^{2},\ i=1,2,\cdots,N\]6  
Then, the computational burden can be reduced. From the recovered parameter estimation \(\hat{\boldsymbol{\theta}}\), the orders and time delays of each input channel can be estimated. The orders \(n_{bi}\) can be determined by the number of elements of each non-zero block. Then, the time delay estimates can be computed from the orders \(n_{bi}\) and the input maximum regression length \(l\). It can be seen from Equation (2) that \(\boldsymbol{\theta}\) contains \((r+1)\) zero blocks. Denote the number of zeros of each zero-block as \(n_{i}\), \(i=1,2,\cdots,r+1\). Then, the time delay estimates can be computed by:  
\[ \begin{split}\hat{d}_{1}&=n_{1},\\ \hat{d}_{i}&=n_{i}-(l-\hat{d}_{i-1}-n_{b_{i}}),\ i= 2,3,\cdots,r\end{split}  \tag{15}\]  
The proposed TH-OMP algorithm can be implemented as the following steps. 1. Set the initial values: let \(k=0\), \(\boldsymbol{r}_{0}=\boldsymbol{y}_{m}\), \(\Lambda_{0}=\varnothing\) and \(\hat{\boldsymbol{\theta}}_{\Lambda_{0}}=\boldsymbol{0}\). Give the error tolerant \(\varepsilon\). 2. Increment \(k\) by one; compute the solution support \(\lambda_{k}\) using Equation (14). 3.