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

Загрузите набор данных смертей от несчастного случая.

load Data_Accidental 
y = Data;
T = length(y);

figure
plot(y/1000)
h1 = gca;
h1.XLim = [0,T];
h1.XTick = 1:12:T;
h1.XTickLabel = datestr(dates(1:12:T),10);
title 'Monthly Accidental Deaths';
ylabel 'Number of Deaths (in thousands)';
hold on

Данные показывают потенциальный квадратичный тренд и сильный сезонный компонент с периодичностью 12.

Шаг 2: соответствуйте квадратичному тренду

Соответствуйте полиному

к наблюдаемому ряду.

t = (1:T)';
X = [ones(T,1) t t.^2];

b = X\y;
tH = X*b;
 
h2 = plot(tH/1000,'r','LineWidth',2);
legend(h2,'Quadratic Trend Estimate')
hold off

Шаг 3. Детрендируйте исходный ряд.

Вычтите подходящую квадратичную линию из исходных данных.

xt = y - tH;

Шаг 4. Оцените сезонные переменные индикатора

Создайте индикатор (макет) переменные в течение каждого месяца. Первый индикатор равен одному для наблюдений в январе и нуля в противном случае. Второй индикатор равен одному для наблюдений в феврале и нуля в противном случае. В общей сложности 12 переменных индикатора создаются в течение этих 12 месяцев. Регресс детрендированный ряд против сезонных индикаторов.

mo = repmat((1:12)',6,1);
sX = dummyvar(mo);
  
bS = sX\xt;
st = sX*bS;

figure
plot(st/1000)
title 'Parametric Estimate of Seasonal Component (Indicators)';
h3 = gca;
h3.XLim = [0,T];
ylabel 'Number of Deaths (in thousands)';
h3.XTick = 1:12:T;
h3.XTickLabel = datestr(dates(1:12:T),10);

В этой регрессии все 12 сезонных индикаторов включены в матрицу проекта. Чтобы предотвратить коллинеарность, термин прерывания не включен (альтернативно, можно включать 11 индикаторов и термин прерывания).

Шаг 5. Deseasonalize исходный ряд

Вычтите предполагаемый сезонный компонент из исходного ряда.

dt = y - st;

figure
plot(dt/1000)
title 'Monthly Accidental Deaths (Deseasonalized)';
h4 = gca;
h4.XLim = [0,T];
ylabel 'Number of Deaths (in thousands)';
h4.XTick = 1:12:T;
h4.XTickLabel = datestr(dates(1:12:T),10);

Квадратичный тренд намного более ясен с сезонным удаленным компонентом.

Шаг 6. Оцените неправильный компонент

Вычтите тренд и сезонные оценки от исходного ряда. Остаток является оценкой неправильного компонента.

bt = y - tH - st;

figure
plot(bt/1000)
title('Irregular Component')
h5 = gca;
h5.XLim = [0,T];
ylabel 'Number of Deaths (in thousands)';
h5.XTick = 1:12:T;
h5.XTickLabel = datestr(dates(1:12:T),10);

Можно опционально смоделировать неправильный компонент с помощью модели стохастического процесса.

Ссылки:

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