How to optimize a bi-exponential signal fitting?

5 次查看（过去 30 天）
Susanna Rampichini
Susanna Rampichini 2024-3-29
回答： Kaustab Pal 2024-8-20
Hi everybody!
I'm trying to solve a fitting curve problem concerning a double-exponential fit.
I have this data set (Data.txt) that behaves like a bi-exponential function following a resting condition. Therefore, since x<0 the signal is almost steady whereas, for x>0 the signal presents a first exponential phase (lasting about 20 sec) followed by a second exponential phase.
The problem is that the resulting fitting does not model the signal behavior in the correct way. Specifically, the first exp phase is too short and/or too low even though correct StartingPoint and Upper and Lower limits are set.
Lastly, I've written two algorithms with and without "Exclude" options. The second one presents a strange pattern at the onset of the second exp function.
I'll appreciate any suggestions.
This is the code.
[x, y] = readvars('prova.txt');
tRest = abs(min(x)); % Resting phase duration (sec)
cdp = 20; % estimated duration (sec) of CardioDynamicPhase (cdp)
cdpEpoch = 5; % epoch duration (sec) to compute mean value of amplitude of cdp
plot(x,y);
% Fitting StartingPoint
yRest = mean( y(1:tRest-1) ); % y value at resting condition before exp onset
A1start = mean( y(1:tRest+cdp-cdpEpoch):y(tRest+cdp) )-yRest; % mean y value of the last 5sec of the cardiodynamic phase
TD1start = 0; % starting value (sec) of the first exponential phase onset
tau1Start = 8; % starting value (sec) of the first exponential phase time constant
A2start = max(y)-A1start-yRest; % starting value of the second exponential phase amplitude
TD2start = 30; % starting value (sec) of the second exponential phase onset
tau2Start = 30; % starting value (sec) of the second exponential phase time constant
% Fitting from exponential onset
g = fittype( @(A1, TD1, tau1, A2, TD2, tau2, yRest, x) yRest+A1.*(1-exp(-(x-TD1)/tau1))+A2.*(1-exp(-(x-TD2)/tau2)), 'problem', 'yRest');
f1 = fit( x(tRest:size(x,1)), y(tRest:size(x,1)), g, 'problem', yRest, 'StartPoint', [A1start, TD1start, tau1Start, A2start, TD2start, tau2Start], 'lower' , [ A1start-100, 0, 5, A2start-200, 10, 10], 'upper', [ A1start+200, 20, 40, max(y), 60, 40], 'MaxIter', 10000);
pg =coeffvalues(f1);
A1 = pg(1);
TD1 = pg(2);
tau1 = pg(3);
A2 = pg(4);
TD2 = pg(5);
tau2 = pg(6);
yFit( 1:(tRest+floor(TD1)) ) = yRest;
yFit( tRest+floor(TD1)+1 : tRest+floor(TD2) ) = yRest + A1.*(1-exp(-( x(tRest+floor(TD1)+1 : tRest+floor(TD2))-TD1 )/tau1));
yFit( tRest+floor(TD2)+1 : size(y,1)) = yFit(tRest+floor(TD2))+ A2.*(1-exp(-( x(tRest+floor(TD2)+1 : size(y,1))-TD2 )/tau2));
% fitting using Exclude option
f2 = fit( x, y, g, 'problem', yRest, 'StartPoint', [A1start, TD1start, tau1Start, A2start, TD2start, tau2Start], 'lower' , [ 390, 0, 5, A2start-200, 10, 10], 'upper', [ 400, 20, 40, max(y), 60, 40], 'Exclude', x<0); %Same fitting algortihm that exclude negative x values
pg2=coeffvalues(f2);
A1e = pg2(1);
TD1e = pg2(2);
tau1e = pg2(3);
A2e = pg2(4);
TD2e = pg2(5);
tau2e = pg2(6);
yFite( 1:(tRest+floor(TD1e)) ) = yRest;
yFite( tRest+floor(TD1e)+1 : tRest+floor(TD2e) ) = yRest + A1e.*(1-exp(-( x(tRest+floor(TD1e)+1 : tRest+floor(TD2e))-TD1e )/tau1e));
yFite( tRest+floor(TD2e)+1 : size(y,1)) = yFit(tRest+floor(TD2e))+ A2e.*(1-exp(-( x(tRest+floor(TD2e)+1 : size(y,1))-TD2e )/tau2e));
% Plot the fitting curves
hold on
plot (x, yFit, 'r');
plot (x, yFite, 'k');
  2 个评论
Umar
Umar 2024-8-9
Hi @Susanna Rampichini,
The code that you have provided is utilizing data from 'prova.txt' but in your comments, you mentioned “I have this data set (Data.txt) that behaves like a bi-exponential function”,
[x, y] = readvars('prova.txt');
Also, you mentioned, “The problem is that the resulting fitting does not model the signal behavior in the correct way”, is this issue occurring due to algorithm with or without "Exclude" options. Please clarify these comments for us to help resolve the problem.

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

Kaustab Pal
Kaustab Pal 2024-8-20
Hi @Susanna Rampichini,
It seems there might be an issue with the way you're calculating A1start. I've generated similar data and adjusted the computation for A1start, which resulted in a fitted curve that closely matches the original data. Below is the revised code for the A1start variable:
A1start = mean(y(tRest+cdp-cdpEpoch:tRest+cdp)) - yRest;
The resulting plot should provide a visual comparison between the original data and the fitted curves, both with and without the exclusion of initial data points.
This approach should help in achieving a more accurate representation of the data.
I hope this helps!

类别

在 Help Center 和 File Exchange 中查找有关 Linear and Nonlinear Regression 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by