www.gusucode.com > econ 案例源码程序 matlab代码 > econ/ForecastObservationsOfKnownTimeInvariantDiffuseStaExample.m

    %% Forecast Observations of Time-Invariant Diffuse State-Space Model
%%
% Suppose that a latent process is a random walk.  Subsequently, the state
% equation is 
%
% $$x_t = x_{t-1} + u_t,$$
%
% where $u_t$ is Gaussian with mean 0 and standard deviation 1.
%%
% Generate a random series of 100 observations from $x_t$, assuming that the
% series starts at 1.5.

% Copyright 2015 The MathWorks, Inc.

T = 100;
x0 = 1.5;
rng(1); % For reproducibility
u = randn(T,1);
x = cumsum([x0;u]);
x = x(2:end);
%%
% Suppose further that the latent process is subject to additive
% measurement error.  Subsequently, the observation equation is
%
% $$y_t = x_t + \varepsilon_t,$$
%
% where $\varepsilon_t$ is Gaussian with mean 0 and standard deviation
% 0.75.  Together, the latent process and observation equations compose a
% state-space model. 

%%
% Use the random latent state process (|x|) and the observation equation to
% generate observations.
y = x + 0.75*randn(T,1);
%%
% Specify the four coefficient matrices.
A = 1;
B = 1;
C = 1;
D = 0.75;
%%
% Create the diffuse state-space model using the coefficient matrices.
% Specify that the inital state distribution is diffuse.
Mdl = dssm(A,B,C,D,'StateType',2)
%%
% |Mdl| is an |dssm| model.  Verify that the model is correctly specified
% using the display in the Command Window.
%%
% Forecast observations 10 periods into the future, and estimate the mean
% squared errors of the forecasts.
numPeriods = 10;
[ForecastedY,YMSE] = forecast(Mdl,numPeriods,y);
%%
% Plot the forecasts with the in-sample responses, and 95% Wald-type
% forecast intervals.
ForecastIntervals(:,1) = ForecastedY - 1.96*sqrt(YMSE);
ForecastIntervals(:,2) = ForecastedY + 1.96*sqrt(YMSE);

figure
plot(T-20:T,y(T-20:T),'-k',T+1:T+numPeriods,ForecastedY,'-.r',...
    T+1:T+numPeriods,ForecastIntervals,'-.b',...
    T:T+1,[y(end)*ones(3,1),[ForecastedY(1);ForecastIntervals(1,:)']],':k',...
    'LineWidth',2)
hold on
title({'Observed Responses and Their Forecasts'})
xlabel('Period')
ylabel('Responses')
legend({'Observations','Forecasted observations','95% forecast intervals'},...
    'Location','Best')
hold off
%%
% The forecast intervals flare out because the process is nonstationary.