How to optimize a bi-exponential signal fitting?

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');

[praguna manvi]

Hi Susanna, please upload the data file

[Umar]

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.

Answer 1

[Kaustab Pal]

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!