% USA_Model2_komplet.m % ************* % PART 3: M U L T I - E Q U A T I O N S M O D E L S % 13.3 Simulation Example % ============================================================== % Robert S. Pindyck & Daniel L. Rubinfeld % Econometric Models & Economic Forecast, 4/e % McGraw-Hill, 1998 % ISBN 0-07-050208-0 % % Data to accompany Econometric Models & Economic Forecast, 4/e % (c) Irwin/McGraw-Hill, 1998 % Disk 1 of 1 % ISBN 0-07-848043-4 % ============================================================== clear all; close all; clc; %data ulozena v textovem souboru USAmacro.dat load USAmacro.dat; maxl = 4; %maximalni zpozdeni v modelu last = 25; %pocet poslednich vynechanych hodnot - odhad do 3/1973 - 1. ropna krize (14*4+2 ctvrtleti) t = USAmacro(maxl+1:end-last,1); % Puvoden t = ctvrtleti od 2Q1950-1Q1988 n = length(t); %endogeni promenne C = USAmacro(maxl+1:end-last,2); % Realna agregatni osobni potreba (v cenach roku 1982) [mld.] I = USAmacro(maxl+1:end-last,3); % Realne hrube domaci investice (v cenach roku 1982) [mld.] R = USAmacro(maxl+1:end-last,4); % Urokova sazba 3 mesicnich Treasury bill [% p.a.] Y = USAmacro(maxl+1:end-last,5); % Realny HDP (ocisten o exporty a importy; v cenach 1982) [mld.] %zpozdene endogenni promenne C_1 = USAmacro(maxl+1-1:end-last-1,2); %Zpozdene C(t-1) Y_1 = USAmacro(maxl+1-1:end-last-1,5); %... dY_1 = USAmacro(maxl+1-1:end-last-1,5)-USAmacro(maxl+1-2:end-last-2,5); %Y(t-1)-Y(t-2) R_4 = USAmacro(maxl+1-4:end-last-4,4); sumR_1 = USAmacro(maxl+1-1:end-last-1,4)+USAmacro(maxl+1-2:end-last-2,4); %R(t-1)+R(t-2) %exogenni promenne konst = ones(n,1); G = USAmacro(maxl+1:end-last,6); % Realne vladni vydaje (v cenach roku 1982) [mld.] dM = USAmacro(maxl+1:end-last,7)-USAmacro(maxl+1-1:end-last-1,7); % Zmena realne penezni zasoby (M1) [mld.] %dalsi endogenni promenne dY = Y-Y_1; %Y(t)-Y(t-1) %pomocne exogenni promenne M = USAmacro(maxl+1:end-last,7); %penezni zasoba - agregat M1 %graficke zobrazeni endogennich promennych figure set(gcf,'Name',' Endogenous variables'); subplot(2,2,1) hndl2=plot(t,C); set(hndl2,'LineWidth',2); title('Real Aggregate Personal Consumption'); ylabel('C [ billions of 1982 dolars ]'); xlabel('Time [ quarter ]'); grid subplot(2,2,2) hndl3=plot(t,I); set(hndl3,'LineWidth',2); title('Real Gross Domestic Investment'); ylabel('I [ billions of 1982 dolars ]'); xlabel('Time [ quarter ]'); grid subplot(2,2,3) hndl4=plot(t,R,'r'); set(hndl4,'LineWidth',2); title('Interest Rate on 3-Month Treasury Bills'); ylabel('R [ percent per year ]'); xlabel('Time [ quarter ]'); grid subplot(2,2,4) hndl5=plot(t,Y,'r'); set(hndl5,'LineWidth',2); title('Real GNP (net export and import)'); ylabel('Y [ billions of 1982 dolars ]'); xlabel('Time [ quarter ]'); grid %graficke zobrazeni exogennich promennych figure set(gcf,'Name',' Exogenous variables'); subplot(2,1,1) hndl6=plot(t,G,'r'); set(hndl6,'LineWidth',2); hold off title('Real Gavernment Spending'); ylabel('G [ billions of 1982 dolars'); xlabel('Time [ quarter ]'); grid subplot(2,1,2) hndl7=plot(t,M,'r'); set(hndl7,'LineWidth',2); hold off title('Real Money Stock M1'); ylabel('M [ billions of 1982 dolars ]'); xlabel('Time [ quarter ]'); grid %odhad spotrebni funkce (13.14) C(t) = a1 + a2*Y(t) + a3*C(t-1) disp('******************'); disp('*Spotrebni funkce*'); disp('******************'); y = C; yendog = [Y]; %endogenni promenne v rovnici xexog = [konst C_1]; %exogenni promenne v rovnici %vsechny exogenni a zpozdenne endogenni promenne v systemu xall = [konst G dM C_1 Y_1 dY_1 R_4 sumR_1]; res2SLS_C = tsls(y,yendog,xexog,xall); resOLS_C = OLS(y,[yendog xexog]); vnames_C = strvcat('Spotreba','Y','konstanta','C_1'); prt_reg(res2SLS_C,vnames_C); prt_reg(resOLS_C,vnames_C); %odhad investicni funkce (13.15) I(t) = b1 + b2*[Y(t-1)-Y(t-2)] + b3*Y(t) + b4*R(t-4) disp('*******************'); disp('*Investicni funkce*'); disp('*******************'); y = I; yendog = [Y]; xexog = [konst dY_1 R_4]; res2SLS_I = tsls(y,yendog,xexog,xall); resOLS_I = OLS(y,[yendog xexog]); vnames_I = strvcat('Investice','Y','konstanta','dY_1','R_4'); prt_reg(res2SLS_I,vnames_I); prt_reg(resOLS_I,vnames_I); %odhad rovnice urokove miry (13.16) R(t) = c1 + c2*Y(t) + c3*[Y(t)-Y(t-1)] + c4*[M(t]-M(t-1)] + c5*[R(t-1)+R(t-2)] disp('**********************'); disp('*Rovnice urokove miry*'); disp('**********************'); y = R; yendog = [Y dY]; xexog = [konst dM sumR_1]; res2SLS_R = tsls(y,yendog,xexog,xall); resOLS_R = OLS(y,[yendog xexog]); vnames_R = strvcat('Ur. mira','Y','dY','konstanta','dM','sumR_1'); prt_reg(res2SLS_R,vnames_R); prt_reg(resOLS_R,vnames_R); neqs = 3; %struktura rovnic y obsahujici zavisle promenne yy(1).eq = C; yy(2).eq = I; yy(3).eq = R; %struktura rovnic Y obshujici endogenni promenne na prave strane YY(1).eq = [Y]; YY(2).eq = [Y]; YY(3).eq = [Y dY]; %struktura rovnic X obsahujici predeterminovane promenne na prave strane XX(1).eq = [konst C_1]; XX(2).eq = [konst dY_1 R_4]; XX(3).eq = [konst dM sumR_1]; res3SLS = thsls(neqs,yy,YY,XX,xall); vnames = strvcat(vnames_C,vnames_I,vnames_R); prt(res3SLS,vnames); %dopocitani Y res2SLS_Y.yhat = res2SLS_C.yhat + res2SLS_I.yhat + G; resOLS_Y.yhat = resOLS_C.yhat + resOLS_I.yhat + G; res3SLS(4).yhat = res3SLS(1).yhat + res3SLS(2).yhat + G; %graficke zobrazeni endogennich promennych figure set(gcf,'Name',' Endogenous variables'); subplot(2,2,1) plot(t,[C res2SLS_C.yhat resOLS_C.yhat res3SLS(1).yhat]); title('Real Aggregate Personal Consumption'); ylabel('C [ billions of 1982 dolars ]'); xlabel('Time [ quarter ]'); legend('actual','predicted 2SLS','predictedOLS','predicted 2SLS'); grid subplot(2,2,2) plot(t,[I res2SLS_I.yhat resOLS_I.yhat res3SLS(2).yhat]); title('Real Gross Domestic Investment'); ylabel('I [ billions of 1982 dolars ]'); xlabel('Time [ quarter ]'); legend('actual','predicted 2SLS','predictedOLS','predicted 2SLS'); grid subplot(2,2,3) plot(t,[R res2SLS_R.yhat resOLS_R.yhat res3SLS(3).yhat]); title('Interest Rate on 3-Month Treasury Bills'); ylabel('R [ percent per year ]'); xlabel('Time [ quarter ]'); legend('actual','predicted 2SLS','predictedOLS','predicted 2SLS'); grid subplot(2,2,4) plot(t,[Y res2SLS_Y.yhat resOLS_Y.yhat res3SLS(4).yhat]); title('Real GNP (net export and import)'); ylabel('Y [ billions of 1982 dolars ]'); xlabel('Time [ quarter ]'); legend('actual','predicted 2SLS','predictedOLS','predicted 2SLS'); grid %simulace modelu simBB*simY = simCC*simX => simYY = inv(simBB)*simCC*simXX %parametry 2SLS a11 = res2SLS_C.beta(2); a12 = res2SLS_C.beta(1); a13 = res2SLS_C.beta(3); a21 = res2SLS_I.beta(2); a22 = res2SLS_I.beta(3); a23 = res2SLS_I.beta(1); a24 = res2SLS_I.beta(4); a31 = res2SLS_R.beta(3); a32 = res2SLS_R.beta(1); a33 = res2SLS_R.beta(2); a34 = res2SLS_R.beta(4); a35 = res2SLS_R.beta(5); simBB =[1 0 0 -a12 0 0 1 0 -a23 0 0 0 1 -a32 -a33 -1 -1 0 1 0 0 0 0 -1 1]; simCC = [a11 0 0 a13 0 0 0 0 a21 0 0 0 a22 a24 0 0 a31 a34 0 0 0 0 a35 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 -1]; sim2SLS_C = zeros(n,1); sim2SLS_I = zeros(n,1); sim2SLS_R = zeros(n,1); sim2SLS_Y = zeros(n,1); sim2SLS_dY= zeros(n,1); for i=1:n, if i<2 sim2SLS_C_1(i) = C_1(i); sim2SLS_Y_1(i) = Y_1(i); else sim2SLS_C_1(i) = sim2SLS_C(i-1); sim2SLS_Y_1(i) = sim2SLS_Y(i-1); end if i<3 sim2SLS_dY_1(i) = dY_1(i); sim2SLS_sumR_1(i) = sumR_1(i); else sim2SLS_dY_1(i) = sim2SLS_dY(i-1); sim2SLS_sumR_1(i) = sim2SLS_R(i-1)+sim2SLS_R(i-2); end if i<5 sim2SLS_R_4(i) = R_4(i); else sim2SLS_R_4(i) = sim2SLS_R(i-4); end % simYY = [simGDP(i);simYPD(i);simTAX(i);simC(i);simdM(i);simINR(i);simIR(i);simINV(i);simRS(i);simRL(i);simWINF(i);simINFL(i);simM(i);simdGDP_C(i)] simXX = [1;dM(i);G(i);sim2SLS_C_1(i);sim2SLS_dY_1(i);sim2SLS_R_4(i);sim2SLS_sumR_1(i);sim2SLS_Y_1(i)]; simYY = inv(simBB)*simCC*simXX; sim2SLS_C(i) = simYY(1); sim2SLS_I(i) = simYY(2); sim2SLS_R(i) = simYY(3); sim2SLS_Y(i) = simYY(4); sim2SLS_dY(i) = simYY(5); end %simulace modelu imBB*simY = simCC*simX => simYY = inv(simBB)*simCC*simXX %parametry OLS a11 = resOLS_C.beta(2); a12 = resOLS_C.beta(1); a13 = resOLS_C.beta(3); a21 = resOLS_I.beta(2); a22 = resOLS_I.beta(3); a23 = resOLS_I.beta(1); a24 = resOLS_I.beta(4); a31 = resOLS_R.beta(3); a32 = resOLS_R.beta(1); a33 = resOLS_R.beta(2); a34 = resOLS_R.beta(4); a35 = resOLS_R.beta(5); simBB =[1 0 0 -a12 0 0 1 0 -a23 0 0 0 1 -a32 -a33 -1 -1 0 1 0 0 0 0 -1 1]; simCC = [a11 0 0 a13 0 0 0 0 a21 0 0 0 a22 a24 0 0 a31 a34 0 0 0 0 a35 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 -1]; simOLS_C = zeros(n,1); simOLS_I = zeros(n,1); simOLS_R = zeros(n,1); simOLS_Y = zeros(n,1); simOLS_dY= zeros(n,1); for i=1:n, if i<2 simOLS_C_1(i) = C_1(i); simOLS_Y_1(i) = Y_1(i); else simOLS_C_1(i) = simOLS_C(i-1); simOLS_Y_1(i) = simOLS_Y(i-1); end if i<3 simOLS_dY_1(i) = dY_1(i); simOLS_sumR_1(i) = sumR_1(i); else simOLS_dY_1(i) = simOLS_dY(i-1); simOLS_sumR_1(i) = simOLS_R(i-1)+simOLS_R(i-2); end if i<5 simOLS_R_4(i) = R_4(i); else simOLS_R_4(i) = simOLS_R(i-4); end % simYY = [simGDP(i);simYPD(i);simTAX(i);simC(i);simdM(i);simINR(i);simIR(i);simINV(i);simRS(i);simRL(i);simWINF(i);simINFL(i);simM(i);simdGDP_C(i)] simXX = [1;dM(i);G(i);simOLS_C_1(i);simOLS_dY_1(i);simOLS_R_4(i);simOLS_sumR_1(i);simOLS_Y_1(i)]; simYY = inv(simBB)*simCC*simXX; simOLS_C(i) = simYY(1); simOLS_I(i) = simYY(2); simOLS_R(i) = simYY(3); simOLS_Y(i) = simYY(4); simOLS_dY(i) = simYY(5); end %simulace modelu imBB*simY = simCC*simX => simYY = inv(simBB)*simCC*simXX %parametry 3SLS a11 = res3SLS(1).beta(2); a12 = res3SLS(1).beta(1); a13 = res3SLS(1).beta(3); a21 = res3SLS(2).beta(2); a22 = res3SLS(2).beta(3); a23 = res3SLS(2).beta(1); a24 = res3SLS(2).beta(4); a31 = res3SLS(3).beta(3); a32 = res3SLS(3).beta(1); a33 = res3SLS(3).beta(2); a34 = res3SLS(3).beta(4); a35 = res3SLS(3).beta(5); simBB =[1 0 0 -a12 0 0 1 0 -a23 0 0 0 1 -a32 -a33 -1 -1 0 1 0 0 0 0 -1 1]; simCC = [a11 0 0 a13 0 0 0 0 a21 0 0 0 a22 a24 0 0 a31 a34 0 0 0 0 a35 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 -1]; sim3SLS_C = zeros(n,1); sim3SLS_I = zeros(n,1); sim3SLS_R = zeros(n,1); sim3SLS_Y = zeros(n,1); sim3SLS_dY= zeros(n,1); for i=1:n, if i<2 sim3SLS_C_1(i) = C_1(i); sim3SLS_Y_1(i) = Y_1(i); else sim3SLS_C_1(i) = sim3SLS_C(i-1); sim3SLS_Y_1(i) = sim3SLS_Y(i-1); end if i<3 sim3SLS_dY_1(i) = dY_1(i); sim3SLS_sumR_1(i) = sumR_1(i); else sim3SLS_dY_1(i) = sim3SLS_dY(i-1); sim3SLS_sumR_1(i) = sim3SLS_R(i-1)+sim3SLS_R(i-2); end if i<5 sim3SLS_R_4(i) = R_4(i); else sim3SLS_R_4(i) = sim3SLS_R(i-4); end % simYY = [simGDP(i);simYPD(i);simTAX(i);simC(i);simdM(i);simINR(i);simIR(i);simINV(i);simRS(i);simRL(i);simWINF(i);simINFL(i);simM(i);simdGDP_C(i)] simXX = [1;dM(i);G(i);sim3SLS_C_1(i);sim3SLS_dY_1(i);sim3SLS_R_4(i);sim3SLS_sumR_1(i);sim3SLS_Y_1(i)]; simYY = inv(simBB)*simCC*simXX; sim3SLS_C(i) = simYY(1); sim3SLS_I(i) = simYY(2); sim3SLS_R(i) = simYY(3); sim3SLS_Y(i) = simYY(4); sim3SLS_dY(i) = simYY(5); end %Summary statistics for historical simulation v_names = ['C' 'I' 'R' 'Y']; sim2SLS(:,1) = sim2SLS_C; sim2SLS(:,2) = sim2SLS_I; sim2SLS(:,3) = sim2SLS_R; sim2SLS(:,4) = sim2SLS_Y; simOLS(:,1) = simOLS_C; simOLS(:,2) = simOLS_I; simOLS(:,3) = simOLS_R; simOLS(:,4) = simOLS_Y; sim3SLS(:,1) = sim3SLS_C; sim3SLS(:,2) = sim3SLS_I; sim3SLS(:,3) = sim3SLS_R; sim3SLS(:,4) = sim3SLS_Y; act(:,1) = C; act(:,2) = I; act(:,3) = R; act(:,4) = Y; for j=1:4 %mean stat2SLS(1,j) = mean(sim2SLS(:,j)); stat3SLS(1,j) = mean(sim3SLS(:,j)); statOLS(1,j) = mean(simOLS(:,j)); %rms error stat2SLS(2,j) = sqrt(mean((sim2SLS(:,j)-act(:,j)).^2)); stat3SLS(2,j) = sqrt(mean((sim3SLS(:,j)-act(:,j)).^2)); statOLS(2,j) = sqrt(mean((simOLS(:,j)-act(:,j)).^2)); %rms percent error stat2SLS(3,j) = sqrt(mean((sim2SLS(:,j)./act(:,j)-1).^2)); stat3SLS(3,j) = sqrt(mean((sim3SLS(:,j)./act(:,j)-1).^2)); statOLS(3,j) = sqrt(mean((simOLS(:,j)./act(:,j)-1).^2)); %mean error stat2SLS(4,j) = mean(sim2SLS(:,j)-act(:,j)); stat3SLS(4,j) = mean(sim3SLS(:,j)-act(:,j)); statOLS(4,j) = mean(simOLS(:,j)-act(:,j)); %mean percent error stat2SLS(5,j) = mean(sim2SLS(:,j)./act(:,j)-1); stat3SLS(5,j) = mean(sim3SLS(:,j)./act(:,j)-1); statOLS(5,j) = mean(simOLS(:,j)./act(:,j)-1); end fprintf('\n \n \n'); fprintf('********Summary statistics for historical simulation 3SLS************ \n'); fprintf(' %6c %6c %6c %6c \n',v_names); fprintf('Mean = %6.3f %6.3f %6.3f %6.3f \n',stat3SLS(1,:)); fprintf('rms error = %6.3f %6.3f %6.3f %6.3f \n',stat3SLS(2,:)); fprintf('rms percent error = %6.3f %6.3f %6.3f %6.3f \n',stat3SLS(3,:)); fprintf('mean error = %6.3f %6.3f %6.3f %6.3f \n',stat3SLS(4,:)); fprintf('mean percent error = %6.3f %6.3f %6.3f %6.3f \n',stat3SLS(5,:)); fprintf('********Summary statistics for historical simulation 2SLS************ \n'); fprintf(' %6c %6c %6c %6c \n',v_names); fprintf('Mean = %6.3f %6.3f %6.3f %6.3f \n',stat2SLS(1,:)); fprintf('rms error = %6.3f %6.3f %6.3f %6.3f \n',stat2SLS(2,:)); fprintf('rms percent error = %6.3f %6.3f %6.3f %6.3f \n',stat2SLS(3,:)); fprintf('mean error = %6.3f %6.3f %6.3f %6.3f \n',stat2SLS(4,:)); fprintf('mean percent error = %6.3f %6.3f %6.3f %6.3f \n',stat2SLS(5,:)); fprintf('\n********Summary statistics for historical simulation OLS************* \n'); fprintf(' %6c %6c %6c %6c \n',v_names); fprintf('Mean = %6.3f %6.3f %6.3f %6.3f \n',statOLS(1,:)); fprintf('rms error = %6.3f %6.3f %6.3f %6.3f \n',statOLS(2,:)); fprintf('rms percent error = %6.3f %6.3f %6.3f %6.3f \n',statOLS(3,:)); fprintf('mean error = %6.3f %6.3f %6.3f %6.3f \n',statOLS(4,:)); fprintf('mean percent error = %6.3f %6.3f %6.3f %6.3f \n',statOLS(5,:)); %ex post predpoved + 3 roky for_t = 12; last = 25; load USAmacro.dat t = USAmacro(maxl+1:end-last+for_t,1); % Puvoden t = ctvrtleti od 2Q1950-1Q1988 n = length(t); %endogeni promenne C = USAmacro(maxl+1:end-last+for_t,2); % Realna agregatni osobni potreba (v cenach roku 1982) [mld.] I = USAmacro(maxl+1:end-last+for_t,3); % Realne hrube domaci investice (v cenach roku 1982) [mld.] R = USAmacro(maxl+1:end-last+for_t,4); % Urokova sazba 3 mesicnich Treasury bill [% p.a.] Y = USAmacro(maxl+1:end-last+for_t,5); % Realny HDP (ocisten o exporty a importy; v cenach 1982) [mld.] %exogenni promenne konst = ones(n,1); G = USAmacro(maxl+1:end-last+for_t,6); % Realne vladni vydaje (v cenach roku 1982) [mld.] dM = USAmacro(maxl+1:end-last+for_t,7)-USAmacro(maxl+1-1:end-last+for_t-1,7); % Zmena realne penezni zasoby (M1) [mld.] %simulace modelu imBB*simY = simCC*simX => simYY = inv(simBB)*simCC*simXX %parametry 2SLS a11 = res2SLS_C.beta(2); a12 = res2SLS_C.beta(1); a13 = res2SLS_C.beta(3); a21 = res2SLS_I.beta(2); a22 = res2SLS_I.beta(3); a23 = res2SLS_I.beta(1); a24 = res2SLS_I.beta(4); a31 = res2SLS_R.beta(3); a32 = res2SLS_R.beta(1); a33 = res2SLS_R.beta(2); a34 = res2SLS_R.beta(4); a35 = res2SLS_R.beta(5); simBB =[1 0 0 -a12 0 0 1 0 -a23 0 0 0 1 -a32 -a33 -1 -1 0 1 0 0 0 0 -1 1]; simCC = [a11 0 0 a13 0 0 0 0 a21 0 0 0 a22 a24 0 0 a31 a34 0 0 0 0 a35 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 -1]; for2SLS_C = zeros(for_t,1); for2SLS_I = zeros(for_t,1); for2SLS_R = zeros(for_t,1); for2SLS_Y = zeros(for_t,1); for2SLS_dY= zeros(for_t,1); for i=1:for_t, if i<2 for2SLS_C_1(i) = C(end-for_t); for2SLS_Y_1(i) = Y(end-for_t); for2SLS_dY_1(i) = Y(end-for_t)-Y(end-for_t-1); for2SLS_sumR_1(i) = R(end-for_t)+R(end-for_t-1); else for2SLS_C_1(i) = for2SLS_C(i-1); for2SLS_Y_1(i) = for2SLS_Y(i-1); end if i==2 for2SLS_dY_1(i) = for2SLS_Y(i-1)-Y(end-for_t); for2SLS_sumR_1(i) = for2SLS_R(i-1)+R(end-for_t); end if i>2 for2SLS_dY_1(i) = for2SLS_dY(i-1); for2SLS_sumR_1(i) = for2SLS_R(i-1)+for2SLS_R(i-2); end if i<5 for2SLS_R_4(i) = R(end-for_t-4+i); else for2SLS_R_4(i) = for2SLS_R(i-4); end % simYY = [simGDP(i);simYPD(i);simTAX(i);simC(i);simdM(i);simINR(i);simIR(i);simINV(i);simRS(i);simRL(i);simWINF(i);simINFL(i);simM(i);simdGDP_C(i)] simXX = [1;dM(end-for_t+i);G(end-for_t+i);for2SLS_C_1(i);for2SLS_dY_1(i);for2SLS_R_4(i);for2SLS_sumR_1(i);for2SLS_Y_1(i)]; simYY = inv(simBB)*simCC*simXX; for2SLS_C(i) = simYY(1); for2SLS_I(i) = simYY(2); for2SLS_R(i) = simYY(3); for2SLS_Y(i) = simYY(4); for2SLS_dY(i) = simYY(5); end %simulace modelu imBB*simY = simCC*simX => simYY = inv(simBB)*simCC*simXX %parametry OLS a11 = resOLS_C.beta(2); a12 = resOLS_C.beta(1); a13 = resOLS_C.beta(3); a21 = resOLS_I.beta(2); a22 = resOLS_I.beta(3); a23 = resOLS_I.beta(1); a24 = resOLS_I.beta(4); a31 = resOLS_R.beta(3); a32 = resOLS_R.beta(1); a33 = resOLS_R.beta(2); a34 = resOLS_R.beta(4); a35 = resOLS_R.beta(5); simBB =[1 0 0 -a12 0 0 1 0 -a23 0 0 0 1 -a32 -a33 -1 -1 0 1 0 0 0 0 -1 1]; simCC = [a11 0 0 a13 0 0 0 0 a21 0 0 0 a22 a24 0 0 a31 a34 0 0 0 0 a35 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 -1]; forOLS_C = zeros(for_t,1); forOLS_I = zeros(for_t,1); forOLS_R = zeros(for_t,1); forOLS_Y = zeros(for_t,1); forOLS_dY= zeros(for_t,1); for i=1:for_t, if i<2 forOLS_C_1(i) = C(end-for_t); forOLS_Y_1(i) = Y(end-for_t); forOLS_dY_1(i) = Y(end-for_t)-Y(end-for_t-1); forOLS_sumR_1(i) = R(end-for_t)+R(end-for_t-1); else forOLS_C_1(i) = forOLS_C(i-1); forOLS_Y_1(i) = forOLS_Y(i-1); end if i==2 forOLS_dY_1(i) = forOLS_Y(i-1)-Y(end-for_t); forOLS_sumR_1(i) = forOLS_R(i-1)+R(end-for_t); end if i>2 forOLS_dY_1(i) = forOLS_dY(i-1); forOLS_sumR_1(i) = forOLS_R(i-1)+forOLS_R(i-2); end if i<5 forOLS_R_4(i) = R(end-for_t-4+i); else forOLS_R_4(i) = forOLS_R(i-4); end % simYY = [simGDP(i);simYPD(i);simTAX(i);simC(i);simdM(i);simINR(i);simIR(i);simINV(i);simRS(i);simRL(i);simWINF(i);simINFL(i);simM(i);simdGDP_C(i)] simXX = [1;dM(end-for_t+i);G(end-for_t+i);forOLS_C_1(i);forOLS_dY_1(i);forOLS_R_4(i);forOLS_sumR_1(i);forOLS_Y_1(i)]; simYY = inv(simBB)*simCC*simXX; forOLS_C(i) = simYY(1); forOLS_I(i) = simYY(2); forOLS_R(i) = simYY(3); forOLS_Y(i) = simYY(4); forOLS_dY(i) = simYY(5); end %simulace modelu imBB*simY = simCC*simX => simYY = inv(simBB)*simCC*simXX %parametry 3SLS a11 = res3SLS(1).beta(2); a12 = res3SLS(1).beta(1); a13 = res3SLS(1).beta(3); a21 = res3SLS(2).beta(2); a22 = res3SLS(2).beta(3); a23 = res3SLS(2).beta(1); a24 = res3SLS(2).beta(4); a31 = res3SLS(3).beta(3); a32 = res3SLS(3).beta(1); a33 = res3SLS(3).beta(2); a34 = res3SLS(3).beta(4); a35 = res3SLS(3).beta(5); simBB =[1 0 0 -a12 0 0 1 0 -a23 0 0 0 1 -a32 -a33 -1 -1 0 1 0 0 0 0 -1 1]; simCC = [a11 0 0 a13 0 0 0 0 a21 0 0 0 a22 a24 0 0 a31 a34 0 0 0 0 a35 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 -1]; for3SLS_C = zeros(for_t,1); for3SLS_I = zeros(for_t,1); for3SLS_R = zeros(for_t,1); for3SLS_Y = zeros(for_t,1); for3SLS_dY= zeros(for_t,1); for i=1:for_t, if i<2 for3SLS_C_1(i) = C(end-for_t); for3SLS_Y_1(i) = Y(end-for_t); for3SLS_dY_1(i) = Y(end-for_t)-Y(end-for_t-1); for3SLS_sumR_1(i) = R(end-for_t)+R(end-for_t-1); else for3SLS_C_1(i) = for3SLS_C(i-1); for3SLS_Y_1(i) = for3SLS_Y(i-1); end if i==2 for3SLS_dY_1(i) = for3SLS_Y(i-1)-Y(end-for_t); for3SLS_sumR_1(i) = for3SLS_R(i-1)+R(end-for_t); end if i>2 for3SLS_dY_1(i) = for3SLS_dY(i-1); for3SLS_sumR_1(i) = for3SLS_R(i-1)+for3SLS_R(i-2); end if i<5 for3SLS_R_4(i) = R(end-for_t-4+i); else for3SLS_R_4(i) = for3SLS_R(i-4); end % simYY = [simGDP(i);simYPD(i);simTAX(i);simC(i);simdM(i);simINR(i);simIR(i);simINV(i);simRS(i);simRL(i);simWINF(i);simINFL(i);simM(i);simdGDP_C(i)] simXX = [1;dM(end-for_t+i);G(end-for_t+i);for3SLS_C_1(i);for3SLS_dY_1(i);for3SLS_R_4(i);for3SLS_sumR_1(i);for3SLS_Y_1(i)]; simYY = inv(simBB)*simCC*simXX; for3SLS_C(i) = simYY(1); for3SLS_I(i) = simYY(2); for3SLS_R(i) = simYY(3); for3SLS_Y(i) = simYY(4); for3SLS_dY(i) = simYY(5); end %graficke zobrazeni nasimulovanych a predpovezenych endogennich promennych figure set(gcf,'Name',' Endogenous variables'); subplot(2,2,1) plot(t(1:end-for_t),[C(1:end-for_t) sim2SLS_C simOLS_C sim3SLS_C]); title('Real Aggregate Personal Consumption'); ylabel('C [ billions of 1982 dolars ]'); xlabel('Time [ quarter ]'); legend('actual','simulated 2SLS','simulated OLS','simulated 3SLS'); grid subplot(2,2,2) plot(t(end-for_t+1:end),[C(end-for_t+1:end) for2SLS_C forOLS_C for3SLS_C]); title('Real Aggregate Personal Consumption'); ylabel('C [ billions of 1982 dolars ]'); xlabel('Time [ quarter ]'); legend('actual','forecast 2SLS','forecast OLS','forecast 3SLS'); grid subplot(2,2,3) plot(t(1:end-for_t),[I(1:end-for_t) sim2SLS_I simOLS_I sim3SLS_I]); title('Real Gross Domestic Investment'); ylabel('I [ billions of 1982 dolars ]'); xlabel('Time [ quarter ]'); legend('actual','simulated 2SLS','simulated OLS','simulated 3SLS'); grid subplot(2,2,4) plot(t(end-for_t+1:end),[I(end-for_t+1:end) for2SLS_I forOLS_I for3SLS_I]); title('Real Gross Domestic Investment'); ylabel('I [ billions of 1982 dolars ]'); xlabel('Time [ quarter ]'); legend('actual','forecast 2SLS','forecast OLS','forecast 3SLS'); grid figure set(gcf,'Name',' Endogenous variables'); subplot(2,2,1) plot(t(1:end-for_t),[R(1:end-for_t) sim2SLS_R simOLS_R sim3SLS_R]); title('Interest Rate on 3-Month Treasury Bills'); ylabel('R [ percent per year ]'); xlabel('Time [ quarter ]'); legend('actual','simulated 2SLS','simulated OLS','simulated 3SLS'); grid subplot(2,2,2) plot(t(end-for_t+1:end),[R(end-for_t+1:end) for2SLS_R forOLS_R for3SLS_R]); title('Interest Rate on 3-Month Treasury Bills'); ylabel('R [ percent per year ]'); xlabel('Time [ quarter ]'); legend('actual','forecast 2SLS','forecast OLS','forecast 3SLS'); grid subplot(2,2,3) plot(t(1:end-for_t),[Y(1:end-for_t) sim2SLS_Y simOLS_Y sim3SLS_Y]); title('Real GNP (net export and import)'); ylabel('Y [ billions of 1982 dolars ]'); xlabel('Time [ quarter ]'); legend('actual','simulated 2SLS','simulated OLS','simulated 3SLS'); grid subplot(2,2,4) plot(t(end-for_t+1:end),[Y(end-for_t+1:end) for2SLS_Y forOLS_Y for3SLS_Y]); title('Real GNP (net export and import)'); ylabel('Y [ billions of 1982 dolars ]'); xlabel('Time [ quarter ]'); legend('actual','forecast 2SLS','forecast OLS','forecast 3SLS'); grid %Summary statistics for ex post forecast 4/1973-4/1976 v_names = ['C' 'I' 'R' 'Y']; act=[]; for2SLS(:,1) = for2SLS_C; for2SLS(:,2) = for2SLS_I; for2SLS(:,3) = for2SLS_R; for2SLS(:,4) = for2SLS_Y; forOLS(:,1) = forOLS_C; forOLS(:,2) = forOLS_I; forOLS(:,3) = forOLS_R; forOLS(:,4) = forOLS_Y; for3SLS(:,1) = for3SLS_C; for3SLS(:,2) = for3SLS_I; for3SLS(:,3) = for3SLS_R; for3SLS(:,4) = for3SLS_Y; act(:,1) = C(end-for_t+1:end); act(:,2) = I(end-for_t+1:end); act(:,3) = R(end-for_t+1:end); act(:,4) = Y(end-for_t+1:end); for j=1:4 %mean stat2SLS(1,j) = mean(for2SLS(:,j)); stat3SLS(1,j) = mean(for3SLS(:,j)); statOLS(1,j) = mean(forOLS(:,j)); %rms error stat2SLS(2,j) = sqrt(mean((for2SLS(:,j)-act(:,j)).^2)); stat3SLS(2,j) = sqrt(mean((for3SLS(:,j)-act(:,j)).^2)); statOLS(2,j) = sqrt(mean((forOLS(:,j)-act(:,j)).^2)); %rms percent error stat2SLS(3,j) = sqrt(mean((for2SLS(:,j)./act(:,j)-1).^2)); stat3SLS(3,j) = sqrt(mean((for3SLS(:,j)./act(:,j)-1).^2)); statOLS(3,j) = sqrt(mean((forOLS(:,j)./act(:,j)-1).^2)); %mean error stat2SLS(4,j) = mean(for2SLS(:,j)-act(:,j)); stat3SLS(4,j) = mean(for3SLS(:,j)-act(:,j)); statOLS(4,j) = mean(forOLS(:,j)-act(:,j)); %mean percent error stat2SLS(5,j) = mean(for2SLS(:,j)./act(:,j)-1); stat3SLS(5,j) = mean(for3SLS(:,j)./act(:,j)-1); statOLS(5,j) = mean(forOLS(:,j)./act(:,j)-1); end fprintf('\n \n \n'); fprintf('********Summary statistics for ex post forecast 3SLS************ \n'); fprintf(' %6c %6c %6c %6c \n',v_names); fprintf('Mean = %6.3f %6.3f %6.3f %6.3f \n',stat3SLS(1,:)); fprintf('rms error = %6.3f %6.3f %6.3f %6.3f \n',stat3SLS(2,:)); fprintf('rms percent error = %6.3f %6.3f %6.3f %6.3f \n',stat3SLS(3,:)); fprintf('mean error = %6.3f %6.3f %6.3f %6.3f \n',stat3SLS(4,:)); fprintf('mean percent error = %6.3f %6.3f %6.3f %6.3f \n',stat3SLS(5,:)); fprintf('\n********Summary statistics for ex post forecast 2SLS************ \n'); fprintf(' %6c %6c %6c %6c \n',v_names); fprintf('Mean = %6.3f %6.3f %6.3f %6.3f \n',stat2SLS(1,:)); fprintf('rms error = %6.3f %6.3f %6.3f %6.3f \n',stat2SLS(2,:)); fprintf('rms percent error = %6.3f %6.3f %6.3f %6.3f \n',stat2SLS(3,:)); fprintf('mean error = %6.3f %6.3f %6.3f %6.3f \n',stat2SLS(4,:)); fprintf('mean percent error = %6.3f %6.3f %6.3f %6.3f \n',stat2SLS(5,:)); fprintf('\n********Summary statistics for ex post forecast OLS************* \n'); fprintf(' %6c %6c %6c %6c \n',v_names); fprintf('Mean = %6.3f %6.3f %6.3f %6.3f \n',statOLS(1,:)); fprintf('rms error = %6.3f %6.3f %6.3f %6.3f \n',statOLS(2,:)); fprintf('rms percent error = %6.3f %6.3f %6.3f %6.3f \n',statOLS(3,:)); fprintf('mean error = %6.3f %6.3f %6.3f %6.3f \n',statOLS(4,:)); fprintf('mean percent error = %6.3f %6.3f %6.3f %6.3f \n',statOLS(5,:));