clear all; close all; clc; %% Setup %Figure grid global nrows ncols fig nrows = 3; % Number of figure rows ncols = 3; % Number of figure columns fig = 0; % Initialize figure counter warning('off','MATLAB:print:FigureTooLargeForPage'); % Figure parameters lw = 1; % line width ms = 250; % marker size fs = 11; % font size %Set up figure %hold on; set(gca,'TickLabelInterpreter','latex'); set(gca,'FontSize',fs); set(gca,'Clipping','off'); set(gca,'Position',[0.10 0.15 0.80 0.70]); set(gca,'TickDir','out'); %Import data RGDPall = importdata('GDPalldata.xls','',1); GDPall = RGDPall.data; lnRGDP = log(GDPall); GDP = GDPall(177:292); %Create vector of GDP values that only include data %from 1991-2019 lnGDP = log(GDP); % Create vector that aligns the dates to the data time = (1991.0:0.25:2019.75)'; T = size(time,1); k = 2; %% Regression constant = ones(T,1); %x = sig*randn(T,1); xmat = [constant,time]; b = inv(xmat'*xmat)*(xmat'*lnGDP); e = lnGDP - xmat*b; %sig = sqrt(e'*e)/(T-k); sig = 0.03; %% Monte Carlo Experiment of Recursive Residuals % Parameters b1 = [1;0.0075]; b2 = [1;0.0074]; r = (1:T)'; % Variables y = zeros(T,1); % Simulate data before the break for tt = (1:1:74) y(tt) = xmat(tt,:)*b1 + sig*randn(1,1); end % Simulate data after the break for tt = (74+1:1:T) y(tt) = xmat(tt,:)*b2 + sig*randn(1,1); end % Calculate recursive residuals rr = zeros(T,1); w = zeros(T,1); for tt = (1:1:T) b = inv(xmat(1:tt-1,:)'*xmat(1:tt-1,:))*(xmat(1:tt-1,:)'*y(1:tt-1)); rr(tt) = y(tt) - xmat(tt,:)*b; s2 = (y(1:tt-1)-xmat(1:tt-1,:)*b)'*(y(1:tt-1)-xmat(1:tt-1,:)*b)/(tt-1-3); w(tt) = rr(tt)/sqrt(s2*(1+xmat(tt,:)*inv(xmat(1:tt-1,:)'*xmat(1:tt-1,:))*xmat(tt,:)')); end % Thresholds lb = -2.64*ones(T,1); ub = 2.64*ones(T,1); % Plot recursive residuals and thresholds newfig; hold on plot(r,lb,'k'); plot(r,ub,'k'); plot(r,w,'k-'); %Create axes labels and title xe = 0.01*(2019.75-1991.0); ye = 0.01*(0.1+0.1); text(120+xe,-7,'Quarters',... 'HorizontalAlignment','left','VerticalAlignment','middle','Interpreter','LaTeX'); text(0,3+ye,'Residuals',... 'HorizontalAlignment','left','VerticalAlignment','bottom','Interpreter','LaTeX'); fs = 13; text(60,3+10*ye,... 'Monte Carlo Experiment of Recursive Residuals of Log of Real GDP from $1991Q1 - 2019Q4$',... 'HorizontalAlignment','center','VerticalAlignment','bottom',... 'Interpreter','LaTeX','FontSize',fs); hold off %Plot recursive residuals and thresholds against time newfig; hold on plot(time,lb,'k'); plot(time,ub,'k'); plot(time,w,'k-'); %Create axes labels and title xe = 0.01*(2019.75-1991.0); ye = 0.01*(0.1+0.1); text(2020.0+xe,-7,'Time',... 'HorizontalAlignment','left','VerticalAlignment','middle','Interpreter','LaTeX'); text(1990.0,3+ye,'Residuals',... 'HorizontalAlignment','left','VerticalAlignment','bottom','Interpreter','LaTeX'); fs = 13; text(2005,3+10*ye,... 'Monte Carlo Experiment of Recursive Residuals of Log of Real GDP from $1991Q1 - 2019Q4$',... 'HorizontalAlignment','center','VerticalAlignment','bottom',... 'Interpreter','LaTeX','FontSize',fs); hold off %% Recursive Residuals %Calculate the standard deviation of the data: sig %sig = sqrt(var(lnGDP)); sig = sqrt(e'*e)/(T-k); %Calculate x, the vector of regressors %x = sig.*rand(T,1); %Calculate the X matrix constant = ones(T,1); xmat = [constant time]; %Calculate recursive residuals %Identify number of coefficients, k k = 2; %Initialize residuals e = zeros(T,1); %Initialize scaled residual w = zeros(T,1); for tt = (1:1:T) %Calculate coefficient estimates, b b = inv(xmat(1:tt-1,:)'*xmat(1:tt-1,:))*(xmat(1:tt-1,:)'*lnGDP(1:tt-1)); %Calculate residuals, e e(tt) = lnGDP(tt)-xmat(tt,:)*b; %Calculate sum of squared errors, s2 s2 = (lnGDP(1:tt-1)-xmat(1:tt-1,:)*b)'*(lnGDP(1:tt-1)-xmat(1:tt-1,:)*b)/(tt-1-k); %Calculate the scaled residual, w w(tt) = e(tt)/sqrt(1+xmat(tt,:)*inv(xmat(1:tt-1,:)'*xmat(1:tt-1,:))*xmat(tt,:)'); end %Calculate the error bounds %Calculation of lower bound, lb lb = (-2.*sqrt(s2).*sqrt(1+xmat(tt,:)*inv(xmat(1:tt-1,:)'*xmat(1:tt-1,:))*xmat(tt,:)'))... *ones(T,1); %Calculation of upper bound, ub ub = (2.*sqrt(s2).*sqrt(1+xmat(tt,:)*inv(xmat(1:tt-1,:)'*xmat(1:tt-1,:))*xmat(tt,:)'))... *ones(T,1); %Plot recursive residuals and error bounds newfig; hold on %Plot lower bound plot(time,lb,'k'); %Plot upper bound plot(time,ub,'k'); %Plot recursive residuals plot(time,w,'k-'); %Create axes labels and title xe = 0.01*(2019.75-1991.0); ye = 0.01*(0.1+0.1); text(2020.0+xe,-0.1,'Time',... 'HorizontalAlignment','left','VerticalAlignment','middle','Interpreter','LaTeX'); text(1990.0,0.1+ye,'Residuals',... 'HorizontalAlignment','left','VerticalAlignment','bottom','Interpreter','LaTeX'); fs = 13; text(2005,0.1+ye,'Recursive Residuals of Log of Real GDP from $1991Q1 - 2019Q4$',... 'HorizontalAlignment','center','VerticalAlignment','bottom',... 'Interpreter','LaTeX','FontSize',fs); hold off %% Recusive Residuals of data from 1947-2020 % Create vector that aligns the dates to the data Rtime = (1947.0:0.25:2020.5)'; RT = size(Rtime,1); %Calculate the standard deviation of the data: sig %sig = sqrt(var(lnGDP)); sig = sqrt(e'*e)/(T-k); %Calculate x, the vector of regressors %x = sig.*rand(RT,1); %Calculate the X matrix constant = ones(RT,1); xmat = [constant Rtime]; %Calculate recursive residuals %Identify number of coefficients, k k = 2; %Initialize residuals e = zeros(RT,1); %Initialize scaled residual w = zeros(RT,1); for tt = (1:1:RT) %Calculate coefficient estimates, b b = inv(xmat(1:tt-1,:)'*xmat(1:tt-1,:))*(xmat(1:tt-1,:)'*lnRGDP(1:tt-1)); %Calculate residuals, e e(tt) = lnRGDP(tt)-xmat(tt,:)*b; %Calculate sum of squared errors, s2 s2 = (lnRGDP(1:tt-1)-xmat(1:tt-1,:)*b)'*(lnRGDP(1:tt-1)-xmat(1:tt-1,:)*b)/(tt-1-k); %Calculate the scaled residual, w w(tt) = e(tt)/sqrt(1+xmat(tt,:)*inv(xmat(1:tt-1,:)'*xmat(1:tt-1,:))*xmat(tt,:)'); end %Calculate the error bounds %Calculation of lower bound, lb lb = (-2.*sqrt(s2).*sqrt(1+xmat(tt,:)*inv(xmat(1:tt-1,:)'*xmat(1:tt-1,:))*xmat(tt,:)'))... *ones(RT,1); %Calculation of upper bound, ub ub = (2.*sqrt(s2).*sqrt(1+xmat(tt,:)*inv(xmat(1:tt-1,:)'*xmat(1:tt-1,:))*xmat(tt,:)'))... *ones(RT,1); %Plot recursive residuals and error bounds newfig; hold on %Plot lower bound plot(Rtime,lb,'k'); %Plot upper bound plot(Rtime,ub,'k'); %Plot recursive residuals plot(Rtime,w,'k-'); %Create axes labels and title xe = 0.01*(2019.75-1991.0); ye = 0.01*(0.1+0.1); text(2030.0+xe,-0.3,'Time',... 'HorizontalAlignment','left','VerticalAlignment','middle','Interpreter','LaTeX'); text(1940.0,0.15+ye,'Residuals',... 'HorizontalAlignment','left','VerticalAlignment','bottom','Interpreter','LaTeX'); fs = 13; text(1985,0.15+10*ye,'Recursive Residuals of Log of Real GDP from $1947Q1 - 2020Q3$',... 'HorizontalAlignment','center','VerticalAlignment','bottom',... 'Interpreter','LaTeX','FontSize',fs); hold off %% Functions %Setup of figure grid function figh = newfig global nrows ncols fig %Get the computer's screensize as a vector screensize = get(groot,'ScreenSize'); screenwidth = screensize(3); screenheight = 0.975*screensize(4); figwidth = screenwidth/ncols; figheight = screenheight/nrows; %Increment the global figure counter fig = fig + 1; figh = figure(fig); set(figh,'OuterPosition', ... [mod(fig-1,ncols)*figwidth, screenheight - ceil(fig/ncols)*figheight, ... figwidth, figheight]); end %savefig.m %Save a figure as pdf, with whitespace around it minimized function savefig(figname) figh = gcf; figh.PaperPositionMode = 'auto'; fig_pos = figh.PaperPosition; figh.PaperSize = [fig_pos(3) fig_pos(4)]; print(figh, figname,'-dpdf'); end