Этот пример показывает, как предсказать мультипликативную сезонную модель ARIMA с помощью forecast
. Ряд ответа является ежемесячными международными числами авиапассажира от 1 949 до 1960.
Загрузите наборы данных рецессий и авиакомпания. Преобразуйте ответ.
load(fullfile(matlabroot,'examples','econ','Data_Airline.mat')) load Data_Recessions y = log(Data);
Создайте предиктор (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
Задайте индексные наборы тот раздел данные в выборки прогноза и оценку.
nSim = 60; % Forecast period
T = length(y);
estInds = 1:(T-nSim);
foreInds = (T-nSim+1):T;
Оцените модель регрессии с мультипликативным сезонным ARIMA} ошибки:
Установите прерывание модели регрессии на 0, поскольку это не идентифицируется в модели с интегрированными ошибками.
Mdl = regARIMA('D',1,'Seasonality',12,'MALags',1,'SMALags',12,... 'Intercept',0); EstMdl = estimate(Mdl,y(estInds),'X',X(estInds));
Regression with ARIMA(0,1,1) Error Model Seasonally Integrated with Seasonal MA(12) (Gaussian Distribution): Value StandardError TStatistic PValue _________ _____________ __________ __________ Intercept 0 0 NaN NaN MA{1} -0.35662 0.10393 -3.4312 0.0006009 SMA{12} -0.67729 0.11294 -5.9972 2.0078e-09 Beta(1) 0.0015098 0.020533 0.073531 0.94138 Variance 0.0015198 0.00021411 7.0983 1.2632e-12
Используйте предполагаемые коэффициенты модели (содержавшийся в EstMdl
), чтобы сгенерировать прогнозы MMSE и соответствующие среднеквадратичные погрешности по 60-месячному горизонту. Используйте наблюдаемый ряд в качестве преддемонстрационных данных. По умолчанию forecast
выводит преддемонстрационные инновации и безусловные воздействия с помощью заданной модели и наблюдений.
[YF,YMSE] = forecast(EstMdl,nSim,'X0',X(estInds),... 'Y0',y(estInds),'XF',X(foreInds)); ForecastInt = [YF,YF] + 1.96*[-sqrt(YMSE), sqrt(YMSE)]; figure h1 = plot(dates,y); title('{\bf Forecasted Monthly Passenger Totals}') hold on h2 = plot(dates(foreInds),YF,'Color','r','LineWidth',2); h3 = plot(dates(foreInds),ForecastInt,'k--','LineWidth',2); datetick legend([h1,h2,h3(1)],'Observations','MMSE Forecasts',... '95% MMSE Forecast Intervals','Location','NorthWest') axis tight hold off
Модель регрессии с ошибками SMA, кажется, предсказывает ряд хорошо, хотя немного переоценивая. Поскольку ошибочный процесс является неустановившимся, интервалы прогноза расширяются, когда время увеличивается.
Сравните прогнозы MMSE с прогнозами Монте-Карло путем симуляции 500 демонстрационных путей от EstMdl
по горизонту прогноза.
[e0,u0] = infer(EstMdl,y(estInds),'X',X(estInds)); rng(5); numPaths = 500; ySim = simulate(EstMdl,nSim,'numPaths',numPaths,... 'E0',e0,'U0',u0,'X',X(foreInds)); meanYSim = mean(ySim,2); ForecastIntMC = [prctile(ySim,2.5,2),prctile(ySim,97.5,2)]; figure h1 = plot(dates(foreInds),y(foreInds)); title('{\bf Forecasted Monthly Passenger Totals}') hold on h2 = plot(dates(foreInds),YF,'Color',[0.85,0.85,0.85],... 'LineWidth',4); h3 = plot(dates(foreInds),ForecastInt,'--','Color',... [0.85,0.85,0.85],'LineWidth',4); h4 = plot(dates(foreInds),meanYSim,'k','LineWidth',2); h5 = plot(dates(foreInds),ForecastIntMC,'k--','LineWidth',2); datetick legend([h1,h2,h3(1),h4,h5(1)],'Observations',... 'MMSE Forecasts','95% MMSE Forecast Intervals',... 'Monte Carlo Forecasts','95% Monte Carlo Forecast Intervals',... 'Location','NorthWest') axis tight hold off
Прогнозы MMSE и средние прогнозы Монте-Карло фактически неразличимы. Однако существуют небольшие несоответствия между теоретическими 95%-ми интервалами прогноза и основанными на симуляции 95%-ми интервалами прогноза.
estimate
| forecast
| regARIMA
| simulate