В этом примере показано, как подогнать регрессионную модель с мультипликативными ошибками ARIMA к данным с помощью estimate.
Загрузка наборов данных авиакомпании и рецессии. Постройте график ежемесячных итоговых записей пассажиров и журнал итоговых записей.
load('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
Преобразование журнала, похоже, линеаризует временной ряд.
Построить предиктор (X), которая заключается в том, находилась ли страна в рецессии в течение выборочного периода. 0 в строке t означает, что страна не находилась в рецессии в месяце t, а 1 в строке t означает, что она находилась в рецессии в месяце t.
X = zeros(numel(dates),1); % Preallocation for j = 1:size(Recessions,1) X(dates >= Recessions(j,1) & dates <= Recessions(j,2)) = 1; end
Подгонка простой модели линейной регрессии
Xtβ + ut
к данным.
Fit = fitlm(X,logY);
Fit является LinearModel содержит оценки наименьших квадратов.
Проверьте отклонения предположений стандартной линейной модели, выведя остатки на график несколькими способами.
figure subplot(2,2,1) plotResiduals(Fit,'caseorder','ResidualType','Standardized',... 'LineStyle','-','MarkerSize',0.5) h = gca; h.FontSize = 8; subplot(2,2,2) plotResiduals(Fit,'lagged','ResidualType','Standardized') h = gca; h.FontSize = 8; subplot(2,2,3) plotResiduals(Fit,'probability','ResidualType','Standardized') h = gca; h.YTick = h.YTick(1:2:end); h.YTickLabel = h.YTickLabel(1:2:end,:); h.FontSize = 8; subplot(2,2,4) plotResiduals(Fit,'histogram','ResidualType','Standardized') h = gca; h.FontSize = 8;
r = Fit.Residuals.Standardized; figure subplot(2,1,1) autocorr(r) h = gca; h.FontSize = 9; subplot(2,1,2) parcorr(r) h = gca; h.FontSize = 9;
Остаточные графики показывают, что безусловные возмущения являются автокоррелированными. График вероятности и гистограмма, по-видимому, указывают на то, что безусловные возмущения являются гауссовыми.
ACF остатков подтверждает, что безусловные нарушения являются автокоррелированными.
Возьмите 1-е различие остатков и постройте график ACF и PACF разностных остатков.
dR = diff(r); figure subplot(2,1,1) autocorr(dR,'NumLags',50) h = gca; h.FontSize = 9; subplot(2,1,2) parcorr(dR,'NumLAgs',50) h = gca; h.FontSize = 9;
ACF показывает, что существуют значительно большие автокорреляции, особенно на каждом 12-м лаге. Это указывает на то, что безусловные нарушения имеют сезонную интеграцию 12-й степени.
Возьмем первое и 12-е различия остатков. Постройте график разностных остатков и их ACF и 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) h = gca; h.FontSize = 7; axis tight subplot(2,2,4) parcorr(diffR) h = gca; h.FontSize = 7; axis tight
Остатки напоминают белый шум (с возможной гетероскедастичностью). Согласно Box and Jenkins (1994), глава 9, ACF и PACF указывают, что безусловные нарушения являются моделью IMA (0,1,1) × (0,1,1) 12.
Укажите регрессионную модель с ошибками IMA (0,1,1) × (0,1,1) 12:
1 + B12L12)
Mdl = regARIMA('MALags',1,'D',1,'Seasonality',12,'SMALags',12)
Mdl =
regARIMA with properties:
Description: "ARIMA(0,1,1) Error Model Seasonally Integrated with Seasonal MA(12) (Gaussian Distribution)"
Distribution: Name = "Gaussian"
Intercept: NaN
Beta: [1×0]
P: 13
D: 1
Q: 13
AR: {}
SAR: {}
MA: {NaN} at lag [1]
SMA: {NaN} at lag [12]
Seasonality: 12
Variance: NaN
Разбейте набор данных на предварительную выборку и выборку оценки, чтобы можно было инициализировать серию. P = Q = 13, так что предварительный пример должен быть не менее 13 периодов.
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
Получение безусловных возмущений предварительной выборки из линейной регрессии данных предварительной выборки.
PreFit = fitlm(preX,preLogY);... % Presample fit for presample residuals EstFit = fitlm(estX,estLogY);... % Estimation sample fit for the intercept U0 = PreFit.Residuals.Raw;
Если модель ошибки интегрирована, то пересечение регрессионной модели не идентифицируется. Набор Intercept к оцененному перехвату из линейной регрессии данных выборки оценки. Оценка регрессионной модели с ошибками IMA.
Mdl.Intercept = EstFit.Coefficients{1,1};
EstMdl = estimate(Mdl,estLogY,'X',estX,'U0',U0);
Regression with ARIMA(0,1,1) Error Model Seasonally Integrated with Seasonal MA(12) (Gaussian Distribution):
Value StandardError TStatistic PValue
_________ _____________ __________ __________
Intercept 5.5722 0 Inf 0
MA{1} -0.025366 0.22197 -0.11427 0.90902
SMA{12} -0.80255 0.052705 -15.227 2.3349e-52
Beta(1) 0.0027588 0.10139 0.02721 0.97829
Variance 0.0072463 0.00015974 45.365 0
MA{1} и Beta1 существенно не отличаются от 0. Можно удалить эти параметры из модели, возможно, добавить другие параметры (например, параметры AR) и сравнить несколько вариантов модели с помощью aicbic. Обратите внимание, что оценка и предварительный пример должны быть одинаковыми по сравнению с конкурирующими моделями.
Ссылки:
Бокс, Г. Э. П., Г. М. Дженкинс и Г. К. Рейнсель. Анализ временных рядов: прогнозирование и контроль. 3-й ред. Энглвуд Клиффс, Нью-Джерси: Прентис Холл, 1994.