%% Simulate Univariate Markov-Switching Dynamic Regression Model
% This example shows how to generate random response and state paths from a two Markov-switching dynamic regression model.
% Suppose that an economy switches between an expansion and a recession, with transition matrix
% \[P = \left[ {\begin{array}{*{20}{c}}
% {0.9}&{0.1}\\
% {0.3}&{0.7}
% \end{array}} \right].\]
% Suppose that y_t is a univariate response process representing an economic measurement that can suggest which state the economy experiences during a period. During an expansion, y_t is this AR(2) model:
% y_t=5+0.3y_{t-1}+0.2y_{t-2}+\varepsilon_{1t},
% where \varepsilon_{1t} is an iid Gaussian process with mean 0 and variance 2. During a recession, y_t is this AR(1) model:
% y_t=-5+0.1y_{t-1}+\varepsilon_{2t},
% where \varepsilon_{2t} is an iid Gaussian process with mean 0 and variance 1.
%
%% Create a Fully Specified Model
% Create the Markov-switching dynamic regression model that describes the state of the economy with respect to y_t (for more details, see docid:econ_ug.mw_51881501-44a1-49a1-93d4-4b1ed234ae61).

% Switching mechanism
P = [0.9 0.1; 0.3 0.7];
mc = dtmc(P,'StateNames',["Expansion" "Recession"]);

% AR submodels 
mdl1 = arima('Constant',5,'AR',[0.3 0.2],...
    'Variance',2,'Description','Expansion State');
mdl2 = arima('Constant',-5,'AR',0.1,...
    'Variance',1,'Description','Recession State');

% Markov-switching model
Mdl = msVAR(mc,[mdl1; mdl2]);

% Mdl is a fully specified msVAR object.
%
%% Simulate One Path
% Simulate one path of responses, innovations, and state indices from the model. Specify a 50-period simulation horizon, that is, a 50-observation path.

rng('default')
[y,e,sp] = simulate(Mdl,50);

% y, e, and sp are 50-by-1 vectors of simulated responses, innovations, and state indices, respectively.
% Plot the simulated observations in separate subplots.

figure
subplot(3,1,1)
plot(y)
ylabel('Response')
grid on

subplot(3,1,2)
plot(e)
ylabel('Innovation') 
grid on

subplot(3,1,3)
plot(sp,'m')
ylabel('State Index')
yticks([1 2])

%
%% Simulate Multiple Paths
% Simulate 3 separate, independent paths of responses, innovations, and state indices from the model. Specify a 5-period simulation horizon.

[Y,E,SP] = simulate(Mdl,5,'NumPaths',3)
Mdl.StateNames(SP(2,1))

% Y, E, and SP are 5-by-3 matrices of simulated responses, innovations, and state indices, respectively.
% The simulated response in the second period of the first path is -4.1065.
% The corresponding simulated innovation in the second period of the first path is 1.4367.
% The corresponding simulated state of the economy in the second period of the first path is a recession.
%
%% Simulate Model Containing Regression Component
% Include a known regression component in each submodel. During an expansion, the ARX(2) submodel is:
% y_t=5+0.3y_{t-1}+0.2y_{t-2}+2x_t+\varepsilon_{1t},
% and, during a recession, the ARX(1) submodel is:
% y_t=-5+0.1y_{t-1}+x_t + \varepsilon_{2t}.
% Specify the regression coefficient values of each submodel by setting the Beta property using dot notation. Then, create a new Markov-switching model from the modified ARX submodels and the existing Markov chain mc.

mdl1.Beta = 2;
mdl2.Beta = 1;
Mdl = msVAR(mc,[mdl1; mdl2]);

% simulate requires exogenous data over the simulation horizon for the regression component. Generate a 50-by-1 matrix of iid standard normal observations for the required exogenous data.

x = randn(50,1);

% Simulate one path of responses, innovations, and state indices from the model. Specify a 50-period simulation horizon. Plot the results.

[y,e,sp] = simulate(Mdl,50,'X',x);

figure
subplot(3,1,1)
plot(y)
ylabel('Response')
grid on

subplot(3,1,2)
plot(e)
ylabel('Innovation') 
grid on

subplot(3,1,3)
plot(sp,'m')
ylabel('State Index')
yticks([1 2])

% Copyright 2019 The MathWorks, Inc.