[Solved] Fitting exponential decay function

I'm trying to fit an exponential decay to a dataset of x and y values (3001 each). Using other software I was able to calculate a k_off around 0.02 however using the fittype and fit to replicate this in MATLAB I get the following results:

Code:

s1 = sprintf('%f*exp(-koff*', y_equil);  %(For y_equil = 0.148356)
s2 = 'x)+plateau'
eq_string = strcat(s1, s2);
f = fittype(eq_string);
f1 = fit(x,y, f)
plot(f1,x,y)
xlabel('Time (s)')
ylabel('Concentration')

------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------

Results:

f1 =

General model:

f1(x) = 0.148356*exp(-koff*x)+plateau

Coefficients (with 95% confidence bounds):

koff = 5.888 (4.858, 6.917)

plateau = 0.1291 (0.1289, 0.1293)

This fit finds the correct value for the plateau but not for the k_off. For some reason the fit is horrible.

How do I correct this?

----------------------------------------------------------------

[Solved]: - I simply modified his code:

Credit to Jing Chen

https://www.mathworks.com/matlabcentral/fileexchange/37321-fitexponential-m

Corrected Code:

yBase           = y_equil;
y               = y - y_equil; 
fh_objective    = @(param) norm(param(2)+(param(3)-param(2))*exp(-param(1)*(x-x(1))) - y, 2);
initGuess(1)    = -(y(2)-y(1))/(x(2)-x(1))/(y(1)-y(end));
initGuess(2)    = y(end);
initGuess(3)    = y(1);
param           = fminsearch(fh_objective,initGuess);
k2               = param(1);
yInf            = param(2) + yBase;
y0              = param(3) + yBase;
yFit2            = yInf + (y0-yInf) * exp(-k2*(x-x(1)));

Answer 1

[Image Analyst]

You can use fitnlm(). See attached demo. I can't use your actual data because you forgot to attach it.

[Knut Andre Skaug]

This model works really good for the type of data I am modeling.

y = a * e^ b x + c (1)

I am however, trying to transform y into a linear shape, as one could if the model was

y = a * e^ b x (2)

If the model were as above, I could approximate a straight line using:

u = log(exp(x));
v = log(y);

P = polyfit (u,v,1)
Yintercept = exp(P(2));
Slope = P(1);

Could you offer some advice into how to carry out a similar transform in order to fit a straight line to the model?

[Walter Roberson]

fit() with model 'exp' will automatically do the log transform internally.

However it you have a linear constant in the equation then you cannot use the log transform, except to get approximate starting values for a more accurate fit (assuming that none of the values are negative.)

[Image Analyst]

Knut, I'm not sure the coefficients would be the same if you did it that way as opposed to using the non-linear fit. Might - I'm just not sure.

It's easy to modify my code to get rid of the offset if you don't want it. Just change the formula and pass in only 2 values in the beta0 vector.

% Uses fitnlm() to fit a non-linear model (an exponential decay curve, Y = a * exp(-b*x)) through noisy data.
% Requires the Statistics and Machine Learning Toolbox, which is where fitnlm() is contained.
% Initialization steps.
clc;    % Clear the command window.
close all;  % Close all figures (except those of imtool.)
clear;  % Erase all existing variables. Or clearvars if you want.
workspace;  % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 20;

% Create the X coordinates from 0 to 20 every 0.5 units.
X = 0 : 0.5 : 20;

% Define coefficients and function that the X values obey.
a = 10 % Arbitrary sample values I picked.
b = 0.4
Y = a * exp(-X * b); % Get a vector.  No noise in this Y yet.

% Add noise to Y.
Y = Y + 0.5 * randn(1, length(Y));

%--------------------------------------------------------------------------------------------------------------------------------------
% Now we have noisy training data that we can send to fitnlm().
% Plot the noisy initial data.
plot(X, Y, 'b*', 'LineWidth', 2, 'MarkerSize', 15);
grid on;

% Convert X and Y into a table, which is the form fitnlm() likes the input data to be in.
% Note: it doesn't matter if X and Y are row vectors or column vectors since we use (:) to get them into column vectors for the table.
tbl = table(X(:), Y(:));
% Define the model as Y = a * exp(-b*x)
% Note how this "x" of modelfun is related to big X and big Y.
% x((:, 1) is actually X and x(:, 2) is actually Y - the first and second columns of the table.
modelfun = @(b,x) b(1) * exp(-b(2)*x(:, 1));  

% Guess values to start with.  Just make your best guess.
aGuessed = 11 % Arbitrary sample values I picked.
bGuessed = 0.5
beta0 = [aGuessed, bGuessed]; % Guess values to start with.  Just make your best guess.

% Now the next line is where the actual model computation is done.
mdl = fitnlm(tbl, modelfun, beta0);
% Now the model creation is done and the coefficients have been determined.
% YAY!!!!

% Extract the coefficient values from the the model object.
% The actual coefficients are in the "Estimate" column of the "Coefficients" table that's part of the mode.
coefficients = mdl.Coefficients{:, 'Estimate'}
% Create smoothed/regressed data using the model:
yFitted = coefficients(1) * exp(-coefficients(2)*X);
% Now we're done and we can plot the smooth model as a red line going through the noisy blue markers.
hold on;
plot(X, yFitted, 'r-', 'LineWidth', 2);
grid on;
title('Exponential Regression with fitnlm()', 'FontSize', fontSize);
xlabel('X', 'FontSize', fontSize);
ylabel('Y', 'FontSize', fontSize);
legendHandle = legend('Noisy Y', 'Fitted Y', 'Location', 'north');
legendHandle.FontSize = 30;

% Place formula text roughly in the middle of the plot.
formulaString = sprintf('Y = %.3f * exp(-%.3f * X)', coefficients(1), coefficients(2))
xl = xlim;
yl = ylim;
xt = xl(1) + abs(xl(2)-xl(1)) * 0.325;
yt = yl(1) + abs(yl(2)-yl(1)) * 0.59;
text(xt, yt, formulaString, 'FontSize', 25, 'FontWeight', 'bold');

% Set up figure properties:
% Enlarge figure to full screen.
set(gcf, 'Units', 'Normalized', 'OuterPosition', [0 0 1 1]);
% Get rid of tool bar and pulldown menus that are along top of figure.
% set(gcf, 'Toolbar', 'none', 'Menu', 'none');
% Give a name to the title bar.
set(gcf, 'Name', 'Demo by ImageAnalyst', 'NumberTitle', 'Off') 

Answer 2

[dpb]

You forced the model to have f(0) = 0.148356 + 0.1291 = 0.2775.

It takes a really large decay constant to make up for that.

The initial value looks like it should be roughly that value of 0.148 so the multiplicating factor should be (0.148-plateau) if it is intended to fix that initial value rather than estimate three parameters instead of only two.

Try

y_equil=0.148356;
s=sprintf('(%f-plateau)*exp(-koff*x)+plateau', y_equil);
f = fittype(s)
f = 
     General model:
     f(koff,plateau,x) = (0.148356-plateau)*exp(-koff*x)+plateau
>> 

and see if don't have better luck

[Grant Huckels]

Thank you for your time and effort! The k that I obtained using your code was around 0.8, and ~0.6 for the 'Robust' on.

However I found an alternative way thanks to Jing Chen's code. I left the updated material for future users that may encounter the same issues I did.