В этом примере показано, как использовать устойчивый сезонный фильтр для deseasonalize временные ряды (использующий аддитивное разложение). Временные ряды являются ежемесячными смертями от несчастного случая в США от 1 973 до 1978 (Броквелл и Дэвис, 2002).
Загрузите набор данных смертей от несчастного случая.
load(fullfile(matlabroot,'examples','econ','Data_Accidental.mat')) 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 (thousands)'; hold on
Данные показывают сильный сезонный компонент с периодичностью 12.
Сглаживайте данные с помощью скользящего среднего значения с 13 терминами. Чтобы предотвратить потерю наблюдения, повторите первые и последние сглаживавшие значения шесть раз. Вычтите сглаживавший ряд из исходного ряда, чтобы детрендировать данные. Добавьте оценку тренда скользящего среднего значения в наблюдаемый график временных рядов.
sW13 = [1/24;repmat(1/12,11,1);1/24]; yS = conv(y,sW13,'same'); yS(1:6) = yS(7); yS(T-5:T) = yS(T-6); xt = y-yS; h = plot(yS/1000,'r','LineWidth',2); legend(h,'13-Term Moving Average') hold off
Детрендированными временными рядами является xt
.
Используя параметр формы 'same'
при вызове conv
возвращает сглаживавший ряд та же длина как исходный ряд.
Создайте массив ячеек, sidx
, сохранить индексы, соответствующие каждому периоду. Данные ежемесячно, с периодичностью 12, таким образом, первый элемент sidx
вектор с элементами 1, 13, 25..., 61 (соответствие наблюдениям в январе). Второй элемент sidx
вектор с элементами 2, 14, 16..., 62 (соответствие наблюдениям в феврале). Это повторяется в течение всех 12 месяцев.
s = 12; sidx = cell(s,1); for i = 1:s sidx{i,1} = i:s:T; end sidx{1:2}
ans = 1×6
1 13 25 37 49 61
ans = 1×6
2 14 26 38 50 62
Используя массив ячеек, чтобы сохранить индексы допускает возможность, что каждый период не происходит то же число раз в промежутке наблюдаемого ряда.
Примените устойчивый сезонный фильтр к детрендированному ряду, xt
. Используя индексы, созданные на Шаге 3, усредните детрендированные данные, соответствующие каждому периоду. Таким образом, среднее значение все январские значения (в индексах 1, 13, 25..., 61), и затем среднее значение все февральские значения (в индексах 2, 14, 26..., 62), и так далее в течение остающихся месяцев. Отложите сглаживавшие значения в один вектор.
Сосредоточьте сезонную оценку, чтобы колебаться вокруг нуля.
sst = cellfun(@(x) mean(xt(x)),sidx); % Put smoothed values back into a vector of length N nc = floor(T/s); % no. complete years rm = mod(T,s); % no. extra months sst = [repmat(sst,nc,1);sst(1:rm)]; % Center the seasonal estimate (additive) sBar = mean(sst); % for centering sst = sst-sBar; figure plot(sst/1000) title 'Stable Seasonal Component'; h2 = gca; h2.XLim = [0 T]; ylabel 'Number of deaths (thousands)'; h2.XTick = 1:12:T; h2.XTickLabel = datestr(dates(1:12:T),10);
Устойчивый сезонный компонент имеет постоянную амплитуду через ряд. Сезонная оценка сосредоточена и колеблется вокруг нуля.
Вычтите предполагаемый сезонный компонент из исходных данных.
dt = y - sst; figure plot(dt/1000) title 'Deseasonalized Series'; ylabel 'Number of deaths (thousands)'; h3 = gca; h3.XLim = [0 T]; h3.XTick = 1:12:T; h3.XTickLabel = datestr(dates(1:12:T),10);
deseasonalized ряд состоит из долгосрочного тренда и неправильных компонентов. Крупномасштабный квадратичный тренд в количестве смертей от несчастного случая ясен с сезонным удаленным компонентом.
Ссылки:
Броквелл, P. J. и Р. А. Дэвис. Введение во Временные ряды и Прогнозирование. 2-й редактор Нью-Йорк, Нью-Йорк: Спрингер, 2002.