# How to Extract/Fit bell curves under a time series...

26 views (last 30 days)

Show older comments

Claudio Iturra
on 12 Apr 2023

Commented: Star Strider
on 13 Apr 2023

##### 2 Comments

### Accepted Answer

Star Strider
on 13 Apr 2023

Edited: Star Strider
on 13 Apr 2023

In order to uniquely identify the peaks and identify their parameters, it is necessary to restrict the range of ‘x’ values for each peak. That is done by identifying the valleys (‘vys’) between them and using those index values to restrict the regression —

% clear all,close all;

x = linspace(0,1,1000);

Pos = [1 2 3 5 7 8]/10;

Hgt = [3 4 4 2 2 3];

Wdt = [2 6 3 3 4 6]/100;

for n = 1:length(Pos)

Gauss(n,:) = Hgt(n)*exp(-((x - Pos(n))/Wdt(n)).^2);

end

% my raw time serie consists of a sum of bell curves + some noise....

PeakSig = sum(Gauss);

subplot(1,2,1)

plot(x,PeakSig,'r')

subplot(1,2,2)

plot(x,Gauss,'--',x,PeakSig)

%--------------------

% Now I am looking for the method to reconstruct the

% bell curves (inverse) from the time series in order to evaluate

% the area under each bell curves.

[pks,plocs] = findpeaks(PeakSig);

[vys,vlocs] = findpeaks(-PeakSig);

vlocsv = [1 vlocs numel(x)];

objfcn = @(b,x) b(1).*exp(-(x-b(2)).^2/b(3));

for k = 1:numel(plocs)

idxrng = vlocsv(k) : vlocsv(k+1);

B0 = [pks(k); x(plocs(k)); 1E-6];

B(:,k) = fminsearch(@(b) norm(PeakSig(idxrng)-objfcn(b,x(idxrng))), B0);

end

Parameters = array2table(B, 'VariableNames',compose('Peak #%d',1:numel(plocs)), 'RowNames',{'Hgt','Pos','Wdt'})

% figure % Plot Individual Peaks

% tiledlayout(3,2)

% for k = 1:size(B,2)

% nexttile

% Bk = B(:,k)

% plot(x, objfcn(B(:,k),x))

% hold on

% xline(B(2,k), '-', string(B(2,k)))

% hold off

% end

figure

plot(x, PeakSig, '-r')

hold on

for k = 1:size(B,2)

Gfit(k,:) = objfcn(B(:,k),x);

plot(x, Gfit(k,:), '--')

AUC(k) = trapz(x, Gfit(k,:));

end

hold off

text(Parameters{2,:}, Parameters{1,:}, "Area = "+AUC, 'Horiz','left', 'Vert','bottom')

EDIT — (13 Apr 2023 at 14:00)

Forgot about the areas. Now included.

.

##### 4 Comments

Star Strider
on 13 Apr 2023

As always, my pleasure!

Filtering isn’t always the complete solution, because there may be peaks that would be passed by the filter yet not be the desired peaks for a specific analysis. That’s when 'MinPeakProminence' really helps!

### More Answers (0)

### See Also

### Categories

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!