3
$\begingroup$

I am trying to implement an Extended Kalman Filter (EKF) and it is becoming harder than I thought.

I have one question. I noticed that the covariance matrix which should get updated over each iteration is not symmetric. I am debugging through MATLAB. I know that P should be symmetric and stay symmetric.

What does it mean if the covariance matrix is not symmetric? Where it could be the error?

EDIT: Any covariance matrix MUST be symmetric no matter what. The symmetry comes from its definition. The covariance tells you how two variables are related and therefore if $x$ is related to $y$ also $y$ will be related to $x$ and of course in the same way.

My problem was that I was checking with the MATLAB function issymmetric(A) if the matrix $A$ was symmetric. Apparently that function checks for exact symmetry and therefore if the involved matrixes are computed numerically and therefore it is not exactly symmetric it will give you false/0 as a result. But, if for the same matrix I checked A-A' I had something of the order of 1e-15.

Thanks in advance.

$\endgroup$
  • $\begingroup$ Unitarily diagonalize $P=V\Lambda V^{*}$, where the adjoint $V^{*}$ is just $V^T$ if your matrix is real. Having done this, do not store $V^{*}$ separately. $\endgroup$ – thb Jun 25 '17 at 17:39
  • $\begingroup$ My comment above does not give enough information, does it? I'll make it a proper answer. $\endgroup$ – thb Jun 25 '17 at 17:42
  • $\begingroup$ Where are you trying to implement an EKF? Microcontroller, PLC ? $\endgroup$ – Daniel Mårtensson Jul 7 '17 at 22:41
0
$\begingroup$

Unitarily diagonalize $P=V\Lambda V^{*}$, where the adjoint $V^{*}$ is just $V^T$ if your matrix is real. Having done this, do not store $V^{*}$ separately.

The matrix $\Lambda$ is to have the eigenvalues of $P$ along its main diagonal and is to be null elsewhere. The columns of $V$ are to be mutually orthogonal unit vectors, whereby $V^{*}V = I$.

The reason this works is that it lends $V$ the useful property that $V^{-1}=V^{*}.$ Therefore, because $\Lambda$ is also trivial to invert, you never ask your software to invoke the Gauss-Jordan algorithm that is desymmetrifying your matrix. Don't let your software decide for itself how to invert anything, but rather work around the problem as here explained by telling your software what the inverses of the various elements are.

EXPLANATION

A matrix represents a linear operation and, if the matrix is invertible, then it has certain properties. Moreover, if the matrix is symmetric, it has other properties. One cannot use a symmetric matrix effectively without exploiting the properties that belong to all symmetric matrices. If you try to treat a symmetric matrix like any old matrix, then, numerically, weird things are likely to happen.

When I said symmetric, I really meant self-adjoint, which is a jargon word but here is an example: $$\left[\begin{array}{cc}2 & 4+i3 \\ 4-i3 & 7\end{array}\right].$$ So, as you see, a self-adjoint matrix is symmetric with complex conjugation. If the matrix is real, then "symmetric" and "self-adjoint" mean the same.

Another word for "self-adjoint" is Hermitian.

The symbol for the transpose is $P^T$ as you know. The symbol for the adjoint is $P^{*}=\text{conjugate}\{P^{T}\}$, but for real matrices $P^{*}=P^T$, so you can read the two symbols as though they were the same.

Now, most matrices are diagonalizable. Significantly, all self-adjoint (real symmetric) matrices are diagonalizable. This means that you can factor the self-adjoint matrix in a special form that resembles $$P = \left[\begin{array}{cc}0.8 & 0.6 \\ -0.6 & 0.8\end{array}\right] \left[\begin{array}{cc}1 & 0 \\ 0 & 2\end{array}\right] \left[\begin{array}{cc}0.8 & -0.6 \\ 0.6 & 0.8\end{array}\right]. $$ For convenience, we give names to the factors, abbreviating the last line as $$ P = V\Lambda V^{*}, $$ the inverse of which conveniently is $$ P^{-1} = V\Lambda^{-1} V^{*}. $$ Writing out the factors, $$ P^{-1} = \left[\begin{array}{cc}0.8 & 0.6 \\ -0.6 & 0.8\end{array}\right] \left[\begin{array}{cc}1 & 0 \\ 0 & 0.5\end{array}\right] \left[\begin{array}{cc}0.8 & -0.6 \\ 0.6 & 0.8\end{array}\right]. $$ So, as matrices go, that's pretty easy.

Interestingly and usefully, $$ P^n = V\Lambda^n V^{*}, $$ wherein, if you think about it, $$ \Lambda^n = \left[\begin{array}{cc}1^n & 0 \\ 0 & 2^n\end{array}\right], $$ so the computer really cannot very well mess that up.

My advice to you is this: don't let the computer directly handle $P$. Instead, make the computer deal in $\Lambda$ and $V$.

One last question is, how do you diagonalize? The answer is that your software almost certainly has a built-in function to diagonalize, but if you wish to learn how to do it, yourself, look up eigenvalue, eigenvector and/or diagonalization in your linear-algebra textbook.

Good luck.

$\endgroup$
  • $\begingroup$ For further information, see sect. 14.11 here. $\endgroup$ – thb Jun 25 '17 at 17:51
  • $\begingroup$ I am sorry but I really have problems understanding what you wrote and therefore what you mean. But thanks for the answer anyway $\endgroup$ – kalmanIsAGameChanger Jun 25 '17 at 20:26
  • 1
    $\begingroup$ Yes, that's always a problem, pitching the answer at the wrong level. I'm an electrical engineer like you (as I suspect), and the majority of the answers on this site fly over my head, but I try. Let me back up, edit the answer, and try to pitch it at a comprehensible level. $\endgroup$ – thb Jun 25 '17 at 20:31
  • $\begingroup$ Thanks a lot. I appreciate it. Take your time $\endgroup$ – kalmanIsAGameChanger Jun 25 '17 at 20:38
  • $\begingroup$ I have updated the answer. As far as I can see, MATLAB's function 'eig()' performs the needed diagonalization. Anyway, either the answer helps or it doesn't. If not, I wish you well finding what you need elsewhere; but my answer tells you what I would do. $\endgroup$ – thb Jun 25 '17 at 21:07
0
$\begingroup$

A quick ad-hoc fix that (in my experience) works great is to simply "symmetrize" the $P$ matrix every time you calculate a new potentially asymmetric value for it by doing:

$P'=\dfrac{P+P^T}{2}$

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.