How to find the valley of a PPG signal? It is getting confused with the dicrotic notch

I used the following code:

[~,index]=findpeaks(PP,'MinPeakDistance',200);

invertedP = - PP;

%findpeaks(invertedP);

[~, indexes] = findpeaks(invertedP,'MinPeakProminence',3,'MinPeakDistance',200);

figure;

plot(PP); hold on; plot(indexes,PP(indexes),'rv','MarkerFaceColor','g');plot(index,PP(index),'rv','MarkerFaceColor','r');

Answer 1

[Star Strider]

It would help to have the data, existing code, (and a more thorough description of the desired result).

In some instances, it is necessary to set two different threshold criteria (for 'MinPeakProminence' or 'MinPeakHeight' or others) and then use the appropriate set operations (such as setdiff) to return only the ones you want.

[Anupama V]

I have attached the dataset file(0127_8min.mat)

data = load('C:\Users\91944\OneDrive\Documents\TBME2013-PPGRR-Benchmark_R3\data\0127_8min.mat');

PP= data.signal.pleth.y;

PP=PP(50000:70000); %first segments

[~,index]=findpeaks(PP,'MinPeakDistance',150);

invertedP = - PP;

%findpeaks(invertedP);

[~, indexes] = findpeaks(invertedP,'MinPeakProminence',3,'MinPeakDistance',150);

figure;

plot(PP); hold on; plot(indexes,PP(indexes),'rv','MarkerFaceColor','g');plot(index,PP(index),'rv','MarkerFaceColor','r');

This is the code for finding peaks and valleys of the ppg signal PP acquired from database.Valley point is not detected here,instead dicrotic notch is detected(As seen in the screenshot attached) Can you please help in correcting the code for finding the valley ?

[Star Strider]

It would help to have some sort of key to whatever are in these structures.

That aside, this is able to distinguish the valleys from the dicrotic dnotches (most of the time) and returns the index locations correctly. It has problems when the diastolic valleys and dicrotic notches are of similar prominences, however adding a highpass filter to eliminate the baseline drift seems to have solved that.

I chose the 1 Hz filter passband empirically based on my experience with physiological signals. It appears to have corrected the problem of misclassified dicriotic notches and diastolic valleys, however I did not rigorously examine the entire signal. (If the empirical approach had not worked, I would have calculated the fft of the signal, and used that result to design the filter.)

Try this —

LD = load(websave('0127_8min','https://www.mathworks.com/matlabcentral/answers/uploaded_files/1183088/0127_8min.mat'))

LD = struct with fields:
    SFresults: [1×1 struct]
       labels: [1×1 struct]
         meta: [1×1 struct]
        param: [1×1 struct]
    reference: [1×1 struct]
       signal: [1×1 struct]

ppg = LD.signal.pleth.y;
L = numel(ppg);
Fs = LD.param.samplingrate.pleth;
t = linspace(0, L-1, L).'/Fs;
% ck = 1/t(2)

ppgf = highpass(ppg, 1, Fs);

[pks,plocs] = findpeaks(ppgf, 'MinPeakProminence',5);                % Systolic Peaks & Indices
[vys,vlocs] = findpeaks(-ppgf, 'MinPeakProminence',5);               % Diastolic Valleys & Indices
[dns,nlocs] = findpeaks(-ppgf, 'MinPeakProminence',1);               % Dicrotic Notehes % Indices

Szm1 = [size(plocs); size(vlocs); size(nlocs)]

Szm1 = 3×2

615 1 615 1 1228 1

Lv = ismember(nlocs,vlocs);                                         % Indices Of Valley Locations In Dicrotic Notch Locations
vlocs = nlocs(Lv);                                                  % 'Pure' Valley Locations
nlocs = nlocs(~Lv);                                                 % 'Pure' Dicrotic Notch Locations
dns = dns(~Lv);                                                     % Dicrotic Notch Values 

Szm2 = [size(plocs); size(vlocs); size(nlocs)]

Szm2 = 3×2

615 1 615 1 613 1

figure
plot(t, ppgf, 'DisplayName','Photoplethysmograph Signal (Highpass Filtered)')
hold on
plot(t(plocs), pks, '^r', 'DisplayName','Systolic Peaks')
plot(t(vlocs), -vys, 'vr', 'DisplayName','Diastolic Valleys')
plot(t(nlocs), -dns, 'sr', 'DisplayName','Dicrotic Notches')
hold off
grid
xlabel('Time (s)')
ylabel('PPG Signal')
legend('location','best')
xlim( [0 10])

.

[Anupama V]

Thank you so much for your help. Can you kindly explain the

Lv = ismember(nlocs,vlocs); vlocs = nlocs(Lv); nlocs = nlocs(~Lv)

dns = dns(~Lv); part?

[Star Strider]

As always, my pleasure!

Sure!

Lv = ismember(nlocs,vlocs); % Indices Of Valley Locations In Dicrotic Notch Locations

vlocs = nlocs(Lv); % 'Pure' Valley Locations

nlocs = nlocs(~Lv); % 'Pure' Dicrotic Notch Locations

dns = dns(~Lv); % Dicrotic Notch Values

The findpeaks call that returns the dicrotic notch locations (‘nlocs’) also returns the diastolic valley locations (‘vlocs’). The ismember call separates them so that it returns a logical vector (‘Lv’) of the values for ‘nlocs’ that are also in ‘vlocs’ however ‘vlocs’ returns only the diastolic valley locations. Thne result of those operations separates them so that after them, ‘vlocs’ contains only the diastolic valley locations and ‘nlocs’ contains only the dicrotic notch locations. The ‘Lv’ variable is then used to return only the dicrotic notch locations in ‘nlocs’ (so they no longer also contain and index into the diastolic valley locations), and the dicrotic notch values in ‘dns’.

It is difficult to describe here, so one way to understand it would be to remove the semicolons from those assignments to see what they produce. That demonstrates that ‘Lv’ returns true values for the locations that are common to both arguments, and false values for locations that are only in ‘nlocs’. Then ‘Lv’ (the true values) are used to identify ‘vlocs’ and its logical complement ‘~Lv’ (the false values) are used to identify ‘nlocs’ and ‘dns’, completely separating them. (In MATLAB, the ~ operator when used with a logical value negates it, so ~true=flase and ~false=true.)

.

[Anupama V]

Got it. Again, thank you , you have helped a lot .

[Star Strider]

As always, my pleasure!

[Anupama V]

In the attached code,all starts of the peak is getting detected. From this how to filter out dicrotic notch and valley? I am attaching my signal below.

LD = load('C:\Users\ANU\Desktop\TBME2013-PPGRR-Benchmark_R3\data\0009_8min.mat.mat')

% remove first trough:

ppg = LD.signal.pleth.y;

ppg(1:10) = [];

% indentify all peaks

[all_peak_value,all_peak_location]= findpeaks(ppg);

% identify all troughs in data

[all_trough_value,all_trough_location]= findpeaks(-ppg);

all_trough_value = -all_trough_value;

% identify main peak (i.e. ones lower than 5)

[main_trough_value,main_trough_location]= findpeaks(ppg,'MinPeakProminence',5);

counter = 0;

for i = 1:length(main_trough_location)

%c=[];

start=main_trough_location(i)-100;

ending=main_trough_location(i);

%c= (all_trough_location>start & all_trough_location<ending);

%c1=c(end);

[r,c] = find(all_trough_location>start & all_trough_location<ending);

end

figure

plot(ppg,'r')

hold on

plot(all_trough_location,all_trough_value,'ko')

[Star Strider]

I do not understand what you want to do.

I do not see any specific dicrotic notches (that would concern me if I knew the corresponding blood pressures), only the individual pulse diastolic minima. The EKG R-deflection (that I added to the plot) approximately corresponds to the diastolic minimum, as I would expect.

There appears to be a periodic baseline variation in the PPG signal that should be straightforward to eliminate with a highpass or bandpass filter.

LD = load(websave('0009_8min.mat','https://www.mathworks.com/matlabcentral/answers/uploaded_files/1251227/0009_8min.mat.mat'))

LD = struct with fields:
    SFresults: [1×1 struct]
       labels: [1×1 struct]
         meta: [1×1 struct]
        param: [1×1 struct]
    reference: [1×1 struct]
       signal: [1×1 struct]

Fs = LD.param.samplingrate.pleth;                   % Sampling Frequency

EKG = LD.signal.ecg.y;
L = numel(EKG);
tv = linspace(0, L-1, L)/Fs;                        % Time Vector
% remove first trough:
ppg = LD.signal.pleth.y;
ppg(1:10) = [];
% indentify all peaks
[all_peak_value,all_peak_location]= findpeaks(ppg);     
% identify all troughs in data
[all_trough_value,all_trough_location]= findpeaks(-ppg);
all_trough_value = -all_trough_value;   
% identify main peak (i.e. ones lower than 5)   
[main_trough_value,main_trough_location]= findpeaks(ppg,'MinPeakProminence',5);
counter = 0;
for i = 1:length(main_trough_location)
    %c=[];
    start=main_trough_location(i)-100;
    ending=main_trough_location(i);
    %c= (all_trough_location>start & all_trough_location<ending);
    %c1=c(end);
    [r,c] = find(all_trough_location>start & all_trough_location<ending);
end
figure
plot(tv(11:end), ppg,'r')
hold on
plot(tv(all_trough_location),all_trough_value,'ko')
plot(tv, EKG)                                           % Add EKG Plot
xlim([0 6.5])

.