MATLAB Answers

0

Exponential decay extrap to y=0

Asked by Jamie Hizzett on 26 Jul 2019
Latest activity Commented on by infinity
on 14 Aug 2019
Hi peeps,
I am trying to fit an exponential decay curve to my data;
x = [0;10;20;30;40;50;60]; %time in minutes
y = [12.5; 7.8;5.1;3.2;2;1.2;0.7];
I am reasonably new to matlab and I am unsure about how to fit an exponential curve to the plot that reaches y=0
Any assistance would be greatly appreciated.
Many thanks,
Jamie

  0 Comments

Sign in to comment.

2 Answers

Answer by infinity
on 26 Jul 2019
 Accepted Answer

Hello,
You could try this
clear
close all
x0 = linspace(0,100,100);
x = [0;10;20;30;40;50;60]; %time in minutes
y = [12.5; 7.8;5.1;3.2;2;1.2;0.7];
f = fittype('exp(-a*x+b)');
[fit1,gof,fitinfo] = fit(x,y,f,'StartPoint',[1 0]); % ,'StartPoint',1
figure
plot(x0,fit1(x0),'r',x,y,'o')

  3 Comments

Hi infinity,
Thank you for your answer, it works really well! How would I modify it to apply the same curve to the same y-axis data, but for x = [0;357;714;892.5;1071;1249.5;1428]?
I have been trying to tinker with it but with no joy. I tried extending the x0 to fit over the new x-axis, but it ruins the extrap
infinity
on 12 Aug 2019
Hello,
You could try this one
clear
close all
x0 = linspace(0,3000,100);
% x = [0;10;20;30;40;50;60]; %time in minutes
x = [0;357;714;892.5;1071;1249.5;1428];
y = [12.5; 7.8;5.1;3.2;2;1.2;0.7];
f = fittype('exp(-a*x+b)');
[fit1,gof,fitinfo] = fit(x,y,f,'StartPoint',[0 1]); % ,'StartPoint',1
figure
plot(x0,fit1(x0),'r',x,y,'o')
Thank you :)

Sign in to comment.


Answer by Image Analyst
on 12 Aug 2019

I don't have the curve fitting toolbox so I did it with fitnlm() from the Statistics and Maching Learning Toolbox.
Do I get the same formula and same zero crossing value as the other answer? It was not posted for some reason.
Here is my code. Don't be afraid that it's long. The meat of the program is only like 4 or 5 lines of code - the rest is just to make a fancy, illustrative plot for you.
% Uses fitnlm() to fit a non-linear model (an exponential decay curve) through noisy data.
% Requires the Statistics and Machine Learning Toolbox, which is where fitnlm() is contained.
% By Image Analyst
% 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.
xTraining = [0;357;714;892.5;1071;1249.5;1428];
yTraining = [12.5; 7.8;5.1;3.2;2;1.2;0.7];
% Now we have noisy training data that we can send to fitnlm().
% Plot the noisy initial data.
plot(xTraining, yTraining, '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.
tbl = table(xTraining, yTraining);
% Define the model as Y = a * exp(-b*x) + c
% 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)) + b(3);
% Guess values to start with. Just make your best guess.
a = 16 % Arbitrary sample values I picked.
b = .001
c = -3
beta0 = [a, b, c];
% 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'}
% Now we're done and we can plot the smooth model as a red line going through the noisy blue markers.
hold on;
% First plot the continuous curve, with hundred of points.
% Create smoothed/regressed data using the model and the training x values:
xContinuous = linspace(min(xTraining), 1800, 2000); % 2000 points (about one pixel per point across the screen).
yFitted = coefficients(1) * exp(-coefficients(2) * xContinuous) + coefficients(3);
plot(xContinuous, yFitted, 'r-', 'LineWidth', 2, 'MarkerSize', 15);
% Next plot markers at just where the training x was.
% Create smoothed/regressed data using the model and the training x values:
yFitted = coefficients(1) * exp(-coefficients(2) * xTraining) + coefficients(3);
plot(xTraining, yFitted, 'r+', 'LineWidth', 2, 'MarkerSize', 15);
% Fancy up the plot
grid on;
title('Exponential Regression with fitnlm()', 'FontSize', fontSize);
xlabel('X', 'FontSize', fontSize);
ylabel('Y', 'FontSize', fontSize);
% Find out where this curve crosses the x axis (at y = 0) by solving the equation.
% Y = a * exp(-b*x) + c = 0
% exp(-b*x) = -c/a;
% -b*x = log(-c/a)
% x = (-1/b) * log(-c/a)
xZero = (-1/coefficients(2)) * log(-coefficients(3) / coefficients(1))
plot(xZero, 0, 'r*', 'MarkerSize', 30, 'LineWidth', 2);
% Make x axis bolder.
line(xlim, [0, 0], 'Color', 'k', 'LineWidth', 2);
legendHandle = legend('Training Y', 'Fitted Y curve', 'Fitted Y at training x', 'xZero', 'x axis', 'Location', 'northeast');
legendHandle.FontSize = 25;
formulaString = sprintf('Y = %.3f * exp(-%.3f * X) + %.3f\nCrosses y = 0 at x = %.2f', ...
coefficients(1), coefficients(2), coefficients(3), xZero)
text(400, 11, formulaString, 'FontSize', 25, 'FontWeight', 'bold');
% Set up figure properties:
% Enlarge figure to full screen.
set(gcf, 'Units', 'Normalized', 'OuterPosition', [0, 0.05, 1, .95]);
% 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')
0000 Screenshot.png

  2 Comments

Thank you very much! It is interesting to see two completely different scripts tackling the same problem. The code written by infinity gets close to 0 but never reaches it! R2 = 0.98.
graph.bmp
infinity
on 14 Aug 2019
Hello,
Depending on how do you choose the type of exponential form. You can get the nearly the same result as @image by modifying the code a bit
clear
close all
x0 = linspace(0,1800,100);
% x = [0;10;20;30;40;50;60]; %time in minutes
x = [0;357;714;892.5;1071;1249.5;1428];
y = [12.5; 7.8;5.1;3.2;2;1.2;0.7];
f = fittype('a*exp(-b*x)+c');
[fit1,gof,fitinfo] = fit(x,y,f,'StartPoint',[0 0 0]); % ,'StartPoint',1
figure
plot(x0,fit1(x0),'r',x,y,'o')
hold on
plot(x0,zeros(1,length(x0)))

Sign in to comment.