Шаг 1. Загрузите данные.

Загрузите набор данных Data_Airline, и постройте естественный журнал ежемесячных подсчётов количества пассажиров.

load('Data_Airline.mat')
dat = log(Data); % Transform to logarithmic scale 
T = size(dat,1);
y = dat(1:103);  % Estimation sample

y является частью dat используется для оценки, а остальная часть dat - выборка удержания для сравнения прогнозов двух моделей.

Шаг 2. Определите и подгоните модель, задающую сезонные лаги.

Создайте ARIMA$$(0,1,1)(0,1,1{)}_{12}$$ модель

где $${\epsilon }_t$$ является независимым и идентично распределенным нормально распределенным рядом со средним 0 и отклонением $${\sigma }^2$$. Использование estimate для подгонки Mdl1 на y.

Mdl1 = arima('Constant',0,'MALags',1,'D',1,...
    'SMALags',12,'Seasonality',12);
EstMdl1 = estimate(Mdl1,y);

 
    ARIMA(0,1,1) Model Seasonally Integrated with Seasonal MA(12) (Gaussian Distribution):
 
                 Value      StandardError    TStatistic      PValue  
                ________    _____________    __________    __________

    Constant           0              0           NaN             NaN
    MA{1}       -0.35732       0.088031        -4.059      4.9286e-05
    SMA{12}     -0.61469       0.096249       -6.3864      1.6985e-10
    Variance    0.001305      0.0001527        8.5467      1.2671e-17

Подобранная модель

где $${\epsilon }_t$$ iid обычно распределенный ряд со средним значением 0 и отклонением 0,0013.

Шаг 3. Определите и подгоните модель с помощью сезонных манекенов.

Создайте модель ARIMAX (0,1,1) с сезонным дифференцированием периода 12 и регрессионым компонентом ,

$$\{x_t;t=1,...,T\}$$ - серия T- векторов-столбцов, имеющая длину 12, которая указывает, в каком месяце наблюдение $$t$$ измеряли. A 1 в строке i $$x_t$$ указывает, что наблюдение было измерено в месяце i, остальные элементы - 0 с.

Обратите внимание, что если вы включаете аддитивную константу в модель, то строки T матрицы проекта X состоят из векторов-строк $$\left[\begin{array}{cc}1& x_t^{\prime }\end{array}\right]$$. Поэтому X является неполный ранг, и один коэффициент регрессии не идентифицируется. Константа остается вне этого примера, чтобы избежать отвлечения от основной цели. Форматируйте матрицу X в выборке

X  = dummyvar(repmat((1:12)',12,1));

% Format the presample X matrix
X0 = [zeros(1,11) 1; dummyvar((1:12)')];
Mdl2 = arima('Constant',0,'MALags',1,'D',1,...
    'Seasonality',12);
EstMdl2   = estimate(Mdl2,y,'X',[X0; X]);

 
    ARIMAX(0,1,1) Model Seasonally Integrated (Gaussian Distribution):
 
                  Value       StandardError    TStatistic      PValue  
                __________    _____________    __________    __________

    Constant             0              0            NaN            NaN
    MA{1}         -0.40711       0.084387        -4.8242     1.4053e-06
    Beta(1)      -0.002577       0.025168       -0.10239        0.91845
    Beta(2)     -0.0057769       0.031885       -0.18118        0.85623
    Beta(3)     -0.0022034       0.030527      -0.072179        0.94246
    Beta(4)     0.00094737       0.019867       0.047686        0.96197
    Beta(5)     -0.0012146       0.017981      -0.067551        0.94614
    Beta(6)        0.00487       0.018374        0.26505        0.79097
    Beta(7)     -0.0087944       0.015285       -0.57535        0.56505
    Beta(8)      0.0048346       0.012484        0.38728        0.69855
    Beta(9)       0.001437       0.018245       0.078758        0.93722
    Beta(10)      0.009274       0.014751        0.62869        0.52955
    Beta(11)     0.0073665         0.0105        0.70158        0.48294
    Beta(12)    0.00098841       0.014295       0.069146        0.94487
    Variance     0.0017715     0.00024657         7.1848     6.7329e-13

Подобранная модель

где $${\epsilon }_t$$ iid обычно распределенный ряд со средним 0 и отклонением 0.0017 и $$\underset{}{\overset{ˆ}{\beta }}$$ - вектор-столбец со значениями Beta1 - Beta12. Обратите внимание, что оценки MA{1} и Variance между Mdl1 и Mdl2 не равны.

Шаг 4. Прогноз с использованием обеих моделей.

Использование forecast прогнозировать обе модели 41 периоды в будущем с июля 1957 года. Постройте график выборки удержания с помощью этих прогнозов.

yF1 = forecast(EstMdl1,41,y);
yF2 = forecast(EstMdl2,41,y,'X0',X(1:103,:),'XF',X(104:end,:));
l1 = plot(100:T,dat(100:end),'k','LineWidth',3);
hold on
l2 = plot(104:144,yF1,'-r','LineWidth',2);
l3 = plot(104:144,yF2,'-b','LineWidth',2);
hold off
title('Passenger Data: Actual vs. Forecasts')
xlabel('Month')
ylabel('Logarithm of Monthly Passenger Data')
legend({'Observations','Polynomial Forecast',...
 			'Regression Forecast'},'Location','NorthWest')

Хотя они переоценивают наблюдения удержания, прогнозы обеих моделей почти эквивалентны. Одно из основных различий моделей заключается в том, что EstMdl1 является более парсимонным, чем EstMdl2.

Ссылки:

Box, G. E. P., G. M. Jenkins, and G. C. Reinsel. Анализ временных рядов: прогнозирование и управление. 3-й эд. Englewood Cliffs, Нью-Джерси: Prentice Hall, 1994.