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

    %% Extract Other Estimates from |Output|
% Estimate a diffuse state-space model, smooth the states, and then
% extract other estimates from the |Output| output argument.
%%
% Consider the diffuse state-space model 
%
% $$\begin{array}{l}
% \left[ {\begin{array}{*{20}{c}}
% {{x_{1,t}}}\\
% {{x_{2,t}}}
% \end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
% \phi &0\\
% 0&1
% \end{array}} \right]\left[ {\begin{array}{*{20}{c}}
% {{x_{1,t - 1}}}\\
% {{x_{2,t - 1}}}
% \end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
% {{\sigma _1}}&0\\
% 0&{{\sigma _2}}
% \end{array}} \right]\left[ {\begin{array}{*{20}{c}}
% {{u_{1,t}}}\\
% {{u_{2,t}}}
% \end{array}} \right]\\
% y_t = \left[ {\begin{array}{*{20}{c}}
% 1&1
% \end{array}} \right]\left[ {\begin{array}{*{20}{c}}
% {{x_{1,t}}}\\
% {{x_{2,t}}}
% \end{array}} \right]
% \end{array}$$
%
% The state variable $x_{1,t}$ is an AR(1) model with autoregressive
% coefficient $\phi$.  $x_{2,t}$ is a random walk. The disturbances
% $u_{1,t}$ and $u_{2,t}$ are independent Gaussian random variables with
% mean 0 and standard deviations $\sigma_1$ and $\sigma_2$, respectively.
% The observation $y_{t}$ is the error-free sum of $x_{1,t}$ and $x_{2,t}$.
%%
% Generate data from the state-space model.  To simulate the data, suppose
% that the sample size $T = 100$, $\phi = 0.6$, $\sigma_1 = 0.2$, $\sigma_2
% = 0.1$, and $x_{1,0} = x_{2,0} = 2$.
rng(1); % For reproducibility
T = 100;
ARMdl = arima('AR',0.6,'Constant',0,'Variance',0.2^2);
x1 = simulate(ARMdl,T,'Y0',2);
u3 = 0.1*randn(T,1);
x3 = cumsum([2;u3]);
x3 = x3(2:end);
y = x1 + x3;
%%
% Specify the coefficient matrices of the state-space model.  To indicate
% unknown parameters, use |NaN| values.
A = [NaN 0; 0 1];
B = [NaN 0; 0 NaN];
C = [1 1];
%%
% Create a diffuse state-space model that describes the model above.
% Specify that $x_{1,t}$ and $x_{2,t}$ have diffuse initial state
% distributions.
StateType = [2 2];
Mdl = dssm(A,B,C,'StateType',StateType);
%%
% Estimate the unknown parameters of |Mdl|. Choose initial parameter
% values for optimization.  Specify that the standard deviations are
% constrained to be positive, but all other parameters are unconstrained
% using the |'lb'| name-value pair argument.
params0 = [0.01 0.1 0.01]; % Initial values chosen arbitrarily
EstMdl = estimate(Mdl,y,params0,'lb',[-Inf 0 0]);
%%
% The parameters are close to their true values.
%%
% Smooth the states of |EstMdl|, and request all other available output.
[X,logL,Output] = smooth(EstMdl,y);
%%
% |X| is a |T|-by-2 matrix of smoothed states, |logL| is the final 
% optimized log-likelihood value, and |Output| is a structure array
% containing various estimates that the Kalman filter requires.  List the
% fields of |output| using |fields|.
fields(Output)
%%
% Convert |Output| to a table.
OutputTbl = struct2table(Output);
OutputTbl(1:10,1:4) % Display first ten rows of first four variables
%%
% The first two rows of the table contain empty cells or zeros.  These
% correspond to the observations required to initialize the diffuse Kalman
% filter.  That is, |SwitchTime| is 2.
SwitchTime = 2;
%%
% Plot the smoothed states and their individual 95% Wald-type confidence
% intervals.
CI = nan(T,2,2);

for j = (SwitchTime + 1):T
    CovX = OutputTbl.SmoothedStatesCov{j};
    CI(j,:,1) = X(j,1) + 1.96*sqrt(CovX(1,1))*[-1 1];
    CI(j,:,2) = X(j,2) + 1.96*sqrt(CovX(2,2))*[-1 1];    
end

figure;
plot(1:T,X(:,1),'k',1:T,CI(:,:,1),'--r');
xlabel('Period');
ylabel('Smoothed states');
title('State 1 Estimates')
legend('Smoothed','95% Individual CIs');
grid on;

figure;
plot(1:T,X(:,2),'k',1:T,CI(:,:,2),'--r');
xlabel('Period');
ylabel('Smoothed states');
title('State 2 Estimates')
legend('Smoothed','95% Individual CIs');
grid on;