% 
% ch4ex2.m simulates the economy from
% mid term 1 then kalman filters for 
% the unobserved states
%

clear all
diary ch4ex2.mout
diary off
delete ch4ex2.mout
diary ch4ex2.mout


% ----- simulate the economy from midterm ----- %

delta = 0.025;
alpha = 0.30;
thetabar = 1.05;

% ----- solutions from midterm ----- %

f1 =[0.4681    0.4957
    0.3569   -0.1313]

f2 =[0.9000         0
    0.2459    0.8858 ]

C = [0.2;0]

% ----- compute states ----- % 

x(1:2,1)=[0;0];

for i = 1:1000
 x(:,i+1) = f2*x(:,i) + C*randn(1);
end

% ----- compute endo vars ----- %

for i = 1:1000
 y(:,i) = f1*x(:,i) ;
end

% ----- add measurement errors to observations ----- %
 
for i = 1:1000
 yobs(:,i) = y(:,i) + 0.5*randn(2,1);
end

% ----- put data into observed data matrix ---- %

Xobs = [yobs(1,100:1000); yobs(2,100:1000)] ;

% ----- keep unobserved x's ----- %

xunobs = [x(1,100:1000); x(2,100:1000)] ;

% ----- initialize filter ----- %

xtt(:,1) = [0;0];
Ptt(:,:) = dlyap1(f1,0.5^2*eye(2));

% ----- start kalman filter where s = t-1 ----- %

for i = 1:900

% - prediction - %

 xts(:,i)    = f2*xtt(:,i);
 Pts         = f2*Ptt*f2' + C*C' ;
 Xtsobs(:,i) = f1*xts(:,i) ;
 fts         = f1*Pts*f1' + 0.5^2*eye(2) ;

% - updating - %

 xtt(:,i+1) = xts(:,i) + Pts*f1'*inv( fts )*(Xobs(:,i)-Xtsobs(:,i)) ;
 Ptt        = Pts - Pts*f1'*pinv( fts )*f1*Pts ;

end

% ----- now plot the two seperate series vs the actual ----- %

figure(1)
 subplot(2,1,1)
 plot(1:300,xtt(1,2:301)',1:300,xunobs(1,1:300)',':')
 legend('filtered','actual')

 subplot(2,1,2)
 plot(1:300,xtt(2,2:301)',1:300,xunobs(2,1:300)',':')
 legend('filtered','actual')