# Calculate area under curve of multiple peaks

1 view (last 30 days)
Rebecca Jones on 17 Jan 2020
Answered: Star Strider on 17 Jan 2020
Hello!
I have a graph that I have used findpeaks to identify multiple peaks in the image (as shown divided by pink). I'd like to measure the area of each peak but am struggling to understand the explanations of cumtrapz on here, which is what I think I need to use.
Can anybody help me out?
Thank you!
My code and the image I'm using is attached.
My code is also here as follows:
% The code will process all TIF files
ALLfiles = dir('*0013.jpg');
% Make empty tables to be filled in the loop
% ScratchSize = table('Size',[0,4],'VariableTypes',{'string' 'string' 'string' 'double'});
% ScratchSize.Properties.VariableNames = {'Condition' 'Time' 'Well' 'WoundArea'};
Streams=strings(length(ALLfiles),6);
% HealOrder = table('Size',[0,4],'VariableTypes',{'string' 'string' 'double' 'double'});
% HealOrder.Properties.VariableNames = {'Condition' 'Time' 'AverageHealedArea' 'StandardDeviation'};
dbstop if error
% Now the code will loop through all the files identified
for i=1:length(ALLfiles)
dbstop if error
%Get the name of the file and remove the extention. Identify the
%condition, the well, and the time.
thisfile=ALLfiles(i).name;
if size(OrigImage,3)==3
GrayScaleScratchFULL=rgb2gray(OrigImage);
else
GrayScaleScratchFULL=OrigImage;
end
%%
greenchannel=OrigImage(:,:,2);
greenchannel(greenchannel>140)=255;
binarytwist=imbinarize(greenchannel);
binaryImage = ~bwareaopen(~binarytwist, 4000);
M = bwareaopen(binaryImage,5);
se = strel('disk',5);
N=imclose(M, se);
%%
texture1 = GrayScaleScratch;
texture1(~N) = 0;
texture2 = GrayScaleScratch;
texture2(N) = 0;
boundary = bwperim(N);
boundary = imdilate(boundary, strel('disk',1));
segmentResults = GrayScaleScratch;
segmentResults(boundary) = 255;
fulldestination = fullfile(destdirectory, thisfile + "_Outline.jpg");
imwrite(segmentResults, fulldestination);
%%
imagebar=sum(binaryImage);
crestbar=1024-imagebar;
%crestbarfig=bar(crestbar);
[peaks,loc]=findpeaks(crestbar, 'MinPeakProminence', 30);
%crestbarfig=bar(crestbar);
findpeaks(crestbar, 'MinPeakProminence', 30);
hold on
text(loc+.02,peaks,num2str((1:numel(peaks))'));
n=length(peaks);
Streams(i,1)=thisfile;
for nn=1:n
Streams(i,nn+1)=peaks(nn);
end
StreamLengths = array2table(Streams);
end
writetable(StreamLengths,destdirectory+'StreamLengths.xlsx');

Cameron B on 17 Jan 2020
Not the most robust code, but it will work for what you have. You will need to check your findpeaks function to make sure it's only getting the points you need.
[lowpeaks,lowloc]=findpeaks(-crestbar, 'MinPeakProminence', 30,'Threshold',2);
cc = 1;
data = [transpose(1:length(crestbar)),transpose(crestbar)];
for ee = 1:length(lowloc) + 1
if ee == 1
intval(cc,1) = trapz(data(1:lowloc(ee),1),data(1:lowloc(ee),2));
elseif ee == length(lowloc) + 1
intval(cc,1) = trapz(data(lowloc(ee-1):end,1),data(lowloc(ee-1):end,2));
else
intval(cc,1) = trapz(data(lowloc(ee-1):lowloc(ee),1),data(lowloc(ee-1):lowloc(ee),2));
end
cc = cc + 1;
end

Star Strider on 17 Jan 2020
Using cumtrapz:
[peaks,loc]=findpeaks(crestbar, 'MinPeakProminence', 30);
[trofs,troflocs] = findpeaks(-crestbar, 'MinPeakProminence', 30, 'MinPeakDistance',10);
troflocs_ext = [1 troflocs numel(crestbar)];
AUC = cumtrapz(crestbar);
for k = 1:numel(troflocs_ext)-1
AUC_Seg(k) = AUC(troflocs_ext(k+1)) - AUC(troflocs_ext(k))
end
The segment areas are:
AUC_Seg =
20341 11820.5 13784.5
where the first one is between 1 and the first trough (at 472), the second segment between that and the second trough (at 542) and the last between that and the end of the vector. Since it is 0 elsewhere, ther is no need to segment the areas outside the peaks.