[R] kalman filter estimation

[Paul Gilbert]

adschai at optonline.net wrote:
> Hi,
> 
> Following convention below:
> y(t) = Ax(t)+Bu(t)+eps(t) # observation eq
> x(t) = Cx(t-1)+Du(t)+eta(t) # state eq
> 
> I modified the following routine (which I copied from: http://www.stat.pitt.edu/stoffer/tsa2/Rcode/Kall.R) to accommodate u(t), an exogenous input to the system. 
> 
> for (i in 2:N){
>  xp[[i]]=C%*%xf[[i-1]]
>  Pp[[i]]=C%*%Pf[[i-1]]%*%t(C)+Q
>    siginv=A[[i]]%*%Pp[[i]]%*%t(A[[i]])+R
>  sig[[i]]=(t(siginv)+siginv)/2     # make sure sig is symmetric
>    siginv=solve(sig[[i]])          # now siginv is sig[[i]]^{-1}
>  K=Pp[[i]]%*%t(A[[i]])%*%siginv
>  innov[[i]]=as.matrix(yobs[i,])-A[[i]]%*%xp[[i]]
>  xf[[i]]=xp[[i]]+K%*%innov[[i]]
>  Pf[[i]]=Pp[[i]]-K%*%A[[i]]%*%Pp[[i]]
>  like= like + log(det(sig[[i]])) + t(innov[[i]])%*%siginv%*%innov[[i]]
>  }
>    like=0.5*like
>    list(xp=xp,Pp=Pp,xf=xf,Pf=Pf,like=like,innov=innov,sig=sig,Kn=K)
> }
> 
> I tried to fit my problem and observe that I got positive log likelihood mainly because the log of determinant of my variance matrix is largely negative. That's not good because they should be positive. Have anyone experience this kind of instability?
> 
The determinant should not be negative (the line "# make sure sig is 
symmetric" should look after that) but the log can be negative.

> Also, I realize that I have about 800 sample points. The above routine when being plugged to optim becomes very slow. Could anyone share a faster way to compute kalman filter? 

The KF recursion you show is not time varying, but the code you show 
looks like it may allow for time varying models. (I'm just guessing, I'm 
not familiar with the code.) If you do not need the time varying aspect 
then this is probably a slow way to implement. Package dse1 in the dse 
bundle implements the recursion the way you show, with an exogenous 
input u(t), so you don't even have to modify the code. It is implemented 
in R for demonstration, but also in fortran for speed. See the dse Users 
Guide for more details, and also ?SS and ?l.SS. One caveat is that the 
way the likelihood is cumulated over i in your code above allows for 
time varying sig, which in theory can be important for a diffuse prior. 
In dse the likelihood is calculated using the residual to get a sample 
estimate of sig, which is faster. I have not found cases where this 
makes much difference (but would be interested to know of one).

> And my last problem is, optim with my defined feasible space does not converge. I have about 20 variables that I need to identify using MLE method. Is there any other way that I can try out? I tried most of the methods available in optim already. They do not converge at all...... Thank you.

This used to be a well known problem for multivariate state space 
models, but seems to be largely forgotten. It does not seriously affect 
univariate models. For a review of the old literature see section 5 of 
Gilbert 1993, available at 
<http://www.bank-banque-canada.ca/pgilbert/research.php#MiscResearch>. 
The bottom line is that you are in trouble trying to use hill climbing 
methods and state-space models unless you have a very good starting 
point, or are luck. This is a feature of the parameter space, not a 
reflection on the optimization routine. One solution is to start by 
estimating a VAR or ARMA model and then convert it to a state space 
model, which was one of the original purposes of dse. (See ?bft for 
example.)

Paul Gilbert

> - adschai
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
====================================================================================

La version française suit le texte anglais.

------------------------------------------------------------------------------------

This email may contain privileged and/or confidential information, and the Bank of
Canada does not waive any related rights. Any distribution, use, or copying of this
email or the information it contains by other than the intended recipient is
unauthorized. If you received this email in error please delete it immediately from
your system and notify the sender promptly by email that you have done so. 

------------------------------------------------------------------------------------

Le présent courriel peut contenir de l'information privilégiée ou confidentielle.
La Banque du Canada ne renonce pas aux droits qui s'y rapportent. Toute diffusion,
utilisation ou copie de ce courriel ou des renseignements qu'il contient par une
personne autre que le ou les destinataires désignés est interdite. Si vous recevez
ce courriel par erreur, veuillez le supprimer immédiatement et envoyer sans délai à
l'expéditeur un message électronique pour l'aviser que vous avez éliminé de votre
ordinateur toute copie du courriel reçu.