diff --git a/matlab/kalman/likelihood/univariate_kalman_filter.m b/matlab/kalman/likelihood/univariate_kalman_filter.m index e2f280a66791fe0a0fadabdbdbdbe1357edeca6c..ccae3665037093853ae3b6b214c2a97c42425d77 100644 --- a/matlab/kalman/likelihood/univariate_kalman_filter.m +++ b/matlab/kalman/likelihood/univariate_kalman_filter.m @@ -132,27 +132,25 @@ while notsteady && t<=last for i=1:rows(z) if Zflag prediction_error = Y(d_index(i),t) - z(i,:)*a; - Fi = z(i,:)*P*z(i,:)' + H(d_index(i)); + PZ = P*z(i,:)'; + Fi = z(i,:)*PZ + H(d_index(i)); else prediction_error = Y(d_index(i),t) - a(z(i)); - Fi = P(z(i),z(i)) + H(d_index(i)); + PZ = P(:,z(i)); + Fi = PZ(z(i)) + H(d_index(i)); end if Fi>kalman_tol - if Zflag - Ki = (P*z(i,:)')/Fi; - else - Ki = P(:,z(i))/Fi; - end + Ki = PZ/Fi; if t>no_more_missing_observations K(:,i) = Ki; end a = a + Ki*prediction_error; - P = P - (Fi*Ki)*Ki'; + P = P - PZ*Ki'; lik(s) = lik(s) + log(Fi) + prediction_error*prediction_error/Fi + l2pi; end end a = T*a; - P = T*P*transpose(T) + QQ; + P = T*P*T' + QQ; if t>=no_more_missing_observations notsteady = max(abs(K(:)-oldK))>riccati_tol; oldK = K(:); diff --git a/matlab/kalman/likelihood/univariate_kalman_filter_ss.m b/matlab/kalman/likelihood/univariate_kalman_filter_ss.m index 78528286323a7bca1333ba3031cbd44f0b136245..24a0cc9b00ea20e4107b90230248e7c079cf66c9 100644 --- a/matlab/kalman/likelihood/univariate_kalman_filter_ss.m +++ b/matlab/kalman/likelihood/univariate_kalman_filter_ss.m @@ -90,19 +90,17 @@ while t<=last for i=1:pp if Zflag prediction_error = Y(i,t) - Z(i,:)*a; - Fi = Z(i,:)*PP*Z(i,:)' + H(i); + PPZ = PP*Z(i,:)'; + Fi = Z(i,:)*PPZ + H(i); else prediction_error = Y(i,t) - a(Z(i)); - Fi = PP(Z(i),Z(i)) + H(i); + PPZ = PP(:,Z(i)); + Fi = PPZ(Z(i)) + H(i); end if Fi>kalman_tol - if Zflag - Ki = (PP*Z(i,:))'/Fi; - else - Ki = PP(:,Z(i))/Fi; - end + Ki = PPZ/Fi; a = a + Ki*prediction_error; - PP = PP - (Fi*Ki)*transpose(Ki); + PP = PP - PPZ*Ki'; likk(s) = likk(s) + log(Fi) + prediction_error*prediction_error/Fi + l2pi; end end