Fit 2 exponential to a curve

Hi,

I am working on fitting algorithm to fit a double exponential to a given curve. These exponentails are RC time constants for two different combinations of R and C, so for R1*C1 the time constant is T1 and for R2*C2 the time constant is T2. Starting from point A in my top image below, the curve follows T1 till Point B, and from Point C (Point C would be much closer to Point B, but for clarity I have shown it to be further away from Point B) till Point D (or even further, but the curve is flat so it doesn't matter) the curve follows T2. The curve is obtained by experiment and I would like to find out Point B, Point C, T1 and T2. I have tried multiple fitting algorithms, nothing works, except for the one where I assume Point B and Point C coinside. The equation that works:

ft = fittype('S + A*((h*(1-exp((T-x)/T1))) + (1-h)*(1-exp((T-x)/T2)))','problem',{'S','A','T'});
fo = fit(x,y,ft,'problem',{Sx,Ax,Tx},'StartPoint',[1,10^-2,10^-2],'Lower',[0,10^-6,10^-6],'Upper',[1,10^-2,10^-2]);

Where:

Sx = Starting Amplitude [Known]

Ax = Final Amplitude [Known]

Tx = Starting Time [Known]

x = Data to fit [Independent Variable]

y = Data to fit [Dependent Variable]

h = coefficient 1 [ Min = 0, Max = 1] = Fraction of amplitude (Ax) that follows time constant T1 [where Point B and Point C are together]

T1 = coefficient 2 [ Min = 10^-6, Max = 10^-2] = Time constant for first combination of RC i.e. R1*C1

T2 = coefficient 3 [ Min = 10^-6, Max = 10^-2] = Time constant for second combination of RC i.e. R2*C2

Original Curve:

Zoomed in image of the fittied Curve, Red curve is Original Signal and Blue curve is the fitted curve. The dashed horizontal line in below image is the value of 'h' in above fitting algorithm that provides the seperation point between the RC time constants.

In the zoomed in image, it is clear that that fitted curve is slightly overestimating the first RC time constant T1 and underestimating the second RC time constant T2. In this particular example of the recorded data, shown above, the fitting algorithm works good and gives me expected values. But, in another example, shown below (again Red curve is original signal, Blue curve is fitted curve and dashed horizontal line is the value of 'h'), the fitting algorithm gives inaccurate results.

If there was a way I could seperate Point B and Point C (shown in first image) in the fitting algorithm, meaning T1 ends at Point B and T2 starts at Point C then I am quite sure I would get exact result in all recorded cases. Any help in this regard is highly aprriciated, thanks in advance!

[darova]

Can't you manually correct your coefficients?

clc,clear

n = 100;
x = linspace(0,10,n);
y = sin(x)+rand(1,n)/5;     % original data
ft = @(a,b,x)a*sin(x)+b;    % function to fit
f = fit(x',y',ft);          % fit data
cc = coeffvalues(f);        % get coefficients
a = cc(1);
b = cc(2)+0.5;              % change b coefficient
plot(x,y,x,ft(a,b,x))

[Chaitanya Jha]

Hi Darova! thanks for the answer, but no, I cannot change the coefficients manually since in real experiment I will not know what values it will have. The curves presented here are from known values that's why I know if the fit is correct or not.

[darova]

So you want to fit data separately?

[Chaitanya Jha]

Yes, but the Point B and Point C is to be determined by the fitting algorithm as it varies with each experiment. What I expect is that between Point B and Point C the two exponentials act together as the curve is in the transition period and so I want to exclude it from my calculations. Point B is the most important point of distinction.

I have tried multiple fitting equations to achieve this seperation of Point B and Point C but nothing works, the equations that I have tried are:

ft = fittype('(x < h1)*(S+A*(1-exp(T-x)/t1)) + (x > h2)*(S+A*(1-exp((T-x)/t2)))','problem',{'S','A','T'});
ft = fittype('(x < h1)*(S+A*(1-exp(T-x)/t1)) + (x > (h1+h2))*(S+A*(1-exp((T-x)/t2)))','problem',{'S','A','T'});
ft = fittype('S + A*(1-exp((T-x)/((h1*t1)+(h2*t2))))','problem',{'S','A','T'});
ft = fittype('S + A*(1-exp((T-x)/((h1*t1)+((1-(h1+h2))*t2))))','problem',{'S','A','T'});

I have also tried to make a different file with my equation to be used in fitting algorithm:

function y = MyFitEqn(x,S,T,A,h1,h2,t1,t2)
if x < h1
    y = S + A*(1-exp((T-x)/t1));
elseif x > (h1+h2)
    y = S + A*(1-exp((T-x)/t2));
end

Where h1 will correspond to Point B and h2 will correspond to Point C.

There are many more variations that I have tried and in almost all the cases decipted here the coefficients are stuck at the lower or upper bounds. I have no idea how to force the fitting algorithm to take the Point B and Point C seperately.

[darova]

Why don't just to separately fit the data?

x1 = x(x<h1);
y1 = y(x<h1);
x2 = x(x>h1+h2);
y2 = y(x>h1+h2);
f1 = fit(x1,y1,func);
f2 = fit(x2,y2,func);

[Chaitanya Jha]

Because I don't know h1 and h2, it depends on Point B and Point C and these are experiment dependent, varies with each data set.

[darova]

Can you attach the data?

[Chaitanya Jha]

I have attached the data in the file 'CurveData.mat', this is a part of the data collected during an experiment and so the time starts at 0.04 seconds, sampling frequency is 250 Kilo Hertz so difference between consecutive x-values is 0.000004 seconds. In this file the data is given as:

M.x => Time data or independent variable

M.y => Amplitude data or dependent variable

M.Sx => Starting amplitude

M.Ax => [Final amplitude where the curve is flat] - M.Sx i.e. the difference in amplitude from starting point to reaching a stable state

M.Tx => Starting Time

M.y is filtered data, original data is quite noisy so I ran it through filter: sgolayfilt(original_y_data,5,101)

[darova]

What about this simple fit function?

load CurveData.mat
x = M.x;
y = M.y;

a = 2e4;
b = 2*400;
c = 0;
ft = @(a,b,c,x) -1./(a*x-b) + c;
plot(x,y,x,ft(a,b,c,x))

x1 = x(x>0.05);
y1 = y(x>0.05);
c = c+mean(y1)-mean(ft(a,b,c,x1));
line(x,ft(a,b,c,x),'color','r')
legend('original data','fit','corrected fit')

[Chaitanya Jha]

No, here the fitted curve and corrected fit curve are both having an amplitude much lower than the original curve and also it doesn't give me any details about Point B and Point C.

[darova]

Sorry it was my best shot. I don't have more ideas

[Chaitanya Jha]

Thanks a lot for helping me! Highly appriciated :)

Answer 1

[Image Analyst]

Are you sure you're picking the right model? Maybe it's the "rate equation" instead of an exponential decay. I attach demos for each type of equation.

[Chaitanya Jha]

Thanks! I tried the method but it is not working for me. I do not have any predictor variables and doing this only with x and y is not producing results rather giving me just errors. I will try again soon, it might take some time. I will update you if it works.

[Image Analyst]

If you need me to do it, then attach your x and y in a .mat file with the paper clip icon.

save('answers.mat', 'x', 'y');