%% Specify Coefficient Starting Values for Simulation Smoother
% Consider the 3-D VAR(4) model of docid:econ_ug.mw_93ab52c2-232a-4877-8152-03c57b96825f. In this case, assume that the prior model is semiconjugate.
% Load the US macroeconomic data set.  Compute the inflation rate, stabilize the unemployment and federal funds rates, and remove missing values (the data includes leading missing values only).

load Data_USEconModel
seriesnames = ["INFL" "UNRATE" "FEDFUNDS"];
DataTable.INFL = 100*[NaN; price2ret(DataTable.CPIAUCSL)];
DataTable.DUNRATE = [NaN; diff(DataTable.UNRATE)];
DataTable.DFEDFUNDS = [NaN; diff(DataTable.FEDFUNDS)];
seriesnames(2:3) = "D" + seriesnames(2:3);
rmDataTable = rmmissing(DataTable); 

% Create a semiconjugate prior model. Specify the response series names.

numseries = numel(seriesnames);
numlags = 4;

PriorMdl = bayesvarm(numseries,numlags,'ModelType','semiconjugate',...
    'SeriesNames',seriesnames);

% To obtain starting values for the coefficients, consider the coefficients from a VAR model fitted to the first 30 observations.
% Define index sets that partition the data into four sets:
% Length p = 4 initialization period for the dynamic component of the model
% n_0 = 30 observations to estimate the VAR for the coefficient initial values
% Length p = 4 initialization period for the dynamic component of the Bayesian VAR model
% Remaining n_1 observations to estimate the posterior

T = size(rmDataTable);
n0 = 30;
idxpre0 = 1:PriorMdl.P;
idxest0 = (idxpre0(end) + 1):(idxpre0(end) + 1 + n0);
idxpre1 = (idxest0(end) + 1 - PriorMdl.P):idxest0(end);
idxest1 = (idxest0(end) + 1):T;
n1 = numel(idxest1);

% Obtain initial coefficient values by fitting a VAR model to the first 34 observations.  Use the first four observations as a presample.

Mdl0 = varm(PriorMdl.NumSeries,PriorMdl.P);
EstMdl0 = estimate(Mdl0,rmDataTable{idxest0,seriesnames},'Y0',rmDataTable{idxpre0,seriesnames});

% estimate returns a model object containing the AR coefficient estimates in matrix form and the constants in a vector. The Bayesian VAR functions require initial coefficient values in a vector.  Format the initial coefficient values.

Coeff0 = [EstMdl0.AR{:} EstMdl0.Constant].';
Coeff0 = Coeff0(:);

% Simulate 1000 coefficients and innovations covariance matrices from the posterior distribution. Specify the remaining observations from which to estimate the posterior. Initialize the VAR model by specifying the last four observations in the previous estimation sample as a presample, and specify the initial coefficient vector to initialize the posterior sample.

rng(1); % For reproducibility
[Coeff,Sigma] = simsmooth(PriorMdl,rmDataTable{idxest1,seriesnames},...
    'Y0',rmDataTable{idxpre1,seriesnames},'Coeff0',Coeff0);

% Coeff is a 39-by-1000 matrix of coefficients and Sigma is a 3-by-3-by-1000 array of innovations covariance matrices. Both Coeff and Sigma are randomly drawn from the posterior distribution.
% Copyright 2020 The MathWorks, Inc.