www.gusucode.com > econ 案例源码程序 matlab代码 > econ/RefineEstimationOfAStateSpaceModelWithARegressionComponeExample.m
%% Refine Estimation of State-Space Model Containing Regression Component % This example shows how to choose a better set of initial parameter values % to pass to |estimate| when you estimate |ssm| models containing a % regression component. %% % State-space models tend to induce likelihood surfaces that have multiple % maxima. Also, the software uses numerical optimization to find the % parameter estimates. Therefore, it is important to use a good initial % value for the parameters to attain the global maximum of the likelihood % function. %% % Suppose that the relationship between the change in the unemployment rate % and the nominal gross national product (nGNP) growth rate is of interest. Suppose % further that the first difference of the unemployment rate is an % ARMA(1,1) series. Symbolically, and in state-space form, the model is % % $$\begin{array}{l} % \left[ {\begin{array}{*{20}{c}} % {{x_{1,t}}}\\ % {{x_{2,t}}} % \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} % \phi &\theta \\ % 0&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\\ % 1&0 % \end{array}} \right]\left[ {\begin{array}{*{20}{c}} % {{u_{1,t}}}\\ % {{u_{2,t}}} % \end{array}} \right]\\ % {y_t} - \beta {Z_t} = {x_{1,t}} + \sigma\varepsilon_t, % \end{array}$$ % % where: % % * $x_{1,t}$ is the change in the unemployment rate at time _t_. % * $x_{2,t}$ is a dummy state for the MA(1) effect. % * $y_{t}$ is the observed change in the unemployment rate being deflated by the % growth rate of nGNP ($Z_t$). % * $u_{1,t}$ is the Gaussian series of state disturbances having mean 0 and % standard deviation 1. % * $\varepsilon_t$ is the Gaussian series of observation innovations having % mean 0 and standard deviation $\sigma$. % %% % Load the Nelson-Plosser data set, which contains the unemployment rate and % nGNP series data. load Data_NelsonPlosser rng(1); % For reproducibility %% % 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); u = DataTable.UR(~isNaN); T = size(gnpn,1); % Sample size Z = [ones(T-1,1) diff(log(gnpn))]; y = diff(u); %% % This example continues excluding |NaN| values from the analysis. % However, using the Kalman filter framework, the software can accommodate % series containing missing values. %% % Specify the coefficient matrices. A = [NaN NaN; 0 0]; B = [1 0;1 0]; C = [1 0]; D = NaN; %% % Specify the state-space model using |ssm|. Mdl = ssm(A,B,C,D); %% % To help you choose a good set of initial parameter values for estimation, % pass |Mdl| and an initial set of parameter values to |refine|. params0 = [0 0 0]; % Chosen arbitrarily Beta0 = [0 0]; Output = refine(Mdl,y,params0,'Predictors',Z,'Beta0',Beta0); %% % |Output| is a 1-by-5 structure array containing the recommended initial % parameter values. %% % Choose the initial parameter values corresponding to the % largest loglikelihood. logL = cell2mat({Output.LogLikelihood})'; [~,maxLogLIndx] = max(logL) refinedParams0 = Output(maxLogLIndx).Parameters Description = Output(maxLogLIndx).Description %% % The algorithm that yields the highest loglikelihood value is |Starting % value perturbation|, which is the fourth struct in the structure array % |Output|. %% % Estimate |Mdl| using the refined initial parameter values % |refinedParams0|. pBeta = numel(Beta0); EstMdl = estimate(Mdl,y,refinedParams0(1:(end - pBeta)),'Predictors',Z,... 'Beta0',refinedParams0((end - pBeta + 1):end),... 'lb',[-Inf,-Inf,0,-Inf,-Inf]); %% % |estimate| returns reasonable parameter estimates and their corresponding % standard errors.