%Macro 202a Homework #2 Solution Anton A Cheremukhin 19 October 2005 close all; clear all; clc; %data [time, gdp, cons, inv1, hours,inv2] = textread('datarbc.txt','%d %f %f %f %f %f'); T=size(cons,1); gdphp=hpfilter(gdp,100); y0=log(gdp./gdphp); conshp=hpfilter(cons,100); c0=log(cons./conshp); hourshp=hpfilter(hours,100); l0=log(hours./hourshp); %c0=log(cons); X=[ones(T,1) time]; beta=inv(X'*X)*(X'*c0); c0=c0-X*beta; %y0=log(gdp); X=[ones(T,1) time]; beta=inv(X'*X)*(X'*y0); y0=y0-X*beta; %l0=log(hours); X=[ones(T,1) time]; beta=inv(X'*X)*(X'*l0); l0=l0-X*beta; [capital] = textread('capital.txt','%f'); capitalhp=hpfilter(capital,100); k0=log(capital./capitalhp); %k0=log(capital); X=[ones(T,1) time]; beta=inv(X'*X)*(X'*k0); k0=k0-X*beta; %predefined parameters beta=0.96; g=0.0018; ro=3; alpha=1/3; fi=0.17; lambda=0.65; delta=0.1; sigmay=0.0004; sigmak=0.003; sigmac=0.0001; sigmal=0.001; vol=0.0002; %set of parameters q0=[log(-1+1/beta), log(-1+1/g), log(ro), log(-1+1/alpha), log(-1+1/fi) ,log(-1+1/lambda), log(-1+1/delta), log(sigmay), log(sigmak), log(sigmac), log(sigmal), log(vol)]; y2=[y0,k0,c0,l0]; H0=0.0001*eye(12); crit=0.00001; nit=1000; % find likelihood maximizer with csminwel [fh,xh,gh,H,itct,fcount,retcodeh] = csminwel(@neolike, q0, H0, [], crit, nit, y2); beta=1/(1+exp(xh(1))); g=1/(1+exp(xh(2))); ro=exp(xh(3)); alpha=1/(1+exp(xh(4))); fi=1/(1+exp(xh(5))); lambda=1/(1+exp(xh(6))); delta=1/(1+exp(xh(7))); sigmay=exp(xh(8)); sigmak=exp(xh(9))+0.00000000001; sigmac=exp(xh(10)); sigmal=exp(xh(11)); vol=exp(xh(12)); %%%%%%%%%%%%%%%%%%% %new parameters - definitions teta=((((1+g)^(1-fi*(1-ro)))/beta)-(1-delta))/alpha;%%=cap/gdp nu=1-(delta+g)/teta; %=cons/gdp labor=1-1/(1+fi*(1-alpha)/(nu*(1-fi))); %steady state values cap=labor/(teta^(1/(1-alpha))); gdp=teta*cap; cons=nu*gdp; %more parameter definitions phi=alpha*beta*gdp/((1+g)*cap); omega=fi*(1-ro)-1; a1=labor/(1-labor); pi=(1-fi)*(1-ro)*a1; %matrices D=[ 0,-phi,-omega,pi,0,phi,omega,-pi; -gdp,(1+g)*cap,cons,0,0,0,0,0; 1,0,0,-(1-alpha),-1,0,0,0; 1,0,-1,-(1+a1),0,0,0,0; 0,0,0,0,1,0,0,0; 1,0,0,0,0,0,0,0; 0,1,0,0,0,0,0,0; 0,0,1,0,0,0,0,0;]; E=[ 0,0,0,0,0,0,0,0; 0,(1-delta)*cap,0,0,0,0,0,0; 0,alpha,0,0,0,0,0,0; 0,0,0,0,0,0,0,0; 0,0,0,0,lambda,0,0,0; 0,0,0,0,0,1,0,0; 0,0,0,0,0,0,1,0; 0,0,0,0,0,0,0,1;]; %use gensys [G,C,impact,fmat,fwt,ywt,gev,eu]=gensys(D,E,zeros(8,1),[0;0;0;0;1;0;0;0],[0,0,0; 0,0,0; 0,0,0; 0,0,0; 0,0,0; 1,0,0; 0,1,0; 0,0,1]); F = G([1,2,3,4,5],[1,2,3,4,5]); %H = [1 0 0 0 0; 0 0 1 0 0; 0 0 0 1 0]; H = [eye(4),zeros(4,1)];%G([1,2,3,4],[2,5]); Q = vol*impact([1,2,3,4,5])*impact([1,2,3,4,5])'; %R = diag([sigmay,sigmac,sigmal]); R = diag([sigmay,sigmak,sigmac,sigmal]); initx = [y0(1);k0(1);c0(1);l0(1);0]; %initx = [y0(1);0;c0(1);l0(1);0]; initV = 0.001*eye(5); %y1=[y0,c0,l0]'; y1=y2'; [xfilt, Vfilt, VVfilt, loglik] = kalman_filter(y1, F, H, Q, R, initx, initV); y1=y2'; figure(1) clf for i=1:4 subplot(2,2,i) hold on plot(time, y1(i,:), 'b*-'); errorbar(time, xfilt(i,:), sqrt(reshape(Vfilt(i, i, :),1,size(Vfilt,3))), 'rx:'); hold off legend('observed', 'filtered', 2) xlabel('time') if i==1 ylabel('output') elseif i==2 ylabel('capital') elseif i==3 ylabel('consumption') elseif i==4 ylabel('labor') end end %[beta, g, ro, alpha, fi ,lambda, delta, sigmay, sigmak, sigmac, sigmal] figure(2) clf hold on plot(time, y1(1,:), 'b*-'); errorbar(time, xfilt(5,:), sqrt(reshape(Vfilt(5, 5, :),1,size(Vfilt,3))), 'rx:'); hold off legend('output', 'productivity filtered', 2) xlabel('time') ylabel('productivity vs output') [beta, g, ro, alpha, fi ,lambda, delta, sigmay, sigmak, sigmac, sigmal, vol]