Consider a Linear Equation,<\/p>\n\n\n\n
$$ y_i = \\sum_{j=1}^n C_{i,j} x_j +v_i,\\quad i=1,2,\u2026$$<\/p>\n\n\n\n
, where C_{i,j}<\/span> are scalars and v_i\\in \\mathbb{R}<\/span> is the measurement noise. The noise is unknown, while we assume it follows certain patterns (the assumptions are due to some statistical properties of the noise). We assume v_i, v_j<\/span> are independent for i\\neq j<\/span>. Properties are mean of zero, and variance equals sigma squared.<\/p>\n\n\n\n $$\\mathbb{E}(v_i)=0$$<\/p>\n\n\n\n $$\\mathbb{E}(v_i^2) = \\sigma_i^2$$<\/p>\n\n\n\n We can rewrite y_i = \\sum_{j=1}^n C_{i,j} x_j +v_i<\/span> as,<\/p>\n\n\n\n $$ \\begin{pmatrix} y_1 \\ y_2 \\ \\vdots\\ y_s\\end{pmatrix} = \\begin{pmatrix} C_{11} & C_{12} & \\cdots & C_{1n} \\ C_{21} & C_{22}& \\cdots & C_{2n} \\ \\vdots & \\vdots & \\cdots & \\vdots \\ C_{s1} & C_{s2} & \\cdots & C_{sn}\\end{pmatrix} \\begin{pmatrix} x_1 \\ x_2 \\ \\vdots\\ x_n\\end{pmatrix} + \\begin{pmatrix} v_1 \\ v_2 \\ \\vdots\\ v_s\\end{pmatrix} $$<\/p>\n\n\n\n , in a matrix form,<\/p>\n\n\n\n $$ \\vec{y} = C \\vec{x} + \\vec{v} $$<\/p>\n\n\n\n , but I would write in a short form,<\/p>\n\n\n\n $$ y= C x +v$$<\/p>\n\n\n\n We solve for the least squared estimator from the optimisation problem, (there is a squared L2 norm)<\/p>\n\n\n\n $$ \\min_x || y-Cx ||_2^2 $$<\/p>\n\n\n\n The classic least squared estimator might not work well when data evolving. So, there emerges a Recursive Least Squared Method to deal with the discrete-time instance. Let’s say, for a discrete-time instance k<\/span>, y_k \\in \\mathbb{R}’<\/span> is within a set of measurements group follows,<\/p>\n\n\n\n $$y_k = C_k x + v_k$$<\/p>\n\n\n\n , where C_k \\in \\mathbb{R}^{l\\times n}<\/span>, and v_k \\in \\mathbb{R}^l<\/span> is the measurement noise vector. We assume that the covariance of the measurement noise is given by,<\/p>\n\n\n\n $$ \\mathbb{E}[v_k v_k^T] = R_k$$<\/p>\n\n\n\n , and<\/p>\n\n\n\n $$\\mathbb{E}[v_k]=0$$<\/p>\n\n\n\n The recursive least squared method has the following form in this section,<\/p>\n\n\n\n $$\\hat{x}k = \\hat{x}{k-1} + K_k (y_k – C_k \\hat{x}_{k-1})$$<\/p>\n\n\n\n , where \\hat{x}k<\/span> and \\hat{x}{k-1}<\/span> are the estimates of the vector x<\/span> at the discrete-time instants k<\/span> and k-1<\/span>, and K_k \\in \\mathbb{R}^{n\\times l}<\/span> is the gain matrix that we need to determine. K_k<\/span> is coined the ‘Gain Matrix’<\/p>\n\n\n\n The above equation updates the estimate of x<\/span> at the time instant k<\/span> on the basis of the estimate \\hat{x}_{k-1}<\/span> at the previous time instant k-1<\/span> and on the basis of the measurement y_k<\/span> obtained at the time instant k<\/span>, as well as on the basis of the gain matrix K_k<\/span> computed at the time instant k<\/span>.<\/p>\n\n\n\n $\\hat{x}$ is the estimate.<\/p>\n\n\n\n $$ \\hat{x}k = \\begin{pmatrix} \\hat{x}{1,k} \\ \\hat{x}{2,k} \\ \\vdots \\\\hat{x}{n,k} \\end{pmatrix} $$<\/p>\n\n\n\n , which is corresponding with the true vector x<\/span>.<\/p>\n\n\n\n $$x = \\begin{pmatrix} x_1 \\ x_2 \\ \\vdots \\ x_n \\end{pmatrix}$$<\/p>\n\n\n\n The estimation error, \\epsilon_{i,k} = x_i – \\hat{x}_{i,k} \\quad i=1,2,\u2026,n<\/span>.<\/p>\n\n\n\n $$\\epsilon_k = \\begin{pmatrix} \\epsilon_{1,k} \\ \\epsilon_{2,k} \\ \\vdots \\\\epsilon_{n,k} \\end{pmatrix} = x – \\hat{x}_k = \\begin{pmatrix} x_1-\\hat{x}_{1,k} \\ x_2 – \\hat{x}_{2,k} \\ \\vdots \\x_n-\\hat{x}_{n,k} \\end{pmatrix} $$<\/p>\n\n\n\n The gain K_k<\/span> is computed by minimising the sum of variances of the estimation errors,<\/p>\n\n\n\n $$ W_k = \\mathbb{E}(\\epsilon_{1,k}^2) + \\mathbb{E}(\\epsilon_{2,k}^2) + \\cdots + \\mathbb{E}(\\epsilon_{n,k}^2) $$<\/p>\n\n\n\n Next, let’s show the cost function could be represented as follows, (tr(.)<\/span> is the trace of a matrix)<\/p>\n\n\n\n $$ W_k = tr(P_k) $$<\/p>\n\n\n\n , and P_k<\/span> is the estimation error covariance matrix defined by<\/p>\n\n\n\n $$ P_k = \\mathbb{E}(\\epsilon_k \\epsilon_k^T )$$<\/p>\n\n\n\n Or, says,<\/p>\n\n\n\n $$ K_k = arg\\min_{K_k} W_k = tr\\bigg( \\mathbb{E}(\\epsilon_k \\epsilon_k^T ) \\bigg)$$<\/p>\n\n\n\n Why is that?<\/p>\n\n\n\n $$\\epsilon_k \\epsilon_k^T = \\begin{pmatrix} \\epsilon_{1,k} \\ \\epsilon_{2,k} \\\\vdots \\ \\epsilon_{n,k} \\end{pmatrix} \\begin{pmatrix} \\epsilon_{1,k} & \\epsilon_{2,k} & \\cdots & \\epsilon_{n,k} \\end{pmatrix}$$<\/p>\n\n\n\n $$ = \\begin{pmatrix} \\epsilon_{1,k}^2 & \\cdots & \\epsilon_{1,k}\\epsilon_{n,k} \\ \\vdots & \\epsilon_{i,k}^2 & \\vdots \\ \\epsilon_{1,k}\\epsilon_{n,k} & \\cdots & \\epsilon_{n,k}^2\\end{pmatrix} $$<\/p>\n\n\n\n So,<\/p>\n\n\n\n $$ P_k = \\mathbb{E}[\\epsilon_k \\epsilon_k^T] $$<\/p>\n\n\n\n $$tr(P_k) = \\mathbb{E}(\\epsilon_{1,k}^2) + \\mathbb{E}(\\epsilon_{2,k}^2) + \\cdots + \\mathbb{E}(\\epsilon_{n,k}^2)$$<\/p>\n\n\n\n $$ K_k = arg\\min_{K_k} W_k = tr\\bigg( \\mathbb{E}(\\epsilon_k \\epsilon_k^T ) \\bigg) = tr(P_k)$$<\/p>\n\n\n\n Let’s derive the optimisation problem.<\/p>\n\n\n\n $$\\epsilon_k = x-\\hat{x}_k$$<\/p>\n\n\n\n $$ =x-\\hat{x}{k-1} – K_k(y_k – C_k \\hat{x}{k-1}) $$<\/p>\n\n\n\n $$ = x- \\hat{x}{k-1} – K_k (C_k x + v_k – C_k \\hat{x}{k-1}) $$<\/p>\n\n\n\n $$ = (I – K_k C_k)(x-\\hat{x}_{k-1}) – K_k v_k $$<\/p>\n\n\n\n $$ =(I-K_k C_k )\\epsilon_{k-1} – K_k v_k $$<\/p>\n\n\n\n Recall y_k = C_k x + v_k<\/span> and \\hat{x}k = \\hat{x}{k-1} + K_k (y_k – C_k \\hat{x}_{k-1})<\/span><\/p>\n\n\n\n So, \\epsilon_k \\epsilon_k^T<\/span> would be,<\/p>\n\n\n\n $$\\epsilon_k \\epsilon_k^T = \\bigg((I-K_k C_k )\\epsilon_{k-1} – K_k v_k\\bigg)\\bigg((I-K_k C_k )\\epsilon_{k-1} – K_k v_k\\bigg)^T$$<\/p>\n\n\n\n $P_k = \\mathbb{E}(\\epsilon_k \\epsilon_k^T)$, and $P_{k-1} = \\mathbb{E}(\\epsilon_{k-1} \\epsilon_{k-1}^T)$.<\/p>\n\n\n\n $\\mathbb{E}(\\epsilon_{k-1} v_k^T) = \\mathbb{E}(\\epsilon_{k-1}) \\mathbb{E}(v_k^T) =0$ by the white noise property of $\\epsilon$ and $v$. However, $\\mathbb{E}(v_k v_k^T) = R_k$. Substituting all those into $P_k$, we would get,<\/p>\n\n\n\n $$P_k = (I – K_k C_k)P_{k-1}(I – K_k C_k)^T + K_k R_k K_k^T$$<\/p>\n\n\n\n $$ P_k = P_{k-1} – P_{k-1} C_k^T K_k^T – K_k C_k P_{k-1} + K_k C_k P_{k-1}C_k^T K_k^T + K_k R_k K_k^T $$<\/p>\n\n\n\n $$W = tr(P_k)= tr(P_{k-1}) – tr(P_{k-1} C_k^T K_k^T) – tr(K_k C_k P_{k-1}) + tr(K_k C_k P_{k-1}C_k^T K_k^T) + tr(K_k R_k K_k^T) $$<\/p>\n\n\n\n We take F.O.C. to solve for K_k = arg\\min_{K_k} W_k = tr\\bigg( \\mathbb{E}(\\epsilon_k \\epsilon_k^T ) \\bigg) = tr(P_k)<\/span>, by letting \\frac{\\partial W_k}{\\partial K_k} = 0<\/span>. See the Matrix Cookbook<\/a> and find how to do derivatives w.r.t. K_k<\/span>.<\/p>\n\n\n\n $$\\frac{\\partial W_k}{\\partial K_k} = -2P_{k-1} C_k^T + 2K_k C_k P_{k-1} C_k^T + 2K_k R_k = 0$$<\/p>\n\n\n\n We solve for K_k<\/span>,<\/p>\n\n\n\n $$ K_k = P_{k-1} C_k^T (R_k + C_k P_{k-1} C_k^T)^{-1}$$<\/p>\n\n\n\n , we let L_k = R_k + C_k P_{k-1} C_k^T<\/span>, and L_k<\/span> has the following property L_k = L_k^T<\/span> and L_k^{-1} = (L_k^{-1})^T<\/span><\/p>\n\n\n\n $$ K_k = P_{k-1} C_k^T L_k^{-1} $$<\/p>\n\n\n\n Plug K_k = P_{k-1} C_k^T K_k^{-1}<\/span> back into P_k<\/span>.<\/p>\n\n\n\n $$ P_k = P_{k-1} – K_kC_k P_{k-1} = (I-K_k C_k)P_{k-1} $$<\/p>\n\n\n\n In the end, the Recursive Least Squared Method could be summarised as the following three equations.<\/p>\n\n\n\n $$ K_k = P_{k-1} C_k^T (R_k + C_k P_{k-1} C_k^T)^{-1}$$<\/p>\n\n\n\n $$\\hat{x}_k = \\hat{x}_{k-1} + K_k (y_k – C_k \\hat{x}_{k-1})$$<\/p>\n\n\n\n (I-K_k C_k)P_{k-1}<\/span><\/p>\n\n\n\n
\n\n\n\nRecursive Least Squared Method<\/h3>\n\n\n\n
Notation<\/h4>\n\n\n\n
\n\n\n\nOptimisation<\/h4>\n\n\n\n
\n\n\n\nSummary<\/h4>\n\n\n\n
Reference<\/h4>\n\n\n\n