www.gusucode.com > econ 案例源码程序 matlab代码 > econ/EstimateaRegressionModelwithMultiplicativeARIMAErrorsExample.m
%% Estimate a Regression Model with Multiplicative ARIMA Errors % This example shows how to fit a regression model with multiplicative ARIMA % errors to data using |estimate|. % Copyright 2015 The MathWorks, Inc. %% % Load the Airline data set from the MATLAB(R) root folder, and load the % Recession data set. Plot the monthly passenger totals and the log of the % totals. load(fullfile(matlabroot,'examples','econ','Data_Airline.mat')) load Data_Recessions y = Data; logY = log(y); figure subplot(2,1,1) plot(y) title('{\bf Monthly Passenger Totals (Jan1949 - Dec1960)}') datetick subplot(2,1,2) plot(log(y)) title('{\bf Monthly Passenger Log-Totals (Jan1949 - Dec1960)}') datetick %% % The log transformation seems to linearize the time series. %% % Construct the predictor (|X|), which is whether the country was in a recession % during the sampled period. A 0 in row _t_ means the country was not in % a recession in month _t_, and a 1 in row _t_ means that it was in a recession % in month _t_. X = zeros(numel(dates),1); % Preallocation for j = 1:size(Recessions,1) X(dates >= Recessions(j,1) & dates <= Recessions(j,2)) = 1; end %% % Fit the simple linear regression model % % $${y_t} = c + {X_t}\beta + {u_t}$$ % % to the data. Fit = fitlm(X,logY); %% % |Fit| is a |LinearModel| that contains the least squares estimates. %% % Check for standard linear model assumption departures by plotting the % residuals several ways. figure subplot(2,2,1) plotResiduals(Fit,'caseorder','ResidualType','Standardized',... 'LineStyle','-','MarkerSize',0.5) subplot(2,2,2) plotResiduals(Fit,'lagged','ResidualType','Standardized') subplot(2,2,3) plotResiduals(Fit,'probability','ResidualType','Standardized') subplot(2,2,4) plotResiduals(Fit,'histogram','ResidualType','Standardized') r = Fit.Residuals.Standardized; figure subplot(2,1,1) autocorr(r) subplot(2,1,2) parcorr(r) %% % The residual plots indicate that the unconditional disturbances are autocorrelated. % The probability plot and histogram seem to indicate that the unconditional % disturbances are Gaussian. %% % The ACF of the residuals confirms that the unconditional disturbances % are autocorrelated. %% % Take the 1st difference of the residuals and plot the ACF and PACF of % the differenced residuals. dR = diff(r); figure subplot(2,1,1) autocorr(dR,50) subplot(2,1,2) parcorr(dR,50) %% % The ACF shows that there are significantly large autocorrelations, particularly % at every 12th lag. This indicates that the unconditional disturbances % have 12th degree seasonal integration. %% % Take the first and 12th differences of the residuals. Plot the differenced % residuals, and their ACF and PACF. DiffPoly = LagOp([1 -1]); SDiffPoly = LagOp([1 -1],'Lags',[0, 12]); diffR = filter(DiffPoly*SDiffPoly,r); figure subplot(2,1,1) plot(diffR) axis tight subplot(2,2,3) autocorr(diffR) axis tight subplot(2,2,4) parcorr(diffR) axis tight %% % The residuals resemble white noise (with possible heteroscedasticity). % According to Box and Jenkins (1994), Chapter 9, the ACF and PACF indicate % that the unconditional disturbances are an % $\rm{IMA}(0,1,1)\times(0,1,1)_{12}$ model. %% % Specify the regression model with $\rm{IMA}(0,1,1)\times(0,1,1)_{12}$ % errors: % % $$\begin{array}{c}{y_t} = {X_t}\beta + {u_t}\\\left( {1 - L} \right)\left( {1 - {L^{12}}} \right){u_t} = \left( {1 + {b_1}L} \right)\left( {1 + {B_{12}}{L^{12}}} \right){\varepsilon _t}.\end{array}$$ % Mdl = regARIMA('MALags',1,'D',1,'Seasonality',12,'SMALags',12) %% % Partition the data set into the presample and estimation sample so that % you can initialize the series. |P| = |Q| = 13, so the presample should % be at least 13 periods long. preLogY = logY(1:13); % Presample responses estLogY = logY(14:end); % Estimation sample responses preX = X(1:13); % Presample predictors estX = X(14:end); % Estimation sample predictors %% % Obtain presample unconditional disturbances from a linear regression of % the presample data. PreFit = fitlm(preX,preLogY);... % Presample fit for presample residuals EstFit = fitlm(estX,estLogY);... % Estimation sample fit for the intercept U0 = PreFit.Residuals.Raw; %% % If the error model is integrated, then the regression model intercept % is not identifiable. Set |Intercept| to the estimated intercept from a % linear regression of the estimation sample data. Estimate the regression % model with IMA errors. Mdl.Intercept = EstFit.Coefficients{1,1}; EstMdl = estimate(Mdl,estLogY,'X',estX,'U0',U0); %% % |MA{1}| and |Beta1| are not significantly different from 0. You can remove % these parameters from the model, possibly add other parameters (e.g., % AR parameters), and compare multiple model fits using |aicbic|. Note that % the estimation and presample should be the same over competing models. %% % References: % % Box, G. E. P., G. M. Jenkins, and G. C. Reinsel. _Time Series Analysis: Forecasting and Control_. 3rd ed. Englewood Cliffs, NJ: Prentice Hall, 1994.