Вейвлет Обнаружения точек изменения

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

Первый пример применяется вейвлетом обнаружения точек изменения к очень старым временным рядам - данным по минимумам по реке Нил за 622 - 1281 годы нашей эры. Минимальные значения уровня реки измерялись на манометре Роды под Каиром. Измерения указаны в метрах.

Загрузите и постройте график данных.

load nileriverminima
years = 622:1284;
figure
plot(years,nileriverminima)
title('Nile River Minimum Levels')
AX = gca;
AX.XTick = 622:100:1222;
grid on
xlabel('Year')
ylabel('Meters')

Конструкция началась на новом измерительном устройстве около 715 года нашей эры. Изучение данных до и после приблизительно 722 года нашей эры, по-видимому, изменило изменчивость данных. Можно использовать вейвлеты, чтобы исследовать гипотезу о том, что на изменчивость измерений повлияло введение нового измерительного прибора.

Получите мультирезолюционный анализ (MRA) данных с помощью вейвлета Haar.

wt = modwt(nileriverminima,'haar',4);
mra = modwtmra(wt,'haar');

Постройте график MRA и фокусируйте внимание на деталях уровня 1 и уровня 2.

figure
subplot(2,1,1)
plot(years,mra(1,:))
title('Level 1 Details')
subplot(2,1,2)
plot(years,mra(2,:))
title('Level 2 Details')
AX = gca;
AX.XTick = 622:100:1222;
xlabel('Years')

Примените тест общего изменения отклонения к вейвлет.

for JJ = 1:5
    pts_Opt = wvarchg(wt(JJ,:),2);
    changepoints{JJ} = pts_Opt;
end
cellfun(@(x) ~isempty(x),changepoints,'uni',0)
ans =

  1x5 cell array

    {[1]}    {[0]}    {[0]}    {[0]}    {[0]}

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

years(cell2mat(changepoints))
ans =

   721

Разделите данные на два сегмента. Первый сегмент включает в себя годы с 622 по 721, когда мелкомасштабные коэффициенты вейвлета указывают на изменение отклонения. Второй сегмент содержит годы с 722 по 1284. Получите объективные оценки отклонения вейвлета по шкале.

tspre = nileriverminima(1:100);
tspost = nileriverminima(101:end);
wpre = modwt(tspre,'haar',4);
wpost = modwt(tspost,'haar',4);
wvarpre = modwtvar(wpre,'haar',0.95,'table')
wvarpost = modwtvar(wpost,'haar',0.95,'table')
wvarpre =

  5x4 table

          NJ     Lower      Variance     Upper 
          __    ________    ________    _______

    D1    99     0.25199     0.36053    0.55846
    D2    97     0.15367     0.25149    0.48477
    D3    93    0.056137     0.11014    0.30622
    D4    85    0.018881    0.047427    0.26453
    S4    85    0.017875      0.0449    0.25044


wvarpost =

  5x4 table

          NJ      Lower      Variance     Upper 
          ___    ________    ________    _______

    D1    562     0.11394     0.13354    0.15869
    D2    560    0.085288     0.10639    0.13648
    D3    556      0.0693    0.094168    0.13539
    D4    548    0.053644    0.081877    0.14024
    S4    548     0.24608     0.37558    0.64329

Сравните результаты.

Vpre = table2array(wvarpre);
Vpost = table2array(wvarpost);
Vpre = Vpre(1:end-1,2:end);
Vpost = Vpost(1:end-1,2:end);

Vpre(:,1) = Vpre(:,2)-Vpre(:,1);
Vpre(:,3) = Vpre(:,3)-Vpre(:,2);

Vpost(:,1) = Vpost(:,2)-Vpost(:,1);
Vpost(:,3) = Vpost(:,3)-Vpost(:,2);

figure
errorbar(1:4,Vpre(:,2),Vpre(:,1),Vpre(:,3),'ko',...
    'MarkerFaceColor',[0 0 0])
hold on
errorbar(1.5:4.5,Vpost(:,2),Vpost(:,1),Vpost(:,3),'b^',...
    'MarkerFaceColor',[0 0 1])
set(gca,'xtick',1.25:4.25)
set(gca,'xticklabel',{'2 Year','4 Years','8 Years','16 Years','32 Years'})
grid on
ylabel('Variance')
title('Wavelet Variance 622-721 and 722-1284 by Scale','fontsize',14)
legend('Years 622-721','Years 722-1284','Location','NorthEast')

Отклонение вейвлета указывает на значительное изменение отклонения между данными 622-721 и 722-1284 в шкалах 2 и 4 лет.

В приведенном выше примере использовался вейвлет Haar только с двумя коэффициентами из-за проблем с граничными эффектами с относительно небольшими временными рядами (100 выборок из 622-721). Если ваши данные являются приблизительно стационарными различиями первого или второго порядка, можно заменить смещенную оценку с помощью контура 'reflection'. Это позволяет вам использовать более длинный вейвлет, не беспокоясь о граничных коэффициентах. Повторите анализ с помощью вейвлета 'sym4' по умолчанию.

wpre = modwt(tspre,4,'reflection');
wpost = modwt(tspost,4,'reflection');
wvarpre = modwtvar(wpre,[],[],'EstimatorType','biased',...
    'Boundary','reflection','table');
wvarpost = modwtvar(wpost,[],[],'EstimatorType','biased',...
    'Boundary','reflection','table');

Постройте график результатов.

Vpre = table2array(wvarpre);
Vpost = table2array(wvarpost);
Vpre = Vpre(1:end-1,2:end);
Vpost = Vpost(1:end-1,2:end);

Vpre(:,1) = Vpre(:,2)-Vpre(:,1);
Vpre(:,3) = Vpre(:,3)-Vpre(:,2);

Vpost(:,1) = Vpost(:,2)-Vpost(:,1);
Vpost(:,3) = Vpost(:,3)-Vpost(:,2);

figure
errorbar(1:4,Vpre(:,2),Vpre(:,1),Vpre(:,3),'ko','MarkerFaceColor',[0 0 0])
hold on
errorbar(1.5:4.5,Vpost(:,2),Vpost(:,1),Vpost(:,3),'b^','MarkerFaceColor',[0 0 1])
set(gca,'xtick',1.25:4.25)
set(gca,'xticklabel',{'2 Years','4 Years', '8 Years', '16 Years','32 Years'})
grid on
ylabel('Variance')
title({'Wavelet Variance 622-721 and 722-1284 by Scale'; ...
    'Biased Estimate -- Reflection Boundary'},'fontsize',14)
legend('622-721','722-1284','Location','NorthEast')
hold off

Вывод подкреплен. Существует значительное различие в отклонении данных по шкалам 2 и 4 года, но не по более длинным шкалам. Можно сделать вывод, что в отклонении процесса что-то изменилось.

В финансовых временных рядах можно использовать вейвлеты для обнаружения изменений волатильности. Чтобы проиллюстрировать это, рассмотрим ежеквартальные взвешенные по цепочке данные о реальном ВВП США для 1974Q1 к 2012Q4. Данные были преобразованы путем сначала взятия естественного логарифма, а затем вычисления различия за год. Получите вейвлет-преобразование (MODWT) реальных данных ВВП до уровня шесть с вейвлетом 'db2'. Исследуйте дисперсию данных и сравните ее с отклонениями по шкале, полученной с помощью MODWT.

load GDPcomponents
realgdpwt = modwt(realgdp,'db2',6,'reflection');
gdpmra = modwtmra(realgdpwt,'db2','reflection');
vardata = var(realgdp,1);
varwt = var(realgdpwt(:,1:numel(realgdp)),1,2);

В vardata у вас есть отклонение для совокупных временных рядов ВВП. В varwt у вас есть отклонение по шкале для MODWT. В varwt семь элементов поскольку вы получили MODWT до шестого уровня, что привело к шести дисперсиям вейвлет-коэффициентов и одной дисперсии масштабирующего коэффициента. Суммируйте отклонения по шкале, чтобы увидеть, что отклонение сохранено. Постройте график отклонений вейвлета с помощью шкалы, игнорирующей дисперсию коэффициента масштабирования.

totalMODWTvar = sum(varwt);
bar(varwt(1:end-1,:))
AX = gca;
AX.XTickLabels = {'[2 4)','[4 8)','[8 16)','[16 32)','[32 64)','[64 128)'};
xlabel('Quarters')
ylabel('Variance')
title('Wavelet Variance by Scale')

Поскольку эти данные являются ежеквартальными, первая шкала захватывает изменения между двумя и четырьмя четвертями, вторая шкала между четырьмя и восемью, третья между 8 и 16 и так далее.

Из MODWT и простой столбиковой диаграммы видно, что циклы в данных от 8 до 32 кварталов составляют наибольшее отклонение в данных ВВП. Если учесть отклонения вейвлет на этих шкалах, они составляют 57% изменчивости в данных ВВП. Это означает, что колебания в ВВП в течение периода от 2 до 8 лет составляют большую часть изменчивости, наблюдаемой во временных рядах.

Постройте график деталей уровня 1, D1. Эти детали захватывают колебания в данных между двумя и четырьмя четвертями в длительности.

helperFinancialDataExample1(gdpmra(1,:),years,...
    'Year over Year Real U.S. GDP - D1')

Исследуя детали уровня один, кажется, что сокращение отклонения начиная с 1980-х годов.

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

pts_Opt = wvarchg(realgdpwt(1,1:numel(realgdp)),2);
years(pts_Opt)
ans =

        1982

Существует точка изменения отклонения, идентифицированная в 1982 году. Этот пример не исправляет задержку, введенную вейвлетом 'db2' на первом уровне. Однако эта задержка составляет только две выборки, поэтому она не оказывает заметного влияния на результаты.

Чтобы оценить изменения в волатильности данных ВВП до и после 1982 года, разделите исходные данные на серии до и после изменений. Получите вейвлет набора наборов данных до и после. В этом случае ряд относительно короток, поэтому используйте wavelet Haar, чтобы минимизировать количество граничных коэффициентов. Вычислите объективные оценки отклонения вейвлета по шкале и постройте график результата.

tspre = realgdp(1:pts_Opt);
tspost = realgdp(pts_Opt+1:end);
wtpre = modwt(tspre,'haar',5);
wtpost = modwt(tspost,'haar',5);
prevar = modwtvar(wtpre,'haar','table');
postvar = modwtvar(wtpost,'haar','table');
xlab = {'[2Q,4Q)','[4Q,8Q)','[8Q,16Q)','[16Q,32Q)','[32Q,64Q)'};
helperFinancialDataExampleVariancePlot(prevar,postvar,'table',xlab)
title('Wavelet Variance By Scale')
legend('Pre 1982 Q2','Post 1982 Q2','Location','NorthWest')

Из предыдущего графика видно, что существуют значительные различия между pre-1982Q2 и post-1982Q2 отклонениями в шкалах от 2 до 16 четвертей.

Поскольку временные ряды в этом примере настолько коротки, может быть полезно использовать смещенные оценки отклонения. Смещенные оценки не удаляют граничные коэффициенты. Используйте вейвлет 'db2' с четырьмя коэффициентами.

wtpre = modwt(tspre,'db2',5,'reflection');
wtpost = modwt(tspost,'db2',5,'reflection');
prevar = modwtvar(wtpre,'db2',0.95,'EstimatorType','biased','table');
postvar = modwtvar(wtpost,'db2',0.95,'EstimatorType','biased','table');
xlab = {'[2Q,4Q)','[4Q,8Q)','[8Q,16Q)','[16Q,32Q)','[32Q,64Q)'};
figure
helperFinancialDataExampleVariancePlot(prevar,postvar,'table',xlab)
title('Wavelet Variance By Scale')
legend('Pre 1982 Q2','Post 1982 Q2','Location','NorthWest')

Результаты подтверждают наше первоначальное открытие о снижении волатильности по шкалам с 2 до 16 четвертей.

Использование вейвлет позволяет вам фокусироваться на шкалах, где изменение волатильности локализовано. Чтобы увидеть это, исследуйте график необработанных данных вместе с деталями вейвлета первого уровня.

subplot(2,1,1)
helperFinancialDataExample1(realgdp,years,...
    'Year over Year Real U.S. GDP -- Raw Data')
subplot(2,1,2)
helperFinancialDataExample1(gdpmra(1,:),years,...
    'Year over Year Real U.S. GDP -- Wavelet Level 1 Details')

Затенённая область упоминается как «Великая умеренность», означающая период пониженной макроэкономической волатильности в США, начиная с середины 1980-х годов.

Исследуя совокупные данные, неясно, что на самом деле в этот период снижается волатильность. Однако вейвлет подробно обнаруживает изменение волатильности.

Для просмотра документации необходимо авторизоваться на сайте