lsqcurvefit not working, large f(x) values during iterations

I am trying to use lsqcurvefit to fit data using a nonlinear function (modelFun in the code), but I just can't get it to work. Do you see anything wrong with my equation or code that would give the figure I've attached here?

t = VarName1; %Time data points in ms

y = VarName2; %Current data points

Diff = 2000; %Diffusion coefficient in um^2/ms

Co = 0.0001; %Initial concentration in uM

%Plot of actual data

plot(t,y,'--b');

xlabel('Time');

ylabel('Current');

%This equation has four parameters: p1, p2, p3, and p4. p1 is k1 (on rate

%constant), p2 is k2 (off rate constant), p3 is the length of the unstirred

%layer (L), and p4 is iN (total current). modelFun = @(p,t) (-169 * (1 - (0.0004 / (p(1) * 0.0001 + 0.0004)) * exp(-((p(1) * 0.0001 + 0.0004) .* t)) + 0.0004 / (p(1) * 0.0001 + 0.0004)));

%Non-linear least squares fit of data using model function

format long

%opts = statset('MaxIter',1000);

options=optimset('disp','iter','LargeScale','off','TolFun',0.000001,'MaxIter',100000,'MaxFunEvals',1000000);

lb = [100]; %k1 is in uM^(-1) ms^(-1), k2 is in ms^(-1), L is in um, and total current, I, is in nA

ub = [400];

startingVals = [296];

coefEsts = lsqcurvefit(modelFun, startingVals, t, y, lb, ub, options);

xgrid = linspace(0,3.11E+5,100);

line(xgrid, modelFun(coefEsts, xgrid), 'Color', 'r');

Answer 1

[John D'Errico]

To quote you:

%This equation has four parameters: p1, p2, p3, and p4. p1 is k1 (on rate
%constant), p2 is k2 (off rate constant), p3 is the length of the unstirred
%layer (L), and p4 is iN (total current).

No. Your model has ONE parameter. P(1).

modelFun = @(p,t) (-169 * (1 - (0.0004 / (p(1) * 0.0001 + 0.0004)) * exp(-((p(1) * 0.0001 + 0.0004) .* t)) + 0.0004 / (p(1) * 0.0001 + 0.0004)));

Everything else besides t is a CONSTANT. Fixed. Immutable. Frozen.

Can you truly be surprised that it did not converge to something reasonable?

Worse, while I lack any data,

modelFun = @(p,t) (-169 * (1 - (0.0004 / (p(1) * 0.0001 + 0.0004)) * exp(-((p(1) * 0.0001 + 0.0004) .* t)) + 0.0004 / (p(1) * 0.0001 + 0.0004)));

t = linspace(0,3,100);
plot(t,modelFun(296,t))

So your model will not even have the chance of a snowball in a very hot place to survive.

[Kelly McGuire]

Ok, but before I had those 4 parameters and I still get the same problem...sorry, I am new to curve fitting, so what exactly is making my equation not fit my data with an exponential, even when I have 4 parameters in the equation? Thanks for the help!

[John D'Errico]

NO. See what I showed in the plot. Not a chance. EVER. Why in the name of god did you decide to fix almost all of the parameters?

In fact, try this:

plot(t,modelFun(100,t))
plot(t,modelFun(296,t))
plot(t,modelFun(400,t))

Do you see any difference in those plots at all?

[John D'Errico]

plot(t',[modelFun(100,t);modelFun(296,t);modelFun(400,t)]')

Which of those curves looks the most like the aim curve? :)

[John D'Errico]

"Ok, but before I had those 4 parameters and I still get the same problem...sorry, I am new to curve fitting, so what exactly is making my equation not fit my data with an exponential, even when I have 4 parameters in the equation?"

You don't have 4 parameters in that equation. I don't know what the original equation was, so I can't tell you why it does not fit.

At the very least, you may have terrible starting values. Or the fundamental shape of that equation may be not at all when your data is.

Again, I don't have your equation, nor do I even have your data to show you how to fit it.

[Kelly McGuire]

Ok, so now I have the two parameters I need not fixed anymore...

modelFun = @(p,t) -169 .* (1 - (p(2) ./ (p(1) .* 100 + p(2))) .* exp(-((p(1) .* 100 + p(2)) .* t)) + p(2) ./ (p(1) .* 100 + p(2)));

[John D'Errico]

Better.

What range do you expect those parameters to live in? I assume you may want this:

 100 <= P(1) <= 400

although I think there are problem with that range.

but anyway, where does P(2) live?

The point is, all I will do at first is plot your function, substituting values of P(1) and P(2).

[John D'Errico]

In fact, as close as I can come to your target curve is this:

plot(t,modelFun([.1,-1],t))

Note that is it not even close to the curve you originally showed. Just because it has an exp in there does not mean it will have the desired shape.

[John D'Errico]

Ok. Its difficult to visualize this without some help. so...

syms t p1 p2
f = -169 .* (1 - (p2 ./ (p1 .* 100 + p2)) .* exp(-((p1 .* 100 + p2) .* t)) + p2 ./ (p1 .* 100 + p2));

pretty(f)
p2 exp(-t (100 p1 + p2)) 169      169 p2
---------------------------- - ----------- - 169
         100 p1 + p2           100 p1 + p2

Now we should be able to see a bit. The problem is, I'm not sure that model can possibly fit your target curve for any values of p(1) and p(2).

So lets look at that question. If we take the second and third terms, and recognize that the model needs to approach roughly 10 from your target curve as t-->inf.

Then we can solve for the value of 100*p(1)+p(2) as approximately:

100*p(1)+p(2) = -169*P(2)/179

In order for this curve to approach an asymptote as t-->inf, we NEED 100*p(1)+p(2) to be a positive number, because it is in the exponential term, and we need that to be a negative exponential.

Therefore, we KNOW that p(2) must be negative. This is true with absolute certainty if your model is to be believed.

But you seem to have indicated that p(2) is a secondary rate parameter. So I'm not at all sure you have the right model.