## 3 The Bias Compensation-Based Parameter and State Estimation Algorithm
### The Parameter Estimation Algorithm
Equation (19) can be rewritten as:  
\[\lim_{k\to\infty}\frac{1}{k}[\mathbf{P}^{-1}(k)(\mathbf{\hat{\vartheta}}_{\text{LS}}( k)-\mathbf{\vartheta})]=-(\mathbf{R}\mathbf{\vartheta}-\mathbf{r}).\]  
or:  
\[\lim_{k\to\infty}\mathbf{\hat{\vartheta}}_{\text{LS}}(k)=\mathbf{\vartheta}-k\mathbf{P}( k)(\mathbf{R}\mathbf{\vartheta}-\mathbf{r}). \tag{25}\]  
It can be seen from Equation (25) that the bias of the least squares estimate \(\mathbf{\hat{\vartheta}}_{\text{LS}}(k)\) can be eliminated by adding a compensation term \(\Delta\mathbf{\vartheta}(k):=k\mathbf{P}(k)(\mathbf{R}\mathbf{\vartheta}-\mathbf{r})\). That is, the estimate \(\mathbf{\hat{\vartheta}}_{\text{C}}(k):=\mathbf{\hat{\vartheta}}_{\text{LS}}(k)+ \Delta\mathbf{\vartheta}(k)\) is an unbiased estimate of the true parameter vector \(\mathbf{\vartheta}\). Thus, the unbiased estimate \(\mathbf{\hat{\vartheta}}_{\text{C}}(k)\) can be computed by the following recursive expressions,  
\[\mathbf{\hat{\vartheta}}_{\text{C}}(k) =\mathbf{\hat{\vartheta}}_{\text{LS}}(k)+k\mathbf{P}(k)[\mathbf{\hat{R}}(k) \mathbf{\hat{\vartheta}}_{\text{C}}(k)  
Let \(J(k):=\sum_{i=1}^{k}\varepsilon_{\text{LS}}^{2}(i)\) be the cost function at time \(k\), \(\hat{\mathbf{Q}}(k)\) and \(\hat{\mathbf{\zeta}}(k)\) be the estimates of \(\mathbf{Q}\) and \(\mathbf{\zeta}\) at instant \(k\), respectively. Then, the estimate of the noise variance \(\delta\) can be computed by:  
\[\hat{\delta}(k) =\frac{\frac{1}{k}J(k)}{\hat{\mathbf{\sigma}}_{\text{C}}^{\text{T}}(k )\hat{\mathbf{Q}}(k)\hat{\mathbf{\sigma}}_{\text{LS}}(k)-\hat{\mathbf{\zeta}}^{\text{T}}(k )[\hat{\mathbf{\sigma}}_{\text{C}}(k)+\hat{\mathbf{\sigma}}_{\text{LS}}(k)]+1+\hat{\mathbf{ \sigma}}_{\text{v}}^{\text{T}}(k)\hat{\mathbf{\sigma}}_{\text{v}}(k)}, \tag{30}\] \[J(k) =J(k-1)+\frac{[y(k)-\mathbf{\phi}_{\text{s}}^{\text{T}}(k)\hat{\mathbf{ \sigma}}_{\text{LS}}(k-1)]^{2}}{1+\mathbf{\phi}_{\text{s}}^{\text{T}}(k)\mathbf{P}(k-1 )\mathbf{\phi}_{\text{s}}(k)},\] (31) \[\hat{\mathbf{\zeta}}(k) =[\hat{R}(n)/\hat{\delta}(k-1),\hat{R}(n-1)/\hat{\delta}(k-1), \cdots,\hat{R}(1)/\hat{\delta}(k-1),0,0\cdots,0]^{\text{T}},\] (32) \[\hat{\mathbf{Q}}(k) =\left[ \begin{array}{cccc}1+\hat{\mathbf{\sigma}}_{\text{v}}^{ \text{T}}(k)\hat{\mathbf{\sigma}}_{\text{v}}(k)&\hat{R}(1)/\hat{\delta}(k-1)&\cdots &\hat{R}(n-1)/\hat{\delta}(k-1)\\ \hat{R}(1)/\hat{\delta}(k-1)&1+\hat{\mathbf{\sigma}}_{\text{v}}^{\text{T}}(k)\hat {\mathbf{\sigma}}_{\text{v}}(k)&\cdots&\hat{R}(n-2)/\hat{\delta}(k-1)\\ \vdots&\vdots&\ddots&\vdots\\ \hat{R}(n-1)/\hat{\delta}(k-1)&\hat{R}(n-2)/\hat{\delta}(k-1)&\cdots&1+\hat{ \mathbf{\sigma}}_{\text{v}}^{\text{T}}(k)\hat{\mathbf{\sigma}}_{\text{v}}(k)\\ &\mathbf{0}\end{array} \right],\] (33) \[\frac{\hat{R}(i)}{\hat{\delta}(k-1)} =\hat{\mathbf{\sigma}}_{\text{v}}^{\text{T}}(i:n_{e})\left[ \begin{array}{ c}1\\ \hat{\mathbf{\sigma}}_{\text{v}}(1:n_{e}-i)\end{array} \right],\;i=1,2,\cdots,n_{e}. \tag{34}\]  
Note that Equations (30)-(34) involve the estimate \(\hat{\mathbf{\sigma}}_{\text{v}}(k)\) of the noise parameter vector \(\mathbf{\phi}_{\text{v}}\), which can be computed by the noise model in Equation (16). Let \(\hat{e}(k)\) and \(\hat{v}(k)\) be the estimates of \(e(k)\) and \(v(k)\), respectively.