Оцените модель регрессии с мультипликативными ошибками ARIMA

В этом примере показано, как подбирать модель регрессии с мультипликативными ошибками 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

Figure contains 2 axes objects. Axes object 1 with title blank M o n t h l y blank P a s s e n g e r blank T o t a l s blank ( J a n 1 9 4 9 blank - blank D e c 1 9 6 0 ) contains an object of type line. Axes object 2 with title blank M o n t h l y blank P a s s e n g e r blank L o g - T o t a l s blank ( J a n 1 9 4 9 blank - blank D e c 1 9 6 0 ) contains an object of type line.

Логарифмическое преобразование, кажется, линеаризует временные ряды.

Создайте предиктор (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

Подбирайте простую модель линейной регрессии

yt=c+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;

Figure contains 4 axes objects. Axes object 1 with title Case order plot of residuals contains 2 objects of type line. Axes object 2 with title Plot of residuals vs. lagged residuals contains 3 objects of type line. Axes object 3 with title Normal probability plot of residuals contains 2 objects of type line. Axes object 4 with title Histogram of residuals contains an object of type patch.

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;

Figure contains 2 axes objects. Axes object 1 with title Sample Autocorrelation Function contains 4 objects of type stem, line. Axes object 2 with title Sample Partial Autocorrelation Function contains 4 objects of type stem, line.

Остаточные графики показывают, что безусловные воздействия автокоррелируются. График вероятности и гистограмма, кажется, указывают, что безусловные воздействия являются Гауссовыми.

ACF остаточных значений подтверждает, что безусловные воздействия автокоррелируются.

Возьмите 1-е различие остаточных значений и постройте ACF и PACF differenced остаточных значений.

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;

Figure contains 2 axes objects. Axes object 1 with title Sample Autocorrelation Function contains 4 objects of type stem, line. Axes object 2 with title Sample Partial Autocorrelation Function contains 4 objects of type stem, line.

ACF показывает, что существуют значительно большие автокорреляции, особенно в каждой 12-й задержке. Это указывает, что безусловные воздействия имеют 12-ю степень сезонное интегрирование.

Возьмите первые и 12-е различия остаточных значений. Постройте differenced остаточные значения, и их 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

Figure contains 3 axes objects. Axes object 1 contains an object of type line. Axes object 2 with title Sample Autocorrelation Function contains 4 objects of type stem, line. Axes object 3 with title Sample Partial Autocorrelation Function contains 4 objects of type stem, line.

Остаточные значения напоминают белый шум (с возможным heteroscedasticity). Согласно Полю и Дженкинсу (1994), Глава 9, ACF и PACF указывают, что безусловные воздействия IMA(0,1,1)×(0,1,1)12 модель.

Задайте модель регрессии с IMA(0,1,1)×(0,1,1)12 ошибки:

yt=Xtβ+ut(1-L)(1-L12)ut=(1+b1L)(1+B12L12)εt.

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. Обратите внимание на то, что оценка и предварительная выборка должны быть тем же самым по конкурирующим моделям.

Ссылки:

Поле, G. E. P. Г. М. Дженкинс и Г. К. Рейнсель. Анализ Временных Рядов: Прогнозирование и Управление. 3-й редактор Englewood Cliffs, NJ: Prentice Hall, 1994.

Смотрите также

| |

Похожие темы