Пиковый анализ

В этом примере показано, как выполнить основной пиковый анализ. Это поможет вам ответить на вопросы, такие как: Как я нахожу peaks в своем сигнале? Как я измеряю расстояние между peaks? Как я измеряю амплитуду peaks сигнала, который затронут трендом? Как я нахожу peaks в сигнале с шумом? Как я нахожу локальные минимумы?

Нахождение максимумов или Peaks

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

load sunspot.dat
year = sunspot(:,1); 
relNums = sunspot(:,2);
findpeaks(relNums,year)
xlabel('Year')
ylabel('Sunspot Number')
title('Find All Peaks')

Figure contains an axes. The axes with title Find All Peaks contains 2 objects of type line.

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

Измерение расстояния между Peaks

Peaks в сигнале, кажется, появляется равномерно. Однако часть peaks очень друг близко к другу. MinPeakProminence свойство может использоваться, отфильтровывают этот peaks. Рассмотрите peaks, который понижается с обеих сторон по крайней мере 40 относительными числами солнечного пятна прежде, чем столкнуться с большим значением.

findpeaks(relNums,year,'MinPeakProminence',40)
xlabel('Year')
ylabel('Sunspot Number')
title('Find Prominent Peaks')

Figure contains an axes. The axes with title Find Prominent Peaks contains 2 objects of type line.

Следующая гистограмма показывает распределение пиковых интервалов в годах:

figure
[pks, locs] = findpeaks(relNums,year,'MinPeakProminence',40);
peakInterval = diff(locs);
histogram(peakInterval)
grid on
xlabel('Year Intervals')
ylabel('Frequency of Occurrence')
title('Histogram of Peak Intervals (years)')

Figure contains an axes. The axes with title Histogram of Peak Intervals (years) contains an object of type histogram.

AverageDistance_Peaks = mean(diff(locs))
AverageDistance_Peaks = 10.9600

Распределение показывает, что большинство пиковых интервалов находится между 10 и 12 годами, указывая, что сигнал имеет циклическую природу. Кроме того, средний интервал 10,96 лет между peaks совпадает с известным циклическим действием солнечного пятна 11 лет.

Нахождение Peaks в отсеченных или влажных сигналах

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

load clippedpeaks.mat

figure

% Show all peaks in the first plot
ax(1) = subplot(2,1,1);
findpeaks(saturatedData)
xlabel('Samples')
ylabel('Amplitude')
title('Detecting Saturated Peaks')

% Specify a minimum excursion in the second plot
ax(2) = subplot(2,1,2);
findpeaks(saturatedData,'threshold',5)
xlabel('Samples')
ylabel('Amplitude')
title('Filtering Out Saturated Peaks')

% link and zoom in to show the changes
linkaxes(ax(1:2),'xy')
axis(ax,[50 70 0 250])

Figure contains 2 axes. Axes 1 with title Detecting Saturated Peaks contains 2 objects of type line. Axes 2 with title Filtering Out Saturated Peaks contains 2 objects of type line.

Первый подграфик показывает, что в случае плоского пика, возрастающее ребро обнаруживается как пик. Второй подграфик показывает, что определение порога может помочь отклонить плоский peaks.

Измерение амплитуд Peaks

Этот пример показывает пиковый анализ в ECG (Электрокардиограмма) сигнал. ECG является мерой электрического действия основы в зависимости от времени. Сигнал измеряется электродами, присоединенными к коже, и чувствителен к воздействиям, таким как интерференция источника питания и шумы из-за артефактов перемещения.

load noisyecg.mat
t = 1:length(noisyECG_withTrend);

figure
plot(t,noisyECG_withTrend)
title('Signal with a Trend')
xlabel('Samples');
ylabel('Voltage(mV)')
legend('Noisy ECG Signal')
grid on

Figure contains an axes. The axes with title Signal with a Trend contains an object of type line. This object represents Noisy ECG Signal.

Вычитание тренда из данных

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

[p,s,mu] = polyfit((1:numel(noisyECG_withTrend))',noisyECG_withTrend,6);
f_y = polyval(p,(1:numel(noisyECG_withTrend))',[],mu);

ECG_data = noisyECG_withTrend - f_y;        % Detrend data

figure
plot(t,ECG_data)
grid on
ax = axis;
axis([ax(1:2) -1.2 1.2])
title('Detrended ECG Signal')
xlabel('Samples')
ylabel('Voltage(mV)')
legend('Detrended ECG Signal')

Figure contains an axes. The axes with title Detrended ECG Signal contains an object of type line. This object represents Detrended ECG Signal.

После удаления тренда найдите комплекс QRS, который является самым видным пиком повторения в сигнале ECG. Комплекс QRS соответствует деполяризации правых и левых желудочков человеческого сердца. Это может использоваться, чтобы определить сердечный уровень пациента или предсказать отклонения в сердечной функции. Следующий рисунок показывает форму комплекса QRS в сигнале ECG.

Пороговая обработка, чтобы найти Peaks интересным

Комплекс QRS состоит из трех главных компонентов: Q волна, R волна, S волна. Волны R могут быть обнаружены peaks пороговой обработки выше 0,5 мВ. Заметьте, что волны R разделяются больше чем 200 выборками. Используйте эту информацию, чтобы удалить нежелательный peaks путем определения 'MinPeakDistance'.

[~,locs_Rwave] = findpeaks(ECG_data,'MinPeakHeight',0.5,...
                                    'MinPeakDistance',200);

Для обнаружения S-волн найдите локальные минимумы в сигнале и примените пороги соответственно.

Нахождение локальных минимумов в сигнале

Локальные минимумы могут быть обнаружены путем нахождения peaks на инвертированной версии исходного сигнала.

ECG_inverted = -ECG_data;
[~,locs_Swave] = findpeaks(ECG_inverted,'MinPeakHeight',0.5,...
                                        'MinPeakDistance',200);

Следующий график показывает волны R и волны S, обнаруженные в сигнале.

figure
hold on 
plot(t,ECG_data)
plot(locs_Rwave,ECG_data(locs_Rwave),'rv','MarkerFaceColor','r')
plot(locs_Swave,ECG_data(locs_Swave),'rs','MarkerFaceColor','b')
axis([0 1850 -1.1 1.1])
grid on
legend('ECG Signal','R waves','S waves')
xlabel('Samples')
ylabel('Voltage(mV)')
title('R wave and S wave in Noisy ECG Signal')

Figure contains an axes. The axes with title R wave and S wave in Noisy ECG Signal contains 3 objects of type line. These objects represent ECG Signal, R waves, S waves.

Затем мы пытаемся определить местоположения волн Q. Пороговая обработка peaks, чтобы определить местоположение результатов волн Q в обнаружении нежелательного peaks как волны Q содержится шум. Мы фильтруем сигнал сначала и затем находим peaks. Фильтрация Savitzky-Golay используется, чтобы удалить шум в сигнале.

smoothECG = sgolayfilt(ECG_data,7,21);

figure
plot(t,ECG_data,'b',t,smoothECG,'r')
grid on
axis tight
xlabel('Samples')
ylabel('Voltage(mV)')
legend('Noisy ECG Signal','Filtered Signal')
title('Filtering Noisy ECG Signal')

Figure contains an axes. The axes with title Filtering Noisy ECG Signal contains 2 objects of type line. These objects represent Noisy ECG Signal, Filtered Signal.

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

[~,min_locs] = findpeaks(-smoothECG,'MinPeakDistance',40);

% Peaks between -0.2mV and -0.5mV
locs_Qwave = min_locs(smoothECG(min_locs)>-0.5 & smoothECG(min_locs)<-0.2);

figure
hold on
plot(t,smoothECG); 
plot(locs_Qwave,smoothECG(locs_Qwave),'rs','MarkerFaceColor','g')
plot(locs_Rwave,smoothECG(locs_Rwave),'rv','MarkerFaceColor','r')
plot(locs_Swave,smoothECG(locs_Swave),'rs','MarkerFaceColor','b')
grid on
title('Thresholding Peaks in Signal')
xlabel('Samples')
ylabel('Voltage(mV)')
ax = axis;
axis([0 1850 -1.1 1.1])
legend('Smooth ECG signal','Q wave','R wave','S wave')

Figure contains an axes. The axes with title Thresholding Peaks in Signal contains 4 objects of type line. These objects represent Smooth ECG signal, Q wave, R wave, S wave.

Вышеупомянутый рисунок показывает, что QRS объединяют успешно обнаруженный в шумном сигнале ECG.

Ошибка между шумным и сглаженным сигналом

Заметьте среднюю разницу между комплексом QRS в сырых данных и детрендированным отфильтрованным сигналом.

% Values of the Extrema
[val_Qwave, val_Rwave, val_Swave] = deal(smoothECG(locs_Qwave), smoothECG(locs_Rwave), smoothECG(locs_Swave));

meanError_Qwave = mean((noisyECG_withTrend(locs_Qwave) - val_Qwave))
meanError_Qwave = 0.2771
meanError_Rwave = mean((noisyECG_withTrend(locs_Rwave) - val_Rwave))
meanError_Rwave = 0.3476
meanError_Swave = mean((noisyECG_withTrend(locs_Swave) - val_Swave))
meanError_Swave = 0.1844

Это демонстрирует, что важно детрендировать сигнал с шумом для эффективного пикового анализа.

Пиковые свойства

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

avg_riseTime = mean(locs_Rwave-locs_Qwave); % Average Rise time
avg_fallTime = mean(locs_Swave-locs_Rwave); % Average Fall time
avg_riseLevel = mean(val_Rwave-val_Qwave);  % Average Rise Level
avg_fallLevel = mean(val_Rwave-val_Swave);  % Average Fall Level

helperPeakAnalysisPlot(t,smoothECG,...
                    locs_Qwave,locs_Rwave,locs_Swave,...
                    val_Qwave,val_Rwave,val_Swave,...
                    avg_riseTime,avg_fallTime,...
                    avg_riseLevel,avg_fallLevel)

Figure contains an axes. The axes with title QRS-complex in an ECG signal contains 11 objects of type line, text. These objects represent QRS-Complex, Peak, Minima.

Смотрите также

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