%% Fit SETAR Model to Simulated Data
% Assess estimation accuracy using simulated data from a known data-generating process (DGP). This example uses arbitrary parameter values.
%
%% Create Model for DGP
% Create a discrete threshold transition at mid-level 1.

ttDGP = threshold(1)

% ttDGP is a threshold object representing the state switching mechanism of the DGP.
% Create the following fully specified self-exciting TAR (SETAR) model for the DGP.
% State 1: y_t =\varepsilon_t.
% State 2: y_t =2+\varepsilon_t.
% \varepsilon_t \sim N\left(0,1\right).
% Specify the submodels by using arima.

mdl1DGP = arima(Constant=0);
mdl2DGP = arima(Constant=2);
mdlDGP = [mdl1DGP mdl2DGP];

% Because the innovations distribution is invariant across states, the tsVAR software ignores the value of the submodel innovations variance (Variance property).
% Create a threshold-switching model for the DGP. Specify the model-wide innovations variance.

MdlDGP = tsVAR(ttDGP,mdlDGP,Covariance=1);

% MdlDGP is a tsVAR object representing the DGP.
%
%% Simulate Response Paths from DGP
% Generate a random response path of length 100 from the DGP. By default, simulate assumes a SETAR model with delay d=1, in other words, the threshold variable is y_{t-1}.

rng(1) % For reproducibiliy
y = simulate(MdlDGP,100);

% y is a 100-by-1 vector of representing the simulated response path.
%
%% Create Model for Estimation
% Create a partially specified threshold-switching model that has the same structure as the data-generating process, but specify the  transition mid-level, submodel coefficients, and model-wide constant as unknown for estimation.

tt = threshold(NaN);
mdl1 = arima('Constant',NaN);
mdl2 = arima('Constant',NaN);
Mdl = tsVAR(tt,[mdl1,mdl2],'Covariance',NaN);

% Mdl is a partially specified tsVAR object representing a template for estimation. NaN-valued elements of the Switch and Submodels properties indicate estimable parameters.
% Mdl is agnostic of the threshold variable; tsVAR object functions enable you to specify threshold variable characteristics or data.
%
%% Create Threshold Transitions Containing Initial Values
% The estimation procedure requires initial values for all estimable threshold transition parameters.
% Fully specify a threshold transition that has the same structure as tt, but set the mid-level to 0.

tt0 = threshold(0);

% tt0 is a fully specified threshold object.
%
%% Estimate Model
% Fit the model to the simulated path. By default, the model is self exciting and the delay of the threshold variable is d=1.

EstMdl = estimate(Mdl,tt0,y)

% EstMdl is a fully specified tsVAR object representing the estimated SETAR model.
% Display an estimation summary of the submodels,.

summarize(EstMdl)

% The estimates are close to their true values.
% Plot the estimated switching mechanism with the threshold data, which is the response data.

figure
ttplot(EstMdl.Switch,'Data',y)

% Copyright 2021 The MathWorks, Inc.