[EMAIL PROTECTED] 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@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. ______________________________________________ R-help@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.