%% Choose State-Space Model Specification Using Backtesting
% This example shows how to choose the state-space model specification with
% the best predictive performance using a rolling window.  A rolling window
% analysis for an explicitly defined state-space model is straightforward,
% so this example focuses on implicitly defined state-space models.
%%
% Suppose that the linear relationship between the change in the observed
% unemployment rate and the nominal gross national product (nGNP) growth
% rate is of interest. Suppose further that you want to choose between an
% AR(1) or an AR(2) model for the first difference of the unemployment rate
% (i.e., the state). That is,
%
% $$\begin{array}{*{20}{l}}
% {{\rm{Model 1:}}\begin{array}{*{20}{c}}
% {{x_t} = {\phi _{11}}{x_{t - 1}} + {\phi _{12}}{x_{t - 1}} + {\varepsilon _t}}\\
% {{y_t} - {\beta _1}{z_t} = {x_t} + {u_t}}
% \end{array}}\\
% {{\rm{Model 2:}}\begin{array}{*{20}{c}}
% {{x_t} = {\phi _{21}}{x_{t - 1}} + {\varepsilon _t}}\\ 
% {{y_t} - {\beta _2}{z_t} = {x_t} + {u_t}}
% \end{array}}
% \end{array}.$$
%
% * $\varepsilon_t$ and $u_t$ are Gaussian process with mean 0 and variance 1.
% * $x_t$ is the true unemployment rate at time $t$.
% * $y_t$ is the observed unemployment.
% * $z_t$ is the nGNP rate.
%
%%
% Declare the functions |rwParamMap.m| and |rwAR2ParamMap.m|, which specify
% how the parameters in |params| map to the state-space model matrices, the
% initial state values, and the type of state. Save them in your working
% folder.
%
% <include>rwParamMap.m</include>
%
% <include>rwAR2ParamMap.m</include>
%
%% 
% Load the Nelson-Plosser data set, which contains the unemployment rate
% and nGNP series, among other things.
load Data_NelsonPlosser
%%
% Preprocess the data by taking the natural logarithm of the nGNP series,
% and the first difference of each.   Also, remove the starting |NaN|
% values from each series.
isNaN = any(ismissing(DataTable),2); % Flag periods containing NaNs
gnpn = DataTable.GNPN(~isNaN);
ur = DataTable.UR(~isNaN);           % Sample size
Z = diff(log(gnpn));
y = diff(ur);
T = size(y,1);
%%
% If you define a state-space model implicitly and the response and
% predictor data (i.e., |y| and |Z|) exist in the MATLAB(R) Workspace, then
% the software creates a link from the parameter-to-matrix mapping function
% those series. If the data do not exist in the MATLAB(R) Workspace, then
% the software creates the model, but you must provide the data using the
% appropriate name-value pair arguments when you, e.g., estimate the
% model.
%%
% Therefore, to conduct a rolling window analysis when the state-space
% model is implicitly defined and there is a regression component, you must
% specify the state-space model indicating the indices of the data to be
% analyzed for each window. Conduct a rolling window analysis of the
% simulated data. Let the rolling window length (|m|) be 40 periods and the
% forecast horizon (|h|) be 10 periods. For this example, assume that the
% time-series are stable (i.e., all parameters are time-invariant).
m = 40; 
N = T - m + 1;      % Number of rolling windows
h = 10;

fError1 = nan(N,h); % Preallocation
fError2 = nan(N,h); 

for j = 1:N;
    idxRW = j:(m + j - h - 1);
		idxFH = (m + j - h):(m + j - 1);
    Mdl1 = ssm(@(c)rwParamMap(c,y(idxRW),Z(idxRW)));
    Mdl2 = ssm(@(c)rwAR2ParamMap(c,y(idxRW),Z(idxRW)));
    [EstMdl1,estParams1] = estimate(Mdl1,y(idxRW),[0.5 -20]',...
        'Display','off');
    [EstMdl2,estParams2] = estimate(Mdl2,y(idxRW),[0.5 0.1 -20]',...
       'Display','off');
    fY1 = forecast(EstMdl1,h,y(idxRW),'Predictors0',Z(idxRW),...
        'PredictorsF',Z(idxFH),'Beta',estParams1(end));
    fY2 = forecast(EstMdl2,h,y(idxRW),'Predictors0',Z(idxRW),...
        'PredictorsF',Z(idxFH),'Beta',estParams2(end));
		fError1(j,:) = y(idxFH) - fY1;
		fError2(j,:) = y(idxFH) - fY2;
end
%%
% Compute the RMSE for each step-ahead forecast, and compare them for each
% model.
fRMSE1 = sqrt(mean(fError1.^2));
fRMSE2 = sqrt(mean(fError2.^2));
fRMSE1 < fRMSE2
%%
% Overall, the predictive performance of the AR(2) model is better than the
% AR(1) model.
%%
% Alternatively, you can compare the predictive performance of the models
% using the Diebold-Mariano test. For more details, see
% <docid:econ_ug.buhn26v [1]>.


%% 
% Copyright 2012 The MathWorks, Inc.