%% Infer Innovations from Model Containing Regression Component
% Consider the model and data in Infer VEC Model Innovations.
% Load the Data_USEconVECModel data set and preprocess the data.

load Data_USEconVECModel
FRED.GDP = 100*log(FRED.GDP);      
FRED.GDPDEF = 100*log(FRED.GDPDEF);
FRED.COE = 100*log(FRED.COE);       
FRED.HOANBS = 100*log(FRED.HOANBS); 
FRED.PCEC = 100*log(FRED.PCEC);     
FRED.GPDI = 100*log(FRED.GPDI);

%
% The Data_Recessions data set contains the beginning and ending serial dates of recessions. Load this data set. Convert the matrix of date serial numbers to a datetime array.

load Data_Recessions
dtrec = datetime(Recessions,'ConvertFrom','datenum');

%
% Create a dummy variable that identifies periods in which the U.S. was in a recession or worse. Specifically, the variable should be 1 if FRED.Time occurs during a recession, and 0 otherwise.

isin = @(x)(any(dtrec(:,1) <= x & x <= dtrec(:,2)));
isrecession = double(arrayfun(isin,FRED.Time));

%
% Create a VEC(1) model using the shorthand syntax. Assume that the appropriate cointegration rank is 4. You do not have to specify the presence of a regression component when creating the model. Specify the variable names.

Mdl = vecm(7,4,1);
Mdl.SeriesNames = FRED.Properties.VariableNames;

%
% Estimate the model using the entire sample. Specify the predictor identifying whether the observation was measured during a recession. Return the standard errors.

EstMdl = estimate(Mdl,FRED.Variables,'X',isrecession);

% Infer innovations from the estimated model. Supply the predictor data. Return the loglikelihood objective function value.
%

[E,logL] = infer(EstMdl,FRED.Variables,'X',isrecession);
logL

% E is a 238-by-7 matrix of inferred innovations.
% Plot the residuals on separate plots. Synchronize the residuals with the dates by removing the first Mdl.P dates.
%

idx = FRED.Time((EstMdl.P + 1):end);

figure;
subplot(2,2,1)
plot(idx,E(:,1));
title('Residuals: GDP');
ylabel('Index (scaled)');
xlabel('Date');
subplot(2,2,2)
plot(idx,E(:,2));
title('Residuals: GDP Deflator');
ylabel('Index (scaled');
xlabel('Date');
subplot(2,2,3)
plot(idx,E(:,3));
title('Residuals: PCE');
ylabel('Billions of $ (scaled)');
xlabel('Date');
subplot(2,2,4)
plot(idx,E(:,4));
title('Residuals: NBSH');
ylabel('Index (scaled)');
xlabel('Date');

figure;
subplot(2,2,1)
plot(idx,E(:,5));
title('Residuals: Federal Funds Rate');
ylabel('Percent');
xlabel('Date');
subplot(2,2,2)
plot(idx,E(:,6));
title('Residuals: COE');
ylabel('Billions of $ (scaled)');
xlabel('Date');
subplot(2,2,3)
plot(idx,E(:,7));
title('Residuals: GPDI');
ylabel('Billions of $ (scaled)');
xlabel('Date');

% The residuals corresponding to the federal funds rate exhibit heteroscedasticity.
% Copyright 2012 The MathWorks, Inc.