www.gusucode.com > econ 案例源码程序 matlab代码 > econ/FilterTimeVaryingStateSpaceModelExample.m
%% Filter Time-Varying State-Space Model % This example shows how to generate data from a known model, fit a % state-space model to the data, and then filter the states. %% % Suppose that a latent process comprises an AR(2) and an MA(1) model. % There are 50 periods, and the MA(1) process drops out of the model for % the final 25 periods. Subsequently, the state equation for the first 25 % periods is % % $$\begin{array}{l} % {x_{1,t}} = 0.7{x_{1,t - 1}} - 0.2{x_{1,t - 2}} + {u_{1,t}}\\ % {x_{2,t}} = {u_{2,t}} + 0.6{u_{2,t - 1}}, % \end{array}$$ % % and for the last 25 periods, it is % % $${x_{1,t}} = 0.7{x_{1,t - 1}} - 0.2{x_{1,t - 2}} + {u_{1,t}},$$ % % where $u_{1,t}$ and $u_{2,t}$ are Gaussian with mean 0 and standard % deviation 1. %% % Assuming that the series starts at 1.5 and 1, respectively, generate a % random series of 50 observations from $x_{1,t}$ and $x_{2,t}$. % Copyright 2015 The MathWorks, Inc. T = 50; ARMdl = arima('AR',{0.7,-0.2},'Constant',0,'Variance',1); MAMdl = arima('MA',0.6,'Constant',0,'Variance',1); x0 = [1.5 1; 1.5 1]; rng(1); x = [simulate(ARMdl,T,'Y0',x0(:,1)),... [simulate(MAMdl,T/2,'Y0',x0(:,2));nan(T/2,1)]]; %% % The last 25 values for the simulated MA(1) data are |NaN| values. %% % Suppose further that the latent processes are measured using % % $${y_t} = 2\left( {{x_{1,t}} + {x_{2,t}}} \right) + {\varepsilon _t},$$ % % for the first 25 periods, and % % $${y_t} = 2{x_{1,t}} + {\varepsilon _t}$$ % % for the last 25 periods, where $\varepsilon_t$ is Gaussian with mean 0 % and standard deviation 1. % % Use the random latent state process (|x|) and the observation equation to % generate observations. y = 2*nansum(x')'+randn(T,1); %% % Together, the latent process and observation equations compose a % state-space model. Supposing that the coefficients are unknown % parameters, the state-space model is % % $$\begin{array}{l} % \left[ {\begin{array}{*{20}{c}} % {{x_{1,t}}}\\ % {{x_{2,t}}}\\ % {{x_{3,t}}}\\ % {{x_{4,t}}} % \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} % {{\phi _1}}&{{\phi _2}}&0&0\\ % 1&0&0&0\\ % 0&0&0&{{\theta _1}}\\ % 0&0&0&0 % \end{array}} \right]\left[ {\begin{array}{*{20}{c}} % {{x_{1,t - 1}}}\\ % {{x_{2,t - 1}}}\\ % {{x_{3,t - 1}}}\\ % {{x_{4,t - 1}}} % \end{array}} \right] + \left[ {\begin{array}{*{20}{c}} % 1&0\\ % 0&0\\ % 0&1\\ % 0&1 % \end{array}} \right]\left[ {\begin{array}{*{20}{c}} % {{u_{1,t}}}\\ % {{u_{2,t}}} % \end{array}} \right]\\ % {y_t} = a({x_{1,t}} + {x_{3,t}}) + {\varepsilon _t} % \end{array}$$ % % for the first 25 periods, % % $$\begin{array}{l} % \left[ {\begin{array}{*{20}{c}} % {{x_{1,t}}}\\ % {{x_{2,t}}} % \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} % {{\phi _1}}&{{\phi _2}}&0&0\\ % 1&0&0&0 % \end{array}} \right]\left[ {\begin{array}{*{20}{c}} % {{x_{1,t - 1}}}\\ % {{x_{2,t - 1}}}\\ % {{x_{3,t - 1}}}\\ % {{x_{4,t - 1}}} % \end{array}} \right] + \left[ {\begin{array}{*{20}{c}} % 1\\ % 0 % \end{array}} \right]{u_{1,t}}\\ % {y_t} = b{x_{1,t}} + {\varepsilon _t} % \end{array}$$ % % for period 26, and % % $$\begin{array}{l} % \left[ {\begin{array}{*{20}{c}} % {{x_{1,t}}}\\ % {{x_{2,t}}} % \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} % {{\phi _1}}&{{\phi _2}}\\ % 1&0 % \end{array}} \right]\left[ {\begin{array}{*{20}{c}} % {{x_{1,t - 1}}}\\ % {{x_{2,t - 1}}} % \end{array}} \right] + \left[ {\begin{array}{*{20}{c}} % 1\\ % 0 % \end{array}} \right]{u_{1,t}}\\ % {y_t} = b{x_{1,t}} + {\varepsilon _t} % \end{array}$$ % % for the last 24 periods. %% % Write a function that specifies how the parameters in |params| map to the % state-space model matrices, the initial state values, and the type of % state. % % <include>AR2MAParamMap.m</include> % %% % Save this code as a file named |AR2MAParamMap| on your MATLAB(R) path. %% % Create the state-space model by passing the function |AR2MAParamMap| as a % function handle to |ssm|. Mdl = ssm(@(params)AR2MAParamMap(params,T)); %% % |ssm| implicitly creates the state-space model. Usually, you % cannot verify an implicitly defined state-space model. %% % Pass the observed responses (|y|) to |estimate| to estimate the % parameters. Specify an arbitrary set of positive initial values for the % unknown parameters. params0 = 0.1*ones(5,1); EstMdl = estimate(Mdl,y,params0); %% % |EstMdl| is an |ssm| model containing the estimated coefficients. % Likelihood surfaces of state-space models might contain local maxima. % Therefore, it is good practice to try several initial parameter values, % or consider using |refine|. %% % Filter the states and obtain state forecasts by passing |EstMdl| and the % observed responses to |filter|. [~,~,Output]= filter(EstMdl,y); %% % |Output| is a |T|-by-1 structure array containing the filtered states and % state forecasts, among other things. %% % Extract the filtered and forecasted states from the cell arrays. Recall % that the two, different states are in positions 1 and 3. The states in % positions 2 and 4 help specify the processes of interest. stateIndx = [1 3]; % State indices of interest FilteredStates = NaN(T,numel(stateIndx)); ForecastedStates = NaN(T,numel(stateIndx)); for t = 1:T maxInd = size(Output(t).FilteredStates,1); mask = stateIndx <= maxInd; FilteredStates(t,mask) = Output(t).FilteredStates(stateIndx(mask),1); ForecastedStates(t,mask) = Output(t).ForecastedStates(stateIndx(mask),1); end %% % Plot the true state values, the filtered states, and the state forecasts % together for each model. figure plot(1:T,x(:,1),'-k',1:T,FilteredStates(:,1),':r',... 1:T,ForecastedStates(:,1),'--g','LineWidth',2); title('AR(2) State Values') xlabel('Period') ylabel('State Value') legend({'True state values','Filtered state values','State forecasts'}); figure plot(1:T,x(:,2),'-k',1:T,FilteredStates(:,2),':r',... 1:T,ForecastedStates(:,2),'--g','LineWidth',2); title('MA(1) State Values') xlabel('Period') ylabel('State Value') legend({'True state values','Filtered state values','State forecasts'});