Мультифрактальный анализ

В этом примере показано, как использовать вейвлеты, чтобы охарактеризовать локальную регулярность сигнала. Способность описать регулярность сигнала важна при контакте с явлениями, которые не имеют никакой характеристической шкалы. Сигналы с динамикой без шкал широко наблюдаются во многих различных областях применения включая биомедицинскую обработку сигналов, геофизику, финансы и интернет-трафик. Каждый раз, когда вы применяете некоторый аналитический метод к своим данным, вы неизменно принимаете что-то о данных. Например, если вы используете автокорреляцию или оценку спектральной плотности мощности (PSD), вы принимаете, что ваши данные являются инвариантом перевода, что означает, что статистические данные сигнала как среднее значение и отклонение не изменяются в зависимости от времени. Сигналы без характеристической шкалы инвариантны к масштабу. Это означает, что статистические данные сигнала не изменяются, если мы расширяем или уменьшаем ось времени. Классическим методам обработки сигналов обычно не удается соответственно описать эти сигналы или показать различия между сигналами с различным поведением масштабирования. В этих случаях фрактальный анализ может обеспечить уникальное понимание. Некоторые следующие примеры используют pwelch для рисунка. Чтобы выполнить тот код, у вас должен быть Signal Processing Toolbox™.

Процессы закона о степени

Важный класс сигналов с динамикой без шкал имеет автокорреляцию или спектральные плотности мощности (PSD), которые следуют закону о степени. Процесс закона степени имеет PSD формы C|ω|-α для некоторой положительной константы, C, и некоторая экспонента α. В некоторых случаях сигнал интереса показывает закон степени PSD. В других случаях сигнал интереса повреждается шумом с законом степени PSD. Эти шумы часто упоминаются, как окрашено. Способность оценить экспоненту от реализации этих процессов имеет важные последствия. Для одного это позволяет вам делать выводы о механизме, генерирующем данные, а также представляющем эмпирические свидетельства, чтобы поддержать или отклонить теоретические предсказания. В случае вмешивающегося шума с законом степени PSD это полезно в разработке эффективных фильтров.

Шум Брауна или Броуновский процесс, является одним таким цветным шумовым процессом с теоретической экспонентой α=2. Один способ оценить экспоненту процесса закона о степени состоит в том, чтобы соответствовать линии наименьших квадратов к графику логарифмического журнала PSD.

load brownnoise;
[Pxx,F] = pwelch(brownnoise,kaiser(1000,10),500,1000,1);
plot(log10(F(2:end)),log10(Pxx(2:end)));
grid on;
xlabel('log10(F)'); ylabel('log10(Pxx)');
title('Log-Log Plot of PSD Estimate')

Figure contains an axes object. The axes object with title Log-Log Plot of PSD Estimate contains an object of type line.

Регресс журнал значения PSD на логарифмических частотах. Обратите внимание, что необходимо проигнорировать нулевую частоту, чтобы не брать журнал нуля.

Xpred = [ones(length(F(2:end)),1) log10(F(2:end))];
b = lscov(Xpred,log10(Pxx(2:end)));
y = b(1)+b(2)*log10(F(2:end));
hold on;
plot(log10(F(2:end)),y,'r--');
title(['Estimated Slope is ' num2str(b(2))]);

Figure contains an axes object. The axes object with title Estimated Slope is -1.8652 contains 2 objects of type line.

В качестве альтернативы можно использовать и дискретные и непрерывные аналитические методы вейвлета, чтобы оценить экспоненту. Отношение между экспонентой Держателя, H, возвращенный dwtleader и wtmm и α этот сценарий, α=2H+1.

[dhbrown,hbrown,cpbrown] = dwtleader(brownnoise);
hexp = wtmm(brownnoise);
fprintf('Wavelet leader estimate is %1.2f\n',-2*cpbrown(1)-1);
Wavelet leader estimate is -1.91
fprintf('WTMM estimate is %1.2f\n',-2*hexp-1);
WTMM estimate is -2.00

В этом случае оценка, полученная путем подбора кривой линии наименьших квадратов к журналу оценки PSD и тех полученные методы вейвлета использования, находится в хорошем соглашении.

Мультифрактальный анализ

Существует много реальных сигналов, которые предоставляют нелинейное поведение закона степени, которое зависит от моментов высшего порядка и шкалы. Мультифрактальный анализ обеспечивает способ описать эти сигналы. Мультифрактальный анализ состоит из определения, существует ли некоторый тип масштабирования закона степени в течение различных статистических моментов в различных шкалах. Если это поведение масштабирования характеризуется одной экспонентой масштабирования, или эквивалентно является линейной функцией моментов, процесс является монофракталом. Если масштабирующееся поведение шкалой является нелинейной функцией моментов, процесс является мультифракталом. Коричневый шум от предыдущего раздела является примером монофрактального процесса, и это продемонстрировано в более позднем разделе.

Чтобы проиллюстрировать, как фрактальный анализ может показать структуру сигнала, не очевидную с более классическими методами обработки сигналов, загрузите RWdata.mat который содержит два временных рядов (Ts1 и Ts2) с 8 000 выборок каждый. Отобразите данные на графике.

load RWdata;
figure;
plot([Ts1 Ts2]); grid on;
legend('Ts1','Ts2','Location','NorthEast');
xlabel('Time'); ylabel('Amplitude');

Figure contains an axes object. The axes object contains 2 objects of type line. These objects represent Ts1, Ts2.

Сигналы имеют очень похожую статистику второго порядка. Если вы смотрите на средние значения, значения RMS и отклонения Ts1 и Ts2, значения почти идентичны. Оценки PSD также очень похожи.

pwelch([Ts1 Ts2],kaiser(1e3,10))

Figure contains an axes object. The axes object with title Welch Power Spectral Density Estimate contains 2 objects of type line.

Последовательности автокорреляции затухают очень медленно и для временных рядов и являются весьма формирующими для дифференциации временных рядов.

[xc1,lags] = xcorr(Ts1,300,'coef');
xc2 = xcorr(Ts2,300,'coef');
subplot(2,1,1)
hs1 = stem(lags(301:end),xc1(301:end));
hs1.Marker = 'none';
title('Autocorrelation Sequence of Ts1');
subplot(2,1,2)
hs2 = stem(lags(301:end),xc2(301:end));
hs2.Marker = 'none';
title('Autocorrelation Sequence of Ts2');
xlabel('Lag')

Figure contains 2 axes objects. Axes object 1 with title Autocorrelation Sequence of Ts1 contains an object of type stem. Axes object 2 with title Autocorrelation Sequence of Ts2 contains an object of type stem.

Даже в задержке 300, автокорреляции 0.94 и 0.96 соответственно.

То, что эти сигналы очень отличаются, показано посредством фрактального анализа. Вычислите и постройте мультифрактальные спектры двух сигналов. В мультифрактальном анализе дискретные методы вейвлета на основе так называемых лидеров вейвлета являются самыми устойчивыми.

[dh1,h1,cp1,tauq1] = dwtleader(Ts1);
[dh2,h2,cp2,tauq2] = dwtleader(Ts2);
figure;
hp = plot(h1,dh1,'b-o',h2,dh2,'b-^');
hp(1).MarkerFaceColor = 'b';
hp(2).MarkerFaceColor = 'r';
grid on;
xlabel('h'); ylabel('D(h)');
legend('Ts1','Ts2','Location','NorthEast');
title('Multifractal Spectrum');

Figure contains an axes object. The axes object with title Multifractal Spectrum contains 2 objects of type line. These objects represent Ts1, Ts2.

Мультифрактальный спектр эффективно показывает распределение масштабирующихся экспонент для сигнала. Эквивалентно, мультифрактальный спектр обеспечивает меру того, насколько локальная регулярность сигнала варьируется вовремя. Сигнал, который является монофрактальными выставками по существу та же регулярность везде вовремя и поэтому имеет мультифрактальный спектр с узкой поддержкой. С другой стороны мультифрактальный сигнал показывает изменения регулярности сигнала в зависимости от времени и имеет мультифрактальный спектр с более широкой поддержкой. От мультифрактальных спектров, показанных здесь, Ts2, кажется, монофрактальный сигнал, охарактеризованный кластером масштабирующихся экспонент приблизительно 0,78. С другой стороны, Ts1, демонстрирует широкий диапазон масштабирующихся экспонент, указывающих, что это - мультифрактал. Отметьте общую область значений масштабирующихся (Держатель) экспоненты для Ts2 всего 0.14, в то время как это является в 4.6 раза более большим для Ts1. Ts2 на самом деле пример монофрактального фракционного броуновского движения (fBm) процесс с экспонентой Держателя 0,8 и Ts1 мультифрактальный случайный обход.

Можно также использовать масштабирующуюся экспоненту выходные параметры от dwtleader наряду с 2-м cumulant, чтобы помочь классифицировать процесс как монофрактал по сравнению с мультифракталом. Вспомните, что монофрактальный процесс имеет линейное правило масштабирования в зависимости от статистических моментов, в то время как мультифрактальный процесс имеет нелинейное правило масштабирования. dwtleader использует область значений моментов от-5 до 5 в оценке этих правил масштабирования. График масштабирующихся экспонент для fBm и процесса мультифрактального случайного обхода (MRW) показывает различие.

plot(-5:5,tauq1,'b-o',-5:5,tauq2,'r-^');
grid on;
xlabel('Q-th Moment'); ylabel('Scaling Exponents');
title('Scaling Exponents');
legend('MRW','fBm','Location','SouthEast');

Figure contains an axes object. The axes object with title Scaling Exponents contains 2 objects of type line. These objects represent MRW, fBm.

Масштабирующиеся экспоненты для fBm процесса являются линейной функцией моментов, в то время как экспоненты для процесса MRW показывают отклонение от линейности. Та же информация получена в итоге 1-м, 2-м, и 3-м cumulants. Первый cumulant является оценкой наклона, другими словами, это получает линейное поведение. Второй cumulant получает первое отклонение от линейности. Можно думать о втором cumulant как о коэффициентах для (квадратичного) термина второго порядка, в то время как третий cumulant характеризует более сложный отъезд масштабирующихся экспонент от линейности. Если вы исследуете 2-й и 3-й cumulants на процесс MRW, они в 6 и 42 раза более большие, чем соответствующий cumulants для fBm данных. В последнем случае 2-е и 3-и cumulants являются почти нулем как ожидалось.

Для сравнения добавьте мультифрактальный спектр для коричневого шума, вычисленного в более раннем примере.

hp = plot(h1,dh1,'b-o',h2,dh2,'b-^',hbrown,dhbrown,'r-v');
hp(1).MarkerFaceColor = 'b';
hp(2).MarkerFaceColor = 'r';
hp(3).MarkerFaceColor = 'k';
grid on;
xlabel('h'); ylabel('D(h)');
legend('Ts1','Ts2','brown noise','Location','SouthEast');
title('Multifractal Spectrum');

Figure contains an axes object. The axes object with title Multifractal Spectrum contains 3 objects of type line. These objects represent Ts1, Ts2, brown noise.

Куда Процесс Идет Затем? Персистентное и Антиперсистентное Поведение

Оба дробный Броуновский процесс (Ts2) и коричневые шумовые ряды являются монофракталом. Однако простой график двух временных рядов показывает, что они кажутся очень отличающимися.

subplot(2,1,1)
plot(brownnoise); title('Brown Noise');
grid on; axis tight;
subplot(2,1,2)
plot(Ts2); title('fBm H=0.8'); grid on;
axis tight;

Figure contains 2 axes objects. Axes object 1 with title Brown Noise contains an object of type line. Axes object 2 with title fBm H=0.8 contains an object of type line.

fBm данные намного более сглаженны, чем коричневый шум. Шум Брауна, также известный как случайный обход, имеет теоретическую экспоненту Держателя 0,5. Это значение формирует контур между процессами с экспонентами Держателя, H, от 0 <H <0.5 и те с экспонентами Держателя в интервале 0.5 <H <1. Первый называется антиперсистентный и показывает короткую память. Последние называются персистентные и показывают хорошую память. В антиперсистентных временных рядах увеличение значения во время t сопровождается с уменьшением в значении во время t+1 с высокой вероятностью. Точно так же уменьшение в значении во время t обычно сопровождается увеличением значения во время t+1. Другими словами, временные ряды имеют тенденцию всегда возвращаться к его среднему значению. В персистентных временных рядах увеличения значения имеют тенденцию сопровождаться последующими увеличениями, в то время как уменьшения в значении имеют тенденцию сопровождаться последующими уменьшениями.

Чтобы видеть некоторые реальные примеры антиперсистентных временных рядов, загрузите и анализируйте ежедневный журнал, возвращается для Взвешенного Тайваня и Сеульские индексы запаса Составного объекта. Ежедневная газета возвращает для обоих покрытий индексов аппроксимированный период с июля 1997 до апреля 2016.

load StockCompositeData;
subplot(2,1,1)
plot(SeoulComposite); title('Seoul Composite Index - 07/1997-04/2016');
ylabel('Log Returns'); grid on;
subplot(2,1,2);
plot(TaiwanWeighted); title('Taiwan Weighted Index - 07/1997-04/2016');
ylabel('Log Returns');
grid on;

Figure contains 2 axes objects. Axes object 1 with title Seoul Composite Index - 07/1997-04/2016 contains an object of type line. Axes object 2 with title Taiwan Weighted Index - 07/1997-04/2016 contains an object of type line.

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

[dhseoul,hseoul,cpseoul] = dwtleader(SeoulComposite);
[dhtaiwan,htaiwan,cptaiwan] = dwtleader(TaiwanWeighted);
figure;
plot(hseoul,dhseoul,'b-o','MarkerFaceColor','b');
hold on;
plot(htaiwan,dhtaiwan,'r-^','MarkerFaceColor','r');
xlabel('h'); ylabel('D(h)'); grid on;
title('Multifractal Spectrum');

Figure contains an axes object. The axes object with title Multifractal Spectrum contains 2 objects of type line.

От мультифрактального спектра ясно, что оба временных рядов являются антиперсистентными. Для сравнения постройте мультифрактальные спектры двух финансовых временных рядов наряду с коричневым шумом и fBm данными, показанными ранее.

plot(hbrown,dhbrown,'k-v','MarkerFaceColor','k');
plot(h2,dh2,'b-*','MarkerFaceColor','b');
legend('Seoul Composite','Taiwan Weighted Index','Brown Noise','FBM',...
    'Location','SouthEast');
hold off;

Figure contains an axes object. The axes object with title Multifractal Spectrum contains 4 objects of type line. These objects represent Seoul Composite, Taiwan Weighted Index, Brown Noise, FBM.

Решение, что процесс является антиперсистентным или персистентный, полезно в предсказании будущего. Например, временные ряды с хорошей памятью, которая увеличивается, как могут ожидать, продолжат увеличиваться. В то время как временные ряды, которые показывают антиперсистентность, как могут ожидать, переместятся в противоположное направление.

Измерение фрактальной динамики изменчивости сердечного ритма

Нормальная изменчивость уровня человеческого сердца, измеренная как интервалы RR, отображает мультифрактальное поведение. Далее, сокращения этого нелинейного поведения масштабирования являются хорошими предикторами сердечной болезни и даже выхода из строя.

Как пример вызванного изменения во фрактальной динамике изменчивости сердечного ритма, полагайте, что пациент администрировал простагландин E1 из-за серьезного гипертонического эпизода. Данные являются частью RHRV, основанного на R пакета программного обеспечения для анализа изменчивости сердечного ритма. Авторы любезно дали разрешение для его использования в этом примере.

Загрузите и отобразите данные на графике. Вертикальная красная линия отмечает начало эффекта простагландина E1 на изменчивости сердечного ритма и сердечного ритма.

load hrvDrug;
plot(hrvDrug); grid on;
hold on;
plot([4642 4642],[min(hrvDrug) max(hrvDrug)],'r','linewidth',2);
hold off;
ylabel('Heart Rate'); xlabel('Sample');

Figure contains an axes object. The axes object contains 2 objects of type line.

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

predrug = hrvDrug(1:4642);
postdrug = hrvDrug(4643:end);
[dhpre,hpre] = dwtleader(predrug);
[dhpost,hpost] = dwtleader(postdrug);
figure;
hl = plot(hpre,dhpre,'b-d',hpost,dhpost,'r-^');
hl(1).MarkerFaceColor = 'b';
hl(2).MarkerFaceColor = 'r';
xlabel('h'); ylabel('D(h)');
grid on;
legend('Predrug','Postdrug');
title('Multifractal Spectrum'); xlabel('h'); ylabel('D(h)');

Figure contains an axes object. The axes object with title Multifractal Spectrum contains 2 objects of type line. These objects represent Predrug, Postdrug.

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

Ссылки

Л. Родригес-Линарес, L., А.Х. Мендес, М.Х. Ладо, Д.Н. Оливиери, Кс.А. Вила и я. Гомес-Конде, "Инструмент с открытым исходным кодом для спектрального анализа изменчивости сердечного ритма", Компьютерные Методы и Программы в Биомедицине, 103 (1):39-50,2011.

Вендт, H. и Abry, P. "Тесты Multifractality с помощью загруженных лидеров вейвлета", Сделка IEEE. Обработка сигналов, издание 55, № 10, стр 4811-4820, 2007.

Вендт, H., Abry, P. и Jaffard, S. "Загрузитесь для эмпирического мультифрактального анализа", Журнал Обработки сигналов IEEE, 24, 4, 38-48, 2007.

Jaffard, S., Lashermes, B. и Abry, P. "Лидеры вейвлета в мультифрактальном анализе". В Т. Цяне, М.И. Вае и С. Юэшэне, редакторах. Анализ вейвлета и Приложения, стр 219-264, Birkhauser, 2006.